MOM_obsolete_params.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!> Methods for testing for, and list of, obsolete run-time parameters.
7
8! This module was first conceived and written by Robert Hallberg, July 2010.
9
10use mom_error_handler, only : mom_error, fatal, warning, is_root_pe
11use mom_file_parser, only : read_param, log_version, param_file_type
12
13implicit none ; private
14
15#include <MOM_memory.h>
16
19
20contains
21
22!> Scans input parameter file for list obsolete parameters.
23subroutine find_obsolete_params(param_file)
24 type(param_file_type), intent(in) :: param_file !< Structure containing parameter file data.
25 ! Local variables
26 character(len=40) :: mdl = "find_obsolete_params" ! This module's name.
27! This include declares and sets the variable "version".
28#include "version_variable.h"
29 integer :: l_seg, nseg
30 logical :: test_logic, split
31 character(len=40) :: temp_string
32
33 if (.not.is_root_pe()) return
34
35 call obsolete_logical(param_file, "BLOCKED_ANALYTIC_FV_PGF", &
36 hint="BLOCKED_ANALYTIC_FV_PGF is no longer available.")
37
38 call obsolete_logical(param_file, "ADD_KV_SLOW", &
39 hint="This option is no longer needed, nor supported.")
40
41 call obsolete_char(param_file, "OBC_CONFIG", &
42 hint="Instead use OBC_USER_CONFIG and use the new segments protocol.")
43 call obsolete_char(param_file, "READ_OBC_ETA", &
44 hint="Instead use OBC_SEGMENT_XXX_DATA.")
45 call obsolete_char(param_file, "READ_OBC_UV", &
46 hint="Instead use OBC_SEGMENT_XXX_DATA.")
47 call obsolete_char(param_file, "READ_OBC_TS", &
48 hint="Instead use OBC_SEGMENT_XXX_DATA.")
49 call obsolete_char(param_file, "EXTEND_OBC_SEGMENTS", &
50 hint="This option is no longer needed, nor supported.")
51 call obsolete_char(param_file, "MEKE_VISCOSITY_COEFF", &
52 hint="This option has been replaced by MEKE_VISCOSITY_COEFF_KU and \n" //&
53 " MEKE_VISCOSITY_COEFF_AU. Please set these parameters instead.")
54 nseg = 0
55 call read_param(param_file, "OBC_NUMBER_OF_SEGMENTS", nseg)
56 do l_seg = 1,nseg
57 write(temp_string(1:22),"('OBC_SEGMENT_',i3.3,'_TNUDGE')") l_seg
58 call obsolete_real(param_file, temp_string, &
59 hint="Instead use OBC_SEGMENT_xxx_VELOCITY_NUDGING_TIMESCALES.")
60 enddo
61
62 call obsolete_logical(param_file, "CONVERT_THICKNESS_UNITS", .true.)
63 call obsolete_logical(param_file, "MASK_MASSLESS_TRACERS", .false.)
64
65 call obsolete_logical(param_file, "SALT_REJECT_BELOW_ML", .false.)
66 call obsolete_logical(param_file, "MLE_USE_MLD_AVE_BUG", .false.)
67 call obsolete_logical(param_file, "CORRECT_DENSITY", .true.)
68 call obsolete_char(param_file, "WINDSTRESS_STAGGER", warning_val="C", &
69 hint="Use WIND_STAGGER instead.")
70
71 call obsolete_char(param_file, "DIAG_REMAP_Z_GRID_DEF", &
72 hint="Use NUM_DIAG_COORDS, DIAG_COORDS and DIAG_COORD_DEF_Z")
73
74 call obsolete_real(param_file, "VSTAR_SCALE_FACTOR", hint="Use EPBL_VEL_SCALE_FACTOR instead.")
75
76 call obsolete_real(param_file, "VSTAR_SCALE_COEF")
77 call obsolete_real(param_file, "ZSTAR_RIGID_SURFACE_THRESHOLD")
78 call obsolete_logical(param_file, "HENYEY_IGW_BACKGROUND_NEW")
79
80 call obsolete_real(param_file, "SLIGHT_DZ_SURFACE")
81 call obsolete_int(param_file, "SLIGHT_NZ_SURFACE_FIXED")
82 call obsolete_real(param_file, "SLIGHT_SURFACE_AVG_DEPTH")
83 call obsolete_real(param_file, "SLIGHT_NLAY_TO_INTERIOR")
84 call obsolete_logical(param_file, "SLIGHT_FIX_HALOCLINES")
85 call obsolete_real(param_file, "HALOCLINE_FILTER_LENGTH")
86 call obsolete_real(param_file, "HALOCLINE_STRAT_TOL")
87
88 ! Test for inconsistent parameter settings.
89 split = .true. ; test_logic = .false.
90 call read_param(param_file,"SPLIT",split)
91 call read_param(param_file,"DYNAMIC_SURFACE_PRESSURE",test_logic)
92 if (test_logic .and. .not.split) call mom_error(fatal, &
93 "find_obsolete_params: #define DYNAMIC_SURFACE_PRESSURE is not yet "//&
94 "implemented without #define SPLIT.")
95 call obsolete_char(param_file, "CONTINUITY_SCHEME", warning_val="PPM", &
96 hint="Only one continuity scheme is available so this need not be specified.")
97 call obsolete_real(param_file, "ETA_TOLERANCE_AUX", only_warn=.true.)
98 call obsolete_real(param_file, "BT_MASS_SOURCE_LIMIT", 0.0)
99 call obsolete_real(param_file, "FIRST_GUESS_SURFACE_LAYER_DEPTH")
100 call obsolete_logical(param_file, "CORRECT_SURFACE_LAYER_AVERAGE")
101 call obsolete_int(param_file, "SEAMOUNT_LENGTH_SCALE", hint="Use SEAMOUNT_X_LENGTH_SCALE instead.")
102 call obsolete_int(param_file, "USE_LATERAL_BOUNDARY_DIFFUSION", &
103 hint="Use USE_HORIZONTAL_BOUNDARY_DIFFUSION instead.")
104
105 call obsolete_logical(param_file, "MSTAR_FIXED", hint="Instead use MSTAR_MODE.")
106 call obsolete_logical(param_file, "USE_VISBECK_SLOPE_BUG", .false.)
107 call obsolete_logical(param_file, "Use_PP81", hint="get_param is case sensitive so use USE_PP81.")
108
109 call obsolete_logical(param_file, "ALLOW_CLOCKS_IN_OMP_LOOPS", .true.)
110 call obsolete_logical(param_file, "LARGE_FILE_SUPPORT", .true.)
111 call obsolete_real(param_file, "MIN_Z_DIAG_INTERVAL")
112 call obsolete_char(param_file, "Z_OUTPUT_GRID_FILE")
113
114 call obsolete_logical(param_file, "CFL_BASED_TRUNCATIONS", .true.)
115 call obsolete_logical(param_file, "KD_BACKGROUND_VIA_KDML_BUG", .false.)
116 call obsolete_logical(param_file, "USE_DIABATIC_TIME_BUG", .false.)
117
118 call read_param(param_file, "INTERPOLATE_SPONGE_TIME_SPACE", test_logic)
119 call obsolete_logical(param_file, "NEW_SPONGES", warning_val=test_logic, &
120 hint="Use INTERPOLATE_SPONGE_TIME_SPACE instead.")
121
122 test_logic = .true. ; call read_param(param_file, "BOUND_KH", test_logic)
123 call obsolete_logical(param_file, "BETTER_BOUND_KH", warning_val=test_logic, hint="Use BOUND_KH alone.")
124 test_logic = .true. ; call read_param(param_file, "BOUND_AH", test_logic)
125 call obsolete_logical(param_file, "BETTER_BOUND_AH", warning_val=test_logic, hint="Use BOUND_AH alone.")
126
127 test_logic = .false. ; call read_param(param_file, "UNSPLIT_DT_VISC_BUG", test_logic)
128 call obsolete_logical(param_file, "FIX_UNSPLIT_DT_VISC_BUG", warning_val=(.not.test_logic), &
129 hint="Use UNSPLIT_DT_VISC_BUG instead, but with the reversed meaning.")
130
131 call obsolete_logical(param_file, "SMOOTH_RI", hint="Instead use N_SMOOTH_RI.")
132
133 call obsolete_logical(param_file, "INTERNAL_TIDE_CORNER_ADVECT", .false.)
134 call obsolete_logical(param_file, "TIDE_USE_SAL_SCALAR", hint="Use SAL_SCALAR_APPROX instead.")
135 call obsolete_logical(param_file, "TIDAL_SAL_SHT", hint="Use SAL_HARMONICS instead.")
136 call obsolete_int(param_file, "TIDAL_SAL_SHT_DEGREE", hint="Use SAL_HARMONICS_DEGREE instead.")
137 call obsolete_real(param_file, "RHO_E", hint="Use RHO_SOLID_EARTH instead.")
138 call obsolete_logical(param_file, "DEFAULT_2018_ANSWERS", hint="Instead use DEFAULT_ANSWER_DATE.")
139
140 call obsolete_logical(param_file, "SURFACE_FORCING_2018_ANSWERS", &
141 hint="Instead use SURFACE_FORCING_ANSWER_DATE.")
142 call obsolete_logical(param_file, "WIND_GYRES_2018_ANSWERS", &
143 hint="Instead use WIND_GYRES_ANSWER_DATE.")
144
145 call obsolete_logical(param_file, "BAROTROPIC_2018_ANSWERS", &
146 hint="Instead use BAROTROPIC_ANSWER_DATE.")
147 call obsolete_logical(param_file, "EPBL_2018_ANSWERS", hint="Instead use EPBL_ANSWER_DATE.")
148 call obsolete_logical(param_file, "HOR_REGRID_2018_ANSWERS", &
149 hint="Instead use HOR_REGRID_ANSWER_DATE.")
150 call obsolete_logical(param_file, "HOR_VISC_2018_ANSWERS", &
151 hint="Instead use HOR_VISC_ANSWER_DATE.")
152 call obsolete_logical(param_file, "IDL_HURR_2018_ANSWERS", &
153 hint="Instead use IDL_HURR_ANSWER_DATE.")
154 call obsolete_logical(param_file, "MEKE_GEOMETRIC_2018_ANSWERS", &
155 hint="Instead use MEKE_GEOMETRIC_ANSWER_DATE.")
156 call obsolete_logical(param_file, "ODA_2018_ANSWERS", hint="Instead use ODA_ANSWER_DATE.")
157 call obsolete_logical(param_file, "OPTICS_2018_ANSWERS", hint="Instead use OPTICS_ANSWER_DATE.")
158 call obsolete_logical(param_file, "REGULARIZE_LAYERS_2018_ANSWERS", &
159 hint="Instead use REGULARIZE_LAYERS_ANSWER_DATE.")
160 call obsolete_logical(param_file, "REMAPPING_2018_ANSWERS", &
161 hint="Instead use REMAPPING_ANSWER_DATE.")
162 call obsolete_logical(param_file, "SET_DIFF_2018_ANSWERS", &
163 hint="Instead use SET_DIFF_ANSWER_DATE.")
164 call obsolete_logical(param_file, "SET_VISC_2018_ANSWERS", &
165 hint="Instead use SET_VISC_ANSWER_DATE.")
166 call obsolete_logical(param_file, "SURFACE_2018_ANSWERS", hint="Instead use SURFACE_ANSWER_DATE.")
167 call obsolete_logical(param_file, "TIDAL_MIXING_2018_ANSWERS", &
168 hint="Instead use TIDAL_MIXING_ANSWER_DATE.")
169 call obsolete_logical(param_file, "VERT_FRICTION_2018_ANSWERS", &
170 hint="Instead use VERT_FRICTION_ANSWER_DATE.")
171
172 call obsolete_logical(param_file, "USE_GRID_SPACE_DIAGNOSTIC_AXES", &
173 hint="Instead use USE_INDEX_DIAGNOSTIC_AXIS.")
174
175 ! Write the file version number to the model log.
176 call log_version(param_file, mdl, version)
177
178end subroutine find_obsolete_params
179
180!> Test for presence of obsolete LOGICAL in parameter file.
181subroutine obsolete_logical(param_file, varname, warning_val, hint)
182 type(param_file_type), intent(in) :: param_file !< Structure containing parameter file data.
183 character(len=*), intent(in) :: varname !< Name of obsolete LOGICAL parameter.
184 logical, optional, intent(in) :: warning_val !< An allowed value that causes a warning instead of an error.
185 character(len=*), optional, intent(in) :: hint !< A hint to the user about what to do.
186 ! Local variables
187 logical :: test_logic, fatal_err
188 logical :: var_is_set ! True if this value was read by read_param.
189 character(len=128) :: hint_msg
190
191 test_logic = .false. ; call read_param(param_file, varname, test_logic, set=var_is_set)
192 fatal_err = .true.
193 if (var_is_set .and. present(warning_val)) fatal_err = (warning_val .neqv. test_logic)
194 hint_msg = " " ; if (present(hint)) hint_msg = hint
195
196 if (var_is_set) then
197 if (fatal_err) then
198 call mom_error(fatal, "MOM_obsolete_params: "//trim(varname)// &
199 " is an obsolete run-time flag, and should not be used. "// &
200 trim(hint_msg))
201 else
202 call mom_error(warning, "MOM_obsolete_params: "//trim(varname)// &
203 " is an obsolete run-time flag. "//trim(hint_msg))
204 endif
205 endif
206
207end subroutine obsolete_logical
208
209!> Test for presence of obsolete STRING in parameter file.
210subroutine obsolete_char(param_file, varname, warning_val, hint)
211 type(param_file_type), intent(in) :: param_file !< Structure containing parameter file data.
212 character(len=*), intent(in) :: varname !< Name of obsolete STRING parameter.
213 character(len=*), optional, intent(in) :: warning_val !< An allowed value that causes a warning instead of an error.
214 character(len=*), optional, intent(in) :: hint !< A hint to the user about what to do.
215 ! Local variables
216 character(len=200) :: test_string, hint_msg
217 logical :: var_is_set ! True if this value was read by read_param.
218 logical :: only_warn
219
220 test_string = '' ; call read_param(param_file, varname, test_string, set=var_is_set)
221 hint_msg = " " ; if (present(hint)) hint_msg = hint
222
223 if (var_is_set) then
224 only_warn = .false.
225 if (present(warning_val)) then ! Check if test_string and warning_val are the same.
226 if (len_trim(warning_val) == len_trim(test_string)) then
227 if (index(trim(test_string), trim(warning_val)) == 1) only_warn = .true.
228 endif
229 endif
230
231 if (only_warn) then
232 call mom_error(warning, &
233 "MOM_obsolete_params: "//trim(varname)// &
234 " is an obsolete run-time flag. "//trim(hint_msg))
235 else
236 call mom_error(fatal, &
237 "MOM_obsolete_params: "//trim(varname)// &
238 " is an obsolete run-time flag, and should not be used. "//trim(hint_msg))
239 endif
240 endif
241end subroutine obsolete_char
242
243!> Test for presence of obsolete REAL in parameter file.
244subroutine obsolete_real(param_file, varname, warning_val, hint, only_warn)
245 type(param_file_type), intent(in) :: param_file !< Structure containing parameter file data.
246 character(len=*), intent(in) :: varname !< Name of obsolete REAL parameter.
247 real, optional, intent(in) :: warning_val !< An allowed value that causes a warning instead of an error.
248 character(len=*), optional, intent(in) :: hint !< A hint to the user about what to do.
249 logical, optional, intent(in) :: only_warn !< If present and true, issue warnings instead of fatal errors.
250
251 ! Local variables
252 real :: test_val, warn_val
253 logical :: var_is_set ! True if this value was read by read_param.
254 logical :: issue_warning
255 character(len=128) :: hint_msg
256
257 test_val = -9e35 ; call read_param(param_file, varname, test_val, set=var_is_set)
258 warn_val = -9e35 ; if (present(warning_val)) warn_val = warning_val
259 hint_msg = " " ; if (present(hint)) hint_msg = hint
260 issue_warning = .false. ; if (present(only_warn)) issue_warning = only_warn
261
262 if (var_is_set) then
263 if ((test_val == warn_val) .or. issue_warning) then
264 call mom_error(warning, "MOM_obsolete_params: "//trim(varname)// &
265 " is an obsolete run-time flag. "//trim(hint_msg))
266 else
267 call mom_error(fatal, "MOM_obsolete_params: "//trim(varname)// &
268 " is an obsolete run-time flag, and should not be used. "// &
269 trim(hint_msg))
270 endif
271 endif
272end subroutine obsolete_real
273
274!> Test for presence of obsolete INTEGER in parameter file.
275subroutine obsolete_int(param_file, varname, warning_val, hint)
276 type(param_file_type), intent(in) :: param_file !< Structure containing parameter file data.
277 character(len=*), intent(in) :: varname !< Name of obsolete INTEGER parameter.
278 integer, optional, intent(in) :: warning_val !< An allowed value that causes a warning instead of an error.
279 character(len=*), optional, intent(in) :: hint !< A hint to the user about what to do.
280 ! Local variables
281 logical :: var_is_set ! True if this value was read by read_param.
282 integer :: test_val, warn_val
283 character(len=128) :: hint_msg
284
285 test_val = -123456788 ; call read_param(param_file, varname, test_val, set=var_is_set)
286 warn_val = -123456788 ; if (present(warning_val)) warn_val = warning_val
287 hint_msg = " " ; if (present(hint)) hint_msg = hint
288
289 if (var_is_set) then
290 if (test_val == warn_val) then
291 call mom_error(warning, "MOM_obsolete_params: "//trim(varname)// &
292 " is an obsolete run-time flag. "//trim(hint_msg))
293 else
294 call mom_error(fatal, "MOM_obsolete_params: "//trim(varname)// &
295 " is an obsolete run-time flag, and should not be used. "// &
296 trim(hint_msg))
297 endif
298 endif
299end subroutine obsolete_int
300
301end module mom_obsolete_params