Changeset 4311

Timestamp:

2013-11-22T11:41:01+01:00 (10 years ago)

Author:

cbricaud

Message:

time split corrections ctd

Location:

branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC

Files:

• DYN/dynnxt.F90 (modified) (6 diffs)
• step.F90 (modified) (2 diffs)

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

@@ -106 +106 @@
       REAL(wp) ::   zue3a, zue3n, zue3b, zuf, zec   ! local scalars
       REAL(wp) ::   zve3a, zve3n, zve3b, zvf        !   -      -
-      REAL(wp), POINTER, DIMENSION(:,:)   ::  zua, zva
+      REAL(wp), POINTER, DIMENSION(:,:)   ::  zua, zva, zhura, zhvra
       REAL(wp), POINTER, DIMENSION(:,:,:) ::  ze3u_f, ze3v_f 
       !!----------------------------------------------------------------------
@@ -113 +113 @@
       !
       CALL wrk_alloc( jpi,jpj,jpk, ze3u_f, ze3v_f )
-      IF ( lk_dynspg_ts ) CALL wrk_alloc( jpi,jpj, zua, zva )
+      IF ( lk_dynspg_ts ) CALL wrk_alloc( jpi,jpj, zua, zva, zhura, zhvra )
       !
       IF( kt == nit000 ) THEN
@@ -162 +162 @@
       zva(:,:) = 0._wp
       IF (lk_vvl) THEN
+         zhura(:,:) = 0._wp
+         zhvra(:,:) = 0._wp
          DO jk = 1, jpkm1
             zua(:,:) = zua(:,:) + fse3u_a(:,:,jk) * ua(:,:,jk) * umask(:,:,jk)
             zva(:,:) = zva(:,:) + fse3v_a(:,:,jk) * va(:,:,jk) * vmask(:,:,jk)    
-         END DO
-         DO jk = 1, jpkm1
-            ua(:,:,jk) = ( ua(:,:,jk) - zua(:,:) * hur_e(:,:) + ua_b(:,:) ) * umask(:,:,jk)
-            va(:,:,jk) = ( va(:,:,jk) - zva(:,:) * hvr_e(:,:) + va_b(:,:) ) * vmask(:,:,jk)
+            zhura(:,:) = zhura(:,:) + fse3u_a(:,:,jk) * umask(:,:,jk)
+            zhvra(:,:) = zhvra(:,:) + fse3v_a(:,:,jk) * vmask(:,:,jk)
+         END DO
+         zhura(:,:) = umask(:,:,1)/(zhura(:,:) + 1._wp - umask(:,:,1)) 
+         zhvra(:,:) = vmask(:,:,1)/(zhvra(:,:) + 1._wp - vmask(:,:,1))  
+
+         DO jk = 1, jpkm1
+            ua(:,:,jk) = ( ua(:,:,jk) - zua(:,:) * zhura(:,:) + ua_b(:,:) ) * umask(:,:,jk)
+            va(:,:,jk) = ( va(:,:,jk) - zva(:,:) * zhvra(:,:) + va_b(:,:) ) * vmask(:,:,jk)
          END DO
       ELSE
@@ -185 +192 @@
          ! prior to asselin filtering     
          ! In the forward case, this is done below after asselin filtering    
+         ! so that asselin contribution is removed at the same time
          DO jk = 1, jpkm1
             un(:,:,jk) = ( un(:,:,jk) - un_adv(:,:) + un_b(:,:) )*umask(:,:,jk)
@@ -333 +341 @@
          !
          IF (lk_dynspg_ts.AND.ln_bt_fw) THEN
-         ! Remove asselin filtering of barotropic velocities if forward time splitting
-         ! note that we replace barotropic velocities by advective velocities       
+         ! Revert "before" velocities to time split estimate
+         ! Doing it here also means that asselin filter contribution is removed  
+         ! 
             zua(:,:) = 0._wp
             zva(:,:) = 0._wp
@@ -360 +369 @@
       ! 
       CALL wrk_dealloc( jpi,jpj,jpk, ze3u_f, ze3v_f )
-      IF ( lk_dynspg_ts ) CALL wrk_dealloc( jpi,jpj, zua, zva )
+      IF ( lk_dynspg_ts ) CALL wrk_dealloc( jpi,jpj, zua, zva , zhura, zhvra)
       !
       IF( nn_timing == 1 )  CALL timing_stop('dyn_nxt')

branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/step.F90

@@ -126 +126 @@
                                   CALL dyn_spg( kstp, indic )  ! surface pressure gradient
 
-                                  hdivb(:,:,:) = hdivn(:,:,:)  ! Store now divergence and rot temporarly, revert to these below 
-                                  rotb(:,:,:)  = rotn(:,:,:)     
-                                  ua_sv(:,:,:) = ua(:,:,:)     ! Save trends (barotropic trend has been fully updated)
+                                  ua_sv(:,:,:) = ua(:,:,:)     ! Save trends (barotropic trend has been fully updated at this stage)
                                   va_sv(:,:,:) = va(:,:,:)
 
@@ -134 +132 @@
       ENDIF
       IF( lk_vvl     )   CALL dom_vvl_sf_nxt( kstp )  ! after vertical scale factors
-                         CALL wzv           ( kstp )  ! now cross-level velocity (original)
+                         CALL wzv           ( kstp )  ! now cross-level velocity
 
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>