[gmx-users] Re: Chloroform (CHCl3) solvent box for G53a5 force field

[Pablo Englebienne]

OK, I started over with the CHCl3 box from scratch. I prepared the following itp file from the CHCL3 parameters in ffG53a5.rtp:

---[chcl3.itp]---
[ moleculetype ]
; Name            nrexcl
CHCL3             3

[ atoms ]
;   nr       type  resnr residue  atom   cgnr     charge       mass
     1       HCHL      1  CHCL3   HChL      1      0.082      1.008   
     2       CCHL      1  CHCL3   CChL      1      0.179     12.011  
     3      CLCHL      1  CHCL3  CLCh1      1     -0.087     35.453 
     4      CLCHL      1  CHCL3  CLCh2      1     -0.087     35.453 
     5      CLCHL      1  CHCL3  CLCh3      1     -0.087     35.453

[ bonds ]
;  ai    aj funct            
    1     2     2    gb_39
    2     3     2    gb_40
    2     4     2    gb_40
    2     5     2    gb_40
    1     3     2    gb_47
    1     4     2    gb_47
    1     5     2    gb_47
    3     4     2    gb_48
    3     5     2    gb_48
    4     5     2    gb_48

[ angles ]
;  ai    aj    ak funct          
    1     2     3     2  ga_43
    1     2     4     2  ga_43
    1     2     5     2  ga_43
    3     2     4     2  ga_44
    3     2     5     2  ga_44
    4     2     5     2  ga_44
---[chcl3.itp]---

I noticed that G53a5 includes 4 types of bond stretching terms specific for CHCl3 (C-Cl, C-H, H-Cl and Cl-Cl), therefore I specified all of them. Should all of these terms be harmonic bonds (function type 2) or some (e.g., the H-Cl and Cl-Cl terms) should be type 6 (as described in section 5.4 of the manual)? I tried with both types and I get the same minimized structure with the following topology and mdp files:

---[topol.top]---
; Include forcefield parameters
#include "ffG53a5.itp"
#include "chcl3.itp"

[ system ]
; Name
Chloroform

[ molecules ]
; Compound        #mols
CHCL3             1
---[topol.top]---

---[em.mdp]---
integrator          =  cg
nsteps              =  50000
;
;	Energy minimizing stuff
;
emtol               =  1e-5
emstep              =  0.01

;
;	Electrostatics
;
coulombtype         =  cut-off
rcoulomb            =  0
ns_type             =  simple
nstlist             =  0
rlist		    =  0

;
;	vdW
;
rvdw                =  0

;
;       PBC
;
pbc		    =  no
---[em.mdp]---

In the minimized structure, not all C-Cl distances are equivalent, although the minimization converges OK.

Any comments?