differential_geometry.dense_ops.dense_vector_divergence_contravariant#

differential_geometry.dense_ops.dense_vector_divergence_contravariant(vector_field: ndarray, Dterm_field: ndarray, *varargs, derivative_field: ndarray | None = None, field_axes: Sequence[int] | None = None, derivative_axes: Sequence[int] | None = None, edge_order: Literal[1, 2] = 2, out: ndarray | None = None, **_) → ndarray[source]#

Compute the divergence of a contravariant vector field in a general coordinate system using the D-terms and raw partial derivatives.

This implements the formula:

\[\nabla_i V^i = D_i V^i + \partial_i V^i\]

where \(D_i = (\partial_i \rho) / \rho\) is the logarithmic derivative of the coordinate density.

Parameters:

  • vector_field (numpy.ndarray) – Contravariant vector field with shape (F_1, ..., F_M, ndim), where the final axis corresponds to coordinate directions. The first m axes are spatial grid axes. Must be broadcast-compatible with Dterm_field.

  • Dterm_field (numpy.ndarray) – D-term array of shape (..., ndim), where the last axis matches the number of coordinate directions. Must broadcast with the spatial axes of vector_field.

  • *varargs –

    Grid spacing for each 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 derivative_axes is provided, then varargs must match its shape. Otherwise, there must be m elements in varargs.

  • derivative_field (numpy.ndarray, optional) – Optional array of precomputed partial derivatives of selected vector components. Must have shape broadcast-compatible with the spatial shape of vector_field, and a final axis indexing the selected derivative components (same length as derivative_axes).

  • field_axes (list of int, optional) – Maps each grid axis (0 to m-1) to a corresponding component index in the vector field. Defaults to identity mapping [0, 1, ..., m-1].

  • derivative_axes (list of int, optional) – Grid axes along which to compute partial derivatives. If not specified, all spatial axes listed in field_axes are used.

  • edge_order ({1, 2}, default 2) – Accuracy order of finite differences used in derivative computation.

  • out (numpy.ndarray, optional) – Optional output buffer for storing the result. Must have shape equal to the broadcast of the grid (non-component) dimensions of vector_field and Dterm_field.

Returns:

A scalar field representing the divergence, with shape equal to the broadcasted grid shape of vector_field and Dterm_field (excluding the final component axis).

Return type:

numpy.ndarray

Raises:

ValueError – If axes length does not match the spatial rank of the input.

Examples

In the most basic case, we only need to supply the D-field and the vector field.

>>> import numpy as np
>>>
>>> # Build a field.
>>> x,y = np.linspace(-1,1,5),np.linspace(-1,1,5)
>>> X,Y = np.meshgrid(x,y,indexing='ij')
>>> V = np.stack([np.ones_like(X),Y],axis=-1)
>>>
>>> # Create cartesian Dfield.
>>> D = np.zeros_like(V)
>>>
>>> # Compute the divergence.
>>> DivV = dense_vector_divergence_contravariant(V,D,x,y)
>>> np.all(DivV == 1.0)
True

In some cases, the D term and the vector field might only be broadcastable. This is perfectly fine, but the fields must be broadcastable already! For example, this will fail:

>>> import numpy as np
>>>
>>> # Build a field.
>>> x,y = np.linspace(-1,1,6),np.linspace(-1,1,5)
>>> X,Y = np.meshgrid(x,y,indexing='ij')
>>> V = np.stack([np.ones_like(y),y],axis=-1)
>>> V.shape
(5, 2)
>>>
>>> # Create cartesian Dfield.
>>> D = np.zeros(x.shape + (2,))
>>>
>>> # Compute the divergence.
>>> DivV = dense_vector_divergence_contravariant(V,D,x,y) 
ValueError: shape mismatch: objects cannot be broadcast to a single shape.  Mismatch is between arg 0 with shape (5,) and arg 1 with shape (6,).

This will succeed:

>>> import numpy as np
>>>
>>> # Build a field.
>>> x,y = np.linspace(-1,1,6),np.linspace(-1,1,5)
>>> X,Y = np.meshgrid(x,y,indexing='ij')
>>> V = np.stack([np.ones_like(y),y],axis=-1)[None,:,:]
>>> V.shape
(1, 5, 2)
>>>
>>> # Create cartesian Dfield.
>>> D = np.zeros(x.shape + (2,))[:,None,:]
>>> D.shape
(6, 1, 2)
>>>
>>> # Compute the divergence.
>>> DivV = dense_vector_divergence_contravariant(V,D,y,derivative_axes=[1])
>>> np.all(DivV == 1.0)
True

The derivative_axes=[1] tells dense_vector_divergence_contravariant() not to compute a derivative over the x axis (its singleton).