GROMACS modification: No
I am attempting a simulation of ~ 200 toluene molecules using gmx. I used LigParGen server to create the .itp file. I used the resulting pdb file with packmol to produce a starting pdb file with 200 toluene molecules. When I attempt to minimize the structure, grompp complains that the system has a net charge and PME is not recommended.
The net charge of the system per molecule from the .itp file is -0.0001 units. However, when grompp processes the pdb file, the net charge is reported as -0.020000000 and PME is flagged as a warning. I wonder if there’s a way to work around this error. We could use coulombtype = cutoff as one alternative, however, I am curious to know if there’s a better fix? I am posting the charge section from the .itp file below.
; Name nrexcl
[ atoms ]
; nr type resnr residue atom cgnr charge mass
1 opls_800 1 TOL C00 1 -0.1325 12.0110
2 opls_801 1 TOL C01 1 -0.1435 12.0110
3 opls_802 1 TOL C02 1 -0.1326 12.0110
4 opls_803 1 TOL C03 1 -0.1354 12.0110
5 opls_804 1 TOL C04 1 -0.0896 12.0110
6 opls_805 1 TOL C05 1 -0.1356 12.0110
7 opls_806 1 TOL H06 1 0.1367 1.0080
8 opls_807 1 TOL H07 1 0.1359 1.0080
9 opls_808 1 TOL H08 1 0.1366 1.0080
10 opls_809 1 TOL H09 1 0.1372 1.0080
11 opls_810 1 TOL H0A 1 0.1368 1.0080
12 opls_811 1 TOL C0B 1 -0.1907 12.0110
13 opls_812 1 TOL H0C 1 0.0922 1.0080
14 opls_813 1 TOL H0D 1 0.0922 1.0080
15 opls_814 1 TOL H0E 1 0.0922 1.0080
-0.0001 is too small, u can just round it off to zero.
Not when the charges are actually reported to four decimal points. It is clear that this is an actual deficiency because the charges sum to -0.0001 when they should not. The charge in the ligand topology needs to be adjusted so it sums to 0. The accumulation of this error when there are 200 molecules of this type (-0.02) is definitely non-negligible.
Hi @jalemkul, can you embellish more on your comment “The charge in the ligand topology needs to be adjusted so it sums to 0”. There are going to be some non-English speaking people that could think there is a program that does/can do this, rather than maybe it is a directive to “do it yourself”. Then what exactly does “adjusted” mean as well. This language can also be misinterpreted. Can we, for instance if we have a total net charge of -0.0003, just remove -.0001 from 3 atoms?
Interestingly, I tested obabel on ribociclib.sdf downloaded from Pubchem to convert to TRIPOS mol2 file and vanilla obabel is even worse than ligpargen and the total charge was -0.0005. This is now seeming like a fundamental issue with small molecule files and charges.
Empirical adjustments to charges are fine, if appropriately validated. For a charge discrepancy that small, it’s clearly a floating-point rounding issue with the software that generated the charges. A lot of force fields don’t bother even going to the ten-thousandths place because that level of precision becomes rather meaningless in MD simulations. So, yes, distribute the charge difference sensibly, but check what you’re doing to make sure your topology is still valid.
Still seems like a software issue. I’ve never had a problem using Antechamber or CGenFF programs leading to nonsense charges like this.
Yes agreed. I created an issue for ligpargen on Github :). No sense in creating an issue for obabel, it won’t get done for a while.