Changeset 4671 for palm/trunk/SOURCE/chemistry_model_mod.f90
- Timestamp:
- Sep 9, 2020 8:27:58 PM (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/chemistry_model_mod.f90
r4637 r4671 26 26 ! ----------------- 27 27 ! $Id$ 28 ! Implementation of downward facing USM and LSM surfaces 29 ! 30 ! 4637 2020-08-07 07:49:33Z suehring 28 31 ! Avoid usage of omp_lib, instead declare omp_get_thread_num explicitly 29 32 ! … … 684 687 DO j = nys, nyn 685 688 match_def = surf_def_h(0)%start_index(j,i) <= surf_def_h(0)%end_index(j,i) 686 match_lsm = surf_lsm_h %start_index(j,i) <= surf_lsm_h%end_index(j,i)687 match_usm = surf_usm_h %start_index(j,i) <= surf_usm_h%end_index(j,i)689 match_lsm = surf_lsm_h(0)%start_index(j,i) <= surf_lsm_h(0)%end_index(j,i) 690 match_usm = surf_usm_h(0)%start_index(j,i) <= surf_usm_h(0)%end_index(j,i) 688 691 689 692 IF ( match_def ) THEN … … 692 695 surf_def_h(0)%cssws(lsp,m) 693 696 ELSEIF ( match_lsm .AND. .NOT. match_usm ) THEN 694 m = surf_lsm_h %end_index(j,i)697 m = surf_lsm_h(0)%end_index(j,i) 695 698 chem_species(lsp)%cssws_av(j,i) = chem_species(lsp)%cssws_av(j,i) + & 696 surf_lsm_h %cssws(lsp,m)699 surf_lsm_h(0)%cssws(lsp,m) 697 700 ELSEIF ( match_usm ) THEN 698 m = surf_usm_h %end_index(j,i)701 m = surf_usm_h(0)%end_index(j,i) 699 702 chem_species(lsp)%cssws_av(j,i) = chem_species(lsp)%cssws_av(j,i) + & 700 surf_usm_h %cssws(lsp,m)703 surf_usm_h(0)%cssws(lsp,m) 701 704 ENDIF 702 705 ENDDO … … 1468 1471 ! 1469 1472 !-- No average output for now. 1470 DO m = 1, surf_lsm_h%ns 1471 local_pf(surf_lsm_h%i(m),surf_lsm_h%j(m),nzb+1) = & 1472 local_pf(surf_lsm_h%i(m),surf_lsm_h%j(m),nzb+1) & 1473 + surf_lsm_h%cssws(lsp,m) 1473 !-- !!! IT NEEDS TO RETHINK - with fully 3D structure, only lower (upper) 1474 !-- !!! upward facing horizontal surfaces should be taken into account here 1475 DO m = 1, surf_lsm_h(0)%ns 1476 local_pf(surf_lsm_h(0)%i(m),surf_lsm_h(0)%j(m),nzb+1) = & 1477 local_pf(surf_lsm_h(0)%i(m),surf_lsm_h(0)%j(m),nzb+1) & 1478 + surf_lsm_h(0)%cssws(lsp,m) 1474 1479 ENDDO 1475 DO m = 1, surf_usm_h %ns1476 local_pf(surf_usm_h %i(m),surf_usm_h%j(m),nzb+1) = &1477 local_pf(surf_usm_h %i(m),surf_usm_h%j(m),nzb+1) &1478 + surf_usm_h %cssws(lsp,m)1480 DO m = 1, surf_usm_h(0)%ns 1481 local_pf(surf_usm_h(0)%i(m),surf_usm_h(0)%j(m),nzb+1) = & 1482 local_pf(surf_usm_h(0)%i(m),surf_usm_h(0)%j(m),nzb+1) & 1483 + surf_usm_h(0)%cssws(lsp,m) 1479 1484 ENDDO 1480 1485 grid = 'zu' … … 1577 1582 ! 1578 1583 !-- no average for now 1579 DO m = 1, surf_usm_h%ns 1580 local_pf(surf_usm_h%i(m),surf_usm_h%j(m),surf_usm_h%k(m)) = & 1581 local_pf(surf_usm_h%i(m),surf_usm_h%j(m),surf_usm_h%k(m)) & 1582 + surf_usm_h%cssws(lsp,m) 1583 ENDDO 1584 DO m = 1, surf_lsm_h%ns 1585 local_pf(surf_lsm_h%i(m),surf_lsm_h%j(m),surf_lsm_h%k(m)) = & 1586 local_pf(surf_lsm_h%i(m),surf_lsm_h%j(m),surf_lsm_h%k(m)) & 1587 + surf_lsm_h%cssws(lsp,m) 1584 DO l = 0, 1 1585 DO m = 1, surf_usm_h(l)%ns 1586 local_pf(surf_usm_h(l)%i(m),surf_usm_h(l)%j(m),surf_usm_h(l)%k(m)) = & 1587 local_pf(surf_usm_h(l)%i(m),surf_usm_h(l)%j(m),surf_usm_h(l)%k(m)) & 1588 + surf_usm_h(l)%cssws(lsp,m) 1589 ENDDO 1590 DO m = 1, surf_lsm_h(l)%ns 1591 local_pf(surf_lsm_h(l)%i(m),surf_lsm_h(l)%j(m),surf_lsm_h(l)%k(m)) = & 1592 local_pf(surf_lsm_h(l)%i(m),surf_lsm_h(l)%j(m),surf_lsm_h(l)%k(m)) & 1593 + surf_lsm_h(l)%cssws(lsp,m) 1594 ENDDO 1588 1595 ENDDO 1589 1596 DO l = 0, 3 … … 2994 3001 surf_def_h(1)%cssws(ilsp,:), & 2995 3002 surf_def_h(2)%cssws(ilsp,:), & 2996 surf_lsm_h%cssws(ilsp,:), & 2997 surf_usm_h%cssws(ilsp,:), & 3003 surf_lsm_h(0)%cssws(ilsp,:), & 3004 surf_lsm_h(1)%cssws(ilsp,:), & 3005 surf_usm_h(0)%cssws(ilsp,:), & 3006 surf_usm_h(1)%cssws(ilsp,:), & 2998 3007 surf_def_v(0)%cssws(ilsp,:), & 2999 3008 surf_def_v(1)%cssws(ilsp,:), & … … 3115 3124 surf_def_h(0)%cssws(ilsp,:), surf_def_h(1)%cssws(ilsp,:), & 3116 3125 surf_def_h(2)%cssws(ilsp,:), & 3117 surf_lsm_h%cssws(ilsp,:), surf_usm_h%cssws(ilsp,:), & 3126 surf_lsm_h(0)%cssws(ilsp,:), surf_lsm_h(1)%cssws(ilsp,:), & 3127 surf_usm_h(0)%cssws(ilsp,:), surf_usm_h(1)%cssws(ilsp,:), & 3118 3128 surf_def_v(0)%cssws(ilsp,:), surf_def_v(1)%cssws(ilsp,:), & 3119 3129 surf_def_v(2)%cssws(ilsp,:), surf_def_v(3)%cssws(ilsp,:), & … … 3273 3283 INTEGER(iwp) :: j !< running index on y-axis 3274 3284 INTEGER(iwp) :: k !< vertical index counter 3285 INTEGER(iwp) :: l !< direction index counter 3275 3286 INTEGER(iwp) :: lpr !< running index chem spcs 3276 3287 INTEGER(iwp) :: m !< running index for surface elements … … 3328 3339 ENDDO 3329 3340 ! 3330 !-- Add surface fluxes 3331 surf_s = surf_def_h(0)%start_index(j,i) 3332 surf_e = surf_def_h(0)%end_index(j,i) 3333 DO m = surf_s, surf_e 3334 k = surf_def_h(0)%k(m) + surf_def_h(0)%koff 3335 sums_tmp(k,tn) = sums_tmp(k,tn) + surf_def_h(0)%cssws(cs_pr_index_fl_sgs(lpr),m) 3341 !-- Add surface fluxes (?Is the order mandatory or could it be done in one cycle?) 3342 DO l = 0, 1 3343 surf_s = surf_def_h(0)%start_index(j,i) 3344 surf_e = surf_def_h(0)%end_index(j,i) 3345 DO m = surf_s, surf_e 3346 k = surf_def_h(0)%k(m) + surf_def_h(0)%koff 3347 sums_tmp(k,tn) = sums_tmp(k,tn) + surf_def_h(0)%cssws(cs_pr_index_fl_sgs(lpr),m) 3348 ENDDO 3336 3349 ENDDO 3337 surf_s = surf_def_h(1)%start_index(j,i) 3338 surf_e = surf_def_h(1)%end_index(j,i) 3339 DO m = surf_s, surf_e 3340 k = surf_def_h(1)%k(m) + surf_def_h(1)%koff 3341 sums_tmp(k,tn) = sums_tmp(k,tn) + surf_def_h(1)%cssws(cs_pr_index_fl_sgs(lpr),m) 3350 DO l = 0, 1 3351 surf_s = surf_lsm_h(0)%start_index(j,i) 3352 surf_e = surf_lsm_h(0)%end_index(j,i) 3353 DO m = surf_s, surf_e 3354 k = surf_lsm_h(0)%k(m) + surf_lsm_h(0)%koff 3355 sums_tmp(k,tn) = sums_tmp(k,tn) + surf_lsm_h(0)%cssws(cs_pr_index_fl_sgs(lpr),m) 3356 ENDDO 3342 3357 ENDDO 3343 surf_s = surf_lsm_h%start_index(j,i) 3344 surf_e = surf_lsm_h%end_index(j,i) 3345 DO m = surf_s, surf_e 3346 k = surf_lsm_h%k(m) + surf_lsm_h%koff 3347 sums_tmp(k,tn) = sums_tmp(k,tn) + surf_lsm_h%cssws(cs_pr_index_fl_sgs(lpr),m) 3348 ENDDO 3349 surf_s = surf_usm_h%start_index(j,i) 3350 surf_e = surf_usm_h%end_index(j,i) 3351 DO m = surf_s, surf_e 3352 k = surf_usm_h%k(m) + surf_usm_h%koff 3353 sums_tmp(k,tn) = sums_tmp(k,tn) + surf_usm_h%cssws(cs_pr_index_fl_sgs(lpr),m) 3358 DO l = 0, 1 3359 surf_s = surf_usm_h(0)%start_index(j,i) 3360 surf_e = surf_usm_h(0)%end_index(j,i) 3361 DO m = surf_s, surf_e 3362 k = surf_usm_h(0)%k(m) + surf_usm_h(0)%koff 3363 sums_tmp(k,tn) = sums_tmp(k,tn) + surf_usm_h(0)%cssws(cs_pr_index_fl_sgs(lpr),m) 3364 ENDDO 3354 3365 ENDDO 3355 3366 ENDDO … … 3450 3461 ! ------------ 3451 3462 !> Subroutine to calculate the deposition of gases and PMs. For now deposition only takes place on 3452 !> lsm and usm horizontal surfaces. Default surfaces are NOT considered. The deposition of particles3453 !> is derived following Zhang et al., 2001, gases are deposited using the DEPAC module3454 !> (van Zanten et al., 2010).3463 !> lsm and usm horizontal upward faceing surfaces. Default surfaces are NOT considered. 3464 !> The deposition of particlesis derived following Zhang et al., 2001, gases are deposited using 3465 !> the DEPAC module (van Zanten et al., 2010). 3455 3466 !> 3456 3467 !> @TODO: Consider deposition on vertical surfaces … … 3768 3779 k = 0 3769 3780 ! 3770 !-- LSM or USM surface present at i,j:3781 !-- LSM or USM horizintal upward facing surface present at i,j: 3771 3782 !-- Default surfaces are NOT considered for deposition 3772 match_lsm = surf_lsm_h %start_index(j,i) <= surf_lsm_h%end_index(j,i)3773 match_usm = surf_usm_h %start_index(j,i) <= surf_usm_h%end_index(j,i)3783 match_lsm = surf_lsm_h(0)%start_index(j,i) <= surf_lsm_h(0)%end_index(j,i) 3784 match_usm = surf_usm_h(0)%start_index(j,i) <= surf_usm_h(0)%end_index(j,i) 3774 3785 ! 3775 3786 !--For LSM surfaces … … 3777 3788 ! 3778 3789 !-- Get surface element information at i,j: 3779 m = surf_lsm_h %start_index(j,i)3780 k = surf_lsm_h %k(m)3790 m = surf_lsm_h(0)%start_index(j,i) 3791 k = surf_lsm_h(0)%k(m) 3781 3792 ! 3782 3793 !-- Get needed variables for surface element m 3783 ustar_surf = surf_lsm_h %us(m)3784 z0h_surf = surf_lsm_h %z0h(m)3785 r_aero_surf = surf_lsm_h %r_a(m)3786 solar_rad = surf_lsm_h %rad_sw_dir(m) + surf_lsm_h%rad_sw_dif(m)3787 3788 lai = surf_lsm_h %lai(m)3794 ustar_surf = surf_lsm_h(0)%us(m) 3795 z0h_surf = surf_lsm_h(0)%z0h(m) 3796 r_aero_surf = surf_lsm_h(0)%r_a(m) 3797 solar_rad = surf_lsm_h(0)%rad_sw_dir(m) + surf_lsm_h(0)%rad_sw_dif(m) 3798 3799 lai = surf_lsm_h(0)%lai(m) 3789 3800 sai = lai + 1 3790 3801 ! … … 3808 3819 ! 3809 3820 !-- Get land use for i,j and assign to DEPAC lu 3810 IF ( surf_lsm_h %frac(m,ind_veg_wall) > 0 ) THEN3811 luv_palm = surf_lsm_h %vegetation_type(m)3821 IF ( surf_lsm_h(0)%frac(m,ind_veg_wall) > 0 ) THEN 3822 luv_palm = surf_lsm_h(0)%vegetation_type(m) 3812 3823 IF ( luv_palm == ind_luv_user ) THEN 3813 3824 message_string = 'No vegetation type defined. Please define vegetation type to enable deposition calculation' … … 3852 3863 ENDIF 3853 3864 3854 IF ( surf_lsm_h %frac(m,ind_pav_green) > 0 ) THEN3855 lup_palm = surf_lsm_h %pavement_type(m)3865 IF ( surf_lsm_h(0)%frac(m,ind_pav_green) > 0 ) THEN 3866 lup_palm = surf_lsm_h(0)%pavement_type(m) 3856 3867 IF ( lup_palm == ind_lup_user ) THEN 3857 3868 message_string = 'No pavement type defined. Please define pavement type to enable deposition calculation' … … 3890 3901 ENDIF 3891 3902 3892 IF ( surf_lsm_h %frac(m,ind_wat_win) > 0 ) THEN3893 luw_palm = surf_lsm_h %water_type(m)3903 IF ( surf_lsm_h(0)%frac(m,ind_wat_win) > 0 ) THEN 3904 luw_palm = surf_lsm_h(0)%water_type(m) 3894 3905 IF ( luw_palm == ind_luw_user ) THEN 3895 3906 message_string = 'No water type defined. Please define water type to enable deposition calculation' … … 3909 3920 ! 3910 3921 !-- Set wetness indicator to dry or wet for lsm vegetation or pavement 3911 IF ( surf_lsm_h %c_liq(m) > 0 ) THEN3922 IF ( surf_lsm_h(0)%c_liq(m) > 0 ) THEN 3912 3923 nwet = 1 3913 3924 ELSE … … 3950 3961 ! 3951 3962 !-- Vegetation 3952 IF ( surf_lsm_h %frac(m,ind_veg_wall) > 0 ) THEN3963 IF ( surf_lsm_h(0)%frac(m,ind_veg_wall) > 0 ) THEN 3953 3964 3954 3965 ! … … 4084 4095 ! 4085 4096 !-- Pavement 4086 IF ( surf_lsm_h %frac(m,ind_pav_green) > 0 ) THEN4097 IF ( surf_lsm_h(0)%frac(m,ind_pav_green) > 0 ) THEN 4087 4098 ! 4088 4099 !-- No vegetation on pavements: … … 4208 4219 ! 4209 4220 !-- Water 4210 IF ( surf_lsm_h %frac(m,ind_wat_win) > 0 ) THEN4221 IF ( surf_lsm_h(0)%frac(m,ind_wat_win) > 0 ) THEN 4211 4222 ! 4212 4223 !-- No vegetation on water: … … 4339 4350 DO lsp = 1, nspec 4340 4351 4341 bud(lsp) = surf_lsm_h %frac(m,ind_veg_wall) * bud_luv(lsp) + &4342 surf_lsm_h %frac(m,ind_pav_green) * bud_lup(lsp) + &4343 surf_lsm_h %frac(m,ind_wat_win) * bud_luw(lsp)4352 bud(lsp) = surf_lsm_h(0)%frac(m,ind_veg_wall) * bud_luv(lsp) + & 4353 surf_lsm_h(0)%frac(m,ind_pav_green) * bud_lup(lsp) + & 4354 surf_lsm_h(0)%frac(m,ind_wat_win) * bud_luw(lsp) 4344 4355 ! 4345 4356 !-- Compute new concentration: … … 4356 4367 ! 4357 4368 !-- Get surface element information at i,j: 4358 m = surf_usm_h %start_index(j,i)4359 k = surf_usm_h %k(m)4369 m = surf_usm_h(0)%start_index(j,i) 4370 k = surf_usm_h(0)%k(m) 4360 4371 ! 4361 4372 !-- Get needed variables for surface element m 4362 ustar_surf = surf_usm_h %us(m)4363 z0h_surf = surf_usm_h %z0h(m)4364 r_aero_surf = surf_usm_h %r_a(m)4365 solar_rad = surf_usm_h %rad_sw_dir(m) + surf_usm_h%rad_sw_dif(m)4366 lai = surf_usm_h %lai(m)4373 ustar_surf = surf_usm_h(0)%us(m) 4374 z0h_surf = surf_usm_h(0)%z0h(m) 4375 r_aero_surf = surf_usm_h(0)%r_a(m) 4376 solar_rad = surf_usm_h(0)%rad_sw_dir(m) + surf_usm_h(0)%rad_sw_dif(m) 4377 lai = surf_usm_h(0)%lai(m) 4367 4378 sai = lai + 1 4368 4379 ! … … 4386 4397 ! 4387 4398 !-- Get land use for i,j and assign to DEPAC lu 4388 IF ( surf_usm_h %frac(m,ind_pav_green) > 0 ) THEN4399 IF ( surf_usm_h(0)%frac(m,ind_pav_green) > 0 ) THEN 4389 4400 ! 4390 4401 !-- For green urban surfaces (e.g. green roofs assume LU short grass … … 4432 4443 ENDIF 4433 4444 4434 IF ( surf_usm_h %frac(m,ind_veg_wall) > 0 ) THEN4445 IF ( surf_usm_h(0)%frac(m,ind_veg_wall) > 0 ) THEN 4435 4446 ! 4436 4447 !-- For walls in USM assume concrete walls/roofs, … … 4473 4484 ENDIF 4474 4485 4475 IF ( surf_usm_h %frac(m,ind_wat_win) > 0 ) THEN4486 IF ( surf_usm_h(0)%frac(m,ind_wat_win) > 0 ) THEN 4476 4487 ! 4477 4488 !-- For windows in USM assume metal as this is as close as we get, assumed LU class desert as … … 4516 4527 !-- @TODO: Activate these lines as soon as new ebsolver branch is merged: 4517 4528 !-- Set wetness indicator to dry or wet for usm vegetation or pavement 4518 !IF ( surf_usm_h %c_liq(m) > 0 ) THEN4529 !IF ( surf_usm_h(0)%c_liq(m) > 0 ) THEN 4519 4530 ! nwet = 1 4520 4531 !ELSE … … 4557 4568 4558 4569 !-- Walls/roofs 4559 IF ( surf_usm_h %frac(m,ind_veg_wall) > 0 ) THEN4570 IF ( surf_usm_h(0)%frac(m,ind_veg_wall) > 0 ) THEN 4560 4571 ! 4561 4572 !-- No vegetation on non-green walls: … … 4684 4695 ! 4685 4696 !-- Green usm surfaces 4686 IF ( surf_usm_h %frac(m,ind_pav_green) > 0 ) THEN4697 IF ( surf_usm_h(0)%frac(m,ind_pav_green) > 0 ) THEN 4687 4698 4688 4699 ! … … 4817 4828 ! 4818 4829 !-- Windows 4819 IF ( surf_usm_h %frac(m,ind_wat_win) > 0 ) THEN4830 IF ( surf_usm_h(0)%frac(m,ind_wat_win) > 0 ) THEN 4820 4831 ! 4821 4832 !-- No vegetation on windows: … … 4947 4958 4948 4959 4949 bud(lsp) = surf_usm_h %frac(m,ind_veg_wall) * bud_luu(lsp) + &4950 surf_usm_h %frac(m,ind_pav_green) * bud_lug(lsp) + &4951 surf_usm_h %frac(m,ind_wat_win) * bud_lud(lsp)4960 bud(lsp) = surf_usm_h(0)%frac(m,ind_veg_wall) * bud_luu(lsp) + & 4961 surf_usm_h(0)%frac(m,ind_pav_green) * bud_lug(lsp) + & 4962 surf_usm_h(0)%frac(m,ind_wat_win) * bud_lud(lsp) 4952 4963 ! 4953 4964 !-- Compute new concentration
Note: See TracChangeset
for help on using the changeset viewer.