#!C:/Python26/python.exe

from simbody.simtk import *

# molecule-specialized simbody System
system = CompoundSystem()

matter = SimbodyMatterSubsystem(system)
decorations = DecorationSubsystem(system)

# molecular force field
forceField = DuMMForceFieldSubsystem(system)
forceField.loadAmber99Parameters()

protein = Protein("SIMTK")

protein.assignBiotypes()
system.adoptCompound(protein)

# finalize multibody system
system.modelCompounds() 

# show me a movie
system.defaultSubsystem().addEventReporter(
    VTKEventReporter(system, 0.020))

# Maintain a constant temperature
system.defaultSubsystem().addEventHandler(
    VelocityRescalingThermostat(system,  293.15, 0.1))

# Instantiate simbody model
system.realizeTopology()
state = system.defaultState()

# Relax the structure before dynamics run
LocalEnergyMinimizer.minimizeEnergy(system, state, 15.0)

# Simulate it.
integ = VerletIntegrator(system)
ts = TimeStepper(system, integ)
ts.initialize(state)
ts.stepTo(2.0)