[gmx-developers] shell polarizable water model

[Ananya Mondal]

Hi Mikhail,
I am using GMX4.5.5 to simulate SWM4-NDP polarizable water model, it
is working fine.
When i am running separately NA+ polarizable ion in SWM4-NDP  water
and Cl- polarizable ion in SWM4-NDP polarizable water model, it is
working fine.
But when I use the both NA+ and CL- (NaCl salt solution) polarizable
model which is in conjuction with SWM4-NDP model, the RMS force is
very large:

can anyone give some idea, what i am missing..
Thanking you
Ananya

here is my itp file
[ defaults ]
LJ      Geometric
[ atomtypes ]
;name    mass         charge ptype    c6        c12
  WO    15.99940       0.0     A     0.0        0.0
  WH     1.00800       0.0     A     0.0        0.0
  WS     0.0           0.0     S     0.0        0.0
  WD     0.0           0.0     D     0.0        0.0
  NAc   22.989769      0.0     A     0.0        0.0
  NAs    0.0000000     0.0     S     0.0        0.0
  CLc   35.45300       0.0     A     0.0        0.0
  CLs    0.0000000     0.0     S     0.0        0.0
[ nonbond_params ]
WO      WO      1      3.67796770E-03   3.83182123E-06
WO      NAc     1      7.84740725E-04   4.51332937E-07
WO      CLc     1      6.24036649E-03   1.83075990E-05
NAc     NAc     1      1.56840324E-04   4.66461074E-08
CLc     CLc     1      9.08923242E-03   6.44595639E-05
NAc     CLc     1      1.46880350E-03   2.62418143E-06
[ moleculetype ]
; molname       nrexcl
SM2             2
[ atoms ]
; id    at type res nr  residu name     at name         cg nr   charge
1       WO      1       SM2             OW1             1      1.71636
2       WH      1       SM2             HW2             1      0.55733
3       WH      1       SM2             HW3             1      0.55733
4       WD      1       SM2             DW              1     -1.11466
5       WS      1       SM2             SW              1     -1.71636
[ polarization ]
1       5       1       0.00097822
[ constraints ]
; i     funct   doh     dhh
1       2       1       0.09572
1       3       1       0.09572
3       2       1       0.15139
[ virtual_sites3 ]
4       1       2       3       2       0.5      0.024034
[ exclusions ]
; iatom excluded from interaction with i
1       2       3       4       5
2       1       3       4       5
3       1       2       4       5
4       1       2       3       5
5       1       2       3       4
[ moleculetype ]
; molname       nrexcl
NA+            1
[ atoms ]
; id    at type res nr  residu name     at name cg nr   charge
1       NAc       1       NA             NAc    1        1.687597
2       NAs       1       NA             NAs    1       -0.687597
[ polarization ]
; See notes above.      alpha (nm^3)
1       2       1       0.000157
[ exclusions ]
; iatom excluded from interaction with i
1       2
2       1

[ moleculetype ]
; molname       nrexcl
Cl-             1

[ atoms ]
; id    at type res nr  residu name     at name cg nr   charge
1       CLc     1       ion             CLc     1        2.457187
2       CLs     1       ion             CLs     1       -3.457187

[ polarization ]
; See notes above.      alpha (nm^3)
1       2       1       0.003969

[ exclusions ]
; iatom excluded from interaction with i
1       2
2       1

here is my g.mdb file

include                  =
define                   =
integrator               = md
tinit                    = 0
dt                       = 0.0005
nsteps                   = 1000000
init_step                = 0
simulation_part          = 1
comm-mode                = Linear
nstcomm                  = 10
comm-grps                = system
bd-fric                  = 0
ld-seed                  = 1993
emtol                    = 0.1
emstep                   = 0.01
niter                    = 30
fcstep                   = 0
nstcgsteep               = 1000
nbfgscorr                = 10
rtpi                     = 0.05
nstxout                  = 10
nstvout                  = 00
nstfout                  = 0
nstlog                   = 50
nstcalcenergy            = -1
nstenergy                = 10
nstxtcout                = 250
xtc-precision            = 1000
xtc-grps                 =
energygrps               = System
nstlist                  = 5
ns-type                  = grid
pbc                      = xyz
rlist                    = 0.75
coulombtype              = PME
rcoulomb-switch          = 0
rcoulomb                 = 0.75
epsilon_r                = 1
epsilon_rf               = 1
vdw-type                 = Cut-off
rvdw-switch              = 0
rvdw                     = 0.75
DispCorr                 = EnerPres
table-extension          = 1
energygrp_table          =
fourierspacing           = 0.12
fourier_nx               = 32
fourier_ny               = 32
fourier_nz               = 32
pme_order                = 4
ewald_rtol               = 1e-05
ewald_geometry           = 3d
epsilon_surface          = 0
optimize_fft             = no
gb_algorithm             = Still
nstgbradii               = 0
rgbradii                 = 0
gb_saltconc              = 0
implicit_solvent         = No
tcoupl                   = nose-hoover
nsttcouple               = -1
nh-chain-length          = 10
tc-grps                  = system
tau-t                    = 0.05
ref-t                    = 300
gen-vel                  = yes
gen-temp                 = 300
gen-seed                 = 173529
constraints              = all-bonds
constraint-algorithm     = shake
continuation             = no
Shake-SOR                = no
shake-tol                = 0.0001
lincs-order              = 4
lincs-iter               = 2
lincs-warnangle          = 30
morse                    = no
energygrp_excl           =
nwall                    = 0
wall_type                = 9-3
wall_r_linpot            = -1
wall_atomtype            =
wall_density             =
wall_ewald_zfac          = 3
pull                     = no
disre                    = No
disre-weighting          = Conservative
disre-mixed              = no
disre-fc                 = 1000
disre-tau                = 0
nstdisreout              = 100
orire                    = no
orire-fc                 = 0
orire-tau                = 0
orire-fitgrp             =
nstorireout              = 100
dihre                    = no
dihre-fc                 = 1000
free-energy              = no
init-lambda              = 0
delta-lambda             = 0
foreign_lambda           =
sc-alpha                 = 0
sc-power                 = 0
sc-sigma                 = 0.3
nstdhdl                  = 10
separate-dhdl-file       = yes
dhdl-derivatives         = yes
dh_hist_size             = 0
dh_hist_spacing          = 0.1
couple-moltype           =
couple-lambda0           = vdw-q
couple-lambda1           = vdw-q
couple-intramol          = no
userreal4                = 0