| GCC Code Coverage Report | |||||||||||||||||||||
|
|||||||||||||||||||||
| Line | Branch | Exec | Source |
1 |
!$Id $ |
||
2 |
|||
3 |
45612576 |
SUBROUTINE cltracrn( itr, dtime,u1lay, v1lay, & |
|
4 |
576 |
cdrag,coef,t,ftsol,pctsrf, & |
|
5 |
576 |
tr,trs,paprs,pplay,delp, & |
|
6 |
masktr,fshtr,hsoltr,tautr,vdeptr, & |
||
7 |
lat,d_tr,d_trs ) |
||
8 |
|||
9 |
USE dimphy |
||
10 |
USE traclmdz_mod, ONLY : id_rn, id_pb |
||
11 |
USE indice_sol_mod |
||
12 |
|||
13 |
IMPLICIT NONE |
||
14 |
!====================================================================== |
||
15 |
! Auteur(s): Alex/LMD) date: fev 99 |
||
16 |
! inspire de clqh + clvent |
||
17 |
! Objet: diffusion verticale de traceurs avec quantite de traceur ds |
||
18 |
! le sol ( reservoir de sol de radon ) |
||
19 |
! |
||
20 |
! note : pour l'instant le traceur dans le sol et le flux sont |
||
21 |
! calcules mais ils ne servent que de diagnostiques |
||
22 |
! seule la tendance sur le traceur est sortie (d_tr) |
||
23 |
!--------------------------------------------------------------------- |
||
24 |
! Arguments: |
||
25 |
! itr......input-R- le type de traceur : id_rn(radon), id_pb(plomb) |
||
26 |
! dtime....input-R- intervalle du temps (en secondes) ~ pdtphys |
||
27 |
! u1lay....input-R- vent u de la premiere couche (m/s) |
||
28 |
! v1lay....input-R- vent v de la premiere couche (m/s) |
||
29 |
! cdrag....input-R- cdrag |
||
30 |
! coef.....input-R- le coefficient d'echange (m**2/s) l>1, valable uniquement pour k entre 2 et klev |
||
31 |
! t........input-R- temperature (K) |
||
32 |
! paprs....input-R- pression a inter-couche (Pa) |
||
33 |
! pplay....input-R- pression au milieu de couche (Pa) |
||
34 |
! delp.....input-R- epaisseur de couche (Pa) |
||
35 |
! ftsol....input-R- temperature du sol (en Kelvin) |
||
36 |
! tr.......input-R- traceurs |
||
37 |
! trs......input-R- traceurs dans le sol |
||
38 |
! masktr...input-R- Masque reservoir de sol traceur (1 = reservoir) |
||
39 |
! fshtr....input-R- Flux surfacique de production dans le sol |
||
40 |
! tautr....input-R- Constante de decroissance du traceur |
||
41 |
! vdeptr...input-R- Vitesse de depot sec dans la couche brownienne |
||
42 |
! hsoltr...input-R- Epaisseur equivalente du reservoir de sol |
||
43 |
! lat......input-R- latitude en degree |
||
44 |
! d_tr.....output-R- le changement de "tr" |
||
45 |
! d_trs....output-R- le changement de "trs" |
||
46 |
!====================================================================== |
||
47 |
include "YOMCST.h" |
||
48 |
! |
||
49 |
!Entrees |
||
50 |
INTEGER,INTENT(IN) :: itr |
||
51 |
REAL,INTENT(IN) :: dtime |
||
52 |
REAL,DIMENSION(klon),INTENT(IN) :: u1lay, v1lay |
||
53 |
REAL,DIMENSION(klon),INTENT(IN) :: cdrag |
||
54 |
REAL,DIMENSION(klon,klev),INTENT(IN) :: coef, t |
||
55 |
REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol, pctsrf |
||
56 |
REAL,DIMENSION(klon,klev),INTENT(IN) :: tr |
||
57 |
REAL,DIMENSION(klon),INTENT(IN) :: trs |
||
58 |
REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs |
||
59 |
REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay, delp |
||
60 |
REAL,DIMENSION(klon),INTENT(IN) :: masktr |
||
61 |
REAL,DIMENSION(klon),INTENT(IN) :: fshtr |
||
62 |
REAL,INTENT(IN) :: hsoltr |
||
63 |
REAL,INTENT(IN) :: tautr |
||
64 |
REAL,INTENT(IN) :: vdeptr |
||
65 |
REAL,DIMENSION(klon),INTENT(IN) :: lat |
||
66 |
|||
67 |
!Sorties |
||
68 |
REAL,DIMENSION(klon,klev),INTENT(OUT) :: d_tr |
||
69 |
REAL,DIMENSION(klon),INTENT(OUT) :: d_trs ! (diagnostic) traceur ds le sol |
||
70 |
|||
71 |
!Locales |
||
72 |
1152 |
REAL,DIMENSION(klon,klev) :: flux_tr ! (diagnostic) flux de traceur |
|
73 |
INTEGER :: i, k, n, l |
||
74 |
1152 |
REAL,DIMENSION(klon) :: rotrhi |
|
75 |
1152 |
REAL,DIMENSION(klon,klev) :: zx_coef |
|
76 |
1152 |
REAL,DIMENSION(klon) :: zx_buf |
|
77 |
1152 |
REAL,DIMENSION(klon,klev) :: zx_ctr |
|
78 |
1152 |
REAL,DIMENSION(klon,klev) :: zx_dtr |
|
79 |
1152 |
REAL,DIMENSION(klon) :: zx_trs |
|
80 |
REAL :: zx_a, zx_b |
||
81 |
|||
82 |
1152 |
REAL,DIMENSION(klon,klev) :: local_tr |
|
83 |
1152 |
REAL,DIMENSION(klon) :: local_trs |
|
84 |
1152 |
REAL,DIMENSION(klon) :: zts ! champ de temperature du sol |
|
85 |
576 |
REAL,DIMENSION(klon) :: zx_alpha1, zx_alpha2 |
|
86 |
!====================================================================== |
||
87 |
!AA Pour l'instant les 4 types de surface ne sont pas pris en compte |
||
88 |
!AA On fabrique avec zts un champ de temperature de sol |
||
89 |
!AA que le pondere par la fraction de nature de sol. |
||
90 |
|||
91 |
✓✓ | 573120 |
DO i = 1,klon |
92 |
573120 |
zts(i) = 0. |
|
93 |
ENDDO |
||
94 |
|||
95 |
✓✓ | 2880 |
DO n=1,nbsrf |
96 |
✓✓ | 2293056 |
DO i = 1,klon |
97 |
2292480 |
zts(i) = zts(i) + ftsol(i,n)*pctsrf(i,n) |
|
98 |
ENDDO |
||
99 |
ENDDO |
||
100 |
|||
101 |
✓✓ | 573120 |
DO i = 1,klon |
102 |
573120 |
rotrhi(i) = RD * zts(i) / hsoltr |
|
103 |
ENDDO |
||
104 |
|||
105 |
✓✓ | 23040 |
DO k = 1, klev |
106 |
✓✓ | 22352256 |
DO i = 1, klon |
107 |
22351680 |
local_tr(i,k) = tr(i,k) |
|
108 |
ENDDO |
||
109 |
ENDDO |
||
110 |
|||
111 |
✓✓ | 573120 |
DO i = 1, klon |
112 |
573120 |
local_trs(i) = trs(i) |
|
113 |
ENDDO |
||
114 |
!====================================================================== |
||
115 |
!AA Attention si dans clmain zx_alf1(i) = 1.0 |
||
116 |
!AA Il doit y avoir coherence (dc la meme chose ici) |
||
117 |
|||
118 |
✓✓ | 573120 |
DO i = 1, klon |
119 |
!AA zx_alpha1(i) = (paprs(i,1)-pplay(i,2))/(pplay(i,1)-pplay(i,2)) |
||
120 |
572544 |
zx_alpha1(i) = 1.0 |
|
121 |
573120 |
zx_alpha2(i) = 1.0 - zx_alpha1(i) |
|
122 |
ENDDO |
||
123 |
!====================================================================== |
||
124 |
✓✓ | 573120 |
DO i = 1, klon |
125 |
zx_coef(i,1) = cdrag(i)*(1.0+SQRT(u1lay(i)**2+v1lay(i)**2)) & |
||
126 |
572544 |
*pplay(i,1)/(RD*t(i,1)) |
|
127 |
573120 |
zx_coef(i,1) = zx_coef(i,1) * dtime*RG |
|
128 |
ENDDO |
||
129 |
|||
130 |
✓✓ | 22464 |
DO k = 2, klev |
131 |
✓✓ | 21779136 |
DO i = 1, klon |
132 |
zx_coef(i,k) = coef(i,k)*RG/(pplay(i,k-1)-pplay(i,k)) & |
||
133 |
21756672 |
*(paprs(i,k)*2/(t(i,k)+t(i,k-1))/RD)**2 |
|
134 |
21778560 |
zx_coef(i,k) = zx_coef(i,k) * dtime*RG |
|
135 |
ENDDO |
||
136 |
ENDDO |
||
137 |
!====================================================================== |
||
138 |
✓✓ | 573120 |
DO i = 1, klon |
139 |
572544 |
zx_buf(i) = delp(i,klev) + zx_coef(i,klev) |
|
140 |
572544 |
zx_ctr(i,klev) = local_tr(i,klev)*delp(i,klev)/zx_buf(i) |
|
141 |
573120 |
zx_dtr(i,klev) = zx_coef(i,klev) / zx_buf(i) |
|
142 |
ENDDO |
||
143 |
|||
144 |
✓✓ | 21888 |
DO l = klev-1, 2 , -1 |
145 |
✓✓ | 21206016 |
DO i = 1, klon |
146 |
zx_buf(i) = delp(i,l)+zx_coef(i,l) & |
||
147 |
21184128 |
+zx_coef(i,l+1)*(1.-zx_dtr(i,l+1)) |
|
148 |
|||
149 |
zx_ctr(i,l) = ( local_tr(i,l)*delp(i,l) & |
||
150 |
21184128 |
+ zx_coef(i,l+1)*zx_ctr(i,l+1) )/zx_buf(i) |
|
151 |
21205440 |
zx_dtr(i,l) = zx_coef(i,l) / zx_buf(i) |
|
152 |
ENDDO |
||
153 |
ENDDO |
||
154 |
|||
155 |
✓✓ | 573120 |
DO i = 1, klon |
156 |
zx_buf(i) = delp(i,1) + zx_coef(i,2)*(1.-zx_dtr(i,2)) & |
||
157 |
+ masktr(i) * zx_coef(i,1) & |
||
158 |
572544 |
*( zx_alpha1(i)+zx_alpha2(i)*zx_dtr(i,2) ) |
|
159 |
|||
160 |
zx_ctr(i,1) = ( local_tr(i,1)*delp(i,1) & |
||
161 |
+ zx_ctr(i,2) & |
||
162 |
*(zx_coef(i,2) & |
||
163 |
- masktr(i) * zx_coef(i,1) & |
||
164 |
572544 |
*zx_alpha2(i) ) ) / zx_buf(i) |
|
165 |
573120 |
zx_dtr(i,1) = masktr(i) * zx_coef(i,1) / zx_buf(i) |
|
166 |
ENDDO |
||
167 |
!====================================================================== |
||
168 |
! Calculer d'abord local_trs nouvelle quantite dans le reservoir |
||
169 |
! de sol |
||
170 |
!===================================================================== |
||
171 |
|||
172 |
✓✓ | 573120 |
DO i = 1, klon |
173 |
!------------------------- |
||
174 |
! Au dessus des continents |
||
175 |
!-- |
||
176 |
! Le pb peut se deposer partout : vdeptr = 10-3 m/s |
||
177 |
! Le Rn est traiter commme une couche Brownienne puisque vdeptr = 0. |
||
178 |
!------------------------------------------------------------------- |
||
179 |
✓✓ | 572544 |
IF ( NINT(masktr(i)) .EQ. 1 ) THEN |
180 |
350208 |
zx_trs(i) = local_trs(i) |
|
181 |
zx_a = zx_trs(i) & |
||
182 |
+fshtr(i)*dtime*rotrhi(i) & |
||
183 |
+rotrhi(i)*masktr(i)*zx_coef(i,1)/RG & |
||
184 |
*(zx_ctr(i,1)*(zx_alpha1(i)+zx_alpha2(i)*zx_dtr(i,2)) & |
||
185 |
350208 |
+zx_alpha2(i)*zx_ctr(i,2)) |
|
186 |
! Pour l'instant, pour aller vite, le depot sec est traite comme une decroissance |
||
187 |
zx_b = 1. + rotrhi(i)*masktr(i)*zx_coef(i,1)/RG & |
||
188 |
* (1.-zx_dtr(i,1) & |
||
189 |
*(zx_alpha1(i)+zx_alpha2(i)*zx_dtr(i,2))) & |
||
190 |
+ dtime / tautr & |
||
191 |
350208 |
+ dtime * vdeptr / hsoltr |
|
192 |
350208 |
zx_trs(i) = zx_a / zx_b |
|
193 |
350208 |
local_trs(i) = zx_trs(i) |
|
194 |
ENDIF |
||
195 |
!-------------------------------------------------------- |
||
196 |
! Si on est entre 60N et 70N on divise par 2 l'emanation |
||
197 |
!-------------------------------------------------------- |
||
198 |
|||
199 |
✓✓✓✗ ✓✓✓✓ ✓✓ |
572544 |
IF ( (itr.eq.id_rn.AND.NINT(masktr(i)).EQ.1.AND.lat(i).GE.60..AND.lat(i).LE.70.).OR. & |
200 |
(itr.eq.id_pb.AND.NINT(masktr(i)).EQ.1.AND.lat(i).GE.60..AND.lat(i).LE.70.) ) THEN |
||
201 |
30816 |
zx_trs(i) = local_trs(i) |
|
202 |
zx_a = zx_trs(i) & |
||
203 |
+(fshtr(i)/2.)*dtime*rotrhi(i) & |
||
204 |
+rotrhi(i)*masktr(i)*zx_coef(i,1)/RG & |
||
205 |
*(zx_ctr(i,1)*(zx_alpha1(i)+zx_alpha2(i)*zx_dtr(i,2)) & |
||
206 |
30816 |
+zx_alpha2(i)*zx_ctr(i,2)) |
|
207 |
! |
||
208 |
zx_b = 1. + rotrhi(i)*masktr(i)*zx_coef(i,1)/RG & |
||
209 |
* (1.-zx_dtr(i,1) & |
||
210 |
*(zx_alpha1(i)+zx_alpha2(i)*zx_dtr(i,2))) & |
||
211 |
+ dtime / tautr & |
||
212 |
30816 |
+ dtime * vdeptr / hsoltr |
|
213 |
! |
||
214 |
30816 |
zx_trs(i) = zx_a / zx_b |
|
215 |
30816 |
local_trs(i) = zx_trs(i) |
|
216 |
ENDIF |
||
217 |
|||
218 |
!---------------------------------------------- |
||
219 |
! Au dessus des oceans et aux hautes latitudes |
||
220 |
!-- |
||
221 |
! au dessous de -60S pas d'emission de radon au dessus |
||
222 |
! des oceans et des continents |
||
223 |
!--------------------------------------------------------------- |
||
224 |
|||
225 |
✓✓✓✓ ✓✓✓✗ ✗✓ |
572544 |
IF ( (itr.EQ.id_rn.AND.NINT(masktr(i)).EQ.0).OR. & |
226 |
(itr.EQ.id_rn.AND.NINT(masktr(i)).EQ.1.AND.lat(i).LT.-60.)) THEN |
||
227 |
222336 |
zx_trs(i) = 0. |
|
228 |
222336 |
local_trs(i) = 0. |
|
229 |
END IF |
||
230 |
!-- |
||
231 |
! au dessus de 70 N pas d'emission de radon au dessus |
||
232 |
! des oceans et des continents |
||
233 |
!-------------------------------------------------------------- |
||
234 |
✓✓✓✓ ✓✓✓✗ ✓✓ |
572544 |
IF ( (itr.EQ.id_rn.AND.NINT(masktr(i)).EQ.0).OR. & |
235 |
(itr.EQ.id_rn.AND.NINT(masktr(i)).EQ.1.AND.lat(i).GT.70.)) THEN |
||
236 |
224352 |
zx_trs(i) = 0. |
|
237 |
224352 |
local_trs(i) = 0. |
|
238 |
END IF |
||
239 |
!--------------------------------------------- |
||
240 |
! Au dessus des oceans la source est nulle |
||
241 |
!-------------------------------------------- |
||
242 |
|||
243 |
✓✓✓✓ |
573120 |
IF (itr.eq.id_rn.AND.NINT(masktr(i)).EQ.0) THEN |
244 |
222336 |
zx_trs(i) = 0. |
|
245 |
222336 |
local_trs(i) = 0. |
|
246 |
END IF |
||
247 |
|||
248 |
ENDDO ! sur le i=1,klon |
||
249 |
! |
||
250 |
!====================================================================== |
||
251 |
! Une fois on a zx_trs, on peut faire l'iteration |
||
252 |
!====================================================================== |
||
253 |
|||
254 |
✓✓ | 573120 |
DO i = 1, klon |
255 |
573120 |
local_tr(i,1) = zx_ctr(i,1)+zx_dtr(i,1)*zx_trs(i) |
|
256 |
ENDDO |
||
257 |
✓✓ | 22464 |
DO l = 2, klev |
258 |
✓✓ | 21779136 |
DO i = 1, klon |
259 |
21778560 |
local_tr(i,l) = zx_ctr(i,l) + zx_dtr(i,l)*local_tr(i,l-1) |
|
260 |
ENDDO |
||
261 |
ENDDO |
||
262 |
!====================================================================== |
||
263 |
! Calcul du flux de traceur (flux_tr): UA/(m**2 s) |
||
264 |
!====================================================================== |
||
265 |
✓✓ | 573120 |
DO i = 1, klon |
266 |
flux_tr(i,1) = masktr(i)*zx_coef(i,1)/RG & |
||
267 |
* (zx_alpha1(i)*local_tr(i,1)+zx_alpha2(i)*local_tr(i,2) & |
||
268 |
573120 |
-zx_trs(i)) / dtime |
|
269 |
ENDDO |
||
270 |
✓✓ | 22464 |
DO l = 2, klev |
271 |
✓✓ | 21779136 |
DO i = 1, klon |
272 |
flux_tr(i,l) = zx_coef(i,l)/RG & |
||
273 |
21778560 |
* (local_tr(i,l)-local_tr(i,l-1)) / dtime |
|
274 |
ENDDO |
||
275 |
ENDDO |
||
276 |
!====================================================================== |
||
277 |
! Calcul des tendances du traceur ds le sol et dans l'atmosphere |
||
278 |
!====================================================================== |
||
279 |
✓✓ | 23040 |
DO l = 1, klev |
280 |
✓✓ | 22352256 |
DO i = 1, klon |
281 |
22351680 |
d_tr(i,l) = local_tr(i,l) - tr(i,l) |
|
282 |
ENDDO |
||
283 |
ENDDO |
||
284 |
✓✓ | 573120 |
DO i = 1, klon |
285 |
573120 |
d_trs(i) = local_trs(i) - trs(i) |
|
286 |
ENDDO |
||
287 |
|||
288 |
576 |
END SUBROUTINE cltracrn |
| Generated by: GCOVR (Version 4.2) |