MOM6
mom_debugging Module Reference

Detailed Description

Provides checksumming functions for debugging.

This module contains subroutines that perform various error checking and debugging functions for MOM6. This routine is similar to it counterpart in the SIS2 code, except for the use of the ocean_grid_type and by keeping them separate we retain the ability to set up MOM6 and SIS2 debugging separately.

Data Types

interface  check_redundant
 Check for consistency between the duplicated points of a C-grid vector. More...
 
interface  check_redundant_b
 Check for consistency between the duplicated points of a B-grid vector or scalar. More...
 
interface  check_redundant_c
 Check for consistency between the duplicated points of a C-grid vector. More...
 
interface  check_redundant_t
 Check for consistency between the duplicated points of an A-grid vector or scalar. More...
 
interface  vec_chksum
 Do checksums on the components of a C-grid vector. More...
 
interface  vec_chksum_a
 Do checksums on the components of an A-grid vector. More...
 
interface  vec_chksum_b
 Do checksums on the components of a B-grid vector. More...
 
interface  vec_chksum_c
 Do checksums on the components of a C-grid vector. More...
 

Functions/Subroutines

subroutine, public mom_debugging_init (param_file)
 MOM_debugging_init initializes the MOM_debugging module, and sets the parameterts that control which checks are active for MOM6. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
subroutine chksum_vec_b2d (mesg, u_comp, v_comp, G, halos, scalars, symmetric)
 
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. More...
 
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. More...
 
real function, public totalstuff (HI, hThick, areaT, stuff)
 This function returns the sum over computational domain of all processors of hThick*stuff, where stuff is a 3-d array at tracer points. More...
 
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 last call. More...
 
logical function, public check_column_integral (nk, field, known_answer)
 Returns false if the column integral of a given quantity is within roundoff. More...
 
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. More...
 

Variables

integer max_redundant_prints = 100
 Maximum number of times to write redundant messages. More...
 
integer, dimension(3) redundant_prints = 0
 Counters for controlling redundant printing. More...
 
logical debug = .false.
 Write out verbose debugging data. More...
 
logical debug_chksums = .true.
 Perform checksums on arrays. More...
 
logical debug_redundant = .true.
 Check redundant values on PE boundaries. More...
 

Function/Subroutine Documentation

◆ check_column_integral()

logical function, public mom_debugging::check_column_integral ( integer, intent(in)  nk,
real, dimension(nk), intent(in)  field,
real, intent(in), optional  known_answer 
)

Returns false if the column integral of a given quantity is within roundoff.

Parameters
[in]nkNumber of levels in column
[in]fieldField to be summed
[in]known_answerIf present is the expected sum, If missing, assumed zero

Definition at line 788 of file MOM_debugging.F90.

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 

◆ check_column_integrals()

logical function, public mom_debugging::check_column_integrals ( integer, intent(in)  nk_1,
real, dimension(nk_1), intent(in)  field_1,
integer, intent(in)  nk_2,
real, dimension(nk_2), intent(in)  field_2,
real, intent(in), optional  missing_value 
)

Returns false if the column integrals of two given quantities are within roundoff of each other.

Parameters
[in]nk_1Number of levels in field 1
[in]nk_2Number of levels in field 2
[in]field_1First field to be summed
[in]field_2Second field to be summed
[in]missing_valueIf column contains missing values, mask them from the sum

Definition at line 823 of file MOM_debugging.F90.

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 

Referenced by mom_ale::ale_offline_inputs().

Here is the caller graph for this function:

◆ check_redundant_sb2d()

subroutine mom_debugging::check_redundant_sb2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsdb:), intent(in)  array,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je 
)
private

Check for consistency between the duplicated points of a 2-D scalar at corner points.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]arrayThe array to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check

Definition at line 236 of file MOM_debugging.F90.

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. &
281  redundant_prints(2) < max_redundant_prints) then
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)
286  redundant_prints(2) = redundant_prints(2) + 1
287  endif
288  enddo ; enddo
289 

References max_redundant_prints, mom_error_handler::mom_error(), redundant_prints, and mom_domains::to_all.

Referenced by check_redundant_sb3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ check_redundant_sb3d()

subroutine mom_debugging::check_redundant_sb3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsdb:,:), intent(in)  array,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je 
)
private

Check for consistency between the duplicated points of a 3-D scalar at corner points.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]arrayThe array to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check

Definition at line 211 of file MOM_debugging.F90.

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

References check_redundant_sb2d().

Here is the call graph for this function:

◆ check_redundant_st2d()

subroutine mom_debugging::check_redundant_st2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isd:,g%jsd:), intent(in)  array,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je 
)
private

Check for consistency between the duplicated points of a 2-D scalar at tracer points.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]arrayThe array to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check

Definition at line 426 of file MOM_debugging.F90.

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. &
457  redundant_prints(1) < max_redundant_prints) then
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)
462  redundant_prints(1) = redundant_prints(1) + 1
463  endif
464  enddo ; enddo
465 

References max_redundant_prints, and redundant_prints.

Referenced by check_redundant_st3d().

Here is the caller graph for this function:

◆ check_redundant_st3d()

subroutine mom_debugging::check_redundant_st3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isd:,g%jsd:,:), intent(in)  array,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je 
)
private

Check for consistency between the duplicated points of a 3-D scalar at tracer points.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]arrayThe array to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check

Definition at line 401 of file MOM_debugging.F90.

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

References check_redundant_st2d().

Here is the call graph for this function:

◆ check_redundant_vb2d()

subroutine mom_debugging::check_redundant_vb2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsdb:), intent(in)  u_comp,
real, dimension(g%isdb:,g%jsdb:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je,
integer, intent(in), optional  direction 
)
private

Check for consistency between the duplicated points of a 2-D B-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector to be checked for consistency
[in]v_compThe v-component of the vector to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check
[in]directionthe direction flag to be passed to pass_vector

Definition at line 325 of file MOM_debugging.F90.

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. &
377  redundant_prints(2) < max_redundant_prints) then
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)
382  redundant_prints(2) = redundant_prints(2) + 1
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. &
387  redundant_prints(2) < max_redundant_prints) then
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)
393  redundant_prints(2) = redundant_prints(2) + 1
394  endif
395  enddo ; enddo
396 

References max_redundant_prints, mom_error_handler::mom_error(), and redundant_prints.

Referenced by check_redundant_vb3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ check_redundant_vb3d()

subroutine mom_debugging::check_redundant_vb3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsdb:,:), intent(in)  u_comp,
real, dimension(g%isdb:,g%jsdb:,:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je,
integer, intent(in), optional  direction 
)
private

Check for consistency between the duplicated points of a 3-D B-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector to be checked for consistency
[in]v_compThe v-component of the vector to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check
[in]directionthe direction flag to be passed to pass_vector

Definition at line 295 of file MOM_debugging.F90.

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

References check_redundant_vb2d().

Here is the call graph for this function:

◆ check_redundant_vc2d()

subroutine mom_debugging::check_redundant_vc2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsd:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsdb:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je,
integer, intent(in), optional  direction 
)
private

Check for consistency between the duplicated points of a 2-D C-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector to be checked for consistency
[in]v_compThe u-component of the vector to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check
[in]directionthe direction flag to be passed to pass_vector

Definition at line 136 of file MOM_debugging.F90.

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. &
187  redundant_prints(3) < max_redundant_prints) then
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)
192  redundant_prints(3) = redundant_prints(3) + 1
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. &
197  redundant_prints(3) < max_redundant_prints) then
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)
203  redundant_prints(3) = redundant_prints(3) + 1
204  endif
205  enddo ; enddo
206 

References max_redundant_prints, mom_error_handler::mom_error(), and redundant_prints.

Referenced by check_redundant_vc3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ check_redundant_vc3d()

subroutine mom_debugging::check_redundant_vc3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsd:,:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsdb:,:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je,
integer, intent(in), optional  direction 
)
private

Check for consistency between the duplicated points of a 3-D C-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector to be checked for consistency
[in]v_compThe u-component of the vector to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check
[in]directionthe direction flag to be passed to pass_vector

Definition at line 106 of file MOM_debugging.F90.

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

References check_redundant_vc2d().

Here is the call graph for this function:

◆ check_redundant_vt2d()

subroutine mom_debugging::check_redundant_vt2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isd:,g%jsd:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsd:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je,
integer, intent(in), optional  direction 
)
private

Check for consistency between the duplicated points of a 2-D A-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector to be checked for consistency
[in]v_compThe v-component of the vector to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check
[in]directionthe direction flag to be passed to pass_vector

Definition at line 501 of file MOM_debugging.F90.

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. &
540  redundant_prints(1) < max_redundant_prints) then
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)
545  redundant_prints(1) = redundant_prints(1) + 1
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. &
550  redundant_prints(1) < max_redundant_prints) then
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)
556  redundant_prints(1) = redundant_prints(1) + 1
557  endif
558  enddo ; enddo
559 

References max_redundant_prints, and redundant_prints.

Referenced by check_redundant_vt3d().

Here is the caller graph for this function:

◆ check_redundant_vt3d()

subroutine mom_debugging::check_redundant_vt3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isd:,g%jsd:,:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsd:,:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  is,
integer, intent(in), optional  ie,
integer, intent(in), optional  js,
integer, intent(in), optional  je,
integer, intent(in), optional  direction 
)
private

Check for consistency between the duplicated points of a 3-D A-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector to be checked for consistency
[in]v_compThe v-component of the vector to be checked for consistency
[in]isThe starting i-index to check
[in]ieThe ending i-index to check
[in]jsThe starting j-index to check
[in]jeThe ending j-index to check
[in]directionthe direction flag to be passed to pass_vector

Definition at line 471 of file MOM_debugging.F90.

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

References check_redundant_vt2d().

Here is the call graph for this function:

◆ chksum_vec_a2d()

subroutine mom_debugging::chksum_vec_a2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isd:,g%jsd:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsd:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  halos,
logical, intent(in), optional  scalars 
)
private

Do a checksum and redundant point check on a 2d C-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector
[in]v_compThe v-component of the vector
[in]halosThe width of halos to check (default 0)
[in]scalarsIf true this is a pair of scalars that are being checked.

Definition at line 699 of file MOM_debugging.F90.

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 

References debug_chksums, debug_redundant, and mom_domains::to_all.

◆ chksum_vec_a3d()

subroutine mom_debugging::chksum_vec_a3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isd:,g%jsd:,:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsd:,:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  halos,
logical, intent(in), optional  scalars 
)
private

Do a checksum and redundant point check on a 3d C-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector
[in]v_compThe v-component of the vector
[in]halosThe width of halos to check (default 0)
[in]scalarsIf true this is a pair of scalars that are being checked.

Definition at line 672 of file MOM_debugging.F90.

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 

References debug_chksums, debug_redundant, and mom_domains::to_all.

◆ chksum_vec_b2d()

subroutine mom_debugging::chksum_vec_b2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsdb:), intent(in)  u_comp,
real, dimension(g%isdb:,g%jsdb:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  halos,
logical, intent(in), optional  scalars,
logical, intent(in), optional  symmetric 
)
private
Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector
[in]v_compThe v-component of the vector
[in]halosThe width of halos to check (default 0)
[in]scalarsIf true this is a pair of scalars that are being checked.
[in]symmetricIf true, do the checksums on the full symmetric computational domain.

Definition at line 643 of file MOM_debugging.F90.

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 

References debug_chksums, debug_redundant, and mom_domains::to_all.

◆ chksum_vec_b3d()

subroutine mom_debugging::chksum_vec_b3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsdb:,:), intent(in)  u_comp,
real, dimension(g%isdb:,g%jsdb:,:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  halos,
logical, intent(in), optional  scalars 
)
private

Do a checksum and redundant point check on a 3d B-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector
[in]v_compThe v-component of the vector
[in]halosThe width of halos to check (default 0)
[in]scalarsIf true this is a pair of scalars that are being checked.

Definition at line 616 of file MOM_debugging.F90.

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 

References debug_chksums, debug_redundant, and mom_domains::to_all.

◆ chksum_vec_c2d()

subroutine mom_debugging::chksum_vec_c2d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsd:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsdb:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  halos,
logical, intent(in), optional  scalars 
)
private

Do a checksum and redundant point check on a 2d C-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector
[in]v_compThe v-component of the vector
[in]halosThe width of halos to check (default 0)
[in]scalarsIf true this is a pair of scalars that are being checked.

Definition at line 590 of file MOM_debugging.F90.

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 

References debug_chksums, debug_redundant, and mom_domains::to_all.

◆ chksum_vec_c3d()

subroutine mom_debugging::chksum_vec_c3d ( character(len=*), intent(in)  mesg,
real, dimension(g%isdb:,g%jsd:,:), intent(in)  u_comp,
real, dimension(g%isd:,g%jsdb:,:), intent(in)  v_comp,
type(ocean_grid_type), intent(inout)  G,
integer, intent(in), optional  halos,
logical, intent(in), optional  scalars 
)
private

Do a checksum and redundant point check on a 3d C-grid vector.

Parameters
[in]mesgAn identifying message
[in,out]gThe ocean's grid structure
[in]u_compThe u-component of the vector
[in]v_compThe v-component of the vector
[in]halosThe width of halos to check (default 0)
[in]scalarsIf true this is a pair of scalars that are being checked.

Definition at line 564 of file MOM_debugging.F90.

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 

References debug_chksums, debug_redundant, and mom_domains::to_all.

◆ mom_debugging_init()

subroutine, public mom_debugging::mom_debugging_init ( type(param_file_type), intent(in)  param_file)

MOM_debugging_init initializes the MOM_debugging module, and sets the parameterts that control which checks are active for MOM6.

Parameters
[in]param_fileA structure to parse for run-time parameters

Definition at line 81 of file MOM_debugging.F90.

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 

References debug, debug_chksums, debug_redundant, and mom_checksums::mom_checksums_init().

Here is the call graph for this function:

◆ totalstuff()

real function, public mom_debugging::totalstuff ( type(hor_index_type), intent(in)  HI,
real, dimension(hi%isd:,hi%jsd:,:), intent(in)  hThick,
real, dimension(hi%isd:,hi%jsd:), intent(in)  areaT,
real, dimension(hi%isd:,hi%jsd:,:), intent(in)  stuff 
)

This function returns the sum over computational domain of all processors of hThick*stuff, where stuff is a 3-d array at tracer points.

Parameters
[in]hiA horizontal index type
[in]hthickThe array of thicknesses to use as weights
[in]areatThe array of cell areas [m2]
[in]stuffThe array of stuff to be summed
Returns
the globally integrated amoutn of stuff

Definition at line 727 of file MOM_debugging.F90.

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 

Referenced by totaltands().

Here is the caller graph for this function:

◆ totaltands()

subroutine, public mom_debugging::totaltands ( type(hor_index_type), intent(in)  HI,
real, dimension(hi%isd:,hi%jsd:,:), intent(in)  hThick,
real, dimension(hi%isd:,hi%jsd:), intent(in)  areaT,
real, dimension(hi%isd:,hi%jsd:,:), intent(in)  temperature,
real, dimension(hi%isd:,hi%jsd:,:), intent(in)  salinity,
character(len=*), intent(in)  mesg 
)

This subroutine display the total thickness, temperature and salinity as well as the change since the last call.

Parameters
[in]hiA horizontal index type
[in]hthickThe array of thicknesses to use as weights
[in]areatThe array of cell areas [m2]
[in]temperatureThe temperature field to sum
[in]salinityThe salinity field to sum
[in]mesgAn identifying message

Definition at line 747 of file MOM_debugging.F90.

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 

References mom_error_handler::is_root_pe(), and totalstuff().

Here is the call graph for this function:

Variable Documentation

◆ debug

logical mom_debugging::debug = .false.
private

Write out verbose debugging data.

Definition at line 72 of file MOM_debugging.F90.

72 logical :: debug = .false. !< Write out verbose debugging data

Referenced by mom_debugging_init().

◆ debug_chksums

logical mom_debugging::debug_chksums = .true.
private

Perform checksums on arrays.

Definition at line 73 of file MOM_debugging.F90.

73 logical :: debug_chksums = .true. !< Perform checksums on arrays

Referenced by chksum_vec_a2d(), chksum_vec_a3d(), chksum_vec_b2d(), chksum_vec_b3d(), chksum_vec_c2d(), chksum_vec_c3d(), and mom_debugging_init().

◆ debug_redundant

logical mom_debugging::debug_redundant = .true.
private

Check redundant values on PE boundaries.

Definition at line 74 of file MOM_debugging.F90.

74 logical :: debug_redundant = .true. !< Check redundant values on PE boundaries

Referenced by chksum_vec_a2d(), chksum_vec_a3d(), chksum_vec_b2d(), chksum_vec_b3d(), chksum_vec_c2d(), chksum_vec_c3d(), and mom_debugging_init().

◆ max_redundant_prints

integer mom_debugging::max_redundant_prints = 100
private

Maximum number of times to write redundant messages.

Definition at line 70 of file MOM_debugging.F90.

70 integer :: max_redundant_prints = 100 !< Maximum number of times to write redundant messages

Referenced by check_redundant_sb2d(), check_redundant_st2d(), check_redundant_vb2d(), check_redundant_vc2d(), and check_redundant_vt2d().

◆ redundant_prints

integer, dimension(3) mom_debugging::redundant_prints = 0
private

Counters for controlling redundant printing.

Definition at line 71 of file MOM_debugging.F90.

71 integer :: redundant_prints(3) = 0 !< Counters for controlling redundant printing

Referenced by check_redundant_sb2d(), check_redundant_st2d(), check_redundant_vb2d(), check_redundant_vc2d(), and check_redundant_vt2d().