MOM6
Private functions | List of all members
mom_coms::reproducing_sum Interface Reference

Detailed Description

Find an accurate and order-invariant sum of a distributed 2d or 3d field.

Definition at line 53 of file MOM_coms.F90.

Private functions

real function reproducing_sum_2d (array, isr, ier, jsr, jer, EFP_sum, reproducing, overflow_check, err, only_on_PE)
 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, EFP_lay_sums, err, only_on_PE)
 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...
 

Detailed Description

Find an accurate and order-invariant sum of a distributed 2d or 3d field.

Definition at line 53 of file MOM_coms.F90.

Functions and subroutines

◆ reproducing_sum_2d()

real function mom_coms::reproducing_sum::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,
logical, intent(in), optional  only_on_PE 
)
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.
[in]only_on_peIf present and true, do not do the sum across processors, only reporting the local sum
Returns
Result

Definition at line 220 of file MOM_coms.F90.

220  real, dimension(:,:), intent(in) :: array !< The array to be summed
221  integer, optional, intent(in) :: isr !< The starting i-index of the sum, noting
222  !! that the array indices starts at 1
223  integer, optional, intent(in) :: ier !< The ending i-index of the sum, noting
224  !! that the array indices starts at 1
225  integer, optional, intent(in) :: jsr !< The starting j-index of the sum, noting
226  !! that the array indices starts at 1
227  integer, optional, intent(in) :: jer !< The ending j-index of the sum, noting
228  !! that the array indices starts at 1
229  type(EFP_type), optional, intent(out) :: EFP_sum !< The result in extended fixed point format
230  logical, optional, intent(in) :: reproducing !< If present and false, do the sum
231  !! using the naive non-reproducing approach
232  logical, optional, intent(in) :: overflow_check !< If present and false, disable
233  !! checking for overflows in incremental results.
234  !! This can speed up calculations if the number
235  !! of values being summed is small enough
236  integer, optional, intent(out) :: err !< If present, return an error code instead of
237  !! triggering any fatal errors directly from
238  !! this routine.
239  logical, optional, intent(in) :: only_on_PE !< If present and true, do not do the sum
240  !! across processors, only reporting the local sum
241  real :: sum !< Result
242 
243  ! This subroutine uses a conversion to an integer representation
244  ! of real numbers to give order-invariant sums that will reproduce
245  ! across PE count. This idea comes from R. Hallberg and A. Adcroft.
246 
247  integer(kind=8), dimension(ni) :: ints_sum
248  integer(kind=8) :: prec_error
249  real :: rsum(1), rs
250  logical :: repro, do_sum_across_PEs
251  character(len=256) :: mesg
252  type(EFP_type) :: EFP_val ! An extended fixed point version of the sum
253  integer :: i, j, n, is, ie, js, je
254 
255  if (num_pes() > max_count_prec) call mom_error(fatal, &
256  "reproducing_sum: Too many processors are being used for the value of "//&
257  "prec. Reduce prec to (2^63-1)/num_PEs.")
258 
259  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
260 
261  is = 1 ; ie = size(array,1) ; js = 1 ; je = size(array,2 )
262  if (present(isr)) then
263  if (isr < is) call mom_error(fatal, "Value of isr too small in reproducing_sum_2d.")
264  is = isr
265  endif
266  if (present(ier)) then
267  if (ier > ie) call mom_error(fatal, "Value of ier too large in reproducing_sum_2d.")
268  ie = ier
269  endif
270  if (present(jsr)) then
271  if (jsr < js) call mom_error(fatal, "Value of jsr too small in reproducing_sum_2d.")
272  js = jsr
273  endif
274  if (present(jer)) then
275  if (jer > je) call mom_error(fatal, "Value of jer too large in reproducing_sum_2d.")
276  je = jer
277  endif
278 
279  repro = .true. ; if (present(reproducing)) repro = reproducing
280  do_sum_across_pes = .true. ; if (present(only_on_pe)) do_sum_across_pes = .not.only_on_pe
281 
282  if (repro) then
283  efp_val = reproducing_efp_sum_2d(array, isr, ier, jsr, jer, overflow_check, err, only_on_pe)
284  sum = ints_to_real(efp_val%v)
285  if (present(efp_sum)) efp_sum = efp_val
286  if (debug) ints_sum(:) = efp_sum%v(:)
287  else
288  rsum(1) = 0.0
289  do j=js,je ; do i=is,ie
290  rsum(1) = rsum(1) + array(i,j)
291  enddo ; enddo
292  if (do_sum_across_pes) call sum_across_pes(rsum,1)
293  sum = rsum(1)
294 
295  if (present(err)) then ; err = 0 ; endif
296 
297  if (debug .or. present(efp_sum)) then
298  overflow_error = .false.
299  ints_sum = real_to_ints(sum, prec_error, overflow_error)
300  if (overflow_error) then
301  if (present(err)) then
302  err = err + 2
303  else
304  write(mesg, '(ES13.5)') sum
305  call mom_error(fatal,"Repro_sum_2d: Overflow in real_to_ints conversion of "//trim(mesg))
306  endif
307  endif
308  endif
309  if (present(efp_sum)) efp_sum%v(:) = ints_sum(:)
310  endif
311 
312  if (debug) then
313  write(mesg,'("2d RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
314  call mom_mesg(mesg, 3)
315  endif
316 

◆ reproducing_sum_3d()

real function mom_coms::reproducing_sum::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,
type(efp_type), dimension(:), intent(out), optional  EFP_lay_sums,
integer, intent(out), optional  err,
logical, intent(in), optional  only_on_PE 
)
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]efp_lay_sumsThe sums by vertical layer in EFP format
[out]errIf present, return an error code instead of triggering any fatal errors directly from this routine.
[in]only_on_peIf present and true, do not do the sum across processors, only reporting the local sum
Returns
Result

Definition at line 325 of file MOM_coms.F90.

325  real, dimension(:,:,:), intent(in) :: array !< The array to be summed
326  integer, optional, intent(in) :: isr !< The starting i-index of the sum, noting
327  !! that the array indices starts at 1
328  integer, optional, intent(in) :: ier !< The ending i-index of the sum, noting
329  !! that the array indices starts at 1
330  integer, optional, intent(in) :: jsr !< The starting j-index of the sum, noting
331  !! that the array indices starts at 1
332  integer, optional, intent(in) :: jer !< The ending j-index of the sum, noting
333  !! that the array indices starts at 1
334  real, dimension(:), optional, intent(out) :: sums !< The sums by vertical layer
335  type(EFP_type), optional, intent(out) :: EFP_sum !< The result in extended fixed point format
336  type(EFP_type), dimension(:), &
337  optional, intent(out) :: EFP_lay_sums !< The sums by vertical layer in EFP format
338  integer, optional, intent(out) :: err !< If present, return an error code instead of
339  !! triggering any fatal errors directly from
340  !! this routine.
341  logical, optional, intent(in) :: only_on_PE !< If present and true, do not do the sum
342  !! across processors, only reporting the local sum
343  real :: sum !< Result
344 
345  ! This subroutine uses a conversion to an integer representation
346  ! of real numbers to give order-invariant sums that will reproduce
347  ! across PE count. This idea comes from R. Hallberg and A. Adcroft.
348 
349  real :: val, max_mag_term
350  integer(kind=8), dimension(ni) :: ints_sum
351  integer(kind=8), dimension(ni,size(array,3)) :: ints_sums
352  integer(kind=8) :: prec_error
353  character(len=256) :: mesg
354  logical :: do_sum_across_PEs
355  integer :: i, j, k, is, ie, js, je, ke, isz, jsz, n
356 
357  if (num_pes() > max_count_prec) call mom_error(fatal, &
358  "reproducing_sum: Too many processors are being used for the value of "//&
359  "prec. Reduce prec to (2^63-1)/num_PEs.")
360 
361  prec_error = (2_8**62 + (2_8**62 - 1)) / num_pes()
362  max_mag_term = 0.0
363 
364  is = 1 ; ie = size(array,1) ; js = 1 ; je = size(array,2) ; ke = size(array,3)
365  if (present(isr)) then
366  if (isr < is) call mom_error(fatal, "Value of isr too small in reproducing_sum(_3d).")
367  is = isr
368  endif
369  if (present(ier)) then
370  if (ier > ie) call mom_error(fatal, "Value of ier too large in reproducing_sum(_3d).")
371  ie = ier
372  endif
373  if (present(jsr)) then
374  if (jsr < js) call mom_error(fatal, "Value of jsr too small in reproducing_sum(_3d).")
375  js = jsr
376  endif
377  if (present(jer)) then
378  if (jer > je) call mom_error(fatal, "Value of jer too large in reproducing_sum(_3d).")
379  je = jer
380  endif
381  jsz = je+1-js; isz = ie+1-is
382 
383  do_sum_across_pes = .true. ; if (present(only_on_pe)) do_sum_across_pes = .not.only_on_pe
384 
385  if (present(sums) .or. present(efp_lay_sums)) then
386  if (present(sums)) then ; if (size(sums) < ke) then
387  call mom_error(fatal, "Sums is smaller than the vertical extent of array in reproducing_sum(_3d).")
388  endif ; endif
389  if (present(efp_lay_sums)) then ; if (size(efp_lay_sums) < ke) then
390  call mom_error(fatal, "Sums is smaller than the vertical extent of array in reproducing_sum(_3d).")
391  endif ; endif
392  ints_sums(:,:) = 0
393  overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
394  if (jsz*isz < max_count_prec) then
395  do k=1,ke
396  do j=js,je ; do i=is,ie
397  call increment_ints_faster(ints_sums(:,k), array(i,j,k), max_mag_term)
398  enddo ; enddo
399  call carry_overflow(ints_sums(:,k), prec_error)
400  enddo
401  elseif (isz < max_count_prec) then
402  do k=1,ke ; do j=js,je
403  do i=is,ie
404  call increment_ints_faster(ints_sums(:,k), array(i,j,k), max_mag_term)
405  enddo
406  call carry_overflow(ints_sums(:,k), prec_error)
407  enddo ; enddo
408  else
409  do k=1,ke ; do j=js,je ; do i=is,ie
410  call increment_ints(ints_sums(:,k), &
411  real_to_ints(array(i,j,k), prec_error), prec_error)
412  enddo ; enddo ; enddo
413  endif
414  if (present(err)) then
415  err = 0
416  if (abs(max_mag_term) >= prec_error*pr(1)) err = err+1
417  if (overflow_error) err = err+2
418  if (nan_error) err = err+2
419  if (err > 0) then ; do k=1,ke ; do n=1,ni ; ints_sums(n,k) = 0 ; enddo ; enddo ; endif
420  else
421  if (nan_error) call mom_error(fatal, "NaN in input field of reproducing_sum(_3d).")
422  if (abs(max_mag_term) >= prec_error*pr(1)) then
423  write(mesg, '(ES13.5)') max_mag_term
424  call mom_error(fatal,"Overflow in reproducing_sum(_3d) conversion of "//trim(mesg))
425  endif
426  if (overflow_error) call mom_error(fatal, "Overflow in reproducing_sum(_3d).")
427  endif
428 
429  if (do_sum_across_pes) call sum_across_pes(ints_sums(:,1:ke), ni*ke)
430 
431  sum = 0.0
432  do k=1,ke
433  call regularize_ints(ints_sums(:,k))
434  val = ints_to_real(ints_sums(:,k))
435  if (present(sums)) sums(k) = val
436  sum = sum + val
437  enddo
438  if (present(efp_lay_sums)) then ; do k=1,ke
439  efp_lay_sums(k)%v(:) = ints_sums(:,k)
440  enddo ; endif
441 
442  if (present(efp_sum)) then
443  efp_sum%v(:) = 0
444  do k=1,ke ; call increment_ints(efp_sum%v(:), ints_sums(:,k)) ; enddo
445  endif
446 
447  if (debug) then
448  do n=1,ni ; ints_sum(n) = 0 ; enddo
449  do k=1,ke ; do n=1,ni ; ints_sum(n) = ints_sum(n) + ints_sums(n,k) ; enddo ; enddo
450  write(mesg,'("3D RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
451  call mom_mesg(mesg, 3)
452  endif
453  else
454  ints_sum(:) = 0
455  overflow_error = .false. ; nan_error = .false. ; max_mag_term = 0.0
456  if (jsz*isz < max_count_prec) then
457  do k=1,ke
458  do j=js,je ; do i=is,ie
459  call increment_ints_faster(ints_sum, array(i,j,k), max_mag_term)
460  enddo ; enddo
461  call carry_overflow(ints_sum, prec_error)
462  enddo
463  elseif (isz < max_count_prec) then
464  do k=1,ke ; do j=js,je
465  do i=is,ie
466  call increment_ints_faster(ints_sum, array(i,j,k), max_mag_term)
467  enddo
468  call carry_overflow(ints_sum, prec_error)
469  enddo ; enddo
470  else
471  do k=1,ke ; do j=js,je ; do i=is,ie
472  call increment_ints(ints_sum, real_to_ints(array(i,j,k), prec_error), &
473  prec_error)
474  enddo ; enddo ; enddo
475  endif
476  if (present(err)) then
477  err = 0
478  if (abs(max_mag_term) >= prec_error*pr(1)) err = err+1
479  if (overflow_error) err = err+2
480  if (nan_error) err = err+2
481  if (err > 0) then ; do n=1,ni ; ints_sum(n) = 0 ; enddo ; endif
482  else
483  if (nan_error) call mom_error(fatal, "NaN in input field of reproducing_sum(_3d).")
484  if (abs(max_mag_term) >= prec_error*pr(1)) then
485  write(mesg, '(ES13.5)') max_mag_term
486  call mom_error(fatal,"Overflow in reproducing_sum(_3d) conversion of "//trim(mesg))
487  endif
488  if (overflow_error) call mom_error(fatal, "Overflow in reproducing_sum(_3d).")
489  endif
490 
491  if (do_sum_across_pes) call sum_across_pes(ints_sum, ni)
492 
493  call regularize_ints(ints_sum)
494  sum = ints_to_real(ints_sum)
495 
496  if (present(efp_sum)) efp_sum%v(:) = ints_sum(:)
497 
498  if (debug) then
499  write(mesg,'("3d RS: ", ES24.16, 6 Z17.16)') sum, ints_sum(1:ni)
500  call mom_mesg(mesg, 3)
501  endif
502  endif
503 
The documentation for this interface was generated from the following file: