In [3]:
# OpenMM Imports
import simtk.openmm as mm
import simtk.openmm.app as app

# ParmEd Imports
from parmed import load_file
from parmed.openmm.reporters import NetCDFReporter
from parmed import unit as u
from parmed import gromacs

In [5]:
gromacs.GROMACS_TOPDIR = "/Users/jan-hendrikprinz/Studium/git/adaptive-sampling/package/examples/ntl9/top"

In [6]:
# Load the Gromacs files
print('Loading Gromacs files...')
top = load_file('files/topol-NTL9.top')
gro = load_file('files/start-NTL9.gro')


Loading Gromacs files...

In [8]:
gro.write_pdb('files/initial.pdb')

In [ ]:
# Create the OpenMM system
print('Creating OpenMM System')
system = top.createSystem(nonbondedMethod=app.NoCutoff,
                          constraints=app.HBonds, implicitSolvent=app.GBn2,
                          implicitSolventSaltConc=0.1*u.moles/u.liter,
)

# Create the integrator to do Langevin dynamics
integrator = mm.LangevinIntegrator(
                        300*u.kelvin,       # Temperature of heat bath
                        1.0/u.picoseconds,  # Friction coefficient
                        2.0*u.femtoseconds, # Time step
)

# Define the platform to use; CUDA, OpenCL, CPU, or Reference. Or do not specify
# the platform to use the default (fastest) platform
platform = mm.Platform.getPlatformByName('CUDA')
prop = dict(CudaPrecision='mixed') # Use mixed single/double precision

# Create the Simulation object
sim = app.Simulation(top.topology, system, integrator, platform, prop)

# Set the particle positions
sim.context.setPositions(gro.positions)

# Minimize the energy
print('Minimizing energy')
sim.minimizeEnergy(maxIterations=500)

# Set up the reporters to report energies and coordinates every 100 steps
sim.reporters.append(
        app.StateDataReporter(sys.stdout, 100, step=True, potentialEnergy=True,
                              kineticEnergy=True, temperature=True)
)
sim.reporters.append(
        NetCDFReporter('dhfr_gb.nc', 100, crds=True)
)

# Run dynamics
print('Running dynamics')
sim.step(10000)

In [3]:
traj = md.formats.GroTrajectoryFile('files/start-NTL9.gro')

In [5]:
traj.topology


Out[5]:
<mdtraj.Topology with 1 chains, 4540 residues, 14100 atoms, 0 bonds at 0x10d5350d0>

In [ ]: