# sgld (c39b2)

Self-Guided Langevin Dynamics (SGLD)

and

Self-Guided Molecular Dynamics (SGMD) Simulation Methods

By Xiongwu Wu and Bernard R. Brooks

SGMD/SGLD is a method to enhance conformational searching and sampling

through enhancing the low frequency motion. The main references for the

SGLD simulation method are:

(1) Wu, X., Brooks, B.R., "Self-guided Langevin dynamics simulation method",

Chem. Phys. Letter, 381, 512-518(2003).

(2) Wu, X., Brooks, B.R., "Toward Canonical Ensemble Distribution from

Self-guided Langevin dynamics simulation", J. Chem. Phys., 134, 134108 (2011).

(3) Wu, X., Brooks, B.R., "Force-momentum based self-guided Langevin dynamics:

a rapid sampling method that approaches the canonical ensemble", J. Chem.

Phys., 135, 204101 (2011)

(4) Wu, X., Damjanovic, A., Brooks, B.R., "Efficient and unbiased sampling

of biomolecular systems in the canonical ensemble: a review of self-guided

Langevin dynamics", Adv. Chem. Phys., Vol.150, 255-326 (2012)

(5) Wu, X., Hodoscek, M., Brooks, B.R., "Replica exchange of self-guided

Langevin dynamics simulation", J. Chem. Phys., 137, 044106 (2012).

* Syntax | Syntax of the SGLD dynamics command

* Background | Description of SGMD/SGLD methods

* Examples | SGLD usage examples

Top

[Syntax SGLD]

SCALar SGWT SET 0.0 ! clean guiding weight array: SGWT

SCALar SGWT SET 1.0 atom-selection ! select atoms to set guiding weights

DYNAmics {LEAP {[LANG SGLD]} } ! SGLD simulation

{ {[HOOVER|CPT] SGMD]} } ! SGMD for a NVT/NPT ensemble

{ {[ SGMD]} sg-spec } ! SGMD for a NVE ensemble

{ } other-dynamic-spec

sg-spec::= {[SGFT real][SGBZ]}[TSGAVG real][TREFLF real] other-sg-spec

{[TEMPSG real] }

other-sg-spec::= [SGFF real] [TSGAVP real] [ISGSTA int] [ISGEND int]

[SGSTOPT] [SGSTOPR] [SGCOM]

Keyword Default Purpose

SGLD false Turn on SGLD simulation (for Langevin dynamics, LANG).

Friction constants, FBETA, must be set at first.

It is equivalent to SGMD for non-Langevin dynamics.

SGMD false Turn on SGMD simulation (for MD simulations).

It is equivalent to SGLD for Langevin dynamics.

SGBZ false Activate SGMDfp/SGLDfp to maintain a Boltzmann conformational

distribution. When SGBZ is on, SGFF is estimated from

simulation to mantain a Boltzmann distribution.

SGFT 0.2(SGMD) Momentum guiding factor. This is the lamda in the SGLD

1.0(SGLD) equation of motion. When TEMPSG=0, SGFT will be

constant throughout a simulation.

TEMPSG 0 K Target self-guiding temperature. Please note that the TEMPSG

definition has been changed since c36. TEMPSG set the

conformatonal search ability comparable to that in a simulation

at temperature TEMPSG. TEMPSG>0 will override SGFT. a TEMPSG

higher than TBATH or TREF will enhance a conformational search

while a TEMPSG lower than TBATH or TREF will slow down a

conformational search. TEMPSG should not be set when SGBZ is set.

When TEMPSG>0, SGFT fluctuates during the simulation to maintain a guiding

temperature of TEMPSG. If a constant SGFT simulation is

desired, one can read the approximate values of SGFT from

the output of a SGMD/SGLD simulation with TEMPSG=2*TREF as a suggested value.

TSGAVG 0.2 ps Local average time. A larger TSGAVG will result in slower

motion to be enhanced. All motions with periods larger than

TSGAVG will be enhanced.

ISGSTA 1 The index of the first atom of the region to apply the

guiding force

ISGEND natom The index of the last atom of the region to apply the

guiding force. If reweighting is to be performed and ISGSTA>1 or

ISGEND<natom, one need turn on analysis with "anal on" before

running dynamics.

SGSTOPT false Stop translation of the center of mass. By default, SGLD

allows the center of mass to have an random motion. By set

SGSTOPT, SGLD stops any translation in the center of mass.

SGSTOPR false Stop rotation around the center of mass. By default, SGLD

allows the center of mass to have an random motion. By set

SGSTOPR, SGLD stops any rotation around the center of mass.

In SGLD simulations, NTRFRQ has no effect if SGSTOPT and

SGSTOPR is set.

SGCOM false Allow net guiding force and torque to the system. By default the net

guiding force and torque on the center of mass is removed. Set SGCOM

allows accelerating the motion of the center of mass.

TSGAVP 10*TSGAVG Convergency time. This time is used to control convergency of

running averages to estimate instaneous parameters like SGFTI,

TEMPSGI, TEMPLFI, etc. A larger TSGAVP will result in a smaller

fluctuation of instaneous parameters.

TREFLF 0 K Reference low frequency temperature in a canonical ensemble.

TREFLF is the low frequency temperature when no guiding force is

applied to a simulation system. TREFLF is required for weighting

factor calculation in SGMD/SGLD simulation or for the guiding

force calculation in SGMDfp/SGLDfp simulations (when SGBZ is

turned on). TREFLF can be determined from a SGMD/SGLD simulation

with SGFT=0 and TEMPSG=0, TREFLF=<TEMPLF>. The average TEMPLF value

can be accessed by querying ?TEMPLF in CHARMM input script.

When TREFLF=0, TREFLF is estimated during a simulation. An actual

value for TREFLF can increase the accuracy in SGMDfp/SGLDfp simulations.

SGFF 0.0 Force guiding factor. Defines how much low frequency forces are included

in the guiding forces. It is used only in SGMDf/SGLDf simulations.

SGWT {1.0} Self-guiding weight array for atoms. Using SCALAR to set SGWT.

SGWT is used to select atoms to apply guiding forces.

[Syntax SGLD]

SCALar SGWT SET 0.0 ! clean guiding weight array: SGWT

SCALar SGWT SET 1.0 atom-selection ! select atoms to set guiding weights

DYNAmics {LEAP {[LANG SGLD]} } ! SGLD simulation

{ {[HOOVER|CPT] SGMD]} } ! SGMD for a NVT/NPT ensemble

{ {[ SGMD]} sg-spec } ! SGMD for a NVE ensemble

{ } other-dynamic-spec

sg-spec::= {[SGFT real][SGBZ]}[TSGAVG real][TREFLF real] other-sg-spec

{[TEMPSG real] }

other-sg-spec::= [SGFF real] [TSGAVP real] [ISGSTA int] [ISGEND int]

[SGSTOPT] [SGSTOPR] [SGCOM]

Keyword Default Purpose

SGLD false Turn on SGLD simulation (for Langevin dynamics, LANG).

Friction constants, FBETA, must be set at first.

It is equivalent to SGMD for non-Langevin dynamics.

SGMD false Turn on SGMD simulation (for MD simulations).

It is equivalent to SGLD for Langevin dynamics.

SGBZ false Activate SGMDfp/SGLDfp to maintain a Boltzmann conformational

distribution. When SGBZ is on, SGFF is estimated from

simulation to mantain a Boltzmann distribution.

SGFT 0.2(SGMD) Momentum guiding factor. This is the lamda in the SGLD

1.0(SGLD) equation of motion. When TEMPSG=0, SGFT will be

constant throughout a simulation.

TEMPSG 0 K Target self-guiding temperature. Please note that the TEMPSG

definition has been changed since c36. TEMPSG set the

conformatonal search ability comparable to that in a simulation

at temperature TEMPSG. TEMPSG>0 will override SGFT. a TEMPSG

higher than TBATH or TREF will enhance a conformational search

while a TEMPSG lower than TBATH or TREF will slow down a

conformational search. TEMPSG should not be set when SGBZ is set.

When TEMPSG>0, SGFT fluctuates during the simulation to maintain a guiding

temperature of TEMPSG. If a constant SGFT simulation is

desired, one can read the approximate values of SGFT from

the output of a SGMD/SGLD simulation with TEMPSG=2*TREF as a suggested value.

TSGAVG 0.2 ps Local average time. A larger TSGAVG will result in slower

motion to be enhanced. All motions with periods larger than

TSGAVG will be enhanced.

ISGSTA 1 The index of the first atom of the region to apply the

guiding force

ISGEND natom The index of the last atom of the region to apply the

guiding force. If reweighting is to be performed and ISGSTA>1 or

ISGEND<natom, one need turn on analysis with "anal on" before

running dynamics.

SGSTOPT false Stop translation of the center of mass. By default, SGLD

allows the center of mass to have an random motion. By set

SGSTOPT, SGLD stops any translation in the center of mass.

SGSTOPR false Stop rotation around the center of mass. By default, SGLD

allows the center of mass to have an random motion. By set

SGSTOPR, SGLD stops any rotation around the center of mass.

In SGLD simulations, NTRFRQ has no effect if SGSTOPT and

SGSTOPR is set.

SGCOM false Allow net guiding force and torque to the system. By default the net

guiding force and torque on the center of mass is removed. Set SGCOM

allows accelerating the motion of the center of mass.

TSGAVP 10*TSGAVG Convergency time. This time is used to control convergency of

running averages to estimate instaneous parameters like SGFTI,

TEMPSGI, TEMPLFI, etc. A larger TSGAVP will result in a smaller

fluctuation of instaneous parameters.

TREFLF 0 K Reference low frequency temperature in a canonical ensemble.

TREFLF is the low frequency temperature when no guiding force is

applied to a simulation system. TREFLF is required for weighting

factor calculation in SGMD/SGLD simulation or for the guiding

force calculation in SGMDfp/SGLDfp simulations (when SGBZ is

turned on). TREFLF can be determined from a SGMD/SGLD simulation

with SGFT=0 and TEMPSG=0, TREFLF=<TEMPLF>. The average TEMPLF value

can be accessed by querying ?TEMPLF in CHARMM input script.

When TREFLF=0, TREFLF is estimated during a simulation. An actual

value for TREFLF can increase the accuracy in SGMDfp/SGLDfp simulations.

SGFF 0.0 Force guiding factor. Defines how much low frequency forces are included

in the guiding forces. It is used only in SGMDf/SGLDf simulations.

SGWT {1.0} Self-guiding weight array for atoms. Using SCALAR to set SGWT.

SGWT is used to select atoms to apply guiding forces.

Top

Backgrounds

SGMD/SGLD simulation enhances conformational search efficiency through

acceleration of low frequency motions in a molecular system. The low frequency

motion represents the slow conformational change that contributes the most

in conformational search. A guiding force resonant with the low freuqency

motion, which is calculated based on a local average of momentums, is introduced

into the equation of motion to enhance the low frequency motion. It has been

demonstrated that within the suggested parameter range, a SGLD simulation

has a dramatically enhanced conformational search efficiency while has little

perturbation on conformational distribution. The conformational distribution of

a SGLD simulation More details are described in

the paper:

The SGLD algorithm is developed based on a previously developed

force-based self-guided molecular dynamics simulation method (SGMD,

J. Phys. Chem. B, 102, 7238(1998) and J. Chem. Phys., 110, 9401(1999)), which

employs a local average of non-bonded interaction forces as the guiding force.

The force-based guiding force correlates with the energy surface that may

bias to the low energy conformations. The SGLD method uses a local

average of velocities to calculate the guiding force, and overcomes the

disadvantage of the old method. When SGLD is used in MD simulations without

a friction, it is referred to as SGMDp (momentum-based self-guided molecular

dynamics).

The self-guided dynamics simulation enhances conformational search and sampling

efficiency by accelerating low frequency motion. The conformational distribution in

a self-guided dynamics simulation can be quantitatively described by the so-called

low frequency properties, like low frequency temperature and low frequency energy.

A weighting factor for each frame is printed out during simulation and ensemble average

properties can be calculated as weighted avearges. The self-guided dynamcs simulation

methods are further improved to directly produce a canonical ensemble distribution

by incooperate the low frequency force and low frequency motions in the guiding

force calculation. This new method is referred to as SGMDfp or SGLDfp, depending on

the presents of the friction forces. This extension makes SGMD/SGLD an excellent tool

for both conformational search and conformational sampling.

In CHARMM input script, only a key word, SGMD or SGLD, in a "DYNA ..."

simulation command line is needed to turn on SGMD or SGLD with default values

for all parameters. For LD or SGLD simulation, the friction constant, FBETA,

need be set. Whether to use SGMD or SGLD is automatically decided depending on

whether the friction constant is zero or not.

SGMD/SGLD has two primary input parameters, the guiding factor, SGFT, and

the local average time, TSGAVG, which define the guiding effect in a SGMD/SGLD

simulation. SGFT set the strength of the guiding force and TSGAVG defines the

slow motion mode to be enhanced. A larger TSGAVG will result in a slower

motion (with period longer than TSGAVG) to be enhanced. A larger SGFT will

introduce stronger guiding forces and result in a larger energy barrier

overcoming ability. When SGFT=0, a SGMD/SGLD simulation reduces to a normal

MD/LD simulation. SGFT should be limited to keep low frequency motion from changing

too much. The low frequency motion is measured by the low frequency temperature,

which is calculated from the local average momentum and depends on TSGAVG and SGFT.

To avoid the difficulty in choosing the value of SGFT, a parameter

called guiding temperature (TEMPSG) is introduced. TEMPSG, as an alternate

to SGFT, specifies the conformational searching ability that is comparable to

a simulation at temperature TEMPSG. When TEMPSG=0, SGFT will be used to calculate

the guiding force. If TEMPSG>0, SGFT will fluctuate during a simulation to

reach a conformational search comparable to a simulation at a temperature of TEMPSG.

TEMPSG>TBATH for a LD simulation or TEMPSG>TREF for a MD simultaion will enhance

conformational search, while TEMPSG<TBATH or TEMPSG<TREF will slow down conformational

search.

The Key word, SGBZ, will turn on the SGMDfp/SGLDfp method to maintain a

canonical ensemble. A input parameter, TREFLF, can be used to input

the reference low frequency temperature, which is the low frequency temperature

when no guiding force is applied. If TREFLF is not input, it will be estimated

during a SGMD/SGLD simulation. When SGBZ is set, TEMPSG should not be used to set

guiding effect, instead, SGFT should be used.

The guiding range, ISGSTA and ISGEND, define the starting and ending atom indexes

to apply the guiding forces. By default, all atoms are included in the guiding range.

Within the range, the guiding weight array, SGWT, can be set with the SCALAR command.

The actual guiding factor is SGFT*SGWT(I) for atom I. By default, SGWT=1.

SGLD can be used for replica exchange, abbreviated as RXSGLD, to obtain

accelerated conformational search while preserves the canonical ensemble. The guiding

effect, SGFT or TEMPSG, is used to define different replicas, with or without

temperature changes. Without temperature changes, RXSGLD can achieve very high

exchange efficiency and can be applied to much larger systems than temperature-based

replica exchange simulations. More details can be found in repdstr.doc.

When SGMD or SGLD is turned on, DYNA simulation output will contain

two lines with the following quantities:

DYNA SGLF> SGFT TEMPSG TEMPLF TREFLF FRCLF EPOTLF SGWT

DYNA SGHF> SGFF SGFD TEMPHF TREFHF FRCHF EPOTHF VIRSG

These quantities are instaneouse values defined as below:

SGFT: Momentum guiding factor

SGFF: Force guiding factor. Adjusted when SGBZ is on

SGFD: Force dumping factor. Adjusted when SGBZ is on.

TEMPSG: Guiding temperature.

SGEXP: Probability exponent. exp(SGEXP) is the weighting factor of current frame.

VIRSG: Virial of the guiding force.

TEMPLF: low frequency temperature

TEMPHF: high frequency temperature

TREFLF: reference low frequency temperature

TREFHF: reference high frequency temperature

FRCLF: low frequency force factor

FRCHF: high frequency force factor

EPOTLF: low frequency potential energy

EPOTHF: high frequency potential energy

The weight of a conformation is calculated by

Weight= exp(SGEXP)

=exp(((FRCLF*TREFLF/TEMPLF-1)*EPOTLF+(FRCHF*TREFHF/TEMPHF-1)*EPOTHF+

VIRSG)/(KBOLTZ*Temp))

The averages of above variables can be queried by the following

parameters:

?SGFT=<SGFTI>, ?SGFF=<SGFFI>, ?TEMPSG=<TEMPSG>, ?SGFD=<SGFD>,?SGEXP=<SGEXP>,

?TEMPLF=<TEMPLF>,?TEMPHF=<TEMPHF>, ?TREFLF=<TREFLF>,?TREFHF=<TREFHF>,

?FRCLF=<FRCLF>,?FRCHF=<FRCHF>,?EPOTLF=<EPOTLF>,?EPOTHF=<EPOTHF>,?VIRSG=<VIRSG>,

In addition to the DYNA output, for each trajectory frame, there is a

SGMD/SGLD output marked by "TRAJ SGLF>" and "TRAJ SGHF>", which can be used for

post-processing ensemble distributions.

SGMD/SGLD simulation results can be weighted to produce canonical ensemble

properties. The weighting information can be extracted from simulation outputs.

A simple command line awk script is provided here to get a weighted average

from a SGMD/SGLD simulation:

%> awk 'BEGIN{n=0;if(N0==0)N0=1000;} \

{if($1=="DYNA>"){ep=$6;n++;} if(n<N0)next; \

if($1=="DYNA"&& $2=="SGLF>"){wt=exp($9);aw+=wt;ae+=ep;aew+=wt*ep;m++;}} \

END{printf("Ep and weighted Ep: %10.4f %10.4f\n",ae/m,aew/aw);}' sgmd-sgld.out

The following awk script will print weighting factors for every frame in

a SGMD/SGLD trajectory:

%> awk '{if($1=="TRAJ"&&$2=="SGLF>")printf(" %10d %10.4f\n",++i,exp($9));}' sgmd-sgld.out

The following awk script will recalculate and print weighting factors for every frame in

a SGMD/SGLD trajectory:

%> cat >sgwt.awk

BEGIN{n=0;if(TEMP==0)TEMP=300;kt=0.001987*TEMP;if(N0==0)N0=10;} \

{if($1=="TRAJ"&&$2=="SGLF>"){templf=$5;treflf=$6;frclf=$7;epotlf=$8;sgwt=$9;} \

if($1=="TRAJ"&&$2=="SGHF>"){temphf=$5;trefhf=$6;frchf=$7;epothf=$8;virsg=$9; \

eplfs[n]=epotlf;ephfs[n]=epothf;virs[n]=virsg;n++;if(n<N0)next; \

avgtlf+=templf;avgrlf+=treflf;avgflf+=frclf;avgelf+=epotlf; \

avgthf+=temphf;avgrhf+=trefhf;avgfhf+=frchf;avgehf+=epothf;avgvir+=virsg;}} \

END{m=n-N0;avgtlf/=m;avgrlf/=m;avgflf/=m;avgelf/=m; \

avgthf/=m;avgrhf/=m;avgfhf/=m;avgehf/=m;avgvir/=m; \

if(TREFLF>0){avgrlf=TREFLF;avgrhf=TEMP-TREFLF;} if(TREFHF==0){avgthf=TEMP-avgtlf;} \

for(i=0;i<n;i++){ \

wt=exp(((avgflf*avgrlf/avgtlf-1)*(eplfs[i]-avgelf)+(avgfhf*avgrhf/avgthf-1)*(ephfs[i]-avgehf)+virs[i]-avgvir)/kt); \

printf(" %6d %10.4f\n", i+1,wt);}}

EOF

Run awk with above script can produce weighting factor for each frame in the trajectory file:

%> awk -f sgwt.awk -v N0=100 -v TREFLF=8.5 -v TEMP=400 sgmd-sgld.out

Here, N0 is equilibrium frame number (default 10), TEMP is simulation temperature (default is 300K),

and TREFLF is the reference low frequency temperature (default is the estimated value from simulation).

Backgrounds

SGMD/SGLD simulation enhances conformational search efficiency through

acceleration of low frequency motions in a molecular system. The low frequency

motion represents the slow conformational change that contributes the most

in conformational search. A guiding force resonant with the low freuqency

motion, which is calculated based on a local average of momentums, is introduced

into the equation of motion to enhance the low frequency motion. It has been

demonstrated that within the suggested parameter range, a SGLD simulation

has a dramatically enhanced conformational search efficiency while has little

perturbation on conformational distribution. The conformational distribution of

a SGLD simulation More details are described in

the paper:

The SGLD algorithm is developed based on a previously developed

force-based self-guided molecular dynamics simulation method (SGMD,

J. Phys. Chem. B, 102, 7238(1998) and J. Chem. Phys., 110, 9401(1999)), which

employs a local average of non-bonded interaction forces as the guiding force.

The force-based guiding force correlates with the energy surface that may

bias to the low energy conformations. The SGLD method uses a local

average of velocities to calculate the guiding force, and overcomes the

disadvantage of the old method. When SGLD is used in MD simulations without

a friction, it is referred to as SGMDp (momentum-based self-guided molecular

dynamics).

The self-guided dynamics simulation enhances conformational search and sampling

efficiency by accelerating low frequency motion. The conformational distribution in

a self-guided dynamics simulation can be quantitatively described by the so-called

low frequency properties, like low frequency temperature and low frequency energy.

A weighting factor for each frame is printed out during simulation and ensemble average

properties can be calculated as weighted avearges. The self-guided dynamcs simulation

methods are further improved to directly produce a canonical ensemble distribution

by incooperate the low frequency force and low frequency motions in the guiding

force calculation. This new method is referred to as SGMDfp or SGLDfp, depending on

the presents of the friction forces. This extension makes SGMD/SGLD an excellent tool

for both conformational search and conformational sampling.

In CHARMM input script, only a key word, SGMD or SGLD, in a "DYNA ..."

simulation command line is needed to turn on SGMD or SGLD with default values

for all parameters. For LD or SGLD simulation, the friction constant, FBETA,

need be set. Whether to use SGMD or SGLD is automatically decided depending on

whether the friction constant is zero or not.

SGMD/SGLD has two primary input parameters, the guiding factor, SGFT, and

the local average time, TSGAVG, which define the guiding effect in a SGMD/SGLD

simulation. SGFT set the strength of the guiding force and TSGAVG defines the

slow motion mode to be enhanced. A larger TSGAVG will result in a slower

motion (with period longer than TSGAVG) to be enhanced. A larger SGFT will

introduce stronger guiding forces and result in a larger energy barrier

overcoming ability. When SGFT=0, a SGMD/SGLD simulation reduces to a normal

MD/LD simulation. SGFT should be limited to keep low frequency motion from changing

too much. The low frequency motion is measured by the low frequency temperature,

which is calculated from the local average momentum and depends on TSGAVG and SGFT.

To avoid the difficulty in choosing the value of SGFT, a parameter

called guiding temperature (TEMPSG) is introduced. TEMPSG, as an alternate

to SGFT, specifies the conformational searching ability that is comparable to

a simulation at temperature TEMPSG. When TEMPSG=0, SGFT will be used to calculate

the guiding force. If TEMPSG>0, SGFT will fluctuate during a simulation to

reach a conformational search comparable to a simulation at a temperature of TEMPSG.

TEMPSG>TBATH for a LD simulation or TEMPSG>TREF for a MD simultaion will enhance

conformational search, while TEMPSG<TBATH or TEMPSG<TREF will slow down conformational

search.

The Key word, SGBZ, will turn on the SGMDfp/SGLDfp method to maintain a

canonical ensemble. A input parameter, TREFLF, can be used to input

the reference low frequency temperature, which is the low frequency temperature

when no guiding force is applied. If TREFLF is not input, it will be estimated

during a SGMD/SGLD simulation. When SGBZ is set, TEMPSG should not be used to set

guiding effect, instead, SGFT should be used.

The guiding range, ISGSTA and ISGEND, define the starting and ending atom indexes

to apply the guiding forces. By default, all atoms are included in the guiding range.

Within the range, the guiding weight array, SGWT, can be set with the SCALAR command.

The actual guiding factor is SGFT*SGWT(I) for atom I. By default, SGWT=1.

SGLD can be used for replica exchange, abbreviated as RXSGLD, to obtain

accelerated conformational search while preserves the canonical ensemble. The guiding

effect, SGFT or TEMPSG, is used to define different replicas, with or without

temperature changes. Without temperature changes, RXSGLD can achieve very high

exchange efficiency and can be applied to much larger systems than temperature-based

replica exchange simulations. More details can be found in repdstr.doc.

When SGMD or SGLD is turned on, DYNA simulation output will contain

two lines with the following quantities:

DYNA SGLF> SGFT TEMPSG TEMPLF TREFLF FRCLF EPOTLF SGWT

DYNA SGHF> SGFF SGFD TEMPHF TREFHF FRCHF EPOTHF VIRSG

These quantities are instaneouse values defined as below:

SGFT: Momentum guiding factor

SGFF: Force guiding factor. Adjusted when SGBZ is on

SGFD: Force dumping factor. Adjusted when SGBZ is on.

TEMPSG: Guiding temperature.

SGEXP: Probability exponent. exp(SGEXP) is the weighting factor of current frame.

VIRSG: Virial of the guiding force.

TEMPLF: low frequency temperature

TEMPHF: high frequency temperature

TREFLF: reference low frequency temperature

TREFHF: reference high frequency temperature

FRCLF: low frequency force factor

FRCHF: high frequency force factor

EPOTLF: low frequency potential energy

EPOTHF: high frequency potential energy

The weight of a conformation is calculated by

Weight= exp(SGEXP)

=exp(((FRCLF*TREFLF/TEMPLF-1)*EPOTLF+(FRCHF*TREFHF/TEMPHF-1)*EPOTHF+

VIRSG)/(KBOLTZ*Temp))

The averages of above variables can be queried by the following

parameters:

?SGFT=<SGFTI>, ?SGFF=<SGFFI>, ?TEMPSG=<TEMPSG>, ?SGFD=<SGFD>,?SGEXP=<SGEXP>,

?TEMPLF=<TEMPLF>,?TEMPHF=<TEMPHF>, ?TREFLF=<TREFLF>,?TREFHF=<TREFHF>,

?FRCLF=<FRCLF>,?FRCHF=<FRCHF>,?EPOTLF=<EPOTLF>,?EPOTHF=<EPOTHF>,?VIRSG=<VIRSG>,

In addition to the DYNA output, for each trajectory frame, there is a

SGMD/SGLD output marked by "TRAJ SGLF>" and "TRAJ SGHF>", which can be used for

post-processing ensemble distributions.

SGMD/SGLD simulation results can be weighted to produce canonical ensemble

properties. The weighting information can be extracted from simulation outputs.

A simple command line awk script is provided here to get a weighted average

from a SGMD/SGLD simulation:

%> awk 'BEGIN{n=0;if(N0==0)N0=1000;} \

{if($1=="DYNA>"){ep=$6;n++;} if(n<N0)next; \

if($1=="DYNA"&& $2=="SGLF>"){wt=exp($9);aw+=wt;ae+=ep;aew+=wt*ep;m++;}} \

END{printf("Ep and weighted Ep: %10.4f %10.4f\n",ae/m,aew/aw);}' sgmd-sgld.out

The following awk script will print weighting factors for every frame in

a SGMD/SGLD trajectory:

%> awk '{if($1=="TRAJ"&&$2=="SGLF>")printf(" %10d %10.4f\n",++i,exp($9));}' sgmd-sgld.out

The following awk script will recalculate and print weighting factors for every frame in

a SGMD/SGLD trajectory:

%> cat >sgwt.awk

BEGIN{n=0;if(TEMP==0)TEMP=300;kt=0.001987*TEMP;if(N0==0)N0=10;} \

{if($1=="TRAJ"&&$2=="SGLF>"){templf=$5;treflf=$6;frclf=$7;epotlf=$8;sgwt=$9;} \

if($1=="TRAJ"&&$2=="SGHF>"){temphf=$5;trefhf=$6;frchf=$7;epothf=$8;virsg=$9; \

eplfs[n]=epotlf;ephfs[n]=epothf;virs[n]=virsg;n++;if(n<N0)next; \

avgtlf+=templf;avgrlf+=treflf;avgflf+=frclf;avgelf+=epotlf; \

avgthf+=temphf;avgrhf+=trefhf;avgfhf+=frchf;avgehf+=epothf;avgvir+=virsg;}} \

END{m=n-N0;avgtlf/=m;avgrlf/=m;avgflf/=m;avgelf/=m; \

avgthf/=m;avgrhf/=m;avgfhf/=m;avgehf/=m;avgvir/=m; \

if(TREFLF>0){avgrlf=TREFLF;avgrhf=TEMP-TREFLF;} if(TREFHF==0){avgthf=TEMP-avgtlf;} \

for(i=0;i<n;i++){ \

wt=exp(((avgflf*avgrlf/avgtlf-1)*(eplfs[i]-avgelf)+(avgfhf*avgrhf/avgthf-1)*(ephfs[i]-avgehf)+virs[i]-avgvir)/kt); \

printf(" %6d %10.4f\n", i+1,wt);}}

EOF

Run awk with above script can produce weighting factor for each frame in the trajectory file:

%> awk -f sgwt.awk -v N0=100 -v TREFLF=8.5 -v TEMP=400 sgmd-sgld.out

Here, N0 is equilibrium frame number (default 10), TEMP is simulation temperature (default is 300K),

and TREFLF is the reference low frequency temperature (default is the estimated value from simulation).

Top

Examples

1) SGMD simulation with default setting.

! Perform SGMD simulation

DYNA LEAP STRT CPT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGMD -

TCON TCOU 0.1 TREF 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

2) SGLD simulation with default setting.

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

! Perform SGLD simulation

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

3) SGLD simulation with a constant guiding factor.

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

! Perform SGLD simulation

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD TSGAVG 0.2 SGFT 1.0 -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

4) SGLD simulation with a guiding temperature of 500 K (to reach a searching ability comparable to

a 500 K high temperature simulation).

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

! clean self-guiding weights.

SCAL SGWT SET 0.0 SELE ALL END

! apply guiding force to only segment 1

SCAL SGWT SET 1.0 SELE iseg 1 end

! Perform SGLD simulation

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD TSGAVG 0.2 TEMPSG 500 -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

5) SGLDfp simulation with a constant guiding factor and to maintain ensemble distribution.

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD SGBZ TSGAVG 0.2 SGFT 1.0 -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

Examples

1) SGMD simulation with default setting.

! Perform SGMD simulation

DYNA LEAP STRT CPT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGMD -

TCON TCOU 0.1 TREF 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

2) SGLD simulation with default setting.

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

! Perform SGLD simulation

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

3) SGLD simulation with a constant guiding factor.

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

! Perform SGLD simulation

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD TSGAVG 0.2 SGFT 1.0 -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

4) SGLD simulation with a guiding temperature of 500 K (to reach a searching ability comparable to

a 500 K high temperature simulation).

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

! clean self-guiding weights.

SCAL SGWT SET 0.0 SELE ALL END

! apply guiding force to only segment 1

SCAL SGWT SET 1.0 SELE iseg 1 end

! Perform SGLD simulation

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD TSGAVG 0.2 TEMPSG 500 -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0

5) SGLDfp simulation with a constant guiding factor and to maintain ensemble distribution.

! Set friction forces

SCAL FBETA SET 1.0 SELE ALL END

DYNA LANG LEAP STRT NSTE 1000000 TIME 0.002 -

IPRFRQ 10000 ISVFRQ 1000 IHTFRQ 0 IEQFRQ 0 INBFRQ 10 IHBFRQ 0 -

IUNREA -1 IUNWRI 31 IUNCRD 32 IUNVEL 33 KUNIT -1 -

NSAVC 1000 NSAVV 00 NPRINT 1000 ISEED 314159 -

SGLD SGBZ TSGAVG 0.2 SGFT 1.0 -

TBATH 300 FIRST 260 -

NBXMOD 5 ATOM CDIEL VATOM VDISTANCE EIPS VIPS PXYZ -

CUTNB 12.0 CTOFNB 10.0 EPS 1.0 E14FAC 1.0 WMIN 1.0