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(ni)	real_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.
subroutine, public	reset_efp_overflow_error ()
	Reset the module's error flag to false.
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.

Variables
integer(kind=8), parameter	prec =2_8**46
	The precision of each integer.
real, parameter	r_prec =2.0**46
	A real version of prec.
real, parameter	i_prec =1.0/(2.0**46)
	The inverse of prec.
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.
integer, parameter	ni =6
	The number of long integers to use to represent a real number.
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.
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.
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.
logical	overflow_error = .false.
	This becomes true if an overflow is encountered.
logical	nan_error = .false.
	This becomes true if a NaN is encountered.
logical	debug = .false.
	Making this true enables debugging output.

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_sum	The array of EFP integers being modified by carries, but without changing value.
[in]	prec_error	The 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.

  integer(kind=8), dimension(ni), intent(inout) :: int_sum  !< The array of EFP integers being
                                              !! modified by carries, but without changing value.
  integer(kind=8),                intent(in)    :: prec_error  !< The 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.
 
  ! This subroutine handles carrying of the overflow.
  integer :: i, num_carry
 
  do i=ni,2,-1 ; if (abs(int_sum(i)) >= prec) then
    num_carry = int(int_sum(i) * i_prec)
    int_sum(i) = int_sum(i) - num_carry*prec
    int_sum(i-1) = int_sum(i-1) + num_carry
  endif ; enddo
  if (abs(int_sum(1)) > prec_error) then
    overflow_error = .true.
  endif
 

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]	efp1	The recipient extended fixed point number
[in]	efp2	The source extended fixed point number

Definition at line 638 of file MOM_coms.F90.

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

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]	efps	The list of extended fixed point numbers
[in]	nval	The number of values being summed.
[out]	errors	A list of error flags for each sum

Definition at line 698 of file MOM_coms.F90.

  type(EFP_type), dimension(:), &
              intent(inout) :: EFPs   !< The list of extended fixed point numbers
                                      !! being summed across PEs.
  integer,    intent(in)    :: nval   !< The number of values being summed.
  logical, dimension(:), &
    optional, intent(out)   :: errors !< A list of error flags for each sum
 
  !   This subroutine does a sum across PEs of a list of EFP variables,
  ! returning the sums in place, with all overflows carried.
 
  integer(kind=8), dimension(ni,nval) :: ints
  integer(kind=8) :: prec_error
  logical :: error_found
  character(len=256) :: mesg
  integer :: i, n
 
  if (num_pes() > max_count_prec) call mom_error(fatal, &
    "reproducing_sum: Too many processors are being used for the value of "//&
    "prec.  Reduce prec to (2^63-1)/num_PEs.")
 
  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
  ! overflow_error is an overflow error flag for the whole module.
  overflow_error = .false. ; error_found = .false.
 
  do i=1,nval ; do n=1,ni ; ints(n,i) = efps(i)%v(n) ; enddo ; enddo
 
  call sum_across_pes(ints(:,:), ni*nval)
 
  if (present(errors)) errors(:) = .false.
  do i=1,nval
    overflow_error = .false.
    call carry_overflow(ints(:,i), prec_error)
    do n=1,ni ; efps(i)%v(n) = ints(n,i) ; enddo
    if (present(errors)) errors(i) = overflow_error
    if (overflow_error) then
      write (mesg,'("EFP_list_sum_across_PEs error at ",i6," val was ",ES12.6, ", prec_error = ",ES12.6)') &
             i, efp_to_real(efps(i)), real(prec_error)
      call mom_error(warning, mesg)
    endif
    error_found = error_found .or. overflow_error
  enddo
  if (error_found .and. .not.(present(errors))) then
    call mom_error(fatal, "Overflow in EFP_list_sum_across_PEs.")
  endif
 

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]	efp1	The first extended fixed point number
[in]	efp2	The extended fixed point number being subtracted from the first extended fixed point number

Definition at line 625 of file MOM_coms.F90.

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

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]	efp1	The first extended fixed point number
[in]	efp2	The second extended fixed point number

Definition at line 614 of file MOM_coms.F90.

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

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]	efp1	The first extended fixed point number
[in]	efp2	The 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.

  type(EFP_type), intent(in) :: EFP1  !< The first extended fixed point number
  type(EFP_type), intent(in) :: EFP2  !< The extended fixed point number being
                        !! subtracted from the first extended fixed point number
  real :: EFP_real_diff !< The real result
 
  type(EFP_type)             :: EFP_diff
 
  efp_diff = efp1 - efp2
  efp_real_diff = efp_to_real(efp_diff)
 

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]

efp1

The extended fixed point number being converted

Definition at line 650 of file MOM_coms.F90.

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

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_sum	The array of EFP integers being incremented
[in]	int2	The array of EFP integers being added
[in]	prec_error	The 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.

  integer(kind=8), dimension(ni), intent(inout) :: int_sum !< The array of EFP integers being incremented
  integer(kind=8), dimension(ni), intent(in)    :: int2    !< The array of EFP integers being added
  integer(kind=8), optional,      intent(in)    :: prec_error !< The 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.
 
  ! This subroutine increments a number with another, both using the integer
  ! representation in real_to_ints.
  integer :: i
 
  do i=ni,2,-1
    int_sum(i) = int_sum(i) + int2(i)
    ! Carry the local overflow.
    if (int_sum(i) > prec) then
      int_sum(i) = int_sum(i) - prec
      int_sum(i-1) = int_sum(i-1) + 1
    elseif (int_sum(i) < -prec) then
      int_sum(i) = int_sum(i) + prec
      int_sum(i-1) = int_sum(i-1) - 1
    endif
  enddo
  int_sum(1) = int_sum(1) + int2(1)
  if (present(prec_error)) then
    if (abs(int_sum(1)) > prec_error) overflow_error = .true.
  else
    if (abs(int_sum(1)) > prec) overflow_error = .true.
  endif
 

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_sum	The array of EFP integers being incremented
[in]	r	The real number being added.
[in,out]	max_mag_term	A running maximum magnitude of the r's.

Definition at line 506 of file MOM_coms.F90.

  integer(kind=8), dimension(ni), intent(inout) :: int_sum  !< The array of EFP integers being incremented
  real,                           intent(in)    :: r        !< The real number being added.
  real,                           intent(inout) :: max_mag_term !< A running maximum magnitude of the r's.
 
  ! This subroutine increments a number with another, both using the integer
  ! representation in real_to_ints, but without doing any carrying of overflow.
  ! The entire operation is embedded in a single call for greater speed.
  real :: rs
  integer(kind=8) :: ival
  integer :: sgn, i
 
  if ((r >= 1e30) .eqv. (r < 1e30)) then ; nan_error = .true. ; return ; endif
  sgn = 1 ; if (r<0.0) sgn = -1
  rs = abs(r)
  if (rs > abs(max_mag_term)) max_mag_term = r
 
  ! Abort if the number has no EFP representation
  if (rs > max_efp_float) then
    overflow_error = .true.
    return
  endif
 
  do i=1,ni
    ival = int(rs*i_pr(i), 8)
    rs = rs - ival*pr(i)
    int_sum(i) = int_sum(i) + sgn*ival
  enddo
 

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]

ints

The array of EFP integers

Definition at line 459 of file MOM_coms.F90.

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

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]	val	The real number being converted
[in,out]	overflow	Returns 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.

  real,              intent(in)    :: val !< The real number being converted
  logical, optional, intent(inout) :: overflow !< Returns true if the conversion is being
                                          !! done on a value that is too large to be represented
  type(EFP_type) :: real_to_EFP
 
  logical :: over
  character(len=80) :: mesg
 
  if (present(overflow)) then
    real_to_efp%v(:) = real_to_ints(val, overflow=overflow)
  else
    over = .false.
    real_to_efp%v(:) = real_to_ints(val, overflow=over)
    if (over) then
      write(mesg, '(ES13.5)') val
      call mom_error(fatal,"Overflow in real_to_EFP conversion of "//trim(mesg))
    endif
  endif
 

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]	r	The real number being converted
[in]	prec_error	The 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]	overflow	Returns 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.

  real,                      intent(in) :: r  !< The real number being converted
  integer(kind=8), optional, intent(in) :: prec_error  !< The 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.
  logical,         optional, intent(inout) :: overflow !< Returns true if the conversion is being
                                              !! done on a value that is too large to be represented
  integer(kind=8), dimension(ni)  :: ints
  !   This subroutine converts a real number to an equivalent representation
  ! using several long integers.
 
  real :: rs
  character(len=80) :: mesg
  integer(kind=8) :: ival, prec_err
  integer :: sgn, i
 
  prec_err = prec ; if (present(prec_error)) prec_err = prec_error
  ints(:) = 0_8
  if ((r >= 1e30) .eqv. (r < 1e30)) then ; nan_error = .true. ; return ; endif
 
  sgn = 1 ; if (r<0.0) sgn = -1
  rs = abs(r)
 
  if (present(overflow)) then
    if (.not.(rs < prec_err*pr(1))) overflow = .true.
    if ((r >= 1e30) .eqv. (r < 1e30)) overflow = .true.
  elseif (.not.(rs < prec_err*pr(1))) then
    write(mesg, '(ES13.5)') r
    call mom_error(fatal,"Overflow in real_to_ints conversion of "//trim(mesg))
  endif
 
  do i=1,ni
    ival = int(rs*i_pr(i), 8)
    rs = rs - ival*pr(i)
    ints(i) = sgn*ival
  enddo
 

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_sum

The array of integers being modified to take a

Definition at line 562 of file MOM_coms.F90.

  integer(kind=8), dimension(ni), &
    intent(inout) :: int_sum !< The array of integers being modified to take a
                             !! regular form with all integers of the same sign,
                             !! but without changing value.
 
  ! This subroutine carries the overflow, and then makes sure that
  ! all integers are of the same sign as the overall value.
  logical :: positive
  integer :: i, num_carry
 
  do i=ni,2,-1 ; if (abs(int_sum(i)) >= prec) then
    num_carry = int(int_sum(i) * i_prec)
    int_sum(i) = int_sum(i) - num_carry*prec
    int_sum(i-1) = int_sum(i-1) + num_carry
  endif ; enddo
 
  ! Determine the sign of the final number.
  positive = .true.
  do i=1,ni
    if (abs(int_sum(i)) > 0) then
      if (int_sum(i) < 0) positive = .false.
      exit
    endif
  enddo
 
  if (positive) then
    do i=ni,2,-1 ; if (int_sum(i) < 0) then
      int_sum(i) = int_sum(i) + prec
      int_sum(i-1) = int_sum(i-1) - 1
    endif ; enddo
  else
    do i=ni,2,-1 ; if (int_sum(i) > 0) then
      int_sum(i) = int_sum(i) - prec
      int_sum(i-1) = int_sum(i-1) + 1
    endif ; enddo
  endif
 

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]	array	The array to be summed
[in]	isr	The starting i-index of the sum, noting that the array indices starts at 1
[in]	ier	The ending i-index of the sum, noting that the array indices starts at 1
[in]	jsr	The starting j-index of the sum, noting that the array indices starts at 1
[in]	jer	The ending j-index of the sum, noting that the array indices starts at 1
[out]	efp_sum	The result in extended fixed point format
[in]	reproducing	If present and false, do the sum using the naive non-reproducing approach
[in]	overflow_check	If 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]	err	If 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.

  real, dimension(:,:),     intent(in)  :: array    !< The array to be summed
  integer,        optional, intent(in)  :: isr     !< The starting i-index of the sum, noting
                                                   !! that the array indices starts at 1
  integer,        optional, intent(in)  :: ier     !< The ending i-index of the sum, noting
                                                   !! that the array indices starts at 1
  integer,        optional, intent(in)  :: jsr     !< The starting j-index of the sum, noting
                                                   !! that the array indices starts at 1
  integer,        optional, intent(in)  :: jer     !< The ending j-index of the sum, noting
                                                   !! that the array indices starts at 1
  type(EFP_type), optional, intent(out) :: EFP_sum  !< The result in extended fixed point format
  logical,        optional, intent(in)  :: reproducing !< If present and false, do the sum
                                                !! using the naive non-reproducing approach
  logical,        optional, intent(in)  :: overflow_check !< If 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
  integer,        optional, intent(out) :: err  !< If present, return an error code instead of
                                                !! triggering any fatal errors directly from
                                                !! this routine.
  real                                  :: sum  !< Result
 
  !   This subroutine uses a conversion to an integer representation
  ! of real numbers to give order-invariant sums that will reproduce
  ! across PE count.  This idea comes from R. Hallberg and A. Adcroft.
 
  integer(kind=8), dimension(ni)  :: ints_sum
  integer(kind=8) :: ival, prec_error
  real    :: rsum(1), rs
  real    :: max_mag_term
  logical :: repro, over_check
  character(len=256) :: mesg
  integer :: i, j, n, is, ie, js, je, sgn
 
  if (num_pes() > max_count_prec) call mom_error(fatal, &
    "reproducing_sum: Too many processors are being used for the value of "//&
    "prec.  Reduce prec to (2^63-1)/num_PEs.")
 
  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
 
  is = 1 ; ie = size(array,1) ; js = 1 ; je = size(array,2 )
  if (present(isr)) then
    if (isr < is) call mom_error(fatal, &
      "Value of isr too small in reproducing_sum_2d.")
    is = isr
  endif
  if (present(ier)) then
    if (ier > ie) call mom_error(fatal, &
      "Value of ier too large in reproducing_sum_2d.")
    ie = ier
  endif
  if (present(jsr)) then
    if (jsr < js) call mom_error(fatal, &
      "Value of jsr too small in reproducing_sum_2d.")
    js = jsr
  endif
  if (present(jer)) then
    if (jer > je) call mom_error(fatal, &
      "Value of jer too large in reproducing_sum_2d.")
    je = jer
  endif
 
  repro = .true. ; if (present(reproducing)) repro = reproducing
  over_check = .true. ; if (present(overflow_check)) over_check = overflow_check
 
  if (repro) then
    overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
    ints_sum(:) = 0
    if (over_check) then
      if ((je+1-js)*(ie+1-is) < max_count_prec) then
        do j=js,je ; do i=is,ie
          call increment_ints_faster(ints_sum, array(i,j), max_mag_term)
        enddo ; enddo
        call carry_overflow(ints_sum, prec_error)
      elseif ((ie+1-is) < max_count_prec) then
        do j=js,je
          do i=is,ie
            call increment_ints_faster(ints_sum, array(i,j), max_mag_term)
          enddo
          call carry_overflow(ints_sum, prec_error)
        enddo
      else
        do j=js,je ; do i=is,ie
          call increment_ints(ints_sum, real_to_ints(array(i,j), prec_error), &
                              prec_error)
        enddo ; enddo
      endif
    else
      do j=js,je ; do i=is,ie
        sgn = 1 ; if (array(i,j)<0.0) sgn = -1
        rs = abs(array(i,j))
        do n=1,ni
          ival = int(rs*i_pr(n), 8)
          rs = rs - ival*pr(n)
          ints_sum(n) = ints_sum(n) + sgn*ival
        enddo
      enddo ; enddo
      call carry_overflow(ints_sum, prec_error)
    endif
 
    if (present(err)) then
      err = 0
      if (overflow_error) &
        err = err+2
      if (nan_error) &
        err = err+4
      if (err > 0) then ; do n=1,ni ; ints_sum(n) = 0 ; enddo ; endif
    else
      if (nan_error) then
        call mom_error(fatal, "NaN in input field of reproducing_sum(_2d).")
      endif
      if (abs(max_mag_term) >= prec_error*pr(1)) then
        write(mesg, '(ES13.5)') max_mag_term
        call mom_error(fatal,"Overflow in reproducing_sum(_2d) conversion of "//trim(mesg))
      endif
      if (overflow_error) then
        call mom_error(fatal, "Overflow in reproducing_sum(_2d).")
      endif
    endif
 
    call sum_across_pes(ints_sum, ni)
 
    call regularize_ints(ints_sum)
    sum = ints_to_real(ints_sum)
  else
    rsum(1) = 0.0
    do j=js,je ; do i=is,ie
      rsum(1) = rsum(1) + array(i,j)
    enddo ; enddo
    call sum_across_pes(rsum,1)
    sum = rsum(1)
 
    if (present(err)) then ; err = 0 ; endif
 
    if (debug .or. present(efp_sum)) then
      overflow_error = .false.
      ints_sum = real_to_ints(sum, prec_error, overflow_error)
      if (overflow_error) then
        if (present(err)) then
          err = err + 2
        else
          write(mesg, '(ES13.5)') sum
          call mom_error(fatal,"Repro_sum_2d: Overflow in real_to_ints conversion of "//trim(mesg))
        endif
      endif
    endif
  endif
 
  if (present(efp_sum)) efp_sum%v(:) = ints_sum(:)
 
  if (debug) then
    write(mesg,'("2d RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
    call mom_mesg(mesg, 3)
  endif
 

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]	array	The array to be summed
[in]	isr	The starting i-index of the sum, noting that the array indices starts at 1
[in]	ier	The ending i-index of the sum, noting that the array indices starts at 1
[in]	jsr	The starting j-index of the sum, noting that the array indices starts at 1
[in]	jer	The ending j-index of the sum, noting that the array indices starts at 1
[out]	sums	The sums by vertical layer
[out]	efp_sum	The result in extended fixed point format
[out]	err	If 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.

  real, dimension(:,:,:),       intent(in)  :: array   !< The array to be summed
  integer,            optional, intent(in)  :: isr     !< The starting i-index of the sum, noting
                                                       !! that the array indices starts at 1
  integer,            optional, intent(in)  :: ier     !< The ending i-index of the sum, noting
                                                       !! that the array indices starts at 1
  integer,            optional, intent(in)  :: jsr     !< The starting j-index of the sum, noting
                                                       !! that the array indices starts at 1
  integer,            optional, intent(in)  :: jer     !< The ending j-index of the sum, noting
                                                       !! that the array indices starts at 1
  real, dimension(:), optional, intent(out) :: sums    !< The sums by vertical layer
  type(EFP_type),     optional, intent(out) :: EFP_sum !< The result in extended fixed point format
  integer,            optional, intent(out) :: err  !< If present, return an error code instead of
                                                    !! triggering any fatal errors directly from
                                                    !! this routine.
  real                                      :: sum  !< Result
 
  !   This subroutine uses a conversion to an integer representation
  ! of real numbers to give order-invariant sums that will reproduce
  ! across PE count.  This idea comes from R. Hallberg and A. Adcroft.
 
  real    :: max_mag_term
  integer(kind=8), dimension(ni)  :: ints_sum
  integer(kind=8), dimension(ni,size(array,3))  :: ints_sums
  integer(kind=8) :: prec_error
  character(len=256) :: mesg
  integer :: i, j, k, is, ie, js, je, ke, isz, jsz, n
 
  if (num_pes() > max_count_prec) call mom_error(fatal, &
    "reproducing_sum: Too many processors are being used for the value of "//&
    "prec.  Reduce prec to (2^63-1)/num_PEs.")
 
  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
  max_mag_term = 0.0
 
  is = 1 ; ie = size(array,1) ; js = 1 ; je = size(array,2) ; ke = size(array,3)
  if (present(isr)) then
    if (isr < is) call mom_error(fatal, &
      "Value of isr too small in reproducing_sum(_3d).")
    is = isr
  endif
  if (present(ier)) then
    if (ier > ie) call mom_error(fatal, &
      "Value of ier too large in reproducing_sum(_3d).")
    ie = ier
  endif
  if (present(jsr)) then
    if (jsr < js) call mom_error(fatal, &
      "Value of jsr too small in reproducing_sum(_3d).")
    js = jsr
  endif
  if (present(jer)) then
    if (jer > je) call mom_error(fatal, &
      "Value of jer too large in reproducing_sum(_3d).")
    je = jer
  endif
  jsz = je+1-js; isz = ie+1-is
 
  if (present(sums)) then
    if (size(sums) > ke) call mom_error(fatal, "Sums is smaller than "//&
      "the vertical extent of array in reproducing_sum(_3d).")
    ints_sums(:,:) = 0
    overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
    if (jsz*isz < max_count_prec) then
      do k=1,ke
        do j=js,je ; do i=is,ie
          call increment_ints_faster(ints_sums(:,k), array(i,j,k), max_mag_term)
        enddo ; enddo
        call carry_overflow(ints_sums(:,k), prec_error)
      enddo
    elseif (isz < max_count_prec) then
      do k=1,ke ; do j=js,je
        do i=is,ie
          call increment_ints_faster(ints_sums(:,k), array(i,j,k), max_mag_term)
        enddo
        call carry_overflow(ints_sums(:,k), prec_error)
      enddo ; enddo
    else
      do k=1,ke ; do j=js,je ; do i=is,ie
        call increment_ints(ints_sums(:,k), &
                            real_to_ints(array(i,j,k), prec_error), prec_error)
      enddo ; enddo ; enddo
    endif
    if (present(err)) then
      err = 0
      if (abs(max_mag_term) >= prec_error*pr(1)) err = err+1
      if (overflow_error) err = err+2
      if (nan_error) err = err+2
      if (err > 0) then ; do k=1,ke ; do n=1,ni ; ints_sums(n,k) = 0 ; enddo ; enddo ; endif
    else
      if (nan_error) call mom_error(fatal, "NaN in input field of reproducing_sum(_3d).")
      if (abs(max_mag_term) >= prec_error*pr(1)) then
        write(mesg, '(ES13.5)') max_mag_term
        call mom_error(fatal,"Overflow in reproducing_sum(_3d) conversion of "//trim(mesg))
      endif
      if (overflow_error) call mom_error(fatal, "Overflow in reproducing_sum(_3d).")
    endif
 
    call sum_across_pes(ints_sums(:,1:ke), ni*ke)
 
    sum = 0.0
    do k=1,ke
      call regularize_ints(ints_sums(:,k))
      sums(k) = ints_to_real(ints_sums(:,k))
      sum = sum + sums(k)
    enddo
 
    if (present(efp_sum)) then
      efp_sum%v(:) = 0
      do k=1,ke ; call increment_ints(efp_sum%v(:), ints_sums(:,k)) ; enddo
    endif
 
    if (debug) then
      do n=1,ni ; ints_sum(n) = 0 ; enddo
      do k=1,ke ; do n=1,ni ; ints_sum(n) = ints_sum(n) + ints_sums(n,k) ; enddo ; enddo
      write(mesg,'("3D RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
      call mom_mesg(mesg, 3)
    endif
  else
    ints_sum(:) = 0
    overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
    if (jsz*isz < max_count_prec) then
      do k=1,ke
        do j=js,je ; do i=is,ie
          call increment_ints_faster(ints_sum, array(i,j,k), max_mag_term)
        enddo ; enddo
        call carry_overflow(ints_sum, prec_error)
      enddo
    elseif (isz < max_count_prec) then
      do k=1,ke ; do j=js,je
        do i=is,ie
          call increment_ints_faster(ints_sum, array(i,j,k), max_mag_term)
        enddo
        call carry_overflow(ints_sum, prec_error)
      enddo ; enddo
    else
      do k=1,ke ; do j=js,je ; do i=is,ie
        call increment_ints(ints_sum, real_to_ints(array(i,j,k), prec_error), &
                            prec_error)
      enddo ; enddo ; enddo
    endif
    if (present(err)) then
      err = 0
      if (abs(max_mag_term) >= prec_error*pr(1)) err = err+1
      if (overflow_error) err = err+2
      if (nan_error) err = err+2
      if (err > 0) then ; do n=1,ni ; ints_sum(n) = 0 ; enddo ; endif
    else
      if (nan_error) call mom_error(fatal, "NaN in input field of reproducing_sum(_3d).")
      if (abs(max_mag_term) >= prec_error*pr(1)) then
        write(mesg, '(ES13.5)') max_mag_term
        call mom_error(fatal,"Overflow in reproducing_sum(_3d) conversion of "//trim(mesg))
      endif
      if (overflow_error) call mom_error(fatal, "Overflow in reproducing_sum(_3d).")
    endif
 
    call sum_across_pes(ints_sum, ni)
 
    call regularize_ints(ints_sum)
    sum = ints_to_real(ints_sum)
 
    if (present(efp_sum)) efp_sum%v(:) = ints_sum(:)
 
    if (debug) then
      write(mesg,'("3d RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
      call mom_mesg(mesg, 3)
    endif
  endif