Changeset 1314 for palm/trunk/SOURCE/lpm_init.f90

Timestamp:

Mar 14, 2014 6:25:17 PM (10 years ago)

Author:

suehring

Message:

Vertical logarithmic interpolation of horizontal particle speed for particles
between roughness height and first vertical grid level.

File:

• palm/trunk/SOURCE/lpm_init.f90 (modified) (3 diffs)

palm/trunk/SOURCE/lpm_init.f90

@@ -20 +20 @@
 ! Current revisions:
 ! -----------------
-! 
+! Vertical logarithmic interpolation of horizontal particle speed for particles 
+! between roughness height and first vertical grid level.
 !
 ! Former revisions:
@@ -110 +111 @@
     IMPLICIT NONE
 
-    INTEGER ::  i, j, n, nn
+    INTEGER ::  i, j, k, n, nn
 #if defined( __parallel )
     INTEGER, DIMENSION(3) ::  blocklengths, displacements, types
 #endif
     LOGICAL ::  uniform_particles_l
-    REAL    ::  pos_x, pos_y, pos_z
+    REAL    ::  height_int, height_p, pos_x, pos_y, pos_z, z_p,            &
+                z0_av_local = 0.0
 
 
@@ -184 +186 @@
                  de_dy(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &
                  de_dz(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+    ENDIF
+
+!
+!-- Allocate array required for logarithmic vertical interpolation of 
+!-- horizontal particle velocities between the surface and the first vertical
+!-- grid level. In order to avoid repeated CPU cost-intensive CALLS of 
+!-- intrinsic FORTRAN procedure LOG(z/z0), LOG(z/z0) is precalculated for 
+!-- several heights. Splitting into 20 sublayers turned out to be sufficient. 
+!-- To obtain exact height levels of particles, linear interpolation is applied
+!-- (see lpm_advec.f90).
+    IF ( prandtl_layer )  THEN
+       
+       ALLOCATE ( log_z_z0(0:number_of_sublayers) ) 
+       z_p         = zu(nzb+1) - zw(nzb)
+
+!
+!--    Calculate horizontal mean value of z0 used for logartihmic
+!--    interpolation. Note: this is not exact for heterogeneous z0. 
+!--    However, sensitivity studies showed that the effect is 
+!--    negligible. 
+       z0_av_local  = SUM( z0(nys:nyn,nxl:nxr) )
+       z0_av_global = 0.0
+
+       CALL MPI_ALLREDUCE(z0_av_local, z0_av_global, 1, MPI_REAL, MPI_SUM, &
+                          comm2d, ierr )
+
+       z0_av_global = z0_av_global  / ( ( ny + 1 ) * ( nx + 1 ) )
+!
+!--    Horizontal wind speed is zero below and at z0
+       log_z_z0(0) = 0.0   
+!
+!--    Calculate vertical depth of the sublayers
+       height_int  = ( z_p - z0_av_global ) / REAL( number_of_sublayers ) 
+!
+!--    Precalculate LOG(z/z0)
+       height_p    = 0.0
+       DO  k = 1, number_of_sublayers
+
+          height_p    = height_p + height_int
+          log_z_z0(k) = LOG( height_p / z0_av_global )
+
+       ENDDO
+
+
     ENDIF