Recon1d_PLM_CWK.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> Piecewise Linear Method 1D reconstruction
6!!
7!! This implementation of PLM follows Colella and Woodward, 1984, except for assuming
8!! uniform cell thicknesses. Cells resort to PCM for extrema including first and last cells in column.
9!! The cell-wise reconstructions are limited so that the edge values (which are also the
10!! extrema in a cell) are bounded by the neighbor cell means. However, this does not yield
11!! monotonic profiles for the whole column.
12!!
13!! Note that internally the edge values, rather than the PLM slope, are stored to ensure
14!! resulting calculations are properly bounded.
15module recon1d_plm_cwk
16
17use recon1d_type, only : testing
18use recon1d_plm_cw, only : plm_cw
19
20implicit none ; private
21
22public plm_cwk, testing
23
24!> PLM reconstruction following Colella and Woodward, 1984
25!!
26!! Implemented by extending recon1d_plm_cw.
27!!
28!! The source for the methods ultimately used by this class are:
29!! - init() -> recon1d_plm_cw.init()
30!! - reconstruct() *locally defined
31!! - average() -> recon1d_plm_cw.average()
32!! - f() -> recon1d_plm_cw.f()
33!! - dfdx() -> recon1d_plm_cw.dfdx()
34!! - x() -> recon1d_plm_cw.x()
35!! - check_reconstruction() -> recon1d_plm_cw.check_reconstruction()
36!! - unit_tests() -> recon1d_plm_cw.unit_tests()
37!! - destroy() -> recon1d_plm_cw.destroy()
38!! - remap_to_sub_grid() -> recon1d_type.remap_to_sub_grid()
39!! - init_parent() -> init()
40!! - reconstruct_parent() -> reconstruct()
41type, extends (plm_cw) :: plm_cwk
42
43contains
44 !> Implementation of the PLM_CWK reconstruction
45 procedure :: reconstruct => reconstruct
46
47end type plm_cwk
48
49contains
50
51!> Calculate a 1D PLM reconstructions based on h(:) and u(:)
52subroutine reconstruct(this, h, u)
53 class(plm_cwk), intent(inout) :: this !< This reconstruction
54 real, intent(in) :: h(*) !< Grid spacing (thickness) [typically H]
55 real, intent(in) :: u(*) !< Cell mean values [A]
56 ! Local variables
57 real :: slp ! The PLM slopes (difference across cell) [A]
58 real :: sigma_l, sigma_c, sigma_r ! Left, central and right slope estimates as
59 ! differences across the cell [A]
60 real :: u_min, u_max ! Minimum and maximum value across cell [A]
61 real :: u_l, u_r, u_c ! Left, right, and center values [A]
62 integer :: k, n
63
64 n = this%n
65
66 ! Loop over all cells
67 do k = 1, n
68 this%u_mean(k) = u(k)
69 enddo
70
71 ! Boundary cells use PCM
72 this%ul(1) = u(1)
73 this%ur(1) = u(1)
74
75 ! Loop over interior cells
76 do k = 2, n-1
77 u_l = u(k-1)
78 u_c = u(k)
79 u_r = u(k+1)
80
81 ! Side differences
82 sigma_r = u_r - u_c
83 sigma_l = u_c - u_l
84
85 ! This is the second order slope given by equation 1.7 of
86 ! Piecewise Parabolic Method, Colella and Woodward (1984),
87 ! but for uniform resolution.
88 sigma_c = 0.5 * ( u_r - u_l )
89
90 ! Limit slope so that reconstructions are bounded by neighbors
91 u_min = min( u_l, u_c, u_r )
92 u_max = max( u_l, u_c, u_r )
93
94 if ( (sigma_l * sigma_r) > 0.0 ) then
95 ! This limits the slope so that the edge values are bounded by the two cell averages spanning the edge
96 slp = sign( min( abs(sigma_c), 2.*min( u_c - u_min, u_max - u_c ) ), sigma_c )
97 else
98 ! Extrema in the mean values require a PCM reconstruction
99 slp = 0.0
100 endif
101
102 ! Left edge
103 u_min = min( u_c, u_l )
104 u_max = max( u_c, u_l )
105 u_l = u_c - 0.5 * slp
106 this%ul(k) = max( min( u_l, u_max), u_min )
107
108 ! Right edge
109 u_min = min( u_c, u_r )
110 u_max = max( u_c, u_r )
111 u_r = u_c + 0.5 * slp
112 this%ur(k) = max( min( u_r, u_max), u_min )
113 enddo
114
115 ! Boundary cells use PCM
116 this%ul(n) = u(n)
117 this%ur(n) = u(n)
118
119end subroutine reconstruct
120
121!> \namespace recon1d_plm_cwk
122!!
123
124end module recon1d_plm_cwk