OpenMM
CharmmPsfFile Class Reference

A chemical structure instantiated from CHARMM files. More...

List of all members.

Public Member Functions

def __init__
 Opens and parses a PSF file, then instantiates a CharmmPsfFile instance from the data.
def loadParameters
 Loads parameters from a parameter set that was loaded via CHARMM RTF, PAR, and STR files.
def setBox
 Sets the periodic box boundary conditions.
def topology
 Create an OpenMM Topology object from the stored bonded network.
def createSystem
 Construct an OpenMM System representing the topology described by the prmtop file.
def system
 Return the cached system class -- it needs to be initialized via "createSystem" first!
def boxLengths
 Return tuple of 3 units.
def boxLengths
def boxVectors
 Return the box vectors.
def boxVectors
 Sets the box vectors.
def deleteCmap
 Deletes the CMAP terms from the CHARMM PSF.

Public Attributes

 residue_list
 atom_list
 bond_list
 angle_list
 dihedral_list
 dihedral_parameter_list
 improper_list
 cmap_list
 donor_list
 acceptor_list
 group_list
 title
 flags
 box_vectors
 urey_bradley_list

Static Public Attributes

int BOND_FORCE_GROUP = 0
int ANGLE_FORCE_GROUP = 1
int DIHEDRAL_FORCE_GROUP = 2
int UREY_BRADLEY_FORCE_GROUP = 3
int IMPROPER_FORCE_GROUP = 4
int CMAP_FORCE_GROUP = 5
int NONBONDED_FORCE_GROUP = 6
int GB_FORCE_GROUP = 6

Detailed Description

A chemical structure instantiated from CHARMM files.

This structure has numerous attributes that are lists of the elements of this structure, including atoms, bonds, torsions, etc. The attributes are

  • residue_list
  • atom_list
  • bond_list
  • angle_list
  • dihedral_list
  • dihedral_parameter_list
  • improper_list
  • cmap_list
  • donor_list # hbonds donors?
  • acceptor_list # hbond acceptors?
  • group_list # list of nonbonded interaction groups

Additional attribute is available if a CharmmParameterSet is loaded into this structure.

  • urey_bradley_list

The lengths of each of these lists gives the pointers (e.g., natom, nres, etc.)


Constructor & Destructor Documentation

def __init__ (   self,
  psf_name 
)

Opens and parses a PSF file, then instantiates a CharmmPsfFile instance from the data.

Parameters:
psf_name(str) Name of the PSF file (it must exist)

Member Function Documentation

def boxLengths (   self)

Return tuple of 3 units.

def boxLengths (   self,
  stuff 
)
def boxVectors (   self)

Return the box vectors.

def boxVectors (   self,
  stuff 
)

Sets the box vectors.

def createSystem (   self,
  params,
  nonbondedMethod = ff.NoCutoff,
  nonbondedCutoff = 1.0*u.nanometer,
  switchDistance = 0.0*u.nanometer,
  constraints = None,
  rigidWater = True,
  implicitSolvent = None,
  implicitSolventKappa = None,
  implicitSolventSaltConc = 0.0*u.moles/u.liter,
  temperature = 298.15*u.kelvin,
  soluteDielectric = 1.0,
  solventDielectric = 78.5,
  removeCMMotion = True,
  hydrogenMass = None,
  ewaldErrorTolerance = 0.0005,
  flexibleConstraints = True,
  verbose = False 
)

Construct an OpenMM System representing the topology described by the prmtop file.

You MUST have loaded a parameter set into this PSF before calling createSystem. If not, AttributeError will be raised. ValueError is raised for illegal input.

Parameters:
params(CharmmParameterSet) The parameter set to use to parametrize this molecule
nonbondedMethod(object=NoCutoff) The method to use for nonbonded interactions. Allowed values are NoCutoff, CutoffNonPeriodic, CutoffPeriodic, Ewald, or PME.
nonbondedCutoff(distance=1*nanometer) The cutoff distance to use for nonbonded interactions.
switchDistance(distance=0*nanometer) The distance at which the switching function is active for nonbonded interactions. If the switchDistance evaluates to boolean False (if it is 0), no switching function will be used. Illegal values will raise a ValueError
constraints(object=None) Specifies which bonds or angles should be implemented with constraints. Allowed values are None, HBonds, AllBonds, or HAngles.
rigidWater(boolean=True) If true, water molecules will be fully rigid regardless of the value passed for the constraints argument
implicitSolvent(object=None) If not None, the implicit solvent model to use. Allowed values are HCT, OBC1, OBC2, or GBn
implicitSolventKappa(float=None) Debye screening parameter to model salt concentrations in GB solvent.
implicitSolventSaltConc(float=0.0*u.moles/u.liter) Salt concentration for GB simulations. Converted to Debye length `kappa'
temperature(float=298.15*u.kelvin) Temperature used in the salt concentration-to-kappa conversion for GB salt concentration term
soluteDielectric(float=1.0) The solute dielectric constant to use in the implicit solvent model.
solventDielectric(float=78.5) The solvent dielectric constant to use in the implicit solvent model.
removeCMMotion(boolean=True) If true, a CMMotionRemover will be added to the System.
hydrogenMass(mass=None) The mass to use for hydrogen atoms bound to heavy atoms. Any mass added to a hydrogen is subtracted from the heavy atom to keep their total mass the same.
ewaldErrorTolerance(float=0.0005) The error tolerance to use if the nonbonded method is Ewald or PME.
flexibleConstraints(bool=True) Are our constraints flexible or not?
verbose(bool=False) Optionally prints out a running progress report
def deleteCmap (   self)

Deletes the CMAP terms from the CHARMM PSF.

def loadParameters (   self,
  parmset 
)

Loads parameters from a parameter set that was loaded via CHARMM RTF, PAR, and STR files.

Parameters:
parmset(CharmmParameterSet) List of all parameters
def setBox (   self,
  a,
  b,
  c,
  alpha = 90.0*u.degrees,
  beta = 90.0*u.degrees,
  gamma = 90.0*u.degrees 
)

Sets the periodic box boundary conditions.

Parameters:
a(length) Lengths of the periodic cell
b(length) Lengths of the periodic cell
c(length) Lengths of the periodic cell
alpha(floats, optional) Angles between the periodic cell vectors.
beta(floats, optional) Angles between the periodic cell vectors.
gamma(floats, optional) Angles between the periodic cell vectors.
def system (   self)

Return the cached system class -- it needs to be initialized via "createSystem" first!

def topology (   self)

Create an OpenMM Topology object from the stored bonded network.


Member Data Documentation

int ANGLE_FORCE_GROUP = 1 [static]
int BOND_FORCE_GROUP = 0 [static]
int CMAP_FORCE_GROUP = 5 [static]
int DIHEDRAL_FORCE_GROUP = 2 [static]
int GB_FORCE_GROUP = 6 [static]
int IMPROPER_FORCE_GROUP = 4 [static]
int NONBONDED_FORCE_GROUP = 6 [static]
int UREY_BRADLEY_FORCE_GROUP = 3 [static]

The documentation for this class was generated from the following file:
 All Classes Functions Variables