MOM_PointAccel.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!> Debug accelerations at a given point
6!!
7!! The two subroutines in this file write out all of the terms
8!! in the u- or v-momentum balance at a given point. Usually
9!! these subroutines are called after the velocities exceed some
10!! threshold, in order to determine which term is culpable.
11!! often this is done for debugging purposes.
12module mom_pointaccel
13
15use mom_domains, only : pe_here
16use mom_error_handler, only : mom_error, note
17use mom_file_parser, only : get_param, log_param, log_version, param_file_type
18use mom_get_input, only : directories
19use mom_grid, only : ocean_grid_type
20use mom_io, only : open_ascii_file, append_file, multiple, single_file
21use mom_time_manager, only : time_type, get_time, get_date, set_date, operator(-)
25
26implicit none ; private
27
28#include <MOM_memory.h>
29
31
32!> The control structure for the MOM_PointAccel module
33type, public :: pointaccel_cs ; private
34 character(len=200) :: u_trunc_file !< The complete path to the file in which a column's worth of
35 !! u-accelerations are written if u-velocity truncations occur.
36 character(len=200) :: v_trunc_file !< The complete path to the file in which a column's worth of
37 !! v-accelerations are written if v-velocity truncations occur.
38 integer :: u_file !< The unit number for an opened u-truncation files, or -1 if it has not yet been opened.
39 integer :: v_file !< The unit number for an opened v-truncation files, or -1 if it has not yet been opened.
40 integer :: cols_written !< The number of columns whose output has been
41 !! written by this PE during the current run.
42 integer :: max_writes !< The maximum number of times any PE can write out
43 !! a column's worth of accelerations during a run.
44 logical :: full_column !< If true, write out the accelerations in all massive layers,
45 !! otherwise just document the ones with large velocities.
46 type(time_type), pointer :: time => null() !< A pointer to the ocean model's clock.
47 type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
48 !! regulate the timing of diagnostic output.
49! The following are pointers to many of the state variables and accelerations
50! that are used to step the physical model forward. They all use the same
51! names as the variables they point to in MOM.F90
52 real, pointer, dimension(:,:,:) :: &
53 u_av => null(), & !< Time average u-velocity [L T-1 ~> m s-1]
54 v_av => null(), & !< Time average velocity [L T-1 ~> m s-1]
55 u_prev => null(), & !< Previous u-velocity [L T-1 ~> m s-1]
56 v_prev => null(), & !< Previous v-velocity [L T-1 ~> m s-1]
57 t => null(), & !< Temperature [C ~> degC]
58 s => null(), & !< Salinity [S ~> ppt]
59 u_accel_bt => null(), & !< Barotropic u-accelerations [L T-2 ~> m s-2]
60 v_accel_bt => null() !< Barotropic v-accelerations [L T-2 ~> m s-2]
61end type pointaccel_cs
62
63contains
64
65!> This subroutine writes to an output file all of the accelerations
66!! that have been applied to a column of zonal velocities over the
67!! previous timestep. This subroutine is called from vertvisc.
68subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, str, a, hv)
69 integer, intent(in) :: i !< The zonal index of the column to be documented.
70 integer, intent(in) :: j !< The meridional index of the column to be documented.
71 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
72 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
73 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
74 real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), &
75 intent(in) :: um !< The new zonal velocity [L T-1 ~> m s-1].
76 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
77 intent(in) :: hin !< The layer thickness [H ~> m or kg m-2].
78 type(accel_diag_ptrs), intent(in) :: adp !< A structure pointing to the various
79 !! accelerations in the momentum equations.
80 type(cont_diag_ptrs), intent(in) :: cdp !< A structure with pointers to various terms
81 !! in the continuity equations.
82 real, intent(in) :: dt !< The ocean dynamics time step [T ~> s].
83 type(pointaccel_cs), pointer :: cs !< The control structure returned by a previous
84 !! call to PointAccel_init.
85 real, intent(in) :: vel_rpt !< The velocity magnitude that triggers a report [L T-1 ~> m s-1]
86 real, optional, intent(in) :: str !< The surface wind stress [R L Z T-2 ~> Pa]
87 real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)+1), &
88 optional, intent(in) :: a !< The layer coupling coefficients from vertvisc
89 !! [H T-1 ~> m s-1 or Pa s m-1]
90 real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), &
91 optional, intent(in) :: hv !< The layer thicknesses at velocity grid points,
92 !! from vertvisc [H ~> m or kg m-2].
93
94 ! Local variables
95 real :: cfl ! The local velocity-based CFL number [nondim]
96 real :: angstrom ! A negligibly small thickness [H ~> m or kg m-2]
97 real :: du ! A velocity change [L T-1 ~> m s-1]
98 real :: inorm(szk_(gv)) ! The inverse of the normalized velocity change [L T-1 ~> m s-1]
99 real :: e(szk_(gv)+1) ! Simple estimates of interface heights based on the sum of thicknesses [m]
100 real :: h_scale ! A scaling factor for thicknesses [m H-1 ~> 1] or [kg m-2 H-1 ~> 1]
101 real :: vel_scale ! A scaling factor for velocities [m T s-1 L-1 ~> 1]
102 real :: uh_scale ! A scaling factor for transport per unit length [m2 T s-1 L-1 H-1 ~> 1]
103 ! or [kg T m-1 s-1 L-1 H-1 ~> 1]
104 real :: temp_scale ! A scaling factor for temperatures [degC C-1 ~> 1]
105 real :: saln_scale ! A scaling factor for salinities [ppt S-1 ~> 1]
106 integer :: yr, mo, day, hr, minute, sec, yearday
107 integer :: k, ks, ke
108 integer :: nz
109 logical :: do_k(szk_(gv)+1)
110 logical :: prev_avail
111 integer :: file
112
113 angstrom = gv%Angstrom_H + gv%H_subroundoff
114 h_scale = gv%H_to_mks ; vel_scale = us%L_T_to_m_s ; uh_scale = h_scale*vel_scale
115 temp_scale = us%C_to_degC ; saln_scale = us%S_to_ppt
116
117! if (.not.associated(CS)) return
118 nz = gv%ke
119 if (cs%cols_written < cs%max_writes) then
120 cs%cols_written = cs%cols_written + 1
121
122 ks = 1 ; ke = nz
123 do_k(:) = .false.
124
125 ! Open up the file for output if this is the first call.
126 if (cs%u_file == -1) then
127 if (len_trim(cs%u_trunc_file) < 1) return
128 call open_ascii_file(cs%u_file, trim(cs%u_trunc_file), action=append_file, &
129 threading=multiple, fileset=single_file)
130 if (cs%u_file == -1) then
131 call mom_error(note, 'Unable to open file '//trim(cs%u_trunc_file)//'.')
132 return
133 endif
134 endif
135 file = cs%u_file
136
137 prev_avail = (associated(cs%u_prev) .and. associated(cs%v_prev))
138
139 ! Determine which layers to write out accelerations for.
140 do k=1,nz
141 if (((max(cs%u_av(i,j,k),um(i,j,k)) >= vel_rpt) .or. &
142 (min(cs%u_av(i,j,k),um(i,j,k)) <= -vel_rpt)) .and. &
143 ((hin(i,j,k) + hin(i+1,j,k)) > 3.0*angstrom)) exit
144 enddo
145 ks = k
146 do k=nz,1,-1
147 if (((max(cs%u_av(i,j,k), um(i,j,k)) >= vel_rpt) .or. &
148 (min(cs%u_av(i,j,k), um(i,j,k)) <= -vel_rpt)) .and. &
149 ((hin(i,j,k) + hin(i+1,j,k)) > 3.0*angstrom)) exit
150 enddo
151 ke = k
152 if (ke < ks) then
153 ks = 1 ; ke = nz ; write(file,'("U: Unable to set ks & ke.")')
154 endif
155 if (cs%full_column) then
156 ks = 1 ; ke = nz
157 endif
158
159 call get_date(cs%Time, yr, mo, day, hr, minute, sec)
160 call get_time((cs%Time - set_date(yr, 1, 1, 0, 0, 0)), sec, yearday)
161 write (file,'(/,"--------------------------")')
162 write (file,'(/,"Time ",I0," ",I0," ",F6.2," U-velocity violation at ",I0,": ",I0,", ",I0, &
163 & " (",F7.2," E ",F7.2," N) Layers ",I0," to ",I0,". dt = ",1PG10.4)') &
164 yr, yearday, (real(sec)/3600.0), pe_here(), i, j, &
165 g%geoLonCu(i,j), g%geoLatCu(i,j), ks, ke, us%T_to_s*dt
166
167 if (ks <= gv%nk_rho_varies) ks = 1
168 do k=ks,ke
169 if ((hin(i,j,k) + hin(i+1,j,k)) > 3.0*angstrom) do_k(k) = .true.
170 enddo
171
172 write(file,'(/,"Layers:")', advance='no')
173 do k=ks,ke ; if (do_k(k)) write(file,'(I10," ")', advance='no') (k) ; enddo
174 write(file,'(/,"u(m): ")', advance='no')
175 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*um(i,j,k)) ; enddo
176 if (prev_avail) then
177 write(file,'(/,"u(mp): ")', advance='no')
178 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%u_prev(i,j,k)) ; enddo
179 endif
180 write(file,'(/,"u(3): ")', advance='no')
181 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%u_av(i,j,k)) ; enddo
182
183 write(file,'(/,"CFL u: ")', advance='no')
184 do k=ks,ke ; if (do_k(k)) then
185 cfl = abs(um(i,j,k)) * dt * g%dy_Cu(i,j)
186 if (um(i,j,k) < 0.0) then ; cfl = cfl * g%IareaT(i+1,j)
187 else ; cfl = cfl * g%IareaT(i,j) ; endif
188 write(file,'(ES10.3," ")', advance='no') cfl
189 endif ; enddo
190 write(file,'(/,"CFL0 u:")', advance='no')
191 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
192 abs(um(i,j,k)) * dt * g%IdxCu(i,j) ; enddo
193
194 if (prev_avail) then
195 write(file,'(/,"du: ")', advance='no')
196 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
197 (vel_scale*(um(i,j,k)-cs%u_prev(i,j,k))) ; enddo
198 endif
199 write(file,'(/,"CAu: ")', advance='no')
200 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*dt*adp%CAu(i,j,k)) ; enddo
201 write(file,'(/,"PFu: ")', advance='no')
202 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*dt*adp%PFu(i,j,k)) ; enddo
203 write(file,'(/,"diffu: ")', advance='no')
204 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*dt*adp%diffu(i,j,k)) ; enddo
205
206 if (associated(adp%gradKEu)) then
207 write(file,'(/,"KEu: ")', advance='no')
208 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
209 (vel_scale*dt*adp%gradKEu(i,j,k)) ; enddo
210 endif
211 if (associated(adp%rv_x_v)) then
212 write(file,'(/,"Coru: ")', advance='no')
213 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
214 vel_scale*dt*(adp%CAu(i,j,k)-adp%rv_x_v(i,j,k)) ; enddo
215 endif
216 if (associated(adp%du_dt_visc)) then
217 write(file,'(/,"ubv: ")', advance='no')
218 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
219 vel_scale*(um(i,j,k) - dt*adp%du_dt_visc(i,j,k)) ; enddo
220 write(file,'(/,"duv: ")', advance='no')
221 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
222 (vel_scale*dt*adp%du_dt_visc(i,j,k)) ; enddo
223 endif
224 if (associated(adp%du_other)) then
225 write(file,'(/,"du_other: ")', advance='no')
226 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
227 (vel_scale*adp%du_other(i,j,k)) ; enddo
228 endif
229 if (present(a)) then
230 write(file,'(/,"a: ",ES10.3," ")', advance='no') h_scale*a(i,j,ks)*dt
231 do k=ks+1,ke+1 ; if (do_k(k-1)) write(file,'(ES10.3," ")', advance='no') (h_scale*a(i,j,k)*dt) ; enddo
232 endif
233 if (present(hv)) then
234 write(file,'(/,"hvel: ")', advance='no')
235 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hv(i,j,k) ; enddo
236 endif
237 if (present(str)) then
238 write(file,'(/,"Stress: ",ES10.3)', advance='no') (uh_scale*gv%RZ_to_H) * (str*dt)
239 endif
240
241 if (associated(cs%u_accel_bt)) then
242 write(file,'(/,"dubt: ")', advance='no')
243 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
244 (vel_scale*dt*cs%u_accel_bt(i,j,k)) ; enddo
245 endif
246 write(file,'(/)')
247
248 write(file,'(/,"h--: ")', advance='no')
249 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (h_scale*hin(i,j-1,k)) ; enddo
250 write(file,'(/,"h+-: ")', advance='no')
251 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (h_scale*hin(i+1,j-1,k)) ; enddo
252 write(file,'(/,"h-0: ")', advance='no')
253 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (h_scale*hin(i,j,k)) ; enddo
254 write(file,'(/,"h+0: ")', advance='no')
255 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (h_scale*hin(i+1,j,k)) ; enddo
256 write(file,'(/,"h-+: ")', advance='no')
257 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (h_scale*hin(i,j+1,k)) ; enddo
258 write(file,'(/,"h++: ")', advance='no')
259 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (h_scale*hin(i+1,j+1,k)) ; enddo
260
261
262 e(nz+1) = -us%Z_to_m*(g%bathyT(i,j) + g%Z_ref)
263 do k=nz,1,-1 ; e(k) = e(k+1) + h_scale*hin(i,j,k) ; enddo
264 write(file,'(/,"e-: ",ES10.3," ")', advance='no') e(ks)
265 do k=ks+1,ke+1 ; if (do_k(k-1)) write(file,'(ES10.3," ")', advance='no') e(k) ; enddo
266
267 e(nz+1) = -us%Z_to_m*(g%bathyT(i+1,j) + g%Z_ref)
268 do k=nz,1,-1 ; e(k) = e(k+1) + h_scale*hin(i+1,j,k) ; enddo
269 write(file,'(/,"e+: ",ES10.3," ")', advance='no') e(ks)
270 do k=ks+1,ke+1 ; if (do_k(k-1)) write(file,'(ES10.3," ")', advance='no') e(k) ; enddo
271 if (associated(cs%T)) then
272 write(file,'(/,"T-: ")', advance='no')
273 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') temp_scale*cs%T(i,j,k) ; enddo
274 write(file,'(/,"T+: ")', advance='no')
275 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') temp_scale*cs%T(i+1,j,k) ; enddo
276 endif
277 if (associated(cs%S)) then
278 write(file,'(/,"S-: ")', advance='no')
279 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') saln_scale*cs%S(i,j,k) ; enddo
280 write(file,'(/,"S+: ")', advance='no')
281 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') saln_scale*cs%S(i+1,j,k) ; enddo
282 endif
283
284 if (prev_avail) then
285 write(file,'(/,"v--: ")', advance='no')
286 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%v_prev(i,j-1,k)) ; enddo
287 write(file,'(/,"v-+: ")', advance='no')
288 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%v_prev(i,j,k)) ; enddo
289 write(file,'(/,"v+-: ")', advance='no')
290 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%v_prev(i+1,j-1,k)) ; enddo
291 write(file,'(/,"v++: ")', advance='no')
292 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%v_prev(i+1,j,k)) ; enddo
293 endif
294
295 write(file,'(/,"vh--: ")', advance='no')
296 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
297 (uh_scale*cdp%vh(i,j-1,k)*g%IdxCv(i,j-1)) ; enddo
298 write(file,'(/," vhC--:")', advance='no')
299 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
300 (0.5*cs%v_av(i,j-1,k)*uh_scale*(hin(i,j-1,k) + hin(i,j,k))) ; enddo
301 if (prev_avail) then
302 write(file,'(/," vhCp--:")', advance='no')
303 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
304 (0.5*cs%v_prev(i,j-1,k)*uh_scale*(hin(i,j-1,k) + hin(i,j,k))) ; enddo
305 endif
306
307 write(file,'(/,"vh-+: ")', advance='no')
308 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
309 (uh_scale*cdp%vh(i,j,k)*g%IdxCv(i,j)) ; enddo
310 write(file,'(/," vhC-+:")', advance='no')
311 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
312 (0.5*cs%v_av(i,j,k)*uh_scale*(hin(i,j,k) + hin(i,j+1,k))) ; enddo
313 if (prev_avail) then
314 write(file,'(/," vhCp-+:")', advance='no')
315 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
316 (0.5*cs%v_prev(i,j,k)*uh_scale*(hin(i,j,k) + hin(i,j+1,k))) ; enddo
317 endif
318
319 write(file,'(/,"vh+-: ")', advance='no')
320 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
321 (uh_scale*cdp%vh(i+1,j-1,k)*g%IdxCv(i+1,j-1)) ; enddo
322 write(file,'(/," vhC+-:")', advance='no')
323 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
324 (0.5*cs%v_av(i+1,j-1,k)*uh_scale*(hin(i+1,j-1,k) + hin(i+1,j,k))) ; enddo
325 if (prev_avail) then
326 write(file,'(/," vhCp+-:")', advance='no')
327 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
328 (0.5*cs%v_prev(i+1,j-1,k)*uh_scale*(hin(i+1,j-1,k) + hin(i+1,j,k))) ; enddo
329 endif
330
331 write(file,'(/,"vh++: ")', advance='no')
332 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
333 (uh_scale*cdp%vh(i+1,j,k)*g%IdxCv(i+1,j)) ; enddo
334 write(file,'(/," vhC++:")', advance='no')
335 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
336 (0.5*cs%v_av(i+1,j,k)*uh_scale*(hin(i+1,j,k) + hin(i+1,j+1,k))) ; enddo
337 if (prev_avail) then
338 write(file,'(/," vhCp++:")', advance='no')
339 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
340 (0.5*cs%v_av(i+1,j,k)*uh_scale*(hin(i+1,j,k) + hin(i+1,j+1,k))) ; enddo
341 endif
342
343 write(file,'(/,"D: ",2(ES10.3))') us%Z_to_m*(g%bathyT(i,j) + g%Z_ref), us%Z_to_m*(g%bathyT(i+1,j) + g%Z_ref)
344
345 ! From here on, the normalized accelerations are written.
346 if (prev_avail) then
347 do k=ks,ke
348 du = um(i,j,k) - cs%u_prev(i,j,k)
349 if (abs(du) < 1.0e-6*us%m_s_to_L_T) du = 1.0e-6*us%m_s_to_L_T
350 inorm(k) = 1.0 / du
351 enddo
352
353 write(file,'(2/,"Norm: ")', advance='no')
354 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') (vel_scale / inorm(k)) ; enddo
355
356 write(file,'(/,"du: ")', advance='no')
357 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
358 ((um(i,j,k)-cs%u_prev(i,j,k)) * inorm(k)) ; enddo
359
360 write(file,'(/,"CAu: ")', advance='no')
361 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
362 (dt*adp%CAu(i,j,k) * inorm(k)) ; enddo
363
364 write(file,'(/,"PFu: ")', advance='no')
365 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
366 (dt*adp%PFu(i,j,k) * inorm(k)) ; enddo
367
368 write(file,'(/,"diffu: ")', advance='no')
369 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
370 (dt*adp%diffu(i,j,k) * inorm(k)) ; enddo
371
372 if (associated(adp%gradKEu)) then
373 write(file,'(/,"KEu: ")', advance='no')
374 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
375 (dt*adp%gradKEu(i,j,k) * inorm(k)) ; enddo
376 endif
377 if (associated(adp%rv_x_v)) then
378 write(file,'(/,"Coru: ")', advance='no')
379 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
380 (dt*(adp%CAu(i,j,k)-adp%rv_x_v(i,j,k)) * inorm(k)) ; enddo
381 endif
382 if (associated(adp%du_dt_visc)) then
383 write(file,'(/,"duv: ")', advance='no')
384 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
385 (dt*adp%du_dt_visc(i,j,k) * inorm(k)) ; enddo
386 endif
387 if (associated(adp%du_other)) then
388 write(file,'(/,"du_other: ")', advance='no')
389 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
390 (adp%du_other(i,j,k) * inorm(k)) ; enddo
391 endif
392 if (associated(cs%u_accel_bt)) then
393 write(file,'(/,"dubt: ")', advance='no')
394 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
395 (dt*cs%u_accel_bt(i,j,k) * inorm(k)) ; enddo
396 endif
397 endif
398
399 write(file,'(2/)')
400
401 flush(file)
402 endif
403
404end subroutine write_u_accel
405
406!> This subroutine writes to an output file all of the accelerations
407!! that have been applied to a column of meridional velocities over
408!! the previous timestep. This subroutine is called from vertvisc.
409subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, str, a, hv)
410 integer, intent(in) :: i !< The zonal index of the column to be documented.
411 integer, intent(in) :: j !< The meridional index of the column to be documented.
412 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
413 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
414 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
415 real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), &
416 intent(in) :: vm !< The new meridional velocity [L T-1 ~> m s-1].
417 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
418 intent(in) :: hin !< The layer thickness [H ~> m or kg m-2].
419 type(accel_diag_ptrs), intent(in) :: adp !< A structure pointing to the various
420 !! accelerations in the momentum equations.
421 type(cont_diag_ptrs), intent(in) :: cdp !< A structure with pointers to various terms in
422 !! the continuity equations.
423 real, intent(in) :: dt !< The ocean dynamics time step [T ~> s].
424 type(pointaccel_cs), pointer :: cs !< The control structure returned by a previous
425 !! call to PointAccel_init.
426 real, intent(in) :: vel_rpt !< The velocity magnitude that triggers a report [L T-1 ~> m s-1]
427 real, optional, intent(in) :: str !< The surface wind stress [R L Z T-2 ~> Pa]
428 real, dimension(SZI_(G),SZJB_(G),SZK_(GV)+1), &
429 optional, intent(in) :: a !< The layer coupling coefficients from vertvisc
430 !! [H T-1 ~> m s-1 or Pa s m-1]
431 real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), &
432 optional, intent(in) :: hv !< The layer thicknesses at velocity grid points,
433 !! from vertvisc [H ~> m or kg m-2].
434
435 ! Local variables
436 real :: cfl ! The local velocity-based CFL number [nondim]
437 real :: angstrom ! A negligibly small thickness [H ~> m or kg m-2]
438 real :: dv ! A velocity change [L T-1 ~> m s-1]
439 real :: inorm(szk_(gv)) ! The inverse of the normalized velocity change [L T-1 ~> m s-1]
440 real :: e(szk_(gv)+1) ! Simple estimates of interface heights based on the sum of thicknesses [m]
441 real :: h_scale ! A scaling factor for thicknesses [m H-1 ~> 1] or [kg m-2 H-1 ~> 1]
442 real :: vel_scale ! A scaling factor for velocities [m T s-1 L-1 ~> 1]
443 real :: uh_scale ! A scaling factor for transport per unit length [m2 T s-1 L-1 H-1 ~> 1]
444 ! or [kg T m-1 s-1 L-1 H-1 ~> 1]
445 real :: temp_scale ! A scaling factor for temperatures [degC C-1 ~> 1]
446 real :: saln_scale ! A scaling factor for salinities [ppt S-1 ~> 1]
447 integer :: yr, mo, day, hr, minute, sec, yearday
448 integer :: k, ks, ke
449 integer :: nz
450 logical :: do_k(szk_(gv)+1)
451 logical :: prev_avail
452 integer :: file
453
454 angstrom = gv%Angstrom_H + gv%H_subroundoff
455 h_scale = gv%H_to_mks ; vel_scale = us%L_T_to_m_s ; uh_scale = h_scale*vel_scale
456 temp_scale = us%C_to_degC ; saln_scale = us%S_to_ppt
457
458! if (.not.associated(CS)) return
459 nz = gv%ke
460 if (cs%cols_written < cs%max_writes) then
461 cs%cols_written = cs%cols_written + 1
462
463 ks = 1 ; ke = nz
464 do_k(:) = .false.
465
466 ! Open up the file for output if this is the first call.
467 if (cs%v_file == -1) then
468 if (len_trim(cs%v_trunc_file) < 1) return
469 call open_ascii_file(cs%v_file, trim(cs%v_trunc_file), action=append_file, &
470 threading=multiple, fileset=single_file)
471 if (cs%v_file == -1) then
472 call mom_error(note, 'Unable to open file '//trim(cs%v_trunc_file)//'.')
473 return
474 endif
475 endif
476 file = cs%v_file
477
478 prev_avail = (associated(cs%u_prev) .and. associated(cs%v_prev))
479
480 do k=1,nz
481 if (((max(cs%v_av(i,j,k), vm(i,j,k)) >= vel_rpt) .or. &
482 (min(cs%v_av(i,j,k), vm(i,j,k)) <= -vel_rpt)) .and. &
483 ((hin(i,j,k) + hin(i,j+1,k)) > 3.0*angstrom)) exit
484 enddo
485 ks = k
486 do k=nz,1,-1
487 if (((max(cs%v_av(i,j,k), vm(i,j,k)) >= vel_rpt) .or. &
488 (min(cs%v_av(i,j,k), vm(i,j,k)) <= -vel_rpt)) .and. &
489 ((hin(i,j,k) + hin(i,j+1,k)) > 3.0*angstrom)) exit
490 enddo
491 ke = k
492 if (ke < ks) then
493 ks = 1 ; ke = nz ; write(file,'("V: Unable to set ks & ke.")')
494 endif
495 if (cs%full_column) then
496 ks = 1 ; ke = nz
497 endif
498
499 call get_date(cs%Time, yr, mo, day, hr, minute, sec)
500 call get_time((cs%Time - set_date(yr, 1, 1, 0, 0, 0)), sec, yearday)
501 write (file,'(/,"--------------------------")')
502 write (file,'(/,"Time ",I0," ",I0," ",F6.2," V-velocity violation at ",I0,": ",I0,", ",I0, &
503 & " (",F7.2," E ",F7.2," N) Layers ",I0," to ",I0,". dt = ",1PG10.4)') &
504 yr, yearday, (real(sec)/3600.0), pe_here(), i, j, &
505 g%geoLonCv(i,j), g%geoLatCv(i,j), ks, ke, us%T_to_s*dt
506
507 if (ks <= gv%nk_rho_varies) ks = 1
508 do k=ks,ke
509 if ((hin(i,j,k) + hin(i,j+1,k)) > 3.0*angstrom) do_k(k) = .true.
510 enddo
511
512 write(file,'(/,"Layers:")', advance='no')
513 do k=ks,ke ; if (do_k(k)) write(file,'(I10," ")', advance='no') (k) ; enddo
514 write(file,'(/,"v(m): ")', advance='no')
515 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*vm(i,j,k)) ; enddo
516
517 if (prev_avail) then
518 write(file,'(/,"v(mp): ")', advance='no')
519 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%v_prev(i,j,k)) ; enddo
520 endif
521
522 write(file,'(/,"v(3): ")', advance='no')
523 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*cs%v_av(i,j,k)) ; enddo
524 write(file,'(/,"CFL v: ")', advance='no')
525 do k=ks,ke ; if (do_k(k)) then
526 cfl = abs(vm(i,j,k)) * dt * g%dx_Cv(i,j)
527 if (vm(i,j,k) < 0.0) then ; cfl = cfl * g%IareaT(i,j+1)
528 else ; cfl = cfl * g%IareaT(i,j) ; endif
529 write(file,'(ES10.3," ")', advance='no') cfl
530 endif ; enddo
531 write(file,'(/,"CFL0 v:")', advance='no')
532 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
533 abs(vm(i,j,k)) * dt * g%IdyCv(i,j) ; enddo
534
535 if (prev_avail) then
536 write(file,'(/,"dv: ")', advance='no')
537 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
538 (vel_scale*(vm(i,j,k)-cs%v_prev(i,j,k))) ; enddo
539 endif
540
541 write(file,'(/,"CAv: ")', advance='no')
542 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*dt*adp%CAv(i,j,k)) ; enddo
543
544 write(file,'(/,"PFv: ")', advance='no')
545 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*dt*adp%PFv(i,j,k)) ; enddo
546
547 write(file,'(/,"diffv: ")', advance='no')
548 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') (vel_scale*dt*adp%diffv(i,j,k)) ; enddo
549
550 if (associated(adp%gradKEv)) then
551 write(file,'(/,"KEv: ")', advance='no')
552 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
553 (vel_scale*dt*adp%gradKEv(i,j,k)) ; enddo
554 endif
555 if (associated(adp%rv_x_u)) then
556 write(file,'(/,"Corv: ")', advance='no')
557 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
558 vel_scale*dt*(adp%CAv(i,j,k)-adp%rv_x_u(i,j,k)) ; enddo
559 endif
560 if (associated(adp%dv_dt_visc)) then
561 write(file,'(/,"vbv: ")', advance='no')
562 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
563 vel_scale*(vm(i,j,k) - dt*adp%dv_dt_visc(i,j,k)) ; enddo
564
565 write(file,'(/,"dvv: ")', advance='no')
566 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
567 (vel_scale*dt*adp%dv_dt_visc(i,j,k)) ; enddo
568 endif
569 if (associated(adp%dv_other)) then
570 write(file,'(/,"dv_other: ")', advance='no')
571 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
572 (vel_scale*adp%dv_other(i,j,k)) ; enddo
573 endif
574 if (present(a)) then
575 write(file,'(/,"a: ",ES10.3," ")', advance='no') h_scale*a(i,j,ks)*dt
576 do k=ks+1,ke+1 ; if (do_k(k-1)) write(file,'(ES10.3," ")', advance='no') (h_scale*a(i,j,k)*dt) ; enddo
577 endif
578 if (present(hv)) then
579 write(file,'(/,"hvel: ")', advance='no')
580 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hv(i,j,k) ; enddo
581 endif
582 if (present(str)) then
583 write(file,'(/,"Stress: ",ES10.3)', advance='no') (uh_scale*gv%RZ_to_H) * (str*dt)
584 endif
585
586 if (associated(cs%v_accel_bt)) then
587 write(file,'("dvbt: ")', advance='no')
588 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
589 (vel_scale*dt*cs%v_accel_bt(i,j,k)) ; enddo
590 endif
591 write(file,'(/)')
592
593 write(file,'("h--: ")', advance='no')
594 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hin(i-1,j,k) ; enddo
595 write(file,'(/,"h0-: ")', advance='no')
596 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hin(i,j,k) ; enddo
597 write(file,'(/,"h+-: ")', advance='no')
598 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hin(i+1,j,k) ; enddo
599 write(file,'(/,"h-+: ")', advance='no')
600 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hin(i-1,j+1,k) ; enddo
601 write(file,'(/,"h0+: ")', advance='no')
602 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hin(i,j+1,k) ; enddo
603 write(file,'(/,"h++: ")', advance='no')
604 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hin(i+1,j+1,k) ; enddo
605
606 e(nz+1) = -us%Z_to_m*(g%bathyT(i,j) + g%Z_ref)
607 do k=nz,1,-1 ; e(k) = e(k+1) + h_scale*hin(i,j,k) ; enddo
608 write(file,'(/,"e-: ",ES10.3," ")', advance='no') e(ks)
609 do k=ks+1,ke+1 ; if (do_k(k-1)) write(file,'(ES10.3," ")', advance='no') e(k) ; enddo
610
611 e(nz+1) = -us%Z_to_m*(g%bathyT(i,j+1) + g%Z_ref)
612 do k=nz,1,-1 ; e(k) = e(k+1) + h_scale*hin(i,j+1,k) ; enddo
613 write(file,'(/,"e+: ",ES10.3," ")', advance='no') e(ks)
614 do k=ks+1,ke+1 ; if (do_k(k-1)) write(file,'(ES10.3," ")', advance='no') e(k) ; enddo
615 if (associated(cs%T)) then
616 write(file,'(/,"T-: ")', advance='no')
617 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') temp_scale*cs%T(i,j,k) ; enddo
618 write(file,'(/,"T+: ")', advance='no')
619 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') temp_scale*cs%T(i,j+1,k) ; enddo
620 endif
621 if (associated(cs%S)) then
622 write(file,'(/,"S-: ")', advance='no')
623 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') saln_scale*cs%S(i,j,k) ; enddo
624 write(file,'(/,"S+: ")', advance='no')
625 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') saln_scale*cs%S(i,j+1,k) ; enddo
626 endif
627
628 if (prev_avail) then
629 write(file,'(/,"u--: ")', advance='no')
630 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') vel_scale*cs%u_prev(i-1,j,k) ; enddo
631 write(file,'(/,"u-+: ")', advance='no')
632 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') vel_scale*cs%u_prev(i-1,j+1,k) ; enddo
633 write(file,'(/,"u+-: ")', advance='no')
634 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') vel_scale*cs%u_prev(i,j,k) ; enddo
635 write(file,'(/,"u++: ")', advance='no')
636 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') vel_scale*cs%u_prev(i,j+1,k) ; enddo
637 endif
638
639 write(file,'(/,"uh--: ")', advance='no')
640 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
641 (uh_scale*cdp%uh(i-1,j,k)*g%IdyCu(i-1,j)) ; enddo
642 write(file,'(/," uhC--: ")', advance='no')
643 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
644 (cs%u_av(i-1,j,k) * uh_scale*0.5*(hin(i-1,j,k) + hin(i,j,k))) ; enddo
645 if (prev_avail) then
646 write(file,'(/," uhCp--:")', advance='no')
647 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
648 (cs%u_prev(i-1,j,k) * uh_scale*0.5*(hin(i-1,j,k) + hin(i,j,k))) ; enddo
649 endif
650
651 write(file,'(/,"uh-+: ")', advance='no')
652 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
653 (uh_scale*cdp%uh(i-1,j+1,k)*g%IdyCu(i-1,j+1)) ; enddo
654 write(file,'(/," uhC-+: ")', advance='no')
655 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
656 (cs%u_av(i-1,j+1,k) * uh_scale*0.5*(hin(i-1,j+1,k) + hin(i,j+1,k))) ; enddo
657 if (prev_avail) then
658 write(file,'(/," uhCp-+:")', advance='no')
659 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
660 (cs%u_prev(i-1,j+1,k) * uh_scale*0.5*(hin(i-1,j+1,k) + hin(i,j+1,k))) ; enddo
661 endif
662
663 write(file,'(/,"uh+-: ")', advance='no')
664 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
665 (uh_scale*cdp%uh(i,j,k)*g%IdyCu(i,j)) ; enddo
666 write(file,'(/," uhC+-: ")', advance='no')
667 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
668 (cs%u_av(i,j,k) * uh_scale*0.5*(hin(i,j,k) + hin(i+1,j,k))) ; enddo
669 if (prev_avail) then
670 write(file,'(/," uhCp+-:")', advance='no')
671 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
672 (cs%u_prev(i,j,k) * uh_scale*0.5*(hin(i,j,k) + hin(i+1,j,k))) ; enddo
673 endif
674
675 write(file,'(/,"uh++: ")', advance='no')
676 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
677 (uh_scale*cdp%uh(i,j+1,k)*g%IdyCu(i,j+1)) ; enddo
678 write(file,'(/," uhC++: ")', advance='no')
679 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
680 (cs%u_av(i,j+1,k) * uh_scale*0.5*(hin(i,j+1,k) + hin(i+1,j+1,k))) ; enddo
681 if (prev_avail) then
682 write(file,'(/," uhCp++:")', advance='no')
683 do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') &
684 (cs%u_prev(i,j+1,k) * uh_scale*0.5*(hin(i,j+1,k) + hin(i+1,j+1,k))) ; enddo
685 endif
686
687 write(file,'(/,"D: ",2(ES10.3))') us%Z_to_m*(g%bathyT(i,j) + g%Z_ref), us%Z_to_m*(g%bathyT(i,j+1) + g%Z_ref)
688
689 ! From here on, the normalized accelerations are written.
690 if (prev_avail) then
691 do k=ks,ke
692 dv = vm(i,j,k) - cs%v_prev(i,j,k)
693 if (abs(dv) < 1.0e-6*us%m_s_to_L_T) dv = 1.0e-6*us%m_s_to_L_T
694 inorm(k) = 1.0 / dv
695 enddo
696
697 write(file,'(2/,"Norm: ")', advance='no')
698 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') (vel_scale / inorm(k)) ; enddo
699 write(file,'(/,"dv: ")', advance='no')
700 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
701 ((vm(i,j,k)-cs%v_prev(i,j,k)) * inorm(k)) ; enddo
702 write(file,'(/,"CAv: ")', advance='no')
703 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
704 (dt*adp%CAv(i,j,k) * inorm(k)) ; enddo
705 write(file,'(/,"PFv: ")', advance='no')
706 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
707 (dt*adp%PFv(i,j,k) * inorm(k)) ; enddo
708 write(file,'(/,"diffv: ")', advance='no')
709 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
710 (dt*adp%diffv(i,j,k) * inorm(k)) ; enddo
711
712 if (associated(adp%gradKEu)) then
713 write(file,'(/,"KEv: ")', advance='no')
714 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
715 (dt*adp%gradKEv(i,j,k) * inorm(k)) ; enddo
716 endif
717 if (associated(adp%rv_x_u)) then
718 write(file,'(/,"Corv: ")', advance='no')
719 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
720 (dt*(adp%CAv(i,j,k)-adp%rv_x_u(i,j,k)) * inorm(k)) ; enddo
721 endif
722 if (associated(adp%dv_dt_visc)) then
723 write(file,'(/,"dvv: ")', advance='no')
724 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
725 (dt*adp%dv_dt_visc(i,j,k) * inorm(k)) ; enddo
726 endif
727 if (associated(adp%dv_other)) then
728 write(file,'(/,"dv_other: ")', advance='no')
729 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
730 (adp%dv_other(i,j,k) * inorm(k)) ; enddo
731 endif
732 if (associated(cs%v_accel_bt)) then
733 write(file,'(/,"dvbt: ")', advance='no')
734 do k=ks,ke ; if (do_k(k)) write(file,'(F10.6," ")', advance='no') &
735 (dt*cs%v_accel_bt(i,j,k) * inorm(k)) ; enddo
736 endif
737 endif
738
739 write(file,'(2/)')
740
741 flush(file)
742 endif
743
744end subroutine write_v_accel
745
746!> This subroutine initializes the parameters regulating how truncations are logged.
747subroutine pointaccel_init(MIS, Time, G, param_file, diag, dirs, CS)
748 type(ocean_internal_state), &
749 target, intent(in) :: mis !< For "MOM Internal State" a set of pointers
750 !! to the fields and accelerations that make
751 !! up the ocean's physical state.
752 type(time_type), target, intent(in) :: time !< The current model time.
753 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
754 type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time
755 !! parameters.
756 type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to regulate
757 !! diagnostic output.
758 type(directories), intent(in) :: dirs !< A structure containing several relevant
759 !! directory paths.
760 type(pointaccel_cs), pointer :: cs !< A pointer that is set to point to the
761 !! control structure for this module.
762 ! This include declares and sets the variable "version".
763# include "version_variable.h"
764 character(len=40) :: mdl = "MOM_PointAccel" ! This module's name.
765
766 if (associated(cs)) return
767 allocate(cs)
768
769 cs%diag => diag ; cs%Time => time
770
771 cs%T => mis%T ; cs%S => mis%S
772 cs%u_accel_bt => mis%u_accel_bt ; cs%v_accel_bt => mis%v_accel_bt
773 cs%u_prev => mis%u_prev ; cs%v_prev => mis%v_prev
774 cs%u_av => mis%u_av ; if (.not.associated(mis%u_av)) cs%u_av => mis%u(:,:,:)
775 cs%v_av => mis%v_av ; if (.not.associated(mis%v_av)) cs%v_av => mis%v(:,:,:)
776
777 ! Read all relevant parameters and write them to the model log.
778 call log_version(param_file, mdl, version, "", debugging=.true.)
779 call get_param(param_file, mdl, "U_TRUNC_FILE", cs%u_trunc_file, &
780 "The absolute path to the file where the accelerations "//&
781 "leading to zonal velocity truncations are written. \n"//&
782 "Leave this empty for efficiency if this diagnostic is "//&
783 "not needed.", default="", debuggingparam=.true.)
784 call get_param(param_file, mdl, "V_TRUNC_FILE", cs%v_trunc_file, &
785 "The absolute path to the file where the accelerations "//&
786 "leading to meridional velocity truncations are written. \n"//&
787 "Leave this empty for efficiency if this diagnostic is "//&
788 "not needed.", default="", debuggingparam=.true.)
789 call get_param(param_file, mdl, "MAX_TRUNC_FILE_SIZE_PER_PE", cs%max_writes, &
790 "The maximum number of columns of truncations that any PE "//&
791 "will write out during a run.", default=50, debuggingparam=.true.)
792 call get_param(param_file, mdl, "DEBUG_FULL_COLUMN", cs%full_column, &
793 "If true, write out the accelerations in all massive layers; otherwise "//&
794 "just document the ones with large velocities.", &
795 default=.false., debuggingparam=.true.)
796
797 if (len_trim(dirs%output_directory) > 0) then
798 if (len_trim(cs%u_trunc_file) > 0) &
799 cs%u_trunc_file = trim(dirs%output_directory)//trim(cs%u_trunc_file)
800 if (len_trim(cs%v_trunc_file) > 0) &
801 cs%v_trunc_file = trim(dirs%output_directory)//trim(cs%v_trunc_file)
802 call log_param(param_file, mdl, "output_dir/U_TRUNC_FILE", cs%u_trunc_file, debuggingparam=.true.)
803 call log_param(param_file, mdl, "output_dir/V_TRUNC_FILE", cs%v_trunc_file, debuggingparam=.true.)
804 endif
805 cs%u_file = -1 ; cs%v_file = -1 ; cs%cols_written = 0
806
807end subroutine pointaccel_init
808
809end module mom_pointaccel