← back to index

src/core/MOM_variables.F90

portedportable, not yet portedexecuted, not portableexecutable, not hit by this run

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
17use MOM_unit_scaling, only : unit_scale_type
18use MOM_verticalGrid, only : verticalGrid_type
19use MOM_tracer_types, only : tracer_type
20
21implicit none ; private
22
23#include <MOM_memory.h>
24
25public allocate_surface_state, deallocate_surface_state, MOM_thermovar_chksum
26public ocean_grid_type, alloc_BT_cont_type, dealloc_BT_cont_type
27public rotate_surface_state
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]
164end type ocean_internal_state
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() !< Instantaneous active mixing layer depth [Z ~> m].
296 real, pointer, dimension(:,:) :: h_ML => NULL() !< Instantaneous active mixing layer thickness [H ~> m or kg m-2].
297 real, pointer, dimension(:,:) :: sfc_buoy_flx => NULL() !< Surface buoyancy flux (derived) [Z2 T-3 ~> m2 s-3].
298 real, pointer, dimension(:,:,:) :: Kd_shear => NULL()
299 !< The shear-driven turbulent diapycnal diffusivity at the interfaces between layers
300 !! in tracer columns [H Z T-1 ~> m2 s-1 or kg m-1 s-1]
301 real, pointer, dimension(:,:,:) :: Kv_shear => NULL()
302 !< The shear-driven turbulent vertical viscosity at the interfaces between layers
303 !! in tracer columns [H Z T-1 ~> m2 s-1 or Pa s]
304 real, pointer, dimension(:,:,:) :: Kv_shear_Bu => NULL()
305 !< The shear-driven turbulent vertical viscosity at the interfaces between layers in
306 !! corner columns [H Z T-1 ~> m2 s-1 or Pa s]
307 real, pointer, dimension(:,:,:) :: Kv_slow => NULL()
308 !< The turbulent vertical viscosity component due to "slow" processes (e.g., tidal,
309 !! background, convection etc) [H Z T-1 ~> m2 s-1 or Pa s]
310 real, pointer, dimension(:,:,:) :: TKE_turb => NULL()
311 !< The turbulent kinetic energy per unit mass at the interfaces [Z2 T-2 ~> m2 s-2].
312 !! This may be at the tracer or corner points
313end type vertvisc_type
314
315!> Container for information about the summed layer transports
316!! and how they will vary as the barotropic velocity is changed.
317type, public :: BT_cont_type
318 real, allocatable :: FA_u_EE(:,:) !< The effective open face area for zonal barotropic transport
319 !! drawing from locations far to the east [H L ~> m2 or kg m-1].
320 real, allocatable :: FA_u_E0(:,:) !< The effective open face area for zonal barotropic transport
321 !! drawing from nearby to the east [H L ~> m2 or kg m-1].
322 real, allocatable :: FA_u_W0(:,:) !< The effective open face area for zonal barotropic transport
323 !! drawing from nearby to the west [H L ~> m2 or kg m-1].
324 real, allocatable :: FA_u_WW(:,:) !< The effective open face area for zonal barotropic transport
325 !! drawing from locations far to the west [H L ~> m2 or kg m-1].
326 real, allocatable :: uBT_WW(:,:) !< uBT_WW is the barotropic velocity [L T-1 ~> m s-1], beyond which the
327 !! marginal open face area is FA_u_WW. uBT_WW must be non-negative.
328 real, allocatable :: uBT_EE(:,:) !< uBT_EE is a barotropic velocity [L T-1 ~> m s-1], beyond which the
329 !! marginal open face area is FA_u_EE. uBT_EE must be non-positive.
330 real, allocatable :: FA_v_NN(:,:) !< The effective open face area for meridional barotropic transport
331 !! drawing from locations far to the north [H L ~> m2 or kg m-1].
332 real, allocatable :: FA_v_N0(:,:) !< The effective open face area for meridional barotropic transport
333 !! drawing from nearby to the north [H L ~> m2 or kg m-1].
334 real, allocatable :: FA_v_S0(:,:) !< The effective open face area for meridional barotropic transport
335 !! drawing from nearby to the south [H L ~> m2 or kg m-1].
336 real, allocatable :: FA_v_SS(:,:) !< The effective open face area for meridional barotropic transport
337 !! drawing from locations far to the south [H L ~> m2 or kg m-1].
338 real, allocatable :: vBT_SS(:,:) !< vBT_SS is the barotropic velocity, [L T-1 ~> m s-1], beyond which the
339 !! marginal open face area is FA_v_SS. vBT_SS must be non-negative.
340 real, allocatable :: vBT_NN(:,:) !< vBT_NN is the barotropic velocity, [L T-1 ~> m s-1], beyond which the
341 !! marginal open face area is FA_v_NN. vBT_NN must be non-positive.
342 real, allocatable :: h_u(:,:,:) !< An effective thickness at zonal faces, taking into account the effects
343 !! of vertical viscosity and fractional open areas [H ~> m or kg m-2].
344 !! This is primarily used as a non-normalized weight in determining
345 !! the depth averaged accelerations for the barotropic solver.
346 real, allocatable :: h_v(:,:,:) !< An effective thickness at meridional faces, taking into account the effects
347 !! of vertical viscosity and fractional open areas [H ~> m or kg m-2].
348 !! This is primarily used as a non-normalized weight in determining
349 !! the depth averaged accelerations for the barotropic solver.
350 type(group_pass_type) :: pass_polarity_BT !< Structure for polarity group halo updates
351 type(group_pass_type) :: pass_FA_uv !< Structure for face area group halo updates
352end type BT_cont_type
353
354!> Container for grids modifying cell metric at porous barriers
355type, public :: porous_barrier_type
356 ! Each of the following fields has nz layers.
357 real, allocatable :: por_face_areaU(:,:,:) !< fractional open area of U-faces [nondim]
358 real, allocatable :: por_face_areaV(:,:,:) !< fractional open area of V-faces [nondim]
359 ! Each of the following fields is found at nz+1 interfaces.
360 real, allocatable :: por_layer_widthU(:,:,:) !< fractional open width of U-faces [nondim]
361 real, allocatable :: por_layer_widthV(:,:,:) !< fractional open width of V-faces [nondim]
362end type porous_barrier_type
363
364contains
365
366!> Allocates the fields for the surface (return) properties of
367!! the ocean model. Unused fields are unallocated.
3681subroutine allocate_surface_state(sfc_state, G, use_temperature, do_integrals, &
369 gas_fields_ocn, use_meltpot, use_iceshelves, &
370 omit_frazil, sfc_state_in, turns, use_marbl_tracers)
371 type(ocean_grid_type), intent(in) :: G !< ocean grid structure
372 type(surface), intent(inout) :: sfc_state !< ocean surface state type to be allocated.
373 logical, optional, intent(in) :: use_temperature !< If true, allocate the space for thermodynamic variables.
374 logical, optional, intent(in) :: do_integrals !< If true, allocate the space for vertically
375 !! integrated fields.
376 type(coupler_1d_bc_type), &
377 optional, intent(in) :: gas_fields_ocn !< If present, this type describes the
378 !! ocean and surface-ice fields that will participate
379 !! in the calculation of additional gas or other
380 !! tracer fluxes, and can be used to spawn related
381 !! internal variables in the ice model.
382 logical, optional, intent(in) :: use_meltpot !< If true, allocate the space for melt potential
383 logical, optional, intent(in) :: use_iceshelves !< If true, allocate the space for the stresses
384 !! under ice shelves.
385 logical, optional, intent(in) :: omit_frazil !< If present and false, do not allocate the space to
386 !! pass frazil fluxes to the coupler
387 type(surface), &
388 optional, intent(in) :: sfc_state_in !< If present and its tr_fields are initialized,
389 !! this type describes the ocean and surface-ice fields that
390 !! will participate in the calculation of additional gas or
391 !! other tracer fluxes, and can be used to spawn related
392 !! internal variables in the ice model. If gas_fields_ocn
393 !! is present, it is used and tr_fields_in is ignored.
394 integer, optional, intent(in) :: turns !< If present, the number of counterclockwise quarter
395 !! turns to use on the new grid.
396 logical, optional, intent(in) :: use_marbl_tracers !< If true, allocate the space for CO2 flux from MARBL
397
398 ! local variables
399 logical :: use_temp, alloc_integ, use_melt_potential, alloc_iceshelves, alloc_frazil, alloc_fco2
400 logical :: even_turns ! True if turns is absent or even
401 integer :: tr_field_i_mem(4), tr_field_j_mem(4)
402 integer :: is, ie, js, je, isd, ied, jsd, jed
403 integer :: isdB, iedB, jsdB, jedB
404
4051 is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec
4061 isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed
4071 isdB = G%isdB ; iedB = G%iedB ; jsdB = G%jsdB ; jedB = G%jedB
408
4091 use_temp = .true. ; if (present(use_temperature)) use_temp = use_temperature
4101 alloc_integ = .true. ; if (present(do_integrals)) alloc_integ = do_integrals
4111 use_melt_potential = .false. ; if (present(use_meltpot)) use_melt_potential = use_meltpot
4121 alloc_iceshelves = .false. ; if (present(use_iceshelves)) alloc_iceshelves = use_iceshelves
4131 alloc_frazil = .true. ; if (present(omit_frazil)) alloc_frazil = .not.omit_frazil
4141 alloc_fco2 = .false. ; if (present(use_marbl_tracers)) alloc_fco2 = use_marbl_tracers
415
4161 if (sfc_state%arrays_allocated) return
417
4181 if (use_temp) then
4198773 allocate(sfc_state%SST(isd:ied,jsd:jed), source=0.0)
4208773 allocate(sfc_state%SSS(isd:ied,jsd:jed), source=0.0)
421 else
4220 allocate(sfc_state%sfc_density(isd:ied,jsd:jed), source=0.0)
423 endif
4241 if (use_temp .and. alloc_frazil) then
4258773 allocate(sfc_state%frazil(isd:ied,jsd:jed), source=0.0)
426 endif
4278773 allocate(sfc_state%sea_lev(isd:ied,jsd:jed), source=0.0)
4288773 allocate(sfc_state%Hml(isd:ied,jsd:jed), source=0.0)
4298841 allocate(sfc_state%u(IsdB:IedB,jsd:jed), source=0.0)
4308902 allocate(sfc_state%v(isd:ied,JsdB:JedB), source=0.0)
431
4321 if (use_melt_potential) then
4330 allocate(sfc_state%melt_potential(isd:ied,jsd:jed), source=0.0)
434 endif
435
4361 if (alloc_integ) then
437 ! Allocate structures for the vertically integrated ocean_mass, ocean_heat, and ocean_salt.
4388773 allocate(sfc_state%ocean_mass(isd:ied,jsd:jed), source=0.0)
4391 if (use_temp) then
4408773 allocate(sfc_state%ocean_heat(isd:ied,jsd:jed), source=0.0)
4418773 allocate(sfc_state%ocean_salt(isd:ied,jsd:jed), source=0.0)
442 endif
443 endif
444
4451 if (alloc_iceshelves) then
4460 allocate(sfc_state%taux_shelf(IsdB:IedB,jsd:jed), source=0.0)
4470 allocate(sfc_state%tauy_shelf(isd:ied,JsdB:JedB), source=0.0)
448 endif
449
450 ! The data fields in the coupler_2d_bc_type are never rotated.
4511 even_turns = .true. ; if (present(turns)) even_turns = (modulo(turns, 2) == 0)
4521 if (even_turns) then
4539 tr_field_i_mem(1:4) = (/is,is,ie,ie/) ; tr_field_j_mem(1:4) = (/js,js,je,je/)
454 else
4550 tr_field_i_mem(1:4) = (/js,js,je,je/) ; tr_field_j_mem(1:4) = (/is,is,ie,ie/)
456 endif
4571 if (present(gas_fields_ocn)) then
458 call coupler_type_spawn(gas_fields_ocn, sfc_state%tr_fields, &
4590 tr_field_i_mem, tr_field_j_mem, as_needed=.true.)
4601 elseif (present(sfc_state_in)) then
4610 if (coupler_type_initialized(sfc_state_in%tr_fields)) then
462 call coupler_type_spawn(sfc_state_in%tr_fields, sfc_state%tr_fields, &
4630 tr_field_i_mem, tr_field_j_mem, as_needed=.true.)
464 endif
465 endif
466
4671 if (alloc_fco2) then
4680 allocate(sfc_state%fco2(isd:ied,jsd:jed), source=0.0)
469 endif
470
4711 sfc_state%arrays_allocated = .true.
472
473end subroutine allocate_surface_state
474
475!> Deallocates the elements of a surface state type.
4760subroutine deallocate_surface_state(sfc_state)
477 type(surface), intent(inout) :: sfc_state !< ocean surface state type to be deallocated here.
478
4790 if (.not.sfc_state%arrays_allocated) return
480
4810 if (allocated(sfc_state%melt_potential)) deallocate(sfc_state%melt_potential)
4820 if (allocated(sfc_state%SST)) deallocate(sfc_state%SST)
4830 if (allocated(sfc_state%SSS)) deallocate(sfc_state%SSS)
4840 if (allocated(sfc_state%sfc_density)) deallocate(sfc_state%sfc_density)
4850 if (allocated(sfc_state%sea_lev)) deallocate(sfc_state%sea_lev)
4860 if (allocated(sfc_state%Hml)) deallocate(sfc_state%Hml)
4870 if (allocated(sfc_state%u)) deallocate(sfc_state%u)
4880 if (allocated(sfc_state%v)) deallocate(sfc_state%v)
4890 if (allocated(sfc_state%ocean_mass)) deallocate(sfc_state%ocean_mass)
4900 if (allocated(sfc_state%ocean_heat)) deallocate(sfc_state%ocean_heat)
4910 if (allocated(sfc_state%ocean_salt)) deallocate(sfc_state%ocean_salt)
4920 if (allocated(sfc_state%fco2)) deallocate(sfc_state%fco2)
4930 call coupler_type_destructor(sfc_state%tr_fields)
494
4950 sfc_state%arrays_allocated = .false.
496
497end subroutine deallocate_surface_state
498
499!> Rotate the surface state fields from the input to the model indices.
5000subroutine rotate_surface_state(sfc_state_in, sfc_state, G, turns)
501 type(surface), intent(in) :: sfc_state_in !< The input unrotated surface state type that is the data source.
502 type(surface), intent(inout) :: sfc_state !< The rotated surface state type whose arrays will be filled in
503 type(ocean_grid_type), intent(in) :: G !< The ocean grid structure
504 integer, intent(in) :: turns !< The number of counterclockwise quarter turns to use on the rotated grid.
505
506 logical :: use_temperature, do_integrals, use_melt_potential, use_iceshelves
507
508 ! NOTE: Many of these are weak tests, since only one is checked
509 use_temperature = allocated(sfc_state_in%SST) &
5100 .and. allocated(sfc_state_in%SSS)
5110 use_melt_potential = allocated(sfc_state_in%melt_potential)
5120 do_integrals = allocated(sfc_state_in%ocean_mass)
513 use_iceshelves = allocated(sfc_state_in%taux_shelf) &
5140 .and. allocated(sfc_state_in%tauy_shelf)
515
5160 if (.not. sfc_state%arrays_allocated) then
517 call allocate_surface_state(sfc_state, G, use_temperature=use_temperature, &
518 do_integrals=do_integrals, use_meltpot=use_melt_potential, &
5190 use_iceshelves=use_iceshelves, sfc_state_in=sfc_state_in, turns=turns)
520 endif
521
5220 if (use_temperature) then
5230 call rotate_array(sfc_state_in%SST, turns, sfc_state%SST)
5240 call rotate_array(sfc_state_in%SSS, turns, sfc_state%SSS)
525 else
5260 call rotate_array(sfc_state_in%sfc_density, turns, sfc_state%sfc_density)
527 endif
528
5290 call rotate_array(sfc_state_in%Hml, turns, sfc_state%Hml)
530 call rotate_vector(sfc_state_in%u, sfc_state_in%v, turns, &
5310 sfc_state%u, sfc_state%v)
5320 call rotate_array(sfc_state_in%sea_lev, turns, sfc_state%sea_lev)
533
5340 if (use_melt_potential) then
5350 call rotate_array(sfc_state_in%melt_potential, turns, sfc_state%melt_potential)
536 endif
537
5380 if (do_integrals) then
5390 call rotate_array(sfc_state_in%ocean_mass, turns, sfc_state%ocean_mass)
5400 if (use_temperature) then
5410 call rotate_array(sfc_state_in%ocean_heat, turns, sfc_state%ocean_heat)
5420 call rotate_array(sfc_state_in%ocean_salt, turns, sfc_state%ocean_salt)
5430 call rotate_array(sfc_state_in%SSS, turns, sfc_state%SSS)
544 endif
545 endif
546
5470 if (use_iceshelves) then
548 call rotate_vector(sfc_state_in%taux_shelf, sfc_state_in%tauy_shelf, turns, &
5490 sfc_state%taux_shelf, sfc_state%tauy_shelf)
550 endif
551
5520 if (use_temperature .and. allocated(sfc_state_in%frazil)) &
5530 call rotate_array(sfc_state_in%frazil, turns, sfc_state%frazil)
554
555 ! Scalar transfers
5560 sfc_state%T_is_conT = sfc_state_in%T_is_conT
5570 sfc_state%S_is_absS = sfc_state_in%S_is_absS
558
559 ! NOTE: Tracer fields are handled by FMS, so are left unrotated. Any
560 ! reads/writes to tr_fields must be appropriately rotated.
5610 if (coupler_type_initialized(sfc_state_in%tr_fields)) then
5620 call coupler_type_copy_data(sfc_state_in%tr_fields, sfc_state%tr_fields)
563 endif
5640end subroutine rotate_surface_state
565
566!> Allocates the arrays contained within a BT_cont_type and initializes them to 0.
5671subroutine alloc_BT_cont_type(BT_cont, G, GV, alloc_faces)
568 type(BT_cont_type), pointer :: BT_cont !< The BT_cont_type whose elements will be allocated
569 type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
570 type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
571 logical, optional, intent(in) :: alloc_faces !< If present and true, allocate
572 !! memory for effective face thicknesses.
573
574 integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, nz
5751 isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed ; nz = GV%ke
5761 IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB
577
5781 if (associated(BT_cont)) call MOM_error(FATAL, &
5790 "alloc_BT_cont_type called with an associated BT_cont_type pointer.")
580
5811 allocate(BT_cont)
582 !$omp target enter data map(to: BT_cont)
5838841 allocate(BT_cont%FA_u_WW(IsdB:IedB,jsd:jed), source=0.0)
5848841 allocate(BT_cont%FA_u_W0(IsdB:IedB,jsd:jed), source=0.0)
5858841 allocate(BT_cont%FA_u_E0(IsdB:IedB,jsd:jed), source=0.0)
5868841 allocate(BT_cont%FA_u_EE(IsdB:IedB,jsd:jed), source=0.0)
5878841 allocate(BT_cont%uBT_WW(IsdB:IedB,jsd:jed), source=0.0)
5888841 allocate(BT_cont%uBT_EE(IsdB:IedB,jsd:jed), source=0.0)
589 !$omp target enter data map(to: BT_cont%FA_u_WW, BT_cont%FA_u_W0)
590 !$omp target enter data map(to: BT_cont%FA_u_E0, BT_cont%FA_u_EE)
591 !$omp target enter data map(to: BT_cont%uBT_WW, BT_cont%uBT_EE)
592
5938902 allocate(BT_cont%FA_v_SS(isd:ied,JsdB:JedB), source=0.0)
5948902 allocate(BT_cont%FA_v_S0(isd:ied,JsdB:JedB), source=0.0)
5958902 allocate(BT_cont%FA_v_N0(isd:ied,JsdB:JedB), source=0.0)
5968902 allocate(BT_cont%FA_v_NN(isd:ied,JsdB:JedB), source=0.0)
5978902 allocate(BT_cont%vBT_SS(isd:ied,JsdB:JedB), source=0.0)
5988902 allocate(BT_cont%vBT_NN(isd:ied,JsdB:JedB), source=0.0)
599 !$omp target enter data map(to: BT_cont%FA_v_SS, BT_cont%FA_v_S0)
600 !$omp target enter data map(to: BT_cont%FA_v_N0, BT_cont%FA_v_NN)
601 !$omp target enter data map(to: BT_cont%vBT_SS, BT_cont%vBT_NN)
602
6031 if (present(alloc_faces)) then ; if (alloc_faces) then
604663076 allocate(BT_cont%h_u(IsdB:IedB,jsd:jed,1:nz), source=0.0)
605667651 allocate(BT_cont%h_v(isd:ied,JsdB:JedB,1:nz), source=0.0)
606 !$omp target enter data map(to: BT_cont%h_u, BT_cont%h_v)
607 endif ; endif
608
6091end subroutine alloc_BT_cont_type
610
611!> Deallocates the arrays contained within a BT_cont_type.
6121subroutine dealloc_BT_cont_type(BT_cont)
613 type(BT_cont_type), pointer :: BT_cont !< The BT_cont_type whose elements will be deallocated.
614
6151 if (.not.associated(BT_cont)) return
616
617 !$omp target exit data map(delete: BT_cont%FA_u_EE, BT_cont%FA_u_E0)
618 !$omp target exit data map(delete: BT_cont%FA_u_W0, BT_cont%FA_u_WW)
619 !$omp target exit data map(delete: BT_cont%uBT_WW, BT_cont%uBT_EE)
6201 deallocate(BT_cont%FA_u_EE) ; deallocate(BT_cont%FA_u_E0)
6211 deallocate(BT_cont%FA_u_W0) ; deallocate(BT_cont%FA_u_WW)
6221 deallocate(BT_cont%uBT_WW) ; deallocate(BT_cont%uBT_EE)
623
624 !$omp target exit data map(delete: BT_cont%FA_v_NN, BT_cont%FA_v_N0)
625 !$omp target exit data map(delete: BT_cont%FA_v_S0, BT_cont%FA_v_SS)
626 !$omp target exit data map(delete: BT_cont%vBT_SS, BT_cont%vBT_NN)
6271 deallocate(BT_cont%FA_v_NN) ; deallocate(BT_cont%FA_v_N0)
6281 deallocate(BT_cont%FA_v_S0) ; deallocate(BT_cont%FA_v_SS)
6291 deallocate(BT_cont%vBT_SS) ; deallocate(BT_cont%vBT_NN)
630
631 ! These are always allocated in pairs.
6321 if (allocated(BT_cont%h_u) .and. allocated(BT_cont%h_v)) then
633 !$omp target exit data map(delete: BT_cont%h_u, BT_cont%h_v)
6341 deallocate(BT_cont%h_u, BT_cont%h_v)
635 endif
636
637 !$omp target exit data map(delete: BT_cont)
6381 deallocate(BT_cont)
639end subroutine dealloc_BT_cont_type
640
641!> Diagnostic checksums on various elements of a thermo_var_ptrs type for debugging.
6420subroutine MOM_thermovar_chksum(mesg, tv, G, US)
643 character(len=*), intent(in) :: mesg !< A message that appears in the checksum lines
644 type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
645 type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
646 type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
647
648 ! Note that for the chksum calls to be useful for reproducing across PE
649 ! counts, there must be no redundant points, so all variables use is..ie
650 ! and js...je as their extent.
6510 if (associated(tv%T)) &
6520 call hchksum(tv%T, mesg//" tv%T", G%HI, unscale=US%C_to_degC)
6530 if (associated(tv%S)) &
6540 call hchksum(tv%S, mesg//" tv%S", G%HI, unscale=US%S_to_ppt)
6550 if (associated(tv%frazil)) &
6560 call hchksum(tv%frazil, mesg//" tv%frazil", G%HI, unscale=US%Q_to_J_kg*US%RZ_to_kg_m2)
6570 if (associated(tv%salt_deficit)) &
6580 call hchksum(tv%salt_deficit, mesg//" tv%salt_deficit", G%HI, unscale=US%RZ_to_kg_m2*US%S_to_ppt)
6590 if (associated(tv%TempxPmE)) &
6600 call hchksum(tv%TempxPmE, mesg//" tv%TempxPmE", G%HI, unscale=US%RZ_to_kg_m2*US%C_to_degC)
6610end subroutine MOM_thermovar_chksum
662
6630end module MOM_variables