From df74f181a53b966a1047db51f1e6fa4f0ed1f0bb Mon Sep 17 00:00:00 2001
From: sebastien masson <smasson@amphitrite.locean-ipsl.upmc.fr>
Date: Mon, 26 Sep 2022 10:02:58 +0200
Subject: [PATCH] minor fix for NP folding reproducibility, #68
---
arch/arch-osx_gfortran.fcm | 4 +-
src/OCE/DYN/dynspg_ts.F90 | 36 +++---
src/OCE/DYN/dynvor.F90 | 244 +++++++++++++++++--------------------
src/OCE/LDF/ldfslp.F90 | 16 +--
4 files changed, 138 insertions(+), 162 deletions(-)
diff --git a/arch/arch-osx_gfortran.fcm b/arch/arch-osx_gfortran.fcm
index 3ce12293a..ece845b94 100644
--- a/arch/arch-osx_gfortran.fcm
+++ b/arch/arch-osx_gfortran.fcm
@@ -35,8 +35,8 @@
%CPP cpp -Dkey_nosignedzero
%FC mpif90
-%PROD_FCFLAGS -fdefault-real-8 -O3 -funroll-all-loops -fcray-pointer -ffree-line-length-none -fallow-argument-mismatch
-%DEBUG_FCFLAGS -fdefault-real-8 -O0 -g -fbacktrace -funroll-all-loops -fcray-pointer -ffree-line-length-none -fcheck=all -finit-real=nan
+%PROD_FCFLAGS -fdefault-real-8 -O3 -funroll-all-loops -fcray-pointer -ffree-line-length-none -fallow-argument-mismatch
+%DEBUG_FCFLAGS -fdefault-real-8 -O0 -g -fbacktrace -funroll-all-loops -fcray-pointer -ffree-line-length-none -fcheck=all -finit-real=nan -fallow-argument-mismatch
%FFLAGS %FCFLAGS
%LD %FC
%LDFLAGS
diff --git a/src/OCE/DYN/dynspg_ts.F90 b/src/OCE/DYN/dynspg_ts.F90
index a1e78f536..bd40bf540 100644
--- a/src/OCE/DYN/dynspg_ts.F90
+++ b/src/OCE/DYN/dynspg_ts.F90
@@ -1247,16 +1247,16 @@ CONTAINS
SELECT CASE( nn_e3f_typ ) !* ff_f/e3 at F-point
CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
DO_2D( 0, 0, 0, 0 )
- zwz(ji,jj) = ( ht(ji,jj+1) + ht(ji+1,jj+1) &
- & + ht(ji,jj ) + ht(ji+1,jj ) ) * 0.25_wp
+ zwz(ji,jj) = ( ( ht(ji,jj+1) + ht(ji+1,jj+1) ) & ! need additional () for reproducibility around NP
+ & + ( ht(ji,jj ) + ht(ji+1,jj ) ) ) * 0.25_wp
IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj)
END_2D
CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
DO_2D( 0, 0, 0, 0 )
- zwz(ji,jj) = ( ht(ji,jj+1) + ht(ji+1,jj+1) &
- & + ht(ji,jj ) + ht(ji+1,jj ) ) &
- & / ( MAX(ssmask(ji,jj+1) + ssmask(ji+1,jj+1) &
- & + ssmask(ji,jj ) + ssmask(ji+1,jj ) , 1._wp ) )
+ zwz(ji,jj) = ( ( ht(ji,jj+1) + ht(ji+1,jj+1) ) & ! need additional () for
+ & + ( ht(ji,jj ) + ht(ji+1,jj ) ) ) & ! reproducibility around NP
+ & / ( MAX( ssmask(ji,jj+1) + ssmask(ji+1,jj+1) &
+ & + ssmask(ji,jj ) + ssmask(ji+1,jj ), 1._wp ) )
IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj)
END_2D
END SELECT
@@ -1327,24 +1327,24 @@ CONTAINS
!
CASE( np_ENS ) ! enstrophy conserving scheme (f-point)
DO_2D( 0, 0, 0, 0 )
- zy1 = r1_8 * ( zhV(ji ,jj-1) + zhV(ji+1,jj-1) &
- & + zhV(ji ,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj)
- zx1 = - r1_8 * ( zhU(ji-1,jj ) + zhU(ji-1,jj+1) &
- & + zhU(ji ,jj ) + zhU(ji ,jj+1) ) * r1_e2v(ji,jj)
+ zy1 = r1_8 * ( ( zhV(ji ,jj-1) + zhV(ji+1,jj-1) ) & ! need additional () for
+ & + ( zhV(ji ,jj ) + zhV(ji+1,jj ) ) ) * r1_e1u(ji,jj) ! reproducibility around NP
+ zx1 = - r1_8 * ( ( zhU(ji-1,jj ) + zhU(ji-1,jj+1) ) & ! need additional () for
+ & + ( zhU(ji ,jj ) + zhU(ji ,jj+1) ) ) * r1_e2v(ji,jj) ! reproducibility around NP
zu_trd(ji,jj) = zy1 * ( zwz(ji ,jj-1) + zwz(ji,jj) )
zv_trd(ji,jj) = zx1 * ( zwz(ji-1,jj ) + zwz(ji,jj) )
END_2D
!
CASE( np_EET , np_EEN ) ! energy & enstrophy scheme (using e3t or e3f)
DO_2D( 0, 0, 0, 0 )
- zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zhV(ji ,jj ) &
- & + ftnw(ji+1,jj) * zhV(ji+1,jj ) &
- & + ftse(ji,jj ) * zhV(ji ,jj-1) &
- & + ftsw(ji+1,jj) * zhV(ji+1,jj-1) )
- zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * ( ftsw(ji,jj+1) * zhU(ji-1,jj+1) &
- & + ftse(ji,jj+1) * zhU(ji ,jj+1) &
- & + ftnw(ji,jj ) * zhU(ji-1,jj ) &
- & + ftne(ji,jj ) * zhU(ji ,jj ) )
+ zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ( ftne(ji,jj ) * zhV(ji ,jj ) & ! need additional () for
+ & + ftnw(ji+1,jj) * zhV(ji+1,jj ) ) & ! reproducibility around NP
+ & + ( ftse(ji,jj ) * zhV(ji ,jj-1) &
+ & + ftsw(ji+1,jj) * zhV(ji+1,jj-1) ) )
+ zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * ( ( ftsw(ji,jj+1) * zhU(ji-1,jj+1) &
+ & + ftse(ji,jj+1) * zhU(ji ,jj+1) ) &
+ & + ( ftnw(ji,jj ) * zhU(ji-1,jj ) &
+ & + ftne(ji,jj ) * zhU(ji ,jj ) ) )
END_2D
!
END SELECT
diff --git a/src/OCE/DYN/dynvor.F90 b/src/OCE/DYN/dynvor.F90
index b7a84124a..4089c0b9f 100644
--- a/src/OCE/DYN/dynvor.F90
+++ b/src/OCE/DYN/dynvor.F90
@@ -361,8 +361,8 @@ CONTAINS
CASE ( np_RVO , np_CRV ) !* relative vorticity at f-point is used
ALLOCATE( zwz(A2D(1)) )
DO_2D( 1, 1, 1, 1 )
- zwz(ji,jj) = ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
+ zwz(ji,jj) = ( ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for
+ & - ( e1u(ji,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj) ! NP repro
END_2D
IF( ln_dynvor_msk ) THEN ! mask relative vorticity
DO_2D( 1, 1, 1, 1 )
@@ -380,8 +380,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( 0, 1, 0, 1 )
- zwt(ji,jj) = r1_4 * ( zwz(ji-1,jj ) + zwz(ji,jj ) &
- & + zwz(ji-1,jj-1) + zwz(ji,jj-1) ) &
+ zwt(ji,jj) = r1_4 * ( ( zwz(ji-1,jj ) + zwz(ji,jj ) ) & ! need additional () for
+ & + ( zwz(ji-1,jj-1) + zwz(ji,jj-1) ) ) & ! reproducibility around NP
& * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm)
END_2D
CASE ( np_MET ) !* metric term
@@ -392,8 +392,8 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( 0, 1, 0, 1 )
- zwt(ji,jj) = ( ff_t(ji,jj) + r1_4 * ( zwz(ji-1,jj ) + zwz(ji,jj ) &
- & + zwz(ji-1,jj-1) + zwz(ji,jj-1) ) ) &
+ zwt(ji,jj) = ( ff_t(ji,jj) + r1_4 * ( ( zwz(ji-1,jj ) + zwz(ji,jj ) ) & ! need additional () for
+ & + ( zwz(ji-1,jj-1) + zwz(ji,jj-1) ) ) ) & ! reproducibility around NP
& * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm)
END_2D
CASE ( np_CME ) !* Coriolis + metric
@@ -473,8 +473,8 @@ CONTAINS
ALLOCATE( zwz(A2D(1),jpk) )
DO jk = 1, jpkm1 ! Horizontal slab
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
- zwz(ji,jj,jk) = ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
+ zwz(ji,jj,jk) = ( ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! need additional () for NP repro
+ & - ( e1u(ji,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask relative vorticity
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
@@ -498,8 +498,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( 0, 1, 0, 1 )
- zwt(ji,jj) = r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) &
- & + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) &
+ zwt(ji,jj) = r1_4 * ( ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) ) & ! need additional () for
+ & + ( zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) ) & ! reproducibility around NP
& * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm)
END_2D
CASE ( np_MET ) !* metric term
@@ -510,8 +510,8 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( 0, 1, 0, 1 )
- zwt(ji,jj) = ( ff_t(ji,jj) + r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) &
- & + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) ) &
+ zwt(ji,jj) = ( ff_t(ji,jj) + r1_4 * ( ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) ) & ! need additional () for
+ & + ( zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) ) ) & ! reproducibility around NP
& * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm)
END_2D
CASE ( np_CME ) !* Coriolis + metric
@@ -596,8 +596,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( 1, 0, 1, 0 )
- zwz(ji,jj) = ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
+ zwz(ji,jj) = ( ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for
+ & - ( e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj) ! NP repro
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
DO_2D( 1, 0, 1, 0 )
@@ -611,8 +611,8 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( 1, 0, 1, 0 )
- zwz(ji,jj) = ff_f(ji,jj) + ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
+ zwz(ji,jj) = ff_f(ji,jj) + ( ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for NP repro
+ & - ( e1u(ji,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity (NOT the Coriolis term)
DO_2D( 1, 0, 1, 0 )
@@ -636,22 +636,22 @@ CONTAINS
SELECT CASE( nn_e3f_typ ) !== potential vorticity ==!
CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
DO_2D( 1, 0, 1, 0 )
- ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) &
- & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * 4._wp / ze3f
ELSE ; zwz(ji,jj) = 0._wp
ENDIF
END_2D
CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
DO_2D( 1, 0, 1, 0 )
- ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) &
- & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) )
- zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
- & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
+ zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
+ & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * zmsk / ze3f
ELSE ; zwz(ji,jj) = 0._wp
ENDIF
@@ -728,8 +728,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( 1, 0, 1, 0 )
- zwz(ji,jj) = ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
+ zwz(ji,jj) = ( ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for NP repro
+ & - ( e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
DO_2D( 1, 0, 1, 0 )
@@ -743,8 +743,8 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( 1, 0, 1, 0 )
- zwz(ji,jj) = ff_f(ji,jj) + ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
+ zwz(ji,jj) = ff_f(ji,jj) + ( ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) &! add () for NP repro
+ & -( e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity (NOT the Coriolis term)
DO_2D( 1, 0, 1, 0 )
@@ -769,22 +769,22 @@ CONTAINS
SELECT CASE( nn_e3f_typ ) !== potential vorticity ==!
CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
DO_2D( 1, 0, 1, 0 )
- ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) &
- & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * 4._wp / ze3f
ELSE ; zwz(ji,jj) = 0._wp
ENDIF
END_2D
CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
DO_2D( 1, 0, 1, 0 )
- ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) &
- & + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) )
- zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
- & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
+ zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
+ & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
IF( ze3f /= 0._wp ) THEN ; zwz(ji,jj) = zwz(ji,jj) * zmsk / ze3f
ELSE ; zwz(ji,jj) = 0._wp
ENDIF
@@ -799,10 +799,10 @@ CONTAINS
!
! !== compute and add the vorticity term trend =!
DO_2D( 0, 0, 0, 0 )
- zuav = r1_8 * r1_e1u(ji,jj) * ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) &
- & + zwy(ji ,jj ) + zwy(ji+1,jj ) )
- zvau =-r1_8 * r1_e2v(ji,jj) * ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) &
- & + zwx(ji ,jj ) + zwx(ji ,jj+1) )
+ zuav = r1_8 * r1_e1u(ji,jj) * ( ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) ) & ! need additional () for
+ & + ( zwy(ji ,jj ) + zwy(ji+1,jj ) ) ) ! reproducibility around NP
+ zvau =-r1_8 * r1_e2v(ji,jj) * ( ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) ) &
+ & + ( zwx(ji ,jj ) + zwx(ji ,jj+1) ) )
pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zuav * ( zwz(ji ,jj-1) + zwz(ji,jj) )
pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zvau * ( zwz(ji-1,jj ) + zwz(ji,jj) )
END_2D
@@ -871,26 +871,22 @@ CONTAINS
SELECT CASE( nn_e3f_typ ) ! == reciprocal of e3 at F-point
CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
DO_2D( 1, 1, 1, 1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- ze3f = ( (e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)) &
- & + (e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk)) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
IF( ze3f /= 0._wp ) THEN ; z1_e3f(ji,jj) = 4._wp / ze3f
ELSE ; z1_e3f(ji,jj) = 0._wp
ENDIF
END_2D
CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
DO_2D( 1, 1, 1, 1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- ze3f = ( (e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)) &
- & + (e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk)) )
- zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
- & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
+ zmsk = (tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
+ & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
IF( ze3f /= 0._wp ) THEN ; z1_e3f(ji,jj) = zmsk / ze3f
ELSE ; z1_e3f(ji,jj) = 0._wp
ENDIF
@@ -906,8 +902,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( 1, 1, 1, 1 )
- zwz(ji,jj) = ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji ,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)*z1_e3f(ji,jj)
+ zwz(ji,jj) = ( ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for NP repro
+ & - ( e1u(ji ,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj)*z1_e3f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
DO_2D( 1, 1, 1, 1 )
@@ -921,12 +917,8 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( 1, 1, 1, 1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- zwz(ji,jj) = ( ff_f(ji,jj) + ( ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & ) & ! bracket for halo 1 - halo 2 compatibility
- & - ( e1u(ji ,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) &
- & ) & ! bracket for halo 1 - halo 2 compatibility
+ zwz(ji,jj) = ( ff_f(ji,jj) + ( ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for
+ & - ( e1u(ji ,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) & ! NP repro
& ) * r1_e1e2f(ji,jj) ) * z1_e3f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
@@ -957,11 +949,11 @@ CONTAINS
ztsw(ji,jj) = zwz(ji ,jj-1) + zwz(ji-1,jj-1) + zwz(ji-1,jj )
END_2D
!
- DO_2D( 0, 0, 0, 0 )
- zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) &
- & + ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) )
- zva = - r1_12 * r1_e2v(ji,jj) * ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) &
- & + ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) )
+ DO_2D( 0, 0, 0, 0 ) ! add () for NP repro
+ zua = + r1_12 * r1_e1u(ji,jj) * ( ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) ) & ! add () for
+ & + ( ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) ) ) ! NP repro
+ zva = - r1_12 * r1_e2v(ji,jj) * ( ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) ) &
+ & + ( ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) ) )
pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zua
pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
END_2D
@@ -1031,26 +1023,22 @@ CONTAINS
SELECT CASE( nn_e3f_typ ) ! == reciprocal of e3 at F-point
CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- ze3f = ( (e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)) &
- & + (e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk)) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
IF( ze3f /= 0._wp ) THEN ; z1_e3f(ji,jj) = 4._wp / ze3f
ELSE ; z1_e3f(ji,jj) = 0._wp
ENDIF
END_2D
CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- ze3f = ( (e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
- & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)) &
- & + (e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
- & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk)) )
- zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
- & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
+ ze3f = ( ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
+ & + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) ) )
+ zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
+ & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
IF( ze3f /= 0._wp ) THEN ; z1_e3f(ji,jj) = zmsk / ze3f
ELSE ; z1_e3f(ji,jj) = 0._wp
ENDIF
@@ -1066,8 +1054,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
- zwz(ji,jj,jk) = ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & - e1u(ji ,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)*z1_e3f(ji,jj)
+ zwz(ji,jj,jk) = ( ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for NP repro
+ & - ( e1u(ji ,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) * r1_e1e2f(ji,jj)*z1_e3f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
@@ -1081,12 +1069,8 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- zwz(ji,jj,jk) = ( ff_f(ji,jj) + ( ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
- & ) & ! bracket for halo 1 - halo 2 compatibility
- & - ( e1u(ji ,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) &
- & ) & ! bracket for halo 1 - halo 2 compatibility
+ zwz(ji,jj,jk) = ( ff_f(ji,jj) + ( ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for
+ & - ( e1u(ji ,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) & ! NP repro
& ) * r1_e1e2f(ji,jj) ) * z1_e3f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
@@ -1131,10 +1115,10 @@ CONTAINS
ztsw_jp1 = zwz(ji ,jj ,jk) + zwz(ji-1,jj ,jk) + zwz(ji-1,jj+1,jk)
ztsw_ip1 = zwz(ji+1,jj-1,jk) + zwz(ji ,jj-1,jk) + zwz(ji ,jj ,jk)
!
- zua = + r1_12 * r1_e1u(ji,jj) * ( ztne * zwy + ztnw_ip1 * zwy_ip1 &
- & + ztse * zwy_jm1 + ztsw_ip1 * zwy_ip1_jm1 )
- zva = - r1_12 * r1_e2v(ji,jj) * ( ztsw_jp1 * zwx_im1_jp1 + ztse_jp1 * zwx_jp1 &
- & + ztnw * zwx_im1 + ztne * zwx )
+ zua = + r1_12 * r1_e1u(ji,jj) * ( ( ztne * zwy + ztnw_ip1 * zwy_ip1 ) & ! need additional () for
+ & + ( ztse * zwy_jm1 + ztsw_ip1 * zwy_ip1_jm1 ) ) ! reproducibility around NP
+ zva = - r1_12 * r1_e2v(ji,jj) * ( ( ztsw_jp1 * zwx_im1_jp1 + ztse_jp1 * zwx_jp1 ) &
+ & + ( ztnw * zwx_im1 + ztne * zwx ) )
pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zua
pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
END_3D
@@ -1156,10 +1140,10 @@ CONTAINS
END_2D
!
DO_2D( 0, 0, 0, 0 )
- zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) &
- & + ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) )
- zva = - r1_12 * r1_e2v(ji,jj) * ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) &
- & + ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) )
+ zua = + r1_12 * r1_e1u(ji,jj) * ( ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) ) & ! add () for
+ & + ( ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) ) ) ! NP repro
+ zva = - r1_12 * r1_e2v(ji,jj) * ( ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) ) &
+ & + ( ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) ) )
pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zua
pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
END_2D
@@ -1223,10 +1207,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( 1, 1, 1, 1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- zwz(ji,jj) = ( (e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk)) &
- & - (e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk)) ) &
+ zwz(ji,jj) = ( ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for
+ & - ( e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) & ! NP reproducibility
& * r1_e1e2f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
@@ -1241,10 +1223,8 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( 1, 1, 1, 1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- zwz(ji,jj) = ( ff_f(ji,jj) + ( (e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk)) &
- & - (e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk)) ) &
+ zwz(ji,jj) = ( ff_f(ji,jj) + ( ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! add () for
+ & - ( e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) & ! NP repro
& * r1_e1e2f(ji,jj) )
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
@@ -1278,10 +1258,10 @@ CONTAINS
END_2D
!
DO_2D( 0, 0, 0, 0 )
- zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) &
- & + ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) )
- zva = - r1_12 * r1_e2v(ji,jj) * ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) &
- & + ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) )
+ zua = + r1_12 * r1_e1u(ji,jj) * ( ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) ) & ! add () for
+ & + ( ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) ) ) ! NP repro
+ zva = - r1_12 * r1_e2v(ji,jj) * ( ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) ) &
+ & + ( ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) ) )
pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zua
pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
END_2D
@@ -1343,10 +1323,8 @@ CONTAINS
END_2D
CASE ( np_RVO ) !* relative vorticity
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- zwz(ji,jj,jk) = ( (e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk)) &
- & - (e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk)) ) &
+ zwz(ji,jj,jk) = ( ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) & ! need additional ()
+ & - ( e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) & ! for NP reproducibility
& * r1_e1e2f(ji,jj)
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
@@ -1361,11 +1339,9 @@ CONTAINS
END_2D
CASE ( np_CRV ) !* Coriolis + relative vorticity
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
- ! round brackets added to fix the order of floating point operations
- ! needed to ensure halo 1 - halo 2 compatibility
- zwz(ji,jj,jk) = ( ff_f(ji,jj) + ( (e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk)) &
- & - (e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk)) ) &
- & * r1_e1e2f(ji,jj) )
+ zwz(ji,jj,jk) = ( ff_f(ji,jj) + ( ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) ) &! add () for
+ & - ( e1u(ji ,jj+1) * pu(ji ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk) ) ) &! NP repro
+ & * r1_e1e2f(ji,jj) )
END_2D
IF( ln_dynvor_msk ) THEN ! mask the relative vorticity
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
@@ -1407,10 +1383,10 @@ CONTAINS
END_2D
!
DO_2D( 0, 0, 0, 0 )
- zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) &
- & + ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) )
- zva = - r1_12 * r1_e2v(ji,jj) * ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) &
- & + ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) )
+ zua = + r1_12 * r1_e1u(ji,jj) * ( ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) ) & ! add () for
+ & + ( ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) ) ) ! NP repro
+ zva = - r1_12 * r1_e2v(ji,jj) * ( ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) ) &
+ & + ( ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) ) )
pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zua
pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
END_2D
@@ -1529,10 +1505,10 @@ CONTAINS
SELECT CASE( nn_e3f_typ )
CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
DO_3D( 0, 0, 0, 0, 1, jpk )
- e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) &
- & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) &
- & + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) &
- & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) * 0.25_wp
+ e3f_0vor(ji,jj,jk) = ( ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) &
+ & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) ) * 0.25_wp
END_3D
CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
DO_3D( 0, 0, 0, 0, 1, jpk )
@@ -1540,10 +1516,10 @@ CONTAINS
& + tmask(ji,jj ,jk) +tmask(ji+1,jj ,jk) )
!
IF( zmsk /= 0._wp ) THEN
- e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) &
- & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) &
- & + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) &
- & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) / zmsk
+ e3f_0vor(ji,jj,jk) = ( ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) & ! need additional () for
+ & + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) ) & ! reproducibility around NP
+ & + ( e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) &
+ & + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) ) / zmsk
ELSE ; e3f_0vor(ji,jj,jk) = 0._wp
ENDIF
END_3D
diff --git a/src/OCE/LDF/ldfslp.F90 b/src/OCE/LDF/ldfslp.F90
index 6231b844d..4dd61e07e 100644
--- a/src/OCE/LDF/ldfslp.F90
+++ b/src/OCE/LDF/ldfslp.F90
@@ -266,10 +266,10 @@ CONTAINS
& + umask(ji-1,jj,jk-1) + umask(ji,jj,jk-1) , zeps ) * e1t(ji,jj)
zcj = MAX( vmask(ji,jj-1,jk ) + vmask(ji,jj,jk-1) &
& + vmask(ji,jj-1,jk-1) + vmask(ji,jj,jk ) , zeps ) * e2t(ji,jj)
- zai = ( zgru (ji-1,jj,jk ) + zgru (ji,jj,jk-1) &
- & + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk ) ) / zci * wmask (ji,jj,jk)
- zaj = ( zgrv (ji,jj-1,jk ) + zgrv (ji,jj,jk-1) &
- & + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk ) ) / zcj * wmask (ji,jj,jk)
+ zai = ( ( zgru (ji-1,jj,jk ) + zgru (ji,jj,jk-1) ) & ! need additional () for reproducibility around NP
+ & + ( zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk ) ) ) / zci * wmask (ji,jj,jk)
+ zaj = ( ( zgrv (ji,jj-1,jk ) + zgrv (ji,jj,jk-1) ) & ! need additional () for reproducibility around NP
+ & + ( zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk ) ) ) / zcj * wmask (ji,jj,jk)
! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0.
! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt)
zbi = MIN( zbw ,- 100._wp* ABS( zai ) , -7.e+3_wp/e3w(ji,jj,jk,Kmm)* ABS( zai ) )
@@ -645,10 +645,10 @@ CONTAINS
& + umask(ji-1,jj,ikm1) + umask(ji,jj,ikm1) , zeps ) * e1t(ji,jj)
zcj = MAX( vmask(ji,jj-1,ik ) + vmask(ji,jj,ik ) &
& + vmask(ji,jj-1,ikm1) + vmask(ji,jj,ikm1) , zeps ) * e2t(ji,jj)
- zai = ( p_gru(ji-1,jj,ik ) + p_gru(ji,jj,ik) &
- & + p_gru(ji-1,jj,ikm1) + p_gru(ji,jj,ikm1 ) ) / zci * tmask(ji,jj,ik)
- zaj = ( p_grv(ji,jj-1,ik ) + p_grv(ji,jj,ik ) &
- & + p_grv(ji,jj-1,ikm1) + p_grv(ji,jj,ikm1) ) / zcj * tmask(ji,jj,ik)
+ zai = ( ( p_gru(ji-1,jj,ik ) + p_gru(ji,jj,ik ) ) & ! need additional () for reproducibility around NP
+ & + ( p_gru(ji-1,jj,ikm1) + p_gru(ji,jj,ikm1) ) ) / zci * tmask(ji,jj,ik)
+ zaj = ( ( p_grv(ji,jj-1,ik ) + p_grv(ji,jj,ik ) ) & ! need additional () for reproducibility around NP
+ & + ( p_grv(ji,jj-1,ikm1) + p_grv(ji,jj,ikm1) ) ) / zcj * tmask(ji,jj,ik)
! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0.
! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt)
zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/e3w(ji,jj,ik,Kmm)* ABS( zai ) )
--
GitLab