MOM6
mom_pressureforce Module Reference

Detailed Description

A thin wrapper for Boussinesq/non-Boussinesq forms of the pressure force calculation.

This thin module provides a branch to two forms of the horizontal accelerations due to pressure gradients. The two options currently available are a Montgomery potential form (used in traditional isopycnal layer models), and the analytic finite volume form.

Data Types

type  pressureforce_cs
 Pressure force control structure. More...
 

Functions/Subroutines

subroutine, public pressureforce (h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm, pbce, eta)
 A thin layer between the model and the Boussinesq and non-Boussinesq pressure force routines. More...
 
subroutine, public pressureforce_init (Time, G, GV, US, param_file, diag, CS, tides_CSp)
 Initialize the pressure force control structure. More...
 
subroutine, public pressureforce_end (CS)
 Deallocate the pressure force control structure. More...
 

Function/Subroutine Documentation

◆ pressureforce()

subroutine, public mom_pressureforce::pressureforce ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(thermo_var_ptrs), intent(in)  tv,
real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke), intent(out)  PFu,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke), intent(out)  PFv,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(pressureforce_cs), pointer  CS,
type(ale_cs), pointer  ALE_CSp,
real, dimension(:,:), optional, pointer  p_atm,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(out), optional  pbce,
real, dimension( g %isd: g %ied, g %jsd: g %jed), intent(out), optional  eta 
)

A thin layer between the model and the Boussinesq and non-Boussinesq pressure force routines.

Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]usA dimensional unit scaling type
[in]hLayer thicknesses [H ~> m or kg m-2]
[in]tvA structure pointing to various thermodynamic variables
[out]pfuZonal pressure force acceleration [L T-2 ~> m s-2]
[out]pfvMeridional pressure force acceleration [L T-2 ~> m s-2]
csPressure force control structure
ale_cspALE control structure
p_atmThe pressure at the ice-ocean or
[out]pbceThe baroclinic pressure anomaly in each layer
[out]etaThe bottom mass used to calculate PFu and PFv,

Definition at line 48 of file MOM_PressureForce.F90.

48  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
49  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
50  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
51  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
52  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
53  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
54  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
55  intent(out) :: PFu !< Zonal pressure force acceleration [L T-2 ~> m s-2]
56  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
57  intent(out) :: PFv !< Meridional pressure force acceleration [L T-2 ~> m s-2]
58  type(PressureForce_CS), pointer :: CS !< Pressure force control structure
59  type(ALE_CS), pointer :: ALE_CSp !< ALE control structure
60  real, dimension(:,:), &
61  optional, pointer :: p_atm !< The pressure at the ice-ocean or
62  !! atmosphere-ocean interface [Pa].
63  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
64  optional, intent(out) :: pbce !< The baroclinic pressure anomaly in each layer
65  !! due to eta anomalies [m2 s-2 H-1 ~> m s-2 or m4 s-2 kg-1].
66  real, dimension(SZI_(G),SZJ_(G)), &
67  optional, intent(out) :: eta !< The bottom mass used to calculate PFu and PFv,
68  !! [H ~> m or kg m-2], with any tidal contributions.
69 
70  if (cs%Analytic_FV_PGF .and. cs%blocked_AFV) then
71  if (gv%Boussinesq) then
72  call pressureforce_blk_afv_bouss(h, tv, pfu, pfv, g, gv, us, &
73  cs%PressureForce_blk_AFV_CSp, ale_csp, p_atm, pbce, eta)
74  else
75  call pressureforce_blk_afv_nonbouss(h, tv, pfu, pfv, g, gv, us, &
76  cs%PressureForce_blk_AFV_CSp, p_atm, pbce, eta)
77  endif
78  elseif (cs%Analytic_FV_PGF) then
79  if (gv%Boussinesq) then
80  call pressureforce_afv_bouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_AFV_CSp, &
81  ale_csp, p_atm, pbce, eta)
82  else
83  call pressureforce_afv_nonbouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_AFV_CSp, &
84  ale_csp, p_atm, pbce, eta)
85  endif
86  else
87  if (gv%Boussinesq) then
88  call pressureforce_mont_bouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_Mont_CSp, &
89  p_atm, pbce, eta)
90  else
91  call pressureforce_mont_nonbouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_Mont_CSp, &
92  p_atm, pbce, eta)
93  endif
94  endif
95 

Referenced by mom_dynamics_split_rk2::step_mom_dyn_split_rk2(), mom_dynamics_unsplit::step_mom_dyn_unsplit(), and mom_dynamics_unsplit_rk2::step_mom_dyn_unsplit_rk2().

Here is the caller graph for this function:

◆ pressureforce_end()

subroutine, public mom_pressureforce::pressureforce_end ( type(pressureforce_cs), pointer  CS)

Deallocate the pressure force control structure.

Parameters
csPressure force control structure

Definition at line 145 of file MOM_PressureForce.F90.

145  type(PressureForce_CS), pointer :: CS !< Pressure force control structure
146 
147  if (cs%Analytic_FV_PGF .and. cs%blocked_AFV) then
148  call pressureforce_blk_afv_end(cs%PressureForce_blk_AFV_CSp)
149  elseif (cs%Analytic_FV_PGF) then
150  call pressureforce_afv_end(cs%PressureForce_AFV_CSp)
151  else
152  call pressureforce_mont_end(cs%PressureForce_Mont_CSp)
153  endif
154 
155  if (associated(cs)) deallocate(cs)

◆ pressureforce_init()

subroutine, public mom_pressureforce::pressureforce_init ( type(time_type), intent(in), target  Time,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(diag_ctrl), intent(inout), target  diag,
type(pressureforce_cs), pointer  CS,
type(tidal_forcing_cs), optional, pointer  tides_CSp 
)

Initialize the pressure force control structure.

Parameters
[in]timeCurrent model time
[in]gOcean grid structure
[in]gvVertical grid structure
[in]usA dimensional unit scaling type
[in]param_fileParameter file handles
[in,out]diagDiagnostics control structure
csPressure force control structure
tides_cspTide control structure

Definition at line 100 of file MOM_PressureForce.F90.

100  type(time_type), target, intent(in) :: Time !< Current model time
101  type(ocean_grid_type), intent(in) :: G !< Ocean grid structure
102  type(verticalGrid_type), intent(in) :: GV !< Vertical grid structure
103  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
104  type(param_file_type), intent(in) :: param_file !< Parameter file handles
105  type(diag_ctrl), target, intent(inout) :: diag !< Diagnostics control structure
106  type(PressureForce_CS), pointer :: CS !< Pressure force control structure
107  type(tidal_forcing_CS), optional, pointer :: tides_CSp !< Tide control structure
108 #include "version_variable.h"
109  character(len=40) :: mdl = "MOM_PressureForce" ! This module's name.
110 
111  if (associated(cs)) then
112  call mom_error(warning, "PressureForce_init called with an associated "// &
113  "control structure.")
114  return
115  else ; allocate(cs) ; endif
116 
117  ! Read all relevant parameters and write them to the model log.
118  call log_version(param_file, mdl, version, "")
119  call get_param(param_file, mdl, "ANALYTIC_FV_PGF", cs%Analytic_FV_PGF, &
120  "If true the pressure gradient forces are calculated "//&
121  "with a finite volume form that analytically integrates "//&
122  "the equations of state in pressure to avoid any "//&
123  "possibility of numerical thermobaric instability, as "//&
124  "described in Adcroft et al., O. Mod. (2008).", default=.true.)
125  call get_param(param_file, mdl, "BLOCKED_ANALYTIC_FV_PGF", cs%blocked_AFV, &
126  "If true, used the blocked version of the ANALYTIC_FV_PGF "//&
127  "code. The value of this parameter should not change answers.", &
128  default=.false., do_not_log=.true., debuggingparam=.true.)
129 
130  if (cs%Analytic_FV_PGF .and. cs%blocked_AFV) then
131  call pressureforce_blk_afv_init(time, g, gv, us, param_file, diag, &
132  cs%PressureForce_blk_AFV_CSp, tides_csp)
133  elseif (cs%Analytic_FV_PGF) then
134  call pressureforce_afv_init(time, g, gv, us, param_file, diag, &
135  cs%PressureForce_AFV_CSp, tides_csp)
136  else
137  call pressureforce_mont_init(time, g, gv, us, param_file, diag, &
138  cs%PressureForce_Mont_CSp, tides_csp)
139  endif
140 

References mom_error_handler::mom_error().

Referenced by mom_dynamics_split_rk2::initialize_dyn_split_rk2(), mom_dynamics_unsplit::initialize_dyn_unsplit(), and mom_dynamics_unsplit_rk2::initialize_dyn_unsplit_rk2().

Here is the call graph for this function:
Here is the caller graph for this function: