2 residues force field

Automatic force field optimization. Given a force field and a set of reference data (e.g. energies and forces from QM calculations, experimental measurements), we tune the force field parameters such that it accurately reproduces the reference data.
POST REPLY
User avatar
Noelia Beatriz Luque
Posts: 16
Joined: Mon Nov 19, 2012 11:12 am

2 residues force field

Post by Noelia Beatriz Luque » Wed Jan 16, 2013 12:28 pm

Dear Lee-Ping,

I am still having some problems in getting good numbers for the parameters, so I have decided to start with a simpler system: Ti4+ as a residue with one atom, and, for the moment, a first shell of 5 or 6 TIP3P water molecules. I decided for this sytem, because I already have 200 configurations computed at DFT level and I just need to parametrize the two vdw parameters of Ti. So, I think it is a easy benchmark case.
I can run gromacs for a simple dynamics with this system, however when I use forcebalance, I get the following error:

FFMolecule = self.FF.FFMolecules[molname]
KeyError: 'TIT'


The .itp file I am using is the one below. I think the problem raises from using two residues in the same .itp file. But I am not sure.

[ defaults ]
1 2 yes 0.5 0.5

[ atomtypes ]
Ti 22 47.86700 4.00000 A 3.17500e-01 7.11756e-02 ; PARM 5 6
OW 8 16.00 -0.834 A 3.15061e-01 6.36386e-01
HW 1 1.008 0.417 A 0.00000e+00 0.00000e+00

[ bondtypes ]
OW HW 1 0.09572 502416.0

[ angletypes ]
HW OW HW 1 104.520 628.020

[ moleculetype ]
TIT 1
SOL 2

[ atoms ]
;id attype res nr resname atname cgnr charge
1 Ti 1 TIT Ti 1 4.00000
2 OW 2 SOL OW 2 -0.834
3 HW 2 SOL HW1 2 0.417
4 HW 2 SOL HW2 2 0.417

[ bonds ]
; i j funct length force_constant
2 3 1
2 4 1

[ angles ]
; i j k funct angle force_constant
3 2 4 1

The topology file is :
include "ti_6h2o.itp"

[ system ]
First intent of Ti in 6 H2O extracted from dft

[ molecules ]
TIT 1
SOL 6

Thanks in advance for your help!
Best regards,
Noelia

User avatar
Lee-Ping Wang
Posts: 102
Joined: Sun Jun 19, 2011 5:14 pm

Re: 2 residues force field

Post by Lee-Ping Wang » Wed Jan 16, 2013 2:36 pm

Hi Noelia,

GROMACS can be a bit confusing with residues and molecules. ForceBalance doesn't look at the residues - it labels atoms internally using the molecule name and the atom number within the molecule (this only affects you if you are parameterizing the charges, since they are properties of the atoms). This is because you can have two of the same residue in a molecule, and that can lead to duplicate labels and other problems.

I think your provided .itp file has a mistake in that you can only have one entry in the [ moleculetype ] section. Did you make sure that this system runs in GROMACS?

There are two ways to set up the GROMACS files: You can either define one big molecule with 19 atoms in the [ atoms ] section (then you may name the residues however you want), or you define two separate molecules in two .itp files. In the latter case, you need to #include both of them in the .top file and specify them in the ForceBalance input file under the "forcefield" option. (I think you can have two [ moleculetype ] sections within one .itp file, but I'm not sure if ForceBalance will parse this correctly.)

The one molecule setup looks like this:
"tic.itp"

Code: Select all

[ defaults ]
1 2 yes 0.5 0.5

[ atomtypes ]
Ti 22 47.86700 4.00000 A 3.17500e-01 7.11756e-02 ; PARM 5 6
OW 8 16.00 -0.834 A 3.15061e-01 6.36386e-01
HW 1 1.008 0.417 A 0.00000e+00 0.00000e+00

[ bondtypes ]
OW HW 1 0.09572 502416.0

[ angletypes ]
HW OW HW 1 104.520 628.020

[ moleculetype ]
TIC 1

[ atoms ]
;id attype res nr resname atname cgnr charge
1 Ti 1 TIT Ti 1 4.00000
2 OW 2 SOL OW 2 -0.834
3 HW 2 SOL HW1 2 0.417
4 HW 2 SOL HW2 2 0.417
5 OW 3 SOL OW 3 -0.834
6 HW 3 SOL HW1 3 0.417
7 HW 3 SOL HW2 3 0.417
(repeat this until you have six molecules)

[ bonds ]
; i j funct length force_constant
2 3 1
2 4 1
5 6 1
5 7 1
(repeat this until you have twelve bonds)

[ angles ]
; i j k funct angle force_constant
3 2 4 1 
6 5 7 1 
(repeat this until you have six angles)
"topol.top"

Code: Select all

include "tic.itp"

[ system ]
First intent of Ti in 6 H2O extracted from dft

[ molecules ]
TIC 1
The two-molecule setup looks something like this:

"tit.itp"

Code: Select all

[ defaults ]
1 2 yes 0.5 0.5

[ atomtypes ]
Ti 22 47.86700 4.00000 A 3.17500e-01 7.11756e-02 ; PARM 5 6
OW 8 16.00 -0.834 A 3.15061e-01 6.36386e-01
HW 1 1.008 0.417 A 0.00000e+00 0.00000e+00

[ bondtypes ]
OW HW 1 0.09572 502416.0

[ angletypes ]
HW OW HW 1 104.520 628.020

[ moleculetype ]
TIT 1

[ atoms ]
;id attype res nr resname atname cgnr charge
1 Ti 1 TIT Ti 1 4.00000
"sol.itp"

Code: Select all

[ moleculetype ]
SOL 2

[ atoms ]
;id attype res nr resname atname cgnr charge
1 OW 1 SOL OW 1 -0.834
2 HW 1 SOL HW1 1 0.417
3 HW 1 SOL HW2 1 0.417

[ bonds ]
; i j funct length force_constant
1 2 1
1 3 1

[ angles ]
; i j k funct angle force_constant
2 1 3 1 
"topol.top"

Code: Select all

#include "tit.itp"
#include "sol.itp"

[ system ]
First intent of Ti in 6 H2O extracted from dft

[ molecules ]
TIT 1
SOL 6
I hope this helps.

- Lee-Ping

POST REPLY