Index: palm/trunk/SOURCE/prognostic_equations.f90
===================================================================
--- palm/trunk/SOURCE/prognostic_equations.f90	(revision 1917)
+++ palm/trunk/SOURCE/prognostic_equations.f90	(revision 1960)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! Separate humidity and passive scalar
 ! 
 ! Former revisions:
@@ -192,12 +192,13 @@
     USE arrays_3d,                                                             &
         ONLY:  diss_l_e, diss_l_nr, diss_l_pt, diss_l_q, diss_l_qr,            &
-               diss_l_sa, diss_s_e, diss_s_nr, diss_s_pt, diss_s_q,            &
-               diss_s_qr, diss_s_sa, e, e_p, flux_s_e, flux_s_nr, flux_s_pt,   &
-               flux_s_q, flux_s_qr, flux_s_sa, flux_l_e, flux_l_nr,            &
-               flux_l_pt, flux_l_q, flux_l_qr, flux_l_sa, nr, nr_p, nrsws,     &
-               nrswst, pt, ptdf_x, ptdf_y, pt_init, pt_p, prho, q, q_init,     &
-               q_p, qsws, qswst, qr, qr_p, qrsws, qrswst, rdf, rdf_sc,         &
-               ref_state, rho, sa, sa_init, sa_p, saswsb, saswst, shf, tend,   &
-               te_m, tnr_m, tpt_m, tq_m, tqr_m, tsa_m, tswst, tu_m, tv_m,      &
+               diss_l_s, diss_l_sa, diss_s_e, diss_s_nr, diss_s_pt, diss_s_q,  &
+               diss_s_qr, diss_s_s, diss_s_sa, e, e_p, flux_s_e, flux_s_nr,    &
+               flux_s_pt, flux_s_q, flux_s_qr, flux_s_s, flux_s_sa, flux_l_e,  &
+               flux_l_nr, flux_l_pt, flux_l_q, flux_l_qr, flux_l_s, flux_l_sa, &
+               nr, nr_p, nrsws, nrswst, pt, ptdf_x, ptdf_y, pt_init, pt_p,     &
+               prho, q, q_init, q_p, qsws, qswst, qr, qr_p, qrsws, qrswst, rdf,&
+               rdf_sc, ref_state, rho, s, s_init, s_p, sa, sa_init, sa_p,      &
+               saswsb, saswst, shf, ssws, sswst, tend,                         &
+               te_m, tnr_m, tpt_m, tq_m, tqr_m, ts_m, tsa_m, tswst, tu_m, tv_m,&
                tw_m, u, ug, u_init, u_p, v, vg, vpt, v_init, v_p, w, w_p
         
@@ -216,5 +217,5 @@
                use_upstream_for_tke, wall_heatflux,                            &
                wall_nrflux, wall_qflux, wall_qflux, wall_qflux, wall_qrflux,   &
-               wall_salinityflux, ws_scheme_mom, ws_scheme_sca
+               wall_salinityflux, wall_sflux, ws_scheme_mom, ws_scheme_sca
 
     USE cpulog,                                                                &
@@ -724,7 +725,6 @@
 
 !
-!--       If required, compute prognostic equation for total water content /
-!--       scalar
-          IF ( humidity  .OR.  passive_scalar )  THEN
+!--       If required, compute prognostic equation for total water content
+          IF ( humidity )  THEN
 
 !
@@ -745,5 +745,5 @@
 
 !
-!--          Sink or source of scalar concentration due to canopy elements
+!--          Sink or source of humidity due to canopy elements
              IF ( plant_canopy )  CALL pcm_tendency( i, j, 5 )
 
@@ -881,5 +881,78 @@
 
           ENDIF
-
+          
+!
+!--       If required, compute prognostic equation for scalar
+          IF ( passive_scalar )  THEN
+!
+!--          Tendency-terms for total water content / scalar
+             tend(:,j,i) = 0.0_wp
+             IF ( timestep_scheme(1:5) == 'runge' ) &
+             THEN
+                IF ( ws_scheme_sca )  THEN
+                   CALL advec_s_ws( i, j, s, 's', flux_s_s, & 
+                                diss_s_s, flux_l_s, diss_l_s, i_omp_start, tn )
+                ELSE 
+                   CALL advec_s_pw( i, j, s )
+                ENDIF
+             ELSE
+                CALL advec_s_up( i, j, s )
+             ENDIF
+             CALL diffusion_s( i, j, s, ssws, sswst, wall_sflux )
+
+!
+!--          Sink or source of scalar concentration due to canopy elements
+             IF ( plant_canopy )  CALL pcm_tendency( i, j, 7 )
+
+!
+!--          Large scale advection, still need to be extended for scalars
+!              IF ( large_scale_forcing )  THEN
+!                 CALL ls_advec( i, j, simulated_time, 's' )
+!              ENDIF
+
+!
+!--          Nudging, still need to be extended for scalars
+!              IF ( nudging )  CALL nudge( i, j, simulated_time, 's' ) 
+
+!
+!--          If required compute influence of large-scale subsidence/ascent.
+!--          Note, the last argument is of no meaning in this case, as it is 
+!--          only used in conjunction with large_scale_forcing, which is to 
+!--          date not implemented for scalars.
+             IF ( large_scale_subsidence  .AND.                                &
+                  .NOT. use_subsidence_tendencies  .AND.                       &
+                  .NOT. large_scale_forcing )  THEN
+                CALL subsidence( i, j, tend, s, s_init, 3 )
+             ENDIF
+
+             CALL user_actions( i, j, 's-tendency' )
+
+!
+!--          Prognostic equation for scalar
+             DO  k = nzb_s_inner(j,i)+1, nzt
+                s_p(k,j,i) = s(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +       &
+                                                  tsc(3) * ts_m(k,j,i) )       &
+                                      - tsc(5) * rdf_sc(k) *                   &
+                                        ( s(k,j,i) - s_init(k) )
+                IF ( s_p(k,j,i) < 0.0_wp )  s_p(k,j,i) = 0.1_wp * s(k,j,i)
+             ENDDO
+
+!
+!--          Calculate tendencies for the next Runge-Kutta step
+             IF ( timestep_scheme(1:5) == 'runge' )  THEN
+                IF ( intermediate_timestep_count == 1 )  THEN
+                   DO  k = nzb_s_inner(j,i)+1, nzt
+                      ts_m(k,j,i) = tend(k,j,i)
+                   ENDDO
+                ELSEIF ( intermediate_timestep_count < &
+                         intermediate_timestep_count_max )  THEN
+                   DO  k = nzb_s_inner(j,i)+1, nzt
+                      ts_m(k,j,i) = -9.5625_wp * tend(k,j,i) + &
+                                     5.3125_wp * ts_m(k,j,i)
+                   ENDDO
+                ENDIF
+             ENDIF
+
+          ENDIF          
 !
 !--       If required, compute prognostic equation for turbulent kinetic 
@@ -1432,8 +1505,8 @@
 
 !
-!-- If required, compute prognostic equation for total water content / scalar
-    IF ( humidity  .OR.  passive_scalar )  THEN
-
-       CALL cpu_log( log_point(29), 'q/s-equation', 'start' )
+!-- If required, compute prognostic equation for total water content
+    IF ( humidity )  THEN
+
+       CALL cpu_log( log_point(29), 'q-equation', 'start' )
 
 !
@@ -1470,5 +1543,5 @@
        
 !
-!--    Sink or source of scalar concentration due to canopy elements
+!--    Sink or source of humidity due to canopy elements
        IF ( plant_canopy ) CALL pcm_tendency( 5 )
 
@@ -1493,5 +1566,5 @@
 
 !
-!--    Prognostic equation for total water content / scalar
+!--    Prognostic equation for total water content
        DO  i = nxl, nxr
           DO  j = nys, nyn
@@ -1529,5 +1602,5 @@
        ENDIF
 
-       CALL cpu_log( log_point(29), 'q/s-equation', 'stop' )
+       CALL cpu_log( log_point(29), 'q-equation', 'stop' )
 
 !
@@ -1686,5 +1759,107 @@
 
     ENDIF
-
+!
+!-- If required, compute prognostic equation for scalar
+    IF ( passive_scalar )  THEN
+
+       CALL cpu_log( log_point(66), 's-equation', 'start' )
+
+!
+!--    Scalar/q-tendency terms with communication
+       sbt = tsc(2)
+       IF ( scalar_advec == 'bc-scheme' )  THEN
+
+          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
+!
+!--          Bott-Chlond scheme always uses Euler time step. Thus:
+             sbt = 1.0_wp
+          ENDIF
+          tend = 0.0_wp
+          CALL advec_s_bc( s, 's' )
+
+       ENDIF
+
+!
+!--    Scalar/q-tendency terms with no communication
+       IF ( scalar_advec /= 'bc-scheme' )  THEN
+          tend = 0.0_wp
+          IF ( timestep_scheme(1:5) == 'runge' )  THEN
+             IF ( ws_scheme_sca )  THEN
+                CALL advec_s_ws( s, 's' )
+             ELSE
+                CALL advec_s_pw( s )
+             ENDIF
+          ELSE
+             CALL advec_s_up( s )
+          ENDIF
+       ENDIF
+
+       CALL diffusion_s( s, ssws, sswst, wall_sflux )
+       
+!
+!--    Sink or source of humidity due to canopy elements
+       IF ( plant_canopy ) CALL pcm_tendency( 7 )
+
+!
+!--    Large scale advection. Not implemented for scalars so far.
+!        IF ( large_scale_forcing )  THEN
+!           CALL ls_advec( simulated_time, 'q' )
+!        ENDIF
+
+!
+!--    Nudging. Not implemented for scalars so far.
+!        IF ( nudging )  CALL nudge( simulated_time, 'q' ) 
+
+!
+!--    If required compute influence of large-scale subsidence/ascent.
+!--    Not implemented for scalars so far.
+       IF ( large_scale_subsidence  .AND.                                      &
+            .NOT. use_subsidence_tendencies  .AND.                             &
+            .NOT. large_scale_forcing )  THEN
+         CALL subsidence( tend, s, s_init, 3 )
+       ENDIF
+
+       CALL user_actions( 's-tendency' )
+
+!
+!--    Prognostic equation for total water content
+       DO  i = nxl, nxr
+          DO  j = nys, nyn
+             DO  k = nzb_s_inner(j,i)+1, nzt
+                s_p(k,j,i) = s(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +          &
+                                                  tsc(3) * ts_m(k,j,i) )       &
+                                      - tsc(5) * rdf_sc(k) *                   &
+                                        ( s(k,j,i) - s_init(k) )
+                IF ( s_p(k,j,i) < 0.0_wp )  s_p(k,j,i) = 0.1_wp * s(k,j,i)
+             ENDDO
+          ENDDO
+       ENDDO
+
+!
+!--    Calculate tendencies for the next Runge-Kutta step
+       IF ( timestep_scheme(1:5) == 'runge' )  THEN
+          IF ( intermediate_timestep_count == 1 )  THEN
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+                   DO  k = nzb_s_inner(j,i)+1, nzt
+                      ts_m(k,j,i) = tend(k,j,i)
+                   ENDDO
+                ENDDO
+             ENDDO
+          ELSEIF ( intermediate_timestep_count < &
+                   intermediate_timestep_count_max )  THEN
+             DO  i = nxl, nxr
+                DO  j = nys, nyn
+                   DO  k = nzb_s_inner(j,i)+1, nzt
+                      ts_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * ts_m(k,j,i)
+                   ENDDO
+                ENDDO
+             ENDDO
+          ENDIF
+       ENDIF
+
+       CALL cpu_log( log_point(66), 's-equation', 'stop' )
+
+    ENDIF
 !
 !-- If required, compute prognostic equation for turbulent kinetic 
@@ -2226,8 +2401,8 @@
 
 !
-!-- If required, compute prognostic equation for total water content / scalar
-    IF ( humidity  .OR.  passive_scalar )  THEN
-
-       CALL cpu_log( log_point(29), 'q/s-equation', 'start' )
+!-- If required, compute prognostic equation for total water content
+    IF ( humidity )  THEN
+
+       CALL cpu_log( log_point(29), 'q-equation', 'start' )
 
 !
@@ -2307,5 +2482,5 @@
        ENDDO
 
-       CALL cpu_log( log_point(29), 'q/s-equation', 'stop' )
+       CALL cpu_log( log_point(29), 'q-equation', 'stop' )
 
 !
@@ -2430,4 +2605,91 @@
     ENDIF
 
+!
+!-- If required, compute prognostic equation for scalar
+    IF ( passive_scalar )  THEN
+
+       CALL cpu_log( log_point(66), 's-equation', 'start' )
+
+!
+!--    Scalar/q-tendency terms with communication
+       sbt = tsc(2)
+       IF ( scalar_advec == 'bc-scheme' )  THEN
+
+          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
+!
+!--          Bott-Chlond scheme always uses Euler time step. Thus:
+             sbt = 1.0_wp
+          ENDIF
+          tend = 0.0_wp
+          CALL advec_s_bc( s, 's' )
+
+       ENDIF
+
+!
+!--    Scalar/q-tendency terms with no communication
+       IF ( scalar_advec /= 'bc-scheme' )  THEN
+          tend = 0.0_wp
+          IF ( timestep_scheme(1:5) == 'runge' )  THEN
+             IF ( ws_scheme_sca )  THEN
+                CALL advec_s_ws( s, 's' )
+             ELSE
+                CALL advec_s_pw( s )
+             ENDIF
+          ELSE
+             CALL advec_s_up( s )
+          ENDIF
+       ENDIF
+
+       CALL diffusion_s( s, ssws, sswst, wall_sflux )
+
+!
+!--    Sink or source of scalar concentration due to canopy elements
+       IF ( plant_canopy ) CALL pcm_tendency( 7 )
+
+!
+!--    Large scale advection. Not implemented so far.
+!        IF ( large_scale_forcing )  THEN
+!           CALL ls_advec( simulated_time, 's' )
+!        ENDIF
+
+!
+!--    Nudging. Not implemented so far.
+!        IF ( nudging )  CALL nudge( simulated_time, 's' ) 
+
+!
+!--    If required compute influence of large-scale subsidence/ascent.
+!--    Not implemented so far.
+       IF ( large_scale_subsidence  .AND.                                      &
+            .NOT. use_subsidence_tendencies  .AND.                             &
+            .NOT. large_scale_forcing )  THEN
+         CALL subsidence( tend, s, s_init, 3 )
+       ENDIF
+
+       CALL user_actions( 's-tendency' )
+
+!
+!--    Prognostic equation for total water content / scalar
+       DO  i = i_left, i_right
+          DO  j = j_south, j_north
+             DO  k = nzb_s_inner(j,i)+1, nzt
+                s_p(k,j,i) = s(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +          &
+                                                  tsc(3) * ts_m(k,j,i) )       &
+                                      - tsc(5) * rdf_sc(k) *                   &
+                                        ( s(k,j,i) - s_init(k) )
+                IF ( s_p(k,j,i) < 0.0_wp )  s_p(k,j,i) = 0.1_wp * s(k,j,i)
+!
+!--             Tendencies for the next Runge-Kutta step
+                IF ( runge_step == 1 )  THEN
+                   ts_m(k,j,i) = tend(k,j,i)
+                ELSEIF ( runge_step == 2 )  THEN
+                   ts_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * ts_m(k,j,i)
+                ENDIF
+             ENDDO
+          ENDDO
+       ENDDO
+
+       CALL cpu_log( log_point(66), 's-equation', 'stop' )
+
+    ENDIF
 !
 !-- If required, compute prognostic equation for turbulent kinetic
