Geometry interface decisions before CADET-Core v6.0

This post is about interface choices that should be settled before the CADET-Core v6.0 release.
The goal is not to prescribe a full internal rewrite.
It is to decide how users and downstream tooling should describe column geometry, dimensionality, and grids before those names and structures become harder to change.

Status quo

The unified column model got the big structural decision right: particle physics is now separated from the column.
ColumnModel1D orchestrates and owns almost no transport numerics.
The particle is a runtime IParticleModel object selected by a factory, its intra-particle discretization is a nested operator inside it (ParticleDiffusionOperatorFV/DG), and the bulk transport operator is a template parameter on the column.
Adding a particle variant no longer means a new column class.
That is a real win and worth keeping intact.

This post looks at the one axis that did not come out as cleanly, and which two open issues are now forcing a decision on: geometry (#510 for the interface, #697 for frustum DG).
My main proposal is to treat this as a release-interface decision first.
We can agree on the input structure and output metadata for v6.0 without deciding, in the same step, how much transport-operator code should be merged internally.

The axes of variation

A column model varies along several independent axes:

  1. Physical geometry: axial cylinder, radial (cylindrical wedge), frustum (spherical conical frustum), …
  2. Dimensionality: how many directions are resolved, as a nested sequence of symmetry assumptions: 3D (fully resolved) to 2D (angular symmetry) to 1D (+ radial symmetry) to 0D (+ axial symmetry, i.e. a CSTR). The physical symmetry is the product of this and the geometry: dimensionality says how far down this sequence you are, geometry says which directions those are and their metric, which is why radial-2D and axial-2D differ.
  3. Discretization method: FV vs DG.
  4. Grid: number of cells/elements, uniform vs non-uniform, optional user-defined faces.
  5. Particle model: equilibrium / homogeneous / general-rate.

These are genuinely orthogonal, so the useful interface question is not “how many classes do we need” but which choices should users express separately.
Each independent user choice should have its own field instead of being hidden inside a unit name or an implementation variant.
A geometry choice should not be communicated by selecting a different unit type; a grid choice should not double as the definition of the physical domain.
For the implementation, separate classes or specialized code paths should be reserved for choices that really change the numerical structure.
Geometry is different from FV/DG in this respect: it changes the shape of the domain, while FV versus DG changes the discretization method.
That distinction is enough to guide the public interface; it does not require every user to know the operator internals.

The restructure already handles most axes this way.
The particle model is selected from the configuration as a runtime object.
Discretization is allowed to specialize internally, which is correct: FV and DG differ in sparsity, basis, and Jacobian structure, so they really are different code.
Geometry is the one that ended up encoded more tightly than the user-facing structure requires.

Where geometry is encoded too tightly

Geometry is currently encoded in two places.

In the transport-operator type.
Geometry is fused with discretization into the single ConvDispOperator template parameter, so the operator set is Axial/Radial x FV/DG, and every new geometry needs a corresponding place in every relevant discretization path.
(Frustum sits outside this today, only on the legacy wrapper path. That is understandable given it is still experimental, and the DG scheme is about to be worked out in [#697]; the question is where it should land when it matures.)

In the unit-type name.
Geometry is part of the unit identity: COLUMN_MODEL_1D versus RADIAL_COLUMN_MODEL_1D.
So a user selects geometry by choosing a differently-named unit, and each new geometry means a new name.
This is what #510 sets out to fix.

The interface conclusion is the same in both cases: users should select geometry as geometry, not by selecting a different unit type.
Internally, each discretization can consume the resulting geometry data in the form it needs.

Why a single geometry source helps

This is not only about avoiding extra unit names.
Several parts of CADET (not just Core) need the same physical map from solution entries to space: setting initial conditions, evaluating space-dependent parameters, reading out solution fields, and exporting or plotting results.
All of them need to know which axes are resolved, what the coordinates mean, and what geometry weights or areas belong to those coordinates.

If each path reconstructs that information from a mix of unit type, dimensionality, grid faces, and geometry-specific special cases, radial, axial, frustum, 0D, 1D, and 2D models can drift apart in small but user-visible ways.
A single geometry object gives those paths one shared source of truth.
It can answer common questions such as: what are the resolved axes, where are the cell or element coordinates, how do relative grid positions map to physical coordinates, and which volume, area, or radial weight belongs to each location?

That is useful even if the numerical operators remain partly separate.
The solver can still use the discretization-specific implementation it needs, while setup, parameter evaluation, output metadata, and plotting all see the same geometry description.

Geometry parameters versus grid parameters

One interface detail is worth separating from the operator discussion.
Geometry should describe the physical domain: length, radius or diameter, inner and outer radii, frustum parameters, and so on.
Discretization should describe how that domain is sampled: method, NCOL/NELEM, polynomial degree, reconstruction, and optional non-uniform grid faces.

This matters for user-defined grids.
If the user provides physical face coordinates, then the geometry parameters and the coordinates both describe the same extent.
That raises avoidable questions: should COL_LENGTH be derived from the last grid face, or should the grid be checked against COL_LENGTH?
The current column implementation takes the second approach.

My recommendation would be to keep physical geometry as the source of truth and let user-defined grids be relative distributions on [0, 1], which are then mapped to the physical domain by the geometry.
This is already the particle precedent: PAR_DISC_VECTOR is relative and is mapped to [PAR_CORERADIUS, PAR_RADIUS].
For columns, the same idea could coexist with the current physical GRID_FACES behavior for backward compatibility.
So this does not have to be a hard requirement; it is mainly a cleaner target interface.

Possible staged proposal

One possible path is to split the work into two steps.
The first is the v6.0 interface decision and can proceed without changing the numerics; the second is a localized refactor of the operator layer if we decide it is worth doing.

Step 1 - parse geometry and dimensionality as fields before v6.0 ([#510]).
Replace the compound unit name with the three fields already proposed in [#510]: unit_type = column, geometry = axial | radial | frustum, dimensionality = 0D | 1D | 2D.
Keep the physical geometry parameters with the geometry and the numerical grid parameters with the discretization.
The geometry object can then derive the metric data used by numerical operators.
It can also publish a resolved-axis descriptor (axis role, coordinates, and metric), giving visualization and output tooling an explicit map of the field structure without adding a redundant symmetry input field.
This is meant as an interface layer over the current code first, not as a demand for an immediate rewrite of all operators.

Step 2 - pass geometry-derived metric data instead of selecting operators by geometry.
Where geometry only changes fixed metric values such as cross sections, face areas, or radial weights, those values can be passed into the operator rather than choosing a different operator solely because geometry changed.
For FV this could collapse Axial/Radial (and eventually frustum) into a single operator that receives those metric values, and gives the existing experimental frustum-FV work a clear path into the unified model.
For DG the metric enters the scheme itself, so how far the operators merge is exactly what the frustum-DG work ([#697]) will clarify; Step 2 does not add frustum-DG, it just avoids growing the Axial/Radial/Frustum x FV/DG list of operator variants while that work happens.
Crucially, Step 1 does not depend on any of this and can land first.

The same split should probably apply to particles.
PAR_GEOM = SPHERE | CYLINDER | SLAB is already a particle-geometry field in spirit, while NCELLS, PAR_DISC_TYPE, and PAR_DISC_VECTOR are particle-discretization fields.
If the column interface is cleaned up, it would be natural to make this explicit for particles as well: particle geometry owns radius/core radius and the geometry factor, particle discretization owns the 1D radial grid.

What to leave alone

Exposing dimensionality as a field does not require collapsing the internal ColumnModel1D / ColumnModel2D class split.
That split is justified by discretization structure, the 2D DG tensor-product operator is genuinely different, not by geometry.
Keep the interface field and the internal class boundary separate: dimensionality selects the class family, geometry supplies metric data within it.

Let’s discuss the field names and allowed inputs

  • UNIT_TYPE → Column_Model, CSTR. Further, last time we agreed to also support the original LRM LRMP GRM names (as a feature, not for bwd compatibility) and cadet throws an error if the particle is not configured accordingly. What about crystallization in a dpfr, then its not a “COLUMN” or “COLUMN_MODEL”?
  • DIMENSIONALITY or SPATIAL_DIMENSIONS or SPATIAL_RESOLUTION ? Here, we have the option to only input 1D, 2D, 3D or have a list like [“AXIAL”, “RADIAL”] ? What about 0D or should the CSTR be a separate unit_type?
  • GEOMETRY or COLUMN_GEOMETRY ? → AXIAL_FLOW_COLUMN, RADIAL_FLOW_COLUMN, FRUSTUM_COLUMN ? or wiothout the column suffix?

I agree to normalize the user-supplied discretization grid to [0,1].

Another discussion we should pick up here again is to rename col_length to bed_length to avoid ambiguities with the different geometries, specifically radial flow

Another discussion we should pick up here again is to rename col_length to bed_length to avoid ambiguities with the different geometries, specifically radial flow

I second that!

Technically, there are two meaningful implementations of 2D flow. Starting from the 3D formulation, one can eliminate either the angular coordinate or the column radius. While eliminating the angular coordinate is more intuitive for a cylindrical column or frustum, both choices are meaningful for the wedge geometry. Any ideas on how best to address this?

Let the user specify “axial” “radial” instead of “2D” ?

For completeness, we should include both 2D options (also in the companion modeling manuscript). We can simply call them 2Da when the angle is removed and 2Db when the other coordinate (column radius in axial flow, column height in radial flow) is removed (first). Since 2Da is the standard case, we can define 2D = 2Da unless specified otherwise.

Giving it a second thought, the alternative 2D model for a cylindrical column or frustum does not make much sense. Eliminating the radius while keeping the angular coordinate would lead to a singularity at the center, where a single point would have different concentrations depending on the angle.

However, for the wedge geometry, there are two different and meaningful 2D versions of the model.

This post pictures the current state of the interface discussion and is always updated once we make a decision. Marked in yellow are the decisions to be made.

  • UNIT_TYPE
    • what should this field do?
      • a.1) remove this field
      • a.2) optional feature field for established units, eg lrm, lrmp, grm, cstr. The feature is that cadet-core checks whether the input fields and parameters are feasible for that model.
      • b) distinguish between chromatography vs non-chromatography and potential future unit types.
      • d) other suggestions?
  • SPATIAL_RESOLUTION
    • What input should we have
      • a single string: 0D or 1D or 2D
      • a list of strings with coordinate direction names: [“axial”], [“radial”], [“axial”, “radial”] ?
    • Field name decision (will be made after the decision above)
      • DIMENSIONALITY
      • SPATIAL_DIMENSIONS
      • SPATIAL_RESOLUTION?
      • other suggestions?
  • GEOMETRY
    • field does not need to be specified / is ignored for 0D
    • Name decision
      • GEOMETRY
      • other suggestions?
    • Input names
      • axial_flow_cylinder, radial_flow_cylinder_annular, axial_flow_frustum, radial_flow_cylinder_shell_wedge
      • other suggestions
  • Geometry related input
    • we want a minimal amount of unified input fields from which the rest is computed
    • what geometric parameter inputs should we have?
      • a) always specify cross_section_area, bed_length and, if needed, further properties
        • axial: cross_section_area and bed_length
        • radial: cross_section_area_inner or outer, bed_length and col_height or a different name for that
        • radial wedge → same as radial plus circle_fraction in (0,1)
        • frustum → cross_section_area_large_end and cross_section_area_small_end or different name, and bed_length
      • other suggestions?
      • What to call new field FLOW_DIRECTION (-1, 1) or FORWARD_FLOW (bool) ? Default value or make it explicit?
  • Decided:
    • User supplied discretization grids are normalized to [0, 1]
    • rename col_length to bed_length
    • remove velocity field, add field FORWARD_FLOW and derive velocity from cross_section_area