The Effective Fragment Potential Method

   The basic idea behind the effective fragment potential 
(EFP) method is to replace the chemically inert part of a 
system by EFPs, while performing a regular ab initio 
calculation on the chemically active part.  Here "inert" 
means that no covalent bond breaking process occurs.  This 
"spectator region" consists of one or more "fragments", 
which interact with the ab initio "active region" through 
non-bonded interactions, and so of course these EFP 
interactions affect the ab initio wavefunction.  The EFP 
particles can be closed shell or open shell (high spin 
ROHF) based potentials.  The "active region" can use nearly 
every kind of wavefunction available in GAMESS.

   A simple example of an active region might be a solute 
molecule, with a surrounding spectator region of solvent 
molecules represented by fragments.  Each discrete solvent 
molecule is represented by a single fragment potential, in 
marked contrast to continuum models for solvation.

   The quantum mechanical part of the system is entered in 
the $DATA group, along with an appropriate basis.  The EFPs 
defining the fragments are input by means of a $EFRAG  
group, and one or more $FRAGNAME groups describing each  
fragment's EFP.  These groups define non-bonded 
interactions between the ab initio system and the 
fragments, and also between the fragments.  The former 
interactions enter via one-electron operators in the ab 
initio Hamiltonian, while the latter interactions are 
treated by analytic functions.  The only electrons 
explicitly treated (with basis functions used to expand 
occupied orbitals) are those in the active region, so there 
are no new two electron terms.  Thus the use of EFPs leads 
to significant time savings, compared to full ab initio 
calculations on the same system. 
   There are two types of EFP available in GAMESS, EFP1 and 
EFP2.  EFP1, the original method, employs a fitted 
repulsive potential.  EFP1 is primarily used to model water 
molecules to study aqueous solvation effects, at the 
RHF/DZP or DFT/DZP (specifically, B3LYP) levels, see 
references 1-3 and 26, respectively.  EFP2 is a more 
general method that is applicable to any species, including 
water, and its repulsive potential is obtained from first 
principles.  EFP2 has been extended to include other 
effects as well, such as charge transfer and dispersion.  
EFP2 forms the basis of the covalent EFP method described 
below for modeling enzymes, see reference 14.

   Parallelization of the EFP1 and EFP2 models is described 
in reference 32.

   MD simulations with EFP are described in reference 31.

   The ab initio/EFP1, or pure EFP system can be wrapped in 
a Polarizable Continuum Model, see references 23, 43, and 

terms in an EFP

   The non-bonded interactions currently implemented are:

1) Coulomb interaction.  The charge distribution of the 
fragments is represented by an arbitrary number of charges, 
dipoles, quadrupoles, and octopoles, which interact with 
the ab initio hamiltonian as well as with multipoles on 
other fragments (see reference 2 and 18).  It is possible 
to use a screening term that accounts for the charge 
penetration (reference 17 and 42).  This screening term is 
automatically included for EFP1.  Typically the multipole 
expansion points are located on atomic nuclei and at bond 

2) Dipole polarizability.  An arbitrary number of dipole 
polarizability tensors can be used to calculate the induced 
dipole on a fragment due to the electric field of the ab 
initio system as well as all the other fragments.  These 
induced dipoles interact with the ab initio system as well 
as the other EFPs, in turn changing their electric fields.  
All induced dipoles are therefore iterated to self-
consistency.  Typically the polarizability tensors are 
located at the centroid of charge of each localized orbital 
of a fragment.  See reference 41.

3) Repulsive potential.  Two different forms are used in 
EFP1: one for ab initio-EFP repulsion and one for EFP-EFP 
repulsion.  The form of the potentials is empirical, and 
consists of distributed Gaussian or exponential functions, 
respectively.  The primary contribution to the repulsion is 
the quantum mechanical exchange repulsion, but the fitting 
technique used to develop this term also includes the 
effects of charge transfer.  Typically these fitted 
potentials are located on each atomic nucleus within the 
fragment (see reference 3).  In EFP2, polarization energies 
can also be augmented by screening terms, analogous to the 
electrostatic screening, to prevent "polarization collapse" 
(MS in preparation)

For EFP2, the third term is divided into separate analytic 
formulae for different physical interactions:
    a) exchange repulsion
    b) dispersion
    c) charge transfer
A summary of EFP2, and its contrast to EFP1 can be found in 
reference 18 and 44.  The repulsive potential for EFP2 is 
based on an overlap expansion using localized molecular 
orbitals, as described in references 5, 6, and 9.  
Dispersion energy is described in reference 34, and charge 
transfer in reference 39 (which supercedes reference 22's 

    EFP2 potentials have no fitted parameters, and can be 
automatically generated during a RUNTYP=MAKEFP job, as 
described below.

constructing an EFP1

   RUNTYP=MOROKUMA assists in the decomposition of inter-
molecular interaction energies into electrostatic, 
polarization, charge transfer, and exchange repulsion 
contributions.  This is very useful in developing EFPs 
since potential problems can be attributed to a particular 
term by comparison to these energy components for a 
particular system.

   A molecular multipole expansion can be obtained using 
$ELMOM.  A distributed multipole expansion can be obtained 
by either a Mulliken-like partitioning of the density 
(using $STONE) or by using localized molecular orbitals 
($LOCAL: DIPDCM and QADDCM).  The dipole polarizability 
tensor can be obtained during a Hessian run ($CPHF), and a 
distributed LMO polarizability expression is also available 

   In EFP1, the repulsive potential is derived by fitting 
the difference between ab initio computed intermolecular 
interaction energies, and the form used for Coulomb and 
polarizability interactions.  This difference is obtained 
at a large number of different interaction geometries, and 
is then fitted.  Thus, the repulsive term is implicitly a 
function of the choices made in representing the Coulomb 
and polarizability terms.  Note that GAMESS currently does 
not provide a way to obtain these EFP1 repulsive potential.

   Since a user cannot generate all of the EFP1 terms 
necessary to define a new $FRAGNAME group using GAMESS, in 
practice the usage of EFP1 is limited to the internally 
stored H2ORHF or H2ODFT potentials mentioned below.

constructing an EFP2

   As noted above, the repulsive potential for EFP2 is 
derived from a localized orbital overlap expansion.  It is 
generally recommended that one use at least a double zeta 
plus diffuse plus polarization basis set, e.g. 6-31++G(d,p) 
to generate the EFP2 repulsive potential.  However, it has 
been observed that 6-31G(d) works reasonably well due to a 
fortuitous cancellation of errors.  The EFP2 potential for 
any molecule can be generated as follows:

(a) Choose a basis set and geometry for the molecule of 
interest.  The geometry is ordinarily optimized at your 
choice of Hartree-Fock/MP2/CCSD(T), with your chosen basis 
set, but this is not a requirement.  It is good to recall, 
however, that EFP internal geometries are fixed, so it is 
important to give some thought to the chosen geometry.

(b) Perform a RUNTYP=MAKEFP run for the chosen molecule 
using the chosen geometry in $DATA and the chosen basis set 
in $BASIS.  This will generate the entire EFP2 potential in 
the run's .efp file.  The only user-defined variable that 
must be filled in is changing the FRAGNAME's group name, to 
$C2H5OH or $DMSO, etc.  This step can use RHF or ROHF to 
describe the electronic structure of the system.

(c) Transfer the entire fragment potential for the molecule 
to any input file in which this fragment is to be used. 
Since the internal geometry of an EFP is fixed, one need 
only specify the first three atoms of any fragment in order 
to position them in $EFRAG.  Coordinates of any other atoms 
in the rigid fragment will be automatically determined by 
the program.

If the EFP contains less than three atoms, you can still 
generate a fragment potential.  After a normal MAKEFP run, 
add dummy atoms (e.g. in the X and/or Y directions) with 
zero nuclear charges, and add corresponding dummy bond 
midpoints too.  Carefully insert zero entries in the 
multipole sections, and in the electrostatic screening 
sections, for each such dummy point, but don't add data to 
any other kind of EFP term such as polarizability.  This 
trick gives the necessary 3 points for use in $EFRAG groups 
to specify "rotational" positions of fragments.
current limitations

1. For EFP1, the energy and energy gradient are programmed, 
which permits RUNTYP=ENERGY, GRADIENT, and numerical 
HESSIAN.  The necessary programing to use the EFP gradients 
to move on the potential surface are programmed for 
gradient based potential surface explorations such as DRC 
are not yet available.  Finally, RUNTYP=PROP is also 

For EFP2, the gradient terms for ab initio-EFP interactions 
have not yet been coded, so geometry optimizations are only 
sensible for a COORD=FRAGONLY run; that is, a run in which 
only EFP2 fragments are present.

2. The ab initio part of the system must be treated with 
RHF, ROHF, UHF, the open shell SCF wavefunctions permitted 
by the GVB code, or MCSCF.  DFT analogs of RHF, ROHF, and 
UHF may also be used.  Correlated methods such as MP2 and 
CI should not be used.

3. EFPs can move relative to the ab initio system and 
relative to each other, but the internal structure of an 
EFP is frozen.

4. The boundary between the ab initio system and EFP1's 
must not be placed across a chemical bond.  However, see 
the discussion below regarding covalent bonds.

5. Calculations must be done in C1 symmetry at present.

6. Reorientation of the fragments and ab initio system is 
not well coordinated.  If you are giving Cartesian 
coordinates for the fragments (COORD=CART in $EFRAG), be 
sure to use $CONTRL's COORD=UNIQUE option so that the ab 
initio molecule is not reoriented.

7. If you need IR intensities, you have to use NVIB=2.  The 
potential surface is usually very soft for EFP motions, and 
double differenced Hessians should usually be obtained.

practical hints for using EFPs

   At the present time, we have only two internally stored 
EFP potentials suitable for general use.  These model 
water, using the fragment name H2ORHF or H2ODFT.  The 
H2ORHF numerical parameters are improved values over the 
values which were presented and used in reference 2, and 
they also include the improved EFP-EFP repulsive term 
defined in reference 3.  The H2ORHF water EFP was derived 
from RHF/DH(d,p) computations on the water dimer system.  
When you use it, therefore, the ab initio part of your 
system should be treated at the SCF level, using a basis 
set of the same quality (ideally DH(d,p), but probably 
other DZP sets such as 6-31G(d,p) will give good results as 
well).  Use of better basis sets than DZP with this water 
EFP has not been tested.  Similarly, H2ODFT was developed 
using B3LYP/DZP water wavefunctions, so this should be used 
(rather than H2ORHF) if you are using DFT to treat the 
solute.  Since H2ODFT water parameters are obtained from a 
correlated calculation, they can also be used when the 
solute is treated by MP2.

   As noted, effective fragments have frozen internal 
geometries, and therefore only translate and rotate with 
respect to the ab initio region.  An EFP's frozen 
coordinates are positioned to the desired location(s) in 
$EFRAG as follows:
  a) the corresponding points are found in $FRAGNAME.
  b) Point -1- in $EFRAG and its FRAGNAME equivalent are 
     made to coincide.
  c) The vector connecting -1- and -2- is aligned with
     the corresponding vector connecting FRAGNAME points.
  d) The plane defined by -1-, -2-, and -3- is made to
     coincide with the corresponding FRAGNAME plane.
Therefore the 3 points in $EFRAG define only the relative 
position of the EFP, and not its internal structure. So, if 
the "internal structure" given by points in $EFRAG differs 
from the true values in $FRAGNAME, then the order in which 
the points are given in $EFRAG can affect the positioning 
of the fragment.  It may be easier to input water EFPs if 
you use the Z-matrix style to define them, because then you 
can ensure you use the actual frozen geometry in your 
$EFRAG.  Note that the H2ORHF EFP uses the frozen geometry 
r(OH)=0.9438636, a(HOH)=106.70327, and the names of its 3 
fragment points are ZO1, ZH2, ZH3.

                         * * *

   Building a large cluster of EFP particles by hand can be 
tedious.  The RUNTYP=GLOBOP program described below has an 
option for constructing dense clusters.  The method tries 
to place particles near the origin, but not colliding with 
other EFP particles already placed there, so that the 
clusters grow outwards from the center.  Here are some 
    a) place 100 water molecules, all with the same coords
       in $EFRAG.  This will build up a droplet of water
       with particles close together, but not on top of
       each other, with various orientations.
    b) place 16 waters (same coords, all first) followed by
       16 methanols (also sharing their same coords, after
       all waters).  A 50-50 mixture of 32 molecules will
       be created, if you choose the default of picking the
       particles randomly from the initial list of 32.
    c) to solvate a solute, add the solute in the $DATA                         
       group at or near the origin.  Add the solvent
       molecules near by (same coords is ok), and run
       the globop run with RNDINI as demonstrated below.
       (optional, add MCTYP=3 to $GLOBOP input)

Example, allowing the random cluster to have 20 geometry 
optimization steps:
 $contrl runtyp=globop coord=fragonly $end
 $globop rndini=.true. riord=rand mcmin=.true.
         mctyp=4 nblock=0 $end
 $statpt nstep=20 $end
O1   -2.8091763203009  -2.1942725073400  -0.2722207394107
H2   -2.3676165499399  -1.6856118830379  -0.9334073942601
H3   -2.1441965467625  -2.5006167998896   0.3234583094693
...repeat this 15 more times...
O1   4.9515153249    .4286994611   1.3368662306
H2   5.3392575544    .1717424606   3.0555957053
C3   6.2191743799   2.5592349960    .4064662379
H4   5.7024200977   2.7548960076  -1.5604873643
H5   5.6658856694   4.2696553371   1.4008542042
H6   8.2588049857   2.3458272252    .5282762681
...repeat 15 more times...
...give a full EFP2 potential for water...
...give a full EFP2 potential for methanol...
Note that the random cluster generation now proceeds into a 
full Monte Carlo simulation.

                         * * *

   The translations and rotations of EFPs with respect to 
the ab initio system and one another are automatically 
quite soft degrees of freedom.  After all, the EFP model is 
meant to handle weak interactions!  Therefore the 
satisfactory location of structures on these flat surfaces 
will require use of a tight convergence on the gradient: 
OPTTOL=0.00001 in the $STATPT group.

   The effect of a bulk continuum surrounding the solute 
plus EFP waters can be obtained by using the PCM model, see 
reference 23 and 43.  To do this, simply add a $PCM group 
to your input, in addition to the $EFRAG.  The simultaneous 
use of EFP and PCM allows for gradients, so geometry 
optimization can be performed.

global optimization

    If there are a large number of particles to move (EFP 
and/or FMO and/or atom groups), it is difficult to locate 
the lowest energy structures by hand.  Typically these are 
numerous, and one would like to have a number of them, not 
just the very lowest energy.  The RUNTYP of GLOBOP contains 
a Monte Carlo procedure to generate a random set of 
starting structures to look for those with the lowest 
energy at a single temperature.  If desired, a simulated 
annealing protocol to cool the temperature may be used.  
These two procedures may be combined with a local minimum 
search, at some or all of the randomly generated 
structures.  The local minimum search is controlled by the 
usual geometry optimizer, namely $STATPT input, and thus 
permits the optimization of any ab initio atoms.

    The Monte Carlo procedure by default uses a Metropolis 
algorithm to move just one of the fragments.  The method of 
Parks to move all fragments simultaneously is also allowed.

    The present program was used to optimize the structure 
of water clusters.  Let us consider the case of the twelve 
water cluster, for which the following ten structures were 
published by Day, Pachter, Gordon, and Merrill:
   1. (D2d)2     -0.170209     6. (D2d)(C2)  -0.167796
   2. (D2d)(S4)  -0.169933     7. S6         -0.167761
   3. (S4)2      -0.169724     8. cage b     -0.167307
   4. D3         -0.168289     9. cage a     -0.167284
   5. (C1c)(Cs)  -0.167930    10. (C1c)(C1c) -0.167261
A test input using Metropolis style Monte Carlo to examine 
300 geometries at each temperature value, using simulated 
annealing cooling from 200 to 50 degrees, and with local 
minimization every 10 structures was run ten times.  Each 
run sampled about 7000 geometries.  One simulation found 
structure 2, while two of the runs found structure 3.  The 
other seven runs located structures with energy values in 
the range -0.163 to -0.164.  In all cases the runs began 
with the same initial geometry, but produced different 
results due to the random number generation used in the 
Monte Carlo.  Clearly one must try a lot of simulations to 
be confident about having found most of the low energy 
structures.  In particular, it is good to try more than one 
initial structure, unlike what was done in this test.

   Ab initio atoms can be addressed using FMO, either in 
multiple fragments, or perhaps a single large fragment.  
Alternatively, ab initio atoms can be put into groups and 
used directly in globop, which for small systems has a 
lower overhead than FMO.  In the case of large molecules                        
separated into multiple fragments, the keywords NPRBND, 
PRSEP, IBNDS, and INDEP are applicable.  These specify the 
atoms in each set of fragments or groups whose bond is cut 
in the fragmentation process.  The paired atoms are 
constrained during the Monte Carlo procedure to ensure that 
the bond is not spacially broken.  In the case where a 
fragment that is being translated or rotated is paired with 
two or more fragments, the movement is repeated on all 
attached fragments, after randomly choosing which pair is 
the starting point.  For example, given a molecule split 
into five fragments such that:


where A,B,C,D,E are the fragments.  If C is chosen for a 
translation, either B or D will be randomly chosen to be 
the starting pair.  When B is chosen as the starting pair, 
C, D, and E will all be translated by the same amount:


which maintains the relative position of C, D, and E.  
Setting INDEP=1 will not propagate the translation:


So that only C is moved.

The same approach is used for rotations.  Since a small 
translation or rotation can result in a significant change 
in the total system, it is advised that case be taken when 
using solvent molecules and to the size of boundary 
conditions.  If a propagated movement moves a fragment 
outside the boundary, a warning will be printed and the 
step will be discarded as a proximity alert.  Also, the 
pair binding is not implemented for RNDINI=.TRUE. To 
initialize a set of solvent molecules around pair bonded 
fragments, include the pair bonded fragments in IFXFMO.

The Metrpolis Monte Carlo procedure involves the movement 
of groups that are internally rigid.  To introduce some 
internal flexibility for FMO and ab initio groups, a 
secondary Monte Carlo search where the entire system is 
held rigid while the atoms in one group are moved is 
implemented.  The secondary Monte Carlo occurs when a FMO 
or ab initio group is translated and occurs for that group.  
The lowest energy internal configuration for the secondary 
Monte Carlo is used when evaluating the step of the primary 
Monte Carlo search.  The temperature at which the secondary 
Monte Carlo is used in the case of simulated annealing is 
set by SMTEMP and the number of steps in each secondary 
search is given by NSMTP. To turn on this feature, set the 
values of SMTEMP and NSMTP to non-zero values.

Monte Carlo references:
  N.Metropolis, A.Rosenbluth, A.Teller
      J.Chem.Phys. 21, 1087(1953).
  G.T.Parks  Nucl.Technol. 89, 233(1990).
Monte Carlo with local minimization:
  Z.Li, H.A.Scheraga
      Proc.Nat.Acad.Sci. USA  84, 6611(1987).
Simulated annealing reference:
  S.Kirkpatrick, C.D.Gelatt, M.P.Vecci
      Science 220, 671(1983).

The present program is described in reference 15.  It is 
patterned on the work of
   D.J.Wales, M.P.Hodges Chem.Phys.Lett. 286, 65-72 (1998).

QM/MM across covalent bonds

    Recent work by Visvaldas Kairys and Jan Jensen has made 
it possible to extend the EFP methodology beyond the simple 
solute/solvent case described above.  When there is a 
covalent bond between the portion of the system to be 
modeled by quantum mechanics, and the portion which is to 
be treated by EFP multipole and polarizability terms, an 
additional layer is needed in the model.  The covalent 
linkage is not so simple as the interactions between closed 
shell solute and solvent molecules.  The "buffer zone" 
between the quantum mechanics and the EFP consists of 
frozen nuclei, and frozen localized orbitals, so that the 
quantum mechanical region sees a orbital representation of 
the closest particles, and multipoles etc. beyond that. 
Since the orbitals in the buffer zone are frozen, it need 
extend only over a few atoms in order to keep the orbitals 
in the fully optimized quantum region within that region.

    The general outline of this kind of computation is as 
    a) a full quantum mechanics computation on a system
       containing the quantum region, the buffer region,
       and a few atoms into the EFP region, to obtain the
       frozen localized orbitals in the buffer zone.
       This is called the "truncation run".
    b) a full quantum mechanics computation on a system
       with all quantum region atoms removed, and with
       the frozen localized orbitals in the buffer zone.
       The necessary multipole and polarizability data
       to construct the EFP that will describes the EFP
       region will be extracted from the wavefunction.
       This is called the "MAKEFP run".  It is possible
       to use several such runs if the total EFP region
       is quite large.
    c) The intended QM/MM run(s), after combining the
       information from these first two types of runs.

    As an example, consider a protonated lysine residue 
which one might want to consider quantum mechanically in a 
protein whose larger parts are to be treated with an EFP.  
The protonated lysine is

  +                             /

The bonds which you see drawn show how the molecule is 
partitioned between the quantum mechanical side chain, a 
CH2CH group in the buffer zone, and eventually two 
different EFPs may be substituted in the area of the NH2 
and COOH groups to form the protein backbone.

   The "truncation run" will be on the entire system as you 
see it, with the 13 atoms in the side chain first in $DATA, 
the 5 atoms in the buffer zone next in $DATA, and the 
simplified EFP region at the end.  This run will compute 
the full quantum wavefunction by RUNTYP=ENERGY, followed by 
the calculation of localized orbitals, and then truncation 
of the localized orbitals that are found in the buffer zone 
so that they contain no contribution from AOs outside the 
buffer zone. The key input groups for this run are
 $truncn doproj=.true. plain=.true. natab=13 natbf=5 $end
This will generate a total of 6 localized molecular 
orbitals in the buffer zone (one CC, three CH, two 1s inner 
shells), expanded in terms of atomic orbitals located only 
on those atoms.

    The truncation run prepares template input files for 
the next run, including adjustments of nuclear charges at 
boundaries, etc.

    The "MAKEFP" run drops all 13 atoms in the quantum 
region, and uses the frozen orbitals just prepared to 
obtain a wavefunction for the EFP region.  The carbon atom 
in the buffer zone that is connected to the now absent QM 
region will have its nuclear charge changed from 6 to 5 to 
account for a missing electron.  The key input for this 
RUNTYP=MAKEFP job is the six orbitals in $VEC, plus the 
 $guess guess=huckel insorb=6 $end
 $mofrz frz=.true. ifrz(1)=1,2,3,4,5,6 $end

which will cause the wavefunction optimization for the 
remaining atoms to optimize orbitals only in the NH2 and 
COOH pieces.  After this wavefunction is found, the run 
extracts the EFP information needed for the QM/MM third 
run(s).  This means running the Stone analysis for 
distributed multipoles, and obtaining a polarizability  
tensor for each localized orbital in the EFP region.

    The QM/MM run might be RUNTYP=OPTIMIZE, etc. depending 
on what you want to do with the quantum atoms, and its 
$DATA group will contain both the 13 fully optimized atoms, 
and the 5 buffer atoms, and a basis set will exist on both 
sets of atoms.  The carbon atom in the buffer zone that 
borders the EFP region will have its nuclear charge set to 
4 since now two bonding electrons to the EFP region are 
lost.  $VEC input will provide the six frozen orbitals in 
the buffer zone.  The EFP atoms are defined in a fragment 
potential group.
    The QM/MM run could use RHF or ROHF wavefunctions, to 
geometry optimize the locations of the quantum atoms (but 
not of course the frozen buffer zone or the EFP piece).  It 
could remove the proton to compute the proton affinity at 
that terminal nitrogen, hunt for transition states, and so 
on.  Presently the gradient for GVB and MCSCF is not quite 
right, so their use is discouraged.

    Input to control the QM/MM preparation is $TRUNCN and 
$MOFRZ groups.  There are a number of other parameters in 
various groups, namely QMMMBUF in $STONE, MOIDON and POLNUM 
in $LOCAL, NBUFFMO in $EFRAG, and INSORB in $GUESS that are 
relevant to this kind of computation.  For RUNTYP=MAKEFP, 
the biggest choices are LOCAL=RUEDENBRG vs. BOYS, and 
POLNUM in $LOCAL, otherwise this is pretty much a standard 
RUNTYP=ENERGY input file.

    Source code distributions of GAMESS contain a directory 
named ~/gamess/tools/efp, which has various tools for EFP 
manipulation in it, described in file readme.1st.  A full 
input file for the protonated lysine molecule is included, 
with instructions about how to proceed to the next steps. 
Tips on more specialized input possibilities are appended 
to the file readme.1st.

Simpler potentials

   Since the EFP model's electrostatics is a set of 
distributed multipoles (monopole to octopole) and 
distributed polarizabilities (dipole), it is possible to 
generate some water potentials found in the literature by 
setting many EFP terms to zero.  It is also necessary to 
provide a Lennard-Jones 6-12 repulsive potential, and then 
make a choice to follow the EFP1 type formula for QM/EFP 
repulsion.  Accordingly, EFP1 type calculations can be made 
with the following water potentials,
The Wikipedia page
defines the first four of these, which are not polarizable 
potentials.  The same web site references the primary 
literature, so that is not repeated here.  POL5P is a 
polarizable potential, with parameters given by
    D.Si and H.Li  J.Chem.Phys. 133, 144112/1-8(2010)

   The first paper is more descriptive, while the second 
presents a very detailed derivation of the EFP1 method.  
Reference 18 is an overview article on EFP2.  Reference 44 
is the most recent review.

   The model development papers are: 1, 2, 3, 5, 6, 9, 14, 
17, 18, 22, 23, 26, 31, 32, 34, 39, 41, 42, 43, 44, 46, 50, 
51, 55, 57, 58.

1. "Effective fragment method for modeling intermolecular 
    hydrogen bonding effects on quantum mechanical 
    J.H.Jensen, P.N.Day, M.S.Gordon, H.Basch, D.Cohen, 
    D.R.Garmer, M.Krauss, W.J.Stevens in "Modeling the
    Hydrogen Bond" (D.A. Smith, ed.) ACS Symposium Series 
    569, 1994, pp 139-151.
2. "An effective fragment method for modeling solvent
    effects in quantum mechanical calculations".
    P.N.Day, J.H.Jensen, M.S.Gordon, S.P.Webb,
    W.J.Stevens, M.Krauss, D.Garmer, H.Basch, D.Cohen
    J.Chem.Phys. 105, 1968-1986(1996).
3. "The effective fragment model for solvation: internal
    rotation in formamide"
    W.Chen, M.S.Gordon, J.Chem.Phys., 105, 11081-90(1996)
4. "Transphosphorylation catalyzed by ribonuclease A: 
    Computational study using ab initio EFPs"
    B.D.Wladkowski, M. Krauss, W.J.Stevens
    J.Am.Chem.Soc. 117, 10537-10545(1995)
5. "Modeling intermolecular exchange integrals between
    nonorthogonal orbitals"
    J.H.Jensen  J.Chem.Phys. 104, 7795-7796(1996)
6. "An approximate formula for the intermolecular Pauli
    repulsion between closed shell molecules"
    J.H.Jensen, M.S.Gordon  Mol.Phys. 89, 1313-1325(1996)
7. "A study of aqueous glutamic acid using the effective
    fragment potential model"
    P.N.Day, R.Pachter  J.Chem.Phys. 107, 2990-9(1997)
8. "Solvation and the excited states of formamide"
    M.Krauss, S.P.Webb  J.Chem.Phys. 107, 5771-5(1997)
9. "An approximate formula for the intermolecular Pauli
    repulsion between closed shell molecules.  Application
    to the effective fragment potential method"
    J.H.Jensen, M.S.Gordon
    J.Chem.Phys. 108, 4772-4782(1998)
10. "Study of small water clusters using the effective
    fragment potential method"
    G.N.Merrill, M.S.Gordon J.Phys.Chem.A 102, 2650-7(1998)
11. "Solvation of the Menshutkin Reaction: A Rigourous
    test of the Effective Fragement Model"
    S.P.Webb, M.S.Gordon  J.Phys.Chem.A  103, 1265-73(1999)
12. "Evaluation of the charge penetration energy between
    nonorthogonal molecular orbitals using the Spherical
    Gaussian Overlap approximation"
    V.Kairys, J.H.Jensen
    Chem.Phys.Lett. 315, 140-144(1999)
13. "Solvation of Sodium Chloride: EFP study of NaCl(H2O)n"
    C.P.Petersen, M.S.Gordon
    J.Phys.Chem.A 103, 4162-6(1999)
14. "QM/MM boundaries across covalent bonds: frozen LMO
    based approach for the Effective Fragment Potential
    V.Kairys, J.H.Jensen  J.Phys.Chem.A  104, 6656-65(2000)
15. "A study of water clusters using the effective fragment
    potential and Monte Carlo simulated annealing"
    P.N.Day, R.Pachter, M.S.Gordon, G.N.Merrill
    J.Chem.Phys. 112, 2063-73(2000)
16. "A combined discrete/continuum solvation model:
    Application to glycine"  P.Bandyopadhyay, M.S.Gordon
    J.Chem.Phys. 113, 1104-9(2000)
17. "Evaluation of charge penetration between distributed
    multipolar expansions"
    M.A.Freitag, M.S.Gordon, J.H.Jensen, W.J.Stevens
    J.Chem.Phys. 112, 7300-7306(2000)
18. "The Effective Fragment Potential Method: a QM-based MM
    approach to modeling environmental effects in
    M.S.Gordon, M.A.Freitag, P.Bandyopadhyay, J.H.Jensen,
    V.Kairys, W.J.Stevens J.Phys.Chem.A  105, 293-307(2001)
19. "Accurate Intraprotein Electrostatics derived from
    first principles: EFP study of proton affinities of
    lysine 55 and tyrosine 20 in Turkey Ovomucoid"
    R.M.Minikis, V.Kairys, J.H.Jensen
    J.Phys.Chem.A  105, 3829-3837(2001)
20. "Active site structure & mechanism of Human Glyoxalase"
    U.Richter, M.Krauss J.Am.Chem.Soc. 123, 6973-6982(2001)
21. "Solvent effect on the global and atomic DFT-based
    reactivity descriptors using the EFP model. Solvation 
    of ammonia."  R.Balawender, B.Safi, P.Geerlings
    J.Phys.Chem.A  105, 6703-6710(2001)
22. "Intermolecular exchange-induction and charge transfer:
    Derivation of approximate formulas using nonorthogonal
    localized molecular orbitals."
    J.H.Jensen J.Chem.Phys. 114, 8775-8783(2001)
23. "An integrated effective fragment-polarizable continuum
    approach to solvation: Theory & application to glycine"
    P.Bandyopadhyay, M.S.Gordon, B.Mennucci, J.Tomasi
    J.Chem.Phys. 116, 5023-5032(2002)
24. "The prediction of protein pKa's using QM/MM: the pKa
    of Lysine 55 in turkey ovomucoid third domain"
    H.Li, A.W.Hains, J.E.Everts, A.D.Robertson, J.H.Jensen
    J.Phys.Chem.B 106, 3486-3494(2002)
25. "Computational studies of aliphatic amine basicity"
    D.C.Caskey, R.Damrauer, D.McGoff
    J.Org.Chem. 67, 5098-5105(2002)
26. "Density Functional Theory based Effective Fragment
    Potential" I.Adamovic, M.A.Freitag, M.S.Gordon
    J.Chem.Phys. 118, 6725-6732(2003)
27. "Intraprotein electrostatics derived from first
    principles: Divid-and-conquer approaches for QM/MM
    calculations"  P.A.Molina, H.Li, J.H.Jensen
    J.Comput.Chem. 24, 1971-1979(2003)
28. "Formation of alkali metal/alkaline earth cation water
    clusters, M(H2O)1-6, M=Li+, K+, Mg+2, Ca+2: an
    effective fragment potential caase study"
    G.N.Merrill, S.P.Webb, D.B.Bivin
    J.Phys.Chem.A  107, 386-396(2003)
29. "Anion-water clusters A-(H2O)1-6, A=OH, F, SH, Cl, and
    Br. An effective fragment potential test case"
    G.N.Merrill, S.P.Webb
    J.Phys.Chem.A  107,7852-7860(2003)
30. "The application of the Effective Fragment Potential to
    molecular anion solvation: a study of ten oxyanion-
    water clusters, A-(H2O)1-4"
    G.N.Merrill, S.P.Webb  J.Phys.Chem.A 108, 833-839(2004)
31. "The effective fragment potential: small clusters and
    radial distribution functions"
    H.M.Netzloff, M.S.Gordon J.Chem.Phys. 121, 2711-4(2004)
32. "Fast fragments: the development of a parallel
    effective fragment potential method"
    H.M.Netzloff, M.S.Gordon
    J.Comput.Chem. 25, 1926-36(2004)
33. "Theoretical investigations of acetylcholine (Ach) and
    acetylthiocholine (ATCh) using ab initio and effective
    fragment potential methods"
    J.Song, M.S.Gordon, C.A.Deakyne, W.Zheng
    J.Phys.Chem.A 108, 11419-11432(2004)
34. "Dynamic polarizability, dispersion coefficient C6, and
    dispersion energy in the effective fragment potential
    I.Adamovic, M.S.Gordon  Mol.Phys. 103, 379-387(2005)
35. "Solvent effects on the SN2 reaction: Application of
    the density functional theory-based effective fragment
    potential method"
    I.Adamovic, M.S.Gordon J.Phys.Chem.A 109, 1629-36(2005)
36. "Theoretical study of the solvation of fluorine and
    chlorine anions by water"
    D.D.Kemp, M.S.Gordon  J.Phys.Chem.A 109, 7688-99(2005)
37. "Modeling styrene-styrene interactions"
    I.Adamovic, H.Li, M.H.Lamm, M.S.Gordon
    J.Phys.Chem.A 110, 519-525(2006)
38. "Methanol-water mixtures: a microsolvation study using
    the Effective Fragment Potential method"
    I.Adamovic, M.S.Gordon
    J.Phys.Chem.A 110, 10267-10273(2006)
39. "Charge transfer interaction in the effective fragment
    potential method"  H.Li, M.S.Gordon, J.H.Jensen
    J.Chem.Phys. 124, 214108/1-16(2006)
40. "Incremental solvation of nonionized and zwitterionic
    C.M.Aikens, M.S.Gordon
    J.Am.Chem.Soc. 128, 12835-12850(2006)
41. "Gradients of the polarization energy in the Effective
    Fragment Potential method"
    H.Li, H.M.Netzloff, M.S.Gordon
    J.Chem.Phys. 125, 194103/1-9(2006)
42. "Electrostatic energy in the Effective Fragment
    Potential method: Theory and application to benzene
    L.V.Slipchenko, M.S.Gordon
    J.Comput.Chem. 28, 276-291(2007)
43. "Polarization energy gradients in combined Quantum
    Mechanics, Effective Fragment Potential, and
    Polarizable Continuum Model Calculations"
    H.Li, M.S.Gordon  J.Chem.Phys. 126, 124112/1-10(2007)
44. "The Effective Fragment Potential: a general method
    for predicting intermolecular interactions"
    M.S.Gordon, L.V.Slipchenko, H.Li, J.H.Jensen
    Annual Reports in Computational Chemistry, Volume 3,
    pp 177-193 (2007).
45. "An Interpretation of the Enhancement of the Water
    Dipole Moment Due to the Presence of Other Water
    D.D.Kemp, M.S.Gordon
    J.Phys.Chem.A  112, 4885-4894(2008)
46. "Solvent effects on optical properties of molecules: a
    combined time-dependent density functional/effective
    fragment potential approach"
    S.Yoo, F.Zahariev, S.Sok, M.S.Gordon
    J.Chem.Phys. 129, 144112/1-8(2008)
47. "Modeling pi-pi interactions with the effective
    fragment potential method: The benzene dimer and
    T.Smith, L.V.Slipchenko, M.S.Gordon
    J.Phys.Chem.A  112, 5286-5294(2008)
48. "Water-benzene interactions: An effective fragment
    potential and correlated quantum chemistry study"
    L.V.Slipchenko, M.S.Gordon
    J.Phys.Chem.A 113, 2092-2102(2009)
49. "Ab initio QM/MM excited-state molecular dynamics study
    of Coumarin 151 in water solution"
    D.Kina, P.Arora, A.Nakayama, T.Noro, M.S.Gordon,
    T.Taketsugu  Int.J.Quantum Chem. 109, 2308-2318(2009)
50. "Damping functions in the effective fragment potential
    L.V.Slipchenko, M.S.Gordon
    Mol.Phys. 197, 999-1016 (2009)
51. "A combined effective fragment potential-fragment
    molecular orbital method. 1. the energy expression"
    T.Nagata, D.G.Fedorov, K.Kitaura, M.S.Gordon
    J.Chem.Phys. 131, 024101/1-12(2009)
52. "Alanine: then there was water"
    J.M.Mullin, M.S.Gordon
    J.Phys.Chem.B 113, 8657-8669(2009)
53. "Water and Alanine: from puddles(32) to ponds(49)"
    J.M.Mullin, M.S.Gordon
    J.Phys.Chem.B 113, 14413-14420(2009)
54. "Structure of large nitrate-water clusters at ambient
    temperatures: simulations with effective fragment
    potentials and force fields with implications for
    atmospheric chemistry"
    Y.Miller, J.L.Thoman, D.D.Kemp, B.J.Finlayson-Pitts,
    M.S.Gordon, D.J.Tobias, R.B.Gerber
    J.Phys.Chem.A  113, 12805-12814(2009)
55. "Quantum mechanical/molecular mechanical/continuum
    style solvation model: linear response theory,
    variational treatment, and nuclear gradients"
    H.Li  J.Chem.Phys. 131, 184103/1-8(2009)
56. "Aqueous solvation of bihalide anions"
    D.D.Kemp, M.S.Gordon
    J.Phys.Chem.A 114, 1298-1303(2010)
57. "Exchange repulsion between effective fragment
    potentials and ab initio molecules"
    D.D.Kemp, J.M.Rintelman, M.S.Gordon, J.H.Jensen
    Theoret.Chem.Acc. 125, 481-491(2010).
58. Modeling Solvent Effects on Electronic Excited States
    A.DeFusco, N.Minezawa, L.V.Slipchenko, F.Zahariev,
    M.S.Gordon  J.Phys.Chem.Lett. 2, 2184-2192(2011)

created on 7/7/2017