DensityProfile¶

Bases: DynamicAnalysisBase

Serial and parallel implementations to calculate the number and charge density profiles \(\rho_i(z)\) and \(\rho_q(z)\) of a constant-volume system along the specified axes.

The microscopic number density profile of species \(i\) is calculated by binning particle positions along an axis \(z\) using

\[\rho_i(z)=\frac{V}{N_\mathrm{bins}}\left\langle \sum_\alpha\delta(z-z_\alpha)\right\rangle\]

where \(V\) is the system volume and \(N_\mathrm{bins}\) is the number of bins. The angular brackets denote an ensemble average.

If the species carry charges, the charge density profile can be obtained using

\[\rho_q(z)=e\sum_i z_i\rho_i(z)\]

where \(z_i\) is the charge number of species \(i\) and \(e\) is the elementary charge.

With the charge density profile, the surface charge density is given by

\[\sigma_q=\frac{1}{L}\left( \varepsilon_0\varepsilon_\mathrm{r}\Delta\Psi -\int_0^L z\rho_q(z)\,\mathrm{d}z\right)\]

where \(\varepsilon_0\) and \(\varepsilon_\mathrm{r}\) are the vacuum and relative permittivities, respectively, \(\Delta\Psi\) is the potential difference, and \(L\) is the system dimension, and the potential profile can be computed by numerically solving Poisson’s equation for electrostatics:

\[\varepsilon_0\varepsilon_\mathrm{r}\nabla^2\Psi(z)=-\rho_q(z)\]

Parameters:

groupsMDAnalysis.AtomGroup or array-like

Groups of atoms for which density profiles are calculated.

groupingsstr or array-like, default: "atoms"

Determines whether the centers of mass are used in lieu of individual atom positions. If groupings is a str, the same value is used for all groups.

Note

If the desired grouping is not "atoms",

the trajectory file should have segments (or chains) containing residues (or molecules) and residues containing atoms, and
residues and segments should be locally unwrapped at the simulation box edges, if not already, using MDAnalysis.transformations.wrap.unwrap, MDAnalysis.core.groups.AtomGroup.unwrap(), or MDAnalysis.lib.mdamath.make_whole().

Valid values:

"atoms": Atom positions (generally or for coarse-grained simulations).
"residues": Residues’ centers of mass (for atomistic simulations).
"segments": Segments’ centers of mass (for atomistic polymer simulations).

axesint, str, or array-like, default: "xyz"

Axes along which to compute the density profiles.

Examples:

2 for the \(z\)-direction.
"xy" for the \(x\)- and \(y\)-directions.
(0, 1) for the \(x\)- and \(y\)-directions.

n_binsint or array-like

Number of histogram bins \(N_\mathrm{bins}\) for each axis in axes. If an int is provided, the same value is used for all axes.

chargesarray-like, openmm.unit.Quantity, or pint.Quantity, keyword-only, optional

Charges \(q_i\) for the entities in the \(N_\mathrm{groups}\) atom groups in groups. If not provided, they will be retrieved from the main MDAnalysis.core.universe.Universe object only if it contains charge information.

Note

Depending on the grouping for a specific atom group, all entities (atoms, residues, or segments) must carry the same charge. Otherwise, the charge density contribution for that atom group would not make sense. If this condition does not hold, change how the atoms are grouped in the atom groups so that all entities share the same charge.

Shape: \((N_\mathrm{groups},)\).

Reference unit: \(\mathrm{e}\).

dimensionsarray-like, openmm.unit.Quantity, or pint.Quantity, keyword-only, optional

System dimensions \((L_x,\,L_y,\,L_z)\). If the MDAnalysis.core.universe.Universe object that the atom groups in groups belong to does not contain dimensionality information, provide it here. Affected by dim_scales.

Shape: \((3,)\).

Reference unit: \(\mathrm{Å}\).

dim_scalesfloat or array-like, keyword-only, optional

Scale factors for the system dimensions. If an int is provided, the same value is used for all axes.

Shape: \((3,)\).

averagebool, keyword-only, default: True

Determines whether the density profiles are averaged over the \(N_\mathrm{frames}\) analysis frames.

recenterint, list, MDAnalysis.AtomGroup, or tuple, keyword-only, optional

Constrains the center of mass of an atom group by adjusting the particle coordinates every analysis frame. Either specify an MDAnalysis.core.groups.AtomGroup, its index within groups, a list of atom groups or their indices, or a tuple containing the aforementioned information and the fixed center of mass coordinates, in that order. To avoid recentering in a specific dimension, set the coordinate to numpy.nan. If the center of mass is not specified, the center of the simulation box is used.

Shape: \((3,)\) for the fixed center of mass.

reducedbool, keyword-only, default: False

Specifies whether the data is in reduced units. Affects results.number_densities, results.charge_densities, etc.

parallelbool, keyword-only, default: False

Determines whether the analysis is performed in parallel.

Note

The Joblib threading backend generally provides the best performance.

verbosebool, keyword-only, default: True

Determines whether detailed progress is shown.

**kwargs

Additional keyword arguments to pass to MDAnalysis.analysis.base.AnalysisBase.

Attributes:

universeMDAnalysis.Universe

MDAnalysis.core.universe.Universe object containing all information describing the simulation system.

axestuple

Axes along which the density profiles are calculated.

results.unitsdict

Reference units for the results. For example, to get the reference units for results.bins, call results.units["bins"].

results.timesnumpy.ndarray

Times \(t\). Only available if average=False.

Shape: \((N_\mathrm{frames},)\).

Reference unit: \(\mathrm{ps}\).

results.binsdict

Bin centers \(z\) corresponding to the density profiles in each dimension. The key is the axis, e.g., results.bins["z"] for the \(z\)-axis.

Shape: Each array has shape \((N_\mathrm{bins},)\).

Reference unit: \(\mathrm{Å}\).

results.bin_edgesdict

Bin edges corresponding to the density profiles in each dimension. The key is the axis, e.g., results.bin_edges["z"] for the \(z\)-axis.

Shape: Each array has shape \((N_\mathrm{bins}+1,)\).

Reference unit: \(\mathrm{Å}\).

results.number_densitiesdict

Number density profiles \(\rho(z)\). The key is the axis, e.g., results.number_densities["z"] for the \(z\)-axis.

Shape: Each array has shape \((N_\mathrm{groups},\,N_\mathrm{bins})\). If average=False, an additional second dimension of length \(N_\mathrm{frames}\) is present.

Reference unit: \(\mathrm{Å}^{-3}\).

results.charge_densitiesdict

Charge density profiles \(\rho_q(z)\). Only available if charge information was found or provided. The key is the axis, e.g., results.charge_densities["z"] for the \(z\)-axis.

Shape: Each array has shape \((N_\mathrm{bins},)\). If average=False, an additional first dimension of length \(N_\mathrm{frames}\) is present.

Reference unit: \(\mathrm{e/Å}^{-3}\).

results.surface_charge_densitiesnumpy.ndarray

Surface charge densities \(\sigma_q\). Only available after running calculate_surface_charge_densities().

Shape: \((N_\mathrm{axes},)\) or \((N_\mathrm{axes},\,N_\mathrm{frames})\).

results.potentialsdict

Potential profiles \(\Psi(z)\). Only available after running calculate_potential_profiles(). The key is the axis, e.g., results.potentials["z"] for the \(z\)-axis.

Shape: Each array has shape \((N_\mathrm{bins},)\). If average=False, an additional second dimension of length \(N_\mathrm{frames}\) is present.

Reference unit: \(\mathrm{V}\).

Methods

`calculate_potential_profiles`	Calculates the potential profiles in the specified dimensions using the charge density profiles by numerically solving Poisson's equation for electrostatics.
`calculate_surface_charge_densities`	Calculates the surface charge densities \(\sigma_q\) for the specified system dimensions using the charge density profiles \(\rho_q(z)\).
`get_supported_backends`	Tuple with backends supported by the core library for a given class.
`run`	Performs the calculation.
`save`	Saves results to a binary or archive file in NumPy format.

Calculates the potential profiles in the specified dimensions using the charge density profiles by numerically solving Poisson’s equation for electrostatics.

Parameters:

axesstr or array-like, optional

Axes along which to compute the potential profiles. If not specified, all axes for which charge density profiles were calculated will be used.

Examples:"xy" or ("x", "y") for the \(x\)- and \(y\)-directions.

dielectricsfloat, optional

Relative permittivities or dielectric constants \(\varepsilon_\mathrm{r}\). Only optional if previously provided to another calculation method in this class.

sigmas_qfloat, array-like, openmm.unit.Quantity, or pint.Quantity, keyword-only, optional

Surface charge densities \(\sigma_q\). Used to ensure that the electric field in the bulk of the solution is zero. If not provided, it is determined using dVs and the charge density profiles, or the average values in the centers of the integrated charge density profiles.

Note

\(\sigma_q\) and \(\Delta\Psi\) should have the same sign.

Shapes: \((N_\mathrm{axes},)\) or \((N_\mathrm{axes},\,N_\mathrm{frames})\).

Reference unit: \(\mathrm{e/Å^2}\).

dVsfloat, array-like, openmm.unit.Quantity, or pint.Quantity, keyword-only, optional

Potential differences \(\Delta\Psi\) across the system dimensions specified in axes. Can be retrieved if previously provided to another calculation method in this class. Has no effect if sigmas_q is provided since this value is used solely to calculate sigmas_q.

Shapes: \((N_\mathrm{axes},)\) or \((N_\mathrm{axes},\,N_\mathrm{frames})\).

Reference unit: \(\mathrm{V}\).

thresholdsfloat or array-like, keyword-only, default: 1e-5

Thresholds for determining the plateau regions of the integrals of the charge density profiles to calculate sigmas_q. Has no effect if sigmas_q is provided, or if sigmas_q can be calculated using dVs and charge_density_profiles.

V0sfloat, array-like, openmm.unit.Quantity, or pint.Quantity, keyword-only, default: 0

Potentials \(\Psi_0\) at the left boundary.

Shapes: \((N_\mathrm{axes},)\) or \((N_\mathrm{axes},\,N_\mathrm{frames})\).

Reference unit: \(\mathrm{V}\).

methodsstr or array-like, keyword-only, default: "integral"

Methods to use to calculate the potential profiles.

Valid values: "integral", "matrix".

pbcsbool, keyword-only, default: False

Specifies whether the system has periodic boundary conditions in each of the axes. Only used when method="matrix".

Calculates the surface charge densities \(\sigma_q\) for the specified system dimensions using the charge density profiles \(\rho_q(z)\).

Parameters:

axesstr or array-like, optional

Axes along which to compute the potential profiles. If not specified, all axes for which charge density profiles were calculated will be used.

Examples:"xy" or ("x", "y") for the \(x\)- and \(y\)-directions.

dielectricsfloat, optional

Relative permittivities or dielectric constants \(\varepsilon_\mathrm{r}\). Only optional if previously provided to another calculation method in this class.

dVsfloat, array-like, openmm.unit.Quantity, or pint.Quantity, keyword-only, optional

Potential differences \(\Delta\Psi\) across the system dimensions specified in axes. Can be retrieved if previously provided to another calculation method in this class.

Shapes: \((N_\mathrm{axes},)\) or \((N_\mathrm{axes},\,N_\mathrm{frames})\).

Reference unit: \(\mathrm{V}\).

classmethod get_supported_backends()¶

Tuple with backends supported by the core library for a given class. User can pass either one of these values as backend=... to run() method, or a custom object that has apply method (see documentation for run()):

‘serial’: no parallelization

‘multiprocessing’: parallelization using multiprocessing.Pool

‘dask’: parallelization using dask.delayed.compute(). Requires installation of mdanalysis[dask]

If you want to add your own backend to an existing class, pass a backends.BackendBase subclass (see its documentation to learn how to implement it properly), and specify unsupported_backend=True.

Returns:

tuple: names of built-in backends that can be used in run(backend=...)()

Added in version 2.8.0: ..

property parallelizable¶

Boolean mark showing that a given class can be parallelizable with split-apply-combine procedure. Namely, if we can safely distribute _single_frame() to multiple workers and then combine them with a proper _conclude() call. If set to False, no backends except for serial are supported.

Note

If you want to check parallelizability of the whole class, without explicitly creating an instance of the class, see _analysis_algorithm_is_parallelizable. Note that you setting it to other value will break things if the algorithm behind the analysis is not trivially parallelizable.

Returns:

bool: if a given AnalysisBase subclass instance is parallelizable with split-apply-combine, or not

Added in version 2.8.0: ..

run(start: int = None, stop: int = None, step: int = None, frames: slice | ndarray[int] = None, verbose: bool = None, **kwargs) → SerialAnalysisBase | ParallelAnalysisBase¶

Performs the calculation.