\hypertarget{namespacemom__density__integrals}{}\section{mom\+\_\+density\+\_\+integrals Module Reference}
\label{namespacemom__density__integrals}\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}


\subsection{Detailed Description}
Provides integrals of density. \subsection*{Functions/\+Subroutines}
\begin{DoxyCompactItemize}
\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_ac36ae5f4af2d02df0a1adf41b762e017}{int\+\_\+density\+\_\+dz}} (T, S, z\+\_\+t, z\+\_\+b, rho\+\_\+ref, rho\+\_\+0, G\+\_\+e, HI, E\+OS, US, dpa, intz\+\_\+dpa, intx\+\_\+dpa, inty\+\_\+dpa, bathyT, dz\+\_\+neglect, use\+Mass\+Wght\+Interp)
\begin{DoxyCompactList}\small\item\em Calls the appropriate subroutine to calculate analytical and nearly-\/analytical integrals in z across layers of pressure anomalies, which are required for calculating the finite-\/volume form pressure accelerations in a Boussinesq model. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_ad0ff54518bfacf7beebbe4dba687a914}{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+pcm}} (T, S, z\+\_\+t, z\+\_\+b, rho\+\_\+ref, rho\+\_\+0, G\+\_\+e, HI, E\+OS, US, dpa, intz\+\_\+dpa, intx\+\_\+dpa, inty\+\_\+dpa, bathyT, dz\+\_\+neglect, use\+Mass\+Wght\+Interp)
\begin{DoxyCompactList}\small\item\em Calculates (by numerical quadrature) integrals of pressure anomalies across layers, which are required for calculating the finite-\/volume form pressure accelerations in a Boussinesq model. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_a454d4d62eb599716cc3389fb8aa90b4b}{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+plm}} (k, tv, T\+\_\+t, T\+\_\+b, S\+\_\+t, S\+\_\+b, e, rho\+\_\+ref, rho\+\_\+0, G\+\_\+e, dz\+\_\+subroundoff, bathyT, HI, GV, E\+OS, US, dpa, intz\+\_\+dpa, intx\+\_\+dpa, inty\+\_\+dpa, use\+Mass\+Wght\+Interp)
\begin{DoxyCompactList}\small\item\em Compute pressure gradient force integrals by quadrature for the case where T and S are linear profiles. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_a3bf090e5cbb58811b3dd0ee4de80f999}{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+ppm}} (k, tv, T\+\_\+t, T\+\_\+b, S\+\_\+t, S\+\_\+b, e, rho\+\_\+ref, rho\+\_\+0, G\+\_\+e, dz\+\_\+subroundoff, bathyT, HI, GV, E\+OS, US, dpa, intz\+\_\+dpa, intx\+\_\+dpa, inty\+\_\+dpa, use\+Mass\+Wght\+Interp)
\begin{DoxyCompactList}\small\item\em Compute pressure gradient force integrals for layer \char`\"{}k\char`\"{} and the case where T and S are parabolic profiles. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_a759c2ae7aec17c59d532050f68a1e518}{int\+\_\+specific\+\_\+vol\+\_\+dp}} (T, S, p\+\_\+t, p\+\_\+b, alpha\+\_\+ref, HI, E\+OS, US, dza, intp\+\_\+dza, intx\+\_\+dza, inty\+\_\+dza, halo\+\_\+size, bathyP, d\+P\+\_\+tiny, use\+Mass\+Wght\+Interp)
\begin{DoxyCompactList}\small\item\em Calls the appropriate subroutine to calculate analytical and nearly-\/analytical integrals in pressure across layers of geopotential anomalies, which are required for calculating the finite-\/volume form pressure accelerations in a non-\/\+Boussinesq model. There are essentially no free assumptions, apart from the use of Boole\textquotesingle{}s rule to do the horizontal integrals, and from a truncation in the series for log(1-\/eps/1+eps) that assumes that $\vert$eps$\vert$ $<$ 0.\+34. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_a86ebdfebaaea2ac0339dcabed482dd11}{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+pcm}} (T, S, p\+\_\+t, p\+\_\+b, alpha\+\_\+ref, HI, E\+OS, US, dza, intp\+\_\+dza, intx\+\_\+dza, inty\+\_\+dza, halo\+\_\+size, bathyP, d\+P\+\_\+neglect, use\+Mass\+Wght\+Interp)
\begin{DoxyCompactList}\small\item\em This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-\/volume form pressure accelerations in a non-\/\+Boussinesq model. There are essentially no free assumptions, apart from the use of Boole\textquotesingle{}s rule quadrature to do the integrals. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_a71a518f80b5f2ecf83fd9c638cad140d}{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+plm}} (T\+\_\+t, T\+\_\+b, S\+\_\+t, S\+\_\+b, p\+\_\+t, p\+\_\+b, alpha\+\_\+ref, d\+P\+\_\+neglect, bathyP, HI, E\+OS, US, dza, intp\+\_\+dza, intx\+\_\+dza, inty\+\_\+dza, use\+Mass\+Wght\+Interp)
\begin{DoxyCompactList}\small\item\em This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-\/volume form pressure accelerations in a non-\/\+Boussinesq model. There are essentially no free assumptions, apart from the use of Boole\textquotesingle{}s rule quadrature to do the integrals. \end{DoxyCompactList}\item 
subroutine, public \mbox{\hyperlink{namespacemom__density__integrals_af9946b0e15d53d8f1fd1ab1c1f81f3a2}{find\+\_\+depth\+\_\+of\+\_\+pressure\+\_\+in\+\_\+cell}} (T\+\_\+t, T\+\_\+b, S\+\_\+t, S\+\_\+b, z\+\_\+t, z\+\_\+b, P\+\_\+t, P\+\_\+tgt, rho\+\_\+ref, G\+\_\+e, E\+OS, US, P\+\_\+b, z\+\_\+out, z\+\_\+tol)
\begin{DoxyCompactList}\small\item\em Find the depth at which the reconstructed pressure matches P\+\_\+tgt. \end{DoxyCompactList}\item 
real function \mbox{\hyperlink{namespacemom__density__integrals_abb8d1f2c24def176a27ef2d30fa373df}{frac\+\_\+dp\+\_\+at\+\_\+pos}} (T\+\_\+t, T\+\_\+b, S\+\_\+t, S\+\_\+b, z\+\_\+t, z\+\_\+b, rho\+\_\+ref, G\+\_\+e, pos, E\+OS)
\begin{DoxyCompactList}\small\item\em Returns change in anomalous pressure change from top to non-\/dimensional position pos between z\+\_\+t and z\+\_\+b. \end{DoxyCompactList}\end{DoxyCompactItemize}


\subsection{Function/\+Subroutine Documentation}
\mbox{\Hypertarget{namespacemom__density__integrals_af9946b0e15d53d8f1fd1ab1c1f81f3a2}\label{namespacemom__density__integrals_af9946b0e15d53d8f1fd1ab1c1f81f3a2}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!find\+\_\+depth\+\_\+of\+\_\+pressure\+\_\+in\+\_\+cell@{find\+\_\+depth\+\_\+of\+\_\+pressure\+\_\+in\+\_\+cell}}
\index{find\+\_\+depth\+\_\+of\+\_\+pressure\+\_\+in\+\_\+cell@{find\+\_\+depth\+\_\+of\+\_\+pressure\+\_\+in\+\_\+cell}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{find\+\_\+depth\+\_\+of\+\_\+pressure\+\_\+in\+\_\+cell()}{find\_depth\_of\_pressure\_in\_cell()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::find\+\_\+depth\+\_\+of\+\_\+pressure\+\_\+in\+\_\+cell (\begin{DoxyParamCaption}\item[{real, intent(in)}]{T\+\_\+t,  }\item[{real, intent(in)}]{T\+\_\+b,  }\item[{real, intent(in)}]{S\+\_\+t,  }\item[{real, intent(in)}]{S\+\_\+b,  }\item[{real, intent(in)}]{z\+\_\+t,  }\item[{real, intent(in)}]{z\+\_\+b,  }\item[{real, intent(in)}]{P\+\_\+t,  }\item[{real, intent(in)}]{P\+\_\+tgt,  }\item[{real, intent(in)}]{rho\+\_\+ref,  }\item[{real, intent(in)}]{G\+\_\+e,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, intent(out)}]{P\+\_\+b,  }\item[{real, intent(out)}]{z\+\_\+out,  }\item[{real, intent(in), optional}]{z\+\_\+tol }\end{DoxyParamCaption})}



Find the depth at which the reconstructed pressure matches P\+\_\+tgt. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em t\+\_\+t} & Potential temperature at the cell top \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em t\+\_\+b} & Potential temperature at the cell bottom \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+t} & Salinity at the cell top \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+b} & Salinity at the cell bottom \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+t} & Absolute height of top of cell \mbox{[}Z $\sim$$>$ m\mbox{]} (Boussinesq ????)\\
\hline
\mbox{\tt in}  & {\em z\+\_\+b} & Absolute height of bottom of cell \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+t} & Anomalous pressure of top of cell, relative to g$\ast$rho\+\_\+ref$\ast$z\+\_\+t \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+tgt} & Target pressure at height z\+\_\+out, relative to g$\ast$rho\+\_\+ref$\ast$z\+\_\+out \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+ref} & Reference density with which calculation are anomalous to \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]}\\
\hline
\mbox{\tt in}  & {\em g\+\_\+e} & Gravitational acceleration \mbox{[}L2 Z-\/1 T-\/2 $\sim$$>$ m s-\/2\mbox{]}\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt out}  & {\em p\+\_\+b} & Pressure at the bottom of the cell \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}\\
\hline
\mbox{\tt out}  & {\em z\+\_\+out} & Absolute depth at which anomalous pressure = p\+\_\+tgt \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+tol} & The tolerance in finding z\+\_\+out \mbox{[}Z $\sim$$>$ m\mbox{]} \\
\hline
\end{DoxyParams}


Definition at line 1527 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
1527   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: T\_t\textcolor{comment}{ !< Potential temperature at the cell top [degC]}
1528   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: T\_b\textcolor{comment}{ !< Potential temperature at the cell bottom [degC]}
1529   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: S\_t\textcolor{comment}{ !< Salinity at the cell top [ppt]}
1530   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: S\_b\textcolor{comment}{ !< Salinity at the cell bottom [ppt]}
1531   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: z\_t\textcolor{comment}{ !< Absolute height of top of cell [Z ~> m]   (Boussinesq ????)}
1532   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: z\_b\textcolor{comment}{ !< Absolute height of bottom of cell [Z ~> m]}
1533   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: P\_t\textcolor{comment}{ !< Anomalous pressure of top of cell, relative}
1534 \textcolor{comment}{                                            !! to g*rho\_ref*z\_t [R L2 T-2 ~> Pa]}
1535   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: P\_tgt\textcolor{comment}{ !< Target pressure at height z\_out, relative}
1536 \textcolor{comment}{                                            !! to g*rho\_ref*z\_out [R L2 T-2 ~> Pa]}
1537   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: rho\_ref\textcolor{comment}{ !< Reference density with which calculation}
1538 \textcolor{comment}{                                            !! are anomalous to [R ~> kg m-3]}
1539   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(in)}  :: G\_e\textcolor{comment}{ !< Gravitational acceleration [L2 Z-1 T-2 ~> m s-2]}
1540   \textcolor{keywordtype}{type}(EOS\_type),        \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
1541   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)}  :: US\textcolor{comment}{ !< A dimensional unit scaling type}
1542   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(out)} :: P\_b\textcolor{comment}{ !< Pressure at the bottom of the cell [R L2 T-2 ~> Pa]}
1543   \textcolor{keywordtype}{real},                  \textcolor{keywordtype}{intent(out)} :: z\_out\textcolor{comment}{ !< Absolute depth at which anomalous pressure = p\_tgt [Z ~>
       m]}
1544   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{optional},        \textcolor{keywordtype}{intent(in)}  :: z\_tol\textcolor{comment}{ !< The tolerance in finding z\_out [Z ~> m]}
1545 
1546   \textcolor{comment}{! Local variables}
1547   \textcolor{keywordtype}{real} :: dp    \textcolor{comment}{! Pressure thickness of the layer [R L2 T-2 ~> Pa]}
1548   \textcolor{keywordtype}{real} :: F\_guess, F\_l, F\_r  \textcolor{comment}{! Fractional positions [nondim]}
1549   \textcolor{keywordtype}{real} :: GxRho \textcolor{comment}{! The product of the gravitational acceleration and reference density [R L2 Z-1 T-2 ~> Pa
       m-1]}
1550   \textcolor{keywordtype}{real} :: Pa, Pa\_left, Pa\_right, Pa\_tol \textcolor{comment}{! Pressure anomalies, P = integral of g*(rho-rho\_ref) dz [R L2 T-2
       ~> Pa]}
1551   \textcolor{keywordtype}{character(len=240)} :: msg
1552 
1553   gxrho = g\_e * rho\_ref
1554 
1555   \textcolor{comment}{! Anomalous pressure difference across whole cell}
1556   dp = frac\_dp\_at\_pos(t\_t, t\_b, s\_t, s\_b, z\_t, z\_b, rho\_ref, g\_e, 1.0, eos)
1557 
1558   p\_b = p\_t + dp \textcolor{comment}{! Anomalous pressure at bottom of cell}
1559 
1560   \textcolor{keywordflow}{if} (p\_tgt <= p\_t ) \textcolor{keywordflow}{then}
1561     z\_out = z\_t
1562     \textcolor{keywordflow}{return}
1563 \textcolor{keywordflow}{  endif}
1564 
1565   \textcolor{keywordflow}{if} (p\_tgt >= p\_b) \textcolor{keywordflow}{then}
1566     z\_out = z\_b
1567     \textcolor{keywordflow}{return}
1568 \textcolor{keywordflow}{  endif}
1569 
1570   f\_l = 0.
1571   pa\_left = p\_t - p\_tgt \textcolor{comment}{! Pa\_left < 0}
1572   f\_r = 1.
1573   pa\_right = p\_b - p\_tgt \textcolor{comment}{! Pa\_right > 0}
1574   pa\_tol = gxrho * 1.0e-5*us%m\_to\_Z
1575   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(z\_tol)) pa\_tol = gxrho * z\_tol
1576 
1577   f\_guess = f\_l - pa\_left / (pa\_right - pa\_left) * (f\_r - f\_l)
1578   pa = pa\_right - pa\_left \textcolor{comment}{! To get into iterative loop}
1579   \textcolor{keywordflow}{do} \textcolor{keywordflow}{while} ( abs(pa) > pa\_tol )
1580 
1581     z\_out = z\_t + ( z\_b - z\_t ) * f\_guess
1582     pa = frac\_dp\_at\_pos(t\_t, t\_b, s\_t, s\_b, z\_t, z\_b, rho\_ref, g\_e, f\_guess, eos) - ( p\_tgt - p\_t )
1583 
1584     \textcolor{keywordflow}{if} (pa<pa\_left) \textcolor{keywordflow}{then}
1585       \textcolor{keyword}{write}(msg,*) pa\_left,pa,pa\_right,p\_t-p\_tgt,p\_b-p\_tgt
1586       \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{'find\_depth\_of\_pressure\_in\_cell out of bounds negative: /n'}//msg)
1587     \textcolor{keywordflow}{elseif} (pa<0.) \textcolor{keywordflow}{then}
1588       pa\_left = pa
1589       f\_l = f\_guess
1590     \textcolor{keywordflow}{elseif} (pa>pa\_right) \textcolor{keywordflow}{then}
1591       \textcolor{keyword}{write}(msg,*) pa\_left,pa,pa\_right,p\_t-p\_tgt,p\_b-p\_tgt
1592       \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{'find\_depth\_of\_pressure\_in\_cell out of bounds positive: /n'}//msg)
1593     \textcolor{keywordflow}{elseif} (pa>0.) \textcolor{keywordflow}{then}
1594       pa\_right = pa
1595       f\_r = f\_guess
1596     \textcolor{keywordflow}{else} \textcolor{comment}{! Pa == 0}
1597       \textcolor{keywordflow}{return}
1598 \textcolor{keywordflow}{    endif}
1599     f\_guess = f\_l - pa\_left / (pa\_right - pa\_left) * (f\_r - f\_l)
1600 
1601 \textcolor{keywordflow}{  enddo}
1602 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_abb8d1f2c24def176a27ef2d30fa373df}\label{namespacemom__density__integrals_abb8d1f2c24def176a27ef2d30fa373df}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!frac\+\_\+dp\+\_\+at\+\_\+pos@{frac\+\_\+dp\+\_\+at\+\_\+pos}}
\index{frac\+\_\+dp\+\_\+at\+\_\+pos@{frac\+\_\+dp\+\_\+at\+\_\+pos}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{frac\+\_\+dp\+\_\+at\+\_\+pos()}{frac\_dp\_at\_pos()}}
{\footnotesize\ttfamily real function mom\+\_\+density\+\_\+integrals\+::frac\+\_\+dp\+\_\+at\+\_\+pos (\begin{DoxyParamCaption}\item[{real, intent(in)}]{T\+\_\+t,  }\item[{real, intent(in)}]{T\+\_\+b,  }\item[{real, intent(in)}]{S\+\_\+t,  }\item[{real, intent(in)}]{S\+\_\+b,  }\item[{real, intent(in)}]{z\+\_\+t,  }\item[{real, intent(in)}]{z\+\_\+b,  }\item[{real, intent(in)}]{rho\+\_\+ref,  }\item[{real, intent(in)}]{G\+\_\+e,  }\item[{real, intent(in)}]{pos,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



Returns change in anomalous pressure change from top to non-\/dimensional position pos between z\+\_\+t and z\+\_\+b. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em t\+\_\+t} & Potential temperature at the cell top \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em t\+\_\+b} & Potential temperature at the cell bottom \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+t} & Salinity at the cell top \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+b} & Salinity at the cell bottom \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+t} & The geometric height at the top of the layer \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+b} & The geometric height at the bottom of the layer \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+ref} & A mean density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]}, that is subtracted out to reduce the magnitude of each of the integrals.\\
\hline
\mbox{\tt in}  & {\em g\+\_\+e} & The Earth\textquotesingle{}s gravitational acceleration \mbox{[}L2 Z-\/1 T-\/2 $\sim$$>$ m s-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em pos} & The fractional vertical position, 0 to 1 \mbox{[}nondim\mbox{]}\\
\hline
 & {\em eos} & Equation of state structure \\
\hline
\end{DoxyParams}


Definition at line 1609 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
1609   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: T\_t\textcolor{comment}{ !< Potential temperature at the cell top [degC]}
1610   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: T\_b\textcolor{comment}{ !< Potential temperature at the cell bottom [degC]}
1611   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: S\_t\textcolor{comment}{ !< Salinity at the cell top [ppt]}
1612   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: S\_b\textcolor{comment}{ !< Salinity at the cell bottom [ppt]}
1613   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: z\_t\textcolor{comment}{ !< The geometric height at the top of the layer [Z ~> m]}
1614   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: z\_b\textcolor{comment}{ !< The geometric height at the bottom of the layer [Z ~> m]}
1615   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: rho\_ref\textcolor{comment}{ !< A mean density [R ~> kg m-3], that is subtracted out to}
1616 \textcolor{comment}{                                     !! reduce the magnitude of each of the integrals.}
1617   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: G\_e\textcolor{comment}{ !< The Earth's gravitational acceleration [L2 Z-1 T-2 ~> m s-2]}
1618   \textcolor{keywordtype}{real},           \textcolor{keywordtype}{intent(in)}  :: pos\textcolor{comment}{ !< The fractional vertical position, 0 to 1 [nondim]}
1619   \textcolor{keywordtype}{type}(EOS\_type), \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
1620   \textcolor{keywordtype}{real}                        :: fract\_dp\_at\_pos\textcolor{comment}{ !< The change in pressure from the layer top to}
1621 \textcolor{comment}{                                     !! fractional position pos [R L2 T-2 ~> Pa]}
1622   \textcolor{comment}{! Local variables}
1623   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_90 = 1.0/90.0  \textcolor{comment}{! A rational constant [nondim]}
1624   \textcolor{keywordtype}{real} :: dz                 \textcolor{comment}{! Distance from the layer top [Z ~> m]}
1625   \textcolor{keywordtype}{real} :: top\_weight, bottom\_weight \textcolor{comment}{! Fractional weights at quadrature points [nondim]}
1626   \textcolor{keywordtype}{real} :: rho\_ave            \textcolor{comment}{! Average density [R ~> kg m-3]}
1627   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(5)} :: T5   \textcolor{comment}{! Temperatures at quadrature points [degC]}
1628   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(5)} :: S5   \textcolor{comment}{! Salinities at quadrature points [ppt]}
1629   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(5)} :: p5   \textcolor{comment}{! Pressures at quadrature points [R L2 T-2 ~> Pa]}
1630   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(5)} :: rho5 \textcolor{comment}{! Densities at quadrature points [R ~> kg m-3]}
1631   \textcolor{keywordtype}{integer} :: n
1632 
1633   \textcolor{keywordflow}{do} n=1,5
1634     \textcolor{comment}{! Evaluate density at five quadrature points}
1635     bottom\_weight = 0.25*\textcolor{keywordtype}{real(n-1)} * pos
1636     top\_weight = 1.0 - bottom\_weight
1637     \textcolor{comment}{! Salinity and temperature points are linearly interpolated}
1638     s5(n) = top\_weight * s\_t + bottom\_weight * s\_b
1639     t5(n) = top\_weight * t\_t + bottom\_weight * t\_b
1640     p5(n) = ( top\_weight * z\_t + bottom\_weight * z\_b ) * ( g\_e * rho\_ref )
1641 \textcolor{keywordflow}{  enddo}
1642   \textcolor{keyword}{call }calculate\_density(t5, s5, p5, rho5, eos)
1643   rho5(:) = rho5(:) \textcolor{comment}{!- rho\_ref ! Work with anomalies relative to rho\_ref}
1644 
1645   \textcolor{comment}{! Use Boole's rule to estimate the average density}
1646   rho\_ave = c1\_90*(7.0*(rho5(1)+rho5(5)) + 32.0*(rho5(2)+rho5(4)) + 12.0*rho5(3))
1647 
1648   dz = ( z\_t - z\_b ) * pos
1649   frac\_dp\_at\_pos = g\_e * dz * rho\_ave
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_ac36ae5f4af2d02df0a1adf41b762e017}\label{namespacemom__density__integrals_ac36ae5f4af2d02df0a1adf41b762e017}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!int\+\_\+density\+\_\+dz@{int\+\_\+density\+\_\+dz}}
\index{int\+\_\+density\+\_\+dz@{int\+\_\+density\+\_\+dz}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{int\+\_\+density\+\_\+dz()}{int\_density\_dz()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::int\+\_\+density\+\_\+dz (\begin{DoxyParamCaption}\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{T,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{S,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{z\+\_\+t,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{z\+\_\+b,  }\item[{real, intent(in)}]{rho\+\_\+ref,  }\item[{real, intent(in)}]{rho\+\_\+0,  }\item[{real, intent(in)}]{G\+\_\+e,  }\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(inout)}]{dpa,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(inout), optional}]{intz\+\_\+dpa,  }\item[{real, dimension(szib\+\_\+(hi),szj\+\_\+(hi)), intent(inout), optional}]{intx\+\_\+dpa,  }\item[{real, dimension(szi\+\_\+(hi),szjb\+\_\+(hi)), intent(inout), optional}]{inty\+\_\+dpa,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in), optional}]{bathyT,  }\item[{real, intent(in), optional}]{dz\+\_\+neglect,  }\item[{logical, intent(in), optional}]{use\+Mass\+Wght\+Interp }\end{DoxyParamCaption})}



Calls the appropriate subroutine to calculate analytical and nearly-\/analytical integrals in z across layers of pressure anomalies, which are required for calculating the finite-\/volume form pressure accelerations in a Boussinesq model. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em hi} & Ocean horizontal index structures for the arrays\\
\hline
\mbox{\tt in}  & {\em t} & Potential temperature referenced to the surface \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s} & Salinity \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+t} & Height at the top of the layer in depth units \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+b} & Height at the bottom of the layer \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+ref} & A mean density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is subtracted out to reduce the magnitude of each of the integrals.\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+0} & A density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is used to calculate the pressure (as p$\sim$=-\/z$\ast$rho\+\_\+0$\ast$\+G\+\_\+e) used in the equation of state.\\
\hline
\mbox{\tt in}  & {\em g\+\_\+e} & The Earth\textquotesingle{}s gravitational acceleration \mbox{[}L2 Z-\/1 T-\/2 $\sim$$>$ m s-\/2\mbox{]} or \mbox{[}m2 Z-\/1 s-\/2 $\sim$$>$ m s-\/2\mbox{]}\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in,out}  & {\em dpa} & The change in the pressure anomaly\\
\hline
\mbox{\tt in,out}  & {\em intz\+\_\+dpa} & The integral through the thickness of the\\
\hline
\mbox{\tt in,out}  & {\em intx\+\_\+dpa} & The integral in x of the difference between\\
\hline
\mbox{\tt in,out}  & {\em inty\+\_\+dpa} & The integral in y of the difference between\\
\hline
\mbox{\tt in}  & {\em bathyt} & The depth of the bathymetry \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em dz\+\_\+neglect} & A minuscule thickness change \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em usemasswghtinterp} & If true, uses mass weighting to interpolate T/S for top and bottom integrals. \\
\hline
\end{DoxyParams}


Definition at line 42 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
42   \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)}  :: HI\textcolor{comment}{  !< Ocean horizontal index structures for the arrays}
43   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
44                         \textcolor{keywordtype}{intent(in)}  :: T\textcolor{comment}{   !< Potential temperature referenced to the surface [degC]}
45   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
46                         \textcolor{keywordtype}{intent(in)}  :: S\textcolor{comment}{   !< Salinity [ppt]}
47   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
48                         \textcolor{keywordtype}{intent(in)}  :: z\_t\textcolor{comment}{ !< Height at the top of the layer in depth units [Z ~> m]}
49   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
50                         \textcolor{keywordtype}{intent(in)}  :: z\_b\textcolor{comment}{ !< Height at the bottom of the layer [Z ~> m]}
51   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_ref\textcolor{comment}{ !< A mean density [R ~> kg m-3] or [kg m-3], that is}
52 \textcolor{comment}{                                           !! subtracted out to reduce the magnitude of each of the}
53 \textcolor{comment}{                                           !! integrals.}
54   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_0\textcolor{comment}{ !< A density [R ~> kg m-3] or [kg m-3], that is used}
55 \textcolor{comment}{                                           !! to calculate the pressure (as p~=-z*rho\_0*G\_e)}
56 \textcolor{comment}{                                           !! used in the equation of state.}
57   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: G\_e\textcolor{comment}{ !< The Earth's gravitational acceleration}
58 \textcolor{comment}{                                           !! [L2 Z-1 T-2 ~> m s-2] or [m2 Z-1 s-2 ~> m s-2]}
59   \textcolor{keywordtype}{type}(EOS\_type),       \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
60   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{  !< A dimensional unit scaling type}
61   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
62                       \textcolor{keywordtype}{intent(inout)} :: dpa\textcolor{comment}{ !< The change in the pressure anomaly}
63 \textcolor{comment}{                                           !! across the layer [R L2 T-2 ~> Pa] or [Pa]}
64   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
65             \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intz\_dpa\textcolor{comment}{ !< The integral through the thickness of the}
66 \textcolor{comment}{                                           !! layer of the pressure anomaly relative to the}
67 \textcolor{comment}{                                           !! anomaly at the top of the layer [R L2 Z T-2 ~> Pa m]}
68   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(HI),SZJ\_(HI))}, &
69             \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intx\_dpa\textcolor{comment}{ !< The integral in x of the difference between}
70 \textcolor{comment}{                                          !! the pressure anomaly at the top and bottom of the}
71 \textcolor{comment}{                                          !! layer divided by the x grid spacing [R L2 T-2 ~> Pa]}
72   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJB\_(HI))}, &
73             \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: inty\_dpa\textcolor{comment}{ !< The integral in y of the difference between}
74 \textcolor{comment}{                                          !! the pressure anomaly at the top and bottom of the}
75 \textcolor{comment}{                                          !! layer divided by the y grid spacing [R L2 T-2 ~> Pa]}
76   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
77               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: bathyT\textcolor{comment}{ !< The depth of the bathymetry [Z ~> m]}
78   \textcolor{keywordtype}{real},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: dz\_neglect\textcolor{comment}{ !< A minuscule thickness change [Z ~> m]}
79   \textcolor{keywordtype}{logical},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: useMassWghtInterp\textcolor{comment}{ !< If true, uses mass weighting to}
80 \textcolor{comment}{                                           !! interpolate T/S for top and bottom integrals.}
81 
82   \textcolor{keywordflow}{if} (eos\_quadrature(eos)) \textcolor{keywordflow}{then}
83     \textcolor{keyword}{call }int\_density\_dz\_generic\_pcm(t, s, z\_t, z\_b, rho\_ref, rho\_0, g\_e, hi, eos, us, dpa, &
84                                     intz\_dpa, intx\_dpa, inty\_dpa, bathyt, dz\_neglect, usemasswghtinterp)
85   \textcolor{keywordflow}{else}
86     \textcolor{keyword}{call }analytic\_int\_density\_dz(t, s, z\_t, z\_b, rho\_ref, rho\_0, g\_e, hi, eos, dpa, &
87                                  intz\_dpa, intx\_dpa, inty\_dpa, bathyt, dz\_neglect, usemasswghtinterp)
88 \textcolor{keywordflow}{  endif}
89 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_ad0ff54518bfacf7beebbe4dba687a914}\label{namespacemom__density__integrals_ad0ff54518bfacf7beebbe4dba687a914}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+pcm@{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+pcm}}
\index{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+pcm@{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+pcm}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+pcm()}{int\_density\_dz\_generic\_pcm()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+pcm (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{T,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{S,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{z\+\_\+t,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{z\+\_\+b,  }\item[{real, intent(in)}]{rho\+\_\+ref,  }\item[{real, intent(in)}]{rho\+\_\+0,  }\item[{real, intent(in)}]{G\+\_\+e,  }\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout)}]{dpa,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intz\+\_\+dpa,  }\item[{real, dimension( hi \%isdb\+: hi \%iedb, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intx\+\_\+dpa,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsdb\+: hi \%jedb), intent(inout), optional}]{inty\+\_\+dpa,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in), optional}]{bathyT,  }\item[{real, intent(in), optional}]{dz\+\_\+neglect,  }\item[{logical, intent(in), optional}]{use\+Mass\+Wght\+Interp }\end{DoxyParamCaption})}



Calculates (by numerical quadrature) integrals of pressure anomalies across layers, which are required for calculating the finite-\/volume form pressure accelerations in a Boussinesq model. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em hi} & Horizontal index type for input variables.\\
\hline
\mbox{\tt in}  & {\em t} & Potential temperature of the layer \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s} & Salinity of the layer \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+t} & Height at the top of the layer in depth units \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em z\+\_\+b} & Height at the bottom of the layer \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+ref} & A mean density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is subtracted out to reduce the magnitude of each of the integrals.\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+0} & A density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is used to calculate the pressure (as p$\sim$=-\/z$\ast$rho\+\_\+0$\ast$\+G\+\_\+e) used in the equation of state.\\
\hline
\mbox{\tt in}  & {\em g\+\_\+e} & The Earth\textquotesingle{}s gravitational acceleration \mbox{[}L2 Z-\/1 T-\/2 $\sim$$>$ m s-\/2\mbox{]} or \mbox{[}m2 Z-\/1 s-\/2 $\sim$$>$ m s-\/2\mbox{]}\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in,out}  & {\em dpa} & The change in the pressure anomaly\\
\hline
\mbox{\tt in,out}  & {\em intz\+\_\+dpa} & The integral through the thickness of the\\
\hline
\mbox{\tt in,out}  & {\em intx\+\_\+dpa} & The integral in x of the difference between\\
\hline
\mbox{\tt in,out}  & {\em inty\+\_\+dpa} & The integral in y of the difference between\\
\hline
\mbox{\tt in}  & {\em bathyt} & The depth of the bathymetry \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em dz\+\_\+neglect} & A minuscule thickness change \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em usemasswghtinterp} & If true, uses mass weighting to interpolate T/S for top and bottom integrals. \\
\hline
\end{DoxyParams}


Definition at line 98 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
98   \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)}  :: HI\textcolor{comment}{  !< Horizontal index type for input variables.}
99   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
100                         \textcolor{keywordtype}{intent(in)}  :: T\textcolor{comment}{  !< Potential temperature of the layer [degC]}
101   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
102                         \textcolor{keywordtype}{intent(in)}  :: S\textcolor{comment}{  !< Salinity of the layer [ppt]}
103   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
104                         \textcolor{keywordtype}{intent(in)}  :: z\_t\textcolor{comment}{ !< Height at the top of the layer in depth units [Z ~> m]}
105   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
106                         \textcolor{keywordtype}{intent(in)}  :: z\_b\textcolor{comment}{ !< Height at the bottom of the layer [Z ~> m]}
107   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_ref\textcolor{comment}{ !< A mean density [R ~> kg m-3] or [kg m-3], that is}
108 \textcolor{comment}{                                          !! subtracted out to reduce the magnitude}
109 \textcolor{comment}{                                          !! of each of the integrals.}
110   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_0\textcolor{comment}{ !< A density [R ~> kg m-3] or [kg m-3], that is used}
111 \textcolor{comment}{                                          !! to calculate the pressure (as p~=-z*rho\_0*G\_e)}
112 \textcolor{comment}{                                          !! used in the equation of state.}
113   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: G\_e\textcolor{comment}{ !< The Earth's gravitational acceleration}
114 \textcolor{comment}{                                          !! [L2 Z-1 T-2 ~> m s-2] or [m2 Z-1 s-2 ~> m s-2]}
115   \textcolor{keywordtype}{type}(EOS\_type),       \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
116   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type}
117   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
118                       \textcolor{keywordtype}{intent(inout)} :: dpa\textcolor{comment}{ !< The change in the pressure anomaly}
119 \textcolor{comment}{                                          !! across the layer [R L2 T-2 ~> Pa] or [Pa]}
120   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
121             \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intz\_dpa\textcolor{comment}{ !< The integral through the thickness of the}
122 \textcolor{comment}{                                          !! layer of the pressure anomaly relative to the}
123 \textcolor{comment}{                                          !! anomaly at the top of the layer [R L2 Z T-2 ~> Pa m]}
124   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(HI),SZJ\_(HI))}, &
125             \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intx\_dpa\textcolor{comment}{ !< The integral in x of the difference between}
126 \textcolor{comment}{                                          !! the pressure anomaly at the top and bottom of the}
127 \textcolor{comment}{                                          !! layer divided by the x grid spacing [R L2 T-2 ~> Pa]}
128   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJB\_(HI))}, &
129             \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: inty\_dpa\textcolor{comment}{ !< The integral in y of the difference between}
130 \textcolor{comment}{                                          !! the pressure anomaly at the top and bottom of the}
131 \textcolor{comment}{                                          !! layer divided by the y grid spacing [R L2 T-2 ~> Pa]}
132   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
133               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: bathyT\textcolor{comment}{ !< The depth of the bathymetry [Z ~> m]}
134   \textcolor{keywordtype}{real},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: dz\_neglect\textcolor{comment}{ !< A minuscule thickness change [Z ~> m]}
135   \textcolor{keywordtype}{logical},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: useMassWghtInterp\textcolor{comment}{ !< If true, uses mass weighting to}
136 \textcolor{comment}{                                          !! interpolate T/S for top and bottom integrals.}
137   \textcolor{comment}{! Local variables}
138   \textcolor{keywordtype}{real} :: T5(5), S5(5) \textcolor{comment}{! Temperatures and salinities at five quadrature points [degC] and [ppt]}
139   \textcolor{keywordtype}{real} :: p5(5)      \textcolor{comment}{! Pressures at five quadrature points, never rescaled from Pa [Pa]}
140   \textcolor{keywordtype}{real} :: r5(5)      \textcolor{comment}{! Densities at five quadrature points [R ~> kg m-3] or [kg m-3]}
141   \textcolor{keywordtype}{real} :: rho\_anom   \textcolor{comment}{! The depth averaged density anomaly [R ~> kg m-3] or [kg m-3]}
142   \textcolor{keywordtype}{real} :: w\_left, w\_right \textcolor{comment}{! Left and right weights [nondim]}
143   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_90 = 1.0/90.0  \textcolor{comment}{! Rational constants.}
144   \textcolor{keywordtype}{real} :: GxRho      \textcolor{comment}{! The gravitational acceleration times density and unit conversion factors [Pa Z-1 ~>
       kg m-2 s-2]}
145   \textcolor{keywordtype}{real} :: I\_Rho      \textcolor{comment}{! The inverse of the Boussinesq density [R-1 ~> m3 kg-1] or [m3 kg-1]}
146   \textcolor{keywordtype}{real} :: rho\_scale  \textcolor{comment}{! A scaling factor for densities from kg m-3 to R [R m3 kg-1 ~> 1]}
147   \textcolor{keywordtype}{real} :: rho\_ref\_mks \textcolor{comment}{! The reference density in MKS units, never rescaled from kg m-3 [kg m-3]}
148   \textcolor{keywordtype}{real} :: dz         \textcolor{comment}{! The layer thickness [Z ~> m]}
149   \textcolor{keywordtype}{real} :: hWght      \textcolor{comment}{! A pressure-thickness below topography [Z ~> m]}
150   \textcolor{keywordtype}{real} :: hL, hR     \textcolor{comment}{! Pressure-thicknesses of the columns to the left and right [Z ~> m]}
151   \textcolor{keywordtype}{real} :: iDenom     \textcolor{comment}{! The inverse of the denominator in the weights [Z-2 ~> m-2]}
152   \textcolor{keywordtype}{real} :: hWt\_LL, hWt\_LR \textcolor{comment}{! hWt\_LA is the weighted influence of A on the left column [nondim]}
153   \textcolor{keywordtype}{real} :: hWt\_RL, hWt\_RR \textcolor{comment}{! hWt\_RA is the weighted influence of A on the right column [nondim]}
154   \textcolor{keywordtype}{real} :: wt\_L, wt\_R \textcolor{comment}{! The linear weights of the left and right columns [nondim]}
155   \textcolor{keywordtype}{real} :: wtT\_L, wtT\_R \textcolor{comment}{! The weights for tracers from the left and right columns [nondim]}
156   \textcolor{keywordtype}{real} :: intz(5)    \textcolor{comment}{! The gravitational acceleration times the integrals of density}
157                      \textcolor{comment}{! with height at the 5 sub-column locations [R L2 T-2 ~> Pa] or [Pa]}
158   \textcolor{keywordtype}{logical} :: do\_massWeight \textcolor{comment}{! Indicates whether to do mass weighting.}
159   \textcolor{keywordtype}{integer} :: is, ie, js, je, Isq, Ieq, Jsq, Jeq, i, j, m, n
160 
161   \textcolor{comment}{! These array bounds work for the indexing convention of the input arrays, but}
162   \textcolor{comment}{! on the computational domain defined for the output arrays.}
163   isq = hi%IscB ; ieq = hi%IecB
164   jsq = hi%JscB ; jeq = hi%JecB
165   is = hi%isc ; ie = hi%iec
166   js = hi%jsc ; je = hi%jec
167 
168   rho\_scale = us%kg\_m3\_to\_R
169   gxrho = us%RL2\_T2\_to\_Pa * g\_e * rho\_0
170   rho\_ref\_mks = rho\_ref * us%R\_to\_kg\_m3
171   i\_rho = 1.0 / rho\_0
172 
173   do\_massweight = .false.
174   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(usemasswghtinterp)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (usemasswghtinterp) \textcolor{keywordflow}{then}
175     do\_massweight = .true.
176     \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{present}(bathyt)) \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"int\_density\_dz\_generic: "}//&
177         \textcolor{stringliteral}{"bathyT must be present if useMassWghtInterp is present and true."})
178     \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{present}(dz\_neglect)) \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"int\_density\_dz\_generic: "}//&
179         \textcolor{stringliteral}{"dz\_neglect must be present if useMassWghtInterp is present and true."})
180 \textcolor{keywordflow}{  endif} ;\textcolor{keywordflow}{ endif}
181 
182   \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1
183     dz = z\_t(i,j) - z\_b(i,j)
184     \textcolor{keywordflow}{do} n=1,5
185       t5(n) = t(i,j) ; s5(n) = s(i,j)
186       p5(n) = -gxrho*(z\_t(i,j) - 0.25*\textcolor{keywordtype}{real}(n-1)*dz)
187 \textcolor{keywordflow}{    enddo}
188     \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
189       \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
190     \textcolor{keywordflow}{else}
191       \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks)
192 \textcolor{keywordflow}{    endif}
193 
194     \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
195     rho\_anom = c1\_90*(7.0*(r5(1)+r5(5)) + 32.0*(r5(2)+r5(4)) + 12.0*r5(3))
196     dpa(i,j) = g\_e*dz*rho\_anom
197     \textcolor{comment}{! Use a Boole's-rule-like fifth-order accurate estimate of the double integral of}
198     \textcolor{comment}{! the pressure anomaly.}
199     \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intz\_dpa)) intz\_dpa(i,j) = 0.5*g\_e*dz**2 * &
200           (rho\_anom - c1\_90*(16.0*(r5(4)-r5(2)) + 7.0*(r5(5)-r5(1))) )
201 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
202 
203   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intx\_dpa)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq
204     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
205     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation of}
206     \textcolor{comment}{! T & S along the top and bottom integrals, akin to thickness weighting.}
207     hwght = 0.0
208     \textcolor{keywordflow}{if} (do\_massweight) &
209       hwght = max(0., -bathyt(i,j)-z\_t(i+1,j), -bathyt(i+1,j)-z\_t(i,j))
210     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
211       hl = (z\_t(i,j) - z\_b(i,j)) + dz\_neglect
212       hr = (z\_t(i+1,j) - z\_b(i+1,j)) + dz\_neglect
213       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
214       idenom = 1.0 / ( hwght*(hr + hl) + hl*hr )
215       hwt\_ll = (hwght*hl + hr*hl) * idenom ; hwt\_lr = (hwght*hr) * idenom
216       hwt\_rr = (hwght*hr + hr*hl) * idenom ; hwt\_rl = (hwght*hl) * idenom
217     \textcolor{keywordflow}{else}
218       hwt\_ll = 1.0 ; hwt\_lr = 0.0 ; hwt\_rr = 1.0 ; hwt\_rl = 0.0
219 \textcolor{keywordflow}{    endif}
220 
221     intz(1) = dpa(i,j) ; intz(5) = dpa(i+1,j)
222     \textcolor{keywordflow}{do} m=2,4
223       \textcolor{comment}{! T, S, and z are interpolated in the horizontal.  The z interpolation}
224       \textcolor{comment}{! is linear, but for T and S it may be thickness weighted.}
225       wt\_l = 0.25*\textcolor{keywordtype}{real(5-m)} ; wt\_R = 1.0-wt\_l
226       wtt\_l = wt\_l*hwt\_ll + wt\_r*hwt\_rl ; wtt\_r = wt\_l*hwt\_lr + wt\_r*hwt\_rr
227       dz = wt\_l*(z\_t(i,j) - z\_b(i,j)) + wt\_r*(z\_t(i+1,j) - z\_b(i+1,j))
228       t5(1) = wtt\_l*t(i,j) + wtt\_r*t(i+1,j)
229       s5(1) = wtt\_l*s(i,j) + wtt\_r*s(i+1,j)
230       p5(1) = -gxrho*(wt\_l*z\_t(i,j) + wt\_r*z\_t(i+1,j))
231       \textcolor{keywordflow}{do} n=2,5
232         t5(n) = t5(1) ; s5(n) = s5(1) ; p5(n) = p5(n-1) + gxrho*0.25*dz
233 \textcolor{keywordflow}{      enddo}
234       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
235         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
236       \textcolor{keywordflow}{else}
237         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks)
238 \textcolor{keywordflow}{      endif}
239 
240     \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
241       intz(m) = g\_e*dz*( c1\_90*(7.0*(r5(1)+r5(5)) + 32.0*(r5(2)+r5(4)) + 12.0*r5(3)))
242 \textcolor{keywordflow}{    enddo}
243     \textcolor{comment}{! Use Boole's rule to integrate the bottom pressure anomaly values in x.}
244     intx\_dpa(i,j) = c1\_90*(7.0*(intz(1)+intz(5)) + 32.0*(intz(2)+intz(4)) + &
245                            12.0*intz(3))
246 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
247 
248   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(inty\_dpa)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie
249     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
250     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation of}
251     \textcolor{comment}{! T & S along the top and bottom integrals, akin to thickness weighting.}
252     hwght = 0.0
253     \textcolor{keywordflow}{if} (do\_massweight) &
254       hwght = max(0., -bathyt(i,j)-z\_t(i,j+1), -bathyt(i,j+1)-z\_t(i,j))
255     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
256       hl = (z\_t(i,j) - z\_b(i,j)) + dz\_neglect
257       hr = (z\_t(i,j+1) - z\_b(i,j+1)) + dz\_neglect
258       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
259       idenom = 1.0 / ( hwght*(hr + hl) + hl*hr )
260       hwt\_ll = (hwght*hl + hr*hl) * idenom ; hwt\_lr = (hwght*hr) * idenom
261       hwt\_rr = (hwght*hr + hr*hl) * idenom ; hwt\_rl = (hwght*hl) * idenom
262     \textcolor{keywordflow}{else}
263       hwt\_ll = 1.0 ; hwt\_lr = 0.0 ; hwt\_rr = 1.0 ; hwt\_rl = 0.0
264 \textcolor{keywordflow}{    endif}
265 
266     intz(1) = dpa(i,j) ; intz(5) = dpa(i,j+1)
267     \textcolor{keywordflow}{do} m=2,4
268       \textcolor{comment}{! T, S, and z are interpolated in the horizontal.  The z interpolation}
269       \textcolor{comment}{! is linear, but for T and S it may be thickness weighted.}
270       wt\_l = 0.25*\textcolor{keywordtype}{real(5-m)} ; wt\_R = 1.0-wt\_l
271       wtt\_l = wt\_l*hwt\_ll + wt\_r*hwt\_rl ; wtt\_r = wt\_l*hwt\_lr + wt\_r*hwt\_rr
272       dz = wt\_l*(z\_t(i,j) - z\_b(i,j)) + wt\_r*(z\_t(i,j+1) - z\_b(i,j+1))
273       t5(1) = wtt\_l*t(i,j) + wtt\_r*t(i,j+1)
274       s5(1) = wtt\_l*s(i,j) + wtt\_r*s(i,j+1)
275       p5(1) = -gxrho*(wt\_l*z\_t(i,j) + wt\_r*z\_t(i,j+1))
276       \textcolor{keywordflow}{do} n=2,5
277         t5(n) = t5(1) ; s5(n) = s5(1)
278         p5(n) = p5(n-1) + gxrho*0.25*dz
279 \textcolor{keywordflow}{      enddo}
280       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
281         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
282       \textcolor{keywordflow}{else}
283         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks)
284 \textcolor{keywordflow}{      endif}
285 
286     \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
287       intz(m) = g\_e*dz*( c1\_90*(7.0*(r5(1)+r5(5)) + 32.0*(r5(2)+r5(4)) + 12.0*r5(3)))
288 \textcolor{keywordflow}{    enddo}
289     \textcolor{comment}{! Use Boole's rule to integrate the values.}
290     inty\_dpa(i,j) = c1\_90*(7.0*(intz(1)+intz(5)) + 32.0*(intz(2)+intz(4)) + &
291                                      12.0*intz(3))
292 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_a454d4d62eb599716cc3389fb8aa90b4b}\label{namespacemom__density__integrals_a454d4d62eb599716cc3389fb8aa90b4b}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+plm@{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+plm}}
\index{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+plm@{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+plm}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+plm()}{int\_density\_dz\_generic\_plm()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+plm (\begin{DoxyParamCaption}\item[{integer, intent(in)}]{k,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed, gv \%ke), intent(in)}]{T\+\_\+t,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed, gv \%ke), intent(in)}]{T\+\_\+b,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed, gv \%ke), intent(in)}]{S\+\_\+t,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed, gv \%ke), intent(in)}]{S\+\_\+b,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed, gv \%ke+1), intent(in)}]{e,  }\item[{real, intent(in)}]{rho\+\_\+ref,  }\item[{real, intent(in)}]{rho\+\_\+0,  }\item[{real, intent(in)}]{G\+\_\+e,  }\item[{real, intent(in)}]{dz\+\_\+subroundoff,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{bathyT,  }\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(inout)}]{dpa,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(inout), optional}]{intz\+\_\+dpa,  }\item[{real, dimension(szib\+\_\+(hi),szj\+\_\+(hi)), intent(inout), optional}]{intx\+\_\+dpa,  }\item[{real, dimension(szi\+\_\+(hi),szjb\+\_\+(hi)), intent(inout), optional}]{inty\+\_\+dpa,  }\item[{logical, intent(in), optional}]{use\+Mass\+Wght\+Interp }\end{DoxyParamCaption})}



Compute pressure gradient force integrals by quadrature for the case where T and S are linear profiles. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em k} & Layer index to calculate integrals for\\
\hline
\mbox{\tt in}  & {\em hi} & Ocean horizontal index structures for the input arrays\\
\hline
\mbox{\tt in}  & {\em gv} & Vertical grid structure\\
\hline
\mbox{\tt in}  & {\em tv} & Thermodynamic variables\\
\hline
\mbox{\tt in}  & {\em t\+\_\+t} & Potential temperature at the cell top \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em t\+\_\+b} & Potential temperature at the cell bottom \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+t} & Salinity at the cell top \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+b} & Salinity at the cell bottom \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em e} & Height of interfaces \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+ref} & A mean density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is subtracted out to reduce the magnitude of each of the integrals.\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+0} & A density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is used to calculate the pressure (as p$\sim$=-\/z$\ast$rho\+\_\+0$\ast$\+G\+\_\+e) used in the equation of state.\\
\hline
\mbox{\tt in}  & {\em g\+\_\+e} & The Earth\textquotesingle{}s gravitational acceleration \mbox{[}L2 Z-\/1 T-\/2 $\sim$$>$ m s-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em dz\+\_\+subroundoff} & A minuscule thickness change \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em bathyt} & The depth of the bathymetry \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in,out}  & {\em dpa} & The change in the pressure anomaly across the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}\\
\hline
\mbox{\tt in,out}  & {\em intz\+\_\+dpa} & The integral through the thickness of the layer of\\
\hline
\mbox{\tt in,out}  & {\em intx\+\_\+dpa} & The integral in x of the difference between the\\
\hline
\mbox{\tt in,out}  & {\em inty\+\_\+dpa} & The integral in y of the difference between the\\
\hline
\mbox{\tt in}  & {\em usemasswghtinterp} & If true, uses mass weighting to interpolate T/S for top and bottom integrals. \\
\hline
\end{DoxyParams}


Definition at line 301 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
301   \textcolor{keywordtype}{integer},              \textcolor{keywordtype}{intent(in)}  :: k\textcolor{comment}{   !< Layer index to calculate integrals for}
302   \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)}  :: HI\textcolor{comment}{  !< Ocean horizontal index structures for the input arrays}
303   \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< Vertical grid structure}
304   \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{  !< Thermodynamic variables}
305   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
306                         \textcolor{keywordtype}{intent(in)}  :: T\_t\textcolor{comment}{ !< Potential temperature at the cell top [degC]}
307   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
308                         \textcolor{keywordtype}{intent(in)}  :: T\_b\textcolor{comment}{ !< Potential temperature at the cell bottom [degC]}
309   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
310                         \textcolor{keywordtype}{intent(in)}  :: S\_t\textcolor{comment}{ !< Salinity at the cell top [ppt]}
311   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
312                         \textcolor{keywordtype}{intent(in)}  :: S\_b\textcolor{comment}{ !< Salinity at the cell bottom [ppt]}
313   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV)+1)}, &
314                         \textcolor{keywordtype}{intent(in)}  :: e\textcolor{comment}{   !< Height of interfaces [Z ~> m]}
315   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_ref\textcolor{comment}{ !< A mean density [R ~> kg m-3] or [kg m-3], that is
       subtracted}
316 \textcolor{comment}{                                           !! out to reduce the magnitude of each of the integrals.}
317   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_0\textcolor{comment}{ !< A density [R ~> kg m-3] or [kg m-3], that is used to
       calculate}
318 \textcolor{comment}{                                           !! the pressure (as p~=-z*rho\_0*G\_e) used in the equation of
       state.}
319   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: G\_e\textcolor{comment}{ !< The Earth's gravitational acceleration [L2 Z-1 T-2 ~> m s-2]}
320   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: dz\_subroundoff\textcolor{comment}{ !< A minuscule thickness change [Z ~> m]}
321   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
322                         \textcolor{keywordtype}{intent(in)}  :: bathyT\textcolor{comment}{ !< The depth of the bathymetry [Z ~> m]}
323   \textcolor{keywordtype}{type}(EOS\_type),       \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
324   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type}
325   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
326                         \textcolor{keywordtype}{intent(inout)} :: dpa\textcolor{comment}{ !< The change in the pressure anomaly across the layer [R L2
       T-2 ~> Pa]}
327   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
328               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intz\_dpa\textcolor{comment}{ !< The integral through the thickness of the layer of}
329 \textcolor{comment}{                                           !! the pressure anomaly relative to the anomaly at the}
330 \textcolor{comment}{                                           !! top of the layer [R L2 Z T-2 ~> Pa Z]}
331   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(HI),SZJ\_(HI))}, &
332               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intx\_dpa\textcolor{comment}{ !< The integral in x of the difference between the}
333 \textcolor{comment}{                                           !! pressure anomaly at the top and bottom of the layer}
334 \textcolor{comment}{                                           !! divided by the x grid spacing [R L2 T-2 ~> Pa]}
335   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJB\_(HI))}, &
336               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: inty\_dpa\textcolor{comment}{ !< The integral in y of the difference between the}
337 \textcolor{comment}{                                           !! pressure anomaly at the top and bottom of the layer}
338 \textcolor{comment}{                                           !! divided by the y grid spacing [R L2 T-2 ~> Pa]}
339   \textcolor{keywordtype}{logical},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: useMassWghtInterp\textcolor{comment}{ !< If true, uses mass weighting to}
340 \textcolor{comment}{                                           !! interpolate T/S for top and bottom integrals.}
341 
342 \textcolor{comment}{! This subroutine calculates (by numerical quadrature) integrals of}
343 \textcolor{comment}{! pressure anomalies across layers, which are required for calculating the}
344 \textcolor{comment}{! finite-volume form pressure accelerations in a Boussinesq model.  The one}
345 \textcolor{comment}{! potentially dodgy assumption here is that rho\_0 is used both in the denominator}
346 \textcolor{comment}{! of the accelerations, and in the pressure used to calculated density (the}
347 \textcolor{comment}{! latter being -z*rho\_0*G\_e).  These two uses could be separated if need be.}
348 \textcolor{comment}{!}
349 \textcolor{comment}{! It is assumed that the salinity and temperature profiles are linear in the}
350 \textcolor{comment}{! vertical. The top and bottom values within each layer are provided and}
351 \textcolor{comment}{! a linear interpolation is used to compute intermediate values.}
352 
353   \textcolor{comment}{! Local variables}
354   \textcolor{keywordtype}{real} :: T5((5*HI%iscB+1):(5*(HI%iecB+2)))  \textcolor{comment}{! Temperatures along a line of subgrid locations [degC]}
355   \textcolor{keywordtype}{real} :: S5((5*HI%iscB+1):(5*(HI%iecB+2)))  \textcolor{comment}{! Salinities along a line of subgrid locations [ppt]}
356   \textcolor{keywordtype}{real} :: T25((5*HI%iscB+1):(5*(HI%iecB+2))) \textcolor{comment}{! SGS temperature variance along a line of subgrid locations
       [degC2]}
357   \textcolor{keywordtype}{real} :: TS5((5*HI%iscB+1):(5*(HI%iecB+2))) \textcolor{comment}{! SGS temp-salt covariance along a line of subgrid locations
       [degC ppt]}
358   \textcolor{keywordtype}{real} :: S25((5*HI%iscB+1):(5*(HI%iecB+2))) \textcolor{comment}{! SGS salinity variance along a line of subgrid locations
       [ppt2]}
359   \textcolor{keywordtype}{real} :: p5((5*HI%iscB+1):(5*(HI%iecB+2)))  \textcolor{comment}{! Pressures along a line of subgrid locations, never}
360                                                \textcolor{comment}{! rescaled from Pa [Pa]}
361   \textcolor{keywordtype}{real} :: r5((5*HI%iscB+1):(5*(HI%iecB+2)))  \textcolor{comment}{! Densities anomalies along a line of subgrid}
362                                                \textcolor{comment}{! locations [R ~> kg m-3] or [kg m-3]}
363   \textcolor{keywordtype}{real} :: T15((15*HI%iscB+1):(15*(HI%iecB+1))) \textcolor{comment}{! Temperatures at an array of subgrid locations [degC]}
364   \textcolor{keywordtype}{real} :: S15((15*HI%iscB+1):(15*(HI%iecB+1))) \textcolor{comment}{! Salinities at an array of subgrid locations [ppt]}
365   \textcolor{keywordtype}{real} :: T215((15*HI%iscB+1):(15*(HI%iecB+1))) \textcolor{comment}{! SGS temperature variance along a line of subgrid
       locations [degC2]}
366   \textcolor{keywordtype}{real} :: TS15((15*HI%iscB+1):(15*(HI%iecB+1))) \textcolor{comment}{! SGS temp-salt covariance along a line of subgrid
       locations [degC ppt]}
367   \textcolor{keywordtype}{real} :: S215((15*HI%iscB+1):(15*(HI%iecB+1))) \textcolor{comment}{! SGS salinity variance along a line of subgrid locations
       [ppt2]}
368   \textcolor{keywordtype}{real} :: p15((15*HI%iscB+1):(15*(HI%iecB+1))) \textcolor{comment}{! Pressures at an array of subgrid locations [Pa]}
369   \textcolor{keywordtype}{real} :: r15((15*HI%iscB+1):(15*(HI%iecB+1))) \textcolor{comment}{! Densities at an array of subgrid locations}
370                                                  \textcolor{comment}{! [R ~> kg m-3] or [kg m-3]}
371   \textcolor{keywordtype}{real} :: wt\_t(5), wt\_b(5)          \textcolor{comment}{! Top and bottom weights [nondim]}
372   \textcolor{keywordtype}{real} :: rho\_anom                  \textcolor{comment}{! A density anomaly [R ~> kg m-3] or [kg m-3]}
373   \textcolor{keywordtype}{real} :: w\_left, w\_right           \textcolor{comment}{! Left and right weights [nondim]}
374   \textcolor{keywordtype}{real} :: intz(5)    \textcolor{comment}{! The gravitational acceleration times the integrals of density}
375                      \textcolor{comment}{! with height at the 5 sub-column locations [R L2 T-2 ~> Pa] or [Pa]}
376   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_90 = 1.0/90.0  \textcolor{comment}{! A rational constant [nondim]}
377   \textcolor{keywordtype}{real} :: GxRho      \textcolor{comment}{! The gravitational acceleration times density and unit conversion factors [Pa Z-1 ~>
       kg m-2 s-2]}
378   \textcolor{keywordtype}{real} :: I\_Rho      \textcolor{comment}{! The inverse of the Boussinesq density [R-1 ~> m3 kg-1] or [m3 kg-1]}
379   \textcolor{keywordtype}{real} :: rho\_scale  \textcolor{comment}{! A scaling factor for densities from kg m-3 to R [R m3 kg-1 ~> 1]}
380   \textcolor{keywordtype}{real} :: rho\_ref\_mks \textcolor{comment}{! The reference density in MKS units, never rescaled from kg m-3 [kg m-3]}
381   \textcolor{keywordtype}{real} :: dz(HI%iscB:HI%iecB+1)   \textcolor{comment}{! Layer thicknesses at tracer points [Z ~> m]}
382   \textcolor{keywordtype}{real} :: dz\_x(5,HI%iscB:HI%iecB) \textcolor{comment}{! Layer thicknesses along an x-line of subrid locations [Z ~> m]}
383   \textcolor{keywordtype}{real} :: dz\_y(5,HI%isc:HI%iec)   \textcolor{comment}{! Layer thicknesses along a y-line of subrid locations [Z ~> m]}
384   \textcolor{keywordtype}{real} :: massWeightToggle          \textcolor{comment}{! A non-dimensional toggle factor (0 or 1) [nondim]}
385   \textcolor{keywordtype}{real} :: Ttl, Tbl, Ttr, Tbr        \textcolor{comment}{! Temperatures at the velocity cell corners [degC]}
386   \textcolor{keywordtype}{real} :: Stl, Sbl, Str, Sbr        \textcolor{comment}{! Salinities at the velocity cell corners [ppt]}
387   \textcolor{keywordtype}{real} :: hWght                     \textcolor{comment}{! A topographically limited thicknes weight [Z ~> m]}
388   \textcolor{keywordtype}{real} :: hL, hR                    \textcolor{comment}{! Thicknesses to the left and right [Z ~> m]}
389   \textcolor{keywordtype}{real} :: iDenom                    \textcolor{comment}{! The denominator of the thickness weight expressions [Z-2 ~> m-2]}
390   \textcolor{keywordtype}{logical} :: use\_stanley\_eos \textcolor{comment}{! True is SGS variance fields exist in tv.}
391   \textcolor{keywordtype}{logical} :: use\_varT, use\_varS, use\_covarTS
392   \textcolor{keywordtype}{integer} :: Isq, Ieq, Jsq, Jeq, i, j, m, n
393   \textcolor{keywordtype}{integer} :: pos
394 
395   isq = hi%IscB ; ieq = hi%IecB ; jsq = hi%JscB ; jeq = hi%JecB
396 
397   rho\_scale = us%kg\_m3\_to\_R
398   gxrho = us%RL2\_T2\_to\_Pa * g\_e * rho\_0
399   rho\_ref\_mks = rho\_ref * us%R\_to\_kg\_m3
400   i\_rho = 1.0 / rho\_0
401   massweighttoggle = 0.
402   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(usemasswghtinterp)) \textcolor{keywordflow}{then}
403     \textcolor{keywordflow}{if} (usemasswghtinterp) massweighttoggle = 1.
404 \textcolor{keywordflow}{  endif}
405 
406   use\_vart = \textcolor{keyword}{associated}(tv%varT)
407   use\_covarts = \textcolor{keyword}{associated}(tv%covarTS)
408   use\_vars = \textcolor{keyword}{associated}(tv%varS)
409   use\_stanley\_eos = use\_vart .or. use\_covarts .or. use\_vars
410   t25(:) = 0.
411   ts5(:) = 0.
412   s25(:) = 0.
413   t215(:) = 0.
414   ts15(:) = 0.
415   s215(:) = 0.
416 
417   \textcolor{keywordflow}{do} n = 1, 5
418     wt\_t(n) = 0.25 * \textcolor{keywordtype}{real}(5-n)
419     wt\_b(n) = 1.0 - wt\_t(n)
420 \textcolor{keywordflow}{  enddo}
421 
422   \textcolor{comment}{! 1. Compute vertical integrals}
423   \textcolor{keywordflow}{do} j=jsq,jeq+1
424     \textcolor{keywordflow}{do} i = isq,ieq+1
425       dz(i) = e(i,j,k) - e(i,j,k+1)
426       \textcolor{keywordflow}{do} n=1,5
427         p5(i*5+n) = -gxrho*(e(i,j,k) - 0.25*\textcolor{keywordtype}{real}(n-1)*dz(i))
428         \textcolor{comment}{! Salinity and temperature points are linearly interpolated}
429         s5(i*5+n) = wt\_t(n) * s\_t(i,j,k) + wt\_b(n) * s\_b(i,j,k)
430         t5(i*5+n) = wt\_t(n) * t\_t(i,j,k) + wt\_b(n) * t\_b(i,j,k)
431 \textcolor{keywordflow}{      enddo}
432       \textcolor{keywordflow}{if} (use\_vart) t25(i*5+1:i*5+5) = tv%varT(i,j,k)
433       \textcolor{keywordflow}{if} (use\_covarts) ts5(i*5+1:i*5+5) = tv%covarTS(i,j,k)
434       \textcolor{keywordflow}{if} (use\_vars) s25(i*5+1:i*5+5) = tv%varS(i,j,k)
435 \textcolor{keywordflow}{    enddo}
436     \textcolor{keywordflow}{if} (use\_stanley\_eos) \textcolor{keywordflow}{then}
437       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
438         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, t25, ts5, s25, r5, 1, (ieq-isq+2)*5, eos, &
439                                rho\_ref=rho\_ref\_mks, scale=rho\_scale)
440       \textcolor{keywordflow}{else}
441         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, t25, ts5, s25, r5, 1, (ieq-isq+2)*5, eos, &
442                                rho\_ref=rho\_ref\_mks)
443 \textcolor{keywordflow}{      endif}
444     \textcolor{keywordflow}{else}
445       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
446         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, (ieq-isq+2)*5, eos, rho\_ref=rho\_ref\_mks, &
447                                scale=rho\_scale)
448       \textcolor{keywordflow}{else}
449         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, (ieq-isq+2)*5, eos, rho\_ref=rho\_ref\_mks)
450 \textcolor{keywordflow}{      endif}
451 \textcolor{keywordflow}{    endif}
452 
453     \textcolor{keywordflow}{do} i=isq,ieq+1
454     \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
455       rho\_anom = c1\_90*(7.0*(r5(i*5+1)+r5(i*5+5)) + 32.0*(r5(i*5+2)+r5(i*5+4)) + 12.0*r5(i*5+3))
456       dpa(i,j) = g\_e*dz(i)*rho\_anom
457       \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intz\_dpa)) \textcolor{keywordflow}{then}
458       \textcolor{comment}{! Use a Boole's-rule-like fifth-order accurate estimate of}
459       \textcolor{comment}{! the double integral of the pressure anomaly.}
460         intz\_dpa(i,j) = 0.5*g\_e*dz(i)**2 * &
461                 (rho\_anom - c1\_90*(16.0*(r5(i*5+4)-r5(i*5+2)) + 7.0*(r5(i*5+5)-r5(i*5+1))) )
462 \textcolor{keywordflow}{      endif}
463 \textcolor{keywordflow}{    enddo}
464 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! end loops on j}
465 
466   \textcolor{comment}{! 2. Compute horizontal integrals in the x direction}
467   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intx\_dpa)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=hi%jsc,hi%jec
468     \textcolor{keywordflow}{do} i=isq,ieq
469       \textcolor{comment}{! Corner values of T and S}
470       \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
471       \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation}
472       \textcolor{comment}{! of T,S along the top and bottom integrals, almost like thickness}
473       \textcolor{comment}{! weighting.}
474       \textcolor{comment}{! Note: To work in terrain following coordinates we could offset}
475       \textcolor{comment}{! this distance by the layer thickness to replicate other models.}
476       hwght = massweighttoggle * &
477               max(0., -bathyt(i,j)-e(i+1,j,k), -bathyt(i+1,j)-e(i,j,k))
478       \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
479         hl = (e(i,j,k) - e(i,j,k+1)) + dz\_subroundoff
480         hr = (e(i+1,j,k) - e(i+1,j,k+1)) + dz\_subroundoff
481         hwght = hwght * ( (hl-hr)/(hl+hr) )**2
482         idenom = 1./( hwght*(hr + hl) + hl*hr )
483         ttl = ( (hwght*hr)*t\_t(i+1,j,k) + (hwght*hl + hr*hl)*t\_t(i,j,k) ) * idenom
484         ttr = ( (hwght*hl)*t\_t(i,j,k) + (hwght*hr + hr*hl)*t\_t(i+1,j,k) ) * idenom
485         tbl = ( (hwght*hr)*t\_b(i+1,j,k) + (hwght*hl + hr*hl)*t\_b(i,j,k) ) * idenom
486         tbr = ( (hwght*hl)*t\_b(i,j,k) + (hwght*hr + hr*hl)*t\_b(i+1,j,k) ) * idenom
487         stl = ( (hwght*hr)*s\_t(i+1,j,k) + (hwght*hl + hr*hl)*s\_t(i,j,k) ) * idenom
488         str = ( (hwght*hl)*s\_t(i,j,k) + (hwght*hr + hr*hl)*s\_t(i+1,j,k) ) * idenom
489         sbl = ( (hwght*hr)*s\_b(i+1,j,k) + (hwght*hl + hr*hl)*s\_b(i,j,k) ) * idenom
490         sbr = ( (hwght*hl)*s\_b(i,j,k) + (hwght*hr + hr*hl)*s\_b(i+1,j,k) ) * idenom
491       \textcolor{keywordflow}{else}
492         ttl = t\_t(i,j,k); tbl = t\_b(i,j,k); ttr = t\_t(i+1,j,k); tbr = t\_b(i+1,j,k)
493         stl = s\_t(i,j,k); sbl = s\_b(i,j,k); str = s\_t(i+1,j,k); sbr = s\_b(i+1,j,k)
494 \textcolor{keywordflow}{      endif}
495 
496       \textcolor{keywordflow}{do} m=2,4
497         w\_left = wt\_t(m) ; w\_right = wt\_b(m)
498         dz\_x(m,i) = w\_left*(e(i,j,k) - e(i,j,k+1)) + w\_right*(e(i+1,j,k) - e(i+1,j,k+1))
499 
500         \textcolor{comment}{! Salinity and temperature points are linearly interpolated in}
501         \textcolor{comment}{! the horizontal. The subscript (1) refers to the top value in}
502         \textcolor{comment}{! the vertical profile while subscript (5) refers to the bottom}
503         \textcolor{comment}{! value in the vertical profile.}
504         pos = i*15+(m-2)*5
505         t15(pos+1) = w\_left*ttl + w\_right*ttr
506         t15(pos+5) = w\_left*tbl + w\_right*tbr
507 
508         s15(pos+1) = w\_left*stl + w\_right*str
509         s15(pos+5) = w\_left*sbl + w\_right*sbr
510 
511         p15(pos+1) = -gxrho*(w\_left*e(i,j,k) + w\_right*e(i+1,j,k))
512 
513         \textcolor{comment}{! Pressure}
514         \textcolor{keywordflow}{do} n=2,5
515           p15(pos+n) = p15(pos+n-1) + gxrho*0.25*dz\_x(m,i)
516 \textcolor{keywordflow}{        enddo}
517 
518         \textcolor{comment}{! Salinity and temperature (linear interpolation in the vertical)}
519         \textcolor{keywordflow}{do} n=2,4
520           s15(pos+n) = wt\_t(n) * s15(pos+1) + wt\_b(n) * s15(pos+5)
521           t15(pos+n) = wt\_t(n) * t15(pos+1) + wt\_b(n) * t15(pos+5)
522 \textcolor{keywordflow}{        enddo}
523         \textcolor{keywordflow}{if} (use\_vart) t215(pos+1:pos+5) = w\_left*tv%varT(i,j,k) + w\_right*tv%varT(i+1,j,k)
524         \textcolor{keywordflow}{if} (use\_covarts) ts15(pos+1:pos+5) = w\_left*tv%covarTS(i,j,k) + w\_right*tv%covarTS(i+1,j,k)
525         \textcolor{keywordflow}{if} (use\_vars) s215(pos+1:pos+5) = w\_left*tv%varS(i,j,k) + w\_right*tv%varS(i+1,j,k)
526 \textcolor{keywordflow}{      enddo}
527 \textcolor{keywordflow}{    enddo}
528 
529     \textcolor{keywordflow}{if} (use\_stanley\_eos) \textcolor{keywordflow}{then}
530       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
531         \textcolor{keyword}{call }calculate\_density(t15, s15, p15, t215, ts15, s215, r15, 1, 15*(ieq-isq+1), eos, &
532                                rho\_ref=rho\_ref\_mks, scale=rho\_scale)
533       \textcolor{keywordflow}{else}
534         \textcolor{keyword}{call }calculate\_density(t15, s15, p15, t215, ts15, s215, r15, 1, 15*(ieq-isq+1), eos, &
535                                rho\_ref=rho\_ref\_mks)
536 \textcolor{keywordflow}{      endif}
537     \textcolor{keywordflow}{else}
538       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
539         \textcolor{keyword}{call }calculate\_density(t15, s15, p15, r15, 1, 15*(ieq-isq+1), eos, rho\_ref=rho\_ref\_mks, &
540                                scale=rho\_scale)
541       \textcolor{keywordflow}{else}
542         \textcolor{keyword}{call }calculate\_density(t15, s15, p15, r15, 1, 15*(ieq-isq+1), eos, rho\_ref=rho\_ref\_mks)
543 \textcolor{keywordflow}{      endif}
544 \textcolor{keywordflow}{    endif}
545 
546     \textcolor{keywordflow}{do} i=isq,ieq
547       intz(1) = dpa(i,j) ; intz(5) = dpa(i+1,j)
548 
549       \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
550       \textcolor{keywordflow}{do} m = 2,4
551         pos = i*15+(m-2)*5
552         intz(m) = g\_e*dz\_x(m,i)*( c1\_90*(7.0*(r15(pos+1)+r15(pos+5)) + 32.0*(r15(pos+2)+r15(pos+4)) + &
553                           12.0*r15(pos+3)))
554 \textcolor{keywordflow}{      enddo}
555       \textcolor{comment}{! Use Boole's rule to integrate the bottom pressure anomaly values in x.}
556       intx\_dpa(i,j) = c1\_90*(7.0*(intz(1)+intz(5)) + 32.0*(intz(2)+intz(4)) + &
557                              12.0*intz(3))
558 \textcolor{keywordflow}{    enddo}
559 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ endif}
560 
561   \textcolor{comment}{! 3. Compute horizontal integrals in the y direction}
562   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(inty\_dpa)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=jsq,jeq
563     \textcolor{keywordflow}{do} i=hi%isc,hi%iec
564     \textcolor{comment}{! Corner values of T and S}
565     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
566     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation}
567     \textcolor{comment}{! of T,S along the top and bottom integrals, almost like thickness}
568     \textcolor{comment}{! weighting.}
569     \textcolor{comment}{! Note: To work in terrain following coordinates we could offset}
570     \textcolor{comment}{! this distance by the layer thickness to replicate other models.}
571       hwght = massweighttoggle * &
572               max(0., -bathyt(i,j)-e(i,j+1,k), -bathyt(i,j+1)-e(i,j,k))
573       \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
574         hl = (e(i,j,k) - e(i,j,k+1)) + dz\_subroundoff
575         hr = (e(i,j+1,k) - e(i,j+1,k+1)) + dz\_subroundoff
576         hwght = hwght * ( (hl-hr)/(hl+hr) )**2
577         idenom = 1./( hwght*(hr + hl) + hl*hr )
578         ttl = ( (hwght*hr)*t\_t(i,j+1,k) + (hwght*hl + hr*hl)*t\_t(i,j,k) ) * idenom
579         ttr = ( (hwght*hl)*t\_t(i,j,k) + (hwght*hr + hr*hl)*t\_t(i,j+1,k) ) * idenom
580         tbl = ( (hwght*hr)*t\_b(i,j+1,k) + (hwght*hl + hr*hl)*t\_b(i,j,k) ) * idenom
581         tbr = ( (hwght*hl)*t\_b(i,j,k) + (hwght*hr + hr*hl)*t\_b(i,j+1,k) ) * idenom
582         stl = ( (hwght*hr)*s\_t(i,j+1,k) + (hwght*hl + hr*hl)*s\_t(i,j,k) ) * idenom
583         str = ( (hwght*hl)*s\_t(i,j,k) + (hwght*hr + hr*hl)*s\_t(i,j+1,k) ) * idenom
584         sbl = ( (hwght*hr)*s\_b(i,j+1,k) + (hwght*hl + hr*hl)*s\_b(i,j,k) ) * idenom
585         sbr = ( (hwght*hl)*s\_b(i,j,k) + (hwght*hr + hr*hl)*s\_b(i,j+1,k) ) * idenom
586       \textcolor{keywordflow}{else}
587         ttl = t\_t(i,j,k); tbl = t\_b(i,j,k); ttr = t\_t(i,j+1,k); tbr = t\_b(i,j+1,k)
588         stl = s\_t(i,j,k); sbl = s\_b(i,j,k); str = s\_t(i,j+1,k); sbr = s\_b(i,j+1,k)
589 \textcolor{keywordflow}{      endif}
590 
591       \textcolor{keywordflow}{do} m=2,4
592         w\_left = wt\_t(m) ; w\_right = wt\_b(m)
593         dz\_y(m,i) = w\_left*(e(i,j,k) - e(i,j,k+1)) + w\_right*(e(i,j+1,k) - e(i,j+1,k+1))
594 
595         \textcolor{comment}{! Salinity and temperature points are linearly interpolated in}
596         \textcolor{comment}{! the horizontal. The subscript (1) refers to the top value in}
597         \textcolor{comment}{! the vertical profile while subscript (5) refers to the bottom}
598         \textcolor{comment}{! value in the vertical profile.}
599         pos = i*15+(m-2)*5
600         t15(pos+1) = w\_left*ttl + w\_right*ttr
601         t15(pos+5) = w\_left*tbl + w\_right*tbr
602 
603         s15(pos+1) = w\_left*stl + w\_right*str
604         s15(pos+5) = w\_left*sbl + w\_right*sbr
605 
606         p15(pos+1) = -gxrho*(w\_left*e(i,j,k) + w\_right*e(i,j+1,k))
607 
608         \textcolor{comment}{! Pressure}
609         \textcolor{keywordflow}{do} n=2,5
610           p15(pos+n) = p15(pos+n-1) + gxrho*0.25*dz\_y(m,i)
611 \textcolor{keywordflow}{        enddo}
612 
613         \textcolor{comment}{! Salinity and temperature (linear interpolation in the vertical)}
614         \textcolor{keywordflow}{do} n=2,4
615           s15(pos+n) = wt\_t(n) * s15(pos+1) + wt\_b(n) * s15(pos+5)
616           t15(pos+n) = wt\_t(n) * t15(pos+1) + wt\_b(n) * t15(pos+5)
617 \textcolor{keywordflow}{        enddo}
618         \textcolor{keywordflow}{if} (use\_vart) t215(pos+1:pos+5) = w\_left*tv%varT(i,j,k) + w\_right*tv%varT(i,j+1,k)
619         \textcolor{keywordflow}{if} (use\_covarts) ts15(pos+1:pos+5) = w\_left*tv%covarTS(i,j,k) + w\_right*tv%covarTS(i,j+1,k)
620         \textcolor{keywordflow}{if} (use\_vars) s215(pos+1:pos+5) = w\_left*tv%varS(i,j,k) + w\_right*tv%varS(i,j+1,k)
621 \textcolor{keywordflow}{      enddo}
622 \textcolor{keywordflow}{    enddo}
623 
624     \textcolor{keywordflow}{if} (use\_stanley\_eos) \textcolor{keywordflow}{then}
625       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
626         \textcolor{keyword}{call }calculate\_density(t15(15*hi%isc+1:), s15(15*hi%isc+1:), p15(15*hi%isc+1:), &
627                                t215(15*hi%isc+1:), ts15(15*hi%isc+1:), s215(15*hi%isc+1:), &
628                                r15(15*hi%isc+1:), 1, 15*(hi%iec-hi%isc+1), eos, &
629                                rho\_ref=rho\_ref\_mks, scale=rho\_scale)
630       \textcolor{keywordflow}{else}
631         \textcolor{keyword}{call }calculate\_density(t15(15*hi%isc+1:), s15(15*hi%isc+1:), p15(15*hi%isc+1:), &
632                                t215(15*hi%isc+1:), ts15(15*hi%isc+1:), s215(15*hi%isc+1:), &
633                                r15(15*hi%isc+1:), 1, 15*(hi%iec-hi%isc+1), eos, rho\_ref=rho\_ref\_mks)
634 \textcolor{keywordflow}{      endif}
635     \textcolor{keywordflow}{else}
636       \textcolor{keywordflow}{if} (rho\_scale /= 1.0) \textcolor{keywordflow}{then}
637         \textcolor{keyword}{call }calculate\_density(t15(15*hi%isc+1:), s15(15*hi%isc+1:), p15(15*hi%isc+1:), &
638                                r15(15*hi%isc+1:), 1, 15*(hi%iec-hi%isc+1), eos, &
639                                rho\_ref=rho\_ref\_mks, scale=rho\_scale)
640       \textcolor{keywordflow}{else}
641         \textcolor{keyword}{call }calculate\_density(t15(15*hi%isc+1:), s15(15*hi%isc+1:), p15(15*hi%isc+1:), &
642                                r15(15*hi%isc+1:), 1, 15*(hi%iec-hi%isc+1), eos, rho\_ref=rho\_ref\_mks)
643 \textcolor{keywordflow}{      endif}
644 \textcolor{keywordflow}{    endif}
645 
646     \textcolor{keywordflow}{do} i=hi%isc,hi%iec
647       intz(1) = dpa(i,j) ; intz(5) = dpa(i,j+1)
648 
649       \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
650       \textcolor{keywordflow}{do} m = 2,4
651         pos = i*15+(m-2)*5
652         intz(m) = g\_e*dz\_y(m,i)*( c1\_90*(7.0*(r15(pos+1)+r15(pos+5)) + &
653                                          32.0*(r15(pos+2)+r15(pos+4)) + &
654                                          12.0*r15(pos+3)))
655 \textcolor{keywordflow}{      enddo}
656       \textcolor{comment}{! Use Boole's rule to integrate the values.}
657       inty\_dpa(i,j) = c1\_90*(7.0*(intz(1)+intz(5)) + 32.0*(intz(2)+intz(4)) + &
658                              12.0*intz(3))
659 \textcolor{keywordflow}{    enddo}
660 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ endif}
661 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_a3bf090e5cbb58811b3dd0ee4de80f999}\label{namespacemom__density__integrals_a3bf090e5cbb58811b3dd0ee4de80f999}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+ppm@{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+ppm}}
\index{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+ppm@{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+ppm}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+ppm()}{int\_density\_dz\_generic\_ppm()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::int\+\_\+density\+\_\+dz\+\_\+generic\+\_\+ppm (\begin{DoxyParamCaption}\item[{integer, intent(in)}]{k,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi),szk\+\_\+(gv)), intent(in)}]{T\+\_\+t,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi),szk\+\_\+(gv)), intent(in)}]{T\+\_\+b,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi),szk\+\_\+(gv)), intent(in)}]{S\+\_\+t,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi),szk\+\_\+(gv)), intent(in)}]{S\+\_\+b,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi),szk\+\_\+(gv)+1), intent(in)}]{e,  }\item[{real, intent(in)}]{rho\+\_\+ref,  }\item[{real, intent(in)}]{rho\+\_\+0,  }\item[{real, intent(in)}]{G\+\_\+e,  }\item[{real, intent(in)}]{dz\+\_\+subroundoff,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{bathyT,  }\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout)}]{dpa,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intz\+\_\+dpa,  }\item[{real, dimension( hi \%isdb\+: hi \%iedb, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intx\+\_\+dpa,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsdb\+: hi \%jedb), intent(inout), optional}]{inty\+\_\+dpa,  }\item[{logical, intent(in), optional}]{use\+Mass\+Wght\+Interp }\end{DoxyParamCaption})}



Compute pressure gradient force integrals for layer \char`\"{}k\char`\"{} and the case where T and S are parabolic profiles. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em k} & Layer index to calculate integrals for\\
\hline
\mbox{\tt in}  & {\em hi} & Ocean horizontal index structures for the input arrays\\
\hline
\mbox{\tt in}  & {\em gv} & Vertical grid structure\\
\hline
\mbox{\tt in}  & {\em tv} & Thermodynamic variables\\
\hline
\mbox{\tt in}  & {\em t\+\_\+t} & Potential temperature at the cell top \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em t\+\_\+b} & Potential temperature at the cell bottom \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+t} & Salinity at the cell top \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+b} & Salinity at the cell bottom \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em e} & Height of interfaces \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+ref} & A mean density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is subtracted out to reduce the magnitude of each of the integrals.\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+0} & A density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]} or \mbox{[}kg m-\/3\mbox{]}, that is used to calculate the pressure (as p$\sim$=-\/z$\ast$rho\+\_\+0$\ast$\+G\+\_\+e) used in the equation of state.\\
\hline
\mbox{\tt in}  & {\em g\+\_\+e} & The Earth\textquotesingle{}s gravitational acceleration \mbox{[}m s-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em dz\+\_\+subroundoff} & A minuscule thickness change \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
\mbox{\tt in}  & {\em bathyt} & The depth of the bathymetry \mbox{[}Z $\sim$$>$ m\mbox{]}\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in,out}  & {\em dpa} & The change in the pressure anomaly across the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}\\
\hline
\mbox{\tt in,out}  & {\em intz\+\_\+dpa} & The integral through the thickness of the layer of\\
\hline
\mbox{\tt in,out}  & {\em intx\+\_\+dpa} & The integral in x of the difference between the\\
\hline
\mbox{\tt in,out}  & {\em inty\+\_\+dpa} & The integral in y of the difference between the\\
\hline
\mbox{\tt in}  & {\em usemasswghtinterp} & If true, uses mass weighting to interpolate T/S for top and bottom integrals. \\
\hline
\end{DoxyParams}


Definition at line 670 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
670   \textcolor{keywordtype}{integer},              \textcolor{keywordtype}{intent(in)}  :: k\textcolor{comment}{   !< Layer index to calculate integrals for}
671   \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)}  :: HI\textcolor{comment}{  !< Ocean horizontal index structures for the input arrays}
672   \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< Vertical grid structure}
673   \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{  !< Thermodynamic variables}
674   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
675                         \textcolor{keywordtype}{intent(in)}  :: T\_t\textcolor{comment}{ !< Potential temperature at the cell top [degC]}
676   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
677                         \textcolor{keywordtype}{intent(in)}  :: T\_b\textcolor{comment}{ !< Potential temperature at the cell bottom [degC]}
678   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
679                         \textcolor{keywordtype}{intent(in)}  :: S\_t\textcolor{comment}{ !< Salinity at the cell top [ppt]}
680   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV))}, &
681                         \textcolor{keywordtype}{intent(in)}  :: S\_b\textcolor{comment}{ !< Salinity at the cell bottom [ppt]}
682   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI),SZK\_(GV)+1)}, &
683                         \textcolor{keywordtype}{intent(in)}  :: e\textcolor{comment}{   !< Height of interfaces [Z ~> m]}
684   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_ref\textcolor{comment}{ !< A mean density [R ~> kg m-3] or [kg m-3], that is}
685 \textcolor{comment}{                                           !! subtracted out to reduce the magnitude of each of the
       integrals.}
686   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: rho\_0\textcolor{comment}{ !< A density [R ~> kg m-3] or [kg m-3], that is used to
       calculate}
687 \textcolor{comment}{                                           !! the pressure (as p~=-z*rho\_0*G\_e) used in the equation of
       state.}
688   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: G\_e\textcolor{comment}{ !< The Earth's gravitational acceleration [m s-2]}
689   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: dz\_subroundoff\textcolor{comment}{ !< A minuscule thickness change [Z ~> m]}
690   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
691                         \textcolor{keywordtype}{intent(in)}  :: bathyT\textcolor{comment}{ !< The depth of the bathymetry [Z ~> m]}
692   \textcolor{keywordtype}{type}(EOS\_type),       \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
693   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{  !< A dimensional unit scaling type}
694   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
695                         \textcolor{keywordtype}{intent(inout)} :: dpa\textcolor{comment}{ !< The change in the pressure anomaly across the layer [R L2
       T-2 ~> Pa]}
696   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
697               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intz\_dpa\textcolor{comment}{ !< The integral through the thickness of the layer of}
698 \textcolor{comment}{                                           !! the pressure anomaly relative to the anomaly at the}
699 \textcolor{comment}{                                           !! top of the layer [R L2 Z T-2 ~> Pa m]}
700   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(HI),SZJ\_(HI))}, &
701               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intx\_dpa\textcolor{comment}{ !< The integral in x of the difference between the}
702 \textcolor{comment}{                                           !! pressure anomaly at the top and bottom of the layer}
703 \textcolor{comment}{                                           !! divided by the x grid spacing [R L2 T-2 ~> Pa]}
704   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJB\_(HI))}, &
705               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: inty\_dpa\textcolor{comment}{ !< The integral in y of the difference between the}
706 \textcolor{comment}{                                           !! pressure anomaly at the top and bottom of the layer}
707 \textcolor{comment}{                                           !! divided by the y grid spacing [R L2 T-2 ~> Pa]}
708   \textcolor{keywordtype}{logical},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: useMassWghtInterp\textcolor{comment}{ !< If true, uses mass weighting to}
709 \textcolor{comment}{                                           !! interpolate T/S for top and bottom integrals.}
710 
711 \textcolor{comment}{! This subroutine calculates (by numerical quadrature) integrals of}
712 \textcolor{comment}{! pressure anomalies across layers, which are required for calculating the}
713 \textcolor{comment}{! finite-volume form pressure accelerations in a Boussinesq model.  The one}
714 \textcolor{comment}{! potentially dodgy assumption here is that rho\_0 is used both in the denominator}
715 \textcolor{comment}{! of the accelerations, and in the pressure used to calculated density (the}
716 \textcolor{comment}{! latter being -z*rho\_0*G\_e).  These two uses could be separated if need be.}
717 \textcolor{comment}{!}
718 \textcolor{comment}{! It is assumed that the salinity and temperature profiles are parabolic in the}
719 \textcolor{comment}{! vertical. The top and bottom values within each layer are provided and}
720 \textcolor{comment}{! a parabolic interpolation is used to compute intermediate values.}
721 
722   \textcolor{comment}{! Local variables}
723   \textcolor{keywordtype}{real} :: T5(5) \textcolor{comment}{! Temperatures along a line of subgrid locations [degC]}
724   \textcolor{keywordtype}{real} :: S5(5) \textcolor{comment}{! Salinities along a line of subgrid locations [ppt]}
725   \textcolor{keywordtype}{real} :: T25(5) \textcolor{comment}{! SGS temperature variance along a line of subgrid locations [degC2]}
726   \textcolor{keywordtype}{real} :: TS5(5) \textcolor{comment}{! SGS temperature-salinity covariance along a line of subgrid locations [degC ppt]}
727   \textcolor{keywordtype}{real} :: S25(5) \textcolor{comment}{! SGS salinity variance along a line of subgrid locations [ppt2]}
728   \textcolor{keywordtype}{real} :: p5(5) \textcolor{comment}{! Pressures at five quadrature points, never rescaled from Pa [Pa]}
729   \textcolor{keywordtype}{real} :: r5(5) \textcolor{comment}{! Density anomalies from rho\_ref at quadrature points [R ~> kg m-3] or [kg m-3]}
730   \textcolor{keywordtype}{real} :: wt\_t(5), wt\_b(5) \textcolor{comment}{! Top and bottom weights [nondim]}
731   \textcolor{keywordtype}{real} :: rho\_anom \textcolor{comment}{! The integrated density anomaly [R ~> kg m-3] or [kg m-3]}
732   \textcolor{keywordtype}{real} :: w\_left, w\_right  \textcolor{comment}{! Left and right weights [nondim]}
733   \textcolor{keywordtype}{real} :: intz(5) \textcolor{comment}{! The gravitational acceleration times the integrals of density}
734                   \textcolor{comment}{! with height at the 5 sub-column locations [R L2 T-2 ~> Pa] or [Pa]}
735   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_90 = 1.0/90.0  \textcolor{comment}{! Rational constants.}
736   \textcolor{keywordtype}{real} :: GxRho \textcolor{comment}{! The gravitational acceleration times density and unit conversion factors [Pa Z-1 ~> kg
       m-2 s-2]}
737   \textcolor{keywordtype}{real} :: I\_Rho \textcolor{comment}{! The inverse of the Boussinesq density [R-1 ~> m3 kg-1] or [m3 kg-1]}
738   \textcolor{keywordtype}{real} :: rho\_scale \textcolor{comment}{! A scaling factor for densities from kg m-3 to R [R m3 kg-1 ~> 1]}
739   \textcolor{keywordtype}{real} :: rho\_ref\_mks \textcolor{comment}{! The reference density in MKS units, never rescaled from kg m-3 [kg m-3]}
740   \textcolor{keywordtype}{real} :: dz \textcolor{comment}{! Layer thicknesses at tracer points [Z ~> m]}
741   \textcolor{keywordtype}{real} :: massWeightToggle \textcolor{comment}{! A non-dimensional toggle factor (0 or 1) [nondim]}
742   \textcolor{keywordtype}{real} :: Ttl, Tbl, Tml, Ttr, Tbr, Tmr \textcolor{comment}{! Temperatures at the velocity cell corners [degC]}
743   \textcolor{keywordtype}{real} :: Stl, Sbl, Sml, Str, Sbr, Smr \textcolor{comment}{! Salinities at the velocity cell corners [ppt]}
744   \textcolor{keywordtype}{real} :: s6 \textcolor{comment}{! PPM curvature coefficient for S [ppt]}
745   \textcolor{keywordtype}{real} :: t6 \textcolor{comment}{! PPM curvature coefficient for T [degC]}
746   \textcolor{keywordtype}{real} :: T\_top, T\_mn, T\_bot \textcolor{comment}{! Left edge, cell mean and right edge values used in PPM reconstructions of T}
747   \textcolor{keywordtype}{real} :: S\_top, S\_mn, S\_bot \textcolor{comment}{! Left edge, cell mean and right edge values used in PPM reconstructions of S}
748   \textcolor{keywordtype}{real} :: hWght  \textcolor{comment}{! A topographically limited thicknes weight [Z ~> m]}
749   \textcolor{keywordtype}{real} :: hL, hR \textcolor{comment}{! Thicknesses to the left and right [Z ~> m]}
750   \textcolor{keywordtype}{real} :: iDenom \textcolor{comment}{! The denominator of the thickness weight expressions [Z-2 ~> m-2]}
751   \textcolor{keywordtype}{integer} :: Isq, Ieq, Jsq, Jeq, i, j, m, n
752   \textcolor{keywordtype}{logical} :: use\_PPM \textcolor{comment}{! If false, assume zero curvature in reconstruction, i.e. PLM}
753   \textcolor{keywordtype}{logical} :: use\_stanley\_eos \textcolor{comment}{! True is SGS variance fields exist in tv.}
754   \textcolor{keywordtype}{logical} :: use\_varT, use\_varS, use\_covarTS
755 
756   isq = hi%IscB ; ieq = hi%IecB ; jsq = hi%JscB ; jeq = hi%JecB
757 
758   rho\_scale = us%kg\_m3\_to\_R
759   gxrho = us%RL2\_T2\_to\_Pa * g\_e * rho\_0
760   rho\_ref\_mks = rho\_ref * us%R\_to\_kg\_m3
761   i\_rho = 1.0 / rho\_0
762   massweighttoggle = 0.
763   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(usemasswghtinterp)) \textcolor{keywordflow}{then}
764     \textcolor{keywordflow}{if} (usemasswghtinterp) massweighttoggle = 1.
765 \textcolor{keywordflow}{  endif}
766 
767   \textcolor{comment}{! In event PPM calculation is bypassed with use\_PPM=False}
768   s6 = 0.
769   t6 = 0.
770   use\_ppm = .true. \textcolor{comment}{! This is a place-holder to allow later re-use of this function}
771 
772   use\_vart = \textcolor{keyword}{associated}(tv%varT)
773   use\_covarts = \textcolor{keyword}{associated}(tv%covarTS)
774   use\_vars = \textcolor{keyword}{associated}(tv%varS)
775   use\_stanley\_eos = use\_vart .or. use\_covarts .or. use\_vars
776   t25(:) = 0.
777   ts5(:) = 0.
778   s25(:) = 0.
779 
780   \textcolor{keywordflow}{do} n = 1, 5
781     wt\_t(n) = 0.25 * \textcolor{keywordtype}{real}(5-n)
782     wt\_b(n) = 1.0 - wt\_t(n)
783 \textcolor{keywordflow}{  enddo}
784 
785   \textcolor{comment}{! 1. Compute vertical integrals}
786   \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1
787     \textcolor{keywordflow}{if} (use\_ppm) \textcolor{keywordflow}{then}
788       \textcolor{comment}{! Curvature coefficient of the parabolas}
789       s6 = 3.0 * ( 2.0*tv%S(i,j,k) - ( s\_t(i,j,k) + s\_b(i,j,k) ) )
790       t6 = 3.0 * ( 2.0*tv%T(i,j,k) - ( t\_t(i,j,k) + t\_b(i,j,k) ) )
791 \textcolor{keywordflow}{    endif}
792     dz = e(i,j,k) - e(i,j,k+1)
793     \textcolor{keywordflow}{do} n=1,5
794       p5(n) = -gxrho*(e(i,j,k) - 0.25*\textcolor{keywordtype}{real}(n-1)*dz)
795       \textcolor{comment}{! Salinity and temperature points are reconstructed with PPM}
796       s5(n) = wt\_t(n) * s\_t(i,j,k) + wt\_b(n) * ( s\_b(i,j,k) + s6 * wt\_t(n) )
797       t5(n) = wt\_t(n) * t\_t(i,j,k) + wt\_b(n) * ( t\_b(i,j,k) + t6 * wt\_t(n) )
798 \textcolor{keywordflow}{    enddo}
799     \textcolor{keywordflow}{if} (use\_stanley\_eos) \textcolor{keywordflow}{then}
800       \textcolor{keywordflow}{if} (use\_vart) t25(:) = tv%varT(i,j,k)
801       \textcolor{keywordflow}{if} (use\_covarts) ts5(:) = tv%covarTS(i,j,k)
802       \textcolor{keywordflow}{if} (use\_vars) s25(:) = tv%varS(i,j,k)
803       \textcolor{keyword}{call }calculate\_density(t5, s5, p5, t25, ts5, s25, r5, &
804                              1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
805     \textcolor{keywordflow}{else}
806       \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
807 \textcolor{keywordflow}{    endif}
808 
809     \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
810     rho\_anom = c1\_90*(7.0*(r5(1)+r5(5)) + 32.0*(r5(2)+r5(4)) + 12.0*r5(3))
811     dpa(i,j) = g\_e*dz*rho\_anom
812     \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intz\_dpa)) \textcolor{keywordflow}{then}
813       \textcolor{comment}{! Use a Boole's-rule-like fifth-order accurate estimate of}
814       \textcolor{comment}{! the double integral of the pressure anomaly.}
815       intz\_dpa(i,j) = 0.5*g\_e*dz**2 * &
816                       (rho\_anom - c1\_90*(16.0*(r5(4)-r5(2)) + 7.0*(r5(5)-r5(1))) )
817 \textcolor{keywordflow}{    endif}
818 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! end loops on j and i}
819 
820   \textcolor{comment}{! 2. Compute horizontal integrals in the x direction}
821   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intx\_dpa)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=hi%jsc,hi%jec ; \textcolor{keywordflow}{do} i=isq,ieq
822     \textcolor{comment}{! Corner values of T and S}
823     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
824     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation}
825     \textcolor{comment}{! of T,S along the top and bottom integrals, almost like thickness}
826     \textcolor{comment}{! weighting.}
827     \textcolor{comment}{! Note: To work in terrain following coordinates we could offset}
828     \textcolor{comment}{! this distance by the layer thickness to replicate other models.}
829     hwght = massweighttoggle * &
830             max(0., -bathyt(i,j)-e(i+1,j,k), -bathyt(i+1,j)-e(i,j,k))
831     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
832       hl = (e(i,j,k) - e(i,j,k+1)) + dz\_subroundoff
833       hr = (e(i+1,j,k) - e(i+1,j,k+1)) + dz\_subroundoff
834       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
835       idenom = 1./( hwght*(hr + hl) + hl*hr )
836       ttl = ( (hwght*hr)*t\_t(i+1,j,k) + (hwght*hl + hr*hl)*t\_t(i,j,k) ) * idenom
837       tbl = ( (hwght*hr)*t\_b(i+1,j,k) + (hwght*hl + hr*hl)*t\_b(i,j,k) ) * idenom
838       tml = ( (hwght*hr)*tv%T(i+1,j,k)+ (hwght*hl + hr*hl)*tv%T(i,j,k) ) * idenom
839       ttr = ( (hwght*hl)*t\_t(i,j,k) + (hwght*hr + hr*hl)*t\_t(i+1,j,k) ) * idenom
840       tbr = ( (hwght*hl)*t\_b(i,j,k) + (hwght*hr + hr*hl)*t\_b(i+1,j,k) ) * idenom
841       tmr = ( (hwght*hl)*tv%T(i,j,k) + (hwght*hr + hr*hl)*tv%T(i+1,j,k) ) * idenom
842       stl = ( (hwght*hr)*s\_t(i+1,j,k) + (hwght*hl + hr*hl)*s\_t(i,j,k) ) * idenom
843       sbl = ( (hwght*hr)*s\_b(i+1,j,k) + (hwght*hl + hr*hl)*s\_b(i,j,k) ) * idenom
844       sml = ( (hwght*hr)*tv%S(i+1,j,k) + (hwght*hl + hr*hl)*tv%S(i,j,k) ) * idenom
845       str = ( (hwght*hl)*s\_t(i,j,k) + (hwght*hr + hr*hl)*s\_t(i+1,j,k) ) * idenom
846       sbr = ( (hwght*hl)*s\_b(i,j,k) + (hwght*hr + hr*hl)*s\_b(i+1,j,k) ) * idenom
847       smr = ( (hwght*hl)*tv%S(i,j,k) + (hwght*hr + hr*hl)*tv%S(i+1,j,k) ) * idenom
848     \textcolor{keywordflow}{else}
849       ttl = t\_t(i,j,k); tbl = t\_b(i,j,k); ttr = t\_t(i+1,j,k); tbr = t\_b(i+1,j,k)
850       tml = tv%T(i,j,k); tmr = tv%T(i+1,j,k)
851       stl = s\_t(i,j,k); sbl = s\_b(i,j,k); str = s\_t(i+1,j,k); sbr = s\_b(i+1,j,k)
852       sml = tv%S(i,j,k); smr = tv%S(i+1,j,k)
853 \textcolor{keywordflow}{    endif}
854 
855     \textcolor{keywordflow}{do} m=2,4
856       w\_left = wt\_t(m) ; w\_right = wt\_b(m)
857 
858       \textcolor{comment}{! Salinity and temperature points are linearly interpolated in}
859       \textcolor{comment}{! the horizontal. The subscript (1) refers to the top value in}
860       \textcolor{comment}{! the vertical profile while subscript (5) refers to the bottom}
861       \textcolor{comment}{! value in the vertical profile.}
862       t\_top = w\_left*ttl + w\_right*ttr
863       t\_mn = w\_left*tml + w\_right*tmr
864       t\_bot = w\_left*tbl + w\_right*tbr
865 
866       s\_top = w\_left*stl + w\_right*str
867       s\_mn = w\_left*sml + w\_right*smr
868       s\_bot = w\_left*sbl + w\_right*sbr
869 
870       \textcolor{comment}{! Pressure}
871       dz = w\_left*(e(i,j,k) - e(i,j,k+1)) + w\_right*(e(i+1,j,k) - e(i+1,j,k+1))
872       p5(1) = -gxrho*(w\_left*e(i,j,k) + w\_right*e(i+1,j,k))
873       \textcolor{keywordflow}{do} n=2,5
874         p5(n) = p5(n-1) + gxrho*0.25*dz
875 \textcolor{keywordflow}{      enddo}
876 
877       \textcolor{comment}{! Parabolic reconstructions in the vertical for T and S}
878       \textcolor{keywordflow}{if} (use\_ppm) \textcolor{keywordflow}{then}
879         \textcolor{comment}{! Coefficients of the parabolas}
880         s6 = 3.0 * ( 2.0*s\_mn - ( s\_top + s\_bot ) )
881         t6 = 3.0 * ( 2.0*t\_mn - ( t\_top + t\_bot ) )
882 \textcolor{keywordflow}{      endif}
883       \textcolor{keywordflow}{do} n=1,5
884         s5(n) = wt\_t(n) * s\_top + wt\_b(n) * ( s\_bot + s6 * wt\_t(n) )
885         t5(n) = wt\_t(n) * t\_top + wt\_b(n) * ( t\_bot + t6 * wt\_t(n) )
886 \textcolor{keywordflow}{      enddo}
887 
888       \textcolor{keywordflow}{if} (use\_stanley\_eos) \textcolor{keywordflow}{then}
889         \textcolor{keywordflow}{if} (use\_vart) t25(:) = w\_left*tv%varT(i,j,k) + w\_right*tv%varT(i+1,j,k)
890         \textcolor{keywordflow}{if} (use\_covarts) ts5(:) = w\_left*tv%covarTS(i,j,k) + w\_right*tv%covarTS(i+1,j,k)
891         \textcolor{keywordflow}{if} (use\_vars) s25(:) = w\_left*tv%varS(i,j,k) + w\_right*tv%varS(i+1,j,k)
892         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, t25, ts5, s25, r5, &
893                                1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
894       \textcolor{keywordflow}{else}
895         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
896 \textcolor{keywordflow}{      endif}
897 
898       \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
899       intz(m) = g\_e*dz*( c1\_90*(7.0*(r5(1)+r5(5)) + 32.0*(r5(2)+r5(4)) + 12.0*r5(3)) )
900 \textcolor{keywordflow}{    enddo} \textcolor{comment}{! m}
901     intz(1) = dpa(i,j) ; intz(5) = dpa(i+1,j)
902 
903     \textcolor{comment}{! Use Boole's rule to integrate the bottom pressure anomaly values in x.}
904     intx\_dpa(i,j) = c1\_90*(7.0*(intz(1)+intz(5)) + 32.0*(intz(2)+intz(4)) + 12.0*intz(3))
905 
906 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
907 
908   \textcolor{comment}{! 3. Compute horizontal integrals in the y direction}
909   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(inty\_dpa)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=hi%isc,hi%iec
910     \textcolor{comment}{! Corner values of T and S}
911     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
912     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation}
913     \textcolor{comment}{! of T,S along the top and bottom integrals, almost like thickness}
914     \textcolor{comment}{! weighting.}
915     \textcolor{comment}{! Note: To work in terrain following coordinates we could offset}
916     \textcolor{comment}{! this distance by the layer thickness to replicate other models.}
917     hwght = massweighttoggle * &
918             max(0., -bathyt(i,j)-e(i,j+1,k), -bathyt(i,j+1)-e(i,j,k))
919     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
920       hl = (e(i,j,k) - e(i,j,k+1)) + dz\_subroundoff
921       hr = (e(i,j+1,k) - e(i,j+1,k+1)) + dz\_subroundoff
922       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
923       idenom = 1./( hwght*(hr + hl) + hl*hr )
924       ttl = ( (hwght*hr)*t\_t(i,j+1,k) + (hwght*hl + hr*hl)*t\_t(i,j,k) ) * idenom
925       tbl = ( (hwght*hr)*t\_b(i,j+1,k) + (hwght*hl + hr*hl)*t\_b(i,j,k) ) * idenom
926       tml = ( (hwght*hr)*tv%T(i,j+1,k)+ (hwght*hl + hr*hl)*tv%T(i,j,k) ) * idenom
927       ttr = ( (hwght*hl)*t\_t(i,j,k) + (hwght*hr + hr*hl)*t\_t(i,j+1,k) ) * idenom
928       tbr = ( (hwght*hl)*t\_b(i,j,k) + (hwght*hr + hr*hl)*t\_b(i,j+1,k) ) * idenom
929       tmr = ( (hwght*hl)*tv%T(i,j,k) + (hwght*hr + hr*hl)*tv%T(i,j+1,k) ) * idenom
930       stl = ( (hwght*hr)*s\_t(i,j+1,k) + (hwght*hl + hr*hl)*s\_t(i,j,k) ) * idenom
931       sbl = ( (hwght*hr)*s\_b(i,j+1,k) + (hwght*hl + hr*hl)*s\_b(i,j,k) ) * idenom
932       sml = ( (hwght*hr)*tv%S(i,j+1,k)+ (hwght*hl + hr*hl)*tv%S(i,j,k) ) * idenom
933       str = ( (hwght*hl)*s\_t(i,j,k) + (hwght*hr + hr*hl)*s\_t(i,j+1,k) ) * idenom
934       sbr = ( (hwght*hl)*s\_b(i,j,k) + (hwght*hr + hr*hl)*s\_b(i,j+1,k) ) * idenom
935       smr = ( (hwght*hl)*tv%S(i,j,k) + (hwght*hr + hr*hl)*tv%S(i,j+1,k) ) * idenom
936     \textcolor{keywordflow}{else}
937       ttl = t\_t(i,j,k); tbl = t\_b(i,j,k); ttr = t\_t(i,j+1,k); tbr = t\_b(i,j+1,k)
938       tml = tv%T(i,j,k); tmr = tv%T(i,j+1,k)
939       stl = s\_t(i,j,k); sbl = s\_b(i,j,k); str = s\_t(i,j+1,k); sbr = s\_b(i,j+1,k)
940       sml = tv%S(i,j,k); smr = tv%S(i,j+1,k)
941 \textcolor{keywordflow}{    endif}
942 
943     \textcolor{keywordflow}{do} m=2,4
944       w\_left = wt\_t(m) ; w\_right = wt\_b(m)
945 
946       \textcolor{comment}{! Salinity and temperature points are linearly interpolated in}
947       \textcolor{comment}{! the horizontal. The subscript (1) refers to the top value in}
948       \textcolor{comment}{! the vertical profile while subscript (5) refers to the bottom}
949       \textcolor{comment}{! value in the vertical profile.}
950       t\_top = w\_left*ttl + w\_right*ttr
951       t\_mn = w\_left*tml + w\_right*tmr
952       t\_bot = w\_left*tbl + w\_right*tbr
953 
954       s\_top = w\_left*stl + w\_right*str
955       s\_mn = w\_left*sml + w\_right*smr
956       s\_bot = w\_left*sbl + w\_right*sbr
957 
958       \textcolor{comment}{! Pressure}
959       dz = w\_left*(e(i,j,k) - e(i,j,k+1)) + w\_right*(e(i,j+1,k) - e(i,j+1,k+1))
960       p5(1) = -gxrho*(w\_left*e(i,j,k) + w\_right*e(i,j+1,k))
961       \textcolor{keywordflow}{do} n=2,5
962         p5(n) = p5(n-1) + gxrho*0.25*dz
963 \textcolor{keywordflow}{      enddo}
964 
965       \textcolor{comment}{! Parabolic reconstructions in the vertical for T and S}
966       \textcolor{keywordflow}{if} (use\_ppm) \textcolor{keywordflow}{then}
967         \textcolor{comment}{! Coefficients of the parabolas}
968         s6 = 3.0 * ( 2.0*s\_mn - ( s\_top + s\_bot ) )
969         t6 = 3.0 * ( 2.0*t\_mn - ( t\_top + t\_bot ) )
970 \textcolor{keywordflow}{      endif}
971       \textcolor{keywordflow}{do} n=1,5
972         s5(n) = wt\_t(n) * s\_top + wt\_b(n) * ( s\_bot + s6 * wt\_t(n) )
973         t5(n) = wt\_t(n) * t\_top + wt\_b(n) * ( t\_bot + t6 * wt\_t(n) )
974 \textcolor{keywordflow}{      enddo}
975 
976       \textcolor{keywordflow}{if} (use\_stanley\_eos) \textcolor{keywordflow}{then}
977         \textcolor{keywordflow}{if} (use\_vart) t25(:) = w\_left*tv%varT(i,j,k) + w\_right*tv%varT(i,j+1,k)
978         \textcolor{keywordflow}{if} (use\_covarts) ts5(:) = w\_left*tv%covarTS(i,j,k) + w\_right*tv%covarTS(i,j+1,k)
979         \textcolor{keywordflow}{if} (use\_vars) s25(:) = w\_left*tv%varS(i,j,k) + w\_right*tv%varS(i,j+1,k)
980         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, t25, ts5, s25, r5, &
981                                1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
982       \textcolor{keywordflow}{else}
983         \textcolor{keyword}{call }calculate\_density(t5, s5, p5, r5, 1, 5, eos, rho\_ref=rho\_ref\_mks, scale=rho\_scale)
984 \textcolor{keywordflow}{      endif}
985 
986       \textcolor{comment}{! Use Boole's rule to estimate the pressure anomaly change.}
987       intz(m) = g\_e*dz*( c1\_90*(7.0*(r5(1)+r5(5)) + 32.0*(r5(2)+r5(4)) + 12.0*r5(3)) )
988 \textcolor{keywordflow}{    enddo} \textcolor{comment}{! m}
989     intz(1) = dpa(i,j) ; intz(5) = dpa(i,j+1)
990 
991     \textcolor{comment}{! Use Boole's rule to integrate the bottom pressure anomaly values in y.}
992     inty\_dpa(i,j) = c1\_90*(7.0*(intz(1)+intz(5)) + 32.0*(intz(2)+intz(4)) + 12.0*intz(3))
993 
994 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
995 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_a86ebdfebaaea2ac0339dcabed482dd11}\label{namespacemom__density__integrals_a86ebdfebaaea2ac0339dcabed482dd11}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+pcm@{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+pcm}}
\index{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+pcm@{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+pcm}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+pcm()}{int\_spec\_vol\_dp\_generic\_pcm()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+pcm (\begin{DoxyParamCaption}\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{T,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{S,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{p\+\_\+t,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{p\+\_\+b,  }\item[{real, intent(in)}]{alpha\+\_\+ref,  }\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout)}]{dza,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intp\+\_\+dza,  }\item[{real, dimension( hi \%isdb\+: hi \%iedb, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intx\+\_\+dza,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsdb\+: hi \%jedb), intent(inout), optional}]{inty\+\_\+dza,  }\item[{integer, intent(in), optional}]{halo\+\_\+size,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in), optional}]{bathyP,  }\item[{real, intent(in), optional}]{d\+P\+\_\+neglect,  }\item[{logical, intent(in), optional}]{use\+Mass\+Wght\+Interp }\end{DoxyParamCaption})}



This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-\/volume form pressure accelerations in a non-\/\+Boussinesq model. There are essentially no free assumptions, apart from the use of Boole\textquotesingle{}s rule quadrature to do the integrals. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em hi} & A horizontal index type structure.\\
\hline
\mbox{\tt in}  & {\em t} & Potential temperature of the layer \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s} & Salinity of the layer \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+t} & Pressure atop the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+b} & Pressure below the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em alpha\+\_\+ref} & A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals \mbox{[}R-\/1 $\sim$$>$ m3 kg-\/1\mbox{]} The calculation is mathematically identical with different values of alpha\+\_\+ref, but alpha\+\_\+ref alters the effects of roundoff, and answers do change.\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in,out}  & {\em dza} & The change in the geopotential anomaly\\
\hline
\mbox{\tt in,out}  & {\em intp\+\_\+dza} & The integral in pressure through the layer of\\
\hline
\mbox{\tt in,out}  & {\em intx\+\_\+dza} & The integral in x of the difference between\\
\hline
\mbox{\tt in,out}  & {\em inty\+\_\+dza} & The integral in y of the difference between\\
\hline
\mbox{\tt in}  & {\em halo\+\_\+size} & The width of halo points on which to calculate dza.\\
\hline
\mbox{\tt in}  & {\em bathyp} & The pressure at the bathymetry \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em dp\+\_\+neglect} & A minuscule pressure change with the same units as p\+\_\+t \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em usemasswghtinterp} & If true, uses mass weighting to interpolate T/S for top and bottom integrals. \\
\hline
\end{DoxyParams}


Definition at line 1065 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
1065   \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)}  :: HI\textcolor{comment}{ !< A horizontal index type structure.}
1066   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1067                         \textcolor{keywordtype}{intent(in)}  :: T\textcolor{comment}{  !< Potential temperature of the layer [degC]}
1068   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1069                         \textcolor{keywordtype}{intent(in)}  :: S\textcolor{comment}{  !< Salinity of the layer [ppt]}
1070   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1071                         \textcolor{keywordtype}{intent(in)}  :: p\_t\textcolor{comment}{ !< Pressure atop the layer [R L2 T-2 ~> Pa] or [Pa]}
1072   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1073                         \textcolor{keywordtype}{intent(in)}  :: p\_b\textcolor{comment}{ !< Pressure below the layer [R L2 T-2 ~> Pa] or [Pa]}
1074   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: alpha\_ref\textcolor{comment}{ !< A mean specific volume that is subtracted out}
1075 \textcolor{comment}{                            !! to reduce the magnitude of each of the integrals [R-1 ~> m3 kg-1]}
1076 \textcolor{comment}{                            !! The calculation is mathematically identical with different values of}
1077 \textcolor{comment}{                            !! alpha\_ref, but alpha\_ref alters the effects of roundoff, and}
1078 \textcolor{comment}{                            !! answers do change.}
1079   \textcolor{keywordtype}{type}(EOS\_type),       \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
1080   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type}
1081   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1082                         \textcolor{keywordtype}{intent(inout)} :: dza\textcolor{comment}{ !< The change in the geopotential anomaly}
1083 \textcolor{comment}{                            !! across the layer [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1084   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1085               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intp\_dza\textcolor{comment}{ !< The integral in pressure through the layer of}
1086 \textcolor{comment}{                            !! the geopotential anomaly relative to the anomaly at the bottom of the}
1087 \textcolor{comment}{                            !! layer [R L4 T-4 ~> Pa m2 s-2] or [Pa m2 s-2]}
1088   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(HI),SZJ\_(HI))}, &
1089               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intx\_dza\textcolor{comment}{  !< The integral in x of the difference between}
1090 \textcolor{comment}{                            !! the geopotential anomaly at the top and bottom of the layer divided}
1091 \textcolor{comment}{                            !! by the x grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1092   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJB\_(HI))}, &
1093               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: inty\_dza\textcolor{comment}{  !< The integral in y of the difference between}
1094 \textcolor{comment}{                            !! the geopotential anomaly at the top and bottom of the layer divided}
1095 \textcolor{comment}{                            !! by the y grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1096   \textcolor{keywordtype}{integer},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: halo\_size\textcolor{comment}{ !< The width of halo points on which to calculate dza.}
1097   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1098               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: bathyP\textcolor{comment}{ !< The pressure at the bathymetry [R L2 T-2 ~> Pa] or [Pa]}
1099   \textcolor{keywordtype}{real},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: dP\_neglect\textcolor{comment}{ !< A minuscule pressure change with}
1100 \textcolor{comment}{                                             !! the same units as p\_t [R L2 T-2 ~> Pa] or [Pa]}
1101   \textcolor{keywordtype}{logical},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: useMassWghtInterp\textcolor{comment}{ !< If true, uses mass weighting}
1102 \textcolor{comment}{                            !! to interpolate T/S for top and bottom integrals.}
1103 
1104 \textcolor{comment}{!   This subroutine calculates analytical and nearly-analytical integrals in}
1105 \textcolor{comment}{! pressure across layers of geopotential anomalies, which are required for}
1106 \textcolor{comment}{! calculating the finite-volume form pressure accelerations in a non-Boussinesq}
1107 \textcolor{comment}{! model.  There are essentially no free assumptions, apart from the use of}
1108 \textcolor{comment}{! Boole's rule to do the horizontal integrals, and from a truncation in the}
1109 \textcolor{comment}{! series for log(1-eps/1+eps) that assumes that |eps| < 0.34.}
1110 
1111   \textcolor{comment}{! Local variables}
1112   \textcolor{keywordtype}{real} :: T5(5)      \textcolor{comment}{! Temperatures at five quadrature points [degC]}
1113   \textcolor{keywordtype}{real} :: S5(5)      \textcolor{comment}{! Salinities at five quadrature points [ppt]}
1114   \textcolor{keywordtype}{real} :: p5(5)      \textcolor{comment}{! Pressures at five quadrature points, scaled back to Pa if necessary [Pa]}
1115   \textcolor{keywordtype}{real} :: a5(5)      \textcolor{comment}{! Specific volumes at five quadrature points [R-1 ~> m3 kg-1] or [m3 kg-1]}
1116   \textcolor{keywordtype}{real} :: alpha\_anom \textcolor{comment}{! The depth averaged specific density anomaly [R-1 ~> m3 kg-1]}
1117   \textcolor{keywordtype}{real} :: dp         \textcolor{comment}{! The pressure change through a layer [R L2 T-2 ~> Pa]}
1118   \textcolor{keywordtype}{real} :: hWght      \textcolor{comment}{! A pressure-thickness below topography [R L2 T-2 ~> Pa]}
1119   \textcolor{keywordtype}{real} :: hL, hR     \textcolor{comment}{! Pressure-thicknesses of the columns to the left and right [R L2 T-2 ~> Pa]}
1120   \textcolor{keywordtype}{real} :: alpha\_ref\_mks \textcolor{comment}{! The reference specific volume in MKS units, never rescaled from m3 kg-1 [m3 kg-1]}
1121   \textcolor{keywordtype}{real} :: iDenom     \textcolor{comment}{! The inverse of the denominator in the weights [T4 R-2 L-4 ~> Pa-2]}
1122   \textcolor{keywordtype}{real} :: hWt\_LL, hWt\_LR \textcolor{comment}{! hWt\_LA is the weighted influence of A on the left column [nondim]}
1123   \textcolor{keywordtype}{real} :: hWt\_RL, hWt\_RR \textcolor{comment}{! hWt\_RA is the weighted influence of A on the right column [nondim]}
1124   \textcolor{keywordtype}{real} :: wt\_L, wt\_R \textcolor{comment}{! The linear weights of the left and right columns [nondim]}
1125   \textcolor{keywordtype}{real} :: wtT\_L, wtT\_R \textcolor{comment}{! The weights for tracers from the left and right columns [nondim]}
1126   \textcolor{keywordtype}{real} :: intp(5)    \textcolor{comment}{! The integrals of specific volume with pressure at the}
1127                      \textcolor{comment}{! 5 sub-column locations [L2 T-2 ~> m2 s-2]}
1128   \textcolor{keywordtype}{real} :: RL2\_T2\_to\_Pa  \textcolor{comment}{! A unit conversion factor from the rescaled units of pressure to Pa [Pa T2 R-1 L-2
       ~> 1]}
1129   \textcolor{keywordtype}{real} :: SV\_scale   \textcolor{comment}{! A multiplicative factor by which to scale specific}
1130                      \textcolor{comment}{! volume from m3 kg-1 to the desired units [kg m-3 R-1 ~> 1]}
1131   \textcolor{keywordtype}{logical} :: do\_massWeight \textcolor{comment}{! Indicates whether to do mass weighting.}
1132   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_90 = 1.0/90.0  \textcolor{comment}{! A rational constant.}
1133   \textcolor{keywordtype}{integer} :: Isq, Ieq, Jsq, Jeq, ish, ieh, jsh, jeh, i, j, m, n, halo
1134 
1135   isq = hi%IscB ; ieq = hi%IecB ; jsq = hi%JscB ; jeq = hi%JecB
1136   halo = 0 ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(halo\_size)) halo = max(halo\_size,0)
1137   ish = hi%isc-halo ; ieh = hi%iec+halo ; jsh = hi%jsc-halo ; jeh = hi%jec+halo
1138   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intx\_dza)) \textcolor{keywordflow}{then} ; ish = min(isq,ish) ; ieh = max(ieq+1,ieh);\textcolor{keywordflow}{ endif}
1139   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(inty\_dza)) \textcolor{keywordflow}{then} ; jsh = min(jsq,jsh) ; jeh = max(jeq+1,jeh);\textcolor{keywordflow}{ endif}
1140 
1141   sv\_scale = us%R\_to\_kg\_m3
1142   rl2\_t2\_to\_pa = us%RL2\_T2\_to\_Pa
1143   alpha\_ref\_mks = alpha\_ref * us%kg\_m3\_to\_R
1144 
1145   do\_massweight = .false.
1146   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(usemasswghtinterp)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (usemasswghtinterp) \textcolor{keywordflow}{then}
1147     do\_massweight = .true.
1148     \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{present}(bathyp)) \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"int\_spec\_vol\_dp\_generic: "}//&
1149         \textcolor{stringliteral}{"bathyP must be present if useMassWghtInterp is present and true."})
1150     \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{present}(dp\_neglect)) \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"int\_spec\_vol\_dp\_generic: "}//&
1151         \textcolor{stringliteral}{"dP\_neglect must be present if useMassWghtInterp is present and true."})
1152 \textcolor{keywordflow}{  endif} ;\textcolor{keywordflow}{ endif}
1153 
1154   \textcolor{keywordflow}{do} j=jsh,jeh ; \textcolor{keywordflow}{do} i=ish,ieh
1155     dp = p\_b(i,j) - p\_t(i,j)
1156     \textcolor{keywordflow}{do} n=1,5
1157       t5(n) = t(i,j) ; s5(n) = s(i,j)
1158       p5(n) = rl2\_t2\_to\_pa * (p\_b(i,j) - 0.25*\textcolor{keywordtype}{real}(n-1)*dp)
1159 \textcolor{keywordflow}{    enddo}
1160 
1161     \textcolor{keywordflow}{if} (sv\_scale /= 1.0) \textcolor{keywordflow}{then}
1162       \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks, scale=sv\_scale)
1163     \textcolor{keywordflow}{else}
1164       \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks)
1165 \textcolor{keywordflow}{    endif}
1166 
1167     \textcolor{comment}{! Use Boole's rule to estimate the interface height anomaly change.}
1168     alpha\_anom = c1\_90*(7.0*(a5(1)+a5(5)) + 32.0*(a5(2)+a5(4)) + 12.0*a5(3))
1169     dza(i,j) = dp*alpha\_anom
1170     \textcolor{comment}{! Use a Boole's-rule-like fifth-order accurate estimate of the double integral of}
1171     \textcolor{comment}{! the interface height anomaly.}
1172     \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intp\_dza)) intp\_dza(i,j) = 0.5*dp**2 * &
1173           (alpha\_anom - c1\_90*(16.0*(a5(4)-a5(2)) + 7.0*(a5(5)-a5(1))) )
1174 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
1175 
1176   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intx\_dza)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=hi%jsc,hi%jec ; \textcolor{keywordflow}{do} i=isq,ieq
1177     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
1178     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation of}
1179     \textcolor{comment}{! T & S along the top and bottom integrals, akin to thickness weighting.}
1180     hwght = 0.0
1181     \textcolor{keywordflow}{if} (do\_massweight) &
1182       hwght = max(0., bathyp(i,j)-p\_t(i+1,j), bathyp(i+1,j)-p\_t(i,j))
1183     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
1184       hl = (p\_b(i,j) - p\_t(i,j)) + dp\_neglect
1185       hr = (p\_b(i+1,j) - p\_t(i+1,j)) + dp\_neglect
1186       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
1187       idenom = 1.0 / ( hwght*(hr + hl) + hl*hr )
1188       hwt\_ll = (hwght*hl + hr*hl) * idenom ; hwt\_lr = (hwght*hr) * idenom
1189       hwt\_rr = (hwght*hr + hr*hl) * idenom ; hwt\_rl = (hwght*hl) * idenom
1190     \textcolor{keywordflow}{else}
1191       hwt\_ll = 1.0 ; hwt\_lr = 0.0 ; hwt\_rr = 1.0 ; hwt\_rl = 0.0
1192 \textcolor{keywordflow}{    endif}
1193 
1194     intp(1) = dza(i,j) ; intp(5) = dza(i+1,j)
1195     \textcolor{keywordflow}{do} m=2,4
1196       wt\_l = 0.25*\textcolor{keywordtype}{real(5-m)} ; wt\_R = 1.0-wt\_l
1197       wtt\_l = wt\_l*hwt\_ll + wt\_r*hwt\_rl ; wtt\_r = wt\_l*hwt\_lr + wt\_r*hwt\_rr
1198 
1199       \textcolor{comment}{! T, S, and p are interpolated in the horizontal.  The p interpolation}
1200       \textcolor{comment}{! is linear, but for T and S it may be thickness weighted.}
1201       p5(1) = rl2\_t2\_to\_pa * (wt\_l*p\_b(i,j) + wt\_r*p\_b(i+1,j))
1202       dp = wt\_l*(p\_b(i,j) - p\_t(i,j)) + wt\_r*(p\_b(i+1,j) - p\_t(i+1,j))
1203       t5(1) = wtt\_l*t(i,j) + wtt\_r*t(i+1,j)
1204       s5(1) = wtt\_l*s(i,j) + wtt\_r*s(i+1,j)
1205 
1206       \textcolor{keywordflow}{do} n=2,5
1207         t5(n) = t5(1) ; s5(n) = s5(1) ; p5(n) = p5(n-1) - rl2\_t2\_to\_pa * 0.25*dp
1208 \textcolor{keywordflow}{      enddo}
1209       \textcolor{keywordflow}{if} (sv\_scale /= 1.0) \textcolor{keywordflow}{then}
1210         \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks, scale=sv\_scale)
1211       \textcolor{keywordflow}{else}
1212         \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks)
1213 \textcolor{keywordflow}{      endif}
1214 
1215     \textcolor{comment}{! Use Boole's rule to estimate the interface height anomaly change.}
1216       intp(m) = dp*( c1\_90*(7.0*(a5(1)+a5(5)) + 32.0*(a5(2)+a5(4)) + &
1217                                 12.0*a5(3)))
1218 \textcolor{keywordflow}{    enddo}
1219     \textcolor{comment}{! Use Boole's rule to integrate the interface height anomaly values in x.}
1220     intx\_dza(i,j) = c1\_90*(7.0*(intp(1)+intp(5)) + 32.0*(intp(2)+intp(4)) + &
1221                            12.0*intp(3))
1222 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
1223 
1224   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(inty\_dza)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=hi%isc,hi%iec
1225     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
1226     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation of}
1227     \textcolor{comment}{! T & S along the top and bottom integrals, akin to thickness weighting.}
1228     hwght = 0.0
1229     \textcolor{keywordflow}{if} (do\_massweight) &
1230       hwght = max(0., bathyp(i,j)-p\_t(i,j+1), bathyp(i,j+1)-p\_t(i,j))
1231     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
1232       hl = (p\_b(i,j) - p\_t(i,j)) + dp\_neglect
1233       hr = (p\_b(i,j+1) - p\_t(i,j+1)) + dp\_neglect
1234       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
1235       idenom = 1.0 / ( hwght*(hr + hl) + hl*hr )
1236       hwt\_ll = (hwght*hl + hr*hl) * idenom ; hwt\_lr = (hwght*hr) * idenom
1237       hwt\_rr = (hwght*hr + hr*hl) * idenom ; hwt\_rl = (hwght*hl) * idenom
1238     \textcolor{keywordflow}{else}
1239       hwt\_ll = 1.0 ; hwt\_lr = 0.0 ; hwt\_rr = 1.0 ; hwt\_rl = 0.0
1240 \textcolor{keywordflow}{    endif}
1241 
1242     intp(1) = dza(i,j) ; intp(5) = dza(i,j+1)
1243     \textcolor{keywordflow}{do} m=2,4
1244       wt\_l = 0.25*\textcolor{keywordtype}{real(5-m)} ; wt\_R = 1.0-wt\_l
1245       wtt\_l = wt\_l*hwt\_ll + wt\_r*hwt\_rl ; wtt\_r = wt\_l*hwt\_lr + wt\_r*hwt\_rr
1246 
1247       \textcolor{comment}{! T, S, and p are interpolated in the horizontal.  The p interpolation}
1248       \textcolor{comment}{! is linear, but for T and S it may be thickness weighted.}
1249       p5(1) = rl2\_t2\_to\_pa * (wt\_l*p\_b(i,j) + wt\_r*p\_b(i,j+1))
1250       dp = wt\_l*(p\_b(i,j) - p\_t(i,j)) + wt\_r*(p\_b(i,j+1) - p\_t(i,j+1))
1251       t5(1) = wtt\_l*t(i,j) + wtt\_r*t(i,j+1)
1252       s5(1) = wtt\_l*s(i,j) + wtt\_r*s(i,j+1)
1253       \textcolor{keywordflow}{do} n=2,5
1254         t5(n) = t5(1) ; s5(n) = s5(1) ; p5(n) = rl2\_t2\_to\_pa * (p5(n-1) - 0.25*dp)
1255 \textcolor{keywordflow}{      enddo}
1256       \textcolor{keywordflow}{if} (sv\_scale /= 1.0) \textcolor{keywordflow}{then}
1257         \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks, scale=sv\_scale)
1258       \textcolor{keywordflow}{else}
1259         \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks)
1260 \textcolor{keywordflow}{      endif}
1261 
1262     \textcolor{comment}{! Use Boole's rule to estimate the interface height anomaly change.}
1263       intp(m) = dp*( c1\_90*(7.0*(a5(1)+a5(5)) + 32.0*(a5(2)+a5(4)) + &
1264                                 12.0*a5(3)))
1265 \textcolor{keywordflow}{    enddo}
1266     \textcolor{comment}{! Use Boole's rule to integrate the interface height anomaly values in y.}
1267     inty\_dza(i,j) = c1\_90*(7.0*(intp(1)+intp(5)) + 32.0*(intp(2)+intp(4)) + &
1268                            12.0*intp(3))
1269 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
1270 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_a71a518f80b5f2ecf83fd9c638cad140d}\label{namespacemom__density__integrals_a71a518f80b5f2ecf83fd9c638cad140d}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+plm@{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+plm}}
\index{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+plm@{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+plm}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+plm()}{int\_spec\_vol\_dp\_generic\_plm()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::int\+\_\+spec\+\_\+vol\+\_\+dp\+\_\+generic\+\_\+plm (\begin{DoxyParamCaption}\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{T\+\_\+t,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{T\+\_\+b,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{S\+\_\+t,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{S\+\_\+b,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{p\+\_\+t,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{p\+\_\+b,  }\item[{real, intent(in)}]{alpha\+\_\+ref,  }\item[{real, intent(in)}]{d\+P\+\_\+neglect,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(in)}]{bathyP,  }\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout)}]{dza,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intp\+\_\+dza,  }\item[{real, dimension( hi \%isdb\+: hi \%iedb, hi \%jsd\+: hi \%jed), intent(inout), optional}]{intx\+\_\+dza,  }\item[{real, dimension( hi \%isd\+: hi \%ied, hi \%jsdb\+: hi \%jedb), intent(inout), optional}]{inty\+\_\+dza,  }\item[{logical, intent(in), optional}]{use\+Mass\+Wght\+Interp }\end{DoxyParamCaption})}



This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-\/volume form pressure accelerations in a non-\/\+Boussinesq model. There are essentially no free assumptions, apart from the use of Boole\textquotesingle{}s rule quadrature to do the integrals. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em hi} & A horizontal index type structure.\\
\hline
\mbox{\tt in}  & {\em t\+\_\+t} & Potential temperature at the top of the layer \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em t\+\_\+b} & Potential temperature at the bottom of the layer \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+t} & Salinity at the top the layer \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s\+\_\+b} & Salinity at the bottom the layer \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+t} & Pressure atop the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+b} & Pressure below the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em alpha\+\_\+ref} & A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals \mbox{[}R-\/1 $\sim$$>$ m3 kg-\/1\mbox{]} The calculation is mathematically identical with different values of alpha\+\_\+ref, but alpha\+\_\+ref alters the effects of roundoff, and answers do change.\\
\hline
\mbox{\tt in}  & {\em dp\+\_\+neglect} & !$<$ A miniscule pressure change with the same units as p\+\_\+t \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em bathyp} & The pressure at the bathymetry \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in,out}  & {\em dza} & The change in the geopotential anomaly\\
\hline
\mbox{\tt in,out}  & {\em intp\+\_\+dza} & The integral in pressure through the layer of\\
\hline
\mbox{\tt in,out}  & {\em intx\+\_\+dza} & The integral in x of the difference between\\
\hline
\mbox{\tt in,out}  & {\em inty\+\_\+dza} & The integral in y of the difference between\\
\hline
\mbox{\tt in}  & {\em usemasswghtinterp} & If true, uses mass weighting to interpolate T/S for top and bottom integrals. \\
\hline
\end{DoxyParams}


Definition at line 1280 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
1280   \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)}  :: HI\textcolor{comment}{ !< A horizontal index type structure.}
1281   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1282                         \textcolor{keywordtype}{intent(in)}  :: T\_t\textcolor{comment}{  !< Potential temperature at the top of the layer [degC]}
1283   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1284                         \textcolor{keywordtype}{intent(in)}  :: T\_b\textcolor{comment}{  !< Potential temperature at the bottom of the layer [degC]}
1285   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1286                         \textcolor{keywordtype}{intent(in)}  :: S\_t\textcolor{comment}{  !< Salinity at the top the layer [ppt]}
1287   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1288                         \textcolor{keywordtype}{intent(in)}  :: S\_b\textcolor{comment}{  !< Salinity at the bottom the layer [ppt]}
1289   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1290                         \textcolor{keywordtype}{intent(in)}  :: p\_t\textcolor{comment}{ !< Pressure atop the layer [R L2 T-2 ~> Pa] or [Pa]}
1291   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1292                         \textcolor{keywordtype}{intent(in)}  :: p\_b\textcolor{comment}{ !< Pressure below the layer [R L2 T-2 ~> Pa] or [Pa]}
1293   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: alpha\_ref\textcolor{comment}{ !< A mean specific volume that is subtracted out}
1294 \textcolor{comment}{                            !! to reduce the magnitude of each of the integrals [R-1 ~> m3 kg-1]}
1295 \textcolor{comment}{                            !! The calculation is mathematically identical with different values of}
1296 \textcolor{comment}{                            !! alpha\_ref, but alpha\_ref alters the effects of roundoff, and}
1297 \textcolor{comment}{                            !! answers do change.}
1298   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: dP\_neglect\textcolor{comment}{ !<!< A miniscule pressure change with}
1299 \textcolor{comment}{                                             !! the same units as p\_t [R L2 T-2 ~> Pa] or [Pa]}
1300   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1301                         \textcolor{keywordtype}{intent(in)}  :: bathyP\textcolor{comment}{ !< The pressure at the bathymetry [R L2 T-2 ~> Pa] or [Pa]}
1302   \textcolor{keywordtype}{type}(EOS\_type),       \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
1303   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type}
1304   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1305                         \textcolor{keywordtype}{intent(inout)} :: dza\textcolor{comment}{ !< The change in the geopotential anomaly}
1306 \textcolor{comment}{                            !! across the layer [L2 T-2 ~> m2 s-2]}
1307   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1308               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intp\_dza\textcolor{comment}{ !< The integral in pressure through the layer of}
1309 \textcolor{comment}{                            !! the geopotential anomaly relative to the anomaly at the bottom of the}
1310 \textcolor{comment}{                            !! layer [R L4 T-4 ~> Pa m2 s-2] or [Pa m2 s-2]}
1311   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(HI),SZJ\_(HI))}, &
1312               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intx\_dza\textcolor{comment}{  !< The integral in x of the difference between}
1313 \textcolor{comment}{                            !! the geopotential anomaly at the top and bottom of the layer divided}
1314 \textcolor{comment}{                            !! by the x grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1315   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJB\_(HI))}, &
1316               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: inty\_dza\textcolor{comment}{  !< The integral in y of the difference between}
1317 \textcolor{comment}{                            !! the geopotential anomaly at the top and bottom of the layer divided}
1318 \textcolor{comment}{                            !! by the y grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1319   \textcolor{keywordtype}{logical},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: useMassWghtInterp\textcolor{comment}{ !< If true, uses mass weighting}
1320 \textcolor{comment}{                            !! to interpolate T/S for top and bottom integrals.}
1321 
1322 \textcolor{comment}{!   This subroutine calculates analytical and nearly-analytical integrals in}
1323 \textcolor{comment}{! pressure across layers of geopotential anomalies, which are required for}
1324 \textcolor{comment}{! calculating the finite-volume form pressure accelerations in a non-Boussinesq}
1325 \textcolor{comment}{! model.  There are essentially no free assumptions, apart from the use of}
1326 \textcolor{comment}{! Boole's rule to do the horizontal integrals, and from a truncation in the}
1327 \textcolor{comment}{! series for log(1-eps/1+eps) that assumes that |eps| < 0.34.}
1328 
1329   \textcolor{keywordtype}{real} :: T5(5)      \textcolor{comment}{! Temperatures at five quadrature points [degC]}
1330   \textcolor{keywordtype}{real} :: S5(5)      \textcolor{comment}{! Salinities at five quadrature points [ppt]}
1331   \textcolor{keywordtype}{real} :: p5(5)      \textcolor{comment}{! Pressures at five quadrature points, scaled back to Pa as necessary [Pa]}
1332   \textcolor{keywordtype}{real} :: a5(5)      \textcolor{comment}{! Specific volumes at five quadrature points [R-1 ~> m3 kg-1] or [m3 kg-1]}
1333   \textcolor{keywordtype}{real} :: T15(15)    \textcolor{comment}{! Temperatures at fifteen interior quadrature points [degC]}
1334   \textcolor{keywordtype}{real} :: S15(15)    \textcolor{comment}{! Salinities at fifteen interior quadrature points [ppt]}
1335   \textcolor{keywordtype}{real} :: p15(15)    \textcolor{comment}{! Pressures at fifteen quadrature points, scaled back to Pa as necessary [Pa]}
1336   \textcolor{keywordtype}{real} :: a15(15)    \textcolor{comment}{! Specific volumes at fifteen quadrature points [R-1 ~> m3 kg-1] or [m3 kg-1]}
1337   \textcolor{keywordtype}{real} :: wt\_t(5), wt\_b(5) \textcolor{comment}{! Weights of top and bottom values at quadrature points [nondim]}
1338   \textcolor{keywordtype}{real} :: T\_top, T\_bot, S\_top, S\_bot, P\_top, P\_bot
1339 
1340   \textcolor{keywordtype}{real} :: alpha\_anom \textcolor{comment}{! The depth averaged specific density anomaly [m3 kg-1]}
1341   \textcolor{keywordtype}{real} :: dp         \textcolor{comment}{! The pressure change through a layer [R L2 T-2 ~> Pa]}
1342   \textcolor{keywordtype}{real} :: dp\_90(2:4) \textcolor{comment}{! The pressure change through a layer divided by 90 [R L2 T-2 ~> Pa]}
1343   \textcolor{keywordtype}{real} :: hWght      \textcolor{comment}{! A pressure-thickness below topography [R L2 T-2 ~> Pa]}
1344   \textcolor{keywordtype}{real} :: hL, hR     \textcolor{comment}{! Pressure-thicknesses of the columns to the left and right [R L2 T-2 ~> Pa]}
1345   \textcolor{keywordtype}{real} :: alpha\_ref\_mks \textcolor{comment}{! The reference specific volume in MKS units, never rescaled from m3 kg-1 [m3 kg-1]}
1346   \textcolor{keywordtype}{real} :: iDenom     \textcolor{comment}{! The inverse of the denominator in the weights [T4 R-2 L-4 ~> Pa-2]}
1347   \textcolor{keywordtype}{real} :: hWt\_LL, hWt\_LR \textcolor{comment}{! hWt\_LA is the weighted influence of A on the left column [nondim]}
1348   \textcolor{keywordtype}{real} :: hWt\_RL, hWt\_RR \textcolor{comment}{! hWt\_RA is the weighted influence of A on the right column [nondim]}
1349   \textcolor{keywordtype}{real} :: wt\_L, wt\_R \textcolor{comment}{! The linear weights of the left and right columns [nondim]}
1350   \textcolor{keywordtype}{real} :: wtT\_L, wtT\_R \textcolor{comment}{! The weights for tracers from the left and right columns [nondim]}
1351   \textcolor{keywordtype}{real} :: intp(5)    \textcolor{comment}{! The integrals of specific volume with pressure at the}
1352                      \textcolor{comment}{! 5 sub-column locations [L2 T-2 ~> m2 s-2]}
1353   \textcolor{keywordtype}{real} :: RL2\_T2\_to\_Pa  \textcolor{comment}{! A unit conversion factor from the rescaled units of pressure to Pa [Pa T2 R-1 L-2
       ~> 1]}
1354   \textcolor{keywordtype}{real} :: SV\_scale   \textcolor{comment}{! A multiplicative factor by which to scale specific}
1355                      \textcolor{comment}{! volume from m3 kg-1 to the desired units [kg m-3 R-1 ~> 1]}
1356   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_90 = 1.0/90.0  \textcolor{comment}{! A rational constant.}
1357   \textcolor{keywordtype}{logical} :: do\_massWeight \textcolor{comment}{! Indicates whether to do mass weighting.}
1358   \textcolor{keywordtype}{integer} :: Isq, Ieq, Jsq, Jeq, i, j, m, n, pos
1359 
1360   isq = hi%IscB ; ieq = hi%IecB ; jsq = hi%JscB ; jeq = hi%JecB
1361 
1362   do\_massweight = .false.
1363   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(usemasswghtinterp)) do\_massweight = usemasswghtinterp
1364 
1365   sv\_scale = us%R\_to\_kg\_m3
1366   rl2\_t2\_to\_pa = us%RL2\_T2\_to\_Pa
1367   alpha\_ref\_mks = alpha\_ref * us%kg\_m3\_to\_R
1368 
1369   \textcolor{keywordflow}{do} n = 1, 5 \textcolor{comment}{! Note that these are reversed from int\_density\_dz.}
1370     wt\_t(n) = 0.25 * \textcolor{keywordtype}{real}(n-1)
1371     wt\_b(n) = 1.0 - wt\_t(n)
1372 \textcolor{keywordflow}{  enddo}
1373 
1374   \textcolor{comment}{! 1. Compute vertical integrals}
1375   \textcolor{keywordflow}{do} j=jsq,jeq+1; \textcolor{keywordflow}{do} i=isq,ieq+1
1376     dp = p\_b(i,j) - p\_t(i,j)
1377     \textcolor{keywordflow}{do} n=1,5 \textcolor{comment}{! T, S and p are linearly interpolated in the vertical.}
1378       p5(n) = rl2\_t2\_to\_pa * (wt\_t(n) * p\_t(i,j) + wt\_b(n) * p\_b(i,j))
1379       s5(n) = wt\_t(n) * s\_t(i,j) + wt\_b(n) * s\_b(i,j)
1380       t5(n) = wt\_t(n) * t\_t(i,j) + wt\_b(n) * t\_b(i,j)
1381 \textcolor{keywordflow}{    enddo}
1382     \textcolor{keywordflow}{if} (sv\_scale /= 1.0) \textcolor{keywordflow}{then}
1383       \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks, scale=sv\_scale)
1384     \textcolor{keywordflow}{else}
1385       \textcolor{keyword}{call }calculate\_spec\_vol(t5, s5, p5, a5, 1, 5, eos, alpha\_ref\_mks)
1386 \textcolor{keywordflow}{    endif}
1387 
1388     \textcolor{comment}{! Use Boole's rule to estimate the interface height anomaly change.}
1389     alpha\_anom = c1\_90*((7.0*(a5(1)+a5(5)) + 32.0*(a5(2)+a5(4))) + 12.0*a5(3))
1390     dza(i,j) = dp*alpha\_anom
1391     \textcolor{comment}{! Use a Boole's-rule-like fifth-order accurate estimate of the double integral of}
1392     \textcolor{comment}{! the interface height anomaly.}
1393     \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intp\_dza)) intp\_dza(i,j) = 0.5*dp**2 * &
1394           (alpha\_anom - c1\_90*(16.0*(a5(4)-a5(2)) + 7.0*(a5(5)-a5(1))) )
1395 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
1396 
1397   \textcolor{comment}{! 2. Compute horizontal integrals in the x direction}
1398   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(intx\_dza)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=hi%jsc,hi%jec ; \textcolor{keywordflow}{do} i=isq,ieq
1399     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
1400     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation}
1401     \textcolor{comment}{! of T,S along the top and bottom integrals, almost like thickness}
1402     \textcolor{comment}{! weighting. Note: To work in terrain following coordinates we could}
1403     \textcolor{comment}{! offset this distance by the layer thickness to replicate other models.}
1404     hwght = 0.0
1405     \textcolor{keywordflow}{if} (do\_massweight) &
1406       hwght =  max(0., bathyp(i,j)-p\_t(i+1,j), bathyp(i+1,j)-p\_t(i,j))
1407     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
1408       hl = (p\_b(i,j) - p\_t(i,j)) + dp\_neglect
1409       hr = (p\_b(i+1,j) - p\_t(i+1,j)) + dp\_neglect
1410       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
1411       idenom = 1.0 / ( hwght*(hr + hl) + hl*hr )
1412       hwt\_ll = (hwght*hl + hr*hl) * idenom ; hwt\_lr = (hwght*hr) * idenom
1413       hwt\_rr = (hwght*hr + hr*hl) * idenom ; hwt\_rl = (hwght*hl) * idenom
1414     \textcolor{keywordflow}{else}
1415       hwt\_ll = 1.0 ; hwt\_lr = 0.0 ; hwt\_rr = 1.0 ; hwt\_rl = 0.0
1416 \textcolor{keywordflow}{    endif}
1417 
1418     \textcolor{keywordflow}{do} m=2,4
1419       wt\_l = 0.25*\textcolor{keywordtype}{real(5-m)} ; wt\_R = 1.0-wt\_l
1420       wtt\_l = wt\_l*hwt\_ll + wt\_r*hwt\_rl ; wtt\_r = wt\_l*hwt\_lr + wt\_r*hwt\_rr
1421 
1422       \textcolor{comment}{! T, S, and p are interpolated in the horizontal.  The p interpolation}
1423       \textcolor{comment}{! is linear, but for T and S it may be thickness weighted.}
1424       p\_top = wt\_l*p\_t(i,j) + wt\_r*p\_t(i+1,j)
1425       p\_bot = wt\_l*p\_b(i,j) + wt\_r*p\_b(i+1,j)
1426       t\_top = wtt\_l*t\_t(i,j) + wtt\_r*t\_t(i+1,j)
1427       t\_bot = wtt\_l*t\_b(i,j) + wtt\_r*t\_b(i+1,j)
1428       s\_top = wtt\_l*s\_t(i,j) + wtt\_r*s\_t(i+1,j)
1429       s\_bot = wtt\_l*s\_b(i,j) + wtt\_r*s\_b(i+1,j)
1430       dp\_90(m) = c1\_90*(p\_bot - p\_top)
1431 
1432       \textcolor{comment}{! Salinity, temperature and pressure with linear interpolation in the vertical.}
1433       pos = (m-2)*5
1434       \textcolor{keywordflow}{do} n=1,5
1435         p15(pos+n) = rl2\_t2\_to\_pa * (wt\_t(n) * p\_top + wt\_b(n) * p\_bot)
1436         s15(pos+n) = wt\_t(n) * s\_top + wt\_b(n) * s\_bot
1437         t15(pos+n) = wt\_t(n) * t\_top + wt\_b(n) * t\_bot
1438 \textcolor{keywordflow}{      enddo}
1439 \textcolor{keywordflow}{    enddo}
1440 
1441     \textcolor{keywordflow}{if} (sv\_scale /= 1.0) \textcolor{keywordflow}{then}
1442       \textcolor{keyword}{call }calculate\_spec\_vol(t15, s15, p15, a15, 1, 15, eos, alpha\_ref\_mks, scale=sv\_scale)
1443     \textcolor{keywordflow}{else}
1444       \textcolor{keyword}{call }calculate\_spec\_vol(t15, s15, p15, a15, 1, 15, eos, alpha\_ref\_mks)
1445 \textcolor{keywordflow}{    endif}
1446 
1447     intp(1) = dza(i,j) ; intp(5) = dza(i+1,j)
1448     \textcolor{keywordflow}{do} m=2,4
1449       \textcolor{comment}{! Use Boole's rule to estimate the interface height anomaly change.}
1450       \textcolor{comment}{! The integrals at the ends of the segment are already known.}
1451       pos = (m-2)*5
1452       intp(m) = dp\_90(m)*((7.0*(a15(pos+1)+a15(pos+5)) + &
1453                           32.0*(a15(pos+2)+a15(pos+4))) + 12.0*a15(pos+3))
1454 \textcolor{keywordflow}{    enddo}
1455     \textcolor{comment}{! Use Boole's rule to integrate the interface height anomaly values in x.}
1456     intx\_dza(i,j) = c1\_90*((7.0*(intp(1)+intp(5)) + 32.0*(intp(2)+intp(4))) + &
1457                            12.0*intp(3))
1458 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
1459 
1460   \textcolor{comment}{! 3. Compute horizontal integrals in the y direction}
1461   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(inty\_dza)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=hi%isc,hi%iec
1462     \textcolor{comment}{! hWght is the distance measure by which the cell is violation of}
1463     \textcolor{comment}{! hydrostatic consistency. For large hWght we bias the interpolation}
1464     \textcolor{comment}{! of T,S along the top and bottom integrals, like thickness weighting.}
1465     hwght = 0.0
1466     \textcolor{keywordflow}{if} (do\_massweight) &
1467       hwght = max(0., bathyp(i,j)-p\_t(i,j+1), bathyp(i,j+1)-p\_t(i,j))
1468     \textcolor{keywordflow}{if} (hwght > 0.) \textcolor{keywordflow}{then}
1469       hl = (p\_b(i,j) - p\_t(i,j)) + dp\_neglect
1470       hr = (p\_b(i,j+1) - p\_t(i,j+1)) + dp\_neglect
1471       hwght = hwght * ( (hl-hr)/(hl+hr) )**2
1472       idenom = 1.0 / ( hwght*(hr + hl) + hl*hr )
1473       hwt\_ll = (hwght*hl + hr*hl) * idenom ; hwt\_lr = (hwght*hr) * idenom
1474       hwt\_rr = (hwght*hr + hr*hl) * idenom ; hwt\_rl = (hwght*hl) * idenom
1475     \textcolor{keywordflow}{else}
1476       hwt\_ll = 1.0 ; hwt\_lr = 0.0 ; hwt\_rr = 1.0 ; hwt\_rl = 0.0
1477 \textcolor{keywordflow}{    endif}
1478 
1479     \textcolor{keywordflow}{do} m=2,4
1480       wt\_l = 0.25*\textcolor{keywordtype}{real(5-m)} ; wt\_R = 1.0-wt\_l
1481       wtt\_l = wt\_l*hwt\_ll + wt\_r*hwt\_rl ; wtt\_r = wt\_l*hwt\_lr + wt\_r*hwt\_rr
1482 
1483       \textcolor{comment}{! T, S, and p are interpolated in the horizontal.  The p interpolation}
1484       \textcolor{comment}{! is linear, but for T and S it may be thickness weighted.}
1485       p\_top = wt\_l*p\_t(i,j) + wt\_r*p\_t(i,j+1)
1486       p\_bot = wt\_l*p\_b(i,j) + wt\_r*p\_b(i,j+1)
1487       t\_top = wtt\_l*t\_t(i,j) + wtt\_r*t\_t(i,j+1)
1488       t\_bot = wtt\_l*t\_b(i,j) + wtt\_r*t\_b(i,j+1)
1489       s\_top = wtt\_l*s\_t(i,j) + wtt\_r*s\_t(i,j+1)
1490       s\_bot = wtt\_l*s\_b(i,j) + wtt\_r*s\_b(i,j+1)
1491       dp\_90(m) = c1\_90*(p\_bot - p\_top)
1492 
1493       \textcolor{comment}{! Salinity, temperature and pressure with linear interpolation in the vertical.}
1494       pos = (m-2)*5
1495       \textcolor{keywordflow}{do} n=1,5
1496         p15(pos+n) = rl2\_t2\_to\_pa * (wt\_t(n) * p\_top + wt\_b(n) * p\_bot)
1497         s15(pos+n) = wt\_t(n) * s\_top + wt\_b(n) * s\_bot
1498         t15(pos+n) = wt\_t(n) * t\_top + wt\_b(n) * t\_bot
1499 \textcolor{keywordflow}{      enddo}
1500 \textcolor{keywordflow}{    enddo}
1501 
1502     \textcolor{keywordflow}{if} (sv\_scale /= 1.0) \textcolor{keywordflow}{then}
1503       \textcolor{keyword}{call }calculate\_spec\_vol(t15, s15, p15, a15, 1, 15, eos, alpha\_ref\_mks, scale=sv\_scale)
1504     \textcolor{keywordflow}{else}
1505       \textcolor{keyword}{call }calculate\_spec\_vol(t15, s15, p15, a15, 1, 15, eos, alpha\_ref\_mks)
1506 \textcolor{keywordflow}{    endif}
1507 
1508     intp(1) = dza(i,j) ; intp(5) = dza(i,j+1)
1509     \textcolor{keywordflow}{do} m=2,4
1510       \textcolor{comment}{! Use Boole's rule to estimate the interface height anomaly change.}
1511       \textcolor{comment}{! The integrals at the ends of the segment are already known.}
1512       pos = (m-2)*5
1513       intp(m) = dp\_90(m) * ((7.0*(a15(pos+1)+a15(pos+5)) + &
1514                             32.0*(a15(pos+2)+a15(pos+4))) + 12.0*a15(pos+3))
1515 \textcolor{keywordflow}{    enddo}
1516     \textcolor{comment}{! Use Boole's rule to integrate the interface height anomaly values in x.}
1517     inty\_dza(i,j) = c1\_90*((7.0*(intp(1)+intp(5)) + 32.0*(intp(2)+intp(4))) + &
1518                            12.0*intp(3))
1519 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
1520 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__density__integrals_a759c2ae7aec17c59d532050f68a1e518}\label{namespacemom__density__integrals_a759c2ae7aec17c59d532050f68a1e518}} 
\index{mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}!int\+\_\+specific\+\_\+vol\+\_\+dp@{int\+\_\+specific\+\_\+vol\+\_\+dp}}
\index{int\+\_\+specific\+\_\+vol\+\_\+dp@{int\+\_\+specific\+\_\+vol\+\_\+dp}!mom\+\_\+density\+\_\+integrals@{mom\+\_\+density\+\_\+integrals}}
\subsubsection{\texorpdfstring{int\+\_\+specific\+\_\+vol\+\_\+dp()}{int\_specific\_vol\_dp()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+density\+\_\+integrals\+::int\+\_\+specific\+\_\+vol\+\_\+dp (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{T,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{S,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{p\+\_\+t,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in)}]{p\+\_\+b,  }\item[{real, intent(in)}]{alpha\+\_\+ref,  }\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI,  }\item[{type(eos\+\_\+type), pointer}]{E\+OS,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(inout)}]{dza,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(inout), optional}]{intp\+\_\+dza,  }\item[{real, dimension(szib\+\_\+(hi),szj\+\_\+(hi)), intent(inout), optional}]{intx\+\_\+dza,  }\item[{real, dimension(szi\+\_\+(hi),szjb\+\_\+(hi)), intent(inout), optional}]{inty\+\_\+dza,  }\item[{integer, intent(in), optional}]{halo\+\_\+size,  }\item[{real, dimension(szi\+\_\+(hi),szj\+\_\+(hi)), intent(in), optional}]{bathyP,  }\item[{real, intent(in), optional}]{d\+P\+\_\+tiny,  }\item[{logical, intent(in), optional}]{use\+Mass\+Wght\+Interp }\end{DoxyParamCaption})}



Calls the appropriate subroutine to calculate analytical and nearly-\/analytical integrals in pressure across layers of geopotential anomalies, which are required for calculating the finite-\/volume form pressure accelerations in a non-\/\+Boussinesq model. There are essentially no free assumptions, apart from the use of Boole\textquotesingle{}s rule to do the horizontal integrals, and from a truncation in the series for log(1-\/eps/1+eps) that assumes that $\vert$eps$\vert$ $<$ 0.\+34. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em hi} & The horizontal index structure\\
\hline
\mbox{\tt in}  & {\em t} & Potential temperature referenced to the surface \mbox{[}degC\mbox{]}\\
\hline
\mbox{\tt in}  & {\em s} & Salinity \mbox{[}ppt\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+t} & Pressure at the top of the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em p\+\_\+b} & Pressure at the bottom of the layer \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em alpha\+\_\+ref} & A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals \mbox{[}R-\/1 $\sim$$>$ m3 kg-\/1\mbox{]} The calculation is mathematically identical with different values of alpha\+\_\+ref, but this reduces the effects of roundoff.\\
\hline
 & {\em eos} & Equation of state structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in,out}  & {\em dza} & The change in the geopotential anomaly across\\
\hline
\mbox{\tt in,out}  & {\em intp\+\_\+dza} & The integral in pressure through the layer of the\\
\hline
\mbox{\tt in,out}  & {\em intx\+\_\+dza} & The integral in x of the difference between the\\
\hline
\mbox{\tt in,out}  & {\em inty\+\_\+dza} & The integral in y of the difference between the\\
\hline
\mbox{\tt in}  & {\em halo\+\_\+size} & The width of halo points on which to calculate dza.\\
\hline
\mbox{\tt in}  & {\em bathyp} & The pressure at the bathymetry \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em dp\+\_\+tiny} & A minuscule pressure change with the same units as p\+\_\+t \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} or \mbox{[}Pa\mbox{]}\\
\hline
\mbox{\tt in}  & {\em usemasswghtinterp} & If true, uses mass weighting to interpolate T/S for top and bottom integrals. \\
\hline
\end{DoxyParams}


Definition at line 1007 of file M\+O\+M\+\_\+density\+\_\+integrals.\+F90.


\begin{DoxyCode}
1007   \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)}  :: HI\textcolor{comment}{  !< The horizontal index structure}
1008   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1009                         \textcolor{keywordtype}{intent(in)}  :: T\textcolor{comment}{   !< Potential temperature referenced to the surface [degC]}
1010   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1011                         \textcolor{keywordtype}{intent(in)}  :: S\textcolor{comment}{   !< Salinity [ppt]}
1012   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1013                         \textcolor{keywordtype}{intent(in)}  :: p\_t\textcolor{comment}{ !< Pressure at the top of the layer [R L2 T-2 ~> Pa] or [Pa]}
1014   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1015                         \textcolor{keywordtype}{intent(in)}  :: p\_b\textcolor{comment}{ !< Pressure at the bottom of the layer [R L2 T-2 ~> Pa] or [Pa]}
1016   \textcolor{keywordtype}{real},                 \textcolor{keywordtype}{intent(in)}  :: alpha\_ref\textcolor{comment}{ !< A mean specific volume that is subtracted out}
1017 \textcolor{comment}{                            !! to reduce the magnitude of each of the integrals [R-1 ~> m3 kg-1]}
1018 \textcolor{comment}{                            !! The calculation is mathematically identical with different values of}
1019 \textcolor{comment}{                            !! alpha\_ref, but this reduces the effects of roundoff.}
1020   \textcolor{keywordtype}{type}(EOS\_type),       \textcolor{keywordtype}{pointer}     :: EOS\textcolor{comment}{ !< Equation of state structure}
1021   \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{  !< A dimensional unit scaling type}
1022   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1023                         \textcolor{keywordtype}{intent(inout)} :: dza\textcolor{comment}{ !< The change in the geopotential anomaly across}
1024 \textcolor{comment}{                            !! the layer [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1025   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1026               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intp\_dza\textcolor{comment}{ !< The integral in pressure through the layer of the}
1027 \textcolor{comment}{                            !! geopotential anomaly relative to the anomaly at the bottom of the}
1028 \textcolor{comment}{                            !! layer [R L4 T-4 ~> Pa m2 s-2] or [Pa m2 s-2]}
1029   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(HI),SZJ\_(HI))}, &
1030               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: intx\_dza\textcolor{comment}{ !< The integral in x of the difference between the}
1031 \textcolor{comment}{                            !! geopotential anomaly at the top and bottom of the layer divided by}
1032 \textcolor{comment}{                            !! the x grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1033   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJB\_(HI))}, &
1034               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: inty\_dza\textcolor{comment}{ !< The integral in y of the difference between the}
1035 \textcolor{comment}{                            !! geopotential anomaly at the top and bottom of the layer divided by}
1036 \textcolor{comment}{                            !! the y grid spacing [L2 T-2 ~> m2 s-2] or [m2 s-2]}
1037   \textcolor{keywordtype}{integer},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: halo\_size\textcolor{comment}{ !< The width of halo points on which to calculate dza.}
1038   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(HI),SZJ\_(HI))}, &
1039               \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: bathyP\textcolor{comment}{  !< The pressure at the bathymetry [R L2 T-2 ~> Pa] or [Pa]}
1040   \textcolor{keywordtype}{real},       \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: dP\_tiny\textcolor{comment}{ !< A minuscule pressure change with}
1041 \textcolor{comment}{                            !! the same units as p\_t [R L2 T-2 ~> Pa] or [Pa]}
1042   \textcolor{keywordtype}{logical},    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}  :: useMassWghtInterp\textcolor{comment}{ !< If true, uses mass weighting}
1043 \textcolor{comment}{                            !! to interpolate T/S for top and bottom integrals.}
1044 
1045   \textcolor{keywordflow}{if} (eos\_quadrature(eos)) \textcolor{keywordflow}{then}
1046     \textcolor{keyword}{call }int\_spec\_vol\_dp\_generic\_pcm(t, s, p\_t, p\_b, alpha\_ref, hi, eos, us, &
1047                                      dza, intp\_dza, intx\_dza, inty\_dza, halo\_size, &
1048                                      bathyp, dp\_tiny, usemasswghtinterp)
1049   \textcolor{keywordflow}{else}
1050     \textcolor{keyword}{call }analytic\_int\_specific\_vol\_dp(t, s, p\_t, p\_b, alpha\_ref, hi, eos, &
1051                                       dza, intp\_dza, intx\_dza, inty\_dza, halo\_size, &
1052                                       bathyp, dp\_tiny, usemasswghtinterp)
1053 \textcolor{keywordflow}{  endif}
1054 
\end{DoxyCode}
