1 | SUBROUTINE co2cloud(ngrid,nlay,ptimestep, |
---|
2 | & pplev,pplay,pdpsrf,pzlay,pt,pdt, |
---|
3 | & pq,pdq,pdqcloudco2,pdtcloudco2, |
---|
4 | & nq,tau,tauscaling,rdust,rice,riceco2,nuice, |
---|
5 | & rsedcloudco2,rhocloudco2,zlev,pdqs_sedco2) |
---|
6 | ! to use 'getin' |
---|
7 | use dimradmars_mod, only: naerkind |
---|
8 | USE comcstfi_h |
---|
9 | USE ioipsl_getincom |
---|
10 | USE updaterad |
---|
11 | use conc_mod, only: mmean |
---|
12 | use tracer_mod, only: nqmx, igcm_co2, igcm_co2_ice, |
---|
13 | & igcm_dust_mass, igcm_dust_number, |
---|
14 | & igcm_dust_mass, igcm_dust_number, |
---|
15 | & igcm_ccnco2_mass, igcm_ccnco2_number, |
---|
16 | & rho_dust, nuiceco2_sed, nuiceco2_ref, |
---|
17 | & rho_ice_co2,r3n_q |
---|
18 | IMPLICIT NONE |
---|
19 | |
---|
20 | |
---|
21 | c======================================================================= |
---|
22 | c CO2 clouds formation |
---|
23 | c |
---|
24 | c There is a time loop specific to cloud formation |
---|
25 | c due to timescales smaller than the GCM integration timestep. |
---|
26 | c microphysics subroutine is improvedCO2clouds.F |
---|
27 | c |
---|
28 | c The co2 clouds tracers (co2_ice, ccn mass and concentration) are |
---|
29 | c sedimented at each microtimestep. pdqs_sedco2 keeps track of the |
---|
30 | c CO2 flux at the surface |
---|
31 | c |
---|
32 | c Authors: 09/2016 Joachim Audouard & Constantino Listowski |
---|
33 | c Adaptation of the water ice clouds scheme (with specific microphysics) |
---|
34 | c of Montmessin, Navarro & al. |
---|
35 | c |
---|
36 | c |
---|
37 | c======================================================================= |
---|
38 | |
---|
39 | c----------------------------------------------------------------------- |
---|
40 | c declarations: |
---|
41 | c ------------- |
---|
42 | |
---|
43 | !#include "dimensions.h" |
---|
44 | !#include "dimphys.h" |
---|
45 | #include "callkeys.h" |
---|
46 | !#include "tracer.h" |
---|
47 | !#include "comgeomfi.h" |
---|
48 | !#include "dimradmars.h" |
---|
49 | ! naerkind is set in scatterers.h (built when compiling with makegcm -s #) |
---|
50 | !#include"scatterers.h" |
---|
51 | #include "microphys.h" |
---|
52 | |
---|
53 | |
---|
54 | c Inputs: |
---|
55 | c ------ |
---|
56 | |
---|
57 | INTEGER ngrid,nlay |
---|
58 | INTEGER nq ! nombre de traceurs |
---|
59 | REAL ptimestep ! pas de temps physique (s) |
---|
60 | REAL pplev(ngrid,nlay+1) ! pression aux inter-couches (Pa) |
---|
61 | REAL pplay(ngrid,nlay) ! pression au milieu des couches (Pa) |
---|
62 | REAL pdpsrf(ngrid) ! tendence surf pressure |
---|
63 | REAL pzlay(ngrid,nlay) ! altitude at the middle of the layers |
---|
64 | REAL pt(ngrid,nlay) ! temperature at the middle of the layers (K) |
---|
65 | REAL pdt(ngrid,nlay) ! tendence temperature des autres param. |
---|
66 | real,intent(in) :: zlev(ngrid,nlay+1) ! altitude at the boundaries of the layers |
---|
67 | |
---|
68 | real pq(ngrid,nlay,nq) ! traceur (kg/kg) |
---|
69 | real pdq(ngrid,nlay,nq) ! tendance avant condensation (kg/kg.s-1) |
---|
70 | |
---|
71 | real rice(ngrid,nlay) ! Water Ice mass mean radius (m) |
---|
72 | ! used for nucleation of CO2 on ice-coated ccns |
---|
73 | |
---|
74 | REAL tau(ngrid,naerkind) ! Column dust optical depth at each point |
---|
75 | REAL tauscaling(ngrid) ! Convertion factor for dust amount |
---|
76 | real rdust(ngrid,nlay) ! Dust geometric mean radius (m) |
---|
77 | |
---|
78 | c Outputs: |
---|
79 | c ------- |
---|
80 | |
---|
81 | real pdqcloudco2(ngrid,nlay,nq) ! tendence de la condensation H2O(kg/kg.s-1) |
---|
82 | REAL pdtcloudco2(ngrid,nlay) ! tendence temperature due |
---|
83 | ! a la chaleur latente |
---|
84 | |
---|
85 | REAL riceco2(ngrid,nlay) ! Ice mass mean radius (m) |
---|
86 | ! (r_c in montmessin_2004) |
---|
87 | REAL nuice(ngrid,nlay) ! Estimated effective variance |
---|
88 | ! of the size distribution |
---|
89 | real rsedcloudco2(ngrid,nlay) ! Cloud sedimentation radius |
---|
90 | real rhocloudco2(ngrid,nlay) ! Cloud density (kg.m-3) |
---|
91 | real rhocloudco2t(ngrid,nlay) ! Cloud density (kg.m-3) |
---|
92 | real pdqs_sedco2(ngrid) ! CO2 flux at the surface |
---|
93 | c local: |
---|
94 | c ------ |
---|
95 | |
---|
96 | ! for ice radius computation |
---|
97 | REAL Mo,No |
---|
98 | REAl ccntyp |
---|
99 | |
---|
100 | ! for time loop |
---|
101 | INTEGER microstep ! time subsampling step variable |
---|
102 | INTEGER imicro ! time subsampling for coupled water microphysics & sedimentation |
---|
103 | SAVE imicro |
---|
104 | REAL microtimestep ! integration timestep for coupled water microphysics & sedimentation |
---|
105 | SAVE microtimestep |
---|
106 | |
---|
107 | ! tendency given by clouds (inside the micro loop) |
---|
108 | REAL subpdqcloudco2(ngrid,nlay,nq) ! cf. pdqcloud |
---|
109 | REAL subpdtcloudco2(ngrid,nlay) ! cf. pdtcloud |
---|
110 | |
---|
111 | ! global tendency (clouds+physics) |
---|
112 | REAL subpdq(ngrid,nlay,nq) ! cf. pdqcloud |
---|
113 | REAL subpdt(ngrid,nlay) ! cf. pdtcloud |
---|
114 | real wq(ngrid,nlay+1) ! ! displaced tracer mass (kg.m-2) during microtimestep because sedim (?/m-2) |
---|
115 | |
---|
116 | REAL satuco2(ngrid,nlay) ! co2 satu ratio for output |
---|
117 | REAL zqsatco2(ngrid,nlay) ! saturation co2 |
---|
118 | |
---|
119 | INTEGER iq,ig,l |
---|
120 | LOGICAL,SAVE :: firstcall=.true. |
---|
121 | DOUBLE PRECISION Nccnco2, Niceco2,mdustJA,ndustJA |
---|
122 | DOUBLE PRECISION Qccnco2 |
---|
123 | real :: beta |
---|
124 | |
---|
125 | real epaisseur (ngrid,nlay) ! Layer thickness (m) |
---|
126 | real masse (ngrid,nlay) ! Layer mass (kg.m-2) |
---|
127 | |
---|
128 | |
---|
129 | real tempo_traceur_t(ngrid,nlay) |
---|
130 | real tempo_traceurs(ngrid,nlay,nq) |
---|
131 | real sav_trac(ngrid,nlay,nq) |
---|
132 | real pdqsed(ngrid,nlay,nq) |
---|
133 | c ** un petit test de coherence |
---|
134 | c -------------------------- |
---|
135 | |
---|
136 | IF (firstcall) THEN |
---|
137 | |
---|
138 | if (nq.gt.nqmx) then |
---|
139 | write(*,*) 'stop in co2cloud (nq.gt.nqmx)!' |
---|
140 | write(*,*) 'nq=',nq,' nqmx=',nqmx |
---|
141 | stop |
---|
142 | endif |
---|
143 | |
---|
144 | write(*,*) "co2cloud: igcm_co2=",igcm_co2 |
---|
145 | write(*,*) " igcm_co2_ice=",igcm_co2_ice |
---|
146 | |
---|
147 | write(*,*) "time subsampling for microphysic ?" |
---|
148 | #ifdef MESOSCALE |
---|
149 | imicro = 2 |
---|
150 | #else |
---|
151 | imicro = 30 |
---|
152 | #endif |
---|
153 | call getin("imicro",imicro) |
---|
154 | imicro=ptimestep/50. |
---|
155 | |
---|
156 | write(*,*)"imicro = ",imicro |
---|
157 | |
---|
158 | microtimestep = ptimestep/real(imicro) |
---|
159 | write(*,*)"Physical timestep is",ptimestep |
---|
160 | write(*,*)"CO2 Microphysics timestep is",microtimestep |
---|
161 | |
---|
162 | firstcall=.false. |
---|
163 | ENDIF ! of IF (firstcall) |
---|
164 | |
---|
165 | c-----Initialization |
---|
166 | beta=0.85 |
---|
167 | subpdq(1:ngrid,1:nlay,1:nq) = 0 |
---|
168 | subpdt(1:ngrid,1:nlay) = 0 |
---|
169 | subpdqcloudco2(1:ngrid,1:nlay,1:nq) = 0 |
---|
170 | subpdtcloudco2(1:ngrid,1:nlay) = 0 |
---|
171 | |
---|
172 | |
---|
173 | wq(:,:)=0 |
---|
174 | ! default value if no ice |
---|
175 | rhocloudco2(1:ngrid,1:nlay) = rho_dust |
---|
176 | rhocloudco2t(1:ngrid,1:nlay) = rho_dust |
---|
177 | epaisseur(1:ngrid,1:nlay)=0 |
---|
178 | masse(1:ngrid,1:nlay)=0 |
---|
179 | |
---|
180 | tempo_traceur_t(1:ngrid,1:nlay)=0 |
---|
181 | tempo_traceurs(1:ngrid,1:nlay,1:nq)=0 |
---|
182 | sav_trac(1:ngrid,1:nlay,1:nq)=0 |
---|
183 | pdqsed(1:ngrid,1:nlay,1:nq)=0 |
---|
184 | |
---|
185 | do l=1,nlay |
---|
186 | do ig=1, ngrid |
---|
187 | masse(ig,l)=(pplev(ig,l) - pplev(ig,l+1)) /g |
---|
188 | epaisseur(ig,l)= zlev(ig,l+1) - zlev(ig,l) |
---|
189 | |
---|
190 | enddo |
---|
191 | enddo |
---|
192 | |
---|
193 | |
---|
194 | |
---|
195 | |
---|
196 | |
---|
197 | |
---|
198 | |
---|
199 | |
---|
200 | c------------------------------------------------------------------- |
---|
201 | c 1. Tendencies: |
---|
202 | c------------------ |
---|
203 | |
---|
204 | |
---|
205 | |
---|
206 | c------------------------------------------------------------------ |
---|
207 | c Time subsampling for microphysics |
---|
208 | c------------------------------------------------------------------ |
---|
209 | DO microstep=1,imicro |
---|
210 | c------ Temperature tendency subpdt |
---|
211 | ! Each microtimestep we give the cloud scheme a stepped entry subpdt instead of pdt |
---|
212 | ! If imicro=1 subpdt is the same as pdt |
---|
213 | DO l=1,nlay |
---|
214 | DO ig=1,ngrid |
---|
215 | c tempo_traceur_t(ig,l)=tempo_traceur_t(ig,l) |
---|
216 | c & + subpdtcloudco2(ig,l) |
---|
217 | !write(*,*) 'T micro= ', tempo_traceur_t(ig,l) |
---|
218 | c tempo_traceurs(ig,l,:)=tempo_traceurs(ig,l,:) |
---|
219 | c & +subpdqcloudco2(ig,l,:) |
---|
220 | |
---|
221 | subpdt(ig,l) = subpdt(ig,l) |
---|
222 | & + pdt(ig,l) ! At each micro timestep we add pdt in order to have a stepped entry |
---|
223 | |
---|
224 | subpdq(ig,l,igcm_dust_mass) = |
---|
225 | & subpdq(ig,l,igcm_dust_mass) |
---|
226 | & + pdq(ig,l,igcm_dust_mass) |
---|
227 | |
---|
228 | subpdq(ig,l,igcm_dust_number) = |
---|
229 | & subpdq(ig,l,igcm_dust_number) |
---|
230 | & + pdq(ig,l,igcm_dust_number) |
---|
231 | |
---|
232 | subpdq(ig,l,igcm_ccnco2_mass) = |
---|
233 | & subpdq(ig,l,igcm_ccnco2_mass) |
---|
234 | & + pdq(ig,l,igcm_ccnco2_mass) |
---|
235 | |
---|
236 | subpdq(ig,l,igcm_ccnco2_number) = |
---|
237 | & subpdq(ig,l,igcm_ccnco2_number) |
---|
238 | & + pdq(ig,l,igcm_ccnco2_number) |
---|
239 | |
---|
240 | subpdq(ig,l,igcm_co2_ice) = |
---|
241 | & subpdq(ig,l,igcm_co2_ice) |
---|
242 | & + pdq(ig,l,igcm_co2_ice) |
---|
243 | subpdq(ig,l,igcm_co2) = |
---|
244 | & subpdq(ig,l,igcm_co2) |
---|
245 | & + pdq(ig,l,igcm_co2) |
---|
246 | |
---|
247 | tempo_traceur_t(ig,l)= pt(ig,l)+subpdt(ig,l)*microtimestep |
---|
248 | tempo_traceurs(ig,l,:)= pq(ig,l,:)+subpdq(ig,l,:) |
---|
249 | & *microtimestep |
---|
250 | !Stepped entry for sedimentation |
---|
251 | ENDDO |
---|
252 | ENDDO |
---|
253 | |
---|
254 | !RSEDCLOUD AND RICECO2 HERE |
---|
255 | |
---|
256 | DO l=1, nlay |
---|
257 | DO ig=1,ngrid |
---|
258 | Niceco2=tempo_traceurs(ig,l,igcm_co2_ice) |
---|
259 | Nccnco2=max(tempo_traceurs(ig,l,igcm_ccnco2_number), |
---|
260 | & 1.e-30) |
---|
261 | Qccnco2=max(tempo_traceurs(ig,l,igcm_ccnco2_mass), |
---|
262 | & 1.e-30) |
---|
263 | mdustJA= tempo_traceurs(ig,l,igcm_dust_mass) |
---|
264 | ndustJA=tempo_traceurs(ig,l,igcm_dust_number) |
---|
265 | if ((ndustJA .lt. tauscaling(ig)) .or. (mdustJA .lt. |
---|
266 | & 1.e-30 *tauscaling(ig))) then |
---|
267 | rdust(ig,l)=1.e-10 |
---|
268 | else |
---|
269 | rdust(ig,l)=(3./4./pi/2500.*mdustJA/ndustJA)**(1./3.) |
---|
270 | rdust(ig,l)=min(rdust(ig,l),5.e-6) |
---|
271 | rdust(ig,l)=max(rdust(ig,l),1.e-9) |
---|
272 | end if |
---|
273 | rhocloudco2t(ig,l) = (Niceco2 *rho_ice_co2 |
---|
274 | & + Qccnco2*rho_dust) |
---|
275 | & / (Niceco2 + Qccnco2) |
---|
276 | riceco2(ig,l)= Niceco2*3.0/ |
---|
277 | & (4.0*rho_ice_co2*pi*Nccnco2) |
---|
278 | & +rdust(ig,l)*rdust(ig,l)*rdust(ig,l) |
---|
279 | riceco2(ig,l)=riceco2(ig,l)**(1.0/3.0) |
---|
280 | write(*,*) "in co2clouds, rice = ",riceco2(ig,l) |
---|
281 | write(*,*) "in co2clouds, rho = ",rhocloudco2t(ig,l) |
---|
282 | |
---|
283 | call updaterice_microCO2(Niceco2,Qccnco2,Nccnco2, |
---|
284 | & tauscaling(ig),riceco2(ig,l),rhocloudco2t(ig,l)) |
---|
285 | write(*,*) "in co2clouds, rice update = ",riceco2(ig,l) |
---|
286 | write(*,*) "in co2clouds, rho update = " |
---|
287 | & ,rhocloudco2t(ig,l) |
---|
288 | |
---|
289 | rsedcloudco2(ig,l)=max(riceco2(ig,l)* |
---|
290 | & (1.+nuiceco2_sed)*(1.+nuiceco2_sed)*(1.+nuiceco2_sed), |
---|
291 | & rdust(ig,l)) |
---|
292 | rsedcloudco2(ig,l)=min(rsedcloudco2(ig,l),5.e-4) |
---|
293 | write(*,*) 'Rsedcloud = ',rsedcloudco2(ig,l) |
---|
294 | !write(*,*) 'Rhocloudco2 = ',rhocloudco2t(ig,l) |
---|
295 | |
---|
296 | ENDDO |
---|
297 | ENDDO |
---|
298 | |
---|
299 | ! Gravitational sedimentation |
---|
300 | |
---|
301 | ! sedimentation computed from radius computed from q in module radii_mod |
---|
302 | sav_trac(:,:,igcm_co2_ice)=tempo_traceurs(:,:,igcm_co2_ice) |
---|
303 | sav_trac(:,:,igcm_ccnco2_mass)= |
---|
304 | & tempo_traceurs(:,:,igcm_ccnco2_mass) |
---|
305 | sav_trac(:,:,igcm_ccnco2_number)= |
---|
306 | & tempo_traceurs(:,:,igcm_ccnco2_number) |
---|
307 | |
---|
308 | call newsedim(ngrid,nlay,ngrid*nlay,ngrid*nlay, |
---|
309 | & microtimestep,pplev,masse,epaisseur,tempo_traceur_t, |
---|
310 | & rsedcloudco2,rhocloudco2t, |
---|
311 | & tempo_traceurs(:,:,igcm_co2_ice),wq,beta) ! 3 traceurs |
---|
312 | |
---|
313 | ! sedim at the surface of co2 ice |
---|
314 | do ig=1,ngrid |
---|
315 | pdqs_sedco2(ig)=pdqs_sedco2(ig)+ wq(ig,1) |
---|
316 | end do |
---|
317 | |
---|
318 | call newsedim(ngrid,nlay,ngrid*nlay,ngrid*nlay, |
---|
319 | & microtimestep,pplev,masse,epaisseur,tempo_traceur_t, |
---|
320 | & rsedcloudco2,rhocloudco2t, |
---|
321 | & tempo_traceurs(:,:,igcm_ccnco2_mass),wq,beta) |
---|
322 | |
---|
323 | call newsedim(ngrid,nlay,ngrid*nlay,ngrid*nlay, |
---|
324 | & microtimestep,pplev,masse,epaisseur,tempo_traceur_t, |
---|
325 | & rsedcloudco2,rhocloudco2t, |
---|
326 | & tempo_traceurs(:,:,igcm_ccnco2_number),wq,beta) |
---|
327 | |
---|
328 | |
---|
329 | DO l = 1, nlay |
---|
330 | DO ig=1,ngrid |
---|
331 | pdqsed(ig,l,igcm_ccnco2_mass)= |
---|
332 | & (tempo_traceurs(ig,l,igcm_ccnco2_mass)- |
---|
333 | & sav_trac(ig,l,igcm_ccnco2_mass))/microtimestep |
---|
334 | pdqsed(ig,l,igcm_ccnco2_number)= |
---|
335 | & (tempo_traceurs(ig,l,igcm_ccnco2_number)- |
---|
336 | & sav_trac(ig,l,igcm_ccnco2_number))/microtimestep |
---|
337 | pdqsed(ig,l,igcm_co2_ice)= |
---|
338 | & (tempo_traceurs(ig,l,igcm_co2_ice)- |
---|
339 | & sav_trac(ig,l,igcm_co2_ice))/microtimestep |
---|
340 | ENDDO |
---|
341 | ENDDO |
---|
342 | !pdqsed est la tendance due a la sedimentation |
---|
343 | |
---|
344 | DO l = 1, nlay |
---|
345 | DO ig=1,ngrid |
---|
346 | pdqsed(ig,l,igcm_ccnco2_mass)= |
---|
347 | & (tempo_traceurs(ig,l,igcm_ccnco2_mass)- |
---|
348 | & sav_trac(ig,l,igcm_ccnco2_mass))/microtimestep |
---|
349 | pdqsed(ig,l,igcm_ccnco2_number)= |
---|
350 | & (tempo_traceurs(ig,l,igcm_ccnco2_number)- |
---|
351 | & sav_trac(ig,l,igcm_ccnco2_number))/microtimestep |
---|
352 | pdqsed(ig,l,igcm_co2_ice)= |
---|
353 | & (tempo_traceurs(ig,l,igcm_co2_ice)- |
---|
354 | & sav_trac(ig,l,igcm_co2_ice))/microtimestep |
---|
355 | ENDDO |
---|
356 | ENDDO |
---|
357 | !pdqsed est la tendance due a la sedimentation |
---|
358 | DO l=1,nlay |
---|
359 | DO ig=1,ngrid |
---|
360 | subpdq(ig,l,igcm_ccnco2_mass) = |
---|
361 | & subpdq(ig,l,igcm_ccnco2_mass) |
---|
362 | & +pdqsed(ig,l,igcm_ccnco2_mass) |
---|
363 | |
---|
364 | subpdq(ig,l,igcm_ccnco2_number) = |
---|
365 | & subpdq(ig,l,igcm_ccnco2_number) |
---|
366 | & +pdqsed(ig,l,igcm_ccnco2_number) |
---|
367 | |
---|
368 | subpdq(ig,l,igcm_co2_ice) = |
---|
369 | & subpdq(ig,l,igcm_co2_ice) |
---|
370 | & +pdqsed(ig,l,igcm_co2_ice) |
---|
371 | ENDDO |
---|
372 | ENDDO |
---|
373 | c------------------------------------------------------------------- |
---|
374 | c 2. Main call to the different cloud schemes: |
---|
375 | c------------------------------------------------ |
---|
376 | IF (microphysco2) THEN |
---|
377 | CALL improvedCO2clouds(ngrid,nlay,microtimestep, |
---|
378 | & pplay,pt,subpdt, |
---|
379 | & pq,subpdq,subpdqcloudco2,subpdtcloudco2, |
---|
380 | & nq,tauscaling) |
---|
381 | |
---|
382 | ELSE |
---|
383 | |
---|
384 | write(*,*) ' no simpleCO2clouds procedure: STOP' ! listo |
---|
385 | STOP |
---|
386 | |
---|
387 | c CALL simpleclouds(ngrid,nlay,microtimestep, ! for water-ice clouds |
---|
388 | c & pplay,pzlay,pt,subpdt, |
---|
389 | c & pq,subpdq,subpdqcloud,subpdtcloud, |
---|
390 | c & nq,tau,riceco2) |
---|
391 | ENDIF |
---|
392 | |
---|
393 | |
---|
394 | c------------------------------------------------------------------- |
---|
395 | c 3. Updating tendencies after cloud scheme: |
---|
396 | c----------------------------------------------- |
---|
397 | |
---|
398 | c IF (microphysco2) THEN |
---|
399 | DO l=1,nlay |
---|
400 | DO ig=1,ngrid |
---|
401 | subpdq(ig,l,igcm_dust_mass) = |
---|
402 | & subpdq(ig,l,igcm_dust_mass) |
---|
403 | & + subpdqcloudco2(ig,l,igcm_dust_mass) |
---|
404 | |
---|
405 | subpdq(ig,l,igcm_dust_number) = |
---|
406 | & subpdq(ig,l,igcm_dust_number) |
---|
407 | & + subpdqcloudco2(ig,l,igcm_dust_number) |
---|
408 | |
---|
409 | subpdq(ig,l,igcm_ccnco2_mass) = |
---|
410 | & subpdq(ig,l,igcm_ccnco2_mass) |
---|
411 | & + subpdqcloudco2(ig,l,igcm_ccnco2_mass) |
---|
412 | c & +pdqsed(ig,l,igcm_ccnco2_mass) |
---|
413 | |
---|
414 | subpdq(ig,l,igcm_ccnco2_number) = |
---|
415 | & subpdq(ig,l,igcm_ccnco2_number) |
---|
416 | & + subpdqcloudco2(ig,l,igcm_ccnco2_number) |
---|
417 | c & +pdqsed(ig,l,igcm_ccnco2_number) |
---|
418 | |
---|
419 | subpdq(ig,l,igcm_co2_ice) = |
---|
420 | & subpdq(ig,l,igcm_co2_ice) |
---|
421 | & + subpdqcloudco2(ig,l,igcm_co2_ice) |
---|
422 | c & +pdqsed(ig,l,igcm_co2_ice) |
---|
423 | |
---|
424 | subpdq(ig,l,igcm_co2) = |
---|
425 | & subpdq(ig,l,igcm_co2) |
---|
426 | & + subpdqcloudco2(ig,l,igcm_co2) |
---|
427 | ENDDO |
---|
428 | ENDDO |
---|
429 | |
---|
430 | |
---|
431 | !ici |
---|
432 | ! call WRITEdiagfi(ngrid,"co2cloud000","co2 traceur","kg/kg",1, |
---|
433 | ! & pq(1,:,igcm_co2_ice) + ptimestep* |
---|
434 | ! & ( subpdq(1,:,igcm_co2_ice))) |
---|
435 | |
---|
436 | |
---|
437 | IF (activice) THEN |
---|
438 | DO l=1,nlay |
---|
439 | DO ig=1,ngrid |
---|
440 | subpdt(ig,l) = |
---|
441 | & subpdt(ig,l) + subpdtcloudco2(ig,l) |
---|
442 | ENDDO |
---|
443 | ENDDO |
---|
444 | ENDIF |
---|
445 | |
---|
446 | |
---|
447 | ENDDO ! of DO microstep=1,imicro |
---|
448 | |
---|
449 | c------------------------------------------------------------------- |
---|
450 | c 6. Compute final tendencies after time loop: |
---|
451 | c------------------------------------------------ |
---|
452 | c CO2 flux at surface (kg.m-2.s-1) |
---|
453 | do ig=1,ngrid |
---|
454 | pdqs_sedco2(ig)=pdqs_sedco2(ig)/ptimestep |
---|
455 | enddo |
---|
456 | |
---|
457 | c------ Temperature tendency pdtcloud |
---|
458 | DO l=1,nlay |
---|
459 | DO ig=1,ngrid |
---|
460 | pdtcloudco2(ig,l) = |
---|
461 | & subpdt(ig,l)/imicro-pdt(ig,l) |
---|
462 | ENDDO |
---|
463 | ENDDO |
---|
464 | |
---|
465 | c------ Tracers tendencies pdqcloud |
---|
466 | DO l=1,nlay |
---|
467 | DO ig=1,ngrid |
---|
468 | |
---|
469 | pdqcloudco2(ig,l,igcm_co2_ice) = |
---|
470 | & subpdq(ig,l,igcm_co2_ice)/imicro |
---|
471 | & - pdq(ig,l,igcm_co2_ice) |
---|
472 | pdqcloudco2(ig,l,igcm_co2) = |
---|
473 | & subpdq(ig,l,igcm_co2)/imicro |
---|
474 | & - pdq(ig,l,igcm_co2) |
---|
475 | ENDDO |
---|
476 | ENDDO |
---|
477 | |
---|
478 | |
---|
479 | ! call WRITEdiagfi(ngrid,"co2cloud00","co2 traceur","kg/kg",1, |
---|
480 | ! & pq(1,:,igcm_co2_ice) + ptimestep* |
---|
481 | ! & (pdq(1,:,igcm_co2_ice) + pdqcloudco2(1,:,igcm_co2_ice))) |
---|
482 | |
---|
483 | |
---|
484 | IF(microphysco2) THEN |
---|
485 | DO l=1,nlay |
---|
486 | DO ig=1,ngrid |
---|
487 | pdqcloudco2(ig,l,igcm_ccnco2_mass) = |
---|
488 | & subpdq(ig,l,igcm_ccnco2_mass)/imicro |
---|
489 | & - pdq(ig,l,igcm_ccnco2_mass) |
---|
490 | pdqcloudco2(ig,l,igcm_ccnco2_number) = |
---|
491 | & subpdq(ig,l,igcm_ccnco2_number)/imicro |
---|
492 | & - pdq(ig,l,igcm_ccnco2_number) |
---|
493 | ENDDO |
---|
494 | ENDDO |
---|
495 | ENDIF |
---|
496 | |
---|
497 | |
---|
498 | IF(scavenging) THEN |
---|
499 | DO l=1,nlay |
---|
500 | DO ig=1,ngrid |
---|
501 | pdqcloudco2(ig,l,igcm_dust_mass) = |
---|
502 | & subpdq(ig,l,igcm_dust_mass)/real(imicro) |
---|
503 | & - pdq(ig,l,igcm_dust_mass) |
---|
504 | pdqcloudco2(ig,l,igcm_dust_number) = |
---|
505 | & subpdq(ig,l,igcm_dust_number)/real(imicro) |
---|
506 | & - pdq(ig,l,igcm_dust_number) |
---|
507 | ENDDO |
---|
508 | ENDDO |
---|
509 | ENDIF |
---|
510 | |
---|
511 | c ENDIF |
---|
512 | c------- Due to stepped entry, other processes tendencies can add up to negative values |
---|
513 | c------- Therefore, enforce positive values and conserve mass |
---|
514 | |
---|
515 | |
---|
516 | IF(microphysco2) THEN |
---|
517 | DO l=1,nlay |
---|
518 | DO ig=1,ngrid |
---|
519 | IF ((pq(ig,l,igcm_ccnco2_number) + |
---|
520 | & ptimestep* (pdq(ig,l,igcm_ccnco2_number) + |
---|
521 | & pdqcloudco2(ig,l,igcm_ccnco2_number)) |
---|
522 | & .lt. 0) |
---|
523 | & .or. (pq(ig,l,igcm_ccnco2_mass) + |
---|
524 | & ptimestep* (pdq(ig,l,igcm_ccnco2_mass) + |
---|
525 | & pdqcloudco2(ig,l,igcm_ccnco2_mass)) |
---|
526 | & .lt. 0)) THEN |
---|
527 | |
---|
528 | pdqcloudco2(ig,l,igcm_ccnco2_number) = |
---|
529 | & - pq(ig,l,igcm_ccnco2_number)/ptimestep |
---|
530 | & - pdq(ig,l,igcm_ccnco2_number) +0 |
---|
531 | |
---|
532 | pdqcloudco2(ig,l,igcm_dust_number) = |
---|
533 | & -pdqcloudco2(ig,l,igcm_ccnco2_number) |
---|
534 | |
---|
535 | pdqcloudco2(ig,l,igcm_ccnco2_mass) = |
---|
536 | & - pq(ig,l,igcm_ccnco2_mass)/ptimestep |
---|
537 | & - pdq(ig,l,igcm_ccnco2_mass)+0 |
---|
538 | |
---|
539 | pdqcloudco2(ig,l,igcm_dust_mass) = |
---|
540 | & -pdqcloudco2(ig,l,igcm_ccnco2_mass) |
---|
541 | |
---|
542 | ENDIF |
---|
543 | ENDDO |
---|
544 | ENDDO |
---|
545 | ENDIF |
---|
546 | |
---|
547 | |
---|
548 | IF(scavenging) THEN |
---|
549 | DO l=1,nlay |
---|
550 | DO ig=1,ngrid |
---|
551 | IF ( (pq(ig,l,igcm_dust_number) + |
---|
552 | & ptimestep* (pdq(ig,l,igcm_dust_number) + |
---|
553 | & pdqcloudco2(ig,l,igcm_dust_number)) .lt. 0.) |
---|
554 | & .or. (pq(ig,l,igcm_dust_mass)+ |
---|
555 | & ptimestep* (pdq(ig,l,igcm_dust_mass) + |
---|
556 | & pdqcloudco2(ig,l,igcm_dust_mass)) |
---|
557 | & .lt. 0.)) then |
---|
558 | |
---|
559 | pdqcloudco2(ig,l,igcm_dust_number) = |
---|
560 | & - pq(ig,l,igcm_dust_number)/ptimestep |
---|
561 | & - pdq(ig,l,igcm_dust_number)+0 |
---|
562 | |
---|
563 | pdqcloudco2(ig,l,igcm_ccnco2_number) = |
---|
564 | & -pdqcloudco2(ig,l,igcm_dust_number) |
---|
565 | |
---|
566 | pdqcloudco2(ig,l,igcm_dust_mass) = |
---|
567 | & - pq(ig,l,igcm_dust_mass)/ptimestep |
---|
568 | & - pdq(ig,l,igcm_dust_mass) +0 |
---|
569 | |
---|
570 | pdqcloudco2(ig,l,igcm_ccnco2_mass) = |
---|
571 | & -pdqcloudco2(ig,l,igcm_dust_mass) |
---|
572 | ENDIF |
---|
573 | ENDDO |
---|
574 | ENDDO |
---|
575 | ENDIF !pq+ptime*(pdq+pdqc)=1 ! pdqc=1-pq/ptime-pdq |
---|
576 | |
---|
577 | |
---|
578 | DO l=1,nlay |
---|
579 | DO ig=1,ngrid |
---|
580 | IF (pq(ig,l,igcm_co2_ice) + ptimestep* |
---|
581 | & (pdq(ig,l,igcm_co2_ice) + pdqcloudco2(ig,l,igcm_co2_ice)) |
---|
582 | & .lt. 1.e-25) THEN |
---|
583 | pdqcloudco2(ig,l,igcm_co2_ice) = |
---|
584 | & - pq(ig,l,igcm_co2_ice)/ptimestep - pdq(ig,l,igcm_co2_ice) |
---|
585 | & +1.e-25 |
---|
586 | pdqcloudco2(ig,l,igcm_co2) = -pdqcloudco2(ig,l,igcm_co2_ice) |
---|
587 | ENDIF |
---|
588 | ENDDO |
---|
589 | ENDDO |
---|
590 | |
---|
591 | |
---|
592 | |
---|
593 | |
---|
594 | c------Update the ice and dust particle size "riceco2" for output or photochemistry |
---|
595 | c------Only rsedcloudco2 is used for the co2 (cloud) cycle |
---|
596 | |
---|
597 | IF(scavenging) THEN |
---|
598 | DO l=1, nlay |
---|
599 | DO ig=1,ngrid |
---|
600 | |
---|
601 | c call updaterdust( |
---|
602 | c & pq(ig,l,igcm_dust_mass) + ! dust mass |
---|
603 | c & (pdq(ig,l,igcm_dust_mass) + ! dust mass |
---|
604 | c & pdqcloudco2(ig,l,igcm_dust_mass))*ptimestep, ! dust mass |
---|
605 | c & pq(ig,l,igcm_dust_number) + ! dust number |
---|
606 | c & (pdq(ig,l,igcm_dust_number) + ! dust number |
---|
607 | c & pdqcloudco2(ig,l,igcm_dust_number))*ptimestep, ! dust number |
---|
608 | c & rdust(ig,l)) |
---|
609 | c write(*,*) "in co2clouds, rdust(ig,l)= ",rdust(ig,l) |
---|
610 | mdustJA= pq(ig,l,igcm_dust_mass) + |
---|
611 | & (pdq(ig,l,igcm_dust_mass) + |
---|
612 | & pdqcloudco2(ig,l,igcm_dust_mass))*ptimestep |
---|
613 | ndustJA=pq(ig,l,igcm_dust_number) + |
---|
614 | & (pdq(ig,l,igcm_dust_number) + |
---|
615 | & pdqcloudco2(ig,l,igcm_dust_number))*ptimestep |
---|
616 | if ((ndustJA .lt. tauscaling(ig)) .or. (mdustJA .lt. |
---|
617 | & 1.e-30 *tauscaling(ig))) then |
---|
618 | rdust(ig,l)=1.e-10 |
---|
619 | else |
---|
620 | rdust(ig,l)=(3./4./pi/2500.*mdustJA/ndustJA)**(1./3.) |
---|
621 | rdust(ig,l)=min(rdust(ig,l),5.e-4) |
---|
622 | rdust(ig,l)=max(rdust(ig,l),1.e-10) |
---|
623 | endif |
---|
624 | ENDDO |
---|
625 | ENDDO |
---|
626 | ENDIF |
---|
627 | |
---|
628 | |
---|
629 | IF(microphysco2) THEN |
---|
630 | |
---|
631 | DO l=1, nlay |
---|
632 | DO ig=1,ngrid |
---|
633 | |
---|
634 | c call updaterice_microco2( |
---|
635 | c & pq(ig,l,igcm_co2_ice) + ! ice mass |
---|
636 | c & (pdq(ig,l,igcm_co2_ice) + ! ice mass |
---|
637 | c & pdqcloudco2(ig,l,igcm_co2_ice))*ptimestep, ! ice mass |
---|
638 | c & pq(ig,l,igcm_ccnco2_mass) + ! ccn mass |
---|
639 | c & (pdq(ig,l,igcm_ccnco2_mass) + ! ccn mass |
---|
640 | c & pdqcloudco2(ig,l,igcm_ccnco2_mass))*ptimestep, ! ccn mass |
---|
641 | c & pq(ig,l,igcm_ccnco2_number) + ! ccn number |
---|
642 | c & (pdq(ig,l,igcm_ccnco2_number) + ! ccn number |
---|
643 | c & pdqcloudco2(ig,l,igcm_ccnco2_number))*ptimestep, ! ccn number |
---|
644 | c & tauscaling(ig),riceco2(ig,l),rhocloudco2(ig,l)) |
---|
645 | c write(*,*) "in co2clouds, riceco2(ig,l)= ",riceco2(ig,l) |
---|
646 | |
---|
647 | |
---|
648 | Niceco2=pq(ig,l,igcm_co2_ice) + |
---|
649 | & (pdq(ig,l,igcm_co2_ice) + |
---|
650 | & pdqcloudco2(ig,l,igcm_co2_ice))*ptimestep |
---|
651 | Nccnco2=max((pq(ig,l,igcm_ccnco2_number) + |
---|
652 | & (pdq(ig,l,igcm_ccnco2_number) + |
---|
653 | & pdqcloudco2(ig,l,igcm_ccnco2_number))*ptimestep)* |
---|
654 | & tauscaling(ig),1.e-30) |
---|
655 | Qccnco2=max((pq(ig,l,igcm_ccnco2_mass) + |
---|
656 | & (pdq(ig,l,igcm_ccnco2_mass) + |
---|
657 | & pdqcloudco2(ig,l,igcm_ccnco2_mass))*ptimestep)* |
---|
658 | & tauscaling(ig),1.e-30) |
---|
659 | rhocloudco2t(ig,l) = (Niceco2 *rho_ice_co2 + Qccnco2*rho_dust) |
---|
660 | & / (Niceco2 + Qccnco2) |
---|
661 | c rhocloudco2(ig,l) = min(max(rhocloudco2t,rho_ice_co2),rho_dust) |
---|
662 | |
---|
663 | c write(*,*) "test, nccnco2 =",nccnco22 |
---|
664 | |
---|
665 | |
---|
666 | riceco2(ig,l)= Niceco2*3.0/ |
---|
667 | & (4.0*rho_ice_co2*pi*Nccnco2) |
---|
668 | & +rdust(ig,l)*rdust(ig,l)*rdust(ig,l) |
---|
669 | |
---|
670 | riceco2(ig,l)=riceco2(ig,l)**(1.0/3.0) |
---|
671 | write(*,*) "In co2cloud, after loop, riceco2 =",riceco2(ig,l) |
---|
672 | write(*,*) "In co2cloud, after loop, rhoco2 =" |
---|
673 | & ,rhocloudco2t(ig,l) |
---|
674 | |
---|
675 | call updaterice_microCO2(Niceco2,Qccnco2,Nccnco2, |
---|
676 | & tauscaling(ig),riceco2(ig,l),rhocloudco2t(ig,l)) |
---|
677 | |
---|
678 | write(*,*) "In co2cloud, after loop and update, riceco2 =" |
---|
679 | & ,riceco2(ig,l) |
---|
680 | write(*,*) "In co2cloud, after loop and update, rhoco2 =" |
---|
681 | & ,rhocloudco2t(ig,l) |
---|
682 | |
---|
683 | if ( Niceco2 |
---|
684 | & .le. 1.e-23 .or. riceco2(ig,l) .le. 1.e-10 .or. |
---|
685 | & riceco2(ig,l) .ge. 4.999e-4) then ! .or. riceco2(ig,l) .gt. 1.e-4 ) then |
---|
686 | riceco2(ig,l)=0. |
---|
687 | |
---|
688 | !NO CLOUD : RESET TRACER AND CONSERVE MASS |
---|
689 | pdqcloudco2(ig,l,igcm_co2)= pq(ig,l,igcm_co2_ice) |
---|
690 | & /ptimestep+pdq(ig,l,igcm_co2_ice) |
---|
691 | |
---|
692 | pdqcloudco2(ig,l,igcm_co2_ice)=-pq(ig,l,igcm_co2_ice) |
---|
693 | & /ptimestep-pdq(ig,l,igcm_co2_ice) |
---|
694 | |
---|
695 | pdqcloudco2(ig,l,igcm_ccnco2_mass)= |
---|
696 | & -pq(ig,l,igcm_ccnco2_mass) |
---|
697 | & /ptimestep-pdq(ig,l,igcm_ccnco2_mass) |
---|
698 | |
---|
699 | pdqcloudco2(ig,l,igcm_ccnco2_number)= |
---|
700 | & -pq(ig,l,igcm_ccnco2_number) |
---|
701 | & /ptimestep-pdq(ig,l,igcm_ccnco2_number) |
---|
702 | |
---|
703 | pdqcloudco2(ig,l,igcm_dust_number)= |
---|
704 | & pq(ig,l,igcm_ccnco2_number) |
---|
705 | & /ptimestep+pdq(ig,l,igcm_ccnco2_number) |
---|
706 | |
---|
707 | pdqcloudco2(ig,l,igcm_dust_mass)= |
---|
708 | & pq(ig,l,igcm_ccnco2_mass) |
---|
709 | & /ptimestep+pdq(ig,l,igcm_ccnco2_mass) |
---|
710 | c$$$ |
---|
711 | |
---|
712 | c$$$ |
---|
713 | endif |
---|
714 | |
---|
715 | c write(*,*) "in co2clouds, riceco2(ig,l)v2= ",riceco2(ig,l) |
---|
716 | |
---|
717 | ENDDO |
---|
718 | ENDDO |
---|
719 | |
---|
720 | ELSE ! no microphys ! not of concern for co2 clouds - listo |
---|
721 | |
---|
722 | ENDIF ! of IF(microphysco2) |
---|
723 | |
---|
724 | |
---|
725 | c TO CHEK for relevancy - listo |
---|
726 | |
---|
727 | c A correction if a lot of subliming CO2 fills the 1st layer FF04/2005 |
---|
728 | c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
---|
729 | c Then that should not affect the ice particle radius |
---|
730 | do ig=1,ngrid |
---|
731 | if(pdpsrf(ig)*ptimestep.gt.0.9*(pplev(ig,1)-pplev(ig,2)))then |
---|
732 | if(pdpsrf(ig)*ptimestep.gt.0.9*(pplev(ig,1)-pplev(ig,3))) |
---|
733 | & riceco2(ig,2)=riceco2(ig,3) |
---|
734 | riceco2(ig,1)=riceco2(ig,2) |
---|
735 | end if |
---|
736 | end do |
---|
737 | |
---|
738 | |
---|
739 | DO l=1,nlay |
---|
740 | DO ig=1,ngrid |
---|
741 | rsedcloudco2(ig,l)=max(riceco2(ig,l)* |
---|
742 | & (1.+nuiceco2_sed)*(1.+nuiceco2_sed)*(1.+nuiceco2_sed), |
---|
743 | & rdust(ig,l)) |
---|
744 | rsedcloudco2(ig,l)=min(rsedcloudco2(ig,l),1.e-4) |
---|
745 | ENDDO |
---|
746 | ENDDO |
---|
747 | |
---|
748 | call co2sat(ngrid*nlay,pt,pplay,zqsatco2) |
---|
749 | do ig=1,ngrid |
---|
750 | do l=1,nlay |
---|
751 | satuco2(ig,l) = pq(ig,l,igcm_co2)* |
---|
752 | & (mmean(ig,l)/44.01)*pplay(ig,l)/zqsatco2(ig,l) |
---|
753 | |
---|
754 | write(*,*) "In CO2 pt,sat ",pt(ig,l),satuco2(ig,l) |
---|
755 | enddo |
---|
756 | enddo |
---|
757 | c call WRITEDIAGFI(ngrid,"satuco2","vap in satu","kg/kg",1, |
---|
758 | c & satuco2) |
---|
759 | c call WRITEdiagfi(ngrid,"riceco2","ice radius","m" |
---|
760 | c & ,1,riceco2) |
---|
761 | ! or output in diagfi.nc (for testphys1d) |
---|
762 | c call WRITEDIAGFI(ngrid,'ps','Surface pressure','Pa',0,ps) |
---|
763 | c call WRITEDIAGFI(ngrid,'temp','Temperature ', |
---|
764 | c & 'K JA',1,pt) |
---|
765 | |
---|
766 | call WRITEdiagfi(ngrid,"rsedcloudco2","rsed co2","m",1, |
---|
767 | & rsedcloudco2) |
---|
768 | |
---|
769 | ! used for rad. transfer calculations |
---|
770 | ! nuice is constant because a lognormal distribution is prescribed |
---|
771 | c nuice(1:ngrid,1:nlay)=nuice_ref |
---|
772 | |
---|
773 | |
---|
774 | |
---|
775 | c======================================================================= |
---|
776 | |
---|
777 | END |
---|
778 | |
---|