Changeset 3896 for branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/DYN

Timestamp:

2013-05-02T14:40:43+02:00 (11 years ago)

Author:

acc

Message:

Branch 2013/dev_r3858_NOC_ZTC, #863. Fix to dynvor.F90 to remove instability when using the EEN vorticity scheme with z-tilde. Also cosmetic changes and better initial settings for zstar/ztilde when using full ice embedding

Location:

branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/DYN

Files:

• dynvor.F90 (modified) (3 diffs)
• sshwzv.F90 (modified) (1 diff)

branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90

@@ -572 +572 @@
       INTEGER  ::   ierr                                          ! local integer
       REAL(wp) ::   zfac12, zua, zva                              ! local scalars
+      REAL(wp) ::   zmsk, ze3                                     ! local scalars
       !                                                           !  3D workspace 
       REAL(wp), POINTER    , DIMENSION(:,:  )         :: zwx, zwy, zwz
@@ -577 +578 @@
 #if defined key_vvl
       REAL(wp), POINTER    , DIMENSION(:,:,:)         :: ze3f     !  3D workspace (lk_vvl=T)
-#endif
-#if ! defined key_vvl
+#else
       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), SAVE   :: ze3f     ! lk_vvl=F, ze3f=1/e3f saved one for all
 #endif
@@ -604 +604 @@
       ENDIF
 
-      IF( kt == nit000 .OR. lk_vvl ) THEN      ! reciprocal of e3 at F-point (masked averaging of e3t)
+      IF( kt == nit000 .OR. lk_vvl ) THEN      ! reciprocal of e3 at F-point (masked averaging of e3t over ocean points)
          DO jk = 1, jpk
             DO jj = 1, jpjm1
                DO ji = 1, jpim1
-                  ze3f(ji,jj,jk) = ( fse3t(ji,jj+1,jk)*tmask(ji,jj+1,jk) + fse3t(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk)   &
-                     &             + fse3t(ji,jj  ,jk)*tmask(ji,jj  ,jk) + fse3t(ji+1,jj  ,jk)*tmask(ji+1,jj  ,jk) ) * 0.25
-                  IF( ze3f(ji,jj,jk) /= 0._wp )   ze3f(ji,jj,jk) = 1._wp / ze3f(ji,jj,jk)
+                  ze3  = ( fse3t(ji,jj+1,jk)*tmask(ji,jj+1,jk) + fse3t(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk)   &
+                     &   + fse3t(ji,jj  ,jk)*tmask(ji,jj  ,jk) + fse3t(ji+1,jj  ,jk)*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( ze3 /= 0._wp )   ze3f(ji,jj,jk) = zmsk / ze3
                END DO
             END DO

branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90

@@ -110 +110 @@
       ! In forward Euler time stepping case, the same formulation as in the leap-frog case can be used
       ! because emp_b field is initialized with the vlaues of emp field. Hence, 0.5 * ( emp + emp_b ) = emp
-      z1_rau0 = 0.5 / rau0
+      z1_rau0 = 0.5_wp * r1_rau0
       ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * tmask(:,:,1)