Changeset 132 for palm/trunk/SOURCE/init_3d_model.f90

Timestamp:

Nov 20, 2007 9:46:11 AM (16 years ago)

Author:

letzel

Message:

Vertical profiles now based on nzb_s_inner; they are divided by
ngp_2dh_s_inner (scalars, procucts of scalars and velocity components) and
ngp_2dh (staggered velocity components and their products), respectively.

Allow new case bc_uv_t = 'dirichlet_0' for channel flow.

File:

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

palm/trunk/SOURCE/init_3d_model.f90

@@ -7 +7 @@
 ! Actual revisions:
 ! -----------------
-! 
+! New counter ngp_2dh_s_inner.
+! Allow new case bc_uv_t = 'dirichlet_0' for channel flow.
+! Corrected calculation of initial volume flow for 'set_1d-model_profiles' and
+! 'set_constant_profiles' in case of buildings in the reference cross-sections.
 !
 ! Former revisions:
@@ -80 +83 @@
     INTEGER, DIMENSION(:), ALLOCATABLE ::  ngp_2dh_l, ngp_3d_inner_l
 
-    INTEGER, DIMENSION(:,:), ALLOCATABLE ::  ngp_2dh_outer_l
+    INTEGER, DIMENSION(:,:), ALLOCATABLE ::  ngp_2dh_outer_l,  &
+         ngp_2dh_s_inner_l
 
     REAL, DIMENSION(1:2) ::  volume_flow_area_l, volume_flow_initial_l
@@ -97 +101 @@
               ngp_2dh_outer(nzb:nzt+1,0:statistic_regions),                 &
               ngp_2dh_outer_l(nzb:nzt+1,0:statistic_regions),               &
+              ngp_2dh_s_inner(nzb:nzt+1,0:statistic_regions),               &
+              ngp_2dh_s_inner_l(nzb:nzt+1,0:statistic_regions),             &
               rmask(nys-1:nyn+1,nxl-1:nxr+1,0:statistic_regions),           &
               sums(nzb:nzt+1,pr_palm+max_pr_user),                          &
@@ -426 +432 @@
                 ENDDO
              ENDDO
+             IF ( conserve_volume_flow )  THEN
+                IF ( nxr == nx )  THEN
+                   DO  j = nys, nyn
+                      DO  k = nzb + 1, nzb_u_inner(j,nx)
+                         u_nzb_p1_for_vfc(j) = u1d(k) * dzu(k)
+                      ENDDO
+                   ENDDO
+                ENDIF
+                IF ( nyn == ny )  THEN
+                   DO  i = nxl, nxr
+                      DO  k = nzb + 1, nzb_v_inner(ny,i)
+                         v_nzb_p1_for_vfc(i) = v1d(k) * dzu(k)
+                      ENDDO
+                   ENDDO
+                ENDIF
+             ENDIF
 !
 !--          WARNING: The extra boundary conditions set after running the
@@ -486 +508 @@
              IF ( nxr == nx )  THEN
                 DO  j = nys, nyn
-                   k = nzb_u_inner(j,nx) + 1
-                   u_nzb_p1_for_vfc(j) = u_init(k) * dzu(k)
+                   DO  k = nzb + 1, nzb_u_inner(j,nx) + 1
+                      u_nzb_p1_for_vfc(j) = u_init(k) * dzu(k)
+                   ENDDO
                 ENDDO
              ENDIF
              IF ( nyn == ny )  THEN
                 DO  i = nxl, nxr
-                   k = nzb_v_inner(ny,i) + 1
-                   v_nzb_p1_for_vfc(i) = v_init(k) * dzu(k)
+                   DO  k = nzb + 1, nzb_v_inner(ny,i) + 1
+                      v_nzb_p1_for_vfc(i) = v_init(k) * dzu(k)
+                   ENDDO
                 ENDDO
              ENDIF
@@ -561 +585 @@
 
 !
+!--    apply channel flow boundary condition
+       IF ( TRIM( bc_uv_t ) == 'dirichlet_0' )  THEN
+
+          u(nzt+1,:,:) = 0.0
+          v(nzt+1,:,:) = 0.0
+
+!--       for the Dirichlet condition to be correctly applied at the top, set
+!--       ug and vg to zero there
+          ug(nzt+1)    = 0.0
+          vg(nzt+1)    = 0.0
+
+       ENDIF
+
+!
 !--    Calculate virtual potential temperature
        IF ( humidity ) vpt = pt * ( 1.0 + 0.61 * q )
@@ -1036 +1074 @@
 !--          total domain
 !-- ngp_3d:  number of grid points of the total domain
-!-- Note: The lower vertical index nzb_s_outer imposes a small error on the 2D 
-!-- ----  averages of staggered variables such as u and v due to the topography
-!--       arrangement on the staggered grid. Maybe revise later.
-    ngp_2dh_outer_l = 0
-    ngp_2dh_outer   = 0
-    ngp_2dh_l       = 0
-    ngp_2dh         = 0
-    ngp_3d_inner_l  = 0
-    ngp_3d_inner    = 0
-    ngp_3d          = 0
-    ngp_sums        = ( nz + 2 ) * ( pr_palm + max_pr_user )
+    ngp_2dh_outer_l   = 0
+    ngp_2dh_outer     = 0
+    ngp_2dh_s_inner_l = 0
+    ngp_2dh_s_inner   = 0
+    ngp_2dh_l         = 0
+    ngp_2dh           = 0
+    ngp_3d_inner_l    = 0
+    ngp_3d_inner      = 0
+    ngp_3d            = 0
+    ngp_sums          = ( nz + 2 ) * ( pr_palm + max_pr_user )
 
     DO  sr = 0, statistic_regions
@@ -1059 +1096 @@
                 DO  k = nzb_s_outer(j,i), nz + 1
                    ngp_2dh_outer_l(k,sr) = ngp_2dh_outer_l(k,sr) + 1
+                ENDDO
+                DO  k = nzb_s_inner(j,i), nz + 1
+                   ngp_2dh_s_inner_l(k,sr) = ngp_2dh_s_inner_l(k,sr) + 1
                 ENDDO
 !
@@ -1075 +1115 @@
     CALL MPI_ALLREDUCE( ngp_2dh_outer_l(0,0), ngp_2dh_outer(0,0), (nz+2)*sr, &
                         MPI_INTEGER, MPI_SUM, comm2d, ierr )
+    CALL MPI_ALLREDUCE( ngp_2dh_s_inner_l(0,0), ngp_2dh_s_inner(0,0),        &
+                        (nz+2)*sr, MPI_INTEGER, MPI_SUM, comm2d, ierr )
     CALL MPI_ALLREDUCE( ngp_3d_inner_l(0), ngp_3d_inner(0), sr, MPI_INTEGER, &
                         MPI_SUM, comm2d, ierr )
 #else
-    ngp_2dh       = ngp_2dh_l
-    ngp_2dh_outer = ngp_2dh_outer_l
-    ngp_3d_inner  = ngp_3d_inner_l
+    ngp_2dh         = ngp_2dh_l
+    ngp_2dh_outer   = ngp_2dh_outer_l
+    ngp_2dh_s_inner = ngp_2dh_s_inner_l
+    ngp_3d_inner    = ngp_3d_inner_l
 #endif