Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
MODULE icealb
!!======================================================================
!! *** MODULE icealb ***
!! Atmospheric forcing: Albedo over sea ice
!!=====================================================================
!! History : 4.0 ! 2017-07 (C. Rousset) Split ice and ocean albedos
!! 4.0 ! 2018 (many people) SI3 [aka Sea Ice cube]
!!----------------------------------------------------------------------
#if defined key_si3
!!----------------------------------------------------------------------
!! 'key_si3' SI3 sea-ice model
!!----------------------------------------------------------------------
!! ice_alb : albedo for ice (clear and overcast skies)
!! ice_alb_init : initialisation of albedo computation
!!----------------------------------------------------------------------
USE phycst ! physical constants
USE dom_oce ! domain: ocean
USE ice, ONLY: jpl ! sea-ice: number of categories
USE icevar ! sea-ice: operations
!
USE in_out_manager ! I/O manager
USE lib_mpp ! MPP library
USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
USE timing ! Timing
IMPLICIT NONE
PRIVATE
PUBLIC ice_alb_init ! called in icestp
PUBLIC ice_alb ! called in icesbc.F90 and iceupdate.F90
REAL(wp), PUBLIC, PARAMETER :: rn_alb_oce = 0.066_wp !: ocean or lead albedo (Pegau and Paulson, Ann. Glac. 2001)
!
! !!* albedo namelist (namalb)
REAL(wp) :: rn_alb_sdry ! dry snow albedo
REAL(wp) :: rn_alb_smlt ! melting snow albedo
REAL(wp) :: rn_alb_idry ! dry ice albedo
REAL(wp) :: rn_alb_imlt ! bare puddled ice albedo
REAL(wp) :: rn_alb_dpnd ! ponded ice albedo
REAL(wp) :: rn_alb_hpiv ! pivotal ice thickness in meters (above which albedo is constant)
!! * Substitutions
# include "do_loop_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/ICE 4.0 , NEMO Consortium (2018)
!! $Id: icealb.F90 15549 2021-11-28 20:00:36Z clem $
!! Software governed by the CeCILL license (see ./LICENSE)
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE ice_alb( pt_su, ph_ice, ph_snw, ld_pnd_alb, pafrac_pnd, ph_pnd, pcloud_fra, palb_ice )
!!----------------------------------------------------------------------
!! *** ROUTINE ice_alb ***
!!
!! ** Purpose : Computation of the albedo of the snow/ice system
!!
!! ** Method : The scheme is "home made" (for cloudy skies) and based on Brandt et al. (J. Climate 2005)
!! and Grenfell & Perovich (JGR 2004)
!! 1) Albedo dependency on ice thickness follows the findings from Brandt et al (2005)
!! which are an update of Allison et al. (JGR 1993) ; Brandt et al. 1999
!! 0-5cm : linear function of ice thickness
!! 5-100cm: log function of ice thickness
!! > 100cm: constant
!! 2) Albedo dependency on snow thickness follows the findings from Grenfell & Perovich (2004)
!! i.e. it increases as -EXP(-snw_thick/0.02) during freezing and -EXP(-snw_thick/0.03) during melting
!! 3) Albedo dependency on clouds is speculated from measurements of Grenfell and Perovich (2004)
!! i.e. cloudy-clear albedo depend on cloudy albedo following a 2d order polynomial law
!! 4) The needed 4 parameters are: dry and melting snow, freezing ice and bare puddled ice
!!
!! compilation of values from literature (reference overcast sky values)
!! rn_alb_sdry = 0.85 ! dry snow
!! rn_alb_smlt = 0.75 ! melting snow
!! rn_alb_idry = 0.60 ! bare frozen ice
!! rn_alb_imlt = 0.50 ! bare puddled ice albedo
!! rn_alb_dpnd = 0.36 ! ponded-ice overcast albedo (Lecomte et al, 2015)
!! ! early melt pnds 0.27, late melt ponds 0.14 Grenfell & Perovich
!! Perovich et al 2002 (Sheba) => the only dataset for which all types of ice/snow were retrieved
!! rn_alb_sdry = 0.85 ! dry snow
!! rn_alb_smlt = 0.72 ! melting snow
!! rn_alb_idry = 0.65 ! bare frozen ice
!! Brandt et al 2005 (East Antarctica)
!! rn_alb_sdry = 0.87 ! dry snow
!! rn_alb_smlt = 0.82 ! melting snow
!! rn_alb_idry = 0.54 ! bare frozen ice
!!
!! ** Note : The old parameterization from Shine & Henderson-Sellers (not here anymore) presented several misconstructions
!! 1) ice albedo when ice thick. tends to 0 is different than ocean albedo
!! 2) for small ice thick. covered with some snow (<3cm?), albedo is larger
!! under melting conditions than under freezing conditions
!! 3) the evolution of ice albedo as a function of ice thickness shows
!! 3 sharp inflexion points (at 5cm, 100cm and 150cm) that look highly unrealistic
!!
!! References : Shine & Henderson-Sellers 1985, JGR, 90(D1), 2243-2250.
!! Brandt et al. 2005, J. Climate, vol 18
!! Grenfell & Perovich 2004, JGR, vol 109
!!----------------------------------------------------------------------
REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: pt_su ! ice surface temperature (Kelvin)
REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: ph_ice ! sea-ice thickness
REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: ph_snw ! snow depth
LOGICAL , INTENT(in ) :: ld_pnd_alb ! effect of melt ponds on albedo
REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: pafrac_pnd ! melt pond relative fraction (per unit ice area)
REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: ph_pnd ! melt pond depth
REAL(wp), INTENT(in ), DIMENSION(:,:) :: pcloud_fra ! cloud fraction
REAL(wp), INTENT( out), DIMENSION(:,:,:) :: palb_ice ! albedo of ice
!
REAL(wp), DIMENSION(jpi,jpj,jpl) :: za_s_fra ! ice fraction covered by snow
INTEGER :: ji, jj, jl ! dummy loop indices
REAL(wp) :: z1_c1, z1_c2,z1_c3, z1_c4 ! local scalar
REAL(wp) :: z1_href_pnd ! inverse of the characteristic length scale (Lecomte et al. 2015)
REAL(wp) :: zalb_pnd, zafrac_pnd ! ponded sea ice albedo & relative pound fraction
REAL(wp) :: zalb_ice, zafrac_ice ! bare sea ice albedo & relative ice fraction
REAL(wp) :: zalb_snw, zafrac_snw ! snow-covered sea ice albedo & relative snow fraction
REAL(wp) :: zalb_cs, zalb_os ! albedo of ice under clear/overcast sky
!!---------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('icealb')
!
z1_href_pnd = 1._wp / 0.05_wp
z1_c1 = 1._wp / ( LOG(rn_alb_hpiv) - LOG(0.05_wp) )
z1_c2 = 1._wp / 0.05_wp
z1_c3 = 1._wp / 0.02_wp
z1_c4 = 1._wp / 0.03_wp
!
CALL ice_var_snwfra( ph_snw, za_s_fra ) ! calculate ice fraction covered by snow
!
DO jl = 1, jpl
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) ! palb_ice used over the full domain in icesbc
!
!---------------------------------------------!
!--- Specific snow, ice and pond fractions ---!
!---------------------------------------------!
zafrac_snw = za_s_fra(ji,jj,jl)
IF( ld_pnd_alb ) THEN
zafrac_pnd = MIN( pafrac_pnd(ji,jj,jl), 1._wp - zafrac_snw ) ! make sure (a_ip_eff + a_s_fra) <= 1
ELSE
zafrac_pnd = 0._wp
ENDIF
zafrac_ice = MAX( 0._wp, 1._wp - zafrac_pnd - zafrac_snw ) ! max for roundoff errors
!
!---------------!
!--- Albedos ---!
!---------------!
! !--- Bare ice albedo (for hi > 100cm)
IF( ld_pnd_alb ) THEN
zalb_ice = rn_alb_idry
ELSE
IF( ph_snw(ji,jj,jl) == 0._wp .AND. pt_su(ji,jj,jl) >= rt0 ) THEN ; zalb_ice = rn_alb_imlt
ELSE ; zalb_ice = rn_alb_idry ; ENDIF
ENDIF
! !--- Bare ice albedo (for hi < 100cm)
IF( 0.05 < ph_ice(ji,jj,jl) .AND. ph_ice(ji,jj,jl) <= rn_alb_hpiv ) THEN ! 5cm < hi < 100cm
zalb_ice = zalb_ice + ( 0.18_wp - zalb_ice ) * z1_c1 * ( LOG(rn_alb_hpiv) - LOG(ph_ice(ji,jj,jl)) )
ELSEIF( ph_ice(ji,jj,jl) <= 0.05_wp ) THEN ! 0cm < hi < 5cm
zalb_ice = rn_alb_oce + ( 0.18_wp - rn_alb_oce ) * z1_c2 * ph_ice(ji,jj,jl)
ENDIF
!
! !--- Snow-covered ice albedo (freezing, melting cases)
IF( pt_su(ji,jj,jl) < rt0 ) THEN
zalb_snw = rn_alb_sdry - ( rn_alb_sdry - zalb_ice ) * EXP( - ph_snw(ji,jj,jl) * z1_c3 )
ELSE
zalb_snw = rn_alb_smlt - ( rn_alb_smlt - zalb_ice ) * EXP( - ph_snw(ji,jj,jl) * z1_c4 )
ENDIF
! !--- Ponded ice albedo
zalb_pnd = rn_alb_dpnd - ( rn_alb_dpnd - zalb_ice ) * EXP( - ph_pnd(ji,jj,jl) * z1_href_pnd )
!
! !--- Surface albedo is weighted mean of snow, ponds and bare ice contributions
zalb_os = ( zafrac_snw * zalb_snw + zafrac_pnd * zalb_pnd + zafrac_ice * zalb_ice ) * tmask(ji,jj,1)
!
zalb_cs = zalb_os - ( - 0.1010_wp * zalb_os * zalb_os &
& + 0.1933_wp * zalb_os - 0.0148_wp ) * tmask(ji,jj,1)
!
! albedo depends on cloud fraction because of non-linear spectral effects
palb_ice(ji,jj,jl) = ( 1._wp - pcloud_fra(ji,jj) ) * zalb_cs + pcloud_fra(ji,jj) * zalb_os
END_2D
END DO
!
!
IF( ln_timing ) CALL timing_stop('icealb')
!
END SUBROUTINE ice_alb
SUBROUTINE ice_alb_init
!!----------------------------------------------------------------------
!! *** ROUTINE alb_init ***
!!
!! ** Purpose : initializations for the albedo parameters
!!
!! ** Method : Read the namelist namalb
!!----------------------------------------------------------------------
INTEGER :: ios ! Local integer output status for namelist read
!!
NAMELIST/namalb/ rn_alb_sdry, rn_alb_smlt, rn_alb_idry, rn_alb_imlt, rn_alb_dpnd, rn_alb_hpiv
!!----------------------------------------------------------------------
!
READ ( numnam_ice_ref, namalb, IOSTAT = ios, ERR = 901)
901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namalb in reference namelist' )
READ ( numnam_ice_cfg, namalb, IOSTAT = ios, ERR = 902 )
902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namalb in configuration namelist' )
IF(lwm) WRITE( numoni, namalb )
!
IF(lwp) THEN ! Control print
WRITE(numout,*)
WRITE(numout,*) 'ice_alb_init: set albedo parameters'
WRITE(numout,*) '~~~~~~~~~~~~'
WRITE(numout,*) ' Namelist namalb:'
WRITE(numout,*) ' albedo of dry snow rn_alb_sdry = ', rn_alb_sdry
WRITE(numout,*) ' albedo of melting snow rn_alb_smlt = ', rn_alb_smlt
WRITE(numout,*) ' albedo of dry ice rn_alb_idry = ', rn_alb_idry
WRITE(numout,*) ' albedo of bare puddled ice rn_alb_imlt = ', rn_alb_imlt
WRITE(numout,*) ' albedo of ponded ice rn_alb_dpnd = ', rn_alb_dpnd
WRITE(numout,*) ' pivotal ice thickness (m) rn_alb_hpiv = ', rn_alb_hpiv
ENDIF
!
END SUBROUTINE ice_alb_init
#else
!!----------------------------------------------------------------------
!! Default option Dummy module NO SI3 sea-ice model
!!----------------------------------------------------------------------
#endif
!!======================================================================
END MODULE icealb