grid.atomgrid module

Module for generating AtomGrid.

class AtomGrid(rgrid, *, degrees=None, sizes=None, center=None, rotate=0, use_spherical=False)

Bases: Grid

Atomic grid construction class for integrating three-dimensional functions.

Atomic grid is composed of a radial grid \(\{(r_i, w_i)\}_{i=1}^{N}\) meant to integrate the radius coordinate in spherical coordinates. Further, each radial point is associated with an Angular (Lebedev or Symmetric spherical t-design) grid \(\{(\theta^i_j, \phi^i_j, w_j^i)\}_{j=1}^{M_i}\) that integrates over a sphere (angles in spherical coordinates). The atomic grid points can also be centered at a given location.

__init__(rgrid, *, degrees=None, sizes=None, center=None, rotate=0, use_spherical=False)

Construct atomic grid for given arguments.

Parameters
  • rgrid (OneDGrid) – The (one-dimensional) radial grid representing the radius of spherical grids.

  • degrees (ndarray(N, dtype=int) or list, keyword-only argument) – Sequence of angular grid degrees used for constructing spherical grids at each radial grid point. If only one degree is given, the specified degree is used for all spherical grids. If the given degree is not supported, the next largest degree is used.

  • sizes (ndarray(N, dtype=int) or list, keyword-only argument) – Sequence of angular grid sizes used for constructing spherical grids at each radial grid point. If only one size is given, the specified size is used for all spherical grids. If the given size is not supported, the next largest size is used. If both degrees and sizes are given, degrees is used for making the spherical grids.

  • center (ndarray(3,), optional, keyword-only argument) – Cartesian coordinates of the grid center. If None, the origin is used.

  • rotate (int, optional) – Integer used as a seed for generating random rotation matrices to rotate the angular spherical grids at each radial grid point. If the integer is zero, then no rotate is used.

  • use_spherical (bool, optional) – If true, loads the symmetric spherical t-design grid rather than the Lebedev-Laikov grid for the angular grid.

property basis

Generate spherical harmonic basis evaluated on atomic grid points.

Type

ndarray(N, 3)

property center

Cartesian coordinates of the grid center.

Type

ndarray(3,)

convert_cartesian_to_spherical(points=None, center=None)

Convert a set of points from Cartesian to spherical coordinates.

The conversion is defined as

\[\begin{split}\begin{align} r &= \sqrt{x^2 + y^2 + z^2}\\ \theta &= arc\tan (\frac{y}{x})\\ \phi &= arc\cos(\frac{z}{r}) \end{align}\end{split}\]

with the canonical choice \(r=0\), then \(\theta,\phi = 0\). If the points attribute is not specified, then atomic grid points are used and the canonical choice when \(r=0\), is the points \((r=0, \theta_j, \phi_j)\) where \((\theta_j, \phi_j)\) come from the Angular grid with the degree at \(r=0\).

Parameters
  • points (ndarray(n, 3), optional) – Points in three-dimensions. Atomic grid points will be used if points is not given

  • center (ndarray(3,), optional) – Center of the atomic grid points. If center is not provided, then the atomic center of this class is used.

Returns

Spherical coordinates of atoms respect to the center (radius \(r\), azimuthal \(\theta\), polar \(\phi\)).

Return type

ndarray(N, 3)

property degrees

Return the degree of each angular grid at each radial point.

Type

ndarray(N,)

classmethod from_preset(rgrid=None, *, atnum, preset, center=None, rotate=0, use_spherical=False)

High level api to construct an atomic grid with preset arguments.

Examples

>>> # construct an atomic grid for H with fine grid setting
>>> atgrid = AtomGrid.from_preset(rgrid, atnum=1, preset="fine")
Parameters
  • rgrid (OneDGrid, optional) – The (1-dimensional) radial grid representing the radius of spherical grids.

  • atnum (int, keyword-only argument) – The atomic number specifying the predefined grid.

  • preset (str, keyword-only argument) – The name of predefined grid specifying the radial sectors and their corresponding number of angular grid points. Supported preset options include: ‘coarse’, ‘medium’, ‘fine’, ‘veryfine’, ‘ultrafine’, and ‘insane’.

  • center (ndarray(3,), optional, keyword-only argument) – Cartesian coordinates of the grid center. If None, the origin is used.

  • rotate (int, optional) – Integer used as a seed for generating random rotation matrices to rotate the angular spherical grids at each radial grid point. If the integer is zero, then no rotate is used.

  • use_spherical (bool, optional) – If true, loads the symmetric spherical t-design grid rather than the Lebedev-Laikov grid for the angular grid.

classmethod from_pruned(rgrid, radius, *_, sectors_r, sectors_degree=None, sectors_size=None, center=None, rotate=0, use_spherical=False)

Initialize AtomGrid class that splits radial sections into sectors which specified degrees.

Given a sequence of radial sectors \(\{a_i\}_{i=1}^Q\), a radius number \(R\) and angular degree sectors \(\{L_i \}_{i=1}^{Q+1}\). This assigned the degrees to the following radial points:

\[\begin{split}\begin{align*} &L_1 \text{ when } r < R a_1 \\ &L_2 \text{ when } R a_1 \leq r < R a_2 \vdots \\ &L_{Q+1} \text{ when } R a_{Q} < r. \end{align*}\end{split}\]

Examples

>>> sectors_r = [0.5, 1., 1.5]
>>> sectors_degree = [3, 7, 5, 3]
# 0 <= r < 0.5 radius, angular grid with degree 3
# 0.5 radius <= r < radius, angular grid with degree 7
# rad <= r < 1.5 radius, angular grid with degree 5
# 1.5 radius <= r, angular grid with degree 3
>>> atgrid = AtomGrid.from_pruned(rgrid, radius, sectors_r, sectors_degree)
Parameters
  • rgrid (OneDGrid) – The (one-dimensional) radial grid representing the radius of spherical grids.

  • radius (float) – The atomic radius to be multiplied with r_sectors (to make them atom specific).

  • sectors_r (ndarray(N,), keyword-only argument) – Sequence of boundary points specifying radial sectors of the pruned grid. The first sector is [0, radius*sectors_r[0]], then [radius*sectors_r[0], radius*sectors_r[1]], and so on.

  • sectors_degree (ndarray(N + 1, dtype=int), keyword-only argument) – Sequence of angular degrees for each radial sector of the pruned grid.

  • sectors_size (ndarray(N + 1, dtype=int), keyword-only argument) – Sequence of angular sizes for each radial sector of the pruned grid. If both sectors_degree and sectors_size are given, sectors_degree is used.

  • center (ndarray(3,), optional, keyword-only argument) – Cartesian coordinates of the grid center. If None, the origin is used.

  • rotate (int, optional) – Integer used as a seed for generating random rotation matrices to rotate the angular spherical grids at each radial grid point. If the integer is zero, then no rotate is used.

  • use_spherical (bool, optional) – If true, loads the symmetric spherical t-design grid rather than the Lebedev-Laikov grid for the angular grid.

Returns

Generated AtomGrid instance for this special init method.

Return type

AtomGrid

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

get_shell_grid(index, r_sq=True)

Get the spherical integral grid at radial point from specified index.

The spherical integration grid has points scaled with the ith radial point and weights multipled by the ith weight of the radial grid.

Note that when \(r=0\) then the Cartesian points are all zeros.

Parameters
  • index (int) – Index of radial points.

  • r_sq (bool, default True) – If true, multiplies the angular grid weights with r**2.

Returns

AngularGrid at given radial index position and weights.

Return type

AngularGrid

property indices

Indices saved for each spherical shell.

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

integrate_angular_coordinates(func_vals)

Integrate the angular coordinates of a sequence of functions.

Given a series of functions \(f_k \in L^2(\mathbb{R}^3)\), this returns the values

\[f_k(r_i) = \int \int f(r_i, \theta, \phi) sin(\theta) d\theta d\phi\]

on each radial point \(r_i\) in the atomic grid.

Parameters

func_vals (ndarray(..., N)) – The function values evaluated on all \(N\) points on the atomic grid for many types of functions. This can also be one-dimensional.

Returns

The function \(f_{...}(r_i)\) on each \(M\) radial points.

Return type

ndarray(…, M)

interpolate(func_vals)

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

Any real-valued function \(f(r, \theta, \phi)\) can be decomposed as

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

A cubic spline is used to interpolate the radial functions \(\sum_i \rho_{ilm}(r)\). This is then multipled by the corresponding spherical harmonics at all \((\theta_j, \phi_j)\) angles and summed to obtain the equation above.

Parameters

func_vals (ndarray(N,)) – The function values evaluated on all \(N\) points on the 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

>>> # First generate a atomic grid with raidal points that have all degree 10.
>>> from grid.basegrid import OneDGrid
>>> radial_grid = OneDGrid(np.linspace(0.01, 10, num=100), np.ones(100), (0, np.inf))
>>> atom_grid = AtomGrid(radial_grid, degrees=[10])
# 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 function values and interpolate them
>>> func_vals = polynomial_func(atom_grid.points)
>>> interpolate_func = atom_grid.interpolate(func_vals)
# To interpolate at new points.
>>> new_pts = np.array([[1.0, 1.0, 1.0], [0.0, 0.0, 0.0]])
>>> 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)
property l_max

Largest angular degree L value in angular grids.

Type

int

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 n_shells

Number of shells in radial points.

Type

int

property points

Cartesian coordinates of the grid points (centered).

Type

ndarray(N, 3)

radial_component_splines(func_vals)

Return spline to interpolate radial components wrt to expansion in real spherical harmonics.

For each pt \(r_i\) of the atomic grid with associated angular degree \(l_i\), the function \(f(r_i, \theta, \phi)\) is projected onto the spherical harmonic expansion:

\[f(r_i, \theta, \phi) \approx \sum_{l=0}^{l_i} \sum_{m=-l}^l \rho^{lm}(r_i) Y^m_l(\theta, \phi)\]

where \(Y^m_l\) is the real Spherical harmonic of degree \(l\) and order \(m\). The radial components \(\rho^{lm}(r_i)\) are calculated via integration on the :math:`i`th Lebedev/angular grid of the atomic grid:

\[\rho^{lm}(r_i) = \int \int f(r_i, \theta, \phi) Y^m_l(\theta, \phi) \sin(\theta) d\theta d\phi,\]

and then interpolated using a cubic spline over all radial points of the atomic grid.

Parameters

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

Returns

A list of size \((l_{max}/2 + 1)^2\) of CubicSpline object for interpolating the coefficients \(\rho^{lm}(r)\). The input of spline is array of \(N\) points on \([0, \infty)\) and the output is \(\{\rho^{lm}(r_i)\}\). The list starts with degrees \(l=0,\cdots l_{max}\), and the for each degree, the zeroth order spline is first, followed by positive orders then negative.

Return type

list[scipy.PPoly]

property rgrid

The radial grid representing the radius of spherical grids.

Type

OneDGrid

property rotate

Integer representing the seed for rotating the angular grid.

Type

int

save(filename)

Save atomic 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

spherical_average(func_vals)

Return spline of the spherical average of a function.

This function takes a function \(f\) evaluated on the atomic grid points and returns the spherical average of it defined as:

\[f_{avg}(r) := \frac{\int \int f(r, \theta, \phi) \sin(\theta) d\theta d\phi}{4 \pi}.\]

The definition is chosen such that \(\int f_{avg}(r) 4\pi r^2 dr\) matches the full integral \(\int \int \int f(x,y,z)dxdydz\).

Parameters

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

Returns

Cubic spline with input r in the positive real axis and output \(f_{avg}(r)\).

Return type

CubicSpline

Examples

>>> # Define a Gaussian function that takes Cartesian coordinates as input
>>> func = lambda cart_pts: np.exp(-np.linalg.norm(cart_pts, axis=1)**2.0)
# Construct atomic grid with degree 10 on a radial grid on [0, \infty)
>>> radial_grid = GaussLaguerre(100, alpha=1.0)
>>> atgrid = AtomGrid(radial_grid, degrees=[10])
# Evaluate func on atmic grid points (which are stored in Cartesian coordinates)
>>> func_vals = func(atgrid.points)
# Compute spherical average spline & evaluate it on a set of (radial) points in [0, \infty)
>>> spherical_avg = atgrid.spherical_average(func_vals)
>>> points = np.arange(0.0, 10.0)
>>> evals = spherical_avg(points)
# the largest error happens at origin because the spline is being extrapolated
>>> assert np.all(abs(evals - np.exp(- points ** 2)) < 1.0e-3)
property use_spherical

True then symmetric spherical t-design is used rather than Lebedev-Laikov grid.

Type

bool

property weights

the weights of each grid point.

Type

np.ndarray(N,)