MOM_check_scaling.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!> This module is used to check the dimensional scaling factors used by the MOM6 ocean model
7
12
13implicit none ; private
14
15! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
16! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
17! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
18! vary with the Boussinesq approximation, the Boussinesq variant is given first.
19
21
22contains
23
24!> Evaluate whether the dimensional scaling factors provide unique tests for all of the combinations
25!! of dimensions that are used in MOM6 (or perhaps widely used), and if they are not unique, explore
26!! whether another combination of scaling factors can be found that is unique or has less common
27!! cases with coinciding scaling.
28subroutine check_mom6_scaling_factors(GV, US)
29 type(verticalgrid_type), pointer :: gv !< The container for vertical grid data
30 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
31
32 ! Local variables
33 integer, parameter :: ndims = 8 ! The number of rescalable dimensional factors.
34 real, dimension(ndims) :: scales ! An array of scaling factors for each of the basic units [various].
35 integer, dimension(ndims) :: scale_pow2 ! The powers of 2 that give each element of scales.
36 character(len=2), dimension(ndims) :: key
37 integer, allocatable :: weights(:)
38 character(len=80), allocatable :: descriptions(:)
39 integer :: n, ns, max_pow
40
41 ! If no scaling is being done, simply return.
42 if ((us%Z_to_m == 1.) .and. (gv%H_to_MKS == 1.) .and. (us%L_to_m == 1.) .and. &
43 (us%T_to_s == 1.) .and. (us%R_to_kg_m3 == 1.) .and. (us%Q_to_J_kg == 1.) .and. &
44 (us%C_to_degC == 1.) .and. (us%S_to_ppt == 1.)) return
45
46 ! Set the names and scaling factors of the dimensions being rescaled.
47 key(:) = ["Z", "H", "L", "T", "R", "Q", "C", "S"]
48 scales(:) = (/ us%Z_to_m, gv%H_to_MKS, us%L_to_m, us%T_to_s, us%R_to_kg_m3, us%Q_to_J_kg, &
49 us%C_to_degC, us%S_to_ppt/)
50 call scales_to_powers(scales, scale_pow2)
51 max_pow = 40 ! 60
52
53 ! The first call is just to find out how many elements are in the list of scaling combinations.
54 call compose_dimension_list(ns, descriptions, weights)
55
56 allocate(descriptions(ns))
57 do n=1,ns ; descriptions(n) = "" ; enddo
58 allocate(weights(ns), source=0)
59 ! This call records all the list of powers, the descriptions, and their weights.
60 call compose_dimension_list(ns, descriptions, weights)
61
62 call check_scaling_uniqueness("MOM6", descriptions, weights, key, scale_pow2, max_pow)
63
64 deallocate(weights)
65 deallocate(descriptions)
66
67end subroutine check_mom6_scaling_factors
68
69
70!> This routine composes a list of the commonly used dimensional scaling factors in the MOM6
71!! code, along with weights reflecting the frequency of their occurrence in the MOM6 code or
72!! other considerations of how likely the variables are be used.
73subroutine compose_dimension_list(ns, des, wts)
74 integer, intent(out) :: ns !< The running sum of valid descriptions
75 character(len=*), allocatable, intent(inout) :: des(:) !< The unit descriptions that have been converted
76 integer, allocatable, intent(inout) :: wts(:) !< A list of the integer weights for each scaling factor,
77 !! perhaps the number of times it occurs in the MOM6 code.
78
79 ns = 0
80 ! Accumulate a list of units used in MOM6, in approximate descending order of frequency of occurrence in
81 ! doxygen comments (i.e., arguments and elements in types), excluding the code in the user, ice_shelf and
82 ! framework directories and the passive tracer packages.
83 call add_scaling(ns, des, wts, "[H ~> m or kg m-2]", 716) ! Layer thicknesses
84 call add_scaling(ns, des, wts, "[L T-1 ~> m s-1]", 264) ! Horizontal velocities
85 call add_scaling(ns, des, wts, "[Z ~> m]", 244) ! Depths and vertical distance
86 call add_scaling(ns, des, wts, "[T ~> s]", 154) ! Time intervals
87 call add_scaling(ns, des, wts, "[S ~> ppt]", 135) ! Salinities
88 call add_scaling(ns, des, wts, "[C ~> degC]", 135) ! Temperatures
89 call add_scaling(ns, des, wts, "[R L2 T-2 ~> Pa]", 133) ! Dynamic pressure
90 ! call add_scaling(ns, des, wts, "[R L2 T-2 ~> J m-3]") ! Energy density
91 call add_scaling(ns, des, wts, "[Z2 T-1 ~> m2 s-1]", 132) ! Vertical viscosities and diffusivities
92 call add_scaling(ns, des, wts, "[R ~> kg m-3]", 122) ! Densities
93
94 call add_scaling(ns, des, wts, "[H L2 T-1 ~> m3 s-1 or kg s-1]", 97) ! Volume or mass transports
95 call add_scaling(ns, des, wts, "[H L2 ~> m3 or kg]", 91) ! Cell volumes or masses
96 call add_scaling(ns, des, wts, "[L T-2 ~> m s-2]", 82) ! Horizontal accelerations
97 call add_scaling(ns, des, wts, "[T-1 ~> s-1]", 67) ! Rates
98 call add_scaling(ns, des, wts, "[Z T-1 ~> m s-1]", 56) ! Friction velocities and viscous coupling
99 call add_scaling(ns, des, wts, "[Q R Z T-1 ~> W m-2]", 42) ! Vertical heat fluxes
100 call add_scaling(ns, des, wts, "[L2 T-1 ~> m2 s-1]", 45) ! Horizontal viscosity or diffusivity
101 call add_scaling(ns, des, wts, "[L2 T-2 ~> m2 s-2]", 37) ! Resolved kinetic energy per unit mass
102 call add_scaling(ns, des, wts, "[L ~> m]", 35) ! Horizontal distances
103 call add_scaling(ns, des, wts, "[T-2 ~> s-2]", 33) ! Squared shears and buoyancy frequency
104
105 call add_scaling(ns, des, wts, "[R Z L T-2 ~> Pa]", 33) ! Wind stresses
106 call add_scaling(ns, des, wts, "[H L ~> m2 or kg m-1]", 32) ! Lateral cell face areas
107 call add_scaling(ns, des, wts, "[L2 ~> m2]", 31) ! Horizontal areas
108 call add_scaling(ns, des, wts, "[R C-1 ~> kg m-3 degC-1]", 26) ! Thermal expansion coefficients
109 call add_scaling(ns, des, wts, "[L2 Z-1 T-2 ~> m s-2]", 26) ! Gravitational acceleration
110 call add_scaling(ns, des, wts, "[R S-1 ~> kg m-3 ppt-1]", 23) ! Haline contraction coefficients
111 call add_scaling(ns, des, wts, "[R Z3 T-3 ~> W m-2]", 23) ! Integrated turbulent kinetic energy sources
112 call add_scaling(ns, des, wts, "[R Z T-1 ~> kg m-2 s-1]", 19) ! Vertical mass fluxes
113 call add_scaling(ns, des, wts, "[C H ~> degC m or degC kg m-2]", 17) ! Heat content
114 call add_scaling(ns, des, wts, "[H-1 ~> m-1 or m2 kg-1]", 17) ! Inverse cell thicknesses
115
116 call add_scaling(ns, des, wts, "[Z-1 ~> m-1]", 14) ! Inverse vertical distances
117 call add_scaling(ns, des, wts, "[R-1 ~> m3 kg-1]", 14) ! Specific volumes
118 call add_scaling(ns, des, wts, "[Z L-1 ~> nondim]", 12) ! Slopes
119 call add_scaling(ns, des, wts, "[L-1 ~> m-1]", 12) ! Inverse horizontal distances
120 call add_scaling(ns, des, wts, "[L2 T-2 H-1 ~> m s-2 or m4 s-2 kg-1]", 12) ! pbce or gtot
121 call add_scaling(ns, des, wts, "[R Z ~> kg m-2]", 11) ! Layer or column mass loads
122 call add_scaling(ns, des, wts, "[Z L2 T-2 ~> m3 s-2]", 11) ! Integrated energy per unit mass
123 call add_scaling(ns, des, wts, "[R Z3 T-2 ~> J m-2]", 11) ! Integrated turbulent kinetic energy density
124 call add_scaling(ns, des, wts, "[H T-1 ~> m s-1 or kg m-2 s-1]", 9) ! Vertical thickness fluxes
125 call add_scaling(ns, des, wts, "[L-1 T-1 ~> m-1 s-1]", 9) ! Laplacian of velocity
126
127 call add_scaling(ns, des, wts, "[Z3 T-3 ~> m3 s-3]", 9) ! Integrated turbulent kinetic energy sources
128 call add_scaling(ns, des, wts, "[S H ~> ppt m or ppt kg m-2]", 8) ! Depth integrated salinity
129 call add_scaling(ns, des, wts, "[Z2 T-2 ~> m2 s-2]", 8) ! Turbulent kinetic energy
130 call add_scaling(ns, des, wts, "[R L2 Z T-2 ~> Pa m]", 7) ! Vertically integrated pressure anomalies
131 call add_scaling(ns, des, wts, "[T2 Z-1 ~> s2 m-1]", 7) ! (TKE_to_Kd)
132 call add_scaling(ns, des, wts, "[L4 T-1 ~> m4 s-1]", 7) ! Biharmonic viscosity
133 call add_scaling(ns, des, wts, "[L3 ~> m3]", 7) ! Metric dependent constants for viscosity
134 call add_scaling(ns, des, wts, "[L2 T-3 ~> m2 s-3]", 7) ! Buoyancy flux or MEKE sources [L2 T-3 ~> W kg-1]
135 call add_scaling(ns, des, wts, "[H2 ~> m2 or kg2 m-4]", 7) ! Squared layer thicknesses
136 call add_scaling(ns, des, wts, "[C H T-1 ~> degC m s-1 or degC kg m-2 s-1]", 7) ! vertical heat fluxes
137
138 call add_scaling(ns, des, wts, "[L-2 ~> m-2]", 6) ! Inverse areas
139 call add_scaling(ns, des, wts, "[R Z L2 T-3 ~> W m-2]", 6) ! Energy sources, including for MEKE
140 call add_scaling(ns, des, wts, "[Z2 T-3 ~> m2 s-3]", 5) ! Certain buoyancy fluxes
141 call add_scaling(ns, des, wts, "[Z2 ~> m2]", 5) ! Squared vertical distances
142 call add_scaling(ns, des, wts, "[S H T-1 ~> ppt m s-1 or ppt kg m-2 s-1]", 5) ! vertical salinity fluxes
143 call add_scaling(ns, des, wts, "[R-1 C-1 ~> m3 kg-1 degC-1]", 5) ! Specific volume temperature gradient
144 call add_scaling(ns, des, wts, "[R-1 S-1 ~> m3 kg-1 ppt-1]", 4) ! Specific volume salnity gradient
145 call add_scaling(ns, des, wts, "[Q R Z ~> J m-2]", 4) ! time-integrated frazil heat flux
146 call add_scaling(ns, des, wts, "[Z C-1 ~> m degC-1]", 4) ! Inverse temperature gradients
147 call add_scaling(ns, des, wts, "[Z S-1 ~> m ppt-1]", 4) ! Inverse salinity gradients
148
149 call add_scaling(ns, des, wts, "[R Z3 T-2 H-1 ~> J m-3 or J kg-1]", 4) ! Partial derivatives of energy
150 call add_scaling(ns, des, wts, "[R Z3 T-2 S-1 ~> J m-2 ppt-1]", 4) ! Sensitity of energy change to salinity
151 call add_scaling(ns, des, wts, "[R Z3 T-2 C-1 ~> J m-2 degC-1]", 4) ! Sensitity of energy change to temperature
152 call add_scaling(ns, des, wts, "[R L4 T-4 ~> Pa m2 s-2]", 4) ! Integral in geopotential of pressure
153 call add_scaling(ns, des, wts, "[Q ~> J kg-1]", 4) ! Latent heats
154 call add_scaling(ns, des, wts, "[Q C-1 ~> J kg-1 degC-1]", 4) ! Heat capacity
155 call add_scaling(ns, des, wts, "[L-3 ~> m-3]", 4) ! Metric dependent constants for viscosity
156 call add_scaling(ns, des, wts, "[L2 Z-2 T-2 ~> s-2]", 4) ! Buoyancy frequency in some params.
157 call add_scaling(ns, des, wts, "[H R ~> kg m-2 or kg2 m-5]", 4) ! Layer-integrated density
158 call add_scaling(ns, des, wts, "[H L T-1 ~> m2 s-1 or kg m-1 s-1]", 4) ! Layer integrated velocities
159
160 call add_scaling(ns, des, wts, "[H T2 L-1 ~> s2 or kg s2 m-3]", 4) ! BT_cont_type face curvature fit
161 call add_scaling(ns, des, wts, "[H L-1 ~> nondim or kg m-3]", 4) ! BT_cont_type face curvature fit
162 call add_scaling(ns, des, wts, "[C2 ~> degC2]", 4) ! Squared temperature anomalies
163 call add_scaling(ns, des, wts, "[S2 ~> ppt2]", 3) ! Squared salinity anomalies
164 call add_scaling(ns, des, wts, "[C S ~> degC ppt]", 3) ! Covariance of temperature and salinity anomalies
165 call add_scaling(ns, des, wts, "[S R Z ~> gSalt m-2]", 3) ! Total ocean column salt
166 call add_scaling(ns, des, wts, "[C R Z ~> degC kg m-2]", 3) ! Total ocean column temperature
167 call add_scaling(ns, des, wts, "[Pa T2 R-1 L-2 ~> 1]", 3) ! Pressure conversions
168 call add_scaling(ns, des, wts, "[Z H-1 ~> nondim or m3 kg-1]", 3) ! Thickness to height conversion
169 call add_scaling(ns, des, wts, "[R Z2 T-2 ~> J m-3]", 3) ! Potential energy height derivatives
170
171 call add_scaling(ns, des, wts, "[H-2 ~> m-2 or m4 kg-2]", 3) ! Mixed layer local work variables
172 call add_scaling(ns, des, wts, "[C S-1 ~> degC ppt-1]", 2) ! T / S gauge transformation
173 call add_scaling(ns, des, wts, "[R S-2 ~> kg m-3 ppt-2]", 2) ! Second derivative of density
174 call add_scaling(ns, des, wts, "[R C-2 ~> kg m-3 degC-2]", 2) ! Second derivative of density
175 call add_scaling(ns, des, wts, "[R S-1 C-1 ~> kg m-3 ppt-1 degC-1]", 2) ! Second derivative of density
176 call add_scaling(ns, des, wts, "[T2 S-1 L-2 ~> kg m-3 ppt-1 Pa-1]", 2) ! Second derivative of density
177 call add_scaling(ns, des, wts, "[T2 C-1 L-2 ~> kg m-3 degC-1 Pa-1]", 2) ! Second derivative of density
178 call add_scaling(ns, des, wts, "[T2 L-2 ~> s2 m-2]", 2) ! Inverse velocities squared
179 call add_scaling(ns, des, wts, "[R Z2 T-3 ~> W m-3]", 2) ! Kinetic energy dissipation rates
180 call add_scaling(ns, des, wts, "[R H-1 ~> kg m-4 or m-1]", 2) ! Vertical density gradients
181
182 call add_scaling(ns, des, wts, "[L4 ~> m4]", 2) ! Metric dependent constants for viscosity
183 call add_scaling(ns, des, wts, "[Z L T-1 ~> m2 s-1]", 2) ! Layer integrated velocities
184 call add_scaling(ns, des, wts, "[C Z ~> degC m]", 2) ! Depth integrated temperature
185 call add_scaling(ns, des, wts, "[S Z ~> ppt m]", 1) ! Layer integrated salinity
186 call add_scaling(ns, des, wts, "[T L4 ~> s m4]", 2) ! Biharmonic metric dependent constant
187 call add_scaling(ns, des, wts, "[L6 ~> m6]", 2) ! Biharmonic Leith metric dependent constant
188 call add_scaling(ns, des, wts, "[L4 Z-1 T-1 ~> m3 s-1]", 2) ! Rigidity of ice
189 call add_scaling(ns, des, wts, "[L4 Z-2 T-1 ~> m2 s-1]", 1) ! Ice rigidity term
190 call add_scaling(ns, des, wts, "[R-1 Z-1 ~> m2 kg-1]", 1) ! Inverse of column mass
191 call add_scaling(ns, des, wts, "[Z-2 ~> m-2]", 1) ! Inverse of denominator in some weighted averages
192
193 call add_scaling(ns, des, wts, "[R Z2 T-1 ~> J s m-3]", 1) ! River mixing term
194 call add_scaling(ns, des, wts, "[R L2 T-2 H-1 ~> Pa m-1 or Pa m2 kg-1]", 1) ! Thickness to pressure conversion
195 call add_scaling(ns, des, wts, "[Z T2 R-1 L-2 ~> m Pa-1]", 1) ! Atmospheric pressure SSH correction
196 call add_scaling(ns, des, wts, "[T Z ~> s m] ", 1) ! Time integrated SSH
197 call add_scaling(ns, des, wts, "[Z-1 T-1 ~> m-1 s-1]", 1) ! barotropic PV
198 call add_scaling(ns, des, wts, "[L2 T ~> m2 s]", 1) ! Greatbatch & Lamb 90 coefficient
199 call add_scaling(ns, des, wts, "[Z L2 T-1 ~> m3 s-1]", 1) ! Overturning (GM) streamfunction
200 call add_scaling(ns, des, wts, "[kg H-1 L-2 ~> kg m-3 or 1]", 1) ! Diagnostic conversions to mass
201 call add_scaling(ns, des, wts, "[S-1 ~> ppt-1]", 1) ! Unscaling salinity
202 call add_scaling(ns, des, wts, "[C-1 ~> degC-1]", 1) ! Unscaling temperature
203
204 call add_scaling(ns, des, wts, "[R Z H-1 ~> kg m-3 or 1] ", 1) ! A unit conversion factor
205 call add_scaling(ns, des, wts, "[H R-1 Z-1 ~> m3 kg-2 or 1]", 1) ! A unit conversion factor
206 call add_scaling(ns, des, wts, "[H Z-1 ~> 1 or kg m-3]", 1) ! A unit conversion factor
207 call add_scaling(ns, des, wts, "[m T s-1 L-1 ~> 1]", 1) ! A unit conversion factor
208
209end subroutine compose_dimension_list
210
211!> Augment the count the valid unit descriptions, and add the provided description and its weight
212!! to the end of the list if that list is allocated.
213subroutine add_scaling(ns, descs, wts, scaling, weight)
214 integer, intent(inout) :: ns !< The running sum of valid descriptions.
215 character(len=*), allocatable, intent(inout) :: descs(:) !< The unit descriptions that have been converted
216 integer, allocatable, intent(inout) :: wts(:) !< A list of the integers for each scaling
217 character(len=*), intent(in) :: scaling !< The unit description that will be converted
218 integer, optional, intent(in) :: weight !< An optional weight or occurrence count
219 !! for this unit description, 1 by default.
220
221 integer :: iend
222
223 iend = index(scaling, "~>")
224 if (iend <= 1) then
225 call mom_mesg("No scaling indicator ~> found for "//trim(scaling))
226 else
227 ! Count and perhaps store this description and its weight.
228 ns = ns + 1
229 if (allocated(descs)) descs(ns) = scaling
230 if (allocated(wts)) then
231 wts(ns) = 1 ; if (present(weight)) wts(ns) = weight
232 endif
233 endif
234
235end subroutine add_scaling
236
237end module mom_check_scaling