pyretis.engines package

Definition of engines.

This package defines engines for PyRETIS. The engines are responsible for carrying out dynamics for a system. This can in principle both be molecular dynamics or Monte Carlo dynamics. Typically, with RETIS, this will be molecular dynamics in some form in order to propagate the equations of motion and obtain new trajectories.

Package structure

Modules

cp2k.py (pyretis.engines.cp2k)
Defines an engine for use with CP2K.
engine.py (pyretis.engines.engine)
Defines the base engine class.
external.py (pyretis.engines.external)
Defines the interface for external engines.
gromacs.py (pyretis.engines.gromacs)
Defines an engine for use with GROMACS.
gromacs2.py (pyretis.engines.gromacs2)
Defines an engine for use with GROMACS. This is an alternative implementation which does not rely on continuously starting and stopping the GROMACS executable.
internal.py (pyretis.engines.internal)
Defines internal PyRETIS engines.
lammps.py (pyretis.engines.lammps)
Defines and engine for use with LAMMPS.
openmm.py (pyretis.engines.openmm)
Defines an engine for use with OpenMM.

Important methods defined here

engine_factory (engine_factory())
A method to create engines from settings.
pyretis.engines.engine_factory(settings)[source]

Create an engine according to the given settings.

This function is included as a convenient way of setting up and selecting an engine. It will return the created engine.

Parameters:settings (dict) – This defines how we set up and select the engine.
Returns:out (object like EngineBase) – The object representing the engine to use in a simulation.

pyretis.engines.cp2k module

A CP2K external MD integrator interface.

This module defines a class for using CP2K as an external engine.

Important classes defined here

CP2KEngine (CP2KEngine)
A class responsible for interfacing CP2K.
class pyretis.engines.cp2k.CP2KEngine(cp2k, input_path, timestep, subcycles, extra_files=None, seed=0)[source]

Bases: pyretis.engines.external.ExternalMDEngine

A class for interfacing CP2K.

This class defines the interface to CP2K.

cp2k

The command for executing CP2K.

Type:string
input_path

The directory where the input files are stored.

Type:string
timestep

The time step used in the CP2K MD simulation.

Type:float
subcycles

The number of steps each CP2K run is composed of.

Type:integer
rgen

An object we use to set seeds for velocity generation.

Type:object like RandomGenerator
extra_files

List of extra files which may be required to run CP2K.

Type:list
__init__(cp2k, input_path, timestep, subcycles, extra_files=None, seed=0)[source]

Set up the CP2K engine.

Parameters:
  • cp2k (string) – The CP2K executable.
  • input_path (string) – The path to where the input files are stored.
  • timestep (float) – The time step used in the CP2K simulation.
  • subcycles (integer) – The number of steps each CP2K run is composed of.
  • extra_files (list) – List of extra files which may be required to run CP2K.
  • seed (integer, optional) – A seed for the random number generator.
  • extra_files (list) – List of extra files which may be required to run CP2K.
_extract_frame(traj_file, idx, out_file)[source]

Extract a frame from a trajectory file.

This method is used by self.dump_config when we are dumping from a trajectory file. It is not used if we are dumping from a single config file.

Parameters:
  • traj_file (string) – The trajectory file to dump from.
  • idx (integer) – The frame number we look for.
  • out_file (string) – The file to dump to.
static _find_backup_files(dirname)[source]

Return backup-files in the given directory.

_name_output(basename)[source]

Return the name of the output file.

_prepare_shooting_point(input_file)[source]

Create initial configuration for a shooting move.

This creates a new initial configuration with random velocities.

Parameters:input_file (string) – The input configuration to generate velocities for.
Returns:
  • output_file (string) – The name of the file created.
  • energy (dict) – The energy terms read from the CP2K energy file.
_propagate_from(name, path, system, order_function, interfaces, msg_file, reverse=False)[source]

Propagate with CP2K from the current system configuration.

Here, we assume that this method is called after the propagate() has been called in the parent. The parent is then responsible for reversing the velocities and also for setting the initial state of the system.

Parameters:
  • name (string) – A name to use for the trajectory we are generating.
  • path (object like PathBase) – This is the path we use to fill in phase-space points.
  • system (object like System) – The system object gives the initial state.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • msg_file (object like FileIO) – An object we use for writing out messages that are useful for inspecting the status of the current propagation.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

static _read_configuration(filename)[source]

Read CP2K output configuration.

This method is used when we calculate the order parameter.

Parameters:filename (string) – The file to read the configuration from.
Returns:
  • box (numpy.array) – The box dimensions if we manage to read it.
  • xyz (numpy.array) – The positions.
  • vel (numpy.array) – The velocities.
  • names (list of strings) – The atom names found in the file.
_reverse_velocities(filename, outfile)[source]

Reverse velocity in a given snapshot.

Parameters:
  • filename (string) – The configuration to reverse velocities in.
  • outfile (string) – The output file for storing the configuration with reversed velocities.
add_input_files(dirname)[source]

Add required input files to a given directory.

Parameters:dirname (string) – The full path to where we want to add the files.
integrate(system, steps, order_function=None, thermo='full')[source]

Propagate several integration steps.

This method will perform several integration steps using CP2K. It will also calculate the order parameter(s) and energy terms if requested.

Parameters:
  • system (object like System) – The system we are integrating.
  • steps (integer) – The number of steps we are going to perform. Note that we do not integrate on the first step (e.g. step 0) but we do obtain the other properties. This is to output the starting configuration.
  • order_function (object like OrderParameter, optional) – An order function can be specified if we want to calculate the order parameter along with the simulation.
  • thermo (string, optional) – Select the thermodynamic properties we are to calculate.
Yields:

results (dict) – The results from a MD step. This contains the state of the system and order parameter(s) and energies (if calculated).

modify_velocities(system, rgen, sigma_v=None, aimless=True, momentum=False, rescale=None)[source]

Modify the velocities of the current state.

This method will modify the velocities of a time slice.

Parameters:
  • system (object like System) – System is used here since we need access to the particle list.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • sigma_v (numpy.array, optional) – These values can be used to set a standard deviation (one for each particle) for the generated velocities.
  • aimless (boolean, optional) – Determines if we should do aimless shooting or not.
  • momentum (boolean, optional) – If True, we reset the linear momentum to zero after generating.
  • rescale (float, optional) – In some NVE simulations, we may wish to re-scale the energy to a fixed value. If rescale is a float > 0, we will re-scale the energy (after modification of the velocities) to match the given float.
Returns:

  • dek (float) – The change in the kinetic energy.
  • kin_new (float) – The new kinetic energy.

run_cp2k(input_file, proj_name)[source]

Run the CP2K executable.

Returns:out (dict) – The files created by the run.
step(system, name)[source]

Perform a single step with CP2K.

Parameters:
  • system (object like System) – The system we are integrating.
  • name (string) – To name the output files from the CP2K step.
Returns:

out (string) – The name of the output configuration, obtained after completing the step.

pyretis.engines.cp2k.write_for_continue(infile, outfile, timestep, subcycles, name='md_continue')[source]

Create input file for a single step.

Note, the single step actually consists of a number of subcycles. But from PyRETIS’ point of view, this is a single step. Here, we make use of restart files named previous.restart and previous.wfn to continue a run.

Parameters:
  • infile (string) – The input template to use.
  • outfile (string) – The file to create.
  • timestep (float) – The time-step to use for the simulation.
  • subcycles (integer) – The number of sub-cycles to perform.
  • name (string, optional) – A name for the CP2K project.
pyretis.engines.cp2k.write_for_genvel(infile, outfile, posfile, seed, name='genvel')[source]

Create input file for velocity generation.

Parameters:
  • infile (string) – The input template to use.
  • outfile (string) – The file to create.
  • posfile (string) – The (base)name for the input file to read positions from.
  • seed (integer) – A seed for generating velocities.
  • name (string, optional) – A name for the CP2K project.
pyretis.engines.cp2k.write_for_integrate(infile, outfile, timestep, subcycles, posfile, name='md_step', print_freq=None)[source]

Create input file for a single step for the integrate method.

Here, we do minimal changes and just set the time step and subcycles and the starting configuration.

Parameters:
  • infile (string) – The input template to use.
  • outfile (string) – The file to create.
  • timestep (float) – The time-step to use for the simulation.
  • subcycles (integer) – The number of sub-cycles to perform.
  • posfile (string) – The (base)name for the input file to read positions from.
  • name (string, optional) – A name for the CP2K project.
  • print_freq (integer, optional) – How often we should print to the trajectory file.
pyretis.engines.cp2k.write_for_step_vel(infile, outfile, timestep, subcycles, posfile, vel, name='md_step', print_freq=None)[source]

Create input file for a single step.

Note, the single step actually consists of a number of subcycles. But from PyRETIS’ point of view, this is a single step. Further, we here assume that we start from a given xyz file and we also explicitly give the velocities here.

Parameters:
  • infile (string) – The input template to use.
  • outfile (string) – The file to create.
  • timestep (float) – The time-step to use for the simulation.
  • subcycles (integer) – The number of sub-cycles to perform.
  • posfile (string) – The (base)name for the input file to read positions from.
  • vel (numpy.array) – The velocities to set in the input.
  • name (string, optional) – A name for the CP2K project.
  • print_freq (integer, optional) – How often we should print to the trajectory file.

pyretis.engines.engine module

Definition of PyRETIS engines.

This module defines the base class for the engines.

Important classes defined here

EngineBase (EngineBase)
The base class for engines.
class pyretis.engines.engine.EngineBase(description)[source]

Bases: object

Abstract base class for engines.

The engines perform molecular dynamics (or Monte Carlo) and they are assumed to act on a system. Typically they will integrate Newtons equation of motion in time for that system.

description

Short string description of the engine. Used for printing information about the integrator.

Type:string
exe_dir

A directory where the engine is going to be executed.

Type:string
engine_type

Describe the type of engine as an “internal” or “external” engine. If this is undefined, this variable is set to None.

Type:string or None
needs_order

Determines if the engine needs an internal order parameter or not. If not, it is assumed that the order parameter is calculated by the engine.

Type:boolean
__init__(description)[source]

Just add the description.

__str__()[source]

Return the string description of the integrator.

static add_to_path(path, phase_point, left, right)[source]

Add a phase point and perform some checks.

This method is intended to be used by the propagate methods.

Parameters:
  • path (object like PathBase) – The path to add to.
  • phase_point (object like py:class:.System) – The phase point to add to the path.
  • left (float) – The left interface.
  • right (float) – The right interface.
abstract calculate_order(order_function, system, xyz=None, vel=None, box=None)[source]

Obtain the order parameter.

classmethod can_use_order_function(order_function)[source]

Fail if the engine can’t be used with an empty order parameter.

abstract clean_up()[source]

Perform clean up after using the engine.

abstract dump_phasepoint(phasepoint, deffnm=None)[source]

Dump phase point to a file.

engine_type = None
property exe_dir

Return the directory we are currently using.

abstract integration_step(system)[source]

Perform one time step of the integration.

abstract kick_across_middle(system, order_function, rgen, middle, tis_settings)[source]

Force a phase point across the middle interface.

abstract modify_velocities(system, rgen, sigma_v=None, aimless=True, momentum=False, rescale=None)[source]

Modify the velocities of the current state.

Parameters:
  • system (object like System) – The system is used here since we need access to the particle list.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • sigma_v (numpy.array, optional) – These values can be used to set a standard deviation (one for each particle) for the generated velocities.
  • aimless (boolean, optional) – Determines if we should do aimless shooting or not.
  • momentum (boolean, optional) – If True, we reset the linear momentum to zero after generating.
  • rescale (float, optional) – In some NVE simulations, we may wish to re-scale the energy to a fixed value. If rescale is a float > 0, we will re-scale the energy (after modification of the velocities) to match the given float.
Returns:

  • dek (float) – The change in the kinetic energy.
  • kin_new (float) – The new kinetic energy.

needs_order = True
abstract propagate(path, system, order_function, interfaces, reverse=False)[source]

Propagate equations of motion.

static snapshot_to_system(system, snapshot)[source]

Convert a snapshot to a system object.

pyretis.engines.external module

Definition of external engines.

This module defines the base class for external MD engines.

Important classes defined here

ExternalMDEngine (ExternalMDEngine)
The base class for external engines which defines the interface to external programs.
class pyretis.engines.external.ExternalMDEngine(description, timestep, subcycles)[source]

Bases: pyretis.engines.engine.EngineBase

Base class for interfacing external MD engines.

This class defines the interface to external programs. The interface will define how we interact with the external programs, how we write input files for them and read output files from them. New engines should inherit from this class and implement the following methods:

description

A string with a description of the external engine. This can, for instance, be what program we are interfacing with. This is used for outputting information to the user.

Type:string
timestep

The time step used for the external engine.

Type:float
subcycles

The number of steps the external step is composed of. That is: each external step is really composed of subcycles number of iterations.

Type:integer
__init__(description, timestep, subcycles)[source]

Set up the external engine.

Here we just set up some common properties which are useful for the execution.

Parameters:
  • description (string) – A string with a description of the external engine. This can, for instance, be what program we are interfacing with. This is used for outputting information to the user.
  • timestep (float) – The time step used in the simulation.
  • subcycles (integer) – The number of sub-cycles each external integration step is composed of.
_abc_impl = <_abc_data object>
static _copyfile(source, dest)[source]

Copy file from source to destination.

abstract _extract_frame(traj_file, idx, out_file)[source]

Extract a frame from a trajectory file.

Parameters:
  • traj_file (string) – The trajectory file to open.
  • idx (integer) – The frame number we look for.
  • out_file (string) – The file to dump to.
_look_for_input_files(input_path, required, optional_files=None)[source]

Check that required files are present.

The required files are needed to run a simulation with the external engine. We will fail if these are missing. There might also be optional files which are not required, but might be passed in here. If these are not present we will not fail, but delete the reference to this file.

Parameters:
  • input_path (string) – The path to the folder where the input files are stored.
  • required (dict of strings) – These are the file names of the required files.
  • optional_files (list of strings, optional) – These are the file names of the optional files.
Returns:

out (dict) – The paths to the required and optional files we found.

static _modify_input(sourcefile, outputfile, settings, delim='=')[source]

Modify input file for external software.

Here we assume that the input file has a syntax consisting of keyword = setting. We will only replace settings for the keywords we find in the file that is also inside the settings dictionary.

Parameters:
  • sourcefile (string) – The path of the file to use for creating the output.
  • outputfile (string) – The path of the file to write.
  • settings (dict) – A dictionary with settings to write.
  • delim (string, optional) – The delimiter used for separation keywords from settings.
static _movefile(source, dest)[source]

Move file from source to destination.

_name_output(basename)[source]

Return the name of the output file.

abstract _propagate_from(name, path, system, order_function, interfaces, msg_file, reverse=False)[source]

Run the actual propagation using the specific engine.

This method is called after propagate(). And we assume that the necessary preparations before the actual propagation (e.g. dumping of the configuration etc.) is handled in that method.

Parameters:
  • name (string) – A name to use for the trajectory we are generating.
  • path (object like PathBase) – This is the path we use to fill in phase-space points.
  • system (object like System) – The system object gives the initial state.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • msg_file (object like FileIO) – An object we use for writing out messages that are useful for inspecting the status of the current propagation.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

abstract _read_configuration(filename)[source]

Read output configuration from external software.

Parameters:filename (string) – The file to open and read a configuration from.
Returns:
  • out[0] (numpy.array) – The dimensions of the simulation box.
  • out[1] (numpy.array) – The positions found in the given filename.
  • out[2] (numpy.array) – The velocities found in the given filename.
static _read_input_settings(sourcefile, delim='=')[source]

Read input settings for simulation input files.

Here we assume that the input file has a syntax consisting of keyword = setting, where = can be any string given in the input parameter delim.

Parameters:
  • sourcefile (string) – The path of the file to use for creating the output.
  • delim (string, optional) – The delimiter used for separation keywords from settings.
Returns:

settings (dict of strings) – The settings found in the file.

Note

Important: We are here assuming that there will ONLY be one keyword per line.

_remove_files(dirname, files)[source]

Remove files from a directory.

Parameters:
  • dirname (string) – Where we are removing.
  • files (list of strings) – A list with files to remove.
static _removefile(filename)[source]

Remove a given file if it exist.

abstract _reverse_velocities(filename, outfile)[source]

Reverse velocities in a given snapshot.

Parameters:
  • filename (string) – Input file with velocities.
  • outfile (string) – File to write with reversed velocities.
calculate_order(order_function, system, xyz=None, vel=None, box=None)[source]

Calculate order parameter from configuration in a file.

Note, if xyz, vel or box are given, we will NOT read positions, velocity and box information from the current configuration file.

Parameters:
  • order_function (object like OrderParameter) – The class used for calculating the order parameter.
  • system (object like System) – The object the order parameter is acting on.
  • xyz (numpy.array, optional) – The positions to use, in case we have already read them somewhere else. We will then not attempt to read them again.
  • vel (numpy.array, optional) – The velocities to use, in case we already have read them.
  • box (numpy.array, optional) – The current box vectors, in case we already have read them.
Returns:

out (list of floats) – The calculated order parameter(s).

clean_up()[source]

Will remove all files from the current directory.

dump_config(config, deffnm='conf')[source]

Extract configuration frame from a system if needed.

Parameters:
  • config (tuple) – The configuration given as (filename, index).
  • deffnm (string, optional) – The base name for the file we dump to.
Returns:

out (string) – The file name we dumped to. If we did not in fact dump, this is because the system contains a single frame and we can use it directly. Then we return simply this file name.

Note

If the velocities should be reversed, this is handled elsewhere.

dump_frame(system, deffnm='conf')[source]

Just dump the frame from a system object.

dump_phasepoint(phasepoint, deffnm='conf')[source]

Just dump the frame from a system object.

engine_type = 'external'
execute_command(cmd, cwd=None, inputs=None)[source]

Execute an external command for the engine.

We are here executing a command and then waiting until it finishes. The standard out and standard error are piped to files during the execution and can be inspected if the command fails. This method returns the return code of the issued command.

Parameters:
  • cmd (list of strings) – The command to execute.
  • cwd (string or None, optional) – The current working directory to set for the command.
  • inputs (bytes or None, optional) – Additional inputs to give to the command. These are not arguments but more akin to keystrokes etc. that the external command may take.
Returns:

out (int) – The return code of the issued command.

integration_step(system)[source]

Perform a single time step of the integration.

For external engines, it does not make much sense to run single steps unless we absolutely have to. We therefore just fail here. I.e. the external engines are not intended for performing pure MD simulations.

If it’s absolutely needed, there is a self.step() method which can be used, for instance in the initialisation.

Parameters:system (object like System) – A system to run the integration step on.
kick_across_middle(system, order_function, rgen, middle, tis_settings)[source]

Force a phase point across the middle interface.

This is accomplished by repeatedly kicking the phase point so that it crosses the middle interface.

Parameters:
  • system (object like System) – This is the system that we are investigating.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • middle (float) – This is the value for the middle interface.
  • tis_settings (dict) – This dictionary contains settings for TIS. Explicitly used here:
    • zero_momentum: boolean, determines if the momentum is zeroed.
    • rescale_energy: boolean, determines if energy is re-scaled.
Returns:

  • out[0] (object like System) – The phase-point just before the crossing the interface.
  • out[1] (object like System) – The phase-point just after crossing the interface.

Note

This function will update the input system state.

propagate(path, initial_state, order_function, interfaces, reverse=False)[source]

Propagate the equations of motion with the external code.

This method will explicitly do the common set-up, before calling more specialised code for doing the actual propagation.

Parameters:
  • path (object like PathBase) – This is the path we use to fill in phase-space points. We are here not returning a new path - this since we want to delegate the creation of the path to the method that is running propagate.
  • initial_state (object like System) – The system object gives the initial state for the integration.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

abstract step(system, name)[source]

Perform a single step with the external engine.

Parameters:
  • system (object like System) – The system we are integrating.
  • name (string) – To name the output files from the external engine.
Returns:

out (string) – The name of the output configuration, obtained after completing the step.

pyretis.engines.gromacs module

A GROMACS external MD integrator interface.

This module defines a class for using GROMACS as an external engine.

Important classes defined here

GromacsEngine (GromacsEngine)
A class responsible for interfacing GROMACS.
class pyretis.engines.gromacs.GromacsEngine(gmx, mdrun, input_path, timestep, subcycles, maxwarn=0, gmx_format='gro', write_vel=True, write_force=False)[source]

Bases: pyretis.engines.external.ExternalMDEngine

A class for interfacing GROMACS.

This class defines the interface to GROMACS.

gmx

The command for executing GROMACS. Note that we are assuming that we are using version 5 (or later) of GROMACS.

Type:string
mdrun

The command for executing GROMACS mdrun. In some cases, this executable can be different from gmx mdrun.

Type:string
mdrun_c

The command for executing GROMACS mdrun when continuing a simulation. This is derived from the mdrun command.

Type:string
input_path

The directory where the input files are stored.

Type:string
input_files

The names of the input files. We expect to find the keys 'conf', 'input' 'topology'.

Type:dict of strings
ext_time

The time to extend simulations by. It is equal to timestep * subcycles.

Type:float
maxwarn

Setting for the GROMACS grompp -maxwarn option.

Type:integer
ext

This string selects the output format for GROMACS.

Type:string
energy_terms

This lists the energy terms we are to extract from GROMACS.

Type:string (binary)
__init__(gmx, mdrun, input_path, timestep, subcycles, maxwarn=0, gmx_format='gro', write_vel=True, write_force=False)[source]

Set up the GROMACS engine.

Parameters:
  • gmx (string) – The GROMACS executable.
  • mdrun (string) – The GROMACS mdrun executable.
  • input_path (string) – The absolute path to where the input files are stored.
  • timestep (float) – The time step used in the GROMACS MD simulation.
  • subcycles (integer) – The number of steps each GROMACS MD run is composed of.
  • maxwarn (integer, optional) – Setting for the GROMACS grompp -maxwarn option.
  • gmx_format (string, optional) – The format used for GROMACS configurations.
  • write_vel (boolean, optional) – Determines if GROMACS should write velocities or not.
  • write_force (boolean, optional) – Determines if GROMACS should write forces or not.
_execute_grompp(mdp_file, config, deffnm)[source]

Execute the GROMACS preprocessor.

Parameters:
  • mdp_file (string) – The path to the mdp file.
  • config (string) – The path to the GROMACS config file to use as input.
  • deffnm (string) – A string used to name the GROMACS files.
Returns:

out_files (dict) – This dict contains files that were created by the GROMACS preprocessor.

_execute_grompp_and_mdrun(config, deffnm)[source]

Execute GROMACS grompp and mdrun.

Here we use the input file given in the input directory.

Parameters:
  • config (string) – The path to the GROMACS config file to use as input.
  • deffnm (string) – A string used to name the GROMACS files.
Returns:

out_files (dict of strings) – The files created by this command.

_execute_mdrun(tprfile, deffnm)[source]

Execute GROMACS mdrun.

This method is intended as the initial gmx mdrun executed. That is, we here assume that we do not continue a simulation.

Parameters:
  • tprfile (string) – The .tpr file to use for executing GROMACS.
  • deffnm (string) – To give the GROMACS simulation a name.
Returns:

out_files (dict) – This dict contains the output files created by mdrun. Note that we here hard code the file names.

_execute_mdrun_continue(tprfile, cptfile, deffnm)[source]

Continue the execution of GROMACS.

Here, we assume that we have already executed gmx mdrun and that we are to append and continue a simulation.

Parameters:
  • tprfile (string) – The .tpr file which defines the simulation.
  • cptfile (string) – The last checkpoint file (.cpt) from the previous run.
  • deffnm (string) – To give the GROMACS simulation a name.
Returns:

out_files (dict) – The output files created/appended by GROMACS when we continue the simulation.

_extend_and_execute_mdrun(tpr_file, cpt_file, deffnm)[source]

Extend GROMACS and execute mdrun.

Parameters:
  • tpr_file (string) – The location of the “current” .tpr file.
  • cpt_file (string) – The last checkpoint file (.cpt) from the previous run.
  • deffnm (string) – To give the GROMACS simulation a name.
Returns:

out_files (dict) – The files created by GROMACS when we extend.

_extend_gromacs(tprfile, time)[source]

Extend a GROMACS simulation.

Parameters:
  • tprfile (string) – The file to read for extending.
  • time (float) – The time (in ps) to extend the simulation by.
Returns:

out_files (dict) – The files created by GROMACS when we extend.

_extract_frame(traj_file, idx, out_file)[source]

Extract a frame from a .trr, .xtc or .trj file.

If the extension is different from .trr, .xtc or .trj, we will basically just copy the given input file.

Parameters:
  • traj_file (string) – The GROMACS file to open.
  • idx (integer) – The frame number we look for.
  • out_file (string) – The file to dump to.

Note

This will only properly work if the frames in the input trajectory are uniformly spaced in time.

_name_output(basename)[source]

Create a file name for the output file.

This method is used when we dump a configuration to add the correct extension for GROMACS (either gro or g96).

Parameters:basename (string) – The base name to give to the file.
Returns:out (string) – A file name with an extension.
_prepare_shooting_point(input_file)[source]

Create the initial configuration for a shooting move.

This creates a new initial configuration with random velocities. Here, the random velocities are obtained by running a zero-step GROMACS simulation.

Parameters:input_file (string) – The input configuration to generate velocities for.
Returns:
  • output_file (string) – The name of the file created.
  • energy (dict) – The energy terms read from the GROMACS .edr file.
_propagate_from(name, path, system, order_function, interfaces, msg_file, reverse=False)[source]

Propagate with GROMACS from the current system configuration.

Here, we assume that this method is called after the propagate() has been called in the parent. The parent is then responsible for reversing the velocities and also for setting the initial state of the system.

Parameters:
  • name (string) – A name to use for the trajectory we are generating.
  • path (object like PathBase) – This is the path we use to fill in phase-space points.
  • system (object like System) – The system object gives the initial state.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • msg_file (object like FileIO) – An object we use for writing out messages that are useful for inspecting the status of the current propagation.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

_read_configuration(filename)[source]

Read output from GROMACS .g96/gro files.

Parameters:filename (string) – The file to read the configuration from.
Returns:
  • box (numpy.array) – The box dimensions.
  • xyz (numpy.array) – The positions.
  • vel (numpy.array) – The velocities.
_remove_gromacs_backup_files(dirname)[source]

Remove files GROMACS has backed up.

These are files starting with a ‘#’

Parameters:dirname (string) – The directory where we are to remove files.
_reverse_velocities(filename, outfile)[source]

Reverse velocity in a given snapshot.

Parameters:
  • filename (string) – The configuration to reverse velocities in.
  • outfile (string) – The output file for storing the configuration with reversed velocities.
get_energies(energy_file, begin=None, end=None)[source]

Return energies from a GROMACS run.

Parameters:
  • energy_file (string) – The file to read energies from.
  • begin (float, optional) – Select the time for the first frame to read.
  • end (float, optional) – Select the time for the last frame to read.
Returns:

energy (dict fo numpy.arrays) – The energies read from the produced GROMACS xvg file.

integrate(system, steps, order_function=None, thermo='full')[source]

Perform several integration steps.

This method will perform several integration steps using GROMACS. It will also calculate order parameter(s) and energy terms if requested.

Parameters:
  • system (object like System) – The system we are integrating.
  • steps (integer) – The number of steps we are going to perform. Note that we do not integrate on the first step (e.g. step 0) but we do obtain the other properties. This is to output the starting configuration.
  • order_function (object like OrderParameter, optional) – An order function can be specified if we want to calculate the order parameter along with the simulation.
  • thermo (string, optional) – Select the thermodynamic properties we are to calculate.
Yields:

results (dict) – The results from a MD step. This contains the state of the system and order parameter(s) and energies (if calculated).

modify_velocities(system, rgen, sigma_v=None, aimless=True, momentum=False, rescale=None)[source]

Modify the velocities of the current state.

This method will modify the velocities of a time slice.

Parameters:
  • system (object like System) – The system is used here since we need access to the particle list.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • sigma_v (numpy.array, optional) – These values can be used to set a standard deviation (one for each particle) for the generated velocities.
  • aimless (boolean, optional) – Determines if we should do aimless shooting or not.
  • momentum (boolean, optional) – If True, we reset the linear momentum to zero after generating.
  • rescale (float, optional) – In some NVE simulations, we may wish to re-scale the energy to a fixed value. If rescale is a float > 0, we will re-scale the energy (after modification of the velocities) to match the given float.
Returns:

  • dek (float) – The change in the kinetic energy.
  • kin_new (float) – The new kinetic energy.

static rename_energies(gmx_energy)[source]

Rename GROMACS energy terms to PyRETIS convention.

static select_energy_terms(terms)[source]

Select energy terms to extract from GROMACS.

Parameters:terms (string) – This string will name the terms to extract. Currently we only allow for two types of output, but this can be customized in the future.
step(system, name)[source]

Perform a single step with GROMACS.

Parameters:
  • system (object like System) – The system we are integrating.
  • name (string) – To name the output files from the GROMACS step.
Returns:

out (string) – The name of the output configuration, obtained after completing the step.

pyretis.engines.gromacs2 module

A GROMACS external MD integrator interface.

This module defines a class for using GROMACS as an external engine.

Important classes defined here

GromacsEngine2
A class responsible for interfacing GROMACS.
class pyretis.engines.gromacs2.GromacsEngine2(gmx, mdrun, input_path, timestep, subcycles, maxwarn=0, gmx_format='g96', write_vel=True, write_force=False)[source]

Bases: pyretis.engines.gromacs.GromacsEngine

A class for interfacing GROMACS.

This class defines an interface to GROMACS. Attributes are similar to GromacsEngine. In this particular interface, GROMACS is executed without starting and stopping and we rely on reading the output TRR file from GROMACS while a simulation is running.

__init__(gmx, mdrun, input_path, timestep, subcycles, maxwarn=0, gmx_format='g96', write_vel=True, write_force=False)[source]

Set up the GROMACS engine.

Parameters:
  • gmx (string) – The GROMACS executable.
  • mdrun (string) – The GROMACS mdrun executable.
  • input_path (string) – The absolute path to where the input files are stored.
  • timestep (float) – The time step used in the GROMACS MD simulation.
  • subcycles (integer) – The number of steps each GROMACS MD run is composed of.
  • maxwarn (integer, optional) – Setting for the GROMACS grompp -maxwarn option.
  • gmx_format (string, optional) – The format used for GROMACS configurations.
  • write_vel (boolean, optional) – Determines if GROMACS should write velocities or not.
  • write_force (boolean, optional) – Determines if GROMACS should write forces or not.
_propagate_from(name, path, system, order_function, interfaces, msg_file, reverse=False)[source]

Propagate with GROMACS from the current system configuration.

Here, we assume that this method is called after the propagate() has been called in the parent. The parent is then responsible for reversing the velocities and also for setting the initial state of the system.

Parameters:
  • name (string) – A name to use for the trajectory we are generating.
  • path (object like pyretis.core.Path.PathBase) – This is the path we use to fill in phase-space points.
  • system (object like System from pyretis.core.system) – The system object gives the initial state.
  • order_function (object like pyretis.orderparameter.OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • msg_file (object like FileIO) – An object we use for writing out messages that are useful for inspecting the status of the current propagation.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

integrate(system, steps, order_function=None, thermo='full')[source]

Perform several integration steps.

This method will perform several integration steps using GROMACS. It will also calculate order parameter(s) and energy terms if requested.

Parameters:
  • system (object like System) – The system we are integrating.
  • steps (integer) – The number of steps we are going to perform. Note that we do not integrate on the first step (e.g. step 0) but we do obtain the other properties. This is to output the starting configuration.
  • order_function (object like OrderParameter, optional) – An order function can be specified if we want to calculate the order parameter along with the simulation.
  • thermo (string, optional) – Select the thermodynamic properties we are to obtain.
Yields:

results (dict) – The results from a MD step. This contains the state of the system and order parameter(s) and energies (if calculated).

class pyretis.engines.gromacs2.GromacsRunner(cmd, trr_file, edr_file, exe_dir)[source]

Bases: object

A helper class for running GROMACS.

This class handles the reading of the TRR on the fly and it is used to decide when to end the GROMACS execution.

cmd

The command for executing GROMACS.

Type:string
trr_file

The GROMACS TRR file we are going to read.

Type:string
edr_file

A .edr file we are going to read.

Type:string
exe_dir

Path to where we are currently running GROMACS.

Type:string
fileh

The current open file object.

Type:file object
running

The process running GROMACS.

Type:None or object like subprocess.Popen
bytes_read

The number of bytes read so far from the TRR file.

Type:integer
ino

The current inode we are using for the file.

Type:integer
stop_read

If this is set to True, we will stop the reading.

Type:boolean
SLEEP

How long we wait after an unsuccessful read before reading again.

Type:float
data_size

The size of the data (x, v, f, box, etc.) in the TRR file.

Type:integer
header_size

The size of the header in the TRR file.

Type:integer
SLEEP = 0.1
__del__()[source]

Just stop execution and close file.

__enter__()[source]

Start running GROMACS, for a context manager.

__exit__(exc_type, exc_val, exc_tb)[source]

Just stop execution and close file for a context manager.

__init__(cmd, trr_file, edr_file, exe_dir)[source]

Set the GROMACS commands and the files we need.

Parameters:
  • cmd (string) – The command for executing GROMACS.
  • trr_file (string) – The GROMACS TRR file we are going to read.
  • edr_file (string) – A .edr file we are going to read.
  • exe_dir (string) – Path to where we are currently running GROMACS.
check_poll()[source]

Check the current status of the running subprocess.

close()[source]

Close the file, in case that is explicitly needed.

get_gromacs_frames()[source]

Read the GROMACS TRR file on-the-fly.

start()[source]

Start execution of GROMACS and wait for output file creation.

stop()[source]

Stop the current GROMACS execution.

pyretis.engines.gromacs2.get_data(fileh, header)[source]

Read data from the TRR file.

Parameters:
  • fileh (file object) – The file we are reading.
  • header (dict) – The previously read header. Contains sizes and what to read.
Returns:

  • data (dict) – The data read from the file.
  • data_size (integer) – The size of the data read.

pyretis.engines.gromacs2.read_remaining_trr(filename, fileh, start)[source]

Read remaining frames from the TRR file.

Parameters:
  • filename (string) – The file we are reading from.
  • fileh (file object) – The file object we are reading from.
  • start (integer) – The current position we are at.
Yields:
  • out[0] (string) – The header read from the file
  • out[1] (dict) – The data read from the file.
  • out[2] (integer) – The size of the data read.
pyretis.engines.gromacs2.reopen_file(filename, fileh, inode, bytes_read)[source]

Reopen a file if the inode has changed.

Parameters:
  • filename (string) – The name of the file we are working with.
  • fileh (file object) – The current open file object.
  • inode (integer) – The current inode we are using.
  • bytes_read (integer) – The position we should start reading at.
Returns:

  • out[0] (file object or None) – The new file object.
  • out[1] (integer or None) – The new inode.

pyretis.engines.internal module

Definition of numerical MD integrators.

These integrators are representations of engines for performing molecular dynamics. Typically they will propagate Newtons equations of motion in time numerically.

Important classes defined here

MDEngine (MDEngine)
Base class for internal MDEngines.
RandomWalk (RandomWalk)
A Random Walk integrator.
Verlet (Verlet)
A Verlet MD integrator.
VelocityVerlet (VelocityVerlet)
A Velocity Verlet MD integrator.
Langevin (Langevin)
A Langevin MD integrator.
class pyretis.engines.internal.MDEngine(timestep, description, dynamics=None)[source]

Bases: pyretis.engines.engine.EngineBase

Base class for internal MD integrators.

This class defines an internal MD integrator. This class of integrators work with the positions and velocities of the system object directly. Further, we make use of the system object in order to update forces etc.

timestep

Time step for the integration.

Type:float
description

Description of the MD integrator.

Type:string
dynamics

A short string to represent the type of dynamics produced by the integrator (NVE, NVT, stochastic, …).

Type:str
__call__(system)[source]

To allow calling MDEngine(system).

Here, we are just calling self.integration_step(system).

Parameters:system (object like System) – The system we are integrating.
Returns:out (None) – Does not return anything, but will update the particles.
__init__(timestep, description, dynamics=None)[source]

Set up the integrator.

Parameters:
  • timestep (float) – The time step for the integrator in internal units.
  • description (string) – A short description of the integrator.
  • dynamics (string or None, optional) – Description of the kind of dynamics the integrator does.
static calculate_order(order_function, system, xyz=None, vel=None, box=None)[source]

Return the order parameter.

This method is just to help to calculate the order parameter in cases where only the engine can do it.

Parameters:
  • order_function (object like OrderParameter) – The object used for calculating the order parameter(s).
  • system (object like System) – The system containing the current positions and velocities.
  • xyz (numpy.array, optional) – The positions to use. Typically for internal engines, this is not needed. It is included here as it can be used for testing and also to be compatible with the generic function defined by the parent.
  • vel (numpy.array, optional) – The velocities to use.
  • box (numpy.array, optional) – The current box vectors.
Returns:

out (list of floats) – The calculated order parameter(s).

clean_up()[source]

Clean up after using the engine.

Currently, this is only included for compatibility with external integrators.

dump_phasepoint(phasepoint, deffnm=None)[source]

For compatibility with external integrators.

engine_type = 'internal'
integrate(system, steps, order_function=None, thermo='full')[source]

Perform several integration steps.

This method will perform several integration steps, but it will also calculate order parameter(s) if requested and energy terms.

Parameters:
  • system (object like System) – The system we are integrating.
  • steps (integer) – The number of steps we are going to perform. Note that we do not integrate on the first step (e.g. step 0) but we do obtain the other properties. This is to output the starting configuration.
  • order_function (object like OrderParameter, optional) – An order function can be specified if we want to calculate the order parameter along with the simulation.
  • thermo (string, optional) – Select the thermodynamic properties we are to calculate.
Yields:

results (dict) – The result of a MD step. This contains the state of the system and also the order parameter(s) (if calculated) and the thermodynamic quantities (if calculated).

integration_step(system)[source]

Perform a single time step of the integration.

Parameters:system (object like System) – The system we are acting on.
Returns:out (None) – Does not return anything, in derived classes it will typically update the given System.
invert_dt()[source]

Invert the time step for the integration.

Returns:out (boolean) – True if the time step is positive, False otherwise.
kick_across_middle(system, order_function, rgen, middle, tis_settings)[source]

Force a phase point across the middle interface.

This is accomplished by repeatedly kicking the phase point so that it crosses the middle interface.

Parameters:
  • system (object like System) – This is the system that contains the particles we are investigating
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • middle (float) – This is the value for the middle interface.
  • tis_settings (dict) – This dictionary contains settings for TIS. Explicitly used here:
    • zero_momentum: boolean, determines if the momentum is zeroed.
    • rescale_energy: boolean, determines if energy is re-scaled.
Returns:

  • out[0] (object like System) – The phase-point just before the interface.
  • out[1] (object like System) – The phase-point just after the interface.

Note

This function will update the system state.

static modify_velocities(system, rgen, sigma_v=None, aimless=True, momentum=False, rescale=None)[source]

Modify the velocities of the current state.

Parameters:
  • system (object like System) – The system is used here since we need access to the particle list.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • sigma_v (numpy.array, optional) – These values can be used to set a standard deviation (one for each particle) for the generated velocities.
  • aimless (boolean, optional) – Determines if we should do aimless shooting or not.
  • momentum (boolean, optional) – If True, we reset the linear momentum to zero after generating.
  • rescale (float, optional) – In some NVE simulations, we may wish to re-scale the energy to a fixed value. If rescale is a float > 0, we will re-scale the energy (after modification of the velocities) to match the given float.
Returns:

  • dek (float) – The change in the kinetic energy.
  • kin_new (float) – The new kinetic energy.

propagate(path, initial_system, order_function, interfaces, reverse=False)[source]

Generate a path by integrating until a criterion is met.

This function will generate a path by calling the function specifying the integration step repeatedly. The integration is carried out until the order parameter has passed the specified interfaces or if we have integrated for more than a specified maximum number of steps. The given system defines the initial state and the system is reset to its initial state when this method is done.

Parameters:
  • path (object like PathBase) – This is the path we use to fill in phase-space points. We are here not returning a new path - this since we want to delegate the creation of the path (type) to the method that is running propagate.
  • initial_system (object like System) – The system object gives the initial state for the integration. The propagation will not alter the state of the system.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

select_thermo_function(thermo='full')[source]

Select function for calculating thermodynamic properties.

Parameters:thermo (string, or None, optional) – String which selects the kind of thermodynamic output.
Returns:thermo_func (callable or None) – The function matching the requested thermodynamic output.
class pyretis.engines.internal.RandomWalk(timestep, rgen=None, seed=0)[source]

Bases: pyretis.engines.internal.MDEngine

A Random Walker integrator.

This class defines a Random walker integrator.

timestep

The length of the step.

Type:float
__init__(timestep, rgen=None, seed=0)[source]

Set up the Random walker integrator.

Parameters:
  • timestep (float) – The time step in internal units.
  • rgen (string, optional) – This string can be used to pick a particular random generator, which is useful for testing.
  • seed (integer, optional) – A seed for the random generator.
integration_step(system)[source]

Random Walker integration, one time step.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but alters the state of the given system.
class pyretis.engines.internal.Verlet(timestep)[source]

Bases: pyretis.engines.internal.MDEngine

The Verlet MD integrator.

This class defines the Verlet MD integrator.

half_idt

Half of the inverse time step: 0.5 / timestep.

Type:float
timestepsq

Squared time step: timestep**2.

Type:float
previous_pos

Stores the previous positions of the particles.

Type:numpy.array
__init__(timestep)[source]

Set up the Verlet MD integrator.

Parameters:timestep (float) – The time step in internal units.
integration_step(system)[source]

Perform one Verlet integration step.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but alters the state of the given system.
set_initial_positions(particles)[source]

Get initial positions for the Verlet integration.

Parameters:particles (object like Particles) – The initial configuration. Positions and velocities are required.
class pyretis.engines.internal.VelocityVerlet(timestep)[source]

Bases: pyretis.engines.internal.MDEngine

The Velocity Verlet MD integrator.

This class defines the Velocity Verlet integrator.

half_timestep

Half of the timestep.

Type:float
__init__(timestep)[source]

Set up the Velocity Verlet integrator.

Parameters:timestep (float) – The time step in internal units.
integration_step(system)[source]

Velocity Verlet integration, one time step.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but alters the state of the given system.
class pyretis.engines.internal.Langevin(timestep, gamma, rgen=None, seed=0, high_friction=False)[source]

Bases: pyretis.engines.internal.MDEngine

The Langevin MD integrator.

This class defines a Langevin integrator.

rgen

This is the class that handles the generation of random numbers.

Type:object like RandomGenerator
gamma

The friction parameter.

Type:float
high_friction

Determines if we are in the high friction limit and should do the over-damped version.

Type:boolean
init_params

If true, we will initiate parameters for the Langevin integrator when integrate_step is invoked.

Type:boolean
param_high

This contains the parameters for the high friction limit. Here we integrate the equations of motion according to: r(t + dt) = r(t) + dt * f(t)/m*gamma + dr. The items in the dict are:

  • sigma : float standard deviation for the positions, used when drawing dr.
  • bddt : numpy.array Equal to dt*gamma/masses, since the masses is a numpy.array this will have the same shape.
Type:dict
param_iner

This dict contains the parameters for the non-high friction limit where we integrate the equations of motion according to: r(t + dt) = r(t) + c1 * dt * v(t) + c2*dt*dt*a(t) + dr and v(r + dt) = c0 * v(t) + (c1-c2)*dt*a(t) + c2*dt*a(t+dt) + dv. The dict contains:

  • c0 : float Corresponds to c0 in the equation above.
  • a1 : float Corresponds to c1*dt in the equation above.
  • a2 : numpy.array Corresponds to c2*dt*dt/mass in the equation above. Here we divide by the masses in order to use the forces rather than the acceleration. Since the masses might be different for different particles, this will result in a numpy.array with shape equal to the shape of the masses.
  • b1 : numpy.array Corresponds to (c1-c2)*dt/mass in the equation above. Here we also divide by the masses, resulting in a numpy.array.
  • b2 : numpy.array Corresponds to c2*dt/mass in the equation above. Here we also divide by the masses, resulting in a numpy.array.
  • mean : numpy.array (2,) The means for the bivariate Gaussian distribution.
  • cov : numpy.array (2,2) This array contains the covariance for the bivariate Gaussian distribution. param_iner[‘mean’] and param_iner[‘cov’] are used as parameters when drawing dr and dv from the bivariate distribution.
Type:dict

Note

Currently, we are using a multi-normal distribution from numpy. Consider replacing this one as it seems somewhat slow.

__init__(timestep, gamma, rgen=None, seed=0, high_friction=False)[source]

Set up the Langevin integrator.

Actually, it is very convenient to set some variables for the different particles. However, to have a uniform initialisation for the different integrators, we postpone this. This initialisation can be done later by calling explicitly the function self._init_parameters(system) or it will be called the first time self.integration_step is invoked.

Parameters:
  • timestep (float) – The time step in internal units.
  • gamma (float) – The gamma parameter for the Langevin integrator.
  • rgen (string, optional) – This string can be used to pick a particular random generator, which is useful for testing.
  • seed (integer, optional) – A seed for the random generator.
  • high_friction (boolean, optional) – Determines if we are in the high_friction limit and should do the over-damped version.
_init_parameters(system)[source]

Extra initialisation of the Langevin integrator.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but updates self.param.
integration_step(system)[source]

Langevin integration, one time step.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but alters the state of the given system.
integration_step_inertia(system)[source]

Langevin integration, one time step.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but alters the state of the given system.
integration_step_overdamped(system)[source]

Over damped Langevin integration, one time step.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but alters the state of the given system.

pyretis.engines.lammps module

This module defines the class for interfacing LAMMPS.

Important classes defined here

LAMMPSEngine (LAMMPSEngine)
The class responsible for interfacing with LAMMPS.
class pyretis.engines.lammps.LAMMPSEngine(lmp, input_path, subcycles, extra_files=None)[source]

Bases: pyretis.engines.external.ExternalMDEngine

A class for interfacing LAMMPS.

lmp

The command for executing LAMMPS

Type:string
input_path

The directory where the input files are stored.

Type:string
input_files

The names of the input files.

Type:dict of strings
__init__(lmp, input_path, subcycles, extra_files=None)[source]

Set up the LAMMPS engine.

Parameters:
  • lmp (string) – The LAMMPS executable.
  • input_path (string) – The absolute path to where the input files are stored.
  • subcycles (integer) – The frequency of output of data by LAMMPS.
  • extra_files (list of strings, optional) – Additional files needed to run the LAMMPS simulation.
_make_order_fix(ordername)[source]

Create a LAMMPS fix for the order parameter.

Note that we here make LAMMPS also output for step zero.

_propagate_from(name, path, system, order_function, interfaces, msg_file, reverse=False)[source]

Propagate with LAMMPS from the current system configuration.

Parameters:
  • name (string) – A name to use for the trajectory we are generating.
  • path (object like PathBase) – This is the path we use to fill in phase-space points.
  • system (object like System) – The system object gives the initial state.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • msg_file (object like FileIO) – An object we use for writing out messages that are useful for inspecting the status of the current propagation.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

_reverse_velocities(filename, outfile)[source]

Reverse the velocities for a given configuration.

add_input_files(dirname)[source]

Add required input files to a given directory.

Parameters:dirname (string) – The path to the directory where we want to add the files.
calculate_order(order_function, system)[source]

Return the last seen order parameter.

dump_phasepoint(phasepoint, deffnm='conf')[source]

Dump a phase point to a new file.

Note

We do not reverse the velocities here.

integrate(system, steps, order_function=None, thermo='full')[source]

Propagate several integration steps with LAMMPS.

modify_velocities(system, rgen, sigma_v=None, aimless=True, momentum=False, rescale=None)[source]

Modify velocities.

needs_order = False
propagate(path, initial_state, order_function, interfaces, reverse=False)[source]

Propagate the equations of motion with the external code.

Parameters:
  • path (object like PathBase) – This is the path we use to fill in phase-space points. We are here not returning a new path - this since we want to delegate the creation of the path to the method that is running propagate.
  • initial_state (object like System) – The system object gives the initial state for the integration.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.

read_energies(name)[source]

Read energies obtained in a LAMMPS run.

Here, we assume that the energies can be found in a log file in the current execute directory.

read_order_parameters(ordername)[source]

Read order parameters as calculated by LAMMPS.

We assume here that these can be found in the current execute directory in a file named ordername, and further that they contain the step number in the first column, followed by the order parameters in following columns.

run_lammps(system, settings, name)[source]

Execute LAMMPS.

This method will handle input files, run LAMMPS and return some data after the run.

Parameters:
  • system (object like System.) – The system defines the initial state we are using.
  • settings (dict) – This dict contains settings for creating the LAMMPS input file.
  • name (string) – This string is used as a base name for some of the LAMMPS input scripts and output files.
Returns:

  • out[0] (numpy.array) – The order parameters obtained during the run.
  • out[1] (dict of numpy.arrays) – The energies, each dictionary key corresponds to a specific energy term.
  • out[2] (string) – The name of the trajectory created by LAMMPS in this run.

step(system, name)[source]

Perform a single step with LAMMPS.

Parameters:
  • system (object like System) – The system we are integrating.
  • name (string) – To name the output files from the LAMMPS step.
Returns:

out (string) – The name of the output configuration, obtained after completing the step.

pyretis.engines.lammps._add_generate_vel(settings)[source]

Add generation of velocities to the LAMMPS input.

pyretis.engines.lammps._add_order_fix(settings)[source]

Add the fix for printing the order parameter for LAMMPS.

pyretis.engines.lammps._add_stopping_condition(settings)[source]

Add the interface stopping condition to the LAMMPS input.

pyretis.engines.lammps._add_traj_dump(settings)[source]

Add the dump commands for storing the LAMMPS trajectory.

pyretis.engines.lammps.add_to_lammps_input(infile, outfile, to_add)[source]

Add PyRETIS specific settings to an input file for LAMMPS.

This method will append settings needed by PyRETIS to an input file for LAMMPS.

Parameters:
  • infile (string) – Path to a file which contains the LAMMPS settings.
  • outfile (string) – Path to the file we should create.
  • to_add (list of strings) – The settings to add for PyRETIS.
pyretis.engines.lammps.create_lammps_md_input(system, infile, outfile, settings)[source]

Create MD input file for LAMMPS.

We will here write to a new file, by appending text to the given template.

Parameters:
  • system (object like System) – The system contains the current particle state, and this determines the initial configuration and if we are to modify velocities in some way (e.g. reversing them before running or drawing new ones).
  • infile (string) – Path to a file which contains the LAMMPS settings.
  • settings (dict) – The settings we are going to use for creating the LAMMPS file.
Returns:

out (list of strings) – The LAMMPS commands written to the output file.

pyretis.engines.lammps.read_lammps_input(filename)[source]

Read a LAMMPS input file.

This will read in a LAMMPS input file and can be used to read the value of particular settings, for instance the time step.

Parameters:filename (string) – The path to the LAMMPS input file.
Returns:out (list of tuples) – The settings found in the LAMMPS input file. The tuples are on the form (keyword, setting).
pyretis.engines.lammps.read_lammps_log(filename)[source]

Read some info from a LAMMPS log file.

In particular, this method is used to read the thermodynamic output from a simulation (e.g. potential and kinetic energies).

Parameters:filename (string) – The path to the LAMMPS log file.
Returns:out (dict) – A dict containing the data we found in the file.
pyretis.engines.lammps.system_to_lammps(system, reverse_velocities, dimension)[source]

Convert a LAMMPS system into LAMMPS commands for loading it.

This method will convert a system created by LAMMPS into LAMMPS commands, so that LAMMPS can read the given snapshot and use it.

Parameters:
  • system (object like System) – The system we are converting.
  • reverse_velocities (boolean) – True if the velocities in the system are to be reversed.
  • dimension (integer) – The number of dimensions used in the LAMMPS simulation.
Returns:

out (list of strings) – The commands needed by LAMMPS to read the configuration.

pyretis.engines.openmm module

Definition of the OpenMM engine.

This module defines the class for the OpenMM MD engine.

Important classes defined here

OpenMMEngine (OpenMMEngine)
The class for running the OpenMM engine.
class pyretis.engines.openmm.OpenMMEngine(simulation, subcycles=1)[source]

Bases: pyretis.engines.engine.EngineBase

A class for interfacing with OpenMM.

This class defines the interface to OpenMM.

simulation

OpenMM simulation object.

Type:OpenMM.simulation
subcyles

Number of OpenMM steps per PyRETIS step.

Type:int
__init__(simulation, subcycles=1)[source]

Set up the OpenMM Engine.

Parameters:
  • simulation (OpenMM.simulation) – The OpenMM simulation object.
  • subcycles (int) – Number of OpenMM integration steps per PyRETIS step.
calculate_order(order_function, system, xyz=None, vel=None, box=None)[source]

Return the order parameter.

This method is just to help to calculate the order parameter in cases where only the engine can do it.

Parameters:
  • order_function (object like OrderParameter) – The object used for calculating the order parameter(s).
  • system (object like System) – The system containing the current positions and velocities.
  • xyz (numpy.array, optional) – The positions to use. Typically for internal engines, this is not needed. It is included here as it can be used for testing and also to be compatible with the generic function defined by the parent.
  • vel (numpy.array, optional) – The velocities to use.
  • box (numpy.array, optional) – The current box vectors.
Returns:

out (list of floats) – The calculated order parameter(s).

clean_up()[source]

Clean up after using the engine.

Currently, this is only included for compatibility with external integrators.

dump_phasepoint(phasepoint, deffnm=None)[source]

For compatibility with external integrators.

engine_type = 'openmm'
integration_step(system)[source]

Perform one integration step of n subcycles.

Parameters:system (object like System) – The system to integrate/act on. Assumed to have a particle list in system.particles.
Returns:out (None) – Does not return anything, but alters the state of the given system.
kick_across_middle(system, order_function, rgen, middle, tis_settings)[source]

Force a phase point across the middle interface.

This is accomplished by repeatedly kicking the phase point so that it crosses the middle interface.

Parameters:
  • system (object like System) – This is the system that contains the particles we are investigating.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • middle (float) – This is the value for the middle interface.
  • tis_settings (dict) – This dictionary contains settings for TIS. Explicitly used here:
    • zero_momentum: boolean, determines if the momentum is zeroed.
    • rescale_energy: boolean, determines if energy is re-scaled.
Returns:

  • out[0] (object like System) – The phase-point just before the interface.
  • out[1] (object like System) – The phase-point just after the interface.

Note

This function will update the system state.

modify_velocities(system, rgen, sigma_v=None, aimless=True, momentum=False, rescale=None)[source]

Modify the velocities of the current state.

This method will modify the velocities of a time slice. And it is part of the integrator since it, conceptually, fits here: we are acting on the system and modifying it.

Parameters:
  • system (object like System) – The system is used here since we need access to the particle list.
  • rgen (object like RandomGenerator) – This is the random generator that will be used.
  • sigma_v (numpy.array, optional) – These values can be used to set a standard deviation (one for each particle) for the generated velocities.
  • aimless (boolean, optional) – Determines if we should do aimless shooting or not.
  • momentum (boolean, optional) – If True, we reset the linear momentum to zero after generating.
  • rescale (float, optional) – In some NVE simulations, we may wish to re-scale the energy to a fixed value. If rescale is a float > 0, we will re-scale the energy (after modification of the velocities) to match the given float.
Returns:

  • dek (float) – The change in the kinetic energy.
  • kin_new (float) – The new kinetic energy.

propagate(path, initial_system, order_function, interfaces, reverse=False)[source]

Generate a path by integrating until a criterion is met.

This function will generate a path by calling the function specifying the integration step repeatedly. The integration is carried out until the order parameter has passed the specified interfaces or if we have integrated for more than a specified maximum number of steps. The given system defines the initial state and the system is reset to its initial state when this method is done.

Parameters:
  • path (object like PathBase) – This is the path we use to fill in phase-space points. We are here not returning a new path - this since we want to delegate the creation of the path (type) to the method that is running propagate.
  • initial_system (object like System) – The system object gives the initial state for the integration. The initial state is stored and the system is reset to the initial state when the integration is done.
  • order_function (object like OrderParameter) – The object used for calculating the order parameter.
  • interfaces (list of floats) – These interfaces define the stopping criterion.
  • reverse (boolean, optional) – If True, the system will be propagated backward in time.
Returns:

  • success (boolean) – This is True if we generated an acceptable path.
  • status (string) – A text description of the current status of the propagation.