Provides five MUSCL reconstructors, one per TVD limiter. The dispatcher in the reconstruction module selects the appropriate procedure pointer at initialisation time based on the limiter name supplied by the user.
Stencil layout (half_stencil=2, stencil_width=3): col 1 = f_{i-1}, col 2 = f_i (centre), col 3 = f_{i+1}
Formula (component-wise, left-biased stencil): ΔL = s(k,2) - s(k,1) ΔR = s(k,3) - s(k,2) r = ΔR / ΔL (when |ΔL| > tiny, else 0 → phi = 0) f_hat(k) = s(k,2) + 0.5 · phi(r) · ΔL
The right-biased (reversed) stencil for F^- is assembled by spatial_discretization in mirror order; the same formula gives the correct downwind reconstruction automatically.
References: [1] B. van Leer, "Towards the Ultimate Conservative Difference Scheme. V. A Second-Order Sequel to Godunov's Method," J. Comput. Phys., 32:101-136, 1979.
Abstract interface for a pure scalar TVD limiter function phi(r).
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | r |
Minmod TVD limiter. Most dissipative of the standard limiters; guaranteed TVD for any CFL within the stability limit.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | r |
Superbee TVD limiter. Least dissipative TVD limiter; may over-sharpen contact discontinuities on smooth problems.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | r |
Monotonized Central (MC) TVD limiter. Good balance between dissipation and resolution; recommended for most shock-capturing applications.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | r |
van Leer smooth TVD limiter.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | r |
Koren TVD limiter. Third-order accurate for smooth data; preserves monotonicity near discontinuities.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | r |
Core MUSCL reconstruction kernel (component-wise, left-biased stencil).
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | f_stencil(:,:) | |||
| real(kind=wp), | intent(out) | :: | f_hat(:) | |||
| procedure(limiter_iface) | :: | lim_fn |
TVD limiter function phi(r) |
MUSCL reconstruction with the Minmod limiter.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | f_stencil(:,:) | |||
| real(kind=wp), | intent(out) | :: | f_hat(:) |
MUSCL reconstruction with the Superbee limiter.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | f_stencil(:,:) | |||
| real(kind=wp), | intent(out) | :: | f_hat(:) |
MUSCL reconstruction with the Monotonized Central (MC) limiter.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | f_stencil(:,:) | |||
| real(kind=wp), | intent(out) | :: | f_hat(:) |
MUSCL reconstruction with the van Leer smooth limiter.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | f_stencil(:,:) | |||
| real(kind=wp), | intent(out) | :: | f_hat(:) |
MUSCL reconstruction with the Koren limiter.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | f_stencil(:,:) | |||
| real(kind=wp), | intent(out) | :: | f_hat(:) |