[gmx-users] Packmol starting .pdb and proper equilibration
Justin Lemkul
jalemkul at vt.edu
Thu Oct 30 02:04:49 CET 2014
On 10/29/14 2:07 PM, Ioanna Styliari wrote:
> Dear Gromacs users,
>
> I am simulating a system that contains a drug nanoparticle in amorphous state, surrounded by water and the water is also surrounded by acetone (imagine 2 concentric spheres for the first two inside one big box for the acetone).
>
> I used to build the system using gromacs editconf and more precisely building one box at a time and running em, nvt and npt equilibration for each one of the extra box-layers I was adding. However this was a very tedious procedure so I have started recently building the system using Packmol. Once I get the pdb from packmol I edit it to add the box dimensions with editconf and then go through an em, nvt, npt (usually 100ps) for the whole system.
> I have already run 4-5 simulations of 50ns and 100ns and they were running fine (they are running on parallel in a cluster). However in my last one I got the following error:
>
> Fatal error:
> 1 particles communicated to PME node 20 are more than 2/3 times the cut-off out of the domain decomposition cell of their charge group in dimension y.
>
> I know this means that the system is not well equilibrated but after the nvt and npt the T and P are converged. I have reduced the time step to 1fs but I am also worried that the equilibration is not the correct one. It is also worth saying that I had similar errors even with my previous building procedure using gromacs.
>
What is your concern regarding the time step? This is a common strategy for
temperamental systems.
> Any suggestions/experiences with using Packmol? I am also attaching my mdp options below just in case.
>
> Thank you very much for your time in advance,
>
> Ioanna Styliari
>
> EM:
> integrator = steep ; Algorithm (steep = steepest descent minimization)
> emtol = 10.0 ; Stop minimization when the maximum force < 10.0 kJ/mol/nm
> emstep = 0.001 ; Energy step size
> nsteps = 500000 ; Maximum number of (minimization) steps to perform
>
> ; Parameters describing how to find the neighbours of each atom and how to calculate the interactions
> Nstlist = 1 ; Frequency to update the neighbour list and long range forces
> ns_type = grid ; Method to determine neighbour list (simple, grid)
> rlist = 1.0 ; Cut-off for making neighbour list (short range forces)
> coulombtype = PME ; Treatment of long range electrostatic interactions
> rcoulomb = 1.0 ; Short-range electrostatic cut-off
> vdw-type = Cut-off
> rvdw = 1.0 ; Short-range Van der Waals cut-off
> pbc = xyz ; Periodic Boundary Conditions
> NVT
> title = pol1-box NVT
> ; Run parameters
> integrator = md ;
> nsteps = 100000 ; 100 ps
> dt = 0.001 ; 1 fs
> ; Bond parameters
> continuation = no ; first dynamics run
> ; Neighborsearching
> ns_type = grid ; search neighboring grid cells
> nstlist = 5 ;
> rlist = 1.0 ; short-range neighborlist cutoff (in nm)
> rcoulomb = 1.0 ; short-range electrostatic cutoff (in nm)
> rvdw = 1.0 ; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
> pme_order = 4 ; cubic interpolation
> fourierspacing = 0.16 ; grid spacing for FFT
> ; Temperature coupling is on
> tcoupl = V-rescale ; modified Berendsen thermostat
> tc-grps = system ;
> tau_t = 0.1 ; time constant, in ps
> ref_t = 300 ; reference temperature, one for each group, in K
> ; Pressure coupling is off
> pcoupl = no ; no pressure coupling in NVT
> ; Periodic boundary conditions
> pbc = xyz ; 3-D PBC
> ; Dispersion correction
> DispCorr = EnerPres ; account for cut-off vdW scheme
> ; Velocity generation
> gen_vel = yes ; assign velocities from Maxwell distribution
> gen_temp = 300 ; temperature for Maxwell distribution
> gen_seed = -1 ; generate a random seed
>
> constraints = none
>
> NPT
> title = pol1-box NPT
> ; Run parameters
> integrator = md
> nsteps = 100000
> dt = 0.001
> ; Bond parameters
> continuation = yes ; Restarting after NVT
> ; Neighborsearching
> ns_type = grid ; search neighbouring grid cells
> nstlist = 5 ; fs
> rlist = 1.0 ; short-range neighbourlist cutoff (in nm)
> rcoulomb = 1.0 ; short-range electrostatic cutoff (in nm)
> rvdw = 1.0 ; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
> pme_order = 4 ; cubic interpolation
> fourierspacing = 0.16 ; grid spacing for FFT
> ; Temperature coupling is on
> tcoupl = V-rescale ; modified Berendsen thermostat
> tc-grps = system ;
> tau_t = 0.1 ; time constant, in ps
> ref_t = 300 ; reference temperature, one for each group, in K
> ; Pressure coupling is on
> pcoupl = berendsen ; Pressure coupling on in NPT
> pcoupltype = isotropic ; uniform scaling of box vectors
> tau_p = 2.0 ; time constant, in ps
> ref_p = 1.0 ; reference pressure, in bar
> compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1
> refcoord_scaling = com
> ; Periodic boundary conditions
> pbc = xyz ; 3-D PBC
> ; Dispersion correction
> DispCorr = EnerPres ; account for cut-off vdW scheme
> ; Velocity generation
> gen_vel = no ; Velocity generation is off
>
>
> and then MD:
>
> title = pol1-box first 2 ns MD
> ; Run parameters
> integrator = md ;
> nsteps = 2000000 ;
> dt = 0.001 ; fs
> ; Bond parameters
> continuation = yes ; Restarting after NPT
> ; Neighborsearching
> ns_type = grid ; search neighboring grid cells
> nstlist = 5 ; fs
> rlist = 1.0 ; short-range neighborlist cutoff (in nm)
> rcoulomb = 1.0 ; short-range electrostatic cutoff (in nm)
> rvdw = 1.0 ; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics
> pme_order = 4 ; cubic interpolation
> fourierspacing = 0.16 ; grid spacing for FFT
> ; Temperature coupling is on
> tcoupl = V-rescale ; modified Berendsen thermostat
> tc-grps = system ; two coupling groups - more accurate
> tau_t = 0.1 ; time constant, in ps
> ref_t = 300 ; reference temperature, one for each group, in K
> ; Pressure coupling is on
> pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT
> pcoupltype = isotropic ; uniform scaling of box vectors
> tau_p = 2.0 ; time constant, in ps
> ref_p = 1.0 ; reference pressure, in bar
> compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1
> ; Periodic boundary conditions
> pbc = xyz ; 3-D PBC
> ; Dispersion correction
> DispCorr = EnerPres ; account for cut-off vdW scheme
> ; Velocity generation
> gen_vel = no ; Velocity generation is off
>
So in reality there is no difference between NPT and the production run, e.g. no
restraints or anything. At what point dos the crash occur? What force field
are you using?
-Justin
--
==================================================
Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow
Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201
jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul
==================================================
More information about the gromacs.org_gmx-users
mailing list