include/pops/core/model/physical_model.hpp Source FileΒΆ

adc_cpp: include/pops/core/model/physical_model.hpp Source File
adc_cpp 0.3.0
Model-free C++23 core for coupled hyperbolic-elliptic systems on adaptive (AMR) meshes, with MPI and GPU (Kokkos) backends
physical_model.hpp
Go to the documentation of this file.
1
17
18#pragma once
19
20#include <pops/core/state/state.hpp> // Aux: the contract fixes the auxiliary to pops::Aux
22#include <pops/core/state/variables.hpp> // Variables: mandatory contract of the hyperbolic model
23
24#include <concepts>
25
26// The contract of the physics layer.
27//
28// A PhysicalModel describes ONE equation: its pointwise formulas. Nothing
29// more. It is the only "what to compute" axis of the architecture, separate from
30// the "where / how to iterate" axis (mesh + dispatch) and the "in what order" axis
31// (integrator + coupler).
32//
33// Everything is a pure function of pointwise states:
34// - flux(U, aux, dir): the physical flux in direction dir
35// - max_wave_speed(U, aux, dir): the largest wave speed (for the CFL
36// and the Riemann solver)
37// - source(U, aux): the pointwise source term
38// - elliptic_rhs(U): the right-hand side of the elliptic equation
39// (charge / mass density depending on the model)
40//
41// flux AND source take aux: this is the point that unifies drift transport
42// (aux in the flux) and the self-gravitating compressible fluid (aux in the
43// source) under one same spatial operator.
44//
45// Aux contract (slice, see milestone 4): the auxiliary is FIXED to pops::Aux (phi, grad
46// phi). This is what load_aux builds and all the spatial operator provides;
47// the concept therefore requires it explicitly (M::Aux == pops::Aux) rather than letting
48// one believe that a model could declare an arbitrary auxiliary that the code would
49// never fill. Generalizing load_aux<Model> to an arbitrary Model::Aux remains
50// possible later; meanwhile the contract says exactly what the code provides.
51
52namespace pops {
53
62// POPS_HD : aux_comps() est evaluee a la compilation (argument de template non-type
63// load_aux<aux_comps<Model>()>) DANS les kernels device (cf. spatial_operator_eb.hpp). Sous nvcc,
64// appeler une constexpr __host__ depuis une fonction __host__ __device__ est refuse (#20013-D) ;
65// la marquer POPS_HD la rend callable des deux cotes. Hors nvcc, POPS_HD est vide -> constexpr pur.
66template <class M>
67POPS_HD constexpr int aux_comps() {
68 if constexpr (requires { M::n_aux; })
69 return M::n_aux;
70 else
71 return kAuxBaseComps;
72}
73
80template <class M>
82 requires(const M m, const typename M::State u, const typename M::Aux a, int dir) {
83 typename M::State;
84 typename M::Aux;
85 requires std::same_as<typename M::Aux, Aux>;
86 { M::n_vars } -> std::convertible_to<int>;
87 { m.flux(u, a, dir) } -> std::same_as<typename M::State>;
88 { m.max_wave_speed(u, a, dir) } -> std::convertible_to<Real>;
89 { m.source(u, a) } -> std::same_as<typename M::State>;
90 { m.elliptic_rhs(u) } -> std::convertible_to<Real>;
91 };
92
93// ---------------------------------------------------------------------------------------------
94// OPTIONAL TIME STEP BOUNDS of the model contract (audit 2026-06, "step_cfl" workstream).
95//
96// Historically, step_cfl knew ONLY the hyperbolic transport bound
97// dt <= cfl * h / max_wave_speed. But a model may impose other bounds: source frequency
98// (collision/reaction, mu = eig(dS/dU), unit 1/time, WITHOUT h), or directly an admissible
99// step (coupled transport-source formula not reducible). These three OPTIONAL traits let the
100// model declare them; a model that declares none keeps STRICTLY the historical
101// behavior (max_wave_speed fallback, bit-identical).
102//
103// SEMANTICS (all bounds apply to the EFFECTIVE SUBSTEP stride*dt/substeps of the block,
104// see SystemStepper::step_cfl):
105// - stability_speed(U, aux, dir): stability speed lambda* [length/time] which REPLACES
106// max_wave_speed in the block CFL reduction (dt <= cfl * h / max_cells(lambda*)). For
107// when the speed relevant for STABILITY is not the physical wave speed (declared
108// conservative bound, speed modified by a coupling...). The Riemann solvers, themselves,
109// keep reading max_wave_speed (accuracy != stability).
110// - source_frequency(U, aux): local frequency mu [1/time] of the local source/coupling;
111// imposes dt <= cfl / max_cells(mu) -- NO h (the source bound has no space
112// dimension). Shortcut for explicit relaxation/collision/reaction.
113// - stability_dt(U, aux): direct ADMISSIBLE step [time] per cell; imposes
114// dt <= min_cells(stability_dt). The cfl is NOT applied (the model already declares an
115// admissible step; applying cfl on top would mix two margins). This is the most general form.
116//
117// STABILITY vs ACCURACY: these traits declare STABILITY bounds. A source treated
118// implicitly (SourceImplicit/IMEX) may no longer impose a stability bound while keeping an
119// ACCURACY constraint: it is up to the model to choose what stability_dt/source_frequency
120// return in that case (or to not declare them). NON-local bounds (multi-block
121// coupling, Schur/Poisson, AMR/scheduler) do NOT go through these cell-by-cell traits:
122// they go through System::add_dt_bound (global host bound, one evaluation per step).
123//
124// GPU/MPI PRODUCTION: like flux/source, these methods must be POPS_HD (evaluated in
125// reduction kernels) -- a per-cell Python callback is not a production path;
126// the DSL compiles them (m.stability_speed(...) / m.stability_dt(...)).
127// ---------------------------------------------------------------------------------------------
128
130template <class M>
131concept HasStabilitySpeed = requires(const M m, const typename M::State u, const Aux a, int dir) {
132 { m.stability_speed(u, a, dir) } -> std::convertible_to<Real>;
133};
134
136template <class M>
137concept HasSourceFrequency = requires(const M m, const typename M::State u, const Aux a) {
138 { m.source_frequency(u, a) } -> std::convertible_to<Real>;
139};
140
142template <class M>
143concept HasStabilityDt = requires(const M m, const typename M::State u, const Aux a) {
144 { m.stability_dt(u, a) } -> std::convertible_to<Real>;
145};
146
153template <class M>
154concept HasPointwiseProjection = requires(const M m, const typename M::State u, const Aux a) {
155 { m.project(u, a) } -> std::same_as<typename M::State>;
156};
157
166template <class M>
168 PhysicalModel<M> && requires(const M m, const typename M::State u, const typename M::Prim p) {
169 typename M::Prim;
170 { m.to_primitive(u) } -> std::same_as<typename M::Prim>;
171 { m.to_conservative(p) } -> std::same_as<typename M::State>;
172 };
173
180template <class M>
182 requires(const M m, const typename M::State u, const typename M::Prim p, const Aux a, int dir) {
183 typename M::State;
184 typename M::Prim;
185 { M::n_vars } -> std::convertible_to<int>;
186 { m.flux(u, a, dir) } -> std::same_as<typename M::State>;
187 { m.max_wave_speed(u, a, dir) } -> std::convertible_to<Real>;
188 { m.to_primitive(u) } -> std::same_as<typename M::Prim>;
189 { m.to_conservative(p) } -> std::same_as<typename M::State>;
190 { M::conservative_vars() } -> std::same_as<VariableSet>;
191 { M::primitive_vars() } -> std::same_as<VariableSet>;
192 };
193
195template <class M>
197
198} // namespace pops
Trait OPTIONNEL : PROJECTION PONCTUELLE post-pas U -> project(U, aux) (ADC-177).
Definition physical_model.hpp:154
OPTIONAL extension of a PhysicalModel: primitive variables + cons<->prim conversions.
Definition physical_model.hpp:167
OPTIONAL trait: local source frequency mu [1/s] (bound dt <= cfl / max mu, without h).
Definition physical_model.hpp:137
OPTIONAL trait: direct admissible step per cell (bound dt <= min stability_dt, without cfl).
Definition physical_model.hpp:143
OPTIONAL trait: stability speed lambda* replacing max_wave_speed in the block CFL.
Definition physical_model.hpp:131
Old name (compat): HyperbolicPhysicalModel used to be HyperbolicModel.
Definition physical_model.hpp:196
Hyperbolic brick of a model: flux + wave speed + variables + cons<->prim conversions.
Definition physical_model.hpp:181
Minimal contract of a physical model.
Definition physical_model.hpp:81
Definition amr_hierarchy.hpp:29
constexpr int kAuxBaseComps
Definition state.hpp:147
POPS_HD constexpr int aux_comps()
Width of the aux channel a model CONSUMES.
Definition physical_model.hpp:67
Pointwise types of the physics layer: StateVec<N> (conserved state) and Aux (auxiliary fields from th...
Base scalar types and the POPS_HD macro (host+device portability).
#define POPS_HD
Definition types.hpp:25
Descriptor of a model's variables (Vars).