Provides init_recon_scheme(), which binds the reconstruction procedure pointers and stencil metadata inside a @p solver_state_t instance. All scheme implementations live in dedicated modules (one per scheme): weno5_js, weno5_z, weno_cu6, weno7_js, weno9_js, weno11_js — WENO family (shared helpers in weno_family) eno3 — Classical ENO3 muscl — MUSCL + TVD limiters mp5 — Monotonicity-Preserving 5th order teno5 — Targeted ENO 5th order upwind1, upwind2, central2 — Low-order baseline schemes
The constants @p max_stencil_width and @p max_coupling_radius are re-exported here for convenience so callers do not need to depend directly on solver_interfaces.
Supported scheme names: 'weno5' — WENO5-JS (Jiang & Shu, 1996) 'weno5z' — WENO5-Z (Borges et al., 2008) 'weno_cu6'— WENO-CU6 (Hu, Wang & Adams, 2010) 'weno7' — WENO7-JS (Balsara & Shu, 2000) 'weno9' — WENO9-JS (Balsara & Shu, 2000) 'weno11' — WENO11-JS (Balsara & Shu, 2000) 'eno3' — Classical ENO3 (Harten et al., 1987) 'muscl' — MUSCL + TVD limiter (van Leer, 1979); limiter selectable via init_recon_scheme 'mp5' — Monotonicity-Preserving 5th order (Suresh & Huynh, 1997) 'teno5' — Targeted ENO 5th order (Fu, Hu & Adams, 2016) 'upwind1' — First-order upwind 'upwind2' — Second-order upwind 'central2'— Second-order central (no dissipation; smooth problems only)
Auto defaults for char_proj: 'weno5', 'weno5z', 'weno_cu6', 'weno7', 'weno9', 'weno11', 'eno3', 'mp5', 'teno5' → .true. 'muscl', 'upwind1', 'upwind2', 'central2' → .false.
References: [1] G.-S. Jiang and C.-W. Shu, "Efficient Implementation of Weighted ENO Schemes," J. Comput. Phys., 126:202-228, 1996. [2] D. S. Balsara and C.-W. Shu, "Monotonicity Preserving Weighted Essentially Non-oscillatory Schemes with Increasingly High Order of Accuracy," J. Comput. Phys., 160:405-452, 2000. [3] R. Borges, M. Carmona, B. Costa, and W. S. Don, "An improved weighted essentially non-oscillatory scheme for hyperbolic conservation laws," J. Comput. Phys., 227:3191-3211, 2008. [4] X. Y. Hu, Q. Wang, and N. A. Adams, "An adaptive central-upwind weighted essentially non-oscillatory scheme," J. Comput. Phys., 229:8952-8965, 2010. [5] A. Harten, B. Engquist, S. Osher, and S. R. Chakravarthy, "Uniformly high-order accurate essentially non-oscillatory schemes, III," J. Comput. Phys., 71:231-303, 1987. [6] 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. [7] A. Suresh and H. T. Huynh, "Accurate Monotonicity-Preserving Schemes with Runge-Kutta Time Stepping," J. Comput. Phys., 136:83-99, 1997. [8] L. Fu, X. Y. Hu, and N. A. Adams, "A family of high-order targeted ENO schemes for compressible-fluid simulations," J. Comput. Phys., 305:333-359, 2016.
Bind the reconstruction procedure pointers and stencil metadata in @p state to the requested scheme, then resolve state%use_char_proj.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(solver_state_t), | intent(inout) | :: | state |
Solver instance; reconstruction pointers, stencil |
||
| character(len=*), | intent(in) | :: | scheme |
Name of the reconstruction scheme. |
||
| character(len=*), | intent(in), | optional | :: | char_proj |
(optional) Characteristic-projection override. |
|
| character(len=*), | intent(in), | optional | :: | limiter |
(optional) TVD limiter for 'muscl'. Valid values: |