!> @file eno3.f90 !> @brief Classical ENO3 (Harten et al. 1987) reconstruction scheme. !! !! Selects the smoothest 3-cell substencil per equation component using !! Newton divided differences. !! !! References: !! [1] 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. module eno3 use precision, only: wp use solver_state, only: neq implicit none private public :: eno3_reconstruct contains ! --------------------------------------------------------------------------- !> Perform classical ENO3 reconstruction on the supplied stencil. !! !! Selects the smoothest 3-cell substencil per equation component using !! Newton divided differences (Harten et al., 1987). !! !! Stencil column mapping (5-point, half_stencil=3): !! col 1 = v_{i-2}, col 2 = v_{i-1}, col 3 = v_i (centre), !! col 4 = v_{i+1}, col 5 = v_{i+2}. !! !! Two-step stencil selection per component: !! Step 1 — compare |d^1[i-1]| vs |d^1[i]| (1st-order divided differences): !! d^1[j] = v_{j+1} - v_j !! Step 2 — compare |d^2[s-1]| vs |d^2[s]| (2nd-order divided differences): !! d^2[j] = (v_{j+2} - 2*v_{j+1} + v_j) / 2 !! !! Selected substencil starts at column s ∈ {1, 2, 3} and uses the same !! reconstruction formulas as the three WENO5 substencils: !! S0 (s=1): f̂ = (1/3)*v1 - (7/6)*v2 + (11/6)*v3 !! S1 (s=2): f̂ = -(1/6)*v1 + (5/6)*v2 + (1/3)*v3 !! S2 (s=3): f̂ = (1/3)*v1 + (5/6)*v2 - (1/6)*v3 !! !! @param[in] f_stencil Stencil at 5 points (neq x 5) !! @param[out] f_hat Reconstructed value at the interface ! --------------------------------------------------------------------------- pure subroutine eno3_reconstruct(f_stencil, f_hat) real(wp), intent(in) :: f_stencil(:, :) real(wp), intent(out) :: f_hat(:) real(wp) :: d1_l, d1_r, d2_l, d2_r real(wp) :: v1, v2, v3 integer :: ieq, s do ieq = 1, neq ! --- Step 1: 1st-order divided differences at i-1 and i --- ! d^1[i-1] = f_i - f_{i-1} = col3 - col2 ! d^1[i] = f_{i+1} - f_i = col4 - col3 d1_l = f_stencil(ieq, 3) - f_stencil(ieq, 2) d1_r = f_stencil(ieq, 4) - f_stencil(ieq, 3) if (abs(d1_l) <= abs(d1_r)) then s = 2 ! 2-cell {i-1, i} = {col2, col3} else s = 3 ! 2-cell {i, i+1} = {col3, col4} end if ! --- Step 2: 2nd-order divided differences for left/right extension --- ! d^2[s-1] = (col(s+1) - 2*col(s) + col(s-1)) / 2 [left extension] ! d^2[s] = (col(s+2) - 2*col(s+1) + col(s)) / 2 [right extension] d2_l = (f_stencil(ieq, s + 1) - 2.0_wp * f_stencil(ieq, s) & & + f_stencil(ieq, s - 1)) * 0.5_wp d2_r = (f_stencil(ieq, s + 2) - 2.0_wp * f_stencil(ieq, s + 1) & & + f_stencil(ieq, s)) * 0.5_wp if (abs(d2_l) <= abs(d2_r)) s = s - 1 ! extend left ! --- Apply selected substencil reconstruction formula --- v1 = f_stencil(ieq, s) v2 = f_stencil(ieq, s + 1) v3 = f_stencil(ieq, s + 2) select case (s) case (1) ! S0 = {i-2, i-1, i} f_hat(ieq) = 1.0_wp / 3.0_wp * v1 - 7.0_wp / 6.0_wp * v2 + 11.0_wp / 6.0_wp * v3 case (2) ! S1 = {i-1, i, i+1} f_hat(ieq) = -1.0_wp / 6.0_wp * v1 + 5.0_wp / 6.0_wp * v2 + 1.0_wp / 3.0_wp * v3 case (3) ! S2 = {i, i+1, i+2} f_hat(ieq) = 1.0_wp / 3.0_wp * v1 + 5.0_wp / 6.0_wp * v2 - 1.0_wp / 6.0_wp * v3 end select end do end subroutine eno3_reconstruct end module eno3