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

source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/albedoice.F90 @ 8373

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

remove most of wrk_alloc

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