differential_geometry.dense_ops.dense_gradient_contravariant_diag#

differential_geometry.dense_ops.dense_gradient_contravariant_diag(tensor_field: ndarray, inverse_metric_field: ndarray, rank: int, ndim: int, *varargs, field_axes: Sequence[int] | None = None, output_indices: Sequence[int] | None = None, out: ndarray | None = None, edge_order: Literal[1, 2] = 2, **_) → ndarray[source]#

Compute the contravariant gradient \(\nabla^\mu T^{\dots}\) of a tensor field using a diagonal inverse metric tensor.

This function computes the covariant partial derivatives over the field dimensions, then raises the newly introduced index using the provided diagonal inverse metric:

\[\nabla^\mu T^{\dots} = g^{\mu\mu} \partial_\mu T^{\dots}\]

No summation is implied since the metric is diagonal.

Warning

This is a low-level routine and does not validate input shapes or consistency. It assumes all inputs are correctly broadcast and aligned.

Parameters:

  • tensor_field (numpy.ndarray) –

    Tensor field of shape (F_1, ..., F_m, N, ...), where the last rank axes are the tensor index dimensions. The partial derivatives are computed over the first m axes.

    Hint

    Because this function is a low-level callable, it does not enforce the density of the elements in the trailing dimensions of the field. This can be used in some cases where the field is not technically a tensor field.

  • inverse_metric_field (numpy.ndarray) –

    Inverse metric tensor with shape (..., ndim, ), where the leading dimensions (denoted here as F1, ..., F_n) must be broadcast-compatible with the grid (spatial) portion of tensor_field.

    Specifically, if tensor_field has shape (S1, ..., S_m, I1, ..., I_rank), then inverse_metric_field must be broadcast-compatible with (S1, ..., S_m). The last dimension must be exactly (ndim,), representing the inverse metric at each grid point.

  • rank (int) – Number of trailing axes that represent tensor indices (i.e., tensor rank). The rank determines the number of identified coordinate axes and therefore determines the shape of the returned array.

  • ndim (int) – The number of total dimensions in the relevant coordinate system. This determines the maximum allowed value for m and the number of elements in the trailing dimension of the output.

  • *varargs –

    Grid spacing for each spatial axis. Follows the same format as numpy.gradient(), and can be:

    • A single scalar (applied to all spatial axes),

    • A list of scalars (one per axis),

    • A list of coordinate arrays (one per axis),

    • A mix of scalars and arrays (broadcast-compatible).

    If field_axes is provided, varargs must match its length. If axes is not provided, then varargs should be tensor_field.ndim - rank in length (or be a scalar).

  • field_axes (list of int, optional) – The spatial axes over which to compute the component-wise partial derivatives. If field_axes is not specified, then all tensor_field.ndim - rank axes are computed.

  • output_indices (list of int, optional) –

    Explicit mapping from each axis listed in field_axes to the indices in the output’s final dimension of size ndim. This parameter allows precise control over the placement of computed gradients within the output array.

    If provided, output_indices must have the same length as field_axes. Each computed gradient along the spatial axis field_axes[i] will be stored at index output_indices[i] in the last dimension of the output.

    If omitted, the default behavior is output_indices = field_axes, i.e., gradients are placed in the output slot corresponding to their source axis.

    Note

    All positions in the final axis of the output array not listed in output_indices are left unmodified (typically zero-filled unless out was pre-populated).

  • edge_order ({1, 2}, optional) – Gradient is calculated using N-th order accurate differences at the boundaries. Default: 1.

  • out (numpy.ndarray, optional) – Buffer in which to store the output to preserve memory. If provided, out must have shape tensor_field.shape + (ndim,). If out is not specified, then it is allocated during the function call.

Returns:

The contravariant gradient of the input tensor field, with shape B + (ndim,), where B is the broadcasted shape of the spatial portion of tensor_field and inverse_metric_field.

This result represents the gradient with the differentiation index raised by the inverse metric, and the final axis of size ndim corresponds to the contravariant coordinate direction of the derivative at each point.

Return type:

numpy.ndarray

See also

numpy.gradient

The computational backend for this operation.

dense_gradient_covariant

Contravariant gradient for full metric.

dense_gradient_contravariant_full

Contravariant gradient for full metric.

dense_gradient

Wrapper for user-facing gradient computations.

Examples

As a more advanced example, let’s look at the co. versus contra. components of a gradient case:

>>> # ------------------------------------------------------------------
>>> # 0)  Set up a toy scalar field  f(r) = sin(r/r_0) cos^2(theta)
>>> # ------------------------------------------------------------------
>>> from scipy.interpolate import RegularGridInterpolator
>>> from pymetric.differential_geometry.dense_ops import (dense_gradient_contravariant_diag,
... dense_gradient_covariant)
>>> import matplotlib.pyplot as plt
>>> r = np.linspace(1e-4,1,1000)
>>> theta = np.linspace(0,np.pi,30)
>>> R,THETA = np.meshgrid(r,theta,indexing='ij')
>>> Z = np.sin(10*R) * np.cos(THETA)**2
>>>
>>> # Build a 1D interpolator for Z.
>>> interpZ = RegularGridInterpolator((r,theta),Z,bounds_error=False)
>>>
>>> # Compute the gradients of Z
>>> gradZ = dense_gradient_covariant(Z,0,2,r,theta)
>>>
>>> # Create the metric
>>> metric = np.zeros(R.shape + (2,))
>>> metric[:,:,0,] = 1
>>> metric[:,:,1,] = 1/R**2
>>>
>>> # Compute the contravariant gradient.
>>> gradZcontra = dense_gradient_contravariant_diag(Z,metric,0,2,r,theta)
>>>
>>> # Build interpolators for the 2 gradient components.
>>> interpZr = RegularGridInterpolator((r,theta),gradZ[...,0],bounds_error=False)
>>> interpZtheta = RegularGridInterpolator((r,theta),gradZ[...,1],bounds_error=False)
>>> interpCZr = RegularGridInterpolator((r,theta),gradZcontra[...,0],bounds_error=False)
>>> interpCZtheta = RegularGridInterpolator((r,theta),gradZcontra[...,1],bounds_error=False)
>>>
>>> # Construct an X/Y grid.
>>> bound = 1/np.sqrt(2)
>>> x,y = np.linspace(-bound,bound,100),np.linspace(-bound,bound,100)
>>> X,Y = np.meshgrid(x,y,indexing='ij')
>>> RG = np.sqrt(X**2+Y**2)
>>> THETAG = np.arccos(Y/RG)
>>> grid_points = np.stack([RG.ravel(),THETAG.ravel()],axis=1)
>>> Zgrid = interpZ(grid_points).reshape(RG.shape)
>>> Zrgrid = interpZr(grid_points).reshape(RG.shape)
>>> Zthetagrid = interpZtheta(grid_points).reshape(RG.shape)
>>> ZCrgrid = interpCZr(grid_points).reshape(RG.shape)
>>> ZCthetagrid = interpCZtheta(grid_points).reshape(RG.shape)

>>> # Setup the figure.
>>> fig,axes = plt.subplots(2,2, sharex=True, sharey=True)
>>> _ = axes[0,0].imshow(Zrgrid.T    ,extent=[-bound,bound,-bound,bound], vmin=-3,vmax=3,cmap='seismic',origin='lower')
>>> _ = axes[0,1].imshow(Zthetagrid.T,extent=[-bound,bound,-bound,bound], vmin=-3,vmax=3,cmap='seismic',origin='lower')
>>> _ = axes[1,0].imshow(ZCrgrid.T    ,extent=[-bound,bound,-bound,bound],vmin=-3,vmax=3,cmap='seismic',origin='lower')
>>> P = axes[1,1].imshow(ZCthetagrid.T,extent=[-bound,bound,-bound,bound],vmin=-3,vmax=3,cmap='seismic',origin='lower')
>>> _ = plt.colorbar(P,ax=axes)
>>> plt.show()

(Source code, png, hires.png, pdf)

../_images/differential_geometry-dense_ops-dense_gradient_contravariant_diag-1.png