dragonfly_sim.core.reconstruction

Green-Gauss gradient reconstruction for the p = 0 discretization.

At p = 0, grad U vanishes inside every cell, which leaves no viscous volume term and no vorticity for the turbulence sources, and gives only first-order face states. As a cell-centred finite-volume code does, a DG0 gradient G is reconstructed by Green-Gauss over each cell’s facets, with facet averages {U} inside and the registered boundary state U_b on the boundary,

|K| G_K = sum_f |f| U_f (x) n_f.

G is used in three places:

• MUSCL: the interior-facet traces of the inviscid flux are reconstructed linearly, U_K + G_K . (x_f - x_K), optionally with a smooth van Albada limiter (see muscl_traces).

• The viscous facet gradient, {G} with its component along the line between the two cell centres replaced by the direct difference (which keeps the two-point coupling that suppresses odd-even modes), or on a boundary facet with the normal component replaced by (U_b - U)/d_n.

• The turbulence sources (vorticity, gradients of the transported scalars).

G is eliminated from the unknowns exactly (utils/gradient_condensation.py), so Newton sees the true Jacobian. For adjoints and shape derivatives G and the cell centres are DerivedFields (utils/derived_fields.py), with exact derivatives with respect to U, alpha and the mesh coordinates.

Attributes

_JIT

Classes

CellCentres	DG0 cell centres, each cell's reference midpoint mapped to the cell. For
GreenGaussGradient	G = M^-1 b(U, alpha, x), as a DerivedField (see the module docstring).
GreenGaussReconstruction	p = 0 gradient reconstruction attached to a CompressibleEulerModel host.
GreenGaussStateOperator	dG/dU = M_G^-1 db/dU as a python-context Mat (G layout x state layout)

Module Contents

class dragonfly_sim.core.reconstruction.CellCentres(mesh_obj)

Bases: dragonfly_sim.utils.derived_fields.DerivedField

DG0 cell centres, each cell’s reference midpoint mapped to the cell. For the degree-1 geometry used here that is the plain average of the cell’s vertices, c = C x, so dc/dx is the constant averaging matrix C.

_mesh_jacobian_loop()

_mesh_jacobian_vectorized()

The same matrix as _mesh_jacobian_loop, as one CSR build: row (cell, component i) holds 1/nv at the cell’s vertex dofs’ component i.

discard(obj)

Free a matrix or vector returned by state_jacobian / alpha_vector / mesh_jacobian once the caller is done with it. A field that caches them keeps its own (override); a python-context operator also frees what its context holds.

mesh_jacobian()

release_cache()

Free whatever the field caches between calls (nothing by default).

update_geometry()

_mesh_jac = None

function

mesh_obj

rec = None

class dragonfly_sim.core.reconstruction.GreenGaussGradient(reconstruction)

Bases: dragonfly_sim.utils.derived_fields.DerivedField

G = M^-1 b(U, alpha, x), as a DerivedField (see the module docstring).

Options (attributes of the reconstruction):

cache_jacobians keep dG/dU, dG/dalpha and dG/dx between calls

while the state, alpha and the mesh are unchanged (several functionals and products of one totals evaluation share them); freed by release_cache

structured_state_jacobian dG/dU as the operator W (x) I plus a small

boundary part (GreenGaussStateOperator) instead of an assembled matrix

_build_alpha_vector()

_build_mesh_jacobian()

_build_state_jacobian()

_cached(name, build)

build() once per (state, alpha, mesh), or every call if caching is off.

_form(key, factory)

_key()

_same_key(key)

alpha_vector()

discard(obj)

Free a matrix or vector returned by state_jacobian / alpha_vector / mesh_jacobian once the caller is done with it. A field that caches them keeps its own (override); a python-context operator also frees what its context holds.

mesh_jacobian()

release_cache()

Free whatever the field caches between calls (nothing by default).

state_jacobian()

update()

update_geometry()

_cache = None

_forms

depends_on_state = True

function

rec

class dragonfly_sim.core.reconstruction.GreenGaussReconstruction(model, muscl=True, limiter=None, limiter_eps=0.01, limiter_scales=None)

p = 0 gradient reconstruction attached to a CompressibleEulerModel host.

install() plugs it in: the MUSCL face traces (if muscl), the derived fields (cell centres and G), and the condensed nonlinear problem. Call it after the function spaces exist and before the boundary conditions are defined; the Green-Gauss forms are built lazily, once every boundary condition has registered its state.

limiter: None (unlimited MUSCL) or “vanalbada”. limiter_scales() returns one magnitude per state component for the limiter’s smoothing parameter.

_weights()

W, the scalar operator of the interior-facet part of dG/dU before M_G^-1: (3 n_cells) x n_cells, assembled once per mesh.

compact_forms(F)

(a_UU, a_UG): the compiled forms of dR_U/dU and dR_U/dG, cached for F.

condensed_jacobian(F)

The exact dR/dU of the residual form F, as a CondensedJacobian.

face_gradient_boundary(U_b, U=None)

The cell’s G with its normal component replaced by (U_b - U)/d_n, d_n the cell centre’s normal distance from the facet.

face_gradient_interior()

green_gauss_rhs()

b(U; tau) with |K| G_K = b: facet averages inside, the registered boundary states on the boundary.

install()

make_problem(F_newton, J_newton)

muscl_traces()

Linear reconstruction of both interior-facet traces from G at the facet’s quadrature points: U_K + G_K . (x_f - x_K) unlimited, or with limiter=”vanalbada” the edge-based van Albada form (Blazek, Computational Fluid Dynamics, ch. 5) generalised to an off-centre facet. For the cell K and its neighbour N across the facet, with s = |x_f - x_K| / |x_N - x_K|,

D+ = s (U_N - U_K) (central increment to the facet) D- = 2 G_K . (x_f - x_K) - D+ (upwind increment) U_f = U_K + Psi(D-, D+), Psi(a, b) = ((a^2 + e) b + (b^2 + e) a) / (a^2 + b^2 + 2e),

componentwise. Psi(a, a) = a, so smooth regions keep the unlimited reconstruction; at a shock D- ~ 0 and the trace falls back to first order. Psi is smooth, so the exact Jacobian survives. e = (limiter_eps * S_k)^2 with S_k the component’s scale.

release_cache()

set_up()

state_operator()

dG/dU at the current state, as a python-context Mat.

update_gradient()

G <- M_G^-1 b(U) at the current state and mesh. Collective.

update_mass()

Refresh M_G^-1 in place after mesh motion.

G

_W = None

_W_key = None

_gg_form = None

cache_jacobians = False

centres

fast_cell_centre_jacobian = False

gg_rhs = None

gradient

inv_mass = None

limiter = None

limiter_eps = 0.01

limiter_scales = None

mesh_obj

model

muscl = True

structured_state_jacobian = False

xc

class dragonfly_sim.core.reconstruction.GreenGaussStateOperator(rec, W, B_bdry)

dG/dU = M_G^-1 db/dU as a python-context Mat (G layout x state layout) with mult and multTranspose, without assembling the 18-by-6 blocks.

The interior-facet part of b is a facet average of U times the area vector, component by component: its derivative is W (x) I_n, with W the scalar operator (3 area-vector components per cell pair, assembled once per mesh) and I_n the identity over the n state components. Only the boundary-facet part, U_b(U) of the registered boundary states, couples components; it is assembled per state (boundary cells only).

destroy(mat=None)

mult(mat, x, y)

multTranspose(mat, v, y)

_t

n

dragonfly_sim.core.reconstruction._JIT