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.
icealb.F90 in branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icealb.F90 @ 8517

Last change on this file since 8517 was 8517, checked in by clem, 7 years ago

changes in style - part6 - one more round

File size: 20.6 KB
RevLine 
[8426]1MODULE icealb
2   !!======================================================================
3   !!                       ***  MODULE  icealb  ***
4   !! Ocean forcing:  bulk thermohaline forcing of the ice
5   !!=====================================================================
[8486]6   !! History :  4.0  ! 2017-07  (C. Rousset) Split ice and ocean albedos
[8426]7   !!----------------------------------------------------------------------
8#if defined key_lim3
9   !!----------------------------------------------------------------------
[8486]10   !!   'key_lim3'                                    LIM 3.0 sea-ice model
[8426]11   !!----------------------------------------------------------------------
[8486]12   !!   ice_alb       : albedo for ice (clear and overcast skies)
13   !!   alb_init      : initialisation of albedo computation
14   !!----------------------------------------------------------------------
[8505]15   USE ice, ONLY: jpl ! sea-ice: number of categories
[8426]16   USE phycst         ! physical constants
17   !
18   USE in_out_manager ! I/O manager
19   USE lib_mpp        ! MPP library
20   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
21
22   IMPLICIT NONE
23   PRIVATE
24
[8505]25   PUBLIC   ice_alb_init   ! called in icestp
26   PUBLIC   ice_alb        ! called in iceforcing.F90 and iceupdate.F90
[8426]27
[8486]28   REAL(wp), PUBLIC, PARAMETER ::   rn_alb_oce = 0.066   !: ocean or lead albedo (Pegau and Paulson, Ann. Glac. 2001)
29
[8498]30   REAL(wp) , PARAMETER ::   rc1    = 0.05    ! snow thickness (only for nn_ice_alb=0)
31   REAL(wp) , PARAMETER ::   rc2    = 0.10    !  "        "
32   REAL(wp) , PARAMETER ::   rcloud = 0.06    ! cloud effect on albedo (only-for nn_ice_alb=0)
33   REAL(wp) , PARAMETER ::   r1_c1 = 1. / rc1
34   REAL(wp) , PARAMETER ::   r1_c2 = 1. / rc2
[8505]35   !
[8512]36   ! ** albedo namelist (namice_alb)
[8505]37   INTEGER  ::   nn_ice_alb       ! type of albedo scheme: 0: Shine & Henderson-Sellers (JGR 1985)
38   !                                      !                         1: "home made" based on Brandt et al. (JClim 2005)
39   !                                      !                             and Grenfell & Perovich (JGR 2004)
40   !                                      !                         2: Same as 1 but with melt ponds
41   REAL(wp) ::   rn_alb_sdry      ! dry snow albedo
42   REAL(wp) ::   rn_alb_smlt      ! melting snow albedo
43   REAL(wp) ::   rn_alb_idry      ! dry ice albedo
44   REAL(wp) ::   rn_alb_imlt      ! bare puddled ice albedo
45   REAL(wp) ::   rn_alb_dpnd      ! ponded ice albedo
[8426]46
47   !!----------------------------------------------------------------------
[8486]48   !! NEMO/ICE 4.0 , NEMO Consortium (2017)
[8426]49   !! $Id: icealb.F90 8268 2017-07-03 15:01:04Z clem $
50   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
51   !!----------------------------------------------------------------------
52CONTAINS
53
54   SUBROUTINE ice_alb( pt_ice, ph_ice, ph_snw, pafrac_pnd, ph_pnd, ld_pnd, pa_ice_cs, pa_ice_os )
55      !!----------------------------------------------------------------------
56      !!               ***  ROUTINE ice_alb  ***
57      !!         
58      !! ** Purpose :   Computation of the albedo of the snow/ice system
59      !!       
60      !! ** Method  :   Two schemes are available (from namelist parameter nn_ice_alb)
61      !!                  0: the scheme is that of Shine & Henderson-Sellers (JGR 1985) for clear-skies
62      !!                  1: the scheme is "home made" (for cloudy skies) and based on Brandt et al. (J. Climate 2005)
63      !!                                                                           and Grenfell & Perovich (JGR 2004)
64      !!                  2: fractional surface-based formulation of scheme 1 (NEW)
65      !!                Description of scheme 1:
66      !!                  1) Albedo dependency on ice thickness follows the findings from Brandt et al (2005)
67      !!                     which are an update of Allison et al. (JGR 1993) ; Brandt et al. 1999
68      !!                     0-5cm  : linear function of ice thickness
69      !!                     5-150cm: log    function of ice thickness
70      !!                     > 150cm: constant
71      !!                  2) Albedo dependency on snow thickness follows the findings from Grenfell & Perovich (2004)
72      !!                     i.e. it increases as -EXP(-snw_thick/0.02) during freezing and -EXP(-snw_thick/0.03) during melting
73      !!                  3) Albedo dependency on clouds is speculated from measurements of Grenfell and Perovich (2004)
74      !!                     i.e. cloudy-clear albedo depend on cloudy albedo following a 2d order polynomial law
75      !!                  4) The needed 4 parameters are: dry and melting snow, freezing ice and bare puddled ice
76      !!
77      !! ** Note    :   The parameterization from Shine & Henderson-Sellers presents several misconstructions:
78      !!                  1) ice albedo when ice thick. tends to 0 is different than ocean albedo
79      !!                  2) for small ice thick. covered with some snow (<3cm?), albedo is larger
80      !!                     under melting conditions than under freezing conditions
81      !!                  3) the evolution of ice albedo as a function of ice thickness shows 
82      !!                     3 sharp inflexion points (at 5cm, 100cm and 150cm) that look highly unrealistic
83      !!
84      !! References :   Shine & Henderson-Sellers 1985, JGR, 90(D1), 2243-2250.
85      !!                Brandt et al. 2005, J. Climate, vol 18
86      !!                Grenfell & Perovich 2004, JGR, vol 109
87      !!----------------------------------------------------------------------
[8486]88      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pt_ice       !  ice surface temperature (Kelvin)
89      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_ice       !  sea-ice thickness
90      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_snw       !  snow depth
91      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pafrac_pnd   !  melt pond relative fraction (per unit ice area)
92      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_pnd       !  melt pond depth
93      LOGICAL , INTENT(in   )                   ::   ld_pnd       !  melt ponds radiatively active or not
94      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   pa_ice_cs    !  albedo of ice under clear    sky
95      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   pa_ice_os    !  albedo of ice under overcast sky
[8426]96      !
[8486]97      INTEGER  ::   ji, jj, jl                ! dummy loop indices
98      REAL(wp) ::   zswitch, z1_c1, z1_c2     ! local scalar
99      REAL(wp) ::   z1_href_pnd               !   -      -                 
[8426]100      REAL(wp) ::   zalb_sm, zalb_sf, zalb_st ! albedo of snow melting, freezing, total
101      !
[8486]102      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb, zalb_it   ! intermediate variable & albedo of ice (snow free)
[8426]103!! MV MP
[8486]104      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb_pnd        ! ponded sea ice albedo
105      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb_ice        ! bare sea ice albedo
106      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb_snw        ! snow-covered sea ice albedo
107      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zafrac_snw      ! relative snow fraction
108      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zafrac_ice      ! relative ice fraction
109      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zafrac_pnd      ! relative ice fraction (effective)
[8426]110      !!---------------------------------------------------------------------
111
112      !-----------------------------------------------------
113      !  Snow-free albedo (no ice thickness dependence yet)
114      !-----------------------------------------------------
115      !
116      ! Part common to nn_ice_alb = 0, 1, 2
117      !
118      IF ( .NOT. ld_pnd ) THEN   !--- No melt ponds OR radiatively inactive melt ponds
119         ! Bare ice albedo is prescribed, with implicit assumption on pond fraction and depth
[8498]120         WHERE     ( ph_snw == 0._wp .AND. pt_ice >= rt0 )   ;   zalb(:,:,:) = rn_alb_imlt
121         ELSEWHERE                                           ;   zalb(:,:,:) = rn_alb_idry
122         END WHERE
[8426]123      ENDIF
124
[8486]125      SELECT CASE ( nn_ice_alb )       ! select a parameterization
126      !
127      !              !------------------------------------------
128      CASE( 0 )      !  Shine and Henderson-Sellers (1985)     !   (based on clear sky values)
129         !           !------------------------------------------
130         !
131         !                       ! Thickness-dependent bare ice albedo
[8498]132         WHERE    ( 1.5  < ph_ice                     )   ;   zalb_it = zalb
133         ELSEWHERE( 1.0  < ph_ice .AND. ph_ice <= 1.5 )   ;   zalb_it = 0.472  + 2.0 * ( zalb - 0.472 ) * ( ph_ice - 1.0 )
134         ELSEWHERE( 0.05 < ph_ice .AND. ph_ice <= 1.0 )   ;   zalb_it = 0.2467 + 0.7049 * ph_ice              &
135            &                                                                  - 0.8608 * ph_ice * ph_ice     &
136            &                                                                  + 0.3812 * ph_ice * ph_ice * ph_ice
137         ELSEWHERE                                        ;   zalb_it = 0.1    + 3.6    * ph_ice
[8426]138         END WHERE
[8486]139         !
140         IF ( ld_pnd ) THEN      ! Depth-dependent ponded ice albedo
141            z1_href_pnd = 1. / 0.05    ! inverse of the characteristic length scale (Lecomte et al. 2015)
142            zalb_pnd  = rn_alb_dpnd - ( rn_alb_dpnd - zalb_it ) * EXP( - ph_pnd * z1_href_pnd ) 
143            !                          ! Snow-free ice albedo is a function of pond fraction
[8498]144            WHERE ( ph_snw == 0._wp .AND. pt_ice >= rt0 )   
[8486]145               zalb_it = zalb_it * ( 1. - pafrac_pnd  ) + zalb_pnd * pafrac_pnd 
146            END WHERE
[8426]147         ENDIF 
[8486]148         !
149!!gm: optimization ( rn_alb_smlt - rn_alb_imlt ) * r1_c2 can be computed one for all
150!!gm  before the DO-loop below
[8426]151         DO jl = 1, jpl
152            DO jj = 1, jpj
153               DO ji = 1, jpi
[8486]154                  !                    ! Freezing snow
[8498]155                  ! no effect of underlying ice layer IF snow thickness > c1. Albedo does not depend on snow thick if > rc2
156                  zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ( ph_snw(ji,jj,jl) - rc1 ) ) )
[8486]157                  zalb_sf = ( 1._wp - zswitch ) * (  zalb_it(ji,jj,jl)  &
158                     &                                   + ph_snw(ji,jj,jl) * ( rn_alb_sdry - zalb_it(ji,jj,jl) ) * r1_c1  )   &
159                     &    +           zswitch   * rn_alb_sdry 
160                     !
161                  !                    ! Melting snow
[8498]162                  ! no effect of underlying ice layer. Albedo does not depend on snow thick IF > rc2
163                  zswitch = MAX( 0._wp , SIGN( 1._wp , ph_snw(ji,jj,jl) - rc2 ) )
[8486]164                  zalb_sm = ( 1._wp - zswitch ) * ( rn_alb_imlt + ph_snw(ji,jj,jl) * ( rn_alb_smlt - rn_alb_imlt ) * r1_c2 )   &
165                     &      +         zswitch   *   rn_alb_smlt 
166                     !
167                  !                    ! Snow albedo
168                  zswitch =  MAX( 0._wp , SIGN( 1._wp , pt_ice(ji,jj,jl) - rt0_snow ) )   
169                  zalb_st =  zswitch * zalb_sm + ( 1._wp - zswitch ) * zalb_sf
[8426]170                  !
[8486]171                  !                    ! Surface albedo
[8426]172                  zswitch             = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ph_snw(ji,jj,jl) ) )
173                  pa_ice_cs(ji,jj,jl) = zswitch * zalb_st + ( 1._wp - zswitch ) * zalb_it(ji,jj,jl)
174                  !
175               END DO
176            END DO
177         END DO
[8486]178         !
[8426]179         pa_ice_os(:,:,:) = pa_ice_cs(:,:,:) + rcloud       ! Oberhuber correction for overcast sky
[8486]180         !
181         !
182         !           !------------------------------------------
183      CASE( 1 )      !  New parameterization (2016)            !    (based on overcast sky values)
184         !           !------------------------------------------
185         !
186         !                    compilation of values from literature (reference overcast sky values)
187         !                          rn_alb_sdry = 0.85      ! dry snow
188         !                          rn_alb_smlt = 0.75      ! melting snow
189         !                          rn_alb_idry = 0.60      ! bare frozen ice
190         !                          rn_alb_dpnd = 0.36      ! ponded-ice overcast albedo (Lecomte et al, 2015)
191         !                                                  ! early melt pnds 0.27, late melt ponds 0.14 Grenfell & Perovich
192         !                    Perovich et al 2002 (Sheba) => the only dataset for which all types of ice/snow were retrieved
193         !                          rn_alb_sdry = 0.85      ! dry snow
194         !                          rn_alb_smlt = 0.72      ! melting snow
195         !                          rn_alb_idry = 0.65      ! bare frozen ice
196         !                    Brandt et al 2005 (East Antarctica)
197         !                          rn_alb_sdry = 0.87      ! dry snow
198         !                          rn_alb_smlt = 0.82      ! melting snow
199         !                          rn_alb_idry = 0.54      ! bare frozen ice
200         !
201         !              !--- Computation of snow-free ice albedo
[8426]202         z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) ) 
203         z1_c2 = 1. / 0.05
204
[8486]205         !              !--- Accounting for the ice-thickness dependency
206         WHERE     ( 1.5  < ph_ice                     )   ;   zalb_it = zalb
207         ELSE WHERE( 0.05 < ph_ice .AND. ph_ice <= 1.5 )   ;   zalb_it = zalb + ( 0.18 - zalb ) * z1_c1 * ( LOG(1.5) - LOG(ph_ice) )
208         ELSE WHERE                                        ;   zalb_it = rn_alb_oce + ( 0.18 - rn_alb_oce ) * z1_c2 * ph_ice
[8426]209         END WHERE
[8486]210         !
211         IF ( ld_pnd ) THEN      ! Depth-dependent ponded ice albedo
212            z1_href_pnd = 0.05        ! inverse of the characteristic length scale (Lecomte et al. 2015)
213            zalb_pnd  = rn_alb_dpnd - ( rn_alb_dpnd - zalb_it ) * EXP( - ph_pnd * z1_href_pnd ) 
214            !
215            !                    ! Snow-free ice albedo is weighted mean of ponded ice and bare ice contributions
[8498]216            WHERE ( ph_snw == 0._wp .AND. pt_ice >= rt0 )
[8486]217               zalb_it = zalb_it * ( 1. - pafrac_pnd  ) + zalb_pnd * pafrac_pnd
218            END WHERE
[8426]219         ENDIF 
[8486]220         !
221         !              !--- Overcast sky surface albedo (accounting for snow, ice melt ponds)
[8426]222         z1_c1 = 1. / 0.02
223         z1_c2 = 1. / 0.03
224         DO jl = 1, jpl
225            DO jj = 1, jpj
226               DO ji = 1, jpi
227                  ! Snow depth dependence of snow albedo
228                  zalb_sf = rn_alb_sdry - ( rn_alb_sdry - zalb_it(ji,jj,jl)) * EXP( - ph_snw(ji,jj,jl) * z1_c1 );
229                  zalb_sm = rn_alb_smlt - ( rn_alb_smlt - zalb_it(ji,jj,jl)) * EXP( - ph_snw(ji,jj,jl) * z1_c2 );
230
231                  ! Snow albedo (MV I guess we could use rt0 instead of rt0_snow)
232                  zswitch = MAX( 0._wp , SIGN( 1._wp , pt_ice(ji,jj,jl) - rt0_snow ) )   
233                  zalb_st = zswitch * zalb_sm + ( 1._wp - zswitch ) * zalb_sf
234
235                  ! Surface albedo   
236                  zswitch             = MAX( 0._wp , SIGN( 1._wp , - ph_snw(ji,jj,jl) ) )
237                  pa_ice_os(ji,jj,jl) = ( 1._wp - zswitch ) * zalb_st + zswitch *  zalb_it(ji,jj,jl)
238
239              END DO
240            END DO
241         END DO
[8486]242         !
243         !              !--- Clear sky surface albedo
[8426]244         pa_ice_cs = pa_ice_os - ( - 0.1010 * pa_ice_os * pa_ice_os + 0.1933 * pa_ice_os - 0.0148 ); 
[8486]245         !
246         !
247         !           !------------------------------------------
248      CASE( 2 )      !  Fractional surface-based parameterization (2017)
249         !           !------------------------------------------
250         !              MV: I propose this formulation that is more elegant, and more easy to expand towards
251         !              varying pond and snow fraction.
252         !              Formulation 1 has issues to handle ponds and snow together that can't easily be fixed.
253         !              This one handles it better, I believe.
254         !
255         !                 !==  Snow, bare ice and ponded ice fractions  ==!
256         !
257         !                       ! Specific fractions (zafrac) refer to relative area covered by snow within each ice category
258         !
259         !                       !--- Effective pond fraction (for now, we prevent melt ponds and snow at the same time)
260         zafrac_pnd = 0._wp
[8498]261         IF ( ld_pnd ) THEN
262            WHERE( ph_snw == 0._wp )   zafrac_pnd = pafrac_pnd   ! Snow fully "shades" melt ponds
[8486]263         ENDIF         
264         !
265         !                       !--- Specific snow fraction (for now, prescribed)
[8498]266         WHERE( ph_snw > 0._wp )   ;   zafrac_snw = 1.
267         ELSEWHERE                 ;   zafrac_snw = 0.
[8486]268         END WHERE
269         !
270         !                       !--- Specific ice fraction
271         zafrac_ice = 1. - zafrac_snw - zafrac_pnd
272         !
273         !                 !==  Snow-covered, pond-covered, and bare ice albedos  ==!
274         !
275         !                       !--- Bare ice albedo
276         z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) ) 
277         z1_c2 = 1. / 0.05
[8498]278         WHERE    ( 1.5  < ph_ice                     )   ;   zalb_ice = zalb
279         ELSEWHERE( 0.05 < ph_ice .AND. ph_ice <= 1.5 )   ;   zalb_ice = zalb       + ( 0.18 - zalb       ) * z1_c1 *  &
[8426]280            &                                                                       ( LOG(1.5) - LOG(ph_ice) )
[8498]281         ELSEWHERE                                        ;   zalb_ice = rn_alb_oce + ( 0.18 - rn_alb_oce ) * z1_c2 * ph_ice
[8486]282         END WHERE
283         !
284         z1_c1 = 1. / 0.02       !--- Snow-covered ice albedo (freezing, melting cases)
285         z1_c2 = 1. / 0.03
286         !
[8498]287         WHERE( pt_ice < rt0_snow )   ;   zalb_snw = rn_alb_sdry - ( rn_alb_sdry - zalb_ice ) * EXP( - ph_snw * z1_c1 )
288         ELSEWHERE                    ;   zalb_snw = rn_alb_smlt - ( rn_alb_smlt - zalb_ice ) * EXP( - ph_snw * z1_c2 )
[8486]289         END WHERE
290         !
291         IF ( ld_pnd ) THEN      !--- Depth-dependent ponded ice albedo
292            z1_href_pnd = 0.05        ! inverse of the characteristic length scale (Lecomte et al. 2015)
293            zalb_pnd  = rn_alb_dpnd - ( rn_alb_dpnd - zalb_ice ) * EXP( - ph_pnd * z1_href_pnd ) 
294         ELSE
295            zalb_pnd  = rn_alb_dpnd
296         ENDIF
297         !                       !--- Surface albedo is weighted mean of snow, ponds and bare ice contributions
298         pa_ice_os = zafrac_snw * zalb_snw  +  zafrac_pnd * zalb_pnd  +  zafrac_ice * zalb_ice
299         pa_ice_cs = pa_ice_os - ( - 0.1010 * pa_ice_os * pa_ice_os + 0.1933 * pa_ice_os - 0.0148 )
[8426]300
301      END SELECT
302      !
303   END SUBROUTINE ice_alb
304
[8486]305
[8505]306   SUBROUTINE ice_alb_init
[8426]307      !!----------------------------------------------------------------------
308      !!                 ***  ROUTINE alb_init  ***
309      !!
310      !! ** Purpose :   initializations for the albedo parameters
311      !!
[8512]312      !! ** Method  :   Read the namelist namice_alb
[8426]313      !!----------------------------------------------------------------------
314      INTEGER  ::   ios                 ! Local integer output status for namelist read
[8486]315      !!
[8512]316      NAMELIST/namice_alb/ nn_ice_alb, rn_alb_sdry, rn_alb_smlt, rn_alb_idry, rn_alb_imlt, rn_alb_dpnd
[8426]317      !!----------------------------------------------------------------------
318      !
[8512]319      REWIND( numnam_ice_ref )              ! Namelist namice_alb in reference namelist : Albedo parameters
320      READ  ( numnam_ice_ref, namice_alb, IOSTAT = ios, ERR = 901)
321901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namice_alb in reference namelist', lwp )
[8426]322
[8517]323      REWIND( numnam_ice_cfg )              ! Namelist namice_alb in configuration namelist : Albedo parameters
[8512]324      READ  ( numnam_ice_cfg, namice_alb, IOSTAT = ios, ERR = 902 )
325902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namice_alb in configuration namelist', lwp )
326      IF(lwm) WRITE ( numoni, namice_alb )
[8426]327      !
328      IF(lwp) THEN                      ! Control print
329         WRITE(numout,*)
[8514]330         WRITE(numout,*) 'ice_alb_init: set albedo parameters'
331         WRITE(numout,*) '~~~~~~~~~~~~'
[8512]332         WRITE(numout,*) '   Namelist namice_alb : albedo '
[8426]333         WRITE(numout,*) '      choose the albedo parameterization   nn_ice_alb  = ', nn_ice_alb
334         WRITE(numout,*) '      albedo of dry snow                   rn_alb_sdry = ', rn_alb_sdry
335         WRITE(numout,*) '      albedo of melting snow               rn_alb_smlt = ', rn_alb_smlt
336         WRITE(numout,*) '      albedo of dry ice                    rn_alb_idry = ', rn_alb_idry
337         WRITE(numout,*) '      albedo of bare puddled ice           rn_alb_imlt = ', rn_alb_imlt
338         WRITE(numout,*) '      albedo of ponded ice                 rn_alb_dpnd = ', rn_alb_dpnd
339      ENDIF
340      !
[8505]341   END SUBROUTINE ice_alb_init
[8426]342
343#else
344   !!----------------------------------------------------------------------
345   !!   Default option           Dummy module      NO LIM 3.0 sea-ice model
346   !!----------------------------------------------------------------------
347#endif
348
349   !!======================================================================
350END MODULE icealb
Note: See TracBrowser for help on using the repository browser.