MOM6
mom_checksums Module Reference

Detailed Description

Routines to calculate checksums of various array and vector types.

Data Types

interface  bchksum
 Checksums an array (2d or 3d) staggered at corner points. More...
 
interface  bchksum_pair
 Checksums a pair of arrays (2d or 3d) staggered at corner points. More...
 
interface  chk_sum_msg
 Write a message with either checksums or numerical statistics of arrays. More...
 
interface  chksum
 This is an older interface for 1-, 2-, or 3-D checksums. More...
 
interface  hchksum
 Checksums an array (2d or 3d) staggered at tracer points. More...
 
interface  hchksum_pair
 Checksums a pair of arrays (2d or 3d) staggered at tracer points. More...
 
interface  is_nan
 Returns .true. if any element of x is a NaN, and .false. otherwise. More...
 
interface  qchksum
 This is an older interface that has been renamed Bchksum. More...
 
interface  uchksum
 Checksums an array (2d or 3d) staggered at C-grid u points. More...
 
interface  uvchksum
 Checksums a pair velocity arrays (2d or 3d) staggered at C-grid locations. More...
 
interface  vchksum
 Checksums an array (2d or 3d) staggered at C-grid v points. More...
 

Functions/Subroutines

subroutine, public chksum0 (scalar, mesg, scale, logunit)
 Checksum a scalar field (consistent with array checksums) More...
 
subroutine, public zchksum (array, mesg, scale, logunit)
 Checksum a 1d array (typically a column). More...
 
subroutine chksum_pair_h_2d (mesg, arrayA, arrayB, HI, haloshift, omit_corners, scale, logunit)
 Checksums on a pair of 2d arrays staggered at tracer points. More...
 
subroutine chksum_pair_h_3d (mesg, arrayA, arrayB, HI, haloshift, omit_corners, scale, logunit)
 Checksums on a pair of 3d arrays staggered at tracer points. More...
 
subroutine chksum_h_2d (array, mesg, HI, haloshift, omit_corners, scale, logunit)
 Checksums a 2d array staggered at tracer points. More...
 
subroutine chksum_pair_b_2d (mesg, arrayA, arrayB, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums on a pair of 2d arrays staggered at q-points. More...
 
subroutine chksum_pair_b_3d (mesg, arrayA, arrayB, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums on a pair of 3d arrays staggered at q-points. More...
 
subroutine chksum_b_2d (array, mesg, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a 2d array staggered at corner points. More...
 
subroutine chksum_uv_2d (mesg, arrayU, arrayV, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a pair of 2d velocity arrays staggered at C-grid locations. More...
 
subroutine chksum_uv_3d (mesg, arrayU, arrayV, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a pair of 3d velocity arrays staggered at C-grid locations. More...
 
subroutine chksum_u_2d (array, mesg, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a 2d array staggered at C-grid u points. More...
 
subroutine chksum_v_2d (array, mesg, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a 2d array staggered at C-grid v points. More...
 
subroutine chksum_h_3d (array, mesg, HI, haloshift, omit_corners, scale, logunit)
 Checksums a 3d array staggered at tracer points. More...
 
subroutine chksum_b_3d (array, mesg, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a 3d array staggered at corner points. More...
 
subroutine chksum_u_3d (array, mesg, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a 3d array staggered at C-grid u points. More...
 
subroutine chksum_v_3d (array, mesg, HI, haloshift, symmetric, omit_corners, scale, logunit)
 Checksums a 3d array staggered at C-grid v points. More...
 
subroutine chksum1d (array, mesg, start_i, end_i, compare_PEs)
 chksum1d does a checksum of a 1-dimensional array. More...
 
subroutine chksum2d (array, mesg)
 chksum2d does a checksum of all data in a 2-d array. More...
 
subroutine chksum3d (array, mesg)
 chksum3d does a checksum of all data in a 2-d array. More...
 
logical function is_nan_0d (x)
 This function returns .true. if x is a NaN, and .false. otherwise. More...
 
logical function is_nan_1d (x, skip_mpp)
 Returns .true. if any element of x is a NaN, and .false. otherwise. More...
 
logical function is_nan_2d (x)
 Returns .true. if any element of x is a NaN, and .false. otherwise. More...
 
logical function is_nan_3d (x)
 Returns .true. if any element of x is a NaN, and .false. otherwise. More...
 
subroutine chk_sum_msg1 (fmsg, bc0, mesg, iounit)
 Write a message including the checksum of the non-shifted array. More...
 
subroutine chk_sum_msg5 (fmsg, bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit)
 Write a message including checksums of non-shifted and diagonally shifted arrays. More...
 
subroutine chk_sum_msg_nsew (fmsg, bc0, bcN, bcS, bcE, bcW, mesg, iounit)
 Write a message including checksums of non-shifted and laterally shifted arrays. More...
 
subroutine chk_sum_msg_s (fmsg, bc0, bcS, mesg, iounit)
 Write a message including checksums of non-shifted and southward shifted arrays. More...
 
subroutine chk_sum_msg_w (fmsg, bc0, bcW, mesg, iounit)
 Write a message including checksums of non-shifted and westward shifted arrays. More...
 
subroutine chk_sum_msg2 (fmsg, bc0, bcSW, mesg, iounit)
 Write a message including checksums of non-shifted and southwestward shifted arrays. More...
 
subroutine chk_sum_msg3 (fmsg, aMean, aMin, aMax, mesg, iounit)
 Write a message including the global mean, maximum and minimum of an array. More...
 
subroutine, public mom_checksums_init (param_file)
 MOM_checksums_init initializes the MOM_checksums module. As it happens, the only thing that it does is to log the version of this module. More...
 
subroutine chksum_error (signal, message)
 A wrapper for MOM_error used in the checksum code. More...
 
integer function bitcount (x)
 Does a bitcount of a number by first casting to an integer and then using BTEST to check bit by bit. More...
 

Variables

integer, parameter bc_modulus = 1000000000
 Modulus of checksum bitcount. More...
 
integer, parameter default_shift =0
 The default array shift. More...
 
logical calculatestatistics =.true.
 If true, report min, max and mean. More...
 
logical writechksums =.true.
 If true, report the bitcount checksum. More...
 
logical checkfornans =.true.
 If true, checks array for NaNs and cause FATAL error is any are found. More...
 

Function/Subroutine Documentation

◆ bitcount()

integer function mom_checksums::bitcount ( real, intent(in)  x)
private

Does a bitcount of a number by first casting to an integer and then using BTEST to check bit by bit.

Parameters
[in]xNumber to be bitcount

Definition at line 1855 of file MOM_checksums.F90.

1855  real, intent(in) :: x !< Number to be bitcount
1856 
1857  integer, parameter :: xk = kind(x) !< Kind type of x
1858 
1859  ! NOTE: Assumes that reals and integers of kind=xk are the same size
1860  bitcount = popcnt(transfer(x, 1_xk))

Referenced by chksum0(), chksum1d(), chksum2d(), chksum3d(), and subchk().

Here is the caller graph for this function:

◆ chk_sum_msg1()

subroutine mom_checksums::chk_sum_msg1 ( character(len=*), intent(in)  fmsg,
integer, intent(in)  bc0,
character(len=*), intent(in)  mesg,
integer, intent(in)  iounit 
)
private

Write a message including the checksum of the non-shifted array.

Parameters
[in]fmsgA checksum code-location specific preamble
[in]mesgAn identifying message supplied by top-level caller
[in]bc0The bitcount of the non-shifted array
[in]iounitChecksum logger IO unit

Definition at line 1741 of file MOM_checksums.F90.

1741  character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble
1742  character(len=*), intent(in) :: mesg !< An identifying message supplied by top-level caller
1743  integer, intent(in) :: bc0 !< The bitcount of the non-shifted array
1744  integer, intent(in) :: iounit !< Checksum logger IO unit
1745 
1746  if (is_root_pe()) &
1747  write(iounit, '(A,1(A,I10,X),A)') fmsg, " c=", bc0, trim(mesg)

References mom_error_handler::is_root_pe().

Here is the call graph for this function:

◆ chk_sum_msg2()

subroutine mom_checksums::chk_sum_msg2 ( character(len=*), intent(in)  fmsg,
integer, intent(in)  bc0,
integer, intent(in)  bcSW,
character(len=*), intent(in)  mesg,
integer, intent(in)  iounit 
)
private

Write a message including checksums of non-shifted and southwestward shifted arrays.

Parameters
[in]fmsgA checksum code-location specific preamble
[in]mesgAn identifying message supplied by top-level caller
[in]bc0The bitcount of the non-shifted array
[in]bcswThe bitcount of the southwest-shifted array
[in]iounitChecksum logger IO unit

Definition at line 1806 of file MOM_checksums.F90.

1806  character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble
1807  character(len=*), intent(in) :: mesg !< An identifying message supplied by top-level caller
1808  integer, intent(in) :: bc0 !< The bitcount of the non-shifted array
1809  integer, intent(in) :: bcSW !< The bitcount of the southwest-shifted array
1810  integer, intent(in) :: iounit !< Checksum logger IO unit
1811 
1812  if (is_root_pe()) write(iounit, '(A,2(A,I9,1X),A)') &
1813  fmsg, " c=", bc0, "s/w=", bcsw, trim(mesg)

References mom_error_handler::is_root_pe().

Here is the call graph for this function:

◆ chk_sum_msg3()

subroutine mom_checksums::chk_sum_msg3 ( character(len=*), intent(in)  fmsg,
real, intent(in)  aMean,
real, intent(in)  aMin,
real, intent(in)  aMax,
character(len=*), intent(in)  mesg,
integer, intent(in)  iounit 
)
private

Write a message including the global mean, maximum and minimum of an array.

Parameters
[in]fmsgA checksum code-location specific preamble
[in]mesgAn identifying message supplied by top-level caller
[in]ameanThe mean value of the array
[in]aminThe minimum value of the array
[in]amaxThe maximum value of the array
[in]iounitChecksum logger IO unit

Definition at line 1818 of file MOM_checksums.F90.

1818  character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble
1819  character(len=*), intent(in) :: mesg !< An identifying message supplied by top-level caller
1820  real, intent(in) :: aMean !< The mean value of the array
1821  real, intent(in) :: aMin !< The minimum value of the array
1822  real, intent(in) :: aMax !< The maximum value of the array
1823  integer, intent(in) :: iounit !< Checksum logger IO unit
1824 
1825  ! NOTE: We add zero to aMin and aMax to remove any negative zeros.
1826  ! This is due to inconsistencies of signed zero in local vs MPI calculations.
1827 
1828  if (is_root_pe()) write(iounit, '(A,3(A,ES25.16,1X),A)') &
1829  fmsg, " mean=", amean, "min=", (0. + amin), "max=", (0. + amax), trim(mesg)

References mom_error_handler::is_root_pe().

Here is the call graph for this function:

◆ chk_sum_msg5()

subroutine mom_checksums::chk_sum_msg5 ( character(len=*), intent(in)  fmsg,
integer, intent(in)  bc0,
integer, intent(in)  bcSW,
integer, intent(in)  bcSE,
integer, intent(in)  bcNW,
integer, intent(in)  bcNE,
character(len=*), intent(in)  mesg,
integer, intent(in)  iounit 
)
private

Write a message including checksums of non-shifted and diagonally shifted arrays.

Parameters
[in]fmsgA checksum code-location specific preamble
[in]mesgAn identifying message supplied by top-level caller
[in]bc0The bitcount of the non-shifted array
[in]bcswThe bitcount for SW shifted array
[in]bcseThe bitcount for SE shifted array
[in]bcnwThe bitcount for NW shifted array
[in]bcneThe bitcount for NE shifted array
[in]iounitChecksum logger IO unit

Definition at line 1752 of file MOM_checksums.F90.

1752  character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble
1753  character(len=*), intent(in) :: mesg !< An identifying message supplied by top-level caller
1754  integer, intent(in) :: bc0 !< The bitcount of the non-shifted array
1755  integer, intent(in) :: bcSW !< The bitcount for SW shifted array
1756  integer, intent(in) :: bcSE !< The bitcount for SE shifted array
1757  integer, intent(in) :: bcNW !< The bitcount for NW shifted array
1758  integer, intent(in) :: bcNE !< The bitcount for NE shifted array
1759  integer, intent(in) :: iounit !< Checksum logger IO unit
1760 
1761  if (is_root_pe()) write(iounit, '(A,5(A,I10,1X),A)') &
1762  fmsg, " c=", bc0, "sw=", bcsw, "se=", bcse, "nw=", bcnw, "ne=", bcne, trim(mesg)

References mom_error_handler::is_root_pe().

Here is the call graph for this function:

◆ chk_sum_msg_nsew()

subroutine mom_checksums::chk_sum_msg_nsew ( character(len=*), intent(in)  fmsg,
integer, intent(in)  bc0,
integer, intent(in)  bcN,
integer, intent(in)  bcS,
integer, intent(in)  bcE,
integer, intent(in)  bcW,
character(len=*), intent(in)  mesg,
integer, intent(in)  iounit 
)
private

Write a message including checksums of non-shifted and laterally shifted arrays.

Parameters
[in]fmsgA checksum code-location specific preamble
[in]mesgAn identifying message supplied by top-level caller
[in]bc0The bitcount of the non-shifted array
[in]bcnThe bitcount for N shifted array
[in]bcsThe bitcount for S shifted array
[in]bceThe bitcount for E shifted array
[in]bcwThe bitcount for W shifted array
[in]iounitChecksum logger IO unit

Definition at line 1767 of file MOM_checksums.F90.

1767  character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble
1768  character(len=*), intent(in) :: mesg !< An identifying message supplied by top-level caller
1769  integer, intent(in) :: bc0 !< The bitcount of the non-shifted array
1770  integer, intent(in) :: bcN !< The bitcount for N shifted array
1771  integer, intent(in) :: bcS !< The bitcount for S shifted array
1772  integer, intent(in) :: bcE !< The bitcount for E shifted array
1773  integer, intent(in) :: bcW !< The bitcount for W shifted array
1774  integer, intent(in) :: iounit !< Checksum logger IO unit
1775 
1776  if (is_root_pe()) write(iounit, '(A,5(A,I10,1X),A)') &
1777  fmsg, " c=", bc0, "N=", bcn, "S=", bcs, "E=", bce, "W=", bcw, trim(mesg)

References mom_error_handler::is_root_pe().

Referenced by chksum_b_2d(), chksum_b_3d(), chksum_h_2d(), chksum_h_3d(), chksum_u_2d(), chksum_u_3d(), chksum_v_2d(), and chksum_v_3d().

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

◆ chk_sum_msg_s()

subroutine mom_checksums::chk_sum_msg_s ( character(len=*), intent(in)  fmsg,
integer, intent(in)  bc0,
integer, intent(in)  bcS,
character(len=*), intent(in)  mesg,
integer, intent(in)  iounit 
)
private

Write a message including checksums of non-shifted and southward shifted arrays.

Parameters
[in]fmsgA checksum code-location specific preamble
[in]mesgAn identifying message supplied by top-level caller
[in]bc0The bitcount of the non-shifted array
[in]bcsThe bitcount of the south-shifted array
[in]iounitChecksum logger IO unit

Definition at line 1782 of file MOM_checksums.F90.

1782  character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble
1783  character(len=*), intent(in) :: mesg !< An identifying message supplied by top-level caller
1784  integer, intent(in) :: bc0 !< The bitcount of the non-shifted array
1785  integer, intent(in) :: bcS !< The bitcount of the south-shifted array
1786  integer, intent(in) :: iounit !< Checksum logger IO unit
1787 
1788  if (is_root_pe()) write(iounit, '(A,2(A,I10,1X),A)') &
1789  fmsg, " c=", bc0, "S=", bcs, trim(mesg)

References mom_error_handler::is_root_pe().

Referenced by chksum_v_2d(), and chksum_v_3d().

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

◆ chk_sum_msg_w()

subroutine mom_checksums::chk_sum_msg_w ( character(len=*), intent(in)  fmsg,
integer, intent(in)  bc0,
integer, intent(in)  bcW,
character(len=*), intent(in)  mesg,
integer, intent(in)  iounit 
)
private

Write a message including checksums of non-shifted and westward shifted arrays.

Parameters
[in]fmsgA checksum code-location specific preamble
[in]mesgAn identifying message supplied by top-level caller
[in]bc0The bitcount of the non-shifted array
[in]bcwThe bitcount of the west-shifted array
[in]iounitChecksum logger IO unit

Definition at line 1794 of file MOM_checksums.F90.

1794  character(len=*), intent(in) :: fmsg !< A checksum code-location specific preamble
1795  character(len=*), intent(in) :: mesg !< An identifying message supplied by top-level caller
1796  integer, intent(in) :: bc0 !< The bitcount of the non-shifted array
1797  integer, intent(in) :: bcW !< The bitcount of the west-shifted array
1798  integer, intent(in) :: iounit !< Checksum logger IO unit
1799 
1800  if (is_root_pe()) write(iounit, '(A,2(A,I10,1X),A)') &
1801  fmsg, " c=", bc0, "W=", bcw, trim(mesg)

References mom_error_handler::is_root_pe().

Referenced by chksum_u_2d(), and chksum_u_3d().

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

◆ chksum0()

subroutine, public mom_checksums::chksum0 ( real, intent(in)  scalar,
character(len=*), intent(in)  mesg,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)

Checksum a scalar field (consistent with array checksums)

Parameters
[in]scalarThe array to be checksummed
[in]mesgAn identifying message
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 88 of file MOM_checksums.F90.

88  real, intent(in) :: scalar !< The array to be checksummed
89  character(len=*), intent(in) :: mesg !< An identifying message
90  real, optional, intent(in) :: scale !< A scaling factor for this array.
91  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
92 
93  real :: scaling !< Explicit rescaling factor
94  integer :: iounit !< Log IO unit
95  real :: rs !< Rescaled scalar
96  integer :: bc !< Scalar bitcount
97 
98  if (checkfornans .and. is_nan(scalar)) &
99  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
100 
101  scaling = 1.0 ; if (present(scale)) scaling = scale
102  iounit = error_unit; if(present(logunit)) iounit = logunit
103 
104  if (calculatestatistics) then
105  rs = scaling * scalar
106  if (is_root_pe()) &
107  call chk_sum_msg(" scalar:", rs, rs, rs, mesg, iounit)
108  endif
109 
110  if (.not. writechksums) return
111 
112  bc = mod(bitcount(abs(scaling * scalar)), bc_modulus)
113  if (is_root_pe()) &
114  call chk_sum_msg(" scalar:", bc, mesg, iounit)
115 

References bc_modulus, bitcount(), calculatestatistics, checkfornans, chksum_error(), mom_error_handler::is_root_pe(), and writechksums.

Here is the call graph for this function:

◆ chksum1d()

subroutine mom_checksums::chksum1d ( real, dimension(:), intent(in)  array,
character(len=*), intent(in)  mesg,
integer, intent(in), optional  start_i,
integer, intent(in), optional  end_i,
logical, intent(in), optional  compare_PEs 
)
private

chksum1d does a checksum of a 1-dimensional array.

Parameters
[in]arrayThe array to be summed (index starts at 1).
[in]mesgAn identifying message.
[in]start_iThe starting index for the sum (default 1)
[in]end_iThe ending index for the sum (default all)
[in]compare_pesIf true, compare across PEs instead of summing and list the root_PE value (default true)

Definition at line 1551 of file MOM_checksums.F90.

1551  real, dimension(:), intent(in) :: array !< The array to be summed (index starts at 1).
1552  character(len=*), intent(in) :: mesg !< An identifying message.
1553  integer, optional, intent(in) :: start_i !< The starting index for the sum (default 1)
1554  integer, optional, intent(in) :: end_i !< The ending index for the sum (default all)
1555  logical, optional, intent(in) :: compare_PEs !< If true, compare across PEs instead of summing
1556  !! and list the root_PE value (default true)
1557 
1558  integer :: is, ie, i, bc, sum1, sum_bc
1559  real :: sum
1560  real, allocatable :: sum_here(:)
1561  logical :: compare
1562  integer :: pe_num ! pe number of the data
1563  integer :: nPEs ! Total number of processsors
1564 
1565  is = lbound(array,1) ; ie = ubound(array,1)
1566  if (present(start_i)) is = start_i
1567  if (present(end_i)) ie = end_i
1568  compare = .true. ; if (present(compare_pes)) compare = compare_pes
1569 
1570  sum = 0.0 ; sum_bc = 0
1571  do i=is,ie
1572  sum = sum + array(i)
1573  bc = bitcount(abs(array(i)))
1574  sum_bc = sum_bc + bc
1575  enddo
1576 
1577  pe_num = pe_here() + 1 - root_pe() ; npes = num_pes()
1578  allocate(sum_here(npes)) ; sum_here(:) = 0.0 ; sum_here(pe_num) = sum
1579  call sum_across_pes(sum_here,npes)
1580 
1581  sum1 = sum_bc
1582  call sum_across_pes(sum1)
1583 
1584  if (.not.compare) then
1585  sum = 0.0
1586  do i=1,npes ; sum = sum + sum_here(i) ; enddo
1587  sum_bc = sum1
1588  elseif (is_root_pe()) then
1589  if (sum1 /= npes*sum_bc) &
1590  write(0, '(A40," bitcounts do not match across PEs: ",I12,1X,I12)') &
1591  mesg, sum1, npes*sum_bc
1592  do i=1,npes ; if (sum /= sum_here(i)) then
1593  write(0, '(A40," PE ",i4," sum mismatches root_PE: ",3(ES22.13,1X))') &
1594  mesg, i, sum_here(i), sum, sum_here(i)-sum
1595  endif ; enddo
1596  endif
1597  deallocate(sum_here)
1598 
1599  if (is_root_pe()) &
1600  write(0,'(A50,1X,ES25.16,1X,I12)') mesg, sum, sum_bc
1601 

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

Here is the call graph for this function:

◆ chksum2d()

subroutine mom_checksums::chksum2d ( real, dimension(:,:)  array,
character(len=*)  mesg 
)
private

chksum2d does a checksum of all data in a 2-d array.

Parameters
arrayThe array to be checksummed
mesgAn identifying message

Definition at line 1609 of file MOM_checksums.F90.

1609 
1610  real, dimension(:,:) :: array !< The array to be checksummed
1611  character(len=*) :: mesg !< An identifying message
1612 
1613  integer :: xs,xe,ys,ye,i,j,sum1,bc
1614  real :: sum
1615 
1616  xs = lbound(array,1) ; xe = ubound(array,1)
1617  ys = lbound(array,2) ; ye = ubound(array,2)
1618 
1619  sum = 0.0 ; sum1 = 0
1620  do i=xs,xe ; do j=ys,ye
1621  bc = bitcount(abs(array(i,j)))
1622  sum1 = sum1 + bc
1623  enddo ; enddo
1624  call sum_across_pes(sum1)
1625 
1626  sum = reproducing_sum(array(:,:))
1627 
1628  if (is_root_pe()) &
1629  write(0,'(A50,1X,ES25.16,1X,I12)') mesg, sum, sum1
1630 ! write(0,'(A40,1X,Z16.16,1X,Z16.16,1X,ES25.16,1X,I12)') &
1631 ! mesg, sum, sum1, sum, sum1
1632 

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

Here is the call graph for this function:

◆ chksum3d()

subroutine mom_checksums::chksum3d ( real, dimension(:,:,:)  array,
character(len=*)  mesg 
)
private

chksum3d does a checksum of all data in a 2-d array.

Parameters
arrayThe array to be checksummed
mesgAn identifying message

Definition at line 1637 of file MOM_checksums.F90.

1637 
1638  real, dimension(:,:,:) :: array !< The array to be checksummed
1639  character(len=*) :: mesg !< An identifying message
1640 
1641  integer :: xs,xe,ys,ye,zs,ze,i,j,k, bc,sum1
1642  real :: sum
1643 
1644  xs = lbound(array,1) ; xe = ubound(array,1)
1645  ys = lbound(array,2) ; ye = ubound(array,2)
1646  zs = lbound(array,3) ; ze = ubound(array,3)
1647 
1648  sum = 0.0 ; sum1 = 0
1649  do i=xs,xe ; do j=ys,ye ; do k=zs,ze
1650  bc = bitcount(abs(array(i,j,k)))
1651  sum1 = sum1 + bc
1652  enddo ; enddo ; enddo
1653 
1654  call sum_across_pes(sum1)
1655  sum = reproducing_sum(array(:,:,:))
1656 
1657  if (is_root_pe()) &
1658  write(0,'(A50,1X,ES25.16,1X,I12)') mesg, sum, sum1
1659 ! write(0,'(A40,1X,Z16.16,1X,Z16.16,1X,ES25.16,1X,I12)') &
1660 ! mesg, sum, sum1, sum, sum1
1661 

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

Here is the call graph for this function:

◆ chksum_b_2d()

subroutine mom_checksums::chksum_b_2d ( real, dimension(hi%isdb:,hi%jsdb:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 2d array staggered at corner points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 439 of file MOM_checksums.F90.

439  type(hor_index_type), intent(in) :: HI !< A horizontal index type
440  real, dimension(HI%IsdB:,HI%JsdB:), &
441  intent(in) :: array !< The array to be checksummed
442  character(len=*), intent(in) :: mesg !< An identifying message
443  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
444  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the
445  !! full symmetric computational domain.
446  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
447  real, optional, intent(in) :: scale !< A scaling factor for this array.
448  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
449 
450  real, allocatable, dimension(:,:) :: rescaled_array
451  real :: scaling
452  integer :: iounit !< Log IO unit
453  integer :: i, j, Is, Js
454  real :: aMean, aMin, aMax
455  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
456  integer :: bcN, bcS, bcE, bcW
457  logical :: do_corners, sym, sym_stats
458 
459  if (checkfornans) then
460  if (is_nan(array(hi%IscB:hi%IecB,hi%JscB:hi%JecB))) &
461  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
462 ! if (is_NaN(array)) &
463 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
464  endif
465 
466  scaling = 1.0 ; if (present(scale)) scaling = scale
467  iounit = error_unit; if(present(logunit)) iounit = logunit
468  sym_stats = .false. ; if (present(symmetric)) sym_stats = symmetric
469  if (present(haloshift)) then ; if (haloshift > 0) sym_stats = .true. ; endif
470 
471  if (calculatestatistics) then
472  if (present(scale)) then
473  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
474  lbound(array,2):ubound(array,2)) )
475  rescaled_array(:,:) = 0.0
476  is = hi%isc ; if (sym_stats) is = hi%isc-1
477  js = hi%jsc ; if (sym_stats) js = hi%jsc-1
478  do j=js,hi%JecB ; do i=is,hi%IecB
479  rescaled_array(i,j) = scale*array(i,j)
480  enddo ; enddo
481  call substats(hi, rescaled_array, sym_stats, amean, amin, amax)
482  deallocate(rescaled_array)
483  else
484  call substats(hi, array, sym_stats, amean, amin, amax)
485  endif
486  if (is_root_pe()) &
487  call chk_sum_msg("B-point:", amean, amin, amax, mesg, iounit)
488  endif
489 
490  if (.not.writechksums) return
491 
492  hshift = default_shift
493  if (present(haloshift)) hshift = haloshift
494  if (hshift<0) hshift = hi%ied-hi%iec
495 
496  if ( hi%iscB-hshift<hi%isdB .or. hi%iecB+hshift>hi%iedB .or. &
497  hi%jscB-hshift<hi%jsdB .or. hi%jecB+hshift>hi%jedB ) then
498  write(0,*) 'chksum_B_2d: haloshift =',hshift
499  write(0,*) 'chksum_B_2d: isd,isc,iec,ied=',hi%isdB,hi%iscB,hi%iecB,hi%iedB
500  write(0,*) 'chksum_B_2d: jsd,jsc,jec,jed=',hi%jsdB,hi%jscB,hi%jecB,hi%jedB
501  call chksum_error(fatal,'Error in chksum_B_2d '//trim(mesg))
502  endif
503 
504  bc0 = subchk(array, hi, 0, 0, scaling)
505 
506  sym = .false. ; if (present(symmetric)) sym = symmetric
507 
508  if ((hshift==0) .and. .not.sym) then
509  if (is_root_pe()) call chk_sum_msg("B-point:", bc0, mesg, iounit)
510  return
511  endif
512 
513  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
514 
515  if (do_corners) then
516  if (sym) then
517  bcsw = subchk(array, hi, -hshift-1, -hshift-1, scaling)
518  bcse = subchk(array, hi, hshift, -hshift-1, scaling)
519  bcnw = subchk(array, hi, -hshift-1, hshift, scaling)
520  else
521  bcsw = subchk(array, hi, -hshift, -hshift, scaling)
522  bcse = subchk(array, hi, hshift, -hshift, scaling)
523  bcnw = subchk(array, hi, -hshift, hshift, scaling)
524  endif
525  bcne = subchk(array, hi, hshift, hshift, scaling)
526 
527  if (is_root_pe()) &
528  call chk_sum_msg("B-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
529  else
530  bcs = subchk(array, hi, 0, -hshift, scaling)
531  bce = subchk(array, hi, hshift, 0, scaling)
532  bcw = subchk(array, hi, -hshift, 0, scaling)
533  bcn = subchk(array, hi, 0, hshift, scaling)
534 
535  if (is_root_pe()) &
536  call chk_sum_msg_nsew("B-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
537  endif
538 
539  contains
540 
541  integer function subchk(array, HI, di, dj, scale)
542  type(hor_index_type), intent(in) :: HI !< A horizontal index type
543  real, dimension(HI%IsdB:,HI%JsdB:), intent(in) :: array !< The array to be checksummed
544  integer, intent(in) :: di !< i- direction array shift for this checksum
545  integer, intent(in) :: dj !< j- direction array shift for this checksum
546  real, intent(in) :: scale !< A scaling factor for this array.
547  integer :: i, j, bc
548  subchk = 0
549  ! This line deliberately uses the h-point computational domain.
550  do j=hi%jsc+dj,hi%jec+dj; do i=hi%isc+di,hi%iec+di
551  bc = bitcount(abs(scale*array(i,j)))
552  subchk = subchk + bc
553  enddo ; enddo
554  call sum_across_pes(subchk)
555  subchk=mod(subchk, bc_modulus)
556  end function subchk
557 
558  subroutine substats(HI, array, sym_stats, aMean, aMin, aMax)
559  type(hor_index_type), intent(in) :: HI !< A horizontal index type
560  real, dimension(HI%IsdB:,HI%JsdB:), intent(in) :: array !< The array to be checksummed
561  logical, intent(in) :: sym_stats !< If true, evaluate the statistics on the
562  !! full symmetric computational domain.
563  real, intent(out) :: aMean, aMin, aMax
564 
565  integer :: i, j, n, IsB, JsB
566 
567  isb = hi%isc ; if (sym_stats) isb = hi%isc-1
568  jsb = hi%jsc ; if (sym_stats) jsb = hi%jsc-1
569 
570  amin = array(hi%isc,hi%jsc) ; amax = amin
571  do j=jsb,hi%JecB ; do i=isb,hi%IecB
572  amin = min(amin, array(i,j))
573  amax = max(amax, array(i,j))
574  enddo ; enddo
575  ! This line deliberately uses the h-point computational domain.
576  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec))
577  n = (1 + hi%jec - hi%jsc) * (1 + hi%iec - hi%isc)
578  call sum_across_pes(n)
579  call min_across_pes(amin)
580  call max_across_pes(amax)
581  amean = amean / real(n)
582  end subroutine substats
583 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_pair_b_2d().

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

◆ chksum_b_3d()

subroutine mom_checksums::chksum_b_3d ( real, dimension(hi%isdb:,hi%jsdb:,:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 3d array staggered at corner points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 1085 of file MOM_checksums.F90.

1085  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1086  real, dimension(HI%IsdB:,HI%JsdB:,:), intent(in) :: array !< The array to be checksummed
1087  character(len=*), intent(in) :: mesg !< An identifying message
1088  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
1089  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
1090  !! symmetric computational domain.
1091  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
1092  real, optional, intent(in) :: scale !< A scaling factor for this array.
1093  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
1094 
1095  real, allocatable, dimension(:,:,:) :: rescaled_array
1096  real :: scaling
1097  integer :: iounit !< Log IO unit
1098  integer :: i, j, k, Is, Js
1099  real :: aMean, aMin, aMax
1100  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
1101  integer :: bcN, bcS, bcE, bcW
1102  logical :: do_corners, sym, sym_stats
1103 
1104  if (checkfornans) then
1105  if (is_nan(array(hi%IscB:hi%IecB,hi%JscB:hi%JecB,:))) &
1106  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
1107 ! if (is_NaN(array)) &
1108 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
1109  endif
1110 
1111  scaling = 1.0 ; if (present(scale)) scaling = scale
1112  iounit = error_unit; if(present(logunit)) iounit = logunit
1113  sym_stats = .false. ; if (present(symmetric)) sym_stats = symmetric
1114  if (present(haloshift)) then ; if (haloshift > 0) sym_stats = .true. ; endif
1115 
1116  if (calculatestatistics) then
1117  if (present(scale)) then
1118  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
1119  lbound(array,2):ubound(array,2), &
1120  lbound(array,3):ubound(array,3)) )
1121  rescaled_array(:,:,:) = 0.0
1122  is = hi%isc ; if (sym_stats) is = hi%isc-1
1123  js = hi%jsc ; if (sym_stats) js = hi%jsc-1
1124  do k=1,size(array,3) ; do j=js,hi%JecB ; do i=is,hi%IecB
1125  rescaled_array(i,j,k) = scale*array(i,j,k)
1126  enddo ; enddo ; enddo
1127  call substats(hi, rescaled_array, sym_stats, amean, amin, amax)
1128  deallocate(rescaled_array)
1129  else
1130  call substats(hi, array, sym_stats, amean, amin, amax)
1131  endif
1132 
1133  if (is_root_pe()) &
1134  call chk_sum_msg("B-point:", amean, amin, amax, mesg, iounit)
1135  endif
1136 
1137  if (.not.writechksums) return
1138 
1139  hshift = default_shift
1140  if (present(haloshift)) hshift = haloshift
1141  if (hshift<0) hshift = hi%ied-hi%iec
1142 
1143  if ( hi%isc-hshift<hi%isd .or. hi%iec+hshift>hi%ied .or. &
1144  hi%jsc-hshift<hi%jsd .or. hi%jec+hshift>hi%jed ) then
1145  write(0,*) 'chksum_B_3d: haloshift =',hshift
1146  write(0,*) 'chksum_B_3d: isd,isc,iec,ied=',hi%isd,hi%isc,hi%iec,hi%ied
1147  write(0,*) 'chksum_B_3d: jsd,jsc,jec,jed=',hi%jsd,hi%jsc,hi%jec,hi%jed
1148  call chksum_error(fatal,'Error in chksum_B_3d '//trim(mesg))
1149  endif
1150 
1151  bc0 = subchk(array, hi, 0, 0, scaling)
1152 
1153  sym = .false. ; if (present(symmetric)) sym = symmetric
1154 
1155  if ((hshift==0) .and. .not.sym) then
1156  if (is_root_pe()) call chk_sum_msg("B-point:", bc0, mesg, iounit)
1157  return
1158  endif
1159 
1160  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
1161 
1162  if (do_corners) then
1163  if (sym) then
1164  bcsw = subchk(array, hi, -hshift-1, -hshift-1, scaling)
1165  bcse = subchk(array, hi, hshift, -hshift-1, scaling)
1166  bcnw = subchk(array, hi, -hshift-1, hshift, scaling)
1167  else
1168  bcsw = subchk(array, hi, -hshift-1, -hshift-1, scaling)
1169  bcse = subchk(array, hi, hshift, -hshift-1, scaling)
1170  bcnw = subchk(array, hi, -hshift-1, hshift, scaling)
1171  endif
1172  bcne = subchk(array, hi, hshift, hshift, scaling)
1173 
1174  if (is_root_pe()) &
1175  call chk_sum_msg("B-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
1176  else
1177  if (sym) then
1178  bcs = subchk(array, hi, 0, -hshift-1, scaling)
1179  bcw = subchk(array, hi, -hshift-1, 0, scaling)
1180  else
1181  bcs = subchk(array, hi, 0, -hshift, scaling)
1182  bcw = subchk(array, hi, -hshift, 0, scaling)
1183  endif
1184  bce = subchk(array, hi, hshift, 0, scaling)
1185  bcn = subchk(array, hi, 0, hshift, scaling)
1186 
1187  if (is_root_pe()) &
1188  call chk_sum_msg_nsew("B-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
1189  endif
1190 
1191  contains
1192 
1193  integer function subchk(array, HI, di, dj, scale)
1194  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1195  real, dimension(HI%IsdB:,HI%JsdB:,:), intent(in) :: array !< The array to be checksummed
1196  integer, intent(in) :: di !< i- direction array shift for this checksum
1197  integer, intent(in) :: dj !< j- direction array shift for this checksum
1198  real, intent(in) :: scale !< A scaling factor for this array.
1199  integer :: i, j, k, bc
1200  subchk = 0
1201  ! This line deliberately uses the h-point computational domain.
1202  do k=lbound(array,3),ubound(array,3) ; do j=hi%jsc+dj,hi%jec+dj ; do i=hi%isc+di,hi%iec+di
1203  bc = bitcount(abs(scale*array(i,j,k)))
1204  subchk = subchk + bc
1205  enddo ; enddo ; enddo
1206  call sum_across_pes(subchk)
1207  subchk=mod(subchk, bc_modulus)
1208  end function subchk
1209 
1210  subroutine substats(HI, array, sym_stats, aMean, aMin, aMax)
1211  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1212  real, dimension(HI%IsdB:,HI%JsdB:,:), intent(in) :: array !< The array to be checksummed
1213  logical, intent(in) :: sym_stats !< If true, evaluate the statistics on the
1214  !! full symmetric computational domain.
1215  real, intent(out) :: aMean, aMin, aMax
1216 
1217  integer :: i, j, k, n, IsB, JsB
1218 
1219  isb = hi%isc ; if (sym_stats) isb = hi%isc-1
1220  jsb = hi%jsc ; if (sym_stats) jsb = hi%jsc-1
1221 
1222  amin = array(hi%isc,hi%jsc,1) ; amax = amin
1223  do k=lbound(array,3),ubound(array,3) ; do j=jsb,hi%JecB ; do i=isb,hi%IecB
1224  amin = min(amin, array(i,j,k))
1225  amax = max(amax, array(i,j,k))
1226  enddo ; enddo ; enddo
1227  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec,:))
1228  n = (1 + hi%jec - hi%jsc) * (1 + hi%iec - hi%isc) * size(array,3)
1229  call sum_across_pes(n)
1230  call min_across_pes(amin)
1231  call max_across_pes(amax)
1232  amean = amean / real(n)
1233  end subroutine substats
1234 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_pair_b_3d().

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

◆ chksum_error()

subroutine mom_checksums::chksum_error ( integer, intent(in)  signal,
character(len=*), intent(in)  message 
)
private

A wrapper for MOM_error used in the checksum code.

Parameters
[in]signalAn error severity level, such as FATAL or WARNING
[in]messageAn error message

Definition at line 1846 of file MOM_checksums.F90.

1846  ! Wrapper for MOM_error to help place specific break points in debuggers
1847  integer, intent(in) :: signal !< An error severity level, such as FATAL or WARNING
1848  character(len=*), intent(in) :: message !< An error message
1849  call mom_error(signal, message)

References mom_error_handler::mom_error().

Referenced by chksum0(), chksum_b_2d(), chksum_b_3d(), chksum_h_2d(), chksum_h_3d(), chksum_u_2d(), chksum_u_3d(), chksum_v_2d(), chksum_v_3d(), and zchksum().

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

◆ chksum_h_2d()

subroutine mom_checksums::chksum_h_2d ( real, dimension(hi%isd:,hi%jsd:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 2d array staggered at tracer points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 247 of file MOM_checksums.F90.

247  type(hor_index_type), intent(in) :: HI !< A horizontal index type
248  real, dimension(HI%isd:,HI%jsd:), intent(in) :: array !< The array to be checksummed
249  character(len=*), intent(in) :: mesg !< An identifying message
250  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
251  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
252  real, optional, intent(in) :: scale !< A scaling factor for this array.
253  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
254 
255  real, allocatable, dimension(:,:) :: rescaled_array
256  real :: scaling
257  integer :: iounit !< Log IO unit
258  integer :: i, j
259  real :: aMean, aMin, aMax
260  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
261  integer :: bcN, bcS, bcE, bcW
262  logical :: do_corners
263 
264  if (checkfornans) then
265  if (is_nan(array(hi%isc:hi%iec,hi%jsc:hi%jec))) &
266  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
267 ! if (is_NaN(array)) &
268 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
269  endif
270 
271  scaling = 1.0 ; if (present(scale)) scaling = scale
272  iounit = error_unit; if(present(logunit)) iounit = logunit
273 
274  if (calculatestatistics) then
275  if (present(scale)) then
276  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
277  lbound(array,2):ubound(array,2)) )
278  rescaled_array(:,:) = 0.0
279  do j=hi%jsc,hi%jec ; do i=hi%isc,hi%iec
280  rescaled_array(i,j) = scale*array(i,j)
281  enddo ; enddo
282  call substats(hi, rescaled_array, amean, amin, amax)
283  deallocate(rescaled_array)
284  else
285  call substats(hi, array, amean, amin, amax)
286  endif
287 
288  if (is_root_pe()) &
289  call chk_sum_msg("h-point:", amean, amin, amax, mesg, iounit)
290  endif
291 
292  if (.not.writechksums) return
293 
294  hshift = default_shift
295  if (present(haloshift)) hshift = haloshift
296  if (hshift<0) hshift = hi%ied-hi%iec
297 
298  if ( hi%isc-hshift<hi%isd .or. hi%iec+hshift>hi%ied .or. &
299  hi%jsc-hshift<hi%jsd .or. hi%jec+hshift>hi%jed ) then
300  write(0,*) 'chksum_h_2d: haloshift =',hshift
301  write(0,*) 'chksum_h_2d: isd,isc,iec,ied=',hi%isd,hi%isc,hi%iec,hi%ied
302  write(0,*) 'chksum_h_2d: jsd,jsc,jec,jed=',hi%jsd,hi%jsc,hi%jec,hi%jed
303  call chksum_error(fatal,'Error in chksum_h_2d '//trim(mesg))
304  endif
305 
306  bc0 = subchk(array, hi, 0, 0, scaling)
307 
308  if (hshift==0) then
309  if (is_root_pe()) call chk_sum_msg("h-point:", bc0, mesg, iounit)
310  return
311  endif
312 
313  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
314 
315  if (do_corners) then
316  bcsw = subchk(array, hi, -hshift, -hshift, scaling)
317  bcse = subchk(array, hi, hshift, -hshift, scaling)
318  bcnw = subchk(array, hi, -hshift, hshift, scaling)
319  bcne = subchk(array, hi, hshift, hshift, scaling)
320 
321  if (is_root_pe()) &
322  call chk_sum_msg("h-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
323  else
324  bcs = subchk(array, hi, 0, -hshift, scaling)
325  bce = subchk(array, hi, hshift, 0, scaling)
326  bcw = subchk(array, hi, -hshift, 0, scaling)
327  bcn = subchk(array, hi, 0, hshift, scaling)
328 
329  if (is_root_pe()) &
330  call chk_sum_msg_nsew("h-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
331  endif
332 
333  contains
334  integer function subchk(array, HI, di, dj, scale)
335  type(hor_index_type), intent(in) :: HI !< A horizontal index type
336  real, dimension(HI%isd:,HI%jsd:), intent(in) :: array !< The array to be checksummed
337  integer, intent(in) :: di !< i- direction array shift for this checksum
338  integer, intent(in) :: dj !< j- direction array shift for this checksum
339  real, intent(in) :: scale !< A scaling factor for this array.
340  integer :: i, j, bc
341  subchk = 0
342  do j=hi%jsc+dj,hi%jec+dj; do i=hi%isc+di,hi%iec+di
343  bc = bitcount(abs(scale*array(i,j)))
344  subchk = subchk + bc
345  enddo ; enddo
346  call sum_across_pes(subchk)
347  subchk=mod(subchk, bc_modulus)
348  end function subchk
349 
350  subroutine substats(HI, array, aMean, aMin, aMax)
351  type(hor_index_type), intent(in) :: HI !< A horizontal index type
352  real, dimension(HI%isd:,HI%jsd:), intent(in) :: array !< The array to be checksummed
353  real, intent(out) :: aMean, aMin, aMax
354 
355  integer :: i, j, n
356 
357  amin = array(hi%isc,hi%jsc)
358  amax = array(hi%isc,hi%jsc)
359  n = 0
360  do j=hi%jsc,hi%jec ; do i=hi%isc,hi%iec
361  amin = min(amin, array(i,j))
362  amax = max(amax, array(i,j))
363  n = n + 1
364  enddo ; enddo
365  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec))
366  call sum_across_pes(n)
367  call min_across_pes(amin)
368  call max_across_pes(amax)
369  amean = amean / real(n)
370  end subroutine substats
371 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_pair_h_2d().

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

◆ chksum_h_3d()

subroutine mom_checksums::chksum_h_3d ( real, dimension(hi%isd:,hi%jsd:,:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 3d array staggered at tracer points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 952 of file MOM_checksums.F90.

952  type(hor_index_type), intent(in) :: HI !< A horizontal index type
953  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: array !< The array to be checksummed
954  character(len=*), intent(in) :: mesg !< An identifying message
955  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
956  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
957  real, optional, intent(in) :: scale !< A scaling factor for this array.
958  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
959 
960  real, allocatable, dimension(:,:,:) :: rescaled_array
961  real :: scaling
962  integer :: iounit !< Log IO unit
963  integer :: i, j, k
964  real :: aMean, aMin, aMax
965  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
966  integer :: bcN, bcS, bcE, bcW
967  logical :: do_corners
968 
969  if (checkfornans) then
970  if (is_nan(array(hi%isc:hi%iec,hi%jsc:hi%jec,:))) &
971  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
972 ! if (is_NaN(array)) &
973 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
974  endif
975 
976  scaling = 1.0 ; if (present(scale)) scaling = scale
977  iounit = error_unit; if(present(logunit)) iounit = logunit
978 
979  if (calculatestatistics) then
980  if (present(scale)) then
981  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
982  lbound(array,2):ubound(array,2), &
983  lbound(array,3):ubound(array,3)) )
984  rescaled_array(:,:,:) = 0.0
985  do k=1,size(array,3) ; do j=hi%jsc,hi%jec ; do i=hi%isc,hi%iec
986  rescaled_array(i,j,k) = scale*array(i,j,k)
987  enddo ; enddo ; enddo
988 
989  call substats(hi, rescaled_array, amean, amin, amax)
990  deallocate(rescaled_array)
991  else
992  call substats(hi, array, amean, amin, amax)
993  endif
994 
995  if (is_root_pe()) &
996  call chk_sum_msg("h-point:", amean, amin, amax, mesg, iounit)
997  endif
998 
999  if (.not.writechksums) return
1000 
1001  hshift = default_shift
1002  if (present(haloshift)) hshift = haloshift
1003  if (hshift<0) hshift = hi%ied-hi%iec
1004 
1005  if ( hi%isc-hshift<hi%isd .or. hi%iec+hshift>hi%ied .or. &
1006  hi%jsc-hshift<hi%jsd .or. hi%jec+hshift>hi%jed ) then
1007  write(0,*) 'chksum_h_3d: haloshift =',hshift
1008  write(0,*) 'chksum_h_3d: isd,isc,iec,ied=',hi%isd,hi%isc,hi%iec,hi%ied
1009  write(0,*) 'chksum_h_3d: jsd,jsc,jec,jed=',hi%jsd,hi%jsc,hi%jec,hi%jed
1010  call chksum_error(fatal,'Error in chksum_h_3d '//trim(mesg))
1011  endif
1012 
1013  bc0 = subchk(array, hi, 0, 0, scaling)
1014 
1015  if (hshift==0) then
1016  if (is_root_pe()) call chk_sum_msg("h-point:", bc0, mesg, iounit)
1017  return
1018  endif
1019 
1020  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
1021 
1022  if (do_corners) then
1023  bcsw = subchk(array, hi, -hshift, -hshift, scaling)
1024  bcse = subchk(array, hi, hshift, -hshift, scaling)
1025  bcnw = subchk(array, hi, -hshift, hshift, scaling)
1026  bcne = subchk(array, hi, hshift, hshift, scaling)
1027 
1028  if (is_root_pe()) &
1029  call chk_sum_msg("h-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
1030  else
1031  bcs = subchk(array, hi, 0, -hshift, scaling)
1032  bce = subchk(array, hi, hshift, 0, scaling)
1033  bcw = subchk(array, hi, -hshift, 0, scaling)
1034  bcn = subchk(array, hi, 0, hshift, scaling)
1035 
1036  if (is_root_pe()) &
1037  call chk_sum_msg_nsew("h-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
1038  endif
1039 
1040  contains
1041 
1042  integer function subchk(array, HI, di, dj, scale)
1043  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1044  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: array !< The array to be checksummed
1045  integer, intent(in) :: di !< i- direction array shift for this checksum
1046  integer, intent(in) :: dj !< j- direction array shift for this checksum
1047  real, intent(in) :: scale !< A scaling factor for this array.
1048  integer :: i, j, k, bc
1049  subchk = 0
1050  do k=lbound(array,3),ubound(array,3) ; do j=hi%jsc+dj,hi%jec+dj ; do i=hi%isc+di,hi%iec+di
1051  bc = bitcount(abs(scale*array(i,j,k)))
1052  subchk = subchk + bc
1053  enddo ; enddo ; enddo
1054  call sum_across_pes(subchk)
1055  subchk=mod(subchk, bc_modulus)
1056  end function subchk
1057 
1058  subroutine substats(HI, array, aMean, aMin, aMax)
1059  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1060  real, dimension(HI%isd:,HI%jsd:,:), intent(in) :: array !< The array to be checksummed
1061  real, intent(out) :: aMean, aMin, aMax
1062 
1063  integer :: i, j, k, n
1064 
1065  amin = array(hi%isc,hi%jsc,1)
1066  amax = array(hi%isc,hi%jsc,1)
1067  n = 0
1068  do k=lbound(array,3),ubound(array,3) ; do j=hi%jsc,hi%jec ; do i=hi%isc,hi%iec
1069  amin = min(amin, array(i,j,k))
1070  amax = max(amax, array(i,j,k))
1071  n = n + 1
1072  enddo ; enddo ; enddo
1073  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec,:))
1074  call sum_across_pes(n)
1075  call min_across_pes(amin)
1076  call max_across_pes(amax)
1077  amean = amean / real(n)
1078  end subroutine substats
1079 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_pair_h_3d().

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

◆ chksum_pair_b_2d()

subroutine mom_checksums::chksum_pair_b_2d ( character(len=*), intent(in)  mesg,
real, dimension(hi%isd:,hi%jsd:), intent(in)  arrayA,
real, dimension(hi%isd:,hi%jsd:), intent(in)  arrayB,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums on a pair of 2d arrays staggered at q-points.

Parameters
[in]mesgIdentifying messages
[in]hiA horizontal index type
[in]arrayaThe first array to be checksummed
[in]arraybThe second array to be checksummed
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]haloshiftThe width of halos to check (default 0)
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 377 of file MOM_checksums.F90.

377  character(len=*), intent(in) :: mesg !< Identifying messages
378  type(hor_index_type), intent(in) :: HI !< A horizontal index type
379  real, dimension(HI%isd:,HI%jsd:), intent(in) :: arrayA !< The first array to be checksummed
380  real, dimension(HI%isd:,HI%jsd:), intent(in) :: arrayB !< The second array to be checksummed
381  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
382  !! symmetric computational domain.
383  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
384  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
385  real, optional, intent(in) :: scale !< A scaling factor for this array.
386  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
387 
388  logical :: sym
389 
390  sym = .false. ; if (present(symmetric)) sym = symmetric
391 
392  if (present(haloshift)) then
393  call chksum_b_2d(arraya, 'x '//mesg, hi, haloshift, symmetric=sym, &
394  omit_corners=omit_corners, scale=scale, logunit=logunit)
395  call chksum_b_2d(arrayb, 'y '//mesg, hi, haloshift, symmetric=sym, &
396  omit_corners=omit_corners, scale=scale, logunit=logunit)
397  else
398  call chksum_b_2d(arraya, 'x '//mesg, hi, symmetric=sym, scale=scale, &
399  logunit=logunit)
400  call chksum_b_2d(arrayb, 'y '//mesg, hi, symmetric=sym, scale=scale, &
401  logunit=logunit)
402  endif
403 

References chksum_b_2d().

Here is the call graph for this function:

◆ chksum_pair_b_3d()

subroutine mom_checksums::chksum_pair_b_3d ( character(len=*), intent(in)  mesg,
real, dimension(hi%isdb:,hi%jsdb:, :), intent(in)  arrayA,
real, dimension(hi%isdb:,hi%jsdb:, :), intent(in)  arrayB,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums on a pair of 3d arrays staggered at q-points.

Parameters
[in]mesgIdentifying messages
[in]hiA horizontal index type
[in]arrayaThe first array to be checksummed
[in]arraybThe second array to be checksummed
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 409 of file MOM_checksums.F90.

409  character(len=*), intent(in) :: mesg !< Identifying messages
410  type(hor_index_type), intent(in) :: HI !< A horizontal index type
411  real, dimension(HI%IsdB:,HI%JsdB:, :), intent(in) :: arrayA !< The first array to be checksummed
412  real, dimension(HI%IsdB:,HI%JsdB:, :), intent(in) :: arrayB !< The second array to be checksummed
413  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
414  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
415  !! symmetric computational domain.
416  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
417  real, optional, intent(in) :: scale !< A scaling factor for this array.
418  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
419 
420  logical :: sym
421 
422  if (present(haloshift)) then
423  call chksum_b_3d(arraya, 'x '//mesg, hi, haloshift, symmetric, &
424  omit_corners, scale=scale, logunit=logunit)
425  call chksum_b_3d(arrayb, 'y '//mesg, hi, haloshift, symmetric, &
426  omit_corners, scale=scale, logunit=logunit)
427  else
428  call chksum_b_3d(arraya, 'x '//mesg, hi, symmetric=symmetric, scale=scale, &
429  logunit=logunit)
430  call chksum_b_3d(arrayb, 'y '//mesg, hi, symmetric=symmetric, scale=scale, &
431  logunit=logunit)
432  endif
433 

References chksum_b_3d().

Here is the call graph for this function:

◆ chksum_pair_h_2d()

subroutine mom_checksums::chksum_pair_h_2d ( character(len=*), intent(in)  mesg,
real, dimension(hi%isd:,hi%jsd:), intent(in)  arrayA,
real, dimension(hi%isd:,hi%jsd:), intent(in)  arrayB,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums on a pair of 2d arrays staggered at tracer points.

Parameters
[in]mesgIdentifying messages
[in]hiA horizontal index type
[in]arrayaThe first array to be checksummed
[in]arraybThe second array to be checksummed
[in]haloshiftThe width of halos to check (default 0)
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 200 of file MOM_checksums.F90.

200  character(len=*), intent(in) :: mesg !< Identifying messages
201  type(hor_index_type), intent(in) :: HI !< A horizontal index type
202  real, dimension(HI%isd:,HI%jsd:), intent(in) :: arrayA !< The first array to be checksummed
203  real, dimension(HI%isd:,HI%jsd:), intent(in) :: arrayB !< The second array to be checksummed
204  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
205  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
206  real, optional, intent(in) :: scale !< A scaling factor for this array.
207  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
208 
209  if (present(haloshift)) then
210  call chksum_h_2d(arraya, 'x '//mesg, hi, haloshift, omit_corners, &
211  scale=scale, logunit=logunit)
212  call chksum_h_2d(arrayb, 'y '//mesg, hi, haloshift, omit_corners, &
213  scale=scale, logunit=logunit)
214  else
215  call chksum_h_2d(arraya, 'x '//mesg, hi, scale=scale, logunit=logunit)
216  call chksum_h_2d(arrayb, 'y '//mesg, hi, scale=scale, logunit=logunit)
217  endif
218 

References chksum_h_2d().

Here is the call graph for this function:

◆ chksum_pair_h_3d()

subroutine mom_checksums::chksum_pair_h_3d ( character(len=*), intent(in)  mesg,
real, dimension(hi%isd:,hi%jsd:, :), intent(in)  arrayA,
real, dimension(hi%isd:,hi%jsd:, :), intent(in)  arrayB,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums on a pair of 3d arrays staggered at tracer points.

Parameters
[in]mesgIdentifying messages
[in]hiA horizontal index type
[in]arrayaThe first array to be checksummed
[in]arraybThe second array to be checksummed
[in]haloshiftThe width of halos to check (default 0)
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 224 of file MOM_checksums.F90.

224  character(len=*), intent(in) :: mesg !< Identifying messages
225  type(hor_index_type), intent(in) :: HI !< A horizontal index type
226  real, dimension(HI%isd:,HI%jsd:, :), intent(in) :: arrayA !< The first array to be checksummed
227  real, dimension(HI%isd:,HI%jsd:, :), intent(in) :: arrayB !< The second array to be checksummed
228  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
229  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
230  real, optional, intent(in) :: scale !< A scaling factor for this array.
231  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
232 
233  if (present(haloshift)) then
234  call chksum_h_3d(arraya, 'x '//mesg, hi, haloshift, omit_corners, &
235  scale=scale, logunit=logunit)
236  call chksum_h_3d(arrayb, 'y '//mesg, hi, haloshift, omit_corners, &
237  scale=scale, logunit=logunit)
238  else
239  call chksum_h_3d(arraya, 'x '//mesg, hi, scale=scale, logunit=logunit)
240  call chksum_h_3d(arrayb, 'y '//mesg, hi, scale=scale, logunit=logunit)
241  endif
242 

References chksum_h_3d().

Here is the call graph for this function:

◆ chksum_u_2d()

subroutine mom_checksums::chksum_u_2d ( real, dimension(hi%isdb:,hi%jsd:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 2d array staggered at C-grid u points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 645 of file MOM_checksums.F90.

645  type(hor_index_type), intent(in) :: HI !< A horizontal index type
646  real, dimension(HI%IsdB:,HI%jsd:), intent(in) :: array !< The array to be checksummed
647  character(len=*), intent(in) :: mesg !< An identifying message
648  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
649  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
650  !! symmetric computational domain.
651  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
652  real, optional, intent(in) :: scale !< A scaling factor for this array.
653  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
654 
655  real, allocatable, dimension(:,:) :: rescaled_array
656  real :: scaling
657  integer :: iounit !< Log IO unit
658  integer :: i, j, Is
659  real :: aMean, aMin, aMax
660  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
661  integer :: bcN, bcS, bcE, bcW
662  logical :: do_corners, sym, sym_stats
663 
664  if (checkfornans) then
665  if (is_nan(array(hi%IscB:hi%IecB,hi%jsc:hi%jec))) &
666  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
667 ! if (is_NaN(array)) &
668 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
669  endif
670 
671  scaling = 1.0 ; if (present(scale)) scaling = scale
672  iounit = error_unit; if(present(logunit)) iounit = logunit
673  sym_stats = .false. ; if (present(symmetric)) sym_stats = symmetric
674  if (present(haloshift)) then ; if (haloshift > 0) sym_stats = .true. ; endif
675 
676  if (calculatestatistics) then
677  if (present(scale)) then
678  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
679  lbound(array,2):ubound(array,2)) )
680  rescaled_array(:,:) = 0.0
681  is = hi%isc ; if (sym_stats) is = hi%isc-1
682  do j=hi%jsc,hi%jec ; do i=is,hi%IecB
683  rescaled_array(i,j) = scale*array(i,j)
684  enddo ; enddo
685  call substats(hi, rescaled_array, sym_stats, amean, amin, amax)
686  deallocate(rescaled_array)
687  else
688  call substats(hi, array, sym_stats, amean, amin, amax)
689  endif
690 
691  if (is_root_pe()) &
692  call chk_sum_msg("u-point:", amean, amin, amax, mesg, iounit)
693  endif
694 
695  if (.not.writechksums) return
696 
697  hshift = default_shift
698  if (present(haloshift)) hshift = haloshift
699  if (hshift<0) hshift = hi%iedB-hi%iecB
700 
701  if ( hi%iscB-hshift<hi%isdB .or. hi%iecB+hshift>hi%iedB .or. &
702  hi%jsc-hshift<hi%jsd .or. hi%jec+hshift>hi%jed ) then
703  write(0,*) 'chksum_u_2d: haloshift =',hshift
704  write(0,*) 'chksum_u_2d: isd,isc,iec,ied=',hi%isdB,hi%iscB,hi%iecB,hi%iedB
705  write(0,*) 'chksum_u_2d: jsd,jsc,jec,jed=',hi%jsd,hi%jsc,hi%jec,hi%jed
706  call chksum_error(fatal,'Error in chksum_u_2d '//trim(mesg))
707  endif
708 
709  bc0 = subchk(array, hi, 0, 0, scaling)
710 
711  sym = .false. ; if (present(symmetric)) sym = symmetric
712 
713  if ((hshift==0) .and. .not.sym) then
714  if (is_root_pe()) call chk_sum_msg("u-point:", bc0, mesg, iounit)
715  return
716  endif
717 
718  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
719 
720  if (hshift==0) then
721  bcw = subchk(array, hi, -hshift-1, 0, scaling)
722  if (is_root_pe()) call chk_sum_msg_w("u-point:", bc0, bcw, mesg, iounit)
723  elseif (do_corners) then
724  if (sym) then
725  bcsw = subchk(array, hi, -hshift-1, -hshift, scaling)
726  bcnw = subchk(array, hi, -hshift-1, hshift, scaling)
727  else
728  bcsw = subchk(array, hi, -hshift, -hshift, scaling)
729  bcnw = subchk(array, hi, -hshift, hshift, scaling)
730  endif
731  bcse = subchk(array, hi, hshift, -hshift, scaling)
732  bcne = subchk(array, hi, hshift, hshift, scaling)
733 
734  if (is_root_pe()) &
735  call chk_sum_msg("u-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
736  else
737  bcs = subchk(array, hi, 0, -hshift, scaling)
738  bce = subchk(array, hi, hshift, 0, scaling)
739  if (sym) then
740  bcw = subchk(array, hi, -hshift-1, 0, scaling)
741  else
742  bcw = subchk(array, hi, -hshift, 0, scaling)
743  endif
744  bcn = subchk(array, hi, 0, hshift, scaling)
745 
746  if (is_root_pe()) &
747  call chk_sum_msg_nsew("u-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
748  endif
749 
750  contains
751 
752  integer function subchk(array, HI, di, dj, scale)
753  type(hor_index_type), intent(in) :: HI !< A horizontal index type
754  real, dimension(HI%IsdB:,HI%jsd:), intent(in) :: array !< The array to be checksummed
755  integer, intent(in) :: di !< i- direction array shift for this checksum
756  integer, intent(in) :: dj !< j- direction array shift for this checksum
757  real, intent(in) :: scale !< A scaling factor for this array.
758  integer :: i, j, bc
759  subchk = 0
760  ! This line deliberately uses the h-point computational domain.
761  do j=hi%jsc+dj,hi%jec+dj; do i=hi%isc+di,hi%iec+di
762  bc = bitcount(abs(scale*array(i,j)))
763  subchk = subchk + bc
764  enddo ; enddo
765  call sum_across_pes(subchk)
766  subchk=mod(subchk, bc_modulus)
767  end function subchk
768 
769  subroutine substats(HI, array, sym_stats, aMean, aMin, aMax)
770  type(hor_index_type), intent(in) :: HI !< A horizontal index type
771  real, dimension(HI%IsdB:,HI%jsd:), intent(in) :: array !< The array to be checksummed
772  logical, intent(in) :: sym_stats !< If true, evaluate the statistics on the
773  !! full symmetric computational domain.
774  real, intent(out) :: aMean, aMin, aMax
775 
776  integer :: i, j, n, IsB
777 
778  isb = hi%isc ; if (sym_stats) isb = hi%isc-1
779 
780  amin = array(hi%isc,hi%jsc) ; amax = amin
781  do j=hi%jsc,hi%jec ; do i=isb,hi%IecB
782  amin = min(amin, array(i,j))
783  amax = max(amax, array(i,j))
784  enddo ; enddo
785  ! This line deliberately uses the h-point computational domain.
786  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec))
787  n = (1 + hi%jec - hi%jsc) * (1 + hi%iec - hi%isc)
788  call sum_across_pes(n)
789  call min_across_pes(amin)
790  call max_across_pes(amax)
791  amean = amean / real(n)
792  end subroutine substats
793 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chk_sum_msg_w(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_uv_2d().

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

◆ chksum_u_3d()

subroutine mom_checksums::chksum_u_3d ( real, dimension(hi%isdb:,hi%jsd:,:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 3d array staggered at C-grid u points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 1240 of file MOM_checksums.F90.

1240  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1241  real, dimension(HI%isdB:,HI%Jsd:,:), intent(in) :: array !< The array to be checksummed
1242  character(len=*), intent(in) :: mesg !< An identifying message
1243  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
1244  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
1245  !! symmetric computational domain.
1246  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
1247  real, optional, intent(in) :: scale !< A scaling factor for this array.
1248  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
1249 
1250  real, allocatable, dimension(:,:,:) :: rescaled_array
1251  real :: scaling
1252  integer :: iounit !< Log IO unit
1253  integer :: i, j, k, Is
1254  real :: aMean, aMin, aMax
1255  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
1256  integer :: bcN, bcS, bcE, bcW
1257  logical :: do_corners, sym, sym_stats
1258 
1259  if (checkfornans) then
1260  if (is_nan(array(hi%IscB:hi%IecB,hi%jsc:hi%jec,:))) &
1261  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
1262 ! if (is_NaN(array)) &
1263 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
1264  endif
1265 
1266  scaling = 1.0 ; if (present(scale)) scaling = scale
1267  iounit = error_unit; if(present(logunit)) iounit = logunit
1268  sym_stats = .false. ; if (present(symmetric)) sym_stats = symmetric
1269  if (present(haloshift)) then ; if (haloshift > 0) sym_stats = .true. ; endif
1270 
1271  if (calculatestatistics) then
1272  if (present(scale)) then
1273  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
1274  lbound(array,2):ubound(array,2), &
1275  lbound(array,3):ubound(array,3)) )
1276  rescaled_array(:,:,:) = 0.0
1277  is = hi%isc ; if (sym_stats) is = hi%isc-1
1278  do k=1,size(array,3) ; do j=hi%jsc,hi%jec ; do i=is,hi%IecB
1279  rescaled_array(i,j,k) = scale*array(i,j,k)
1280  enddo ; enddo ; enddo
1281  call substats(hi, rescaled_array, sym_stats, amean, amin, amax)
1282  deallocate(rescaled_array)
1283  else
1284  call substats(hi, array, sym_stats, amean, amin, amax)
1285  endif
1286  if (is_root_pe()) &
1287  call chk_sum_msg("u-point:", amean, amin, amax, mesg, iounit)
1288  endif
1289 
1290  if (.not.writechksums) return
1291 
1292  hshift = default_shift
1293  if (present(haloshift)) hshift = haloshift
1294  if (hshift<0) hshift = hi%ied-hi%iec
1295 
1296  if ( hi%isc-hshift<hi%isd .or. hi%iec+hshift>hi%ied .or. &
1297  hi%jsc-hshift<hi%jsd .or. hi%jec+hshift>hi%jed ) then
1298  write(0,*) 'chksum_u_3d: haloshift =',hshift
1299  write(0,*) 'chksum_u_3d: isd,isc,iec,ied=',hi%isd,hi%isc,hi%iec,hi%ied
1300  write(0,*) 'chksum_u_3d: jsd,jsc,jec,jed=',hi%jsd,hi%jsc,hi%jec,hi%jed
1301  call chksum_error(fatal,'Error in chksum_u_3d '//trim(mesg))
1302  endif
1303 
1304  bc0 = subchk(array, hi, 0, 0, scaling)
1305 
1306  sym = .false. ; if (present(symmetric)) sym = symmetric
1307 
1308  if ((hshift==0) .and. .not.sym) then
1309  if (is_root_pe()) call chk_sum_msg("u-point:", bc0, mesg, iounit)
1310  return
1311  endif
1312 
1313  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
1314 
1315  if (hshift==0) then
1316  bcw = subchk(array, hi, -hshift-1, 0, scaling)
1317  if (is_root_pe()) call chk_sum_msg_w("u-point:", bc0, bcw, mesg, iounit)
1318  elseif (do_corners) then
1319  if (sym) then
1320  bcsw = subchk(array, hi, -hshift-1, -hshift, scaling)
1321  bcnw = subchk(array, hi, -hshift-1, hshift, scaling)
1322  else
1323  bcsw = subchk(array, hi, -hshift, -hshift, scaling)
1324  bcnw = subchk(array, hi, -hshift, hshift, scaling)
1325  endif
1326  bcse = subchk(array, hi, hshift, -hshift, scaling)
1327  bcne = subchk(array, hi, hshift, hshift, scaling)
1328 
1329  if (is_root_pe()) &
1330  call chk_sum_msg("u-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
1331  else
1332  bcs = subchk(array, hi, 0, -hshift, scaling)
1333  bce = subchk(array, hi, hshift, 0, scaling)
1334  if (sym) then
1335  bcw = subchk(array, hi, -hshift-1, 0, scaling)
1336  else
1337  bcw = subchk(array, hi, -hshift, 0, scaling)
1338  endif
1339  bcn = subchk(array, hi, 0, hshift, scaling)
1340 
1341  if (is_root_pe()) &
1342  call chk_sum_msg_nsew("u-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
1343  endif
1344 
1345  contains
1346 
1347  integer function subchk(array, HI, di, dj, scale)
1348  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1349  real, dimension(HI%IsdB:,HI%jsd:,:), intent(in) :: array !< The array to be checksummed
1350  integer, intent(in) :: di !< i- direction array shift for this checksum
1351  integer, intent(in) :: dj !< j- direction array shift for this checksum
1352  real, intent(in) :: scale !< A scaling factor for this array.
1353  integer :: i, j, k, bc
1354  subchk = 0
1355  ! This line deliberately uses the h-point computational domain.
1356  do k=lbound(array,3),ubound(array,3) ; do j=hi%jsc+dj,hi%jec+dj ; do i=hi%isc+di,hi%iec+di
1357  bc = bitcount(abs(scale*array(i,j,k)))
1358  subchk = subchk + bc
1359  enddo ; enddo ; enddo
1360  call sum_across_pes(subchk)
1361  subchk=mod(subchk, bc_modulus)
1362  end function subchk
1363 
1364  subroutine substats(HI, array, sym_stats, aMean, aMin, aMax)
1365  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1366  real, dimension(HI%IsdB:,HI%jsd:,:), intent(in) :: array !< The array to be checksummed
1367  logical, intent(in) :: sym_stats !< If true, evaluate the statistics on the
1368  !! full symmetric computational domain.
1369  real, intent(out) :: aMean, aMin, aMax
1370 
1371  integer :: i, j, k, n, IsB
1372 
1373  isb = hi%isc ; if (sym_stats) isb = hi%isc-1
1374 
1375  amin = array(hi%isc,hi%jsc,1) ; amax = amin
1376  do k=lbound(array,3),ubound(array,3) ; do j=hi%jsc,hi%jec ; do i=isb,hi%IecB
1377  amin = min(amin, array(i,j,k))
1378  amax = max(amax, array(i,j,k))
1379  enddo ; enddo ; enddo
1380  ! This line deliberately uses the h-point computational domain.
1381  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec,:))
1382  n = (1 + hi%jec - hi%jsc) * (1 + hi%iec - hi%isc) * size(array,3)
1383  call sum_across_pes(n)
1384  call min_across_pes(amin)
1385  call max_across_pes(amax)
1386  amean = amean / real(n)
1387  end subroutine substats
1388 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chk_sum_msg_w(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_uv_3d().

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

◆ chksum_uv_2d()

subroutine mom_checksums::chksum_uv_2d ( character(len=*), intent(in)  mesg,
real, dimension(hi%isdb:,hi%jsd:), intent(in)  arrayU,
real, dimension(hi%isd:,hi%jsdb:), intent(in)  arrayV,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a pair of 2d velocity arrays staggered at C-grid locations.

Parameters
[in]mesgIdentifying messages
[in]hiA horizontal index type
[in]arrayuThe u-component array to be checksummed
[in]arrayvThe v-component array to be checksummed
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for these arrays.
[in]logunitIO unit for checksum logging

Definition at line 589 of file MOM_checksums.F90.

589  character(len=*), intent(in) :: mesg !< Identifying messages
590  type(hor_index_type), intent(in) :: HI !< A horizontal index type
591  real, dimension(HI%IsdB:,HI%jsd:), intent(in) :: arrayU !< The u-component array to be checksummed
592  real, dimension(HI%isd:,HI%JsdB:), intent(in) :: arrayV !< The v-component array to be checksummed
593  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
594  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
595  !! symmetric computational domain.
596  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
597  real, optional, intent(in) :: scale !< A scaling factor for these arrays.
598  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
599 
600  if (present(haloshift)) then
601  call chksum_u_2d(arrayu, 'u '//mesg, hi, haloshift, symmetric, &
602  omit_corners, scale, logunit=logunit)
603  call chksum_v_2d(arrayv, 'v '//mesg, hi, haloshift, symmetric, &
604  omit_corners, scale, logunit=logunit)
605  else
606  call chksum_u_2d(arrayu, 'u '//mesg, hi, symmetric=symmetric, &
607  logunit=logunit)
608  call chksum_v_2d(arrayv, 'v '//mesg, hi, symmetric=symmetric, &
609  logunit=logunit)
610  endif
611 

References chksum_u_2d(), and chksum_v_2d().

Here is the call graph for this function:

◆ chksum_uv_3d()

subroutine mom_checksums::chksum_uv_3d ( character(len=*), intent(in)  mesg,
real, dimension(hi%isdb:,hi%jsd:,:), intent(in)  arrayU,
real, dimension(hi%isd:,hi%jsdb:,:), intent(in)  arrayV,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a pair of 3d velocity arrays staggered at C-grid locations.

Parameters
[in]mesgIdentifying messages
[in]hiA horizontal index type
[in]arrayuThe u-component array to be checksummed
[in]arrayvThe v-component array to be checksummed
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for these arrays.
[in]logunitIO unit for checksum logging

Definition at line 617 of file MOM_checksums.F90.

617  character(len=*), intent(in) :: mesg !< Identifying messages
618  type(hor_index_type), intent(in) :: HI !< A horizontal index type
619  real, dimension(HI%IsdB:,HI%jsd:,:), intent(in) :: arrayU !< The u-component array to be checksummed
620  real, dimension(HI%isd:,HI%JsdB:,:), intent(in) :: arrayV !< The v-component array to be checksummed
621  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
622  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
623  !! symmetric computational domain.
624  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
625  real, optional, intent(in) :: scale !< A scaling factor for these arrays.
626  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
627 
628  if (present(haloshift)) then
629  call chksum_u_3d(arrayu, 'u '//mesg, hi, haloshift, symmetric, &
630  omit_corners, scale, logunit=logunit)
631  call chksum_v_3d(arrayv, 'v '//mesg, hi, haloshift, symmetric, &
632  omit_corners, scale, logunit=logunit)
633  else
634  call chksum_u_3d(arrayu, 'u '//mesg, hi, symmetric=symmetric, &
635  logunit=logunit)
636  call chksum_v_3d(arrayv, 'v '//mesg, hi, symmetric=symmetric, &
637  logunit=logunit)
638  endif
639 

References chksum_u_3d(), and chksum_v_3d().

Here is the call graph for this function:

◆ chksum_v_2d()

subroutine mom_checksums::chksum_v_2d ( real, dimension(hi%isd:,hi%jsdb:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 2d array staggered at C-grid v points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 799 of file MOM_checksums.F90.

799  type(hor_index_type), intent(in) :: HI !< A horizontal index type
800  real, dimension(HI%isd:,HI%JsdB:), intent(in) :: array !< The array to be checksummed
801  character(len=*), intent(in) :: mesg !< An identifying message
802  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
803  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
804  !! symmetric computational domain.
805  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
806  real, optional, intent(in) :: scale !< A scaling factor for this array.
807  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
808 
809  real, allocatable, dimension(:,:) :: rescaled_array
810  real :: scaling
811  integer :: iounit !< Log IO unit
812  integer :: i, j, Js
813  real :: aMean, aMin, aMax
814  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
815  integer :: bcN, bcS, bcE, bcW
816  logical :: do_corners, sym, sym_stats
817 
818  if (checkfornans) then
819  if (is_nan(array(hi%isc:hi%iec,hi%JscB:hi%JecB))) &
820  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
821 ! if (is_NaN(array)) &
822 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
823  endif
824 
825  scaling = 1.0 ; if (present(scale)) scaling = scale
826  iounit = error_unit; if(present(logunit)) iounit = logunit
827  sym_stats = .false. ; if (present(symmetric)) sym_stats = symmetric
828  if (present(haloshift)) then ; if (haloshift > 0) sym_stats = .true. ; endif
829 
830  if (calculatestatistics) then
831  if (present(scale)) then
832  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
833  lbound(array,2):ubound(array,2)) )
834  rescaled_array(:,:) = 0.0
835  js = hi%jsc ; if (sym_stats) js = hi%jsc-1
836  do j=js,hi%JecB ; do i=hi%isc,hi%iec
837  rescaled_array(i,j) = scale*array(i,j)
838  enddo ; enddo
839  call substats(hi, rescaled_array, sym_stats, amean, amin, amax)
840  deallocate(rescaled_array)
841  else
842  call substats(hi, array, sym_stats, amean, amin, amax)
843  endif
844 
845  if (is_root_pe()) &
846  call chk_sum_msg("v-point:", amean, amin, amax, mesg, iounit)
847  endif
848 
849  if (.not.writechksums) return
850 
851  hshift = default_shift
852  if (present(haloshift)) hshift = haloshift
853  if (hshift<0) hshift = hi%ied-hi%iec
854 
855  if ( hi%isc-hshift<hi%isd .or. hi%iec+hshift>hi%ied .or. &
856  hi%jscB-hshift<hi%jsdB .or. hi%jecB+hshift>hi%jedB ) then
857  write(0,*) 'chksum_v_2d: haloshift =',hshift
858  write(0,*) 'chksum_v_2d: isd,isc,iec,ied=',hi%isd,hi%isc,hi%iec,hi%ied
859  write(0,*) 'chksum_v_2d: jsd,jsc,jec,jed=',hi%jsdB,hi%jscB,hi%jecB,hi%jedB
860  call chksum_error(fatal,'Error in chksum_v_2d '//trim(mesg))
861  endif
862 
863  bc0 = subchk(array, hi, 0, 0, scaling)
864 
865  sym = .false. ; if (present(symmetric)) sym = symmetric
866 
867  if ((hshift==0) .and. .not.sym) then
868  if (is_root_pe()) call chk_sum_msg("v-point:", bc0, mesg, iounit)
869  return
870  endif
871 
872  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
873 
874  if (hshift==0) then
875  bcs = subchk(array, hi, 0, -hshift-1, scaling)
876  if (is_root_pe()) call chk_sum_msg_s("v-point:", bc0, bcs, mesg, iounit)
877  elseif (do_corners) then
878  if (sym) then
879  bcsw = subchk(array, hi, -hshift, -hshift-1, scaling)
880  bcse = subchk(array, hi, hshift, -hshift-1, scaling)
881  else
882  bcsw = subchk(array, hi, -hshift, -hshift, scaling)
883  bcse = subchk(array, hi, hshift, -hshift, scaling)
884  endif
885  bcnw = subchk(array, hi, -hshift, hshift, scaling)
886  bcne = subchk(array, hi, hshift, hshift, scaling)
887 
888  if (is_root_pe()) &
889  call chk_sum_msg("v-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
890  else
891  if (sym) then
892  bcs = subchk(array, hi, 0, -hshift-1, scaling)
893  else
894  bcs = subchk(array, hi, 0, -hshift, scaling)
895  endif
896  bce = subchk(array, hi, hshift, 0, scaling)
897  bcw = subchk(array, hi, -hshift, 0, scaling)
898  bcn = subchk(array, hi, 0, hshift, scaling)
899 
900  if (is_root_pe()) &
901  call chk_sum_msg_nsew("v-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
902  endif
903 
904  contains
905 
906  integer function subchk(array, HI, di, dj, scale)
907  type(hor_index_type), intent(in) :: HI !< A horizontal index type
908  real, dimension(HI%isd:,HI%JsdB:), intent(in) :: array !< The array to be checksummed
909  integer, intent(in) :: di !< i- direction array shift for this checksum
910  integer, intent(in) :: dj !< j- direction array shift for this checksum
911  real, intent(in) :: scale !< A scaling factor for this array.
912  integer :: i, j, bc
913  subchk = 0
914  ! This line deliberately uses the h-point computational domain.
915  do j=hi%jsc+dj,hi%jec+dj; do i=hi%isc+di,hi%iec+di
916  bc = bitcount(abs(scale*array(i,j)))
917  subchk = subchk + bc
918  enddo ; enddo
919  call sum_across_pes(subchk)
920  subchk=mod(subchk, bc_modulus)
921  end function subchk
922 
923  subroutine substats(HI, array, sym_stats, aMean, aMin, aMax)
924  type(hor_index_type), intent(in) :: HI !< A horizontal index type
925  real, dimension(HI%isd:,HI%JsdB:), intent(in) :: array !< The array to be checksummed
926  logical, intent(in) :: sym_stats !< If true, evaluate the statistics on the
927  !! full symmetric computational domain.
928  real, intent(out) :: aMean, aMin, aMax
929 
930  integer :: i, j, n, JsB
931 
932  jsb = hi%jsc ; if (sym_stats) jsb = hi%jsc-1
933 
934  amin = array(hi%isc,hi%jsc) ; amax = amin
935  do j=jsb,hi%JecB ; do i=hi%isc,hi%iec
936  amin = min(amin, array(i,j))
937  amax = max(amax, array(i,j))
938  enddo ; enddo
939  ! This line deliberately uses the h-computational domain.
940  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec))
941  n = (1 + hi%jec - hi%jsc) * (1 + hi%iec - hi%isc)
942  call sum_across_pes(n)
943  call min_across_pes(amin)
944  call max_across_pes(amax)
945  amean = amean / real(n)
946  end subroutine substats
947 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chk_sum_msg_s(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_uv_2d().

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

◆ chksum_v_3d()

subroutine mom_checksums::chksum_v_3d ( real, dimension(hi%isd:,hi%jsdb:,:), intent(in)  array,
character(len=*), intent(in)  mesg,
type(hor_index_type), intent(in)  HI,
integer, intent(in), optional  haloshift,
logical, intent(in), optional  symmetric,
logical, intent(in), optional  omit_corners,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)
private

Checksums a 3d array staggered at C-grid v points.

Parameters
[in]hiA horizontal index type
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]haloshiftThe width of halos to check (default 0)
[in]symmetricIf true, do the checksums on the full symmetric computational domain.
[in]omit_cornersIf true, avoid checking diagonal shifts
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 1394 of file MOM_checksums.F90.

1394  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1395  real, dimension(HI%isd:,HI%JsdB:,:), intent(in) :: array !< The array to be checksummed
1396  character(len=*), intent(in) :: mesg !< An identifying message
1397  integer, optional, intent(in) :: haloshift !< The width of halos to check (default 0)
1398  logical, optional, intent(in) :: symmetric !< If true, do the checksums on the full
1399  !! symmetric computational domain.
1400  logical, optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
1401  real, optional, intent(in) :: scale !< A scaling factor for this array.
1402  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
1403 
1404  real, allocatable, dimension(:,:,:) :: rescaled_array
1405  real :: scaling
1406  integer :: iounit !< Log IO unit
1407  integer :: i, j, k, Js
1408  integer :: bc0, bcSW, bcSE, bcNW, bcNE, hshift
1409  integer :: bcN, bcS, bcE, bcW
1410  real :: aMean, aMin, aMax
1411  logical :: do_corners, sym, sym_stats
1412 
1413  if (checkfornans) then
1414  if (is_nan(array(hi%isc:hi%iec,hi%JscB:hi%JecB,:))) &
1415  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
1416 ! if (is_NaN(array)) &
1417 ! call chksum_error(FATAL, 'NaN detected in halo: '//trim(mesg))
1418  endif
1419 
1420  scaling = 1.0 ; if (present(scale)) scaling = scale
1421  iounit = error_unit; if(present(logunit)) iounit = logunit
1422  sym_stats = .false. ; if (present(symmetric)) sym_stats = symmetric
1423  if (present(haloshift)) then ; if (haloshift > 0) sym_stats = .true. ; endif
1424 
1425  if (calculatestatistics) then
1426  if (present(scale)) then
1427  allocate( rescaled_array(lbound(array,1):ubound(array,1), &
1428  lbound(array,2):ubound(array,2), &
1429  lbound(array,3):ubound(array,3)) )
1430  rescaled_array(:,:,:) = 0.0
1431  js = hi%jsc ; if (sym_stats) js = hi%jsc-1
1432  do k=1,size(array,3) ; do j=js,hi%JecB ; do i=hi%isc,hi%iec
1433  rescaled_array(i,j,k) = scale*array(i,j,k)
1434  enddo ; enddo ; enddo
1435  call substats(hi, rescaled_array, sym_stats, amean, amin, amax)
1436  deallocate(rescaled_array)
1437  else
1438  call substats(hi, array, sym_stats, amean, amin, amax)
1439  endif
1440  if (is_root_pe()) &
1441  call chk_sum_msg("v-point:", amean, amin, amax, mesg, iounit)
1442  endif
1443 
1444  if (.not.writechksums) return
1445 
1446  hshift = default_shift
1447  if (present(haloshift)) hshift = haloshift
1448  if (hshift<0) hshift = hi%ied-hi%iec
1449 
1450  if ( hi%isc-hshift<hi%isd .or. hi%iec+hshift>hi%ied .or. &
1451  hi%jsc-hshift<hi%jsd .or. hi%jec+hshift>hi%jed ) then
1452  write(0,*) 'chksum_v_3d: haloshift =',hshift
1453  write(0,*) 'chksum_v_3d: isd,isc,iec,ied=',hi%isd,hi%isc,hi%iec,hi%ied
1454  write(0,*) 'chksum_v_3d: jsd,jsc,jec,jed=',hi%jsd,hi%jsc,hi%jec,hi%jed
1455  call chksum_error(fatal,'Error in chksum_v_3d '//trim(mesg))
1456  endif
1457 
1458  bc0 = subchk(array, hi, 0, 0, scaling)
1459 
1460  sym = .false. ; if (present(symmetric)) sym = symmetric
1461 
1462  if ((hshift==0) .and. .not.sym) then
1463  if (is_root_pe()) call chk_sum_msg("v-point:", bc0, mesg, iounit)
1464  return
1465  endif
1466 
1467  do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners
1468 
1469  if (hshift==0) then
1470  bcs = subchk(array, hi, 0, -hshift-1, scaling)
1471  if (is_root_pe()) call chk_sum_msg_s("v-point:", bc0, bcs, mesg, iounit)
1472  elseif (do_corners) then
1473  if (sym) then
1474  bcsw = subchk(array, hi, -hshift, -hshift-1, scaling)
1475  bcse = subchk(array, hi, hshift, -hshift-1, scaling)
1476  else
1477  bcsw = subchk(array, hi, -hshift, -hshift, scaling)
1478  bcse = subchk(array, hi, hshift, -hshift, scaling)
1479  endif
1480  bcnw = subchk(array, hi, -hshift, hshift, scaling)
1481  bcne = subchk(array, hi, hshift, hshift, scaling)
1482 
1483  if (is_root_pe()) &
1484  call chk_sum_msg("v-point:", bc0, bcsw, bcse, bcnw, bcne, mesg, iounit)
1485  else
1486  if (sym) then
1487  bcs = subchk(array, hi, 0, -hshift-1, scaling)
1488  else
1489  bcs = subchk(array, hi, 0, -hshift, scaling)
1490  endif
1491  bce = subchk(array, hi, hshift, 0, scaling)
1492  bcw = subchk(array, hi, -hshift, 0, scaling)
1493  bcn = subchk(array, hi, 0, hshift, scaling)
1494 
1495  if (is_root_pe()) &
1496  call chk_sum_msg_nsew("v-point:", bc0, bcn, bcs, bce, bcw, mesg, iounit)
1497  endif
1498 
1499  contains
1500 
1501  integer function subchk(array, HI, di, dj, scale)
1502  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1503  real, dimension(HI%isd:,HI%JsdB:,:), intent(in) :: array !< The array to be checksummed
1504  integer, intent(in) :: di !< i- direction array shift for this checksum
1505  integer, intent(in) :: dj !< j- direction array shift for this checksum
1506  real, intent(in) :: scale !< A scaling factor for this array.
1507  integer :: i, j, k, bc
1508  subchk = 0
1509  ! This line deliberately uses the h-point computational domain.
1510  do k=lbound(array,3),ubound(array,3) ; do j=hi%jsc+dj,hi%jec+dj ; do i=hi%isc+di,hi%iec+di
1511  bc = bitcount(abs(scale*array(i,j,k)))
1512  subchk = subchk + bc
1513  enddo ; enddo ; enddo
1514  call sum_across_pes(subchk)
1515  subchk=mod(subchk, bc_modulus)
1516  end function subchk
1517 
1518  !subroutine subStats(HI, array, mesg, sym_stats)
1519  subroutine substats(HI, array, sym_stats, aMean, aMin, aMax)
1520  type(hor_index_type), intent(in) :: HI !< A horizontal index type
1521  real, dimension(HI%isd:,HI%JsdB:,:), intent(in) :: array !< The array to be checksummed
1522  logical, intent(in) :: sym_stats !< If true, evaluate the statistics on the
1523  !! full symmetric computational domain.
1524  real, intent(out) :: aMean, aMin, aMax !< Mean/min/max of array over domain
1525 
1526  integer :: i, j, k, n, JsB
1527 
1528  jsb = hi%jsc ; if (sym_stats) jsb = hi%jsc-1
1529 
1530  amin = array(hi%isc,hi%jsc,1) ; amax = amin
1531  do k=lbound(array,3),ubound(array,3) ; do j=jsb,hi%JecB ; do i=hi%isc,hi%iec
1532  amin = min(amin, array(i,j,k))
1533  amax = max(amax, array(i,j,k))
1534  enddo ; enddo ; enddo
1535  ! This line deliberately uses the h-point computational domain.
1536  amean = reproducing_sum(array(hi%isc:hi%iec,hi%jsc:hi%jec,:))
1537  n = (1 + hi%jec - hi%jsc) * (1 + hi%iec - hi%isc) * size(array,3)
1538  call sum_across_pes(n)
1539  call min_across_pes(amin)
1540  call max_across_pes(amax)
1541  amean = amean / real(n)
1542  end subroutine substats
1543 

References calculatestatistics, checkfornans, chk_sum_msg_nsew(), chk_sum_msg_s(), chksum_error(), default_shift, mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Referenced by chksum_uv_3d().

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

◆ is_nan_0d()

logical function mom_checksums::is_nan_0d ( real, intent(in)  x)
private

This function returns .true. if x is a NaN, and .false. otherwise.

Parameters
[in]xThe value to be checked for NaNs.

Definition at line 1666 of file MOM_checksums.F90.

1666  real, intent(in) :: x !< The value to be checked for NaNs.
1667  logical :: is_NaN_0d
1668 
1669  !is_NaN_0d = (((x < 0.0) .and. (x >= 0.0)) .or. &
1670  ! (.not.(x < 0.0) .and. .not.(x >= 0.0)))
1671  if (((x < 0.0) .and. (x >= 0.0)) .or. &
1672  (.not.(x < 0.0) .and. .not.(x >= 0.0))) then
1673  is_nan_0d = .true.
1674  else
1675  is_nan_0d = .false.
1676  endif
1677 

Referenced by is_nan_1d(), is_nan_2d(), and is_nan_3d().

Here is the caller graph for this function:

◆ is_nan_1d()

logical function mom_checksums::is_nan_1d ( real, dimension(:), intent(in)  x,
logical, intent(in), optional  skip_mpp 
)
private

Returns .true. if any element of x is a NaN, and .false. otherwise.

Parameters
[in]xThe array to be checked for NaNs.
[in]skip_mppIf true, only check this array only on the local PE (default false).

Definition at line 1682 of file MOM_checksums.F90.

1682  real, dimension(:), intent(in) :: x !< The array to be checked for NaNs.
1683  logical, optional, intent(in) :: skip_mpp !< If true, only check this array only
1684  !! on the local PE (default false).
1685  logical :: is_NaN_1d
1686 
1687  integer :: i, n
1688  logical :: call_mpp
1689 
1690  n = 0
1691  do i = lbound(x,1), ubound(x,1)
1692  if (is_nan_0d(x(i))) n = n + 1
1693  enddo
1694  call_mpp = .true.
1695  if (present(skip_mpp)) call_mpp = .not.skip_mpp
1696 
1697  if (call_mpp) call sum_across_pes(n)
1698  is_nan_1d = .false.
1699  if (n>0) is_nan_1d = .true.
1700 

References is_nan_0d().

Here is the call graph for this function:

◆ is_nan_2d()

logical function mom_checksums::is_nan_2d ( real, dimension(:,:), intent(in)  x)
private

Returns .true. if any element of x is a NaN, and .false. otherwise.

Parameters
[in]xThe array to be checked for NaNs.

Definition at line 1705 of file MOM_checksums.F90.

1705  real, dimension(:,:), intent(in) :: x !< The array to be checked for NaNs.
1706  logical :: is_NaN_2d
1707 
1708  integer :: i, j, n
1709 
1710  n = 0
1711  do j = lbound(x,2), ubound(x,2) ; do i = lbound(x,1), ubound(x,1)
1712  if (is_nan_0d(x(i,j))) n = n + 1
1713  enddo ; enddo
1714  call sum_across_pes(n)
1715  is_nan_2d = .false.
1716  if (n>0) is_nan_2d = .true.
1717 

References is_nan_0d().

Here is the call graph for this function:

◆ is_nan_3d()

logical function mom_checksums::is_nan_3d ( real, dimension(:,:,:), intent(in)  x)
private

Returns .true. if any element of x is a NaN, and .false. otherwise.

Parameters
[in]xThe array to be checked for NaNs.

Definition at line 1722 of file MOM_checksums.F90.

1722  real, dimension(:,:,:), intent(in) :: x !< The array to be checked for NaNs.
1723  logical :: is_NaN_3d
1724 
1725  integer :: i, j, k, n
1726 
1727  n = 0
1728  do k = lbound(x,3), ubound(x,3)
1729  do j = lbound(x,2), ubound(x,2) ; do i = lbound(x,1), ubound(x,1)
1730  if (is_nan_0d(x(i,j,k))) n = n + 1
1731  enddo ; enddo
1732  enddo
1733  call sum_across_pes(n)
1734  is_nan_3d = .false.
1735  if (n>0) is_nan_3d = .true.
1736 

References is_nan_0d().

Here is the call graph for this function:

◆ mom_checksums_init()

subroutine, public mom_checksums::mom_checksums_init ( type(param_file_type), intent(in)  param_file)

MOM_checksums_init initializes the MOM_checksums module. As it happens, the only thing that it does is to log the version of this module.

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

Definition at line 1835 of file MOM_checksums.F90.

1835  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
1836 ! This include declares and sets the variable "version".
1837 #include "version_variable.h"
1838  character(len=40) :: mdl = "MOM_checksums" ! This module's name.
1839 
1840  call log_version(param_file, mdl, version)
1841 

Referenced by mom_debugging::mom_debugging_init().

Here is the caller graph for this function:

◆ zchksum()

subroutine, public mom_checksums::zchksum ( real, dimension(:), intent(in)  array,
character(len=*), intent(in)  mesg,
real, intent(in), optional  scale,
integer, intent(in), optional  logunit 
)

Checksum a 1d array (typically a column).

Parameters
[in]arrayThe array to be checksummed
[in]mesgAn identifying message
[in]scaleA scaling factor for this array.
[in]logunitIO unit for checksum logging

Definition at line 121 of file MOM_checksums.F90.

121  real, dimension(:), intent(in) :: array !< The array to be checksummed
122  character(len=*), intent(in) :: mesg !< An identifying message
123  real, optional, intent(in) :: scale !< A scaling factor for this array.
124  integer, optional, intent(in) :: logunit !< IO unit for checksum logging
125 
126  real, allocatable, dimension(:) :: rescaled_array
127  real :: scaling
128  integer :: iounit !< Log IO unit
129  integer :: k
130  real :: aMean, aMin, aMax
131  integer :: bc0
132 
133  if (checkfornans) then
134  if (is_nan(array(:))) &
135  call chksum_error(fatal, 'NaN detected: '//trim(mesg))
136  endif
137 
138  scaling = 1.0 ; if (present(scale)) scaling = scale
139  iounit = error_unit; if(present(logunit)) iounit = logunit
140 
141  if (calculatestatistics) then
142  if (present(scale)) then
143  allocate(rescaled_array(lbound(array,1):ubound(array,1)))
144  rescaled_array(:) = 0.0
145  do k=1, size(array, 1)
146  rescaled_array(k) = scale * array(k)
147  enddo
148 
149  call substats(rescaled_array, amean, amin, amax)
150  deallocate(rescaled_array)
151  else
152  call substats(array, amean, amin, amax)
153  endif
154 
155  if (is_root_pe()) &
156  call chk_sum_msg(" column:", amean, amin, amax, mesg, iounit)
157  endif
158 
159  if (.not. writechksums) return
160 
161  bc0 = subchk(array, scaling)
162  if (is_root_pe()) call chk_sum_msg(" column:", bc0, mesg, iounit)
163 
164  contains
165 
166  integer function subchk(array, scale)
167  real, dimension(:), intent(in) :: array !< The array to be checksummed
168  real, intent(in) :: scale !< A scaling factor for this array.
169  integer :: k, bc
170  subchk = 0
171  do k=lbound(array, 1), ubound(array, 1)
172  bc = bitcount(abs(scale * array(k)))
173  subchk = subchk + bc
174  enddo
175  subchk=mod(subchk, bc_modulus)
176  end function subchk
177 
178  subroutine substats(array, aMean, aMin, aMax)
179  real, dimension(:), intent(in) :: array !< The array to be checksummed
180  real, intent(out) :: aMean, aMin, aMax
181 
182  integer :: k, n
183 
184  amin = array(1)
185  amax = array(1)
186  n = 0
187  do k=lbound(array,1), ubound(array,1)
188  amin = min(amin, array(k))
189  amax = max(amax, array(k))
190  n = n + 1
191  enddo
192  amean = sum(array(:)) / real(n)
193  end subroutine substats
194 

References calculatestatistics, checkfornans, chksum_error(), mom_error_handler::is_root_pe(), subchk(), substats(), and writechksums.

Here is the call graph for this function:

Variable Documentation

◆ bc_modulus

integer, parameter mom_checksums::bc_modulus = 1000000000
private

Modulus of checksum bitcount.

Definition at line 77 of file MOM_checksums.F90.

77 integer, parameter :: bc_modulus = 1000000000 !< Modulus of checksum bitcount

Referenced by chksum0(), and subchk().

◆ calculatestatistics

logical mom_checksums::calculatestatistics =.true.
private

If true, report min, max and mean.

Definition at line 79 of file MOM_checksums.F90.

79 logical :: calculateStatistics=.true. !< If true, report min, max and mean.

Referenced by chksum0(), chksum_b_2d(), chksum_b_3d(), chksum_h_2d(), chksum_h_3d(), chksum_u_2d(), chksum_u_3d(), chksum_v_2d(), chksum_v_3d(), and zchksum().

◆ checkfornans

logical mom_checksums::checkfornans =.true.
private

If true, checks array for NaNs and cause FATAL error is any are found.

Definition at line 81 of file MOM_checksums.F90.

81 logical :: checkForNaNs=.true. !< If true, checks array for NaNs and cause

Referenced by chksum0(), chksum_b_2d(), chksum_b_3d(), chksum_h_2d(), chksum_h_3d(), chksum_u_2d(), chksum_u_3d(), chksum_v_2d(), chksum_v_3d(), and zchksum().

◆ default_shift

integer, parameter mom_checksums::default_shift =0
private

The default array shift.

Definition at line 78 of file MOM_checksums.F90.

78 integer, parameter :: default_shift=0 !< The default array shift

Referenced by chksum_b_2d(), chksum_b_3d(), chksum_h_2d(), chksum_h_3d(), chksum_u_2d(), chksum_u_3d(), chksum_v_2d(), and chksum_v_3d().

◆ writechksums

logical mom_checksums::writechksums =.true.
private

If true, report the bitcount checksum.

Definition at line 80 of file MOM_checksums.F90.

80 logical :: writeChksums=.true. !< If true, report the bitcount checksum

Referenced by chksum0(), chksum_b_2d(), chksum_b_3d(), chksum_h_2d(), chksum_h_3d(), chksum_u_2d(), chksum_u_3d(), chksum_v_2d(), chksum_v_3d(), and zchksum().

substats
subroutine substats(array, aMean, aMin, aMax)
Definition: MOM_checksums.F90:179
subchk
integer function subchk(array, scale)
Definition: MOM_checksums.F90:167