\hypertarget{namespaceplm__functions}{}\section{plm\+\_\+functions Module Reference}
\label{namespaceplm__functions}\index{plm\+\_\+functions@{plm\+\_\+functions}}


\subsection{Detailed Description}
Piecewise linear reconstruction functions. 

Date of creation\+: 2008.\+06.\+06 L. White.

This module contains routines that handle one-\/dimensionnal finite volume reconstruction using the piecewise linear method (P\+LM). \subsection*{Functions/\+Subroutines}
\begin{DoxyCompactItemize}
\item 
real elemental pure function, public \mbox{\hyperlink{namespaceplm__functions_a072affa78922591148b954fa63872246}{plm\+\_\+slope\+\_\+wa}} (h\+\_\+l, h\+\_\+c, h\+\_\+r, h\+\_\+neglect, u\+\_\+l, u\+\_\+c, u\+\_\+r)
\begin{DoxyCompactList}\small\item\em Returns a limited P\+LM slope following White and Adcroft, 2008. \mbox{[}units of u\mbox{]} Note that this is not the same as the Colella and Woodward method. \end{DoxyCompactList}\item 
real elemental pure function, public \mbox{\hyperlink{namespaceplm__functions_a29d2baec39200c5aa448c79181a16f4f}{plm\+\_\+slope\+\_\+cw}} (h\+\_\+l, h\+\_\+c, h\+\_\+r, h\+\_\+neglect, u\+\_\+l, u\+\_\+c, u\+\_\+r)
\begin{DoxyCompactList}\small\item\em Returns a limited P\+LM slope following Colella and Woodward 1984. \end{DoxyCompactList}\item 
real elemental pure function, public \mbox{\hyperlink{namespaceplm__functions_a497e5e73108e08afcb4e2186710ae094}{plm\+\_\+monotonized\+\_\+slope}} (u\+\_\+l, u\+\_\+c, u\+\_\+r, s\+\_\+l, s\+\_\+c, s\+\_\+r)
\begin{DoxyCompactList}\small\item\em Returns a limited P\+LM slope following Colella and Woodward 1984. \end{DoxyCompactList}\item 
real elemental pure function, public \mbox{\hyperlink{namespaceplm__functions_a42fbf62545902eeb6c9e035763496b07}{plm\+\_\+extrapolate\+\_\+slope}} (h\+\_\+l, h\+\_\+c, h\+\_\+neglect, u\+\_\+l, u\+\_\+c)
\begin{DoxyCompactList}\small\item\em Returns a P\+LM slope using h2 extrapolation from a cell to the left. Use the negative to extrapolate from the a cell to the right. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespaceplm__functions_a0e3a6bedfb3064ed8fae8c8449649054}{plm\+\_\+reconstruction}} (N, h, u, edge\+\_\+values, ppoly\+\_\+coef, h\+\_\+neglect)
\begin{DoxyCompactList}\small\item\em Reconstruction by linear polynomials within each cell. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespaceplm__functions_a85e83fd8ad314548d861d5048dad8f72}{plm\+\_\+boundary\+\_\+extrapolation}} (N, h, u, edge\+\_\+values, ppoly\+\_\+coef, h\+\_\+neglect)
\begin{DoxyCompactList}\small\item\em Reconstruction by linear polynomials within boundary cells. \end{DoxyCompactList}\end{DoxyCompactItemize}
\subsection*{Variables}
\begin{DoxyCompactItemize}
\item 
\mbox{\Hypertarget{namespaceplm__functions_abc6a0167b8f90d1d1a95c0dec6921053}\label{namespaceplm__functions_abc6a0167b8f90d1d1a95c0dec6921053}} 
real, parameter \mbox{\hyperlink{namespaceplm__functions_abc6a0167b8f90d1d1a95c0dec6921053}{hneglect\+\_\+dflt}} = 1.E-\/30
\begin{DoxyCompactList}\small\item\em Default negligible cell thickness. \end{DoxyCompactList}\end{DoxyCompactItemize}


\subsection{Function/\+Subroutine Documentation}
\mbox{\Hypertarget{namespaceplm__functions_a85e83fd8ad314548d861d5048dad8f72}\label{namespaceplm__functions_a85e83fd8ad314548d861d5048dad8f72}} 
\index{plm\+\_\+functions@{plm\+\_\+functions}!plm\+\_\+boundary\+\_\+extrapolation@{plm\+\_\+boundary\+\_\+extrapolation}}
\index{plm\+\_\+boundary\+\_\+extrapolation@{plm\+\_\+boundary\+\_\+extrapolation}!plm\+\_\+functions@{plm\+\_\+functions}}
\subsubsection{\texorpdfstring{plm\+\_\+boundary\+\_\+extrapolation()}{plm\_boundary\_extrapolation()}}
{\footnotesize\ttfamily subroutine, public plm\+\_\+functions\+::plm\+\_\+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\+\_\+coef,  }\item[{real, intent(in), optional}]{h\+\_\+neglect }\end{DoxyParamCaption})}



Reconstruction by linear polynomials within boundary cells. 

The left and right edge values in the left and right boundary cells, respectively, are estimated using a linear extrapolation within the cells.

This extrapolation is E\+X\+A\+CT when the underlying profile is linear.

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{\tt in}  & {\em n} & Number of cells\\
\hline
\mbox{\tt in}  & {\em h} & cell widths (size N)\\
\hline
\mbox{\tt in}  & {\em u} & cell averages (size N)\\
\hline
\mbox{\tt in,out}  & {\em edge\+\_\+values} & edge values of piecewise polynomials, with the same units as u.\\
\hline
\mbox{\tt in,out}  & {\em ppoly\+\_\+coef} & coefficients of piecewise polynomials, mainly with the same units as u.\\
\hline
\mbox{\tt in}  & {\em h\+\_\+neglect} & A negligibly small width for the purpose of cell reconstructions in the same units as h \\
\hline
\end{DoxyParams}


Definition at line 269 of file P\+L\+M\+\_\+functions.\+F90.


\begin{DoxyCode}
269   \textcolor{keywordtype}{integer},              \textcolor{keywordtype}{intent(in)}    :: N\textcolor{comment}{ !< Number of cells}
270   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:)},   \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{ !< cell widths (size N)}
271   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:)},   \textcolor{keywordtype}{intent(in)}    :: u\textcolor{comment}{ !< cell averages (size N)}
272   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: edge\_values\textcolor{comment}{ !< edge values of piecewise polynomials,}
273 \textcolor{comment}{                                           !! with the same units as u.}
274   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_coef\textcolor{comment}{ !< coefficients of piecewise polynomials, mainly}
275 \textcolor{comment}{                                           !! with the same units as u.}
276   \textcolor{keywordtype}{real},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}    :: h\_neglect\textcolor{comment}{ !< A negligibly small width for}
277 \textcolor{comment}{                                           !! the purpose of cell reconstructions}
278 \textcolor{comment}{                                           !! in the same units as h}
279   \textcolor{comment}{! Local variables}
280   \textcolor{keywordtype}{real}    :: slope                \textcolor{comment}{! retained PLM slope}
281   \textcolor{keywordtype}{real}    :: hNeglect
282 
283   hneglect = hneglect\_dflt ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(h\_neglect)) hneglect = h\_neglect
284 
285   \textcolor{comment}{! Extrapolate from 2 to 1 to estimate slope}
286   slope = - plm\_extrapolate\_slope( h(2), h(1), hneglect, u(2), u(1) )
287 
288   edge\_values(1,1) = u(1) - 0.5 * slope
289   edge\_values(1,2) = u(1) + 0.5 * slope
290 
291   ppoly\_coef(1,1) = edge\_values(1,1)
292   ppoly\_coef(1,2) = edge\_values(1,2) - edge\_values(1,1)
293 
294   \textcolor{comment}{! Extrapolate from N-1 to N to estimate slope}
295   slope = plm\_extrapolate\_slope( h(n-1), h(n), hneglect, u(n-1), u(n) )
296 
297   edge\_values(n,1) = u(n) - 0.5 * slope
298   edge\_values(n,2) = u(n) + 0.5 * slope
299 
300   ppoly\_coef(n,1) = edge\_values(n,1)
301   ppoly\_coef(n,2) = edge\_values(n,2) - edge\_values(n,1)
302 
\end{DoxyCode}
\mbox{\Hypertarget{namespaceplm__functions_a42fbf62545902eeb6c9e035763496b07}\label{namespaceplm__functions_a42fbf62545902eeb6c9e035763496b07}} 
\index{plm\+\_\+functions@{plm\+\_\+functions}!plm\+\_\+extrapolate\+\_\+slope@{plm\+\_\+extrapolate\+\_\+slope}}
\index{plm\+\_\+extrapolate\+\_\+slope@{plm\+\_\+extrapolate\+\_\+slope}!plm\+\_\+functions@{plm\+\_\+functions}}
\subsubsection{\texorpdfstring{plm\+\_\+extrapolate\+\_\+slope()}{plm\_extrapolate\_slope()}}
{\footnotesize\ttfamily real elemental pure function, public plm\+\_\+functions\+::plm\+\_\+extrapolate\+\_\+slope (\begin{DoxyParamCaption}\item[{real, intent(in)}]{h\+\_\+l,  }\item[{real, intent(in)}]{h\+\_\+c,  }\item[{real, intent(in)}]{h\+\_\+neglect,  }\item[{real, intent(in)}]{u\+\_\+l,  }\item[{real, intent(in)}]{u\+\_\+c }\end{DoxyParamCaption})}



Returns a P\+LM slope using h2 extrapolation from a cell to the left. Use the negative to extrapolate from the a cell to the right. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em h\+\_\+l} & Thickness of left cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+c} & Thickness of center cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+neglect} & A negligible thickness \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+l} & Value of left cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+c} & Value of center cell \mbox{[}units of u\mbox{]} \\
\hline
\end{DoxyParams}


Definition at line 161 of file P\+L\+M\+\_\+functions.\+F90.


\begin{DoxyCode}
161   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_l\textcolor{comment}{ !< Thickness of left cell [units of grid thickness]}
162   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_c\textcolor{comment}{ !< Thickness of center cell [units of grid thickness]}
163   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_neglect\textcolor{comment}{ !< A negligible thickness [units of grid thickness]}
164   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_l\textcolor{comment}{ !< Value of left cell [units of u]}
165   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_c\textcolor{comment}{ !< Value of center cell [units of u]}
166   \textcolor{comment}{! Local variables}
167   \textcolor{keywordtype}{real} :: left\_edge \textcolor{comment}{! Left edge value [units of u]}
168   \textcolor{keywordtype}{real} :: hl, hc \textcolor{comment}{! Left and central cell thicknesses [units of grid thickness]}
169 
170   \textcolor{comment}{! Avoid division by zero for vanished cells}
171    hl = h\_l + h\_neglect
172    hc = h\_c + h\_neglect
173 
174   \textcolor{comment}{! The h2 scheme is used to compute the left edge value}
175   left\_edge = (u\_l*hc + u\_c*hl) / (hl + hc)
176 
177   plm\_extrapolate\_slope = 2.0 * ( u\_c - left\_edge )
178 
\end{DoxyCode}
\mbox{\Hypertarget{namespaceplm__functions_a497e5e73108e08afcb4e2186710ae094}\label{namespaceplm__functions_a497e5e73108e08afcb4e2186710ae094}} 
\index{plm\+\_\+functions@{plm\+\_\+functions}!plm\+\_\+monotonized\+\_\+slope@{plm\+\_\+monotonized\+\_\+slope}}
\index{plm\+\_\+monotonized\+\_\+slope@{plm\+\_\+monotonized\+\_\+slope}!plm\+\_\+functions@{plm\+\_\+functions}}
\subsubsection{\texorpdfstring{plm\+\_\+monotonized\+\_\+slope()}{plm\_monotonized\_slope()}}
{\footnotesize\ttfamily real elemental pure function, public plm\+\_\+functions\+::plm\+\_\+monotonized\+\_\+slope (\begin{DoxyParamCaption}\item[{real, intent(in)}]{u\+\_\+l,  }\item[{real, intent(in)}]{u\+\_\+c,  }\item[{real, intent(in)}]{u\+\_\+r,  }\item[{real, intent(in)}]{s\+\_\+l,  }\item[{real, intent(in)}]{s\+\_\+c,  }\item[{real, intent(in)}]{s\+\_\+r }\end{DoxyParamCaption})}



Returns a limited P\+LM slope following Colella and Woodward 1984. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em u\+\_\+l} & Value of left cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+c} & Value of center cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+r} & Value of right cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+l} & P\+LM slope of left cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+c} & P\+LM slope of center cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+r} & P\+LM slope of right cell \mbox{[}units of u\mbox{]} \\
\hline
\end{DoxyParams}


Definition at line 122 of file P\+L\+M\+\_\+functions.\+F90.


\begin{DoxyCode}
122   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_l\textcolor{comment}{ !< Value of left cell [units of u]}
123   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_c\textcolor{comment}{ !< Value of center cell [units of u]}
124   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_r\textcolor{comment}{ !< Value of right cell [units of u]}
125   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: s\_l\textcolor{comment}{ !< PLM slope of left cell [units of u]}
126   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: s\_c\textcolor{comment}{ !< PLM slope of center cell [units of u]}
127   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: s\_r\textcolor{comment}{ !< PLM slope of right cell [units of u]}
128   \textcolor{comment}{! Local variables}
129   \textcolor{keywordtype}{real} :: e\_r, e\_l, edge \textcolor{comment}{! Right, left and temporary edge values [units of u]}
130   \textcolor{keywordtype}{real} :: almost\_two \textcolor{comment}{! The number 2, almost.}
131   \textcolor{keywordtype}{real} :: slp \textcolor{comment}{! Magnitude of PLM central slope [units of u]}
132 
133   almost\_two = 2. * ( 1. - epsilon(s\_c) )
134 
135   \textcolor{comment}{! Edge values of neighbors abutting this cell}
136   e\_r = u\_l + 0.5*s\_l
137   e\_l = u\_r - 0.5*s\_r
138   slp = abs(s\_c)
139 
140   \textcolor{comment}{! Check that left edge is between right edge of cell to the left and this cell mean}
141   edge = u\_c - 0.5 * s\_c
142   \textcolor{keywordflow}{if} ( ( edge - e\_r ) * ( u\_c - edge ) < 0. ) \textcolor{keywordflow}{then}
143     edge = 0.5 * ( edge + e\_r )
144     slp = min( slp, abs( edge - u\_c ) * almost\_two )
145 \textcolor{keywordflow}{  endif}
146 
147   \textcolor{comment}{! Check that right edge is between left edge of cell to the right and this cell mean}
148   edge = u\_c + 0.5 * s\_c
149   \textcolor{keywordflow}{if} ( ( edge - u\_c ) * ( e\_l - edge ) < 0. ) \textcolor{keywordflow}{then}
150     edge = 0.5 * ( edge + e\_l )
151     slp = min( slp, abs( edge - u\_c ) * almost\_two )
152 \textcolor{keywordflow}{  endif}
153 
154   plm\_monotonized\_slope = sign( slp, s\_c )
155 
\end{DoxyCode}
\mbox{\Hypertarget{namespaceplm__functions_a0e3a6bedfb3064ed8fae8c8449649054}\label{namespaceplm__functions_a0e3a6bedfb3064ed8fae8c8449649054}} 
\index{plm\+\_\+functions@{plm\+\_\+functions}!plm\+\_\+reconstruction@{plm\+\_\+reconstruction}}
\index{plm\+\_\+reconstruction@{plm\+\_\+reconstruction}!plm\+\_\+functions@{plm\+\_\+functions}}
\subsubsection{\texorpdfstring{plm\+\_\+reconstruction()}{plm\_reconstruction()}}
{\footnotesize\ttfamily subroutine, public plm\+\_\+functions\+::plm\+\_\+reconstruction (\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\+\_\+coef,  }\item[{real, intent(in), optional}]{h\+\_\+neglect }\end{DoxyParamCaption})}



Reconstruction by linear polynomials within each cell. 

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{\tt in}  & {\em n} & Number of cells\\
\hline
\mbox{\tt in}  & {\em h} & cell widths (size N)\\
\hline
\mbox{\tt in}  & {\em u} & cell averages (size N)\\
\hline
\mbox{\tt in,out}  & {\em edge\+\_\+values} & edge values of piecewise polynomials, with the same units as u.\\
\hline
\mbox{\tt in,out}  & {\em ppoly\+\_\+coef} & coefficients of piecewise polynomials, mainly with the same units as u.\\
\hline
\mbox{\tt in}  & {\em h\+\_\+neglect} & A negligibly small width for the purpose of cell reconstructions in the same units as h \\
\hline
\end{DoxyParams}


Definition at line 187 of file P\+L\+M\+\_\+functions.\+F90.


\begin{DoxyCode}
187   \textcolor{keywordtype}{integer},              \textcolor{keywordtype}{intent(in)}    :: N\textcolor{comment}{ !< Number of cells}
188   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:)},   \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{ !< cell widths (size N)}
189   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:)},   \textcolor{keywordtype}{intent(in)}    :: u\textcolor{comment}{ !< cell averages (size N)}
190   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: edge\_values\textcolor{comment}{ !< edge values of piecewise polynomials,}
191 \textcolor{comment}{                                           !! with the same units as u.}
192   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{intent(inout)} :: ppoly\_coef\textcolor{comment}{ !< coefficients of piecewise polynomials, mainly}
193 \textcolor{comment}{                                           !! with the same units as u.}
194   \textcolor{keywordtype}{real},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}    :: h\_neglect\textcolor{comment}{ !< A negligibly small width for}
195 \textcolor{comment}{                                           !! the purpose of cell reconstructions}
196 \textcolor{comment}{                                           !! in the same units as h}
197 
198   \textcolor{comment}{! Local variables}
199   \textcolor{keywordtype}{integer}       :: k                    \textcolor{comment}{! loop index}
200   \textcolor{keywordtype}{real}          :: u\_l, u\_c, u\_r        \textcolor{comment}{! left, center and right cell averages}
201   \textcolor{keywordtype}{real}          :: h\_l, h\_c, h\_r, h\_cn  \textcolor{comment}{! left, center and right cell widths}
202   \textcolor{keywordtype}{real}          :: slope                \textcolor{comment}{! retained PLM slope}
203   \textcolor{keywordtype}{real}          :: a, b                 \textcolor{comment}{! auxiliary variables}
204   \textcolor{keywordtype}{real}          :: u\_min, u\_max, e\_l, e\_r, edge
205   \textcolor{keywordtype}{real}          :: almost\_one
206   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(N)} :: slp, mslp
207   \textcolor{keywordtype}{real}    :: hNeglect
208 
209   hneglect = hneglect\_dflt ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(h\_neglect)) hneglect = h\_neglect
210 
211   almost\_one = 1. - epsilon(slope)
212 
213   \textcolor{comment}{! Loop on interior cells}
214   \textcolor{keywordflow}{do} k = 2,n-1
215     slp(k) = plm\_slope\_wa(h(k-1), h(k), h(k+1), hneglect, u(k-1), u(k), u(k+1))
216 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! end loop on interior cells}
217 
218   \textcolor{comment}{! Boundary cells use PCM. Extrapolation is handled after monotonization.}
219   slp(1) = 0.
220   slp(n) = 0.
221 
222   \textcolor{comment}{! This loop adjusts the slope so that edge values are monotonic.}
223   \textcolor{keywordflow}{do} k = 2, n-1
224     mslp(k) = plm\_monotonized\_slope( u(k-1), u(k), u(k+1), slp(k-1), slp(k), slp(k+1) )
225 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! end loop on interior cells}
226   mslp(1) = 0.
227   mslp(n) = 0.
228 
229   \textcolor{comment}{! Store and return edge values and polynomial coefficients.}
230   edge\_values(1,1) = u(1)
231   edge\_values(1,2) = u(1)
232   ppoly\_coef(1,1) = u(1)
233   ppoly\_coef(1,2) = 0.
234   \textcolor{keywordflow}{do} k = 2, n-1
235     slope = mslp(k)
236     u\_l = u(k) - 0.5 * slope \textcolor{comment}{! Left edge value of cell k}
237     u\_r = u(k) + 0.5 * slope \textcolor{comment}{! Right edge value of cell k}
238 
239     edge\_values(k,1) = u\_l
240     edge\_values(k,2) = u\_r
241     ppoly\_coef(k,1) = u\_l
242     ppoly\_coef(k,2) = ( u\_r - u\_l )
243     \textcolor{comment}{! Check to see if this evaluation of the polynomial at x=1 would be}
244     \textcolor{comment}{! monotonic w.r.t. the next cell's edge value. If not, scale back!}
245     edge = ppoly\_coef(k,2) + ppoly\_coef(k,1)
246     e\_r = u(k+1) - 0.5 * sign( mslp(k+1), slp(k+1) )
247     \textcolor{keywordflow}{if} ( (edge-u(k))*(e\_r-edge)<0.) \textcolor{keywordflow}{then}
248       ppoly\_coef(k,2) = ppoly\_coef(k,2) * almost\_one
249 \textcolor{keywordflow}{    endif}
250 \textcolor{keywordflow}{  enddo}
251   edge\_values(n,1) = u(n)
252   edge\_values(n,2) = u(n)
253   ppoly\_coef(n,1) = u(n)
254   ppoly\_coef(n,2) = 0.
255 
\end{DoxyCode}
\mbox{\Hypertarget{namespaceplm__functions_a29d2baec39200c5aa448c79181a16f4f}\label{namespaceplm__functions_a29d2baec39200c5aa448c79181a16f4f}} 
\index{plm\+\_\+functions@{plm\+\_\+functions}!plm\+\_\+slope\+\_\+cw@{plm\+\_\+slope\+\_\+cw}}
\index{plm\+\_\+slope\+\_\+cw@{plm\+\_\+slope\+\_\+cw}!plm\+\_\+functions@{plm\+\_\+functions}}
\subsubsection{\texorpdfstring{plm\+\_\+slope\+\_\+cw()}{plm\_slope\_cw()}}
{\footnotesize\ttfamily real elemental pure function, public plm\+\_\+functions\+::plm\+\_\+slope\+\_\+cw (\begin{DoxyParamCaption}\item[{real, intent(in)}]{h\+\_\+l,  }\item[{real, intent(in)}]{h\+\_\+c,  }\item[{real, intent(in)}]{h\+\_\+r,  }\item[{real, intent(in)}]{h\+\_\+neglect,  }\item[{real, intent(in)}]{u\+\_\+l,  }\item[{real, intent(in)}]{u\+\_\+c,  }\item[{real, intent(in)}]{u\+\_\+r }\end{DoxyParamCaption})}



Returns a limited P\+LM slope following Colella and Woodward 1984. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em h\+\_\+l} & Thickness of left cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+c} & Thickness of center cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+r} & Thickness of right cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+neglect} & A negligible thickness \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+l} & Value of left cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+c} & Value of center cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+r} & Value of right cell \mbox{[}units of u\mbox{]} \\
\hline
\end{DoxyParams}


Definition at line 68 of file P\+L\+M\+\_\+functions.\+F90.


\begin{DoxyCode}
68   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_l\textcolor{comment}{ !< Thickness of left cell [units of grid thickness]}
69   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_c\textcolor{comment}{ !< Thickness of center cell [units of grid thickness]}
70   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_r\textcolor{comment}{ !< Thickness of right cell [units of grid thickness]}
71   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_neglect\textcolor{comment}{ !< A negligible thickness [units of grid thickness]}
72   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_l\textcolor{comment}{ !< Value of left cell [units of u]}
73   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_c\textcolor{comment}{ !< Value of center cell [units of u]}
74   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_r\textcolor{comment}{ !< Value of right cell [units of u]}
75   \textcolor{comment}{! Local variables}
76   \textcolor{keywordtype}{real} :: sigma\_l, sigma\_c, sigma\_r \textcolor{comment}{! Left, central and right slope estimates as}
77                                     \textcolor{comment}{! differences across the cell [units of u]}
78   \textcolor{keywordtype}{real} :: u\_min, u\_max \textcolor{comment}{! Minimum and maximum value across cell [units of u]}
79   \textcolor{keywordtype}{real} :: h\_cn \textcolor{comment}{! Thickness of center cell [units of grid thickness]}
80 
81   h\_cn = h\_c + h\_neglect
82 
83   \textcolor{comment}{! Side differences}
84   sigma\_r = u\_r - u\_c
85   sigma\_l = u\_c - u\_l
86 
87   \textcolor{comment}{! This is the second order slope given by equation 1.7 of}
88   \textcolor{comment}{! Piecewise Parabolic Method, Colella and Woodward (1984),}
89   \textcolor{comment}{! http://dx.doi.org/10.1016/0021-991(84)90143-8.}
90   \textcolor{comment}{! For uniform resolution it simplifies to ( u\_r - u\_l )/2 .}
91   sigma\_c = ( h\_c / ( h\_cn + ( h\_l + h\_r ) ) ) * ( &
92                 ( 2.*h\_l + h\_c ) / ( h\_r + h\_cn ) * sigma\_r &
93               + ( 2.*h\_r + h\_c ) / ( h\_l + h\_cn ) * sigma\_l )
94 
95   \textcolor{comment}{! Limit slope so that reconstructions are bounded by neighbors}
96   u\_min = min( u\_l, u\_c, u\_r )
97   u\_max = max( u\_l, u\_c, u\_r )
98   \textcolor{keywordflow}{if} ( (sigma\_l * sigma\_r) > 0.0 ) \textcolor{keywordflow}{then}
99     \textcolor{comment}{! This limits the slope so that the edge values are bounded by the}
100     \textcolor{comment}{! two cell averages spanning the edge.}
101     plm\_slope\_cw = sign( min( abs(sigma\_c), 2.*min( u\_c - u\_min, u\_max - u\_c ) ), sigma\_c )
102   \textcolor{keywordflow}{else}
103     \textcolor{comment}{! Extrema in the mean values require a PCM reconstruction avoid generating}
104     \textcolor{comment}{! larger extreme values.}
105     plm\_slope\_cw = 0.0
106 \textcolor{keywordflow}{  endif}
107 
108   \textcolor{comment}{! This block tests to see if roundoff causes edge values to be out of bounds}
109   \textcolor{keywordflow}{if} (u\_c - 0.5*abs(plm\_slope\_cw) < u\_min .or.  u\_c + 0.5*abs(plm\_slope\_cw) > u\_max) \textcolor{keywordflow}{then}
110     plm\_slope\_cw = plm\_slope\_cw * ( 1. - epsilon(plm\_slope\_cw) )
111 \textcolor{keywordflow}{  endif}
112 
113   \textcolor{comment}{! An attempt to avoid inconsistency when the values become unrepresentable.}
114   \textcolor{comment}{! ### The following 1.E-140 is dimensionally inconsistent. A newer version of}
115   \textcolor{comment}{! PLM is progress that will avoid the need for such rounding.}
116   \textcolor{keywordflow}{if} (abs(plm\_slope\_cw) < 1.e-140) plm\_slope\_cw = 0.
117 
\end{DoxyCode}
\mbox{\Hypertarget{namespaceplm__functions_a072affa78922591148b954fa63872246}\label{namespaceplm__functions_a072affa78922591148b954fa63872246}} 
\index{plm\+\_\+functions@{plm\+\_\+functions}!plm\+\_\+slope\+\_\+wa@{plm\+\_\+slope\+\_\+wa}}
\index{plm\+\_\+slope\+\_\+wa@{plm\+\_\+slope\+\_\+wa}!plm\+\_\+functions@{plm\+\_\+functions}}
\subsubsection{\texorpdfstring{plm\+\_\+slope\+\_\+wa()}{plm\_slope\_wa()}}
{\footnotesize\ttfamily real elemental pure function, public plm\+\_\+functions\+::plm\+\_\+slope\+\_\+wa (\begin{DoxyParamCaption}\item[{real, intent(in)}]{h\+\_\+l,  }\item[{real, intent(in)}]{h\+\_\+c,  }\item[{real, intent(in)}]{h\+\_\+r,  }\item[{real, intent(in)}]{h\+\_\+neglect,  }\item[{real, intent(in)}]{u\+\_\+l,  }\item[{real, intent(in)}]{u\+\_\+c,  }\item[{real, intent(in)}]{u\+\_\+r }\end{DoxyParamCaption})}



Returns a limited P\+LM slope following White and Adcroft, 2008. \mbox{[}units of u\mbox{]} Note that this is not the same as the Colella and Woodward method. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em h\+\_\+l} & Thickness of left cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+c} & Thickness of center cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+r} & Thickness of right cell \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h\+\_\+neglect} & A negligible thickness \mbox{[}units of grid thickness\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+l} & Value of left cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+c} & Value of center cell \mbox{[}units of u\mbox{]}\\
\hline
\mbox{\tt in}  & {\em u\+\_\+r} & Value of right cell \mbox{[}units of u\mbox{]} \\
\hline
\end{DoxyParams}


Definition at line 22 of file P\+L\+M\+\_\+functions.\+F90.


\begin{DoxyCode}
22   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_l\textcolor{comment}{ !< Thickness of left cell [units of grid thickness]}
23   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_c\textcolor{comment}{ !< Thickness of center cell [units of grid thickness]}
24   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_r\textcolor{comment}{ !< Thickness of right cell [units of grid thickness]}
25   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: h\_neglect\textcolor{comment}{ !< A negligible thickness [units of grid thickness]}
26   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_l\textcolor{comment}{ !< Value of left cell [units of u]}
27   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_c\textcolor{comment}{ !< Value of center cell [units of u]}
28   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: u\_r\textcolor{comment}{ !< Value of right cell [units of u]}
29   \textcolor{comment}{! Local variables}
30   \textcolor{keywordtype}{real} :: sigma\_l, sigma\_c, sigma\_r \textcolor{comment}{! Left, central and right slope estimates as}
31                                     \textcolor{comment}{! differences across the cell [units of u]}
32   \textcolor{keywordtype}{real} :: u\_min, u\_max \textcolor{comment}{! Minimum and maximum value across cell [units of u]}
33 
34   \textcolor{comment}{! Side differences}
35   sigma\_r = u\_r - u\_c
36   sigma\_l = u\_c - u\_l
37 
38   \textcolor{comment}{! Quasi-second order difference}
39   sigma\_c = 2.0 * ( u\_r - u\_l ) * ( h\_c / ( h\_l + 2.0*h\_c + h\_r + h\_neglect) )
40 
41   \textcolor{comment}{! Limit slope so that reconstructions are bounded by neighbors}
42   u\_min = min( u\_l, u\_c, u\_r )
43   u\_max = max( u\_l, u\_c, u\_r )
44   \textcolor{keywordflow}{if} ( (sigma\_l * sigma\_r) > 0.0 ) \textcolor{keywordflow}{then}
45     \textcolor{comment}{! This limits the slope so that the edge values are bounded by the}
46     \textcolor{comment}{! two cell averages spanning the edge.}
47     plm\_slope\_wa = sign( min( abs(sigma\_c), 2.*min( u\_c - u\_min, u\_max - u\_c ) ), sigma\_c )
48   \textcolor{keywordflow}{else}
49     \textcolor{comment}{! Extrema in the mean values require a PCM reconstruction avoid generating}
50     \textcolor{comment}{! larger extreme values.}
51     plm\_slope\_wa = 0.0
52 \textcolor{keywordflow}{  endif}
53 
54   \textcolor{comment}{! This block tests to see if roundoff causes edge values to be out of bounds}
55   \textcolor{keywordflow}{if} (u\_c - 0.5*abs(plm\_slope\_wa) < u\_min .or.  u\_c + 0.5*abs(plm\_slope\_wa) > u\_max) \textcolor{keywordflow}{then}
56     plm\_slope\_wa = plm\_slope\_wa * ( 1. - epsilon(plm\_slope\_wa) )
57 \textcolor{keywordflow}{  endif}
58 
59   \textcolor{comment}{! An attempt to avoid inconsistency when the values become unrepresentable.}
60   \textcolor{comment}{! ### The following 1.E-140 is dimensionally inconsistent. A newer version of}
61   \textcolor{comment}{! PLM is progress that will avoid the need for such rounding.}
62   \textcolor{keywordflow}{if} (abs(plm\_slope\_wa) < 1.e-140) plm\_slope\_wa = 0.
63 
\end{DoxyCode}
