/[lmdze]/trunk/phylmd/Interface_surf/interfsurf_hq.f
ViewVC logotype

Annotation of /trunk/phylmd/Interface_surf/interfsurf_hq.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 98 - (hide annotations)
Tue May 13 17:23:16 2014 UTC (9 years, 11 months ago) by guez
File size: 24544 byte(s) 
Split inter_barxy.f : one procedure per module, one module per
file. Grouped the files into a directory.

Split orbite.f.

Value of raz_date read from the namelist is taken into account
(resetting the step counter) even if annee_ref == anneeref and day_ref
== dayref. raz_date is no longer modified by gcm main unit. (Following
LMDZ.)

Removed argument klon of interfsur_lim. Renamed arguments lmt_alb,
lmt_rug to alb_new, z0_new (same name as corresponding actual
arguments in interfsurf_hq).

Removed argument klon of interfsurf_hq.

Removed arguments qs and d_qs of diagetpq. Were always
zero. Downgraded arguments d_qw, d_ql of diagetpq to local variables,
they were not used in physiq. Removed all computations for solid water
in diagetpq, was just zero.


Downgraded arguments fs_bound, fq_bound of diagphy to local variables,
they were not used in physiq. Encapsulated in a test on iprt all
computations in diagphy.

Removed parameter nbtr of module dimphy. Replaced it everywhere in the
program by nqmx - 2.

Removed parameter rnpb of procedure physiq. Kept the true case in
physiq and phytrac. Could not work with false case anyway.

Removed arguments klon, llm, airephy of qcheck. Removed argument ftsol
of initrrnpb, was not used.
1 guez 54 module interfsurf_hq_m
2    
3     implicit none
4    
5     contains
6    
7 guez 98 SUBROUTINE interfsurf_hq(itime, dtime, jour, rmu0, iim, jjm, nisurf, knon, &
8     knindex, pctsrf, rlat, debut, ok_veget, soil_model, nsoilmx, tsoil, &
9     qsol, u1_lay, v1_lay, temp_air, spechum, tq_cdrag, petAcoef, peqAcoef, &
10     petBcoef, peqBcoef, precip_rain, precip_snow, fder, rugos, rugoro, &
11     snow, qsurf, tsurf, p1lay, ps, radsol, ocean, evap, fluxsens, fluxlat, &
12     dflux_l, dflux_s, tsurf_new, alb_new, alblw, z0_new, pctsrf_new, &
13     agesno, fqcalving, ffonte, run_off_lic_0, flux_o, flux_g, tslab, seaice)
14 guez 54
15     ! Cette routine sert d'aiguillage entre l'atmosphère et la surface
16     ! en général (sols continentaux, océans, glaces) pour les flux de
17 guez 98 ! chaleur et d'humidité.
18 guez 54
19 guez 98 ! Laurent Fairhead, 02/2000
20 guez 54
21 guez 72 USE abort_gcm_m, ONLY: abort_gcm
22     USE albsno_m, ONLY: albsno
23     USE calcul_fluxs_m, ONLY: calcul_fluxs
24 guez 98 USE dimphy, ONLY: klon
25 guez 72 USE fonte_neige_m, ONLY: fonte_neige
26     USE gath_cpl, ONLY: gath2cpl
27     USE indicesol, ONLY: epsfra, is_lic, is_oce, is_sic, is_ter, nbsrf
28     USE interface_surf, ONLY: coastalflow, riverflow, run_off, &
29     run_off_lic, conf_interface, tmp_rcoa, tmp_rlic, tmp_rriv
30     USE interfoce_lim_m, ONLY: interfoce_lim
31     USE interfoce_slab_m, ONLY: interfoce_slab
32     USE interfsur_lim_m, ONLY: interfsur_lim
33     USE suphec_m, ONLY: rcpd, rlstt, rlvtt, rtt
34 guez 54
35     ! Parametres d'entree
36     ! input:
37     ! iim, jjm nbres de pts de grille
38     ! dtime pas de temps de la physique (en s)
39     ! jour jour dans l'annee en cours,
40     ! rmu0 cosinus de l'angle solaire zenithal
41     ! nisurf index de la surface a traiter (1 = sol continental)
42     ! knon nombre de points de la surface a traiter
43     ! knindex index des points de la surface a traiter
44     ! pctsrf tableau des pourcentages de surface de chaque maille
45     ! rlat latitudes
46     ! debut logical: 1er appel a la physique
47     ! ok_veget logical: appel ou non au schema de surface continental
48     ! (si false calcul simplifie des fluxs sur les continents)
49     ! u1_lay vitesse u 1ere couche
50     ! v1_lay vitesse v 1ere couche
51     ! temp_air temperature de l'air 1ere couche
52     ! spechum humidite specifique 1ere couche
53     ! tq_cdrag cdrag
54     ! petAcoef coeff. A de la resolution de la CL pour t
55     ! peqAcoef coeff. A de la resolution de la CL pour q
56     ! petBcoef coeff. B de la resolution de la CL pour t
57     ! peqBcoef coeff. B de la resolution de la CL pour q
58     ! precip_rain precipitation liquide
59     ! precip_snow precipitation solide
60     ! tsurf temperature de surface
61     ! tslab temperature slab ocean
62     ! pctsrf_slab pourcentages (0-1) des sous-surfaces dans le slab
63     ! tmp_pctsrf_slab = pctsrf_slab
64     ! p1lay pression 1er niveau (milieu de couche)
65     ! ps pression au sol
66     ! radsol rayonnement net aus sol (LW + SW)
67     ! ocean type d'ocean utilise ("force" ou "slab" mais pas "couple")
68     ! fder derivee des flux (pour le couplage)
69     ! rugos rugosite
70     ! rugoro rugosite orographique
71     ! run_off_lic_0 runoff glacier du pas de temps precedent
72 guez 72 integer, intent(IN):: itime ! numero du pas de temps
73     integer, intent(IN):: iim, jjm
74     real, intent(IN):: dtime
75     integer, intent(IN):: jour
76     real, intent(IN):: rmu0(klon)
77     integer, intent(IN):: nisurf
78     integer, intent(IN):: knon
79     integer, dimension(klon), intent(in):: knindex
80 guez 62 real, intent(IN):: pctsrf(klon, nbsrf)
81 guez 72 logical, intent(IN):: debut, ok_veget
82     real, dimension(klon), intent(IN):: rlat
83     real, dimension(klon), intent(INOUT):: tq_cdrag
84     real, dimension(klon), intent(IN):: u1_lay, v1_lay
85     real, dimension(klon), intent(IN):: temp_air, spechum
86     real, dimension(klon), intent(IN):: petAcoef, peqAcoef
87     real, dimension(klon), intent(IN):: petBcoef, peqBcoef
88     real, dimension(klon), intent(IN):: precip_rain, precip_snow
89     real, dimension(klon), intent(IN):: ps
90     real, dimension(klon), intent(IN):: tsurf, p1lay
91 guez 54 !IM: "slab" ocean
92 guez 72 real, dimension(klon), intent(INOUT):: tslab
93     real, allocatable, dimension(:), save:: tmp_tslab
94     real, dimension(klon), intent(OUT):: flux_o, flux_g
95     real, dimension(klon), intent(INOUT):: seaice ! glace de mer (kg/m2)
96     REAL, DIMENSION(klon), INTENT(INOUT):: radsol, fder
97     real, dimension(klon), intent(IN):: rugos, rugoro
98 guez 54 character(len=*), intent(IN):: ocean
99 guez 72 real, dimension(klon), intent(INOUT):: evap, snow, qsurf
100 guez 54 !! PB ajout pour soil
101     logical, intent(in):: soil_model
102 guez 72 integer:: nsoilmx
103     REAL, DIMENSION(klon, nsoilmx):: tsoil
104     REAL, dimension(klon), intent(INOUT):: qsol
105     REAL, dimension(klon):: soilcap
106     REAL, dimension(klon):: soilflux
107 guez 54
108     ! Parametres de sortie
109     ! output:
110     ! evap evaporation totale
111     ! fluxsens flux de chaleur sensible
112     ! fluxlat flux de chaleur latente
113     ! tsurf_new temperature au sol
114     ! alb_new albedo
115     ! z0_new surface roughness
116     ! pctsrf_new nouvelle repartition des surfaces
117     real, dimension(klon), intent(OUT):: fluxsens, fluxlat
118 guez 72 real, dimension(klon), intent(OUT):: tsurf_new, alb_new
119 guez 54 real, dimension(klon), intent(OUT):: alblw
120 guez 72 real, dimension(klon), intent(OUT):: z0_new
121 guez 54 real, dimension(klon), intent(OUT):: dflux_l, dflux_s
122 guez 72 real, dimension(klon, nbsrf), intent(OUT):: pctsrf_new
123 guez 54 real, dimension(klon), intent(INOUT):: agesno
124     real, dimension(klon), intent(INOUT):: run_off_lic_0
125    
126     ! Flux thermique utiliser pour fondre la neige
127     !jld a rajouter real, dimension(klon), intent(INOUT):: ffonte
128     real, dimension(klon), intent(INOUT):: ffonte
129     ! Flux d'eau "perdue" par la surface et nécessaire pour que limiter la
130     ! hauteur de neige, en kg/m2/s
131     !jld a rajouter real, dimension(klon), intent(INOUT):: fqcalving
132     real, dimension(klon), intent(INOUT):: fqcalving
133     !IM: "slab" ocean - Local
134 guez 72 real, parameter:: t_grnd=271.35
135     real, dimension(klon):: zx_sl
136 guez 54 integer i
137 guez 72 real, allocatable, dimension(:), save:: tmp_flux_o, tmp_flux_g
138     real, allocatable, dimension(:), save:: tmp_radsol
139     real, allocatable, dimension(:, :), save:: tmp_pctsrf_slab
140     real, allocatable, dimension(:), save:: tmp_seaice
141 guez 54
142     ! Local
143 guez 72 character (len = 20), save:: modname = 'interfsurf_hq'
144     character (len = 80):: abort_message
145     logical, save:: first_call = .true.
146     integer, save:: error
147     integer:: ii
148     logical, save:: check = .false.
149 guez 54 real, dimension(klon):: cal, beta, dif_grnd, capsol
150 guez 72 real, parameter:: calice=1.0/(5.1444e+06*0.15), tau_gl=86400.*5.
151     real, parameter:: calsno=1./(2.3867e+06*.15)
152 guez 54 real, dimension(klon):: tsurf_temp
153     real, dimension(klon):: alb_neig, alb_eau
154     real, DIMENSION(klon):: zfra
155 guez 72 logical:: cumul = .false.
156     INTEGER, dimension(1):: iloc
157 guez 54 real, dimension(klon):: fder_prev
158 guez 72 REAL, dimension(klon):: bidule
159 guez 54
160     !-------------------------------------------------------------
161    
162     if (check) write(*, *) 'Entree ', modname
163    
164     ! On doit commencer par appeler les schemas de surfaces continentales
165     ! car l'ocean a besoin du ruissellement qui est y calcule
166    
167     if (first_call) then
168 guez 72 call conf_interface
169 guez 54 if (nisurf /= is_ter .and. klon > 1) then
170     write(*, *)' *** Warning ***'
171     write(*, *)' nisurf = ', nisurf, ' /= is_ter = ', is_ter
172     write(*, *)'or on doit commencer par les surfaces continentales'
173     abort_message='voir ci-dessus'
174     call abort_gcm(modname, abort_message, 1)
175     endif
176     if (ocean /= 'slab' .and. ocean /= 'force') then
177     write(*, *)' *** Warning ***'
178     write(*, *)'Option couplage pour l''ocean = ', ocean
179     abort_message='option pour l''ocean non valable'
180     call abort_gcm(modname, abort_message, 1)
181     endif
182     if ( is_oce > is_sic ) then
183     write(*, *)' *** Warning ***'
184     write(*, *)' Pour des raisons de sequencement dans le code'
185     write(*, *)' l''ocean doit etre traite avant la banquise'
186     write(*, *)' or is_oce = ', is_oce, '> is_sic = ', is_sic
187     abort_message='voir ci-dessus'
188     call abort_gcm(modname, abort_message, 1)
189     endif
190     endif
191     first_call = .false.
192    
193     ! Initialisations diverses
194    
195     ffonte(1:knon)=0.
196     fqcalving(1:knon)=0.
197    
198     cal = 999999.
199     beta = 999999.
200     dif_grnd = 999999.
201     capsol = 999999.
202     alb_new = 999999.
203     z0_new = 999999.
204     alb_neig = 999999.
205     tsurf_new = 999999.
206     alblw = 999999.
207    
208     !IM: "slab" ocean; initialisations
209     flux_o = 0.
210     flux_g = 0.
211    
212     if (.not. allocated(tmp_flux_o)) then
213     allocate(tmp_flux_o(klon), stat = error)
214     DO i=1, knon
215     tmp_flux_o(knindex(i))=flux_o(i)
216     ENDDO
217     if (error /= 0) then
218     abort_message='Pb allocation tmp_flux_o'
219     call abort_gcm(modname, abort_message, 1)
220     endif
221     endif
222     if (.not. allocated(tmp_flux_g)) then
223     allocate(tmp_flux_g(klon), stat = error)
224     DO i=1, knon
225     tmp_flux_g(knindex(i))=flux_g(i)
226     ENDDO
227     if (error /= 0) then
228     abort_message='Pb allocation tmp_flux_g'
229     call abort_gcm(modname, abort_message, 1)
230     endif
231     endif
232     if (.not. allocated(tmp_radsol)) then
233     allocate(tmp_radsol(klon), stat = error)
234     if (error /= 0) then
235     abort_message='Pb allocation tmp_radsol'
236     call abort_gcm(modname, abort_message, 1)
237     endif
238     endif
239     DO i=1, knon
240     tmp_radsol(knindex(i))=radsol(i)
241     ENDDO
242     if (.not. allocated(tmp_pctsrf_slab)) then
243     allocate(tmp_pctsrf_slab(klon, nbsrf), stat = error)
244     if (error /= 0) then
245     abort_message='Pb allocation tmp_pctsrf_slab'
246     call abort_gcm(modname, abort_message, 1)
247     endif
248     DO i=1, klon
249     tmp_pctsrf_slab(i, 1:nbsrf)=pctsrf(i, 1:nbsrf)
250     ENDDO
251     endif
252    
253     if (.not. allocated(tmp_seaice)) then
254     allocate(tmp_seaice(klon), stat = error)
255     if (error /= 0) then
256     abort_message='Pb allocation tmp_seaice'
257     call abort_gcm(modname, abort_message, 1)
258     endif
259     DO i=1, klon
260     tmp_seaice(i)=seaice(i)
261     ENDDO
262     endif
263    
264     if (.not. allocated(tmp_tslab)) then
265     allocate(tmp_tslab(klon), stat = error)
266     if (error /= 0) then
267     abort_message='Pb allocation tmp_tslab'
268     call abort_gcm(modname, abort_message, 1)
269     endif
270     endif
271     DO i=1, klon
272     tmp_tslab(i)=tslab(i)
273     ENDDO
274    
275     ! Aiguillage vers les differents schemas de surface
276    
277     if (nisurf == is_ter) then
278     ! Surface "terre" appel a l'interface avec les sols continentaux
279    
280     ! allocation du run-off
281     if (.not. allocated(coastalflow)) then
282     allocate(coastalflow(knon), stat = error)
283     if (error /= 0) then
284     abort_message='Pb allocation coastalflow'
285     call abort_gcm(modname, abort_message, 1)
286     endif
287     allocate(riverflow(knon), stat = error)
288     if (error /= 0) then
289     abort_message='Pb allocation riverflow'
290     call abort_gcm(modname, abort_message, 1)
291     endif
292     allocate(run_off(knon), stat = error)
293     if (error /= 0) then
294     abort_message='Pb allocation run_off'
295     call abort_gcm(modname, abort_message, 1)
296     endif
297     !cym
298     run_off=0.0
299     !cym
300    
301     ALLOCATE (tmp_rriv(iim, jjm+1), stat=error)
302     if (error /= 0) then
303     abort_message='Pb allocation tmp_rriv'
304     call abort_gcm(modname, abort_message, 1)
305     endif
306     ALLOCATE (tmp_rcoa(iim, jjm+1), stat=error)
307     if (error /= 0) then
308     abort_message='Pb allocation tmp_rcoa'
309     call abort_gcm(modname, abort_message, 1)
310     endif
311     ALLOCATE (tmp_rlic(iim, jjm+1), stat=error)
312     if (error /= 0) then
313     abort_message='Pb allocation tmp_rlic'
314     call abort_gcm(modname, abort_message, 1)
315     endif
316     tmp_rriv = 0.0
317     tmp_rcoa = 0.0
318     tmp_rlic = 0.0
319     else if (size(coastalflow) /= knon) then
320     write(*, *)'Bizarre, le nombre de points continentaux'
321     write(*, *)'a change entre deux appels. J''arrete ...'
322     abort_message='voir ci-dessus'
323     call abort_gcm(modname, abort_message, 1)
324     endif
325     coastalflow = 0.
326     riverflow = 0.
327    
328     ! Calcul age de la neige
329    
330     if (.not. ok_veget) then
331     ! calcul albedo: lecture albedo fichier boundary conditions
332     ! puis ajout albedo neige
333 guez 98 call interfsur_lim(itime, dtime, jour, nisurf, knon, knindex, &
334 guez 54 debut, alb_new, z0_new)
335    
336     ! calcul snow et qsurf, hydrol adapté
337     CALL calbeta(dtime, nisurf, knon, snow, qsol, beta, capsol, dif_grnd)
338    
339     IF (soil_model) THEN
340     CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, &
341     soilflux)
342     cal(1:knon) = RCPD / soilcap(1:knon)
343     radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
344     ELSE
345     cal = RCPD * capsol
346     ENDIF
347     CALL calcul_fluxs( klon, knon, nisurf, dtime, &
348     tsurf, p1lay, cal, beta, tq_cdrag, ps, &
349     precip_rain, precip_snow, snow, qsurf, &
350     radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
351     petAcoef, peqAcoef, petBcoef, peqBcoef, &
352     tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
353    
354     CALL fonte_neige( klon, knon, nisurf, dtime, &
355     tsurf, p1lay, cal, beta, tq_cdrag, ps, &
356     precip_rain, precip_snow, snow, qsol, &
357     radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
358     petAcoef, peqAcoef, petBcoef, peqBcoef, &
359     tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
360     fqcalving, ffonte, run_off_lic_0)
361    
362     call albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)
363     where (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
364     zfra(1:knon) = max(0.0, min(1.0, snow(1:knon)/(snow(1:knon)+10.0)))
365     alb_new(1 : knon) = alb_neig(1 : knon) *zfra(1:knon) + &
366     alb_new(1 : knon)*(1.0-zfra(1:knon))
367     z0_new = sqrt(z0_new**2+rugoro**2)
368     alblw(1 : knon) = alb_new(1 : knon)
369     endif
370    
371     ! Remplissage des pourcentages de surface
372     pctsrf_new(:, nisurf) = pctsrf(:, nisurf)
373     else if (nisurf == is_oce) then
374     ! Surface "ocean" appel a l'interface avec l'ocean
375     if (ocean == 'slab') then
376     tsurf_new = tsurf
377     pctsrf_new = tmp_pctsrf_slab
378     else
379     ! lecture conditions limites
380     call interfoce_lim(itime, dtime, jour, klon, nisurf, knon, knindex, &
381     debut, tsurf_new, pctsrf_new)
382     endif
383    
384     tsurf_temp = tsurf_new
385     cal = 0.
386     beta = 1.
387     dif_grnd = 0.
388     alb_neig = 0.
389     agesno = 0.
390    
391     call calcul_fluxs( klon, knon, nisurf, dtime, &
392     tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
393     precip_rain, precip_snow, snow, qsurf, &
394     radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
395     petAcoef, peqAcoef, petBcoef, peqBcoef, &
396     tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
397    
398     fder_prev = fder
399     fder = fder_prev + dflux_s + dflux_l
400    
401     iloc = maxloc(fder(1:klon))
402     if (check.and.fder(iloc(1))> 0.) then
403     WRITE(*, *)'**** Debug fder****'
404     WRITE(*, *)'max fder(', iloc(1), ') = ', fder(iloc(1))
405     WRITE(*, *)'fder_prev, dflux_s, dflux_l', fder_prev(iloc(1)), &
406     dflux_s(iloc(1)), dflux_l(iloc(1))
407     endif
408    
409     !IM: flux ocean-atmosphere utile pour le "slab" ocean
410     DO i=1, knon
411     zx_sl(i) = RLVTT
412     if (tsurf_new(i) .LT. RTT) zx_sl(i) = RLSTT
413     flux_o(i) = fluxsens(i)-evap(i)*zx_sl(i)
414     tmp_flux_o(knindex(i)) = flux_o(i)
415     tmp_radsol(knindex(i))=radsol(i)
416     ENDDO
417    
418     ! 2eme appel a interfoce pour le cumul des champs (en particulier
419     ! fluxsens et fluxlat calcules dans calcul_fluxs)
420    
421     if (ocean == 'slab ') then
422     seaice=tmp_seaice
423     cumul = .true.
424     call interfoce_slab(klon, debut, itime, dtime, jour, &
425     tmp_radsol, tmp_flux_o, tmp_flux_g, pctsrf, &
426     tslab, seaice, pctsrf_new)
427    
428     tmp_pctsrf_slab=pctsrf_new
429     DO i=1, knon
430     tsurf_new(i)=tslab(knindex(i))
431     ENDDO
432     endif
433    
434     ! calcul albedo
435     if ( minval(rmu0) == maxval(rmu0) .and. minval(rmu0) == -999.999 ) then
436     CALL alboc(FLOAT(jour), rlat, alb_eau)
437     else ! cycle diurne
438     CALL alboc_cd(rmu0, alb_eau)
439     endif
440     DO ii =1, knon
441     alb_new(ii) = alb_eau(knindex(ii))
442     enddo
443    
444     z0_new = sqrt(rugos**2 + rugoro**2)
445     alblw(1:knon) = alb_new(1:knon)
446     else if (nisurf == is_sic) then
447     if (check) write(*, *)'sea ice, nisurf = ', nisurf
448    
449     ! Surface "glace de mer" appel a l'interface avec l'ocean
450    
451    
452     if (ocean == 'slab ') then
453     pctsrf_new=tmp_pctsrf_slab
454    
455     DO ii = 1, knon
456     tsurf_new(ii) = tsurf(ii)
457     IF (pctsrf_new(knindex(ii), nisurf) < EPSFRA) then
458     snow(ii) = 0.0
459     tsurf_new(ii) = RTT - 1.8
460     IF (soil_model) tsoil(ii, :) = RTT -1.8
461     ENDIF
462     ENDDO
463    
464     CALL calbeta(dtime, nisurf, knon, snow, qsol, beta, capsol, dif_grnd)
465    
466     IF (soil_model) THEN
467 guez 72 CALL soil(dtime, nisurf, knon, snow, tsurf_new, tsoil, soilcap, &
468     soilflux)
469 guez 54 cal(1:knon) = RCPD / soilcap(1:knon)
470     radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
471     ELSE
472     dif_grnd = 1.0 / tau_gl
473     cal = RCPD * calice
474     WHERE (snow > 0.0) cal = RCPD * calsno
475     ENDIF
476     tsurf_temp = tsurf_new
477     beta = 1.0
478    
479     ELSE
480     ! ! lecture conditions limites
481     CALL interfoce_lim(itime, dtime, jour, &
482     klon, nisurf, knon, knindex, &
483     debut, &
484     tsurf_new, pctsrf_new)
485    
486     !IM cf LF
487     DO ii = 1, knon
488     tsurf_new(ii) = tsurf(ii)
489     !IMbad IF (pctsrf_new(ii, nisurf) < EPSFRA) then
490     IF (pctsrf_new(knindex(ii), nisurf) < EPSFRA) then
491     snow(ii) = 0.0
492     !IM cf LF/JLD tsurf(ii) = RTT - 1.8
493     tsurf_new(ii) = RTT - 1.8
494     IF (soil_model) tsoil(ii, :) = RTT -1.8
495     endif
496     enddo
497    
498     CALL calbeta(dtime, nisurf, knon, snow, qsol, beta, capsol, dif_grnd)
499    
500     IF (soil_model) THEN
501 guez 72 CALL soil(dtime, nisurf, knon, snow, tsurf_new, tsoil, soilcap, &
502     soilflux)
503 guez 54 cal(1:knon) = RCPD / soilcap(1:knon)
504     radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
505     dif_grnd = 0.
506     ELSE
507     dif_grnd = 1.0 / tau_gl
508     cal = RCPD * calice
509     WHERE (snow > 0.0) cal = RCPD * calsno
510     ENDIF
511     !IMbadtsurf_temp = tsurf
512     tsurf_temp = tsurf_new
513     beta = 1.0
514     ENDIF
515    
516     CALL calcul_fluxs( klon, knon, nisurf, dtime, &
517     tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
518     precip_rain, precip_snow, snow, qsurf, &
519     radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
520     petAcoef, peqAcoef, petBcoef, peqBcoef, &
521     tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
522    
523     !IM: flux entre l'ocean et la glace de mer pour le "slab" ocean
524     DO i = 1, knon
525     flux_g(i) = 0.0
526    
527     !IM: faire dependre le coefficient de conduction de la glace de mer
528     ! de l'epaisseur de la glace de mer, dans l'hypothese ou le coeff.
529     ! actuel correspond a 3m de glace de mer, cf. L.Li
530    
531     ! IF(1.EQ.0) THEN
532     ! IF(siceh(i).GT.0.) THEN
533     ! new_dif_grnd(i) = dif_grnd(i)*3./siceh(i)
534     ! ELSE
535     ! new_dif_grnd(i) = 0.
536     ! ENDIF
537     ! ENDIF !(1.EQ.0) THEN
538    
539     IF (cal(i).GT.1.0e-15) flux_g(i)=(tsurf_new(i)-t_grnd) &
540     * dif_grnd(i) *RCPD/cal(i)
541     ! & * new_dif_grnd(i) *RCPD/cal(i)
542     tmp_flux_g(knindex(i))=flux_g(i)
543     tmp_radsol(knindex(i))=radsol(i)
544     ENDDO
545    
546     CALL fonte_neige( klon, knon, nisurf, dtime, &
547     tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
548     precip_rain, precip_snow, snow, qsol, &
549     radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
550     petAcoef, peqAcoef, petBcoef, peqBcoef, &
551     tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
552     fqcalving, ffonte, run_off_lic_0)
553    
554     ! calcul albedo
555    
556     CALL albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)
557     WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
558     zfra(1:knon) = MAX(0.0, MIN(1.0, snow(1:knon)/(snow(1:knon)+10.0)))
559     alb_new(1 : knon) = alb_neig(1 : knon) *zfra(1:knon) + &
560     0.6 * (1.0-zfra(1:knon))
561    
562     fder_prev = fder
563     fder = fder_prev + dflux_s + dflux_l
564    
565     iloc = maxloc(fder(1:klon))
566     if (check.and.fder(iloc(1))> 0.) then
567     WRITE(*, *)'**** Debug fder ****'
568     WRITE(*, *)'max fder(', iloc(1), ') = ', fder(iloc(1))
569     WRITE(*, *)'fder_prev, dflux_s, dflux_l', fder_prev(iloc(1)), &
570     dflux_s(iloc(1)), dflux_l(iloc(1))
571     endif
572    
573    
574     ! 2eme appel a interfoce pour le cumul et le passage des flux a l'ocean
575    
576     z0_new = 0.002
577     z0_new = SQRT(z0_new**2+rugoro**2)
578     alblw(1:knon) = alb_new(1:knon)
579    
580     else if (nisurf == is_lic) then
581    
582     if (check) write(*, *)'glacier, nisurf = ', nisurf
583    
584     if (.not. allocated(run_off_lic)) then
585     allocate(run_off_lic(knon), stat = error)
586     if (error /= 0) then
587     abort_message='Pb allocation run_off_lic'
588     call abort_gcm(modname, abort_message, 1)
589     endif
590     run_off_lic = 0.
591     endif
592    
593     ! Surface "glacier continentaux" appel a l'interface avec le sol
594    
595     IF (soil_model) THEN
596     CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, soilflux)
597     cal(1:knon) = RCPD / soilcap(1:knon)
598     radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
599     ELSE
600     cal = RCPD * calice
601     WHERE (snow > 0.0) cal = RCPD * calsno
602     ENDIF
603     beta = 1.0
604     dif_grnd = 0.0
605    
606     call calcul_fluxs( klon, knon, nisurf, dtime, &
607     tsurf, p1lay, cal, beta, tq_cdrag, ps, &
608     precip_rain, precip_snow, snow, qsurf, &
609     radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
610     petAcoef, peqAcoef, petBcoef, peqBcoef, &
611     tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
612    
613     call fonte_neige( klon, knon, nisurf, dtime, &
614     tsurf, p1lay, cal, beta, tq_cdrag, ps, &
615     precip_rain, precip_snow, snow, qsol, &
616     radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
617     petAcoef, peqAcoef, petBcoef, peqBcoef, &
618     tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
619     fqcalving, ffonte, run_off_lic_0)
620    
621     ! passage du run-off des glaciers calcule dans fonte_neige au coupleur
622     bidule=0.
623     bidule(1:knon)= run_off_lic(1:knon)
624     call gath2cpl(bidule, tmp_rlic, klon, knon, iim, jjm, knindex)
625    
626     ! calcul albedo
627    
628     CALL albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)
629     WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
630     zfra(1:knon) = MAX(0.0, MIN(1.0, snow(1:knon)/(snow(1:knon)+10.0)))
631     alb_new(1 : knon) = alb_neig(1 : knon)*zfra(1:knon) + &
632     0.6 * (1.0-zfra(1:knon))
633    
634     !IM: plusieurs choix/tests sur l'albedo des "glaciers continentaux"
635     ! alb_new(1 : knon) = 0.6 !IM cf FH/GK
636     ! alb_new(1 : knon) = 0.82
637     ! alb_new(1 : knon) = 0.77 !211003 Ksta0.77
638     ! alb_new(1 : knon) = 0.8 !KstaTER0.8 & LMD_ARMIP5
639     !IM: KstaTER0.77 & LMD_ARMIP6
640     alb_new(1 : knon) = 0.77
641    
642    
643     ! Rugosite
644    
645     z0_new = rugoro
646    
647     ! Remplissage des pourcentages de surface
648    
649     pctsrf_new(:, nisurf) = pctsrf(:, nisurf)
650    
651     alblw(1:knon) = alb_new(1:knon)
652     else
653     write(*, *)'Index surface = ', nisurf
654     abort_message = 'Index surface non valable'
655     call abort_gcm(modname, abort_message, 1)
656     endif
657    
658     END SUBROUTINE interfsurf_hq
659    
660     end module interfsurf_hq_m
  ViewVC Help
Powered by ViewVC 1.1.21