MDMC.MD.engine_facades package
Submodules
MDMC.MD.engine_facades.dlpoly_engine module
Facade for DL_POLY MD engine
This is a facade to the DL_POLY MD engine and the Python wrapper dlpoly-py that can interface with it.
- class MDMC.MD.engine_facades.dlpoly_engine.DLPOLYAttribute(dlpoly: DLPoly = None, control: Control = None, config: Config = None, field: Field = None, statis: str = None, output: str = None, dest_config: str = None, rdf: str = None, workdir: str = None)[source]
Bases:
ABC
A class which has a
dlpoly-py
object as an attribute and possesses attributes and methods relating to it.- Parameters:
dlpoly (dlpoly-py, optional) – Set the
dlpoly
attribute to adlpoly-py
object. Default is None, which results in a newdlpoly-py
object being initialised.
- dlpoly
The
dlpoly-py
object owned by this class- Type:
dlpoly-py
- class MDMC.MD.engine_facades.dlpoly_engine.DLPOLYEngine(dlpoly: DLPoly = None, control: Control = None, config: Config = None, field: Field = None, statis: str = None, output: str = None, dest_config: str = None, rdf: str = None, workdir: str = None)[source]
Bases:
DLPOLYAttribute
,MDEngine
Facade for DL_POLY
- HANDLED_PARAMS = {'barostat', 'ensemble', 'field', 'p_damp', 'pressure', 'rescale_step', 't_damp', 't_fraction', 't_window', 'temperature', 'thermostat', 'timestep'}
- property barostat: str
Get or set the str which specifies the barostat
- Returns:
The
barostat
name- Return type:
str
- clear() None [source]
Deletes all atoms of the MD engine, restores all settings to their default values, and frees all memory.
- convert_trajectory(start: int = 0, stop: int = None, step: int = 1, **settings: dict) CompactTrajectory [source]
Parses the trajectory from the
DL_POLY
format into MDMC format.
- property ensemble: DLPOLYEnsemble
Get or set the ensemble object which applies a
thermostat
and/orbarostat
to DL_POLY- Returns:
The simulation thermodynamic ensemble
- Return type:
- minimize(n_steps: int, minimize_every: int = 10, output_log: str = None, work_dir: str = None, **settings: dict) None [source]
Minimizes the simulation energy
- Parameters:
n_steps (int) – Maximum number of steps for the MD run.
minimize_every (int, optional, default 10) – Number of MD steps between two consecutive minimizations.
output_log (str, optional, default None) – file where the output goes.
work_dir (str, optional, default None) – folder where the run happens
**settings – The majority of these are generic but some are specific to the
MDEngine
that is being used.etol (float, energy tolerance criteria for energy minimisation)
ftol (float, force tolerance criteria for force minimisation, active only if non-zero)
maxiter (int, not used in this facade)
maxeval (int, not used in this facade)
- property pressure
Get or set the pressure of the simulation in
katm
- Returns:
Pressure in
katm
- Return type:
float
- run(n_steps: int, equilibration=False, output_log: str = None, work_dir: str = None, **settings: dict) None [source]
Runs a simulation. Must follow a call to
setup_universe()
andsetup_simulation()
.- Parameters:
n_steps (int) – Number of steps for the time integrator.
equilibration (bool (optional, default False)) – If True, just runs an equilibration which does not store the
trajectory
.output_log (str, optional, default None) – file where the output goes.
work_dir (str, optional, default None) – folder where the run happens
**settings – The majority of these are generic but some are specific to the
MDEngine
that is being used.time_equilibration (int, number of steps to run equilibrarion for,) – typically 0 in production stage
traj_calculate (on/off prints trajectory or not)
traj_start (int, at which the trajectory printing shall start)
traj_key (pos,pos-vel,pos-vel-force, level of details for trajectory) – prints positions, positions and velocities, positions, velocities and forces
numprocs (int, number of mpi processes to start to run dl_poly)
- property saved_config: Config
Get the saved configuration of the atomic positions
- Returns:
The atomic positions
- Return type:
Configuration
- setup_simulation(**settings: dict) None [source]
Sets the options required to perform a simulation on a setup
Universe
. Must follow a call tosetup_universe()
.- Parameters:
**settings – The majority of these are generic but some are specific to the
MDEngine
that is being used.
- setup_universe(universe: Universe, **settings: dict) None [source]
Creates the simulation box, the atomic configuration, and the topology in DL_POLY
- Parameters:
universe (Universe) – The MDMC
Universe
which will be setup in DL_POLY.**settings – The majority of these are generic but some are specific to the
MDEngine
that is being used.
- property temperature: float
Get or set the temperature of the simulation in
K
- Returns:
Temperature in
K
- Return type:
float
- class MDMC.MD.engine_facades.dlpoly_engine.DLPOLYEnsemble(dlpoly: DLPoly, temperature: str = None, pressure: float = None, thermostat: str = None, barostat: str = None, **settings: dict)[source]
Bases:
DLPOLYAttribute
A thermodynamic ensemble determined by applying a thermostat and/or barostat
- Parameters:
dlpoly (dlpoly-py) – Set the
dlpoly
attribute to adlpoly-py
object.temperature (float, optional) – Thermostat temperature. Default is None, which is only valid if a
thermostat
is also None.pressure (float, optional) – Barostat pressure. Default is None, which is only valid if a
barostat
is also None.thermostat (str) – Name of a thermostat to be applied.
barostat (str) – Name of a barostat to be applied.
**settings –
ensemble
(str)time_step
(float)t_damp
(int)p_damp
(int)t_window
(float)t_fraction
(float)rescale_step
(int)
- rescale_step
Number of steps between applying temperature rescaling. This only applies to rescale thermostats.
- Type:
- property barostat: str
Get or set the str which specifies the barostat
- Raises
If
self.pressure
has not been set.
- property thermostat: str
Get or set the str which specifies the thermostat
- Raises:
AttributeError – If
self.temperature
has not been set.
- class MDMC.MD.engine_facades.dlpoly_engine.DLPOLYSimulation(universe: Universe, traj_step: int, time_step: float = 1.0, dlpoly=None, **settings: dict)[source]
Bases:
DLPOLYAttribute
The attributes and methods related running a simulation in DL_POLY using a
DLPOLYUniverse
object- Parameters:
universe (Universe) – The MDMC
Universe
used to create theDLPOLYUniverse
.dlpoly (dlpoly-py, optional) – Set the
dlpoly
attribute to adlpoly-py
object. Default is None, which results in a newdlpoly-py
object being initialised.traj_step (int) – How many steps the simulation should take between dumping each
CompactTrajectory
frametime_step (float, optional) – Simulation timestep in
fs
, default is1.
**settings – The majority of these are generic but some are specific to the
MDEngine
that is being used.
- traj_step
Number of simulation steps that elapse between the
CompactTrajectory
being stored.- Type:
- ensemble
Simulation ensemble, which applies a
thermostat
andbarostat
.- Type:
- \*\*settings
The majority of these are generic but some are specific to the
MDEngine
that is being used.
- property barostat: str
Get or set the string which specifies the barostat
- Returns:
The barostat name
- Return type:
str
- property pressure: float
Get or set the pressure of the simulation in
katm
- Returns:
Pressure in
atm
- Return type:
float
- property temperature: float
Get or set the temperature of the simulation in
K
- Returns:
Temperature in
K
- Return type:
float
- class MDMC.MD.engine_facades.dlpoly_engine.DLPOLYUniverse(universe: Universe, dlpoly: DLPoly = None, **settings: dict)[source]
Bases:
DLPOLYAttribute
A class with what would be the equivalent in DL_Poly to the MDMC universe (i.e.the configuration and topology)
- Parameters:
universe (Universe) – The MDMC
Universe
used to create theDLPOLYUniverse
dlpoly (dlpoly-py, optional) – Set the
dlpoly
attribute to adlpoly-py
object. Default is None, which results in a newdlpoly-py
object being initialised.**settings – The majority of these are generic but some are specific to the
MDEngine
that is being used.
- There might be a lot of Attributes needed (see DL_POLYUniverse for example)
- MDMC.MD.engine_facades.dlpoly_engine.convert_unit(value: ndarray | float, unit: Unit = None, to_dlpoly: bool = True)[source]
Converts between MDMC units and DL_POLY real units
- Parameters:
value (array_like or float_like) – The value of the physical property to be converted, in MDMC units. Must derive from either
ndarray
or float.unit (Unit, optional) – The unit of the
value
. If None, thevalue
must possess aunit
attribute i.e. derive fromUnitFloat
orUnitArray
. Default is None.to_dlpoly (bool, optional) – If True the conversion is from MDMC units to DL_POLY units. Default is True.
- Returns:
Value in DL_POLY units if
to_dlpoly
is True, otherwise value in MDMC units. Return type is same asvalue
type.- Return type:
float or numpy.ndarray
MDMC.MD.engine_facades.facade module
Module containing an abstract base class for MD engine facades
- class MDMC.MD.engine_facades.facade.MDEngine[source]
Bases:
ABC
Abstract base class for MD engine facades
- clear() None [source]
Deletes all atoms of the MD engine, restores all settings to their default values, and frees all memory.
- abstract convert_trajectory(start: int = 0, stop: int = None, step: int = 1, **settings: dict) CompactTrajectory [source]
Parses the trajectory from the
MDEngine
format into MDMC format- Parameters:
start (int) – The index of the first trajectory, inclusive.
stop (int) – The index of the last trajectory, exclusive.
step (int) – The step size between trajectories.
**settings –
scaled_positions
(bool)If the
trajectory_file
has scaledpositions
atom_IDs
(list)LAMMPS
ID
of the atoms which should be included. If not passed then all atoms are included in the converted trajectory.
- Returns:
The
CompactTrajectory
from the most recent production simulation- Return type:
CompactTrajectory
- abstract minimize(n_steps: int, minimize_every: int = 10, **settings: dict) None [source]
Minimizes the simulation energy
- property parent_simulation: Simulation
Get or set the simulation that created this engine facade
- Returns:
The Simulation object that created this engine facade
- Return type:
Simulation
- abstract reset_config() None [source]
Resets the configuration of the simulation to that in
saved_config
- abstract run(n_steps: int, equilibration: bool) None [source]
Runs a simulation. Must follow a call to
setup_universe()
andsetup_simulation()
.
- abstract property saved_config: Configuration
Get the saved configuration of the atomic positions
- Returns:
The atomic positions
- Return type:
Configuration
- abstract setup_simulation(**settings: dict) None [source]
Sets the options required to perform a simulation on a setup
Universe
. Must follow a call tosetup_universe()
.- Parameters:
universe (Universe) – A molecular dynamics
Universe
which will be simulated in theMDEngine
.settings** – The majority of these are generic but some are specific to the
MDEngine
that is being used.
- abstract setup_universe(universe: str, **settings: dict) None [source]
Creates a
Universe.configuration
and populates withStructure
- Parameters:
universe (Universe) – A molecular dynamics
Universe
which will be setup in theMDEngine
.**settings – The majority of these are generic but some are specific to the
MDEngine
that is being used.
- property time_step: float
Get the simulation time step in
fs
from the parent simulation- Returns:
Simulation time step in
fs
- Return type:
float
MDMC.MD.engine_facades.facade_factory module
Factory class for generating MD engine facades
- class MDMC.MD.engine_facades.facade_factory.MDEngineFacadeFactory[source]
Bases:
object
Provides a factory for creating facades to
MDEngine
. Any facade within theengine_facades
folder can be created with a str of the classname
, as long as it is a subclass ofMDEngine
.- static create_facade(module_name: str) MDEngine [source]
- Parameters:
module_name (str) – A module name in
engine_facades
. Aliases to these module names are also valid.- Returns:
The specified
MDEngine
, as determined by themodule_name
- Return type:
MDEngine
- static standardise_alias(alias: str) ModuleType [source]
Converts an
alias
into a module name
MDMC.MD.engine_facades.lammps_engine module
Facade for LAMMPS MD engine
This is a facade to PyLammps
(added in 30th-Jul-2016 version), a convenience
wrapper for the LAMMPS Python interface i.e. where Python is extended with
LAMMPS.
Defining all interaction types requires that LAMMPS was built with the MOLECULE package.
Notes
When variables are either passed to or from PyLammps
, the ctypes conversion
can mean that they are unnecessarily cast, particularly from float to int.
This can cause issues as LAMMPS requires certain variables, e.g. number of
steps, to be int. Therefore it is always a good idea to be cast these
variables when they are read from PyLammps
e.g.
int(lmp.variables['steps'].value)
.
A minor bug in LAMMPS (Dec 2018 version) means that nangletypes
returned
by PyLammps
is incorrectly set to ndihedraltypes
. This was corrected in
Mar 2020 release.
- class MDMC.MD.engine_facades.lammps_engine.LAMMPSEngine[source]
Bases:
PyLammpsAttribute
,MDEngine
Facade for LAMMPS
One notable issue with creating the LAMMPS facade is that some
pair_styles
are combinations of both dispersion and coulombic interactions. For example, although LAMMPS haslj/cut
,coul/cut
, andcoul/long
, there is nolj/long
pair_style; it only exists in combination withcoul/long
(i.e.lj/long/coul/long
). This means that the facade has to not just parse each nonbonded interaction, but also has to determine if thatpair_style
should be combined with anotherpair_style
(e.g.lj/long
andcoul/long
have to be combined before passing topair_style
andpair_coeff
). An additional complication of this is that when setting the coefficients (Parameters
) of these interactions, they again have to be set together. So in theory a coulombic interaction needs to know the dispersion interaction it has been paired with an also pass those coefficients when callingpair_coeff
. In practice a simplification can be made, at least in the case of currently implemented pair styles: As coulombicpair_styles
do no take any coefficients, any coulombicpair_styles
that comprise a combinedpair_style
can be ignored for the purposes of setting the coulombicpair_style
; this can be taken care of whenpair_coeff
is called from the correspoding dispersion interaction, as this possesses the dispersion coefficients that also need to be passed when the coefficients of this pair style are set.- property barostat: str
Get or set the str which specifies the barostat
- Returns:
The
barostat
name- Return type:
str
- clear() None [source]
Deletes all atoms of the MD engine, restores all settings to their default values, and frees all memory in LAMMPS.
- convert_trajectory(start: int = 0, stop: int = None, step: int = 1, **settings: dict) CompactTrajectory [source]
Converts between a LAMMPS trajectory dump and an MDMC
CompactTrajectory
The LAMMPS dump must include at least
id
,atom_type
, andxyz
positions
. Thexyz
positions
must be consecutive and in that order. The same is true of thexyz
components of thevelocity
, if they are provided.- Parameters:
start (int) – The index of the first trajectory, inclusive.
stop (int) – The index of the last trajectory, exclusive.
step (int) – The step size between trajectories.
**settings –
scaled_positions
(bool)If the
trajectory_file
has scaledpositions
atom_IDs
(list)LAMMPS
ID
of the atoms which should be included. If not passed then all atoms are included in the converted trajectory.
- Returns:
The MDMC
CompactTrajectory
corresponding to the LAMMPStrajectory_file
.- Return type:
CompactTrajectory
- Raises:
AssertionError – If
universe
is passed, and the number of atoms in thetrajectory_file
is not the same as in theuniverse
.TypeError – If
trajectory_file
describes a triclinic universe.
- property ensemble: LAMMPSEnsemble
Get or set the ensemble object which applies a
thermostat
and/orbarostat
to LAMMPS- Returns:
The simulation thermodynamic ensemble
- Return type:
- minimize(n_steps: int, minimize_every: int = 10, output_log: str = None, work_dir: str = None, **settings: dict) None [source]
Moves the atoms towards a potential energy minimum, by performing an MD simulation interrupted periodically by a structure relaxation. In the end, the configuration with the lowest potential energy reached during the run is kept.
- Parameters:
n_steps (int) – Number of MD simulation steps
minimize_every (int, optional, default 10) – The structure relaxation will be performed every ‘minimize_every’ steps of the MD simulation
output_log (str, optional) – Not used in this facade
work_dir (str, optional) – Not used in this facade
**settings
etol (float, energy tolerance criteria for energy minimisation)
ftol (float, force tolerance criteria for force minimisation, active only if non-zero)
maxiter (int, maximum number of steps in a single structure relaxation)
maxeval (int, maximum number of force calculations in a single structure relaxation)
- property pressure: float
Get or set the pressure of the simulation in
atm
- Returns:
Pressure in
atm
- Return type:
float
- run(n_steps: int, equilibration=False, output_log: str = None, work_dir: str = None, **settings: dict)[source]
Runs a simulation. Must follow a call to
setup_universe()
andsetup_simulation()
.
- property saved_config: ndarray
Get the saved configuration of the atomic positions
- Returns:
The configuration from the start of the run. Each row of the
array
corresponds to the LAMMPS atomID - 1
(offset is due to array zero indexing) and the columns of thearray
are thex
,y
,z
components of theposition
, themass
and thecharge
of eachAtom
.- Return type:
- setup_simulation(**settings: dict) None [source]
Sets simulation parameters in LAMMPS, such as the thermodynamic variables, thermostat/barostat parameters and trajectory settings
- Parameters:
**settings – Passed to
LAMMPSSimulation
- setup_universe(universe: Universe, **settings: dict) None [source]
Creates the simulation box, the atomic configuration, and the topology in LAMMPS
- Parameters:
universe (Universe) – The MDMC
Universe
which will be setup in LAMMPS.**settings –
atom_style
(str)A LAMMPS
atom_style
str. The default setting ofreal
will generally be appropriate.
- property temperature: float
Get or set the temperature of the simulation in
K
- Returns:
Temperature in
K
- Return type:
float
- class MDMC.MD.engine_facades.lammps_engine.LAMMPSEnsemble(lmp, temperature: float = None, pressure: float = None, thermostat: str = None, barostat: str = None, **settings: dict)[source]
Bases:
PyLammpsAttribute
A thermodynamic ensemble determined by applying a thermostat and/or barostat
- Parameters:
lmp (PyLammps) – Set the
lmp
attribute to aPyLammps
object.temperature (float, optional) – Thermostat temperature. Default is None, which is only valid if a
thermostat
is also None.pressure (float, optional) – Barostat pressure. Default is None, which is only valid if a
barostat
is also None.thermostat (str) – Name of a thermostat to be applied.
barostat (str) – Name of a barostat to be applied.
**settings –
time_step
(float)t_damp
(int)p_damp
(int)t_window
(float)t_fraction
(float)rescale_step
(int)
- rescale_step
Number of steps between applying temperature rescaling. This only applies to rescale thermostats.
- Type:
- apply_ensemble_fixes() None [source]
Passes the required LAMMPS
fixes
to apply a specific thermodynamic ensemble to the simulationRemoves all pre-existing thermostat and barostat fixes
- property barostat: str
Get or set the str which specifies the barostat
- Raises:
AttributeError – If
self.pressure
has not been set.
- property p_damp
- remove_ensemble_fixes() None [source]
Removes all LAMMPS
fixes
relating to the ensemble i.e. removes all thermostats and barostatsThis must be done before thermostat and barostat fixes are added, so that there is no conflict with existing thermostat and barostats
fixes
. It is also required for Shake and Rattlefixes
which cannot be added after barostat fixes have been applied.
- property t_damp
- property t_fraction: float
Get or set the fraction by which the
temperature
is rescaled to the target temperatureThis is required for the rescale
thermostat
.- Returns:
Fraction (i.e. between 0.0 and 1.0 inclusive) by which the
temperature
is rescaled- Return type:
float
- Raises:
ValueError – If set to a value outside of 0.0 and 1.0 inclusive.
- property t_window: float
Get or set the
temperature
range inK
in which thetemperature
is not rescaledThis only applies to rescale
thermostats
.- Returns:
temperature range in
K
- Return type:
float
- property thermostat: str
Get or set the str which specifies the thermostat
- Raises:
AttributeError – If
self.temperature
has not been set.
- class MDMC.MD.engine_facades.lammps_engine.LAMMPSSimulation(universe, traj_step: int, time_step: float = 1.0, lmp=None, **settings: dict)[source]
Bases:
PyLammpsAttribute
The attributes and methods related running a simulation in LAMMPS using a
LAMMPSUniverse
object- Parameters:
universe (Universe) – The MDMC
Universe
used to create theLAMMPSUniverse
.traj_step (int) – How many steps the simulation should take between dumping each
CompactTrajectory
frametime_step (float, optional) – Simulation timestep in
fs
, default is1.
lmp (PyLammps, optional) – Set the
lmp
attribute to aPyLammps
object. Default is None, which results in a newPyLammps
object being initialised.**settings –
temperature
(float)skin
(float)neighbor_steps
(int)remove_linear_momentum
(int)remove_angular_momentum
(int)velocity_seed
(int): The seed to be used by LAMMPS to create velocities.
- traj_step
Number of simulation steps that elapse between the
CompactTrajectory
being stored.- Type:
- ensemble
Simulation ensemble, which applies a
thermostat
andbarostat
.- Type:
- property ang_momentum_steps: int
Get or set the number of steps between resetting the angular momentum
- Returns:
Number of steps between the angular momentum being removed
- Return type:
int
- property barostat: str
Get or set the string which specifies the barostat
- Returns:
The barostat name
- Return type:
str
- property lin_momentum_steps: int
Get or set the number of steps between resetting the linear momentum
- Returns:
Number of steps between the linear momentum being removed
- Return type:
int
- property neighbor_steps: int
Get or set the number of steps between neighbor list updates
- Returns:
Number of steps between neighbor list updates
- Return type:
int
- property pressure: float
Get or set the pressure of the simulation in
atm
- Returns:
Pressure in
atm
- Return type:
float
- property skin: float
Get or set the skin distance in
Ang
- Returns:
The skin distance in
Ang
. This distance plus the forcecutoff
distance determines which atom pairs are stored in the neighbor list.- Return type:
float
- property temperature: float
Get or set the temperature of the simulation in
K
- Returns:
Temperature in
K
- Return type:
float
- class MDMC.MD.engine_facades.lammps_engine.LAMMPSUniverse(universe: Universe, lmp: PyLammps = None, **settings: dict)[source]
Bases:
PyLammpsAttribute
A class with what would be the equivalent in LAMMPS to the universe (i.e. the configuration and topology)
- Parameters:
universe (Universe) – The MDMC
Universe
used to create theLAMMPSUniverse
lmp (PyLammps, optional) – Set the
lmp
attribute to aPyLammps
object. Default is None, which results in a newPyLammps
object being initialised.**settings –
atom_style
(str)A LAMMPS
atom_style
string. The default setting ofreal
will generally be appropriate.nonbonded_mix
(str)The name of the formula which determines non-bonded mixing
- atom_dict
A dict with {
MDMC_atom
:LAMMPS_atom_id
}, whereMDMC_atom
is an MDMCAtom
andLAMMPS_atom
is the corresponding LAMMPSAtom
.- Type:
- atom_types
A dict with {
type_ID
:MDMC_atom_group
}, where thetype_ID
is a unique int andMDMC_atom_group
is a list ofAtom
which are identical in terms of element and interactions.- Type:
- atom_type_properties
Each tuple is (
symbol
,mass
) for allatom_types
(ordered) byatom_type
, wheresymbol
is a str specifying the element of theatom_type
andmass
is a float specifying themass
of theatom_type
.- Type:
list of tuples
- proper_ID
A dict of {
proper
:ID pairs
} relating eachDihedralAngle
(withimproper == False
) to a LAMMPSID
.- Type:
- improper_ID
A dict of {
improper
:ID pairs
} relating eachDihedralAngle
(withimproper == True
) to a LAMMPSID
.- Type:
- apply_constraints() None [source]
Adds a constraint
fix
to LAMMPS for all bonds and bond angles which are constrained
- property nonbonded_mix: str
Get or set the formula used to calculate nonbonded interactions between different
atom_types
Options are
geometric
,arithmetic
andsixthpower
, which are defined in the LAMMPS documentation.- Returns:
The name of the formula which determines non-bonded mixing
- Return type:
str
- Raises:
ValueError – str specifies an unsupported mix name
- set_config(config: ndarray) None [source]
Changes the positions of all of the atoms in the LAMMPS wrapper
- Parameters:
config (numpy.ndarray) – The
positions
,mass
andcharge
of theAtom
objects, used to set the LAMMPS configuration. Each row of the array must correspond to the LAMMPSatom ID - 1
(offset is due toarray
zero indexing) and the columns of thearray
must be thex
,y
,z
components of theposition
, themass
and thecharge
of eachAtom
.- Raises:
IndexError – If
config
does not contain the same number of atoms as LAMMPS possesses.
- class MDMC.MD.engine_facades.lammps_engine.PyLammpsAttribute(lmp: PyLammps = None, atom_style: str = 'full')[source]
Bases:
object
A class which has a
PyLammps
object as an attributeIt possesses attributes and methods relating to the
PyLammps
object- Parameters:
lmp (PyLammps, optional) – Set the
lmp
attribute to aPyLammps
object. Default is None, which results in a newPyLammps
object being initialised.atom_style (str, optional) – The LAMMPS
atom_style
, which determines the properties that can be associated which the atoms (e.g.charge
,bonds
). Default isfull
.
- lmp
The
PyLammps
object owned by this class- Type:
PyLammps
- property dumps: list
Get the PyLammps wrapper list of dumps
Dumps are LAMMPS commands which write atom quantities to file for specified timesteps
- property fix_names: list[str]
Get the names of the
fixes
applied in LAMMPS- Returns:
The names of the
fixes
- Return type:
list of str
- property fix_styles: list[str]
Get the styles of the
fixes
applied in LAMMPS- Returns:
The styles of the
fixes
- Return type:
list of str
- property fixes: list[dict]
Get the
PyLammps
wrapper list offixes
- Returns:
Each dict states the
group
,name
andstyle
of a LAMMPS ``fix` which is applied- Return type:
list of dict
- property system_state: namedtuple
Get the
PyLammps
wrapper systemstate
dict- Returns:
Contains the properties of the simulation box.
- Return type:
System
- MDMC.MD.engine_facades.lammps_engine.convert_unit(value: ndarray | float, unit: Unit = None, to_lammps: bool = True)[source]
Converts between MDMC units and LAMMPS real units
- Parameters:
value (array_like or float_like) – The value of the physical property to be converted, in MDMC units. Must derive from either
ndarray
or float.unit (Unit, optional) – The unit of the
value
. If None, thevalue
must possess aunit
attribute i.e. derive fromUnitFloat
orUnitArray
. Default is None.to_lammps (bool, optional) – If True the conversion is from MDMC units to LAMMPS units. Default is True.
- Returns:
Value in LAMMPS units if
to_lammps
is True, otherwise value in MDMC units. Return type is same asvalue
type.- Return type:
float or numpy.ndarray
- MDMC.MD.engine_facades.lammps_engine.parse_all_nonbonded_styles(interactions: NonBondedInteraction) dict[tuple, list[str]] [source]
Converts all
NonBondedInteractions
to LAMMPS pair stylesThis is required because LAMMPS frequently treats
Coulombic
andDisperion
interactions together; these cases need to be dealt with to generate the correct input to LAMMPSpair_styles
. For example, while thepair_styles
buck
,lj/cut
,coul/cut
andcoul/long
can all be passed separately, the dispersive and coulombic styles are combined (dispersive always preceding coulombic) as this is dealt with more efficiently in the LAMMPS engine.buck/long
andlj/long
only exist as part of other pair styles, such asbuck/long/coul/long
andlj/long/coul/long
, so must be combined.- Parameters:
interactions (list of NonBondedInteractions) –
NonBondedInteractions
to be parsed into LAMMPSpair_styles
.- Returns:
A dict with {
pair_styles
:pair_modifiers
} wherepair_styles
is a tuple of all of the combined LAMMPS pair styles corresponding tointeractions
. Each tuple contains the combinedpair_style
as one element, and thecutoffs
as the second element e.g.('lj/cut/coul/cut', 5.0 10.0)
. If additional arguments are required when combining styles, they will be added and returned in the tuple as the third and middle element e.g.('buck/long/coul/long', 'long long', '10.0'
).pair_modifiers
is a list of str for settingpair_modify
for the correspondingpair_style
.- Return type:
dict
- Raises:
ValueError – If
Dispersion
andCoulombic
interactions have different long rangecutoff
, which is not implemented in LAMMPS.ValueError – If long range
Dispersion
not defined in conjunction with a long rangeCoulombic
.
- MDMC.MD.engine_facades.lammps_engine.parse_bonded_coefficients(interaction: BondedInteraction) list[str] [source]
Orders MDMC
Parameter
objects for input to LAMMPSbond_coeff
andangle_coeff
- Parameters:
interaction (BondedInteraction) –
BondedInteraction
where its style and parameters will be parsed.- Returns:
Style and parameters converted to the input format for LAMMPS
bond_coeff
andangle_coeff
- Return type:
list of str
- Raises:
NotImplementedError – If
interaction
has afunction
that has not been implemented in the LAMMPS facade.
- MDMC.MD.engine_facades.lammps_engine.parse_bonded_styles(interaction: BondedInteraction) str [source]
Converts MDMC
InteractionFunction
names forBondedInteractions
to LAMMP bond styles- Parameters:
interaction (BondedInteraction) –
BondedInteraction
to be parsed into LAMMPS bond style.- Returns:
LAMMPS bond style corresponding to
interaction
- Return type:
str
- Raises:
NotImplementedError – If
interaction
has afunction
that has not been implemented in the LAMMPS facade.
- MDMC.MD.engine_facades.lammps_engine.parse_constraint(constraint_algorithm: ConstraintAlgorithm, bonds: list[Bond] = None, bond_ID_dict: dict = None, angles: list[BondAngle] = None, angle_ID_dict: dict = None) list [source]
Converts an MDMC
ConstraintAlgorithm
for input to LAMMPS fixAt least one of bonds and angles must be passed.
- Parameters:
constraint_algorithm (ConstraintAlgorithm) – An object that derives from
ConstraintAlgorithm
to be parsed.bonds (list of Bonds, optional) – Constrained
Bond
interactions.bond_ID_dict (dict, optional) – dict with {
bond
:ID pairs
} relating eachBond
object to a LAMMPSID
.angles (list of BondAngles, optional) – Constrained
BondAngle
interactions.angle_ID_dict (dict, optional) – dict with {
angle
:ID pairs
} relating eachBondAngle
object to a LAMMPSID
.
- Returns:
Input parameters for LAMMPS
fix
, not including the first two terms (fix ID
,group-ID
). The output list has a maximum length of 7, where the last four entries are optional but a minimum of two is required:[algorithm name, accuracy, max iterations, 'b', bond IDs, 'a', angle IDs]
- Return type:
list
- Raises:
TypeError – If there is not at least one
constrained
Interaction
passed.NotImplementedError – If
constraint_algorithm
not been implemented in the LAMMPS facade.
- MDMC.MD.engine_facades.lammps_engine.parse_dispersion_coefficients(interactions: list[NonBondedInteraction], nonbonded_styles: list[tuple] = None) list[str] [source]
Orders MDMC
Parameter
objects for input to LAMMPSpair_coeff
- Parameters:
interactions (list of NonBondedInteraction) –
NonBondedInteractions
where their style and parameters will be parsed.nonbonded_styles (list of tuple) – The parsed
NonBondedInteractions
, where the combinedpair_styles
have been created by theparse_all_nonbonded_styles
function. If None (default), theinteractions
are parsed manually.
- Returns:
Each element is a partial
pair_coeff
command, containing apair_style
, the orderedDispersion
Parameter
objects converted to the correct input format for LAMMPSpair_coeff
, and thecutoff
. For example, if interactions contains aBuckingham
withParameter
A, B, C = 4.184, 1.0, 4.184
, andcutoff = 20.0
, and aCoulombic
with acutoff
of10.0
, the str element of the list will be:'buck/coul/cut 1.0 1.0 1.0 20.0 10.0'
- Return type:
list of str
- Raises:
NotImplementedError – If
interaction
has afunction
that has not been implemented in the LAMMPS facade.ValueError – If the LAMMPS Buckingham parameter rho is less than or equal to zero.
- MDMC.MD.engine_facades.lammps_engine.parse_kspace_solver(solver: KSpaceSolver) list [source]
Converts an MDMC
KSpaceSolver
for input to LAMMPSkspace_style
- Parameters:
solver (KSpaceSolver) – A
KSpaceSolver
to be parsed.- Returns:
Style and parameters for input to LAMMPS
kspace_style
- Return type:
list
- Raises:
NotImplementedError – If
solver
type has has not been implemented in the LAMMPS facade.
- MDMC.MD.engine_facades.lammps_engine.parse_nonbonded_modifications(interaction: Interaction) list[str] [source]
Parses MDMC
Interaction
attributes into list that can be used with LAMMPSpair_modify
command- MDMC.MD.engine_facades.lammps_engine.interaction
The
Interaction
for which the attributes are parsed- Type:
- Returns:
A list of str which are a valid
pair_modify
command, or an empty list if noInteraction
attributes require settingpair_modify
- Return type:
list
- MDMC.MD.engine_facades.lammps_engine.parse_nonbonded_styles(interaction: NonBondedInteraction) tuple [source]
Converts MDMC
InteractionFunction
names forNonBondedInteractions
to LAMMPS pair styles- Parameters:
interaction (NonBondedInteraction) –
NonBondedInteraction
to be parsed into LAMMPS pair style.- Returns:
(
styles
,mods
) wherestyles
is a list of str of LAMMPS pair style corresponding to theinteraction
, andmods
is a list of str of LAMMPS pair modifications corresponding to the interaction- Return type:
tuple
- Raises:
NotImplementedError – If
interaction
has afunction
that has not been implemented in the LAMMPS facade.
Module contents
Facades for molecular dynamics packages
Contents
ff force_field_factory lammps_engine (requires external module lammps.py)