Use ADER CFL Correction Factors When Determining FV Patch Resolution
The ADER-DG time step size does not only depend on the DG CFL number that
scales with 1/(2*N+1)
but by an additional order-dependent correction factor C_ADER(N)
.
For order 5, e.g., this factor is roughly 0.5.
We currently use 2*N+1 FV cells per FV patch in our implementation of the a posteriori limiting ADER-DG method. The motivation is to match the time step size of the coupled ADER-DG and the FV scheme. This does not take the correction factor into account, i.e. we add unneccessary numerical dissipation as we use unnecessary small time step sizes in the FV cells.
Using the same ADER-DG time step size, we can increase the FV cell sizes by
a factor 1/C_ADER(N)
. For order 5, e.g., we can double the size of the FV patches.
The following table corrects the patch sizes. I round down to not go below the ADER-DG time step size with the FV time step sizes.
Order N | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 |
---|---|---|---|---|---|---|---|---|---|---|
C_ADER | 100.0 % | 99.0 % | 85.0 % | 70.0 % | 62.1 % | 49.5 % | 49.4 % | 45.0 % | 34.0 % | 28.5 % |
Patch Scaling | 1.00 | 1.01 | 1.18 | 1.43 | 1.61 | 2.02 | 2.02 | 2.22 | 2.94 | 3.51 |
Corr. Patch Size (Orig. Patch Size) | 1 | 3 | 5 | 10 (7) | 14 (9) | 22 (11) | 26 (13) | 33 (15) | 50 (17) | 66 (19) |
Tasks
Proposed switch "Limiting-ADER-DG" solver specification:
"patch_size" : ("original" | "max" | number)
-
max
would look up the corrected patch size in a lookup table. -
original
would use the original 2*order+1. -
number
would use a user prescribed number. Control code must check if the FV time step size is smaller or equal to the ADER-DG time step size. Otherwise, an error must be printed out OR the ADER-DG method's CFL condition must be adjusted.If a number is specified, then the ADER-DG to FV projectors have to be constructed manually. For "max" and "original", they can be precomputed once.
-
Decide if a generic implementation of the patch size feasible or if the corrected patch sizes should be hardcoded?. A: Eventually both but first the hardcoded versions.
-
Modify toolkit's solver controllers and templates: -
exahype-specfile.schema.json
: Extendpatch_size
option in spec file schema -
exahype/toolkit/solverController.py
: Takepatch_size
parameter into account. Found issue: Defaults are not set. -
Generic kernels: -
Check that the CFL factor is picked up properly in the ADER-DG stableTimeStepSize kernel -
GenericLimiterSolver(Header|Implementation).template
:-
Update arguments passed to limiter kernels -
Allocate Limiter projection matrices -
Dimensions:
uh2lim[n] = new double[basisSize*basisSizeLim]();
uh2lob[n] = new double[basisSize*basisSize]();
lim2uh[n] = new double[basisSizeLim*basisSize]();
-
Look up the FV to DG least-squares projection -
Main reference: M. Dumbser et al., A posteriori subcell limiting of the discontinuous Galerkin finite element method for hyperbolic conservation laws, Journal of Computational Physics, vol. 278, pp. 47–75, Dec. 2014.
-
Least-squares reference: M. Dumbser and M. Käser, Arbitrary high order non-oscillatory finite volume schemes on unstructured meshes for linear hyperbolic systems, Journal of Computational Physics, vol. 221, no. 2, pp. 693–723, Feb. 2007.
-
Objective function is
\min_{u_h} \| u_h(x) - v_h(x) \|_{L^2(K)} = ([A] \{u\} - \{v\}) \cdot ([A] \{u\} - \{v\})
s.t. mass conservation constraint.
[A]
is the DG to FV projection matrix. The objective function evaluates the DG solutionu_h
at the FV cell centersv_h
. -
Mass conservation constraint:
\int_K u_h(x) dx = \int_K v_h(x)\,\text{d}x
\sum_i u_i \int_{\bar K} \bar\phi_i (\bar x) \text{d}\bar x = \sum_j v_j \cdot \frac{\bar V}{\bar K}
\sum_i u_i w_i = \sum_j v_j \cdot \frac{1}{\text{volumes per patch}}
[C] \{u\} = [R] \{v\}
Above steps exploited properties of DG basis functions (Gauss-Legendre support points).
-
Constraint coupled via Lagrange multiplier. Updated objective function:
f = ([A] \{u\} - \{v\}) \cdot ([A] \{u\} - \{v\}) + \lambda \cdot ([C] \{u\} - [R] \{v\})
-
Optimality condition (degrees of freedom/parameters:
\{u\}
and\lambda
):\begin{pmatrix} 2 [A] [A] & -[C]\\ [C] & 0 \end{pmatrix} \cdot \begin{pmatrix} \{u\}\\ \lambda \end{pmatrix} = \begin{pmatrix} 2\,[A]\,\{v\}\\ [R]\,\{v\} \end{pmatrix}
To obtain solution
\{u\}
and Lagrange multipliers\lambda
, the equations must be inverted.C
is sparse due to the tensor product structure.[A]
is dense. As each nodal DG DOF is associated with a basis functions, the above procedure yields a FV to DG projector. -
I can reuse JM's existing code. Essential is that I replace
[A]
, which I do anyway. Convenient.
-
-
-
-
kernels/limiter/generic/Limiter.(h|cpp)
: Update signatures to take limiter patch size into account -
LimiterProjectionMatrices.(h|cpp)
: Do not use original functions anymore
-
-
Optimised kernels: -
Check if changes are necessary in OptimisedLimiterSolver(Header|Implementation).template
-
Make sure dg2fv
,fv2dg
andnDofLim
,nDofLimPad
,getBasisSizeLimiter()`, and
getBasisSizeLimiterPadded()`` are defined and used correctly:-
Usage: -
CodeGenerator/codegenerator/templates/configurationParameters_cpph.template (20 matches) -
CodeGenerator/codegenerator/templates/limiter_cpp.template (51 matches) -
CodeGenerator/codegenerator/templates/Quadrature_cpp.template (2 matches)
-
-
Definition: -
CodeGenerator/codegenerator/templates/controller.py (5 matches) -
CodeGenerator/codegenerator/models/limiterModel.py (29 matches) -
CodeGenerator/codegenerator/models/quadratureModel.py (4 matches)
-
-
-
-