MOM_unit_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!> Provides a transparent unit rescaling type to facilitate dimensional consistency testing
7
8use mom_error_handler, only : mom_error, mom_mesg, fatal
9use mom_file_parser, only : get_param, log_param, log_version, param_file_type
10
11implicit none ; private
12
13! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
14! consistency testing. These are noted in comments with units like Z, H, L, T, R and Q, along with
15! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the rescaled
16! combination is a nondimensional variable, the notation would be "a slope [Z L-1 ~> nondim]",
17! but if (as the case for the variables here), the rescaled combination is exactly 1, the right
18! notation would be something like "a dimensional scaling factor [Z m-1 ~> 1]". If the units
19! vary with the Boussinesq approximation, the Boussinesq variant is given first.
20
22
23!> Describes various unit conversion factors
24type, public :: unit_scale_type
25 real :: m_to_z !< A constant that translates distances in meters to the units of depth [Z m-1 ~> 1]
26 real :: z_to_m !< A constant that translates distances in the units of depth to meters [m Z-1 ~> 1]
27 real :: m_to_l !< A constant that translates lengths in meters to the units of horizontal lengths [L m-1 ~> 1]
28 real :: l_to_m !< A constant that translates lengths in the units of horizontal lengths to meters [m L-1 ~> 1]
29 real :: s_to_t !< A constant that translates time intervals in seconds to the units of time [T s-1 ~> 1]
30 real :: t_to_s !< A constant that translates the units of time to seconds [s T-1 ~> 1]
31 real :: r_to_kg_m3 !< A constant that translates the units of density to kilograms per meter cubed [kg m-3 R-1 ~> 1]
32 real :: kg_m3_to_r !< A constant that translates kilograms per meter cubed to the units of density [R m3 kg-1 ~> 1]
33 real :: q_to_j_kg !< A constant that translates the units of enthalpy to Joules per kilogram [J kg-1 Q-1 ~> 1]
34 real :: j_kg_to_q !< A constant that translates Joules per kilogram to the units of enthalpy [Q kg J-1 ~> 1]
35 real :: c_to_degc !< A constant that translates the units of temperature to degrees Celsius [degC C-1 ~> 1]
36 real :: degc_to_c !< A constant that translates degrees Celsius to the units of temperature [C degC-1 ~> 1]
37 real :: s_to_ppt !< A constant that translates the units of salinity to parts per thousand [ppt S-1 ~> 1]
38 real :: ppt_to_s !< A constant that translates parts per thousand to the units of salinity [S ppt-1 ~> 1]
39
40 ! These are useful combinations of the fundamental scale conversion factors above.
41 real :: z_to_l !< Convert vertical distances to lateral lengths [L Z-1 ~> 1]
42 real :: l_to_z !< Convert lateral lengths to vertical distances [Z L-1 ~> 1]
43 real :: l_t_to_m_s !< Convert lateral velocities from L T-1 to m s-1 [T m L-1 s-1 ~> 1]
44 real :: m_s_to_l_t !< Convert lateral velocities from m s-1 to L T-1 [L s T-1 m-1 ~> 1]
45 real :: l_t2_to_m_s2 !< Convert lateral accelerations from L T-2 to m s-2 [L s2 T-2 m-1 ~> 1]
46 real :: z2_t_to_m2_s !< Convert vertical diffusivities from Z2 T-1 to m2 s-1 [T m2 Z-2 s-1 ~> 1]
47 real :: m2_s_to_z2_t !< Convert vertical diffusivities from m2 s-1 to Z2 T-1 [Z2 s T-1 m-2 ~> 1]
48 real :: w_m2_to_qrz_t !< Convert heat fluxes from W m-2 to Q R Z T-1 [Q R Z m2 T-1 W-1 ~> 1]
49 real :: qrz_t_to_w_m2 !< Convert heat fluxes from Q R Z T-1 to W m-2 [W T Q-1 R-1 Z-1 m-2 ~> 1]
50 ! Not used enough: real :: kg_m2_to_RZ !< Convert mass loads from kg m-2 to R Z [R Z m2 kg-1 ~> 1]
51 real :: rz_to_kg_m2 !< Convert mass loads from R Z to kg m-2 [kg R-1 Z-1 m-2 ~> 1]
52 real :: rzl2_to_kg !< Convert masses from R Z L2 to kg [kg R-1 Z-1 L-2 ~> 1]
53 real :: kg_m2s_to_rz_t !< Convert mass fluxes from kg m-2 s-1 to R Z T-1 [R Z m2 s T-1 kg-1 ~> 1]
54 real :: rz_t_to_kg_m2s !< Convert mass fluxes from R Z T-1 to kg m-2 s-1 [T kg R-1 Z-1 m-2 s-1 ~> 1]
55 real :: rz3_t3_to_w_m2 !< Convert turbulent kinetic energy fluxes from R Z3 T-3 to W m-2 [W T3 R-1 Z-3 m-2 ~> 1]
56 real :: w_m2_to_rz3_t3 !< Convert turbulent kinetic energy fluxes from W m-2 to R Z3 T-3 [R Z3 m2 T-3 W-1 ~> 1]
57 real :: rl2_t2_to_pa !< Convert pressures from R L2 T-2 to Pa [Pa T2 R-1 L-2 ~> 1]
58 real :: rlz_t2_to_pa !< Convert wind stresses from R L Z T-2 to Pa [Pa T2 R-1 L-1 Z-1 ~> 1]
59 real :: pa_to_rl2_t2 !< Convert pressures from Pa to R L2 T-2 [R L2 T-2 Pa-1 ~> 1]
60 real :: pa_to_rlz_t2 !< Convert wind stresses from Pa to R L Z T-2 [R L Z T-2 Pa-1 ~> 1]
61
62 ! These are no longer used for changing scaling across restarts.
63 real :: m_to_z_restart = 1.0 !< A copy of the m_to_Z that is used in restart files.
64 real :: m_to_l_restart = 1.0 !< A copy of the m_to_L that is used in restart files.
65 real :: s_to_t_restart = 1.0 !< A copy of the s_to_T that is used in restart files.
66 real :: kg_m3_to_r_restart = 1.0 !< A copy of the kg_m3_to_R that is used in restart files.
67 real :: j_kg_to_q_restart = 1.0 !< A copy of the J_kg_to_Q that is used in restart files.
68end type unit_scale_type
69
70contains
71
72!> Allocates and initializes the ocean model unit scaling type
73subroutine unit_scaling_init( param_file, US )
74 type(param_file_type), intent(in) :: param_file !< Parameter file handle/type
75 type(unit_scale_type), pointer :: us !< A dimensional unit scaling type
76
77 ! This routine initializes a unit_scale_type structure (US).
78
79 ! Local variables
80 integer :: z_power, l_power, t_power, r_power, q_power, c_power, s_power
81 real :: z_rescale_factor, l_rescale_factor, t_rescale_factor, r_rescale_factor, q_rescale_factor
82 real :: c_rescale_factor, s_rescale_factor
83 ! This include declares and sets the variable "version".
84# include "version_variable.h"
85 character(len=16) :: mdl = "MOM_unit_scaling"
86
87 if (associated(us)) call mom_error(fatal, &
88 'unit_scaling_init: called with an associated US pointer.')
89 allocate(us)
90
91 ! Read all relevant parameters and write them to the model log.
92 call log_version(param_file, mdl, version, &
93 "Parameters for doing unit scaling of variables.", debugging=.true.)
94 call get_param(param_file, mdl, "Z_RESCALE_POWER", z_power, &
95 "An integer power of 2 that is used to rescale the model's "//&
96 "internal units of depths and heights. Valid values range from -300 to 300.", &
97 default=0, debuggingparam=.true.)
98 call get_param(param_file, mdl, "L_RESCALE_POWER", l_power, &
99 "An integer power of 2 that is used to rescale the model's "//&
100 "internal units of lateral distances. Valid values range from -300 to 300.", &
101 default=0, debuggingparam=.true.)
102 call get_param(param_file, mdl, "T_RESCALE_POWER", t_power, &
103 "An integer power of 2 that is used to rescale the model's "//&
104 "internal units of time. Valid values range from -300 to 300.", &
105 default=0, debuggingparam=.true.)
106 call get_param(param_file, mdl, "R_RESCALE_POWER", r_power, &
107 "An integer power of 2 that is used to rescale the model's "//&
108 "internal units of density. Valid values range from -300 to 300.", &
109 default=0, debuggingparam=.true.)
110 call get_param(param_file, mdl, "Q_RESCALE_POWER", q_power, &
111 "An integer power of 2 that is used to rescale the model's "//&
112 "internal units of heat content. Valid values range from -300 to 300.", &
113 default=0, debuggingparam=.true.)
114 call get_param(param_file, mdl, "C_RESCALE_POWER", c_power, &
115 "An integer power of 2 that is used to rescale the model's "//&
116 "internal units of temperature. Valid values range from -300 to 300.", &
117 default=0, debuggingparam=.true.)
118 call get_param(param_file, mdl, "S_RESCALE_POWER", s_power, &
119 "An integer power of 2 that is used to rescale the model's "//&
120 "internal units of salinity. Valid values range from -300 to 300.", &
121 default=0, debuggingparam=.true.)
122
123 if (abs(z_power) > 300) call mom_error(fatal, "unit_scaling_init: "//&
124 "Z_RESCALE_POWER is outside of the valid range of -300 to 300.")
125 if (abs(l_power) > 300) call mom_error(fatal, "unit_scaling_init: "//&
126 "L_RESCALE_POWER is outside of the valid range of -300 to 300.")
127 if (abs(t_power) > 300) call mom_error(fatal, "unit_scaling_init: "//&
128 "T_RESCALE_POWER is outside of the valid range of -300 to 300.")
129 if (abs(r_power) > 300) call mom_error(fatal, "unit_scaling_init: "//&
130 "R_RESCALE_POWER is outside of the valid range of -300 to 300.")
131 if (abs(q_power) > 300) call mom_error(fatal, "unit_scaling_init: "//&
132 "Q_RESCALE_POWER is outside of the valid range of -300 to 300.")
133 if (abs(c_power) > 300) call mom_error(fatal, "unit_scaling_init: "//&
134 "C_RESCALE_POWER is outside of the valid range of -300 to 300.")
135 if (abs(s_power) > 300) call mom_error(fatal, "unit_scaling_init: "//&
136 "S_RESCALE_POWER is outside of the valid range of -300 to 300.")
137
138 z_rescale_factor = 1.0
139 if (z_power /= 0) z_rescale_factor = 2.0**z_power
140 us%Z_to_m = 1.0 * z_rescale_factor
141 us%m_to_Z = 1.0 / z_rescale_factor
142
143 l_rescale_factor = 1.0
144 if (l_power /= 0) l_rescale_factor = 2.0**l_power
145 us%L_to_m = 1.0 * l_rescale_factor
146 us%m_to_L = 1.0 / l_rescale_factor
147
148 t_rescale_factor = 1.0
149 if (t_power /= 0) t_rescale_factor = 2.0**t_power
150 us%T_to_s = 1.0 * t_rescale_factor
151 us%s_to_T = 1.0 / t_rescale_factor
152
153 r_rescale_factor = 1.0
154 if (r_power /= 0) r_rescale_factor = 2.0**r_power
155 us%R_to_kg_m3 = 1.0 * r_rescale_factor
156 us%kg_m3_to_R = 1.0 / r_rescale_factor
157
158 q_rescale_factor = 1.0
159 if (q_power /= 0) q_rescale_factor = 2.0**q_power
160 us%Q_to_J_kg = 1.0 * q_rescale_factor
161 us%J_kg_to_Q = 1.0 / q_rescale_factor
162
163 c_rescale_factor = 1.0
164 if (c_power /= 0) c_rescale_factor = 2.0**c_power
165 us%C_to_degC = 1.0 * c_rescale_factor
166 us%degC_to_C = 1.0 / c_rescale_factor
167
168 s_rescale_factor = 1.0
169 if (s_power /= 0) s_rescale_factor = 2.0**s_power
170 us%S_to_ppt = 1.0 * s_rescale_factor
171 us%ppt_to_S = 1.0 / s_rescale_factor
172
174end subroutine unit_scaling_init
175
176!> Allocates and initializes the ocean model unit scaling type to unscaled values.
177subroutine unit_no_scaling_init(US)
178 type(unit_scale_type), pointer :: us !< A dimensional unit scaling type
179
180 if (associated(us)) call mom_error(fatal, &
181 'unit_scaling_init: called with an associated US pointer.')
182 allocate(us)
183
184 us%Z_to_m = 1.0 ; us%m_to_Z = 1.0
185 us%L_to_m = 1.0 ; us%m_to_L = 1.0
186 us%T_to_s = 1.0 ; us%s_to_T = 1.0
187 us%R_to_kg_m3 = 1.0 ; us%kg_m3_to_R = 1.0
188 us%Q_to_J_kg = 1.0 ; us%J_kg_to_Q = 1.0
189 us%C_to_degC = 1.0 ; us%degC_to_C = 1.0
190 us%S_to_ppt = 1.0 ; us%ppt_to_S = 1.0
191
193end subroutine unit_no_scaling_init
194
195!> This subroutine sets useful combinations of the fundamental scale conversion factors
196!! in the unit scaling type.
197subroutine set_unit_scaling_combos(US)
198 type(unit_scale_type), intent(inout) :: US !< A dimensional unit scaling type
199
200 ! Convert vertical to horizontal length scales and the reverse:
201 us%Z_to_L = us%Z_to_m * us%m_to_L
202 us%L_to_Z = us%L_to_m * us%m_to_Z
203 ! Horizontal velocities:
204 us%L_T_to_m_s = us%L_to_m * us%s_to_T
205 us%m_s_to_L_T = us%m_to_L * us%T_to_s
206 ! Horizontal accelerations:
207 us%L_T2_to_m_s2 = us%L_to_m * us%s_to_T**2
208 ! It does not look like US%m_s2_to_L_T2 would be used, so it does not exist.
209 ! Vertical diffusivities and viscosities:
210 us%Z2_T_to_m2_s = us%Z_to_m**2 * us%s_to_T
211 us%m2_s_to_Z2_T = us%m_to_Z**2 * us%T_to_s
212 ! Column mass loads:
213 us%RZ_to_kg_m2 = us%R_to_kg_m3 * us%Z_to_m
214 ! It does not seem like US%kg_m2_to_RZ would be used enough in MOM6 to justify its existence.
215 ! Vertical mass fluxes:
216 us%kg_m2s_to_RZ_T = us%kg_m3_to_R * us%m_to_Z * us%T_to_s
217 us%RZ_T_to_kg_m2s = us%R_to_kg_m3 * us%Z_to_m * us%s_to_T
218 ! Turbulent kinetic energy vertical fluxes:
219 us%RZ3_T3_to_W_m2 = us%R_to_kg_m3 * us%Z_to_m**3 * us%s_to_T**3
220 us%W_m2_to_RZ3_T3 = us%kg_m3_to_R * us%m_to_Z**3 * us%T_to_s**3
221 ! Vertical heat fluxes:
222 us%W_m2_to_QRZ_T = us%J_kg_to_Q * us%kg_m3_to_R * us%m_to_Z * us%T_to_s
223 us%QRZ_T_to_W_m2 = us%Q_to_J_kg * us%R_to_kg_m3 * us%Z_to_m * us%s_to_T
224 ! Pressures:
225 us%RL2_T2_to_Pa = us%R_to_kg_m3 * us%L_T_to_m_s**2
226 us%Pa_to_RL2_T2 = us%kg_m3_to_R * us%m_s_to_L_T**2
227 ! Wind stresses:
228 us%RLZ_T2_to_Pa = us%R_to_kg_m3 * us%L_T_to_m_s**2 * us%Z_to_L
229 us%Pa_to_RLZ_T2 = us%kg_m3_to_R * us%m_s_to_L_T**2 * us%L_to_Z
230 ! Masses:
231 us%RZL2_to_kg = us%R_to_kg_m3 * us%Z_to_m * us%L_to_m**2
232
233end subroutine set_unit_scaling_combos
234
235!> Set the unit scaling factors for output to restart files to the unit scaling
236!! factors for this run.
237subroutine fix_restart_unit_scaling(US, unscaled)
238 type(unit_scale_type), intent(inout) :: us !< A dimensional unit scaling type
239 logical, optional, intent(in) :: unscaled !< If true, set the restart factors as though the
240 !! model would be unscaled, which is appropriate if the
241 !! scaling is undone when writing a restart file.
242
243 us%m_to_Z_restart = 1.0 ! US%m_to_Z
244 us%m_to_L_restart = 1.0 ! US%m_to_L
245 us%s_to_T_restart = 1.0 ! US%s_to_T
246 us%kg_m3_to_R_restart = 1.0 ! US%kg_m3_to_R
247 us%J_kg_to_Q_restart = 1.0 ! US%J_kg_to_Q
248
249 if (present(unscaled)) then ; if (unscaled) then
250 us%m_to_Z_restart = 1.0
251 us%m_to_L_restart = 1.0
252 us%s_to_T_restart = 1.0
253 us%kg_m3_to_R_restart = 1.0
254 us%J_kg_to_Q_restart = 1.0
255 endif ; endif
256
257end subroutine fix_restart_unit_scaling
258
259!> Deallocates a unit scaling structure.
260subroutine unit_scaling_end( US )
261 type(unit_scale_type), pointer :: us !< A dimensional unit scaling type
262
263 deallocate( us )
264
265end subroutine unit_scaling_end
266
267end module mom_unit_scaling