MOM_variables.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!> Provides transparent structures with groups of MOM6 variables and supporting routines
6module mom_variables
7
8use mom_array_transform, only : rotate_array, rotate_vector
9use mom_coupler_types, only : coupler_1d_bc_type, coupler_2d_bc_type
10use mom_coupler_types, only : coupler_type_spawn, coupler_type_destructor, coupler_type_initialized
11use mom_coupler_types, only : coupler_type_copy_data
12use mom_debugging, only : hchksum
13use mom_domains, only : mom_domain_type, get_domain_extent, group_pass_type
14use mom_eos, only : eos_type
15use mom_error_handler, only : mom_error, fatal
16use mom_grid, only : ocean_grid_type
20
21implicit none ; private
22
23#include <MOM_memory.h>
24
28
29! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
30! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
31! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
32! vary with the Boussinesq approximation, the Boussinesq variant is given first.
33
34!> A structure for creating arrays of pointers to 3D arrays
35type, public :: p3d
36 real, dimension(:,:,:), pointer :: p => null() !< A pointer to a 3D array [various]
37end type p3d
38!> A structure for creating arrays of pointers to 2D arrays
39type, public :: p2d
40 real, dimension(:,:), pointer :: p => null() !< A pointer to a 2D array [various]
41end type p2d
42
43!> Pointers to various fields which may be used describe the surface state of MOM, and which
44!! will be returned to the calling program
45type, public :: surface
46 real, allocatable, dimension(:,:) :: &
47 sst, & !< The sea surface temperature [C ~> degC].
48 sss, & !< The sea surface salinity [S ~> psu or gSalt/kg].
49 sfc_density, & !< The mixed layer density [R ~> kg m-3].
50 hml, & !< The mixed layer depth [Z ~> m].
51 u, & !< The mixed layer zonal velocity [L T-1 ~> m s-1].
52 v, & !< The mixed layer meridional velocity [L T-1 ~> m s-1].
53 sea_lev, & !< The sea level [Z ~> m]. If a reduced surface gravity is
54 !! used, that is compensated for in sea_lev.
55 frazil, & !< The energy needed to heat the ocean column to the freezing point during
56 !! the call to step_MOM [Q R Z ~> J m-2].
57 melt_potential, & !< Instantaneous amount of heat that can be used to melt sea ice [Q R Z ~> J m-2].
58 !! This is computed w.r.t. surface freezing temperature.
59 ocean_mass, & !< The total mass of the ocean [R Z ~> kg m-2].
60 ocean_heat, & !< The total heat content of the ocean in [C R Z ~> degC kg m-2].
61 ocean_salt, & !< The total salt content of the ocean in [1e-3 S R Z ~> kgSalt m-2].
62 taux_shelf, & !< The zonal stresses on the ocean under shelves [R L Z T-2 ~> Pa].
63 tauy_shelf, & !< The meridional stresses on the ocean under shelves [R L Z T-2 ~> Pa].
64 fco2 !< CO2 flux from the ocean to the atmosphere [R Z T-1 ~> kgCO2 m-2 s-1]
65 logical :: t_is_cont = .false. !< If true, the temperature variable SST is actually the
66 !! conservative temperature in [C ~> degC].
67 logical :: s_is_abss = .false. !< If true, the salinity variable SSS is actually the
68 !! absolute salinity in [S ~> gSalt kg-1].
69 type(coupler_2d_bc_type) :: tr_fields !< A structure that may contain an
70 !! array of named fields describing tracer-related quantities.
71 !### NOTE: ALL OF THE ARRAYS IN TR_FIELDS USE THE COUPLER'S INDEXING CONVENTION AND HAVE NO
72 !### HALOS! THIS IS DONE TO CONFORM TO THE TREATMENT IN MOM4, BUT I DON'T LIKE IT! -RWH
73 logical :: arrays_allocated = .false. !< A flag that indicates whether the surface type
74 !! has had its memory allocated.
75end type surface
76
77!> Pointers to an assortment of thermodynamic fields that may be available, including
78!! potential temperature, salinity, heat capacity, and the equation of state control structure.
79type, public :: thermo_var_ptrs
80 ! If allocated, the following variables have nz layers.
81 real, pointer :: t(:,:,:) => null() !< Potential temperature [C ~> degC].
82 real, pointer :: s(:,:,:) => null() !< Salinity [PSU] or [gSalt/kg], generically [S ~> ppt].
83 real, pointer :: p_surf(:,:) => null() !< Ocean surface pressure used in equation of state
84 !! calculations [R L2 T-2 ~> Pa]
85 type(eos_type), pointer :: eqn_of_state => null() !< Type that indicates the
86 !! equation of state to use.
87 real :: p_ref !< The coordinate-density reference pressure [R L2 T-2 ~> Pa].
88 !! This is the pressure used to calculate Rml from
89 !! T and S when eqn_of_state is associated.
90 real :: c_p !< The heat capacity of seawater [Q C-1 ~> J degC-1 kg-1].
91 !! When conservative temperature is used, this is
92 !! constant and exactly 3991.86795711963 J degC-1 kg-1.
93 logical :: t_is_cont = .false. !< If true, the temperature variable tv%T is
94 !! actually the conservative temperature [C ~> degC].
95 logical :: s_is_abss = .false. !< If true, the salinity variable tv%S is
96 !! actually the absolute salinity in units of [S ~> gSalt kg-1].
97 real :: min_salinity !< The minimum value of salinity when BOUND_SALINITY=True [S ~> ppt].
98 real, allocatable, dimension(:,:,:) :: spv_avg
99 !< The layer averaged in situ specific volume [R-1 ~> m3 kg-1].
100 integer :: valid_spv_halo = -1 !< If positive, the valid halo size for SpV_avg, or if negative
101 !! SpV_avg is not currently set.
102
103 ! These arrays are accumulated fluxes for communication with other components.
104 real, dimension(:,:), pointer :: frazil => null()
105 !< The energy needed to heat the ocean column to the
106 !! freezing point since calculate_surface_state was
107 !! last called [Q Z R ~> J m-2].
108 logical :: frazil_was_reset !< If true, frazil has not accumulated since it was last reset.
109 real, dimension(:,:), pointer :: salt_deficit => null()
110 !< The salt needed to maintain the ocean column
111 !! at a minimum salinity of MIN_SALINITY since the last time
112 !! that calculate_surface_state was called, [S R Z ~> gSalt m-2].
113 real, dimension(:,:), pointer :: tempxpme => null()
114 !< The net inflow of water into the ocean times the
115 !! temperature at which this inflow occurs since the
116 !! last call to calculate_surface_state [C R Z ~> degC kg m-2].
117 !! This should be prescribed in the forcing fields, but
118 !! as it often is not, this is a useful heat budget diagnostic.
119 real, dimension(:,:), pointer :: internal_heat => null()
120 !< Any internal or geothermal heat sources that
121 !! have been applied to the ocean since the last call to
122 !! calculate_surface_state [C R Z ~> degC kg m-2].
123 ! The following variables are most normally not used but when they are they
124 ! will be either set by parameterizations or prognostic.
125 real, pointer :: vart(:,:,:) => null() !< SGS variance of potential temperature [C2 ~> degC2].
126 real, pointer :: vars(:,:,:) => null() !< SGS variance of salinity [S2 ~> ppt2].
127 real, pointer :: covarts(:,:,:) => null() !< SGS covariance of salinity and potential
128 !! temperature [C S ~> degC ppt].
129 type(tracer_type), pointer :: tr_t => null() !< pointer to temp in tracer registry
130 type(tracer_type), pointer :: tr_s => null() !< pointer to salinty in tracer registry
131end type thermo_var_ptrs
132
133!> Pointers to all of the prognostic variables allocated in MOM_variables.F90 and MOM.F90.
134!!
135!! It is useful for sending these variables for diagnostics, and in preparation for ensembles
136!! later on. All variables have the same names as the local (public) variables
137!! they refer to in MOM.F90.
138type, public :: ocean_internal_state
139 real, pointer, dimension(:,:,:) :: &
140 t => null(), & !< Pointer to the temperature state variable [C ~> degC]
141 s => null(), & !< Pointer to the salinity state variable [S ~> ppt] (i.e., PSU or g/kg)
142 u => null(), & !< Pointer to the zonal velocity [L T-1 ~> m s-1]
143 v => null(), & !< Pointer to the meridional velocity [L T-1 ~> m s-1]
144 h => null() !< Pointer to the layer thicknesses [H ~> m or kg m-2]
145 real, pointer, dimension(:,:,:) :: &
146 uh => null(), & !< Pointer to zonal transports [H L2 T-1 ~> m3 s-1 or kg s-1]
147 vh => null() !< Pointer to meridional transports [H L2 T-1 ~> m3 s-1 or kg s-1]
148 real, pointer, dimension(:,:,:) :: &
149 cau => null(), & !< Pointer to the zonal Coriolis and Advective acceleration [L T-2 ~> m s-2]
150 cav => null(), & !< Pointer to the meridional Coriolis and Advective acceleration [L T-2 ~> m s-2]
151 pfu => null(), & !< Pointer to the zonal Pressure force acceleration [L T-2 ~> m s-2]
152 pfv => null(), & !< Pointer to the meridional Pressure force acceleration [L T-2 ~> m s-2]
153 diffu => null(), & !< Pointer to the zonal acceleration due to lateral viscosity [L T-2 ~> m s-2]
154 diffv => null(), & !< Pointer to the meridional acceleration due to lateral viscosity [L T-2 ~> m s-2]
155 pbce => null(), & !< Pointer to the baroclinic pressure force dependency on free surface movement
156 !! [L2 T-2 H-1 ~> m s-2 or m4 kg-1 s-2]
157 u_accel_bt => null(), & !< Pointer to the zonal barotropic-solver acceleration [L T-2 ~> m s-2]
158 v_accel_bt => null() !< Pointer to the meridional barotropic-solver acceleration [L T-2 ~> m s-2]
159 real, pointer, dimension(:,:,:) :: &
160 u_av => null(), & !< Pointer to zonal velocity averaged over the timestep [L T-1 ~> m s-1]
161 v_av => null(), & !< Pointer to meridional velocity averaged over the timestep [L T-1 ~> m s-1]
162 u_prev => null(), & !< Pointer to zonal velocity at the end of the last timestep [L T-1 ~> m s-1]
163 v_prev => null() !< Pointer to meridional velocity at the end of the last timestep [L T-1 ~> m s-1]
165
166!> Pointers to arrays with accelerations, which can later be used for derived diagnostics, like energy balances.
167type, public :: accel_diag_ptrs
168
169 ! Each of the following fields has nz layers.
170 real, pointer, dimension(:,:,:) :: &
171 diffu => null(), & !< Zonal acceleration due to along isopycnal viscosity [L T-2 ~> m s-2]
172 diffv => null(), & !< Meridional acceleration due to along isopycnal viscosity [L T-2 ~> m s-2]
173 cau => null(), & !< Zonal Coriolis and momentum advection accelerations [L T-2 ~> m s-2]
174 cav => null(), & !< Meridional Coriolis and momentum advection accelerations [L T-2 ~> m s-2]
175 pfu => null(), & !< Zonal acceleration due to pressure forces [L T-2 ~> m s-2]
176 pfv => null(), & !< Meridional acceleration due to pressure forces [L T-2 ~> m s-2]
177 du_dt_visc => null(), &!< Zonal acceleration due to vertical viscosity [L T-2 ~> m s-2]
178 dv_dt_visc => null(), &!< Meridional acceleration due to vertical viscosity [L T-2 ~> m s-2]
179 du_dt_visc_gl90 => null(), &!< Zonal acceleration due to GL90 vertical viscosity
180 ! (is included in du_dt_visc) [L T-2 ~> m s-2]
181 dv_dt_visc_gl90 => null(), &!< Meridional acceleration due to GL90 vertical viscosity
182 ! (is included in dv_dt_visc) [L T-2 ~> m s-2]
183 du_dt_str => null(), & !< Zonal acceleration due to the surface stress (included
184 !! in du_dt_visc) [L T-2 ~> m s-2]
185 dv_dt_str => null(), & !< Meridional acceleration due to the surface stress (included
186 !! in dv_dt_visc) [L T-2 ~> m s-2]
187 du_dt_dia => null(), & !< Zonal acceleration due to diapycnal mixing [L T-2 ~> m s-2]
188 dv_dt_dia => null(), & !< Meridional acceleration due to diapycnal mixing [L T-2 ~> m s-2]
189 u_accel_bt => null(), &!< Pointer to the zonal barotropic-solver acceleration [L T-2 ~> m s-2]
190 v_accel_bt => null(), &!< Pointer to the meridional barotropic-solver acceleration [L T-2 ~> m s-2]
191
192 ! sal_[uv] and tide_[uv] are 3D fields because of their baroclinic component in Boussinesq mode.
193 sal_u => null(), & !< Zonal acceleration due to self-attraction and loading [L T-2 ~> m s-2]
194 sal_v => null(), & !< Meridional acceleration due to self-attraction and loading [L T-2 ~> m s-2]
195 tides_u => null(), & !< Zonal acceleration due to astronomical tidal forcing [L T-2 ~> m s-2]
196 tides_v => null() !< Meridional acceleration due to astronomical tidal forcing [L T-2 ~> m s-2]
197 real, pointer, dimension(:,:,:) :: du_other => null()
198 !< Zonal velocity changes due to any other processes that are
199 !! not due to any explicit accelerations [L T-1 ~> m s-1].
200 real, pointer, dimension(:,:,:) :: dv_other => null()
201 !< Meridional velocity changes due to any other processes that are
202 !! not due to any explicit accelerations [L T-1 ~> m s-1].
203
204 ! Sub-terms of [uv]_accel_bt
205 real, pointer :: bt_pgf_u(:,:,:) => null() !< Zonal acceleration due to anomalous pressure gradient from
206 !! barotropic solver, a 3D component of u_accel_bt that includes both
207 !! PFuBT and the offset term for central differencing timestepping
208 !! [L T-2 ~> m s-2]
209 real, pointer :: bt_pgf_v(:,:,:) => null() !< Meridional acceleration due to anomalous pressure gradient from
210 !! barotropic solver, a 3D component of v_accel_bt that includes both
211 !! PFvBT and the offset term for central differencing timestepping
212 !! [L T-2 ~> m s-2]
213 real, pointer :: bt_cor_u(:,:) => null() !< Zonal acceleration due to anomalous Coriolis force from barotropic
214 !! solver, a 2D component of u_accel_bt [L T-2 ~> m s-2]
215 real, pointer :: bt_cor_v(:,:) => null() !< Meridional acceleration due to anomalous Coriolis force from barotropic
216 !! solver, a 2D component of v_accel_bt [L T-2 ~> m s-2]
217 real, pointer :: bt_lwd_u(:,:) => null() !< Zonal acceleration due to linear wave drag from barotropic solver,
218 !! a 2D component of u_accel_bt [L T-2 ~> m s-2]
219 real, pointer :: bt_lwd_v(:,:) => null() !< Meridional acceleration due to linear wave drag from barotropic solver,
220 !! a 2D component of v_accel_bt [L T-2 ~> m s-2]
221
222 ! These accelerations are sub-terms included in the accelerations above.
223 real, pointer :: gradkeu(:,:,:) => null() !< gradKEu = - d/dx(u2) [L T-2 ~> m s-2]
224 real, pointer :: gradkev(:,:,:) => null() !< gradKEv = - d/dy(u2) [L T-2 ~> m s-2]
225 real, pointer :: rv_x_v(:,:,:) => null() !< rv_x_v = rv * v at u [L T-2 ~> m s-2]
226 real, pointer :: rv_x_u(:,:,:) => null() !< rv_x_u = rv * u at v [L T-2 ~> m s-2]
227
228 real, pointer :: diag_hfrac_u(:,:,:) => null() !< Fractional layer thickness at u points [nondim]
229 real, pointer :: diag_hfrac_v(:,:,:) => null() !< Fractional layer thickness at v points [nondim]
230 real, pointer :: diag_hu(:,:,:) => null() !< layer thickness at u points, modulated by the viscous
231 !! remnant and fractional open areas [H ~> m or kg m-2]
232 real, pointer :: diag_hv(:,:,:) => null() !< layer thickness at v points, modulated by the viscous
233 !! remnant and fractional open areas [H ~> m or kg m-2]
234
235 real, pointer :: visc_rem_u(:,:,:) => null() !< viscous remnant at u points [nondim]
236 real, pointer :: visc_rem_v(:,:,:) => null() !< viscous remnant at v points [nondim]
237
238end type accel_diag_ptrs
239
240!> Pointers to arrays with transports, which can later be used for derived diagnostics, like energy balances.
241type, public :: cont_diag_ptrs
242
243! Each of the following fields has nz layers.
244 real, pointer, dimension(:,:,:) :: &
245 uh => null(), & !< Resolved zonal layer thickness fluxes, [H L2 T-1 ~> m3 s-1 or kg s-1]
246 vh => null(), & !< Resolved meridional layer thickness fluxes, [H L2 T-1 ~> m3 s-1 or kg s-1]
247 uh_smooth => null(), & !< Interface height smoothing induced zonal volume fluxes [H L2 T-1 ~> m3 s-1 or kg s-1]
248 vh_smooth => null(), & !< Interface height smoothing induced meridional volume fluxes [H L2 T-1 ~> m3 s-1 or kg s-1]
249 uhgm => null(), & !< Isopycnal height diffusion induced zonal volume fluxes [H L2 T-1 ~> m3 s-1 or kg s-1]
250 vhgm => null() !< Isopycnal height diffusion induced meridional volume fluxes [H L2 T-1 ~> m3 s-1 or kg s-1]
251
252! Each of the following fields is found at nz+1 interfaces.
253 real, pointer :: diapyc_vel(:,:,:) => null() !< The net diapycnal velocity [H T-1 ~> m s-1 or kg m-2 s-1]
254
255end type cont_diag_ptrs
256
257!> Vertical viscosities, drag coefficients, and related fields.
258type, public :: vertvisc_type
259 real, allocatable, dimension(:,:) :: &
260 bbl_thick_u, & !< The bottom boundary layer thickness at the u-points [Z ~> m].
261 bbl_thick_v, & !< The bottom boundary layer thickness at the v-points [Z ~> m].
262 kv_bbl_u, & !< The bottom boundary layer viscosity at the u-points [H Z T-1 ~> m2 s-1 or Pa s]
263 kv_bbl_v, & !< The bottom boundary layer viscosity at the v-points [H Z T-1 ~> m2 s-1 or Pa s]
264 ustar_bbl, & !< The turbulence velocity in the bottom boundary layer at
265 !! h points [H T-1 ~> m s-1 or kg m-2 s-1].
266 bbl_meanke_loss, & !< The viscous loss of mean kinetic energy in the bottom boundary layer
267 !! [H L2 T-3 ~> m3 s-3 or W m-2].
268 bbl_meanke_loss_sqrtcd, & !< The viscous loss of mean kinetic energy in the bottom boundary layer
269 !! divided by the square root of the drag coefficient [H L2 T-3 ~> m3 s-3 or W m-2].
270 !! This is being set only to retain old answers, and should be phased out.
271 taux_shelf, & !< The zonal stresses on the ocean under shelves [R Z L T-2 ~> Pa].
272 tauy_shelf !< The meridional stresses on the ocean under shelves [R Z L T-2 ~> Pa].
273 real, allocatable, dimension(:,:) :: tbl_thick_shelf_u
274 !< Thickness of the viscous top boundary layer under ice shelves at u-points [Z ~> m].
275 real, allocatable, dimension(:,:) :: tbl_thick_shelf_v
276 !< Thickness of the viscous top boundary layer under ice shelves at v-points [Z ~> m].
277 real, allocatable, dimension(:,:) :: kv_tbl_shelf_u
278 !< Viscosity in the viscous top boundary layer under ice shelves at
279 !! u-points [H Z T-1 ~> m2 s-1 or Pa s]
280 real, allocatable, dimension(:,:) :: kv_tbl_shelf_v
281 !< Viscosity in the viscous top boundary layer under ice shelves at
282 !! v-points [H Z T-1 ~> m2 s-1 or Pa s]
283 real, allocatable, dimension(:,:) :: nkml_visc_u
284 !< The number of layers in the viscous surface mixed layer at u-points [nondim].
285 !! This is not an integer because there may be fractional layers, and it is stored in
286 !! terms of layers, not depth, to facilitate the movement of the viscous boundary layer
287 !! with the flow.
288 real, allocatable, dimension(:,:) :: nkml_visc_v
289 !< The number of layers in the viscous surface mixed layer at v-points [nondim].
290 real, allocatable, dimension(:,:,:) :: &
291 ray_u, & !< The Rayleigh drag velocity to be applied to each layer at u-points [H T-1 ~> m s-1 or Pa s m-1].
292 ray_v !< The Rayleigh drag velocity to be applied to each layer at v-points [H T-1 ~> m s-1 or Pa s m-1].
293
294 ! The following elements are pointers so they can be used as targets for pointers in the restart registry.
295 real, pointer, dimension(:,:) :: mld => null()
296 !< Instantaneous active mixing layer depth as used by the mixed layer restratification
297 !! parameterization [Z ~> m].
298 real, pointer, dimension(:,:) :: h_ml => null()
299 !< Instantaneous active mixing layer thickness as used by the mixed layer restratification
300 !! parameterization [H ~> m or kg m-2]
301 real, pointer, dimension(:,:) :: mld_param => null()
302 !< Instantaneous active mixed or mixing layer depth as used by the brine plume parameterization.
303 !! It could be coordinated with MLD above but we may want the ability to use different scales
304 !! in different parameterizations [Z ~> m].
305 real, pointer, dimension(:,:) :: h_ml_param => null()
306 !< Instantaneous active mixed or mixing layer thickness as used by the brine plume
307 !! parameterization [H ~> m or kg m-2].
308 real, pointer, dimension(:,:) :: sfc_buoy_flx => null() !< Surface buoyancy flux (derived) [Z2 T-3 ~> m2 s-3].
309 real, pointer, dimension(:,:,:) :: kd_shear => null()
310 !< The shear-driven turbulent diapycnal diffusivity at the interfaces between layers
311 !! in tracer columns [H Z T-1 ~> m2 s-1 or kg m-1 s-1]
312 real, pointer, dimension(:,:,:) :: kv_shear => null()
313 !< The shear-driven turbulent vertical viscosity at the interfaces between layers
314 !! in tracer columns [H Z T-1 ~> m2 s-1 or Pa s]
315 real, pointer, dimension(:,:,:) :: kv_shear_bu => null()
316 !< The shear-driven turbulent vertical viscosity at the interfaces between layers in
317 !! corner columns [H Z T-1 ~> m2 s-1 or Pa s]
318 real, pointer, dimension(:,:,:) :: kv_slow => null()
319 !< The turbulent vertical viscosity component due to "slow" processes (e.g., tidal,
320 !! background, convection etc) [H Z T-1 ~> m2 s-1 or Pa s]
321 real, pointer, dimension(:,:,:) :: tke_turb => null()
322 !< The turbulent kinetic energy per unit mass at the interfaces [Z2 T-2 ~> m2 s-2].
323 !! This may be at the tracer or corner points
324end type vertvisc_type
325
326!> Container for information about the summed layer transports
327!! and how they will vary as the barotropic velocity is changed.
328type, public :: bt_cont_type
329 real, allocatable :: fa_u_ee(:,:) !< The effective open face area for zonal barotropic transport
330 !! drawing from locations far to the east [H L ~> m2 or kg m-1].
331 real, allocatable :: fa_u_e0(:,:) !< The effective open face area for zonal barotropic transport
332 !! drawing from nearby to the east [H L ~> m2 or kg m-1].
333 real, allocatable :: fa_u_w0(:,:) !< The effective open face area for zonal barotropic transport
334 !! drawing from nearby to the west [H L ~> m2 or kg m-1].
335 real, allocatable :: fa_u_ww(:,:) !< The effective open face area for zonal barotropic transport
336 !! drawing from locations far to the west [H L ~> m2 or kg m-1].
337 real, allocatable :: ubt_ww(:,:) !< uBT_WW is the barotropic velocity [L T-1 ~> m s-1], beyond which the
338 !! marginal open face area is FA_u_WW. uBT_WW must be non-negative.
339 real, allocatable :: ubt_ee(:,:) !< uBT_EE is a barotropic velocity [L T-1 ~> m s-1], beyond which the
340 !! marginal open face area is FA_u_EE. uBT_EE must be non-positive.
341 real, allocatable :: fa_v_nn(:,:) !< The effective open face area for meridional barotropic transport
342 !! drawing from locations far to the north [H L ~> m2 or kg m-1].
343 real, allocatable :: fa_v_n0(:,:) !< The effective open face area for meridional barotropic transport
344 !! drawing from nearby to the north [H L ~> m2 or kg m-1].
345 real, allocatable :: fa_v_s0(:,:) !< The effective open face area for meridional barotropic transport
346 !! drawing from nearby to the south [H L ~> m2 or kg m-1].
347 real, allocatable :: fa_v_ss(:,:) !< The effective open face area for meridional barotropic transport
348 !! drawing from locations far to the south [H L ~> m2 or kg m-1].
349 real, allocatable :: vbt_ss(:,:) !< vBT_SS is the barotropic velocity, [L T-1 ~> m s-1], beyond which the
350 !! marginal open face area is FA_v_SS. vBT_SS must be non-negative.
351 real, allocatable :: vbt_nn(:,:) !< vBT_NN is the barotropic velocity, [L T-1 ~> m s-1], beyond which the
352 !! marginal open face area is FA_v_NN. vBT_NN must be non-positive.
353 real, allocatable :: h_u(:,:,:) !< An effective thickness at zonal faces, taking into account the effects
354 !! of vertical viscosity and fractional open areas [H ~> m or kg m-2].
355 !! This is primarily used as a non-normalized weight in determining
356 !! the depth averaged accelerations for the barotropic solver.
357 real, allocatable :: h_v(:,:,:) !< An effective thickness at meridional faces, taking into account the effects
358 !! of vertical viscosity and fractional open areas [H ~> m or kg m-2].
359 !! This is primarily used as a non-normalized weight in determining
360 !! the depth averaged accelerations for the barotropic solver.
361 type(group_pass_type) :: pass_polarity_bt !< Structure for polarity group halo updates
362 type(group_pass_type) :: pass_fa_uv !< Structure for face area group halo updates
363end type bt_cont_type
364
365!> Container for grids modifying cell metric at porous barriers
366type, public :: porous_barrier_type
367 ! Each of the following fields has nz layers.
368 real, allocatable :: por_face_areau(:,:,:) !< fractional open area of U-faces [nondim]
369 real, allocatable :: por_face_areav(:,:,:) !< fractional open area of V-faces [nondim]
370 ! Each of the following fields is found at nz+1 interfaces.
371 real, allocatable :: por_layer_widthu(:,:,:) !< fractional open width of U-faces [nondim]
372 real, allocatable :: por_layer_widthv(:,:,:) !< fractional open width of V-faces [nondim]
373end type porous_barrier_type
374
375contains
376
377!> Allocates the fields for the surface (return) properties of
378!! the ocean model. Unused fields are unallocated.
379subroutine allocate_surface_state(sfc_state, G, use_temperature, do_integrals, &
380 gas_fields_ocn, use_meltpot, use_iceshelves, &
381 omit_frazil, sfc_state_in, turns, use_marbl_tracers)
382 type(ocean_grid_type), intent(in) :: g !< ocean grid structure
383 type(surface), intent(inout) :: sfc_state !< ocean surface state type to be allocated.
384 logical, optional, intent(in) :: use_temperature !< If true, allocate the space for thermodynamic variables.
385 logical, optional, intent(in) :: do_integrals !< If true, allocate the space for vertically
386 !! integrated fields.
387 type(coupler_1d_bc_type), &
388 optional, intent(in) :: gas_fields_ocn !< If present, this type describes the
389 !! ocean and surface-ice fields that will participate
390 !! in the calculation of additional gas or other
391 !! tracer fluxes, and can be used to spawn related
392 !! internal variables in the ice model.
393 logical, optional, intent(in) :: use_meltpot !< If true, allocate the space for melt potential
394 logical, optional, intent(in) :: use_iceshelves !< If true, allocate the space for the stresses
395 !! under ice shelves.
396 logical, optional, intent(in) :: omit_frazil !< If present and false, do not allocate the space to
397 !! pass frazil fluxes to the coupler
398 type(surface), &
399 optional, intent(in) :: sfc_state_in !< If present and its tr_fields are initialized,
400 !! this type describes the ocean and surface-ice fields that
401 !! will participate in the calculation of additional gas or
402 !! other tracer fluxes, and can be used to spawn related
403 !! internal variables in the ice model. If gas_fields_ocn
404 !! is present, it is used and tr_fields_in is ignored.
405 integer, optional, intent(in) :: turns !< If present, the number of counterclockwise quarter
406 !! turns to use on the new grid.
407 logical, optional, intent(in) :: use_marbl_tracers !< If true, allocate the space for CO2 flux from MARBL
408
409 ! local variables
410 logical :: use_temp, alloc_integ, use_melt_potential, alloc_iceshelves, alloc_frazil, alloc_fco2
411 logical :: even_turns ! True if turns is absent or even
412 integer :: tr_field_i_mem(4), tr_field_j_mem(4)
413 integer :: is, ie, js, je, isd, ied, jsd, jed
414 integer :: isdb, iedb, jsdb, jedb
415
416 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
417 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
418 isdb = g%isdB ; iedb = g%iedB ; jsdb = g%jsdB ; jedb = g%jedB
419
420 use_temp = .true. ; if (present(use_temperature)) use_temp = use_temperature
421 alloc_integ = .true. ; if (present(do_integrals)) alloc_integ = do_integrals
422 use_melt_potential = .false. ; if (present(use_meltpot)) use_melt_potential = use_meltpot
423 alloc_iceshelves = .false. ; if (present(use_iceshelves)) alloc_iceshelves = use_iceshelves
424 alloc_frazil = .true. ; if (present(omit_frazil)) alloc_frazil = .not.omit_frazil
425 alloc_fco2 = .false. ; if (present(use_marbl_tracers)) alloc_fco2 = use_marbl_tracers
426
427 if (sfc_state%arrays_allocated) return
428
429 if (use_temp) then
430 allocate(sfc_state%SST(isd:ied,jsd:jed), source=0.0)
431 allocate(sfc_state%SSS(isd:ied,jsd:jed), source=0.0)
432 else
433 allocate(sfc_state%sfc_density(isd:ied,jsd:jed), source=0.0)
434 endif
435 if (use_temp .and. alloc_frazil) then
436 allocate(sfc_state%frazil(isd:ied,jsd:jed), source=0.0)
437 endif
438 allocate(sfc_state%sea_lev(isd:ied,jsd:jed), source=0.0)
439 allocate(sfc_state%Hml(isd:ied,jsd:jed), source=0.0)
440 allocate(sfc_state%u(isdb:iedb,jsd:jed), source=0.0)
441 allocate(sfc_state%v(isd:ied,jsdb:jedb), source=0.0)
442
443 if (use_melt_potential) then
444 allocate(sfc_state%melt_potential(isd:ied,jsd:jed), source=0.0)
445 endif
446
447 if (alloc_integ) then
448 ! Allocate structures for the vertically integrated ocean_mass, ocean_heat, and ocean_salt.
449 allocate(sfc_state%ocean_mass(isd:ied,jsd:jed), source=0.0)
450 if (use_temp) then
451 allocate(sfc_state%ocean_heat(isd:ied,jsd:jed), source=0.0)
452 allocate(sfc_state%ocean_salt(isd:ied,jsd:jed), source=0.0)
453 endif
454 endif
455
456 if (alloc_iceshelves) then
457 allocate(sfc_state%taux_shelf(isdb:iedb,jsd:jed), source=0.0)
458 allocate(sfc_state%tauy_shelf(isd:ied,jsdb:jedb), source=0.0)
459 endif
460
461 ! The data fields in the coupler_2d_bc_type are never rotated.
462 even_turns = .true. ; if (present(turns)) even_turns = (modulo(turns, 2) == 0)
463 if (even_turns) then
464 tr_field_i_mem(1:4) = (/is,is,ie,ie/) ; tr_field_j_mem(1:4) = (/js,js,je,je/)
465 else
466 tr_field_i_mem(1:4) = (/js,js,je,je/) ; tr_field_j_mem(1:4) = (/is,is,ie,ie/)
467 endif
468 if (present(gas_fields_ocn)) then
469 call coupler_type_spawn(gas_fields_ocn, sfc_state%tr_fields, &
470 tr_field_i_mem, tr_field_j_mem, as_needed=.true.)
471 elseif (present(sfc_state_in)) then
472 if (coupler_type_initialized(sfc_state_in%tr_fields)) then
473 call coupler_type_spawn(sfc_state_in%tr_fields, sfc_state%tr_fields, &
474 tr_field_i_mem, tr_field_j_mem, as_needed=.true.)
475 endif
476 endif
477
478 if (alloc_fco2) then
479 allocate(sfc_state%fco2(isd:ied,jsd:jed), source=0.0)
480 endif
481
482 sfc_state%arrays_allocated = .true.
483
484end subroutine allocate_surface_state
485
486!> Deallocates the elements of a surface state type.
487subroutine deallocate_surface_state(sfc_state)
488 type(surface), intent(inout) :: sfc_state !< ocean surface state type to be deallocated here.
489
490 if (.not.sfc_state%arrays_allocated) return
491
492 if (allocated(sfc_state%melt_potential)) deallocate(sfc_state%melt_potential)
493 if (allocated(sfc_state%SST)) deallocate(sfc_state%SST)
494 if (allocated(sfc_state%SSS)) deallocate(sfc_state%SSS)
495 if (allocated(sfc_state%sfc_density)) deallocate(sfc_state%sfc_density)
496 if (allocated(sfc_state%sea_lev)) deallocate(sfc_state%sea_lev)
497 if (allocated(sfc_state%Hml)) deallocate(sfc_state%Hml)
498 if (allocated(sfc_state%u)) deallocate(sfc_state%u)
499 if (allocated(sfc_state%v)) deallocate(sfc_state%v)
500 if (allocated(sfc_state%ocean_mass)) deallocate(sfc_state%ocean_mass)
501 if (allocated(sfc_state%ocean_heat)) deallocate(sfc_state%ocean_heat)
502 if (allocated(sfc_state%ocean_salt)) deallocate(sfc_state%ocean_salt)
503 if (allocated(sfc_state%fco2)) deallocate(sfc_state%fco2)
504 call coupler_type_destructor(sfc_state%tr_fields)
505
506 sfc_state%arrays_allocated = .false.
507
508end subroutine deallocate_surface_state
509
510!> Rotate the surface state fields from the input to the model indices.
511subroutine rotate_surface_state(sfc_state_in, sfc_state, G, turns)
512 type(surface), intent(in) :: sfc_state_in !< The input unrotated surface state type that is the data source.
513 type(surface), intent(inout) :: sfc_state !< The rotated surface state type whose arrays will be filled in
514 type(ocean_grid_type), intent(in) :: g !< The ocean grid structure
515 integer, intent(in) :: turns !< The number of counterclockwise quarter turns to use on the rotated grid.
516
517 logical :: use_temperature, do_integrals, use_melt_potential, use_iceshelves
518
519 ! NOTE: Many of these are weak tests, since only one is checked
520 use_temperature = allocated(sfc_state_in%SST) &
521 .and. allocated(sfc_state_in%SSS)
522 use_melt_potential = allocated(sfc_state_in%melt_potential)
523 do_integrals = allocated(sfc_state_in%ocean_mass)
524 use_iceshelves = allocated(sfc_state_in%taux_shelf) &
525 .and. allocated(sfc_state_in%tauy_shelf)
526
527 if (.not. sfc_state%arrays_allocated) then
528 call allocate_surface_state(sfc_state, g, use_temperature=use_temperature, &
529 do_integrals=do_integrals, use_meltpot=use_melt_potential, &
530 use_iceshelves=use_iceshelves, sfc_state_in=sfc_state_in, turns=turns)
531 endif
532
533 if (use_temperature) then
534 call rotate_array(sfc_state_in%SST, turns, sfc_state%SST)
535 call rotate_array(sfc_state_in%SSS, turns, sfc_state%SSS)
536 else
537 call rotate_array(sfc_state_in%sfc_density, turns, sfc_state%sfc_density)
538 endif
539
540 call rotate_array(sfc_state_in%Hml, turns, sfc_state%Hml)
541 call rotate_vector(sfc_state_in%u, sfc_state_in%v, turns, &
542 sfc_state%u, sfc_state%v)
543 call rotate_array(sfc_state_in%sea_lev, turns, sfc_state%sea_lev)
544
545 if (use_melt_potential) then
546 call rotate_array(sfc_state_in%melt_potential, turns, sfc_state%melt_potential)
547 endif
548
549 if (do_integrals) then
550 call rotate_array(sfc_state_in%ocean_mass, turns, sfc_state%ocean_mass)
551 if (use_temperature) then
552 call rotate_array(sfc_state_in%ocean_heat, turns, sfc_state%ocean_heat)
553 call rotate_array(sfc_state_in%ocean_salt, turns, sfc_state%ocean_salt)
554 call rotate_array(sfc_state_in%SSS, turns, sfc_state%SSS)
555 endif
556 endif
557
558 if (use_iceshelves) then
559 call rotate_vector(sfc_state_in%taux_shelf, sfc_state_in%tauy_shelf, turns, &
560 sfc_state%taux_shelf, sfc_state%tauy_shelf)
561 endif
562
563 if (use_temperature .and. allocated(sfc_state_in%frazil)) &
564 call rotate_array(sfc_state_in%frazil, turns, sfc_state%frazil)
565
566 ! Scalar transfers
567 sfc_state%T_is_conT = sfc_state_in%T_is_conT
568 sfc_state%S_is_absS = sfc_state_in%S_is_absS
569
570 ! NOTE: Tracer fields are handled by FMS, so are left unrotated. Any
571 ! reads/writes to tr_fields must be appropriately rotated.
572 if (coupler_type_initialized(sfc_state_in%tr_fields)) then
573 call coupler_type_copy_data(sfc_state_in%tr_fields, sfc_state%tr_fields)
574 endif
575end subroutine rotate_surface_state
576
577!> Allocates the arrays contained within a BT_cont_type and initializes them to 0.
578subroutine alloc_bt_cont_type(BT_cont, G, GV, alloc_faces)
579 type(bt_cont_type), pointer :: bt_cont !< The BT_cont_type whose elements will be allocated
580 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
581 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
582 logical, optional, intent(in) :: alloc_faces !< If present and true, allocate
583 !! memory for effective face thicknesses.
584
585 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb, nz
586 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = gv%ke
587 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
588
589 if (associated(bt_cont)) call mom_error(fatal, &
590 "alloc_BT_cont_type called with an associated BT_cont_type pointer.")
591
592 allocate(bt_cont)
593 allocate(bt_cont%FA_u_WW(isdb:iedb,jsd:jed), source=0.0)
594 allocate(bt_cont%FA_u_W0(isdb:iedb,jsd:jed), source=0.0)
595 allocate(bt_cont%FA_u_E0(isdb:iedb,jsd:jed), source=0.0)
596 allocate(bt_cont%FA_u_EE(isdb:iedb,jsd:jed), source=0.0)
597 allocate(bt_cont%uBT_WW(isdb:iedb,jsd:jed), source=0.0)
598 allocate(bt_cont%uBT_EE(isdb:iedb,jsd:jed), source=0.0)
599
600 allocate(bt_cont%FA_v_SS(isd:ied,jsdb:jedb), source=0.0)
601 allocate(bt_cont%FA_v_S0(isd:ied,jsdb:jedb), source=0.0)
602 allocate(bt_cont%FA_v_N0(isd:ied,jsdb:jedb), source=0.0)
603 allocate(bt_cont%FA_v_NN(isd:ied,jsdb:jedb), source=0.0)
604 allocate(bt_cont%vBT_SS(isd:ied,jsdb:jedb), source=0.0)
605 allocate(bt_cont%vBT_NN(isd:ied,jsdb:jedb), source=0.0)
606
607 if (present(alloc_faces)) then ; if (alloc_faces) then
608 allocate(bt_cont%h_u(isdb:iedb,jsd:jed,1:nz), source=0.0)
609 allocate(bt_cont%h_v(isd:ied,jsdb:jedb,1:nz), source=0.0)
610 endif ; endif
611
612end subroutine alloc_bt_cont_type
613
614!> Deallocates the arrays contained within a BT_cont_type.
615subroutine dealloc_bt_cont_type(BT_cont)
616 type(bt_cont_type), pointer :: bt_cont !< The BT_cont_type whose elements will be deallocated.
617
618 if (.not.associated(bt_cont)) return
619
620 deallocate(bt_cont%FA_u_WW) ; deallocate(bt_cont%FA_u_W0)
621 deallocate(bt_cont%FA_u_E0) ; deallocate(bt_cont%FA_u_EE)
622 deallocate(bt_cont%uBT_WW) ; deallocate(bt_cont%uBT_EE)
623
624 deallocate(bt_cont%FA_v_SS) ; deallocate(bt_cont%FA_v_S0)
625 deallocate(bt_cont%FA_v_N0) ; deallocate(bt_cont%FA_v_NN)
626 deallocate(bt_cont%vBT_SS) ; deallocate(bt_cont%vBT_NN)
627
628 if (allocated(bt_cont%h_u)) deallocate(bt_cont%h_u)
629 if (allocated(bt_cont%h_v)) deallocate(bt_cont%h_v)
630
631 deallocate(bt_cont)
632
633end subroutine dealloc_bt_cont_type
634
635!> Diagnostic checksums on various elements of a thermo_var_ptrs type for debugging.
636subroutine mom_thermovar_chksum(mesg, tv, G, US)
637 character(len=*), intent(in) :: mesg !< A message that appears in the checksum lines
638 type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
639 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
640 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
641
642 ! Note that for the chksum calls to be useful for reproducing across PE
643 ! counts, there must be no redundant points, so all variables use is..ie
644 ! and js...je as their extent.
645 if (associated(tv%T)) &
646 call hchksum(tv%T, mesg//" tv%T", g%HI, unscale=us%C_to_degC)
647 if (associated(tv%S)) &
648 call hchksum(tv%S, mesg//" tv%S", g%HI, unscale=us%S_to_ppt)
649 if (associated(tv%frazil)) &
650 call hchksum(tv%frazil, mesg//" tv%frazil", g%HI, unscale=us%Q_to_J_kg*us%RZ_to_kg_m2)
651 if (associated(tv%salt_deficit)) &
652 call hchksum(tv%salt_deficit, mesg//" tv%salt_deficit", g%HI, unscale=us%RZ_to_kg_m2*us%S_to_ppt)
653 if (associated(tv%TempxPmE)) &
654 call hchksum(tv%TempxPmE, mesg//" tv%TempxPmE", g%HI, unscale=us%RZ_to_kg_m2*us%C_to_degC)
655end subroutine mom_thermovar_chksum
656
657end module mom_variables