From e5aa2707e14622ca983f876e474c4922d2898ae7 Mon Sep 17 00:00:00 2001 From: Alistair Adcroft Date: Thu, 8 Aug 2024 12:17:33 -0400 Subject: [PATCH] Replaced lr_edge() and inv_f() with f() and dfdx() --- src/ALE/MOM_remapping.F90 | 1 + src/ALE/Recon1d_EMPLM_WA.F90 | 36 +++---- src/ALE/Recon1d_EMPLM_WA_poly.F90 | 37 +++----- src/ALE/Recon1d_MPLM_WA.F90 | 36 +++---- src/ALE/Recon1d_MPLM_WA_poly.F90 | 41 ++++---- src/ALE/Recon1d_PCM.F90 | 79 ++++++---------- src/ALE/Recon1d_PLM_CW.F90 | 130 +++++++++----------------- src/ALE/Recon1d_PPM_H4_2019.F90 | 150 ++++++++++++------------------ src/ALE/Recon1d_type.F90 | 60 +++++------- 9 files changed, 223 insertions(+), 347 deletions(-) diff --git a/src/ALE/MOM_remapping.F90 b/src/ALE/MOM_remapping.F90 index 95333806e1..47718630f0 100644 --- a/src/ALE/MOM_remapping.F90 +++ b/src/ALE/MOM_remapping.F90 @@ -1723,6 +1723,7 @@ logical function remapping_unit_tests(verbose, num_comp_samp) type(PPM_H4_2019) :: PPM_H4_2019 call test%set( verbose=verbose ) ! Sets the verbosity flag in test +! call test%set( stop_instantly=.true. ) ! While debugging answer_date = 20190101 ! 20181231 h_neglect = hNeglect_dflt diff --git a/src/ALE/Recon1d_EMPLM_WA.F90 b/src/ALE/Recon1d_EMPLM_WA.F90 index 991c505b49..69c272f9bb 100644 --- a/src/ALE/Recon1d_EMPLM_WA.F90 +++ b/src/ALE/Recon1d_EMPLM_WA.F90 @@ -17,13 +17,12 @@ module Recon1d_EMPLM_WA !! The source for the methods ultimately used by this class are: !! init() -> MPLM_WA%init() -> PLM_CW%init() !! reconstruct() *locally defined -!! lr_edge() -> MPLM_WA%lr_edge() -> PLM_CW%lr_edge() !! average() -> MPLM_WA%average() -> PLM_CW%average() -!! inv_f() -> MPLM_WA%inv_f() -> PLM_CW%inv_f() +!! f() -> MPLM_WA%f() -> PLM_CW%f() +!! dfdx() -> MPLM_WA%dfdx() -> PLM_CW%dfdx() !! check_reconstruction() -> MPLM_WA%check_reconstruction() !! unit_tests() *locally defined !! destroy() -> MPLM_WA%destroy() -> PLM_CW%destroy() -!! cell_mean() -> MPLM_WA%cell_mean() -> PLM_CW%cell_mean() -> Recon1d%cell_mean() !! remap_to_sub_grid() -> MPLM_WA%remap_to_sub_grid() -> PLM_CW%remap_to_sub_grd() -> Recon1d%remap_to_sub_grid() !! init_parent() -> MPLM_WA%init_parent() -> PLM_CW%init() !! reconstruct_parent() -> MPLM_WA%reconstruct() @@ -140,35 +139,30 @@ logical function unit_tests(this, verbose, stdout, stderr) call this%reconstruct( (/2.,2.,2./), (/1.,3.,5./) ) call test%real_arr(3, this%u_mean, (/1.,3.,5./), 'Setting cell values') + do k = 1, 3 - um(k) = this%cell_mean(k) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo - call test%real_arr(3, um, (/1.,3.,5./), 'Return cell mean') + call test%real_arr(3, ul, (/0.,2.,4./), 'Evaluation on left edge') + call test%real_arr(3, um, (/1.,3.,5./), 'Evaluation in center') + call test%real_arr(3, ur, (/2.,4.,6./), 'Evaluation on right edge') do k = 1, 3 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%dfdx(k, 0.) + um(k) = this%dfdx(k, 0.5) + ur(k) = this%dfdx(k, 1.) enddo - call test%real_arr(3, ul, (/0.,2.,4./), 'Return left edge') - call test%real_arr(3, ur, (/2.,4.,6./), 'Return right edge') + call test%real_arr(3, ul, (/2.,2.,2./), 'dfdx on left edge') + call test%real_arr(3, um, (/2.,2.,2./), 'dfdx in center') + call test%real_arr(3, ur, (/2.,2.,2./), 'dfdx on right edge') do k = 1, 3 um(k) = this%average(k, 0.5, 0.75) ! Average from x=0.25 to 0.75 in each cell enddo call test%real_arr(3, um, (/1.25,3.25,5.25/), 'Return interval average') - do k = 1, 3 - ull(k) = this%inv_f(k, real(2*k-2)) - ul(k) = this%inv_f(k, real(2*k-1)-0.5) - um(k) = this%inv_f(k, real(2*k-1)) - ur(k) = this%inv_f(k, real(2*k-1)+0.5) - urr(k) = this%inv_f(k, real(2*k-0)) - enddo - call test%real_arr(3, ull, (/0.,0.,0./), 'Return position of f<<') - call test%real_arr(3, ul, (/0.25,0.25,0.25/), 'Return position of f<') - call test%real_arr(3, um, (/0.5,0.5,0.5/), 'Return position of f=') - call test%real_arr(3, ur, (/0.75,0.75,0.75/), 'Return position of f>') - call test%real_arr(3, urr, (/1.,1.,1./), 'Return position of f>>') - unit_tests = test%summarize('EMPLM_WA:unit_tests') end function unit_tests diff --git a/src/ALE/Recon1d_EMPLM_WA_poly.F90 b/src/ALE/Recon1d_EMPLM_WA_poly.F90 index 75f71cb511..68a34c9cb9 100644 --- a/src/ALE/Recon1d_EMPLM_WA_poly.F90 +++ b/src/ALE/Recon1d_EMPLM_WA_poly.F90 @@ -20,14 +20,12 @@ module Recon1d_EMPLM_WA_poly !! The source for the methods ultimately used by this class are: !! init() -> MPLM_WA_poly%init() !! reconstruct() -> MPLM_WA_poly%reconstruct() -!! lr_edge() -> MPLM_WA_poly%lr_edge() -> MPLM_WA_poly%lr_edge() -> PLM_CW%lr_edge() !! average() -> MPLM_WA_poly%average() -!! inv_f() -> MPLM_WA_poly%inv_f() -> MPLM_WA_poly%inv_f() -> PLM_CW%inv_f() +!! f() -> MPLM_WA_poly%f() -> MPLM_WA_poly%f() -> PLM_CW%f() +!! dfdx() -> MPLM_WA_poly%dfdx() -> MPLM_WA_poly%dfdx() -> PLM_CW%dfdx() !! check_reconstruction() -> MPLM_WA_poly%check_reconstruction() !! unit_tests() *locally defined !! destroy() -> MPLM_WA_poly%destroy() -> MPLM_WA%destroy() -> PLM_CW%destroy() -!! cell_mean() -> MPLM_WA_poly%cell_mean() -> MPLM_WA%cell_mean() -> PLM_CW%cell_mean() -!! -> Recon1d%cell_mean() !! remap_to_sub_grid() *locally defined !! init_parent() -> MPLM_WA_poly%init() !! reconstruct_parent() -> MPLM_WA_poly%reconstruct() @@ -120,35 +118,30 @@ logical function unit_tests(this, verbose, stdout, stderr) call this%reconstruct( (/2.,2.,2./), (/1.,3.,5./) ) call test%real_arr(3, this%u_mean, (/1.,3.,5./), 'Setting cell values') + do k = 1, 3 - um(k) = this%cell_mean(k) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo - call test%real_arr(3, um, (/1.,3.,5./), 'Return cell mean') + call test%real_arr(3, ul, (/0.,2.,4./), 'Evaluation on left edge') + call test%real_arr(3, um, (/1.,3.,5./), 'Evaluation in center') + call test%real_arr(3, ur, (/2.,4.,6./), 'Evaluation on right edge') do k = 1, 3 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%dfdx(k, 0.) + um(k) = this%dfdx(k, 0.5) + ur(k) = this%dfdx(k, 1.) enddo - call test%real_arr(3, ul, (/0.,2.,4./), 'Return left edge') - call test%real_arr(3, ur, (/2.,4.,6./), 'Return right edge') + call test%real_arr(3, ul, (/2.,2.,2./), 'dfdx on left edge') + call test%real_arr(3, um, (/2.,2.,2./), 'dfdx in center') + call test%real_arr(3, ur, (/2.,2.,2./), 'dfdx on right edge') do k = 1, 3 um(k) = this%average(k, 0.5, 0.75) ! Average from x=0.25 to 0.75 in each cell enddo call test%real_arr(3, um, (/1.25,3.25,5.25/), 'Return interval average') - do k = 1, 3 - ull(k) = this%inv_f(k, real(2*k-2)) - ul(k) = this%inv_f(k, real(2*k-1)-0.5) - um(k) = this%inv_f(k, real(2*k-1)) - ur(k) = this%inv_f(k, real(2*k-1)+0.5) - urr(k) = this%inv_f(k, real(2*k-0)) - enddo - call test%real_arr(3, ull, (/0.,0.,0./), 'Return position of f<<') - call test%real_arr(3, ul, (/0.25,0.25,0.25/), 'Return position of f<') - call test%real_arr(3, um, (/0.5,0.5,0.5/), 'Return position of f=') - call test%real_arr(3, ur, (/0.75,0.75,0.75/), 'Return position of f>') - call test%real_arr(3, urr, (/1.,1.,1./), 'Return position of f>>') - unit_tests = test%summarize('EMPLM_WA_poly:unit_tests') end function unit_tests diff --git a/src/ALE/Recon1d_MPLM_WA.F90 b/src/ALE/Recon1d_MPLM_WA.F90 index d6c969e7e7..064c79b3dd 100644 --- a/src/ALE/Recon1d_MPLM_WA.F90 +++ b/src/ALE/Recon1d_MPLM_WA.F90 @@ -18,13 +18,12 @@ module Recon1d_MPLM_WA !! The source for the methods ultimately used by this class are: !! init() -> PLM_CW%init() !! reconstruct() *locally defined -!! lr_edge() -> PLM_CW%lr_edge() !! average() -> PLM_CW%average() -!! inv_f() -> PLM_CW%inv_f() +!! f() -> PLM_CW%f() +!! dfdx() -> PLM_CW%dfdx() !! check_reconstruction() *locally defined !! unit_tests() *locally defined !! destroy() -> PLM_CW%destroy() -!! cell_mean() -> PLM_CW%cell_mean() -> Recon1d%cell_mean() !! remap_to_sub_grid() -> PLM_CW%remap_to_sub_grd() -> Recon1d%remap_to_sub_grid() !! init_parent() -> PLM_CW%init() !! reconstruct_parent() -> reconstruct() @@ -234,35 +233,30 @@ logical function unit_tests(this, verbose, stdout, stderr) call this%reconstruct( (/2.,2.,2./), (/1.,3.,5./) ) call test%real_arr(3, this%u_mean, (/1.,3.,5./), 'Setting cell values') + do k = 1, 3 - um(k) = this%cell_mean(k) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo - call test%real_arr(3, um, (/1.,3.,5./), 'Return cell mean') + call test%real_arr(3, ul, (/1.,2.,5./), 'Evaluation on left edge') + call test%real_arr(3, um, (/1.,3.,5./), 'Evaluation in center') + call test%real_arr(3, ur, (/1.,4.,5./), 'Evaluation on right edge') do k = 1, 3 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%dfdx(k, 0.) + um(k) = this%dfdx(k, 0.5) + ur(k) = this%dfdx(k, 1.) enddo - call test%real_arr(3, ul, (/1.,2.,5./), 'Return left edge') - call test%real_arr(3, ur, (/1.,4.,5./), 'Return right edge') + call test%real_arr(3, ul, (/0.,2.,0./), 'dfdx on left edge') + call test%real_arr(3, um, (/0.,2.,0./), 'dfdx in center') + call test%real_arr(3, ur, (/0.,2.,0./), 'dfdx on right edge') do k = 1, 3 um(k) = this%average(k, 0.5, 0.75) ! Average from x=0.25 to 0.75 in each cell enddo call test%real_arr(3, um, (/1.,3.25,5./), 'Return interval average') - do k = 1, 3 - ull(k) = this%inv_f(k, real(2*k-2)) - ul(k) = this%inv_f(k, real(2*k-1)-0.5) - um(k) = this%inv_f(k, real(2*k-1)) - ur(k) = this%inv_f(k, real(2*k-1)+0.5) - urr(k) = this%inv_f(k, real(2*k-0)) - enddo - call test%real_arr(3, ull, (/0.,0.,0./), 'Return position of f<<') - call test%real_arr(3, ul, (/0.,0.25,0./), 'Return position of f<') - call test%real_arr(3, um, (/0.5,0.5,0.5/), 'Return position of f=') - call test%real_arr(3, ur, (/1.,0.75,1./), 'Return position of f>') - call test%real_arr(3, urr, (/1.,1.,1./), 'Return position of f>>') - unit_tests = test%summarize('MPLM_WA:unit_tests') end function unit_tests diff --git a/src/ALE/Recon1d_MPLM_WA_poly.F90 b/src/ALE/Recon1d_MPLM_WA_poly.F90 index 5fd152553b..cd77746365 100644 --- a/src/ALE/Recon1d_MPLM_WA_poly.F90 +++ b/src/ALE/Recon1d_MPLM_WA_poly.F90 @@ -21,13 +21,12 @@ module Recon1d_MPLM_WA_poly !! The source for the methods ultimately used by this class are: !! init() *locally defined !! reconstruct() *locally defined -!! lr_edge() -> MPLM_WA%lr_edge() -> PLM_CW%lr_edge() !! average() *locally defined -!! inv_f() -> MPLM_WA%inv_f() -> PLM_CW%inv_f() +!! f() -> MPLM_WA%f() -> PLM_CW%f() +!! dfdx() -> MPLM_WA%dfdx() -> PLM_CW%dfdx() !! check_reconstruction() *locally defined !! unit_tests() *locally defined !! destroy() -> MPLM_WA%destroy() -> PLM_CW%destroy() -!! cell_mean() -> MPLM_WA%cell_mean() -> PLM_CW%cell_mean() -> Recon1d%cell_mean() !! remap_to_sub_grid() *locally defined !! init_parent() -> init() !! reconstruct_parent() -> reconstruct() @@ -89,7 +88,7 @@ subroutine init(this, n, h_neglect, check) end subroutine init -!> Calculate a 1D PLM reconstructions based on h(:) and u(:) +!> Calculate a 1D MPLM_WA_poly reconstructions based on h(:) and u(:) subroutine reconstruct(this, h, u) class(MPLM_WA_poly), intent(inout) :: this !< This reconstruction real, intent(in) :: h(*) !< Grid spacing (thickness) [typically H] @@ -294,7 +293,8 @@ subroutine remap_to_sub_grid(this, h0, u0, n1, h_sub, & i0_last_thick_cell = 0 do i0 = 1, n0 - call this%lr_edge(i0, ul, ur) + ul = this%ul(i0) + ur = this%ur(i0) u0_min(i0) = min(ul, ur) u0_max(i0) = max(ul, ur) if (h0(i0)>0.) i0_last_thick_cell = i0 @@ -451,35 +451,30 @@ logical function unit_tests(this, verbose, stdout, stderr) call this%reconstruct( (/2.,2.,2./), (/1.,3.,5./) ) call test%real_arr(3, this%u_mean, (/1.,3.,5./), 'Setting cell values') + do k = 1, 3 - um(k) = this%cell_mean(k) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo - call test%real_arr(3, um, (/1.,3.,5./), 'Return cell mean') + call test%real_arr(3, ul, (/1.,2.,5./), 'Evaluation on left edge') + call test%real_arr(3, um, (/1.,3.,5./), 'Evaluation in center') + call test%real_arr(3, ur, (/1.,4.,5./), 'Evaluation on right edge') do k = 1, 3 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%dfdx(k, 0.) + um(k) = this%dfdx(k, 0.5) + ur(k) = this%dfdx(k, 1.) enddo - call test%real_arr(3, ul, (/1.,2.,5./), 'Return left edge') - call test%real_arr(3, ur, (/1.,4.,5./), 'Return right edge') + call test%real_arr(3, ul, (/0.,2.,0./), 'dfdx on left edge') + call test%real_arr(3, um, (/0.,2.,0./), 'dfdx in center') + call test%real_arr(3, ur, (/0.,2.,0./), 'dfdx on right edge') do k = 1, 3 um(k) = this%average(k, 0.5, 0.75) ! Average from x=0.25 to 0.75 in each cell enddo call test%real_arr(3, um, (/1.,3.25,5./), 'Return interval average') - do k = 1, 3 - ull(k) = this%inv_f(k, real(2*k-2)) - ul(k) = this%inv_f(k, real(2*k-1)-0.5) - um(k) = this%inv_f(k, real(2*k-1)) - ur(k) = this%inv_f(k, real(2*k-1)+0.5) - urr(k) = this%inv_f(k, real(2*k-0)) - enddo - call test%real_arr(3, ull, (/0.,0.,0./), 'Return position of f<<') - call test%real_arr(3, ul, (/0.,0.25,0./), 'Return position of f<') - call test%real_arr(3, um, (/0.5,0.5,0.5/), 'Return position of f=') - call test%real_arr(3, ur, (/1.,0.75,1./), 'Return position of f>') - call test%real_arr(3, urr, (/1.,1.,1./), 'Return position of f>>') - unit_tests = test%summarize('MPLM_WA_poly:unit_tests') end function unit_tests diff --git a/src/ALE/Recon1d_PCM.F90 b/src/ALE/Recon1d_PCM.F90 index 764bbdfad1..3b64844983 100644 --- a/src/ALE/Recon1d_PCM.F90 +++ b/src/ALE/Recon1d_PCM.F90 @@ -14,13 +14,12 @@ module Recon1d_PCM !! The source for the methods ultimately used by this class are: !! init() *locally defined !! reconstruct() *locally defined -!! lr_edge() *locally defined !! average() *locally defined -!! inv_f() *locally defined +!! f() *locally defined +!! dfdx() *locally defined !! check_reconstruction() *locally defined !! unit_tests() *locally defined !! destroy() *locally defined -!! cell_mean() -> Recon1d%cell_mean() !! remap_to_sub_grid() -> Recon1d%remap_to_sub_grid() !! init_parent() -> init() !! reconstruct_parent() -> parent() @@ -31,12 +30,12 @@ module Recon1d_PCM procedure :: init => init !> Implementation of the PCM reconstruction procedure :: reconstruct => reconstruct - !> Implementation of function returning the PCM edge values - procedure :: lr_edge => lr_edge !> Implementation of the PCM average over an interval [A] procedure :: average => average - !> Implementation of finding the PCM position of a value - procedure :: inv_f => inv_f + !> Implementation of evaluating the PCM reconstruction at a point [A] + procedure :: f => f + !> Implementation of the derivative of the PCM reconstruction at a point [A] + procedure :: dfdx => dfdx !> Implementation of deallocation for PCM procedure :: destroy => destroy !> Implementation of check reconstruction for the PCM reconstruction @@ -86,40 +85,25 @@ subroutine reconstruct(this, h, u) end subroutine reconstruct -!> Returns the left/right edge values in cell k of a 1D PCM reconstruction -subroutine lr_edge(this, k, u_left, u_right) - class(PCM), intent(in) :: this !< This reconstruction - integer, intent(in) :: k !< Cell number - real, intent(out) :: u_left !< Left edge value [A] - real, intent(out) :: u_right !< Right edge value [A] +!> Value of PCM reconstruction at a point in cell k [A] +real function f(this, k, x) + class(PCM), intent(in) :: this !< This reconstruction + integer, intent(in) :: k !< Cell number + real, intent(in) :: x !< Non-dimensional position within element [nondim] - u_left = this%u_mean(k) - u_right = this%u_mean(k) + f = this%u_mean(k) -end subroutine lr_edge +end function f -!> Position at which 1D PCM reconstruction has a value f in cell k [0..1] -real function inv_f(this, k, f) +!> Derivative of PCM reconstruction at a point in cell k [A] +real function dfdx(this, k, x) class(PCM), intent(in) :: this !< This reconstruction integer, intent(in) :: k !< Cell number - real, intent(in) :: f !< Value of reconstruction to solve for [A] - real :: dfdx_sign ! sign of dfdx, i.e. -1 or +1 [nondim] - real :: df ! Difference in values relative to direction of slope [A] + real, intent(in) :: x !< Non-dimensional position within element [nondim] - ! This is the sign of dfdx across cell k - dfdx_sign = sign( 1., this%u_mean( min(k+1,this%n) ) - this%u_mean( max(k-1,1) ) ) - ! This is the difference between searched for value and cell mean, with possible sign change - df = ( f - this%u_mean(k) ) * dfdx_sign + dfdx = 0. - inv_f = 0.5 ! For PCM, every value between x=0 and x=1 has the same value - ! and x=0.5 is the average position - if ( df < 0. ) then - inv_f = 0. - elseif ( df > 0. ) then - inv_f = 1. - endif - -end function inv_f +end function dfdx !> Average between xa and xb for cell k of a 1D PCM reconstruction [A] real function average(this, k, xa, xb) @@ -178,31 +162,30 @@ logical function unit_tests(this, verbose, stdout, stderr) call this%reconstruct( (/2.,2.,2./), (/1.,3.,5./) ) call test%real_arr(3, this%u_mean, (/1.,3.,5./), 'Setting cell values') + do k = 1, 3 - um(k) = this%cell_mean(k) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo - call test%real_arr(3, um, (/1.,3.,5./), 'Return cell mean') + call test%real_arr(3, ul, (/1.,3.,5./), 'Evaluation on left edge') + call test%real_arr(3, um, (/1.,3.,5./), 'Evaluation in center') + call test%real_arr(3, ur, (/1.,3.,5./), 'Evaluation on right edge') do k = 1, 3 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%dfdx(k, 0.) + um(k) = this%dfdx(k, 0.5) + ur(k) = this%dfdx(k, 1.) enddo - call test%real_arr(3, ul, (/1.,3.,5./), 'Return left edge') - call test%real_arr(3, ur, (/1.,3.,5./), 'Return right edge') + call test%real_arr(3, ul, (/0.,0.,0./), 'dfdx on left edge') + call test%real_arr(3, um, (/0.,0.,0./), 'dfdx in center') + call test%real_arr(3, ur, (/0.,0.,0./), 'dfdx on right edge') do k = 1, 3 um(k) = this%average(k, 0.5, 0.75) enddo call test%real_arr(3, um, (/1.,3.,5./), 'Return interval average') - do k = 1, 3 - ul(k) = this%inv_f(k, real(2*k-2)) - um(k) = this%inv_f(k, real(2*k-1)) - ur(k) = this%inv_f(k, real(2*k-0)) - enddo - call test%real_arr(3, ul, (/0.,0.,0./), 'Return position of f<') - call test%real_arr(3, um, (/0.5,0.5,0.5/), 'Return position of f=') - call test%real_arr(3, ur, (/1.,1.,1./), 'Return position of f>') - unit_tests = test%summarize('PCM:unit_tests') end function unit_tests diff --git a/src/ALE/Recon1d_PLM_CW.F90 b/src/ALE/Recon1d_PLM_CW.F90 index 5b6e71b6a6..7b70ce0006 100644 --- a/src/ALE/Recon1d_PLM_CW.F90 +++ b/src/ALE/Recon1d_PLM_CW.F90 @@ -21,13 +21,12 @@ module Recon1d_PLM_CW !! The source for the methods ultimately used by this class are: !! init() *locally defined !! reconstruct() *locally defined -!! lr_edge() *locally defined !! average() *locally defined -!! inv_f() *locally defined +!! f() *locally defined +!! dfdx() *locally defined !! check_reconstruction() *locally defined !! unit_tests() *locally defined !! destroy() *locally defined -!! cell_mean() -> Recon1d%cell_mean() !! remap_to_sub_grid() -> Recon1d%remap_to_sub_grid() !! init_parent() -> init() !! reconstruct_parent() -> reconstruct() @@ -41,12 +40,12 @@ module Recon1d_PLM_CW procedure :: init => init !> Implementation of the PLM_CW reconstruction procedure :: reconstruct => reconstruct - !> Implementation of function returning the PLM_CW edge values - procedure :: lr_edge => lr_edge !> Implementation of the PLM_CW average over an interval [A] procedure :: average => average - !> Implementation of finding the PLM_CW position of a value - procedure :: inv_f => inv_f + !> Implementation of evaluating the PLM_CW reconstruction at a point [A] + procedure :: f => f + !> Implementation of the derivative of the PLM_CW reconstruction at a point [A] + procedure :: dfdx => dfdx !> Implementation of deallocation for PLM_CW procedure :: destroy => destroy !> Implementation of check reconstruction for the PLM_CW reconstruction @@ -164,52 +163,35 @@ subroutine reconstruct(this, h, u) end subroutine reconstruct -!> Returns the left/right edge values in cell k of a 1D PLM reconstruction -subroutine lr_edge(this, k, u_left, u_right) - class(PLM_CW), intent(in) :: this !< This reconstruction - integer, intent(in) :: k !< Cell number - real, intent(out) :: u_left !< Left edge value [A] - real, intent(out) :: u_right !< Right edge value [A] +!> Value of PLM_CW reconstruction at a point in cell k [A] +real function f(this, k, x) + class(PLM_CW), intent(in) :: this !< This reconstruction + integer, intent(in) :: k !< Cell number + real, intent(in) :: x !< Non-dimensional position within element [nondim] + real :: du ! Difference across cell [A] + real :: u_a, u_b ! Two estimate of f [A] - u_left = this%ul(k) - u_right = this%ur(k) + du = this%ur(k) - this%ul(k) -end subroutine lr_edge + ! This expression for u_a can overshoot u_r but is good for x<<1 + u_a = this%ul(k) + du * x + ! This expression for u_b can overshoot u_l but is good for 1-x<<1 + u_b = this%ur(k) + du * ( x - 1. ) -!> Position at which 1D PLM reconstruction has a value f in cell k [0..1] -real function inv_f(this, k, f) + ! Since u_a and u_b are both bounded, this will perserve uniformity + f = 0.5 * ( u_a + u_b ) + +end function f + +!> Derivative of PLM_CW reconstruction at a point in cell k [A] +real function dfdx(this, k, x) class(PLM_CW), intent(in) :: this !< This reconstruction integer, intent(in) :: k !< Cell number - real, intent(in) :: f !< Value of reconstruction to solve for [A] - real :: dfdx_sign ! sign of dfdx, i.e. -1 or +1 [nondim] - real :: sf ! Values with potentially flipped sign [A] - real :: du ! Difference across cell [A] + real, intent(in) :: x !< Non-dimensional position within element [nondim] - du = this%ur(k) - this%ul(k) + dfdx = this%ur(k) - this%ul(k) - ! This is the sign of dfdx across cell k - if (abs(du) <= 0.) then - dfdx_sign = sign( 1., this%u_mean(min(this%n,k+1)) - this%u_mean(max(1,k-1)) ) - else - dfdx_sign = sign( 1., du ) - endif - sf = f * dfdx_sign - - inv_f = 0.5 ! Avoid compiler warnings - - if ( sf < this%ul(k) * dfdx_sign ) then - inv_f = 0. - elseif ( sf > this%ur(k) * dfdx_sign ) then - inv_f = 1. - elseif ( du == 0. ) then - inv_f = 0.5 ! Same as for PCM - else ! ul < f < ur or ul > f > ur - ! Note that if ur=ul (i.e. ur-ul=0 then we would meet one of the previous - ! conditions that avoid division by zero - inv_f = ( f - this%ul(k) ) / du - endif - -end function inv_f +end function dfdx !> Average between xa and xb for cell k of a 1D PLM reconstruction [A] real function average(this, k, xa, xb) @@ -314,53 +296,30 @@ logical function unit_tests(this, verbose, stdout, stderr) call this%reconstruct( (/2.,2.,2./), (/1.,3.,5./) ) call test%real_arr(3, this%u_mean, (/1.,3.,5./), 'Setting cell values') + do k = 1, 3 - um(k) = this%cell_mean(k) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo - call test%real_arr(3, um, (/1.,3.,5./), 'Return cell mean') + call test%real_arr(3, ul, (/1.,2.,5./), 'Evaluation on left edge') + call test%real_arr(3, um, (/1.,3.,5./), 'Evaluation in center') + call test%real_arr(3, ur, (/1.,4.,5./), 'Evaluation on right edge') do k = 1, 3 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%dfdx(k, 0.) + um(k) = this%dfdx(k, 0.5) + ur(k) = this%dfdx(k, 1.) enddo - call test%real_arr(3, ul, (/1.,2.,5./), 'Return left edge') - call test%real_arr(3, ur, (/1.,4.,5./), 'Return right edge') + call test%real_arr(3, ul, (/0.,2.,0./), 'dfdx on left edge') + call test%real_arr(3, um, (/0.,2.,0./), 'dfdx in center') + call test%real_arr(3, ur, (/0.,2.,0./), 'dfdx on right edge') do k = 1, 3 um(k) = this%average(k, 0.5, 0.75) ! Average from x=0.25 to 0.75 in each cell enddo call test%real_arr(3, um, (/1.,3.25,5./), 'Return interval average') - ! With f(x) = 2*k-1, spacing h=2, find cell-relative positions - do k = 1, 3 - ull(k) = this%inv_f(k, real(2*k-1)-1.0) - ul(k) = this%inv_f(k, real(2*k-1)-0.5) - um(k) = this%inv_f(k, real(2*k-1)) - ur(k) = this%inv_f(k, real(2*k-1)+0.5) - urr(k) = this%inv_f(k, real(2*k-1)+1.0) - enddo - call test%real_arr(3, ull, (/0.,0.,0./), 'Return position of f<<') - call test%real_arr(3, ul, (/0.,0.25,0./), 'Return position of f<') - call test%real_arr(3, um, (/0.5,0.5,0.5/), 'Return position of f=') - call test%real_arr(3, ur, (/1.,0.75,1./), 'Return position of f>') - call test%real_arr(3, urr, (/1.,1.,1./), 'Return position of f>>') - - ! With f(x) = -(2*k-1), find cell-relative positions - call this%reconstruct( (/2.,2.,2./), (/-1.,-3.,-5./) ) - call test%real_arr(3, this%ul, (/-1.,-2.,-5./), 'Return left edge -f') - call test%real_arr(3, this%ur, (/-1.,-4.,-5./), 'Return right edge -f') - do k = 1, 3 - ull(k) = this%inv_f(k, -real(2*k-1)+1.0) - ul(k) = this%inv_f(k, -real(2*k-1)+0.5) - um(k) = this%inv_f(k, -real(2*k-1)) - ur(k) = this%inv_f(k, -real(2*k-1)-0.5) - urr(k) = this%inv_f(k, -real(2*k-1)-1.0) - enddo - call test%real_arr(3, ull, (/0.,0.,0./), 'Return position of f<< -f') - call test%real_arr(3, ul, (/0.,0.25,0./), 'Return position of f< -f') - call test%real_arr(3, um, (/0.5,0.5,0.5/), 'Return position of f= -f') - call test%real_arr(3, ur, (/1.,0.75,1./), 'Return position of f> -f') - call test%real_arr(3, urr, (/1.,1.,1./), 'Return position of f>> -f') - call this%destroy() deallocate( um, ul, ur, ull, urr ) @@ -375,10 +334,11 @@ logical function unit_tests(this, verbose, stdout, stderr) ! are bounded by neighbors but ur(2) and ul(3) are out-of-order. call this%reconstruct( (/1.,1.,1.,1./), (/0.,3.,4.,7./) ) do k = 1, 4 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%f(k, 0.) + ur(k) = this%f(k, 1.) enddo - call test%real_arr(4, ul, (/0.,2.,3.,7./), 'Return left edge') - call test%real_arr(4, ur, (/0.,4.,5.,7./), 'Return right edge') + call test%real_arr(4, ul, (/0.,2.,3.,7./), 'Evaluation on left edge') + call test%real_arr(4, ur, (/0.,4.,5.,7./), 'Evaluation on right edge') deallocate( um, ul, ur ) diff --git a/src/ALE/Recon1d_PPM_H4_2019.F90 b/src/ALE/Recon1d_PPM_H4_2019.F90 index 74aac3e0f9..e031fd5d01 100644 --- a/src/ALE/Recon1d_PPM_H4_2019.F90 +++ b/src/ALE/Recon1d_PPM_H4_2019.F90 @@ -21,13 +21,12 @@ module Recon1d_PPM_H4_2019 !! The source for the methods ultimately used by this class are: !! init() *locally defined !! reconstruct() *locally defined -!! lr_edge() *locally defined !! average() *locally defined -!! inv_f() *locally defined +!! f() *locally defined +!! dfdx() *locally defined !! check_reconstruction() *locally defined !! unit_tests() *locally defined !! destroy() *locally defined -!! cell_mean() -> Recon1d%cell_mean() !! remap_to_sub_grid() -> Recon1d%remap_to_sub_grid() !! init_parent() -> init() !! reconstruct_parent() -> reconstruct() @@ -41,12 +40,12 @@ module Recon1d_PPM_H4_2019 procedure :: init => init !> Implementation of the PPM_H4_2019 reconstruction procedure :: reconstruct => reconstruct - !> Implementation of function returning the PPM_H4_2019 edge values - procedure :: lr_edge => lr_edge !> Implementation of the PPM_H4_2019 average over an interval [A] procedure :: average => average - !> Implementation of finding the PPM_H4_2019 position of a value - procedure :: inv_f => inv_f + !> Implementation of evaluating the PPM_H4_2019 reconstruction at a point [A] + procedure :: f => f + !> Implementation of the derivative of the PPM_H4_2019 reconstruction at a point [A] + procedure :: dfdx => dfdx !> Implementation of deallocation for PPM_H4_2019 procedure :: destroy => destroy !> Implementation of check reconstruction for the PPM_H4_2019 reconstruction @@ -115,69 +114,46 @@ subroutine reconstruct(this, h, u) end subroutine reconstruct -!> Returns the left/right edge values in cell k of a 1D PPM reconstruction -subroutine lr_edge(this, k, u_left, u_right) - class(PPM_H4_2019), intent(in) :: this !< This reconstruction - integer, intent(in) :: k !< Cell number - real, intent(out) :: u_left !< Left edge value [A] - real, intent(out) :: u_right !< Right edge value [A] +!> Value of PPM_H4_2019 reconstruction at a point in cell k [A] +real function f(this, k, x) + class(PPM_H4_2019), intent(in) :: this !< This reconstruction + integer, intent(in) :: k !< Cell number + real, intent(in) :: x !< Non-dimensional position within element [nondim] + real :: du ! Difference across cell [A] + real :: a6 ! Collela and Woordward curvature parameter [A] + real :: u_a, u_b ! Two estimate of f [A] + real :: lmx ! 1 - x [nondim] + real :: wb ! Weight based on x [nondim] + + du = this%ur(k) - this%ul(k) + a6 = 3.0 * ( ( this%u_mean(k) - this%ul(k) ) + ( this%u_mean(k) - this%ur(k) ) ) + lmx = 1.0 - x - u_left = this%ul(k) - u_right = this%ur(k) + ! This expression for u_a can overshoot u_r but is good for x<<1 + u_a = this%ul(k) + x * ( du + a6 * x ) + ! This expression for u_b can overshoot u_l but is good for 1-x<<1 + u_b = this%ur(k) + lmx * ( - du + a6 * lmx ) -end subroutine lr_edge + ! Since u_a and u_b are both side-bounded, using weights=0 or 1 will preserve uniformity + wb = 0.5 + sign(0.5, x - 0.5 ) ! = 1 @ x=0, = 0 @ x=1 + f = ( ( 1. - wb ) * u_a ) + ( wb * u_b ) -!> Position at which 1D PPM_H4_2019 reconstruction has a value f in cell k [0..1] -real function inv_f(this, k, f) +end function f + +!> Derivative of PPM_H4_2019 reconstruction at a point in cell k [A] +real function dfdx(this, k, x) class(PPM_H4_2019), intent(in) :: this !< This reconstruction - integer, intent(in) :: k !< Cell number - real, intent(in) :: f !< Value of reconstruction to solve for [A] - real :: dfdx_sign ! sign of dfdx, i.e. -1 or +1 [nondim] - real :: sf ! Values with potentially flipped sign [A] + integer, intent(in) :: k !< Cell number + real, intent(in) :: x !< Non-dimensional position within element [nondim] real :: du ! Difference across cell [A] - real :: a6 ! The Colella and Woodward PPM curviture parameter [A] - real :: qA, qB, qC ! Respective multiplier of x^2, x, 1 in quadratic [nondim] - real :: discr ! Discriminant in quadratic problem [nondim] + real :: a6 ! Collela and Woordward curvature parameter [A] du = this%ur(k) - this%ul(k) + a6 = 3.0 * ( ( this%u_mean(k) - this%ul(k) ) + ( this%u_mean(k) - this%ur(k) ) ) - ! This is the sign of dfdx across cell k - if (abs(du) <= 0.) then - dfdx_sign = sign( 1., this%u_mean(min(this%n,k+1)) - this%u_mean(max(1,k-1)) ) - else - dfdx_sign = sign( 1., du ) - endif - sf = f * dfdx_sign - - inv_f = 0.5 ! Avoid compiler warnings - - if ( sf < this%ul(k) * dfdx_sign ) then - inv_f = 0. - elseif ( sf > this%ur(k) * dfdx_sign ) then - inv_f = 1. - elseif ( du == 0. ) then - inv_f = 0.5 ! Same as for PCM - else ! ul < f < ur or ul > f > ur - ! Invert PPM since we know f is between ul and ur - a6 = 3.0 * ( ( this%u_mean(k) - this%ul(k) ) + ( this%u_mean(k) - this%ur(k) ) ) - if ( abs(a6) > 0. ) then ! u_mean not equal to (ul+ur)/2 - qA = - a6 - qB = du + a6 - qC = this%ul(k) - f - discr = sqrt( (qB * qB) - 4. * ( qA * qC ) ) - if ( qB >= 0. ) then - inv_f = - 0.5 * ( qB + discr ) / qA - inv_f = - 2.0 * qC / ( qB + discr ) - else ! qB < 0 - inv_f = - 0.5 * ( qB + discr ) / qA - inv_f = 2.0 * qC / ( - qB + discr ) - endif - else ! a6==0 - inv_f = ( f - this%ul(k) ) / du - endif - endif + dfdx = du + a6 * ( 2.0 * x - 1.0 ) -end function inv_f +end function dfdx !> Average between xa and xb for cell k of a 1D PPM reconstruction [A] real function average(this, k, xa, xb) @@ -274,37 +250,38 @@ logical function unit_tests(this, verbose, stdout, stderr) call test%test( this%n /= 5, 'Setting number of levels') allocate( um(5), ul(5), ur(5), ull(5), urr(5) ) + ! Straight line, f(x) = x , or f(K) = 2*K call this%reconstruct( (/2.,2.,2.,2.,2./), (/1.,3.,5.,7.,9./) ) call test%real_arr(5, this%u_mean, (/1.,3.,5.,7.,9./), 'Setting cell values') + call test%real_arr(5, this%ul, (/1.,2.,4.,6.,9./), 'Left edge values', robits=2) + call test%real_arr(5, this%ur, (/1.,4.,6.,8.,9./), 'Right edge values') do k = 1, 5 - um(k) = this%cell_mean(k) + um(k) = this%u_mean(k) enddo call test%real_arr(5, um, (/1.,3.,5.,7.,9./), 'Return cell mean') do k = 1, 5 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo - ! Values of 2 at k=2 differs by roundoff - call test%real_arr(5, ul, (/1.,2.,4.,6.,9./), 'Return left edge', robits=2) - call test%real_arr(5, ur, (/1.,4.,6.,8.,9./), 'Return right edge') + call test%real_arr(5, ul, this%ul, 'Evaluation on left edge') + call test%real_arr(5, um, (/1.,3.,5.,7.,9./), 'Evaluation in center') + call test%real_arr(5, ur, this%ur, 'Evaluation on right edge') do k = 1, 5 - um(k) = this%average(k, 0.5, 0.75) ! Average from x=0.25 to 0.75 in each cell + ul(k) = this%dfdx(k, 0.) + um(k) = this%dfdx(k, 0.5) + ur(k) = this%dfdx(k, 1.) enddo - call test%real_arr(5, um, (/1.,3.25,5.25,7.25,9./), 'Return interval average') + call test%real_arr(5, ul, (/0.,2.,2.,2.,0./), 'dfdx on left edge', robits=3) + call test%real_arr(5, um, (/0.,2.,2.,2.,0./), 'dfdx in center', robits=2) + call test%real_arr(5, ur, (/0.,2.,2.,2.,0./), 'dfdx on right edge', robits=6) do k = 1, 5 - ull(k) = this%inv_f(k, real(2*k-1)-1.0) - ul(k) = this%inv_f(k, real(2*k-1)-0.5) - um(k) = this%inv_f(k, real(2*k-1)) - ur(k) = this%inv_f(k, real(2*k-1)+0.5) - urr(k) = this%inv_f(k, real(2*k-1)+1.0) + um(k) = this%average(k, 0.5, 0.75) ! Average from x=0.25 to 0.75 in each cell enddo - call test%real_arr(5, ull, (/0.,0.,0.,0.,0./), 'Return position of f<<', tol=4e-16) - call test%real_arr(5, ul, (/0.,0.25,0.25,0.25,0./), 'Return position of f<', robits=4) - call test%real_arr(5, um, (/0.5,0.5,0.5,0.5,0.5/), 'Return position of f=', robits=1) - call test%real_arr(5, ur, (/1.,0.75,0.75,0.75,1./), 'Return position of f>', robits=1) - call test%real_arr(5, urr, (/1.,1.,1.,1.,1./), 'Return position of f>>') + call test%real_arr(5, um, (/1.,3.25,5.25,7.25,9./), 'Return interval average') if (verbose) write(stdout,'(a)') 'PPM_H4_2019:unit_tests testing with parabola' @@ -315,24 +292,13 @@ logical function unit_tests(this, verbose, stdout, stderr) ! edges: 0, 3, 12, 27, 48, 75 call this%reconstruct( (/3.,3.,3.,3.,3./), (/1.,7.,19.,37.,61./) ) do k = 1, 5 - call this%lr_edge(k, ul(k), ur(k)) + ul(k) = this%f(k, 0.) + um(k) = this%f(k, 0.5) + ur(k) = this%f(k, 1.) enddo call test%real_arr(5, ul, (/1.,3.,12.,27.,61./), 'Return left edge', robits=2) call test%real_arr(5, ur, (/1.,12.,27.,48.,61./), 'Return right edge', robits=1) -! do k = 1, 5 -! ull(k) = this%inv_f(k, 3.*(real(k)-1.0)**2 ) -! ul(k) = this%inv_f(k, 3.*(real(k)-0.75)**2 ) -! um(k) = this%inv_f(k, 3.*(real(k)-0.5)**2 ) -! ur(k) = this%inv_f(k, (real(k)-0.5)**2 ) -! urr(k) = this%inv_f(k, (real(k)-0.5)**2 ) -! enddo -! call test%real_arr(5, ull, (/0.,0.,0.,0.,0./), 'Return position of f<<', tol=5e-16) -! call test%real_arr(5, ul, (/0.,0.25,0.25,0.25,0./), 'Return position of f<', robits=2) -! call test%real_arr(5, um, (/0.5,0.5,0.5,0.5,0.5/), 'Return position of f=', robits=1) -! call test%real_arr(5, ur, (/1.,0.75,0.75,0.75,1./), 'Return position of f>', robits=1) -! call test%real_arr(5, urr, (/1.,1.,1.,1.,1./), 'Return position of f>>') - call this%destroy() deallocate( um, ul, ur, ull, urr ) diff --git a/src/ALE/Recon1d_type.F90 b/src/ALE/Recon1d_type.F90 index f581d018c4..0fd7d49075 100644 --- a/src/ALE/Recon1d_type.F90 +++ b/src/ALE/Recon1d_type.F90 @@ -25,12 +25,12 @@ module Recon1d_type procedure(i_init), deferred :: init !> Deferred implementation of reconstruction function procedure(i_reconstruct), deferred :: reconstruct - !> Deferred implementation of function returning the edge values - procedure(i_lr_edge), deferred :: lr_edge !> Deferred implementation of the average over an interval procedure(i_average), deferred :: average - !> Deferred implementation of the finding position of a value - procedure(i_inv_f), deferred :: inv_f + !> Deferred implementation of evaluating the reconstruction at a point + procedure(i_f), deferred :: f + !> Deferred implementation of the derivative of the reconstruction at a point + procedure(i_dfdx), deferred :: dfdx !> Deferred implementation of check_reconstruction procedure(i_check_reconstruction), deferred :: check_reconstruction !> Deferred implementation of unit tests for the reconstruction @@ -41,8 +41,6 @@ module Recon1d_type ! The following functions/subroutines are shared across all reconstructions and provided by this module ! unless replaced for the purpose of optimization - !> Cell mean of cell k [A] - procedure :: cell_mean => a_cell_mean !> Remaps the column to subgrid h_sub procedure :: remap_to_sub_grid => remap_to_sub_grid !> Set debugging @@ -58,7 +56,7 @@ module Recon1d_type end type Recon1d -!> Class to assist in unit tests +!> Class to assist in unit tests, not to be used outside of Recon1d types type :: testing private !> True if any fail has been encountered since this instance of "testing" was created @@ -110,15 +108,6 @@ subroutine i_reconstruct(this, h, u) real, intent(in) :: u(*) !< Cell mean values [A] end subroutine i_reconstruct - !> Returns the left/right edge values in cell k of a 1D reconstruction - subroutine i_lr_edge(this, k, u_left, u_right) - import :: Recon1d - class(Recon1d), intent(in) :: this !< This reconstruction - integer, intent(in) :: k !< Cell number - real, intent(out) :: u_left !< Left edge value [A] - real, intent(out) :: u_right !< Right edge value [A] - end subroutine i_lr_edge - !> Average between xa and xb for cell k of a 1D reconstruction [A] !! !! It is assumed that 0<=xa<=1, 0<=xb<=1, and xa<=xb @@ -130,15 +119,25 @@ real function i_average(this, k, xa, xb) real, intent(in) :: xb !< End of averaging interval on element (0 to 1) end function i_average - !> Position at which 1D reconstruction has a value f [0..1] + !> Point-wise value of reconstruction [A] + !! + !! THe function is only valid for 0 <= x <= 1. x is effectively clipped to this range. + real function i_f(this, k, x) + import :: Recon1d + class(Recon1d), intent(in) :: this !< This reconstruction + integer, intent(in) :: k !< Cell number + real, intent(in) :: x !< Non-dimensional position within element [nondim] + end function i_f + + !> Point-wise value of derivative reconstruction [A] !! - !! It is assumed the reconstruction is not decreasing with cell index. - real function i_inv_f(this, k, f) + !! THe function is only valid for 0 <= x <= 1. x is effectively clipped to this range. + real function i_dfdx(this, k, x) import :: Recon1d class(Recon1d), intent(in) :: this !< This reconstruction integer, intent(in) :: k !< Cell number - real, intent(in) :: f !< Value of reconstruction to solve for [A] - end function i_inv_f + real, intent(in) :: x !< Non-dimensional position within element [nondim] + end function i_dfdx !> Returns true if some inconsistency is detected, false otherwise !! @@ -189,21 +188,12 @@ end function i_unit_tests contains -!> Mean value for cell k [A] -real function a_cell_mean(this, k) - class(Recon1d), intent(in) :: this !< This reconstruction - integer, intent(in) :: k !< Cell number - - a_cell_mean = this%u_mean(k) - -end function a_cell_mean - !> Remaps the column to subgrid h_sub !! !! It is assumed that h_sub is a perfect sub-grid of h0, meaning each h0 cell !! can be constructed by joining a contiguous set of h_sub cells. The integer !! indices isrc_start, isrc_end, isub_src provide this mapping, and are -!! calcualted in MOM_remapping +!! calculated in MOM_remapping subroutine remap_to_sub_grid(this, h0, u0, n1, h_sub, & isrc_start, isrc_end, isrc_max, isub_src, & u_sub, uh_sub, u02_err) @@ -338,7 +328,7 @@ subroutine test(this, state, label, ignore) class(testing), intent(inout) :: this !< This testing class logical, intent(in) :: state !< True to indicate a fail, false otherwise character(len=*), intent(in) :: label !< Message - logical, optional, intent(in) :: ignore !< If presnt and true, ignore a fail + logical, optional, intent(in) :: ignore !< If present and true, ignore a fail ! Local variables logical :: ignore_this_fail @@ -420,7 +410,7 @@ subroutine real_arr(this, n, u_test, u_true, label, tol, robits, ignore) character(len=*), intent(in) :: label !< Message real, optional, intent(in) :: tol !< The tolerance for differences between u and u_true [A] integer, optional, intent(in) :: robits !< Number of bits of round-off to allow - logical, optional, intent(in) :: ignore !< If presnt and true, ignore a fail + logical, optional, intent(in) :: ignore !< If present and true, ignore a fail ! Local variables integer :: k logical :: this_test, ignore_this_fail @@ -473,7 +463,7 @@ subroutine int_arr(this, n, i_test, i_true, label, ignore) integer, dimension(n), intent(in) :: i_test !< Values to test [A] integer, dimension(n), intent(in) :: i_true !< Values to test against (correct answer) [A] character(len=*), intent(in) :: label !< Message - logical, optional, intent(in) :: ignore !< If presnt and true, ignore a fail + logical, optional, intent(in) :: ignore !< If present and true, ignore a fail ! Local variables integer :: k logical :: this_test, ignore_this_fail @@ -524,5 +514,5 @@ end subroutine int_arr !! but reconstruc_parent() is not redefined so that it still is defined by Recon1d_XYZ !! !! The module also defines a type "testing" which is used to abbreviate the unit_tests() member. -!! This should no be used outside of these derived reconstructions. +!! This should no be used outside of these reconstruction types. end module Recon1d_type