\hypertarget{namespacep3m__functions}{}\doxysection{p3m\+\_\+functions Module Reference} \label{namespacep3m__functions}\index{p3m\_functions@{p3m\_functions}} \doxysubsection{Detailed Description} Cubic interpolation functions. Date of creation\+: 2008.\+06.\+09 L. White. This module contains p3m interpolation routines. p3m interpolation is performed by estimating the edge values and slopes and constructing a cubic polynomial. We then make sure that the edge values are bounded and continuous and we then modify the slopes to get a monotonic cubic curve. \doxysubsection*{Functions/\+Subroutines} \begin{DoxyCompactItemize} \item subroutine, public \mbox{\hyperlink{namespacep3m__functions_a685d3ef292536dae810b2059ccaa5819}{p3m\+\_\+interpolation}} (N, h, u, edge\+\_\+values, ppoly\+\_\+S, ppoly\+\_\+coef, h\+\_\+neglect, answers\+\_\+2018) \begin{DoxyCompactList}\small\item\em Set up a piecewise cubic interpolation from cell averages and estimated edge slopes and values. \end{DoxyCompactList}\item subroutine \mbox{\hyperlink{namespacep3m__functions_a44db59cb5df7f139e05b745746342fcf}{p3m\+\_\+limiter}} (N, h, u, edge\+\_\+values, ppoly\+\_\+S, ppoly\+\_\+coef, h\+\_\+neglect, answers\+\_\+2018) \begin{DoxyCompactList}\small\item\em Adust a piecewise cubic reconstruction with a limiter that adjusts the edge values and slopes. \end{DoxyCompactList}\item subroutine, public \mbox{\hyperlink{namespacep3m__functions_ab70eb9e69fc6586e1d6a371d2eeb44d1}{p3m\+\_\+boundary\+\_\+extrapolation}} (N, h, u, edge\+\_\+values, ppoly\+\_\+S, ppoly\+\_\+coef, h\+\_\+neglect, h\+\_\+neglect\+\_\+edge) \begin{DoxyCompactList}\small\item\em Calculate the edge values and slopes at boundary cells as part of building a piecewise cubic sub-\/grid scale profiles. \end{DoxyCompactList}\item subroutine \mbox{\hyperlink{namespacep3m__functions_a82c244f63ffd696c15d58a029e8e24f7}{build\+\_\+cubic\+\_\+interpolant}} (h, k, edge\+\_\+values, ppoly\+\_\+S, ppoly\+\_\+coef) \begin{DoxyCompactList}\small\item\em Build cubic interpolant in cell k. \end{DoxyCompactList}\item logical function \mbox{\hyperlink{namespacep3m__functions_a3c461688b7c3ae5b2b4a7cf3311b8a69}{is\+\_\+cubic\+\_\+monotonic}} (ppoly\+\_\+coef, k) \begin{DoxyCompactList}\small\item\em Check whether the cubic reconstruction in cell k is monotonic. \end{DoxyCompactList}\item subroutine \mbox{\hyperlink{namespacep3m__functions_adb96651fe725f11e90dec2b8509989b0}{monotonize\+\_\+cubic}} (h, u0\+\_\+l, u0\+\_\+r, sigma\+\_\+l, sigma\+\_\+r, slope, u1\+\_\+l, u1\+\_\+r) \begin{DoxyCompactList}\small\item\em Monotonize a cubic curve by modifying the edge slopes. \end{DoxyCompactList}\end{DoxyCompactItemize} \doxysubsection*{Variables} \begin{DoxyCompactItemize} \item \mbox{\Hypertarget{namespacep3m__functions_aedfe0a0f65453ca61b90d4c601747c90}\label{namespacep3m__functions_aedfe0a0f65453ca61b90d4c601747c90}} real, parameter \mbox{\hyperlink{namespacep3m__functions_aedfe0a0f65453ca61b90d4c601747c90}{hneglect\+\_\+dflt}} = 1.E-\/30 \begin{DoxyCompactList}\small\item\em Default value of a negligible cell thickness. \end{DoxyCompactList}\item \mbox{\Hypertarget{namespacep3m__functions_abdb3a0dd6c75e81bbcb5f5c982c31e7b}\label{namespacep3m__functions_abdb3a0dd6c75e81bbcb5f5c982c31e7b}} real, parameter \mbox{\hyperlink{namespacep3m__functions_abdb3a0dd6c75e81bbcb5f5c982c31e7b}{hneglect\+\_\+edge\+\_\+dflt}} = 1.E-\/10 \begin{DoxyCompactList}\small\item\em Default value of a negligible edge thickness. \end{DoxyCompactList}\end{DoxyCompactItemize} \doxysubsection{Function/\+Subroutine Documentation} \mbox{\Hypertarget{namespacep3m__functions_a82c244f63ffd696c15d58a029e8e24f7}\label{namespacep3m__functions_a82c244f63ffd696c15d58a029e8e24f7}} \index{p3m\_functions@{p3m\_functions}!build\_cubic\_interpolant@{build\_cubic\_interpolant}} \index{build\_cubic\_interpolant@{build\_cubic\_interpolant}!p3m\_functions@{p3m\_functions}} \doxysubsubsection{\texorpdfstring{build\_cubic\_interpolant()}{build\_cubic\_interpolant()}} {\footnotesize\ttfamily subroutine p3m\+\_\+functions\+::build\+\_\+cubic\+\_\+interpolant (\begin{DoxyParamCaption}\item[{real, dimension(\+:), intent(in)}]{h, }\item[{integer, intent(in)}]{k, }\item[{real, dimension(\+:,\+:), intent(in)}]{edge\+\_\+values, }\item[{real, dimension(\+:,\+:), intent(in)}]{ppoly\+\_\+S, }\item[{real, dimension(\+:,\+:), intent(inout)}]{ppoly\+\_\+coef }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} Build cubic interpolant in cell k. Given edge values and edge slopes, compute coefficients of cubic in cell k. N\+O\+TE\+: edge values and slopes M\+U\+ST have been properly calculated prior to calling this routine. \begin{DoxyParams}[1]{Parameters} \mbox{\texttt{ in}} & {\em h} & cell widths (size N) \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em k} & The index of the cell to work on \\ \hline \mbox{\texttt{ in}} & {\em edge\+\_\+values} & Edge value of polynomial in arbitrary units \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em ppoly\+\_\+s} & Edge slope of polynomial \mbox{[}A H-\/1\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em ppoly\+\_\+coef} & Coefficients of polynomial \mbox{[}A\mbox{]} \\ \hline \end{DoxyParams} Definition at line 343 of file P3\+M\+\_\+functions.\+F90. \begin{DoxyCode}{0} \DoxyCodeLine{344 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< cell widths (size N) [H]}} \DoxyCodeLine{345 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: k\textcolor{comment}{ !< The index of the cell to work on}} \DoxyCodeLine{346 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(in)} :: edge\_values\textcolor{comment}{ !< Edge value of polynomial in arbitrary units [A]}} \DoxyCodeLine{347 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(in)} :: ppoly\_S\textcolor{comment}{ !< Edge slope of polynomial [A H-\/1]}} \DoxyCodeLine{348 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_coef\textcolor{comment}{ !< Coefficients of polynomial [A]}} \DoxyCodeLine{349 } \DoxyCodeLine{350 \textcolor{comment}{! Local variables}} \DoxyCodeLine{351 \textcolor{keywordtype}{ real} :: u0\_l, u0\_r \textcolor{comment}{! edge values [A]}} \DoxyCodeLine{352 \textcolor{keywordtype}{ real} :: u1\_l, u1\_r \textcolor{comment}{! edge slopes times the cell width [A]}} \DoxyCodeLine{353 \textcolor{keywordtype}{ real} :: h\_c \textcolor{comment}{! cell width [H]}} \DoxyCodeLine{354 \textcolor{keywordtype}{ real} :: a0, a1, a2, a3 \textcolor{comment}{! cubic coefficients [A]}} \DoxyCodeLine{355 } \DoxyCodeLine{356 h\_c = h(k)} \DoxyCodeLine{357 } \DoxyCodeLine{358 u0\_l = edge\_values(k,1)} \DoxyCodeLine{359 u0\_r = edge\_values(k,2)} \DoxyCodeLine{360 } \DoxyCodeLine{361 u1\_l = ppoly\_s(k,1) * h\_c} \DoxyCodeLine{362 u1\_r = ppoly\_s(k,2) * h\_c} \DoxyCodeLine{363 } \DoxyCodeLine{364 a0 = u0\_l} \DoxyCodeLine{365 a1 = u1\_l} \DoxyCodeLine{366 a2 = 3.0 * ( u0\_r -\/ u0\_l ) -\/ u1\_r -\/ 2.0 * u1\_l} \DoxyCodeLine{367 a3 = u1\_r + u1\_l + 2.0 * ( u0\_l -\/ u0\_r )} \DoxyCodeLine{368 } \DoxyCodeLine{369 ppoly\_coef(k,1) = a0} \DoxyCodeLine{370 ppoly\_coef(k,2) = a1} \DoxyCodeLine{371 ppoly\_coef(k,3) = a2} \DoxyCodeLine{372 ppoly\_coef(k,4) = a3} \DoxyCodeLine{373 } \end{DoxyCode} \mbox{\Hypertarget{namespacep3m__functions_a3c461688b7c3ae5b2b4a7cf3311b8a69}\label{namespacep3m__functions_a3c461688b7c3ae5b2b4a7cf3311b8a69}} \index{p3m\_functions@{p3m\_functions}!is\_cubic\_monotonic@{is\_cubic\_monotonic}} \index{is\_cubic\_monotonic@{is\_cubic\_monotonic}!p3m\_functions@{p3m\_functions}} \doxysubsubsection{\texorpdfstring{is\_cubic\_monotonic()}{is\_cubic\_monotonic()}} {\footnotesize\ttfamily logical function p3m\+\_\+functions\+::is\+\_\+cubic\+\_\+monotonic (\begin{DoxyParamCaption}\item[{real, dimension(\+:,\+:), intent(in)}]{ppoly\+\_\+coef, }\item[{integer, intent(in)}]{k }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} Check whether the cubic reconstruction in cell k is monotonic. This function checks whether the cubic curve in cell k is monotonic. If so, returns 1. Otherwise, returns 0. The cubic is monotonic if the first derivative is single-\/signed in (0,1). Hence, we check whether the roots (if any) lie inside this interval. If there is no root or if both roots lie outside this interval, the cubic is monotonic. \begin{DoxyParams}[1]{Parameters} \mbox{\texttt{ in}} & {\em ppoly\+\_\+coef} & Coefficients of cubic polynomial in arbitary units \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em k} & The index of the cell to work on \\ \hline \end{DoxyParams} Definition at line 385 of file P3\+M\+\_\+functions.\+F90. \begin{DoxyCode}{0} \DoxyCodeLine{386 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(in)} :: ppoly\_coef\textcolor{comment}{ !< Coefficients of cubic polynomial in arbitary units [A]}} \DoxyCodeLine{387 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: k\textcolor{comment}{ !< The index of the cell to work on}} \DoxyCodeLine{388 \textcolor{comment}{! Local variables}} \DoxyCodeLine{389 \textcolor{keywordtype}{ real} :: a, b, c \textcolor{comment}{! Coefficients of the first derivative of the cubic [A]}} \DoxyCodeLine{390 } \DoxyCodeLine{391 a = ppoly\_coef(k,2)} \DoxyCodeLine{392 b = 2.0 * ppoly\_coef(k,3)} \DoxyCodeLine{393 c = 3.0 * ppoly\_coef(k,4)} \DoxyCodeLine{394 } \DoxyCodeLine{395 \textcolor{comment}{! Look for real roots of the quadratic derivative equation, c*x**2 + b*x + a = 0, in (0, 1)}} \DoxyCodeLine{396 \textcolor{keywordflow}{if} (b*b -\/ 4.0*a*c <= 0.0) \textcolor{keywordflow}{then} \textcolor{comment}{! The cubic is monotonic everywhere.}} \DoxyCodeLine{397 is\_cubic\_monotonic = .true.} \DoxyCodeLine{398 \textcolor{keywordflow}{elseif} (a * (a + (b + c)) < 0.0) \textcolor{keywordflow}{then} \textcolor{comment}{! The derivative changes sign between the endpoints of (0, 1)}} \DoxyCodeLine{399 is\_cubic\_monotonic = .false.} \DoxyCodeLine{400 \textcolor{keywordflow}{elseif} (b * (b + 2.0*c) < 0.0) \textcolor{keywordflow}{then} \textcolor{comment}{! The second derivative changes sign inside of (0, 1)}} \DoxyCodeLine{401 is\_cubic\_monotonic = .false.} \DoxyCodeLine{402 \textcolor{keywordflow}{else}} \DoxyCodeLine{403 is\_cubic\_monotonic = .true.} \DoxyCodeLine{404 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{405 } \end{DoxyCode} \mbox{\Hypertarget{namespacep3m__functions_adb96651fe725f11e90dec2b8509989b0}\label{namespacep3m__functions_adb96651fe725f11e90dec2b8509989b0}} \index{p3m\_functions@{p3m\_functions}!monotonize\_cubic@{monotonize\_cubic}} \index{monotonize\_cubic@{monotonize\_cubic}!p3m\_functions@{p3m\_functions}} \doxysubsubsection{\texorpdfstring{monotonize\_cubic()}{monotonize\_cubic()}} {\footnotesize\ttfamily subroutine p3m\+\_\+functions\+::monotonize\+\_\+cubic (\begin{DoxyParamCaption}\item[{real, intent(in)}]{h, }\item[{real, intent(in)}]{u0\+\_\+l, }\item[{real, intent(in)}]{u0\+\_\+r, }\item[{real, intent(in)}]{sigma\+\_\+l, }\item[{real, intent(in)}]{sigma\+\_\+r, }\item[{real, intent(in)}]{slope, }\item[{real, intent(inout)}]{u1\+\_\+l, }\item[{real, intent(inout)}]{u1\+\_\+r }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} Monotonize a cubic curve by modifying the edge slopes. This routine takes care of monotonizing a cubic on \mbox{[}0,1\mbox{]} by modifying the edge slopes. The edge values are N\+OT modified. The cubic is entirely determined by the four degrees of freedom u0\+\_\+l, u0\+\_\+r, u1\+\_\+l and u1\+\_\+r. u1\+\_\+l and u1\+\_\+r are the edge slopes expressed in the G\+L\+O\+B\+AL coordinate system. The monotonization occurs as follows. \begin{DoxyParams}[1]{Parameters} \mbox{\texttt{ in}} & {\em h} & cell width \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em u0\+\_\+l} & left edge value in arbitrary units \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em u0\+\_\+r} & right edge value \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em sigma\+\_\+l} & left 2nd-\/order slopes \mbox{[}A H-\/1\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em sigma\+\_\+r} & right 2nd-\/order slopes \mbox{[}A H-\/1\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em slope} & limited P\+LM slope \mbox{[}A H-\/1\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em u1\+\_\+l} & left edge slopes \mbox{[}A H-\/1\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em u1\+\_\+r} & right edge slopes \mbox{[}A H-\/1\mbox{]} \\ \hline \end{DoxyParams} Definition at line 435 of file P3\+M\+\_\+functions.\+F90. \begin{DoxyCode}{0} \DoxyCodeLine{436 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< cell width [H]}} \DoxyCodeLine{437 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(in)} :: u0\_l\textcolor{comment}{ !< left edge value in arbitrary units [A]}} \DoxyCodeLine{438 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(in)} :: u0\_r\textcolor{comment}{ !< right edge value [A]}} \DoxyCodeLine{439 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(in)} :: sigma\_l\textcolor{comment}{ !< left 2nd-\/order slopes [A H-\/1]}} \DoxyCodeLine{440 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(in)} :: sigma\_r\textcolor{comment}{ !< right 2nd-\/order slopes [A H-\/1]}} \DoxyCodeLine{441 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(in)} :: slope\textcolor{comment}{ !< limited PLM slope [A H-\/1]}} \DoxyCodeLine{442 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(inout)} :: u1\_l\textcolor{comment}{ !< left edge slopes [A H-\/1]}} \DoxyCodeLine{443 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{intent(inout)} :: u1\_r\textcolor{comment}{ !< right edge slopes [A H-\/1]}} \DoxyCodeLine{444 \textcolor{comment}{! Local variables}} \DoxyCodeLine{445 \textcolor{keywordtype}{logical} :: found\_ip} \DoxyCodeLine{446 \textcolor{keywordtype}{logical} :: inflexion\_l \textcolor{comment}{! bool telling if inflex. pt must be on left}} \DoxyCodeLine{447 \textcolor{keywordtype}{logical} :: inflexion\_r \textcolor{comment}{! bool telling if inflex. pt must be on right}} \DoxyCodeLine{448 \textcolor{keywordtype}{ real} :: a1, a2, a3 \textcolor{comment}{! Temporary slopes times the cell width [A]}} \DoxyCodeLine{449 \textcolor{keywordtype}{ real} :: u1\_l\_tmp \textcolor{comment}{! trial left edge slope [A H-\/1]}} \DoxyCodeLine{450 \textcolor{keywordtype}{ real} :: u1\_r\_tmp \textcolor{comment}{! trial right edge slope [A H-\/1]}} \DoxyCodeLine{451 \textcolor{keywordtype}{ real} :: xi\_ip \textcolor{comment}{! location of inflexion point in cell coordinates (0,1) [nondim]}} \DoxyCodeLine{452 \textcolor{keywordtype}{ real} :: slope\_ip \textcolor{comment}{! slope at inflexion point times cell width [A]}} \DoxyCodeLine{453 } \DoxyCodeLine{454 found\_ip = .false.} \DoxyCodeLine{455 inflexion\_l = .false.} \DoxyCodeLine{456 inflexion\_r = .false.} \DoxyCodeLine{457 } \DoxyCodeLine{458 \textcolor{comment}{! If the edge slopes are inconsistent w.r.t. the limited PLM slope,}} \DoxyCodeLine{459 \textcolor{comment}{! set them to zero}} \DoxyCodeLine{460 \textcolor{keywordflow}{if} ( u1\_l*slope <= 0.0 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{461 u1\_l = 0.0} \DoxyCodeLine{462 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{463 } \DoxyCodeLine{464 \textcolor{keywordflow}{if} ( u1\_r*slope <= 0.0 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{465 u1\_r = 0.0} \DoxyCodeLine{466 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{467 } \DoxyCodeLine{468 \textcolor{comment}{! Compute the location of the inflexion point, which is the root}} \DoxyCodeLine{469 \textcolor{comment}{! of the second derivative}} \DoxyCodeLine{470 a1 = h * u1\_l} \DoxyCodeLine{471 a2 = 3.0 * ( u0\_r -\/ u0\_l ) -\/ h*(u1\_r + 2.0*u1\_l)} \DoxyCodeLine{472 a3 = h*(u1\_r + u1\_l) + 2.0*(u0\_l -\/ u0\_r)} \DoxyCodeLine{473 } \DoxyCodeLine{474 \textcolor{comment}{! There is a possible root (and inflexion point) only if a3 is nonzero.}} \DoxyCodeLine{475 \textcolor{comment}{! When a3 is zero, the second derivative of the cubic is constant (the}} \DoxyCodeLine{476 \textcolor{comment}{! cubic degenerates into a parabola) and no inflexion point exists.}} \DoxyCodeLine{477 \textcolor{keywordflow}{if} ( a3 /= 0.0 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{478 \textcolor{comment}{! Location of inflexion point}} \DoxyCodeLine{479 xi\_ip = -\/ a2 / (3.0 * a3)} \DoxyCodeLine{480 } \DoxyCodeLine{481 \textcolor{comment}{! If the inflexion point lies in [0,1], change boolean value}} \DoxyCodeLine{482 \textcolor{keywordflow}{if} ( (xi\_ip >= 0.0) .AND. (xi\_ip <= 1.0) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{483 found\_ip = .true.} \DoxyCodeLine{484 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{485 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{486 } \DoxyCodeLine{487 \textcolor{comment}{! When there is an inflexion point within [0,1], check the slope}} \DoxyCodeLine{488 \textcolor{comment}{! to see if it is consistent with the limited PLM slope. If not,}} \DoxyCodeLine{489 \textcolor{comment}{! decide on which side we want to collapse the inflexion point.}} \DoxyCodeLine{490 \textcolor{comment}{! If the inflexion point lies on one of the edges, the cubic is}} \DoxyCodeLine{491 \textcolor{comment}{! guaranteed to be monotonic}} \DoxyCodeLine{492 \textcolor{keywordflow}{if} ( found\_ip ) \textcolor{keywordflow}{then}} \DoxyCodeLine{493 slope\_ip = a1 + 2.0*a2*xi\_ip + 3.0*a3*xi\_ip*xi\_ip} \DoxyCodeLine{494 } \DoxyCodeLine{495 \textcolor{comment}{! Check whether slope is consistent}} \DoxyCodeLine{496 \textcolor{keywordflow}{if} ( slope\_ip*slope < 0.0 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{497 \textcolor{keywordflow}{if} ( abs(sigma\_l) < abs(sigma\_r) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{498 inflexion\_l = .true.} \DoxyCodeLine{499 \textcolor{keywordflow}{else}} \DoxyCodeLine{500 inflexion\_r = .true.} \DoxyCodeLine{501 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{502 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{503 \textcolor{keywordflow}{ endif} \textcolor{comment}{! found\_ip}} \DoxyCodeLine{504 } \DoxyCodeLine{505 \textcolor{comment}{! At this point, if the cubic is not monotonic, we know where the}} \DoxyCodeLine{506 \textcolor{comment}{! inflexion point should lie. When the cubic is monotonic, both}} \DoxyCodeLine{507 \textcolor{comment}{! 'inflexion\_l' and 'inflexion\_r' are false and nothing is to be done.}} \DoxyCodeLine{508 } \DoxyCodeLine{509 \textcolor{comment}{! Move inflexion point on the left}} \DoxyCodeLine{510 \textcolor{keywordflow}{if} ( inflexion\_l ) \textcolor{keywordflow}{then}} \DoxyCodeLine{511 } \DoxyCodeLine{512 u1\_l\_tmp = 1.5*(u0\_r-\/u0\_l)/h -\/ 0.5*u1\_r} \DoxyCodeLine{513 u1\_r\_tmp = 3.0*(u0\_r-\/u0\_l)/h -\/ 2.0*u1\_l} \DoxyCodeLine{514 } \DoxyCodeLine{515 \textcolor{keywordflow}{if} ( (u1\_l\_tmp*slope < 0.0) .AND. (u1\_r\_tmp*slope < 0.0) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{516 } \DoxyCodeLine{517 u1\_l = 0.0} \DoxyCodeLine{518 u1\_r = 3.0 * (u0\_r -\/ u0\_l) / h} \DoxyCodeLine{519 } \DoxyCodeLine{520 \textcolor{keywordflow}{elseif} (u1\_l\_tmp*slope < 0.0) \textcolor{keywordflow}{then}} \DoxyCodeLine{521 } \DoxyCodeLine{522 u1\_r = u1\_r\_tmp} \DoxyCodeLine{523 u1\_l = 1.5*(u0\_r -\/ u0\_l)/h -\/ 0.5*u1\_r} \DoxyCodeLine{524 } \DoxyCodeLine{525 \textcolor{keywordflow}{elseif} (u1\_r\_tmp*slope < 0.0) \textcolor{keywordflow}{then}} \DoxyCodeLine{526 } \DoxyCodeLine{527 u1\_l = u1\_l\_tmp} \DoxyCodeLine{528 u1\_r = 3.0*(u0\_r -\/ u0\_l)/h -\/ 2.0*u1\_l} \DoxyCodeLine{529 } \DoxyCodeLine{530 \textcolor{keywordflow}{else}} \DoxyCodeLine{531 } \DoxyCodeLine{532 u1\_l = u1\_l\_tmp} \DoxyCodeLine{533 u1\_r = u1\_r\_tmp} \DoxyCodeLine{534 } \DoxyCodeLine{535 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{536 } \DoxyCodeLine{537 \textcolor{keywordflow}{ endif} \textcolor{comment}{! end treating case with inflexion point on the left}} \DoxyCodeLine{538 } \DoxyCodeLine{539 \textcolor{comment}{! Move inflexion point on the right}} \DoxyCodeLine{540 \textcolor{keywordflow}{if} ( inflexion\_r ) \textcolor{keywordflow}{then}} \DoxyCodeLine{541 } \DoxyCodeLine{542 u1\_l\_tmp = 3.0*(u0\_r-\/u0\_l)/h -\/ 2.0*u1\_r} \DoxyCodeLine{543 u1\_r\_tmp = 1.5*(u0\_r-\/u0\_l)/h -\/ 0.5*u1\_l} \DoxyCodeLine{544 } \DoxyCodeLine{545 \textcolor{keywordflow}{if} ( (u1\_l\_tmp*slope < 0.0) .AND. (u1\_r\_tmp*slope < 0.0) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{546 } \DoxyCodeLine{547 u1\_l = 3.0 * (u0\_r -\/ u0\_l) / h} \DoxyCodeLine{548 u1\_r = 0.0} \DoxyCodeLine{549 } \DoxyCodeLine{550 \textcolor{keywordflow}{elseif} (u1\_l\_tmp*slope < 0.0) \textcolor{keywordflow}{then}} \DoxyCodeLine{551 } \DoxyCodeLine{552 u1\_r = u1\_r\_tmp} \DoxyCodeLine{553 u1\_l = 3.0*(u0\_r -\/ u0\_l)/h -\/ 2.0*u1\_r} \DoxyCodeLine{554 } \DoxyCodeLine{555 \textcolor{keywordflow}{elseif} (u1\_r\_tmp*slope < 0.0) \textcolor{keywordflow}{then}} \DoxyCodeLine{556 } \DoxyCodeLine{557 u1\_l = u1\_l\_tmp} \DoxyCodeLine{558 u1\_r = 1.5*(u0\_r -\/ u0\_l)/h -\/ 0.5*u1\_l} \DoxyCodeLine{559 } \DoxyCodeLine{560 \textcolor{keywordflow}{else}} \DoxyCodeLine{561 } \DoxyCodeLine{562 u1\_l = u1\_l\_tmp} \DoxyCodeLine{563 u1\_r = u1\_r\_tmp} \DoxyCodeLine{564 } \DoxyCodeLine{565 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{566 } \DoxyCodeLine{567 \textcolor{keywordflow}{ endif} \textcolor{comment}{! end treating case with inflexion point on the right}} \DoxyCodeLine{568 } \DoxyCodeLine{569 \textcolor{comment}{! Zero out negligibly small slopes.}} \DoxyCodeLine{570 \textcolor{keywordflow}{if} ( abs(u1\_l*h) < epsilon(u0\_l) * (abs(u0\_l) + abs(u0\_r)) ) u1\_l = 0.0} \DoxyCodeLine{571 \textcolor{keywordflow}{if} ( abs(u1\_r*h) < epsilon(u0\_l) * (abs(u0\_l) + abs(u0\_r)) ) u1\_r = 0.0} \DoxyCodeLine{572 } \end{DoxyCode} \mbox{\Hypertarget{namespacep3m__functions_ab70eb9e69fc6586e1d6a371d2eeb44d1}\label{namespacep3m__functions_ab70eb9e69fc6586e1d6a371d2eeb44d1}} \index{p3m\_functions@{p3m\_functions}!p3m\_boundary\_extrapolation@{p3m\_boundary\_extrapolation}} \index{p3m\_boundary\_extrapolation@{p3m\_boundary\_extrapolation}!p3m\_functions@{p3m\_functions}} \doxysubsubsection{\texorpdfstring{p3m\_boundary\_extrapolation()}{p3m\_boundary\_extrapolation()}} {\footnotesize\ttfamily subroutine, public p3m\+\_\+functions\+::p3m\+\_\+boundary\+\_\+extrapolation (\begin{DoxyParamCaption}\item[{integer, intent(in)}]{N, }\item[{real, dimension(\+:), intent(in)}]{h, }\item[{real, dimension(\+:), intent(in)}]{u, }\item[{real, dimension(\+:,\+:), intent(inout)}]{edge\+\_\+values, }\item[{real, dimension(\+:,\+:), intent(inout)}]{ppoly\+\_\+S, }\item[{real, dimension(\+:,\+:), intent(inout)}]{ppoly\+\_\+coef, }\item[{real, intent(in), optional}]{h\+\_\+neglect, }\item[{real, intent(in), optional}]{h\+\_\+neglect\+\_\+edge }\end{DoxyParamCaption})} Calculate the edge values and slopes at boundary cells as part of building a piecewise cubic sub-\/grid scale profiles. The following explanations apply to the left boundary cell. The same reasoning holds for the right boundary cell. A cubic needs to be built in the cell and requires four degrees of freedom, which are the edge values and slopes. The right edge values and slopes are taken to be that of the neighboring cell (i.\+e., the left edge value and slope of the neighboring cell). The left edge value and slope are determined by computing the parabola based on the cell average and the right edge value and slope. The resulting cubic is not necessarily monotonic and the slopes are subsequently modified to yield a monotonic cubic. \begin{DoxyParams}[1]{Parameters} \mbox{\texttt{ in}} & {\em n} & Number of cells \\ \hline \mbox{\texttt{ in}} & {\em h} & cell widths (size N) \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em u} & cell averages (size N) in arbitrary units \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em edge\+\_\+values} & Edge value of polynomial \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em ppoly\+\_\+s} & Edge slope of polynomial \mbox{[}A H-\/1\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em ppoly\+\_\+coef} & Coefficients of polynomial \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em h\+\_\+neglect} & A negligibly small width for the purpose of cell reconstructions \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em h\+\_\+neglect\+\_\+edge} & A negligibly small width for the purpose of finding edge values \mbox{[}H\mbox{]} \\ \hline \end{DoxyParams} Definition at line 191 of file P3\+M\+\_\+functions.\+F90. \begin{DoxyCode}{0} \DoxyCodeLine{193 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: N\textcolor{comment}{ !< Number of cells}} \DoxyCodeLine{194 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< cell widths (size N) [H]}} \DoxyCodeLine{195 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: u\textcolor{comment}{ !< cell averages (size N) in arbitrary units [A]}} \DoxyCodeLine{196 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: edge\_values\textcolor{comment}{ !< Edge value of polynomial [A]}} \DoxyCodeLine{197 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_S\textcolor{comment}{ !< Edge slope of polynomial [A H-\/1]}} \DoxyCodeLine{198 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_coef\textcolor{comment}{ !< Coefficients of polynomial [A]}} \DoxyCodeLine{199 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: h\_neglect\textcolor{comment}{ !< A negligibly small width for the}} \DoxyCodeLine{200 \textcolor{comment}{ !! purpose of cell reconstructions [H]}} \DoxyCodeLine{201 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: h\_neglect\_edge\textcolor{comment}{ !< A negligibly small width}} \DoxyCodeLine{202 \textcolor{comment}{ !! for the purpose of finding edge values [H]}} \DoxyCodeLine{203 \textcolor{comment}{! Local variables}} \DoxyCodeLine{204 \textcolor{keywordtype}{integer} :: i0, i1} \DoxyCodeLine{205 \textcolor{keywordtype}{logical} :: monotonic \textcolor{comment}{! boolean indicating whether the cubic is monotonic}} \DoxyCodeLine{206 \textcolor{keywordtype}{ real} :: u0, u1 \textcolor{comment}{! Values of u in two adjacent cells [A]}} \DoxyCodeLine{207 \textcolor{keywordtype}{ real} :: h0, h1 \textcolor{comment}{! Values of h in two adjacent cells, plus a smal increment [H]}} \DoxyCodeLine{208 \textcolor{keywordtype}{ real} :: b, c, d \textcolor{comment}{! Temporary variables [A]}} \DoxyCodeLine{209 \textcolor{keywordtype}{ real} :: u0\_l, u0\_r \textcolor{comment}{! Left and right edge values [A]}} \DoxyCodeLine{210 \textcolor{keywordtype}{ real} :: u1\_l, u1\_r \textcolor{comment}{! Left and right edge slopes [A H-\/1]}} \DoxyCodeLine{211 \textcolor{keywordtype}{ real} :: slope \textcolor{comment}{! The cell center slope [A H-\/1]}} \DoxyCodeLine{212 \textcolor{keywordtype}{ real} :: hNeglect, hNeglect\_edge \textcolor{comment}{! Negligibly small thickness [H]}} \DoxyCodeLine{213 } \DoxyCodeLine{214 hneglect = hneglect\_dflt ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(h\_neglect)) hneglect = h\_neglect} \DoxyCodeLine{215 hneglect\_edge = hneglect\_edge\_dflt ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(h\_neglect\_edge)) hneglect\_edge = h\_neglect\_edge} \DoxyCodeLine{216 } \DoxyCodeLine{217 \textcolor{comment}{! -\/-\/-\/-\/-\/ Left boundary -\/-\/-\/-\/-\/}} \DoxyCodeLine{218 i0 = 1} \DoxyCodeLine{219 i1 = 2} \DoxyCodeLine{220 h0 = h(i0) + hneglect\_edge} \DoxyCodeLine{221 h1 = h(i1) + hneglect\_edge} \DoxyCodeLine{222 u0 = u(i0)} \DoxyCodeLine{223 u1 = u(i1)} \DoxyCodeLine{224 } \DoxyCodeLine{225 \textcolor{comment}{! Compute the left edge slope in neighboring cell and express it in}} \DoxyCodeLine{226 \textcolor{comment}{! the global coordinate system}} \DoxyCodeLine{227 b = ppoly\_coef(i1,2)} \DoxyCodeLine{228 u1\_r = b / h1 \textcolor{comment}{! derivative evaluated at xi = 0.0, expressed w.r.t. x}} \DoxyCodeLine{229 } \DoxyCodeLine{230 \textcolor{comment}{! Limit the right slope by the PLM limited slope}} \DoxyCodeLine{231 slope = 2.0 * ( u1 -\/ u0 ) / ( h0 + hneglect )} \DoxyCodeLine{232 \textcolor{keywordflow}{if} ( abs(u1\_r) > abs(slope) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{233 u1\_r = slope} \DoxyCodeLine{234 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{235 } \DoxyCodeLine{236 \textcolor{comment}{! The right edge value in the boundary cell is taken to be the left}} \DoxyCodeLine{237 \textcolor{comment}{! edge value in the neighboring cell}} \DoxyCodeLine{238 u0\_r = edge\_values(i1,1)} \DoxyCodeLine{239 } \DoxyCodeLine{240 \textcolor{comment}{! Given the right edge value and slope, we determine the left}} \DoxyCodeLine{241 \textcolor{comment}{! edge value and slope by computing the parabola as determined by}} \DoxyCodeLine{242 \textcolor{comment}{! the right edge value and slope and the boundary cell average}} \DoxyCodeLine{243 u0\_l = 3.0 * u0 + 0.5 * h0*u1\_r -\/ 2.0 * u0\_r} \DoxyCodeLine{244 u1\_l = ( -\/ 6.0 * u0 -\/ 2.0 * h0*u1\_r + 6.0 * u0\_r) / ( h0 + hneglect )} \DoxyCodeLine{245 } \DoxyCodeLine{246 \textcolor{comment}{! Check whether the edge values are monotonic. For example, if the left edge}} \DoxyCodeLine{247 \textcolor{comment}{! value is larger than the right edge value while the slope is positive, the}} \DoxyCodeLine{248 \textcolor{comment}{! edge values are inconsistent and we need to modify the left edge value}} \DoxyCodeLine{249 \textcolor{keywordflow}{if} ( (u0\_r-\/u0\_l) * slope < 0.0 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{250 u0\_l = u0\_r} \DoxyCodeLine{251 u1\_l = 0.0} \DoxyCodeLine{252 u1\_r = 0.0} \DoxyCodeLine{253 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{254 } \DoxyCodeLine{255 \textcolor{comment}{! Store edge values and slope, build cubic and check monotonicity}} \DoxyCodeLine{256 edge\_values(i0,1) = u0\_l} \DoxyCodeLine{257 edge\_values(i0,2) = u0\_r} \DoxyCodeLine{258 ppoly\_s(i0,1) = u1\_l} \DoxyCodeLine{259 ppoly\_s(i0,2) = u1\_r} \DoxyCodeLine{260 } \DoxyCodeLine{261 \textcolor{comment}{! Store edge values and slope, build cubic and check monotonicity}} \DoxyCodeLine{262 \textcolor{keyword}{call }build\_cubic\_interpolant( h, i0, edge\_values, ppoly\_s, ppoly\_coef )} \DoxyCodeLine{263 monotonic = is\_cubic\_monotonic( ppoly\_coef, i0 )} \DoxyCodeLine{264 } \DoxyCodeLine{265 \textcolor{keywordflow}{if} ( .not.monotonic ) \textcolor{keywordflow}{then}} \DoxyCodeLine{266 \textcolor{keyword}{call }monotonize\_cubic( h0, u0\_l, u0\_r, 0.0, slope, slope, u1\_l, u1\_r )} \DoxyCodeLine{267 } \DoxyCodeLine{268 \textcolor{comment}{! Rebuild cubic after monotonization}} \DoxyCodeLine{269 ppoly\_s(i0,1) = u1\_l} \DoxyCodeLine{270 ppoly\_s(i0,2) = u1\_r} \DoxyCodeLine{271 \textcolor{keyword}{call }build\_cubic\_interpolant( h, i0, edge\_values, ppoly\_s, ppoly\_coef )} \DoxyCodeLine{272 } \DoxyCodeLine{273 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{274 } \DoxyCodeLine{275 \textcolor{comment}{! -\/-\/-\/-\/-\/ Right boundary -\/-\/-\/-\/-\/}} \DoxyCodeLine{276 i0 = n-\/1} \DoxyCodeLine{277 i1 = n} \DoxyCodeLine{278 h0 = h(i0) + hneglect\_edge} \DoxyCodeLine{279 h1 = h(i1) + hneglect\_edge} \DoxyCodeLine{280 u0 = u(i0)} \DoxyCodeLine{281 u1 = u(i1)} \DoxyCodeLine{282 } \DoxyCodeLine{283 \textcolor{comment}{! Compute the right edge slope in neighboring cell and express it in}} \DoxyCodeLine{284 \textcolor{comment}{! the global coordinate system}} \DoxyCodeLine{285 b = ppoly\_coef(i0,2)} \DoxyCodeLine{286 c = ppoly\_coef(i0,3)} \DoxyCodeLine{287 d = ppoly\_coef(i0,4)} \DoxyCodeLine{288 u1\_l = (b + 2*c + 3*d) / ( h0 + hneglect ) \textcolor{comment}{! derivative evaluated at xi = 1.0}} \DoxyCodeLine{289 } \DoxyCodeLine{290 \textcolor{comment}{! Limit the left slope by the PLM limited slope}} \DoxyCodeLine{291 slope = 2.0 * ( u1 -\/ u0 ) / ( h1 + hneglect )} \DoxyCodeLine{292 \textcolor{keywordflow}{if} ( abs(u1\_l) > abs(slope) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{293 u1\_l = slope} \DoxyCodeLine{294 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{295 } \DoxyCodeLine{296 \textcolor{comment}{! The left edge value in the boundary cell is taken to be the right}} \DoxyCodeLine{297 \textcolor{comment}{! edge value in the neighboring cell}} \DoxyCodeLine{298 u0\_l = edge\_values(i0,2)} \DoxyCodeLine{299 } \DoxyCodeLine{300 \textcolor{comment}{! Given the left edge value and slope, we determine the right}} \DoxyCodeLine{301 \textcolor{comment}{! edge value and slope by computing the parabola as determined by}} \DoxyCodeLine{302 \textcolor{comment}{! the left edge value and slope and the boundary cell average}} \DoxyCodeLine{303 u0\_r = 3.0 * u1 -\/ 0.5 * h1*u1\_l -\/ 2.0 * u0\_l} \DoxyCodeLine{304 u1\_r = ( 6.0 * u1 -\/ 2.0 * h1*u1\_l -\/ 6.0 * u0\_l) / ( h1 + hneglect )} \DoxyCodeLine{305 } \DoxyCodeLine{306 \textcolor{comment}{! Check whether the edge values are monotonic. For example, if the right edge}} \DoxyCodeLine{307 \textcolor{comment}{! value is smaller than the left edge value while the slope is positive, the}} \DoxyCodeLine{308 \textcolor{comment}{! edge values are inconsistent and we need to modify the right edge value}} \DoxyCodeLine{309 \textcolor{keywordflow}{if} ( (u0\_r-\/u0\_l) * slope < 0.0 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{310 u0\_r = u0\_l} \DoxyCodeLine{311 u1\_l = 0.0} \DoxyCodeLine{312 u1\_r = 0.0} \DoxyCodeLine{313 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{314 } \DoxyCodeLine{315 \textcolor{comment}{! Store edge values and slope, build cubic and check monotonicity}} \DoxyCodeLine{316 edge\_values(i1,1) = u0\_l} \DoxyCodeLine{317 edge\_values(i1,2) = u0\_r} \DoxyCodeLine{318 ppoly\_s(i1,1) = u1\_l} \DoxyCodeLine{319 ppoly\_s(i1,2) = u1\_r} \DoxyCodeLine{320 } \DoxyCodeLine{321 \textcolor{keyword}{call }build\_cubic\_interpolant( h, i1, edge\_values, ppoly\_s, ppoly\_coef )} \DoxyCodeLine{322 monotonic = is\_cubic\_monotonic( ppoly\_coef, i1 )} \DoxyCodeLine{323 } \DoxyCodeLine{324 \textcolor{keywordflow}{if} ( .not.monotonic ) \textcolor{keywordflow}{then}} \DoxyCodeLine{325 \textcolor{keyword}{call }monotonize\_cubic( h1, u0\_l, u0\_r, slope, 0.0, slope, u1\_l, u1\_r )} \DoxyCodeLine{326 } \DoxyCodeLine{327 \textcolor{comment}{! Rebuild cubic after monotonization}} \DoxyCodeLine{328 ppoly\_s(i1,1) = u1\_l} \DoxyCodeLine{329 ppoly\_s(i1,2) = u1\_r} \DoxyCodeLine{330 \textcolor{keyword}{call }build\_cubic\_interpolant( h, i1, edge\_values, ppoly\_s, ppoly\_coef )} \DoxyCodeLine{331 } \DoxyCodeLine{332 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{333 } \end{DoxyCode} \mbox{\Hypertarget{namespacep3m__functions_a685d3ef292536dae810b2059ccaa5819}\label{namespacep3m__functions_a685d3ef292536dae810b2059ccaa5819}} \index{p3m\_functions@{p3m\_functions}!p3m\_interpolation@{p3m\_interpolation}} \index{p3m\_interpolation@{p3m\_interpolation}!p3m\_functions@{p3m\_functions}} \doxysubsubsection{\texorpdfstring{p3m\_interpolation()}{p3m\_interpolation()}} {\footnotesize\ttfamily subroutine, public p3m\+\_\+functions\+::p3m\+\_\+interpolation (\begin{DoxyParamCaption}\item[{integer, intent(in)}]{N, }\item[{real, dimension(\+:), intent(in)}]{h, }\item[{real, dimension(\+:), intent(in)}]{u, }\item[{real, dimension(\+:,\+:), intent(inout)}]{edge\+\_\+values, }\item[{real, dimension(\+:,\+:), intent(inout)}]{ppoly\+\_\+S, }\item[{real, dimension(\+:,\+:), intent(inout)}]{ppoly\+\_\+coef, }\item[{real, intent(in), optional}]{h\+\_\+neglect, }\item[{logical, intent(in), optional}]{answers\+\_\+2018 }\end{DoxyParamCaption})} Set up a piecewise cubic interpolation from cell averages and estimated edge slopes and values. Cubic interpolation between edges. The edge values and slopes M\+U\+ST have been estimated prior to calling this routine. It is assumed that the size of the array \textquotesingle{}u\textquotesingle{} is equal to the number of cells defining \textquotesingle{}grid\textquotesingle{} and \textquotesingle{}ppoly\textquotesingle{}. No consistency check is performed here. \begin{DoxyParams}[1]{Parameters} \mbox{\texttt{ in}} & {\em n} & Number of cells \\ \hline \mbox{\texttt{ in}} & {\em h} & cell widths (size N) \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em u} & cell averages (size N) in arbitrary units \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em edge\+\_\+values} & Edge value of polynomial \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em ppoly\+\_\+s} & Edge slope of polynomial \mbox{[}A H-\/1\mbox{]}. \\ \hline \mbox{\texttt{ in,out}} & {\em ppoly\+\_\+coef} & Coefficients of polynomial \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em h\+\_\+neglect} & A negligibly small width for the purpose of cell reconstructions \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em answers\+\_\+2018} & If true use older, less acccurate expressions. \\ \hline \end{DoxyParams} Definition at line 28 of file P3\+M\+\_\+functions.\+F90. \begin{DoxyCode}{0} \DoxyCodeLine{29 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: N\textcolor{comment}{ !< Number of cells}} \DoxyCodeLine{30 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< cell widths (size N) [H]}} \DoxyCodeLine{31 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: u\textcolor{comment}{ !< cell averages (size N) in arbitrary units [A]}} \DoxyCodeLine{32 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: edge\_values\textcolor{comment}{ !< Edge value of polynomial [A]}} \DoxyCodeLine{33 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_S\textcolor{comment}{ !< Edge slope of polynomial [A H-\/1].}} \DoxyCodeLine{34 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_coef\textcolor{comment}{ !< Coefficients of polynomial [A]}} \DoxyCodeLine{35 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: h\_neglect\textcolor{comment}{ !< A negligibly small width for the}} \DoxyCodeLine{36 \textcolor{comment}{ !! purpose of cell reconstructions [H]}} \DoxyCodeLine{37 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: answers\_2018\textcolor{comment}{ !< If true use older, less acccurate expressions.}} \DoxyCodeLine{38 } \DoxyCodeLine{39 \textcolor{comment}{! Call the limiter for p3m, which takes care of everything from}} \DoxyCodeLine{40 \textcolor{comment}{! computing the coefficients of the cubic to monotonizing it.}} \DoxyCodeLine{41 \textcolor{comment}{! This routine could be called directly instead of having to call}} \DoxyCodeLine{42 \textcolor{comment}{! 'P3M\_interpolation' first but we do that to provide an homogeneous}} \DoxyCodeLine{43 \textcolor{comment}{! interface.}} \DoxyCodeLine{44 \textcolor{keyword}{call }p3m\_limiter( n, h, u, edge\_values, ppoly\_s, ppoly\_coef, h\_neglect, answers\_2018 )} \DoxyCodeLine{45 } \end{DoxyCode} \mbox{\Hypertarget{namespacep3m__functions_a44db59cb5df7f139e05b745746342fcf}\label{namespacep3m__functions_a44db59cb5df7f139e05b745746342fcf}} \index{p3m\_functions@{p3m\_functions}!p3m\_limiter@{p3m\_limiter}} \index{p3m\_limiter@{p3m\_limiter}!p3m\_functions@{p3m\_functions}} \doxysubsubsection{\texorpdfstring{p3m\_limiter()}{p3m\_limiter()}} {\footnotesize\ttfamily subroutine p3m\+\_\+functions\+::p3m\+\_\+limiter (\begin{DoxyParamCaption}\item[{integer, intent(in)}]{N, }\item[{real, dimension(\+:), intent(in)}]{h, }\item[{real, dimension(\+:), intent(in)}]{u, }\item[{real, dimension(\+:,\+:), intent(inout)}]{edge\+\_\+values, }\item[{real, dimension(\+:,\+:), intent(inout)}]{ppoly\+\_\+S, }\item[{real, dimension(\+:,\+:), intent(inout)}]{ppoly\+\_\+coef, }\item[{real, intent(in), optional}]{h\+\_\+neglect, }\item[{logical, intent(in), optional}]{answers\+\_\+2018 }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} Adust a piecewise cubic reconstruction with a limiter that adjusts the edge values and slopes. The p3m limiter operates as follows\+: \begin{DoxyEnumerate} \item Edge values are bounded \item Discontinuous edge values are systematically averaged \item Loop on cells and do the following a. Build cubic curve b. Check if cubic curve is monotonic c. If not, monotonize cubic curve and rebuild it \end{DoxyEnumerate} Step 3 of the monotonization process leaves all edge values unchanged. \begin{DoxyParams}[1]{Parameters} \mbox{\texttt{ in}} & {\em n} & Number of cells \\ \hline \mbox{\texttt{ in}} & {\em h} & cell widths (size N) \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em u} & cell averages (size N) in arbitrary units \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em edge\+\_\+values} & Edge value of polynomial \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em ppoly\+\_\+s} & Edge slope of polynomial \mbox{[}A H-\/1\mbox{]} \\ \hline \mbox{\texttt{ in,out}} & {\em ppoly\+\_\+coef} & Coefficients of polynomial \mbox{[}A\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em h\+\_\+neglect} & A negligibly small width for the purpose of cell reconstructions \mbox{[}H\mbox{]} \\ \hline \mbox{\texttt{ in}} & {\em answers\+\_\+2018} & If true use older, less acccurate expressions. \\ \hline \end{DoxyParams} Definition at line 61 of file P3\+M\+\_\+functions.\+F90. \begin{DoxyCode}{0} \DoxyCodeLine{62 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: N\textcolor{comment}{ !< Number of cells}} \DoxyCodeLine{63 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< cell widths (size N) [H]}} \DoxyCodeLine{64 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:)}, \textcolor{keywordtype}{intent(in)} :: u\textcolor{comment}{ !< cell averages (size N) in arbitrary units [A]}} \DoxyCodeLine{65 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: edge\_values\textcolor{comment}{ !< Edge value of polynomial [A]}} \DoxyCodeLine{66 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_S\textcolor{comment}{ !< Edge slope of polynomial [A H-\/1]}} \DoxyCodeLine{67 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_coef\textcolor{comment}{ !< Coefficients of polynomial [A]}} \DoxyCodeLine{68 \textcolor{keywordtype}{ real}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: h\_neglect\textcolor{comment}{ !< A negligibly small width for}} \DoxyCodeLine{69 \textcolor{comment}{ !! the purpose of cell reconstructions [H]}} \DoxyCodeLine{70 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: answers\_2018\textcolor{comment}{ !< If true use older, less acccurate expressions.}} \DoxyCodeLine{71 } \DoxyCodeLine{72 \textcolor{comment}{! Local variables}} \DoxyCodeLine{73 \textcolor{keywordtype}{integer} :: k \textcolor{comment}{! loop index}} \DoxyCodeLine{74 \textcolor{keywordtype}{logical} :: monotonic \textcolor{comment}{! boolean indicating whether the cubic is monotonic}} \DoxyCodeLine{75 \textcolor{keywordtype}{ real} :: u0\_l, u0\_r \textcolor{comment}{! edge values [A]}} \DoxyCodeLine{76 \textcolor{keywordtype}{ real} :: u1\_l, u1\_r \textcolor{comment}{! edge slopes [A H-\/1]}} \DoxyCodeLine{77 \textcolor{keywordtype}{ real} :: u\_l, u\_c, u\_r \textcolor{comment}{! left, center and right cell averages [A]}} \DoxyCodeLine{78 \textcolor{keywordtype}{ real} :: h\_l, h\_c, h\_r \textcolor{comment}{! left, center and right cell widths [H]}} \DoxyCodeLine{79 \textcolor{keywordtype}{ real} :: sigma\_l, sigma\_c, sigma\_r \textcolor{comment}{! left, center and right van Leer slopes [A H-\/1]}} \DoxyCodeLine{80 \textcolor{keywordtype}{ real} :: slope \textcolor{comment}{! retained PLM slope [A H-\/1]}} \DoxyCodeLine{81 \textcolor{keywordtype}{ real} :: eps} \DoxyCodeLine{82 \textcolor{keywordtype}{ real} :: hNeglect \textcolor{comment}{! A negligibly small thickness [H]}} \DoxyCodeLine{83 } \DoxyCodeLine{84 hneglect = hneglect\_dflt ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(h\_neglect)) hneglect = h\_neglect} \DoxyCodeLine{85 } \DoxyCodeLine{86 eps = 1e-\/10} \DoxyCodeLine{87 } \DoxyCodeLine{88 \textcolor{comment}{! 1. Bound edge values (boundary cells are assumed to be local extrema)}} \DoxyCodeLine{89 \textcolor{keyword}{call }bound\_edge\_values( n, h, u, edge\_values, hneglect, answers\_2018 )} \DoxyCodeLine{90 } \DoxyCodeLine{91 \textcolor{comment}{! 2. Systematically average discontinuous edge values}} \DoxyCodeLine{92 \textcolor{keyword}{call }average\_discontinuous\_edge\_values( n, edge\_values )} \DoxyCodeLine{93 } \DoxyCodeLine{94 } \DoxyCodeLine{95 \textcolor{comment}{! 3. Loop on cells and do the following}} \DoxyCodeLine{96 \textcolor{comment}{! (a) Build cubic curve}} \DoxyCodeLine{97 \textcolor{comment}{! (b) Check if cubic curve is monotonic}} \DoxyCodeLine{98 \textcolor{comment}{! (c) If not, monotonize cubic curve and rebuild it}} \DoxyCodeLine{99 \textcolor{keywordflow}{do} k = 1,n} \DoxyCodeLine{100 } \DoxyCodeLine{101 \textcolor{comment}{! Get edge values, edge slopes and cell width}} \DoxyCodeLine{102 u0\_l = edge\_values(k,1)} \DoxyCodeLine{103 u0\_r = edge\_values(k,2)} \DoxyCodeLine{104 u1\_l = ppoly\_s(k,1)} \DoxyCodeLine{105 u1\_r = ppoly\_s(k,2)} \DoxyCodeLine{106 } \DoxyCodeLine{107 \textcolor{comment}{! Get cell widths and cell averages (boundary cells are assumed to}} \DoxyCodeLine{108 \textcolor{comment}{! be local extrema for the sake of slopes)}} \DoxyCodeLine{109 u\_c = u(k)} \DoxyCodeLine{110 h\_c = h(k)} \DoxyCodeLine{111 } \DoxyCodeLine{112 \textcolor{keywordflow}{if} ( k == 1 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{113 h\_l = h(k)} \DoxyCodeLine{114 u\_l = u(k)} \DoxyCodeLine{115 \textcolor{keywordflow}{else}} \DoxyCodeLine{116 h\_l = h(k-\/1)} \DoxyCodeLine{117 u\_l = u(k-\/1)} \DoxyCodeLine{118 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{119 } \DoxyCodeLine{120 \textcolor{keywordflow}{if} ( k == n ) \textcolor{keywordflow}{then}} \DoxyCodeLine{121 h\_r = h(k)} \DoxyCodeLine{122 u\_r = u(k)} \DoxyCodeLine{123 \textcolor{keywordflow}{else}} \DoxyCodeLine{124 h\_r = h(k+1)} \DoxyCodeLine{125 u\_r = u(k+1)} \DoxyCodeLine{126 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{127 } \DoxyCodeLine{128 \textcolor{comment}{! Compute limited slope}} \DoxyCodeLine{129 sigma\_l = 2.0 * ( u\_c -\/ u\_l ) / ( h\_c + hneglect )} \DoxyCodeLine{130 sigma\_c = 2.0 * ( u\_r -\/ u\_l ) / ( h\_l + 2.0*h\_c + h\_r + hneglect )} \DoxyCodeLine{131 sigma\_r = 2.0 * ( u\_r -\/ u\_c ) / ( h\_c + hneglect )} \DoxyCodeLine{132 } \DoxyCodeLine{133 \textcolor{keywordflow}{if} ( (sigma\_l * sigma\_r) > 0.0 ) \textcolor{keywordflow}{then}} \DoxyCodeLine{134 slope = sign( min(abs(sigma\_l),abs(sigma\_c),abs(sigma\_r)), sigma\_c )} \DoxyCodeLine{135 \textcolor{keywordflow}{else}} \DoxyCodeLine{136 slope = 0.0} \DoxyCodeLine{137 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{138 } \DoxyCodeLine{139 \textcolor{comment}{! If the slopes are small, set them to zero to prevent asymmetric representation near extrema.}} \DoxyCodeLine{140 \textcolor{keywordflow}{if} ( abs(u1\_l*h\_c) < epsilon(u\_c)*abs(u\_c) ) u1\_l = 0.0} \DoxyCodeLine{141 \textcolor{keywordflow}{if} ( abs(u1\_r*h\_c) < epsilon(u\_c)*abs(u\_c) ) u1\_r = 0.0} \DoxyCodeLine{142 } \DoxyCodeLine{143 \textcolor{comment}{! The edge slopes are limited from above by the respective}} \DoxyCodeLine{144 \textcolor{comment}{! one-\/sided slopes}} \DoxyCodeLine{145 \textcolor{keywordflow}{if} ( abs(u1\_l) > abs(sigma\_l) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{146 u1\_l = sigma\_l} \DoxyCodeLine{147 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{148 } \DoxyCodeLine{149 \textcolor{keywordflow}{if} ( abs(u1\_r) > abs(sigma\_r) ) \textcolor{keywordflow}{then}} \DoxyCodeLine{150 u1\_r = sigma\_r} \DoxyCodeLine{151 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{152 } \DoxyCodeLine{153 \textcolor{comment}{! Build cubic interpolant (compute the coefficients)}} \DoxyCodeLine{154 \textcolor{keyword}{call }build\_cubic\_interpolant( h, k, edge\_values, ppoly\_s, ppoly\_coef )} \DoxyCodeLine{155 } \DoxyCodeLine{156 \textcolor{comment}{! Check whether cubic is monotonic}} \DoxyCodeLine{157 monotonic = is\_cubic\_monotonic( ppoly\_coef, k )} \DoxyCodeLine{158 } \DoxyCodeLine{159 \textcolor{comment}{! If cubic is not monotonic, monotonize it by modifiying the}} \DoxyCodeLine{160 \textcolor{comment}{! edge slopes, store the new edge slopes and recompute the}} \DoxyCodeLine{161 \textcolor{comment}{! cubic coefficients}} \DoxyCodeLine{162 \textcolor{keywordflow}{if} ( .not.monotonic ) \textcolor{keywordflow}{then}} \DoxyCodeLine{163 \textcolor{keyword}{call }monotonize\_cubic( h\_c, u0\_l, u0\_r, sigma\_l, sigma\_r, slope, u1\_l, u1\_r )} \DoxyCodeLine{164 \textcolor{keywordflow}{ endif}} \DoxyCodeLine{165 } \DoxyCodeLine{166 \textcolor{comment}{! Store edge slopes}} \DoxyCodeLine{167 ppoly\_s(k,1) = u1\_l} \DoxyCodeLine{168 ppoly\_s(k,2) = u1\_r} \DoxyCodeLine{169 } \DoxyCodeLine{170 \textcolor{comment}{! Recompute coefficients of cubic}} \DoxyCodeLine{171 \textcolor{keyword}{call }build\_cubic\_interpolant( h, k, edge\_values, ppoly\_s, ppoly\_coef )} \DoxyCodeLine{172 } \DoxyCodeLine{173 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! loop on cells}} \DoxyCodeLine{174 } \end{DoxyCode}