\hypertarget{namespacemom__pressureforce__fv}{}\section{mom\+\_\+pressureforce\+\_\+fv Module Reference} \label{namespacemom__pressureforce__fv}\index{mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}} \subsection{Detailed Description} Finite volume pressure gradient (integrated by quadrature or analytically) Provides the Boussinesq and non-\/\+Boussinesq forms of horizontal accelerations due to pressure gradients using a vertically integrated finite volume form, as described by Adcroft et al., 2008. Integration in the vertical is made either by quadrature or analytically. This form eliminates the thermobaric instabilities that had been a problem with previous forms of the pressure gradient force calculation, as described by Hallberg, 2005. Adcroft, A., R. Hallberg, and M. Harrison, 2008\+: A finite volume discretization of the pressure gradient force using analytic integration. Ocean Modelling, 22, 106-\/113. \href{http://doi.org/10.1016/j.ocemod.2008.02.001}{\tt http\+://doi.\+org/10.\+1016/j.\+ocemod.\+2008.\+02.\+001} Hallberg, 2005\+: A thermobaric instability of Lagrangian vertical coordinate ocean models. Ocean Modelling, 8, 279-\/300. \href{http://dx.doi.org/10.1016/j.ocemod.2004.01.001}{\tt http\+://dx.\+doi.\+org/10.\+1016/j.\+ocemod.\+2004.\+01.\+001} \subsection*{Data Types} \begin{DoxyCompactItemize} \item type \mbox{\hyperlink{structmom__pressureforce__fv_1_1pressureforce__fv__cs}{pressureforce\+\_\+fv\+\_\+cs}} \begin{DoxyCompactList}\small\item\em Finite volume pressure gradient control structure. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection*{Functions/\+Subroutines} \begin{DoxyCompactItemize} \item subroutine, public \mbox{\hyperlink{namespacemom__pressureforce__fv_a50c4a61827e473e643f3f330adf62872}{pressureforce\+\_\+fv\+\_\+nonbouss}} (h, tv, P\+Fu, P\+Fv, G, GV, US, CS, A\+L\+E\+\_\+\+C\+Sp, p\+\_\+atm, pbce, eta) \begin{DoxyCompactList}\small\item\em Non-\/\+Boussinesq analytically-\/integrated finite volume form of pressure gradient. \end{DoxyCompactList}\item subroutine, public \mbox{\hyperlink{namespacemom__pressureforce__fv_a08e9a212bae769c475093a6e0fd6b12b}{pressureforce\+\_\+fv\+\_\+bouss}} (h, tv, P\+Fu, P\+Fv, G, GV, US, CS, A\+L\+E\+\_\+\+C\+Sp, p\+\_\+atm, pbce, eta) \begin{DoxyCompactList}\small\item\em Boussinesq analytically-\/integrated finite volume form of pressure gradient. \end{DoxyCompactList}\item subroutine, public \mbox{\hyperlink{namespacemom__pressureforce__fv_a59da32e314b0a8a8feb6e8cd63ef88e0}{pressureforce\+\_\+fv\+\_\+init}} (Time, G, GV, US, param\+\_\+file, diag, CS, tides\+\_\+\+C\+Sp) \begin{DoxyCompactList}\small\item\em Initializes the finite volume pressure gradient control structure. \end{DoxyCompactList}\item subroutine, public \mbox{\hyperlink{namespacemom__pressureforce__fv_ae4ced0f46de081e5cbd55134d878e802}{pressureforce\+\_\+fv\+\_\+end}} (CS) \begin{DoxyCompactList}\small\item\em Deallocates the finite volume pressure gradient control structure. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection{Function/\+Subroutine Documentation} \mbox{\Hypertarget{namespacemom__pressureforce__fv_a08e9a212bae769c475093a6e0fd6b12b}\label{namespacemom__pressureforce__fv_a08e9a212bae769c475093a6e0fd6b12b}} \index{mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}!pressureforce\+\_\+fv\+\_\+bouss@{pressureforce\+\_\+fv\+\_\+bouss}} \index{pressureforce\+\_\+fv\+\_\+bouss@{pressureforce\+\_\+fv\+\_\+bouss}!mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}} \subsubsection{\texorpdfstring{pressureforce\+\_\+fv\+\_\+bouss()}{pressureforce\_fv\_bouss()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+fv\+::pressureforce\+\_\+fv\+\_\+bouss (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{h, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{real, dimension(szib\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(out)}]{P\+Fu, }\item[{real, dimension(szi\+\_\+(g),szjb\+\_\+(g),szk\+\_\+(g)), intent(out)}]{P\+Fv, }\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__pressureforce__fv_1_1pressureforce__fv__cs}{pressureforce\+\_\+fv\+\_\+cs}}), pointer}]{CS, }\item[{type(ale\+\_\+cs), pointer}]{A\+L\+E\+\_\+\+C\+Sp, }\item[{real, dimension(\+:,\+:), optional, pointer}]{p\+\_\+atm, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(out), optional}]{pbce, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g)), intent(out), optional}]{eta }\end{DoxyParamCaption})} Boussinesq analytically-\/integrated finite volume form of pressure gradient. Determines the acceleration due to hydrostatic pressure forces, using the finite volume form of the terms and analytic integrals in depth. To work, the following fields must be set outside of the usual (is\+:ie,js\+:je) range before this subroutine is called\+: h(isB\+:ie+1,jsB\+:je+1), T(isB\+:ie+1,jsB\+:je+1), and S(isB\+:ie+1,jsB\+:je+1). \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & Ocean grid structure\\ \hline \mbox{\tt in} & {\em gv} & Vertical grid structure\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em h} & Layer thickness \mbox{[}H $\sim$$>$ m\mbox{]}\\ \hline \mbox{\tt in} & {\em tv} & Thermodynamic variables\\ \hline \mbox{\tt out} & {\em pfu} & Zonal acceleration \mbox{[}L T-\/2 $\sim$$>$ m s-\/2\mbox{]}\\ \hline \mbox{\tt out} & {\em pfv} & Meridional acceleration \mbox{[}L T-\/2 $\sim$$>$ m s-\/2\mbox{]}\\ \hline & {\em cs} & Finite volume P\+GF control structure\\ \hline & {\em ale\+\_\+csp} & A\+LE control structure\\ \hline & {\em p\+\_\+atm} & The pressure at the ice-\/ocean or atmosphere-\/ocean interface \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}.\\ \hline \mbox{\tt out} & {\em pbce} & The baroclinic pressure anomaly in each layer due to eta anomalies \mbox{[}L2 T-\/2 H-\/1 $\sim$$>$ m s-\/2 or m4 s-\/2 kg-\/1\mbox{]}.\\ \hline \mbox{\tt out} & {\em eta} & The bottom mass used to calculate P\+Fu and P\+Fv \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}, with any tidal contributions or compressibility compensation. \\ \hline \end{DoxyParams} Definition at line 417 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+F\+V.\+F90. \begin{DoxyCode} 417 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{ !< Ocean grid structure} 418 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< Vertical grid structure} 419 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 420 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thickness [H ~> m]} 421 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< Thermodynamic variables} 422 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: PFu\textcolor{comment}{ !< Zonal acceleration [L T-2 ~> m s-2]} 423 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: PFv\textcolor{comment}{ !< Meridional acceleration [L T-2 ~> m s-2]} 424 \textcolor{keywordtype}{type}(PressureForce\_FV\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< Finite volume PGF control structure} 425 \textcolor{keywordtype}{type}(ALE\_CS), \textcolor{keywordtype}{pointer} :: ALE\_CSp\textcolor{comment}{ !< ALE control structure} 426 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: p\_atm\textcolor{comment}{ !< The pressure at the ice-ocean} 427 \textcolor{comment}{ !! or atmosphere-ocean interface [R L2 T-2 ~> Pa].} 428 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: pbce\textcolor{comment}{ !< The baroclinic pressure} 429 \textcolor{comment}{ !! anomaly in each layer due to eta anomalies} 430 \textcolor{comment}{ !! [L2 T-2 H-1 ~> m s-2 or m4 s-2 kg-1].} 431 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: eta\textcolor{comment}{ !< The bottom mass used to} 432 \textcolor{comment}{ !! calculate PFu and PFv [H ~> m or kg m-2], with any} 433 \textcolor{comment}{ !! tidal contributions or compressibility compensation.} 434 \textcolor{comment}{! Local variables} 435 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)} :: e \textcolor{comment}{! Interface height in depth units [Z ~> m].} 436 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))} :: & 437 e\_tidal, & \textcolor{comment}{! The bottom geopotential anomaly due to tidal forces from} 438 \textcolor{comment}{! astronomical sources and self-attraction and loading [Z ~> m].} 439 dm \textcolor{comment}{! The barotropic adjustment to the Montgomery potential to} 440 \textcolor{comment}{! account for a reduced gravity model [L2 T-2 ~> m2 s-2].} 441 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: & 442 Rho\_cv\_BL \textcolor{comment}{! The coordinate potential density in the deepest variable} 443 \textcolor{comment}{! density near-surface layer [R ~> kg m-3].} 444 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))} :: & 445 dz\_geo, & \textcolor{comment}{! The change in geopotential thickness through a layer [L2 T-2 ~> m2 s-2].} 446 pa, & \textcolor{comment}{! The pressure anomaly (i.e. pressure + g*RHO\_0*e) at the} 447 \textcolor{comment}{! the interface atop a layer [R L2 T-2 ~> Pa].} 448 dpa, & \textcolor{comment}{! The change in pressure anomaly between the top and bottom} 449 \textcolor{comment}{! of a layer [R L2 T-2 ~> Pa].} 450 intz\_dpa \textcolor{comment}{! The vertical integral in depth of the pressure anomaly less the} 451 \textcolor{comment}{! pressure anomaly at the top of the layer [H R L2 T-2 ~> m Pa or kg m-2 Pa].} 452 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G))} :: & 453 intx\_pa, & \textcolor{comment}{! The zonal integral of the pressure anomaly along the interface} 454 \textcolor{comment}{! atop a layer, divided by the grid spacing [R L2 T-2 ~> Pa].} 455 intx\_dpa \textcolor{comment}{! The change in intx\_pa through a layer [R L2 T-2 ~> Pa].} 456 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G))} :: & 457 inty\_pa, & \textcolor{comment}{! The meridional integral of the pressure anomaly along the} 458 \textcolor{comment}{! interface atop a layer, divided by the grid spacing [R L2 T-2 ~> Pa].} 459 inty\_dpa \textcolor{comment}{! The change in inty\_pa through a layer [R L2 T-2 ~> Pa].} 460 461 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{target} :: & 462 T\_tmp, & \textcolor{comment}{! Temporary array of temperatures where layers that are lighter} 463 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [degC].} 464 s\_tmp \textcolor{comment}{! Temporary array of salinities where layers that are lighter} 465 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [ppt].} 466 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))} :: & 467 S\_t, S\_b, T\_t, T\_b \textcolor{comment}{! Top and bottom edge values for linear reconstructions} 468 \textcolor{comment}{! of salinity and temperature within each layer.} 469 \textcolor{keywordtype}{real} :: rho\_in\_situ(SZI\_(G)) \textcolor{comment}{! The in situ density [R ~> kg m-3].} 470 \textcolor{keywordtype}{real} :: p\_ref(SZI\_(G)) \textcolor{comment}{! The pressure used to calculate the coordinate} 471 \textcolor{comment}{! density, [R L2 T-2 ~> Pa] (usually 2e7 Pa = 2000 dbar).} 472 \textcolor{keywordtype}{real} :: p0(SZI\_(G)) \textcolor{comment}{! An array of zeros to use for pressure [R L2 T-2 ~> Pa].} 473 \textcolor{keywordtype}{real} :: h\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 474 \textcolor{comment}{! in roundoff and can be neglected [H ~> m].} 475 \textcolor{keywordtype}{real} :: I\_Rho0 \textcolor{comment}{! The inverse of the Boussinesq reference density [R-1 ~> m3 kg-1].} 476 \textcolor{keywordtype}{real} :: G\_Rho0 \textcolor{comment}{! G\_Earth / Rho0 in [L2 Z-1 T-2 R-1 ~> m4 s-2 kg-1].} 477 \textcolor{keywordtype}{real} :: rho\_ref \textcolor{comment}{! The reference density [R ~> kg m-3].} 478 \textcolor{keywordtype}{real} :: dz\_neglect \textcolor{comment}{! A minimal thickness [Z ~> m], like e.} 479 \textcolor{keywordtype}{logical} :: use\_p\_atm \textcolor{comment}{! If true, use the atmospheric pressure.} 480 \textcolor{keywordtype}{logical} :: use\_ALE \textcolor{comment}{! If true, use an ALE pressure reconstruction.} 481 \textcolor{keywordtype}{logical} :: use\_EOS \textcolor{comment}{! If true, density is calculated from T & S using an equation of state.} 482 \textcolor{keywordtype}{type}(thermo\_var\_ptrs) :: tv\_tmp\textcolor{comment}{! A structure of temporary T & S.} 483 \textcolor{keywordtype}{real} :: Tl(5) \textcolor{comment}{! copy and T in local stencil [degC]} 484 \textcolor{keywordtype}{real} :: mn\_T \textcolor{comment}{! mean of T in local stencil [degC]} 485 \textcolor{keywordtype}{real} :: mn\_T2 \textcolor{comment}{! mean of T**2 in local stencil [degC]} 486 \textcolor{keywordtype}{real} :: hl(5) \textcolor{comment}{! Copy of local stencil of H [H ~> m]} 487 \textcolor{keywordtype}{real} :: r\_sm\_H \textcolor{comment}{! Reciprocal of sum of H in local stencil [H-1 ~> m-1]} 488 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_6 = 1.0/6.0 489 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: EOSdom \textcolor{comment}{! The i-computational domain for the equation of state} 490 \textcolor{keywordtype}{integer} :: is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz, nkmb 491 \textcolor{keywordtype}{integer} :: i, j, k 492 493 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke 494 nkmb=gv%nk\_rho\_varies 495 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB 496 eosdom(1) = isq - (g%isd-1) ; eosdom(2) = g%iec+1 - (g%isd-1) 497 498 \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal, & 499 \textcolor{stringliteral}{"MOM\_PressureForce\_FV\_Bouss: Module must be initialized before it is used."}) 500 501 use\_p\_atm = .false. 502 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(p\_atm)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(p\_atm)) use\_p\_atm = .true. ;\textcolor{keywordflow}{ endif} 503 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 504 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p0(i) = 0.0 ;\textcolor{keywordflow}{ enddo} 505 use\_ale = .false. 506 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(ale\_csp)) use\_ale = cs%reconstruct .and. use\_eos 507 508 \textcolor{keywordflow}{if} (cs%Stanley\_T2\_det\_coeff>=0.) \textcolor{keywordflow}{then} 509 \textcolor{keywordflow}{if} (.not. \textcolor{keyword}{associated}(tv%varT)) \textcolor{keyword}{call }safe\_alloc\_ptr(tv%varT, g%isd, g%ied, g%jsd, g%jed, gv%ke) 510 \textcolor{keywordflow}{do} k=1, nz ; \textcolor{keywordflow}{do} j=js-1,je+1 ; \textcolor{keywordflow}{do} i=is-1,ie+1 511 \textcolor{comment}{! Strictly speaking we should estimate the *horizontal* grid-scale variance} 512 \textcolor{comment}{! but neither of the following blocks make a rotation to the horizontal} 513 \textcolor{comment}{! and instead work along coordinate.} 514 515 \textcolor{comment}{! This block calculates a simple |delta T| along coordinates and does} 516 \textcolor{comment}{! not allow vanishing layer thicknesses or layers tracking topography} 517 \textcolor{comment}{!! SGS variance in i-direction [degC2]} 518 \textcolor{comment}{!dTdi2 = ( ( G%mask2dCu(I ,j) * G%IdxCu(I ,j) * ( tv%T(i+1,j,k) - tv%T(i,j,k) ) &} 519 \textcolor{comment}{! + G%mask2dCu(I-1,j) * G%IdxCu(I-1,j) * ( tv%T(i,j,k) - tv%T(i-1,j,k) ) &} 520 \textcolor{comment}{! ) * G%dxT(i,j) * 0.5 )**2} 521 \textcolor{comment}{!! SGS variance in j-direction [degC2]} 522 \textcolor{comment}{!dTdj2 = ( ( G%mask2dCv(i,J ) * G%IdyCv(i,J ) * ( tv%T(i,j+1,k) - tv%T(i,j,k) ) &} 523 \textcolor{comment}{! + G%mask2dCv(i,J-1) * G%IdyCv(i,J-1) * ( tv%T(i,j,k) - tv%T(i,j-1,k) ) &} 524 \textcolor{comment}{! ) * G%dyT(i,j) * 0.5 )**2} 525 \textcolor{comment}{!tv%varT(i,j,k) = CS%Stanley\_T2\_det\_coeff * 0.5 * ( dTdi2 + dTdj2 )} 526 527 \textcolor{comment}{! This block does a thickness weighted variance calculation and helps control for} 528 \textcolor{comment}{! extreme gradients along layers which are vanished against topography. It is} 529 \textcolor{comment}{! still a poor approximation in the interior when coordinates are strongly tilted.} 530 hl(1) = h(i,j,k) * g%mask2dT(i,j) 531 hl(2) = h(i-1,j,k) * g%mask2dCu(i-1,j) 532 hl(3) = h(i+1,j,k) * g%mask2dCu(i,j) 533 hl(4) = h(i,j-1,k) * g%mask2dCv(i,j-1) 534 hl(5) = h(i,j+1,k) * g%mask2dCv(i,j) 535 r\_sm\_h = 1. / ( ( hl(1) + ( ( hl(2) + hl(3) ) + ( hl(4) + hl(5) ) ) ) + gv%H\_subroundoff ) 536 \textcolor{comment}{! Mean of T} 537 tl(1) = tv%T(i,j,k) ; tl(2) = tv%T(i-1,j,k) ; tl(3) = tv%T(i+1,j,k) 538 tl(4) = tv%T(i,j-1,k) ; tl(5) = tv%T(i,j+1,k) 539 mn\_t = ( hl(1)*tl(1) + ( ( hl(2)*tl(2) + hl(3)*tl(3) ) + ( hl(4)*tl(4) + hl(5)*tl(5) ) ) ) * r\_sm\_h 540 \textcolor{comment}{! Adjust T vectors to have zero mean} 541 tl(:) = tl(:) - mn\_t ; mn\_t = 0. 542 \textcolor{comment}{! Variance of T} 543 mn\_t2 = ( hl(1)*tl(1)*tl(1) + ( ( hl(2)*tl(2)*tl(2) + hl(3)*tl(3)*tl(3) ) & 544 + ( hl(4)*tl(4)*tl(4) + hl(5)*tl(5)*tl(5) ) ) ) * r\_sm\_h 545 \textcolor{comment}{! Variance should be positive but round-off can violate this. Calculating} 546 \textcolor{comment}{! variance directly would fix this but requires more operations.} 547 tv%varT(i,j,k) = cs%Stanley\_T2\_det\_coeff * max(0., mn\_t2) 548 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 549 \textcolor{keywordflow}{ endif} 550 551 h\_neglect = gv%H\_subroundoff 552 dz\_neglect = gv%H\_subroundoff * gv%H\_to\_Z 553 i\_rho0 = 1.0 / gv%Rho0 554 g\_rho0 = gv%g\_Earth / gv%Rho0 555 rho\_ref = cs%Rho0 556 557 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 558 \textcolor{comment}{! Determine the surface height anomaly for calculating self attraction} 559 \textcolor{comment}{! and loading. This should really be based on bottom pressure anomalies,} 560 \textcolor{comment}{! but that is not yet implemented, and the current form is correct for} 561 \textcolor{comment}{! barotropic tides.} 562 \textcolor{comment}{!$OMP parallel do default(shared)} 563 \textcolor{keywordflow}{do} j=jsq,jeq+1 564 \textcolor{keywordflow}{do} i=isq,ieq+1 565 e(i,j,1) = -g%bathyT(i,j) 566 \textcolor{keywordflow}{ enddo} 567 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 568 e(i,j,1) = e(i,j,1) + h(i,j,k)*gv%H\_to\_Z 569 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 570 \textcolor{keywordflow}{ enddo} 571 \textcolor{keyword}{call }calc\_tidal\_forcing(cs%Time, e(:,:,1), e\_tidal, g, cs%tides\_CSp, m\_to\_z=us%m\_to\_Z) 572 \textcolor{keywordflow}{ endif} 573 574 \textcolor{comment}{! Here layer interface heights, e, are calculated.} 575 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 576 \textcolor{comment}{!$OMP parallel do default(shared)} 577 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 578 e(i,j,nz+1) = -(g%bathyT(i,j) + e\_tidal(i,j)) 579 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 580 \textcolor{keywordflow}{else} 581 \textcolor{comment}{!$OMP parallel do default(shared)} 582 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 583 e(i,j,nz+1) = -g%bathyT(i,j) 584 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 585 \textcolor{keywordflow}{ endif} 586 \textcolor{comment}{!$OMP parallel do default(shared)} 587 \textcolor{keywordflow}{do} j=jsq,jeq+1; \textcolor{keywordflow}{do} k=nz,1,-1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 588 e(i,j,k) = e(i,j,k+1) + h(i,j,k)*gv%H\_to\_Z 589 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 590 591 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 592 \textcolor{comment}{! With a bulk mixed layer, replace the T & S of any layers that are} 593 \textcolor{comment}{! lighter than the the buffer layer with the properties of the buffer} 594 \textcolor{comment}{! layer. These layers will be massless anyway, and it avoids any} 595 \textcolor{comment}{! formal calculations with hydrostatically unstable profiles.} 596 597 \textcolor{keywordflow}{if} (nkmb>0) \textcolor{keywordflow}{then} 598 tv\_tmp%T => t\_tmp ; tv\_tmp%S => s\_tmp 599 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 600 601 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p\_ref(i) = tv%P\_Ref ;\textcolor{keywordflow}{ enddo} 602 \textcolor{comment}{!$OMP parallel do default(shared) private(Rho\_cv\_BL)} 603 \textcolor{keywordflow}{do} j=jsq,jeq+1 604 \textcolor{keywordflow}{do} k=1,nkmb ; \textcolor{keywordflow}{do} i=isq,ieq+1 605 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 606 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 607 \textcolor{keyword}{call }calculate\_density(tv%T(:,j,nkmb), tv%S(:,j,nkmb), p\_ref, rho\_cv\_bl(:), & 608 tv%eqn\_of\_state, eosdom) 609 610 \textcolor{keywordflow}{do} k=nkmb+1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 611 \textcolor{keywordflow}{if} (gv%Rlay(k) < rho\_cv\_bl(i)) \textcolor{keywordflow}{then} 612 tv\_tmp%T(i,j,k) = tv%T(i,j,nkmb) ; tv\_tmp%S(i,j,k) = tv%S(i,j,nkmb) 613 \textcolor{keywordflow}{else} 614 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 615 \textcolor{keywordflow}{ endif} 616 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 617 \textcolor{keywordflow}{ enddo} 618 \textcolor{keywordflow}{else} 619 tv\_tmp%T => tv%T ; tv\_tmp%S => tv%S 620 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 621 \textcolor{keywordflow}{ endif} 622 \textcolor{keywordflow}{ endif} 623 624 \textcolor{keywordflow}{if} (cs%GFS\_scale < 1.0) \textcolor{keywordflow}{then} 625 \textcolor{comment}{! Adjust the Montgomery potential to make this a reduced gravity model.} 626 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 627 \textcolor{comment}{!$OMP parallel do default(shared)} 628 \textcolor{keywordflow}{do} j=jsq,jeq+1 629 \textcolor{keywordflow}{if} (use\_p\_atm) \textcolor{keywordflow}{then} 630 \textcolor{keyword}{call }calculate\_density(tv\_tmp%T(:,j,1), tv\_tmp%S(:,j,1), p\_atm(:,j), rho\_in\_situ, & 631 tv%eqn\_of\_state, eosdom) 632 \textcolor{keywordflow}{else} 633 \textcolor{keyword}{call }calculate\_density(tv\_tmp%T(:,j,1), tv\_tmp%S(:,j,1), p0, rho\_in\_situ, & 634 tv%eqn\_of\_state, eosdom) 635 \textcolor{keywordflow}{ endif} 636 \textcolor{keywordflow}{do} i=isq,ieq+1 637 dm(i,j) = (cs%GFS\_scale - 1.0) * (g\_rho0 * rho\_in\_situ(i)) * e(i,j,1) 638 \textcolor{keywordflow}{ enddo} 639 \textcolor{keywordflow}{ enddo} 640 \textcolor{keywordflow}{else} 641 \textcolor{comment}{!$OMP parallel do default(shared)} 642 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 643 dm(i,j) = (cs%GFS\_scale - 1.0) * (g\_rho0 * gv%Rlay(1)) * e(i,j,1) 644 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 645 \textcolor{keywordflow}{ endif} 646 \textcolor{keywordflow}{ endif} 647 \textcolor{comment}{! I have checked that rho\_0 drops out and that the 1-layer case is right. RWH.} 648 649 \textcolor{comment}{! If regridding is activated, do a linear reconstruction of salinity} 650 \textcolor{comment}{! and temperature across each layer. The subscripts 't' and 'b' refer} 651 \textcolor{comment}{! to top and bottom values within each layer (these are the only degrees} 652 \textcolor{comment}{! of freedeom needed to know the linear profile).} 653 \textcolor{keywordflow}{if} ( use\_ale ) \textcolor{keywordflow}{then} 654 \textcolor{keywordflow}{if} ( cs%Recon\_Scheme == 1 ) \textcolor{keywordflow}{then} 655 \textcolor{keyword}{call }ts\_plm\_edge\_values(ale\_csp, s\_t, s\_b, t\_t, t\_b, g, gv, tv, h, cs%boundary\_extrap) 656 \textcolor{keywordflow}{elseif} ( cs%Recon\_Scheme == 2 ) \textcolor{keywordflow}{then} 657 \textcolor{keyword}{call }ts\_ppm\_edge\_values(ale\_csp, s\_t, s\_b, t\_t, t\_b, g, gv, tv, h, cs%boundary\_extrap) 658 \textcolor{keywordflow}{ endif} 659 \textcolor{keywordflow}{ endif} 660 661 \textcolor{comment}{! Set the surface boundary conditions on pressure anomaly and its horizontal} 662 \textcolor{comment}{! integrals, assuming that the surface pressure anomaly varies linearly} 663 \textcolor{comment}{! in x and y.} 664 \textcolor{keywordflow}{if} (use\_p\_atm) \textcolor{keywordflow}{then} 665 \textcolor{comment}{!$OMP parallel do default(shared)} 666 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 667 pa(i,j) = (rho\_ref*gv%g\_Earth)*e(i,j,1) + p\_atm(i,j) 668 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 669 \textcolor{keywordflow}{else} 670 \textcolor{comment}{!$OMP parallel do default(shared)} 671 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 672 pa(i,j) = (rho\_ref*gv%g\_Earth)*e(i,j,1) 673 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 674 \textcolor{keywordflow}{ endif} 675 \textcolor{comment}{!$OMP parallel do default(shared)} 676 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 677 intx\_pa(i,j) = 0.5*(pa(i,j) + pa(i+1,j)) 678 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 679 \textcolor{comment}{!$OMP parallel do default(shared)} 680 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 681 inty\_pa(i,j) = 0.5*(pa(i,j) + pa(i,j+1)) 682 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 683 684 \textcolor{keywordflow}{do} k=1,nz 685 \textcolor{comment}{! Calculate 4 integrals through the layer that are required in the} 686 \textcolor{comment}{! subsequent calculation.} 687 688 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 689 \textcolor{comment}{! The following routine computes the integrals that are needed to} 690 \textcolor{comment}{! calculate the pressure gradient force. Linear profiles for T and S are} 691 \textcolor{comment}{! assumed when regridding is activated. Otherwise, the previous version} 692 \textcolor{comment}{! is used, whereby densities within each layer are constant no matter} 693 \textcolor{comment}{! where the layers are located.} 694 \textcolor{keywordflow}{if} ( use\_ale ) \textcolor{keywordflow}{then} 695 \textcolor{keywordflow}{if} ( cs%Recon\_Scheme == 1 ) \textcolor{keywordflow}{then} 696 \textcolor{keyword}{call }int\_density\_dz\_generic\_plm(k, tv, t\_t, t\_b, s\_t, s\_b, e, & 697 rho\_ref, cs%Rho0, gv%g\_Earth, dz\_neglect, g%bathyT, & 698 g%HI, gv, tv%eqn\_of\_state, us, dpa, intz\_dpa, intx\_dpa, inty\_dpa, & 699 usemasswghtinterp=cs%useMassWghtInterp) 700 \textcolor{keywordflow}{elseif} ( cs%Recon\_Scheme == 2 ) \textcolor{keywordflow}{then} 701 \textcolor{keyword}{call }int\_density\_dz\_generic\_ppm(k, tv, t\_t, t\_b, s\_t, s\_b, e, & 702 rho\_ref, cs%Rho0, gv%g\_Earth, dz\_neglect, g%bathyT, & 703 g%HI, gv, tv%eqn\_of\_state, us, dpa, intz\_dpa, intx\_dpa, inty\_dpa, & 704 usemasswghtinterp=cs%useMassWghtInterp) 705 \textcolor{keywordflow}{ endif} 706 \textcolor{keywordflow}{else} 707 \textcolor{keyword}{call }int\_density\_dz(tv\_tmp%T(:,:,k), tv\_tmp%S(:,:,k), e(:,:,k), e(:,:,k+1), & 708 rho\_ref, cs%Rho0, gv%g\_Earth, g%HI, tv%eqn\_of\_state, us, dpa, & 709 intz\_dpa, intx\_dpa, inty\_dpa, g%bathyT, dz\_neglect, cs%useMassWghtInterp) 710 \textcolor{keywordflow}{ endif} 711 \textcolor{comment}{!$OMP parallel do default(shared)} 712 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 713 intz\_dpa(i,j) = intz\_dpa(i,j)*gv%Z\_to\_H 714 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 715 \textcolor{keywordflow}{else} 716 \textcolor{comment}{!$OMP parallel do default(shared)} 717 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 718 dz\_geo(i,j) = gv%g\_Earth * gv%H\_to\_Z*h(i,j,k) 719 dpa(i,j) = (gv%Rlay(k) - rho\_ref) * dz\_geo(i,j) 720 intz\_dpa(i,j) = 0.5*(gv%Rlay(k) - rho\_ref) * dz\_geo(i,j)*h(i,j,k) 721 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 722 \textcolor{comment}{!$OMP parallel do default(shared)} 723 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 724 intx\_dpa(i,j) = 0.5*(gv%Rlay(k) - rho\_ref) * (dz\_geo(i,j) + dz\_geo(i+1,j)) 725 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 726 \textcolor{comment}{!$OMP parallel do default(shared)} 727 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 728 inty\_dpa(i,j) = 0.5*(gv%Rlay(k) - rho\_ref) * (dz\_geo(i,j) + dz\_geo(i,j+1)) 729 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 730 \textcolor{keywordflow}{ endif} 731 732 \textcolor{comment}{! Compute pressure gradient in x direction} 733 \textcolor{comment}{!$OMP parallel do default(shared)} 734 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 735 pfu(i,j,k) = (((pa(i,j)*h(i,j,k) + intz\_dpa(i,j)) - & 736 (pa(i+1,j)*h(i+1,j,k) + intz\_dpa(i+1,j))) + & 737 ((h(i+1,j,k) - h(i,j,k)) * intx\_pa(i,j) - & 738 (e(i+1,j,k+1) - e(i,j,k+1)) * intx\_dpa(i,j) * gv%Z\_to\_H)) * & 739 ((2.0*i\_rho0*g%IdxCu(i,j)) / & 740 ((h(i,j,k) + h(i+1,j,k)) + h\_neglect)) 741 intx\_pa(i,j) = intx\_pa(i,j) + intx\_dpa(i,j) 742 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 743 \textcolor{comment}{! Compute pressure gradient in y direction} 744 \textcolor{comment}{!$OMP parallel do default(shared)} 745 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 746 pfv(i,j,k) = (((pa(i,j)*h(i,j,k) + intz\_dpa(i,j)) - & 747 (pa(i,j+1)*h(i,j+1,k) + intz\_dpa(i,j+1))) + & 748 ((h(i,j+1,k) - h(i,j,k)) * inty\_pa(i,j) - & 749 (e(i,j+1,k+1) - e(i,j,k+1)) * inty\_dpa(i,j) * gv%Z\_to\_H)) * & 750 ((2.0*i\_rho0*g%IdyCv(i,j)) / & 751 ((h(i,j,k) + h(i,j+1,k)) + h\_neglect)) 752 inty\_pa(i,j) = inty\_pa(i,j) + inty\_dpa(i,j) 753 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 754 \textcolor{comment}{!$OMP parallel do default(shared)} 755 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 756 pa(i,j) = pa(i,j) + dpa(i,j) 757 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 758 \textcolor{keywordflow}{ enddo} 759 760 \textcolor{keywordflow}{if} (cs%GFS\_scale < 1.0) \textcolor{keywordflow}{then} 761 \textcolor{keywordflow}{do} k=1,nz 762 \textcolor{comment}{!$OMP parallel do default(shared)} 763 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 764 pfu(i,j,k) = pfu(i,j,k) - (dm(i+1,j) - dm(i,j)) * g%IdxCu(i,j) 765 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 766 \textcolor{comment}{!$OMP parallel do default(shared)} 767 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 768 pfv(i,j,k) = pfv(i,j,k) - (dm(i,j+1) - dm(i,j)) * g%IdyCv(i,j) 769 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 770 \textcolor{keywordflow}{ enddo} 771 \textcolor{keywordflow}{ endif} 772 773 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(pbce)) \textcolor{keywordflow}{then} 774 \textcolor{keyword}{call }set\_pbce\_bouss(e, tv\_tmp, g, gv, us, cs%Rho0, cs%GFS\_scale, pbce) 775 \textcolor{keywordflow}{ endif} 776 777 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(eta)) \textcolor{keywordflow}{then} 778 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 779 \textcolor{comment}{! eta is the sea surface height relative to a time-invariant geoid, for comparison with} 780 \textcolor{comment}{! what is used for eta in btstep. See how e was calculated about 200 lines above.} 781 \textcolor{comment}{!$OMP parallel do default(shared)} 782 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 783 eta(i,j) = e(i,j,1)*gv%Z\_to\_H + e\_tidal(i,j)*gv%Z\_to\_H 784 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 785 \textcolor{keywordflow}{else} 786 \textcolor{comment}{!$OMP parallel do default(shared)} 787 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 788 eta(i,j) = e(i,j,1)*gv%Z\_to\_H 789 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 790 \textcolor{keywordflow}{ endif} 791 \textcolor{keywordflow}{ endif} 792 793 \textcolor{keywordflow}{if} (cs%id\_e\_tidal>0) \textcolor{keyword}{call }post\_data(cs%id\_e\_tidal, e\_tidal, cs%diag) 794 \textcolor{keywordflow}{if} (cs%id\_tvar\_sgs>0) \textcolor{keyword}{call }post\_data(cs%id\_tvar\_sgs, tv%varT, cs%diag) 795 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__fv_ae4ced0f46de081e5cbd55134d878e802}\label{namespacemom__pressureforce__fv_ae4ced0f46de081e5cbd55134d878e802}} \index{mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}!pressureforce\+\_\+fv\+\_\+end@{pressureforce\+\_\+fv\+\_\+end}} \index{pressureforce\+\_\+fv\+\_\+end@{pressureforce\+\_\+fv\+\_\+end}!mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}} \subsubsection{\texorpdfstring{pressureforce\+\_\+fv\+\_\+end()}{pressureforce\_fv\_end()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+fv\+::pressureforce\+\_\+fv\+\_\+end (\begin{DoxyParamCaption}\item[{type(\mbox{\hyperlink{structmom__pressureforce__fv_1_1pressureforce__fv__cs}{pressureforce\+\_\+fv\+\_\+cs}}), pointer}]{CS }\end{DoxyParamCaption})} Deallocates the finite volume pressure gradient control structure. \begin{DoxyParams}{Parameters} {\em cs} & Finite volume pressure control structure that will be deallocated in this subroutine. \\ \hline \end{DoxyParams} Definition at line 881 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+F\+V.\+F90. \begin{DoxyCode} 881 \textcolor{keywordtype}{type}(PressureForce\_FV\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< Finite volume pressure control structure that} 882 \textcolor{comment}{ !! will be deallocated in this subroutine.} 883 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{deallocate}(cs) \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__fv_a59da32e314b0a8a8feb6e8cd63ef88e0}\label{namespacemom__pressureforce__fv_a59da32e314b0a8a8feb6e8cd63ef88e0}} \index{mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}!pressureforce\+\_\+fv\+\_\+init@{pressureforce\+\_\+fv\+\_\+init}} \index{pressureforce\+\_\+fv\+\_\+init@{pressureforce\+\_\+fv\+\_\+init}!mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}} \subsubsection{\texorpdfstring{pressureforce\+\_\+fv\+\_\+init()}{pressureforce\_fv\_init()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+fv\+::pressureforce\+\_\+fv\+\_\+init (\begin{DoxyParamCaption}\item[{type(time\+\_\+type), intent(in), target}]{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__pressureforce__fv_1_1pressureforce__fv__cs}{pressureforce\+\_\+fv\+\_\+cs}}), pointer}]{CS, }\item[{type(tidal\+\_\+forcing\+\_\+cs), optional, pointer}]{tides\+\_\+\+C\+Sp }\end{DoxyParamCaption})} Initializes the finite volume pressure gradient control structure. \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em time} & Current model time\\ \hline \mbox{\tt in} & {\em g} & Ocean grid structure\\ \hline \mbox{\tt in} & {\em gv} & Vertical grid structure\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em param\+\_\+file} & Parameter file handles\\ \hline \mbox{\tt in,out} & {\em diag} & Diagnostics control structure\\ \hline & {\em cs} & Finite volume P\+GF control structure\\ \hline & {\em tides\+\_\+csp} & Tides control structure \\ \hline \end{DoxyParams} Definition at line 800 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+F\+V.\+F90. \begin{DoxyCode} 800 \textcolor{keywordtype}{type}(time\_type), \textcolor{keywordtype}{target}, \textcolor{keywordtype}{intent(in)} :: Time\textcolor{comment}{ !< Current model time} 801 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{ !< Ocean grid structure} 802 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< Vertical grid structure} 803 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 804 \textcolor{keywordtype}{type}(param\_file\_type), \textcolor{keywordtype}{intent(in)} :: param\_file\textcolor{comment}{ !< Parameter file handles} 805 \textcolor{keywordtype}{type}(diag\_ctrl), \textcolor{keywordtype}{target}, \textcolor{keywordtype}{intent(inout)} :: diag\textcolor{comment}{ !< Diagnostics control structure} 806 \textcolor{keywordtype}{type}(PressureForce\_FV\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< Finite volume PGF control structure} 807 \textcolor{keywordtype}{type}(tidal\_forcing\_CS), \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: tides\_CSp\textcolor{comment}{ !< Tides control structure} 808 \textcolor{comment}{! This include declares and sets the variable "version".} 809 \textcolor{preprocessor}{# include "version\_variable.h"} 810 \textcolor{preprocessor}{} \textcolor{keywordtype}{character(len=40)} :: mdl \textcolor{comment}{! This module's name.} 811 \textcolor{keywordtype}{logical} :: use\_ALE 812 813 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keywordflow}{then} 814 \textcolor{keyword}{call }mom\_error(warning, \textcolor{stringliteral}{"PressureForce\_init called with an associated "}// & 815 \textcolor{stringliteral}{"control structure."}) 816 \textcolor{keywordflow}{return} 817 \textcolor{keywordflow}{else} ; \textcolor{keyword}{allocate}(cs) ;\textcolor{keywordflow}{ endif} 818 819 cs%diag => diag ; cs%Time => time 820 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(tides\_csp)) \textcolor{keywordflow}{then} 821 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tides\_csp)) cs%tides\_CSp => tides\_csp 822 \textcolor{keywordflow}{ endif} 823 824 mdl = \textcolor{stringliteral}{"MOM\_PressureForce\_FV"} 825 \textcolor{keyword}{call }log\_version(param\_file, mdl, version, \textcolor{stringliteral}{""}) 826 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"RHO\_0"}, cs%Rho0, & 827 \textcolor{stringliteral}{"The mean ocean density used with BOUSSINESQ true to "}//& 828 \textcolor{stringliteral}{"calculate accelerations and the mass for conservation "}//& 829 \textcolor{stringliteral}{"properties, or with BOUSSINSEQ false to convert some "}//& 830 \textcolor{stringliteral}{"parameters from vertical units of m to kg m-2."}, & 831 units=\textcolor{stringliteral}{"kg m-3"}, default=1035.0, scale=us%kg\_m3\_to\_R) 832 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TIDES"}, cs%tides, & 833 \textcolor{stringliteral}{"If true, apply tidal momentum forcing."}, default=.false.) 834 \textcolor{keyword}{call }get\_param(param\_file, \textcolor{stringliteral}{"MOM"}, \textcolor{stringliteral}{"USE\_REGRIDDING"}, use\_ale, & 835 \textcolor{stringliteral}{"If True, use the ALE algorithm (regridding/remapping). "}//& 836 \textcolor{stringliteral}{"If False, use the layered isopycnal algorithm."}, default=.false. ) 837 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"MASS\_WEIGHT\_IN\_PRESSURE\_GRADIENT"}, cs%useMassWghtInterp, & 838 \textcolor{stringliteral}{"If true, use mass weighting when interpolating T/S for "}//& 839 \textcolor{stringliteral}{"integrals near the bathymetry in FV pressure gradient "}//& 840 \textcolor{stringliteral}{"calculations."}, default=.false.) 841 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"RECONSTRUCT\_FOR\_PRESSURE"}, cs%reconstruct, & 842 \textcolor{stringliteral}{"If True, use vertical reconstruction of T & S within "}//& 843 \textcolor{stringliteral}{"the integrals of the FV pressure gradient calculation. "}//& 844 \textcolor{stringliteral}{"If False, use the constant-by-layer algorithm. "}//& 845 \textcolor{stringliteral}{"The default is set by USE\_REGRIDDING."}, & 846 default=use\_ale ) 847 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"PRESSURE\_RECONSTRUCTION\_SCHEME"}, cs%Recon\_Scheme, & 848 \textcolor{stringliteral}{"Order of vertical reconstruction of T/S to use in the "}//& 849 \textcolor{stringliteral}{"integrals within the FV pressure gradient calculation.\(\backslash\)n"}//& 850 \textcolor{stringliteral}{" 0: PCM or no reconstruction.\(\backslash\)n"}//& 851 \textcolor{stringliteral}{" 1: PLM reconstruction.\(\backslash\)n"}//& 852 \textcolor{stringliteral}{" 2: PPM reconstruction."}, default=1) 853 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"BOUNDARY\_EXTRAPOLATION\_PRESSURE"}, cs%boundary\_extrap, & 854 \textcolor{stringliteral}{"If true, the reconstruction of T & S for pressure in "}//& 855 \textcolor{stringliteral}{"boundary cells is extrapolated, rather than using PCM "}//& 856 \textcolor{stringliteral}{"in these cells. If true, the same order polynomial is "}//& 857 \textcolor{stringliteral}{"used as is used for the interior cells."}, default=.true.) 858 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"PGF\_STANLEY\_T2\_DET\_COEFF"}, cs%Stanley\_T2\_det\_coeff, & 859 \textcolor{stringliteral}{"The coefficient correlating SGS temperature variance with "}// & 860 \textcolor{stringliteral}{"the mean temperature gradient in the deterministic part of "}// & 861 \textcolor{stringliteral}{"the Stanley form of the Brankart correction. "}// & 862 \textcolor{stringliteral}{"Negative values disable the scheme."}, units=\textcolor{stringliteral}{"nondim"}, default=-1.0) 863 \textcolor{keywordflow}{if} (cs%Stanley\_T2\_det\_coeff>=0.) \textcolor{keywordflow}{then} 864 cs%id\_tvar\_sgs = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'tvar\_sgs\_pgf'}, diag%axesTL, & 865 time, \textcolor{stringliteral}{'SGS temperature variance used in PGF'}, \textcolor{stringliteral}{'degC2'}) 866 \textcolor{keywordflow}{ endif} 867 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 868 cs%id\_e\_tidal = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'e\_tidal'}, diag%axesT1, & 869 time, \textcolor{stringliteral}{'Tidal Forcing Astronomical and SAL Height Anomaly'}, \textcolor{stringliteral}{'meter'}, conversion=us%Z\_to\_m) 870 \textcolor{keywordflow}{ endif} 871 872 cs%GFS\_scale = 1.0 873 \textcolor{keywordflow}{if} (gv%g\_prime(1) /= gv%g\_Earth) cs%GFS\_scale = gv%g\_prime(1) / gv%g\_Earth 874 875 \textcolor{keyword}{call }log\_param(param\_file, mdl, \textcolor{stringliteral}{"GFS / G\_EARTH"}, cs%GFS\_scale) 876 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__fv_a50c4a61827e473e643f3f330adf62872}\label{namespacemom__pressureforce__fv_a50c4a61827e473e643f3f330adf62872}} \index{mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}!pressureforce\+\_\+fv\+\_\+nonbouss@{pressureforce\+\_\+fv\+\_\+nonbouss}} \index{pressureforce\+\_\+fv\+\_\+nonbouss@{pressureforce\+\_\+fv\+\_\+nonbouss}!mom\+\_\+pressureforce\+\_\+fv@{mom\+\_\+pressureforce\+\_\+fv}} \subsubsection{\texorpdfstring{pressureforce\+\_\+fv\+\_\+nonbouss()}{pressureforce\_fv\_nonbouss()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+fv\+::pressureforce\+\_\+fv\+\_\+nonbouss (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in)}]{h, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\item[{real, dimension(szib\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(out)}]{P\+Fu, }\item[{real, dimension(szi\+\_\+(g),szjb\+\_\+(g),szk\+\_\+(g)), intent(out)}]{P\+Fv, }\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__pressureforce__fv_1_1pressureforce__fv__cs}{pressureforce\+\_\+fv\+\_\+cs}}), pointer}]{CS, }\item[{type(ale\+\_\+cs), pointer}]{A\+L\+E\+\_\+\+C\+Sp, }\item[{real, dimension(\+:,\+:), optional, pointer}]{p\+\_\+atm, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(out), optional}]{pbce, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g)), intent(out), optional}]{eta }\end{DoxyParamCaption})} Non-\/\+Boussinesq analytically-\/integrated finite volume form of pressure gradient. Determines the acceleration due to hydrostatic pressure forces, using the analytic finite volume form of the Pressure gradient, and does not make the Boussinesq approximation. To work, the following fields must be set outside of the usual (is\+:ie,js\+:je) range before this subroutine is called\+: h(isB\+:ie+1,jsB\+:je+1), T(isB\+:ie+1,jsB\+:je+1), and S(isB\+:ie+1,jsB\+:je+1). \begin{DoxyParams}[1]{Parameters} \mbox{\tt in} & {\em g} & Ocean grid structure\\ \hline \mbox{\tt in} & {\em gv} & Vertical grid structure\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em h} & Layer thickness \mbox{[}H $\sim$$>$ kg/m2\mbox{]}\\ \hline \mbox{\tt in} & {\em tv} & Thermodynamic variables\\ \hline \mbox{\tt out} & {\em pfu} & Zonal acceleration \mbox{[}L T-\/2 $\sim$$>$ m s-\/2\mbox{]}\\ \hline \mbox{\tt out} & {\em pfv} & Meridional acceleration \mbox{[}L T-\/2 $\sim$$>$ m s-\/2\mbox{]}\\ \hline & {\em cs} & Finite volume P\+GF control structure\\ \hline & {\em ale\+\_\+csp} & A\+LE control structure\\ \hline & {\em p\+\_\+atm} & The pressure at the ice-\/ocean or atmosphere-\/ocean interface \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}.\\ \hline \mbox{\tt out} & {\em pbce} & The baroclinic pressure anomaly in each layer due to eta anomalies \mbox{[}L2 T-\/2 H-\/1 $\sim$$>$ m s-\/2 or m4 s-\/2 kg-\/1\mbox{]}.\\ \hline \mbox{\tt out} & {\em eta} & The bottom mass used to calculate P\+Fu and P\+Fv \mbox{[}H $\sim$$>$ m or kg m-\/2\mbox{]}, with any tidal contributions or compressibility compensation. \\ \hline \end{DoxyParams} Definition at line 77 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+F\+V.\+F90. \begin{DoxyCode} 77 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: G\textcolor{comment}{ !< Ocean grid structure} 78 \textcolor{keywordtype}{type}(verticalGrid\_type), \textcolor{keywordtype}{intent(in)} :: GV\textcolor{comment}{ !< Vertical grid structure} 79 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: US\textcolor{comment}{ !< A dimensional unit scaling type} 80 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thickness [H ~> kg/m2]} 81 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< Thermodynamic variables} 82 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: PFu\textcolor{comment}{ !< Zonal acceleration [L T-2 ~> m s-2]} 83 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: PFv\textcolor{comment}{ !< Meridional acceleration [L T-2 ~> m s-2]} 84 \textcolor{keywordtype}{type}(PressureForce\_FV\_CS), \textcolor{keywordtype}{pointer} :: CS\textcolor{comment}{ !< Finite volume PGF control structure} 85 \textcolor{keywordtype}{type}(ALE\_CS), \textcolor{keywordtype}{pointer} :: ALE\_CSp\textcolor{comment}{ !< ALE control structure} 86 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: p\_atm\textcolor{comment}{ !< The pressure at the ice-ocean} 87 \textcolor{comment}{ !! or atmosphere-ocean interface [R L2 T-2 ~> Pa].} 88 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: pbce\textcolor{comment}{ !< The baroclinic pressure} 89 \textcolor{comment}{ !! anomaly in each layer due to eta anomalies} 90 \textcolor{comment}{ !! [L2 T-2 H-1 ~> m s-2 or m4 s-2 kg-1].} 91 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: eta\textcolor{comment}{ !< The bottom mass used to} 92 \textcolor{comment}{ !! calculate PFu and PFv [H ~> m or kg m-2], with any tidal} 93 \textcolor{comment}{ !! contributions or compressibility compensation.} 94 \textcolor{comment}{! Local variables} 95 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)} :: p \textcolor{comment}{! Interface pressure [R L2 T-2 ~> Pa].} 96 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{target} :: & 97 T\_tmp, & \textcolor{comment}{! Temporary array of temperatures where layers that are lighter} 98 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [degC].} 99 s\_tmp \textcolor{comment}{! Temporary array of salinities where layers that are lighter} 100 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [ppt].} 101 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))} :: & 102 S\_t, & \textcolor{comment}{! Top and bottom edge values for linear reconstructions} 103 S\_b, & \textcolor{comment}{! of salinity within each layer [ppt].} 104 T\_t, & \textcolor{comment}{! Top and bottom edge values for linear reconstructions} 105 T\_b \textcolor{comment}{! of temperature within each layer [degC].} 106 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))} :: & 107 dza, & \textcolor{comment}{! The change in geopotential anomaly between the top and bottom} 108 \textcolor{comment}{! of a layer [L2 T-2 ~> m2 s-2].} 109 intp\_dza \textcolor{comment}{! The vertical integral in depth of the pressure anomaly less} 110 \textcolor{comment}{! the pressure anomaly at the top of the layer [R L4 Z-4 ~> Pa m2 s-2].} 111 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))} :: & 112 dp, & \textcolor{comment}{! The (positive) change in pressure across a layer [R L2 Z-2 ~> Pa].} 113 SSH, & \textcolor{comment}{! The sea surface height anomaly, in depth units [Z ~> m].} 114 e\_tidal, & \textcolor{comment}{! The bottom geopotential anomaly due to tidal forces from} 115 \textcolor{comment}{! astronomical sources and self-attraction and loading [Z ~> m].} 116 dm, & \textcolor{comment}{! The barotropic adjustment to the Montgomery potential to} 117 \textcolor{comment}{! account for a reduced gravity model [L2 T-2 ~> m2 s-2].} 118 za \textcolor{comment}{! The geopotential anomaly (i.e. g*e + alpha\_0*pressure) at the} 119 \textcolor{comment}{! interface atop a layer [L2 T-2 ~> m2 s-2].} 120 121 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: Rho\_cv\_BL \textcolor{comment}{! The coordinate potential density in the deepest variable} 122 \textcolor{comment}{! density near-surface layer [R ~> kg m-3].} 123 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G))} :: & 124 intx\_za \textcolor{comment}{! The zonal integral of the geopotential anomaly along the} 125 \textcolor{comment}{! interface below a layer, divided by the grid spacing [L2 T-2 ~> m2 s-2].} 126 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))} :: & 127 intx\_dza \textcolor{comment}{! The change in intx\_za through a layer [L2 T-2 ~> m2 s-2].} 128 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G))} :: & 129 inty\_za \textcolor{comment}{! The meridional integral of the geopotential anomaly along the} 130 \textcolor{comment}{! interface below a layer, divided by the grid spacing [L2 T-2 ~> m2 s-2].} 131 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))} :: & 132 inty\_dza \textcolor{comment}{! The change in inty\_za through a layer [L2 T-2 ~> m2 s-2].} 133 \textcolor{keywordtype}{real} :: p\_ref(SZI\_(G)) \textcolor{comment}{! The pressure used to calculate the coordinate} 134 \textcolor{comment}{! density, [R L2 T-2 ~> Pa] (usually 2e7 Pa = 2000 dbar).} 135 136 \textcolor{keywordtype}{real} :: dp\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 137 \textcolor{comment}{! in roundoff and can be neglected [R L2 T-2 ~> Pa].} 138 \textcolor{keywordtype}{real} :: I\_gEarth \textcolor{comment}{! The inverse of GV%g\_Earth [L2 Z L-2 ~> s2 m-1]} 139 \textcolor{keywordtype}{real} :: alpha\_anom \textcolor{comment}{! The in-situ specific volume, averaged over a} 140 \textcolor{comment}{! layer, less alpha\_ref [R-1 ~> m3 kg-1].} 141 \textcolor{keywordtype}{logical} :: use\_p\_atm \textcolor{comment}{! If true, use the atmospheric pressure.} 142 \textcolor{keywordtype}{logical} :: use\_ALE \textcolor{comment}{! If true, use an ALE pressure reconstruction.} 143 \textcolor{keywordtype}{logical} :: use\_EOS \textcolor{comment}{! If true, density is calculated from T & S using an equation of state.} 144 \textcolor{keywordtype}{type}(thermo\_var\_ptrs) :: tv\_tmp\textcolor{comment}{! A structure of temporary T & S.} 145 146 \textcolor{keywordtype}{real} :: alpha\_ref \textcolor{comment}{! A reference specific volume [R-1 ~> m3 kg-1] that is used} 147 \textcolor{comment}{! to reduce the impact of truncation errors.} 148 \textcolor{keywordtype}{real} :: rho\_in\_situ(SZI\_(G)) \textcolor{comment}{! The in situ density [R ~> kg m-3].} 149 \textcolor{keywordtype}{real} :: Pa\_to\_H \textcolor{comment}{! A factor to convert from Pa to the thicknesss units (H) [H T2 R-1 L-2 ~> H Pa-1].} 150 \textcolor{keywordtype}{real} :: H\_to\_RL2\_T2 \textcolor{comment}{! A factor to convert from thicknesss units (H) to pressure units [R L2 T-2 H-1 ~> Pa H-1].} 151 \textcolor{comment}{! real :: oneatm = 101325.0 ! 1 atm in [Pa] = [kg m-1 s-2]} 152 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{parameter} :: C1\_6 = 1.0/6.0 153 \textcolor{keywordtype}{integer} :: is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz, nkmb 154 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: EOSdom \textcolor{comment}{! The i-computational domain for the equation of state} 155 \textcolor{keywordtype}{integer} :: i, j, k 156 157 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke 158 nkmb=gv%nk\_rho\_varies 159 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB 160 eosdom(1) = isq - (g%isd-1) ; eosdom(2) = g%iec+1 - (g%isd-1) 161 162 \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal, & 163 \textcolor{stringliteral}{"MOM\_PressureForce\_FV\_nonBouss: Module must be initialized before it is used."}) 164 \textcolor{keywordflow}{if} (cs%Stanley\_T2\_det\_coeff>=0.) \textcolor{keyword}{call }mom\_error(fatal, & 165 \textcolor{stringliteral}{"MOM\_PressureForce\_FV\_nonBouss: The Stanley parameterization is not yet"}//& 166 \textcolor{stringliteral}{"implemented in non-Boussinesq mode."}) 167 168 use\_p\_atm = .false. 169 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(p\_atm)) \textcolor{keywordflow}{then} ; \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(p\_atm)) use\_p\_atm = .true. ;\textcolor{keywordflow}{ endif} 170 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 171 use\_ale = .false. 172 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(ale\_csp)) use\_ale = cs%reconstruct .and. use\_eos 173 174 h\_to\_rl2\_t2 = gv%g\_Earth*gv%H\_to\_RZ 175 dp\_neglect = gv%g\_Earth*gv%H\_to\_RZ * gv%H\_subroundoff 176 alpha\_ref = 1.0 / cs%Rho0 177 i\_gearth = 1.0 / gv%g\_Earth 178 179 \textcolor{keywordflow}{if} (use\_p\_atm) \textcolor{keywordflow}{then} 180 \textcolor{comment}{!$OMP parallel do default(shared)} 181 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 182 p(i,j,1) = p\_atm(i,j) 183 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 184 \textcolor{keywordflow}{else} 185 \textcolor{comment}{!$OMP parallel do default(shared)} 186 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 187 p(i,j,1) = 0.0 \textcolor{comment}{! or oneatm} 188 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 189 \textcolor{keywordflow}{ endif} 190 \textcolor{comment}{!$OMP parallel do default(shared)} 191 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} k=2,nz+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 192 p(i,j,k) = p(i,j,k-1) + h\_to\_rl2\_t2 * h(i,j,k-1) 193 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 194 195 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 196 \textcolor{comment}{! With a bulk mixed layer, replace the T & S of any layers that are} 197 \textcolor{comment}{! lighter than the the buffer layer with the properties of the buffer} 198 \textcolor{comment}{! layer. These layers will be massless anyway, and it avoids any} 199 \textcolor{comment}{! formal calculations with hydrostatically unstable profiles.} 200 \textcolor{keywordflow}{if} (nkmb>0) \textcolor{keywordflow}{then} 201 tv\_tmp%T => t\_tmp ; tv\_tmp%S => s\_tmp 202 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 203 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p\_ref(i) = tv%P\_Ref ;\textcolor{keywordflow}{ enddo} 204 \textcolor{comment}{!$OMP parallel do default(shared) private(Rho\_cv\_BL)} 205 \textcolor{keywordflow}{do} j=jsq,jeq+1 206 \textcolor{keywordflow}{do} k=1,nkmb ; \textcolor{keywordflow}{do} i=isq,ieq+1 207 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 208 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 209 \textcolor{keyword}{call }calculate\_density(tv%T(:,j,nkmb), tv%S(:,j,nkmb), p\_ref, rho\_cv\_bl(:), & 210 tv%eqn\_of\_state, eosdom) 211 \textcolor{keywordflow}{do} k=nkmb+1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 212 \textcolor{keywordflow}{if} (gv%Rlay(k) < rho\_cv\_bl(i)) \textcolor{keywordflow}{then} 213 tv\_tmp%T(i,j,k) = tv%T(i,j,nkmb) ; tv\_tmp%S(i,j,k) = tv%S(i,j,nkmb) 214 \textcolor{keywordflow}{else} 215 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 216 \textcolor{keywordflow}{ endif} 217 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 218 \textcolor{keywordflow}{ enddo} 219 \textcolor{keywordflow}{else} 220 tv\_tmp%T => tv%T ; tv\_tmp%S => tv%S 221 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 222 \textcolor{keywordflow}{ endif} 223 \textcolor{keywordflow}{ endif} 224 225 \textcolor{comment}{! If regridding is activated, do a linear reconstruction of salinity} 226 \textcolor{comment}{! and temperature across each layer. The subscripts 't' and 'b' refer} 227 \textcolor{comment}{! to top and bottom values within each layer (these are the only degrees} 228 \textcolor{comment}{! of freedeom needed to know the linear profile).} 229 \textcolor{keywordflow}{if} ( use\_ale ) \textcolor{keywordflow}{then} 230 \textcolor{keywordflow}{if} ( cs%Recon\_Scheme == 1 ) \textcolor{keywordflow}{then} 231 \textcolor{keyword}{call }ts\_plm\_edge\_values(ale\_csp, s\_t, s\_b, t\_t, t\_b, g, gv, tv, h, cs%boundary\_extrap) 232 \textcolor{keywordflow}{elseif} ( cs%Recon\_Scheme == 2) \textcolor{keywordflow}{then} 233 \textcolor{keyword}{call }ts\_ppm\_edge\_values(ale\_csp, s\_t, s\_b, t\_t, t\_b, g, gv, tv, h, cs%boundary\_extrap) 234 \textcolor{keywordflow}{ endif} 235 \textcolor{keywordflow}{ endif} 236 237 \textcolor{comment}{!$OMP parallel do default(shared) private(alpha\_anom,dp)} 238 \textcolor{keywordflow}{do} k=1,nz 239 \textcolor{comment}{! Calculate 4 integrals through the layer that are required in the} 240 \textcolor{comment}{! subsequent calculation.} 241 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 242 \textcolor{keywordflow}{if} ( use\_ale ) \textcolor{keywordflow}{then} 243 \textcolor{keywordflow}{if} ( cs%Recon\_Scheme == 1 ) \textcolor{keywordflow}{then} 244 \textcolor{keyword}{call }int\_spec\_vol\_dp\_generic\_plm( t\_t(:,:,k), t\_b(:,:,k), s\_t(:,:,k), s\_b(:,:,k), & 245 p(:,:,k), p(:,:,k+1), alpha\_ref, dp\_neglect, p(:,:,nz+1), g%HI, & 246 tv%eqn\_of\_state, us, dza(:,:,k), intp\_dza(:,:,k), intx\_dza(:,:,k), inty\_dza(:,:,k), & 247 usemasswghtinterp=cs%useMassWghtInterp) 248 \textcolor{keywordflow}{elseif} ( cs%Recon\_Scheme == 2 ) \textcolor{keywordflow}{then} 249 \textcolor{keyword}{call }mom\_error(fatal, \textcolor{stringliteral}{"PressureForce\_FV\_nonBouss: "}//& 250 \textcolor{stringliteral}{"int\_spec\_vol\_dp\_generic\_ppm does not exist yet."}) 251 \textcolor{comment}{! call int\_spec\_vol\_dp\_generic\_ppm ( tv%T(:,:,k), T\_t(:,:,k), T\_b(:,:,k), &} 252 \textcolor{comment}{! tv%S(:,:,k), S\_t(:,:,k), S\_b(:,:,k), p(:,:,K), p(:,:,K+1), &} 253 \textcolor{comment}{! alpha\_ref, G%HI, tv%eqn\_of\_state, dza(:,:,k), intp\_dza(:,:,k), &} 254 \textcolor{comment}{! intx\_dza(:,:,k), inty\_dza(:,:,k))} 255 \textcolor{keywordflow}{ endif} 256 \textcolor{keywordflow}{else} 257 \textcolor{keyword}{call }int\_specific\_vol\_dp(tv\_tmp%T(:,:,k), tv\_tmp%S(:,:,k), p(:,:,k), & 258 p(:,:,k+1), alpha\_ref, g%HI, tv%eqn\_of\_state, & 259 us, dza(:,:,k), intp\_dza(:,:,k), intx\_dza(:,:,k), & 260 inty\_dza(:,:,k), bathyp=p(:,:,nz+1), dp\_tiny=dp\_neglect, & 261 usemasswghtinterp=cs%useMassWghtInterp) 262 \textcolor{keywordflow}{ endif} 263 \textcolor{keywordflow}{else} 264 alpha\_anom = 1.0 / gv%Rlay(k) - alpha\_ref 265 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 266 dp(i,j) = h\_to\_rl2\_t2 * h(i,j,k) 267 dza(i,j,k) = alpha\_anom * dp(i,j) 268 intp\_dza(i,j,k) = 0.5 * alpha\_anom * dp(i,j)**2 269 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 270 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 271 intx\_dza(i,j,k) = 0.5 * alpha\_anom * (dp(i,j)+dp(i+1,j)) 272 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 273 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 274 inty\_dza(i,j,k) = 0.5 * alpha\_anom * (dp(i,j)+dp(i,j+1)) 275 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 276 \textcolor{keywordflow}{ endif} 277 \textcolor{keywordflow}{ enddo} 278 279 \textcolor{comment}{! The bottom geopotential anomaly is calculated first so that the increments} 280 \textcolor{comment}{! to the geopotential anomalies can be reused. Alternately, the surface} 281 \textcolor{comment}{! geopotential could be calculated directly with separate calls to} 282 \textcolor{comment}{! int\_specific\_vol\_dp with alpha\_ref=0, and the anomalies used going} 283 \textcolor{comment}{! downward, which would relieve the need for dza, intp\_dza, intx\_dza, and} 284 \textcolor{comment}{! inty\_dza to be 3-D arrays.} 285 286 \textcolor{comment}{! Sum vertically to determine the surface geopotential anomaly.} 287 \textcolor{comment}{!$OMP parallel do default(shared)} 288 \textcolor{keywordflow}{do} j=jsq,jeq+1 289 \textcolor{keywordflow}{do} i=isq,ieq+1 290 za(i,j) = alpha\_ref*p(i,j,nz+1) - gv%g\_Earth*g%bathyT(i,j) 291 \textcolor{keywordflow}{ enddo} 292 \textcolor{keywordflow}{do} k=nz,1,-1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 293 za(i,j) = za(i,j) + dza(i,j,k) 294 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 295 \textcolor{keywordflow}{ enddo} 296 297 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 298 \textcolor{comment}{! Find and add the tidal geopotential anomaly.} 299 \textcolor{comment}{!$OMP parallel do default(shared)} 300 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 301 ssh(i,j) = (za(i,j) - alpha\_ref*p(i,j,1)) * i\_gearth 302 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 303 \textcolor{keyword}{call }calc\_tidal\_forcing(cs%Time, ssh, e\_tidal, g, cs%tides\_CSp, m\_to\_z=us%m\_to\_Z) 304 \textcolor{comment}{!$OMP parallel do default(shared)} 305 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 306 za(i,j) = za(i,j) - gv%g\_Earth * e\_tidal(i,j) 307 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 308 \textcolor{keywordflow}{ endif} 309 310 \textcolor{keywordflow}{if} (cs%GFS\_scale < 1.0) \textcolor{keywordflow}{then} 311 \textcolor{comment}{! Adjust the Montgomery potential to make this a reduced gravity model.} 312 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 313 \textcolor{comment}{!$OMP parallel do default(shared) private(rho\_in\_situ)} 314 \textcolor{keywordflow}{do} j=jsq,jeq+1 315 \textcolor{keyword}{call }calculate\_density(tv\_tmp%T(:,j,1), tv\_tmp%S(:,j,1), p(:,j,1), rho\_in\_situ, & 316 tv%eqn\_of\_state, eosdom) 317 318 \textcolor{keywordflow}{do} i=isq,ieq+1 319 dm(i,j) = (cs%GFS\_scale - 1.0) * (p(i,j,1)*(1.0/rho\_in\_situ(i) - alpha\_ref) + za(i,j)) 320 \textcolor{keywordflow}{ enddo} 321 \textcolor{keywordflow}{ enddo} 322 \textcolor{keywordflow}{else} 323 \textcolor{comment}{!$OMP parallel do default(shared)} 324 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 325 dm(i,j) = (cs%GFS\_scale - 1.0) * (p(i,j,1)*(1.0/gv%Rlay(1) - alpha\_ref) + za(i,j)) 326 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 327 \textcolor{keywordflow}{ endif} 328 \textcolor{comment}{! else} 329 \textcolor{comment}{! do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 ; dM(i,j) = 0.0 ; enddo ; enddo} 330 \textcolor{keywordflow}{ endif} 331 332 \textcolor{comment}{! This order of integrating upward and then downward again is necessary with} 333 \textcolor{comment}{! a nonlinear equation of state, so that the surface geopotentials will go} 334 \textcolor{comment}{! linearly between the values at thickness points, but the bottom} 335 \textcolor{comment}{! geopotentials will not now be linear at the sub-grid-scale. Doing this} 336 \textcolor{comment}{! ensures no motion with flat isopycnals, even with a nonlinear equation of state.} 337 \textcolor{comment}{!$OMP parallel do default(shared)} 338 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 339 intx\_za(i,j) = 0.5*(za(i,j) + za(i+1,j)) 340 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 341 \textcolor{comment}{!$OMP parallel do default(shared)} 342 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 343 inty\_za(i,j) = 0.5*(za(i,j) + za(i,j+1)) 344 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 345 \textcolor{keywordflow}{do} k=1,nz 346 \textcolor{comment}{! These expressions for the acceleration have been carefully checked in} 347 \textcolor{comment}{! a set of idealized cases, and should be bug-free.} 348 \textcolor{comment}{!$OMP parallel do default(shared)} 349 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 350 dp(i,j) = h\_to\_rl2\_t2 * h(i,j,k) 351 za(i,j) = za(i,j) - dza(i,j,k) 352 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 353 \textcolor{comment}{!$OMP parallel do default(shared)} 354 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 355 intx\_za(i,j) = intx\_za(i,j) - intx\_dza(i,j,k) 356 pfu(i,j,k) = ( ((za(i,j)*dp(i,j) + intp\_dza(i,j,k)) - & 357 (za(i+1,j)*dp(i+1,j) + intp\_dza(i+1,j,k))) + & 358 ((dp(i+1,j) - dp(i,j)) * intx\_za(i,j) - & 359 (p(i+1,j,k) - p(i,j,k)) * intx\_dza(i,j,k)) ) * & 360 (2.0*g%IdxCu(i,j) / ((dp(i,j) + dp(i+1,j)) + dp\_neglect)) 361 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 362 \textcolor{comment}{!$OMP parallel do default(shared)} 363 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 364 inty\_za(i,j) = inty\_za(i,j) - inty\_dza(i,j,k) 365 pfv(i,j,k) = (((za(i,j)*dp(i,j) + intp\_dza(i,j,k)) - & 366 (za(i,j+1)*dp(i,j+1) + intp\_dza(i,j+1,k))) + & 367 ((dp(i,j+1) - dp(i,j)) * inty\_za(i,j) - & 368 (p(i,j+1,k) - p(i,j,k)) * inty\_dza(i,j,k))) * & 369 (2.0*g%IdyCv(i,j) / ((dp(i,j) + dp(i,j+1)) + dp\_neglect)) 370 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 371 372 \textcolor{keywordflow}{if} (cs%GFS\_scale < 1.0) \textcolor{keywordflow}{then} 373 \textcolor{comment}{! Adjust the Montgomery potential to make this a reduced gravity model.} 374 \textcolor{comment}{!$OMP parallel do default(shared)} 375 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 376 pfu(i,j,k) = pfu(i,j,k) - (dm(i+1,j) - dm(i,j)) * g%IdxCu(i,j) 377 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 378 \textcolor{comment}{!$OMP parallel do default(shared)} 379 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 380 pfv(i,j,k) = pfv(i,j,k) - (dm(i,j+1) - dm(i,j)) * g%IdyCv(i,j) 381 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 382 \textcolor{keywordflow}{ endif} 383 \textcolor{keywordflow}{ enddo} 384 385 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(pbce)) \textcolor{keywordflow}{then} 386 \textcolor{keyword}{call }set\_pbce\_nonbouss(p, tv\_tmp, g, gv, us, cs%GFS\_scale, pbce) 387 \textcolor{keywordflow}{ endif} 388 389 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(eta)) \textcolor{keywordflow}{then} 390 pa\_to\_h = 1.0 / (gv%g\_Earth * gv%H\_to\_RZ) 391 \textcolor{keywordflow}{if} (use\_p\_atm) \textcolor{keywordflow}{then} 392 \textcolor{comment}{!$OMP parallel do default(shared)} 393 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 394 eta(i,j) = (p(i,j,nz+1) - p\_atm(i,j))*pa\_to\_h \textcolor{comment}{! eta has the same units as h.} 395 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 396 \textcolor{keywordflow}{else} 397 \textcolor{comment}{!$OMP parallel do default(shared)} 398 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 399 eta(i,j) = p(i,j,nz+1)*pa\_to\_h \textcolor{comment}{! eta has the same units as h.} 400 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 401 \textcolor{keywordflow}{ endif} 402 \textcolor{keywordflow}{ endif} 403 404 \textcolor{keywordflow}{if} (cs%id\_e\_tidal>0) \textcolor{keyword}{call }post\_data(cs%id\_e\_tidal, e\_tidal, cs%diag) 405 \end{DoxyCode}