grid.molgrid module

Molecular grid class.

class MolGrid(atnums, atgrids, aim_weights, store=False)

Bases: Grid

Molecular grid class for integration of three-dimensional functions.

Molecular grid is defined here to be a weighted average of \(M\) atomic grids (see AtomGrid). This is defined by a atom in molecule weights (or nuclear weight functions) \(w_n(r)\) for each center n such that \(\sum^M_n w_n(r) = 1\) for all points \(r\in\mathbb{R}^3.\)

Examples

There are multiple methods of specifiying molecular grids. This example chooses Becke weights as the atom in molecule/nuclear weights and the radial grid is the same for all atoms. Two atoms are considered with charges [1, 2], respectively.

>>> from grid.becke BeckeWeights
>>> from grid.onedgrid import GaussLaguerre
>>> becke = BeckeWeights(order=3)
>>> rgrid = GaussLaguerre(100)
>>> charges = [1, 2]
>>> coords = np.array([[-1.0, 0.0, 0.0], [1.0, 0.0, 0.0]])

The default method is based on explicitly specifing the atomic grids (AtomGrids) for each atom.

>>> from grid.atomgrid import AtomGrid
>>> atgrid1 = AtomGrid(rgrid, degrees=5, center=coords[0])
>>> atgrid2 = AtomGrid(rgrid, degrees=10, center=coords[1])
>>> molgrid = MolGrid(charges, [atgrid1, atgrid2], aim_weights=becke)

The from_size method constructs AtomGrids with degree_size specified from integer size.

>>> size = 100  # Number of angular points used in each shell in the atomic grid.
>>> molgrid = MolGrid.from_size(charges, coords, rgrid, size=5, aim_weights=becke)

The from_pruned method is based on AtomGrid.from_pruned method on the idea of spliting radial grid points into sectors that have the same angular degrees.

>>> sectors_r = [[0.5, 1., 1.5], [0.25, 0.5]]
>>> sectors_deg = [[3, 7, 5, 3], [3, 2, 2]]
>>> radius = 1.0
>>> mol_grid = MolGrid.from_pruned(charges, coords, rgrid, radius, becke,
>>>                                sectors_r=sectors_r, sectors_degree=sectors_deg)

The from_preset method is based on AtomGrid.from_preset method based on a string specifying the size of each Levedev grid at each radial points.

>>> preset = "fine"  # Many choices available.
>>> molgrid = MolGrid.from_preset(charges, coords, rgrid, preset, aim_weights=becke)

The general way to integrate is the following.

>>> integrand = integrand_func(molgrid.points)
>>> integrate = molgrid.integrate(integrand)

References

1

Becke, Axel D. “A multicenter numerical integration scheme for polyatomic molecules.” The Journal of chemical physics 88.4 (1988): 2547-2553.

__init__(atnums, atgrids, aim_weights, store=False)

Initialize class.

Parameters

  • atnums (np.ndarray(M, 3)) – Atomic number of \(M\) atoms in molecule.

  • atgrids (list[AtomGrid]) – List of atomic grids of size \(M\) for each atom in molecule.

  • aim_weights (Callable or np.ndarray(sum^M_n N_n,)) – Atoms in molecule weights \({ {w_n(r_k)}_k^{N_i}}_n^{M}\), where \(N_i\) is the number of points in the ith atomic grid.

  • store (bool) – If true, then the atomic grids atgrids are stored as attribute atgrids.

property aim_weights

Get the atom in molecules/nuclear weight function of the molecular grid.

Returns

List of atomic weights where \(K = \sum_n N_i\) and \(N_i\) is the number of points in the ith atomic grid.

Return type

ndarray(K,)

property atcoords

Center/Atomic coordinates.

Type

ndarray(M, 3)

property atgrids

List of atomic grids for each center. Returns None if store false.

Type

List[AtomGrid]

property atweights

Get the atomic grid weights for all points, i.e. without atom in molecule weights.

Returns

List of weights of atomic grid for each point,where \(K = \sum_n N_i\) and \(N_i\) is the number of points in the ith atomic grid.

Return type

ndarray(K,)

classmethod from_preset(atnums, atcoords, rgrid, preset, aim_weights, *_, rotate=37, store=False)

Construct molecular grid wih preset parameters.

Parameters

  • atnums (np.ndarray(M,)) – Array of atomic numbers.

  • atcoords (np.ndarray(M, 3)) – Atomic coordinates of atoms.

  • rgrid ((OneDGrid, list[OneDGrid], dict[int: OneDGrid])) – One dimensional radial grid. If of type OneDGrid then this radial grid is used for all atoms. If a list is provided,then ith grid correspond to the ith atom. If dictionary is provided, then the keys correspond to the `atnums[i]`attribute.

  • preset ((str, list[str], dict[int: str])) – Preset grid accuracy scheme, support “coarse”, “medium”, “fine”, “veryfine”, “ultrafine”, “insane”. If string is provided ,then preset is used for all atoms, either it is specified by a list, or a dictionary whose keys are from atnums.

  • aim_weights (Callable or np.ndarray(K,)) – Atoms in molecule weights.

  • rotate (bool or int, optional) – Flag to set auto rotation for atomic grid, if given int, the number will be used as a seed to generate rantom matrix.

  • store (bool, optional) – Store atomic grid as a class attribute.

classmethod from_pruned(atnums, atcoords, rgrid, radius, aim_weights, sectors_r, sectors_degree=None, sectors_size=None, rotate=37, store=False)

Initialize a MolGrid instance with pruned method from AtomGrid.

Parameters

  • atnums (ndarray(M,)) – List of atomic numbers for each atom.

  • atcoords (np.ndarray(M, 3)) – Cartesian coordinates for each atoms

  • rgrid (OneDGrid or List[OneDGrid] or Dict[int: OneDGrid]) – One dimensional grid for the radial component. If a list is provided,then ith grid correspond to the ith atom. If dictionary is provided, then the keys are correspond to the atnums[i] attribute.

  • radius (float, List[float]) – The atomic radius to be multiplied with r_sectors (to make them atom specific). If float, then the same atomic radius is used for all atoms, else a list specifies it for each atom.

  • aim_weights (Callable or np.ndarray(sum^M_n N_n,)) – Atoms in molecule/nuclear weights \({ {w_n(r_k)}_k^{N_i}}_n^{M}\), where \(N_i\) is the number of points in the ith atomic grid.

  • sectors_r (List[List], keyword-only argument) – Each row is a sequence of boundary points specifying radial sectors of the pruned grid for the m`th atom. The first sector is ``[0, radius*sectors_r[0]]`, then [radius*sectors_r[0], radius*sectors_r[1]], and so on.

  • sectors_degree (List[List], keyword-only argument) – Each row is a sequence of Lebedev/angular degrees for each radial sector of the pruned grid for the m`th atom. If both `sectors_degree and sectors_size are given, sectors_degree is used.

  • sectors_size (List[List], keyword-only argument) – Each row is a sequence of Lebedev sizes for each radial sector of the pruned grid for the m`th atom. If both `sectors_degree and sectors_size are given, sectors_degree is used.

  • rotate (bool or int , optional) – Flag to set auto rotation for atomic grid, if given int, the number will be used as a seed to generate rantom matrix.

  • store (bool, optional) – Flag to store each original atomic grid information.

Returns

Molecular grid class for integration.

Return type

MolGrid

classmethod from_size(atnums, atcoords, rgrid, size, aim_weights, rotate=37, store=False)

Initialize a MolGrid instance with Horton Style input.

Example

>>> onedg = UniformInteger(100) # number of points, oned grid before TF.
>>> rgrid = ExpRTransform(1e-5, 2e1).generate_radial(onedg) # radial grid
>>> becke = BeckeWeights(order=3)
>>> molgrid = MolGrid.from_size(atnums, atcoords, rgrid, 110, becke)
Parameters

  • atnums (np.ndarray(M, 3)) – Atomic number of \(M\) atoms in molecule.

  • atcoords (np.ndarray(N, 3)) – Cartesian coordinates for each atoms

  • rgrid (OneDGrid) – one dimension grid to construct spherical grid

  • size (int) – Num of points on each shell of angular grid

  • aim_weights (Callable or np.ndarray(K,)) – Atoms in molecule weights.

  • rotate (bool or int , optional) – Flag to set auto rotation for atomic grid, if given int, the number will be used as a seed to generate rantom matrix.

  • store (bool, optional) – Flag to store each original atomic grid information

Returns

MolGrid instance with specified grid property

Return type

MolGrid

get_atomic_grid(index)

Get atomic grid corresponding to the given atomic index.

Parameters

index (int) – index of atom starting from 0 to \(N-1\) where \(N\) is the number of atoms in molecule.

Returns

If store=True, the AtomGrid instance used is returned. If store=False, the LocalGrid containing points and weights of AtomGrid is returned.

Return type

AtomGrid or LocalGrid

get_localgrid(center, radius)

Create a grid contain points within the given radius of center.

Parameters

  • center (float or np.array(M,)) – Cartesian coordinates of the center of the local grid.

  • radius (float) – Radius of sphere around the center. When equal to np.inf, the local grid coincides with the whole grid, which can be useful for debugging.

Returns

Instance of LocalGrid.

Return type

LocalGrid

property indices

Get the indices of the molecualr grid.

Returns

List of indices \([i_0, i_1, \cdots]\) that whose indices range [i_k, i_{k+1}] specify which points in points correspond to kth atomic grid.

Return type

ndarray(M + 1,)

integrate(*value_arrays)

Integrate over the whole grid for given multiple value arrays.

Product of all value_arrays will be computed element-wise then integrated on the grid with its weights:

\[\int w(x) \prod_i f_i(x) dx.\]
Parameters

*value_arrays (np.ndarray(N, )) – One or multiple value array to integrate.

Returns

The calculated integral over given integrand or function

Return type

float

interpolate(func_vals)

Return function that interpolates (and its derivatives) from function values.

Consider a real-valued function \(f(r, \theta, \phi) = \sum_A f_A(r, \theta, \phi)\) written as a sum of atomic centers \(f_A(r, \theta, \phi)\). Each of these functions can be further decomposed based on the atom grid centered at \(A\):

\[f_A(r, \theta, \phi) = \sum_l \sum_{m=-l}^l \sum_i \rho^A_{ilm}(r) Y_{lm}(\theta, \phi)\]

A cubic spline is used to interpolate the radial functions \(\sum_i \rho^A_{ilm}(r)\) based on the atomic grid centered at \(A\). This is then multipled by the corresponding spherical harmonics at all \((\theta_j, \phi_j)\) angles and summed to obtain approximation to \(f_A\). This is then further summed over all centers to get \(f\).

Parameters

func_vals (ndarray(sum_i N_i,)) – The function values evaluated on all \(N_i\) points on the :math:`i`th atomic grid.

Returns

Callable function that interpolates the function and its derivative provided. The function takes the following attributes:

pointsndarray(N, 3)

Cartesian coordinates of \(N\) points to evaluate the splines on.

derivint, optional

If deriv is zero, then only returns function values. If it is one, then returns the first derivative of the interpolated function with respect to either Cartesian or spherical coordinates. Only higher-order derivatives (`deriv`=2,3) are supported for the derivatives wrt to radial components.

deriv_sphericalbool

If True, then returns the derivatives with respect to spherical coordinates \((r, \theta, \phi)\). Default False.

only_radial_derivbool

If true, then the derivative wrt to radius \(r\) is returned.

This function returns the following.

ndarray(M,…):

The interpolated function values or its derivatives with respect to Cartesian \((x,y,z)\) or if deriv_spherical then \((r, \theta, \phi)\) or if only_radial_derivs then derivative wrt to \(r\) is only returned.

Return type

Callable[[ndarray(M, 3), int] -> ndarray(M)]

Examples

>>> # Consider the function (3x^2 + 4y^2 + 5z^2)
>>> def polynomial_func(pts) :
>>>     return 3.0 * points[:, 0]**2.0 + 4.0 * points[:, 1]**2.0 + 5.0 * points[:, 2]**2.0
>>> # Evaluate the polynomial over the molecular grid points and interpolate
>>> polynomial_vals = polynomial_func(molgrid.points)
>>> interpolate_func = molgrid.interpolate(polynomial_vals)
>>> # Use it to interpolate at new points.
>>> interpolate_vals = interpolate_func(new_pts)
>>> # Can calculate first derivative wrt to Cartesian or spherical
>>> interpolate_derivs = interpolate_func(new_pts, deriv=1)
>>> interpolate_derivs_sph = interpolate_func(new_pts, deriv=1, deriv_spherical=True)
>>> # Only higher-order derivatives are supported for the radius coordinate r.
>>> interpolated_derivs_radial = interpolate_func(new_pts, deriv=2, only_radial_derivs=True)
moments(orders, centers, func_vals, type_mom='cartesian', return_orders=False)

Compute the multipole moment integral of a function over centers.

The Cartesian type moments are:

\[m_{n_x, n_y, n_z} = \int (x - X_c)^{n_x} (y - Y_c)^{n_y} (z - Z_c)^{n_z} f(r) dr,\]

where \(\textbf{R}_c = (X_c, Y_c, Z_c)\) is the center of the moment, \(\f(r)\) is the density, and \((n_x, n_y, n_z)\) are the Cartesian orders.

The spherical/pure moments with \((l, m)\) parameter are:

\[m_{lm} = \int | \textbf{r} - \textbf{R}_c|^l S_l^m(\theta, \phi) f(\textbf{r}) d\textbf{r},\]

where \(S_l^m\) is a regular, real solid harmonic.

The radial moments with \(n\) parameter are:

\[m_n = \int | \textbf{r} - \textbf{R}_c|^{n} f(\textbf{r}) d\textbf{r}\]

The radial combined with spherical/pure moments \((n, l, m)\) are:

\[m_{nlm} = \int | \textbf{r} - \textbf{R}_c|^{n+1} S_l^m(\theta, \phi) f(\textbf{r}) d\textbf{r}\]
Parameters

  • orders (int) – Generates all orders with Horton order depending on the type of the multipole moment type_mom.

  • centers (ndarray(M, 3)) – The centers \(\textbf{R}_c\) of the moments to compute from.

  • func_vals (ndarray(N,)) – The function \(f\) values evaluated on all \(N\) points on the integration grid.

  • type_mom (str) – The type of multipole moments: “cartesian”, “pure”, “radial” and “pure-radial”.

  • return_orders (bool) – If true, it will also return a list of size \(L\) of the orders corresponding to each integral/row of the output.

Returns

Computes the moment integral of the function on the m`th center for all orders. If `return_orders is true, then this also returns a list that describes what each row/order is, e.g. for Cartesian, [(0, 0, 0), (1, 0, 0) ,…].

Return type

ndarray(L, M), or (ndarray(L, M), list)

property points

Positions of the grid points.

Type

np.ndarray(N,) or np.ndarray(N, M)

save(filename)

Save molecular grid attributes as a npz file.

Parameters

filename (str) – The path/name of the .npz file.

property size

the total number of points on the grid.

Type

int

property weights

the weights of each grid point.

Type

np.ndarray(N,)