| 1 | SUBROUTINE phytrac (firstcall,lastcall, |
|---|
| 2 | . nqmax,nmicro,ptimestep,appkim,dtkim, |
|---|
| 3 | . pplev,pplay,delp,ptemp,pmu0,pfract,pdecli, |
|---|
| 4 | . lonsol,tr_seri,d_tr_mph,d_tr_kim) |
|---|
| 5 | |
|---|
| 6 | USE infotrac |
|---|
| 7 | |
|---|
| 8 | IMPLICIT none |
|---|
| 9 | |
|---|
| 10 | c====================================================================== |
|---|
| 11 | c S. Lebonnois, mai 2008 |
|---|
| 12 | c |
|---|
| 13 | c Arguments: |
|---|
| 14 | c |
|---|
| 15 | c firstcall----input-L-variable logique indiquant le premier passage |
|---|
| 16 | c lastcall-----input-L-variable logique indiquant le dernier passage |
|---|
| 17 | c nqmax--------input-I-nombre de traceurs (total) |
|---|
| 18 | c nmicro-------input-I-nombre de traceurs microphysiques !! doivent etre toujours en premiers!! |
|---|
| 19 | c ptimestep----input-R-pas d'integration pour la physique (seconde) |
|---|
| 20 | c appkim-------input-I-appel a la chimie |
|---|
| 21 | c dtkim--------input-R-pas de temps chimique (seconde) |
|---|
| 22 | c pplev--------input-R-pression pour chaque inter-couche (en Pa) |
|---|
| 23 | c pplay--------input-R-pression pour chaque couche (en Pa) |
|---|
| 24 | c delp---------input-R-epaisseur d'une couche (en Pa) |
|---|
| 25 | c ptemp--------input-R-temperature (K) |
|---|
| 26 | c pmu0---------input-R-cos angle zenithal |
|---|
| 27 | c pfract-------input-R-fractional day |
|---|
| 28 | c pdecli-------input-R-declinaison en radian |
|---|
| 29 | c lonsol-------input-R-longitude solaire en radian |
|---|
| 30 | c tr_seri------input-R-mass mixing ratio traceurs (kg/kg) |
|---|
| 31 | c d_tr_mph----output-R-tendance microphysique de "qx" (kg/kg/s) |
|---|
| 32 | c d_tr_kim----output-R-tendance chimique de "qx" (kg/kg/s) |
|---|
| 33 | c====================================================================== |
|---|
| 34 | #include "dimensions.h" |
|---|
| 35 | #include "dimphy.h" |
|---|
| 36 | #include "clesphys.h" |
|---|
| 37 | #include "YOMCST.h" |
|---|
| 38 | |
|---|
| 39 | c====================================================================== |
|---|
| 40 | c Variables argument: |
|---|
| 41 | c |
|---|
| 42 | LOGICAL firstcall,lastcall |
|---|
| 43 | INTEGER nqmax,nmicro,nlat,appkim |
|---|
| 44 | REAL ptimestep,dtkim |
|---|
| 45 | REAL pplev(klon,klev+1),pplay(klon,klev+1),delp(klon,klev) |
|---|
| 46 | REAL ptemp(klon,klev) |
|---|
| 47 | REAL pmu0(klon), pfract(klon), pdecli, lonsol |
|---|
| 48 | REAL tr_seri(klon,klev,nqmax) |
|---|
| 49 | REAL d_tr_mph(klon,klev,nqmax),d_tr_kim(klon,klev,nqmax) |
|---|
| 50 | |
|---|
| 51 | c====================================================================== |
|---|
| 52 | c Local variables |
|---|
| 53 | |
|---|
| 54 | * common relatifs aux aerosols |
|---|
| 55 | REAL qaer(klon,klev,nqtot) |
|---|
| 56 | common/traceurs/qaer |
|---|
| 57 | |
|---|
| 58 | c grandeurs en moyennes zonales |
|---|
| 59 | REAL zplev(jjm+1,klev+1),zplay(jjm+1,klev) |
|---|
| 60 | REAL ztemp(jjm+1,klev),zmu0(jjm+1),zfract(jjm+1) |
|---|
| 61 | real temp_eq(klev),press_eq(klev) |
|---|
| 62 | REAL zqaer(jjm+1,klev,nqmax) ! et non nmicro... Permet nmicro=0. |
|---|
| 63 | REAL zqaer0(jjm+1,klev,nqmax) |
|---|
| 64 | REAL pdqmfi(jjm+1,klev,nqmax) |
|---|
| 65 | REAL ychim(jjm+1,klev,nqmax-nmicro) |
|---|
| 66 | REAL qysat(klev,nqtot) ! dim nqtot, mais en fait nqmax-nmicro (save...) |
|---|
| 67 | REAL pdyfi(jjm+1,klev,nqtot) ! dim nqtot, mais en fait nqmax-nmicro (save...) |
|---|
| 68 | character*10 nomqy(nqmax-nmicro+1) |
|---|
| 69 | integer i,j,l,iq,ig0 |
|---|
| 70 | |
|---|
| 71 | c La saturation n est calculee qu une seule fois: sauvegarde qysat |
|---|
| 72 | c La chimie n est pas calculee tous les pas, il faut donc |
|---|
| 73 | c sauvegarder les sorties de la chimie |
|---|
| 74 | |
|---|
| 75 | SAVE pdyfi,qysat |
|---|
| 76 | |
|---|
| 77 | c====================================================================== |
|---|
| 78 | c====================================================================== |
|---|
| 79 | |
|---|
| 80 | c----------------------------------------------------------------------- |
|---|
| 81 | c convertion moyennes zonales et changement d unites pour microphy |
|---|
| 82 | c --------------------------------- |
|---|
| 83 | |
|---|
| 84 | c print*,'CONVERSION 2D ET CHANGEMENT UNITES (PHYTRAC)' |
|---|
| 85 | |
|---|
| 86 | zplev = 0.0 |
|---|
| 87 | zplay = 0.0 |
|---|
| 88 | ztemp = 0.0 |
|---|
| 89 | zqaer = 0.0 |
|---|
| 90 | ychim = 0.0 |
|---|
| 91 | zmu0 = 0.0 |
|---|
| 92 | zfract= 0.0 |
|---|
| 93 | |
|---|
| 94 | do l=1,llm+1 |
|---|
| 95 | zplev(1,l) = pplev(1,l) |
|---|
| 96 | do j=2,jjm |
|---|
| 97 | ig0=1+(j-2)*iim |
|---|
| 98 | do i=1,iim |
|---|
| 99 | zplev(j,l) = zplev(j,l) + pplev(ig0+i,l)/iim |
|---|
| 100 | enddo |
|---|
| 101 | enddo |
|---|
| 102 | zplev(jjm+1,l) = pplev(klon,l) |
|---|
| 103 | enddo |
|---|
| 104 | |
|---|
| 105 | do l=1,llm |
|---|
| 106 | ztemp(1,l) = ptemp(1,l) |
|---|
| 107 | zplay(1,l) = pplay(1,l) |
|---|
| 108 | do j=2,jjm |
|---|
| 109 | ig0=1+(j-2)*iim |
|---|
| 110 | do i=1,iim |
|---|
| 111 | ztemp(j,l) = ztemp(j,l) + ptemp(ig0+i,l)/iim |
|---|
| 112 | zplay(j,l) = zplay(j,l) + pplay(ig0+i,l)/iim |
|---|
| 113 | enddo |
|---|
| 114 | enddo |
|---|
| 115 | ztemp(jjm+1,l) = ptemp(klon,l) |
|---|
| 116 | zplay(jjm+1,l) = pplay(klon,l) |
|---|
| 117 | temp_eq = ztemp((jjm+1)/2,:) |
|---|
| 118 | press_eq = zplay((jjm+1)/2,:)/100. ! en mbar |
|---|
| 119 | enddo |
|---|
| 120 | |
|---|
| 121 | zmu0(1) = pmu0(1) |
|---|
| 122 | zfract(1) = pfract(1) |
|---|
| 123 | do j=2,jjm |
|---|
| 124 | ig0=1+(j-2)*iim |
|---|
| 125 | do i=1,iim |
|---|
| 126 | zmu0(j) = zmu0(j) + pmu0(ig0+i)/iim |
|---|
| 127 | zfract(j) = zfract(j) + pfract(ig0+i)/iim |
|---|
| 128 | enddo |
|---|
| 129 | enddo |
|---|
| 130 | zmu0(jjm+1) = pmu0(klon) |
|---|
| 131 | zfract(jjm+1) = pfract(klon) |
|---|
| 132 | |
|---|
| 133 | c TRACEURS MICROPHYSIQUES |
|---|
| 134 | |
|---|
| 135 | if (microfi.eq.1) then |
|---|
| 136 | |
|---|
| 137 | do iq=1,nmicro |
|---|
| 138 | c print*,tname(iq) |
|---|
| 139 | |
|---|
| 140 | c Traceurs microphysiques: passage en extensif: n/kg --> n/m^2 |
|---|
| 141 | DO l=1,llm |
|---|
| 142 | DO i = 1, klon |
|---|
| 143 | qaer(i,l,iq) = tr_seri(i,l,iq)*delp(i,l)/RG |
|---|
| 144 | ENDDO |
|---|
| 145 | ENDDO |
|---|
| 146 | |
|---|
| 147 | do l=1,llm |
|---|
| 148 | zqaer(1,l,iq) = qaer(1,l,iq) |
|---|
| 149 | do j=2,jjm |
|---|
| 150 | ig0=1+(j-2)*iim |
|---|
| 151 | do i=1,iim |
|---|
| 152 | zqaer(j,l,iq) = zqaer(j,l,iq) + qaer(ig0+i,l,iq)/iim |
|---|
| 153 | enddo |
|---|
| 154 | enddo |
|---|
| 155 | zqaer(jjm+1,l,iq) = qaer(klon,l,iq) |
|---|
| 156 | enddo |
|---|
| 157 | enddo |
|---|
| 158 | |
|---|
| 159 | endif ! microfi |
|---|
| 160 | |
|---|
| 161 | c AUTRES TRACEURS |
|---|
| 162 | |
|---|
| 163 | if (nqmax.gt.nmicro) then |
|---|
| 164 | do iq=nmicro+1,nqmax |
|---|
| 165 | do l=1,llm |
|---|
| 166 | ychim(1,l,iq-nmicro) = tr_seri(1,l,iq) |
|---|
| 167 | do j=2,jjm |
|---|
| 168 | ig0=1+(j-2)*iim |
|---|
| 169 | do i=1,iim |
|---|
| 170 | ychim(j,l,iq-nmicro) = ychim(j,l,iq-nmicro) |
|---|
| 171 | . + tr_seri(ig0+i,l,iq)/iim |
|---|
| 172 | enddo |
|---|
| 173 | enddo |
|---|
| 174 | ychim(jjm+1,l,iq-nmicro) = tr_seri(klon,l,iq) |
|---|
| 175 | enddo |
|---|
| 176 | nomqy(iq-nmicro) = tname(iq) |
|---|
| 177 | |
|---|
| 178 | c print*,iq-nmicro,nomqy(iq-nmicro) |
|---|
| 179 | |
|---|
| 180 | enddo |
|---|
| 181 | nomqy(nqmax-nmicro+1) = "HV" |
|---|
| 182 | endif |
|---|
| 183 | |
|---|
| 184 | c----------------------------------------------------------------------- |
|---|
| 185 | c initialisation des qysat au premier appel: |
|---|
| 186 | c --------------------------------- |
|---|
| 187 | |
|---|
| 188 | c!! ATTENTION, qysat pris uniquement a l'equateur |
|---|
| 189 | c!! justifie puisque dans cette region, les var de t et p sont faibles... |
|---|
| 190 | |
|---|
| 191 | if (firstcall .and. chimi .and.(nqmax.gt.nmicro)) then |
|---|
| 192 | call inicondens(nqmax-nmicro,press_eq,temp_eq,nomqy,qysat) |
|---|
| 193 | endif |
|---|
| 194 | |
|---|
| 195 | c----------------------------------------------------------------------- |
|---|
| 196 | c Appel de la microphysique |
|---|
| 197 | c -------------------------- |
|---|
| 198 | |
|---|
| 199 | pdqmfi = 0. |
|---|
| 200 | zqaer0 = zqaer |
|---|
| 201 | |
|---|
| 202 | IF(firstcall) THEN |
|---|
| 203 | print*,'MICROPHYSIQUE ',MICROFI |
|---|
| 204 | ENDIF |
|---|
| 205 | |
|---|
| 206 | IF(MICROFI.eq.1) THEN |
|---|
| 207 | |
|---|
| 208 | IF(firstcall) THEN |
|---|
| 209 | print*,'OH! On passe dans la microphysique' |
|---|
| 210 | ENDIF |
|---|
| 211 | |
|---|
| 212 | CALL pg2(zplev,ztemp,zqaer,pdqmfi ! tendances 2D, /s |
|---|
| 213 | & ,nmicro,ptimestep,zmu0,zfract,lastcall) |
|---|
| 214 | |
|---|
| 215 | ELSE |
|---|
| 216 | |
|---|
| 217 | IF(firstcall) THEN |
|---|
| 218 | print*,'MICROPHYSIQUE OFF-LINE',MICROFI |
|---|
| 219 | if (nmicro.gt.0) then |
|---|
| 220 | CALL pg2(zplev,ztemp,zqaer,pdqmfi |
|---|
| 221 | & ,nmicro,ptimestep,zmu0,zfract,lastcall) |
|---|
| 222 | endif |
|---|
| 223 | ENDIF |
|---|
| 224 | |
|---|
| 225 | ENDIF |
|---|
| 226 | |
|---|
| 227 | zqaer = zqaer+pdqmfi*ptimestep |
|---|
| 228 | pdqmfi = (zqaer-zqaer0)/ptimestep |
|---|
| 229 | |
|---|
| 230 | c----------------------------------------------------------------------- |
|---|
| 231 | c Condensation |
|---|
| 232 | c ------------- |
|---|
| 233 | |
|---|
| 234 | if ((chimi).and.(nqmax.gt.nmicro)) then |
|---|
| 235 | |
|---|
| 236 | c tendance (en /s) passee sur pdqmfi(nmicro+1 a nqmax) |
|---|
| 237 | c print*,'Condensation' |
|---|
| 238 | |
|---|
| 239 | do iq=1,nqmax-nmicro |
|---|
| 240 | do l=1,llm |
|---|
| 241 | do j=1,jjm+1 |
|---|
| 242 | if (ychim(j,l,iq).gt.qysat(l,iq)) then |
|---|
| 243 | pdqmfi(j,l,nmicro+iq)= (-ychim(j,l,iq)+qysat(l,iq)) !delta y |
|---|
| 244 | . / ptimestep ! / dt |
|---|
| 245 | endif |
|---|
| 246 | enddo |
|---|
| 247 | enddo |
|---|
| 248 | enddo |
|---|
| 249 | |
|---|
| 250 | ENDIF |
|---|
| 251 | |
|---|
| 252 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 253 | c eventuellement, modif initiale de la compo |
|---|
| 254 | c |
|---|
| 255 | c tendance (en /s) passee sur pdqmfi(nmicro+1 a nqmax) |
|---|
| 256 | c |
|---|
| 257 | c if (firstcall .and. chimi .and.(nqmax.gt.nmicro)) then |
|---|
| 258 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 259 | c!!!remise de CH4 a 1.5%!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 260 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 261 | c do iq=1,nqmax-nmicro |
|---|
| 262 | c if (nomqy(iq).eq."CH4") then |
|---|
| 263 | c do l=1,llm |
|---|
| 264 | c do j=1,jjm+1 |
|---|
| 265 | c if (ychim(j,l,iq).le.0.015) then |
|---|
| 266 | c pdqmfi(j,l,nmicro+iq)= (-ychim(j,l,iq)+0.015) !delta y |
|---|
| 267 | c . / ptimestep ! / dt |
|---|
| 268 | c endif |
|---|
| 269 | c enddo |
|---|
| 270 | c enddo |
|---|
| 271 | c endif |
|---|
| 272 | c enddo |
|---|
| 273 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 274 | c |
|---|
| 275 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 276 | c!!!remise de C2H2 a 1.e-5 max !!!!!!!!!!!!! |
|---|
| 277 | c!!!remise de C2H6 a 3.e-5 max !!!!!!!!!!!!! |
|---|
| 278 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 279 | c do iq=1,nqmax-nmicro |
|---|
| 280 | c if (nomqy(iq).eq."C2H2") then |
|---|
| 281 | c do l=1,llm |
|---|
| 282 | c do j=1,jjm+1 |
|---|
| 283 | c if (ychim(j,l,iq).gt.1.e-5) then |
|---|
| 284 | c pdqmfi(j,l,nmicro+iq)= (-ychim(j,l,iq)+1.e-5) !delta y |
|---|
| 285 | c . / ptimestep ! / dt |
|---|
| 286 | c endif |
|---|
| 287 | c enddo |
|---|
| 288 | c enddo |
|---|
| 289 | c endif |
|---|
| 290 | c if (nomqy(iq).eq."C2H6") then |
|---|
| 291 | c do l=1,llm |
|---|
| 292 | c do j=1,jjm+1 |
|---|
| 293 | c if (ychim(j,l,iq).gt.3.e-5) then |
|---|
| 294 | c pdqmfi(j,l,nmicro+iq)= (-ychim(j,l,iq)+3.e-5) !delta y |
|---|
| 295 | c . / ptimestep ! / dt |
|---|
| 296 | c endif |
|---|
| 297 | c enddo |
|---|
| 298 | c enddo |
|---|
| 299 | c endif |
|---|
| 300 | c enddo |
|---|
| 301 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 302 | c endif |
|---|
| 303 | c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 304 | |
|---|
| 305 | c----------------------------------------------------------------------- |
|---|
| 306 | c Appel de la chimie |
|---|
| 307 | c -------------------------- |
|---|
| 308 | |
|---|
| 309 | if((appkim.eq.1).and.(chimi)) then |
|---|
| 310 | print*,'On passe dans la CHIMIE' |
|---|
| 311 | |
|---|
| 312 | c do iq=1,nqmax-nmicro |
|---|
| 313 | c if (nomqy(iq).eq."C2H2") then |
|---|
| 314 | c print*,"C2H2top=",ychim(:,klev,iq) |
|---|
| 315 | c endif |
|---|
| 316 | c enddo |
|---|
| 317 | |
|---|
| 318 | c Appel Chimie |
|---|
| 319 | c ------------ |
|---|
| 320 | CALL calchim(nqmax-nmicro,ychim,nomqy,pdecli,lonsol,dtkim, |
|---|
| 321 | . ztemp,zplay,zplev, |
|---|
| 322 | . pdyfi) |
|---|
| 323 | c ychim ne doit pas etre modifie, pdyfi en /s |
|---|
| 324 | |
|---|
| 325 | endif |
|---|
| 326 | |
|---|
| 327 | c----------------------------------------------------------------------- |
|---|
| 328 | c retour des tendances vers 3D |
|---|
| 329 | c --------------------------------- |
|---|
| 330 | |
|---|
| 331 | c TRACEURS MICROPHYSIQUES |
|---|
| 332 | |
|---|
| 333 | if (microfi.eq.1) then |
|---|
| 334 | |
|---|
| 335 | DO iq=1,nmicro |
|---|
| 336 | DO l=1,llm |
|---|
| 337 | d_tr_mph(1,l,iq) = pdqmfi(1,l,iq) |
|---|
| 338 | ig0 = 2 |
|---|
| 339 | DO j=2,jjm |
|---|
| 340 | DO i = 1, iim |
|---|
| 341 | d_tr_mph(ig0,l,iq) = pdqmfi(j,l,iq) |
|---|
| 342 | & *qaer(ig0,l,iq)/zqaer0(j,l,iq) |
|---|
| 343 | ig0 = ig0 + 1 |
|---|
| 344 | ENDDO |
|---|
| 345 | ENDDO |
|---|
| 346 | d_tr_mph(ig0,l,iq) = pdqmfi(jjm+1,l,iq) |
|---|
| 347 | ENDDO |
|---|
| 348 | ENDDO |
|---|
| 349 | |
|---|
| 350 | do iq=1,nmicro |
|---|
| 351 | DO l=1,llm |
|---|
| 352 | DO i = 1, klon |
|---|
| 353 | c incrementation de la tendance sur qaer (pour sorties dans physiq.F) |
|---|
| 354 | qaer(i,l,iq) = qaer(i,l,iq) + d_tr_mph(i,l,iq)*ptimestep |
|---|
| 355 | c Traceurs microphysiques: passage en intensif: n/m^2 --> n/kg |
|---|
| 356 | d_tr_mph(i,l,iq) = d_tr_mph(i,l,iq)*RG/delp(i,l) |
|---|
| 357 | ENDDO |
|---|
| 358 | ENDDO |
|---|
| 359 | enddo |
|---|
| 360 | |
|---|
| 361 | endif ! microfi |
|---|
| 362 | |
|---|
| 363 | c AUTRES TRACEURS |
|---|
| 364 | |
|---|
| 365 | if ((chimi).and.(nqmax.gt.nmicro)) then |
|---|
| 366 | c on passe de pdyfi (tendance chimique en /s calculee quand chimie appelee) |
|---|
| 367 | c a d_tr_kim (tendance chimique 3D en /s, passee a physiq) |
|---|
| 368 | c et de pdqmfi a d_tr_mph (tendance condensation 3D en /s passee a physiq) |
|---|
| 369 | |
|---|
| 370 | DO iq=nmicro+1,nqmax |
|---|
| 371 | DO l=1,llm |
|---|
| 372 | d_tr_kim(1,l,iq) = pdyfi(1,l,iq-nmicro) |
|---|
| 373 | d_tr_mph(1,l,iq) = pdqmfi(1,l,iq) |
|---|
| 374 | ig0 = 2 |
|---|
| 375 | DO j=2,jjm |
|---|
| 376 | DO i = 1, iim |
|---|
| 377 | d_tr_kim(ig0,l,iq) = pdyfi(j,l,iq-nmicro) |
|---|
| 378 | & *tr_seri(ig0,l,iq)/ychim(j,l,iq-nmicro) |
|---|
| 379 | d_tr_mph(ig0,l,iq) = pdqmfi(j,l,iq) |
|---|
| 380 | & *tr_seri(ig0,l,iq)/ychim(j,l,iq-nmicro) |
|---|
| 381 | ig0 = ig0 + 1 |
|---|
| 382 | ENDDO |
|---|
| 383 | ENDDO |
|---|
| 384 | d_tr_kim(ig0,l,iq) = pdyfi(jjm+1,l,iq-nmicro) |
|---|
| 385 | d_tr_mph(ig0,l,iq) = pdqmfi(jjm+1,l,iq) |
|---|
| 386 | ENDDO |
|---|
| 387 | ENDDO |
|---|
| 388 | |
|---|
| 389 | endif ! chimi |
|---|
| 390 | |
|---|
| 391 | RETURN |
|---|
| 392 | END |
|---|
