# Selecting which parameters to fit

For complex simulations, there may be a lot of parameters; a simple argon simulation has two parameters, and later on in this guide we will see a water simulation with a total of 8 - this only gets bigger as the molecules and forces become more complex.

Thus, a user may only want to refine a subset of all their parameters. In MDMC this can be done in two ways:

## Fixing a Parameter

Parameter objects can be fixed. This can either be set when they are initialised (created) or changed for an existing Parameter. By default, an initialised Parameter is not fixed.

[1]:

from MDMC.MD.parameters import Parameter
charge = Parameter(value=0.5, name='charge', fixed=True, unit='e')


Attempting to change a fixed parameter produces a warning and does not change the parameter.

[2]:

charge.value = 12
print(f"The value of the charge parameter is: {charge.value}")

The value of the charge parameter is: 0.5 e


For parameters which already exist (for example, a parameter created by a potential), the parameter can be fixed by setting their fixed attribute to True.

[3]:

sigma = Parameter(1.0, name='sigma', unit='Ang')
print('Is sigma fixed: {}'.format(sigma.fixed))
sigma.fixed = True
print('Is sigma fixed: {}'.format(sigma.fixed))

Is sigma fixed: False
Is sigma fixed: True


## Filtering parameters

Only parameters that are passed to Control (as fit_parameters) will be refined. While it is simplest to pass all parameters in a Universe to Control, it is also possible to filter out a subset. To demonstrate this, we use a universe filled with water under a SPCE force field, as used in running a simulation.

[4]:

from MDMC.MD import *

universe = Universe(dimensions=21.75, constraint_algorithm=Shake(1e-4, 100), electrostatic_solver=PPPM(accuracy=1e-5))
H1 = Atom('H')
H2 = Atom('H', position=(0., 1.63298, 0.))
O = Atom('O', position=(0., 0.81649, 0.57736))
H_coulombic = Coulombic(atoms=[H1, H2], cutoff=10.)
O_coulombic = Coulombic(atoms=O, cutoff=10.)
water_mol = Molecule(position=(0, 0, 0),
velocity=(0, 0, 0),
atoms=[H1, H2, O],
interactions=[Bond((H1, O), (H2, O), constrained=True),
BondAngle(H1, O, H2, constrained=True)],
name='water')
universe.fill(water_mol, num_density=0.03356718472021752)
O_dispersion = Dispersion(universe, [O.atom_type, O.atom_type], cutoff=10., vdw_tail_correction=True)

Universe created with:
Dimensions [21.75 21.75 21.75]


There are 8 parameters in the universe:

[5]:

print(universe.parameters)

{'equilibrium_state (#23)': <Parameter
{ID: 23,
type: 'equilibrium_state',
value: 1.0 Ang,
unit: 'Ang',
fixed: False,
constraints: None,
interactions_name: 'Bond',
functions_name: 'HarmonicPotential',
tied: False}>, 'potential_strength (#24)': <Parameter
{ID: 24,
type: 'potential_strength',
value: 4637.0 kJ / mol Ang ^ 2,
unit: 'kJ / mol Ang ^ 2',
fixed: False,
constraints: None,
interactions_name: 'Bond',
functions_name: 'HarmonicPotential',
tied: False}>, 'equilibrium_state (#17)': <Parameter
{ID: 17,
type: 'equilibrium_state',
value: 109.47 deg,
unit: 'deg',
fixed: False,
constraints: None,
interactions_name: 'BondAngle',
functions_name: 'HarmonicPotential',
tied: False}>, 'potential_strength (#18)': <Parameter
{ID: 18,
type: 'potential_strength',
value: 383.0 kJ / mol rad ^ 2,
unit: 'kJ / mol rad ^ 2',
fixed: False,
constraints: None,
interactions_name: 'BondAngle',
functions_name: 'HarmonicPotential',
tied: False}>, 'epsilon (#37)': <Parameter
{ID: 37,
type: 'epsilon',
value: 0.6502 kJ / mol,
unit: 'kJ / mol',
fixed: False,
constraints: None,
interactions_name: 'Dispersion',
functions_name: 'LennardJones',
tied: False}>, 'sigma (#38)': <Parameter
{ID: 38,
type: 'sigma',
value: 3.166 Ang,
unit: 'Ang',
fixed: False,
constraints: None,
interactions_name: 'Dispersion',
functions_name: 'LennardJones',
tied: False}>, 'charge (#4)': <Parameter
{ID: 4,
type: 'charge',
value: 0.4238 e,
unit: 'e',
fixed: False,
constraints: None,
interactions_name: 'Coulombic',
functions_name: 'Coulomb',
tied: False}>, 'charge (#27)': <Parameter
{ID: 27,
type: 'charge',
value: -0.8476 e,
unit: 'e',
fixed: False,
constraints: None,
interactions_name: 'Coulombic',
functions_name: 'Coulomb',
tied: False}>}


So while all 8 parameters can be passed when initiliasing Control, they can also be filtered. Parameters objects have a number of convenience methods to assist with this:

[6]:

parameters = universe.parameters
help(parameters)

Help on Parameters in module MDMC.MD.parameters object:

class Parameters(builtins.dict)
|  Parameters(init_parameters: "Optional[Union[Parameter, 'list[Parameter]']]" = None)
|
|  A dict-like object where every element is a Parameter indexed by name,
|  which contains a number of helper methods for filtering.
|
|  Although Parameters is a dict, it should be treated like a list when writing to it;
|  i.e. initialise it using a list and use append to add to it. These parameters can
|  then be accessed by their name as a key.
|
|  In short; Parameters writes like a list and reads like a dict.
|
|  Parameters
|  ----------
|  init_parameters: Parameter or list of Parameters, optional, default None
|      The initial Parameter objects that the Parameters object contains.
|
|  Attributes
|  ----------
|  array: np.ndarray
|      An alphabetically-sorted numpy array of the Parameters stored in this object.
|
|  Method resolution order:
|      Parameters
|      builtins.dict
|      builtins.object
|
|  Methods defined here:
|
|  __getitem__(self, key: 'str') -> "Union[Parameter, 'list[Parameter]']"
|      Return self[key].
|
|  __init__(self, init_parameters: "Optional[Union[Parameter, 'list[Parameter]']]" = None)
|      Initialize self.  See help(type(self)) for accurate signature.
|
|  __setitem__(self, key: 'str', value: 'Parameter') -> 'NoReturn'
|      Set self[key] to value.
|
|  append(self, parameters: "Union['list[Parameter]', Parameter]") -> 'None'
|      Appends a Parameter or list of Parameters to the dict,
|      with the parameter name as its key.
|
|      Parameters
|      ----------
|      parameters: Parameter or list of Parameters
|          The parameter(s) to be added to the dict.
|
|  filter(self, predicate: 'Callable[[Parameter], bool]') -> 'Parameters'
|      Filters using a predicate
|
|      Parameters
|      ----------
|      predicate : function
|          A function that returns a bool which takes a Parameter as an
|          argument.
|
|      Returns
|      -------
|      Parameters
|          The Parameter objects which meet the condition of the predicate
|
|  filter_atom_attribute(self, attribute: 'str', value: 'Union[str, float]') -> 'Parameters'
|      Filters based on the attribute of Atom objects which have each
|      Parameter applied to them
|
|
|      Parameters
|      ----------
|      attribute : str
|          An attribute of an Atom. Attributes to match to must be either
|          float or str.
|      value : str, float
|          The value of the Atom attribute.
|
|      Returns
|      -------
|      Parameters
|          The Parameter objects which are applied to an Atom object
|          which has the specified value of the specified attribute
|
|  filter_function(self, function_name: 'str') -> 'Parameters'
|      Filters based on the name of the InteractionFunction of each
|      Parameter
|
|      Parameters
|      ----------
|      function_name : str
|          The name of the InteractionFunction of Parameter objects to
|          return, for example 'LennardJones' or 'HarmonicPotential'.
|
|      Returns
|      -------
|      Parameters
|          The Parameter objects which have a function with the
|          specified function_name
|
|  filter_interaction(self, interaction_name: 'str') -> 'Parameters'
|      Filters based on the name of the Interaction of each Parameter
|
|      Parameters
|      ----------
|      interaction_name : str
|          The name of the Interaction of Parameter objects to return,
|          for example 'Bond'.
|
|      Returns
|      -------
|      Parameters
|          The Parameter objects which have an Interaction with the
|          specified interaction_name
|
|  filter_name(self, name: 'str') -> 'Parameters'
|      Filters by name
|
|      Parameters
|      ----------
|      name : str
|          The name of the Parameter objects to return.
|
|      Returns
|      -------
|      Parameters
|          The Parameter objects with name
|
|  filter_structure(self, structure_name: 'str') -> 'Parameters'
|      Filters based on the name of the Structure to which each
|      Parameter applies
|
|      Parameters
|      ----------
|      structure_name : str
|          The name of a Structure.
|
|      Returns
|      -------
|      Parameters
|          The Parameter objects which are applied to a Structure
|          which has the specified structure_name
|
|  filter_value(self, comparison: 'str', value: 'float') -> 'Parameters'
|      Filters by value
|
|      Parameters
|      ----------
|      comparison : str
|          A str representing a comparison operator, '>', '<',
|          '>=', '<=', '==', '!='.
|      value : float
|          A float with which Parameter values are compared, using the
|          comparison operator.
|
|      Returns
|      -------
|      Parameters
|          The Parameter objects which return a True when their values
|          are compared with value using the comparison operator
|
|  log_parameters(self) -> 'None'
|      Logs all Parameters by ID
|
|  ----------------------------------------------------------------------
|
|  as_array
|      The parameters in the object as a sorted numpy array.
|
|      Returns
|      -------
|      np.ndarray
|          An alphabetically-sorted array of parameter values in the object.
|
|  ----------------------------------------------------------------------
|  Data descriptors defined here:
|
|  __dict__
|      dictionary for instance variables
|
|  __weakref__
|      list of weak references to the object
|
|  ----------------------------------------------------------------------
|  Methods inherited from builtins.dict:
|
|  __contains__(self, key, /)
|      True if the dictionary has the specified key, else False.
|
|  __delitem__(self, key, /)
|      Delete self[key].
|
|  __eq__(self, value, /)
|      Return self==value.
|
|  __ge__(self, value, /)
|      Return self>=value.
|
|  __getattribute__(self, name, /)
|      Return getattr(self, name).
|
|  __gt__(self, value, /)
|      Return self>value.
|
|  __ior__(self, value, /)
|      Return self|=value.
|
|  __iter__(self, /)
|      Implement iter(self).
|
|  __le__(self, value, /)
|      Return self<=value.
|
|  __len__(self, /)
|      Return len(self).
|
|  __lt__(self, value, /)
|      Return self<value.
|
|  __ne__(self, value, /)
|      Return self!=value.
|
|  __or__(self, value, /)
|      Return self|value.
|
|  __repr__(self, /)
|      Return repr(self).
|
|  __reversed__(self, /)
|      Return a reverse iterator over the dict keys.
|
|  __ror__(self, value, /)
|      Return value|self.
|
|  __sizeof__(...)
|      D.__sizeof__() -> size of D in memory, in bytes
|
|  clear(...)
|      D.clear() -> None.  Remove all items from D.
|
|  copy(...)
|      D.copy() -> a shallow copy of D
|
|  get(self, key, default=None, /)
|      Return the value for key if key is in the dictionary, else default.
|
|  items(...)
|      D.items() -> a set-like object providing a view on D's items
|
|  keys(...)
|      D.keys() -> a set-like object providing a view on D's keys
|
|  pop(...)
|      D.pop(k[,d]) -> v, remove specified key and return the corresponding value.
|
|      If the key is not found, return the default if given; otherwise,
|      raise a KeyError.
|
|  popitem(self, /)
|      Remove and return a (key, value) pair as a 2-tuple.
|
|      Pairs are returned in LIFO (last-in, first-out) order.
|      Raises KeyError if the dict is empty.
|
|  setdefault(self, key, default=None, /)
|      Insert key with a value of default if key is not in the dictionary.
|
|      Return the value for key if key is in the dictionary, else default.
|
|  update(...)
|      D.update([E, ]**F) -> None.  Update D from dict/iterable E and F.
|      If E is present and has a .keys() method, then does:  for k in E: D[k] = E[k]
|      If E is present and lacks a .keys() method, then does:  for k, v in E: D[k] = v
|      In either case, this is followed by: for k in F:  D[k] = F[k]
|
|  values(...)
|      D.values() -> an object providing a view on D's values
|
|  ----------------------------------------------------------------------
|  Class methods inherited from builtins.dict:
|
|  __class_getitem__(...)
|      See PEP 585
|
|  fromkeys(iterable, value=None, /)
|      Create a new dictionary with keys from iterable and values set to value.
|
|  ----------------------------------------------------------------------
|  Static methods inherited from builtins.dict:
|
|  __new__(*args, **kwargs) class method of builtins.dict
|      Create and return a new object.  See help(type) for accurate signature.
|
|  ----------------------------------------------------------------------
|  Data and other attributes inherited from builtins.dict:
|
|  __hash__ = None



For example, if only the charge parameters should be refined, parameters could be filtered by name:

[7]:

charges = parameters.filter_name('charge')
print(charges)

# charges is a Parameters object
print('\nThe class of charges is: {}'.format(type(charges)))

{'charge (#4)': <Parameter
{ID: 4,
type: 'charge',
value: 0.4238 e,
unit: 'e',
fixed: False,
constraints: None,
interactions_name: 'Coulombic',
functions_name: 'Coulomb',
tied: False}>, 'charge (#27)': <Parameter
{ID: 27,
type: 'charge',
value: -0.8476 e,
unit: 'e',
fixed: False,
constraints: None,
interactions_name: 'Coulombic',
functions_name: 'Coulomb',
tied: False}>}

The class of charges is: <class 'MDMC.MD.parameters.Parameters'>


As each filter returns a Parameters object, filters can be chained together. For example, to find the potential strengths of all bonds:

[8]:

bond_potential_strengths = parameters.filter_name('potential_strength').filter_interaction('Bond')
print(bond_potential_strengths)

# These operations are commutative
print('\nThe order these methods are'
' applied does not matter: {}'.format(bond_potential_strengths
== parameters.filter_interaction('Bond').filter_name('potential_strength')))

{'potential_strength (#24)': <Parameter
{ID: 24,
type: 'potential_strength',
value: 4637.0 kJ / mol Ang ^ 2,
unit: 'kJ / mol Ang ^ 2',
fixed: False,
constraints: None,
interactions_name: 'Bond',
functions_name: 'HarmonicPotential',
tied: False}>}

The order these methods are applied does not matter: True


It is also possible to filter parameters based on the properties of the atoms to which they apply. For instance, we can filter the SPCE parameters so that only parameters of interactions on H atoms are shown:

[9]:

H_parameters = parameters.filter_atom_attribute('name', 'H')
print(H_parameters)

{'equilibrium_state (#23)': <Parameter
{ID: 23,
type: 'equilibrium_state',
value: 1.0 Ang,
unit: 'Ang',
fixed: False,
constraints: None,
interactions_name: 'Bond',
functions_name: 'HarmonicPotential',
tied: False}>, 'potential_strength (#24)': <Parameter
{ID: 24,
type: 'potential_strength',
value: 4637.0 kJ / mol Ang ^ 2,
unit: 'kJ / mol Ang ^ 2',
fixed: False,
constraints: None,
interactions_name: 'Bond',
functions_name: 'HarmonicPotential',
tied: False}>, 'equilibrium_state (#17)': <Parameter
{ID: 17,
type: 'equilibrium_state',
value: 109.47 deg,
unit: 'deg',
fixed: False,
constraints: None,
interactions_name: 'BondAngle',
functions_name: 'HarmonicPotential',
tied: False}>, 'potential_strength (#18)': <Parameter
{ID: 18,
type: 'potential_strength',
value: 383.0 kJ / mol rad ^ 2,
unit: 'kJ / mol rad ^ 2',
fixed: False,
constraints: None,
interactions_name: 'BondAngle',
functions_name: 'HarmonicPotential',
tied: False}>, 'charge (#4)': <Parameter
{ID: 4,
type: 'charge',
value: 0.4238 e,
unit: 'e',
fixed: False,
constraints: None,
interactions_name: 'Coulombic',
functions_name: 'Coulomb',
tied: False}>}


Finally there is also a more flexible method (Parameters.filter) which can be used in conjunction with any function to filter the parameters by some predicate function (a function which takes a parameter and returns True or False according to what is ‘filtered out’). This is most easily done with a lambda function but a regular function can also be used.

[10]:

def is_length(parameter):
"""A function which tells you whether the parameter measures length or not."""
return parameter.unit == 'Ang'

length_parameters = parameters.filter(is_length)
print(length_parameters)

# For those more familiar with Python, this can also be done using a lambda
lambda_length_parameters = parameters.filter(lambda p: p.unit == 'Ang')
print('\nThe same filter can be achieved using lambdas: {}'.format(lambda_length_parameters == length_parameters))

{'equilibrium_state (#23)': <Parameter
{ID: 23,
type: 'equilibrium_state',
value: 1.0 Ang,
unit: 'Ang',
fixed: False,
constraints: None,
interactions_name: 'Bond',
functions_name: 'HarmonicPotential',
tied: False}>, 'sigma (#38)': <Parameter
{ID: 38,
type: 'sigma',
value: 3.166 Ang,
unit: 'Ang',
fixed: False,
constraints: None,
interactions_name: 'Dispersion',
functions_name: 'LennardJones',
tied: False}>}

The same filter can be achieved using lambdas: True