Re: Problems with temperature control.

From: Josep Ivan Balaguer Molins (jobamo1_at_etsid.upv.es)
Date: Tue Jul 31 2018 - 06:37:49 CDT

Hell again Brian, I think the problem is using the reassign
temperature control. I would like to try doing the same with Langvine
if posible:

    for { set TEMP 600 } { $TEMP >= 300 } { incr TEMP -10 } {

  langevin on ;# do langevin dynamics
  langevinDamping 2 ;# damping coefficient (gamma) of 1/ps
  langevinTemp TEMP
  langevinHydrogen no ;# don't couple langevin bath to hydrogens
  langevinPiston on
  langevinPistonTarget 1.01325 ;# in bar -> 1 atm
  langevinPistonPeriod 200
  langevinPistonDecay 150
  langevinPistonTemp TEMP
}

This doesn't work, how should I code it? I don't see any example on
internet nor in the manual.

Best regards,

Ivan
Quoting Josep Ivan Balaguer Molins <jobamo1_at_etsid.upv.es>:

> Thank you very much Brian. I've turned off the pressure control. And
> I've changed the temperature increment from -10 to -0.1 and even to
> -0.001:
>
> run 1000
>
> for { set TEMP 600 } { $TEMP >= 300 } { incr TEMP -0.1 } {
> reassignTemp $TEMP
>
> run 1000
> }
>
> reassignFreq 100
>
> In any case I get the warning "low global exclusion count" as you
> expected before. And the simulation stops:
>
> ERROR: Constraint failure in RATTLE algorithm for atom 1347!
> ERROR: Constraint failure; simulation has become unstable.
> ERROR: Constraint failure in RATTLE algorithm for atom 9897!
> ERROR: Constraint failure; simulation has become unstable.
> ERROR: Constraint failure in RATTLE algorithm for atom 12147!
> ERROR: Constraint failure; simulation has become unstable.
> ERROR: Constraint failure in RATTLE algorithm for atom 8997!
> ERROR: Constraint failure; simulation has become unstable.
>
> I've read this:
> http://www.ks.uiuc.edu/Research/namd/wiki/index.cgi?NamdTroubleshooting
>
> I've increased the margin from 4 to even 20 without succes. The same
> atoms become unestable. I've decreased the reassingFreq but it
> doesn't solve the problem either.
>
> Any idea?
>
> Thank you very much.
>
> Ivan
>
>
>
>
> Quoting Brian Radak <brian.radak_at_gmail.com>:
>
>> Sorry, I meant reassignFreq. In any event, that is not what is apparently
>> causing the error.
>>
>> I'm actually not sure why you would get a low global bond count (I might
>> have expected a low global exclusion count). My guess is the system is
>> cooling too quickly and thus pressure control is causing issues. In
>> retrospect, the pressure control is probably not a good idea in this case,
>> since your temperature change is very large.
>>
>> On Mon, Jul 30, 2018 at 10:41 AM, Josep Ivan Balaguer Molins <
>> jobamo1_at_etsid.upv.es> wrote:
>>
>>> Sorry, Brian but I don't understand... Here is the code were I get the
>>> error.
>>>
>>>
>>> #############################################################
>>> ## ADJUSTABLE PARAMETERS ##
>>> #############################################################
>>> gromacs on
>>> grotopfile 32x50LA_crystal.top
>>> grocoorfile 32x50LA_crystal.gro
>>> bincoordinates testD2.coor
>>>
>>> paraTypeCharmm on
>>> exclude scaled1-4
>>> 1-4scaling 0.5
>>>
>>> switching on
>>> switchdist 8 # at 8Å we start to smooth electrostatic to 0
>>> cutoff 12.0
>>> pairlistdist 26.0
>>> margin 4
>>>
>>> # Integrator Parameters
>>>
>>>
>>>
>>> rigidBonds all ;# needed for 2fs steps
>>> nonbondedFreq 1
>>> vdwGeometricSigma yes
>>> fullElectFrequency 2
>>> stepspercycle 20
>>> pairlistsperCycle 2
>>>
>>> wrapAll on
>>>
>>> timestep 1.0 ;# 2fs/step 2 en simulación, 1 minimizacion
>>>
>>> outputEnergies 100
>>> outputtiming 100
>>> outputPressure 100
>>> binaryoutput yes
>>> outputname testD3
>>> dcdfreq 200
>>> restartfreq 1000
>>> restartname rest_TestD3
>>>
>>> binvelocities testD2.vel
>>>
>>> reassignFreq 100
>>> run 80
>>> for { set TEMP 550 } { $TEMP >= 300 } { incr TEMP -10 } {
>>> reassignTemp $TEMP
>>> run 100
>>> }
>>>
>>>
>>>
>>>
>>> # Constant Temperature Control
>>>
>>> langevin off ;# do langevin dynamics
>>> langevinDamping 2 ;# damping coefficient (gamma) of 1/ps
>>> langevinTemp 300
>>> langevinHydrogen no ;# don't couple langevin bath to hydrogens
>>> langevinPiston on
>>> langevinPistonTarget 1.01325 ;# in bar -> 1 atm
>>> langevinPistonPeriod 10000
>>> langevinPistonDecay 15000
>>> langevinPistonTemp 300
>>>
>>> useGroupPressure yes ;# needed for rigidBonds
>>> useFlexibleCell no
>>> useConstantArea no
>>>
>>> extendedSystem testD2.xsc
>>>
>>> run 300000
>>>
>>>
>>> I don't have a reassignIncr here... or do you mean I have to make it
>>> smaller on the previous minimization and simulation?
>>> If I turn langevin off, tha same error appears.
>>>
>>> Thank you very much for your answer.
>>>
>>> Best regards
>>>
>>> Ivan
>>>
>>>
>>> Quoting Brian Radak <brian.radak_at_gmail.com>:
>>>
>>> Note that the Langevin and reassign thermostats are completely decoupled
>>>> code paths. My guess is that the reassignIncr is too large compared to the
>>>> relaxation rate dictated by Langevin damping and the two are redundant -
>>>> personally I would only use one or the other (you can also script changes
>>>> to langevinTemp). You could also try using a smaller damping coefficient,
>>>> but this really isn't any different from not using Langevin at all.
>>>>
>>>>
>>>> On Mon, Jul 30, 2018 at 6:21 AM, Josep Ivan Balaguer Molins <
>>>> jobamo1_at_etsid.upv.es> wrote:
>>>>
>>>> Hi everyone!
>>>>>
>>>>> I have a system I've minimized and during minimization, it's been heat up
>>>>> to 600K using this code lines:
>>>>>
>>>>> temperature 0
>>>>>
>>>>> reassignFreq 1000
>>>>> reassignIncr 10
>>>>> reassignHold 600 # Final temperature
>>>>>
>>>>> langevin on ;# do langevin dynamics
>>>>> langevinDamping 2 ;# damping coefficient (gamma) of 1/ps
>>>>> langevinTemp 600 # Final temperature
>>>>> langevinHydrogen no ;# don't couple langevin bath to hydrogens
>>>>>
>>>>> langevinPiston on
>>>>> langevinPistonTarget 1.01325 ;# in bar -> 1 atm
>>>>> langevinPistonPeriod 200
>>>>> langevinPistonDecay 100
>>>>> langevinPistonTemp 600 # Final temperature
>>>>>
>>>>> It does work, I can see the temperature converging to 600K in VMD.
>>>>>
>>>>> After this, I perform an Isothermal simulation to see if the systems
>>>>> keeps
>>>>> the stability and works too:
>>>>>
>>>>> binvelocities testD1.vel # Use velocities from minimization output.
>>>>>
>>>>> langevin on ;# do langevin dynamics
>>>>> langevinDamping 2 ;# damping coefficient (gamma) of 1/ps
>>>>>
>>>>> langevinTemp 600 # I keep the 600K
>>>>>
>>>>> langevinHydrogen no ;# don't couple langevin bath to hydrogens
>>>>>
>>>>> useGroupPressure yes ;# needed for rigidBonds
>>>>> useFlexibleCell no
>>>>> useConstantArea no
>>>>>
>>>>> langevinPiston on
>>>>> langevinPistonTarget 1.01325 ;# in bar -> 1 atm
>>>>> langevinPistonPeriod 200
>>>>> langevinPistonDecay 100
>>>>> langevinPistonTemp 600 # I keep the 600K
>>>>>
>>>>> extendedSystem testD1.xsc
>>>>>
>>>>> run 400000
>>>>>
>>>>> Seeing the result in VMD, the code works well and the simulation has no
>>>>> errors, the temperature is stable at 600K.
>>>>>
>>>>> Finally to perform the final simulation, I have to cool the system down
>>>>> to
>>>>> 300K. So I use the output files and with Langevin I am able to cool down
>>>>> the sistem without problems to 300K. But seeing in VMD the temperature
>>>>> graph, the temperature goes down so fast. In 5000 integration steps goes
>>>>> from 600K to 300K. So I need to change the temperature graph. To do so, I
>>>>> use this code:
>>>>>
>>>>> binvelocities testD2.vel
>>>>>
>>>>> reassignFreq 100
>>>>>
>>>>> run 1000
>>>>>
>>>>> for { set TEMP 550 } { $TEMP >= 300 } { incr TEMP -10 } {
>>>>> reassignTemp $TEMP
>>>>>
>>>>> run 100
>>>>> }
>>>>>
>>>>> And I get the following error: FATAL ERROR: Bad global bond count! (14395
>>>>> vs 14400)
>>>>>
>>>>> How is that possible taking into account that before everything worked?
>>>>>
>>>>> Changing cutoff and pairlistdist parameters, I got it down to 14397
>>>>> instead of 14395. But I can get any better. And I don't understand why is
>>>>> this happening.
>>>>>
>>>>> Thank you very much.
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>
>>>

This archive was generated by hypermail 2.1.6 : Mon Dec 31 2018 - 23:21:20 CST