Changeset 4079 for palm


Ignore:
Timestamp:
Jul 9, 2019 6:04:41 PM (5 years ago)
Author:
suehring
Message:

Implementation of a monotonic flux limiter for vertical advection term in Wicker-Skamarock scheme. The flux limiter is currently only applied for passive scalars (passive scalar, chemical species, aerosols) within the region up to the highest topography, in order to avoid the built-up of large concentrations within poorly resolved cavities in urban environments. To enable the limiter monotonic_limiter_z = .T. must be set. Note, the limiter is currently only implemented for the cache-optimized version of advec_ws. Further changes in offline nesting: Set boundary condition for w at nzt+1 at all lateral boundaries (even though these won't enter the numerical solution), in order to avoid high vertical velocities in the run-control file which might built-up due to the mass-conservation; bugfix in offline nesting for chemical species

Location:
palm/trunk/SOURCE
Files:
7 edited

Legend:

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

    r3873 r4079  
    2525! -----------------
    2626! $Id$
     27! Implementation of a flux limiter according to Skamarock (2006) for the
     28! vertical scalar advection. Please note, this is only implemented for the
     29! cache-optimized version at the moment. Implementation for the vector-
     30! optimized version will follow after critical issues concerning
     31! vectorization are fixed.
     32!
     33! 3873 2019-04-08 15:44:30Z knoop
    2734! Moved ocean_mode specific code to ocean_mod
    2835!
     
    273280!> the divergence is not sufficiently reduced, resulting in erroneous fluxes
    274281!> and could lead to numerical instabilities.
     282!>
     283!> @todo Implement monotonic flux limiter also for vector version.
    275284!------------------------------------------------------------------------------!
    276285 MODULE advec_ws
     
    296305               bc_dirichlet_s, bc_radiation_l, bc_radiation_n,                 &
    297306               bc_radiation_r, bc_radiation_s, intermediate_timestep_count,    &
    298                u_gtrans, v_gtrans
     307               u_gtrans, v_gtrans, dt_3d
    299308
    300309    USE indices,                                                               &
     
    11161125    SUBROUTINE advec_s_ws_ij( i, j, sk, sk_char, swap_flux_y_local,            &
    11171126                              swap_diss_y_local, swap_flux_x_local,            &
    1118                               swap_diss_x_local, i_omp, tn )
     1127                              swap_diss_x_local, i_omp, tn, flux_limitation )
    11191128
    11201129
     
    11261135       INTEGER(iwp) ::  k         !< grid index along z-direction
    11271136       INTEGER(iwp) ::  k_mm      !< k-2 index in disretization, can be modified to avoid segmentation faults
     1137       INTEGER(iwp) ::  k_mmm     !< k-3 index in disretization, can be modified to avoid segmentation faults
    11281138       INTEGER(iwp) ::  k_pp      !< k+2 index in disretization, can be modified to avoid segmentation faults
    11291139       INTEGER(iwp) ::  k_ppp     !< k+3 index in disretization, can be modified to avoid segmentation faults
    11301140       INTEGER(iwp) ::  nzb_max_l !< index indicating upper bound for order degradation of horizontal advection terms 
    11311141       INTEGER(iwp) ::  tn        !< number of OpenMP thread
    1132        
    1133        REAL(wp)     ::  ibit0  !< flag indicating 1st-order scheme along x-direction
    1134        REAL(wp)     ::  ibit1  !< flag indicating 3rd-order scheme along x-direction
    1135        REAL(wp)     ::  ibit2  !< flag indicating 5th-order scheme along x-direction
    1136        REAL(wp)     ::  ibit3  !< flag indicating 1st-order scheme along y-direction
    1137        REAL(wp)     ::  ibit4  !< flag indicating 3rd-order scheme along y-direction
    1138        REAL(wp)     ::  ibit5  !< flag indicating 5th-order scheme along y-direction
    1139        REAL(wp)     ::  ibit6  !< flag indicating 1st-order scheme along z-direction
    1140        REAL(wp)     ::  ibit7  !< flag indicating 3rd-order scheme along z-direction
    1141        REAL(wp)     ::  ibit8  !< flag indicating 5th-order scheme along z-direction
    1142        REAL(wp)     ::  diss_d !< artificial dissipation term at grid box bottom
    1143        REAL(wp)     ::  div    !< diverence on scalar grid
    1144        REAL(wp)     ::  flux_d !< 6th-order flux at grid box bottom
    1145        REAL(wp)     ::  u_comp !< advection velocity along x-direction
    1146        REAL(wp)     ::  v_comp !< advection velocity along y-direction
     1142
     1143       LOGICAL, OPTIONAL ::  flux_limitation !< flag indicating flux limitation of the vertical advection
     1144       LOGICAL           ::  limiter         !< control flag indicating the application of flux limitation
     1145
     1146       REAL(wp) ::  diss_d        !< artificial dissipation term at grid box bottom
     1147       REAL(wp) ::  div           !< velocity diverence on scalar grid
     1148       REAL(wp) ::  div_in        !< vertical flux divergence of ingoing fluxes
     1149       REAL(wp) ::  div_out       !< vertical flux divergence of outgoing fluxes     
     1150       REAL(wp) ::  f_corr_t      !< correction flux at grid-cell top, i.e. the difference between high and low-order flux
     1151       REAL(wp) ::  f_corr_d      !< correction flux at grid-cell bottom, i.e. the difference between high and low-order flux
     1152       REAL(wp) ::  f_corr_t_in   !< correction flux of ingoing flux part at grid-cell top
     1153       REAL(wp) ::  f_corr_d_in   !< correction flux of ingoing flux part at grid-cell bottom
     1154       REAL(wp) ::  f_corr_t_out  !< correction flux of outgoing flux part at grid-cell top
     1155       REAL(wp) ::  f_corr_d_out  !< correction flux of outgoing flux part at grid-cell bottom
     1156       REAL(wp) ::  fac_correction!< factor to limit the in- and outgoing fluxes
     1157       REAL(wp) ::  flux_d        !< 6th-order flux at grid box bottom
     1158       REAL(wp) ::  ibit0         !< flag indicating 1st-order scheme along x-direction
     1159       REAL(wp) ::  ibit1         !< flag indicating 3rd-order scheme along x-direction
     1160       REAL(wp) ::  ibit2         !< flag indicating 5th-order scheme along x-direction
     1161       REAL(wp) ::  ibit3         !< flag indicating 1st-order scheme along y-direction
     1162       REAL(wp) ::  ibit4         !< flag indicating 3rd-order scheme along y-direction
     1163       REAL(wp) ::  ibit5         !< flag indicating 5th-order scheme along y-direction
     1164       REAL(wp) ::  ibit6         !< flag indicating 1st-order scheme along z-direction
     1165       REAL(wp) ::  ibit7         !< flag indicating 3rd-order scheme along z-direction
     1166       REAL(wp) ::  ibit8         !< flag indicating 5th-order scheme along z-direction
     1167       REAL(wp) ::  max_val       !< maximum value of the quanitity along the numerical stencil (in vertical direction)
     1168       REAL(wp) ::  min_val       !< maximum value of the quanitity along the numerical stencil (in vertical direction)
     1169       REAL(wp) ::  mon           !< monotone solution of the advection equation using 1st-order fluxes
     1170       REAL(wp) ::  u_comp        !< advection velocity along x-direction
     1171       REAL(wp) ::  v_comp        !< advection velocity along y-direction
    11471172!       
    11481173!--    sk is an array from parameter list. It should not be a pointer, because
     
    11531178       REAL(wp), INTENT(IN),DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  sk !<  advected scalar
    11541179
    1155        REAL(wp), DIMENSION(nzb:nzt+1)         ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
    1156        REAL(wp), DIMENSION(nzb:nzt+1)         ::  diss_r !< discretized artificial dissipation at rightward-side of the grid box
    1157        REAL(wp), DIMENSION(nzb:nzt+1)         ::  diss_t !< discretized artificial dissipation at top of the grid box
    1158        REAL(wp), DIMENSION(nzb:nzt+1)         ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
    1159        REAL(wp), DIMENSION(nzb:nzt+1)         ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
    1160        REAL(wp), DIMENSION(nzb:nzt+1)         ::  flux_t !< discretized 6th-order flux at top of the grid box
     1180       REAL(wp), DIMENSION(nzb:nzt+1)         ::  diss_n     !< discretized artificial dissipation at northward-side of the grid box
     1181       REAL(wp), DIMENSION(nzb:nzt+1)         ::  diss_r     !< discretized artificial dissipation at rightward-side of the grid box
     1182       REAL(wp), DIMENSION(nzb:nzt+1)         ::  diss_t     !< discretized artificial dissipation at top of the grid box
     1183       REAL(wp), DIMENSION(nzb:nzt+1)         ::  flux_n     !< discretized 6th-order flux at northward-side of the grid box
     1184       REAL(wp), DIMENSION(nzb:nzt+1)         ::  flux_r     !< discretized 6th-order flux at rightward-side of the grid box
     1185       REAL(wp), DIMENSION(nzb:nzt+1)         ::  flux_t     !< discretized 6th-order flux at top of the grid box
     1186       REAL(wp), DIMENSION(nzb:nzt+1)         ::  flux_t_1st !< discretized 1st-order flux at top of the grid box
    11611187       
    11621188       REAL(wp), DIMENSION(nzb+1:nzt,0:threads_per_task-1) ::  swap_diss_y_local !< discretized artificial dissipation at southward-side of the grid box
     
    11781204          nzb_max_l = nzb_max
    11791205       END IF
     1206!
     1207!--    Set control flag for flux limiter
     1208       limiter = .FALSE.
     1209       IF ( PRESENT( flux_limitation) )  limiter = flux_limitation
    11801210!
    11811211!--    Compute southside fluxes of the respective PE bounds.
     
    12311261                                  +           ( sk(k,j+2,i) + sk(k,j-3,i) )   &
    12321262                                                ) * adv_sca_5
    1233               swap_diss_y_local(k,tn) = -ABS( v_comp ) * (                    &
     1263             swap_diss_y_local(k,tn) = -ABS( v_comp ) * (                     &
    12341264                                    10.0_wp * ( sk(k,j,i)   - sk(k,j-1,i) )   &
    12351265                                  -  5.0_wp * ( sk(k,j+1,i) - sk(k,j-2,i) )   &
    12361266                                  +             sk(k,j+2,i) - sk(k,j-3,i)     &
    1237                                                          ) * adv_sca_5
     1267                                                        ) * adv_sca_5
    12381268
    12391269          ENDDO
     
    12661296                                                  )
    12671297
    1268               swap_diss_x_local(k,j,tn) = -ABS( u_comp ) * (                  &
     1298             swap_diss_x_local(k,j,tn) = -ABS( u_comp ) * (                   &
    12691299                                               ( 10.0_wp * ibit2 * adv_sca_5  &
    12701300                                            +     3.0_wp * ibit1 * adv_sca_3  &
     
    12791309                                               ) *                            &
    12801310                                            ( sk(k,j,i+2) - sk(k,j,i-3)    )  &
    1281                                                            )
     1311                                                          )
    12821312
    12831313          ENDDO
     
    15351565                             ( sk(k_ppp,j,i) - sk(k_mm,j,i) )                 &
    15361566                                                         )
    1537        ENDDO
     1567       ENDDO
     1568
     1569       IF ( limiter )  THEN
     1570!
     1571!--       Compute monotone first-order fluxes which are required for mononte
     1572!--       flux limitation.
     1573          flux_t_1st(nzb) = 0.0_wp
     1574          DO  k = nzb+1, nzb_max_l
     1575             flux_t_1st(k) = ( w(k,j,i)   * ( sk(k+1,j,i)  + sk(k,j,i) )       &
     1576                       -  ABS( w(k,j,i) ) * ( sk(k+1,j,i)  - sk(k,j,i) ) )     &
     1577                           * rho_air_zw(k) * adv_sca_1
     1578!
     1579!--          In flux limitation the total flux will be corrected. For the sake
     1580!--          of cleariness the higher-order advective and disspative fluxes
     1581!--          will be merged onto flux_t.
     1582             flux_t(k) = flux_t(k) + diss_t(k)
     1583             diss_t(k) = 0.0_wp
     1584          ENDDO
     1585!
     1586!--       Flux limitation of vertical fluxes according to Skamarock (2006).
     1587!--       Please note, as flux limitation implies linear dependencies of fluxes,
     1588!--       flux limitation is only made for the vertical advection term. Limitation
     1589!--       of the horizontal terms cannot be parallelized.
     1590!--       Due to the linear dependency, the following loop will not be vectorized.
     1591!--       Further, note that the flux limiter is only applied within the urban
     1592!--       layer, i.e up to the topography top. 
     1593          DO  k = nzb+1, nzb_max_l
     1594!
     1595!--          Compute one-dimensional divergence along the vertical direction,
     1596!--          which is used to correct the advection discretization. This is
     1597!--          necessary as in one-dimensional space the advection velocity
     1598!--          should actually be constant.
     1599             div = ( w(k,j,i)   * rho_air_zw(k)                                &
     1600                   - w(k-1,j,i) * rho_air_zw(k-1)                              &     
     1601                   ) * drho_air(k) * ddzw(k)
     1602!
     1603!--          Compute monotone solution of the advection equation from
     1604!--          1st-order fluxes. Please note, the advection equation is corrected
     1605!--          by the divergence term (in 1D the advective flow should be divergence
     1606!--          free). Moreover, please note, as time-increment the full timestep
     1607!--          is used, even though a Runge-Kutta scheme will be used. However,
     1608!--          the length of the actual time increment is not important at all
     1609!--          since it cancels out later when the fluxes are limited.   
     1610             mon = sk(k,j,i) + ( - ( flux_t_1st(k) - flux_t_1st(k-1) )         &
     1611                             * drho_air(k) * ddzw(k)                           &
     1612                             + div * sk(k,j,i)                                 &
     1613                               ) * dt_3d
     1614!
     1615!--          Determine minimum and maximum values along the numerical stencil.
     1616             k_mmm = MAX( k - 3, nzb + 1 )
     1617             k_ppp = MIN( k + 3, nzt + 1 )
     1618
     1619             min_val = MINVAL( sk(k_mmm:k_ppp,j,i) )
     1620             max_val = MAXVAL( sk(k_mmm:k_ppp,j,i) )
     1621!
     1622!--          Compute difference between high- and low-order fluxes, which may
     1623!--          act as correction fluxes
     1624             f_corr_t = flux_t(k)   - flux_t_1st(k)
     1625             f_corr_d = flux_t(k-1) - flux_t_1st(k-1)
     1626!
     1627!--          Determine outgoing fluxes, i.e. the part of the fluxes which can
     1628!--          decrease the value within the grid box
     1629             f_corr_t_out = MAX( 0.0_wp, f_corr_t )
     1630             f_corr_d_out = MIN( 0.0_wp, f_corr_d )
     1631!
     1632!--          Determine ingoing fluxes, i.e. the part of the fluxes which can
     1633!--          increase the value within the grid box
     1634             f_corr_t_in = MIN( 0.0_wp, f_corr_t)
     1635             f_corr_d_in = MAX( 0.0_wp, f_corr_d)
     1636!
     1637!--          Compute divergence of outgoing correction fluxes
     1638             div_out = - ( f_corr_t_out - f_corr_d_out ) * drho_air(k)         &
     1639                                                         * ddzw(k) * dt_3d
     1640!
     1641!--          Compute divergence of ingoing correction fluxes
     1642             div_in = - ( f_corr_t_in - f_corr_d_in )    * drho_air(k)         &
     1643                                                         * ddzw(k) * dt_3d
     1644!
     1645!--          Check if outgoing fluxes can lead to undershoots, i.e. values smaller
     1646!--          than the minimum value within the numerical stencil. If so, limit
     1647!--          them.
     1648             IF ( mon - min_val < - div_out  .AND.  ABS( div_out ) > 0.0_wp )  &
     1649             THEN
     1650                fac_correction = ( mon - min_val ) / ( - div_out )
     1651                f_corr_t_out = f_corr_t_out * fac_correction
     1652                f_corr_d_out = f_corr_d_out * fac_correction
     1653             ENDIF
     1654!
     1655!--          Check if ingoing fluxes can lead to overshoots, i.e. values larger
     1656!--          than the maximum value within the numerical stencil. If so, limit
     1657!--          them.
     1658             IF ( mon - max_val > - div_in  .AND.  ABS( div_in ) > 0.0_wp )    &
     1659             THEN
     1660                fac_correction = ( mon - max_val ) / ( - div_in )
     1661                f_corr_t_in = f_corr_t_in * fac_correction
     1662                f_corr_d_in = f_corr_d_in * fac_correction
     1663             ENDIF
     1664!
     1665!--          Finally add the limited fluxes to the original ones. If no
     1666!--          flux limitation was done, the fluxes equal the original ones.
     1667             flux_t(k)   = flux_t_1st(k)   + f_corr_t_out + f_corr_t_in
     1668             flux_t(k-1) = flux_t_1st(k-1) + f_corr_d_out + f_corr_d_in
     1669          ENDDO
     1670       ENDIF
    15381671
    15391672       DO  k = nzb+1, nzb_max_l
  • palm/trunk/SOURCE/chemistry_model_mod.f90

    r4069 r4079  
    2727! -----------------
    2828! $Id$
     29! Application of monotonic flux limiter for the vertical scalar advection
     30! up to the topography top (only for the cache-optimized version at the
     31! moment).
     32!
     33! 4069 2019-07-01 14:05:51Z Giersch
    2934! Masked output running index mid has been introduced as a local variable to
    3035! avoid runtime error (Loop variable has been modified) in time_integration
     
    354359                dt_3d, humidity, initializing_actions, message_string,                             &
    355360                omega, tsc, intermediate_timestep_count, intermediate_timestep_count_max,          &
    356          max_pr_user, timestep_scheme, use_prescribed_profile_data, ws_scheme_sca, air_chemistry
     361                max_pr_user,                                                                       &
     362                monotonic_limiter_z,                                                               &
     363                timestep_scheme, use_prescribed_profile_data, ws_scheme_sca, air_chemistry
    357364
    358365    USE arrays_3d,                                                                                 &
     
    27502757             CALL advec_s_ws( i, j, chem_species(ilsp)%conc, 'kc', chem_species(ilsp)%flux_s_cs,   &
    27512758                              chem_species(ilsp)%diss_s_cs, chem_species(ilsp)%flux_l_cs,          &
    2752                               chem_species(ilsp)%diss_l_cs, i_omp_start, tn )
     2759                              chem_species(ilsp)%diss_l_cs, i_omp_start, tn, monotonic_limiter_z )
    27532760          ELSE
    27542761             CALL advec_s_pw( i, j, chem_species(ilsp)%conc )
     
    27822789               BTEST( wall_flags_0(k,j,i), 0 )                                                     &
    27832790               )
    2784 
    2785           IF ( chem_species(ilsp)%conc_p(k,j,i) < 0.0_wp )  THEN
    2786              chem_species(ilsp)%conc_p(k,j,i) = 0.1_wp * chem_species(ilsp)%conc(k,j,i)    !FKS6
    2787           ENDIF
    27882791       ENDDO
    27892792!
     
    33153318       r_aero_surf = surf_lsm_h%r_a(m)
    33163319       solar_rad   = surf_lsm_h%rad_sw_dir(m) + surf_lsm_h%rad_sw_dif(m)
     3320       
    33173321       lai = surf_lsm_h%lai(m)
    33183322       sai = lai + 1
  • palm/trunk/SOURCE/modules.f90

    r4069 r4079  
    2525! -----------------
    2626! $Id$
     27! + monotonic_limiter_z
     28!
     29! 4069 2019-07-01 14:05:51Z Giersch
    2730! Masked output running index mid has been introduced as a local variable to
    2831! avoid runtime error (Loop variable has been modified) in time_integration
     
    13601363    LOGICAL ::  masking_method = .FALSE.                         !< namelist parameter
    13611364    LOGICAL ::  mg_switch_to_pe0 = .FALSE.                       !< internal multigrid switch for steering the ghost point exchange in case that data has been collected on PE0
     1365    LOGICAL ::  monotonic_limiter_z = .FALSE.                    !< use monotonic flux limiter for vertical scalar advection
    13621366    LOGICAL ::  nesting_offline = .FALSE.                        !< flag controlling offline nesting in COSMO model 
    13631367    LOGICAL ::  neutral = .FALSE.                                !< namelist parameter
  • palm/trunk/SOURCE/nesting_offl_mod.f90

    r4022 r4079  
    2525! -----------------
    2626! $Id$
     27! - Set boundary condition for w at nzt+1 at the lateral boundaries, even
     28!   though these won't enter the numerical solution. However, due to the mass
     29!   conservation these values might some up to very large values which will
     30!   occur in the run-control file
     31! - Bugfix in offline nesting of chemical species
     32! - Do not set Neumann conditions for TKE and passive scalar
     33!
     34! 4022 2019-06-12 11:52:39Z suehring
    2735! Detection of boundary-layer depth in stable boundary layer on basis of
    2836! boundary data improved
     
    353361                                   BTEST( wall_flags_0(k,j,-1), 3 ) )
    354362             ENDDO
     363             w(nzt,j,-1) = w(nzt-1,j,-1)
    355364          ENDDO
    356365
     
    427436                                    BTEST( wall_flags_0(k,j,nxr+1), 3 ) )
    428437             ENDDO
     438             w(nzt,j,nxr+1) = w(nzt-1,j,nxr+1)
    429439          ENDDO
    430440
     
    504514                                  BTEST( wall_flags_0(k,-1,i), 3 ) )
    505515             ENDDO
     516             w(nzt,-1,i) = w(nzt-1,-1,i)
    506517          ENDDO
    507518
     
    580591                                      BTEST( wall_flags_0(k,nyn+1,i), 3 ) )
    581592             ENDDO
     593             w(nzt,nyn+1,i) = w(nzt-1,nyn+1,i)
    582594          ENDDO
    583595
     
    710722                   DO  j = nys, nyn
    711723                      chem_species(n)%conc(nzt+1,j,i) = interpolate_in_time(   &
    712                                               nest_offl%chem_north(0,j,i,n),   &
    713                                               nest_offl%chem_north(1,j,i,n),   &
     724                                              nest_offl%chem_top(0,j,i,n),   &
     725                                              nest_offl%chem_top(1,j,i,n),   &
    714726                                              fac_dt                       )
    715727                   ENDDO
     
    727739          IF (  bc_dirichlet_n )  diss(:,nyn+1,:) = diss(:,nyn,:)
    728740       ENDIF
    729        IF ( .NOT. constant_diffusion )  THEN
    730           IF (  bc_dirichlet_l )  e(:,:,nxl-1) = e(:,:,nxl)
    731           IF (  bc_dirichlet_r )  e(:,:,nxr+1) = e(:,:,nxr)
    732           IF (  bc_dirichlet_s )  e(:,nys-1,:) = e(:,nys,:)
    733           IF (  bc_dirichlet_n )  e(:,nyn+1,:) = e(:,nyn,:)
    734           e(nzt+1,:,:) = e(nzt,:,:)
    735        ENDIF
    736        IF ( passive_scalar )  THEN
    737           IF (  bc_dirichlet_l )  s(:,:,nxl-1) = s(:,:,nxl)
    738           IF (  bc_dirichlet_r )  s(:,:,nxr+1) = s(:,:,nxr)
    739           IF (  bc_dirichlet_s )  s(:,nys-1,:) = s(:,nys,:)
    740           IF (  bc_dirichlet_n )  s(:,nyn+1,:) = s(:,nyn,:)
    741        ENDIF
     741!        IF ( .NOT. constant_diffusion )  THEN
     742!           IF (  bc_dirichlet_l )  e(:,:,nxl-1) = e(:,:,nxl)
     743!           IF (  bc_dirichlet_r )  e(:,:,nxr+1) = e(:,:,nxr)
     744!           IF (  bc_dirichlet_s )  e(:,nys-1,:) = e(:,nys,:)
     745!           IF (  bc_dirichlet_n )  e(:,nyn+1,:) = e(:,nyn,:)
     746!           e(nzt+1,:,:) = e(nzt,:,:)
     747!        ENDIF
     748!        IF ( passive_scalar )  THEN
     749!           IF (  bc_dirichlet_l )  s(:,:,nxl-1) = s(:,:,nxl)
     750!           IF (  bc_dirichlet_r )  s(:,:,nxr+1) = s(:,:,nxr)
     751!           IF (  bc_dirichlet_s )  s(:,nys-1,:) = s(:,nys,:)
     752!           IF (  bc_dirichlet_n )  s(:,nyn+1,:) = s(:,nyn,:)
     753!        ENDIF
    742754
    743755       CALL exchange_horiz( u, nbgp )
     
    755767          ENDDO
    756768       ENDIF
    757        
     769!
     770!--    Set top boundary condition at all horizontal grid points, also at the
     771!--    lateral boundary grid points.
     772       w(nzt+1,:,:) = w(nzt,:,:)       
    758773!
    759774!--    In case of Rayleigh damping, where the profiles u_init, v_init
  • palm/trunk/SOURCE/parin.f90

    r4022 r4079  
    2525! -----------------
    2626! $Id$
     27! +monotonic_limiter_z
     28!
     29! 4022 2019-06-12 11:52:39Z suehring
    2730! Change default top boundary condition for pressure to Neumann in offline
    2831! nesting case
     
    590593             loop_optimization, lsf_exception, masking_method, mg_cycles,      &
    591594             mg_switch_to_pe0_level, mixing_length_1d, momentum_advec,         &
     595             monotonic_limiter_z,                                              &
    592596             netcdf_precision, neutral, ngsrb,                                 &
    593597             nsor, nsor_ini, nudging, nx, ny, nz, ocean_mode, omega,           &
     
    660664             loop_optimization, lsf_exception, masking_method, mg_cycles,      &
    661665             mg_switch_to_pe0_level, mixing_length_1d, momentum_advec,         &
     666             monotonic_limiter_z,                                              &
    662667             netcdf_precision, neutral, ngsrb,                                 &
    663668             nsor, nsor_ini, nudging, nx, ny, nz, ocean_mode, omega,           &
  • palm/trunk/SOURCE/prognostic_equations.f90

    r4048 r4079  
    2525! -----------------
    2626! $Id$
     27! Application of monotonic flux limiter for the vertical scalar advection
     28! up to the topography top (only for the cache-optimized version at the
     29! moment). Please note, at the moment the limiter is only applied for passive
     30! scalars.
     31!
     32! 4048 2019-06-21 21:00:21Z knoop
    2733! Moved tcm_prognostic_equations to module_interface
    2834!
     
    419425               humidity, intermediate_timestep_count,                          &
    420426               intermediate_timestep_count_max, large_scale_forcing,           &
    421                large_scale_subsidence, neutral, nudging,                       &
     427               large_scale_subsidence,                                         &
     428               monotonic_limiter_z,                                            &
     429               neutral, nudging,                                               &
    422430               ocean_mode, passive_scalar, plant_canopy, pt_reference,         &
    423431               scalar_advec, scalar_advec, simulated_time, sloping_surface,    &
     
    978986             THEN
    979987                IF ( ws_scheme_sca )  THEN
     988!
     989!--                For scalar advection apply monotonic flux limiter near
     990!--                topography.
    980991                   CALL advec_s_ws( i, j, s, 's', flux_s_s, &
    981                                 diss_s_s, flux_l_s, diss_l_s, i_omp_start, tn )
     992                                    diss_s_s, flux_l_s, diss_l_s, i_omp_start, &
     993                                    tn, monotonic_limiter_z )
    982994                ELSE
    983995                   CALL advec_s_pw( i, j, s )
  • palm/trunk/SOURCE/salsa_mod.f90

    r4069 r4079  
    2626! -----------------
    2727! $Id$
     28! Application of monotonic flux limiter for the vertical scalar advection
     29! up to the topography top (only for the cache-optimized version at the
     30! moment).
     31!
     32! 4069 2019-07-01 14:05:51Z Giersch
    2833! Masked output running index mid has been introduced as a local variable to
    2934! avoid runtime error (Loop variable has been modified) in time_integration
     
    77737778    IF ( timestep_scheme(1:5) == 'runge' )  THEN
    77747779       IF ( ws_scheme_sca )  THEN
    7775           CALL advec_s_ws( i, j, rs, id, flux_s, diss_s, flux_l, diss_l, i_omp_start, tn )
     7780          CALL advec_s_ws( i, j, rs, id, flux_s, diss_s, flux_l, diss_l, i_omp_start, tn,          &
     7781                           monotonic_limiter_z )
    77767782       ELSE
    77777783          CALL advec_s_pw( i, j, rs )
Note: See TracChangeset for help on using the changeset viewer.