MOM6
mom_checksums::hchksum Interface Reference

Detailed Description

Checksums an array (2d or 3d) staggered at tracer points.

Definition at line 48 of file MOM_checksums.F90.

Private functions

subroutine chksum_h_2d (array, mesg, HI, haloshift, omit_corners, scale, logunit)
 Checksums a 2d array staggered at tracer points. More...
 
subroutine chksum_h_3d (array, mesg, HI, haloshift, omit_corners, scale, logunit)
 Checksums a 3d array staggered at tracer points. More...
 

Functions and subroutines

◆ chksum_h_2d()

subroutine mom_checksums::hchksum::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 

◆ chksum_h_3d()

subroutine mom_checksums::hchksum::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 

The documentation for this interface was generated from the following file:
substats
subroutine substats(array, aMean, aMin, aMax)
Definition: MOM_checksums.F90:179
subchk
integer function subchk(array, scale)
Definition: MOM_checksums.F90:167