Commit 79e76a4e authored by ChayaSt's avatar ChayaSt
Browse files

added a test for nbfix generator

parent 8a95c7c7
...@@ -11,6 +11,7 @@ try: ...@@ -11,6 +11,7 @@ try:
except ImportError: except ImportError:
from io import StringIO from io import StringIO
import os import os
import warnings
class TestForceField(unittest.TestCase): class TestForceField(unittest.TestCase):
"""Test the ForceField.createSystem() method.""" """Test the ForceField.createSystem() method."""
...@@ -535,6 +536,66 @@ class AmoebaTestForceField(unittest.TestCase): ...@@ -535,6 +536,66 @@ class AmoebaTestForceField(unittest.TestCase):
diff = norm(f1-f2) diff = norm(f1-f2)
self.assertTrue(diff < 0.1 or diff/norm(f1) < 1e-3) self.assertTrue(diff < 0.1 or diff/norm(f1) < 1e-3)
def test_nbfix_generator(self):
""" Test the nbfix generator"""
warnings.filterwarnings('ignore', category=CharmmPSFWarning)
psf = CharmmPsfFile('systems/ions.psf')
pdb = PDBFile('systems/ions.pdb')
params = CharmmParameterSet('systems/toppar_water_ions.str')
# Box dimensions (found from bounding box)
psf.setBox(12.009*angstroms, 12.338*angstroms, 11.510*angstroms)
# Turn off charges so we only test the Lennard-Jones energies
for a in psf.atom_list:
a.charge = 0.0
# Now compute the full energy
plat = Platform.getPlatformByName('Reference')
system = psf.createSystem(params, nonbondedMethod=PME,
nonbondedCutoff=5*angstroms)
con = Context(system, VerletIntegrator(2*femtoseconds), plat)
con.setPositions(pdb.positions)
# Now set up stystem from ffxml. Setting chareges to 0 so we only test the Lennard-Jones energies
xml = """
<ForceField>
<AtomTypes>
<Type name="SOD" class="SOD" element="Na" mass="22.98977"/>
<Type name="CLA" class="CLA" element="Cl" mass="35.45"/>
</AtomTypes>
<Residues>
<Residue name="CLA">
<Atom name="CLA" type="CLA" charge="0"/>
</Residue>
<Residue name="SOD">
<Atom name="SOD" type="SOD" charge="0"/>
</Residue>
</Residues>
<NonbondedForce coulomb14scale="1.0" lj14scale="1.0">
<UseAttributeFromResidue name="charge"/>
<Atom type="SOD" sigma="0.251367073323" epsilon="0.1962296"/>
<Atom type="CLA" sigma="0.404468018036" epsilon="0.6276"/>
</NonbondedForce>
<NBFixForce>
<NBFix type1="CLA" type2="SOD" emin="0.350933" rmin="0.3731"/>
</NBFixForce>
</ForceField> """
ff = ForceField(StringIO(xml))
system2 = ff.createSystem(pdb.topology, nonbondedMethod=PME,
nonbondedCutoff=5*angstroms)
con2 = Context(system2, VerletIntegrator(2*femtoseconds), plat)
con2.setPositions(pdb.positions)
state = con.getState(getEnergy=True, enforcePeriodicBox=True)
ene = state.getPotentialEnergy().value_in_unit(kilocalories_per_mole)
state2 = con2.getState(getEnergy=True, enforcePeriodicBox=True)
ene2 = state2.getPotentialEnergy().value_in_unit(kilocalories_per_mole)
self.assertAlmostEqual(ene, ene2)
def test_Wildcard(self): def test_Wildcard(self):
"""Test that PeriodicTorsionForces using wildcard ('') for atom types / classes in the ffxml are correctly registered""" """Test that PeriodicTorsionForces using wildcard ('') for atom types / classes in the ffxml are correctly registered"""
......
CRYST1 12.009 12.338 11.510 90.00 90.00 90.00 P 1
ATOM 1 SOD SOD I 1 -5.844 1.432 3.239 1.00 0.00 ION NA
ATOM 2 CLA CLA I 2 -5.605 -3.153 1.213 1.00 0.00 ION CL
END
\ No newline at end of file
PSF
4 !NTITLE
REMARKS original generated structure x-plor psf file
REMARKS topology ../../../../../Charmm36_FF/toppar/top_all36_cgenff.rtf
REMARKS topology ../../../../../Charmm36_FF/toppar/toppar_water_ions.str
REMARKS segment ION { first NONE; last NONE; auto angles dihedrals }
2 !NATOM
1 ION 1 SOD SOD SOD 1.000000 22.9898 0
2 ION 2 CLA CLA CLA -1.000000 35.4500 0
0 !NBOND: bonds
0 !NTHETA: angles
0 !NPHI: dihedrals
0 !NIMPHI: impropers
0 !NDON: donors
0 !NACC: acceptors
0 !NNB
0 0
1 0 !NGRP
0 0 0
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment