!> @file weno11_js.f90 !> @brief WENO11-JS finite-difference reconstruction scheme. !! !! Classical Jiang-Shu WENO on an 11-point stencil using six 6-point !! substencils. As with WENO9, these are the finite-difference !! flux-reconstruction coefficients rather than plain pointwise interpolants. module weno11_js use precision, only: wp use weno_family, only: apply_weno_js_generic implicit none private real(wp), parameter :: candidate_coeffs(6, 6) = reshape([ & -1.0_wp / 6.0_wp, 31.0_wp / 30.0_wp, -163.0_wp / 60.0_wp, 79.0_wp & / 20.0_wp, & -71.0_wp / 20.0_wp, 49.0_wp / 20.0_wp, 1.0_wp / 30.0_wp, -13.0_wp & / 60.0_wp, & 37.0_wp / 60.0_wp, -21.0_wp / 20.0_wp, 29.0_wp / 20.0_wp, 1.0_wp & / 6.0_wp, & -1.0_wp / 60.0_wp, 7.0_wp / 60.0_wp, -23.0_wp / 60.0_wp, 19.0_wp & / 20.0_wp, & 11.0_wp / 30.0_wp, -1.0_wp / 30.0_wp, 1.0_wp / 60.0_wp, -2.0_wp & / 15.0_wp, & 37.0_wp / 60.0_wp, 37.0_wp / 60.0_wp, -2.0_wp / 15.0_wp, 1.0_wp & / 60.0_wp, & -1.0_wp / 30.0_wp, 11.0_wp / 30.0_wp, 19.0_wp / 20.0_wp, -23.0_wp & / 60.0_wp, & 7.0_wp / 60.0_wp, -1.0_wp / 60.0_wp, 1.0_wp / 6.0_wp, 29.0_wp / 20.0_wp, & -21.0_wp / 20.0_wp, 37.0_wp / 60.0_wp, -13.0_wp / 60.0_wp, 1.0_wp & / 30.0_wp & ], [6, 6]) real(wp), parameter :: linear_weights(6) = [ & 1.0_wp / 462.0_wp, 5.0_wp / 77.0_wp, 25.0_wp / 77.0_wp, 100.0_wp / 231.0_wp, & 25.0_wp / 154.0_wp, 1.0_wp / 77.0_wp] real(wp), parameter :: beta_mats(6, 6, 6) = reshape([ & 384187.0_wp / 40320.0_wp, -539591.0_wp / 5040.0_wp, 1840141.0_wp / & 7560.0_wp, -314063.0_wp / 1120.0_wp, & 661145.0_wp / 4032.0_wp, -235637.0_wp / 6048.0_wp, -539591.0_wp / 5040.0_wp, & 12160229.0_wp / 40320.0_wp, & -41615261.0_wp / 30240.0_wp, 2674951.0_wp / 1680.0_wp, -1048211.0_wp / & 1120.0_wp, 2706017.0_wp / 12096.0_wp, & 1840141.0_wp / 7560.0_wp, -41615261.0_wp / 30240.0_wp, 47689393.0_wp / & 30240.0_wp, -6937561.0_wp / 1890.0_wp, & 32862709.0_wp / 15120.0_wp, -15848531.0_wp / 30240.0_wp, -314063.0_wp / & 1120.0_wp, 2674951.0_wp / 1680.0_wp, & -6937561.0_wp / 1890.0_wp, 21703781.0_wp / 10080.0_wp, -25980937.0_wp / & 10080.0_wp, 4762921.0_wp / 7560.0_wp, & 661145.0_wp / 4032.0_wp, -1048211.0_wp / 1120.0_wp, 32862709.0_wp / & 15120.0_wp, -25980937.0_wp / 10080.0_wp, & 31617079.0_wp / 40320.0_wp, -2966279.0_wp / 7560.0_wp, -235637.0_wp / & 6048.0_wp, 2706017.0_wp / 12096.0_wp, & -15848531.0_wp / 30240.0_wp, 4762921.0_wp / 7560.0_wp, -2966279.0_wp / & 7560.0_wp, 6150211.0_wp / 120960.0_wp, & 90593.0_wp / 40320.0_wp, -188483.0_wp / 7560.0_wp, 139471.0_wp / 2520.0_wp, & -68601.0_wp / 1120.0_wp, & 2033509.0_wp / 60480.0_wp, -73379.0_wp / 10080.0_wp, -188483.0_wp / & 7560.0_wp, 8449957.0_wp / 120960.0_wp, & -9478331.0_wp / 30240.0_wp, 5300629.0_wp / 15120.0_wp, -5877617.0_wp / & 30240.0_wp, 2567287.0_wp / 60480.0_wp, & 139471.0_wp / 2520.0_wp, -9478331.0_wp / 30240.0_wp, 1197047.0_wp / & 3360.0_wp, -169859.0_wp / 210.0_wp, & 6881719.0_wp / 15120.0_wp, -1015303.0_wp / 10080.0_wp, -68601.0_wp / & 1120.0_wp, 5300629.0_wp / 15120.0_wp, & -169859.0_wp / 210.0_wp, 4721851.0_wp / 10080.0_wp, -16306061.0_wp / & 30240.0_wp, 61427.0_wp / 504.0_wp, & 2033509.0_wp / 60480.0_wp, -5877617.0_wp / 30240.0_wp, 6881719.0_wp / & 15120.0_wp, -16306061.0_wp / 30240.0_wp, & 19365967.0_wp / 120960.0_wp, -1139749.0_wp / 15120.0_wp, -73379.0_wp / & 10080.0_wp, 2567287.0_wp / 60480.0_wp, & -1015303.0_wp / 10080.0_wp, 61427.0_wp / 504.0_wp, -1139749.0_wp / & 15120.0_wp, 384187.0_wp / 40320.0_wp, & 139633.0_wp / 120960.0_wp, -178747.0_wp / 15120.0_wp, 178999.0_wp / & 7560.0_wp, -139633.0_wp / 6048.0_wp, & 662503.0_wp / 60480.0_wp, -12281.0_wp / 6048.0_wp, -178747.0_wp / & 15120.0_wp, 141661.0_wp / 4480.0_wp, & -1323367.0_wp / 10080.0_wp, 1991239.0_wp / 15120.0_wp, -643999.0_wp / & 10080.0_wp, 243127.0_wp / 20160.0_wp, & 178999.0_wp / 7560.0_wp, -1323367.0_wp / 10080.0_wp, 159219.0_wp / & 1120.0_wp, -559651.0_wp / 1890.0_wp, & 248681.0_wp / 1680.0_wp, -288521.0_wp / 10080.0_wp, -139633.0_wp / & 6048.0_wp, 1991239.0_wp / 15120.0_wp, & -559651.0_wp / 1890.0_wp, 4877743.0_wp / 30240.0_wp, -5106971.0_wp / & 30240.0_wp, 255397.0_wp / 7560.0_wp, & 662503.0_wp / 60480.0_wp, -643999.0_wp / 10080.0_wp, 248681.0_wp / & 1680.0_wp, -5106971.0_wp / 30240.0_wp, & 1884439.0_wp / 40320.0_wp, -1240.0_wp / 63.0_wp, -12281.0_wp / 6048.0_wp, & 243127.0_wp / 20160.0_wp, & -288521.0_wp / 10080.0_wp, 255397.0_wp / 7560.0_wp, -1240.0_wp / 63.0_wp, & 90593.0_wp / 40320.0_wp, & 90593.0_wp / 40320.0_wp, -1240.0_wp / 63.0_wp, 255397.0_wp / 7560.0_wp, & -288521.0_wp / 10080.0_wp, & 243127.0_wp / 20160.0_wp, -12281.0_wp / 6048.0_wp, -1240.0_wp / 63.0_wp, & 1884439.0_wp / 40320.0_wp, & -5106971.0_wp / 30240.0_wp, 248681.0_wp / 1680.0_wp, -643999.0_wp / & 10080.0_wp, 662503.0_wp / 60480.0_wp, & 255397.0_wp / 7560.0_wp, -5106971.0_wp / 30240.0_wp, 4877743.0_wp / & 30240.0_wp, -559651.0_wp / 1890.0_wp, & 1991239.0_wp / 15120.0_wp, -139633.0_wp / 6048.0_wp, -288521.0_wp / & 10080.0_wp, 248681.0_wp / 1680.0_wp, & -559651.0_wp / 1890.0_wp, 159219.0_wp / 1120.0_wp, -1323367.0_wp / & 10080.0_wp, 178999.0_wp / 7560.0_wp, & 243127.0_wp / 20160.0_wp, -643999.0_wp / 10080.0_wp, 1991239.0_wp / & 15120.0_wp, -1323367.0_wp / 10080.0_wp, & 141661.0_wp / 4480.0_wp, -178747.0_wp / 15120.0_wp, -12281.0_wp / 6048.0_wp, & 662503.0_wp / 60480.0_wp, & -139633.0_wp / 6048.0_wp, 178999.0_wp / 7560.0_wp, -178747.0_wp / & 15120.0_wp, 139633.0_wp / 120960.0_wp, & 384187.0_wp / 40320.0_wp, -1139749.0_wp / 15120.0_wp, 61427.0_wp / 504.0_wp, & -1015303.0_wp / 10080.0_wp, & 2567287.0_wp / 60480.0_wp, -73379.0_wp / 10080.0_wp, -1139749.0_wp / & 15120.0_wp, 19365967.0_wp / 120960.0_wp, & -16306061.0_wp / 30240.0_wp, 6881719.0_wp / 15120.0_wp, -5877617.0_wp / & 30240.0_wp, 2033509.0_wp / 60480.0_wp, & 61427.0_wp / 504.0_wp, -16306061.0_wp / 30240.0_wp, 4721851.0_wp / & 10080.0_wp, -169859.0_wp / 210.0_wp, & 5300629.0_wp / 15120.0_wp, -68601.0_wp / 1120.0_wp, -1015303.0_wp / & 10080.0_wp, 6881719.0_wp / 15120.0_wp, & -169859.0_wp / 210.0_wp, 1197047.0_wp / 3360.0_wp, -9478331.0_wp / & 30240.0_wp, 139471.0_wp / 2520.0_wp, & 2567287.0_wp / 60480.0_wp, -5877617.0_wp / 30240.0_wp, 5300629.0_wp / & 15120.0_wp, -9478331.0_wp / 30240.0_wp, & 8449957.0_wp / 120960.0_wp, -188483.0_wp / 7560.0_wp, -73379.0_wp / & 10080.0_wp, 2033509.0_wp / 60480.0_wp, & -68601.0_wp / 1120.0_wp, 139471.0_wp / 2520.0_wp, -188483.0_wp / 7560.0_wp, & 90593.0_wp / 40320.0_wp, & 6150211.0_wp / 120960.0_wp, -2966279.0_wp / 7560.0_wp, 4762921.0_wp / & 7560.0_wp, -15848531.0_wp / 30240.0_wp, & 2706017.0_wp / 12096.0_wp, -235637.0_wp / 6048.0_wp, -2966279.0_wp / & 7560.0_wp, 31617079.0_wp / 40320.0_wp, & -25980937.0_wp / 10080.0_wp, 32862709.0_wp / 15120.0_wp, -1048211.0_wp / & 1120.0_wp, 661145.0_wp / 4032.0_wp, & 4762921.0_wp / 7560.0_wp, -25980937.0_wp / 10080.0_wp, 21703781.0_wp / & 10080.0_wp, -6937561.0_wp / 1890.0_wp, & 2674951.0_wp / 1680.0_wp, -314063.0_wp / 1120.0_wp, -15848531.0_wp / & 30240.0_wp, 32862709.0_wp / 15120.0_wp, & -6937561.0_wp / 1890.0_wp, 47689393.0_wp / 30240.0_wp, -41615261.0_wp / & 30240.0_wp, 1840141.0_wp / 7560.0_wp, & 2706017.0_wp / 12096.0_wp, -1048211.0_wp / 1120.0_wp, 2674951.0_wp / & 1680.0_wp, -41615261.0_wp / 30240.0_wp, & 12160229.0_wp / 40320.0_wp, -539591.0_wp / 5040.0_wp, -235637.0_wp / & 6048.0_wp, 661145.0_wp / 4032.0_wp, & -314063.0_wp / 1120.0_wp, 1840141.0_wp / 7560.0_wp, -539591.0_wp / & 5040.0_wp, 384187.0_wp / 40320.0_wp & ], [6, 6, 6]) public :: weno11_js_reconstruct contains pure subroutine weno11_js_reconstruct(f_stencil, f_hat) real(wp), intent(in) :: f_stencil(:, :) real(wp), intent(out) :: f_hat(:) call apply_weno_js_generic(f_stencil, candidate_coeffs, beta_mats, linear_weights, f_hat) end subroutine weno11_js_reconstruct end module weno11_js