\hypertarget{namespacemom__set__visc}{}\section{mom\+\_\+set\+\_\+visc Module Reference} \label{namespacemom__set__visc}\index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsection{Detailed Description} Calculates various values related to the bottom boundary layer, such as the viscosity and thickness of the B\+BL (set\+\_\+viscous\+\_\+\+B\+BL). This would also be the module in which other viscous quantities that are flow-\/independent might be set. This information is transmitted to other modules via a vertvisc type structure. The same code is used for the two velocity components, by indirectly referencing the velocities and defining a handful of direction-\/specific defined variables. \subsection*{Data Types} \begin{DoxyCompactItemize} \item type \hyperlink{structmom__set__visc_1_1set__visc__cs}{set\+\_\+visc\+\_\+cs} \begin{DoxyCompactList}\small\item\em Control structure for M\+O\+M\+\_\+set\+\_\+visc. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection*{Functions/\+Subroutines} \begin{DoxyCompactItemize} \item subroutine, public \hyperlink{namespacemom__set__visc_a9865113fe07928e7a240c2868ed45e5f}{set\+\_\+viscous\+\_\+bbl} (u, v, h, tv, visc, G, GV, US, CS, symmetrize) \begin{DoxyCompactList}\small\item\em Calculates the thickness of the bottom boundary layer and the viscosity within that layer. \end{DoxyCompactList}\item real function \hyperlink{namespacemom__set__visc_a0356a4e81cca9f7f31bbf87c717a6600}{set\+\_\+v\+\_\+at\+\_\+u} (v, h, G, i, j, k, mask2d\+Cv, O\+BC) \begin{DoxyCompactList}\small\item\em This subroutine finds a thickness-\/weighted value of v at the u-\/points. \end{DoxyCompactList}\item real function \hyperlink{namespacemom__set__visc_a46583b82467e74d8654c3c0a037a25cd}{set\+\_\+u\+\_\+at\+\_\+v} (u, h, G, i, j, k, mask2d\+Cu, O\+BC) \begin{DoxyCompactList}\small\item\em This subroutine finds a thickness-\/weighted value of u at the v-\/points. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__set__visc_aba41cd4f8baa1cda9036d97087ce8a22}{set\+\_\+viscous\+\_\+ml} (u, v, h, tv, forces, visc, dt, G, GV, US, CS, symmetrize) \begin{DoxyCompactList}\small\item\em Calculates the thickness of the surface boundary layer for applying an elevated viscosity. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__set__visc_ae2d9d9f74c1e9aec56257cfad372b0fd}{set\+\_\+visc\+\_\+register\+\_\+restarts} (HI, GV, param\+\_\+file, visc, restart\+\_\+\+CS) \begin{DoxyCompactList}\small\item\em Register any fields associated with the vertvisc\+\_\+type. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__set__visc_a2ec9be1a61c7d4a062aab4ff03e3ca29}{set\+\_\+visc\+\_\+init} (Time, G, GV, US, param\+\_\+file, diag, visc, CS, restart\+\_\+\+CS, O\+BC) \begin{DoxyCompactList}\small\item\em Initializes the M\+O\+M\+\_\+set\+\_\+visc control structure. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__set__visc_abbac7a7410557759c74444787de415a9}{set\+\_\+visc\+\_\+end} (visc, CS) \begin{DoxyCompactList}\small\item\em This subroutine dellocates any memory in the set\+\_\+visc control structure. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection{Function/\+Subroutine Documentation} \mbox{\Hypertarget{namespacemom__set__visc_a46583b82467e74d8654c3c0a037a25cd}\label{namespacemom__set__visc_a46583b82467e74d8654c3c0a037a25cd}} \index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}!set\+\_\+u\+\_\+at\+\_\+v@{set\+\_\+u\+\_\+at\+\_\+v}} \index{set\+\_\+u\+\_\+at\+\_\+v@{set\+\_\+u\+\_\+at\+\_\+v}!mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsubsection{\texorpdfstring{set\+\_\+u\+\_\+at\+\_\+v()}{set\_u\_at\_v()}} {\footnotesize\ttfamily real function mom\+\_\+set\+\_\+visc\+::set\+\_\+u\+\_\+at\+\_\+v (\begin{DoxyParamCaption}\item[{real, dimension( g \%isdb\+: g \%iedb, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{u, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h, }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G, }\item[{integer, intent(in)}]{i, }\item[{integer, intent(in)}]{j, }\item[{integer, intent(in)}]{k, }\item[{real, dimension( g \%isdb\+: g \%iedb, g \%jsd\+: g \%jed), intent(in)}]{mask2d\+Cu, }\item[{type(ocean\+\_\+obc\+\_\+type), pointer}]{O\+BC }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} This subroutine finds a thickness-\/weighted value of u at the v-\/points. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & The ocean\textquotesingle{}s grid structure\\ \hline \mbox{\tt in} & {\em u} & The zonal velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]} or other units.\\ \hline \mbox{\tt in} & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\ \hline \mbox{\tt in} & {\em i} & The i-\/index of the u-\/location to work on.\\ \hline \mbox{\tt in} & {\em j} & The j-\/index of the u-\/location to work on.\\ \hline \mbox{\tt in} & {\em k} & The k-\/index of the u-\/location to work on.\\ \hline \mbox{\tt in} & {\em mask2dcu} & A multiplicative mask of the u-\/points\\ \hline & {\em obc} & A pointer to an open boundary condition structure\\ \hline \end{DoxyParams} \begin{DoxyReturn}{Returns} The return value of u at v points in the same units as u, i.\+e. \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]} or other units. \end{DoxyReturn} Definition at line 1160 of file M\+O\+M\+\_\+set\+\_\+viscosity.\+F90. \begin{DoxyCode} 1160 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< The ocean's grid structure} 1161 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, & 1162 \textcolor{keywordtype}{intent(in)} :: u\textcolor{comment}{ !< The zonal velocity [L T-1 ~> m s-1] or other units.} 1163 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, & 1164 \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2]} 1165 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: i\textcolor{comment}{ !< The i-index of the u-location to work on.} 1166 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: j\textcolor{comment}{ !< The j-index of the u-location to work on.} 1167 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: k\textcolor{comment}{ !< The k-index of the u-location to work on.} 1168 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G))}, & 1169 \textcolor{keywordtype}{intent(in)} :: mask2dcu\textcolor{comment}{ !< A multiplicative mask of the u-points} 1170 \textcolor{keywordtype}{type}(ocean\_obc\_type), \textcolor{keywordtype}{pointer} :: obc\textcolor{comment}{ !< A pointer to an open boundary condition structure} 1171 \textcolor{keywordtype}{real} :: set\_u\_at\_v\textcolor{comment}{ !< The return value of u at v points in the} 1172 \textcolor{comment}{ !! same units as u, i.e. [L T-1 ~> m s-1] or other units.} 1173 1174 \textcolor{comment}{! This subroutine finds a thickness-weighted value of u at the v-points.} 1175 \textcolor{keywordtype}{real} :: hwt(-1:0,0:1) \textcolor{comment}{! Masked weights used to average u onto v [H ~> m or kg m-2].} 1176 \textcolor{keywordtype}{real} :: hwt\_tot \textcolor{comment}{! The sum of the masked thicknesses [H ~> m or kg m-2].} 1177 \textcolor{keywordtype}{integer} :: i0, j0, i1, j1 1178 1179 \textcolor{keywordflow}{do} j0 = 0,1 ; \textcolor{keywordflow}{do} i0 = -1,0 ; i1 = i+i0 ; j1 = j+j0 1180 hwt(i0,j0) = (h(i1,j1,k) + h(i1+1,j1,k)) * mask2dcu(i1,j1) 1181 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1182 1183 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (obc%number\_of\_segments > 0) \textcolor{keywordflow}{then} 1184 \textcolor{keywordflow}{do} j0 = 0,1 ; \textcolor{keywordflow}{do} i0 = -1,0 ; \textcolor{keywordflow}{if} ((obc%segnum\_u(i+i0,j+j0) /= obc\_none)) \textcolor{keywordflow}{then} 1185 i1 = i+i0 ; j1 = j+j0 1186 \textcolor{keywordflow}{if} (obc%segment(obc%segnum\_u(i1,j1))%direction == obc\_direction\_e) \textcolor{keywordflow}{then} 1187 hwt(i0,j0) = 2.0 * h(i1,j1,k) * mask2dcu(i1,j1) 1188 \textcolor{keywordflow}{elseif} (obc%segment(obc%segnum\_u(i1,j1))%direction == obc\_direction\_w) \textcolor{keywordflow}{then} 1189 hwt(i0,j0) = 2.0 * h(i1+1,j1,k) * mask2dcu(i1,j1) 1190 \textcolor{keywordflow}{ endif} 1191 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1192 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 1193 1194 hwt\_tot = (hwt(-1,0) + hwt(0,1)) + (hwt(0,0) + hwt(-1,1)) 1195 set\_u\_at\_v = 0.0 1196 \textcolor{keywordflow}{if} (hwt\_tot > 0.0) set\_u\_at\_v = & 1197 ((hwt(0,0) * u(i,j,k) + hwt(-1,1) * u(i-1,j+1,k)) + & 1198 (hwt(-1,0) * u(i-1,j,k) + hwt(0,1) * u(i,j+1,k))) / hwt\_tot 1199 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__set__visc_a0356a4e81cca9f7f31bbf87c717a6600}\label{namespacemom__set__visc_a0356a4e81cca9f7f31bbf87c717a6600}} \index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}!set\+\_\+v\+\_\+at\+\_\+u@{set\+\_\+v\+\_\+at\+\_\+u}} \index{set\+\_\+v\+\_\+at\+\_\+u@{set\+\_\+v\+\_\+at\+\_\+u}!mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsubsection{\texorpdfstring{set\+\_\+v\+\_\+at\+\_\+u()}{set\_v\_at\_u()}} {\footnotesize\ttfamily real function mom\+\_\+set\+\_\+visc\+::set\+\_\+v\+\_\+at\+\_\+u (\begin{DoxyParamCaption}\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(ocean\+\_\+grid\+\_\+type), intent(in)}]{G, }\item[{integer, intent(in)}]{i, }\item[{integer, intent(in)}]{j, }\item[{integer, intent(in)}]{k, }\item[{real, dimension(szi\+\_\+(g),szjb\+\_\+(g)), intent(in)}]{mask2d\+Cv, }\item[{type(ocean\+\_\+obc\+\_\+type), pointer}]{O\+BC }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} This subroutine finds a thickness-\/weighted value of v at the u-\/points. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & The ocean\textquotesingle{}s grid structure\\ \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 i} & The i-\/index of the u-\/location to work on.\\ \hline \mbox{\tt in} & {\em j} & The j-\/index of the u-\/location to work on.\\ \hline \mbox{\tt in} & {\em k} & The k-\/index of the u-\/location to work on.\\ \hline \mbox{\tt in} & {\em mask2dcv} & A multiplicative mask of the v-\/points\\ \hline & {\em obc} & A pointer to an open boundary condition structure\\ \hline \end{DoxyParams} \begin{DoxyReturn}{Returns} The return value of v at u points points in the same units as u, i.\+e. \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]} or other units. \end{DoxyReturn} Definition at line 1116 of file M\+O\+M\+\_\+set\+\_\+viscosity.\+F90. \begin{DoxyCode} 1116 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< The ocean's grid structure} 1117 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, & 1118 \textcolor{keywordtype}{intent(in)} :: v\textcolor{comment}{ !< The meridional velocity [L T-1 ~> m s-1]} 1119 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, & 1120 \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2]} 1121 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: i\textcolor{comment}{ !< The i-index of the u-location to work on.} 1122 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: j\textcolor{comment}{ !< The j-index of the u-location to work on.} 1123 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: k\textcolor{comment}{ !< The k-index of the u-location to work on.} 1124 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G))},& 1125 \textcolor{keywordtype}{intent(in)} :: mask2dcv\textcolor{comment}{ !< A multiplicative mask of the v-points} 1126 \textcolor{keywordtype}{type}(ocean\_obc\_type), \textcolor{keywordtype}{pointer} :: obc\textcolor{comment}{ !< A pointer to an open boundary condition structure} 1127 \textcolor{keywordtype}{real} :: set\_v\_at\_u\textcolor{comment}{ !< The return value of v at u points points in the} 1128 \textcolor{comment}{ !! same units as u, i.e. [L T-1 ~> m s-1] or other units.} 1129 1130 \textcolor{comment}{! This subroutine finds a thickness-weighted value of v at the u-points.} 1131 \textcolor{keywordtype}{real} :: hwt(0:1,-1:0) \textcolor{comment}{! Masked weights used to average u onto v [H ~> m or kg m-2].} 1132 \textcolor{keywordtype}{real} :: hwt\_tot \textcolor{comment}{! The sum of the masked thicknesses [H ~> m or kg m-2].} 1133 \textcolor{keywordtype}{integer} :: i0, j0, i1, j1 1134 1135 \textcolor{keywordflow}{do} j0 = -1,0 ; \textcolor{keywordflow}{do} i0 = 0,1 ; i1 = i+i0 ; j1 = j+j0 1136 hwt(i0,j0) = (h(i1,j1,k) + h(i1,j1+1,k)) * mask2dcv(i1,j1) 1137 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1138 1139 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (obc%number\_of\_segments > 0) \textcolor{keywordflow}{then} 1140 \textcolor{keywordflow}{do} j0 = -1,0 ; \textcolor{keywordflow}{do} i0 = 0,1 ; \textcolor{keywordflow}{if} ((obc%segnum\_v(i+i0,j+j0) /= obc\_none)) \textcolor{keywordflow}{then} 1141 i1 = i+i0 ; j1 = j+j0 1142 \textcolor{keywordflow}{if} (obc%segment(obc%segnum\_v(i1,j1))%direction == obc\_direction\_n) \textcolor{keywordflow}{then} 1143 hwt(i0,j0) = 2.0 * h(i1,j1,k) * mask2dcv(i1,j1) 1144 \textcolor{keywordflow}{elseif} (obc%segment(obc%segnum\_v(i1,j1))%direction == obc\_direction\_s) \textcolor{keywordflow}{then} 1145 hwt(i0,j0) = 2.0 * h(i1,j1+1,k) * mask2dcv(i1,j1) 1146 \textcolor{keywordflow}{ endif} 1147 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1148 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 1149 1150 hwt\_tot = (hwt(0,-1) + hwt(1,0)) + (hwt(1,-1) + hwt(0,0)) 1151 set\_v\_at\_u = 0.0 1152 \textcolor{keywordflow}{if} (hwt\_tot > 0.0) set\_v\_at\_u = & 1153 ((hwt(0,0) * v(i,j,k) + hwt(1,-1) * v(i+1,j-1,k)) + & 1154 (hwt(1,0) * v(i+1,j,k) + hwt(0,-1) * v(i,j-1,k))) / hwt\_tot 1155 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__set__visc_abbac7a7410557759c74444787de415a9}\label{namespacemom__set__visc_abbac7a7410557759c74444787de415a9}} \index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}!set\+\_\+visc\+\_\+end@{set\+\_\+visc\+\_\+end}} \index{set\+\_\+visc\+\_\+end@{set\+\_\+visc\+\_\+end}!mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsubsection{\texorpdfstring{set\+\_\+visc\+\_\+end()}{set\_visc\_end()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+visc\+::set\+\_\+visc\+\_\+end (\begin{DoxyParamCaption}\item[{type(vertvisc\+\_\+type), intent(inout)}]{visc, }\item[{type(\hyperlink{structmom__set__visc_1_1set__visc__cs}{set\+\_\+visc\+\_\+cs}), pointer}]{CS }\end{DoxyParamCaption})} This subroutine dellocates any memory in the set\+\_\+visc control structure. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in,out} & {\em visc} & A structure containing vertical viscosities and related fields. Elements are deallocated here.\\ \hline & {\em cs} & The control structure returned by a previous call to set\+\_\+visc\+\_\+init. \\ \hline \end{DoxyParams} Definition at line 2325 of file M\+O\+M\+\_\+set\+\_\+viscosity.\+F90. \begin{DoxyCode} 2325 \textcolor{keywordtype}{type}(vertvisc\_type), \textcolor{keywordtype}{intent(inout)} :: visc\textcolor{comment}{ !< A structure containing vertical viscosities and} 2326 \textcolor{comment}{ !! related fields. Elements are deallocated here.} 2327 \textcolor{keywordtype}{type}(set\_visc\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< The control structure returned by a previous} 2328 \textcolor{comment}{ !! call to set\_visc\_init.} 2329 \textcolor{keywordflow}{if} (cs%bottomdraglaw) \textcolor{keywordflow}{then} 2330 \textcolor{keyword}{deallocate}(visc%bbl\_thick\_u) ; \textcolor{keyword}{deallocate}(visc%bbl\_thick\_v) 2331 \textcolor{keyword}{deallocate}(visc%kv\_bbl\_u) ; \textcolor{keyword}{deallocate}(visc%kv\_bbl\_v) 2332 \textcolor{keywordflow}{if} (\textcolor{keyword}{allocated}(cs%bbl\_u)) \textcolor{keyword}{deallocate}(cs%bbl\_u) 2333 \textcolor{keywordflow}{if} (\textcolor{keyword}{allocated}(cs%bbl\_v)) \textcolor{keyword}{deallocate}(cs%bbl\_v) 2334 \textcolor{keywordflow}{ endif} 2335 \textcolor{keywordflow}{if} (cs%Channel\_drag) \textcolor{keywordflow}{then} 2336 \textcolor{keyword}{deallocate}(visc%Ray\_u) ; \textcolor{keyword}{deallocate}(visc%Ray\_v) 2337 \textcolor{keywordflow}{ endif} 2338 \textcolor{keywordflow}{if} (cs%dynamic\_viscous\_ML) \textcolor{keywordflow}{then} 2339 \textcolor{keyword}{deallocate}(visc%nkml\_visc\_u) ; \textcolor{keyword}{deallocate}(visc%nkml\_visc\_v) 2340 \textcolor{keywordflow}{ endif} 2341 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kd\_shear)) \textcolor{keyword}{deallocate}(visc%Kd\_shear) 2342 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_slow)) \textcolor{keyword}{deallocate}(visc%Kv\_slow) 2343 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%TKE\_turb)) \textcolor{keyword}{deallocate}(visc%TKE\_turb) 2344 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_shear)) \textcolor{keyword}{deallocate}(visc%Kv\_shear) 2345 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_shear\_Bu)) \textcolor{keyword}{deallocate}(visc%Kv\_shear\_Bu) 2346 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%ustar\_bbl)) \textcolor{keyword}{deallocate}(visc%ustar\_bbl) 2347 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%TKE\_bbl)) \textcolor{keyword}{deallocate}(visc%TKE\_bbl) 2348 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%taux\_shelf)) \textcolor{keyword}{deallocate}(visc%taux\_shelf) 2349 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%tauy\_shelf)) \textcolor{keyword}{deallocate}(visc%tauy\_shelf) 2350 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%tbl\_thick\_shelf\_u)) \textcolor{keyword}{deallocate}(visc%tbl\_thick\_shelf\_u) 2351 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%tbl\_thick\_shelf\_v)) \textcolor{keyword}{deallocate}(visc%tbl\_thick\_shelf\_v) 2352 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%kv\_tbl\_shelf\_u)) \textcolor{keyword}{deallocate}(visc%kv\_tbl\_shelf\_u) 2353 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%kv\_tbl\_shelf\_v)) \textcolor{keyword}{deallocate}(visc%kv\_tbl\_shelf\_v) 2354 2355 \textcolor{keyword}{deallocate}(cs) \end{DoxyCode} \mbox{\Hypertarget{namespacemom__set__visc_a2ec9be1a61c7d4a062aab4ff03e3ca29}\label{namespacemom__set__visc_a2ec9be1a61c7d4a062aab4ff03e3ca29}} \index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}!set\+\_\+visc\+\_\+init@{set\+\_\+visc\+\_\+init}} \index{set\+\_\+visc\+\_\+init@{set\+\_\+visc\+\_\+init}!mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsubsection{\texorpdfstring{set\+\_\+visc\+\_\+init()}{set\_visc\_init()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+visc\+::set\+\_\+visc\+\_\+init (\begin{DoxyParamCaption}\item[{type(time\+\_\+type), intent(in), target}]{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(vertvisc\+\_\+type), intent(inout)}]{visc, }\item[{type(\hyperlink{structmom__set__visc_1_1set__visc__cs}{set\+\_\+visc\+\_\+cs}), pointer}]{CS, }\item[{type(mom\+\_\+restart\+\_\+cs), pointer}]{restart\+\_\+\+CS, }\item[{type(ocean\+\_\+obc\+\_\+type), pointer}]{O\+BC }\end{DoxyParamCaption})} Initializes the M\+O\+M\+\_\+set\+\_\+visc control structure. \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 that is used to regulate diagnostic output.\\ \hline \mbox{\tt in,out} & {\em visc} & A structure containing vertical viscosities and related fields. Allocated here.\\ \hline & {\em cs} & A pointer that is set to point to the control structure for this module\\ \hline & {\em restart\+\_\+cs} & A pointer to the restart control structure.\\ \hline & {\em obc} & A pointer to an open boundary condition structure \\ \hline \end{DoxyParams} Definition at line 1973 of file M\+O\+M\+\_\+set\+\_\+viscosity.\+F90. \begin{DoxyCode} 1973 \textcolor{keywordtype}{type}(time\_type), \textcolor{keywordtype}{target}, \textcolor{keywordtype}{intent(in)} :: time\textcolor{comment}{ !< The current model time.} 1974 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(inout)} :: g\textcolor{comment}{ !< The ocean's grid structure.} 1975 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< The ocean's vertical grid structure.} 1976 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 1977 \textcolor{keywordtype}{type}(param\_file\_type), \textcolor{keywordtype}{intent(in)} :: param\_file\textcolor{comment}{ !< A structure to parse for run-time} 1978 \textcolor{comment}{ !! parameters.} 1979 \textcolor{keywordtype}{type}(diag\_ctrl), \textcolor{keywordtype}{target}, \textcolor{keywordtype}{intent(inout)} :: diag\textcolor{comment}{ !< A structure that is used to regulate diagnostic} 1980 \textcolor{comment}{ !! output.} 1981 \textcolor{keywordtype}{type}(vertvisc\_type), \textcolor{keywordtype}{intent(inout)} :: visc\textcolor{comment}{ !< A structure containing vertical viscosities and} 1982 \textcolor{comment}{ !! related fields. Allocated here.} 1983 \textcolor{keywordtype}{type}(set\_visc\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< A pointer that is set to point to the control} 1984 \textcolor{comment}{ !! structure for this module} 1985 \textcolor{keywordtype}{type}(mom\_restart\_cs), \textcolor{keywordtype}{pointer} :: restart\_cs\textcolor{comment}{ !< A pointer to the restart control structure.} 1986 \textcolor{keywordtype}{type}(ocean\_obc\_type), \textcolor{keywordtype}{pointer} :: obc\textcolor{comment}{ !< A pointer to an open boundary condition structure} 1987 1988 \textcolor{comment}{! Local variables} 1989 \textcolor{keywordtype}{real} :: csmag\_chan\_dflt, smag\_const1, tke\_decay\_dflt, bulk\_ri\_ml\_dflt 1990 \textcolor{keywordtype}{real} :: kv\_background 1991 \textcolor{keywordtype}{real} :: omega\_frac\_dflt 1992 \textcolor{keywordtype}{real} :: z\_rescale \textcolor{comment}{! A rescaling factor for heights from the representation in} 1993 \textcolor{comment}{! a restart file to the internal representation in this run.} 1994 \textcolor{keywordtype}{real} :: i\_t\_rescale \textcolor{comment}{! A rescaling factor for time from the internal representation in this run} 1995 \textcolor{comment}{! to the representation in a restart file.} 1996 \textcolor{keywordtype}{real} :: z2\_t\_rescale \textcolor{comment}{! A rescaling factor for vertical diffusivities and viscosities from the} 1997 \textcolor{comment}{! representation in a restart file to the internal representation in this run.} 1998 \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, n 1999 \textcolor{keywordtype}{integer} :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb, nz 2000 \textcolor{keywordtype}{logical} :: default\_2018\_answers 2001 \textcolor{keywordtype}{logical} :: use\_kappa\_shear, adiabatic, use\_omega, mle\_use\_pbl\_mld 2002 \textcolor{keywordtype}{logical} :: use\_cvmix\_ddiff, differential\_diffusion, use\_kpp 2003 \textcolor{keywordtype}{logical} :: use\_regridding 2004 \textcolor{keywordtype}{character(len=200)} :: filename, tideamp\_file 2005 \textcolor{keywordtype}{type}(obc\_segment\_type), \textcolor{keywordtype}{pointer} :: segment => null() \textcolor{comment}{! pointer to OBC segment type} 2006 \textcolor{comment}{! This include declares and sets the variable "version".} 2007 \textcolor{preprocessor}{# include "version\_variable.h"} 2008 \textcolor{preprocessor}{} \textcolor{keywordtype}{character(len=40)} :: mdl = \textcolor{stringliteral}{"MOM\_set\_visc"} \textcolor{comment}{! This module's name.} 2009 2010 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keywordflow}{then} 2011 \textcolor{keyword}{call }mom\_error(warning, \textcolor{stringliteral}{"set\_visc\_init called with an associated "}// & 2012 \textcolor{stringliteral}{"control structure."}) 2013 \textcolor{keywordflow}{return} 2014 \textcolor{keywordflow}{ endif} 2015 \textcolor{keyword}{allocate}(cs) 2016 2017 cs%OBC => obc 2018 2019 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec 2020 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = gv%ke 2021 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB 2022 2023 cs%diag => diag 2024 2025 \textcolor{comment}{! Set default, read and log parameters} 2026 \textcolor{keyword}{call }log\_version(param\_file, mdl, version, \textcolor{stringliteral}{""}) 2027 cs%RiNo\_mix = .false. ; use\_cvmix\_ddiff = .false. 2028 differential\_diffusion = .false. 2029 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"INPUTDIR"}, cs%inputdir, default=\textcolor{stringliteral}{"."}) 2030 cs%inputdir = slasher(cs%inputdir) 2031 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DEFAULT\_2018\_ANSWERS"}, default\_2018\_answers, & 2032 \textcolor{stringliteral}{"This sets the default value for the various \_2018\_ANSWERS parameters."}, & 2033 default=.false.) 2034 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"SET\_VISC\_2018\_ANSWERS"}, cs%answers\_2018, & 2035 \textcolor{stringliteral}{"If true, use the order of arithmetic and expressions that recover the "}//& 2036 \textcolor{stringliteral}{"answers from the end of 2018. Otherwise, use updated and more robust "}//& 2037 \textcolor{stringliteral}{"forms of the same expressions."}, default=default\_2018\_answers) 2038 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BOTTOMDRAGLAW"}, cs%bottomdraglaw, & 2039 \textcolor{stringliteral}{"If true, the bottom stress is calculated with a drag "}//& 2040 \textcolor{stringliteral}{"law of the form c\_drag*|u|*u. The velocity magnitude "}//& 2041 \textcolor{stringliteral}{"may be an assumed value or it may be based on the "}//& 2042 \textcolor{stringliteral}{"actual velocity in the bottommost HBBL, depending on "}//& 2043 \textcolor{stringliteral}{"LINEAR\_DRAG."}, default=.true.) 2044 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"CHANNEL\_DRAG"}, cs%Channel\_drag, & 2045 \textcolor{stringliteral}{"If true, the bottom drag is exerted directly on each "}//& 2046 \textcolor{stringliteral}{"layer proportional to the fraction of the bottom it "}//& 2047 \textcolor{stringliteral}{"overlies."}, default=.false.) 2048 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"LINEAR\_DRAG"}, cs%linear\_drag, & 2049 \textcolor{stringliteral}{"If LINEAR\_DRAG and BOTTOMDRAGLAW are defined the drag "}//& 2050 \textcolor{stringliteral}{"law is cdrag*DRAG\_BG\_VEL*u."}, default=.false.) 2051 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ADIABATIC"}, adiabatic, default=.false., & 2052 do\_not\_log=.true.) 2053 \textcolor{keywordflow}{if} (adiabatic) \textcolor{keywordflow}{then} 2054 \textcolor{keyword}{call }log\_param(param\_file, mdl, \textcolor{stringliteral}{"ADIABATIC"},adiabatic, & 2055 \textcolor{stringliteral}{"There are no diapycnal mass fluxes if ADIABATIC is "}//& 2056 \textcolor{stringliteral}{"true. This assumes that KD = KDML = 0.0 and that "}//& 2057 \textcolor{stringliteral}{"there is no buoyancy forcing, but makes the model "}//& 2058 \textcolor{stringliteral}{"faster by eliminating subroutine calls."}, default=.false.) 2059 \textcolor{keywordflow}{ endif} 2060 2061 \textcolor{keywordflow}{if} (.not.adiabatic) \textcolor{keywordflow}{then} 2062 cs%RiNo\_mix = kappa\_shear\_is\_used(param\_file) 2063 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DOUBLE\_DIFFUSION"}, differential\_diffusion, & 2064 \textcolor{stringliteral}{"If true, increase diffusivites for temperature or salt "}//& 2065 \textcolor{stringliteral}{"based on double-diffusive parameterization from MOM4/KPP."}, & 2066 default=.false.) 2067 use\_cvmix\_ddiff = cvmix\_ddiff\_is\_used(param\_file) 2068 \textcolor{keywordflow}{ endif} 2069 2070 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"PRANDTL\_TURB"}, visc%Prandtl\_turb, & 2071 \textcolor{stringliteral}{"The turbulent Prandtl number applied to shear "}//& 2072 \textcolor{stringliteral}{"instability."}, units=\textcolor{stringliteral}{"nondim"}, default=1.0) 2073 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DEBUG"}, cs%debug, default=.false.) 2074 2075 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DYNAMIC\_VISCOUS\_ML"}, cs%dynamic\_viscous\_ML, & 2076 \textcolor{stringliteral}{"If true, use a bulk Richardson number criterion to "}//& 2077 \textcolor{stringliteral}{"determine the mixed layer thickness for viscosity."}, & 2078 default=.false.) 2079 \textcolor{keywordflow}{if} (cs%dynamic\_viscous\_ML) \textcolor{keywordflow}{then} 2080 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BULK\_RI\_ML"}, bulk\_ri\_ml\_dflt, default=0.0) 2081 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BULK\_RI\_ML\_VISC"}, cs%bulk\_Ri\_ML, & 2082 \textcolor{stringliteral}{"The efficiency with which mean kinetic energy released "}//& 2083 \textcolor{stringliteral}{"by mechanically forced entrainment of the mixed layer "}//& 2084 \textcolor{stringliteral}{"is converted to turbulent kinetic energy. By default, "}//& 2085 \textcolor{stringliteral}{"BULK\_RI\_ML\_VISC = BULK\_RI\_ML or 0."}, units=\textcolor{stringliteral}{"nondim"}, & 2086 default=bulk\_ri\_ml\_dflt) 2087 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TKE\_DECAY"}, tke\_decay\_dflt, default=0.0) 2088 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TKE\_DECAY\_VISC"}, cs%TKE\_decay, & 2089 \textcolor{stringliteral}{"TKE\_DECAY\_VISC relates the vertical rate of decay of "}//& 2090 \textcolor{stringliteral}{"the TKE available for mechanical entrainment to the "}//& 2091 \textcolor{stringliteral}{"natural Ekman depth for use in calculating the dynamic "}//& 2092 \textcolor{stringliteral}{"mixed layer viscosity. By default, "}//& 2093 \textcolor{stringliteral}{"TKE\_DECAY\_VISC = TKE\_DECAY or 0."}, units=\textcolor{stringliteral}{"nondim"}, & 2094 default=tke\_decay\_dflt) 2095 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_USE\_OMEGA"}, use\_omega, & 2096 \textcolor{stringliteral}{"If true, use the absolute rotation rate instead of the "}//& 2097 \textcolor{stringliteral}{"vertical component of rotation when setting the decay "}//& 2098 \textcolor{stringliteral}{"scale for turbulence."}, default=.false., do\_not\_log=.true.) 2099 omega\_frac\_dflt = 0.0 2100 \textcolor{keywordflow}{if} (use\_omega) \textcolor{keywordflow}{then} 2101 \textcolor{keyword}{call }mom\_error(warning, \textcolor{stringliteral}{"ML\_USE\_OMEGA is depricated; use ML\_OMEGA\_FRAC=1.0 instead."}) 2102 omega\_frac\_dflt = 1.0 2103 \textcolor{keywordflow}{ endif} 2104 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ML\_OMEGA\_FRAC"}, cs%omega\_frac, & 2105 \textcolor{stringliteral}{"When setting the decay scale for turbulence, use this "}//& 2106 \textcolor{stringliteral}{"fraction of the absolute rotation rate blended with the "}//& 2107 \textcolor{stringliteral}{"local value of f, as sqrt((1-of)*f^2 + of*4*omega^2)."}, & 2108 units=\textcolor{stringliteral}{"nondim"}, default=omega\_frac\_dflt) 2109 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"OMEGA"}, cs%omega, & 2110 \textcolor{stringliteral}{"The rotation rate of the earth."}, units=\textcolor{stringliteral}{"s-1"}, & 2111 default=7.2921e-5, scale=us%T\_to\_s) 2112 \textcolor{comment}{! This give a minimum decay scale that is typically much less than Angstrom.} 2113 cs%ustar\_min = 2e-4*cs%omega*(gv%Angstrom\_Z + gv%H\_to\_Z*gv%H\_subroundoff) 2114 \textcolor{keywordflow}{else} 2115 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"OMEGA"}, cs%omega, & 2116 \textcolor{stringliteral}{"The rotation rate of the earth."}, units=\textcolor{stringliteral}{"s-1"}, & 2117 default=7.2921e-5, scale=us%T\_to\_s) 2118 \textcolor{keywordflow}{ endif} 2119 2120 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"HBBL"}, cs%Hbbl, & 2121 \textcolor{stringliteral}{"The thickness of a bottom boundary layer with a "}//& 2122 \textcolor{stringliteral}{"viscosity of KVBBL if BOTTOMDRAGLAW is not defined, or "}//& 2123 \textcolor{stringliteral}{"the thickness over which near-bottom velocities are "}//& 2124 \textcolor{stringliteral}{"averaged for the drag law if BOTTOMDRAGLAW is defined "}//& 2125 \textcolor{stringliteral}{"but LINEAR\_DRAG is not."}, units=\textcolor{stringliteral}{"m"}, fail\_if\_missing=.true.) \textcolor{comment}{! Rescaled later} 2126 \textcolor{keywordflow}{if} (cs%bottomdraglaw) \textcolor{keywordflow}{then} 2127 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"CDRAG"}, cs%cdrag, & 2128 \textcolor{stringliteral}{"CDRAG is the drag coefficient relating the magnitude of "}//& 2129 \textcolor{stringliteral}{"the velocity field to the bottom stress. CDRAG is only "}//& 2130 \textcolor{stringliteral}{"used if BOTTOMDRAGLAW is defined."}, units=\textcolor{stringliteral}{"nondim"}, & 2131 default=0.003) 2132 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BBL\_USE\_TIDAL\_BG"}, cs%BBL\_use\_tidal\_bg, & 2133 \textcolor{stringliteral}{"Flag to use the tidal RMS amplitude in place of constant "}//& 2134 \textcolor{stringliteral}{"background velocity for computing u* in the BBL. "}//& 2135 \textcolor{stringliteral}{"This flag is only used when BOTTOMDRAGLAW is true and "}//& 2136 \textcolor{stringliteral}{"LINEAR\_DRAG is false."}, default=.false.) 2137 \textcolor{keywordflow}{if} (cs%BBL\_use\_tidal\_bg) \textcolor{keywordflow}{then} 2138 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TIDEAMP\_FILE"}, tideamp\_file, & 2139 \textcolor{stringliteral}{"The path to the file containing the spatially varying "}//& 2140 \textcolor{stringliteral}{"tidal amplitudes with INT\_TIDE\_DISSIPATION."}, default=\textcolor{stringliteral}{"tideamp.nc"}) 2141 \textcolor{keywordflow}{else} 2142 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DRAG\_BG\_VEL"}, cs%drag\_bg\_vel, & 2143 \textcolor{stringliteral}{"DRAG\_BG\_VEL is either the assumed bottom velocity (with "}//& 2144 \textcolor{stringliteral}{"LINEAR\_DRAG) or an unresolved velocity that is "}//& 2145 \textcolor{stringliteral}{"combined with the resolved velocity to estimate the "}//& 2146 \textcolor{stringliteral}{"velocity magnitude. DRAG\_BG\_VEL is only used when "}//& 2147 \textcolor{stringliteral}{"BOTTOMDRAGLAW is defined."}, units=\textcolor{stringliteral}{"m s-1"}, default=0.0, scale=us%m\_s\_to\_L\_T) 2148 \textcolor{keywordflow}{ endif} 2149 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_REGRIDDING"}, use\_regridding, & 2150 do\_not\_log = .true., default = .false. ) 2151 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BBL\_USE\_EOS"}, cs%BBL\_use\_EOS, & 2152 \textcolor{stringliteral}{"If true, use the equation of state in determining the "}//& 2153 \textcolor{stringliteral}{"properties of the bottom boundary layer. Otherwise use "}//& 2154 \textcolor{stringliteral}{"the layer target potential densities. The default of "}//& 2155 \textcolor{stringliteral}{"this is determined by USE\_REGRIDDING."}, default=use\_regridding) 2156 \textcolor{keywordflow}{if} (use\_regridding .and. (.not. cs%BBL\_use\_EOS)) & 2157 \textcolor{keyword}{call }mom\_error(fatal,\textcolor{stringliteral}{"When using MOM6 in ALE mode it is required to "}//& 2158 \textcolor{stringliteral}{"set BBL\_USE\_EOS to True"}) 2159 \textcolor{keywordflow}{ endif} 2160 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BBL\_THICK\_MIN"}, cs%BBL\_thick\_min, & 2161 \textcolor{stringliteral}{"The minimum bottom boundary layer thickness that can be "}//& 2162 \textcolor{stringliteral}{"used with BOTTOMDRAGLAW. This might be "}//& 2163 \textcolor{stringliteral}{"Kv/(cdrag*drag\_bg\_vel) to give Kv as the minimum "}//& 2164 \textcolor{stringliteral}{"near-bottom viscosity."}, units=\textcolor{stringliteral}{"m"}, default=0.0) \textcolor{comment}{! Rescaled later} 2165 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"HTBL\_SHELF\_MIN"}, cs%Htbl\_shelf\_min, & 2166 \textcolor{stringliteral}{"The minimum top boundary layer thickness that can be "}//& 2167 \textcolor{stringliteral}{"used with BOTTOMDRAGLAW. This might be "}//& 2168 \textcolor{stringliteral}{"Kv/(cdrag*drag\_bg\_vel) to give Kv as the minimum "}//& 2169 \textcolor{stringliteral}{"near-top viscosity."}, units=\textcolor{stringliteral}{"m"}, default=cs%BBL\_thick\_min, scale=gv%m\_to\_H) 2170 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"HTBL\_SHELF"}, cs%Htbl\_shelf, & 2171 \textcolor{stringliteral}{"The thickness over which near-surface velocities are "}//& 2172 \textcolor{stringliteral}{"averaged for the drag law under an ice shelf. By "}//& 2173 \textcolor{stringliteral}{"default this is the same as HBBL"}, units=\textcolor{stringliteral}{"m"}, default=cs%Hbbl, scale=gv%m\_to\_H) 2174 \textcolor{comment}{! These unit conversions are out outside the get\_param calls because the are also defaults.} 2175 cs%Hbbl = cs%Hbbl * gv%m\_to\_H \textcolor{comment}{! Rescale} 2176 cs%BBL\_thick\_min = cs%BBL\_thick\_min * gv%m\_to\_H \textcolor{comment}{! Rescale} 2177 2178 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KV"}, kv\_background, & 2179 \textcolor{stringliteral}{"The background kinematic viscosity in the interior. "}//& 2180 \textcolor{stringliteral}{"The molecular value, ~1e-6 m2 s-1, may be used."}, & 2181 units=\textcolor{stringliteral}{"m2 s-1"}, fail\_if\_missing=.true.) 2182 2183 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_KPP"}, use\_kpp, & 2184 \textcolor{stringliteral}{"If true, turns on the [CVMix] KPP scheme of Large et al., 1994, "}//& 2185 \textcolor{stringliteral}{"to calculate diffusivities and non-local transport in the OBL."}, & 2186 do\_not\_log=.true., default=.false.) 2187 2188 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KV\_BBL\_MIN"}, cs%KV\_BBL\_min, & 2189 \textcolor{stringliteral}{"The minimum viscosities in the bottom boundary layer."}, & 2190 units=\textcolor{stringliteral}{"m2 s-1"}, default=kv\_background, scale=us%m2\_s\_to\_Z2\_T) 2191 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KV\_TBL\_MIN"}, cs%KV\_TBL\_min, & 2192 \textcolor{stringliteral}{"The minimum viscosities in the top boundary layer."}, & 2193 units=\textcolor{stringliteral}{"m2 s-1"}, default=kv\_background, scale=us%m2\_s\_to\_Z2\_T) 2194 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"CORRECT\_BBL\_BOUNDS"}, cs%correct\_BBL\_bounds, & 2195 \textcolor{stringliteral}{"If true, uses the correct bounds on the BBL thickness and "}//& 2196 \textcolor{stringliteral}{"viscosity so that the bottom layer feels the intended drag."}, & 2197 default=.false.) 2198 2199 \textcolor{keywordflow}{if} (cs%Channel\_drag) \textcolor{keywordflow}{then} 2200 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"SMAG\_LAP\_CONST"}, smag\_const1, default=-1.0) 2201 2202 csmag\_chan\_dflt = 0.15 2203 \textcolor{keywordflow}{if} (smag\_const1 >= 0.0) csmag\_chan\_dflt = smag\_const1 2204 2205 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"SMAG\_CONST\_CHANNEL"}, cs%c\_Smag, & 2206 \textcolor{stringliteral}{"The nondimensional Laplacian Smagorinsky constant used "}//& 2207 \textcolor{stringliteral}{"in calculating the channel drag if it is enabled. The "}//& 2208 \textcolor{stringliteral}{"default is to use the same value as SMAG\_LAP\_CONST if "}//& 2209 \textcolor{stringliteral}{"it is defined, or 0.15 if it is not. The value used is "}//& 2210 \textcolor{stringliteral}{"also 0.15 if the specified value is negative."}, & 2211 units=\textcolor{stringliteral}{"nondim"}, default=csmag\_chan\_dflt) 2212 \textcolor{keywordflow}{if} (cs%c\_Smag < 0.0) cs%c\_Smag = 0.15 2213 \textcolor{keywordflow}{ endif} 2214 2215 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MLE\_USE\_PBL\_MLD"}, mle\_use\_pbl\_mld, & 2216 default=.false., do\_not\_log=.true.) 2217 2218 \textcolor{keywordflow}{if} (cs%RiNo\_mix .and. kappa\_shear\_at\_vertex(param\_file)) \textcolor{keywordflow}{then} 2219 \textcolor{comment}{! This is necessary for reproduciblity across restarts in non-symmetric mode.} 2220 \textcolor{keyword}{call }pass\_var(visc%Kv\_shear\_Bu, g%Domain, position=corner, complete=.true.) 2221 \textcolor{keywordflow}{ endif} 2222 2223 \textcolor{keywordflow}{if} (cs%bottomdraglaw) \textcolor{keywordflow}{then} 2224 \textcolor{keyword}{allocate}(visc%bbl\_thick\_u(isdb:iedb,jsd:jed)) ; visc%bbl\_thick\_u(:,:) = 0.0 2225 \textcolor{keyword}{allocate}(visc%kv\_bbl\_u(isdb:iedb,jsd:jed)) ; visc%kv\_bbl\_u(:,:) = 0.0 2226 \textcolor{keyword}{allocate}(visc%bbl\_thick\_v(isd:ied,jsdb:jedb)) ; visc%bbl\_thick\_v(:,:) = 0.0 2227 \textcolor{keyword}{allocate}(visc%kv\_bbl\_v(isd:ied,jsdb:jedb)) ; visc%kv\_bbl\_v(:,:) = 0.0 2228 \textcolor{keyword}{allocate}(visc%ustar\_bbl(isd:ied,jsd:jed)) ; visc%ustar\_bbl(:,:) = 0.0 2229 \textcolor{keyword}{allocate}(visc%TKE\_bbl(isd:ied,jsd:jed)) ; visc%TKE\_bbl(:,:) = 0.0 2230 2231 cs%id\_bbl\_thick\_u = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'bbl\_thick\_u'}, & 2232 diag%axesCu1, time, \textcolor{stringliteral}{'BBL thickness at u points'}, \textcolor{stringliteral}{'m'}, conversion=us%Z\_to\_m) 2233 cs%id\_kv\_bbl\_u = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'kv\_bbl\_u'}, diag%axesCu1, & 2234 time, \textcolor{stringliteral}{'BBL viscosity at u points'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s) 2235 cs%id\_bbl\_u = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'bbl\_u'}, diag%axesCu1, & 2236 time, \textcolor{stringliteral}{'BBL mean u current'}, \textcolor{stringliteral}{'m s-1'}, conversion=us%L\_T\_to\_m\_s) 2237 \textcolor{keywordflow}{if} (cs%id\_bbl\_u>0) \textcolor{keywordflow}{then} 2238 \textcolor{keyword}{allocate}(cs%bbl\_u(isdb:iedb,jsd:jed)) ; cs%bbl\_u(:,:) = 0.0 2239 \textcolor{keywordflow}{ endif} 2240 cs%id\_bbl\_thick\_v = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'bbl\_thick\_v'}, & 2241 diag%axesCv1, time, \textcolor{stringliteral}{'BBL thickness at v points'}, \textcolor{stringliteral}{'m'}, conversion=us%Z\_to\_m) 2242 cs%id\_kv\_bbl\_v = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'kv\_bbl\_v'}, diag%axesCv1, & 2243 time, \textcolor{stringliteral}{'BBL viscosity at v points'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s) 2244 cs%id\_bbl\_v = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'bbl\_v'}, diag%axesCv1, & 2245 time, \textcolor{stringliteral}{'BBL mean v current'}, \textcolor{stringliteral}{'m s-1'}, conversion=us%L\_T\_to\_m\_s) 2246 \textcolor{keywordflow}{if} (cs%id\_bbl\_v>0) \textcolor{keywordflow}{then} 2247 \textcolor{keyword}{allocate}(cs%bbl\_v(isd:ied,jsdb:jedb)) ; cs%bbl\_v(:,:) = 0.0 2248 \textcolor{keywordflow}{ endif} 2249 \textcolor{keywordflow}{if} (cs%BBL\_use\_tidal\_bg) \textcolor{keywordflow}{then} 2250 \textcolor{keyword}{allocate}(cs%tideamp(isd:ied,jsd:jed)) ; cs%tideamp(:,:) = 0.0 2251 filename = trim(cs%inputdir) // trim(tideamp\_file) 2252 \textcolor{keyword}{call }log\_param(param\_file, mdl, \textcolor{stringliteral}{"INPUTDIR/TIDEAMP\_FILE"}, filename) 2253 \textcolor{keyword}{call }mom\_read\_data(filename, \textcolor{stringliteral}{'tideamp'}, cs%tideamp, g%domain, timelevel=1, scale=us%m\_to\_Z*us%T\_to\_s) 2254 \textcolor{keyword}{call }pass\_var(cs%tideamp,g%domain) 2255 \textcolor{keywordflow}{ endif} 2256 \textcolor{keywordflow}{ endif} 2257 \textcolor{keywordflow}{if} (cs%Channel\_drag) \textcolor{keywordflow}{then} 2258 \textcolor{keyword}{allocate}(visc%Ray\_u(isdb:iedb,jsd:jed,nz)) ; visc%Ray\_u(:,:,:) = 0.0 2259 \textcolor{keyword}{allocate}(visc%Ray\_v(isd:ied,jsdb:jedb,nz)) ; visc%Ray\_v(:,:,:) = 0.0 2260 cs%id\_Ray\_u = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Rayleigh\_u'}, diag%axesCuL, & 2261 time, \textcolor{stringliteral}{'Rayleigh drag velocity at u points'}, \textcolor{stringliteral}{'m s-1'}, conversion=us%Z\_to\_m*us%s\_to\_T) 2262 cs%id\_Ray\_v = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'Rayleigh\_v'}, diag%axesCvL, & 2263 time, \textcolor{stringliteral}{'Rayleigh drag velocity at v points'}, \textcolor{stringliteral}{'m s-1'}, conversion=us%Z\_to\_m*us%s\_to\_T) 2264 \textcolor{keywordflow}{ endif} 2265 2266 \textcolor{keywordflow}{if} (use\_cvmix\_ddiff .or. differential\_diffusion) \textcolor{keywordflow}{then} 2267 \textcolor{keyword}{allocate}(visc%Kd\_extra\_T(isd:ied,jsd:jed,nz+1)) ; visc%Kd\_extra\_T(:,:,:) = 0.0 2268 \textcolor{keyword}{allocate}(visc%Kd\_extra\_S(isd:ied,jsd:jed,nz+1)) ; visc%Kd\_extra\_S(:,:,:) = 0.0 2269 \textcolor{keywordflow}{ endif} 2270 2271 \textcolor{keywordflow}{if} (cs%dynamic\_viscous\_ML) \textcolor{keywordflow}{then} 2272 \textcolor{keyword}{allocate}(visc%nkml\_visc\_u(isdb:iedb,jsd:jed)) ; visc%nkml\_visc\_u(:,:) = 0.0 2273 \textcolor{keyword}{allocate}(visc%nkml\_visc\_v(isd:ied,jsdb:jedb)) ; visc%nkml\_visc\_v(:,:) = 0.0 2274 cs%id\_nkml\_visc\_u = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'nkml\_visc\_u'}, & 2275 diag%axesCu1, time, \textcolor{stringliteral}{'Number of layers in viscous mixed layer at u points'}, \textcolor{stringliteral}{'m'}) 2276 cs%id\_nkml\_visc\_v = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'nkml\_visc\_v'}, & 2277 diag%axesCv1, time, \textcolor{stringliteral}{'Number of layers in viscous mixed layer at v points'}, \textcolor{stringliteral}{'m'}) 2278 \textcolor{keywordflow}{ endif} 2279 2280 \textcolor{keyword}{call }register\_restart\_field\_as\_obsolete(\textcolor{stringliteral}{'Kd\_turb'},\textcolor{stringliteral}{'Kd\_shear'}, restart\_cs) 2281 \textcolor{keyword}{call }register\_restart\_field\_as\_obsolete(\textcolor{stringliteral}{'Kv\_turb'},\textcolor{stringliteral}{'Kv\_shear'}, restart\_cs) 2282 2283 \textcolor{comment}{! Account for possible changes in dimensional scaling for variables that have been} 2284 \textcolor{comment}{! read from a restart file.} 2285 z\_rescale = 1.0 2286 \textcolor{keywordflow}{if} ((us%m\_to\_Z\_restart /= 0.0) .and. (us%m\_to\_Z\_restart /= us%m\_to\_Z)) & 2287 z\_rescale = us%m\_to\_Z / us%m\_to\_Z\_restart 2288 i\_t\_rescale = 1.0 2289 \textcolor{keywordflow}{if} ((us%s\_to\_T\_restart /= 0.0) .and. (us%s\_to\_T\_restart /= us%s\_to\_T)) & 2290 i\_t\_rescale = us%s\_to\_T\_restart / us%s\_to\_T 2291 z2\_t\_rescale = z\_rescale**2*i\_t\_rescale 2292 2293 \textcolor{keywordflow}{if} (z2\_t\_rescale /= 1.0) \textcolor{keywordflow}{then} 2294 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kd\_shear)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (query\_initialized(visc%Kd\_shear, \textcolor{stringliteral}{"Kd\_shear"}, restart\_cs)) \textcolor{keywordflow}{ then} 2295 \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=is,ie 2296 visc%Kd\_shear(i,j,k) = z2\_t\_rescale * visc%Kd\_shear(i,j,k) 2297 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 2298 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 2299 2300 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_shear)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (query\_initialized(visc%Kv\_shear, \textcolor{stringliteral}{"Kv\_shear"}, restart\_cs)) \textcolor{keywordflow}{ then} 2301 \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=is,ie 2302 visc%Kv\_shear(i,j,k) = z2\_t\_rescale * visc%Kv\_shear(i,j,k) 2303 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 2304 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 2305 2306 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Kv\_shear\_Bu)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (query\_initialized(visc%Kv\_shear\_Bu, \textcolor{stringliteral}{"Kv\_shear\_Bu"}, restart\_cs)) \textcolor{keywordflow}{then} 2307 \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} j=js-1,je ; \textcolor{keywordflow}{do} i=is-1,ie 2308 visc%Kv\_shear\_Bu(i,j,k) = z2\_t\_rescale * visc%Kv\_shear\_Bu(i,j,k) 2309 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 2310 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 2311 \textcolor{keywordflow}{ endif} 2312 2313 \textcolor{keywordflow}{if} (mle\_use\_pbl\_mld .and. (z\_rescale /= 1.0)) \textcolor{keywordflow}{then} 2314 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%MLD)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (query\_initialized(visc%MLD, \textcolor{stringliteral}{"MLD"}, restart\_cs)) \textcolor{keywordflow}{then} 2315 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=is,ie 2316 visc%MLD(i,j) = z\_rescale * visc%MLD(i,j) 2317 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 2318 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 2319 \textcolor{keywordflow}{ endif} 2320 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__set__visc_ae2d9d9f74c1e9aec56257cfad372b0fd}\label{namespacemom__set__visc_ae2d9d9f74c1e9aec56257cfad372b0fd}} \index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}!set\+\_\+visc\+\_\+register\+\_\+restarts@{set\+\_\+visc\+\_\+register\+\_\+restarts}} \index{set\+\_\+visc\+\_\+register\+\_\+restarts@{set\+\_\+visc\+\_\+register\+\_\+restarts}!mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsubsection{\texorpdfstring{set\+\_\+visc\+\_\+register\+\_\+restarts()}{set\_visc\_register\_restarts()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+visc\+::set\+\_\+visc\+\_\+register\+\_\+restarts (\begin{DoxyParamCaption}\item[{type(hor\+\_\+index\+\_\+type), intent(in)}]{HI, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file, }\item[{type(vertvisc\+\_\+type), intent(inout)}]{visc, }\item[{type(mom\+\_\+restart\+\_\+cs), pointer}]{restart\+\_\+\+CS }\end{DoxyParamCaption})} Register any fields associated with the vertvisc\+\_\+type. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em hi} & A horizontal index type structure.\\ \hline \mbox{\tt in} & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\ \hline \mbox{\tt in} & {\em param\+\_\+file} & A structure to parse for run-\/time parameters.\\ \hline \mbox{\tt in,out} & {\em visc} & A structure containing vertical viscosities and related fields. Allocated here.\\ \hline & {\em restart\+\_\+cs} & A pointer to the restart control structure. \\ \hline \end{DoxyParams} Definition at line 1887 of file M\+O\+M\+\_\+set\+\_\+viscosity.\+F90. \begin{DoxyCode} 1887 \textcolor{keywordtype}{type}(hor\_index\_type), \textcolor{keywordtype}{intent(in)} :: hi\textcolor{comment}{ !< A horizontal index type structure.} 1888 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< The ocean's vertical grid structure.} 1889 \textcolor{keywordtype}{type}(param\_file\_type), \textcolor{keywordtype}{intent(in)} :: param\_file\textcolor{comment}{ !< A structure to parse for run-time} 1890 \textcolor{comment}{ !! parameters.} 1891 \textcolor{keywordtype}{type}(vertvisc\_type), \textcolor{keywordtype}{intent(inout)} :: visc\textcolor{comment}{ !< A structure containing vertical} 1892 \textcolor{comment}{ !! viscosities and related fields.} 1893 \textcolor{comment}{ !! Allocated here.} 1894 \textcolor{keywordtype}{type}(mom\_restart\_cs), \textcolor{keywordtype}{pointer} :: restart\_cs\textcolor{comment}{ !< A pointer to the restart control structure.} 1895 \textcolor{comment}{! Local variables} 1896 \textcolor{keywordtype}{logical} :: use\_kappa\_shear, ks\_at\_vertex 1897 \textcolor{keywordtype}{logical} :: adiabatic, usekpp, useepbl 1898 \textcolor{keywordtype}{logical} :: use\_cvmix\_shear, mle\_use\_pbl\_mld, use\_cvmix\_conv 1899 \textcolor{keywordtype}{integer} :: isd, ied, jsd, jed, nz 1900 \textcolor{keywordtype}{real} :: hfreeze\textcolor{comment}{ !< If hfreeze > 0 [m], melt potential will be computed.} 1901 \textcolor{keywordtype}{character(len=40)} :: mdl = \textcolor{stringliteral}{"MOM\_set\_visc"} \textcolor{comment}{! This module's name.} 1902 isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke 1903 1904 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ADIABATIC"}, adiabatic, default=.false., & 1905 do\_not\_log=.true.) 1906 1907 use\_kappa\_shear = .false. ; ks\_at\_vertex = .false. ; use\_cvmix\_shear = .false. 1908 usekpp = .false. ; useepbl = .false. ; use\_cvmix\_conv = .false. 1909 1910 \textcolor{keywordflow}{if} (.not.adiabatic) \textcolor{keywordflow}{then} 1911 use\_kappa\_shear = kappa\_shear\_is\_used(param\_file) 1912 ks\_at\_vertex = kappa\_shear\_at\_vertex(param\_file) 1913 use\_cvmix\_shear = cvmix\_shear\_is\_used(param\_file) 1914 use\_cvmix\_conv = cvmix\_conv\_is\_used(param\_file) 1915 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_KPP"}, usekpp, & 1916 \textcolor{stringliteral}{"If true, turns on the [CVMix] KPP scheme of Large et al., 1994, "}//& 1917 \textcolor{stringliteral}{"to calculate diffusivities and non-local transport in the OBL."}, & 1918 default=.false., do\_not\_log=.true.) 1919 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"ENERGETICS\_SFC\_PBL"}, useepbl, & 1920 \textcolor{stringliteral}{"If true, use an implied energetics planetary boundary "}//& 1921 \textcolor{stringliteral}{"layer scheme to determine the diffusivity and viscosity "}//& 1922 \textcolor{stringliteral}{"in the surface boundary layer."}, default=.false., do\_not\_log=.true.) 1923 \textcolor{keywordflow}{ endif} 1924 1925 \textcolor{keywordflow}{if} (use\_kappa\_shear .or. usekpp .or. useepbl .or. use\_cvmix\_shear .or. use\_cvmix\_conv) \textcolor{keywordflow}{then} 1926 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%Kd\_shear, isd, ied, jsd, jed, nz+1) 1927 \textcolor{keyword}{call }register\_restart\_field(visc%Kd\_shear, \textcolor{stringliteral}{"Kd\_shear"}, .false., restart\_cs, & 1928 \textcolor{stringliteral}{"Shear-driven turbulent diffusivity at interfaces"}, \textcolor{stringliteral}{"m2 s-1"}, z\_grid=\textcolor{stringliteral}{'i'}) 1929 \textcolor{keywordflow}{ endif} 1930 \textcolor{keywordflow}{if} (usekpp .or. useepbl .or. use\_cvmix\_shear .or. use\_cvmix\_conv .or. & 1931 (use\_kappa\_shear .and. .not.ks\_at\_vertex )) \textcolor{keywordflow}{then} 1932 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%Kv\_shear, isd, ied, jsd, jed, nz+1) 1933 \textcolor{keyword}{call }register\_restart\_field(visc%Kv\_shear, \textcolor{stringliteral}{"Kv\_shear"}, .false., restart\_cs, & 1934 \textcolor{stringliteral}{"Shear-driven turbulent viscosity at interfaces"}, \textcolor{stringliteral}{"m2 s-1"}, z\_grid=\textcolor{stringliteral}{'i'}) 1935 \textcolor{keywordflow}{ endif} 1936 \textcolor{keywordflow}{if} (use\_kappa\_shear .and. ks\_at\_vertex) \textcolor{keywordflow}{then} 1937 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%TKE\_turb, hi%IsdB, hi%IedB, hi%JsdB, hi%JedB, nz+1) 1938 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%Kv\_shear\_Bu, hi%IsdB, hi%IedB, hi%JsdB, hi%JedB, nz+1) 1939 \textcolor{keyword}{call }register\_restart\_field(visc%Kv\_shear\_Bu, \textcolor{stringliteral}{"Kv\_shear\_Bu"}, .false., restart\_cs, & 1940 \textcolor{stringliteral}{"Shear-driven turbulent viscosity at vertex interfaces"}, \textcolor{stringliteral}{"m2 s-1"}, & 1941 hor\_grid=\textcolor{stringliteral}{"Bu"}, z\_grid=\textcolor{stringliteral}{'i'}) 1942 \textcolor{keywordflow}{elseif} (use\_kappa\_shear) \textcolor{keywordflow}{then} 1943 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%TKE\_turb, isd, ied, jsd, jed, nz+1) 1944 \textcolor{keywordflow}{ endif} 1945 1946 \textcolor{keywordflow}{if} (usekpp) \textcolor{keywordflow}{then} 1947 \textcolor{comment}{! MOM\_bkgnd\_mixing uses Kv\_slow when KPP is defined.} 1948 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%Kv\_slow, isd, ied, jsd, jed, nz+1) 1949 \textcolor{keywordflow}{ endif} 1950 1951 \textcolor{comment}{! visc%MLD is used to communicate the state of the (e)PBL or KPP to the rest of the model} 1952 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MLE\_USE\_PBL\_MLD"}, mle\_use\_pbl\_mld, & 1953 default=.false., do\_not\_log=.true.) 1954 \textcolor{comment}{! visc%MLD needs to be allocated when melt potential is computed (HFREEZE>0)} 1955 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"HFREEZE"}, hfreeze, & 1956 default=-1.0, do\_not\_log=.true.) 1957 1958 \textcolor{keywordflow}{if} (mle\_use\_pbl\_mld) \textcolor{keywordflow}{then} 1959 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%MLD, isd, ied, jsd, jed) 1960 \textcolor{keyword}{call }register\_restart\_field(visc%MLD, \textcolor{stringliteral}{"MLD"}, .false., restart\_cs, & 1961 \textcolor{stringliteral}{"Instantaneous active mixing layer depth"}, \textcolor{stringliteral}{"m"}) 1962 \textcolor{keywordflow}{ endif} 1963 1964 \textcolor{keywordflow}{if} (hfreeze >= 0.0 .and. .not.mle\_use\_pbl\_mld) \textcolor{keywordflow}{then} 1965 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%MLD, isd, ied, jsd, jed) 1966 \textcolor{keywordflow}{ endif} 1967 1968 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__set__visc_a9865113fe07928e7a240c2868ed45e5f}\label{namespacemom__set__visc_a9865113fe07928e7a240c2868ed45e5f}} \index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}!set\+\_\+viscous\+\_\+bbl@{set\+\_\+viscous\+\_\+bbl}} \index{set\+\_\+viscous\+\_\+bbl@{set\+\_\+viscous\+\_\+bbl}!mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsubsection{\texorpdfstring{set\+\_\+viscous\+\_\+bbl()}{set\_viscous\_bbl()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+visc\+::set\+\_\+viscous\+\_\+bbl (\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(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{type(vertvisc\+\_\+type), intent(inout)}]{visc, }\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(\hyperlink{structmom__set__visc_1_1set__visc__cs}{set\+\_\+visc\+\_\+cs}), pointer}]{CS, }\item[{logical, intent(in), optional}]{symmetrize }\end{DoxyParamCaption})} Calculates the thickness of the bottom boundary layer and the viscosity within that layer. A drag law is used, either linearized about an assumed bottom velocity or using the actual near-\/bottom velocities combined with an assumed unresolved velocity. The bottom boundary layer thickness is limited by a combination of stratification and rotation, as in the paper of Killworth and Edwards, J\+PO 1999. It is not necessary to calculate the thickness and viscosity every time step; instead previous values may be used.\hypertarget{namespacemom__set__visc_set_viscous_BBL}{}\subsection{Viscous Bottom Boundary Layer}\label{namespacemom__set__visc_set_viscous_BBL} If \hyperlink{structmom__set__visc_1_1set__visc__cs_aca87adb9636e724f7ccb78e7c640f399}{set\+\_\+visc\+\_\+cs.\+bottomdraglaw} is True then a bottom boundary layer viscosity and thickness are calculated so that the bottom stress is \[ \mathbf{\tau}_b = C_d | U_{bbl} | \mathbf{u}_{bbl} \] If \hyperlink{structmom__set__visc_1_1set__visc__cs_aca87adb9636e724f7ccb78e7c640f399}{set\+\_\+visc\+\_\+cs.\+bottomdraglaw} is True then the term $|U_{bbl}|$ is set equal to the value in \hyperlink{structmom__set__visc_1_1set__visc__cs_a8e5f72291c7e5c5e602c9b4765d63ae0}{set\+\_\+visc\+\_\+cs.\+drag\+\_\+bg\+\_\+vel} so that $C_d |U_{bbl}|$ becomes a Rayleigh bottom drag. Otherwise $|U_{bbl}|$ is found by averaging the flow over the bottom \hyperlink{structmom__set__visc_1_1set__visc__cs_a20768c8872b3eab25e59dfe7c5ff62a0}{set\+\_\+visc\+\_\+cs.\+hbbl} of the model, adding the amplitude of tides \hyperlink{structmom__set__visc_1_1set__visc__cs_a765c0eec288ee86e83afc501b757a821}{set\+\_\+visc\+\_\+cs.\+tideamp} and a constant \hyperlink{structmom__set__visc_1_1set__visc__cs_a8e5f72291c7e5c5e602c9b4765d63ae0}{set\+\_\+visc\+\_\+cs.\+drag\+\_\+bg\+\_\+vel}. For these calculations the vertical grid at the velocity component locations is found by \[ \begin{array}{ll} \frac{2 h^- h^+}{h^- + h^+} & u \left( h^+ - h^-\right) >= 0 \\ \frac{1}{2} \left( h^- + h^+ \right) & u \left( h^+ - h^-\right) < 0 \end{array} \] which biases towards the thin cell if the thin cell is upwind. Biasing the grid toward thin upwind cells helps increase the effect of viscosity and inhibits flow out of these thin cells. After diagnosing $|U_{bbl}|$ over a fixed depth an active viscous boundary layer thickness is found using the ideas of Killworth and Edwards, 1999 (hereafter K\+W99). K\+W99 solve the equation \[ \left( \frac{h_{bbl}}{h_f} \right)^2 + \frac{h_{bbl}}{h_N} = 1 \] for the boundary layer depth $h_{bbl}$. Here \[ h_f = \frac{C_n u_*}{f} \] is the rotation controlled boundary layer depth in the absence of stratification. $u_*$ is the surface friction speed given by \[ u_*^2 = C_d |U_{bbl}|^2 \] and is a function of near bottom model flow. \[ h_N = \frac{C_i u_*}{N} = \frac{ (C_i u_* )^2 }{g^\prime} \] is the stratification controlled boundary layer depth. The non-\/dimensional parameters $C_n=0.5$ and $C_i=20$ are suggested by Zilitinkevich and Mironov, 1996. If a Richardson number dependent mixing scheme is being used, as indicated by \hyperlink{structmom__set__visc_1_1set__visc__cs_abf6d0c52b4c7697848e11790fed715dc}{set\+\_\+visc\+\_\+cs.\+rino\+\_\+mix}, then the boundary layer thickness is bounded to be no larger than a half of \hyperlink{structmom__set__visc_1_1set__visc__cs_a20768c8872b3eab25e59dfe7c5ff62a0}{set\+\_\+visc\+\_\+cs.\+hbbl} . \begin{DoxyRefDesc}{Todo} \item[\hyperlink{todo__todo000007}{Todo}]Channel drag needs to be explained\end{DoxyRefDesc} A B\+BL viscosity is calculated so that the no-\/slip boundary condition in the vertical viscosity solver implies the stress $\mathbf{\tau}_b$.\hypertarget{namespacemom__set__visc_set_viscous_BBL_ref}{}\subsubsection{References}\label{namespacemom__set__visc_set_viscous_BBL_ref} \begin{DoxyItemize} \item Killworth, P. D., and N. R. Edwards, 1999\+: A Turbulent Bottom Boundary Layer Code for Use in Numerical Ocean Models. J. Phys. Oceanogr., 29, 1221-\/1238, \href{https://doi.org/10.1175/1520-0485(1999)029<1221:ATBBLC>2.0.CO;2}{\tt doi\+:10.\+1175/1520-\/0485(1999)029$<$1221\+:A\+T\+B\+B\+LC$>$2.\+0.\+CO;2} \item Zilitinkevich, S., Mironov, D.\+V., 1996\+: A multi-\/limit formulation for the equilibrium depth of a stably stratified boundary layer. Boundary-\/\+Layer Meteorology 81, 325-\/351. \href{https://doi.org/10.1007/BF02430334}{\tt doi\+:10.\+1007/\+B\+F02430334}\end{DoxyItemize} \begin{DoxyParams}[1]{Parameters} \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 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} & A structure containing pointers to any available thermodynamic fields. Absent fields have N\+U\+LL ptrs..\\ \hline \mbox{\tt in,out} & {\em visc} & A structure containing vertical viscosities and related fields.\\ \hline & {\em cs} & The control structure returned by a previous call to set\+\_\+visc\+\_\+init.\\ \hline \mbox{\tt in} & {\em symmetrize} & If present and true, do extra calculations of those values in visc that would be calculated with symmetric memory. \\ \hline \end{DoxyParams} Definition at line 193 of file M\+O\+M\+\_\+set\+\_\+viscosity.\+F90. \begin{DoxyCode} 193 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(inout)} :: g\textcolor{comment}{ !< The ocean's grid structure.} 194 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< The ocean's vertical grid structure.} 195 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 196 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, & 197 \textcolor{keywordtype}{intent(in)} :: u\textcolor{comment}{ !< The zonal velocity [L T-1 ~> m s-1].} 198 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, & 199 \textcolor{keywordtype}{intent(in)} :: v\textcolor{comment}{ !< The meridional velocity [L T-1 ~> m s-1].} 200 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, & 201 \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2].} 202 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< A structure containing pointers to any} 203 \textcolor{comment}{ !! available thermodynamic fields. Absent fields} 204 \textcolor{comment}{ !! have NULL ptrs..} 205 \textcolor{keywordtype}{type}(vertvisc\_type), \textcolor{keywordtype}{intent(inout)} :: visc\textcolor{comment}{ !< A structure containing vertical viscosities and} 206 \textcolor{comment}{ !! related fields.} 207 \textcolor{keywordtype}{type}(set\_visc\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< The control structure returned by a previous} 208 \textcolor{comment}{ !! call to set\_visc\_init.} 209 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: symmetrize\textcolor{comment}{ !< If present and true, do extra calculations} 210 \textcolor{comment}{ !! of those values in visc that would be} 211 \textcolor{comment}{ !! calculated with symmetric memory.} 212 213 \textcolor{comment}{! Local variables} 214 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G))} :: & 215 ustar, & \textcolor{comment}{! The bottom friction velocity [Z T-1 ~> m s-1].} 216 t\_eos, & \textcolor{comment}{! The temperature used to calculate the partial derivatives} 217 \textcolor{comment}{! of density with T and S [degC].} 218 s\_eos, & \textcolor{comment}{! The salinity used to calculate the partial derivatives} 219 \textcolor{comment}{! of density with T and S [ppt].} 220 dr\_dt, & \textcolor{comment}{! Partial derivative of the density in the bottom boundary} 221 \textcolor{comment}{! layer with temperature [R degC-1 ~> kg m-3 degC-1].} 222 dr\_ds, & \textcolor{comment}{! Partial derivative of the density in the bottom boundary} 223 \textcolor{comment}{! layer with salinity [R ppt-1 ~> kg m-3 ppt-1].} 224 press \textcolor{comment}{! The pressure at which dR\_dT and dR\_dS are evaluated [R L2 T-2 ~> Pa].} 225 \textcolor{keywordtype}{real} :: htot \textcolor{comment}{! Sum of the layer thicknesses up to some point [H ~> m or kg m-2].} 226 \textcolor{keywordtype}{real} :: htot\_vel \textcolor{comment}{! Sum of the layer thicknesses up to some point [H ~> m or kg m-2].} 227 228 \textcolor{keywordtype}{real} :: rhtot \textcolor{comment}{! Running sum of thicknesses times the layer potential} 229 \textcolor{comment}{! densities [H R ~> kg m-2 or kg2 m-5].} 230 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G))} :: & 231 d\_u, & \textcolor{comment}{! Bottom depth interpolated to u points [Z ~> m].} 232 mask\_u \textcolor{comment}{! A mask that disables any contributions from u points that} 233 \textcolor{comment}{! are land or past open boundary conditions [nondim], 0 or 1.} 234 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G))} :: & 235 d\_v, & \textcolor{comment}{! Bottom depth interpolated to v points [Z ~> m].} 236 mask\_v \textcolor{comment}{! A mask that disables any contributions from v points that} 237 \textcolor{comment}{! are land or past open boundary conditions [nondim], 0 or 1.} 238 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZK\_(G))} :: & 239 h\_at\_vel, & \textcolor{comment}{! Layer thickness at a velocity point, using an upwind-biased} 240 \textcolor{comment}{! second order accurate estimate based on the previous velocity} 241 \textcolor{comment}{! direction [H ~> m or kg m-2].} 242 h\_vel, & \textcolor{comment}{! Arithmetic mean of the layer thicknesses adjacent to a} 243 \textcolor{comment}{! velocity point [H ~> m or kg m-2].} 244 t\_vel, & \textcolor{comment}{! Arithmetic mean of the layer temperatures adjacent to a} 245 \textcolor{comment}{! velocity point [degC].} 246 s\_vel, & \textcolor{comment}{! Arithmetic mean of the layer salinities adjacent to a} 247 \textcolor{comment}{! velocity point [ppt].} 248 rml\_vel \textcolor{comment}{! Arithmetic mean of the layer coordinate densities adjacent} 249 \textcolor{comment}{! to a velocity point [R ~> kg m-3].} 250 251 \textcolor{keywordtype}{real} :: h\_vel\_pos \textcolor{comment}{! The arithmetic mean thickness at a velocity point} 252 \textcolor{comment}{! plus H\_neglect to avoid 0 values [H ~> m or kg m-2].} 253 \textcolor{keywordtype}{real} :: ustarsq \textcolor{comment}{! 400 times the square of ustar, times} 254 \textcolor{comment}{! Rho0 divided by G\_Earth and the conversion} 255 \textcolor{comment}{! from m to thickness units [H R ~> kg m-2 or kg2 m-5].} 256 \textcolor{keywordtype}{real} :: cdrag\_sqrt\_z \textcolor{comment}{! Square root of the drag coefficient, times a unit conversion} 257 \textcolor{comment}{! factor from lateral lengths to vertical depths [Z L-1 ~> 1].} 258 \textcolor{keywordtype}{real} :: cdrag\_sqrt \textcolor{comment}{! Square root of the drag coefficient [nondim].} 259 \textcolor{keywordtype}{real} :: oldfn \textcolor{comment}{! The integrated energy required to} 260 \textcolor{comment}{! entrain up to the bottom of the layer,} 261 \textcolor{comment}{! divided by G\_Earth [H R ~> kg m-2 or kg2 m-5].} 262 \textcolor{keywordtype}{real} :: dfn \textcolor{comment}{! The increment in oldfn for entraining} 263 \textcolor{comment}{! the layer [H R ~> kg m-2 or kg2 m-5].} 264 \textcolor{keywordtype}{real} :: dh \textcolor{comment}{! The increment in layer thickness from} 265 \textcolor{comment}{! the present layer [H ~> m or kg m-2].} 266 \textcolor{keywordtype}{real} :: bbl\_thick \textcolor{comment}{! The thickness of the bottom boundary layer [H ~> m or kg m-2].} 267 \textcolor{keywordtype}{real} :: bbl\_thick\_z \textcolor{comment}{! The thickness of the bottom boundary layer [Z ~> m].} 268 \textcolor{keywordtype}{real} :: kv\_bbl \textcolor{comment}{! The bottom boundary layer viscosity [Z2 T-1 ~> m2 s-1].} 269 \textcolor{keywordtype}{real} :: c2f \textcolor{comment}{! C2f = 2*f at velocity points [T-1 ~> s-1].} 270 271 \textcolor{keywordtype}{real} :: u\_bg\_sq \textcolor{comment}{! The square of an assumed background} 272 \textcolor{comment}{! velocity, for calculating the mean} 273 \textcolor{comment}{! magnitude near the bottom for use in the} 274 \textcolor{comment}{! quadratic bottom drag [L2 T-2 ~> m2 s-2].} 275 \textcolor{keywordtype}{real} :: hwtot \textcolor{comment}{! Sum of the thicknesses used to calculate} 276 \textcolor{comment}{! the near-bottom velocity magnitude [H ~> m or kg m-2].} 277 \textcolor{keywordtype}{real} :: hutot \textcolor{comment}{! Running sum of thicknesses times the} 278 \textcolor{comment}{! velocity magnitudes [H T T-1 ~> m2 s-1 or kg m-1 s-1].} 279 \textcolor{keywordtype}{real} :: thtot \textcolor{comment}{! Running sum of thickness times temperature [degC H ~> degC m or degC kg m-2].} 280 \textcolor{keywordtype}{real} :: shtot \textcolor{comment}{! Running sum of thickness times salinity [ppt H ~> ppt m or ppt kg m-2].} 281 \textcolor{keywordtype}{real} :: hweight \textcolor{comment}{! The thickness of a layer that is within Hbbl} 282 \textcolor{comment}{! of the bottom [H ~> m or kg m-2].} 283 \textcolor{keywordtype}{real} :: v\_at\_u, u\_at\_v \textcolor{comment}{! v at a u point or vice versa [L T-1 ~> m s-1].} 284 \textcolor{keywordtype}{real} :: rho0x400\_g \textcolor{comment}{! 400*Rho0/G\_Earth, times unit conversion factors} 285 \textcolor{comment}{! [R T2 H Z-2 ~> kg s2 m-4 or kg2 s2 m-7].} 286 \textcolor{comment}{! The 400 is a constant proposed by Killworth and Edwards, 1999.} 287 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),max(GV%nk\_rho\_varies,1))} :: & 288 rml \textcolor{comment}{! The mixed layer coordinate density [R ~> kg m-3].} 289 \textcolor{keywordtype}{real} :: p\_ref(szi\_(g)) \textcolor{comment}{! The pressure used to calculate the coordinate} 290 \textcolor{comment}{! density [R L2 T-2 ~> Pa] (usually set to 2e7 Pa = 2000 dbar).} 291 292 \textcolor{keywordtype}{real} :: d\_vel \textcolor{comment}{! The bottom depth at a velocity point [H ~> m or kg m-2].} 293 \textcolor{keywordtype}{real} :: dp, dm \textcolor{comment}{! The depths at the edges of a velocity cell [H ~> m or kg m-2].} 294 \textcolor{keywordtype}{real} :: a \textcolor{comment}{! a is the curvature of the bottom depth across a} 295 \textcolor{comment}{! cell, times the cell width squared [H ~> m or kg m-2].} 296 \textcolor{keywordtype}{real} :: a\_3, a\_12 \textcolor{comment}{! a/3 and a/12 [H ~> m or kg m-2].} 297 \textcolor{keywordtype}{real} :: c24\_a \textcolor{comment}{! 24/a [H-1 ~> m-1 or m2 kg-1].} 298 \textcolor{keywordtype}{real} :: slope \textcolor{comment}{! The absolute value of the bottom depth slope across} 299 \textcolor{comment}{! a cell times the cell width [H ~> m or kg m-2].} 300 \textcolor{keywordtype}{real} :: apb\_4a, ax2\_3apb \textcolor{comment}{! Various nondimensional ratios of a and slope.} 301 \textcolor{keywordtype}{real} :: a2x48\_apb3, iapb, ibma\_2 \textcolor{comment}{! Combinations of a and slope [H-1 ~> m-1 or m2 kg-1].} 302 \textcolor{comment}{! All of the following "volumes" have units of thickness because they are normalized} 303 \textcolor{comment}{! by the full horizontal area of a velocity cell.} 304 \textcolor{keywordtype}{real} :: vol\_open \textcolor{comment}{! The cell volume above which it is open [H ~> m or kg m-2].} 305 \textcolor{keywordtype}{real} :: vol\_direct \textcolor{comment}{! With less than Vol\_direct [H ~> m or kg m-2], there is a direct} 306 \textcolor{comment}{! solution of a cubic equation for L.} 307 \textcolor{keywordtype}{real} :: vol\_2\_reg \textcolor{comment}{! The cell volume above which there are two separate} 308 \textcolor{comment}{! open areas that must be integrated [H ~> m or kg m-2].} 309 \textcolor{keywordtype}{real} :: vol \textcolor{comment}{! The volume below the interface whose normalized} 310 \textcolor{comment}{! width is being sought [H ~> m or kg m-2].} 311 \textcolor{keywordtype}{real} :: vol\_below \textcolor{comment}{! The volume below the interface below the one that} 312 \textcolor{comment}{! is currently under consideration [H ~> m or kg m-2].} 313 \textcolor{keywordtype}{real} :: vol\_err \textcolor{comment}{! The error in the volume with the latest estimate of} 314 \textcolor{comment}{! L, or the error for the interface below [H ~> m or kg m-2].} 315 \textcolor{keywordtype}{real} :: vol\_quit \textcolor{comment}{! The volume error below which to quit iterating [H ~> m or kg m-2].} 316 \textcolor{keywordtype}{real} :: vol\_tol \textcolor{comment}{! A volume error tolerance [H ~> m or kg m-2].} 317 \textcolor{keywordtype}{real} :: l(szk\_(g)+1) \textcolor{comment}{! The fraction of the full cell width that is open at} 318 \textcolor{comment}{! the depth of each interface [nondim].} 319 \textcolor{keywordtype}{real} :: l\_direct \textcolor{comment}{! The value of L above volume Vol\_direct [nondim].} 320 \textcolor{keywordtype}{real} :: l\_max, l\_min \textcolor{comment}{! Upper and lower bounds on the correct value for L [nondim].} 321 \textcolor{keywordtype}{real} :: vol\_err\_max \textcolor{comment}{! The volume errors for the upper and lower bounds on} 322 \textcolor{keywordtype}{real} :: vol\_err\_min \textcolor{comment}{! the correct value for L [H ~> m or kg m-2].} 323 \textcolor{keywordtype}{real} :: vol\_0 \textcolor{comment}{! A deeper volume with known width L0 [H ~> m or kg m-2].} 324 \textcolor{keywordtype}{real} :: l0 \textcolor{comment}{! The value of L above volume Vol\_0 [nondim].} 325 \textcolor{keywordtype}{real} :: dvol \textcolor{comment}{! vol - Vol\_0 [H ~> m or kg m-2].} 326 \textcolor{keywordtype}{real} :: dv\_dl2 \textcolor{comment}{! The partial derivative of volume with L squared} 327 \textcolor{comment}{! evaluated at L=L0 [H ~> m or kg m-2].} 328 \textcolor{keywordtype}{real} :: h\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 329 \textcolor{comment}{! in roundoff and can be neglected [H ~> m or kg m-2].} 330 \textcolor{keywordtype}{real} :: usth \textcolor{comment}{! ustar converted to units of H T-1 [H T-1 ~> m s-1 or kg m-2 s-1].} 331 \textcolor{keywordtype}{real} :: root \textcolor{comment}{! A temporary variable [H T-1 ~> m s-1 or kg m-2 s-1].} 332 333 \textcolor{keywordtype}{real} :: cell\_width \textcolor{comment}{! The transverse width of the velocity cell [L ~> m].} 334 \textcolor{keywordtype}{real} :: rayleigh \textcolor{comment}{! A nondimensional value that is multiplied by the layer's} 335 \textcolor{comment}{! velocity magnitude to give the Rayleigh drag velocity, times} 336 \textcolor{comment}{! a lateral to vertical distance conversion factor [Z L-1 ~> 1].} 337 \textcolor{keywordtype}{real} :: gam \textcolor{comment}{! The ratio of the change in the open interface width} 338 \textcolor{comment}{! to the open interface width atop a cell [nondim].} 339 \textcolor{keywordtype}{real} :: bbl\_frac \textcolor{comment}{! The fraction of a layer's drag that goes into the} 340 \textcolor{comment}{! viscous bottom boundary layer [nondim].} 341 \textcolor{keywordtype}{real} :: bbl\_visc\_frac \textcolor{comment}{! The fraction of all the drag that is expressed as} 342 \textcolor{comment}{! a viscous bottom boundary layer [nondim].} 343 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: c1\_3 = 1.0/3.0, c1\_6 = 1.0/6.0, c1\_12 = 1.0/12.0 344 \textcolor{keywordtype}{real} :: c2pi\_3 \textcolor{comment}{! An irrational constant, 2/3 pi.} 345 \textcolor{keywordtype}{real} :: tmp \textcolor{comment}{! A temporary variable.} 346 \textcolor{keywordtype}{real} :: tmp\_val\_m1\_to\_p1 347 \textcolor{keywordtype}{real} :: curv\_tol \textcolor{comment}{! Numerator of curvature cubed, used to estimate} 348 \textcolor{comment}{! accuracy of a single L(:) Newton iteration} 349 \textcolor{keywordtype}{logical} :: use\_l0, do\_one\_l\_iter \textcolor{comment}{! Control flags for L(:) Newton iteration} 350 \textcolor{keywordtype}{logical} :: use\_bbl\_eos, do\_i(szib\_(g)) 351 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The computational domain for the equation of state} 352 \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, isq, ieq, jsq, jeq, nz, m, n, k2, nkmb, nkml 353 \textcolor{keywordtype}{integer} :: itt, maxitt=20 354 \textcolor{keywordtype}{type}(ocean\_obc\_type), \textcolor{keywordtype}{pointer} :: obc => null() 355 356 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke 357 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB 358 nkmb = gv%nk\_rho\_varies ; nkml = gv%nkml 359 h\_neglect = gv%H\_subroundoff 360 rho0x400\_g = 400.0*(gv%Rho0 / (us%L\_to\_Z**2 * gv%g\_Earth)) * gv%Z\_to\_H 361 vol\_quit = 0.9*gv%Angstrom\_H + h\_neglect 362 c2pi\_3 = 8.0*atan(1.0)/3.0 363 364 \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal,\textcolor{stringliteral}{"MOM\_set\_viscosity(BBL): "}//& 365 \textcolor{stringliteral}{"Module must be initialized before it is used."}) 366 \textcolor{keywordflow}{if} (.not.cs%bottomdraglaw) \textcolor{keywordflow}{return} 367 368 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(symmetrize)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (symmetrize) \textcolor{keywordflow}{then} 369 jsq = js-1 ; isq = is-1 370 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 371 372 \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then} 373 \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"Start set\_viscous\_BBL [uv]"}, u, v, g%HI, haloshift=1, scale=us%L\_T\_to\_m\_s) 374 \textcolor{keyword}{call }hchksum(h,\textcolor{stringliteral}{"Start set\_viscous\_BBL h"}, g%HI, haloshift=1, scale=gv%H\_to\_m) 375 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%T)) \textcolor{keyword}{call }hchksum(tv%T, \textcolor{stringliteral}{"Start set\_viscous\_BBL T"}, g%HI, haloshift=1) 376 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%S)) \textcolor{keyword}{call }hchksum(tv%S, \textcolor{stringliteral}{"Start set\_viscous\_BBL S"}, g%HI, haloshift=1) 377 \textcolor{keywordflow}{ endif} 378 379 use\_bbl\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) .and. cs%BBL\_use\_EOS 380 obc => cs%OBC 381 382 u\_bg\_sq = cs%drag\_bg\_vel * cs%drag\_bg\_vel 383 cdrag\_sqrt = sqrt(cs%cdrag) 384 cdrag\_sqrt\_z = us%L\_to\_Z * sqrt(cs%cdrag) 385 k2 = max(nkmb+1, 2) 386 387 \textcolor{comment}{! With a linear drag law, the friction velocity is already known.} 388 \textcolor{comment}{! if (CS%linear\_drag) ustar(:) = cdrag\_sqrt\_Z*CS%drag\_bg\_vel} 389 390 \textcolor{keywordflow}{if} ((nkml>0) .and. .not.use\_bbl\_eos) \textcolor{keywordflow}{then} 391 eosdom(1) = isq - (g%isd-1) ; eosdom(2) = g%iec+1 - (g%isd-1) 392 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p\_ref(i) = tv%P\_Ref ;\textcolor{keywordflow}{ enddo} 393 \textcolor{comment}{!$OMP parallel do default(shared)} 394 \textcolor{keywordflow}{do} k=1,nkmb ; \textcolor{keywordflow}{do} j=jsq,jeq+1 395 \textcolor{keyword}{call }calculate\_density(tv%T(:,j,k), tv%S(:,j,k), p\_ref, rml(:,j,k), tv%eqn\_of\_state, & 396 eosdom) 397 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 398 \textcolor{keywordflow}{ endif} 399 400 \textcolor{comment}{!$OMP parallel do default(shared)} 401 \textcolor{keywordflow}{do} j=js-1,je ; \textcolor{keywordflow}{do} i=is-1,ie+1 402 d\_v(i,j) = 0.5*(g%bathyT(i,j) + g%bathyT(i,j+1)) 403 mask\_v(i,j) = g%mask2dCv(i,j) 404 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 405 \textcolor{comment}{!$OMP parallel do default(shared)} 406 \textcolor{keywordflow}{do} j=js-1,je+1 ; \textcolor{keywordflow}{do} i=is-1,ie 407 d\_u(i,j) = 0.5*(g%bathyT(i,j) + g%bathyT(i+1,j)) 408 mask\_u(i,j) = g%mask2dCu(i,j) 409 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 410 411 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} n=1,obc%number\_of\_segments 412 \textcolor{keywordflow}{if} (.not. obc%segment(n)%on\_pe) cycle 413 \textcolor{comment}{! Use a one-sided projection of bottom depths at OBC points.} 414 i = obc%segment(n)%HI%IsdB ; j = obc%segment(n)%HI%JsdB 415 \textcolor{keywordflow}{if} (obc%segment(n)%is\_N\_or\_S .and. (j >= js-1) .and. (j <= je)) \textcolor{keywordflow}{then} 416 \textcolor{keywordflow}{do} i = max(is-1,obc%segment(n)%HI%isd), min(ie+1,obc%segment(n)%HI%ied) 417 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_n) d\_v(i,j) = g%bathyT(i,j) 418 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_s) d\_v(i,j) = g%bathyT(i,j+1) 419 \textcolor{keywordflow}{ enddo} 420 \textcolor{keywordflow}{elseif} (obc%segment(n)%is\_E\_or\_W .and. (i >= is-1) .and. (i <= ie)) \textcolor{keywordflow}{then} 421 \textcolor{keywordflow}{do} j = max(js-1,obc%segment(n)%HI%jsd), min(je+1,obc%segment(n)%HI%jed) 422 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_e) d\_u(i,j) = g%bathyT(i,j) 423 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_w) d\_u(i,j) = g%bathyT(i+1,j) 424 \textcolor{keywordflow}{ enddo} 425 \textcolor{keywordflow}{ endif} 426 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 427 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} n=1,obc%number\_of\_segments 428 \textcolor{comment}{! Now project bottom depths across cell-corner points in the OBCs. The two} 429 \textcolor{comment}{! projections have to occur in sequence and can not be combined easily.} 430 \textcolor{keywordflow}{if} (.not. obc%segment(n)%on\_pe) cycle 431 \textcolor{comment}{! Use a one-sided projection of bottom depths at OBC points.} 432 i = obc%segment(n)%HI%IsdB ; j = obc%segment(n)%HI%JsdB 433 \textcolor{keywordflow}{if} (obc%segment(n)%is\_N\_or\_S .and. (j >= js-1) .and. (j <= je)) \textcolor{keywordflow}{then} 434 \textcolor{keywordflow}{do} i = max(is-1,obc%segment(n)%HI%IsdB), min(ie,obc%segment(n)%HI%IedB) 435 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_n) \textcolor{keywordflow}{then} 436 d\_u(i,j+1) = d\_u(i,j) ; mask\_u(i,j+1) = 0.0 437 \textcolor{keywordflow}{elseif} (obc%segment(n)%direction == obc\_direction\_s) \textcolor{keywordflow}{then} 438 d\_u(i,j) = d\_u(i,j+1) ; mask\_u(i,j) = 0.0 439 \textcolor{keywordflow}{ endif} 440 \textcolor{keywordflow}{ enddo} 441 \textcolor{keywordflow}{elseif} (obc%segment(n)%is\_E\_or\_W .and. (i >= is-1) .and. (i <= ie)) \textcolor{keywordflow}{then} 442 \textcolor{keywordflow}{do} j = max(js-1,obc%segment(n)%HI%JsdB), min(je,obc%segment(n)%HI%JedB) 443 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_e) \textcolor{keywordflow}{then} 444 d\_v(i+1,j) = d\_v(i,j) ; mask\_v(i+1,j) = 0.0 445 \textcolor{keywordflow}{elseif} (obc%segment(n)%direction == obc\_direction\_w) \textcolor{keywordflow}{then} 446 d\_v(i,j) = d\_v(i+1,j) ; mask\_v(i,j) = 0.0 447 \textcolor{keywordflow}{ endif} 448 \textcolor{keywordflow}{ enddo} 449 \textcolor{keywordflow}{ endif} 450 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 451 452 \textcolor{keywordflow}{if} (.not.use\_bbl\_eos) rml\_vel(:,:) = 0.0 453 454 \textcolor{comment}{!$OMP parallel do default(private) shared(u,v,h,tv,visc,G,GV,US,CS,Rml,nz,nkmb, &} 455 \textcolor{comment}{!$OMP nkml,Isq,Ieq,Jsq,Jeq,h\_neglect,Rho0x400\_G,C2pi\_3, &} 456 \textcolor{comment}{!$OMP U\_bg\_sq,cdrag\_sqrt\_Z,cdrag\_sqrt,K2,use\_BBL\_EOS, &} 457 \textcolor{comment}{!$OMP OBC,maxitt,D\_u,D\_v,mask\_u,mask\_v) &} 458 \textcolor{comment}{!$OMP firstprivate(Vol\_quit)} 459 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} m=1,2 460 461 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} 462 \textcolor{comment}{! m=1 refers to u-points} 463 \textcolor{keywordflow}{if} (j 0) do\_i(i) = .true. 468 \textcolor{keywordflow}{ enddo} 469 \textcolor{keywordflow}{else} 470 \textcolor{comment}{! m=2 refers to v-points} 471 is = g%isc ; ie = g%iec 472 \textcolor{keywordflow}{do} i=is,ie 473 do\_i(i) = .false. 474 \textcolor{keywordflow}{if} (g%mask2dCv(i,j) > 0) do\_i(i) = .true. 475 \textcolor{keywordflow}{ enddo} 476 \textcolor{keywordflow}{ endif} 477 478 \textcolor{comment}{! Calculate thickness at velocity points (u or v depending on value of m).} 479 \textcolor{comment}{! Also interpolate the ML density or T/S properties.} 480 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} \textcolor{comment}{! u-points} 481 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 482 \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 483 \textcolor{keywordflow}{if} (u(i,j,k) *(h(i+1,j,k) - h(i,j,k)) >= 0) \textcolor{keywordflow}{then} 484 \textcolor{comment}{! If the flow is from thin to thick then bias towards the thinner thickness} 485 h\_at\_vel(i,k) = 2.0*h(i,j,k)*h(i+1,j,k) / & 486 (h(i,j,k) + h(i+1,j,k) + h\_neglect) 487 \textcolor{keywordflow}{else} 488 \textcolor{comment}{! If the flow is from thick to thin then use the simple average thickness} 489 h\_at\_vel(i,k) = 0.5 * (h(i,j,k) + h(i+1,j,k)) 490 \textcolor{keywordflow}{ endif} 491 \textcolor{keywordflow}{ endif} 492 h\_vel(i,k) = 0.5 * (h(i,j,k) + h(i+1,j,k)) 493 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 494 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 495 \textcolor{comment}{! Perhaps these should be thickness weighted.} 496 t\_vel(i,k) = 0.5 * (tv%T(i,j,k) + tv%T(i+1,j,k)) 497 s\_vel(i,k) = 0.5 * (tv%S(i,j,k) + tv%S(i+1,j,k)) 498 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ; \textcolor{keywordflow}{else} ; \textcolor{keywordflow}{do} k=1,nkmb ; \textcolor{keywordflow}{do} i=is,ie 499 rml\_vel(i,k) = 0.5 * (rml(i,j,k) + rml(i+1,j,k)) 500 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 501 \textcolor{keywordflow}{else} \textcolor{comment}{! v-points} 502 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 503 \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 504 \textcolor{keywordflow}{if} (v(i,j,k) * (h(i,j+1,k) - h(i,j,k)) >= 0) \textcolor{keywordflow}{then} 505 \textcolor{comment}{! If the flow is from thin to thick then bias towards the thinner thickness} 506 h\_at\_vel(i,k) = 2.0*h(i,j,k)*h(i,j+1,k) / & 507 (h(i,j,k) + h(i,j+1,k) + h\_neglect) 508 \textcolor{keywordflow}{else} 509 \textcolor{comment}{! If the flow is from thick to thin then use the simple average thickness} 510 h\_at\_vel(i,k) = 0.5 * (h(i,j,k) + h(i,j+1,k)) 511 \textcolor{keywordflow}{ endif} 512 \textcolor{keywordflow}{ endif} 513 h\_vel(i,k) = 0.5 * (h(i,j,k) + h(i,j+1,k)) 514 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 515 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 516 t\_vel(i,k) = 0.5 * (tv%T(i,j,k) + tv%T(i,j+1,k)) 517 s\_vel(i,k) = 0.5 * (tv%S(i,j,k) + tv%S(i,j+1,k)) 518 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ; \textcolor{keywordflow}{else} ; \textcolor{keywordflow}{do} k=1,nkmb ; \textcolor{keywordflow}{do} i=is,ie 519 rml\_vel(i,k) = 0.5 * (rml(i,j,k) + rml(i,j+1,k)) 520 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 521 \textcolor{keywordflow}{ endif} 522 523 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (obc%number\_of\_segments > 0) \textcolor{keywordflow}{then} 524 \textcolor{comment}{! Apply a zero gradient projection of thickness across OBC points.} 525 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} 526 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i) .and. (obc%segnum\_u(i,j) /= obc\_none)) \textcolor{keywordflow}{then} 527 \textcolor{keywordflow}{if} (obc%segment(obc%segnum\_u(i,j))%direction == obc\_direction\_e) \textcolor{keywordflow}{then} 528 \textcolor{keywordflow}{do} k=1,nz 529 h\_at\_vel(i,k) = h(i,j,k) ; h\_vel(i,k) = h(i,j,k) 530 \textcolor{keywordflow}{ enddo} 531 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} 532 \textcolor{keywordflow}{do} k=1,nz 533 t\_vel(i,k) = tv%T(i,j,k) ; s\_vel(i,k) = tv%S(i,j,k) 534 \textcolor{keywordflow}{ enddo} 535 \textcolor{keywordflow}{else} 536 \textcolor{keywordflow}{do} k=1,nkmb 537 rml\_vel(i,k) = rml(i,j,k) 538 \textcolor{keywordflow}{ enddo} 539 \textcolor{keywordflow}{ endif} 540 \textcolor{keywordflow}{elseif} (obc%segment(obc%segnum\_u(i,j))%direction == obc\_direction\_w) \textcolor{keywordflow}{then} 541 \textcolor{keywordflow}{do} k=1,nz 542 h\_at\_vel(i,k) = h(i+1,j,k) ; h\_vel(i,k) = h(i+1,j,k) 543 \textcolor{keywordflow}{ enddo} 544 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} 545 \textcolor{keywordflow}{do} k=1,nz 546 t\_vel(i,k) = tv%T(i+1,j,k) ; s\_vel(i,k) = tv%S(i+1,j,k) 547 \textcolor{keywordflow}{ enddo} 548 \textcolor{keywordflow}{else} 549 \textcolor{keywordflow}{do} k=1,nkmb 550 rml\_vel(i,k) = rml(i+1,j,k) 551 \textcolor{keywordflow}{ enddo} 552 \textcolor{keywordflow}{ endif} 553 \textcolor{keywordflow}{ endif} 554 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 555 \textcolor{keywordflow}{else} 556 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i) .and. (obc%segnum\_v(i,j) /= obc\_none)) \textcolor{keywordflow}{then} 557 \textcolor{keywordflow}{if} (obc%segment(obc%segnum\_v(i,j))%direction == obc\_direction\_n) \textcolor{keywordflow}{then} 558 \textcolor{keywordflow}{do} k=1,nz 559 h\_at\_vel(i,k) = h(i,j,k) ; h\_vel(i,k) = h(i,j,k) 560 \textcolor{keywordflow}{ enddo} 561 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} 562 \textcolor{keywordflow}{do} k=1,nz 563 t\_vel(i,k) = tv%T(i,j,k) ; s\_vel(i,k) = tv%S(i,j,k) 564 \textcolor{keywordflow}{ enddo} 565 \textcolor{keywordflow}{else} 566 \textcolor{keywordflow}{do} k=1,nkmb 567 rml\_vel(i,k) = rml(i,j,k) 568 \textcolor{keywordflow}{ enddo} 569 \textcolor{keywordflow}{ endif} 570 \textcolor{keywordflow}{elseif} (obc%segment(obc%segnum\_v(i,j))%direction == obc\_direction\_s) \textcolor{keywordflow}{then} 571 \textcolor{keywordflow}{do} k=1,nz 572 h\_at\_vel(i,k) = h(i,j+1,k) ; h\_vel(i,k) = h(i,j+1,k) 573 \textcolor{keywordflow}{ enddo} 574 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} 575 \textcolor{keywordflow}{do} k=1,nz 576 t\_vel(i,k) = tv%T(i,j+1,k) ; s\_vel(i,k) = tv%S(i,j+1,k) 577 \textcolor{keywordflow}{ enddo} 578 \textcolor{keywordflow}{else} 579 \textcolor{keywordflow}{do} k=1,nkmb 580 rml\_vel(i,k) = rml(i,j+1,k) 581 \textcolor{keywordflow}{ enddo} 582 \textcolor{keywordflow}{ endif} 583 \textcolor{keywordflow}{ endif} 584 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 585 \textcolor{keywordflow}{ endif} 586 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 587 588 \textcolor{keywordflow}{if} (use\_bbl\_eos .or. .not.cs%linear\_drag) \textcolor{keywordflow}{then} 589 \textcolor{comment}{! Calculate the mean velocity magnitude over the bottommost CS%Hbbl of} 590 \textcolor{comment}{! the water column for determining the quadratic bottom drag.} 591 \textcolor{comment}{! Used in ustar(i)} 592 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 593 htot\_vel = 0.0 ; hwtot = 0.0 ; hutot = 0.0 594 thtot = 0.0 ; shtot = 0.0 595 \textcolor{keywordflow}{do} k=nz,1,-1 596 597 \textcolor{keywordflow}{if} (htot\_vel>=cs%Hbbl) \textcolor{keywordflow}{exit} \textcolor{comment}{! terminate the k loop} 598 599 hweight = min(cs%Hbbl - htot\_vel, h\_at\_vel(i,k)) 600 \textcolor{keywordflow}{if} (hweight < 1.5*gv%Angstrom\_H + h\_neglect) cycle 601 602 htot\_vel = htot\_vel + h\_at\_vel(i,k) 603 hwtot = hwtot + hweight 604 605 \textcolor{keywordflow}{if} ((.not.cs%linear\_drag) .and. (hweight >= 0.0)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} 606 v\_at\_u = set\_v\_at\_u(v, h, g, i, j, k, mask\_v, obc) 607 \textcolor{keywordflow}{if} (cs%BBL\_use\_tidal\_bg) \textcolor{keywordflow}{then} 608 u\_bg\_sq = 0.5*( g%mask2dT(i,j)*(cs%tideamp(i,j)*cs%tideamp(i,j))+ & 609 g%mask2dT(i+1,j)*(cs%tideamp(i+1,j)*cs%tideamp(i+1,j)) ) 610 \textcolor{keywordflow}{ endif} 611 hutot = hutot + hweight * sqrt(u(i,j,k)*u(i,j,k) + & 612 v\_at\_u*v\_at\_u + u\_bg\_sq) 613 \textcolor{keywordflow}{else} 614 u\_at\_v = set\_u\_at\_v(u, h, g, i, j, k, mask\_u, obc) 615 \textcolor{keywordflow}{if} (cs%BBL\_use\_tidal\_bg) \textcolor{keywordflow}{then} 616 u\_bg\_sq = 0.5*( g%mask2dT(i,j)*(cs%tideamp(i,j)*cs%tideamp(i,j))+ & 617 g%mask2dT(i,j+1)*(cs%tideamp(i,j+1)*cs%tideamp(i,j+1)) ) 618 \textcolor{keywordflow}{ endif} 619 hutot = hutot + hweight * sqrt(v(i,j,k)*v(i,j,k) + & 620 u\_at\_v*u\_at\_v + u\_bg\_sq) 621 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 622 623 \textcolor{keywordflow}{if} (use\_bbl\_eos .and. (hweight >= 0.0)) \textcolor{keywordflow}{then} 624 thtot = thtot + hweight * t\_vel(i,k) 625 shtot = shtot + hweight * s\_vel(i,k) 626 \textcolor{keywordflow}{ endif} 627 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of k loop} 628 629 \textcolor{comment}{! Set u* based on u*^2 = Cdrag u\_bbl^2} 630 \textcolor{keywordflow}{if} (.not.cs%linear\_drag .and. (hwtot > 0.0)) \textcolor{keywordflow}{then} 631 ustar(i) = cdrag\_sqrt\_z*hutot/hwtot 632 \textcolor{keywordflow}{else} 633 ustar(i) = cdrag\_sqrt\_z*cs%drag\_bg\_vel 634 \textcolor{keywordflow}{ endif} 635 636 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (hwtot > 0.0) \textcolor{keywordflow}{then} 637 t\_eos(i) = thtot/hwtot ; s\_eos(i) = shtot/hwtot 638 \textcolor{keywordflow}{else} 639 t\_eos(i) = 0.0 ; s\_eos(i) = 0.0 640 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 641 642 \textcolor{comment}{! Diagnostic BBL flow speed at u- and v-points.} 643 \textcolor{keywordflow}{if} (cs%id\_bbl\_u>0 .and. m==1) \textcolor{keywordflow}{then} 644 \textcolor{keywordflow}{if} (hwtot > 0.0) cs%bbl\_u(i,j) = hutot/hwtot 645 \textcolor{keywordflow}{elseif} (cs%id\_bbl\_v>0 .and. m==2) \textcolor{keywordflow}{then} 646 \textcolor{keywordflow}{if} (hwtot > 0.0) cs%bbl\_v(i,j) = hutot/hwtot 647 \textcolor{keywordflow}{ endif} 648 649 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 650 \textcolor{keywordflow}{else} 651 \textcolor{keywordflow}{do} i=is,ie ; ustar(i) = cdrag\_sqrt\_z*cs%drag\_bg\_vel ;\textcolor{keywordflow}{ enddo} 652 \textcolor{keywordflow}{ endif} \textcolor{comment}{! Not linear\_drag} 653 654 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} 655 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%p\_surf)) \textcolor{keywordflow}{then} 656 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} i=is,ie ; press(i) = 0.5*(tv%p\_surf(i,j) + tv%p\_surf(i+1,j)) ;\textcolor{keywordflow}{ enddo} 657 \textcolor{keywordflow}{else} ; \textcolor{keywordflow}{do} i=is,ie ; press(i) = 0.5*(tv%p\_surf(i,j) + tv%p\_surf(i,j+1)) ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 658 \textcolor{keywordflow}{else} 659 \textcolor{keywordflow}{do} i=is,ie ; press(i) = 0.0 ;\textcolor{keywordflow}{ enddo} 660 \textcolor{keywordflow}{ endif} 661 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (.not.do\_i(i)) \textcolor{keywordflow}{then} ; t\_eos(i) = 0.0 ; s\_eos(i) = 0.0 ;\textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 662 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 663 press(i) = press(i) + (gv%H\_to\_RZ*gv%g\_Earth) * h\_vel(i,k) 664 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 665 \textcolor{keyword}{call }calculate\_density\_derivs(t\_eos, s\_eos, press, dr\_dt, dr\_ds, tv%eqn\_of\_state, & 666 (/is-g%IsdB+1,ie-g%IsdB+1/) ) 667 \textcolor{keywordflow}{ endif} 668 669 \textcolor{comment}{! Find a BBL thickness given by equation 2.20 of Killworth and Edwards, 1999:} 670 \textcolor{comment}{! ( f h / Cn u* )^2 + ( N h / Ci u* ) = 1} 671 \textcolor{comment}{! where Cn=0.5 and Ci=20 (constants suggested by Zilitinkevich and Mironov, 1996).} 672 \textcolor{comment}{! Eq. 2.20 can be expressed in terms of boundary layer thicknesses limited by} 673 \textcolor{comment}{! rotation (h\_f) and stratification (h\_N):} 674 \textcolor{comment}{! ( h / h\_f )^2 + ( h / h\_N ) = 1} 675 \textcolor{comment}{! When stratification dominates h\_N<0)} 683 \textcolor{comment}{! For layer mode, N^2 = g'/h. Since (Ci u*)^2 = (h\_N N)^2 = h\_N g' then} 684 \textcolor{comment}{! h\_N = (Ci u*)^2 / g' (KW99, eq, 2.22)} 685 \textcolor{comment}{! Starting from the bottom, integrate the stratification upward until h\_N N balances Ci u*} 686 \textcolor{comment}{! or in layer mode} 687 \textcolor{comment}{! h\_N Delta rho ~ (Ci u*)^2 rho0 / g} 688 \textcolor{comment}{! where the rhs is stored in variable ustarsq.} 689 \textcolor{comment}{! The method was described in Stephens and Hallberg 2000 (unpublished and lost manuscript).} 690 \textcolor{keywordflow}{if} (use\_bbl\_eos) \textcolor{keywordflow}{then} 691 thtot = 0.0 ; shtot = 0.0 ; oldfn = 0.0 692 \textcolor{keywordflow}{do} k=nz,2,-1 693 \textcolor{keywordflow}{if} (h\_at\_vel(i,k) <= 0.0) cycle 694 695 \textcolor{comment}{! Delta rho * h\_bbl assuming everything below is homogenized} 696 oldfn = dr\_dt(i)*(thtot - t\_vel(i,k)*htot) + & 697 dr\_ds(i)*(shtot - s\_vel(i,k)*htot) 698 \textcolor{keywordflow}{if} (oldfn >= ustarsq) \textcolor{keywordflow}{exit} 699 700 \textcolor{comment}{! Local Delta rho * h\_bbl at interface} 701 dfn = (dr\_dt(i)*(t\_vel(i,k) - t\_vel(i,k-1)) + & 702 dr\_ds(i)*(s\_vel(i,k) - s\_vel(i,k-1))) * & 703 (h\_at\_vel(i,k) + htot) 704 705 \textcolor{keywordflow}{if} ((oldfn + dfn) <= ustarsq) \textcolor{keywordflow}{then} 706 \textcolor{comment}{! Use whole layer} 707 dh = h\_at\_vel(i,k) 708 \textcolor{keywordflow}{else} 709 \textcolor{comment}{! Use only part of the layer} 710 dh = h\_at\_vel(i,k) * sqrt((ustarsq-oldfn) / (dfn)) 711 \textcolor{keywordflow}{ endif} 712 713 \textcolor{comment}{! Increment total BBL thickness and cumulative T and S} 714 htot = htot + dh 715 thtot = thtot + t\_vel(i,k)*dh ; shtot = shtot + s\_vel(i,k)*dh 716 \textcolor{keywordflow}{ enddo} 717 \textcolor{keywordflow}{if} ((oldfn < ustarsq) .and. h\_at\_vel(i,1) > 0.0) \textcolor{keywordflow}{then} 718 \textcolor{comment}{! Layer 1 might be part of the BBL.} 719 \textcolor{keywordflow}{if} (dr\_dt(i) * (thtot - t\_vel(i,1)*htot) + & 720 dr\_ds(i) * (shtot - s\_vel(i,1)*htot) < ustarsq) & 721 htot = htot + h\_at\_vel(i,1) 722 \textcolor{keywordflow}{ endif} \textcolor{comment}{! Examination of layer 1.} 723 \textcolor{keywordflow}{else} \textcolor{comment}{! Use Rlay and/or the coordinate density as density variables.} 724 rhtot = 0.0 725 \textcolor{keywordflow}{do} k=nz,k2,-1 726 oldfn = rhtot - gv%Rlay(k)*htot 727 dfn = (gv%Rlay(k) - gv%Rlay(k-1))*(h\_at\_vel(i,k)+htot) 728 729 \textcolor{keywordflow}{if} (oldfn >= ustarsq) \textcolor{keywordflow}{then} 730 cycle 731 \textcolor{keywordflow}{elseif} ((oldfn + dfn) <= ustarsq) \textcolor{keywordflow}{then} 732 dh = h\_at\_vel(i,k) 733 \textcolor{keywordflow}{else} 734 dh = h\_at\_vel(i,k) * sqrt((ustarsq-oldfn) / (dfn)) 735 \textcolor{keywordflow}{ endif} 736 737 htot = htot + dh 738 rhtot = rhtot + gv%Rlay(k)*dh 739 \textcolor{keywordflow}{ enddo} 740 \textcolor{keywordflow}{if} (nkml>0) \textcolor{keywordflow}{then} 741 \textcolor{keywordflow}{do} k=nkmb,2,-1 742 oldfn = rhtot - rml\_vel(i,k)*htot 743 dfn = (rml\_vel(i,k) - rml\_vel(i,k-1)) * (h\_at\_vel(i,k)+htot) 744 745 \textcolor{keywordflow}{if} (oldfn >= ustarsq) \textcolor{keywordflow}{then} 746 cycle 747 \textcolor{keywordflow}{elseif} ((oldfn + dfn) <= ustarsq) \textcolor{keywordflow}{then} 748 dh = h\_at\_vel(i,k) 749 \textcolor{keywordflow}{else} 750 dh = h\_at\_vel(i,k) * sqrt((ustarsq-oldfn) / (dfn)) 751 \textcolor{keywordflow}{ endif} 752 753 htot = htot + dh 754 rhtot = rhtot + rml\_vel(i,k)*dh 755 \textcolor{keywordflow}{ enddo} 756 \textcolor{keywordflow}{if} (rhtot - rml\_vel(i,1)*htot < ustarsq) htot = htot + h\_at\_vel(i,1) 757 \textcolor{keywordflow}{else} 758 \textcolor{keywordflow}{if} (rhtot - gv%Rlay(1)*htot < ustarsq) htot = htot + h\_at\_vel(i,1) 759 \textcolor{keywordflow}{ endif} 760 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_BBL\_EOS} 761 762 \textcolor{comment}{! Value of 2*f at u- or v-points.} 763 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} ; c2f = g%CoriolisBu(i,j-1) + g%CoriolisBu(i,j) 764 \textcolor{keywordflow}{else} ; c2f = g%CoriolisBu(i-1,j) + g%CoriolisBu(i,j) ;\textcolor{keywordflow}{ endif} 765 766 \textcolor{comment}{! The thickness of a rotation limited BBL ignoring stratification is} 767 \textcolor{comment}{! h\_f ~ Cn u* / f (limit of KW99 eq. 2.20 for N->0).} 768 \textcolor{comment}{! The buoyancy limit of BBL thickness (h\_N) is already in the variable htot from above.} 769 \textcolor{comment}{! Substituting x = h\_N/h into KW99 eq. 2.20 yields the quadratic} 770 \textcolor{comment}{! x^2 - x = (h\_N / h\_f)^2} 771 \textcolor{comment}{! for which the positive root is} 772 \textcolor{comment}{! xp = 1/2 + sqrt( 1/4 + (h\_N/h\_f)^2 )} 773 \textcolor{comment}{! and thus h\_bbl = h\_N / xp . Since h\_f = Cn u*/f and Cn=0.5} 774 \textcolor{comment}{! xp = 1/2 + sqrt( 1/4 + (2 f h\_N/u*)^2 )} 775 \textcolor{comment}{! To avoid dividing by zero if u*=0 then} 776 \textcolor{comment}{! xp u* = 1/2 u* + sqrt( 1/4 u*^2 + (2 f h\_N)^2 )} 777 \textcolor{keywordflow}{if} (cs%cdrag * u\_bg\_sq <= 0.0) \textcolor{keywordflow}{then} 778 \textcolor{comment}{! This avoids NaNs and overflows, and could be used in all cases,} 779 \textcolor{comment}{! but is not bitwise identical to the current code.} 780 usth = ustar(i)*gv%Z\_to\_H ; root = sqrt(0.25*usth**2 + (htot*c2f)**2) 781 \textcolor{keywordflow}{if} (htot*usth <= (cs%BBL\_thick\_min+h\_neglect) * (0.5*usth + root)) \textcolor{keywordflow}{then} 782 bbl\_thick = cs%BBL\_thick\_min 783 \textcolor{keywordflow}{else} 784 \textcolor{comment}{! The following expression reads} 785 \textcolor{comment}{! h\_bbl = h\_N u* / ( 1/2 u* + sqrt( 1/4 u*^2 + ( 2 f h\_N )^2 ) )} 786 \textcolor{comment}{! which is h\_bbl = h\_N u*/(xp u*) as described above.} 787 bbl\_thick = (htot * usth) / (0.5*usth + root) 788 \textcolor{keywordflow}{ endif} 789 \textcolor{keywordflow}{else} 790 \textcolor{comment}{! The following expression reads} 791 \textcolor{comment}{! h\_bbl = h\_N / ( 1/2 + sqrt( 1/4 + ( 2 f h\_N / u* )^2 ) )} 792 \textcolor{comment}{! which is h\_bbl = h\_N/xp as described above.} 793 bbl\_thick = htot / (0.5 + sqrt(0.25 + htot*htot*c2f*c2f/ & 794 ((ustar(i)*ustar(i)) * (gv%Z\_to\_H**2)) ) ) 795 796 \textcolor{keywordflow}{if} (bbl\_thick < cs%BBL\_thick\_min) bbl\_thick = cs%BBL\_thick\_min 797 \textcolor{keywordflow}{ endif} 798 799 \textcolor{comment}{! If there is Richardson number dependent mixing, that determines} 800 \textcolor{comment}{! the vertical extent of the bottom boundary layer, and there is no} 801 \textcolor{comment}{! need to set that scale here. In fact, viscously reducing the} 802 \textcolor{comment}{! shears over an excessively large region reduces the efficacy of} 803 \textcolor{comment}{! the Richardson number dependent mixing.} 804 \textcolor{comment}{! In other words, if using RiNo\_mix then CS%Hbbl acts as an upper bound on} 805 \textcolor{comment}{! bbl\_thick.} 806 \textcolor{keywordflow}{if} ((bbl\_thick > 0.5*cs%Hbbl) .and. (cs%RiNo\_mix)) bbl\_thick = 0.5*cs%Hbbl 807 808 \textcolor{keywordflow}{if} (cs%Channel\_drag) \textcolor{keywordflow}{then} 809 \textcolor{comment}{! The drag within the bottommost bbl\_thick is applied as a part of} 810 \textcolor{comment}{! an enhanced bottom viscosity, while above this the drag is applied} 811 \textcolor{comment}{! directly to the layers in question as a Rayleigh drag term.} 812 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} 813 d\_vel = d\_u(i,j) 814 tmp = g%mask2dCu(i,j+1) * d\_u(i,j+1) 815 dp = 2.0 * d\_vel * tmp / (d\_vel + tmp) 816 tmp = g%mask2dCu(i,j-1) * d\_u(i,j-1) 817 dm = 2.0 * d\_vel * tmp / (d\_vel + tmp) 818 \textcolor{keywordflow}{else} 819 d\_vel = d\_v(i,j) 820 tmp = g%mask2dCv(i+1,j) * d\_v(i+1,j) 821 dp = 2.0 * d\_vel * tmp / (d\_vel + tmp) 822 tmp = g%mask2dCv(i-1,j) * d\_v(i-1,j) 823 dm = 2.0 * d\_vel * tmp / (d\_vel + tmp) 824 \textcolor{keywordflow}{ endif} 825 \textcolor{keywordflow}{if} (dm > dp) \textcolor{keywordflow}{then} ; tmp = dp ; dp = dm ; dm = tmp ;\textcolor{keywordflow}{ endif} 826 827 \textcolor{comment}{! Convert the D's to the units of thickness.} 828 dp = gv%Z\_to\_H*dp ; dm = gv%Z\_to\_H*dm ; d\_vel = gv%Z\_to\_H*d\_vel 829 830 a\_3 = (dp + dm - 2.0*d\_vel) ; a = 3.0*a\_3 ; a\_12 = 0.25*a\_3 831 slope = dp - dm 832 \textcolor{comment}{! If the curvature is small enough, there is no reason not to assume} 833 \textcolor{comment}{! a uniformly sloping or flat bottom.} 834 \textcolor{keywordflow}{if} (abs(a) < 1e-2*(slope + cs%BBL\_thick\_min)) a = 0.0 835 \textcolor{comment}{! Each cell extends from x=-1/2 to 1/2, and has a topography} 836 \textcolor{comment}{! given by D(x) = a*x^2 + b*x + D - a/12.} 837 838 \textcolor{comment}{! Calculate the volume above which the entire cell is open and the} 839 \textcolor{comment}{! other volumes at which the equation that is solved for L changes.} 840 \textcolor{keywordflow}{if} (a > 0.0) \textcolor{keywordflow}{then} 841 \textcolor{keywordflow}{if} (slope >= a) \textcolor{keywordflow}{then} 842 vol\_open = d\_vel - dm ; vol\_2\_reg = vol\_open 843 \textcolor{keywordflow}{else} 844 tmp = slope/a 845 vol\_open = 0.25*slope*tmp + c1\_12*a 846 vol\_2\_reg = 0.5*tmp**2 * (a - c1\_3*slope) 847 \textcolor{keywordflow}{ endif} 848 \textcolor{comment}{! Define some combinations of a & b for later use.} 849 c24\_a = 24.0/a ; iapb = 1.0/(a+slope) 850 apb\_4a = (slope+a)/(4.0*a) ; a2x48\_apb3 = (48.0*(a*a))*(iapb**3) 851 ax2\_3apb = 2.0*c1\_3*a*iapb 852 \textcolor{keywordflow}{elseif} (a == 0.0) \textcolor{keywordflow}{then} 853 vol\_open = 0.5*slope 854 \textcolor{keywordflow}{if} (slope > 0) iapb = 1.0/slope 855 \textcolor{keywordflow}{else} \textcolor{comment}{! a < 0.0} 856 vol\_open = d\_vel - dm 857 \textcolor{keywordflow}{if} (slope >= -a) \textcolor{keywordflow}{then} 858 iapb = 1.0e30 ; \textcolor{keywordflow}{if} (slope+a /= 0.0) iapb = 1.0/(a+slope) 859 vol\_direct = 0.0 ; l\_direct = 0.0 ; c24\_a = 0.0 860 \textcolor{keywordflow}{else} 861 c24\_a = 24.0/a ; iapb = 1.0/(a+slope) 862 l\_direct = 1.0 + slope/a \textcolor{comment}{! L\_direct < 1 because a < 0} 863 vol\_direct = -c1\_6*a*l\_direct**3 864 \textcolor{keywordflow}{ endif} 865 ibma\_2 = 2.0 / (slope - a) 866 \textcolor{keywordflow}{ endif} 867 868 l(nz+1) = 0.0 ; vol = 0.0 ; vol\_err = 0.0 ; bbl\_visc\_frac = 0.0 869 \textcolor{comment}{! Determine the normalized open length at each interface.} 870 \textcolor{keywordflow}{do} k=nz,1,-1 871 vol\_below = vol 872 873 vol = vol + h\_vel(i,k) 874 h\_vel\_pos = h\_vel(i,k) + h\_neglect 875 876 \textcolor{keywordflow}{if} (vol >= vol\_open) \textcolor{keywordflow}{then} ; l(k) = 1.0 877 \textcolor{keywordflow}{elseif} (a == 0) \textcolor{keywordflow}{then} \textcolor{comment}{! The bottom has no curvature.} 878 l(k) = sqrt(2.0*vol*iapb) 879 \textcolor{keywordflow}{elseif} (a > 0) \textcolor{keywordflow}{then} 880 \textcolor{comment}{! There may be a minimum depth, and there are} 881 \textcolor{comment}{! analytic expressions for L for all cases.} 882 \textcolor{keywordflow}{if} (vol < vol\_2\_reg) \textcolor{keywordflow}{then} 883 \textcolor{comment}{! In this case, there is a contiguous open region and} 884 \textcolor{comment}{! vol = 0.5*L^2*(slope + a/3*(3-4L)).} 885 \textcolor{keywordflow}{if} (a2x48\_apb3*vol < 1e-8) \textcolor{keywordflow}{then} \textcolor{comment}{! Could be 1e-7?} 886 \textcolor{comment}{! There is a very good approximation here for massless layers.} 887 l0 = sqrt(2.0*vol*iapb) ; l(k) = l0*(1.0 + ax2\_3apb*l0) 888 \textcolor{keywordflow}{else} 889 l(k) = apb\_4a * (1.0 - & 890 2.0 * cos(c1\_3*acos(a2x48\_apb3*vol - 1.0) - c2pi\_3)) 891 \textcolor{keywordflow}{ endif} 892 \textcolor{comment}{! To check the answers.} 893 \textcolor{comment}{! Vol\_err = 0.5*(L(K)*L(K))*(slope + a\_3*(3.0-4.0*L(K))) - vol} 894 \textcolor{keywordflow}{else} \textcolor{comment}{! There are two separate open regions.} 895 \textcolor{comment}{! vol = slope^2/4a + a/12 - (a/12)*(1-L)^2*(1+2L)} 896 \textcolor{comment}{! At the deepest volume, L = slope/a, at the top L = 1.} 897 \textcolor{comment}{!L(K) = 0.5 - cos(C1\_3*acos(1.0 - C24\_a*(Vol\_open - vol)) - C2pi\_3)} 898 tmp\_val\_m1\_to\_p1 = 1.0 - c24\_a*(vol\_open - vol) 899 tmp\_val\_m1\_to\_p1 = max(-1., min(1., tmp\_val\_m1\_to\_p1)) 900 l(k) = 0.5 - cos(c1\_3*acos(tmp\_val\_m1\_to\_p1) - c2pi\_3) 901 \textcolor{comment}{! To check the answers.} 902 \textcolor{comment}{! Vol\_err = Vol\_open - a\_12*(1.0+2.0*L(K)) * (1.0-L(K))**2 - vol} 903 \textcolor{keywordflow}{ endif} 904 \textcolor{keywordflow}{else} \textcolor{comment}{! a < 0.} 905 \textcolor{keywordflow}{if} (vol <= vol\_direct) \textcolor{keywordflow}{then} 906 \textcolor{comment}{! Both edges of the cell are bounded by walls.} 907 l(k) = (-0.25*c24\_a*vol)**c1\_3 908 \textcolor{keywordflow}{else} 909 \textcolor{comment}{! x\_R is at 1/2 but x\_L is in the interior & L is found by solving} 910 \textcolor{comment}{! vol = 0.5*L^2*(slope + a/3*(3-4L))} 911 912 \textcolor{comment}{! Vol\_err = 0.5*(L(K+1)*L(K+1))*(slope + a\_3*(3.0-4.0*L(K+1))) - vol\_below} 913 \textcolor{comment}{! Change to ...} 914 \textcolor{comment}{! if (min(Vol\_below + Vol\_err, vol) <= Vol\_direct) then ?} 915 \textcolor{keywordflow}{if} (vol\_below + vol\_err <= vol\_direct) \textcolor{keywordflow}{then} 916 l0 = l\_direct ; vol\_0 = vol\_direct 917 \textcolor{keywordflow}{else} 918 l0 = l(k+1) ; vol\_0 = vol\_below + vol\_err 919 \textcolor{comment}{! Change to Vol\_0 = min(Vol\_below + Vol\_err, vol) ?} 920 \textcolor{keywordflow}{ endif} 921 922 \textcolor{comment}{! Try a relatively simple solution that usually works well} 923 \textcolor{comment}{! for massless layers.} 924 dv\_dl2 = 0.5*(slope+a) - a*l0 ; dvol = (vol-vol\_0) 925 \textcolor{comment}{! dV\_dL2 = 0.5*(slope+a) - a*L0 ; dVol = max(vol-Vol\_0, 0.0)} 926 927 use\_l0 = .false. 928 do\_one\_l\_iter = .false. 929 \textcolor{keywordflow}{if} (cs%answers\_2018) \textcolor{keywordflow}{then} 930 curv\_tol = gv%Angstrom\_H*dv\_dl2**2 & 931 * (0.25 * dv\_dl2 * gv%Angstrom\_H - a * l0 * dvol) 932 do\_one\_l\_iter = (a * a * dvol**3) < curv\_tol 933 \textcolor{keywordflow}{else} 934 \textcolor{comment}{! The following code is more robust when GV%Angstrom\_H=0, but} 935 \textcolor{comment}{! it changes answers.} 936 use\_l0 = (dvol <= 0.) 937 938 vol\_tol = max(0.5 * gv%Angstrom\_H + gv%H\_subroundoff, 1e-14 * vol) 939 vol\_quit = max(0.9 * gv%Angstrom\_H + gv%H\_subroundoff, 1e-14 * vol) 940 941 curv\_tol = vol\_tol * dv\_dl2**2 & 942 * (dv\_dl2 * vol\_tol - 2.0 * a * l0 * dvol) 943 do\_one\_l\_iter = (a * a * dvol**3) < curv\_tol 944 \textcolor{keywordflow}{ endif} 945 946 \textcolor{keywordflow}{if} (use\_l0) \textcolor{keywordflow}{then} 947 l(k) = l0 948 vol\_err = 0.5*(l(k)*l(k))*(slope + a\_3*(3.0-4.0*l(k))) - vol 949 \textcolor{keywordflow}{elseif} (do\_one\_l\_iter) \textcolor{keywordflow}{then} 950 \textcolor{comment}{! One iteration of Newton's method should give an estimate} 951 \textcolor{comment}{! that is accurate to within Vol\_tol.} 952 l(k) = sqrt(l0*l0 + dvol / dv\_dl2) 953 vol\_err = 0.5*(l(k)*l(k))*(slope + a\_3*(3.0-4.0*l(k))) - vol 954 \textcolor{keywordflow}{else} 955 \textcolor{keywordflow}{if} (dv\_dl2*(1.0-l0*l0) < dvol + & 956 dv\_dl2 * (vol\_open - vol)*ibma\_2) \textcolor{keywordflow}{then} 957 l\_max = sqrt(1.0 - (vol\_open - vol)*ibma\_2) 958 \textcolor{keywordflow}{else} 959 l\_max = sqrt(l0*l0 + dvol / dv\_dl2) 960 \textcolor{keywordflow}{ endif} 961 l\_min = sqrt(l0*l0 + dvol / (0.5*(slope+a) - a*l\_max)) 962 963 vol\_err\_min = 0.5*(l\_min**2)*(slope + a\_3*(3.0-4.0*l\_min)) - vol 964 vol\_err\_max = 0.5*(l\_max**2)*(slope + a\_3*(3.0-4.0*l\_max)) - vol 965 \textcolor{comment}{! if ((abs(Vol\_err\_min) <= Vol\_quit) .or. (Vol\_err\_min >= Vol\_err\_max)) then} 966 \textcolor{keywordflow}{if} (abs(vol\_err\_min) <= vol\_quit) \textcolor{keywordflow}{then} 967 l(k) = l\_min ; vol\_err = vol\_err\_min 968 \textcolor{keywordflow}{else} 969 l(k) = sqrt((l\_min**2*vol\_err\_max - l\_max**2*vol\_err\_min) / & 970 (vol\_err\_max - vol\_err\_min)) 971 \textcolor{keywordflow}{do} itt=1,maxitt 972 vol\_err = 0.5*(l(k)*l(k))*(slope + a\_3*(3.0-4.0*l(k))) - vol 973 \textcolor{keywordflow}{if} (abs(vol\_err) <= vol\_quit) \textcolor{keywordflow}{exit} 974 \textcolor{comment}{! Take a Newton's method iteration. This equation has proven} 975 \textcolor{comment}{! robust enough not to need bracketing.} 976 l(k) = l(k) - vol\_err / (l(k)* (slope + a - 2.0*a*l(k))) 977 \textcolor{comment}{! This would be a Newton's method iteration for L^2:} 978 \textcolor{comment}{! L(K) = sqrt(L(K)*L(K) - Vol\_err / (0.5*(slope+a) - a*L(K)))} 979 \textcolor{keywordflow}{ enddo} 980 \textcolor{keywordflow}{ endif} \textcolor{comment}{! end of iterative solver} 981 \textcolor{keywordflow}{ endif} \textcolor{comment}{! end of 1-boundary alternatives.} 982 \textcolor{keywordflow}{ endif} \textcolor{comment}{! end of a<0 cases.} 983 \textcolor{keywordflow}{ endif} 984 985 \textcolor{comment}{! Determine the drag contributing to the bottom boundary layer} 986 \textcolor{comment}{! and the Raleigh drag that acts on each layer.} 987 \textcolor{keywordflow}{if} (l(k) > l(k+1)) \textcolor{keywordflow}{then} 988 \textcolor{keywordflow}{if} (vol\_below < bbl\_thick) \textcolor{keywordflow}{then} 989 bbl\_frac = (1.0-vol\_below/bbl\_thick)**2 990 bbl\_visc\_frac = bbl\_visc\_frac + bbl\_frac*(l(k) - l(k+1)) 991 \textcolor{keywordflow}{else} 992 bbl\_frac = 0.0 993 \textcolor{keywordflow}{ endif} 994 995 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} ; cell\_width = g%dy\_Cu(i,j) 996 \textcolor{keywordflow}{else} ; cell\_width = g%dx\_Cv(i,j) ;\textcolor{keywordflow}{ endif} 997 gam = 1.0 - l(k+1)/l(k) 998 rayleigh = us%L\_to\_Z * cs%cdrag * (l(k)-l(k+1)) * (1.0-bbl\_frac) * & 999 (12.0*cs%c\_Smag*h\_vel\_pos) / (12.0*cs%c\_Smag*h\_vel\_pos + & 1000 us%L\_to\_Z*gv%Z\_to\_H * cs%cdrag * gam*(1.0-gam)*(1.0-1.5*gam) * l(k)**2 * cell\_width) 1001 \textcolor{keywordflow}{else} \textcolor{comment}{! This layer feels no drag.} 1002 rayleigh = 0.0 1003 \textcolor{keywordflow}{ endif} 1004 1005 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} 1006 \textcolor{keywordflow}{if} (rayleigh > 0.0) \textcolor{keywordflow}{then} 1007 v\_at\_u = set\_v\_at\_u(v, h, g, i, j, k, mask\_v, obc) 1008 visc%Ray\_u(i,j,k) = rayleigh*sqrt(u(i,j,k)*u(i,j,k) + & 1009 v\_at\_u*v\_at\_u + u\_bg\_sq) 1010 \textcolor{keywordflow}{else} ; visc%Ray\_u(i,j,k) = 0.0 ;\textcolor{keywordflow}{ endif} 1011 \textcolor{keywordflow}{else} 1012 \textcolor{keywordflow}{if} (rayleigh > 0.0) \textcolor{keywordflow}{then} 1013 u\_at\_v = set\_u\_at\_v(u, h, g, i, j, k, mask\_u, obc) 1014 visc%Ray\_v(i,j,k) = rayleigh*sqrt(v(i,j,k)*v(i,j,k) + & 1015 u\_at\_v*u\_at\_v + u\_bg\_sq) 1016 \textcolor{keywordflow}{else} ; visc%Ray\_v(i,j,k) = 0.0 ;\textcolor{keywordflow}{ endif} 1017 \textcolor{keywordflow}{ endif} 1018 1019 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! k loop to determine L(K).} 1020 1021 \textcolor{comment}{! Set the near-bottom viscosity to a value which will give} 1022 \textcolor{comment}{! the correct stress when the shear occurs over bbl\_thick.} 1023 \textcolor{comment}{! See next block for explanation.} 1024 bbl\_thick\_z = bbl\_thick * gv%H\_to\_Z 1025 \textcolor{keywordflow}{if} (cs%correct\_BBL\_bounds .and. & 1026 cdrag\_sqrt*ustar(i)*bbl\_thick\_z*bbl\_visc\_frac <= cs%Kv\_BBL\_min) \textcolor{keywordflow}{then} 1027 \textcolor{comment}{! If the bottom stress implies less viscosity than Kv\_BBL\_min then} 1028 \textcolor{comment}{! set kv\_bbl to the bound and recompute bbl\_thick to be consistent} 1029 \textcolor{comment}{! but with a ridiculously large upper bound on thickness (for Cd u*=0)} 1030 kv\_bbl = cs%Kv\_BBL\_min 1031 \textcolor{keywordflow}{if} (cdrag\_sqrt*ustar(i)*bbl\_visc\_frac*g%Rad\_Earth*us%m\_to\_Z > kv\_bbl) \textcolor{keywordflow}{then} 1032 bbl\_thick\_z = kv\_bbl / ( cdrag\_sqrt*ustar(i)*bbl\_visc\_frac ) 1033 \textcolor{keywordflow}{else} 1034 bbl\_thick\_z = g%Rad\_Earth * us%m\_to\_Z 1035 \textcolor{keywordflow}{ endif} 1036 \textcolor{keywordflow}{else} 1037 kv\_bbl = cdrag\_sqrt*ustar(i)*bbl\_thick\_z*bbl\_visc\_frac 1038 \textcolor{keywordflow}{ endif} 1039 1040 \textcolor{keywordflow}{else} \textcolor{comment}{! Not Channel\_drag.} 1041 \textcolor{comment}{! Set the near-bottom viscosity to a value which will give} 1042 \textcolor{comment}{! the correct stress when the shear occurs over bbl\_thick.} 1043 \textcolor{comment}{! - The bottom stress is tau\_b = Cdrag * u\_bbl^2} 1044 \textcolor{comment}{! - u\_bbl was calculated by averaging flow over CS%Hbbl} 1045 \textcolor{comment}{! (and includes unresolved tidal components)} 1046 \textcolor{comment}{! - u\_bbl is embedded in u* since u*^2 = Cdrag u\_bbl^2} 1047 \textcolor{comment}{! - The average shear in the BBL is du/dz = 2 * u\_bbl / h\_bbl} 1048 \textcolor{comment}{! (which assumes a linear profile, hence the "2")} 1049 \textcolor{comment}{! - bbl\_thick was bounded to <= 0.5 * CS%Hbbl} 1050 \textcolor{comment}{! - The viscous stress kv\_bbl du/dz should balance tau\_b} 1051 \textcolor{comment}{! Cdrag u\_bbl^2 = kv\_bbl du/dz} 1052 \textcolor{comment}{! = 2 kv\_bbl u\_bbl} 1053 \textcolor{comment}{! so} 1054 \textcolor{comment}{! kv\_bbl = 0.5 h\_bbl Cdrag u\_bbl} 1055 \textcolor{comment}{! = 0.5 h\_bbl sqrt(Cdrag) u*} 1056 bbl\_thick\_z = bbl\_thick * gv%H\_to\_Z 1057 \textcolor{keywordflow}{if} (cs%correct\_BBL\_bounds .and. & 1058 cdrag\_sqrt*ustar(i)*bbl\_thick\_z <= cs%Kv\_BBL\_min) \textcolor{keywordflow}{then} 1059 \textcolor{comment}{! If the bottom stress implies less viscosity than Kv\_BBL\_min then} 1060 \textcolor{comment}{! set kv\_bbl to the bound and recompute bbl\_thick to be consistent} 1061 \textcolor{comment}{! but with a ridiculously large upper bound on thickness (for Cd u*=0)} 1062 kv\_bbl = cs%Kv\_BBL\_min 1063 \textcolor{keywordflow}{if} (cdrag\_sqrt*ustar(i)*g%Rad\_Earth*us%m\_to\_Z > kv\_bbl) \textcolor{keywordflow}{then} 1064 bbl\_thick\_z = kv\_bbl / ( cdrag\_sqrt*ustar(i) ) 1065 \textcolor{keywordflow}{else} 1066 bbl\_thick\_z = g%Rad\_Earth * us%m\_to\_Z 1067 \textcolor{keywordflow}{ endif} 1068 \textcolor{keywordflow}{else} 1069 kv\_bbl = cdrag\_sqrt*ustar(i)*bbl\_thick\_z 1070 \textcolor{keywordflow}{ endif} 1071 \textcolor{keywordflow}{ endif} 1072 kv\_bbl = max(cs%Kv\_BBL\_min, kv\_bbl) 1073 \textcolor{keywordflow}{if} (m==1) \textcolor{keywordflow}{then} 1074 visc%Kv\_bbl\_u(i,j) = kv\_bbl 1075 visc%bbl\_thick\_u(i,j) = bbl\_thick\_z 1076 \textcolor{keywordflow}{else} 1077 visc%Kv\_bbl\_v(i,j) = kv\_bbl 1078 visc%bbl\_thick\_v(i,j) = bbl\_thick\_z 1079 \textcolor{keywordflow}{ endif} 1080 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of i loop} 1081 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of m & j loops} 1082 1083 \textcolor{comment}{! Offer diagnostics for averaging} 1084 \textcolor{keywordflow}{if} (cs%id\_bbl\_thick\_u > 0) & 1085 \textcolor{keyword}{call }post\_data(cs%id\_bbl\_thick\_u, visc%bbl\_thick\_u, cs%diag) 1086 \textcolor{keywordflow}{if} (cs%id\_kv\_bbl\_u > 0) & 1087 \textcolor{keyword}{call }post\_data(cs%id\_kv\_bbl\_u, visc%kv\_bbl\_u, cs%diag) 1088 \textcolor{keywordflow}{if} (cs%id\_bbl\_u > 0) & 1089 \textcolor{keyword}{call }post\_data(cs%id\_bbl\_u, cs%bbl\_u, cs%diag) 1090 \textcolor{keywordflow}{if} (cs%id\_bbl\_thick\_v > 0) & 1091 \textcolor{keyword}{call }post\_data(cs%id\_bbl\_thick\_v, visc%bbl\_thick\_v, cs%diag) 1092 \textcolor{keywordflow}{if} (cs%id\_kv\_bbl\_v > 0) & 1093 \textcolor{keyword}{call }post\_data(cs%id\_kv\_bbl\_v, visc%kv\_bbl\_v, cs%diag) 1094 \textcolor{keywordflow}{if} (cs%id\_bbl\_v > 0) & 1095 \textcolor{keyword}{call }post\_data(cs%id\_bbl\_v, cs%bbl\_v, cs%diag) 1096 \textcolor{keywordflow}{if} (cs%id\_Ray\_u > 0) & 1097 \textcolor{keyword}{call }post\_data(cs%id\_Ray\_u, visc%Ray\_u, cs%diag) 1098 \textcolor{keywordflow}{if} (cs%id\_Ray\_v > 0) & 1099 \textcolor{keyword}{call }post\_data(cs%id\_Ray\_v, visc%Ray\_v, cs%diag) 1100 1101 \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then} 1102 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%Ray\_u) .and. \textcolor{keyword}{associated}(visc%Ray\_v)) & 1103 \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"Ray [uv]"}, visc%Ray\_u, visc%Ray\_v, g%HI, haloshift=0, scale=us%Z\_to\_m*us%s\_to\_T) 1104 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%kv\_bbl\_u) .and. \textcolor{keyword}{associated}(visc%kv\_bbl\_v)) & 1105 \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"kv\_bbl\_[uv]"}, visc%kv\_bbl\_u, visc%kv\_bbl\_v, g%HI, & 1106 haloshift=0, scale=us%Z2\_T\_to\_m2\_s, scalar\_pair=.true.) 1107 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%bbl\_thick\_u) .and. \textcolor{keyword}{associated}(visc%bbl\_thick\_v)) & 1108 \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"bbl\_thick\_[uv]"}, visc%bbl\_thick\_u, visc%bbl\_thick\_v, & 1109 g%HI, haloshift=0, scale=us%Z\_to\_m, scalar\_pair=.true.) 1110 \textcolor{keywordflow}{ endif} 1111 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__set__visc_aba41cd4f8baa1cda9036d97087ce8a22}\label{namespacemom__set__visc_aba41cd4f8baa1cda9036d97087ce8a22}} \index{mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}!set\+\_\+viscous\+\_\+ml@{set\+\_\+viscous\+\_\+ml}} \index{set\+\_\+viscous\+\_\+ml@{set\+\_\+viscous\+\_\+ml}!mom\+\_\+set\+\_\+visc@{mom\+\_\+set\+\_\+visc}} \subsubsection{\texorpdfstring{set\+\_\+viscous\+\_\+ml()}{set\_viscous\_ml()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+set\+\_\+visc\+::set\+\_\+viscous\+\_\+ml (\begin{DoxyParamCaption}\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[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{type(mech\+\_\+forcing), intent(in)}]{forces, }\item[{type(vertvisc\+\_\+type), intent(inout)}]{visc, }\item[{real, intent(in)}]{dt, }\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(\hyperlink{structmom__set__visc_1_1set__visc__cs}{set\+\_\+visc\+\_\+cs}), pointer}]{CS, }\item[{logical, intent(in), optional}]{symmetrize }\end{DoxyParamCaption})} Calculates the thickness of the surface boundary layer for applying an elevated viscosity. A bulk Richardson criterion or the thickness of the topmost N\+K\+ML layers (with a bulk mixed layer) are currently used. The thicknesses are given in terms of fractional layers, so that this thickness will move as the thickness of the topmost layers change. \begin{DoxyParams}[1]{Parameters} \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 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} & A structure containing pointers to any available thermodynamic fields. Absent fields have N\+U\+LL ptrs.\\ \hline \mbox{\tt in} & {\em forces} & A structure with the driving mechanical forces\\ \hline \mbox{\tt in,out} & {\em visc} & A structure containing vertical viscosities and related fields.\\ \hline \mbox{\tt in} & {\em dt} & Time increment \mbox{[}T $\sim$$>$ s\mbox{]}.\\ \hline & {\em cs} & The control structure returned by a previous call to set\+\_\+visc\+\_\+init.\\ \hline \mbox{\tt in} & {\em symmetrize} & If present and true, do extra calculations of those values in visc that would be calculated with symmetric memory. \\ \hline \end{DoxyParams} Definition at line 1208 of file M\+O\+M\+\_\+set\+\_\+viscosity.\+F90. \begin{DoxyCode} 1208 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(inout)} :: g\textcolor{comment}{ !< The ocean's grid structure.} 1209 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< The ocean's vertical grid structure.} 1210 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 1211 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, & 1212 \textcolor{keywordtype}{intent(in)} :: u\textcolor{comment}{ !< The zonal velocity [L T-1 ~> m s-1].} 1213 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, & 1214 \textcolor{keywordtype}{intent(in)} :: v\textcolor{comment}{ !< The meridional velocity [L T-1 ~> m s-1].} 1215 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, & 1216 \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2].} 1217 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< A structure containing pointers to any available} 1218 \textcolor{comment}{ !! thermodynamic fields. Absent fields have} 1219 \textcolor{comment}{ !! NULL ptrs.} 1220 \textcolor{keywordtype}{type}(mech\_forcing), \textcolor{keywordtype}{intent(in)} :: forces\textcolor{comment}{ !< A structure with the driving mechanical forces} 1221 \textcolor{keywordtype}{type}(vertvisc\_type), \textcolor{keywordtype}{intent(inout)} :: visc\textcolor{comment}{ !< A structure containing vertical viscosities and} 1222 \textcolor{comment}{ !! related fields.} 1223 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: dt\textcolor{comment}{ !< Time increment [T ~> s].} 1224 \textcolor{keywordtype}{type}(set\_visc\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< The control structure returned by a previous} 1225 \textcolor{comment}{ !! call to set\_visc\_init.} 1226 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: symmetrize\textcolor{comment}{ !< If present and true, do extra calculations} 1227 \textcolor{comment}{ !! of those values in visc that would be} 1228 \textcolor{comment}{ !! calculated with symmetric memory.} 1229 1230 \textcolor{comment}{! Local variables} 1231 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G))} :: & 1232 htot, & \textcolor{comment}{! The total depth of the layers being that are within the} 1233 \textcolor{comment}{! surface mixed layer [H ~> m or kg m-2].} 1234 thtot, & \textcolor{comment}{! The integrated temperature of layers that are within the} 1235 \textcolor{comment}{! surface mixed layer [H degC ~> m degC or kg degC m-2].} 1236 shtot, & \textcolor{comment}{! The integrated salt of layers that are within the} 1237 \textcolor{comment}{! surface mixed layer [H ppt ~> m ppt or kg ppt m-2].} 1238 rhtot, & \textcolor{comment}{! The integrated density of layers that are within the surface mixed layer} 1239 \textcolor{comment}{! [H R ~> kg m-2 or kg2 m-5]. Rhtot is only used if no} 1240 \textcolor{comment}{! equation of state is used.} 1241 uhtot, & \textcolor{comment}{! The depth integrated zonal and meridional velocities within} 1242 vhtot, & \textcolor{comment}{! the surface mixed layer [H L T-1 ~> m2 s-1 or kg m-1 s-1].} 1243 idecay\_len\_tke, & \textcolor{comment}{! The inverse of a turbulence decay length scale [H-1 ~> m-1 or m2 kg-1].} 1244 dr\_dt, & \textcolor{comment}{! Partial derivative of the density at the base of layer nkml} 1245 \textcolor{comment}{! (roughly the base of the mixed layer) with temperature [R degC-1 ~> kg m-3 degC-1].} 1246 dr\_ds, & \textcolor{comment}{! Partial derivative of the density at the base of layer nkml} 1247 \textcolor{comment}{! (roughly the base of the mixed layer) with salinity [R ppt-1 ~> kg m-3 ppt-1].} 1248 ustar, & \textcolor{comment}{! The surface friction velocity under ice shelves [Z T-1 ~> m s-1].} 1249 press, & \textcolor{comment}{! The pressure at which dR\_dT and dR\_dS are evaluated [R L2 T-2 ~> Pa].} 1250 t\_eos, & \textcolor{comment}{! The potential temperature at which dR\_dT and dR\_dS are evaluated [degC]} 1251 s\_eos \textcolor{comment}{! The salinity at which dR\_dT and dR\_dS are evaluated [ppt].} 1252 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G))} :: & 1253 mask\_u \textcolor{comment}{! A mask that disables any contributions from u points that} 1254 \textcolor{comment}{! are land or past open boundary conditions [nondim], 0 or 1.} 1255 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G))} :: & 1256 mask\_v \textcolor{comment}{! A mask that disables any contributions from v points that} 1257 \textcolor{comment}{! are land or past open boundary conditions [nondim], 0 or 1.} 1258 \textcolor{keywordtype}{real} :: h\_at\_vel(szib\_(g),szk\_(g))\textcolor{comment}{! Layer thickness at velocity points,} 1259 \textcolor{comment}{! using an upwind-biased second order accurate estimate based} 1260 \textcolor{comment}{! on the previous velocity direction [H ~> m or kg m-2].} 1261 \textcolor{keywordtype}{integer} :: k\_massive(szib\_(g)) \textcolor{comment}{! The k-index of the deepest layer yet found} 1262 \textcolor{comment}{! that has more than h\_tiny thickness and will be in the} 1263 \textcolor{comment}{! viscous mixed layer.} 1264 \textcolor{keywordtype}{real} :: uh2 \textcolor{comment}{! The squared magnitude of the difference between the velocity} 1265 \textcolor{comment}{! integrated through the mixed layer and the velocity of the} 1266 \textcolor{comment}{! interior layer layer times the depth of the the mixed layer} 1267 \textcolor{comment}{! [H2 Z2 T-2 ~> m4 s-2 or kg2 m-2 s-2].} 1268 \textcolor{keywordtype}{real} :: htot\_vel \textcolor{comment}{! Sum of the layer thicknesses up to some point [H ~> m or kg m-2].} 1269 \textcolor{keywordtype}{real} :: hwtot \textcolor{comment}{! Sum of the thicknesses used to calculate} 1270 \textcolor{comment}{! the near-bottom velocity magnitude [H ~> m or kg m-2].} 1271 \textcolor{keywordtype}{real} :: hutot \textcolor{comment}{! Running sum of thicknesses times the} 1272 \textcolor{comment}{! velocity magnitudes [H L T-1 ~> m2 s-1 or kg m-1 s-1].} 1273 \textcolor{keywordtype}{real} :: hweight \textcolor{comment}{! The thickness of a layer that is within Hbbl} 1274 \textcolor{comment}{! of the bottom [H ~> m or kg m-2].} 1275 \textcolor{keywordtype}{real} :: tbl\_thick\_z \textcolor{comment}{! The thickness of the top boundary layer [Z ~> m].} 1276 1277 \textcolor{keywordtype}{real} :: hlay \textcolor{comment}{! The layer thickness at velocity points [H ~> m or kg m-2].} 1278 \textcolor{keywordtype}{real} :: i\_2hlay \textcolor{comment}{! 1 / 2*hlay [H-1 ~> m-1 or m2 kg-1].} 1279 \textcolor{keywordtype}{real} :: t\_lay \textcolor{comment}{! The layer temperature at velocity points [degC].} 1280 \textcolor{keywordtype}{real} :: s\_lay \textcolor{comment}{! The layer salinity at velocity points [ppt].} 1281 \textcolor{keywordtype}{real} :: rlay \textcolor{comment}{! The layer potential density at velocity points [R ~> kg m-3].} 1282 \textcolor{keywordtype}{real} :: rlb \textcolor{comment}{! The potential density of the layer below [R ~> kg m-3].} 1283 \textcolor{keywordtype}{real} :: v\_at\_u \textcolor{comment}{! The meridonal velocity at a zonal velocity point [L T-1 ~> m s-1].} 1284 \textcolor{keywordtype}{real} :: u\_at\_v \textcolor{comment}{! The zonal velocity at a meridonal velocity point [L T-1 ~> m s-1].} 1285 \textcolor{keywordtype}{real} :: ghprime \textcolor{comment}{! The mixed-layer internal gravity wave speed squared, based} 1286 \textcolor{comment}{! on the mixed layer thickness and density difference across} 1287 \textcolor{comment}{! the base of the mixed layer [L2 T-2 ~> m2 s-2].} 1288 \textcolor{keywordtype}{real} :: ribulk \textcolor{comment}{! The bulk Richardson number below which water is in the} 1289 \textcolor{comment}{! viscous mixed layer, including reduction for turbulent} 1290 \textcolor{comment}{! decay. Nondimensional.} 1291 \textcolor{keywordtype}{real} :: dt\_rho0 \textcolor{comment}{! The time step divided by the conversion from the layer} 1292 \textcolor{comment}{! thickness to layer mass [T H Z-1 R-1 ~> s m3 kg-1 or s].} 1293 \textcolor{keywordtype}{real} :: g\_h\_rho0 \textcolor{comment}{! The gravitational acceleration times the conversion from H to m divided} 1294 \textcolor{comment}{! by the mean density [L2 T-2 H-1 R-1 ~> m4 s-2 kg-1 or m7 s-2 kg-2].} 1295 \textcolor{keywordtype}{real} :: ustarsq \textcolor{comment}{! 400 times the square of ustar, times} 1296 \textcolor{comment}{! Rho0 divided by G\_Earth and the conversion} 1297 \textcolor{comment}{! from m to thickness units [H R ~> kg m-2 or kg2 m-5].} 1298 \textcolor{keywordtype}{real} :: cdrag\_sqrt\_z \textcolor{comment}{! Square root of the drag coefficient, times a unit conversion} 1299 \textcolor{comment}{! factor from lateral lengths to vertical depths [Z L-1 ~> 1].} 1300 \textcolor{keywordtype}{real} :: cdrag\_sqrt \textcolor{comment}{! Square root of the drag coefficient [nondim].} 1301 \textcolor{keywordtype}{real} :: oldfn \textcolor{comment}{! The integrated energy required to} 1302 \textcolor{comment}{! entrain up to the bottom of the layer,} 1303 \textcolor{comment}{! divided by G\_Earth [H R ~> kg m-2 or kg2 m-5].} 1304 \textcolor{keywordtype}{real} :: dfn \textcolor{comment}{! The increment in oldfn for entraining} 1305 \textcolor{comment}{! the layer [H R ~> kg m-2 or kg2 m-5].} 1306 \textcolor{keywordtype}{real} :: dh \textcolor{comment}{! The increment in layer thickness from} 1307 \textcolor{comment}{! the present layer [H ~> m or kg m-2].} 1308 \textcolor{keywordtype}{real} :: u\_bg\_sq \textcolor{comment}{! The square of an assumed background velocity, for} 1309 \textcolor{comment}{! calculating the mean magnitude near the top for use in} 1310 \textcolor{comment}{! the quadratic surface drag [L2 T-2 ~> m2 s-2].} 1311 \textcolor{keywordtype}{real} :: h\_tiny \textcolor{comment}{! A very small thickness [H ~> m or kg m-2]. Layers that are less than} 1312 \textcolor{comment}{! h\_tiny can not be the deepest in the viscous mixed layer.} 1313 \textcolor{keywordtype}{real} :: absf \textcolor{comment}{! The absolute value of f averaged to velocity points [T-1 ~> s-1].} 1314 \textcolor{keywordtype}{real} :: u\_star \textcolor{comment}{! The friction velocity at velocity points [Z T-1 ~> m s-1].} 1315 \textcolor{keywordtype}{real} :: h\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 1316 \textcolor{comment}{! in roundoff and can be neglected [H ~> m or kg m-2].} 1317 \textcolor{keywordtype}{real} :: rho0x400\_g \textcolor{comment}{! 400*Rho0/G\_Earth, times unit conversion factors} 1318 \textcolor{comment}{! [R T2 H Z-2 ~> kg s2 m-4 or kg2 s2 m-7].} 1319 \textcolor{comment}{! The 400 is a constant proposed by Killworth and Edwards, 1999.} 1320 \textcolor{keywordtype}{real} :: ustar1 \textcolor{comment}{! ustar [H T-1 ~> m s-1 or kg m-2 s-1]} 1321 \textcolor{keywordtype}{real} :: h2f2 \textcolor{comment}{! (h*2*f)^2 [H2 T-2 ~> m2 s-2 or kg2 m-4 s-2]} 1322 \textcolor{keywordtype}{logical} :: use\_eos, do\_any, do\_any\_shelf, do\_i(szib\_(g)) 1323 \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, isq, ieq, jsq, jeq, nz, k2, nkmb, nkml, n 1324 \textcolor{keywordtype}{type}(ocean\_obc\_type), \textcolor{keywordtype}{pointer} :: obc => null() 1325 1326 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke 1327 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB 1328 nkmb = gv%nk\_rho\_varies ; nkml = gv%nkml 1329 1330 \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal,\textcolor{stringliteral}{"MOM\_set\_viscosity(visc\_ML): "}//& 1331 \textcolor{stringliteral}{"Module must be initialized before it is used."}) 1332 \textcolor{keywordflow}{if} (.not.(cs%dynamic\_viscous\_ML .or. \textcolor{keyword}{associated}(forces%frac\_shelf\_u) .or. & 1333 \textcolor{keyword}{associated}(forces%frac\_shelf\_v)) ) \textcolor{keywordflow}{return} 1334 1335 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(symmetrize)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (symmetrize) \textcolor{keywordflow}{then} 1336 jsq = js-1 ; isq = is-1 1337 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 1338 1339 rho0x400\_g = 400.0*(gv%Rho0/(us%L\_to\_Z**2 * gv%g\_Earth)) * gv%Z\_to\_H 1340 u\_bg\_sq = cs%drag\_bg\_vel * cs%drag\_bg\_vel 1341 cdrag\_sqrt = sqrt(cs%cdrag) 1342 cdrag\_sqrt\_z = us%L\_to\_Z * sqrt(cs%cdrag) 1343 1344 obc => cs%OBC 1345 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 1346 dt\_rho0 = dt / gv%H\_to\_RZ 1347 h\_neglect = gv%H\_subroundoff 1348 h\_tiny = 2.0*gv%Angstrom\_H + h\_neglect 1349 g\_h\_rho0 = (gv%g\_Earth*gv%H\_to\_Z) / (gv%Rho0) 1350 1351 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(forces%frac\_shelf\_u) .neqv. \textcolor{keyword}{associated}(forces%frac\_shelf\_v)) & 1352 \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"set\_viscous\_ML: one of forces%frac\_shelf\_u and "}//& 1353 \textcolor{stringliteral}{"forces%frac\_shelf\_v is associated, but the other is not."}) 1354 1355 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(forces%frac\_shelf\_u)) \textcolor{keywordflow}{then} 1356 \textcolor{comment}{! This configuration has ice shelves, and the appropriate variables need to be} 1357 \textcolor{comment}{! allocated. If the arrays have already been allocated, these calls do nothing.} 1358 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%tauy\_shelf, g%isd, g%ied, g%JsdB, g%JedB) 1359 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%tbl\_thick\_shelf\_u, g%IsdB, g%IedB, g%jsd, g%jed) 1360 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%tbl\_thick\_shelf\_v, g%isd, g%ied, g%JsdB, g%JedB) 1361 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%kv\_tbl\_shelf\_u, g%IsdB, g%IedB, g%jsd, g%jed) 1362 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%kv\_tbl\_shelf\_v, g%isd, g%ied, g%JsdB, g%JedB) 1363 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%taux\_shelf, g%IsdB, g%IedB, g%jsd, g%jed) 1364 \textcolor{keyword}{call }safe\_alloc\_ptr(visc%tauy\_shelf, g%isd, g%ied, g%JsdB, g%JedB) 1365 1366 \textcolor{comment}{! With a linear drag law under shelves, the friction velocity is already known.} 1367 \textcolor{comment}{! if (CS%linear\_drag) ustar(:) = cdrag\_sqrt\_Z*CS%drag\_bg\_vel} 1368 \textcolor{keywordflow}{ endif} 1369 1370 \textcolor{comment}{!$OMP parallel do default(shared)} 1371 \textcolor{keywordflow}{do} j=js-1,je ; \textcolor{keywordflow}{do} i=is-1,ie+1 1372 mask\_v(i,j) = g%mask2dCv(i,j) 1373 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1374 \textcolor{comment}{!$OMP parallel do default(shared)} 1375 \textcolor{keywordflow}{do} j=js-1,je+1 ; \textcolor{keywordflow}{do} i=is-1,ie 1376 mask\_u(i,j) = g%mask2dCu(i,j) 1377 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1378 1379 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} n=1,obc%number\_of\_segments 1380 \textcolor{comment}{! Now project bottom depths across cell-corner points in the OBCs. The two} 1381 \textcolor{comment}{! projections have to occur in sequence and can not be combined easily.} 1382 \textcolor{keywordflow}{if} (.not. obc%segment(n)%on\_pe) cycle 1383 \textcolor{comment}{! Use a one-sided projection of bottom depths at OBC points.} 1384 i = obc%segment(n)%HI%IsdB ; j = obc%segment(n)%HI%JsdB 1385 \textcolor{keywordflow}{if} (obc%segment(n)%is\_N\_or\_S .and. (j >= js-1) .and. (j <= je)) \textcolor{keywordflow}{then} 1386 \textcolor{keywordflow}{do} i = max(is-1,obc%segment(n)%HI%IsdB), min(ie,obc%segment(n)%HI%IedB) 1387 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_n) mask\_u(i,j+1) = 0.0 1388 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_s) mask\_u(i,j) = 0.0 1389 \textcolor{keywordflow}{ enddo} 1390 \textcolor{keywordflow}{elseif} (obc%segment(n)%is\_E\_or\_W .and. (i >= is-1) .and. (i <= je)) \textcolor{keywordflow}{then} 1391 \textcolor{keywordflow}{do} j = max(js-1,obc%segment(n)%HI%JsdB), min(je,obc%segment(n)%HI%JedB) 1392 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_e) mask\_v(i+1,j) = 0.0 1393 \textcolor{keywordflow}{if} (obc%segment(n)%direction == obc\_direction\_w) mask\_v(i,j) = 0.0 1394 \textcolor{keywordflow}{ enddo} 1395 \textcolor{keywordflow}{ endif} 1396 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1397 1398 \textcolor{comment}{!$OMP parallel do default(private) shared(u,v,h,tv,forces,visc,dt,G,GV,US,CS,use\_EOS,dt\_Rho0, &} 1399 \textcolor{comment}{!$OMP h\_neglect,h\_tiny,g\_H\_Rho0,js,je,OBC,Isq,Ieq,nz, &} 1400 \textcolor{comment}{!$OMP U\_bg\_sq,mask\_v,cdrag\_sqrt,cdrag\_sqrt\_Z,Rho0x400\_G,nkml)} 1401 \textcolor{keywordflow}{do} j=js,je \textcolor{comment}{! u-point loop} 1402 \textcolor{keywordflow}{if} (cs%dynamic\_viscous\_ML) \textcolor{keywordflow}{then} 1403 do\_any = .false. 1404 \textcolor{keywordflow}{do} i=isq,ieq 1405 htot(i) = 0.0 1406 \textcolor{keywordflow}{if} (g%mask2dCu(i,j) < 0.5) \textcolor{keywordflow}{then} 1407 do\_i(i) = .false. ; visc%nkml\_visc\_u(i,j) = nkml 1408 \textcolor{keywordflow}{else} 1409 do\_i(i) = .true. ; do\_any = .true. 1410 k\_massive(i) = nkml 1411 thtot(i) = 0.0 ; shtot(i) = 0.0 ; rhtot(i) = 0.0 1412 uhtot(i) = dt\_rho0 * forces%taux(i,j) 1413 vhtot(i) = 0.25 * dt\_rho0 * ((forces%tauy(i,j) + forces%tauy(i+1,j-1)) + & 1414 (forces%tauy(i,j-1) + forces%tauy(i+1,j))) 1415 1416 \textcolor{keywordflow}{if} (cs%omega\_frac >= 1.0) \textcolor{keywordflow}{then} ; absf = 2.0*cs%omega ; \textcolor{keywordflow}{else} 1417 absf = 0.5*(abs(g%CoriolisBu(i,j)) + abs(g%CoriolisBu(i,j-1))) 1418 \textcolor{keywordflow}{if} (cs%omega\_frac > 0.0) & 1419 absf = sqrt(cs%omega\_frac*4.0*cs%omega**2 + (1.0-cs%omega\_frac)*absf**2) 1420 \textcolor{keywordflow}{ endif} 1421 u\_star = max(cs%ustar\_min, 0.5 * (forces%ustar(i,j) + forces%ustar(i+1,j))) 1422 idecay\_len\_tke(i) = ((absf / u\_star) * cs%TKE\_decay) * gv%H\_to\_Z 1423 \textcolor{keywordflow}{ endif} 1424 \textcolor{keywordflow}{ enddo} 1425 1426 \textcolor{keywordflow}{if} (do\_any) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz 1427 \textcolor{keywordflow}{if} (k > nkml) \textcolor{keywordflow}{then} 1428 do\_any = .false. 1429 \textcolor{keywordflow}{if} (use\_eos .and. (k==nkml+1)) \textcolor{keywordflow}{then} 1430 \textcolor{comment}{! Find dRho/dT and dRho\_dS.} 1431 \textcolor{keywordflow}{do} i=isq,ieq 1432 press(i) = (gv%H\_to\_RZ*gv%g\_Earth) * htot(i) 1433 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%p\_surf)) press(i) = press(i) + 0.5*(tv%p\_surf(i,j)+tv%p\_surf(i+1,j)) 1434 k2 = max(1,nkml) 1435 i\_2hlay = 1.0 / (h(i,j,k2) + h(i+1,j,k2) + h\_neglect) 1436 t\_eos(i) = (h(i,j,k2)*tv%T(i,j,k2) + h(i+1,j,k2)*tv%T(i+1,j,k2)) * i\_2hlay 1437 s\_eos(i) = (h(i,j,k2)*tv%S(i,j,k2) + h(i+1,j,k2)*tv%S(i+1,j,k2)) * i\_2hlay 1438 \textcolor{keywordflow}{ enddo} 1439 \textcolor{keyword}{call }calculate\_density\_derivs(t\_eos, s\_eos, press, dr\_dt, dr\_ds, tv%eqn\_of\_state, & 1440 (/isq-g%IsdB+1,ieq-g%IsdB+1/) ) 1441 \textcolor{keywordflow}{ endif} 1442 1443 \textcolor{keywordflow}{do} i=isq,ieq ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1444 1445 hlay = 0.5*(h(i,j,k) + h(i+1,j,k)) 1446 \textcolor{keywordflow}{if} (hlay > h\_tiny) \textcolor{keywordflow}{then} \textcolor{comment}{! Only consider non-vanished layers.} 1447 i\_2hlay = 1.0 / (h(i,j,k) + h(i+1,j,k)) 1448 v\_at\_u = 0.5 * (h(i,j,k) * (v(i,j,k) + v(i,j-1,k)) + & 1449 h(i+1,j,k) * (v(i+1,j,k) + v(i+1,j-1,k))) * i\_2hlay 1450 uh2 = ((uhtot(i) - htot(i)*u(i,j,k))**2 + (vhtot(i) - htot(i)*v\_at\_u)**2) 1451 1452 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1453 t\_lay = (h(i,j,k)*tv%T(i,j,k) + h(i+1,j,k)*tv%T(i+1,j,k)) * i\_2hlay 1454 s\_lay = (h(i,j,k)*tv%S(i,j,k) + h(i+1,j,k)*tv%S(i+1,j,k)) * i\_2hlay 1455 ghprime = g\_h\_rho0 * (dr\_dt(i) * (t\_lay*htot(i) - thtot(i)) + & 1456 dr\_ds(i) * (s\_lay*htot(i) - shtot(i))) 1457 \textcolor{keywordflow}{else} 1458 ghprime = g\_h\_rho0 * (gv%Rlay(k)*htot(i) - rhtot(i)) 1459 \textcolor{keywordflow}{ endif} 1460 1461 \textcolor{keywordflow}{if} (ghprime > 0.0) \textcolor{keywordflow}{then} 1462 ribulk = cs%bulk\_Ri\_ML * exp(-htot(i) * idecay\_len\_tke(i)) 1463 \textcolor{keywordflow}{if} (ribulk * uh2 <= (htot(i)**2) * ghprime) \textcolor{keywordflow}{then} 1464 visc%nkml\_visc\_u(i,j) = \textcolor{keywordtype}{real}(k\_massive(i)) 1465 do\_i(i) = .false. 1466 \textcolor{keywordflow}{elseif} (ribulk * uh2 <= (htot(i) + hlay)**2 * ghprime) \textcolor{keywordflow}{then} 1467 visc%nkml\_visc\_u(i,j) = \textcolor{keywordtype}{real(k-1)} + & 1468 ( sqrt(ribulk * uh2 / ghprime) - htot(i) ) / hlay 1469 do\_i(i) = .false. 1470 \textcolor{keywordflow}{ endif} 1471 \textcolor{keywordflow}{ endif} 1472 k\_massive(i) = k 1473 \textcolor{keywordflow}{ endif} \textcolor{comment}{! hlay > h\_tiny} 1474 1475 \textcolor{keywordflow}{if} (do\_i(i)) do\_any = .true. 1476 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 1477 1478 \textcolor{keywordflow}{if} (.not.do\_any) \textcolor{keywordflow}{exit} \textcolor{comment}{! All columns are done.} 1479 \textcolor{keywordflow}{ endif} 1480 1481 \textcolor{keywordflow}{do} i=isq,ieq ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1482 htot(i) = htot(i) + 0.5 * (h(i,j,k) + h(i+1,j,k)) 1483 uhtot(i) = uhtot(i) + 0.5 * (h(i,j,k) + h(i+1,j,k)) * u(i,j,k) 1484 vhtot(i) = vhtot(i) + 0.25 * (h(i,j,k) * (v(i,j,k) + v(i,j-1,k)) + & 1485 h(i+1,j,k) * (v(i+1,j,k) + v(i+1,j-1,k))) 1486 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1487 thtot(i) = thtot(i) + 0.5 * (h(i,j,k)*tv%T(i,j,k) + h(i+1,j,k)*tv%T(i+1,j,k)) 1488 shtot(i) = shtot(i) + 0.5 * (h(i,j,k)*tv%S(i,j,k) + h(i+1,j,k)*tv%S(i+1,j,k)) 1489 \textcolor{keywordflow}{else} 1490 rhtot(i) = rhtot(i) + 0.5 * (h(i,j,k) + h(i+1,j,k)) * gv%Rlay(k) 1491 \textcolor{keywordflow}{ endif} 1492 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 1493 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1494 1495 \textcolor{keywordflow}{if} (do\_any) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} i=isq,ieq ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1496 visc%nkml\_visc\_u(i,j) = k\_massive(i) 1497 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1498 \textcolor{keywordflow}{ endif} \textcolor{comment}{! dynamic\_viscous\_ML} 1499 1500 do\_any\_shelf = .false. 1501 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(forces%frac\_shelf\_u)) \textcolor{keywordflow}{then} 1502 \textcolor{keywordflow}{do} i=isq,ieq 1503 \textcolor{keywordflow}{if} (forces%frac\_shelf\_u(i,j)*g%mask2dCu(i,j) == 0.0) \textcolor{keywordflow}{then} 1504 do\_i(i) = .false. 1505 visc%tbl\_thick\_shelf\_u(i,j) = 0.0 ; visc%kv\_tbl\_shelf\_u(i,j) = 0.0 1506 \textcolor{keywordflow}{else} 1507 do\_i(i) = .true. ; do\_any\_shelf = .true. 1508 \textcolor{keywordflow}{ endif} 1509 \textcolor{keywordflow}{ enddo} 1510 \textcolor{keywordflow}{ endif} 1511 1512 \textcolor{keywordflow}{if} (do\_any\_shelf) \textcolor{keywordflow}{then} 1513 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isq,ieq ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1514 \textcolor{keywordflow}{if} (u(i,j,k) *(h(i+1,j,k) - h(i,j,k)) >= 0) \textcolor{keywordflow}{then} 1515 h\_at\_vel(i,k) = 2.0*h(i,j,k)*h(i+1,j,k) / & 1516 (h(i,j,k) + h(i+1,j,k) + h\_neglect) 1517 \textcolor{keywordflow}{else} 1518 h\_at\_vel(i,k) = 0.5 * (h(i,j,k) + h(i+1,j,k)) 1519 \textcolor{keywordflow}{ endif} 1520 \textcolor{keywordflow}{else} 1521 h\_at\_vel(i,k) = 0.0 ; ustar(i) = 0.0 1522 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1523 1524 \textcolor{keywordflow}{do} i=isq,ieq ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1525 htot\_vel = 0.0 ; hwtot = 0.0 ; hutot = 0.0 1526 thtot(i) = 0.0 ; shtot(i) = 0.0 1527 \textcolor{keywordflow}{if} (use\_eos .or. .not.cs%linear\_drag) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz 1528 \textcolor{keywordflow}{if} (htot\_vel>=cs%Htbl\_shelf) \textcolor{keywordflow}{exit} \textcolor{comment}{! terminate the k loop} 1529 hweight = min(cs%Htbl\_shelf - htot\_vel, h\_at\_vel(i,k)) 1530 \textcolor{keywordflow}{if} (hweight <= 1.5*gv%Angstrom\_H + h\_neglect) cycle 1531 1532 htot\_vel = htot\_vel + h\_at\_vel(i,k) 1533 hwtot = hwtot + hweight 1534 1535 \textcolor{keywordflow}{if} (.not.cs%linear\_drag) \textcolor{keywordflow}{then} 1536 v\_at\_u = set\_v\_at\_u(v, h, g, i, j, k, mask\_v, obc) 1537 hutot = hutot + hweight * sqrt(u(i,j,k)**2 + & 1538 v\_at\_u**2 + u\_bg\_sq) 1539 \textcolor{keywordflow}{ endif} 1540 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1541 thtot(i) = thtot(i) + hweight * 0.5 * (tv%T(i,j,k) + tv%T(i+1,j,k)) 1542 shtot(i) = shtot(i) + hweight * 0.5 * (tv%S(i,j,k) + tv%S(i+1,j,k)) 1543 \textcolor{keywordflow}{ endif} 1544 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1545 1546 \textcolor{keywordflow}{if} ((.not.cs%linear\_drag) .and. (hwtot > 0.0)) \textcolor{keywordflow}{then} 1547 ustar(i) = cdrag\_sqrt\_z * hutot/hwtot 1548 \textcolor{keywordflow}{else} 1549 ustar(i) = cdrag\_sqrt\_z * cs%drag\_bg\_vel 1550 \textcolor{keywordflow}{ endif} 1551 1552 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (hwtot > 0.0) \textcolor{keywordflow}{then} 1553 t\_eos(i) = thtot(i)/hwtot ; s\_eos(i) = shtot(i)/hwtot 1554 \textcolor{keywordflow}{else} 1555 t\_eos(i) = 0.0 ; s\_eos(i) = 0.0 1556 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 1557 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! I-loop} 1558 1559 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1560 \textcolor{keyword}{call }calculate\_density\_derivs(t\_eos, s\_eos, forces%p\_surf(:,j), dr\_dt, dr\_ds, & 1561 tv%eqn\_of\_state, (/isq-g%IsdB+1,ieq-g%IsdB+1/) ) 1562 \textcolor{keywordflow}{ endif} 1563 1564 \textcolor{keywordflow}{do} i=isq,ieq ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1565 \textcolor{comment}{! The 400.0 in this expression is the square of a constant proposed} 1566 \textcolor{comment}{! by Killworth and Edwards, 1999, in equation (2.20).} 1567 ustarsq = rho0x400\_g * ustar(i)**2 1568 htot(i) = 0.0 1569 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1570 thtot(i) = 0.0 ; shtot(i) = 0.0 1571 \textcolor{keywordflow}{do} k=1,nz-1 1572 \textcolor{keywordflow}{if} (h\_at\_vel(i,k) <= 0.0) cycle 1573 t\_lay = 0.5 * (tv%T(i,j,k) + tv%T(i+1,j,k)) 1574 s\_lay = 0.5 * (tv%S(i,j,k) + tv%S(i+1,j,k)) 1575 oldfn = dr\_dt(i)*(t\_lay*htot(i) - thtot(i)) + dr\_ds(i)*(s\_lay*htot(i) - shtot(i)) 1576 \textcolor{keywordflow}{if} (oldfn >= ustarsq) \textcolor{keywordflow}{exit} 1577 1578 dfn = (dr\_dt(i)*(0.5*(tv%T(i,j,k+1)+tv%T(i+1,j,k+1)) - t\_lay) + & 1579 dr\_ds(i)*(0.5*(tv%S(i,j,k+1)+tv%S(i+1,j,k+1)) - s\_lay)) * & 1580 (h\_at\_vel(i,k)+htot(i)) 1581 \textcolor{keywordflow}{if} ((oldfn + dfn) <= ustarsq) \textcolor{keywordflow}{then} 1582 dh = h\_at\_vel(i,k) 1583 \textcolor{keywordflow}{else} 1584 dh = h\_at\_vel(i,k) * sqrt((ustarsq-oldfn) / (dfn)) 1585 \textcolor{keywordflow}{ endif} 1586 1587 htot(i) = htot(i) + dh 1588 thtot(i) = thtot(i) + t\_lay*dh ; shtot(i) = shtot(i) + s\_lay*dh 1589 \textcolor{keywordflow}{ enddo} 1590 \textcolor{keywordflow}{if} ((oldfn < ustarsq) .and. (h\_at\_vel(i,nz) > 0.0)) \textcolor{keywordflow}{then} 1591 t\_lay = 0.5*(tv%T(i,j,nz) + tv%T(i+1,j,nz)) 1592 s\_lay = 0.5*(tv%S(i,j,nz) + tv%S(i+1,j,nz)) 1593 \textcolor{keywordflow}{if} (dr\_dt(i)*(t\_lay*htot(i) - thtot(i)) + & 1594 dr\_ds(i)*(s\_lay*htot(i) - shtot(i)) < ustarsq) & 1595 htot(i) = htot(i) + h\_at\_vel(i,nz) 1596 \textcolor{keywordflow}{ endif} \textcolor{comment}{! Examination of layer nz.} 1597 \textcolor{keywordflow}{else} \textcolor{comment}{! Use Rlay as the density variable.} 1598 rhtot = 0.0 1599 \textcolor{keywordflow}{do} k=1,nz-1 1600 rlay = gv%Rlay(k) ; rlb = gv%Rlay(k+1) 1601 1602 oldfn = rlay*htot(i) - rhtot(i) 1603 \textcolor{keywordflow}{if} (oldfn >= ustarsq) \textcolor{keywordflow}{exit} 1604 1605 dfn = (rlb - rlay)*(h\_at\_vel(i,k)+htot(i)) 1606 \textcolor{keywordflow}{if} ((oldfn + dfn) <= ustarsq) \textcolor{keywordflow}{then} 1607 dh = h\_at\_vel(i,k) 1608 \textcolor{keywordflow}{else} 1609 dh = h\_at\_vel(i,k) * sqrt((ustarsq-oldfn) / (dfn)) 1610 \textcolor{keywordflow}{ endif} 1611 1612 htot(i) = htot(i) + dh 1613 rhtot(i) = rhtot(i) + rlay*dh 1614 \textcolor{keywordflow}{ enddo} 1615 \textcolor{keywordflow}{if} (gv%Rlay(nz)*htot(i) - rhtot(i) < ustarsq) & 1616 htot(i) = htot(i) + h\_at\_vel(i,nz) 1617 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 1618 1619 \textcolor{comment}{!visc%tbl\_thick\_shelf\_u(I,j) = GV%H\_to\_Z * max(CS%Htbl\_shelf\_min, &} 1620 \textcolor{comment}{! htot(I) / (0.5 + sqrt(0.25 + &} 1621 \textcolor{comment}{! (htot(i)*(G%CoriolisBu(I,J-1)+G%CoriolisBu(I,J)))**2 / &} 1622 \textcolor{comment}{! (ustar(i)*GV%Z\_to\_H)**2 )) )} 1623 ustar1 = ustar(i)*gv%Z\_to\_H 1624 h2f2 = (htot(i)*(g%CoriolisBu(i,j-1)+g%CoriolisBu(i,j)) + h\_neglect*cs%omega)**2 1625 tbl\_thick\_z = gv%H\_to\_Z * max(cs%Htbl\_shelf\_min, & 1626 ( htot(i)*ustar1 ) / ( 0.5*ustar1 + sqrt((0.5*ustar1)**2 + h2f2 ) ) ) 1627 visc%tbl\_thick\_shelf\_u(i,j) = tbl\_thick\_z 1628 visc%Kv\_tbl\_shelf\_u(i,j) = max(cs%Kv\_TBL\_min, cdrag\_sqrt*ustar(i)*tbl\_thick\_z) 1629 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! I-loop} 1630 \textcolor{keywordflow}{ endif} \textcolor{comment}{! do\_any\_shelf} 1631 1632 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! j-loop at u-points} 1633 1634 \textcolor{comment}{!$OMP parallel do default(private) shared(u,v,h,tv,forces,visc,dt,G,GV,US,CS,use\_EOS,dt\_Rho0, &} 1635 \textcolor{comment}{!$OMP h\_neglect,h\_tiny,g\_H\_Rho0,is,ie,OBC,Jsq,Jeq,nz, &} 1636 \textcolor{comment}{!$OMP U\_bg\_sq,cdrag\_sqrt,cdrag\_sqrt\_Z,Rho0x400\_G,nkml,mask\_u)} 1637 \textcolor{keywordflow}{do} j=jsq,jeq \textcolor{comment}{! v-point loop} 1638 \textcolor{keywordflow}{if} (cs%dynamic\_viscous\_ML) \textcolor{keywordflow}{then} 1639 do\_any = .false. 1640 \textcolor{keywordflow}{do} i=is,ie 1641 htot(i) = 0.0 1642 \textcolor{keywordflow}{if} (g%mask2dCv(i,j) < 0.5) \textcolor{keywordflow}{then} 1643 do\_i(i) = .false. ; visc%nkml\_visc\_v(i,j) = nkml 1644 \textcolor{keywordflow}{else} 1645 do\_i(i) = .true. ; do\_any = .true. 1646 k\_massive(i) = nkml 1647 thtot(i) = 0.0 ; shtot(i) = 0.0 ; rhtot(i) = 0.0 1648 vhtot(i) = dt\_rho0 * forces%tauy(i,j) 1649 uhtot(i) = 0.25 * dt\_rho0 * ((forces%taux(i,j) + forces%taux(i-1,j+1)) + & 1650 (forces%taux(i-1,j) + forces%taux(i,j+1))) 1651 1652 \textcolor{keywordflow}{if} (cs%omega\_frac >= 1.0) \textcolor{keywordflow}{then} ; absf = 2.0*cs%omega ; \textcolor{keywordflow}{else} 1653 absf = 0.5*(abs(g%CoriolisBu(i-1,j)) + abs(g%CoriolisBu(i,j))) 1654 \textcolor{keywordflow}{if} (cs%omega\_frac > 0.0) & 1655 absf = sqrt(cs%omega\_frac*4.0*cs%omega**2 + (1.0-cs%omega\_frac)*absf**2) 1656 \textcolor{keywordflow}{ endif} 1657 1658 u\_star = max(cs%ustar\_min, 0.5 * (forces%ustar(i,j) + forces%ustar(i,j+1))) 1659 idecay\_len\_tke(i) = ((absf / u\_star) * cs%TKE\_decay) * gv%H\_to\_Z 1660 1661 \textcolor{keywordflow}{ endif} 1662 \textcolor{keywordflow}{ enddo} 1663 1664 \textcolor{keywordflow}{if} (do\_any) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz 1665 \textcolor{keywordflow}{if} (k > nkml) \textcolor{keywordflow}{then} 1666 do\_any = .false. 1667 \textcolor{keywordflow}{if} (use\_eos .and. (k==nkml+1)) \textcolor{keywordflow}{then} 1668 \textcolor{comment}{! Find dRho/dT and dRho\_dS.} 1669 \textcolor{keywordflow}{do} i=is,ie 1670 press(i) = (gv%H\_to\_RZ * gv%g\_Earth) * htot(i) 1671 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%p\_surf)) press(i) = press(i) + 0.5*(tv%p\_surf(i,j)+tv%p\_surf(i,j+1)) 1672 k2 = max(1,nkml) 1673 i\_2hlay = 1.0 / (h(i,j,k2) + h(i,j+1,k2) + h\_neglect) 1674 t\_eos(i) = (h(i,j,k2)*tv%T(i,j,k2) + h(i,j+1,k2)*tv%T(i,j+1,k2)) * i\_2hlay 1675 s\_eos(i) = (h(i,j,k2)*tv%S(i,j,k2) + h(i,j+1,k2)*tv%S(i,j+1,k2)) * i\_2hlay 1676 \textcolor{keywordflow}{ enddo} 1677 \textcolor{keyword}{call }calculate\_density\_derivs(t\_eos, s\_eos, press, dr\_dt, dr\_ds, & 1678 tv%eqn\_of\_state, (/is-g%IsdB+1,ie-g%IsdB+1/) ) 1679 \textcolor{keywordflow}{ endif} 1680 1681 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1682 1683 hlay = 0.5*(h(i,j,k) + h(i,j+1,k)) 1684 \textcolor{keywordflow}{if} (hlay > h\_tiny) \textcolor{keywordflow}{then} \textcolor{comment}{! Only consider non-vanished layers.} 1685 i\_2hlay = 1.0 / (h(i,j,k) + h(i,j+1,k)) 1686 u\_at\_v = 0.5 * (h(i,j,k) * (u(i-1,j,k) + u(i,j,k)) + & 1687 h(i,j+1,k) * (u(i-1,j+1,k) + u(i,j+1,k))) * i\_2hlay 1688 uh2 = ((uhtot(i) - htot(i)*u\_at\_v)**2 + (vhtot(i) - htot(i)*v(i,j,k))**2) 1689 1690 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1691 t\_lay = (h(i,j,k)*tv%T(i,j,k) + h(i,j+1,k)*tv%T(i,j+1,k)) * i\_2hlay 1692 s\_lay = (h(i,j,k)*tv%S(i,j,k) + h(i,j+1,k)*tv%S(i,j+1,k)) * i\_2hlay 1693 ghprime = g\_h\_rho0 * (dr\_dt(i) * (t\_lay*htot(i) - thtot(i)) + & 1694 dr\_ds(i) * (s\_lay*htot(i) - shtot(i))) 1695 \textcolor{keywordflow}{else} 1696 ghprime = g\_h\_rho0 * (gv%Rlay(k)*htot(i) - rhtot(i)) 1697 \textcolor{keywordflow}{ endif} 1698 1699 \textcolor{keywordflow}{if} (ghprime > 0.0) \textcolor{keywordflow}{then} 1700 ribulk = cs%bulk\_Ri\_ML * exp(-htot(i) * idecay\_len\_tke(i)) 1701 \textcolor{keywordflow}{if} (ribulk * uh2 <= htot(i)**2 * ghprime) \textcolor{keywordflow}{then} 1702 visc%nkml\_visc\_v(i,j) = \textcolor{keywordtype}{real}(k\_massive(i)) 1703 do\_i(i) = .false. 1704 \textcolor{keywordflow}{elseif} (ribulk * uh2 <= (htot(i) + hlay)**2 * ghprime) \textcolor{keywordflow}{then} 1705 visc%nkml\_visc\_v(i,j) = \textcolor{keywordtype}{real(k-1)} + & 1706 ( sqrt(ribulk * uh2 / ghprime) - htot(i) ) / hlay 1707 do\_i(i) = .false. 1708 \textcolor{keywordflow}{ endif} 1709 \textcolor{keywordflow}{ endif} 1710 k\_massive(i) = k 1711 \textcolor{keywordflow}{ endif} \textcolor{comment}{! hlay > h\_tiny} 1712 1713 \textcolor{keywordflow}{if} (do\_i(i)) do\_any = .true. 1714 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 1715 1716 \textcolor{keywordflow}{if} (.not.do\_any) \textcolor{keywordflow}{exit} \textcolor{comment}{! All columns are done.} 1717 \textcolor{keywordflow}{ endif} 1718 1719 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1720 htot(i) = htot(i) + 0.5 * (h(i,j,k) + h(i,j+1,k)) 1721 vhtot(i) = vhtot(i) + 0.5 * (h(i,j,k) + h(i,j+1,k)) * v(i,j,k) 1722 uhtot(i) = uhtot(i) + 0.25 * (h(i,j,k) * (u(i-1,j,k) + u(i,j,k)) + & 1723 h(i,j+1,k) * (u(i-1,j+1,k) + u(i,j+1,k))) 1724 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1725 thtot(i) = thtot(i) + 0.5 * (h(i,j,k)*tv%T(i,j,k) + h(i,j+1,k)*tv%T(i,j+1,k)) 1726 shtot(i) = shtot(i) + 0.5 * (h(i,j,k)*tv%S(i,j,k) + h(i,j+1,k)*tv%S(i,j+1,k)) 1727 \textcolor{keywordflow}{else} 1728 rhtot(i) = rhtot(i) + 0.5 * (h(i,j,k) + h(i,j+1,k)) * gv%Rlay(k) 1729 \textcolor{keywordflow}{ endif} 1730 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 1731 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1732 1733 \textcolor{keywordflow}{if} (do\_any) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1734 visc%nkml\_visc\_v(i,j) = k\_massive(i) 1735 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1736 \textcolor{keywordflow}{ endif} \textcolor{comment}{! dynamic\_viscous\_ML} 1737 1738 do\_any\_shelf = .false. 1739 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(forces%frac\_shelf\_v)) \textcolor{keywordflow}{then} 1740 \textcolor{keywordflow}{do} i=is,ie 1741 \textcolor{keywordflow}{if} (forces%frac\_shelf\_v(i,j)*g%mask2dCv(i,j) == 0.0) \textcolor{keywordflow}{then} 1742 do\_i(i) = .false. 1743 visc%tbl\_thick\_shelf\_v(i,j) = 0.0 ; visc%kv\_tbl\_shelf\_v(i,j) = 0.0 1744 \textcolor{keywordflow}{else} 1745 do\_i(i) = .true. ; do\_any\_shelf = .true. 1746 \textcolor{keywordflow}{ endif} 1747 \textcolor{keywordflow}{ enddo} 1748 \textcolor{keywordflow}{ endif} 1749 1750 \textcolor{keywordflow}{if} (do\_any\_shelf) \textcolor{keywordflow}{then} 1751 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1752 \textcolor{keywordflow}{if} (v(i,j,k) * (h(i,j+1,k) - h(i,j,k)) >= 0) \textcolor{keywordflow}{then} 1753 h\_at\_vel(i,k) = 2.0*h(i,j,k)*h(i,j+1,k) / & 1754 (h(i,j,k) + h(i,j+1,k) + h\_neglect) 1755 \textcolor{keywordflow}{else} 1756 h\_at\_vel(i,k) = 0.5 * (h(i,j,k) + h(i,j+1,k)) 1757 \textcolor{keywordflow}{ endif} 1758 \textcolor{keywordflow}{else} 1759 h\_at\_vel(i,k) = 0.0 ; ustar(i) = 0.0 1760 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 1761 1762 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1763 htot\_vel = 0.0 ; hwtot = 0.0 ; hutot = 0.0 1764 thtot(i) = 0.0 ; shtot(i) = 0.0 1765 \textcolor{keywordflow}{if} (use\_eos .or. .not.cs%linear\_drag) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz 1766 \textcolor{keywordflow}{if} (htot\_vel>=cs%Htbl\_shelf) \textcolor{keywordflow}{exit} \textcolor{comment}{! terminate the k loop} 1767 hweight = min(cs%Htbl\_shelf - htot\_vel, h\_at\_vel(i,k)) 1768 \textcolor{keywordflow}{if} (hweight <= 1.5*gv%Angstrom\_H + h\_neglect) cycle 1769 1770 htot\_vel = htot\_vel + h\_at\_vel(i,k) 1771 hwtot = hwtot + hweight 1772 1773 \textcolor{keywordflow}{if} (.not.cs%linear\_drag) \textcolor{keywordflow}{then} 1774 u\_at\_v = set\_u\_at\_v(u, h, g, i, j, k, mask\_u, obc) 1775 hutot = hutot + hweight * sqrt(v(i,j,k)**2 + & 1776 u\_at\_v**2 + u\_bg\_sq) 1777 \textcolor{keywordflow}{ endif} 1778 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1779 thtot(i) = thtot(i) + hweight * 0.5 * (tv%T(i,j,k) + tv%T(i,j+1,k)) 1780 shtot(i) = shtot(i) + hweight * 0.5 * (tv%S(i,j,k) + tv%S(i,j+1,k)) 1781 \textcolor{keywordflow}{ endif} 1782 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1783 1784 \textcolor{keywordflow}{if} (.not.cs%linear\_drag) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (hwtot > 0.0) \textcolor{keywordflow}{then} 1785 ustar(i) = cdrag\_sqrt\_z * hutot/hwtot 1786 \textcolor{keywordflow}{else} 1787 ustar(i) = cdrag\_sqrt\_z * cs%drag\_bg\_vel 1788 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 1789 1790 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (hwtot > 0.0) \textcolor{keywordflow}{then} 1791 t\_eos(i) = thtot(i)/hwtot ; s\_eos(i) = shtot(i)/hwtot 1792 \textcolor{keywordflow}{else} 1793 t\_eos(i) = 0.0 ; s\_eos(i) = 0.0 1794 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 1795 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! I-loop} 1796 1797 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1798 \textcolor{keyword}{call }calculate\_density\_derivs(t\_eos, s\_eos, forces%p\_surf(:,j), dr\_dt, dr\_ds, & 1799 tv%eqn\_of\_state, (/is-g%IsdB+1,ie-g%IsdB+1/) ) 1800 \textcolor{keywordflow}{ endif} 1801 1802 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (do\_i(i)) \textcolor{keywordflow}{then} 1803 \textcolor{comment}{! The 400.0 in this expression is the square of a constant proposed} 1804 \textcolor{comment}{! by Killworth and Edwards, 1999, in equation (2.20).} 1805 ustarsq = rho0x400\_g * ustar(i)**2 1806 htot(i) = 0.0 1807 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 1808 thtot(i) = 0.0 ; shtot(i) = 0.0 1809 \textcolor{keywordflow}{do} k=1,nz-1 1810 \textcolor{keywordflow}{if} (h\_at\_vel(i,k) <= 0.0) cycle 1811 t\_lay = 0.5 * (tv%T(i,j,k) + tv%T(i,j+1,k)) 1812 s\_lay = 0.5 * (tv%S(i,j,k) + tv%S(i,j+1,k)) 1813 oldfn = dr\_dt(i)*(t\_lay*htot(i) - thtot(i)) + dr\_ds(i)*(s\_lay*htot(i) - shtot(i)) 1814 \textcolor{keywordflow}{if} (oldfn >= ustarsq) \textcolor{keywordflow}{exit} 1815 1816 dfn = (dr\_dt(i)*(0.5*(tv%T(i,j,k+1)+tv%T(i,j+1,k+1)) - t\_lay) + & 1817 dr\_ds(i)*(0.5*(tv%S(i,j,k+1)+tv%S(i,j+1,k+1)) - s\_lay)) * & 1818 (h\_at\_vel(i,k)+htot(i)) 1819 \textcolor{keywordflow}{if} ((oldfn + dfn) <= ustarsq) \textcolor{keywordflow}{then} 1820 dh = h\_at\_vel(i,k) 1821 \textcolor{keywordflow}{else} 1822 dh = h\_at\_vel(i,k) * sqrt((ustarsq-oldfn) / (dfn)) 1823 \textcolor{keywordflow}{ endif} 1824 1825 htot(i) = htot(i) + dh 1826 thtot(i) = thtot(i) + t\_lay*dh ; shtot(i) = shtot(i) + s\_lay*dh 1827 \textcolor{keywordflow}{ enddo} 1828 \textcolor{keywordflow}{if} ((oldfn < ustarsq) .and. (h\_at\_vel(i,nz) > 0.0)) \textcolor{keywordflow}{then} 1829 t\_lay = 0.5*(tv%T(i,j,nz) + tv%T(i,j+1,nz)) 1830 s\_lay = 0.5*(tv%S(i,j,nz) + tv%S(i,j+1,nz)) 1831 \textcolor{keywordflow}{if} (dr\_dt(i)*(t\_lay*htot(i) - thtot(i)) + & 1832 dr\_ds(i)*(s\_lay*htot(i) - shtot(i)) < ustarsq) & 1833 htot(i) = htot(i) + h\_at\_vel(i,nz) 1834 \textcolor{keywordflow}{ endif} \textcolor{comment}{! Examination of layer nz.} 1835 \textcolor{keywordflow}{else} \textcolor{comment}{! Use Rlay as the density variable.} 1836 rhtot = 0.0 1837 \textcolor{keywordflow}{do} k=1,nz-1 1838 rlay = gv%Rlay(k) ; rlb = gv%Rlay(k+1) 1839 1840 oldfn = rlay*htot(i) - rhtot(i) 1841 \textcolor{keywordflow}{if} (oldfn >= ustarsq) \textcolor{keywordflow}{exit} 1842 1843 dfn = (rlb - rlay)*(h\_at\_vel(i,k)+htot(i)) 1844 \textcolor{keywordflow}{if} ((oldfn + dfn) <= ustarsq) \textcolor{keywordflow}{then} 1845 dh = h\_at\_vel(i,k) 1846 \textcolor{keywordflow}{else} 1847 dh = h\_at\_vel(i,k) * sqrt((ustarsq-oldfn) / (dfn)) 1848 \textcolor{keywordflow}{ endif} 1849 1850 htot(i) = htot(i) + dh 1851 rhtot = rhtot + rlay*dh 1852 \textcolor{keywordflow}{ enddo} 1853 \textcolor{keywordflow}{if} (gv%Rlay(nz)*htot(i) - rhtot(i) < ustarsq) & 1854 htot(i) = htot(i) + h\_at\_vel(i,nz) 1855 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 1856 1857 \textcolor{comment}{!visc%tbl\_thick\_shelf\_v(i,J) = GV%H\_to\_Z * max(CS%Htbl\_shelf\_min, &} 1858 \textcolor{comment}{! htot(i) / (0.5 + sqrt(0.25 + &} 1859 \textcolor{comment}{! (htot(i)*(G%CoriolisBu(I-1,J)+G%CoriolisBu(I,J)))**2 / &} 1860 \textcolor{comment}{! (ustar(i)*GV%Z\_to\_H)**2 )) )} 1861 ustar1 = ustar(i)*gv%Z\_to\_H 1862 h2f2 = (htot(i)*(g%CoriolisBu(i-1,j)+g%CoriolisBu(i,j)) + h\_neglect*cs%omega)**2 1863 tbl\_thick\_z = gv%H\_to\_Z * max(cs%Htbl\_shelf\_min, & 1864 ( htot(i)*ustar1 ) / ( 0.5*ustar1 + sqrt((0.5*ustar1)**2 + h2f2 ) ) ) 1865 visc%tbl\_thick\_shelf\_v(i,j) = tbl\_thick\_z 1866 visc%Kv\_tbl\_shelf\_v(i,j) = max(cs%Kv\_TBL\_min, cdrag\_sqrt*ustar(i)*tbl\_thick\_z) 1867 1868 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! i-loop} 1869 \textcolor{keywordflow}{ endif} \textcolor{comment}{! do\_any\_shelf} 1870 1871 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! J-loop at v-points} 1872 1873 \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then} 1874 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(visc%nkml\_visc\_u) .and. \textcolor{keyword}{associated}(visc%nkml\_visc\_v)) & 1875 \textcolor{keyword}{call }uvchksum(\textcolor{stringliteral}{"nkml\_visc\_[uv]"}, visc%nkml\_visc\_u, visc%nkml\_visc\_v, & 1876 g%HI, haloshift=0, scalar\_pair=.true.) 1877 \textcolor{keywordflow}{ endif} 1878 \textcolor{keywordflow}{if} (cs%id\_nkml\_visc\_u > 0) & 1879 \textcolor{keyword}{call }post\_data(cs%id\_nkml\_visc\_u, visc%nkml\_visc\_u, cs%diag) 1880 \textcolor{keywordflow}{if} (cs%id\_nkml\_visc\_v > 0) & 1881 \textcolor{keyword}{call }post\_data(cs%id\_nkml\_visc\_v, visc%nkml\_visc\_v, cs%diag) 1882 \end{DoxyCode}