=====================================
I     IMC                           I
I     (Inverse Monte Carlo)         I
I     v.2.1     from 18 May 1999    I
===================================================== 
	Author: Alexander Lyubartsev
		Division of Physical Chemistry
		Arrhenius Lab.
		Stockholm University
		S-10691 Stockholm Sweden

		e-mail: sasha@physc.su.se
======================================================

The program IMC is intended for reconstruction of pair interaction 
potentials for monoatomic multicomponent system if the radial distribution 
functions (RDF) are known. 

The program performs one iteration of the inverse (reverse) Monte 
Carlo procedure originally suggested in  A.P.Lyubartsev, A.Laaksonen, 
Phys.Rev.E, v.52(4) 3730 (1995) (relevant references are also:
A.P.Lyubartesv, A.Laaksonen, Comp.Phys.Comm,v.121-122, 57 (1999) or
A.P.Lyubartsev, A.Laaksonen, J.Chem.Phys., v.111, 11207 (1999))

The program runs MC simulation using some trial potential (at the first
iteration, the potential of mean force is used), compares resulting RDF
with the reference (input) RDF, and produces an improved potential.  
This improved potential is used at the next iteration, and the procedure 
is repeated untill convergency. 


Installation
------------
System requirements: 

- Any operating system with standard Fortran-77 compiler.
(The program has been tested on PC-Linux with g77 or Portland Group 
Fortran, and on AIX operating system, xlf fortran compiler).  
- Linear algebra library (e.g, LAPACK+BLAS, 
available from ftp://ftp.netlib.org/lapack;  LINPACK or ESSL libraries may
be also used)

Before compilation, check array boundaries in file invmc.h

You should also check subroutines which solve linear equation systems. 
The subroutines are taken from an external linear algebra library. 
Check a peace of code after commentary line 
*   Linear equations solution
(about 400 lines from the top)
Default is LAPACK library. Lines which use ESSL library are commented out.

Another place to check - normalisation of RDF 
for likewise particles. It might be either V/(N**2/2) or
V/(N*(N-1)/2), depending on definition in the input RDF. 
The program uses V/(N**2/2)

Compile using e.g:

f77 -O3 -o imc imc.f -llapack -lblas

Input files for the program:
----------------------------
1) File with input RDF
Format:
1-st line contains 4 numbers:
number of species (NTYP)
number of grid points on interval [Rmin,Rmax] for which RDF and potential 
are defined (NA);
Rmin (usually 0)
Rmax (usually cut-off distance)
Then follows NTYP lines with species names, one line per each
Then follows lines with RDF(s). Each line have 4 parameters:
1st - r
2nd - RDF(r)
3rd - species type number of the 1st particle
4th - species type number of the 2nd particle

Everything in free format.
An example of input RDF file is nacl.rdf which contains 3 RDFs obtained
in simulation of NaCl solution.

2) File with trial potential (optional)
At first iteration, the potential of mean force is used as a trial interaction
potential. On the next iteractions, output of the previous iteration is used
as input for the potential. File nacl.pot1 is an example.

3) The main input file (NAMELIST format) with parameters of the simulation.
Meaning of parameters are: 

 &INPUT
 IPRINT=5,                  ! level of output ( 5 is a good value)
 NTYP=2,                    ! num. of particle's species
 NA=200,                    ! num. of slices of [Rmin;Rmax]
 NSPEC=18,18,               ! num. of particles of each specii
 CH=1.,-1.,                 ! charges (for out-cutoff electrostatics) 
 NMKS=6100000,              ! num. of MC-steps
 NMKS0=100000,              ! num of MC steps for equilibration
 LPOT=.f.,                  ! .t. : trial potential from the input file,
                               otherwise mean force potential used
 LDMP=.t.,                  ! dump restart file at the end
 LRST=.f.,                  ! .t. - start from the restart file
			  !   otherwise - new run
 FILRDF='nacl.rdf',       ! file with reference RDF
 FILPOT='nacl.pot',       ! input file with "trial" potential 
                          ! (not needed if LPOT=.f.)
 FOUT='nacl.pot1',        ! improved potential (output), to be used as
                          ! input at next iteration        
 FDMP='nacl.dmp2',        ! restart file
 AF=3.,                ! Ewald parameter, keep this. (erfc(AF) must be small) 
                       ! Put zero if no electrostatics forces out cut-off
 FQ=9.,                ! Defines k-cut-off in reciprocal Ewald 
		       ! (exp(-FQ) must be small)
 EPS=78.,              ! dielectric permittivity
 TEMP=300.,            ! temperature
                       ! Note: Ewald parameters, Eps and temperature define
                       !  only out cut-off corrections of electrostatic 
                       !  interactions
 BOXL=19.4,            !  box size ()
 DR=5.,                !  max particle displacement at each step
 IOUT=1000,            !  output parameter
 IAV=50,               !  how often < SaSb > evaluated
 REGP=1.,              !  regularization parameter - between  0 and 1
 DPOTM=5.,             !  maximum change of the potential at this iteration
 RTM=10.,              ! keep this
 NRS=51,               ! random seed
 LPOUT=.f.,            !   not really used
 /

A sample input file is included (imc.in)
The potential is given in the reduce units (i.e., potential/kT )

Running the program:
-------------------
Run the program with the supplied input file:

imc < imc.in  [ > imc.out ]

The program produces output file "nacl.pot1" which contains improved 
interaction potential. On the next iteration, this potential may be used 
as an input potential (change parameters   LPOT=.t., FILPOT="nacl.pot1", 
FOUT="nacl.pot2" for the second iteration).

Number of iteration required for convergency is differ from case to case.
For ionic solution it is typically 3-5. If you see that the process 
disconverge, repeat specifying lower regularisation parameters (if you have
a restart file, it is possible to recalculate new potential with another 
regularisation parameter without running the simulation again)

Files:
------


README 		- this file
imc.f  		- Fortran source code
invmc.h		- include file
imc.in 		- sample input file
nacl.rdf	- sample input RDF file
nacl.pot1	- sample file with interaction potential which is obtained as
		  output and then used as an input at the next iteration

