Changeset 318 for palm/trunk

Timestamp:

May 15, 2009 4:25:33 PM (16 years ago)

Author:

letzel

Message:

• Modified timestep criterion in case of simulations with plant canopy (timestep; modules rolled back to revision 305)

Location:

palm/trunk/SOURCE

Files:

• CURRENT_MODIFICATIONS (modified) (1 diff)
• modules.f90 (modified) (2 diffs)
• timestep.f90 (modified) (2 diffs)

palm/trunk/SOURCE/CURRENT_MODIFICATIONS

@@ -1 +1 @@
 New:
 ---
-Additional timestep criterion in case of simulations with plant canopy (modules, timestep)
+Additional timestep criterion in case of simulations with plant canopy (timestep)
 
 Check for illegal entries in section_xy|xz|yz that exceed nz+1|ny+1|nx+1 

palm/trunk/SOURCE/modules.f90

@@ -5 +5 @@
 ! Actual revisions:
 ! -----------------
-! +cdc_max, lad_u_max, lad_v_max, lad_w_max, 
-! cdc_max_ijk, lad_u_max_ijk, lad_v_max_ijk, lad_w_max_ijk
 ! typographical error in dots_unit fixed
 ! +clip_dvrp_*, cluster_size, color_interval, dvrpsize_interval, dvrp_overlap,
@@ -1219 +1217 @@
     CHARACTER (LEN=40) ::  region(0:9)
     INTEGER ::  pr_palm = 80, statistic_regions = 0, var_ts = 100
-    INTEGER ::  cdc_max_ijk(3), lad_u_max_ijk(3), lad_v_max_ijk(3),  &
-                lad_w_max_ijk(3), u_max_ijk(3), v_max_ijk(3), w_max_ijk(3)
+    INTEGER ::  u_max_ijk(3), v_max_ijk(3), w_max_ijk(3)
     LOGICAL ::  flow_statistics_called = .FALSE.
-    REAL ::     cdc_max, lad_u_max, lad_v_max, lad_w_max, u_max, v_max, w_max
+    REAL ::     u_max, v_max, w_max
     REAL, DIMENSION(:), ALLOCATABLE       ::  sums_divnew_l, sums_divold_l
     REAL, DIMENSION(:,:), ALLOCATABLE     ::  sums, sums_wsts_bc_l, ts_value

palm/trunk/SOURCE/timestep.f90

@@ -47 +47 @@
     INTEGER ::  i, j, k
 
-    REAL ::  div, dt_diff, dt_diff_l, dt_plant_canopy,  & 
-             dt_u, dt_v, dt_w, lad_max, percent_change, &
+    REAL ::  div, dt_diff, dt_diff_l, dt_plant_canopy,                 &
+             dt_plant_canopy_l,                                        &
+             dt_plant_canopy_u, dt_plant_canopy_v, dt_plant_canopy_w,  & 
+             dt_u, dt_v, dt_w, lad_max, percent_change,                &
              u_gtrans_l, vabs_max, value, v_gtrans_l
 
@@ -182 +184 @@
 !--    it is not allowed to extract more than the available momentum
        IF ( plant_canopy ) THEN
-!
-!--       Determine the maxima of the leaf drag coefficient and the leaf 
-!--       area density.
-          CALL global_min_max( nzb, nzt+1, nys-1, nyn+1, nxl-1, nxr+1, cdc, &
-                               'max', cdc_max, cdc_max_ijk )
-          CALL global_min_max( nzb, nzt+1, nys-1, nyn+1, nxl-1, nxr+1, lad_u, &
-                               'max', lad_u_max, lad_u_max_ijk )
-          CALL global_min_max( nzb, nzt+1, nys-1, nyn+1, nxl-1, nxr+1, lad_v, &
-                               'max', lad_v_max, lad_v_max_ijk )
-          CALL global_min_max( nzb, nzt+1, nys-1, nyn+1, nxl-1, nxr+1, lad_w, &
-                               'max', lad_w_max, lad_w_max_ijk )
-!
-!--       Determine the potential maximum of the absolute value of the velocity
-          vabs_max = SQRT( u_max**2 + v_max**2 + w_max**2 )
-          lad_max  = MAX( lad_u_max, lad_v_max, lad_w_max )
-
-          dt_plant_canopy = 1.0 / ( cdc_max * lad_max * vabs_max )
+
+          dt_plant_canopy_l = 0.0
+          DO  i = nxl, nxr
+             DO  j = nys, nyn
+                DO k = nzb+1, nzt
+                   dt_plant_canopy_u = cdc(k,j,i) * lad_u(k,j,i) *  &
+                                       SQRT(     u(k,j,i)**2     +  &
+                                             ( ( v(k,j,i-1)      +  &
+                                                 v(k,j,i)        +  &
+                                                 v(k,j+1,i)      +  &
+                                                 v(k,j+1,i-1) )     &
+                                               / 4.0 )**2        +  &
+                                             ( ( w(k-1,j,i-1)    +  &
+                                                 w(k-1,j,i)      +  &
+                                                 w(k,j,i-1)      +  &
+                                                 w(k,j,i) )         &
+                                                 / 4.0 )**2 ) 
+                   IF ( dt_plant_canopy_u > dt_plant_canopy_l ) THEN
+                      dt_plant_canopy_l = dt_plant_canopy_u  
+                   ENDIF
+                   dt_plant_canopy_v = cdc(k,j,i) * lad_v(k,j,i) *  &
+                                       SQRT( ( ( u(k,j-1,i)      +  &
+                                                 u(k,j-1,i+1)    +  &
+                                                 u(k,j,i)        +  &
+                                                 u(k,j,i+1) )       &
+                                               / 4.0 )**2        +  &
+                                                 v(k,j,i)**2     +  &
+                                             ( ( w(k-1,j-1,i)    +  &
+                                                 w(k-1,j,i)      +  &
+                                                 w(k,j-1,i)      +  &
+                                                 w(k,j,i) )         &
+                                                 / 4.0 )**2 ) 
+                   IF ( dt_plant_canopy_v > dt_plant_canopy_l ) THEN
+                      dt_plant_canopy_l = dt_plant_canopy_v
+                   ENDIF                   
+                   dt_plant_canopy_w = cdc(k,j,i) * lad_w(k,j,i) *  &
+                                       SQRT( ( ( u(k,j,i)        +  &
+                                                 u(k,j,i+1)      +  &
+                                                 u(k+1,j,i)      +  &
+                                                 u(k+1,j,i+1) )     &
+                                               / 4.0 )**2        +  &
+                                             ( ( v(k,j,i)        +  &
+                                                 v(k,j+1,i)      +  &
+                                                 v(k+1,j,i)      +  &
+                                                 v(k+1,j+1,i) )     &
+                                               / 4.0 )**2        +  &
+                                                 w(k,j,i)**2 )     
+                   IF ( dt_plant_canopy_w > dt_plant_canopy_l ) THEN
+                      dt_plant_canopy_l = dt_plant_canopy_w
+                   ENDIF
+                ENDDO
+             ENDDO
+          ENDDO 
+
+          IF ( dt_plant_canopy_l > 0.0 ) THEN
+             dt_plant_canopy_l = 0.1 / dt_plant_canopy_l
+          ENDIF
+!
+!--       Determine the global minumum
+#if defined( __parallel )
+          CALL MPI_ALLREDUCE( dt_plant_canopy_l, dt_plant_canopy, 1, MPI_REAL,  &
+                              MPI_MIN, comm2d, ierr )
+#else
+          dt_plant_canopy = dt_plant_canopy_l
+#endif
 
        ELSE