[gmx-users] BAR calculations with tabulated non-bonded potentials

[Rustad, James R]

Dear GMX users,

I cannot figure out from the manual whether it is possible to use tabulated non-bonded potentials in a BAR calculation.

I want to estimate the free energy for changing Fe2+ to Fe3+.  These have different charges, of course, and also different non-bonded interactions.  These non-bonded interactions are neither LJ nor Buckingham, so I am using a lookup table (table_O_Fe2.xvg and table_O_Fe3.xvg) for the "vdw" part (vdw_type=User).

I think I should be able to define Fe2 and Fe3 atoms types, and associate the tables (table_O_Fe2.xvg and table_O_Fe3.xvg) with those types in the mdp file:

energygrps               = O Si Na Ca Al Fe2 Fe3
energygrp_table           = O Al  O Ca  O Fe2  O Fe3  O Na  O O  O Si  Si Al  Si Ca  Si Fe2  Si Fe3  Si Na  Si Si  Na Al  Na Ca  Na Fe2  Na Fe3  Na Na  Ca Al  Ca Ca  Ca Fe2  Ca Fe3  Al Al  Al Fe2  Al Fe3  Fe2 Fe2  Fe2 Fe3  Fe3 Fe3

and in the .top file:

[ atomtypes ]
;name       mass          charge    ptype   A            B
 O      15.99940     -1.20000   A       0.0          0.0
 Na     22.99000      0.60000   A       0.0          0.0
 Si     28.08000      2.40000  A        0.0          0.0
 Ca     40.07800      1.20000   A       0.0          0.0
 Al     26.98200      1.80000   A       0.0          0.0
 Fe2    55.30500      1.20000   A       0.0          0.0
 Fe3    55.30500      1.80000   A       0.0          0.0

[ nonbond_params ]
; i     j       func   C6                      C12
 O    O     1     4.0904959800            2.12268E-09
 O    Si    1     32.858468400            9.64853E-11
 O    Na    1     2.2541869930            4.82427E-10
 O    Ca    1     2.9149186000            4.82427E-10
 O    Al    1     34.887265640            8.68368E-11
 O    Fe3   1     40.425524240            1.92800E-10
 O    Fe2   1      7.542                  1.92800E-10

...

Also in the .top file, I define a moleculetype Fe, which is Fe2 in state A and Fe3 in state B:



[ moleculetype ]
; name nrexcl
Fe     0
[ atoms ]
; nr   type   resnr  residue   atom   cgnr   charge  mass     typeB chargeB massB
1       Fe2    1     Fe        Fe      2       1.2   55.3050   Fe3    1.8    55.305


I expect that I could specify this atom in the gro file like this:

...
 1764Na      Na 1764   3.448   3.335   0.801  0.9202  0.1009 -0.6229
 1765Fe      Fe 1765   1.123   1.300   0.562  0.3254 -1.0409 -0.6704
 1766Ca      Ca 1766   2.536   3.320   1.379 -0.5737 -0.4492 -1.8149
...

The problem is that grompp fails, telling me
"
Fatal error:
Group Fe2 not found in indexfile
"

So, while I want to associate the tables with atom types, it seems GMX is expecting that I am associating the tables (through the energy groups) with a particular atom, which, in my case, moves between two atom types in the BAR calculation.

Can anyone tell me whether I am just suffering from a misconception about how the tables (or energy groups) are used, or whether it is in fact just not possible to do BAR (or other free energy) calculations in this way?

Thanks

Jim Rustad




-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20120518/5b6b3264/attachment.html>