\hypertarget{namespacemom__pressureforce__mont}{}\section{mom\+\_\+pressureforce\+\_\+mont Module Reference} \label{namespacemom__pressureforce__mont}\index{mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}} \subsection{Detailed Description} Provides the Montgomery potential form of pressure gradient. Provides the Boussunesq and non-\/\+Boussinesq forms of the horizontal accelerations due to pressure gradients using the Montgomery potential. A second-\/order accurate, centered scheme is used. If a split time stepping scheme is used, the vertical decomposition into barotropic and baroclinic contributions described by Hallberg (J Comp Phys 1997) is used. With a nonlinear equation of state, compressibility is added along the lines proposed by Sun et al. (J\+PO 1999), but with compressibility coefficients based on a fit to a user-\/provided reference profile. \subsection*{Data Types} \begin{DoxyCompactItemize} \item type \hyperlink{structmom__pressureforce__mont_1_1pressureforce__mont__cs}{pressureforce\+\_\+mont\+\_\+cs} \begin{DoxyCompactList}\small\item\em Control structure for the Montgomery potential form of pressure gradient. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection*{Functions/\+Subroutines} \begin{DoxyCompactItemize} \item subroutine, public \hyperlink{namespacemom__pressureforce__mont_a6880a913a82b65eb65a728abb487ef91}{pressureforce\+\_\+mont\+\_\+nonbouss} (h, tv, P\+Fu, P\+Fv, G, GV, US, CS, p\+\_\+atm, pbce, eta) \begin{DoxyCompactList}\small\item\em Non-\/\+Boussinesq Montgomery-\/potential form of pressure gradient. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__pressureforce__mont_a5f96546655ff2d3fc9090ba2467ef521}{pressureforce\+\_\+mont\+\_\+bouss} (h, tv, P\+Fu, P\+Fv, G, GV, US, CS, p\+\_\+atm, pbce, eta) \begin{DoxyCompactList}\small\item\em Boussinesq Montgomery-\/potential form of pressure gradient. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__pressureforce__mont_a0779efd30a447c2bc20294c46aeea180}{set\+\_\+pbce\+\_\+bouss} (e, tv, G, GV, US, Rho0, G\+F\+S\+\_\+scale, pbce, rho\+\_\+star) \begin{DoxyCompactList}\small\item\em Determines the partial derivative of the acceleration due to pressure forces with the free surface height. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__pressureforce__mont_a297cdf6e4eb83d250444c5f527a8a232}{set\+\_\+pbce\+\_\+nonbouss} (p, tv, G, GV, US, G\+F\+S\+\_\+scale, pbce, alpha\+\_\+star) \begin{DoxyCompactList}\small\item\em Determines the partial derivative of the acceleration due to pressure forces with the column mass. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__pressureforce__mont_a55f5cfcd322f96938285ff4eb7646d14}{pressureforce\+\_\+mont\+\_\+init} (Time, G, GV, US, param\+\_\+file, diag, CS, tides\+\_\+\+C\+Sp) \begin{DoxyCompactList}\small\item\em Initialize the Montgomery-\/potential form of P\+GF control structure. \end{DoxyCompactList}\item subroutine, public \hyperlink{namespacemom__pressureforce__mont_a1aee187502d5f5b78db6c1bfb4fe72a6}{pressureforce\+\_\+mont\+\_\+end} (CS) \begin{DoxyCompactList}\small\item\em Deallocates the Montgomery-\/potential form of P\+GF control structure. \end{DoxyCompactList}\end{DoxyCompactItemize} \subsection{Function/\+Subroutine Documentation} \mbox{\Hypertarget{namespacemom__pressureforce__mont_a5f96546655ff2d3fc9090ba2467ef521}\label{namespacemom__pressureforce__mont_a5f96546655ff2d3fc9090ba2467ef521}} \index{mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}!pressureforce\+\_\+mont\+\_\+bouss@{pressureforce\+\_\+mont\+\_\+bouss}} \index{pressureforce\+\_\+mont\+\_\+bouss@{pressureforce\+\_\+mont\+\_\+bouss}!mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}} \subsubsection{\texorpdfstring{pressureforce\+\_\+mont\+\_\+bouss()}{pressureforce\_mont\_bouss()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+mont\+::pressureforce\+\_\+mont\+\_\+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(\hyperlink{structmom__pressureforce__mont_1_1pressureforce__mont__cs}{pressureforce\+\_\+mont\+\_\+cs}), pointer}]{CS, }\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 Montgomery-\/potential form of pressure gradient. Determines the acceleration due to pressure forces. 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 due to pressure gradients (equal to -\/d\+M/dx) \mbox{[}L T-\/2 $\sim$$>$ m s-\/2\mbox{]}.\\ \hline \mbox{\tt out} & {\em pfv} & Meridional acceleration due to pressure gradients (equal to -\/d\+M/dy) \mbox{[}L T-\/2 $\sim$$>$ m s2\mbox{]}.\\ \hline & {\em cs} & Control structure for Montgomery potential P\+GF\\ \hline & {\em p\+\_\+atm} & The pressure at the ice-\/ocean or atmosphere-\/ocean \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}.\\ \hline \mbox{\tt out} & {\em pbce} & The baroclinic pressure anomaly in each layer due to free surface height anomalies \mbox{[}L2 T-\/2 H-\/1 $\sim$$>$ m s-\/2\mbox{]}.\\ \hline \mbox{\tt out} & {\em eta} & Free surface height \mbox{[}H $\sim$$>$ m\mbox{]}. \\ \hline \end{DoxyParams} Definition at line 361 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+Montgomery.\+F90. \begin{DoxyCode} 361 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< Ocean grid structure.} 362 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< Vertical grid structure.} 363 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 364 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thickness [H ~> m].} 365 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< Thermodynamic variables.} 366 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: pfu\textcolor{comment}{ !< Zonal acceleration due to pressure gradients} 367 \textcolor{comment}{ !! (equal to -dM/dx) [L T-2 ~> m s-2].} 368 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: pfv\textcolor{comment}{ !< Meridional acceleration due to pressure gradients} 369 \textcolor{comment}{ !! (equal to -dM/dy) [L T-2 ~> m s2].} 370 \textcolor{keywordtype}{type}(pressureforce\_mont\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< Control structure for Montgomery potential PGF} 371 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: p\_atm\textcolor{comment}{ !< The pressure at the ice-ocean or} 372 \textcolor{comment}{ !! atmosphere-ocean [R L2 T-2 ~> Pa].} 373 \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 anomaly in} 374 \textcolor{comment}{ !! each layer due to free surface height anomalies} 375 \textcolor{comment}{ !! [L2 T-2 H-1 ~> m s-2].} 376 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: eta\textcolor{comment}{ !< Free surface height [H ~> m].} 377 \textcolor{comment}{! Local variables} 378 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))} :: & 379 m, & \textcolor{comment}{! The Montgomery potential, M = (p/rho + gz) [L2 T-2 ~> m2 s-2].} 380 rho\_star \textcolor{comment}{! In-situ density divided by the derivative with depth of the} 381 \textcolor{comment}{! corrected e times (G\_Earth/Rho0) [m2 Z-1 s-2 ~> m s-2].} 382 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)} :: e \textcolor{comment}{! Interface height in m.} 383 \textcolor{comment}{! e may be adjusted (with a nonlinear equation of state) so that} 384 \textcolor{comment}{! its derivative compensates for the adiabatic compressibility} 385 \textcolor{comment}{! in seawater, but e will still be close to the interface depth.} 386 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{target} :: & 387 t\_tmp, & \textcolor{comment}{! Temporary array of temperatures where layers that are lighter} 388 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [degC].} 389 s\_tmp \textcolor{comment}{! Temporary array of salinities where layers that are lighter} 390 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [ppt].} 391 392 \textcolor{keywordtype}{real} :: rho\_cv\_bl(szi\_(g)) \textcolor{comment}{! The coordinate potential density in} 393 \textcolor{comment}{! the deepest variable density near-surface layer [R ~> kg m-3].} 394 \textcolor{keywordtype}{real} :: h\_star(szi\_(g),szj\_(g)) \textcolor{comment}{! Layer thickness after compensation} 395 \textcolor{comment}{! for compressibility [Z ~> m].} 396 \textcolor{keywordtype}{real} :: e\_tidal(szi\_(g),szj\_(g)) \textcolor{comment}{! Bottom geopotential anomaly due to tidal} 397 \textcolor{comment}{! forces from astronomical sources and self-} 398 \textcolor{comment}{! attraction and loading, in depth units [Z ~> m].} 399 \textcolor{keywordtype}{real} :: p\_ref(szi\_(g)) \textcolor{comment}{! The pressure used to calculate the coordinate} 400 \textcolor{comment}{! density [R L2 T-2 ~> Pa] (usually 2e7 Pa = 2000 dbar).} 401 \textcolor{keywordtype}{real} :: i\_rho0 \textcolor{comment}{! 1/Rho0 [R-1 ~> m3 kg-1].} 402 \textcolor{keywordtype}{real} :: g\_rho0 \textcolor{comment}{! G\_Earth / Rho0 [L2 Z-1 T-2 R-1 ~> m4 s-2 kg-1].} 403 \textcolor{keywordtype}{real} :: pfu\_bc, pfv\_bc \textcolor{comment}{! The pressure gradient force due to along-layer} 404 \textcolor{comment}{! compensated density gradients [L T-2 ~> m s-2]} 405 \textcolor{keywordtype}{real} :: h\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 406 \textcolor{comment}{! in roundoff and can be neglected [Z ~> m].} 407 \textcolor{keywordtype}{logical} :: use\_p\_atm \textcolor{comment}{! If true, use the atmospheric pressure.} 408 \textcolor{keywordtype}{logical} :: use\_eos \textcolor{comment}{! If true, density is calculated from T & S using} 409 \textcolor{comment}{! an equation of state.} 410 \textcolor{keywordtype}{logical} :: is\_split \textcolor{comment}{! A flag indicating whether the pressure} 411 \textcolor{comment}{! gradient terms are to be split into} 412 \textcolor{comment}{! barotropic and baroclinic pieces.} 413 \textcolor{keywordtype}{type}(thermo\_var\_ptrs) :: tv\_tmp\textcolor{comment}{! A structure of temporary T & S.} 414 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The computational domain for the equation of state} 415 \textcolor{keywordtype}{integer} :: is, ie, js, je, isq, ieq, jsq, jeq, nz, nkmb 416 \textcolor{keywordtype}{integer} :: i, j, k 417 418 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke 419 nkmb=gv%nk\_rho\_varies 420 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB 421 eosdom(1) = isq - (g%isd-1) ; eosdom(2) = g%iec+1 - (g%isd-1) 422 423 use\_p\_atm = .false. 424 \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} 425 is\_split = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(pbce)) is\_split = .true. 426 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 427 428 \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal, & 429 \textcolor{stringliteral}{"MOM\_PressureForce\_Mont: Module must be initialized before it is used."}) 430 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 431 \textcolor{keywordflow}{if} (query\_compressible(tv%eqn\_of\_state)) \textcolor{keyword}{call }mom\_error(fatal, & 432 \textcolor{stringliteral}{"PressureForce\_Mont\_Bouss: The Montgomery form of the pressure force "}//& 433 \textcolor{stringliteral}{"can no longer be used with a compressible EOS. Use #define ANALYTIC\_FV\_PGF."}) 434 \textcolor{keywordflow}{ endif} 435 436 h\_neglect = gv%H\_subroundoff * gv%H\_to\_Z 437 i\_rho0 = 1.0/cs%Rho0 438 g\_rho0 = gv%g\_Earth / gv%Rho0 439 440 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 441 \textcolor{comment}{! Determine the surface height anomaly for calculating self attraction} 442 \textcolor{comment}{! and loading. This should really be based on bottom pressure anomalies,} 443 \textcolor{comment}{! but that is not yet implemented, and the current form is correct for} 444 \textcolor{comment}{! barotropic tides.} 445 \textcolor{comment}{!$OMP parallel do default(shared)} 446 \textcolor{keywordflow}{do} j=jsq,jeq+1 447 \textcolor{keywordflow}{do} i=isq,ieq+1 ; e(i,j,1) = -g%bathyT(i,j) ;\textcolor{keywordflow}{ enddo} 448 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 449 e(i,j,1) = e(i,j,1) + h(i,j,k)*gv%H\_to\_Z 450 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 451 \textcolor{keywordflow}{ enddo} 452 \textcolor{keyword}{call }calc\_tidal\_forcing(cs%Time, e(:,:,1), e\_tidal, g, cs%tides\_CSp, m\_to\_z=us%m\_to\_Z) 453 \textcolor{keywordflow}{ endif} 454 455 \textcolor{comment}{! Here layer interface heights, e, are calculated.} 456 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 457 \textcolor{comment}{!$OMP parallel do default(shared)} 458 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 459 e(i,j,nz+1) = -(g%bathyT(i,j) + e\_tidal(i,j)) 460 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 461 \textcolor{keywordflow}{else} 462 \textcolor{comment}{!$OMP parallel do default(shared)} 463 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 464 e(i,j,nz+1) = -g%bathyT(i,j) 465 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 466 \textcolor{keywordflow}{ endif} 467 \textcolor{comment}{!$OMP parallel do default(shared)} 468 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} k=nz,1,-1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 469 e(i,j,k) = e(i,j,k+1) + h(i,j,k)*gv%H\_to\_Z 470 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 471 472 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 473 \textcolor{comment}{! Calculate in-situ densities (rho\_star).} 474 475 \textcolor{comment}{! With a bulk mixed layer, replace the T & S of any layers that are} 476 \textcolor{comment}{! lighter than the the buffer layer with the properties of the buffer} 477 \textcolor{comment}{! layer. These layers will be massless anyway, and it avoids any} 478 \textcolor{comment}{! formal calculations with hydrostatically unstable profiles.} 479 480 \textcolor{keywordflow}{if} (nkmb>0) \textcolor{keywordflow}{then} 481 tv\_tmp%T => t\_tmp ; tv\_tmp%S => s\_tmp 482 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 483 484 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p\_ref(i) = tv%P\_Ref ;\textcolor{keywordflow}{ enddo} 485 \textcolor{comment}{!$OMP parallel do default(shared) private(Rho\_cv\_BL)} 486 \textcolor{keywordflow}{do} j=jsq,jeq+1 487 \textcolor{keywordflow}{do} k=1,nkmb ; \textcolor{keywordflow}{do} i=isq,ieq+1 488 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 489 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 490 \textcolor{keyword}{call }calculate\_density(tv%T(:,j,nkmb), tv%S(:,j,nkmb), p\_ref, rho\_cv\_bl(:), & 491 tv%eqn\_of\_state, eosdom) 492 493 \textcolor{keywordflow}{do} k=nkmb+1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 494 \textcolor{keywordflow}{if} (gv%Rlay(k) < rho\_cv\_bl(i)) \textcolor{keywordflow}{then} 495 tv\_tmp%T(i,j,k) = tv%T(i,j,nkmb) ; tv\_tmp%S(i,j,k) = tv%S(i,j,nkmb) 496 \textcolor{keywordflow}{else} 497 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 498 \textcolor{keywordflow}{ endif} 499 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 500 \textcolor{keywordflow}{ enddo} 501 \textcolor{keywordflow}{else} 502 tv\_tmp%T => tv%T ; tv\_tmp%S => tv%S 503 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 504 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p\_ref(i) = 0.0 ;\textcolor{keywordflow}{ enddo} 505 \textcolor{keywordflow}{ endif} 506 507 \textcolor{comment}{! This no longer includes any pressure dependency, since this routine} 508 \textcolor{comment}{! will come down with a fatal error if there is any compressibility.} 509 \textcolor{comment}{!$OMP parallel do default(shared)} 510 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} j=jsq,jeq+1 511 \textcolor{keyword}{call }calculate\_density(tv\_tmp%T(:,j,k), tv\_tmp%S(:,j,k), p\_ref, rho\_star(:,j,k), & 512 tv%eqn\_of\_state, eosdom) 513 \textcolor{keywordflow}{do} i=isq,ieq+1 ; rho\_star(i,j,k) = g\_rho0*rho\_star(i,j,k) ;\textcolor{keywordflow}{ enddo} 514 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 515 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 516 517 \textcolor{comment}{! Here the layer Montgomery potentials, M, are calculated.} 518 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 519 \textcolor{comment}{!$OMP parallel do default(shared)} 520 \textcolor{keywordflow}{do} j=jsq,jeq+1 521 \textcolor{keywordflow}{do} i=isq,ieq+1 522 m(i,j,1) = cs%GFS\_scale * (rho\_star(i,j,1) * e(i,j,1)) 523 \textcolor{keywordflow}{if} (use\_p\_atm) m(i,j,1) = m(i,j,1) + p\_atm(i,j) * i\_rho0 524 \textcolor{keywordflow}{ enddo} 525 \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 526 m(i,j,k) = m(i,j,k-1) + (rho\_star(i,j,k) - rho\_star(i,j,k-1)) * e(i,j,k) 527 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 528 \textcolor{keywordflow}{ enddo} 529 \textcolor{keywordflow}{else} \textcolor{comment}{! not use\_EOS} 530 \textcolor{comment}{!$OMP parallel do default(shared)} 531 \textcolor{keywordflow}{do} j=jsq,jeq+1 532 \textcolor{keywordflow}{do} i=isq,ieq+1 533 m(i,j,1) = gv%g\_prime(1) * e(i,j,1) 534 \textcolor{keywordflow}{if} (use\_p\_atm) m(i,j,1) = m(i,j,1) + p\_atm(i,j) * i\_rho0 535 \textcolor{keywordflow}{ enddo} 536 \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 537 m(i,j,k) = m(i,j,k-1) + gv%g\_prime(k) * e(i,j,k) 538 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 539 \textcolor{keywordflow}{ enddo} 540 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 541 542 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(pbce)) \textcolor{keywordflow}{then} 543 \textcolor{keyword}{call }set\_pbce\_bouss(e, tv\_tmp, g, gv, us, cs%Rho0, cs%GFS\_scale, pbce, rho\_star) 544 \textcolor{keywordflow}{ endif} 545 546 \textcolor{comment}{! Calculate the pressure force. On a Cartesian grid,} 547 \textcolor{comment}{! PFu = - dM/dx and PFv = - dM/dy.} 548 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 549 \textcolor{comment}{!$OMP parallel do default(shared) private(h\_star,PFu\_bc,PFv\_bc)} 550 \textcolor{keywordflow}{do} k=1,nz 551 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 552 h\_star(i,j) = (e(i,j,k) - e(i,j,k+1)) + h\_neglect 553 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 554 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 555 pfu\_bc = -1.0*(rho\_star(i+1,j,k) - rho\_star(i,j,k)) * (g%IdxCu(i,j) * & 556 ((h\_star(i,j) * h\_star(i+1,j) - (e(i,j,k) * h\_star(i+1,j) + & 557 e(i+1,j,k) * h\_star(i,j))) / (h\_star(i,j) + h\_star(i+1,j)))) 558 pfu(i,j,k) = -(m(i+1,j,k) - m(i,j,k)) * g%IdxCu(i,j) + pfu\_bc 559 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs%PFu\_bc)) cs%PFu\_bc(i,j,k) = pfu\_bc 560 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 561 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 562 pfv\_bc = -1.0*(rho\_star(i,j+1,k) - rho\_star(i,j,k)) * (g%IdyCv(i,j) * & 563 ((h\_star(i,j) * h\_star(i,j+1) - (e(i,j,k) * h\_star(i,j+1) + & 564 e(i,j+1,k) * h\_star(i,j))) / (h\_star(i,j) + h\_star(i,j+1)))) 565 pfv(i,j,k) = -(m(i,j+1,k) - m(i,j,k)) * g%IdyCv(i,j) + pfv\_bc 566 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs%PFv\_bc)) cs%PFv\_bc(i,j,k) = pfv\_bc 567 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 568 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! k-loop} 569 \textcolor{keywordflow}{else} \textcolor{comment}{! .not. use\_EOS} 570 \textcolor{comment}{!$OMP parallel do default(shared)} 571 \textcolor{keywordflow}{do} k=1,nz 572 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 573 pfu(i,j,k) = -(m(i+1,j,k) - m(i,j,k)) * g%IdxCu(i,j) 574 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 575 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 576 pfv(i,j,k) = -(m(i,j+1,k) - m(i,j,k)) * g%IdyCv(i,j) 577 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 578 \textcolor{keywordflow}{ enddo} 579 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 580 581 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(eta)) \textcolor{keywordflow}{then} 582 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 583 \textcolor{comment}{! eta is the sea surface height relative to a time-invariant geoid, for} 584 \textcolor{comment}{! comparison with what is used for eta in btstep. See how e was calculated} 585 \textcolor{comment}{! about 200 lines above.} 586 \textcolor{comment}{!$OMP parallel do default(shared)} 587 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 588 eta(i,j) = e(i,j,1)*gv%Z\_to\_H + e\_tidal(i,j)*gv%Z\_to\_H 589 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 590 \textcolor{keywordflow}{else} 591 \textcolor{comment}{!$OMP parallel do default(shared)} 592 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 593 eta(i,j) = e(i,j,1)*gv%Z\_to\_H 594 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 595 \textcolor{keywordflow}{ endif} 596 \textcolor{keywordflow}{ endif} 597 598 \textcolor{keywordflow}{if} (cs%id\_PFu\_bc>0) \textcolor{keyword}{call }post\_data(cs%id\_PFu\_bc, cs%PFu\_bc, cs%diag) 599 \textcolor{keywordflow}{if} (cs%id\_PFv\_bc>0) \textcolor{keyword}{call }post\_data(cs%id\_PFv\_bc, cs%PFv\_bc, cs%diag) 600 \textcolor{keywordflow}{if} (cs%id\_e\_tidal>0) \textcolor{keyword}{call }post\_data(cs%id\_e\_tidal, e\_tidal, cs%diag) 601 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__mont_a1aee187502d5f5b78db6c1bfb4fe72a6}\label{namespacemom__pressureforce__mont_a1aee187502d5f5b78db6c1bfb4fe72a6}} \index{mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}!pressureforce\+\_\+mont\+\_\+end@{pressureforce\+\_\+mont\+\_\+end}} \index{pressureforce\+\_\+mont\+\_\+end@{pressureforce\+\_\+mont\+\_\+end}!mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}} \subsubsection{\texorpdfstring{pressureforce\+\_\+mont\+\_\+end()}{pressureforce\_mont\_end()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+mont\+::pressureforce\+\_\+mont\+\_\+end (\begin{DoxyParamCaption}\item[{type(\hyperlink{structmom__pressureforce__mont_1_1pressureforce__mont__cs}{pressureforce\+\_\+mont\+\_\+cs}), pointer}]{CS }\end{DoxyParamCaption})} Deallocates the Montgomery-\/potential form of P\+GF control structure. \begin{DoxyParams}{Parameters} {\em cs} & Control structure for Montgomery potential P\+GF \\ \hline \end{DoxyParams} Definition at line 890 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+Montgomery.\+F90. \begin{DoxyCode} 890 \textcolor{keywordtype}{type}(pressureforce\_mont\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< Control structure for Montgomery potential PGF} 891 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{deallocate}(cs) \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__mont_a55f5cfcd322f96938285ff4eb7646d14}\label{namespacemom__pressureforce__mont_a55f5cfcd322f96938285ff4eb7646d14}} \index{mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}!pressureforce\+\_\+mont\+\_\+init@{pressureforce\+\_\+mont\+\_\+init}} \index{pressureforce\+\_\+mont\+\_\+init@{pressureforce\+\_\+mont\+\_\+init}!mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}} \subsubsection{\texorpdfstring{pressureforce\+\_\+mont\+\_\+init()}{pressureforce\_mont\_init()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+mont\+::pressureforce\+\_\+mont\+\_\+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(\hyperlink{structmom__pressureforce__mont_1_1pressureforce__mont__cs}{pressureforce\+\_\+mont\+\_\+cs}), pointer}]{CS, }\item[{type(tidal\+\_\+forcing\+\_\+cs), optional, pointer}]{tides\+\_\+\+C\+Sp }\end{DoxyParamCaption})} Initialize the Montgomery-\/potential form of P\+GF 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} & Montgomery P\+GF control structure\\ \hline & {\em tides\+\_\+csp} & Tides control structure \\ \hline \end{DoxyParams} Definition at line 822 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+Montgomery.\+F90. \begin{DoxyCode} 822 \textcolor{keywordtype}{type}(time\_type), \textcolor{keywordtype}{target}, \textcolor{keywordtype}{intent(in)} :: time\textcolor{comment}{ !< Current model time} 823 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< ocean grid structure} 824 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< Vertical grid structure} 825 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 826 \textcolor{keywordtype}{type}(param\_file\_type), \textcolor{keywordtype}{intent(in)} :: param\_file\textcolor{comment}{ !< Parameter file handles} 827 \textcolor{keywordtype}{type}(diag\_ctrl), \textcolor{keywordtype}{target}, \textcolor{keywordtype}{intent(inout)} :: diag\textcolor{comment}{ !< Diagnostics control structure} 828 \textcolor{keywordtype}{type}(pressureforce\_mont\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< Montgomery PGF control structure} 829 \textcolor{keywordtype}{type}(tidal\_forcing\_cs), \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: tides\_csp\textcolor{comment}{ !< Tides control structure} 830 831 \textcolor{comment}{! Local variables} 832 \textcolor{keywordtype}{logical} :: use\_temperature, use\_eos 833 \textcolor{comment}{! This include declares and sets the variable "version".} 834 \textcolor{preprocessor}{# include "version\_variable.h"} 835 \textcolor{preprocessor}{} \textcolor{keywordtype}{character(len=40)} :: mdl \textcolor{comment}{! This module's name.} 836 837 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs)) \textcolor{keywordflow}{then} 838 \textcolor{keyword}{call }mom\_error(warning, \textcolor{stringliteral}{"PressureForce\_init called with an associated "}// & 839 \textcolor{stringliteral}{"control structure."}) 840 \textcolor{keywordflow}{return} 841 \textcolor{keywordflow}{else} ; \textcolor{keyword}{allocate}(cs) ;\textcolor{keywordflow}{ endif} 842 843 cs%diag => diag ; cs%Time => time 844 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(tides\_csp)) \textcolor{keywordflow}{then} 845 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(tides\_csp)) cs%tides\_CSp => tides\_csp 846 \textcolor{keywordflow}{ endif} 847 848 mdl = \textcolor{stringliteral}{"MOM\_PressureForce\_Mont"} 849 \textcolor{keyword}{call }log\_version(param\_file, mdl, version, \textcolor{stringliteral}{""}) 850 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"RHO\_0"}, cs%Rho0, & 851 \textcolor{stringliteral}{"The mean ocean density used with BOUSSINESQ true to "}//& 852 \textcolor{stringliteral}{"calculate accelerations and the mass for conservation "}//& 853 \textcolor{stringliteral}{"properties, or with BOUSSINSEQ false to convert some "}//& 854 \textcolor{stringliteral}{"parameters from vertical units of m to kg m-2."}, & 855 units=\textcolor{stringliteral}{"kg m-3"}, default=1035.0, scale=us%R\_to\_kg\_m3) 856 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"TIDES"}, cs%tides, & 857 \textcolor{stringliteral}{"If true, apply tidal momentum forcing."}, default=.false.) 858 \textcolor{keyword}{call }get\_param(param\_file, mdl, \textcolor{stringliteral}{"USE\_EOS"}, use\_eos, default=.true., & 859 do\_not\_log=.true.) \textcolor{comment}{! Input for diagnostic use only.} 860 861 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 862 cs%id\_PFu\_bc = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'PFu\_bc'}, diag%axesCuL, time, & 863 \textcolor{stringliteral}{'Density Gradient Zonal Pressure Force Accel.'}, \textcolor{stringliteral}{"meter second-2"}, conversion=us%L\_T2\_to\_m\_s2) 864 cs%id\_PFv\_bc = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'PFv\_bc'}, diag%axesCvL, time, & 865 \textcolor{stringliteral}{'Density Gradient Meridional Pressure Force Accel.'}, \textcolor{stringliteral}{"meter second-2"}, conversion=us%L\_T2\_to\_m\_s2) 866 \textcolor{keywordflow}{if} (cs%id\_PFu\_bc > 0) \textcolor{keywordflow}{then} 867 \textcolor{keyword}{call }safe\_alloc\_ptr(cs%PFu\_bc,g%IsdB,g%IedB,g%jsd,g%jed,g%ke) 868 cs%PFu\_bc(:,:,:) = 0.0 869 \textcolor{keywordflow}{ endif} 870 \textcolor{keywordflow}{if} (cs%id\_PFv\_bc > 0) \textcolor{keywordflow}{then} 871 \textcolor{keyword}{call }safe\_alloc\_ptr(cs%PFv\_bc,g%isd,g%ied,g%JsdB,g%JedB,g%ke) 872 cs%PFv\_bc(:,:,:) = 0.0 873 \textcolor{keywordflow}{ endif} 874 \textcolor{keywordflow}{ endif} 875 876 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 877 cs%id\_e\_tidal = register\_diag\_field(\textcolor{stringliteral}{'ocean\_model'}, \textcolor{stringliteral}{'e\_tidal'}, diag%axesT1, & 878 time, \textcolor{stringliteral}{'Tidal Forcing Astronomical and SAL Height Anomaly'}, \textcolor{stringliteral}{'meter'}, conversion=us%Z\_to\_m) 879 \textcolor{keywordflow}{ endif} 880 881 cs%GFS\_scale = 1.0 882 \textcolor{keywordflow}{if} (gv%g\_prime(1) /= gv%g\_Earth) cs%GFS\_scale = gv%g\_prime(1) / gv%g\_Earth 883 884 \textcolor{keyword}{call }log\_param(param\_file, mdl, \textcolor{stringliteral}{"GFS / G\_EARTH"}, cs%GFS\_scale) 885 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__mont_a6880a913a82b65eb65a728abb487ef91}\label{namespacemom__pressureforce__mont_a6880a913a82b65eb65a728abb487ef91}} \index{mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}!pressureforce\+\_\+mont\+\_\+nonbouss@{pressureforce\+\_\+mont\+\_\+nonbouss}} \index{pressureforce\+\_\+mont\+\_\+nonbouss@{pressureforce\+\_\+mont\+\_\+nonbouss}!mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}} \subsubsection{\texorpdfstring{pressureforce\+\_\+mont\+\_\+nonbouss()}{pressureforce\_mont\_nonbouss()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+mont\+::pressureforce\+\_\+mont\+\_\+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(\hyperlink{structmom__pressureforce__mont_1_1pressureforce__mont__cs}{pressureforce\+\_\+mont\+\_\+cs}), pointer}]{CS, }\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 Montgomery-\/potential form of pressure gradient. Determines the acceleration due to pressure forces in a non-\/\+Boussinesq fluid using the compressibility compensated (if appropriate) Montgomery-\/potential form described in Hallberg (Ocean Mod., 2005). 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 m-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em tv} & Thermodynamic variables.\\ \hline \mbox{\tt out} & {\em pfu} & Zonal acceleration due to pressure gradients (equal to -\/d\+M/dx) \mbox{[}L T-\/2 $\sim$$>$ m s-\/2\mbox{]}.\\ \hline \mbox{\tt out} & {\em pfv} & Meridional acceleration due to pressure gradients (equal to -\/d\+M/dy) \mbox{[}L T-\/2 $\sim$$>$ m s-\/2\mbox{]}.\\ \hline & {\em cs} & Control structure for Montgomery potential P\+GF\\ \hline & {\em p\+\_\+atm} & The pressure at the ice-\/ocean or atmosphere-\/ocean \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}.\\ \hline \mbox{\tt out} & {\em pbce} & The baroclinic pressure anomaly in\\ \hline \mbox{\tt out} & {\em eta} & Free surface height \mbox{[}H $\sim$$>$ kg m-\/1\mbox{]}. \\ \hline \end{DoxyParams} Definition at line 64 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+Montgomery.\+F90. \begin{DoxyCode} 64 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< Ocean grid structure.} 65 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< Vertical grid structure.} 66 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 67 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(in)} :: h\textcolor{comment}{ !< Layer thickness, [H ~> kg m-2].} 68 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< Thermodynamic variables.} 69 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZIB\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: pfu\textcolor{comment}{ !< Zonal acceleration due to pressure gradients} 70 \textcolor{comment}{ !! (equal to -dM/dx) [L T-2 ~> m s-2].} 71 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJB\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: pfv\textcolor{comment}{ !< Meridional acceleration due to pressure gradients} 72 \textcolor{comment}{ !! (equal to -dM/dy) [L T-2 ~> m s-2].} 73 \textcolor{keywordtype}{type}(pressureforce\_mont\_cs), \textcolor{keywordtype}{pointer} :: cs\textcolor{comment}{ !< Control structure for Montgomery potential PGF} 74 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(:,:)}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{pointer} :: p\_atm\textcolor{comment}{ !< The pressure at the ice-ocean or} 75 \textcolor{comment}{ !! atmosphere-ocean [R L2 T-2 ~> Pa].} 76 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, & 77 \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: pbce\textcolor{comment}{ !< The baroclinic pressure anomaly in} 78 \textcolor{comment}{ !! each layer due to free surface height anomalies,} 79 \textcolor{comment}{ !! [L2 T-2 H-1 ~> m s-2 or m4 kg-1 s-2].} 80 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(out)} :: eta\textcolor{comment}{ !< Free surface height [H ~> kg m-1].} 81 82 \textcolor{comment}{! Local variables} 83 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))} :: & 84 m, & \textcolor{comment}{! The Montgomery potential, M = (p/rho + gz) [L2 T-2 ~> m2 s-2].} 85 alpha\_star, & \textcolor{comment}{! Compression adjusted specific volume [R-1 ~> m3 kg-1].} 86 dz\_geo \textcolor{comment}{! The change in geopotential across a layer [L2 T-2 ~> m2 s-2].} 87 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)} :: p \textcolor{comment}{! Interface pressure [R L2 T-2 ~> Pa].} 88 \textcolor{comment}{! p may be adjusted (with a nonlinear equation of state) so that} 89 \textcolor{comment}{! its derivative compensates for the adiabatic compressibility} 90 \textcolor{comment}{! in seawater, but p will still be close to the pressure.} 91 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{target} :: & 92 t\_tmp, & \textcolor{comment}{! Temporary array of temperatures where layers that are lighter} 93 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [degC].} 94 s\_tmp \textcolor{comment}{! Temporary array of salinities where layers that are lighter} 95 \textcolor{comment}{! than the mixed layer have the mixed layer's properties [ppt].} 96 97 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G))} :: rho\_cv\_bl \textcolor{comment}{! The coordinate potential density in the} 98 \textcolor{comment}{! deepest variable density near-surface layer [R ~> kg m-3].} 99 100 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))} :: & 101 dm, & \textcolor{comment}{! A barotropic correction to the Montgomery potentials to enable the use} 102 \textcolor{comment}{! of a reduced gravity form of the equations [L2 T-2 ~> m2 s-2].} 103 dp\_star, & \textcolor{comment}{! Layer thickness after compensation for compressibility [R L2 T-2 ~> Pa].} 104 ssh, & \textcolor{comment}{! The sea surface height anomaly, in depth units [Z ~> m].} 105 e\_tidal, & \textcolor{comment}{! Bottom geopotential anomaly due to tidal forces from} 106 \textcolor{comment}{! astronomical sources and self-attraction and loading [Z ~> m].} 107 geopot\_bot \textcolor{comment}{! Bottom geopotential relative to time-mean sea level,} 108 \textcolor{comment}{! including any tidal contributions [L2 T-2 ~> m2 s-2].} 109 \textcolor{keywordtype}{real} :: p\_ref(szi\_(g)) \textcolor{comment}{! The pressure used to calculate the coordinate} 110 \textcolor{comment}{! density [R L2 T-2 ~> Pa] (usually 2e7 Pa = 2000 dbar).} 111 \textcolor{keywordtype}{real} :: rho\_in\_situ(szi\_(g)) \textcolor{comment}{!In-situ density of a layer [R ~> kg m-3].} 112 \textcolor{keywordtype}{real} :: pfu\_bc, pfv\_bc \textcolor{comment}{! The pressure gradient force due to along-layer} 113 \textcolor{comment}{! compensated density gradients [L T-2 ~> m s-2]} 114 \textcolor{keywordtype}{real} :: dp\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 115 \textcolor{comment}{! in roundoff and can be neglected [R L2 T-2 ~> Pa].} 116 \textcolor{keywordtype}{logical} :: use\_p\_atm \textcolor{comment}{! If true, use the atmospheric pressure.} 117 \textcolor{keywordtype}{logical} :: use\_eos \textcolor{comment}{! If true, density is calculated from T & S using an equation of state.} 118 \textcolor{keywordtype}{logical} :: is\_split \textcolor{comment}{! A flag indicating whether the pressure gradient terms are to be} 119 \textcolor{comment}{! split into barotropic and baroclinic pieces.} 120 \textcolor{keywordtype}{type}(thermo\_var\_ptrs) :: tv\_tmp\textcolor{comment}{! A structure of temporary T & S.} 121 122 \textcolor{keywordtype}{real} :: i\_gearth \textcolor{comment}{! The inverse of g\_Earth [T2 Z L-2 ~> s2 m-1]} 123 \textcolor{comment}{! real :: dalpha} 124 \textcolor{keywordtype}{real} :: pa\_to\_h \textcolor{comment}{! A factor to convert from R L2 T-2 to the thickness units (H).} 125 \textcolor{keywordtype}{real} :: alpha\_lay(szk\_(g)) \textcolor{comment}{! The specific volume of each layer [R-1 ~> m3 kg-1].} 126 \textcolor{keywordtype}{real} :: dalpha\_int(szk\_(g)+1) \textcolor{comment}{! The change in specific volume across each} 127 \textcolor{comment}{! interface [R-1 ~> m3 kg-1].} 128 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The computational domain for the equation of state} 129 \textcolor{keywordtype}{integer} :: is, ie, js, je, isq, ieq, jsq, jeq, nz, nkmb 130 \textcolor{keywordtype}{integer} :: i, j, k 131 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke 132 nkmb=gv%nk\_rho\_varies 133 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB 134 eosdom(1) = isq - (g%isd-1) ; eosdom(2) = g%iec+1 - (g%isd-1) 135 136 use\_p\_atm = .false. 137 \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} 138 is\_split = .false. ; \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(pbce)) is\_split = .true. 139 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 140 141 \textcolor{keywordflow}{if} (.not.\textcolor{keyword}{associated}(cs)) \textcolor{keyword}{call }mom\_error(fatal, & 142 \textcolor{stringliteral}{"MOM\_PressureForce\_Mont: Module must be initialized before it is used."}) 143 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 144 \textcolor{keywordflow}{if} (query\_compressible(tv%eqn\_of\_state)) \textcolor{keyword}{call }mom\_error(fatal, & 145 \textcolor{stringliteral}{"PressureForce\_Mont\_nonBouss: The Montgomery form of the pressure force "}//& 146 \textcolor{stringliteral}{"can no longer be used with a compressible EOS. Use #define ANALYTIC\_FV\_PGF."}) 147 \textcolor{keywordflow}{ endif} 148 149 i\_gearth = 1.0 / gv%g\_Earth 150 dp\_neglect = gv%g\_Earth * gv%H\_to\_RZ * gv%H\_subroundoff 151 \textcolor{keywordflow}{do} k=1,nz ; alpha\_lay(k) = 1.0 / (gv%Rlay(k)) ;\textcolor{keywordflow}{ enddo} 152 \textcolor{keywordflow}{do} k=2,nz ; dalpha\_int(k) = alpha\_lay(k-1) - alpha\_lay(k) ;\textcolor{keywordflow}{ enddo} 153 154 \textcolor{keywordflow}{if} (use\_p\_atm) \textcolor{keywordflow}{then} 155 \textcolor{comment}{!$OMP parallel do default(shared)} 156 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 ; p(i,j,1) = p\_atm(i,j) ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 157 \textcolor{keywordflow}{else} 158 \textcolor{comment}{!$OMP parallel do default(shared)} 159 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 ; p(i,j,1) = 0.0 ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 160 \textcolor{keywordflow}{ endif} 161 \textcolor{comment}{!$OMP parallel do default(shared)} 162 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 163 p(i,j,k+1) = p(i,j,k) + gv%g\_Earth * gv%H\_to\_RZ * h(i,j,k) 164 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 165 166 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(eta)) \textcolor{keywordflow}{then} 167 pa\_to\_h = 1.0 / (gv%g\_Earth * gv%H\_to\_RZ) 168 \textcolor{keywordflow}{if} (use\_p\_atm) \textcolor{keywordflow}{then} 169 \textcolor{comment}{!$OMP parallel do default(shared)} 170 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 171 eta(i,j) = (p(i,j,nz+1) - p\_atm(i,j)) * pa\_to\_h \textcolor{comment}{! eta has the same units as h.} 172 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 173 \textcolor{keywordflow}{else} 174 \textcolor{comment}{!$OMP parallel do default(shared)} 175 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 176 eta(i,j) = p(i,j,nz+1) * pa\_to\_h \textcolor{comment}{! eta has the same units as h.} 177 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 178 \textcolor{keywordflow}{ endif} 179 \textcolor{keywordflow}{ endif} 180 181 \textcolor{keywordflow}{if} (cs%tides) \textcolor{keywordflow}{then} 182 \textcolor{comment}{! Determine the sea surface height anomalies, to enable the calculation} 183 \textcolor{comment}{! of self-attraction and loading.} 184 \textcolor{comment}{!$OMP parallel do default(shared)} 185 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 186 ssh(i,j) = -g%bathyT(i,j) 187 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 188 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 189 \textcolor{comment}{!$OMP parallel do default(shared)} 190 \textcolor{keywordflow}{do} k=1,nz 191 \textcolor{keyword}{call }int\_specific\_vol\_dp(tv%T(:,:,k), tv%S(:,:,k), p(:,:,k), p(:,:,k+1), & 192 0.0, g%HI, tv%eqn\_of\_state, us, dz\_geo(:,:,k), halo\_size=1) 193 \textcolor{keywordflow}{ enddo} 194 \textcolor{comment}{!$OMP parallel do default(shared)} 195 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 196 ssh(i,j) = ssh(i,j) + i\_gearth * dz\_geo(i,j,k) 197 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 198 \textcolor{keywordflow}{else} 199 \textcolor{comment}{!$OMP parallel do default(shared)} 200 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 201 ssh(i,j) = ssh(i,j) + gv%H\_to\_RZ * h(i,j,k) * alpha\_lay(k) 202 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 203 \textcolor{keywordflow}{ endif} 204 205 \textcolor{keyword}{call }calc\_tidal\_forcing(cs%Time, ssh, e\_tidal, g, cs%tides\_CSp, m\_to\_z=us%m\_to\_Z) 206 \textcolor{comment}{!$OMP parallel do default(shared)} 207 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 208 geopot\_bot(i,j) = -gv%g\_Earth*(e\_tidal(i,j) + g%bathyT(i,j)) 209 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 210 \textcolor{keywordflow}{else} 211 \textcolor{comment}{!$OMP parallel do default(shared)} 212 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 213 geopot\_bot(i,j) = -gv%g\_Earth*g%bathyT(i,j) 214 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 215 \textcolor{keywordflow}{ endif} 216 217 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 218 \textcolor{comment}{! Calculate in-situ specific volumes (alpha\_star).} 219 220 \textcolor{comment}{! With a bulk mixed layer, replace the T & S of any layers that are} 221 \textcolor{comment}{! lighter than the the buffer layer with the properties of the buffer} 222 \textcolor{comment}{! layer. These layers will be massless anyway, and it avoids any} 223 \textcolor{comment}{! formal calculations with hydrostatically unstable profiles.} 224 \textcolor{keywordflow}{if} (nkmb>0) \textcolor{keywordflow}{then} 225 tv\_tmp%T => t\_tmp ; tv\_tmp%S => s\_tmp 226 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 227 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p\_ref(i) = tv%P\_Ref ;\textcolor{keywordflow}{ enddo} 228 \textcolor{comment}{!$OMP parallel do default(shared) private(Rho\_cv\_BL)} 229 \textcolor{keywordflow}{do} j=jsq,jeq+1 230 \textcolor{keywordflow}{do} k=1,nkmb ; \textcolor{keywordflow}{do} i=isq,ieq+1 231 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 232 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 233 \textcolor{keyword}{call }calculate\_density(tv%T(:,j,nkmb), tv%S(:,j,nkmb), p\_ref, rho\_cv\_bl(:), & 234 tv%eqn\_of\_state, eosdom) 235 \textcolor{keywordflow}{do} k=nkmb+1,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 236 \textcolor{keywordflow}{if} (gv%Rlay(k) < rho\_cv\_bl(i)) \textcolor{keywordflow}{then} 237 tv\_tmp%T(i,j,k) = tv%T(i,j,nkmb) ; tv\_tmp%S(i,j,k) = tv%S(i,j,nkmb) 238 \textcolor{keywordflow}{else} 239 tv\_tmp%T(i,j,k) = tv%T(i,j,k) ; tv\_tmp%S(i,j,k) = tv%S(i,j,k) 240 \textcolor{keywordflow}{ endif} 241 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 242 \textcolor{keywordflow}{ enddo} 243 \textcolor{keywordflow}{else} 244 tv\_tmp%T => tv%T ; tv\_tmp%S => tv%S 245 tv\_tmp%eqn\_of\_state => tv%eqn\_of\_state 246 \textcolor{keywordflow}{do} i=isq,ieq+1 ; p\_ref(i) = 0 ;\textcolor{keywordflow}{ enddo} 247 \textcolor{keywordflow}{ endif} 248 \textcolor{comment}{!$OMP parallel do default(shared) private(rho\_in\_situ)} 249 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} j=jsq,jeq+1 250 \textcolor{keyword}{call }calculate\_density(tv\_tmp%T(:,j,k), tv\_tmp%S(:,j,k), p\_ref, rho\_in\_situ, & 251 tv%eqn\_of\_state, eosdom) 252 \textcolor{keywordflow}{do} i=isq,ieq+1 ; alpha\_star(i,j,k) = 1.0 / rho\_in\_situ(i) ;\textcolor{keywordflow}{ enddo} 253 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 254 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 255 256 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 257 \textcolor{comment}{!$OMP parallel do default(shared)} 258 \textcolor{keywordflow}{do} j=jsq,jeq+1 259 \textcolor{keywordflow}{do} i=isq,ieq+1 260 m(i,j,nz) = geopot\_bot(i,j) + p(i,j,nz+1) * alpha\_star(i,j,nz) 261 \textcolor{keywordflow}{ enddo} 262 \textcolor{keywordflow}{do} k=nz-1,1,-1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 263 m(i,j,k) = m(i,j,k+1) + p(i,j,k+1) * (alpha\_star(i,j,k) - alpha\_star(i,j,k+1)) 264 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 265 \textcolor{keywordflow}{ enddo} 266 \textcolor{keywordflow}{else} \textcolor{comment}{! not use\_EOS} 267 \textcolor{comment}{!$OMP parallel do default(shared)} 268 \textcolor{keywordflow}{do} j=jsq,jeq+1 269 \textcolor{keywordflow}{do} i=isq,ieq+1 270 m(i,j,nz) = geopot\_bot(i,j) + p(i,j,nz+1) * alpha\_lay(nz) 271 \textcolor{keywordflow}{ enddo} 272 \textcolor{keywordflow}{do} k=nz-1,1,-1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 273 m(i,j,k) = m(i,j,k+1) + p(i,j,k+1) * dalpha\_int(k+1) 274 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 275 \textcolor{keywordflow}{ enddo} 276 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 277 278 \textcolor{keywordflow}{if} (cs%GFS\_scale < 1.0) \textcolor{keywordflow}{then} 279 \textcolor{comment}{! Adjust the Montgomery potential to make this a reduced gravity model.} 280 \textcolor{comment}{!$OMP parallel do default(shared)} 281 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 282 dm(i,j) = (cs%GFS\_scale - 1.0) * m(i,j,1) 283 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 284 \textcolor{comment}{!$OMP parallel do default(shared)} 285 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 286 m(i,j,k) = m(i,j,k) + dm(i,j) 287 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 288 289 \textcolor{comment}{! Could instead do the following, to avoid taking small differences} 290 \textcolor{comment}{! of large numbers...} 291 \textcolor{comment}{! do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1} 292 \textcolor{comment}{! M(i,j,1) = CS%GFS\_scale * M(i,j,1)} 293 \textcolor{comment}{! enddo ; enddo} 294 \textcolor{comment}{! if (use\_EOS) then} 295 \textcolor{comment}{! do k=2,nz ; do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1} 296 \textcolor{comment}{! M(i,j,k) = M(i,j,k-1) - p(i,j,K) * (alpha\_star(i,j,k-1) - alpha\_star(i,j,k))} 297 \textcolor{comment}{! enddo ; enddo ; enddo} 298 \textcolor{comment}{! else ! not use\_EOS} 299 \textcolor{comment}{! do k=2,nz ; do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1} 300 \textcolor{comment}{! M(i,j,k) = M(i,j,k-1) - p(i,j,K) * dalpha\_int(K)} 301 \textcolor{comment}{! enddo ; enddo ; enddo} 302 \textcolor{comment}{! endif ! use\_EOS} 303 304 \textcolor{keywordflow}{ endif} 305 306 \textcolor{comment}{! Note that ddM/dPb = alpha\_star(i,j,1)} 307 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(pbce)) \textcolor{keywordflow}{then} 308 \textcolor{keyword}{call }set\_pbce\_nonbouss(p, tv\_tmp, g, gv, us, cs%GFS\_scale, pbce, alpha\_star) 309 \textcolor{keywordflow}{ endif} 310 311 \textcolor{comment}{! Calculate the pressure force. On a Cartesian grid,} 312 \textcolor{comment}{! PFu = - dM/dx and PFv = - dM/dy.} 313 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 314 \textcolor{comment}{!$OMP parallel do default(shared) private(dp\_star,PFu\_bc,PFv\_bc)} 315 \textcolor{keywordflow}{do} k=1,nz 316 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 317 dp\_star(i,j) = (p(i,j,k+1) - p(i,j,k)) + dp\_neglect 318 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 319 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 320 \textcolor{comment}{! PFu\_bc = p* grad alpha*} 321 pfu\_bc = (alpha\_star(i+1,j,k) - alpha\_star(i,j,k)) * (g%IdxCu(i,j) * & 322 ((dp\_star(i,j)*dp\_star(i+1,j) + (p(i,j,k)*dp\_star(i+1,j) + p(i+1,j,k)*dp\_star(i,j))) / & 323 (dp\_star(i,j) + dp\_star(i+1,j)))) 324 pfu(i,j,k) = -(m(i+1,j,k) - m(i,j,k)) * g%IdxCu(i,j) + pfu\_bc 325 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs%PFu\_bc)) cs%PFu\_bc(i,j,k) = pfu\_bc 326 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 327 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 328 pfv\_bc = (alpha\_star(i,j+1,k) - alpha\_star(i,j,k)) * (g%IdyCv(i,j) * & 329 ((dp\_star(i,j)*dp\_star(i,j+1) + (p(i,j,k)*dp\_star(i,j+1) + p(i,j+1,k)*dp\_star(i,j))) / & 330 (dp\_star(i,j) + dp\_star(i,j+1)))) 331 pfv(i,j,k) = -(m(i,j+1,k) - m(i,j,k)) * g%IdyCv(i,j) + pfv\_bc 332 \textcolor{keywordflow}{if} (\textcolor{keyword}{associated}(cs%PFv\_bc)) cs%PFv\_bc(i,j,k) = pfv\_bc 333 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 334 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! k-loop} 335 \textcolor{keywordflow}{else} \textcolor{comment}{! .not. use\_EOS} 336 \textcolor{comment}{!$OMP parallel do default(shared)} 337 \textcolor{keywordflow}{do} k=1,nz 338 \textcolor{keywordflow}{do} j=js,je ; \textcolor{keywordflow}{do} i=isq,ieq 339 pfu(i,j,k) = -(m(i+1,j,k) - m(i,j,k)) * g%IdxCu(i,j) 340 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 341 \textcolor{keywordflow}{do} j=jsq,jeq ; \textcolor{keywordflow}{do} i=is,ie 342 pfv(i,j,k) = -(m(i,j+1,k) - m(i,j,k)) * g%IdyCv(i,j) 343 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 344 \textcolor{keywordflow}{ enddo} 345 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 346 347 \textcolor{keywordflow}{if} (cs%id\_PFu\_bc>0) \textcolor{keyword}{call }post\_data(cs%id\_PFu\_bc, cs%PFu\_bc, cs%diag) 348 \textcolor{keywordflow}{if} (cs%id\_PFv\_bc>0) \textcolor{keyword}{call }post\_data(cs%id\_PFv\_bc, cs%PFv\_bc, cs%diag) 349 \textcolor{keywordflow}{if} (cs%id\_e\_tidal>0) \textcolor{keyword}{call }post\_data(cs%id\_e\_tidal, e\_tidal, cs%diag) 350 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__mont_a0779efd30a447c2bc20294c46aeea180}\label{namespacemom__pressureforce__mont_a0779efd30a447c2bc20294c46aeea180}} \index{mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}!set\+\_\+pbce\+\_\+bouss@{set\+\_\+pbce\+\_\+bouss}} \index{set\+\_\+pbce\+\_\+bouss@{set\+\_\+pbce\+\_\+bouss}!mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}} \subsubsection{\texorpdfstring{set\+\_\+pbce\+\_\+bouss()}{set\_pbce\_bouss()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+mont\+::set\+\_\+pbce\+\_\+bouss (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)+1), intent(in)}]{e, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\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, intent(in)}]{Rho0, }\item[{real, intent(in)}]{G\+F\+S\+\_\+scale, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(out)}]{pbce, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in), optional}]{rho\+\_\+star }\end{DoxyParamCaption})} Determines the partial derivative of the acceleration due to pressure forces with the free surface height. \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 e} & Interface height \mbox{[}Z $\sim$$>$ m\mbox{]}.\\ \hline \mbox{\tt in} & {\em tv} & Thermodynamic variables\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em rho0} & The \char`\"{}\+Boussinesq\char`\"{} ocean density \mbox{[}R $\sim$$>$ kg m-\/3\mbox{]}.\\ \hline \mbox{\tt in} & {\em gfs\+\_\+scale} & Ratio between gravity applied to top interface and the gravitational acceleration of the planet \mbox{[}nondim\mbox{]}. Usually this ratio is 1.\\ \hline \mbox{\tt out} & {\em pbce} & The baroclinic pressure anomaly in each layer due\\ \hline \mbox{\tt in} & {\em rho\+\_\+star} & The layer densities (maybe compressibility \\ \hline \end{DoxyParams} Definition at line 607 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+Montgomery.\+F90. \begin{DoxyCode} 607 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< Ocean grid structure} 608 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< Vertical grid structure} 609 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)}, \textcolor{keywordtype}{intent(in)} :: e\textcolor{comment}{ !< Interface height [Z ~> m].} 610 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< Thermodynamic variables} 611 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 612 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: rho0\textcolor{comment}{ !< The "Boussinesq" ocean density [R ~> kg m-3].} 613 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: gfs\_scale\textcolor{comment}{ !< Ratio between gravity applied to top} 614 \textcolor{comment}{ !! interface and the gravitational acceleration of} 615 \textcolor{comment}{ !! the planet [nondim]. Usually this ratio is 1.} 616 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, & 617 \textcolor{keywordtype}{intent(out)} :: pbce\textcolor{comment}{ !< The baroclinic pressure anomaly in each layer due} 618 \textcolor{comment}{ !! to free surface height anomalies} 619 \textcolor{comment}{ !! [L2 T-2 H-1 ~> m s-2].} 620 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, & 621 \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: rho\_star\textcolor{comment}{ !< The layer densities (maybe compressibility} 622 \textcolor{comment}{ !! compensated), times g/rho\_0 [L2 Z-1 T-2 ~> m s-2].} 623 624 \textcolor{comment}{! Local variables} 625 \textcolor{keywordtype}{real} :: ihtot(szi\_(g)) \textcolor{comment}{! The inverse of the sum of the layer thicknesses [H-1 ~> m-1 or m2 kg-1].} 626 \textcolor{keywordtype}{real} :: press(szi\_(g)) \textcolor{comment}{! Interface pressure [R L2 T-2 ~> Pa].} 627 \textcolor{keywordtype}{real} :: t\_int(szi\_(g)) \textcolor{comment}{! Interface temperature [degC].} 628 \textcolor{keywordtype}{real} :: s\_int(szi\_(g)) \textcolor{comment}{! Interface salinity [ppt].} 629 \textcolor{keywordtype}{real} :: dr\_dt(szi\_(g)) \textcolor{comment}{! Partial derivative of density with temperature [R degC-1 ~> kg m-3 degC-1].} 630 \textcolor{keywordtype}{real} :: dr\_ds(szi\_(g)) \textcolor{comment}{! Partial derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1].} 631 \textcolor{keywordtype}{real} :: rho\_in\_situ(szi\_(g)) \textcolor{comment}{! In-situ density at the top of a layer [R ~> kg m-3].} 632 \textcolor{keywordtype}{real} :: g\_rho0 \textcolor{comment}{! A scaled version of g\_Earth / Rho0 [L2 Z-1 T-2 R-1 ~> m4 s-2 kg-1]} 633 \textcolor{keywordtype}{real} :: rho0xg \textcolor{comment}{! g\_Earth * Rho0 [kg s-2 m-1 Z-1 ~> kg s-2 m-2]} 634 \textcolor{keywordtype}{logical} :: use\_eos \textcolor{comment}{! If true, density is calculated from T & S using} 635 \textcolor{comment}{! an equation of state.} 636 \textcolor{keywordtype}{real} :: z\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 637 \textcolor{comment}{! in roundoff and can be neglected [Z ~> m].} 638 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The computational domain for the equation of state} 639 \textcolor{keywordtype}{integer} :: isq, ieq, jsq, jeq, nz, i, j, k 640 641 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB ; nz = g%ke 642 eosdom(1) = isq - (g%isd-1) ; eosdom(2) = g%iec+1 - (g%isd-1) 643 644 rho0xg = rho0 * gv%g\_Earth 645 g\_rho0 = gv%g\_Earth / gv%Rho0 646 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 647 z\_neglect = gv%H\_subroundoff*gv%H\_to\_Z 648 649 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 650 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(rho\_star)) \textcolor{keywordflow}{then} 651 \textcolor{comment}{!$OMP parallel do default(shared) private(Ihtot)} 652 \textcolor{keywordflow}{do} j=jsq,jeq+1 653 \textcolor{keywordflow}{do} i=isq,ieq+1 654 ihtot(i) = gv%H\_to\_Z / ((e(i,j,1)-e(i,j,nz+1)) + z\_neglect) 655 pbce(i,j,1) = gfs\_scale * rho\_star(i,j,1) * gv%H\_to\_Z 656 \textcolor{keywordflow}{ enddo} 657 \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 658 pbce(i,j,k) = pbce(i,j,k-1) + (rho\_star(i,j,k)-rho\_star(i,j,k-1)) * & 659 ((e(i,j,k) - e(i,j,nz+1)) * ihtot(i)) 660 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 661 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of j loop} 662 \textcolor{keywordflow}{else} 663 \textcolor{comment}{!$OMP parallel do default(shared) private(Ihtot,press,rho\_in\_situ,T\_int,S\_int,dR\_dT,dR\_dS)} 664 \textcolor{keywordflow}{do} j=jsq,jeq+1 665 \textcolor{keywordflow}{do} i=isq,ieq+1 666 ihtot(i) = gv%H\_to\_Z / ((e(i,j,1)-e(i,j,nz+1)) + z\_neglect) 667 press(i) = -rho0xg*e(i,j,1) 668 \textcolor{keywordflow}{ enddo} 669 \textcolor{keyword}{call }calculate\_density(tv%T(:,j,1), tv%S(:,j,1), press, rho\_in\_situ, & 670 tv%eqn\_of\_state, eosdom) 671 \textcolor{keywordflow}{do} i=isq,ieq+1 672 pbce(i,j,1) = g\_rho0*(gfs\_scale * rho\_in\_situ(i)) * gv%H\_to\_Z 673 \textcolor{keywordflow}{ enddo} 674 \textcolor{keywordflow}{do} k=2,nz 675 \textcolor{keywordflow}{do} i=isq,ieq+1 676 press(i) = -rho0xg*e(i,j,k) 677 t\_int(i) = 0.5*(tv%T(i,j,k-1)+tv%T(i,j,k)) 678 s\_int(i) = 0.5*(tv%S(i,j,k-1)+tv%S(i,j,k)) 679 \textcolor{keywordflow}{ enddo} 680 \textcolor{keyword}{call }calculate\_density\_derivs(t\_int, s\_int, press, dr\_dt, dr\_ds, & 681 tv%eqn\_of\_state, eosdom) 682 \textcolor{keywordflow}{do} i=isq,ieq+1 683 pbce(i,j,k) = pbce(i,j,k-1) + g\_rho0 * & 684 ((e(i,j,k) - e(i,j,nz+1)) * ihtot(i)) * & 685 (dr\_dt(i)*(tv%T(i,j,k)-tv%T(i,j,k-1)) + & 686 dr\_ds(i)*(tv%S(i,j,k)-tv%S(i,j,k-1))) 687 \textcolor{keywordflow}{ enddo} 688 \textcolor{keywordflow}{ enddo} 689 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of j loop} 690 \textcolor{keywordflow}{ endif} 691 \textcolor{keywordflow}{else} \textcolor{comment}{! not use\_EOS} 692 \textcolor{comment}{!$OMP parallel do default(shared) private(Ihtot)} 693 \textcolor{keywordflow}{do} j=jsq,jeq+1 694 \textcolor{keywordflow}{do} i=isq,ieq+1 695 ihtot(i) = 1.0 / ((e(i,j,1)-e(i,j,nz+1)) + z\_neglect) 696 pbce(i,j,1) = gv%g\_prime(1) * gv%H\_to\_Z 697 \textcolor{keywordflow}{ enddo} 698 \textcolor{keywordflow}{do} k=2,nz ; \textcolor{keywordflow}{do} i=isq,ieq+1 699 pbce(i,j,k) = pbce(i,j,k-1) + & 700 (gv%g\_prime(k)*gv%H\_to\_Z) * ((e(i,j,k) - e(i,j,nz+1)) * ihtot(i)) 701 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 702 \textcolor{keywordflow}{ enddo} \textcolor{comment}{! end of j loop} 703 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 704 \end{DoxyCode} \mbox{\Hypertarget{namespacemom__pressureforce__mont_a297cdf6e4eb83d250444c5f527a8a232}\label{namespacemom__pressureforce__mont_a297cdf6e4eb83d250444c5f527a8a232}} \index{mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}!set\+\_\+pbce\+\_\+nonbouss@{set\+\_\+pbce\+\_\+nonbouss}} \index{set\+\_\+pbce\+\_\+nonbouss@{set\+\_\+pbce\+\_\+nonbouss}!mom\+\_\+pressureforce\+\_\+mont@{mom\+\_\+pressureforce\+\_\+mont}} \subsubsection{\texorpdfstring{set\+\_\+pbce\+\_\+nonbouss()}{set\_pbce\_nonbouss()}} {\footnotesize\ttfamily subroutine, public mom\+\_\+pressureforce\+\_\+mont\+::set\+\_\+pbce\+\_\+nonbouss (\begin{DoxyParamCaption}\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)+1), intent(in)}]{p, }\item[{type(thermo\+\_\+var\+\_\+ptrs), intent(in)}]{tv, }\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, intent(in)}]{G\+F\+S\+\_\+scale, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(out)}]{pbce, }\item[{real, dimension(szi\+\_\+(g),szj\+\_\+(g),szk\+\_\+(g)), intent(in), optional}]{alpha\+\_\+star }\end{DoxyParamCaption})} Determines the partial derivative of the acceleration due to pressure forces with the column mass. \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 p} & Interface pressures \mbox{[}R L2 T-\/2 $\sim$$>$ Pa\mbox{]}.\\ \hline \mbox{\tt in} & {\em tv} & Thermodynamic variables\\ \hline \mbox{\tt in} & {\em us} & A dimensional unit scaling type\\ \hline \mbox{\tt in} & {\em gfs\+\_\+scale} & Ratio between gravity applied to top interface and the gravitational acceleration of the planet \mbox{[}nondim\mbox{]}. Usually this ratio is 1.\\ \hline \mbox{\tt out} & {\em pbce} & The baroclinic pressure anomaly in each layer due to free surface height anomalies \mbox{[}L2 H-\/1 T-\/2 $\sim$$>$ m4 kg-\/1 s-\/2\mbox{]}.\\ \hline \mbox{\tt in} & {\em alpha\+\_\+star} & The layer specific volumes (maybe compressibility compensated) \mbox{[}R-\/1 $\sim$$>$ m3 kg-\/1\mbox{]}. \\ \hline \end{DoxyParams} Definition at line 710 of file M\+O\+M\+\_\+\+Pressure\+Force\+\_\+\+Montgomery.\+F90. \begin{DoxyCode} 710 \textcolor{keywordtype}{type}(ocean\_grid\_type), \textcolor{keywordtype}{intent(in)} :: g\textcolor{comment}{ !< Ocean grid structure} 711 \textcolor{keywordtype}{type}(verticalgrid\_type), \textcolor{keywordtype}{intent(in)} :: gv\textcolor{comment}{ !< Vertical grid structure} 712 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G)+1)}, \textcolor{keywordtype}{intent(in)} :: p\textcolor{comment}{ !< Interface pressures [R L2 T-2 ~> Pa].} 713 \textcolor{keywordtype}{type}(thermo\_var\_ptrs), \textcolor{keywordtype}{intent(in)} :: tv\textcolor{comment}{ !< Thermodynamic variables} 714 \textcolor{keywordtype}{type}(unit\_scale\_type), \textcolor{keywordtype}{intent(in)} :: us\textcolor{comment}{ !< A dimensional unit scaling type} 715 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{intent(in)} :: gfs\_scale\textcolor{comment}{ !< Ratio between gravity applied to top} 716 \textcolor{comment}{ !! interface and the gravitational acceleration of} 717 \textcolor{comment}{ !! the planet [nondim]. Usually this ratio is 1.} 718 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{intent(out)} :: pbce\textcolor{comment}{ !< The baroclinic pressure anomaly in each} 719 \textcolor{comment}{ !! layer due to free surface height anomalies} 720 \textcolor{comment}{ !! [L2 H-1 T-2 ~> m4 kg-1 s-2].} 721 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G),SZK\_(G))}, \textcolor{keywordtype}{optional}, \textcolor{keywordtype}{intent(in)} :: alpha\_star\textcolor{comment}{ !< The layer specific volumes} 722 \textcolor{comment}{ !! (maybe compressibility compensated) [R-1 ~> m3 kg-1].} 723 \textcolor{comment}{! Local variables} 724 \textcolor{keywordtype}{real}, \textcolor{keywordtype}{dimension(SZI\_(G),SZJ\_(G))} :: & 725 dpbce, & \textcolor{comment}{! A barotropic correction to the pbce to enable the use of} 726 \textcolor{comment}{! a reduced gravity form of the equations [L2 H-1 T-2 ~> m4 kg-1 s-2].} 727 c\_htot \textcolor{comment}{! dP\_dH divided by the total ocean pressure [H-1 ~> m2 kg-1].} 728 \textcolor{keywordtype}{real} :: t\_int(szi\_(g)) \textcolor{comment}{! Interface temperature [degC].} 729 \textcolor{keywordtype}{real} :: s\_int(szi\_(g)) \textcolor{comment}{! Interface salinity [ppt].} 730 \textcolor{keywordtype}{real} :: dr\_dt(szi\_(g)) \textcolor{comment}{! Partial derivative of density with temperature [R degC-1 ~> kg m-3 degC-1].} 731 \textcolor{keywordtype}{real} :: dr\_ds(szi\_(g)) \textcolor{comment}{! Partial derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1].} 732 \textcolor{keywordtype}{real} :: rho\_in\_situ(szi\_(g)) \textcolor{comment}{! In-situ density at an interface [R ~> kg m-3].} 733 \textcolor{keywordtype}{real} :: alpha\_lay(szk\_(g)) \textcolor{comment}{! The specific volume of each layer [R-1 ~> m3 kg-1].} 734 \textcolor{keywordtype}{real} :: dalpha\_int(szk\_(g)+1) \textcolor{comment}{! The change in specific volume across each interface [R-1 ~> m3 kg-1].} 735 \textcolor{keywordtype}{real} :: dp\_dh \textcolor{comment}{! A factor that converts from thickness to pressure times other dimensional} 736 \textcolor{comment}{! conversion factors [R L2 T-2 H-1 ~> Pa m2 kg-1].} 737 \textcolor{keywordtype}{real} :: dp\_neglect \textcolor{comment}{! A thickness that is so small it is usually lost} 738 \textcolor{comment}{! in roundoff and can be neglected [R L2 T-2 ~> Pa].} 739 \textcolor{keywordtype}{logical} :: use\_eos \textcolor{comment}{! If true, density is calculated from T & S using an equation of state.} 740 \textcolor{keywordtype}{integer}, \textcolor{keywordtype}{dimension(2)} :: eosdom \textcolor{comment}{! The computational domain for the equation of state} 741 \textcolor{keywordtype}{integer} :: isq, ieq, jsq, jeq, nz, i, j, k 742 743 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB ; nz = g%ke 744 eosdom(1) = isq - (g%isd-1) ; eosdom(2) = g%iec+1 - (g%isd-1) 745 746 use\_eos = \textcolor{keyword}{associated}(tv%eqn\_of\_state) 747 748 dp\_dh = gv%g\_Earth * gv%H\_to\_RZ 749 dp\_neglect = gv%g\_Earth * gv%H\_to\_RZ * gv%H\_subroundoff 750 751 \textcolor{keywordflow}{if} (use\_eos) \textcolor{keywordflow}{then} 752 \textcolor{keywordflow}{if} (\textcolor{keyword}{present}(alpha\_star)) \textcolor{keywordflow}{then} 753 \textcolor{comment}{!$OMP parallel do default(shared)} 754 \textcolor{keywordflow}{do} j=jsq,jeq+1 755 \textcolor{keywordflow}{do} i=isq,ieq+1 756 c\_htot(i,j) = dp\_dh / ((p(i,j,nz+1)-p(i,j,1)) + dp\_neglect) 757 pbce(i,j,nz) = dp\_dh * alpha\_star(i,j,nz) 758 \textcolor{keywordflow}{ enddo} 759 \textcolor{keywordflow}{do} k=nz-1,1,-1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 760 pbce(i,j,k) = pbce(i,j,k+1) + ((p(i,j,k+1)-p(i,j,1)) * c\_htot(i,j)) * & 761 (alpha\_star(i,j,k) - alpha\_star(i,j,k+1)) 762 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 763 \textcolor{keywordflow}{ enddo} 764 \textcolor{keywordflow}{else} 765 \textcolor{comment}{!$OMP parallel do default(shared) private(T\_int,S\_int,dR\_dT,dR\_dS,rho\_in\_situ)} 766 \textcolor{keywordflow}{do} j=jsq,jeq+1 767 \textcolor{keyword}{call }calculate\_density(tv%T(:,j,nz), tv%S(:,j,nz), p(:,j,nz+1), rho\_in\_situ, & 768 tv%eqn\_of\_state, eosdom) 769 \textcolor{keywordflow}{do} i=isq,ieq+1 770 c\_htot(i,j) = dp\_dh / ((p(i,j,nz+1)-p(i,j,1)) + dp\_neglect) 771 pbce(i,j,nz) = dp\_dh / (rho\_in\_situ(i)) 772 \textcolor{keywordflow}{ enddo} 773 \textcolor{keywordflow}{do} k=nz-1,1,-1 774 \textcolor{keywordflow}{do} i=isq,ieq+1 775 t\_int(i) = 0.5*(tv%T(i,j,k)+tv%T(i,j,k+1)) 776 s\_int(i) = 0.5*(tv%S(i,j,k)+tv%S(i,j,k+1)) 777 \textcolor{keywordflow}{ enddo} 778 \textcolor{keyword}{call }calculate\_density(t\_int, s\_int, p(:,j,k+1), rho\_in\_situ, tv%eqn\_of\_state, eosdom) 779 \textcolor{keyword}{call }calculate\_density\_derivs(t\_int, s\_int, p(:,j,k+1), dr\_dt, dr\_ds, & 780 tv%eqn\_of\_state, eosdom) 781 \textcolor{keywordflow}{do} i=isq,ieq+1 782 pbce(i,j,k) = pbce(i,j,k+1) + ((p(i,j,k+1)-p(i,j,1))*c\_htot(i,j)) * & 783 ((dr\_dt(i)*(tv%T(i,j,k+1)-tv%T(i,j,k)) + & 784 dr\_ds(i)*(tv%S(i,j,k+1)-tv%S(i,j,k))) / (rho\_in\_situ(i)**2)) 785 \textcolor{keywordflow}{ enddo} 786 \textcolor{keywordflow}{ enddo} 787 \textcolor{keywordflow}{ enddo} 788 \textcolor{keywordflow}{ endif} 789 \textcolor{keywordflow}{else} \textcolor{comment}{! not use\_EOS} 790 791 \textcolor{keywordflow}{do} k=1,nz ; alpha\_lay(k) = 1.0 / (gv%Rlay(k)) ;\textcolor{keywordflow}{ enddo} 792 \textcolor{keywordflow}{do} k=2,nz ; dalpha\_int(k) = alpha\_lay(k-1) - alpha\_lay(k) ;\textcolor{keywordflow}{ enddo} 793 794 \textcolor{comment}{!$OMP parallel do default(shared)} 795 \textcolor{keywordflow}{do} j=jsq,jeq+1 796 \textcolor{keywordflow}{do} i=isq,ieq+1 797 c\_htot(i,j) = dp\_dh / ((p(i,j,nz+1)-p(i,j,1)) + dp\_neglect) 798 pbce(i,j,nz) = dp\_dh * alpha\_lay(nz) 799 \textcolor{keywordflow}{ enddo} 800 \textcolor{keywordflow}{do} k=nz-1,1,-1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 801 pbce(i,j,k) = pbce(i,j,k+1) + ((p(i,j,k+1)-p(i,j,1))*c\_htot(i,j)) * dalpha\_int(k+1) 802 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 803 \textcolor{keywordflow}{ enddo} 804 \textcolor{keywordflow}{ endif} \textcolor{comment}{! use\_EOS} 805 806 \textcolor{keywordflow}{if} (gfs\_scale < 1.0) \textcolor{keywordflow}{then} 807 \textcolor{comment}{! Adjust the Montgomery potential to make this a reduced gravity model.} 808 \textcolor{comment}{!$OMP parallel do default(shared)} 809 \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 810 dpbce(i,j) = (gfs\_scale - 1.0) * pbce(i,j,1) 811 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 812 \textcolor{comment}{!$OMP parallel do default(shared)} 813 \textcolor{keywordflow}{do} k=1,nz ; \textcolor{keywordflow}{do} j=jsq,jeq+1 ; \textcolor{keywordflow}{do} i=isq,ieq+1 814 pbce(i,j,k) = pbce(i,j,k) + dpbce(i,j) 815 \textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} ;\textcolor{keywordflow}{ enddo} 816 \textcolor{keywordflow}{ endif} 817 \end{DoxyCode}