MOM_PressureForce.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> A thin wrapper for Boussinesq/non-Boussinesq forms of the pressure force calculation.
7
8use mom_diag_mediator, only : diag_ctrl, time_type
9use mom_error_handler, only : mom_error, mom_mesg, fatal, warning, is_root_pe
10use mom_file_parser, only : get_param, log_version, param_file_type
11use mom_grid, only : ocean_grid_type
18use mom_self_attr_load, only : sal_cs
23use mom_ale, only: ale_cs
24implicit none ; private
25
26#include <MOM_memory.h>
27
29
30!> Pressure force control structure
31type, 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 !> Control structure for the analytically integrated finite volume pressure force
35 type(pressureforce_fv_cs) :: pressureforce_fv
36 !> Control structure for the Montgomery potential form of pressure force
37 type(pressureforce_mont_cs) :: pressureforce_mont
38end type pressureforce_cs
39
40contains
41
42!> A thin layer between the model and the Boussinesq and non-Boussinesq pressure force routines.
43subroutine pressureforce(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, ADp, p_atm, pbce, eta)
44 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
45 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
46 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
47 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
48 intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
49 type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
50 real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), &
51 intent(out) :: pfu !< Zonal pressure force acceleration [L T-2 ~> m s-2]
52 real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), &
53 intent(out) :: pfv !< Meridional pressure force acceleration [L T-2 ~> m s-2]
54 type(pressureforce_cs), intent(inout) :: cs !< Pressure force control structure
55 type(ale_cs), pointer :: ale_csp !< ALE control structure
56 type(accel_diag_ptrs), pointer :: adp !< Acceleration diagnostic pointers
57 real, dimension(:,:), pointer :: p_atm !< The pressure at the ice-ocean or
58 !! atmosphere-ocean interface [R L2 T-2 ~> Pa].
59 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
60 optional, intent(out) :: pbce !< The baroclinic pressure anomaly in each layer
61 !! due to eta anomalies [L2 T-2 H-1 ~> m s-2 or m4 s-2 kg-1].
62 real, dimension(SZI_(G),SZJ_(G)), &
63 optional, intent(out) :: eta !< The bottom mass used to calculate PFu and PFv,
64 !! [H ~> m or kg m-2], with any tidal contributions.
65
66 if (cs%Analytic_FV_PGF) then
67 if (gv%Boussinesq) then
68 call pressureforce_fv_bouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_FV, &
69 ale_csp, adp, p_atm, pbce, eta)
70 else
71 call pressureforce_fv_nonbouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_FV, &
72 ale_csp, adp, p_atm, pbce, eta)
73 endif
74 else
75 if (gv%Boussinesq) then
76 call pressureforce_mont_bouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_Mont, &
77 p_atm, pbce, eta)
78 else
79 call pressureforce_mont_nonbouss(h, tv, pfu, pfv, g, gv, us, cs%PressureForce_Mont, &
80 p_atm, pbce, eta)
81 endif
82 endif
83
84end subroutine pressureforce
85
86!> Initialize the pressure force control structure
87subroutine pressureforce_init(Time, G, GV, US, param_file, diag, CS, ADp, SAL_CSp, tides_CSp)
88 type(time_type), target, intent(in) :: time !< Current model time
89 type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
90 type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure
91 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
92 type(param_file_type), intent(in) :: param_file !< Parameter file handles
93 type(diag_ctrl), target, intent(inout) :: diag !< Diagnostics control structure
94 type(pressureforce_cs), intent(inout) :: cs !< Pressure force control structure
95 type(accel_diag_ptrs), pointer :: adp !< Acceleration diagnostic pointers
96 type(sal_cs), intent(in), optional :: sal_csp !< SAL control structure
97 type(tidal_forcing_cs), intent(in), optional :: tides_csp !< Tide control structure
98#include "version_variable.h"
99 character(len=40) :: mdl = "MOM_PressureForce" ! This module's name.
100
101 ! Read all relevant parameters and write them to the model log.
102 call log_version(param_file, mdl, version, "")
103 call get_param(param_file, mdl, "ANALYTIC_FV_PGF", cs%Analytic_FV_PGF, &
104 "If true the pressure gradient forces are calculated "//&
105 "with a finite volume form that analytically integrates "//&
106 "the equations of state in pressure to avoid any "//&
107 "possibility of numerical thermobaric instability, as "//&
108 "described in Adcroft et al., O. Mod. (2008).", default=.true.)
109
110 if (cs%Analytic_FV_PGF) then
111 call pressureforce_fv_init(time, g, gv, us, param_file, diag, &
112 cs%PressureForce_FV, adp, sal_csp, tides_csp)
113 else
114 call pressureforce_mont_init(time, g, gv, us, param_file, diag, &
115 cs%PressureForce_Mont, sal_csp, tides_csp)
116 endif
117end subroutine pressureforce_init
118
119!> \namespace mom_pressureforce
120!!
121!! This thin module provides a branch to two forms of the horizontal accelerations
122!! due to pressure gradients. The two options currently available are a
123!! Montgomery potential form (used in traditional isopycnal layer models), and the
124!! analytic finite volume form.
125
126end module mom_pressureforce