grid.utils module¶
Utils function module.
- convert_cart_to_sph(points, center=None)¶
Convert a set of points from cartesian to spherical coordinates.
- Parameters
points (np.ndarray(n, 3)) – The (3-dimensional) Cartesian coordinates of points.
center (np.ndarray(3,), list, optional) – Cartesian coordinates of the center of spherical coordinate system. If None, origin is used.
- Returns
Spherical coordinates of atoms respect to the center [radius \(r\), azumuthal \(\theta\), polar \(\phi\)]
- Return type
np.ndarray(N, 3)
- convert_derivative_from_spherical_to_cartesian(deriv_r, deriv_theta, deriv_phi, r, theta, phi)¶
Convert the derivative form spherical coordinates to Cartesian.
If \(r\) is zero, then the derivative wrt to theta and phi are set to zero. If \(\phi\) is zero, then the derivative wrt to theta is set to zero.
- Parameters
deriv_r (float, float, float) – Derivative wrt to spherical coordinates
deriv_theta (float, float, float) – Derivative wrt to spherical coordinates
deriv_phi (float, float, float) – Derivative wrt to spherical coordinates
r (float, float, float) – The spherical points where the derivative is taken.
theta (float, float, float) – The spherical points where the derivative is taken.
phi (float, float, float) – The spherical points where the derivative is taken.
- Returns
The derivative wrt to (x, y, z) Cartesian coordinates.
- Return type
ndarray(3,)
- dipole_moment_of_molecule(grid, density, coords, charges)¶
Calculate the dipole of a molecule.
This is defined as the observable form a wavefunction \(\Psi\):
\[\vec{\mu} = \int \Psi \hat{\mu} \Psi \vec{r}\]which results in
\[\vec{\mu} = \sum_{i=1}^{N_{atoms}} Z_i (\vec{R_i} - \vec{R_c}) - \int ((\vec{r} - \vec{R_c})) \rho(\vec{r}) dr,\]where \(N_{atoms}\) is the number of atoms, \(Z_i\) is the atomic charge of the ith atom, \(\vec{R_i}\) is the ith coordinate of the atom, \(\vec{R_c}\) is the center of the molecule and \(\rho\) is the electron density of the molecule.
- Parameters
grid (Grid) – The grid used to perform the integration.
density (ndarray) – The electron density evaluated on points on grid object.
coords (ndarray[M, 3]) – The coordinates of the atoms.
charges (ndarray[M]) – The atomic charge of each atom.
- Returns
Returns array of size three corresponding to the dipole moment of a molecule.
- Return type
ndarray
- generate_derivative_real_spherical_harmonics(l_max, theta, phi)¶
Generate derivatives of real spherical harmonics.
If \(\phi\) is zero, then the first component of the derivative wrt to \(phi\) is set to zero.
- Parameters
l_max (int) – Largest angular degree of the spherical harmonics.
theta (np.ndarray(N,)) – Azimuthal angle \(\theta \in [0, 2\pi]\) that are being evaluated on. If this angle is outside of bounds, then periodicity is used.
phi (np.ndarray(N,)) – Polar angle \(\phi \in [0, \pi]\) that are being evaluated on. If this angle is outside of bounds, then periodicity is used.
- Returns
Derivative of spherical harmonics, (theta first, then phi) of all degrees up to \(l_{max}\) and orders \(m\) in Horton 2 order, i.e. \(m=0, 1, -1, \cdots, l, -l\).
- Return type
ndarray(2, (l_max^2 + 1)^2, M)
Notes
The derivative with respect to \(\phi\) is
\[\frac{\partial Y_l^m(\theta, \phi)}{\partial\phi} = -m Y_l^{-m}(\theta, \phi),\]The derivative with respect to \(\theta\) is
\[\frac{\partial Y_l^m(\theta, \phi)}{\partial\phi} = |m| \cot(\phi) \Re (Y_l^{|m|}(\theta, \phi)) + \sqrt{(n - |m|)(n + |m| + 1)} Re(e^{-i\theta} Y_l^{|m| + 1}(\theta, \phi))\]for positive \(m\) and for negative \(m\), the real projection \(\Re\) is replaced with the imaginary projection \(\Im\).
- generate_orders_horton_order(order, type_ord, dim=3)¶
Generate all orders from an integer \(l\).
For Cartesian, the orders are \((n_x, n_y, n_z)\) such that they sum to order. If dim=1,2, then it generates Cartesian orders \((n_x)\), \((n_x, n_y)\), respectively, such that they sum to order.
For radial, the orders is just the order \(l\).
For pure, the orders \((l, m)\) following the order \([(l, 0), (l, 1), (l, -1), \cdots, (l, l), (l, -l)]\).
- Parameters
order (int) – The order \(l\).
type_ord (str, optional) – The type of the order, it is either “cartesian”, “radial”, “pure” or “pure-radial”.
dim (int, optional) – The dimension of the orders for only Cartesian.
- Returns
Each row is a list of D integers (e.g. \((n_x, n_y, n_z)\) or \((l, m)\)). The output is in Horton 2 order. e.g. order=2 it’s [[2, 0, 0], [1, 1, 0], [1, 0, 1], [0, 2, 0], ….]
- Return type
ndarray(3 * order , D)
- generate_real_spherical_harmonics(l_max, theta, phi)¶
Compute the real spherical harmonics recursively up to a maximum angular degree l.
\[Y_l^m(\theta, \phi) = \frac{(2l + 1) (l - m)!}{4\pi (l + m)!} f(m, \theta) P_l^m(\cos(\phi)),\]where \(l\) is the angular degree, \(m\) is the order and \(f(m, \theta) = \sqrt{2} \cos(m \theta)\) when \(m>0\) otherwise \(f(m, \theta) = \sqrt{2} \sin(m\theta)\) when \(m<0\), and equal to one when \(m= 0\). \(P_l^m\) is the associated Legendre polynomial without the Conway phase factor. The angle \(\theta \in [0, 2\pi]\) is the azimuthal angle and \(\phi \in [0, \pi]\) is the polar angle.
- Parameters
l_max (int) – Largest angular degree of the spherical harmonics.
theta (np.ndarray(N,)) – Azimuthal angle \(\theta \in [0, 2\pi]\) that are being evaluated on. If this angle is outside of bounds, then periodicity is used.
phi (np.ndarray(N,)) – Polar angle \(\phi \in [0, \pi]\) that are being evaluated on. If this angle is outside of bounds, then periodicity is used.
- Returns
Value of real spherical harmonics of all orders \(m\),and degree \(l\) spherical harmonics. For each degree \(l\), the orders \(m\) are in Horton 2 order, i.e. \(m=0, 1, -1, 2, -2, \cdots, l, -l\).
- Return type
ndarray((l_max + 1)**2, N)
Notes
The associated Legendre polynomials are computed using the forward recursion: \(P_l^m(\phi) = \frac{2l + 1}{l - m + 1}\cos(\phi) P_{l-1}^m(x) - \frac{(l + m)}{l - m + 1} P_{l-1}^m(x)\), and \(P_l^l(\phi) = (2l + 1) \sin(\phi) P_{l-1}^{l-1}(\phi)\).
For higher maximum degree \(l_{max} > 1900\) with double precision the computation of spherical harmonic will underflow within the range \(20^\circ \leq \phi \leq 160^\circ\). This code does not implement the modified recursion formulas in [2] and instead opts for higher precision defined by the computer system and np.longdouble.
The mapping from \((l, m)\) to the correct row index in the spherical harmonic array is given by \((l + 1)^2 + 2 * m - 1\) if \(m > 0\), else it is \((l+1)^2 + 2 |m|\).
References
- 1
Colombo, Oscar L. Numerical methods for harmonic analysis on the sphere. Ohio State Univ Columbus Dept of Geodetic Science And Surveying, 1981.
- 2
Holmes, Simon A., and Will E. Featherstone. “A unified approach to the Clenshaw summation and the recursive computation of very high degree and order normalised associated Legendre functions.” Journal of Geodesy 76.5 (2002): 279-299.
- generate_real_spherical_harmonics_scipy(l_max, theta, phi)¶
Generate real spherical harmonics up to degree \(l\) and for all orders \(m\).
The real spherical harmonics are defined as a function \(Y_{lm} : L^2(S^2) \rightarrow \mathbb{R}\) such that
\[\begin{split}Y_{lm} = \begin{cases} \frac{i}{\sqrt{2}} (Y^m_l - (-1)^m Y_l^{-m} & \text{if } < m < 0 \\ Y_l^0 & \text{if } m = 0 \\ \frac{1}{\sqrt{2}} (Y^{-|m|}_{l} + (-1)^m Y_l^m) & \text{if } m > 0 \end{cases},\end{split}\]where \(l \in \mathbb{N}\), \(m \in \{-l_{max}, \cdots, l_{max} \}\) and \(Y^m_l\) is the complex spherical harmonic.
- Parameters
l_max (int) – Largest angular degree of the spherical harmonics.
theta (np.ndarray(N,)) – Azimuthal angle \(\theta \in [0, 2\pi]\) that are being evaluated on. If this angle is outside of bounds, then periodicity is used.
phi (np.ndarray(N,)) – Polar angle \(\phi \in [0, \pi]\) that are being evaluated on. If this angle is outside of bounds, then periodicity is used.
- Returns
Value of real spherical harmonics of all orders \(m\),and degree \(l\) spherical harmonics. For each degree, the zeroth order is stored, followed by positive orders then negative orders,ordered as: \((-l, -l + 1, \cdots -1)\).
- Return type
ndarray((l_max + 1)**2, N)
Examples
To obtain the l-th degree for all orders >>> spherical_harmonic = generate_real_spherical_harmonics(5, theta, phi) To obtain specific degrees, e.g. l=2 >>> desired_degree = 2 >>> spherical_harmonic[(desired_degree)**2: (desired_degree + 1)**2, :]
- get_cov_radii(atnums, type='bragg')¶
Get the covalent radii for given atomic number(s).
- Parameters
atnums (int or np.ndarray) – atomic number of interested
type (str, default to bragg) – types of covalent radii for elements. “bragg”: Bragg-Slater empirically measured covalent radii “cambridge”: Covalent radii from analysis of the Cambridge Database” “alvarez”: Covalent radii from https://doi.org/10.1039/B801115J
- Returns
covalent radii of desired atom(s)
- Return type
np.ndarray
- Raises
ValueError – Invalid covalent type, or input atomic number is 0
- solid_harmonics(l_max, sph_pts)¶
Generate the solid harmonics from zero to a maximum angular degree.
\[R_l^m(r, \theta, \phi) = \sqrt{\frac{4\pi}{2l + 1}} r^l Y_l^m(\theta, \phi)\]- Parameters
l_max (int) – Largest angular degree of the spherical harmonics.
sph_pts (ndarray(M, 3)) – Three-dimensional points in spherical coordinates \((r, \theta, \phi)\), where \(r\geq 0, \theta \in [0, 2\pi]\) and \(\phi \in [0, \pi]\).
- Returns
The solid harmonics from degree \(l=0=, \cdots, l_{max}\) and for all orders \(m\), ordered as \(m=0, 1, -1, 2, -2, \cdots, l, -l\) evaluated on \(M\) points.
- Return type
ndarray((l_max + 1)^2, M)