Changes in src/conccalc.f90 [4c64400:9287c01] in flexpart.git
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src/conccalc.f90
r4c64400 r9287c01 21 21 22 22 subroutine conccalc(itime,weight) 23 ! i i23 ! i i 24 24 !***************************************************************************** 25 25 ! * … … 59 59 real :: xl,yl,wx,wy,w 60 60 real,parameter :: factor=.596831, hxmax=6.0, hymax=4.0, hzmax=150. 61 61 ! integer xscav_count 62 62 63 63 ! For forward simulations, make a loop over the number of species; … … 65 65 ! releasepoints 66 66 !*************************************************************************** 67 67 ! xscav_count=0 68 68 do i=1,numpart 69 69 if (itra1(i).ne.itime) goto 20 … … 76 76 33 continue 77 77 78 78 ! if (xscav_frac1(i,1).lt.0) xscav_count=xscav_count+1 79 79 80 ! For special runs, interpolate the air density to the particle position 80 81 !************************************************************************ … … 126 127 !***************************************************************************** 127 128 do ind=indz,indzp 128 rhoprof(ind-indz+1)=p1*rho(ix ,jy ,ind, memind(2)) &129 rhoprof(ind-indz+1)=p1*rho(ix ,jy ,ind,2) & 129 130 +p2*rho(ixp,jy ,ind,2) & 130 131 +p3*rho(ix ,jyp,ind,2) & … … 172 173 if (yl.lt.0.) jy=jy-1 173 174 174 ! if (i.eq.10000) write(*,*) itime,xtra1(i),ytra1(i),ztra1(i),xl,yl175 175 176 176 … … 181 181 !***************************************************************************** 182 182 183 if ( lnokernel.or.(itage.lt.10800).or.(xl.lt.0.5).or.(yl.lt.0.5).or. &183 if (((itage.lt.10800).or.(xl.lt.0.5).or.(yl.lt.0.5).or. & 184 184 (xl.gt.real(numxgrid-1)-0.5).or. & 185 (yl.gt.real(numygrid-1)-0.5)) then ! no kernel, direct attribution to grid cell185 (yl.gt.real(numygrid-1)-0.5)).or.(.not.usekernel)) then ! no kernel, direct attribution to grid cell 186 186 if ((ix.ge.0).and.(jy.ge.0).and.(ix.le.numxgrid-1).and. & 187 187 (jy.le.numygrid-1)) then 188 do ks=1,nspec 189 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 188 if (DRYBKDEP.or.WETBKDEP) then 189 do ks=1,nspec 190 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 191 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 192 xmass1(i,ks)/rhoi*weight*max(xscav_frac1(i,ks),0.0) 193 end do 194 else 195 do ks=1,nspec 196 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 190 197 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 191 198 xmass1(i,ks)/rhoi*weight 192 end do 199 end do 200 endif 193 201 endif 194 202 195 else ! attribution via uniform kernel 203 else ! attribution via uniform kernel 196 204 197 205 ddx=xl-real(ix) ! distance to left cell border … … 220 228 if ((jy.ge.0).and.(jy.le.numygrid-1)) then 221 229 w=wx*wy 222 do ks=1,nspec 223 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 230 if (DRYBKDEP.or.WETBKDEP) then 231 do ks=1,nspec 232 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 233 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 234 xmass1(i,ks)/rhoi*w*weight*max(xscav_frac1(i,ks),0.0) 235 end do 236 else 237 do ks=1,nspec 238 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 224 239 gridunc(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 225 240 xmass1(i,ks)/rhoi*weight*w 226 end do 241 end do 242 endif 227 243 endif 228 244 229 245 if ((jyp.ge.0).and.(jyp.le.numygrid-1)) then 230 246 w=wx*(1.-wy) 231 do ks=1,nspec 232 gridunc(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 247 if (DRYBKDEP.or.WETBKDEP) then 248 do ks=1,nspec 249 gridunc(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 250 gridunc(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 251 xmass1(i,ks)/rhoi*weight*w*max(xscav_frac1(i,ks),0.0) 252 end do 253 else 254 do ks=1,nspec 255 gridunc(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 233 256 gridunc(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 234 257 xmass1(i,ks)/rhoi*weight*w 235 end do 236 endif 237 endif 258 end do 259 endif 260 endif 261 endif !ix ge 0 238 262 239 263 … … 241 265 if ((jyp.ge.0).and.(jyp.le.numygrid-1)) then 242 266 w=(1.-wx)*(1.-wy) 243 do ks=1,nspec 244 gridunc(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 267 if (DRYBKDEP.or.WETBKDEP) then 268 do ks=1,nspec 269 gridunc(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 270 gridunc(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 271 xmass1(i,ks)/rhoi*w*weight*max(xscav_frac1(i,ks),0.0) 272 end do 273 else 274 do ks=1,nspec 275 gridunc(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 245 276 gridunc(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 246 277 xmass1(i,ks)/rhoi*weight*w 247 end do 278 end do 279 endif 248 280 endif 249 281 250 282 if ((jy.ge.0).and.(jy.le.numygrid-1)) then 251 283 w=(1.-wx)*wy 252 do ks=1,nspec 253 gridunc(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)= & 284 if (DRYBKDEP.or.WETBKDEP) then 285 do ks=1,nspec 286 gridunc(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)= & 287 gridunc(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 288 xmass1(i,ks)/rhoi*weight*w*max(xscav_frac1(i,ks),0.0) 289 end do 290 else 291 do ks=1,nspec 292 gridunc(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)= & 254 293 gridunc(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 255 294 xmass1(i,ks)/rhoi*weight*w 256 end do 257 endif 258 endif 259 endif 260 261 295 end do 296 endif 297 endif 298 endif !ixp ge 0 299 endif 262 300 263 301 !************************************ … … 282 320 if ((itage.lt.10800).or.(xl.lt.0.5).or.(yl.lt.0.5).or. & 283 321 (xl.gt.real(numxgridn-1)-0.5).or. & 284 (yl.gt.real(numygridn-1)-0.5) ) then ! no kernel, direct attribution to grid cell322 (yl.gt.real(numygridn-1)-0.5).or.(.not.usekernel)) then ! no kernel, direct attribution to grid cell 285 323 if ((ix.ge.0).and.(jy.ge.0).and.(ix.le.numxgridn-1).and. & 286 324 (jy.le.numygridn-1)) then 287 do ks=1,nspec 288 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 325 if (DRYBKDEP.or.WETBKDEP) then 326 do ks=1,nspec 327 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 328 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 329 xmass1(i,ks)/rhoi*weight*max(xscav_frac1(i,ks),0.0) 330 end do 331 else 332 do ks=1,nspec 333 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 289 334 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 290 335 xmass1(i,ks)/rhoi*weight 291 end do 336 end do 337 endif 292 338 endif 293 339 … … 319 365 if ((jy.ge.0).and.(jy.le.numygridn-1)) then 320 366 w=wx*wy 321 do ks=1,nspec 322 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 367 if (DRYBKDEP.or.WETBKDEP) then 368 do ks=1,nspec 369 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 370 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 371 xmass1(i,ks)/rhoi*weight*w*max(xscav_frac1(i,ks),0.0) 372 end do 373 else 374 do ks=1,nspec 375 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)= & 323 376 griduncn(ix,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 324 377 xmass1(i,ks)/rhoi*weight*w 325 end do 378 end do 379 endif 326 380 endif 327 381 328 382 if ((jyp.ge.0).and.(jyp.le.numygridn-1)) then 329 383 w=wx*(1.-wy) 330 do ks=1,nspec 331 griduncn(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 384 if (DRYBKDEP.or.WETBKDEP) then 385 do ks=1,nspec 386 griduncn(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 387 griduncn(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 388 xmass1(i,ks)/rhoi*weight*w*max(xscav_frac1(i,ks),0.0) 389 end do 390 else 391 do ks=1,nspec 392 griduncn(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 332 393 griduncn(ix,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 333 394 xmass1(i,ks)/rhoi*weight*w 334 end do 395 end do 396 endif 335 397 endif 336 398 endif … … 340 402 if ((jyp.ge.0).and.(jyp.le.numygridn-1)) then 341 403 w=(1.-wx)*(1.-wy) 342 do ks=1,nspec 343 griduncn(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 404 if (DRYBKDEP.or.WETBKDEP) then 405 do ks=1,nspec 406 griduncn(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 407 griduncn(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 408 xmass1(i,ks)/rhoi*weight*w*max(xscav_frac1(i,ks),0.0) 409 end do 410 else 411 do ks=1,nspec 412 griduncn(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)= & 344 413 griduncn(ixp,jyp,kz,ks,nrelpointer,nclass(i),nage)+ & 345 414 xmass1(i,ks)/rhoi*weight*w 346 end do 415 end do 416 endif 347 417 endif 348 418 349 419 if ((jy.ge.0).and.(jy.le.numygridn-1)) then 350 420 w=(1.-wx)*wy 351 do ks=1,nspec 352 griduncn(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)= & 421 if (DRYBKDEP.or.WETBKDEP) then 422 do ks=1,nspec 423 griduncn(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)= & 424 griduncn(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 425 xmass1(i,ks)/rhoi*weight*w*max(xscav_frac1(i,ks),0.0) 426 end do 427 else 428 do ks=1,nspec 429 griduncn(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)= & 353 430 griduncn(ixp,jy,kz,ks,nrelpointer,nclass(i),nage)+ & 354 431 xmass1(i,ks)/rhoi*weight*w 355 end do 432 end do 433 endif 356 434 endif 357 435 endif 358 436 endif 359 360 437 endif 361 438 endif 362 439 20 continue 363 440 end do 441 ! write(*,*) 'xscav count:',xscav_count 364 442 365 443 !***********************************************************************
Note: See TracChangeset
for help on using the changeset viewer.