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

Timestamp:

Nov 28, 2007 10:03:58 AM (16 years ago)

Author:

letzel

Message:

Plant canopy model of Watanabe (2004,BLM 112,307-341) added.

File:

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

palm/trunk/SOURCE/init_3d_model.f90

@@ -235 +235 @@
 
 !
+!-- 3D-arrays for the leaf area density and the canopy drag coefficient
+    IF ( plant_canopy ) THEN
+       ALLOCATE ( lad_s(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1),  &
+                  lad_u(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1),  &
+                  lad_v(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1),  &
+                  lad_w(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1),  &
+                  cdc(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
+    ENDIF
+
+!
 !-- 4D-array for storing the Rif-values at vertical walls
     IF ( topography /= 'flat' )  THEN
@@ -923 +933 @@
                                                     'problem'
        CALL local_stop
+    ENDIF
+
+!
+!-- Initialization of the leaf area density
+    IF ( plant_canopy ) THEN
+ 
+       SELECT CASE ( TRIM( canopy_mode ) )
+
+          CASE( 'block' )
+
+             DO  i = nxl-1, nxr+1
+                DO  j = nys-1, nyn+1
+                   lad_s(:,j,i) = lad(:)
+                   cdc(:,j,i)   = drag_coefficient
+                ENDDO
+             ENDDO
+
+          CASE DEFAULT
+
+!
+!--          The DEFAULT case is reached either if the parameter 
+!--          canopy mode contains a wrong character string or if the 
+!--          user has coded a special case in the user interface. 
+!--          There, the subroutine user_init_plant_canopy checks 
+!--          which of these two conditions applies.
+             CALL user_init_plant_canopy
+ 
+          END SELECT
+
+       CALL exchange_horiz( lad_s )
+       CALL exchange_horiz( cdc )
+
+       DO  i = nxl, nxr
+          DO  j = nys, nyn
+             DO  k = nzb, nzt+1 
+                lad_u(k,j,i) = 0.5 * ( lad_s(k,j,i-1) + lad_s(k,j,i) )
+                lad_v(k,j,i) = 0.5 * ( lad_s(k,j-1,i) + lad_s(k,j,i) )
+             ENDDO
+             DO  k = nzb, nzt
+                lad_w(k,j,i) = 0.5 * ( lad_s(k+1,j,i) + lad_s(k,j,i) )
+             ENDDO
+          ENDDO
+       ENDDO
+
+       lad_w(nzt+1,:,:) = lad_w(nzt,:,:)
+
+       CALL exchange_horiz( lad_u )
+       CALL exchange_horiz( lad_v )
+       CALL exchange_horiz( lad_w )
+  
     ENDIF