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