/[lmdze]/trunk/Sources/phylmd/clqh.f
ViewVC logotype

Annotation of /trunk/Sources/phylmd/clqh.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 99 - (hide annotations)
Wed Jul 2 18:39:15 2014 UTC (9 years, 9 months ago) by guez
Original Path: trunk/phylmd/clqh.f
File size: 11501 byte(s)
Created procedure test_disvert (following LMDZ). Added procedures
hybrid and funcd in module disvert_m. Upgraded compute_ab from
internal procedure of disvert to module procedure. Added variables y,
ya in module disvert_m. Upgraded s from local variable of procedure
disvert to module variable.

Renamed allowed value of variable vert_sampling in procedure disvert
from "read" to "read_hybrid". Added possibility to read pressure
values, value "read_pressure". Replaced vertical distribution for
value "param" by the distribution "strato_correct" from LMDZ (but kept
the value "param"). In case "tropo", replaced 1 by dsigmin (following
LMDZ). In case "strato", replaced 0.3 by dsigmin (following LMDZ).

Changed computation of bp in procedure compute_ab.

Removed debugindex case in clmain. Removed useless argument rlon of
procedure clmain. Removed useless variables ytaux, ytauy of procedure
clmain.

Removed intermediary variables tsol, qsol, tsolsrf, tslab in procedure
etat0.

Removed variable ok_veget:. coupling with the model Orchid is not
possible. Removed variable ocean: modeling an ocean slab is not
possible.

Removed useless variables tmp_rriv and tmp_rcoa from module
interface_surf.

Moved initialization of variables da, mp, phi in procedure physiq to
to inside the test iflag_con >= 3.

1 guez 49 module clqh_m
2 guez 3
3 guez 38 IMPLICIT none
4 guez 3
5 guez 49 contains
6 guez 3
7 guez 72 SUBROUTINE clqh(dtime, itime, jour, debut, rlat, knon, nisurf, knindex, &
8 guez 99 pctsrf, soil_model, tsoil, qsol, rmu0, co2_ppm, rugos, rugoro, u1lay, &
9     v1lay, coef, t, q, ts, paprs, pplay, delp, radsol, albedo, alblw, &
10     snow, qsurf, precip_rain, precip_snow, fder, swnet, fluxlat, &
11     pctsrf_new, agesno, d_t, d_q, d_ts, z0_new, flux_t, flux_q, dflux_s, &
12     dflux_l, fqcalving, ffonte, run_off_lic_0, flux_o, flux_g, tslab, seaice)
13 guez 3
14 guez 62 ! Author: Z. X. Li (LMD/CNRS)
15 guez 49 ! Date: 1993/08/18
16     ! Objet : diffusion verticale de "q" et de "h"
17 guez 3
18 guez 49 USE conf_phys_m, ONLY : iflag_pbl
19     USE dimens_m, ONLY : iim, jjm
20     USE dimphy, ONLY : klev, klon, zmasq
21     USE dimsoil, ONLY : nsoilmx
22     USE indicesol, ONLY : is_ter, nbsrf
23 guez 54 USE interfsurf_hq_m, ONLY : interfsurf_hq
24 guez 49 USE suphec_m, ONLY : rcpd, rd, rg, rkappa
25 guez 38
26 guez 49 ! Arguments:
27 guez 70 INTEGER, intent(in):: knon
28 guez 49 REAL, intent(in):: dtime ! intervalle du temps (s)
29     REAL u1lay(klon) ! vitesse u de la 1ere couche (m/s)
30     REAL v1lay(klon) ! vitesse v de la 1ere couche (m/s)
31 guez 70
32     REAL, intent(in):: coef(:, :) ! (knon, klev)
33     ! Le coefficient d'echange (m**2/s) multiplie par le cisaillement
34     ! du vent (dV/dz). La premiere valeur indique la valeur de Cdrag
35     ! (sans unite).
36    
37 guez 49 REAL t(klon, klev) ! temperature (K)
38     REAL q(klon, klev) ! humidite specifique (kg/kg)
39     REAL ts(klon) ! temperature du sol (K)
40     REAL evap(klon) ! evaporation au sol
41     REAL paprs(klon, klev+1) ! pression a inter-couche (Pa)
42     REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
43     REAL delp(klon, klev) ! epaisseur de couche en pression (Pa)
44     REAL radsol(klon) ! ray. net au sol (Solaire+IR) W/m2
45     REAL albedo(klon) ! albedo de la surface
46     REAL alblw(klon)
47     REAL snow(klon) ! hauteur de neige
48     REAL qsurf(klon) ! humidite de l'air au dessus de la surface
49     real precip_rain(klon), precip_snow(klon)
50     REAL agesno(klon)
51     REAL rugoro(klon)
52     REAL run_off_lic_0(klon)! runof glacier au pas de temps precedent
53 guez 62 integer, intent(in):: jour ! jour de l'annee en cours
54     real, intent(in):: rmu0(klon) ! cosinus de l'angle solaire zenithal
55 guez 49 real rugos(klon) ! rugosite
56 guez 99 integer, intent(in):: knindex(klon)
57 guez 62 real, intent(in):: pctsrf(klon, nbsrf)
58 guez 72 real, intent(in):: rlat(klon)
59 guez 99 REAL, intent(in):: co2_ppm ! taux CO2 atmosphere
60 guez 3
61 guez 49 REAL d_t(klon, klev) ! incrementation de "t"
62     REAL d_q(klon, klev) ! incrementation de "q"
63     REAL d_ts(klon) ! incrementation de "ts"
64     REAL flux_t(klon, klev) ! (diagnostic) flux de la chaleur
65     ! sensible, flux de Cp*T, positif vers
66     ! le bas: j/(m**2 s) c.a.d.: W/m2
67     REAL flux_q(klon, klev) ! flux de la vapeur d'eau:kg/(m**2 s)
68     REAL dflux_s(klon) ! derivee du flux sensible dF/dTs
69     REAL dflux_l(klon) ! derivee du flux latent dF/dTs
70     !IM cf JLD
71     ! Flux thermique utiliser pour fondre la neige
72     REAL ffonte(klon)
73     ! Flux d'eau "perdue" par la surface et nĂ©cessaire pour que limiter la
74     ! hauteur de neige, en kg/m2/s
75     REAL fqcalving(klon)
76     !IM "slab" ocean
77     REAL tslab(klon) !temperature du slab ocean (K) (OCEAN='slab ')
78     REAL seaice(klon) ! glace de mer en kg/m2
79     REAL flux_o(klon) ! flux entre l'ocean et l'atmosphere W/m2
80     REAL flux_g(klon) ! flux entre l'ocean et la glace de mer W/m2
81 guez 3
82 guez 49 INTEGER i, k
83     REAL zx_cq(klon, klev)
84     REAL zx_dq(klon, klev)
85     REAL zx_ch(klon, klev)
86     REAL zx_dh(klon, klev)
87     REAL zx_buf1(klon)
88     REAL zx_buf2(klon)
89     REAL zx_coef(klon, klev)
90     REAL local_h(klon, klev) ! enthalpie potentielle
91     REAL local_q(klon, klev)
92     REAL local_ts(klon)
93     REAL psref(klon) ! pression de reference pour temperature potent.
94     REAL zx_pkh(klon, klev), zx_pkf(klon, klev)
95 guez 3
96 guez 49 ! contre-gradient pour la vapeur d'eau: (kg/kg)/metre
97     REAL gamq(klon, 2:klev)
98     ! contre-gradient pour la chaleur sensible: Kelvin/metre
99     REAL gamt(klon, 2:klev)
100     REAL z_gamaq(klon, 2:klev), z_gamah(klon, 2:klev)
101     REAL zdelz
102 guez 3
103 guez 49 ! Rajout pour l'interface
104     integer, intent(in):: itime
105     integer nisurf
106     logical, intent(in):: debut
107     real zlev1(klon)
108 guez 72 real fder(klon)
109 guez 49 real temp_air(klon), spechum(klon)
110     real epot_air(klon), ccanopy(klon)
111     real tq_cdrag(klon), petAcoef(klon), peqAcoef(klon)
112     real petBcoef(klon), peqBcoef(klon)
113 guez 72 real swnet(klon), swdown(klon)
114 guez 49 real p1lay(klon)
115     !$$$C PB ajout pour soil
116     LOGICAL, intent(in):: soil_model
117     REAL tsoil(klon, nsoilmx)
118     REAL qsol(klon)
119 guez 3
120 guez 49 ! Parametres de sortie
121     real fluxsens(klon), fluxlat(klon)
122 guez 72 real tsurf_new(klon), alb_new(klon)
123     real z0_new(klon)
124 guez 49 real pctsrf_new(klon, nbsrf)
125     ! JLD
126     real zzpk
127 guez 3
128 guez 49 character (len = 20) :: modname = 'Debut clqh'
129     LOGICAL check
130     PARAMETER (check=.false.)
131 guez 38
132 guez 49 !----------------------------------------------------------------
133 guez 3
134 guez 49 if (check) THEN
135     write(*, *) modname, ' nisurf=', nisurf
136     !C call flush(6)
137     endif
138    
139     if (check) THEN
140     WRITE(*, *)' qsurf (min, max)' &
141     , minval(qsurf(1:knon)), maxval(qsurf(1:knon))
142     !C call flush(6)
143     ENDIF
144    
145     if (iflag_pbl.eq.1) then
146     do k = 3, klev
147     do i = 1, knon
148     gamq(i, k)= 0.0
149     gamt(i, k)= -1.0e-03
150     enddo
151     enddo
152     do i = 1, knon
153     gamq(i, 2) = 0.0
154     gamt(i, 2) = -2.5e-03
155     enddo
156     else
157     do k = 2, klev
158     do i = 1, knon
159     gamq(i, k) = 0.0
160     gamt(i, k) = 0.0
161     enddo
162     enddo
163     endif
164    
165     DO i = 1, knon
166     psref(i) = paprs(i, 1) !pression de reference est celle au sol
167     local_ts(i) = ts(i)
168     ENDDO
169     DO k = 1, klev
170     DO i = 1, knon
171     zx_pkh(i, k) = (psref(i)/paprs(i, k))**RKAPPA
172     zx_pkf(i, k) = (psref(i)/pplay(i, k))**RKAPPA
173     local_h(i, k) = RCPD * t(i, k) * zx_pkf(i, k)
174     local_q(i, k) = q(i, k)
175     ENDDO
176     ENDDO
177    
178     ! Convertir les coefficients en variables convenables au calcul:
179    
180     DO k = 2, klev
181     DO i = 1, knon
182     zx_coef(i, k) = coef(i, k)*RG/(pplay(i, k-1)-pplay(i, k)) &
183     *(paprs(i, k)*2/(t(i, k)+t(i, k-1))/RD)**2
184     zx_coef(i, k) = zx_coef(i, k) * dtime*RG
185     ENDDO
186     ENDDO
187    
188     ! Preparer les flux lies aux contre-gardients
189    
190     DO k = 2, klev
191     DO i = 1, knon
192     zdelz = RD * (t(i, k-1)+t(i, k))/2.0 / RG /paprs(i, k) &
193     *(pplay(i, k-1)-pplay(i, k))
194     z_gamaq(i, k) = gamq(i, k) * zdelz
195     z_gamah(i, k) = gamt(i, k) * zdelz *RCPD * zx_pkh(i, k)
196     ENDDO
197     ENDDO
198     DO i = 1, knon
199     zx_buf1(i) = zx_coef(i, klev) + delp(i, klev)
200     zx_cq(i, klev) = (local_q(i, klev)*delp(i, klev) &
201     -zx_coef(i, klev)*z_gamaq(i, klev))/zx_buf1(i)
202     zx_dq(i, klev) = zx_coef(i, klev) / zx_buf1(i)
203    
204     zzpk=(pplay(i, klev)/psref(i))**RKAPPA
205     zx_buf2(i) = zzpk*delp(i, klev) + zx_coef(i, klev)
206     zx_ch(i, klev) = (local_h(i, klev)*zzpk*delp(i, klev) &
207     -zx_coef(i, klev)*z_gamah(i, klev))/zx_buf2(i)
208     zx_dh(i, klev) = zx_coef(i, klev) / zx_buf2(i)
209     ENDDO
210     DO k = klev-1, 2 , -1
211     DO i = 1, knon
212     zx_buf1(i) = delp(i, k)+zx_coef(i, k) &
213     +zx_coef(i, k+1)*(1.-zx_dq(i, k+1))
214     zx_cq(i, k) = (local_q(i, k)*delp(i, k) &
215     +zx_coef(i, k+1)*zx_cq(i, k+1) &
216     +zx_coef(i, k+1)*z_gamaq(i, k+1) &
217     -zx_coef(i, k)*z_gamaq(i, k))/zx_buf1(i)
218     zx_dq(i, k) = zx_coef(i, k) / zx_buf1(i)
219    
220     zzpk=(pplay(i, k)/psref(i))**RKAPPA
221     zx_buf2(i) = zzpk*delp(i, k)+zx_coef(i, k) &
222     +zx_coef(i, k+1)*(1.-zx_dh(i, k+1))
223     zx_ch(i, k) = (local_h(i, k)*zzpk*delp(i, k) &
224     +zx_coef(i, k+1)*zx_ch(i, k+1) &
225     +zx_coef(i, k+1)*z_gamah(i, k+1) &
226     -zx_coef(i, k)*z_gamah(i, k))/zx_buf2(i)
227     zx_dh(i, k) = zx_coef(i, k) / zx_buf2(i)
228     ENDDO
229     ENDDO
230    
231     DO i = 1, knon
232     zx_buf1(i) = delp(i, 1) + zx_coef(i, 2)*(1.-zx_dq(i, 2))
233     zx_cq(i, 1) = (local_q(i, 1)*delp(i, 1) &
234     +zx_coef(i, 2)*(z_gamaq(i, 2)+zx_cq(i, 2))) &
235     /zx_buf1(i)
236     zx_dq(i, 1) = -1. * RG / zx_buf1(i)
237    
238     zzpk=(pplay(i, 1)/psref(i))**RKAPPA
239     zx_buf2(i) = zzpk*delp(i, 1) + zx_coef(i, 2)*(1.-zx_dh(i, 2))
240     zx_ch(i, 1) = (local_h(i, 1)*zzpk*delp(i, 1) &
241     +zx_coef(i, 2)*(z_gamah(i, 2)+zx_ch(i, 2))) &
242     /zx_buf2(i)
243     zx_dh(i, 1) = -1. * RG / zx_buf2(i)
244     ENDDO
245    
246     ! Appel a interfsurf (appel generique) routine d'interface avec la surface
247    
248     ! initialisation
249     petAcoef =0.
250     peqAcoef = 0.
251     petBcoef =0.
252     peqBcoef = 0.
253     p1lay =0.
254    
255     petAcoef(1:knon) = zx_ch(1:knon, 1)
256     peqAcoef(1:knon) = zx_cq(1:knon, 1)
257     petBcoef(1:knon) = zx_dh(1:knon, 1)
258     peqBcoef(1:knon) = zx_dq(1:knon, 1)
259 guez 70 tq_cdrag(1:knon) =coef(:knon, 1)
260 guez 49 temp_air(1:knon) =t(1:knon, 1)
261     epot_air(1:knon) =local_h(1:knon, 1)
262     spechum(1:knon)=q(1:knon, 1)
263     p1lay(1:knon) = pplay(1:knon, 1)
264     zlev1(1:knon) = delp(1:knon, 1)
265    
266     if(nisurf.eq.is_ter) THEN
267     swdown(1:knon) = swnet(1:knon)/(1-albedo(1:knon))
268     else
269     swdown(1:knon) = swnet(1:knon)
270     endif
271     ccanopy = co2_ppm
272    
273 guez 99 CALL interfsurf_hq(itime, dtime, jour, rmu0, iim, jjm, nisurf, knon, &
274     knindex, pctsrf, rlat, debut, soil_model, nsoilmx, tsoil, qsol, &
275     u1lay, v1lay, temp_air, spechum, tq_cdrag, petAcoef, peqAcoef, &
276     petBcoef, peqBcoef, precip_rain, precip_snow, fder, rugos, rugoro, &
277     snow, qsurf, ts, p1lay, psref, radsol, evap, fluxsens, &
278     fluxlat, dflux_l, dflux_s, tsurf_new, alb_new, alblw, z0_new, &
279     pctsrf_new, agesno, fqcalving, ffonte, run_off_lic_0, flux_o, &
280     flux_g, tslab, seaice)
281 guez 49
282     do i = 1, knon
283     flux_t(i, 1) = fluxsens(i)
284     flux_q(i, 1) = - evap(i)
285     d_ts(i) = tsurf_new(i) - ts(i)
286     albedo(i) = alb_new(i)
287     enddo
288    
289     !==== une fois on a zx_h_ts, on peut faire l'iteration ========
290     DO i = 1, knon
291     local_h(i, 1) = zx_ch(i, 1) + zx_dh(i, 1)*flux_t(i, 1)*dtime
292     local_q(i, 1) = zx_cq(i, 1) + zx_dq(i, 1)*flux_q(i, 1)*dtime
293     ENDDO
294     DO k = 2, klev
295     DO i = 1, knon
296     local_q(i, k) = zx_cq(i, k) + zx_dq(i, k)*local_q(i, k-1)
297     local_h(i, k) = zx_ch(i, k) + zx_dh(i, k)*local_h(i, k-1)
298     ENDDO
299     ENDDO
300     !======================================================================
301     !== flux_q est le flux de vapeur d'eau: kg/(m**2 s) positive vers bas
302     !== flux_t est le flux de cpt (energie sensible): j/(m**2 s)
303     DO k = 2, klev
304     DO i = 1, knon
305     flux_q(i, k) = (zx_coef(i, k)/RG/dtime) &
306     * (local_q(i, k)-local_q(i, k-1)+z_gamaq(i, k))
307     flux_t(i, k) = (zx_coef(i, k)/RG/dtime) &
308     * (local_h(i, k)-local_h(i, k-1)+z_gamah(i, k)) &
309     / zx_pkh(i, k)
310     ENDDO
311     ENDDO
312     !======================================================================
313     ! Calcul tendances
314     DO k = 1, klev
315     DO i = 1, knon
316     d_t(i, k) = local_h(i, k)/zx_pkf(i, k)/RCPD - t(i, k)
317     d_q(i, k) = local_q(i, k) - q(i, k)
318     ENDDO
319     ENDDO
320    
321     END SUBROUTINE clqh
322    
323     end module clqh_m

  ViewVC Help
Powered by ViewVC 1.1.21