Re: Problems with temperature control.

From: Josep Ivan Balaguer Molins (jobamo1_at_etsid.upv.es)
Date: Tue Jul 31 2018 - 10:54:40 CDT

Hi Brian, sorry for my stupidity and thank you very much for being so
helpful. This is how I've implemented the rescalevels, the same error
of unestable system appears. I don't understand why...

Is it correct? Or is it wrong?

Thank you very much.

binvelocities testD2.vel

    run 100

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

    reassignTemp $TEMP

    rescalevels sqrt($TEMP+300)

    run 10
    }

reassignFreq 40

# Constant Temperature Control

  langevin off ;# do langevin dynamics

Quoting Brian Radak <brian.radak_at_gmail.com>:

> You should only need to adjust langevinTemp in order to accomplish
> temperature control. I might also recommend rescaling the velocities
> (rescalevels <scale factor>) at each interval by the sqrt ratio of the old
> and new temperatures as is done for simulated tempering, but this is
> probably overkill for a simple annealing run.
>
> Again, pressure control is not your friend when the temperature changes
> quickly. The exclusion errors you are seeing almost always come up when the
> simulation cell changes rapidly for any reason.
>
> On Tue, Jul 31, 2018, 7:37 AM Josep Ivan Balaguer Molins <
> jobamo1_at_etsid.upv.es> wrote:
>
>> 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

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