Changeset 1683

Timestamp:

2009-10-28T11:24:41+01:00 (14 years ago)

Author:

rblod

Message:

Correct time stepping for volume variable and filtered free surface, see #568

File:

• trunk/NEMO/OPA_SRC/DYN/dynspg_flt.F90 (modified) (3 diffs)

trunk/NEMO/OPA_SRC/DYN/dynspg_flt.F90

@@ -29 +29 @@
    USE phycst          ! physical constants
    USE domvvl          ! variable volume
+   USE dynadv          ! advection 
    USE solmat          ! matrix construction for elliptic solvers
    USE solver          ! solver initialization
@@ -136 +137 @@
       zraur  = 1. / rauw
 
-      !! Explicit physics with thickness weighted updates
-      IF( lk_vvl ) THEN          ! variable volume 
-
-         ! Evaluate the masked next velocity (effect of the additional force not included)
-         ! -------------------   (thickness weighted velocity, surface pressure gradient already included in dyn_hpg)
-         DO jk = 1, jpkm1
-            DO jj = 2, jpjm1
-               DO ji = fs_2, fs_jpim1   ! vector opt.
-                  ua(ji,jj,jk) = (        ub(ji,jj,jk) * fse3u_b(ji,jj,jk)      &
-                     &           + z2dt * ua(ji,jj,jk) * fse3u_n(ji,jj,jk)  )   &
-                     &         / fse3u_a(ji,jj,jk) * umask(ji,jj,jk)
-                  va(ji,jj,jk) = (        vb(ji,jj,jk) * fse3v_b(ji,jj,jk)      &
-                     &           + z2dt * va(ji,jj,jk) * fse3v_n(ji,jj,jk)  )   &
-                     &         / fse3v_a(ji,jj,jk) * vmask(ji,jj,jk)
+      ! Evaluate the masked next velocity (effect of the additional force not included)
+      ! ---------------------------------  
+      IF( lk_vvl ) THEN          ! variable volume  (surface pressure gradient already included in dyn_hpg)
+         !
+         IF( ln_dynadv_vec ) THEN      ! vector form : applied on velocity
+            DO jk = 1, jpkm1
+               DO jj = 2, jpjm1
+                  DO ji = fs_2, fs_jpim1   ! vector opt.
+                     ua(ji,jj,jk) = (  ub(ji,jj,jk) + z2dt * ua(ji,jj,jk)  ) * umask(ji,jj,jk)
+                     va(ji,jj,jk) = (  vb(ji,jj,jk) + z2dt * va(ji,jj,jk)  ) * vmask(ji,jj,jk)
+                  END DO
                END DO
             END DO
-         END DO
-
-      ELSE                       ! fixed volume 
-
-         ! Surface pressure gradient (now)
-         DO jj = 2, jpjm1
+            !
+         ELSE                          ! flux form : applied on thickness weighted velocity
+            DO jk = 1, jpkm1
+               DO jj = 2, jpjm1
+                  DO ji = fs_2, fs_jpim1   ! vector opt.
+                     ua(ji,jj,jk) = (        ub(ji,jj,jk) * fse3u_b(ji,jj,jk)      &
+                        &           + z2dt * ua(ji,jj,jk) * fse3u_n(ji,jj,jk)  )   &
+                        &         / fse3u_a(ji,jj,jk) * umask(ji,jj,jk)
+                     va(ji,jj,jk) = (        vb(ji,jj,jk) * fse3v_b(ji,jj,jk)      &
+                        &           + z2dt * va(ji,jj,jk) * fse3v_n(ji,jj,jk)  )   &
+                        &         / fse3v_a(ji,jj,jk) * vmask(ji,jj,jk)
+                 END DO
+               END DO
+            END DO
+            !
+         ENDIF
+         !
+      ELSE                       ! fixed volume  (add the surface pressure gradient + unweighted time stepping)
+         !
+         DO jj = 2, jpjm1              ! Surface pressure gradient (now)
             DO ji = fs_2, fs_jpim1   ! vector opt.
                spgu(ji,jj) = - grav * ( sshn(ji+1,jj) - sshn(ji,jj) ) / e1u(ji,jj)
@@ -163 +175 @@
             END DO 
          END DO 
-         !
-         ! add the surface pressure trend to the general trend and
-         ! evaluate the masked next velocity (effect of the additional force not included)
-         DO jk = 1, jpkm1
+         DO jk = 1, jpkm1              ! unweighted time stepping 
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.