Index: palm/trunk/SOURCE/advec_ws.f90
===================================================================
--- palm/trunk/SOURCE/advec_ws.f90	(revision 672)
+++ palm/trunk/SOURCE/advec_ws.f90	(revision 673)
@@ -5,5 +5,8 @@
 ! Current revisions:
 ! -----------------
-! 
+! Allocation of weight_substep moved to init_3d_model. 
+! Declaration of ws_scheme_sca and ws_scheme_mom moved to check_parameters.
+! Setting bc for the horizontal velocity variances added (moved from 
+! flow_statistics)
 !
 ! Former revisions:
@@ -81,8 +84,4 @@
        USE statistics
 
-!       
-!--    Set the LOGICALS to enhance the performance.
-       IF ( momentum_advec == 'ws-scheme' )  ws_scheme_mom = .TRUE.
-       IF ( scalar_advec == 'ws-scheme' )  ws_scheme_sca = .TRUE.
 !
 !--    Allocate arrays needed for dissipation control.
@@ -101,15 +100,5 @@
        adv_mom_5 = 0.0083333333333333
        adv_mom_3 = 0.041666666666666
-!
-!--    Allocate arrays needed for statistical evaluation of fluxes when 
-!--    ws-scheme is used.
-!--    The following array contains the weighting factors of the RK3 average
-!--    of tendecies. 
-       ALLOCATE (weight_substep(1:intermediate_timestep_count_max) )
-       IF ( intermediate_timestep_count_max == 3)  THEN   !RK3
-          weight_substep(1) = 0.166666666666666
-          weight_substep(2) = 0.3
-          weight_substep(3) = 0.533333333333333
-       ENDIF
+
 !         
 !--    Arrays needed for statical evaluation of fluxes.
@@ -197,4 +186,6 @@
        
     END SUBROUTINE ws_init
+    
+
     
     SUBROUTINE ws_statistics
@@ -821,5 +812,8 @@
              / (abs(u_comp(k) - gu) + 1.0E-20)         )                      &
              * weight_substep(intermediate_timestep_count) * rmask(j,i,:)
-       ENDDO    
+       ENDDO 
+       sums_us2_ws_l(nzb_u_inner(j,i),:) =                                   &
+                                           sums_us2_ws_l(nzb_u_inner(j,i)+1,:)
+
 !
 !--    Vertical advection, degradation of order near surface and top. 
@@ -1153,4 +1147,6 @@
 
        ENDDO
+       sums_vs2_ws_l(nzb_v_inner(j,i),:) =                                   &  
+                                           sums_vs2_ws_l(nzb_v_inner(j,i)+1,:)
 !
 !--    Vertical advection, degradation of order near surface and top. 
@@ -2208,4 +2204,6 @@
                  * weight_substep(intermediate_timestep_count) * rmask(j,i,:)
             ENDDO
+            sums_us2_ws_l(nzb_u_inner(j,i),:) =                               &
+                                           sums_us2_ws_l(nzb_u_inner(j,i)+1,:)
          ENDDO
        ENDDO
@@ -2570,4 +2568,6 @@
                   * weight_substep(intermediate_timestep_count) * rmask(j,i,:)
             ENDDO
+            sums_vs2_ws_l(nzb_v_inner(j,i),:) =                              &  
+                                           sums_vs2_ws_l(nzb_v_inner(j,i)+1,:)
          ENDDO
        ENDDO
Index: palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- palm/trunk/SOURCE/check_parameters.f90	(revision 672)
+++ palm/trunk/SOURCE/check_parameters.f90	(revision 673)
@@ -4,5 +4,6 @@
 ! Current revisions:
 ! -----------------
-!
+! Declaration of ws_scheme_sca and ws_scheme_mom added (moved from 
+! advec_ws.f90).
 !
 ! Former revisions:
@@ -570,4 +571,9 @@
 !
 !-- Advection schemes:
+!       
+!-- Set the LOGICALS to enhance the performance.
+    IF ( momentum_advec == 'ws-scheme' )    ws_scheme_mom = .TRUE.
+    IF ( scalar_advec   == 'ws-scheme'   )  ws_scheme_sca = .TRUE.
+    
     IF ( momentum_advec /= 'pw-scheme' .AND. momentum_advec /= 'ws-scheme' .AND. &
          momentum_advec /= 'ups-scheme' ) THEN
Index: palm/trunk/SOURCE/data_output_2d.f90
===================================================================
--- palm/trunk/SOURCE/data_output_2d.f90	(revision 672)
+++ palm/trunk/SOURCE/data_output_2d.f90	(revision 673)
@@ -4,5 +4,6 @@
 ! Current revisions:
 ! -----------------
-!
+! When using Multigrid or SOR solver an additional CALL exchange_horiz is 
+! is needed for pressure output.
 ! Former revisions:
 ! -----------------
@@ -239,4 +240,6 @@
 
              CASE ( 'p_xy', 'p_xz', 'p_yz' )
+                IF ( psolver == 'multigrid' .OR. psolver == 'sor' )           &
+                   CALL exchange_horiz( p, nbgp )
                 IF ( av == 0 )  THEN
                    to_be_resorted => p
Index: palm/trunk/SOURCE/data_output_3d.f90
===================================================================
--- palm/trunk/SOURCE/data_output_3d.f90	(revision 672)
+++ palm/trunk/SOURCE/data_output_3d.f90	(revision 673)
@@ -4,5 +4,6 @@
 ! Current revisions:
 ! -----------------
-!
+! When using Multigrid or SOR solver an additional CALL exchange_horiz is 
+! is needed for pressure output.
 ! Former revisions:
 ! -----------------
@@ -152,4 +153,6 @@
 
           CASE ( 'p' )
+             IF ( psolver == 'multigrid' .OR. psolver == 'sor' )              &
+                CALL exchange_horiz( p, nbgp )
              IF ( av == 0 )  THEN
                 to_be_resorted => p
Index: palm/trunk/SOURCE/flow_statistics.f90
===================================================================
--- palm/trunk/SOURCE/flow_statistics.f90	(revision 672)
+++ palm/trunk/SOURCE/flow_statistics.f90	(revision 673)
@@ -4,14 +4,6 @@
 ! Current revisions:
 ! -----------------
-! When advection is computed with ws-scheme, turbulent fluxes are already
-! computed in the respective advection routines and buffered in arrays
-! sums_xx_ws_l(). This is due to a consistent treatment of statistics with the
-! numerics and to avoid unphysical kinks near the surface. 
-! So some if requests has to be done to dicern between fluxes from ws-scheme
-! other advection schemes. 
-! Furthermore the computation of z_i is only done if the heat flux exceeds a
-! minimum value. This affects only simulations of a neutral boundary layer and
-! is due to reasons of computations in the advection scheme.
-!
+! Top bc for the horizontal velocity variances added for ocean runs.
+! Setting the corresponding bottom bc moved to advec_ws.
 !
 ! Former revisions:
@@ -157,14 +149,10 @@
        IF ( ws_scheme_mom )  THEN
 !       
-!--       Boundary condition for u'u' and v'v', because below the surface no
-!--       computation for these quantities is done.
-          DO  i = nxl, nxr
-             DO  j =  nys, nyn
-                sums_us2_ws_l(nzb_u_inner(j,i),sr) =                          & 
-                    sums_us2_ws_l(nzb_u_inner(j,i)+1,sr)
-                sums_vs2_ws_l(nzb_v_inner(j,i),sr) =                          &  
-                    sums_vs2_ws_l(nzb_v_inner(j,i)+1,sr)
-             ENDDO
-          ENDDO
+!--       According to the Neumann bc for the horizontal velocity components,
+!--       the corresponding fluxes has to satisfiy the same bc.
+          IF ( ocean )  THEN
+             sums_us2_ws_l(nzt+1,sr) = sums_us2_ws_l(nzt,sr)
+             sums_vs2_ws_l(nzt+1,sr) = sums_vs2_ws_l(nzt,sr)    
+          ENDIF
 !          
 !--       Swap the turbulent quantities evaluated in advec_ws.
@@ -184,5 +172,5 @@
        ENDIF
        IF ( ws_scheme_sca )  THEN
-          sums_l(:,17,0) = sums_wspts_ws_l(:,sr)      ! w*pts* from advec_s_ws
+          sums_l(:,17,0) = sums_wspts_ws_l(:,sr)      ! w*pt* from advec_s_ws
           IF ( ocean ) sums_l(:,66,0) = sums_wssas_ws_l(:,sr) ! w*sa*
           IF ( humidity  .OR.  passive_scalar ) sums_l(:,49,0) =              &
Index: palm/trunk/SOURCE/init_3d_model.f90
===================================================================
--- palm/trunk/SOURCE/init_3d_model.f90	(revision 672)
+++ palm/trunk/SOURCE/init_3d_model.f90	(revision 673)
@@ -7,5 +7,7 @@
 ! Current revisions:
 ! -----------------
-!
+! weight_substep (moved from advec_ws) and weight_pres added.
+! Allocate p_sub when using Multigrid or SOR solver.
+! Call of ws_init moved behind the if requests.
 ! Former revisions:
 ! -----------------
@@ -208,4 +210,10 @@
               w_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg),  &
               w_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+!
+!-- Following array is required to buffer the perturbation pressure during 
+!-- Runge-Kutta 3rd order time integration.
+    IF ( psolver == 'multigrid' .OR. psolver == 'sor' )  THEN
+       ALLOCATE( p_sub(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+    ENDIF
 
     IF ( timestep_scheme(1:5) /= 'runge' )  THEN
@@ -425,5 +433,16 @@
 
     ENDIF
-
+    
+!
+!-- Allocate arrays containing the RK coefficient for right evaluation of 
+!-- perturbation pressure and turbulent fluxes. At this point it is needed
+!-- for right pressure correction during initialization. Further below 
+!-- the real values will be set.
+    ALLOCATE (weight_substep(1:intermediate_timestep_count_max),           &
+              weight_pres(1:intermediate_timestep_count_max) )
+    weight_substep = 1.
+    weight_pres = 1.
+    intermediate_timestep_count = 1 ! needed for simulated_time=0
+       
 !
 !-- Initialize model variables
@@ -823,4 +842,8 @@
           IF ( precipitation )  precipitation_amount = 0.0
        ENDIF
+       
+!
+!--    Initialize quantities for special advections schemes
+       CALL init_advec
 
 !
@@ -1304,8 +1327,28 @@
 
 !
-!-- Initialize quantities for special advections schemes
-    CALL init_advec
-    IF ( momentum_advec == 'ws-scheme' .OR.  &
-         scalar_advec == 'ws-scheme' ) CALL ws_init
+!-- Initialize the ws-scheme.    
+    IF ( ws_scheme_sca .OR. ws_scheme_mom )  CALL ws_init        
+
+!
+!-- Setting weighting factors for right evaluation of perturbation pressure 
+!-- and turbulent quantities during the RK substeps.                
+    IF ( TRIM(timestep_scheme) == 'runge-kutta-3' )  THEN    ! RK3       
+       weight_substep(1) = 0.166666666666666
+       weight_substep(2) = 0.3
+       weight_substep(3) = 0.533333333333333
+          
+       weight_pres(1) = 0.333333333333333
+       weight_pres(2) = 0.416666666666666
+       weight_pres(3) = 0.25          
+    ELSEIF ( TRIM(timestep_scheme) == 'runge-kutta-2' )  THEN  ! RK2       
+       weight_substep(1) = 0.5
+       weight_substep(2) = 0.5
+          
+       weight_pres(1) = 0.5
+       weight_pres(2) = 0.5          
+    ELSE                                                  ! Euler and Leapfrog       
+       weight_substep(1) = 1.0      
+       weight_pres(1) = 1.0                   
+    ENDIF
 
 !
Index: palm/trunk/SOURCE/modules.f90
===================================================================
--- palm/trunk/SOURCE/modules.f90	(revision 672)
+++ palm/trunk/SOURCE/modules.f90	(revision 673)
@@ -5,5 +5,7 @@
 ! Current revisions:
 ! -----------------
-!
+! +weight_pres to weight the respective contribution of the Runge-Kutta 
+! substeps. +p_sub to buffer the intermediate contributions for Multigrid and 
+! SOR.
 ! Former revisions:
 ! -----------------
@@ -252,6 +254,6 @@
     REAL, DIMENSION(:,:,:), ALLOCATABLE ::                                     &
           canopy_heat_flux, cdc, d, diss, lad_s, lad_u, lad_v, lad_w, lai,     &
-          l_wall, sec, sls, tend, u_m_l, u_m_n, u_m_r, u_m_s, v_m_l, v_m_n,    &
-          v_m_r, v_m_s, w_m_l, w_m_n, w_m_r, w_m_s
+          l_wall, p_sub, sec, sls, tend, u_m_l, u_m_n, u_m_r, u_m_s, v_m_l,    &
+          v_m_n, v_m_r, v_m_s, w_m_l, w_m_n, w_m_r, w_m_s
 
     REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET ::                             &
@@ -1382,5 +1384,5 @@
     REAL ::     u_max, v_max, w_max
     REAL, DIMENSION(:), ALLOCATABLE       ::  sums_divnew_l, sums_divold_l, &
-                                              weight_substep
+                                              weight_substep, weight_pres
     REAL, DIMENSION(:,:), ALLOCATABLE     ::  sums, sums_wsts_bc_l,        &
                                               sums_wsus_ws_l, sums_wsvs_ws_l,&
Index: palm/trunk/SOURCE/pres.f90
===================================================================
--- palm/trunk/SOURCE/pres.f90	(revision 672)
+++ palm/trunk/SOURCE/pres.f90	(revision 673)
@@ -4,5 +4,6 @@
 ! Current revisions:
 ! -----------------
-!
+! Weighting coefficients added for right computation of the pressure during
+! Runge-Kutta substeps. 
 ! Former revisions:
 ! -----------------
@@ -72,5 +73,5 @@
     INTEGER ::  i, j, k, sr
 
-    REAL    ::  ddt_3d, localsum, threadsum
+    REAL    ::  ddt_3d, localsum, threadsum, d_weight_pres
 
     REAL, DIMENSION(1:2) ::  volume_flow_l, volume_flow_offset
@@ -82,4 +83,5 @@
 
     ddt_3d = 1.0 / dt_3d
+    d_weight_pres = 1. / weight_pres(intermediate_timestep_count)
 
 !
@@ -253,5 +255,6 @@
              d(k,j,i) = ( ( u(k,j,i+1) - u(k,j,i) ) * ddx + &
                           ( v(k,j+1,i) - v(k,j,i) ) * ddy + &
-                          ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d
+                          ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d      &
+                                                                * d_weight_pres 
           ENDDO
 !
@@ -268,5 +271,5 @@
                                        - g * ( pt(k+1,j,i) - sums(k+1,4) ) / &
                                          sums(k+1,4)                         &
-                                       ) * ddzw(k+1) * ddt_3d
+                                       ) * ddzw(k+1) * ddt_3d * d_weight_pres 
           ENDIF
 
@@ -291,5 +294,6 @@
              d(k,j,i) = ( ( u(k,j,i+1) - u(k,j,i) ) * ddx + &
                           ( v(k,j+1,i) - v(k,j,i) ) * ddy + &
-                          ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d
+                          ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d      &
+                                                                * d_weight_pres 
              ENDDO
           ENDDO
@@ -307,5 +311,6 @@
                            - g * ( pt(k+1,j,i) - sums(k+1,4) ) / &
                              sums(k+1,4)                         &
-                                        ) * ddzw(k+1) * ddt_3d
+                                        ) * ddzw(k+1) * ddt_3d   &
+                                          * d_weight_pres 
           ENDDO
        ENDDO
@@ -321,5 +326,6 @@
                 d(k,j,i) = ( ( u(k,j,i+1) - u(k,j,i) ) * ddx + &
                           ( v(k,j+1,i) - v(k,j,i) ) * ddy + &
-                          ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d
+                          ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d      &
+                                                                * d_weight_pres 
              ENDDO
           ENDDO
@@ -340,5 +346,6 @@
        ENDDO
     ENDDO
-    localsum = ( localsum + threadsum ) * dt_3d
+    localsum = ( localsum + threadsum ) * dt_3d                               &
+               * weight_pres(intermediate_timestep_count)
     !$OMP END PARALLEL
 #endif
@@ -496,11 +503,5 @@
           DEALLOCATE( p )
           ALLOCATE( p(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
-          DO  i = nxl, nxr
-             DO  j = nys, nyn
-                DO  k = nzb_s_inner(j,i), nzt
-                   p(k,j,i) = tend(k,j,i)
-                ENDDO
-             ENDDO
-          ENDDO
+          
        ENDIF
 
@@ -510,15 +511,37 @@
 !-- Store perturbation pressure on array p, used in the momentum equations
     IF ( psolver(1:7) == 'poisfft' )  THEN
-!
-!--    Here, only the values from the left and right boundaries are copied
-!--    The remaining values are copied in the following loop due to speed
-!--    optimization
-       !$OMP PARALLEL DO
-       DO  j = nysg, nyng
-          DO  k = nzb, nzt+1
-             p(k,j,nxlg:nxl-1) = tend(k,j,nxlg:nxl-1)
-             p(k,j,nxr+1:nxrg) = tend(k,j,nxr+1:nxrg)
-          ENDDO
-       ENDDO
+
+       IF ( intermediate_timestep_count == 1 )  THEN
+          !$OMP PARALLEL PRIVATE (i,j,k)
+          !$OMP DO
+          DO  i = nxlg, nxrg
+             DO  j = nysg, nyng
+                DO  k = nzb, nzt+1
+                   p(k,j,i) = tend(k,j,i)                                     &
+                            * weight_substep(intermediate_timestep_count)
+                   p(k,j,i) = tend(k,j,i)                                     &
+                            * weight_substep(intermediate_timestep_count)
+                ENDDO
+             ENDDO
+          ENDDO
+          !$OMP END PARALLEL
+          
+       ELSE 
+          !$OMP PARALLEL PRIVATE (i,j,k)
+          !$OMP DO
+          DO  i = nxlg, nxrg
+             DO  j = nysg, nyng
+                DO  k = nzb, nzt+1
+                   p(k,j,i) = p(k,j,i) + tend(k,j,i)                          &
+                            * weight_substep(intermediate_timestep_count)
+                   p(k,j,i) = p(k,j,i) + tend(k,j,i)                          &
+                            * weight_substep(intermediate_timestep_count)
+                ENDDO
+             ENDDO
+          ENDDO
+          !$OMP END PARALLEL
+          
+       ENDIF
+       
     ENDIF
 
@@ -533,24 +556,21 @@
     !$OMP PARALLEL PRIVATE (i,j,k)
     !$OMP DO
-    DO  i = nxl, nxr
-       IF ( psolver(1:7) == 'poisfft' )  THEN
-          DO  j = nysg, nyng
-             DO  k = nzb, nzt+1
-                p(k,j,i) = tend(k,j,i)
-             ENDDO
-          ENDDO
-       ENDIF
+    DO  i = nxl, nxr   
        DO  j = nys, nyn
           DO  k = nzb_w_inner(j,i)+1, nzt
-             w(k,j,i) = w(k,j,i) - dt_3d * &
-                                   ( tend(k+1,j,i) - tend(k,j,i) ) * ddzu(k+1)
+             w(k,j,i) = w(k,j,i) - dt_3d *                                    &
+                        ( tend(k+1,j,i) - tend(k,j,i) ) * ddzu(k+1)           &
+                        * weight_pres(intermediate_timestep_count)
           ENDDO
           DO  k = nzb_u_inner(j,i)+1, nzt
-             u(k,j,i) = u(k,j,i) - dt_3d * ( tend(k,j,i) - tend(k,j,i-1) ) * ddx
+             u(k,j,i) = u(k,j,i) - dt_3d *                                 &
+                        ( tend(k,j,i) - tend(k,j,i-1) ) * ddx              &
+                        * weight_pres(intermediate_timestep_count)
           ENDDO
           DO  k = nzb_v_inner(j,i)+1, nzt
-             v(k,j,i) = v(k,j,i) - dt_3d * ( tend(k,j,i) - tend(k,j-1,i) ) * ddy
-          ENDDO
-
+             v(k,j,i) = v(k,j,i) - dt_3d *                                 &
+                        ( tend(k,j,i) - tend(k,j-1,i) ) * ddy              &
+                        * weight_pres(intermediate_timestep_count)
+          ENDDO                                                         
 !
 !--       Sum up the volume flow through the right and north boundary
@@ -575,5 +595,47 @@
     ENDDO
     !$OMP END PARALLEL
-
+    
+    IF ( psolver == 'multigrid' .OR. psolver == 'sor' )  THEN       
+       IF ( intermediate_timestep_count == 1 .OR. simulated_time == 0)  THEN 
+          !$OMP PARALLEL PRIVATE (i,j,k)
+          !$OMP DO
+          DO  i = nxl, nxr
+             DO  j = nys, nyn
+                DO  k = nzb, nzt+1
+                   p_sub(k,j,i) = tend(k,j,i)                                 &
+                                * weight_substep(intermediate_timestep_count)
+                ENDDO
+             ENDDO
+          ENDDO
+          !$OMP END PARALLEL
+       ELSE
+          !$OMP PARALLEL PRIVATE (i,j,k)
+          !$OMP DO
+          DO  i = nxl, nxr
+             DO  j = nys, nyn
+                DO  k = nzb, nzt+1
+                   p_sub(k,j,i) = p_sub(k,j,i) + tend(k,j,i)                  &
+                                * weight_substep(intermediate_timestep_count)
+                ENDDO
+             ENDDO
+          ENDDO
+          !$OMP END PARALLEL
+       ENDIF
+          
+       IF ( intermediate_timestep_count == intermediate_timestep_count_max )  &
+          THEN
+          !$OMP PARALLEL PRIVATE (i,j,k)
+          !$OMP DO
+          DO  i = nxl, nxr
+             DO  j = nys, nyn
+                DO  k = nzb, nzt+1
+                   p(k,j,i) = p_sub(k,j,i) 
+                ENDDO
+             ENDDO
+          ENDDO
+          !$OMP END PARALLEL
+       ENDIF
+    ENDIF
+ 
 !
 !-- Resize tend to its normal size in case of multigrid and ws-scheme.
Index: palm/trunk/SOURCE/prognostic_equations.f90
===================================================================
--- palm/trunk/SOURCE/prognostic_equations.f90	(revision 672)
+++ palm/trunk/SOURCE/prognostic_equations.f90	(revision 673)
@@ -4,5 +4,6 @@
 ! Current revisions:
 ! -----------------
-!
+! Consideration of the pressure gradient (steered by tsc(4)) during the time 
+! integration removed.
 ! Former revisions:
 ! -----------------
@@ -212,5 +213,4 @@
                           dt_3d * (                                            &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) &
-                                 - tsc(4) * ( p(k,j,i) - p(k,j,i-1) ) * ddx    &
                                   ) -                                          &
                            tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
@@ -297,5 +297,4 @@
                           dt_3d * (                                            &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) &
-                                 - tsc(4) * ( p(k,j,i) - p(k,j-1,i) ) * ddy    &
                                   ) -                                          &
                           tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
@@ -383,5 +382,4 @@
                           dt_3d * (                                            &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
-                              - tsc(4) * ( p(k+1,j,i) - p(k,j,i) ) * ddzu(k+1) &
                                   ) -                                          &
                           tsc(5) * rdf(k) * w(k,j,i)
@@ -899,5 +897,5 @@
 !------------------------------------------------------------------------------!
 ! Version with one optimized loop over all equations. It is only allowed to
-! be called for the standard Piascek-Williams advection scheme.
+! be called for the Wicker and Skamarock or Piascek-Williams advection scheme.
 !
 ! Here the calls of most subroutines are embedded in two DO loops over i and j,
@@ -983,5 +981,4 @@
                              dt_3d * (                                         &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) &
-                                 - tsc(4) * ( p(k,j,i) - p(k,j,i-1) ) * ddx    &
                                      ) -                                       &
                              tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
@@ -1050,5 +1047,4 @@
                              dt_3d * (                                         &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) &
-                                 - tsc(4) * ( p(k,j,i) - p(k,j-1,i) ) * ddy    &
                                      ) -                                       &
                              tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
@@ -1116,5 +1112,4 @@
                           dt_3d * (                                         &
                                 tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
-                           - tsc(4) * ( p(k+1,j,i) - p(k,j,i) ) * ddzu(k+1) &
                                   ) -                                       &
                           tsc(5) * rdf(k) * w(k,j,i)
@@ -1525,5 +1520,4 @@
                           dt_3d * (                                            &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) &
-                                 - tsc(4) * ( p(k,j,i) - p(k,j,i-1) ) * ddx    &
                                   ) -                                          &
                           tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
@@ -1617,5 +1611,4 @@
                           dt_3d * (                                            &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) &
-                                 - tsc(4) * ( p(k,j,i) - p(k,j-1,i) ) * ddy    &
                                   ) -                                          &
                           tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
@@ -1706,5 +1699,4 @@
                           dt_3d * (                                            &
                                    tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
-                              - tsc(4) * ( p(k+1,j,i) - p(k,j,i) ) * ddzu(k+1) &
                                   ) -                                          &
                           tsc(5) * rdf(k) * w(k,j,i)
Index: palm/trunk/SOURCE/timestep_scheme_steering.f90
===================================================================
--- palm/trunk/SOURCE/timestep_scheme_steering.f90	(revision 672)
+++ palm/trunk/SOURCE/timestep_scheme_steering.f90	(revision 673)
@@ -4,5 +4,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! No pressure term during time integration (tsc(4)=0.0).
 !
 ! Former revisions:
@@ -60,7 +60,7 @@
              IF ( timestep_scheme == 'leapfrog+euler'  .OR. &
                   timestep_scheme == 'euler' .OR. simulated_time == 0.0 )  THEN
-                tsc(1:5) = (/ 1.0, 1.0, 0.0, 1.0, 1.0 /)
+                tsc(1:5) = (/ 1.0, 1.0, 0.0, 0.0, 1.0 /)
              ELSE
-                tsc(1:5) = (/ 0.0, 2.0, 0.0, 1.0, 2.0 /)
+                tsc(1:5) = (/ 0.0, 2.0, 0.0, 0.0, 2.0 /)
              ENDIF
           ELSE
@@ -69,7 +69,7 @@
 !--          user.
              IF ( timestep_scheme == 'euler' )  THEN
-                tsc(1:5) = (/ 1.0, 1.0, 0.0, 1.0, 1.0 /)
+                tsc(1:5) = (/ 1.0, 1.0, 0.0, 0.0, 1.0 /)
              ELSE
-                tsc(1:5) = (/ 0.0, 2.0, 0.0, 1.0, 2.0 /)
+                tsc(1:5) = (/ 0.0, 2.0, 0.0, 0.0, 2.0 /)
              ENDIF
           ENDIF
@@ -84,11 +84,11 @@
           IF ( simulated_time == 0.0 )  THEN
              dt_changed = .TRUE.
-             tsc(1:5) = (/ 1.0, 1.0, 0.0, 1.0, 1.0 /)
+             tsc(1:5) = (/ 1.0, 1.0, 0.0, 0.0, 1.0 /)
           ELSE
              dt_changed = .FALSE.
              IF ( timestep_scheme == 'euler' )  THEN
-                tsc(1:5) = (/ 1.0, 1.0, 0.0, 1.0, 1.0 /)
+                tsc(1:5) = (/ 1.0, 1.0, 0.0, 0.0, 1.0 /)
              ELSE
-                tsc(1:5) = (/ 0.0, 2.0, 0.0, 1.0, 2.0 /)
+                tsc(1:5) = (/ 0.0, 2.0, 0.0, 0.0, 2.0 /)
              ENDIF
           ENDIF
