\hypertarget{namespacemom__kappa__shear}{}\section{mom\+\_\+kappa\+\_\+shear Module Reference} \label{namespacemom__kappa__shear}\index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsection{Detailed Description} Shear-\/dependent mixing following Jackson et al. 2008. By Laura Jackson and Robert Hallberg, 2006-\/2008. This file contains the subroutines that determine the diapycnal diffusivity driven by resolved shears, as specified by the parameterizations described in Jackson and Hallberg (J\+PO, 2008). The technique by which the 6 equations (for kappa, T\+KE, u, v, T, and S) are solved simultaneously has been dramatically revised from the previous version. The previous version was not converging in some cases, especially near the surface mixed layer, while the revised version does. The revised version solves for kappa and T\+KE with shear and stratification fixed, then marches the density and velocities forward with an adaptive (and aggressive) time step in a predictor-\/corrector-\/corrector emulation of a trapezoidal scheme. Run-\/time-\/settable parameters determine the tolerence to which the kappa and T\+KE equations are solved and the minimum time step that can be taken. \subsection*{Data Types} \begin{DoxyCompactItemize} \item type \mbox{\hyperlink{structmom__kappa__shear_1_1kappa__shear__cs}{kappa\+\_\+shear\+\_\+cs}} \begin{DoxyCompactList}\small\item\em This control structure holds the parameters that regulate shear mixing. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection*{Functions/\+Subroutines} \begin{DoxyCompactItemize} \item subroutine, public \mbox{\hyperlink{namespacemom__kappa__shear_a3f00b08e1174690d40c0c2065fa9a8b1}{calculate\+\_\+kappa\+\_\+shear}} (u\+\_\+in, v\+\_\+in, h, tv, p\+\_\+surf, kappa\+\_\+io, tke\+\_\+io, kv\+\_\+io, dt, G, GV, US, CS, initialize\+\_\+all) \begin{DoxyCompactList}\small\item\em Subroutine for calculating shear-\/driven diffusivity and T\+KE in tracer columns. \end{DoxyCompactList}\item subroutine, public \mbox{\hyperlink{namespacemom__kappa__shear_a2d8e291656bab5f83179523c4bea4d85}{calc\+\_\+kappa\+\_\+shear\+\_\+vertex}} (u\+\_\+in, v\+\_\+in, h, T\+\_\+in, S\+\_\+in, tv, p\+\_\+surf, kappa\+\_\+io, tke\+\_\+io, kv\+\_\+io, dt, G, GV, US, CS, initialize\+\_\+all) \begin{DoxyCompactList}\small\item\em Subroutine for calculating shear-\/driven diffusivity and T\+KE in corner columns. \end{DoxyCompactList}\item subroutine \mbox{\hyperlink{namespacemom__kappa__shear_a26cc5bb15545f04cfaf07e53410e09ec}{kappa\+\_\+shear\+\_\+column}} (kappa, tke, dt, nzc, f2, surface\+\_\+pres, dz, u0xdz, v0xdz, T0xdz, S0xdz, kappa\+\_\+avg, tke\+\_\+avg, tv, CS, GV, US, I\+\_\+\+Ld2\+\_\+1d, dz\+\_\+\+Int\+\_\+1d) \begin{DoxyCompactList}\small\item\em This subroutine calculates shear-\/driven diffusivity and T\+KE in a single column. \end{DoxyCompactList}\item subroutine \mbox{\hyperlink{namespacemom__kappa__shear_a0b931b0b834d887e321eb6eb1924fa9a}{calculate\+\_\+projected\+\_\+state}} (kappa, u0, v0, T0, S0, dt, nz, dz, I\+\_\+dz\+\_\+int, dbuoy\+\_\+dT, dbuoy\+\_\+dS, u, v, T, Sal, GV, US, N2, S2, ks\+\_\+int, ke\+\_\+int, vel\+\_\+underflow) \begin{DoxyCompactList}\small\item\em This subroutine calculates the velocities, temperature and salinity that the water column will have after mixing for dt with diffusivities kappa. It may also calculate the projected buoyancy frequency and shear. \end{DoxyCompactList}\item subroutine \mbox{\hyperlink{namespacemom__kappa__shear_a351d44e4fe5cfb5852d019a0c1e66100}{find\+\_\+kappa\+\_\+tke}} (N2, S2, kappa\+\_\+in, Idz, dz\+\_\+\+Int, I\+\_\+\+L2\+\_\+bdry, f2, nz, CS, GV, US, K\+\_\+Q, tke, kappa, kappa\+\_\+src, local\+\_\+src) \begin{DoxyCompactList}\small\item\em This subroutine calculates new, consistent estimates of T\+KE and kappa. \end{DoxyCompactList}\item logical function, public \mbox{\hyperlink{namespacemom__kappa__shear_a82428168ba463cf7276ae96d3e32475c}{kappa\+\_\+shear\+\_\+init}} (Time, G, GV, US, param\+\_\+file, diag, CS) \begin{DoxyCompactList}\small\item\em This subroutine initializes the parameters that regulate shear-\/driven mixing. \end{DoxyCompactList}\item logical function, public \mbox{\hyperlink{namespacemom__kappa__shear_ac7859c609e462000ca8fd763d68d141e}{kappa\+\_\+shear\+\_\+is\+\_\+used}} (param\+\_\+file) \begin{DoxyCompactList}\small\item\em This function indicates to other modules whether the Jackson et al shear mixing parameterization will be used without needing to duplicate the log entry. \end{DoxyCompactList}\item logical function, public \mbox{\hyperlink{namespacemom__kappa__shear_ad4d87b0928aea195213e682b493eb555}{kappa\+\_\+shear\+\_\+at\+\_\+vertex}} (param\+\_\+file) \begin{DoxyCompactList}\small\item\em This function indicates to other modules whether the Jackson et al shear mixing parameterization will be used at the vertices without needing to duplicate the log entry. It returns false if the Jackson et al scheme is not used or if it is used via calculations at the tracer points. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection{Function/\+Subroutine Documentation} \mbox{\Hypertarget{namespacemom__kappa__shear_a2d8e291656bab5f83179523c4bea4d85}\label{namespacemom__kappa__shear_a2d8e291656bab5f83179523c4bea4d85}} \index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}!calc\+\_\+kappa\+\_\+shear\+\_\+vertex@{calc\+\_\+kappa\+\_\+shear\+\_\+vertex}} \index{calc\+\_\+kappa\+\_\+shear\+\_\+vertex@{calc\+\_\+kappa\+\_\+shear\+\_\+vertex}!mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsubsection{\texorpdfstring{calc\+\_\+kappa\+\_\+shear\+\_\+vertex()}{calc\_kappa\_shear\_vertex()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+kappa\+\_\+shear\+::calc\+\_\+kappa\+\_\+shear\+\_\+vertex (\begin{DoxyParamCaption}\item[{real, dimension(szib\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{u\+\_\+in, }\item[{real, dimension(szi\+\_\+(g),szjb\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{v\+\_\+in, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{h, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{T\+\_\+in, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{S\+\_\+in, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{real, dimension(\+:,\+:), pointer}]{p\+\_\+surf, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)+1), intent(out)}]{kappa\+\_\+io, }\item[{real, dimension(szib\+\_\+(g),szjb\+\_\+(g),szk\+\_\+(gv)+1), intent(out)}]{tke\+\_\+io, }\item[{real, dimension(szib\+\_\+(g),szjb\+\_\+(g),szk\+\_\+(gv)+1), intent(inout)}]{kv\+\_\+io, }\item[{real, intent(in)}]{dt, }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US, }\item[{type(\mbox{\hyperlink{structmom__kappa__shear_1_1kappa__shear__cs}{kappa\+\_\+shear\+\_\+cs}}), pointer}]{CS, }\item[{logical, intent(in), optional}]{initialize\+\_\+all }\end{DoxyParamCaption})} Subroutine for calculating shear-\/driven diffusivity and T\+KE in corner columns. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & The ocean\textquotesingle{}s grid structure.\\ \hline \mbox{\tt in} & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em u\+\_\+in} & Initial zonal velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em v\+\_\+in} & Initial 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 t\+\_\+in} & Layer potential temperatures \mbox{[}degC\mbox{]}\\ \hline \mbox{\tt in} & {\em s\+\_\+in} & Layer salinities in ppt.\\ \hline \mbox{\tt in} & {\em tv} & A structure containing pointers to any available thermodynamic fields. Absent fields have N\+U\+LL ptrs.\\ \hline & {\em p\+\_\+surf} & The pressure at the ocean surface \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} (or N\+U\+LL).\\ \hline \mbox{\tt out} & {\em kappa\+\_\+io} & The diapycnal diffusivity at each interface\\ \hline \mbox{\tt out} & {\em tke\+\_\+io} & The turbulent kinetic energy per unit mass at\\ \hline \mbox{\tt in,out} & {\em kv\+\_\+io} & The vertical viscosity at each interface \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \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 kappa\+\_\+shear\+\_\+init.\\ \hline \mbox{\tt in} & {\em initialize\+\_\+all} & If present and false, the previous value of kappa is used to start the iterations \\ \hline \end{DoxyParams} Definition at line 344 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 344 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{ !< The ocean's grid structure.} 345 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< The ocean's vertical grid structure.} 346 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 347 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(GV))}, & 348 \textcolor{keywordtype}{intent(in)} :: u\_in\textcolor{comment}{ !< Initial zonal velocity [L T-1 ~> m s-1].} 349 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(GV))}, & 350 \textcolor{keywordtype}{intent(in)} :: v\_in\textcolor{comment}{ !< Initial meridional velocity [L T-1 ~> m s-1].} 351 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))}, & 352 \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2].} 353 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))}, & 354 \textcolor{keywordtype}{intent(in)} :: T\_in\textcolor{comment}{ !< Layer potential temperatures [degC]} 355 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))}, & 356 \textcolor{keywordtype}{intent(in)} :: S\_in\textcolor{comment}{ !< Layer salinities in ppt.} 357 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< A structure containing pointers to any} 358 \textcolor{comment}{ !! available thermodynamic fields. Absent fields} 359 \textcolor{comment}{ !! have NULL ptrs.} 360 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{pointer} :: p\_surf\textcolor{comment}{ !< The pressure at the ocean surface [R L2 T-2 ~> Pa]} 361 \textcolor{comment}{ !! (or NULL).} 362 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV)+1)}, & 363 \textcolor{keywordtype}{intent(out)} :: kappa\_io\textcolor{comment}{ !< The diapycnal diffusivity at each interface} 364 \textcolor{comment}{ !! (not layer!) [Z2 T-1 ~> m2 s-1].} 365 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJB\_(G),SZK\_(GV)+1)}, & 366 \textcolor{keywordtype}{intent(out)} :: tke\_io\textcolor{comment}{ !< The turbulent kinetic energy per unit mass at} 367 \textcolor{comment}{ !! each interface (not layer!) [Z2 T-2 ~> m2 s-2].} 368 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJB\_(G),SZK\_(GV)+1)}, & 369 \textcolor{keywordtype}{intent(inout)} :: kv\_io\textcolor{comment}{ !< The vertical viscosity at each interface [Z2 T-1 ~> m2 s-1].} 370 \textcolor{comment}{ !! The previous value is used to initialize kappa} 371 \textcolor{comment}{ !! in the vertex columes as Kappa = Kv/Prandtl} 372 \textcolor{comment}{ !! to accelerate the iteration toward covergence.} 373 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: dt\textcolor{comment}{ !< Time increment [T ~> s].} 374 \textcolor{keywordtype}{type}(Kappa\_shear\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< The control structure returned by a previous} 375 \textcolor{comment}{ !! call to kappa\_shear\_init.} 376 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: initialize\_all\textcolor{comment}{ !< If present and false, the previous} 377 \textcolor{comment}{ !! value of kappa is used to start the iterations} 378 379 \textcolor{comment}{! Local variables} 380 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZK\_(GV))} :: & 381 h\_2d, & \textcolor{comment}{! A 2-D version of h, but converted to [Z ~> m].} 382 u\_2d, v\_2d, & \textcolor{comment}{! 2-D versions of u\_in and v\_in, converted to [L T-1 ~> m s-1].} 383 T\_2d, S\_2d, rho\_2d \textcolor{comment}{! 2-D versions of T [degC], S [ppt], and rho [R ~> kg m-3].} 384 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZK\_(GV)+1,2)} :: & 385 kappa\_2d \textcolor{comment}{! Quasi 2-D versions of kappa\_io [Z2 T-1 ~> m2 s-1].} 386 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZK\_(GV)+1)} :: & 387 tke\_2d \textcolor{comment}{! 2-D version tke\_io [Z2 T-2 ~> m2 s-2].} 388 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))} :: & 389 Idz, & \textcolor{comment}{! The inverse of the distance between TKE points [Z-1 ~> m-1].} 390 dz, & \textcolor{comment}{! The layer thickness [Z ~> m].} 391 u0xdz, & \textcolor{comment}{! The initial zonal velocity times dz [L Z T-1 ~> m2 s-1].} 392 v0xdz, & \textcolor{comment}{! The initial meridional velocity times dz [L Z T-1 ~> m2 s-1].} 393 T0xdz, & \textcolor{comment}{! The initial temperature times dz [degC Z ~> degC m].} 394 S0xdz \textcolor{comment}{! The initial salinity times dz [ppt Z ~> ppt m].} 395 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: & 396 kappa, & \textcolor{comment}{! The shear-driven diapycnal diffusivity at an interface [Z2 T-1 ~> m2 s-1].} 397 tke, & \textcolor{comment}{! The Turbulent Kinetic Energy per unit mass at an interface [Z2 T-2 ~> m2 s-2].} 398 kappa\_avg, & \textcolor{comment}{! The time-weighted average of kappa [Z2 T-1 ~> m2 s-1].} 399 tke\_avg \textcolor{comment}{! The time-weighted average of TKE [Z2 T-2 ~> m2 s-2].} 400 \textcolor{keywordtype}{real} :: f2 \textcolor{comment}{! The squared Coriolis parameter of each column [T-2 ~> s-2].} 401 \textcolor{keywordtype}{real} :: surface\_pres \textcolor{comment}{! The top surface pressure [R L2 T-2 ~> Pa].} 402 403 \textcolor{keywordtype}{real} :: dz\_in\_lay \textcolor{comment}{! The running sum of the thickness in a layer [Z ~> m].} 404 \textcolor{keywordtype}{real} :: k0dt \textcolor{comment}{! The background diffusivity times the timestep [Z2 ~> m2].} 405 \textcolor{keywordtype}{real} :: dz\_massless \textcolor{comment}{! A layer thickness that is considered massless [Z ~> m].} 406 \textcolor{keywordtype}{real} :: I\_hwt \textcolor{comment}{! The inverse of the masked thickness weights [H-1 ~> m-1 or m2 kg-1].} 407 \textcolor{keywordtype}{real} :: I\_Prandtl \textcolor{comment}{! The inverse of the turbulent Prandtl number [nondim].} 408 \textcolor{keywordtype}{logical} :: use\_temperature \textcolor{comment}{! If true, temperature and salinity have been} 409 \textcolor{comment}{! allocated and are being used as state variables.} 410 \textcolor{keywordtype}{logical} :: new\_kappa = .true. \textcolor{comment}{! If true, ignore the value of kappa from the} 411 \textcolor{comment}{! last call to this subroutine.} 412 \textcolor{keywordtype}{logical} :: do\_i \textcolor{comment}{! If true, work on this column.} 413 414 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: kc \textcolor{comment}{! The index map between the original} 415 \textcolor{comment}{! interfaces and the interfaces with massless layers} 416 \textcolor{comment}{! merged into nearby massive layers.} 417 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: kf \textcolor{comment}{! The fractional weight of interface kc+1 for} 418 \textcolor{comment}{! interpolating back to the original index space [nondim].} 419 \textcolor{keywordtype}{integer} :: IsB, IeB, JsB, JeB, i, j, k, nz, nzc, J2, J2m1 420 421 \textcolor{comment}{! Diagnostics that should be deleted?} 422 isb = g%isc-1 ; ieb = g%iecB ; jsb = g%jsc-1 ; jeb = g%jecB ; nz = gv%ke 423 424 use\_temperature = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%T)) use\_temperature = .true. 425 new\_kappa = .true. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(initialize\_all)) new\_kappa = initialize\_all 426 427 k0dt = dt*cs%kappa\_0 428 dz\_massless = 0.1*sqrt(k0dt) 429 i\_prandtl = 0.0 ; \textcolor{keywordflow}{if} (cs%Prandtl\_turb > 0.0) i\_prandtl = 1.0 / cs%Prandtl\_turb 430 431 \textcolor{comment}{!$OMP parallel do default(private) shared(jsB,jeB,isB,ieB,nz,h,u\_in,v\_in,use\_temperature,new\_kappa, &} 432 \textcolor{comment}{!$OMP tv,G,GV,US,CS,kappa\_io,dz\_massless,k0dt,p\_surf,dt,tke\_io,kv\_io,I\_Prandtl)} 433 \textcolor{keywordflow}{do} j=jsb,jeb 434 j2 = mod(j,2)+1 ; j2m1 = 3-j2 \textcolor{comment}{! = mod(J-1,2)+1} 435 436 \textcolor{comment}{! Interpolate the various quantities to the corners, using masks.} 437 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isb,ieb 438 u\_2d(i,k) = (u\_in(i,j,k) * (g%mask2dCu(i,j) * (h(i,j,k) + h(i+1,j,k))) + & 439 u\_in(i,j+1,k) * (g%mask2dCu(i,j+1) * (h(i,j+1,k) + h(i+1,j+1,k))) ) / & 440 ((g%mask2dCu(i,j) * (h(i,j,k) + h(i+1,j,k)) + & 441 g%mask2dCu(i,j+1) * (h(i,j+1,k) + h(i+1,j+1,k))) + gv%H\_subroundoff) 442 v\_2d(i,k) = (v\_in(i,j,k) * (g%mask2dCv(i,j) * (h(i,j,k) + h(i,j+1,k))) + & 443 v\_in(i+1,j,k) * (g%mask2dCv(i+1,j) * (h(i+1,j,k) + h(i+1,j+1,k))) ) / & 444 ((g%mask2dCv(i,j) * (h(i,j,k) + h(i,j+1,k)) + & 445 g%mask2dCv(i+1,j) * (h(i+1,j,k) + h(i+1,j+1,k))) + gv%H\_subroundoff) 446 i\_hwt = 1.0 / (((g%mask2dT(i,j) * h(i,j,k) + g%mask2dT(i+1,j+1) * h(i+1,j+1,k)) + & 447 (g%mask2dT(i+1,j) * h(i+1,j,k) + g%mask2dT(i,j+1) * h(i,j+1,k))) + & 448 gv%H\_subroundoff) 449 \textcolor{keywordflow}{if} (use\_temperature) \textcolor{keywordflow}{then} 450 t\_2d(i,k) = ( ((g%mask2dT(i,j) * h(i,j,k)) * t\_in(i,j,k) + & 451 (g%mask2dT(i+1,j+1) * h(i+1,j+1,k)) * t\_in(i+1,j+1,k)) + & 452 ((g%mask2dT(i+1,j) * h(i+1,j,k)) * t\_in(i+1,j,k) + & 453 (g%mask2dT(i,j+1) * h(i,j+1,k)) * t\_in(i,j+1,k)) ) * i\_hwt 454 s\_2d(i,k) = ( ((g%mask2dT(i,j) * h(i,j,k)) * s\_in(i,j,k) + & 455 (g%mask2dT(i+1,j+1) * h(i+1,j+1,k)) * s\_in(i+1,j+1,k)) + & 456 ((g%mask2dT(i+1,j) * h(i+1,j,k)) * s\_in(i+1,j,k) + & 457 (g%mask2dT(i,j+1) * h(i,j+1,k)) * s\_in(i,j+1,k)) ) * i\_hwt 458 \textcolor{keywordflow}{ endif} 459 h\_2d(i,k) = gv%H\_to\_Z * ((g%mask2dT(i,j) * h(i,j,k) + g%mask2dT(i+1,j+1) * h(i+1,j+1,k)) + & 460 (g%mask2dT(i+1,j) * h(i+1,j,k) + g%mask2dT(i,j+1) * h(i,j+1,k)) ) / & 461 ((g%mask2dT(i,j) + g%mask2dT(i+1,j+1)) + & 462 (g%mask2dT(i+1,j) + g%mask2dT(i,j+1)) + 1.0e-36 ) 463 \textcolor{comment}{! h\_2d(I,k) = 0.25*((h(i,j,k) + h(i+1,j+1,k)) + (h(i+1,j,k) + h(i,j+1,k)))*GV%H\_to\_Z} 464 \textcolor{comment}{! h\_2d(I,k) = ((h(i,j,k)**2 + h(i+1,j+1,k)**2) + &} 465 \textcolor{comment}{! (h(i+1,j,k)**2 + h(i,j+1,k)**2))*GV%H\_to\_Z * I\_hwt} 466 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 467 \textcolor{keywordflow}{if} (.not.use\_temperature) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isb,ieb 468 rho\_2d(i,k) = gv%Rlay(k) 469 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 470 \textcolor{keywordflow}{if} (.not.new\_kappa) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=isb,ieb 471 kappa\_2d(i,k,j2) = kv\_io(i,j,k) * i\_prandtl 472 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 473 474 \textcolor{comment}{!---------------------------------------} 475 \textcolor{comment}{! Work on each column.} 476 \textcolor{comment}{!---------------------------------------} 477 \textcolor{keywordflow}{do} i=isb,ieb ; \textcolor{keywordflow}{if} ((g%mask2dCu(i,j) + g%mask2dCu(i,j+1)) + & 478 (g%mask2dCv(i,j) + g%mask2dCv(i+1,j)) > 0.0) \textcolor{keywordflow}{then} 479 \textcolor{comment}{! call cpu\_clock\_begin(Id\_clock\_setup)} 480 \textcolor{comment}{! Store a transposed version of the initial arrays.} 481 \textcolor{comment}{! Any elimination of massless layers would occur here.} 482 \textcolor{keywordflow}{if} (cs%eliminate\_massless) \textcolor{keywordflow}{then} 483 nzc = 1 484 \textcolor{keywordflow}{do} k=1,nz 485 \textcolor{comment}{! Zero out the thicknesses of all layers, even if they are unused.} 486 dz(k) = 0.0 ; u0xdz(k) = 0.0 ; v0xdz(k) = 0.0 487 t0xdz(k) = 0.0 ; s0xdz(k) = 0.0 488 489 \textcolor{comment}{! Add a new layer if this one has mass.} 490 \textcolor{comment}{! if ((dz(nzc) > 0.0) .and. (h\_2d(I,k) > dz\_massless)) nzc = nzc+1} 491 \textcolor{keywordflow}{if} ((k>cs%nkml) .and. (dz(nzc) > 0.0) .and. & 492 (h\_2d(i,k) > dz\_massless)) nzc = nzc+1 493 494 \textcolor{comment}{! Only merge clusters of massless layers.} 495 \textcolor{comment}{! if ((dz(nzc) > dz\_massless) .or. &} 496 \textcolor{comment}{! ((dz(nzc) > 0.0) .and. (h\_2d(I,k) > dz\_massless))) nzc = nzc+1} 497 498 kc(k) = nzc 499 dz(nzc) = dz(nzc) + h\_2d(i,k) 500 u0xdz(nzc) = u0xdz(nzc) + u\_2d(i,k)*h\_2d(i,k) 501 v0xdz(nzc) = v0xdz(nzc) + v\_2d(i,k)*h\_2d(i,k) 502 \textcolor{keywordflow}{if} (use\_temperature) \textcolor{keywordflow}{then} 503 t0xdz(nzc) = t0xdz(nzc) + t\_2d(i,k)*h\_2d(i,k) 504 s0xdz(nzc) = s0xdz(nzc) + s\_2d(i,k)*h\_2d(i,k) 505 \textcolor{keywordflow}{else} 506 t0xdz(nzc) = t0xdz(nzc) + rho\_2d(i,k)*h\_2d(i,k) 507 s0xdz(nzc) = s0xdz(nzc) + rho\_2d(i,k)*h\_2d(i,k) 508 \textcolor{keywordflow}{ endif} 509 \textcolor{keywordflow}{ enddo} 510 kc(nz+1) = nzc+1 511 512 \textcolor{comment}{! Set up Idz as the inverse of layer thicknesses.} 513 \textcolor{keywordflow}{do} k=1,nzc ; idz(k) = 1.0 / dz(k) ;\textcolor{keywordflow}{ enddo} 514 515 \textcolor{comment}{! Now determine kf, the fractional weight of interface kc when} 516 \textcolor{comment}{! interpolating between interfaces kc and kc+1.} 517 kf(1) = 0.0 ; dz\_in\_lay = h\_2d(i,1) 518 \textcolor{keywordflow}{do} k=2,nz 519 \textcolor{keywordflow}{if} (kc(k) > kc(k-1)) \textcolor{keywordflow}{then} 520 kf(k) = 0.0 ; dz\_in\_lay = h\_2d(i,k) 521 \textcolor{keywordflow}{else} 522 kf(k) = dz\_in\_lay*idz(kc(k)) ; dz\_in\_lay = dz\_in\_lay + h\_2d(i,k) 523 \textcolor{keywordflow}{ endif} 524 \textcolor{keywordflow}{ enddo} 525 kf(nz+1) = 0.0 526 \textcolor{keywordflow}{else} 527 \textcolor{keywordflow}{do} k=1,nz 528 dz(k) = h\_2d(i,k) 529 u0xdz(k) = u\_2d(i,k)*dz(k) ; v0xdz(k) = v\_2d(i,k)*dz(k) 530 \textcolor{keywordflow}{ enddo} 531 \textcolor{keywordflow}{if} (use\_temperature) \textcolor{keywordflow}{then} 532 \textcolor{keywordflow}{do} k=1,nz 533 t0xdz(k) = t\_2d(i,k)*dz(k) ; s0xdz(k) = s\_2d(i,k)*dz(k) 534 \textcolor{keywordflow}{ enddo} 535 \textcolor{keywordflow}{else} 536 \textcolor{keywordflow}{do} k=1,nz 537 t0xdz(k) = rho\_2d(i,k)*dz(k) ; s0xdz(k) = rho\_2d(i,k)*dz(k) 538 \textcolor{keywordflow}{ enddo} 539 \textcolor{keywordflow}{ endif} 540 nzc = nz 541 \textcolor{keywordflow}{do} k=1,nzc+1 ; kc(k) = k ; kf(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 542 \textcolor{keywordflow}{ endif} 543 f2 = g%CoriolisBu(i,j)**2 544 surface\_pres = 0.0 545 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(p\_surf)) \textcolor{keywordflow}{then} 546 \textcolor{keywordflow}{if} (cs%psurf\_bug) \textcolor{keywordflow}{then} 547 \textcolor{comment}{! This is wrong because it is averaging values from land in some places.} 548 surface\_pres = 0.25 * ((p\_surf(i,j) + p\_surf(i+1,j+1)) + & 549 (p\_surf(i+1,j) + p\_surf(i,j+1))) 550 \textcolor{keywordflow}{else} 551 surface\_pres = ((g%mask2dT(i,j) * p\_surf(i,j) + g%mask2dT(i+1,j+1) * p\_surf(i+1,j+1)) + & 552 (g%mask2dT(i+1,j) * p\_surf(i+1,j) + g%mask2dT(i,j+1) * p\_surf(i,j+1)) ) / & 553 ((g%mask2dT(i,j) + g%mask2dT(i+1,j+1)) + & 554 (g%mask2dT(i+1,j) + g%mask2dT(i,j+1)) + 1.0e-36 ) 555 \textcolor{keywordflow}{ endif} 556 \textcolor{keywordflow}{ endif} 557 558 \textcolor{comment}{! ----------------------------------------------------} 559 \textcolor{comment}{! Set the initial guess for kappa, here defined at interfaces.} 560 \textcolor{comment}{! ----------------------------------------------------} 561 \textcolor{keywordflow}{if} (new\_kappa) \textcolor{keywordflow}{then} 562 \textcolor{keywordflow}{do} k=1,nzc+1 ; kappa(k) = us%m2\_s\_to\_Z2\_T*1.0 ;\textcolor{keywordflow}{ enddo} 563 \textcolor{keywordflow}{else} 564 \textcolor{keywordflow}{do} k=1,nzc+1 ; kappa(k) = kappa\_2d(i,k,j2) ;\textcolor{keywordflow}{ enddo} 565 \textcolor{keywordflow}{ endif} 566 567 \textcolor{keyword}{call }kappa\_shear\_column(kappa, tke, dt, nzc, f2, surface\_pres, & 568 dz, u0xdz, v0xdz, t0xdz, s0xdz, kappa\_avg, & 569 tke\_avg, tv, cs, gv, us) 570 \textcolor{comment}{! call cpu\_clock\_begin(Id\_clock\_setup)} 571 \textcolor{comment}{! Extrapolate from the vertically reduced grid back to the original layers.} 572 \textcolor{keywordflow}{if} (nz == nzc) \textcolor{keywordflow}{then} 573 \textcolor{keywordflow}{do} k=1,nz+1 574 kappa\_2d(i,k,j2) = kappa\_avg(k) 575 \textcolor{keywordflow}{if} (cs%all\_layer\_TKE\_bug) \textcolor{keywordflow}{then} 576 tke\_2d(i,k) = tke(k) 577 \textcolor{keywordflow}{else} 578 tke\_2d(i,k) = tke\_avg(k) 579 \textcolor{keywordflow}{ endif} 580 \textcolor{keywordflow}{ enddo} 581 \textcolor{keywordflow}{else} 582 \textcolor{keywordflow}{do} k=1,nz+1 583 \textcolor{keywordflow}{if} (kf(k) == 0.0) \textcolor{keywordflow}{then} 584 kappa\_2d(i,k,j2) = kappa\_avg(kc(k)) 585 tke\_2d(i,k) = tke\_avg(kc(k)) 586 \textcolor{keywordflow}{else} 587 kappa\_2d(i,k,j2) = (1.0-kf(k)) * kappa\_avg(kc(k)) + kf(k) * kappa\_avg(kc(k)+1) 588 tke\_2d(i,k) = (1.0-kf(k)) * tke\_avg(kc(k)) + kf(k) * tke\_avg(kc(k)+1) 589 \textcolor{keywordflow}{ endif} 590 \textcolor{keywordflow}{ enddo} 591 \textcolor{keywordflow}{ endif} 592 \textcolor{comment}{! call cpu\_clock\_end(Id\_clock\_setup)} 593 \textcolor{keywordflow}{else} \textcolor{comment}{! Land points, still inside the i-loop.} 594 \textcolor{keywordflow}{do} k=1,nz+1 595 kappa\_2d(i,k,j2) = 0.0 ; tke\_2d(i,k) = 0.0 596 \textcolor{keywordflow}{ enddo} 597 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! i-loop} 598 599 \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=isb,ieb 600 tke\_io(i,j,k) = g%mask2dBu(i,j) * tke\_2d(i,k) 601 kv\_io(i,j,k) = ( g%mask2dBu(i,j) * kappa\_2d(i,k,j2) ) * cs%Prandtl\_turb 602 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 603 \textcolor{keywordflow}{if} (j>=g%jsc) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=g%isc,g%iec 604 \textcolor{comment}{! Set the diffusivities in tracer columns from the values at vertices.} 605 kappa\_io(i,j,k) = g%mask2dT(i,j) * 0.25 * & 606 ((kappa\_2d(i-1,k,j2m1) + kappa\_2d(i,k,j2)) + & 607 (kappa\_2d(i-1,k,j2) + kappa\_2d(i,k,j2m1))) 608 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 609 610 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of J-loop} 611 612 \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then} 613 \textcolor{keyword}{call }hchksum(kappa\_io, \textcolor{stringliteral}{"kappa"}, g%HI, scale=us%Z2\_T\_to\_m2\_s) 614 \textcolor{keyword}{call }bchksum(tke\_io, \textcolor{stringliteral}{"tke"}, g%HI, scale=us%Z\_to\_m**2*us%s\_to\_T**2) 615 \textcolor{keywordflow}{ endif} 616 617 \textcolor{keywordflow}{if} (cs%id\_Kd\_shear > 0) \textcolor{keyword}{call }post\_data(cs%id\_Kd\_shear, kappa\_io, cs%diag) 618 \textcolor{keywordflow}{if} (cs%id\_TKE > 0) \textcolor{keyword}{call }post\_data(cs%id\_TKE, tke\_io, cs%diag) 619 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__kappa__shear_a3f00b08e1174690d40c0c2065fa9a8b1}\label{namespacemom__kappa__shear_a3f00b08e1174690d40c0c2065fa9a8b1}} \index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}!calculate\+\_\+kappa\+\_\+shear@{calculate\+\_\+kappa\+\_\+shear}} \index{calculate\+\_\+kappa\+\_\+shear@{calculate\+\_\+kappa\+\_\+shear}!mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsubsection{\texorpdfstring{calculate\+\_\+kappa\+\_\+shear()}{calculate\_kappa\_shear()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+kappa\+\_\+shear\+::calculate\+\_\+kappa\+\_\+shear (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{u\+\_\+in, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{v\+\_\+in, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)), intent(in)}]{h, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{real, dimension(\+:,\+:), pointer}]{p\+\_\+surf, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)+1), intent(inout)}]{kappa\+\_\+io, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)+1), intent(out)}]{tke\+\_\+io, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(gv)+1), intent(inout)}]{kv\+\_\+io, }\item[{real, intent(in)}]{dt, }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US, }\item[{type(\mbox{\hyperlink{structmom__kappa__shear_1_1kappa__shear__cs}{kappa\+\_\+shear\+\_\+cs}}), pointer}]{CS, }\item[{logical, intent(in), optional}]{initialize\+\_\+all }\end{DoxyParamCaption})} Subroutine for calculating shear-\/driven diffusivity and T\+KE in tracer columns. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & The ocean\textquotesingle{}s grid structure.\\ \hline \mbox{\tt in} & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em u\+\_\+in} & Initial zonal velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em v\+\_\+in} & Initial 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 & {\em p\+\_\+surf} & The pressure at the ocean surface \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]} (or N\+U\+LL).\\ \hline \mbox{\tt in,out} & {\em kappa\+\_\+io} & The diapycnal diffusivity at each interface\\ \hline \mbox{\tt out} & {\em tke\+\_\+io} & The turbulent kinetic energy per unit mass at\\ \hline \mbox{\tt in,out} & {\em kv\+\_\+io} & The vertical viscosity at each interface\\ \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 kappa\+\_\+shear\+\_\+init.\\ \hline \mbox{\tt in} & {\em initialize\+\_\+all} & If present and false, the previous value of kappa is used to start the iterations \\ \hline \end{DoxyParams} Definition at line 111 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 111 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{ !< The ocean's grid structure.} 112 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< The ocean's vertical grid structure.} 113 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 114 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))}, & 115 \textcolor{keywordtype}{intent(in)} :: u\_in\textcolor{comment}{ !< Initial zonal velocity [L T-1 ~> m s-1].} 116 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))}, & 117 \textcolor{keywordtype}{intent(in)} :: v\_in\textcolor{comment}{ !< Initial meridional velocity [L T-1 ~> m s-1].} 118 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV))}, & 119 \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2].} 120 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< A structure containing pointers to any} 121 \textcolor{comment}{ !! available thermodynamic fields. Absent fields} 122 \textcolor{comment}{ !! have NULL ptrs.} 123 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{pointer} :: p\_surf\textcolor{comment}{ !< The pressure at the ocean surface [R L2 T-2 ~> Pa] (or NULL).} 124 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV)+1)}, & 125 \textcolor{keywordtype}{intent(inout)} :: kappa\_io\textcolor{comment}{ !< The diapycnal diffusivity at each interface} 126 \textcolor{comment}{ !! (not layer!) [Z2 T-1 ~> m2 s-1]. Initially this is the} 127 \textcolor{comment}{ !! value from the previous timestep, which may} 128 \textcolor{comment}{ !! accelerate the iteration toward convergence.} 129 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV)+1)}, & 130 \textcolor{keywordtype}{intent(out)} :: tke\_io\textcolor{comment}{ !< The turbulent kinetic energy per unit mass at} 131 \textcolor{comment}{ !! each interface (not layer!) [Z2 T-2 ~> m2 s-2].} 132 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(GV)+1)}, & 133 \textcolor{keywordtype}{intent(inout)} :: kv\_io\textcolor{comment}{ !< The vertical viscosity at each interface} 134 \textcolor{comment}{ !! (not layer!) [Z2 T-1 ~> m2 s-1]. This discards any} 135 \textcolor{comment}{ !! previous value (i.e. it is intent out) and} 136 \textcolor{comment}{ !! simply sets Kv = Prandtl * Kd\_shear} 137 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: dt\textcolor{comment}{ !< Time increment [T ~> s].} 138 \textcolor{keywordtype}{type}(Kappa\_shear\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< The control structure returned by a previous} 139 \textcolor{comment}{ !! call to kappa\_shear\_init.} 140 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: initialize\_all\textcolor{comment}{ !< If present and false, the previous} 141 \textcolor{comment}{ !! value of kappa is used to start the iterations} 142 143 \textcolor{comment}{! Local variables} 144 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(GV))} :: & 145 h\_2d, & \textcolor{comment}{! A 2-D version of h, but converted to [Z ~> m].} 146 u\_2d, v\_2d, & \textcolor{comment}{! 2-D versions of u\_in and v\_in, converted to [L T-1 ~> m s-1].} 147 T\_2d, S\_2d, rho\_2d \textcolor{comment}{! 2-D versions of T [degC], S [ppt], and rho [R ~> kg m-3].} 148 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZK\_(GV)+1)} :: & 149 kappa\_2d, & \textcolor{comment}{! 2-D version of kappa\_io [Z2 T-1 ~> m2 s-1].} 150 tke\_2d \textcolor{comment}{! 2-D version tke\_io [Z2 T-2 ~> m2 s-2].} 151 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))} :: & 152 Idz, & \textcolor{comment}{! The inverse of the distance between TKE points [Z-1 ~> m-1].} 153 dz, & \textcolor{comment}{! The layer thickness [Z ~> m].} 154 u0xdz, & \textcolor{comment}{! The initial zonal velocity times dz [Z L T-1 ~> m2 s-1].} 155 v0xdz, & \textcolor{comment}{! The initial meridional velocity times dz [Z L T-1 ~> m2 s-1].} 156 T0xdz, & \textcolor{comment}{! The initial temperature times dz [degC Z ~> degC m].} 157 S0xdz \textcolor{comment}{! The initial salinity times dz [ppt Z ~> ppt m].} 158 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: & 159 kappa, & \textcolor{comment}{! The shear-driven diapycnal diffusivity at an interface [Z2 T-1 ~> m2 s-1].} 160 tke, & \textcolor{comment}{! The Turbulent Kinetic Energy per unit mass at an interface [Z2 T-2 ~> m2 s-2].} 161 kappa\_avg, & \textcolor{comment}{! The time-weighted average of kappa [Z2 T-1 ~> m2 s-1].} 162 tke\_avg \textcolor{comment}{! The time-weighted average of TKE [Z2 T-2 ~> m2 s-2].} 163 \textcolor{keywordtype}{real} :: f2 \textcolor{comment}{! The squared Coriolis parameter of each column [T-2 ~> s-2].} 164 \textcolor{keywordtype}{real} :: surface\_pres \textcolor{comment}{! The top surface pressure [R L2 T-2 ~> Pa].} 165 166 \textcolor{keywordtype}{real} :: dz\_in\_lay \textcolor{comment}{! The running sum of the thickness in a layer [Z ~> m].} 167 \textcolor{keywordtype}{real} :: k0dt \textcolor{comment}{! The background diffusivity times the timestep [Z2 ~> m2].} 168 \textcolor{keywordtype}{real} :: dz\_massless \textcolor{comment}{! A layer thickness that is considered massless [Z ~> m].} 169 \textcolor{keywordtype}{logical} :: use\_temperature \textcolor{comment}{! If true, temperature and salinity have been} 170 \textcolor{comment}{! allocated and are being used as state variables.} 171 \textcolor{keywordtype}{logical} :: new\_kappa = .true. \textcolor{comment}{! If true, ignore the value of kappa from the} 172 \textcolor{comment}{! last call to this subroutine.} 173 174 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: kc \textcolor{comment}{! The index map between the original} 175 \textcolor{comment}{! interfaces and the interfaces with massless layers} 176 \textcolor{comment}{! merged into nearby massive layers.} 177 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)} :: kf \textcolor{comment}{! The fractional weight of interface kc+1 for} 178 \textcolor{comment}{! interpolating back to the original index space [nondim].} 179 \textcolor{keywordtype}{integer} :: is, ie, js, je, i, j, k, nz, nzc 180 181 is = g%isc ; ie = g%iec; js = g%jsc ; je = g%jec ; nz = gv%ke 182 183 use\_temperature = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%T)) use\_temperature = .true. 184 new\_kappa = .true. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(initialize\_all)) new\_kappa = initialize\_all 185 186 k0dt = dt*cs%kappa\_0 187 dz\_massless = 0.1*sqrt(k0dt) 188 189 \textcolor{comment}{!$OMP parallel do default(private) shared(js,je,is,ie,nz,h,u\_in,v\_in,use\_temperature,new\_kappa, &} 190 \textcolor{comment}{!$OMP tv,G,GV,US,CS,kappa\_io,dz\_massless,k0dt,p\_surf,dt,tke\_io,kv\_io)} 191 \textcolor{keywordflow}{do} j=js,je 192 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 193 h\_2d(i,k) = h(i,j,k)*gv%H\_to\_Z 194 u\_2d(i,k) = u\_in(i,j,k) ; v\_2d(i,k) = v\_in(i,j,k) 195 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 196 \textcolor{keywordflow}{if} (use\_temperature) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 197 t\_2d(i,k) = tv%T(i,j,k) ; s\_2d(i,k) = tv%S(i,j,k) 198 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ; \textcolor{keywordflow}{else} ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is,ie 199 rho\_2d(i,k) = gv%Rlay(k) \textcolor{comment}{! Could be tv%Rho(i,j,k) ?} 200 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 201 \textcolor{keywordflow}{if} (.not.new\_kappa) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie 202 kappa\_2d(i,k) = kappa\_io(i,j,k) 203 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 204 205 \textcolor{comment}{!---------------------------------------} 206 \textcolor{comment}{! Work on each column.} 207 \textcolor{comment}{!---------------------------------------} 208 \textcolor{keywordflow}{do} i=is,ie ; \textcolor{keywordflow}{if} (g%mask2dT(i,j) > 0.5) \textcolor{keywordflow}{then} 209 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_setup)} 210 \textcolor{comment}{! Store a transposed version of the initial arrays.} 211 \textcolor{comment}{! Any elimination of massless layers would occur here.} 212 \textcolor{keywordflow}{if} (cs%eliminate\_massless) \textcolor{keywordflow}{then} 213 nzc = 1 214 \textcolor{keywordflow}{do} k=1,nz 215 \textcolor{comment}{! Zero out the thicknesses of all layers, even if they are unused.} 216 dz(k) = 0.0 ; u0xdz(k) = 0.0 ; v0xdz(k) = 0.0 217 t0xdz(k) = 0.0 ; s0xdz(k) = 0.0 218 219 \textcolor{comment}{! Add a new layer if this one has mass.} 220 \textcolor{comment}{! if ((dz(nzc) > 0.0) .and. (h\_2d(i,k) > dz\_massless)) nzc = nzc+1} 221 \textcolor{keywordflow}{if} ((k>cs%nkml) .and. (dz(nzc) > 0.0) .and. & 222 (h\_2d(i,k) > dz\_massless)) nzc = nzc+1 223 224 \textcolor{comment}{! Only merge clusters of massless layers.} 225 \textcolor{comment}{! if ((dz(nzc) > dz\_massless) .or. &} 226 \textcolor{comment}{! ((dz(nzc) > 0.0) .and. (h\_2d(i,k) > dz\_massless))) nzc = nzc+1} 227 228 kc(k) = nzc 229 dz(nzc) = dz(nzc) + h\_2d(i,k) 230 u0xdz(nzc) = u0xdz(nzc) + u\_2d(i,k)*h\_2d(i,k) 231 v0xdz(nzc) = v0xdz(nzc) + v\_2d(i,k)*h\_2d(i,k) 232 \textcolor{keywordflow}{if} (use\_temperature) \textcolor{keywordflow}{then} 233 t0xdz(nzc) = t0xdz(nzc) + t\_2d(i,k)*h\_2d(i,k) 234 s0xdz(nzc) = s0xdz(nzc) + s\_2d(i,k)*h\_2d(i,k) 235 \textcolor{keywordflow}{else} 236 t0xdz(nzc) = t0xdz(nzc) + rho\_2d(i,k)*h\_2d(i,k) 237 s0xdz(nzc) = s0xdz(nzc) + rho\_2d(i,k)*h\_2d(i,k) 238 \textcolor{keywordflow}{ endif} 239 \textcolor{keywordflow}{ enddo} 240 kc(nz+1) = nzc+1 241 242 \textcolor{comment}{! Set up Idz as the inverse of layer thicknesses.} 243 \textcolor{keywordflow}{do} k=1,nzc ; idz(k) = 1.0 / dz(k) ;\textcolor{keywordflow}{ enddo} 244 245 \textcolor{comment}{! Now determine kf, the fractional weight of interface kc when} 246 \textcolor{comment}{! interpolating between interfaces kc and kc+1.} 247 kf(1) = 0.0 ; dz\_in\_lay = h\_2d(i,1) 248 \textcolor{keywordflow}{do} k=2,nz 249 \textcolor{keywordflow}{if} (kc(k) > kc(k-1)) \textcolor{keywordflow}{then} 250 kf(k) = 0.0 ; dz\_in\_lay = h\_2d(i,k) 251 \textcolor{keywordflow}{else} 252 kf(k) = dz\_in\_lay*idz(kc(k)) ; dz\_in\_lay = dz\_in\_lay + h\_2d(i,k) 253 \textcolor{keywordflow}{ endif} 254 \textcolor{keywordflow}{ enddo} 255 kf(nz+1) = 0.0 256 \textcolor{keywordflow}{else} 257 \textcolor{keywordflow}{do} k=1,nz 258 dz(k) = h\_2d(i,k) 259 u0xdz(k) = u\_2d(i,k)*dz(k) ; v0xdz(k) = v\_2d(i,k)*dz(k) 260 \textcolor{keywordflow}{ enddo} 261 \textcolor{keywordflow}{if} (use\_temperature) \textcolor{keywordflow}{then} 262 \textcolor{keywordflow}{do} k=1,nz 263 t0xdz(k) = t\_2d(i,k)*dz(k) ; s0xdz(k) = s\_2d(i,k)*dz(k) 264 \textcolor{keywordflow}{ enddo} 265 \textcolor{keywordflow}{else} 266 \textcolor{keywordflow}{do} k=1,nz 267 t0xdz(k) = rho\_2d(i,k)*dz(k) ; s0xdz(k) = rho\_2d(i,k)*dz(k) 268 \textcolor{keywordflow}{ enddo} 269 \textcolor{keywordflow}{ endif} 270 nzc = nz 271 \textcolor{keywordflow}{do} k=1,nzc+1 ; kc(k) = k ; kf(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 272 \textcolor{keywordflow}{ endif} 273 f2 = 0.25 * ((g%CoriolisBu(i,j)**2 + g%CoriolisBu(i-1,j-1)**2) + & 274 (g%CoriolisBu(i,j-1)**2 + g%CoriolisBu(i-1,j)**2)) 275 surface\_pres = 0.0 ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(p\_surf)) surface\_pres = p\_surf(i,j) 276 277 \textcolor{comment}{! ---------------------------------------------------- I\_Ld2\_1d, dz\_Int\_1d} 278 279 \textcolor{comment}{! Set the initial guess for kappa, here defined at interfaces.} 280 \textcolor{comment}{! ----------------------------------------------------} 281 \textcolor{keywordflow}{if} (new\_kappa) \textcolor{keywordflow}{then} 282 \textcolor{keywordflow}{do} k=1,nzc+1 ; kappa(k) = us%m2\_s\_to\_Z2\_T*1.0 ;\textcolor{keywordflow}{ enddo} 283 \textcolor{keywordflow}{else} 284 \textcolor{keywordflow}{do} k=1,nzc+1 ; kappa(k) = kappa\_2d(i,k) ;\textcolor{keywordflow}{ enddo} 285 \textcolor{keywordflow}{ endif} 286 287 \textcolor{keyword}{call }kappa\_shear\_column(kappa, tke, dt, nzc, f2, surface\_pres, & 288 dz, u0xdz, v0xdz, t0xdz, s0xdz, kappa\_avg, & 289 tke\_avg, tv, cs, gv, us) 290 291 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_setup)} 292 \textcolor{comment}{! Extrapolate from the vertically reduced grid back to the original layers.} 293 \textcolor{keywordflow}{if} (nz == nzc) \textcolor{keywordflow}{then} 294 \textcolor{keywordflow}{do} k=1,nz+1 295 kappa\_2d(i,k) = kappa\_avg(k) 296 \textcolor{keywordflow}{if} (cs%all\_layer\_TKE\_bug) \textcolor{keywordflow}{then} 297 tke\_2d(i,k) = tke(k) 298 \textcolor{keywordflow}{else} 299 tke\_2d(i,k) = tke\_avg(k) 300 \textcolor{keywordflow}{ endif} 301 \textcolor{keywordflow}{ enddo} 302 \textcolor{keywordflow}{else} 303 \textcolor{keywordflow}{do} k=1,nz+1 304 \textcolor{keywordflow}{if} (kf(k) == 0.0) \textcolor{keywordflow}{then} 305 kappa\_2d(i,k) = kappa\_avg(kc(k)) 306 tke\_2d(i,k) = tke\_avg(kc(k)) 307 \textcolor{keywordflow}{else} 308 kappa\_2d(i,k) = (1.0-kf(k)) * kappa\_avg(kc(k)) + & 309 kf(k) * kappa\_avg(kc(k)+1) 310 tke\_2d(i,k) = (1.0-kf(k)) * tke\_avg(kc(k)) + & 311 kf(k) * tke\_avg(kc(k)+1) 312 \textcolor{keywordflow}{ endif} 313 \textcolor{keywordflow}{ enddo} 314 \textcolor{keywordflow}{ endif} 315 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_setup)} 316 \textcolor{keywordflow}{else} \textcolor{comment}{! Land points, still inside the i-loop.} 317 \textcolor{keywordflow}{do} k=1,nz+1 318 kappa\_2d(i,k) = 0.0 ; tke\_2d(i,k) = 0.0 319 \textcolor{keywordflow}{ enddo} 320 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! i-loop} 321 322 \textcolor{keywordflow}{do} k=1,nz+1 ; \textcolor{keywordflow}{do} i=is,ie 323 kappa\_io(i,j,k) = g%mask2dT(i,j) * kappa\_2d(i,k) 324 tke\_io(i,j,k) = g%mask2dT(i,j) * tke\_2d(i,k) 325 kv\_io(i,j,k) = ( g%mask2dT(i,j) * kappa\_2d(i,k) ) * cs%Prandtl\_turb 326 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 327 328 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of j-loop} 329 330 \textcolor{keywordflow}{if} (cs%debug) \textcolor{keywordflow}{then} 331 \textcolor{keyword}{call }hchksum(kappa\_io, \textcolor{stringliteral}{"kappa"}, g%HI, scale=us%Z2\_T\_to\_m2\_s) 332 \textcolor{keyword}{call }hchksum(tke\_io, \textcolor{stringliteral}{"tke"}, g%HI, scale=us%Z\_to\_m**2*us%s\_to\_T**2) 333 \textcolor{keywordflow}{ endif} 334 335 \textcolor{keywordflow}{if} (cs%id\_Kd\_shear > 0) \textcolor{keyword}{call }post\_data(cs%id\_Kd\_shear, kappa\_io, cs%diag) 336 \textcolor{keywordflow}{if} (cs%id\_TKE > 0) \textcolor{keyword}{call }post\_data(cs%id\_TKE, tke\_io, cs%diag) 337 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__kappa__shear_a0b931b0b834d887e321eb6eb1924fa9a}\label{namespacemom__kappa__shear_a0b931b0b834d887e321eb6eb1924fa9a}} \index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}!calculate\+\_\+projected\+\_\+state@{calculate\+\_\+projected\+\_\+state}} \index{calculate\+\_\+projected\+\_\+state@{calculate\+\_\+projected\+\_\+state}!mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsubsection{\texorpdfstring{calculate\+\_\+projected\+\_\+state()}{calculate\_projected\_state()}} {\footnotesize\ttfamily subroutine mom\+\_\+kappa\+\_\+shear\+::calculate\+\_\+projected\+\_\+state (\begin{DoxyParamCaption}\item[{real, dimension(nz+1), intent(in)}]{kappa, }\item[{real, dimension(nz), intent(in)}]{u0, }\item[{real, dimension(nz), intent(in)}]{v0, }\item[{real, dimension(nz), intent(in)}]{T0, }\item[{real, dimension(nz), intent(in)}]{S0, }\item[{real, intent(in)}]{dt, }\item[{integer, intent(in)}]{nz, }\item[{real, dimension(nz), intent(in)}]{dz, }\item[{real, dimension(nz+1), intent(in)}]{I\+\_\+dz\+\_\+int, }\item[{real, dimension(nz+1), intent(in)}]{dbuoy\+\_\+dT, }\item[{real, dimension(nz+1), intent(in)}]{dbuoy\+\_\+dS, }\item[{real, dimension(nz), intent(inout)}]{u, }\item[{real, dimension(nz), intent(inout)}]{v, }\item[{real, dimension(nz), intent(inout)}]{T, }\item[{real, dimension(nz), intent(inout)}]{Sal, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US, }\item[{real, dimension(nz+1), intent(inout), optional}]{N2, }\item[{real, dimension(nz+1), intent(inout), optional}]{S2, }\item[{integer, intent(in), optional}]{ks\+\_\+int, }\item[{integer, intent(in), optional}]{ke\+\_\+int, }\item[{real, intent(in), optional}]{vel\+\_\+underflow }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} This subroutine calculates the velocities, temperature and salinity that the water column will have after mixing for dt with diffusivities kappa. It may also calculate the projected buoyancy frequency and shear. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em nz} & The number of layers (after eliminating massless layers?).\\ \hline \mbox{\tt in} & {\em kappa} & The diapycnal diffusivity at interfaces, \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em u0} & The initial zonal velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em v0} & The initial meridional velocity \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em t0} & The initial temperature \mbox{[}degC\mbox{]}.\\ \hline \mbox{\tt in} & {\em s0} & The initial salinity \mbox{[}ppt\mbox{]}.\\ \hline \mbox{\tt in} & {\em dz} & The grid spacing of layers \mbox{[}Z $\sim$$>$ m\mbox{]}.\\ \hline \mbox{\tt in} & {\em i\+\_\+dz\+\_\+int} & The inverse of the layer\textquotesingle{}s thicknesses \mbox{[}Z-\/1 $\sim$$>$ m-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em dbuoy\+\_\+dt} & The partial derivative of buoyancy with temperature \mbox{[}Z T-\/2 deg\+C-\/1 $\sim$$>$ m s-\/2 deg\+C-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em dbuoy\+\_\+ds} & The partial derivative of buoyancy with salinity \mbox{[}Z T-\/2 ppt-\/1 $\sim$$>$ m s-\/2 ppt-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em dt} & The time step \mbox{[}T $\sim$$>$ s\mbox{]}.\\ \hline \mbox{\tt in,out} & {\em u} & The zonal velocity after dt \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\ \hline \mbox{\tt in,out} & {\em v} & The meridional velocity after dt \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}.\\ \hline \mbox{\tt in,out} & {\em t} & The temperature after dt \mbox{[}degC\mbox{]}.\\ \hline \mbox{\tt in,out} & {\em sal} & The salinity after dt \mbox{[}ppt\mbox{]}.\\ \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,out} & {\em n2} & The buoyancy frequency squared at interfaces \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\ \hline \mbox{\tt in,out} & {\em s2} & The squared shear at interfaces \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em ks\+\_\+int} & The topmost k-\/index with a non-\/zero diffusivity.\\ \hline \mbox{\tt in} & {\em ke\+\_\+int} & The bottommost k-\/index with a non-\/zero diffusivity.\\ \hline \mbox{\tt in} & {\em vel\+\_\+underflow} & If present and true, any velocities that are smaller in magnitude than this value are set to 0 \mbox{[}L T-\/1 $\sim$$>$ m s-\/1\mbox{]}. \\ \hline \end{DoxyParams} Definition at line 1082 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 1082 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: nz\textcolor{comment}{ !< The number of layers (after eliminating massless} 1083 \textcolor{comment}{ !! layers?).} 1084 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: kappa\textcolor{comment}{ !< The diapycnal diffusivity at interfaces,} 1085 \textcolor{comment}{ !! [Z2 T-1 ~> m2 s-1].} 1086 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(in)} :: u0\textcolor{comment}{ !< The initial zonal velocity [L T-1 ~> m s-1].} 1087 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(in)} :: v0\textcolor{comment}{ !< The initial meridional velocity [L T-1 ~> m s-1].} 1088 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(in)} :: T0\textcolor{comment}{ !< The initial temperature [degC].} 1089 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(in)} :: S0\textcolor{comment}{ !< The initial salinity [ppt].} 1090 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(in)} :: dz\textcolor{comment}{ !< The grid spacing of layers [Z ~> m].} 1091 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: I\_dz\_int\textcolor{comment}{ !< The inverse of the layer's thicknesses} 1092 \textcolor{comment}{ !! [Z-1 ~> m-1].} 1093 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: dbuoy\_dT\textcolor{comment}{ !< The partial derivative of buoyancy with} 1094 \textcolor{comment}{ !! temperature [Z T-2 degC-1 ~> m s-2 degC-1].} 1095 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: dbuoy\_dS\textcolor{comment}{ !< The partial derivative of buoyancy with} 1096 \textcolor{comment}{ !! salinity [Z T-2 ppt-1 ~> m s-2 ppt-1].} 1097 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: dt\textcolor{comment}{ !< The time step [T ~> s].} 1098 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(inout)} :: u\textcolor{comment}{ !< The zonal velocity after dt [L T-1 ~> m s-1].} 1099 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(inout)} :: v\textcolor{comment}{ !< The meridional velocity after dt [L T-1 ~> m s-1].} 1100 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(inout)} :: T\textcolor{comment}{ !< The temperature after dt [degC].} 1101 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(inout)} :: Sal\textcolor{comment}{ !< The salinity after dt [ppt].} 1102 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< The ocean's vertical grid structure.} 1103 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 1104 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{optional}, & 1105 \textcolor{keywordtype}{intent(inout)} :: N2\textcolor{comment}{ !< The buoyancy frequency squared at interfaces [T-2 ~> s-2].} 1106 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{optional}, & 1107 \textcolor{keywordtype}{intent(inout)} :: S2\textcolor{comment}{ !< The squared shear at interfaces [T-2 ~> s-2].} 1108 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: ks\_int\textcolor{comment}{ !< The topmost k-index with a non-zero diffusivity.} 1109 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: ke\_int\textcolor{comment}{ !< The bottommost k-index with a non-zero} 1110 \textcolor{comment}{ !! diffusivity.} 1111 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: vel\_underflow\textcolor{comment}{ !< If present and true, any velocities that} 1112 \textcolor{comment}{ !! are smaller in magnitude than this value are} 1113 \textcolor{comment}{ !! set to 0 [L T-1 ~> m s-1].} 1114 1115 \textcolor{comment}{! Local variables} 1116 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)} :: c1 1117 \textcolor{keywordtype}{real} :: L2\_to\_Z2 \textcolor{comment}{! A conversion factor from horizontal length units to vertical depth} 1118 \textcolor{comment}{! units squared [Z2 s2 T-2 m-2 ~> 1].} 1119 \textcolor{keywordtype}{real} :: underflow\_vel \textcolor{comment}{! Velocities smaller in magnitude than underflow\_vel are set to 0 [L T-1 ~> m s-1].} 1120 \textcolor{keywordtype}{real} :: a\_a, a\_b, b1, d1, bd1, b1nz\_0 1121 \textcolor{keywordtype}{integer} :: k, ks, ke 1122 1123 ks = 1 ; ke = nz 1124 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(ks\_int)) ks = max(ks\_int-1,1) 1125 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(ke\_int)) ke = min(ke\_int,nz) 1126 underflow\_vel = 0.0 ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(vel\_underflow)) underflow\_vel = vel\_underflow 1127 1128 \textcolor{keywordflow}{if} (ks > ke) \textcolor{keywordflow}{return} 1129 1130 \textcolor{keywordflow}{if} (dt > 0.0) \textcolor{keywordflow}{then} 1131 a\_b = dt*(kappa(ks+1)*i\_dz\_int(ks+1)) 1132 b1 = 1.0 / (dz(ks) + a\_b) 1133 c1(ks+1) = a\_b * b1 ; d1 = dz(ks) * b1 \textcolor{comment}{! = 1 - c1} 1134 1135 u(ks) = (b1 * dz(ks))*u0(ks) ; v(ks) = (b1 * dz(ks))*v0(ks) 1136 t(ks) = (b1 * dz(ks))*t0(ks) ; sal(ks) = (b1 * dz(ks))*s0(ks) 1137 \textcolor{keywordflow}{do} k=ks+1,ke-1 1138 a\_a = a\_b 1139 a\_b = dt*(kappa(k+1)*i\_dz\_int(k+1)) 1140 bd1 = dz(k) + d1*a\_a 1141 b1 = 1.0 / (bd1 + a\_b) 1142 c1(k+1) = a\_b * b1 ; d1 = bd1 * b1 \textcolor{comment}{! d1 = 1 - c1} 1143 1144 u(k) = b1 * (dz(k)*u0(k) + a\_a*u(k-1)) 1145 v(k) = b1 * (dz(k)*v0(k) + a\_a*v(k-1)) 1146 t(k) = b1 * (dz(k)*t0(k) + a\_a*t(k-1)) 1147 sal(k) = b1 * (dz(k)*s0(k) + a\_a*sal(k-1)) 1148 \textcolor{keywordflow}{ enddo} 1149 \textcolor{comment}{! T and S have insulating boundary conditions, u & v use no-slip} 1150 \textcolor{comment}{! bottom boundary conditions at the solid bottom.} 1151 1152 \textcolor{comment}{! For insulating boundary conditions or mixing simply stopping, use...} 1153 a\_a = a\_b 1154 b1 = 1.0 / (dz(ke) + d1*a\_a) 1155 t(ke) = b1 * (dz(ke)*t0(ke) + a\_a*t(ke-1)) 1156 sal(ke) = b1 * (dz(ke)*s0(ke) + a\_a*sal(ke-1)) 1157 1158 \textcolor{comment}{! There is no distinction between the effective boundary conditions for} 1159 \textcolor{comment}{! tracers and velocities if the mixing is separated from the bottom, but if} 1160 \textcolor{comment}{! the mixing goes all the way to the bottom, use no-slip BCs for velocities.} 1161 \textcolor{keywordflow}{if} (ke == nz) \textcolor{keywordflow}{then} 1162 a\_b = dt*(kappa(nz+1)*i\_dz\_int(nz+1)) 1163 b1nz\_0 = 1.0 / ((dz(nz) + d1*a\_a) + a\_b) 1164 \textcolor{keywordflow}{else} 1165 b1nz\_0 = b1 1166 \textcolor{keywordflow}{ endif} 1167 u(ke) = b1nz\_0 * (dz(ke)*u0(ke) + a\_a*u(ke-1)) 1168 v(ke) = b1nz\_0 * (dz(ke)*v0(ke) + a\_a*v(ke-1)) 1169 \textcolor{keywordflow}{if} (abs(u(ke)) < underflow\_vel) u(ke) = 0.0 1170 \textcolor{keywordflow}{if} (abs(v(ke)) < underflow\_vel) v(ke) = 0.0 1171 1172 \textcolor{keywordflow}{do} k=ke-1,ks,-1 1173 u(k) = u(k) + c1(k+1)*u(k+1) 1174 v(k) = v(k) + c1(k+1)*v(k+1) 1175 \textcolor{keywordflow}{if} (abs(u(k)) < underflow\_vel) u(k) = 0.0 1176 \textcolor{keywordflow}{if} (abs(v(k)) < underflow\_vel) v(k) = 0.0 1177 t(k) = t(k) + c1(k+1)*t(k+1) 1178 sal(k) = sal(k) + c1(k+1)*sal(k+1) 1179 \textcolor{keywordflow}{ enddo} 1180 \textcolor{keywordflow}{else} \textcolor{comment}{! dt <= 0.0} 1181 \textcolor{keywordflow}{do} k=1,nz 1182 u(k) = u0(k) ; v(k) = v0(k) ; t(k) = t0(k) ; sal(k) = s0(k) 1183 \textcolor{keywordflow}{if} (abs(u(k)) < underflow\_vel) u(k) = 0.0 1184 \textcolor{keywordflow}{if} (abs(v(k)) < underflow\_vel) v(k) = 0.0 1185 \textcolor{keywordflow}{ enddo} 1186 \textcolor{keywordflow}{ endif} 1187 1188 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(s2)) \textcolor{keywordflow}{then} 1189 \textcolor{comment}{! L2\_to\_Z2 = US%m\_to\_Z**2 * US%T\_to\_s**2} 1190 l2\_to\_z2 = us%L\_to\_Z**2 1191 s2(1) = 0.0 ; s2(nz+1) = 0.0 1192 \textcolor{keywordflow}{if} (ks > 1) & 1193 s2(ks) = ((u(ks)-u0(ks-1))**2 + (v(ks)-v0(ks-1))**2) * (l2\_to\_z2*i\_dz\_int(ks)**2) 1194 \textcolor{keywordflow}{do} k=ks+1,ke 1195 s2(k) = ((u(k)-u(k-1))**2 + (v(k)-v(k-1))**2) * (l2\_to\_z2*i\_dz\_int(k)**2) 1196 \textcolor{keywordflow}{ enddo} 1197 \textcolor{keywordflow}{if} (ke 1) & 1204 n2(ks) = max(0.0, i\_dz\_int(ks) * & 1205 (dbuoy\_dt(ks) * (t0(ks-1)-t(ks)) + dbuoy\_ds(ks) * (s0(ks-1)-sal(ks)))) 1206 \textcolor{keywordflow}{do} k=ks+1,ke 1207 n2(k) = max(0.0, i\_dz\_int(k) * & 1208 (dbuoy\_dt(k) * (t(k-1)-t(k)) + dbuoy\_ds(k) * (sal(k-1)-sal(k)))) 1209 \textcolor{keywordflow}{ enddo} 1210 \textcolor{keywordflow}{if} (ke$ s-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em s2} & The squared shear at interfaces \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em kappa\+\_\+in} & The initial guess at the diffusivity \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em dz\+\_\+int} & The thicknesses associated with interfaces \mbox{[}Z-\/1 $\sim$$>$ m-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em i\+\_\+l2\+\_\+bdry} & The inverse of the squared distance to boundaries \mbox{[}Z-\/2 $\sim$$>$ m-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em idz} & The inverse grid spacing of layers \mbox{[}Z-\/1 $\sim$$>$ m-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em f2} & The squared Coriolis parameter \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\ \hline & {\em cs} & A pointer to this module\textquotesingle{}s control 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,out} & {\em k\+\_\+q} & The shear-\/driven diapycnal diffusivity divided by the turbulent kinetic energy per unit mass at interfaces \mbox{[}T $\sim$$>$ s\mbox{]}.\\ \hline \mbox{\tt out} & {\em tke} & The turbulent kinetic energy per unit mass at interfaces \mbox{[}Z2 T-\/2 $\sim$$>$ m2 s-\/2\mbox{]}.\\ \hline \mbox{\tt out} & {\em kappa} & The diapycnal diffusivity at interfaces \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \hline \mbox{\tt out} & {\em kappa\+\_\+src} & The source term for kappa \mbox{[}T-\/1 $\sim$$>$ s-\/1\mbox{]}.\\ \hline \mbox{\tt out} & {\em local\+\_\+src} & The sum of all local sources for kappa, \\ \hline \end{DoxyParams} Definition at line 1220 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 1220 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: nz\textcolor{comment}{ !< The number of layers to work on.} 1221 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: N2\textcolor{comment}{ !< The buoyancy frequency squared at interfaces [T-2 ~> s-2].} 1222 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: S2\textcolor{comment}{ !< The squared shear at interfaces [T-2 ~> s-2].} 1223 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: kappa\_in\textcolor{comment}{ !< The initial guess at the diffusivity} 1224 \textcolor{comment}{ !! [Z2 T-1 ~> m2 s-1].} 1225 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: dz\_Int\textcolor{comment}{ !< The thicknesses associated with interfaces} 1226 \textcolor{comment}{ !! [Z-1 ~> m-1].} 1227 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(in)} :: I\_L2\_bdry\textcolor{comment}{ !< The inverse of the squared distance to} 1228 \textcolor{comment}{ !! boundaries [Z-2 ~> m-2].} 1229 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)}, \textcolor{keywordtype}{intent(in)} :: Idz\textcolor{comment}{ !< The inverse grid spacing of layers [Z-1 ~> m-1].} 1230 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: f2\textcolor{comment}{ !< The squared Coriolis parameter [T-2 ~> s-2].} 1231 \textcolor{keywordtype}{type}(Kappa\_shear\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< A pointer to this module's control structure.} 1232 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< The ocean's vertical grid structure.} 1233 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 1234 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(inout)} :: K\_Q\textcolor{comment}{ !< The shear-driven diapycnal diffusivity divided by} 1235 \textcolor{comment}{ !! the turbulent kinetic energy per unit mass at} 1236 \textcolor{comment}{ !! interfaces [T ~> s].} 1237 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(out)} :: tke\textcolor{comment}{ !< The turbulent kinetic energy per unit mass at} 1238 \textcolor{comment}{ !! interfaces [Z2 T-2 ~> m2 s-2].} 1239 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{intent(out)} :: kappa\textcolor{comment}{ !< The diapycnal diffusivity at interfaces} 1240 \textcolor{comment}{ !! [Z2 T-1 ~> m2 s-1].} 1241 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{optional}, & 1242 \textcolor{keywordtype}{intent(out)} :: kappa\_src\textcolor{comment}{ !< The source term for kappa [T-1 ~> s-1].} 1243 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)}, \textcolor{keywordtype}{optional}, & 1244 \textcolor{keywordtype}{intent(out)} :: local\_src\textcolor{comment}{ !< The sum of all local sources for kappa,} 1245 \textcolor{comment}{ !! [T-1 ~> s-1].} 1246 \textcolor{comment}{! This subroutine calculates new, consistent estimates of TKE and kappa.} 1247 1248 \textcolor{comment}{! Local variables} 1249 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz)} :: & 1250 aQ, & \textcolor{comment}{! aQ is the coupling between adjacent interfaces in the TKE equations [Z T-1 ~> m s-1].} 1251 dQdz \textcolor{comment}{! Half the partial derivative of TKE with depth [Z T-2 ~> m s-2].} 1252 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)} :: & 1253 dK, & \textcolor{comment}{! The change in kappa [Z2 T-1 ~> m2 s-1].} 1254 dQ, & \textcolor{comment}{! The change in TKE [Z2 T-2 ~> m2 s-2].} 1255 cQ, cK, & \textcolor{comment}{! cQ and cK are the upward influences in the tridiagonal and} 1256 \textcolor{comment}{! hexadiagonal solvers for the TKE and kappa equations [nondim].} 1257 i\_ld2, & \textcolor{comment}{! 1/Ld^2, where Ld is the effective decay length scale for kappa [Z-2 ~> m-2].} 1258 tke\_decay, & \textcolor{comment}{! The local TKE decay rate [T-1 ~> s-1].} 1259 k\_src, & \textcolor{comment}{! The source term in the kappa equation [T-1 ~> s-1].} 1260 dqmdk, & \textcolor{comment}{! With Newton's method the change in dQ(k-1) due to dK(k) [T ~> s].} 1261 dkdq, & \textcolor{comment}{! With Newton's method the change in dK(k) due to dQ(k) [T-1 ~> s-1].} 1262 e1 \textcolor{comment}{! The fractional change in a layer TKE due to a change in the} 1263 \textcolor{comment}{! TKE of the layer above when all the kappas below are 0.} 1264 \textcolor{comment}{! e1 is nondimensional, and 0 < e1 < 1.} 1265 \textcolor{keywordtype}{real} :: tke\_src \textcolor{comment}{! The net source of TKE due to mixing against the shear} 1266 \textcolor{comment}{! and stratification [Z2 T-3 ~> m2 s-3]. (For convenience,} 1267 \textcolor{comment}{! a term involving the non-dissipation of q0 is also} 1268 \textcolor{comment}{! included here.)} 1269 \textcolor{keywordtype}{real} :: bQ \textcolor{comment}{! The inverse of the pivot in the tridiagonal equations [T Z-1 ~> s m-1].} 1270 \textcolor{keywordtype}{real} :: bK \textcolor{comment}{! The inverse of the pivot in the tridiagonal equations [Z-1 ~> m-1].} 1271 \textcolor{keywordtype}{real} :: bQd1 \textcolor{comment}{! A term in the denominator of bQ [Z T-1 ~> m s-1].} 1272 \textcolor{keywordtype}{real} :: bKd1 \textcolor{comment}{! A term in the denominator of bK [Z ~> m].} 1273 \textcolor{keywordtype}{real} :: cQcomp, cKcomp \textcolor{comment}{! 1 - cQ or 1 - cK in the tridiagonal equations.} 1274 \textcolor{keywordtype}{real} :: c\_s2 \textcolor{comment}{! The coefficient for the decay of TKE due to} 1275 \textcolor{comment}{! shear (i.e. proportional to |S|*tke), nondimensional.} 1276 \textcolor{keywordtype}{real} :: c\_n2 \textcolor{comment}{! The coefficient for the decay of TKE due to} 1277 \textcolor{comment}{! stratification (i.e. proportional to N*tke) [nondim].} 1278 \textcolor{keywordtype}{real} :: Ri\_crit \textcolor{comment}{! The critical shear Richardson number for shear-} 1279 \textcolor{comment}{! driven mixing. The theoretical value is 0.25.} 1280 \textcolor{keywordtype}{real} :: q0 \textcolor{comment}{! The background level of TKE [Z2 T-2 ~> m2 s-2].} 1281 \textcolor{keywordtype}{real} :: Ilambda2 \textcolor{comment}{! 1.0 / CS%lambda**2 [nondim]} 1282 \textcolor{keywordtype}{real} :: TKE\_min \textcolor{comment}{! The minimum value of shear-driven TKE that can be} 1283 \textcolor{comment}{! solved for [Z2 T-2 ~> m2 s-2].} 1284 \textcolor{keywordtype}{real} :: kappa0 \textcolor{comment}{! The background diapycnal diffusivity [Z2 T-1 ~> m2 s-1].} 1285 \textcolor{keywordtype}{real} :: kappa\_trunc \textcolor{comment}{! Diffusivities smaller than this are rounded to 0 [Z2 T-1 ~> m2 s-1].} 1286 1287 \textcolor{keywordtype}{real} :: eden1, eden2, I\_eden, ome \textcolor{comment}{! Variables used in calculating e1.} 1288 \textcolor{keywordtype}{real} :: diffusive\_src \textcolor{comment}{! The diffusive source in the kappa equation [Z T-1 ~> m s-1].} 1289 \textcolor{keywordtype}{real} :: chg\_by\_k0 \textcolor{comment}{! The value of k\_src that leads to an increase of} 1290 \textcolor{comment}{! kappa\_0 if only the diffusive term is a sink [T-1 ~> s-1].} 1291 1292 \textcolor{keywordtype}{real} :: kappa\_mean \textcolor{comment}{! A mean value of kappa [Z2 T-1 ~> m2 s-1].} 1293 \textcolor{keywordtype}{real} :: Newton\_test \textcolor{comment}{! The value of relative error that will cause the next} 1294 \textcolor{comment}{! iteration to use Newton's method.} 1295 \textcolor{comment}{! Temporary variables used in the Newton's method iterations.} 1296 \textcolor{keywordtype}{real} :: decay\_term\_k \textcolor{comment}{! The decay term in the diffusivity equation} 1297 \textcolor{keywordtype}{real} :: decay\_term\_Q \textcolor{comment}{! The decay term in the TKE equation - proportional to [T-1 ~> s-1]} 1298 \textcolor{keywordtype}{real} :: I\_Q \textcolor{comment}{! The inverse of TKE [T2 Z-2 ~> s2 m-2]} 1299 \textcolor{keywordtype}{real} :: kap\_src 1300 \textcolor{keywordtype}{real} :: v1 \textcolor{comment}{! A temporary variable proportional to [T-1 ~> s-1]} 1301 \textcolor{keywordtype}{real} :: v2 1302 \textcolor{keywordtype}{real} :: tol\_err \textcolor{comment}{! The tolerance for max\_err that determines when to} 1303 \textcolor{comment}{! stop iterating.} 1304 \textcolor{keywordtype}{real} :: Newton\_err \textcolor{comment}{! The tolerance for max\_err that determines when to} 1305 \textcolor{comment}{! start using Newton's method. Empirically, an initial} 1306 \textcolor{comment}{! value of about 0.2 seems to be most efficient.} 1307 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: roundoff = 1.0e-16 \textcolor{comment}{! A negligible fractional change in TKE.} 1308 \textcolor{comment}{! This could be larger but performance gains are small.} 1309 1310 \textcolor{keywordtype}{logical} :: tke\_noflux\_bottom\_BC = .false. \textcolor{comment}{! Specify the boundary conditions} 1311 \textcolor{keywordtype}{logical} :: tke\_noflux\_top\_BC = .false. \textcolor{comment}{! that are applied to the TKE eqns.} 1312 \textcolor{keywordtype}{logical} :: do\_Newton \textcolor{comment}{! If .true., use Newton's method for the next iteration.} 1313 \textcolor{keywordtype}{logical} :: abort\_Newton \textcolor{comment}{! If .true., an Newton's method has encountered a 0} 1314 \textcolor{comment}{! pivot, and should not have been used.} 1315 \textcolor{keywordtype}{logical} :: was\_Newton \textcolor{comment}{! The value of do\_Newton before checking convergence.} 1316 \textcolor{keywordtype}{logical} :: within\_tolerance \textcolor{comment}{! If .true., all points are within tolerance to} 1317 \textcolor{comment}{! enable this subroutine to return.} 1318 \textcolor{keywordtype}{integer} :: ks\_src, ke\_src \textcolor{comment}{! The range indices that have nonzero k\_src.} 1319 \textcolor{keywordtype}{integer} :: ks\_kappa, ke\_kappa, ke\_tke \textcolor{comment}{! The ranges of k-indices that are or} 1320 \textcolor{keywordtype}{integer} :: ks\_kappa\_prev, ke\_kappa\_prev \textcolor{comment}{! were being worked on.} 1321 \textcolor{keywordtype}{integer} :: itt, k, k2 1322 1323 \textcolor{comment}{! These variables are used only for debugging.} 1324 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{parameter} :: debug\_soln = .false. 1325 \textcolor{keywordtype}{real} :: K\_err\_lin, Q\_err\_lin, TKE\_src\_norm 1326 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nz+1)} :: & 1327 I\_Ld2\_debug, & \textcolor{comment}{! A separate version of I\_Ld2 for debugging [Z-2 ~> m-2].} 1328 kappa\_prev, & \textcolor{comment}{! The value of kappa at the start of the current iteration [Z2 T-1 ~> m2 s-1].} 1329 TKE\_prev \textcolor{comment}{! The value of TKE at the start of the current iteration [Z2 T-2 ~> m2 s-2].} 1330 1331 c\_n2 = cs%C\_N**2 ; c\_s2 = cs%C\_S**2 1332 q0 = cs%TKE\_bg ; kappa0 = cs%kappa\_0 1333 tke\_min = max(cs%TKE\_bg, 1.0e-20*us%m\_to\_Z**2*us%T\_to\_s**2) 1334 ri\_crit = cs%Rino\_crit 1335 ilambda2 = 1.0 / cs%lambda**2 1336 kappa\_trunc = cs%kappa\_trunc 1337 do\_newton = .false. ; abort\_newton = .false. 1338 tol\_err = cs%kappa\_tol\_err 1339 newton\_err = 0.2 \textcolor{comment}{! This initial value may be automatically reduced later.} 1340 1341 ks\_kappa = 2 ; ke\_kappa = nz ; ks\_kappa\_prev = 2 ; ke\_kappa\_prev = nz 1342 1343 ke\_src = 0 ; ks\_src = nz+1 1344 \textcolor{keywordflow}{do} k=2,nz 1345 \textcolor{keywordflow}{if} (n2(k) < ri\_crit * s2(k)) \textcolor{keywordflow}{then} \textcolor{comment}{! Equivalent to Ri < Ri\_crit.} 1346 \textcolor{comment}{! Ri = N2(K) / S2(K)} 1347 \textcolor{comment}{! k\_src(K) = (2.0 * CS%Shearmix\_rate * sqrt(S2(K))) * &} 1348 \textcolor{comment}{! ((Ri\_crit - Ri) / (Ri\_crit + CS%FRi\_curvature*Ri))} 1349 k\_src(k) = (2.0 * cs%Shearmix\_rate * sqrt(s2(k))) * & 1350 ((ri\_crit*s2(k) - n2(k)) / (ri\_crit*s2(k) + cs%FRi\_curvature*n2(k))) 1351 ke\_src = k 1352 \textcolor{keywordflow}{if} (ks\_src > k) ks\_src = k 1353 \textcolor{keywordflow}{else} 1354 k\_src(k) = 0.0 1355 \textcolor{keywordflow}{ endif} 1356 \textcolor{keywordflow}{ enddo} 1357 1358 \textcolor{comment}{! If there is no source anywhere, return kappa(K) = 0.} 1359 \textcolor{keywordflow}{if} (ks\_src > ke\_src) \textcolor{keywordflow}{then} 1360 \textcolor{keywordflow}{do} k=1,nz+1 1361 kappa(k) = 0.0 ; k\_q(k) = 0.0 ; tke(k) = tke\_min 1362 \textcolor{keywordflow}{ enddo} 1363 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kappa\_src)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; kappa\_src(k) = 0.0 ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1364 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(local\_src)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; local\_src(k) = 0.0 ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1365 \textcolor{keywordflow}{return} 1366 \textcolor{keywordflow}{ endif} 1367 1368 \textcolor{keywordflow}{do} k=1,nz+1 1369 kappa(k) = kappa\_in(k) 1370 \textcolor{comment}{! TKE\_decay(K) = c\_n*sqrt(N2(K)) + c\_s*sqrt(S2(K)) ! The expression in JHL.} 1371 tke\_decay(k) = sqrt(c\_n2*n2(k) + c\_s2*s2(k)) 1372 \textcolor{keywordflow}{if} ((kappa(k) > 0.0) .and. (k\_q(k) > 0.0)) \textcolor{keywordflow}{then} 1373 tke(k) = kappa(k) / k\_q(k) \textcolor{comment}{! Perhaps take the max with TKE\_min} 1374 \textcolor{keywordflow}{else} 1375 tke(k) = tke\_min 1376 \textcolor{keywordflow}{ endif} 1377 \textcolor{keywordflow}{ enddo} 1378 \textcolor{comment}{! Apply boundary conditions to kappa.} 1379 kappa(1) = 0.0 ; kappa(nz+1) = 0.0 1380 1381 \textcolor{comment}{! Calculate the term (e1) that allows changes in TKE to be calculated quickly} 1382 \textcolor{comment}{! below the deepest nonzero value of kappa. If kappa = 0, below interface} 1383 \textcolor{comment}{! k-1, the final changes in TKE are related by dQ(K+1) = e1(K+1)*dQ(K).} 1384 eden2 = kappa0 * idz(nz) 1385 \textcolor{keywordflow}{if} (tke\_noflux\_bottom\_bc) \textcolor{keywordflow}{then} 1386 eden1 = dz\_int(nz+1)*tke\_decay(nz+1) 1387 i\_eden = 1.0 / (eden2 + eden1) 1388 e1(nz+1) = eden2 * i\_eden ; ome = eden1 * i\_eden 1389 \textcolor{keywordflow}{else} 1390 e1(nz+1) = 0.0 ; ome = 1.0 1391 \textcolor{keywordflow}{ endif} 1392 \textcolor{keywordflow}{do} k=nz,2,-1 1393 eden1 = dz\_int(k)*tke\_decay(k) + ome * eden2 1394 eden2 = kappa0 * idz(k-1) 1395 i\_eden = 1.0 / (eden2 + eden1) 1396 e1(k) = eden2 * i\_eden ; ome = eden1 * i\_eden \textcolor{comment}{! = 1-e1} 1397 \textcolor{keywordflow}{ enddo} 1398 e1(1) = 0.0 1399 1400 1401 \textcolor{comment}{! Iterate here to convergence to within some tolerance of order tol\_err.} 1402 \textcolor{keywordflow}{do} itt=1,cs%max\_RiNo\_it 1403 1404 \textcolor{comment}{! ----------------------------------------------------} 1405 \textcolor{comment}{! Calculate TKE} 1406 \textcolor{comment}{! ----------------------------------------------------} 1407 1408 \textcolor{keywordflow}{if} (debug\_soln) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=1,nz+1 ; kappa\_prev(k) = kappa(k) ; tke\_prev(k) = tke(k) ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1409 1410 \textcolor{keywordflow}{if} (.not.do\_newton) \textcolor{keywordflow}{then} 1411 \textcolor{comment}{! Use separate steps of the TKE and kappa equations, that are} 1412 \textcolor{comment}{! explicit in the nonlinear source terms, implicit in a linearized} 1413 \textcolor{comment}{! version of the nonlinear sink terms, and implicit in the linear} 1414 \textcolor{comment}{! terms.} 1415 1416 ke\_tke = max(ke\_kappa,ke\_kappa\_prev)+1 1417 \textcolor{comment}{! aQ is the coupling between adjacent interfaces [Z T-1 ~> m s-1].} 1418 \textcolor{keywordflow}{do} k=1,min(ke\_tke,nz) 1419 aq(k) = (0.5*(kappa(k)+kappa(k+1)) + kappa0) * idz(k) 1420 \textcolor{keywordflow}{ enddo} 1421 dq(1) = -tke(1) 1422 \textcolor{keywordflow}{if} (tke\_noflux\_top\_bc) \textcolor{keywordflow}{then} 1423 tke\_src = kappa0*s2(1) + q0 * tke\_decay(1) \textcolor{comment}{! Uses that kappa(1) = 0} 1424 bqd1 = dz\_int(1) * tke\_decay(1) 1425 bq = 1.0 / (bqd1 + aq(1)) 1426 tke(1) = bq * (dz\_int(1)*tke\_src) 1427 cq(2) = aq(1) * bq ; cqcomp = bqd1 * bq \textcolor{comment}{! = 1 - cQ} 1428 \textcolor{keywordflow}{else} 1429 tke(1) = q0 ; cq(2) = 0.0 ; cqcomp = 1.0 1430 \textcolor{keywordflow}{ endif} 1431 \textcolor{keywordflow}{do} k=2,ke\_tke-1 1432 dq(k) = -tke(k) 1433 tke\_src = (kappa(k) + kappa0)*s2(k) + q0*tke\_decay(k) 1434 bqd1 = dz\_int(k)*(tke\_decay(k) + n2(k)*k\_q(k)) + cqcomp*aq(k-1) 1435 bq = 1.0 / (bqd1 + aq(k)) 1436 tke(k) = bq * (dz\_int(k)*tke\_src + aq(k-1)*tke(k-1)) 1437 cq(k+1) = aq(k) * bq ; cqcomp = bqd1 * bq \textcolor{comment}{! = 1 - cQ} 1438 \textcolor{keywordflow}{ enddo} 1439 \textcolor{keywordflow}{if} ((ke\_tke == nz+1) .and. .not.(tke\_noflux\_bottom\_bc)) \textcolor{keywordflow}{then} 1440 tke(nz+1) = tke\_min 1441 dq(nz+1) = 0.0 1442 \textcolor{keywordflow}{else} 1443 k = ke\_tke 1444 tke\_src = kappa0*s2(k) + q0*tke\_decay(k) \textcolor{comment}{! Uses that kappa(ke\_tke) = 0} 1445 \textcolor{keywordflow}{if} (k == nz+1) \textcolor{keywordflow}{then} 1446 dq(k) = -tke(k) 1447 bq = 1.0 / (dz\_int(k)*tke\_decay(k) + cqcomp*aq(k-1)) 1448 tke(k) = max(tke\_min, bq * (dz\_int(k)*tke\_src + aq(k-1)*tke(k-1))) 1449 dq(k) = tke(k) + dq(k) 1450 \textcolor{keywordflow}{else} 1451 bq = 1.0 / ((dz\_int(k)*tke\_decay(k) + cqcomp*aq(k-1)) + aq(k)) 1452 cq(k+1) = aq(k) * bq 1453 \textcolor{comment}{! Account for all changes deeper in the water column.} 1454 dq(k) = -tke(k) 1455 tke(k) = max((bq * (dz\_int(k)*tke\_src + aq(k-1)*tke(k-1)) + & 1456 cq(k+1)*(tke(k+1) - e1(k+1)*tke(k))) / (1.0 - cq(k+1)*e1(k+1)), tke\_min) 1457 dq(k) = tke(k) + dq(k) 1458 1459 \textcolor{comment}{! Adjust TKE deeper in the water column in case ke\_tke increases.} 1460 \textcolor{comment}{! This might not be strictly necessary?} 1461 \textcolor{keywordflow}{do} k=ke\_tke+1,nz+1 1462 dq(k) = e1(k)*dq(k-1) 1463 tke(k) = max(tke(k) + dq(k), tke\_min) 1464 \textcolor{keywordflow}{if} (abs(dq(k)) < roundoff*tke(k)) \textcolor{keywordflow}{exit} 1465 \textcolor{keywordflow}{ enddo} 1466 \textcolor{keywordflow}{do} k2=k+1,nz 1467 \textcolor{keywordflow}{if} (dq(k2) == 0.0) \textcolor{keywordflow}{exit} 1468 dq(k2) = 0.0 1469 \textcolor{keywordflow}{ enddo} 1470 \textcolor{keywordflow}{ endif} 1471 \textcolor{keywordflow}{ endif} 1472 \textcolor{keywordflow}{do} k=ke\_tke-1,1,-1 1473 tke(k) = max(tke(k) + cq(k+1)*tke(k+1), tke\_min) 1474 dq(k) = tke(k) + dq(k) 1475 \textcolor{keywordflow}{ enddo} 1476 1477 \textcolor{comment}{! ----------------------------------------------------} 1478 \textcolor{comment}{! Calculate kappa, here defined at interfaces.} 1479 \textcolor{comment}{! ----------------------------------------------------} 1480 1481 ke\_kappa\_prev = ke\_kappa ; ks\_kappa\_prev = ks\_kappa 1482 1483 dk(1) = 0.0 \textcolor{comment}{! kappa takes boundary values of 0.} 1484 ck(2) = 0.0 ; ckcomp = 1.0 1485 \textcolor{keywordflow}{if} (itt == 1) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{dO} k=2,nz 1486 i\_ld2(k) = (n2(k)*ilambda2 + f2) / tke(k) + i\_l2\_bdry(k) 1487 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1488 \textcolor{keywordflow}{do} k=2,nz 1489 dk(k) = -kappa(k) 1490 \textcolor{keywordflow}{if} (itt>1) & 1491 i\_ld2(k) = (n2(k)*ilambda2 + f2) / tke(k) + i\_l2\_bdry(k) 1492 bkd1 = dz\_int(k)*i\_ld2(k) + ckcomp*idz(k-1) 1493 bk = 1.0 / (bkd1 + idz(k)) 1494 1495 kappa(k) = bk * (idz(k-1)*kappa(k-1) + dz\_int(k) * k\_src(k)) 1496 ck(k+1) = idz(k) * bk ; ckcomp = bkd1 * bk \textcolor{comment}{! = 1 - cK(K+1)} 1497 1498 \textcolor{comment}{! Neglect values that are smaller than kappa\_trunc.} 1499 \textcolor{keywordflow}{if} (kappa(k) < ckcomp*kappa\_trunc) \textcolor{keywordflow}{then} 1500 kappa(k) = 0.0 1501 \textcolor{keywordflow}{if} (k > ke\_src) \textcolor{keywordflow}{then} ; ke\_kappa = k-1 ; k\_q(k) = 0.0 ; \textcolor{keywordflow}{exit} ;\textcolor{keywordflow}{ endif} 1502 \textcolor{keywordflow}{elseif} (kappa(k) < 2.0*ckcomp*kappa\_trunc) \textcolor{keywordflow}{then} 1503 kappa(k) = 2.0 * (kappa(k) - ckcomp*kappa\_trunc) 1504 \textcolor{keywordflow}{ endif} 1505 \textcolor{keywordflow}{ enddo} 1506 k\_q(ke\_kappa) = kappa(ke\_kappa) / tke(ke\_kappa) 1507 dk(ke\_kappa) = dk(ke\_kappa) + kappa(ke\_kappa) 1508 \textcolor{keywordflow}{do} k=ke\_kappa+2,ke\_kappa\_prev 1509 dk(k) = -kappa(k) ; kappa(k) = 0.0 ; k\_q(k) = 0.0 1510 \textcolor{keywordflow}{ enddo} 1511 \textcolor{keywordflow}{do} k=ke\_kappa-1,2,-1 1512 kappa(k) = kappa(k) + ck(k+1)*kappa(k+1) 1513 \textcolor{comment}{! Neglect values that are smaller than kappa\_trunc.} 1514 \textcolor{keywordflow}{if} (kappa(k) <= kappa\_trunc) \textcolor{keywordflow}{then} 1515 kappa(k) = 0.0 1516 \textcolor{keywordflow}{if} (k < ks\_src) \textcolor{keywordflow}{then} ; ks\_kappa = k+1 ; k\_q(k) = 0.0 ; \textcolor{keywordflow}{exit} ;\textcolor{keywordflow}{ endif} 1517 \textcolor{keywordflow}{elseif} (kappa(k) < 2.0*kappa\_trunc) \textcolor{keywordflow}{then} 1518 kappa(k) = 2.0 * (kappa(k) - kappa\_trunc) 1519 \textcolor{keywordflow}{ endif} 1520 1521 dk(k) = dk(k) + kappa(k) 1522 k\_q(k) = kappa(k) / tke(k) 1523 \textcolor{keywordflow}{ enddo} 1524 \textcolor{keywordflow}{do} k=ks\_kappa\_prev,ks\_kappa-2 ; kappa(k) = 0.0 ; k\_q(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 1525 1526 \textcolor{keywordflow}{else} \textcolor{comment}{! do\_Newton is .true.} 1527 \textcolor{comment}{! Once the solutions are close enough, use a Newton's method solver of the} 1528 \textcolor{comment}{! whole system to accelerate convergence.} 1529 ks\_kappa\_prev = ks\_kappa ; ke\_kappa\_prev = ke\_kappa ; ke\_kappa = nz 1530 ks\_kappa = 2 1531 dk(1) = 0.0 ; ck(2) = 0.0 ; ckcomp = 1.0 ; dkdq(1) = 0.0 1532 aq(1) = (0.5*(kappa(1)+kappa(2))+kappa0) * idz(1) 1533 dqdz(1) = 0.5*(tke(1) - tke(2))*idz(1) 1534 \textcolor{keywordflow}{if} (tke\_noflux\_top\_bc) \textcolor{keywordflow}{then} 1535 tke\_src = dz\_int(1) * (kappa0*s2(1) - (tke(1) - q0)*tke\_decay(1)) - & 1536 aq(1) * (tke(1) - tke(2)) 1537 1538 bq = 1.0 / (aq(1) + dz\_int(1)*tke\_decay(1)) 1539 cq(2) = aq(1) * bq 1540 cqcomp = (dz\_int(1)*tke\_decay(1)) * bq \textcolor{comment}{! = 1 - cQ(2)} 1541 dqmdk(2) = -dqdz(1) * bq 1542 dq(1) = bq * tke\_src 1543 \textcolor{keywordflow}{else} 1544 dq(1) = 0.0 ; cq(2) = 0.0 ; cqcomp = 1.0 ; dqmdk(2) = 0.0 1545 \textcolor{keywordflow}{ endif} 1546 \textcolor{keywordflow}{do} k=2,nz 1547 i\_q = 1.0 / tke(k) 1548 i\_ld2(k) = (n2(k)*ilambda2 + f2) * i\_q + i\_l2\_bdry(k) 1549 1550 kap\_src = dz\_int(k) * (k\_src(k) - i\_ld2(k)*kappa(k)) + & 1551 idz(k-1)*(kappa(k-1)-kappa(k)) - idz(k)*(kappa(k)-kappa(k+1)) 1552 1553 \textcolor{comment}{! Ensure that the pivot is always positive, and that 0 <= cK <= 1.} 1554 \textcolor{comment}{! Otherwise do not use Newton's method.} 1555 decay\_term\_k = -idz(k-1)*dqmdk(k)*dkdq(k-1) + dz\_int(k)*i\_ld2(k) 1556 \textcolor{keywordflow}{if} (decay\_term\_k < 0.0) \textcolor{keywordflow}{then} ; abort\_newton = .true. ; \textcolor{keywordflow}{exit} ;\textcolor{keywordflow}{ endif} 1557 bk = 1.0 / (idz(k) + idz(k-1)*ckcomp + decay\_term\_k) 1558 1559 ck(k+1) = bk * idz(k) 1560 ckcomp = bk * (idz(k-1)*ckcomp + decay\_term\_k) \textcolor{comment}{! = 1-cK(K+1)} 1561 \textcolor{keywordflow}{if} (cs%dKdQ\_iteration\_bug) \textcolor{keywordflow}{then} 1562 dkdq(k) = bk * (idz(k-1)*dkdq(k-1)*cq(k) + & 1563 us%m\_to\_Z*(n2(k)*ilambda2 + f2) * i\_q**2 * kappa(k) ) 1564 \textcolor{keywordflow}{else} 1565 dkdq(k) = bk * (idz(k-1)*dkdq(k-1)*cq(k) + & 1566 dz\_int(k)*(n2(k)*ilambda2 + f2) * i\_q**2 * kappa(k) ) 1567 \textcolor{keywordflow}{ endif} 1568 dk(k) = bk * (kap\_src + idz(k-1)*dk(k-1) + idz(k-1)*dkdq(k-1)*dq(k-1)) 1569 1570 \textcolor{comment}{! Truncate away negligibly small values of kappa.} 1571 \textcolor{keywordflow}{if} (dk(k) <= ckcomp*(kappa\_trunc - kappa(k))) \textcolor{keywordflow}{then} 1572 dk(k) = -ckcomp*kappa(k) 1573 \textcolor{comment}{! if (K > ke\_src) then ; ke\_kappa = k-1 ; K\_Q(K) = 0.0 ; exit ; endif} 1574 \textcolor{keywordflow}{elseif} (dk(k) < ckcomp*(2.0*kappa\_trunc - kappa(k))) \textcolor{keywordflow}{then} 1575 dk(k) = 2.0 * dk(k) - ckcomp*(2.0*kappa\_trunc - kappa(k)) 1576 \textcolor{keywordflow}{ endif} 1577 1578 \textcolor{comment}{! Solve for dQ(K)...} 1579 aq(k) = (0.5*(kappa(k)+kappa(k+1))+kappa0) * idz(k) 1580 dqdz(k) = 0.5*(tke(k) - tke(k+1))*idz(k) 1581 tke\_src = dz\_int(k) * (((kappa(k) + kappa0)*s2(k) - kappa(k)*n2(k)) - & 1582 (tke(k) - q0)*tke\_decay(k)) - & 1583 (aq(k) * (tke(k) - tke(k+1)) - aq(k-1) * (tke(k-1) - tke(k))) 1584 v1 = aq(k-1) + dqdz(k-1)*dkdq(k-1) 1585 v2 = (v1*dqmdk(k) + dqdz(k-1)*ck(k)) + & 1586 ((dqdz(k-1) - dqdz(k)) + dz\_int(k)*(s2(k) - n2(k))) 1587 1588 \textcolor{comment}{! Ensure that the pivot is always positive, and that 0 <= cQ <= 1.} 1589 \textcolor{comment}{! Otherwise do not use Newton's method.} 1590 decay\_term\_q = dz\_int(k)*tke\_decay(k) - dqdz(k-1)*dkdq(k-1)*cq(k) - v2*dkdq(k) 1591 \textcolor{keywordflow}{if} (decay\_term\_q < 0.0) \textcolor{keywordflow}{then} ; abort\_newton = .true. ; \textcolor{keywordflow}{exit} ;\textcolor{keywordflow}{ endif} 1592 bq = 1.0 / (aq(k) + (cqcomp*aq(k-1) + decay\_term\_q)) 1593 1594 cq(k+1) = aq(k) * bq 1595 cqcomp = (cqcomp*aq(k-1) + decay\_term\_q) * bq 1596 dqmdk(k+1) = (v2 * ck(k+1) - dqdz(k)) * bq 1597 1598 \textcolor{comment}{! Ensure that TKE+dQ will not drop below 0.5*TKE.} 1599 dq(k) = max(bq * ((v1 * dq(k-1) + dqdz(k-1)*dk(k-1)) + & 1600 (v2 * dk(k) + tke\_src)), cqcomp*(-0.5*tke(k))) 1601 1602 \textcolor{comment}{! Check whether the next layer will be affected by any nonzero kappas.} 1603 \textcolor{keywordflow}{if} ((itt > 1) .and. (k > ke\_src) .and. (dk(k) == 0.0) .and. & 1604 ((kappa(k) + kappa(k+1)) == 0.0)) \textcolor{keywordflow}{then} 1605 \textcolor{comment}{! Could also do .and. (bQ*abs(tke\_src) < roundoff*TKE(K)) then} 1606 ke\_kappa = k-1 ; \textcolor{keywordflow}{exit} 1607 \textcolor{keywordflow}{ endif} 1608 \textcolor{keywordflow}{ enddo} 1609 \textcolor{keywordflow}{if} ((ke\_kappa == nz) .and. (.not. abort\_newton)) \textcolor{keywordflow}{then} 1610 dk(nz+1) = 0.0 ; dkdq(nz+1) = 0.0 1611 \textcolor{keywordflow}{if} (tke\_noflux\_bottom\_bc) \textcolor{keywordflow}{then} 1612 k = nz+1 1613 tke\_src = dz\_int(k) * (kappa0*s2(k) - (tke(k) - q0)*tke\_decay(k)) + & 1614 aq(k-1) * (tke(k-1) - tke(k)) 1615 1616 v1 = aq(k-1) + dqdz(k-1)*dkdq(k-1) 1617 decay\_term\_q = max(0.0, dz\_int(k)*tke\_decay(k) - dqdz(k-1)*dkdq(k-1)*cq(k)) 1618 \textcolor{keywordflow}{if} (decay\_term\_q < 0.0) \textcolor{keywordflow}{then} 1619 abort\_newton = .true. 1620 \textcolor{keywordflow}{else} 1621 bq = 1.0 / (aq(k) + (cqcomp*aq(k-1) + decay\_term\_q)) 1622 \textcolor{comment}{! Ensure that TKE+dQ will not drop below 0.5*TKE.} 1623 dq(k) = max(bq * ((v1 * dq(k-1) + dqdz(k-1)*dk(k-1)) + tke\_src), -0.5*tke(k)) 1624 tke(k) = max(tke(k) + dq(k), tke\_min) 1625 \textcolor{keywordflow}{ endif} 1626 \textcolor{keywordflow}{else} 1627 dq(nz+1) = 0.0 1628 \textcolor{keywordflow}{ endif} 1629 \textcolor{keywordflow}{elseif} (.not. abort\_newton) \textcolor{keywordflow}{then} 1630 \textcolor{comment}{! Alter the first-guess determination of dQ(K).} 1631 dq(ke\_kappa+1) = dq(ke\_kappa+1) / (1.0 - cq(ke\_kappa+2)*e1(ke\_kappa+2)) 1632 tke(ke\_kappa+1) = max(tke(ke\_kappa+1) + dq(ke\_kappa+1), tke\_min) 1633 \textcolor{keywordflow}{do} k=ke\_kappa+2,nz+1 1634 \textcolor{keywordflow}{if} (debug\_soln .and. (k < nz+1)) \textcolor{keywordflow}{then} 1635 \textcolor{comment}{! Ignore this source?} 1636 aq(k) = (0.5*(kappa(k)+kappa(k+1))+kappa0) * idz(k) 1637 \textcolor{comment}{! tke\_src\_norm = (dz\_Int(K) * (kappa0*S2(K) - (TKE(K)-q0)*TKE\_decay(K)) - &} 1638 \textcolor{comment}{! (aQ(k) * (TKE(K) - TKE(K+1)) - aQ(k-1) * (TKE(K-1) - TKE(K))) ) / &} 1639 \textcolor{comment}{! (aQ(k) + (aQ(k-1) + dz\_Int(K)*TKE\_decay(K)))} 1640 \textcolor{keywordflow}{ endif} 1641 dk(k) = 0.0 1642 \textcolor{comment}{! Ensure that TKE+dQ will not drop below 0.5*TKE.} 1643 dq(k) = max(e1(k)*dq(k-1),-0.5*tke(k)) 1644 tke(k) = max(tke(k) + dq(k), tke\_min) 1645 \textcolor{keywordflow}{if} (abs(dq(k)) < roundoff*tke(k)) \textcolor{keywordflow}{exit} 1646 \textcolor{keywordflow}{ enddo} 1647 \textcolor{keywordflow}{if} (debug\_soln) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k2=k+1,nz+1 ; dq(k2) = 0.0 ; dk(k2) = 0.0 ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1648 \textcolor{keywordflow}{ endif} 1649 \textcolor{keywordflow}{if} (.not. abort\_newton) \textcolor{keywordflow}{then} 1650 \textcolor{keywordflow}{do} k=ke\_kappa,2,-1 1651 \textcolor{comment}{! Ensure that TKE+dQ will not drop below 0.5*TKE.} 1652 dq(k) = max(dq(k) + (cq(k+1)*dq(k+1) + dqmdk(k+1) * dk(k+1)), -0.5*tke(k)) 1653 tke(k) = max(tke(k) + dq(k), tke\_min) 1654 dk(k) = dk(k) + (ck(k+1)*dk(k+1) + dkdq(k) * dq(k)) 1655 \textcolor{comment}{! Truncate away negligibly small values of kappa.} 1656 \textcolor{keywordflow}{if} (dk(k) <= kappa\_trunc - kappa(k)) \textcolor{keywordflow}{then} 1657 dk(k) = -kappa(k) 1658 kappa(k) = 0.0 1659 \textcolor{keywordflow}{if} ((k < ks\_src) .and. (k+1 > ks\_kappa)) ks\_kappa = k+1 1660 \textcolor{keywordflow}{elseif} (dk(k) < 2.0*kappa\_trunc - kappa(k)) \textcolor{keywordflow}{then} 1661 dk(k) = 2.0*dk(k) - (2.0*kappa\_trunc - kappa(k)) 1662 kappa(k) = max(kappa(k) + dk(k), 0.0) \textcolor{comment}{! The max is for paranoia.} 1663 \textcolor{keywordflow}{if} (k<=ks\_kappa) ks\_kappa = 2 1664 \textcolor{keywordflow}{else} 1665 kappa(k) = kappa(k) + dk(k) 1666 \textcolor{keywordflow}{if} (k<=ks\_kappa) ks\_kappa = 2 1667 \textcolor{keywordflow}{ endif} 1668 \textcolor{keywordflow}{ enddo} 1669 dq(1) = max(dq(1) + cq(2)*dq(2) + dqmdk(2) * dk(2), tke\_min - tke(1)) 1670 tke(1) = max(tke(1) + dq(1), tke\_min) 1671 dk(1) = 0.0 1672 \textcolor{keywordflow}{ endif} 1673 1674 \textcolor{comment}{! Check these solutions for consistency.} 1675 \textcolor{comment}{! The unit conversions here have not been carefully tested.} 1676 \textcolor{keywordflow}{if} (debug\_soln) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{do} k=2,nz 1677 \textcolor{comment}{! In these equations, K\_err\_lin and Q\_err\_lin should be at round-off levels} 1678 \textcolor{comment}{! compared with the dominant terms, perhaps, dz\_Int*I\_Ld2*kappa and} 1679 \textcolor{comment}{! dz\_Int*TKE\_decay*TKE. The exception is where, either 1) the decay term has been} 1680 \textcolor{comment}{! been increased to ensure a positive pivot, or 2) negative TKEs have been} 1681 \textcolor{comment}{! truncated, or 3) small or negative kappas have been rounded toward 0.} 1682 i\_q = 1.0 / tke(k) 1683 i\_ld2\_debug(k) = (n2(k)*ilambda2 + f2) * i\_q + i\_l2\_bdry(k) 1684 1685 kap\_src = dz\_int(k) * (k\_src(k) - i\_ld2(k)*kappa\_prev(k)) + & 1686 (idz(k-1)*(kappa\_prev(k-1)-kappa\_prev(k)) - & 1687 idz(k)*(kappa\_prev(k)-kappa\_prev(k+1))) 1688 k\_err\_lin = -idz(k-1)*(dk(k-1)-dk(k)) + idz(k)*(dk(k)-dk(k+1)) + & 1689 dz\_int(k)*i\_ld2\_debug(k)*dk(k) - kap\_src - & 1690 dz\_int(k)*(n2(k)*ilambda2 + f2)*i\_q**2*kappa\_prev(k) * dq(k) 1691 1692 tke\_src = dz\_int(k) * ((kappa\_prev(k) + kappa0)*s2(k) - & 1693 kappa\_prev(k)*n2(k) - (tke\_prev(k) - q0)*tke\_decay(k)) - & 1694 (aq(k) * (tke\_prev(k) - tke\_prev(k+1)) - aq(k-1) * (tke\_prev(k-1) - tke\_prev(k))) 1695 q\_err\_lin = tke\_src + (aq(k-1) * (dq(k-1)-dq(k)) - aq(k) * (dq(k)-dq(k+1))) - & 1696 0.5*(tke\_prev(k)-tke\_prev(k+1))*idz(k) * (dk(k) + dk(k+1)) - & 1697 0.5*(tke\_prev(k)-tke\_prev(k-1))*idz(k-1)* (dk(k-1) + dk(k)) + & 1698 dz\_int(k) * (dk(k) * (s2(k) - n2(k)) - dq(k)*tke\_decay(k)) 1699 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ endif} 1700 1701 \textcolor{keywordflow}{ endif} \textcolor{comment}{! End of the Newton's method solver.} 1702 1703 \textcolor{comment}{! Test kappa for convergence...} 1704 \textcolor{keywordflow}{if} ((tol\_err < newton\_err) .and. (.not.abort\_newton)) \textcolor{keywordflow}{then} 1705 \textcolor{comment}{! A lower tolerance is used to switch to Newton's method than to switch back.} 1706 newton\_test = newton\_err ; \textcolor{keywordflow}{if} (do\_newton) newton\_test = 2.0*newton\_err 1707 was\_newton = do\_newton 1708 within\_tolerance = .true. ; do\_newton = .true. 1709 \textcolor{keywordflow}{do} k=min(ks\_kappa,ks\_kappa\_prev),max(ke\_kappa,ke\_kappa\_prev) 1710 kappa\_mean = kappa0 + (kappa(k) - 0.5*dk(k)) 1711 \textcolor{keywordflow}{if} (abs(dk(k)) > newton\_test * kappa\_mean) \textcolor{keywordflow}{then} 1712 \textcolor{keywordflow}{if} (do\_newton) abort\_newton = .true. 1713 within\_tolerance = .false. ; do\_newton = .false. ; \textcolor{keywordflow}{exit} 1714 \textcolor{keywordflow}{elseif} (abs(dk(k)) > tol\_err * kappa\_mean) \textcolor{keywordflow}{then} 1715 within\_tolerance = .false. ; \textcolor{keywordflow}{if} (.not.do\_newton) \textcolor{keywordflow}{exit} 1716 \textcolor{keywordflow}{ endif} 1717 \textcolor{keywordflow}{if} (abs(dq(k)) > newton\_test*(tke(k) - 0.5*dq(k))) \textcolor{keywordflow}{then} 1718 \textcolor{keywordflow}{if} (do\_newton) abort\_newton = .true. 1719 do\_newton = .false. ; \textcolor{keywordflow}{if} (.not.within\_tolerance) \textcolor{keywordflow}{exit} 1720 \textcolor{keywordflow}{ endif} 1721 \textcolor{keywordflow}{ enddo} 1722 1723 \textcolor{keywordflow}{else} \textcolor{comment}{! Newton's method will not be used again, so no need to check.} 1724 within\_tolerance = .true. 1725 \textcolor{keywordflow}{do} k=min(ks\_kappa,ks\_kappa\_prev),max(ke\_kappa,ke\_kappa\_prev) 1726 \textcolor{keywordflow}{if} (abs(dk(k)) > tol\_err * (kappa0 + (kappa(k) - 0.5*dk(k)))) \textcolor{keywordflow}{then} 1727 within\_tolerance = .false. ; \textcolor{keywordflow}{exit} 1728 \textcolor{keywordflow}{ endif} 1729 \textcolor{keywordflow}{ enddo} 1730 \textcolor{keywordflow}{ endif} 1731 1732 \textcolor{keywordflow}{if} (abort\_newton) \textcolor{keywordflow}{then} 1733 do\_newton = .false. ; abort\_newton = .false. 1734 \textcolor{comment}{! We went to Newton too quickly last time, so restrict the tolerance.} 1735 newton\_err = 0.5*newton\_err 1736 ke\_kappa\_prev = nz 1737 \textcolor{keywordflow}{do} k=2,nz ; k\_q(k) = kappa(k) / max(tke(k), tke\_min) ;\textcolor{keywordflow}{ enddo} 1738 \textcolor{keywordflow}{ endif} 1739 1740 \textcolor{keywordflow}{if} (within\_tolerance) \textcolor{keywordflow}{exit} 1741 1742 \textcolor{keywordflow}{ enddo} 1743 1744 \textcolor{keywordflow}{if} (do\_newton) \textcolor{keywordflow}{then} \textcolor{comment}{! K\_Q needs to be calculated.} 1745 \textcolor{keywordflow}{do} k=1,ks\_kappa-1 ; k\_q(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 1746 \textcolor{keywordflow}{do} k=ks\_kappa,ke\_kappa ; k\_q(k) = kappa(k) / tke(k) ;\textcolor{keywordflow}{ enddo} 1747 \textcolor{keywordflow}{do} k=ke\_kappa+1,nz+1 ; k\_q(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 1748 \textcolor{keywordflow}{ endif} 1749 1750 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(local\_src)) \textcolor{keywordflow}{then} 1751 local\_src(1) = 0.0 ; local\_src(nz+1) = 0.0 1752 \textcolor{keywordflow}{do} k=2,nz 1753 diffusive\_src = idz(k-1)*(kappa(k-1)-kappa(k)) + idz(k)*(kappa(k+1)-kappa(k)) 1754 chg\_by\_k0 = kappa0 * ((idz(k-1)+idz(k)) / dz\_int(k) + i\_ld2(k)) 1755 \textcolor{keywordflow}{if} (diffusive\_src <= 0.0) \textcolor{keywordflow}{then} 1756 local\_src(k) = k\_src(k) + chg\_by\_k0 1757 \textcolor{keywordflow}{else} 1758 local\_src(k) = (k\_src(k) + chg\_by\_k0) + diffusive\_src / dz\_int(k) 1759 \textcolor{keywordflow}{ endif} 1760 \textcolor{keywordflow}{ enddo} 1761 \textcolor{keywordflow}{ endif} 1762 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(kappa\_src)) \textcolor{keywordflow}{then} 1763 kappa\_src(1) = 0.0 ; kappa\_src(nz+1) = 0.0 1764 \textcolor{keywordflow}{do} k=2,nz 1765 kappa\_src(k) = k\_src(k) 1766 \textcolor{keywordflow}{ enddo} 1767 \textcolor{keywordflow}{ endif} 1768 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__kappa__shear_ad4d87b0928aea195213e682b493eb555}\label{namespacemom__kappa__shear_ad4d87b0928aea195213e682b493eb555}} \index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}!kappa\+\_\+shear\+\_\+at\+\_\+vertex@{kappa\+\_\+shear\+\_\+at\+\_\+vertex}} \index{kappa\+\_\+shear\+\_\+at\+\_\+vertex@{kappa\+\_\+shear\+\_\+at\+\_\+vertex}!mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsubsection{\texorpdfstring{kappa\+\_\+shear\+\_\+at\+\_\+vertex()}{kappa\_shear\_at\_vertex()}} {\footnotesize\ttfamily logical function, public mom\+\_\+kappa\+\_\+shear\+::kappa\+\_\+shear\+\_\+at\+\_\+vertex (\begin{DoxyParamCaption}\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file }\end{DoxyParamCaption})} This function indicates to other modules whether the Jackson et al shear mixing parameterization will be used at the vertices without needing to duplicate the log entry. It returns false if the Jackson et al scheme is not used or if it is used via calculations at the tracer points. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em param\+\_\+file} & A structure to parse for run-\/time parameters \\ \hline \end{DoxyParams} Definition at line 1975 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 1975 \textcolor{keywordtype}{type}(param\_file\_type), \textcolor{keywordtype}{intent(in)} :: param\_file\textcolor{comment}{ !< A structure to parse for run-time parameters} 1976 1977 \textcolor{comment}{! Local variables} 1978 \textcolor{keywordtype}{character(len=40)} :: mdl = \textcolor{stringliteral}{"MOM\_kappa\_shear"} \textcolor{comment}{! This module's name.} 1979 \textcolor{keywordtype}{logical} :: do\_kappa\_shear 1980 \textcolor{comment}{! This function returns true only if the parameters "USE\_JACKSON\_PARAM" and "VERTEX\_SHEAR" are both true.} 1981 1982 kappa\_shear\_at\_vertex = .false. 1983 1984 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_JACKSON\_PARAM"}, do\_kappa\_shear, & 1985 default=.false., do\_not\_log=.true.) 1986 \textcolor{keywordflow}{if} (do\_kappa\_shear) & 1987 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"VERTEX\_SHEAR"}, kappa\_shear\_at\_vertex, & 1988 \textcolor{stringliteral}{"If true, do the calculations of the shear-driven mixing "}//& 1989 \textcolor{stringliteral}{"at the cell vertices (i.e., the vorticity points)."}, & 1990 default=.false., do\_not\_log=.true.) 1991 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__kappa__shear_a26cc5bb15545f04cfaf07e53410e09ec}\label{namespacemom__kappa__shear_a26cc5bb15545f04cfaf07e53410e09ec}} \index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}!kappa\+\_\+shear\+\_\+column@{kappa\+\_\+shear\+\_\+column}} \index{kappa\+\_\+shear\+\_\+column@{kappa\+\_\+shear\+\_\+column}!mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsubsection{\texorpdfstring{kappa\+\_\+shear\+\_\+column()}{kappa\_shear\_column()}} {\footnotesize\ttfamily subroutine mom\+\_\+kappa\+\_\+shear\+::kappa\+\_\+shear\+\_\+column (\begin{DoxyParamCaption}\item[{real, dimension( gv \%ke+1), intent(inout)}]{kappa, }\item[{real, dimension( gv \%ke+1), intent(out)}]{tke, }\item[{real, intent(in)}]{dt, }\item[{integer, intent(in)}]{nzc, }\item[{real, intent(in)}]{f2, }\item[{real, intent(in)}]{surface\+\_\+pres, }\item[{real, dimension( gv \%ke), intent(in)}]{dz, }\item[{real, dimension( gv \%ke), intent(in)}]{u0xdz, }\item[{real, dimension( gv \%ke), intent(in)}]{v0xdz, }\item[{real, dimension( gv \%ke), intent(in)}]{T0xdz, }\item[{real, dimension( gv \%ke), intent(in)}]{S0xdz, }\item[{real, dimension( gv \%ke+1), intent(out)}]{kappa\+\_\+avg, }\item[{real, dimension( gv \%ke+1), intent(out)}]{tke\+\_\+avg, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{type(\mbox{\hyperlink{structmom__kappa__shear_1_1kappa__shear__cs}{kappa\+\_\+shear\+\_\+cs}}), pointer}]{CS, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US, }\item[{real, dimension( gv \%ke+1), intent(out), optional}]{I\+\_\+\+Ld2\+\_\+1d, }\item[{real, dimension( gv \%ke+1), intent(out), optional}]{dz\+\_\+\+Int\+\_\+1d }\end{DoxyParamCaption})\hspace{0.3cm}{\ttfamily [private]}} This subroutine calculates shear-\/driven diffusivity and T\+KE in a single column. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\ \hline \mbox{\tt in,out} & {\em kappa} & The time-\/weighted average of kappa \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \hline \mbox{\tt out} & {\em tke} & The Turbulent Kinetic Energy per unit mass at\\ \hline \mbox{\tt in} & {\em nzc} & The number of active layers in the column.\\ \hline \mbox{\tt in} & {\em f2} & The square of the Coriolis parameter \mbox{[}T-\/2 $\sim$$>$ s-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em surface\+\_\+pres} & The surface pressure \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}.\\ \hline \mbox{\tt in} & {\em dz} & The layer thickness \mbox{[}Z $\sim$$>$ m\mbox{]}.\\ \hline \mbox{\tt in} & {\em u0xdz} & The initial zonal velocity times dz \mbox{[}Z L T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em v0xdz} & The initial meridional velocity times dz \mbox{[}Z L T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \hline \mbox{\tt in} & {\em t0xdz} & The initial temperature times dz \mbox{[}degC Z $\sim$$>$ degC m\mbox{]}.\\ \hline \mbox{\tt in} & {\em s0xdz} & The initial salinity times dz \mbox{[}ppt Z $\sim$$>$ ppt m\mbox{]}.\\ \hline \mbox{\tt out} & {\em kappa\+\_\+avg} & The time-\/weighted average of kappa \mbox{[}Z2 T-\/1 $\sim$$>$ m2 s-\/1\mbox{]}.\\ \hline \mbox{\tt out} & {\em tke\+\_\+avg} & The time-\/weighted average of T\+KE \mbox{[}Z2 T-\/2 $\sim$$>$ m2 s-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em dt} & Time increment \mbox{[}T $\sim$$>$ s\mbox{]}.\\ \hline \mbox{\tt in} & {\em tv} & A structure containing pointers to any available thermodynamic fields. Absent fields have N\+U\+LL ptrs.\\ \hline & {\em cs} & The control structure returned by a previous call to kappa\+\_\+shear\+\_\+init.\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt out} & {\em i\+\_\+ld2\+\_\+1d} & The inverse of the squared mixing length \mbox{[}Z-\/2 $\sim$$>$ m-\/2\mbox{]}.\\ \hline \mbox{\tt out} & {\em dz\+\_\+int\+\_\+1d} & The extent of a finite-\/volume space surrounding an interface, \\ \hline \end{DoxyParams} Definition at line 627 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 627 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< The ocean's vertical grid structure.} 628 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)}, & 629 \textcolor{keywordtype}{intent(inout)} :: kappa\textcolor{comment}{ !< The time-weighted average of kappa [Z2 T-1 ~> m2 s-1].} 630 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)}, & 631 \textcolor{keywordtype}{intent(out)} :: tke\textcolor{comment}{ !< The Turbulent Kinetic Energy per unit mass at} 632 \textcolor{comment}{ !! an interface [Z2 T-2 ~> m2 s-2].} 633 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{intent(in)} :: nzc\textcolor{comment}{ !< The number of active layers in the column.} 634 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: f2\textcolor{comment}{ !< The square of the Coriolis parameter [T-2 ~> s-2].} 635 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: surface\_pres\textcolor{comment}{ !< The surface pressure [R L2 T-2 ~> Pa].} 636 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))}, & 637 \textcolor{keywordtype}{intent(in)} :: dz\textcolor{comment}{ !< The layer thickness [Z ~> m].} 638 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))}, & 639 \textcolor{keywordtype}{intent(in)} :: u0xdz\textcolor{comment}{ !< The initial zonal velocity times dz [Z L T-1 ~> m2 s-1].} 640 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))}, & 641 \textcolor{keywordtype}{intent(in)} :: v0xdz\textcolor{comment}{ !< The initial meridional velocity times dz [Z L T-1 ~> m2 s-1].} 642 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))}, & 643 \textcolor{keywordtype}{intent(in)} :: T0xdz\textcolor{comment}{ !< The initial temperature times dz [degC Z ~> degC m].} 644 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV))}, & 645 \textcolor{keywordtype}{intent(in)} :: S0xdz\textcolor{comment}{ !< The initial salinity times dz [ppt Z ~> ppt m].} 646 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)}, & 647 \textcolor{keywordtype}{intent(out)} :: kappa\_avg\textcolor{comment}{ !< The time-weighted average of kappa [Z2 T-1 ~> m2 s-1].} 648 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)}, & 649 \textcolor{keywordtype}{intent(out)} :: tke\_avg\textcolor{comment}{ !< The time-weighted average of TKE [Z2 T-2 ~> m2 s-2].} 650 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: dt\textcolor{comment}{ !< Time increment [T ~> s].} 651 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< A structure containing pointers to any} 652 \textcolor{comment}{ !! available thermodynamic fields. Absent fields} 653 \textcolor{comment}{ !! have NULL ptrs.} 654 \textcolor{keywordtype}{type}(Kappa\_shear\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< The control structure returned by a previous} 655 \textcolor{comment}{ !! call to kappa\_shear\_init.} 656 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 657 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)}, & 658 \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: I\_Ld2\_1d\textcolor{comment}{ !< The inverse of the squared mixing length [Z-2 ~> m-2].} 659 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZK\_(GV)+1)}, & 660 \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: dz\_Int\_1d\textcolor{comment}{ !< The extent of a finite-volume space surrounding an interface,} 661 \textcolor{comment}{ !! as used in calculating kappa and TKE [Z ~> m].} 662 663 \textcolor{comment}{! Local variables} 664 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nzc)} :: & 665 u, & \textcolor{comment}{! The zonal velocity after a timestep of mixing [L T-1 ~> m s-1].} 666 v, & \textcolor{comment}{! The meridional velocity after a timestep of mixing [L T-1 ~> m s-1].} 667 Idz, & \textcolor{comment}{! The inverse of the distance between TKE points [Z-1 ~> m-1].} 668 T, & \textcolor{comment}{! The potential temperature after a timestep of mixing [degC].} 669 Sal, & \textcolor{comment}{! The salinity after a timestep of mixing [ppt].} 670 u\_test, v\_test, & \textcolor{comment}{! Temporary velocities [L T-1 ~> m s-1].} 671 T\_test, S\_test \textcolor{comment}{! Temporary temperatures [degC] and salinities [ppt].} 672 673 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(nzc+1)} :: & 674 N2, & \textcolor{comment}{! The squared buoyancy frequency at an interface [T-2 ~> s-2].} 675 dz\_Int, & \textcolor{comment}{! The extent of a finite-volume space surrounding an interface,} 676 \textcolor{comment}{! as used in calculating kappa and TKE [Z ~> m].} 677 i\_dz\_int, & \textcolor{comment}{! The inverse of the distance between velocity & density points} 678 \textcolor{comment}{! above and below an interface [Z-1 ~> m-1]. This is used to} 679 \textcolor{comment}{! calculate N2, shear, and fluxes, and it might differ from} 680 \textcolor{comment}{! 1/dz\_Int, as they have different uses.} 681 s2, & \textcolor{comment}{! The squared shear at an interface [T-2 ~> s-2].} 682 a1, & \textcolor{comment}{! a1 is the coupling between adjacent interfaces in the TKE,} 683 \textcolor{comment}{! velocity, and density equations [Z s-1 ~> m s-1] or [Z ~> m]} 684 c1, & \textcolor{comment}{! c1 is used in the tridiagonal (and similar) solvers.} 685 k\_src, & \textcolor{comment}{! The shear-dependent source term in the kappa equation [T-1 ~> s-1].} 686 kappa\_src, & \textcolor{comment}{! The shear-dependent source term in the kappa equation [T-1 ~> s-1].} 687 kappa\_out, & \textcolor{comment}{! The kappa that results from the kappa equation [Z2 T-1 ~> m2 s-1].} 688 kappa\_mid, & \textcolor{comment}{! The average of the initial and predictor estimates of kappa [Z2 T-1 ~> m2 s-1].} 689 tke\_pred, & \textcolor{comment}{! The value of TKE from a predictor step [Z2 T-2 ~> m2 s-2].} 690 kappa\_pred, & \textcolor{comment}{! The value of kappa from a predictor step [Z2 T-1 ~> m2 s-1].} 691 pressure, & \textcolor{comment}{! The pressure at an interface [R L2 T-2 ~> Pa].} 692 t\_int, & \textcolor{comment}{! The temperature interpolated to an interface [degC].} 693 sal\_int, & \textcolor{comment}{! The salinity interpolated to an interface [ppt].} 694 dbuoy\_dt, & \textcolor{comment}{! The partial derivatives of buoyancy with changes in temperature} 695 dbuoy\_ds, & \textcolor{comment}{! and salinity, [Z T-2 degC-1 ~> m s-2 degC-1] and [Z T-2 ppt-1 ~> m s-2 ppt-1].} 696 i\_l2\_bdry, & \textcolor{comment}{! The inverse of the square of twice the harmonic mean} 697 \textcolor{comment}{! distance to the top and bottom boundaries [Z-2 ~> m-2].} 698 k\_q, & \textcolor{comment}{! Diffusivity divided by TKE [T ~> s].} 699 k\_q\_tmp, & \textcolor{comment}{! A temporary copy of diffusivity divided by TKE [T ~> s].} 700 local\_src\_avg, & \textcolor{comment}{! The time-integral of the local source [nondim].} 701 tol\_min, & \textcolor{comment}{! Minimum tolerated ksrc for the corrector step [T-1 ~> s-1].} 702 tol\_max, & \textcolor{comment}{! Maximum tolerated ksrc for the corrector step [T-1 ~> s-1].} 703 tol\_chg, & \textcolor{comment}{! The tolerated kappa change integrated over a timestep [nondim].} 704 dist\_from\_top, & \textcolor{comment}{! The distance from the top surface [Z ~> m].} 705 local\_src \textcolor{comment}{! The sum of all sources of kappa, including kappa\_src and} 706 \textcolor{comment}{! sources from the elliptic term [T-1 ~> s-1].} 707 708 \textcolor{keywordtype}{real} :: dist\_from\_bot \textcolor{comment}{! The distance from the bottom surface [Z ~> m].} 709 \textcolor{keywordtype}{real} :: b1 \textcolor{comment}{! The inverse of the pivot in the tridiagonal equations.} 710 \textcolor{keywordtype}{real} :: bd1 \textcolor{comment}{! A term in the denominator of b1.} 711 \textcolor{keywordtype}{real} :: d1 \textcolor{comment}{! 1 - c1 in the tridiagonal equations.} 712 \textcolor{keywordtype}{real} :: gR0 \textcolor{comment}{! A conversion factor from Z to pressure, given by Rho\_0 times g} 713 \textcolor{comment}{! [R L2 T-2 Z-1 ~> kg m-2 s-2].} 714 \textcolor{keywordtype}{real} :: g\_R0 \textcolor{comment}{! g\_R0 is a rescaled version of g/Rho [Z R-1 T-2 ~> m4 kg-1 s-2].} 715 \textcolor{keywordtype}{real} :: Norm \textcolor{comment}{! A factor that normalizes two weights to 1 [Z-2 ~> m-2].} 716 \textcolor{keywordtype}{real} :: tol\_dksrc \textcolor{comment}{! Tolerance for the change in the kappa source within an iteration} 717 \textcolor{comment}{! relative to the local source [nondim]. This must be greater than 1.} 718 \textcolor{keywordtype}{real} :: tol2 \textcolor{comment}{! The tolerance for the change in the kappa source within an iteration} 719 \textcolor{comment}{! relative to the average local source over previous iterations [nondim].} 720 \textcolor{keywordtype}{real} :: tol\_dksrc\_low \textcolor{comment}{! The tolerance for the fractional decrease in ksrc} 721 \textcolor{comment}{! within an iteration [nondim]. 0 < tol\_dksrc\_low < 1.} 722 \textcolor{keywordtype}{real} :: Ri\_crit \textcolor{comment}{! The critical shear Richardson number for shear-} 723 \textcolor{comment}{! driven mixing [nondim]. The theoretical value is 0.25.} 724 \textcolor{keywordtype}{real} :: dt\_rem \textcolor{comment}{! The remaining time to advance the solution [T ~> s].} 725 \textcolor{keywordtype}{real} :: dt\_now \textcolor{comment}{! The time step used in the current iteration [T ~> s].} 726 \textcolor{keywordtype}{real} :: dt\_wt \textcolor{comment}{! The fractional weight of the current iteration [nondim].} 727 \textcolor{keywordtype}{real} :: dt\_test \textcolor{comment}{! A time-step that is being tested for whether it} 728 \textcolor{comment}{! gives acceptably small changes in k\_src [T ~> s].} 729 \textcolor{keywordtype}{real} :: Idtt \textcolor{comment}{! Idtt = 1 / dt\_test [T-1 ~> s-1].} 730 \textcolor{keywordtype}{real} :: dt\_inc \textcolor{comment}{! An increment to dt\_test that is being tested [T ~> s].} 731 732 \textcolor{keywordtype}{real} :: k0dt \textcolor{comment}{! The background diffusivity times the timestep [Z2 ~> m2].} 733 \textcolor{keywordtype}{logical} :: valid\_dt \textcolor{comment}{! If true, all levels so far exhibit acceptably small changes in k\_src.} 734 \textcolor{keywordtype}{logical} :: use\_temperature \textcolor{comment}{! If true, temperature and salinity have been} 735 \textcolor{comment}{! allocated and are being used as state variables.} 736 \textcolor{keywordtype}{integer} :: ks\_kappa, ke\_kappa \textcolor{comment}{! The k-range with nonzero kappas.} 737 \textcolor{keywordtype}{integer} :: dt\_halvings \textcolor{comment}{! The number of times that the time-step is halved} 738 \textcolor{comment}{! in seeking an acceptable timestep. If none is} 739 \textcolor{comment}{! found, dt\_rem*0.5^dt\_halvings is used.} 740 \textcolor{keywordtype}{integer} :: dt\_refinements \textcolor{comment}{! The number of 2-fold refinements that will be used} 741 \textcolor{comment}{! to estimate the maximum permitted time step. I.e.,} 742 \textcolor{comment}{! the resolution is 1/2^dt\_refinements.} 743 \textcolor{keywordtype}{integer} :: k, itt, itt\_dt 744 745 \textcolor{comment}{! This calculation of N2 is for debugging only.} 746 \textcolor{comment}{! real, dimension(SZK\_(GV)+1) :: &} 747 \textcolor{comment}{! N2\_debug, & ! A version of N2 for debugging [T-2 ~> s-2]} 748 749 ri\_crit = cs%Rino\_crit 750 gr0 = gv%Rho0 * gv%g\_Earth 751 g\_r0 = (us%L\_to\_Z**2 * gv%g\_Earth) / (gv%Rho0) 752 k0dt = dt*cs%kappa\_0 753 754 tol\_dksrc = cs%kappa\_src\_max\_chg 755 \textcolor{keywordflow}{if} (tol\_dksrc == 10.0) \textcolor{keywordflow}{then} 756 \textcolor{comment}{! This is equivalent to the expression below, but avoids changes at roundoff for the default value.} 757 tol\_dksrc\_low = 0.95 758 \textcolor{keywordflow}{else} 759 tol\_dksrc\_low = (tol\_dksrc - 0.5)/tol\_dksrc 760 \textcolor{keywordflow}{ endif} 761 tol2 = 2.0*cs%kappa\_tol\_err 762 dt\_refinements = 5 \textcolor{comment}{! Selected so that 1/2^dt\_refinements < 1-tol\_dksrc\_low} 763 use\_temperature = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%T)) use\_temperature = .true. 764 765 766 \textcolor{comment}{! Set up Idz as the inverse of layer thicknesses.} 767 \textcolor{keywordflow}{do} k=1,nzc ; idz(k) = 1.0 / dz(k) ;\textcolor{keywordflow}{ enddo} 768 \textcolor{comment}{! Set up I\_dz\_int as the inverse of the distance between} 769 \textcolor{comment}{! adjacent layer centers.} 770 i\_dz\_int(1) = 2.0 / dz(1) 771 dist\_from\_top(1) = 0.0 772 \textcolor{keywordflow}{do} k=2,nzc 773 i\_dz\_int(k) = 2.0 / (dz(k-1) + dz(k)) 774 dist\_from\_top(k) = dist\_from\_top(k-1) + dz(k-1) 775 \textcolor{keywordflow}{ enddo} 776 i\_dz\_int(nzc+1) = 2.0 / dz(nzc) 777 778 \textcolor{comment}{! Determine the velocities and thicknesses after eliminating massless} 779 \textcolor{comment}{! layers and applying a time-step of background diffusion.} 780 \textcolor{keywordflow}{if} (nzc > 1) \textcolor{keywordflow}{then} 781 a1(2) = k0dt*i\_dz\_int(2) 782 b1 = 1.0 / (dz(1) + a1(2)) 783 u(1) = b1 * u0xdz(1) ; v(1) = b1 * v0xdz(1) 784 t(1) = b1 * t0xdz(1) ; sal(1) = b1 * s0xdz(1) 785 c1(2) = a1(2) * b1 ; d1 = dz(1) * b1 \textcolor{comment}{! = 1 - c1} 786 \textcolor{keywordflow}{do} k=2,nzc-1 787 bd1 = dz(k) + d1*a1(k) 788 a1(k+1) = k0dt*i\_dz\_int(k+1) 789 b1 = 1.0 / (bd1 + a1(k+1)) 790 u(k) = b1 * (u0xdz(k) + a1(k)*u(k-1)) 791 v(k) = b1 * (v0xdz(k) + a1(k)*v(k-1)) 792 t(k) = b1 * (t0xdz(k) + a1(k)*t(k-1)) 793 sal(k) = b1 * (s0xdz(k) + a1(k)*sal(k-1)) 794 c1(k+1) = a1(k+1) * b1 ; d1 = bd1 * b1 \textcolor{comment}{! d1 = 1 - c1} 795 \textcolor{keywordflow}{ enddo} 796 \textcolor{comment}{! rho or T and S have insulating boundary conditions, u & v use no-slip} 797 \textcolor{comment}{! bottom boundary conditions (if kappa0 > 0).} 798 \textcolor{comment}{! For no-slip bottom boundary conditions} 799 b1 = 1.0 / ((dz(nzc) + d1*a1(nzc)) + k0dt*i\_dz\_int(nzc+1)) 800 u(nzc) = b1 * (u0xdz(nzc) + a1(nzc)*u(nzc-1)) 801 v(nzc) = b1 * (v0xdz(nzc) + a1(nzc)*v(nzc-1)) 802 \textcolor{comment}{! For insulating boundary conditions} 803 b1 = 1.0 / (dz(nzc) + d1*a1(nzc)) 804 t(nzc) = b1 * (t0xdz(nzc) + a1(nzc)*t(nzc-1)) 805 sal(nzc) = b1 * (s0xdz(nzc) + a1(nzc)*sal(nzc-1)) 806 \textcolor{keywordflow}{do} k=nzc-1,1,-1 807 u(k) = u(k) + c1(k+1)*u(k+1) ; v(k) = v(k) + c1(k+1)*v(k+1) 808 t(k) = t(k) + c1(k+1)*t(k+1) ; sal(k) = sal(k) + c1(k+1)*sal(k+1) 809 \textcolor{keywordflow}{ enddo} 810 \textcolor{keywordflow}{else} 811 \textcolor{comment}{! This is correct, but probably unnecessary.} 812 b1 = 1.0 / (dz(1) + k0dt*i\_dz\_int(2)) 813 u(1) = b1 * u0xdz(1) ; v(1) = b1 * v0xdz(1) 814 b1 = 1.0 / dz(1) 815 t(1) = b1 * t0xdz(1) ; sal(1) = b1 * s0xdz(1) 816 \textcolor{keywordflow}{ endif} 817 818 \textcolor{comment}{! This uses half the harmonic mean of thicknesses to provide two estimates} 819 \textcolor{comment}{! of the boundary between cells, and the inverse of the harmonic mean to} 820 \textcolor{comment}{! weight the two estimates. The net effect is that interfaces around thin} 821 \textcolor{comment}{! layers have thin cells, and the total thickness adds up properly.} 822 \textcolor{comment}{! The top- and bottom- interfaces have zero thickness, consistent with} 823 \textcolor{comment}{! adding additional zero thickness layers.} 824 dz\_int(1) = 0.0 ; dz\_int(2) = dz(1) 825 \textcolor{keywordflow}{do} k=2,nzc-1 826 norm = 1.0 / (dz(k)*(dz(k-1)+dz(k+1)) + 2.0*dz(k-1)*dz(k+1)) 827 dz\_int(k) = dz\_int(k) + dz(k) * ( ((dz(k)+dz(k+1)) * dz(k-1)) * norm) 828 dz\_int(k+1) = dz(k) * ( ((dz(k-1)+dz(k)) * dz(k+1)) * norm) 829 \textcolor{keywordflow}{ enddo} 830 dz\_int(nzc) = dz\_int(nzc) + dz(nzc) ; dz\_int(nzc+1) = 0.0 831 832 dist\_from\_bot = 0.0 833 \textcolor{keywordflow}{do} k=nzc,2,-1 834 dist\_from\_bot = dist\_from\_bot + dz(k) 835 i\_l2\_bdry(k) = (dist\_from\_top(k) + dist\_from\_bot)**2 / & 836 (dist\_from\_top(k) * dist\_from\_bot)**2 837 \textcolor{keywordflow}{ enddo} 838 839 \textcolor{comment}{! Calculate thermodynamic coefficients and an initial estimate of N2.} 840 \textcolor{keywordflow}{if} (use\_temperature) \textcolor{keywordflow}{then} 841 pressure(1) = surface\_pres 842 \textcolor{keywordflow}{do} k=2,nzc 843 pressure(k) = pressure(k-1) + gr0*dz(k-1) 844 t\_int(k) = 0.5*(t(k-1) + t(k)) 845 sal\_int(k) = 0.5*(sal(k-1) + sal(k)) 846 \textcolor{keywordflow}{ enddo} 847 \textcolor{keyword}{call }calculate\_density\_derivs(t\_int, sal\_int, pressure, dbuoy\_dt, dbuoy\_ds, & 848 tv%eqn\_of\_state, (/2,nzc/), scale=-g\_r0 ) 849 \textcolor{keywordflow}{else} 850 \textcolor{keywordflow}{do} k=1,nzc+1 ; dbuoy\_dt(k) = -g\_r0 ; dbuoy\_ds(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 851 \textcolor{keywordflow}{ endif} 852 853 \textcolor{comment}{! N2\_debug(1) = 0.0 ; N2\_debug(nzc+1) = 0.0} 854 \textcolor{comment}{! do K=2,nzc} 855 \textcolor{comment}{! N2\_debug(K) = max((dbuoy\_dT(K) * (T0xdz(k-1)*Idz(k-1) - T0xdz(k)*Idz(k)) + &} 856 \textcolor{comment}{! dbuoy\_dS(K) * (S0xdz(k-1)*Idz(k-1) - S0xdz(k)*Idz(k))) * &} 857 \textcolor{comment}{! I\_dz\_int(K), 0.0)} 858 \textcolor{comment}{! enddo} 859 860 \textcolor{comment}{! This call just calculates N2 and S2.} 861 \textcolor{keyword}{call }calculate\_projected\_state(kappa, u, v, t, sal, 0.0, nzc, dz, i\_dz\_int, & 862 dbuoy\_dt, dbuoy\_ds, u, v, t, sal, gv, us, & 863 n2=n2, s2=s2, vel\_underflow=cs%vel\_underflow) 864 \textcolor{comment}{! ----------------------------------------------------} 865 \textcolor{comment}{! Iterate} 866 \textcolor{comment}{! ----------------------------------------------------} 867 dt\_rem = dt 868 \textcolor{keywordflow}{do} k=1,nzc+1 869 k\_q(k) = 0.0 870 kappa\_avg(k) = 0.0 ; tke\_avg(k) = 0.0 871 local\_src\_avg(k) = 0.0 872 \textcolor{comment}{! Use the grid spacings to scale errors in the source.} 873 \textcolor{keywordflow}{if} ( dz\_int(k) > 0.0 ) & 874 local\_src\_avg(k) = 0.1 * k0dt * i\_dz\_int(k) / dz\_int(k) 875 \textcolor{keywordflow}{ enddo} 876 877 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_setup)} 878 879 \textcolor{comment}{! do itt=1,CS%max\_RiNo\_it} 880 \textcolor{keywordflow}{do} itt=1,cs%max\_KS\_it 881 882 \textcolor{comment}{! ----------------------------------------------------} 883 \textcolor{comment}{! Calculate new values of u, v, rho, N^2 and S.} 884 \textcolor{comment}{! ----------------------------------------------------} 885 886 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_KQ)} 887 \textcolor{keyword}{call }find\_kappa\_tke(n2, s2, kappa, idz, dz\_int, i\_l2\_bdry, f2, & 888 nzc, cs, gv, us, k\_q, tke, kappa\_out, kappa\_src, local\_src) 889 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_KQ)} 890 891 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_avg)} 892 \textcolor{comment}{! Determine the range of non-zero values of kappa\_out.} 893 ks\_kappa = gv%ke+1 ; ke\_kappa = 0 894 \textcolor{keywordflow}{do} k=2,nzc ; \textcolor{keywordflow}{if} (kappa\_out(k) > 0.0) \textcolor{keywordflow}{then} 895 ks\_kappa = k ; \textcolor{keywordflow}{exit} 896 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 897 \textcolor{keywordflow}{do} k=nzc,ks\_kappa,-1 ; \textcolor{keywordflow}{if} (kappa\_out(k) > 0.0) \textcolor{keywordflow}{then} 898 ke\_kappa = k ; \textcolor{keywordflow}{exit} 899 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ enddo} 900 \textcolor{keywordflow}{if} (ke\_kappa == nzc) kappa\_out(nzc+1) = 0.0 901 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_avg)} 902 903 \textcolor{comment}{! Determine how long to use this value of kappa (dt\_now).} 904 905 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_project)} 906 \textcolor{keywordflow}{if} ((ke\_kappa < ks\_kappa) .or. (itt==cs%max\_KS\_it)) \textcolor{keywordflow}{then} 907 dt\_now = dt\_rem 908 \textcolor{keywordflow}{else} 909 \textcolor{comment}{! Limit dt\_now so that |k\_src(k)-kappa\_src(k)| < tol * local\_src(k)} 910 dt\_test = dt\_rem 911 \textcolor{keywordflow}{do} k=2,nzc 912 tol\_max(k) = kappa\_src(k) + tol\_dksrc * local\_src(k) 913 tol\_min(k) = kappa\_src(k) - tol\_dksrc\_low * local\_src(k) 914 tol\_chg(k) = tol2 * local\_src\_avg(k) 915 \textcolor{keywordflow}{ enddo} 916 917 \textcolor{keywordflow}{do} itt\_dt=1,(cs%max\_KS\_it+1-itt)/2 918 \textcolor{comment}{! The maximum number of times that the time-step is halved in} 919 \textcolor{comment}{! seeking an acceptable timestep is reduced with each iteration,} 920 \textcolor{comment}{! so that as the maximum number of iterations is approached, the} 921 \textcolor{comment}{! whole remaining timestep is used. Typically, an acceptable} 922 \textcolor{comment}{! timestep is found long before the minimum is reached, so the} 923 \textcolor{comment}{! value of max\_KS\_it may be unimportant, especially if it is large} 924 \textcolor{comment}{! enough.} 925 \textcolor{keyword}{call }calculate\_projected\_state(kappa\_out, u, v, t, sal, 0.5*dt\_test, nzc, dz, i\_dz\_int, & 926 dbuoy\_dt, dbuoy\_ds, u\_test, v\_test, t\_test, s\_test, & 927 gv, us, n2, s2, ks\_int=ks\_kappa, ke\_int=ke\_kappa, & 928 vel\_underflow=cs%vel\_underflow) 929 valid\_dt = .true. 930 idtt = 1.0 / dt\_test 931 \textcolor{keywordflow}{do} k=max(ks\_kappa-1,2),min(ke\_kappa+1,nzc) 932 \textcolor{keywordflow}{if} (n2(k) < ri\_crit * s2(k)) \textcolor{keywordflow}{then} \textcolor{comment}{! Equivalent to Ri < Ri\_crit.} 933 k\_src(k) = (2.0 * cs%Shearmix\_rate * sqrt(s2(k))) * & 934 ((ri\_crit*s2(k) - n2(k)) / (ri\_crit*s2(k) + cs%FRi\_curvature*n2(k))) 935 \textcolor{keywordflow}{if} (cs%restrictive\_tolerance\_check) \textcolor{keywordflow}{then} 936 \textcolor{keywordflow}{if} ((k\_src(k) > min(tol\_max(k), kappa\_src(k) + idtt*tol\_chg(k))) .or. & 937 (k\_src(k) < max(tol\_min(k), kappa\_src(k) - idtt*tol\_chg(k)))) \textcolor{keywordflow}{then} 938 valid\_dt = .false. ; \textcolor{keywordflow}{exit} 939 \textcolor{keywordflow}{ endif} 940 \textcolor{keywordflow}{else} 941 \textcolor{keywordflow}{if} ((k\_src(k) > max(tol\_max(k), kappa\_src(k) + idtt*tol\_chg(k))) .or. & 942 (k\_src(k) < min(tol\_min(k), kappa\_src(k) - idtt*tol\_chg(k)))) \textcolor{keywordflow}{then} 943 valid\_dt = .false. ; \textcolor{keywordflow}{exit} 944 \textcolor{keywordflow}{ endif} 945 \textcolor{keywordflow}{ endif} 946 \textcolor{keywordflow}{else} 947 \textcolor{keywordflow}{if} (0.0 < min(tol\_min(k), kappa\_src(k) - idtt*tol\_chg(k))) \textcolor{keywordflow}{then} 948 valid\_dt = .false. ; k\_src(k) = 0.0 ; \textcolor{keywordflow}{exit} 949 \textcolor{keywordflow}{ endif} 950 \textcolor{keywordflow}{ endif} 951 \textcolor{keywordflow}{ enddo} 952 953 \textcolor{keywordflow}{if} (valid\_dt) \textcolor{keywordflow}{exit} 954 dt\_test = 0.5*dt\_test 955 \textcolor{keywordflow}{ enddo} 956 \textcolor{keywordflow}{if} ((dt\_test < dt\_rem) .and. valid\_dt) \textcolor{keywordflow}{then} 957 dt\_inc = 0.5*dt\_test 958 \textcolor{keywordflow}{do} itt\_dt=1,dt\_refinements 959 \textcolor{keyword}{call }calculate\_projected\_state(kappa\_out, u, v, t, sal, 0.5*(dt\_test+dt\_inc), & 960 nzc, dz, i\_dz\_int, dbuoy\_dt, dbuoy\_ds, u\_test, v\_test, t\_test, s\_test, & 961 gv, us, n2, s2, ks\_int=ks\_kappa, ke\_int=ke\_kappa, vel\_underflow=cs%vel\_underflow) 962 valid\_dt = .true. 963 idtt = 1.0 / (dt\_test+dt\_inc) 964 \textcolor{keywordflow}{do} k=max(ks\_kappa-1,2),min(ke\_kappa+1,nzc) 965 \textcolor{keywordflow}{if} (n2(k) < ri\_crit * s2(k)) \textcolor{keywordflow}{then} \textcolor{comment}{! Equivalent to Ri < Ri\_crit.} 966 k\_src(k) = (2.0 * cs%Shearmix\_rate * sqrt(s2(k))) * & 967 ((ri\_crit*s2(k) - n2(k)) / (ri\_crit*s2(k) + cs%FRi\_curvature*n2(k))) 968 \textcolor{keywordflow}{if} ((k\_src(k) > max(tol\_max(k), kappa\_src(k) + idtt*tol\_chg(k))) .or. & 969 (k\_src(k) < min(tol\_min(k), kappa\_src(k) - idtt*tol\_chg(k)))) \textcolor{keywordflow}{then} 970 valid\_dt = .false. ; \textcolor{keywordflow}{exit} 971 \textcolor{keywordflow}{ endif} 972 \textcolor{keywordflow}{else} 973 \textcolor{keywordflow}{if} (0.0 < min(tol\_min(k), kappa\_src(k) - idtt*tol\_chg(k))) \textcolor{keywordflow}{then} 974 valid\_dt = .false. ; k\_src(k) = 0.0 ; \textcolor{keywordflow}{exit} 975 \textcolor{keywordflow}{ endif} 976 \textcolor{keywordflow}{ endif} 977 \textcolor{keywordflow}{ enddo} 978 979 \textcolor{keywordflow}{if} (valid\_dt) dt\_test = dt\_test + dt\_inc 980 dt\_inc = 0.5*dt\_inc 981 \textcolor{keywordflow}{ enddo} 982 \textcolor{keywordflow}{else} 983 dt\_inc = 0.0 984 \textcolor{keywordflow}{ endif} 985 986 dt\_now = min(dt\_test*(1.0+cs%kappa\_tol\_err)+dt\_inc, dt\_rem) 987 \textcolor{keywordflow}{do} k=2,nzc 988 local\_src\_avg(k) = local\_src\_avg(k) + dt\_now * local\_src(k) 989 \textcolor{keywordflow}{ enddo} 990 \textcolor{keywordflow}{ endif} \textcolor{comment}{! Are all the values of kappa\_out 0?} 991 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_project)} 992 993 \textcolor{comment}{! The state has already been projected forward. Now find new values of kappa.} 994 995 \textcolor{keywordflow}{if} (ke\_kappa < ks\_kappa) \textcolor{keywordflow}{then} 996 \textcolor{comment}{! There is no mixing now, and will not be again.} 997 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_avg)} 998 dt\_wt = dt\_rem / dt ; dt\_rem = 0.0 999 \textcolor{keywordflow}{do} k=1,nzc+1 1000 kappa\_mid(k) = 0.0 1001 \textcolor{comment}{! This would be here but does nothing.} 1002 \textcolor{comment}{! kappa\_avg(K) = kappa\_avg(K) + kappa\_mid(K)*dt\_wt} 1003 tke\_avg(k) = tke\_avg(k) + dt\_wt*tke(k) 1004 \textcolor{keywordflow}{ enddo} 1005 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_avg)} 1006 \textcolor{keywordflow}{else} 1007 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_project)} 1008 \textcolor{keyword}{call }calculate\_projected\_state(kappa\_out, u, v, t, sal, dt\_now, nzc, dz, i\_dz\_int, & 1009 dbuoy\_dt, dbuoy\_ds, u\_test, v\_test, t\_test, s\_test, & 1010 gv, us, n2=n2, s2=s2, ks\_int=ks\_kappa, ke\_int=ke\_kappa, & 1011 vel\_underflow=cs%vel\_underflow) 1012 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_project)} 1013 1014 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_KQ)} 1015 \textcolor{keywordflow}{do} k=1,nzc+1 ; k\_q\_tmp(k) = k\_q(k) ;\textcolor{keywordflow}{ enddo} 1016 \textcolor{keyword}{call }find\_kappa\_tke(n2, s2, kappa\_out, idz, dz\_int, i\_l2\_bdry, f2, & 1017 nzc, cs, gv, us, k\_q\_tmp, tke\_pred, kappa\_pred) 1018 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_KQ)} 1019 1020 ks\_kappa = gv%ke+1 ; ke\_kappa = 0 1021 \textcolor{keywordflow}{do} k=1,nzc+1 1022 kappa\_mid(k) = 0.5*(kappa\_out(k) + kappa\_pred(k)) 1023 \textcolor{keywordflow}{if} ((kappa\_mid(k) > 0.0) .and. (k 0.0) ke\_kappa = k 1025 \textcolor{keywordflow}{ enddo} 1026 1027 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_project)} 1028 \textcolor{keyword}{call }calculate\_projected\_state(kappa\_mid, u, v, t, sal, dt\_now, nzc, dz, i\_dz\_int, & 1029 dbuoy\_dt, dbuoy\_ds, u\_test, v\_test, t\_test, s\_test, & 1030 gv, us, n2=n2, s2=s2, ks\_int=ks\_kappa, ke\_int=ke\_kappa, & 1031 vel\_underflow=cs%vel\_underflow) 1032 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_project)} 1033 1034 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_KQ)} 1035 \textcolor{keyword}{call }find\_kappa\_tke(n2, s2, kappa\_out, idz, dz\_int, i\_l2\_bdry, f2, & 1036 nzc, cs, gv, us, k\_q, tke\_pred, kappa\_pred) 1037 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_KQ)} 1038 1039 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_avg)} 1040 dt\_wt = dt\_now / dt ; dt\_rem = dt\_rem - dt\_now 1041 \textcolor{keywordflow}{do} k=1,nzc+1 1042 kappa\_mid(k) = 0.5*(kappa\_out(k) + kappa\_pred(k)) 1043 kappa\_avg(k) = kappa\_avg(k) + kappa\_mid(k)*dt\_wt 1044 tke\_avg(k) = tke\_avg(k) + dt\_wt*0.5*(tke\_pred(k) + tke(k)) 1045 kappa(k) = kappa\_pred(k) \textcolor{comment}{! First guess for the next iteration.} 1046 \textcolor{keywordflow}{ enddo} 1047 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_avg)} 1048 \textcolor{keywordflow}{ endif} 1049 1050 \textcolor{keywordflow}{if} (dt\_rem > 0.0) \textcolor{keywordflow}{then} 1051 \textcolor{comment}{! Update the values of u, v, T, Sal, N2, and S2 for the next iteration.} 1052 \textcolor{comment}{! call cpu\_clock\_begin(id\_clock\_project)} 1053 \textcolor{keyword}{call }calculate\_projected\_state(kappa\_mid, u, v, t, sal, dt\_now, nzc, & 1054 dz, i\_dz\_int, dbuoy\_dt, dbuoy\_ds, u, v, t, sal, & 1055 gv, us, n2, s2, vel\_underflow=cs%vel\_underflow) 1056 \textcolor{comment}{! call cpu\_clock\_end(id\_clock\_project)} 1057 \textcolor{keywordflow}{ endif} 1058 1059 \textcolor{keywordflow}{if} (dt\_rem <= 0.0) \textcolor{keywordflow}{exit} 1060 1061 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! end itt loop} 1062 1063 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(i\_ld2\_1d)) \textcolor{keywordflow}{then} 1064 \textcolor{keywordflow}{do} k=1,gv%ke+1 ; i\_ld2\_1d(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 1065 \textcolor{keywordflow}{do} k=2,nzc ; \textcolor{keywordflow}{if} (tke(k) > 0.0) & 1066 i\_ld2\_1d(k) = i\_l2\_bdry(k) + (n2(k) / cs%lambda**2 + f2) / tke(k) 1067 \textcolor{keywordflow}{ enddo} 1068 \textcolor{keywordflow}{ endif} 1069 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(dz\_int\_1d)) \textcolor{keywordflow}{then} 1070 \textcolor{keywordflow}{do} k=1,nzc+1 ; dz\_int\_1d(k) = dz\_int(k) ;\textcolor{keywordflow}{ enddo} 1071 \textcolor{keywordflow}{do} k=nzc+2,gv%ke ; dz\_int\_1d(k) = 0.0 ;\textcolor{keywordflow}{ enddo} 1072 \textcolor{keywordflow}{ endif} 1073 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__kappa__shear_a82428168ba463cf7276ae96d3e32475c}\label{namespacemom__kappa__shear_a82428168ba463cf7276ae96d3e32475c}} \index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}!kappa\+\_\+shear\+\_\+init@{kappa\+\_\+shear\+\_\+init}} \index{kappa\+\_\+shear\+\_\+init@{kappa\+\_\+shear\+\_\+init}!mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsubsection{\texorpdfstring{kappa\+\_\+shear\+\_\+init()}{kappa\_shear\_init()}} {\footnotesize\ttfamily logical function, public mom\+\_\+kappa\+\_\+shear\+::kappa\+\_\+shear\+\_\+init (\begin{DoxyParamCaption}\item[{type(time\+\_\+type), intent(in)}]{Time, }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US, }\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file, }\item[{type(diag\+\_\+ctrl), intent(inout), target}]{diag, }\item[{type(\mbox{\hyperlink{structmom__kappa__shear_1_1kappa__shear__cs}{kappa\+\_\+shear\+\_\+cs}}), pointer}]{CS }\end{DoxyParamCaption})} This subroutine initializes the parameters that regulate shear-\/driven mixing. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em time} & The current model time.\\ \hline \mbox{\tt in} & {\em g} & The ocean\textquotesingle{}s grid structure.\\ \hline \mbox{\tt in} & {\em gv} & The ocean\textquotesingle{}s vertical grid structure.\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em 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 & {\em cs} & A pointer that is set to point to the control structure for this module\\ \hline \end{DoxyParams} \begin{DoxyReturn}{Returns} True if module is to be used, False otherwise \end{DoxyReturn} Definition at line 1773 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 1773 \textcolor{keywordtype}{type}(time\_type), \textcolor{keywordtype}{intent(in)} :: Time\textcolor{comment}{ !< The current model time.} 1774 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{ !< The ocean's grid structure.} 1775 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< The ocean's vertical grid structure.} 1776 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 1777 \textcolor{keywordtype}{type}(param\_file\_type), \textcolor{keywordtype}{intent(in)} :: param\_file\textcolor{comment}{ !< A structure to parse for run-time} 1778 \textcolor{comment}{ !! parameters.} 1779 \textcolor{keywordtype}{type}(diag\_ctrl), \textcolor{keywordtype}{target}, \textcolor{keywordtype}{intent(inout)} :: diag\textcolor{comment}{ !< A structure that is used to regulate diagnostic} 1780 \textcolor{comment}{ !! output.} 1781 \textcolor{keywordtype}{type}(Kappa\_shear\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< A pointer that is set to point to the control} 1782 \textcolor{comment}{ !! structure for this module} 1783 \textcolor{keywordtype}{logical} :: kappa\_shear\_init\textcolor{comment}{ !< True if module is to be used, False otherwise} 1784 1785 \textcolor{comment}{! Local variables} 1786 \textcolor{keywordtype}{logical} :: merge\_mixedlayer 1787 \textcolor{keywordtype}{logical} :: debug\_shear 1788 \textcolor{keywordtype}{logical} :: just\_read \textcolor{comment}{! If true, this module is not used, so only read the parameters.} 1789 \textcolor{comment}{! This include declares and sets the variable "version".} 1790 \textcolor{preprocessor}{# include "version\_variable.h"} 1791 \textcolor{preprocessor}{} \textcolor{keywordtype}{character(len=40)} :: mdl = \textcolor{stringliteral}{"MOM\_kappa\_shear"} \textcolor{comment}{! This module's name.} 1792 \textcolor{keywordtype}{real} :: kappa\_0\_unscaled \textcolor{comment}{! The value of kappa\_0 in MKS units [m2 s-1]} 1793 \textcolor{keywordtype}{real} :: KD\_normal \textcolor{comment}{! The KD of the main model, read here only as a parameter} 1794 \textcolor{comment}{! for setting the default of KD\_SMOOTH in MKS units [m2 s-1]} 1795 1796 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keywordflow}{then} 1797 \textcolor{keyword}{call }mom\_error(warning, \textcolor{stringliteral}{"kappa\_shear\_init called with an associated "}// & 1798 \textcolor{stringliteral}{"control structure."}) 1799 \textcolor{keywordflow}{return} 1800 \textcolor{keywordflow}{ endif} 1801 \textcolor{keyword}{allocate}(cs) 1802 1803 \textcolor{comment}{! The Jackson-Hallberg-Legg shear mixing parameterization uses the following} 1804 \textcolor{comment}{! 6 nondimensional coefficients. That paper gives 3 best fit parameter sets.} 1805 \textcolor{comment}{! Ri\_Crit Rate FRi\_Curv K\_buoy TKE\_N TKE\_Shear} 1806 \textcolor{comment}{! p1: 0.25 0.089 -0.97 0.82 0.24 0.14} 1807 \textcolor{comment}{! p2: 0.30 0.085 -0.94 0.86 0.26 0.13} 1808 \textcolor{comment}{! p3: 0.35 0.088 -0.89 0.81 0.28 0.12} 1809 \textcolor{comment}{! Future research will reveal how these should be modified to take} 1810 \textcolor{comment}{! subgridscale inhomogeneity into account.} 1811 1812 \textcolor{comment}{! Set default, read and log parameters} 1813 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_JACKSON\_PARAM"}, kappa\_shear\_init, default=.false., do\_not\_log=.true. ) 1814 \textcolor{keyword}{call }log\_version(param\_file, mdl, version, & 1815 \textcolor{stringliteral}{"Parameterization of shear-driven turbulence following Jackson, Hallberg and Legg, JPO 2008"}, & 1816 log\_to\_all=.true., debugging=kappa\_shear\_init, all\_default=.not.kappa\_shear\_init) 1817 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_JACKSON\_PARAM"}, kappa\_shear\_init, & 1818 \textcolor{stringliteral}{"If true, use the Jackson-Hallberg-Legg (JPO 2008) "}//& 1819 \textcolor{stringliteral}{"shear mixing parameterization."}, default=.false.) 1820 just\_read = .not.kappa\_shear\_init 1821 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"VERTEX\_SHEAR"}, cs%KS\_at\_vertex, & 1822 \textcolor{stringliteral}{"If true, do the calculations of the shear-driven mixing "}//& 1823 \textcolor{stringliteral}{"at the cell vertices (i.e., the vorticity points)."}, & 1824 default=.false., do\_not\_log=just\_read) 1825 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"RINO\_CRIT"}, cs%RiNo\_crit, & 1826 \textcolor{stringliteral}{"The critical Richardson number for shear mixing."}, & 1827 units=\textcolor{stringliteral}{"nondim"}, default=0.25, do\_not\_log=just\_read) 1828 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"SHEARMIX\_RATE"}, cs%Shearmix\_rate, & 1829 \textcolor{stringliteral}{"A nondimensional rate scale for shear-driven entrainment. "}//& 1830 \textcolor{stringliteral}{"Jackson et al find values in the range of 0.085-0.089."}, & 1831 units=\textcolor{stringliteral}{"nondim"}, default=0.089, do\_not\_log=just\_read) 1832 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MAX\_RINO\_IT"}, cs%max\_RiNo\_it, & 1833 \textcolor{stringliteral}{"The maximum number of iterations that may be used to "}//& 1834 \textcolor{stringliteral}{"estimate the Richardson number driven mixing."}, & 1835 units=\textcolor{stringliteral}{"nondim"}, default=50, do\_not\_log=just\_read) 1836 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD"}, kd\_normal, default=0.0, do\_not\_log=.true.) 1837 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD\_KAPPA\_SHEAR\_0"}, cs%kappa\_0, & 1838 \textcolor{stringliteral}{"The background diffusivity that is used to smooth the "}//& 1839 \textcolor{stringliteral}{"density and shear profiles before solving for the "}//& 1840 \textcolor{stringliteral}{"diffusivities. The default is the greater of KD and 1e-7 m2 s-1."}, & 1841 units=\textcolor{stringliteral}{"m2 s-1"}, default=max(kd\_normal, 1.0e-7), scale=us%m2\_s\_to\_Z2\_T, & 1842 unscaled=kappa\_0\_unscaled, do\_not\_log=just\_read) 1843 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KD\_TRUNC\_KAPPA\_SHEAR"}, cs%kappa\_trunc, & 1844 \textcolor{stringliteral}{"The value of shear-driven diffusivity that is considered negligible "}//& 1845 \textcolor{stringliteral}{"and is rounded down to 0. The default is 1% of KD\_KAPPA\_SHEAR\_0."}, & 1846 units=\textcolor{stringliteral}{"m2 s-1"}, default=0.01*kappa\_0\_unscaled, scale=us%m2\_s\_to\_Z2\_T, do\_not\_log=just\_read ) 1847 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"FRI\_CURVATURE"}, cs%FRi\_curvature, & 1848 \textcolor{stringliteral}{"The nondimensional curvature of the function of the "}//& 1849 \textcolor{stringliteral}{"Richardson number in the kappa source term in the "}//& 1850 \textcolor{stringliteral}{"Jackson et al. scheme."}, units=\textcolor{stringliteral}{"nondim"}, default=-0.97, do\_not\_log=just\_read) 1851 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TKE\_N\_DECAY\_CONST"}, cs%C\_N, & 1852 \textcolor{stringliteral}{"The coefficient for the decay of TKE due to "}//& 1853 \textcolor{stringliteral}{"stratification (i.e. proportional to N*tke). "}//& 1854 \textcolor{stringliteral}{"The values found by Jackson et al. are 0.24-0.28."}, & 1855 units=\textcolor{stringliteral}{"nondim"}, default=0.24, do\_not\_log=just\_read) 1856 \textcolor{comment}{! call get\_param(param\_file, mdl, "LAYER\_KAPPA\_STAGGER", CS%layer\_stagger, &} 1857 \textcolor{comment}{! default=.false., do\_not\_log=just\_read)} 1858 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TKE\_SHEAR\_DECAY\_CONST"}, cs%C\_S, & 1859 \textcolor{stringliteral}{"The coefficient for the decay of TKE due to shear (i.e. "}//& 1860 \textcolor{stringliteral}{"proportional to |S|*tke). The values found by Jackson "}//& 1861 \textcolor{stringliteral}{"et al. are 0.14-0.12."}, units=\textcolor{stringliteral}{"nondim"}, default=0.14, do\_not\_log=just\_read) 1862 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_BUOY\_SCALE\_COEF"}, cs%lambda, & 1863 \textcolor{stringliteral}{"The coefficient for the buoyancy length scale in the "}//& 1864 \textcolor{stringliteral}{"kappa equation. The values found by Jackson et al. are "}//& 1865 \textcolor{stringliteral}{"in the range of 0.81-0.86."}, units=\textcolor{stringliteral}{"nondim"}, default=0.82, do\_not\_log=just\_read) 1866 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_N\_OVER\_S\_SCALE\_COEF2"}, cs%lambda2\_N\_S, & 1867 \textcolor{stringliteral}{"The square of the ratio of the coefficients of the "}//& 1868 \textcolor{stringliteral}{"buoyancy and shear scales in the diffusivity equation, "}//& 1869 \textcolor{stringliteral}{"Set this to 0 (the default) to eliminate the shear scale. "}//& 1870 \textcolor{stringliteral}{"This is only used if USE\_JACKSON\_PARAM is true."}, & 1871 units=\textcolor{stringliteral}{"nondim"}, default=0.0, do\_not\_log=just\_read) 1872 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_SHEAR\_TOL\_ERR"}, cs%kappa\_tol\_err, & 1873 \textcolor{stringliteral}{"The fractional error in kappa that is tolerated. "}//& 1874 \textcolor{stringliteral}{"Iteration stops when changes between subsequent "}//& 1875 \textcolor{stringliteral}{"iterations are smaller than this everywhere in a "}//& 1876 \textcolor{stringliteral}{"column. The peak diffusivities usually converge most "}//& 1877 \textcolor{stringliteral}{"rapidly, and have much smaller errors than this."}, & 1878 units=\textcolor{stringliteral}{"nondim"}, default=0.1, do\_not\_log=just\_read) 1879 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TKE\_BACKGROUND"}, cs%TKE\_bg, & 1880 \textcolor{stringliteral}{"A background level of TKE used in the first iteration "}//& 1881 \textcolor{stringliteral}{"of the kappa equation. TKE\_BACKGROUND could be 0."}, & 1882 units=\textcolor{stringliteral}{"m2 s-2"}, default=0.0, scale=us%m\_to\_Z**2*us%T\_to\_s**2) 1883 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_SHEAR\_ELIM\_MASSLESS"}, cs%eliminate\_massless, & 1884 \textcolor{stringliteral}{"If true, massless layers are merged with neighboring "}//& 1885 \textcolor{stringliteral}{"massive layers in this calculation. The default is "}//& 1886 \textcolor{stringliteral}{"true and I can think of no good reason why it should "}//& 1887 \textcolor{stringliteral}{"be false. This is only used if USE\_JACKSON\_PARAM is true."}, & 1888 default=.true., do\_not\_log=just\_read) 1889 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MAX\_KAPPA\_SHEAR\_IT"}, cs%max\_KS\_it, & 1890 \textcolor{stringliteral}{"The maximum number of iterations that may be used to "}//& 1891 \textcolor{stringliteral}{"estimate the time-averaged diffusivity."}, units=\textcolor{stringliteral}{"nondim"}, & 1892 default=13, do\_not\_log=just\_read) 1893 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"PRANDTL\_TURB"}, cs%Prandtl\_turb, & 1894 \textcolor{stringliteral}{"The turbulent Prandtl number applied to shear instability."}, & 1895 units=\textcolor{stringliteral}{"nondim"}, default=1.0, do\_not\_log=.true.) 1896 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"VEL\_UNDERFLOW"}, cs%vel\_underflow, & 1897 \textcolor{stringliteral}{"A negligibly small velocity magnitude below which velocity components are set "}//& 1898 \textcolor{stringliteral}{"to 0. A reasonable value might be 1e-30 m/s, which is less than an "}//& 1899 \textcolor{stringliteral}{"Angstrom divided by the age of the universe."}, & 1900 units=\textcolor{stringliteral}{"m s-1"}, default=0.0, scale=us%m\_s\_to\_L\_T, do\_not\_log=just\_read) 1901 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_SHEAR\_MAX\_KAP\_SRC\_CHG"}, cs%kappa\_src\_max\_chg, & 1902 \textcolor{stringliteral}{"The maximum permitted increase in the kappa source within an iteration relative "}//& 1903 \textcolor{stringliteral}{"to the local source; this must be greater than 1. The lower limit for the "}//& 1904 \textcolor{stringliteral}{"permitted fractional decrease is (1 - 0.5/kappa\_src\_max\_chg). These limits "}//& 1905 \textcolor{stringliteral}{"could perhaps be made dynamic with an improved iterative solver."}, & 1906 default=10.0, units=\textcolor{stringliteral}{"nondim"}, do\_not\_log=just\_read) 1907 1908 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DEBUG"}, cs%debug, & 1909 \textcolor{stringliteral}{"If true, write out verbose debugging data."}, & 1910 default=.false., debuggingparam=.true., do\_not\_log=just\_read) 1911 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"DEBUG\_KAPPA\_SHEAR"}, debug\_shear, & 1912 \textcolor{stringliteral}{"If true, write debugging data for the kappa-shear code."}, & 1913 default=.false., debuggingparam=.true., do\_not\_log=.true.) 1914 \textcolor{keywordflow}{if} (debug\_shear) cs%debug = .true. 1915 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_SHEAR\_VERTEX\_PSURF\_BUG"}, cs%psurf\_bug, & 1916 \textcolor{stringliteral}{"If true, do a simple average of the cell surface pressures to get a pressure "}//& 1917 \textcolor{stringliteral}{"at the corner if VERTEX\_SHEAR=True. Otherwise mask out any land points in "}//& 1918 \textcolor{stringliteral}{"the average."}, default=.true., do\_not\_log=(just\_read .or. (.not.cs%KS\_at\_vertex))) 1919 1920 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_SHEAR\_ITER\_BUG"}, cs%dKdQ\_iteration\_bug, & 1921 \textcolor{stringliteral}{"If true, use an older, dimensionally inconsistent estimate of the "}//& 1922 \textcolor{stringliteral}{"derivative of diffusivity with energy in the Newton's method iteration. "}//& 1923 \textcolor{stringliteral}{"The bug causes undercorrections when dz > 1 m."}, default=.false., do\_not\_log=just\_read) 1924 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_SHEAR\_ALL\_LAYER\_TKE\_BUG"}, cs%all\_layer\_TKE\_bug, & 1925 \textcolor{stringliteral}{"If true, report back the latest estimate of TKE instead of the time average "}//& 1926 \textcolor{stringliteral}{"TKE when there is mass in all layers. Otherwise always report the time "}//& 1927 \textcolor{stringliteral}{"averaged TKE, as is currently done when there are some massless layers."}, & 1928 default=.false., do\_not\_log=just\_read) 1929 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_RESTRICTIVE\_TOLERANCE\_CHECK"}, cs%restrictive\_tolerance\_check, & 1930 \textcolor{stringliteral}{"If true, uses the more restrictive tolerance check to determine if a timestep "}//& 1931 \textcolor{stringliteral}{"is acceptable for the KS\_it outer iteration loop. False uses the original less "}//& 1932 \textcolor{stringliteral}{"restrictive check."}, default=.false., do\_not\_log=just\_read) 1933 \textcolor{comment}{! id\_clock\_KQ = cpu\_clock\_id('Ocean KS kappa\_shear', grain=CLOCK\_ROUTINE)} 1934 \textcolor{comment}{! id\_clock\_avg = cpu\_clock\_id('Ocean KS avg', grain=CLOCK\_ROUTINE)} 1935 \textcolor{comment}{! id\_clock\_project = cpu\_clock\_id('Ocean KS project', grain=CLOCK\_ROUTINE)} 1936 \textcolor{comment}{! id\_clock\_setup = cpu\_clock\_id('Ocean KS setup', grain=CLOCK\_ROUTINE)} 1937 1938 cs%nkml = 1 1939 \textcolor{keywordflow}{if} (gv%nkml>0) \textcolor{keywordflow}{then} 1940 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"KAPPA\_SHEAR\_MERGE\_ML"},merge\_mixedlayer, & 1941 \textcolor{stringliteral}{"If true, combine the mixed layers together before solving the "}//& 1942 \textcolor{stringliteral}{"kappa-shear equations."}, default=.true., do\_not\_log=just\_read) 1943 \textcolor{keywordflow}{if} (merge\_mixedlayer) cs%nkml = gv%nkml 1944 \textcolor{keywordflow}{ endif} 1945 1946 \textcolor{comment}{! Forego remainder of initialization if not using this scheme} 1947 \textcolor{keywordflow}{if} (.not. kappa\_shear\_init) \textcolor{keywordflow}{return} 1948 1949 cs%diag => diag 1950 1951 cs%id\_Kd\_shear = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'},\textcolor{stringliteral}{'Kd\_shear'}, diag%axesTi, time, & 1952 \textcolor{stringliteral}{'Shear-driven Diapycnal Diffusivity'}, \textcolor{stringliteral}{'m2 s-1'}, conversion=us%Z2\_T\_to\_m2\_s) 1953 cs%id\_TKE = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'},\textcolor{stringliteral}{'TKE\_shear'}, diag%axesTi, time, & 1954 \textcolor{stringliteral}{'Shear-driven Turbulent Kinetic Energy'}, \textcolor{stringliteral}{'m2 s-2'}, conversion=us%Z\_to\_m**2*us%s\_to\_T**2) 1955 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__kappa__shear_ac7859c609e462000ca8fd763d68d141e}\label{namespacemom__kappa__shear_ac7859c609e462000ca8fd763d68d141e}} \index{mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}!kappa\+\_\+shear\+\_\+is\+\_\+used@{kappa\+\_\+shear\+\_\+is\+\_\+used}} \index{kappa\+\_\+shear\+\_\+is\+\_\+used@{kappa\+\_\+shear\+\_\+is\+\_\+used}!mom\+\_\+kappa\+\_\+shear@{mom\+\_\+kappa\+\_\+shear}} \subsubsection{\texorpdfstring{kappa\+\_\+shear\+\_\+is\+\_\+used()}{kappa\_shear\_is\_used()}} {\footnotesize\ttfamily logical function, public mom\+\_\+kappa\+\_\+shear\+::kappa\+\_\+shear\+\_\+is\+\_\+used (\begin{DoxyParamCaption}\item[{type(param\+\_\+file\+\_\+type), intent(in)}]{param\+\_\+file }\end{DoxyParamCaption})} This function indicates to other modules whether the Jackson et al shear mixing parameterization will be used without needing to duplicate the log entry. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em param\+\_\+file} & A structure to parse for run-\/time parameters \\ \hline \end{DoxyParams} Definition at line 1961 of file M\+O\+M\+\_\+kappa\+\_\+shear.\+F90. \begin{DoxyCode} 1961 \textcolor{keywordtype}{type}(param\_file\_type), \textcolor{keywordtype}{intent(in)} :: param\_file\textcolor{comment}{ !< A structure to parse for run-time parameters} 1962 1963 \textcolor{comment}{! Local variables} 1964 \textcolor{keywordtype}{character(len=40)} :: mdl = \textcolor{stringliteral}{"MOM\_kappa\_shear"} \textcolor{comment}{! This module's name.} 1965 \textcolor{comment}{! This function reads the parameter "USE\_JACKSON\_PARAM" and returns its value.} 1966 1967 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_JACKSON\_PARAM"}, kappa\_shear\_is\_used, & 1968 default=.false., do\_not\_log=.true.) \end{DoxyCode}