MOM6
Functions/Subroutines | Variables
ppm_functions Module Reference

Detailed Description

Provides functions used with the Piecewise-Parabolic-Method in the vertical ALE algorithm.

Functions/Subroutines

subroutine, public ppm_reconstruction (N, h, u, edge_values, ppoly_coef, h_neglect, answers_2018)
 Builds quadratic polynomials coefficients from cell mean and edge values. More...
 
subroutine ppm_limiter_standard (N, h, u, edge_values, h_neglect, answers_2018)
 Adjusts edge values using the standard PPM limiter (Colella & Woodward, JCP 1984) after first checking that the edge values are bounded by neighbors cell averages and that the edge values are monotonic between cell averages. More...
 
subroutine, public ppm_boundary_extrapolation (N, h, u, edge_values, ppoly_coef, h_neglect)
 Reconstruction by parabolas within boundary cells. More...
 

Variables

real, parameter hneglect_dflt = 1.E-30
 A tiny width that is so small that adding it to cell widths does not change the value due to a computational representation. It is used to avoid division by zero. More...
 

Function/Subroutine Documentation

◆ ppm_boundary_extrapolation()

subroutine, public ppm_functions::ppm_boundary_extrapolation ( integer, intent(in)  N,
real, dimension(:), intent(in)  h,
real, dimension(:), intent(in)  u,
real, dimension(:,:), intent(inout)  edge_values,
real, dimension(:,:), intent(inout)  ppoly_coef,
real, intent(in), optional  h_neglect 
)

Reconstruction by parabolas within boundary cells.

Parameters
[in]nNumber of cells
[in]hcell widths (size N) [H]
[in]ucell averages (size N) [A]
[in,out]edge_valuesedge values of piecewise polynomials [A]
[in,out]ppoly_coefcoefficients of piecewise polynomials, mainly [A]
[in]h_neglectA negligibly small width for the purpose of cell reconstructions [H]

Definition at line 134 of file PPM_functions.F90.

134 !------------------------------------------------------------------------------
135 ! Reconstruction by parabolas within boundary cells.
136 !
137 ! The following explanations apply to the left boundary cell. The same
138 ! reasoning holds for the right boundary cell.
139 !
140 ! A parabola needs to be built in the cell and requires three degrees of
141 ! freedom, which are the right edge value and slope and the cell average.
142 ! The right edge values and slopes are taken to be that of the neighboring
143 ! cell (i.e., the left edge value and slope of the neighboring cell).
144 ! The resulting parabola is not necessarily monotonic and the traditional
145 ! PPM limiter is used to modify one of the edge values in order to yield
146 ! a monotonic parabola.
147 !
148 ! N: number of cells in grid
149 ! h: thicknesses of grid cells
150 ! u: cell averages to use in constructing piecewise polynomials
151 ! edge_values : edge values of piecewise polynomials
152 ! ppoly_coef : coefficients of piecewise polynomials
153 !
154 ! It is assumed that the size of the array 'u' is equal to the number of cells
155 ! defining 'grid' and 'ppoly'. No consistency check is performed here.
156 !------------------------------------------------------------------------------
157 
158  ! Arguments
159  integer, intent(in) :: n !< Number of cells
160  real, dimension(:), intent(in) :: h !< cell widths (size N) [H]
161  real, dimension(:), intent(in) :: u !< cell averages (size N) [A]
162  real, dimension(:,:), intent(inout) :: edge_values !< edge values of piecewise polynomials [A]
163  real, dimension(:,:), intent(inout) :: ppoly_coef !< coefficients of piecewise polynomials, mainly [A]
164  real, optional, intent(in) :: h_neglect !< A negligibly small width for
165  !! the purpose of cell reconstructions [H]
166 
167  ! Local variables
168  integer :: i0, i1
169  real :: u0, u1
170  real :: h0, h1
171  real :: a, b, c
172  real :: u0_l, u0_r
173  real :: u1_l, u1_r
174  real :: slope
175  real :: exp1, exp2
176  real :: hneglect
177 
178  hneglect = hneglect_dflt ; if (present(h_neglect)) hneglect = h_neglect
179 
180  ! ----- Left boundary -----
181  i0 = 1
182  i1 = 2
183  h0 = h(i0)
184  h1 = h(i1)
185  u0 = u(i0)
186  u1 = u(i1)
187 
188  ! Compute the left edge slope in neighboring cell and express it in
189  ! the global coordinate system
190  b = ppoly_coef(i1,2)
191  u1_r = b *((h0+hneglect)/(h1+hneglect)) ! derivative evaluated at xi = 0.0,
192  ! expressed w.r.t. xi (local coord. system)
193 
194  ! Limit the right slope by the PLM limited slope
195  slope = 2.0 * ( u1 - u0 )
196  if ( abs(u1_r) > abs(slope) ) then
197  u1_r = slope
198  endif
199 
200  ! The right edge value in the boundary cell is taken to be the left
201  ! edge value in the neighboring cell
202  u0_r = edge_values(i1,1)
203 
204  ! Given the right edge value and slope, we determine the left
205  ! edge value and slope by computing the parabola as determined by
206  ! the right edge value and slope and the boundary cell average
207  u0_l = 3.0 * u0 + 0.5 * u1_r - 2.0 * u0_r
208 
209  ! Apply the traditional PPM limiter
210  exp1 = (u0_r - u0_l) * (u0 - 0.5*(u0_l+u0_r))
211  exp2 = (u0_r - u0_l) * (u0_r - u0_l) / 6.0
212 
213  if ( exp1 > exp2 ) then
214  u0_l = 3.0 * u0 - 2.0 * u0_r
215  endif
216 
217  if ( exp1 < -exp2 ) then
218  u0_r = 3.0 * u0 - 2.0 * u0_l
219  endif
220 
221  edge_values(i0,1) = u0_l
222  edge_values(i0,2) = u0_r
223 
224  a = u0_l
225  b = 6.0 * u0 - 4.0 * u0_l - 2.0 * u0_r
226  c = 3.0 * ( u0_r + u0_l - 2.0 * u0 )
227 
228  ppoly_coef(i0,1) = a
229  ppoly_coef(i0,2) = b
230  ppoly_coef(i0,3) = c
231 
232  ! ----- Right boundary -----
233  i0 = n-1
234  i1 = n
235  h0 = h(i0)
236  h1 = h(i1)
237  u0 = u(i0)
238  u1 = u(i1)
239 
240  ! Compute the right edge slope in neighboring cell and express it in
241  ! the global coordinate system
242  b = ppoly_coef(i0,2)
243  c = ppoly_coef(i0,3)
244  u1_l = (b + 2*c) ! derivative evaluated at xi = 1.0
245  u1_l = u1_l * ((h1+hneglect)/(h0+hneglect))
246 
247  ! Limit the left slope by the PLM limited slope
248  slope = 2.0 * ( u1 - u0 )
249  if ( abs(u1_l) > abs(slope) ) then
250  u1_l = slope
251  endif
252 
253  ! The left edge value in the boundary cell is taken to be the right
254  ! edge value in the neighboring cell
255  u0_l = edge_values(i0,2)
256 
257  ! Given the left edge value and slope, we determine the right
258  ! edge value and slope by computing the parabola as determined by
259  ! the left edge value and slope and the boundary cell average
260  u0_r = 3.0 * u1 - 0.5 * u1_l - 2.0 * u0_l
261 
262  ! Apply the traditional PPM limiter
263  exp1 = (u0_r - u0_l) * (u1 - 0.5*(u0_l+u0_r))
264  exp2 = (u0_r - u0_l) * (u0_r - u0_l) / 6.0
265 
266  if ( exp1 > exp2 ) then
267  u0_l = 3.0 * u1 - 2.0 * u0_r
268  endif
269 
270  if ( exp1 < -exp2 ) then
271  u0_r = 3.0 * u1 - 2.0 * u0_l
272  endif
273 
274  edge_values(i1,1) = u0_l
275  edge_values(i1,2) = u0_r
276 
277  a = u0_l
278  b = 6.0 * u1 - 4.0 * u0_l - 2.0 * u0_r
279  c = 3.0 * ( u0_r + u0_l - 2.0 * u1 )
280 
281  ppoly_coef(i1,1) = a
282  ppoly_coef(i1,2) = b
283  ppoly_coef(i1,3) = c
284 

◆ ppm_limiter_standard()

subroutine ppm_functions::ppm_limiter_standard ( integer, intent(in)  N,
real, dimension(:), intent(in)  h,
real, dimension(:), intent(in)  u,
real, dimension(:,:), intent(inout)  edge_values,
real, intent(in), optional  h_neglect,
logical, intent(in), optional  answers_2018 
)
private

Adjusts edge values using the standard PPM limiter (Colella & Woodward, JCP 1984) after first checking that the edge values are bounded by neighbors cell averages and that the edge values are monotonic between cell averages.

Parameters
[in]nNumber of cells
[in]hcell widths (size N) [H]
[in]ucell average properties (size N) [A]
[in,out]edge_valuesPotentially modified edge values [A]
[in]h_neglectA negligibly small width [H]
[in]answers_2018If true use older, less acccurate expressions.

Definition at line 63 of file PPM_functions.F90.

63  integer, intent(in) :: n !< Number of cells
64  real, dimension(:), intent(in) :: h !< cell widths (size N) [H]
65  real, dimension(:), intent(in) :: u !< cell average properties (size N) [A]
66  real, dimension(:,:), intent(inout) :: edge_values !< Potentially modified edge values [A]
67  real, optional, intent(in) :: h_neglect !< A negligibly small width [H]
68  logical, optional, intent(in) :: answers_2018 !< If true use older, less acccurate expressions.
69 
70  ! Local variables
71  integer :: k ! Loop index
72  real :: u_l, u_c, u_r ! Cell averages (left, center and right)
73  real :: edge_l, edge_r ! Edge values (left and right)
74  real :: expr1, expr2
75 
76  ! Bound edge values
77  call bound_edge_values( n, h, u, edge_values, h_neglect, answers_2018 )
78 
79  ! Make discontinuous edge values monotonic
80  call check_discontinuous_edge_values( n, u, edge_values )
81 
82  ! Loop on interior cells to apply the standard
83  ! PPM limiter (Colella & Woodward, JCP 84)
84  do k = 2,n-1
85 
86  ! Get cell averages
87  u_l = u(k-1)
88  u_c = u(k)
89  u_r = u(k+1)
90 
91  edge_l = edge_values(k,1)
92  edge_r = edge_values(k,2)
93 
94  if ( (u_r - u_c)*(u_c - u_l) <= 0.0) then
95  ! Flatten extremum
96  edge_l = u_c
97  edge_r = u_c
98  else
99  expr1 = 3.0 * (edge_r - edge_l) * ( (u_c - edge_l) + (u_c - edge_r))
100  expr2 = (edge_r - edge_l) * (edge_r - edge_l)
101  if ( expr1 > expr2 ) then
102  ! Place extremum at right edge of cell by adjusting left edge value
103  edge_l = u_c + 2.0 * ( u_c - edge_r )
104  edge_l = max( min( edge_l, max(u_l, u_c) ), min(u_l, u_c) ) ! In case of round off
105  elseif ( expr1 < -expr2 ) then
106  ! Place extremum at left edge of cell by adjusting right edge value
107  edge_r = u_c + 2.0 * ( u_c - edge_l )
108  edge_r = max( min( edge_r, max(u_r, u_c) ), min(u_r, u_c) ) ! In case of round off
109  endif
110  endif
111  ! This checks that the difference in edge values is representable
112  ! and avoids overshoot problems due to round off.
113  !### The 1.e-60 needs to have units of [A], so this dimensionally inconsisent.
114  if ( abs( edge_r - edge_l )<max(1.e-60,epsilon(u_c)*abs(u_c)) ) then
115  edge_l = u_c
116  edge_r = u_c
117  endif
118 
119  edge_values(k,1) = edge_l
120  edge_values(k,2) = edge_r
121 
122  enddo ! end loop on interior cells
123 
124  ! PCM within boundary cells
125  edge_values(1,:) = u(1)
126  edge_values(n,:) = u(n)
127 

◆ ppm_reconstruction()

subroutine, public ppm_functions::ppm_reconstruction ( integer, intent(in)  N,
real, dimension(n), intent(in)  h,
real, dimension(n), intent(in)  u,
real, dimension(n,2), intent(inout)  edge_values,
real, dimension(n,3), intent(inout)  ppoly_coef,
real, intent(in), optional  h_neglect,
logical, intent(in), optional  answers_2018 
)

Builds quadratic polynomials coefficients from cell mean and edge values.

Parameters
[in]nNumber of cells
[in]hCell widths [H]
[in]uCell averages [A]
[in,out]edge_valuesEdge values [A]
[in,out]ppoly_coefPolynomial coefficients, mainly [A]
[in]h_neglectA negligibly small width [H]
[in]answers_2018If true use older, less acccurate expressions.

Definition at line 29 of file PPM_functions.F90.

29  integer, intent(in) :: n !< Number of cells
30  real, dimension(N), intent(in) :: h !< Cell widths [H]
31  real, dimension(N), intent(in) :: u !< Cell averages [A]
32  real, dimension(N,2), intent(inout) :: edge_values !< Edge values [A]
33  real, dimension(N,3), intent(inout) :: ppoly_coef !< Polynomial coefficients, mainly [A]
34  real, optional, intent(in) :: h_neglect !< A negligibly small width [H]
35  logical, optional, intent(in) :: answers_2018 !< If true use older, less acccurate expressions.
36 
37  ! Local variables
38  integer :: k ! Loop index
39  real :: edge_l, edge_r ! Edge values (left and right)
40 
41  ! PPM limiter
42  call ppm_limiter_standard( n, h, u, edge_values, h_neglect, answers_2018 )
43 
44  ! Loop over all cells
45  do k = 1,n
46 
47  edge_l = edge_values(k,1)
48  edge_r = edge_values(k,2)
49 
50  ! Store polynomial coefficients
51  ppoly_coef(k,1) = edge_l
52  ppoly_coef(k,2) = 4.0 * ( u(k) - edge_l ) + 2.0 * ( u(k) - edge_r )
53  ppoly_coef(k,3) = 3.0 * ( ( edge_r - u(k) ) + ( edge_l - u(k) ) )
54 
55  enddo
56 

Variable Documentation

◆ hneglect_dflt

real, parameter ppm_functions::hneglect_dflt = 1.E-30
private

A tiny width that is so small that adding it to cell widths does not change the value due to a computational representation. It is used to avoid division by zero.

Note
This is a dimensional parameter and should really include a unit conversion.

Definition at line 23 of file PPM_functions.F90.

23 real, parameter :: hneglect_dflt = 1.e-30