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 4905 for branches/2014/dev_CNRS_2014 – NEMO

Ignore:
Timestamp:
2014-11-27T17:32:55+01:00 (9 years ago)
Author:
cetlod
Message:

2014/dev_CNRS_2014 :minor correction on sbcice_lim.F90, see ticket #1415

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_CNRS_2014/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90

    r4902 r4905  
    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  
     
    8080   !!---------------------------------------------------------------------- 
    8181CONTAINS 
    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 
    11182 
    11283   !!====================================================================== 
     
    133104      !!--------------------------------------------------------------------- 
    134105      INTEGER, INTENT(in) ::   kt      ! ocean time step 
    135       INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE) 
     106      INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE, =5 COUPLED) 
    136107      !! 
    137       INTEGER  ::   ji, jj, jl, jk      ! dummy loop index 
     108      INTEGER  ::   jl      ! dummy loop index 
    138109      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   ! 
     110      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_os, zalb_cs  ! ice albedo under overcast/clear sky 
     111      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice          ! mean ice albedo (for coupled) 
    152112      !!---------------------------------------------------------------------- 
    153113 
    154       !- O.M. : why do we allocate all these arrays even when MOD( kt-1, nn_fsbc ) /= 0 ????? 
    155  
    156114      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,   & 
    163             &                             z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all) 
    164       ENDIF 
    165115 
    166116      IF( kt == nit000 ) THEN 
     
    186136         v_oce(:,:) = ssv_m(:,:)                     ! (C-grid dynamics :  U- & V-points as the ocean) 
    187137         ! 
    188          t_bo(:,:) = ( eos_fzp( sss_m ) +  rt0 ) * tmask(:,:,1) + rt0 * ( 1. - tmask(:,:,1) )   ! masked sea surface freezing temperature [Kelvin] 
    189          !                                                                                      ! (set to rt0 over land) 
    190          CALL albedo_ice( t_su, ht_i, ht_s, zalb_ice_cs, zalb_ice_os )  ! ... ice albedo 
    191  
     138         t_bo(:,:) = ( eos_fzp( sss_m ) +  rt0 ) * tmask(:,:,1) + rt0 * ( 1. - tmask(:,:,1) )  ! masked sea surface freezing temperature [Kelvin] 
     139         !                                                                                  ! (set to rt0 over land) 
     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( jp_core , jp_cpl )   ! 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] 
     
    214167         ! 
    215168         SELECT CASE( kblk ) 
    216          CASE( 3 )                                       ! CLIO bulk formulation 
    217             CALL blk_ice_clio( t_su , zalb_ice_cs, zalb_ice_os,                           & 
     169         CASE( jp_clio )                                       ! CLIO bulk formulation 
     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            !          
    223          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 
     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 
     179         CASE( jp_core )                                       ! CORE bulk formulation 
    232180            CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice   ,               & 
    233181               &                      utau_ice  , vtau_ice  , qns_ice    , qsr_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            ! 
    238          CASE ( 5 ) 
    239             zalb_ice (:,:,:) = 0.5 * ( zalb_ice_cs (:,:,:) +  zalb_ice_os (:,:,:) ) 
     189         CASE ( jp_cpl ) 
    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  
    245             ! Latent heat flux is forced to 0 in coupled : 
    246             !  it is included in qns (non-solar heat flux) 
    247             qla_ice  (:,:,:) = 0.0e0_wp 
    248             dqla_ice (:,:,:) = 0.0e0_wp 
     193            ! MV -> seb  
     194!           CALL sbc_cpl_ice_flx( p_frld=ato_i, palbi=zalb_ice, psst=sst_m, pist=t_su    ) 
     195 
     196!           IF( nn_limflx == 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
     197!              &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
     198!           ! Latent heat flux is forced to 0 in coupled : 
     199!           !  it is included in qns (non-solar heat flux) 
     200!           qla_ice  (:,:,:) = 0._wp 
     201!           dqla_ice (:,:,:) = 0._wp 
     202            ! END MV -> seb 
    249203            ! 
    250204         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 
     205          
    284206         !                                           !----------------------! 
    285207         !                                           ! LIM-3  time-stepping ! 
     
    388310                          pfrld(:,:)   = 1._wp - at_i(:,:) 
    389311                          phicif(:,:)  = vt_i(:,:) 
     312 
     313                          ! MV -> seb 
     314                          SELECT CASE( kblk ) 
     315                             CASE ( jp_cpl ) 
     316                             CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su    ) 
     317                             IF( nn_limflx == 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
     318                          &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
     319                           ! Latent heat flux is forced to 0 in coupled : 
     320                           !  it is included in qns (non-solar heat flux) 
     321                             qla_ice  (:,:,:) = 0._wp 
     322                             dqla_ice (:,:,:) = 0._wp 
     323                          END SELECT 
     324                          ! END MV -> seb 
    390325                          ! 
    391326                          CALL lim_var_bv                 ! bulk brine volume (diag) 
     
    419354         IF( ln_nicep )   CALL lim_ctl( kt )              ! alerts in case of model crash 
    420355         ! 
     356         CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice ) 
     357         ! 
    421358      ENDIF                                    ! End sea-ice time step only 
    422359 
     
    428365      !                                                   ! otherwise the atm.-ocean stresses are used everywhere 
    429366      IF( ln_limdyn )     CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) )  ! using before instantaneous surf. currents 
    430        
    431367!!gm   remark, the ocean-ice stress is not saved in ice diag call above .....  find a solution!!! 
    432       CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs, zalb_ice ) 
    433  
    434       IF( lk_cpl ) THEN 
    435          IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) & 
    436             &    CALL wrk_dealloc( jpi, jpj, ztem_ice_all , zalb_ice_all , z_qsr_ice_all, z_qns_ice_all,   & 
    437             &                                z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all) 
    438       ENDIF 
     368 
    439369      ! 
    440370      IF( nn_timing == 1 )  CALL timing_stop('sbc_ice_lim') 
    441371      ! 
    442372   END SUBROUTINE sbc_ice_lim 
    443  
    444  
     373    
     374    
     375      SUBROUTINE ice_lim_flx( ptn_ice, palb_ice, pqns_ice, pqsr_ice,   & 
     376         &                          pdqn_ice, pqla_ice, pdql_ice, k_limflx ) 
     377      !!--------------------------------------------------------------------- 
     378      !!                  ***  ROUTINE sbc_ice_lim  *** 
     379      !!                    
     380      !! ** Purpose :   update the ice surface boundary condition by averaging and / or 
     381      !!                redistributing fluxes on ice categories                    
     382      !! 
     383      !! ** Method  :   average then redistribute  
     384      !! 
     385      !! ** Action  :    
     386      !!--------------------------------------------------------------------- 
     387      INTEGER                   , INTENT(in   ) ::   k_limflx   ! =-1 do nothing; =0 average ;  
     388                                                                ! =1 average and redistribute ; =2 redistribute 
     389      REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   ptn_ice    ! ice surface temperature  
     390      REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   palb_ice   ! ice albedo 
     391      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqns_ice   ! non solar flux 
     392      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqsr_ice   ! net solar flux 
     393      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdqn_ice   ! non solar flux sensitivity 
     394      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqla_ice   ! latent heat flux 
     395      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdql_ice   ! latent heat flux sensitivity 
     396      ! 
     397      INTEGER  ::   jl      ! dummy loop index 
     398      ! 
     399      REAL(wp), POINTER, DIMENSION(:,:) :: zalb_m    ! Mean albedo over all categories 
     400      REAL(wp), POINTER, DIMENSION(:,:) :: ztem_m    ! Mean temperature over all categories 
     401      ! 
     402      REAL(wp), POINTER, DIMENSION(:,:) :: z_qsr_m   ! Mean solar heat flux over all categories 
     403      REAL(wp), POINTER, DIMENSION(:,:) :: z_qns_m   ! Mean non solar heat flux over all categories 
     404      REAL(wp), POINTER, DIMENSION(:,:) :: z_qla_m   ! Mean latent heat flux over all categories 
     405      REAL(wp), POINTER, DIMENSION(:,:) :: z_dqn_m   ! Mean d(qns)/dT over all categories 
     406      REAL(wp), POINTER, DIMENSION(:,:) :: z_dql_m   ! Mean d(qla)/dT over all categories 
     407      !!---------------------------------------------------------------------- 
     408 
     409      IF( nn_timing == 1 )  CALL timing_start('ice_lim_flx') 
     410      ! 
     411      ! 
     412      SELECT CASE( k_limflx )                              !==  averaged on all ice categories  ==! 
     413      CASE( 0 , 1 ) 
     414         CALL wrk_alloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m) 
     415         ! 
     416         z_qns_m(:,:) = fice_ice_ave ( pqns_ice (:,:,:) ) 
     417         z_qsr_m(:,:) = fice_ice_ave ( pqsr_ice (:,:,:) ) 
     418         z_dqn_m(:,:) = fice_ice_ave ( pdqn_ice (:,:,:) ) 
     419         z_qla_m(:,:) = fice_ice_ave ( pqla_ice (:,:,:) ) 
     420         z_dql_m(:,:) = fice_ice_ave ( pdql_ice (:,:,:) ) 
     421         DO jl = 1, jpl 
     422            pdqn_ice(:,:,jl) = z_dqn_m(:,:) 
     423            pdql_ice(:,:,jl) = z_dql_m(:,:) 
     424         END DO 
     425         ! 
     426         DO jl = 1, jpl 
     427            pqns_ice(:,:,jl) = z_qns_m(:,:) 
     428            pqsr_ice(:,:,jl) = z_qsr_m(:,:) 
     429            pqla_ice(:,:,jl) = z_qla_m(:,:) 
     430         END DO 
     431         ! 
     432         CALL wrk_dealloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m) 
     433      END SELECT 
     434 
     435      SELECT CASE( k_limflx )                              !==  redistribution on all ice categories  ==! 
     436      CASE( 1 , 2 ) 
     437         CALL wrk_alloc( jpi,jpj, zalb_m, ztem_m ) 
     438         ! 
     439         zalb_m(:,:) = fice_ice_ave ( palb_ice (:,:,:) )  
     440         ztem_m(:,:) = fice_ice_ave ( ptn_ice  (:,:,:) )  
     441         DO jl = 1, jpl 
     442            pqns_ice(:,:,jl) = pqns_ice(:,:,jl) + pdqn_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:)) 
     443            pqla_ice(:,:,jl) = pqla_ice(:,:,jl) + pdql_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:)) 
     444            pqsr_ice(:,:,jl) = pqsr_ice(:,:,jl) * ( 1._wp - palb_ice(:,:,jl) ) / ( 1._wp - zalb_m(:,:) )  
     445         END DO 
     446         ! 
     447         CALL wrk_dealloc( jpi,jpj, zalb_m, ztem_m ) 
     448      END SELECT 
     449      ! 
     450      IF( nn_timing == 1 )  CALL timing_stop('ice_lim_flx') 
     451      ! 
     452   END SUBROUTINE ice_lim_flx 
     453    
     454    
    445455   SUBROUTINE lim_ctl( kt ) 
    446456      !!----------------------------------------------------------------------- 
     
    854864      ! 
    855865   END SUBROUTINE lim_prt_state 
     866    
     867      
     868   FUNCTION fice_cell_ave ( ptab ) 
     869      !!-------------------------------------------------------------------------- 
     870      !! * Compute average over categories, for grid cell (ice covered and free ocean) 
     871      !!-------------------------------------------------------------------------- 
     872      REAL (wp), DIMENSION (jpi,jpj) :: fice_cell_ave 
     873      REAL (wp), DIMENSION (jpi,jpj,jpl), INTENT (in) :: ptab 
     874      INTEGER :: jl ! Dummy loop index 
     875       
     876      fice_cell_ave (:,:) = 0.0_wp 
     877       
     878      DO jl = 1, jpl 
     879         fice_cell_ave (:,:) = fice_cell_ave (:,:) & 
     880            &                  + a_i (:,:,jl) * ptab (:,:,jl) 
     881      END DO 
     882       
     883   END FUNCTION fice_cell_ave 
     884    
     885    
     886   FUNCTION fice_ice_ave ( ptab ) 
     887      !!-------------------------------------------------------------------------- 
     888      !! * Compute average over categories, for ice covered part of grid cell 
     889      !!-------------------------------------------------------------------------- 
     890      REAL (kind=wp), DIMENSION (jpi,jpj) :: fice_ice_ave 
     891      REAL (kind=wp), DIMENSION (jpi,jpj,jpl), INTENT(in) :: ptab 
     892 
     893      fice_ice_ave (:,:) = 0.0_wp 
     894      WHERE ( at_i (:,:) .GT. 0.0_wp ) fice_ice_ave (:,:) = fice_cell_ave ( ptab (:,:,:)) / at_i (:,:) 
     895 
     896   END FUNCTION fice_ice_ave 
     897 
    856898 
    857899#else 
Note: See TracChangeset for help on using the changeset viewer.