1 |
|
|
! |
2 |
|
|
MODULE surf_landice_mod |
3 |
|
|
|
4 |
|
|
IMPLICIT NONE |
5 |
|
|
|
6 |
|
|
CONTAINS |
7 |
|
|
! |
8 |
|
|
!**************************************************************************************** |
9 |
|
|
! |
10 |
|
1152 |
SUBROUTINE surf_landice(itime, dtime, knon, knindex, & |
11 |
|
|
rlon, rlat, debut, lafin, & |
12 |
|
|
rmu0, lwdownm, albedo, pphi1, & |
13 |
|
|
swnet, lwnet, tsurf, p1lay, & |
14 |
|
|
cdragh, cdragm, precip_rain, precip_snow, precip_bs, temp_air, spechum, & |
15 |
|
|
AcoefH, AcoefQ, BcoefH, BcoefQ, & |
16 |
|
|
AcoefU, AcoefV, BcoefU, BcoefV, & |
17 |
|
|
AcoefQBS, BcoefQBS, & |
18 |
|
288 |
ps, u1, v1, gustiness, rugoro, pctsrf, & |
19 |
|
|
snow, qsurf, qsol, qbs1, agesno, & |
20 |
|
288 |
tsoil, z0m, z0h, SFRWL, alb_dir, alb_dif, evap, fluxsens, fluxlat, fluxbs, & |
21 |
|
|
tsurf_new, dflux_s, dflux_l, & |
22 |
|
|
alt, slope, cloudf, & |
23 |
|
|
snowhgt, qsnow, to_ice, sissnow, & |
24 |
|
288 |
alb3, runoff, & |
25 |
|
|
flux_u1, flux_v1) |
26 |
|
|
|
27 |
|
|
USE dimphy |
28 |
|
|
USE geometry_mod, ONLY : longitude,latitude |
29 |
|
|
USE surface_data, ONLY : type_ocean, calice, calsno, landice_opt, iflag_albcalc |
30 |
|
|
USE fonte_neige_mod, ONLY : fonte_neige,run_off_lic,fqcalving_global,ffonte_global,fqfonte_global,runofflic_global |
31 |
|
|
USE cpl_mod, ONLY : cpl_send_landice_fields |
32 |
|
|
USE calcul_fluxs_mod |
33 |
|
|
USE phys_local_var_mod, ONLY : zxrhoslic, zxustartlic |
34 |
|
|
USE phys_output_var_mod, ONLY : snow_o,zfra_o |
35 |
|
|
!FC |
36 |
|
|
USE ioipsl_getin_p_mod, ONLY : getin_p |
37 |
|
|
USE blowing_snow_ini_mod, ONLY : zeta_bs, pbst_bs, prt_bs, iflag_saltation_bs |
38 |
|
|
|
39 |
|
|
#ifdef CPP_INLANDSIS |
40 |
|
|
USE surf_inlandsis_mod, ONLY : surf_inlandsis |
41 |
|
|
#endif |
42 |
|
|
|
43 |
|
|
USE indice_sol_mod |
44 |
|
|
|
45 |
|
|
! INCLUDE "indicesol.h" |
46 |
|
|
INCLUDE "dimsoil.h" |
47 |
|
|
INCLUDE "YOMCST.h" |
48 |
|
|
INCLUDE "clesphys.h" |
49 |
|
|
|
50 |
|
|
! Input variables |
51 |
|
|
!**************************************************************************************** |
52 |
|
|
INTEGER, INTENT(IN) :: itime, knon |
53 |
|
|
INTEGER, DIMENSION(klon), INTENT(in) :: knindex |
54 |
|
|
REAL, INTENT(in) :: dtime |
55 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: swnet ! net shortwave radiance |
56 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: lwnet ! net longwave radiance |
57 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: tsurf |
58 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: p1lay |
59 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: cdragh, cdragm |
60 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: precip_rain, precip_snow, precip_bs |
61 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: temp_air, spechum |
62 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: AcoefH, AcoefQ |
63 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: BcoefH, BcoefQ |
64 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: AcoefU, AcoefV, BcoefU, BcoefV |
65 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: AcoefQBS, BcoefQBS |
66 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: ps |
67 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: u1, v1, gustiness, qbs1 |
68 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: rugoro |
69 |
|
|
REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf |
70 |
|
|
|
71 |
|
|
LOGICAL, INTENT(IN) :: debut !true if first step |
72 |
|
|
LOGICAL, INTENT(IN) :: lafin !true if last step |
73 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: rlon, rlat |
74 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: rmu0 |
75 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: lwdownm !ylwdown |
76 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: albedo !mean albedo |
77 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: pphi1 |
78 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: alt !mean altitude of the grid box |
79 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: slope !mean slope in grid box |
80 |
|
|
REAL, DIMENSION(klon), INTENT(IN) :: cloudf !total cloud fraction |
81 |
|
|
|
82 |
|
|
! In/Output variables |
83 |
|
|
!**************************************************************************************** |
84 |
|
|
REAL, DIMENSION(klon), INTENT(INOUT) :: snow, qsol |
85 |
|
|
REAL, DIMENSION(klon), INTENT(INOUT) :: agesno |
86 |
|
|
REAL, DIMENSION(klon, nsoilmx), INTENT(INOUT) :: tsoil |
87 |
|
|
|
88 |
|
|
! Output variables |
89 |
|
|
!**************************************************************************************** |
90 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: qsurf |
91 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: z0m, z0h |
92 |
|
|
!albedo SB >>> |
93 |
|
|
! REAL, DIMENSION(klon), INTENT(OUT) :: alb1 ! new albedo in visible SW interval |
94 |
|
|
! REAL, DIMENSION(klon), INTENT(OUT) :: alb2 ! new albedo in near IR interval |
95 |
|
|
REAL, DIMENSION(6), INTENT(IN) :: SFRWL |
96 |
|
|
REAL, DIMENSION(klon,nsw), INTENT(OUT) :: alb_dir,alb_dif |
97 |
|
|
!albedo SB <<< |
98 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: evap, fluxsens, fluxlat |
99 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: fluxbs |
100 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: tsurf_new |
101 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: dflux_s, dflux_l |
102 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: flux_u1, flux_v1 |
103 |
|
|
|
104 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: alb3 |
105 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: qsnow !column water in snow [kg/m2] |
106 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: snowhgt !Snow height (m) |
107 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: to_ice |
108 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: sissnow |
109 |
|
|
REAL, DIMENSION(klon), INTENT(OUT) :: runoff !Land ice runoff |
110 |
|
|
|
111 |
|
|
|
112 |
|
|
! Local variables |
113 |
|
|
!**************************************************************************************** |
114 |
|
576 |
REAL, DIMENSION(klon) :: soilcap, soilflux |
115 |
|
576 |
REAL, DIMENSION(klon) :: cal, beta, dif_grnd |
116 |
|
576 |
REAL, DIMENSION(klon) :: zfra, alb_neig |
117 |
|
576 |
REAL, DIMENSION(klon) :: radsol |
118 |
|
576 |
REAL, DIMENSION(klon) :: u0, v0, u1_lay, v1_lay, ustar |
119 |
|
|
INTEGER :: i,j,nt |
120 |
|
|
REAL, DIMENSION(klon) :: fqfonte,ffonte |
121 |
|
576 |
REAL, DIMENSION(klon) :: run_off_lic_frac |
122 |
|
|
REAL, DIMENSION(klon) :: emis_new !Emissivity |
123 |
|
|
REAL, DIMENSION(klon) :: swdown,lwdown |
124 |
|
|
REAL, DIMENSION(klon) :: precip_snow_adv, snow_adv !Snow Drift precip./advection (not used in inlandsis) |
125 |
|
|
REAL, DIMENSION(klon) :: erod !erosion of surface snow (flux, kg/m2/s like evap) |
126 |
|
|
REAL, DIMENSION(klon) :: zsl_height, wind_velo !surface layer height, wind spd |
127 |
|
|
REAL, DIMENSION(klon) :: dens_air, snow_cont_air !air density; snow content air |
128 |
|
|
REAL, DIMENSION(klon) :: alb_soil !albedo of underlying ice |
129 |
|
|
REAL, DIMENSION(klon) :: pexner !Exner potential |
130 |
|
|
REAL :: pref |
131 |
|
|
REAL, DIMENSION(klon,nsoilmx) :: tsoil0 !modif |
132 |
|
|
REAL :: dtis ! subtimestep |
133 |
|
|
LOGICAL :: debut_is, lafin_is ! debut and lafin for inlandsis |
134 |
|
|
|
135 |
|
|
CHARACTER (len = 20) :: modname = 'surf_landice' |
136 |
|
|
CHARACTER (len = 80) :: abort_message |
137 |
|
|
|
138 |
|
|
|
139 |
|
576 |
REAL,DIMENSION(klon) :: alb1,alb2 |
140 |
|
576 |
REAL,DIMENSION(klon) :: precip_totsnow, evap_totsnow |
141 |
|
576 |
REAL, DIMENSION (klon,6) :: alb6 |
142 |
|
|
REAL :: rho0, rhoice, ustart0, hsalt, esalt, rhod |
143 |
|
|
REAL :: lambdasalt,fluxsalt, csalt, nunu, aa, bb, cc |
144 |
|
|
REAL :: tau_dens, tau_dens0, tau_densmin, rhomax, rhohard |
145 |
|
576 |
REAL, DIMENSION(klon) :: ws1, rhos, ustart, qsalt |
146 |
|
|
! End definition |
147 |
|
|
!**************************************************************************************** |
148 |
|
|
!FC |
149 |
|
|
!FC |
150 |
|
|
REAL,SAVE :: alb_vis_sno_lic |
151 |
|
|
!$OMP THREADPRIVATE(alb_vis_sno_lic) |
152 |
|
|
REAL,SAVE :: alb_nir_sno_lic |
153 |
|
|
!$OMP THREADPRIVATE(alb_nir_sno_lic) |
154 |
|
|
LOGICAL, SAVE :: firstcall = .TRUE. |
155 |
|
|
!$OMP THREADPRIVATE(firstcall) |
156 |
|
|
|
157 |
|
|
|
158 |
|
|
!FC firtscall initializations |
159 |
|
|
!****************************************************************************************** |
160 |
✓✓ |
288 |
IF (firstcall) THEN |
161 |
|
1 |
alb_vis_sno_lic=0.77 |
162 |
|
1 |
CALL getin_p('alb_vis_sno_lic',alb_vis_sno_lic) |
163 |
|
1 |
PRINT*, 'alb_vis_sno_lic',alb_vis_sno_lic |
164 |
|
1 |
alb_nir_sno_lic=0.77 |
165 |
|
1 |
CALL getin_p('alb_nir_sno_lic',alb_nir_sno_lic) |
166 |
|
1 |
PRINT*, 'alb_nir_sno_lic',alb_nir_sno_lic |
167 |
|
|
|
168 |
|
1 |
firstcall=.false. |
169 |
|
|
ENDIF |
170 |
|
|
!****************************************************************************************** |
171 |
|
|
|
172 |
|
|
! Initialize output variables |
173 |
✓✓ |
286560 |
alb3(:) = 999999. |
174 |
✓✓ |
286560 |
alb2(:) = 999999. |
175 |
✓✓ |
286560 |
alb1(:) = 999999. |
176 |
✓✓ |
286560 |
fluxbs(:)=0. |
177 |
✓✓ |
286560 |
runoff(:) = 0. |
178 |
|
|
!**************************************************************************************** |
179 |
|
|
! Calculate total absorbed radiance at surface |
180 |
|
|
! |
181 |
|
|
!**************************************************************************************** |
182 |
✓✓ |
286560 |
radsol(:) = 0.0 |
183 |
✓✓ |
44064 |
radsol(1:knon) = swnet(1:knon) + lwnet(1:knon) |
184 |
|
|
|
185 |
|
|
!**************************************************************************************** |
186 |
|
|
|
187 |
|
|
!**************************************************************************************** |
188 |
|
|
! landice_opt = 0 : soil_model, calcul_flux, fonte_neige, ... |
189 |
|
|
! landice_opt = 1 : prepare and call INterace Lmdz SISvat (INLANDSIS) |
190 |
|
|
!**************************************************************************************** |
191 |
|
|
|
192 |
|
|
|
193 |
✗✓ |
288 |
IF (landice_opt .EQ. 1) THEN |
194 |
|
|
|
195 |
|
|
!**************************************************************************************** |
196 |
|
|
! CALL to INLANDSIS interface |
197 |
|
|
!**************************************************************************************** |
198 |
|
|
#ifdef CPP_INLANDSIS |
199 |
|
|
|
200 |
|
|
debut_is=debut |
201 |
|
|
lafin_is=.false. |
202 |
|
|
! Suppose zero surface speed |
203 |
|
|
u0(:) = 0.0 |
204 |
|
|
v0(:) = 0.0 |
205 |
|
|
|
206 |
|
|
|
207 |
|
|
CALL calcul_flux_wind(knon, dtime, & |
208 |
|
|
u0, v0, u1, v1, gustiness, cdragm, & |
209 |
|
|
AcoefU, AcoefV, BcoefU, BcoefV, & |
210 |
|
|
p1lay, temp_air, & |
211 |
|
|
flux_u1, flux_v1) |
212 |
|
|
|
213 |
|
|
|
214 |
|
|
! Set constants and compute some input for SISVAT |
215 |
|
|
! = 1000 hPa |
216 |
|
|
! and calculate incoming flux for SW and LW interval: swdown, lwdown |
217 |
|
|
swdown(:) = 0.0 |
218 |
|
|
lwdown(:) = 0.0 |
219 |
|
|
snow_cont_air(:) = 0. ! the snow content in air is not a prognostic variable of the model |
220 |
|
|
alb_soil(:) = 0.4 ! before albedo(:) but here it is the ice albedo that we have to set |
221 |
|
|
ustar(:) = 0. |
222 |
|
|
pref = 100000. |
223 |
|
|
DO i = 1, knon |
224 |
|
|
swdown(i) = swnet(i)/(1-albedo(i)) |
225 |
|
|
lwdown(i) = lwdownm(i) |
226 |
|
|
wind_velo(i) = u1(i)**2 + v1(i)**2 |
227 |
|
|
wind_velo(i) = wind_velo(i)**0.5 |
228 |
|
|
pexner(i) = (p1lay(i)/pref)**(RD/RCPD) |
229 |
|
|
dens_air(i) = p1lay(i)/RD/temp_air(i) ! dry air density |
230 |
|
|
zsl_height(i) = pphi1(i)/RG |
231 |
|
|
tsoil0(i,:) = tsoil(i,:) |
232 |
|
|
ustar(i)= (cdragm(i)*(wind_velo(i)**2))**0.5 |
233 |
|
|
END DO |
234 |
|
|
|
235 |
|
|
|
236 |
|
|
|
237 |
|
|
dtis=dtime |
238 |
|
|
|
239 |
|
|
IF (lafin) THEN |
240 |
|
|
lafin_is=.true. |
241 |
|
|
END IF |
242 |
|
|
|
243 |
|
|
CALL surf_inlandsis(knon, rlon, rlat, knindex, itime, dtis, debut_is, lafin_is,& |
244 |
|
|
rmu0, swdown, lwdown, albedo, pexner, ps, p1lay, precip_rain, precip_snow, & |
245 |
|
|
zsl_height, wind_velo, ustar, temp_air, dens_air, spechum, tsurf,& |
246 |
|
|
rugoro, snow_cont_air, alb_soil, alt, slope, cloudf, & |
247 |
|
|
radsol, qsol, tsoil0, snow, zfra, snowhgt, qsnow, to_ice, sissnow,agesno, & |
248 |
|
|
AcoefH, AcoefQ, BcoefH, BcoefQ, cdragm, cdragh, & |
249 |
|
|
run_off_lic, fqfonte, ffonte, evap, erod, fluxsens, fluxlat,dflux_s, dflux_l, & |
250 |
|
|
tsurf_new, alb1, alb2, alb3, alb6, & |
251 |
|
|
emis_new, z0m, z0h, qsurf) |
252 |
|
|
|
253 |
|
|
debut_is=.false. |
254 |
|
|
|
255 |
|
|
|
256 |
|
|
! Treatment of snow melting and calving |
257 |
|
|
|
258 |
|
|
! for consistency with standard LMDZ, add calving to run_off_lic |
259 |
|
|
run_off_lic(:)=run_off_lic(:) + to_ice(:) |
260 |
|
|
|
261 |
|
|
DO i = 1, knon |
262 |
|
|
ffonte_global(knindex(i),is_lic) = ffonte(i) |
263 |
|
|
fqfonte_global(knindex(i),is_lic) = fqfonte(i)! net melting= melting - refreezing |
264 |
|
|
fqcalving_global(knindex(i),is_lic) = to_ice(i) ! flux |
265 |
|
|
runofflic_global(knindex(i)) = run_off_lic(i) |
266 |
|
|
ENDDO |
267 |
|
|
! Here, we assume that the calving term is equal to the to_ice term |
268 |
|
|
! (no ice accumulation) |
269 |
|
|
|
270 |
|
|
|
271 |
|
|
#else |
272 |
|
|
abort_message='Pb de coherence: landice_opt = 1 mais CPP_INLANDSIS = .false.' |
273 |
|
|
CALL abort_physic(modname,abort_message,1) |
274 |
|
|
#endif |
275 |
|
|
|
276 |
|
|
|
277 |
|
|
ELSE |
278 |
|
|
|
279 |
|
|
!**************************************************************************************** |
280 |
|
|
! Soil calculations |
281 |
|
|
! |
282 |
|
|
!**************************************************************************************** |
283 |
|
|
|
284 |
|
|
! EV: use calbeta |
285 |
|
288 |
CALL calbeta(dtime, is_lic, knon, snow, qsol, beta, cal, dif_grnd) |
286 |
|
|
|
287 |
|
|
|
288 |
|
|
! use soil model and recalculate properly cal |
289 |
✓✗ |
288 |
IF (soil_model) THEN |
290 |
|
|
CALL soil(dtime, is_lic, knon, snow, tsurf, qsol, & |
291 |
✓✗✓✓ ✓✓ |
87840 |
& longitude(knindex(1:knon)), latitude(knindex(1:knon)), tsoil, soilcap, soilflux) |
292 |
✓✓ |
44064 |
cal(1:knon) = RCPD / soilcap(1:knon) |
293 |
✓✓ |
44064 |
radsol(1:knon) = radsol(1:knon) + soilflux(1:knon) |
294 |
|
|
ELSE |
295 |
|
|
cal = RCPD * calice |
296 |
|
|
WHERE (snow > 0.0) cal = RCPD * calsno |
297 |
|
|
ENDIF |
298 |
|
|
|
299 |
|
|
|
300 |
|
|
!**************************************************************************************** |
301 |
|
|
! Calulate fluxes |
302 |
|
|
! |
303 |
|
|
!**************************************************************************************** |
304 |
|
|
! beta(:) = 1.0 |
305 |
|
|
! dif_grnd(:) = 0.0 |
306 |
|
|
|
307 |
|
|
! Suppose zero surface speed |
308 |
✓✓ |
286560 |
u0(:)=0.0 |
309 |
✓✓ |
286560 |
v0(:)=0.0 |
310 |
✓✓ |
286560 |
u1_lay(:) = u1(:) - u0(:) |
311 |
✓✓ |
286560 |
v1_lay(:) = v1(:) - v0(:) |
312 |
|
|
|
313 |
|
|
CALL calcul_fluxs(knon, is_lic, dtime, & |
314 |
|
|
tsurf, p1lay, cal, beta, cdragh, cdragh, ps, & |
315 |
|
|
precip_rain, precip_snow, snow, qsurf, & |
316 |
|
|
radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, gustiness, & |
317 |
|
|
1.,AcoefH, AcoefQ, BcoefH, BcoefQ, & |
318 |
|
288 |
tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l) |
319 |
|
|
|
320 |
|
|
CALL calcul_flux_wind(knon, dtime, & |
321 |
|
|
u0, v0, u1, v1, gustiness, cdragm, & |
322 |
|
|
AcoefU, AcoefV, BcoefU, BcoefV, & |
323 |
|
|
p1lay, temp_air, & |
324 |
|
288 |
flux_u1, flux_v1) |
325 |
|
|
|
326 |
|
|
|
327 |
|
|
!**************************************************************************************** |
328 |
|
|
! Calculate albedo |
329 |
|
|
! |
330 |
|
|
!**************************************************************************************** |
331 |
|
|
|
332 |
|
288 |
CALL albsno(klon,knon,dtime,agesno(:),alb_neig(:), precip_snow(:)) |
333 |
|
|
|
334 |
|
|
|
335 |
|
|
! EV: following lines are obsolete since we set alb1 and alb2 to constant values |
336 |
|
|
! I therefore comment them |
337 |
|
|
! alb1(1:knon) = alb_neig(1:knon)*zfra(1:knon) + & |
338 |
|
|
! 0.6 * (1.0-zfra(1:knon)) |
339 |
|
|
! |
340 |
|
|
!IM: plusieurs choix/tests sur l'albedo des "glaciers continentaux" |
341 |
|
|
! alb1(1 : knon) = 0.6 !IM cf FH/GK |
342 |
|
|
! alb1(1 : knon) = 0.82 |
343 |
|
|
! alb1(1 : knon) = 0.77 !211003 Ksta0.77 |
344 |
|
|
! alb1(1 : knon) = 0.8 !KstaTER0.8 & LMD_ARMIP5 |
345 |
|
|
!IM: KstaTER0.77 & LMD_ARMIP6 |
346 |
|
|
|
347 |
|
|
! Attantion: alb1 and alb2 are not the same! |
348 |
✓✓ |
44064 |
alb1(1:knon) = alb_vis_sno_lic |
349 |
✓✓ |
44064 |
alb2(1:knon) = alb_nir_sno_lic |
350 |
|
|
|
351 |
|
|
|
352 |
|
|
!**************************************************************************************** |
353 |
|
|
! Rugosity |
354 |
|
|
! |
355 |
|
|
!**************************************************************************************** |
356 |
✓✓ |
286560 |
z0m = z0m_landice |
357 |
✓✓ |
286560 |
z0h = z0h_landice |
358 |
|
|
!z0m = SQRT(z0m**2+rugoro**2) |
359 |
|
|
|
360 |
|
|
|
361 |
|
|
|
362 |
|
|
! Simple blowing snow param |
363 |
✗✓ |
288 |
if (ok_bs) then |
364 |
|
|
ustart0 = 0.211 |
365 |
|
|
rhoice = 920.0 |
366 |
|
|
rho0 = 200.0 |
367 |
|
|
rhomax=450.0 |
368 |
|
|
rhohard=400.0 |
369 |
|
|
tau_dens0=86400.0*10. ! 10 days by default, in s |
370 |
|
|
tau_densmin=86400.0 ! 1 days according to in situ obs by C. Amory |
371 |
|
|
|
372 |
|
|
! computation of threshold friction velocity |
373 |
|
|
! which depends on surface snow density |
374 |
|
|
do i = 1, knon |
375 |
|
|
! estimation of snow density |
376 |
|
|
! snow density increases with snow age and |
377 |
|
|
! increases even faster in case of sedimentation of blowing snow |
378 |
|
|
tau_dens=max(tau_densmin, tau_dens0*exp(-abs(precip_bs(i))/pbst_bs-abs(precip_rain(i))/prt_bs)) |
379 |
|
|
rhos(i)=rho0+(rhohard-rho0)*(1.-exp(-agesno(i)*86400.0/tau_dens)) |
380 |
|
|
! blowing snow flux formula used in MAR |
381 |
|
|
ws1(i)=(u1(i)**2+v1(i)**2)**0.5 |
382 |
|
|
ustar(i)=(cdragm(i)*(u1(i)**2+v1(i)**2))**0.5 |
383 |
|
|
ustart(i)=ustart0*exp(max(rhoice/rho0-rhoice/rhos(i),0.))*exp(max(0.,rhos(i)-rhomax)) |
384 |
|
|
! we have multiplied by exp to prevent erosion when rhos>rhomax (usefull till |
385 |
|
|
! rhohard<450) |
386 |
|
|
enddo |
387 |
|
|
|
388 |
|
|
! computation of qbs at the top of the saltation layer |
389 |
|
|
! two formulations possible |
390 |
|
|
! pay attention that qbs is a mixing ratio and has to be converted |
391 |
|
|
! to specific content |
392 |
|
|
|
393 |
|
|
if (iflag_saltation_bs .eq. 1) then |
394 |
|
|
! expression from CRYOWRF (Sharma et al. 2022) |
395 |
|
|
aa=2.6 |
396 |
|
|
bb=2.5 |
397 |
|
|
cc=2.0 |
398 |
|
|
lambdasalt=0.45 |
399 |
|
|
do i =1, knon |
400 |
|
|
rhod=p1lay(i)/RD/temp_air(i) |
401 |
|
|
nunu=max(ustar(i)/ustart(i),1.e-3) |
402 |
|
|
fluxsalt=rhod/RG*(ustar(i)**3)*(1.-nunu**(-2)) * & |
403 |
|
|
(aa+bb*nunu**(-2)+cc*nunu**(-1)) |
404 |
|
|
csalt=fluxsalt/(2.8*ustart(i)) |
405 |
|
|
hsalt=0.08436*ustar(i)**1.27 |
406 |
|
|
qsalt(i)=1./rhod*csalt*lambdasalt*RG/(max(ustar(i)**2,1E-6)) & |
407 |
|
|
* exp(-lambdasalt*RG*hsalt/max(ustar(i)**2,1E-6)) |
408 |
|
|
qsalt(i)=max(qsalt(i),0.) |
409 |
|
|
enddo |
410 |
|
|
|
411 |
|
|
|
412 |
|
|
else |
413 |
|
|
! default formulation from MAR model (Amory et al. 2021, Gallee et al. 2001) |
414 |
|
|
do i=1, knon |
415 |
|
|
esalt=1./(3.25*max(ustar(i),0.001)) |
416 |
|
|
hsalt=0.08436*ustar(i)**1.27 |
417 |
|
|
qsalt(i)=(max(ustar(i)**2-ustart(i)**2,0.))/(RG*hsalt)*esalt |
418 |
|
|
!ep=qsalt*cdragm(i)*sqrt(u1(i)**2+v1(i)**2) |
419 |
|
|
enddo |
420 |
|
|
endif |
421 |
|
|
|
422 |
|
|
! calculation of erosion (emission flux towards the first atmospheric level) |
423 |
|
|
! consistent with implicit resolution of turbulent mixing equation |
424 |
|
|
do i=1, knon |
425 |
|
|
rhod=p1lay(i)/RD/temp_air(i) |
426 |
|
|
fluxbs(i)=rhod*ws1(i)*cdragm(i)*zeta_bs*(AcoefQBS(i)-qsalt(i)) & |
427 |
|
|
/ (1.-rhod*ws1(i)*zeta_bs*cdragm(i)*BcoefQBS(i)*dtime) |
428 |
|
|
!fluxbs(i)= zeta_bs*rhod*ws1(i)*cdragm(i)*(qbs1(i)-qsalt(i)) |
429 |
|
|
enddo |
430 |
|
|
|
431 |
|
|
! for outputs |
432 |
|
|
do j = 1, knon |
433 |
|
|
i = knindex(j) |
434 |
|
|
zxustartlic(i) = ustart(j) |
435 |
|
|
zxrhoslic(i) = rhos(j) |
436 |
|
|
enddo |
437 |
|
|
|
438 |
|
|
endif |
439 |
|
|
|
440 |
|
|
|
441 |
|
|
|
442 |
|
|
!**************************************************************************************** |
443 |
|
|
! Calculate surface snow amount |
444 |
|
|
! |
445 |
|
|
!**************************************************************************************** |
446 |
✗✓ |
288 |
IF (ok_bs) THEN |
447 |
|
|
precip_totsnow(:)=precip_snow(:)+precip_bs(:) |
448 |
|
|
evap_totsnow(:)=evap(:)-fluxbs(:) ! flux bs is positive towards the surface (snow erosion) |
449 |
|
|
ELSE |
450 |
✓✓ |
286560 |
precip_totsnow(:)=precip_snow(:) |
451 |
✓✓ |
286560 |
evap_totsnow(:)=evap(:) |
452 |
|
|
ENDIF |
453 |
|
|
|
454 |
|
|
CALL fonte_neige(knon, is_lic, knindex, dtime, & |
455 |
|
|
tsurf, precip_rain, precip_totsnow, & |
456 |
|
288 |
snow, qsol, tsurf_new, evap_totsnow) |
457 |
|
|
|
458 |
✓✓✓✓
|
44064 |
WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0. |
459 |
✓✓ |
44064 |
zfra(1:knon) = MAX(0.0,MIN(1.0,snow(1:knon)/(snow(1:knon)+10.0))) |
460 |
|
|
|
461 |
|
|
|
462 |
|
|
END IF ! landice_opt |
463 |
|
|
|
464 |
|
|
|
465 |
|
|
!**************************************************************************************** |
466 |
|
|
! Send run-off on land-ice to coupler if coupled ocean. |
467 |
|
|
! run_off_lic has been calculated in fonte_neige or surf_inlandsis |
468 |
|
|
! If landice_opt>=2, corresponding call is done from surf_land_orchidee |
469 |
|
|
!**************************************************************************************** |
470 |
✗✓✗✗
|
288 |
IF (type_ocean=='couple' .AND. landice_opt .LT. 2) THEN |
471 |
|
|
! Compress fraction where run_off_lic is active (here all pctsrf(is_lic)) |
472 |
|
|
run_off_lic_frac(:)=0.0 |
473 |
|
|
DO j = 1, knon |
474 |
|
|
i = knindex(j) |
475 |
|
|
run_off_lic_frac(j) = pctsrf(i,is_lic) |
476 |
|
|
ENDDO |
477 |
|
|
|
478 |
|
|
CALL cpl_send_landice_fields(itime, knon, knindex, run_off_lic, run_off_lic_frac) |
479 |
|
|
ENDIF |
480 |
|
|
|
481 |
|
|
! transfer runoff rate [kg/m2/s](!) to physiq for output |
482 |
✓✓ |
44064 |
runoff(1:knon)=run_off_lic(1:knon)/dtime |
483 |
|
|
|
484 |
✓✓ |
286560 |
snow_o=0. |
485 |
✓✓ |
286560 |
zfra_o = 0. |
486 |
✓✓ |
44064 |
DO j = 1, knon |
487 |
|
43776 |
i = knindex(j) |
488 |
|
43776 |
snow_o(i) = snow(j) |
489 |
|
44064 |
zfra_o(i) = zfra(j) |
490 |
|
|
ENDDO |
491 |
|
|
|
492 |
|
|
|
493 |
|
|
!albedo SB >>> |
494 |
|
|
select case(NSW) |
495 |
|
|
case(2) |
496 |
|
|
alb_dir(1:knon,1)=alb1(1:knon) |
497 |
|
|
alb_dir(1:knon,2)=alb2(1:knon) |
498 |
|
|
case(4) |
499 |
|
|
alb_dir(1:knon,1)=alb1(1:knon) |
500 |
|
|
alb_dir(1:knon,2)=alb2(1:knon) |
501 |
|
|
alb_dir(1:knon,3)=alb2(1:knon) |
502 |
|
|
alb_dir(1:knon,4)=alb2(1:knon) |
503 |
|
|
case(6) |
504 |
✓✓ |
44064 |
alb_dir(1:knon,1)=alb1(1:knon) |
505 |
✓✓ |
44064 |
alb_dir(1:knon,2)=alb1(1:knon) |
506 |
✓✓ |
44064 |
alb_dir(1:knon,3)=alb1(1:knon) |
507 |
✓✓ |
44064 |
alb_dir(1:knon,4)=alb2(1:knon) |
508 |
✓✓ |
44064 |
alb_dir(1:knon,5)=alb2(1:knon) |
509 |
✓✓ |
44064 |
alb_dir(1:knon,6)=alb2(1:knon) |
510 |
|
|
|
511 |
✗✗✗✓ ✓✗✗✗
|
576 |
IF ((landice_opt .EQ. 1) .AND. (iflag_albcalc .EQ. 2)) THEN |
512 |
|
|
alb_dir(1:knon,1)=alb6(1:knon,1) |
513 |
|
|
alb_dir(1:knon,2)=alb6(1:knon,2) |
514 |
|
|
alb_dir(1:knon,3)=alb6(1:knon,3) |
515 |
|
|
alb_dir(1:knon,4)=alb6(1:knon,4) |
516 |
|
|
alb_dir(1:knon,5)=alb6(1:knon,5) |
517 |
|
|
alb_dir(1:knon,6)=alb6(1:knon,6) |
518 |
|
|
ENDIF |
519 |
|
|
|
520 |
|
|
end select |
521 |
✓✓✓✓
|
1719648 |
alb_dif=alb_dir |
522 |
|
|
!albedo SB <<< |
523 |
|
|
|
524 |
|
|
|
525 |
|
288 |
END SUBROUTINE surf_landice |
526 |
|
|
! |
527 |
|
|
!**************************************************************************************** |
528 |
|
|
! |
529 |
|
|
END MODULE surf_landice_mod |
530 |
|
|
|
531 |
|
|
|
532 |
|
|
|