MOM6
MOM_debugging.F90
Go to the documentation of this file.
1 !> Provides checksumming functions for debugging
2 !!
3 !! This module contains subroutines that perform various error checking and
4 !! debugging functions for MOM6. This routine is similar to it counterpart in
5 !! the SIS2 code, except for the use of the ocean_grid_type and by keeping them
6 !! separate we retain the ability to set up MOM6 and SIS2 debugging separately.
8 
9 ! This file is part of MOM6. See LICENSE.md for the license.
10 
13 use mom_coms, only : pe_here, root_pe, num_pes, sum_across_pes
14 use mom_coms, only : min_across_pes, max_across_pes, reproducing_sum
15 use mom_domains, only : pass_vector, pass_var, pe_here
16 use mom_domains, only : bgrid_ne, agrid, to_all, scalar_pair
17 use mom_error_handler, only : mom_error, fatal, warning, is_root_pe
19 use mom_grid, only : ocean_grid_type
20 use mom_hor_index, only : hor_index_type
21 
22 implicit none ; private
23 
28 
29 ! These interfaces come from MOM_checksums.
31 
32 !> Check for consistency between the duplicated points of a C-grid vector
33 interface check_redundant
35 end interface check_redundant
36 !> Check for consistency between the duplicated points of a C-grid vector
39 end interface check_redundant_c
40 !> Check for consistency between the duplicated points of a B-grid vector or scalar
44 end interface check_redundant_b
45 !> Check for consistency between the duplicated points of an A-grid vector or scalar
49 end interface check_redundant_t
50 
51 !> Do checksums on the components of a C-grid vector
52 interface vec_chksum
53  module procedure chksum_vec_c3d, chksum_vec_c2d
54 end interface vec_chksum
55 !> Do checksums on the components of a C-grid vector
56 interface vec_chksum_c
57  module procedure chksum_vec_c3d, chksum_vec_c2d
58 end interface vec_chksum_c
59 !> Do checksums on the components of a B-grid vector
60 interface vec_chksum_b
61  module procedure chksum_vec_b3d, chksum_vec_b2d
62 end interface vec_chksum_b
63 !> Do checksums on the components of an A-grid vector
64 interface vec_chksum_a
65  module procedure chksum_vec_a3d, chksum_vec_a2d
66 end interface vec_chksum_a
67 
68 ! Note: these parameters are module data but ONLY used when debugging and
69 ! so can violate the thread-safe requirement of no module/global data.
70 integer :: max_redundant_prints = 100 !< Maximum number of times to write redundant messages
71 integer :: redundant_prints(3) = 0 !< Counters for controlling redundant printing
72 logical :: debug = .false. !< Write out verbose debugging data
73 logical :: debug_chksums = .true. !< Perform checksums on arrays
74 logical :: debug_redundant = .true. !< Check redundant values on PE boundaries
75 
76 contains
77 
78 !> MOM_debugging_init initializes the MOM_debugging module, and sets
79 !! the parameterts that control which checks are active for MOM6.
80 subroutine mom_debugging_init(param_file)
81  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
82 ! This include declares and sets the variable "version".
83 #include "version_variable.h"
84  character(len=40) :: mdl = "MOM_debugging" ! This module's name.
85 
86  call log_version(param_file, mdl, version)
87  call get_param(param_file, mdl, "DEBUG", debug, &
88  "If true, write out verbose debugging data.", &
89  default=.false., debuggingparam=.true.)
90  call get_param(param_file, mdl, "DEBUG_CHKSUMS", debug_chksums, &
91  "If true, checksums are performed on arrays in the "//&
92  "various vec_chksum routines.", default=debug, &
93  debuggingparam=.true.)
94  call get_param(param_file, mdl, "DEBUG_REDUNDANT", debug_redundant, &
95  "If true, debug redundant data points during calls to "//&
96  "the various vec_chksum routines.", default=debug, &
97  debuggingparam=.true.)
98 
99  call mom_checksums_init(param_file)
100 
101 end subroutine mom_debugging_init
102 
103 !> Check for consistency between the duplicated points of a 3-D C-grid vector
104 subroutine check_redundant_vc3d(mesg, u_comp, v_comp, G, is, ie, js, je, &
105  direction)
106  character(len=*), intent(in) :: mesg !< An identifying message
107  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
108  real, dimension(G%IsdB:,G%jsd:,:), intent(in) :: u_comp !< The u-component of the vector
109  !! to be checked for consistency
110  real, dimension(G%isd:,G%JsdB:,:), intent(in) :: v_comp !< The u-component of the vector
111  !! to be checked for consistency
112  integer, optional, intent(in) :: is !< The starting i-index to check
113  integer, optional, intent(in) :: ie !< The ending i-index to check
114  integer, optional, intent(in) :: js !< The starting j-index to check
115  integer, optional, intent(in) :: je !< The ending j-index to check
116  integer, optional, intent(in) :: direction !< the direction flag to be
117  !! passed to pass_vector
118  ! Local variables
119  character(len=24) :: mesg_k
120  integer :: k
121 
122  do k=1,size(u_comp,3)
123  if (k < 10) then ; write(mesg_k,'(" Layer",i2," ")') k
124  elseif (k < 100) then ; write(mesg_k,'(" Layer",i3," ")') k
125  elseif (k < 1000) then ; write(mesg_k,'(" Layer",i4," ")') k
126  else ; write(mesg_k,'(" Layer",i9," ")') k ; endif
127 
128  call check_redundant_vc2d(trim(mesg)//trim(mesg_k), u_comp(:,:,k), &
129  v_comp(:,:,k), g, is, ie, js, je, direction)
130  enddo
131 end subroutine check_redundant_vc3d
132 
133 !> Check for consistency between the duplicated points of a 2-D C-grid vector
134 subroutine check_redundant_vc2d(mesg, u_comp, v_comp, G, is, ie, js, je, &
135  direction)
136  character(len=*), intent(in) :: mesg !< An identifying message
137  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
138  real, dimension(G%IsdB:,G%jsd:), intent(in) :: u_comp !< The u-component of the vector
139  !! to be checked for consistency
140  real, dimension(G%isd:,G%JsdB:), intent(in) :: v_comp !< The u-component of the vector
141  !! to be checked for consistency
142  integer, optional, intent(in) :: is !< The starting i-index to check
143  integer, optional, intent(in) :: ie !< The ending i-index to check
144  integer, optional, intent(in) :: js !< The starting j-index to check
145  integer, optional, intent(in) :: je !< The ending j-index to check
146  integer, optional, intent(in) :: direction !< the direction flag to be
147  !! passed to pass_vector
148  ! Local variables
149  real :: u_nonsym(G%isd:G%ied,G%jsd:G%jed)
150  real :: v_nonsym(G%isd:G%ied,G%jsd:G%jed)
151  real :: u_resym(G%IsdB:G%IedB,G%jsd:G%jed)
152  real :: v_resym(G%isd:G%ied,G%JsdB:G%JedB)
153  character(len=128) :: mesg2
154  integer :: i, j, is_ch, ie_ch, js_ch, je_ch
155  integer :: Isq, Ieq, Jsq, Jeq, isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
156 
157  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
158  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
159  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
160 
161  if (.not.(present(is) .or. present(ie) .or. present(js) .or. present(je))) then
162  ! This only works with symmetric memory, so otherwise return.
163  if ((isd == isdb) .and. (jsd == jsdb)) return
164  endif
165 
166  do i=isd,ied ; do j=jsd,jed
167  u_nonsym(i,j) = u_comp(i,j) ; v_nonsym(i,j) = v_comp(i,j)
168  enddo ; enddo
169 
170  if (.not.associated(g%Domain_aux)) call mom_error(fatal," check_redundant"//&
171  " called with a non-associated auxiliary domain the grid type.")
172  call pass_vector(u_nonsym, v_nonsym, g%Domain_aux, direction)
173 
174  do i=isdb,iedb ; do j=jsd,jed ; u_resym(i,j) = u_comp(i,j) ; enddo ; enddo
175  do i=isd,ied ; do j=jsdb,jedb ; v_resym(i,j) = v_comp(i,j) ; enddo ; enddo
176  do i=isd,ied ; do j=jsd,jed
177  u_resym(i,j) = u_nonsym(i,j) ; v_resym(i,j) = v_nonsym(i,j)
178  enddo ; enddo
179  call pass_vector(u_resym, v_resym, g%Domain, direction)
180 
181  is_ch = isq ; ie_ch = ieq ; js_ch = jsq ; je_ch = jeq
182  if (present(is)) is_ch = is ; if (present(ie)) ie_ch = ie
183  if (present(js)) js_ch = js ; if (present(js)) je_ch = je
184 
185  do i=is_ch,ie_ch ; do j=js_ch+1,je_ch
186  if (u_resym(i,j) /= u_comp(i,j) .and. &
188  write(mesg2,'(" redundant u-components",2(1pe12.4)," differ by ", &
189  & 1pe12.4," at i,j = ",2i4," on pe ",i4)') &
190  u_comp(i,j), u_resym(i,j),u_comp(i,j)-u_resym(i,j),i,j,pe_here()
191  write(0,'(A130)') trim(mesg)//trim(mesg2)
193  endif
194  enddo ; enddo
195  do i=is_ch+1,ie_ch ; do j=js_ch,je_ch
196  if (v_resym(i,j) /= v_comp(i,j) .and. &
198  write(mesg2,'(" redundant v-comps",2(1pe12.4)," differ by ", &
199  & 1pe12.4," at i,j = ",2i4," x,y = ",2(1pe12.4)" on pe ",i4)') &
200  v_comp(i,j), v_resym(i,j),v_comp(i,j)-v_resym(i,j),i,j, &
201  g%geoLonBu(i,j), g%geoLatBu(i,j), pe_here()
202  write(0,'(A155)') trim(mesg)//trim(mesg2)
204  endif
205  enddo ; enddo
206 
207 end subroutine check_redundant_vc2d
208 
209 !> Check for consistency between the duplicated points of a 3-D scalar at corner points
210 subroutine check_redundant_sb3d(mesg, array, G, is, ie, js, je)
211  character(len=*), intent(in) :: mesg !< An identifying message
212  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
213  real, dimension(G%IsdB:,G%JsdB:,:), intent(in) :: array !< The array to be checked for consistency
214  integer, optional, intent(in) :: is !< The starting i-index to check
215  integer, optional, intent(in) :: ie !< The ending i-index to check
216  integer, optional, intent(in) :: js !< The starting j-index to check
217  integer, optional, intent(in) :: je !< The ending j-index to check
218 
219  ! Local variables
220  character(len=24) :: mesg_k
221  integer :: k
222 
223  do k=1,size(array,3)
224  if (k < 10) then ; write(mesg_k,'(" Layer",i2," ")') k
225  elseif (k < 100) then ; write(mesg_k,'(" Layer",i3," ")') k
226  elseif (k < 1000) then ; write(mesg_k,'(" Layer",i4," ")') k
227  else ; write(mesg_k,'(" Layer",i9," ")') k ; endif
228 
229  call check_redundant_sb2d(trim(mesg)//trim(mesg_k), array(:,:,k), &
230  g, is, ie, js, je)
231  enddo
232 end subroutine check_redundant_sb3d
233 
234 !> Check for consistency between the duplicated points of a 2-D scalar at corner points
235 subroutine check_redundant_sb2d(mesg, array, G, is, ie, js, je)
236  character(len=*), intent(in) :: mesg !< An identifying message
237  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
238  real, dimension(G%IsdB:,G%JsdB:), intent(in) :: array !< The array to be checked for consistency
239  integer, optional, intent(in) :: is !< The starting i-index to check
240  integer, optional, intent(in) :: ie !< The ending i-index to check
241  integer, optional, intent(in) :: js !< The starting j-index to check
242  integer, optional, intent(in) :: je !< The ending j-index to check
243  ! Local variables
244  real :: a_nonsym(G%isd:G%ied,G%jsd:G%jed)
245  real :: a_resym(G%IsdB:G%IedB,G%JsdB:G%JedB)
246  character(len=128) :: mesg2
247  integer :: i, j, is_ch, ie_ch, js_ch, je_ch
248  integer :: Isq, Ieq, Jsq, Jeq, isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
249 
250  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
251  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
252  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
253 
254  if (.not.(present(is) .or. present(ie) .or. present(js) .or. present(je))) then
255  ! This only works with symmetric memory, so otherwise return.
256  if ((isd == isdb) .and. (jsd == jsdb)) return
257  endif
258 
259  do i=isd,ied ; do j=jsd,jed
260  a_nonsym(i,j) = array(i,j)
261  enddo ; enddo
262 
263  if (.not.associated(g%Domain_aux)) call mom_error(fatal," check_redundant"//&
264  " called with a non-associated auxiliary domain the grid type.")
265  call pass_vector(a_nonsym, a_nonsym, g%Domain_aux, &
266  direction=to_all+scalar_pair, stagger=bgrid_ne)
267 
268  do i=isdb,iedb ; do j=jsdb,jedb ; a_resym(i,j) = array(i,j) ; enddo ; enddo
269  do i=isd,ied ; do j=jsd,jed
270  a_resym(i,j) = a_nonsym(i,j)
271  enddo ; enddo
272  call pass_vector(a_resym, a_resym, g%Domain, direction=to_all+scalar_pair, &
273  stagger=bgrid_ne)
274 
275  is_ch = isq ; ie_ch = ieq ; js_ch = jsq ; je_ch = jeq
276  if (present(is)) is_ch = is ; if (present(ie)) ie_ch = ie
277  if (present(js)) js_ch = js ; if (present(js)) je_ch = je
278 
279  do i=is_ch,ie_ch ; do j=js_ch,je_ch
280  if (a_resym(i,j) /= array(i,j) .and. &
282  write(mesg2,'(" Redundant points",2(1pe12.4)," differ by ", &
283  & 1pe12.4," at i,j = ",2i4," on pe ",i4)') &
284  array(i,j), a_resym(i,j),array(i,j)-a_resym(i,j),i,j,pe_here()
285  write(0,'(A130)') trim(mesg)//trim(mesg2)
287  endif
288  enddo ; enddo
289 
290 end subroutine check_redundant_sb2d
291 
292 !> Check for consistency between the duplicated points of a 3-D B-grid vector
293 subroutine check_redundant_vb3d(mesg, u_comp, v_comp, G, is, ie, js, je, &
294  direction)
295  character(len=*), intent(in) :: mesg !< An identifying message
296  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
297  real, dimension(G%IsdB:,G%JsdB:,:), intent(in) :: u_comp !< The u-component of the vector
298  !! to be checked for consistency
299  real, dimension(G%IsdB:,G%JsdB:,:), intent(in) :: v_comp !< The v-component of the vector
300  !! to be checked for consistency
301  integer, optional, intent(in) :: is !< The starting i-index to check
302  integer, optional, intent(in) :: ie !< The ending i-index to check
303  integer, optional, intent(in) :: js !< The starting j-index to check
304  integer, optional, intent(in) :: je !< The ending j-index to check
305  integer, optional, intent(in) :: direction !< the direction flag to be
306  !! passed to pass_vector
307  ! Local variables
308  character(len=24) :: mesg_k
309  integer :: k
310 
311  do k=1,size(u_comp,3)
312  if (k < 10) then ; write(mesg_k,'(" Layer",i2," ")') k
313  elseif (k < 100) then ; write(mesg_k,'(" Layer",i3," ")') k
314  elseif (k < 1000) then ; write(mesg_k,'(" Layer",i4," ")') k
315  else ; write(mesg_k,'(" Layer",i9," ")') k ; endif
316 
317  call check_redundant_vb2d(trim(mesg)//trim(mesg_k), u_comp(:,:,k), &
318  v_comp(:,:,k), g, is, ie, js, je, direction)
319  enddo
320 end subroutine check_redundant_vb3d
321 
322 !> Check for consistency between the duplicated points of a 2-D B-grid vector
323 subroutine check_redundant_vb2d(mesg, u_comp, v_comp, G, is, ie, js, je, &
324  direction)
325  character(len=*), intent(in) :: mesg !< An identifying message
326  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
327  real, dimension(G%IsdB:,G%JsdB:), intent(in) :: u_comp !< The u-component of the vector
328  !! to be checked for consistency
329  real, dimension(G%IsdB:,G%JsdB:), intent(in) :: v_comp !< The v-component of the vector
330  !! to be checked for consistency
331  integer, optional, intent(in) :: is !< The starting i-index to check
332  integer, optional, intent(in) :: ie !< The ending i-index to check
333  integer, optional, intent(in) :: js !< The starting j-index to check
334  integer, optional, intent(in) :: je !< The ending j-index to check
335  integer, optional, intent(in) :: direction !< the direction flag to be
336  !! passed to pass_vector
337  ! Local variables
338  real :: u_nonsym(G%isd:G%ied,G%jsd:G%jed)
339  real :: v_nonsym(G%isd:G%ied,G%jsd:G%jed)
340  real :: u_resym(G%IsdB:G%IedB,G%JsdB:G%JedB)
341  real :: v_resym(G%IsdB:G%IedB,G%JsdB:G%JedB)
342  character(len=128) :: mesg2
343  integer :: i, j, is_ch, ie_ch, js_ch, je_ch
344  integer :: Isq, Ieq, Jsq, Jeq, isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
345 
346  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
347  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
348  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
349 
350  if (.not.(present(is) .or. present(ie) .or. present(js) .or. present(je))) then
351  ! This only works with symmetric memory, so otherwise return.
352  if ((isd == isdb) .and. (jsd == jsdb)) return
353  endif
354 
355  do i=isd,ied ; do j=jsd,jed
356  u_nonsym(i,j) = u_comp(i,j) ; v_nonsym(i,j) = v_comp(i,j)
357  enddo ; enddo
358 
359  if (.not.associated(g%Domain_aux)) call mom_error(fatal," check_redundant"//&
360  " called with a non-associated auxiliary domain the grid type.")
361  call pass_vector(u_nonsym, v_nonsym, g%Domain_aux, direction, stagger=bgrid_ne)
362 
363  do i=isdb,iedb ; do j=jsdb,jedb
364  u_resym(i,j) = u_comp(i,j) ; v_resym(i,j) = v_comp(i,j)
365  enddo ; enddo
366  do i=isd,ied ; do j=jsd,jed
367  u_resym(i,j) = u_nonsym(i,j) ; v_resym(i,j) = v_nonsym(i,j)
368  enddo ; enddo
369  call pass_vector(u_resym, v_resym, g%Domain, direction, stagger=bgrid_ne)
370 
371  is_ch = isq ; ie_ch = ieq ; js_ch = jsq ; je_ch = jeq
372  if (present(is)) is_ch = is ; if (present(ie)) ie_ch = ie
373  if (present(js)) js_ch = js ; if (present(js)) je_ch = je
374 
375  do i=is_ch,ie_ch ; do j=js_ch,je_ch
376  if (u_resym(i,j) /= u_comp(i,j) .and. &
378  write(mesg2,'(" redundant u-components",2(1pe12.4)," differ by ", &
379  & 1pe12.4," at i,j = ",2i4," on pe ",i4)') &
380  u_comp(i,j), u_resym(i,j),u_comp(i,j)-u_resym(i,j),i,j,pe_here()
381  write(0,'(A130)') trim(mesg)//trim(mesg2)
383  endif
384  enddo ; enddo
385  do i=is_ch,ie_ch ; do j=js_ch,je_ch
386  if (v_resym(i,j) /= v_comp(i,j) .and. &
388  write(mesg2,'(" redundant v-comps",2(1pe12.4)," differ by ", &
389  & 1pe12.4," at i,j = ",2i4," x,y = ",2(1pe12.4)" on pe ",i4)') &
390  v_comp(i,j), v_resym(i,j),v_comp(i,j)-v_resym(i,j),i,j, &
391  g%geoLonBu(i,j), g%geoLatBu(i,j), pe_here()
392  write(0,'(A155)') trim(mesg)//trim(mesg2)
394  endif
395  enddo ; enddo
396 
397 end subroutine check_redundant_vb2d
398 
399 !> Check for consistency between the duplicated points of a 3-D scalar at tracer points
400 subroutine check_redundant_st3d(mesg, array, G, is, ie, js, je)
401  character(len=*), intent(in) :: mesg !< An identifying message
402  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
403  real, dimension(G%isd:,G%jsd:,:), intent(in) :: array !< The array to be checked for consistency
404  integer, optional, intent(in) :: is !< The starting i-index to check
405  integer, optional, intent(in) :: ie !< The ending i-index to check
406  integer, optional, intent(in) :: js !< The starting j-index to check
407  integer, optional, intent(in) :: je !< The ending j-index to check
408  ! Local variables
409  character(len=24) :: mesg_k
410  integer :: k
411 
412  do k=1,size(array,3)
413  if (k < 10) then ; write(mesg_k,'(" Layer",i2," ")') k
414  elseif (k < 100) then ; write(mesg_k,'(" Layer",i3," ")') k
415  elseif (k < 1000) then ; write(mesg_k,'(" Layer",i4," ")') k
416  else ; write(mesg_k,'(" Layer",i9," ")') k ; endif
417 
418  call check_redundant_st2d(trim(mesg)//trim(mesg_k), array(:,:,k), &
419  g, is, ie, js, je)
420  enddo
421 end subroutine check_redundant_st3d
422 
423 
424 !> Check for consistency between the duplicated points of a 2-D scalar at tracer points
425 subroutine check_redundant_st2d(mesg, array, G, is, ie, js, je)
426  character(len=*), intent(in) :: mesg !< An identifying message
427  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
428  real, dimension(G%isd:,G%jsd:), intent(in) :: array !< The array to be checked for consistency
429  integer, optional, intent(in) :: is !< The starting i-index to check
430  integer, optional, intent(in) :: ie !< The ending i-index to check
431  integer, optional, intent(in) :: js !< The starting j-index to check
432  integer, optional, intent(in) :: je !< The ending j-index to check
433  ! Local variables
434  real :: a_nonsym(G%isd:G%ied,G%jsd:G%jed)
435  character(len=128) :: mesg2
436 
437  integer :: i, j, is_ch, ie_ch, js_ch, je_ch
438  integer :: Isq, Ieq, Jsq, Jeq, isd, ied, jsd, jed
439  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
440 
441  is_ch = g%isc ; ie_ch = g%iec ; js_ch = g%jsc ; je_ch = g%jec
442  if (present(is)) is_ch = is ; if (present(ie)) ie_ch = ie
443  if (present(js)) js_ch = js ; if (present(js)) je_ch = je
444 
445  ! This only works on points outside of the standard computational domain.
446  if ((is_ch == g%isc) .and. (ie_ch == g%iec) .and. &
447  (js_ch == g%jsc) .and. (je_ch == g%jec)) return
448 
449  do i=isd,ied ; do j=jsd,jed
450  a_nonsym(i,j) = array(i,j)
451  enddo ; enddo
452 
453  call pass_var(a_nonsym, g%Domain)
454 
455  do i=is_ch,ie_ch ; do j=js_ch,je_ch
456  if (a_nonsym(i,j) /= array(i,j) .and. &
458  write(mesg2,'(" Redundant points",2(1pe12.4)," differ by ", &
459  & 1pe12.4," at i,j = ",2i4," on pe ",i4)') &
460  array(i,j), a_nonsym(i,j),array(i,j)-a_nonsym(i,j),i,j,pe_here()
461  write(0,'(A130)') trim(mesg)//trim(mesg2)
463  endif
464  enddo ; enddo
465 
466 end subroutine check_redundant_st2d
467 
468 !> Check for consistency between the duplicated points of a 3-D A-grid vector
469 subroutine check_redundant_vt3d(mesg, u_comp, v_comp, G, is, ie, js, je, &
470  direction)
471  character(len=*), intent(in) :: mesg !< An identifying message
472  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
473  real, dimension(G%isd:,G%jsd:,:), intent(in) :: u_comp !< The u-component of the vector
474  !! to be checked for consistency
475  real, dimension(G%isd:,G%jsd:,:), intent(in) :: v_comp !< The v-component of the vector
476  !! to be checked for consistency
477  integer, optional, intent(in) :: is !< The starting i-index to check
478  integer, optional, intent(in) :: ie !< The ending i-index to check
479  integer, optional, intent(in) :: js !< The starting j-index to check
480  integer, optional, intent(in) :: je !< The ending j-index to check
481  integer, optional, intent(in) :: direction !< the direction flag to be
482  !! passed to pass_vector
483  ! Local variables
484  character(len=24) :: mesg_k
485  integer :: k
486 
487  do k=1,size(u_comp,3)
488  if (k < 10) then ; write(mesg_k,'(" Layer",i2," ")') k
489  elseif (k < 100) then ; write(mesg_k,'(" Layer",i3," ")') k
490  elseif (k < 1000) then ; write(mesg_k,'(" Layer",i4," ")') k
491  else ; write(mesg_k,'(" Layer",i9," ")') k ; endif
492 
493  call check_redundant_vt2d(trim(mesg)//trim(mesg_k), u_comp(:,:,k), &
494  v_comp(:,:,k), g, is, ie, js, je, direction)
495  enddo
496 end subroutine check_redundant_vt3d
497 
498 !> Check for consistency between the duplicated points of a 2-D A-grid vector
499 subroutine check_redundant_vt2d(mesg, u_comp, v_comp, G, is, ie, js, je, &
500  direction)
501  character(len=*), intent(in) :: mesg !< An identifying message
502  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
503  real, dimension(G%isd:,G%jsd:), intent(in) :: u_comp !< The u-component of the vector
504  !! to be checked for consistency
505  real, dimension(G%isd:,G%jsd:), intent(in) :: v_comp !< The v-component of the vector
506  !! to be checked for consistency
507  integer, optional, intent(in) :: is !< The starting i-index to check
508  integer, optional, intent(in) :: ie !< The ending i-index to check
509  integer, optional, intent(in) :: js !< The starting j-index to check
510  integer, optional, intent(in) :: je !< The ending j-index to check
511  integer, optional, intent(in) :: direction !< the direction flag to be
512  !! passed to pass_vector
513  ! Local variables
514  real :: u_nonsym(G%isd:G%ied,G%jsd:G%jed)
515  real :: v_nonsym(G%isd:G%ied,G%jsd:G%jed)
516  character(len=128) :: mesg2
517 
518  integer :: i, j, is_ch, ie_ch, js_ch, je_ch
519  integer :: Isq, Ieq, Jsq, Jeq, isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
520  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
521  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
522  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
523 
524  is_ch = g%isc ; ie_ch = g%iec ; js_ch = g%jsc ; je_ch = g%jec
525  if (present(is)) is_ch = is ; if (present(ie)) ie_ch = ie
526  if (present(js)) js_ch = js ; if (present(js)) je_ch = je
527 
528  ! This only works on points outside of the standard computational domain.
529  if ((is_ch == g%isc) .and. (ie_ch == g%iec) .and. &
530  (js_ch == g%jsc) .and. (je_ch == g%jec)) return
531 
532  do i=isd,ied ; do j=jsd,jed
533  u_nonsym(i,j) = u_comp(i,j) ; v_nonsym(i,j) = v_comp(i,j)
534  enddo ; enddo
535 
536  call pass_vector(u_nonsym, v_nonsym, g%Domain, direction, stagger=agrid)
537 
538  do i=is_ch,ie_ch ; do j=js_ch+1,je_ch
539  if (u_nonsym(i,j) /= u_comp(i,j) .and. &
541  write(mesg2,'(" redundant u-components",2(1pe12.4)," differ by ", &
542  & 1pe12.4," at i,j = ",2i4," on pe ",i4)') &
543  u_comp(i,j), u_nonsym(i,j),u_comp(i,j)-u_nonsym(i,j),i,j,pe_here()
544  write(0,'(A130)') trim(mesg)//trim(mesg2)
546  endif
547  enddo ; enddo
548  do i=is_ch+1,ie_ch ; do j=js_ch,je_ch
549  if (v_nonsym(i,j) /= v_comp(i,j) .and. &
551  write(mesg2,'(" redundant v-comps",2(1pe12.4)," differ by ", &
552  & 1pe12.4," at i,j = ",2i4," x,y = ",2(1pe12.4)" on pe ",i4)') &
553  v_comp(i,j), v_nonsym(i,j),v_comp(i,j)-v_nonsym(i,j),i,j, &
554  g%geoLonBu(i,j), g%geoLatBu(i,j), pe_here()
555  write(0,'(A155)') trim(mesg)//trim(mesg2)
557  endif
558  enddo ; enddo
559 
560 end subroutine check_redundant_vt2d
561 
562 !> Do a checksum and redundant point check on a 3d C-grid vector.
563 subroutine chksum_vec_c3d(mesg, u_comp, v_comp, G, halos, scalars)
564  character(len=*), intent(in) :: mesg !< An identifying message
565  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
566  real, dimension(G%IsdB:,G%jsd:,:), intent(in) :: u_comp !< The u-component of the vector
567  real, dimension(G%isd:,G%JsdB:,:), intent(in) :: v_comp !< The v-component of the vector
568  integer, optional, intent(in) :: halos !< The width of halos to check (default 0)
569  logical, optional, intent(in) :: scalars !< If true this is a pair of
570  !! scalars that are being checked.
571  ! Local variables
572  logical :: are_scalars
573  are_scalars = .false. ; if (present(scalars)) are_scalars = scalars
574 
575  if (debug_chksums) then
576  call uvchksum(mesg, u_comp, v_comp, g%HI, halos)
577  endif
578  if (debug_redundant) then
579  if (are_scalars) then
580  call check_redundant_c(mesg, u_comp, v_comp, g, direction=to_all+scalar_pair)
581  else
582  call check_redundant_c(mesg, u_comp, v_comp, g)
583  endif
584  endif
585 
586 end subroutine chksum_vec_c3d
587 
588 !> Do a checksum and redundant point check on a 2d C-grid vector.
589 subroutine chksum_vec_c2d(mesg, u_comp, v_comp, G, halos, scalars)
590  character(len=*), intent(in) :: mesg !< An identifying message
591  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
592  real, dimension(G%IsdB:,G%jsd:), intent(in) :: u_comp !< The u-component of the vector
593  real, dimension(G%isd:,G%JsdB:), intent(in) :: v_comp !< The v-component of the vector
594  integer, optional, intent(in) :: halos !< The width of halos to check (default 0)
595  logical, optional, intent(in) :: scalars !< If true this is a pair of
596  !! scalars that are being checked.
597  ! Local variables
598  logical :: are_scalars
599  are_scalars = .false. ; if (present(scalars)) are_scalars = scalars
600 
601  if (debug_chksums) then
602  call uvchksum(mesg, u_comp, v_comp, g%HI, halos)
603  endif
604  if (debug_redundant) then
605  if (are_scalars) then
606  call check_redundant_c(mesg, u_comp, v_comp, g, direction=to_all+scalar_pair)
607  else
608  call check_redundant_c(mesg, u_comp, v_comp, g)
609  endif
610  endif
611 
612 end subroutine chksum_vec_c2d
613 
614 !> Do a checksum and redundant point check on a 3d B-grid vector.
615 subroutine chksum_vec_b3d(mesg, u_comp, v_comp, G, halos, scalars)
616  character(len=*), intent(in) :: mesg !< An identifying message
617  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
618  real, dimension(G%IsdB:,G%JsdB:,:), intent(in) :: u_comp !< The u-component of the vector
619  real, dimension(G%IsdB:,G%JsdB:,:), intent(in) :: v_comp !< The v-component of the vector
620  integer, optional, intent(in) :: halos !< The width of halos to check (default 0)
621  logical, optional, intent(in) :: scalars !< If true this is a pair of
622  !! scalars that are being checked.
623  ! Local variables
624  logical :: are_scalars
625  are_scalars = .false. ; if (present(scalars)) are_scalars = scalars
626 
627  if (debug_chksums) then
628  call bchksum(u_comp, mesg//"(u)", g%HI, halos)
629  call bchksum(v_comp, mesg//"(v)", g%HI, halos)
630  endif
631  if (debug_redundant) then
632  if (are_scalars) then
633  call check_redundant_b(mesg, u_comp, v_comp, g, direction=to_all+scalar_pair)
634  else
635  call check_redundant_b(mesg, u_comp, v_comp, g)
636  endif
637  endif
638 
639 end subroutine chksum_vec_b3d
640 
641 ! Do a checksum and redundant point check on a 2d B-grid vector.
642 subroutine chksum_vec_b2d(mesg, u_comp, v_comp, G, halos, scalars, symmetric)
643  character(len=*), intent(in) :: mesg !< An identifying message
644  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
645  real, dimension(G%IsdB:,G%JsdB:), intent(in) :: u_comp !< The u-component of the vector
646  real, dimension(G%IsdB:,G%JsdB:), intent(in) :: v_comp !< The v-component of the vector
647  integer, optional, intent(in) :: halos !< The width of halos to check (default 0)
648  logical, optional, intent(in) :: scalars !< If true this is a pair of
649  !! scalars that are being checked.
650  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the
651  !! full symmetric computational domain.
652  ! Local variables
653  logical :: are_scalars
654  are_scalars = .false. ; if (present(scalars)) are_scalars = scalars
655 
656  if (debug_chksums) then
657  call bchksum(u_comp, mesg//"(u)", g%HI, halos, symmetric=symmetric)
658  call bchksum(v_comp, mesg//"(v)", g%HI, halos, symmetric=symmetric)
659  endif
660  if (debug_redundant) then
661  if (are_scalars) then
662  call check_redundant_b(mesg, u_comp, v_comp, g, direction=to_all+scalar_pair)
663  else
664  call check_redundant_b(mesg, u_comp, v_comp, g)
665  endif
666  endif
667 
668 end subroutine chksum_vec_b2d
669 
670 !> Do a checksum and redundant point check on a 3d C-grid vector.
671 subroutine chksum_vec_a3d(mesg, u_comp, v_comp, G, halos, scalars)
672  character(len=*), intent(in) :: mesg !< An identifying message
673  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
674  real, dimension(G%isd:,G%jsd:,:), intent(in) :: u_comp !< The u-component of the vector
675  real, dimension(G%isd:,G%jsd:,:), intent(in) :: v_comp !< The v-component of the vector
676  integer, optional, intent(in) :: halos !< The width of halos to check (default 0)
677  logical, optional, intent(in) :: scalars !< If true this is a pair of
678  !! scalars that are being checked.
679  ! Local variables
680  logical :: are_scalars
681  are_scalars = .false. ; if (present(scalars)) are_scalars = scalars
682 
683  if (debug_chksums) then
684  call hchksum(u_comp, mesg//"(u)", g%HI, halos)
685  call hchksum(v_comp, mesg//"(v)", g%HI, halos)
686  endif
687  if (debug_redundant) then
688  if (are_scalars) then
689  call check_redundant_t(mesg, u_comp, v_comp, g, direction=to_all+scalar_pair)
690  else
691  call check_redundant_t(mesg, u_comp, v_comp, g)
692  endif
693  endif
694 
695 end subroutine chksum_vec_a3d
696 
697 !> Do a checksum and redundant point check on a 2d C-grid vector.
698 subroutine chksum_vec_a2d(mesg, u_comp, v_comp, G, halos, scalars)
699  character(len=*), intent(in) :: mesg !< An identifying message
700  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
701  real, dimension(G%isd:,G%jsd:), intent(in) :: u_comp !< The u-component of the vector
702  real, dimension(G%isd:,G%jsd:), intent(in) :: v_comp !< The v-component of the vector
703  integer, optional, intent(in) :: halos !< The width of halos to check (default 0)
704  logical, optional, intent(in) :: scalars !< If true this is a pair of
705  !! scalars that are being checked.
706  ! Local variables
707  logical :: are_scalars
708  are_scalars = .false. ; if (present(scalars)) are_scalars = scalars
709 
710  if (debug_chksums) then
711  call hchksum(u_comp, mesg//"(u)", g%HI, halos)
712  call hchksum(v_comp, mesg//"(v)", g%HI, halos)
713  endif
714  if (debug_redundant) then
715  if (are_scalars) then
716  call check_redundant_t(mesg, u_comp, v_comp, g, direction=to_all+scalar_pair)
717  else
718  call check_redundant_t(mesg, u_comp, v_comp, g)
719  endif
720  endif
721 
722 end subroutine chksum_vec_a2d
723 
724 !> This function returns the sum over computational domain of all
725 !! processors of hThick*stuff, where stuff is a 3-d array at tracer points.
726 function totalstuff(HI, hThick, areaT, stuff)
727  type(hor_index_type), intent(in) :: hi !< A horizontal index type
728  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: hthick !< The array of thicknesses to use as weights
729  real, dimension(HI%isd:,HI%jsd:), intent(in) :: areat !< The array of cell areas [m2]
730  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: stuff !< The array of stuff to be summed
731  real :: totalstuff !< the globally integrated amoutn of stuff
732  ! Local variables
733  integer :: i, j, k, nz
734 
735  nz = size(hthick,3)
736  totalstuff = 0.
737  do k = 1, nz ; do j = hi%jsc, hi%jec ; do i = hi%isc, hi%iec
738  totalstuff = totalstuff + hthick(i,j,k) * stuff(i,j,k) * areat(i,j)
739  enddo ; enddo ; enddo
740  call sum_across_pes(totalstuff)
741 
742 end function totalstuff
743 
744 !> This subroutine display the total thickness, temperature and salinity
745 !! as well as the change since the last call.
746 subroutine totaltands(HI, hThick, areaT, temperature, salinity, mesg)
747  type(hor_index_type), intent(in) :: hi !< A horizontal index type
748  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: hthick !< The array of thicknesses to use as weights
749  real, dimension(HI%isd:,HI%jsd:), intent(in) :: areat !< The array of cell areas [m2]
750  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: temperature !< The temperature field to sum
751  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: salinity !< The salinity field to sum
752  character(len=*), intent(in) :: mesg !< An identifying message
753  ! NOTE: This subroutine uses "save" data which is not thread safe and is purely for
754  ! extreme debugging without a proper debugger.
755  real, save :: totalh = 0., totalt = 0., totals = 0.
756  ! Local variables
757  logical, save :: firstcall = .true.
758  real :: thish, thist, thiss, delh, delt, dels
759  integer :: i, j, k, nz
760 
761  nz = size(hthick,3)
762  thish = 0.
763  do k = 1, nz ; do j = hi%jsc, hi%jec ; do i = hi%isc, hi%iec
764  thish = thish + hthick(i,j,k) * areat(i,j)
765  enddo ; enddo ; enddo
766  call sum_across_pes(thish)
767  thist = totalstuff(hi, hthick, areat, temperature)
768  thiss = totalstuff(hi, hthick, areat, salinity)
769 
770  if (is_root_pe()) then
771  if (firstcall) then
772  totalh = thish ; totalt = thist ; totals = thiss
773  write(0,*) 'Totals H,T,S:',thish,thist,thiss,' ',mesg
774  firstcall = .false.
775  else
776  delh = thish - totalh
777  delt = thist - totalt
778  dels = thiss - totals
779  totalh = thish ; totalt = thist ; totals = thiss
780  write(0,*) 'Tot/del H,T,S:',thish,thist,thiss,delh,delt,dels,' ',mesg
781  endif
782  endif
783 
784 end subroutine totaltands
785 
786 !> Returns false if the column integral of a given quantity is within roundoff
787 logical function check_column_integral(nk, field, known_answer)
788  integer, intent(in) :: nk !< Number of levels in column
789  real, dimension(nk), intent(in) :: field !< Field to be summed
790  real, optional, intent(in) :: known_answer !< If present is the expected sum,
791  !! If missing, assumed zero
792  ! Local variables
793  real :: u_sum, error, expected
794  integer :: k
795 
796  u_sum = field(1)
797  error = 0.
798 
799  ! Reintegrate and sum roundoff errors
800  do k=2,nk
801  u_sum = u_sum + field(k)
802  error = error + epsilon(u_sum)*max(abs(u_sum),abs(field(k)))
803  enddo
804 
805  ! Assign expected answer to either the optional input or 0
806  if (present(known_answer)) then
807  expected = known_answer
808  else
809  expected = 0.
810  endif
811 
812  ! Compare the column integrals against calculated roundoff error
813  if (abs(u_sum-expected) > error) then
814  check_column_integral = .true.
815  else
816  check_column_integral = .false.
817  endif
818 
819 end function check_column_integral
820 
821 !> Returns false if the column integrals of two given quantities are within roundoff of each other
822 logical function check_column_integrals(nk_1, field_1, nk_2, field_2, missing_value)
823  integer, intent(in) :: nk_1 !< Number of levels in field 1
824  integer, intent(in) :: nk_2 !< Number of levels in field 2
825  real, dimension(nk_1), intent(in) :: field_1 !< First field to be summed
826  real, dimension(nk_2), intent(in) :: field_2 !< Second field to be summed
827  real, optional, intent(in) :: missing_value !< If column contains missing values,
828  !! mask them from the sum
829  ! Local variables
830  real :: u1_sum, error1, u2_sum, error2, misval
831  integer :: k
832 
833  ! Assign missing value
834  if (present(missing_value)) then
835  misval = missing_value
836  else
837  misval = 0.
838  endif
839 
840  u1_sum = field_1(1)
841  error1 = 0.
842 
843  ! Reintegrate and sum roundoff errors
844  do k=2,nk_1
845  if (field_1(k)/=misval) then
846  u1_sum = u1_sum + field_1(k)
847  error1 = error1 + epsilon(u1_sum)*max(abs(u1_sum),abs(field_1(k)))
848  endif
849  enddo
850 
851  u2_sum = field_2(1)
852  error2 = 0.
853 
854  ! Reintegrate and sum roundoff errors
855  do k=2,nk_2
856  if (field_2(k)/=misval) then
857  u2_sum = u2_sum + field_2(k)
858  error2 = error2 + epsilon(u2_sum)*max(abs(u2_sum),abs(field_2(k)))
859  endif
860  enddo
861 
862  ! Compare the column integrals against calculated roundoff error
863  if (abs(u1_sum-u2_sum) > (error1+error2)) then
864  check_column_integrals = .true.
865  else
866  check_column_integrals = .false.
867  endif
868 
869 end function check_column_integrals
870 
871 end module mom_debugging
mom_debugging::check_redundant_vb2d
subroutine check_redundant_vb2d(mesg, u_comp, v_comp, G, is, ie, js, je, direction)
Check for consistency between the duplicated points of a 2-D B-grid vector.
Definition: MOM_debugging.F90:325
mom_debugging::vec_chksum_b
Do checksums on the components of a B-grid vector.
Definition: MOM_debugging.F90:60
mom_debugging::check_redundant_vt3d
subroutine check_redundant_vt3d(mesg, u_comp, v_comp, G, is, ie, js, je, direction)
Check for consistency between the duplicated points of a 3-D A-grid vector.
Definition: MOM_debugging.F90:471
mom_debugging::check_redundant_t
Check for consistency between the duplicated points of an A-grid vector or scalar.
Definition: MOM_debugging.F90:46
mom_debugging::vec_chksum
Do checksums on the components of a C-grid vector.
Definition: MOM_debugging.F90:52
mom_debugging::check_column_integral
logical function, public check_column_integral(nk, field, known_answer)
Returns false if the column integral of a given quantity is within roundoff.
Definition: MOM_debugging.F90:788
mom_debugging::check_redundant_vb3d
subroutine check_redundant_vb3d(mesg, u_comp, v_comp, G, is, ie, js, je, direction)
Check for consistency between the duplicated points of a 3-D B-grid vector.
Definition: MOM_debugging.F90:295
mom_checksums::bchksum
Checksums an array (2d or 3d) staggered at corner points.
Definition: MOM_checksums.F90:53
mom_checksums::mom_checksums_init
subroutine, public mom_checksums_init(param_file)
MOM_checksums_init initializes the MOM_checksums module. As it happens, the only thing that it does i...
Definition: MOM_checksums.F90:1835
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_checksums::is_nan
Returns .true. if any element of x is a NaN, and .false. otherwise.
Definition: MOM_checksums.F90:73
mom_debugging::max_redundant_prints
integer max_redundant_prints
Maximum number of times to write redundant messages.
Definition: MOM_debugging.F90:70
mom_checksums::qchksum
This is an older interface that has been renamed Bchksum.
Definition: MOM_checksums.F90:58
mom_debugging::check_redundant_b
Check for consistency between the duplicated points of a B-grid vector or scalar.
Definition: MOM_debugging.F90:41
mom_checksums::uvchksum
Checksums a pair velocity arrays (2d or 3d) staggered at C-grid locations.
Definition: MOM_checksums.F90:28
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
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_checksums::chksum
This is an older interface for 1-, 2-, or 3-D checksums.
Definition: MOM_checksums.F90:63
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_hor_index
Defines the horizontal index type (hor_index_type) used for providing index ranges.
Definition: MOM_hor_index.F90:2
mom_debugging::check_column_integrals
logical function, public check_column_integrals(nk_1, field_1, nk_2, field_2, missing_value)
Returns false if the column integrals of two given quantities are within roundoff of each other.
Definition: MOM_debugging.F90:823
mom_coms
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Definition: MOM_coms.F90:3
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_debugging::mom_debugging_init
subroutine, public mom_debugging_init(param_file)
MOM_debugging_init initializes the MOM_debugging module, and sets the parameterts that control which ...
Definition: MOM_debugging.F90:81
mom_domains::to_all
integer, parameter, public to_all
A flag for passing in all directions.
Definition: MOM_domains.F90:132
mom_debugging::totaltands
subroutine, public totaltands(HI, hThick, areaT, temperature, salinity, mesg)
This subroutine display the total thickness, temperature and salinity as well as the change since the...
Definition: MOM_debugging.F90:747
mom_debugging::check_redundant_vc2d
subroutine check_redundant_vc2d(mesg, u_comp, v_comp, G, is, ie, js, je, direction)
Check for consistency between the duplicated points of a 2-D C-grid vector.
Definition: MOM_debugging.F90:136
mom_debugging::vec_chksum_c
Do checksums on the components of a C-grid vector.
Definition: MOM_debugging.F90:56
mom_checksums
Routines to calculate checksums of various array and vector types.
Definition: MOM_checksums.F90:2
mom_debugging::debug
logical debug
Write out verbose debugging data.
Definition: MOM_debugging.F90:72
mom_debugging::check_redundant_st3d
subroutine check_redundant_st3d(mesg, array, G, is, ie, js, je)
Check for consistency between the duplicated points of a 3-D scalar at tracer points.
Definition: MOM_debugging.F90:401
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_debugging::check_redundant_st2d
subroutine check_redundant_st2d(mesg, array, G, is, ie, js, je)
Check for consistency between the duplicated points of a 2-D scalar at tracer points.
Definition: MOM_debugging.F90:426
mom_debugging::debug_chksums
logical debug_chksums
Perform checksums on arrays.
Definition: MOM_debugging.F90:73
mom_debugging::debug_redundant
logical debug_redundant
Check redundant values on PE boundaries.
Definition: MOM_debugging.F90:74
mom_debugging::chksum_vec_a2d
subroutine chksum_vec_a2d(mesg, u_comp, v_comp, G, halos, scalars)
Do a checksum and redundant point check on a 2d C-grid vector.
Definition: MOM_debugging.F90:699
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_debugging::check_redundant_sb2d
subroutine check_redundant_sb2d(mesg, array, G, is, ie, js, je)
Check for consistency between the duplicated points of a 2-D scalar at corner points.
Definition: MOM_debugging.F90:236
mom_debugging::check_redundant_c
Check for consistency between the duplicated points of a C-grid vector.
Definition: MOM_debugging.F90:37
mom_debugging::chksum_vec_b2d
subroutine chksum_vec_b2d(mesg, u_comp, v_comp, G, halos, scalars, symmetric)
Definition: MOM_debugging.F90:643
mom_debugging::redundant_prints
integer, dimension(3) redundant_prints
Counters for controlling redundant printing.
Definition: MOM_debugging.F90:71
mom_debugging::chksum_vec_a3d
subroutine chksum_vec_a3d(mesg, u_comp, v_comp, G, halos, scalars)
Do a checksum and redundant point check on a 3d C-grid vector.
Definition: MOM_debugging.F90:672
mom_hor_index::hor_index_type
Container for horizontal index ranges for data, computational and global domains.
Definition: MOM_hor_index.F90:15
mom_debugging::chksum_vec_b3d
subroutine chksum_vec_b3d(mesg, u_comp, v_comp, G, halos, scalars)
Do a checksum and redundant point check on a 3d B-grid vector.
Definition: MOM_debugging.F90:616
mom_debugging::check_redundant_vc3d
subroutine check_redundant_vc3d(mesg, u_comp, v_comp, G, is, ie, js, je, direction)
Check for consistency between the duplicated points of a 3-D C-grid vector.
Definition: MOM_debugging.F90:106
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_debugging::check_redundant
Check for consistency between the duplicated points of a C-grid vector.
Definition: MOM_debugging.F90:33
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
mom_debugging::chksum_vec_c2d
subroutine chksum_vec_c2d(mesg, u_comp, v_comp, G, halos, scalars)
Do a checksum and redundant point check on a 2d C-grid vector.
Definition: MOM_debugging.F90:590
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_debugging::vec_chksum_a
Do checksums on the components of an A-grid vector.
Definition: MOM_debugging.F90:64
mom_debugging::totalstuff
real function, public totalstuff(HI, hThick, areaT, stuff)
This function returns the sum over computational domain of all processors of hThick*stuff,...
Definition: MOM_debugging.F90:727
mom_coms::reproducing_sum
Find an accurate and order-invariant sum of distributed 2d or 3d fields.
Definition: MOM_coms.F90:54
mom_checksums::hchksum
Checksums an array (2d or 3d) staggered at tracer points.
Definition: MOM_checksums.F90:48
mom_debugging::chksum_vec_c3d
subroutine chksum_vec_c3d(mesg, u_comp, v_comp, G, halos, scalars)
Do a checksum and redundant point check on a 3d C-grid vector.
Definition: MOM_debugging.F90:564
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_debugging::check_redundant_sb3d
subroutine check_redundant_sb3d(mesg, array, G, is, ie, js, je)
Check for consistency between the duplicated points of a 3-D scalar at corner points.
Definition: MOM_debugging.F90:211
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_debugging::check_redundant_vt2d
subroutine check_redundant_vt2d(mesg, u_comp, v_comp, G, is, ie, js, je, direction)
Check for consistency between the duplicated points of a 2-D A-grid vector.
Definition: MOM_debugging.F90:501