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

From: faride badalkhani (farideh.khamseh_at_gmail.com)
Date: Tue Aug 23 2016 - 14:16:28 CDT

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?

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?

Your help is really appreciated.

Regards,
Farideh

On Aug 23, 2016 11:29 PM, "Brian Radak" <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> 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>
>> 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> 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> 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> 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-dJEEx0SBl
>>>>>> U2Fuekma?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> 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
>>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> 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
>>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> 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
>>>>>
>>>>
>>>> --
>>>> 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
>>>>
>>>
>> --
>> 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
>>
>
>
> --
> 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 : Tue Dec 27 2016 - 23:22:24 CST