diff --git a/cime_config/config_component.xml b/cime_config/config_component.xml index bad82b988..a1e9e681b 100644 --- a/cime_config/config_component.xml +++ b/cime_config/config_component.xml @@ -1107,6 +1107,21 @@ pes or cores used relative to MAX_MPITASKS_PER_NODE for accounting (0 means TOTALPES is valid) + + logical + TRUE,FALSE + FALSE + run_flags + env_run.xml + + If TRUE, perform some calculations (in CMEPS and elsewhere) such as global sums + using a less efficient algorithm that is independent of processor count. This allows + for Bit-For-Bit (BFB) reproducibility when changing PE layouts, for many compsets. + However, setting this flag does NOT guarantee that all calculations throughout all + model configurations are BFB independent of processor count. + + + @@ -2515,21 +2530,6 @@ add aoflux calculation to runseq - - - - - - - logical - TRUE,FALSE - FALSE - run_flags - env_run.xml - turns on coupler bit-for-bit reproducibility with varying pe counts - - - ========================================= Notes: diff --git a/cime_config/namelist_definition_drv.xml b/cime_config/namelist_definition_drv.xml index 7eae8def6..b5ba5f863 100644 --- a/cime_config/namelist_definition_drv.xml +++ b/cime_config/namelist_definition_drv.xml @@ -2676,6 +2676,26 @@ + + + logical + reprosum + ALLCOMP_attributes + + If true, perform some calculations (in CMEPS and elsewhere) such as global sums + using a less efficient algorithm that is independent of processor count. This allows + for Bit-For-Bit (BFB) reproducibility when changing PE layouts, for many compsets. + However, setting this flag does NOT guarantee that all calculations throughout all + model configurations are BFB independent of processor count. + + + $BFBFLAG + + + diff --git a/mediator/CMakeLists.txt b/mediator/CMakeLists.txt index 4b35db94b..0a65b69e9 100644 --- a/mediator/CMakeLists.txt +++ b/mediator/CMakeLists.txt @@ -9,6 +9,7 @@ set(SRCFILES esmFldsExchange_cesm_mod.F90 med_fraction_mod.F90 esmFldsExchange_ufs_mod.F90 med_io_mod.F90 med_phases_history_mod.F90 med_phases_prep_ocn_mod.F90 med_utils_mod.F90 esmFlds.F90 med_kind_mod.F90 + med_global_sums_mod.F90 med_phases_prep_rof_mod.F90 med_constants_mod.F90 med_map_mod.F90 med_phases_prep_atm_mod.F90 med_phases_prep_wav_mod.F90 diff --git a/mediator/med_global_sums_mod.F90 b/mediator/med_global_sums_mod.F90 new file mode 100644 index 000000000..d611d064a --- /dev/null +++ b/mediator/med_global_sums_mod.F90 @@ -0,0 +1,157 @@ +module med_global_sums_mod + + !----------------------------------------------------------------------------- + ! Computation of global sums, optionally in a manner that is independent of processor count. + ! + ! NOTE: This is deliberately duplicated from dshr_global_sums in + ! cdeps/dshr/dshr_global_sums_mod.F90. The duplication is accepted because neither CMEPS nor + ! CDEPS is available in all of the builds that need this functionality, and putting it in the + ! CESM share code would require a separate, nearly complete reimplementation for builds that + ! do not have access to that code. If you change this module, make the corresponding change + ! in CDEPS. + !----------------------------------------------------------------------------- + + use med_kind_mod, only : r8=>SHR_KIND_R8, cs=>SHR_KIND_CS + use ESMF , only : ESMF_GridComp, ESMF_GridCompGet, ESMF_VM, ESMF_VMGet + use ESMF , only : ESMF_VMAllreduce, ESMF_REDUCE_SUM, ESMF_SUCCESS + use NUOPC , only : NUOPC_CompAttributeGet + use shr_log_mod, only : shr_log_error + use med_utils_mod, only : chkerr => med_utils_ChkErr +#ifdef CESMCOUPLED + use shr_reprosum_mod, only : shr_reprosum_calc +#endif + + implicit none + private + + public :: med_global_sums + + character(*), parameter :: u_FILE_u = & + __FILE__ + + ! med_global_sums is generic over the number of fields being summed: the 1fld version takes + ! summands for a single field and returns a scalar sum; the nflds version takes summands for + ! an arbitrary number of fields and returns an array of sums. + interface med_global_sums + module procedure med_global_sums_1fld + module procedure med_global_sums_nflds + end interface med_global_sums + +!=============================================================================== +contains +!=============================================================================== + + subroutine med_global_sums_nflds(gcomp, local_summands, global_sums, rc) + + ! Compute global sums of the given local summands. + ! + ! local_summands has dimensions [number of local grid cells, number of fields]; + ! global_sums has dimensions [number of fields]. + ! + ! If the bfbflag attribute is true, the sums are computed in a manner that is + ! independent of processor count; otherwise a cheaper, processor-count-dependent + ! algorithm is used. Note that, for the sums to truly be independent of processor + ! count, the caller needs to provide the individual summands here rather than + ! pre-summing them locally. + + ! input/output variables + type(ESMF_GridComp), intent(in) :: gcomp + real(r8), intent(in) :: local_summands(:,:) + real(r8), intent(out) :: global_sums(:) + integer , intent(out) :: rc + + ! local variables + type(ESMF_VM) :: vm + integer :: nsummands ! number of local grid cells + integer :: nflds ! number of fields + integer :: n, nf +#ifdef CESMCOUPLED + integer :: mpicom ! MPI communicator of this component +#endif + character(len=CS) :: cvalue + logical :: isPresent, isSet + logical :: bfbflag ! value of the bfbflag attribute + real(r8) :: local_sums(size(local_summands, 2)) ! local sums over grid cells, for each field + + character(len=*), parameter :: subname = '(med_global_sums_nflds)' + !--------------------------------------------------------------- + + rc = ESMF_SUCCESS + + call NUOPC_CompAttributeGet(gcomp, name='bfbflag', value=cvalue, & + isPresent=isPresent, isSet=isSet, rc=rc) + if (chkerr(rc,__LINE__,u_FILE_u)) return + bfbflag = .false. + if (isPresent .and. isSet) then + read(cvalue,*) bfbflag + end if + + nsummands = size(local_summands, 1) + nflds = size(local_summands, 2) + + if (size(global_sums) /= nflds) then + call shr_log_error(subname//' ERROR: size of global_sums must agree with the '// & + 'second dimension of local_summands', line=__LINE__, file=u_FILE_u, rc=rc) + return + end if + + call ESMF_GridCompGet(gcomp, vm=vm, rc=rc) + if (chkerr(rc,__LINE__,u_FILE_u)) return + + if (bfbflag) then +#ifdef CESMCOUPLED + call ESMF_VMGet(vm, mpiCommunicator=mpicom, rc=rc) + if (chkerr(rc,__LINE__,u_FILE_u)) return + call shr_reprosum_calc(arr=local_summands, arr_gsum=global_sums, & + nsummands=nsummands, dsummands=nsummands, nflds=nflds, commid=mpicom) +#else + call shr_log_error(subname//' ERROR: bfbflag is true, but this build does not have '// & + 'access to the shr_reprosum module needed for bfbflag to be true', & + line=__LINE__, file=u_FILE_u, rc=rc) + return +#endif + else + do nf = 1, nflds + local_sums(nf) = 0.0_r8 + do n = 1, nsummands + local_sums(nf) = local_sums(nf) + local_summands(n,nf) + end do + end do + call ESMF_VMAllreduce(vm, senddata=local_sums, recvdata=global_sums, count=nflds, & + reduceflag=ESMF_REDUCE_SUM, rc=rc) + if (chkerr(rc,__LINE__,u_FILE_u)) return + end if + + end subroutine med_global_sums_nflds + +!=============================================================================== + + subroutine med_global_sums_1fld(gcomp, local_summands, global_sum, rc) + + ! Compute the global sum of the given local summands, for a single field. + ! + ! This is a convenience wrapper around med_global_sums_nflds for the common case of a + ! single field; see that routine for details. + + ! input/output variables + type(ESMF_GridComp), intent(in) :: gcomp + real(r8) , intent(in) :: local_summands(:) + real(r8) , intent(out) :: global_sum + integer , intent(out) :: rc + + ! local variables + real(r8), allocatable :: summands_2d(:,:) + real(r8) :: global_sums(1) + !--------------------------------------------------------------- + + rc = ESMF_SUCCESS + + allocate(summands_2d(size(local_summands), 1)) + summands_2d(:,1) = local_summands(:) + call med_global_sums_nflds(gcomp, summands_2d, global_sums, rc) + if (chkerr(rc,__LINE__,u_FILE_u)) return + global_sum = global_sums(1) + + end subroutine med_global_sums_1fld + +end module med_global_sums_mod diff --git a/mediator/med_phases_post_rof_mod.F90 b/mediator/med_phases_post_rof_mod.F90 index a809bceee..7e3774db9 100644 --- a/mediator/med_phases_post_rof_mod.F90 +++ b/mediator/med_phases_post_rof_mod.F90 @@ -6,16 +6,16 @@ module med_phases_post_rof_mod use NUOPC , only : NUOPC_CompAttributeGet use ESMF , only : ESMF_Clock, ESMF_ClockIsCreated use ESMF , only : ESMF_LogWrite, ESMF_LOGMSG_INFO, ESMF_SUCCESS - use ESMF , only : ESMF_GridComp, ESMF_GridCompGet + use ESMF , only : ESMF_GridComp use ESMF , only : ESMF_Mesh, ESMF_MESHLOC_ELEMENT, ESMF_TYPEKIND_R8 use ESMF , only : ESMF_Field use ESMF , only : ESMF_FieldBundle, ESMF_FieldBundleCreate use ESMF , only : ESMF_FieldBundleGet, ESMF_FieldBundleAdd - use ESMF , only : ESMF_VM, ESMF_VMAllreduce, ESMF_REDUCE_SUM use med_kind_mod , only : CX=>SHR_KIND_CX, CS=>SHR_KIND_CS, CL=>SHR_KIND_CL, R8=>SHR_KIND_R8 use med_internalstate_mod , only : complnd, compocn, compice, comprof use med_internalstate_mod , only : InternalState, maintask, logunit use med_utils_mod , only : chkerr => med_utils_ChkErr + use med_global_sums_mod , only : med_global_sums use med_constants_mod , only : dbug_flag => med_constants_dbug_flag use med_phases_history_mod, only : med_phases_history_write_comp use med_map_mod , only : med_map_field_packed @@ -323,15 +323,18 @@ subroutine med_phases_post_rof_remove_negative_runoff(gcomp, field_name, rc) ! local variables type(InternalState) :: is_local - type(ESMF_VM) :: vm real(r8), pointer :: runoff_flux(:) ! temporary 1d pointer real(r8), pointer :: areas(:) - real(r8) :: local_positive(1), global_positive(1) - real(r8) :: local_negative(1), global_negative(1) + ! The local_* variables below hold per-grid-cell contributions to the corresponding global + ! sums. These are kept per grid cell, rather than being summed up locally here, so that the + ! global sums can be computed in a manner that is independent of processor count if bfbflag + ! is set. + real(r8), allocatable :: local_positive(:), local_negative(:) + real(r8) :: global_positive, global_negative real(r8) :: global_sum real(r8) :: multiplier - real(r8) :: local_positive_final(1), global_positive_final(1) - real(r8) :: local_negative_final(1), global_negative_final(1) + real(r8), allocatable :: local_positive_final(:), local_negative_final(:) + real(r8) :: global_positive_final, global_negative_final real(r8) :: global_sum_final integer :: n @@ -358,36 +361,36 @@ subroutine med_phases_post_rof_remove_negative_runoff(gcomp, field_name, rc) call fldbun_getdata1d(FBrof_r, trim(field_name), runoff_flux, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return - local_positive(1) = 0.0_r8 - local_negative(1) = 0.0_r8 + allocate(local_positive(size(runoff_flux))) + allocate(local_negative(size(runoff_flux))) + local_positive(:) = 0.0_r8 + local_negative(:) = 0.0_r8 do n = 1, size(runoff_flux) if (runoff_flux(n) >= 0.0_r8) then - local_positive(1) = local_positive(1) + areas(n) * runoff_flux(n) + local_positive(n) = areas(n) * runoff_flux(n) else - local_negative(1) = local_negative(1) + areas(n) * runoff_flux(n) + local_negative(n) = areas(n) * runoff_flux(n) end if end do - call ESMF_GridCompGet(gcomp, vm=vm, rc=rc) + call med_global_sums(gcomp, local_positive, global_positive, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return - call ESMF_VMAllreduce(vm, senddata=local_positive, recvdata=global_positive, count=1, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) + call med_global_sums(gcomp, local_negative, global_negative, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return - call ESMF_VMAllreduce(vm, senddata=local_negative, recvdata=global_negative, count=1, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - global_sum = global_positive(1) + global_negative(1) + deallocate(local_positive) + deallocate(local_negative) + global_sum = global_positive + global_negative if (maintask .and. dbug_flag > dbug_threshold) then write(logunit,'(a)') subname//' Before correction: '//trim(field_name) - write(logunit,'(a,e27.17)') subname//' global_positive = ', global_positive(1) - write(logunit,'(a,e27.17)') subname//' global_negative = ', global_negative(1) + write(logunit,'(a,e27.17)') subname//' global_positive = ', global_positive + write(logunit,'(a,e27.17)') subname//' global_negative = ', global_negative write(logunit,'(a,e27.17)') subname//' global_sum = ', global_sum end if if (global_sum > 0.0_r8) then ! There is enough positive runoff to absorb all of the negative runoff; so set ! negative runoff to 0 and downweight positive runoff to conserve. - multiplier = global_sum/global_positive(1) + multiplier = global_sum/global_positive do n = 1, size(runoff_flux) if (runoff_flux(n) > 0.0_r8) then runoff_flux(n) = runoff_flux(n) * multiplier @@ -399,7 +402,7 @@ subroutine med_phases_post_rof_remove_negative_runoff(gcomp, field_name, rc) ! There is more negative than positive runoff. Hopefully this happens rarely, if ! ever; so set positive runoff to 0 and downweight negative runoff to minimize ! negative runoff and conserve. - multiplier = global_sum/global_negative(1) + multiplier = global_sum/global_negative do n = 1, size(runoff_flux) if (runoff_flux(n) < 0.0_r8) then runoff_flux(n) = runoff_flux(n) * multiplier @@ -417,26 +420,28 @@ subroutine med_phases_post_rof_remove_negative_runoff(gcomp, field_name, rc) if (dbug_flag > dbug_threshold) then ! Recompute positives, negatives and total sum for output diagnostic purposes - local_positive_final(1) = 0.0_r8 - local_negative_final(1) = 0.0_r8 + allocate(local_positive_final(size(runoff_flux))) + allocate(local_negative_final(size(runoff_flux))) + local_positive_final(:) = 0.0_r8 + local_negative_final(:) = 0.0_r8 do n = 1, size(runoff_flux) if (runoff_flux(n) >= 0.0_r8) then - local_positive_final(1) = local_positive_final(1) + areas(n) * runoff_flux(n) + local_positive_final(n) = areas(n) * runoff_flux(n) else - local_negative_final(1) = local_negative_final(1) + areas(n) * runoff_flux(n) + local_negative_final(n) = areas(n) * runoff_flux(n) end if end do - call ESMF_VMAllreduce(vm, senddata=local_positive_final, recvdata=global_positive_final, count=1, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) + call med_global_sums(gcomp, local_positive_final, global_positive_final, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return - call ESMF_VMAllreduce(vm, senddata=local_negative_final, recvdata=global_negative_final, count=1, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) + call med_global_sums(gcomp, local_negative_final, global_negative_final, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return - global_sum_final = global_positive_final(1) + global_negative_final(1) + deallocate(local_positive_final) + deallocate(local_negative_final) + global_sum_final = global_positive_final + global_negative_final if (maintask) then write(logunit,'(a)') subname//' After correction: '//trim(field_name) - write(logunit,'(a,e27.17)') subname//' global_positive_final = ', global_positive_final(1) - write(logunit,'(a,e27.17)') subname//' global_negative_final = ', global_negative_final(1) + write(logunit,'(a,e27.17)') subname//' global_positive_final = ', global_positive_final + write(logunit,'(a,e27.17)') subname//' global_negative_final = ', global_negative_final write(logunit,'(a,e27.17)') subname//' global_sum_final = ', global_sum_final end if end if diff --git a/mediator/med_phases_prep_atm_mod.F90 b/mediator/med_phases_prep_atm_mod.F90 index b230ade31..b38147f8a 100644 --- a/mediator/med_phases_prep_atm_mod.F90 +++ b/mediator/med_phases_prep_atm_mod.F90 @@ -7,10 +7,11 @@ module med_phases_prep_atm_mod use med_kind_mod , only : CX=>SHR_KIND_CX, CS=>SHR_KIND_CS, CL=>SHR_KIND_CL, R8=>SHR_KIND_R8 use ESMF , only : ESMF_LogWrite, ESMF_LOGMSG_INFO, ESMF_SUCCESS use ESMF , only : ESMF_Field, ESMF_FieldGet, ESMF_FieldBundleGet - use ESMF , only : ESMF_GridComp, ESMF_GridCompGet + use ESMF , only : ESMF_GridComp use med_constants_mod , only : dbug_flag => med_constants_dbug_flag use med_utils_mod , only : memcheck => med_memcheck use med_utils_mod , only : chkerr => med_utils_ChkErr + use med_global_sums_mod , only : med_global_sums use med_methods_mod , only : FB_diagnose => med_methods_FB_diagnose use med_methods_mod , only : FB_fldchk => med_methods_FB_FldChk use med_methods_mod , only : FB_getfldptr=> med_methods_FB_GetFldPtr @@ -270,9 +271,6 @@ subroutine med_phases_prep_atm_enthalpy_correction (gcomp, hcorr, rc) ! 'Foxx_hrofl','Foxx_rofi','Foxx_hrofi','Foxx_rofl_glc', ! 'Foxx_hrofl_glc','Foxx_rofi_glc','Foxx_hrofi_glc' - use ESMF , only : ESMF_VMAllreduce, ESMF_GridCompGet, ESMF_REDUCE_SUM - use ESMF , only : ESMF_VM - ! input/output variables type(ESMF_GridComp) , intent(in) :: gcomp real(r8) , intent(in) :: hcorr(:) @@ -280,9 +278,6 @@ subroutine med_phases_prep_atm_enthalpy_correction (gcomp, hcorr, rc) ! local variables type(InternalState) :: is_local - integer :: n - real(r8) :: local_htot_corr(1) - type(ESMF_VM) :: vm !--------------------------------------- rc = ESMF_SUCCESS @@ -291,15 +286,7 @@ subroutine med_phases_prep_atm_enthalpy_correction (gcomp, hcorr, rc) call ESMF_GridCompGetInternalState(gcomp, is_local, rc) if (chkErr(rc,__LINE__,u_FILE_u)) return - ! Determine sum of enthalpy correction for each hcorr index locally - local_htot_corr(1) = 0._r8 - do n = 1,size(hcorr) - local_htot_corr(1) = local_htot_corr(1) + hcorr(n) - end do - call ESMF_GridCompGet(gcomp, vm=vm, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - call ESMF_VMAllreduce(vm, senddata=local_htot_corr, recvdata=global_htot_corr, count=1, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) + call med_global_sums(gcomp, hcorr, global_htot_corr(1), rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return end subroutine med_phases_prep_atm_enthalpy_correction diff --git a/mediator/med_phases_prep_glc_mod.F90 b/mediator/med_phases_prep_glc_mod.F90 index 0d28197be..8c6be77b6 100644 --- a/mediator/med_phases_prep_glc_mod.F90 +++ b/mediator/med_phases_prep_glc_mod.F90 @@ -8,13 +8,12 @@ module med_phases_prep_glc_mod use NUOPC , only : NUOPC_CompAttributeGet use NUOPC_Model , only : NUOPC_ModelGet use ESMF , only : ESMF_LogWrite, ESMF_LOGMSG_INFO, ESMF_SUCCESS, ESMF_FAILURE - use ESMF , only : ESMF_VM, ESMF_VMGet, ESMF_VMAllReduce, ESMF_REDUCE_SUM, ESMF_REDUCE_MAX use ESMF , only : ESMF_Clock, ESMF_ClockCreate, ESMF_ClockIsCreated use ESMF , only : ESMF_ClockGetAlarm, ESMF_ClockAdvance, ESMF_ClockGet use ESMF , only : ESMF_Time, ESMF_TimeGet use ESMF , only : ESMF_Alarm, ESMF_AlarmCreate, ESMF_AlarmSet, ESMF_AlarmGet use ESMF , only : ESMF_AlarmIsRinging, ESMF_AlarmRingerOff - use ESMF , only : ESMF_GridComp, ESMF_GridCompGet + use ESMF , only : ESMF_GridComp use ESMF , only : ESMF_FieldBundle, ESMF_FieldBundleGet, ESMF_FieldBundleAdd use ESMF , only : ESMF_FieldBundleCreate, ESMF_FieldBundleIsCreated use ESMF , only : ESMF_Field, ESMF_FieldGet, ESMF_FieldCreate @@ -44,6 +43,7 @@ module med_phases_prep_glc_mod use med_methods_mod , only : med_methods_FB_check_wtracers use med_field_info_mod , only : med_field_info_type, med_field_info_array_from_state use med_utils_mod , only : chkerr => med_utils_ChkErr + use med_global_sums_mod , only : med_global_sums use nuopc_shr_methods , only : alarmInit use glc_elevclass_mod , only : glc_get_num_elevation_classes use glc_elevclass_mod , only : glc_get_elevation_classes @@ -63,8 +63,8 @@ module med_phases_prep_glc_mod private :: med_phases_prep_glc_map_lnd2glc private :: med_phases_prep_glc_renormalize_smb - private :: renormalize_smb_accumulate_sums_l - private :: renormalize_smb_accumulate_sums_g + private :: renormalize_smb_compute_summands_l + private :: renormalize_smb_compute_summands_g private :: renormalize_smb_do_renormalization ! ----------------- @@ -1114,7 +1114,6 @@ subroutine med_phases_prep_glc_renormalize_smb(gcomp, ns, rc) ! local variables type(InternalState) :: is_local - type(ESMF_VM) :: vm real(r8) , pointer :: qice_g(:) ! SMB (Flgl_qice) on glc grid without elev classes real(r8) , pointer :: qice_l_ec(:,:) ! SMB (Flgl_qice_elev) on land grid with elev classes real(r8) , pointer :: qice_g_wtracers(:,:) ! SMB water tracers (Flgl_qice_wtracers) on glc grid without elev classes @@ -1129,12 +1128,19 @@ subroutine med_phases_prep_glc_renormalize_smb(gcomp, ns, rc) real(r8) , pointer :: dataptr1d(:) ! temporary 1d pointer integer :: n, t - ! local and global sums of accumulation and ablation; used to compute renormalization factors - ! the first element of each of these is for bulk water; the remaining elements are for water tracers - real(r8) :: local_accum_lnd(1+num_wtracers), global_accum_lnd(1+num_wtracers) - real(r8) :: local_accum_glc(1+num_wtracers), global_accum_glc(1+num_wtracers) - real(r8) :: local_ablat_lnd(1+num_wtracers), global_ablat_lnd(1+num_wtracers) - real(r8) :: local_ablat_glc(1+num_wtracers), global_ablat_glc(1+num_wtracers) + ! local summands and global sums of accumulation and ablation; used to compute + ! renormalization factors. + ! + ! For the local variables, the first dimension is over local grid cells (on the land or glc + ! grid, as appropriate) and the second dimension is over bulk water plus water tracers; the + ! global variables have just the latter dimension. For the latter dimension, the first + ! element is for bulk water; the remaining elements are for water tracers. + real(r8), allocatable :: local_accum_lnd(:,:) + real(r8), allocatable :: local_ablat_lnd(:,:) + real(r8), allocatable :: local_accum_glc(:,:) + real(r8), allocatable :: local_ablat_glc(:,:) + real(r8) :: global_accum_lnd(1+num_wtracers), global_ablat_lnd(1+num_wtracers) + real(r8) :: global_accum_glc(1+num_wtracers), global_ablat_glc(1+num_wtracers) ! areas real(r8), allocatable :: effective_area_l(:) ! effective areas on the land grid: grid cell area multiplied by min(lndfrac,icemask_l). @@ -1226,40 +1232,40 @@ subroutine med_phases_prep_glc_renormalize_smb(gcomp, ns, rc) effective_area_l(n) = min(lndfrac(n), icemask_l(n)) * is_local%wrap%mesh_info(complnd)%areas(n) end do - ! determine accumulation and ablation sums for qice on the land grid + ! determine accumulation and ablation summands for qice on the land grid + allocate(local_accum_lnd(size(lndfrac), 1+num_wtracers)) + allocate(local_ablat_lnd(size(lndfrac), 1+num_wtracers)) call fldbun_getdata2d(FBlndAccum2glc_l, trim(qice_elev_fieldname), qice_l_ec, rc) if (chkErr(rc,__LINE__,u_FILE_u)) return - call renormalize_smb_accumulate_sums_l(qice_l_ec, effective_area_l, frac_l_ec, & - local_accum_lnd(1), local_ablat_lnd(1)) + call renormalize_smb_compute_summands_l(qice_l_ec, effective_area_l, frac_l_ec, & + local_accum_lnd(:,1), local_ablat_lnd(:,1)) - ! and, similarly, determine accumulation and ablation sums for qice water tracers on the land grid + ! and, similarly, determine accumulation and ablation summands for qice water tracers on the land grid if (has_wtracers) then call fldbun_getdata3d(FBlndAccum2glc_l, trim(qice_elev_wtracers_fieldname), qice_l_ec_wtracers, rc) if (chkErr(rc,__LINE__,u_FILE_u)) return do t = 1, num_wtracers - call renormalize_smb_accumulate_sums_l(qice_l_ec_wtracers(:,t,:), effective_area_l, frac_l_ec, & - local_accum_lnd(1+t), local_ablat_lnd(1+t)) + call renormalize_smb_compute_summands_l(qice_l_ec_wtracers(:,t,:), effective_area_l, frac_l_ec, & + local_accum_lnd(:,1+t), local_ablat_lnd(:,1+t)) end do end if deallocate(effective_area_l) ! determine global accum/ablat on the land grid - call ESMF_GridCompGet(gcomp, vm=vm, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - call ESMF_VMAllreduce(vm, senddata=local_accum_lnd, recvdata=global_accum_lnd, count=1+num_wtracers, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) + call med_global_sums(gcomp, local_accum_lnd, global_accum_lnd, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return - call ESMF_VMAllreduce(vm, senddata=local_ablat_lnd, recvdata=global_ablat_lnd, count=1+num_wtracers, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) + call med_global_sums(gcomp, local_ablat_lnd, global_ablat_lnd, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return + deallocate(local_accum_lnd) + deallocate(local_ablat_lnd) if (maintask) then - write(logunit,'(a,d21.10)') trim(subname)//'global_accum_lnd = ', global_accum_lnd(1) - write(logunit,'(a,d21.10)') trim(subname)//'global_ablat_lnd = ', global_ablat_lnd(1) + write(logunit,'(a,d28.17)') trim(subname)//'global_accum_lnd = ', global_accum_lnd(1) + write(logunit,'(a,d28.17)') trim(subname)//'global_ablat_lnd = ', global_ablat_lnd(1) do t = 1, num_wtracers - write(logunit,'(a,i0,a,d21.10)') trim(subname)//'tracer #', t, & + write(logunit,'(a,i0,a,d28.17)') trim(subname)//'tracer #', t, & ': global_accum_lnd = ', global_accum_lnd(1+t) - write(logunit,'(a,i0,a,d21.10)') trim(subname)//'tracer #', t, & + write(logunit,'(a,i0,a,d28.17)') trim(subname)//'tracer #', t, & ': global_ablat_lnd = ', global_ablat_lnd(1+t) end do endif @@ -1272,34 +1278,38 @@ subroutine med_phases_prep_glc_renormalize_smb(gcomp, ns, rc) call fldbun_getdata1d(is_local%wrap%FBImp(compglc(ns),compglc(ns)), 'Sg_area', area_g, rc) if (chkerr(rc,__LINE__,u_FILE_u)) return - ! determine accumulation and ablation sums for qice on the glc grid + ! determine accumulation and ablation summands for qice on the glc grid call fldbun_getdata1d(is_local%wrap%FBExp(compglc(ns)), qice_fieldname, qice_g, rc) if (chkerr(rc,__LINE__,u_FILE_u)) return - call renormalize_smb_accumulate_sums_g(qice_g, icemask_g, area_g, & - local_accum_glc(1), local_ablat_glc(1)) + allocate(local_accum_glc(size(qice_g), 1+num_wtracers)) + allocate(local_ablat_glc(size(qice_g), 1+num_wtracers)) + call renormalize_smb_compute_summands_g(qice_g, icemask_g, area_g, & + local_accum_glc(:,1), local_ablat_glc(:,1)) - ! and, similarly, determine accumulation and ablation sums for qice water tracers on the glc grid + ! and, similarly, determine accumulation and ablation summands for qice water tracers on the glc grid if (has_wtracers) then call fldbun_getdata2d(is_local%wrap%FBExp(compglc(ns)), qice_wtracers_fieldname, qice_g_wtracers, rc) if (chkerr(rc,__LINE__,u_FILE_u)) return do t = 1, num_wtracers - call renormalize_smb_accumulate_sums_g(qice_g_wtracers(t,:), icemask_g, area_g, & - local_accum_glc(1+t), local_ablat_glc(1+t)) + call renormalize_smb_compute_summands_g(qice_g_wtracers(t,:), icemask_g, area_g, & + local_accum_glc(:,1+t), local_ablat_glc(:,1+t)) end do end if ! determine global accum/ablat on the glc grid - call ESMF_VMAllreduce(vm, senddata=local_accum_glc, recvdata=global_accum_glc, count=1+num_wtracers, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) - call ESMF_VMAllreduce(vm, senddata=local_ablat_glc, recvdata=global_ablat_glc, count=1+num_wtracers, & - reduceflag=ESMF_REDUCE_SUM, rc=rc) + call med_global_sums(gcomp, local_accum_glc, global_accum_glc, rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + call med_global_sums(gcomp, local_ablat_glc, global_ablat_glc, rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + deallocate(local_accum_glc) + deallocate(local_ablat_glc) if (maintask) then - write(logunit,'(a,d21.10)') trim(subname)//'global_accum_glc = ', global_accum_glc(1) - write(logunit,'(a,d21.10)') trim(subname)//'global_ablat_glc = ', global_ablat_glc(1) + write(logunit,'(a,d28.17)') trim(subname)//'global_accum_glc = ', global_accum_glc(1) + write(logunit,'(a,d28.17)') trim(subname)//'global_ablat_glc = ', global_ablat_glc(1) do t = 1, num_wtracers - write(logunit,'(a,i0,a,d21.10)') trim(subname)//'tracer #', t, & + write(logunit,'(a,i0,a,d28.17)') trim(subname)//'tracer #', t, & ': global_accum_glc = ', global_accum_glc(1+t) - write(logunit,'(a,i0,a,d21.10)') trim(subname)//'tracer #', t, & + write(logunit,'(a,i0,a,d28.17)') trim(subname)//'tracer #', t, & ': global_ablat_glc = ', global_ablat_glc(1+t) end do endif @@ -1318,66 +1328,73 @@ subroutine med_phases_prep_glc_renormalize_smb(gcomp, ns, rc) end subroutine med_phases_prep_glc_renormalize_smb !================================================================================================ - subroutine renormalize_smb_accumulate_sums_l(qice_l_ec, effective_area_l, frac_l_ec, & + subroutine renormalize_smb_compute_summands_l(qice_l_ec, effective_area_l, frac_l_ec, & accum_l, ablat_l) - ! Compute local accumulation and ablation sums on land grid + ! Compute the per-grid-cell accumulation and ablation summands on the land grid. These are + ! the summands of the global accumulation and ablation sums; they are returned per grid cell + ! (rather than being summed here) so that the caller can compute the global sums in a + ! processor-count-independent manner if desired. Note that the sum over elevation classes + ! for a given grid cell is done here, since that dimension is not distributed. ! input/output variables real(r8), intent(in) :: qice_l_ec(:,:) real(r8), intent(in) :: effective_area_l(:) real(r8), intent(in) :: frac_l_ec(:,:) - real(r8), intent(out) :: accum_l - real(r8), intent(out) :: ablat_l + real(r8), intent(out) :: accum_l(:) + real(r8), intent(out) :: ablat_l(:) ! local variables integer :: n, ec !--------------------------------------------------------------- - accum_l = 0.0_r8 - ablat_l = 0.0_r8 + accum_l(:) = 0.0_r8 + ablat_l(:) = 0.0_r8 do n = 1, size(effective_area_l) if (effective_area_l(n) > 0.0_r8) then do ec = 1, ungriddedCount if (qice_l_ec(ec,n) >= 0.0_r8) then - accum_l = accum_l + effective_area_l(n) * frac_l_ec(ec,n) * qice_l_ec(ec,n) + accum_l(n) = accum_l(n) + effective_area_l(n) * frac_l_ec(ec,n) * qice_l_ec(ec,n) else - ablat_l = ablat_l + effective_area_l(n) * frac_l_ec(ec,n) * qice_l_ec(ec,n) + ablat_l(n) = ablat_l(n) + effective_area_l(n) * frac_l_ec(ec,n) * qice_l_ec(ec,n) end if end do end if end do - end subroutine renormalize_smb_accumulate_sums_l + end subroutine renormalize_smb_compute_summands_l !================================================================================================ - subroutine renormalize_smb_accumulate_sums_g(qice_g, icemask_g, area_g, & + subroutine renormalize_smb_compute_summands_g(qice_g, icemask_g, area_g, & accum_g, ablat_g) - ! Compute local accumulation and ablation sums on glc grid + ! Compute the per-grid-cell accumulation and ablation summands on the glc grid. These are + ! the summands of the global accumulation and ablation sums; they are returned per grid cell + ! (rather than being summed here) so that the caller can compute the global sums in a + ! processor-count-independent manner if desired. ! input/output variables real(r8), intent(in) :: qice_g(:) real(r8), intent(in) :: icemask_g(:) real(r8), intent(in) :: area_g(:) - real(r8), intent(out) :: accum_g - real(r8), intent(out) :: ablat_g + real(r8), intent(out) :: accum_g(:) + real(r8), intent(out) :: ablat_g(:) ! local variables integer :: n !--------------------------------------------------------------- - accum_g = 0.0_r8 - ablat_g = 0.0_r8 + accum_g(:) = 0.0_r8 + ablat_g(:) = 0.0_r8 do n = 1, size(qice_g) if (qice_g(n) >= 0.0_r8) then - accum_g = accum_g + icemask_g(n) * area_g(n) * qice_g(n) + accum_g(n) = icemask_g(n) * area_g(n) * qice_g(n) else - ablat_g = ablat_g + icemask_g(n) * area_g(n) * qice_g(n) + ablat_g(n) = icemask_g(n) * area_g(n) * qice_g(n) end if end do - end subroutine renormalize_smb_accumulate_sums_g + end subroutine renormalize_smb_compute_summands_g !================================================================================================ subroutine renormalize_smb_do_renormalization(global_accum_lnd, global_ablat_lnd, & @@ -1415,13 +1432,13 @@ subroutine renormalize_smb_do_renormalization(global_accum_lnd, global_ablat_lnd endif if (maintask) then if (present(tracer_num)) then - write(logunit,'(a,i0,a,d21.10)') trim(caller_subname)//'tracer #', tracer_num, & + write(logunit,'(a,i0,a,d28.17)') trim(caller_subname)//'tracer #', tracer_num, & ': accum_renorm_factor = ', accum_renorm_factor - write(logunit,'(a,i0,a,d21.10)') trim(caller_subname)//'tracer #', tracer_num, & + write(logunit,'(a,i0,a,d28.17)') trim(caller_subname)//'tracer #', tracer_num, & ': ablat_renorm_factor = ', ablat_renorm_factor else - write(logunit,'(a,d21.10)') trim(caller_subname)//'accum_renorm_factor = ', accum_renorm_factor - write(logunit,'(a,d21.10)') trim(caller_subname)//'ablat_renorm_factor = ', ablat_renorm_factor + write(logunit,'(a,d28.17)') trim(caller_subname)//'accum_renorm_factor = ', accum_renorm_factor + write(logunit,'(a,d28.17)') trim(caller_subname)//'ablat_renorm_factor = ', ablat_renorm_factor endif endif