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

Annotation of /trunk/phylmd/clqh.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 207 - (hide annotations)
Thu Sep 1 10:30:53 2016 UTC (7 years, 7 months ago) by guez
Original Path: trunk/Sources/phylmd/clqh.f
File size: 9556 byte(s)
New philosophy on compiler options.

Removed source code for thermcep = f. (Not used in LMDZ either.)

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 206 SUBROUTINE clqh(dtime, jour, debut, rlat, nisurf, knindex, tsoil, qsol, &
8     rmu0, rugos, rugoro, u1lay, v1lay, coef, t, q, ts, paprs, pplay, delp, &
9     radsol, albedo, snow, qsurf, precip_rain, precip_snow, fder, fluxlat, &
10     pctsrf_new_sic, agesno, d_t, d_q, d_ts, z0_new, flux_t, flux_q, &
11     dflux_s, dflux_l, fqcalving, ffonte, run_off_lic_0)
12 guez 3
13 guez 62 ! Author: Z. X. Li (LMD/CNRS)
14 guez 49 ! Date: 1993/08/18
15     ! Objet : diffusion verticale de "q" et de "h"
16 guez 3
17 guez 154 USE conf_phys_m, ONLY: iflag_pbl
18     USE dimphy, ONLY: klev, klon
19     USE dimsoil, ONLY: nsoilmx
20     USE interfsurf_hq_m, ONLY: interfsurf_hq
21     USE suphec_m, ONLY: rcpd, rd, rg, rkappa
22 guez 38
23 guez 155 REAL, intent(in):: dtime ! intervalle du temps (s)
24     integer, intent(in):: jour ! jour de l'annee en cours
25 guez 154 logical, intent(in):: debut
26     real, intent(in):: rlat(klon)
27 guez 175 integer, intent(in):: nisurf
28 guez 154 integer, intent(in):: knindex(:) ! (knon)
29 guez 202
30 guez 154 REAL tsoil(klon, nsoilmx)
31    
32     REAL, intent(inout):: qsol(klon)
33     ! column-density of water in soil, in kg m-2
34    
35 guez 155 real, intent(in):: rmu0(klon) ! cosinus de l'angle solaire zenithal
36     real rugos(klon) ! rugosite
37 guez 154 REAL rugoro(klon)
38 guez 155 REAL u1lay(klon) ! vitesse u de la 1ere couche (m / s)
39     REAL v1lay(klon) ! vitesse v de la 1ere couche (m / s)
40 guez 70
41     REAL, intent(in):: coef(:, :) ! (knon, klev)
42 guez 155 ! Le coefficient d'echange (m**2 / s) multiplie par le cisaillement
43     ! du vent (dV / dz). La premiere valeur indique la valeur de Cdrag
44 guez 70 ! (sans unite).
45    
46 guez 155 REAL t(klon, klev) ! temperature (K)
47     REAL q(klon, klev) ! humidite specifique (kg / kg)
48 guez 207 REAL, intent(in):: ts(:) ! (knon) temperature du sol (K)
49 guez 49 REAL paprs(klon, klev+1) ! pression a inter-couche (Pa)
50 guez 155 REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
51     REAL delp(klon, klev) ! epaisseur de couche en pression (Pa)
52     REAL radsol(klon) ! ray. net au sol (Solaire+IR) W / m2
53     REAL, intent(inout):: albedo(:) ! (knon) albedo de la surface
54 guez 191 REAL, intent(inout):: snow(klon) ! hauteur de neige
55 guez 155 REAL qsurf(klon) ! humidite de l'air au dessus de la surface
56 guez 101
57     real, intent(in):: precip_rain(klon)
58 guez 155 ! liquid water mass flux (kg / m2 / s), positive down
59 guez 101
60     real, intent(in):: precip_snow(klon)
61 guez 155 ! solid water mass flux (kg / m2 / s), positive down
62 guez 101
63 guez 154 real, intent(inout):: fder(klon)
64     real fluxlat(klon)
65 guez 202 real, intent(in):: pctsrf_new_sic(:) ! (klon)
66 guez 175 REAL, intent(inout):: agesno(:) ! (knon)
67 guez 155 REAL d_t(klon, klev) ! incrementation de "t"
68     REAL d_q(klon, klev) ! incrementation de "q"
69 guez 207 REAL, intent(out):: d_ts(:) ! (knon) incr\'ementation de "ts"
70 guez 154 real z0_new(klon)
71 guez 206
72     REAL, intent(out):: flux_t(:) ! (knon)
73     ! (diagnostic) flux de chaleur sensible (Cp T) à la surface,
74     ! positif vers le bas, W / m2
75    
76     REAL, intent(out):: flux_q(:) ! (knon)
77     ! flux de la vapeur d'eau à la surface, en kg / (m**2 s)
78    
79 guez 155 REAL dflux_s(klon) ! derivee du flux sensible dF / dTs
80     REAL dflux_l(klon) ! derivee du flux latent dF / dTs
81 guez 154
82 guez 150 ! Flux d'eau "perdue" par la surface et n\'ecessaire pour que limiter la
83 guez 155 ! hauteur de neige, en kg / m2 / s
84 guez 49 REAL fqcalving(klon)
85 guez 101
86 guez 154 ! Flux thermique utiliser pour fondre la neige
87     REAL ffonte(klon)
88    
89     REAL run_off_lic_0(klon)! runof glacier au pas de temps precedent
90    
91     ! Local:
92    
93 guez 206 INTEGER knon
94     REAL evap(size(knindex)) ! (knon) evaporation au sol
95 guez 154
96 guez 49 INTEGER i, k
97     REAL zx_cq(klon, klev)
98     REAL zx_dq(klon, klev)
99     REAL zx_ch(klon, klev)
100     REAL zx_dh(klon, klev)
101     REAL zx_buf1(klon)
102     REAL zx_buf2(klon)
103     REAL zx_coef(klon, klev)
104     REAL local_h(klon, klev) ! enthalpie potentielle
105     REAL local_q(klon, klev)
106     REAL psref(klon) ! pression de reference pour temperature potent.
107     REAL zx_pkh(klon, klev), zx_pkf(klon, klev)
108 guez 3
109 guez 155 ! contre-gradient pour la vapeur d'eau: (kg / kg) / metre
110 guez 49 REAL gamq(klon, 2:klev)
111 guez 155 ! contre-gradient pour la chaleur sensible: Kelvin / metre
112 guez 49 REAL gamt(klon, 2:klev)
113     REAL z_gamaq(klon, 2:klev), z_gamah(klon, 2:klev)
114     REAL zdelz
115 guez 3
116 guez 49 real temp_air(klon), spechum(klon)
117     real tq_cdrag(klon), petAcoef(klon), peqAcoef(klon)
118     real petBcoef(klon), peqBcoef(klon)
119     real p1lay(klon)
120 guez 3
121 guez 206 real tsurf_new(size(knindex)) ! (knon)
122 guez 49 real zzpk
123 guez 3
124 guez 49 !----------------------------------------------------------------
125 guez 3
126 guez 206 knon = size(knindex)
127    
128 guez 155 if (iflag_pbl == 1) then
129 guez 49 do k = 3, klev
130     do i = 1, knon
131     gamq(i, k)= 0.0
132 guez 155 gamt(i, k)= - 1.0e-03
133 guez 49 enddo
134     enddo
135     do i = 1, knon
136     gamq(i, 2) = 0.0
137 guez 155 gamt(i, 2) = - 2.5e-03
138 guez 49 enddo
139     else
140     do k = 2, klev
141     do i = 1, knon
142     gamq(i, k) = 0.0
143     gamt(i, k) = 0.0
144     enddo
145     enddo
146     endif
147    
148     DO i = 1, knon
149     psref(i) = paprs(i, 1) !pression de reference est celle au sol
150     ENDDO
151     DO k = 1, klev
152     DO i = 1, knon
153 guez 155 zx_pkh(i, k) = (psref(i) / paprs(i, k))**RKAPPA
154     zx_pkf(i, k) = (psref(i) / pplay(i, k))**RKAPPA
155 guez 49 local_h(i, k) = RCPD * t(i, k) * zx_pkf(i, k)
156     local_q(i, k) = q(i, k)
157     ENDDO
158     ENDDO
159    
160     ! Convertir les coefficients en variables convenables au calcul:
161    
162     DO k = 2, klev
163     DO i = 1, knon
164 guez 206 zx_coef(i, k) = coef(i, k) * RG / (pplay(i, k - 1) - pplay(i, k)) &
165     * (paprs(i, k) * 2 / (t(i, k)+t(i, k - 1)) / RD)**2
166     zx_coef(i, k) = zx_coef(i, k) * dtime * RG
167 guez 49 ENDDO
168     ENDDO
169    
170     ! Preparer les flux lies aux contre-gardients
171    
172     DO k = 2, klev
173     DO i = 1, knon
174 guez 155 zdelz = RD * (t(i, k - 1)+t(i, k)) / 2.0 / RG / paprs(i, k) &
175 guez 206 * (pplay(i, k - 1) - pplay(i, k))
176 guez 49 z_gamaq(i, k) = gamq(i, k) * zdelz
177 guez 206 z_gamah(i, k) = gamt(i, k) * zdelz * RCPD * zx_pkh(i, k)
178 guez 49 ENDDO
179     ENDDO
180     DO i = 1, knon
181     zx_buf1(i) = zx_coef(i, klev) + delp(i, klev)
182 guez 206 zx_cq(i, klev) = (local_q(i, klev) * delp(i, klev) &
183     - zx_coef(i, klev) * z_gamaq(i, klev)) / zx_buf1(i)
184 guez 49 zx_dq(i, klev) = zx_coef(i, klev) / zx_buf1(i)
185    
186 guez 155 zzpk=(pplay(i, klev) / psref(i))**RKAPPA
187 guez 206 zx_buf2(i) = zzpk * delp(i, klev) + zx_coef(i, klev)
188     zx_ch(i, klev) = (local_h(i, klev) * zzpk * delp(i, klev) &
189     - zx_coef(i, klev) * z_gamah(i, klev)) / zx_buf2(i)
190 guez 49 zx_dh(i, klev) = zx_coef(i, klev) / zx_buf2(i)
191     ENDDO
192 guez 155 DO k = klev - 1, 2, - 1
193 guez 49 DO i = 1, knon
194     zx_buf1(i) = delp(i, k)+zx_coef(i, k) &
195 guez 206 +zx_coef(i, k+1) * (1. - zx_dq(i, k+1))
196     zx_cq(i, k) = (local_q(i, k) * delp(i, k) &
197     +zx_coef(i, k+1) * zx_cq(i, k+1) &
198     +zx_coef(i, k+1) * z_gamaq(i, k+1) &
199     - zx_coef(i, k) * z_gamaq(i, k)) / zx_buf1(i)
200 guez 49 zx_dq(i, k) = zx_coef(i, k) / zx_buf1(i)
201    
202 guez 155 zzpk=(pplay(i, k) / psref(i))**RKAPPA
203 guez 206 zx_buf2(i) = zzpk * delp(i, k)+zx_coef(i, k) &
204     +zx_coef(i, k+1) * (1. - zx_dh(i, k+1))
205     zx_ch(i, k) = (local_h(i, k) * zzpk * delp(i, k) &
206     +zx_coef(i, k+1) * zx_ch(i, k+1) &
207     +zx_coef(i, k+1) * z_gamah(i, k+1) &
208     - zx_coef(i, k) * z_gamah(i, k)) / zx_buf2(i)
209 guez 49 zx_dh(i, k) = zx_coef(i, k) / zx_buf2(i)
210     ENDDO
211     ENDDO
212    
213     DO i = 1, knon
214 guez 206 zx_buf1(i) = delp(i, 1) + zx_coef(i, 2) * (1. - zx_dq(i, 2))
215     zx_cq(i, 1) = (local_q(i, 1) * delp(i, 1) &
216     +zx_coef(i, 2) * (z_gamaq(i, 2)+zx_cq(i, 2))) &
217 guez 155 / zx_buf1(i)
218     zx_dq(i, 1) = - 1. * RG / zx_buf1(i)
219 guez 49
220 guez 155 zzpk=(pplay(i, 1) / psref(i))**RKAPPA
221 guez 206 zx_buf2(i) = zzpk * delp(i, 1) + zx_coef(i, 2) * (1. - zx_dh(i, 2))
222     zx_ch(i, 1) = (local_h(i, 1) * zzpk * delp(i, 1) &
223     +zx_coef(i, 2) * (z_gamah(i, 2)+zx_ch(i, 2))) &
224 guez 155 / zx_buf2(i)
225     zx_dh(i, 1) = - 1. * RG / zx_buf2(i)
226 guez 49 ENDDO
227    
228     ! Appel a interfsurf (appel generique) routine d'interface avec la surface
229    
230     ! initialisation
231     petAcoef =0.
232     peqAcoef = 0.
233     petBcoef =0.
234     peqBcoef = 0.
235     p1lay =0.
236    
237     petAcoef(1:knon) = zx_ch(1:knon, 1)
238     peqAcoef(1:knon) = zx_cq(1:knon, 1)
239 guez 155 petBcoef(1:knon) = zx_dh(1:knon, 1)
240 guez 49 peqBcoef(1:knon) = zx_dq(1:knon, 1)
241 guez 70 tq_cdrag(1:knon) =coef(:knon, 1)
242 guez 49 temp_air(1:knon) =t(1:knon, 1)
243     spechum(1:knon)=q(1:knon, 1)
244     p1lay(1:knon) = pplay(1:knon, 1)
245    
246 guez 202 CALL interfsurf_hq(dtime, jour, rmu0, nisurf, knon, knindex, rlat, debut, &
247     nsoilmx, tsoil, qsol, u1lay, v1lay, temp_air, spechum, tq_cdrag, &
248     petAcoef, peqAcoef, petBcoef, peqBcoef, precip_rain, precip_snow, &
249 guez 207 fder, rugos, rugoro, snow, qsurf, ts, p1lay, psref, radsol, &
250 guez 206 evap, flux_t, fluxlat, dflux_l, dflux_s, tsurf_new, albedo, &
251 guez 202 z0_new, pctsrf_new_sic, agesno, fqcalving, ffonte, run_off_lic_0)
252 guez 49
253 guez 206 flux_q = - evap
254 guez 207 d_ts = tsurf_new - ts
255 guez 49
256     !==== une fois on a zx_h_ts, on peut faire l'iteration ========
257     DO i = 1, knon
258 guez 206 local_h(i, 1) = zx_ch(i, 1) + zx_dh(i, 1) * flux_t(i) * dtime
259     local_q(i, 1) = zx_cq(i, 1) + zx_dq(i, 1) * flux_q(i) * dtime
260 guez 49 ENDDO
261     DO k = 2, klev
262     DO i = 1, knon
263 guez 206 local_q(i, k) = zx_cq(i, k) + zx_dq(i, k) * local_q(i, k - 1)
264     local_h(i, k) = zx_ch(i, k) + zx_dh(i, k) * local_h(i, k - 1)
265 guez 49 ENDDO
266     ENDDO
267 guez 155
268 guez 49 ! Calcul tendances
269     DO k = 1, klev
270     DO i = 1, knon
271 guez 155 d_t(i, k) = local_h(i, k) / zx_pkf(i, k) / RCPD - t(i, k)
272 guez 49 d_q(i, k) = local_q(i, k) - q(i, k)
273     ENDDO
274     ENDDO
275    
276     END SUBROUTINE clqh
277    
278     end module clqh_m

  ViewVC Help
Powered by ViewVC 1.1.21