Changeset 4109 for palm/trunk/SOURCE/chemistry_model_mod.f90
- Timestamp:
- Jul 22, 2019 5:00:34 PM (5 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/chemistry_model_mod.f90 (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/chemistry_model_mod.f90
r4102 r4109 22 22 ! Current revisions: 23 23 ! ----------------- 24 ! 24 ! - Decycling boundary conditions are only set at the ghost points not on the 25 ! prognostic grid points 26 ! - Allocation and initialization of special advection flags cs_advc_flags_s 27 ! used for chemistry. These are exclusively used for chemical species in 28 ! order to distinguish from the usually-used flags which might be different 29 ! when decycling is applied in combination with cyclic boundary conditions. 30 ! Moreover, cs_advc_flags_s considers extended zones around buildings where 31 ! first-order upwind scheme is applied for the horizontal advection terms, 32 ! in order to overcome high concentration peaks due to stationary numerical 33 ! oscillations caused by horizontal advection discretization. 25 34 ! 26 35 ! Former revisions: … … 344 353 345 354 USE advec_ws, & 346 ONLY: advec_s_ws 355 ONLY: advec_s_ws, ws_init_flags_scalar 347 356 348 357 USE diffusion_s_mod, & … … 353 362 354 363 USE indices, & 355 ONLY: nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, nz, nzb, nzt, wall_flags_0 364 ONLY: advc_flags_s, & 365 nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, nz, nzb, nzt, wall_flags_0 356 366 357 367 USE pegrid, & … … 363 373 USE control_parameters, & 364 374 ONLY: bc_lr_cyc, bc_ns_cyc, & 375 bc_dirichlet_l, & 376 bc_dirichlet_n, & 377 bc_dirichlet_r, & 378 bc_dirichlet_s, & 379 bc_radiation_l, & 380 bc_radiation_n, & 381 bc_radiation_r, & 382 bc_radiation_s, & 365 383 debug_output, & 366 384 dt_3d, humidity, initializing_actions, message_string, & … … 368 386 max_pr_user, & 369 387 monotonic_limiter_z, & 388 scalar_advec, & 370 389 timestep_scheme, use_prescribed_profile_data, ws_scheme_sca, air_chemistry 371 390 … … 962 981 IF ( boundary == 1 .AND. nxl == 0 ) THEN 963 982 ss = nxlg 964 ee = nxl +2983 ee = nxl-1 965 984 ELSEIF ( boundary == 2 .AND. nxr == nx ) THEN 966 ss = nxr -2985 ss = nxr+1 967 986 ee = nxrg 968 987 ENDIF … … 1026 1045 IF ( boundary == 3 .AND. nys == 0 ) THEN 1027 1046 ss = nysg 1028 ee = nys +21047 ee = nys-1 1029 1048 ELSEIF ( boundary == 4 .AND. nyn == ny ) THEN 1030 ss = nyn -21049 ss = nyn+1 1031 1050 ee = nyng 1032 1051 ENDIF … … 1854 1873 !------------------------------------------------------------------------------! 1855 1874 SUBROUTINE chem_init_internal 1856 1875 1857 1876 USE pegrid 1858 1877 … … 1921 1940 1922 1941 ENDDO 1942 ! 1943 !-- Set control flags for decycling only at lateral boundary cores, within the 1944 !-- inner cores the decycle flag is set to .False.. Even though it does not 1945 !-- affect the setting of chemistry boundary conditions, this flag is used to 1946 !-- set advection control flags appropriately. 1947 decycle_chem_lr = MERGE( decycle_chem_lr, .FALSE., & 1948 nxl == 0 .OR. nxr == nx ) 1949 decycle_chem_ns = MERGE( decycle_chem_ns, .FALSE., & 1950 nys == 0 .OR. nyn == ny ) 1951 ! 1952 !-- For some passive scalars decycling may be enabled. This case, the lateral 1953 !-- boundary conditions are non-cyclic for these scalars (chemical species 1954 !-- and aerosols), while the other scalars may have 1955 !-- cyclic boundary conditions. However, large gradients near the boundaries 1956 !-- may produce stationary numerical oscillations near the lateral boundaries 1957 !-- when a higher-order scheme is applied near these boundaries. 1958 !-- To get rid-off this, set-up additional flags that control the order of the 1959 !-- scalar advection scheme near the lateral boundaries for passive scalars 1960 !-- with decycling. 1961 IF ( scalar_advec == 'ws-scheme' ) THEN 1962 ALLOCATE( cs_advc_flags_s(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 1963 ! 1964 !-- In case of decyling, set Neumann boundary conditions for wall_flags_0 1965 !-- bit 31 instead of cyclic boundary conditions. 1966 !-- Bit 31 is used to identify extended degradation zones (please see 1967 !-- following comment). 1968 !-- Note, since several also other modules like Salsa or other future 1969 !-- one may access this bit but may have other boundary conditions, the 1970 !-- original value of wall_flags_0 bit 31 must not be modified. Hence, 1971 !-- store the boundary conditions directly on cs_advc_flags_s. 1972 !-- cs_advc_flags_s will be later overwritten in ws_init_flags_scalar and 1973 !-- bit 31 won't be used to control the numerical order. 1974 !-- Initialize with flag 31 only. 1975 cs_advc_flags_s = 0 1976 cs_advc_flags_s = MERGE( IBSET( cs_advc_flags_s, 31 ), 0, & 1977 BTEST( wall_flags_0, 31 ) ) 1978 1979 IF ( decycle_chem_ns ) THEN 1980 IF ( nys == 0 ) THEN 1981 DO i = 1, nbgp 1982 cs_advc_flags_s(:,nys-i,:) = MERGE( & 1983 IBSET( cs_advc_flags_s(:,nys,:), 31 ), & 1984 IBCLR( cs_advc_flags_s(:,nys,:), 31 ), & 1985 BTEST( cs_advc_flags_s(:,nys,:), 31 ) & 1986 ) 1987 ENDDO 1988 ENDIF 1989 IF ( nyn == ny ) THEN 1990 DO i = 1, nbgp 1991 cs_advc_flags_s(:,nyn+i,:) = MERGE( & 1992 IBSET( cs_advc_flags_s(:,nyn,:), 31 ), & 1993 IBCLR( cs_advc_flags_s(:,nyn,:), 31 ), & 1994 BTEST( cs_advc_flags_s(:,nyn,:), 31 ) & 1995 ) 1996 ENDDO 1997 ENDIF 1998 ENDIF 1999 IF ( .NOT. decycle_chem_lr ) THEN 2000 IF ( nxl == 0 ) THEN 2001 DO i = 1, nbgp 2002 cs_advc_flags_s(:,:,nxl-i) = MERGE( & 2003 IBSET( cs_advc_flags_s(:,:,nxl), 31 ), & 2004 IBCLR( cs_advc_flags_s(:,:,nxl), 31 ), & 2005 BTEST( cs_advc_flags_s(:,:,nxl), 31 ) & 2006 ) 2007 ENDDO 2008 ENDIF 2009 IF ( nxr == nx ) THEN 2010 DO i = 1, nbgp 2011 cs_advc_flags_s(:,:,nxr+i) = MERGE( & 2012 IBSET( cs_advc_flags_s(:,:,nxr), 31 ), & 2013 IBCLR( cs_advc_flags_s(:,:,nxr), 31 ), & 2014 BTEST( cs_advc_flags_s(:,:,nxr), 31 ) & 2015 ) 2016 ENDDO 2017 ENDIF 2018 ENDIF 2019 ! 2020 !-- To initialize advection flags appropriately, pass the boundary flags. 2021 !-- The last argument indicates that a passive scalar is treated, where 2022 !-- the horizontal advection terms are degraded already 2 grid points before 2023 !-- the lateral boundary to avoid stationary oscillations at large-gradients. 2024 !-- Also, extended degradation zones are applied, where horizontal advection of 2025 !-- passive scalars is discretized by first-order scheme at all grid points 2026 !-- that in the vicinity of buildings (<= 3 grid points). Even though no 2027 !-- building is within the numerical stencil, first-order scheme is used. 2028 !-- At fourth and fifth grid point the order of the horizontal advection scheme 2029 !-- is successively upgraded. 2030 !-- These extended degradation zones are used to avoid stationary numerical 2031 !-- oscillations, which are responsible for high concentration maxima that may 2032 !-- appear under shear-free stable conditions. 2033 CALL ws_init_flags_scalar( & 2034 bc_dirichlet_l .OR. bc_radiation_l .OR. decycle_chem_lr, & 2035 bc_dirichlet_n .OR. bc_radiation_n .OR. decycle_chem_ns, & 2036 bc_dirichlet_r .OR. bc_radiation_r .OR. decycle_chem_lr, & 2037 bc_dirichlet_s .OR. bc_radiation_s .OR. decycle_chem_ns, & 2038 cs_advc_flags_s, .TRUE. ) 2039 ENDIF 1923 2040 ! 1924 2041 !-- Initial concentration of profiles is prescribed by parameters cs_profile … … 2650 2767 IF ( timestep_scheme(1:5) == 'runge' ) THEN 2651 2768 IF ( ws_scheme_sca ) THEN 2652 CALL advec_s_ws( chem_species(ilsp)%conc, 'kc' ) 2769 CALL advec_s_ws( cs_advc_flags_s, chem_species(ilsp)%conc, 'kc', & 2770 bc_dirichlet_l .OR. bc_radiation_l .OR. decycle_chem_lr, & 2771 bc_dirichlet_n .OR. bc_radiation_n .OR. decycle_chem_ns, & 2772 bc_dirichlet_r .OR. bc_radiation_r .OR. decycle_chem_lr, & 2773 bc_dirichlet_s .OR. bc_radiation_s .OR. decycle_chem_ns ) 2653 2774 ELSE 2654 2775 CALL advec_s_pw( chem_species(ilsp)%conc ) … … 2755 2876 IF ( timestep_scheme(1:5) == 'runge' ) THEN 2756 2877 IF ( ws_scheme_sca ) THEN 2757 CALL advec_s_ws( i, j, chem_species(ilsp)%conc, 'kc', chem_species(ilsp)%flux_s_cs, & 2758 chem_species(ilsp)%diss_s_cs, chem_species(ilsp)%flux_l_cs, & 2759 chem_species(ilsp)%diss_l_cs, i_omp_start, tn, monotonic_limiter_z ) 2878 CALL advec_s_ws( cs_advc_flags_s, & 2879 i, & 2880 j, & 2881 chem_species(ilsp)%conc, & 2882 'kc', & 2883 chem_species(ilsp)%flux_s_cs, & 2884 chem_species(ilsp)%diss_s_cs, & 2885 chem_species(ilsp)%flux_l_cs, & 2886 chem_species(ilsp)%diss_l_cs, & 2887 i_omp_start, & 2888 tn, & 2889 bc_dirichlet_l .OR. bc_radiation_l .OR. decycle_chem_lr, & 2890 bc_dirichlet_n .OR. bc_radiation_n .OR. decycle_chem_ns, & 2891 bc_dirichlet_r .OR. bc_radiation_r .OR. decycle_chem_lr, & 2892 bc_dirichlet_s .OR. bc_radiation_s .OR. decycle_chem_ns, & 2893 monotonic_limiter_z ) 2760 2894 ELSE 2761 2895 CALL advec_s_pw( i, j, chem_species(ilsp)%conc )
Note: See TracChangeset
for help on using the changeset viewer.
