Changeset 2068 for branches/DEV_r1837_MLF/NEMO/OPA_SRC/DYN/sshwzv.F90

Timestamp:

2010-09-06T17:56:51+02:00 (14 years ago)

Author:

mlelod

Message:

ticket: #663 ensuring restartability and conservation

File:

• branches/DEV_r1837_MLF/NEMO/OPA_SRC/DYN/sshwzv.F90 (modified) (7 diffs)

branches/DEV_r1837_MLF/NEMO/OPA_SRC/DYN/sshwzv.F90

@@ -104 +104 @@
             DO jj = 1, jpjm1
                DO ji = 1, jpim1      ! NO Vector Opt.
-                  sshf_b(ji,jj) = 0.5  * umask(ji,jj,1) * umask(ji,jj+1,1)                   &
-                     &                 / ( e1f(ji,jj  ) * e2f(ji,jj  ) )                     &
-                     &                 * ( e1u(ji,jj  ) * e2u(ji,jj  ) * sshu_b(ji,jj  )     &
-                     &                   + e1u(ji,jj+1) * e2u(ji,jj+1) * sshu_b(ji,jj+1) )
-               END DO
-            END DO
-            DO jj = 1, jpjm1
-               DO ji = 1, jpim1      ! NO Vector Opt.
                   sshf_n(ji,jj) = 0.5  * umask(ji,jj,1) * umask(ji,jj+1,1)                   &
                        &               / ( e1f(ji,jj  ) * e2f(ji,jj  ) )                     &
@@ -118 +110 @@
                END DO
             END DO
-            CALL lbc_lnk( sshf_b, 'F', 1. )   ;   CALL lbc_lnk( sshf_n, 'F', 1. )
+            CALL lbc_lnk( sshf_n, 'F', 1. )
          ENDIF
          !
@@ -157 +149 @@
       !                                           !   After Sea Surface Height   !
       !                                           !------------------------------!
-      z1_rau0 = 0.5 / rau0
-      !
       zhdiv(:,:) = 0.e0
       DO jk = 1, jpkm1                                 ! Horizontal divergence of barotropic transports
@@ -166 +156 @@
       ! 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
-      ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * tmask(:,:,1)
+      IF( neuler == 0 .AND. kt == nit000 ) THEN
+         z1_rau0 = 1.e0 / rau0
+         ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 *   emp(:,:)                + zhdiv(:,:) )  ) * tmask(:,:,1)
+      ELSE
+         z1_rau0 = 0.5 / rau0
+         ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * tmask(:,:,1)
+      ENDIF
 
 #if defined key_obc
@@ -194 +190 @@
       !                                           !----------------------------------------!
       ! Needed for Robert-Asselin time filter and for Brown & Campana semi implicit hydrostatic presure gradient
-      fse3t_m(:,:,:) =          fse3t_b(:,:,:)   &
-         &             - 2.e0 * fse3t_n(:,:,:)   &
-         &             +        fse3t_a(:,:,:)
+      DO jk = 1, jpk
+         fse3t_d(:,:,jk) =          fse3t_b(:,:,jk)   &
+            &              - 2.e0 * fse3t_n(:,:,jk)   &
+            &              +        fse3t_a(:,:,jk)
+      ENDDO
       !                                           !------------------------------!
       !                                           !     Now Vertical Velocity    !
@@ -251 +249 @@
       !!
       INTEGER  ::   ji, jj   ! dummy loop indices
-      REAL(wp) ::   zec      ! temporary scalare
+      REAL(wp) ::   zec      ! temporary scalar
       !!----------------------------------------------------------------------
 
@@ -281 +279 @@
             CALL lbc_lnk( sshf_n, 'F', 1. )
          ELSE                                           ! Leap-Frog time-stepping: Asselin filter + swap
+            zec = atfp * rdt / rau0
             DO jj = 1, jpj
                DO ji = 1, jpi                                ! before <-- now filtered