zoomy_core.model.models.basisfunctions module

zoomy_core.model.models.basisfunctions module#

Module zoomy_core.model.models.basisfunctions.

class zoomy_core.model.models.basisfunctions.Basisfunction(level=0, symbol='phi', **kwargs)#

Bases: object

Basisfunction. (class).

Each instance carries:

self.basis — list of concrete polynomial expressions in z (used by eval() and numeric reconstruction).
self.phi — list of opaque sympy Function subclasses phi_0, phi_1, …, phi_L (one per index k). Each subclass carries class attributes _basis = self and _index = k. When an Integral’s integrand contains atoms whose func._basis points back to a Basisfunction, EvaluateIntegrals routes the integral to that basis’s evaluate_integral() for resolution against the concrete polynomial — without the model author ever substituting the polynomial form themselves.

The Function subclass __name__ is "{symbol}_{k}" (e.g. phi_0, eta_2, mu_1), so different bases coexist in the same equation as distinct sympy classes — mixing is automatic.

name = 'Basisfunction'#

bounds()#: Bounds.

basis_definition()#: Basis definition.

weight(z)#: Weight.

mean_coefficients()#

Return coefficients c_k such that sum(c_k * phi_k(z)) = 1.

Used for computing the depth-averaged velocity:: u_mean = sum(c_k * alpha_k)

Default: c = [1, 0, 0, …] (assumes phi_0 = 1). Override for bases where the constant function is a non-trivial combination (e.g., B-splines: c = [1, 1, 1, …]).

weight_eval(z)#: Weight eval.

resolve_atoms(expr)#

Replace every phi_fn(k_concrete, arg) atom of this basis with the concrete polynomial self.eval(int(k_concrete), arg).

Calls with a symbolic k (i.e. atoms inside an unevaluated sp.Sum whose summation index hasn’t been substituted) are left untouched — they’ll resolve once the Sum is .doit()-ed and the index becomes a concrete integer.

Atoms of other bases (different _basis back-reference) match different Function classes and are left untouched, so multi-basis expressions resolve correctly by calling resolve_atoms once per basis present.

get(k)#: Get.

eval(k, _z)#: Eval.

eval_psi(k, _z)#: Eval psi.

get_lambda(k)#: Get lambda.

plot(ax)#: Plot.

reconstruct_velocity_profile(alpha, N=100)#: Reconstruct velocity profile.

reconstruct_velocity_profile_at(alpha, z)#: Reconstruct velocity profile at.

reconstruct_alpha(velocities, z)#: Reconstruct alpha.

project_onto_basis(Y)#: Project onto basis.

get_diff_basis()#: Get diff basis.

class zoomy_core.model.models.basisfunctions.Monomials(level=0, symbol='phi', **kwargs)#

Bases: Basisfunction

Monomials. (class).

name = 'Monomials'#

class zoomy_core.model.models.basisfunctions.LayeredBasis(inner_cls=None, inner_level=0, interfaces=None, n_layers=None, **inner_kwargs)#

Bases: Basisfunction

Composable multi-layer basis: an inner basis rescaled onto each of N sub-intervals of [0, 1] (or physical z-space).

For an inner basis with m = inner_level + 1 functions and a partition into N layers, LayeredBasis has N * m basis functions. Basis function j = i*m + k is the inner function phi_k rescaled onto layer i’s local coordinate zeta_i = (z - lower_i) / h_i and zero outside the layer.

This subsumes several useful cases:

Piecewise constant (multi-layer SWE): inner_cls=Monomials, inner_level=0. m = 1, phi_inner_0 = 1, so each basis function is an indicator on its layer and basis.alpha lists the layer-average velocities.
Multi-layer SME: inner_cls=Legendre_shifted, inner_level = L. Each layer carries L + 1 moments; basis.alpha is flat-indexed as alpha_{i*m + k}.
Single-layer smooth basis: n_layers=1 recovers the plain inner basis on [0, 1].

Parameters:

inner_cls (Basisfunction subclass, default Monomials) – The per-layer basis. Default produces piecewise-constant.
inner_level (int, default 0) – Polynomial level of the inner basis.
interfaces (list of sympy, optional) – Explicit breakpoints (len N + 1). May live in either ζ-space ([0, 1/2, 1]) or physical z-space ([state.b, z_1, state.eta]).
n_layers (int, optional) – If interfaces isn’t given, create n_layers + 1 equal breakpoints on [0, 1].
**inner_kwargs – Forwarded to inner_cls(level=inner_level, **inner_kwargs).
Usage:: –
# Two-layer SWE (piecewise constant on physical interfaces): bf = LayeredBasis(interfaces=[state.b, z_1, state.eta])

# Two-layer level-1 SME (Legendre inside each layer): bf = LayeredBasis(Legendre_shifted, inner_level=1,

interfaces=[state.b, z_1, state.eta])

# Three equal layers in ζ-space, linear inner: bf = LayeredBasis(Legendre_shifted, inner_level=1, n_layers=3)

name = 'LayeredBasis'#

basis_definition()#: Basis definition.

bounds()#: Bounds.

weight(z)#: Weight.

mean_coefficients()#

c_{i*m+k} such that sum(c_j * phi_j) = 1 where defined.

For a partition, exactly one layer is active at each point. Inside layer i, the inner basis reconstructs the constant 1 via its own mean coefficients. Flatten them into the layered ordering.

flat_index(layer_idx, moment_idx)#: (layer_idx, moment_idx) → flat basis index.

layer_zeta(layer_idx, z_value)#: Inner coordinate ζ_i(z_value) = (z_value - lower_i) / h_i.

class zoomy_core.model.models.basisfunctions.Legendre_shifted(level=0, symbol='phi', **kwargs)#

Bases: Basisfunction

Legendre shifted. (class).

name = 'Legendre_shifted'#

basis_definition()#

Standard shifted Legendre on σ ∈ [0, 1].

phi_k(σ) = legendre(k, 2σ − 1) — the canonical convention:

phi_k(σ=0) = (−1)^k (BOTTOM of the layer)
phi_k(σ=1) = 1 (TOP / FREE SURFACE)

Removed the historical (−1)^k pre-factor that swapped these endpoint values — it had inverted the physical interpretation of odd modes (q_1 > 0 meant u(bottom) > u(top), opposite to the intuitive “shear” direction). The new convention matches K&T 2019 §3 directly and makes the multi-layer upwind formulas u_ℓ(σ=1) = Σ q_ℓ_k / h_ℓ (top of layer = sum of modes) used in MLSME / MLVAM correct without any further changes.

analytical_weighted_integral(poly_expr, var)#

Compute int_0^1 poly_expr(z) * 1 dz exactly via antiderivative.

Since the Legendre weight is 1, this is just the polynomial antiderivative evaluated at the bounds. ~100x faster than sympy.integrate for polynomial integrands.

class zoomy_core.model.models.basisfunctions.Chebyshevu(level=0, symbol='phi', **kwargs)#

Bases: Basisfunction

Chebyshevu. (class).

name = 'Chebyshevu'#

bounds()#: Bounds.

weight(z)#: Weight.

basis_definition()#: Basis definition.

class zoomy_core.model.models.basisfunctions.Chebyshevu_shifted(level=0, symbol='phi', **kwargs)#

Bases: Basisfunction

Chebyshev U polynomials shifted to [0,1] with phi_0 = 1.

phi_k(z) = U_k(2z - 1) (unnormalized Chebyshev of the second kind) weight = sqrt(z * (1 - z)) domain = [0, 1]

Properties: - phi_0 = 1 (recovers SWE at level 0) - Orthogonal: M[i,j] = (pi/8) * delta_{ij} - Analytical quadrature nodes: z_k = (1 + cos(k*pi/(n+1))) / 2

name = 'Chebyshevu_shifted'#

bounds()#: Bounds.

weight(z)#: Weight.

basis_definition()#: Basis definition.

analytical_weighted_integral(poly_expr, var)#

Compute int_0^1 poly_expr(z) * sqrt(z*(1-z)) dz exactly.

Uses Chebyshev U orthogonality: expand poly_expr in the U_k(2z-1) basis, then int U_k * w dz = pi/8 * delta_{k,0}.

Returns c_0 * pi/8 where c_0 is the U_0 expansion coefficient. Returns None if poly_expr is not a polynomial.

quadrature_nodes(n=None)#

Gauss-Chebyshev U quadrature on [0,1].

Integrates integral(f(z) * sqrt(z*(1-z)), 0, 1) exactly for polynomial f of degree <= 2*n + 1.

Uses n+1 nodes. Default n = 2*(level+1) for safety with triple products.

class zoomy_core.model.models.basisfunctions.Legendre_DN(level=0, symbol='phi', **kwargs)#

Bases: Basisfunction

Legendre DN. (class).

name = 'Legendre_DN - satifying no-slip and no-stress. This is a non-SWE basis'#

bounds()#: Bounds.

basis_definition()#: Basis definition.

class zoomy_core.model.models.basisfunctions.GalerkinBasis(level=0, parent='legendre', bc_bottom='slip', bc_top='nostress', slip_length=None, **kwargs)#

Bases: Basisfunction

General BC-aware basis via Shen-type recombination.

Constructs basis functions from a parent polynomial family where: - phi_0 = 1 (constant, for mass conservation / SWE limit) - phi_k (k >= 1) are linear combinations of parent polynomials that

automatically satisfy user-specified linear boundary conditions

Supported BC types at bottom (z=bounds[0]) and top (z=bounds[1]): - ‘free’: no constraint (free surface, du/dz = 0 implied by no BC) - ‘noslip’: u(z_bc) = 0 - ‘nostress’: du/dz(z_bc) = 0 - ‘slip’: du/dz(z_bc) = u(z_bc) / slip_length

The BC-satisfying basis is built by solving a small linear system for each k, combining parent polynomials P_k, P_{k+1}, P_{k+2} with coefficients that enforce the two BCs.

Usage:: basis = GalerkinBasis(level=3, parent=’legendre’, bc_bottom=’noslip’, bc_top=’nostress’) basis = GalerkinBasis(level=2, parent=’chebyshev’, bc_bottom=’slip’, bc_top=’free’,

slip_length=0.5)

name = 'GalerkinBasis'#

bounds()#: Bounds.

weight(z)#: Weight.

basis_definition()#: Basis definition.

class zoomy_core.model.models.basisfunctions.SplineBasis(level=0, degree=1, **kwargs)#

Bases: Basisfunction

B-spline basis on [0,1] using raw B-splines (nodal DOFs).

phi_k = B_k (k-th B-spline hat function) - phi_0 peaks at z=0 (bottom): phi_0(0) = 1 - phi_{n-1} peaks at z=1 (top): phi_{n-1}(1) = 1 - Interior phi_k peak at internal knots

The partition of unity (sum B_k = 1) means: - Depth-averaged velocity = weighted sum of nodal velocities - Mass flux = h * sum(alpha_k * integral(B_k)) - Bottom velocity = alpha_0 (direct nodal access)

Provides get_knot_spans() for piecewise integration.

name = 'SplineBasis'#

bounds()#: Bounds.

weight(z)#: Weight.

basis_definition()#: Basis definition.

get_knot_spans()#: Return list of (a, b) intervals between consecutive knot values.

mean_coefficients()#: B-splines form a partition of unity: sum(B_k) = 1, so c_k = 1 for all k.

class zoomy_core.model.models.basisfunctions.Spline(level=0, symbol='phi', **kwargs)#

Bases: Basisfunction

Spline. (class). Legacy — use SplineBasis instead.

name = 'Spline'#

basis_definition(degree=1, knots=[0, 0, 0.001, 1, 1])#: Basis definition.

class zoomy_core.model.models.basisfunctions.OrthogonalSplineWithConstant(level=0, symbol='phi', **kwargs)#

Bases: Basisfunction

OrthogonalSplineWithConstant. (class).

name = 'OrthogonalSplineWithConstant'#

basis_definition(degree=1, knots=[0, 0, 0.5, 1, 1])#: Basis definition.