README file for running isomolar semigrand ensemble molecular dynamics simulations (iSGMD) in ProtoMol2 **************

For an indepth description of the iSG ensemble and the equations of motion used by ProtoMol2, please see
T.I. Morrow and E.J. Maginn, "Isomolar semigrand ensemble molecular dynamics: Development and application to
liquid-liquid equilibria", Journal of Chemical Physics, 122, 054504 (2005). 

The iSG ensemble is a constant pressure, constant temperature, constant total # of molecules, and constant chemical potential
difference ensemble.  This ensemble is only valid for mixtures.  One mixture component must be chosen to be the reference component.
The chemical potential difference between component i and the reference component, Mu[i] - Mu[1], is commonly expressed using a quantity called
the fugacity fraction.  The chemical potential difference is related to the fugacity fraction by:

Mu[i] - Mu[1] = k * T * ln (X[i] / X[1]) + ideal gas terms,

where k is Boltzmann's constant, T is the temperature, X[i] is the fugacity fraction of component i, and the ideal gas terms are temperature-dependent
constants.  For an n-component system you must specify n fugacity fractions.  The fugacity fractions can vary only between 0 and 1 and must sum to one.
So for a 3 component system, for example, if X[1] = 0.3 and X[2] = 0.5, then X[3] must be equal to 0.2.

In order to maintain the mixture at the specified fugacity fractions, the simulation will attempt to transform molecules between different identities using
a dynamical transformation scheme.  As the iSGMD simulation progresses you will see the mixture composition change with time in order to reach a
chemical equilibrium with the specified fugacity fractions.  The mixture's composition, volume, and energy are written to disk in the iSGProperties file
as the simulation progresses.  This information can be used in conjuction with Gibbs-Duhem integration or histogram reweighiting to determine
phase coexistence.

To run an iSGMD simulation in Protomol2, you will need to supply the following input files:

1) an xyz coordinate file for all the atoms in the system.  This can be in either XYZ format or PDB format
2) either an initial temperature or an xyz velocity file for all the atoms in the system.  The velocity file can be in either XYZ format or PDB format
3) a structure file in the Protein Structure Format (PSF).  Please see the file ethane_methane_300_Pure_Trappe.psf for an example.  VERY 
   IMPORTANT!  You must specify the identity of each molecule/atom with an integer number in the far right column of the !NATOM section.  Number your
   molecule identities starting from zero (i.e. component 1 is 0, component 2 is 1, component 3 is 2, etc.).  Below is an example of this:
*******************
PSF

       5 !NTITLE
 REMARKS FILENAME=ethane_methane_300.psf by tmorrow
 REMARKS ProtoMol (built on Mar 14 2005 at 12:49:33)
 REMARKS This .psf file was created by PSFWriter
 REMARKS It was not manually assembled
 REMARKS Time : 150000, step : 150000.

     600 !NATOM
       1 MTET 1    BTW  U1   U1    0               15.0344           0
       2 MTET 1    BTW  U2   U2    0               15.0344           0
       3 MTET 2    BTW  U1   U1    0               15.0344           0
       4 MTET 2    BTW  U2   U2    0               15.0344           0
       5 MTET 3    BTW  U1   U1    0               16.0423           1
       6 MTET 3    BTW  U2   U2    0               15.0344           1
       7 MTET 4    BTW  U1   U1    0               15.0344           0
       8 MTET 4    BTW  U2   U2    0               15.0344           0
       9 MTET 5    BTW  U1   U1    0               15.0344           0
      10 MTET 5    BTW  U2   U2    0               15.0344           0
      11 MTET 6    BTW  U1   U1    0               16.0423           1
      12 MTET 6    BTW  U2   U2    0               15.0344           1
      13 MTET 7    BTW  U1   U1    0               15.0344           0
      14 MTET 7    BTW  U2   U2    0               15.0344           0
      15 MTET 8    BTW  U1   U1    0               16.0423           1
      16 MTET 8    BTW  U2   U2    0               15.0344           1
      17 MTET 9    BTW  U1   U1    0               15.0344           0
      18 MTET 9    BTW  U2   U2    0               15.0344           0
      19 MTET 10   BTW  U1   U1    0               15.0344           0
      20 MTET 10   BTW  U2   U2    0               15.0344           0
      .
      .
      .

      ^   ^   ^     ^   ^     ^    ^                 ^               ^
      a   b   c     d   e     f    g                 h               i
  a) Atom ID#
  b) molecule name (not very important)
  c) Molecule ID#
  d) residue name
  e) atom name
  f) atom force field type
  g) atomic charge (in this methane/ethane example none of the atoms are charged)
  h) atomic mass
  i) atom's identity, in this example ethane = 0 and methane = 1, so molecules 3, 6, and 8 are methane molecules and the rest are ethane molecules.
********************

4) If you desire to start a new simulation using the output of an old simulation as the starting point, you can specify an XSC file.  The XSC file
   contains the exact values of the simulation volume, volume velocity, thermostats and their velocities, chemostat and its velocity, and the ID# of the
   last molecule that was being transformed at the end of the simulation.  The XSC file will be created for you by ProtoMol2 when an iSGMD simulation is 
   finished, so you would normally not need to edit it yourself.
5) A ProtoMol2 configuration (.conf) file.  Please see ethane_methane_300.conf for the details.
6) A forcefield file (PAR file).  For each atom type in the PAR file you must specify one force parameter for each identity.  This means that for
   a 3-component system, for example, you must supply 3 Lennard-Jones parameters for each atom type, 3 force constants for every bond type, etc.  You must
   also specify the particular molecule identity to which each force parameter applies in the far right column of each section.  Please see
   the example file ethane_methane.par, a section of which is shown below:

****************
BONDS
!
!V(bond) = Kb(b - b0)**2
!
!Kb: kcal/mole/A**2
!b0: A
!
!atom type Kb          b0    identity
!
!-- single-topology bonds for iSGMD simulations
!-- for the UA model
U1   U2     500.0       1.540    0   ! UA alkane C-C bond -- ethane identity   (TraPPe model)
U1   U2     500.0       1.540    1   ! dummy-C bond -- methane identity        (TraPPe model)
.
.
.
NONBONDED nbxmod  5 atom cdiel shift vatom vdistance vswitch -
cutnb 14.0 ctofnb 12.0 ctonnb 10.0 eps 1.0 e14fac 1.0 wmin 1.5
                !adm jr., 5/08/91, suggested cutoff scheme
!
!V(Lennard-Jones) = Eps,i,j[(Rmin,i,j/ri,j)**12 - 2(Rmin,i,j/ri,j)**6]
!
!epsilon: kcal/mole, Eps,i,j = sqrt(eps,i * eps,j)
!Rmin/2: A, Rmin,i,j = Rmin/2,i + Rmin/2,j
!
!atom  ignored   epsilon      Rmin/2   ignored   eps,1-4       Rmin/2,1-4  identity
!
!-- single topology LJ parameters for iSGMD simulations
!-- for the UA model
U1     0.0       -0.19460      2.105      0.0     -0.0000        1.000      0 ! ethane methyl CH3      -- ethane identity  (TraPPe model)
U1     0.0       -0.29390      2.093      0.0     -0.0000        1.000      1 ! methane atom           -- methane identity (TraPPe model)
U2     0.0       -0.19460      2.105      0.0     -0.0000        1.000      0 ! ethane methyl CH3      -- ethane identity  (TraPPe model)
U2     0.0       -0.00000      2.105      0.0     -0.0000        1.000      1 ! dummy atom             -- methane identity (TraPPe model)
.
.
.

  In this example the bond is exactly the same in both the methane and ethane identities.  Even though it doesn't change you still have to list
  the force constant twice, once for methane and a second time for ethane.  In the Lennard-Jones section you can see that atom type U1 has ethane
  parameters for identity 0 and methane parameters for identity 1, while atom type U2 has ethane parameters for identity 0 and is a dummy atom (epsilon
  = 0) in identity 1.
********************

7) A transformation pathway file (TRANS file).  This file will define for the computer the exact pathway you want to use to accomplish the molecular
   transformation, including any intermediate states if necessary.  A detailed description of how to fill out each section of the TRANS file is given in
   the file ethane_methane.trans, and will not be repeated here.
