Changeset 9987 for branches/UKMO/dev_r5518_obs_oper_update_icethick/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

Timestamp:

2018-07-23T11:33:03+02:00 (6 years ago)

Author:

emmafiedler

Message:

Merge with GO6 FOAMv14 package branch r9288

File:

• branches/UKMO/dev_r5518_obs_oper_update_icethick/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 (modified) (11 diffs)

branches/UKMO/dev_r5518_obs_oper_update_icethick/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

@@ -39 +39 @@
    USE zdfmxl          ! mixed layer
    USE dianam          ! build name of file (routine)
+   USE zdftke          ! vertical physics: one-equation scheme 
    USE zdfddm          ! vertical  physics: double diffusion
    USE diahth          ! thermocline diagnostics
@@ -46 +47 @@
    USE iom
    USE ioipsl
-   USE dynspg_oce, ONLY: un_adv, vn_adv ! barotropic velocities     
-
+   USE dynspg_oce, ONLY: un_adv, vn_adv ! barotropic velocities    
+   USE insitu_tem, ONLY: insitu_t, theta2t
 #if defined key_lim2
    USE limwri_2 
@@ -145 +146 @@
       ENDIF
 
-      IF( .NOT.lk_vvl ) THEN
-         CALL iom_put( "e3t" , fse3t_n(:,:,:) )
-         CALL iom_put( "e3u" , fse3u_n(:,:,:) )
-         CALL iom_put( "e3v" , fse3v_n(:,:,:) )
-         CALL iom_put( "e3w" , fse3w_n(:,:,:) )
-      ENDIF
+      ! Output of initial vertical scale factor
+      CALL iom_put("e3t_0", e3t_0(:,:,:) )
+      CALL iom_put("e3u_0", e3t_0(:,:,:) )
+      CALL iom_put("e3v_0", e3t_0(:,:,:) )
+      !
+      CALL iom_put( "e3t" , fse3t_n(:,:,:) )
+      CALL iom_put( "e3u" , fse3u_n(:,:,:) )
+      CALL iom_put( "e3v" , fse3v_n(:,:,:) )
+      CALL iom_put( "e3w" , fse3w_n(:,:,:) )
+      IF( iom_use("e3tdef") )   &
+         CALL iom_put( "e3tdef"  , ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 )
+      CALL iom_put("tpt_dep", fsdept_n(:,:,:) )
+      CALL iom_put("wpt_dep", fsdepw_n(:,:,:) )
+
 
       CALL iom_put( "ssh" , sshn )                 ! sea surface height
-      if( iom_use('ssh2') )   CALL iom_put( "ssh2", sshn(:,:) * sshn(:,:) )   ! square of sea surface height
       
       CALL iom_put( "toce", tsn(:,:,:,jp_tem) )    ! 3D temperature
+      CALL theta2t ! in-situ temperature conversion
+      CALL iom_put( "tinsitu", insitu_t(:,:,:))    ! in-situ temperature
       CALL iom_put(  "sst", tsn(:,:,1,jp_tem) )    ! surface temperature
       IF ( iom_use("sbt") ) THEN
@@ -194 +204 @@
          CALL iom_put( "taubot", z2d )           
       ENDIF
-         
+      
       CALL iom_put( "uoce", un(:,:,:)         )    ! 3D i-current
       CALL iom_put(  "ssu", un(:,:,1)         )    ! surface i-current
@@ -242 +252 @@
       CALL iom_put( "avt" , avt                        )    ! T vert. eddy diff. coef.
       CALL iom_put( "avm" , avmu                       )    ! T vert. eddy visc. coef.
+      IF( lk_zdftke ) THEN   
+         CALL iom_put( "tke"      , en                               )    ! TKE budget: Turbulent Kinetic Energy   
+         CALL iom_put( "tke_niw"  , e_niw                            )    ! TKE budget: Near-inertial waves   
+      ENDIF 
       CALL iom_put( "avs" , fsavs(:,:,:)               )    ! S vert. eddy diff. coef. (useful only with key_zdfddm)
+                                                            ! Log of eddy diff coef
+      IF( iom_use('logavt') )   CALL iom_put( "logavt", LOG( MAX( 1.e-20_wp, avt  (:,:,:) ) ) )
+      IF( iom_use('logavs') )   CALL iom_put( "logavs", LOG( MAX( 1.e-20_wp, fsavs(:,:,:) ) ) )
 
       IF ( iom_use("sstgrad") .OR. iom_use("sstgrad2") ) THEN
@@ -307 +324 @@
          CALL iom_put( "eken", rke )           
       ENDIF
-         
-      IF( iom_use("u_masstr") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN
+      !
+      CALL iom_put( "hdiv", hdivn )                  ! Horizontal divergence
+      !
+      IF( iom_use("u_masstr") .OR. iom_use("u_masstr_vint") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN
          z3d(:,:,jpk) = 0.e0
+         z2d(:,:) = 0.e0
          DO jk = 1, jpkm1
             z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * fse3u(:,:,jk) * umask(:,:,jk)
+            z2d(:,:) = z2d(:,:) + z3d(:,:,jk)
          END DO
          CALL iom_put( "u_masstr", z3d )                  ! mass transport in i-direction
+         CALL iom_put( "u_masstr_vint", z2d )             ! mass transport in i-direction vertical sum
       ENDIF
       
@@ -376 +398 @@
          CALL iom_put( "v_salttr", 0.5 * z2d )            !  heat transport in j-direction
       ENDIF
+
+      ! Vertical integral of temperature
+      IF( iom_use("tosmint") ) THEN
+         z2d(:,:)=0._wp
+         DO jk = 1, jpkm1
+            DO jj = 2, jpjm1
+               DO ji = fs_2, fs_jpim1   ! vector opt.
+                  z2d(ji,jj) = z2d(ji,jj) + rau0 * fse3t(ji,jj,jk) *  tsn(ji,jj,jk,jp_tem)
+               END DO
+            END DO
+         END DO
+         CALL lbc_lnk( z2d, 'T', -1. )
+         CALL iom_put( "tosmint", z2d ) 
+      ENDIF
+
+      ! Vertical integral of salinity
+      IF( iom_use("somint") ) THEN
+         z2d(:,:)=0._wp
+         DO jk = 1, jpkm1
+            DO jj = 2, jpjm1
+               DO ji = fs_2, fs_jpim1   ! vector opt.
+                  z2d(ji,jj) = z2d(ji,jj) + rau0 * fse3t(ji,jj,jk) * tsn(ji,jj,jk,jp_sal)
+               END DO
+            END DO
+         END DO
+         CALL lbc_lnk( z2d, 'T', -1. )
+         CALL iom_put( "somint", z2d ) 
+      ENDIF
+
+      CALL iom_put( "bn2", rn2 )  !Brunt-Vaisala buoyancy frequency (N^2)
       !
       CALL wrk_dealloc( jpi , jpj      , z2d )
@@ -438 +490 @@
       zdt = rdt
       IF( nacc == 1 ) zdt = rdtmin
-      IF( ln_mskland )   THEN   ;   clop = "only(x)"   ! put 1.e+20 on land (very expensive!!)
-      ELSE                      ;   clop = "x"         ! no use of the mask value (require less cpu time)
-      ENDIF
+      clop = "x"         ! no use of the mask value (require less cpu time, and otherwise the model crashes)
 #if defined key_diainstant
       zsto = nwrite * zdt
@@ -1020 +1070 @@
          CALL histdef( id_i, "vovvldep", "T point depth"         , "m"      ,   &   ! t-point depth
             &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
+         CALL histdef( id_i, "vovvle3t", "T point thickness"         , "m"      ,   &   ! t-point depth
+            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
       END IF
 
@@ -1050 +1102 @@
       CALL histwrite( id_i, "sozotaux", kt, utau             , jpi*jpj    , idex )    ! i-wind stress
       CALL histwrite( id_i, "sometauy", kt, vtau             , jpi*jpj    , idex )    ! j-wind stress
+      IF( lk_vvl ) THEN
+         CALL histwrite( id_i, "vovvldep", kt, fsdept_n(:,:,:), jpi*jpj*jpk, idex )!  T-cell depth       
+         CALL histwrite( id_i, "vovvle3t", kt, fse3t_n (:,:,:), jpi*jpj*jpk, idex )!  T-cell thickness  
+      END IF
 
       ! 3. Close the file
@@ -1062 +1118 @@
       ENDIF
 #endif
+
+      IF (cdfile_name == "output.abort") THEN
+         CALL ctl_stop('MPPSTOP', 'NEMO abort from dia_wri_state')
+      END IF
        
 !     IF( nn_timing == 1 )   CALL timing_stop('dia_wri_state') ! not sure this works for routines not called in first timestep