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

Detailed Description

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

Definition at line 54 of file MOM_coms.F90.

Private functions

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...
 

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 
)
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 

◆ 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,
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 
The documentation for this interface was generated from the following file: