libMobility

Classes

`PSE`(self, periodicityX, periodicityY, ...)	This module computes the RPY mobility in triply periodic boundaries using Ewald splitting with the Positively Split Ewald method [1].
`SelfMobility`(self, periodicityX, ...)	This module ignores hydrodynamic interactions, i.e. the mobility matrix is \(\frac{1}{6\pi\eta a}\mathbb{I}\).
`NBody`(self, periodicityX, periodicityY, ...)	This module computes hydrodynamic interactions using an \(O(N^2)\) algorithm.
`DPStokes`(self, periodicityX, periodicityY, ...)	In the Doubly periodic Stokes geometry (DPStokes), an incompressible fluid exists in a domain which is periodic in the plane and open (or walled) in the third direction.

class libMobility.PSE(self, periodicityX, periodicityY, periodicityZ)

Bases: object

This module computes the RPY mobility in triply periodic boundaries using Ewald splitting with the Positively Split Ewald method [1].

This module will only accept periodic boundary conditions in the three directions.

References

[1] Andrew M. Fiore, Florencio Balboa Usabiaga, Aleksandar Donev, James W. Swan; Rapid sampling of stochastic displacements in Brownian dynamics simulations. J. Chem. Phys. 28 March 2017; 146 (12): 124116. https://doi.org/10.1063/1.4978242

LangevinVelocities(self, dt, kbt, forces, torques)

Computes the hydrodynamic (deterministic and stochastic) velocities according to the Langevin equation,

\[\begin{split}\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix} + \sqrt{\frac{2k_BT}{\Delta t}}\boldsymbol{\mathcal{M}}^{1/2} \boldsymbol{W} + k_BT \partial_{\boldsymbol{q}}\cdot \boldsymbol{\mathcal{M}}.\end{split}\]

If forces and torques are omitted then only the stochastic part is computed. Calling this function is equivalent to calling Mdot, sqrtMdotW and divM in sequence and applying their respective scalar coefficients, but this can be done more efficiently in a combined fashion in some solvers. By default, this is equvialent to an Euler-Maruyama scheme that uses a random finite difference to compute the divergence term, if necessary, which may not be accurate enough for some applications.

Parameters:

dt (float) – Time step \(\Delta t\) for the Langevin equation.
kbt (float) – Boltzmann constant times temperature, \(k_B T\), in units of energy.
forces (array_like, optional) – Forces acting on the particles.
torques (array_like, optional) – Torques acting on the particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

Mdot(self, forces, torques)

Computes the product of the Mobility matrix with a group of forces and/or torques, \(\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix}\).

It is required that setPositions has been called before calling this function.

At least one of the forces or torques must be provided.

Parameters:

forces (array_like, optional) – Forces acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles.
torques (array_like, optional) – Torques acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The linear displacements. The result will have the same format as the forces (or torques if forces are unspecified).
array_like – The angular displacements. The result will have the same format as the torques (or forces if torques are unspecified). This will be None if the solver was initialized with includeAngular=False.

clean(self): Frees any memory allocated by the module.

divM(self, prefactor=1)

Computes the divergence term, \(\boldsymbol{\partial}_\boldsymbol{x}\cdot \boldsymbol{\mathcal{M}}\). It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

initialize(self, viscosity, hydrodynamicRadius, includeAngular = False, tolerance = 0.0001, seed, lanczos_callback)

Initialize the module with a given set of parameters.

Warning

setParameters must be called before this function.

Parameters:

viscosity (float) – Viscosity of the fluid.
hydrodynamicRadius (float) – Hydrodynamic radius of the particles.
includeAngular (bool, optional) – Whether the solver will produce angular velocities. Needed if torques are given. Default is false.
tolerance (float, optional) – Tolerance, used for approximate methods and also for Lanczos (default fluctuation computation). Default is 1e-4.
lanczos_callback (callable, optional) – Callback function to be called during the Lanczos process. It should take two arguments: the current iteration (int) and the current error (float). Default is None, which means no callback will be used.

setParameters(self, psi, Lx, Ly, Lz, shearStrain)

Set the parameters for the PSE solver.

Parameters:

psi (float) – The Splitting parameter.
Lx (float) – The box size in the x direction.
Ly (float) – The box size in the y direction.
Lz (float) – The box size in the z direction.
shearStrain (float) – The shear strain.

setPositions(self, positions): The module will compute the mobility according to this set of positions.

sqrtMdotW(self, prefactor=1.0)

Computes the stochastic contribution, \({\mathcal{M}}^{1/2} \boldsymbol{W}\), where \(\boldsymbol{\mathcal{M}}\) is the grand mobility matrix and \(\boldsymbol{W}\) is a standard normal Gaussian process.

It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear fluctuations. Returned shape will be the same shape as the positions.
array_like – The resulting angular fluctuations. Returned shape will be the same shape as the positions.

class libMobility.SelfMobility(self, periodicityX, periodicityY, periodicityZ)

Bases: object

This module ignores hydrodynamic interactions, i.e. the mobility matrix is \(\frac{1}{6\pi\eta a}\mathbb{I}\). This module will only accept open boundary conditions in the three directions.

LangevinVelocities(self, dt, kbt, forces, torques)

Computes the hydrodynamic (deterministic and stochastic) velocities according to the Langevin equation,

\[\begin{split}\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix} + \sqrt{\frac{2k_BT}{\Delta t}}\boldsymbol{\mathcal{M}}^{1/2} \boldsymbol{W} + k_BT \partial_{\boldsymbol{q}}\cdot \boldsymbol{\mathcal{M}}.\end{split}\]

If forces and torques are omitted then only the stochastic part is computed. Calling this function is equivalent to calling Mdot, sqrtMdotW and divM in sequence and applying their respective scalar coefficients, but this can be done more efficiently in a combined fashion in some solvers. By default, this is equvialent to an Euler-Maruyama scheme that uses a random finite difference to compute the divergence term, if necessary, which may not be accurate enough for some applications.

Parameters:

dt (float) – Time step \(\Delta t\) for the Langevin equation.
kbt (float) – Boltzmann constant times temperature, \(k_B T\), in units of energy.
forces (array_like, optional) – Forces acting on the particles.
torques (array_like, optional) – Torques acting on the particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

Mdot(self, forces, torques)

Computes the product of the Mobility matrix with a group of forces and/or torques, \(\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix}\).

It is required that setPositions has been called before calling this function.

At least one of the forces or torques must be provided.

Parameters:

forces (array_like, optional) – Forces acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles.
torques (array_like, optional) – Torques acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The linear displacements. The result will have the same format as the forces (or torques if forces are unspecified).
array_like – The angular displacements. The result will have the same format as the torques (or forces if torques are unspecified). This will be None if the solver was initialized with includeAngular=False.

clean(self): Frees any memory allocated by the module.

divM(self, prefactor=1)

Computes the divergence term, \(\boldsymbol{\partial}_\boldsymbol{x}\cdot \boldsymbol{\mathcal{M}}\). It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

initialize(self, viscosity, hydrodynamicRadius, includeAngular = False, tolerance = 0.0001, seed, lanczos_callback)

Initialize the module with a given set of parameters.

Warning

setParameters must be called before this function.

Parameters:

viscosity (float) – Viscosity of the fluid.
hydrodynamicRadius (float) – Hydrodynamic radius of the particles.
includeAngular (bool, optional) – Whether the solver will produce angular velocities. Needed if torques are given. Default is false.
tolerance (float, optional) – Tolerance, used for approximate methods and also for Lanczos (default fluctuation computation). Default is 1e-4.
lanczos_callback (callable, optional) – Callback function to be called during the Lanczos process. It should take two arguments: the current iteration (int) and the current error (float). Default is None, which means no callback will be used.

setParameters(self, parameter)

Fixes the parameters of the SelfMobility module. For this module, this is just an example on how to set parameters and it is ignored.

Parameters:: parameter (float) – An example parameter that is not used in this module.

setPositions(self, positions): The module will compute the mobility according to this set of positions.

sqrtMdotW(self, prefactor=1.0)

Computes the stochastic contribution, \({\mathcal{M}}^{1/2} \boldsymbol{W}\), where \(\boldsymbol{\mathcal{M}}\) is the grand mobility matrix and \(\boldsymbol{W}\) is a standard normal Gaussian process.

It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear fluctuations. Returned shape will be the same shape as the positions.
array_like – The resulting angular fluctuations. Returned shape will be the same shape as the positions.

class libMobility.NBody(self, periodicityX, periodicityY, periodicityZ)

Bases: object

This module computes hydrodynamic interactions using an \(O(N^2)\) algorithm. Different hydrodynamic kernels can be chosen depending on the periodicity.

This module only accepts open boundaries in the X and Y directions. The Z direction can be one of:

open: The Rotne-Prager-Yamakawa mobility is used.
single_wall: The Rotne-Prager-Blake mobility is used, with a single wall at the bottom of the simulation box (see setParameters).

LangevinVelocities(self, dt, kbt, forces, torques)

Computes the hydrodynamic (deterministic and stochastic) velocities according to the Langevin equation,

\[\begin{split}\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix} + \sqrt{\frac{2k_BT}{\Delta t}}\boldsymbol{\mathcal{M}}^{1/2} \boldsymbol{W} + k_BT \partial_{\boldsymbol{q}}\cdot \boldsymbol{\mathcal{M}}.\end{split}\]

If forces and torques are omitted then only the stochastic part is computed. Calling this function is equivalent to calling Mdot, sqrtMdotW and divM in sequence and applying their respective scalar coefficients, but this can be done more efficiently in a combined fashion in some solvers. By default, this is equvialent to an Euler-Maruyama scheme that uses a random finite difference to compute the divergence term, if necessary, which may not be accurate enough for some applications.

Parameters:

dt (float) – Time step \(\Delta t\) for the Langevin equation.
kbt (float) – Boltzmann constant times temperature, \(k_B T\), in units of energy.
forces (array_like, optional) – Forces acting on the particles.
torques (array_like, optional) – Torques acting on the particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

Mdot(self, forces, torques)

Computes the product of the Mobility matrix with a group of forces and/or torques, \(\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix}\).

It is required that setPositions has been called before calling this function.

At least one of the forces or torques must be provided.

Parameters:

forces (array_like, optional) – Forces acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles.
torques (array_like, optional) – Torques acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The linear displacements. The result will have the same format as the forces (or torques if forces are unspecified).
array_like – The angular displacements. The result will have the same format as the torques (or forces if torques are unspecified). This will be None if the solver was initialized with includeAngular=False.

clean(self): Frees any memory allocated by the module.

divM(self, prefactor=1)

Computes the divergence term, \(\boldsymbol{\partial}_\boldsymbol{x}\cdot \boldsymbol{\mathcal{M}}\). It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

initialize(self, viscosity, hydrodynamicRadius, includeAngular = False, tolerance = 0.0001, seed, lanczos_callback)

Initialize the module with a given set of parameters.

Warning

setParameters must be called before this function.

Parameters:

viscosity (float) – Viscosity of the fluid.
hydrodynamicRadius (float) – Hydrodynamic radius of the particles.
includeAngular (bool, optional) – Whether the solver will produce angular velocities. Needed if torques are given. Default is false.
tolerance (float, optional) – Tolerance, used for approximate methods and also for Lanczos (default fluctuation computation). Default is 1e-4.
lanczos_callback (callable, optional) – Callback function to be called during the Lanczos process. It should take two arguments: the current iteration (int) and the current error (float). Default is None, which means no callback will be used.

setParameters(self, algorithm = 'advise', Nbatch = -1, NperBatch = -1, wallHeight, delta = 0.001)

Set the parameters for the NBody solver.

Parameters:

algorithm (str) – The algorithm to use. Options are “naive”, “fast”, “block” and “advise”. Default is “advise”.
NBatch (int) – The number of batches to use. If -1 (default), the number of batches is automatically determined.
NperBatch (int) – The number of particles per batch. If -1 (default), the number of particles per batch is automatically determined.
wallHeight (float) – The height of the wall. Only valid if periodicityZ is single_wall.
delta (float) – The finite difference step size for random finite differences. Specified in units of hydrodynamicRadius. Default is 1e-3.

setPositions(self, positions): The module will compute the mobility according to this set of positions.

sqrtMdotW(self, prefactor=1.0)

Computes the stochastic contribution, \({\mathcal{M}}^{1/2} \boldsymbol{W}\), where \(\boldsymbol{\mathcal{M}}\) is the grand mobility matrix and \(\boldsymbol{W}\) is a standard normal Gaussian process.

It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear fluctuations. Returned shape will be the same shape as the positions.
array_like – The resulting angular fluctuations. Returned shape will be the same shape as the positions.

class libMobility.DPStokes(self, periodicityX, periodicityY, periodicityZ)

Bases: object

In the Doubly periodic Stokes geometry (DPStokes), an incompressible fluid exists in a domain which is periodic in the plane and open (or walled) in the third direction. The algorithm is described in [1].

The periodicity must be set to periodic in the X and Y directions. The Z periodicity can be set to open, single_wall, or two_walls. The open option allows for an open boundary condition in the Z direction, while single_wall and two_walls add walls at the bottom and/or top of the simulation box.

References

[1] Aref Hashemi, Raúl P. Peláez, Sachin Natesh, Brennan Sprinkle, Ondrej Maxian, Zecheng Gan, Aleksandar Donev; Computing hydrodynamic interactions in confined doubly periodic geometries in linear time. J. Chem. Phys. 21 April 2023; 158 (15): 154101. https://doi.org/10.1063/5.0141371

LangevinVelocities(self, dt, kbt, forces, torques)

Computes the hydrodynamic (deterministic and stochastic) velocities according to the Langevin equation,

\[\begin{split}\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix} + \sqrt{\frac{2k_BT}{\Delta t}}\boldsymbol{\mathcal{M}}^{1/2} \boldsymbol{W} + k_BT \partial_{\boldsymbol{q}}\cdot \boldsymbol{\mathcal{M}}.\end{split}\]

If forces and torques are omitted then only the stochastic part is computed. Calling this function is equivalent to calling Mdot, sqrtMdotW and divM in sequence and applying their respective scalar coefficients, but this can be done more efficiently in a combined fashion in some solvers. By default, this is equvialent to an Euler-Maruyama scheme that uses a random finite difference to compute the divergence term, if necessary, which may not be accurate enough for some applications.

Parameters:

dt (float) – Time step \(\Delta t\) for the Langevin equation.
kbt (float) – Boltzmann constant times temperature, \(k_B T\), in units of energy.
forces (array_like, optional) – Forces acting on the particles.
torques (array_like, optional) – Torques acting on the particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

Mdot(self, forces, torques)

Computes the product of the Mobility matrix with a group of forces and/or torques, \(\boldsymbol{\mathcal{M}}\begin{bmatrix}\boldsymbol{F}\\\boldsymbol{T}\end{bmatrix}\).

It is required that setPositions has been called before calling this function.

At least one of the forces or torques must be provided.

Parameters:

forces (array_like, optional) – Forces acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles.
torques (array_like, optional) – Torques acting on the particles. Must have shape (N, 3) or (3N), where N is the number of particles. The solver must have been initialized with includeAngular=True.

Returns:

array_like – The linear displacements. The result will have the same format as the forces (or torques if forces are unspecified).
array_like – The angular displacements. The result will have the same format as the torques (or forces if torques are unspecified). This will be None if the solver was initialized with includeAngular=False.

clean(self): Frees any memory allocated by the module.

divM(self, prefactor=1)

Computes the divergence term, \(\boldsymbol{\partial}_\boldsymbol{x}\cdot \boldsymbol{\mathcal{M}}\). It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear displacements. Returned shape will be the same as the input forces if given, or the positions if no forces are given.
array_like – The resulting angular displacements. Returned shape will be the same as the input torques if given, or the forces/positions if no torques are given. This array will be None if the solver was initialized with includeAngular=False.

initialize(self, viscosity, hydrodynamicRadius, includeAngular = False, tolerance = 0.0001, seed, lanczos_callback)

Initialize the module with a given set of parameters.

Warning

setParameters must be called before this function.

Parameters:

viscosity (float) – Viscosity of the fluid.
hydrodynamicRadius (float) – Hydrodynamic radius of the particles.
includeAngular (bool, optional) – Whether the solver will produce angular velocities. Needed if torques are given. Default is false.
tolerance (float, optional) – Tolerance, used for approximate methods and also for Lanczos (default fluctuation computation). Default is 1e-4.
lanczos_callback (callable, optional) – Callback function to be called during the Lanczos process. It should take two arguments: the current iteration (int) and the current error (float). Default is None, which means no callback will be used.

setParameters(self, Lx, Ly, zmin, zmax, allowChangingBoxSize=False, delta=0.001, allowUnsafeForces=False)

When the periodicity is set to single_wall a wall in the bottom of the domain is added. When the periodicity is set to two_walls a wall in the bottom and top of the domain is added.

Even in open mode (Z periodicity set to open) the values of zmin and zmax are still required. The algorithm needs to define a grid in the z direction, and these values define the extents of that grid. The code will fail if a position outside of these extents is used.

Parameters:

Lx (float) – The box size in the x direction.
Ly (float) – The box size in the y direction.
zmin (float) – The minimum value of the z coordinate. This is the position of the bottom wall if the Z periodicity is single_wall or two_walls.
zmax (float) – The maximum value of the z coordinate. This is the position of the top wall if the Z periodicity is two_walls.
allowChangingBoxSize (bool) – Whether the periodic extents Lx & Ly can be modified during parameter selection. Default: false.
delta (float) – The finite difference step size for random finite differences. Specified in units of hydrodynamicRadius. Default: 1e-3.
allowUnsafeForces (bool) – When the Z periodicity is open the net force/torque on the domain must be zero. Numerically, this is checked to some tolerance. Enabling this flag will by bypass the check but can lead to unphysical results if used incorrectly. Default: false.

setPositions(self, positions): The module will compute the mobility according to this set of positions.

sqrtMdotW(self, prefactor=1.0)

Computes the stochastic contribution, \({\mathcal{M}}^{1/2} \boldsymbol{W}\), where \(\boldsymbol{\mathcal{M}}\) is the grand mobility matrix and \(\boldsymbol{W}\) is a standard normal Gaussian process.

It is required that setPositions has been called before calling this function.

Parameters:

prefactor (float, optional) – Prefactor to multiply the result by. Default is 1.0.

Returns:

array_like – The resulting linear fluctuations. Returned shape will be the same shape as the positions.
array_like – The resulting angular fluctuations. Returned shape will be the same shape as the positions.