Re: periodic boundary condition in XY, hard wall in Z - errors

From: Olya Kravchenko (ovkrav_at_gmail.com)
Date: Fri Mar 11 2016 - 14:57:55 CST

Maxim,

thank you for replying. I am trying to figure out MSM right now and
see if it makes difference, I'll try setting exact values for
cellbasisvectors too.

As for the hard wall: I expected to see expansion in the beginning
that does not change further, as it happens when I use PBC in 3D. What
I saw is oscillating continuous expansion, i.e. a couple of times
during simulation the whole box shrinks and then expands more than it
did previous time; the change in proportion (elongation in z) is
continuous.

I ran simulations in NVT today just to see if it helps me keep the
size of the box fixed and I still see the change in proportions and
increase in Z, actually the whole box becomes much bigger, and
proportions change as well, although not as much as it happens when I
run NPT.

Here is my configuration file for NVT:

# Minimization and Equilibration of argon with PBC in XY

structure ar.psf
coordinates ar.pdb
set outputname ar_force_out

firsttimestep 0

set temperature 310

# Input
paraTypeCharmm on
parameters My_parameters.rtf
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 2.0
rigidBonds none
nonbondedFreq 1
fullElectFrequency 2
stepspercycle 10

# Constant Temperature Control
langevin on
langevinDamping 5
langevinTemp $temperature
langevinHydrogen off

# Periodic Boundary Conditions
cellBasisVector1 62.0 0. 0.0
cellBasisVector2 0.0 62.0 0.0
#cellBasisVector3 0.0 0 0.0
cellOrigin 30.0 30.0 32.5

wrapAll on

margin 30

# Constant Pressure Control (variable volume)
useGroupPressure no ;# needed for rigidBonds
useFlexibleCell no
useConstantArea no

langevinPiston on
langevinPistonTarget 1.01325 ;# in bar -> 1 atm
langevinPistonPeriod 100.0
langevinPistonDecay 50.0
langevinPistonTemp $temperature

tclForces on
tclForcesScript my_force.tcl

outputName $outputname

restartfreq 500
dcdfreq 250
xstFreq 250
outputEnergies 100
outputPressure 100

minimize 1000
run 50000

On Fri, Mar 11, 2016 at 3:34 PM, Maxim Belkin <mbelkin_at_ks.uiuc.edu> wrote:
> It is better to set cellbasisvector[1,2] to exact values.
> Because you are not using MSM (and you can’t use PME with non-periodic systems), interactions between atoms in the system aren’t right (they are cut off at 12 AA). You should look into MSM first, it should not be that difficult.
>
> Now, the fact that your simulation crash ~150,000 steps might indicate a problem that emerges due to multiple time-stepping. Try: "stepspercycle 2".
>
> Finally, what do you expect to get from an NPT simulation with vacuum along z? How do you expect langevinPiston to work in this situation?
>
>
>
>> On Mar 10, 2016, at 20:01, Olya Kravchenko <ovkrav_at_gmail.com> wrote:
>>
>> Thank you! I will study MSM method.
>>
>> I have changed the size of cell vectors back to 62, that is the actual
>> size rounded up, and I commented out cellBasisVector3. I ran a number
>> of simulations today and I keep seeing the same error when I increase
>> the number of timesteps. The margin size does seem to make the job run
>> longer, at least I don't see any other correlation. My last run quit
>> near 150 000 timesteps (I set it up for 500K timesteps).
>>
>> Apart from that, I see in vmd that my simulation box actually becomes
>> longer in Z-direction, it is a cube in the beginning and then the
>> proportion changes. I am not sure how to interpret it, but it looks
>> like the size in z does not hold, certainly no hard wall. I wonder if
>> there is a way to set up the hard wall in a simple way?
>>
>> Olga
>>
>>
>>
>> On Thu, Mar 10, 2016 at 4:58 PM, Maxim Belkin <mbelkin_at_ks.uiuc.edu> wrote:
>>> Olya,
>>>
>>> If your system size is 60A in x and y, you should use vectors of length 60 as your "cellbasisvector[1,2]". You are currently using 70, which means you have vacuum in x and y directions and this is not what you want. To simulate a system that is non-periodic along Z with NAMD you have to use Multi-level summation method (look for MSM option in NAMD User’s Guide). I believe you have to remove "cellbasisvector3" instead of specifying 0 0 0 but I’m not 100% positive on that.
>>>
>>> Maxim
>>>
>>>
>>>> On Mar 10, 2016, at 11:21, Olya Kravchenko <ovkrav_at_gmail.com> wrote:
>>>>
>>>> My system consists of 200 argon atoms in a box, PBC in X and Y, I
>>>> would like hard wall in Z and my goal is to measure
>>>> density/concentration. So I think that NPT is reasonable in my case,
>>>> is that correct?
>>>>
>>>> I commented out PME, because I read that it will speed up
>>>> calculations, also because, if I understand it correctly, it requires
>>>> all three cell basis vectors to be defined.
>>>>
>>>> linaccel is commented because I defined my own force in a separate
>>>> file for the force.
>>>>
>>>> Here is my configuration file:
>>>>
>>>> #############################################################
>>>> ## JOB DESCRIPTION ##
>>>> #############################################################
>>>>
>>>> # Equilibration of argon atoms with PBC in XY direction
>>>>
>>>> #############################################################
>>>> ## ADJUSTABLE PARAMETERS ##
>>>> #############################################################
>>>>
>>>> structure ar.psf
>>>> coordinates ar.pdb
>>>>
>>>> set temperature 310
>>>> set outputname ar_force_out
>>>>
>>>> firsttimestep 0
>>>>
>>>>
>>>> #############################################################
>>>> ## SIMULATION PARAMETERS ##
>>>> #############################################################
>>>>
>>>> # Input
>>>> paraTypeCharmm on
>>>> parameters My_parameters.rtf
>>>> 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 2.0 ;# 2fs/step
>>>> rigidBonds all ;# needed for 2fs steps
>>>> nonbondedFreq 1
>>>> fullElectFrequency 2
>>>> stepspercycle 10
>>>>
>>>>
>>>> # Constant Temperature Control
>>>> langevin on ;# do langevin dynamics
>>>> langevinDamping 1 ;# damping coefficient (gamma) of 1/ps
>>>> langevinTemp $temperature
>>>> langevinHydrogen off ;# don't couple langevin bath to hydrogens
>>>>
>>>>
>>>> # Periodic Boundary Conditions
>>>> cellBasisVector1 62.0 0. 0.0
>>>> cellBasisVector2 0.0 62.0 0.0
>>>> cellBasisVector3 0.0 0 0.0
>>>> cellOrigin 30.0 30.0 32.5
>>>>
>>>> wrapAll on
>>>>
>>>> margin 10
>>>>
>>>> # PME (for full-system periodic electrostatics)
>>>> #PME yes
>>>> #PMEGridSpacing 1.0
>>>>
>>>>
>>>> # Constant Pressure Control (variable volume)
>>>> useGroupPressure yes ;# needed for rigidBonds
>>>> useFlexibleCell no
>>>> useConstantArea no
>>>>
>>>> langevinPiston on
>>>> langevinPistonTarget 1.01325 ;# in bar -> 1 atm
>>>> langevinPistonPeriod 100.0
>>>> langevinPistonDecay 50.0
>>>> langevinPistonTemp $temperature
>>>>
>>>> ##################
>>>> ###DEFINE FORCES##
>>>> ##################
>>>>
>>>> tclForces on
>>>> tclForcesScript my_force.tcl
>>>> #set linaccel "30 0 0"
>>>>
>>>> # Output
>>>> outputName $outputname
>>>>
>>>> restartfreq 500 ;# 500steps = every 1ps
>>>> dcdfreq 250
>>>> xstFreq 250
>>>> outputEnergies 100
>>>> outputPressure 100
>>>>
>>>>
>>>> # Minimization
>>>> minimize 100
>>>> reinitvels $temperature
>>>>
>>>> run 100000
>>>>
>>>> Thanks!!
>>>>
>>>> Olga
>>>>
>>>> On Thu, Mar 10, 2016 at 2:05 AM, Norman Geist
>>>> <norman.geist_at_uni-greifswald.de> wrote:
>>>>> We need to see your full input script, since it looks you're doing npt?
>>>>>
>>>>> Norman Geist
>>>>>
>>>>>> -----Ursprüngliche Nachricht-----
>>>>>> Von: owner-namd-l_at_ks.uiuc.edu [mailto:owner-namd-l_at_ks.uiuc.edu] Im
>>>>>> Auftrag von Olya Kravchenko
>>>>>> Gesendet: Donnerstag, 10. März 2016 04:23
>>>>>> An: NAMD <namd-l_at_ks.uiuc.edu>
>>>>>> Betreff: namd-l: periodic boundary condition in XY, hard wall in Z - errors
>>>>>>
>>>>>> Hi all,
>>>>>>
>>>>>> I would like to apply periodic boundary conditions only in X and Y
>>>>>> directions. In my config file I used the following entry:
>>>>>>
>>>>>> # Periodic Boundary Conditions
>>>>>> cellBasisVector1 70.0 0. 0.0
>>>>>> cellBasisVector2 0.0 70.0 0.0
>>>>>> cellBasisVector3 0.0 0 0.0
>>>>>> cellOrigin 30.0 30.0 32.5
>>>>>>
>>>>>> wrapAll on
>>>>>> margin 5
>>>>>>
>>>>>>
>>>>>> When I submit the job it runs for a while until I see the following error:
>>>>>>
>>>>>> WRITING COORDINATES TO DCD FILE ar_force_out.dcd AT STEP 20750
>>>>>> FATAL ERROR: Periodic cell has become too small for original patch grid!
>>>>>> Possible solutions are to restart from a recent checkpoint,
>>>>>> increase margin, or disable useFlexibleCell for liquid simulation.
>>>>>>
>>>>>> How can I address this error? Is this the right way to set up PBC in 2
>>>>>> dimensions?
>>>>>>
>>>>>> My simulation cell is actually around (~61 ~61 ~62), as measured by
>>>>>> minmax. I played with numbers and noticed that my simulation runs
>>>>>> longer if I increase both the margin and cell basis vectors.
>>>>>>
>>>>>> I saw some of the older threads in the archives where the margin was
>>>>>> discussed but it's not clear how large should it be. If I want to run
>>>>>> my simulation for 500 000 steps, what margin should I set up?
>>>>>>
>>>>>> Also, when choosing cell basis vector, how much bigger should it be
>>>>>> than what is measured by "measure minmax"?
>>>>>
>>>
>

This archive was generated by hypermail 2.1.6 : Tue Dec 27 2016 - 23:21:52 CST