Table Of Contents

Previous topic

7. Core modules

Next topic

7.2. Fundamental building blocks — MDAnalysis.core.AtomGroup

This Page

7.1. Core functions of MDAnalysis

The basic class is an AtomGroup; the whole simulation is called the Universe. Selections are computed on an AtomGroup and return another AtomGroup.

Timeseries are a convenient way to analyse trajectories.

To get started, load the Universe:

u = Universe(psffilename,dcdfilename)

A simple selection of all water oxygens within 4 A of the protein:

water_shell = u.selectAtoms('name OH2 and around 4.0 protein')
water_shell.numberOfAtoms()   # how many waters were selected
water_shell.totalMass()       # their total mass

AtomGroup instances have various methods that allow calculation of simple properties. For more complicated analysis, obtain the coordinates as a numpy array

coords = water_shell.coordinates()

and write your own Python code.

7.1.1. Flags

(This is an advanced topic and can probably be skipped by most people.)

There are a number flags that influence how MDAnalysis behaves. They are accessible through the pseudo-dictionary MDAnalysis.core.flags.

The entries appear as ‘name’-‘value’ pairs. Flags check values and illegal ones raise a ValueError. Documentation on all flags can be obtained with

print MDAnalysis.core.flags.doc()

7.1.1.1. List of MDAnalysis flags with default values

class MDAnalysis.core.__init__.flagsDocs[source]

use_periodic_selections = True

Determines if distance selections (AROUND, POINT) respect periodicity.

>>> flags['use_periodic_selections'] = value
Values of flag:
  • True - periodicity is taken into account if supported
  • False - periodicity is ignored

The MDAnalysis preset of this flag is True.

Note that KD-tree based distance selections always ignore this flag. (For details see the docs for the ‘use_KDTree_routines’ flag.)

use_KDTree_routines = ‘fast’

Determines which KDTree routines are used for distance selections

>>> flags['use_KDTree_routines'] = value

Values for flag:

  • True, ‘fast’ - only use KDTree routines that are typically faster than others - POINT uses distance matrix routines (with periodicity) - AROUND uses KDTree routines (always ignores periodicity)
  • ‘always’ - always use KDTree routines where available (eg for benchmarking)
  • False, ‘never’ - always use alternatives

The preset value for MDAnalysis is ‘fast’.

MDAnalysis.KDTree routines are significantly faster for some distance selections. However, they cannot deal with periodic boxes and thus ignore periodicity; if periodicity is crucial, disable KDTree routines with

>>> MDAnalysis.core.flags['use_KDTree_routines'] = False

length_unit = ‘Angstrom’

Base unit for lengths (in particular coordinates in trajectories)

>>> flags['length_unit'] = value

Warning

Do not change, only Angstrom fully supported.

use_pbc = False

Choose whether to consider periodic boundary conditions when performing many MDAnalysis.core.AtomGroup.AtomGroup methods. This is set to False by default but can be enabled with:

>>> MDAnalysis.core.flags['use_pbc'] = True

Values for flag:

  • True - Move all atoms within the primary unit cell before calculation
  • False - Use coordinates as supplied

Warning

Changing this to True changes the default behaviour of commonly used MDAnalysis.core.AtomGroup.AtomGroup methods such as MDAnalysis.core.AtomGroup.AtomGroup.centerOfMass() and MDAnalysis.core.AtomGroup.AtomGroup.centerOfGeometry()!

permissive_pdb_reader = True

Select the default reader for PDB Brookhaven databank files.

>>> flags['permissive_pdb_reader'] = value

The Bio.PDB reader (value=``False``) can deal with ‘proper’ PDB files from the Protein Databank that contain special PDB features such as insertion codes and it can auto-correct some common mistakes; see Bio.PDB for details. However, Bio.PDB has been known to read some simulation system PDB files incompletely; a sure sign of problems is a warning that an atom has appeared twice in a residue.

Therefore, the default for the PDB reader is True, which selects the “primitive” (or “permissive”) reader MDAnalysis.coordinates.PDB.PrimitivePDBReader, which essentially just reads ATOM and HETATM lines and puts atoms in a list.

One can manually switch between the two by providing the permissive keyword to MDAnalysis.Universe.

convert_lengths = True

Determine if trajectory reader and writer converts length units between native and MDAnalysis default.

>>> flags['convert_lengths'] = value

Some trajectories are in other length units than the MDAnalysis default (see the flag length_unit); e.g. Gromacs XTC and TRR trajectories are in nm. If True then coordinates are automatically converted, with False the coordinate values are presented as read from the trajectories.

Note

The conversion of lengths also affects conversion of velocities.

force_unit = ‘kJ/(mol*Angstrom)’

Base unit for forces (in particular forces in trajectories)

>>> flags['force_unit'] = value

time_unit = ‘ps’

Base unit for times (in particular time steps in trajectories)

>>> flags['time_unit'] = value

speed_unit = ‘Angstrom/ps’

Base unit for speed (in particular velocities in trajectories)

>>> flags['speed_unit'] = value

charge_unit = ‘e’

Base unit for charge

>>> flags['charge_unit'] = value

x.__init__(...) initializes x; see help(type(x)) for signature

7.1.1.2. Classes

MDAnalysis.core.__init__.flags
class MDAnalysis.core.__init__.Flags(*args)[source]

Global registry of flags. Acts like a dict for item access.

There are a number flags defined that influence how MDAnalysis behaves. They are accessible through the pseudo-dictionary

MDAnalysis.core.flags

The entries appear as ‘name’-‘value’ pairs. Flags check values and illegal ones raise a ValueError. Documentation on all flags can be obtained with

print MDAnalysis.core.flags.__doc__

New flags are added with the Flags.register() method which takes a new Flag instance as an argument.

For DEVELOPERS: Initialize Flags registry with a list of Flag instances.

doc()[source]

Shows doc strings for all flags.

register(flag)[source]

Register a new Flag instance with the Flags registry.

update(*flags)[source]

Update Flags registry with a list of Flag instances.

class MDAnalysis.core.__init__.Flag(name, default, mapping=None, doc=None)[source]

A Flag, essentially a variable that knows its default and legal values.

Create a new flag which will be registered with FLags.

newflag = Flag(name,default,mapping,doc)
Arguments:
name

name of the flag, must be a legal python name

default

default value

mapping

dict that maps allowed input values to canonical values; if None then no argument checking will be performed and all values are directly set.

doc

doc string; may contain string interpolation mappings for:

%%(name)s        name of the flag
%%(default)r     default value
%%(value)r       current value
%%(mapping)r     mapping

Doc strings are generated dynamically and reflect the current state.

prop()[source]

Use this for property(**flag.prop())