First step towards FD methods and surface integrals/terms #617
Conversation
This could be considered to be breaking since it changes a documented constructor of our `DG` type. On the other hand, the documentation of this constructor was just wrong (since it contained the basis twice but not the mortar). Hence, this is more like a bug fix, I think ;).
Codecov Report
@@ Coverage Diff @@
## main #617 +/- ##
==========================================
- Coverage 93.36% 93.32% -0.05%
==========================================
Files 155 156 +1
Lines 15353 15467 +114
==========================================
+ Hits 14334 14434 +100
- Misses 1019 1033 +14
Flags with carried forward coverage won't be shown. Click here to find out more.
Continue to review full report at Codecov.
|
src/solvers/dg_tree/dg.jl
Outdated
| DGSEM(; RealT=Float64, polydeg::Integer, | ||
| surface_flux=flux_central, | ||
| volume_integral::AbstractVolumeIntegral=VolumeIntegralWeakForm(), | ||
| mortar=MortarL2(basis)) |
There was a problem hiding this comment.
This docstring is redundant with the lines above?
There was a problem hiding this comment.
Also, shouldnt it be
| DGSEM(; RealT=Float64, polydeg::Integer, | |
| surface_flux=flux_central, | |
| volume_integral::AbstractVolumeIntegral=VolumeIntegralWeakForm(), | |
| mortar=MortarL2(basis)) | |
| DGSEM(; RealT=Float64, polydeg::Integer, | |
| surface_integral=SurfaceIntegralWeakForm(), | |
| volume_integral::AbstractVolumeIntegral=VolumeIntegralWeakForm(), | |
| mortar=MortarL2(basis)) |
for consistency? Also, why the dispatch on AbstractVolumeIntegral?
There was a problem hiding this comment.
I removed the restriction to AbstractVolumeIntegral and the other docstring in 05801bf.
DGSEM uses the SurfaceIntegralWeakForm by definition. We could allow more general surface terms, but we don't have a good reason to do so right now. We need to keep the keyword argument surface_flux for backward compatibility.
There was a problem hiding this comment.
Besides the weak/strong argument, an additional reason I see for the abstraction with the SurfaceIntegral type is to pave the way for alternative surface terms such as SAT (simultaneous approximation terms) that sometimes do not follow our "numerical flux" structure.
There was a problem hiding this comment.
Yeah, that's something I had in mind, too.
ranocha
changed the title
|
Initial draft of SBP FD methods which needs modifications of SummatioByPartsOperators.jl that I will push later to upstream. Convergence test for linear advection: julia> using Trixi
julia> convergence_test("examples/tree_2d_fd/elixir_advection_extended.jl", 5, initial_refinement_level=0, abstol=1.0e-10, reltol=1.0e-10, N=32)
[...]
####################################################################################################
l2
scalar
error EOC
7.68e-05 -
1.33e-05 2.53
1.84e-06 2.86
2.29e-07 3.00
2.86e-08 3.00
mean 2.85
----------------------------------------------------------------------------------------------------
linf
scalar
error EOC
2.71e-04 -
3.59e-05 2.92
5.22e-06 2.78
6.51e-07 3.00
8.67e-08 2.91
mean 2.90
----------------------------------------------------------------------------------------------------
Dict{Symbol, Any} with 3 entries:
:variables => ("scalar",)
:l2 => [2.84821]
:linf => [2.90261] |
|
@sloede @andrewwinters5000 @gregorgassner: I implemented a first draft of FD methods. I would like to get feedback before polishing this since I have no motivation to work on this if it will be rejected anyway since people have some strong opinions on coding/naming conventions (cf. #616). |
|
This looks very nice to me. I also looked through the changed files, looks all good to me so far. (i) Did you test the multi-element? Is your EOC test done by increasing the SBP nodes inside the element, or by increasing the number of elements? (ii) General question: Would it make sense to have Strong Form (Volume, Surface) also for the DGSEM operator? (iii) The current volume integral implementation is with mul!() as far as I can tell. It seems to me, that this uses the full (non-sparse) D matrix to do the multiplication. I assume, it will be relative straight forward in future, to specify the volume integral further with the specific optimized implementations for a given SBP D matrix? |
|
Thanks for the encouraging feedback 👍
I tested multiple elements connected via strong form surface terms/integral. The convergence test increases the number of elements since our
It's definitely an option. However, we don't have a use case for it right now and the weak form should be more efficient. Thus, I don't object implementing it for DGSEM but I think we shouldn't do it unless we have a use case.
I have specialized Some numbers: Using
That's not too bad for an un-optimized implementation of the FD methods, I think. |
Oh, wow! |
|
If you promise to finish the rest as well, I'd suggest to consider merging this with the current feature set and doing restructuring/cleanup/renaming in a separate PR |
|
This PR should be in a good shape and more or less ready to merge after another round of reviews. Other features and extensions can be added in new PRs later. |
|
Tests only fail because of the following warnings on Windows: How shall we handle that?
|
|
I'll merge this since the failure on windows-mpi is spurious and not related to this PR. |
This PR basically collects two parts.
Part 1: Surface integrals/terms for additional dispatch (basically finished)
Summary: I introduced a surface integral/term type which is passed to every surface computation, mimicking the approach we have for our volume integrals/terms.
This could be considered to be breaking since it changes a documented constructor of our
DGtype. On the other hand, the documentation of this constructor was just wrong (since it contained the basis twice but not the mortar). Hence, this is more like a bug fix, I think ;).Currently, the
SurfaceintegralWeakFormis used everywhere; theSurfaceIntegralStrongFormwill be used for FD methods later. This will also be useful to get more general methods of @jlchan. When we merge this, we need to remember to notify others.Part 2: First step towards FD methods in Trixi (finished)
This is a rough draft adapting the 2D
TreeMeshDG solver. TODO:TreeMeshstuff (and other mesh types) accordinglyFixadd them to the list in Finite difference SBP methods #607TODO: FDorAdd referencesListed in Finite difference SBP methods #607Add docsListed in Finite difference SBP methods #607We need to consolidate terms/integrals (Naming convention: volume/surface *integral* → *terms* #616) before merging thisListed in Finite difference SBP methods #607Move code to appropriate places (cf. Consolidate mesh infrastructure after adding structured and unstructured ones #542). We need more levels in our solver hierarchy.Listed in Finite difference SBP methods #607DG: Current sub-methods are FDSBP and DGSEMdg_common.jland adgsem_common.jl.tree_2d_dgsem,tree_2d_fdsbpetc.