Re: Autocorrelation function of Rg goes down without fluctuation around zero

From: Brian Radak (brian.radak.accts_at_gmail.com)
Date: Tue Aug 23 2016 - 14:28:36 CDT

On 08/23/2016 02:16 PM, faride badalkhani wrote:

> Hi,
> I mean did I include all the necessary parameters in the NVT
> configuration file?
> As you said I performed a NVT (for 300 ps) then NPT (for 300 ps) and
> then NPT for 20 ns (it is running now) and I have some questions:
> 1) when I performed NVT simulation water box shape changed to a semi
> spherical one but after performing NPT for 300 ps it got the cubic
> shape again. Could you tell me if what happened in NVT should not have
> happened?
>
I could not say. This all depends on your boundary conditions and the
Boltzmann factor of the initial configuration (is it a reasonable
structure for these parameters?).

> 2) what is the difference between config files we execute for first
> NPT (the short one) and the final NPT run? For example, in GROMACS for
> NVT and NPT we impose some forces on atoms and then for production run
> we remove all the forces and atoms can move. Is there something
> similar in NAMD, too?
>
That all depends. If you want to include a restraint or constraint
protocol, that is separate from NpT. Generally constraints can be
confusing for NpT because the virial must be modified - in my opinion
restraints are almost always better.
>
> Your help is really appreciated.
>
> Regards,
> Farideh
>
>
> On Aug 23, 2016 11:29 PM, "Brian Radak" <brian.radak.accts_at_gmail.com
> <mailto:brian.radak.accts_at_gmail.com>> wrote:
>
> How do you mean "correct"? Does the file execute successfully? It
> can be very useful to try this and then read the INFO header lines.
>
> NpT generally requires /more/ options than NVT. So long as you set
> "langevin on" (or some other proper thermostat) and do not set
> "langevinPiston on", you should get NVT results.
>
>
> On 08/22/2016 07:34 AM, faride badalkhani wrote:
>> Dear NAMD users,
>>
>> I want to perform a sequence of equilibration at NVT, then NPT,
>> and then production at NPT. I have some questions:
>> 1- I have never performed NVT simulations before, so could you
>> tell me if the following configuration file is correct?
>>
>> #############################################################
>> ## JOB DESCRIPTION ##
>> #############################################################
>>
>> # Minimization and Equilibration of
>> # G5 OH terminated dendrimer in a Water Box
>> # NVT ensemble
>>
>>
>> #############################################################
>> ## ADJUSTABLE PARAMETERS ##
>> #############################################################
>>
>> structure ./OH-OH/input/OH_wb.psf
>> coordinates ./OH-OH/input/OH_wb.pdb
>>
>> set temperature 300
>> set outputname ./OH-OH/output/OH_eq
>>
>> firsttimestep 0
>>
>>
>> #############################################################
>> ## SIMULATION PARAMETERS ##
>> #############################################################
>>
>> # Input
>> paraTypeCharmm on
>> parameters ./OH-OH/input/par_all36_OH.prm
>> temperature $temperature
>>
>>
>> # Force-Field Parameters
>> exclude scaled1-4
>> 1-4scaling 1.0
>> cutoff 12.0
>> switching on
>> switchdist 10.0
>> pairlistdist 14.0
>>
>>
>> # Integrator Parameters
>> timestep 1.0 ;# 1 fs/step
>> rigidBonds water ;# needed for 1 fs steps
>> nonbondedFreq 1
>> fullElectFrequency 2
>> stepspercycle 10
>>
>>
>> # Constant Temperature Control
>> langevin on ;# do langevin dynamics
>> langevinDamping 5 ;# damping coefficient (gamma) of 5/ps
>> langevinTemp $temperature
>> langevinHydrogen off ;# don't couple langevin bath to hydrogens
>>
>>
>> # Periodic Boundary Conditions
>> cellBasisVector1 73.5 0.0 0.0
>> cellBasisVector2 0.0 86.5 0.0
>> cellBasisVector3 0.0 0.0 88.4
>> cellOrigin -12.7 -12.5 -26.3
>>
>> wrapAll on
>>
>>
>> # PME (for full-system periodic electrostatics)
>> PME yes
>>
>>
>> # let NAMD determine grid
>> PMEGridSpacing 1.0
>>
>>
>> useGroupPressure yes ;# needed for rigidBonds
>> useFlexibleCell no
>> useConstantArea no
>>
>>
>> # Output
>> outputName $outputname
>>
>>
>> restartfreq 5000 ;# 5000 steps = every 5.0 ps
>> dcdfreq 5000
>> xstFreq 5000
>> outputEnergies 5000
>> outputPressure 5000
>>
>>
>> #############################################################
>> ## EXTRA PARAMETERS ##
>> #############################################################
>>
>>
>> #############################################################
>> ## EXECUTION SCRIPT ##
>> #############################################################
>>
>> # Minimization
>> minimize 20000
>> reinitvels $temperature
>>
>> run 300000 ;# 300 ps
>>
>>
>> 2- It would be necessary to set "rigidBonds water " when time
>> step is 1 fs?
>>
>> 3- How many time steps are standard values for the first NVT and
>> NPT simulations? (for example in Gromacs it is recommended to run
>> NVT and NPT for 300 ps if someone wants to run the priduction run
>> for 30 ns).
>>
>> Regards,
>> Farideh
>>
>> On Thu, Aug 18, 2016 at 9:07 PM, Brian Radak <bradak_at_anl.gov
>> <mailto:bradak_at_anl.gov>> wrote:
>>
>> I would think the standard sequence is equilibrate at NpT,
>> then NVT, and then production at NVT. Some people do things
>> differently, and that's fine. There's nothing wrong with that
>> - it just might be a bit slower (you could even check this
>> during equilibration). It's up to you.
>>
>> Restarting probably shouldn't affect results. It might if you
>> were trying to get very precise Newtonian dynamics, but there
>> are plenty of other concerns in that case also.
>>
>>
>> On 08/18/2016 11:31 AM, faride badalkhani wrote:
>>>
>>> And I have one more questions! Could you tell me if the
>>> restarting the simulations every 5 ns affects the results or
>>> not?
>>>
>>> Regards,
>>> Farideh
>>>
>>>
>>> On Aug 18, 2016 8:57 PM, "faride badalkhani"
>>> <farideh.khamseh_at_gmail.com
>>> <mailto:farideh.khamseh_at_gmail.com>> wrote:
>>>
>>> So, do you think it is more reasonable to perform
>>> minimization and then NVT simulation for a few time
>>> steps before production run?
>>>
>>> Regards,
>>> Farideh
>>>
>>>
>>> On Aug 18, 2016 8:50 PM, "Brian Radak" <bradak_at_anl.gov
>>> <mailto:bradak_at_anl.gov>> wrote:
>>>
>>> After equilibration to a desired density, pressure
>>> control is usually just an added cost (maybe 5-10%?)
>>> with not much value. I usually do everything at NVT.
>>> For equilibration at the force field density, these
>>> are usually indistinguishable (unless you want
>>> temperature or pressure dependent properties).
>>>
>>> On 08/18/2016 11:16 AM, faride badalkhani wrote:
>>>>
>>>> Usually, the autocorrelation function of Rg is used
>>>> to find if the structure is well equilibrated or
>>>> not, and to find relaxation time for them.
>>>> What do you mean with "It's also not common to use
>>>> pressure control for these kinds of calculations"?
>>>> Which method do you suggest for perfotming
>>>> simulations?
>>>>
>>>> Regards,
>>>> Farideh
>>>>
>>>>
>>>> On Aug 18, 2016 8:30 PM, "Brian Radak"
>>>> <bradak_at_anl.gov <mailto:bradak_at_anl.gov>> wrote:
>>>>
>>>> I've not used either of "measure rgyr" or
>>>> "g_analyze", so I can only assume they are
>>>> correctly implemented.
>>>>
>>>> I have pretty significant reservations about
>>>> including your equilibration components in the
>>>> calculation. It's also not common to use
>>>> pressure control for these kinds of
>>>> calculations, but I don't have a specific
>>>> reason as to why that's bad other than that it
>>>> changes the physical meaning of the result (I
>>>> don't know what you want C(t) for, I assume
>>>> this is not a problem).
>>>>
>>>> My guess is that you just have very poor
>>>> statistics due to a long characteristic time.
>>>> This is not unexpected for a (presumably large)
>>>> system like a dendrimer.
>>>>
>>>> Cheers,
>>>>
>>>> Brian
>>>>
>>>>
>>>> On 08/18/2016 10:31 AM, faride badalkhani wrote:
>>>>> This system is a Dendrimer, which is a
>>>>> hyperbranched polymer. I used the following
>>>>> script to calculate the radius of gyration:
>>>>>
>>>>> set outfile [open rg.dat w]
>>>>> set mol [molinfo top]
>>>>> set sel [atomselect top "all not waters"]
>>>>> set frames [molinfo $mol get numframes]
>>>>> for {set i 0} {$i < $frames} {incr i} {
>>>>> $sel frame $i
>>>>> $sel update
>>>>> puts $outfile "$i [measure rgyr $sel]"
>>>>> }
>>>>> close $outfile
>>>>> $sel delete
>>>>>
>>>>> For calculating the C(t) I used the g_analyze
>>>>> in gromacs and the Rg computed using the above
>>>>> script. for plotting the data I used the whole
>>>>> simulation time (40 ns). I loaded all the
>>>>> trajectory files in VMD and then calculated
>>>>> the Rg.
>>>>>
>>>>> Regards,
>>>>> Farideh
>>>>>
>>>>> On Thu, Aug 18, 2016 at 7:51 PM, Brian Radak
>>>>> <bradak_at_anl.gov <mailto:bradak_at_anl.gov>> wrote:
>>>>>
>>>>> What system is this? How did you compute
>>>>> Rg and C(t)? Which parts of the trajectory
>>>>> were used?
>>>>>
>>>>> Correlation functions can be very slow to
>>>>> converge, 5 ns might not be enough.
>>>>>
>>>>>
>>>>> On 08/18/2016 10:16 AM, faride badalkhani
>>>>> wrote:
>>>>>> Thank you for your answer and time. At
>>>>>> first I performed minimization for 20 ps
>>>>>> and then NPT simulation for 10 ns. After
>>>>>> that, I restarted the simulation for 5
>>>>>> more ns and continued simulation in 5 ns
>>>>>> steps till 40 ns. I did not get any error
>>>>>> during the simulations, but the
>>>>>> autocorrelation function does not
>>>>>> fluctuate around zero. I have shared the
>>>>>> Rg, RMSD, and C(t) ans also .top, .par,
>>>>>> .conf, .pdb, and the .log file for the
>>>>>> first 10 ns in the link below:
>>>>>>
>>>>>> https://www.dropbox.com/sh/mcnkxt9pefnqgrz/AADRRw4-dJEEx0SBlU2Fuekma?dl=0
>>>>>> <https://www.dropbox.com/sh/mcnkxt9pefnqgrz/AADRRw4-dJEEx0SBlU2Fuekma?dl=0>
>>>>>>
>>>>>> p.s. Bond distances, angles, dihedrals
>>>>>> and impropers were taken form CHARMM
>>>>>> GENERAL FORCE FIELD and the charges were
>>>>>> assigned using CGen FF on paramchem.
>>>>>>
>>>>>> Thank you so much for your help!
>>>>>> Regards,
>>>>>> Farideh
>>>>>>
>>>>>> On Thu, Aug 18, 2016 at 7:17 PM, Brian
>>>>>> Radak <bradak_at_anl.gov
>>>>>> <mailto:bradak_at_anl.gov>> wrote:
>>>>>>
>>>>>> You either made a mistake in your
>>>>>> calculation or else your statistics
>>>>>> are bad. Some details might help
>>>>>> distinguish which of those is the case.
>>>>>>
>>>>>> HTH,
>>>>>> Brian
>>>>>>
>>>>>>
>>>>>> On 08/18/2016 04:21 AM, faride
>>>>>> badalkhani wrote:
>>>>>>
>>>>>> Dear NAMD users,
>>>>>>
>>>>>> I have plotted the
>>>>>> autocorrelation function of
>>>>>> radius of gyration as a function
>>>>>> of time to investigate the
>>>>>> equilibration of a system of
>>>>>> hyperbranched polymer. However,
>>>>>> the plot goes down and does not
>>>>>> fluctuate around zero.
>>>>>>
>>>>>> Could anybody tell me what is the
>>>>>> reason and what should I do?
>>>>>> Any help will be appreciated.
>>>>>>
>>>>>> Regards,
>>>>>> Farideh
>>>>>>
>>>>>>
>>>>>> --
>>>>>> Brian Radak
>>>>>> Postdoctoral Appointee
>>>>>> Leadership Computing Facility
>>>>>> Argonne National Laboratory
>>>>>>
>>>>>> 9700 South Cass Avenue, Bldg. 240
>>>>>> Argonne, IL 60439-4854
>>>>>> (630) 252-8643
>>>>>> brian.radak_at_anl.gov
>>>>>> <mailto:brian.radak_at_anl.gov>
>>>>>>
>>>>>>
>>>>>
>>>>> --
>>>>> Brian Radak
>>>>> Postdoctoral Appointee
>>>>> Leadership Computing Facility
>>>>> Argonne National Laboratory
>>>>>
>>>>> 9700 South Cass Avenue, Bldg. 240
>>>>> Argonne, IL 60439-4854
>>>>> (630) 252-8643
>>>>> brian.radak_at_anl.gov
>>>>> <mailto:brian.radak_at_anl.gov>
>>>>>
>>>>>
>>>>
>>>> --
>>>> Brian Radak
>>>> Postdoctoral Appointee
>>>> Leadership Computing Facility
>>>> Argonne National Laboratory
>>>>
>>>> 9700 South Cass Avenue, Bldg. 240
>>>> Argonne, IL 60439-4854
>>>> (630) 252-8643
>>>> brian.radak_at_anl.gov <mailto:brian.radak_at_anl.gov>
>>>>
>>>
>>> --
>>> Brian Radak
>>> Postdoctoral Appointee
>>> Leadership Computing Facility
>>> Argonne National Laboratory
>>>
>>> 9700 South Cass Avenue, Bldg. 240
>>> Argonne, IL 60439-4854
>>> (630) 252-8643
>>> brian.radak_at_anl.gov <mailto:brian.radak_at_anl.gov>
>>>
>>
>> --
>> Brian Radak
>> Postdoctoral Appointee
>> Leadership Computing Facility
>> Argonne National Laboratory
>>
>> 9700 South Cass Avenue, Bldg. 240
>> Argonne, IL 60439-4854
>> (630) 252-8643
>> brian.radak_at_anl.gov <mailto:brian.radak_at_anl.gov>
>>
>>
>
> --
> Brian Radak
> Postdoctoral Scholar
> Gordon Center for Integrative Science, W323A
> Department of Biochemistry & Molecular Biology
> University of Chicago
> 929 E. 57th St.
> Chicago, IL 60637-1454
> Tel: 773/834-2812
> email: radak_at_uchicago.edu <mailto:radak_at_uchicago.edu>
>

-- 
Brian Radak
Postdoctoral Scholar
Gordon Center for Integrative Science, W323A
Department of Biochemistry & Molecular Biology
University of Chicago
929 E. 57th St.
Chicago, IL 60637-1454
Tel: 773/834-2812
email: radak_at_uchicago.edu

This archive was generated by hypermail 2.1.6 : Sun Dec 31 2017 - 23:20:40 CST