Open boundary conditions using SplitExplicitFreeSurface#5351
Open boundary conditions using SplitExplicitFreeSurface#5351simone-silvestri wants to merge 70 commits intomainfrom
SplitExplicitFreeSurface#5351Conversation
@simone-silvestri that was fast! I will see if I can get to this at the end of the week, or @rafmudaf may be able to test first with our current Med config. |
| BoundaryCondition, | ||
| FluxBoundaryCondition, ValueBoundaryCondition, GradientBoundaryCondition, OpenBoundaryCondition, | ||
| PerturbationAdvection, | ||
| PerturbationAdvection, Flather, Radiation, |
There was a problem hiding this comment.
What's the difference between these three? Isn't Radiation a subset of PerturbationAdvection?
There was a problem hiding this comment.
I think it's the opposite. PerturbationAdvection requires specifying a phase velocity at the boundary while radiation computes it (from what I understood of perturbation advection)
There was a problem hiding this comment.
I think you could make the phase velocity computation the condition and then use perturbation advection?
There was a problem hiding this comment.
Heh, I mean in either case --- whether Radiation is a subset of PerturbationAdvection, or whether PerturbationAdvection is a subset of Radiation --- in either of those cases, we should not introduce a new type but instead reuse the same type to express both schemes.
There was a problem hiding this comment.
we should not introduce a new type but instead reuse the same type to express both schemes.
I'm not sure. Radiation schemes can differ widely in their internals. Having different scheme types for each one may be useful just for the sake of simplicity. That said, the name Radiation is pretty generic and I think most people would call PertubationAdvection, Orlanski, the oblique scheme, and even Flather, radiation schemes.
@simone-silvestri from talking to you at the conference it seems like what you're implementing here for the baroclinic flow is an Orkanski scheme, no?
There was a problem hiding this comment.
I'd say just replace the existing periodic one with this. We don't need both.
…ith `SplitExplicitFreeSurface` (#5378) * Validate any bc on nothing-dimension fields * Add U and V to field location map * Fix missing u, v reference * Add missing b_bg helper function for sponge layer * Maintain barotropic velocities for free surface Co-authored-by: Simone Silvestri <[email protected]> * Add sea-surface height to field location map Co-authored-by: Simone Silvestri <[email protected]> * Remove extra validation dispatch --------- Co-authored-by: Simone Silvestri <[email protected]>
…ananigans.jl into ss/prepare-for-open-boundaries
…ananigans.jl into ss/prepare-for-open-boundaries
XLA folds Δt + zero(Δt) → Δt, so the aliasing between last_Δt and last_stage_Δt persists. Use Reactant.Ops.optimization_barrier which XLA cannot fold, ensuring distinct buffers for the while loop carry. Co-Authored-By: Claude Opus 4.6 (1M context) <[email protected]>
Add `materialize_clock!(clock, timestepper)` called in HFSM and NH model constructors. For QuasiAdamsBashforth2TimeStepper, this sets `clock.last_Δt = clock.last_stage_Δt` to ensure they are distinct objects. This is needed for Reactant, where aliased ConcreteRNumber fields cause XLA buffer donation errors in compiled loops. With materialize_clock! breaking the alias at construction time, the optimization_barrier workaround in the Reactant tick! override is no longer needed and is removed. Co-Authored-By: Claude Opus 4.6 (1M context) <[email protected]>
materialize_clock! in the constructor doesn't prevent MLIR-level aliasing that occurs at runtime inside tick! when both last_Δt and last_stage_Δt are assigned the same Δt value. The optimization_barrier is still needed until Reactant provides a proper buffer copy op. See #5389 (comment) Co-Authored-By: Claude Opus 4.6 (1M context) <[email protected]>
materialize_clock! in the model constructor ensures last_Δt and last_stage_Δt are distinct ConcreteRNumber objects. Assigning the same value to both in tick! via .mlir_data does not re-alias the objects. Co-Authored-By: Claude Opus 4.6 (1M context) <[email protected]>
For QAB2, last_Δt and last_stage_Δt always hold the same value after tick!. In the Reactant extension, materialize_clock! aliases them via setfield! (bypassing the ReactantClock setproperty! override) so that Reactant's tracer sees one buffer, avoiding XLA buffer donation errors. The src/ QAB2 method is removed since aliasing is only needed for Reactant's ConcreteRNumber fields (normal Float64 fields are immutable values with no aliasing concern). Co-Authored-By: Claude Opus 4.6 (1M context) <[email protected]>
…ananigans.jl into ss/prepare-for-open-boundaries
…tion-update-state
…liMA/Oceananigans.jl into glw/remove-initialization-update-state
…-state' into ss/prepare-for-open-boundaries
…-state' into ss/prepare-for-open-boundaries
…-state' into ss/prepare-for-open-boundaries
…into ss/open-boundary-conditions
This PR implements open boundary conditions for the hydrostatic free surface model compatible with the split explicit free surface. The boundary conditions follow what discussed in #5229, meaning
Flatherboundary conditions for the barotropic component and "Orlanski" (Radiation) boundary conditions for the three-dimensional evolution.I need to clean up a bit the PR (for example Flather now stores the external values inside its own struct, but it would be probably better to store it in the
conditionof theboundary_conditionstruct.Since the flow of the algorithm has to change, because the barotropic correction needs to be done after filling up the boundary conditions for the 3d velocities, such that the barotropic components of the boundaries is consistent with the 2d evolution, some refactor is necessary and I would like to open a PR specifically for this (this is actually a positive because it means we can compute the entire
wvelocity in theupdate_statein distributed computations).There is an example that showcases the use of these boundary conditions which we can remove later on.
A part from all the changes that need to be done to clean up the code, this feature should be already fully functional and testable, so if you want to test it @vtamsitt you can give it a try
cc @tomchor