Skip to content

Commit

Permalink
Merge branch 'maltrud/ocean/sss-restoring-reprosum' (PR #6456)
Browse files Browse the repository at this point in the history
Replace global sum with full-precision reproducible sum

The global sum in the surface salinity restoring routine was forced to
be BFB using a truncation calculation that resulted in reduced precision
and some unexplained behavior. This PR replaces this calculation with an
mpas standard reproducible sum with full precision.

Testing was done using a G-case with the IcoswISC30E3r5 grid for one
year, then calculating the global average SSS restoring tendency term,
which should be effectively zero. Also confirmed that the new code
resulted in BFB solution when processor layout was changed from 10 nodes
on chrysalis to 15 nodes.

Fixes #6446
[non-BFB] for C- and G-cases with salinity restoring
  • Loading branch information
jonbob committed Jul 10, 2024
2 parents cb73e3a + 7ca7088 commit 429dc36
Showing 1 changed file with 29 additions and 13 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@ module ocn_tracer_surface_restoring
use mpas_timekeeping
use mpas_forcing
use mpas_stream_manager
use mpas_global_sum_mod
use ocn_constants
use ocn_config
use ocn_framework_forcing
Expand Down Expand Up @@ -250,7 +251,11 @@ subroutine ocn_get_surfaceSalinityData( streamManager, &
integer :: ierr

integer, parameter :: nSums = 2
real (kind=RKIND), dimension(nSums) :: reductions, sums
real (kind=RKIND), dimension(nSums) :: reductions

integer, dimension(2) :: indexForReproSum ! min, max indices for 1 dimensional sums
real (kind=RKIND), dimension(:,:), allocatable :: &
localArrayForReproSum

! initialize monthly forcing to be read from file

Expand Down Expand Up @@ -339,6 +344,9 @@ subroutine ocn_get_surfaceSalinityData( streamManager, &
call mpas_pool_get_array(forcingPool, 'landIceMask', landIceMask)
call mpas_pool_get_array(meshPool, 'areaCell', areaCell)

allocate (localArrayForReproSum(nCellsSolve,2))
localArrayForReproSum(:,:) = 0.0_RKIND

! This is not in a threaded region, so no openMP pragmas are needed.
if ( associated(landIceMask)) then

Expand Down Expand Up @@ -416,11 +424,17 @@ subroutine ocn_get_surfaceSalinityData( streamManager, &

endif

! For reproducible sums like mpas_global_sum_mod(), the range must be specified in min/max pairs for each array dimension.
! In this case, the array is getting packed in advance so the min is always 1 but the max can vary.
! Both need to be provided for the interface.
indexForReproSum(1) = 1
indexForReproSum(2) = 0
do iCell=1,nCellsSolve
deltaS = activeTracersSurfaceRestoringValue(indexSalinitySurfaceRestoringValue,iCell)
if (deltaS .ne. 0.0_RKIND) then
sumAreaDeltaS = sumAreaDeltaS + deltaS*areaCell(iCell)
sumArea = sumArea + areaCell(iCell)
indexForReproSum(2) = indexForReproSum(2) + 1
localArrayForReproSum(indexForReproSum(2),1) = deltaS
localArrayForReproSum(indexForReproSum(2),2) = areaCell(iCell)
endif
enddo

Expand All @@ -429,16 +443,18 @@ subroutine ocn_get_surfaceSalinityData( streamManager, &

! Global sum to subtract global mean of deltaS
dminfo = domain % dminfo
sums(1) = sumAreaDeltaS
sums(2) = sumArea
call mpas_dmpar_sum_real_array(dminfo, nSums, sums(1:nSums), reductions(1:nSums))
sumAreaDeltaSGlobal = reductions(1)
sumAreaGlobal = reductions(2)
avgDeltaS1 = sumAreaDeltaSGlobal/(sumAreaGlobal + 1.e-20_RKIND)
!LPV: There is no guarantee that a global sum of deltaS is BFB across
!decompositions, the next line rounds each value to 10 decimal places to
!ensure BFB results
avgDeltaS = float (int(avgDeltaS1 * 1.0E10_RKIND + 0.5_RKIND)) / 1.0E10_RKIND
! first do 1 field sum for area
reductions(1) = mpas_global_sum(localArrayForReproSum(:,2), &
domain%dminfo%comm, indexForReproSum)
! now do sum of 2 mupltiplied fields for area*deltaS
reductions(2) = mpas_global_sum(localArrayForReproSum(:,1), localArrayForReproSum(:,2), &
domain%dminfo%comm, indexForReproSum)

sumAreaDeltaSGlobal = reductions(2)
sumAreaGlobal = reductions(1)
avgDeltaS = sumAreaDeltaSGlobal/(sumAreaGlobal + 1.e-20_RKIND)

deallocate (localArrayForReproSum)

block => domain % blocklist
do while (associated(block))
Expand Down

0 comments on commit 429dc36

Please sign in to comment.