diff --git a/src/framework/MOM_checksums.F90 b/src/framework/MOM_checksums.F90 index 4e349e0fb7..7f0dfd5d9b 100644 --- a/src/framework/MOM_checksums.F90 +++ b/src/framework/MOM_checksums.F90 @@ -11,6 +11,7 @@ module MOM_checksums use MOM_error_handler, only : MOM_error, FATAL, is_root_pe use MOM_file_parser, only : log_version, param_file_type use MOM_hor_index, only : hor_index_type, rotate_hor_index +use MOM_murmur_hash, only : murmur_hash use iso_fortran_env, only : error_unit, int32, int64 @@ -97,7 +98,9 @@ module MOM_checksums logical :: calculateStatistics=.true. !< If true, report min, max and mean. logical :: writeChksums=.true. !< If true, report the bitcount checksum logical :: checkForNaNs=.true. !< If true, checks array for NaNs and cause - !! FATAL error is any are found + !! FATAL error if any are found +logical :: writeHash = .false. !< If true, report the murmur hash + !! NOTE: writeHash is currently disabled due to non-compliant diagnostics. contains @@ -144,6 +147,9 @@ subroutine chksum0(scalar, mesg, scale, logunit, unscale) if (is_root_pe()) & call chk_sum_msg(" scalar:", bc, mesg, iounit) + if (writeHash .and. is_root_pe()) & + write(iounit, '(" scalar: hash=", z8, 1x, a)') & + murmur_hash(scaling * scalar), mesg end subroutine chksum0 @@ -204,6 +210,10 @@ subroutine zchksum(array, mesg, scale, logunit, unscale) bc0 = subchk(array, scaling) if (is_root_pe()) call chk_sum_msg(" column:", bc0, mesg, iounit) + if (writeHash .and. is_root_pe()) & + write(iounit, '(" column: hash=", z8, 1x, a)') & + murmur_hash(scaling * array), mesg + contains integer function subchk(array, unscale) @@ -381,6 +391,7 @@ subroutine chksum_h_2d(array_m, mesg, HI_m, haloshift, omit_corners, scale, logu ! for checksums and output real, pointer :: array(:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -391,6 +402,7 @@ subroutine chksum_h_2d(array_m, mesg, HI_m, haloshift, omit_corners, scale, logu logical :: do_corners integer :: turns ! Quarter turns from input to model grid + ! Rotate array to the input grid turns = HI_m%turns if (modulo(turns, 4) /= 0) then @@ -451,27 +463,36 @@ subroutine chksum_h_2d(array_m, mesg, HI_m, haloshift, omit_corners, scale, logu if (hshift==0) then if (is_root_pe()) call chk_sum_msg("h-point:", bc0, mesg, iounit) - return - endif + else + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners + + if (do_corners) then + bcSW = subchk(array, HI, -hshift, -hshift, scaling) + bcSE = subchk(array, HI, hshift, -hshift, scaling) + bcNW = subchk(array, HI, -hshift, hshift, scaling) + bcNE = subchk(array, HI, hshift, hshift, scaling) - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners + if (is_root_pe()) & + call chk_sum_msg("h-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) + else + bcS = subchk(array, HI, 0, -hshift, scaling) + bcE = subchk(array, HI, hshift, 0, scaling) + bcW = subchk(array, HI, -hshift, 0, scaling) + bcN = subchk(array, HI, 0, hshift, scaling) - if (do_corners) then - bcSW = subchk(array, HI, -hshift, -hshift, scaling) - bcSE = subchk(array, HI, hshift, -hshift, scaling) - bcNW = subchk(array, HI, -hshift, hshift, scaling) - bcNE = subchk(array, HI, hshift, hshift, scaling) + if (is_root_pe()) & + call chk_sum_msg_NSEW("h-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + endif + endif - if (is_root_pe()) & - call chk_sum_msg("h-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) - else - bcS = subchk(array, HI, 0, -hshift, scaling) - bcE = subchk(array, HI, hshift, 0, scaling) - bcW = subchk(array, HI, -hshift, 0, scaling) - bcN = subchk(array, HI, 0, hshift, scaling) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec)) + hash_array(:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec) - if (is_root_pe()) & - call chk_sum_msg_NSEW("h-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + write(iounit, '("h-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains @@ -675,6 +696,7 @@ subroutine chksum_B_2d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, ! for checksums and output real, pointer :: array(:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -750,33 +772,42 @@ subroutine chksum_B_2d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, if ((hshift==0) .and. .not.sym) then if (is_root_pe()) call chk_sum_msg("B-point:", bc0, mesg, iounit) - return - endif - - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners - - if (do_corners) then - if (sym) then - bcSW = subchk(array, HI, -hshift-1, -hshift-1, scaling) - bcSE = subchk(array, HI, hshift, -hshift-1, scaling) - bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + else + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners + + if (do_corners) then + if (sym) then + bcSW = subchk(array, HI, -hshift-1, -hshift-1, scaling) + bcSE = subchk(array, HI, hshift, -hshift-1, scaling) + bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + else + bcSW = subchk(array, HI, -hshift, -hshift, scaling) + bcSE = subchk(array, HI, hshift, -hshift, scaling) + bcNW = subchk(array, HI, -hshift, hshift, scaling) + endif + bcNE = subchk(array, HI, hshift, hshift, scaling) + + if (is_root_pe()) & + call chk_sum_msg("B-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) else - bcSW = subchk(array, HI, -hshift, -hshift, scaling) - bcSE = subchk(array, HI, hshift, -hshift, scaling) - bcNW = subchk(array, HI, -hshift, hshift, scaling) + bcS = subchk(array, HI, 0, -hshift, scaling) + bcE = subchk(array, HI, hshift, 0, scaling) + bcW = subchk(array, HI, -hshift, 0, scaling) + bcN = subchk(array, HI, 0, hshift, scaling) + + if (is_root_pe()) & + call chk_sum_msg_NSEW("B-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) endif - bcNE = subchk(array, HI, hshift, hshift, scaling) + endif - if (is_root_pe()) & - call chk_sum_msg("B-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) - else - bcS = subchk(array, HI, 0, -hshift, scaling) - bcE = subchk(array, HI, hshift, 0, scaling) - bcW = subchk(array, HI, -hshift, 0, scaling) - bcN = subchk(array, HI, 0, hshift, scaling) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec)) + hash_array(:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec) - if (is_root_pe()) & - call chk_sum_msg_NSEW("B-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + write(iounit, '("B-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains @@ -981,6 +1012,7 @@ subroutine chksum_u_2d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, ! for checksums and output real, pointer :: array(:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -1065,39 +1097,48 @@ subroutine chksum_u_2d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, if ((hshift==0) .and. .not.sym) then if (is_root_pe()) call chk_sum_msg("u-point:", bc0, mesg, iounit) - return - endif + else + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners + if (hshift==0) then + bcW = subchk(array, HI, -hshift-1, 0, scaling) + if (is_root_pe()) call chk_sum_msg_W("u-point:", bc0, bcW, mesg, iounit) + elseif (do_corners) then + if (sym) then + bcSW = subchk(array, HI, -hshift-1, -hshift, scaling) + bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + else + bcSW = subchk(array, HI, -hshift, -hshift, scaling) + bcNW = subchk(array, HI, -hshift, hshift, scaling) + endif + bcSE = subchk(array, HI, hshift, -hshift, scaling) + bcNE = subchk(array, HI, hshift, hshift, scaling) - if (hshift==0) then - bcW = subchk(array, HI, -hshift-1, 0, scaling) - if (is_root_pe()) call chk_sum_msg_W("u-point:", bc0, bcW, mesg, iounit) - elseif (do_corners) then - if (sym) then - bcSW = subchk(array, HI, -hshift-1, -hshift, scaling) - bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + if (is_root_pe()) & + call chk_sum_msg("u-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) else - bcSW = subchk(array, HI, -hshift, -hshift, scaling) - bcNW = subchk(array, HI, -hshift, hshift, scaling) + bcS = subchk(array, HI, 0, -hshift, scaling) + bcE = subchk(array, HI, hshift, 0, scaling) + if (sym) then + bcW = subchk(array, HI, -hshift-1, 0, scaling) + else + bcW = subchk(array, HI, -hshift, 0, scaling) + endif + bcN = subchk(array, HI, 0, hshift, scaling) + + if (is_root_pe()) & + call chk_sum_msg_NSEW("u-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) endif - bcSE = subchk(array, HI, hshift, -hshift, scaling) - bcNE = subchk(array, HI, hshift, hshift, scaling) + endif - if (is_root_pe()) & - call chk_sum_msg("u-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) - else - bcS = subchk(array, HI, 0, -hshift, scaling) - bcE = subchk(array, HI, hshift, 0, scaling) - if (sym) then - bcW = subchk(array, HI, -hshift-1, 0, scaling) - else - bcW = subchk(array, HI, -hshift, 0, scaling) - endif - bcN = subchk(array, HI, 0, hshift, scaling) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec)) + hash_array(:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec) - if (is_root_pe()) & - call chk_sum_msg_NSEW("u-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + write(iounit, '("u-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains @@ -1175,6 +1216,7 @@ subroutine chksum_v_2d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, ! for checksums and output real, pointer :: array(:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -1259,39 +1301,48 @@ subroutine chksum_v_2d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, if ((hshift==0) .and. .not.sym) then if (is_root_pe()) call chk_sum_msg("v-point:", bc0, mesg, iounit) - return - endif + else + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners + if (hshift==0) then + bcS = subchk(array, HI, 0, -hshift-1, scaling) + if (is_root_pe()) call chk_sum_msg_S("v-point:", bc0, bcS, mesg, iounit) + elseif (do_corners) then + if (sym) then + bcSW = subchk(array, HI, -hshift, -hshift-1, scaling) + bcSE = subchk(array, HI, hshift, -hshift-1, scaling) + else + bcSW = subchk(array, HI, -hshift, -hshift, scaling) + bcSE = subchk(array, HI, hshift, -hshift, scaling) + endif + bcNW = subchk(array, HI, -hshift, hshift, scaling) + bcNE = subchk(array, HI, hshift, hshift, scaling) - if (hshift==0) then - bcS = subchk(array, HI, 0, -hshift-1, scaling) - if (is_root_pe()) call chk_sum_msg_S("v-point:", bc0, bcS, mesg, iounit) - elseif (do_corners) then - if (sym) then - bcSW = subchk(array, HI, -hshift, -hshift-1, scaling) - bcSE = subchk(array, HI, hshift, -hshift-1, scaling) + if (is_root_pe()) & + call chk_sum_msg("v-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) else - bcSW = subchk(array, HI, -hshift, -hshift, scaling) - bcSE = subchk(array, HI, hshift, -hshift, scaling) - endif - bcNW = subchk(array, HI, -hshift, hshift, scaling) - bcNE = subchk(array, HI, hshift, hshift, scaling) + if (sym) then + bcS = subchk(array, HI, 0, -hshift-1, scaling) + else + bcS = subchk(array, HI, 0, -hshift, scaling) + endif + bcE = subchk(array, HI, hshift, 0, scaling) + bcW = subchk(array, HI, -hshift, 0, scaling) + bcN = subchk(array, HI, 0, hshift, scaling) - if (is_root_pe()) & - call chk_sum_msg("v-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) - else - if (sym) then - bcS = subchk(array, HI, 0, -hshift-1, scaling) - else - bcS = subchk(array, HI, 0, -hshift, scaling) + if (is_root_pe()) & + call chk_sum_msg_NSEW("v-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) endif - bcE = subchk(array, HI, hshift, 0, scaling) - bcW = subchk(array, HI, -hshift, 0, scaling) - bcN = subchk(array, HI, 0, hshift, scaling) + endif - if (is_root_pe()) & - call chk_sum_msg_NSEW("v-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec)) + hash_array(:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec) + + write(iounit, '("v-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains @@ -1366,6 +1417,7 @@ subroutine chksum_h_3d(array_m, mesg, HI_m, haloshift, omit_corners, scale, logu ! for checksums and output real, pointer :: array(:,:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -1438,27 +1490,36 @@ subroutine chksum_h_3d(array_m, mesg, HI_m, haloshift, omit_corners, scale, logu if (hshift==0) then if (is_root_pe()) call chk_sum_msg("h-point:", bc0, mesg, iounit) - return - endif + else + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners + + if (do_corners) then + bcSW = subchk(array, HI, -hshift, -hshift, scaling) + bcSE = subchk(array, HI, hshift, -hshift, scaling) + bcNW = subchk(array, HI, -hshift, hshift, scaling) + bcNE = subchk(array, HI, hshift, hshift, scaling) - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners + if (is_root_pe()) & + call chk_sum_msg("h-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) + else + bcS = subchk(array, HI, 0, -hshift, scaling) + bcE = subchk(array, HI, hshift, 0, scaling) + bcW = subchk(array, HI, -hshift, 0, scaling) + bcN = subchk(array, HI, 0, hshift, scaling) - if (do_corners) then - bcSW = subchk(array, HI, -hshift, -hshift, scaling) - bcSE = subchk(array, HI, hshift, -hshift, scaling) - bcNW = subchk(array, HI, -hshift, hshift, scaling) - bcNE = subchk(array, HI, hshift, hshift, scaling) + if (is_root_pe()) & + call chk_sum_msg_NSEW("h-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + endif + endif - if (is_root_pe()) & - call chk_sum_msg("h-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) - else - bcS = subchk(array, HI, 0, -hshift, scaling) - bcE = subchk(array, HI, hshift, 0, scaling) - bcW = subchk(array, HI, -hshift, 0, scaling) - bcN = subchk(array, HI, 0, hshift, scaling) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec, size(array, 3))) + hash_array(:,:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec, :) - if (is_root_pe()) & - call chk_sum_msg_NSEW("h-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + write(iounit, '("h-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains @@ -1532,6 +1593,7 @@ subroutine chksum_B_3d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, ! for checksums and output real, pointer :: array(:,:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -1609,38 +1671,47 @@ subroutine chksum_B_3d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, if ((hshift==0) .and. .not.sym) then if (is_root_pe()) call chk_sum_msg("B-point:", bc0, mesg, iounit) - return - endif - - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners - - if (do_corners) then - if (sym) then - bcSW = subchk(array, HI, -hshift-1, -hshift-1, scaling) - bcSE = subchk(array, HI, hshift, -hshift-1, scaling) - bcNW = subchk(array, HI, -hshift-1, hshift, scaling) - else - bcSW = subchk(array, HI, -hshift-1, -hshift-1, scaling) - bcSE = subchk(array, HI, hshift, -hshift-1, scaling) - bcNW = subchk(array, HI, -hshift-1, hshift, scaling) - endif - bcNE = subchk(array, HI, hshift, hshift, scaling) - - if (is_root_pe()) & - call chk_sum_msg("B-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) else - if (sym) then - bcS = subchk(array, HI, 0, -hshift-1, scaling) - bcW = subchk(array, HI, -hshift-1, 0, scaling) + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners + + if (do_corners) then + if (sym) then + bcSW = subchk(array, HI, -hshift-1, -hshift-1, scaling) + bcSE = subchk(array, HI, hshift, -hshift-1, scaling) + bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + else + bcSW = subchk(array, HI, -hshift-1, -hshift-1, scaling) + bcSE = subchk(array, HI, hshift, -hshift-1, scaling) + bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + endif + bcNE = subchk(array, HI, hshift, hshift, scaling) + + if (is_root_pe()) & + call chk_sum_msg("B-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) else - bcS = subchk(array, HI, 0, -hshift, scaling) - bcW = subchk(array, HI, -hshift, 0, scaling) + if (sym) then + bcS = subchk(array, HI, 0, -hshift-1, scaling) + bcW = subchk(array, HI, -hshift-1, 0, scaling) + else + bcS = subchk(array, HI, 0, -hshift, scaling) + bcW = subchk(array, HI, -hshift, 0, scaling) + endif + bcE = subchk(array, HI, hshift, 0, scaling) + bcN = subchk(array, HI, 0, hshift, scaling) + + if (is_root_pe()) & + call chk_sum_msg_NSEW("B-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) endif - bcE = subchk(array, HI, hshift, 0, scaling) - bcN = subchk(array, HI, 0, hshift, scaling) + endif - if (is_root_pe()) & - call chk_sum_msg_NSEW("B-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec, size(array, 3))) + hash_array(:,:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec, :) + + write(iounit, '("B-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains @@ -1718,6 +1789,7 @@ subroutine chksum_u_3d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, ! for checksums and output real, pointer :: array(:,:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -1802,39 +1874,48 @@ subroutine chksum_u_3d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, if ((hshift==0) .and. .not.sym) then if (is_root_pe()) call chk_sum_msg("u-point:", bc0, mesg, iounit) - return - endif + else + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners + if (hshift==0) then + bcW = subchk(array, HI, -hshift-1, 0, scaling) + if (is_root_pe()) call chk_sum_msg_W("u-point:", bc0, bcW, mesg, iounit) + elseif (do_corners) then + if (sym) then + bcSW = subchk(array, HI, -hshift-1, -hshift, scaling) + bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + else + bcSW = subchk(array, HI, -hshift, -hshift, scaling) + bcNW = subchk(array, HI, -hshift, hshift, scaling) + endif + bcSE = subchk(array, HI, hshift, -hshift, scaling) + bcNE = subchk(array, HI, hshift, hshift, scaling) - if (hshift==0) then - bcW = subchk(array, HI, -hshift-1, 0, scaling) - if (is_root_pe()) call chk_sum_msg_W("u-point:", bc0, bcW, mesg, iounit) - elseif (do_corners) then - if (sym) then - bcSW = subchk(array, HI, -hshift-1, -hshift, scaling) - bcNW = subchk(array, HI, -hshift-1, hshift, scaling) + if (is_root_pe()) & + call chk_sum_msg("u-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) else - bcSW = subchk(array, HI, -hshift, -hshift, scaling) - bcNW = subchk(array, HI, -hshift, hshift, scaling) + bcS = subchk(array, HI, 0, -hshift, scaling) + bcE = subchk(array, HI, hshift, 0, scaling) + if (sym) then + bcW = subchk(array, HI, -hshift-1, 0, scaling) + else + bcW = subchk(array, HI, -hshift, 0, scaling) + endif + bcN = subchk(array, HI, 0, hshift, scaling) + + if (is_root_pe()) & + call chk_sum_msg_NSEW("u-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) endif - bcSE = subchk(array, HI, hshift, -hshift, scaling) - bcNE = subchk(array, HI, hshift, hshift, scaling) + endif - if (is_root_pe()) & - call chk_sum_msg("u-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) - else - bcS = subchk(array, HI, 0, -hshift, scaling) - bcE = subchk(array, HI, hshift, 0, scaling) - if (sym) then - bcW = subchk(array, HI, -hshift-1, 0, scaling) - else - bcW = subchk(array, HI, -hshift, 0, scaling) - endif - bcN = subchk(array, HI, 0, hshift, scaling) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec, size(array, 3))) + hash_array(:,:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec, :) - if (is_root_pe()) & - call chk_sum_msg_NSEW("u-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + write(iounit, '("u-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains @@ -1912,6 +1993,7 @@ subroutine chksum_v_3d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, ! for checksums and output real, pointer :: array(:,:,:) ! Field array on the input grid [A ~> a] real, allocatable, dimension(:,:,:) :: rescaled_array ! The array with scaling undone [a] + real, allocatable :: hash_array(:,:,:) ! Subarray used to compute hash [a] type(hor_index_type), pointer :: HI ! Horizontal index bounds of the input grid real :: scaling ! Explicit rescaling factor [a A-1 ~> 1] integer :: iounit !< Log IO unit @@ -1996,39 +2078,48 @@ subroutine chksum_v_3d(array_m, mesg, HI_m, haloshift, symmetric, omit_corners, if ((hshift==0) .and. .not.sym) then if (is_root_pe()) call chk_sum_msg("v-point:", bc0, mesg, iounit) - return - endif + else + do_corners = .true. + if (present(omit_corners)) do_corners = .not. omit_corners - do_corners = .true. ; if (present(omit_corners)) do_corners = .not.omit_corners + if (hshift==0) then + bcS = subchk(array, HI, 0, -hshift-1, scaling) + if (is_root_pe()) call chk_sum_msg_S("v-point:", bc0, bcS, mesg, iounit) + elseif (do_corners) then + if (sym) then + bcSW = subchk(array, HI, -hshift, -hshift-1, scaling) + bcSE = subchk(array, HI, hshift, -hshift-1, scaling) + else + bcSW = subchk(array, HI, -hshift, -hshift, scaling) + bcSE = subchk(array, HI, hshift, -hshift, scaling) + endif + bcNW = subchk(array, HI, -hshift, hshift, scaling) + bcNE = subchk(array, HI, hshift, hshift, scaling) - if (hshift==0) then - bcS = subchk(array, HI, 0, -hshift-1, scaling) - if (is_root_pe()) call chk_sum_msg_S("v-point:", bc0, bcS, mesg, iounit) - elseif (do_corners) then - if (sym) then - bcSW = subchk(array, HI, -hshift, -hshift-1, scaling) - bcSE = subchk(array, HI, hshift, -hshift-1, scaling) + if (is_root_pe()) & + call chk_sum_msg("v-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) else - bcSW = subchk(array, HI, -hshift, -hshift, scaling) - bcSE = subchk(array, HI, hshift, -hshift, scaling) - endif - bcNW = subchk(array, HI, -hshift, hshift, scaling) - bcNE = subchk(array, HI, hshift, hshift, scaling) + if (sym) then + bcS = subchk(array, HI, 0, -hshift-1, scaling) + else + bcS = subchk(array, HI, 0, -hshift, scaling) + endif + bcE = subchk(array, HI, hshift, 0, scaling) + bcW = subchk(array, HI, -hshift, 0, scaling) + bcN = subchk(array, HI, 0, hshift, scaling) - if (is_root_pe()) & - call chk_sum_msg("v-point:", bc0, bcSW, bcSE, bcNW, bcNE, mesg, iounit) - else - if (sym) then - bcS = subchk(array, HI, 0, -hshift-1, scaling) - else - bcS = subchk(array, HI, 0, -hshift, scaling) + if (is_root_pe()) & + call chk_sum_msg_NSEW("v-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) endif - bcE = subchk(array, HI, hshift, 0, scaling) - bcW = subchk(array, HI, -hshift, 0, scaling) - bcN = subchk(array, HI, 0, hshift, scaling) + endif - if (is_root_pe()) & - call chk_sum_msg_NSEW("v-point:", bc0, bcN, bcS, bcE, bcW, mesg, iounit) + if (writeHash .and. is_root_pe()) then + allocate(hash_array(HI%isc:HI%iec, HI%jsc:HI%jec, size(array, 3))) + hash_array(:,:,:) = scaling * array(HI%isc:HI%iec, HI%jsc:HI%jec, :) + + write(iounit, '("v-point: hash=", z8, 1x, a)') & + murmur_hash(hash_array), mesg + deallocate(hash_array) endif contains diff --git a/src/framework/MOM_murmur_hash.F90 b/src/framework/MOM_murmur_hash.F90 new file mode 100644 index 0000000000..16283f61e3 --- /dev/null +++ b/src/framework/MOM_murmur_hash.F90 @@ -0,0 +1,251 @@ +!> MurmurHash is a non-cryptographic hash function developed by Austin Appleby. +!! +!! This module provides an implementation of the 32-bit MurmurHash3 algorithm. +!! It is used in MOM6 to generate unique hashes of field arrays. The hash is +!! sensitive to order of elements and can detect changes that would otherwise +!! be missed by the mean/min/max/bitcount tests. +!! +!! Sensitivity to order means that it must be used with care for tests such as +!! processor layout. +!! +!! This implementation assumes data sizes of either 32 or 64 bits. It cannot +!! be used for smaller types such as strings. +!! +!! https://github.com/aappleby/smhasher +module MOM_murmur_hash + +use, intrinsic :: iso_fortran_env, only : int32, int64, real32, real64 + +implicit none ; private + +public :: murmur_hash + +!> Return the murmur3 hash of an array. +interface murmur_hash + procedure murmurhash3_i32 + procedure murmurhash3_i64 + procedure murmurhash3_r32 + procedure murmurhash3_r32_1d + procedure murmurhash3_r32_2d + procedure murmurhash3_r32_3d + procedure murmurhash3_r32_4d + procedure murmurhash3_r64 + procedure murmurhash3_r64_1d + procedure murmurhash3_r64_2d + procedure murmurhash3_r64_3d + procedure murmurhash3_r64_4d +end interface murmur_hash + +contains + +!> Return the murmur3 hash for a 32-bit integer array. +function murmurhash3_i32(key, seed) result(hash) + integer(int32), intent(in) :: key(:) + !< Input array + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32), parameter :: c1 = int(z'cc9e2d51', kind=int32) + integer(int32), parameter :: c2 = int(z'1b873593', kind=int32) + integer(int32), parameter :: c3 = int(z'e6546b64', kind=int32) + + integer(int32), parameter :: c4 = int(z'85ebca6b', kind=int32) + integer(int32), parameter :: c5 = int(z'c2b2ae35', kind=int32) + + integer :: i + integer(int32) :: k + + hash = 0 + if (present(seed)) hash = seed + + do i = 1, size(key) + k = key(i) + k = k * c1 + k = ishftc(k, 15) + k = k * c2 + + hash = ieor(hash, k) + hash = ishftc(hash, 13) + hash = 5 * hash + c3 + enddo + + ! NOTE: This is the point where the algorithm would handle trailing bytes. + ! Since our arrays are comprised of 4 or 8 byte elements, we skip this part. + + hash = ieor(hash, 4*size(key)) + + hash = ieor(hash, ishft(hash, -16)) + hash = hash * c4 + hash = ieor(hash, ishft(hash, -13)) + hash = hash * c5 + hash = ieor(hash, ishft(hash, -16)) +end function murmurhash3_i32 + + +!> Return the murmur3 hash for a 64-bit integer array. +function murmurhash3_i64(key, seed) result(hash) + integer(int64), intent(in) :: key(:) + !< Input array + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(2*size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_i64 + + +!> Return the murmur3 hash for a 32-bit real array. +function murmurhash3_r32(key, seed) result(hash) + real(real32), intent(in) :: key + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(1) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r32 + + +!> Return the murmur3 hash for a 32-bit real array. +function murmurhash3_r32_1d(key, seed) result(hash) + real(real32), intent(in) :: key(:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r32_1d + + +!> Return the murmur3 hash for a 32-bit real 2D array. +function murmurhash3_r32_2d(key, seed) result(hash) + real(real32), intent(in) :: key(:,:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r32_2d + + +!> Return the murmur3 hash for a 32-bit real 3D array. +function murmurhash3_r32_3d(key, seed) result(hash) + real(real32), intent(in) :: key(:,:,:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r32_3d + + +!> Return the murmur3 hash for a 32-bit real 4D array. +function murmurhash3_r32_4d(key, seed) result(hash) + real(real32), intent(in) :: key(:,:,:,:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r32_4d + + +!> Return the murmur3 hash for a 64-bit real array. +function murmurhash3_r64(key, seed) result(hash) + real(real64), intent(in) :: key + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(2) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r64 + + +!> Return the murmur3 hash for a 64-bit real array. +function murmurhash3_r64_1d(key, seed) result(hash) + real(real64), intent(in) :: key(:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(2*size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r64_1d + + +!> Return the murmur3 hash for a 64-bit real 2D array. +function murmurhash3_r64_2d(key, seed) result(hash) + real(real64), intent(in) :: key(:,:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(2*size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r64_2d + + +!> Return the murmur3 hash for a 64-bit real 3D array. +function murmurhash3_r64_3d(key, seed) result(hash) + real(real64), intent(in) :: key(:,:,:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(2*size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r64_3d + + +!> Return the murmur3 hash for a 64-bit real 4D array. +function murmurhash3_r64_4d(key, seed) result(hash) + real(real64), intent(in) :: key(:,:,:,:) + !< Input array [arbitrary] + integer(int32), intent(in), optional :: seed + !< Hash seed + integer(int32) :: hash + !< Murmur hash of array + + integer(int32) :: ikey(2*size(key)) + + hash = murmur_hash(transfer(key, ikey), seed=seed) +end function murmurhash3_r64_4d + +end module MOM_murmur_hash