quickpaver.make_rlg_spatial_gradient_matrices#

quickpaver.make_rlg_spatial_gradient_matrices(grid: RectilinearGrid, sub_selection: ndarray[tuple[Any, ...], dtype[int64]] | None = None, which: Literal['forward', 'backward', 'both'] = 'both') → Tuple[csc_array, csc_array, csc_array][source]#

Build sparse spatial-gradient matrices for a rectilinear grid.

The returned matrices approximate spatial gradients along the x-, y-, and z-axes of a flattened grid-cell field. Each matrix is built from neighbouring cell pairs along one axis and uses the finite-volume ratio between interface area and cell volume.

If field is a one-dimensional array of shape (grid.n_grid_cells,), then the directional gradient-like quantities can be obtained with:

grad_x = gx @ field

grad_y = gy @ field

grad_z = gz @ field

where gx, gy, and gz are the matrices returned by this function.

Parameters:

  • grid (RectilinearGrid) – Rectilinear grid defining the cell dimensions, shape, indexing convention, and total number of grid cells.

  • sub_selection (Optional[NDArrayInt], optional) – Optional flattened grid-cell indices to include in the gradient computation. Neighbour pairs are kept only when both the owner cell and the neighbouring cell belong to this selection. If None, all grid cells are used. By default None.

  • which (Literal["forward", "backward", "both"], optional) –

    Difference scheme to assemble along each axis:

    • "forward" builds contributions from each cell to its forward neighbour.

    • "backward" builds contributions from each cell to its backward neighbour.

    • "both" combines forward and backward contributions.

    By default "both".

Returns:

Sparse gradient matrices (gx, gy, gz) in CSC format, where:

  • gx acts along the x-axis,

  • gy acts along the y-axis,

  • gz acts along the z-axis.

Each matrix has shape (grid.n_grid_cells, grid.n_grid_cells).

Return type:

Tuple[csc_array, csc_array, csc_array]

Notes

The flattened grid-cell numbering convention is:

node_number = ix + iy * grid.nx + iz * grid.ny * grid.nx

The coefficient magnitude along each axis is computed as:

interface_area / cell_volume

which gives:

  • x-axis: grid.gamma_ij_x_m2 / grid.grid_cell_volume_m3

  • y-axis: grid.gamma_ij_y_m2 / grid.grid_cell_volume_m3

  • z-axis: grid.gamma_ij_z_m2 / grid.grid_cell_volume_m3

Boundary cells without a valid neighbour along a given direction do not receive a contribution for that direction.