Ignore:
Timestamp:
May 30, 2017 5:47:52 PM (4 years ago)
Author:
suehring
Message:

Adjustments according new topography and surface-modelling concept implemented

File:
1 edited

Legend:

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

    r2214 r2232  
    2020! Current revisions:
    2121! -----------------
    22 !
     22! Adjustments to new topography concept
    2323!
    2424! Former revisions:
     
    118118 
    119119    USE arrays_3d,                                                             &
    120         ONLY:  dzu, dzw, e, q, s, shf, tend, u, v, w, zu, zw
     120        ONLY:  dzu, dzw, e, q, s, tend, u, v, w, zu, zw
    121121
    122122    USE indices,                                                               &
    123123        ONLY:  nbgp, nxl, nxlg, nxlu, nxr, nxrg, nyn, nyng, nys, nysg, nysv,   &
    124                nz, nzb, nzb_s_inner, nzb_u_inner, nzb_v_inner, nzb_w_inner, nzt
     124               nz, nzb, nzb_max, nzt, wall_flags_0
    125125
    126126    USE kinds
     
    549549       USE control_parameters,                                                 &
    550550           ONLY:  coupling_char, dz, humidity, io_blocks, io_group,            &
    551                   message_string, ocean, passive_scalar 
    552 
    553        USE control_parameters,                                                 &
    554            ONLY:  urban_surface
     551                  message_string, ocean, passive_scalar, urban_surface
     552
     553       USE surface_mod,                                                        &
     554           ONLY: surf_def_h, surf_lsm_h, surf_usm_h
    555555
    556556       IMPLICIT NONE
     
    562562       INTEGER(iwp) ::  j   !< running index
    563563       INTEGER(iwp) ::  k   !< running index
     564       INTEGER(iwp) ::  m   !< running index
    564565
    565566       REAL(wp) ::  int_bpdf        !< vertical integral for lad-profile construction
     
    639640          canopy_height = pch_index * dz
    640641
    641           DO k = nzb, pch_index
     642          DO k = 0, pch_index
    642643             int_bpdf = int_bpdf +                                             &
    643644                      ( ( ( zw(k) / canopy_height )**( alpha_lad-1.0_wp ) ) *  &
     
    649650!
    650651!--       Preliminary lad profile (defined on w-grid)
    651           DO k = nzb, pch_index
     652          DO k = 0, pch_index
    652653             pre_lad(k) =  lai_beta *                                          &
    653654                        ( ( ( zw(k) / canopy_height )**( alpha_lad-1.0_wp ) )  &
     
    662663!--       when calculating the canopy tendencies)
    663664          lad(0) = pre_lad(0)
    664           DO k = nzb+1, pch_index
     665          DO k = 1, pch_index
    665666             lad(k) = 0.5 * ( pre_lad(k-1) + pre_lad(k) )
    666667          ENDDO         
     
    761762          ENDDO
    762763
     764!           
     765!--       In areas with canopy the surface value of the canopy heat
     766!--       flux distribution overrides the surface heat flux (shf)
     767!--       Start with default surface type
     768          DO  m = 1, surf_def_h(0)%ns
     769             k = surf_def_h(0)%k(m)
     770             IF ( cum_lai_hf(0,j,i) /= 0.0_wp )                                &
     771                surf_def_h(0)%shf(m) = cthf * exp( -ext_coef * cum_lai_hf(0,j,i) )
     772          ENDDO
     773!
     774!--       Natural surfaces
     775          DO  m = 1, surf_lsm_h%ns
     776             k = surf_lsm_h%k(m)
     777             IF ( cum_lai_hf(0,j,i) /= 0.0_wp )                                &
     778                surf_lsm_h%shf(m) = cthf * exp( -ext_coef * cum_lai_hf(0,j,i) )
     779          ENDDO
     780!
     781!--       Urban surfaces
     782          DO  m = 1, surf_usm_h%ns
     783             k = surf_usm_h%k(m)
     784             IF ( cum_lai_hf(0,j,i) /= 0.0_wp )                                &
     785                surf_usm_h%shf(m) = cthf * exp( -ext_coef * cum_lai_hf(0,j,i) )
     786          ENDDO
     787!
    763788!
    764789!--       Calculation of the heating rate (K/s) within the different layers of
    765 !--       the plant canopy
     790!--       the plant canopy. Calculation is only necessary in areas covered with
     791!--       canopy.
     792!--       Within the different canopy layers the plant-canopy heating
     793!--       rate (pc_heating_rate) is calculated as the vertical
     794!--       divergence of the canopy heat fluxes at the top and bottom
     795!--       of the respective layer
    766796          DO  i = nxlg, nxrg
    767797             DO  j = nysg, nyng
    768 !
    769 !--             Calculation only necessary in areas covered with canopy
    770                 IF ( cum_lai_hf(0,j,i) /= 0.0_wp )  THEN
    771 !--             
    772 !--                In areas with canopy the surface value of the canopy heat
    773 !--                flux distribution overrides the surface heat flux (shf)
    774                    shf(j,i) = cthf * exp( -ext_coef * cum_lai_hf(0,j,i) )
    775 !
    776 !--                Within the different canopy layers the plant-canopy heating
    777 !--                rate (pc_heating_rate) is calculated as the vertical
    778 !--                divergence of the canopy heat fluxes at the top and bottom
    779 !--                of the respective layer
    780                    DO  k = 1, pch_index
    781                       pc_heating_rate(k,j,i) = cthf *                         &
    782                                                 ( exp(-ext_coef*cum_lai_hf(k,j,i)) -  &
    783                                                   exp(-ext_coef*cum_lai_hf(k-1,j,i)) ) / dzw(k)
    784                    ENDDO
    785                 ENDIF
     798                DO  k = 1, pch_index
     799                   IF ( cum_lai_hf(0,j,i) /= 0.0_wp )  THEN
     800                      pc_heating_rate(k,j,i) = cthf *                             &
     801                                ( exp(-ext_coef*cum_lai_hf(k,j,i)) -              &
     802                                  exp(-ext_coef*cum_lai_hf(k-1,j,i) ) ) / dzw(k)
     803                   ENDIF
     804                ENDDO
    786805             ENDDO
    787806          ENDDO
     
    966985       INTEGER(iwp) ::  j         !< running index
    967986       INTEGER(iwp) ::  k         !< running index
     987       INTEGER(iwp) ::  k_wall    !< vertical index of topography top
    968988       INTEGER(iwp) ::  kk        !< running index for flat lad arrays
    969989
     
    9871007             DO  i = nxlu, nxr
    9881008                DO  j = nys, nyn
    989                    DO  k = nzb_u_inner(j,i)+1, nzb_u_inner(j,i)+pch_index
    990 
    991                       kk = k - nzb_u_inner(j,i)  !- lad arrays are defined flat
     1009!
     1010!--                Determine topography-top index on u-grid
     1011                   k_wall = MAXLOC(                                            &
     1012                          MERGE( 1, 0,                                         &
     1013                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 14 )    &
     1014                               ), DIM = 1                                      &
     1015                                  ) - 1
     1016                   DO  k = k_wall+1, k_wall+pch_index
     1017
     1018                      kk = k - k_wall   !- lad arrays are defined flat
    9921019!
    9931020!--                   In order to create sharp boundaries of the plant canopy,
     
    10481075             DO  i = nxl, nxr
    10491076                DO  j = nysv, nyn
    1050                    DO  k = nzb_v_inner(j,i)+1, nzb_v_inner(j,i)+pch_index
    1051 
    1052                       kk = k - nzb_v_inner(j,i)  !- lad arrays are defined flat
     1077!
     1078!--                Determine topography-top index on v-grid
     1079                   k_wall = MAXLOC(                                            &
     1080                          MERGE( 1, 0,                                         &
     1081                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 16 )    &
     1082                               ), DIM = 1                                      &
     1083                                  ) - 1
     1084                   DO  k = k_wall+1, k_wall+pch_index
     1085
     1086                      kk = k - k_wall   !- lad arrays are defined flat
    10531087!
    10541088!--                   In order to create sharp boundaries of the plant canopy,
     
    11091143             DO  i = nxl, nxr
    11101144                DO  j = nys, nyn
    1111                    DO  k = nzb_w_inner(j,i)+1, nzb_w_inner(j,i)+pch_index-1
    1112 
    1113                       kk = k - nzb_w_inner(j,i)  !- lad arrays are defined flat
     1145!
     1146!--                Determine topography-top index on w-grid
     1147                   k_wall = MAXLOC(                                            &
     1148                          MERGE( 1, 0,                                         &
     1149                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 18 )    &
     1150                               ), DIM = 1                                      &
     1151                                  ) - 1
     1152                   DO  k = k_wall+1, k_wall+pch_index-1
     1153
     1154                      kk = k - k_wall   !- lad arrays are defined flat
    11141155
    11151156                      pre_tend = 0.0_wp
     
    11571198             DO  i = nxl, nxr
    11581199                DO  j = nys, nyn
    1159                    DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1160                       kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1200!
     1201!--                Determine topography-top index on scalar-grid
     1202                   k_wall = MAXLOC(                                            &
     1203                          MERGE( 1, 0,                                         &
     1204                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1205                               ), DIM = 1                                      &
     1206                                  ) - 1
     1207                   DO  k = k_wall+1, k_wall+pch_index
     1208
     1209                      kk = k - k_wall   !- lad arrays are defined flat
    11611210                      tend(k,j,i) = tend(k,j,i) + pc_heating_rate(kk,j,i)
    11621211                   ENDDO
     
    11691218             DO  i = nxl, nxr
    11701219                DO  j = nys, nyn
    1171                    DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1172                       kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1220!
     1221!--                Determine topography-top index on scalar-grid
     1222                   k_wall = MAXLOC(                                            &
     1223                          MERGE( 1, 0,                                         &
     1224                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1225                               ), DIM = 1                                      &
     1226                                  ) - 1
     1227                   DO  k = k_wall+1, k_wall+pch_index
     1228
     1229                      kk = k - k_wall   !- lad arrays are defined flat
    11731230                      tend(k,j,i) = tend(k,j,i) -                              &
    11741231                                       lsec *                                  &
     
    11941251             DO  i = nxl, nxr
    11951252                DO  j = nys, nyn
    1196                    DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1197                       kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1253!
     1254!--                Determine topography-top index on scalar-grid
     1255                   k_wall = MAXLOC(                                            &
     1256                          MERGE( 1, 0,                                         &
     1257                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1258                               ), DIM = 1                                      &
     1259                                  ) - 1
     1260                   DO  k = k_wall+1, k_wall+pch_index
     1261
     1262                      kk = k - k_wall   !- lad arrays are defined flat
    11981263                      tend(k,j,i) = tend(k,j,i) -                              &
    11991264                                       2.0_wp * cdc *                          &
     
    12181283             DO  i = nxl, nxr
    12191284                DO  j = nys, nyn
    1220                    DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1221                       kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1285!
     1286!--                Determine topography-top index on scalar-grid
     1287                   k_wall = MAXLOC(                                            &
     1288                          MERGE( 1, 0,                                         &
     1289                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1290                               ), DIM = 1                                      &
     1291                                  ) - 1
     1292                   DO  k = k_wall+1, k_wall+pch_index
     1293
     1294                      kk = k - k_wall   !- lad arrays are defined flat
    12221295                      tend(k,j,i) = tend(k,j,i) -                              &
    12231296                                       lsec *                                  &
     
    12881361       INTEGER(iwp) ::  j         !< running index
    12891362       INTEGER(iwp) ::  k         !< running index
     1363       INTEGER(iwp) ::  k_wall    !< vertical index of topography top
    12901364       INTEGER(iwp) ::  kk        !< running index for flat lad arrays
    12911365
     
    13071381!--       u-component
    13081382          CASE ( 1 )
    1309              DO  k = nzb_u_inner(j,i)+1, nzb_u_inner(j,i)+pch_index
    1310 
    1311                 kk = k - nzb_u_inner(j,i)  !- lad arrays are defined flat
     1383!
     1384!--          Determine topography-top index on u-grid
     1385             k_wall = MAXLOC(                                                  &
     1386                          MERGE( 1, 0,                                         &
     1387                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 14 )    &
     1388                               ), DIM = 1                                      &
     1389                            ) - 1
     1390             DO  k = k_wall+1, k_wall+pch_index
     1391
     1392                kk = k - k_wall  !- lad arrays are defined flat
    13121393!
    13131394!--             In order to create sharp boundaries of the plant canopy,
     
    13631444!--       v-component
    13641445          CASE ( 2 )
    1365              DO  k = nzb_v_inner(j,i)+1, nzb_v_inner(j,i)+pch_index
    1366 
    1367                 kk = k - nzb_v_inner(j,i)  !- lad arrays are defined flat
     1446!
     1447!--          Determine topography-top index on v-grid
     1448             k_wall = MAXLOC(                                                  &
     1449                          MERGE( 1, 0,                                         &
     1450                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 16 )    &
     1451                               ), DIM = 1                                      &
     1452                            ) - 1
     1453             DO  k = k_wall+1, k_wall+pch_index
     1454
     1455                kk = k - k_wall  !- lad arrays are defined flat
    13681456!
    13691457!--             In order to create sharp boundaries of the plant canopy,
     
    14191507!--       w-component
    14201508          CASE ( 3 )
    1421              DO  k = nzb_w_inner(j,i)+1, nzb_w_inner(j,i)+pch_index-1
    1422 
    1423                 kk = k - nzb_w_inner(j,i)  !- lad arrays are defined flat
     1509!
     1510!--          Determine topography-top index on w-grid
     1511             k_wall = MAXLOC(                                                  &
     1512                          MERGE( 1, 0,                                         &
     1513                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 18 )    &
     1514                               ), DIM = 1                                      &
     1515                            ) - 1
     1516             DO  k = k_wall+1, k_wall+pch_index-1
     1517
     1518                kk = k - k_wall  !- lad arrays are defined flat
    14241519
    14251520                pre_tend = 0.0_wp
     
    14621557!--       potential temperature
    14631558          CASE ( 4 )
    1464              DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1465                 kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1559!
     1560!--          Determine topography-top index on scalar grid
     1561             k_wall = MAXLOC(                                                  &
     1562                          MERGE( 1, 0,                                         &
     1563                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1564                               ), DIM = 1                                      &
     1565                            ) - 1
     1566             DO  k = k_wall+1, k_wall+pch_index
     1567                kk = k - k_wall  !- lad arrays are defined flat
    14661568                tend(k,j,i) = tend(k,j,i) + pc_heating_rate(kk,j,i)
    14671569             ENDDO
     
    14711573!--       humidity
    14721574          CASE ( 5 )
    1473              DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1474                 kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1575!
     1576!--          Determine topography-top index on scalar grid
     1577             k_wall = MAXLOC(                                                  &
     1578                          MERGE( 1, 0,                                         &
     1579                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1580                               ), DIM = 1                                      &
     1581                            ) - 1
     1582             DO  k = k_wall+1, k_wall+pch_index
     1583
     1584                kk = k - k_wall
    14751585                tend(k,j,i) = tend(k,j,i) -                                    &
    14761586                                 lsec *                                        &
     
    14921602!--       sgs-tke
    14931603          CASE ( 6 )
    1494              DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1495                 kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1604!
     1605!--          Determine topography-top index on scalar grid
     1606             k_wall = MAXLOC(                                                  &
     1607                          MERGE( 1, 0,                                         &
     1608                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1609                               ), DIM = 1                                      &
     1610                            ) - 1
     1611             DO  k = k_wall+1, k_wall+pch_index
     1612
     1613                kk = k - k_wall
    14961614                tend(k,j,i) = tend(k,j,i) -                                    &
    14971615                                 2.0_wp * cdc *                                &
     
    15131631!--       scalar concentration
    15141632          CASE ( 7 )
    1515              DO  k = nzb_s_inner(j,i)+1, nzb_s_inner(j,i)+pch_index
    1516                 kk = k - nzb_s_inner(j,i)  !- lad arrays are defined flat
     1633!
     1634!--          Determine topography-top index on scalar grid
     1635             k_wall = MAXLOC(                                                  &
     1636                          MERGE( 1, 0,                                         &
     1637                                 BTEST( wall_flags_0(nzb:nzb_max,j,i), 12 )    &
     1638                               ), DIM = 1                                      &
     1639                            ) - 1
     1640             DO  k = k_wall+1, k_wall+pch_index
     1641
     1642                kk = k - k_wall
    15171643                tend(k,j,i) = tend(k,j,i) -                                    &
    15181644                                 lsec *                                        &
Note: See TracChangeset for help on using the changeset viewer.