update to USER-FEP docs and examples
This commit is contained in:
@ -119,7 +119,7 @@ USER-CG-CMM, coarse-graining model, Axel Kohlmeyer (Temple U), "pair_style lj/sd
|
||||
USER-COLVARS, collective variables, Fiorin & Henin & Kohlmeyer (3), "fix colvars"_fix_colvars.html, USER/colvars, "colvars"_colvars, lib/colvars
|
||||
USER-CUDA, NVIDIA GPU styles, Christian Trott (U Tech Ilmenau), "Section accelerate"_accelerate_cuda.html, USER/cuda, -, lib/cuda
|
||||
USER-EFF, electron force field, Andres Jaramillo-Botero (Caltech), "pair_style eff/cut"_pair_eff.html, USER/eff, "eff"_eff, -
|
||||
USER-FEP, free energy perturbation, Agilio Padua (U Blaise Pascal Clermont-Ferrand), "fix adapt/fep"_fix_adapt.html, USER/fep, -, -
|
||||
USER-FEP, free energy perturbation, Agilio Padua (U Blaise Pascal Clermont-Ferrand), "compute fep"_compute_fep.html, USER/fep, -, -
|
||||
USER-INTEL, Vectorized CPU and Intel(R) coprocessor styles, W. Michael Brown (Intel), "Section accelerate"_accelerate_intel.html, examples/intel, -, -
|
||||
USER-LB, Lattice Boltzmann fluid, Colin Denniston (U Western Ontario), "fix lb/fluid"_fix_lb_fluid.html, USER/lb, -, -
|
||||
USER-MISC, single-file contributions, USER-MISC/README, USER-MISC/README, -, -, -
|
||||
|
||||
@ -23,7 +23,7 @@ attribute = {pair} or {atom} :l
|
||||
I,J = type pair(s) to set parameter for
|
||||
v_delta = variable with perturbation to apply (in the units of the parameter)
|
||||
{atom} args = aparam I v_delta
|
||||
aparam = parameter to adapt over time
|
||||
aparam = parameter to perturb
|
||||
I = type to set parameter for
|
||||
v_delta = variable with perturbation to apply (in the units of the parameter) :pre
|
||||
zero or more keyword/value pairs may be appended :l
|
||||
|
||||
@ -106,12 +106,12 @@ pair_coeff 1 1 1.0 9.5 :pre
|
||||
The {lj/cut/soft} style and substyles compute the 12/6 Lennard-Jones
|
||||
and Coulomb potential modified by a soft core, in order to avoid
|
||||
singularities during free energy calculations when sites are created
|
||||
or anihilated "(Beutler)"_#Beutler.
|
||||
or anihilated "(Beutler)"_#Beutler,
|
||||
|
||||
:c,image(Eqs/pair_lj_soft.jpg)
|
||||
|
||||
Coulomb interactions are also damped with a soft core at short
|
||||
distance as
|
||||
distance,
|
||||
|
||||
:c,image(Eqs/pair_coul_soft.jpg)
|
||||
|
||||
@ -120,11 +120,11 @@ are the charges on the 2 atoms, and epsilon is the dielectric constant
|
||||
which can be set by the "dielectric"_dielectric.html command.
|
||||
|
||||
The coefficient lambda is an activation parameter. When lambda = 1 the
|
||||
pair potentiel is identical to a Lennard-Jones term plus a Coulomb
|
||||
term. When lambda = 0 the pair interactions are deactivated. The
|
||||
transition between these two extrema is smoothed by a repulsive soft
|
||||
core in order to avoid singularities in potential energy and forces
|
||||
when sites are created or anihilated and can overlap
|
||||
pair potentiel is identical to a Lennard-Jones term or a Coulomb term
|
||||
or a combination of both. When lambda = 0 the interactions are
|
||||
deactivated. The transition between these two extrema is smoothed by a
|
||||
soft repulsive core in order to avoid singularities in potential
|
||||
energy and forces when sites are created or anihilated and can overlap
|
||||
"(Beutler)"_#Beutler.
|
||||
|
||||
The paratemers n, alpha_LJ and alpha_C are set in the
|
||||
@ -132,11 +132,34 @@ The paratemers n, alpha_LJ and alpha_C are set in the
|
||||
choices for the exponent are n = 2 or n = 1. For the remaining
|
||||
coefficients alpha_LJ = 0.5 and alpha_C = 10 Angstrom^2 are
|
||||
appropriate choices. Plots of the LJ and Coulomb terms are shown
|
||||
below, for lambda ranging from 1 ro 0 every 0.1.
|
||||
below, for lambda ranging from 1 to 0 every 0.1.
|
||||
|
||||
:c,image(JPG/lj_soft.jpg),image(JPG/coul_soft.jpg)
|
||||
|
||||
The {lj/cut/tip4p/long/soft} implements a soft-core version of the
|
||||
For the {lj/cut/coul/cut/soft} or {lj/cut/coul/long/soft} pair styles,
|
||||
the following coefficients must be defined for each pair of atoms
|
||||
types via the "pair_coeff"_pair_coeff.html command as in the examples
|
||||
above, or in the data file or restart files read by the
|
||||
"read_data"_read_data.html or "read_restart"_read_restart.html
|
||||
commands, or by mixing as described below:
|
||||
|
||||
epsilon (energy units)
|
||||
sigma (distance units)
|
||||
lambda (activation parameter between 0 and 1)
|
||||
cutoff1 (distance units)
|
||||
cutoff2 (distance units) :ul
|
||||
|
||||
The latter two coefficients are optional. If not specified, the global
|
||||
LJ and Coulombic cutoffs specified in the pair_style command are used.
|
||||
If only one cutoff is specified, it is used as the cutoff for both LJ
|
||||
and Coulombic interactions for this type pair. If both coefficients
|
||||
are specified, they are used as the LJ and Coulombic cutoffs for this
|
||||
type pair. You cannot specify 2 cutoffs for style {lj/cut/soft},
|
||||
since it has no Coulombic terms. For the {coul/cut/soft} and
|
||||
{coul/long/soft} only lambda and the optional cutoff2 are to be
|
||||
specified.
|
||||
|
||||
Style {lj/cut/tip4p/long/soft} implements a soft-core version of the
|
||||
TIP4P water model. The usage of this pair style is documented in the
|
||||
"pair_lj"_pair_lj.html styles. The soft-core version introduces the
|
||||
lambda parameter to the list of arguments, after epsilon and sigma in
|
||||
@ -165,28 +188,18 @@ occur. These substyles are suitable to represent charges embedded in
|
||||
the Lennard-Jones radius of another site (for example hydrogen atoms
|
||||
in several water models).
|
||||
|
||||
For the {lj/cut/coul/cut/soft} or {lj/cut/coul/long/soft} pair styles,
|
||||
the following coefficients must be defined for each pair of atoms
|
||||
types via the "pair_coeff"_pair_coeff.html command as in the examples
|
||||
above, or in the data file or restart files read by the
|
||||
"read_data"_read_data.html or "read_restart"_read_restart.html
|
||||
commands, or by mixing as described below:
|
||||
|
||||
epsilon (energy units)
|
||||
sigma (distance units)
|
||||
lambda (activation parameter between 0 and 1)
|
||||
cutoff1 (distance units)
|
||||
cutoff2 (distance units) :ul
|
||||
|
||||
The latter two coefficients are optional. If not specified, the global
|
||||
LJ and Coulombic cutoffs specified in the pair_style command are used.
|
||||
If only one cutoff is specified, it is used as the cutoff for both LJ
|
||||
and Coulombic interactions for this type pair. If both coefficients
|
||||
are specified, they are used as the LJ and Coulombic cutoffs for this
|
||||
type pair. You cannot specify 2 cutoffs for style {lj/cut/soft},
|
||||
since it has no Coulombic terms. For the {coul/cut/soft} and
|
||||
{coul/long/soft} only lambda and the optional cutoff2 are to be
|
||||
specified.
|
||||
IMPORTANT NOTES: When using the core-softed Coulomb potentials with
|
||||
long-range solvers ({coul/long/soft}, {lj/cut/coul/long/soft}, etc.)
|
||||
in a free energy calculation in which sites holding electrostatic
|
||||
charges are being created or anihilated (using
|
||||
"fix_adapt/fep"_fix_adapt_fep.html and "compute_fep"_compute_fep.html)
|
||||
it is important to adapt both the lambda activation parameter (from 0
|
||||
to 1, or the reverse) and the value of the charge (from 0 to its final
|
||||
value, or the reverse). This ensures that long-range electrostatic
|
||||
terms (kspace) are correct. It is not necessary to use core-softed
|
||||
Coulomb potentials if the van der Waals site is present during the
|
||||
free-energy route, thus avoiding overlap of the charges. Examples are
|
||||
provided in the LAMMPS source directory tree, under examples/USER/fep.
|
||||
|
||||
:line
|
||||
|
||||
@ -236,9 +249,14 @@ to be specified in an input script that reads a restart file.
|
||||
|
||||
[Restrictions:]
|
||||
|
||||
To avoid division by zero do not set sigma = 0. When sites do not
|
||||
interact though the Lennard-Jones term use epsilon = 0 and sigma = 1
|
||||
for example, or else use the {coul/long/soft} or similar substyle.
|
||||
To avoid division by zero do not set sigma = 0; use the lambda
|
||||
parameter instead to activate/deactivate interactions, or use
|
||||
epsilon = 0 and sigma = 1. Alternatively, when sites do not
|
||||
interact though the Lennard-Jones term the {coul/long/soft} or
|
||||
similar substyle can be used via the
|
||||
"pair_style hybrid/overlay"_pair_hybrid.html command.
|
||||
|
||||
:line
|
||||
|
||||
All of the plain {soft} pair styles are part of the USER-FEP package.
|
||||
The {long} styles also requires the KSPACE package to be installed.
|
||||
|
||||
@ -1,22 +1,22 @@
|
||||
# Time-averaged data for fix fFEP
|
||||
# TimeStep c_cFEP[1] c_cFEP[2]
|
||||
125000 0.000490778 0.999186
|
||||
250000 0.00619786 0.989679
|
||||
375000 0.0128045 0.978813
|
||||
500000 0.0203102 0.96667
|
||||
625000 0.0272161 0.955808
|
||||
750000 0.0215735 0.965087
|
||||
875000 0.0148509 0.975821
|
||||
1000000 0.00936378 0.984752
|
||||
1125000 0.00490088 0.992007
|
||||
1250000 0.00388501 0.99366
|
||||
1375000 0.00122887 0.998043
|
||||
1500000 0.000498063 0.999247
|
||||
1625000 -0.000468616 1.00085
|
||||
1750000 -0.00131731 1.00226
|
||||
1875000 -0.00247035 1.00419
|
||||
2000000 -0.00306799 1.00519
|
||||
2125000 -0.00355565 1.006
|
||||
2250000 -0.00434289 1.00733
|
||||
2375000 -0.00497634 1.00839
|
||||
2500000 -0.00576373 1.00972
|
||||
# Time-averaged data for fix FEP
|
||||
# TimeStep c_FEP[1] c_FEP[2]
|
||||
100000 0.000289416 0.999521
|
||||
200000 0.00590502 0.990163
|
||||
300000 0.0115179 0.980934
|
||||
400000 0.0216052 0.96457
|
||||
500000 0.0222451 0.963797
|
||||
600000 0.0200038 0.967603
|
||||
700000 0.0152292 0.975176
|
||||
800000 0.00896315 0.985384
|
||||
900000 0.00585213 0.990434
|
||||
1000000 0.00327599 0.99467
|
||||
1100000 0.00159845 0.997437
|
||||
1200000 0.000171108 0.999804
|
||||
1300000 -0.000313183 1.00059
|
||||
1400000 -0.00148 1.00254
|
||||
1500000 -0.00297976 1.00504
|
||||
1600000 -0.00287926 1.00487
|
||||
1700000 -0.00430671 1.00727
|
||||
1800000 -0.00453729 1.00765
|
||||
1900000 -0.00507356 1.00856
|
||||
2000000 -0.00586954 1.0099
|
||||
|
||||
@ -30,51 +30,54 @@ pair_coeff 1 4 lj/cut/tip4p/long/soft 0.0000 1.0000 0.0 # C4H Hw
|
||||
pair_coeff 2 3 lj/cut/tip4p/long/soft 0.0699 2.8126 0.0 # H Ow
|
||||
pair_coeff 2 4 lj/cut/tip4p/long/soft 0.0000 1.0000 0.0 # H Hw
|
||||
|
||||
variable nsteps equal 2500000
|
||||
variable nprint equal ${nsteps}/500
|
||||
variable ndump equal ${nsteps}/100
|
||||
variable TK equal 300.0
|
||||
variable PBAR equal 1.0
|
||||
|
||||
variable temp equal 300.0
|
||||
variable press equal 1.0
|
||||
|
||||
fix fSHAKE all shake 0.0001 20 ${nprint} b 2 a 2
|
||||
fix SHAKE all shake 0.0001 20 0 b 2 a 2
|
||||
|
||||
neighbor 2.0 bin
|
||||
|
||||
timestep 2.0
|
||||
timestep 1.0
|
||||
|
||||
velocity all create ${temp} 12345
|
||||
velocity all create ${TK} 12345
|
||||
|
||||
thermo_style multi
|
||||
thermo ${nprint}
|
||||
thermo 5000
|
||||
|
||||
fix fNPT all npt temp ${temp} ${temp} 100 iso ${press} ${press} 1000
|
||||
fix TPSTAT all npt temp ${TK} ${TK} 100 iso ${PBAR} ${PBAR} 1000
|
||||
|
||||
run 250000
|
||||
set type 1*2 charge 0.0
|
||||
|
||||
run 100000
|
||||
reset_timestep 0
|
||||
|
||||
variable lambda equal ramp(0.0,1.0)
|
||||
variable q1 equal -0.24*v_lambda
|
||||
variable q2 equal 0.06*v_lambda
|
||||
|
||||
fix fADAPT all adapt/fep 125000 pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda after yes
|
||||
fix ADAPT all adapt/fep 100000 &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda &
|
||||
atom charge 1 v_q1 &
|
||||
atom charge 2 v_q2 &
|
||||
after yes
|
||||
|
||||
fix PRINT all print 100000 "adapt lambda = ${lambda} q1 = ${q1} q2 = ${q2}"
|
||||
|
||||
variable dlambda equal 0.002
|
||||
variable dq1 equal -0.24*v_dlambda
|
||||
variable dq2 equal 0.06*v_dlambda
|
||||
|
||||
compute cFEP all fep ${temp} pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda
|
||||
compute FEP all fep ${TK} &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda &
|
||||
atom charge 1 v_dq1 &
|
||||
atom charge 2 v_dq2
|
||||
|
||||
fix fFEP all ave/time 25 4000 125000 c_cFEP[1] c_cFEP[2] file fdti01.lmp
|
||||
fix FEP all ave/time 20 4000 100000 c_FEP[1] c_FEP[2] file fdti01.lmp
|
||||
|
||||
# compute cRDF all rdf 100 3 3 3 4 4 4
|
||||
# fix fRDF all ave/time 2000 500 ${nsteps} c_cRDF file rdf.lammps mode vector
|
||||
dump TRAJ all custom 20000 dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify TRAJ element C H O H
|
||||
|
||||
group owater type 3
|
||||
compute cMSD owater msd
|
||||
fix fMSD owater ave/time 1 1 ${ndump} c_cMSD[1] c_cMSD[2] c_cMSD[3] c_cMSD[4] file msd.lammps
|
||||
|
||||
dump dCONF all custom ${ndump} dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify dCONF element C H O H
|
||||
|
||||
run ${nsteps}
|
||||
run 2000000
|
||||
|
||||
# write_restart restart.*.lmp
|
||||
write_data data.*.lmp
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
@ -1,22 +1,22 @@
|
||||
# Time-averaged data for fix fFEP
|
||||
# TimeStep c_cFEP[1] c_cFEP[2]
|
||||
125000 0.00652155 0.989123
|
||||
250000 0.00572718 0.990445
|
||||
375000 0.00506788 0.991544
|
||||
500000 0.00448303 0.992522
|
||||
625000 0.00370796 0.993818
|
||||
750000 0.00277535 0.995385
|
||||
875000 0.00230833 0.996171
|
||||
1000000 0.000982551 0.998409
|
||||
1125000 0.000314678 0.99954
|
||||
1250000 -0.000839025 1.0015
|
||||
1375000 -0.00155615 1.00273
|
||||
1500000 -0.00320598 1.00553
|
||||
1625000 -0.00587133 1.01011
|
||||
1750000 -0.00867306 1.01495
|
||||
1875000 -0.0134768 1.02327
|
||||
2000000 -0.0229314 1.03995
|
||||
2125000 -0.0282333 1.04896
|
||||
2250000 -0.0202762 1.03477
|
||||
2375000 -0.0123279 1.02096
|
||||
2500000 -0.00546203 1.00923
|
||||
# Time-averaged data for fix FEP
|
||||
# TimeStep c_FEP[1] c_FEP[2]
|
||||
100000 0.00671766 0.988799
|
||||
200000 0.00591179 0.990139
|
||||
300000 0.00495535 0.991733
|
||||
400000 0.00428164 0.992859
|
||||
500000 0.00367253 0.993879
|
||||
600000 0.00334988 0.994424
|
||||
700000 0.00246906 0.995908
|
||||
800000 0.00139088 0.997725
|
||||
900000 0.000491084 0.99924
|
||||
1000000 0.000615407 0.999044
|
||||
1100000 -0.00172714 1.00302
|
||||
1200000 -0.00333962 1.00579
|
||||
1300000 -0.0040424 1.00699
|
||||
1400000 -0.0102385 1.01763
|
||||
1500000 -0.0173611 1.03006
|
||||
1600000 -0.0243837 1.04234
|
||||
1700000 -0.0211518 1.03654
|
||||
1800000 -0.0220119 1.03776
|
||||
1900000 -0.013304 1.02261
|
||||
2000000 -0.00648381 1.01095
|
||||
|
||||
@ -30,51 +30,53 @@ pair_coeff 1 4 lj/cut/tip4p/long/soft 0.0000 1.0000 1.0 # C4H Hw
|
||||
pair_coeff 2 3 lj/cut/tip4p/long/soft 0.0699 2.8126 1.0 # H Ow
|
||||
pair_coeff 2 4 lj/cut/tip4p/long/soft 0.0000 1.0000 1.0 # H Hw
|
||||
|
||||
variable nsteps equal 2500000
|
||||
variable nprint equal ${nsteps}/500
|
||||
variable ndump equal ${nsteps}/100
|
||||
variable TK equal 300.0
|
||||
variable PBAR equal 1.0
|
||||
|
||||
variable temp equal 300.0
|
||||
variable press equal 1.0
|
||||
|
||||
fix fSHAKE all shake 0.0001 20 ${nprint} b 2 a 2
|
||||
fix SHAKE all shake 0.0001 20 0 b 2 a 2
|
||||
|
||||
neighbor 2.0 bin
|
||||
|
||||
timestep 2.0
|
||||
timestep 1.0
|
||||
|
||||
velocity all create ${temp} 12345
|
||||
velocity all create ${TK} 12345
|
||||
|
||||
thermo_style multi
|
||||
thermo ${nprint}
|
||||
thermo 5000
|
||||
|
||||
fix fNPT all npt temp ${temp} ${temp} 100 iso ${press} ${press} 1000
|
||||
fix TPSTAT all npt temp ${TK} ${TK} 100 iso ${PBAR} ${PBAR} 1000
|
||||
|
||||
run 250000
|
||||
run 100000
|
||||
|
||||
reset_timestep 0
|
||||
|
||||
variable lambda equal ramp(1.0,0.0)
|
||||
variable q1 equal -0.24*v_lambda
|
||||
variable q2 equal 0.06*v_lambda
|
||||
|
||||
fix fADAPT all adapt/fep 125000 pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda after yes
|
||||
fix ADAPT all adapt/fep 100000 &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda &
|
||||
atom charge 1 v_q1 &
|
||||
atom charge 2 v_q2 &
|
||||
after yes
|
||||
|
||||
fix PRINT all print 100000 "adapt lambda = ${lambda} q1 = ${q1} q2 = ${q2}"
|
||||
|
||||
variable dlambda equal -0.002
|
||||
variable dq1 equal -0.24*v_dlambda
|
||||
variable dq2 equal 0.06*v_dlambda
|
||||
|
||||
compute cFEP all fep ${temp} pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda
|
||||
compute FEP all fep ${TK} &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda &
|
||||
atom charge 1 v_dq1 &
|
||||
atom charge 2 v_dq2
|
||||
|
||||
fix fFEP all ave/time 25 4000 125000 c_cFEP[1] c_cFEP[2] file fdti10.lmp
|
||||
fix FEP all ave/time 20 4000 100000 c_FEP[1] c_FEP[2] file fdti10.lmp
|
||||
|
||||
# compute cRDF all rdf 100 3 3 3 4 4 4
|
||||
# fix fRDF all ave/time 2000 500 ${nsteps} c_cRDF file rdf.lammps mode vector
|
||||
dump TRAJ all custom 20000 dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify TRAJ element C H O H
|
||||
|
||||
group owater type 3
|
||||
compute cMSD owater msd
|
||||
fix fMSD owater ave/time 1 1 ${ndump} c_cMSD[1] c_cMSD[2] c_cMSD[3] c_cMSD[4] file msd.lammps
|
||||
|
||||
dump dCONF all custom ${ndump} dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify dCONF element C H O H
|
||||
|
||||
run ${nsteps}
|
||||
run 2000000
|
||||
|
||||
# write_restart restart.*.lmp
|
||||
write_data data.*.lmp
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
@ -1,22 +1,22 @@
|
||||
# Time-averaged data for fix fFEP
|
||||
# TimeStep c_cFEP[1] c_cFEP[2]
|
||||
125000 0.0806308 0.881044
|
||||
250000 0.251798 0.66945
|
||||
375000 0.448135 0.502774
|
||||
500000 0.670275 0.384034
|
||||
625000 0.873864 0.343148
|
||||
750000 0.689533 0.484145
|
||||
875000 0.474153 0.600305
|
||||
1000000 0.303566 0.740249
|
||||
1125000 0.165693 0.855143
|
||||
1250000 0.130949 0.88524
|
||||
1375000 0.0518948 0.976114
|
||||
1500000 0.0283335 1.0025
|
||||
1625000 -0.000790228 1.04259
|
||||
1750000 -0.0263791 1.0736
|
||||
1875000 -0.0587947 1.12568
|
||||
2000000 -0.0764999 1.15472
|
||||
2125000 -0.091264 1.17869
|
||||
2250000 -0.11294 1.21778
|
||||
2375000 -0.130721 1.25162
|
||||
2500000 -0.151808 1.29392
|
||||
# Time-averaged data for fix FEP
|
||||
# TimeStep c_FEP[1] c_FEP[2]
|
||||
100000 0.0735182 0.889583
|
||||
200000 0.241868 0.679931
|
||||
300000 0.407677 0.542008
|
||||
400000 0.709112 0.360902
|
||||
500000 0.718538 0.428553
|
||||
600000 0.639674 0.516854
|
||||
700000 0.482835 0.586307
|
||||
800000 0.289216 0.746055
|
||||
900000 0.192641 0.823932
|
||||
1000000 0.113029 0.908737
|
||||
1100000 0.0619301 0.96572
|
||||
1200000 0.0197356 1.01976
|
||||
1300000 0.00310596 1.03223
|
||||
1400000 -0.0300484 1.08295
|
||||
1500000 -0.0714914 1.14599
|
||||
1600000 -0.0712604 1.14573
|
||||
1700000 -0.109089 1.2123
|
||||
1800000 -0.117256 1.22671
|
||||
1900000 -0.132337 1.25534
|
||||
2000000 -0.153557 1.29745
|
||||
|
||||
@ -30,51 +30,54 @@ pair_coeff 1 4 lj/cut/tip4p/long/soft 0.0000 1.0000 0.0 # C4H Hw
|
||||
pair_coeff 2 3 lj/cut/tip4p/long/soft 0.0699 2.8126 0.0 # H Ow
|
||||
pair_coeff 2 4 lj/cut/tip4p/long/soft 0.0000 1.0000 0.0 # H Hw
|
||||
|
||||
variable nsteps equal 2500000
|
||||
variable nprint equal ${nsteps}/500
|
||||
variable ndump equal ${nsteps}/100
|
||||
variable TK equal 300.0
|
||||
variable PBAR equal 1.0
|
||||
|
||||
variable temp equal 300.0
|
||||
variable press equal 1.0
|
||||
|
||||
fix fSHAKE all shake 0.0001 20 ${nprint} b 2 a 2
|
||||
fix SHAKE all shake 0.0001 20 0 b 2 a 2
|
||||
|
||||
neighbor 2.0 bin
|
||||
|
||||
timestep 2.0
|
||||
timestep 1.0
|
||||
|
||||
velocity all create ${temp} 12345
|
||||
velocity all create ${TK} 12345
|
||||
|
||||
thermo_style multi
|
||||
thermo ${nprint}
|
||||
thermo 5000
|
||||
|
||||
fix fNPT all npt temp ${temp} ${temp} 100 iso ${press} ${press} 1000
|
||||
fix TPSTAT all npt temp ${TK} ${TK} 100 iso ${PBAR} ${PBAR} 1000
|
||||
|
||||
run 250000
|
||||
set type 1*2 charge 0.0
|
||||
|
||||
run 100000
|
||||
reset_timestep 0
|
||||
|
||||
variable lambda equal ramp(0.0,1.0)
|
||||
variable q1 equal -0.24*v_lambda
|
||||
variable q2 equal 0.06*v_lambda
|
||||
|
||||
fix fADAPT all adapt/fep 125000 pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda after yes
|
||||
fix ADAPT all adapt/fep 100000 &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda &
|
||||
atom charge 1 v_q1 &
|
||||
atom charge 2 v_q2 &
|
||||
after yes
|
||||
|
||||
fix PRINT all print 100000 "adapt lambda = ${lambda} q1 = ${q1} q2 = ${q2}"
|
||||
|
||||
variable dlambda equal 0.05
|
||||
variable dq1 equal -0.24*v_dlambda
|
||||
variable dq2 equal 0.06*v_dlambda
|
||||
|
||||
compute cFEP all fep ${temp} pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda
|
||||
compute FEP all fep ${TK} &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda &
|
||||
atom charge 1 v_dq1 &
|
||||
atom charge 2 v_dq2
|
||||
|
||||
fix fFEP all ave/time 25 4000 125000 c_cFEP[1] c_cFEP[2] file fep01.lmp
|
||||
fix FEP all ave/time 20 4000 100000 c_FEP[1] c_FEP[2] file fep01.lmp
|
||||
|
||||
# compute cRDF all rdf 100 3 3 3 4 4 4
|
||||
# fix fRDF all ave/time 2000 500 ${nsteps} c_cRDF file rdf.lammps mode vector
|
||||
dump TRAJ all custom 20000 dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify TRAJ element C H O H
|
||||
|
||||
group owater type 3
|
||||
compute cMSD owater msd
|
||||
fix fMSD owater ave/time 1 1 ${ndump} c_cMSD[1] c_cMSD[2] c_cMSD[3] c_cMSD[4] file msd.lammps
|
||||
|
||||
dump dCONF all custom ${ndump} dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify dCONF element C H O H
|
||||
|
||||
run ${nsteps}
|
||||
run 2000000
|
||||
|
||||
# write_restart restart.*.lmp
|
||||
write_data data.*.lmp
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
@ -1,22 +1,22 @@
|
||||
# Time-averaged data for fix fFEP
|
||||
# TimeStep c_cFEP[1] c_cFEP[2]
|
||||
125000 0.156203 0.771285
|
||||
250000 0.136448 0.799233
|
||||
375000 0.12198 0.820426
|
||||
500000 0.109453 0.840255
|
||||
625000 0.0921141 0.868102
|
||||
750000 0.0713015 0.905735
|
||||
875000 0.0622431 0.922352
|
||||
1000000 0.0330546 0.981942
|
||||
1125000 0.0197866 1.01134
|
||||
1250000 -0.00412826 1.06765
|
||||
1375000 -0.0174302 1.11343
|
||||
1500000 -0.0511646 1.18453
|
||||
1625000 -0.105676 1.36325
|
||||
1750000 -0.161616 1.56478
|
||||
1875000 -0.258263 1.93825
|
||||
2000000 -0.447487 3.0595
|
||||
2125000 -0.549377 3.01163
|
||||
2250000 -0.379009 2.01933
|
||||
2375000 -0.20934 1.45678
|
||||
2500000 -0.0644056 1.1252
|
||||
# Time-averaged data for fix FEP
|
||||
# TimeStep c_FEP[1] c_FEP[2]
|
||||
100000 0.160455 0.766272
|
||||
200000 0.141912 0.791902
|
||||
300000 0.119585 0.824744
|
||||
400000 0.104694 0.847809
|
||||
500000 0.0917124 0.869318
|
||||
600000 0.0859541 0.881974
|
||||
700000 0.0666319 0.920268
|
||||
800000 0.0428654 0.965893
|
||||
900000 0.0240106 0.999889
|
||||
1000000 0.0294147 0.997671
|
||||
1100000 -0.0214947 1.11583
|
||||
1200000 -0.0543642 1.23078
|
||||
1300000 -0.0665772 1.25986
|
||||
1400000 -0.196675 1.66521
|
||||
1500000 -0.341037 2.39926
|
||||
1600000 -0.480217 2.97287
|
||||
1700000 -0.405599 2.36311
|
||||
1800000 -0.415165 2.12414
|
||||
1900000 -0.229669 1.49644
|
||||
2000000 -0.0838847 1.15743
|
||||
|
||||
@ -30,51 +30,53 @@ pair_coeff 1 4 lj/cut/tip4p/long/soft 0.0000 1.0000 1.0 # C4H Hw
|
||||
pair_coeff 2 3 lj/cut/tip4p/long/soft 0.0699 2.8126 1.0 # H Ow
|
||||
pair_coeff 2 4 lj/cut/tip4p/long/soft 0.0000 1.0000 1.0 # H Hw
|
||||
|
||||
variable nsteps equal 2500000
|
||||
variable nprint equal ${nsteps}/500
|
||||
variable ndump equal ${nsteps}/100
|
||||
variable TK equal 300.0
|
||||
variable PBAR equal 1.0
|
||||
|
||||
variable temp equal 300.0
|
||||
variable press equal 1.0
|
||||
|
||||
fix fSHAKE all shake 0.0001 20 ${nprint} b 2 a 2
|
||||
fix SHAKE all shake 0.0001 20 0 b 2 a 2
|
||||
|
||||
neighbor 2.0 bin
|
||||
|
||||
timestep 2.0
|
||||
timestep 1.0
|
||||
|
||||
velocity all create ${temp} 12345
|
||||
velocity all create ${TK} 12345
|
||||
|
||||
thermo_style multi
|
||||
thermo ${nprint}
|
||||
thermo 5000
|
||||
|
||||
fix fNPT all npt temp ${temp} ${temp} 100 iso ${press} ${press} 1000
|
||||
fix TPSTAT all npt temp ${TK} ${TK} 100 iso ${PBAR} ${PBAR} 1000
|
||||
|
||||
run 250000
|
||||
run 100000
|
||||
|
||||
reset_timestep 0
|
||||
|
||||
variable lambda equal ramp(1.0,0.0)
|
||||
variable q1 equal -0.24*v_lambda
|
||||
variable q2 equal 0.06*v_lambda
|
||||
|
||||
fix fADAPT all adapt/fep 125000 pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda after yes
|
||||
fix ADAPT all adapt/fep 100000 &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_lambda &
|
||||
atom charge 1 v_q1 &
|
||||
atom charge 2 v_q2 &
|
||||
after yes
|
||||
|
||||
fix PRINT all print 100000 "adapt lambda = ${lambda} q1 = ${q1} q2 = ${q2}"
|
||||
|
||||
variable dlambda equal -0.05
|
||||
variable dq1 equal -0.24*v_dlambda
|
||||
variable dq2 equal 0.06*v_dlambda
|
||||
|
||||
compute cFEP all fep ${temp} pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda
|
||||
compute FEP all fep ${TK} &
|
||||
pair lj/cut/tip4p/long/soft lambda 1*2 3*4 v_dlambda &
|
||||
atom charge 1 v_dq1 &
|
||||
atom charge 2 v_dq2
|
||||
|
||||
fix fFEP all ave/time 25 4000 125000 c_cFEP[1] c_cFEP[2] file fep10.lmp
|
||||
fix FEP all ave/time 20 4000 100000 c_FEP[1] c_FEP[2] file fep10.lmp
|
||||
|
||||
# compute cRDF all rdf 100 3 3 3 4 4 4
|
||||
# fix fRDF all ave/time 2000 500 ${nsteps} c_cRDF file rdf.lammps mode vector
|
||||
dump TRAJ all custom 20000 dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify TRAJ element C H O H
|
||||
|
||||
group owater type 3
|
||||
compute cMSD owater msd
|
||||
fix fMSD owater ave/time 1 1 ${ndump} c_cMSD[1] c_cMSD[2] c_cMSD[3] c_cMSD[4] file msd.lammps
|
||||
|
||||
dump dCONF all custom ${ndump} dump.lammpstrj id mol type element x y z ix iy iz
|
||||
dump_modify dCONF element C H O H
|
||||
|
||||
run ${nsteps}
|
||||
run 2000000
|
||||
|
||||
# write_restart restart.*.lmp
|
||||
write_data data.*.lmp
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
Reference in New Issue
Block a user