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 5621 for branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90 – NEMO

Ignore:
Timestamp:
2015-07-21T13:25:36+02:00 (9 years ago)
Author:
mathiot
Message:

UKMO_ISF: upgrade to NEMO_3.6_STABLE (r5554)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90

    r5146 r5621  
    3737   USE limdyn          ! Ice dynamics 
    3838   USE limtrp          ! Ice transport 
     39   USE limhdf          ! Ice horizontal diffusion 
    3940   USE limthd          ! Ice thermodynamics 
    4041   USE limitd_me       ! Mechanics on ice thickness distribution 
     
    110111      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_os, zalb_cs  ! ice albedo under overcast/clear sky 
    111112      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice          ! mean ice albedo (for coupled) 
     113      REAL(wp), POINTER, DIMENSION(:,:  )   ::   zutau_ice, zvtau_ice  
    112114      !!---------------------------------------------------------------------- 
    113115 
    114116      IF( nn_timing == 1 )  CALL timing_start('sbc_ice_lim') 
    115117 
    116       IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !  Ice time-step only 
    117          !-----------------------!                                            
    118          ! --- Bulk Formulae --- !                                            
    119          !-----------------------! 
    120          u_oce(:,:) = ssu_m(:,:) * umask(:,:,1)      ! mean surface ocean current at ice velocity point 
    121          v_oce(:,:) = ssv_m(:,:) * vmask(:,:,1)      ! (C-grid dynamics :  U- & V-points as the ocean) 
     118      !-----------------------! 
     119      ! --- Ice time step --- ! 
     120      !-----------------------! 
     121      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
     122 
     123         ! mean surface ocean current at ice velocity point (C-grid dynamics :  U- & V-points as the ocean) 
     124         u_oce(:,:) = ssu_m(:,:) * umask(:,:,1) 
     125         v_oce(:,:) = ssv_m(:,:) * vmask(:,:,1) 
    122126          
    123127         ! masked sea surface freezing temperature [Kelvin] (set to rt0 over land) 
    124          t_bo(:,:) = ( eos_fzp( sss_m ) + rt0 ) * tmask(:,:,1) + rt0 * ( 1._wp - tmask(:,:,1) )   
    125          !                                                                                       
    126          ! Ice albedo 
    127          CALL wrk_alloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice ) 
    128          CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos 
    129  
    130          ! CORE and COUPLED bulk formulations 
    131          SELECT CASE( kblk ) 
    132          CASE( jp_core , jp_cpl ) 
    133  
    134             ! albedo depends on cloud fraction because of non-linear spectral effects 
    135             zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
    136             ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo  
    137             ! (zalb_ice) is computed within the bulk routine 
    138              
    139          END SELECT 
    140           
     128         CALL eos_fzp( sss_m(:,:) , t_bo(:,:) ) 
     129         t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1) + rt0 * ( 1._wp - tmask(:,:,1) ) 
     130           
    141131         ! Mask sea ice surface temperature (set to rt0 over land) 
    142132         DO jl = 1, jpl 
    143133            t_su(:,:,jl) = t_su(:,:,jl) * tmask(:,:,1) + rt0 * ( 1._wp - tmask(:,:,1) ) 
    144          END DO 
    145       
    146          ! Bulk formulae  - provides the following fields: 
    147          ! utau_ice, vtau_ice : surface ice stress                     (U- & V-points)   [N/m2] 
     134         END DO      
     135         ! 
     136         !------------------------------------------------!                                            
     137         ! --- Dynamical coupling with the atmosphere --- !                                            
     138         !------------------------------------------------! 
     139         ! It provides the following fields: 
     140         ! utau_ice, vtau_ice : surface ice stress (U- & V-points)   [N/m2] 
     141         !----------------------------------------------------------------- 
     142         SELECT CASE( kblk ) 
     143         CASE( jp_clio    )   ;   CALL blk_ice_clio_tau                         ! CLIO bulk formulation             
     144         CASE( jp_core    )   ;   CALL blk_ice_core_tau                         ! CORE bulk formulation 
     145         CASE( jp_purecpl )   ;   CALL sbc_cpl_ice_tau( utau_ice , vtau_ice )   ! Coupled   formulation 
     146         END SELECT 
     147          
     148         IF( ln_mixcpl) THEN   ! Case of a mixed Bulk/Coupled formulation 
     149            CALL wrk_alloc( jpi,jpj    , zutau_ice, zvtau_ice) 
     150            CALL sbc_cpl_ice_tau( zutau_ice , zvtau_ice ) 
     151            utau_ice(:,:) = utau_ice(:,:) * xcplmask(:,:,0) + zutau_ice(:,:) * ( 1. - xcplmask(:,:,0) ) 
     152            vtau_ice(:,:) = vtau_ice(:,:) * xcplmask(:,:,0) + zvtau_ice(:,:) * ( 1. - xcplmask(:,:,0) ) 
     153            CALL wrk_dealloc( jpi,jpj  , zutau_ice, zvtau_ice) 
     154         ENDIF 
     155 
     156         !-------------------------------------------------------! 
     157         ! --- ice dynamics and transport (except in 1D case) ---! 
     158         !-------------------------------------------------------! 
     159         numit = numit + nn_fsbc                  ! Ice model time step 
     160         !                                                    
     161         CALL sbc_lim_bef                         ! Store previous ice values 
     162         CALL sbc_lim_diag0                       ! set diag of mass, heat and salt fluxes to 0 
     163         CALL lim_rst_opn( kt )                   ! Open Ice restart file 
     164         ! 
     165         IF( .NOT. lk_c1d ) THEN 
     166            ! 
     167            CALL lim_dyn( kt )                    ! Ice dynamics    ( rheology/dynamics )    
     168            ! 
     169            CALL lim_trp( kt )                    ! Ice transport   ( Advection/diffusion ) 
     170            ! 
     171            IF( nn_monocat /= 2 ) CALL lim_itd_me ! Mechanical redistribution ! (ridging/rafting) 
     172            ! 
     173#if defined key_bdy 
     174            CALL bdy_ice_lim( kt )                ! bdy ice thermo  
     175            IF( ln_icectl )       CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermo bdy - ' ) 
     176#endif 
     177            ! 
     178            CALL lim_update1( kt )                ! Corrections 
     179            ! 
     180         ENDIF 
     181          
     182         ! previous lead fraction and ice volume for flux calculations 
     183         CALL sbc_lim_bef                         
     184         CALL lim_var_glo2eqv                     ! ht_i and ht_s for ice albedo calculation 
     185         CALL lim_var_agg(1)                      ! at_i for coupling (via pfrld)  
     186         pfrld(:,:)   = 1._wp - at_i(:,:) 
     187         phicif(:,:)  = vt_i(:,:) 
     188          
     189         !------------------------------------------------------!                                            
     190         ! --- Thermodynamical coupling with the atmosphere --- !                                            
     191         !------------------------------------------------------! 
     192         ! It provides the following fields: 
    148193         ! qsr_ice , qns_ice  : solar & non solar heat flux over ice   (T-point)         [W/m2] 
    149194         ! qla_ice            : latent heat flux over ice              (T-point)         [W/m2] 
     
    151196         ! tprecip , sprecip  : total & solid precipitation            (T-point)         [Kg/m2/s] 
    152197         ! fr1_i0  , fr2_i0   : 1sr & 2nd fraction of qsr penetration in ice             [%] 
    153          ! 
     198         !---------------------------------------------------------------------------------------- 
     199         CALL wrk_alloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice ) 
     200         CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos 
     201 
    154202         SELECT CASE( kblk ) 
    155203         CASE( jp_clio )                                       ! CLIO bulk formulation 
    156             CALL blk_ice_clio( t_su , zalb_cs    , zalb_os    , zalb_ice  ,               & 
    157                &                      utau_ice   , vtau_ice   , qns_ice   , qsr_ice   ,   & 
    158                &                      qla_ice    , dqns_ice   , dqla_ice  ,               & 
    159                &                      tprecip    , sprecip    ,                           & 
    160                &                      fr1_i0     , fr2_i0     , cp_ice_msh, jpl  ) 
    161             !          
    162             IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
    163                &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
    164  
     204            ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo  
     205            ! (zalb_ice) is computed within the bulk routine 
     206            CALL blk_ice_clio_flx( t_su, zalb_cs, zalb_os, zalb_ice ) 
     207            IF( ln_mixcpl      ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su ) 
     208            IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
    165209         CASE( jp_core )                                       ! CORE bulk formulation 
    166             CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice   ,               & 
    167                &                      utau_ice  , vtau_ice  , qns_ice    , qsr_ice   ,   & 
    168                &                      qla_ice   , dqns_ice  , dqla_ice   ,               & 
    169                &                      tprecip   , sprecip   ,                            & 
    170                &                      fr1_i0    , fr2_i0    , cp_ice_msh, jpl ) 
    171                ! 
    172             IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
    173                &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
    174             ! 
    175          CASE ( jp_cpl ) 
    176              
    177             CALL sbc_cpl_ice_tau( utau_ice , vtau_ice ) 
    178  
     210            ! albedo depends on cloud fraction because of non-linear spectral effects 
     211            zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     212            CALL blk_ice_core_flx( t_su, zalb_ice ) 
     213            IF( ln_mixcpl      ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su ) 
     214            IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
     215         CASE ( jp_purecpl ) 
     216            ! albedo depends on cloud fraction because of non-linear spectral effects 
     217            zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     218                                 CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su ) 
     219            ! clem: evap_ice is forced to 0 in coupled mode for now  
     220            !       but it needs to be changed (along with modif in limthd_dh) once heat flux from evap will be avail. from atm. models 
     221            evap_ice  (:,:,:) = 0._wp   ;   devap_ice (:,:,:) = 0._wp 
     222            IF( nn_limflx == 2 ) CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
    179223         END SELECT 
    180           
    181          !------------------------------! 
    182          ! --- LIM-3 main time-step --- ! 
    183          !------------------------------! 
    184          numit = numit + nn_fsbc                     ! Ice model time step 
    185          !                                                    
    186          CALL sbc_lim_update                ! Store previous ice values 
    187  
    188          CALL sbc_lim_diag0                 ! set diag of mass, heat and salt fluxes to 0 
    189           
    190          CALL lim_rst_opn( kt )             ! Open Ice restart file 
    191          ! 
    192          ! ---------------------------------------------- 
    193          ! ice dynamics and transport (except in 1D case) 
    194          ! ---------------------------------------------- 
    195          IF( .NOT. lk_c1d ) THEN 
    196              
    197             CALL lim_dyn( kt )              ! Ice dynamics    ( rheology/dynamics ) 
    198              
    199             CALL lim_trp( kt )              ! Ice transport   ( Advection/diffusion ) 
    200              
    201             IF( nn_monocat /= 2 ) CALL lim_itd_me  ! Mechanical redistribution ! (ridging/rafting) 
    202  
    203 #if defined key_bdy 
    204             CALL lim_var_glo2eqv 
    205             CALL bdy_ice_lim( kt )         ! bdy ice thermo  
    206             CALL lim_var_zapsmall 
    207             CALL lim_var_agg(1) 
    208             IF( ln_icectl )   CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermo bdy - ' ) 
    209 #endif 
    210             CALL lim_update1( kt ) 
    211              
    212          ENDIF 
    213           
    214          CALL sbc_lim_update                ! Store previous ice values 
    215   
    216          ! ---------------------------------------------- 
    217          ! ice thermodynamics 
    218          ! ---------------------------------------------- 
    219          CALL lim_var_glo2eqv 
    220          CALL lim_var_agg(1) 
    221           
    222          ! previous lead fraction and ice volume for flux calculations 
    223          pfrld(:,:)   = 1._wp - at_i(:,:) 
    224          phicif(:,:)  = vt_i(:,:) 
    225           
    226          SELECT CASE( kblk ) 
    227          CASE ( jp_cpl ) 
    228             CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su    ) 
    229             IF( nn_limflx == 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   & 
    230                &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx ) 
    231             ! Latent heat flux is forced to 0 in coupled: it is included in qns (non-solar heat flux) 
    232             qla_ice  (:,:,:) = 0._wp 
    233             dqla_ice (:,:,:) = 0._wp 
    234          END SELECT 
    235          ! 
    236          CALL lim_thd( kt )                         ! Ice thermodynamics  
    237           
     224         CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice ) 
     225 
     226         !----------------------------! 
     227         ! --- ice thermodynamics --- ! 
     228         !----------------------------! 
     229         CALL lim_thd( kt )                         ! Ice thermodynamics       
     230         ! 
    238231         CALL lim_update2( kt )                     ! Corrections 
    239232         ! 
     
    241234         ! 
    242235         IF(ln_limdiaout) CALL lim_diahsb           ! Diagnostics and outputs  
    243           
     236         ! 
    244237         CALL lim_wri( 1 )                          ! Ice outputs  
    245           
     238         ! 
    246239         IF( kt == nit000 .AND. ln_rstart )   & 
    247240            &             CALL iom_close( numrir )  ! close input ice restart file 
    248241         ! 
    249242         IF( lrst_ice )   CALL lim_rst_write( kt )  ! Ice restart file  
    250          CALL lim_var_glo2eqv                       ! ??? 
    251          ! 
    252          IF( ln_icectl )   CALL lim_ctl( kt )        ! alerts in case of model crash 
    253          ! 
    254          CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice ) 
     243         ! 
     244         IF( ln_icectl )  CALL lim_ctl( kt )        ! alerts in case of model crash 
    255245         ! 
    256246      ENDIF   ! End sea-ice time step only 
    257247 
    258       !--------------------------------! 
    259       ! --- at all ocean time step --- ! 
    260       !--------------------------------! 
    261       ! Update surface ocean stresses (only in ice-dynamic case) 
    262       !    otherwise the atm.-ocean stresses are used everywhere 
     248      !-------------------------! 
     249      ! --- Ocean time step --- ! 
     250      !-------------------------! 
     251      ! Update surface ocean stresses (only in ice-dynamic case) otherwise the atm.-ocean stresses are used everywhere 
    263252      IF( ln_limdyn )     CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) )  ! using before instantaneous surf. currents 
    264253!!gm   remark, the ocean-ice stress is not saved in ice diag call above .....  find a solution!!! 
    265254      ! 
    266       IF( nn_timing == 1 )  CALL timing_stop('sbc_ice_lim') 
     255      IF( nn_timing == 1 ) CALL timing_stop('sbc_ice_lim') 
    267256      ! 
    268257   END SUBROUTINE sbc_ice_lim 
     
    305294      ! 
    306295      CALL lim_itd_init                ! ice thickness distribution initialization 
     296      ! 
     297      CALL lim_hdf_init                ! set ice horizontal diffusion computation parameters 
    307298      ! 
    308299      CALL lim_thd_init                ! set ice thermodynics parameters 
     
    351342      !!------------------------------------------------------------------- 
    352343      INTEGER  ::   ios                 ! Local integer output status for namelist read 
    353       NAMELIST/namicerun/ jpl, nlay_i, nlay_s, cn_icerst_in, cn_icerst_out,   & 
     344      NAMELIST/namicerun/ jpl, nlay_i, nlay_s, cn_icerst_in, cn_icerst_indir, cn_icerst_out, cn_icerst_outdir,  & 
    354345         &                ln_limdyn, rn_amax, ln_limdiahsb, ln_limdiaout, ln_icectl, iiceprt, jiceprt   
    355346      !!------------------------------------------------------------------- 
     
    389380      r1_nlay_s = 1._wp / REAL( nlay_s, wp ) 
    390381      ! 
     382#if defined key_bdy 
     383      IF( lwp .AND. ln_limdiahsb )  CALL ctl_warn('online conservation check activated but it does not work with BDY') 
     384#endif 
     385      ! 
    391386   END SUBROUTINE ice_run 
    392387 
     
    476471 
    477472    
    478    SUBROUTINE ice_lim_flx( ptn_ice, palb_ice, pqns_ice, pqsr_ice,   & 
    479          &                          pdqn_ice, pqla_ice, pdql_ice, k_limflx ) 
     473   SUBROUTINE ice_lim_flx( ptn_ice, palb_ice, pqns_ice, pqsr_ice, pdqn_ice, pevap_ice, pdevap_ice, k_limflx ) 
    480474      !!--------------------------------------------------------------------- 
    481475      !!                  ***  ROUTINE ice_lim_flx  *** 
     
    495489      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqsr_ice   ! net solar flux 
    496490      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdqn_ice   ! non solar flux sensitivity 
    497       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqla_ice   ! latent heat flux 
    498       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdql_ice   ! latent heat flux sensitivity 
     491      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pevap_ice  ! sublimation 
     492      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdevap_ice ! sublimation sensitivity 
    499493      ! 
    500494      INTEGER  ::   jl      ! dummy loop index 
     
    505499      REAL(wp), POINTER, DIMENSION(:,:) :: z_qsr_m   ! Mean solar heat flux over all categories 
    506500      REAL(wp), POINTER, DIMENSION(:,:) :: z_qns_m   ! Mean non solar heat flux over all categories 
    507       REAL(wp), POINTER, DIMENSION(:,:) :: z_qla_m   ! Mean latent heat flux over all categories 
     501      REAL(wp), POINTER, DIMENSION(:,:) :: z_evap_m  ! Mean sublimation over all categories 
    508502      REAL(wp), POINTER, DIMENSION(:,:) :: z_dqn_m   ! Mean d(qns)/dT over all categories 
    509       REAL(wp), POINTER, DIMENSION(:,:) :: z_dql_m   ! Mean d(qla)/dT over all categories 
     503      REAL(wp), POINTER, DIMENSION(:,:) :: z_devap_m ! Mean d(evap)/dT over all categories 
    510504      !!---------------------------------------------------------------------- 
    511505 
     
    515509      SELECT CASE( k_limflx )                              !==  averaged on all ice categories  ==! 
    516510      CASE( 0 , 1 ) 
    517          CALL wrk_alloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m) 
    518          ! 
    519          z_qns_m(:,:) = fice_ice_ave ( pqns_ice (:,:,:) ) 
    520          z_qsr_m(:,:) = fice_ice_ave ( pqsr_ice (:,:,:) ) 
    521          z_dqn_m(:,:) = fice_ice_ave ( pdqn_ice (:,:,:) ) 
    522          z_qla_m(:,:) = fice_ice_ave ( pqla_ice (:,:,:) ) 
    523          z_dql_m(:,:) = fice_ice_ave ( pdql_ice (:,:,:) ) 
     511         CALL wrk_alloc( jpi,jpj, z_qsr_m, z_qns_m, z_evap_m, z_dqn_m, z_devap_m) 
     512         ! 
     513         z_qns_m  (:,:) = fice_ice_ave ( pqns_ice (:,:,:) ) 
     514         z_qsr_m  (:,:) = fice_ice_ave ( pqsr_ice (:,:,:) ) 
     515         z_dqn_m  (:,:) = fice_ice_ave ( pdqn_ice (:,:,:) ) 
     516         z_evap_m (:,:) = fice_ice_ave ( pevap_ice (:,:,:) ) 
     517         z_devap_m(:,:) = fice_ice_ave ( pdevap_ice (:,:,:) ) 
    524518         DO jl = 1, jpl 
    525             pdqn_ice(:,:,jl) = z_dqn_m(:,:) 
    526             pdql_ice(:,:,jl) = z_dql_m(:,:) 
     519            pdqn_ice  (:,:,jl) = z_dqn_m(:,:) 
     520            pdevap_ice(:,:,jl) = z_devap_m(:,:) 
    527521         END DO 
    528522         ! 
    529523         DO jl = 1, jpl 
    530             pqns_ice(:,:,jl) = z_qns_m(:,:) 
    531             pqsr_ice(:,:,jl) = z_qsr_m(:,:) 
    532             pqla_ice(:,:,jl) = z_qla_m(:,:) 
     524            pqns_ice (:,:,jl) = z_qns_m(:,:) 
     525            pqsr_ice (:,:,jl) = z_qsr_m(:,:) 
     526            pevap_ice(:,:,jl) = z_evap_m(:,:) 
    533527         END DO 
    534528         ! 
    535          CALL wrk_dealloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m) 
     529         CALL wrk_dealloc( jpi,jpj, z_qsr_m, z_qns_m, z_evap_m, z_dqn_m, z_devap_m) 
    536530      END SELECT 
    537531 
     
    543537         ztem_m(:,:) = fice_ice_ave ( ptn_ice  (:,:,:) )  
    544538         DO jl = 1, jpl 
    545             pqns_ice(:,:,jl) = pqns_ice(:,:,jl) + pdqn_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:)) 
    546             pqla_ice(:,:,jl) = pqla_ice(:,:,jl) + pdql_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:)) 
    547             pqsr_ice(:,:,jl) = pqsr_ice(:,:,jl) * ( 1._wp - palb_ice(:,:,jl) ) / ( 1._wp - zalb_m(:,:) )  
     539            pqns_ice (:,:,jl) = pqns_ice (:,:,jl) + pdqn_ice  (:,:,jl) * ( ptn_ice(:,:,jl) - ztem_m(:,:) ) 
     540            pevap_ice(:,:,jl) = pevap_ice(:,:,jl) + pdevap_ice(:,:,jl) * ( ptn_ice(:,:,jl) - ztem_m(:,:) ) 
     541            pqsr_ice (:,:,jl) = pqsr_ice (:,:,jl) * ( 1._wp - palb_ice(:,:,jl) ) / ( 1._wp - zalb_m(:,:) )  
    548542         END DO 
    549543         ! 
     
    555549   END SUBROUTINE ice_lim_flx 
    556550 
    557    SUBROUTINE sbc_lim_update 
    558       !!---------------------------------------------------------------------- 
    559       !!                  ***  ROUTINE sbc_lim_update  *** 
     551   SUBROUTINE sbc_lim_bef 
     552      !!---------------------------------------------------------------------- 
     553      !!                  ***  ROUTINE sbc_lim_bef  *** 
    560554      !! 
    561555      !! ** purpose :  store ice variables at "before" time step  
     
    571565      v_ice_b(:,:)     = v_ice(:,:) 
    572566       
    573    END SUBROUTINE sbc_lim_update 
     567   END SUBROUTINE sbc_lim_bef 
    574568 
    575569   SUBROUTINE sbc_lim_diag0 
     
    594588      wfx_spr(:,:) = 0._wp   ;    
    595589       
    596       hfx_in (:,:) = 0._wp   ;   hfx_out(:,:) = 0._wp 
    597590      hfx_thd(:,:) = 0._wp   ;    
    598591      hfx_snw(:,:) = 0._wp   ;   hfx_opw(:,:) = 0._wp 
     
    602595      hfx_spr(:,:) = 0._wp   ;   hfx_dif(:,:) = 0._wp  
    603596      hfx_err(:,:) = 0._wp   ;   hfx_err_rem(:,:) = 0._wp 
    604       hfx_err_dif(:,:) = 0._wp 
     597      hfx_err_dif(:,:) = 0._wp   ; 
    605598 
    606599      afx_tot(:,:) = 0._wp   ; 
    607600      afx_dyn(:,:) = 0._wp   ;   afx_thd(:,:) = 0._wp 
    608601 
    609       diag_heat_dhc(:,:) = 0._wp ; 
     602      diag_heat(:,:) = 0._wp ;   diag_smvi(:,:) = 0._wp ; 
     603      diag_vice(:,:) = 0._wp ;   diag_vsnw(:,:) = 0._wp ; 
    610604       
    611605   END SUBROUTINE sbc_lim_diag0 
    612        
     606 
     607      
    613608   FUNCTION fice_cell_ave ( ptab ) 
    614609      !!-------------------------------------------------------------------------- 
     
    620615       
    621616      fice_cell_ave (:,:) = 0.0_wp 
    622        
    623617      DO jl = 1, jpl 
    624618         fice_cell_ave (:,:) = fice_cell_ave (:,:) + a_i (:,:,jl) * ptab (:,:,jl) 
     
    636630 
    637631      fice_ice_ave (:,:) = 0.0_wp 
    638       WHERE ( at_i (:,:) .GT. 0.0_wp ) fice_ice_ave (:,:) = fice_cell_ave ( ptab (:,:,:)) / at_i (:,:) 
     632      WHERE ( at_i (:,:) > 0.0_wp ) fice_ice_ave (:,:) = fice_cell_ave ( ptab (:,:,:)) / at_i (:,:) 
    639633 
    640634   END FUNCTION fice_ice_ave 
Note: See TracChangeset for help on using the changeset viewer.