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

Timestamp:

Aug 19, 2009 4:56:44 PM (15 years ago)

Author:

letzel

Message:

• Bugfix for extended NetCDF files: In order to avoid 'data mode' errors if updated attributes are larger than their original size, NF90_PUT_ATT is called in 'define mode' enclosed by NF90_REDEF and NF90_ENDDEF calls. This implies a possible performance loss; an alternative strategy would be to ensure equal attribute size in a job chain. (netcdf)
• Changed: For extended NetCDF files, the updated title attribute includes an update of time_average_text where appropriate. (netcdf)
• Bugfix: correction of initial volume flow for non-flat topography (init_3d_model)
• Bugfix: zero initialization of arrays within buildings for 'cyclic_fill' (init_3d_model)
• New: MPI debug option (mrun)
• Documentation updates:
  • topography file format (chapter 3.5.4)
  • NC identifier range expanded to 444 (appendix_a)

File:

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

palm/trunk/SOURCE/init_3d_model.f90

@@ -7 +7 @@
 ! Current revisions:
 ! -----------------
+! bugfix: correction of initial volume flow for non-flat topography
+! bugfix: zero initialization of arrays within buildings for 'cyclic_fill'
 ! bugfix: avoid that ngp_2dh_s_inner becomes zero
 ! initializing_actions='read_data_for_recycling' renamed to 'cyclic_fill', now
@@ -485 +487 @@
                    DO  j = nys, nyn
                       DO  k = nzb + 1, nzb_u_inner(j,nx)
-                         u_nzb_p1_for_vfc(j) = u1d(k) * dzu(k)
+                         u_nzb_p1_for_vfc(j) = u_nzb_p1_for_vfc(j) + &
+                                               u1d(k) * dzu(k)
                       ENDDO
                    ENDDO
@@ -492 +495 @@
                    DO  i = nxl, nxr
                       DO  k = nzb + 1, nzb_v_inner(ny,i)
-                         v_nzb_p1_for_vfc(i) = v1d(k) * dzu(k)
+                         v_nzb_p1_for_vfc(i) = v_nzb_p1_for_vfc(i) + &
+                                               v1d(k) * dzu(k)
                       ENDDO
                    ENDDO
@@ -558 +562 @@
                 DO  j = nys, nyn
                    DO  k = nzb + 1, nzb_u_inner(j,nx) + 1
-                      u_nzb_p1_for_vfc(j) = u_init(k) * dzu(k)
+                      u_nzb_p1_for_vfc(j) = u_nzb_p1_for_vfc(j) + &
+                                            u_init(k) * dzu(k)
                    ENDDO
                 ENDDO
@@ -565 +570 @@
                 DO  i = nxl, nxr
                    DO  k = nzb + 1, nzb_v_inner(ny,i) + 1
-                      v_nzb_p1_for_vfc(i) = v_init(k) * dzu(k)
+                      v_nzb_p1_for_vfc(i) = v_nzb_p1_for_vfc(i) + &
+                                            v_init(k) * dzu(k)
                    ENDDO
                 ENDDO
@@ -1028 +1034 @@
 
 !
+!--    Inside buildings set velocities and TKE back to zero
+       IF ( TRIM( initializing_actions ) == 'cyclic_fill' .AND.  &
+            topography /= 'flat' )  THEN
+!
+!--       Correction of initial volume flow
+          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) = u_nzb_p1_for_vfc(j) + &
+                                            u(k,j,nx) * 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) = v_nzb_p1_for_vfc(i) + &
+                                            v(k,ny,i) * dzu(k)
+                   ENDDO
+                ENDDO
+             ENDIF
+          ENDIF
+
+!
+!--       Inside buildings set velocities and TKE back to zero.
+!--       Other scalars (pt, q, s, km, kh, p, sa, ...) are ignored at present,
+!--       maybe revise later.
+          IF ( timestep_scheme(1:5) == 'runge' )  THEN
+             DO  i = nxl-1, nxr+1
+                DO  j = nys-1, nyn+1
+                   u  (nzb:nzb_u_inner(j,i),j,i) = 0.0
+                   v  (nzb:nzb_v_inner(j,i),j,i) = 0.0
+                   w  (nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   e  (nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   u_m(nzb:nzb_u_inner(j,i),j,i) = 0.0
+                   v_m(nzb:nzb_v_inner(j,i),j,i) = 0.0
+                   w_m(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   e_m(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   tu_m(nzb:nzb_u_inner(j,i),j,i) = 0.0
+                   tv_m(nzb:nzb_v_inner(j,i),j,i) = 0.0
+                   tw_m(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   te_m(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   tpt_m(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                ENDDO
+             ENDDO
+          ELSE
+             DO  i = nxl-1, nxr+1
+                DO  j = nys-1, nyn+1
+                   u  (nzb:nzb_u_inner(j,i),j,i) = 0.0
+                   v  (nzb:nzb_v_inner(j,i),j,i) = 0.0
+                   w  (nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   e  (nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   u_m(nzb:nzb_u_inner(j,i),j,i) = 0.0
+                   v_m(nzb:nzb_v_inner(j,i),j,i) = 0.0
+                   w_m(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   e_m(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   u_p(nzb:nzb_u_inner(j,i),j,i) = 0.0
+                   v_p(nzb:nzb_v_inner(j,i),j,i) = 0.0
+                   w_p(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                   e_p(nzb:nzb_w_inner(j,i),j,i) = 0.0
+                ENDDO
+             ENDDO
+          ENDIF
+
+       ENDIF
+
+!
 !--    Calculate the initial volume flow at the right and north boundary
        IF ( conserve_volume_flow  .AND.  &
@@ -1043 +1117 @@
                 ENDDO
 !
-!--             Correction if velocity at nzb+1 has been set zero further above
-!--             Note: at present, u_nzb_p1_for_vfc is zero (maybe revise later)
+!--             Correction if velocity inside buildings has been set to zero
+!--             further above
                 volume_flow_initial_l(1) = volume_flow_initial_l(1) + &
                                            u_nzb_p1_for_vfc(j)
@@ -1058 +1132 @@
                 ENDDO
 !
-!--             Correction if velocity at nzb+1 has been set zero further above
-!--             Note: at present, v_nzb_p1_for_vfc is zero (maybe revise later)
+!--             Correction if velocity inside buildings has been set to zero
+!--             further above
                 volume_flow_initial_l(2) = volume_flow_initial_l(2) + &
                                            v_nzb_p1_for_vfc(i)