MOM6
mom_coms Module Reference

Detailed Description

Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums facility.

Data Types

interface  assignment(=)
 Copy the value of one extended-fixed-point number into another. More...
 
type  efp_type
 The Extended Fixed Point (EFP) type provides a public interface for doing sums and taking differences with this type. More...
 
interface  operator(+)
 Add two extended-fixed-point numbers. More...
 
interface  operator(-)
 Subtract one extended-fixed-point number from another. More...
 
interface  reproducing_sum
 Find an accurate and order-invariant sum of distributed 2d or 3d fields. More...
 

Functions/Subroutines

real function reproducing_sum_2d (array, isr, ier, jsr, jer, EFP_sum, reproducing, overflow_check, err)
 This subroutine uses a conversion to an integer representation of real numbers to give an order-invariant sum of distributed 2-D arrays that reproduces across domain decomposition. This technique is described in Hallberg & Adcroft, 2014, Parallel Computing, doi:10.1016/j.parco.2014.04.007. More...
 
real function reproducing_sum_3d (array, isr, ier, jsr, jer, sums, EFP_sum, err)
 This subroutine uses a conversion to an integer representation of real numbers to give an order-invariant sum of distributed 3-D arrays that reproduces across domain decomposition. This technique is described in Hallberg & Adcroft, 2014, Parallel Computing, doi:10.1016/j.parco.2014.04.007. More...
 
integer(kind=8) function, dimension(nireal_to_ints (r, prec_error, overflow)
 Convert a real number into the array of integers constitute its extended-fixed-point representation. More...
 
real function ints_to_real (ints)
 Convert the array of integers that constitute an extended-fixed-point representation into a real number. More...
 
subroutine increment_ints (int_sum, int2, prec_error)
 Increment an array of integers that constitutes an extended-fixed-point representation with a another EFP number. More...
 
subroutine increment_ints_faster (int_sum, r, max_mag_term)
 Increment an EFP number with a real number without doing any carrying of of overflows and using only minimal error checking. More...
 
subroutine carry_overflow (int_sum, prec_error)
 This subroutine handles carrying of the overflow. More...
 
subroutine regularize_ints (int_sum)
 This subroutine carries the overflow, and then makes sure that all integers are of the same sign as the overall value. More...
 
logical function, public query_efp_overflow_error ()
 Returns the status of the module's error flag. More...
 
subroutine, public reset_efp_overflow_error ()
 Reset the module's error flag to false. More...
 
type(efp_type) function, public efp_plus (EFP1, EFP2)
 Add two extended-fixed-point numbers. More...
 
type(efp_type) function, public efp_minus (EFP1, EFP2)
 Subract one extended-fixed-point number from another. More...
 
subroutine efp_assign (EFP1, EFP2)
 Copy one extended-fixed-point number into another. More...
 
real function, public efp_to_real (EFP1)
 Return the real number that an extended-fixed-point number corresponds with. More...
 
real function, public efp_real_diff (EFP1, EFP2)
 Take the difference between two extended-fixed-point numbers (EFP1 - EFP2) and return the result as a real number. More...
 
type(efp_type) function, public real_to_efp (val, overflow)
 Return the extended-fixed-point number that a real number corresponds with. More...
 
subroutine, public efp_list_sum_across_pes (EFPs, nval, errors)
 This subroutine does a sum across PEs of a list of EFP variables, returning the sums in place, with all overflows carried. More...
 
subroutine, public mom_infra_end
 This subroutine carries out all of the calls required to close out the infrastructure cleanly. This should only be called in ocean-only runs, as the coupler takes care of this in coupled runs. More...
 

Variables

integer(kind=8), parameter prec =2_8**46
 The precision of each integer. More...
 
real, parameter r_prec =2.0**46
 A real version of prec. More...
 
real, parameter i_prec =1.0/(2.0**46)
 The inverse of prec. More...
 
integer, parameter max_count_prec =2**(63-46)-1
 The number of values that can be added together with the current value of prec before there will be roundoff problems. More...
 
integer, parameter ni =6
 The number of long integers to use to represent a real number. More...
 
real, dimension(ni), parameter pr = (/ r_prec**2, r_prec, 1.0, 1.0/r_prec, 1.0/r_prec**2, 1.0/r_prec**3 /)
 An array of the real precision of each of the integers. More...
 
real, dimension(ni), parameter i_pr = (/ 1.0/r_prec**2, 1.0/r_prec, 1.0, r_prec, r_prec**2, r_prec**3 /)
 An array of the inverse of the real precision of each of the integers. More...
 
real, parameter max_efp_float = pr(1) * (2.**63 - 1.)
 The largest float with an EFP representation. NOTE: Only the first bin can exceed precision, but is bounded by the largest signed integer. More...
 
logical overflow_error = .false.
 This becomes true if an overflow is encountered. More...
 
logical nan_error = .false.
 This becomes true if a NaN is encountered. More...
 
logical debug = .false.
 Making this true enables debugging output. More...
 

Function/Subroutine Documentation

◆ carry_overflow()

subroutine mom_coms::carry_overflow ( integer(kind=8), dimension(ni), intent(inout)  int_sum,
integer(kind=8), intent(in)  prec_error 
)
private

This subroutine handles carrying of the overflow.

Parameters
[in,out]int_sumThe array of EFP integers being modified by carries, but without changing value.
[in]prec_errorThe PE-count dependent precision of the integers that is safe from overflows during global sums. This will be larger than the compile-time precision parameter, and is used to detect overflows.

Definition at line 538 of file MOM_coms.F90.

538  integer(kind=8), dimension(ni), intent(inout) :: int_sum !< The array of EFP integers being
539  !! modified by carries, but without changing value.
540  integer(kind=8), intent(in) :: prec_error !< The PE-count dependent precision of the
541  !! integers that is safe from overflows during global
542  !! sums. This will be larger than the compile-time
543  !! precision parameter, and is used to detect overflows.
544 
545  ! This subroutine handles carrying of the overflow.
546  integer :: i, num_carry
547 
548  do i=ni,2,-1 ; if (abs(int_sum(i)) >= prec) then
549  num_carry = int(int_sum(i) * i_prec)
550  int_sum(i) = int_sum(i) - num_carry*prec
551  int_sum(i-1) = int_sum(i-1) + num_carry
552  endif ; enddo
553  if (abs(int_sum(1)) > prec_error) then
554  overflow_error = .true.
555  endif
556 

References i_prec, ni, overflow_error, and prec.

Referenced by efp_list_sum_across_pes(), reproducing_sum_2d(), and reproducing_sum_3d().

Here is the caller graph for this function:

◆ efp_assign()

subroutine mom_coms::efp_assign ( type(efp_type), intent(out)  EFP1,
type(efp_type), intent(in)  EFP2 
)
private

Copy one extended-fixed-point number into another.

Parameters
[out]efp1The recipient extended fixed point number
[in]efp2The source extended fixed point number

Definition at line 638 of file MOM_coms.F90.

638  type(EFP_type), intent(out) :: EFP1 !< The recipient extended fixed point number
639  type(EFP_type), intent(in) :: EFP2 !< The source extended fixed point number
640  integer i
641  ! This subroutine assigns all components of the extended fixed point type
642  ! variable on the RHS (EFP2) to the components of the variable on the LHS
643  ! (EFP1).
644 
645  do i=1,ni ; efp1%v(i) = efp2%v(i) ; enddo

References ni.

◆ efp_list_sum_across_pes()

subroutine, public mom_coms::efp_list_sum_across_pes ( type(efp_type), dimension(:), intent(inout)  EFPs,
integer, intent(in)  nval,
logical, dimension(:), intent(out), optional  errors 
)

This subroutine does a sum across PEs of a list of EFP variables, returning the sums in place, with all overflows carried.

Parameters
[in,out]efpsThe list of extended fixed point numbers
[in]nvalThe number of values being summed.
[out]errorsA list of error flags for each sum

Definition at line 698 of file MOM_coms.F90.

698  type(EFP_type), dimension(:), &
699  intent(inout) :: EFPs !< The list of extended fixed point numbers
700  !! being summed across PEs.
701  integer, intent(in) :: nval !< The number of values being summed.
702  logical, dimension(:), &
703  optional, intent(out) :: errors !< A list of error flags for each sum
704 
705  ! This subroutine does a sum across PEs of a list of EFP variables,
706  ! returning the sums in place, with all overflows carried.
707 
708  integer(kind=8), dimension(ni,nval) :: ints
709  integer(kind=8) :: prec_error
710  logical :: error_found
711  character(len=256) :: mesg
712  integer :: i, n
713 
714  if (num_pes() > max_count_prec) call mom_error(fatal, &
715  "reproducing_sum: Too many processors are being used for the value of "//&
716  "prec. Reduce prec to (2^63-1)/num_PEs.")
717 
718  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
719  ! overflow_error is an overflow error flag for the whole module.
720  overflow_error = .false. ; error_found = .false.
721 
722  do i=1,nval ; do n=1,ni ; ints(n,i) = efps(i)%v(n) ; enddo ; enddo
723 
724  call sum_across_pes(ints(:,:), ni*nval)
725 
726  if (present(errors)) errors(:) = .false.
727  do i=1,nval
728  overflow_error = .false.
729  call carry_overflow(ints(:,i), prec_error)
730  do n=1,ni ; efps(i)%v(n) = ints(n,i) ; enddo
731  if (present(errors)) errors(i) = overflow_error
732  if (overflow_error) then
733  write (mesg,'("EFP_list_sum_across_PEs error at ",i6," val was ",ES12.6, ", prec_error = ",ES12.6)') &
734  i, efp_to_real(efps(i)), real(prec_error)
735  call mom_error(warning, mesg)
736  endif
737  error_found = error_found .or. overflow_error
738  enddo
739  if (error_found .and. .not.(present(errors))) then
740  call mom_error(fatal, "Overflow in EFP_list_sum_across_PEs.")
741  endif
742 

References carry_overflow(), efp_to_real(), max_count_prec, mom_error_handler::mom_error(), ni, and overflow_error.

Here is the call graph for this function:

◆ efp_minus()

type(efp_type) function, public mom_coms::efp_minus ( type(efp_type), intent(in)  EFP1,
type(efp_type), intent(in)  EFP2 
)

Subract one extended-fixed-point number from another.

Returns
The result in extended fixed point format
Parameters
[in]efp1The first extended fixed point number
[in]efp2The extended fixed point number being subtracted from the first extended fixed point number

Definition at line 625 of file MOM_coms.F90.

625  type(EFP_type) :: EFP_minus !< The result in extended fixed point format
626  type(EFP_type), intent(in) :: EFP1 !< The first extended fixed point number
627  type(EFP_type), intent(in) :: EFP2 !< The extended fixed point number being
628  !! subtracted from the first extended fixed point number
629  integer :: i
630 
631  do i=1,ni ; efp_minus%v(i) = -1*efp2%v(i) ; enddo
632 
633  call increment_ints(efp_minus%v(:), efp1%v(:))

References increment_ints(), and ni.

Here is the call graph for this function:

◆ efp_plus()

type(efp_type) function, public mom_coms::efp_plus ( type(efp_type), intent(in)  EFP1,
type(efp_type), intent(in)  EFP2 
)

Add two extended-fixed-point numbers.

Returns
The result in extended fixed point format
Parameters
[in]efp1The first extended fixed point number
[in]efp2The second extended fixed point number

Definition at line 614 of file MOM_coms.F90.

614  type(EFP_type) :: EFP_plus !< The result in extended fixed point format
615  type(EFP_type), intent(in) :: EFP1 !< The first extended fixed point number
616  type(EFP_type), intent(in) :: EFP2 !< The second extended fixed point number
617 
618  efp_plus = efp1
619 
620  call increment_ints(efp_plus%v(:), efp2%v(:))

References increment_ints().

Here is the call graph for this function:

◆ efp_real_diff()

real function, public mom_coms::efp_real_diff ( type(efp_type), intent(in)  EFP1,
type(efp_type), intent(in)  EFP2 
)

Take the difference between two extended-fixed-point numbers (EFP1 - EFP2) and return the result as a real number.

Parameters
[in]efp1The first extended fixed point number
[in]efp2The extended fixed point number being subtracted from the first extended fixed point number
Returns
The real result

Definition at line 660 of file MOM_coms.F90.

660  type(EFP_type), intent(in) :: EFP1 !< The first extended fixed point number
661  type(EFP_type), intent(in) :: EFP2 !< The extended fixed point number being
662  !! subtracted from the first extended fixed point number
663  real :: EFP_real_diff !< The real result
664 
665  type(EFP_type) :: EFP_diff
666 
667  efp_diff = efp1 - efp2
668  efp_real_diff = efp_to_real(efp_diff)
669 

References efp_to_real().

Here is the call graph for this function:

◆ efp_to_real()

real function, public mom_coms::efp_to_real ( type(efp_type), intent(inout)  EFP1)

Return the real number that an extended-fixed-point number corresponds with.

Parameters
[in,out]efp1The extended fixed point number being converted

Definition at line 650 of file MOM_coms.F90.

650  type(EFP_type), intent(inout) :: EFP1 !< The extended fixed point number being converted
651  real :: EFP_to_real
652 
653  call regularize_ints(efp1%v)
654  efp_to_real = ints_to_real(efp1%v)

References ints_to_real(), and regularize_ints().

Referenced by efp_list_sum_across_pes(), and efp_real_diff().

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

◆ increment_ints()

subroutine mom_coms::increment_ints ( integer(kind=8), dimension(ni), intent(inout)  int_sum,
integer(kind=8), dimension(ni), intent(in)  int2,
integer(kind=8), intent(in), optional  prec_error 
)
private

Increment an array of integers that constitutes an extended-fixed-point representation with a another EFP number.

Parameters
[in,out]int_sumThe array of EFP integers being incremented
[in]int2The array of EFP integers being added
[in]prec_errorThe PE-count dependent precision of the integers that is safe from overflows during global sums. This will be larger than the compile-time precision parameter, and is used to detect overflows.

Definition at line 472 of file MOM_coms.F90.

472  integer(kind=8), dimension(ni), intent(inout) :: int_sum !< The array of EFP integers being incremented
473  integer(kind=8), dimension(ni), intent(in) :: int2 !< The array of EFP integers being added
474  integer(kind=8), optional, intent(in) :: prec_error !< The PE-count dependent precision of the
475  !! integers that is safe from overflows during global
476  !! sums. This will be larger than the compile-time
477  !! precision parameter, and is used to detect overflows.
478 
479  ! This subroutine increments a number with another, both using the integer
480  ! representation in real_to_ints.
481  integer :: i
482 
483  do i=ni,2,-1
484  int_sum(i) = int_sum(i) + int2(i)
485  ! Carry the local overflow.
486  if (int_sum(i) > prec) then
487  int_sum(i) = int_sum(i) - prec
488  int_sum(i-1) = int_sum(i-1) + 1
489  elseif (int_sum(i) < -prec) then
490  int_sum(i) = int_sum(i) + prec
491  int_sum(i-1) = int_sum(i-1) - 1
492  endif
493  enddo
494  int_sum(1) = int_sum(1) + int2(1)
495  if (present(prec_error)) then
496  if (abs(int_sum(1)) > prec_error) overflow_error = .true.
497  else
498  if (abs(int_sum(1)) > prec) overflow_error = .true.
499  endif
500 

References ni, overflow_error, and prec.

Referenced by efp_minus(), efp_plus(), reproducing_sum_2d(), and reproducing_sum_3d().

Here is the caller graph for this function:

◆ increment_ints_faster()

subroutine mom_coms::increment_ints_faster ( integer(kind=8), dimension(ni), intent(inout)  int_sum,
real, intent(in)  r,
real, intent(inout)  max_mag_term 
)
private

Increment an EFP number with a real number without doing any carrying of of overflows and using only minimal error checking.

Parameters
[in,out]int_sumThe array of EFP integers being incremented
[in]rThe real number being added.
[in,out]max_mag_termA running maximum magnitude of the r's.

Definition at line 506 of file MOM_coms.F90.

506  integer(kind=8), dimension(ni), intent(inout) :: int_sum !< The array of EFP integers being incremented
507  real, intent(in) :: r !< The real number being added.
508  real, intent(inout) :: max_mag_term !< A running maximum magnitude of the r's.
509 
510  ! This subroutine increments a number with another, both using the integer
511  ! representation in real_to_ints, but without doing any carrying of overflow.
512  ! The entire operation is embedded in a single call for greater speed.
513  real :: rs
514  integer(kind=8) :: ival
515  integer :: sgn, i
516 
517  if ((r >= 1e30) .eqv. (r < 1e30)) then ; nan_error = .true. ; return ; endif
518  sgn = 1 ; if (r<0.0) sgn = -1
519  rs = abs(r)
520  if (rs > abs(max_mag_term)) max_mag_term = r
521 
522  ! Abort if the number has no EFP representation
523  if (rs > max_efp_float) then
524  overflow_error = .true.
525  return
526  endif
527 
528  do i=1,ni
529  ival = int(rs*i_pr(i), 8)
530  rs = rs - ival*pr(i)
531  int_sum(i) = int_sum(i) + sgn*ival
532  enddo
533 

References i_pr, max_efp_float, nan_error, ni, overflow_error, and pr.

Referenced by reproducing_sum_2d(), and reproducing_sum_3d().

Here is the caller graph for this function:

◆ ints_to_real()

real function mom_coms::ints_to_real ( integer(kind=8), dimension(ni), intent(in)  ints)
private

Convert the array of integers that constitute an extended-fixed-point representation into a real number.

Parameters
[in]intsThe array of EFP integers

Definition at line 459 of file MOM_coms.F90.

459  integer(kind=8), dimension(ni), intent(in) :: ints !< The array of EFP integers
460  real :: r
461  ! This subroutine reverses the conversion in real_to_ints.
462 
463  integer :: i
464 
465  r = 0.0
466  do i=1,ni ; r = r + pr(i)*ints(i) ; enddo

References ni, and pr.

Referenced by efp_to_real(), reproducing_sum_2d(), and reproducing_sum_3d().

Here is the caller graph for this function:

◆ mom_infra_end()

subroutine, public mom_coms::mom_infra_end ( )

This subroutine carries out all of the calls required to close out the infrastructure cleanly. This should only be called in ocean-only runs, as the coupler takes care of this in coupled runs.

Definition at line 748 of file MOM_coms.F90.

748  call print_memuse_stats( 'Memory HiWaterMark', always=.true. )
749  call fms_end

◆ query_efp_overflow_error()

logical function, public mom_coms::query_efp_overflow_error ( )

Returns the status of the module's error flag.

Definition at line 603 of file MOM_coms.F90.

603  logical :: query_EFP_overflow_error
604  query_efp_overflow_error = overflow_error

References overflow_error.

Referenced by mom_spatial_means::global_i_mean(), and mom_spatial_means::global_j_mean().

Here is the caller graph for this function:

◆ real_to_efp()

type(efp_type) function, public mom_coms::real_to_efp ( real, intent(in)  val,
logical, intent(inout), optional  overflow 
)

Return the extended-fixed-point number that a real number corresponds with.

Parameters
[in]valThe real number being converted
[in,out]overflowReturns true if the conversion is being done on a value that is too large to be represented

Definition at line 674 of file MOM_coms.F90.

674  real, intent(in) :: val !< The real number being converted
675  logical, optional, intent(inout) :: overflow !< Returns true if the conversion is being
676  !! done on a value that is too large to be represented
677  type(EFP_type) :: real_to_EFP
678 
679  logical :: over
680  character(len=80) :: mesg
681 
682  if (present(overflow)) then
683  real_to_efp%v(:) = real_to_ints(val, overflow=overflow)
684  else
685  over = .false.
686  real_to_efp%v(:) = real_to_ints(val, overflow=over)
687  if (over) then
688  write(mesg, '(ES13.5)') val
689  call mom_error(fatal,"Overflow in real_to_EFP conversion of "//trim(mesg))
690  endif
691  endif
692 

References mom_error_handler::mom_error(), and real_to_ints().

Here is the call graph for this function:

◆ real_to_ints()

integer(kind=8) function, dimension(ni) mom_coms::real_to_ints ( real, intent(in)  r,
integer(kind=8), intent(in), optional  prec_error,
logical, intent(inout), optional  overflow 
)
private

Convert a real number into the array of integers constitute its extended-fixed-point representation.

Parameters
[in]rThe real number being converted
[in]prec_errorThe PE-count dependent precision of the integers that is safe from overflows during global sums. This will be larger than the compile-time precision parameter, and is used to detect overflows.
[in,out]overflowReturns true if the conversion is being done on a value that is too large to be represented

Definition at line 417 of file MOM_coms.F90.

417  real, intent(in) :: r !< The real number being converted
418  integer(kind=8), optional, intent(in) :: prec_error !< The PE-count dependent precision of the
419  !! integers that is safe from overflows during global
420  !! sums. This will be larger than the compile-time
421  !! precision parameter, and is used to detect overflows.
422  logical, optional, intent(inout) :: overflow !< Returns true if the conversion is being
423  !! done on a value that is too large to be represented
424  integer(kind=8), dimension(ni) :: ints
425  ! This subroutine converts a real number to an equivalent representation
426  ! using several long integers.
427 
428  real :: rs
429  character(len=80) :: mesg
430  integer(kind=8) :: ival, prec_err
431  integer :: sgn, i
432 
433  prec_err = prec ; if (present(prec_error)) prec_err = prec_error
434  ints(:) = 0_8
435  if ((r >= 1e30) .eqv. (r < 1e30)) then ; nan_error = .true. ; return ; endif
436 
437  sgn = 1 ; if (r<0.0) sgn = -1
438  rs = abs(r)
439 
440  if (present(overflow)) then
441  if (.not.(rs < prec_err*pr(1))) overflow = .true.
442  if ((r >= 1e30) .eqv. (r < 1e30)) overflow = .true.
443  elseif (.not.(rs < prec_err*pr(1))) then
444  write(mesg, '(ES13.5)') r
445  call mom_error(fatal,"Overflow in real_to_ints conversion of "//trim(mesg))
446  endif
447 
448  do i=1,ni
449  ival = int(rs*i_pr(i), 8)
450  rs = rs - ival*pr(i)
451  ints(i) = sgn*ival
452  enddo
453 

References i_pr, mom_error_handler::mom_error(), nan_error, ni, pr, and prec.

Referenced by real_to_efp(), reproducing_sum_2d(), and reproducing_sum_3d().

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

◆ regularize_ints()

subroutine mom_coms::regularize_ints ( integer(kind=8), dimension(ni), intent(inout)  int_sum)
private

This subroutine carries the overflow, and then makes sure that all integers are of the same sign as the overall value.

Parameters
[in,out]int_sumThe array of integers being modified to take a

Definition at line 562 of file MOM_coms.F90.

562  integer(kind=8), dimension(ni), &
563  intent(inout) :: int_sum !< The array of integers being modified to take a
564  !! regular form with all integers of the same sign,
565  !! but without changing value.
566 
567  ! This subroutine carries the overflow, and then makes sure that
568  ! all integers are of the same sign as the overall value.
569  logical :: positive
570  integer :: i, num_carry
571 
572  do i=ni,2,-1 ; if (abs(int_sum(i)) >= prec) then
573  num_carry = int(int_sum(i) * i_prec)
574  int_sum(i) = int_sum(i) - num_carry*prec
575  int_sum(i-1) = int_sum(i-1) + num_carry
576  endif ; enddo
577 
578  ! Determine the sign of the final number.
579  positive = .true.
580  do i=1,ni
581  if (abs(int_sum(i)) > 0) then
582  if (int_sum(i) < 0) positive = .false.
583  exit
584  endif
585  enddo
586 
587  if (positive) then
588  do i=ni,2,-1 ; if (int_sum(i) < 0) then
589  int_sum(i) = int_sum(i) + prec
590  int_sum(i-1) = int_sum(i-1) - 1
591  endif ; enddo
592  else
593  do i=ni,2,-1 ; if (int_sum(i) > 0) then
594  int_sum(i) = int_sum(i) - prec
595  int_sum(i-1) = int_sum(i-1) + 1
596  endif ; enddo
597  endif
598 

References i_prec, ni, and prec.

Referenced by efp_to_real(), reproducing_sum_2d(), and reproducing_sum_3d().

Here is the caller graph for this function:

◆ reproducing_sum_2d()

real function mom_coms::reproducing_sum_2d ( real, dimension(:,:), intent(in)  array,
integer, intent(in), optional  isr,
integer, intent(in), optional  ier,
integer, intent(in), optional  jsr,
integer, intent(in), optional  jer,
type(efp_type), intent(out), optional  EFP_sum,
logical, intent(in), optional  reproducing,
logical, intent(in), optional  overflow_check,
integer, intent(out), optional  err 
)
private

This subroutine uses a conversion to an integer representation of real numbers to give an order-invariant sum of distributed 2-D arrays that reproduces across domain decomposition. This technique is described in Hallberg & Adcroft, 2014, Parallel Computing, doi:10.1016/j.parco.2014.04.007.

Parameters
[in]arrayThe array to be summed
[in]isrThe starting i-index of the sum, noting that the array indices starts at 1
[in]ierThe ending i-index of the sum, noting that the array indices starts at 1
[in]jsrThe starting j-index of the sum, noting that the array indices starts at 1
[in]jerThe ending j-index of the sum, noting that the array indices starts at 1
[out]efp_sumThe result in extended fixed point format
[in]reproducingIf present and false, do the sum using the naive non-reproducing approach
[in]overflow_checkIf present and false, disable checking for overflows in incremental results. This can speed up calculations if the number of values being summed is small enough
[out]errIf present, return an error code instead of triggering any fatal errors directly from this routine.
Returns
Result

Definition at line 83 of file MOM_coms.F90.

83  real, dimension(:,:), intent(in) :: array !< The array to be summed
84  integer, optional, intent(in) :: isr !< The starting i-index of the sum, noting
85  !! that the array indices starts at 1
86  integer, optional, intent(in) :: ier !< The ending i-index of the sum, noting
87  !! that the array indices starts at 1
88  integer, optional, intent(in) :: jsr !< The starting j-index of the sum, noting
89  !! that the array indices starts at 1
90  integer, optional, intent(in) :: jer !< The ending j-index of the sum, noting
91  !! that the array indices starts at 1
92  type(EFP_type), optional, intent(out) :: EFP_sum !< The result in extended fixed point format
93  logical, optional, intent(in) :: reproducing !< If present and false, do the sum
94  !! using the naive non-reproducing approach
95  logical, optional, intent(in) :: overflow_check !< If present and false, disable
96  !! checking for overflows in incremental results.
97  !! This can speed up calculations if the number
98  !! of values being summed is small enough
99  integer, optional, intent(out) :: err !< If present, return an error code instead of
100  !! triggering any fatal errors directly from
101  !! this routine.
102  real :: sum !< Result
103 
104  ! This subroutine uses a conversion to an integer representation
105  ! of real numbers to give order-invariant sums that will reproduce
106  ! across PE count. This idea comes from R. Hallberg and A. Adcroft.
107 
108  integer(kind=8), dimension(ni) :: ints_sum
109  integer(kind=8) :: ival, prec_error
110  real :: rsum(1), rs
111  real :: max_mag_term
112  logical :: repro, over_check
113  character(len=256) :: mesg
114  integer :: i, j, n, is, ie, js, je, sgn
115 
116  if (num_pes() > max_count_prec) call mom_error(fatal, &
117  "reproducing_sum: Too many processors are being used for the value of "//&
118  "prec. Reduce prec to (2^63-1)/num_PEs.")
119 
120  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
121 
122  is = 1 ; ie = size(array,1) ; js = 1 ; je = size(array,2 )
123  if (present(isr)) then
124  if (isr < is) call mom_error(fatal, &
125  "Value of isr too small in reproducing_sum_2d.")
126  is = isr
127  endif
128  if (present(ier)) then
129  if (ier > ie) call mom_error(fatal, &
130  "Value of ier too large in reproducing_sum_2d.")
131  ie = ier
132  endif
133  if (present(jsr)) then
134  if (jsr < js) call mom_error(fatal, &
135  "Value of jsr too small in reproducing_sum_2d.")
136  js = jsr
137  endif
138  if (present(jer)) then
139  if (jer > je) call mom_error(fatal, &
140  "Value of jer too large in reproducing_sum_2d.")
141  je = jer
142  endif
143 
144  repro = .true. ; if (present(reproducing)) repro = reproducing
145  over_check = .true. ; if (present(overflow_check)) over_check = overflow_check
146 
147  if (repro) then
148  overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
149  ints_sum(:) = 0
150  if (over_check) then
151  if ((je+1-js)*(ie+1-is) < max_count_prec) then
152  do j=js,je ; do i=is,ie
153  call increment_ints_faster(ints_sum, array(i,j), max_mag_term)
154  enddo ; enddo
155  call carry_overflow(ints_sum, prec_error)
156  elseif ((ie+1-is) < max_count_prec) then
157  do j=js,je
158  do i=is,ie
159  call increment_ints_faster(ints_sum, array(i,j), max_mag_term)
160  enddo
161  call carry_overflow(ints_sum, prec_error)
162  enddo
163  else
164  do j=js,je ; do i=is,ie
165  call increment_ints(ints_sum, real_to_ints(array(i,j), prec_error), &
166  prec_error)
167  enddo ; enddo
168  endif
169  else
170  do j=js,je ; do i=is,ie
171  sgn = 1 ; if (array(i,j)<0.0) sgn = -1
172  rs = abs(array(i,j))
173  do n=1,ni
174  ival = int(rs*i_pr(n), 8)
175  rs = rs - ival*pr(n)
176  ints_sum(n) = ints_sum(n) + sgn*ival
177  enddo
178  enddo ; enddo
179  call carry_overflow(ints_sum, prec_error)
180  endif
181 
182  if (present(err)) then
183  err = 0
184  if (overflow_error) &
185  err = err+2
186  if (nan_error) &
187  err = err+4
188  if (err > 0) then ; do n=1,ni ; ints_sum(n) = 0 ; enddo ; endif
189  else
190  if (nan_error) then
191  call mom_error(fatal, "NaN in input field of reproducing_sum(_2d).")
192  endif
193  if (abs(max_mag_term) >= prec_error*pr(1)) then
194  write(mesg, '(ES13.5)') max_mag_term
195  call mom_error(fatal,"Overflow in reproducing_sum(_2d) conversion of "//trim(mesg))
196  endif
197  if (overflow_error) then
198  call mom_error(fatal, "Overflow in reproducing_sum(_2d).")
199  endif
200  endif
201 
202  call sum_across_pes(ints_sum, ni)
203 
204  call regularize_ints(ints_sum)
205  sum = ints_to_real(ints_sum)
206  else
207  rsum(1) = 0.0
208  do j=js,je ; do i=is,ie
209  rsum(1) = rsum(1) + array(i,j)
210  enddo ; enddo
211  call sum_across_pes(rsum,1)
212  sum = rsum(1)
213 
214  if (present(err)) then ; err = 0 ; endif
215 
216  if (debug .or. present(efp_sum)) then
217  overflow_error = .false.
218  ints_sum = real_to_ints(sum, prec_error, overflow_error)
219  if (overflow_error) then
220  if (present(err)) then
221  err = err + 2
222  else
223  write(mesg, '(ES13.5)') sum
224  call mom_error(fatal,"Repro_sum_2d: Overflow in real_to_ints conversion of "//trim(mesg))
225  endif
226  endif
227  endif
228  endif
229 
230  if (present(efp_sum)) efp_sum%v(:) = ints_sum(:)
231 
232  if (debug) then
233  write(mesg,'("2d RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
234  call mom_mesg(mesg, 3)
235  endif
236 

References carry_overflow(), debug, i_pr, increment_ints(), increment_ints_faster(), ints_to_real(), max_count_prec, mom_error_handler::mom_error(), mom_error_handler::mom_mesg(), nan_error, ni, overflow_error, pr, real_to_ints(), and regularize_ints().

Here is the call graph for this function:

◆ reproducing_sum_3d()

real function mom_coms::reproducing_sum_3d ( real, dimension(:,:,:), intent(in)  array,
integer, intent(in), optional  isr,
integer, intent(in), optional  ier,
integer, intent(in), optional  jsr,
integer, intent(in), optional  jer,
real, dimension(:), intent(out), optional  sums,
type(efp_type), intent(out), optional  EFP_sum,
integer, intent(out), optional  err 
)
private

This subroutine uses a conversion to an integer representation of real numbers to give an order-invariant sum of distributed 3-D arrays that reproduces across domain decomposition. This technique is described in Hallberg & Adcroft, 2014, Parallel Computing, doi:10.1016/j.parco.2014.04.007.

Parameters
[in]arrayThe array to be summed
[in]isrThe starting i-index of the sum, noting that the array indices starts at 1
[in]ierThe ending i-index of the sum, noting that the array indices starts at 1
[in]jsrThe starting j-index of the sum, noting that the array indices starts at 1
[in]jerThe ending j-index of the sum, noting that the array indices starts at 1
[out]sumsThe sums by vertical layer
[out]efp_sumThe result in extended fixed point format
[out]errIf present, return an error code instead of triggering any fatal errors directly from this routine.
Returns
Result

Definition at line 245 of file MOM_coms.F90.

245  real, dimension(:,:,:), intent(in) :: array !< The array to be summed
246  integer, optional, intent(in) :: isr !< The starting i-index of the sum, noting
247  !! that the array indices starts at 1
248  integer, optional, intent(in) :: ier !< The ending i-index of the sum, noting
249  !! that the array indices starts at 1
250  integer, optional, intent(in) :: jsr !< The starting j-index of the sum, noting
251  !! that the array indices starts at 1
252  integer, optional, intent(in) :: jer !< The ending j-index of the sum, noting
253  !! that the array indices starts at 1
254  real, dimension(:), optional, intent(out) :: sums !< The sums by vertical layer
255  type(EFP_type), optional, intent(out) :: EFP_sum !< The result in extended fixed point format
256  integer, optional, intent(out) :: err !< If present, return an error code instead of
257  !! triggering any fatal errors directly from
258  !! this routine.
259  real :: sum !< Result
260 
261  ! This subroutine uses a conversion to an integer representation
262  ! of real numbers to give order-invariant sums that will reproduce
263  ! across PE count. This idea comes from R. Hallberg and A. Adcroft.
264 
265  real :: max_mag_term
266  integer(kind=8), dimension(ni) :: ints_sum
267  integer(kind=8), dimension(ni,size(array,3)) :: ints_sums
268  integer(kind=8) :: prec_error
269  character(len=256) :: mesg
270  integer :: i, j, k, is, ie, js, je, ke, isz, jsz, n
271 
272  if (num_pes() > max_count_prec) call mom_error(fatal, &
273  "reproducing_sum: Too many processors are being used for the value of "//&
274  "prec. Reduce prec to (2^63-1)/num_PEs.")
275 
276  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
277  max_mag_term = 0.0
278 
279  is = 1 ; ie = size(array,1) ; js = 1 ; je = size(array,2) ; ke = size(array,3)
280  if (present(isr)) then
281  if (isr < is) call mom_error(fatal, &
282  "Value of isr too small in reproducing_sum(_3d).")
283  is = isr
284  endif
285  if (present(ier)) then
286  if (ier > ie) call mom_error(fatal, &
287  "Value of ier too large in reproducing_sum(_3d).")
288  ie = ier
289  endif
290  if (present(jsr)) then
291  if (jsr < js) call mom_error(fatal, &
292  "Value of jsr too small in reproducing_sum(_3d).")
293  js = jsr
294  endif
295  if (present(jer)) then
296  if (jer > je) call mom_error(fatal, &
297  "Value of jer too large in reproducing_sum(_3d).")
298  je = jer
299  endif
300  jsz = je+1-js; isz = ie+1-is
301 
302  if (present(sums)) then
303  if (size(sums) > ke) call mom_error(fatal, "Sums is smaller than "//&
304  "the vertical extent of array in reproducing_sum(_3d).")
305  ints_sums(:,:) = 0
306  overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
307  if (jsz*isz < max_count_prec) then
308  do k=1,ke
309  do j=js,je ; do i=is,ie
310  call increment_ints_faster(ints_sums(:,k), array(i,j,k), max_mag_term)
311  enddo ; enddo
312  call carry_overflow(ints_sums(:,k), prec_error)
313  enddo
314  elseif (isz < max_count_prec) then
315  do k=1,ke ; do j=js,je
316  do i=is,ie
317  call increment_ints_faster(ints_sums(:,k), array(i,j,k), max_mag_term)
318  enddo
319  call carry_overflow(ints_sums(:,k), prec_error)
320  enddo ; enddo
321  else
322  do k=1,ke ; do j=js,je ; do i=is,ie
323  call increment_ints(ints_sums(:,k), &
324  real_to_ints(array(i,j,k), prec_error), prec_error)
325  enddo ; enddo ; enddo
326  endif
327  if (present(err)) then
328  err = 0
329  if (abs(max_mag_term) >= prec_error*pr(1)) err = err+1
330  if (overflow_error) err = err+2
331  if (nan_error) err = err+2
332  if (err > 0) then ; do k=1,ke ; do n=1,ni ; ints_sums(n,k) = 0 ; enddo ; enddo ; endif
333  else
334  if (nan_error) call mom_error(fatal, "NaN in input field of reproducing_sum(_3d).")
335  if (abs(max_mag_term) >= prec_error*pr(1)) then
336  write(mesg, '(ES13.5)') max_mag_term
337  call mom_error(fatal,"Overflow in reproducing_sum(_3d) conversion of "//trim(mesg))
338  endif
339  if (overflow_error) call mom_error(fatal, "Overflow in reproducing_sum(_3d).")
340  endif
341 
342  call sum_across_pes(ints_sums(:,1:ke), ni*ke)
343 
344  sum = 0.0
345  do k=1,ke
346  call regularize_ints(ints_sums(:,k))
347  sums(k) = ints_to_real(ints_sums(:,k))
348  sum = sum + sums(k)
349  enddo
350 
351  if (present(efp_sum)) then
352  efp_sum%v(:) = 0
353  do k=1,ke ; call increment_ints(efp_sum%v(:), ints_sums(:,k)) ; enddo
354  endif
355 
356  if (debug) then
357  do n=1,ni ; ints_sum(n) = 0 ; enddo
358  do k=1,ke ; do n=1,ni ; ints_sum(n) = ints_sum(n) + ints_sums(n,k) ; enddo ; enddo
359  write(mesg,'("3D RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
360  call mom_mesg(mesg, 3)
361  endif
362  else
363  ints_sum(:) = 0
364  overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
365  if (jsz*isz < max_count_prec) then
366  do k=1,ke
367  do j=js,je ; do i=is,ie
368  call increment_ints_faster(ints_sum, array(i,j,k), max_mag_term)
369  enddo ; enddo
370  call carry_overflow(ints_sum, prec_error)
371  enddo
372  elseif (isz < max_count_prec) then
373  do k=1,ke ; do j=js,je
374  do i=is,ie
375  call increment_ints_faster(ints_sum, array(i,j,k), max_mag_term)
376  enddo
377  call carry_overflow(ints_sum, prec_error)
378  enddo ; enddo
379  else
380  do k=1,ke ; do j=js,je ; do i=is,ie
381  call increment_ints(ints_sum, real_to_ints(array(i,j,k), prec_error), &
382  prec_error)
383  enddo ; enddo ; enddo
384  endif
385  if (present(err)) then
386  err = 0
387  if (abs(max_mag_term) >= prec_error*pr(1)) err = err+1
388  if (overflow_error) err = err+2
389  if (nan_error) err = err+2
390  if (err > 0) then ; do n=1,ni ; ints_sum(n) = 0 ; enddo ; endif
391  else
392  if (nan_error) call mom_error(fatal, "NaN in input field of reproducing_sum(_3d).")
393  if (abs(max_mag_term) >= prec_error*pr(1)) then
394  write(mesg, '(ES13.5)') max_mag_term
395  call mom_error(fatal,"Overflow in reproducing_sum(_3d) conversion of "//trim(mesg))
396  endif
397  if (overflow_error) call mom_error(fatal, "Overflow in reproducing_sum(_3d).")
398  endif
399 
400  call sum_across_pes(ints_sum, ni)
401 
402  call regularize_ints(ints_sum)
403  sum = ints_to_real(ints_sum)
404 
405  if (present(efp_sum)) efp_sum%v(:) = ints_sum(:)
406 
407  if (debug) then
408  write(mesg,'("3d RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
409  call mom_mesg(mesg, 3)
410  endif
411  endif
412 

References carry_overflow(), debug, increment_ints(), increment_ints_faster(), ints_to_real(), max_count_prec, mom_error_handler::mom_error(), mom_error_handler::mom_mesg(), nan_error, ni, overflow_error, pr, real_to_ints(), and regularize_ints().

Here is the call graph for this function:

◆ reset_efp_overflow_error()

subroutine, public mom_coms::reset_efp_overflow_error ( )

Reset the module's error flag to false.

Definition at line 609 of file MOM_coms.F90.

609  overflow_error = .false.

References overflow_error.

Referenced by mom_spatial_means::global_i_mean(), and mom_spatial_means::global_j_mean().

Here is the caller graph for this function:

Variable Documentation

◆ debug

logical mom_coms::debug = .false.
private

Making this true enables debugging output.

Definition at line 51 of file MOM_coms.F90.

51 logical :: debug = .false. !< Making this true enables debugging output.

Referenced by reproducing_sum_2d(), and reproducing_sum_3d().

◆ i_pr

real, dimension(ni), parameter mom_coms::i_pr = (/ 1.0/r_prec**2, 1.0/r_prec, 1.0, r_prec, r_prec**2, r_prec**3 /)
private

An array of the inverse of the real precision of each of the integers.

Definition at line 41 of file MOM_coms.F90.

41 real, parameter, dimension(ni) :: &
42  I_pr = (/ 1.0/r_prec**2, 1.0/r_prec, 1.0, r_prec, r_prec**2, r_prec**3 /)

Referenced by increment_ints_faster(), real_to_ints(), and reproducing_sum_2d().

◆ i_prec

real, parameter mom_coms::i_prec =1.0/(2.0**46)
private

The inverse of prec.

Definition at line 30 of file MOM_coms.F90.

30 real, parameter :: I_prec=1.0/(2.0**46) !< The inverse of prec.

Referenced by carry_overflow(), and regularize_ints().

◆ max_count_prec

integer, parameter mom_coms::max_count_prec =2**(63-46)-1
private

The number of values that can be added together with the current value of prec before there will be roundoff problems.

Definition at line 31 of file MOM_coms.F90.

31 integer, parameter :: max_count_prec=2**(63-46)-1

Referenced by efp_list_sum_across_pes(), reproducing_sum_2d(), and reproducing_sum_3d().

◆ max_efp_float

real, parameter mom_coms::max_efp_float = pr(1) * (2.**63 - 1.)
private

The largest float with an EFP representation. NOTE: Only the first bin can exceed precision, but is bounded by the largest signed integer.

Definition at line 44 of file MOM_coms.F90.

44 real, parameter :: max_efp_float = pr(1) * (2.**63 - 1.)

Referenced by increment_ints_faster().

◆ nan_error

logical mom_coms::nan_error = .false.
private

This becomes true if a NaN is encountered.

Definition at line 50 of file MOM_coms.F90.

50 logical :: NaN_error = .false. !< This becomes true if a NaN is encountered.

Referenced by increment_ints_faster(), real_to_ints(), reproducing_sum_2d(), and reproducing_sum_3d().

◆ ni

integer, parameter mom_coms::ni =6
private

The number of long integers to use to represent a real number.

Definition at line 36 of file MOM_coms.F90.

36 integer, parameter :: ni=6 !< The number of long integers to use to represent

Referenced by carry_overflow(), efp_assign(), efp_list_sum_across_pes(), efp_minus(), increment_ints(), increment_ints_faster(), ints_to_real(), real_to_ints(), regularize_ints(), reproducing_sum_2d(), and reproducing_sum_3d().

◆ overflow_error

logical mom_coms::overflow_error = .false.
private

This becomes true if an overflow is encountered.

Definition at line 49 of file MOM_coms.F90.

49 logical :: overflow_error = .false. !< This becomes true if an overflow is encountered.

Referenced by carry_overflow(), efp_list_sum_across_pes(), increment_ints(), increment_ints_faster(), query_efp_overflow_error(), reproducing_sum_2d(), reproducing_sum_3d(), and reset_efp_overflow_error().

◆ pr

real, dimension(ni), parameter mom_coms::pr = (/ r_prec**2, r_prec, 1.0, 1.0/r_prec, 1.0/r_prec**2, 1.0/r_prec**3 /)
private

An array of the real precision of each of the integers.

Definition at line 38 of file MOM_coms.F90.

38 real, parameter, dimension(ni) :: &
39  pr = (/ r_prec**2, r_prec, 1.0, 1.0/r_prec, 1.0/r_prec**2, 1.0/r_prec**3 /)

Referenced by increment_ints_faster(), ints_to_real(), real_to_ints(), reproducing_sum_2d(), and reproducing_sum_3d().

◆ prec

integer(kind=8), parameter mom_coms::prec =2_8**46
private

The precision of each integer.

Definition at line 28 of file MOM_coms.F90.

28 integer(kind=8), parameter :: prec=2_8**46 !< The precision of each integer.

Referenced by carry_overflow(), increment_ints(), real_to_ints(), and regularize_ints().

◆ r_prec

real, parameter mom_coms::r_prec =2.0**46
private

A real version of prec.

Definition at line 29 of file MOM_coms.F90.

29 real, parameter :: r_prec=2.0**46 !< A real version of prec.