Repository navigation
PLIC: mass and momentum fluxes consistent with the volume fluxes - #256
Draft
ThomasHowarth wants to merge 7 commits into
Draft
ThomasHowarth wants to merge 7 commits into
ThomasHowarth wants to merge 7 commits into
Conversation
Add a PLIC directory, the first part of geometric volume-of-fluid transport with a piecewise-linear interface calculation. The interface in each mixed cell is a plane, fixed by the volume fraction of the cell through the analytic relations of Scardovelli & Zaleski (JCP 164, 2000). - hydro_plic_K.H: the volume <-> plane relations and the plane fractions used to set up problems. - hydro_plic_normals_K.H: the Youngs, Mixed Youngs-Centered and height-function normals, each a struct passed as a template parameter, as they run inside device kernels. - PLIC::FitPlanes fits the plane of each mixed cell, with the reconstruction selected by name (youngs, myc or myc_height), as advection_type is in HydroUtils. PLIC::nGrowF gives the ghost cells it reads, and PLIC::CheckMesh rejects non-Cartesian coordinates and embedded boundaries. The geometric kernels derive from Kynema-SGF (formerly AMR-Wind), under its BSD 3-Clause License; the headers keep its copyright notice. Tests/PLIC checks to round-off that the volume <-> plane relations invert each other, that each reconstruction's plane holds the volume of its cell, that the myc_height normal is exact for a planar interface, and that FitPlanes on a level of several boxes fits the mixed cells, valid and in the first ghost layer. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Add the two operators of the PLIC library that couple the interface to the solver's boundaries and projection. - PLIC::FillWallGhosts extends the interior interface through the domain faces where the volume fraction is extrapolated: a free contact line at a wall, a smooth continuation at an outflow. The faces are selected by the BCRec of the volume fraction, as for any other scalar: foextrap and hoextrap faces are extended, ext_dir (inflow) and reflect_even (symmetry) ghost cells are kept, and reflect_odd aborts. A ghost value depends on cells further along the wall than one box holds, so it works on a whole level, and its result does not depend on the box layout. - PLIC::SharpFaceDensity gives the face density of a balanced-force projection (Francois et al., JCP 213, 2006): for a planar interface normal to gravity, rho_face g is a discrete gradient. Its inverse is the face coefficient of Hydro::MacProjector. Tests/PLIC checks to round-off that the wall ghost cells continue a planar interface, that ext_dir and reflect_even ghost cells are left as filled, and that the curl of rho_face g vanishes for a planar interface. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Add PLIC::ComputeEdgeStateAndFluxes, the fluxes and edge states of the
volume fraction on the faces of one box, from the unsplit geometric fluxes
of Owkes & Desjardins (JCP 270, 2014). It is called in the solver's own
MFIter loop, as Godunov::ComputeEdgeState is, and its fluxes follow the
convention of HydroUtils::ComputeFluxesOnBoxFromState: edge state times
face velocity, optionally times the face area, for an update
F^{n+1} = F^n - dt div(flux).
The flux through a face is the liquid inside its flux volume, bounded by
the face and its corners traced back along the flow, divided by dt. The
update conserves the liquid volume to round-off and keeps F in [0,1] up to
round-off, without clipping. The edge state is the liquid fraction of the
flux volume.
Each flux depends only on data around its face: the nodes are traced from
the face velocities within one ghost layer, and the flux volume is cut in
coordinates relative to the face. Boxes sharing a face, and periodic
images, get the same flux to the bit, so no synchronization of the faces
is needed and the solver's usual average_down_faces or FluxRegister
handles the levels.
Tests/PLIC becomes a small AmrCore solver that advances F with it in its
own loop: the single vortex, Zalesak's disk, a contact line in the vortex
and LeVeque's 3D deformation. The volume is conserved to round-off per step
and the results do not depend on the box layout or the number of ranks.
L1 shape errors on uniform grids of 64, 128 and 256 cells per direction:
- single vortex: 1.87e-3, 4.68e-4, 1.39e-4
- Zalesak's disk: 2.96e-3, 1.03e-3, 4.68e-4
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01TwkdHV9vdS67quGfmN9xjG
HydroUtils::ComputeFluxesOnBoxFromState advects a volume fraction with advection_type = "PLIC", as any other scalar: the edge states come from PLIC::ComputeEdgeStateAndFluxes with the myc_height reconstruction, and HydroUtils::ComputeFluxes forms the fluxes from them, so u_flux, knownFaceState, area weighting and ComputeConvectiveTerm work unchanged. The boundary conditions are carried by the ghost cells of F, filled with PLIC::FillWallGhosts. "PLIC" aborts with embedded boundaries, for velocities and for more than one component, and ExtrapVelToFaces does not accept it. PLIC::ComputeEdgeStateAndFluxes accepts empty flux arrays, for the edge states only. Tests/PLIC advects through HydroUtils with test.use_hydroutils = 1, with the same shape errors as the direct path and the volume conserved to round-off. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01TwkdHV9vdS67quGfmN9xjG
Add PLIC::ComputeMassMomentumFluxes, per box like ComputeEdgeStateAndFluxes, for momentum advected in conservative form with a large density ratio (Rudman, 1998). It gives the mass and momentum fluxes of the same flux volumes: near the interface each piece carries the density of its phase and the velocity of the cell it lies in, and away from it the flux carries a face velocity given by the solver. The density from the mass fluxes is then the mixture density of the new F to round-off, so that a uniform velocity stays uniform across the interface. Add PLIC::mixture. The setup of the interface band, planes and node displacements is shared by both functions, so that they trace the nodes the same way. Tests/PLIC/inputs.momentum: a drop in a uniform and in a vortex flow with a density ratio of 1000, on one periodic level of several boxes. The uniform velocity stays exactly uniform, the density matches the mixture density of F to round-off and the momentum is conserved. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01TwkdHV9vdS67quGfmN9xjG
The cuts of the flux volumes are most of the cost of the fluxes, and ComputeEdgeStateAndFluxes and ComputeMassMomentumFluxes each made their own. A solver that advances the volume fraction before it has the face velocity of its momentum scheme had to call both, and cut twice. - An overload of PLIC::ComputeEdgeStateAndFluxes also gives the mass fluxes, and the momentum fluxes near the interface, from the same cuts. The momentum flux of the faces away from the interface, which carries the solver's face velocity, is left unset. - PLIC::FarFieldMomentumFluxes sets it later, as mflux * vel_face. It finds those faces again from F, as those with neither neighbouring cell in the interface band, so no mask is kept between the two calls. The fluxes are then those of the separate calls, to the bit. All the flux functions share one implementation, with the mass and momentum fluxes compiled in only where they are asked for. Tests/PLIC/inputs.momentum checks at every step that the fluxes from one set of cuts equal those of the separate calls, to the bit, in 2D and 3D, on one and four ranks. In the PeleLMeX dam break the cuts take 0.26 s instead of 0.51 s. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
nvcc rejects an extended __device__ lambda that first captures a variable inside an if constexpr branch, so the GPU build of hydro_plic_owkes_desjardins.cpp failed. The flux kernel now names the variables of its MassMomentum branches once at the top, which changes nothing on the CPU. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This branch has not been deployed
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
This PR adds mass and momentum fluxes for momentum advected in conservative form at a large density ratio (Rudman, JCP 1998). They come from the same flux volumes as the volume fraction. With rho^n = PLIC::mixture(F^n, rho_l, rho_g), the rho^{n+1} from the mass fluxes is the mixture density of F^{n+1} to round-off, so a uniform velocity stays uniform across the interface whatever the density ratio.
New API
The single-cut path matches the separate calls to the bit. All the flux functions share one implementation, and the mass and momentum terms are compiled in only where they're asked for.
Tests
Tests/PLIC/inputs.momentum runs a drop in a uniform flow and in a vortex flow, with a density ratio of 1000. It checks, to round-off, that the density matches the mixture density of the new F, that a uniform velocity stays uniform and that momentum is conserved. It also checks that the single-cut fluxes equal those of the separate calls, to the bit. Docs/source/PLIC.rst and Tests/PLIC/README.md are updated.
🤖 Generated with Claude Code (https://claude.com/claude-code)