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 @ 8592

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

rewrite ice albedo routine
File size: 12.5 KB
Line 
1MODULE icealb
2   !!======================================================================
3   !!                       ***  MODULE  icealb  ***
4   !! Atmospheric forcing:  Albedo over sea ice
5   !!=====================================================================
6   !! History :  4.0  ! 2017-07  (C. Rousset) Split ice and ocean albedos
7   !!----------------------------------------------------------------------
8#if defined key_lim3
9   !!----------------------------------------------------------------------
10   !!   'key_lim3'                                       ESIM sea-ice model
11   !!----------------------------------------------------------------------
12   !!   ice_alb        : albedo for ice (clear and overcast skies)
13   !!   ice_alb_init   : initialisation of albedo computation
14   !!----------------------------------------------------------------------
15   USE ice, ONLY: jpl ! sea-ice: number of categories
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   USE timing         ! Timing
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC   ice_alb_init   ! called in icestp
27   PUBLIC   ice_alb        ! called in iceforcing.F90 and iceupdate.F90
28
29   REAL(wp), PUBLIC, PARAMETER ::   rn_alb_oce = 0.066   !: ocean or lead albedo (Pegau and Paulson, Ann. Glac. 2001)
30
31   REAL(wp) , PARAMETER ::   ppc1    = 0.05    ! snow thickness (only for nn_ice_alb=0)
32   REAL(wp) , PARAMETER ::   ppc2    = 0.10    !  "        "
33   REAL(wp) , PARAMETER ::   ppcloud = 0.06    ! cloud effect on albedo (only-for nn_ice_alb=0)
34   REAL(wp) , PARAMETER ::   pp1_c1 = 1. / ppc1
35   REAL(wp) , PARAMETER ::   pp1_c2 = 1. / ppc2
36   !
37   ! ** albedo namelist (namalb)
38   REAL(wp) ::   rn_alb_sdry      ! dry snow albedo
39   REAL(wp) ::   rn_alb_smlt      ! melting snow albedo
40   REAL(wp) ::   rn_alb_idry      ! dry ice albedo
41   REAL(wp) ::   rn_alb_imlt      ! bare puddled ice albedo
42   REAL(wp) ::   rn_alb_dpnd      ! ponded ice albedo
43
44   !!----------------------------------------------------------------------
45   !! NEMO/ICE 4.0 , NEMO Consortium (2017)
46   !! $Id: icealb.F90 8268 2017-07-03 15:01:04Z clem $
47   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
48   !!----------------------------------------------------------------------
49CONTAINS
50
51   SUBROUTINE ice_alb( pt_su, ph_ice, ph_snw, pafrac_pnd, ph_pnd, palb_cs, palb_os )
52      !!----------------------------------------------------------------------
53      !!               ***  ROUTINE ice_alb  ***
54      !!         
55      !! ** Purpose :   Computation of the albedo of the snow/ice system
56      !!       
57      !! ** Method  :   The scheme is "home made" (for cloudy skies) and based on Brandt et al. (J. Climate 2005)
58      !!                                                                      and Grenfell & Perovich (JGR 2004)
59      !!                  1) Albedo dependency on ice thickness follows the findings from Brandt et al (2005)
60      !!                     which are an update of Allison et al. (JGR 1993) ; Brandt et al. 1999
61      !!                     0-5cm  : linear function of ice thickness
62      !!                     5-150cm: log    function of ice thickness
63      !!                     > 150cm: constant
64      !!                  2) Albedo dependency on snow thickness follows the findings from Grenfell & Perovich (2004)
65      !!                     i.e. it increases as -EXP(-snw_thick/0.02) during freezing and -EXP(-snw_thick/0.03) during melting
66      !!                  3) Albedo dependency on clouds is speculated from measurements of Grenfell and Perovich (2004)
67      !!                     i.e. cloudy-clear albedo depend on cloudy albedo following a 2d order polynomial law
68      !!                  4) The needed 4 parameters are: dry and melting snow, freezing ice and bare puddled ice
69      !!
70      !!                     compilation of values from literature (reference overcast sky values)
71      !!                        rn_alb_sdry = 0.85      ! dry snow
72      !!                        rn_alb_smlt = 0.75      ! melting snow
73      !!                        rn_alb_idry = 0.60      ! bare frozen ice
74      !!                        rn_alb_imlt = 0.50      ! bare puddled ice albedo
75      !!                        rn_alb_dpnd = 0.36      ! ponded-ice overcast albedo (Lecomte et al, 2015)
76      !!                                                ! early melt pnds 0.27, late melt ponds 0.14 Grenfell & Perovich
77      !!                     Perovich et al 2002 (Sheba) => the only dataset for which all types of ice/snow were retrieved
78      !!                        rn_alb_sdry = 0.85      ! dry snow
79      !!                        rn_alb_smlt = 0.72      ! melting snow
80      !!                        rn_alb_idry = 0.65      ! bare frozen ice
81      !!                     Brandt et al 2005 (East Antarctica)
82      !!                        rn_alb_sdry = 0.87      ! dry snow
83      !!                        rn_alb_smlt = 0.82      ! melting snow
84      !!                        rn_alb_idry = 0.54      ! bare frozen ice
85      !!
86      !! ** Note    :   The old parameterization from Shine & Henderson-Sellers (not here anymore) presented several misconstructions
87      !!                  1) ice albedo when ice thick. tends to 0 is different than ocean albedo
88      !!                  2) for small ice thick. covered with some snow (<3cm?), albedo is larger
89      !!                     under melting conditions than under freezing conditions
90      !!                  3) the evolution of ice albedo as a function of ice thickness shows 
91      !!                     3 sharp inflexion points (at 5cm, 100cm and 150cm) that look highly unrealistic
92      !!
93      !! References :   Shine & Henderson-Sellers 1985, JGR, 90(D1), 2243-2250.
94      !!                Brandt et al. 2005, J. Climate, vol 18
95      !!                Grenfell & Perovich 2004, JGR, vol 109
96      !!----------------------------------------------------------------------
97      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pt_su        !  ice surface temperature (Kelvin)
98      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_ice       !  sea-ice thickness
99      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_snw       !  snow depth
100      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pafrac_pnd   !  melt pond relative fraction (per unit ice area)
101      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_pnd       !  melt pond depth
102      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   palb_cs      !  albedo of ice under clear    sky
103      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   palb_os      !  albedo of ice under overcast sky
104      !
105      INTEGER  ::   ji, jj, jl                ! dummy loop indices
106      REAL(wp) ::   z1_c1, z1_c2,z1_c3, z1_c4 ! local scalar
107      REAL(wp) ::   z1_href_pnd               ! inverse of the characteristic length scale (Lecomte et al. 2015)
108      REAL(wp) ::   zalb_pnd, zafrac_pnd      ! ponded sea ice albedo & relative pound fraction
109      REAL(wp) ::   zalb_ice, zafrac_ice      ! bare sea ice albedo & relative ice fraction
110      REAL(wp) ::   zalb_snw, zafrac_snw      ! snow-covered sea ice albedo & relative snow fraction
111      !!---------------------------------------------------------------------
112      !
113      IF( nn_timing == 1 )   CALL timing_start('icealb')
114      !
115      z1_href_pnd = 0.05
116      z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) ) 
117      z1_c2 = 1. / 0.05
118      z1_c3 = 1. / 0.02
119      z1_c4 = 1. / 0.03
120      !
121      DO jl = 1, jpl
122         DO jj = 1, jpj
123            DO ji = 1, jpi
124               !                       !--- Specific snow, ice and pond fractions (for now, we prevent melt ponds and snow at the same time)
125               IF( ph_snw(ji,jj,jl) == 0._wp ) THEN
126                  zafrac_snw = 0._wp
127                  zafrac_pnd = pafrac_pnd(ji,jj,jl)
128                  zafrac_ice = 1._wp - zafrac_pnd
129               ELSE
130                  zafrac_snw = 1._wp      ! Snow fully "shades" melt ponds and ice
131                  zafrac_pnd = 0._wp
132                  zafrac_ice = 0._wp
133               ENDIF
134               !
135               !                       !--- Bare ice albedo (for hi > 150cm)
136               IF( zafrac_pnd > 0._wp ) THEN
137                  zalb_ice = rn_alb_idry
138               ELSE
139                  IF( ph_snw(ji,jj,jl) == 0._wp .AND. pt_su(ji,jj,jl) >= rt0 ) THEN  ;   zalb_ice = rn_alb_imlt
140                  ELSE                                                               ;   zalb_ice = rn_alb_idry   ;   ENDIF
141               ENDIF
142               !                       !--- Bare ice albedo (for hi < 150cm)
143               IF( 0.05 < ph_ice(ji,jj,jl) .AND. ph_ice(ji,jj,jl) <= 1.5 ) THEN      ! 5cm < hi < 150cm
144                  zalb_ice = zalb_ice    + ( 0.18 - zalb_ice   ) * z1_c1 * ( LOG(1.5) - LOG(ph_ice(ji,jj,jl)) )
145               ELSEIF( ph_ice(ji,jj,jl) <= 0.05 ) THEN                               ! 0cm < hi < 5cm
146                  zalb_ice = rn_alb_oce  + ( 0.18 - rn_alb_oce ) * z1_c2 * ph_ice(ji,jj,jl)
147               ENDIF
148               !
149               !                       !--- Snow-covered ice albedo (freezing, melting cases)
150               IF( pt_su(ji,jj,jl) < rt0_snow ) THEN
151                  zalb_snw = rn_alb_sdry - ( rn_alb_sdry - zalb_ice ) * EXP( - ph_snw(ji,jj,jl) * z1_c3 )
152               ELSE
153                  zalb_snw = rn_alb_smlt - ( rn_alb_smlt - zalb_ice ) * EXP( - ph_snw(ji,jj,jl) * z1_c4 )
154               ENDIF
155               !                       !--- Ponded ice albedo
156               IF( zafrac_pnd > 0._wp ) THEN
157                  zalb_pnd = rn_alb_dpnd - ( rn_alb_dpnd - zalb_ice ) * EXP( - ph_pnd(ji,jj,jl) * z1_href_pnd ) 
158               ELSE
159                  zalb_pnd = rn_alb_dpnd
160               ENDIF
161               !                       !--- Surface albedo is weighted mean of snow, ponds and bare ice contributions
162               palb_os(ji,jj,jl) = zafrac_snw * zalb_snw + zafrac_pnd * zalb_pnd + zafrac_ice * zalb_ice
163               !
164            END DO
165         END DO
166      END DO
167      !
168      palb_cs(:,:,:) = palb_os(:,:,:) - ( - 0.1010 * palb_os(:,:,:) * palb_os(:,:,:) + 0.1933 * palb_os(:,:,:) - 0.0148 )
169      !
170      IF( nn_timing == 1 )   CALL timing_stop('icealb')
171      !
172   END SUBROUTINE ice_alb
173
174
175   SUBROUTINE ice_alb_init
176      !!----------------------------------------------------------------------
177      !!                 ***  ROUTINE alb_init  ***
178      !!
179      !! ** Purpose :   initializations for the albedo parameters
180      !!
181      !! ** Method  :   Read the namelist namalb
182      !!----------------------------------------------------------------------
183      INTEGER  ::   ios                 ! Local integer output status for namelist read
184      !!
185      NAMELIST/namalb/ rn_alb_sdry, rn_alb_smlt, rn_alb_idry, rn_alb_imlt, rn_alb_dpnd
186      !!----------------------------------------------------------------------
187      !
188      REWIND( numnam_ice_ref )              ! Namelist namalb in reference namelist : Albedo parameters
189      READ  ( numnam_ice_ref, namalb, IOSTAT = ios, ERR = 901)
190901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namalb in reference namelist', lwp )
191
192      REWIND( numnam_ice_cfg )              ! Namelist namalb in configuration namelist : Albedo parameters
193      READ  ( numnam_ice_cfg, namalb, IOSTAT = ios, ERR = 902 )
194902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namalb in configuration namelist', lwp )
195      IF(lwm) WRITE ( numoni, namalb )
196      !
197      IF(lwp) THEN                      ! Control print
198         WRITE(numout,*)
199         WRITE(numout,*) 'ice_alb_init: set albedo parameters'
200         WRITE(numout,*) '~~~~~~~~~~~~'
201         WRITE(numout,*) '   Namelist namalb:'
202         WRITE(numout,*) '      albedo of dry snow                   rn_alb_sdry = ', rn_alb_sdry
203         WRITE(numout,*) '      albedo of melting snow               rn_alb_smlt = ', rn_alb_smlt
204         WRITE(numout,*) '      albedo of dry ice                    rn_alb_idry = ', rn_alb_idry
205         WRITE(numout,*) '      albedo of bare puddled ice           rn_alb_imlt = ', rn_alb_imlt
206         WRITE(numout,*) '      albedo of ponded ice                 rn_alb_dpnd = ', rn_alb_dpnd
207      ENDIF
208      !
209   END SUBROUTINE ice_alb_init
210
211#else
212   !!----------------------------------------------------------------------
213   !!   Default option           Dummy module         NO ESIM sea-ice model
214   !!----------------------------------------------------------------------
215#endif
216
217   !!======================================================================
218END MODULE icealb
Note: See TracBrowser for help on using the repository browser.