[gmx-users] pbc problem

[jkrieger at mrc-lmb.cam.ac.uk]

I usually use -pbc nojump for my protein simulations and this works every
time.

> Dear gromacs users
>
> My system contains DOPC + CHOLESTEROLO + WATER + drug molecules in a
> rectangular box.
>
> I put drug molecule in 2 position: a) drug in the center of bilayer
> membrane, b) drug inside water molecules in top leaflet.
>
> For both positions, I did energy minimization successfully with following
> mdp file.
> --------------------------------------------------------------------------------------
> ; Parameters describing what to do, when to stop and what to save
> integrator    = steep        ; Algorithm (steep = steepest descent
> minimization)
> emtol        = 1000.0      ; Stop minimization when the maximum force <
> 1000.0 kJ/mol/nm
> emstep      = 0.01      ; Energy step size
> nsteps        = 50000          ; Maximum number of (minimization) steps to
> perform
>
> ; Parameters describing how to find the neighbors of each atom
> nstlist        = 1            ; Frequency to update the neighbor list and
> long range forces
> ns_type        = grid        ; Method to determine neighbor list (simple,
> grid)
> rlist        = 1.2        ; Cut-off for making neighbor list (short range
> forces)
> coulombtype    = PME        ; Treatment of long range electrostatic
> interactions
> rcoulomb    = 1.2        ; Short-range electrostatic cut-off
> rvdw        = 1.2        ; Short-range Van der Waals cut-off
> pbc            = xyz         ; Periodic Boundary Conditions
> ---------------------------------------------------------------------------------------
> After energy minimization, I saw obtained file (em.gro) by VMD. All things
> were true and intact.
>
> For both positions, I did equilibration in NPT ensemble with following mdp
> file.
> ---------------------------------------------------------------------------------------
> ; Run parameters
> integrator    = md        ; leap-frog integrator
> nsteps        = 250000    ; 2 * 500000 = 1000 ps (1 ns)
> dt            = 0.002        ; 2 fs
> ; Output control
> nstxout        = 100        ; save coordinates every 0.2 ps
> nstvout        = 100        ; save velocities every 0.2 ps
> nstxtcout   = 100        ; xtc compressed trajectory output every 2 ps
> nstenergy    = 100        ; save energies every 0.2 ps
> nstlog        = 100        ; update log file every 0.2 ps
> energygrps  = CHOL DOPC drg SOL
> ; Bond parameters
> continuation    = no            ; Restarting after NVT
> constraint_algorithm = lincs    ; holonomic constraints
> constraints    = all-bonds            ; all bonds (even heavy atom-H
> bonds)
> constrained
> lincs_iter    = 1                    ; accuracy of LINCS
> lincs_order    = 4                    ; also related to accuracy
> ; Neighborsearching
> ns_type        = grid        ; search neighboring grid cels
> nstlist        = 5            ; 10 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            ; More accurate thermostat
> tc-grps        = CHOL_DOPC    drg SOL    ; three coupling groups - more
> accurate
> tau_t        = 0.5    0.5       0.5       ; time constant, in ps
> ref_t        = 323     323       323     ; reference temperature, one for
> each group, in K
> ; Pressure coupling is on
> pcoupl        = Parrinello-Rahman        ; Pressure coupling on in NPT
> pcoupltype    = semiisotropic            ; uniform scaling of x-y box
> vectors, independent z
> tau_p        = 5.0                    ; time constant, in ps
> ref_p        = 1.0    1.0                ; reference pressure, x-y, z (in
> bar)
> compressibility = 4.5e-5    4.5e-5    ; isothermal compressibility, bar^-1
> ; 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    = 323        ; temperature for Maxwell distribution
> gen_seed    = -1        ; generate a random seed
> ; COM motion removal
> ; These options remove motion of the protein/bilayer relative to the
> solvent/ions
> nstcomm         = 1
> comm-mode       = Linear
> comm-grps       = CHOL_DOPC_drg  SOL
> ; Scale COM of reference coordinates
> refcoord_scaling = com
>
> ---------------------------------------------------------------------------------------
> For 2 positions, I chechked tempreture and pressure fluctuation and box
> dimention during equilibration. All things were good. When I saw
> trajectory
> by VMD (npt.gro and npt xtc), I had pbc problem (some atoms leave box and
> enter the box in opposit direction).
>
> For position (a): I corrected pbc problem by
>
> trjconv -f npt.xtc -s npt.tpr -n index.ndx -o 2npt.xtc -pbc mol -center
>
> I selected CHOL_DOPC-drg group for centering. So problem was solved,
> approximately.
>
> For position (b) in which drug molecule is not in the center of lipid, I
> can not use -center or I can not use group containing drug molecule for
> centering.
>
> How to fix this problem.
>
> Any help will highly appreciated.
> --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
> * Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>