Changeset 4630 for palm/trunk/SOURCE/urban_surface_mod.f90

Timestamp:

Jul 30, 2020 2:54:34 PM (4 years ago)

Author:

suehring

Message:

Land-surface model: bugfix in level 3 initialization of root-area-density; Avoid double classifiation of vertical walls (at surfaces that are alo covered by buildings); Land/urban surface: bugfix in resistance calculation - avoid potential divisions by zero; init_grid: in case of ASCII topography flag grid points as terrain and building to allow application of land/urban-surface model

File:

• palm/trunk/SOURCE/urban_surface_mod.f90 (modified) (5 diffs)

palm/trunk/SOURCE/urban_surface_mod.f90

@@ -27 +27 @@
 ! -----------------
 ! $Id$
+! - Bugfix in resistance calculation - avoid potential divisions by zero
+! - Minor formatting adjustment
+!
+! 4602 2020-07-14 14:49:45Z suehring
 ! Add missing initialization of albedo type with values given from static input file
 ! 
@@ -463 +467 @@
     INTEGER(iwp) ::  wall_category = 2               !< default category for wall surface over pedestrian zone
 
+    REAL(wp)     ::  d_roughness_concrete            !< inverse roughness length of average concrete surface
     REAL(wp)     ::  roughness_concrete = 0.001_wp   !< roughness length of average concrete surface
 
@@ -3583 +3588 @@
     ENDIF
 !
+!-- Calculate constant values
+    d_roughness_concrete = 1.0_wp / roughness_concrete
+!
 !-- Flag surface elements belonging to the ground floor level. Therefore, use terrain height array
 !-- from file, if available. This flag is later used to control initialization of surface attributes.
@@ -7539 +7547 @@
                   q_s,                     &  !< saturation specific humidity
                   rho_lv,                  &  !< frequently used parameter for green layers
-                  tend                      !< tendency
+                  tend,                    &  !< tendency
+                  ueff                        !< limited near-surface wind speed - used for calculation of resistance
 
 
@@ -8046 +8055 @@
 !--        normal component is obtained by simple linear interpolation. (An alternative would be an
 !--        logarithmic interpolation.) Parameter roughness_concrete (default value = 0.001) is used
-!--        to calculation of roughness relative to concrete
-           surf_usm_v(l)%r_a(m) = rho_cp / ( surf_usm_v(l)%z0(m) / roughness_concrete              &
-                                  * ( 11.8_wp + 4.2_wp                                             &
-                                      * SQRT( MAX( ( ( u(k,j,i) + u(k,j,i+1) ) * 0.5_wp )**2 +     &
-                                          ( ( v(k,j,i) + v(k,j+1,i) ) * 0.5_wp )**2 +              &
-                                          ( ( w(k,j,i) + w(k-1,j,i) ) * 0.5_wp )**2, 0.01_wp ) )   &
-                                    )  - 4.0_wp  )
+!--        to calculation of roughness relative to concrete. Note, wind velocity is limited
+!--        to avoid division by zero. The nominator can become <= 0.0 for values z0 < 3*10E-4.
+           ueff        = MAX ( SQRT( ( ( u(k,j,i) + u(k,j,i+1) ) * 0.5_wp )**2 +                   &
+                                     ( ( v(k,j,i) + v(k,j+1,i) ) * 0.5_wp )**2 +                   &
+                                     ( ( w(k,j,i) + w(k-1,j,i) ) * 0.5_wp )**2 ),                  &
+                               1.0_wp / 4.2_wp                                                     &
+                               * ( 4.0_wp / ( surf_usm_v(l)%z0(m) * d_roughness_concrete )         &
+                               - 11.8_wp ),                                                        &
+                               0.1_wp                                                              &
+                              )
+
+           surf_usm_v(l)%r_a(m) = rho_cp / ( surf_usm_v(l)%z0(m) * d_roughness_concrete            &
+                                  * ( 11.8_wp + 4.2_wp * ueff )  - 4.0_wp  )
 !
 !--        Limit aerodynamic resistance