MOM6
MOM_PressureForce.F90
Go to the documentation of this file.
1 !> A thin wrapper for Boussinesq/non-Boussinesq forms of the pressure force calculation.
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_diag_mediator, only : diag_ctrl, time_type
7 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning, is_root_pe
9 use mom_grid, only : ocean_grid_type
23 use mom_ale, only: ale_cs
24 implicit none ; private
25 
26 #include <MOM_memory.h>
27 
29 
30 !> Pressure force control structure
31 type, public :: pressureforce_cs ; private
32  logical :: analytic_fv_pgf !< If true, use the analytic finite volume form
33  !! (Adcroft et al., Ocean Mod. 2008) of the PGF.
34  logical :: blocked_afv !< If true, used the blocked version of the ANALYTIC_FV_PGF
35  !! code. The value of this parameter should not change answers.
36  !> Control structure for the analytically integrated finite volume pressure force
37  type(pressureforce_afv_cs), pointer :: pressureforce_afv_csp => null()
38  !> Control structure for the analytically integrated finite volume pressure force
39  type(pressureforce_blk_afv_cs), pointer :: pressureforce_blk_afv_csp => null()
40  !> Control structure for the Montgomery potential form of pressure force
41  type(pressureforce_mont_cs), pointer :: pressureforce_mont_csp => null()
42 end type pressureforce_cs
43 
44 contains
45 
46 !> A thin layer between the model and the Boussinesq and non-Boussinesq pressure force routines.
47 subroutine pressureforce(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm, pbce, eta)
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 
96 end subroutine pressureforce
97 
98 !> Initialize the pressure force control structure
99 subroutine pressureforce_init(Time, G, GV, US, param_file, diag, CS, tides_CSp)
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 
141 end subroutine pressureforce_init
142 
143 !> Deallocate the pressure force control structure
144 subroutine pressureforce_end(CS)
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)
156 end subroutine pressureforce_end
157 
158 !> \namespace mom_pressureforce
159 !!
160 !! This thin module provides a branch to two forms of the horizontal accelerations
161 !! due to pressure gradients. The two options currently available are a
162 !! Montgomery potential form (used in traditional isopycnal layer models), and the
163 !! analytic finite volume form.
164 
165 end module mom_pressureforce
mom_pressureforce_blk_afv
Analytically integrated finite volume pressure gradient.
Definition: MOM_PressureForce_blocked_AFV.F90:2
mom_pressureforce_blk_afv::pressureforce_blk_afv_end
subroutine, public pressureforce_blk_afv_end(CS)
Deallocates the finite volume pressure gradient control structure.
Definition: MOM_PressureForce_blocked_AFV.F90:856
mom_pressureforce_afv::pressureforce_afv_nonbouss
subroutine, public pressureforce_afv_nonbouss(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm, pbce, eta)
Non-Boussinesq analytically-integrated finite volume form of pressure gradient.
Definition: MOM_PressureForce_analytic_FV.F90:103
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_pressureforce_blk_afv::pressureforce_blk_afv_cs
Finite volume pressure gradient control structure.
Definition: MOM_PressureForce_blocked_AFV.F90:36
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_pressureforce::pressureforce_init
subroutine, public pressureforce_init(Time, G, GV, US, param_file, diag, CS, tides_CSp)
Initialize the pressure force control structure.
Definition: MOM_PressureForce.F90:100
mom_pressureforce::pressureforce_end
subroutine, public pressureforce_end(CS)
Deallocate the pressure force control structure.
Definition: MOM_PressureForce.F90:145
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_pressureforce_blk_afv::pressureforce_blk_afv_bouss
subroutine, public pressureforce_blk_afv_bouss(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm, pbce, eta)
Boussinesq analytically-integrated finite volume form of pressure gradient.
Definition: MOM_PressureForce_blocked_AFV.F90:427
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_pressureforce_blk_afv::pressureforce_blk_afv_nonbouss
subroutine, public pressureforce_blk_afv_nonbouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pbce, eta)
Non-Boussinesq analytically-integrated finite volume form of pressure gradient.
Definition: MOM_PressureForce_blocked_AFV.F90:103
mom_pressureforce_mont::pressureforce_mont_nonbouss
subroutine, public pressureforce_mont_nonbouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pbce, eta)
Non-Boussinesq Montgomery-potential form of pressure gradient.
Definition: MOM_PressureForce_Montgomery.F90:65
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_tidal_forcing
Tidal contributions to geopotential.
Definition: MOM_tidal_forcing.F90:2
mom_tidal_forcing::tidal_forcing_cs
The control structure for the MOM_tidal_forcing module.
Definition: MOM_tidal_forcing.F90:26
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_pressureforce_afv::pressureforce_afv_cs
Finite volume pressure gradient control structure.
Definition: MOM_PressureForce_analytic_FV.F90:36
mom_ale
This module contains the main regridding routines.
Definition: MOM_ALE.F90:9
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_pressureforce_afv::pressureforce_afv_end
subroutine, public pressureforce_afv_end(CS)
Deallocates the finite volume pressure gradient control structure.
Definition: MOM_PressureForce_analytic_FV.F90:864
mom_pressureforce_mont::pressureforce_mont_bouss
subroutine, public pressureforce_mont_bouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pbce, eta)
Boussinesq Montgomery-potential form of pressure gradient.
Definition: MOM_PressureForce_Montgomery.F90:366
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_ale::ale_cs
ALE control structure.
Definition: MOM_ALE.F90:65
mom_pressureforce_afv
Analytically integrated finite volume pressure gradient.
Definition: MOM_PressureForce_analytic_FV.F90:2
mom_pressureforce_mont::pressureforce_mont_init
subroutine, public pressureforce_mont_init(Time, G, GV, US, param_file, diag, CS, tides_CSp)
Initialize the Montgomery-potential form of PGF control structure.
Definition: MOM_PressureForce_Montgomery.F90:823
mom_pressureforce_mont
Provides the Montgomery potential form of pressure gradient.
Definition: MOM_PressureForce_Montgomery.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_pressureforce_mont::pressureforce_mont_cs
Control structure for the Montgomery potential form of pressure gradient.
Definition: MOM_PressureForce_Montgomery.F90:31
mom_pressureforce_mont::pressureforce_mont_end
subroutine, public pressureforce_mont_end(CS)
Deallocates the Montgomery-potential form of PGF control structure.
Definition: MOM_PressureForce_Montgomery.F90:891
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_pressureforce::pressureforce
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.
Definition: MOM_PressureForce.F90:48
mom_pressureforce_afv::pressureforce_afv_init
subroutine, public pressureforce_afv_init(Time, G, GV, US, param_file, diag, CS, tides_CSp)
Initializes the finite volume pressure gradient control structure.
Definition: MOM_PressureForce_analytic_FV.F90:791
mom_pressureforce::pressureforce_cs
Pressure force control structure.
Definition: MOM_PressureForce.F90:31
mom_pressureforce_blk_afv::pressureforce_blk_afv_init
subroutine, public pressureforce_blk_afv_init(Time, G, GV, US, param_file, diag, CS, tides_CSp)
Initializes the finite volume pressure gradient control structure.
Definition: MOM_PressureForce_blocked_AFV.F90:783
mom_pressureforce
A thin wrapper for Boussinesq/non-Boussinesq forms of the pressure force calculation.
Definition: MOM_PressureForce.F90:2
mom_pressureforce_afv::pressureforce_afv_bouss
subroutine, public pressureforce_afv_bouss(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm, pbce, eta)
Boussinesq analytically-integrated finite volume form of pressure gradient.
Definition: MOM_PressureForce_analytic_FV.F90:446
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239