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 8150 for branches/2016/v3_6_CMIP6_ice_diagnostics/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90 – NEMO

Ignore:
Timestamp:
2017-06-07T16:37:36+02:00 (7 years ago)
Author:
vancop
Message:

SIMIP outputs, phase 2

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2016/v3_6_CMIP6_ice_diagnostics/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90

    r7536 r8150  
    8080      INTEGER  :: Nu, Nv                                                       ! passage size 
    8181 
    82       INTEGER, PARAMETER :: i_grid = 2                                         ! grid type (eORCA1 = 1, ORCA2 = 2) 
     82      INTEGER, DIMENSION(4)  :: ji0, ji1, jj0, jj1 
    8383       
    8484      !!------------------------------------------------------------------- 
     
    142142      ! 
    143143      IF ( iom_use( "miceage" ) )       CALL iom_put( "miceage"     , om_i * zswi * z1_365   )  ! mean ice age 
    144       IF ( iom_use( "icethic_cea" ) )   CALL iom_put( "icethic_cea" , htm_i * zswi           )  ! ice thickness mean 
    145       IF ( iom_use( "snowthic_cea" ) )  CALL iom_put( "snowthic_cea", htm_s * zswi           )  ! snow thickness mean 
    146144      IF ( iom_use( "micet" ) )         CALL iom_put( "micet"       , ( tm_i  - rt0 ) * zswi )  ! ice mean    temperature 
    147145      IF ( iom_use( "icest" ) )         CALL iom_put( "icest"       , ( tm_su - rt0 ) * zswi )  ! ice surface temperature 
     
    228226      CALL iom_put ('hfxspr'     , hfx_spr(:,:)         )   ! Heat content of snow precip  
    229227 
    230       !-------------------------------- 
    231       ! Add-ons for SIMIP 
    232       !-------------------------------- 
    233       zrho1 = ( rau0 - rhoic ) / rau0; zrho2 = rhosn / rau0 
    234  
    235       IF  ( iom_use( "icethic"  ) ) CALL iom_put( "icethic"     , htm_i * zswi               )          ! ice thickness  
    236       IF  ( iom_use( "icepres"  ) ) CALL iom_put( "icepres"     , zswi                       )          ! ice presence (1 or 0)  
    237       IF  ( iom_use( "snowthic" ) ) CALL iom_put( "snowthic"    , htm_s * zswi               )          ! snow thickness        
    238       IF  ( iom_use( "icemass"  ) ) CALL iom_put( "icemass"     , rhoic * vt_i(:,:) * zswi   )          ! ice mass per cell area  
    239       IF  ( iom_use( "snomass"  ) ) CALL iom_put( "snomass"     , rhosn * vt_s(:,:) * zswi   )          ! snow mass per cell area 
    240       IF  ( iom_use( "icesnt"   ) ) CALL iom_put( "icesnt"      , ( tm_si - rt0 ) * zswi     )          ! snow-ice interface temperature 
    241       IF  ( iom_use( "icebot"   ) ) CALL iom_put( "icebot"      , ( t_bo  - rt0 ) * zswi     )          ! ice bottom temperature 
    242       IF  ( iom_use( "icesmass" ) ) CALL iom_put( "icesmass"    , smt_i * vt_i * rhoic / 1000. * zswi ) ! mass of salt in sea ice per cell area 
    243       IF  ( iom_use( "icefb"    ) ) CALL iom_put( "icefb"       , ( zrho1 * htm_i(:,:) - zrho2 * htm_s(:,:) ) * zswi ) ! mass of salt in sea ice per cell area 
    244  
    245       IF  ( iom_use( "wfxsum"   ) ) CALL iom_put( "wfxsum"      ,   wfx_sum                  )          ! Freshwater flux from sea-ice surface 
    246       IF  ( iom_use( "dmithd"   ) ) CALL iom_put( "dmithd"      , - wfx_bog - wfx_bom - wfx_sum   &     ! Sea-ice mass change from thermodynamics 
    247               &                     - wfx_sni - wfx_opw - wfx_res ) 
    248       IF  ( iom_use( "dmidyn"   ) ) CALL iom_put( "dmidyn"      ,   diag_dmi_dyn             )          ! Sea-ice mass change from dynamics 
    249       IF  ( iom_use( "dmiopw"   ) ) CALL iom_put( "dmiopw"      , - wfx_opw                  )          ! Sea-ice mass change through growth in open water 
    250       IF  ( iom_use( "dmibog"   ) ) CALL iom_put( "dmibog"      , - wfx_bog                  )          ! Sea-ice mass change through basal growth 
    251       IF  ( iom_use( "dmisni"   ) ) CALL iom_put( "dmisni"      , - wfx_sni                  )          ! Sea-ice mass change through snow-to-ice conversion 
    252       IF  ( iom_use( "dmisum"   ) ) CALL iom_put( "dmisum"      , - wfx_sum                  )          ! Sea-ice mass change through surface melting 
    253       IF  ( iom_use( "dmibom"   ) ) CALL iom_put( "dmibom"      , - wfx_bom                  )          ! Sea-ice mass change through bottom melting 
    254       IF  ( iom_use( "dmtsub"   ) ) CALL iom_put( "dmtsub"      , - wfx_sub                  )          ! Sea-ice mass change through evaporation and sublimation 
    255       IF  ( iom_use( "dmsspr"   ) ) CALL iom_put( "dmsspr"      , - wfx_spr                  )          ! snow mass change through snow fall 
    256       IF  ( iom_use( "dmsssi"   ) ) CALL iom_put( "dmsssi"      ,   wfx_sni*rhosn/rhoic      )          ! snow mass change through snow-to-ice conversion 
    257  
    258       IF  ( iom_use( "dmsmel"   ) ) CALL iom_put( "dmsmel"      ,   diag_dms_mel             )          ! snow mass change through melt 
    259       IF  ( iom_use( "dmsdyn"   ) ) CALL iom_put( "dmsdyn"      ,   diag_dms_dyn             )          ! snow mass change through dynamics 
    260  
    261       IF  ( iom_use( "hfxconbo" ) ) CALL iom_put( "hfxconbo"    ,   diag_fc_bo               )          ! bottom conduction flux 
    262       IF  ( iom_use( "hfxconsu" ) ) CALL iom_put( "hfxconsu"    ,   diag_fc_su               )          ! surface conduction flux 
    263  
    264       IF  ( iom_use( "wfxtot"   ) ) CALL iom_put( "wfxtot"      ,   -wfx_ice                 )          ! total freshwater flux from sea ice 
    265  
    266       IF  ( iom_use( "dmtxdyn"  ) ) CALL iom_put( "dmtxdyn"     ,   diag_dmtx_dyn            )          ! X-component of sea-ice mass transport 
    267       IF  ( iom_use( "dmtydyn"  ) ) CALL iom_put( "dmtydyn"     ,   diag_dmty_dyn            )          ! Y-component of sea-ice mass transport 
    268  
    269       IF  ( iom_use( "utau_oce" ) ) CALL iom_put( "utau_oce"    ,   diag_utau_oi*zswi        )          ! X-component of ocean stress on sea ice 
    270       IF  ( iom_use( "vtau_oce" ) ) CALL iom_put( "vtau_oce"    ,   diag_vtau_oi*zswi        )          ! Y-component of ocean stress on sea ice 
    271  
    272       IF  ( iom_use( "dssh_dx"  ) ) CALL iom_put( "dssh_dx"     ,   diag_dssh_dx*zswi        )          ! Sea-surface tilt term in force balance (x-component) 
    273       IF  ( iom_use( "dssh_dy"  ) ) CALL iom_put( "dssh_dy"     ,   diag_dssh_dy*zswi        )          ! Sea-surface tilt term in force balance (y-component) 
    274  
    275       IF  ( iom_use( "corstrx"  ) ) CALL iom_put( "corstrx"     ,   diag_corstrx*zswi        )          ! Coriolis force term in force balance (x-component) 
    276       IF  ( iom_use( "corstry"  ) ) CALL iom_put( "corstry"     ,   diag_corstry*zswi        )          ! Coriolis force term in force balance (y-component) 
    277  
    278       IF  ( iom_use( "intstrx"  ) ) CALL iom_put( "intstrx"     ,   diag_intstrx*zswi        )          ! Internal force term in force balance (x-component) 
    279       IF  ( iom_use( "intstry"  ) ) CALL iom_put( "intstry"     ,   diag_intstry*zswi        )          ! Internal force term in force balance (y-component) 
    280  
    281       IF  ( iom_use( "normstr"  ) ) CALL iom_put( "normstr"     ,   diag_sig1   *zswi        )          ! Normal stress 
    282       IF  ( iom_use( "sheastr"  ) ) CALL iom_put( "sheastr"     ,   diag_sig2   *zswi        )          ! Shear stress 
    283  
    284       !-------------------------------- 
    285       ! Global ice diagnostics (SIMIP) 
    286       !-------------------------------- 
    287       ! 
    288       ! What follows could be a separate routine 
    289       ! 
    290  
    291        IF ( iom_use( "NH_icearea" ) .OR. iom_use( "NH_icevolu" ) ) THEN    ! NH ice area  
    292   
    293          WHERE( fcor > 0 ); zswi(:,:) = 1.0e-12;  
    294          ELSEWHERE        ; zswi(:,:) = 0. 
    295          END WHERE  
    296  
    297            IF ( iom_use( "NH_icearea" ) ) THEN 
    298               zdiag_area_nh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) ) 
    299               CALL iom_put( "NH_icearea", zdiag_area_nh ) 
    300            ENDIF 
    301            IF ( iom_use( "NH_icevolu" ) ) THEN 
    302               zdiag_volu_nh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) ) 
    303               CALL iom_put( "NH_icevolu", zdiag_volu_nh ) 
    304            ENDIF 
    305  
    306       ENDIF  
    307  
    308       IF ( iom_use( "NH_iceextt" ) ) THEN                                   ! NH ice extt  
    309  
    310          WHERE( fcor > 0 .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12;  
    311          ELSEWHERE                          ; zswi(:,:) = 0. 
    312          END WHERE  
    313  
    314          zdiag_extt_nh = glob_sum( zswi(:,:) * e12t(:,:) ) 
    315          CALL iom_put( "NH_iceextt", zdiag_extt_nh ) 
    316  
    317       ENDIF  
    318  
    319       IF ( iom_use( "SH_icearea" ) .OR. iom_use( "SH_icevolu" ) ) THEN    ! SH ice area  / volume 
    320  
    321          WHERE( fcor < 0 ); zswi(:,:) = 1.0e-12;  
    322          ELSEWHERE        ; zswi(:,:) = 0. 
    323          END WHERE  
    324  
    325          IF ( iom_use( "SH_icearea" ) ) THEN 
    326             zdiag_area_sh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) )  
    327             CALL iom_put( "SH_icearea", zdiag_area_sh ) 
    328          ENDIF 
    329          IF ( iom_use( "SH_icevolu" ) ) THEN 
    330             zdiag_volu_sh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) ) 
    331             CALL iom_put( "SH_icevolu", zdiag_volu_sh ) 
    332          ENDIF 
    333  
    334       ENDIF  
    335  
    336       IF ( iom_use( "SH_iceextt" ) ) THEN                                   ! SH ice extt  
    337  
    338          WHERE( fcor < 0 .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12;  
    339          ELSEWHERE                          ; zswi(:,:) = 0. 
    340          END WHERE  
    341  
    342          zdiag_extt_sh = glob_sum( zswi(:,:) * e12t(:,:) ) 
    343  
    344          CALL iom_put( "SH_iceextt", zdiag_extt_sh ) 
    345  
    346       ENDIF  
    347  
    348       !------------------------ 
    349       ! Fluxes through straits 
    350       !------------------------ 
    351       ! 
    352       ! This piece of code is quite awful and should probably be a separate routine 
    353       ! 
    354       ! See Notz et al 2016 for definitions 
    355       ! 4 Arctic passages are considered (Fram, CAA, Barents, Bering) 
    356       ! 
    357  
    358       IF ( iom_use("strait_arfl") .OR. iom_use("strait_mifl") .OR. iom_use("strait_msfl") ) THEN 
    359  
    360          zdiag_area_strait(:) = 0._wp 
    361          zdiag_mice_strait(:) = 0._wp 
    362          zdiag_msno_strait(:) = 0._wp 
    363     
    364          ! === Fram Strait === 
    365          IF ( i_grid == 2 ) THEN     ! ORCA2 
    366             Nv  = 4 
    367             zvi(1:Nv) = (/ 133, 134, 135, 136 /) 
    368             zvj(1:Nv) = (/ 136, 136, 136, 136 /) 
    369          ENDIF 
    370     
    371          IF ( i_grid == 1 ) THEN    ! eORCA1 
    372             Nv  = 10 
    373             zvi(1:Nv) = (/ 268,269,270,271,272,273,274,275,276,277 /) 
    374             zvj(1:Nv) = (/ 311,311,311,311,311,311,311,311,311,311 /) 
    375          ENDIF 
    376     
    377          DO ii = 1, Nv 
    378             ji = zvi(ii)     
    379             jj = zvj(ii) 
    380             zdiag_area_strait(1) = zdiag_area_strait(1)                               & ! --- ice area flux --- 
    381                 &                + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         & ! northwards (positive) flow 
    382                 &                + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 )           ! southwards (negative) flow 
    383     
    384             zdiag_mice_strait(1) = zdiag_mice_strait(1) + rhoic *                     & ! --- ice mass flux  --- 
    385                 &                ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         & ! 
    386                 &                + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) )         ! 
    387     
    388             zdiag_msno_strait(1) = zdiag_msno_strait(1) + rhosn *                     & ! --- snow mass flux --- 
    389                 &                ( vt_s(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         &  
    390                 &                + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) )   
    391          END DO 
    392     
    393          ! === Bering Strait === 
    394          IF ( i_grid == 1 ) THEN ! eORCA1 
    395             Nv  = 3 
    396             zvi(1:Nv) = (/ 113,114,115 /) 
    397             zvj(1:Nv) = (/ 285,285,285 /) 
    398          ENDIF 
    399     
    400          IF ( i_grid == 2 ) THEN ! ORCA2 
    401             Nv  = 2 
    402             zvi(1:Nv) = (/ 55 , 56 /) 
    403             zvj(1:Nv) = (/ 122,122 /) 
    404          ENDIF 
    405     
    406          DO ii = 1, Nv 
    407             ji = zvi(ii)      
    408             jj = zvj(ii) 
    409             zdiag_area_strait(4) = zdiag_area_strait(4)                               & ! --- ice area flux --- 
    410                 &                + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         & ! northwards (positive) flow 
    411                 &                + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 )           ! southwards (negative) flow 
    412     
    413             zdiag_mice_strait(4) = zdiag_mice_strait(4) + rhoic *                     & ! --- ice mass flux  --- 
    414                 &                ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )            & ! 
    415                 &                + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) )            ! 
    416     
    417             zdiag_msno_strait(4) = zdiag_msno_strait(4) + rhosn *                     & ! --- snow mass flux --- 
    418                 &                ( vt_s(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )            &  
    419                 &                + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) )   
    420          END DO 
    421     
    422          ! === Barents throughflow (eORCA1) 
    423  
    424          ! U-flow 
    425          IF ( i_grid == 1 ) THEN ! 'eORCA1' 
    426             Nu = 11 
    427             zui(1:Nu) = (/ 282,283,284,285,286,286,287,288,289,290,292/) 
    428             zuj(1:Nu) = (/ 308,307,306,305,304,303,302,301,300,299,298/) 
    429          ENDIF 
    430          IF ( i_grid == 2 ) THEN ! 'ORCA2' 
    431             Nu = 5 
    432             zui(1:Nu) = (/ 141,142,142,143,144 /) 
    433             zuj(1:Nu) = (/ 134,133,132,131,130 /) 
    434          ENDIF 
    435          zfarea_u = 0._wp 
    436          zfmice_u = 0._wp 
    437          zfmsno_u = 0._wp 
    438     
    439          DO ii = 1, Nu 
    440             ji = zui(ii)      
    441             jj = zuj(ii) 
    442             zfarea_u          = zfarea_u                                           & ! --- ice area zonal flux --- 
    443                 &             + at_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         & ! --- eastward 
    444                 &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 )           ! --- westward 
    445             zfmice_u          = zfmice_u + rhoic *                                 & ! --- ice mass zonal flux ---  
    446                 &             ( vt_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
    447                 &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
    448             zfmsno_u          = zfmsno_u + rhosn *                                 & ! --- snow mass zonal flux ---  
    449                 &             ( vt_s(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
    450                 &             + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
    451          END DO 
    452     
    453          ! V-flow 
    454          IF ( i_grid == 1 ) THEN ! 'eORCA1' 
    455             Nv = 9 
    456             zvi(1:Nv) = (/ 282,283,284,285,286,287,288,289,290/) 
    457             zvj(1:Nv) = (/ 308,307,306,305,303,302,301,300,299/) 
    458          ENDIF 
    459          IF ( i_grid == 2 ) THEN ! 'ORCA2' 
    460             Nv = 4 
    461             zvi(1:Nv) = (/ 140,141,142,143 /) 
    462             zvj(1:Nv) = (/ 135,134,132,131 /) 
    463          ENDIF 
    464          zfarea_v = 0._wp 
    465          zfmice_v = 0._wp 
    466          zfmsno_v = 0._wp 
    467     
    468          DO ii  = 1, Nv   
    469             ji = zvi(ii)      
    470             jj = zvj(ii) 
    471             zfarea_v          = zfarea_v                                           & ! --- ice area meridian flux --- 
    472                 &             + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! --- eastward 
    473                 &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 )           ! --- westward 
    474             zfmice_v          = zfmice_v + rhoic *                                 & ! --- ice mass meridian flux --- 
    475                 &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
    476                 &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
    477             zfmsno_v          = zfmsno_v + rhosn *                                 & ! --- snow mass meridian flux --- 
    478                 &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
    479                 &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
    480          END DO 
    481     
    482          ! Sum U/V contributions 
    483          zdiag_area_strait(3) = zfarea_u + zfarea_v  
    484          zdiag_mice_strait(3) = zfmice_u + zfmice_v 
    485          zdiag_msno_strait(3) = zfmsno_u + zfmsno_v 
    486     
    487          ! === CAA throughflow === 
    488          ! U-flow through Queen Elisabeth Islands and McClure straits 
    489          IF ( i_grid == 1 ) THEN  ! eORCA1 
    490             Nu = 8  
    491             zui(1:Nu) = (/ 231,231,231,  132,132,132,132,132  /) 
    492             zuj(1:Nu) = (/ 328,329,330,  318,319,320,321,322  /) 
    493             zfarea_u = 0._wp 
    494             zfmice_u = 0._wp 
    495             zfmsno_u = 0._wp 
    496        
    497             DO ii = 1, Nu 
    498                ji = zui(ii)      
    499                jj = zuj(ii) 
    500                zfarea_u          = zfarea_u                                           & ! --- ice area zonal flux --- 
    501                    &             + at_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         & ! --- eastward 
    502                    &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 )           ! --- westward 
    503                zfmice_u          = zfmice_u + rhoic *                                 & ! --- ice mass zonal flux ---  
    504                    &             ( vt_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
    505                    &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
    506                zfmsno_u          = zfmsno_u + rhosn *                                 & ! --- snow mass zonal flux ---  
    507                    &             ( vt_s(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
    508                    &             + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
    509             END DO 
    510          ENDIF 
    511  
    512          IF ( i_grid == 2 ) THEN   ! ORCA2 
    513             zfarea_u = 0._wp       ! QEI and McClure straits are not resolved in ORCA2 
    514             zfmice_u = 0._wp 
    515             zfmsno_u = 0._wp 
    516          ENDIF 
    517     
    518          ! V-flow through Nares Strait 
    519          IF ( i_grid == 1 ) THEN   ! eORCA1 
    520             Nv = 4 
    521             zvi(1:Nv) = (/ 254,255,256,257 /) 
    522             zvj(1:Nv) = (/ 317,317,317,317 /) 
    523          ENDIF 
    524          IF ( i_grid == 2 ) THEN   ! ORCA2 
    525             Nv = 2 
    526             zvi(1:Nv) = (/ 117,118 /) 
    527             zvj(1:Nv) = (/ 145,145 /) 
    528          ENDIF 
    529          zfarea_v = 0._wp 
    530          zfmice_v = 0._wp 
    531          zfmsno_v = 0._wp 
    532     
    533          DO ii = 1, Nv   
    534             ji = zvi(ii)      
    535             jj = zvj(ii) 
    536             zfarea_v          = zfarea_v                                           & ! --- ice area meridian flux --- 
    537                 &             + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! --- eastward 
    538                 &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 )           ! --- westward 
    539             zfmice_v          = zfmice_v + rhoic *                                 & ! --- ice mass meridian flux --- 
    540                 &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
    541                 &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
    542             zfmsno_v          = zfmsno_v + rhosn *                                 & ! --- snow mass meridian flux --- 
    543                 &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
    544                 &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
    545          END DO 
    546  
    547          ! Sum U/V contributions 
    548          zdiag_area_strait(2) = zfarea_u + zfarea_v  
    549          zdiag_mice_strait(2) = zfmice_u + zfmice_v 
    550          zdiag_msno_strait(2) = zfmsno_u + zfmsno_v 
    551     
    552          ! === Write in file 
    553          IF ( iom_use("strait_arfl") ) CALL iom_put( "strait_arfl", zdiag_area_strait ) 
    554          IF ( iom_use("strait_mifl") ) CALL iom_put( "strait_mifl", zdiag_mice_strait ) 
    555          IF ( iom_use("strait_msfl") ) CALL iom_put( "strait_msfl", zdiag_msno_strait )  
    556  
    557          WRITE(numout,*) " area flx ", zdiag_area_strait(:) 
    558          WRITE(numout,*) " mice flx ", zdiag_mice_strait(:) 
    559          WRITE(numout,*) " msno flx ", zdiag_msno_strait(:) 
    560  
    561       ENDIF 
    562  
    563228      !---------------------------------- 
    564229      ! Output category-dependent fields 
     
    576241      ! brine volume 
    577242      IF ( iom_use( "brinevol_cat" ) )  CALL iom_put( "brinevol_cat", bv_i * 100. * zswi2 ) 
     243 
     244      !-------------------------------- 
     245      ! Add-ons for SIMIP 
     246      !-------------------------------- 
     247      zrho1 = ( rau0 - rhoic ) / rau0; zrho2 = rhosn / rau0 
     248 
     249      IF  ( iom_use( "icethic"  ) ) CALL iom_put( "icethic"     , htm_i * zswi               )          ! ice thickness  
     250      IF  ( iom_use( "icepres"  ) ) CALL iom_put( "icepres"     , zswi                       )          ! ice presence (1 or 0)  
     251      IF  ( iom_use( "snowthic" ) ) CALL iom_put( "snowthic"    , htm_s * zswi               )          ! snow thickness        
     252      IF  ( iom_use( "icemass"  ) ) CALL iom_put( "icemass"     , rhoic * vt_i(:,:) * zswi   )          ! ice mass per cell area  
     253      IF  ( iom_use( "snomass"  ) ) CALL iom_put( "snomass"     , rhosn * vt_s(:,:) * zswi   )          ! snow mass per cell area 
     254      IF  ( iom_use( "icesnt"   ) ) CALL iom_put( "icesnt"      , ( tm_si - rt0 ) * zswi     )          ! snow-ice interface temperature 
     255      IF  ( iom_use( "icebot"   ) ) CALL iom_put( "icebot"      , ( t_bo  - rt0 ) * zswi     )          ! ice bottom temperature 
     256      IF  ( iom_use( "icesmass" ) ) CALL iom_put( "icesmass"    , SUM( smv_i, DIM = 3 ) * rhoic / 1000. * zswi )         ! mass of salt in sea ice per cell area 
     257      IF  ( iom_use( "icefb"    ) ) CALL iom_put( "icefb"       , ( zrho1 * htm_i(:,:) - zrho2 * htm_s(:,:) ) * zswi )   ! ice freeboard 
     258 
     259      IF  ( iom_use( "wfxsum"   ) ) CALL iom_put( "wfxsum"      ,   wfx_sum                  )          ! Freshwater flux from sea-ice surface 
     260      IF  ( iom_use( "dmithd"   ) ) CALL iom_put( "dmithd"      , - wfx_bog - wfx_bom - wfx_sum   &     ! Sea-ice mass change from thermodynamics 
     261              &                     - wfx_sni - wfx_opw - wfx_res ) 
     262      IF  ( iom_use( "dmidyn"   ) ) CALL iom_put( "dmidyn"      ,   diag_dmi_dyn             )          ! Sea-ice mass change from dynamics 
     263      IF  ( iom_use( "dmiopw"   ) ) CALL iom_put( "dmiopw"      , - wfx_opw                  )          ! Sea-ice mass change through growth in open water 
     264      IF  ( iom_use( "dmibog"   ) ) CALL iom_put( "dmibog"      , - wfx_bog                  )          ! Sea-ice mass change through basal growth 
     265      IF  ( iom_use( "dmisni"   ) ) CALL iom_put( "dmisni"      , - wfx_sni                  )          ! Sea-ice mass change through snow-to-ice conversion 
     266      IF  ( iom_use( "dmisum"   ) ) CALL iom_put( "dmisum"      , - wfx_sum                  )          ! Sea-ice mass change through surface melting 
     267      IF  ( iom_use( "dmibom"   ) ) CALL iom_put( "dmibom"      , - wfx_bom                  )          ! Sea-ice mass change through bottom melting 
     268      IF  ( iom_use( "dmtsub"   ) ) CALL iom_put( "dmtsub"      , - wfx_sub                  )          ! Sea-ice mass change through evaporation and sublimation 
     269      IF  ( iom_use( "dmsspr"   ) ) CALL iom_put( "dmsspr"      , - wfx_spr                  )          ! snow mass change through snow fall 
     270      IF  ( iom_use( "dmsssi"   ) ) CALL iom_put( "dmsssi"      ,   wfx_sni*rhosn/rhoic      )          ! snow mass change through snow-to-ice conversion 
     271 
     272      IF  ( iom_use( "dmsmel"   ) ) CALL iom_put( "dmsmel"      , - wfx_snw_sum              )          ! snow mass change through melt 
     273      IF  ( iom_use( "dmsdyn"   ) ) CALL iom_put( "dmsdyn"      ,   diag_dms_dyn             )          ! snow mass change through dynamics 
     274 
     275      IF  ( iom_use( "hfxconbo" ) ) CALL iom_put( "hfxconbo"    ,   diag_fc_bo               )          ! bottom conduction flux 
     276      IF  ( iom_use( "hfxconsu" ) ) CALL iom_put( "hfxconsu"    ,   diag_fc_su               )          ! surface conduction flux 
     277 
     278      IF  ( iom_use( "wfxtot"   ) ) CALL iom_put( "wfxtot"      ,   -wfx_ice                 )          ! total freshwater flux from sea ice 
     279 
     280      IF  ( iom_use( "dmtxdyn"  ) ) CALL iom_put( "dmtxdyn"     ,   diag_dmtx_dyn            )          ! X-component of sea-ice mass transport 
     281      IF  ( iom_use( "dmtydyn"  ) ) CALL iom_put( "dmtydyn"     ,   diag_dmty_dyn            )          ! Y-component of sea-ice mass transport 
     282 
     283      IF  ( iom_use( "utau_oi"  ) ) CALL iom_put( "utau_oi"     ,   diag_utau_oi*zswi        )          ! X-component of ocean stress on sea ice 
     284      IF  ( iom_use( "vtau_oi"  ) ) CALL iom_put( "vtau_oi"     ,   diag_vtau_oi*zswi        )          ! Y-component of ocean stress on sea ice 
     285 
     286      IF  ( iom_use( "dssh_dx"  ) ) CALL iom_put( "dssh_dx"     ,   diag_dssh_dx*zswi        )          ! Sea-surface tilt term in force balance (x-component) 
     287      IF  ( iom_use( "dssh_dy"  ) ) CALL iom_put( "dssh_dy"     ,   diag_dssh_dy*zswi        )          ! Sea-surface tilt term in force balance (y-component) 
     288 
     289      IF  ( iom_use( "corstrx"  ) ) CALL iom_put( "corstrx"     ,   diag_corstrx*zswi        )          ! Coriolis force term in force balance (x-component) 
     290      IF  ( iom_use( "corstry"  ) ) CALL iom_put( "corstry"     ,   diag_corstry*zswi        )          ! Coriolis force term in force balance (y-component) 
     291 
     292      IF  ( iom_use( "intstrx"  ) ) CALL iom_put( "intstrx"     ,   diag_intstrx*zswi        )          ! Internal force term in force balance (x-component) 
     293      IF  ( iom_use( "intstry"  ) ) CALL iom_put( "intstry"     ,   diag_intstry*zswi        )          ! Internal force term in force balance (y-component) 
     294 
     295      IF  ( iom_use( "normstr"  ) ) CALL iom_put( "normstr"     ,   diag_sig1   *zswi        )          ! Normal stress 
     296      IF  ( iom_use( "sheastr"  ) ) CALL iom_put( "sheastr"     ,   diag_sig2   *zswi        )          ! Shear stress 
     297 
     298      !-------------------------------- 
     299      ! Global ice diagnostics (SIMIP) 
     300      !-------------------------------- 
     301 
     302      IF ( iom_use( "NH_icearea" ) .OR. iom_use( "NH_icevolu" ) .OR. iom_use( "NH_iceextt" ) )   THEN   ! NH integrated diagnostics 
     303  
     304         WHERE( fcor > 0._wp ); zswi(:,:) = 1.0e-12 
     305         ELSEWHERE            ; zswi(:,:) = 0. 
     306         END WHERE  
     307 
     308         zdiag_area_nh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) ) 
     309         zdiag_volu_nh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) ) 
     310 
     311         WHERE( fcor > 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12 
     312         ELSEWHERE                              ; zswi(:,:) = 0. 
     313         END WHERE  
     314 
     315         zdiag_extt_nh = glob_sum( zswi(:,:) * e12t(:,:) ) 
     316 
     317         IF ( iom_use( "NH_icearea" ) ) CALL iom_put( "NH_icearea" ,  zdiag_area_nh  ) 
     318         IF ( iom_use( "NH_icevolu" ) ) CALL iom_put( "NH_icevolu" ,  zdiag_volu_nh  ) 
     319         IF ( iom_use( "NH_iceextt" ) ) CALL iom_put( "NH_iceextt" ,  zdiag_extt_nh  ) 
     320 
     321      ENDIF 
     322 
     323      IF ( iom_use( "SH_icearea" ) .OR. iom_use( "SH_icevolu" ) .OR. iom_use( "SH_iceextt" ) )   THEN   ! SH integrated diagnostics 
     324 
     325         WHERE( fcor < 0._wp ); zswi(:,:) = 1.0e-12;  
     326         ELSEWHERE            ; zswi(:,:) = 0. 
     327         END WHERE  
     328 
     329         zdiag_area_sh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) )  
     330         zdiag_volu_sh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) ) 
     331 
     332         WHERE( fcor < 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12 
     333         ELSEWHERE                              ; zswi(:,:) = 0. 
     334         END WHERE  
     335 
     336         zdiag_extt_sh = glob_sum( zswi(:,:) * e12t(:,:) ) 
     337 
     338         IF ( iom_use( "SH_icearea" ) ) CALL iom_put( "SH_icearea", zdiag_area_sh ) 
     339         IF ( iom_use( "SH_icevolu" ) ) CALL iom_put( "SH_icevolu", zdiag_volu_sh ) 
     340         IF ( iom_use( "SH_iceextt" ) ) CALL iom_put( "SH_iceextt", zdiag_extt_sh ) 
     341 
     342      ENDIF  
     343 
     344      !-------------------------------- 
     345      ! Fluxes through straits (SIMIP) 
     346      !-------------------------------- 
     347      ! 
     348      ! Valid only for ORCA-like grids 
     349      ! 
     350      ! 4 Arctic passages are considered (Fram, CAA, Barents, Bering; see Notz et al (GMD 2016) for definitions) 
     351      ! 
     352      ! Fram and Bering  straits are easy because they follow parallels 
     353      ! Barents and Canadian Arctic Archipelago are less easy because they do not, which is why they look so awful. 
     354      !  
     355 
     356      IF ( iom_use( "strait_arfl" ) .OR. iom_use( "strait_mifl" ) .OR. iom_use( "strait_msfl" ) .AND. cp_cfg == "orca" ) THEN 
     357 
     358         zdiag_area_strait(:) = 0._wp   ;   zdiag_mice_strait(:) = 0._wp   ;   zdiag_msno_strait(:) = 0._wp 
     359    
     360         !------------------------------ 
     361         ! === Fram & Bering Straits === 
     362         !------------------------------ 
     363           
     364         SELECT CASE ( jp_cfg )  
     365           
     366         CASE ( 2 )   ! --- ORCA2 
     367           
     368            ! Fram Strait   (i_strait = 1) 
     369            ji0(1) = 133   ;   ji1(1) = 136 
     370            jj0(1) = 136  
     371 
     372            ! Bering Strait (i_strait = 4) 
     373            ji0(4) = 55    ;   ji1(4) = 56 
     374            jj0(4) = 122 
     375           
     376         CASE ( 1 )   ! --- eORCA1 
     377           
     378            ! Fram Strait 
     379            ji0(1) = 268   ;   ji1(1) = 277 
     380            jj0(1) = 311 
     381 
     382            ! Bering Strait 
     383            ji0(4) = 113   ;   jj1(4) = 115 
     384            jj0(4) = 285 
     385            
     386         END SELECT 
     387           
     388         DO i_strait = 1, 4, 3 
     389 
     390            DO ji = mi0(ji0), mi1(ji1) 
     391               jj = mj0(jj0) 
     392    
     393               zdiag_area_strait(i_strait) = zdiag_area_strait(i_strait)                                 &     ! --- ice area flux --- 
     394                   &                + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         &     ! northwards (positive) flow 
     395                   &                + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 )               ! southwards (negative) flow 
     396       
     397               zdiag_mice_strait(i_strait) = zdiag_mice_strait(i_strait) + rhoic *                       &     ! --- ice mass flux  --- 
     398                   &                ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         &  
     399                   &                + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) )          
     400       
     401               zdiag_msno_strait(i_strait) = zdiag_msno_strait(i_strait) + rhosn *                       &     ! --- snow mass flux --- 
     402                   &                ( vt_s(ji,jj-1) * e12t(ji,jj-1) * MAX( v_ice(ji,jj-1), 0.0 )         &  
     403                   &                + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( v_ice(ji,jj-1), 0.0 ) )   
     404    
     405            END DO 
     406 
     407         END DO 
     408 
     409         !--------------------- 
     410         ! === Barents opening 
     411         !--------------------- 
     412    
     413         SELECT CASE ( jp_cfg )  
     414 
     415            CASE ( 1 )   ! 'eORCA1' 
     416 
     417               Nu = 11   ! U-Flow 
     418               zui(1:Nu) = (/ 282,283,284,285,286,286,287,288,289,290,292/) 
     419               zuj(1:Nu) = (/ 308,307,306,305,304,303,302,301,300,299,298/) 
     420 
     421               Nv = 9    ! V-Flow 
     422               zvi(1:Nv) = (/ 282,283,284,285,286,287,288,289,290/) 
     423               zvj(1:Nv) = (/ 308,307,306,305,303,302,301,300,299/) 
     424 
     425            CASE ( 2 )   ! 'ORCA2' 
     426 
     427               Nu = 5    ! U-Flow 
     428               zui(1:Nu) = (/ 141,142,142,143,144 /) 
     429               zuj(1:Nu) = (/ 134,133,132,131,130 /) 
     430 
     431               Nv = 4    ! V-Flow 
     432               zvi(1:Nv) = (/ 140,141,142,143 /) 
     433               zvj(1:Nv) = (/ 135,134,132,131 /) 
     434 
     435         END SELECT 
     436 
     437         ! Barents U-flow 
     438         zfarea_u = 0._wp   ;   zfmice_u = 0._wp   ;   zfmsno_u = 0._wp 
     439    
     440         DO ii = 1, Nu 
     441 
     442            ji = mi0(zui(ii)) 
     443            jj = mj0(zuj(ii)) 
     444 
     445            zfarea_u          = zfarea_u                                           & ! --- ice area zonal flux --- 
     446                &             + at_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         & ! --- northward 
     447                &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 )           ! --- southward 
     448            zfmice_u          = zfmice_u + rhoic *                                 & ! --- ice mass zonal flux ---  
     449                &             ( vt_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
     450                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
     451            zfmsno_u          = zfmsno_u + rhosn *                                 & ! --- snow mass zonal flux ---  
     452                &             ( vt_s(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
     453                &             + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
     454         END DO 
     455    
     456         ! Barents V-flow 
     457         zfarea_v = 0._wp   ;   zfmice_v = 0._wp   ;   zfmsno_v = 0._wp 
     458 
     459         DO ii  = 1, Nv   
     460 
     461            ji = mi0(zvi(ii)) 
     462            jj = mj0(zvj(ii)) 
     463 
     464            zfarea_v          = zfarea_v                                           & ! --- ice area meridian flux --- 
     465                &             + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! --- eastward 
     466                &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 )           ! --- westward 
     467            zfmice_v          = zfmice_v + rhoic *                                 & ! --- ice mass meridian flux --- 
     468                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
     469                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
     470            zfmsno_v          = zfmsno_v + rhosn *                                 & ! --- snow mass meridian flux --- 
     471                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
     472                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
     473         END DO 
     474    
     475         ! Sum Barents U-/V- contributions 
     476         zdiag_area_strait(3) = zfarea_u + zfarea_v  
     477         zdiag_mice_strait(3) = zfmice_u + zfmice_v 
     478         zdiag_msno_strait(3) = zfmsno_u + zfmsno_v 
     479 
     480         !--------------------- 
     481         ! === CAA throughflow 
     482         !--------------------- 
     483    
     484         SELECT CASE ( jp_cfg )  
     485 
     486            CASE ( 1 )   ! eORCA1 
     487 
     488               ! V-flow through Nares Strait 
     489               Nv = 4 
     490               zvi(1:Nv) = (/ 254,255,256,257 /) 
     491               zvj(1:Nv) = (/ 317,317,317,317 /) 
     492 
     493               ! U-flow through Queen Elisabeth Islands and McClure straits 
     494               Nu = 8  
     495               zui(1:Nu) = (/ 231,231,231,  132,132,132,132,132  /) 
     496               zuj(1:Nu) = (/ 328,329,330,  318,319,320,321,322  /) 
     497 
     498               zfarea_u = 0._wp   ;   zfmice_u = 0._wp   ;   zfmsno_u = 0._wp 
     499 
     500               DO ii = 1, Nu 
     501 
     502                  ji = mi0(zui(ii)) 
     503                  jj = mj0(zuj(ii)) 
     504 
     505                  zfarea_u          = zfarea_u                                                           & ! --- ice area zonal flux --- 
     506                      &             + at_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         & ! --- eastward 
     507                      &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 )           ! --- westward 
     508                  zfmice_u          = zfmice_u + rhoic *                                                 & ! --- ice mass zonal flux ---  
     509                      &             ( vt_i(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
     510                      &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
     511                  zfmsno_u          = zfmsno_u + rhosn *                                                 & ! --- snow mass zonal flux ---  
     512                      &             ( vt_s(ji-1,jj) * e12t(ji-1,jj) * MAX( u_ice(ji-1,jj), 0.0 )         &    
     513                      &             + vt_s(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji-1,jj), 0.0 ) )           
     514 
     515               END DO 
     516 
     517 
     518            CASE ( 2 )   ! ORCA2 
     519 
     520               ! V-flow through Nares Strait 
     521               Nv = 2 
     522               zvi(1:Nv) = (/ 117,118 /) 
     523               zvj(1:Nv) = (/ 145,145 /) 
     524 
     525               ! U-flow through Queen Elisabeth Islands and McClure straits (not resolved in ORCA2) 
     526               zfarea_u = 0._wp   ;   zfmice_u = 0._wp   ;   zfmsno_u = 0._wp 
     527 
     528            END SELECT 
     529    
     530         ! V-flow through Nares Strait 
     531         zfarea_v = 0._wp   ;   zfmice_v = 0._wp   ;   zfmsno_v = 0._wp 
     532    
     533         DO ii = 1, Nv   
     534 
     535            ji = mi0(zvi(ii)) 
     536            jj = mj0(zvj(ii)) 
     537 
     538            zfarea_v          = zfarea_v                                           & ! --- ice area meridian flux --- 
     539                &             + at_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! --- eastward 
     540                &             + at_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 )           ! --- westward 
     541            zfmice_v          = zfmice_v + rhoic *                                 & ! --- ice mass meridian flux --- 
     542                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
     543                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
     544            zfmsno_v          = zfmsno_v + rhosn *                                 & ! --- snow mass meridian flux --- 
     545                &             ( vt_i(ji,jj-1) * e12t(ji,jj-1) * MAX( u_ice(ji,jj-1), 0.0 )         & ! 
     546                &             + vt_i(ji,jj  ) * e12t(ji,jj)   * MIN( u_ice(ji,jj-1), 0.0 ) )         ! 
     547 
     548         END DO 
     549 
     550         ! Sum U/V contributions 
     551         zdiag_area_strait(2) = zfarea_u + zfarea_v  
     552         zdiag_mice_strait(2) = zfmice_u + zfmice_v 
     553         zdiag_msno_strait(2) = zfmsno_u + zfmsno_v 
     554    
     555         ! === Ncdf output 
     556         IF ( iom_use("strait_arfl") ) CALL iom_put( "strait_arfl", zdiag_area_strait ) 
     557         IF ( iom_use("strait_mifl") ) CALL iom_put( "strait_mifl", zdiag_mice_strait ) 
     558         IF ( iom_use("strait_msfl") ) CALL iom_put( "strait_msfl", zdiag_msno_strait )  
     559 
     560         WRITE(numout,*) " area flx ", zdiag_area_strait(:) 
     561         WRITE(numout,*) " mice flx ", zdiag_mice_strait(:) 
     562         WRITE(numout,*) " msno flx ", zdiag_msno_strait(:) 
     563 
     564      ENDIF 
    578565 
    579566      !     !  Create an output files (output.lim.abort.nc) if S < 0 or u > 20 m/s 
Note: See TracChangeset for help on using the changeset viewer.