WIP: Well-balanced mortars

[andrewwinters5000]

This seeks to rewrite the mortar projection routines specifically for the shallow water equations to preserve well-balancedness with non-conforming approximations. That is, the shallow water solver could utilize the AMR capabilities of TreeMesh or P4estMesh. The implementation might become delicate at wet/dry transition regions as one needs to take care in what solution values and/or fluxes are projected to and from the mortar while maintaining positivity.

• Decide on file names and structure in the src/solvers folder
• Implement and test on TreeMesh in 2D
• Implement and test on P4estMesh in 2D

[codecov]

Codecov Report

Attention: Patch coverage is 0% with 351 lines in your changes missing coverage. Please review.

Project coverage is 87.62%. Comparing base (cbc57dc) to head (c33375c).
Report is 2 commits behind head on main.

Files	Patch %	Lines
src/solvers/scratch_p4est.jl	0.00%	147 Missing ⚠️
src/solvers/scratch_tree.jl	0.00%	110 Missing ⚠️
..._dgsem/elixir_shallowwater_perturbation_wet_dry.jl	0.00%	29 Missing ⚠️
...dgsem/elixir_shallowwater_well_balanced_wet_dry.jl	0.00%	26 Missing ⚠️
...4est_2d_dgsem/elixir_shallowwater_well_balanced.jl	0.00%	25 Missing ⚠️
...p4est_2d_dgsem/elixir_shallowwater_perturbation.jl	0.00%	14 Missing ⚠️

Additional details and impacted files

@@             Coverage Diff             @@
##             main      #45       +/-   ##
===========================================
- Coverage   99.27%   87.62%   -11.66%     
===========================================
  Files          61       67        +6     
  Lines        2638     2989      +351     
===========================================
  Hits         2619     2619               
- Misses         19      370      +351     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[andrewwinters5000]

The current version provides a well-balanced approximation for a fully wet configuration. The examples/p4est_2d_dgsem/elixir_shallowwater_well_balanced_wet_dry.jl crashes in the first time step. So this requires further investigation. The main issues with the fully wet version are:

1. There are still significant allocations in the computation of the mortar fluxes.
2. More pressing is that we require a slight modification to the P4estMortarContainer within Trixi.jl to add additional storage of the unprojected large element solution in order to have this information available to subtract off a superfluous contribution of the local physical flux in the large element. I want to avoid this because this specailized mortar treatment should only need to live within the new package.

The modification below in Trixi's src/solvers/dgsem_p4est/containers.jl is slight, but incredibly hacky. I added a third "slot" to store the large element information that I needed as given in the code snippet below. There is probably a better way to precompute and store the parent flux contributions on each mortar. I am just not familiar enough with the P4estMesh implementations and its various caches to see where one could store such information.

# Create mortar container and initialize mortar data.
function init_mortars(mesh::Union{P4estMesh, T8codeMesh}, equations, basis, elements)
    NDIMS = ndims(elements)
    uEltype = eltype(elements)

    # Initialize container
    n_mortars = count_required_surfaces(mesh).mortars

    _u = Vector{uEltype}(undef,
                         3 * nvariables(equations) * 2^(NDIMS - 1) *
                         nnodes(basis)^(NDIMS - 1) * n_mortars)
    u = unsafe_wrap(Array, pointer(_u),
                    (3, nvariables(equations), 2^(NDIMS - 1),
                     ntuple(_ -> nnodes(basis), NDIMS - 1)..., n_mortars))

    _neighbor_ids = Vector{Int}(undef, (2^(NDIMS - 1) + 1) * n_mortars)
    neighbor_ids = unsafe_wrap(Array, pointer(_neighbor_ids),
                               (2^(NDIMS - 1) + 1, n_mortars))

    _node_indices = Vector{NTuple{NDIMS, Symbol}}(undef, 2 * n_mortars)
    node_indices = unsafe_wrap(Array, pointer(_node_indices), (2, n_mortars))

    mortars = P4estMortarContainer{NDIMS, uEltype, NDIMS + 1, NDIMS + 3}(u,
                                                                         neighbor_ids,
                                                                         node_indices,
                                                                         _u,
                                                                         _neighbor_ids,
                                                                         _node_indices)

    if n_mortars > 0
        init_mortars!(mortars, mesh)
    end

    return mortars
end