Dear all,
I am trying to perform a Single Point Calculation (SPC) (to compute the
potential energy) of a series of molecular configurations stored in a
trajectory file.
To do that I am doing a geometry optimization with 0 steps in each
structure, but I am not sure it is working. The input file I am using is:
spc.mdin:
&cntrl
imin = 5,
maxcyc = 0,
cut = 999.0, ! The system is just an isolated molecule
rgbmax = 999.0,
ifqnt = 1
ntb=0
/
&qmmm
qmcut = 9999.9,
qmmask = '@*',
qmcharge = 0,
qm_theory = PM6-DH+,
qmshake = 1,
/
&wt type='END'
/
To check that this is indeed doing what I want, I applied it to a
trajectory generated with the same potential energy function, using the
following input file:
remd.mdin.001:
&cntrl
ntx=1,
nstlim=200, dt=0.001,
numexchg=20000,
irest=0, ntt=3, gamma_ln=1.0,
temp0=300, ig=963,
ntc=2, ntf=2, nscm=100,
ntb=0, igb=0,
cut=999.0, rgbmax=999.0,
ntpr=100, ntwx=100, ntwr=10000,
ioutfm=1,
ntc=2
ifqnt=1
/
&qmmm
qmcut = 9999.9,
qmmask = '.*',
qmcharge = 0,
qm_theory = PM6-DH+,
dftb_disper = 0,
qmshake = 1,
/
&wt type='END'
/
Which as I mentioned it has the same potential energy function. However,
the energies printed from the trajectory don't quite match the ones
computed using the SPC.
This is the remd.mdout.001 file were I get the energies from the simulation
4. RESULTS
--------------------------------------------------------------------------------
NSTEP = 0 TIME(PS) = 0.000 TEMP(K) = 0.00 PRESS
= 0.0
Etot = 55.8629 EKtot = 0.0000 EPtot = 55.8629
BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
1-4 NB = 0.0000 1-4 EEL = 0.0000 VDWAALS =
0.0000
EELEC = 0.0000 EHBOND = 0.0000 RESTRAINT =
0.0000
PM6ESCF= 55.8629
TEMP0 = 300.0000 REPNUM = 1 EXCHANGE#
= 0
------------------------------------------------------------------------------
NMR restraints: Bond = 0.000 Angle = 0.000 Torsion = 0.000
===============================================================================
| RE_POSITION Moving by 0.041055 0.028284 -0.019654
NetCDF error: NetCDF: Variable not found
at write replica mytargettemp
NSTEP = 100 TIME(PS) = 0.100 TEMP(K) = 55.15 PRESS
= 0.0
Etot = 60.8920 EKtot = 5.5893 EPtot = 55.3027
BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
1-4 NB = 0.0000 1-4 EEL = 0.0000 VDWAALS =
0.0000
EELEC = 0.0000 EHBOND = 0.0000 RESTRAINT =
0.0000
PM6ESCF= 55.3027
TEMP0 = 300.0000 REPNUM = 1 EXCHANGE#
= 1
------------------------------------------------------------------------------
NMR restraints: Bond = 0.000 Angle = 0.000 Torsion = 0.000
===============================================================================
| RE_POSITION Moving by -0.035491 -0.048501 0.026206
NetCDF error: NetCDF: Variable not found
at write replica mytargettemp
NSTEP = 200 TIME(PS) = 0.200 TEMP(K) = 87.63 PRESS
= 0.0
Etot = 65.2281 EKtot = 8.8813 EPtot = 56.3467
BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
1-4 NB = 0.0000 1-4 EEL = 0.0000 VDWAALS =
0.0000
EELEC = 0.0000 EHBOND = 0.0000 RESTRAINT =
0.0000
PM6ESCF= 56.3467
TEMP0 = 300.0000 REPNUM = 1 EXCHANGE#
= 1
------------------------------------------------------------------------------
NMR restraints: Bond = 0.000 Angle = 0.000 Torsion = 0.000
===============================================================================
| RE_POSITION Moving by 0.021387 -0.059176 0.025757
NetCDF error: NetCDF: Variable not found
at write replica mytargettemp
NSTEP = 300 TIME(PS) = 0.300 TEMP(K) = 105.27 PRESS
= 0.0
Etot = 70.0486 EKtot = 10.6692 EPtot = 59.3794
BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
1-4 NB = 0.0000 1-4 EEL = 0.0000 VDWAALS =
0.0000
EELEC = 0.0000 EHBOND = 0.0000 RESTRAINT =
0.0000
PM6ESCF= 59.3794
TEMP0 = 300.0000 REPNUM = 1 EXCHANGE#
= 2
------------------------------------------------------------------------------
And this is the output from the SPC in the same first 3 structures from
the trajectory (corresponding to the previous output):
4. RESULTS
--------------------------------------------------------------------------------
Maximum number of minimization cycles reached.
FINAL RESULTS
NSTEP ENERGY RMS GMAX NAME NUMBER
1 5.5525E+01 7.3313E+00 1.9797E+01 C7 14
BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
VDWAALS = 0.0000 EEL = 0.0000 HBOND = 0.0000
1-4 VDW = 0.0000 1-4 EEL = 0.0000 RESTRAINT =
0.0000
PM6ESCF = 55.5251
minimization completed, ENE= 0.55525091E+02 RMS= 0.733131E+01
minimizing coord set # 2
Maximum number of minimization cycles reached.
FINAL RESULTS
NSTEP ENERGY RMS GMAX NAME NUMBER
1 5.6441E+01 7.4006E+00 2.0826E+01 C11 28
BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
VDWAALS = 0.0000 EEL = 0.0000 HBOND = 0.0000
1-4 VDW = 0.0000 1-4 EEL = 0.0000 RESTRAINT =
0.0000
PM6ESCF = 56.4414
minimization completed, ENE= 0.56441442E+02 RMS= 0.740056E+01
minimizing coord set # 3
Maximum number of minimization cycles reached.
FINAL RESULTS
NSTEP ENERGY RMS GMAX NAME NUMBER
1 5.9277E+01 9.9923E+00 4.2688E+01 C3 4
BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
VDWAALS = 0.0000 EEL = 0.0000 HBOND = 0.0000
1-4 VDW = 0.0000 1-4 EEL = 0.0000 RESTRAINT =
0.0000
PM6ESCF = 59.2766
minimization completed, ENE= 0.59276627E+02 RMS= 0.999230E+01
minimizing coord set # 4
Maximum number of minimization cycles reached.
The energies then are:
From simulation From SPC
55.5251 55.3027
56.4414 56.3467
59.2766 59.3794
63.4846 63.6437
63.5524 63.3793
63.1787 62.8946
65.2079 65.0181
74.7678 74.1447
67.0064 68.0420
70.3709 70.7692
They are very correlated, but still different. Does anybody know what is
going on???
Thanks a lot for your help, I deeply appreciate it!!
All the best,
Raimon
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Fri Jul 13 2018 - 06:30:02 PDT
