\hypertarget{namespacep3m__functions}{}\section{p3m\+\_\+functions Module Reference}
\label{namespacep3m__functions}\index{p3m\_functions@{p3m\_functions}}


\subsection{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. \subsection*{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}
\subsection*{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}


\subsection{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}}
\subsubsection{\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 344 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}}
\subsubsection{\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 386 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}}
\subsubsection{\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 436 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}}
\subsubsection{\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 193 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}}
\subsubsection{\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 29 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}}
\subsubsection{\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 62 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}
