New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 13463 for NEMO/branches/2019/dev_r11351_fldread_with_XIOS/tests/CANAL/MY_SRC/diawri.F90 – NEMO

Ignore:
Timestamp:
2020-09-14T17:40:34+02:00 (4 years ago)
Author:
andmirek
Message:

Ticket #2195:update to trunk 13461

Location:
NEMO/branches/2019/dev_r11351_fldread_with_XIOS
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11351_fldread_with_XIOS

    • Property svn:externals
      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
         8 
         9# SETTE 
         10^/utils/CI/sette@13382        sette 
  • NEMO/branches/2019/dev_r11351_fldread_with_XIOS/tests/CANAL/MY_SRC/diawri.F90

    r10425 r13463  
    2626   !!---------------------------------------------------------------------- 
    2727   USE oce            ! ocean dynamics and tracers  
     28   USE isf_oce 
     29   USE isfcpl 
     30   USE abl            ! abl variables in case ln_abl = .true. 
    2831   USE dom_oce        ! ocean space and time domain 
    2932   USE phycst         ! physical constants 
     
    4649   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    4750   USE in_out_manager ! I/O manager 
    48    USE diatmb         ! Top,middle,bottom output 
    4951   USE dia25h         ! 25h Mean output 
    5052   USE iom            !  
     
    5759   USE lib_mpp         ! MPP library 
    5860   USE timing          ! preformance summary 
    59    USE diurnal_bulk    ! diurnal warm layer 
    60    USE cool_skin       ! Cool skin 
     61   USE diu_bulk        ! diurnal warm layer 
     62   USE diu_coolskin    ! Cool skin 
    6163 
    6264   IMPLICIT NONE 
     
    6668   PUBLIC   dia_wri_state 
    6769   PUBLIC   dia_wri_alloc           ! Called by nemogcm module 
    68  
     70#if ! defined key_iomput    
     71   PUBLIC   dia_wri_alloc_abl       ! Called by sbcabl  module (if ln_abl = .true.) 
     72#endif 
    6973   INTEGER ::   nid_T, nz_T, nh_T, ndim_T, ndim_hT   ! grid_T file 
    7074   INTEGER ::          nb_T              , ndim_bT   ! grid_T file 
     
    7276   INTEGER ::   nid_V, nz_V, nh_V, ndim_V, ndim_hV   ! grid_V file 
    7377   INTEGER ::   nid_W, nz_W, nh_W                    ! grid_W file 
     78   INTEGER ::   nid_A, nz_A, nh_A, ndim_A, ndim_hA   ! grid_ABL file    
    7479   INTEGER ::   ndex(1)                              ! ??? 
    7580   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV 
     81   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hA, ndex_A ! ABL 
    7682   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V 
    7783   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_bT 
    7884 
    7985   !! * Substitutions 
    80 #  include "vectopt_loop_substitute.h90" 
     86#  include "do_loop_substitute.h90" 
    8187   !!---------------------------------------------------------------------- 
    8288   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    97103 
    98104    
    99    SUBROUTINE dia_wri( kt ) 
     105   SUBROUTINE dia_wri( kt, Kmm ) 
    100106      !!--------------------------------------------------------------------- 
    101107      !!                  ***  ROUTINE dia_wri  *** 
     
    107113      !!---------------------------------------------------------------------- 
    108114      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
     115      INTEGER, INTENT( in ) ::   Kmm     ! ocean time level index 
    109116      !! 
    110117      INTEGER ::   ji, jj, jk       ! dummy loop indices 
     
    115122      REAL(wp), DIMENSION(jpi,jpj)     ::   z2d   ! 2D workspace 
    116123      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   z3d   ! 3D workspace 
    117       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   bu, bv   ! volume of u- and v-boxes 
    118       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   r1_bt    ! inverse of t-box volume 
    119124      !!---------------------------------------------------------------------- 
    120125      !  
     
    123128      ! Output the initial state and forcings 
    124129      IF( ninist == 1 ) THEN                        
    125          CALL dia_wri_state( 'output.init' ) 
     130         CALL dia_wri_state( Kmm, 'output.init' ) 
    126131         ninist = 0 
    127132      ENDIF 
     
    132137      CALL iom_put("e3v_0", e3v_0(:,:,:) ) 
    133138      ! 
    134       CALL iom_put( "e3t" , e3t_n(:,:,:) ) 
    135       CALL iom_put( "e3u" , e3u_n(:,:,:) ) 
    136       CALL iom_put( "e3v" , e3v_n(:,:,:) ) 
    137       CALL iom_put( "e3w" , e3w_n(:,:,:) ) 
     139      CALL iom_put( "e3t" , e3t(:,:,:,Kmm) ) 
     140      CALL iom_put( "e3u" , e3u(:,:,:,Kmm) ) 
     141      CALL iom_put( "e3v" , e3v(:,:,:,Kmm) ) 
     142      CALL iom_put( "e3w" , e3w(:,:,:,Kmm) ) 
    138143      IF( iom_use("e3tdef") )   & 
    139          CALL iom_put( "e3tdef"  , ( ( e3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 ) 
     144         CALL iom_put( "e3tdef"  , ( ( e3t(:,:,:,Kmm) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 ) 
    140145 
    141146      IF( ll_wd ) THEN 
    142          CALL iom_put( "ssh" , (sshn+ssh_ref)*tmask(:,:,1) )   ! sea surface height (brought back to the reference used for wetting and drying) 
     147         CALL iom_put( "ssh" , (ssh(:,:,Kmm)+ssh_ref)*tmask(:,:,1) )   ! sea surface height (brought back to the reference used for wetting and drying) 
    143148      ELSE 
    144          CALL iom_put( "ssh" , sshn )              ! sea surface height 
     149         CALL iom_put( "ssh" , ssh(:,:,Kmm) )              ! sea surface height 
    145150      ENDIF 
    146151 
    147152      IF( iom_use("wetdep") )   &                  ! wet depth 
    148          CALL iom_put( "wetdep" , ht_0(:,:) + sshn(:,:) ) 
     153         CALL iom_put( "wetdep" , ht_0(:,:) + ssh(:,:,Kmm) ) 
    149154       
    150       CALL iom_put( "toce", tsn(:,:,:,jp_tem) )    ! 3D temperature 
    151       CALL iom_put(  "sst", tsn(:,:,1,jp_tem) )    ! surface temperature 
     155      CALL iom_put( "toce", ts(:,:,:,jp_tem,Kmm) )    ! 3D temperature 
     156      CALL iom_put(  "sst", ts(:,:,1,jp_tem,Kmm) )    ! surface temperature 
    152157      IF ( iom_use("sbt") ) THEN 
    153          DO jj = 1, jpj 
    154             DO ji = 1, jpi 
    155                ikbot = mbkt(ji,jj) 
    156                z2d(ji,jj) = tsn(ji,jj,ikbot,jp_tem) 
    157             END DO 
    158          END DO 
     158         DO_2D( 1, 1, 1, 1 ) 
     159            ikbot = mbkt(ji,jj) 
     160            z2d(ji,jj) = ts(ji,jj,ikbot,jp_tem,Kmm) 
     161         END_2D 
    159162         CALL iom_put( "sbt", z2d )                ! bottom temperature 
    160163      ENDIF 
    161164       
    162       CALL iom_put( "soce", tsn(:,:,:,jp_sal) )    ! 3D salinity 
    163       CALL iom_put(  "sss", tsn(:,:,1,jp_sal) )    ! surface salinity 
     165      CALL iom_put( "soce", ts(:,:,:,jp_sal,Kmm) )    ! 3D salinity 
     166      CALL iom_put(  "sss", ts(:,:,1,jp_sal,Kmm) )    ! surface salinity 
    164167      IF ( iom_use("sbs") ) THEN 
    165          DO jj = 1, jpj 
    166             DO ji = 1, jpi 
    167                ikbot = mbkt(ji,jj) 
    168                z2d(ji,jj) = tsn(ji,jj,ikbot,jp_sal) 
    169             END DO 
    170          END DO 
     168         DO_2D( 1, 1, 1, 1 ) 
     169            ikbot = mbkt(ji,jj) 
     170            z2d(ji,jj) = ts(ji,jj,ikbot,jp_sal,Kmm) 
     171         END_2D 
    171172         CALL iom_put( "sbs", z2d )                ! bottom salinity 
    172173      ENDIF 
    173174 
    174175      IF ( iom_use("taubot") ) THEN                ! bottom stress 
    175          zztmp = rau0 * 0.25 
     176         zztmp = rho0 * 0.25 
    176177         z2d(:,:) = 0._wp 
    177          DO jj = 2, jpjm1 
    178             DO ji = fs_2, fs_jpim1   ! vector opt. 
    179                zztmp2 = (  ( rCdU_bot(ji+1,jj)+rCdU_bot(ji  ,jj) ) * un(ji  ,jj,mbku(ji  ,jj))  )**2   & 
    180                   &   + (  ( rCdU_bot(ji  ,jj)+rCdU_bot(ji-1,jj) ) * un(ji-1,jj,mbku(ji-1,jj))  )**2   & 
    181                   &   + (  ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj  ) ) * vn(ji,jj  ,mbkv(ji,jj  ))  )**2   & 
    182                   &   + (  ( rCdU_bot(ji,jj  )+rCdU_bot(ji,jj-1) ) * vn(ji,jj-1,mbkv(ji,jj-1))  )**2 
    183                z2d(ji,jj) = zztmp * SQRT( zztmp2 ) * tmask(ji,jj,1)  
    184                ! 
    185             END DO 
    186          END DO 
     178         DO_2D( 0, 0, 0, 0 ) 
     179            zztmp2 = (  ( rCdU_bot(ji+1,jj)+rCdU_bot(ji  ,jj) ) * uu(ji  ,jj,mbku(ji  ,jj),Kmm)  )**2   & 
     180               &   + (  ( rCdU_bot(ji  ,jj)+rCdU_bot(ji-1,jj) ) * uu(ji-1,jj,mbku(ji-1,jj),Kmm)  )**2   & 
     181               &   + (  ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj  ) ) * vv(ji,jj  ,mbkv(ji,jj  ),Kmm)  )**2   & 
     182               &   + (  ( rCdU_bot(ji,jj  )+rCdU_bot(ji,jj-1) ) * vv(ji,jj-1,mbkv(ji,jj-1),Kmm)  )**2 
     183            z2d(ji,jj) = zztmp * SQRT( zztmp2 ) * tmask(ji,jj,1)  
     184            ! 
     185         END_2D 
    187186         CALL lbc_lnk( 'diawri', z2d, 'T', 1. ) 
    188187         CALL iom_put( "taubot", z2d )            
    189188      ENDIF 
    190189          
    191       CALL iom_put( "uoce", un(:,:,:) )            ! 3D i-current 
    192       CALL iom_put(  "ssu", un(:,:,1) )            ! surface i-current 
     190      CALL iom_put( "uoce", uu(:,:,:,Kmm) )            ! 3D i-current 
     191      CALL iom_put(  "ssu", uu(:,:,1,Kmm) )            ! surface i-current 
    193192      IF ( iom_use("sbu") ) THEN 
    194          DO jj = 1, jpj 
    195             DO ji = 1, jpi 
    196                ikbot = mbku(ji,jj) 
    197                z2d(ji,jj) = un(ji,jj,ikbot) 
    198             END DO 
    199          END DO 
     193         DO_2D( 1, 1, 1, 1 ) 
     194            ikbot = mbku(ji,jj) 
     195            z2d(ji,jj) = uu(ji,jj,ikbot,Kmm) 
     196         END_2D 
    200197         CALL iom_put( "sbu", z2d )                ! bottom i-current 
    201198      ENDIF 
    202199       
    203       CALL iom_put( "voce", vn(:,:,:) )            ! 3D j-current 
    204       CALL iom_put(  "ssv", vn(:,:,1) )            ! surface j-current 
     200      CALL iom_put( "voce", vv(:,:,:,Kmm) )            ! 3D j-current 
     201      CALL iom_put(  "ssv", vv(:,:,1,Kmm) )            ! surface j-current 
    205202      IF ( iom_use("sbv") ) THEN 
    206          DO jj = 1, jpj 
    207             DO ji = 1, jpi 
    208                ikbot = mbkv(ji,jj) 
    209                z2d(ji,jj) = vn(ji,jj,ikbot) 
    210             END DO 
    211          END DO 
     203         DO_2D( 1, 1, 1, 1 ) 
     204            ikbot = mbkv(ji,jj) 
     205            z2d(ji,jj) = vv(ji,jj,ikbot,Kmm) 
     206         END_2D 
    212207         CALL iom_put( "sbv", z2d )                ! bottom j-current 
    213208      ENDIF 
    214209 
    215       CALL iom_put( "woce", wn )                   ! vertical velocity 
     210      IF( ln_zad_Aimp ) ww = ww + wi               ! Recombine explicit and implicit parts of vertical velocity for diagnostic output 
     211      ! 
     212      CALL iom_put( "woce", ww )                   ! vertical velocity 
    216213      IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN   ! vertical mass transport & its square value 
    217214         ! Caution: in the VVL case, it only correponds to the baroclinic mass transport. 
    218          z2d(:,:) = rau0 * e1e2t(:,:) 
     215         z2d(:,:) = rho0 * e1e2t(:,:) 
    219216         DO jk = 1, jpk 
    220             z3d(:,:,jk) = wn(:,:,jk) * z2d(:,:) 
     217            z3d(:,:,jk) = ww(:,:,jk) * z2d(:,:) 
    221218         END DO 
    222219         CALL iom_put( "w_masstr" , z3d )   
    223220         IF( iom_use('w_masstr2') )   CALL iom_put( "w_masstr2", z3d(:,:,:) * z3d(:,:,:) ) 
    224221      ENDIF 
     222      ! 
     223      IF( ln_zad_Aimp ) ww = ww - wi               ! Remove implicit part of vertical velocity that was added for diagnostic output 
    225224 
    226225      CALL iom_put( "avt" , avt )                  ! T vert. eddy diff. coef. 
     
    231230      IF( iom_use('logavs') )   CALL iom_put( "logavs", LOG( MAX( 1.e-20_wp, avs(:,:,:) ) ) ) 
    232231 
    233       IF ( iom_use("salgrad") .OR. iom_use("salgrad2") ) THEN 
    234          z3d(:,:,jpk) = 0. 
    235          DO jk = 1, jpkm1 
    236             DO jj = 2, jpjm1                                    ! sal gradient 
    237                DO ji = fs_2, fs_jpim1   ! vector opt. 
    238                   zztmp  = tsn(ji,jj,jk,jp_sal) 
    239                   zztmpx = ( tsn(ji+1,jj,jk,jp_sal) - zztmp ) * r1_e1u(ji,jj) + ( zztmp - tsn(ji-1,jj  ,jk,jp_sal) ) * r1_e1u(ji-1,jj) 
    240                   zztmpy = ( tsn(ji,jj+1,jk,jp_sal) - zztmp ) * r1_e2v(ji,jj) + ( zztmp - tsn(ji  ,jj-1,jk,jp_sal) ) * r1_e2v(ji,jj-1) 
    241                   z3d(ji,jj,jk) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy )   & 
    242                      &                 * umask(ji,jj,jk) * umask(ji-1,jj,jk) * vmask(ji,jj,jk) * umask(ji,jj-1,jk) 
    243                END DO 
    244             END DO 
    245          END DO 
    246          CALL lbc_lnk( 'diawri', z3d, 'T', 1. ) 
    247          CALL iom_put( "salgrad2",  z3d )          ! square of module of sal gradient 
    248          z3d(:,:,:) = SQRT( z3d(:,:,:) ) 
    249          CALL iom_put( "salgrad" ,  z3d )          ! module of sal gradient 
    250       ENDIF 
    251           
    252232      IF ( iom_use("sstgrad") .OR. iom_use("sstgrad2") ) THEN 
    253          DO jj = 2, jpjm1                                    ! sst gradient 
    254             DO ji = fs_2, fs_jpim1   ! vector opt. 
    255                zztmp  = tsn(ji,jj,1,jp_tem) 
    256                zztmpx = ( tsn(ji+1,jj,1,jp_tem) - zztmp ) * r1_e1u(ji,jj) + ( zztmp - tsn(ji-1,jj  ,1,jp_tem) ) * r1_e1u(ji-1,jj) 
    257                zztmpy = ( tsn(ji,jj+1,1,jp_tem) - zztmp ) * r1_e2v(ji,jj) + ( zztmp - tsn(ji  ,jj-1,1,jp_tem) ) * r1_e2v(ji,jj-1) 
    258                z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy )   & 
    259                   &              * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1) 
    260             END DO 
    261          END DO 
     233         DO_2D( 0, 0, 0, 0 ) 
     234            zztmp  = ts(ji,jj,1,jp_tem,Kmm) 
     235            zztmpx = ( ts(ji+1,jj,1,jp_tem,Kmm) - zztmp ) * r1_e1u(ji,jj) + ( zztmp - ts(ji-1,jj  ,1,jp_tem,Kmm) ) * r1_e1u(ji-1,jj) 
     236            zztmpy = ( ts(ji,jj+1,1,jp_tem,Kmm) - zztmp ) * r1_e2v(ji,jj) + ( zztmp - ts(ji  ,jj-1,1,jp_tem,Kmm) ) * r1_e2v(ji,jj-1) 
     237            z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy )   & 
     238               &              * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1) 
     239         END_2D 
    262240         CALL lbc_lnk( 'diawri', z2d, 'T', 1. ) 
    263241         CALL iom_put( "sstgrad2",  z2d )          ! square of module of sst gradient 
     
    269247      IF( iom_use("heatc") ) THEN 
    270248         z2d(:,:)  = 0._wp  
    271          DO jk = 1, jpkm1 
    272             DO jj = 1, jpj 
    273                DO ji = 1, jpi 
    274                   z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) 
    275                END DO 
    276             END DO 
    277          END DO 
    278          CALL iom_put( "heatc", rau0_rcp * z2d )   ! vertically integrated heat content (J/m2) 
     249         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     250            z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_tem,Kmm) * tmask(ji,jj,jk) 
     251         END_3D 
     252         CALL iom_put( "heatc", rho0_rcp * z2d )   ! vertically integrated heat content (J/m2) 
    279253      ENDIF 
    280254 
    281255      IF( iom_use("saltc") ) THEN 
    282256         z2d(:,:)  = 0._wp  
    283          DO jk = 1, jpkm1 
    284             DO jj = 1, jpj 
    285                DO ji = 1, jpi 
    286                   z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk) 
    287                END DO 
    288             END DO 
    289          END DO 
    290          CALL iom_put( "saltc", rau0 * z2d )          ! vertically integrated salt content (PSU*kg/m2) 
     257         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     258            z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) * tmask(ji,jj,jk) 
     259         END_3D 
     260         CALL iom_put( "saltc", rho0 * z2d )          ! vertically integrated salt content (PSU*kg/m2) 
    291261      ENDIF 
    292262      ! 
    293263      IF( iom_use("salt2c") ) THEN 
    294264         z2d(:,:)  = 0._wp  
    295          DO jk = 1, jpkm1 
    296             DO jj = 1, jpj 
    297                DO ji = 1, jpi 
    298                   z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) * tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk) 
    299                END DO 
    300             END DO 
    301          END DO 
    302          CALL iom_put( "salt2c", rau0 * z2d )          ! vertically integrated salt content (PSU*kg/m2) 
     265         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     266            z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) * tmask(ji,jj,jk) 
     267         END_3D 
     268         CALL iom_put( "salt2c", rho0 * z2d )          ! vertically integrated salt content (PSU*kg/m2) 
    303269      ENDIF 
    304270      ! 
    305271      IF ( iom_use("eken") ) THEN 
    306272         z3d(:,:,jpk) = 0._wp  
    307          DO jk = 1, jpkm1 
    308             DO jj = 2, jpj 
    309                DO ji = 2, jpi 
    310                   zztmpx = 0.5 * ( un(ji-1,jj  ,jk) + un(ji,jj,jk) ) 
    311                   zztmpy = 0.5 * ( vn(ji  ,jj-1,jk) + vn(ji,jj,jk) ) 
    312                   z3d(ji,jj,jk) = 0.5 * ( zztmpx*zztmpx + zztmpy*zztmpy ) 
    313                END DO 
    314             END DO 
    315          END DO 
     273         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     274            zztmp  = 0.25_wp * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     275            z3d(ji,jj,jk) = zztmp * (  uu(ji-1,jj,jk,Kmm)**2 * e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm)   & 
     276               &                     + uu(ji  ,jj,jk,Kmm)**2 * e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm)   & 
     277               &                     + vv(ji,jj-1,jk,Kmm)**2 * e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm)   & 
     278               &                     + vv(ji,jj  ,jk,Kmm)**2 * e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm)   ) 
     279         END_3D 
    316280         CALL lbc_lnk( 'diawri', z3d, 'T', 1. ) 
    317281         CALL iom_put( "eken", z3d )                 ! kinetic energy 
     
    323287         z3d(1,:, : ) = 0._wp 
    324288         z3d(:,1, : ) = 0._wp 
    325          DO jk = 1, jpkm1 
    326             DO jj = 2, jpj 
    327                DO ji = 2, jpi 
    328                   z3d(ji,jj,jk) = 0.25_wp * ( un(ji  ,jj,jk) * un(ji  ,jj,jk) * e1e2u(ji  ,jj) * e3u_n(ji  ,jj,jk)  & 
    329                      &                      + un(ji-1,jj,jk) * un(ji-1,jj,jk) * e1e2u(ji-1,jj) * e3u_n(ji-1,jj,jk)  & 
    330                      &                      + vn(ji,jj  ,jk) * vn(ji,jj  ,jk) * e1e2v(ji,jj  ) * e3v_n(ji,jj  ,jk)  & 
    331                      &                      + vn(ji,jj-1,jk) * vn(ji,jj-1,jk) * e1e2v(ji,jj-1) * e3v_n(ji,jj-1,jk)  )  & 
    332                      &                    * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) * tmask(ji,jj,jk) 
    333                END DO 
    334             END DO 
    335          END DO 
    336           
     289         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     290            z3d(ji,jj,jk) = 0.25_wp * ( uu(ji  ,jj,jk,Kmm) * uu(ji  ,jj,jk,Kmm) * e1e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm)  & 
     291               &                      + uu(ji-1,jj,jk,Kmm) * uu(ji-1,jj,jk,Kmm) * e1e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm)  & 
     292               &                      + vv(ji,jj  ,jk,Kmm) * vv(ji,jj  ,jk,Kmm) * e1e2v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm)  & 
     293               &                      + vv(ji,jj-1,jk,Kmm) * vv(ji,jj-1,jk,Kmm) * e1e2v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm)  )  & 
     294               &                    * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) 
     295         END_3D 
    337296         CALL lbc_lnk( 'diawri', z3d, 'T', 1. ) 
    338297         CALL iom_put( "ke", z3d ) ! kinetic energy 
    339298 
    340299         z2d(:,:)  = 0._wp  
    341          DO jk = 1, jpkm1 
    342             DO jj = 1, jpj 
    343                DO ji = 1, jpi 
    344                   z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * z3d(ji,jj,jk) * tmask(ji,jj,jk) 
    345                END DO 
    346             END DO 
    347          END DO 
     300         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     301            z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * z3d(ji,jj,jk) * tmask(ji,jj,jk) 
     302         END_3D 
    348303         CALL iom_put( "ke_zint", z2d )   ! vertically integrated kinetic energy 
    349304 
    350305      ENDIF 
    351306      ! 
    352       CALL iom_put( "hdiv", hdivn )                  ! Horizontal divergence 
     307      CALL iom_put( "hdiv", hdiv )                  ! Horizontal divergence 
    353308 
    354309      IF ( iom_use("relvor") .OR. iom_use("absvor") .OR. iom_use("potvor") ) THEN 
    355310          
    356311         z3d(:,:,jpk) = 0._wp  
    357          DO jk = 1, jpkm1 
    358             DO jj = 1, jpjm1 
    359                DO ji = 1, fs_jpim1   ! vector opt. 
    360                   z3d(ji,jj,jk) = (  e2v(ji+1,jj  ) * vn(ji+1,jj  ,jk) - e2v(ji,jj) * vn(ji,jj,jk)    & 
    361                      &              - e1u(ji  ,jj+1) * un(ji  ,jj+1,jk) + e1u(ji,jj) * un(ji,jj,jk)  ) * r1_e1e2f(ji,jj) 
    362                END DO 
    363             END DO 
    364          END DO 
     312         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     313            z3d(ji,jj,jk) = (   e2v(ji+1,jj  ) * vv(ji+1,jj  ,jk,Kmm) - e2v(ji,jj) * vv(ji,jj,jk,Kmm)    & 
     314               &              - e1u(ji  ,jj+1) * uu(ji  ,jj+1,jk,Kmm) + e1u(ji,jj) * uu(ji,jj,jk,Kmm)  ) * r1_e1e2f(ji,jj) 
     315         END_3D 
    365316         CALL lbc_lnk( 'diawri', z3d, 'F', 1. ) 
    366317         CALL iom_put( "relvor", z3d )                  ! relative vorticity 
    367318 
    368          DO jk = 1, jpkm1 
    369             DO jj = 1, jpj 
    370                DO ji = 1, jpi 
    371                   z3d(ji,jj,jk) = ff_f(ji,jj) + z3d(ji,jj,jk)  
    372                END DO 
    373             END DO 
    374          END DO 
     319         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     320            z3d(ji,jj,jk) = ff_f(ji,jj) + z3d(ji,jj,jk)  
     321         END_3D 
    375322         CALL iom_put( "absvor", z3d )                  ! absolute vorticity 
    376323 
    377          DO jk = 1, jpkm1 
    378             DO jj = 1, jpjm1 
    379                DO ji = 1, fs_jpim1   ! vector opt. 
    380                   ze3  = (  e3t_n(ji,jj+1,jk)*tmask(ji,jj+1,jk) + e3t_n(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk)   & 
    381                      &    + e3t_n(ji,jj  ,jk)*tmask(ji,jj  ,jk) + e3t_n(ji+1,jj  ,jk)*tmask(ji+1,jj  ,jk)  ) 
    382                   IF( ze3 /= 0._wp ) THEN   ;   ze3 = 4._wp / ze3 
    383                   ELSE                      ;   ze3 = 0._wp 
    384                   ENDIF 
    385                   z3d(ji,jj,jk) = ze3 * z3d(ji,jj,jk)  
    386                END DO 
    387             END DO 
    388          END DO 
     324         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     325            ze3  = (  e3t(ji,jj+1,jk,Kmm)*tmask(ji,jj+1,jk) + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)   & 
     326               &    + e3t(ji,jj  ,jk,Kmm)*tmask(ji,jj  ,jk) + e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk)  ) 
     327            IF( ze3 /= 0._wp ) THEN   ;   ze3 = 4._wp / ze3 
     328            ELSE                      ;   ze3 = 0._wp 
     329            ENDIF 
     330            z3d(ji,jj,jk) = ze3 * z3d(ji,jj,jk)  
     331         END_3D 
    389332         CALL lbc_lnk( 'diawri', z3d, 'F', 1. ) 
    390333         CALL iom_put( "potvor", z3d )                  ! potential vorticity 
    391334 
    392335      ENDIF 
    393     
    394336      ! 
    395337      IF( iom_use("u_masstr") .OR. iom_use("u_masstr_vint") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN 
     
    397339         z2d(:,:) = 0.e0 
    398340         DO jk = 1, jpkm1 
    399             z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * e3u_n(:,:,jk) * umask(:,:,jk) 
     341            z3d(:,:,jk) = rho0 * uu(:,:,jk,Kmm) * e2u(:,:) * e3u(:,:,jk,Kmm) * umask(:,:,jk) 
    400342            z2d(:,:) = z2d(:,:) + z3d(:,:,jk) 
    401343         END DO 
     
    406348      IF( iom_use("u_heattr") ) THEN 
    407349         z2d(:,:) = 0._wp  
    408          DO jk = 1, jpkm1 
    409             DO jj = 2, jpjm1 
    410                DO ji = fs_2, fs_jpim1   ! vector opt. 
    411                   z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji+1,jj,jk,jp_tem) ) 
    412                END DO 
    413             END DO 
    414          END DO 
     350         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     351            z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_tem,Kmm) + ts(ji+1,jj,jk,jp_tem,Kmm) ) 
     352         END_3D 
    415353         CALL lbc_lnk( 'diawri', z2d, 'U', -1. ) 
    416354         CALL iom_put( "u_heattr", 0.5*rcp * z2d )    ! heat transport in i-direction 
     
    419357      IF( iom_use("u_salttr") ) THEN 
    420358         z2d(:,:) = 0.e0  
    421          DO jk = 1, jpkm1 
    422             DO jj = 2, jpjm1 
    423                DO ji = fs_2, fs_jpim1   ! vector opt. 
    424                   z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji+1,jj,jk,jp_sal) ) 
    425                END DO 
    426             END DO 
    427          END DO 
     359         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     360            z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_sal,Kmm) + ts(ji+1,jj,jk,jp_sal,Kmm) ) 
     361         END_3D 
    428362         CALL lbc_lnk( 'diawri', z2d, 'U', -1. ) 
    429363         CALL iom_put( "u_salttr", 0.5 * z2d )        ! heat transport in i-direction 
     
    434368         z3d(:,:,jpk) = 0.e0 
    435369         DO jk = 1, jpkm1 
    436             z3d(:,:,jk) = rau0 * vn(:,:,jk) * e1v(:,:) * e3v_n(:,:,jk) * vmask(:,:,jk) 
     370            z3d(:,:,jk) = rho0 * vv(:,:,jk,Kmm) * e1v(:,:) * e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
    437371         END DO 
    438372         CALL iom_put( "v_masstr", z3d )              ! mass transport in j-direction 
     
    441375      IF( iom_use("v_heattr") ) THEN 
    442376         z2d(:,:) = 0.e0  
    443          DO jk = 1, jpkm1 
    444             DO jj = 2, jpjm1 
    445                DO ji = fs_2, fs_jpim1   ! vector opt. 
    446                   z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_tem) + tsn(ji,jj+1,jk,jp_tem) ) 
    447                END DO 
    448             END DO 
    449          END DO 
     377         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     378            z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_tem,Kmm) + ts(ji,jj+1,jk,jp_tem,Kmm) ) 
     379         END_3D 
    450380         CALL lbc_lnk( 'diawri', z2d, 'V', -1. ) 
    451381         CALL iom_put( "v_heattr", 0.5*rcp * z2d )    !  heat transport in j-direction 
     
    454384      IF( iom_use("v_salttr") ) THEN 
    455385         z2d(:,:) = 0._wp  
    456          DO jk = 1, jpkm1 
    457             DO jj = 2, jpjm1 
    458                DO ji = fs_2, fs_jpim1   ! vector opt. 
    459                   z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( tsn(ji,jj,jk,jp_sal) + tsn(ji,jj+1,jk,jp_sal) ) 
    460                END DO 
    461             END DO 
    462          END DO 
     386         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     387            z2d(ji,jj) = z2d(ji,jj) + z3d(ji,jj,jk) * ( ts(ji,jj,jk,jp_sal,Kmm) + ts(ji,jj+1,jk,jp_sal,Kmm) ) 
     388         END_3D 
    463389         CALL lbc_lnk( 'diawri', z2d, 'V', -1. ) 
    464390         CALL iom_put( "v_salttr", 0.5 * z2d )        !  heat transport in j-direction 
     
    467393      IF( iom_use("tosmint") ) THEN 
    468394         z2d(:,:) = 0._wp 
    469          DO jk = 1, jpkm1 
    470             DO jj = 2, jpjm1 
    471                DO ji = fs_2, fs_jpim1   ! vector opt. 
    472                   z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) *  tsn(ji,jj,jk,jp_tem) 
    473                END DO 
    474             END DO 
    475          END DO 
     395         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     396            z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) *  ts(ji,jj,jk,jp_tem,Kmm) 
     397         END_3D 
    476398         CALL lbc_lnk( 'diawri', z2d, 'T', -1. ) 
    477          CALL iom_put( "tosmint", rau0 * z2d )        ! Vertical integral of temperature 
     399         CALL iom_put( "tosmint", rho0 * z2d )        ! Vertical integral of temperature 
    478400      ENDIF 
    479401      IF( iom_use("somint") ) THEN 
    480402         z2d(:,:)=0._wp 
    481          DO jk = 1, jpkm1 
    482             DO jj = 2, jpjm1 
    483                DO ji = fs_2, fs_jpim1   ! vector opt. 
    484                   z2d(ji,jj) = z2d(ji,jj) + e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) 
    485                END DO 
    486             END DO 
    487          END DO 
     403         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     404            z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) 
     405         END_3D 
    488406         CALL lbc_lnk( 'diawri', z2d, 'T', -1. ) 
    489          CALL iom_put( "somint", rau0 * z2d )         ! Vertical integral of salinity 
     407         CALL iom_put( "somint", rho0 * z2d )         ! Vertical integral of salinity 
    490408      ENDIF 
    491409 
    492410      CALL iom_put( "bn2", rn2 )                      ! Brunt-Vaisala buoyancy frequency (N^2) 
    493411      ! 
    494  
    495       IF (ln_diatmb)   CALL dia_tmb                   ! tmb values  
    496            
    497       IF (ln_dia25h)   CALL dia_25h( kt )             ! 25h averaging 
     412       
     413      IF (ln_dia25h)   CALL dia_25h( kt, Kmm )        ! 25h averaging 
    498414 
    499415      IF( ln_timing )   CALL timing_stop('dia_wri') 
     
    510426      INTEGER, DIMENSION(2) :: ierr 
    511427      !!---------------------------------------------------------------------- 
    512       ierr = 0 
    513       ALLOCATE( ndex_hT(jpi*jpj) , ndex_T(jpi*jpj*jpk) ,     & 
    514          &      ndex_hU(jpi*jpj) , ndex_U(jpi*jpj*jpk) ,     & 
    515          &      ndex_hV(jpi*jpj) , ndex_V(jpi*jpj*jpk) , STAT=ierr(1) ) 
     428      IF( nn_write == -1 ) THEN 
     429         dia_wri_alloc = 0 
     430      ELSE     
     431         ierr = 0 
     432         ALLOCATE( ndex_hT(jpi*jpj) , ndex_T(jpi*jpj*jpk) ,     & 
     433            &      ndex_hU(jpi*jpj) , ndex_U(jpi*jpj*jpk) ,     & 
     434            &      ndex_hV(jpi*jpj) , ndex_V(jpi*jpj*jpk) , STAT=ierr(1) ) 
    516435         ! 
    517       dia_wri_alloc = MAXVAL(ierr) 
    518       CALL mpp_sum( 'diawri', dia_wri_alloc ) 
     436         dia_wri_alloc = MAXVAL(ierr) 
     437         CALL mpp_sum( 'diawri', dia_wri_alloc ) 
     438         ! 
     439      ENDIF 
    519440      ! 
    520441   END FUNCTION dia_wri_alloc 
     442  
     443   INTEGER FUNCTION dia_wri_alloc_abl() 
     444      !!---------------------------------------------------------------------- 
     445     ALLOCATE(   ndex_hA(jpi*jpj), ndex_A (jpi*jpj*jpkam1), STAT=dia_wri_alloc_abl) 
     446      CALL mpp_sum( 'diawri', dia_wri_alloc_abl ) 
     447      ! 
     448   END FUNCTION dia_wri_alloc_abl 
    521449 
    522450    
    523    SUBROUTINE dia_wri( kt ) 
     451   SUBROUTINE dia_wri( kt, Kmm ) 
    524452      !!--------------------------------------------------------------------- 
    525453      !!                  ***  ROUTINE dia_wri  *** 
     
    531459      !!      define all the NETCDF files and fields 
    532460      !!      At each time step call histdef to compute the mean if ncessary 
    533       !!      Each nwrite time step, output the instantaneous or mean fields 
     461      !!      Each nn_write time step, output the instantaneous or mean fields 
    534462      !!---------------------------------------------------------------------- 
    535463      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index 
     464      INTEGER, INTENT( in ) ::   Kmm  ! ocean time level index 
    536465      ! 
    537466      LOGICAL ::   ll_print = .FALSE.                        ! =T print and flush numout 
     
    541470      INTEGER  ::   ierr                                     ! error code return from allocation 
    542471      INTEGER  ::   iimi, iima, ipk, it, itmod, ijmi, ijma   ! local integers 
     472      INTEGER  ::   ipka                                     ! ABL 
    543473      INTEGER  ::   jn, ierror                               ! local integers 
    544474      REAL(wp) ::   zsto, zout, zmax, zjulian                ! local scalars 
     
    546476      REAL(wp), DIMENSION(jpi,jpj)   :: zw2d       ! 2D workspace 
    547477      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zw3d       ! 3D workspace 
     478      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d_abl   ! ABL 3D workspace 
    548479      !!---------------------------------------------------------------------- 
     480      ! 
     481      IF( ninist == 1 ) THEN     !==  Output the initial state and forcings  ==! 
     482         CALL dia_wri_state( Kmm, 'output.init' ) 
     483         ninist = 0 
     484      ENDIF 
     485      ! 
     486      IF( nn_write == -1 )   RETURN   ! we will never do any output 
    549487      !  
    550488      IF( ln_timing )   CALL timing_start('dia_wri') 
    551       ! 
    552       IF( ninist == 1 ) THEN     !==  Output the initial state and forcings  ==! 
    553          CALL dia_wri_state( 'output.init' ) 
    554          ninist = 0 
    555       ENDIF 
    556489      ! 
    557490      ! 0. Initialisation 
     
    564497      clop = "x"         ! no use of the mask value (require less cpu time and otherwise the model crashes) 
    565498#if defined key_diainstant 
    566       zsto = nwrite * rdt 
     499      zsto = nn_write * rn_Dt 
    567500      clop = "inst("//TRIM(clop)//")" 
    568501#else 
    569       zsto=rdt 
     502      zsto=rn_Dt 
    570503      clop = "ave("//TRIM(clop)//")" 
    571504#endif 
    572       zout = nwrite * rdt 
    573       zmax = ( nitend - nit000 + 1 ) * rdt 
     505      zout = nn_write * rn_Dt 
     506      zmax = ( nitend - nit000 + 1 ) * rn_Dt 
    574507 
    575508      ! Define indices of the horizontal output zoom and vertical limit storage 
     
    577510      ijmi = 1      ;      ijma = jpj 
    578511      ipk = jpk 
     512      IF(ln_abl) ipka = jpkam1 
    579513 
    580514      ! define time axis 
     
    591525 
    592526         ! Compute julian date from starting date of the run 
    593          CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian ) 
     527         CALL ymds2ju( nyear, nmonth, nday, rn_Dt, zjulian ) 
    594528         zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment 
    595529         IF(lwp)WRITE(numout,*) 
     
    601535         ! WRITE root name in date.file for use by postpro 
    602536         IF(lwp) THEN 
    603             CALL dia_nam( clhstnam, nwrite,' ' ) 
     537            CALL dia_nam( clhstnam, nn_write,' ' ) 
    604538            CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
    605539            WRITE(inum,*) clhstnam 
     
    609543         ! Define the T grid FILE ( nid_T ) 
    610544 
    611          CALL dia_nam( clhstnam, nwrite, 'grid_T' ) 
     545         CALL dia_nam( clhstnam, nn_write, 'grid_T' ) 
    612546         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename 
    613547         CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit 
    614548            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       & 
    615             &          nit000-1, zjulian, rdt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set ) 
     549            &          nit000-1, zjulian, rn_Dt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set ) 
    616550         CALL histvert( nid_T, "deptht", "Vertical T levels",      &  ! Vertical grid: gdept 
    617551            &           "m", ipk, gdept_1d, nz_T, "down" ) 
     
    645579         ! Define the U grid FILE ( nid_U ) 
    646580 
    647          CALL dia_nam( clhstnam, nwrite, 'grid_U' ) 
     581         CALL dia_nam( clhstnam, nn_write, 'grid_U' ) 
    648582         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename 
    649583         CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu,           &  ! Horizontal grid: glamu and gphiu 
    650584            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       & 
    651             &          nit000-1, zjulian, rdt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set ) 
     585            &          nit000-1, zjulian, rn_Dt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set ) 
    652586         CALL histvert( nid_U, "depthu", "Vertical U levels",      &  ! Vertical grid: gdept 
    653587            &           "m", ipk, gdept_1d, nz_U, "down" ) 
     
    658592         ! Define the V grid FILE ( nid_V ) 
    659593 
    660          CALL dia_nam( clhstnam, nwrite, 'grid_V' )                   ! filename 
     594         CALL dia_nam( clhstnam, nn_write, 'grid_V' )                   ! filename 
    661595         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam 
    662596         CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv,           &  ! Horizontal grid: glamv and gphiv 
    663597            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       & 
    664             &          nit000-1, zjulian, rdt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set ) 
     598            &          nit000-1, zjulian, rn_Dt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set ) 
    665599         CALL histvert( nid_V, "depthv", "Vertical V levels",      &  ! Vertical grid : gdept 
    666600            &          "m", ipk, gdept_1d, nz_V, "down" ) 
     
    671605         ! Define the W grid FILE ( nid_W ) 
    672606 
    673          CALL dia_nam( clhstnam, nwrite, 'grid_W' )                   ! filename 
     607         CALL dia_nam( clhstnam, nn_write, 'grid_W' )                   ! filename 
    674608         IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam 
    675609         CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit 
    676610            &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       & 
    677             &          nit000-1, zjulian, rdt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set ) 
     611            &          nit000-1, zjulian, rn_Dt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set ) 
    678612         CALL histvert( nid_W, "depthw", "Vertical W levels",      &  ! Vertical grid: gdepw 
    679613            &          "m", ipk, gdepw_1d, nz_W, "down" ) 
    680614 
     615         IF( ln_abl ) THEN  
     616         ! Define the ABL grid FILE ( nid_A ) 
     617            CALL dia_nam( clhstnam, nn_write, 'grid_ABL' ) 
     618            IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam    ! filename 
     619            CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,           &  ! Horizontal grid: glamt and gphit 
     620               &          iimi, iima-iimi+1, ijmi, ijma-ijmi+1,       & 
     621               &          nit000-1, zjulian, rn_Dt, nh_A, nid_A, domain_id=nidom, snc4chunks=snc4set ) 
     622            CALL histvert( nid_A, "ght_abl", "Vertical T levels",      &  ! Vertical grid: gdept 
     623               &           "m", ipka, ght_abl(2:jpka), nz_A, "up" ) 
     624            !                                                            ! Index of ocean points 
     625         ALLOCATE( zw3d_abl(jpi,jpj,ipka) )  
     626         zw3d_abl(:,:,:) = 1._wp  
     627         CALL wheneq( jpi*jpj*ipka, zw3d_abl, 1, 1., ndex_A , ndim_A  )      ! volume 
     628            CALL wheneq( jpi*jpj     , zw3d_abl, 1, 1., ndex_hA, ndim_hA )      ! surface 
     629         DEALLOCATE(zw3d_abl) 
     630         ENDIF 
    681631 
    682632         ! Declare all the output fields as NETCDF variables 
     
    688638            &          jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    689639         IF(  .NOT.ln_linssh  ) THEN 
    690             CALL histdef( nid_T, "vovvle3t", "Level thickness"                    , "m"      ,&  ! e3t_n 
     640            CALL histdef( nid_T, "vovvle3t", "Level thickness"                    , "m"      ,&  ! e3t(:,:,:,Kmm) 
    691641            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    692             CALL histdef( nid_T, "vovvldep", "T point depth"                      , "m"      ,&  ! e3t_n 
     642            CALL histdef( nid_T, "vovvldep", "T point depth"                      , "m"      ,&  ! e3t(:,:,:,Kmm) 
    693643            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    694             CALL histdef( nid_T, "vovvldef", "Squared level deformation"          , "%^2"    ,&  ! e3t_n 
     644            CALL histdef( nid_T, "vovvldef", "Squared level deformation"          , "%^2"    ,&  ! e3t(:,:,:,Kmm) 
    695645            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    696646         ENDIF 
     
    709659            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    710660         IF(  ln_linssh  ) THEN 
    711             CALL histdef( nid_T, "sosst_cd", "Concentration/Dilution term on temperature"     &  ! emp * tsn(:,:,1,jp_tem) 
     661            CALL histdef( nid_T, "sosst_cd", "Concentration/Dilution term on temperature"     &  ! emp * ts(:,:,1,jp_tem,Kmm) 
    712662            &                                                                  , "KgC/m2/s",  &  ! sosst_cd 
    713663            &             jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    714             CALL histdef( nid_T, "sosss_cd", "Concentration/Dilution term on salinity"        &  ! emp * tsn(:,:,1,jp_sal) 
     664            CALL histdef( nid_T, "sosss_cd", "Concentration/Dilution term on salinity"        &  ! emp * ts(:,:,1,jp_sal,Kmm) 
    715665            &                                                                  , "KgPSU/m2/s",&  ! sosss_cd 
    716666            &             jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
     
    728678         CALL histdef( nid_T, "sowindsp", "wind speed at 10m"                  , "m/s"    ,   &  ! wndm 
    729679            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    730 ! 
     680         ! 
     681         IF( ln_abl ) THEN 
     682            CALL histdef( nid_A, "t_abl", "Potential Temperature"     , "K"        ,       &  ! t_abl 
     683               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 
     684            CALL histdef( nid_A, "q_abl", "Humidity"                  , "kg/kg"    ,       &  ! q_abl 
     685               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )  
     686            CALL histdef( nid_A, "u_abl", "Atmospheric U-wind   "     , "m/s"        ,     &  ! u_abl 
     687               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 
     688            CALL histdef( nid_A, "v_abl", "Atmospheric V-wind   "     , "m/s"    ,         &  ! v_abl 
     689               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )  
     690            CALL histdef( nid_A, "tke_abl", "Atmospheric TKE   "     , "m2/s2"    ,        &  ! tke_abl 
     691               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )  
     692            CALL histdef( nid_A, "avm_abl", "Atmospheric turbulent viscosity", "m2/s"   ,  &  ! avm_abl 
     693               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )  
     694            CALL histdef( nid_A, "avt_abl", "Atmospheric turbulent diffusivity", "m2/s2",  &  ! avt_abl 
     695               &          jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout )  
     696            CALL histdef( nid_A, "pblh", "Atmospheric boundary layer height "  , "m",      &  ! pblh 
     697               &          jpi, jpj, nh_A,  1  , 1, 1   , -99 , 32, clop, zsto, zout )                  
     698#if defined key_si3 
     699            CALL histdef( nid_A, "oce_frac", "Fraction of open ocean"  , " ",      &  ! ato_i 
     700               &          jpi, jpj, nh_A,  1  , 1, 1   , -99 , 32, clop, zsto, zout ) 
     701#endif 
     702            CALL histend( nid_A, snc4chunks=snc4set ) 
     703         ENDIF 
     704         ! 
    731705         IF( ln_icebergs ) THEN 
    732706            CALL histdef( nid_T, "calving"             , "calving mass input"                       , "kg/s"   , & 
     
    762736         ENDIF 
    763737 
    764          IF( .NOT. ln_cpl ) THEN 
     738         IF( ln_ssr ) THEN 
    765739            CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp 
    766740               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
     
    770744               &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    771745         ENDIF 
    772  
    773          IF( ln_cpl .AND. nn_ice <= 1 ) THEN 
    774             CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp 
    775                &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    776             CALL histdef( nid_T, "sowafldp", "Surface Water Flux: Damping"        , "Kg/m2/s",   &  ! erp 
    777                &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    778             CALL histdef( nid_T, "sosafldp", "Surface salt flux: Damping"         , "Kg/m2/s",   &  ! erp * sn 
    779                &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    780          ENDIF 
    781           
     746        
    782747         clmx ="l_max(only(x))"    ! max index on a period 
    783748!         CALL histdef( nid_T, "sobowlin", "Bowl Index"                         , "W-point",   &  ! bowl INDEX  
     
    797762 
    798763         !                                                                                      !!! nid_U : 3D 
    799          CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! un 
     764         CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! uu(:,:,:,Kmm) 
    800765            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 
    801766         IF( ln_wave .AND. ln_sdw) THEN 
     
    810775 
    811776         !                                                                                      !!! nid_V : 3D 
    812          CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vn 
     777         CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vv(:,:,:,Kmm) 
    813778            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 
    814779         IF( ln_wave .AND. ln_sdw) THEN 
     
    823788 
    824789         !                                                                                      !!! nid_W : 3D 
    825          CALL histdef( nid_W, "vovecrtz", "Vertical Velocity"                  , "m/s"    ,   &  ! wn 
     790         CALL histdef( nid_W, "vovecrtz", "Vertical Velocity"                  , "m/s"    ,   &  ! ww 
    826791            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
    827792         CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity"          , "m2/s"   ,   &  ! avt 
     
    855820      ! donne le nombre d'elements, et ndex la liste des indices a sortir 
    856821 
    857       IF( lwp .AND. MOD( itmod, nwrite ) == 0 ) THEN  
     822      IF( lwp .AND. MOD( itmod, nn_write ) == 0 ) THEN  
    858823         WRITE(numout,*) 'dia_wri : write model outputs in NetCDF files at ', kt, 'time-step' 
    859824         WRITE(numout,*) '~~~~~~ ' 
     
    861826 
    862827      IF( .NOT.ln_linssh ) THEN 
    863          CALL histwrite( nid_T, "votemper", it, tsn(:,:,:,jp_tem) * e3t_n(:,:,:) , ndim_T , ndex_T  )   ! heat content 
    864          CALL histwrite( nid_T, "vosaline", it, tsn(:,:,:,jp_sal) * e3t_n(:,:,:) , ndim_T , ndex_T  )   ! salt content 
    865          CALL histwrite( nid_T, "sosstsst", it, tsn(:,:,1,jp_tem) * e3t_n(:,:,1) , ndim_hT, ndex_hT )   ! sea surface heat content 
    866          CALL histwrite( nid_T, "sosaline", it, tsn(:,:,1,jp_sal) * e3t_n(:,:,1) , ndim_hT, ndex_hT )   ! sea surface salinity content 
     828         CALL histwrite( nid_T, "votemper", it, ts(:,:,:,jp_tem,Kmm) * e3t(:,:,:,Kmm) , ndim_T , ndex_T  )   ! heat content 
     829         CALL histwrite( nid_T, "vosaline", it, ts(:,:,:,jp_sal,Kmm) * e3t(:,:,:,Kmm) , ndim_T , ndex_T  )   ! salt content 
     830         CALL histwrite( nid_T, "sosstsst", it, ts(:,:,1,jp_tem,Kmm) * e3t(:,:,1,Kmm) , ndim_hT, ndex_hT )   ! sea surface heat content 
     831         CALL histwrite( nid_T, "sosaline", it, ts(:,:,1,jp_sal,Kmm) * e3t(:,:,1,Kmm) , ndim_hT, ndex_hT )   ! sea surface salinity content 
    867832      ELSE 
    868          CALL histwrite( nid_T, "votemper", it, tsn(:,:,:,jp_tem) , ndim_T , ndex_T  )   ! temperature 
    869          CALL histwrite( nid_T, "vosaline", it, tsn(:,:,:,jp_sal) , ndim_T , ndex_T  )   ! salinity 
    870          CALL histwrite( nid_T, "sosstsst", it, tsn(:,:,1,jp_tem) , ndim_hT, ndex_hT )   ! sea surface temperature 
    871          CALL histwrite( nid_T, "sosaline", it, tsn(:,:,1,jp_sal) , ndim_hT, ndex_hT )   ! sea surface salinity 
     833         CALL histwrite( nid_T, "votemper", it, ts(:,:,:,jp_tem,Kmm) , ndim_T , ndex_T  )   ! temperature 
     834         CALL histwrite( nid_T, "vosaline", it, ts(:,:,:,jp_sal,Kmm) , ndim_T , ndex_T  )   ! salinity 
     835         CALL histwrite( nid_T, "sosstsst", it, ts(:,:,1,jp_tem,Kmm) , ndim_hT, ndex_hT )   ! sea surface temperature 
     836         CALL histwrite( nid_T, "sosaline", it, ts(:,:,1,jp_sal,Kmm) , ndim_hT, ndex_hT )   ! sea surface salinity 
    872837      ENDIF 
    873838      IF( .NOT.ln_linssh ) THEN 
    874          zw3d(:,:,:) = ( ( e3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 
    875          CALL histwrite( nid_T, "vovvle3t", it, e3t_n (:,:,:) , ndim_T , ndex_T  )   ! level thickness 
    876          CALL histwrite( nid_T, "vovvldep", it, gdept_n(:,:,:) , ndim_T , ndex_T  )   ! t-point depth 
     839         zw3d(:,:,:) = ( ( e3t(:,:,:,Kmm) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 
     840         CALL histwrite( nid_T, "vovvle3t", it, e3t (:,:,:,Kmm) , ndim_T , ndex_T  )   ! level thickness 
     841         CALL histwrite( nid_T, "vovvldep", it, gdept(:,:,:,Kmm) , ndim_T , ndex_T  )   ! t-point depth 
    877842         CALL histwrite( nid_T, "vovvldef", it, zw3d             , ndim_T , ndex_T  )   ! level thickness deformation 
    878843      ENDIF 
    879       CALL histwrite( nid_T, "sossheig", it, sshn          , ndim_hT, ndex_hT )   ! sea surface height 
     844      CALL histwrite( nid_T, "sossheig", it, ssh(:,:,Kmm)          , ndim_hT, ndex_hT )   ! sea surface height 
    880845      CALL histwrite( nid_T, "sowaflup", it, ( emp-rnf )   , ndim_hT, ndex_hT )   ! upward water flux 
    881846      CALL histwrite( nid_T, "sorunoff", it, rnf           , ndim_hT, ndex_hT )   ! river runoffs 
     
    884849                                                                                  ! in linear free surface case) 
    885850      IF( ln_linssh ) THEN 
    886          zw2d(:,:) = emp (:,:) * tsn(:,:,1,jp_tem) 
     851         zw2d(:,:) = emp (:,:) * ts(:,:,1,jp_tem,Kmm) 
    887852         CALL histwrite( nid_T, "sosst_cd", it, zw2d, ndim_hT, ndex_hT )          ! c/d term on sst 
    888          zw2d(:,:) = emp (:,:) * tsn(:,:,1,jp_sal) 
     853         zw2d(:,:) = emp (:,:) * ts(:,:,1,jp_sal,Kmm) 
    889854         CALL histwrite( nid_T, "sosss_cd", it, zw2d, ndim_hT, ndex_hT )          ! c/d term on sss 
    890855      ENDIF 
     
    895860      CALL histwrite( nid_T, "soicecov", it, fr_i          , ndim_hT, ndex_hT )   ! ice fraction    
    896861      CALL histwrite( nid_T, "sowindsp", it, wndm          , ndim_hT, ndex_hT )   ! wind speed    
    897 ! 
     862      ! 
     863      IF( ln_abl ) THEN  
     864         ALLOCATE( zw3d_abl(jpi,jpj,jpka) ) 
     865         IF( ln_mskland )   THEN  
     866            DO jk=1,jpka 
     867               zw3d_abl(:,:,jk) = tmask(:,:,1) 
     868            END DO        
     869         ELSE 
     870            zw3d_abl(:,:,:) = 1._wp      
     871         ENDIF        
     872         CALL histwrite( nid_A,  "pblh"   , it, pblh(:,:)                  *zw3d_abl(:,:,1     ), ndim_hA, ndex_hA )   ! pblh  
     873         CALL histwrite( nid_A,  "u_abl"  , it, u_abl   (:,:,2:jpka,nt_n  )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! u_abl 
     874         CALL histwrite( nid_A,  "v_abl"  , it, v_abl   (:,:,2:jpka,nt_n  )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! v_abl 
     875         CALL histwrite( nid_A,  "t_abl"  , it, tq_abl  (:,:,2:jpka,nt_n,1)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! t_abl 
     876         CALL histwrite( nid_A,  "q_abl"  , it, tq_abl  (:,:,2:jpka,nt_n,2)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! q_abl        
     877         CALL histwrite( nid_A,  "tke_abl", it, tke_abl (:,:,2:jpka,nt_n  )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! tke_abl 
     878         CALL histwrite( nid_A,  "avm_abl", it, avm_abl (:,:,2:jpka       )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! avm_abl 
     879         CALL histwrite( nid_A,  "avt_abl", it, avt_abl (:,:,2:jpka       )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A  )   ! avt_abl  
     880#if defined key_si3 
     881         CALL histwrite( nid_A,  "oce_frac"   , it, ato_i(:,:)                                  , ndim_hA, ndex_hA )   ! ato_i 
     882#endif 
     883         DEALLOCATE(zw3d_abl) 
     884      ENDIF 
     885      ! 
    898886      IF( ln_icebergs ) THEN 
    899887         ! 
     
    919907      ENDIF 
    920908 
    921       IF( .NOT. ln_cpl ) THEN 
     909      IF( ln_ssr ) THEN 
    922910         CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping 
    923911         CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping 
    924          IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1) 
    925          CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping 
    926       ENDIF 
    927       IF( ln_cpl .AND. nn_ice <= 1 ) THEN 
    928          CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping 
    929          CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping 
    930          IF( ln_ssr ) zw2d(:,:) = erp(:,:) * tsn(:,:,1,jp_sal) * tmask(:,:,1) 
     912         zw2d(:,:) = erp(:,:) * ts(:,:,1,jp_sal,Kmm) * tmask(:,:,1) 
    931913         CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping 
    932914      ENDIF 
     
    941923#endif 
    942924 
    943       CALL histwrite( nid_U, "vozocrtx", it, un            , ndim_U , ndex_U )    ! i-current 
     925      CALL histwrite( nid_U, "vozocrtx", it, uu(:,:,:,Kmm)            , ndim_U , ndex_U )    ! i-current 
    944926      CALL histwrite( nid_U, "sozotaux", it, utau          , ndim_hU, ndex_hU )   ! i-wind stress 
    945927 
    946       CALL histwrite( nid_V, "vomecrty", it, vn            , ndim_V , ndex_V  )   ! j-current 
     928      CALL histwrite( nid_V, "vomecrty", it, vv(:,:,:,Kmm)            , ndim_V , ndex_V  )   ! j-current 
    947929      CALL histwrite( nid_V, "sometauy", it, vtau          , ndim_hV, ndex_hV )   ! j-wind stress 
    948930 
    949       CALL histwrite( nid_W, "vovecrtz", it, wn             , ndim_T, ndex_T )    ! vert. current 
     931      IF( ln_zad_Aimp ) THEN 
     932         CALL histwrite( nid_W, "vovecrtz", it, ww + wi     , ndim_T, ndex_T )    ! vert. current 
     933      ELSE 
     934         CALL histwrite( nid_W, "vovecrtz", it, ww          , ndim_T, ndex_T )    ! vert. current 
     935      ENDIF 
    950936      CALL histwrite( nid_W, "votkeavt", it, avt            , ndim_T, ndex_T )    ! T vert. eddy diff. coef. 
    951937      CALL histwrite( nid_W, "votkeavm", it, avm            , ndim_T, ndex_T )    ! T vert. eddy visc. coef. 
     
    967953         CALL histclo( nid_V ) 
    968954         CALL histclo( nid_W ) 
     955         IF(ln_abl) CALL histclo( nid_A ) 
    969956      ENDIF 
    970957      ! 
     
    974961#endif 
    975962 
    976    SUBROUTINE dia_wri_state( cdfile_name ) 
     963   SUBROUTINE dia_wri_state( Kmm, cdfile_name ) 
    977964      !!--------------------------------------------------------------------- 
    978965      !!                 ***  ROUTINE dia_wri_state  *** 
     
    987974      !!      File 'output.abort.nc' is created in case of abnormal job end 
    988975      !!---------------------------------------------------------------------- 
     976      INTEGER           , INTENT( in ) ::   Kmm              ! time level index 
    989977      CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created 
    990978      !! 
    991       INTEGER :: inum 
     979      INTEGER :: inum, jk 
    992980      !!---------------------------------------------------------------------- 
    993981      !  
     
    996984      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created ' 
    997985      IF(lwp) WRITE(numout,*) '                and named :', cdfile_name, '...nc' 
    998  
    999 #if defined key_si3 
    1000      CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE., kdlev = jpl ) 
    1001 #else 
    1002      CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE. ) 
    1003 #endif 
    1004  
    1005       CALL iom_rstput( 0, 0, inum, 'votemper', tsn(:,:,:,jp_tem) )    ! now temperature 
    1006       CALL iom_rstput( 0, 0, inum, 'vosaline', tsn(:,:,:,jp_sal) )    ! now salinity 
    1007       CALL iom_rstput( 0, 0, inum, 'sossheig', sshn              )    ! sea surface height 
    1008       CALL iom_rstput( 0, 0, inum, 'vozocrtx', un                )    ! now i-velocity 
    1009       CALL iom_rstput( 0, 0, inum, 'vomecrty', vn                )    ! now j-velocity 
    1010       CALL iom_rstput( 0, 0, inum, 'vovecrtz', wn                )    ! now k-velocity 
     986      ! 
     987      CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE. ) 
     988      ! 
     989      CALL iom_rstput( 0, 0, inum, 'votemper', ts(:,:,:,jp_tem,Kmm) )    ! now temperature 
     990      CALL iom_rstput( 0, 0, inum, 'vosaline', ts(:,:,:,jp_sal,Kmm) )    ! now salinity 
     991      CALL iom_rstput( 0, 0, inum, 'sossheig', ssh(:,:,Kmm)              )    ! sea surface height 
     992      CALL iom_rstput( 0, 0, inum, 'vozocrtx', uu(:,:,:,Kmm)                )    ! now i-velocity 
     993      CALL iom_rstput( 0, 0, inum, 'vomecrty', vv(:,:,:,Kmm)                )    ! now j-velocity 
     994      IF( ln_zad_Aimp ) THEN 
     995         CALL iom_rstput( 0, 0, inum, 'vovecrtz', ww + wi        )    ! now k-velocity 
     996      ELSE 
     997         CALL iom_rstput( 0, 0, inum, 'vovecrtz', ww             )    ! now k-velocity 
     998      ENDIF 
     999      CALL iom_rstput( 0, 0, inum, 'risfdep', risfdep            )    ! now k-velocity 
     1000      CALL iom_rstput( 0, 0, inum, 'ht'     , ht                 )    ! now water column height 
     1001      ! 
     1002      IF ( ln_isf ) THEN 
     1003         IF (ln_isfcav_mlt) THEN 
     1004            CALL iom_rstput( 0, 0, inum, 'fwfisf_cav', fwfisf_cav          )    ! now k-velocity 
     1005            CALL iom_rstput( 0, 0, inum, 'rhisf_cav_tbl', rhisf_tbl_cav    )    ! now k-velocity 
     1006            CALL iom_rstput( 0, 0, inum, 'rfrac_cav_tbl', rfrac_tbl_cav    )    ! now k-velocity 
     1007            CALL iom_rstput( 0, 0, inum, 'misfkb_cav', REAL(misfkb_cav,wp) )    ! now k-velocity 
     1008            CALL iom_rstput( 0, 0, inum, 'misfkt_cav', REAL(misfkt_cav,wp) )    ! now k-velocity 
     1009            CALL iom_rstput( 0, 0, inum, 'mskisf_cav', REAL(mskisf_cav,wp), ktype = jp_i1 ) 
     1010         END IF 
     1011         IF (ln_isfpar_mlt) THEN 
     1012            CALL iom_rstput( 0, 0, inum, 'isfmsk_par', REAL(mskisf_par,wp) )    ! now k-velocity 
     1013            CALL iom_rstput( 0, 0, inum, 'fwfisf_par', fwfisf_par          )    ! now k-velocity 
     1014            CALL iom_rstput( 0, 0, inum, 'rhisf_par_tbl', rhisf_tbl_par    )    ! now k-velocity 
     1015            CALL iom_rstput( 0, 0, inum, 'rfrac_par_tbl', rfrac_tbl_par    )    ! now k-velocity 
     1016            CALL iom_rstput( 0, 0, inum, 'misfkb_par', REAL(misfkb_par,wp) )    ! now k-velocity 
     1017            CALL iom_rstput( 0, 0, inum, 'misfkt_par', REAL(misfkt_par,wp) )    ! now k-velocity 
     1018            CALL iom_rstput( 0, 0, inum, 'mskisf_par', REAL(mskisf_par,wp), ktype = jp_i1 ) 
     1019         END IF 
     1020      END IF 
     1021      ! 
    10111022      IF( ALLOCATED(ahtu) ) THEN 
    10121023         CALL iom_rstput( 0, 0, inum,  'ahtu', ahtu              )    ! aht at u-point 
     
    10241035      CALL iom_rstput( 0, 0, inum, 'sometauy', vtau              )    ! j-wind stress 
    10251036      IF(  .NOT.ln_linssh  ) THEN              
    1026          CALL iom_rstput( 0, 0, inum, 'vovvldep', gdept_n        )    !  T-cell depth  
    1027          CALL iom_rstput( 0, 0, inum, 'vovvle3t', e3t_n          )    !  T-cell thickness   
     1037         CALL iom_rstput( 0, 0, inum, 'vovvldep', gdept(:,:,:,Kmm)        )    !  T-cell depth  
     1038         CALL iom_rstput( 0, 0, inum, 'vovvle3t', e3t(:,:,:,Kmm)          )    !  T-cell thickness   
    10281039      END IF 
    10291040      IF( ln_wave .AND. ln_sdw ) THEN 
     
    10321043         CALL iom_rstput( 0, 0, inum, 'sdvecrtz', wsd            )    ! now StokesDrift k-velocity 
    10331044      ENDIF 
    1034   
     1045      IF ( ln_abl ) THEN 
     1046         CALL iom_rstput ( 0, 0, inum, "uz1_abl",   u_abl(:,:,2,nt_a  ) )   ! now first level i-wind 
     1047         CALL iom_rstput ( 0, 0, inum, "vz1_abl",   v_abl(:,:,2,nt_a  ) )   ! now first level j-wind 
     1048         CALL iom_rstput ( 0, 0, inum, "tz1_abl",  tq_abl(:,:,2,nt_a,1) )   ! now first level temperature 
     1049         CALL iom_rstput ( 0, 0, inum, "qz1_abl",  tq_abl(:,:,2,nt_a,2) )   ! now first level humidity 
     1050      ENDIF 
     1051      ! 
     1052      CALL iom_close( inum ) 
     1053      !  
    10351054#if defined key_si3 
    10361055      IF( nn_ice == 2 ) THEN   ! condition needed in case agrif + ice-model but no-ice in child grid 
     1056         CALL iom_open( TRIM(cdfile_name)//'_ice', inum, ldwrt = .TRUE., kdlev = jpl, cdcomp = 'ICE' ) 
    10371057         CALL ice_wri_state( inum ) 
     1058         CALL iom_close( inum ) 
    10381059      ENDIF 
    10391060#endif 
    1040       ! 
    1041       CALL iom_close( inum ) 
    1042       !  
     1061 
    10431062   END SUBROUTINE dia_wri_state 
    10441063 
Note: See TracChangeset for help on using the changeset viewer.