20679838146

Committed 03 Jan 2026 04:24PM UTC coverage: 96.993%. Remained the same

Build # 20679838146

Build Type

push

github

Committed by

ranocha

Commit Message

set development version to v0.13.22-DEV

Coverage Stats

42292 of 43603 relevant lines covered (96.99%)

104628030.91 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

100.0

/src/solvers/dgsem_tree/subcell_finite_volume_O2.jl

"""
    reconstruction_constant(u_ll, u_lr, u_rl, u_rr,
                            sc_interface_coords,
                            node_index, limiter, dg)

Returns the constant "reconstructed" values `u_lr, u_rl` at the interface `sc_interface_coords[node_index - 1]`.
Supposed to be used in conjunction with [`VolumeIntegralPureLGLFiniteVolumeO2`](@ref).
Formally first order accurate.
If a first-order finite volume scheme is desired, [`VolumeIntegralPureLGLFiniteVolume`](@ref) is an
equivalent, but more efficient choice.
"""
@inline function reconstruction_constant(u_ll, u_lr, u_rl, u_rr,
                                         sc_interface_coords, node_index,
                                         limiter, dg)
    return u_lr, u_rl
end

# Helper functions for reconstructions below
@inline function reconstruction_linear(u_lr, u_rl, s_l, s_r,
                                       x_lr, x_rl, sc_interface_coords, node_index)
    # Linear reconstruction at the interface
    u_lr = u_lr + s_l * (sc_interface_coords[node_index - 1] - x_lr)
    u_rl = u_rl + s_r * (sc_interface_coords[node_index - 1] - x_rl)

    return u_lr, u_rl
end

#             Reference element:             
#  -1 ------------------0------------------ 1 -> x
# Gauss-Lobatto-Legendre nodes (schematic for k = 3):
#   .          .                  .         .
#   ^          ^                  ^         ^
# Node indices:
#   1          2                  3         4
# The inner subcell boundaries are governed by the
# cumulative sum of the quadrature weights - 1 .
#  -1 ------------------0------------------ 1 -> x
#        w1-1      (w1+w2)-1   (w1+w2+w3)-1
#   |     |             |             |     |
# Note that only the inner boundaries are stored.
# Subcell interface indices, loop only over 2 -> nnodes(dg) = 4
#   1     2             3             4     5
#
# In general a four-point stencil is required, since we reconstruct the
# piecewise linear solution in both subcells next to the subcell interface.
# Since these subcell boundaries are not aligned with the DG nodes,
# on each neighboring subcell two linear solutions are reconstructed => 4 point stencil.
# For the outer interfaces the stencil shrinks since we do not consider values 
# outside the element (volume integral).
# 
# The left subcell node values are labelled `_ll` (left-left) and `_lr` (left-right), while
# the right subcell node values are labelled `_rl` (right-left) and `_rr` (right-right).

"""
    reconstruction_O2_full(u_ll, u_lr, u_rl, u_rr,
                           sc_interface_coords, node_index,
                           limiter, dg::DGSEM)

Returns the reconstructed values `u_lr, u_rl` at the interface `sc_interface_coords[node_index - 1]`.
Computes limited (linear) slopes on the subcells for a DGSEM element.
Supposed to be used in conjunction with [`VolumeIntegralPureLGLFiniteVolumeO2`](@ref).

The supplied `limiter` governs the choice of slopes given the nodal values
`u_ll`, `u_lr`, `u_rl`, and `u_rr` at the (Gauss-Lobatto Legendre) nodes.
Total-Variation-Diminishing (TVD) choices for the limiter are
    1) [`minmod`](@ref)
    2) [`monotonized_central`](@ref)
    3) [`superbee`](@ref)
    4) [`vanLeer`](@ref)

The reconstructed slopes are for `reconstruction_O2_full` not limited at the cell boundaries.
Formally second order accurate when used without a limiter, i.e., `limiter = `[`central_slope`](@ref).
This approach corresponds to equation (79) described in
- Rueda-Ramírez, Hennemann, Hindenlang, Winters, & Gassner (2021).
  "An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations.
   Part II: Subcell finite volume shock capturing"
  [JCP: 2021.110580](https://doi.org/10.1016/j.jcp.2021.110580)
"""
@inline function reconstruction_O2_full(u_ll, u_lr, u_rl, u_rr,
                                        sc_interface_coords, node_index,
                                        limiter, dg::DGSEM)
    @unpack nodes = dg.basis
    x_lr = nodes[node_index - 1]
    x_rl = nodes[node_index]

    # Slope between "middle" nodes
    s_m = (u_rl - u_lr) / (x_rl - x_lr)

    if node_index == 2 # Catch case ll == lr
        s_l = s_m # Use unlimited "central" slope
    else
        x_ll = nodes[node_index - 2]
        # Slope between "left" nodes
        s_lr = (u_lr - u_ll) / (x_lr - x_ll)
        # Select slope between extrapolated (left) and crossing (middle) slope
        s_l = limiter.(s_lr, s_m)
    end

    if node_index == nnodes(dg) # Catch case rl == rr
        s_r = s_m # Use unlimited "central" slope
    else
        x_rr = nodes[node_index + 1]
        # Slope between "right" nodes
        s_rl = (u_rr - u_rl) / (x_rr - x_rl)
        # Select slope between crossing (middle) and extrapolated (right) slope
        s_r = limiter.(s_m, s_rl)
    end

    return reconstruction_linear(u_lr, u_rl, s_l, s_r,
                                 x_lr, x_rl, sc_interface_coords, node_index)
end

"""
    reconstruction_O2_inner(u_ll, u_lr, u_rl, u_rr,
                            sc_interface_coords, node_index,
                            limiter, dg::DGSEM)

Returns the reconstructed values `u_lr, u_rl` at the interface `sc_interface_coords[node_index - 1]`.
Computes limited (linear) slopes on the *inner* subcells for a DGSEM element.
Supposed to be used in conjunction with [`VolumeIntegralPureLGLFiniteVolumeO2`](@ref).

The supplied `limiter` governs the choice of slopes given the nodal values
`u_ll`, `u_lr`, `u_rl`, and `u_rr` at the (Gauss-Lobatto Legendre) nodes.
Total-Variation-Diminishing (TVD) choices for the limiter are
    1) [`minmod`](@ref)
    2) [`monotonized_central`](@ref)
    3) [`superbee`](@ref)
    4) [`vanLeer`](@ref)

For the outer, i.e., boundary subcells, constant values are used, i.e, no reconstruction.
This reduces the order of the scheme below 2.
This approach corresponds to equation (78) described in
- Rueda-Ramírez, Hennemann, Hindenlang, Winters, & Gassner (2021).
  "An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations. 
   Part II: Subcell finite volume shock capturing"
  [JCP: 2021.110580](https://doi.org/10.1016/j.jcp.2021.110580)
"""
@inline function reconstruction_O2_inner(u_ll, u_lr, u_rl, u_rr,
                                         sc_interface_coords, node_index,
                                         limiter, dg::DGSEM)
    @unpack nodes = dg.basis
    x_lr = nodes[node_index - 1]
    x_rl = nodes[node_index]

    # Slope between "middle" nodes
    s_m = (u_rl - u_lr) / (x_rl - x_lr)

    if node_index == 2 # Catch case ll == lr
        # Do not reconstruct at the boundary
        s_l = zero(s_m)
    else
        x_ll = nodes[node_index - 2]
        # Slope between "left" nodes
        s_lr = (u_lr - u_ll) / (x_lr - x_ll)
        # Select slope between extrapolated (left) and crossing (middle) slope
        s_l = limiter.(s_lr, s_m)
    end

    if node_index == nnodes(dg) # Catch case rl == rr
        # Do not reconstruct at the boundary
        s_r = zero(s_m)
    else
        x_rr = nodes[node_index + 1]
        # Slope between "right" nodes
        s_rl = (u_rr - u_rl) / (x_rr - x_rl)
        # Select slope between crossing (middle) and extrapolated (right) slope
        s_r = limiter.(s_m, s_rl)
    end

    return reconstruction_linear(u_lr, u_rl, s_l, s_r,
                                 x_lr, x_rl, sc_interface_coords, node_index)
end

"""
    central_slope(sl, sr)

Central, non-TVD reconstruction given left and right slopes `sl` and `sr`.
Gives formally full order of accuracy at the expense of sacrificed nonlinear stability.
Similar in spirit to [`flux_central`](@ref).
"""
@inline function central_slope(sl, sr)
    return 0.5f0 * (sl + sr)
end

"""
    minmod(sl, sr)

Classic minmod limiter function for a TVD reconstruction given left and right slopes `sl` and `sr`.
There are many different ways how the minmod limiter can be implemented.
For reference, see for instance Eq. (6.27) in

- Randall J. LeVeque (2002)
  Finite Volume Methods for Hyperbolic Problems
  [DOI: 10.1017/CBO9780511791253](https://doi.org/10.1017/CBO9780511791253)
"""
@inline function minmod(sl, sr)
    return 0.5f0 * (sign(sl) + sign(sr)) * min(abs(sl), abs(sr))
end

"""
    monotonized_central(sl, sr)

Monotonized central limiter function for a TVD reconstruction given left and right slopes `sl` and `sr`.
There are many different ways how the monotonized central limiter can be implemented.
For reference, see for instance Eq. (6.29) in

- Randall J. LeVeque (2002)
  Finite Volume Methods for Hyperbolic Problems
  [DOI: 10.1017/CBO9780511791253](https://doi.org/10.1017/CBO9780511791253)
"""
@inline function monotonized_central(sl, sr)
    # Use recursive property of minmod function
    return minmod(0.5f0 * (sl + sr), minmod(2 * sl, 2 * sr))
end

"""
    superbee(sl, sr)

Superbee limiter function for a TVD reconstruction given left and right slopes `sl` and `sr`.
There are many different ways how the superbee limiter can be implemented.
For reference, see for instance Eq. (6.28) in

- Randall J. LeVeque (2002)
  Finite Volume Methods for Hyperbolic Problems
  [DOI: 10.1017/CBO9780511791253](https://doi.org/10.1017/CBO9780511791253)
"""
@inline function superbee(sl, sr)
    return maxmod(minmod(sl, 2 * sr), minmod(2 * sl, sr))
end

"""
    vanLeer(sl, sr)

Symmetric limiter by van Leer.
See for reference page 70 in 

- Siddhartha Mishra, Ulrik Skre Fjordholm and Rémi Abgrall
  Numerical methods for conservation laws and related equations.
  [Link](https://metaphor.ethz.ch/x/2019/hs/401-4671-00L/literature/mishra_hyperbolic_pdes.pdf)
"""
@inline function vanLeer(sl, sr)
    if abs(sl) + abs(sr) > zero(sl)
        return (abs(sr) * sl + abs(sl) * sr) / (abs(sl) + abs(sr))
    else
        return zero(sl)
    end
end

1	"""
2	reconstruction_constant(u_ll, u_lr, u_rl, u_rr,
3	sc_interface_coords,
4	node_index, limiter, dg)
5
6	Returns the constant "reconstructed" values `u_lr, u_rl` at the interface `sc_interface_coords[node_index - 1]`.
7	Supposed to be used in conjunction with [`VolumeIntegralPureLGLFiniteVolumeO2`](@ref).
8	Formally first order accurate.
9	If a first-order finite volume scheme is desired, [`VolumeIntegralPureLGLFiniteVolume`](@ref) is an
10	equivalent, but more efficient choice.
11	"""
12	@inline function reconstruction_constant(u_ll, u_lr, u_rl, u_rr,	42,469,777✔
13	sc_interface_coords, node_index,
14	limiter, dg)
15	return u_lr, u_rl	42,469,777✔
16	end
17
18	# Helper functions for reconstructions below
19	@inline function reconstruction_linear(u_lr, u_rl, s_l, s_r,	241,955,856✔
20	x_lr, x_rl, sc_interface_coords, node_index)
21	# Linear reconstruction at the interface
22	u_lr = u_lr + s_l * (sc_interface_coords[node_index - 1] - x_lr)	241,955,856✔
23	u_rl = u_rl + s_r * (sc_interface_coords[node_index - 1] - x_rl)	241,955,856✔
24
25	return u_lr, u_rl	241,955,856✔
26	end
27
28	# Reference element:
29	# -1 ------------------0------------------ 1 -> x
30	# Gauss-Lobatto-Legendre nodes (schematic for k = 3):
31	# . . . .
32	# ^ ^ ^ ^
33	# Node indices:
34	# 1 2 3 4
35	# The inner subcell boundaries are governed by the
36	# cumulative sum of the quadrature weights - 1 .
37	# -1 ------------------0------------------ 1 -> x
38	# w1-1 (w1+w2)-1 (w1+w2+w3)-1
39	# \| \| \| \| \|
40	# Note that only the inner boundaries are stored.
41	# Subcell interface indices, loop only over 2 -> nnodes(dg) = 4
42	# 1 2 3 4 5
43	#
44	# In general a four-point stencil is required, since we reconstruct the
45	# piecewise linear solution in both subcells next to the subcell interface.
46	# Since these subcell boundaries are not aligned with the DG nodes,
47	# on each neighboring subcell two linear solutions are reconstructed => 4 point stencil.
48	# For the outer interfaces the stencil shrinks since we do not consider values
49	# outside the element (volume integral).
50	#
51	# The left subcell node values are labelled `_ll` (left-left) and `_lr` (left-right), while
52	# the right subcell node values are labelled `_rl` (right-left) and `_rr` (right-right).
53
54	"""
55	reconstruction_O2_full(u_ll, u_lr, u_rl, u_rr,
56	sc_interface_coords, node_index,
57	limiter, dg::DGSEM)
58
59	Returns the reconstructed values `u_lr, u_rl` at the interface `sc_interface_coords[node_index - 1]`.
60	Computes limited (linear) slopes on the subcells for a DGSEM element.
61	Supposed to be used in conjunction with [`VolumeIntegralPureLGLFiniteVolumeO2`](@ref).
62
63	The supplied `limiter` governs the choice of slopes given the nodal values
64	`u_ll`, `u_lr`, `u_rl`, and `u_rr` at the (Gauss-Lobatto Legendre) nodes.
65	Total-Variation-Diminishing (TVD) choices for the limiter are
66	1) [`minmod`](@ref)
67	2) [`monotonized_central`](@ref)
68	3) [`superbee`](@ref)
69	4) [`vanLeer`](@ref)
70
71	The reconstructed slopes are for `reconstruction_O2_full` not limited at the cell boundaries.
72	Formally second order accurate when used without a limiter, i.e., `limiter = `[`central_slope`](@ref).
73	This approach corresponds to equation (79) described in
74	- Rueda-Ramírez, Hennemann, Hindenlang, Winters, & Gassner (2021).
75	"An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations.
76	Part II: Subcell finite volume shock capturing"
77	[JCP: 2021.110580](https://doi.org/10.1016/j.jcp.2021.110580)
78	"""
79	@inline function reconstruction_O2_full(u_ll, u_lr, u_rl, u_rr,	196,385,616✔
80	sc_interface_coords, node_index,
81	limiter, dg::DGSEM)
82	@unpack nodes = dg.basis	196,385,616✔
83	x_lr = nodes[node_index - 1]	196,385,616✔
84	x_rl = nodes[node_index]	196,385,616✔
85
86	# Slope between "middle" nodes
87	s_m = (u_rl - u_lr) / (x_rl - x_lr)	196,385,616✔
88
89	if node_index == 2 # Catch case ll == lr	196,385,616✔
90	s_l = s_m # Use unlimited "central" slope	65,461,872✔
91	else
92	x_ll = nodes[node_index - 2]	130,923,744✔
93	# Slope between "left" nodes
94	s_lr = (u_lr - u_ll) / (x_lr - x_ll)	130,923,744✔
95	# Select slope between extrapolated (left) and crossing (middle) slope
96	s_l = limiter.(s_lr, s_m)	130,923,744✔
97	end
98
99	if node_index == nnodes(dg) # Catch case rl == rr	196,385,616✔
100	s_r = s_m # Use unlimited "central" slope	65,461,872✔
101	else
102	x_rr = nodes[node_index + 1]	130,923,744✔
103	# Slope between "right" nodes
104	s_rl = (u_rr - u_rl) / (x_rr - x_rl)	130,923,744✔
105	# Select slope between crossing (middle) and extrapolated (right) slope
106	s_r = limiter.(s_m, s_rl)	130,923,744✔
107	end
108
109	return reconstruction_linear(u_lr, u_rl, s_l, s_r,	196,385,616✔
110	x_lr, x_rl, sc_interface_coords, node_index)
111	end
112
113	"""
114	reconstruction_O2_inner(u_ll, u_lr, u_rl, u_rr,
115	sc_interface_coords, node_index,
116	limiter, dg::DGSEM)
117
118	Returns the reconstructed values `u_lr, u_rl` at the interface `sc_interface_coords[node_index - 1]`.
119	Computes limited (linear) slopes on the inner subcells for a DGSEM element.
120	Supposed to be used in conjunction with [`VolumeIntegralPureLGLFiniteVolumeO2`](@ref).
121
122	The supplied `limiter` governs the choice of slopes given the nodal values
123	`u_ll`, `u_lr`, `u_rl`, and `u_rr` at the (Gauss-Lobatto Legendre) nodes.
124	Total-Variation-Diminishing (TVD) choices for the limiter are
125	1) [`minmod`](@ref)
126	2) [`monotonized_central`](@ref)
127	3) [`superbee`](@ref)
128	4) [`vanLeer`](@ref)
129
130	For the outer, i.e., boundary subcells, constant values are used, i.e, no reconstruction.
131	This reduces the order of the scheme below 2.
132	This approach corresponds to equation (78) described in
133	- Rueda-Ramírez, Hennemann, Hindenlang, Winters, & Gassner (2021).
134	"An entropy stable nodal discontinuous Galerkin method for the resistive MHD equations.
135	Part II: Subcell finite volume shock capturing"
136	[JCP: 2021.110580](https://doi.org/10.1016/j.jcp.2021.110580)
137	"""
138	@inline function reconstruction_O2_inner(u_ll, u_lr, u_rl, u_rr,	45,570,240✔
139	sc_interface_coords, node_index,
140	limiter, dg::DGSEM)
141	@unpack nodes = dg.basis	45,570,240✔
142	x_lr = nodes[node_index - 1]	45,570,240✔
143	x_rl = nodes[node_index]	45,570,240✔
144
145	# Slope between "middle" nodes
146	s_m = (u_rl - u_lr) / (x_rl - x_lr)	45,570,240✔
147
148	if node_index == 2 # Catch case ll == lr	45,570,240✔
149	# Do not reconstruct at the boundary
150	s_l = zero(s_m)	14,849,880✔
151	else
152	x_ll = nodes[node_index - 2]	30,720,360✔
153	# Slope between "left" nodes
154	s_lr = (u_lr - u_ll) / (x_lr - x_ll)	30,720,360✔
155	# Select slope between extrapolated (left) and crossing (middle) slope
156	s_l = limiter.(s_lr, s_m)	30,720,360✔
157	end
158
159	if node_index == nnodes(dg) # Catch case rl == rr	45,570,240✔
160	# Do not reconstruct at the boundary
161	s_r = zero(s_m)	14,849,880✔
162	else
163	x_rr = nodes[node_index + 1]	30,720,360✔
164	# Slope between "right" nodes
165	s_rl = (u_rr - u_rl) / (x_rr - x_rl)	30,720,360✔
166	# Select slope between crossing (middle) and extrapolated (right) slope
167	s_r = limiter.(s_m, s_rl)	30,720,360✔
168	end
169
170	return reconstruction_linear(u_lr, u_rl, s_l, s_r,	45,570,240✔
171	x_lr, x_rl, sc_interface_coords, node_index)
172	end
173
174	"""
175	central_slope(sl, sr)
176
177	Central, non-TVD reconstruction given left and right slopes `sl` and `sr`.
178	Gives formally full order of accuracy at the expense of sacrificed nonlinear stability.
179	Similar in spirit to [`flux_central`](@ref).
180	"""
181	@inline function central_slope(sl, sr)	1✔
182	return 0.5f0 * (sl + sr)	1✔
183	end
184
185	"""
186	minmod(sl, sr)
187
188	Classic minmod limiter function for a TVD reconstruction given left and right slopes `sl` and `sr`.
189	There are many different ways how the minmod limiter can be implemented.
190	For reference, see for instance Eq. (6.27) in
191
192	- Randall J. LeVeque (2002)
193	Finite Volume Methods for Hyperbolic Problems
194	[DOI: 10.1017/CBO9780511791253](https://doi.org/10.1017/CBO9780511791253)
195	"""
196	@inline function minmod(sl, sr)	267,585,438✔
197	return 0.5f0 * (sign(sl) + sign(sr)) * min(abs(sl), abs(sr))	267,585,438✔
198	end
199
200	"""
201	monotonized_central(sl, sr)
202
203	Monotonized central limiter function for a TVD reconstruction given left and right slopes `sl` and `sr`.
204	There are many different ways how the monotonized central limiter can be implemented.
205	For reference, see for instance Eq. (6.29) in
206
207	- Randall J. LeVeque (2002)
208	Finite Volume Methods for Hyperbolic Problems
209	[DOI: 10.1017/CBO9780511791253](https://doi.org/10.1017/CBO9780511791253)
210	"""
211	@inline function monotonized_central(sl, sr)	133,792,710✔
212	# Use recursive property of minmod function
213	return minmod(0.5f0 * (sl + sr), minmod(2 * sl, 2 * sr))	133,792,710✔
214	end
215
216	"""
217	superbee(sl, sr)
218
219	Superbee limiter function for a TVD reconstruction given left and right slopes `sl` and `sr`.
220	There are many different ways how the superbee limiter can be implemented.
221	For reference, see for instance Eq. (6.28) in
222
223	- Randall J. LeVeque (2002)
224	Finite Volume Methods for Hyperbolic Problems
225	[DOI: 10.1017/CBO9780511791253](https://doi.org/10.1017/CBO9780511791253)
226	"""
227	@inline function superbee(sl, sr)	6✔
228	return maxmod(minmod(sl, 2 * sr), minmod(2 * sl, sr))	6✔
229	end
230
231	"""
232	vanLeer(sl, sr)
233
234	Symmetric limiter by van Leer.
235	See for reference page 70 in
236
237	- Siddhartha Mishra, Ulrik Skre Fjordholm and Rémi Abgrall
238	Numerical methods for conservation laws and related equations.
239	[Link](https://metaphor.ethz.ch/x/2019/hs/401-4671-00L/literature/mishra_hyperbolic_pdes.pdf)
240	"""
241	@inline function vanLeer(sl, sr)	1,118,717,623✔
242	if abs(sl) + abs(sr) > zero(sl)	1,118,717,623✔
243	return (abs(sr) * sl + abs(sl) * sr) / (abs(sl) + abs(sr))	1,110,416,946✔
244	else
245	return zero(sl)	8,300,677✔
246	end
247	end

trixi-framework / Trixi.jl / 20679838146

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous