Pde

Level-set fire front propagation PDE: Hamilton-Jacobi equation ∂ψ/∂t + S(n̂)|∇ψ| = 0, where S(n̂) is the Rothermel/Mandel (2011) direction-dependent spread rate projected onto the fire-front normal. Migrated from WildlandFire.jl src/level_set_fire_spread.jl:73. v0.8.0: a DISCRETIZATION-AGNOSTIC library PDE leaf. (1) Spatial derivatives are bare `D` ops (op:D, wrt:x / wrt:y) instead of the `grad` sugar. (2) This leaf imports NO discretization of its own — the scheme is chosen by the ASSEMBLER that mounts this component, which carries the rule imports on its model-ref edge (esm-spec §9.7.10 form A: the edge's `expression_template_imports` register into this component's own template scope). The reference assembler (wildlandfire.esm/level_set_hj_2d.esm) mounts the first-order Godunov (Rouy-Tourin / Osher-Sethian) rule godunov_norm_D2_periodic, which fires by compound match on the gradient-magnitude term sqrt(D(psi,x)² + D(psi,y)²), plus the central rules central_D1x/D1y for the individual ∂ψ/∂x, ∂ψ/∂y in the wind/slope normal projections — but any grid + scheme satisfying the same D-op / grid-geometry contract composes here unchanged. (3) The direction-dependent right-hand side is authored as READABLE INTERMEDIATE array-shaped observed fields over the [x,y] index sets — psi_x, psi_y (central ∂ψ), grad_mag (Godunov |∇ψ|), grad_safe, U_n, tanphi_n (wind/slope normal projections), phi_W, phi_S, and S_n (the spread rate) — so the state equation is simply D(psi,t) = −S_n·grad_mag. Both toolkits evaluate these per cell: the Julia evaluator folds the elementwise ones into the state RHS (leaving the producer stencils psi_x/grad_mag for its index beta-reduction), and the Rust array runtime materializes each as a companion buffer in dependency order — identical numerics either way (front radii match to 4 places at NX=NY=32). The producer observeds (psi_x, psi_y, grad_mag) carry the bare D ops the mounted discretization lowers; grad_mag's exponent is the integer 2 on the bare D ops so the Godunov compound matches. Only the scalar grid observeds (dx, dy) remain scalar. (4) Grid geometry (index sets x,y; spacing dx=Lx/NX, dy=Ly/NY; extent x0,y0) is the consumer-supplied contract of grids/cartesian_uniform_2d/grid.esm (metaparameters NX,NY bound at load; extent Lx,Ly parameters), and BCs come from the mounted discretization. The initial perimeter is supplied by the ASSEMBLER (a §11.4.1 dedicated IC model), not baked here. The D ops assume Cartesian (x,y in meters) coordinates; geographic (lon/lat) usage needs metric-tensor scaling not yet expressible in ESM (other:P3-A-metric-tensor-coord-transforms gap).

The level-set propagation equation of fire_behavior (physics/level_set_mod.F90 Calc_tend_ls; Munoz-Esparza et al. 2018) as a PDE on a uniform 2-D grid: d psi/dt = -R |grad psi| + nu |R| (dx d2psi/dx2 + dy d2psi/dy2), the Hamilton-Jacobi front propagation with the rate of spread R normal to the front plus a grid-scaled artificial viscosity. |grad psi| is written as the compound sqrt(D(psi,x)^2 + D(psi,y)^2) and lowered by EarthSciDiscretizations' grids/cartesian_uniform_2d/rules/weno5_eno1_norm_D2_extrapolate_bc, which is fire_behavior's default scheme (fire_upwinding = 9): fifth-order WENO within fire_lsm_band_ngp = 4 cells of the front with the upwind side chosen by the mid-flame wind, first-order ENO (minmod of the one-sided differences) elsewhere and within BDY_ENO1 = 10 cells of the domain boundary. The second derivatives of the viscosity term are lowered by central_D2x_extrapolate_bc and central_D2y_extrapolate_bc. All three rules close the walls as fire_behavior does (Extrapol_var_at_bdys): each wall cell sees a ghost extrapolated linearly from the two innermost cells, max(2a - b, a, b), so the domain edge carries no prescribed level-set value and does not act as a fire front. R is supplied as an input field (it depends on the upwind normal, tested in ROSWRFFire); the propagation term and viscosity are what this component asserts. The viscosity coefficient is the single value fire_viscosity: fire_behavior blends it toward fire_viscosity_bg near the front, and test7 sets the two equal. Tests: the ls_tend dump of Runge-Kutta stage 1 of step 20 (t = 10 s) of test7 with the default scheme, checked at wall, corner, ENO1-band, WENO-band and far-field cells, real32.

Full 1D Saint-Venant PDE system for surface overland flow: mass conservation ∂h̃/∂t = −∂q/∂l + (P − I) and momentum conservation ∂q/∂t = −∂F/∂l + g·h̃·(S₀ − Sᶠ), where the momentum flux F = q²/h̃ + g·h̃²/2 is an auxiliary (observed) PDE field. Migrated from EnvironmentalTransport.jl@a48559c4 src/surface_runoff.jl:117 (SaintVenantPDE builder). This is a DISCRETIZATION-AGNOSTIC leaf (esm-spec §9.7.10): the spatial fluxes are written with the rewrite-target spatial derivative op `D(·, wrt: l)` and the component imports NO discretization of its own, so it can be composed under different schemes. Its inline tests and examples inject a discretization (form C) — here the sibling adapter `saint_venant_upwind1_inflow.esm` (first-order upwind, per-variable inflow BC) — for that run only; the persisted leaf keeps its spatial `D` intact. The 1D spatial variable l is the along-slope coordinate; the sibling SurfaceRunoff component (surface_runoff_component.esm) is the single-node ODE reduction of this same PDE system with the spatial-divergence terms exposed as input parameters.

The YSU (Yonsei University) nonlocal boundary-layer scheme of Hong, Noh and Dudhia (2006), with the stable-regime (Hong 2010) and top-down (theta_li) boundary-layer-top diagnostics of WRF 4.8, written as a single-column PDE component on the EarthSciDiscretizations grids/column_nonuniform_1d grid (layer index `lev`, layer 1 at the surface; interface index `lev_nodes`, interface 1 the surface). The physics is transcribed from the CCPP kernel bl_ysu_run in phys/physics_mmm/bl_ysu.F90 (NCAR/MMM-physics, the WRF 4.8 submodule). Prognostic column fields theta (potential temperature), u, v, qv, qc and qi obey the density-weighted flux-form vertical diffusion d(c)/dt = (1/rho) d/dz[ rho K_c dc/dz - Phi_c ], where K_c is the interface eddy diffusivity of the field (K_h for theta, K_m for u and v, K_q for water) and Phi_c is the prescribed (non-gradient) part of the turbulent flux: the surface flux on interface 1 (-hfx/cp for theta, -qfx for qv, rho_s u_*^2 u_1/|U| for momentum -- the sign convention of the ESD flux rule, Phi = rho K dc/dz on the face) and, inside the boundary layer, the YSU nonlocal (counter-gradient) flux rho K_c gamma_c plus the entrainment flux (w'c')_h ((z - z_1)/(h - z_1))^3, zero on the top interface. The gradient term is lowered by the ESD rule varcoeff_face_laplacian_lev_flux_bc applied to (rho K) on the interfaces (its own prescribed face fluxes kdudz_bot/kdudz_top are zero here because every prescribed flux is carried by Phi so that the six prognostic fields can share one rule instance), the prescribed-flux term by face_flux_D_lev_supplied_faces. The density weighting is exactly the hydrostatic sigma-coordinate form of the Fortran: rho on an interface is (p[k] - p[k+1])/(g (zc[k+1] - zc[k])) and in a layer (p_int[k] - p_int[k+1])/(g dz[k]), so with WRF's pressures and heights the lowered operator reproduces bl_ysu_run's tridiagonal coefficients term by term (dsig xkzh rdz / del). Boundary-layer height: the bulk-Richardson search of components/atmospheric_dynamics/ysu/bulk_richardson.esm is applied (1) with the lowest-level theta_v, (2) with the thermal excess of the surface heat and moisture fluxes when the surface layer is unstable and the buoyancy flux positive, (3) the top-down theta_li scan that raises the top to just above the highest theta_li-unstable level (ysu_topdown_pblmix; it starts at the previous kpbl and evaluates the bulk Richardson number from theta_li at every level it visits, level 1 included -- brli1 -- rather than seeding level 1 with the surface-layer br), and (4) the stable-regime search with the Hong (2010) critical value (0.25 over land, Rossby-number dependent over water) when the surface layer is stable and no boundary layer was found; each search is an aggregate over the level index and never a recurrence (see bulk_richardson.esm). Diffusivities: inside the boundary layer the K-profile K_m = k w_s z (1 - z/h)^p with the surface-layer velocity scale w_s = (u_*^3 + phifac k w_*^3 (1 - z/h))^(1/3) (or u_*/phi_m(z) when stable), K_h = K_m / Pr with the height-tapered Prandtl number of Hong et al. (2006, their eq. A10), K_q likewise; above the boundary layer the local Ri-dependent mixing-length diffusivity of Louis type with the in-cloud moist Richardson number; in the entrainment zone ((z - h)/delta)^2 < 4.6 the geometric mean of the local value and the entrainment-flux value -w_e dz_h exp(-((z - h)/delta)^2). Every diffusivity is floored by the background values and capped by xkzmax. Consumer-supplied inputs: the column fields (p, p_int, exner, ze and the initial theta, u, v, qv, qc, qi) are REWRITE-TARGET ops `input_<name>` (esm-spec 4.2 open tier) that the consumer lowers -- in the inline tests through a per-test `expression_template_imports` library (tests/<test id>_inputs.esm) because esm-spec 6.6 `parameter_overrides` admit scalars only; surface-layer scalars (hfx, qfx, u_*, psi_m, psi_h, br, |U|, U10, land mask, z0, psfc) are ordinary parameters. NOT transcribed (inactive in every test regime): the BEP/BEP+BEM urban forcing terms, the optional ocean-current (uox/vox) surface drag, the topo_wind corrections (ctopo = ctopo2 = 1), the extra passive-tracer mixing (nmix = 0), and the WRF diagnostic outputs u10/v10 (pass-through), dusfc/dvsfc/dtsfc/dqsfc (column sums) and the TKE/hybrid-PBLH diagnostics. The cloud-top radiative-entrainment branch (ysu_topdown_pblmix, WRF cloudflg: a cloud layer at kpbl-1 whose integrated radiative cooling radsum drives top-down thermals wstar3_2, extra entrainment we_rad, a counter-gradient term hgamt2 and a second K-profile term wscalek2) is transcribed (exercised by synthetic_cloudtop_call120, and switched off by topdown = 0 in synthetic_notopdown_cloudtop_call120); radsum is the cumulative column integral of the radiative cooling. Fractional powers of quantities that vanish identically outside their regime (w_*^3 in a stable column, wstar3_2 below a cloud-free top) are gated with ifelse, so that a Jacobian of the right-hand side never evaluates the derivative of a root at zero. Tendencies are observed variables (`dthdt`, `dudt`, `dvdt`, `dqvdt`, `dqcdt`, `dqidt`) so that the inline tests assert the instantaneous right-hand side at time 0 against the real64 kernel's tendency extrapolated to dt -> 0 (Richardson extrapolation of the implicit step over dt = 0.04, 0.02, 0.01 s; the Fortran solves the diffusion implicitly, so its finite-step tendency differs from the derivative by O(dt K / dz^2)). Test provenance: the inputs and expected values come from a WRF v4.8.0 em_scm_xy single-column run of the instrumented fork ctessum-claude/WRF (branch earthsciml-instrumented, module_esm_dump; physics_mmm at ctessum-claude/MMM-physics 52cfbd7) and from replays of the dumped kernel calls through that fork's real64 kernels/ysu_driver; the synthetic_* tests replay a dumped column after the modification their description states. The numbers were transcribed by tools/fill_tests.py of EarthSciML/EqWeFiC. The dumps, the driver and the generator are not part of this repository, so the values are recorded data rather than regenerable output.