Dear Feng,
Indeed, I’m running the example and I was surprised to observe such huge differences only by changing the output_freq. Here are the input files and a part of the output files referring to the NFE section for:
_harm=10 kcal/mol and output_freq=1000:
Malonaldehyde: NCSU steered MD: LCOD rxn coord
&cntrl
imin = 0,
irest = 1,
ntx = 5,
ntb = 0,
cut = 999.0,
tempi = 300.0,
temp0 = 300.0,
ntt = 3,
gamma_ln = 1.0,
ntf = 2, ntc = 2, tol = 0.00001,
dt = 0.0005,
ntpr = 50, ntwr=100, ntwx = 100,
nstlim = 100000,
ifqnt = 1, infe=1
/
&qmmm
qmmask=':*',
qmcharge=0,
qm_theory='DFTB',
qmshake=0,
writepdb=1,
/
&smd
output_file = 'smd_10.txt'
output_freq = 1000
cv_file='cv.in'
/
NFE : output_file = smd_10.txt
NFE : output_freq = 1000 (0.5000 ps)
NFE :
NFE : CV #1
NFE :
NFE : <> path = (0.7000, -0.7000)
NFE : <> path_mode = SPLINE
NFE : <> harm = (10.0000)
NFE : <> harm_mode = SPLINE
NFE :
NFE : type = 'LCOD' (Linear Combination Of Distances)
NFE : 1.000 * (8 [O2] <=> 9 [H4])
NFE : -1.000 * (7 [O1] <=> 9 [H4])
_harm=10 kcal/mol and output_freq=500:
Malonaldehyde: NCSU steered MD: LCOD rxn coord
&cntrl
imin = 0,
irest = 1,
ntx = 5,
ntb = 0,
cut = 999.0,
tempi = 300.0,
temp0 = 300.0,
ntt = 3,
gamma_ln = 1.0,
ntf = 2, ntc = 2, tol = 0.00001,
dt = 0.0005,
ntpr = 50, ntwr=100, ntwx = 100,
nstlim = 100000,
ifqnt = 1, infe=1
/
&qmmm
qmmask=':*',
qmcharge=0,
qm_theory='PM3',
qmshake=0,
writepdb=1,
/
&smd
output_file = 'smd_10-500.txt'
output_freq = 500
cv_file='cv.in'
/
NFE : output_file = smd_10-500.txt
NFE : output_freq = 500 (0.2500 ps)
NFE :
NFE : CV #1
NFE :
NFE : <> path = (0.7000, -0.7000)
NFE : <> path_mode = SPLINE
NFE : <> harm = (10.0000)
NFE : <> harm_mode = SPLINE
NFE :
NFE : type = 'LCOD' (Linear Combination Of Distances)
NFE : 1.000 * (8 [O2] <=> 9 [H4])
NFE : -1.000 * (7 [O1] <=> 9 [H4])
I’ve plotted the CV vs work and as you point out the CV doesn’t follow very well the steering path when harm=10 kcal/mol.. That is when I started changing other parameters.
As I haven’t set the random seed I guess it was different for each run.. On this matter, when doing smd for my protein and its ligand, should I run the simulation with a random seed? If yes, how many replicas should I do in order to obtain “reliable” results?
And, which would be the “ideal” value for the string constant in order not to have such huge differences in the total work? I’ve tried the example with harm=100 kcal/mol and the CV vs work plot is not as bad as for harm=10 kcal/mol, but is 100 kcal/mol a reasonable value (not too high)? Also, I’ve just got the output files of a SMD of my protein performed with harm= 1000 kcal/mol and:
_During many steps the Self Consistency is not achieved
_The water molecule involved in a proton transfer I wish to observe completely blows as one of its hydrogens is 18 A° away from the oxygen atom and yet the proton transfer doesn’t occur.. But, I’should confirm this results again as I think I’ve made a mistake on the steering path..
I really appreciate the time you spend in order to enlighten me on this subject.
Thank you for any forward suggestion/enlightenment,
Sincerely,
Ruth
> On 24 Nov 2016, at 19:55, Feng Pan <fpan3.ncsu.edu> wrote:
>
> Hi, Ruth
>
> The output_freq only affects the output so it should not affect the total
> work. Are you sure you ran both cases with exact same conditions? The
> random seed may also affect the
> results of SMD.
> Are you running the example? I think the harm strength 10kcal/mol is too
> low because the CV path does not follow the Steering path very well. Try a
> larger harm strength, could be better.
> If the strength is too low, the steering does not take many effects, so the
> total work may differ a lot from different runs.
>
> Feng
>
> On Thu, Nov 24, 2016 at 6:34 AM, Ruth Helena Tichauer <rhtichau.laas.fr>
> wrote:
>
>> Dear Amber users,
>>
>> I’ve been following tutorial A10 in order to perform SMD for my protein
>> and its ligand using the LCOD method. In the tutorial, the spring constant
>> for the proton transfer to occur is set to 1000 kcal/mol whose value,
>> according to a few related questions on the mailing list, is too high.
>>
>> In order to find a reasonable value, I’ve run the simulation several times
>> and the proton transfer takes place for “low” values such as 10 kcal/mol if
>> the output_freq is set to 1000, like in the tutorial. But, if the
>> output_freq is lowered (I’ve tried 500 and 50) the proton transfer doesn’t
>> occur anymore (for harm=10kcal/mol).
>>
>> Moreover, in the output_file .txt the CV, handle_position and work are not
>> the same for a given time during the MD, resulting in different values for
>> the total work done.
>>
>> For harm=10 kcal/mol and output_freq=1000 I obtain:
>>
>> #
>> # MD time (ps), CV, handle_position, spring_constant, work
>> #
>> 925.0000 0.50139387 0.70000000 10.00000000
>> 0.00000000
>> 925.5000 0.89660018 0.69958280 10.00000000
>> 0.00030631
>> 926.0000 0.83105576 0.69834240 10.00000000
>> 0.00081418
>> 926.5000 0.66819394 0.69629560 10.00000000
>> 0.00204283
>> 927.0000 0.60709575 0.69345920 10.00000000
>> 0.00467988
>> 927.5000 0.93546672 0.68985000 10.00000000
>> 0.00719652
>> 928.0000 0.56761201 0.68548480 10.00000000
>> 0.01073920
>> 928.5000 0.62737834 0.68038040 10.00000000
>> 0.01357336
>> 929.0000 0.92623955 0.67455360 10.00000000
>> 0.01711340
>> 929.5000 0.66861857 0.66802120 10.00000000
>> 0.02039487
>> 930.0000 0.91611501 0.66080000 10.00000000
>> 0.02901486
>> ...
>> #
>> # <> total work done: 0.7690942380
>>
>> For harm=10 kcal/mol and output_freq=500 I obtain:
>>
>> #
>> # MD time (ps), CV, handle_position, spring_constant, work
>> #
>> 925.0000 0.50139387 0.70000000 10.00000000
>> 0.00000000
>> 925.2500 0.80711710 0.69989535 10.00000000
>> 0.00022605
>> 925.5000 1.00173702 0.69958280 10.00000000
>> 0.00094337
>> 925.7500 0.81672099 0.69906445 10.00000000
>> 0.00198729
>> 926.0000 0.81074280 0.69834240 10.00000000
>> 0.00416796
>> 926.2500 0.93641120 0.69741875 10.00000000
>> 0.00658384
>> 926.5000 0.79972556 0.69629560 10.00000000
>> 0.00872512
>> 926.7500 0.88288670 0.69497505 10.00000000
>> 0.01112337
>> 927.0000 0.99659450 0.69345920 10.00000000
>> 0.01380751
>> 927.2500 0.91748709 0.69175015 10.00000000
>> 0.01664516
>> 927.5000 0.98411987 0.68985000 10.00000000
>> 0.01972750
>> 927.7500 0.76115131 0.68776085 10.00000000
>> 0.02316773
>> 928.0000 0.82942938 0.68548480 10.00000000
>> 0.02737620
>> 928.2500 0.74668613 0.68302395 10.00000000
>> 0.03189408
>> 928.5000 0.95082231 0.68038040 10.00000000
>> 0.03673629
>> 928.7500 0.79738517 0.67755625 10.00000000
>> 0.04177707
>> 929.0000 1.03951245 0.67455360 10.00000000
>> 0.04712890
>> 929.2500 0.71139449 0.67137455 10.00000000
>> 0.05324293
>> 929.5000 0.85144838 0.66802120 10.00000000
>> 0.05961067
>> 929.7500 0.92964077 0.66449565 10.00000000
>> 0.06619512
>> 930.0000 0.84766718 0.66080000 10.00000000
>> 0.07329916
>> ...
>> #
>> # <> total work done: 11.2413859667
>>
>> So I wonder what exactly the output_freq achieves during the SMD? As I
>> want to run the same kind of simulation for my protein and its ligand, I’d
>> like to know how to choose this value.
>>
>> Thank you for any insight on this subject,
>>
>> Sincerely,
>>
>> Ruth
>> _______________________________________________
>> AMBER mailing list
>> AMBER.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber
>>
>
>
>
> --
> Feng Pan
> Ph.D. Candidate
> North Carolina State University
> Department of Physics
> Email: fpan3.ncsu.edu
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Fri Nov 25 2016 - 04:00:03 PST
