[gmx-users] simulations with wall potential

[Beibei Wang]

Hi everyone,

I am trying to simulate a polymer in a water box with Martini force field, and with a constant thickness. To realize this, I am testing simulations with wall potential. My mdp file is as following:

integrator               = md
dt                       = 0.002
nsteps                   = 500000
nstcomm                  = 10

nstxout                  = 0
nstvout                  = 0
nstfout                  = 0
nstlog                   = 1000
nstenergy                = 5000

nstxtcout                = 5000
xtc_precision            = 100
xtc-grps                 =
energygrps               = Sol W

cutoff-scheme            = group
nstlist                  = 10
ns_type                  = grid
pbc                      = xy
periodic-molecules       = no
rlist                    = 1.4

coulombtype              = Shift
rcoulomb_switch          = 0.0
rcoulomb                 = 1.2
epsilon_r                = 15
vdw_type                 = Shift
rvdw_switch              = 0.9
rvdw                     = 1.2

;walls
nwall                     = 2
wall-atomtype             = P5 P5 ;P5 is a Martini particle type
wall-type                 = 10-4
wall-r-linpot             = 1 1
wall-density             = 1.2 1.2
wall-ewald-zfac           = 3

tcoupl                   = berendsen
tc-grps                  = Sytem ;W
tau_t                    = 1.0 ;1.0
ref_t                    = 310 ;310
pcoupl                   = no ; no pressure coupling in NVT

gen_vel                  = yes
gen_temp                 = 310
gen_seed                 = -1


constraints              = none
constraint_algorithm     = Lincs
unconstrained_start      = no
lincs_order              = 4
lincs_warnangle          = 30

The system box was 7.0 by 7.0 by 7.0 nm, was placed at zero point by editconf initially. Here is the head of my system:

    1GLF     B4    1   4.487   3.844   3.538
    1GLF     B5    2   4.360   4.133   3.770
    1GLF     B6    3   4.286   3.687   3.477
    2GLF     B1    4   4.409   4.396   3.647
    2GLF     B2    5   4.119   4.561   3.704
    2GLF     B3    6   4.533   4.499   3.435
    2GLF     B4    7   4.022   4.860   3.578
    2GLF     B5    8   3.961   4.977   3.209

But soon as the simulations started, some water molecules moved across the walls, and far away in Z direction. The main box was moved in Z direction as well:

    1GLF     B4    1   3.175   4.176-309.166  0.0008  0.0012 -0.3582
    1GLF     B5    2   3.082   4.558-309.239  0.0010 -0.0008 -0.3564
    1GLF     B6    3   3.188   3.981-309.342 -0.0042  0.0017 -0.3578
    2GLF     B1    4   2.936   4.683-309.448  0.0004 -0.0040 -0.3567
    2GLF     B2    5   2.602   4.708-309.427  0.0010  0.0007 -0.3586
    2GLF     B3    6   2.990   4.709-309.710  0.0004  0.0001 -0.3614
    2GLF     B4    7   2.375   4.788-309.664  0.0026 -0.0004 -0.3591
2GLF     B5    8   2.150   4.640-309.969 -0.0006 -0.0007 -0.3620

As I used NVT ensemble, the box size shown in the gro file was still 7.0 by 7.0 by 7.0 nm. But I am sure in this case how the wall potential was applied. Where is the two walls? Are they still at z=0 and box_z?

Could anyone help to check if there are anything wrong in my input files? I have tried to change wall-atomtype, wall-type, wall-r-linpot, and wall-density, with higher wall-density (such as 12), the box was separated into parts with a huge z range (see the following figure). The figure is hard to see, the light blue "particles" you could see are my system.

[cid:image001.jpg at 01D27491.8FEB2D00]

Thanks very much.