Changeset 10336

Timestamp:

2018-11-19T15:34:28+01:00 (6 years ago)

Author:

smasson

Message:

dev_r10164_HPC09_ESIWACE_PREP_MERGE: action 7: zdfiwm.F90 and zdftke.F90 optimization from ATOS, see #2133

Location:

NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ZDF

Files:

• zdfiwm.F90 (modified) (4 diffs)
• zdftke.F90 (modified) (2 diffs)

NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ZDF/zdfiwm.F90

@@ -122 +122 @@
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
-      REAL(wp) ::   zztmp        ! scalar workspace
+      REAL(wp) ::   zztmp, ztmp1, ztmp2        ! scalar workspace
       REAL(wp), DIMENSION(jpi,jpj)     ::   zfact       ! Used for vertical structure
       REAL(wp), DIMENSION(jpi,jpj)     ::   zhdep       ! Ocean depth
@@ -157 +157 @@
 !!gm gde3w ==>>>  check for ssh taken into account.... seem OK gde3w_n=gdept_n - sshn
       DO jk = 2, jpkm1              ! complete with the level-dependent part
-         zemx_iwm(:,:,jk) = zfact(:,:) * (  EXP( ( gde3w_n(:,:,jk  ) - zhdep(:,:) ) / hcri_iwm(:,:) )                      &
-            &                             - EXP( ( gde3w_n(:,:,jk-1) - zhdep(:,:) ) / hcri_iwm(:,:) )  ) * wmask(:,:,jk)   &
-            &                          / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) )
-
+         DO jj = 1, jpj             
+            DO ji = 1, jpi
+               IF ( zfact(ji,jj) == 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization
+                  zemx_iwm(ji,jj,jk) = 0._wp
+               ELSE
+                  zemx_iwm(ji,jj,jk) = zfact(ji,jj) * (  EXP( ( gde3w_n(ji,jj,jk  ) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) )     &
+                       &                               - EXP( ( gde3w_n(ji,jj,jk-1) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) ) )   &
+                       &                            / ( gde3w_n(ji,jj,jk) - gde3w_n(ji,jj,jk-1) )
+               ENDIF
+            END DO
+         END DO
 !!gm delta(gde3w_n) = e3t_n  !!  Please verify the grid-point position w versus t-point
 !!gm it seems to me that only 1/hcri_iwm  is used ==>  compute it one for all
@@ -234 +241 @@
       zwkb(:,:,1) = zhdep(:,:) * wmask(:,:,1)
       !
-      zweight(:,:,:) = 0._wp
-      DO jk = 2, jpkm1
-         zweight(:,:,jk) = MAX( 0._wp, rn2(:,:,jk) ) * hbot_iwm(:,:) * wmask(:,:,jk)                    &
-            &   * (  EXP( -zwkb(:,:,jk) / hbot_iwm(:,:) ) - EXP( -zwkb(:,:,jk-1) / hbot_iwm(:,:) )  )
+      DO jk = 2, jpkm1
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF ( rn2(ji,jj,jk) <= 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization
+                  zweight(ji,jj,jk) = 0._wp
+               ELSE
+                  zweight(ji,jj,jk) = rn2(ji,jj,jk) * hbot_iwm(ji,jj)    &
+                     &   * (  EXP( -zwkb(ji,jj,jk) / hbot_iwm(ji,jj) ) - EXP( -zwkb(ji,jj,jk-1) / hbot_iwm(ji,jj) )  )
+               ENDIF
+            END DO
+         END DO
       END DO
       !
@@ -322 +336 @@
       !      
       IF( ln_tsdiff ) THEN          !* Option for differential mixing of salinity and temperature
+         ztmp1 = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10( 1.e-20_wp ) - 0.60_wp ) )
          DO jk = 2, jpkm1              ! Calculate S/T diffusivity ratio as a function of Reb
             DO jj = 1, jpj
                DO ji = 1, jpi
-                  zav_ratio(ji,jj,jk) = ( 0.505_wp + 0.495_wp *                                                                  &
-                      &   TANH(    0.92_wp * (   LOG10(  MAX( 1.e-20_wp, zReb(ji,jj,jk) * 5._wp * r1_6 )  ) - 0.60_wp   )    )   &
-                      &                 ) * wmask(ji,jj,jk)
+                  ztmp2 = zReb(ji,jj,jk) * 5._wp * r1_6
+                  IF ( ztmp2 > 1.e-20_wp .AND. wmask(ji,jj,jk) == 1._wp ) THEN
+                     zav_ratio(ji,jj,jk) = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10(ztmp2) - 0.60_wp ) )
+                  ELSE
+                     zav_ratio(ji,jj,jk) = ztmp1 * wmask(ji,jj,jk)
+                  ENDIF
                END DO
             END DO

NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ZDF/zdftke.F90

@@ -202 +202 @@
       REAL(wp) ::   zzd_up, zzd_lw             !   -         -
       INTEGER , DIMENSION(jpi,jpj)     ::   imlc
-      REAL(wp), DIMENSION(jpi,jpj)     ::   zhlc
+      REAL(wp), DIMENSION(jpi,jpj)     ::   zhlc, zfr_i
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zpelc, zdiag, zd_up, zd_lw
       !!--------------------------------------------------------------------
@@ -290 +290 @@
          END DO
          zcof = 0.016 / SQRT( zrhoa * zcdrag )
+         DO jj = 2, jpjm1
+            DO ji = fs_2, fs_jpim1   ! vector opt.
+               zus  = zcof * SQRT( taum(ji,jj) )           ! Stokes drift
+               zfr_i(ji,jj) = ( 1._wp - 4._wp * fr_i(ji,jj) ) * zus * zus * zus * tmask(ji,jj,1) ! zus > 0. ok
+               IF (zfr_i(ji,jj) < 0. ) zfr_i(ji,jj) = 0.
+            END DO
+         END DO         
          DO jk = 2, jpkm1         !* TKE Langmuir circulation source term added to en
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
-                  zus  = zcof * SQRT( taum(ji,jj) )           ! Stokes drift
-                  !                                           ! vertical velocity due to LC
-                  zind = 0.5 - SIGN( 0.5, pdepw(ji,jj,jk) - zhlc(ji,jj) )
-                  zwlc = zind * rn_lc * zus * SIN( rpi * pdepw(ji,jj,jk) / zhlc(ji,jj) )
-                  !                                           ! TKE Langmuir circulation source term
-                  en(ji,jj,jk) = en(ji,jj,jk) + rdt * MAX(0.,1._wp - 4.*fr_i(ji,jj) ) * ( zwlc * zwlc * zwlc )   &
-                     &                              / zhlc(ji,jj) * wmask(ji,jj,jk) * tmask(ji,jj,1)
+                  IF ( zfr_i(ji,jj) /= 0. ) THEN               
+                     ! vertical velocity due to LC   
+                     IF ( pdepw(ji,jj,jk) - zhlc(ji,jj) < 0 .AND. wmask(ji,jj,jk) /= 0. ) THEN
+                        !                                           ! vertical velocity due to LC
+                        zwlc = rn_lc * SIN( rpi * pdepw(ji,jj,jk) / zhlc(ji,jj) )   ! warning: optimization: zus^3 is in zfr_i
+                        !                                           ! TKE Langmuir circulation source term
+                        en(ji,jj,jk) = en(ji,jj,jk) + rdt * zfr_i(ji,jj) * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj)
+                     ENDIF
+                  ENDIF
                END DO
             END DO