This subroutine allocates space for the static variables used by this module. The metrics may be effectively 0, 1, or 2-D arrays, while fields like the background viscosities are 2-D arrays. ALLOC is a macro defined in MOM_memory.h for allocate or nothing with static memory.
67 type(time_type),
intent(in) :: Time
68 type(ocean_grid_type),
intent(inout) :: G
69 type(param_file_type),
intent(in) :: param_file
70 type(tidal_forcing_CS),
pointer :: CS
73 real,
dimension(SZI_(G), SZJ_(G)) :: &
77 real,
dimension(MAX_CONSTITUENTS) :: freq_def, phase0_def, amp_def, love_def
79 logical :: use_M2, use_S2, use_N2, use_K2, use_K1, use_O1, use_P1, use_Q1
80 logical :: use_MF, use_MM
82 logical :: FAIL_IF_MISSING = .true.
84 #include "version_variable.h"
85 character(len=40) :: mdl =
"MOM_tidal_forcing"
86 character(len=128) :: mesg
87 character(len=200) :: tidal_input_files(4*MAX_CONSTITUENTS)
88 integer :: i, j, c, is, ie, js, je, isd, ied, jsd, jed, nc
89 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
90 isd = g%isd ; ied = g%ied ; jsd = g%jsd; jed = g%jed
92 if (
associated(cs))
then
93 call mom_error(warning,
"tidal_forcing_init called with an associated "// &
99 call log_version(param_file, mdl, version,
"")
100 call get_param(param_file, mdl,
"TIDES", tides, &
101 "If true, apply tidal momentum forcing.", default=.false.)
103 if (.not.tides)
return
109 allocate(cs%sin_struct(isd:ied,jsd:jed,3)) ; cs%sin_struct(:,:,:) = 0.0
110 allocate(cs%cos_struct(isd:ied,jsd:jed,3)) ; cs%cos_struct(:,:,:) = 0.0
111 deg_to_rad = 4.0*atan(1.0)/180.0
112 do j=js-1,je+1 ;
do i=is-1,ie+1
113 lat_rad(i,j) = g%geoLatT(i,j)*deg_to_rad
114 lon_rad(i,j) = g%geoLonT(i,j)*deg_to_rad
116 do j=js-1,je+1 ;
do i=is-1,ie+1
117 cs%sin_struct(i,j,1) = -sin(2.0*lat_rad(i,j)) * sin(lon_rad(i,j))
118 cs%cos_struct(i,j,1) = sin(2.0*lat_rad(i,j)) * cos(lon_rad(i,j))
119 cs%sin_struct(i,j,2) = -cos(lat_rad(i,j))**2 * sin(2.0*lon_rad(i,j))
120 cs%cos_struct(i,j,2) = cos(lat_rad(i,j))**2 * cos(2.0*lon_rad(i,j))
121 cs%sin_struct(i,j,3) = 0.0
122 cs%cos_struct(i,j,3) = (0.5-1.5*sin(lat_rad(i,j))**2)
125 call get_param(param_file, mdl,
"TIDE_M2", use_m2, &
126 "If true, apply tidal momentum forcing at the M2 "//&
127 "frequency. This is only used if TIDES is true.", &
129 call get_param(param_file, mdl,
"TIDE_S2", use_s2, &
130 "If true, apply tidal momentum forcing at the S2 "//&
131 "frequency. This is only used if TIDES is true.", &
133 call get_param(param_file, mdl,
"TIDE_N2", use_n2, &
134 "If true, apply tidal momentum forcing at the N2 "//&
135 "frequency. This is only used if TIDES is true.", &
137 call get_param(param_file, mdl,
"TIDE_K2", use_k2, &
138 "If true, apply tidal momentum forcing at the K2 "//&
139 "frequency. This is only used if TIDES is true.", &
141 call get_param(param_file, mdl,
"TIDE_K1", use_k1, &
142 "If true, apply tidal momentum forcing at the K1 "//&
143 "frequency. This is only used if TIDES is true.", &
145 call get_param(param_file, mdl,
"TIDE_O1", use_o1, &
146 "If true, apply tidal momentum forcing at the O1 "//&
147 "frequency. This is only used if TIDES is true.", &
149 call get_param(param_file, mdl,
"TIDE_P1", use_p1, &
150 "If true, apply tidal momentum forcing at the P1 "//&
151 "frequency. This is only used if TIDES is true.", &
153 call get_param(param_file, mdl,
"TIDE_Q1", use_q1, &
154 "If true, apply tidal momentum forcing at the Q1 "//&
155 "frequency. This is only used if TIDES is true.", &
157 call get_param(param_file, mdl,
"TIDE_MF", use_mf, &
158 "If true, apply tidal momentum forcing at the MF "//&
159 "frequency. This is only used if TIDES is true.", &
161 call get_param(param_file, mdl,
"TIDE_MM", use_mm, &
162 "If true, apply tidal momentum forcing at the MM "//&
163 "frequency. This is only used if TIDES is true.", &
168 if (use_m2) nc=nc+1 ;
if (use_s2) nc=nc+1
169 if (use_n2) nc=nc+1 ;
if (use_k2) nc=nc+1
170 if (use_k1) nc=nc+1 ;
if (use_o1) nc=nc+1
171 if (use_p1) nc=nc+1 ;
if (use_q1) nc=nc+1
172 if (use_mf) nc=nc+1 ;
if (use_mm) nc=nc+1
176 call mom_error(fatal,
"tidal_forcing_init: "// &
177 "TIDES are defined, but no tidal constituents are used.")
181 call get_param(param_file, mdl,
"TIDAL_SAL_FROM_FILE", cs%tidal_sal_from_file, &
182 "If true, read the tidal self-attraction and loading "//&
183 "from input files, specified by TIDAL_INPUT_FILE. "//&
184 "This is only used if TIDES is true.", default=.false.)
185 call get_param(param_file, mdl,
"USE_PREVIOUS_TIDES", cs%use_prev_tides, &
186 "If true, use the SAL from the previous iteration of the "//&
187 "tides to facilitate convergent iteration. "//&
188 "This is only used if TIDES is true.", default=.false.)
189 call get_param(param_file, mdl,
"TIDE_USE_SAL_SCALAR", cs%use_sal_scalar, &
190 "If true and TIDES is true, use the scalar approximation "//&
191 "when calculating self-attraction and loading.", &
192 default=.not.cs%tidal_sal_from_file)
194 if (cs%use_sal_scalar .or. cs%use_prev_tides) &
195 call get_param(param_file, mdl,
"TIDE_SAL_SCALAR_VALUE", cs%sal_scalar, &
196 "The constant of proportionality between sea surface "//&
197 "height (really it should be bottom pressure) anomalies "//&
198 "and bottom geopotential anomalies. This is only used if "//&
199 "TIDES and TIDE_USE_SAL_SCALAR are true.", units=
"m m-1", &
200 fail_if_missing=.true.)
202 if (nc > max_constituents)
then
203 write(mesg,
'("Increase MAX_CONSTITUENTS in MOM_tidal_forcing.F90 to at least",I3, &
204 &"to accommodate all the registered tidal constituents.")') nc
205 call mom_error(fatal,
"MOM_tidal_forcing"//mesg)
208 do c=1,4*max_constituents ; tidal_input_files(c) =
"" ;
enddo
210 if (cs%tidal_sal_from_file .or. cs%use_prev_tides)
then
211 call get_param(param_file, mdl,
"TIDAL_INPUT_FILE", tidal_input_files, &
212 "A list of input files for tidal information.", &
213 default =
"", fail_if_missing=.true.)
219 c=c+1 ; cs%const_name(c) =
"M2" ; cs%freq(c) = 1.4051890e-4 ; cs%struct(c) = 2
220 cs%love_no(c) = 0.693 ; cs%amp(c) = 0.242334 ; cs%phase0(c) = 0.0
221 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
222 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
226 c=c+1 ; cs%const_name(c) =
"S2" ; cs%freq(c) = 1.4544410e-4 ; cs%struct(c) = 2
227 cs%love_no(c) = 0.693 ; cs%amp(c) = 0.112743 ; cs%phase0(c) = 0.0
228 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
229 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
233 c=c+1 ; cs%const_name(c) =
"N2" ; cs%freq(c) = 1.3787970e-4 ; cs%struct(c) = 2
234 cs%love_no(c) = 0.693 ; cs%amp(c) = 0.046397 ; cs%phase0(c) = 0.0
235 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
236 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
240 c=c+1 ; cs%const_name(c) =
"K2" ; cs%freq(c) = 1.4584234e-4 ; cs%struct(c) = 2
241 cs%love_no(c) = 0.693 ; cs%amp(c) = 0.030684 ; cs%phase0(c) = 0.0
242 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
243 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
247 c=c+1 ; cs%const_name(c) =
"K1" ; cs%freq(c) = 0.7292117e-4 ; cs%struct(c) = 1
248 cs%love_no(c) = 0.736 ; cs%amp(c) = 0.141565 ; cs%phase0(c) = 0.0
249 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
250 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
254 c=c+1 ; cs%const_name(c) =
"O1" ; cs%freq(c) = 0.6759774e-4 ; cs%struct(c) = 1
255 cs%love_no(c) = 0.695 ; cs%amp(c) = 0.100661 ; cs%phase0(c) = 0.0
256 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
257 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
261 c=c+1 ; cs%const_name(c) =
"P1" ; cs%freq(c) = 0.7252295e-4 ; cs%struct(c) = 1
262 cs%love_no(c) = 0.706 ; cs%amp(c) = 0.046848 ; cs%phase0(c) = 0.0
263 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
264 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
268 c=c+1 ; cs%const_name(c) =
"Q1" ; cs%freq(c) = 0.6495854e-4 ; cs%struct(c) = 1
269 cs%love_no(c) = 0.695 ; cs%amp(c) = 0.019273 ; cs%phase0(c) = 0.0
270 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
271 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
275 c=c+1 ; cs%const_name(c) =
"MF" ; cs%freq(c) = 0.053234e-4 ; cs%struct(c) = 3
276 cs%love_no(c) = 0.693 ; cs%amp(c) = 0.042041 ; cs%phase0(c) = 0.0
277 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
278 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
282 c=c+1 ; cs%const_name(c) =
"MM" ; cs%freq(c) = 0.026392e-4 ; cs%struct(c) = 3
283 cs%love_no(c) = 0.693 ; cs%amp(c) = 0.022191 ; cs%phase0(c) = 0.0
284 freq_def(c) = cs%freq(c) ; love_def(c) = cs%love_no(c)
285 amp_def(c) = cs%amp(c) ; phase0_def(c) = cs%phase0(c)
292 call get_param(param_file, mdl,
"TIDE_"//trim(cs%const_name(c))//
"_FREQ", cs%freq(c), &
293 "Frequency of the "//trim(cs%const_name(c))//
" tidal constituent. "//&
294 "This is only used if TIDES and TIDE_"//trim(cs%const_name(c))// &
295 " are true.", units=
"s-1", default=freq_def(c))
296 call get_param(param_file, mdl,
"TIDE_"//trim(cs%const_name(c))//
"_AMP", cs%amp(c), &
297 "Amplitude of the "//trim(cs%const_name(c))//
" tidal constituent. "//&
298 "This is only used if TIDES and TIDE_"//trim(cs%const_name(c))// &
299 " are true.", units=
"m", default=amp_def(c))
300 call get_param(param_file, mdl,
"TIDE_"//trim(cs%const_name(c))//
"_PHASE_T0", cs%phase0(c), &
301 "Phase of the "//trim(cs%const_name(c))//
" tidal constituent at time 0. "//&
302 "This is only used if TIDES and TIDE_"//trim(cs%const_name(c))// &
303 " are true.", units=
"radians", default=phase0_def(c))
306 if (cs%tidal_sal_from_file)
then
307 allocate(cs%cosphasesal(isd:ied,jsd:jed,nc))
308 allocate(cs%sinphasesal(isd:ied,jsd:jed,nc))
309 allocate(cs%ampsal(isd:ied,jsd:jed,nc))
312 call find_in_files(tidal_input_files,
"PHASE_SAL_"//trim(cs%const_name(c)),phase,g)
313 call find_in_files(tidal_input_files,
"AMP_SAL_"//trim(cs%const_name(c)),cs%ampsal(:,:,c),g)
314 call pass_var(phase, g%domain,complete=.false.)
315 call pass_var(cs%ampsal(:,:,c),g%domain,complete=.true.)
316 do j=js-1,je+1 ;
do i=is-1,ie+1
317 cs%cosphasesal(i,j,c) = cos(phase(i,j)*deg_to_rad)
318 cs%sinphasesal(i,j,c) = sin(phase(i,j)*deg_to_rad)
323 if (cs%USE_PREV_TIDES)
then
324 allocate(cs%cosphase_prev(isd:ied,jsd:jed,nc))
325 allocate(cs%sinphase_prev(isd:ied,jsd:jed,nc))
326 allocate(cs%amp_prev(isd:ied,jsd:jed,nc))
329 call find_in_files(tidal_input_files,
"PHASE_PREV_"//trim(cs%const_name(c)),phase,g)
330 call find_in_files(tidal_input_files,
"AMP_PREV_"//trim(cs%const_name(c)),cs%amp_prev(:,:,c),g)
331 call pass_var(phase, g%domain,complete=.false.)
332 call pass_var(cs%amp_prev(:,:,c),g%domain,complete=.true.)
333 do j=js-1,je+1 ;
do i=is-1,ie+1
334 cs%cosphase_prev(i,j,c) = cos(phase(i,j)*deg_to_rad)
335 cs%sinphase_prev(i,j,c) = sin(phase(i,j)*deg_to_rad)
340 id_clock_tides = cpu_clock_id(
'(Ocean tides)', grain=clock_module)