Changeset 4314


Ignore:
Timestamp:
Nov 29, 2019 10:29:20 AM (5 years ago)
Author:
suehring
Message:

plant-canopy: Bugfix, plant canopy was still considered at building edges on for the u- and v-component; Relax restriction of LAD on building tops. LAD is only omitted at locations where building grid points emerged artificially by the topography filtering

Location:
palm/trunk/SOURCE
Files:
2 edited

Legend:

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

    r4294 r4314  
    2525! -----------------
    2626! $Id$
     27! Set additional topography flag 4 to mark topography grid points emerged
     28! from the filtering process.
     29!
     30! 4294 2019-11-13 18:34:16Z suehring
    2731! Bugfix, always set bit 5 and 6 of wall_flags, indicating terrain- and
    2832! building surfaces in all  cases, in order to enable terrain-following output
     
    16001604                      num_hole_l     = num_hole_l + 1
    16011605!
    1602 !--                   Clear flag 0 and set special flag ( bit 3) to indicate
     1606!--                   Clear flag 0 and set special flag ( bit 4) to indicate
    16031607!--                   that new topography point is a result of filtering process.
    16041608                      topo_3d(k,j,i) = IBCLR( topo_3d(k,j,i), 0 )
    1605                       topo_3d(k,j,i) = IBSET( topo_3d(k,j,i), 3 )
     1609                      topo_3d(k,j,i) = IBSET( topo_3d(k,j,i), 4 )
    16061610!
    16071611!--                   If filled grid point is occupied by a building, classify
     
    21802184          CALL exchange_horiz_int( topo, nys, nyn, nxl, nxr, nzt, nbgp )
    21812185!
    2182 !--       Set lateral boundary conditions for topography on all ghost layers         
     2186!--       Set lateral boundary conditions for topography on all ghost layers
    21832187          IF ( .NOT. bc_ns_cyc )  THEN
    21842188             IF ( nys == 0  )  THEN
     
    26012605    ENDDO
    26022606!
     2607!-- Set flag 4, indicating new topography grid points due to filtering.
     2608    DO i = nxl, nxr
     2609       DO j = nys, nyn
     2610          DO k = nzb, nzt+1
     2611             IF ( BTEST( topo(k,j,i), 4 ) )                                    &
     2612                wall_flags_0(k,j,i) = IBSET( wall_flags_0(k,j,i), 4 )
     2613          ENDDO
     2614       ENDDO
     2615    ENDDO
     2616!
    26032617!-- Exchange ghost points for wall flags
    26042618    CALL exchange_horiz_int( wall_flags_0, nys, nyn, nxl, nxr, nzt, nbgp )
  • palm/trunk/SOURCE/plant_canopy_model_mod.f90

    r4309 r4314  
    2727! -----------------
    2828! $Id$
     29! - Bugfix, plant canopy was still considered at building edges on for the u-
     30!   and v-component.
     31! - Relax restriction of LAD on building tops. LAD is only omitted at
     32!   locations where building grid points emerged artificially by the
     33!   topography filtering.
     34!
     35! 4309 2019-11-26 18:49:59Z suehring
    2936! Typo
    3037!
     
    114121    USE indices,                                                               &
    115122        ONLY:  nbgp, nxl, nxlg, nxlu, nxr, nxrg, nyn, nyng, nys, nysg, nysv,   &
    116                nz, nzb, nzt, topo_top_ind
     123               nz, nzb, nzt, topo_top_ind, wall_flags_0
    117124
    118125    USE kinds
     
    938945           ONLY: message_string, ocean_mode
    939946
    940        USE indices,                                                            &
    941            ONLY:  wall_flags_0
    942 
    943947       USE netcdf_data_input_mod,                                              &
    944948           ONLY:  leaf_area_density_f
     
    11081112!
    11091113!--                   Check if resolved vegetation is mapped onto buildings.
    1110 !--                   Even if the user defines LAD only for non-building covered
    1111 !--                   grid points, the filter process of topography may cause
    1112 !--                   that LAD overlaps with buildings. This case, an
    1113 !--                   informative message is given and the LAD field is set to
    1114 !--                   zero at these
     1114!--                   In general, this is allowed and also meaningful, e.g.
     1115!--                   when trees carry across roofs. However, due to the
     1116!--                   topography filtering, new building grid points can emerge
     1117!--                   at location where also plant canopy is defined. As a
     1118!--                   result, plant canopy is mapped on top of roofs, with
     1119!--                   siginficant impact on the downstream flow field and the
     1120!--                   nearby surface radiation. In order to avoid that
     1121!--                   plant canopy is mistakenly mapped onto building roofs,
     1122!--                   check for building grid points (bit 6) that emerge from
     1123!--                   the filtering (bit 4) and set LAD to zero at these
     1124!--                   artificially  created building grid points. This case,
     1125!--                   an informative message is given.
    11151126                      IF ( ANY( lad_s(:,j,i) /= 0.0_wp )  .AND.                &
    1116                            ANY( BTEST( wall_flags_0(:,j,i), 6 ) ) )  THEN
     1127                           ANY( BTEST( wall_flags_0(:,j,i), 6 ) ) .AND.        &
     1128                           ANY( BTEST( wall_flags_0(:,j,i), 4 ) ) )  THEN
    11171129                         lad_s(:,j,i) = 0.0_wp
    11181130                         WRITE( message_string, * )                            &
    1119                                           'Resolved plant-canopy is ' //       &
    1120                                           'defined on top of a building ' //   &
    1121                                           '- lad is omitted at this grid '   //&
    1122                                           'point: (i,j) = ', i, j
     1131                                        'Resolved plant-canopy is ' //         &
     1132                                        'defined on top of an artificially '// &
     1133                                        'created building grid point ' //      &
     1134                                        '(emerged from the filering) - ' //    &
     1135                                        'LAD profile is omitted at this ' //   &
     1136                                        'grid point: (i,j) = ', i, j
    11231137                         CALL message( 'pcm_init', 'PA0313', 0, 0, myid, 6, 0 )
    11241138                      ENDIF
     
    19181932       INTEGER(iwp) ::  kk        !< running index for flat lad arrays
    19191933
     1934       LOGICAL ::  building_edge_e !< control flag indicating an eastward-facing building edge
     1935       LOGICAL ::  building_edge_n !< control flag indicating a north-facing building edge
     1936       LOGICAL ::  building_edge_s !< control flag indicating a south-facing building edge
     1937       LOGICAL ::  building_edge_w !< control flag indicating a westward-facing building edge
     1938
    19201939       REAL(wp) ::  ddt_3d    !< inverse of the LES timestep (dt_3d)
    19211940       REAL(wp) ::  lad_local !< local lad value
     
    19351954          CASE ( 1 )
    19361955!
     1956!--          Set control flags indicating east- and westward-orientated
     1957!--          building edges. Note, building_egde_w is set from the perspective
     1958!--          of the potential rooftop grid point, while building_edge_e is
     1959!--          is set from the perspective of the non-building grid point.
     1960             building_edge_w = ANY( BTEST( wall_flags_0(:,j,i),   6 ) )  .AND. &
     1961                         .NOT. ANY( BTEST( wall_flags_0(:,j,i-1), 6 ) )
     1962             building_edge_e = ANY( BTEST( wall_flags_0(:,j,i-1), 6 ) )  .AND. &
     1963                         .NOT. ANY( BTEST( wall_flags_0(:,j,i),   6 ) )
     1964!
    19371965!--          Determine topography-top index on u-grid
    19381966             k_wall = topo_top_ind(j,i,1)
     
    19481976!--             than inside of the canopy.
    19491977!--             For the same reason, the lad at the rightmost(i+1)canopy
    1950 !--             boundary on the u-grid equals lad_s(k,j,i).
     1978!--             boundary on the u-grid equals lad_s(k,j,i), which is considered
     1979!--             in the next if-statement. Note, at left-sided building edges
     1980!--             this is not applied, here the LAD is equals the LAD at grid
     1981!--             point (k,j,i), in order to avoid that LAD is mistakenly mapped
     1982!--             on top of a roof where (usually) is no LAD is defined.
    19511983                lad_local = lad_s(kk,j,i)
    1952                 IF ( lad_local == 0.0_wp .AND. lad_s(kk,j,i-1) > 0.0_wp )  THEN
    1953                    lad_local = lad_s(kk,j,i-1)
    1954                 ENDIF
     1984                IF ( lad_local == 0.0_wp .AND. lad_s(kk,j,i-1) > 0.0_wp  .AND. &
     1985                     .NOT. building_edge_w )  lad_local = lad_s(kk,j,i-1)
     1986!
     1987!--             In order to avoid that LAD is mistakenly considered at right-
     1988!--             sided building edges (here the topography-top index for the
     1989!--             u-component at index j,i is still on the building while the
     1990!--             topography top for the scalar isn't), LAD is taken from grid
     1991!--             point (j,i-1).
     1992                IF ( lad_local > 0.0_wp .AND. lad_s(kk,j,i-1) == 0.0_wp  .AND. &
     1993                     building_edge_e )  lad_local = lad_s(kk,j,i-1)
    19551994
    19561995                pre_tend = 0.0_wp
     
    19952034          CASE ( 2 )
    19962035!
     2036!--          Set control flags indicating north- and southward-orientated
     2037!--          building edges. Note, building_egde_s is set from the perspective
     2038!--          of the potential rooftop grid point, while building_edge_n is
     2039!--          is set from the perspective of the non-building grid point.
     2040             building_edge_s = ANY( BTEST( wall_flags_0(:,j,i),   6 ) )  .AND. &
     2041                         .NOT. ANY( BTEST( wall_flags_0(:,j-1,i), 6 ) )
     2042             building_edge_n = ANY( BTEST( wall_flags_0(:,j-1,i), 6 ) )  .AND. &
     2043                         .NOT. ANY( BTEST( wall_flags_0(:,j,i),   6 ) )
     2044!
    19972045!--          Determine topography-top index on v-grid
    19982046             k_wall = topo_top_ind(j,i,2)
    1999 
    20002047             DO  k = k_wall + 1, k_wall + pch_index_ji(j,i)
    20012048
     
    20082055!--             than inside of the canopy.
    20092056!--             For the same reason, the lad at the northmost(j+1)canopy
    2010 !--             boundary on the v-grid equals lad_s(k,j,i).
     2057!--             boundary on the v-grid equals lad_s(k,j,i), which is considered
     2058!--             in the next if-statement. Note, at left-sided building edges
     2059!--             this is not applied, here the LAD is equals the LAD at grid
     2060!--             point (k,j,i), in order to avoid that LAD is mistakenly mapped
     2061!--             on top of a roof where (usually) is no LAD is defined.
    20112062                lad_local = lad_s(kk,j,i)
    2012                 IF ( lad_local == 0.0_wp .AND. lad_s(kk,j-1,i) > 0.0_wp )  THEN
    2013                    lad_local = lad_s(kk,j-1,i)
    2014                 ENDIF
     2063                IF ( lad_local == 0.0_wp .AND. lad_s(kk,j-1,i) > 0.0_wp  .AND. &
     2064                     .NOT. building_edge_s )  lad_local = lad_s(kk,j-1,i)
     2065!
     2066!--             In order to avoid that LAD is mistakenly considered at right-
     2067!--             sided building edges (here the topography-top index for the
     2068!--             u-component at index j,i is still on the building while the
     2069!--             topography top for the scalar isn't), LAD is taken from grid
     2070!--             point (j,i-1).
     2071                IF ( lad_local > 0.0_wp .AND. lad_s(kk,j-1,i) == 0.0_wp  .AND. &
     2072                     building_edge_n )  lad_local = lad_s(kk,j-1,i)
    20152073
    20162074                pre_tend = 0.0_wp
Note: See TracChangeset for help on using the changeset viewer.