[gmx-users] Shell molecular dynamics

David spoel at xray.bmc.uu.se
Thu Mar 20 19:36:18 CET 2003


On Thu, 2003-03-20 at 15:43, Martin Mucha wrote:
> Hi all,
> 
> I would like to use polarizability model in my simulations, 
> but I do not know how. There are several paragraphs about the topic in the 
> manual, however
> I did not find any mention about how to actually tell gromacs that I am using 
> the shell model, e.g. How to set up topology, mdp and gro files. I would be 
> grateful for any idea, or sample input files of a working shell MD.

see attached itp file for our polarizable water model. This is not a
perfect model! See JPCB 105 (2001) 2618-2626

You can use the mkyaw program in contrib to transform an SPC water box
into a SW box.
 



> 
>   Martin Mucha
> 
> 
> -- 
> address: UFCH JH AV CR, m. 224
>              Dolejskova 3, 18223 Praha 8, Czech Republic
> tel: +420266053495, +420736175269 
> -----------------------------------------------------
> 
> 
> _______________________________________________
> gmx-users mailing list
> gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please don't post (un)subscribe requests to the list. Use the 
> www interface or send it to gmx-users-request at gromacs.org.
-- 
Groeten, David.
________________________________________________________________________
Dr. David van der Spoel, 	Dept. of Cell and Molecular Biology
Husargatan 3, Box 596,  	75124 Uppsala, Sweden
phone:	46 18 471 4205		fax: 46 18 511 755
spoel at xray.bmc.uu.se	spoel at gromacs.org   http://xray.bmc.uu.se/~spoel
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-------------- next part --------------
; 
; Topology file for SW
;
; Paul van Maaren and David van der Spoel
; Molecular Dynamics Simulations of Water with Novel Shell Model Potentials
; J. Phys. Chem. B. 105 (2618-2626), 2001
;
; Force constants for the shell are given by:
;
; k = qs^2/(4 pi eps0 alpha)
;
; Some data: mu (water) = 1.8546 D ( 0.0386116 e nm)
;            1/(4 pi eps0 alpha) = 94513.94
;
; Alpha-X = 1.415   kx = 608069
; Alpha-Y = 1.528   ky = 563101
; Alpha-Z = 1.468   kz = 586116
;
; Alpha   = 1.470   k  = 585318
;
; Bonding parameters from (but without cubic term):
; D. M. Ferguson: 
; Parametrization and Evaluation of a Flexible Water Model 
; J. Comp. Chem. 16(4), 501-511 (1995)
;
; Possible defines:
; -DANISOTROPIC Select anisotropic polarizibility.
; -DRIGID       Rigid model (flexible is default)
; -DPOSRES      Position restrain oxygen atoms
;

[ defaults ]
LJ	Arithmetic

[ 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

[ nonbond_params ]
#ifdef RIGID
#ifdef ANISOTROPIC
WH      WH      1       4.0e-5          4.0e-8
WS      WO      1       1.0e-6          1.0e-12
WS      WH      1       4.0e-5          2.766e-08
WO      WO      1       2.0e-3          1.174e-06
#else
WH      WH      1       4.0e-5          4.0e-8
WS      WO      1       1.0e-6          1.0e-12
WS      WH      1       4.0e-5          2.769e-08
WO      WO      1       2.0e-3          1.176e-06
#endif
#else
#ifdef ANISOTROPIC
WH      WH      1       4.0e-5          4.0e-8
WS      WO      1       1.0e-6          1.0e-12
WS      WH      1       4.0e-5          2.910e-08
WO      WO      1       2.0e-3          1.189e-06
#else
WH      WH      1       4.0e-5          4.0e-8
WS      WO      1       1.0e-6          1.0e-12
WS      WH      1       4.0e-5          2.937e-08
WO      WO      1       2.0e-3          1.187e-06
#endif
#endif

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; This is a the 'classical YAW' model, in which we do have the dummy.
;; The shell is attached to the dummy, in this case the gas-phase
;; quadrupole is correct. Water_pol routine can be used for this
;; model. This has four interaction sites.
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
[ moleculetype ]
; molname	nrexcl
SW		2

[ atoms ]
; id	at type	res nr 	residu name	at name		cg nr	charge
1	WO	1	SM2		OW1		1	1.24588
2	WH	1	SM2		HW2		1	0.62134
3	WH	1	SM2		HW3		1	0.62134
4	WD	1	SM2		DW		1 	0.0
5	WS	1	SM2		SW		1	-2.48856

[ position_restraints ]
; See notes above.
; i	funct	kx	ky	kz	rOH	rHH	rOD
#ifdef ANISOTROPIC
1	2 	608069	563101	586116	0.09572	0.15139 0.0137408
#else
1	2 	585318 	585318 	585318 	0.09572 0.15139 0.0137408
#endif

#ifdef RIGID
[ settles ]
; i	funct	dOH	dHH
1	1	0.09572	0.15139

#else

[ bonds ]
1	2	1   0.09572     458148.
1	3	1   0.09572     458148.

[ angles ]
; i	j	k
2	1	3    1   104.52     417.6
#endif

[ dummies3 ]
; The position of the dummies is computed as follows:
;
;		O
;  	      
;	    	D
;	  
;	H		H
;
; 2 * b = distance (OD) / [ cos (angle(DOH)) 	* distance (OH) ]
;	  0.0137408 nm	/ [ cos (104.52 / 2 deg) * 0.09572 nm	]
;         0.01557 nm 
; Dummy pos x4 = x1 + a*(x2-x1) + b*(x3-X1)
;
; Dummy from			funct	a		b
4	1	2	3	1	0.117265878	0.117265878

[ 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


#ifdef POSRES
; Restrain the oxygen...
[ position_restraints ]
; iatom type	fx	fy	fz
1	1	100	100	100
#endif



More information about the gromacs.org_gmx-users mailing list