Hi,
I was originally having this problem as in #3415 and I think I've found that it's an issue similar to #3815.
I can successfully minimize and start a production simulation for the following system (solvated but waters aren't shown because the file is large):
MODEL 1
ATOM 1 N GLY 1 11.370 20.100 6.000 1.00 0.00 N
ATOM 2 H1 GLY 1 11.610 20.830 5.300 1.00 0.00 H
ATOM 3 H2 GLY 1 12.190 19.460 5.920 1.00 0.00 H
ATOM 4 H3 GLY 1 10.480 19.610 5.780 1.00 0.00 H
ATOM 5 CA GLY 1 11.290 20.770 7.360 1.00 0.00 C
ATOM 6 HA1 GLY 1 12.240 21.260 7.530 1.00 0.00 H
ATOM 7 HA2 GLY 1 10.480 21.480 7.250 1.00 0.00 H
ATOM 8 C GLY 1 11.000 19.870 8.550 1.00 0.00 C
ATOM 9 O GLY 1 10.840 18.700 8.270 1.00 0.00 O
ATOM 10 N GLY 2 11.010 20.400 9.810 1.00 0.00 N
ATOM 11 HN GLY 2 11.260 21.360 9.920 1.00 0.00 H
ATOM 12 CA GLY 2 10.670 19.710 11.040 1.00 0.00 C
ATOM 13 HA1 GLY 2 11.570 19.210 11.350 1.00 0.00 H
ATOM 14 HA2 GLY 2 10.430 20.520 11.710 1.00 0.00 H
ATOM 15 C GLY 2 9.490 18.840 11.090 1.00 0.00 C
ATOM 16 O GLY 2 8.450 19.230 10.640 1.00 0.00 O
ATOM 17 N GLY 3 9.650 17.670 11.660 1.00 0.00 N
ATOM 18 HN GLY 3 10.620 17.520 11.810 1.00 0.00 H
ATOM 19 CA GLY 3 8.750 16.580 11.890 1.00 0.00 C
ATOM 20 HA1 GLY 3 8.500 16.020 11.000 1.00 0.00 H
ATOM 21 HA2 GLY 3 9.140 16.020 12.730 1.00 0.00 H
ATOM 22 C GLY 3 7.380 16.940 12.400 1.00 0.00 C
ATOM 23 OT1 GLY 3 6.390 16.770 11.770 1.00 0.00 O
ATOM 24 OT2 GLY 3 7.360 17.440 13.660 1.00 0.00 O
ATOM 25 HT GLY 3 6.440 17.680 13.780 1.00 0.00 H
TER
ENDMDL
Changing the final line atom type from HXT to HT got rid of the template error I described at the bottom of #3415. However, early in a production simulation I get this stack trace:
Traceback (most recent call last):
File "/ifs/groups/urbancMRIGrp/andrews/GGG_CHARMM_Drude_2019_SWM4-NDP_water/openmm_setup.py", line 74, in <module>
simulation.step(600000000)
File "/ifs/groups/urbancMRIGrp/opt/miniconda/envs/openMM/lib/python3.9/site-packages/openmm/app/simulation.py", line 141, in step
self._simulate(endStep=self.currentStep+steps)
File "/ifs/groups/urbancMRIGrp/opt/miniconda/envs/openMM/lib/python3.9/site-packages/openmm/app/simulation.py", line 206, in _simul$
self.integrator.step(10) # Only take 10 steps at a time, to give Python more chances to respond to a control-c.
File "/ifs/groups/urbancMRIGrp/opt/miniconda/envs/openMM/lib/python3.9/site-packages/openmm/openmm.py", line 5009, in step
return _openmm.VerletIntegrator_step(self, steps)
openmm.OpenMMException: Particle coordinate is nan
I did some testing and finally tried a simulation with charged C-terminus by deleting the hydrogen HT from the above pdb. Everything worked fine after the removal of the COOH hydrogen. I went into the CHARMM drude xml file to look at the patches and I think I found something similar to #3815. The patch options for C-terminal glycine are:
<AllowPatch name="GNTE"/>
<AllowPatch name="ACE"/>
<AllowPatch name="CTEG"/>
<AllowPatch name="CNEG"/>
<AllowPatch name="CT1G"/>
<AllowPatch name="CT2G"/>
<AllowPatch name="CT3"/>
The only option of those that puts a hydrogen on one of the oxygens is here :
<Patch name="CNEG">
<ChangeAtom charge="2.0363207668819214" name="N" type="CNEG-N"/>
<ChangeAtom charge="-2.4183207668819215" name="DN" type="Drude-CNEG-N"/>
<ChangeAtom charge="0.272" name="HN" type="CNEG-HN"/>
<ChangeAtom charge="1.8202974717429852" name="CA" type="CNEG-CA"/>
<ChangeAtom charge="-1.7002974717429853" name="DCA" type="Drude-CNEG-CA"/>
<ChangeAtom charge="-0.017" name="HA1" type="CNEG-HA1"/>
<ChangeAtom charge="-0.017" name="HA2" type="CNEG-HA2"/>
<ChangeAtom charge="2.764526857311924" name="C" type="CNEG-C"/>
<ChangeAtom charge="-1.906526857311924" name="DC" type="Drude-CNEG-C"/>
<AddAtom charge="1.6663059805572338" name="OT1" type="CNEG-OT1"/>
<AddAtom charge="-1.6663059805572338" name="DOT1" type="Drude-CNEG-OT1"/>
<AddAtom charge="1.9633344060265052" name="OT2" type="CNEG-OT2"/>
<AddAtom charge="-1.9633344060265052" name="DOT2" type="Drude-CNEG-OT2"/>
<AddAtom charge="0.374" name="HT" type="CNEG-HT"/>
<AddAtom charge="-0.319" name="LPT1" type="CNEG-LPT1"/>
<AddAtom charge="-0.319" name="LPT2" type="CNEG-LPT2"/>
<AddAtom charge="-0.285" name="LPT3" type="CNEG-LPT3"/>
<AddAtom charge="-0.285" name="LPT4" type="CNEG-LPT4"/>
<RemoveAtom name="O"/>
<RemoveAtom name="DO"/>
<RemoveAtom name="LPOA"/>
<RemoveAtom name="LPOB"/>
<AddBond atomName1="C" atomName2="OT1"/>
<AddBond atomName1="C" atomName2="OT2"/>
<RemoveBond atomName1="C" atomName2="O"/>
<RemoveExternalBond atomName="C"/>
<VirtualSite atomName1="OT1" atomName2="C" atomName3="OT2" p1="-0.0005235721931185042" p2="-0.029995430854691738" p3="0" siteName="LPT1" type="localCoords" wo1="1" wo2="0" wo3="0" wx1="-1.0" wx2="0.0" wx3="1.0" wy1="0" wy2="-1" wy3="1"/>
<VirtualSite atomName1="OT1" atomName2="C" atomName3="OT2" p1="-0.0005235721931185042" p2="0.029995430854691738" p3="0" siteName="LPT2" type="localCoords" wo1="1" wo2="0" wo3="0" wx1="-1.0" wx2="0.0" wx3="1.0" wy1="0" wy2="-1" wy3="1"/>
<VirtualSite atomName1="OT2" atomName2="C" atomName3="HT" p1="-0.012485829976186883" p2="0.0005706440663482743" p3="0.03269217666591343" siteName="LPT3" type="localCoords" wo1="1" wo2="0" wo3="0" wx1="-1.0" wx2="0.0" wx3="1.0" wy1="0" wy2="-1" wy3="1"/>
<VirtualSite atomName1="OT2" atomName2="C" atomName3="HT" p1="-0.012485829976186883" p2="0.0005706440663482743" p3="-0.03269217666591343" siteName="LPT4" type="localCoords" wo1="1" wo2="0" wo3="0" wx1="-1.0" wx2="0.0" wx3="1.0" wy1="0" wy2="-1" wy3="1"/>
</Patch>
But there is no bond added for the hydrogen
<RemoveAtom name="O"/>
<RemoveAtom name="DO"/>
<RemoveAtom name="LPOA"/>
<RemoveAtom name="LPOB"/>
<AddBond atomName1="C" atomName2="OT1"/>
<AddBond atomName1="C" atomName2="OT2"/>
<RemoveBond atomName1="C" atomName2="O"/>
<RemoveExternalBond atomName="C"/>
So I guess the particle acts like an unconstrainted charge in the simulation? Either way, I saw that in #3815 they manually added the bond. So I can also do that. I just wanted to report this as also an issue within the drude force field and also, more importantly, ask if manually changing this entry would affect updates within conda when new versions release.
Hi,
I was originally having this problem as in #3415 and I think I've found that it's an issue similar to #3815.
I can successfully minimize and start a production simulation for the following system (solvated but waters aren't shown because the file is large):
Changing the final line atom type from
HXTtoHTgot rid of the template error I described at the bottom of #3415. However, early in a production simulation I get this stack trace:I did some testing and finally tried a simulation with charged C-terminus by deleting the hydrogen HT from the above pdb. Everything worked fine after the removal of the COOH hydrogen. I went into the CHARMM drude xml file to look at the patches and I think I found something similar to #3815. The patch options for C-terminal glycine are:
The only option of those that puts a hydrogen on one of the oxygens is here :
But there is no bond added for the hydrogen
So I guess the particle acts like an unconstrainted charge in the simulation? Either way, I saw that in #3815 they manually added the bond. So I can also do that. I just wanted to report this as also an issue within the drude force field and also, more importantly, ask if manually changing this entry would affect updates within conda when new versions release.