Changeset 4511 for palm


Ignore:
Timestamp:
Apr 30, 2020 12:20:40 PM (7 months ago)
Author:
raasch
Message:

chemistry decycling replaced by explicit setting of lateral boundary conditions

Location:
palm/trunk/SOURCE
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/check_parameters.f90

    r4495 r4511  
    2525! -----------------
    2626! $Id$
     27! call of chem_boundary_conds removed (respective settings are now done in the chemistry module)
     28!
     29! 4495 2020-04-13 20:11:20Z raasch
    2730! check new restart_data_format parameters
    2831!
     
    131134
    132135    USE chem_modules
    133 
    134     USE chemistry_model_mod,                                                   &
    135         ONLY:  chem_boundary_conds
    136136
    137137    USE control_parameters
     
    15541554
    15551555!
    1556 !-- Boundary conditions for chemical species
    1557     IF ( air_chemistry )  CALL chem_boundary_conds( 'init' )
    1558 
    1559 !
    15601556!-- Boundary conditions for horizontal components of wind speed
    15611557    IF ( bc_uv_b == 'dirichlet' )  THEN
  • palm/trunk/SOURCE/chem_modules.f90

    r4481 r4511  
    2222! Current revisions:
    2323! -----------------
    24 !
     24! 
    2525!
    2626! Former revisions:
    2727! -----------------
    2828! $Id$
     29! new variables for explicit settings of lateral boundary conditions introduced
     30!
     31! 4481 2020-03-31 18:55:54Z maronga
    2932! added namelist flag 'emiss_read_legacy_mode' to allow concurrent
    3033! functioning of new emission read mode under development (ECC)
     
    9295
    9396    CHARACTER (LEN=20) ::  bc_cs_b        = 'dirichlet'         !< namelist parameter: surface boundary condition for concentration
     97    CHARACTER (LEN=20) ::  bc_cs_l        = 'undefined'         !< left boundary condition
     98    CHARACTER (LEN=20) ::  bc_cs_n        = 'undefined'         !< north boundary condition
     99    CHARACTER (LEN=20) ::  bc_cs_r        = 'undefined'         !< right boundary condition
     100    CHARACTER (LEN=20) ::  bc_cs_s        = 'undefined'         !< south boundary condition
    94101    CHARACTER (LEN=20) ::  bc_cs_t        = 'initial_gradient'  !< namelist parameter: top boudary condition for concentration
    95102    CHARACTER (LEN=30) ::  chem_mechanism = 'phstatp'           !< namelist parameter: chemistry mechanism
     
    105112    CHARACTER (LEN=11), DIMENSION(99) ::  data_output_pr_cs   = 'novalue'  !< namelist parameter: tbc...???
    106113    CHARACTER (LEN=11), DIMENSION(99) ::  surface_csflux_name = 'novalue'  !< namelist parameter: tbc...???
     114
     115    INTEGER(iwp) ::  communicator_chem      !< stores the number of the MPI communicator to be used
     116                                            !< for ghost layer data exchange
     117                                            !< 1: cyclic, 2: cyclic along x, 3: cyclic along y,
     118                                            !< 4: non-cyclic
    107119
    108120    INTEGER(iwp) ::  cs_pr_count                           = 0      !< counter for chemical species profiles
     
    131143                                                                     !< chemical species near walls and lateral boundaries
    132144
     145    LOGICAL ::  bc_dirichlet_cs_l         = .FALSE.  !< flag for indicating a dirichlet condition at
     146                                                     !< the left boundary
     147    LOGICAL ::  bc_dirichlet_cs_n         = .FALSE.  !< flag for indicating a dirichlet condition at
     148                                                     !< the north boundary
     149    LOGICAL ::  bc_dirichlet_cs_r         = .FALSE.  !< flag for indicating a dirichlet condition at
     150                                                     !< the right boundary
     151    LOGICAL ::  bc_dirichlet_cs_s         = .FALSE.  !< flag for indicating a dirichlet condition at
     152                                                     !< the south boundary
     153    LOGICAL ::  bc_radiation_cs_l         = .FALSE.  !< flag for indicating a radiation/neumann
     154                                                     !< condition at the left boundary
     155    LOGICAL ::  bc_radiation_cs_n         = .FALSE.  !< flag for indicating a radiation/neumann
     156                                                     !< condition at the north boundary
     157    LOGICAL ::  bc_radiation_cs_r         = .FALSE.  !< flag for indicating a radiation/neumann
     158                                                     !< condition at the right boundary
     159    LOGICAL ::  bc_radiation_cs_s         = .FALSE.  !< flag for indicating a radiation/neumann
     160                                                     !< condition at the south boundary
    133161    LOGICAL ::  constant_top_csflux(99)   = .TRUE.   !< internal flag, set to .FALSE. if no top_csflux is prescribed
    134162    LOGICAL ::  constant_csflux(99)       = .TRUE.   !< internal flag, set to .FALSE. if no surface_csflux is prescribed
  • palm/trunk/SOURCE/chemistry_model_mod.f90

    r4487 r4511  
    2727! -----------------
    2828! $Id$
     29! decycling replaced by explicit setting of lateral boundary conditions
     30!
     31! 4487 2020-04-03 09:38:20Z raasch
    2932! bugfix for subroutine calls that contain the decycle_chem switches as arguments
    30 ! 
     33!
    3134! 4481 2020-03-31 18:55:54Z maronga
    3235! use statement for exchange horiz added,
     
    347350    REAL(kind=wp), DIMENSION(nkppctrl)                     ::  rcntrl       !< 20 real parameters for fine tuning of KPP code
    348351                                                                            !< (e.g starting internal timestep of solver)
    349 !
    350 !-- Decycling of chemistry variables: Dirichlet BCs with cyclic is frequently not
    351 !-- approproate for chemicals compounds since they may accumulate too much.
    352 !-- If no proper boundary conditions from a DYNAMIC input file are available,
    353 !-- de-cycling applies the initial profiles at the inflow boundaries for
    354 !-- Dirichlet boundary conditions
    355     LOGICAL ::  decycle_chem_lr    = .FALSE.    !< switch for setting decycling in left-right direction
    356     LOGICAL ::  decycle_chem_ns    = .FALSE.    !< switch for setting decycling in south-norht direction
    357     CHARACTER (LEN=20), DIMENSION(4) ::  decycle_method = &
    358          (/'dirichlet','dirichlet','dirichlet','dirichlet'/)
    359                               !< Decycling method at horizontal boundaries,
    360                               !< 1=left, 2=right, 3=south, 4=north
    361                               !< dirichlet = initial size distribution and
    362                               !< chemical composition set for the ghost and
    363                               !< first three layers
    364                               !< neumann = zero gradient
    365352
    366353    REAL(kind=wp), PUBLIC ::  cs_time_step = 0._wp
     
    438425    END INTERFACE chem_3d_data_averaging
    439426
    440     INTERFACE chem_boundary_conds
    441        MODULE PROCEDURE chem_boundary_conds
    442        MODULE PROCEDURE chem_boundary_conds_decycle
    443     END INTERFACE chem_boundary_conds
    444 
    445427    INTERFACE chem_boundary_conditions
    446428       MODULE PROCEDURE chem_boundary_conditions
     
    590572    END INTERFACE rc_rctot
    591573
    592 !    INTERFACE rc_comp_point_rc_eff
    593 !       MODULE PROCEDURE rc_comp_point_rc_eff
    594 !    END INTERFACE rc_comp_point_rc_eff
    595 
    596574    INTERFACE drydepo_aero_zhang_vd
    597575       MODULE PROCEDURE drydepo_aero_zhang_vd
     
    604582
    605583
    606     PUBLIC chem_3d_data_averaging, chem_boundary_conds, chem_boundary_conditions, &
    607             chem_boundary_conds_decycle, chem_check_data_output,              &
    608          chem_check_data_output_pr, chem_check_parameters,                    &
    609          chem_data_output_2d, chem_data_output_3d, chem_data_output_mask,     &
    610          chem_define_netcdf_grid, chem_header, chem_init, chem_init_arrays,   &
    611          chem_init_profiles, chem_integrate, chem_parin,                      &
    612          chem_actions, chem_prognostic_equations, chem_rrd_local,             &
    613          chem_statistics, chem_swap_timelevel, chem_wrd_local, chem_depo,     &
    614          chem_non_advective_processes, chem_exchange_horiz_bounds
     584    PUBLIC chem_3d_data_averaging, chem_boundary_conditions, chem_check_data_output,               &
     585           chem_check_data_output_pr, chem_check_parameters, chem_data_output_2d,                  &
     586           chem_data_output_3d, chem_data_output_mask, chem_define_netcdf_grid, chem_header,       &
     587           chem_init, chem_init_arrays, chem_init_profiles, chem_integrate, chem_parin,            &
     588           chem_actions, chem_prognostic_equations, chem_rrd_local, chem_statistics,               &
     589           chem_swap_timelevel, chem_wrd_local, chem_depo, chem_non_advective_processes,           &
     590           chem_exchange_horiz_bounds
    615591
    616592 CONTAINS
     
    742718!> Subroutine to initialize and set all boundary conditions for chemical species
    743719!------------------------------------------------------------------------------!
    744  SUBROUTINE chem_boundary_conds( mode )                                           
    745 
    746     USE control_parameters,                                                    & 
    747         ONLY:  bc_radiation_l, bc_radiation_n, bc_radiation_r, bc_radiation_s
    748 
    749     USE arrays_3d,                                                             &     
    750         ONLY:  dzu                                               
     720 SUBROUTINE chem_boundary_conditions( horizontal_conditions_only )
     721
     722    USE arrays_3d,                                                             &
     723        ONLY:  dzu
    751724
    752725    USE surface_mod,                                                           &
    753726        ONLY:  bc_h                                                           
    754727
    755     CHARACTER (LEN=*), INTENT(IN) ::  mode
    756728    INTEGER(iwp) ::  i                            !< grid index x direction.
    757729    INTEGER(iwp) ::  j                            !< grid index y direction.
     
    761733    INTEGER(iwp) ::  lsp                          !< running index for chem spcs.
    762734
    763 
    764     SELECT CASE ( TRIM( mode ) )       
    765        CASE ( 'init' )
    766 
    767           IF ( bc_cs_b == 'dirichlet' )  THEN
    768              ibc_cs_b = 0
    769           ELSEIF ( bc_cs_b == 'neumann' )  THEN
    770              ibc_cs_b = 1
    771           ELSE
    772              message_string = 'unknown boundary condition: bc_cs_b ="' // TRIM( bc_cs_b ) // '"' 
    773              CALL message( 'chem_boundary_conds', 'CM0429', 1, 2, 0, 6, 0 )
    774           ENDIF                                                                 
    775 !
    776 !--       Set Integer flags and check for possible erroneous settings for top
    777 !--       boundary condition.
    778           IF ( bc_cs_t == 'dirichlet' )  THEN
    779              ibc_cs_t = 0
    780           ELSEIF ( bc_cs_t == 'neumann' )  THEN
    781              ibc_cs_t = 1
    782           ELSEIF ( bc_cs_t == 'initial_gradient' )  THEN
    783              ibc_cs_t = 2
    784           ELSEIF ( bc_cs_t == 'nested' )  THEN         
    785              ibc_cs_t = 3
    786           ELSE
    787              message_string = 'unknown boundary condition: bc_c_t ="' // TRIM( bc_cs_t ) // '"'     
    788              CALL message( 'check_parameters', 'CM0430', 1, 2, 0, 6, 0 )
    789           ENDIF
    790 
    791        CASE ( 'set_bc_bottomtop' )                   
    792 !
    793 !--       Boundary condtions for chemical species at horizontal walls     
    794           DO  lsp = 1, nspec                                                     
    795              IF ( ibc_cs_b == 0 )  THEN
    796                 DO  l = 0, 1
    797                    !$OMP PARALLEL DO PRIVATE( i, j, k )
    798                    DO  m = 1, bc_h(l)%ns
    799                        i = bc_h(l)%i(m)           
    800                        j = bc_h(l)%j(m)
    801                        k = bc_h(l)%k(m)
    802                       chem_species(lsp)%conc_p(k+bc_h(l)%koff,j,i) =           &
     735    LOGICAL, OPTIONAL ::  horizontal_conditions_only  !< switch to set horizontal bc only
     736
     737
     738    IF ( .NOT. PRESENT( horizontal_conditions_only ) )  THEN
     739!
     740!--    Boundary condtions for chemical species at horizontal walls
     741       DO  lsp = 1, nspec
     742
     743          IF ( ibc_cs_b == 0 )  THEN
     744             DO  l = 0, 1
     745                !$OMP PARALLEL DO PRIVATE( i, j, k )
     746                DO  m = 1, bc_h(l)%ns
     747                    i = bc_h(l)%i(m)
     748                    j = bc_h(l)%j(m)
     749                    k = bc_h(l)%k(m)
     750                   chem_species(lsp)%conc_p(k+bc_h(l)%koff,j,i) =           &
    803751                                      chem_species(lsp)%conc(k+bc_h(l)%koff,j,i)
    804                    ENDDO                                       
    805                 ENDDO                                       
    806 
    807              ELSEIF ( ibc_cs_b == 1 )  THEN
    808 !
    809 !--             in boundary_conds there is som extra loop over m here for passive tracer
    810                 DO  l = 0, 1
    811                    !$OMP PARALLEL DO PRIVATE( i, j, k )                                           
    812                    DO m = 1, bc_h(l)%ns
    813                       i = bc_h(l)%i(m)           
    814                       j = bc_h(l)%j(m)
    815                       k = bc_h(l)%k(m)
    816                       chem_species(lsp)%conc_p(k+bc_h(l)%koff,j,i) =           &
     752                ENDDO
     753             ENDDO
     754
     755          ELSEIF ( ibc_cs_b == 1 )  THEN
     756!
     757!--          In boundary_conds there is som extra loop over m here for passive tracer
     758!>           TODO: clarify the meaning of the above comment. Explain in more detail or remove it. (Siggi)
     759             DO  l = 0, 1
     760                !$OMP PARALLEL DO PRIVATE( i, j, k )
     761                DO m = 1, bc_h(l)%ns
     762                   i = bc_h(l)%i(m)
     763                   j = bc_h(l)%j(m)
     764                   k = bc_h(l)%k(m)
     765                   chem_species(lsp)%conc_p(k+bc_h(l)%koff,j,i) =           &
    817766                                         chem_species(lsp)%conc_p(k,j,i)
    818 
    819                    ENDDO
    820767                ENDDO
    821              ENDIF
    822        ENDDO    ! end lsp loop 
    823 !
    824 !--    Top boundary conditions for chemical species - Should this not be done for all species?
    825           IF ( ibc_cs_t == 0 )  THEN                   
    826              DO  lsp = 1, nspec
    827                 chem_species(lsp)%conc_p(nzt+1,:,:) = chem_species(lsp)%conc(nzt+1,:,:)       
    828              ENDDO
    829           ELSEIF ( ibc_cs_t == 1 )  THEN
    830              DO  lsp = 1, nspec
    831                 chem_species(lsp)%conc_p(nzt+1,:,:) = chem_species(lsp)%conc_p(nzt,:,:)
    832              ENDDO
    833           ELSEIF ( ibc_cs_t == 2 )  THEN
    834              DO  lsp = 1, nspec
    835                 chem_species(lsp)%conc_p(nzt+1,:,:) = chem_species(lsp)%conc_p(nzt,:,:) + bc_cs_t_val(lsp) * dzu(nzt+1)
    836768             ENDDO
    837769          ENDIF
    838770
    839        CASE ( 'set_bc_lateral' )                       
    840 !
    841 !--             Lateral boundary conditions for chem species at inflow boundary
    842 !--             are automatically set when chem_species concentration is
    843 !--             initialized. The initially set value at the inflow boundary is not
    844 !--             touched during time integration, hence, this boundary value remains
    845 !--             at a constant value, which is the concentration that flows into the
    846 !--             domain.                                                           
    847 !--             Lateral boundary conditions for chem species at outflow boundary
    848 
    849           IF ( bc_radiation_s )  THEN
    850              DO  lsp = 1, nspec
    851                 chem_species(lsp)%conc_p(:,nys-1,:) = chem_species(lsp)%conc_p(:,nys,:)
    852              ENDDO
    853           ELSEIF ( bc_radiation_n )  THEN
    854              DO  lsp = 1, nspec
    855                 chem_species(lsp)%conc_p(:,nyn+1,:) = chem_species(lsp)%conc_p(:,nyn,:)
    856              ENDDO
    857           ELSEIF ( bc_radiation_l )  THEN
    858              DO  lsp = 1, nspec
    859                 chem_species(lsp)%conc_p(:,:,nxl-1) = chem_species(lsp)%conc_p(:,:,nxl)
    860              ENDDO
    861           ELSEIF ( bc_radiation_r )  THEN
    862              DO  lsp = 1, nspec
    863                 chem_species(lsp)%conc_p(:,:,nxr+1) = chem_species(lsp)%conc_p(:,:,nxr)
    864              ENDDO
    865           ENDIF
    866 
    867     END SELECT
    868 
    869  END SUBROUTINE chem_boundary_conds
    870 
    871 !------------------------------------------------------------------------------!
    872 ! Description:
    873 ! ------------
    874 !> Subroutine for boundary conditions
    875 !------------------------------------------------------------------------------!
    876  SUBROUTINE chem_boundary_conditions
    877 
    878     IMPLICIT NONE
    879 
    880     INTEGER(iwp) ::  lsp             !<
    881     INTEGER(iwp) ::  lsp_usr         !<
    882 
    883 !
    884 !-- Top/bottom boundary conditions for chemical species
    885     CALL chem_boundary_conds( 'set_bc_bottomtop' )
    886 !
    887 !-- Lateral boundary conditions for chemical species
    888     CALL chem_boundary_conds( 'set_bc_lateral' )
    889 
    890 !
    891 !--  Boundary conditions for prognostic quantitites of other modules:
    892 !--  Here, only decycling is carried out
    893      DO  lsp = 1, nvar
    894         lsp_usr = 1
    895         DO WHILE ( TRIM( cs_name( lsp_usr ) ) /= 'novalue' )
    896            IF ( TRIM( chem_species(lsp)%name ) == TRIM( cs_name(lsp_usr) ) )  THEN
    897               CALL chem_boundary_conds( chem_species(lsp)%conc_p,                          &
    898                                         chem_species(lsp)%conc_pr_init )
    899            ENDIF
    900            lsp_usr = lsp_usr + 1
    901         ENDDO
    902      ENDDO
    903 
    904 
    905  END SUBROUTINE chem_boundary_conditions
    906 
    907 !------------------------------------------------------------------------------!
    908 ! Description:
    909 ! ------------
    910 !> Boundary conditions for prognostic variables in chemistry: decycling in the
    911 !> x-direction-
    912 !> Decycling of chemistry variables: Dirichlet BCs with cyclic is frequently not
    913 !> approproate for chemicals compounds since they may accumulate too much.
    914 !> If no proper boundary conditions from a DYNAMIC input file are available,
    915 !> de-cycling applies the initial profiles at the inflow boundaries for
    916 !> Dirichlet boundary conditions
    917 !------------------------------------------------------------------------------!
    918  SUBROUTINE chem_boundary_conds_decycle( cs_3d, cs_pr_init )
    919 
    920     USE control_parameters,                                                                         &
    921         ONLY:  nesting_offline
    922 
    923     INTEGER(iwp) ::  boundary  !<
    924     INTEGER(iwp) ::  ee        !<
    925     INTEGER(iwp) ::  copied    !<
    926     INTEGER(iwp) ::  i         !<
    927     INTEGER(iwp) ::  j         !<
    928     INTEGER(iwp) ::  k         !<
    929     INTEGER(iwp) ::  ss        !<
    930 
    931     REAL(wp), DIMENSION(nzb:nzt+1) ::  cs_pr_init
    932     REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  cs_3d
    933     REAL(wp) ::  flag          !< flag to mask topography grid points
    934 
    935 
    936     flag = 0.0_wp
    937 !
    938 !-- Skip input if forcing from a larger-scale model is applied
    939     IF ( nesting_offline  .AND.  nesting_offline_chem )  RETURN
    940 !
    941 !-- Left and right boundaries
    942     IF ( decycle_chem_lr  .AND.  bc_lr_cyc )  THEN
    943 
    944        DO  boundary = 1, 2
    945 
    946           IF ( decycle_method(boundary) == 'dirichlet' )  THEN
    947 !
    948 !--          Initial profile is copied to ghost and first three layers
    949              ss = 1
    950              ee = 0
    951              IF ( boundary == 1  .AND.  nxl == 0 )  THEN
    952                 ss = nxlg
    953                 ee = nxl-1
    954              ELSEIF ( boundary == 2  .AND.  nxr == nx )  THEN
    955                 ss = nxr+1
    956                 ee = nxrg
    957              ENDIF
    958 
    959              DO  i = ss, ee
    960                 DO  j = nysg, nyng
    961                    DO  k = nzb+1, nzt
    962                       flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
    963                       cs_3d(k,j,i) = cs_pr_init(k) * flag
    964                    ENDDO
    965                 ENDDO
    966              ENDDO
    967 
    968           ELSEIF ( decycle_method(boundary) == 'neumann' )  THEN
    969 !
    970 !--          The value at the boundary is copied to the ghost layers to simulate
    971 !--          an outlet with zero gradient
    972              ss = 1
    973              ee = 0
    974              IF ( boundary == 1  .AND.  nxl == 0 )  THEN
    975                 ss = nxlg
    976                 ee = nxl-1
    977                 copied = nxl
    978              ELSEIF ( boundary == 2  .AND.  nxr == nx )  THEN
    979                 ss = nxr+1
    980                 ee = nxrg
    981                 copied = nxr
    982              ENDIF
    983 
    984              DO  i = ss, ee
    985                 DO  j = nysg, nyng
    986                    DO  k = nzb+1, nzt
    987                       flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
    988                       cs_3d(k,j,i) = cs_3d(k,j,copied) * flag
    989                    ENDDO
    990                 ENDDO
    991              ENDDO
    992 
    993           ELSE
    994              WRITE(message_string,*)                                           &
    995                   'unknown decycling method: decycle_method (', &
    996                   boundary, ') ="' // TRIM( decycle_method(boundary) ) // '"'
    997              CALL message( 'chem_boundary_conds_decycle', 'CM0431',           &
    998                   1, 2, 0, 6, 0 )
     771       ENDDO    ! end lsp loop 
     772
     773!
     774!--    Top boundary conditions for chemical species - Should this not be done for all species?
     775!>     TODO: This question also needs to be clarified. I guess it can be removed because the loop
     776!>           already runs over all species? (Siggi)
     777       DO  lsp = 1, nspec
     778          IF ( ibc_cs_t == 0 )  THEN                   
     779             chem_species(lsp)%conc_p(nzt+1,:,:) = chem_species(lsp)%conc(nzt+1,:,:)
     780          ELSEIF ( ibc_cs_t == 1 )  THEN
     781             chem_species(lsp)%conc_p(nzt+1,:,:) = chem_species(lsp)%conc_p(nzt,:,:)
     782          ELSEIF ( ibc_cs_t == 2 )  THEN
     783             chem_species(lsp)%conc_p(nzt+1,:,:) = chem_species(lsp)%conc_p(nzt,:,:) +             &
     784                                                   bc_cs_t_val(lsp) * dzu(nzt+1)
    999785          ENDIF
    1000786       ENDDO
     787
     788!
     789!--    Lateral boundary conditions.
     790!--    Dirichlet conditions have been already set when chem_species concentration is initialized.
     791!--    The initially set value is not touched during time integration, hence, this boundary value
     792!--    remains at a constant value.
     793!--    Neumann conditions:
     794       IF ( bc_radiation_cs_s )  THEN
     795          DO  lsp = 1, nspec
     796             chem_species(lsp)%conc_p(:,nys-1,:) = chem_species(lsp)%conc_p(:,nys,:)
     797          ENDDO
     798       ENDIF
     799       IF ( bc_radiation_cs_n )  THEN
     800          DO  lsp = 1, nspec
     801             chem_species(lsp)%conc_p(:,nyn+1,:) = chem_species(lsp)%conc_p(:,nyn,:)
     802          ENDDO
     803       ENDIF
     804       IF ( bc_radiation_cs_l )  THEN
     805          DO  lsp = 1, nspec
     806             chem_species(lsp)%conc_p(:,:,nxl-1) = chem_species(lsp)%conc_p(:,:,nxl)
     807          ENDDO
     808       ENDIF
     809       IF ( bc_radiation_cs_r )  THEN
     810          DO  lsp = 1, nspec
     811             chem_species(lsp)%conc_p(:,:,nxr+1) = chem_species(lsp)%conc_p(:,:,nxr)
     812          ENDDO
     813       ENDIF
     814
     815!--    For testing: set Dirichlet conditions for all boundaries
     816!      IF ( bc_dirichlet_cs_s )  THEN
     817!         DO  lsp = 1, nspec
     818!            DO  i = nxlg, nxrg
     819!               DO  j = nysg, nys-1
     820!                  DO  k = nzb, nzt
     821!                     IF ( k /= nzb )  THEN
     822!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     823!                     ELSE
     824!                        flag = 1.0
     825!                     ENDIF
     826!                     chem_species(lsp)%conc_p(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     827!                  ENDDO
     828!               ENDDO
     829!            ENDDO
     830!         ENDDO
     831!      ENDIF
     832!      IF ( bc_dirichlet_cs_n )  THEN
     833!         DO  lsp = 1, nspec
     834!            DO  i = nxlg, nxrg
     835!               DO  j = nyn+1, nyng
     836!                  DO  k = nzb, nzt
     837!                     IF ( k /= nzb )  THEN
     838!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     839!                     ELSE
     840!                        flag = 1.0
     841!                     ENDIF
     842!                     chem_species(lsp)%conc_p(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     843!                  ENDDO
     844!               ENDDO
     845!            ENDDO
     846!         ENDDO
     847!      ENDIF
     848!      IF ( bc_dirichlet_cs_l )  THEN
     849!         DO  lsp = 1, nspec
     850!            DO  i = nxlg, nxl-1
     851!               DO  j = nysg, nyng
     852!                  DO  k = nzb, nzt
     853!                     IF ( k /= nzb )  THEN
     854!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     855!                     ELSE
     856!                        flag = 1.0
     857!                     ENDIF
     858!                     chem_species(lsp)%conc_p(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     859!                  ENDDO
     860!               ENDDO
     861!            ENDDO
     862!         ENDDO
     863!      ENDIF
     864!      IF ( bc_dirichlet_cs_r )  THEN
     865!         DO  lsp = 1, nspec
     866!            DO  i = nxr+1, nxrg
     867!               DO  j = nysg, nyng
     868!                  DO  k = nzb, nzt
     869!                     IF ( k /= nzb )  THEN
     870!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     871!                     ELSE
     872!                        flag = 1.0
     873!                     ENDIF
     874!                     chem_species(lsp)%conc_p(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     875!                  ENDDO
     876!               ENDDO
     877!            ENDDO
     878!         ENDDO
     879!      ENDIF
     880
     881    ELSE
     882!
     883!--    Lateral Neumann booundary conditions for timelevel t.
     884!--    This branch is executed when routine is called after the non-advective processes /
     885!--    before the prognostic equations.
     886       IF ( bc_radiation_cs_s )  THEN
     887          DO  lsp = 1, nspec
     888             chem_species(lsp)%conc(:,nys-1,:) = chem_species(lsp)%conc(:,nys,:)
     889          ENDDO
     890       ENDIF
     891       IF ( bc_radiation_cs_n )  THEN
     892          DO  lsp = 1, nspec
     893             chem_species(lsp)%conc(:,nyn+1,:) = chem_species(lsp)%conc(:,nyn,:)
     894          ENDDO
     895       ENDIF
     896       IF ( bc_radiation_cs_l )  THEN
     897          DO  lsp = 1, nspec
     898             chem_species(lsp)%conc(:,:,nxl-1) = chem_species(lsp)%conc(:,:,nxl)
     899          ENDDO
     900       ENDIF
     901       IF ( bc_radiation_cs_r )  THEN
     902          DO  lsp = 1, nspec
     903             chem_species(lsp)%conc(:,:,nxr+1) = chem_species(lsp)%conc(:,:,nxr)
     904          ENDDO
     905       ENDIF
     906
     907!--    For testing: set Dirichlet conditions for all boundaries
     908!      IF ( bc_dirichlet_cs_s )  THEN
     909!         DO  lsp = 1, nspec
     910!            DO  i = nxlg, nxrg
     911!               DO  j = nysg, nys-1
     912!                  DO  k = nzb, nzt
     913!                     IF ( k /= nzb )  THEN
     914!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     915!                     ELSE
     916!                        flag = 1.0
     917!                     ENDIF
     918!                     chem_species(lsp)%conc(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     919!                  ENDDO
     920!               ENDDO
     921!            ENDDO
     922!         ENDDO
     923!      ENDIF
     924!      IF ( bc_dirichlet_cs_n )  THEN
     925!         DO  lsp = 1, nspec
     926!            DO  i = nxlg, nxrg
     927!               DO  j = nyn+1, nyng
     928!                  DO  k = nzb, nzt
     929!                     IF ( k /= nzb )  THEN
     930!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     931!                     ELSE
     932!                        flag = 1.0
     933!                     ENDIF
     934!                     chem_species(lsp)%conc(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     935!                  ENDDO
     936!               ENDDO
     937!            ENDDO
     938!         ENDDO
     939!      ENDIF
     940!      IF ( bc_dirichlet_cs_l )  THEN
     941!         DO  lsp = 1, nspec
     942!            DO  i = nxlg, nxl-1
     943!               DO  j = nysg, nyng
     944!                  DO  k = nzb, nzt
     945!                     IF ( k /= nzb )  THEN
     946!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     947!                     ELSE
     948!                        flag = 1.0
     949!                     ENDIF
     950!                     chem_species(lsp)%conc(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     951!                  ENDDO
     952!               ENDDO
     953!            ENDDO
     954!         ENDDO
     955!      ENDIF
     956!      IF ( bc_dirichlet_cs_r )  THEN
     957!         DO  lsp = 1, nspec
     958!            DO  i = nxr+1, nxrg
     959!               DO  j = nysg, nyng
     960!                  DO  k = nzb, nzt
     961!                     IF ( k /= nzb )  THEN
     962!                        flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
     963!                     ELSE
     964!                        flag = 1.0
     965!                     ENDIF
     966!                     chem_species(lsp)%conc(k,j,i) = chem_species(lsp)%conc_pr_init(k) * flag
     967!                  ENDDO
     968!               ENDDO
     969!            ENDDO
     970!         ENDDO
     971!      ENDIF
     972
    1001973    ENDIF
    1002 !
    1003 !-- South and north boundaries
    1004     IF ( decycle_chem_ns  .AND.  bc_ns_cyc )  THEN
    1005 
    1006        DO  boundary = 3, 4
    1007 
    1008           IF ( decycle_method(boundary) == 'dirichlet' )  THEN
    1009 !
    1010 !--          Initial profile is copied to ghost and first three layers
    1011              ss = 1
    1012              ee = 0
    1013              IF ( boundary == 3  .AND.  nys == 0 )  THEN
    1014                 ss = nysg
    1015                 ee = nys-1
    1016              ELSEIF ( boundary == 4  .AND.  nyn == ny )  THEN
    1017                 ss = nyn+1
    1018                 ee = nyng
    1019              ENDIF
    1020 
    1021              DO  i = nxlg, nxrg
    1022                 DO  j = ss, ee
    1023                    DO  k = nzb+1, nzt
    1024                       flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
    1025                       cs_3d(k,j,i) = cs_pr_init(k) * flag
    1026                    ENDDO
    1027                 ENDDO
    1028              ENDDO
    1029 
    1030 
    1031           ELSEIF ( decycle_method(boundary) == 'neumann' )  THEN
    1032 !
    1033 !--          The value at the boundary is copied to the ghost layers to simulate
    1034 !--          an outlet with zero gradient
    1035              ss = 1
    1036              ee = 0
    1037              IF ( boundary == 3  .AND.  nys == 0 )  THEN
    1038                 ss = nysg
    1039                 ee = nys-1
    1040                 copied = nys
    1041              ELSEIF ( boundary == 4  .AND.  nyn == ny )  THEN
    1042                 ss = nyn+1
    1043                 ee = nyng
    1044                 copied = nyn
    1045              ENDIF
    1046 
    1047              DO  i = nxlg, nxrg
    1048                 DO  j = ss, ee
    1049                    DO  k = nzb+1, nzt
    1050                       flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
    1051                       cs_3d(k,j,i) = cs_3d(k,copied,i) * flag
    1052                    ENDDO
    1053                 ENDDO
    1054              ENDDO
    1055 
    1056           ELSE
    1057              WRITE(message_string,*)                                           &
    1058                   'unknown decycling method: decycle_method (', &
    1059                   boundary, ') ="' // TRIM( decycle_method(boundary) ) // '"'
    1060              CALL message( 'chem_boundary_conds_decycle', 'CM0432',           &
    1061                   1, 2, 0, 6, 0 )
    1062           ENDIF
    1063        ENDDO
    1064     ENDIF
    1065 
    1066 
    1067  END SUBROUTINE chem_boundary_conds_decycle
     974
     975 END SUBROUTINE chem_boundary_conditions
    1068976
    1069977
     
    11941102 SUBROUTINE chem_check_parameters
    11951103
     1104    USE control_parameters,                                                                        &
     1105        ONLY:  bc_lr, bc_ns
    11961106
    11971107    LOGICAL  ::  found
     
    12231133    ENDIF
    12241134!
    1225 !-- check for  decycling of chem species
    1226     IF ( (.NOT. any(decycle_method == 'neumann') ) .AND. (.NOT. any(decycle_method == 'dirichlet') ) )  THEN
    1227        message_string = 'Incorrect boundary conditions. Only neumann or '   &
    1228                 // 'dirichlet &available for decycling chemical species '
    1229        CALL message( 'chem_check_parameters', 'CM0426', 1, 2, 0, 6, 0 )
    1230     ENDIF
    1231 !
    12321135!-- check for chemical mechanism used
    12331136    CALL get_mechanism_name
     
    12361139       CALL message( 'chem_check_parameters', 'CM0462', 1, 2, 0, 6, 0 )
    12371140    ENDIF
     1141
     1142!
     1143!-- Check bottom boundary condition and set internal steering parameter
     1144    IF ( bc_cs_b == 'dirichlet' )  THEN
     1145       ibc_cs_b = 0
     1146    ELSEIF ( bc_cs_b == 'neumann' )  THEN
     1147       ibc_cs_b = 1
     1148    ELSE
     1149       message_string = 'unknown boundary condition: bc_cs_b ="' // TRIM( bc_cs_b ) // '"'
     1150       CALL message( 'chem_boundary_conds', 'CM0429', 1, 2, 0, 6, 0 )
     1151    ENDIF
     1152!
     1153!-- Check top boundary condition and set internal steering parameter
     1154    IF ( bc_cs_t == 'dirichlet' )  THEN
     1155       ibc_cs_t = 0
     1156    ELSEIF ( bc_cs_t == 'neumann' )  THEN
     1157       ibc_cs_t = 1
     1158    ELSEIF ( bc_cs_t == 'initial_gradient' )  THEN
     1159       ibc_cs_t = 2
     1160    ELSEIF ( bc_cs_t == 'nested' )  THEN
     1161       ibc_cs_t = 3
     1162    ELSE
     1163       message_string = 'unknown boundary condition: bc_c_t ="' // TRIM( bc_cs_t ) // '"'
     1164       CALL message( 'check_parameters', 'CM0430', 1, 2, 0, 6, 0 )
     1165    ENDIF
     1166
    12381167!
    12391168!-- If nesting_chem = .F., set top boundary condition to its default value
     
    12421171       bc_cs_t = 'initial_gradient'
    12431172    ENDIF
     1173
     1174!
     1175!-- Check left and right boundary conditions. First set default value if not set by user.
     1176    IF ( bc_cs_l == 'undefined' )  THEN
     1177       IF ( bc_lr == 'cyclic' )  THEN
     1178          bc_cs_l = 'cyclic'
     1179       ELSEIF ( bc_lr == 'dirichlet/radiation' )  THEN
     1180          bc_cs_l = 'dirichlet'
     1181       ELSEIF ( bc_lr == 'radiation/dirichlet' )  THEN
     1182          bc_cs_l = 'neumann'
     1183       ENDIF
     1184    ENDIF
     1185    IF ( bc_cs_r == 'undefined' )  THEN
     1186       IF ( bc_lr == 'cyclic' )  THEN
     1187          bc_cs_r = 'cyclic'
     1188       ELSEIF ( bc_lr == 'dirichlet/radiation' )  THEN
     1189          bc_cs_r = 'neumann'
     1190       ELSEIF ( bc_lr == 'radiation/dirichlet' )  THEN
     1191          bc_cs_r = 'dirichlet'
     1192       ENDIF
     1193    ENDIF
     1194    IF ( bc_cs_l /= 'dirichlet'  .AND.  bc_cs_l /= 'neumann'  .AND.  bc_cs_l /= 'cyclic' )  THEN
     1195       message_string = 'unknown boundary condition: bc_cs_l = "' // TRIM( bc_cs_l ) // '"'
     1196       CALL message( 'chem_check_parameters','PA0505', 1, 2, 0, 6, 0 )
     1197    ENDIF
     1198    IF ( bc_cs_r /= 'dirichlet'  .AND.  bc_cs_r /= 'neumann'  .AND.  bc_cs_r /= 'cyclic' )  THEN
     1199       message_string = 'unknown boundary condition: bc_cs_r = "' // TRIM( bc_cs_r ) // '"'
     1200       CALL message( 'chem_check_parameters','PA0551', 1, 2, 0, 6, 0 )
     1201    ENDIF
     1202!
     1203!-- Check north and south boundary conditions. First set default value if not set by user.
     1204    IF ( bc_cs_n == 'undefined' )  THEN
     1205       IF ( bc_ns == 'cyclic' )  THEN
     1206          bc_cs_n = 'cyclic'
     1207       ELSEIF ( bc_ns == 'dirichlet/radiation' )  THEN
     1208          bc_cs_n = 'dirichlet'
     1209       ELSEIF ( bc_ns == 'radiation/dirichlet' )  THEN
     1210          bc_cs_n = 'neumann'
     1211       ENDIF
     1212    ENDIF
     1213    IF ( bc_cs_s == 'undefined' )  THEN
     1214       IF ( bc_ns == 'cyclic' )  THEN
     1215          bc_cs_s = 'cyclic'
     1216       ELSEIF ( bc_ns == 'dirichlet/radiation' )  THEN
     1217          bc_cs_s = 'neumann'
     1218       ELSEIF ( bc_ns == 'radiation/dirichlet' )  THEN
     1219          bc_cs_s = 'dirichlet'
     1220       ENDIF
     1221    ENDIF
     1222    IF ( bc_cs_n /= 'dirichlet'  .AND.  bc_cs_n /= 'neumann'  .AND.  bc_cs_n /= 'cyclic' )  THEN
     1223       message_string = 'unknown boundary condition: bc_cs_n = "' // TRIM( bc_cs_n ) // '"'
     1224       CALL message( 'chem_check_parameters','PA0714', 1, 2, 0, 6, 0 )
     1225    ENDIF
     1226    IF ( bc_cs_s /= 'dirichlet'  .AND.  bc_cs_s /= 'neumann'  .AND.  bc_cs_s /= 'cyclic' )  THEN
     1227       message_string = 'unknown boundary condition: bc_cs_s = "' // TRIM( bc_cs_s ) // '"'
     1228       CALL message( 'chem_check_parameters','PA0715', 1, 2, 0, 6, 0 )
     1229    ENDIF
     1230!
     1231!-- Cyclic conditions must be set identically at opposing boundaries
     1232    IF ( ( bc_cs_l == 'cyclic' .AND. bc_cs_r /= 'cyclic' )  .OR.                                 &
     1233         ( bc_cs_r == 'cyclic' .AND. bc_cs_l /= 'cyclic' ) )  THEN
     1234       message_string = 'boundary conditions bc_cs_l and bc_cs_r must both be cyclic or non-cyclic'
     1235       CALL message( 'chem_check_parameters','PA0716', 1, 2, 0, 6, 0 )
     1236    ENDIF
     1237    IF ( ( bc_cs_n == 'cyclic' .AND. bc_cs_s /= 'cyclic' )  .OR.                                 &
     1238         ( bc_cs_s == 'cyclic' .AND. bc_cs_n /= 'cyclic' ) )  THEN
     1239       message_string = 'boundary conditions bc_cs_n and bc_cs_s must both be cyclic or non-cyclic'
     1240       CALL message( 'chem_check_parameters','PA0717', 1, 2, 0, 6, 0 )
     1241    ENDIF
     1242!
     1243!-- Set the switches that control application of horizontal boundary conditions at the boundaries
     1244!-- of the total domain
     1245    IF ( bc_cs_n == 'dirichlet'  .AND.  nyn == ny )  bc_dirichlet_cs_n = .TRUE.
     1246    IF ( bc_cs_n == 'neumann'    .AND.  nyn == ny )  bc_radiation_cs_n = .TRUE.
     1247    IF ( bc_cs_s == 'dirichlet'  .AND.  nys ==  0 )  bc_dirichlet_cs_s = .TRUE.
     1248    IF ( bc_cs_s == 'neumann'    .AND.  nys ==  0 )  bc_radiation_cs_s = .TRUE.
     1249    IF ( bc_cs_l == 'dirichlet'  .AND.  nxl ==  0 )  bc_dirichlet_cs_l = .TRUE.
     1250    IF ( bc_cs_l == 'neumann'    .AND.  nxl ==  0 )  bc_radiation_cs_l = .TRUE.
     1251    IF ( bc_cs_r == 'dirichlet'  .AND.  nxr == nx )  bc_dirichlet_cs_r = .TRUE.
     1252    IF ( bc_cs_r == 'neumann'    .AND.  nxr == nx )  bc_radiation_cs_r = .TRUE.
     1253
     1254!
     1255!-- Set the communicator to be used for ghost layer data exchange
     1256!-- 1: cyclic, 2: cyclic along x, 3: cyclic along y, 4: non-cyclic
     1257    IF ( bc_cs_l == 'cyclic' )  THEN
     1258       IF ( bc_cs_s == 'cyclic' )  THEN
     1259          communicator_chem = 1
     1260       ELSE
     1261          communicator_chem = 2
     1262       ENDIF
     1263    ELSE
     1264       IF ( bc_cs_s == 'cyclic' )  THEN
     1265          communicator_chem = 3
     1266       ELSE
     1267          communicator_chem = 4
     1268       ENDIF
     1269    ENDIF
     1270
    12441271!
    12451272!-- chem_check_parameters is called before the array chem_species is allocated!
    12461273!-- temporary switch of this part of the check
    12471274!    RETURN                !bK commented
    1248 
     1275!> TODO: this workaround definitely needs to be removed from here!!!
    12491276    CALL chem_init_internal
    12501277!
     
    17321759!
    17331760!-- initializatoin of Surface and profile chemical species
    1734 
    17351761    lsp = 1
    17361762    docsinit_chr ='Chemical species for initial surface and profile emissions: '
     
    17501776    IF ( nesting_chem )  WRITE( io, 12 )  nesting_chem
    17511777    IF ( nesting_offline_chem )  WRITE( io, 13 )  nesting_offline_chem
     1778
     1779    WRITE( io, 14 )  TRIM( bc_cs_b ), TRIM( bc_cs_t ), TRIM( bc_cs_s ), TRIM( bc_cs_n ),           &
     1780                     TRIM( bc_cs_l ), TRIM( bc_cs_r )
     1781
    17521782!
    17531783!-- number of variable and fix chemical species and number of reactions
    17541784    cs_fixed = nspec - nvar
    1755 
    17561785    WRITE ( io, * ) '   --> Chemical Mechanism        : ', cs_mech
    17571786    WRITE ( io, * ) '   --> Chemical species, variable: ', nvar
     
    1774180312  FORMAT (/'   Nesting for chemistry variables: ', L1 )
    1775180413  FORMAT (/'   Offline nesting for chemistry variables: ', L1 )
    1776 !
    1777 !
     180514  FORMAT (/'   Boundary conditions for chemical species:', /                                     &
     1806             '      bottom/top:   ',A10,' / ',A10, /                                               &
     1807             '      north/south:  ',A10,' / ',A10, /                                               &
     1808             '      left/right:   ',A10,' / ',A10)
     1809
    17781810 END SUBROUTINE chem_header
    17791811
     
    18841916    INTEGER(iwp) ::  lpr_lev           !< running index for chem spcs profile level
    18851917
     1918    REAL(wp)     ::  flag              !< flag for masking topography/building grid points
    18861919!
    18871920!-- 20200203 ECC
     
    19481981
    19491982    ENDDO
    1950 !
    1951 !-- Set control flags for decycling only at lateral boundary cores, within the
    1952 !-- inner cores the decycle flag is set to .False.. Even though it does not
    1953 !-- affect the setting of chemistry boundary conditions, this flag is used to
    1954 !-- set advection control flags appropriately.
    1955     decycle_chem_lr = MERGE( decycle_chem_lr, .FALSE.,                         &
    1956                              nxl == 0  .OR.  nxr == nx )
    1957     decycle_chem_ns = MERGE( decycle_chem_ns, .FALSE.,                         &
    1958                              nys == 0  .OR.  nyn == ny )
     1983
    19591984!
    19601985!-- For some passive scalars decycling may be enabled. This case, the lateral
     
    19822007!--    Initialize with flag 31 only.
    19832008       cs_advc_flags_s = 0
    1984        cs_advc_flags_s = MERGE( IBSET( cs_advc_flags_s, 31 ), 0,               &
    1985                                 BTEST( wall_flags_total_0, 31 ) )
    1986 
    1987        IF ( decycle_chem_ns )  THEN
     2009       cs_advc_flags_s = MERGE( IBSET( cs_advc_flags_s, 31 ), 0, BTEST( wall_flags_total_0, 31 ) )
     2010
     2011       IF ( bc_dirichlet_cs_n .OR. bc_dirichlet_cs_s )  THEN
    19882012          IF ( nys == 0  )  THEN
    1989              DO  i = 1, nbgp     
     2013             DO  i = 1, nbgp
    19902014                cs_advc_flags_s(:,nys-i,:) = MERGE(                            &
    19912015                                        IBSET( cs_advc_flags_s(:,nys,:), 31 ), &
     
    19962020          ENDIF
    19972021          IF ( nyn == ny )  THEN
    1998              DO  i = 1, nbgp 
     2022             DO  i = 1, nbgp
    19992023                cs_advc_flags_s(:,nyn+i,:) = MERGE(                            &
    20002024                                        IBSET( cs_advc_flags_s(:,nyn,:), 31 ), &
     
    20052029          ENDIF
    20062030       ENDIF
    2007        IF ( decycle_chem_lr )  THEN
     2031       IF ( bc_dirichlet_cs_l .OR. bc_dirichlet_cs_r )  THEN
    20082032          IF ( nxl == 0  )  THEN
    2009              DO  i = 1, nbgp   
     2033             DO  i = 1, nbgp
    20102034                cs_advc_flags_s(:,:,nxl-i) = MERGE(                            &
    20112035                                        IBSET( cs_advc_flags_s(:,:,nxl), 31 ), &
     
    20152039             ENDDO
    20162040          ENDIF
    2017           IF ( nxr == nx )  THEN 
    2018              DO  i = 1, nbgp   
     2041          IF ( nxr == nx )  THEN
     2042             DO  i = 1, nbgp
    20192043                cs_advc_flags_s(:,:,nxr+i) = MERGE(                            &
    20202044                                        IBSET( cs_advc_flags_s(:,:,nxr), 31 ), &
     
    20232047                                                  )
    20242048             ENDDO
    2025           ENDIF     
     2049          ENDIF
    20262050       ENDIF
    20272051!
     
    20392063!--    oscillations, which are responsible for high concentration maxima that may
    20402064!--    appear under shear-free stable conditions.
    2041        CALL ws_init_flags_scalar(                                              &
    2042                   bc_dirichlet_l  .OR.  bc_radiation_l  .OR.  ( decycle_chem_lr .AND. nxl == 0  ), &
    2043                   bc_dirichlet_n  .OR.  bc_radiation_n  .OR.  ( decycle_chem_ns .AND. nyn == ny ), &
    2044                   bc_dirichlet_r  .OR.  bc_radiation_r  .OR.  ( decycle_chem_lr .AND. nxr == nx ), &
    2045                   bc_dirichlet_s  .OR.  bc_radiation_s  .OR.  ( decycle_chem_ns .AND. nys == 0  ), &
    2046                   cs_advc_flags_s, .TRUE. )
     2065       CALL ws_init_flags_scalar( bc_dirichlet_cs_l  .OR.  bc_radiation_cs_l,                      &
     2066                                  bc_dirichlet_cs_n  .OR.  bc_radiation_cs_n,                      &
     2067                                  bc_dirichlet_cs_r  .OR.  bc_radiation_cs_r,                      &
     2068                                  bc_dirichlet_cs_s  .OR.  bc_radiation_cs_s,                      &
     2069                                  cs_advc_flags_s, .TRUE. )
    20472070    ENDIF
    20482071!
     
    20742097!
    20752098!--       Transfer initial profiles to the arrays of the 3D model
     2099!--       Concentrations within buildings are set to zero.
    20762100          DO  lsp = 1, nspec
    20772101             DO  i = nxlg, nxrg
    20782102                DO  j = nysg, nyng
    20792103                   DO lpr_lev = 1, nz + 1
    2080                       chem_species(lsp)%conc(lpr_lev,j,i) = chem_species(lsp)%conc_pr_init(lpr_lev)
     2104                      flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(lpr_lev,j,i), 0 ) )
     2105                      chem_species(lsp)%conc(lpr_lev,j,i) = chem_species(lsp)%conc_pr_init(lpr_lev)&
     2106                                                            * flag
    20812107                   ENDDO
    20822108                ENDDO
     
    20902116             DO  i = nxlg, nxrg
    20912117                DO  j = nysg, nyng
    2092                    chem_species(lsp)%conc(:,j,i) = chem_species(lsp)%conc_pr_init   
     2118                   DO  lpr_lev = nzb, nz+1
     2119                      flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(lpr_lev,j,i), 0 ) )
     2120                      chem_species(lsp)%conc(lpr_lev,j,i) = chem_species(lsp)%conc_pr_init(lpr_lev)&
     2121                                                            * flag
     2122                   ENDDO
    20932123                ENDDO
    20942124             ENDDO
     
    21002130       IF ( cs_surface_initial_change(1) /= 0.0_wp )  THEN           
    21012131          DO  lsp = 1, nspec
    2102              chem_species(lsp)%conc(nzb,:,:) = chem_species(lsp)%conc(nzb,:,:) +  &
    2103                                                cs_surface_initial_change(lsp)
     2132             DO  i = nxlg, nxrg
     2133                DO  j = nysg, nyng
     2134                   flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(nzb,j,i), 0 ) )
     2135                   chem_species(lsp)%conc(nzb,j,i) = chem_species(lsp)%conc(nzb,j,i) +  &
     2136                                                     cs_surface_initial_change(lsp) * flag
     2137                ENDDO
     2138             ENDDO
    21042139          ENDDO
    21052140       ENDIF
     
    23462381
    23472382    NAMELIST /chemistry_parameters/  bc_cs_b,                          &
     2383         bc_cs_l,                          &
     2384         bc_cs_n,                          &
     2385         bc_cs_r,                          &
     2386         bc_cs_s,                          &
    23482387         bc_cs_t,                          &
    23492388         call_chem_at_all_substeps,        &
     
    23602399         cs_vertical_gradient_level,       &
    23612400         daytype_mdh,                      &
    2362          decycle_chem_lr,                  &
    2363          decycle_chem_ns,                  &           
    2364          decycle_method,                   &
    23652401         deposition_dry,                   &
    23662402         emissions_anthropogenic,          &
     
    25902626       solver_type = 'Rang3'
    25912627    ELSE
    2592        message_string = 'illegal solver type'
     2628       message_string = 'illegal Rosenbrock-solver type'
    25932629       CALL message( 'chem_parin', 'PA0506', 1, 2, 0, 6, 0 )
    25942630    END IF
     
    27742810
    27752811   INTEGER(iwp) ::  lsp       !<
    2776    INTEGER(iwp) ::  lsp_usr   !<
    27772812 
    27782813!
     
    27802815       CALL cpu_log( log_point_s(84), 'chem.exch-horiz', 'start' )
    27812816       DO  lsp = 1, nvar
    2782           CALL exchange_horiz( chem_species(lsp)%conc, nbgp )   
    2783           lsp_usr = 1 
    2784           DO WHILE ( TRIM( cs_name( lsp_usr ) ) /= 'novalue' )
    2785              IF ( TRIM(chem_species(lsp)%name) == TRIM(cs_name(lsp_usr)) )  THEN
    2786 !
    2787 !--             As chem_exchange_horiz_bounds is called at the beginning
    2788 !--             of prognostic_equations, boundary conditions are set on
    2789 !--             %conc.
    2790                 CALL chem_boundary_conds( chem_species(lsp)%conc,              &
    2791                                           chem_species(lsp)%conc_pr_init )
    2792              
    2793              ENDIF
    2794              lsp_usr = lsp_usr + 1
    2795           ENDDO
    2796 
    2797 
     2817          CALL exchange_horiz( chem_species(lsp)%conc, nbgp,                                       &
     2818                               alternative_communicator = communicator_chem )
    27982819       ENDDO
     2820
     2821       CALL chem_boundary_conditions( horizontal_conditions_only = .TRUE. )
     2822
    27992823       CALL cpu_log( log_point_s(84), 'chem.exch-horiz', 'stop' )
    28002824
     
    28322856          IF ( ws_scheme_sca )  THEN
    28332857             CALL advec_s_ws( cs_advc_flags_s, chem_species(ilsp)%conc, 'kc',                      &
    2834                   bc_dirichlet_l  .OR.  bc_radiation_l  .OR.  ( decycle_chem_lr .AND. nxl == 0  ), &
    2835                   bc_dirichlet_n  .OR.  bc_radiation_n  .OR.  ( decycle_chem_ns .AND. nyn == ny ), &
    2836                   bc_dirichlet_r  .OR.  bc_radiation_r  .OR.  ( decycle_chem_lr .AND. nxr == nx ), &
    2837                   bc_dirichlet_s  .OR.  bc_radiation_s  .OR.  ( decycle_chem_ns .AND. nys == 0  ) )
     2858                              bc_dirichlet_cs_l  .OR.  bc_radiation_cs_l,                          &
     2859                              bc_dirichlet_cs_n  .OR.  bc_radiation_cs_n,                          &
     2860                              bc_dirichlet_cs_r  .OR.  bc_radiation_cs_r,                          &
     2861                              bc_dirichlet_cs_s  .OR.  bc_radiation_cs_s )
    28382862          ELSE
    28392863             CALL advec_s_pw( chem_species(ilsp)%conc )
     
    29422966       IF ( timestep_scheme(1:5) == 'runge' )  THEN
    29432967          IF ( ws_scheme_sca )  THEN
    2944              CALL advec_s_ws( cs_advc_flags_s, i, j, chem_species(ilsp)%conc, 'kc',                &
     2968             CALL advec_s_ws( cs_advc_flags_s,                                                     &
     2969                              i,                                                                   &
     2970                              j,                                                                   &
     2971                              chem_species(ilsp)%conc,                                             &
     2972                              'kc',                                                                &
    29452973                              chem_species(ilsp)%flux_s_cs,                                        &
    29462974                              chem_species(ilsp)%diss_s_cs,                                        &
     
    29492977                              i_omp_start,                                                         &
    29502978                              tn,                                                                  &
    2951                   bc_dirichlet_l  .OR.  bc_radiation_l  .OR.  ( decycle_chem_lr .AND. nxl == 0  ), &
    2952                   bc_dirichlet_n  .OR.  bc_radiation_n  .OR.  ( decycle_chem_ns .AND. nyn == ny ), &
    2953                   bc_dirichlet_r  .OR.  bc_radiation_r  .OR.  ( decycle_chem_lr .AND. nxr == nx ), &
    2954                   bc_dirichlet_s  .OR.  bc_radiation_s  .OR.  ( decycle_chem_ns .AND. nys == 0  ), &
    2955                   monotonic_limiter_z )
     2979                              bc_dirichlet_cs_l  .OR.  bc_radiation_cs_l,                          &
     2980                              bc_dirichlet_cs_n  .OR.  bc_radiation_cs_n,                          &
     2981                              bc_dirichlet_cs_r  .OR.  bc_radiation_cs_r,                          &
     2982                              bc_dirichlet_cs_s  .OR.  bc_radiation_cs_s,                          &
     2983                              monotonic_limiter_z )
    29562984          ELSE
    29572985             CALL advec_s_pw( i, j, chem_species(ilsp)%conc )
  • palm/trunk/SOURCE/time_integration.f90

    r4508 r4511  
    2525! -----------------
    2626! $Id$
     27! chemistry decycling replaced by explicit setting of lateral boundary conditions
     28!
     29! 4508 2020-04-24 13:32:20Z raasch
    2730! salsa decycling replaced by explicit setting of lateral boundary conditions
    2831!
     
    242245
    243246    USE chem_modules,                                                                              &
    244         ONLY:  bc_cs_t_val, chem_species, emissions_anthropogenic, emiss_read_legacy_mode,         &
    245                n_matched_vars
    246 
    247 #if defined( __parallel )
    248     USE chem_modules,                                                                              &
    249         ONLY:  cs_name
    250 #endif
    251 
    252     USE chemistry_model_mod,                                                                       &
    253         ONLY:  chem_boundary_conds
     247        ONLY:  bc_cs_t_val, chem_species, communicator_chem, emissions_anthropogenic,              &
     248               emiss_read_legacy_mode, n_matched_vars
    254249
    255250    USE control_parameters,                                                                        &
     
    469464    INTEGER(iwp) ::  icc                 !< additional index for aerosol mass bins
    470465    INTEGER(iwp) ::  ig                  !< index for salsa gases
    471     INTEGER(iwp) ::  lsp                 !<
    472466    INTEGER(iwp) ::  mid                 !< masked output running index
    473 #if defined( __parallel )
    474     INTEGER(iwp) ::  lsp_usr             !<
    475467    INTEGER(iwp) ::  n                   !< loop counter for chemistry species
    476 #endif
    477468
    478469    REAL(wp) ::  dt_3d_old  !< temporary storage of timestep to be used for
     
    675666           bc_q_t_val  = ( q_init(nzt+1) - q_init(nzt) ) / dzu(nzt+1)
    676667           IF ( air_chemistry )  THEN
    677               DO  lsp = 1, nvar
    678                  bc_cs_t_val = (  chem_species(lsp)%conc_pr_init(nzt+1)                            &
    679                                 - chem_species(lsp)%conc_pr_init(nzt) )                            &
     668              DO  n = 1, nvar
     669                 bc_cs_t_val = (  chem_species(n)%conc_pr_init(nzt+1)                              &
     670                                - chem_species(n)%conc_pr_init(nzt) )                              &
    680671                               / dzu(nzt+1)
    681672              ENDDO
     
    824815          IF ( passive_scalar )  CALL exchange_horiz( s_p, nbgp )
    825816          IF ( air_chemistry )  THEN
    826              DO  lsp = 1, nvar
    827                 CALL exchange_horiz( chem_species(lsp)%conc_p, nbgp )
     817             DO  n = 1, nvar
     818                CALL exchange_horiz( chem_species(n)%conc_p, nbgp,                                 &
     819                                     alternative_communicator = communicator_chem )
    828820             ENDDO
    829821          ENDIF
     
    932924                IF ( air_chemistry )  THEN
    933925                   DO  n = 1, nvar
    934                       CALL exchange_horiz( chem_species(n)%conc, nbgp )
     926                      CALL exchange_horiz( chem_species(n)%conc, nbgp,                             &
     927                                           alternative_communicator = communicator_chem )
    935928                   ENDDO
    936929                ENDIF
     
    960953!--          Set boundary conditions again after interpolation and anterpolation.
    961954             CALL pmci_boundary_conds
    962 
    963 !
    964 !--          Set chemistry boundary conditions (decycling)
    965              IF ( air_chemistry )  THEN
    966                 DO  lsp = 1, nvar
    967                    lsp_usr = 1
    968                    DO WHILE ( TRIM( cs_name( lsp_usr ) ) /= 'novalue' )
    969                       IF ( TRIM( chem_species(lsp)%name ) == TRIM( cs_name(lsp_usr) ) )  THEN
    970                          CALL chem_boundary_conds( chem_species(lsp)%conc,                         &
    971                                                    chem_species(lsp)%conc_pr_init )
    972                       ENDIF
    973                       lsp_usr = lsp_usr + 1
    974                    ENDDO
    975                 ENDDO
    976              ENDIF
    977955
    978956             CALL cpu_log( log_point(60), 'nesting', 'stop' )
Note: See TracChangeset for help on using the changeset viewer.