Changeset 2005 for branches/DEV_r1837_MLF/NEMO/OPA_SRC/DYN/dynnxt.F90

Timestamp:

2010-07-09T15:07:02+02:00 (14 years ago)

Author:

mlelod

Message:

ticket: #663 MLF: second part (local compatibility essentially)

File:

• branches/DEV_r1837_MLF/NEMO/OPA_SRC/DYN/dynnxt.F90 (modified) (5 diffs)

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

@@ -87 +87 @@
       !!               un,vn   now horizontal velocity of next time-step
       !!----------------------------------------------------------------------
+#if defined key_vvl
+      USE oce, ONLY :   ze3u_f => ta   ! use ta as 3D workspace
+      USE oce, ONLY :   ze3v_f => sa   ! use sa as 3D workspace
+#endif
       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
       !!
@@ -95 +99 @@
       REAL(wp) ::   zue3a , zue3n , zue3b    ! temporary scalar
       REAL(wp) ::   zve3a , zve3n , zve3b    !    -         -
-      REAL(wp) ::   ze3u_b, ze3u_n, ze3u_a   !    -         -
-      REAL(wp) ::   ze3v_b, ze3v_n, ze3v_a   !    -         - 
       REAL(wp) ::   zuf   , zvf              !    -         - 
+      REAL(wp) ::   zec                      !    -         - 
+      REAL(wp) ::   zv_t_ij  , zv_t_ip1j     !     -        -
+      REAL(wp) ::   zv_t_ijp1                !     -        -
+      REAL(wp), DIMENSION(jpi,jpj) ::  zs_t, zs_u_1, zs_v_1      ! temporary 2D workspace
       !!----------------------------------------------------------------------
 
@@ -216 +222 @@
                DO jj = 1, jpj
                   DO ji = 1, jpi    
-                     zuf = atfp * ( ub(ji,jj,jk) + ua(ji,jj,jk) ) + atfp1 * un(ji,jj,jk)
-                     zvf = atfp * ( vb(ji,jj,jk) + va(ji,jj,jk) ) + atfp1 * vn(ji,jj,jk)
+                     zuf = un(ji,jj,jk) + atfp * ( ub(ji,jj,jk) - 2.e0 * un(ji,jj,jk) + ua(ji,jj,jk) )
+                     zvf = vn(ji,jj,jk) + atfp * ( vb(ji,jj,jk) - 2.e0 * vn(ji,jj,jk) + va(ji,jj,jk) )
                      !
                      ub(ji,jj,jk) = zuf                      ! ub <-- filtered velocity
@@ -227 +233 @@
             END DO
          ELSE                                                ! applied on thickness weighted velocity
+            ! Before scale factors at (t/u/v)-points (actually "now filtered" and futur "before")
+            ! ======================================
+            ! Scale factor at t-points
+            ! ------------------------
+            fse3t_b(:,:,:) = fse3t_n(:,:,:) + atfp * fse3t_m(:,:,:)
+            ! Add volume filter correction only at the first level of t-point scale factors
+            zec = atfp * rdt / rau0
+            fse3t_b(:,:,1) = fse3t_b(:,:,1) - zec * ( emp_b(:,:) - emp(:,:) ) * tmask(:,:,1)
+            ! Scale factor at (u/v)-points
+            ! ------------------------
+            ! surface at t-points and inverse surface at (u/v)-points used in surface averaging computations
+            zs_t  (:,:) =       e1t(:,:) * e2t(:,:)
+            zs_u_1(:,:) = 0.5 / e1u(:,:) * e2u(:,:)
+            zs_v_1(:,:) = 0.5 / e1v(:,:) * e2v(:,:)
+            ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points
+            DO jk = 1, jpkm1
+               DO jj = 1, jpjm1
+                  DO ji = 1, jpim1
+                     zv_t_ij          = zs_t(ji  ,jj  ) * fse3t_b(ji  ,jj  ,jk)
+                     zv_t_ip1j        = zs_t(ji+1,jj  ) * fse3t_b(ji+1,jj  ,jk)
+                     zv_t_ijp1        = zs_t(ji  ,jj+1) * fse3t_b(ji  ,jj+1,jk)
+                     ze3u_f(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) )
+                     ze3v_f(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) )
+                  END DO
+               END DO
+            END DO
+            CALL lbc_lnk( ze3u_f, 'U', 1. )               ! lateral boundary conditions
+            CALL lbc_lnk( ze3v_f, 'U', 1. )
+            ! Add initial scale factor to scale factor anomaly
+            ze3u_f(:,:,:) = ze3u_f(:,:,:) + fse3u_0(:,:,:)
+            ze3v_f(:,:,:) = ze3v_f(:,:,:) + fse3v_0(:,:,:)
+            
             DO jk = 1, jpkm1
                DO jj = 1, jpj
-                  DO ji = 1, jpi
-                     ze3u_a = fse3u_a(ji,jj,jk)
-                     ze3v_a = fse3v_a(ji,jj,jk)
-                     ze3u_n = fse3u_n(ji,jj,jk)
-                     ze3v_n = fse3v_n(ji,jj,jk)
-                     ze3u_b = fse3u_b(ji,jj,jk)
-                     ze3v_b = fse3v_b(ji,jj,jk)
+                  DO ji = 1, jpim
+                     zue3a = ua(ji,jj,jk) * fse3u_a(ji,jj,jk)
+                     zve3a = va(ji,jj,jk) * fse3v_a(ji,jj,jk)
+                     zue3n = un(ji,jj,jk) * fse3u_n(ji,jj,jk)
+                     zve3n = vn(ji,jj,jk) * fse3v_n(ji,jj,jk)
+                     zue3b = ub(ji,jj,jk) * fse3u_b(ji,jj,jk)
+                     zve3b = vb(ji,jj,jk) * fse3v_b(ji,jj,jk)
                      !
-                     zue3a = ua(ji,jj,jk) * ze3u_a
-                     zve3a = va(ji,jj,jk) * ze3v_a
-                     zue3n = un(ji,jj,jk) * ze3u_n
-                     zve3n = vn(ji,jj,jk) * ze3v_n
-                     zue3b = ub(ji,jj,jk) * ze3u_b
-                     zve3b = vb(ji,jj,jk) * ze3v_b
-                     !
-                     zuf  = ( atfp * ( zue3b  + zue3a  ) + atfp1 * zue3n  )   &
-                        & / ( atfp * ( ze3u_b + ze3u_a ) + atfp1 * ze3u_n ) * umask(ji,jj,jk)
-                     zvf  = ( atfp * ( zve3b  + zve3a  ) + atfp1 * zve3n  )   &
-                        & / ( atfp * ( ze3v_b + ze3v_a ) + atfp1 * ze3v_n ) * vmask(ji,jj,jk)
+                     zuf  = ( zue3n + atfp * ( zue3b  - 2.e0 * zue3n  + zue3a  ) ) / ze3u_f(ji,jj,jk)
+                     zvf  = ( zve3n + atfp * ( zve3b  - 2.e0 * zve3n  + zve3a  ) ) / ze3v_f(ji,jj,jk)
                      !
                      ub(ji,jj,jk) = zuf                      ! ub <-- filtered velocity
                      vb(ji,jj,jk) = zvf
+                     fse3u_b(ji,jj,jk) = ze3u_f(ji,jj,jk)    ! e3u_b <-- filtered scale factor
+                     fse3v_b(ji,jj,jk) = ze3v_f(ji,jj,jk)
                      un(ji,jj,jk) = ua(ji,jj,jk)             ! un <-- ua
                      vn(ji,jj,jk) = va(ji,jj,jk)
@@ -256 +287 @@
                END DO
             END DO
+            CALL lbc_lnk( ub, 'U', -1. )         ! local domain boundaries
+            CALL lbc_lnk( vb, 'V', -1. ) 
          ENDIF
       ENDIF