Changeset 8102 for branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90

Timestamp:

2017-05-31T11:08:52+02:00 (8 years ago)

Author:

davestorkey

Message:

NEMO 3.6_stable: tracer trends diagnostics. See ticket #1877 for more details.

1. Correct bugs in calculation of total trends and trends due to vertical diffusion.
2. Output component trends every second timestep so that sum of component trends plus Asselin filter trend equals total trend.
3. Layer-integrated versions of trends (as per CMIP6 definition) available in field_def.xml.

File:

• branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90 (modified) (8 diffs)

branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90

@@ -129 +129 @@
       IF( l_trdtra )   THEN                    
          CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds )
-         ztrdt(:,:,jk) = 0._wp
-         ztrds(:,:,jk) = 0._wp
+         ztrdt(:,:,jpk) = 0._wp
+         ztrds(:,:,jpk) = 0._wp
          IF( ln_traldf_iso ) THEN              ! diagnose the "pure" Kz diffusive trend 
             CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdfp, ztrdt )
             CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdfp, ztrds )
          ENDIF
-         ! total trend for the non-time-filtered variables. 
-         DO jk = 1, jpkm1
-            zfact = 1.0 / rdttra(jk)
-            ztrdt(:,:,jk) = ( tsa(:,:,jk,jp_tem) - tsn(:,:,jk,jp_tem) ) * zfact 
-            ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal) - tsn(:,:,jk,jp_sal) ) * zfact 
-         END DO
+         ! total trend for the non-time-filtered variables.
+         ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3t*T)/e3tn; e3tn cancel from tsn terms
+         IF( lk_vvl ) THEN
+            DO jk = 1, jpkm1
+               zfact = 1.0 / rdttra(jk)
+               ztrdt(:,:,jk) = ( tsa(:,:,jk,jp_tem)*fse3t_a(:,:,jk) / fse3t_n(:,:,jk) - tsn(:,:,jk,jp_tem)) * zfact
+               ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal)*fse3t_a(:,:,jk) / fse3t_n(:,:,jk) - tsn(:,:,jk,jp_sal)) * zfact
+            END DO
+         ELSE
+            DO jk = 1, jpkm1
+               zfact = 1.0 / rdttra(jk)
+               ztrdt(:,:,jk) = ( tsa(:,:,jk,jp_tem) - tsn(:,:,jk,jp_tem) ) * zfact 
+               ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal) - tsn(:,:,jk,jp_sal) ) * zfact 
+            END DO
+         END IF
          CALL trd_tra( kt, 'TRA', jp_tem, jptra_tot, ztrdt )
          CALL trd_tra( kt, 'TRA', jp_sal, jptra_tot, ztrds )
-         ! Store now fields before applying the Asselin filter 
-         ! in order to calculate Asselin filter trend later.
-         ztrdt(:,:,:) = tsn(:,:,:,jp_tem) 
-         ztrds(:,:,:) = tsn(:,:,:,jp_sal)
+         IF( .NOT.lk_vvl )  THEN
+            ! Store now fields before applying the Asselin filter 
+            ! in order to calculate Asselin filter trend later.
+            ztrdt(:,:,:) = tsn(:,:,:,jp_tem) 
+            ztrds(:,:,:) = tsn(:,:,:,jp_sal)
+         END IF
       ENDIF
 
@@ -155 +166 @@
             END DO
          END DO
+         IF (l_trdtra.AND.lk_vvl) THEN      ! Zero Asselin filter contribution must be explicitly written out since for vvl
+                                            ! Asselin filter is output by tra_nxt_vvl that is not called on this time step
+            ztrdt(:,:,:) = 0._wp
+            ztrds(:,:,:) = 0._wp
+            CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt )
+            CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds )
+         END IF
       ELSE                                            ! Leap-Frog + Asselin filter time stepping
          !
@@ -164 +182 @@
       !
      ! trends computation
-      IF( l_trdtra ) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt     
+      IF( l_trdtra.AND..NOT.lk_vvl) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt     
          DO jk = 1, jpkm1
             zfact = 1._wp / r2dtra(jk)             
@@ -172 +190 @@
          CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt )
          CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds )
-         CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds )
       END IF
+      IF( l_trdtra) CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds )
       !
       !                        ! control print
@@ -289 +307 @@
       LOGICAL  ::   ll_tra_hpg, ll_traqsr, ll_rnf, ll_isf   ! local logical
       INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices
-      REAL(wp) ::   zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar
+      REAL(wp) ::   zfact, zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar
       REAL(wp) ::   zfact2, ze3t_b, ze3t_n, ze3t_a, ze3t_f, ze3t_d   !   -      -
+      REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztrd_atf
       !!----------------------------------------------------------------------
       !
@@ -315 +334 @@
       ENDIF
       !
+      IF( l_trdtra )   THEN
+         CALL wrk_alloc( jpi, jpj, jpk, kjpt, ztrd_atf )
+         ztrd_atf(:,:,:,:) = 0.0_wp
+      ENDIF
       DO jn = 1, kjpt      
          DO jk = 1, jpkm1
+            zfact = 1._wp / r2dtra(jk)
             zfact1 = atfp * p2dt(jk)
             zfact2 = zfact1 / rau0
@@ -371 +395 @@
                      pta(ji,jj,jk,jn) = ze3t_d * ( ztc_n  + rbcp * ztc_d  )   ! ta <-- Brown & Campana average
                   ENDIF
+                  IF( l_trdtra ) THEN
+                     ztrd_atf(ji,jj,jk,jn) = (ztc_f - ztc_n) * zfact/ze3t_n
+                  ENDIF
                END DO
             END DO
@@ -377 +404 @@
       END DO
       !
+      IF( l_trdtra ) THEN
+         CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrd_atf(:,:,:,jp_tem) )
+         CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrd_atf(:,:,:,jp_sal) )
+         CALL wrk_dealloc( jpi, jpj, jpk, kjpt, ztrd_atf )
+      ENDIF
+
    END SUBROUTINE tra_nxt_vvl