\hypertarget{namespacemom__isopycnal__slopes}{}\section{mom\+\_\+isopycnal\+\_\+slopes Module Reference} \label{namespacemom__isopycnal__slopes}\index{mom\+\_\+isopycnal\+\_\+slopes@{mom\+\_\+isopycnal\+\_\+slopes}} \subsection{Detailed Description} Calculations of isoneutral slopes and stratification. \subsection*{Functions/\+Subroutines} \begin{DoxyCompactItemize} \item subroutine, public \hyperlink{namespacemom__isopycnal__slopes_a9a4cbf819be46d9babab62f7f09734c8}{calc\+\_\+isoneutral\+\_\+slopes} (G, GV, US, h, e, tv, dt\+\_\+kappa\+\_\+smooth, slope\+\_\+x, slope\+\_\+y, N2\+\_\+u, N2\+\_\+v, halo, O\+BC) \begin{DoxyCompactList}\small\item\em Calculate isopycnal slopes, and optionally return N2 used in calculation. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__isopycnal__slopes_a34691482caaff356da3c5182657dba0d}{vert\+\_\+fill\+\_\+ts} (h, T\+\_\+in, S\+\_\+in, kappa\+\_\+dt, T\+\_\+f, S\+\_\+f, G, GV, halo\+\_\+here, larger\+\_\+h\+\_\+denom) \begin{DoxyCompactList}\small\item\em Returns tracer arrays (nominally T and S) with massless layers filled with sensible values, by diffusing vertically with a small but constant diffusivity. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection{Function/\+Subroutine Documentation} \mbox{\Hypertarget{namespacemom__isopycnal__slopes_a9a4cbf819be46d9babab62f7f09734c8}\label{namespacemom__isopycnal__slopes_a9a4cbf819be46d9babab62f7f09734c8}} \index{mom\+\_\+isopycnal\+\_\+slopes@{mom\+\_\+isopycnal\+\_\+slopes}!calc\+\_\+isoneutral\+\_\+slopes@{calc\+\_\+isoneutral\+\_\+slopes}} \index{calc\+\_\+isoneutral\+\_\+slopes@{calc\+\_\+isoneutral\+\_\+slopes}!mom\+\_\+isopycnal\+\_\+slopes@{mom\+\_\+isopycnal\+\_\+slopes}} \subsubsection{\texorpdfstring{calc\+\_\+isoneutral\+\_\+slopes()}{calc\_isoneutral\_slopes()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+isopycnal\+\_\+slopes\+::calc\+\_\+isoneutral\+\_\+slopes (\begin{DoxyParamCaption}\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{type(unit\+\_\+scale\+\_\+type), intent(in)}]{US, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke+1), intent(in)}]{e, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{real, intent(in)}]{dt\+\_\+kappa\+\_\+smooth, }\item[{real, dimension( g \%isdb\+: g \%iedb, g \%jsd\+: g \%jed, g \%ke+1), intent(inout)}]{slope\+\_\+x, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsdb\+: g \%jedb, g \%ke+1), intent(inout)}]{slope\+\_\+y, }\item[{real, dimension( g \%isdb\+: g \%iedb, g \%jsd\+: g \%jed, g \%ke+1), intent(inout), optional}]{N2\+\_\+u, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsdb\+: g \%jedb, g \%ke+1), intent(inout), optional}]{N2\+\_\+v, }\item[{integer, intent(in), optional}]{halo, }\item[{type(ocean\+\_\+obc\+\_\+type), optional, pointer}]{O\+BC }\end{DoxyParamCaption})} Calculate isopycnal slopes, and optionally return N2 used in calculation. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & The ocean\textquotesingle{}s grid structure\\ \hline \mbox{\tt in} & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\ \hline \mbox{\tt in} & {\em e} & Interface heights \mbox{[}Z $\sim$$>$ m\mbox{]} or units given by 1/eta\+\_\+to\+\_\+m)\\ \hline \mbox{\tt in} & {\em tv} & A structure pointing to various thermodynamic variables\\ \hline \mbox{\tt in} & {\em dt\+\_\+kappa\+\_\+smooth} & A smoothing vertical diffusivity times a smoothing timescale \mbox{[}Z2 $\sim$$>$ m2\mbox{]}.\\ \hline \mbox{\tt in,out} & {\em slope\+\_\+x} & Isopycnal slope in i-\/direction \mbox{[}nondim\mbox{]}\\ \hline \mbox{\tt in,out} & {\em slope\+\_\+y} & Isopycnal slope in j-\/direction \mbox{[}nondim\mbox{]}\\ \hline \mbox{\tt in,out} & {\em n2\+\_\+u} & Brunt-\/\+Vaisala frequency squared at\\ \hline \mbox{\tt in,out} & {\em n2\+\_\+v} & Brunt-\/\+Vaisala frequency squared at\\ \hline \mbox{\tt in} & {\em halo} & Halo width over which to compute\\ \hline & {\em obc} & Open boundaries control structure. \\ \hline \end{DoxyParams} Definition at line 30 of file M\+O\+M\+\_\+isopycnal\+\_\+slopes.\+F90. \begin{DoxyCode} 30 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< The ocean's grid structure} 31 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< The ocean's vertical grid structure} 32 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 33 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2]} 34 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)}, \textcolor{keywordtype}{intent(in)} :: e\textcolor{comment}{ !< Interface heights [Z ~> m] or units} 35 \textcolor{comment}{ !! given by 1/eta\_to\_m)} 36 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< A structure pointing to various} 37 \textcolor{comment}{ !! thermodynamic variables} 38 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: dt\_kappa\_smooth\textcolor{comment}{ !< A smoothing vertical diffusivity} 39 \textcolor{comment}{ !! times a smoothing timescale [Z2 ~> m2].} 40 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G)+1)}, \textcolor{keywordtype}{intent(inout)} :: slope\_x\textcolor{comment}{ !< Isopycnal slope in i-direction [nondim]} 41 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G)+1)}, \textcolor{keywordtype}{intent(inout)} :: slope\_y\textcolor{comment}{ !< Isopycnal slope in j-direction [nondim]} 42 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G)+1)}, & 43 \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: n2\_u\textcolor{comment}{ !< Brunt-Vaisala frequency squared at} 44 \textcolor{comment}{ !! interfaces between u-points [T-2 ~> s-2]} 45 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G)+1)}, & 46 \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(inout)} :: n2\_v\textcolor{comment}{ !< Brunt-Vaisala frequency squared at} 47 \textcolor{comment}{ !! interfaces between u-points [T-2 ~> s-2]} 48 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: halo\textcolor{comment}{ !< Halo width over which to compute} 49 \textcolor{keywordtype}{type}(ocean\_obc\_type), \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: obc\textcolor{comment}{ !< Open boundaries control structure.} 50 51 \textcolor{comment}{! real, optional, intent(in) :: eta\_to\_m !< The conversion factor from the units} 52 \textcolor{comment}{! (This argument has been tested but for now serves no purpose.) !! of eta to m; US%Z\_to\_m by default.} 53 \textcolor{comment}{! Local variables} 54 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G), SZJ\_(G), SZK\_(G))} :: & 55 t, & \textcolor{comment}{! The temperature [degC], with the values in} 56 \textcolor{comment}{! in massless layers filled vertically by diffusion.} 57 s \textcolor{comment}{!, & ! The filled salinity [ppt], with the values in} 58 \textcolor{comment}{! in massless layers filled vertically by diffusion.} 59 \textcolor{comment}{! Rho ! Density itself, when a nonlinear equation of state is not in use [R ~> kg m-3].} 60 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G), SZJ\_(G), SZK\_(G)+1)} :: & 61 pres \textcolor{comment}{! The pressure at an interface [R L2 T-2 ~> Pa].} 62 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G))} :: & 63 drho\_dt\_u, & \textcolor{comment}{! The derivative of density with temperature at u points [R degC-1 ~> kg m-3 degC-1].} 64 drho\_ds\_u \textcolor{comment}{! The derivative of density with salinity at u points [R ppt-1 ~> kg m-3 ppt-1].} 65 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: & 66 drho\_dt\_v, & \textcolor{comment}{! The derivative of density with temperature at v points [R degC-1 ~> kg m-3 degC-1].} 67 drho\_ds\_v \textcolor{comment}{! The derivative of density with salinity at v points [R ppt-1 ~> kg m-3 ppt-1].} 68 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G))} :: & 69 t\_u, & \textcolor{comment}{! Temperature on the interface at the u-point [degC].} 70 s\_u, & \textcolor{comment}{! Salinity on the interface at the u-point [ppt].} 71 pres\_u \textcolor{comment}{! Pressure on the interface at the u-point [R L2 T-2 ~> Pa].} 72 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: & 73 t\_v, & \textcolor{comment}{! Temperature on the interface at the v-point [degC].} 74 s\_v, & \textcolor{comment}{! Salinity on the interface at the v-point [ppt].} 75 pres\_v \textcolor{comment}{! Pressure on the interface at the v-point [R L2 T-2 ~> Pa].} 76 \textcolor{keywordtype}{real} :: drdia, drdib \textcolor{comment}{! Along layer zonal- and meridional- potential density} 77 \textcolor{keywordtype}{real} :: drdja, drdjb \textcolor{comment}{! gradients in the layers above (A) and below (B) the} 78 \textcolor{comment}{! interface times the grid spacing [R ~> kg m-3].} 79 \textcolor{keywordtype}{real} :: drdkl, drdkr \textcolor{comment}{! Vertical density differences across an interface [R ~> kg m-3].} 80 \textcolor{keywordtype}{real} :: hg2a, hg2b \textcolor{comment}{! Squares of geometric mean thicknesses [H2 ~> m2 or kg2 m-4].} 81 \textcolor{keywordtype}{real} :: hg2l, hg2r \textcolor{comment}{! Squares of geometric mean thicknesses [H2 ~> m2 or kg2 m-4].} 82 \textcolor{keywordtype}{real} :: haa, hab, hal, har \textcolor{comment}{! Arithmetic mean thicknesses [H ~> m or kg m-2].} 83 \textcolor{keywordtype}{real} :: dzal, dzar \textcolor{comment}{! Temporary thicknesses in eta units [Z ~> m].} 84 \textcolor{keywordtype}{real} :: wta, wtb, wtl, wtr \textcolor{comment}{! Unscaled weights, with various units.} 85 \textcolor{keywordtype}{real} :: drdx, drdy \textcolor{comment}{! Zonal and meridional density gradients [R L-1 ~> kg m-4].} 86 \textcolor{keywordtype}{real} :: drdz \textcolor{comment}{! Vertical density gradient [R Z-1 ~> kg m-4].} 87 \textcolor{keywordtype}{real} :: slope \textcolor{comment}{! The slope of density surfaces, calculated in a way} 88 \textcolor{comment}{! that is always between -1 and 1.} 89 \textcolor{keywordtype}{real} :: mag\_grad2 \textcolor{comment}{! The squared magnitude of the 3-d density gradient [R2 L-2 ~> kg2 m-8].} 90 \textcolor{keywordtype}{real} :: slope2\_ratio \textcolor{comment}{! The ratio of the slope squared to slope\_max squared.} 91 \textcolor{keywordtype}{real} :: h\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 92 \textcolor{comment}{! in roundoff and can be neglected [H ~> m or kg m-2].} 93 \textcolor{keywordtype}{real} :: h\_neglect2 \textcolor{comment}{! h\_neglect^2 [H2 ~> m2 or kg2 m-4].} 94 \textcolor{keywordtype}{real} :: dz\_neglect \textcolor{comment}{! A change in interface heighs that is so small it is usually lost} 95 \textcolor{comment}{! in roundoff and can be neglected [Z ~> m].} 96 \textcolor{keywordtype}{logical} :: use\_eos \textcolor{comment}{! If true, density is calculated from T & S using an equation of state.} 97 \textcolor{keywordtype}{real} :: g\_rho0 \textcolor{comment}{! The gravitational acceleration divided by density [Z2 T-2 R-1 ~> m5 kg-2 s-2]} 98 \textcolor{keywordtype}{real} :: z\_to\_l \textcolor{comment}{! A conversion factor between from units for e to the} 99 \textcolor{comment}{! units for lateral distances.} 100 \textcolor{keywordtype}{real} :: l\_to\_z \textcolor{comment}{! A conversion factor between from units for lateral distances} 101 \textcolor{comment}{! to the units for e.} 102 \textcolor{keywordtype}{real} :: h\_to\_z \textcolor{comment}{! A conversion factor from thickness units to the units of e.} 103 104 \textcolor{keywordtype}{logical} :: present\_n2\_u, present\_n2\_v 105 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom\_u, eosdom\_v \textcolor{comment}{! Domains for the equation of state calculations at u and v points} 106 \textcolor{keywordtype}{integer} :: is, ie, js, je, nz, isdb 107 \textcolor{keywordtype}{integer} :: i, j, k 108 \textcolor{keywordtype}{integer} :: l\_seg 109 \textcolor{keywordtype}{logical} :: local\_open\_u\_bc, local\_open\_v\_bc 110 111 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(halo)) \textcolor{keywordflow}{then} 112 is = g%isc-halo ; ie = g%iec+halo ; js = g%jsc-halo ; je = g%jec+halo 113 \textcolor{keywordflow}{else} 114 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec 115 \textcolor{keywordflow}{ endif} 116 nz = g%ke ; isdb = g%IsdB 117 118 h\_neglect = gv%H\_subroundoff ; h\_neglect2 = h\_neglect**2 119 z\_to\_l = us%Z\_to\_L ; h\_to\_z = gv%H\_to\_Z 120 \textcolor{comment}{! if (present(eta\_to\_m)) then} 121 \textcolor{comment}{! Z\_to\_L = eta\_to\_m*US%m\_to\_L ; H\_to\_Z = GV%H\_to\_m / eta\_to\_m} 122 \textcolor{comment}{! endif} 123 l\_to\_z = 1.0 / z\_to\_l 124 dz\_neglect = gv%H\_subroundoff * h\_to\_z 125 126 local\_open\_u\_bc = .false. 127 local\_open\_v\_bc = .false. 128 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(obc)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(obc)) \textcolor{keywordflow}{then} 129 local\_open\_u\_bc = obc%open\_u\_BCs\_exist\_globally 130 local\_open\_v\_bc = obc%open\_v\_BCs\_exist\_globally 131 \textcolor{keywordflow}{ endif} ;\textcolor{keywordflow}{ endif} 132 133 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 134 135 present\_n2\_u = \textcolor{keyword}{PRESENT}(n2\_u) 136 present\_n2\_v = \textcolor{keyword}{PRESENT}(n2\_v) 137 g\_rho0 = (us%L\_to\_Z*l\_to\_z*gv%g\_Earth) / gv%Rho0 138 \textcolor{keywordflow}{if} (present\_n2\_u) \textcolor{keywordflow}{then} 139 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=is-1,ie 140 n2\_u(i,j,1) = 0. 141 n2\_u(i,j,nz+1) = 0. 142 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 143 \textcolor{keywordflow}{ endif} 144 \textcolor{keywordflow}{if} (present\_n2\_v) \textcolor{keywordflow}{then} 145 \textcolor{keywordflow}{do} j=js-1,je ; \textcolor{keywordflow}{do} i=is,ie 146 n2\_v(i,j,1) = 0. 147 n2\_v(i,j,nz+1) = 0. 148 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 149 \textcolor{keywordflow}{ endif} 150 151 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 152 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(halo)) \textcolor{keywordflow}{then} 153 \textcolor{keyword}{call }vert\_fill\_ts(h, tv%T, tv%S, dt\_kappa\_smooth, t, s, g, gv, halo+1) 154 \textcolor{keywordflow}{else} 155 \textcolor{keyword}{call }vert\_fill\_ts(h, tv%T, tv%S, dt\_kappa\_smooth, t, s, g, gv, 1) 156 \textcolor{keywordflow}{ endif} 157 \textcolor{keywordflow}{ endif} 158 159 \textcolor{comment}{! Find the maximum and minimum permitted streamfunction.} 160 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tv%p\_surf)) \textcolor{keywordflow}{then} 161 \textcolor{comment}{!$OMP parallel do default(shared)} 162 \textcolor{keywordflow}{do} j=js-1,je+1 ; \textcolor{keywordflow}{do} i=is-1,ie+1 163 pres(i,j,1) = tv%p\_surf(i,j) 164 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 165 \textcolor{keywordflow}{else} 166 \textcolor{comment}{!$OMP parallel do default(shared)} 167 \textcolor{keywordflow}{do} j=js-1,je+1 ; \textcolor{keywordflow}{do} i=is-1,ie+1 168 pres(i,j,1) = 0.0 169 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 170 \textcolor{keywordflow}{ endif} 171 \textcolor{comment}{!$OMP parallel do default(shared)} 172 \textcolor{keywordflow}{do} j=js-1,je+1 173 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=is-1,ie+1 174 pres(i,j,k+1) = pres(i,j,k) + gv%g\_Earth * gv%H\_to\_RZ * h(i,j,k) 175 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 176 \textcolor{keywordflow}{ enddo} 177 178 eosdom\_u(1) = is-1 - (g%IsdB-1) ; eosdom\_u(2) = ie - (g%IsdB-1) 179 180 \textcolor{comment}{!$OMP parallel do default(none) shared(nz,is,ie,js,je,IsdB,use\_EOS,G,GV,US,pres,T,S,tv,h,e, &} 181 \textcolor{comment}{!$OMP h\_neglect,dz\_neglect,Z\_to\_L,L\_to\_Z,H\_to\_Z,h\_neglect2, &} 182 \textcolor{comment}{!$OMP present\_N2\_u,G\_Rho0,N2\_u,slope\_x,EOSdom\_u,local\_open\_u\_BC, &} 183 \textcolor{comment}{!$OMP OBC) &} 184 \textcolor{comment}{!$OMP private(drdiA,drdiB,drdkL,drdkR,pres\_u,T\_u,S\_u, &} 185 \textcolor{comment}{!$OMP drho\_dT\_u,drho\_dS\_u,hg2A,hg2B,hg2L,hg2R,haA, &} 186 \textcolor{comment}{!$OMP haB,haL,haR,dzaL,dzaR,wtA,wtB,wtL,wtR,drdz, &} 187 \textcolor{comment}{!$OMP drdx,mag\_grad2,Slope,slope2\_Ratio,l\_seg)} 188 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} k=nz,2,-1 189 \textcolor{keywordflow}{if} (.not.(use\_eos)) \textcolor{keywordflow}{then} 190 drdia = 0.0 ; drdib = 0.0 191 drdkl = gv%Rlay(k)-gv%Rlay(k-1) ; drdkr = gv%Rlay(k)-gv%Rlay(k-1) 192 \textcolor{keywordflow}{ endif} 193 194 \textcolor{comment}{! Calculate the zonal isopycnal slope.} 195 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 196 \textcolor{keywordflow}{do} i=is-1,ie 197 pres\_u(i) = 0.5*(pres(i,j,k) + pres(i+1,j,k)) 198 t\_u(i) = 0.25*((t(i,j,k) + t(i+1,j,k)) + (t(i,j,k-1) + t(i+1,j,k-1))) 199 s\_u(i) = 0.25*((s(i,j,k) + s(i+1,j,k)) + (s(i,j,k-1) + s(i+1,j,k-1))) 200 \textcolor{keywordflow}{ enddo} 201 \textcolor{keyword}{call }calculate\_density\_derivs(t\_u, s\_u, pres\_u, drho\_dt\_u, drho\_ds\_u, & 202 tv%eqn\_of\_state, eosdom\_u) 203 \textcolor{keywordflow}{ endif} 204 205 \textcolor{keywordflow}{do} i=is-1,ie 206 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 207 \textcolor{comment}{! Estimate the horizontal density gradients along layers.} 208 drdia = drho\_dt\_u(i) * (t(i+1,j,k-1)-t(i,j,k-1)) + & 209 drho\_ds\_u(i) * (s(i+1,j,k-1)-s(i,j,k-1)) 210 drdib = drho\_dt\_u(i) * (t(i+1,j,k)-t(i,j,k)) + & 211 drho\_ds\_u(i) * (s(i+1,j,k)-s(i,j,k)) 212 213 \textcolor{comment}{! Estimate the vertical density gradients times the grid spacing.} 214 drdkl = (drho\_dt\_u(i) * (t(i,j,k)-t(i,j,k-1)) + & 215 drho\_ds\_u(i) * (s(i,j,k)-s(i,j,k-1))) 216 drdkr = (drho\_dt\_u(i) * (t(i+1,j,k)-t(i+1,j,k-1)) + & 217 drho\_ds\_u(i) * (s(i+1,j,k)-s(i+1,j,k-1))) 218 \textcolor{keywordflow}{ endif} 219 220 221 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 222 hg2a = h(i,j,k-1)*h(i+1,j,k-1) + h\_neglect2 223 hg2b = h(i,j,k)*h(i+1,j,k) + h\_neglect2 224 hg2l = h(i,j,k-1)*h(i,j,k) + h\_neglect2 225 hg2r = h(i+1,j,k-1)*h(i+1,j,k) + h\_neglect2 226 haa = 0.5*(h(i,j,k-1) + h(i+1,j,k-1)) 227 hab = 0.5*(h(i,j,k) + h(i+1,j,k)) + h\_neglect 228 hal = 0.5*(h(i,j,k-1) + h(i,j,k)) + h\_neglect 229 har = 0.5*(h(i+1,j,k-1) + h(i+1,j,k)) + h\_neglect 230 \textcolor{keywordflow}{if} (gv%Boussinesq) \textcolor{keywordflow}{then} 231 dzal = hal * h\_to\_z ; dzar = har * h\_to\_z 232 \textcolor{keywordflow}{else} 233 dzal = 0.5*(e(i,j,k-1) - e(i,j,k+1)) + dz\_neglect 234 dzar = 0.5*(e(i+1,j,k-1) - e(i+1,j,k+1)) + dz\_neglect 235 \textcolor{keywordflow}{ endif} 236 \textcolor{comment}{! Use the harmonic mean thicknesses to weight the horizontal gradients.} 237 \textcolor{comment}{! These unnormalized weights have been rearranged to minimize divisions.} 238 wta = hg2a*hab ; wtb = hg2b*haa 239 wtl = hg2l*(har*dzar) ; wtr = hg2r*(hal*dzal) 240 241 drdz = (wtl * drdkl + wtr * drdkr) / (dzal*wtl + dzar*wtr) 242 \textcolor{comment}{! The expression for drdz above is mathematically equivalent to:} 243 \textcolor{comment}{! drdz = ((hg2L/haL) * drdkL/dzaL + (hg2R/haR) * drdkR/dzaR) / &} 244 \textcolor{comment}{! ((hg2L/haL) + (hg2R/haR))} 245 \textcolor{comment}{! This is the gradient of density along geopotentials.} 246 drdx = ((wta * drdia + wtb * drdib) / (wta + wtb) - & 247 drdz * (e(i,j,k)-e(i+1,j,k))) * g%IdxCu(i,j) 248 249 \textcolor{comment}{! This estimate of slope is accurate for small slopes, but bounded} 250 \textcolor{comment}{! to be between -1 and 1.} 251 mag\_grad2 = drdx**2 + (l\_to\_z*drdz)**2 252 \textcolor{keywordflow}{if} (mag\_grad2 > 0.0) \textcolor{keywordflow}{then} 253 slope\_x(i,j,k) = drdx / sqrt(mag\_grad2) 254 \textcolor{keywordflow}{else} \textcolor{comment}{! Just in case mag\_grad2 = 0 ever.} 255 slope\_x(i,j,k) = 0.0 256 \textcolor{keywordflow}{ endif} 257 258 \textcolor{keywordflow}{if} (present\_n2\_u) n2\_u(i,j,k) = g\_rho0 * drdz * g%mask2dCu(i,j) \textcolor{comment}{! Square of buoyancy frequency [T-2 ~> s-2]} 259 260 \textcolor{keywordflow}{else} \textcolor{comment}{! With .not.use\_EOS, the layers are constant density.} 261 slope\_x(i,j,k) = (z\_to\_l*(e(i,j,k)-e(i+1,j,k))) * g%IdxCu(i,j) 262 \textcolor{keywordflow}{ endif} 263 \textcolor{keywordflow}{if} (local\_open\_u\_bc) \textcolor{keywordflow}{then} 264 l\_seg = obc%segnum\_u(i,j) 265 \textcolor{keywordflow}{if} (l\_seg /= obc\_none) \textcolor{keywordflow}{then} 266 \textcolor{keywordflow}{if} (obc%segment(l\_seg)%open) \textcolor{keywordflow}{then} 267 slope\_x(i,j,k) = 0. 268 \textcolor{comment}{! This and/or the masking code below is to make slopes match inside} 269 \textcolor{comment}{! land mask. Might not be necessary except for DEBUG output.} 270 \textcolor{comment}{! if (OBC%segment(OBC%segnum\_u(I,j))%direction == OBC\_DIRECTION\_E) then} 271 \textcolor{comment}{! slope\_x(I+1,j,K) = 0.} 272 \textcolor{comment}{! else} 273 \textcolor{comment}{! slope\_x(I-1,j,K) = 0.} 274 \textcolor{comment}{! endif} 275 \textcolor{keywordflow}{ endif} 276 \textcolor{keywordflow}{ endif} 277 slope\_x(i,j,k) = slope\_x(i,j,k) * max(g%mask2dT(i,j),g%mask2dT(i+1,j)) 278 \textcolor{keywordflow}{ endif} 279 280 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! I} 281 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of j-loop} 282 283 eosdom\_v(1) = is - (g%isd-1) ; eosdom\_v(2) = ie - (g%isd-1) 284 285 \textcolor{comment}{! Calculate the meridional isopycnal slope.} 286 \textcolor{comment}{!$OMP parallel do default(none) shared(nz,is,ie,js,je,IsdB,use\_EOS,G,GV,US,pres,T,S,tv, &} 287 \textcolor{comment}{!$OMP h,h\_neglect,e,dz\_neglect,Z\_to\_L,L\_to\_Z,H\_to\_Z, &} 288 \textcolor{comment}{!$OMP h\_neglect2,present\_N2\_v,G\_Rho0,N2\_v,slope\_y,EOSdom\_v, &} 289 \textcolor{comment}{!$OMP local\_open\_v\_BC,OBC) &} 290 \textcolor{comment}{!$OMP private(drdjA,drdjB,drdkL,drdkR,pres\_v,T\_v,S\_v, &} 291 \textcolor{comment}{!$OMP drho\_dT\_v,drho\_dS\_v,hg2A,hg2B,hg2L,hg2R,haA, &} 292 \textcolor{comment}{!$OMP haB,haL,haR,dzaL,dzaR,wtA,wtB,wtL,wtR,drdz, &} 293 \textcolor{comment}{!$OMP drdy,mag\_grad2,Slope,slope2\_Ratio,l\_seg)} 294 \textcolor{keywordflow}{do} j=js-1,je ; \textcolor{keywordflow}{do} k=nz,2,-1 295 \textcolor{keywordflow}{if} (.not.(use\_eos)) \textcolor{keywordflow}{then} 296 drdja = 0.0 ; drdjb = 0.0 297 drdkl = gv%Rlay(k)-gv%Rlay(k-1) ; drdkr = gv%Rlay(k)-gv%Rlay(k-1) 298 \textcolor{keywordflow}{ endif} 299 300 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 301 \textcolor{keywordflow}{do} i=is,ie 302 pres\_v(i) = 0.5*(pres(i,j,k) + pres(i,j+1,k)) 303 t\_v(i) = 0.25*((t(i,j,k) + t(i,j+1,k)) + (t(i,j,k-1) + t(i,j+1,k-1))) 304 s\_v(i) = 0.25*((s(i,j,k) + s(i,j+1,k)) + (s(i,j,k-1) + s(i,j+1,k-1))) 305 \textcolor{keywordflow}{ enddo} 306 \textcolor{keyword}{call }calculate\_density\_derivs(t\_v, s\_v, pres\_v, drho\_dt\_v, drho\_ds\_v, tv%eqn\_of\_state, & 307 eosdom\_v) 308 \textcolor{keywordflow}{ endif} 309 \textcolor{keywordflow}{do} i=is,ie 310 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 311 \textcolor{comment}{! Estimate the horizontal density gradients along layers.} 312 drdja = drho\_dt\_v(i) * (t(i,j+1,k-1)-t(i,j,k-1)) + & 313 drho\_ds\_v(i) * (s(i,j+1,k-1)-s(i,j,k-1)) 314 drdjb = drho\_dt\_v(i) * (t(i,j+1,k)-t(i,j,k)) + & 315 drho\_ds\_v(i) * (s(i,j+1,k)-s(i,j,k)) 316 317 \textcolor{comment}{! Estimate the vertical density gradients times the grid spacing.} 318 drdkl = (drho\_dt\_v(i) * (t(i,j,k)-t(i,j,k-1)) + & 319 drho\_ds\_v(i) * (s(i,j,k)-s(i,j,k-1))) 320 drdkr = (drho\_dt\_v(i) * (t(i,j+1,k)-t(i,j+1,k-1)) + & 321 drho\_ds\_v(i) * (s(i,j+1,k)-s(i,j+1,k-1))) 322 \textcolor{keywordflow}{ endif} 323 324 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 325 hg2a = h(i,j,k-1)*h(i,j+1,k-1) + h\_neglect2 326 hg2b = h(i,j,k)*h(i,j+1,k) + h\_neglect2 327 hg2l = h(i,j,k-1)*h(i,j,k) + h\_neglect2 328 hg2r = h(i,j+1,k-1)*h(i,j+1,k) + h\_neglect2 329 haa = 0.5*(h(i,j,k-1) + h(i,j+1,k-1)) + h\_neglect 330 hab = 0.5*(h(i,j,k) + h(i,j+1,k)) + h\_neglect 331 hal = 0.5*(h(i,j,k-1) + h(i,j,k)) + h\_neglect 332 har = 0.5*(h(i,j+1,k-1) + h(i,j+1,k)) + h\_neglect 333 \textcolor{keywordflow}{if} (gv%Boussinesq) \textcolor{keywordflow}{then} 334 dzal = hal * h\_to\_z ; dzar = har * h\_to\_z 335 \textcolor{keywordflow}{else} 336 dzal = 0.5*(e(i,j,k-1) - e(i,j,k+1)) + dz\_neglect 337 dzar = 0.5*(e(i,j+1,k-1) - e(i,j+1,k+1)) + dz\_neglect 338 \textcolor{keywordflow}{ endif} 339 \textcolor{comment}{! Use the harmonic mean thicknesses to weight the horizontal gradients.} 340 \textcolor{comment}{! These unnormalized weights have been rearranged to minimize divisions.} 341 wta = hg2a*hab ; wtb = hg2b*haa 342 wtl = hg2l*(har*dzar) ; wtr = hg2r*(hal*dzal) 343 344 drdz = (wtl * drdkl + wtr * drdkr) / (dzal*wtl + dzar*wtr) 345 \textcolor{comment}{! The expression for drdz above is mathematically equivalent to:} 346 \textcolor{comment}{! drdz = ((hg2L/haL) * drdkL/dzaL + (hg2R/haR) * drdkR/dzaR) / &} 347 \textcolor{comment}{! ((hg2L/haL) + (hg2R/haR))} 348 \textcolor{comment}{! This is the gradient of density along geopotentials.} 349 drdy = ((wta * drdja + wtb * drdjb) / (wta + wtb) - & 350 drdz * (e(i,j,k)-e(i,j+1,k))) * g%IdyCv(i,j) 351 352 \textcolor{comment}{! This estimate of slope is accurate for small slopes, but bounded} 353 \textcolor{comment}{! to be between -1 and 1.} 354 mag\_grad2 = drdy**2 + (l\_to\_z*drdz)**2 355 \textcolor{keywordflow}{if} (mag\_grad2 > 0.0) \textcolor{keywordflow}{then} 356 slope\_y(i,j,k) = drdy / sqrt(mag\_grad2) 357 \textcolor{keywordflow}{else} \textcolor{comment}{! Just in case mag\_grad2 = 0 ever.} 358 slope\_y(i,j,k) = 0.0 359 \textcolor{keywordflow}{ endif} 360 361 \textcolor{keywordflow}{if} (present\_n2\_v) n2\_v(i,j,k) = g\_rho0 * drdz * g%mask2dCv(i,j) \textcolor{comment}{! Square of buoyancy frequency [T-2 ~> s-2]} 362 363 \textcolor{keywordflow}{else} \textcolor{comment}{! With .not.use\_EOS, the layers are constant density.} 364 slope\_y(i,j,k) = (z\_to\_l*(e(i,j,k)-e(i,j+1,k))) * g%IdyCv(i,j) 365 \textcolor{keywordflow}{ endif} 366 \textcolor{keywordflow}{if} (local\_open\_v\_bc) \textcolor{keywordflow}{then} 367 l\_seg = obc%segnum\_v(i,j) 368 \textcolor{keywordflow}{if} (l\_seg /= obc\_none) \textcolor{keywordflow}{then} 369 \textcolor{keywordflow}{if} (obc%segment(l\_seg)%open) \textcolor{keywordflow}{then} 370 slope\_y(i,j,k) = 0. 371 \textcolor{comment}{! This and/or the masking code below is to make slopes match inside} 372 \textcolor{comment}{! land mask. Might not be necessary except for DEBUG output.} 373 \textcolor{comment}{! if (OBC%segment(OBC%segnum\_v(i,J))%direction == OBC\_DIRECTION\_N) then} 374 \textcolor{comment}{! slope\_y(i,J+1,K) = 0.} 375 \textcolor{comment}{! else} 376 \textcolor{comment}{! slope\_y(i,J-1,K) = 0.} 377 \textcolor{comment}{! endif} 378 \textcolor{keywordflow}{ endif} 379 \textcolor{keywordflow}{ endif} 380 slope\_y(i,j,k) = slope\_y(i,j,k) * max(g%mask2dT(i,j),g%mask2dT(i,j+1)) 381 \textcolor{keywordflow}{ endif} 382 383 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! i} 384 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of j-loop} 385 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__isopycnal__slopes_a34691482caaff356da3c5182657dba0d}\label{namespacemom__isopycnal__slopes_a34691482caaff356da3c5182657dba0d}} \index{mom\+\_\+isopycnal\+\_\+slopes@{mom\+\_\+isopycnal\+\_\+slopes}!vert\+\_\+fill\+\_\+ts@{vert\+\_\+fill\+\_\+ts}} \index{vert\+\_\+fill\+\_\+ts@{vert\+\_\+fill\+\_\+ts}!mom\+\_\+isopycnal\+\_\+slopes@{mom\+\_\+isopycnal\+\_\+slopes}} \subsubsection{\texorpdfstring{vert\+\_\+fill\+\_\+ts()}{vert\_fill\_ts()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+isopycnal\+\_\+slopes\+::vert\+\_\+fill\+\_\+ts (\begin{DoxyParamCaption}\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{h, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{T\+\_\+in, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(in)}]{S\+\_\+in, }\item[{real, intent(in)}]{kappa\+\_\+dt, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(out)}]{T\+\_\+f, }\item[{real, dimension( g \%isd\+: g \%ied, g \%jsd\+: g \%jed, g \%ke), intent(out)}]{S\+\_\+f, }\item[{type(ocean\+\_\+grid\+\_\+type), intent(in)}]{G, }\item[{type(verticalgrid\+\_\+type), intent(in)}]{GV, }\item[{integer, intent(in), optional}]{halo\+\_\+here, }\item[{logical, intent(in), optional}]{larger\+\_\+h\+\_\+denom }\end{DoxyParamCaption})} Returns tracer arrays (nominally T and S) with massless layers filled with sensible values, by diffusing vertically with a small but constant diffusivity. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & The ocean\textquotesingle{}s grid structure\\ \hline \mbox{\tt in} & {\em gv} & The ocean\textquotesingle{}s vertical grid structure\\ \hline \mbox{\tt in} & {\em h} & Layer thicknesses \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}\\ \hline \mbox{\tt in} & {\em t\+\_\+in} & Input temperature \mbox{[}degC\mbox{]}\\ \hline \mbox{\tt in} & {\em s\+\_\+in} & Input salinity \mbox{[}ppt\mbox{]}\\ \hline \mbox{\tt in} & {\em kappa\+\_\+dt} & A vertical diffusivity to use for smoothing times a smoothing timescale \mbox{[}Z2 $\sim$$>$ m2\mbox{]}.\\ \hline \mbox{\tt out} & {\em t\+\_\+f} & Filled temperature \mbox{[}degC\mbox{]}\\ \hline \mbox{\tt out} & {\em s\+\_\+f} & Filled salinity \mbox{[}ppt\mbox{]}\\ \hline \mbox{\tt in} & {\em halo\+\_\+here} & Number of halo points to work on, 0 by default\\ \hline \mbox{\tt in} & {\em larger\+\_\+h\+\_\+denom} & Present and true, add a large enough minimal thickness in the denominator of the flux calculations so that the fluxes are never so large as eliminate the transmission of information across groups of massless layers. \\ \hline \end{DoxyParams} Definition at line 391 of file M\+O\+M\+\_\+isopycnal\+\_\+slopes.\+F90. \begin{DoxyCode} 391 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< The ocean's grid structure} 392 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< The ocean's vertical grid structure} 393 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thicknesses [H ~> m or kg m-2]} 394 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: t\_in\textcolor{comment}{ !< Input temperature [degC]} 395 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: s\_in\textcolor{comment}{ !< Input salinity [ppt]} 396 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: kappa\_dt\textcolor{comment}{ !< A vertical diffusivity to use for smoothing} 397 \textcolor{comment}{ !! times a smoothing timescale [Z2 ~> m2].} 398 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: t\_f\textcolor{comment}{ !< Filled temperature [degC]} 399 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: s\_f\textcolor{comment}{ !< Filled salinity [ppt]} 400 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: halo\_here\textcolor{comment}{ !< Number of halo points to work on,} 401 \textcolor{comment}{ !! 0 by default} 402 \textcolor{keywordtype}{logical}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: larger\_h\_denom\textcolor{comment}{ !< Present and true, add a large} 403 \textcolor{comment}{ !! enough minimal thickness in the denominator of} 404 \textcolor{comment}{ !! the flux calculations so that the fluxes are} 405 \textcolor{comment}{ !! never so large as eliminate the transmission} 406 \textcolor{comment}{ !! of information across groups of massless layers.} 407 \textcolor{comment}{! Local variables} 408 \textcolor{keywordtype}{real} :: ent(szi\_(g),szk\_(g)+1) \textcolor{comment}{! The diffusive entrainment (kappa*dt)/dz} 409 \textcolor{comment}{! between layers in a timestep [H ~> m or kg m-2].} 410 \textcolor{keywordtype}{real} :: b1(szi\_(g)), d1(szi\_(g)) \textcolor{comment}{! b1, c1, and d1 are variables used by the} 411 \textcolor{keywordtype}{real} :: c1(szi\_(g),szk\_(g)) \textcolor{comment}{! tridiagonal solver.} 412 \textcolor{keywordtype}{real} :: kap\_dt\_x2 \textcolor{comment}{! The 2*kappa\_dt converted to H units [H2 ~> m2 or kg2 m-4].} 413 \textcolor{keywordtype}{real} :: h\_neglect \textcolor{comment}{! A negligible thickness [H ~> m or kg m-2], to allow for zero thicknesses.} 414 \textcolor{keywordtype}{real} :: h0 \textcolor{comment}{! A negligible thickness to allow for zero thickness layers without} 415 \textcolor{comment}{! completely decouping groups of layers [H ~> m or kg m-2].} 416 \textcolor{comment}{! Often 0 < h\_neglect << h0.} 417 \textcolor{keywordtype}{real} :: h\_tr \textcolor{comment}{! h\_tr is h at tracer points with a tiny thickness} 418 \textcolor{comment}{! added to ensure positive definiteness [H ~> m or kg m-2].} 419 \textcolor{keywordtype}{integer} :: i, j, k, is, ie, js, je, nz, halo 420 421 halo=0 ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(halo\_here)) halo = max(halo\_here,0) 422 423 is = g%isc-halo ; ie = g%iec+halo ; js = g%jsc-halo ; je = g%jec+halo ; nz = gv%ke 424 425 h\_neglect = gv%H\_subroundoff 426 kap\_dt\_x2 = (2.0*kappa\_dt)*gv%Z\_to\_H**2 427 h0 = h\_neglect 428 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(larger\_h\_denom)) \textcolor{keywordflow}{then} 429 \textcolor{keywordflow}{if} (larger\_h\_denom) h0 = 1.0e-16*sqrt(kappa\_dt)*gv%Z\_to\_H 430 \textcolor{keywordflow}{ endif} 431 432 \textcolor{keywordflow}{if} (kap\_dt\_x2 <= 0.0) \textcolor{keywordflow}{then} 433 \textcolor{comment}{!$OMP parallel do default(shared)} 434 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=is,ie 435 t\_f(i,j,k) = t\_in(i,j,k) ; s\_f(i,j,k) = s\_in(i,j,k) 436 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 437 \textcolor{keywordflow}{else} 438 \textcolor{comment}{!$OMP parallel do default(shared) private(ent,b1,d1,c1,h\_tr)} 439 \textcolor{keywordflow}{do} j=js,je 440 \textcolor{keywordflow}{do} i=is,ie 441 ent(i,2) = kap\_dt\_x2 / ((h(i,j,1)+h(i,j,2)) + h0) 442 h\_tr = h(i,j,1) + h\_neglect 443 b1(i) = 1.0 / (h\_tr + ent(i,2)) 444 d1(i) = b1(i) * h\_tr 445 t\_f(i,j,1) = (b1(i)*h\_tr)*t\_in(i,j,1) 446 s\_f(i,j,1) = (b1(i)*h\_tr)*s\_in(i,j,1) 447 \textcolor{keywordflow}{ enddo} 448 \textcolor{keywordflow}{do} k=2,nz-1 ; \textcolor{keywordflow}{do} i=is,ie 449 ent(i,k+1) = kap\_dt\_x2 / ((h(i,j,k)+h(i,j,k+1)) + h0) 450 h\_tr = h(i,j,k) + h\_neglect 451 c1(i,k) = ent(i,k) * b1(i) 452 b1(i) = 1.0 / ((h\_tr + d1(i)*ent(i,k)) + ent(i,k+1)) 453 d1(i) = b1(i) * (h\_tr + d1(i)*ent(i,k)) 454 t\_f(i,j,k) = b1(i) * (h\_tr*t\_in(i,j,k) + ent(i,k)*t\_f(i,j,k-1)) 455 s\_f(i,j,k) = b1(i) * (h\_tr*s\_in(i,j,k) + ent(i,k)*s\_f(i,j,k-1)) 456 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 457 \textcolor{keywordflow}{do} i=is,ie 458 c1(i,nz) = ent(i,nz) * b1(i) 459 h\_tr = h(i,j,nz) + h\_neglect 460 b1(i) = 1.0 / (h\_tr + d1(i)*ent(i,nz)) 461 t\_f(i,j,nz) = b1(i) * (h\_tr*t\_in(i,j,nz) + ent(i,nz)*t\_f(i,j,nz-1)) 462 s\_f(i,j,nz) = b1(i) * (h\_tr*s\_in(i,j,nz) + ent(i,nz)*s\_f(i,j,nz-1)) 463 \textcolor{keywordflow}{ enddo} 464 \textcolor{keywordflow}{do} k=nz-1,1,-1 ; \textcolor{keywordflow}{do} i=is,ie 465 t\_f(i,j,k) = t\_f(i,j,k) + c1(i,k+1)*t\_f(i,j,k+1) 466 s\_f(i,j,k) = s\_f(i,j,k) + c1(i,k+1)*s\_f(i,j,k+1) 467 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 468 \textcolor{keywordflow}{ enddo} 469 \textcolor{keywordflow}{ endif} 470 \end{DoxyCode}