[gmx-users] Changes in the simulation box after the production run

[Mahsa E]

Dear gmx-users,

I did a 200 ns production md run in NVT ensemble for a simulation box of
polymer chains. Before this step, I did the energy minimisation, NVT and
NPT equilibration on the system. The problem is after the production run, I
don't get the initial equilibrated packed box of polymer and it seems more
like a circular shape with some parts of the chains out of the box. What is
the reason for getting this result?
For the MD run I used the mdp file below:

; 7.3.2 Preprocessing

;define                  =       ; defines to pass to the preprocessor


; 7.3.3 Run Control

integrator              = md                    ; md integrator

tinit                   = 0                     ; [ps] starting time for run

dt                      = 0.002                 ; [ps] time step for
integration

nsteps                  = 100000000                ; maximum number of
steps to integrate, 0.002 * 100000000 = 200000 ps

comm_mode               = Linear                ; remove center of mass
translation

nstcomm                 = 100                     ; [steps] frequency of
mass motion removal

;comm_grps               = Protein Non-Protein   ; group(s) for center of
mass motion removal


; 7.3.8 Output Control

nstxout                 = 0         ; [steps] freq to write coordinates to
trajectory

nstvout                 = 0         ; [steps] freq to write velocities to
trajectory

nstfout                 = 0         ; [steps] freq to write forces to
trajectory

nstlog                  = 1000           ; [steps] freq to write energies
to log file

nstenergy               = 1000           ; [steps] freq to write energies
to energy file

nstxtcout               = 1000           ; [steps] freq to write
coordinates to xtc trajectory

xtc_precision           = 1000          ; [real] precision to write xtc
trajectory

xtc_grps                = System        ; group(s) to write to xtc
trajectory

energygrps              = System        ; group(s) to write to energy file

cutoff-scheme            = verlet

; 7.3.9 Neighbor Searching

nstlist                 = 20             ; [steps] freq to update neighbor
list

ns_type                 = grid          ; method of updating neighbor list

pbc                     = xyz           ; periodic boundary conditions in
all directions

rlist                   = 0.8           ; [nm] cut-off distance for the
short-range neighbor list


; 7.3.10 Electrostatics

coulombtype             = PME           ; Particle-Mesh Ewald electrostatics

rcoulomb                = 1.2           ; [nm] distance for Coulomb cut-off


; 7.3.11 VdW

vdwtype                 = cut-off       ; twin-range cut-off with rlist
where rvdw >= rlist

vdw-modifier             = potential-switch

rvdw-switch              = 1.1

rvdw                    = 1.2           ; [nm] distance for LJ cut-off

DispCorr                = EnerPres      ; apply long range dispersion
corrections


; 7.3.13 Ewald

fourierspacing          = 0.12          ; [nm] grid spacing for FFT grid
when using PME

pme_order               = 4             ; interpolation order for PME, 4 =
cubic

ewald_rtol              = 1e-6          ; relative strength of
Ewald-shifted potential at rcoulomb

ewald-rtol-lj            = 0.001

lj-pme-comb-rule         = Geometric

ewald-geometry           = 3d

epsilon_surface          = 0


; 7.3.14 Temperature Coupling

tcoupl                  = v-rescale                    ; temperature
coupling with

tc_grps                 = system        ; groups to couple seperately to
temperature bath

tau_t                   = 0.1                        ; [ps] time constant
for coupling

ref_t                   = 303                        ; [K] reference
temperature for coupling


; 7.3.17 Velocity Generation

gen_vel                 = no           ; generate velocities according to
Maxwell distribution of temperature

gen_temp                = 303           ; [K] temperature for Maxwell
distribution

gen_seed                = -1            ; [integer] used to initialize
random generator for random velocities


; 7.3.18 Bonds

constraints             = h-bonds

constraint_algorithm    = LINCS         ; LINear Constraint Solver

continuation            = no            ; no = apply constraints to the
start configuration

lincs_order             = 4             ; highest order in the expansion of
the contraint coupling matrix

lincs_iter              = 1             ; number of iterations to correct
for rotational lengthening

lincs_warnangle         = 30            ; [degrees] maximum angle that a
bond can rotate before LINCS will complain

Best regards,

Mahsa