regrid_edge_values module reference

Edge value estimation for high-order resconstruction.

More…

Functions/Subroutines

bound_edge_values()

Bound edge values by neighboring cell averages.

average_discontinuous_edge_values()

Replace discontinuous collocated edge values with their average.

check_discontinuous_edge_values()

Check discontinuous edge values and replace them with their average if not monotonic.

edge_values_explicit_h2()

Compute h2 edge values (explicit second order accurate) in the same units as u.

edge_values_explicit_h4()

Compute h4 edge values (explicit fourth order accurate) in the same units as u.

edge_values_implicit_h4()

Compute ih4 edge values (implicit fourth order accurate) in the same units as u.

end_value_h4()

Determine a one-sided 4th order polynomial fit of u to the data points for the purposes of specifying edge values, as described in the appendix of White and Adcroft JCP 2008.

edge_slopes_implicit_h3()

Compute ih3 edge slopes (implicit third order accurate) in the same units as h.

edge_slopes_implicit_h5()

Compute ih5 edge slopes (implicit fifth order accurate)

edge_values_implicit_h6()

Compute ih6 edge values (implicit sixth order accurate) in the same units as u.

solve_diag_dominant_tridiag()

Solve the tridiagonal system AX = R.

linear_solver()

Solve the linear system AX = R by Gaussian elimination.

test_line()

Test that A*C = R to within a tolerance, issuing a fatal error with an explanatory message if they do not.

Detailed Description

Edge value estimation for high-order resconstruction.

Function/Subroutine Documentation

subroutine regrid_edge_values/bound_edge_values(N, h, u, edge_val, h_neglect, answers_2018)

Bound edge values by neighboring cell averages.

In this routine, we loop on all cells to bound their left and right edge values by the cell averages. That is, the left edge value must lie between the left cell average and the central cell average. A similar reasoning applies to the right edge values.

Both boundary edge values are set equal to the boundary cell averages. Any extrapolation scheme is applied after this routine has been called. Therefore, boundary cells are treated as if they were local extrama.

Parameters

  • n :: [in] Number of cells

  • h :: [in] cell widths [H]

  • u :: [in] cell average properties in arbitrary units [A]

  • edge_val :: [inout] Potentially modified edge values [A]; the second index is for the two edges of each cell.

  • h_neglect :: [in] A negligibly small width [H]

  • answers_2018 :: [in] If true use older, less acccurate expressions.

Call to

hneglect_dflt

Called from

p1m_functions::p1m_interpolation p3m_functions::p3m_limiter ppm_functions::ppm_limiter_standard pqm_functions::pqm_limiter

subroutine regrid_edge_values/average_discontinuous_edge_values(N, edge_val)

Replace discontinuous collocated edge values with their average.

For each interior edge, check whether the edge values are discontinuous. If so, compute the average and replace the edge values by the average.

Parameters

  • n :: [in] Number of cells

  • edge_val :: [inout] Edge values that may be modified [A]; the second index is for the two edges of each cell.

Called from

p1m_functions::p1m_interpolation p3m_functions::p3m_limiter

subroutine regrid_edge_values/check_discontinuous_edge_values(N, u, edge_val)

Check discontinuous edge values and replace them with their average if not monotonic.

For each interior edge, check whether the edge values are discontinuous. If so and if they are not monotonic, replace each edge value by their average.

Parameters

  • n :: [in] Number of cells

  • u :: [in] cell averages in arbitrary units [A]

  • edge_val :: [inout] Cell edge values [A]; the second index is for the two edges of each cell.

Called from

ppm_functions::ppm_limiter_standard pqm_functions::pqm_limiter

subroutine regrid_edge_values/edge_values_explicit_h2(N, h, u, edge_val)

Compute h2 edge values (explicit second order accurate) in the same units as u.

Parameters

  • n :: [in] Number of cells

  • h :: [in] cell widths [H]

  • u :: [in] cell average properties in arbitrary units [A]

  • edge_val :: [inout] Returned edge values [A]; the second index is for the two edges of each cell.

subroutine regrid_edge_values/edge_values_explicit_h4(N, h, u, edge_val, h_neglect, answers_2018)

Compute h4 edge values (explicit fourth order accurate) in the same units as u.

Compute edge values based on fourth-order explicit estimates. These estimates are based on a cubic interpolant spanning four cells and evaluated at the location of the middle edge. An interpolant spanning cells i-2, i-1, i and i+1 is evaluated at edge i-1/2. The estimate for each edge is unique.

i-2    i-1     i     i+1

..ooooo.. i-1/2

The first two edge values are estimated by evaluating the first available cubic interpolant, i.e., the interpolant spanning cells 1, 2, 3 and 4. Similarly, the last two edge values are estimated by evaluating the last available interpolant.

For this fourth-order scheme, at least four cells must exist.

Parameters

  • n :: [in] Number of cells

  • h :: [in] cell widths [H]

  • u :: [in] cell average properties in arbitrary units [A]

  • edge_val :: [inout] Returned edge values [A]; the second index is for the two edges of each cell.

  • h_neglect :: [in] A negligibly small width [H]

  • answers_2018 :: [in] If true use older, less acccurate expressions.

Call to

end_value_h4 polynomial_functions::evaluation_polynomial hminfrac hneglect_dflt hneglect_edge_dflt regrid_solvers::solve_linear_system

subroutine regrid_edge_values/edge_values_implicit_h4(N, h, u, edge_val, h_neglect, answers_2018)

Compute ih4 edge values (implicit fourth order accurate) in the same units as u.

Compute edge values based on fourth-order implicit estimates.

Fourth-order implicit estimates of edge values are based on a two-cell stencil. A tridiagonal system is set up and is based on expressing the edge values in terms of neighboring cell averages. The generic relationship is

\[\alpha u_{i-1/2} + u_{i+1/2} + \beta u_{i+3/2} = a \bar{u}_i + b \bar{u}_{i+1}\]

and the stencil looks like this

i     i+1

..ooo.. i-1/2 i+1/2 i+3/2

In this routine, the coefficients \(\alpha\) , \(\beta\) , \(a\) and \(b\) are computed, the tridiagonal system is built, boundary conditions are prescribed and the system is solved to yield edge-value estimates.

There are N+1 unknowns and we are able to write N-1 equations. The boundary conditions close the system.

Parameters

  • n :: [in] Number of cells

  • h :: [in] cell widths [H]

  • u :: [in] cell average properties in arbitrary units [A]

  • edge_val :: [inout] Returned edge values [A]; the second index is for the two edges of each cell.

  • h_neglect :: [in] A negligibly small width [H]

  • answers_2018 :: [in] If true use older, less acccurate expressions.

Call to

end_value_h4 polynomial_functions::evaluation_polynomial hminfrac hneglect_dflt hneglect_edge_dflt solve_diag_dominant_tridiag regrid_solvers::solve_linear_system regrid_solvers::solve_tridiagonal_system

Called from

mom_ale::ts_ppm_edge_values

subroutine regrid_edge_values/end_value_h4(dz, u, Csys)

Determine a one-sided 4th order polynomial fit of u to the data points for the purposes of specifying edge values, as described in the appendix of White and Adcroft JCP 2008.

Parameters

  • dz :: [in] The thicknesses of 4 layers, starting at the edge [H]. The values of dz must be positive.

  • u :: [in] The average properties of 4 layers, starting at the edge [A]

  • csys :: [out] The four coefficients of a 4th order polynomial fit of u as a function of z [A H-(n-1)]

Called from

edge_slopes_implicit_h3 edge_values_explicit_h4 edge_values_implicit_h4

subroutine regrid_edge_values/edge_slopes_implicit_h3(N, h, u, edge_slopes, h_neglect, answers_2018)

Compute ih3 edge slopes (implicit third order accurate) in the same units as h.

Compute edge slopes based on third-order implicit estimates. Note that the estimates are fourth-order accurate on uniform grids

Third-order implicit estimates of edge slopes are based on a two-cell stencil. A tridiagonal system is set up and is based on expressing the edge slopes in terms of neighboring cell averages. The generic relationship is

\[\alpha u'_{i-1/2} + u'_{i+1/2} + \beta u'_{i+3/2} = a \bar{u}_i + b \bar{u}_{i+1}\]

and the stencil looks like this

i     i+1

..ooo.. i-1/2 i+1/2 i+3/2

In this routine, the coefficients \(\alpha\) , \(\beta\) , a and b are computed, the tridiagonal system is built, boundary conditions are prescribed and the system is solved to yield edge-slope estimates.

There are N+1 unknowns and we are able to write N-1 equations. The boundary conditions close the system.

Parameters

  • n :: [in] Number of cells

  • h :: [in] cell widths [H]

  • u :: [in] cell average properties in arbitrary units [A]

  • edge_slopes :: [inout] Returned edge slopes [A H-1]; the second index is for the two edges of each cell.

  • h_neglect :: [in] A negligibly small width [H]

  • answers_2018 :: [in] If true use older, less acccurate expressions.

Call to

end_value_h4 polynomial_functions::evaluation_polynomial hneglect_dflt solve_diag_dominant_tridiag regrid_solvers::solve_linear_system regrid_solvers::solve_tridiagonal_system

Called from

mom_remapping::build_reconstructions_1d regrid_interp::regridding_set_ppolys

subroutine regrid_edge_values/edge_slopes_implicit_h5(N, h, u, edge_slopes, h_neglect, answers_2018)

Compute ih5 edge slopes (implicit fifth order accurate)

Parameters

  • n :: [in] Number of cells

  • h :: [in] cell widths [H]

  • u :: [in] cell average properties in arbitrary units [A]

  • edge_slopes :: [inout] Returned edge slopes [A H-1]; the second index is for the two edges of each cell.

  • h_neglect :: [in] A negligibly small width [H]

  • answers_2018 :: [in] If true use older, less acccurate expressions.

Call to

hneglect_dflt linear_solver regrid_solvers::solve_tridiagonal_system

Called from

mom_remapping::build_reconstructions_1d regrid_interp::regridding_set_ppolys

subroutine regrid_edge_values/edge_values_implicit_h6(N, h, u, edge_val, h_neglect, answers_2018)

Compute ih6 edge values (implicit sixth order accurate) in the same units as u.

Sixth-order implicit estimates of edge values are based on a four-cell, three-edge stencil. A tridiagonal system is set up and is based on expressing the edge values in terms of neighboring cell averages.

The generic relationship is

\[\alpha u_{i-1/2} + u_{i+1/2} + \beta u_{i+3/2} = a \bar{u}_{i-1} + b \bar{u}_i + c \bar{u}_{i+1} + d \bar{u}_{i+2}\]

and the stencil looks like this

i-1     i     i+1    i+2

..ooooo.. i-1/2 i+1/2 i+3/2

In this routine, the coefficients \(\alpha\) , \(\beta\) , a, b, c and d are computed, the tridiagonal system is built, boundary conditions are prescribed and the system is solved to yield edge-value estimates. This scheme is described in detail by White and Adcroft, 2009, J. Comp. Phys, https://doi.org/10.1016/j.jcp.2008.04.026

Note that the centered stencil only applies to edges 3 to N-1 (edges are numbered 1 to n+1), which yields N-3 equations for N+1 unknowns. Two other equations are written by using a right-biased stencil for edge 2 and a left-biased stencil for edge N. The prescription of boundary conditions (using sixth-order polynomials) closes the system.

CAUTION: For each edge, in order to determine the coefficients of the implicit expression, a 6x6 linear system is solved. This may become computationally expensive if regridding is carried out often. Figuring out closed-form expressions for these coefficients on nonuniform meshes turned out to be intractable.

Parameters

  • n :: [in] Number of cells

  • h :: [in] cell widths [H]

  • u :: [in] cell average properties (size N) in arbitrary units [A]

  • edge_val :: [inout] Returned edge values [A]; the second index is for the two edges of each cell.

  • h_neglect :: [in] A negligibly small width [H]

  • answers_2018 :: [in] If true use older, less acccurate expressions.

Call to

polynomial_functions::evaluation_polynomial hminfrac hneglect_edge_dflt linear_solver regrid_solvers::solve_tridiagonal_system

subroutine regrid_edge_values/solve_diag_dominant_tridiag(Al, Ac, Au, R, X, N)

Solve the tridiagonal system AX = R.

This routine uses a variant of Thomas’s algorithm to solve the tridiagonal system AX = R, in a form that is guaranteed to avoid dividing by a zero pivot. The matrix A is made up of lower (Al) and upper diagonals (Au) and a central diagonal Ad = Ac+Al+Au, where Al, Au, and Ac are all positive (or negative) definite. However when Ac is smaller than roundoff compared with (Al+Au), the answers are prone to inaccuracy.

Parameters

  • n :: [in] The size of the system

  • ac :: [in] Matrix center diagonal offset from Al + Au

  • al :: [in] Matrix lower diagonal

  • au :: [in] Matrix upper diagonal

  • r :: [in] system right-hand side

  • x :: [out] solution vector

Called from

edge_slopes_implicit_h3 edge_values_implicit_h4

subroutine regrid_edge_values/linear_solver(N, A, R, X)

Solve the linear system AX = R by Gaussian elimination.

This routine uses Gauss’s algorithm to transform the system’s original matrix into an upper triangular matrix. Back substitution then yields the answer. The matrix A must be square, with the first index varing along the row.

Parameters

  • n :: [in] The size of the system

  • a :: [inout] The matrix being inverted [nondim]

  • r :: [inout] system right-hand side [A]

  • x :: [inout] solution vector [A]

Call to

mom_error_handler::mom_error

Called from

edge_slopes_implicit_h5 edge_values_implicit_h6

subroutine regrid_edge_values/test_line(msg, N, A, C, R, mag, tol)

Test that A*C = R to within a tolerance, issuing a fatal error with an explanatory message if they do not.

Parameters

  • mag :: [in] The magnitude of leading order terms in this line

  • n :: [in] The number of points in the system

  • a :: [in] One of the two vectors being multiplied

  • c :: [in] One of the two vectors being multiplied

  • r :: [in] The expected solution of the equation

  • msg :: [in] An identifying message for this test

  • tol :: [in] The fractional tolerance for the two solutions

Call to

mom_error_handler::mom_error