MOM6
P1M_functions.F90
Go to the documentation of this file.
1 !> Linear interpolation functions
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
7 
8 implicit none ; private
9 
10 ! The following routines are visible to the outside world
12 
13 contains
14 
15 !> Linearly interpolate between edge values
16 !!
17 !! The resulting piecewise interpolant is stored in 'ppoly'.
18 !! See 'ppoly.F90' for a definition of this structure.
19 !!
20 !! The edge values MUST have been estimated prior to calling this routine.
21 !!
22 !! The estimated edge values must be limited to ensure monotonicity of the
23 !! interpolant. We also make sure that edge values are NOT discontinuous.
24 !!
25 !! It is assumed that the size of the array 'u' is equal to the number of cells
26 !! defining 'grid' and 'ppoly'. No consistency check is performed here.
27 subroutine p1m_interpolation( N, h, u, ppoly_E, ppoly_coef, h_neglect )
28  integer, intent(in) :: n !< Number of cells
29  real, dimension(:), intent(in) :: h !< cell widths (size N)
30  real, dimension(:), intent(in) :: u !< cell average properties (size N)
31  real, dimension(:,:), intent(inout) :: ppoly_e !< Potentially modified edge values,
32  !! with the same units as u.
33  real, dimension(:,:), intent(inout) :: ppoly_coef !< Potentially modified
34  !! piecewise polynomial coefficients, mainly
35  !! with the same units as u.
36  real, optional, intent(in) :: h_neglect !< A negligibly small width
37  !! in the same units as h.
38  ! Local variables
39  integer :: k ! loop index
40  real :: u0_l, u0_r ! edge values (left and right)
41 
42  ! Bound edge values (routine found in 'edge_values.F90')
43  call bound_edge_values( n, h, u, ppoly_e, h_neglect )
44 
45  ! Systematically average discontinuous edge values (routine found in
46  ! 'edge_values.F90')
47  call average_discontinuous_edge_values( n, ppoly_e )
48 
49  ! Loop on interior cells to build interpolants
50  do k = 1,n
51 
52  u0_l = ppoly_e(k,1)
53  u0_r = ppoly_e(k,2)
54 
55  ppoly_coef(k,1) = u0_l
56  ppoly_coef(k,2) = u0_r - u0_l
57 
58  enddo ! end loop on interior cells
59 
60 end subroutine p1m_interpolation
61 
62 !> Interpolation by linear polynomials within boundary cells
63 !!
64 !! The left and right edge values in the left and right boundary cells,
65 !! respectively, are estimated using a linear extrapolation within the cells.
66 !!
67 !! It is assumed that the size of the array 'u' is equal to the number of cells
68 !! defining 'grid' and 'ppoly'. No consistency check is performed here.
69 subroutine p1m_boundary_extrapolation( N, h, u, ppoly_E, ppoly_coef )
70  ! Arguments
71  integer, intent(in) :: n !< Number of cells
72  real, dimension(:), intent(in) :: h !< cell widths (size N)
73  real, dimension(:), intent(in) :: u !< cell averages (size N)
74  real, dimension(:,:), intent(inout) :: ppoly_e !< edge values of piecewise polynomials,
75  !! with the same units as u.
76  real, dimension(:,:), intent(inout) :: ppoly_coef !< coefficients of piecewise polynomials, mainly
77  !! with the same units as u.
78  ! Local variables
79  real :: u0, u1 ! cell averages
80  real :: h0, h1 ! corresponding cell widths
81  real :: slope ! retained PLM slope
82  real :: u0_l, u0_r ! edge values
83 
84  ! -----------------------------------------
85  ! Left edge value in the left boundary cell
86  ! -----------------------------------------
87  h0 = h(1)
88  h1 = h(2)
89 
90  u0 = u(1)
91  u1 = u(2)
92 
93  ! The standard PLM slope is computed as a first estimate for the
94  ! interpolation within the cell
95  slope = 2.0 * ( u1 - u0 )
96 
97  ! The right edge value is then computed and we check whether this
98  ! right edge value is consistent: it cannot be larger than the edge
99  ! value in the neighboring cell if the data set is increasing.
100  ! If the right value is found to too large, the slope is further limited
101  ! by using the edge value in the neighboring cell.
102  u0_r = u0 + 0.5 * slope
103 
104  if ( (u1 - u0) * (ppoly_e(2,1) - u0_r) < 0.0 ) then
105  slope = 2.0 * ( ppoly_e(2,1) - u0 )
106  endif
107 
108  ! Using the limited slope, the left edge value is reevaluated and
109  ! the interpolant coefficients recomputed
110  if ( h0 /= 0.0 ) then
111  ppoly_e(1,1) = u0 - 0.5 * slope
112  else
113  ppoly_e(1,1) = u0
114  endif
115 
116  ppoly_coef(1,1) = ppoly_e(1,1)
117  ppoly_coef(1,2) = ppoly_e(1,2) - ppoly_e(1,1)
118 
119  ! ------------------------------------------
120  ! Right edge value in the left boundary cell
121  ! ------------------------------------------
122  h0 = h(n-1)
123  h1 = h(n)
124 
125  u0 = u(n-1)
126  u1 = u(n)
127 
128  slope = 2.0 * ( u1 - u0 )
129 
130  u0_l = u1 - 0.5 * slope
131 
132  if ( (u1 - u0) * (u0_l - ppoly_e(n-1,2)) < 0.0 ) then
133  slope = 2.0 * ( u1 - ppoly_e(n-1,2) )
134  endif
135 
136  if ( h1 /= 0.0 ) then
137  ppoly_e(n,2) = u1 + 0.5 * slope
138  else
139  ppoly_e(n,2) = u1
140  endif
141 
142  ppoly_coef(n,1) = ppoly_e(n,1)
143  ppoly_coef(n,2) = ppoly_e(n,2) - ppoly_e(n,1)
144 
145 end subroutine p1m_boundary_extrapolation
146 
147 !> \namespace p1m_functions
148 !!
149 !! Date of creation: 2008.06.09
150 !! L. White
151 !!
152 !! This module contains p1m (linear) interpolation routines.
153 !!
154 !! p1m interpolation is performed by estimating the edge values and
155 !! linearly interpolating between them.
156 !
157 !! Once the edge values are estimated, the limiting process takes care of
158 !! ensuring that (1) edge values are bounded by neighoring cell averages
159 !! and (2) discontinuous edge values are averaged in order to provide a
160 !! fully continuous interpolant throughout the domain. This last step is
161 !! essential for the regridding problem to yield a unique solution.
162 !! Also, a routine is provided that takes care of linear extrapolation
163 !! within the boundary cells.
164 
165 end module p1m_functions
regrid_edge_values::average_discontinuous_edge_values
subroutine, public average_discontinuous_edge_values(N, edge_val)
Replace discontinuous collocated edge values with their average.
Definition: regrid_edge_values.F90:141
regrid_edge_values
Edge value estimation for high-order resconstruction.
Definition: regrid_edge_values.F90:2
p1m_functions
Linear interpolation functions.
Definition: P1M_functions.F90:2
p1m_functions::p1m_boundary_extrapolation
subroutine, public p1m_boundary_extrapolation(N, h, u, ppoly_E, ppoly_coef)
Interpolation by linear polynomials within boundary cells.
Definition: P1M_functions.F90:70
regrid_edge_values::bound_edge_values
subroutine, public bound_edge_values(N, h, u, edge_val, h_neglect)
Bound edge values by neighboring cell averages.
Definition: regrid_edge_values.F90:48
p1m_functions::p1m_interpolation
subroutine, public p1m_interpolation(N, h, u, ppoly_E, ppoly_coef, h_neglect)
Linearly interpolate between edge values.
Definition: P1M_functions.F90:28