WIP: Add flexibility to nonconservative BCs

[MarcoArtiano]

This PR allows to have more flexibility when defining the boundary conditions for a set of equations with nonconservative terms. This issues has been discussed in #2175. Now the user can define his own boundary condition, considering that the function accepts both surface flux and the nonconservative flux. An example can be found in elixir_mhd_reflective_wall

[github-actions]

Review checklist

This checklist is meant to assist creators of PRs (to let them know what reviewers will typically look for) and reviewers (to guide them in a structured review process). Items do not need to be checked explicitly for a PR to be eligible for merging.

Purpose and scope

• The PR has a single goal that is clear from the PR title and/or description.
• All code changes represent a single set of modifications that logically belong together.
• No more than 500 lines of code are changed or there is no obvious way to split the PR into multiple PRs.

Code quality

• The code can be understood easily.
• Newly introduced names for variables etc. are self-descriptive and consistent with existing naming conventions.
• There are no redundancies that can be removed by simple modularization/refactoring.
• There are no leftover debug statements or commented code sections.
• The code adheres to our conventions and style guide, and to the Julia guidelines.

Documentation

• New functions and types are documented with a docstring or top-level comment.
• Relevant publications are referenced in docstrings (see example for formatting).
• Inline comments are used to document longer or unusual code sections.
• Comments describe intent ("why?") and not just functionality ("what?").
• If the PR introduces a significant change or new feature, it is documented in NEWS.md with its PR number.

Testing

• The PR passes all tests.
• New or modified lines of code are covered by tests.
• New or modified tests run in less then 10 seconds.

Performance

• There are no type instabilities or memory allocations in performance-critical parts.
• If the PR intent is to improve performance, before/after time measurements are posted in the PR.

Verification

• The correctness of the code was verified using appropriate tests.
• If new equations/methods are added, a convergence test has been run and the results
  are posted in the PR.

Created with ❤️ by the Trixi.jl community.

[codecov]

Codecov Report

Attention: Patch coverage is 93.02326% with 3 lines in your changes missing coverage. Please review.

Project coverage is 89.92%. Comparing base (60d5c32) to head (c2b44ae).
Report is 16 commits behind head on main.

Files with missing lines	Patch %	Lines
src/basic_types.jl	50.00%	3 Missing ⚠️

❗ There is a different number of reports uploaded between BASE (60d5c32) and HEAD (c2b44ae). Click for more details.

HEAD has 1 upload less than BASE

Flag	BASE (60d5c32)	HEAD (c2b44ae)
unittests	26	25

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #2200      +/-   ##
==========================================
- Coverage   96.15%   89.92%   -6.23%     
==========================================
  Files         483      483              
  Lines       38353    38364      +11     
==========================================
- Hits        36876    34497    -2379     
- Misses       1477     3867    +2390     

Flag	Coverage Δ
unittests	89.92% <93.02%> (-6.23%)	⬇️

Flags with carried forward coverage won't be shown. Click here to find out more.

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

[MarcoArtiano]

According to the way this feature has been implemented, the functions boundary_condition_slip_wall in shallow_water_wet_dry_1d.jl and shallow_water_wet_dry_2d.jl should have been modified accordingly.

CC @andrewwinters5000

[andrewwinters5000]

Thanks for tackling this @MarcoArtiano! The functionality works as discussed where it is on the user to put in the boundary conditions for the conservative and nonconservative terms and then combine them appropriately, i.e., scaling noncons with 0.5f0. This is the only aspect I am not a huge fan of, as it is related to how the nonconservative terms are implemented at interfaces and has nothing to do with boundary conditions.

Because of this, we should maybe document where this factor of 0.5f0 magically comes from in the new boundary condition routines (as it was documented beofre inside the compute boundary flux call).

What do you think @patrickersing ?

[andrewwinters5000]

This comment would possibly need added to where the factor of 0.5 scaling occurs.

[patrickersing]

I am also not very content that we have to expose this to the user in order to specify the boundary condition, but I don't see a better way to accomplish this. I definitely think that we should comment this to give some explanation where the factor 0.5 comes from and that it is not specific to any boundary condition.

[MarcoArtiano]

I agree with all of you. I'm also not a huge fan of the factor 0.5. We may alternatively return a tuple from the boundary conditions for non-conservative systems and multiply by 0.5 inside the solver. As an example:

flux_, noncons_ = boundary_condition(u_inner, normal_direction, x, t, surface_integral.surface_flux, equations)

  # Copy flux to element storage in the correct orientation
  for v in eachvariable(equations)
      # Note the factor 0.5 necessary for the nonconservative fluxes based on
      # the interpretation of global SBP operators coupled discontinuously via
      # central fluxes/SATs
      surface_flux_values[v, node_index, direction_index, element_index] = flux_[v] + 0.5f0 * noncons_[v]
  end

What do you think?

[patrickersing]

I like this idea, since it doesn't require knowledge about the implementation aspect of adding the 0.5 to create the boundary condition. It should work well for the BCs that are currently implemented, but I am not sure what this would look like for the new type of boundary condition that you plan to add. Would you then only set the flux_ part and set the noncons_ part to zero?
I think it would be helpful to include a specific example with such a new BC in the PR to better understand and evaluate what the implementation should look like.

[andrewwinters5000]

I also like this idea with the documentation of how it would work in practice from Patrick's comment above. That is, for something like the standard shallow water equations the jump in the bottom topography is zero at the physical boundary (typically) so one's new boundary condition could compute the conservative flux pieces and then have

   return flux_, SVector(0,0,0,0) # flux, noncons