From 26dd003121531e7621e89e0eca4abd75b619fa58 Mon Sep 17 00:00:00 2001 From: Philippe Blain Date: Tue, 12 Jul 2022 11:47:59 -0400 Subject: [PATCH] ice_dyn_vp: use 'global_sum_prod' for bit-for-bit log reproducibility In the previous commits we ensured bit-for-bit reproducibility of the outputs when using the VP solver. Some global norms computed during the nonlinear iteration still use the same non-reproducible pattern of summing over blocks locally before performing the reduction. However, these norms are used only to monitor the convergence in the log file, as well as to exit the iteration when the required convergence level is reached ('nlres_norm'). Only 'nlres_norm' could (in theory) influence the output, but it is unlikely that a difference due to floating point errors would influence the 'if (nlres_norm < tol_nl)' condition used to exist the nonlinear iteration. Change these remaining cases to also use 'global_sum_prod' to compute the global norm, leading to bit-for-bit log reproducibility. --- cicecore/cicedynB/dynamics/ice_dyn_vp.F90 | 23 ++++++----------------- 1 file changed, 6 insertions(+), 17 deletions(-) diff --git a/cicecore/cicedynB/dynamics/ice_dyn_vp.F90 b/cicecore/cicedynB/dynamics/ice_dyn_vp.F90 index 791d8c1f2..982fa253f 100644 --- a/cicecore/cicedynB/dynamics/ice_dyn_vp.F90 +++ b/cicecore/cicedynB/dynamics/ice_dyn_vp.F90 @@ -901,7 +901,8 @@ subroutine anderson_solver (icellt , icellu , & L2norm (iblk)) enddo !$OMP END PARALLEL DO - nlres_norm = sqrt(global_sum(sum(L2norm), distrb_info)) + nlres_norm = sqrt(global_sum_prod(Fx(:,:,:), Fx(:,:,:), distrb_info, field_loc_NEcorner) + & + global_sum_prod(Fy(:,:,:), Fy(:,:,:), distrb_info, field_loc_NEcorner)) if (my_task == master_task .and. monitor_nonlin) then write(nu_diag, '(a,i4,a,d26.16)') "monitor_nonlin: iter_nonlin= ", it_nl, & " nonlin_res_L2norm= ", nlres_norm @@ -991,16 +992,8 @@ subroutine anderson_solver (icellt , icellu , & indxui (:,:), indxuj(:,:) , & res (:), & fpresx (:,:,:), fpresy (:,:,:)) - !$OMP PARALLEL DO PRIVATE(iblk) - do iblk = 1, nblocks - call calc_L2norm_squared (nx_block , ny_block , & - icellu (iblk), & - indxui (:,iblk), indxuj (:,iblk), & - fpresx(:,:,iblk), fpresy(:,:,iblk), & - L2norm (iblk)) - enddo - !$OMP END PARALLEL DO - fpres_norm = sqrt(global_sum(sum(L2norm), distrb_info)) + fpres_norm = sqrt(global_sum_prod(fpresx(:,:,:), fpresx(:,:,:), distrb_info, field_loc_NEcorner) + & + global_sum_prod(fpresy(:,:,:), fpresy(:,:,:), distrb_info, field_loc_NEcorner)) #endif if (my_task == master_task .and. monitor_nonlin) then write(nu_diag, '(a,i4,a,d26.16)') "monitor_nonlin: iter_nonlin= ", it_nl, & @@ -1128,14 +1121,10 @@ subroutine anderson_solver (icellt , icellu , & do iblk = 1, nblocks fpresx(:,:,iblk) = uvel(:,:,iblk) - uprev_k(:,:,iblk) fpresy(:,:,iblk) = vvel(:,:,iblk) - vprev_k(:,:,iblk) - call calc_L2norm_squared (nx_block , ny_block , & - icellu (iblk), & - indxui (:,iblk), indxuj (:,iblk), & - fpresx(:,:,iblk), fpresy(:,:,iblk), & - L2norm (iblk)) enddo !$OMP END PARALLEL DO - prog_norm = sqrt(global_sum(sum(L2norm), distrb_info)) + prog_norm = sqrt(global_sum_prod(fpresx(:,:,:), fpresx(:,:,:), distrb_info, field_loc_NEcorner) + & + global_sum_prod(fpresy(:,:,:), fpresy(:,:,:), distrb_info, field_loc_NEcorner)) if (my_task == master_task .and. monitor_nonlin) then write(nu_diag, '(a,i4,a,d26.16)') "monitor_nonlin: iter_nonlin= ", it_nl, & " progress_res_L2norm= ", prog_norm