\hypertarget{namespacemom__set__diffusivity}{}\section{mom\+\_\+set\+\_\+diffusivity Module Reference}
\label{namespacemom__set__diffusivity}\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}


\subsection{Detailed Description}
Calculate vertical diffusivity from all mixing processes. \subsection*{Data Types}
\begin{DoxyCompactItemize}
\item 
type \hyperlink{structmom__set__diffusivity_1_1diffusivity__diags}{diffusivity\+\_\+diags}
\begin{DoxyCompactList}\small\item\em This structure has memory for used in calculating diagnostics of diffusivity. \end{DoxyCompactList}\item 
type \hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}
\begin{DoxyCompactList}\small\item\em This control structure contains parameters for M\+O\+M\+\_\+set\+\_\+diffusivity. \end{DoxyCompactList}\end{DoxyCompactItemize}
\subsection*{Functions/\+Subroutines}
\begin{DoxyCompactItemize}
\item 
subroutine, public \hyperlink{namespacemom__set__diffusivity_a7c293162d6c8efb882c8b04b4ea5241d}{set\+\_\+diffusivity} (u, v, h, u\+\_\+h, v\+\_\+h, tv, fluxes, optics, visc, dt, G, GV, US, CS, Kd\+\_\+lay, Kd\+\_\+int)
\begin{DoxyCompactList}\small\item\em Sets the interior vertical diffusion of scalars due to the following processes\+: \end{DoxyCompactList}\item 
subroutine \hyperlink{namespacemom__set__diffusivity_a07c0ab3f141f8f9e057be3150a940a94}{find\+\_\+tke\+\_\+to\+\_\+kd} (h, tv, d\+Rho\+\_\+int, N2\+\_\+lay, j, dt, G, GV, US, CS, T\+K\+E\+\_\+to\+\_\+\+Kd, max\+T\+KE, kb)
\begin{DoxyCompactList}\small\item\em Convert turbulent kinetic energy to diffusivity. \end{DoxyCompactList}\item 
subroutine \hyperlink{namespacemom__set__diffusivity_afef80c2221be24f63f1ca96c2abe6fa9}{find\+\_\+n2} (h, tv, T\+\_\+f, S\+\_\+f, fluxes, j, G, GV, US, CS, d\+Rho\+\_\+int, N2\+\_\+lay, N2\+\_\+int, N2\+\_\+bot)
\begin{DoxyCompactList}\small\item\em Calculate Brunt-\/\+Vaisala frequency, N$^\wedge$2. \end{DoxyCompactList}\item 
subroutine \hyperlink{namespacemom__set__diffusivity_a99d0eb7701f8e04d856b75117fe7b83c}{double\+\_\+diffusion} (tv, h, T\+\_\+f, S\+\_\+f, j, G, GV, US, CS, Kd\+\_\+\+T\+\_\+dd, Kd\+\_\+\+S\+\_\+dd)
\begin{DoxyCompactList}\small\item\em This subroutine sets the additional diffusivities of temperature and salinity due to double diffusion, using the same functional form as is used in M\+O\+M4.\+1, and taken from an N\+C\+AR technical note (R\+EF?) that updates what was in Large et al. (1994). All the coefficients here should probably be made run-\/time variables rather than hard-\/coded constants. \end{DoxyCompactList}\item 
subroutine \hyperlink{namespacemom__set__diffusivity_ac48033315c2bbdb5551b272a235c16e5}{add\+\_\+drag\+\_\+diffusivity} (h, u, v, tv, fluxes, visc, j, T\+K\+E\+\_\+to\+\_\+\+Kd, max\+T\+KE, kb, G, GV, US, CS, Kd\+\_\+lay, Kd\+\_\+int, Kd\+\_\+\+B\+BL)
\begin{DoxyCompactList}\small\item\em This routine adds diffusion sustained by flow energy extracted by bottom drag. \end{DoxyCompactList}\item 
subroutine \hyperlink{namespacemom__set__diffusivity_a968690a1c9efb859ee67965f9032c7d7}{add\+\_\+lotw\+\_\+bbl\+\_\+diffusivity} (h, u, v, tv, fluxes, visc, j, N2\+\_\+int, G, GV, US, CS, Kd\+\_\+\+B\+BL, Kd\+\_\+lay, Kd\+\_\+int)
\begin{DoxyCompactList}\small\item\em Calculates a B\+BL diffusivity use a Prandtl number 1 diffusivity with a law of the wall turbulent viscosity, up to a B\+BL height where the energy used for mixing has consumed the mechanical T\+KE input. \end{DoxyCompactList}\item 
subroutine \hyperlink{namespacemom__set__diffusivity_a21162cb5173d52ce92dc65beee9d0ef1}{add\+\_\+mlrad\+\_\+diffusivity} (h, fluxes, j, G, GV, US, CS, T\+K\+E\+\_\+to\+\_\+\+Kd, Kd\+\_\+lay, Kd\+\_\+int)
\begin{DoxyCompactList}\small\item\em This routine adds effects of mixed layer radiation to the layer diffusivities. \end{DoxyCompactList}\item 
subroutine, public \hyperlink{namespacemom__set__diffusivity_a66f77b7e2f9c0c8254da4a0acd5a9996}{set\+\_\+bbl\+\_\+tke} (u, v, h, fluxes, visc, G, GV, US, CS, O\+BC)
\begin{DoxyCompactList}\small\item\em This subroutine calculates several properties related to bottom boundary layer turbulence. \end{DoxyCompactList}\item 
subroutine \hyperlink{namespacemom__set__diffusivity_a5ba8a3be6234304aa5f1dfd0b831078a}{set\+\_\+density\+\_\+ratios} (h, tv, kb, G, GV, US, CS, j, ds\+\_\+dsp1, rho\+\_\+0)
\item 
subroutine, public \hyperlink{namespacemom__set__diffusivity_a2296f20ef1f3a4c2390897a7376f88e7}{set\+\_\+diffusivity\+\_\+init} (Time, G, GV, US, param\+\_\+file, diag, CS, int\+\_\+tide\+\_\+\+C\+Sp, halo\+\_\+\+TS)
\item 
subroutine, public \hyperlink{namespacemom__set__diffusivity_ace82f133d3cee42aa36ec10bcce79e75}{set\+\_\+diffusivity\+\_\+end} (CS)
\begin{DoxyCompactList}\small\item\em Clear pointers and dealocate memory. \end{DoxyCompactList}\end{DoxyCompactItemize}
\subsection*{Variables}
\textbf{ }\par
\begin{DoxyCompactItemize}
\item 
\mbox{\Hypertarget{namespacemom__set__diffusivity_a6b2144f31cb69fe2481815545253d7c9}\label{namespacemom__set__diffusivity_a6b2144f31cb69fe2481815545253d7c9}} 
integer \hyperlink{namespacemom__set__diffusivity_a6b2144f31cb69fe2481815545253d7c9}{id\+\_\+clock\+\_\+kappashear}
\begin{DoxyCompactList}\small\item\em C\+PU time clocks. \end{DoxyCompactList}\item 
\mbox{\Hypertarget{namespacemom__set__diffusivity_a7ceffee2f68b88087dafa02f6235c4fc}\label{namespacemom__set__diffusivity_a7ceffee2f68b88087dafa02f6235c4fc}} 
integer \hyperlink{namespacemom__set__diffusivity_a7ceffee2f68b88087dafa02f6235c4fc}{id\+\_\+clock\+\_\+cvmix\+\_\+ddiff}
\begin{DoxyCompactList}\small\item\em C\+PU time clocks. \end{DoxyCompactList}\end{DoxyCompactItemize}



\subsection{Function/\+Subroutine Documentation}
\mbox{\Hypertarget{namespacemom__set__diffusivity_ac48033315c2bbdb5551b272a235c16e5}\label{namespacemom__set__diffusivity_ac48033315c2bbdb5551b272a235c16e5}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!add\+\_\+drag\+\_\+diffusivity@{add\+\_\+drag\+\_\+diffusivity}}
\index{add\+\_\+drag\+\_\+diffusivity@{add\+\_\+drag\+\_\+diffusivity}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{add\+\_\+drag\+\_\+diffusivity()}{add\_drag\_diffusivity()}}
{\footnotesize\ttfamily subroutine mom\+\_\+set\+\_\+diffusivity\+::add\+\_\+drag\+\_\+diffusivity (\begin{DoxyParamCaption}\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h,  }\item[{real, dimension( g \%isdb\+: g \%iedb, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{u,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsdb\+: g \%jedb, g \%ke), intent(in)}]{v,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{type(forcing), intent(in)}]{fluxes,  }\item[{type(vertvisc\+\_\+type), intent(in)}]{visc,  }\item[{integer, intent(in)}]{j,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(in)}]{T\+K\+E\+\_\+to\+\_\+\+Kd,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(in)}]{max\+T\+KE,  }\item[{integer, dimension( g \%isd\+: g \%ied), intent(in)}]{kb,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(inout)}]{Kd\+\_\+lay,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(inout), optional}]{Kd\+\_\+int,  }\item[{real, dimension(\+:,\+:,\+:), pointer}]{Kd\+\_\+\+B\+BL }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



This routine adds diffusion sustained by flow energy extracted by bottom drag. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em u} & The zonal velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em v} & The meridional velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em tv} & Structure containing pointers to any available thermodynamic fields.\\
\hline
\mbox{\tt in}  & {\em fluxes} & A structure of thermodynamic surface fluxes\\
\hline
\mbox{\tt in}  & {\em visc} & Structure containing vertical viscosities, bottom boundary layer properies, and related fields\\
\hline
\mbox{\tt in}  & {\em j} & j-\/index of row to work on\\
\hline
\mbox{\tt in}  & {\em tke\+\_\+to\+\_\+kd} & The conversion rate between the T\+KE T\+KE dissipated within a layer and the diapycnal diffusivity witin that layer, usually ($\sim$\+Rho\+\_\+0 / (G\+\_\+\+Earth $\ast$ d\+Rho\+\_\+lay)) \mbox{[}Z2 T-\/1 / Z3 T-\/3 = T2 Z-\/1 $\sim$$>$ s2 m-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em maxtke} & The energy required to for a layer to entrain to its maximum-\/realizable thickness \mbox{[}Z3 T-\/3 $\sim$$>$ m3 s-\/3\mbox{]}\\
\hline
\mbox{\tt in}  & {\em kb} & Index of lightest layer denser than the buffer layer, or -\/1 without a bulk mixed layer\\
\hline
 & {\em cs} & Diffusivity control structure\\
\hline
\mbox{\tt in,out}  & {\em kd\+\_\+lay} & The diapycnal diffusivity in layers, \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\
\hline
\mbox{\tt in,out}  & {\em kd\+\_\+int} & The diapycnal diffusivity at interfaces,\\
\hline
 & {\em kd\+\_\+bbl} & Interface B\+BL diffusivity \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}. \\
\hline
\end{DoxyParams}


Definition at line 1165 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
1165   \textcolor{keywordtype}{type}(ocean\_grid\_type),            \textcolor{keywordtype}{intent(in)}    :: g\textcolor{comment}{    !< The ocean's grid structure}
1166   \textcolor{keywordtype}{type}(verticalgrid\_type),          \textcolor{keywordtype}{intent(in)}    :: gv\textcolor{comment}{   !< The ocean's vertical grid structure}
1167   \textcolor{keywordtype}{type}(unit\_scale\_type),            \textcolor{keywordtype}{intent(in)}    :: us\textcolor{comment}{   !< A dimensional unit scaling type}
1168   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, &
1169                                     \textcolor{keywordtype}{intent(in)}    :: u\textcolor{comment}{    !< The zonal velocity [L T-1 ~> m s-1]}
1170   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, &
1171                                     \textcolor{keywordtype}{intent(in)}    :: v\textcolor{comment}{    !< The meridional velocity [L T-1 ~> m s-1]}
1172   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1173                                     \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{    !< Layer thicknesses [H ~> m or kg m-2]}
1174   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),            \textcolor{keywordtype}{intent(in)}    :: tv\textcolor{comment}{   !< Structure containing pointers to any available}
1175 \textcolor{comment}{                                                          !! thermodynamic fields.}
1176   \textcolor{keywordtype}{type}(forcing),                    \textcolor{keywordtype}{intent(in)}    :: fluxes\textcolor{comment}{ !< A structure of thermodynamic surface fluxes}
1177   \textcolor{keywordtype}{type}(vertvisc\_type),              \textcolor{keywordtype}{intent(in)}    :: visc\textcolor{comment}{ !< Structure containing vertical viscosities,
       bottom}
1178 \textcolor{comment}{                                                          !! boundary layer properies, and related fields}
1179   \textcolor{keywordtype}{integer},                          \textcolor{keywordtype}{intent(in)}    :: j\textcolor{comment}{    !< j-index of row to work on}
1180   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)}    :: tke\_to\_kd\textcolor{comment}{ !< The conversion rate between the TKE}
1181 \textcolor{comment}{                                                          !! TKE dissipated within  a layer and the}
1182 \textcolor{comment}{                                                          !! diapycnal diffusivity witin that layer,}
1183 \textcolor{comment}{                                                          !! usually (~Rho\_0 / (G\_Earth * dRho\_lay))}
1184 \textcolor{comment}{                                                          !! [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]}
1185   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)}    :: maxtke\textcolor{comment}{ !< The energy required to for a layer to
       entrain}
1186 \textcolor{comment}{                                                          !! to its maximum-realizable thickness [Z3 T-3 ~>
       m3 s-3]}
1187   \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(SZI\_(G))},      \textcolor{keywordtype}{intent(in)}    :: kb\textcolor{comment}{   !< Index of lightest layer denser than the buffer}
1188 \textcolor{comment}{                                                          !! layer, or -1 without a bulk mixed layer}
1189   \textcolor{keywordtype}{type}(set\_diffusivity\_cs),         \textcolor{keywordtype}{pointer}       :: cs\textcolor{comment}{   !< Diffusivity control structure}
1190   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(inout)} :: kd\_lay\textcolor{comment}{ !< The diapycnal diffusivity in layers,}
1191 \textcolor{comment}{                                                            !! [Z2 T-1 ~> m2 s-1].}
1192   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)}, &
1193                           \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: kd\_int\textcolor{comment}{ !< The diapycnal diffusivity at interfaces,}
1194 \textcolor{comment}{                                                            !! [Z2 T-1 ~> m2 s-1].}
1195   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:,:)},           \textcolor{keywordtype}{pointer}       :: kd\_bbl\textcolor{comment}{ !< Interface BBL diffusivity [Z2 T-1 ~> m2
       s-1].}
1196 
1197 \textcolor{comment}{! This routine adds diffusion sustained by flow energy extracted by bottom drag.}
1198 
1199   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(G)+1)} :: &
1200     rint          \textcolor{comment}{! coordinate density of an interface [R ~> kg m-3]}
1201   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: &
1202     htot, &       \textcolor{comment}{! total thickness above or below a layer, or the}
1203                   \textcolor{comment}{! integrated thickness in the BBL [Z ~> m].}
1204     rho\_htot, &   \textcolor{comment}{! running integral with depth of density [Z R ~> kg m-2]}
1205     gh\_sum\_top, & \textcolor{comment}{! BBL value of g'h that can be supported by}
1206                   \textcolor{comment}{! the local ustar, times R0\_g [R ~> kg m-2]}
1207     rho\_top, &    \textcolor{comment}{! density at top of the BBL [R ~> kg m-3]}
1208     tke, &        \textcolor{comment}{! turbulent kinetic energy available to drive}
1209                   \textcolor{comment}{! bottom-boundary layer mixing in a layer [Z3 T-3 ~> m3 s-3]}
1210     i2decay       \textcolor{comment}{! inverse of twice the TKE decay scale [Z-1 ~> m-1].}
1211 
1212   \textcolor{keywordtype}{real}    :: tke\_to\_layer   \textcolor{comment}{! TKE used to drive mixing in a layer [Z3 T-3 ~> m3 s-3]}
1213   \textcolor{keywordtype}{real}    :: tke\_ray        \textcolor{comment}{! TKE from layer Rayleigh drag used to drive mixing in layer [Z3 T-3 ~> m3 s-3]}
1214   \textcolor{keywordtype}{real}    :: tke\_here       \textcolor{comment}{! TKE that goes into mixing in this layer [Z3 T-3 ~> m3 s-3]}
1215   \textcolor{keywordtype}{real}    :: drl, drbot     \textcolor{comment}{! temporaries holding density differences [R ~> kg m-3]}
1216   \textcolor{keywordtype}{real}    :: cdrag\_sqrt     \textcolor{comment}{! square root of the drag coefficient [nondim]}
1217   \textcolor{keywordtype}{real}    :: ustar\_h        \textcolor{comment}{! value of ustar at a thickness point [Z T-1 ~> m s-1].}
1218   \textcolor{keywordtype}{real}    :: absf           \textcolor{comment}{! average absolute Coriolis parameter around a thickness point [T-1 ~> s-1]}
1219   \textcolor{keywordtype}{real}    :: r0\_g           \textcolor{comment}{! Rho0 / G\_Earth [R T2 Z-1 m-1 ~> kg s2 m-5]}
1220   \textcolor{keywordtype}{real}    :: i\_rho0         \textcolor{comment}{! 1 / RHO0 [R-1 ~> m3 kg-1]}
1221   \textcolor{keywordtype}{real}    :: delta\_kd       \textcolor{comment}{! increment to Kd from the bottom boundary layer mixing [Z2 T-1 ~> m2 s-1].}
1222   \textcolor{keywordtype}{logical} :: rayleigh\_drag  \textcolor{comment}{! Set to true if Rayleigh drag velocities}
1223                             \textcolor{comment}{! defined in visc, on the assumption that this}
1224                             \textcolor{comment}{! extracted energy also drives diapycnal mixing.}
1225 
1226   \textcolor{keywordtype}{logical} :: domore, do\_i(szi\_(g))
1227   \textcolor{keywordtype}{logical} :: do\_diag\_kd\_bbl
1228 
1229   \textcolor{keywordtype}{integer} :: i, k, is, ie, nz, i\_rem, kb\_min
1230   is = g%isc ; ie = g%iec ; nz = g%ke
1231 
1232   do\_diag\_kd\_bbl = \textcolor{keyword}{associated}(kd\_bbl)
1233 
1234   \textcolor{keywordflow}{if} (.not.(cs%bottomdraglaw .and. (cs%BBL\_effic>0.0))) \textcolor{keywordflow}{return}
1235 
1236   cdrag\_sqrt = sqrt(cs%cdrag)
1237   tke\_ray = 0.0 ; rayleigh\_drag = .false.
1238   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Ray\_u) .and. \textcolor{keyword}{associated}(visc%Ray\_v)) rayleigh\_drag = .true.
1239 
1240   i\_rho0 = 1.0 / (gv%Rho0)
1241   r0\_g = gv%Rho0 / (us%L\_to\_Z**2 * gv%g\_Earth)
1242 
1243   \textcolor{keywordflow}{do} k=2,nz ; rint(k) = 0.5*(gv%Rlay(k-1)+gv%Rlay(k)) ;\textcolor{keywordflow}{ enddo}
1244 
1245   kb\_min = max(gv%nk\_rho\_varies+1,2)
1246 
1247   \textcolor{comment}{! The turbulence decay scale is 0.5*ustar/f from K&E & MOM\_vertvisc.F90}
1248   \textcolor{comment}{! Any turbulence that makes it into the mixed layers is assumed}
1249   \textcolor{comment}{! to be relatively small and is discarded.}
1250   \textcolor{keywordflow}{do} i=is,ie
1251     ustar\_h = visc%ustar\_BBL(i,j)
1252     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(fluxes%ustar\_tidal)) &
1253       ustar\_h = ustar\_h + fluxes%ustar\_tidal(i,j)
1254     absf = 0.25 * ((abs(g%CoriolisBu(i-1,j-1)) + abs(g%CoriolisBu(i,j))) + &
1255                    (abs(g%CoriolisBu(i-1,j)) + abs(g%CoriolisBu(i,j-1))))
1256     \textcolor{keywordflow}{if} ((ustar\_h > 0.0) .and. (absf > 0.5*cs%IMax\_decay*ustar\_h))  \textcolor{keywordflow}{then}
1257       i2decay(i) = absf / ustar\_h
1258     \textcolor{keywordflow}{else}
1259       \textcolor{comment}{! The maximum decay scale should be something of order 200 m.}
1260       \textcolor{comment}{! If ustar\_h = 0, this is land so this value doesn't matter.}
1261       i2decay(i) = 0.5*cs%IMax\_decay
1262 \textcolor{keywordflow}{    endif}
1263     tke(i) = ((cs%BBL\_effic * cdrag\_sqrt) * exp(-i2decay(i)*(gv%H\_to\_Z*h(i,j,nz))) ) * &
1264              visc%TKE\_BBL(i,j)
1265 
1266     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(fluxes%TKE\_tidal)) &
1267       tke(i) = tke(i) + fluxes%TKE\_tidal(i,j) * i\_rho0 * &
1268            (cs%BBL\_effic * exp(-i2decay(i)*(gv%H\_to\_Z*h(i,j,nz))))
1269 
1270     \textcolor{comment}{! Distribute the work over a BBL of depth 20^2 ustar^2 / g' following}
1271     \textcolor{comment}{! Killworth & Edwards (1999) and Zilitikevich & Mironov (1996).}
1272     \textcolor{comment}{! Rho\_top is determined by finding the density where}
1273     \textcolor{comment}{! integral(bottom, Z) (rho(z') - rho(Z)) dz' = rho\_0 400 ustar^2 / g}
1274 
1275     gh\_sum\_top(i) = r0\_g * 400.0 * ustar\_h**2
1276 
1277     do\_i(i) = (g%mask2dT(i,j) > 0.5)
1278     htot(i) = gv%H\_to\_Z*h(i,j,nz)
1279     rho\_htot(i) = gv%Rlay(nz)*(gv%H\_to\_Z*h(i,j,nz))
1280     rho\_top(i) = gv%Rlay(1)
1281     \textcolor{keywordflow}{if} (cs%bulkmixedlayer .and. do\_i(i)) rho\_top(i) = gv%Rlay(kb(i)-1)
1282 \textcolor{keywordflow}{  enddo}
1283 
1284   \textcolor{keywordflow}{do} k=nz-1,2,-1 ; domore = .false.
1285     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1286       htot(i) = htot(i) + gv%H\_to\_Z*h(i,j,k)
1287       rho\_htot(i) = rho\_htot(i) + gv%Rlay(k)*(gv%H\_to\_Z*h(i,j,k))
1288       \textcolor{keywordflow}{if} (htot(i)*gv%Rlay(k-1) <= (rho\_htot(i) - gh\_sum\_top(i))) \textcolor{keywordflow}{then}
1289         \textcolor{comment}{! The top of the mixing is in the interface atop the current layer.}
1290         rho\_top(i) = (rho\_htot(i) - gh\_sum\_top(i)) / htot(i)
1291         do\_i(i) = .false.
1292       \textcolor{keywordflow}{elseif} (k <= kb(i)) \textcolor{keywordflow}{then} ; do\_i(i) = .false.
1293       \textcolor{keywordflow}{else} ; domore = .true. ;\textcolor{keywordflow}{ endif}
1294 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1295     \textcolor{keywordflow}{if} (.not.domore) \textcolor{keywordflow}{exit}
1296 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! k-loop}
1297 
1298   \textcolor{keywordflow}{do} i=is,ie ; do\_i(i) = (g%mask2dT(i,j) > 0.5) ;\textcolor{keywordflow}{ enddo}
1299   \textcolor{keywordflow}{do} k=nz-1,kb\_min,-1
1300     i\_rem = 0
1301     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1302       \textcolor{keywordflow}{if} (k<kb(i)) \textcolor{keywordflow}{then} ; do\_i(i) = .false. ; cycle ;\textcolor{keywordflow}{ endif}
1303       i\_rem = i\_rem + 1  \textcolor{comment}{! Count the i-rows that are still being worked on.}
1304       \textcolor{comment}{!   Apply vertical decay of the turbulent energy.  This energy is}
1305       \textcolor{comment}{! simply lost.}
1306       tke(i) = tke(i) * exp(-i2decay(i) * (gv%H\_to\_Z*(h(i,j,k) + h(i,j,k+1))))
1307 
1308 \textcolor{comment}{!      if (maxEnt(i,k) <= 0.0) cycle}
1309       \textcolor{keywordflow}{if} (maxtke(i,k) <= 0.0) cycle
1310 
1311   \textcolor{comment}{! This is an analytic integral where diffusity is a quadratic function of}
1312   \textcolor{comment}{! rho that goes asymptotically to 0 at Rho\_top (vaguely following KPP?).}
1313       \textcolor{keywordflow}{if} (tke(i) > 0.0) \textcolor{keywordflow}{then}
1314         \textcolor{keywordflow}{if} (rint(k) <= rho\_top(i)) \textcolor{keywordflow}{then}
1315           tke\_to\_layer = tke(i)
1316         \textcolor{keywordflow}{else}
1317           drl = rint(k+1) - rint(k) ; drbot = rint(k+1) - rho\_top(i)
1318           tke\_to\_layer = tke(i) * drl * &
1319               (3.0*drbot*(rint(k) - rho\_top(i)) + drl**2) / (drbot**3)
1320 \textcolor{keywordflow}{        endif}
1321       \textcolor{keywordflow}{else} ; tke\_to\_layer = 0.0 ;\textcolor{keywordflow}{ endif}
1322 
1323       \textcolor{comment}{! TKE\_Ray has been initialized to 0 above.}
1324       \textcolor{keywordflow}{if} (rayleigh\_drag) tke\_ray = 0.5*cs%BBL\_effic * us%L\_to\_Z**2 * g%IareaT(i,j) * &
1325             ((g%areaCu(i-1,j) * visc%Ray\_u(i-1,j,k) * u(i-1,j,k)**2 + &
1326               g%areaCu(i,j)   * visc%Ray\_u(i,j,k)   * u(i,j,k)**2) + &
1327              (g%areaCv(i,j-1) * visc%Ray\_v(i,j-1,k) * v(i,j-1,k)**2 + &
1328               g%areaCv(i,j)   * visc%Ray\_v(i,j,k)   * v(i,j,k)**2))
1329 
1330       \textcolor{keywordflow}{if} (tke\_to\_layer + tke\_ray > 0.0) \textcolor{keywordflow}{then}
1331         \textcolor{keywordflow}{if} (cs%BBL\_mixing\_as\_max) \textcolor{keywordflow}{then}
1332           \textcolor{keywordflow}{if} (tke\_to\_layer + tke\_ray > maxtke(i,k)) &
1333               tke\_to\_layer = maxtke(i,k) - tke\_ray
1334 
1335           tke(i) = tke(i) - tke\_to\_layer
1336 
1337           \textcolor{keywordflow}{if} (kd\_lay(i,k) < (tke\_to\_layer + tke\_ray) * tke\_to\_kd(i,k)) \textcolor{keywordflow}{then}
1338             delta\_kd = (tke\_to\_layer + tke\_ray) * tke\_to\_kd(i,k) - kd\_lay(i,k)
1339             \textcolor{keywordflow}{if} ((cs%Kd\_max >= 0.0) .and. (delta\_kd > cs%Kd\_max)) \textcolor{keywordflow}{then}
1340               delta\_kd = cs%Kd\_max
1341               kd\_lay(i,k) = kd\_lay(i,k) + delta\_kd
1342             \textcolor{keywordflow}{else}
1343               kd\_lay(i,k) = (tke\_to\_layer + tke\_ray) * tke\_to\_kd(i,k)
1344 \textcolor{keywordflow}{            endif}
1345             \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) \textcolor{keywordflow}{then}
1346               kd\_int(i,k)   = kd\_int(i,k)   + 0.5 * delta\_kd
1347               kd\_int(i,k+1) = kd\_int(i,k+1) + 0.5 * delta\_kd
1348 \textcolor{keywordflow}{            endif}
1349             \textcolor{keywordflow}{if} (do\_diag\_kd\_bbl) \textcolor{keywordflow}{then}
1350               kd\_bbl(i,j,k) = kd\_bbl(i,j,k) + 0.5 * delta\_kd
1351               kd\_bbl(i,j,k+1) = kd\_bbl(i,j,k+1) + 0.5 * delta\_kd
1352 \textcolor{keywordflow}{            endif}
1353 \textcolor{keywordflow}{          endif}
1354         \textcolor{keywordflow}{else}
1355           \textcolor{keywordflow}{if} (kd\_lay(i,k) >= maxtke(i,k) * tke\_to\_kd(i,k)) \textcolor{keywordflow}{then}
1356             tke\_here = 0.0
1357             tke(i) = tke(i) + tke\_ray
1358           \textcolor{keywordflow}{elseif} (kd\_lay(i,k) + (tke\_to\_layer + tke\_ray) * tke\_to\_kd(i,k) > &
1359                   maxtke(i,k) * tke\_to\_kd(i,k)) \textcolor{keywordflow}{then}
1360             tke\_here = ((tke\_to\_layer + tke\_ray) + kd\_lay(i,k) / tke\_to\_kd(i,k)) - maxtke(i,k)
1361             tke(i) = (tke(i) - tke\_here) + tke\_ray
1362           \textcolor{keywordflow}{else}
1363             tke\_here = tke\_to\_layer + tke\_ray
1364             tke(i) = tke(i) - tke\_to\_layer
1365 \textcolor{keywordflow}{          endif}
1366           \textcolor{keywordflow}{if} (tke(i) < 0.0) tke(i) = 0.0 \textcolor{comment}{! This should be unnecessary?}
1367 
1368           \textcolor{keywordflow}{if} (tke\_here > 0.0) \textcolor{keywordflow}{then}
1369             delta\_kd = tke\_here * tke\_to\_kd(i,k)
1370             \textcolor{keywordflow}{if} (cs%Kd\_max >= 0.0) delta\_kd = min(delta\_kd, cs%Kd\_max)
1371             kd\_lay(i,k) = kd\_lay(i,k) + delta\_kd
1372             \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) \textcolor{keywordflow}{then}
1373               kd\_int(i,k)   = kd\_int(i,k)   + 0.5 * delta\_kd
1374               kd\_int(i,k+1) = kd\_int(i,k+1) + 0.5 * delta\_kd
1375 \textcolor{keywordflow}{            endif}
1376             \textcolor{keywordflow}{if} (do\_diag\_kd\_bbl) \textcolor{keywordflow}{then}
1377               kd\_bbl(i,j,k) = kd\_bbl(i,j,k) + 0.5 * delta\_kd
1378               kd\_bbl(i,j,k+1) = kd\_bbl(i,j,k+1) + 0.5 * delta\_kd
1379 \textcolor{keywordflow}{            endif}
1380 \textcolor{keywordflow}{          endif}
1381 \textcolor{keywordflow}{        endif}
1382 \textcolor{keywordflow}{      endif}
1383 
1384       \textcolor{comment}{! This may be risky - in the case that there are exactly zero}
1385       \textcolor{comment}{! velocities at 4 neighboring points, but nonzero velocities}
1386       \textcolor{comment}{! above the iterations would stop too soon. I don't see how this}
1387       \textcolor{comment}{! could happen in practice. RWH}
1388       \textcolor{keywordflow}{if} ((tke(i)<= 0.0) .and. (tke\_ray == 0.0)) \textcolor{keywordflow}{then}
1389         do\_i(i) = .false. ; i\_rem = i\_rem - 1
1390 \textcolor{keywordflow}{      endif}
1391 
1392 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1393     \textcolor{keywordflow}{if} (i\_rem == 0) \textcolor{keywordflow}{exit}
1394 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! k-loop}
1395 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a968690a1c9efb859ee67965f9032c7d7}\label{namespacemom__set__diffusivity_a968690a1c9efb859ee67965f9032c7d7}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!add\+\_\+lotw\+\_\+bbl\+\_\+diffusivity@{add\+\_\+lotw\+\_\+bbl\+\_\+diffusivity}}
\index{add\+\_\+lotw\+\_\+bbl\+\_\+diffusivity@{add\+\_\+lotw\+\_\+bbl\+\_\+diffusivity}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{add\+\_\+lotw\+\_\+bbl\+\_\+diffusivity()}{add\_lotw\_bbl\_diffusivity()}}
{\footnotesize\ttfamily subroutine mom\+\_\+set\+\_\+diffusivity\+::add\+\_\+lotw\+\_\+bbl\+\_\+diffusivity (\begin{DoxyParamCaption}\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h,  }\item[{real, dimension( g \%isdb\+: g \%iedb, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{u,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsdb\+: g \%jedb, g \%ke), intent(in)}]{v,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{type(forcing), intent(in)}]{fluxes,  }\item[{type(vertvisc\+\_\+type), intent(in)}]{visc,  }\item[{integer, intent(in)}]{j,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(in)}]{N2\+\_\+int,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{real, dimension(\+:,\+:,\+:), pointer}]{Kd\+\_\+\+B\+BL,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(inout), optional}]{Kd\+\_\+lay,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(inout), optional}]{Kd\+\_\+int }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



Calculates a B\+BL diffusivity use a Prandtl number 1 diffusivity with a law of the wall turbulent viscosity, up to a B\+BL height where the energy used for mixing has consumed the mechanical T\+KE input. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & Grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & Vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em u} & u component of flow \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em v} & v component of flow \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h} & Layer thickness \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em tv} & Structure containing pointers to any available thermodynamic fields.\\
\hline
\mbox{\tt in}  & {\em fluxes} & Surface fluxes structure\\
\hline
\mbox{\tt in}  & {\em visc} & Structure containing vertical viscosities, bottom boundary layer properies, and related fields.\\
\hline
\mbox{\tt in}  & {\em j} & j-\/index of row to work on\\
\hline
\mbox{\tt in}  & {\em n2\+\_\+int} & Square of Brunt-\/\+Vaisala at interfaces \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}\\
\hline
 & {\em cs} & Diffusivity control structure\\
\hline
 & {\em kd\+\_\+bbl} & Interface B\+BL diffusivity \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}\\
\hline
\mbox{\tt in,out}  & {\em kd\+\_\+lay} & Layer net diffusivity \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}\\
\hline
\mbox{\tt in,out}  & {\em kd\+\_\+int} & Interface net diffusivity \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]} \\
\hline
\end{DoxyParams}


Definition at line 1403 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
1403   \textcolor{keywordtype}{type}(ocean\_grid\_type),    \textcolor{keywordtype}{intent(in)}    :: g\textcolor{comment}{  !< Grid structure}
1404   \textcolor{keywordtype}{type}(verticalgrid\_type),  \textcolor{keywordtype}{intent(in)}    :: gv\textcolor{comment}{ !< Vertical grid structure}
1405   \textcolor{keywordtype}{type}(unit\_scale\_type),    \textcolor{keywordtype}{intent(in)}    :: us\textcolor{comment}{ !< A dimensional unit scaling type}
1406   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, &
1407                             \textcolor{keywordtype}{intent(in)}    :: u\textcolor{comment}{  !< u component of flow [L T-1 ~> m s-1]}
1408   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, &
1409                             \textcolor{keywordtype}{intent(in)}    :: v\textcolor{comment}{  !< v component of flow [L T-1 ~> m s-1]}
1410   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1411                             \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{  !< Layer thickness [H ~> m or kg m-2]}
1412   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),    \textcolor{keywordtype}{intent(in)}    :: tv\textcolor{comment}{ !< Structure containing pointers to any available}
1413 \textcolor{comment}{                                                !! thermodynamic fields.}
1414   \textcolor{keywordtype}{type}(forcing),            \textcolor{keywordtype}{intent(in)}    :: fluxes\textcolor{comment}{ !< Surface fluxes structure}
1415   \textcolor{keywordtype}{type}(vertvisc\_type),      \textcolor{keywordtype}{intent(in)}    :: visc\textcolor{comment}{ !< Structure containing vertical viscosities, bottom}
1416 \textcolor{comment}{                                                  !! boundary layer properies, and related fields.}
1417   \textcolor{keywordtype}{integer},                  \textcolor{keywordtype}{intent(in)}    :: j\textcolor{comment}{  !< j-index of row to work on}
1418   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)}, &
1419                             \textcolor{keywordtype}{intent(in)}    :: n2\_int\textcolor{comment}{ !< Square of Brunt-Vaisala at interfaces [T-2 ~> s-2]}
1420   \textcolor{keywordtype}{type}(set\_diffusivity\_cs), \textcolor{keywordtype}{pointer}       :: cs\textcolor{comment}{ !< Diffusivity control structure}
1421   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:,:)},   \textcolor{keywordtype}{pointer}       :: kd\_bbl\textcolor{comment}{ !< Interface BBL diffusivity [Z2 T-1 ~> m2 s-1]}
1422   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, &
1423                   \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: kd\_lay\textcolor{comment}{ !< Layer net diffusivity [Z2 T-1 ~> m2 s-1]}
1424   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)}, &
1425                   \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: kd\_int\textcolor{comment}{ !< Interface net diffusivity [Z2 T-1 ~> m2 s-1]}
1426 
1427   \textcolor{comment}{! Local variables}
1428   \textcolor{keywordtype}{real} :: tke\_column       \textcolor{comment}{! net TKE input into the column [Z3 T-3 ~> m3 s-3]}
1429   \textcolor{keywordtype}{real} :: tke\_remaining    \textcolor{comment}{! remaining TKE available for mixing in this layer and above [Z3 T-3 ~> m3 s-3]}
1430   \textcolor{keywordtype}{real} :: tke\_consumed     \textcolor{comment}{! TKE used for mixing in this layer [Z3 T-3 ~> m3 s-3]}
1431   \textcolor{keywordtype}{real} :: tke\_kd\_wall      \textcolor{comment}{! TKE associated with unlimited law of the wall mixing [Z3 T-3 ~> m3 s-3]}
1432   \textcolor{keywordtype}{real} :: cdrag\_sqrt       \textcolor{comment}{! square root of the drag coefficient [nondim]}
1433   \textcolor{keywordtype}{real} :: ustar            \textcolor{comment}{! value of ustar at a thickness point [Z T-1 ~> m s-1].}
1434   \textcolor{keywordtype}{real} :: ustar2           \textcolor{comment}{! square of ustar, for convenience [Z2 T-2 ~> m2 s-2]}
1435   \textcolor{keywordtype}{real} :: absf             \textcolor{comment}{! average absolute value of Coriolis parameter around a thickness point [T-1 ~>
       s-1]}
1436   \textcolor{keywordtype}{real} :: dh, dhm1         \textcolor{comment}{! thickness of layers k and k-1, respecitvely [Z ~> m].}
1437   \textcolor{keywordtype}{real} :: z\_bot            \textcolor{comment}{! distance to interface k from bottom [Z ~> m].}
1438   \textcolor{keywordtype}{real} :: d\_minus\_z        \textcolor{comment}{! distance to interface k from surface [Z ~> m].}
1439   \textcolor{keywordtype}{real} :: total\_thickness  \textcolor{comment}{! total thickness of water column [Z ~> m].}
1440   \textcolor{keywordtype}{real} :: idecay           \textcolor{comment}{! inverse of decay scale used for "Joule heating" loss of TKE with height [Z-1
       ~> m-1].}
1441   \textcolor{keywordtype}{real} :: kd\_wall          \textcolor{comment}{! Law of the wall diffusivity [Z2 T-1 ~> m2 s-1].}
1442   \textcolor{keywordtype}{real} :: kd\_lower         \textcolor{comment}{! diffusivity for lower interface [Z2 T-1 ~> m2 s-1]}
1443   \textcolor{keywordtype}{real} :: ustar\_d          \textcolor{comment}{! u* x D  [Z2 T-1 ~> m2 s-1].}
1444   \textcolor{keywordtype}{real} :: i\_rho0           \textcolor{comment}{! 1 / rho0 [R-1  ~> m3 kg-1]}
1445   \textcolor{keywordtype}{real} :: n2\_min           \textcolor{comment}{! Minimum value of N2 to use in calculation of TKE\_Kd\_wall [T-2 ~> s-2]}
1446   \textcolor{keywordtype}{logical} :: rayleigh\_drag \textcolor{comment}{! Set to true if there are Rayleigh drag velocities defined in visc, on}
1447                            \textcolor{comment}{! the assumption that this extracted energy also drives diapycnal mixing.}
1448   \textcolor{keywordtype}{integer} :: i, k, km1
1449   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: von\_karm = 0.41 \textcolor{comment}{! Von Karman constant
       (http://en.wikipedia.org/wiki/Von\_Karman\_constant)}
1450   \textcolor{keywordtype}{logical} :: do\_diag\_kd\_bbl
1451 
1452   \textcolor{keywordflow}{if} (.not.(cs%bottomdraglaw .and. (cs%BBL\_effic>0.0))) \textcolor{keywordflow}{return}
1453   do\_diag\_kd\_bbl = \textcolor{keyword}{associated}(kd\_bbl)
1454 
1455   n2\_min = 0.
1456   \textcolor{keywordflow}{if} (cs%LOTW\_BBL\_use\_omega) n2\_min = cs%omega**2
1457 
1458   \textcolor{comment}{! Determine whether to add Rayleigh drag contribution to TKE}
1459   rayleigh\_drag = .false.
1460   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Ray\_u) .and. \textcolor{keyword}{associated}(visc%Ray\_v)) rayleigh\_drag = .true.
1461   i\_rho0 = 1.0 / (gv%Rho0)
1462   cdrag\_sqrt = sqrt(cs%cdrag)
1463 
1464   \textcolor{keywordflow}{do} i=g%isc,g%iec \textcolor{comment}{! Developed in single-column mode}
1465 
1466     \textcolor{comment}{! Column-wise parameters.}
1467     absf = 0.25 * ((abs(g%CoriolisBu(i-1,j-1)) + abs(g%CoriolisBu(i,j))) + &
1468                    (abs(g%CoriolisBu(i-1,j)) + abs(g%CoriolisBu(i,j-1)))) \textcolor{comment}{! Non-zero on equator!}
1469 
1470     \textcolor{comment}{! u* at the bottom [Z T-1 ~> m s-1].}
1471     ustar = visc%ustar\_BBL(i,j)
1472     ustar2 = ustar**2
1473     \textcolor{comment}{! In add\_drag\_diffusivity(), fluxes%ustar\_tidal is added in. This might be double counting}
1474     \textcolor{comment}{! since ustar\_BBL should already include all contributions to u*? -AJA}
1475     \textcolor{comment}{!### Examine the question of whether there is double counting of fluxes%ustar\_tidal.}
1476     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(fluxes%ustar\_tidal)) ustar = ustar + fluxes%ustar\_tidal(i,j)
1477 
1478     \textcolor{comment}{! The maximum decay scale should be something of order 200 m. We use the smaller of u*/f and}
1479     \textcolor{comment}{! (IMax\_decay)^-1 as the decay scale. If ustar = 0, this is land so this value doesn't matter.}
1480     idecay = cs%IMax\_decay
1481     \textcolor{keywordflow}{if} ((ustar > 0.0) .and. (absf > cs%IMax\_decay * ustar)) idecay = absf / ustar
1482 
1483     \textcolor{comment}{! Energy input at the bottom [Z3 T-3 ~> m3 s-3].}
1484     \textcolor{comment}{! (Note that visc%TKE\_BBL is in [Z3 T-3 ~> m3 s-3], set in set\_BBL\_TKE().)}
1485     \textcolor{comment}{! I am still unsure about sqrt(cdrag) in this expressions - AJA}
1486     tke\_column = cdrag\_sqrt * visc%TKE\_BBL(i,j)
1487     \textcolor{comment}{! Add in tidal dissipation energy at the bottom [R Z3 T-3 ~> m3 s-3].}
1488     \textcolor{comment}{! Note that TKE\_tidal is in [R Z3 T-3 ~> W m-2].}
1489     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(fluxes%TKE\_tidal)) &
1490       tke\_column = tke\_column + fluxes%TKE\_tidal(i,j) * i\_rho0
1491     tke\_column = cs%BBL\_effic * tke\_column \textcolor{comment}{! Only use a fraction of the mechanical dissipation for mixing.}
1492 
1493     tke\_remaining = tke\_column
1494     total\_thickness = ( sum(h(i,j,:)) + gv%H\_subroundoff )* gv%H\_to\_Z \textcolor{comment}{! Total column thickness [Z ~> m].}
1495     ustar\_d = ustar * total\_thickness
1496     z\_bot = 0.
1497     kd\_lower = 0. \textcolor{comment}{! Diffusivity on bottom boundary.}
1498 
1499     \textcolor{comment}{! Work upwards from the bottom, accumulating work used until it exceeds the available TKE input}
1500     \textcolor{comment}{! at the bottom.}
1501     \textcolor{keywordflow}{do} k=g%ke,2,-1
1502       dh = gv%H\_to\_Z * h(i,j,k) \textcolor{comment}{! Thickness of this level [Z ~> m].}
1503       km1 = max(k-1, 1)
1504       dhm1 = gv%H\_to\_Z * h(i,j,km1) \textcolor{comment}{! Thickness of level above [Z ~> m].}
1505 
1506       \textcolor{comment}{! Add in additional energy input from bottom-drag against slopes (sides)}
1507       \textcolor{keywordflow}{if} (rayleigh\_drag) tke\_remaining = tke\_remaining + &
1508             0.5*cs%BBL\_effic * us%L\_to\_Z**2 * g%IareaT(i,j) * &
1509             ((g%areaCu(i-1,j) * visc%Ray\_u(i-1,j,k) * u(i-1,j,k)**2 + &
1510               g%areaCu(i,j)   * visc%Ray\_u(i,j,k)   * u(i,j,k)**2) + &
1511              (g%areaCv(i,j-1) * visc%Ray\_v(i,j-1,k) * v(i,j-1,k)**2 + &
1512               g%areaCv(i,j)   * visc%Ray\_v(i,j,k)   * v(i,j,k)**2))
1513 
1514       \textcolor{comment}{! Exponentially decay TKE across the thickness of the layer.}
1515       \textcolor{comment}{! This is energy loss in addition to work done as mixing, apparently to Joule heating.}
1516       tke\_remaining = exp(-idecay*dh) * tke\_remaining
1517 
1518       z\_bot = z\_bot + h(i,j,k)*gv%H\_to\_Z \textcolor{comment}{! Distance between upper interface of layer and the bottom [Z ~>
       m].}
1519       d\_minus\_z = max(total\_thickness - z\_bot, 0.) \textcolor{comment}{! Thickness above layer [Z ~> m].}
1520 
1521       \textcolor{comment}{! Diffusivity using law of the wall, limited by rotation, at height z [Z2 T-1 ~> m2 s-1].}
1522       \textcolor{comment}{! This calculation is at the upper interface of the layer}
1523       \textcolor{keywordflow}{if} ( ustar\_d + absf * ( z\_bot * d\_minus\_z ) == 0.) \textcolor{keywordflow}{then}
1524         kd\_wall = 0.
1525       \textcolor{keywordflow}{else}
1526         kd\_wall = ((von\_karm * ustar2) * (z\_bot * d\_minus\_z)) &
1527                   / (ustar\_d + absf * (z\_bot * d\_minus\_z))
1528 \textcolor{keywordflow}{      endif}
1529 
1530       \textcolor{comment}{! TKE associated with Kd\_wall [Z3 T-3 ~> m3 s-3].}
1531       \textcolor{comment}{! This calculation if for the volume spanning the interface.}
1532       tke\_kd\_wall = kd\_wall * 0.5 * (dh + dhm1) * max(n2\_int(i,k), n2\_min)
1533 
1534       \textcolor{comment}{! Now bound Kd such that the associated TKE is no greater than available TKE for mixing.}
1535       \textcolor{keywordflow}{if} (tke\_kd\_wall > 0.) \textcolor{keywordflow}{then}
1536         tke\_consumed = min(tke\_kd\_wall, tke\_remaining)
1537         kd\_wall = (tke\_consumed / tke\_kd\_wall) * kd\_wall \textcolor{comment}{! Scale Kd so that only TKE\_consumed is used.}
1538       \textcolor{keywordflow}{else}
1539         \textcolor{comment}{! Either N2=0 or dh = 0.}
1540         \textcolor{keywordflow}{if} (tke\_remaining > 0.) \textcolor{keywordflow}{then}
1541           kd\_wall = cs%Kd\_max
1542         \textcolor{keywordflow}{else}
1543           kd\_wall = 0.
1544 \textcolor{keywordflow}{        endif}
1545         tke\_consumed = 0.
1546 \textcolor{keywordflow}{      endif}
1547 
1548       \textcolor{comment}{! Now use up the appropriate about of TKE associated with the diffusivity chosen}
1549       tke\_remaining = tke\_remaining - tke\_consumed \textcolor{comment}{! Note this will be non-negative}
1550 
1551       \textcolor{comment}{! Add this BBL diffusivity to the model net diffusivity.}
1552       \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) kd\_int(i,k) = kd\_int(i,k) + kd\_wall
1553       \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_lay)) kd\_lay(i,k) = kd\_lay(i,k) + 0.5 * (kd\_wall + kd\_lower)
1554       kd\_lower = kd\_wall \textcolor{comment}{! Store for next layer up.}
1555       \textcolor{keywordflow}{if} (do\_diag\_kd\_bbl) kd\_bbl(i,j,k) = kd\_wall
1556 \textcolor{keywordflow}{    enddo} \textcolor{comment}{! k}
1557 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! i}
1558 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a21162cb5173d52ce92dc65beee9d0ef1}\label{namespacemom__set__diffusivity_a21162cb5173d52ce92dc65beee9d0ef1}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!add\+\_\+mlrad\+\_\+diffusivity@{add\+\_\+mlrad\+\_\+diffusivity}}
\index{add\+\_\+mlrad\+\_\+diffusivity@{add\+\_\+mlrad\+\_\+diffusivity}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{add\+\_\+mlrad\+\_\+diffusivity()}{add\_mlrad\_diffusivity()}}
{\footnotesize\ttfamily subroutine mom\+\_\+set\+\_\+diffusivity\+::add\+\_\+mlrad\+\_\+diffusivity (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{h,  }\item[{type(forcing), intent(in)}]{fluxes,  }\item[{integer, intent(in)}]{j,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{real, dimension(szi\+\_\+(g),szk\+\_\+(g)), intent(in)}]{T\+K\+E\+\_\+to\+\_\+\+Kd,  }\item[{real, dimension(szi\+\_\+(g),szk\+\_\+(g)), intent(inout), optional}]{Kd\+\_\+lay,  }\item[{real, dimension(szi\+\_\+(g),szk\+\_\+(g)+1), intent(inout), optional}]{Kd\+\_\+int }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



This routine adds effects of mixed layer radiation to the layer diffusivities. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em fluxes} & Surface fluxes structure\\
\hline
\mbox{\tt in}  & {\em j} & The j-\/index to work on\\
\hline
 & {\em cs} & Diffusivity control structure\\
\hline
\mbox{\tt in}  & {\em tke\+\_\+to\+\_\+kd} & The conversion rate between the T\+KE T\+KE dissipated within a layer and the diapycnal diffusivity witin that layer, usually ($\sim$\+Rho\+\_\+0 / (G\+\_\+\+Earth $\ast$ d\+Rho\+\_\+lay)) \mbox{[}Z2 T-\/1 / Z3 T-\/3 = T2 Z-\/1 $\sim$$>$ s2 m-\/1\mbox{]}\\
\hline
\mbox{\tt in,out}  & {\em kd\+\_\+lay} & The diapycnal diffusivity in layers \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\
\hline
\mbox{\tt in,out}  & {\em kd\+\_\+int} & The diapycnal diffusivity at interfaces \\
\hline
\end{DoxyParams}


Definition at line 1563 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
1563   \textcolor{keywordtype}{type}(ocean\_grid\_type),            \textcolor{keywordtype}{intent(in)}    :: g\textcolor{comment}{      !< The ocean's grid structure}
1564   \textcolor{keywordtype}{type}(verticalgrid\_type),          \textcolor{keywordtype}{intent(in)}    :: gv\textcolor{comment}{     !< The ocean's vertical grid structure}
1565   \textcolor{keywordtype}{type}(unit\_scale\_type),            \textcolor{keywordtype}{intent(in)}    :: us\textcolor{comment}{     !< A dimensional unit scaling type}
1566   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1567                                     \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{      !< Layer thicknesses [H ~> m or kg m-2]}
1568   \textcolor{keywordtype}{type}(forcing),                    \textcolor{keywordtype}{intent(in)}    :: fluxes\textcolor{comment}{ !< Surface fluxes structure}
1569   \textcolor{keywordtype}{integer},                          \textcolor{keywordtype}{intent(in)}    :: j\textcolor{comment}{      !< The j-index to work on}
1570   \textcolor{keywordtype}{type}(set\_diffusivity\_cs),         \textcolor{keywordtype}{pointer}       :: cs\textcolor{comment}{     !< Diffusivity control structure}
1571   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)}    :: tke\_to\_kd\textcolor{comment}{ !< The conversion rate between the TKE}
1572 \textcolor{comment}{                                                            !! TKE dissipated within  a layer and the}
1573 \textcolor{comment}{                                                            !! diapycnal diffusivity witin that layer,}
1574 \textcolor{comment}{                                                            !! usually (~Rho\_0 / (G\_Earth * dRho\_lay))}
1575 \textcolor{comment}{                                                            !! [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]}
1576   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, &
1577                           \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: kd\_lay\textcolor{comment}{ !< The diapycnal diffusivity in layers [Z2 T-1
       ~> m2 s-1].}
1578   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)}, &
1579                           \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: kd\_int\textcolor{comment}{ !< The diapycnal diffusivity at interfaces}
1580 \textcolor{comment}{                                                            !! [Z2 T-1 ~> m2 s-1].}
1581 
1582 \textcolor{comment}{! This routine adds effects of mixed layer radiation to the layer diffusivities.}
1583 
1584   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: h\_ml  \textcolor{comment}{! Mixed layer thickness [Z ~> m].}
1585   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: tke\_ml\_flux \textcolor{comment}{! Mixed layer TKE flux [Z3 T-3 ~> m3 s-3]}
1586   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: i\_decay \textcolor{comment}{! A decay rate [Z-1 ~> m-1].}
1587   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: kd\_mlr\_ml \textcolor{comment}{! Diffusivities associated with mixed layer radiation [Z2 T-1 ~> m2
       s-1].}
1588 
1589   \textcolor{keywordtype}{real} :: f\_sq              \textcolor{comment}{! The square of the local Coriolis parameter or a related variable [T-2 ~>
       s-2].}
1590   \textcolor{keywordtype}{real} :: h\_ml\_sq           \textcolor{comment}{! The square of the mixed layer thickness [Z2 ~> m2].}
1591   \textcolor{keywordtype}{real} :: ustar\_sq          \textcolor{comment}{! ustar squared [Z2 T-2 ~> m2 s-2]}
1592   \textcolor{keywordtype}{real} :: kd\_mlr            \textcolor{comment}{! A diffusivity associated with mixed layer turbulence radiation [Z2 T-1 ~> m2
       s-1].}
1593   \textcolor{keywordtype}{real} :: c1\_6              \textcolor{comment}{! 1/6}
1594   \textcolor{keywordtype}{real} :: omega2            \textcolor{comment}{! rotation rate squared [T-2 ~> s-2].}
1595   \textcolor{keywordtype}{real} :: z1                \textcolor{comment}{! layer thickness times I\_decay [nondim]}
1596   \textcolor{keywordtype}{real} :: dzl               \textcolor{comment}{! thickness converted to heights [Z ~> m].}
1597   \textcolor{keywordtype}{real} :: i\_decay\_len2\_tke  \textcolor{comment}{! squared inverse decay lengthscale for}
1598                             \textcolor{comment}{! TKE, as used in the mixed layer code [Z-2 ~> m-2].}
1599   \textcolor{keywordtype}{real} :: h\_neglect         \textcolor{comment}{! negligibly small thickness [Z ~> m].}
1600 
1601   \textcolor{keywordtype}{logical} :: do\_any, do\_i(szi\_(g))
1602   \textcolor{keywordtype}{integer} :: i, k, is, ie, nz, kml
1603   is = g%isc ; ie = g%iec ; nz = g%ke
1604 
1605   omega2    = cs%omega**2
1606   c1\_6      = 1.0 / 6.0
1607   kml       = gv%nkml
1608   h\_neglect = gv%H\_subroundoff*gv%H\_to\_Z
1609 
1610   \textcolor{keywordflow}{if} (.not.cs%ML\_radiation) \textcolor{keywordflow}{return}
1611 
1612   \textcolor{keywordflow}{do} i=is,ie ; h\_ml(i) = 0.0 ; do\_i(i) = (g%mask2dT(i,j) > 0.5) ;\textcolor{keywordflow}{ enddo}
1613   \textcolor{keywordflow}{do} k=1,kml ; \textcolor{keywordflow}{do} i=is,ie ; h\_ml(i) = h\_ml(i) + gv%H\_to\_Z*h(i,j,k) ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo}
1614 
1615   \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1616     \textcolor{keywordflow}{if} (cs%ML\_omega\_frac >= 1.0) \textcolor{keywordflow}{then}
1617       f\_sq = 4.0 * omega2
1618     \textcolor{keywordflow}{else}
1619       f\_sq = 0.25 * ((g%CoriolisBu(i,j)**2 + g%CoriolisBu(i-1,j-1)**2) + &
1620                      (g%CoriolisBu(i,j-1)**2 + g%CoriolisBu(i-1,j)**2))
1621       \textcolor{keywordflow}{if} (cs%ML\_omega\_frac > 0.0) &
1622         f\_sq = cs%ML\_omega\_frac * 4.0 * omega2 + (1.0 - cs%ML\_omega\_frac) * f\_sq
1623 \textcolor{keywordflow}{    endif}
1624 
1625     ustar\_sq = max(fluxes%ustar(i,j), cs%ustar\_min)**2
1626 
1627     tke\_ml\_flux(i) = (cs%mstar * cs%ML\_rad\_coeff) * (ustar\_sq * (fluxes%ustar(i,j)))
1628     i\_decay\_len2\_tke = cs%TKE\_decay**2 * (f\_sq / ustar\_sq)
1629 
1630     \textcolor{keywordflow}{if} (cs%ML\_rad\_TKE\_decay) &
1631       tke\_ml\_flux(i) = tke\_ml\_flux(i) * exp(-h\_ml(i) * sqrt(i\_decay\_len2\_tke))
1632 
1633     \textcolor{comment}{! Calculate the inverse decay scale}
1634     h\_ml\_sq = (cs%ML\_rad\_efold\_coeff * (h\_ml(i)+h\_neglect))**2
1635     i\_decay(i) = sqrt((i\_decay\_len2\_tke * h\_ml\_sq + 1.0) / h\_ml\_sq)
1636 
1637     \textcolor{comment}{! Average the dissipation layer kml+1, using}
1638     \textcolor{comment}{! a more accurate Taylor series approximations for very thin layers.}
1639     z1 = (gv%H\_to\_Z*h(i,j,kml+1)) * i\_decay(i)
1640     \textcolor{keywordflow}{if} (z1 > 1e-5) \textcolor{keywordflow}{then}
1641       kd\_mlr = tke\_ml\_flux(i) * tke\_to\_kd(i,kml+1) * (1.0 - exp(-z1))
1642     \textcolor{keywordflow}{else}
1643       kd\_mlr = tke\_ml\_flux(i) * tke\_to\_kd(i,kml+1) * (z1 * (1.0 - z1 * (0.5 - c1\_6 * z1)))
1644 \textcolor{keywordflow}{    endif}
1645     kd\_mlr\_ml(i) = min(kd\_mlr, cs%ML\_rad\_kd\_max)
1646     tke\_ml\_flux(i) = tke\_ml\_flux(i) * exp(-z1)
1647 \textcolor{keywordflow}{  endif} ;\textcolor{keywordflow}{ enddo}
1648 
1649   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_lay)) \textcolor{keywordflow}{then}
1650     \textcolor{keywordflow}{do} k=1,kml+1 ; \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1651       kd\_lay(i,k) = kd\_lay(i,k) + kd\_mlr\_ml(i)
1652 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo}
1653 \textcolor{keywordflow}{  endif}
1654   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) \textcolor{keywordflow}{then}
1655     \textcolor{keywordflow}{do} k=2,kml+1 ; \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1656       kd\_int(i,k) = kd\_int(i,k) + kd\_mlr\_ml(i)
1657 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo}
1658     \textcolor{keywordflow}{if} (kml<=nz-1) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1659       kd\_int(i,kml+2) = kd\_int(i,kml+2) + 0.5 * kd\_mlr\_ml(i)
1660 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
1661 \textcolor{keywordflow}{  endif}
1662 
1663   \textcolor{keywordflow}{do} k=kml+2,nz-1
1664     do\_any = .false.
1665     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1666       dzl = gv%H\_to\_Z*h(i,j,k) ;  z1 = dzl*i\_decay(i)
1667       \textcolor{keywordflow}{if} (cs%ML\_Rad\_bug) \textcolor{keywordflow}{then}
1668         \textcolor{comment}{! These expresssions are dimensionally inconsistent. -RWH}
1669         \textcolor{comment}{! This is supposed to be the integrated energy deposited in the layer,}
1670         \textcolor{comment}{! not the average over the layer as in these expressions.}
1671         \textcolor{keywordflow}{if} (z1 > 1e-5) \textcolor{keywordflow}{then}
1672           kd\_mlr = (tke\_ml\_flux(i) * tke\_to\_kd(i,k)) * & \textcolor{comment}{! Units of Z2 T-1}
1673                    us%m\_to\_Z * ((1.0 - exp(-z1)) / dzl)  \textcolor{comment}{! Units of m-1}
1674         \textcolor{keywordflow}{else}
1675           kd\_mlr = (tke\_ml\_flux(i) * tke\_to\_kd(i,k)) * &  \textcolor{comment}{! Units of Z2 T-1}
1676                    us%m\_to\_Z * (i\_decay(i) * (1.0 - z1 * (0.5 - c1\_6*z1))) \textcolor{comment}{! Units of m-1}
1677 \textcolor{keywordflow}{        endif}
1678       \textcolor{keywordflow}{else}
1679         \textcolor{keywordflow}{if} (z1 > 1e-5) \textcolor{keywordflow}{then}
1680           kd\_mlr = (tke\_ml\_flux(i) * tke\_to\_kd(i,k)) * (1.0 - exp(-z1))
1681         \textcolor{keywordflow}{else}
1682           kd\_mlr = (tke\_ml\_flux(i) * tke\_to\_kd(i,k)) * (z1 * (1.0 - z1 * (0.5 - c1\_6*z1)))
1683 \textcolor{keywordflow}{        endif}
1684 \textcolor{keywordflow}{      endif}
1685       kd\_mlr = min(kd\_mlr, cs%ML\_rad\_kd\_max)
1686       \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_lay)) \textcolor{keywordflow}{then}
1687         kd\_lay(i,k) = kd\_lay(i,k) + kd\_mlr
1688 \textcolor{keywordflow}{      endif}
1689       \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) \textcolor{keywordflow}{then}
1690         kd\_int(i,k)   = kd\_int(i,k)   + 0.5 * kd\_mlr
1691         kd\_int(i,k+1) = kd\_int(i,k+1) + 0.5 * kd\_mlr
1692 \textcolor{keywordflow}{      endif}
1693 
1694       tke\_ml\_flux(i) = tke\_ml\_flux(i) * exp(-z1)
1695       \textcolor{keywordflow}{if} (tke\_ml\_flux(i) * i\_decay(i) < 0.1 * cs%Kd\_min * omega2) \textcolor{keywordflow}{then}
1696         do\_i(i) = .false.
1697       \textcolor{keywordflow}{else} ; do\_any = .true. ;\textcolor{keywordflow}{ endif}
1698 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1699     \textcolor{keywordflow}{if} (.not.do\_any) \textcolor{keywordflow}{exit}
1700 \textcolor{keywordflow}{  enddo}
1701 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a99d0eb7701f8e04d856b75117fe7b83c}\label{namespacemom__set__diffusivity_a99d0eb7701f8e04d856b75117fe7b83c}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!double\+\_\+diffusion@{double\+\_\+diffusion}}
\index{double\+\_\+diffusion@{double\+\_\+diffusion}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{double\+\_\+diffusion()}{double\_diffusion()}}
{\footnotesize\ttfamily subroutine mom\+\_\+set\+\_\+diffusivity\+::double\+\_\+diffusion (\begin{DoxyParamCaption}\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{T\+\_\+f,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{S\+\_\+f,  }\item[{integer, intent(in)}]{j,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(out)}]{Kd\+\_\+\+T\+\_\+dd,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(out)}]{Kd\+\_\+\+S\+\_\+dd }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



This subroutine sets the additional diffusivities of temperature and salinity due to double diffusion, using the same functional form as is used in M\+O\+M4.\+1, and taken from an N\+C\+AR technical note (R\+EF?) that updates what was in Large et al. (1994). All the coefficients here should probably be made run-\/time variables rather than hard-\/coded constants. 

\begin{DoxyRefDesc}{Todo}
\item[\hyperlink{todo__todo000006}{Todo}]Find reference for N\+C\+AR tech note above. \end{DoxyRefDesc}



\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure.\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em tv} & Structure containing pointers to any available thermodynamic fields; absent fields have N\+U\+LL ptrs.\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em t\+\_\+f} & layer temperatures with the values in massless layers\\
\hline
\mbox{\tt in}  & {\em s\+\_\+f} & Layer salinities with values in massless\\
\hline
\mbox{\tt in}  & {\em j} & Meridional index upon which to work.\\
\hline
 & {\em cs} & Module control structure.\\
\hline
\mbox{\tt out}  & {\em kd\+\_\+t\+\_\+dd} & Interface double diffusion diapycnal\\
\hline
\mbox{\tt out}  & {\em kd\+\_\+s\+\_\+dd} & Interface double diffusion diapycnal \\
\hline
\end{DoxyParams}


Definition at line 1077 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
1077   \textcolor{keywordtype}{type}(ocean\_grid\_type),    \textcolor{keywordtype}{intent(in)}  :: g\textcolor{comment}{   !< The ocean's grid structure.}
1078   \textcolor{keywordtype}{type}(verticalgrid\_type),  \textcolor{keywordtype}{intent(in)}  :: gv\textcolor{comment}{  !< The ocean's vertical grid structure.}
1079   \textcolor{keywordtype}{type}(unit\_scale\_type),    \textcolor{keywordtype}{intent(in)}  :: us\textcolor{comment}{  !< A dimensional unit scaling type}
1080   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),    \textcolor{keywordtype}{intent(in)}  :: tv\textcolor{comment}{  !< Structure containing pointers to any available}
1081 \textcolor{comment}{                                               !! thermodynamic fields; absent fields have NULL ptrs.}
1082   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1083                             \textcolor{keywordtype}{intent(in)}  :: h\textcolor{comment}{   !< Layer thicknesses [H ~> m or kg m-2].}
1084   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1085                             \textcolor{keywordtype}{intent(in)}  :: t\_f\textcolor{comment}{ !< layer temperatures with the values in massless layers}
1086 \textcolor{comment}{                                               !! filled vertically by diffusion [degC].}
1087   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1088                             \textcolor{keywordtype}{intent(in)}  :: s\_f\textcolor{comment}{ !< Layer salinities with values in massless}
1089 \textcolor{comment}{                                               !! layers filled vertically by diffusion [ppt].}
1090   \textcolor{keywordtype}{integer},                  \textcolor{keywordtype}{intent(in)}  :: j\textcolor{comment}{   !< Meridional index upon which to work.}
1091   \textcolor{keywordtype}{type}(set\_diffusivity\_cs), \textcolor{keywordtype}{pointer}     :: cs\textcolor{comment}{  !< Module control structure.}
1092   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)},       &
1093                             \textcolor{keywordtype}{intent(out)} :: kd\_t\_dd\textcolor{comment}{ !< Interface double diffusion diapycnal}
1094 \textcolor{comment}{                                               !! diffusivity for temp [Z2 T-1 ~> m2 s-1].}
1095   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)},       &
1096                             \textcolor{keywordtype}{intent(out)} :: kd\_s\_dd\textcolor{comment}{ !< Interface double diffusion diapycnal}
1097 \textcolor{comment}{                                               !! diffusivity for saln [Z2 T-1 ~> m2 s-1].}
1098 
1099   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: &
1100     drho\_dt,  &    \textcolor{comment}{! partial derivatives of density wrt temp [R degC-1 ~> kg m-3 degC-1]}
1101     drho\_ds,  &    \textcolor{comment}{! partial derivatives of density wrt saln [R ppt-1 ~> kg m-3 ppt-1]}
1102     pres,     &    \textcolor{comment}{! pressure at each interface [R L2 T-2 ~> Pa]}
1103     temp\_int, &    \textcolor{comment}{! temperature at interfaces [degC]}
1104     salin\_int      \textcolor{comment}{! Salinity at interfaces [ppt]}
1105 
1106   \textcolor{keywordtype}{real} ::  alpha\_dt \textcolor{comment}{! density difference between layers due to temp diffs [R ~> kg m-3]}
1107   \textcolor{keywordtype}{real} ::  beta\_ds  \textcolor{comment}{! density difference between layers due to saln diffs [R ~> kg m-3]}
1108 
1109   \textcolor{keywordtype}{real} :: rrho    \textcolor{comment}{! vertical density ratio [nondim]}
1110   \textcolor{keywordtype}{real} :: diff\_dd \textcolor{comment}{! factor for double-diffusion [nondim]}
1111   \textcolor{keywordtype}{real} :: kd\_dd   \textcolor{comment}{! The dominant double diffusive diffusivity [Z2 T-1 ~> m2 s-1]}
1112   \textcolor{keywordtype}{real} :: prandtl \textcolor{comment}{! flux ratio for diffusive convection regime}
1113 
1114   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: rrho0  = 1.9 \textcolor{comment}{! limit for double-diffusive density ratio [nondim]}
1115 
1116   \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The i-computational domain for the equation of state}
1117   \textcolor{keywordtype}{integer} :: i, k, is, ie, nz
1118   is = g%isc ; ie = g%iec ; nz = g%ke
1119 
1120   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%eqn\_of\_state)) \textcolor{keywordflow}{then}
1121     \textcolor{keywordflow}{do} i=is,ie
1122       pres(i) = 0.0 ; kd\_t\_dd(i,1) = 0.0 ; kd\_s\_dd(i,1) = 0.0
1123       kd\_t\_dd(i,nz+1) = 0.0 ; kd\_s\_dd(i,nz+1) = 0.0
1124 \textcolor{keywordflow}{    enddo}
1125     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%p\_surf)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} i=is,ie ; pres(i) = tv%p\_surf(i,j) ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
1126     eosdom(:) = eos\_domain(g%HI)
1127     \textcolor{keywordflow}{do} k=2,nz
1128       \textcolor{keywordflow}{do} i=is,ie
1129         pres(i) = pres(i) + (gv%g\_Earth*gv%H\_to\_RZ)*h(i,j,k-1)
1130         temp\_int(i) = 0.5 * (t\_f(i,j,k-1) + t\_f(i,j,k))
1131         salin\_int(i) = 0.5 * (s\_f(i,j,k-1) + s\_f(i,j,k))
1132 \textcolor{keywordflow}{      enddo}
1133       \textcolor{keyword}{call }calculate\_density\_derivs(temp\_int, salin\_int, pres, drho\_dt, drho\_ds, &
1134                                     tv%eqn\_of\_state, eosdom)
1135 
1136       \textcolor{keywordflow}{do} i=is,ie
1137         alpha\_dt = -1.0*drho\_dt(i) * (t\_f(i,j,k-1) - t\_f(i,j,k))
1138         beta\_ds  = drho\_ds(i) * (s\_f(i,j,k-1) - s\_f(i,j,k))
1139 
1140         \textcolor{keywordflow}{if} ((alpha\_dt > beta\_ds) .and. (beta\_ds > 0.0)) \textcolor{keywordflow}{then}  \textcolor{comment}{! salt finger case}
1141           rrho = min(alpha\_dt / beta\_ds, rrho0)
1142           diff\_dd = 1.0 - ((rrho-1.0)/(rrho0-1.0))
1143           kd\_dd = cs%Max\_salt\_diff\_salt\_fingers * diff\_dd*diff\_dd*diff\_dd
1144           kd\_t\_dd(i,k) = 0.7 * kd\_dd
1145           kd\_s\_dd(i,k) = kd\_dd
1146         \textcolor{keywordflow}{elseif} ((alpha\_dt < 0.) .and. (beta\_ds < 0.) .and. (alpha\_dt > beta\_ds)) \textcolor{keywordflow}{then} \textcolor{comment}{! diffusive
       convection}
1147           rrho = alpha\_dt / beta\_ds
1148           kd\_dd = cs%Kv\_molecular * 0.909 * exp(4.6 * exp(-0.54 * (1/rrho - 1)))
1149           prandtl = 0.15*rrho
1150           \textcolor{keywordflow}{if} (rrho > 0.5) prandtl = (1.85-0.85/rrho)*rrho
1151           kd\_t\_dd(i,k) = kd\_dd
1152           kd\_s\_dd(i,k) = prandtl * kd\_dd
1153         \textcolor{keywordflow}{else}
1154           kd\_t\_dd(i,k) = 0.0 ; kd\_s\_dd(i,k) = 0.0
1155 \textcolor{keywordflow}{        endif}
1156 \textcolor{keywordflow}{      enddo}
1157 \textcolor{keywordflow}{    enddo}
1158 \textcolor{keywordflow}{  endif}
1159 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_afef80c2221be24f63f1ca96c2abe6fa9}\label{namespacemom__set__diffusivity_afef80c2221be24f63f1ca96c2abe6fa9}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!find\+\_\+n2@{find\+\_\+n2}}
\index{find\+\_\+n2@{find\+\_\+n2}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{find\+\_\+n2()}{find\_n2()}}
{\footnotesize\ttfamily subroutine mom\+\_\+set\+\_\+diffusivity\+::find\+\_\+n2 (\begin{DoxyParamCaption}\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{T\+\_\+f,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{S\+\_\+f,  }\item[{type(forcing), intent(in)}]{fluxes,  }\item[{integer, intent(in)}]{j,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(out)}]{d\+Rho\+\_\+int,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(out)}]{N2\+\_\+lay,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(out)}]{N2\+\_\+int,  }\item[{real, dimension( g \%isd\+: g \%ied), intent(out)}]{N2\+\_\+bot }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



Calculate Brunt-\/\+Vaisala frequency, N$^\wedge$2. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em tv} & Structure containing pointers to any available thermodynamic fields.\\
\hline
\mbox{\tt in}  & {\em t\+\_\+f} & layer temperature with the values in massless layers\\
\hline
\mbox{\tt in}  & {\em s\+\_\+f} & Layer salinities with values in massless\\
\hline
\mbox{\tt in}  & {\em fluxes} & A structure of thermodynamic surface fluxes\\
\hline
\mbox{\tt in}  & {\em j} & j-\/index of row to work on\\
\hline
 & {\em cs} & Diffusivity control structure\\
\hline
\mbox{\tt out}  & {\em drho\+\_\+int} & Change in locally referenced potential density\\
\hline
\mbox{\tt out}  & {\em n2\+\_\+int} & The squared buoyancy frequency at the interfaces \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\
\hline
\mbox{\tt out}  & {\em n2\+\_\+lay} & The squared buoyancy frequency of the layers \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\
\hline
\mbox{\tt out}  & {\em n2\+\_\+bot} & The near-\/bottom squared buoyancy frequency \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}. \\
\hline
\end{DoxyParams}


Definition at line 906 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
906   \textcolor{keywordtype}{type}(ocean\_grid\_type),    \textcolor{keywordtype}{intent(in)}  :: g\textcolor{comment}{    !< The ocean's grid structure}
907   \textcolor{keywordtype}{type}(verticalgrid\_type),  \textcolor{keywordtype}{intent(in)}  :: gv\textcolor{comment}{   !< The ocean's vertical grid structure}
908   \textcolor{keywordtype}{type}(unit\_scale\_type),    \textcolor{keywordtype}{intent(in)}  :: us\textcolor{comment}{   !< A dimensional unit scaling type}
909   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
910                             \textcolor{keywordtype}{intent(in)}  :: h\textcolor{comment}{    !< Layer thicknesses [H ~> m or kg m-2]}
911   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),    \textcolor{keywordtype}{intent(in)}  :: tv\textcolor{comment}{   !< Structure containing pointers to any available}
912 \textcolor{comment}{                                                !! thermodynamic fields.}
913   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
914                             \textcolor{keywordtype}{intent(in)}  :: t\_f\textcolor{comment}{  !< layer temperature with the values in massless layers}
915 \textcolor{comment}{                                                !! filled vertically by diffusion [degC].}
916   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
917                             \textcolor{keywordtype}{intent(in)}  :: s\_f\textcolor{comment}{  !< Layer salinities with values in massless}
918 \textcolor{comment}{                                                !! layers filled vertically by diffusion [ppt].}
919   \textcolor{keywordtype}{type}(forcing),            \textcolor{keywordtype}{intent(in)}  :: fluxes\textcolor{comment}{ !< A structure of thermodynamic surface fluxes}
920   \textcolor{keywordtype}{integer},                  \textcolor{keywordtype}{intent(in)}  :: j\textcolor{comment}{    !< j-index of row to work on}
921   \textcolor{keywordtype}{type}(set\_diffusivity\_cs), \textcolor{keywordtype}{pointer}     :: cs\textcolor{comment}{   !< Diffusivity control structure}
922   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)}, &
923                             \textcolor{keywordtype}{intent(out)} :: drho\_int\textcolor{comment}{ !< Change in locally referenced potential density}
924 \textcolor{comment}{                                                !! across each interface [R ~> kg m-3].}
925   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)}, &
926                             \textcolor{keywordtype}{intent(out)} :: n2\_int\textcolor{comment}{ !< The squared buoyancy frequency at the interfaces [T-2
       ~> s-2].}
927   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, &
928                             \textcolor{keywordtype}{intent(out)} :: n2\_lay\textcolor{comment}{ !< The squared buoyancy frequency of the layers [T-2 ~>
       s-2].}
929   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))}, \textcolor{keywordtype}{intent(out)} :: n2\_bot\textcolor{comment}{ !< The near-bottom squared buoyancy frequency [T-2 ~>
       s-2].}
930   \textcolor{comment}{! Local variables}
931   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)} :: &
932     drho\_int\_unfilt, & \textcolor{comment}{! unfiltered density differences across interfaces [R ~> kg m-3]}
933     drho\_dt,         & \textcolor{comment}{! partial derivative of density wrt temp [R degC-1 ~> kg m-3 degC-1]}
934     drho\_ds            \textcolor{comment}{! partial derivative of density wrt saln [R ppt-1 ~> kg m-3 ppt-1]}
935 
936   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: &
937     pres,      &  \textcolor{comment}{! pressure at each interface [R L2 T-2 ~> Pa]}
938     temp\_int,  &  \textcolor{comment}{! temperature at each interface [degC]}
939     salin\_int, &  \textcolor{comment}{! salinity at each interface [ppt]}
940     drho\_bot,  &  \textcolor{comment}{! A density difference [R ~> kg m-3]}
941     h\_amp,     &  \textcolor{comment}{! The topographic roughness amplitude [Z ~> m].}
942     hb,        &  \textcolor{comment}{! The thickness of the bottom layer [Z ~> m].}
943     z\_from\_bot    \textcolor{comment}{! The hieght above the bottom [Z ~> m].}
944 
945   \textcolor{keywordtype}{real} :: rml\_base  \textcolor{comment}{! density of the deepest variable density layer}
946   \textcolor{keywordtype}{real} :: dz\_int    \textcolor{comment}{! thickness associated with an interface [Z ~> m].}
947   \textcolor{keywordtype}{real} :: g\_rho0    \textcolor{comment}{! gravitation acceleration divided by Bouss reference density}
948                     \textcolor{comment}{! times some unit conversion factors [Z T-2 R-1 ~> m4 s-2 kg-1].}
949   \textcolor{keywordtype}{real} :: h\_neglect \textcolor{comment}{! negligibly small thickness, in the same units as h.}
950 
951   \textcolor{keywordtype}{logical} :: do\_i(szi\_(g)), do\_any
952   \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The i-computational domain for the equation of state}
953   \textcolor{keywordtype}{integer} :: i, k, is, ie, nz
954 
955   is = g%isc ; ie = g%iec ; nz = g%ke
956   g\_rho0    = (us%L\_to\_Z**2 * gv%g\_Earth) / (gv%Rho0)
957   h\_neglect = gv%H\_subroundoff
958 
959   \textcolor{comment}{! Find the (limited) density jump across each interface.}
960   \textcolor{keywordflow}{do} i=is,ie
961     drho\_int(i,1) = 0.0 ; drho\_int(i,nz+1) = 0.0
962     drho\_int\_unfilt(i,1) = 0.0 ; drho\_int\_unfilt(i,nz+1) = 0.0
963 \textcolor{keywordflow}{  enddo}
964   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%eqn\_of\_state)) \textcolor{keywordflow}{then}
965     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(fluxes%p\_surf)) \textcolor{keywordflow}{then}
966       \textcolor{keywordflow}{do} i=is,ie ; pres(i) = fluxes%p\_surf(i,j) ;\textcolor{keywordflow}{ enddo}
967     \textcolor{keywordflow}{else}
968       \textcolor{keywordflow}{do} i=is,ie ; pres(i) = 0.0 ;\textcolor{keywordflow}{ enddo}
969 \textcolor{keywordflow}{    endif}
970     eosdom(:) = eos\_domain(g%HI)
971     \textcolor{keywordflow}{do} k=2,nz
972       \textcolor{keywordflow}{do} i=is,ie
973         pres(i) = pres(i) + (gv%g\_Earth*gv%H\_to\_RZ)*h(i,j,k-1)
974         temp\_int(i) = 0.5 * (t\_f(i,j,k) + t\_f(i,j,k-1))
975         salin\_int(i) = 0.5 * (s\_f(i,j,k) + s\_f(i,j,k-1))
976 \textcolor{keywordflow}{      enddo}
977       \textcolor{keyword}{call }calculate\_density\_derivs(temp\_int, salin\_int, pres, drho\_dt(:,k), drho\_ds(:,k), &
978                                     tv%eqn\_of\_state, eosdom)
979       \textcolor{keywordflow}{do} i=is,ie
980         drho\_int(i,k) = max(drho\_dt(i,k)*(t\_f(i,j,k) - t\_f(i,j,k-1)) + &
981                             drho\_ds(i,k)*(s\_f(i,j,k) - s\_f(i,j,k-1)), 0.0)
982         drho\_int\_unfilt(i,k) = max(drho\_dt(i,k)*(tv%T(i,j,k) - tv%T(i,j,k-1)) + &
983                             drho\_ds(i,k)*(tv%S(i,j,k) - tv%S(i,j,k-1)), 0.0)
984 \textcolor{keywordflow}{      enddo}
985 \textcolor{keywordflow}{    enddo}
986   \textcolor{keywordflow}{else}
987     \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=is,ie
988       drho\_int(i,k) = gv%Rlay(k) - gv%Rlay(k-1)
989 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo}
990 \textcolor{keywordflow}{  endif}
991 
992   \textcolor{comment}{! Set the buoyancy frequencies.}
993   \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
994     n2\_lay(i,k) = g\_rho0 * 0.5*(drho\_int(i,k) + drho\_int(i,k+1)) / &
995                   (gv%H\_to\_Z*(h(i,j,k) + h\_neglect))
996 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
997   \textcolor{keywordflow}{do} i=is,ie ; n2\_int(i,1) = 0.0 ; n2\_int(i,nz+1) = 0.0 ;\textcolor{keywordflow}{ enddo}
998   \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=is,ie
999     n2\_int(i,k) = g\_rho0 * drho\_int(i,k) / &
1000                   (0.5*gv%H\_to\_Z*(h(i,j,k-1) + h(i,j,k) + h\_neglect))
1001 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
1002 
1003   \textcolor{comment}{! Find the bottom boundary layer stratification, and use this in the deepest layers.}
1004   \textcolor{keywordflow}{do} i=is,ie
1005     hb(i) = 0.0 ; drho\_bot(i) = 0.0 ; h\_amp(i) = 0.0
1006     z\_from\_bot(i) = 0.5*gv%H\_to\_Z*h(i,j,nz)
1007     do\_i(i) = (g%mask2dT(i,j) > 0.5)
1008 \textcolor{keywordflow}{  enddo}
1009   \textcolor{keywordflow}{if} (cs%use\_tidal\_mixing) \textcolor{keyword}{call }tidal\_mixing\_h\_amp(h\_amp, g, j, cs%tidal\_mixing\_CSp)
1010 
1011   \textcolor{keywordflow}{do} k=nz,2,-1
1012     do\_any = .false.
1013     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1014       dz\_int = 0.5*gv%H\_to\_Z*(h(i,j,k) + h(i,j,k-1))
1015       z\_from\_bot(i) = z\_from\_bot(i) + dz\_int \textcolor{comment}{! middle of the layer above}
1016 
1017       hb(i) = hb(i) + dz\_int
1018       drho\_bot(i) = drho\_bot(i) + drho\_int(i,k)
1019 
1020       \textcolor{keywordflow}{if} (z\_from\_bot(i) > h\_amp(i)) \textcolor{keywordflow}{then}
1021         \textcolor{keywordflow}{if} (k>2) \textcolor{keywordflow}{then}
1022           \textcolor{comment}{! Always include at least one full layer.}
1023           hb(i) = hb(i) + 0.5*gv%H\_to\_Z*(h(i,j,k-1) + h(i,j,k-2))
1024           drho\_bot(i) = drho\_bot(i) + drho\_int(i,k-1)
1025 \textcolor{keywordflow}{        endif}
1026         do\_i(i) = .false.
1027       \textcolor{keywordflow}{else}
1028         do\_any = .true.
1029 \textcolor{keywordflow}{      endif}
1030 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1031     \textcolor{keywordflow}{if} (.not.do\_any) \textcolor{keywordflow}{exit}
1032 \textcolor{keywordflow}{  enddo}
1033 
1034   \textcolor{keywordflow}{do} i=is,ie
1035     \textcolor{keywordflow}{if} (hb(i) > 0.0) \textcolor{keywordflow}{then}
1036       n2\_bot(i) = (g\_rho0 * drho\_bot(i)) / hb(i)
1037     \textcolor{keywordflow}{else} ;  n2\_bot(i) = 0.0 ;\textcolor{keywordflow}{ endif}
1038     z\_from\_bot(i) = 0.5*gv%H\_to\_Z*h(i,j,nz)
1039     do\_i(i) = (g%mask2dT(i,j) > 0.5)
1040 \textcolor{keywordflow}{  enddo}
1041 
1042   \textcolor{keywordflow}{do} k=nz,2,-1
1043     do\_any = .false.
1044     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1045       dz\_int = 0.5*gv%H\_to\_Z*(h(i,j,k) + h(i,j,k-1))
1046       z\_from\_bot(i) = z\_from\_bot(i) + dz\_int \textcolor{comment}{! middle of the layer above}
1047 
1048       n2\_int(i,k) = n2\_bot(i)
1049       \textcolor{keywordflow}{if} (k>2) n2\_lay(i,k-1) = n2\_bot(i)
1050 
1051       \textcolor{keywordflow}{if} (z\_from\_bot(i) > h\_amp(i)) \textcolor{keywordflow}{then}
1052         \textcolor{keywordflow}{if} (k>2) n2\_int(i,k-1) = n2\_bot(i)
1053         do\_i(i) = .false.
1054       \textcolor{keywordflow}{else}
1055         do\_any = .true.
1056 \textcolor{keywordflow}{      endif}
1057 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1058     \textcolor{keywordflow}{if} (.not.do\_any) \textcolor{keywordflow}{exit}
1059 \textcolor{keywordflow}{  enddo}
1060 
1061   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%eqn\_of\_state)) \textcolor{keywordflow}{then}
1062     \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
1063       drho\_int(i,k) = drho\_int\_unfilt(i,k)
1064 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo}
1065 \textcolor{keywordflow}{  endif}
1066 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a07c0ab3f141f8f9e057be3150a940a94}\label{namespacemom__set__diffusivity_a07c0ab3f141f8f9e057be3150a940a94}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!find\+\_\+tke\+\_\+to\+\_\+kd@{find\+\_\+tke\+\_\+to\+\_\+kd}}
\index{find\+\_\+tke\+\_\+to\+\_\+kd@{find\+\_\+tke\+\_\+to\+\_\+kd}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{find\+\_\+tke\+\_\+to\+\_\+kd()}{find\_tke\_to\_kd()}}
{\footnotesize\ttfamily subroutine mom\+\_\+set\+\_\+diffusivity\+::find\+\_\+tke\+\_\+to\+\_\+kd (\begin{DoxyParamCaption}\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke+1), intent(in)}]{d\+Rho\+\_\+int,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(in)}]{N2\+\_\+lay,  }\item[{integer, intent(in)}]{j,  }\item[{real, intent(in)}]{dt,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(out)}]{T\+K\+E\+\_\+to\+\_\+\+Kd,  }\item[{real, dimension( g \%isd\+: g \%ied, g \%ke), intent(out)}]{max\+T\+KE,  }\item[{integer, dimension( g \%isd\+: g \%ied), intent(out)}]{kb }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}



Convert turbulent kinetic energy to diffusivity. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em tv} & Structure containing pointers to any available thermodynamic fields.\\
\hline
\mbox{\tt in}  & {\em drho\+\_\+int} & Change in locally referenced potential density across each interface \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em n2\+\_\+lay} & The squared buoyancy frequency of the layers \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em j} & j-\/index of row to work on\\
\hline
\mbox{\tt in}  & {\em dt} & Time increment \mbox{[}T $\sim$$>$ s\mbox{]}.\\
\hline
 & {\em cs} & Diffusivity control structure\\
\hline
\mbox{\tt out}  & {\em tke\+\_\+to\+\_\+kd} & The conversion rate between the T\+KE dissipated within a layer and the diapycnal diffusivity witin that layer, usually ($\sim$\+Rho\+\_\+0 / (G\+\_\+\+Earth $\ast$ d\+Rho\+\_\+lay)) \mbox{[}Z2 T-\/1 / Z3 T-\/3 = T2 Z-\/1 $\sim$$>$ s2 m-\/1\mbox{]}\\
\hline
\mbox{\tt out}  & {\em maxtke} & The energy required to for a layer to entrain to its maximum realizable thickness \mbox{[}Z3 T-\/3 $\sim$$>$ m3 s-\/3\mbox{]}\\
\hline
\mbox{\tt out}  & {\em kb} & Index of lightest layer denser than the buffer layer, or -\/1 without a bulk mixed layer. \\
\hline
\end{DoxyParams}


Definition at line 692 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
692   \textcolor{keywordtype}{type}(ocean\_grid\_type),            \textcolor{keywordtype}{intent(in)}    :: g\textcolor{comment}{    !< The ocean's grid structure}
693   \textcolor{keywordtype}{type}(verticalgrid\_type),          \textcolor{keywordtype}{intent(in)}    :: gv\textcolor{comment}{   !< The ocean's vertical grid structure}
694   \textcolor{keywordtype}{type}(unit\_scale\_type),            \textcolor{keywordtype}{intent(in)}    :: us\textcolor{comment}{   !< A dimensional unit scaling type}
695   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
696                                     \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{    !< Layer thicknesses [H ~> m or kg m-2]}
697   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),            \textcolor{keywordtype}{intent(in)}    :: tv\textcolor{comment}{   !< Structure containing pointers to any available}
698 \textcolor{comment}{                                                          !! thermodynamic fields.}
699   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)}, \textcolor{keywordtype}{intent(in)}  :: drho\_int\textcolor{comment}{ !< Change in locally referenced potential
       density}
700 \textcolor{comment}{                                                          !! across each interface [R ~> kg m-3].}
701   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)}    :: n2\_lay\textcolor{comment}{ !< The squared buoyancy frequency of the}
702 \textcolor{comment}{                                                          !! layers [T-2 ~> s-2].}
703   \textcolor{keywordtype}{integer},                          \textcolor{keywordtype}{intent(in)}    :: j\textcolor{comment}{    !< j-index of row to work on}
704   \textcolor{keywordtype}{real},                             \textcolor{keywordtype}{intent(in)}    :: dt\textcolor{comment}{   !< Time increment [T ~> s].}
705   \textcolor{keywordtype}{type}(set\_diffusivity\_cs),         \textcolor{keywordtype}{pointer}       :: cs\textcolor{comment}{   !< Diffusivity control structure}
706   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)}   :: tke\_to\_kd\textcolor{comment}{ !< The conversion rate between the}
707 \textcolor{comment}{                                                          !! TKE dissipated within a layer and the}
708 \textcolor{comment}{                                                          !! diapycnal diffusivity witin that layer,}
709 \textcolor{comment}{                                                          !! usually (~Rho\_0 / (G\_Earth * dRho\_lay))}
710 \textcolor{comment}{                                                          !! [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]}
711   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)}   :: maxtke\textcolor{comment}{ !< The energy required to for a layer to
       entrain}
712 \textcolor{comment}{                                                          !! to its maximum realizable thickness [Z3 T-3 ~>
       m3 s-3]}
713   \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(SZI\_(G))},      \textcolor{keywordtype}{intent(out)}   :: kb\textcolor{comment}{   !< Index of lightest layer denser than the buffer}
714 \textcolor{comment}{                                                          !! layer, or -1 without a bulk mixed layer.}
715   \textcolor{comment}{! Local variables}
716   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))} :: &
717     ds\_dsp1, &    \textcolor{comment}{! coordinate variable (sigma-2) difference across an}
718                   \textcolor{comment}{! interface divided by the difference across the interface}
719                   \textcolor{comment}{! below it [nondim]}
720     dsp1\_ds, &    \textcolor{comment}{! inverse coordinate variable (sigma-2) difference}
721                   \textcolor{comment}{! across an interface times the difference across the}
722                   \textcolor{comment}{! interface above it [nondim]}
723     rho\_0,   &    \textcolor{comment}{! Layer potential densities relative to surface pressure [R ~> kg m-3]}
724     maxent        \textcolor{comment}{! maxEnt is the maximum value of entrainment from below (with}
725                   \textcolor{comment}{! compensating entrainment from above to keep the layer}
726                   \textcolor{comment}{! density from changing) that will not deplete all of the}
727                   \textcolor{comment}{! layers above or below a layer within a timestep [Z ~> m].}
728   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: &
729     htot,    &    \textcolor{comment}{! total thickness above or below a layer, or the}
730                   \textcolor{comment}{! integrated thickness in the BBL [Z ~> m].}
731     mfkb,    &    \textcolor{comment}{! total thickness in the mixed and buffer layers times ds\_dsp1 [Z ~> m].}
732     p\_ref,   &    \textcolor{comment}{! array of tv%P\_Ref pressures [R L2 T-2 ~> Pa]}
733     rcv\_kmb, &    \textcolor{comment}{! coordinate density in the lowest buffer layer [R ~> kg m-3]}
734     p\_0           \textcolor{comment}{! An array of 0 pressures [R L2 T-2 ~> Pa]}
735 
736   \textcolor{keywordtype}{real} :: dh\_max      \textcolor{comment}{! maximum amount of entrainment a layer could}
737                       \textcolor{comment}{! undergo before entraining all fluid in the layers}
738                       \textcolor{comment}{! above or below [Z ~> m].}
739   \textcolor{keywordtype}{real} :: drho\_lay    \textcolor{comment}{! density change across a layer [R ~> kg m-3]}
740   \textcolor{keywordtype}{real} :: omega2      \textcolor{comment}{! rotation rate squared [T-2 ~> s-2]}
741   \textcolor{keywordtype}{real} :: g\_rho0      \textcolor{comment}{! gravitation accel divided by Bouss ref density [Z T-2 R-1 ~> m4 s-2 kg-1]}
742   \textcolor{keywordtype}{real} :: g\_irho0     \textcolor{comment}{! Alternate calculation of G\_Rho0 for reproducibility [Z T-2 R-1 ~> m4 s-2 kg-1]}
743   \textcolor{keywordtype}{real} :: i\_rho0      \textcolor{comment}{! inverse of Boussinesq reference density [R-1 ~> m3 kg-1]}
744   \textcolor{keywordtype}{real} :: i\_dt        \textcolor{comment}{! 1/dt [T-1 ~> s-1]}
745   \textcolor{keywordtype}{real} :: h\_neglect   \textcolor{comment}{! negligibly small thickness [H ~> m or kg m-2]}
746   \textcolor{keywordtype}{real} :: hn2po2      \textcolor{comment}{! h (N^2 + Omega^2), in [Z T-2 ~> m s-2].}
747   \textcolor{keywordtype}{logical} :: do\_i(szi\_(g))
748 
749   \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The i-computational domain for the equation of state}
750   \textcolor{keywordtype}{integer} :: i, k, is, ie, nz, i\_rem, kmb, kb\_min
751   is = g%isc ; ie = g%iec ; nz = g%ke
752 
753   i\_dt      = 1.0 / dt
754   omega2    = cs%omega**2
755   h\_neglect = gv%H\_subroundoff
756   g\_rho0    = (us%L\_to\_Z**2 * gv%g\_Earth) / (gv%Rho0)
757   \textcolor{keywordflow}{if} (cs%answers\_2018) \textcolor{keywordflow}{then}
758     i\_rho0    = 1.0 / (gv%Rho0)
759     g\_irho0 = (us%L\_to\_Z**2 * gv%g\_Earth) * i\_rho0
760   \textcolor{keywordflow}{else}
761     g\_irho0 = g\_rho0
762 \textcolor{keywordflow}{  endif}
763 
764   \textcolor{comment}{! Simple but coordinate-independent estimate of Kd/TKE}
765   \textcolor{keywordflow}{if} (cs%simple\_TKE\_to\_Kd) \textcolor{keywordflow}{then}
766     \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
767       hn2po2 = (gv%H\_to\_Z * h(i,j,k)) * (n2\_lay(i,k) + omega2) \textcolor{comment}{! Units of Z T-2.}
768       \textcolor{keywordflow}{if} (hn2po2>0.) \textcolor{keywordflow}{then}
769         tke\_to\_kd(i,k) = 1.0 / hn2po2 \textcolor{comment}{! Units of T2 Z-1.}
770       else; tke\_to\_kd(i,k) = 0.;\textcolor{keywordflow}{ endif}
771       \textcolor{comment}{! The maximum TKE conversion we allow is really a statement}
772       \textcolor{comment}{! about the upper diffusivity we allow. Kd\_max must be set.}
773       maxtke(i,k) = hn2po2 * cs%Kd\_max \textcolor{comment}{! Units of Z3 T-3.}
774 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo}
775     kb(is:ie) = -1 \textcolor{comment}{! kb should not be used by any code in non-layered mode -AJA}
776     \textcolor{keywordflow}{return}
777 \textcolor{keywordflow}{  endif}
778 
779   \textcolor{comment}{! Determine kb - the index of the shallowest active interior layer.}
780   \textcolor{keywordflow}{if} (cs%bulkmixedlayer) \textcolor{keywordflow}{then}
781     kmb = gv%nk\_rho\_varies
782     \textcolor{keywordflow}{do} i=is,ie ; p\_0(i) = 0.0 ; p\_ref(i) = tv%P\_Ref ;\textcolor{keywordflow}{ enddo}
783     eosdom(:) = eos\_domain(g%HI)
784     \textcolor{keywordflow}{do} k=1,nz
785       \textcolor{keyword}{call }calculate\_density(tv%T(:,j,k), tv%S(:,j,k), p\_0, rho\_0(:,k), tv%eqn\_of\_state, eosdom)
786 \textcolor{keywordflow}{    enddo}
787     \textcolor{keyword}{call }calculate\_density(tv%T(:,j,kmb), tv%S(:,j,kmb), p\_ref, rcv\_kmb, tv%eqn\_of\_state, eosdom)
788 
789     kb\_min = kmb+1
790     \textcolor{keywordflow}{do} i=is,ie
791       \textcolor{comment}{!   Determine the next denser layer than the buffer layer in the}
792       \textcolor{comment}{! coordinate density (sigma-2).}
793       \textcolor{keywordflow}{do} k=kmb+1,nz-1 ; \textcolor{keywordflow}{if} (rcv\_kmb(i) <= gv%Rlay(k)) \textcolor{keywordflow}{exit} ;\textcolor{keywordflow}{ enddo}
794       kb(i) = k
795 
796     \textcolor{comment}{!   Backtrack, in case there are massive layers above that are stable}
797     \textcolor{comment}{! in sigma-0.}
798       \textcolor{keywordflow}{do} k=kb(i)-1,kmb+1,-1
799         \textcolor{keywordflow}{if} (rho\_0(i,kmb) > rho\_0(i,k)) \textcolor{keywordflow}{exit}
800         \textcolor{keywordflow}{if} (h(i,j,k)>2.0*gv%Angstrom\_H) kb(i) = k
801 \textcolor{keywordflow}{      enddo}
802 \textcolor{keywordflow}{    enddo}
803 
804     \textcolor{keyword}{call }set\_density\_ratios(h, tv, kb, g, gv, us, cs, j, ds\_dsp1, rho\_0)
805   \textcolor{keywordflow}{else} \textcolor{comment}{! not bulkmixedlayer}
806     kb\_min = 2 ; kmb = 0
807     \textcolor{keywordflow}{do} i=is,ie ; kb(i) = 1 ;\textcolor{keywordflow}{ enddo}
808     \textcolor{keyword}{call }set\_density\_ratios(h, tv, kb, g, gv, us, cs, j, ds\_dsp1)
809 \textcolor{keywordflow}{  endif}
810 
811   \textcolor{comment}{! Determine maxEnt - the maximum permitted entrainment from below by each}
812   \textcolor{comment}{! interior layer.}
813   \textcolor{keywordflow}{do} k=2,nz-1 ; \textcolor{keywordflow}{do} i=is,ie
814     dsp1\_ds(i,k) = 1.0 / ds\_dsp1(i,k)
815 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
816   \textcolor{keywordflow}{do} i=is,ie ; dsp1\_ds(i,nz) = 0.0 ;\textcolor{keywordflow}{ enddo}
817 
818   \textcolor{keywordflow}{if} (cs%bulkmixedlayer) \textcolor{keywordflow}{then}
819     kmb = gv%nk\_rho\_varies
820     \textcolor{keywordflow}{do} i=is,ie
821       htot(i) = gv%H\_to\_Z*h(i,j,kmb)
822       mfkb(i) = 0.0
823       \textcolor{keywordflow}{if} (kb(i) < nz) &
824         mfkb(i) = ds\_dsp1(i,kb(i)) * (gv%H\_to\_Z*(h(i,j,kmb) - gv%Angstrom\_H))
825 \textcolor{keywordflow}{    enddo}
826     \textcolor{keywordflow}{do} k=1,kmb-1 ; \textcolor{keywordflow}{do} i=is,ie
827       htot(i) = htot(i) + gv%H\_to\_Z*h(i,j,k)
828       mfkb(i) = mfkb(i) + ds\_dsp1(i,k+1)*(gv%H\_to\_Z*(h(i,j,k) - gv%Angstrom\_H))
829 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo}
830   \textcolor{keywordflow}{else}
831     \textcolor{keywordflow}{do} i=is,i
832       maxent(i,1) = 0.0 ; htot(i) = gv%H\_to\_Z*(h(i,j,1) - gv%Angstrom\_H)
833 \textcolor{keywordflow}{    enddo}
834 \textcolor{keywordflow}{  endif}
835   \textcolor{keywordflow}{do} k=kb\_min,nz-1 ; \textcolor{keywordflow}{do} i=is,ie
836     \textcolor{keywordflow}{if} (k == kb(i)) \textcolor{keywordflow}{then}
837       maxent(i,kb(i)) = mfkb(i)
838     \textcolor{keywordflow}{elseif} (k > kb(i)) \textcolor{keywordflow}{then}
839       \textcolor{keywordflow}{if} (cs%answers\_2018) \textcolor{keywordflow}{then}
840         maxent(i,k) = (1.0/dsp1\_ds(i,k))*(maxent(i,k-1) + htot(i))
841       \textcolor{keywordflow}{else}
842         maxent(i,k) = ds\_dsp1(i,k)*(maxent(i,k-1) + htot(i))
843 \textcolor{keywordflow}{      endif}
844       htot(i) = htot(i) + gv%H\_to\_Z*(h(i,j,k) - gv%Angstrom\_H)
845 \textcolor{keywordflow}{    endif}
846 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
847 
848   \textcolor{keywordflow}{do} i=is,ie
849     htot(i) = gv%H\_to\_Z*(h(i,j,nz) - gv%Angstrom\_H) ; maxent(i,nz) = 0.0
850     do\_i(i) = (g%mask2dT(i,j) > 0.5)
851 \textcolor{keywordflow}{  enddo}
852   \textcolor{keywordflow}{do} k=nz-1,kb\_min,-1
853     i\_rem = 0
854     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
855       \textcolor{keywordflow}{if} (k<kb(i)) \textcolor{keywordflow}{then} ; do\_i(i) = .false. ; cycle ;\textcolor{keywordflow}{ endif}
856       i\_rem = i\_rem + 1  \textcolor{comment}{! Count the i-rows that are still being worked on.}
857       maxent(i,k) = min(maxent(i,k),dsp1\_ds(i,k+1)*maxent(i,k+1) + htot(i))
858       htot(i) = htot(i) + gv%H\_to\_Z*(h(i,j,k) - gv%Angstrom\_H)
859 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
860     \textcolor{keywordflow}{if} (i\_rem == 0) \textcolor{keywordflow}{exit}
861 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! k-loop}
862 
863   \textcolor{comment}{! Now set maxTKE and TKE\_to\_Kd.}
864   \textcolor{keywordflow}{do} i=is,ie
865     maxtke(i,1) = 0.0 ; tke\_to\_kd(i,1) = 0.0
866     maxtke(i,nz) = 0.0 ; tke\_to\_kd(i,nz) = 0.0
867 \textcolor{keywordflow}{  enddo}
868   \textcolor{keywordflow}{do} k=2,kmb ; \textcolor{keywordflow}{do} i=is,ie
869     maxtke(i,k) = 0.0
870     tke\_to\_kd(i,k) = 1.0 / ((n2\_lay(i,k) + omega2) * &
871                             (gv%H\_to\_Z*(h(i,j,k) + h\_neglect)))
872 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
873   \textcolor{keywordflow}{do} k=kmb+1,kb\_min-1 ; \textcolor{keywordflow}{do} i=is,ie
874     \textcolor{comment}{!   These are the properties in the deeper mixed and buffer layers, and}
875     \textcolor{comment}{! should perhaps be revisited.}
876     maxtke(i,k) = 0.0 ; tke\_to\_kd(i,k) = 0.0
877 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
878   \textcolor{keywordflow}{do} k=kb\_min,nz-1 ; \textcolor{keywordflow}{do} i=is,ie
879     \textcolor{keywordflow}{if} (k<kb(i)) \textcolor{keywordflow}{then}
880       maxtke(i,k) = 0.0
881       tke\_to\_kd(i,k) = 0.0
882     \textcolor{keywordflow}{else}
883       \textcolor{comment}{! maxTKE is found by determining the kappa that gives maxEnt.}
884       \textcolor{comment}{!  kappa\_max = I\_dt * dRho\_int(i,K+1) * maxEnt(i,k) * &}
885       \textcolor{comment}{!             (GV%H\_to\_Z*h(i,j,k) + dh\_max) / dRho\_lay}
886       \textcolor{comment}{!  maxTKE(i,k) = (GV%g\_Earth*US%L\_to\_Z**2) * dRho\_lay * kappa\_max}
887       \textcolor{comment}{! dRho\_int should already be non-negative, so the max is redundant?}
888       dh\_max = maxent(i,k) * (1.0 + dsp1\_ds(i,k))
889       drho\_lay = 0.5 * max(drho\_int(i,k) + drho\_int(i,k+1), 0.0)
890       maxtke(i,k) = i\_dt * (g\_irho0 * &
891           (0.5*max(drho\_int(i,k+1) + dsp1\_ds(i,k)*drho\_int(i,k), 0.0))) * &
892            ((gv%H\_to\_Z*h(i,j,k) + dh\_max) * maxent(i,k))
893       \textcolor{comment}{! TKE\_to\_Kd should be rho\_InSitu / G\_Earth * (delta rho\_InSitu)}
894       \textcolor{comment}{! The omega^2 term in TKE\_to\_Kd is due to a rescaling of the efficiency of turbulent}
895       \textcolor{comment}{! mixing by a factor of N^2 / (N^2 + Omega^2), as proposed by Melet et al., 2013?}
896       tke\_to\_kd(i,k) = 1.0 / (g\_rho0 * drho\_lay + &
897                               cs%omega**2 * gv%H\_to\_Z*(h(i,j,k) + h\_neglect))
898 \textcolor{keywordflow}{    endif}
899 \textcolor{keywordflow}{  enddo} ;\textcolor{keywordflow}{ enddo}
900 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a66f77b7e2f9c0c8254da4a0acd5a9996}\label{namespacemom__set__diffusivity_a66f77b7e2f9c0c8254da4a0acd5a9996}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!set\+\_\+bbl\+\_\+tke@{set\+\_\+bbl\+\_\+tke}}
\index{set\+\_\+bbl\+\_\+tke@{set\+\_\+bbl\+\_\+tke}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{set\+\_\+bbl\+\_\+tke()}{set\_bbl\_tke()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+diffusivity\+::set\+\_\+bbl\+\_\+tke (\begin{DoxyParamCaption}\item[{real, dimension(szib\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{u,  }\item[{real, dimension(szi\+\_\+(g),szjb\+\_\+(g),szk\+\_\+(g)), intent(in)}]{v,  }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{h,  }\item[{type(forcing), intent(in)}]{fluxes,  }\item[{type(vertvisc\+\_\+type), intent(in)}]{visc,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{type(ocean\+\_\+obc\+\_\+type), optional, pointer}]{O\+BC }\end{DoxyParamCaption})}



This subroutine calculates several properties related to bottom boundary layer turbulence. 


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em u} & The zonal velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em v} & The meridional velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\
\hline
\mbox{\tt in}  & {\em fluxes} & A structure of thermodynamic surface fluxes\\
\hline
\mbox{\tt in}  & {\em visc} & Structure containing vertical viscosities, bottom boundary layer properies, and related fields.\\
\hline
 & {\em cs} & Diffusivity control structure\\
\hline
 & {\em obc} & Open boundaries control structure. \\
\hline
\end{DoxyParams}


Definition at line 1707 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
1707   \textcolor{keywordtype}{type}(ocean\_grid\_type),    \textcolor{keywordtype}{intent(in)}    :: g\textcolor{comment}{    !< The ocean's grid structure}
1708   \textcolor{keywordtype}{type}(verticalgrid\_type),  \textcolor{keywordtype}{intent(in)}    :: gv\textcolor{comment}{   !< The ocean's vertical grid structure}
1709   \textcolor{keywordtype}{type}(unit\_scale\_type),    \textcolor{keywordtype}{intent(in)}    :: us\textcolor{comment}{   !< A dimensional unit scaling type}
1710   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, &
1711                             \textcolor{keywordtype}{intent(in)}    :: u\textcolor{comment}{    !< The zonal velocity [L T-1 ~> m s-1]}
1712   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, &
1713                             \textcolor{keywordtype}{intent(in)}    :: v\textcolor{comment}{    !< The meridional velocity [L T-1 ~> m s-1]}
1714   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1715                             \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{    !< Layer thicknesses [H ~> m or kg m-2]}
1716   \textcolor{keywordtype}{type}(forcing),            \textcolor{keywordtype}{intent(in)}    :: fluxes\textcolor{comment}{ !< A structure of thermodynamic surface fluxes}
1717   \textcolor{keywordtype}{type}(vertvisc\_type),      \textcolor{keywordtype}{intent(in)}    :: visc\textcolor{comment}{ !< Structure containing vertical viscosities, bottom}
1718 \textcolor{comment}{                                                  !! boundary layer properies, and related fields.}
1719   \textcolor{keywordtype}{type}(set\_diffusivity\_cs), \textcolor{keywordtype}{pointer}       :: cs\textcolor{comment}{   !< Diffusivity control structure}
1720   \textcolor{keywordtype}{type}(ocean\_obc\_type), \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: obc\textcolor{comment}{  !< Open boundaries control structure.}
1721 
1722   \textcolor{comment}{! This subroutine calculates several properties related to bottom}
1723   \textcolor{comment}{! boundary layer turbulence.}
1724 
1725   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: &
1726     htot          \textcolor{comment}{! total thickness above or below a layer, or the}
1727                   \textcolor{comment}{! integrated thickness in the BBL [Z ~> m].}
1728 
1729   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G))} :: &
1730     uhtot, &      \textcolor{comment}{! running integral of u in the BBL [Z L T-1 ~> m2 s-1]}
1731     ustar, &      \textcolor{comment}{! bottom boundary layer turbulence speed [Z T-1 ~> m s-1].}
1732     u2\_bbl        \textcolor{comment}{! square of the mean zonal velocity in the BBL [L2 T-2 ~> m2 s-2]}
1733 
1734   \textcolor{keywordtype}{real} :: vhtot(szi\_(g)) \textcolor{comment}{! running integral of v in the BBL [Z L T-1 ~> m2 s-1]}
1735 
1736   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G))} :: &
1737     vstar, & \textcolor{comment}{! ustar at at v-points [Z T-1 ~> m s-1].}
1738     v2\_bbl   \textcolor{comment}{! square of average meridional velocity in BBL [L2 T-2 ~> m2 s-2]}
1739 
1740   \textcolor{keywordtype}{real} :: cdrag\_sqrt  \textcolor{comment}{! square root of the drag coefficient [nondim]}
1741   \textcolor{keywordtype}{real} :: hvel        \textcolor{comment}{! thickness at velocity points [Z ~> m].}
1742 
1743   \textcolor{keywordtype}{logical} :: domore, do\_i(szi\_(g))
1744   \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, nz
1745   \textcolor{keywordtype}{integer} :: l\_seg
1746   \textcolor{keywordtype}{logical} :: local\_open\_u\_bc, local\_open\_v\_bc
1747   \textcolor{keywordtype}{logical} :: has\_obc
1748 
1749   local\_open\_u\_bc = .false.
1750   local\_open\_v\_bc = .false.
1751   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then}
1752     local\_open\_u\_bc = obc%open\_u\_BCs\_exist\_globally
1753     local\_open\_v\_bc = obc%open\_v\_BCs\_exist\_globally
1754 \textcolor{keywordflow}{  endif} ;\textcolor{keywordflow}{ endif}
1755 
1756   is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
1757 
1758   \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal,\textcolor{stringliteral}{"set\_BBL\_TKE: "}//&
1759          \textcolor{stringliteral}{"Module must be initialized before it is used."})
1760 
1761   \textcolor{keywordflow}{if} (.not.cs%bottomdraglaw .or. (cs%BBL\_effic<=0.0)) \textcolor{keywordflow}{then}
1762     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%ustar\_BBL)) \textcolor{keywordflow}{then}
1763       \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=is,ie ; visc%ustar\_BBL(i,j) = 0.0 ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo}
1764 \textcolor{keywordflow}{    endif}
1765     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%TKE\_BBL)) \textcolor{keywordflow}{then}
1766       \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=is,ie ; visc%TKE\_BBL(i,j) = 0.0 ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo}
1767 \textcolor{keywordflow}{    endif}
1768     \textcolor{keywordflow}{return}
1769 \textcolor{keywordflow}{  endif}
1770 
1771   cdrag\_sqrt = sqrt(cs%cdrag)
1772 
1773   \textcolor{comment}{!$OMP parallel default(shared) private(do\_i,vhtot,htot,domore,hvel,uhtot,ustar,u2\_bbl)}
1774   \textcolor{comment}{!$OMP do}
1775   \textcolor{keywordflow}{do} j=js-1,je
1776     \textcolor{comment}{! Determine ustar and the square magnitude of the velocity in the}
1777     \textcolor{comment}{! bottom boundary layer. Together these give the TKE source and}
1778     \textcolor{comment}{! vertical decay scale.}
1779     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} ((g%mask2dCv(i,j) > 0.5) .and. (cdrag\_sqrt*visc%bbl\_thick\_v(i,j) > 0.0)) \textcolor{keywordflow}{then}
1780       do\_i(i) = .true. ; vhtot(i) = 0.0 ; htot(i) = 0.0
1781       vstar(i,j) = visc%Kv\_bbl\_v(i,j) / (cdrag\_sqrt*visc%bbl\_thick\_v(i,j))
1782     \textcolor{keywordflow}{else}
1783       do\_i(i) = .false. ; vstar(i,j) = 0.0 ; htot(i) = 0.0
1784 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1785     \textcolor{keywordflow}{do} k=nz,1,-1
1786       domore = .false.
1787       \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1788         \textcolor{comment}{! Determine if grid point is an OBC}
1789         has\_obc = .false.
1790         \textcolor{keywordflow}{if} (local\_open\_v\_bc) \textcolor{keywordflow}{then}
1791           l\_seg = obc%segnum\_v(i,j)
1792           \textcolor{keywordflow}{if} (l\_seg /= obc\_none) \textcolor{keywordflow}{then}
1793             has\_obc = obc%segment(l\_seg)%open
1794 \textcolor{keywordflow}{          endif}
1795 \textcolor{keywordflow}{        endif}
1796 
1797         \textcolor{comment}{! Compute h based on OBC state}
1798         \textcolor{keywordflow}{if} (has\_obc) \textcolor{keywordflow}{then}
1799           \textcolor{keywordflow}{if} (obc%segment(l\_seg)%direction == obc\_direction\_n) \textcolor{keywordflow}{then}
1800             hvel = gv%H\_to\_Z*h(i,j,k)
1801           \textcolor{keywordflow}{else}
1802             hvel = gv%H\_to\_Z*h(i,j+1,k)
1803 \textcolor{keywordflow}{          endif}
1804         \textcolor{keywordflow}{else}
1805           hvel = 0.5*gv%H\_to\_Z*(h(i,j,k) + h(i,j+1,k))
1806 \textcolor{keywordflow}{        endif}
1807 
1808         \textcolor{keywordflow}{if} ((htot(i) + hvel) >= visc%bbl\_thick\_v(i,j)) \textcolor{keywordflow}{then}
1809           vhtot(i) = vhtot(i) + (visc%bbl\_thick\_v(i,j) - htot(i))*v(i,j,k)
1810           htot(i) = visc%bbl\_thick\_v(i,j)
1811           do\_i(i) = .false.
1812         \textcolor{keywordflow}{else}
1813           vhtot(i) = vhtot(i) + hvel*v(i,j,k)
1814           htot(i) = htot(i) + hvel
1815           domore = .true.
1816 \textcolor{keywordflow}{        endif}
1817 \textcolor{keywordflow}{      endif} ;\textcolor{keywordflow}{ enddo}
1818       \textcolor{keywordflow}{if} (.not.domore) \textcolor{keywordflow}{exit}
1819 \textcolor{keywordflow}{    enddo}
1820     \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} ((g%mask2dCv(i,j) > 0.5) .and. (htot(i) > 0.0)) \textcolor{keywordflow}{then}
1821       v2\_bbl(i,j) = (vhtot(i)*vhtot(i))/(htot(i)*htot(i))
1822     \textcolor{keywordflow}{else}
1823       v2\_bbl(i,j) = 0.0
1824 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1825 \textcolor{keywordflow}{  enddo}
1826   \textcolor{comment}{!$OMP do}
1827   \textcolor{keywordflow}{do} j=js,je
1828     \textcolor{keywordflow}{do} i=is-1,ie ; \textcolor{keywordflow}{if} ((g%mask2dCu(i,j) > 0.5) .and. (cdrag\_sqrt*visc%bbl\_thick\_u(i,j) > 0.0))  \textcolor{keywordflow}{then}
1829       do\_i(i) = .true. ; uhtot(i) = 0.0 ; htot(i) = 0.0
1830       ustar(i) = visc%Kv\_bbl\_u(i,j) / (cdrag\_sqrt*visc%bbl\_thick\_u(i,j))
1831     \textcolor{keywordflow}{else}
1832       do\_i(i) = .false. ; ustar(i) = 0.0 ; htot(i) = 0.0
1833 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1834     \textcolor{keywordflow}{do} k=nz,1,-1 ; domore = .false.
1835       \textcolor{keywordflow}{do} i=is-1,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then}
1836         \textcolor{comment}{! Determine if grid point is an OBC}
1837         has\_obc = .false.
1838         \textcolor{keywordflow}{if} (local\_open\_u\_bc) \textcolor{keywordflow}{then}
1839           l\_seg = obc%segnum\_u(i,j)
1840           \textcolor{keywordflow}{if} (l\_seg /= obc\_none) \textcolor{keywordflow}{then}
1841             has\_obc = obc%segment(l\_seg)%open
1842 \textcolor{keywordflow}{          endif}
1843 \textcolor{keywordflow}{        endif}
1844 
1845         \textcolor{comment}{! Compute h based on OBC state}
1846         \textcolor{keywordflow}{if} (has\_obc) \textcolor{keywordflow}{then}
1847           \textcolor{keywordflow}{if} (obc%segment(l\_seg)%direction == obc\_direction\_e) \textcolor{keywordflow}{then}
1848             hvel = gv%H\_to\_Z*h(i,j,k)
1849           \textcolor{keywordflow}{else} \textcolor{comment}{! OBC\_DIRECTION\_W}
1850             hvel = gv%H\_to\_Z*h(i+1,j,k)
1851 \textcolor{keywordflow}{          endif}
1852         \textcolor{keywordflow}{else}
1853           hvel = 0.5*gv%H\_to\_Z*(h(i,j,k) + h(i+1,j,k))
1854 \textcolor{keywordflow}{        endif}
1855 
1856         \textcolor{keywordflow}{if} ((htot(i) + hvel) >= visc%bbl\_thick\_u(i,j)) \textcolor{keywordflow}{then}
1857           uhtot(i) = uhtot(i) + (visc%bbl\_thick\_u(i,j) - htot(i))*u(i,j,k)
1858           htot(i) = visc%bbl\_thick\_u(i,j)
1859           do\_i(i) = .false.
1860         \textcolor{keywordflow}{else}
1861           uhtot(i) = uhtot(i) + hvel*u(i,j,k)
1862           htot(i) = htot(i) + hvel
1863           domore = .true.
1864 \textcolor{keywordflow}{        endif}
1865 \textcolor{keywordflow}{      endif} ;\textcolor{keywordflow}{ enddo}
1866       \textcolor{keywordflow}{if} (.not.domore) \textcolor{keywordflow}{exit}
1867 \textcolor{keywordflow}{    enddo}
1868     \textcolor{keywordflow}{do} i=is-1,ie ; \textcolor{keywordflow}{if} ((g%mask2dCu(i,j) > 0.5) .and. (htot(i) > 0.0)) \textcolor{keywordflow}{then}
1869       u2\_bbl(i) = (uhtot(i)*uhtot(i))/(htot(i)*htot(i))
1870     \textcolor{keywordflow}{else}
1871       u2\_bbl(i) = 0.0
1872 \textcolor{keywordflow}{    endif} ;\textcolor{keywordflow}{ enddo}
1873 
1874     \textcolor{keywordflow}{do} i=is,ie
1875       visc%ustar\_BBL(i,j) = sqrt(0.5*g%IareaT(i,j) * &
1876                 ((g%areaCu(i-1,j)*(ustar(i-1)*ustar(i-1)) + &
1877                   g%areaCu(i,j)*(ustar(i)*ustar(i))) + &
1878                  (g%areaCv(i,j-1)*(vstar(i,j-1)*vstar(i,j-1)) + &
1879                   g%areaCv(i,j)*(vstar(i,j)*vstar(i,j))) ) )
1880       visc%TKE\_BBL(i,j) = us%L\_to\_Z**2 * &
1881                  (((g%areaCu(i-1,j)*(ustar(i-1)*u2\_bbl(i-1)) + &
1882                     g%areaCu(i,j) * (ustar(i)*u2\_bbl(i))) + &
1883                    (g%areaCv(i,j-1)*(vstar(i,j-1)*v2\_bbl(i,j-1)) + &
1884                     g%areaCv(i,j) * (vstar(i,j)*v2\_bbl(i,j))) )*g%IareaT(i,j))
1885 \textcolor{keywordflow}{    enddo}
1886 \textcolor{keywordflow}{  enddo}
1887   \textcolor{comment}{!$OMP end parallel}
1888 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a5ba8a3be6234304aa5f1dfd0b831078a}\label{namespacemom__set__diffusivity_a5ba8a3be6234304aa5f1dfd0b831078a}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!set\+\_\+density\+\_\+ratios@{set\+\_\+density\+\_\+ratios}}
\index{set\+\_\+density\+\_\+ratios@{set\+\_\+density\+\_\+ratios}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{set\+\_\+density\+\_\+ratios()}{set\_density\_ratios()}}
{\footnotesize\ttfamily subroutine mom\+\_\+set\+\_\+diffusivity\+::set\+\_\+density\+\_\+ratios (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{h,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv,  }\item[{integer, dimension(szi\+\_\+(g)), intent(in)}]{kb,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{integer, intent(in)}]{j,  }\item[{real, dimension(szi\+\_\+(g),szk\+\_\+(g)), intent(out)}]{ds\+\_\+dsp1,  }\item[{real, dimension(szi\+\_\+(g),szk\+\_\+(g)), intent(in), optional}]{rho\+\_\+0 }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}}


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure.\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em tv} & Structure containing pointers to any available thermodynamic fields; absent fields have N\+U\+LL ptrs.\\
\hline
\mbox{\tt in}  & {\em kb} & Index of lightest layer denser than the buffer layer, or -\/1 without a bulk mixed layer.\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
 & {\em cs} & Control structure returned by previous call to diabatic\+\_\+entrain\+\_\+init.\\
\hline
\mbox{\tt in}  & {\em j} & Meridional index upon which to work.\\
\hline
\mbox{\tt out}  & {\em ds\+\_\+dsp1} & Coordinate variable (sigma-\/2) difference across an interface divided by the difference across the interface below it \mbox{[}nondim\mbox{]}\\
\hline
\mbox{\tt in}  & {\em rho\+\_\+0} & Layer potential densities relative to \\
\hline
\end{DoxyParams}


Definition at line 1892 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
1892   \textcolor{keywordtype}{type}(ocean\_grid\_type),            \textcolor{keywordtype}{intent(in)}   :: g\textcolor{comment}{  !< The ocean's grid structure.}
1893   \textcolor{keywordtype}{type}(verticalgrid\_type),          \textcolor{keywordtype}{intent(in)}   :: gv\textcolor{comment}{ !< The ocean's vertical grid structure.}
1894   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
1895                                     \textcolor{keywordtype}{intent(in)}   :: h\textcolor{comment}{  !< Layer thicknesses [H ~> m or kg m-2].}
1896   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),            \textcolor{keywordtype}{intent(in)}   :: tv\textcolor{comment}{ !< Structure containing pointers to any}
1897 \textcolor{comment}{                                                       !! available thermodynamic fields; absent}
1898 \textcolor{comment}{                                                       !! fields have NULL ptrs.}
1899   \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(SZI\_(G))},      \textcolor{keywordtype}{intent(in)}   :: kb\textcolor{comment}{ !< Index of lightest layer denser than the buffer}
1900 \textcolor{comment}{                                                       !! layer, or -1 without a bulk mixed layer.}
1901   \textcolor{keywordtype}{type}(unit\_scale\_type),            \textcolor{keywordtype}{intent(in)}   :: us\textcolor{comment}{ !< A dimensional unit scaling type}
1902   \textcolor{keywordtype}{type}(set\_diffusivity\_cs),         \textcolor{keywordtype}{pointer}      :: cs\textcolor{comment}{ !< Control structure returned by previous}
1903 \textcolor{comment}{                                                       !! call to diabatic\_entrain\_init.}
1904   \textcolor{keywordtype}{integer},                          \textcolor{keywordtype}{intent(in)}   :: j\textcolor{comment}{  !< Meridional index upon which to work.}
1905   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)}  :: ds\_dsp1\textcolor{comment}{ !< Coordinate variable (sigma-2)}
1906 \textcolor{comment}{                                                       !! difference across an interface divided by}
1907 \textcolor{comment}{                                                       !! the difference across the interface below}
1908 \textcolor{comment}{                                                       !! it [nondim]}
1909   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))}, &
1910                           \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)}   :: rho\_0\textcolor{comment}{ !< Layer potential densities relative to}
1911 \textcolor{comment}{                                                       !! surface press [R ~> kg m-3].}
1912 
1913   \textcolor{comment}{! Local variables}
1914   \textcolor{keywordtype}{real} :: g\_r0                     \textcolor{comment}{! g\_R0 is a rescaled version of g/Rho [L2 Z-1 R-1 T-2 ~> m4 kg-1 s-2]}
1915   \textcolor{keywordtype}{real} :: eps, tmp                 \textcolor{comment}{! nondimensional temproray variables}
1916   \textcolor{keywordtype}{real} :: a(szk\_(g)), a\_0(szk\_(g)) \textcolor{comment}{! nondimensional temporary variables}
1917   \textcolor{keywordtype}{real} :: p\_ref(szi\_(g))           \textcolor{comment}{! an array of tv%P\_Ref pressures [R L2 T-2 ~> Pa]}
1918   \textcolor{keywordtype}{real} :: rcv(szi\_(g),szk\_(g))     \textcolor{comment}{! coordinate density in the mixed and buffer layers [R ~> kg m-3]}
1919   \textcolor{keywordtype}{real} :: i\_drho                   \textcolor{comment}{! temporary variable [R-1 ~> m3 kg-1]}
1920 
1921   \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The i-computational domain for the equation of state}
1922   \textcolor{keywordtype}{integer} :: i, k, k3, is, ie, nz, kmb
1923   is = g%isc ; ie = g%iec ; nz = g%ke
1924 
1925   \textcolor{keywordflow}{do} k=2,nz-1
1926     \textcolor{keywordflow}{if} (gv%g\_prime(k+1) /= 0.0) \textcolor{keywordflow}{then}
1927       \textcolor{keywordflow}{do} i=is,ie
1928         ds\_dsp1(i,k) = gv%g\_prime(k) / gv%g\_prime(k+1)
1929 \textcolor{keywordflow}{      enddo}
1930     \textcolor{keywordflow}{else}
1931       \textcolor{keywordflow}{do} i=is,ie
1932         ds\_dsp1(i,k) = 1.
1933 \textcolor{keywordflow}{      enddo}
1934 \textcolor{keywordflow}{    endif}
1935 \textcolor{keywordflow}{  enddo}
1936 
1937   \textcolor{keywordflow}{if} (cs%bulkmixedlayer) \textcolor{keywordflow}{then}
1938     g\_r0 = gv%g\_Earth / (gv%Rho0)
1939     kmb = gv%nk\_rho\_varies
1940     eps = 0.1
1941     \textcolor{keywordflow}{do} i=is,ie ; p\_ref(i) = tv%P\_Ref ;\textcolor{keywordflow}{ enddo}
1942     eosdom(:) = eos\_domain(g%HI)
1943     \textcolor{keywordflow}{do} k=1,kmb
1944       \textcolor{keyword}{call }calculate\_density(tv%T(:,j,k), tv%S(:,j,k), p\_ref, rcv(:,k), tv%eqn\_of\_state, eosdom)
1945 \textcolor{keywordflow}{    enddo}
1946     \textcolor{keywordflow}{do} i=is,ie
1947       \textcolor{keywordflow}{if} (kb(i) <= nz-1) \textcolor{keywordflow}{then}
1948 \textcolor{comment}{!   Set up appropriately limited ratios of the reduced gravities of the}
1949 \textcolor{comment}{! interfaces above and below the buffer layer and the next denser layer.}
1950         k = kb(i)
1951 
1952         i\_drho = g\_r0 / gv%g\_prime(k+1)
1953         \textcolor{comment}{! The indexing convention for a is appropriate for the interfaces.}
1954         \textcolor{keywordflow}{do} k3=1,kmb
1955           a(k3+1) = (gv%Rlay(k) - rcv(i,k3)) * i\_drho
1956 \textcolor{keywordflow}{        enddo}
1957         \textcolor{keywordflow}{if} ((\textcolor{keyword}{present}(rho\_0)) .and. (a(kmb+1) < 2.0*eps*ds\_dsp1(i,k))) \textcolor{keywordflow}{then}
1958 \textcolor{comment}{!   If the buffer layer nearly matches the density of the layer below in the}
1959 \textcolor{comment}{! coordinate variable (sigma-2), use the sigma-0-based density ratio if it is}
1960 \textcolor{comment}{! greater (and stable).}
1961           \textcolor{keywordflow}{if} ((rho\_0(i,k) > rho\_0(i,kmb)) .and. &
1962               (rho\_0(i,k+1) > rho\_0(i,k))) \textcolor{keywordflow}{then}
1963             i\_drho = 1.0 / (rho\_0(i,k+1)-rho\_0(i,k))
1964             a\_0(kmb+1) = min((rho\_0(i,k)-rho\_0(i,kmb)) * i\_drho, ds\_dsp1(i,k))
1965             \textcolor{keywordflow}{if} (a\_0(kmb+1) > a(kmb+1)) \textcolor{keywordflow}{then}
1966               \textcolor{keywordflow}{do} k3=2,kmb
1967                 a\_0(k3) = a\_0(kmb+1) + (rho\_0(i,kmb)-rho\_0(i,k3-1)) * i\_drho
1968 \textcolor{keywordflow}{              enddo}
1969               \textcolor{keywordflow}{if} (a(kmb+1) <= eps*ds\_dsp1(i,k)) \textcolor{keywordflow}{then}
1970                 \textcolor{keywordflow}{do} k3=2,kmb+1 ; a(k3) = a\_0(k3) ;\textcolor{keywordflow}{ enddo}
1971               \textcolor{keywordflow}{else}
1972 \textcolor{comment}{! Alternative...  tmp = 0.5*(1.0 - cos(PI*(a(K2+1)/(eps*ds\_dsp1(i,k)) - 1.0)) )}
1973                 tmp = a(kmb+1)/(eps*ds\_dsp1(i,k)) - 1.0
1974                 \textcolor{keywordflow}{do} k3=2,kmb+1 ; a(k3) = tmp*a(k3) + (1.0-tmp)*a\_0(k3) ;\textcolor{keywordflow}{ enddo}
1975 \textcolor{keywordflow}{              endif}
1976 \textcolor{keywordflow}{            endif}
1977 \textcolor{keywordflow}{          endif}
1978 \textcolor{keywordflow}{        endif}
1979 
1980         ds\_dsp1(i,k) = max(a(kmb+1),1e-5)
1981 
1982         \textcolor{keywordflow}{do} k3=2,kmb
1983 \textcolor{comment}{!           ds\_dsp1(i,k3) = MAX(a(k3),1e-5)}
1984           \textcolor{comment}{! Deliberately treat convective instabilies of the upper mixed}
1985           \textcolor{comment}{! and buffer layers with respect to the deepest buffer layer as}
1986           \textcolor{comment}{! though they don't exist.  They will be eliminated by the upcoming}
1987           \textcolor{comment}{! call to the mixedlayer code anyway.}
1988           \textcolor{comment}{! The indexing convention is appropriate for the interfaces.}
1989           ds\_dsp1(i,k3) = max(a(k3),ds\_dsp1(i,k))
1990 \textcolor{keywordflow}{        enddo}
1991 \textcolor{keywordflow}{      endif} \textcolor{comment}{! (kb(i) <= nz-1)}
1992 \textcolor{keywordflow}{    enddo} \textcolor{comment}{! I-loop.}
1993 \textcolor{keywordflow}{  endif} \textcolor{comment}{! bulkmixedlayer}
1994 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a7c293162d6c8efb882c8b04b4ea5241d}\label{namespacemom__set__diffusivity_a7c293162d6c8efb882c8b04b4ea5241d}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!set\+\_\+diffusivity@{set\+\_\+diffusivity}}
\index{set\+\_\+diffusivity@{set\+\_\+diffusivity}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{set\+\_\+diffusivity()}{set\_diffusivity()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+diffusivity\+::set\+\_\+diffusivity (\begin{DoxyParamCaption}\item[{real, dimension(szib\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{u,  }\item[{real, dimension(szi\+\_\+(g),szjb\+\_\+(g),szk\+\_\+(g)), intent(in)}]{v,  }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{h,  }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{u\+\_\+h,  }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{v\+\_\+h,  }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(inout)}]{tv,  }\item[{type(forcing), intent(in)}]{fluxes,  }\item[{type(optics\+\_\+type), pointer}]{optics,  }\item[{type(vertvisc\+\_\+type), intent(inout)}]{visc,  }\item[{real, intent(in)}]{dt,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(out), optional}]{Kd\+\_\+lay,  }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)+1), intent(out), optional}]{Kd\+\_\+int }\end{DoxyParamCaption})}



Sets the interior vertical diffusion of scalars due to the following processes\+: 


\begin{DoxyEnumerate}
\item Shear-\/driven mixing\+: two options, Jackson et at. and K\+PP interior;
\item Background mixing via C\+V\+Mix (Bryan-\/\+Lewis profile) or the scheme described by Harrison \& Hallberg, J\+PO 2008;
\item Double-\/diffusion, old method and new method via C\+V\+Mix;
\item Tidal mixing\+: many options available, see \hyperlink{MOM__tidal__mixing_8F90_source}{M\+O\+M\+\_\+tidal\+\_\+mixing.\+F90}; In addition, this subroutine has the option to set the interior vertical viscosity associated with processes 1,2 and 4 listed above, which is stored in viscKv\+\_\+slow. Vertical viscosity due to shear-\/driven mixing is passed via viscKv\+\_\+shear
\end{DoxyEnumerate}


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em g} & The ocean\textquotesingle{}s grid structure.\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em u} & The zonal velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em v} & The meridional velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em u\+\_\+h} & Zonal velocity interpolated to h points \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\
\hline
\mbox{\tt in}  & {\em v\+\_\+h} & Meridional velocity interpolated to h points \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\
\hline
\mbox{\tt in,out}  & {\em tv} & Structure with pointers to thermodynamic fields. Out is for tvTempx\+PmE.\\
\hline
\mbox{\tt in}  & {\em fluxes} & A structure of thermodynamic surface fluxes\\
\hline
 & {\em optics} & A structure describing the optical properties of the ocean.\\
\hline
\mbox{\tt in,out}  & {\em visc} & Structure containing vertical viscosities, bottom boundary layer properies, and related fields.\\
\hline
\mbox{\tt in}  & {\em dt} & Time increment \mbox{[}T $\sim$$>$ s\mbox{]}.\\
\hline
 & {\em cs} & Module control structure.\\
\hline
\mbox{\tt out}  & {\em kd\+\_\+lay} & Diapycnal diffusivity of each layer \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\
\hline
\mbox{\tt out}  & {\em kd\+\_\+int} & Diapycnal diffusivity at each interface \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}. \\
\hline
\end{DoxyParams}


Definition at line 213 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
213   \textcolor{keywordtype}{type}(ocean\_grid\_type),     \textcolor{keywordtype}{intent(in)}    :: g\textcolor{comment}{    !< The ocean's grid structure.}
214   \textcolor{keywordtype}{type}(verticalgrid\_type),   \textcolor{keywordtype}{intent(in)}    :: gv\textcolor{comment}{   !< The ocean's vertical grid structure.}
215   \textcolor{keywordtype}{type}(unit\_scale\_type),     \textcolor{keywordtype}{intent(in)}    :: us\textcolor{comment}{   !< A dimensional unit scaling type}
216   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, &
217                              \textcolor{keywordtype}{intent(in)}    :: u\textcolor{comment}{    !< The zonal velocity [L T-1 ~> m s-1].}
218   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, &
219                              \textcolor{keywordtype}{intent(in)}    :: v\textcolor{comment}{    !< The meridional velocity [L T-1 ~> m s-1].}
220   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))},  &
221                              \textcolor{keywordtype}{intent(in)}    :: h\textcolor{comment}{    !< Layer thicknesses [H ~> m or kg m-2].}
222   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))},  &
223                              \textcolor{keywordtype}{intent(in)}    :: u\_h\textcolor{comment}{  !< Zonal velocity interpolated to h points [L T-1 ~> m
       s-1].}
224   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))},  &
225                              \textcolor{keywordtype}{intent(in)}    :: v\_h\textcolor{comment}{  !< Meridional velocity interpolated to h points [L T-1
       ~> m s-1].}
226   \textcolor{keywordtype}{type}(thermo\_var\_ptrs),     \textcolor{keywordtype}{intent(inout)} :: tv\textcolor{comment}{   !< Structure with pointers to thermodynamic}
227 \textcolor{comment}{                                                   !! fields. Out is for tv%TempxPmE.}
228   \textcolor{keywordtype}{type}(forcing),             \textcolor{keywordtype}{intent(in)}    :: fluxes\textcolor{comment}{ !< A structure of thermodynamic surface fluxes}
229   \textcolor{keywordtype}{type}(optics\_type),         \textcolor{keywordtype}{pointer}       :: optics\textcolor{comment}{ !< A structure describing the optical}
230 \textcolor{comment}{                                                   !!  properties of the ocean.}
231   \textcolor{keywordtype}{type}(vertvisc\_type),       \textcolor{keywordtype}{intent(inout)} :: visc\textcolor{comment}{ !< Structure containing vertical viscosities, bottom}
232 \textcolor{comment}{                                                   !! boundary layer properies, and related fields.}
233   \textcolor{keywordtype}{real},                      \textcolor{keywordtype}{intent(in)}    :: dt\textcolor{comment}{   !< Time increment [T ~> s].}
234   \textcolor{keywordtype}{type}(set\_diffusivity\_cs),  \textcolor{keywordtype}{pointer}       :: cs\textcolor{comment}{   !< Module control structure.}
235   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, &
236                    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)}   :: kd\_lay\textcolor{comment}{ !< Diapycnal diffusivity of each layer [Z2 T-1 ~> m2
       s-1].}
237   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)}, &
238                    \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)}   :: kd\_int\textcolor{comment}{ !< Diapycnal diffusivity at each interface [Z2 T-1 ~>
       m2 s-1].}
239 
240   \textcolor{comment}{! local variables}
241   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: &
242     n2\_bot        \textcolor{comment}{! bottom squared buoyancy frequency [T-2 ~> s-2]}
243 
244   \textcolor{keywordtype}{type}(diffusivity\_diags)  :: dd \textcolor{comment}{! structure with arrays of available diags}
245 
246   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))} :: &
247     t\_f, s\_f      \textcolor{comment}{! Temperature and salinity [degC] and [ppt] with properties in massless layers}
248                   \textcolor{comment}{! filled vertically by diffusion or the properties after full convective adjustment.}
249 
250   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G))} :: &
251     n2\_lay, &     !< Squared buoyancy frequency associated with layers [T-2 ~> s-2]
252     kd\_lay\_2d, &  !< The layer diffusivities [Z2 T-1 ~> m2 s-1]
253     maxtke, &     !< Energy required to entrain to h\_max [Z3 T-3 ~> m3 s-3]
254     tke\_to\_kd\textcolor{comment}{     !< Conversion rate (~1.0 / (G\_Earth + dRho\_lay)) between}
255 \textcolor{comment}{                  !< TKE dissipated within a layer and Kd in that layer}
256 \textcolor{comment}{                  !< [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]}
257 
258   \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(G)+1)} :: &
259     n2\_int,   &   !< squared buoyancy frequency associated at interfaces [T-2 ~> s-2]
260     kd\_int\_2d, &  !< The interface diffusivities [Z2 T-1 ~> m2 s-1]
261     kv\_bkgnd, &   !< The background diffusion related interface viscosities [Z2 T-1 ~> m2 s-1]
262     drho\_int, &   !< Locally referenced potential density difference across interfaces [R ~> kg m-3]
263     kt\_extra, &   !< Double difusion diffusivity of temperature [Z2 T-1 ~> m2 s-1]
264     ks\_extra\textcolor{comment}{      !< Double difusion diffusivity of salinity [Z2 T-1 ~> m2 s-1]}
265 
266   \textcolor{keywordtype}{real} :: dissip        \textcolor{comment}{! local variable for dissipation calculations [Z2 R T-3 ~> W m-3]}
267   \textcolor{keywordtype}{real} :: omega2        \textcolor{comment}{! squared absolute rotation rate [T-2 ~> s-2]}
268 
269   \textcolor{keywordtype}{logical}   :: use\_eos      \textcolor{comment}{! If true, compute density from T/S using equation of state.}
270   \textcolor{keywordtype}{logical}   :: tke\_to\_kd\_used \textcolor{comment}{! If true, TKE\_to\_Kd and maxTKE need to be calculated.}
271   \textcolor{keywordtype}{integer}   :: kb(szi\_(g))  \textcolor{comment}{! The index of the lightest layer denser than the}
272                             \textcolor{comment}{! buffer layer, or -1 without a bulk mixed layer.}
273   \textcolor{keywordtype}{logical}   :: showcalltree \textcolor{comment}{! If true, show the call tree.}
274 
275   \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, nz, isd, ied, jsd, jed
276 
277   \textcolor{keywordtype}{real}      :: kappa\_dt\_fill \textcolor{comment}{! diffusivity times a timestep used to fill massless layers [Z2 ~> m2]}
278 
279   is  = g%isc ; ie  = g%iec ; js  = g%jsc ; je  = g%jec ; nz = g%ke
280   isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
281   showcalltree = calltree\_showquery()
282   \textcolor{keywordflow}{if} (showcalltree) \textcolor{keyword}{call }calltree\_enter(\textcolor{stringliteral}{"set\_diffusivity(), MOM\_set\_diffusivity.F90"})
283 
284   \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal,\textcolor{stringliteral}{"set\_diffusivity: "}//&
285          \textcolor{stringliteral}{"Module must be initialized before it is used."})
286 
287   \textcolor{keywordflow}{if} (cs%answers\_2018) \textcolor{keywordflow}{then}
288     \textcolor{comment}{! These hard-coded dimensional parameters are being replaced.}
289     kappa\_dt\_fill = us%m\_to\_Z**2 * 1.e-3 * 7200.
290   \textcolor{keywordflow}{else}
291     kappa\_dt\_fill = cs%Kd\_smooth * dt
292 \textcolor{keywordflow}{  endif}
293   omega2 = cs%omega * cs%omega
294 
295   use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state)
296 
297   \textcolor{keywordflow}{if} ((cs%use\_CVMix\_ddiff .or. cs%double\_diffusion) .and. .not. &
298      (\textcolor{keyword}{associated}(visc%Kd\_extra\_T) .and. \textcolor{keyword}{associated}(visc%Kd\_extra\_S))) &
299      \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"set\_diffusivity: both visc%Kd\_extra\_T and "}//&
300          \textcolor{stringliteral}{"visc%Kd\_extra\_S must be associated when USE\_CVMIX\_DDIFF or DOUBLE\_DIFFUSION are true."})
301 
302   tke\_to\_kd\_used = (cs%use\_tidal\_mixing .or. cs%ML\_radiation .or. &
303                    (cs%bottomdraglaw .and. .not.cs%use\_LOTW\_BBL\_diffusivity))
304 
305   \textcolor{comment}{! Set Kd\_lay, Kd\_int and Kv\_slow to constant values, mostly to fill the halos.}
306   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_lay)) kd\_lay(:,:,:) = cs%Kd
307   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) kd\_int(:,:,:) = cs%Kd
308   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_slow)) visc%Kv\_slow(:,:,:) = cs%Kv
309 
310   \textcolor{comment}{! Set up arrays for diagnostics.}
311 
312   \textcolor{keywordflow}{if} (cs%id\_N2 > 0) \textcolor{keywordflow}{then}
313     \textcolor{keyword}{allocate}(dd%N2\_3d(isd:ied,jsd:jed,nz+1)) ; dd%N2\_3d(:,:,:) = 0.0
314 \textcolor{keywordflow}{  endif}
315   \textcolor{keywordflow}{if} (cs%id\_Kd\_user > 0) \textcolor{keywordflow}{then}
316     \textcolor{keyword}{allocate}(dd%Kd\_user(isd:ied,jsd:jed,nz+1)) ; dd%Kd\_user(:,:,:) = 0.0
317 \textcolor{keywordflow}{  endif}
318   \textcolor{keywordflow}{if} (cs%id\_Kd\_work > 0) \textcolor{keywordflow}{then}
319     \textcolor{keyword}{allocate}(dd%Kd\_work(isd:ied,jsd:jed,nz)) ; dd%Kd\_work(:,:,:) = 0.0
320 \textcolor{keywordflow}{  endif}
321   \textcolor{keywordflow}{if} (cs%id\_maxTKE > 0) \textcolor{keywordflow}{then}
322     \textcolor{keyword}{allocate}(dd%maxTKE(isd:ied,jsd:jed,nz)) ; dd%maxTKE(:,:,:) = 0.0
323 \textcolor{keywordflow}{  endif}
324   \textcolor{keywordflow}{if} (cs%id\_TKE\_to\_Kd > 0) \textcolor{keywordflow}{then}
325     \textcolor{keyword}{allocate}(dd%TKE\_to\_Kd(isd:ied,jsd:jed,nz)) ; dd%TKE\_to\_Kd(:,:,:) = 0.0
326 \textcolor{keywordflow}{  endif}
327   \textcolor{keywordflow}{if} ((cs%double\_diffusion) .and. (cs%id\_KT\_extra > 0)) \textcolor{keywordflow}{then}
328     \textcolor{keyword}{allocate}(dd%KT\_extra(isd:ied,jsd:jed,nz+1)) ; dd%KT\_extra(:,:,:) = 0.0
329 \textcolor{keywordflow}{  endif}
330   \textcolor{keywordflow}{if} ((cs%double\_diffusion) .and. (cs%id\_KS\_extra > 0)) \textcolor{keywordflow}{then}
331     \textcolor{keyword}{allocate}(dd%KS\_extra(isd:ied,jsd:jed,nz+1)) ; dd%KS\_extra(:,:,:) = 0.0
332 \textcolor{keywordflow}{  endif}
333   \textcolor{keywordflow}{if} (cs%id\_R\_rho > 0) \textcolor{keywordflow}{then}
334     \textcolor{keyword}{allocate}(dd%drho\_rat(isd:ied,jsd:jed,nz+1)) ; dd%drho\_rat(:,:,:) = 0.0
335 \textcolor{keywordflow}{  endif}
336   \textcolor{keywordflow}{if} (cs%id\_Kd\_BBL > 0) \textcolor{keywordflow}{then}
337     \textcolor{keyword}{allocate}(dd%Kd\_BBL(isd:ied,jsd:jed,nz+1)) ; dd%Kd\_BBL(:,:,:) = 0.0
338 \textcolor{keywordflow}{  endif}
339 
340   \textcolor{keywordflow}{if} (cs%id\_Kd\_bkgnd > 0) \textcolor{keywordflow}{then}
341     \textcolor{keyword}{allocate}(dd%Kd\_bkgnd(isd:ied,jsd:jed,nz+1)) ; dd%Kd\_bkgnd(:,:,:) = 0.
342 \textcolor{keywordflow}{  endif}
343   \textcolor{keywordflow}{if} (cs%id\_Kv\_bkgnd > 0) \textcolor{keywordflow}{then}
344     \textcolor{keyword}{allocate}(dd%Kv\_bkgnd(isd:ied,jsd:jed,nz+1)) ; dd%Kv\_bkgnd(:,:,:) = 0.
345 \textcolor{keywordflow}{  endif}
346 
347   \textcolor{comment}{! set up arrays for tidal mixing diagnostics}
348   \textcolor{keyword}{call }setup\_tidal\_diagnostics(g, cs%tidal\_mixing\_CSp)
349 
350   \textcolor{keywordflow}{if} (cs%useKappaShear) \textcolor{keywordflow}{then}
351     \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then}
352       \textcolor{keyword}{call }hchksum\_pair(\textcolor{stringliteral}{"before calc\_KS [uv]\_h"}, u\_h, v\_h, g%HI, scale=us%L\_T\_to\_m\_s)
353 \textcolor{keywordflow}{    endif}
354     \textcolor{keyword}{call }cpu\_clock\_begin(id\_clock\_kappashear)
355     \textcolor{keywordflow}{if} (cs%Vertex\_shear) \textcolor{keywordflow}{then}
356       \textcolor{keyword}{call }full\_convection(g, gv, us, h, tv, t\_f, s\_f, fluxes%p\_surf, &
357                            (gv%Z\_to\_H**2)*kappa\_dt\_fill, halo=1)
358 
359       \textcolor{keyword}{call }calc\_kappa\_shear\_vertex(u, v, h, t\_f, s\_f, tv, fluxes%p\_surf, visc%Kd\_shear, &
360                                    visc%TKE\_turb, visc%Kv\_shear\_Bu, dt, g, gv, us, cs%kappaShear\_CSp)
361       \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_shear)) visc%Kv\_shear(:,:,:) = 0.0 \textcolor{comment}{! needed for other parameterizations}
362       \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then}
363         \textcolor{keyword}{call }hchksum(visc%Kd\_shear, \textcolor{stringliteral}{"after calc\_KS\_vert visc%Kd\_shear"}, g%HI, scale=us%Z2\_T\_to\_m2\_s)
364         \textcolor{keyword}{call }bchksum(visc%Kv\_shear\_Bu, \textcolor{stringliteral}{"after calc\_KS\_vert visc%Kv\_shear\_Bu"}, g%HI, scale=us%Z2\_T\_to\_m2\_s)
365         \textcolor{keyword}{call }bchksum(visc%TKE\_turb, \textcolor{stringliteral}{"after calc\_KS\_vert visc%TKE\_turb"}, g%HI, scale=us%Z\_to\_m**2*us%s\_to\_T*
      *2)
366 \textcolor{keywordflow}{      endif}
367     \textcolor{keywordflow}{else}
368       \textcolor{comment}{! Changes: visc%Kd\_shear ;  Sets: visc%Kv\_shear and visc%TKE\_turb}
369       \textcolor{keyword}{call }calculate\_kappa\_shear(u\_h, v\_h, h, tv, fluxes%p\_surf, visc%Kd\_shear, visc%TKE\_turb, &
370                                  visc%Kv\_shear, dt, g, gv, us, cs%kappaShear\_CSp)
371       \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then}
372         \textcolor{keyword}{call }hchksum(visc%Kd\_shear, \textcolor{stringliteral}{"after calc\_KS visc%Kd\_shear"}, g%HI, scale=us%Z2\_T\_to\_m2\_s)
373         \textcolor{keyword}{call }hchksum(visc%Kv\_shear, \textcolor{stringliteral}{"after calc\_KS visc%Kv\_shear"}, g%HI, scale=us%Z2\_T\_to\_m2\_s)
374         \textcolor{keyword}{call }hchksum(visc%TKE\_turb, \textcolor{stringliteral}{"after calc\_KS visc%TKE\_turb"}, g%HI, scale=us%Z\_to\_m**2*us%s\_to\_T**2)
375 \textcolor{keywordflow}{      endif}
376 \textcolor{keywordflow}{    endif}
377     \textcolor{keyword}{call }cpu\_clock\_end(id\_clock\_kappashear)
378     \textcolor{keywordflow}{if} (showcalltree) \textcolor{keyword}{call }calltree\_waypoint(\textcolor{stringliteral}{"done with calculate\_kappa\_shear (set\_diffusivity)"})
379   \textcolor{keywordflow}{elseif} (cs%use\_CVMix\_shear) \textcolor{keywordflow}{then}
380     \textcolor{comment}{!NOTE\{BGR\}: this needs to be cleaned up.  It works in 1D case, but has not been tested outside.}
381     \textcolor{keyword}{call }calculate\_cvmix\_shear(u\_h, v\_h, h, tv, visc%Kd\_shear, visc%Kv\_shear, g, gv, us, cs%CVMix\_shear\_CSp
      )
382     \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then}
383       \textcolor{keyword}{call }hchksum(visc%Kd\_shear, \textcolor{stringliteral}{"after CVMix\_shear visc%Kd\_shear"}, g%HI, scale=us%Z2\_T\_to\_m2\_s)
384       \textcolor{keyword}{call }hchksum(visc%Kv\_shear, \textcolor{stringliteral}{"after CVMix\_shear visc%Kv\_shear"}, g%HI, scale=us%Z2\_T\_to\_m2\_s)
385 \textcolor{keywordflow}{    endif}
386   \textcolor{keywordflow}{elseif} (\textcolor{keyword}{associated}(visc%Kv\_shear)) \textcolor{keywordflow}{then}
387     visc%Kv\_shear(:,:,:) = 0.0 \textcolor{comment}{! needed if calculate\_kappa\_shear is not enabled}
388 \textcolor{keywordflow}{  endif}
389 
390   \textcolor{comment}{! Smooth the properties through massless layers.}
391   \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then}
392     \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then}
393       \textcolor{keyword}{call }hchksum(tv%T, \textcolor{stringliteral}{"before vert\_fill\_TS tv%T"},g%HI)
394       \textcolor{keyword}{call }hchksum(tv%S, \textcolor{stringliteral}{"before vert\_fill\_TS tv%S"},g%HI)
395       \textcolor{keyword}{call }hchksum(h, \textcolor{stringliteral}{"before vert\_fill\_TS h"},g%HI, scale=gv%H\_to\_m)
396 \textcolor{keywordflow}{    endif}
397     \textcolor{keyword}{call }vert\_fill\_ts(h, tv%T, tv%S, kappa\_dt\_fill, t\_f, s\_f, g, gv, larger\_h\_denom=.true.)
398     \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then}
399       \textcolor{keyword}{call }hchksum(tv%T, \textcolor{stringliteral}{"after vert\_fill\_TS tv%T"},g%HI)
400       \textcolor{keyword}{call }hchksum(tv%S, \textcolor{stringliteral}{"after vert\_fill\_TS tv%S"},g%HI)
401       \textcolor{keyword}{call }hchksum(h, \textcolor{stringliteral}{"after vert\_fill\_TS h"},g%HI, scale=gv%H\_to\_m)
402 \textcolor{keywordflow}{    endif}
403 \textcolor{keywordflow}{  endif}
404 
405   \textcolor{comment}{!   Calculate the diffusivities, Kd\_lay and Kd\_int, for each layer and interface.  This would}
406   \textcolor{comment}{! be an appropriate place to add a depth-dependent parameterization or another explicit}
407   \textcolor{comment}{! parameterization of Kd.}
408 
409   \textcolor{comment}{!$OMP parallel do default(shared) private(dRho\_int,N2\_lay,Kd\_lay\_2d,Kd\_int\_2d,Kv\_bkgnd,N2\_int,&}
410   \textcolor{comment}{!$OMP                                     N2\_bot,KT\_extra,KS\_extra,TKE\_to\_Kd,maxTKE,dissip,kb)}
411   \textcolor{keywordflow}{do} j=js,je
412 
413     \textcolor{comment}{! Set up variables related to the stratification.}
414     \textcolor{keyword}{call }find\_n2(h, tv, t\_f, s\_f, fluxes, j, g, gv, us, cs, drho\_int, n2\_lay, n2\_int, n2\_bot)
415 
416     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%N2\_3d)) \textcolor{keywordflow}{then}
417       \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie ; dd%N2\_3d(i,j,k) = n2\_int(i,k) ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo}
418 \textcolor{keywordflow}{    endif}
419 
420     \textcolor{comment}{! Add background mixing}
421     \textcolor{keyword}{call }calculate\_bkgnd\_mixing(h, tv, n2\_lay, kd\_lay\_2d, kd\_int\_2d, kv\_bkgnd, j, g, gv, us, cs
      %bkgnd\_mixing\_csp)
422     \textcolor{comment}{! Update Kv and 3-d diffusivity diagnostics.}
423     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_slow)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
424       visc%Kv\_slow(i,j,k) = visc%Kv\_slow(i,j,k) + kv\_bkgnd(i,k)
425 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
426     \textcolor{keywordflow}{if} (cs%id\_Kv\_bkgnd > 0) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
427       dd%Kv\_bkgnd(i,j,k) = kv\_bkgnd(i,k)
428 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
429     \textcolor{keywordflow}{if} (cs%id\_Kd\_bkgnd > 0) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
430       dd%Kd\_bkgnd(i,j,k) = kd\_int\_2d(i,k)
431 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
432 
433     \textcolor{comment}{! Double-diffusion (old method)}
434     \textcolor{keywordflow}{if} (cs%double\_diffusion) \textcolor{keywordflow}{then}
435       \textcolor{keyword}{call }double\_diffusion(tv, h, t\_f, s\_f, j, g, gv, us, cs, kt\_extra, ks\_extra)
436       \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=is,ie
437         \textcolor{keywordflow}{if} (ks\_extra(i,k) > kt\_extra(i,k)) \textcolor{keywordflow}{then} \textcolor{comment}{! salt fingering}
438           kd\_lay\_2d(i,k-1) = kd\_lay\_2d(i,k-1) + 0.5 * kt\_extra(i,k)
439           kd\_lay\_2d(i,k)   = kd\_lay\_2d(i,k)   + 0.5 * kt\_extra(i,k)
440           visc%Kd\_extra\_S(i,j,k) = (ks\_extra(i,k) - kt\_extra(i,k))
441           visc%Kd\_extra\_T(i,j,k) = 0.0
442         \textcolor{keywordflow}{elseif} (kt\_extra(i,k) > 0.0) \textcolor{keywordflow}{then} \textcolor{comment}{! double-diffusive convection}
443           kd\_lay\_2d(i,k-1) = kd\_lay\_2d(i,k-1) + 0.5 * ks\_extra(i,k)
444           kd\_lay\_2d(i,k)   = kd\_lay\_2d(i,k)   + 0.5 * ks\_extra(i,k)
445           visc%Kd\_extra\_T(i,j,k) = (kt\_extra(i,k) - ks\_extra(i,k))
446           visc%Kd\_extra\_S(i,j,k) = 0.0
447         \textcolor{keywordflow}{else} \textcolor{comment}{! There is no double diffusion at this interface.}
448           visc%Kd\_extra\_T(i,j,k) = 0.0
449           visc%Kd\_extra\_S(i,j,k) = 0.0
450 \textcolor{keywordflow}{        endif}
451 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
452       \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%KT\_extra)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
453         dd%KT\_extra(i,j,k) = kt\_extra(i,k)
454 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
455 
456       \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%KS\_extra)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
457         dd%KS\_extra(i,j,k) = ks\_extra(i,k)
458 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
459 \textcolor{keywordflow}{    endif}
460 
461     \textcolor{comment}{! Apply double diffusion via CVMix}
462     \textcolor{comment}{! GMM, we need to pass HBL to compute\_ddiff\_coeffs, but it is not yet available.}
463     \textcolor{keywordflow}{if} (cs%use\_CVMix\_ddiff) \textcolor{keywordflow}{then}
464       \textcolor{keyword}{call }cpu\_clock\_begin(id\_clock\_cvmix\_ddiff)
465       \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%drho\_rat)) \textcolor{keywordflow}{then}
466         \textcolor{keyword}{call }compute\_ddiff\_coeffs(h, tv, g, gv, us, j, visc%Kd\_extra\_T, visc%Kd\_extra\_S, &
467                                   cs%CVMix\_ddiff\_csp, dd%drho\_rat)
468       \textcolor{keywordflow}{else}
469         \textcolor{keyword}{call }compute\_ddiff\_coeffs(h, tv, g, gv, us, j, visc%Kd\_extra\_T, visc%Kd\_extra\_S, cs%CVMix\_ddiff\_csp
      )
470 \textcolor{keywordflow}{      endif}
471       \textcolor{keyword}{call }cpu\_clock\_end(id\_clock\_cvmix\_ddiff)
472 \textcolor{keywordflow}{    endif}
473 
474     \textcolor{comment}{! Calculate conversion ratios from TKE to layer diffusivities.}
475     \textcolor{keywordflow}{if} (tke\_to\_kd\_used) \textcolor{keywordflow}{then}
476       \textcolor{keyword}{call }find\_tke\_to\_kd(h, tv, drho\_int, n2\_lay, j, dt, g, gv, us, cs, tke\_to\_kd, maxtke, kb)
477       \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%maxTKE)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
478         dd%maxTKE(i,j,k) = maxtke(i,k)
479 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
480       \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%TKE\_to\_Kd)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
481         dd%TKE\_to\_Kd(i,j,k) = tke\_to\_kd(i,k)
482 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
483 \textcolor{keywordflow}{    endif}
484 
485     \textcolor{comment}{! Add the input turbulent diffusivity.}
486     \textcolor{keywordflow}{if} (cs%useKappaShear .or. cs%use\_CVMix\_shear) \textcolor{keywordflow}{then}
487       \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) \textcolor{keywordflow}{then}
488         \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=is,ie
489           kd\_int\_2d(i,k) = visc%Kd\_shear(i,j,k) + 0.5 * (kd\_lay\_2d(i,k-1) + kd\_lay\_2d(i,k))
490 \textcolor{keywordflow}{        enddo} ;\textcolor{keywordflow}{ enddo}
491         \textcolor{keywordflow}{do} i=is,ie
492           kd\_int\_2d(i,1) = visc%Kd\_shear(i,j,1) \textcolor{comment}{! This isn't actually used. It could be 0.}
493           kd\_int\_2d(i,nz+1) = 0.0
494 \textcolor{keywordflow}{        enddo}
495 \textcolor{keywordflow}{      endif}
496       \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
497         kd\_lay\_2d(i,k) = kd\_lay\_2d(i,k) + 0.5 * (visc%Kd\_shear(i,j,k) + visc%Kd\_shear(i,j,k+1))
498 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
499     \textcolor{keywordflow}{else}
500       \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) \textcolor{keywordflow}{then}
501         \textcolor{keywordflow}{do} i=is,ie
502           kd\_int\_2d(i,1) = kd\_lay\_2d(i,1) ; kd\_int\_2d(i,nz+1) = 0.0
503 \textcolor{keywordflow}{        enddo}
504         \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=is,ie
505           kd\_int\_2d(i,k) = 0.5 * (kd\_lay\_2d(i,k-1) + kd\_lay\_2d(i,k))
506 \textcolor{keywordflow}{        enddo} ;\textcolor{keywordflow}{ enddo}
507 \textcolor{keywordflow}{      endif}
508 \textcolor{keywordflow}{    endif}
509 
510     \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_int)) \textcolor{keywordflow}{then}
511       \textcolor{comment}{! Add the ML\_Rad diffusivity.}
512       \textcolor{keywordflow}{if} (cs%ML\_radiation) &
513         \textcolor{keyword}{call }add\_mlrad\_diffusivity(h, fluxes, j, g, gv, us, cs, tke\_to\_kd, kd\_lay\_2d, kd\_int\_2d)
514 
515       \textcolor{comment}{! Add the Nikurashin and / or tidal bottom-driven mixing}
516       \textcolor{keywordflow}{if} (cs%use\_tidal\_mixing) &
517         \textcolor{keyword}{call }calculate\_tidal\_mixing(h, n2\_bot, j, tke\_to\_kd, maxtke, g, gv, us, cs%tidal\_mixing\_CSp, &
518                                     n2\_lay, n2\_int, kd\_lay\_2d, kd\_int\_2d, cs%Kd\_max, visc%Kv\_slow)
519 
520       \textcolor{comment}{! This adds the diffusion sustained by the energy extracted from the flow by the bottom drag.}
521       \textcolor{keywordflow}{if} (cs%bottomdraglaw .and. (cs%BBL\_effic>0.0)) \textcolor{keywordflow}{then}
522         \textcolor{keywordflow}{if} (cs%use\_LOTW\_BBL\_diffusivity) \textcolor{keywordflow}{then}
523           \textcolor{keyword}{call }add\_lotw\_bbl\_diffusivity(h, u, v, tv, fluxes, visc, j, n2\_int, g, gv, us, cs,  &
524                                         dd%Kd\_BBL, kd\_lay\_2d, kd\_int\_2d)
525         \textcolor{keywordflow}{else}
526           \textcolor{keyword}{call }add\_drag\_diffusivity(h, u, v,  tv, fluxes, visc, j, tke\_to\_kd, &
527                                     maxtke, kb, g, gv, us, cs, kd\_lay\_2d, kd\_int\_2d, dd%Kd\_BBL)
528 \textcolor{keywordflow}{        endif}
529 \textcolor{keywordflow}{      endif}
530 
531       \textcolor{keywordflow}{if} (cs%limit\_dissipation) \textcolor{keywordflow}{then}
532         \textcolor{comment}{! This calculates the dissipation ONLY from Kd calculated in this routine}
533         \textcolor{comment}{! dissip has units of W/m3 (= kg/m3 * m2/s * 1/s2)}
534         \textcolor{comment}{!   1) a global constant,}
535         \textcolor{comment}{!   2) a dissipation proportional to N (aka Gargett) and}
536         \textcolor{comment}{!   3) dissipation corresponding to a (nearly) constant diffusivity.}
537         \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=is,ie
538           dissip = max( cs%dissip\_min, &   \textcolor{comment}{! Const. floor on dissip.}
539                         cs%dissip\_N0 + cs%dissip\_N1 * sqrt(n2\_int(i,k)), & \textcolor{comment}{! Floor aka Gargett}
540                         cs%dissip\_N2 * n2\_int(i,k)) \textcolor{comment}{! Floor of Kd\_min*rho0/F\_Ri}
541           kd\_int\_2d(i,k) = max(kd\_int\_2d(i,k) , &  \textcolor{comment}{! Apply floor to Kd}
542                               dissip * (cs%FluxRi\_max / (gv%Rho0 * (n2\_int(i,k) + omega2))))
543 \textcolor{keywordflow}{        enddo} ;\textcolor{keywordflow}{ enddo}
544 \textcolor{keywordflow}{      endif}
545 
546       \textcolor{comment}{! Optionally add a uniform diffusivity at the interfaces.}
547       \textcolor{keywordflow}{if} (cs%Kd\_add > 0.0) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
548         kd\_int\_2d(i,k) = kd\_int\_2d(i,k) + cs%Kd\_add
549 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
550 
551       \textcolor{comment}{! Copy the 2-d slices into the 3-d array that is exported.}
552       \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie
553         kd\_int(i,j,k) = kd\_int\_2d(i,k)
554 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
555 
556     \textcolor{keywordflow}{else} \textcolor{comment}{! Kd\_int is not present.}
557 
558       \textcolor{comment}{! Add the ML\_Rad diffusivity.}
559       \textcolor{keywordflow}{if} (cs%ML\_radiation) &
560         \textcolor{keyword}{call }add\_mlrad\_diffusivity(h, fluxes, j, g, gv, us, cs, tke\_to\_kd, kd\_lay\_2d)
561 
562       \textcolor{comment}{! Add the Nikurashin and / or tidal bottom-driven mixing}
563       \textcolor{keywordflow}{if} (cs%use\_tidal\_mixing) &
564         \textcolor{keyword}{call }calculate\_tidal\_mixing(h, n2\_bot, j, tke\_to\_kd, maxtke, g, gv, us, cs%tidal\_mixing\_CSp, &
565                                     n2\_lay, n2\_int, kd\_lay\_2d, kd\_max=cs%Kd\_max, kv=visc%Kv\_slow)
566 
567       \textcolor{comment}{! This adds the diffusion sustained by the energy extracted from the flow by the bottom drag.}
568       \textcolor{keywordflow}{if} (cs%bottomdraglaw .and. (cs%BBL\_effic>0.0)) \textcolor{keywordflow}{then}
569         \textcolor{keywordflow}{if} (cs%use\_LOTW\_BBL\_diffusivity) \textcolor{keywordflow}{then}
570           \textcolor{keyword}{call }add\_lotw\_bbl\_diffusivity(h, u, v, tv, fluxes, visc, j, n2\_int, g, gv, us, cs,  &
571                                         dd%Kd\_BBL, kd\_lay\_2d)
572         \textcolor{keywordflow}{else}
573           \textcolor{keyword}{call }add\_drag\_diffusivity(h, u, v,  tv, fluxes, visc, j, tke\_to\_kd, &
574                                     maxtke, kb, g, gv, us, cs, kd\_lay\_2d, kd\_bbl=dd%Kd\_BBL)
575 \textcolor{keywordflow}{        endif}
576 \textcolor{keywordflow}{      endif}
577 
578 \textcolor{keywordflow}{    endif}
579 
580     \textcolor{keywordflow}{if} (cs%limit\_dissipation) \textcolor{keywordflow}{then}
581       \textcolor{comment}{! This calculates the layer dissipation ONLY from Kd calculated in this routine}
582       \textcolor{comment}{! dissip has units of W/m3 (= kg/m3 * m2/s * 1/s2)}
583       \textcolor{comment}{!   1) a global constant,}
584       \textcolor{comment}{!   2) a dissipation proportional to N (aka Gargett) and}
585       \textcolor{comment}{!   3) dissipation corresponding to a (nearly) constant diffusivity.}
586       \textcolor{keywordflow}{do} k=2,nz-1 ; \textcolor{keywordflow}{do} i=is,ie
587         dissip = max( cs%dissip\_min, &   \textcolor{comment}{! Const. floor on dissip.}
588                       cs%dissip\_N0 + cs%dissip\_N1 * sqrt(n2\_lay(i,k)), & \textcolor{comment}{! Floor aka Gargett}
589                       cs%dissip\_N2 * n2\_lay(i,k)) \textcolor{comment}{! Floor of Kd\_min*rho0/F\_Ri}
590         kd\_lay\_2d(i,k) = max(kd\_lay\_2d(i,k) , &  \textcolor{comment}{! Apply floor to Kd}
591                             dissip * (cs%FluxRi\_max / (gv%Rho0 * (n2\_lay(i,k) + omega2))))
592 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
593 \textcolor{keywordflow}{    endif}
594 
595     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%Kd\_work)) \textcolor{keywordflow}{then}
596       \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
597         dd%Kd\_Work(i,j,k) = gv%Rho0 * kd\_lay\_2d(i,k) * n2\_lay(i,k) * &
598                             gv%H\_to\_Z*h(i,j,k)  \textcolor{comment}{! Watt m-2 s = kg s-3}
599 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
600 \textcolor{keywordflow}{    endif}
601 
602     \textcolor{comment}{! Optionally add a uniform diffusivity to the layers.}
603     \textcolor{keywordflow}{if} ((cs%Kd\_add > 0.0) .and. (\textcolor{keyword}{present}(kd\_lay))) \textcolor{keywordflow}{then}
604       \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
605         kd\_lay\_2d(i,k) = kd\_lay\_2d(i,k) + cs%Kd\_add
606 \textcolor{keywordflow}{      enddo} ;\textcolor{keywordflow}{ enddo}
607 \textcolor{keywordflow}{    endif}
608 
609     \textcolor{comment}{! Copy the 2-d slices into the 3-d array that is exported; this was done above for Kd\_int.}
610     \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_lay)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie
611       kd\_lay(i,j,k) = kd\_lay\_2d(i,k)
612 \textcolor{keywordflow}{    enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif}
613 \textcolor{keywordflow}{  enddo} \textcolor{comment}{! j-loop}
614 
615   \textcolor{keywordflow}{if} (cs%user\_change\_diff) \textcolor{keywordflow}{then}
616     \textcolor{keyword}{call }user\_change\_diff(h, tv, g, gv, us, cs%user\_change\_diff\_CSp, kd\_lay, kd\_int, &
617                           t\_f, s\_f, dd%Kd\_user)
618 \textcolor{keywordflow}{  endif}
619 
620   \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then}
621     \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_lay)) \textcolor{keyword}{call }hchksum(kd\_lay, \textcolor{stringliteral}{"Kd\_lay"}, g%HI, haloshift=0, scale=us%Z2\_T\_to\_m2\_s)
622 
623     \textcolor{keywordflow}{if} (cs%useKappaShear) \textcolor{keyword}{call }hchksum(visc%Kd\_shear, \textcolor{stringliteral}{"Turbulent Kd"}, g%HI, haloshift=0, scale=us
      %Z2\_T\_to\_m2\_s)
624 
625     \textcolor{keywordflow}{if} (cs%use\_CVMix\_ddiff) \textcolor{keywordflow}{then}
626       \textcolor{keyword}{call }hchksum(visc%Kd\_extra\_T, \textcolor{stringliteral}{"MOM\_set\_diffusivity: Kd\_extra\_T"}, g%HI, haloshift=0, scale=us
      %Z2\_T\_to\_m2\_s)
627       \textcolor{keyword}{call }hchksum(visc%Kd\_extra\_S, \textcolor{stringliteral}{"MOM\_set\_diffusivity: Kd\_extra\_S"}, g%HI, haloshift=0, scale=us
      %Z2\_T\_to\_m2\_s)
628 \textcolor{keywordflow}{    endif}
629 
630     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%kv\_bbl\_u) .and. \textcolor{keyword}{associated}(visc%kv\_bbl\_v)) \textcolor{keywordflow}{then}
631       \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"BBL Kv\_bbl\_[uv]"}, visc%kv\_bbl\_u, visc%kv\_bbl\_v, g%HI, &
632                     haloshift=0, symmetric=.true., scale=us%Z2\_T\_to\_m2\_s, &
633                     scalar\_pair=.true.)
634 \textcolor{keywordflow}{    endif}
635 
636     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%bbl\_thick\_u) .and. \textcolor{keyword}{associated}(visc%bbl\_thick\_v)) \textcolor{keywordflow}{then}
637       \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"BBL bbl\_thick\_[uv]"}, visc%bbl\_thick\_u, visc%bbl\_thick\_v, &
638                     g%HI, haloshift=0, symmetric=.true., scale=us%Z\_to\_m, &
639                     scalar\_pair=.true.)
640 \textcolor{keywordflow}{    endif}
641 
642     \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Ray\_u) .and. \textcolor{keyword}{associated}(visc%Ray\_v)) \textcolor{keywordflow}{then}
643       \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"Ray\_[uv]"}, visc%Ray\_u, visc%Ray\_v, g%HI, 0, symmetric=.true., scale=us%Z\_to\_m*us
      %s\_to\_T)
644 \textcolor{keywordflow}{    endif}
645 
646 \textcolor{keywordflow}{  endif}
647 
648   \textcolor{comment}{! post diagnostics}
649   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kd\_lay) .and. (cs%id\_Kd\_layer > 0)) \textcolor{keyword}{call }post\_data(cs%id\_Kd\_layer, kd\_lay, cs%diag)
650 
651   \textcolor{comment}{! background mixing}
652   \textcolor{keywordflow}{if} (cs%id\_Kd\_bkgnd > 0) \textcolor{keyword}{call }post\_data(cs%id\_Kd\_bkgnd, dd%Kd\_bkgnd, cs%diag)
653   \textcolor{keywordflow}{if} (cs%id\_Kv\_bkgnd > 0) \textcolor{keyword}{call }post\_data(cs%id\_Kv\_bkgnd, dd%Kv\_bkgnd, cs%diag)
654 
655   \textcolor{comment}{! tidal mixing}
656   \textcolor{keyword}{call }post\_tidal\_diagnostics(g, gv, h, cs%tidal\_mixing\_CSp)
657   \textcolor{keywordflow}{if} (cs%id\_N2 > 0)         \textcolor{keyword}{call }post\_data(cs%id\_N2,        dd%N2\_3d,     cs%diag)
658   \textcolor{keywordflow}{if} (cs%id\_Kd\_Work > 0)    \textcolor{keyword}{call }post\_data(cs%id\_Kd\_Work,   dd%Kd\_Work,   cs%diag)
659   \textcolor{keywordflow}{if} (cs%id\_maxTKE > 0)     \textcolor{keyword}{call }post\_data(cs%id\_maxTKE,    dd%maxTKE,    cs%diag)
660   \textcolor{keywordflow}{if} (cs%id\_TKE\_to\_Kd > 0)  \textcolor{keyword}{call }post\_data(cs%id\_TKE\_to\_Kd, dd%TKE\_to\_Kd, cs%diag)
661 
662   \textcolor{keywordflow}{if} (cs%id\_Kd\_user > 0)    \textcolor{keyword}{call }post\_data(cs%id\_Kd\_user,   dd%Kd\_user,   cs%diag)
663 
664   \textcolor{comment}{! double diffusive mixing}
665   \textcolor{keywordflow}{if} (cs%double\_diffusion) \textcolor{keywordflow}{then}
666     \textcolor{keywordflow}{if} (cs%id\_KT\_extra > 0) \textcolor{keyword}{call }post\_data(cs%id\_KT\_extra, dd%KT\_extra, cs%diag)
667     \textcolor{keywordflow}{if} (cs%id\_KS\_extra > 0) \textcolor{keyword}{call }post\_data(cs%id\_KS\_extra, dd%KS\_extra, cs%diag)
668   \textcolor{keywordflow}{else}
669     \textcolor{keywordflow}{if} (cs%id\_KT\_extra > 0) \textcolor{keyword}{call }post\_data(cs%id\_KT\_extra, visc%Kd\_extra\_T, cs%diag)
670     \textcolor{keywordflow}{if} (cs%id\_KS\_extra > 0) \textcolor{keyword}{call }post\_data(cs%id\_KS\_extra, visc%Kd\_extra\_S, cs%diag)
671     \textcolor{keywordflow}{if} (cs%id\_R\_rho > 0) \textcolor{keyword}{call }post\_data(cs%id\_R\_rho, dd%drho\_rat, cs%diag)
672 \textcolor{keywordflow}{  endif}
673   \textcolor{keywordflow}{if} (cs%id\_Kd\_BBL > 0)   \textcolor{keyword}{call }post\_data(cs%id\_Kd\_BBL, dd%Kd\_BBL, cs%diag)
674 
675   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%N2\_3d)) \textcolor{keyword}{deallocate}(dd%N2\_3d)
676   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%Kd\_work)) \textcolor{keyword}{deallocate}(dd%Kd\_work)
677   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%Kd\_user)) \textcolor{keyword}{deallocate}(dd%Kd\_user)
678   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%maxTKE)) \textcolor{keyword}{deallocate}(dd%maxTKE)
679   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%TKE\_to\_Kd)) \textcolor{keyword}{deallocate}(dd%TKE\_to\_Kd)
680   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%KT\_extra)) \textcolor{keyword}{deallocate}(dd%KT\_extra)
681   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%KS\_extra)) \textcolor{keyword}{deallocate}(dd%KS\_extra)
682   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%drho\_rat)) \textcolor{keyword}{deallocate}(dd%drho\_rat)
683   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(dd%Kd\_BBL)) \textcolor{keyword}{deallocate}(dd%Kd\_BBL)
684 
685   \textcolor{keywordflow}{if} (showcalltree) \textcolor{keyword}{call }calltree\_leave(\textcolor{stringliteral}{"set\_diffusivity()"})
686 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_ace82f133d3cee42aa36ec10bcce79e75}\label{namespacemom__set__diffusivity_ace82f133d3cee42aa36ec10bcce79e75}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!set\+\_\+diffusivity\+\_\+end@{set\+\_\+diffusivity\+\_\+end}}
\index{set\+\_\+diffusivity\+\_\+end@{set\+\_\+diffusivity\+\_\+end}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{set\+\_\+diffusivity\+\_\+end()}{set\_diffusivity\_end()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+diffusivity\+::set\+\_\+diffusivity\+\_\+end (\begin{DoxyParamCaption}\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS }\end{DoxyParamCaption})}



Clear pointers and dealocate memory. 


\begin{DoxyParams}{Parameters}
{\em cs} & Control structure for this module \\
\hline
\end{DoxyParams}


Definition at line 2322 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
2322   \textcolor{keywordtype}{type}(set\_diffusivity\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< Control structure for this module}
2323 
2324   \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keywordflow}{return}
2325 
2326   \textcolor{keyword}{call }bkgnd\_mixing\_end(cs%bkgnd\_mixing\_csp)
2327 
2328   \textcolor{keywordflow}{if} (cs%use\_tidal\_mixing) \textcolor{keyword}{call }tidal\_mixing\_end(cs%tidal\_mixing\_CSp)
2329 
2330   \textcolor{keywordflow}{if} (cs%user\_change\_diff) \textcolor{keyword}{call }user\_change\_diff\_end(cs%user\_change\_diff\_CSp)
2331 
2332   \textcolor{keywordflow}{if} (cs%use\_CVMix\_shear)  \textcolor{keyword}{call }cvmix\_shear\_end(cs%CVMix\_shear\_csp)
2333 
2334   \textcolor{keywordflow}{if} (cs%use\_CVMix\_ddiff)  \textcolor{keyword}{call }cvmix\_ddiff\_end(cs%CVMix\_ddiff\_csp)
2335 
2336   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{deallocate}(cs)
2337 
\end{DoxyCode}
\mbox{\Hypertarget{namespacemom__set__diffusivity_a2296f20ef1f3a4c2390897a7376f88e7}\label{namespacemom__set__diffusivity_a2296f20ef1f3a4c2390897a7376f88e7}} 
\index{mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}!set\+\_\+diffusivity\+\_\+init@{set\+\_\+diffusivity\+\_\+init}}
\index{set\+\_\+diffusivity\+\_\+init@{set\+\_\+diffusivity\+\_\+init}!mom\+\_\+set\+\_\+diffusivity@{mom\+\_\+set\+\_\+diffusivity}}
\subsubsection{\texorpdfstring{set\+\_\+diffusivity\+\_\+init()}{set\_diffusivity\_init()}}
{\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+diffusivity\+::set\+\_\+diffusivity\+\_\+init (\begin{DoxyParamCaption}\item[{type(time\+\_\+type), intent(in)}]{Time,  }\item[{type(ocean\+\_\+grid\+\_\+type), intent(inout)}]{G,  }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV,  }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US,  }\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file,  }\item[{type(diag\+\_\+ctrl), intent(inout), target}]{diag,  }\item[{type(\hyperlink{structmom__set__diffusivity_1_1set__diffusivity__cs}{set\+\_\+diffusivity\+\_\+cs}), pointer}]{CS,  }\item[{type(int\+\_\+tide\+\_\+cs), pointer}]{int\+\_\+tide\+\_\+\+C\+Sp,  }\item[{integer, intent(out), optional}]{halo\+\_\+\+TS }\end{DoxyParamCaption})}


\begin{DoxyParams}[1]{Parameters}
\mbox{\tt in}  & {\em time} & The current model time\\
\hline
\mbox{\tt in,out}  & {\em g} & The ocean\textquotesingle{}s grid structure.\\
\hline
\mbox{\tt in}  & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\
\hline
\mbox{\tt in}  & {\em us} & A dimensional unit scaling type\\
\hline
\mbox{\tt in}  & {\em param\+\_\+file} & A structure to parse for run-\/time parameters.\\
\hline
\mbox{\tt in,out}  & {\em diag} & A structure used to regulate diagnostic output.\\
\hline
 & {\em cs} & pointer set to point to the module control structure.\\
\hline
 & {\em int\+\_\+tide\+\_\+csp} & A pointer to the internal tides control structure\\
\hline
\mbox{\tt out}  & {\em halo\+\_\+ts} & The halo size of tracer points that must be valid for the calculations in set\+\_\+diffusivity. \\
\hline
\end{DoxyParams}


Definition at line 1998 of file M\+O\+M\+\_\+set\+\_\+diffusivity.\+F90.


\begin{DoxyCode}
1998   \textcolor{keywordtype}{type}(time\_type),          \textcolor{keywordtype}{intent(in)}    :: time\textcolor{comment}{ !< The current model time}
1999   \textcolor{keywordtype}{type}(ocean\_grid\_type),    \textcolor{keywordtype}{intent(inout)} :: g\textcolor{comment}{    !< The ocean's grid structure.}
2000   \textcolor{keywordtype}{type}(verticalgrid\_type),  \textcolor{keywordtype}{intent(in)}    :: gv\textcolor{comment}{   !< The ocean's vertical grid structure.}
2001   \textcolor{keywordtype}{type}(unit\_scale\_type),    \textcolor{keywordtype}{intent(in)}    :: us\textcolor{comment}{   !< A dimensional unit scaling type}
2002   \textcolor{keywordtype}{type}(param\_file\_type),    \textcolor{keywordtype}{intent(in)}    :: param\_file\textcolor{comment}{ !< A structure to parse for run-time}
2003 \textcolor{comment}{                                                  !! parameters.}
2004   \textcolor{keywordtype}{type}(diag\_ctrl), \textcolor{keywordtype}{target},  \textcolor{keywordtype}{intent(inout)} :: diag\textcolor{comment}{ !< A structure used to regulate diagnostic output.}
2005   \textcolor{keywordtype}{type}(set\_diffusivity\_cs), \textcolor{keywordtype}{pointer}       :: cs\textcolor{comment}{   !< pointer set to point to the module control}
2006 \textcolor{comment}{                                                  !! structure.}
2007   \textcolor{keywordtype}{type}(int\_tide\_cs),        \textcolor{keywordtype}{pointer}       :: int\_tide\_csp\textcolor{comment}{ !< A pointer to the internal tides control}
2008 \textcolor{comment}{                                                  !! structure}
2009   \textcolor{keywordtype}{integer},        \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)}   :: halo\_ts\textcolor{comment}{ !< The halo size of tracer points that must be}
2010 \textcolor{comment}{                                                  !! valid for the calculations in set\_diffusivity.}
2011 
2012   \textcolor{comment}{! Local variables}
2013   \textcolor{keywordtype}{real} :: decay\_length
2014   \textcolor{keywordtype}{logical} :: ml\_use\_omega
2015   \textcolor{keywordtype}{logical} :: default\_2018\_answers
2016   \textcolor{comment}{! This include declares and sets the variable "version".}
2017 \textcolor{preprocessor}{# include "version\_variable.h"}
2018 \textcolor{preprocessor}{}  \textcolor{keywordtype}{character(len=40)}  :: mdl = \textcolor{stringliteral}{"MOM\_set\_diffusivity"}  \textcolor{comment}{! This module's name.}
2019   \textcolor{keywordtype}{real} :: omega\_frac\_dflt
2020   \textcolor{keywordtype}{logical} :: bryan\_lewis\_diffusivity \textcolor{comment}{! If true, the background diapycnal diffusivity uses}
2021                                      \textcolor{comment}{! the Bryan-Lewis (1979) style tanh profile.}
2022   \textcolor{keywordtype}{integer} :: i, j, is, ie, js, je
2023   \textcolor{keywordtype}{integer} :: isd, ied, jsd, jed
2024 
2025   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keywordflow}{then}
2026     \textcolor{keyword}{call }mom\_error(warning, \textcolor{stringliteral}{"diabatic\_entrain\_init called with an associated "}// &
2027                             \textcolor{stringliteral}{"control structure."})
2028     \textcolor{keywordflow}{return}
2029 \textcolor{keywordflow}{  endif}
2030   \textcolor{keyword}{allocate}(cs)
2031 
2032   is  = g%isc ; ie  = g%iec ; js  = g%jsc ; je  = g%jec
2033   isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
2034 
2035   cs%diag => diag
2036   \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(int\_tide\_csp))  cs%int\_tide\_CSp  => int\_tide\_csp
2037 
2038   \textcolor{comment}{! These default values always need to be set.}
2039   cs%BBL\_mixing\_as\_max = .true.
2040   cs%cdrag = 0.003 ; cs%BBL\_effic = 0.0
2041   cs%bulkmixedlayer = (gv%nkml > 0)
2042 
2043   \textcolor{comment}{! Read all relevant parameters and write them to the model log.}
2044   \textcolor{keyword}{call }log\_version(param\_file, mdl, version, \textcolor{stringliteral}{""})
2045 
2046   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"INPUTDIR"}, cs%inputdir, default=\textcolor{stringliteral}{"."})
2047   cs%inputdir = slasher(cs%inputdir)
2048   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"FLUX\_RI\_MAX"}, cs%FluxRi\_max, &
2049                  \textcolor{stringliteral}{"The flux Richardson number where the stratification is "}//&
2050                  \textcolor{stringliteral}{"large enough that N2 > omega2.  The full expression for "}//&
2051                  \textcolor{stringliteral}{"the Flux Richardson number is usually "}//&
2052                  \textcolor{stringliteral}{"FLUX\_RI\_MAX*N2/(N2+OMEGA2)."}, units=\textcolor{stringliteral}{"nondim"}, default=0.2)
2053   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"OMEGA"}, cs%omega, &
2054                  \textcolor{stringliteral}{"The rotation rate of the earth."}, units=\textcolor{stringliteral}{"s-1"}, default=7.2921e-5, scale=us%T\_to\_s)
2055 
2056   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DEFAULT\_2018\_ANSWERS"}, default\_2018\_answers, &
2057                  \textcolor{stringliteral}{"This sets the default value for the various \_2018\_ANSWERS parameters."}, &
2058                  default=.false.)
2059   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"SET\_DIFF\_2018\_ANSWERS"}, cs%answers\_2018, &
2060                  \textcolor{stringliteral}{"If true, use the order of arithmetic and expressions that recover the "}//&
2061                  \textcolor{stringliteral}{"answers from the end of 2018.  Otherwise, use updated and more robust "}//&
2062                  \textcolor{stringliteral}{"forms of the same expressions."}, default=default\_2018\_answers)
2063 
2064   \textcolor{comment}{! CS%use\_tidal\_mixing is set to True if an internal tidal dissipation scheme is to be used.}
2065   cs%use\_tidal\_mixing = tidal\_mixing\_init(time, g, gv, us, param\_file, diag, &
2066                                           cs%tidal\_mixing\_CSp)
2067 
2068   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_RADIATION"}, cs%ML\_radiation, &
2069                  \textcolor{stringliteral}{"If true, allow a fraction of TKE available from wind "}//&
2070                  \textcolor{stringliteral}{"work to penetrate below the base of the mixed layer "}//&
2071                  \textcolor{stringliteral}{"with a vertical decay scale determined by the minimum "}//&
2072                  \textcolor{stringliteral}{"of: (1) The depth of the mixed layer, (2) an Ekman "}//&
2073                  \textcolor{stringliteral}{"length scale."}, default=.false.)
2074   \textcolor{keywordflow}{if} (cs%ML\_radiation) \textcolor{keywordflow}{then}
2075     \textcolor{comment}{! This give a minimum decay scale that is typically much less than Angstrom.}
2076     cs%ustar\_min = 2e-4 * cs%omega * (gv%Angstrom\_Z + gv%H\_subroundoff*gv%H\_to\_Z)
2077 
2078     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_RAD\_EFOLD\_COEFF"}, cs%ML\_rad\_efold\_coeff, &
2079                  \textcolor{stringliteral}{"A coefficient that is used to scale the penetration "}//&
2080                  \textcolor{stringliteral}{"depth for turbulence below the base of the mixed layer. "}//&
2081                  \textcolor{stringliteral}{"This is only used if ML\_RADIATION is true."}, units=\textcolor{stringliteral}{"nondim"}, default=0.2)
2082     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_RAD\_BUG"}, cs%ML\_rad\_bug, &
2083                  \textcolor{stringliteral}{"If true use code with a bug that reduces the energy available "}//&
2084                  \textcolor{stringliteral}{"in the transition layer by a factor of the inverse of the energy "}//&
2085                  \textcolor{stringliteral}{"deposition lenthscale (in m)."}, default=.false.)
2086     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_RAD\_KD\_MAX"}, cs%ML\_rad\_kd\_max, &
2087                  \textcolor{stringliteral}{"The maximum diapycnal diffusivity due to turbulence "}//&
2088                  \textcolor{stringliteral}{"radiated from the base of the mixed layer. "}//&
2089                  \textcolor{stringliteral}{"This is only used if ML\_RADIATION is true."}, &
2090                  units=\textcolor{stringliteral}{"m2 s-1"}, default=1.0e-3, scale=us%m2\_s\_to\_Z2\_T)
2091     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_RAD\_COEFF"}, cs%ML\_rad\_coeff, &
2092                  \textcolor{stringliteral}{"The coefficient which scales MSTAR*USTAR^3 to obtain "}//&
2093                  \textcolor{stringliteral}{"the energy available for mixing below the base of the "}//&
2094                  \textcolor{stringliteral}{"mixed layer. This is only used if ML\_RADIATION is true."}, &
2095                  units=\textcolor{stringliteral}{"nondim"}, default=0.2)
2096     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_RAD\_APPLY\_TKE\_DECAY"}, cs%ML\_rad\_TKE\_decay, &
2097                  \textcolor{stringliteral}{"If true, apply the same exponential decay to ML\_rad as "}//&
2098                  \textcolor{stringliteral}{"is applied to the other surface sources of TKE in the "}//&
2099                  \textcolor{stringliteral}{"mixed layer code. This is only used if ML\_RADIATION is true."}, default=.true.)
2100     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MSTAR"}, cs%mstar, &
2101                  \textcolor{stringliteral}{"The ratio of the friction velocity cubed to the TKE "}//&
2102                  \textcolor{stringliteral}{"input to the mixed layer."}, units=\textcolor{stringliteral}{"nondim"}, default=1.2)
2103     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TKE\_DECAY"}, cs%TKE\_decay, &
2104                  \textcolor{stringliteral}{"The ratio of the natural Ekman depth to the TKE decay scale."}, &
2105                  units=\textcolor{stringliteral}{"nondim"}, default=2.5)
2106     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_USE\_OMEGA"}, ml\_use\_omega, &
2107                  \textcolor{stringliteral}{"If true, use the absolute rotation rate instead of the "}//&
2108                  \textcolor{stringliteral}{"vertical component of rotation when setting the decay "}//&
2109                  \textcolor{stringliteral}{"scale for turbulence."}, default=.false., do\_not\_log=.true.)
2110     omega\_frac\_dflt = 0.0
2111     \textcolor{keywordflow}{if} (ml\_use\_omega) \textcolor{keywordflow}{then}
2112       \textcolor{keyword}{call }mom\_error(warning, \textcolor{stringliteral}{"ML\_USE\_OMEGA is depricated; use ML\_OMEGA\_FRAC=1.0 instead."})
2113       omega\_frac\_dflt = 1.0
2114 \textcolor{keywordflow}{    endif}
2115     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_OMEGA\_FRAC"}, cs%ML\_omega\_frac, &
2116                    \textcolor{stringliteral}{"When setting the decay scale for turbulence, use this "}//&
2117                    \textcolor{stringliteral}{"fraction of the absolute rotation rate blended with the "}//&
2118                    \textcolor{stringliteral}{"local value of f, as sqrt((1-of)*f^2 + of*4*omega^2)."}, &
2119                    units=\textcolor{stringliteral}{"nondim"}, default=omega\_frac\_dflt)
2120 \textcolor{keywordflow}{  endif}
2121 
2122   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BOTTOMDRAGLAW"}, cs%bottomdraglaw, &
2123                  \textcolor{stringliteral}{"If true, the bottom stress is calculated with a drag "}//&
2124                  \textcolor{stringliteral}{"law of the form c\_drag*|u|*u. The velocity magnitude "}//&
2125                  \textcolor{stringliteral}{"may be an assumed value or it may be based on the actual "}//&
2126                  \textcolor{stringliteral}{"velocity in the bottommost HBBL, depending on LINEAR\_DRAG."}, default=.true.)
2127   \textcolor{keywordflow}{if}  (cs%bottomdraglaw) \textcolor{keywordflow}{then}
2128     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"CDRAG"}, cs%cdrag, &
2129                  \textcolor{stringliteral}{"The drag coefficient relating the magnitude of the "}//&
2130                  \textcolor{stringliteral}{"velocity field to the bottom stress. CDRAG is only used "}//&
2131                  \textcolor{stringliteral}{"if BOTTOMDRAGLAW is true."}, units=\textcolor{stringliteral}{"nondim"}, default=0.003)
2132     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BBL\_EFFIC"}, cs%BBL\_effic, &
2133                  \textcolor{stringliteral}{"The efficiency with which the energy extracted by "}//&
2134                  \textcolor{stringliteral}{"bottom drag drives BBL diffusion.  This is only "}//&
2135                  \textcolor{stringliteral}{"used if BOTTOMDRAGLAW is true."}, units=\textcolor{stringliteral}{"nondim"}, default=0.20)
2136     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BBL\_MIXING\_MAX\_DECAY"}, decay\_length, &
2137                  \textcolor{stringliteral}{"The maximum decay scale for the BBL diffusion, or 0 to allow the mixing "}//&
2138                  \textcolor{stringliteral}{"to penetrate as far as stratification and rotation permit.  The default "}//&
2139                  \textcolor{stringliteral}{"for now is 200 m. This is only used if BOTTOMDRAGLAW is true."}, &
2140                  units=\textcolor{stringliteral}{"m"}, default=200.0, scale=us%m\_to\_Z)
2141 
2142     cs%IMax\_decay = 0.0
2143     \textcolor{keywordflow}{if} (decay\_length > 0.0) cs%IMax\_decay = 1.0/decay\_length
2144     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BBL\_MIXING\_AS\_MAX"}, cs%BBL\_mixing\_as\_max, &
2145                  \textcolor{stringliteral}{"If true, take the maximum of the diffusivity from the "}//&
2146                  \textcolor{stringliteral}{"BBL mixing and the other diffusivities. Otherwise, "}//&
2147                  \textcolor{stringliteral}{"diffusivity from the BBL\_mixing is simply added."}, &
2148                  default=.true.)
2149     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_LOTW\_BBL\_DIFFUSIVITY"}, cs%use\_LOTW\_BBL\_diffusivity, &
2150                  \textcolor{stringliteral}{"If true, uses a simple, imprecise but non-coordinate dependent, model "}//&
2151                  \textcolor{stringliteral}{"of BBL mixing diffusivity based on Law of the Wall. Otherwise, uses "}//&
2152                  \textcolor{stringliteral}{"the original BBL scheme."}, default=.false.)
2153     \textcolor{keywordflow}{if} (cs%use\_LOTW\_BBL\_diffusivity) \textcolor{keywordflow}{then}
2154       \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"LOTW\_BBL\_USE\_OMEGA"}, cs%LOTW\_BBL\_use\_omega, &
2155                  \textcolor{stringliteral}{"If true, use the maximum of Omega and N for the TKE to diffusion "}//&
2156                  \textcolor{stringliteral}{"calculation. Otherwise, N is N."}, default=.true.)
2157 \textcolor{keywordflow}{    endif}
2158   \textcolor{keywordflow}{else}
2159     cs%use\_LOTW\_BBL\_diffusivity = .false. \textcolor{comment}{! This parameterization depends on a u* from viscous BBL}
2160 \textcolor{keywordflow}{  endif}
2161   cs%id\_Kd\_BBL = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Kd\_BBL'}, diag%axesTi, time, &
2162                  \textcolor{stringliteral}{'Bottom Boundary Layer Diffusivity'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s)
2163   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"SIMPLE\_TKE\_TO\_KD"}, cs%simple\_TKE\_to\_Kd, &
2164                  \textcolor{stringliteral}{"If true, uses a simple estimate of Kd/TKE that will "}//&
2165                  \textcolor{stringliteral}{"work for arbitrary vertical coordinates. If false, "}//&
2166                  \textcolor{stringliteral}{"calculates Kd/TKE and bounds based on exact energetics "}//&
2167                  \textcolor{stringliteral}{"for an isopycnal layer-formulation."}, default=.false.)
2168 
2169   \textcolor{comment}{! set params related to the background mixing}
2170   \textcolor{keyword}{call }bkgnd\_mixing\_init(time, g, gv, us, param\_file, cs%diag, cs%bkgnd\_mixing\_csp)
2171 
2172   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KV"}, cs%Kv, &
2173                  \textcolor{stringliteral}{"The background kinematic viscosity in the interior. "}//&
2174                  \textcolor{stringliteral}{"The molecular value, ~1e-6 m2 s-1, may be used."}, &
2175                  units=\textcolor{stringliteral}{"m2 s-1"}, scale=us%m2\_s\_to\_Z2\_T, fail\_if\_missing=.true.)
2176 
2177   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD"}, cs%Kd, &
2178                  \textcolor{stringliteral}{"The background diapycnal diffusivity of density in the "}//&
2179                  \textcolor{stringliteral}{"interior. Zero or the molecular value, ~1e-7 m2 s-1, "}//&
2180                  \textcolor{stringliteral}{"may be used."}, default=0.0, units=\textcolor{stringliteral}{"m2 s-1"}, scale=us%m2\_s\_to\_Z2\_T)
2181   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD\_MIN"}, cs%Kd\_min, &
2182                  \textcolor{stringliteral}{"The minimum diapycnal diffusivity."}, &
2183                  units=\textcolor{stringliteral}{"m2 s-1"}, default=0.01*cs%Kd*us%Z2\_T\_to\_m2\_s, scale=us%m2\_s\_to\_Z2\_T)
2184   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD\_MAX"}, cs%Kd\_max, &
2185                  \textcolor{stringliteral}{"The maximum permitted increment for the diapycnal "}//&
2186                  \textcolor{stringliteral}{"diffusivity from TKE-based parameterizations, or a negative "}//&
2187                  \textcolor{stringliteral}{"value for no limit."}, units=\textcolor{stringliteral}{"m2 s-1"}, default=-1.0, scale=us%m2\_s\_to\_Z2\_T)
2188   \textcolor{keywordflow}{if} (cs%simple\_TKE\_to\_Kd .and. cs%Kd\_max<=0.) \textcolor{keyword}{call }mom\_error(fatal, &
2189          \textcolor{stringliteral}{"set\_diffusivity\_init: To use SIMPLE\_TKE\_TO\_KD, KD\_MAX must be set to >0."})
2190   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD\_ADD"}, cs%Kd\_add, &
2191                  \textcolor{stringliteral}{"A uniform diapycnal diffusivity that is added "}//&
2192                  \textcolor{stringliteral}{"everywhere without any filtering or scaling."}, &
2193                  units=\textcolor{stringliteral}{"m2 s-1"}, default=0.0, scale=us%m2\_s\_to\_Z2\_T)
2194   \textcolor{keywordflow}{if} (cs%use\_LOTW\_BBL\_diffusivity .and. cs%Kd\_max<=0.) \textcolor{keyword}{call }mom\_error(fatal, &
2195                  \textcolor{stringliteral}{"set\_diffusivity\_init: KD\_MAX must be set (positive) when "}// &
2196                  \textcolor{stringliteral}{"USE\_LOTW\_BBL\_DIFFUSIVITY=True."})
2197   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD\_SMOOTH"}, cs%Kd\_smooth, &
2198                  \textcolor{stringliteral}{"A diapycnal diffusivity that is used to interpolate "}//&
2199                  \textcolor{stringliteral}{"more sensible values of T & S into thin layers."}, &
2200                  units=\textcolor{stringliteral}{"m2 s-1"}, default=1.0e-6, scale=us%m2\_s\_to\_Z2\_T)
2201 
2202   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DEBUG"}, cs%debug, &
2203                  \textcolor{stringliteral}{"If true, write out verbose debugging data."}, &
2204                  default=.false., debuggingparam=.true.)
2205 
2206   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USER\_CHANGE\_DIFFUSIVITY"}, cs%user\_change\_diff, &
2207                  \textcolor{stringliteral}{"If true, call user-defined code to change the diffusivity."}, default=.false.)
2208 
2209   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DISSIPATION\_MIN"}, cs%dissip\_min, &
2210                  \textcolor{stringliteral}{"The minimum dissipation by which to determine a lower "}//&
2211                  \textcolor{stringliteral}{"bound of Kd (a floor)."}, &
2212                  units=\textcolor{stringliteral}{"W m-3"}, default=0.0, scale=us%W\_m2\_to\_RZ3\_T3*us%Z\_to\_m)
2213   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DISSIPATION\_N0"}, cs%dissip\_N0, &
2214                  \textcolor{stringliteral}{"The intercept when N=0 of the N-dependent expression "}//&
2215                  \textcolor{stringliteral}{"used to set a minimum dissipation by which to determine "}//&
2216                  \textcolor{stringliteral}{"a lower bound of Kd (a floor): A in eps\_min = A + B*N."}, &
2217                  units=\textcolor{stringliteral}{"W m-3"}, default=0.0, scale=us%W\_m2\_to\_RZ3\_T3*us%Z\_to\_m)
2218   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DISSIPATION\_N1"}, cs%dissip\_N1, &
2219                  \textcolor{stringliteral}{"The coefficient multiplying N, following Gargett, used to "}//&
2220                  \textcolor{stringliteral}{"set a minimum dissipation by which to determine a lower "}//&
2221                  \textcolor{stringliteral}{"bound of Kd (a floor): B in eps\_min = A + B*N"}, &
2222                  units=\textcolor{stringliteral}{"J m-3"}, default=0.0, scale=us%W\_m2\_to\_RZ3\_T3*us%Z\_to\_m*us%s\_to\_T)
2223   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DISSIPATION\_KD\_MIN"}, cs%dissip\_Kd\_min, &
2224                  \textcolor{stringliteral}{"The minimum vertical diffusivity applied as a floor."}, &
2225                  units=\textcolor{stringliteral}{"m2 s-1"}, default=0.0, scale=us%m2\_s\_to\_Z2\_T)
2226 
2227   cs%limit\_dissipation = (cs%dissip\_min>0.) .or. (cs%dissip\_N1>0.) .or. &
2228                          (cs%dissip\_N0>0.) .or. (cs%dissip\_Kd\_min>0.)
2229   cs%dissip\_N2 = 0.0
2230   \textcolor{keywordflow}{if} (cs%FluxRi\_max > 0.0) &
2231     cs%dissip\_N2 = cs%dissip\_Kd\_min * gv%Rho0 / cs%FluxRi\_max
2232 
2233   cs%id\_Kd\_bkgnd = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Kd\_bkgnd'}, diag%axesTi, time, &
2234       \textcolor{stringliteral}{'Background diffusivity added by MOM\_bkgnd\_mixing module'}, \textcolor{stringliteral}{'m2/s'}, conversion=us%Z2\_T\_to\_m2\_s)
2235   cs%id\_Kv\_bkgnd = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Kv\_bkgnd'}, diag%axesTi, time, &
2236       \textcolor{stringliteral}{'Background viscosity added by MOM\_bkgnd\_mixing module'}, \textcolor{stringliteral}{'m2/s'}, conversion=us%Z2\_T\_to\_m2\_s)
2237 
2238   cs%id\_Kd\_layer = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Kd\_layer'}, diag%axesTL, time, &
2239       \textcolor{stringliteral}{'Diapycnal diffusivity of layers (as set)'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s)
2240 
2241   \textcolor{keywordflow}{if} (cs%use\_tidal\_mixing) \textcolor{keywordflow}{then}
2242     cs%id\_Kd\_Work = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Kd\_Work'}, diag%axesTL, time, &
2243          \textcolor{stringliteral}{'Work done by Diapycnal Mixing'}, \textcolor{stringliteral}{'W m-2'}, conversion=us%RZ3\_T3\_to\_W\_m2)
2244     cs%id\_maxTKE = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'maxTKE'}, diag%axesTL, time, &
2245            \textcolor{stringliteral}{'Maximum layer TKE'}, \textcolor{stringliteral}{'m3 s-3'}, conversion=(us%Z\_to\_m**3*us%s\_to\_T**3))
2246     cs%id\_TKE\_to\_Kd = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'TKE\_to\_Kd'}, diag%axesTL, time, &
2247            \textcolor{stringliteral}{'Convert TKE to Kd'}, \textcolor{stringliteral}{'s2 m'}, conversion=us%Z2\_T\_to\_m2\_s*(us%m\_to\_Z**3*us%T\_to\_s**3))
2248     cs%id\_N2 = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'N2'}, diag%axesTi, time, &
2249          \textcolor{stringliteral}{'Buoyancy frequency squared'}, \textcolor{stringliteral}{'s-2'}, conversion=us%s\_to\_T**2, cmor\_field\_name=\textcolor{stringliteral}{'obvfsq'}, &
2250           cmor\_long\_name=\textcolor{stringliteral}{'Square of seawater buoyancy frequency'}, &
2251           cmor\_standard\_name=\textcolor{stringliteral}{'square\_of\_brunt\_vaisala\_frequency\_in\_sea\_water'})
2252 \textcolor{keywordflow}{  endif}
2253 
2254   \textcolor{keywordflow}{if} (cs%user\_change\_diff) &
2255     cs%id\_Kd\_user = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Kd\_user'}, diag%axesTi, time, &
2256          \textcolor{stringliteral}{'User-specified Extra Diffusivity'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s)
2257 
2258   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DOUBLE\_DIFFUSION"}, cs%double\_diffusion, &
2259                  \textcolor{stringliteral}{"If true, increase diffusivites for temperature or salt "}//&
2260                  \textcolor{stringliteral}{"based on double-diffusive parameterization from MOM4/KPP."}, &
2261                  default=.false.)
2262 
2263   \textcolor{keywordflow}{if} (cs%double\_diffusion) \textcolor{keywordflow}{then}
2264     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MAX\_RRHO\_SALT\_FINGERS"}, cs%Max\_Rrho\_salt\_fingers, &
2265                  \textcolor{stringliteral}{"Maximum density ratio for salt fingering regime."}, &
2266                  default=2.55, units=\textcolor{stringliteral}{"nondim"})
2267     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MAX\_SALT\_DIFF\_SALT\_FINGERS"}, cs%Max\_salt\_diff\_salt\_fingers, &
2268                  \textcolor{stringliteral}{"Maximum salt diffusivity for salt fingering regime."}, &
2269                  default=1.e-4, units=\textcolor{stringliteral}{"m2 s-1"}, scale=us%m2\_s\_to\_Z2\_T)
2270     \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KV\_MOLECULAR"}, cs%Kv\_molecular, &
2271                  \textcolor{stringliteral}{"Molecular viscosity for calculation of fluxes under double-diffusive "}//&
2272                  \textcolor{stringliteral}{"convection."}, default=1.5e-6, units=\textcolor{stringliteral}{"m2 s-1"}, scale=us%m2\_s\_to\_Z2\_T)
2273     \textcolor{comment}{! The default molecular viscosity follows the CCSM4.0 and MOM4p1 defaults.}
2274 \textcolor{keywordflow}{  endif} \textcolor{comment}{! old double-diffusion}
2275 
2276   \textcolor{keywordflow}{if} (cs%user\_change\_diff) \textcolor{keywordflow}{then}
2277     \textcolor{keyword}{call }user\_change\_diff\_init(time, g, gv, us, param\_file, diag, cs%user\_change\_diff\_CSp)
2278 \textcolor{keywordflow}{  endif}
2279 
2280   \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BRYAN\_LEWIS\_DIFFUSIVITY"}, bryan\_lewis\_diffusivity, &
2281                  \textcolor{stringliteral}{"If true, use a Bryan & Lewis (JGR 1979) like tanh "}//&
2282                  \textcolor{stringliteral}{"profile of background diapycnal diffusivity with depth. "}//&
2283                  \textcolor{stringliteral}{"This is done via CVMix."}, default=.false., do\_not\_log=.true.)
2284   \textcolor{keywordflow}{if} (cs%use\_tidal\_mixing .and. bryan\_lewis\_diffusivity) &
2285     \textcolor{keyword}{call }mom\_error(fatal,\textcolor{stringliteral}{"MOM\_Set\_Diffusivity: "}// &
2286          \textcolor{stringliteral}{"Bryan-Lewis and internal tidal dissipation are both enabled. Choose one."})
2287 
2288   cs%useKappaShear = kappa\_shear\_init(time, g, gv, us, param\_file, cs%diag, cs%kappaShear\_CSp)
2289   cs%Vertex\_Shear = kappa\_shear\_at\_vertex(param\_file)
2290 
2291   \textcolor{keywordflow}{if} (cs%useKappaShear) &
2292     id\_clock\_kappashear = cpu\_clock\_id(\textcolor{stringliteral}{'(Ocean kappa\_shear)'}, grain=clock\_module)
2293 
2294   \textcolor{comment}{! CVMix shear-driven mixing}
2295   cs%use\_CVMix\_shear = cvmix\_shear\_init(time, g, gv, us, param\_file, cs%diag, cs%CVMix\_shear\_csp)
2296 
2297   \textcolor{comment}{! CVMix double diffusion mixing}
2298   cs%use\_CVMix\_ddiff = cvmix\_ddiff\_init(time, g, gv, us, param\_file, cs%diag, cs%CVMix\_ddiff\_csp)
2299   \textcolor{keywordflow}{if} (cs%use\_CVMix\_ddiff) &
2300     id\_clock\_cvmix\_ddiff = cpu\_clock\_id(\textcolor{stringliteral}{'(Double diffusion via CVMix)'}, grain=clock\_module)
2301 
2302   \textcolor{keywordflow}{if} (cs%double\_diffusion .or. cs%use\_CVMix\_ddiff) \textcolor{keywordflow}{then}
2303     cs%id\_KT\_extra = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'KT\_extra'}, diag%axesTi, time, &
2304          \textcolor{stringliteral}{'Double-diffusive diffusivity for temperature'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s)
2305     cs%id\_KS\_extra = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'KS\_extra'}, diag%axesTi, time, &
2306          \textcolor{stringliteral}{'Double-diffusive diffusivity for salinity'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s)
2307 \textcolor{keywordflow}{  endif}
2308   \textcolor{keywordflow}{if} (cs%use\_CVMix\_ddiff) \textcolor{keywordflow}{then}
2309     cs%id\_R\_rho = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'R\_rho'}, diag%axesTi, time, &
2310          \textcolor{stringliteral}{'Double-diffusion density ratio'}, \textcolor{stringliteral}{'nondim'})
2311 \textcolor{keywordflow}{  endif}
2312 
2313   \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(halo\_ts)) \textcolor{keywordflow}{then}
2314     halo\_ts = 0
2315     \textcolor{keywordflow}{if} (cs%Vertex\_Shear) halo\_ts = 1
2316 \textcolor{keywordflow}{  endif}
2317 
\end{DoxyCode}
