Using TOPO from Python

Most users drive TOPO entirely through md.ini and the topo-mdrun command (Simulation control options). But every step is also a documented Python function, so you can build a model, inspect its forces, script a custom protocol, or post- process trajectories from your own code. This page is a tour of the public API.

The most useful names are re-exported at the top level of the topo package:

import topo
# model building
topo.models.buildCoarseGrainModel
topo.system
# config + reporter
topo.read_simulation_config, topo.SimulationConfig, topo.topoReporter
# multi-copy + trajectory splitting
topo.make_noninteracting_copies, topo.split_chains
# native-contact (Q) analysis
topo.build_native_contacts, topo.fraction_native_contacts, topo.load_domains

Build a model from a structure

topo.models.buildCoarseGrainModel() is the single entry point. It reads a structure, keeps the Cα atoms, and assembles the full force field (bonds, angles, torsions, Yukawa electrostatics, and the structure-based contacts). It returns a topo.core.system.system object holding the OpenMM System, Topology, and positions.

import topo

cg = topo.models.buildCoarseGrainModel(
    "P0CX28_clean.pdb",
    domain_def="domain.yaml",       # optional: per-domain contact scaling
    stride_output_file=None,        # None -> run STRIDE automatically, cache it
    constraints="AllBonds",         # rigid bonds (default); None -> flexible
    box_dimension=None,             # None -> no PBC; float or [x,y,z] -> PBC
    minimize=False,                 # native PDB is already the minimum
)

print(cg.n_chains, "chains,", cg.n_atoms, "CA beads")
print(type(cg.system))              # openmm.System with all forces added

Key arguments (full reference in topo.core.models.models):

Argument

Default

Meaning

structure_file

required

Input PDB/CIF; defines topology, force field, and native contacts.

domain_def

None

domain.yaml for per-domain/interface sidechain-contact scaling (Domain definition file (domain.yaml)).

stride_output_file

None

Precomputed STRIDE file; None runs STRIDE and caches <prefix>_stride.dat.

constraints

'AllBonds'

'AllBonds' = rigid (constraints, 15 fs step); None = flexible harmonic bonds.

box_dimension

None

None = no PBC; a float = cubic box (nm); [x,y,z] = rectangular box.

minimize

False

Energy-minimize the input geometry if large forces are found.

check_forces

True

Run the build-time energy/large-force check (set False on restart).

Useful attributes and methods of the returned object:

  • cg.system / cg.topology / cg.positions — the OpenMM objects to feed into a Simulation.

  • cg.n_atoms / cg.n_chains / cg.n_bonds / cg.n_angles / cg.n_torsions — geometry counts.

  • cg.rmin_matrix / cg.energy_matrix — the contact \(R_{ij}\) and \(\varepsilon_{ij}\) matrices (nm, kJ/mol).

  • cg.forceGroups — ordered map of force name → force object (one per log column).

  • cg.dumpTopology("model.psf") — write the CA topology as PSF.

  • cg.dumpStructure("model.pdb") — write the current CA coordinates as PDB.

  • cg.dumpForceFieldData("ff.txt") — dump all force-field parameters to text.

  • cg.reportEnergy(simulation, header=...) — print the total and per-force-group energy of a simulation’s current state.

Run dynamics yourself

If you want full control over the integrator and protocol, drive OpenMM directly. This reproduces the core of what topo.mdrun does:

import openmm as mm
from openmm import unit, app
import topo

cg = topo.models.buildCoarseGrainModel("P0CX28_clean.pdb", domain_def="domain.yaml")

integrator = mm.LangevinIntegrator(300*unit.kelvin, 0.01/unit.picosecond,
                                   0.015*unit.picoseconds)
sim = app.Simulation(cg.topology, cg.system, integrator)
sim.context.setPositions(cg.positions)
sim.context.setVelocitiesToTemperature(300*unit.kelvin)

# TOPO's fixed-width log with one energy column per force group:
sim.reporters.append(topo.topoReporter("traj.log", 1000, sbmObject=cg,
                                       step=True, time=True,
                                       potentialEnergy=True, temperature=True))
sim.reporters.append(app.DCDReporter("traj.dcd", 1000))
sim.step(100000)

Passing sbmObject=cg to topoReporter adds one energy column per force group, so you can watch the contact energy, angle energy, etc., separately.

In practice you usually don’t need this — topo.read_simulation_config() plus the topo.engine / topo.mdrun helpers already wrap build → set-up → protocol → finalize. Reach for the manual route only for a non-standard protocol the config file can’t express.

Read a control file

topo.read_simulation_config() parses an md.ini into a SimulationConfig dataclass with OpenMM units already applied — handy for scripting parameter sweeps:

import topo
cfg = topo.read_simulation_config("md.ini")
print(cfg.md_steps, cfg.ref_t, cfg.dt)
kwargs = cfg.build_kwargs()          # ready to pass to buildCoarseGrainModel
cg = topo.models.buildCoarseGrainModel(cfg.pdb_file, **kwargs)

SimulationConfig also resolves output paths (cfg.output_path('.dcd')), the checkpoint path (cfg.checkpoint_path()), the OpenMM platform (cfg.make_platform()), and the annealing step counts (cfg.quench_steps(), cfg.total_steps()).

Multi-copy replication and splitting

Pack independent copies of a chain into one system, then split the combined trajectory afterwards (Tutorial 4 — Many copies in one run (better GPU utilization)):

import topo
from openmm import unit

cg = topo.models.buildCoarseGrainModel("P0CX28_clean.pdb", domain_def="domain.yaml")

system, topology, positions = topo.make_noninteracting_copies(
    cg.system, cg.topology, cg.positions,
    n_copies=10, shift=5.0*unit.nanometer)

# ... run, producing a combined traj.dcd ...

# split back into one DCD per copy (memory-bounded streaming):
topo.split_chains("traj.dcd",
                  [f"traj_{k}.dcd" for k in range(10)],
                  center=True)

The copies are guaranteed non-interacting: bonded terms are duplicated per copy and each CustomNonbondedForce is restricted to intra-copy interaction groups, so the total energy is exactly n_copies × the single-chain energy.

Score native contacts (Q)

The Native-contact analysis (the Q score) analysis is fully scriptable; see that page for a complete example using topo.build_native_contacts() and topo.fraction_native_contacts().

Parse a log

topo.reporter.topo_reporter.readOpenMMReporterFile() reads a TOPO .log into a {column: [values]} dict (see Output files and the run log):

from topo.reporter.topo_reporter import readOpenMMReporterFile
data = readOpenMMReporterFile("traj/traj.log")
pe = data["Potential Energy (kJ/mole)"]

See also