Changeset 4159 for palm/trunk/SOURCE/indoor_model_mod.f90

Timestamp:

Aug 15, 2019 1:31:35 PM (5 years ago)

Author:

suehring

Message:

Revision of topography processing to have a consistent treatment of 2D and 3D buildings (init_grid, surface_mod); Bugfix in indoor model in case of non grid-resolved buildings

File:

• palm/trunk/SOURCE/indoor_model_mod.f90 (modified) (9 diffs)

palm/trunk/SOURCE/indoor_model_mod.f90

@@ -26 +26 @@
 ! -----------------
 ! $Id: indoor_model_mod.f90 
+! Bugfix in case of non grid-resolved buildings. Further, vertical grid spacing
+! is now considered at the correct level.  
+! 
+! 
 ! - change calculation of a_m and c_m
 ! - change calculation of u-values (use h_es in building array)
@@ -449 +453 @@
     REAL(wp) ::  du_tmp                                    !< 1/u_tmp
     REAL(wp) ::  du_win_tmp                                !< 1/building(nb)%u_value_win
+    REAL(wp) ::  facade_area_v                             !< dummy to compute the total facade area from vertical walls
 
     REAL(wp), DIMENSION(:), ALLOCATABLE ::  volume         !< total building volume at each discrete height level
@@ -588 +593 @@
              nb = MINLOC( ABS( buildings(:)%id - building_id_f%var(j,i) ),    &
                          DIM = 1 )
-             DO  k = nzb+1, nzt+1
+             DO  k = nzb, nzt+1
 !
 !--             Check if grid point belongs to a building. 
@@ -622 +627 @@
        DO  k = buildings(nb)%kb_min, buildings(nb)%kb_max
           buildings(nb)%building_height = buildings(nb)%building_height        &
-                                        + dzw(k)
+                                        + dzw(k+1)
        ENDDO
     ENDDO
@@ -650 +655 @@
                    DO  k = buildings(nb)%kb_min, buildings(nb)%kb_max
                       IF ( building_id_f%var(j,i) /= building_id_f%fill )      &
-                         volume_l(k) = volume_l(k) + dx * dy * dzw(k)
+                         volume_l(k) = volume_l(k) + dx * dy * dzw(k+1)
                    ENDDO
                 ENDDO
@@ -720 +725 @@
 !--    Determine building index and check whether building is on PE
        nb = MINLOC( ABS( buildings(:)%id - building_id_f%var(j,i) ), DIM = 1 )
+
        IF ( buildings(nb)%on_pe )  THEN
 !
@@ -846 +852 @@
        IF ( buildings(nb)%on_pe )  THEN
          
+          facade_area_v = 0.0_wp          
+          DO  k = buildings(nb)%kb_min, buildings(nb)%kb_max
+             facade_area_v = facade_area_v + buildings(nb)%num_facade_v(k)     &
+                             * dzw(k+1) * dx
+          ENDDO
+         
           buildings(nb)%vpf = buildings(nb)%vol_tot /                          &
-               ( buildings(nb)%num_facades_per_building_h * dx * dy +          &
-                 SUM( buildings(nb)%num_facade_v(buildings(nb)%kb_min:         &
-                                                 buildings(nb)%kb_max)         &
-                    * dzw(buildings(nb)%kb_min:buildings(nb)%kb_max) ) * dx )
+                        ( buildings(nb)%num_facades_per_building_h * dx * dy + &
+                          facade_area_v )
        ENDIF
     ENDDO
@@ -1130 +1140 @@
 !--          Building geometries --> not time-dependent
              facade_element_area          = dx * dy                                                   !< [m2] surface area per facade element    
-             floor_area_per_facade        = buildings(nb)%vpf(kk) * ddzw(kk)                          !< [m2/m2] floor area per facade area
+             floor_area_per_facade        = buildings(nb)%vpf(kk) * ddzw(kk+1)                        !< [m2/m2] floor area per facade area
              indoor_volume_per_facade     = buildings(nb)%vpf(kk)                                     !< [m3/m2] indoor air volume per facade area 
              buildings(nb)%area_facade    = facade_element_area *                                   &
@@ -1337 +1347 @@
 !--          EXCEPT facade_element_area, EVERYTHING IS CALCULATED EQUALLY)
 !--          Building geometries  --> not time-dependent
-             IF ( l == 0  .OR. l == 1 ) facade_element_area = dx * dzw(kk)                            !< [m2] surface area per facade element
-             IF ( l == 2  .OR. l == 3 ) facade_element_area = dy * dzw(kk)                            !< [m2] surface area per facade element
-             floor_area_per_facade        = buildings(nb)%vpf(kk) * ddzw(kk)                          !< [m2/m2] floor area per facade area
+             IF ( l == 0  .OR. l == 1 ) facade_element_area = dx * dzw(kk+1)                          !< [m2] surface area per facade element
+             IF ( l == 2  .OR. l == 3 ) facade_element_area = dy * dzw(kk+1)                          !< [m2] surface area per facade element
+             floor_area_per_facade        = buildings(nb)%vpf(kk) * ddzw(kk+1)                        !< [m2/m2] floor area per facade area
              indoor_volume_per_facade     = buildings(nb)%vpf(kk)                                     !< [m3/m2] indoor air volume per facade area 
              buildings(nb)%area_facade    = facade_element_area *                                   &