Changeset 4730


Ignore:
Timestamp:
2014-07-25T09:39:23+02:00 (6 years ago)
Author:
vancop
Message:

coupled interface modifications for LIM3

Location:
branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/CONFIG/SHARED/namelist_ref

    r4699 r4730  
    234234   nn_ice      = 2         !  =0 no ice boundary condition   , 
    235235                           !  =1 use observed ice-cover      , 
    236                            !  =2 ice-model used                         ("key_lim3" or "key_lim2) 
     236                           !  =2 ice-model used                         ("key_lim3" or "key_lim2") 
    237237   nn_ice_embd = 1         !  =0 levitating ice (no mass exchange, concentration/dilution effect) 
    238238                           !  =1 levitating ice with mass and salt exchange but no presure effect 
     
    250250   nn_lsm  = 0             !  =0 land/sea mask for input fields is not applied (keep empty land/sea mask filename field) , 
    251251                           !  =1:n number of iterations of land/sea mask application for input fields (fill land/sea mask filename field) 
    252    cn_iceflx = 'linear'    !  redistribution of solar input into ice categories during coupling ice/atm. 
     252   nn_limflx = -1          !  LIM3 Multi-category heat flux formulation (use -1 if LIM3 is not used) 
     253                           !  =-1  Use per-category fluxes, bypass redistributor, forced mode only, not yet implemented coupled 
     254                           !  = 0  Average per-category fluxes (forced and coupled mode) 
     255                           !  = 1  Average and redistribute per-category fluxes, forced mode only, not yet implemented coupled 
     256                           !  = 2  Redistribute a single flux over categories (coupled mode only) 
    253257/ 
    254258!----------------------------------------------------------------------- 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/LIM_SRC_2/limsbc_2.F90

    r4306 r4730  
    101101      !!              - fr_i    : ice fraction 
    102102      !!              - tn_ice  : sea-ice surface temperature 
    103       !!              - alb_ice : sea-ice alberdo (lk_cpl=T) 
     103      !!              - alb_ice : sea-ice albedo (lk_cpl=T) 
    104104      !! 
    105105      !! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90. 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90

    r4688 r4730  
    115115 
    116116      CALL lim_istate_init     !  reading the initials parameters of the ice 
    117  
    118 # if defined key_coupled 
    119       albege(:,:)   = 0.8 * tms(:,:) 
    120 # endif 
    121117 
    122118      ! surface temperature 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

    r4688 r4730  
    102102      !! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90. 
    103103      !!              Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 
     104      !!              These refs are now obsolete since everything has been revised 
     105      !!              The ref should be Rousset et al., 2015? 
    104106      !!--------------------------------------------------------------------- 
    105       INTEGER, INTENT(in) ::   kt    ! number of iteration 
    106       ! 
    107       INTEGER  ::   ji, jj, jl, jk           ! dummy loop indices 
    108       REAL(wp) ::   zinda, zemp      ! local scalars 
    109       REAL(wp) ::   zf_mass         ! Heat flux associated with mass exchange ice->ocean (W.m-2) 
    110       REAL(wp) ::   zfcm1           ! New solar flux received by the ocean 
    111       REAL(wp), POINTER, DIMENSION(:,:,:) ::   zalb, zalbp     ! 2D/3D workspace 
     107      INTEGER, INTENT(in) ::   kt                                   ! number of iteration 
     108      ! 
     109      INTEGER  ::   ji, jj, jl, jk                                  ! dummy loop indices 
     110      ! 
     111      REAL(wp) ::   zinda, zemp                                     !  local scalars 
     112      REAL(wp) ::   zf_mass                                         ! Heat flux associated with mass exchange ice->ocean (W.m-2) 
     113      REAL(wp) ::   zfcm1                                           ! New solar flux received by the ocean 
     114      ! 
     115      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zalb_cs, zalb_os     ! 2D/3D workspace 
    112116      !!--------------------------------------------------------------------- 
    113        
    114       IF( lk_cpl )   CALL wrk_alloc( jpi, jpj, jpl, zalb, zalbp ) 
    115117 
    116118      ! make calls for heat fluxes before it is modified 
     
    134136            ! Solar heat flux reaching the ocean = zfcm1 (W.m-2)  
    135137            !--------------------------------------------------- 
    136             IF( lk_cpl ) THEN ! be carfeful: not been tested yet 
     138            IF( lk_cpl ) THEN ! be careful: not been tested yet 
    137139               ! original line 
    138140               zfcm1 = qsr_tot(ji,jj) 
    139                !!!zfcm1 = qsr_tot(ji,jj) + ftr_ice(ji,jj) * ( 1._wp - pfrld(ji,jj) ) / ( 1._wp - zinda + zinda * iatte(ji,jj) ) 
     141               !!! LIM2 version zqsr = qsr_tot(ji,jj) + ( fstric(ji,jj) - qsr_ice(ji,jj,1) ) * ( 1.0 - pfrld(ji,jj) ) 
    140142               DO jl = 1, jpl 
    141                   zfcm1 = zfcm1 - ( qsr_ice(ji,jj,jl) - ftr_ice(ji,jj,jl) ) * old_a_i(ji,jj,jl) 
     143                  zfcm1 = zfcm1 + ( ftr_ice(ji,jj,jl) - qsr_ice(ji,jj,jl) ) * old_a_i(ji,jj,jl) 
    142144               END DO 
    143145            ELSE 
    144                !!!zfcm1   = pfrld(ji,jj) * qsr(ji,jj)  + & 
    145                !!!     &    ( 1._wp - pfrld(ji,jj) ) * ftr_ice(ji,jj) / ( 1._wp - zinda + zinda * iatte(ji,jj) ) 
     146               !!! LIM2 version zqsr = pfrld(ji,jj) * qsr(ji,jj)  + ( 1.  - pfrld(ji,jj) ) * fstric(ji,jj) 
    146147               zfcm1   = pfrld(ji,jj) * qsr(ji,jj) 
    147148               DO jl = 1, jpl 
     
    215216 
    216217      !------------------------------------------------! 
    217       !    Computation of snow/ice and ocean albedo    ! 
     218      !    Snow/ice albedo (only if sent to coupler)   ! 
    218219      !------------------------------------------------! 
    219220      IF( lk_cpl ) THEN          ! coupled case 
    220          CALL albedo_ice( t_su, ht_i, ht_s, zalbp, zalb )                  ! snow/ice albedo 
    221          alb_ice(:,:,:) =  0.5_wp * zalbp(:,:,:) + 0.5_wp * zalb (:,:,:)   ! Ice albedo (mean clear and overcast skys) 
     221 
     222            CALL wrk_alloc( jpi, jpj, jpl, zalb_cs, zalb_os ) 
     223 
     224            CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos 
     225 
     226            alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     227 
     228            CALL wrk_dealloc( jpi, jpj, jpl, zalb_cs, zalb_os ) 
     229 
    222230      ENDIF 
    223231 
     
    229237         CALL prt_ctl( tab3d_1=tn_ice, clinfo1=' lim_sbc: tn_ice : ', kdim=jpl ) 
    230238      ENDIF 
    231       ! 
    232       IF( lk_cpl )   CALL wrk_dealloc( jpi, jpj, jpl, zalb, zalbp ) 
    233       !  
     239 
    234240   END SUBROUTINE lim_sbc_flx 
    235241 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90

    r4306 r4730  
    5656 
    5757#if defined key_lim3 || defined key_lim2  
    58    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qns_ice       !: non solar heat flux over ice                  [W/m2] 
    59    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qsr_ice       !: solar heat flux over ice                      [W/m2] 
    60    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qsr_ice_mean  !: dauly mean solar heat flux over ice       [W/m2] 
    61    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qla_ice       !: latent flux over ice                          [W/m2] 
    62    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   dqla_ice      !: latent sensibility over ice                 [W/m2/K] 
    63    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   dqns_ice      !: non solar heat flux over ice (LW+SEN+LA)    [W/m2/K] 
    64    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tn_ice        !: ice surface temperature                          [K] 
    65    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   alb_ice       !: albedo of ice 
     58   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qns_ice            !: non solar heat flux over ice                  [W/m2] 
     59   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qsr_ice            !: solar heat flux over ice                      [W/m2] 
     60   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qsr_ice_mean       !: daily mean solar heat flux over ice       [W/m2] 
     61   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qla_ice            !: latent flux over ice                          [W/m2] 
     62   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   dqla_ice           !: latent sensibility over ice                 [W/m2/K] 
     63   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   dqns_ice           !: non solar heat flux over ice (LW+SEN+LA)    [W/m2/K] 
     64   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tn_ice             !: ice surface temperature                          [K] 
     65   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   alb_ice            !: ice albedo                                       [-] 
    6666 
    67    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   utau_ice  !: atmos-ice u-stress. VP: I-pt ; EVP: U,V-pts   [N/m2] 
    68    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   vtau_ice  !: atmos-ice v-stress. VP: I-pt ; EVP: U,V-pts   [N/m2] 
    69    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr1_i0    !: 1st Qsr fraction penetrating inside ice cover    [-] 
    70    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr2_i0    !: 2nd Qsr fraction penetrating inside ice cover    [-] 
    71    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   emp_ice   !: sublimation-snow budget over ice             [kg/m2] 
     67   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   utau_ice           !: atmos-ice u-stress. VP: I-pt ; EVP: U,V-pts   [N/m2] 
     68   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   vtau_ice           !: atmos-ice v-stress. VP: I-pt ; EVP: U,V-pts   [N/m2] 
     69   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr1_i0             !: Solar surface transmission parameter, thick ice  [-] 
     70   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr2_i0             !: Solar surface transmission parameter, thin ice   [-] 
     71   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   emp_ice            !: sublimation-snow budget over ice             [kg/m2] 
    7272 
    7373# if defined key_lim3 
    74    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   tatm_ice  !: air temperature 
     74   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   tatm_ice           !: air temperature [K] 
     75   REAL(wp), PUBLIC, SAVE                                ::   cldf_ice = 0.81    !: cloud fraction over sea ice, summer CLIO value   [-] 
    7576# endif 
    7677 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90

    r4306 r4730  
    4545   !                                             !: =1 levitating ice with mass and salt exchange but no presure effect 
    4646   !                                             !: =2 embedded sea-ice (full salt and mass exchanges and pressure) 
     47   INTEGER , PUBLIC :: nn_limflx        !: LIM3 Multi-category heat flux formulation 
     48   !                                             !: =-1  Use of per-category fluxes 
     49   !                                             !: = 0  Average per-category fluxes 
     50   !                                             !: = 1  Average then redistribute per-category fluxes 
     51   !                                             !: = 2  Redistribute a single flux over categories 
    4752   INTEGER , PUBLIC ::   nn_fwb         !: FreshWater Budget:  
    4853   !                                             !:  = 0 unchecked  
     
    5560   LOGICAL , PUBLIC ::   ln_icebergs    !: Icebergs 
    5661   ! 
    57    CHARACTER (len=8), PUBLIC :: cn_iceflx  !: Flux handling over ice categories 
    58    LOGICAL, PUBLIC :: ln_iceflx_ave     ! Average heat fluxes over all ice categories 
    59    LOGICAL, PUBLIC :: ln_iceflx_linear  ! Redistribute mean heat fluxes over all ice categories, using ice temperature and albedo 
    60    ! 
    61    INTEGER , PUBLIC ::   nn_lsm        !: Number of iteration if seaoverland is applied 
     62   INTEGER , PUBLIC ::   nn_lsm         !: Number of iteration if seaoverland is applied 
    6263   !!---------------------------------------------------------------------- 
    6364   !!              Ocean Surface Boundary Condition fields 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_clio.F90

    r4624 r4730  
    398398 
    399399 
    400    SUBROUTINE blk_ice_clio(  pst   , palb_cs, palb_os ,       & 
     400   SUBROUTINE blk_ice_clio(  pst   , palb_cs, palb_os , palb  & 
    401401      &                      p_taui, p_tauj, p_qns , p_qsr,   & 
    402402      &                      p_qla , p_dqns, p_dqla,          & 
     
    427427      !!---------------------------------------------------------------------- 
    428428      REAL(wp), INTENT(in   ), DIMENSION(:,:,:)   ::   pst      ! ice surface temperature                   [Kelvin] 
    429       REAL(wp), INTENT(in   ), DIMENSION(:,:,:)   ::   palb_cs  ! ice albedo (clear    sky) (alb_ice_cs)         [%] 
    430       REAL(wp), INTENT(in   ), DIMENSION(:,:,:)   ::   palb_os  ! ice albedo (overcast sky) (alb_ice_os)         [%] 
     429      REAL(wp), INTENT(in   ), DIMENSION(:,:,:)   ::   palb_cs  ! ice albedo (clear    sky) (alb_ice_cs)         [-] 
     430      REAL(wp), INTENT(in   ), DIMENSION(:,:,:)   ::   palb_os  ! ice albedo (overcast sky) (alb_ice_os)         [-] 
     431      REAL(wp), INTENT(  out), DIMENSION(:,:,:)   ::   palb     ! ice albedo (actual value)                      [-] 
    431432      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_taui   ! surface ice stress at I-point (i-component) [N/m2] 
    432433      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_tauj   ! surface ice stress at I-point (j-component) [N/m2] 
     
    438439      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_tpr    ! total precipitation          (T-point)   [Kg/m2/s] 
    439440      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_spr    ! solid precipitation          (T-point)   [Kg/m2/s] 
    440       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_fr1    ! 1sr fraction of qsr penetration in ice         [%] 
    441       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_fr2    ! 2nd fraction of qsr penetration in ice         [%] 
     441      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_fr1    ! 1sr fraction of qsr penetration in ice         [-] 
     442      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_fr2    ! 2nd fraction of qsr penetration in ice         [-] 
    442443      CHARACTER(len=1), INTENT(in   )             ::   cd_grid  ! type of sea-ice grid ("C" or "B" grid) 
    443444      INTEGER, INTENT(in   )                      ::   pdim     ! number of ice categories 
     
    542543      !-----------------------------------------------------------! 
    543544      CALL blk_clio_qsr_ice( palb_cs, palb_os, p_qsr ) 
     545       
     546      DO jl = 1, jpl 
     547         palb(:,:,jl) = ( palb_cs(:,:,jl) * ( 1.e0 - sf(jp_ccov)%fnow(ji,jj,1) )   & 
     548            &         +   palb_os(:,:,jl) * sf(jp_ccov)%fnow(ji,jj,1) ) 
     549      END DO 
    544550 
    545551      !                                     ! ========================== ! 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90

    r4689 r4730  
    538538      REAL(wp) ::   zcoef_wnorm, zcoef_wnorm2, zcoef_dqlw, zcoef_dqla, zcoef_dqsb 
    539539      REAL(wp) ::   zztmp                                        ! temporary variable 
    540       REAL(wp) ::   zcoef_frca                                   ! fractional cloud amount 
    541540      REAL(wp) ::   zwnorm_f, zwndi_f , zwndj_f                  ! relative wind module and components at F-point 
    542541      REAL(wp) ::             zwndi_t , zwndj_t                  ! relative wind components at T-point 
     
    562561      zcoef_dqla   = -Ls * Cice * 11637800. * (-5897.8) 
    563562      zcoef_dqsb   = rhoa * cpa * Cice 
    564       zcoef_frca   = 1.0  - 0.3 
    565       ! MV 2014 the proper cloud fraction (mean summer months from the CLIO climato, NH+SH) is 0.19 
    566       zcoef_frca   = 1.0  - 0.19 
    567563 
    568564!!gm brutal.... 
     
    694690     
    695691!CDIR COLLAPSE 
    696       p_fr1(:,:) = ( 0.18 * ( 1.0 - zcoef_frca ) + 0.35 * zcoef_frca ) 
    697 !CDIR COLLAPSE 
    698       p_fr2(:,:) = ( 0.82 * ( 1.0 - zcoef_frca ) + 0.65 * zcoef_frca ) 
     692      p_fr1(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 
     693!CDIR COLLAPSE 
     694      p_fr2(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 
    699695        
    700696!CDIR COLLAPSE 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r4664 r4730  
    210210      !!             ***  ROUTINE sbc_cpl_init  *** 
    211211      !! 
    212       !! ** Purpose :   Initialisation of send and recieved information from 
     212      !! ** Purpose :   Initialisation of send and received information from 
    213213      !!                the atmospheric component 
    214214      !! 
     
    13141314      END SELECT 
    13151315 
    1316       !    Ice Qsr penetration used (only?)in lim2 or lim3  
    1317       ! fraction of net shortwave radiation which is not absorbed in the thin surface layer  
    1318       ! and penetrates inside the ice cover ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 ) 
     1316      ! Surface transimission parameter io (Maykut Untersteiner , 1971 ; Ebert and Curry, 1993 ) 
     1317      ! Used for LIM2 and LIM3 
    13191318      ! Coupled case: since cloud cover is not received from atmosphere  
    1320       !               ===> defined as constant value -> definition done in sbc_cpl_init 
    1321       fr1_i0(:,:) = 0.18 
    1322       fr2_i0(:,:) = 0.82 
    1323  
     1319      !               ===> used prescribed cloud fraction representative for polar oceans in summer (0.81) 
     1320      fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 
     1321      fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 
    13241322 
    13251323      CALL wrk_dealloc( jpi,jpj, zcptn, ztmp, zicefr ) 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90

    r4689 r4730  
    1212   !!            3.4  ! 2011-01  (A Porter)  dynamical allocation 
    1313   !!             -   ! 2012-10  (C. Rousset)  add lim_diahsb 
     14   !!            3.6  ! 2014-07  (M. Vancoppenolle, G. Madec, O. Marti) revise coupled interface 
    1415   !!---------------------------------------------------------------------- 
    1516#if defined key_lim3 
     
    5960   USE prtctl          ! Print control 
    6061   USE lib_fortran     !  
    61    USE cpl_oasis3, ONLY : lk_cpl 
    6262 
    6363#if defined key_bdy  
     
    6969 
    7070   PUBLIC sbc_ice_lim  ! routine called by sbcmod.F90 
    71    PUBLIC lim_prt_state 
    7271    
    7372   !! * Substitutions 
     
    8079   !!---------------------------------------------------------------------- 
    8180CONTAINS 
    82  
    83    FUNCTION fice_cell_ave ( ptab) 
    84       !!-------------------------------------------------------------------------- 
    85       !! * Compute average over categories, for grid cell (ice covered and free ocean) 
    86       !!-------------------------------------------------------------------------- 
    87       REAL (wp), DIMENSION (jpi,jpj) :: fice_cell_ave 
    88       REAL (wp), DIMENSION (jpi,jpj,jpl), INTENT (in) :: ptab 
    89       INTEGER :: jl ! Dummy loop index 
    90        
    91       fice_cell_ave (:,:) = 0.0_wp 
    92        
    93       DO jl = 1, jpl 
    94          fice_cell_ave (:,:) = fice_cell_ave (:,:) & 
    95             &                  + a_i (:,:,jl) * ptab (:,:,jl) 
    96       END DO 
    97        
    98    END FUNCTION fice_cell_ave 
    99     
    100    FUNCTION fice_ice_ave ( ptab) 
    101       !!-------------------------------------------------------------------------- 
    102       !! * Compute average over categories, for ice covered part of grid cell 
    103       !!-------------------------------------------------------------------------- 
    104       REAL (kind=wp), DIMENSION (jpi,jpj) :: fice_ice_ave 
    105       REAL (kind=wp), DIMENSION (jpi,jpj,jpl), INTENT(in) :: ptab 
    106  
    107       fice_ice_ave (:,:) = 0.0_wp 
    108       WHERE ( at_i (:,:) .GT. 0.0_wp ) fice_ice_ave (:,:) = fice_cell_ave ( ptab (:,:,:)) / at_i (:,:) 
    109  
    110    END FUNCTION fice_ice_ave 
    11181 
    11282   !!====================================================================== 
     
    133103      !!--------------------------------------------------------------------- 
    134104      INTEGER, INTENT(in) ::   kt      ! ocean time step 
    135       INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE) 
     105      INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE, =5 COUPLED) 
    136106      !! 
    137       INTEGER  ::   ji, jj, jl, jk      ! dummy loop index 
     107      INTEGER  ::   jl      ! dummy loop index 
    138108      REAL(wp) ::   zcoef   ! local scalar 
    139       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice_os, zalb_ice_cs  ! albedo of the ice under overcast/clear sky 
    140       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice      ! mean albedo of ice (for coupled) 
    141  
    142       REAL(wp), POINTER, DIMENSION(:,:) :: zalb_ice_all    ! Mean albedo over all categories 
    143       REAL(wp), POINTER, DIMENSION(:,:) :: ztem_ice_all    ! Mean temperature over all categories 
    144        
    145       REAL(wp), POINTER, DIMENSION(:,:) :: z_qsr_ice_all   ! Mean solar heat flux over all categories 
    146       REAL(wp), POINTER, DIMENSION(:,:) :: z_qns_ice_all   ! Mean non solar heat flux over all categories 
    147       REAL(wp), POINTER, DIMENSION(:,:) :: z_qla_ice_all   ! Mean latent heat flux over all categories 
    148       REAL(wp), POINTER, DIMENSION(:,:) :: z_dqns_ice_all  ! Mean d(qns)/dT over all categories 
    149       REAL(wp), POINTER, DIMENSION(:,:) :: z_dqla_ice_all  ! Mean d(qla)/dT over all categories 
    150       REAL(wp) ::   ztmelts           ! clem 2014: for HC diags 
    151       REAL(wp) ::   epsi20 = 1.e-20   ! 
     109      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_os, zalb_cs  ! ice albedo under overcast/clear sky 
     110      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice      ! mean ice albedo of ice (for coupled) 
    152111      !!---------------------------------------------------------------------- 
    153112 
    154       !- O.M. : why do we allocate all these arrays even when MOD( kt-1, nn_fsbc ) /= 0 ????? 
    155  
    156113      IF( nn_timing == 1 )  CALL timing_start('sbc_ice_lim') 
    157  
    158       CALL wrk_alloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs, zalb_ice ) 
    159  
    160       IF( lk_cpl ) THEN 
    161          IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) & 
    162             &   CALL wrk_alloc( jpi,jpj, ztem_ice_all, zalb_ice_all, z_qsr_ice_all, z_qns_ice_all, z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all) 
    163       ENDIF 
    164114 
    165115      IF( kt == nit000 ) THEN 
     
    171121         ! 
    172122         IF( ln_nicep ) THEN      ! control print at a given point 
    173             jiindx = 15    ;   jjindx =  44 
     123            jiindx = 177   ;   jjindx = 112 
    174124            IF(lwp) WRITE(numout,*) ' The debugging point is : jiindx : ',jiindx, ' jjindx : ',jjindx 
    175125         ENDIF 
     
    184134         u_oce(:,:) = ssu_m(:,:)                     ! mean surface ocean current at ice velocity point 
    185135         v_oce(:,:) = ssv_m(:,:)                     ! (C-grid dynamics :  U- & V-points as the ocean) 
    186  
    187          ! masked sea surface freezing temperature [Kelvin] 
    188          t_bo(:,:) = ( tfreez( sss_m ) +  rt0 ) * tmask(:,:,1) + rt0 * ( 1. - tmask(:,:,1) ) 
    189  
    190          CALL albedo_ice( t_su, ht_i, ht_s, zalb_ice_cs, zalb_ice_os )  ! ... ice albedo 
    191  
     136         ! 
     137         t_bo(:,:) = eos_fzp( sss_m ) +  rt0         ! masked sea surface freezing temperature [Kelvin] 
     138         !                                           ! (set to rt0 over land) 
     139 
     140         !                                           ! Ice albedo 
     141         CALL wrk_alloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice )       
     142 
     143         CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos 
     144 
     145         SELECT CASE( kblk ) 
     146         CASE( 4 , 5 )   ! CORE and COUPLED bulk formulations 
     147 
     148            ! albedo depends on cloud fraction because of non-linear spectral effects 
     149            zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     150            ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo  
     151            ! zalb_ice is computed within the bulk routine 
     152             
     153         END SELECT 
     154          
     155         !                                           ! Mask sea ice surface temperature 
    192156         DO jl = 1, jpl 
    193157            t_su(:,:,jl) = t_su(:,:,jl) +  rt0 * ( 1. - tmask(:,:,1) ) 
    194158         END DO 
    195  
    196          IF ( ln_cpl ) zalb_ice (:,:,:) = 0.5 * ( zalb_ice_cs (:,:,:) +  zalb_ice_os (:,:,:) ) 
    197           
    198          IF( lk_cpl ) THEN 
    199             IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) THEN 
    200                ! 
    201                ! Compute mean albedo and temperature 
    202                zalb_ice_all (:,:) = fice_ice_ave ( zalb_ice (:,:,:) )  
    203                ztem_ice_all (:,:) = fice_ice_ave ( tn_ice   (:,:,:) )  
    204                ! 
    205             ENDIF 
    206          ENDIF 
    207                                                ! Bulk formulea - provides the following fields: 
     159      
     160         ! Bulk formulae  - provides the following fields: 
    208161         ! utau_ice, vtau_ice : surface ice stress                     (U- & V-points)   [N/m2] 
    209162         ! qsr_ice , qns_ice  : solar & non solar heat flux over ice   (T-point)         [W/m2] 
     
    215168         SELECT CASE( kblk ) 
    216169         CASE( 3 )                                       ! CLIO bulk formulation 
    217             CALL blk_ice_clio( t_su , zalb_ice_cs, zalb_ice_os,                           & 
     170            CALL blk_ice_clio( t_su , zalb_cs, zalb_os, zalb_ice                  & 
    218171               &                      utau_ice   , vtau_ice   , qns_ice   , qsr_ice   ,   & 
    219172               &                      qla_ice    , dqns_ice   , dqla_ice  ,               & 
     
    221174               &                      fr1_i0     , fr2_i0     , cp_ice_msh, jpl  ) 
    222175            !          
     176            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
     177               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
     178 
    223179         CASE( 4 )                                       ! CORE bulk formulation 
    224             ! MV 2014 
    225             ! We must account for cloud fraction in the computation of the albedo 
    226             ! The present ref just uses the clear sky value 
    227             ! The overcast sky value is 0.06 higher, and polar skies are mostly overcast 
    228             ! CORE has no cloud fraction, hence we must prescribe it 
    229             ! Mean summer cloud fraction computed from CLIO = 0.81 
    230             zalb_ice(:,:,:) = 0.19 * zalb_ice_cs(:,:,:) + 0.81 * zalb_ice_os(:,:,:) 
    231             ! Following line, we replace zalb_ice_cs by simply zalb_ice 
    232             CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice   ,               & 
     180            CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice,                  & 
    233181               &                      utau_ice  , vtau_ice  , qns_ice    , qsr_ice   ,   & 
    234182               &                      qla_ice   , dqns_ice  , dqla_ice   ,               & 
    235183               &                      tprecip   , sprecip   ,                            & 
    236184               &                      fr1_i0    , fr2_i0    , cp_ice_msh, jpl  ) 
     185               ! 
     186            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
     187               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
    237188            ! 
    238189         CASE ( 5 ) 
    239             zalb_ice (:,:,:) = 0.5 * ( zalb_ice_cs (:,:,:) +  zalb_ice_os (:,:,:) ) 
    240190             
    241191            CALL sbc_cpl_ice_tau( utau_ice , vtau_ice ) 
    242192 
    243             CALL sbc_cpl_ice_flx( p_frld=ato_i, palbi=zalb_ice, psst=sst_m, pist=tn_ice    ) 
    244  
     193            CALL sbc_cpl_ice_flx( p_frld=ato_i, palbi=zalb_ice, psst=sst_m, pist=t_su    ) 
     194 
     195            IF( nn_limflx == 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
     196               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
    245197            ! Latent heat flux is forced to 0 in coupled : 
    246198            !  it is included in qns (non-solar heat flux) 
    247             qla_ice  (:,:,:) = 0.0e0_wp 
    248             dqla_ice (:,:,:) = 0.0e0_wp 
     199            qla_ice  (:,:,:) = 0._wp 
     200            dqla_ice (:,:,:) = 0._wp 
    249201            ! 
    250202         END SELECT 
    251  
    252          ! Average over all categories 
    253          IF( lk_cpl ) THEN 
    254          IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) THEN 
    255  
    256             z_qns_ice_all  (:,:) = fice_ice_ave ( qns_ice  (:,:,:) ) 
    257             z_qsr_ice_all  (:,:) = fice_ice_ave ( qsr_ice  (:,:,:) ) 
    258             z_dqns_ice_all (:,:) = fice_ice_ave ( dqns_ice (:,:,:) ) 
    259             z_qla_ice_all  (:,:) = fice_ice_ave ( qla_ice  (:,:,:) ) 
    260             z_dqla_ice_all (:,:) = fice_ice_ave ( dqla_ice (:,:,:) ) 
    261  
    262             DO jl = 1, jpl 
    263                dqns_ice (:,:,jl) = z_dqns_ice_all (:,:) 
    264                dqla_ice (:,:,jl) = z_dqla_ice_all (:,:) 
    265             END DO 
    266             ! 
    267             IF ( ln_iceflx_ave ) THEN 
    268                DO jl = 1, jpl 
    269                   qns_ice  (:,:,jl) = z_qns_ice_all  (:,:) 
    270                   qsr_ice  (:,:,jl) = z_qsr_ice_all  (:,:) 
    271                   qla_ice  (:,:,jl) = z_qla_ice_all  (:,:) 
    272                END DO 
    273             END IF 
    274             ! 
    275             IF ( ln_iceflx_linear ) THEN 
    276                DO jl = 1, jpl 
    277                   qns_ice  (:,:,jl) = z_qns_ice_all(:,:) + z_dqns_ice_all(:,:) * (tn_ice(:,:,jl) - ztem_ice_all(:,:)) 
    278                   qla_ice  (:,:,jl) = z_qla_ice_all(:,:) + z_dqla_ice_all(:,:) * (tn_ice(:,:,jl) - ztem_ice_all(:,:)) 
    279                   qsr_ice  (:,:,jl) = (1.0e0_wp-zalb_ice(:,:,jl)) / (1.0e0_wp-zalb_ice_all(:,:)) * z_qsr_ice_all(:,:) 
    280                END DO 
    281             END IF 
    282          END IF 
    283          ENDIF 
     203          
     204         CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice ) 
     205 
    284206         !                                           !----------------------! 
    285207         !                                           ! LIM-3  time-stepping ! 
     
    297219         old_smv_i(:,:,:)   = smv_i(:,:,:)     ! salt content 
    298220         old_oa_i (:,:,:)   = oa_i (:,:,:)     ! areal age content 
    299          old_u_ice(:,:)     = u_ice(:,:) 
    300          old_v_ice(:,:)     = v_ice(:,:) 
    301  
    302          ! trends    !!gm is it truly necessary ??? 
     221         ! 
     222         old_u_ice(:,:) = u_ice(:,:) 
     223         old_v_ice(:,:) = v_ice(:,:) 
     224         !                                           ! intialisation to zero    !!gm is it truly necessary ??? 
    303225         d_a_i_thd  (:,:,:)   = 0._wp   ;   d_a_i_trp  (:,:,:)   = 0._wp 
    304226         d_v_i_thd  (:,:,:)   = 0._wp   ;   d_v_i_trp  (:,:,:)   = 0._wp 
     
    308230         d_smv_i_thd(:,:,:)   = 0._wp   ;   d_smv_i_trp(:,:,:)   = 0._wp 
    309231         d_oa_i_thd (:,:,:)   = 0._wp   ;   d_oa_i_trp (:,:,:)   = 0._wp 
    310          d_u_ice_dyn(:,:)     = 0._wp   ;   d_v_ice_dyn(:,:)     = 0._wp 
    311  
    312          ! salt, heat and mass fluxes 
    313          sfx    (:,:) = 0._wp   ; 
    314          sfx_bri(:,:) = 0._wp   ;   sfx_dyn(:,:) = 0._wp  
    315          sfx_sni(:,:) = 0._wp   ;   sfx_opw(:,:) = 0._wp 
    316          sfx_bog(:,:) = 0._wp   ;   sfx_dyn(:,:) = 0._wp 
    317          sfx_bom(:,:) = 0._wp   ;   sfx_sum(:,:) = 0._wp 
    318          sfx_res(:,:) = 0._wp 
    319  
    320          wfx_snw(:,:) = 0._wp   ;   wfx_ice(:,:) = 0._wp 
    321          wfx_sni(:,:) = 0._wp   ;   wfx_opw(:,:) = 0._wp 
    322          wfx_bog(:,:) = 0._wp   ;   wfx_dyn(:,:) = 0._wp 
    323          wfx_bom(:,:) = 0._wp   ;   wfx_sum(:,:) = 0._wp 
    324          wfx_res(:,:) = 0._wp   ;   wfx_sub(:,:) = 0._wp 
    325          wfx_spr(:,:) = 0._wp   ;    
    326  
    327          hfx_in (:,:) = 0._wp   ;   hfx_out(:,:) = 0._wp 
    328          hfx_thd(:,:) = 0._wp   ;    
    329          hfx_snw(:,:) = 0._wp   ;   hfx_opw(:,:) = 0._wp 
    330          hfx_bog(:,:) = 0._wp   ;   hfx_dyn(:,:) = 0._wp 
    331          hfx_bom(:,:) = 0._wp   ;   hfx_sum(:,:) = 0._wp 
    332          hfx_res(:,:) = 0._wp   ;   hfx_sub(:,:) = 0._wp 
    333          hfx_spr(:,:) = 0._wp   ;   hfx_dif(:,:) = 0._wp  
    334          hfx_err(:,:) = 0._wp   ;   hfx_err_rem(:,:) = 0._wp 
    335  
    336          ! 
    337          fhld  (:,:)    = 0._wp  
    338          fmmflx(:,:)    = 0._wp      
    339          ! part of solar radiation transmitted through the ice 
    340          ftr_ice(:,:,:) = 0._wp 
    341  
    342          ! diags 
    343          diag_trp_vi  (:,:) = 0._wp  ; diag_trp_vs(:,:) = 0._wp  ;  diag_trp_ei(:,:) = 0._wp  ;  diag_trp_es(:,:) = 0._wp 
    344          diag_heat_dhc(:,:) = 0._wp   
    345  
     232         ! 
     233         d_u_ice_dyn(:,:) = 0._wp 
     234         d_v_ice_dyn(:,:) = 0._wp 
     235         ! 
     236         sfx    (:,:) = 0._wp   ;   sfx_thd  (:,:) = 0._wp 
     237         sfx_bri(:,:) = 0._wp   ;   sfx_mec  (:,:) = 0._wp   ;   sfx_res  (:,:) = 0._wp 
     238         fhbri  (:,:) = 0._wp   ;   fheat_mec(:,:) = 0._wp   ;   fheat_res(:,:) = 0._wp 
     239         fhmec  (:,:) = 0._wp   ;    
     240         fmmec  (:,:) = 0._wp 
     241         fmmflx (:,:) = 0._wp      
     242         focea2D(:,:) = 0._wp 
     243         fsup2D (:,:) = 0._wp 
     244 
     245         ! used in limthd.F90 
     246         rdvosif(:,:) = 0._wp   ! variation of ice volume at surface 
     247         rdvobif(:,:) = 0._wp   ! variation of ice volume at bottom 
     248         fdvolif(:,:) = 0._wp   ! total variation of ice volume 
     249         rdvonif(:,:) = 0._wp   ! lateral variation of ice volume 
     250         fstric (:,:) = 0._wp   ! part of solar radiation transmitted through the ice 
     251         ffltbif(:,:) = 0._wp   ! linked with fstric 
     252         qfvbq  (:,:) = 0._wp   ! linked with fstric 
     253         rdm_snw(:,:) = 0._wp   ! variation of snow mass per unit area 
     254         rdm_ice(:,:) = 0._wp   ! variation of ice mass per unit area 
     255         hicifp (:,:) = 0._wp   ! daily thermodynamic ice production.  
     256         ! 
     257         diag_sni_gr(:,:) = 0._wp   ;   diag_lat_gr(:,:) = 0._wp 
     258         diag_bot_gr(:,:) = 0._wp   ;   diag_dyn_gr(:,:) = 0._wp 
     259         diag_bot_me(:,:) = 0._wp   ;   diag_sur_me(:,:) = 0._wp 
     260         diag_res_pr(:,:) = 0._wp   ;   diag_trp_vi(:,:) = 0._wp 
    346261         ! dynamical invariants 
    347262         delta_i(:,:) = 0._wp       ;   divu_i(:,:) = 0._wp       ;   shear_i(:,:) = 0._wp 
     
    394309                          zcoef = rdt_ice /rday           !  Ice natural aging 
    395310                          oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * zcoef 
     311                          CALL lim_var_glo2eqv            ! this CALL is maybe not necessary (Martin) 
    396312         IF( ln_nicep )   CALL lim_prt_state( kt, jiindx, jjindx, 1, ' - ice thermodyn. - ' )   ! control print 
    397313                          CALL lim_itd_th( kt )           !  Remap ice categories, lateral accretion  ! 
     
    409325         !                                           ! Diagnostics and outputs  
    410326         IF (ln_limdiaout) CALL lim_diahsb 
    411  
     327!clem # if ! defined key_iomput 
    412328                          CALL lim_wri( 1  )              ! Ice outputs  
    413  
     329!clem # endif 
    414330         IF( kt == nit000 .AND. ln_rstart )   & 
    415331            &             CALL iom_close( numrir )        ! clem: close input ice restart file 
     
    431347       
    432348!!gm   remark, the ocean-ice stress is not saved in ice diag call above .....  find a solution!!! 
    433       CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs, zalb_ice ) 
    434  
    435       IF( lk_cpl ) THEN 
    436          IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) & 
    437             &    CALL wrk_dealloc( jpi,jpj, ztem_ice_all, zalb_ice_all, z_qsr_ice_all, z_qns_ice_all, z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all) 
    438       ENDIF 
     349 
    439350      ! 
    440351      IF( nn_timing == 1 )  CALL timing_stop('sbc_ice_lim') 
    441352      ! 
    442353   END SUBROUTINE sbc_ice_lim 
    443  
    444  
     354    
     355    
     356      SUBROUTINE ice_lim_flx( ptn_ice, palb_ice, pqns_ice, pqsr_ice,   & 
     357         &                          pdqn_ice, pqla_ice, pdql_ice, k_limflx ) 
     358      !!--------------------------------------------------------------------- 
     359      !!                  ***  ROUTINE sbc_ice_lim  *** 
     360      !!                    
     361      !! ** Purpose :   update the ice surface boundary condition by averaging and / or 
     362     !!                redistributing fluxes on ice categories                    
     363      !! 
     364      !! ** Method  :   average then redistribute  
     365      !! 
     366      !! ** Action  :    
     367      !!--------------------------------------------------------------------- 
     368      INTEGER                   , INTENT(in   ) ::   k_limflx   ! =-1 do nothing; =0 average ;  
     369                                                                ! =1 average and redistribute ; =2 redistribute 
     370      REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   ptn_ice    ! ice surface temperature  
     371      REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   palb_ice   ! ice albedo 
     372      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqns_ice   ! non solar flux 
     373      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqsr_ice   ! net solar flux 
     374      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdqn_ice   ! non solar flux sensitivity 
     375      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqla_ice   ! latent heat flux 
     376      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdql_ice   ! latent heat flux sensitivity 
     377      ! 
     378      INTEGER  ::   jl      ! dummy loop index 
     379      ! 
     380      REAL(wp), POINTER, DIMENSION(:,:) :: zalb_m    ! Mean albedo over all categories 
     381      REAL(wp), POINTER, DIMENSION(:,:) :: ztem_m    ! Mean temperature over all categories 
     382      ! 
     383      REAL(wp), POINTER, DIMENSION(:,:) :: z_qsr_m   ! Mean solar heat flux over all categories 
     384      REAL(wp), POINTER, DIMENSION(:,:) :: z_qns_m   ! Mean non solar heat flux over all categories 
     385      REAL(wp), POINTER, DIMENSION(:,:) :: z_qla_m   ! Mean latent heat flux over all categories 
     386      REAL(wp), POINTER, DIMENSION(:,:) :: z_dqn_m   ! Mean d(qns)/dT over all categories 
     387      REAL(wp), POINTER, DIMENSION(:,:) :: z_dql_m   ! Mean d(qla)/dT over all categories 
     388      !!---------------------------------------------------------------------- 
     389 
     390      IF( nn_timing == 1 )  CALL timing_start('ice_lim_flx') 
     391      ! 
     392      ! 
     393      SELECT CASE( k_limflx )                              !==  averaged on all ice categories  ==! 
     394      CASE( 0 , 1 ) 
     395         CALL wrk_alloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m) 
     396         ! 
     397         z_qns_m(:,:) = fice_ice_ave ( pqns_ice (:,:,:) ) 
     398         z_qsr_m(:,:) = fice_ice_ave ( pqsr_ice (:,:,:) ) 
     399         z_dqn_m(:,:) = fice_ice_ave ( pdqn_ice (:,:,:) ) 
     400         z_qla_m(:,:) = fice_ice_ave ( pqla_ice (:,:,:) ) 
     401         z_dql_m(:,:) = fice_ice_ave ( pdql_ice (:,:,:) ) 
     402         DO jl = 1, jpl 
     403            pdqns_ice(:,:,jl) = z_dqn_m(:,:) 
     404            pdqla_ice(:,:,jl) = z_dql_m(:,:) 
     405         END DO 
     406         ! 
     407         DO jl = 1, jpl 
     408            pqns_ice(:,:,jl) = z_qns_m(:,:) 
     409            pqsr_ice(:,:,jl) = z_qsr_m(:,:) 
     410            pqla_ice(:,:,jl) = z_qla_m(:,:) 
     411         END DO 
     412         ! 
     413         CALL wrk_dealloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m) 
     414      END SELECT 
     415 
     416      SELECT CASE( k_limflx )                              !==  redistribution on all ice categories  ==! 
     417      CASE( 1 , 2 ) 
     418         CALL wrk_alloc( jpi,jpj, zalb_m, ztem_m ) 
     419         ! 
     420         zalb_m(:,:) = fice_ice_ave ( palb_ice (:,:,:) )  
     421         ztem_m(:,:) = fice_ice_ave ( ptn_ice  (:,:,:) )  
     422         DO jl = 1, jpl 
     423            pqns_ice(:,:,jl) = pqns_ice(:,:,jl) + pdqn_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:)) 
     424            pqla_ice(:,:,jl) = pqla_ice(:,:,jl) + pdql_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:)) 
     425            pqsr_ice(:,:,jl) = pqsr_ice(:,:,jl) * ( 1._wp - palb_ice(:,:,jl) ) / ( 1._wp - zalb_m(:,:) )  
     426         END DO 
     427         ! 
     428         CALL wrk_dealloc( jpi,jpj, zalb_m, ztem_m ) 
     429      END SELECT 
     430      ! 
     431      IF( nn_timing == 1 )  CALL timing_stop('ice_lim_flx') 
     432      ! 
     433   END SUBROUTINE ice_lim_flx 
     434    
     435    
    445436   SUBROUTINE lim_ctl( kt ) 
    446437      !!----------------------------------------------------------------------- 
     
    550541!                 WRITE(numout,*) ' sst                  : ', sst_m(ji,jj) 
    551542!                 WRITE(numout,*) ' sss                  : ', sss_m(ji,jj) 
     543!                 WRITE(numout,*) ' s_i_newice           : ', s_i_newice(ji,jj,1:jpl) 
    552544!                 WRITE(numout,*)  
    553545                  inb_alp(ialert_id) = inb_alp(ialert_id) + 1 
     
    606598               !WRITE(numout,*) ' sst       : ', sst_m(ji,jj) 
    607599               !WRITE(numout,*) ' sss       : ', sss_m(ji,jj) 
     600               !WRITE(numout,*) ' qcmif     : ', qcmif(ji,jj) 
     601               !WRITE(numout,*) ' qldif     : ', qldif(ji,jj) 
     602               !WRITE(numout,*) ' qcmif     : ', qcmif(ji,jj) / rdt_ice 
     603               !WRITE(numout,*) ' qldif     : ', qldif(ji,jj) / rdt_ice 
     604               !WRITE(numout,*) ' qfvbq     : ', qfvbq(ji,jj) 
     605               !WRITE(numout,*) ' qdtcn     : ', qdtcn(ji,jj) 
     606               !WRITE(numout,*) ' qfvbq / dt: ', qfvbq(ji,jj) / rdt_ice 
     607               !WRITE(numout,*) ' qdtcn / dt: ', qdtcn(ji,jj) / rdt_ice 
     608               !WRITE(numout,*) ' fdtcn     : ', fdtcn(ji,jj)  
     609               !WRITE(numout,*) ' fhmec     : ', fhmec(ji,jj)  
     610               !WRITE(numout,*) ' fheat_mec : ', fheat_mec(ji,jj)  
     611               !WRITE(numout,*) ' fheat_res : ', fheat_res(ji,jj)  
     612               !WRITE(numout,*) ' fhbri     : ', fhbri(ji,jj)  
    608613               ! 
    609614               !CALL lim_prt_state( kt, ji, jj, 2, '   ') 
     
    674679      !!                n : number of the option 
    675680      !!------------------------------------------------------------------- 
    676       INTEGER         , INTENT(in) ::   kt      ! ocean time step 
     681      INTEGER         , INTENT(in) ::   kt            ! ocean time step 
    677682      INTEGER         , INTENT(in) ::   ki, kj, kn    ! ocean gridpoint indices 
    678683      CHARACTER(len=*), INTENT(in) ::   cd1           ! 
     
    792797               WRITE(numout,*) ' - Heat / FW fluxes ' 
    793798               WRITE(numout,*) '   ~~~~~~~~~~~~~~~~ ' 
    794                WRITE(numout,*) ' - Heat fluxes in and out the ice ***' 
    795                WRITE(numout,*) ' qsr_ini       : ', pfrld(ji,jj) * qsr(ji,jj) + SUM( old_a_i(ji,jj,:) * qsr_ice(ji,jj,:) ) 
    796                WRITE(numout,*) ' qns_ini       : ', pfrld(ji,jj) * qns(ji,jj) + SUM( old_a_i(ji,jj,:) * qns_ice(ji,jj,:) ) 
    797                WRITE(numout,*) 
     799               WRITE(numout,*) ' emp        : ', emp      (ji,jj) 
     800               WRITE(numout,*) ' sfx        : ', sfx      (ji,jj) 
     801               WRITE(numout,*) ' sfx_thd    : ', sfx_thd(ji,jj) 
     802               WRITE(numout,*) ' sfx_bri    : ', sfx_bri  (ji,jj) 
     803               WRITE(numout,*) ' sfx_mec    : ', sfx_mec  (ji,jj) 
     804               WRITE(numout,*) ' sfx_res    : ', sfx_res(ji,jj) 
     805               WRITE(numout,*) ' fmmec      : ', fmmec    (ji,jj) 
     806               WRITE(numout,*) ' fhmec      : ', fhmec    (ji,jj) 
     807               WRITE(numout,*) ' fhbri      : ', fhbri    (ji,jj) 
     808               WRITE(numout,*) ' fheat_mec  : ', fheat_mec(ji,jj) 
    798809               WRITE(numout,*)  
    799810               WRITE(numout,*) ' sst        : ', sst_m(ji,jj)   
     
    825836               WRITE(numout,*) ' qsr       : ', qsr(ji,jj) 
    826837               WRITE(numout,*) ' qns       : ', qns(ji,jj) 
    827                WRITE(numout,*) 
    828                WRITE(numout,*) ' hfx_mass     : ', hfx_thd(ji,jj) + hfx_dyn(ji,jj) + hfx_snw(ji,jj) + hfx_res(ji,jj) 
    829                WRITE(numout,*) ' hfx_in       : ', hfx_in(ji,jj) 
    830                WRITE(numout,*) ' hfx_out      : ', hfx_out(ji,jj) 
    831                WRITE(numout,*) ' dhc          : ', diag_heat_dhc(ji,jj)               
    832                WRITE(numout,*) 
    833                WRITE(numout,*) ' hfx_dyn      : ', hfx_dyn(ji,jj) 
    834                WRITE(numout,*) ' hfx_thd      : ', hfx_thd(ji,jj) 
    835                WRITE(numout,*) ' hfx_res      : ', hfx_res(ji,jj) 
    836                WRITE(numout,*) ' fhtur        : ', fhtur(ji,jj)  
    837                WRITE(numout,*) ' qlead        : ', qlead(ji,jj) * r1_rdtice 
     838               WRITE(numout,*) ' fdtcn     : ', fdtcn(ji,jj) 
     839               WRITE(numout,*) ' qcmif     : ', qcmif(ji,jj) * r1_rdtice 
     840               WRITE(numout,*) ' qldif     : ', qldif(ji,jj) * r1_rdtice 
    838841               WRITE(numout,*) 
    839842               WRITE(numout,*) ' - Salt fluxes at bottom interface ***' 
    840843               WRITE(numout,*) ' emp       : ', emp    (ji,jj) 
     844               WRITE(numout,*) ' sfx_bri   : ', sfx_bri(ji,jj) 
    841845               WRITE(numout,*) ' sfx       : ', sfx    (ji,jj) 
    842846               WRITE(numout,*) ' sfx_res   : ', sfx_res(ji,jj) 
    843                WRITE(numout,*) ' sfx_bri   : ', sfx_bri(ji,jj) 
    844                WRITE(numout,*) ' sfx_dyn   : ', sfx_dyn(ji,jj) 
     847               WRITE(numout,*) ' sfx_mec   : ', sfx_mec(ji,jj) 
     848               WRITE(numout,*) ' - Heat fluxes at bottom interface ***' 
     849               WRITE(numout,*) ' fheat_res : ', fheat_res(ji,jj) 
    845850               WRITE(numout,*) 
    846851               WRITE(numout,*) ' - Momentum fluxes ' 
    847852               WRITE(numout,*) ' utau      : ', utau(ji,jj)  
    848853               WRITE(numout,*) ' vtau      : ', vtau(ji,jj) 
    849             ENDIF  
     854            ENDIF 
    850855            WRITE(numout,*) ' ' 
    851856            ! 
    852857         END DO 
    853858      END DO 
    854  
     859      ! 
    855860   END SUBROUTINE lim_prt_state 
     861    
     862      
     863   FUNCTION fice_cell_ave ( ptab ) 
     864      !!-------------------------------------------------------------------------- 
     865      !! * Compute average over categories, for grid cell (ice covered and free ocean) 
     866      !!-------------------------------------------------------------------------- 
     867      REAL (wp), DIMENSION (jpi,jpj) :: fice_cell_ave 
     868      REAL (wp), DIMENSION (jpi,jpj,jpl), INTENT (in) :: ptab 
     869      INTEGER :: jl ! Dummy loop index 
     870       
     871      fice_cell_ave (:,:) = 0.0_wp 
     872       
     873      DO jl = 1, jpl 
     874         fice_cell_ave (:,:) = fice_cell_ave (:,:) & 
     875            &                  + a_i (:,:,jl) * ptab (:,:,jl) 
     876      END DO 
     877       
     878   END FUNCTION fice_cell_ave 
     879    
     880    
     881   FUNCTION fice_ice_ave ( ptab ) 
     882      !!-------------------------------------------------------------------------- 
     883      !! * Compute average over categories, for ice covered part of grid cell 
     884      !!-------------------------------------------------------------------------- 
     885      REAL (kind=wp), DIMENSION (jpi,jpj) :: fice_ice_ave 
     886      REAL (kind=wp), DIMENSION (jpi,jpj,jpl), INTENT(in) :: ptab 
     887 
     888      fice_ice_ave (:,:) = 0.0_wp 
     889      WHERE ( at_i (:,:) .GT. 0.0_wp ) fice_ice_ave (:,:) = fice_cell_ave ( ptab (:,:,:)) / at_i (:,:) 
     890 
     891   END FUNCTION fice_ice_ave 
     892 
    856893 
    857894#else 
  • branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r4624 r4730  
    8484      NAMELIST/namsbc/ nn_fsbc   , ln_ana    , ln_flx,  ln_blk_clio, ln_blk_core, ln_cpl,   & 
    8585         &             ln_blk_mfs, ln_apr_dyn, nn_ice,  nn_ice_embd, ln_dm2dc   , ln_rnf,   & 
    86          &             ln_ssr    , nn_fwb    , ln_cdgw , ln_wave , ln_sdw, nn_lsm, cn_iceflx 
     86         &             ln_ssr    , nn_fwb    , ln_cdgw , ln_wave , ln_sdw, nn_lsm, nn_iceflx 
    8787      INTEGER  ::   ios 
    8888      !!---------------------------------------------------------------------- 
     
    124124         WRITE(numout,*) '              MFS  bulk  formulation                     ln_blk_mfs  = ', ln_blk_mfs 
    125125         WRITE(numout,*) '              coupled    formulation (T if key_sbc_cpl)  ln_cpl      = ', ln_cpl 
    126          WRITE(numout,*) '              Flux handling over ice categories          cn_iceflx   = ', TRIM (cn_iceflx) 
     126         WRITE(numout,*) '              Multicategory heat flux formulation (LIM3) nn_limflx   = ', nn_limflx 
    127127         WRITE(numout,*) '           Misc. options of sbc : ' 
    128128         WRITE(numout,*) '              Patm gradient added in ocean & ice Eqs.    ln_apr_dyn  = ', ln_apr_dyn 
     
    137137      ENDIF 
    138138 
    139       !   Flux handling over ice categories 
    140 #if defined key_coupled  
    141       SELECT CASE ( TRIM (cn_iceflx)) 
    142       CASE ('ave') 
    143          ln_iceflx_ave    = .TRUE. 
    144          ln_iceflx_linear = .FALSE. 
    145       CASE ('linear') 
    146          ln_iceflx_ave    = .FALSE. 
    147          ln_iceflx_linear = .TRUE. 
    148       CASE default 
    149          ln_iceflx_ave    = .FALSE. 
    150          ln_iceflx_linear = .FALSE. 
     139      ! LIM3 Multi-category heat flux formulation 
     140      SELECT CASE ( nn_limflx) 
     141      CASE ( -1 ) 
     142         IF(lwp) THEN WRITE(numout,*) '              Use of per-category fluxes (nn_limflx = -1) ' 
     143      CASE ( 0  ) 
     144         IF(lwp) THEN WRITE(numout,*) '              Average per-category fluxes (nn_limflx = 0) '  
     145      CASE ( 1  ) 
     146         IF(lwp) THEN WRITE(numout,*) '              Average then redistribute per-category fluxes (nn_limflx = 1) ' 
     147      CASE ( 2  ) 
     148         IF(lwp) THEN WRITE(numout,*) '              Redistribute a single flux over categories (nn_limflx = 2) ' 
    151149      END SELECT 
    152       IF(lwp) WRITE(numout,*) '              Fluxes averaged over all ice categories         ln_iceflx_ave    = ', ln_iceflx_ave 
    153       IF(lwp) WRITE(numout,*) '              Fluxes distributed linearly over ice categories ln_iceflx_linear = ', ln_iceflx_linear 
    154 #endif 
    155150      ! 
    156151#if defined key_top && ! defined key_offline 
     
    206201      IF( ( nn_ice == 3 .OR. nn_ice == 4 ) .AND. nn_ice_embd == 0 )   & 
    207202         &   CALL ctl_stop( 'LIM3 and CICE sea-ice models require nn_ice_embd = 1 or 2' ) 
    208 #if defined key_coupled 
    209       IF( ln_iceflx_ave .AND. ln_iceflx_linear ) & 
    210          &   CALL ctl_stop( ' ln_iceflx_ave and ln_iceflx_linear options are not compatible' ) 
    211       IF( ( nn_ice ==3 .AND. lk_cpl) .AND. .NOT. ( ln_iceflx_ave .OR. ln_iceflx_linear ) ) & 
    212          &   CALL ctl_stop( ' With lim3 coupled, either ln_iceflx_ave or ln_iceflx_linear must be set to .TRUE.' ) 
    213 #endif       
     203      IF( ( nn_ice /= 3 ) .AND. ( nn_limflx >= 0 ) )   & 
     204         &   WRITE(numout,*) 'The nn_limflx>=0 option has no effect if sea ice model is not LIM3' 
     205      IF( ( nn_ice == 3 ) .AND. ( lk_cpl ) .AND. ( ( nn_limflx == -1 ) .OR. ( nn_limflx == 1 ) ) )   & 
     206         &   CALL ctl_stop( 'The chosen nn_limflx for LIM3 in coupled mode must be 0 or 2' ) 
     207      IF( ( nn_ice == 3 ) .AND. ( .NOT. lk_cpl ) .AND. ( nn_limflx == 2 ) )   & 
     208         &   CALL ctl_stop( 'The chosen nn_limflx for LIM3 in forced mode cannot be 2' ) 
     209 
    214210      IF( ln_dm2dc )   nday_qsr = -1   ! initialisation flag 
    215211 
Note: See TracChangeset for help on using the changeset viewer.