Index: /palm/trunk/SOURCE/advec_ws.f90
===================================================================
--- /palm/trunk/SOURCE/advec_ws.f90	(revision 1114)
+++ /palm/trunk/SOURCE/advec_ws.f90	(revision 1115)
@@ -20,4 +20,5 @@
 ! Current revisions:
 ! ------------------
+! calculation of qr and nr is restricted to precipitation
 !
 ! Former revisions:
@@ -210,10 +211,11 @@
           sums_wspts_ws_l = 0.0
 
-          IF ( humidity .OR. passive_scalar )  THEN
+          IF ( humidity  .OR.  passive_scalar )  THEN
              ALLOCATE( sums_wsqs_ws_l(nzb:nzt+1,0:threads_per_task-1) )
              sums_wsqs_ws_l = 0.0
           ENDIF
 
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+          IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+               precipitation )  THEN
              ALLOCATE( sums_wsqrs_ws_l(nzb:nzt+1,0:threads_per_task-1) )
              ALLOCATE( sums_wsnrs_ws_l(nzb:nzt+1,0:threads_per_task-1) )
@@ -264,5 +266,5 @@
                        diss_l_e(nzb+1:nzt,nys:nyn,0:threads_per_task-1) )
 
-             IF ( humidity .OR. passive_scalar )  THEN
+             IF ( humidity  .OR.  passive_scalar )  THEN
                 ALLOCATE( flux_s_q(nzb+1:nzt,0:threads_per_task-1),          &
                           diss_s_q(nzb+1:nzt,0:threads_per_task-1) )
@@ -271,5 +273,6 @@
              ENDIF
 
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.            &
+                  precipitation )  THEN
                 ALLOCATE( flux_s_qr(nzb+1:nzt,0:threads_per_task-1),         &
                           diss_s_qr(nzb+1:nzt,0:threads_per_task-1),         &
@@ -280,5 +283,5 @@
                           flux_l_nr(nzb+1:nzt,nys:nyn,0:threads_per_task-1), &
                           diss_l_nr(nzb+1:nzt,nys:nyn,0:threads_per_task-1) ) 
-             END IF
+             ENDIF
 
              IF ( ocean )  THEN
@@ -297,5 +300,5 @@
 
 !------------------------------------------------------------------------------!
-! Initialize variables used for storing statistic qauntities (fluxes, variances)
+! Initialize variables used for storing statistic quantities (fluxes, variances)
 !------------------------------------------------------------------------------!
     SUBROUTINE ws_statistics
@@ -319,6 +322,7 @@
        IF ( ws_scheme_sca )  THEN
           sums_wspts_ws_l = 0.0
-          IF ( humidity .OR. passive_scalar )  sums_wsqs_ws_l = 0.0
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+          IF ( humidity  .OR.  passive_scalar )  sums_wsqs_ws_l = 0.0
+          IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+               precipitation )  THEN
              sums_wsqrs_ws_l = 0.0
              sums_wsnrs_ws_l = 0.0
@@ -770,5 +774,4 @@
                                  * weight_substep(intermediate_timestep_count)
              ENDDO
-
 
          END SELECT
Index: /palm/trunk/SOURCE/average_3d_data.f90
===================================================================
--- /palm/trunk/SOURCE/average_3d_data.f90	(revision 1114)
+++ /palm/trunk/SOURCE/average_3d_data.f90	(revision 1115)
@@ -20,4 +20,5 @@
 ! Current revisions:
 ! -----------------
+! +qc
 !
 ! Former revisions:
@@ -181,4 +182,13 @@
              ENDDO
 
+          CASE ( 'qc' )
+             DO  i = nxlg, nxrg
+                DO  j = nysg, nyng
+                   DO  k = nzb, nzt+1
+                      qc_av(k,j,i) = qc_av(k,j,i) / REAL( average_count_3d )
+                   ENDDO
+                ENDDO
+             ENDDO
+
           CASE ( 'ql' )
              DO  i = nxlg, nxrg
Index: /palm/trunk/SOURCE/boundary_conds.f90
===================================================================
--- /palm/trunk/SOURCE/boundary_conds.f90	(revision 1114)
+++ /palm/trunk/SOURCE/boundary_conds.f90	(revision 1115)
@@ -20,5 +20,6 @@
 ! Current revisions:
 ! -----------------
-! 
+! boundary conditions of two-moment cloud scheme are restricted to Neumann-
+! boundary-conditions
 !  
 ! Former revisions:
@@ -246,43 +247,20 @@
        q_p(nzt+1,:,:) = q_p(nzt,:,:)   + bc_q_t_val * dzu(nzt+1)
 
-       IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+       IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+            precipitation )  THEN
 !             
-!--       Surface conditions for constant_humidity_flux
-          IF ( ibc_qr_b == 0 ) THEN
-             DO  i = nxlg, nxrg
-                DO  j = nysg, nyng
-                   qr_p(nzb_s_inner(j,i),j,i) = qr(nzb_s_inner(j,i),j,i)
-                ENDDO
-             ENDDO
-          ELSE
-             DO  i = nxlg, nxrg
-                DO  j = nysg, nyng
-                   qr_p(nzb_s_inner(j,i),j,i) = qr_p(nzb_s_inner(j,i)+1,j,i)
-                ENDDO
-             ENDDO
-          ENDIF
-!
-!--       Top boundary
-          qr_p(nzt+1,:,:) = qr_p(nzt,:,:) + bc_qr_t_val * dzu(nzt+1)
-!             
-!--       Surface conditions for constant_humidity_flux
-          IF ( ibc_nr_b == 0 ) THEN
-             DO  i = nxlg, nxrg
-                DO  j = nysg, nyng
-                   nr_p(nzb_s_inner(j,i),j,i) = nr(nzb_s_inner(j,i),j,i)
-                ENDDO
-             ENDDO
-          ELSE
-             DO  i = nxlg, nxrg
-                DO  j = nysg, nyng
-                   nr_p(nzb_s_inner(j,i),j,i) = nr_p(nzb_s_inner(j,i)+1,j,i)
-                ENDDO
-             ENDDO
-          ENDIF
-!
-!--       Top boundary
-          nr_p(nzt+1,:,:) = nr_p(nzt,:,:) + bc_nr_t_val * dzu(nzt+1)
+!--       Surface conditions rain water (Neumann)
+          DO  i = nxlg, nxrg
+             DO  j = nysg, nyng
+                qr_p(nzb_s_inner(j,i),j,i) = qr_p(nzb_s_inner(j,i)+1,j,i)
+                nr_p(nzb_s_inner(j,i),j,i) = nr_p(nzb_s_inner(j,i)+1,j,i)
+             ENDDO
+          ENDDO
+!
+!--       Top boundary condition for rain water (Neumann)
+          qr_p(nzt+1,:,:) = qr_p(nzt,:,:)
+          nr_p(nzt+1,:,:) = nr_p(nzt,:,:)
+           
        ENDIF
-
 !
 !--    In case of inflow at the south boundary the boundary for v is at nys
@@ -308,7 +286,8 @@
           pt_p(:,nys-1,:)     = pt_p(:,nys,:)
           IF ( .NOT. constant_diffusion     )  e_p(:,nys-1,:) = e_p(:,nys,:)
-          IF ( humidity .OR. passive_scalar )  THEN
+          IF ( humidity  .OR.  passive_scalar )  THEN
              q_p(:,nys-1,:) = q_p(:,nys,:)
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation)  THEN
                 qr_p(:,nys-1,:) = qr_p(:,nys,:)
                 nr_p(:,nys-1,:) = nr_p(:,nys,:)
@@ -318,7 +297,8 @@
           pt_p(:,nyn+1,:)     = pt_p(:,nyn,:)
           IF ( .NOT. constant_diffusion     )  e_p(:,nyn+1,:) = e_p(:,nyn,:)
-          IF ( humidity .OR. passive_scalar )  THEN
+          IF ( humidity  .OR.  passive_scalar )  THEN
              q_p(:,nyn+1,:) = q_p(:,nyn,:)
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation )  THEN
                 qr_p(:,nyn+1,:) = qr_p(:,nyn,:)
                 nr_p(:,nyn+1,:) = nr_p(:,nyn,:)
@@ -328,7 +308,8 @@
           pt_p(:,:,nxl-1)     = pt_p(:,:,nxl)
           IF ( .NOT. constant_diffusion     )  e_p(:,:,nxl-1) = e_p(:,:,nxl)
-          IF ( humidity .OR. passive_scalar )  THEN
+          IF ( humidity  .OR.  passive_scalar )  THEN
              q_p(:,:,nxl-1) = q_p(:,:,nxl)
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation )  THEN
                 qr_p(:,:,nxl-1) = qr_p(:,:,nxl)
                 nr_p(:,:,nxl-1) = nr_p(:,:,nxl)
@@ -340,5 +321,5 @@
           IF ( humidity .OR. passive_scalar )  THEN
              q_p(:,:,nxr+1) = q_p(:,:,nxr)
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  precipitation )  THEN
                 qr_p(:,:,nxr+1) = qr_p(:,:,nxr)
                 nr_p(:,:,nxr+1) = nr_p(:,:,nxr)
Index: /palm/trunk/SOURCE/calc_liquid_water_content.f90
===================================================================
--- /palm/trunk/SOURCE/calc_liquid_water_content.f90	(revision 1114)
+++ /palm/trunk/SOURCE/calc_liquid_water_content.f90	(revision 1115)
@@ -20,4 +20,5 @@
 ! Current revisions:
 ! -----------------
+! drizzle can be used independently from precipitation
 !
 ! Former revisions:
@@ -63,4 +64,5 @@
     USE cloud_parameters
     USE constants
+    USE control_parameters
     USE grid_variables
     USE indices
@@ -75,5 +77,5 @@
     DO  i = nxlg, nxrg
        DO  j = nysg, nyng
-          DO  k = nzb_2d(j,i)+1, nzt
+          DO  k = nzb_s_inner(j,i)+1, nzt
 
 !
@@ -101,10 +103,27 @@
 !
 !--          Compute the liquid water content
-             IF ( ( q(k,j,i) - q_s ) > 0.0 ) THEN
-                ql(k,j,i) = q(k,j,i) - q_s
+             IF ( icloud_scheme == 0  .AND.  precipitation)  THEN
+                IF ( ( q(k,j,i) - q_s - qr(k,j,i) ) > 0.0 ) THEN
+                   qc(k,j,i) = q(k,j,i) - q_s - qr(k,j,i)
+                   ql(k,j,i) = qc(k,j,i) + qr(k,j,i)
+                ELSE
+                   qc(k,j,i) = 0.0 
+                   ql(k,j,i) = qr(k,j,i)
+                ENDIF
+             ELSEIF ( icloud_scheme == 0  .AND.  .NOT. precipitation )  THEN
+                IF ( ( q(k,j,i) - q_s ) > 0.0 ) THEN
+                   qc(k,j,i) = q(k,j,i) - q_s
+                   ql(k,j,i) = qc(k,j,i)
+                ELSE
+                   qc(k,j,i) = 0.0 
+                   ql(k,j,i) = 0.0
+                ENDIF
              ELSE
-                ql(k,j,i) = 0.0 
+                IF ( ( q(k,j,i) - q_s ) > 0.0 ) THEN
+                   ql(k,j,i) = q(k,j,i) - q_s
+                ELSE
+                   ql(k,j,i) = 0.0 
+                ENDIF
              ENDIF
-
           ENDDO
        ENDDO
Index: /palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- /palm/trunk/SOURCE/check_parameters.f90	(revision 1114)
+++ /palm/trunk/SOURCE/check_parameters.f90	(revision 1115)
@@ -20,5 +20,6 @@
 ! Current revisions:
 ! -----------------
-! 
+! unused variables removed
+! drizzle can be used without precipitation
 !
 ! Former revisions:
@@ -569,6 +570,9 @@
           WRITE( action, '(A)' )  'galilei_transformation = .TRUE.'
        ENDIF
+       IF ( cloud_physics )  THEN
+          WRITE( action, '(A)' )  'cloud_physics = .TRUE.'
+       ENDIF
        IF ( cloud_droplets )  THEN
-          WRITE( action, '(A)' )  'cloud_droplets = .TRUE. (has not been tested)'
+          WRITE( action, '(A)' )  'cloud_droplets = .TRUE.'
        ENDIF
        IF ( .NOT. prandtl_layer )  THEN
@@ -899,10 +903,10 @@
     ENDIF 
 
-    IF ( cloud_physics  .AND. icloud_scheme == 0  .AND.                      &
-         .NOT. precipitation ) THEN
-       message_string = 'cloud_scheme = seifert_beheng requires ' // &
-                        'precipitation = .TRUE.'
-       CALL message( 'check_parameters', 'PA0363', 1, 2, 0, 6, 0 )
-    ENDIF 
+!    IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+!         .NOT. precipitation  .AND.  .NOT. drizzle ) THEN
+!       message_string = 'cloud_scheme = seifert_beheng requires ' // &
+!                        'precipitation = .TRUE. or drizzle = .TRUE.'
+!       CALL message( 'check_parameters', 'PA0363', 1, 2, 0, 6, 0 )
+!    ENDIF 
 
 !
@@ -916,11 +920,4 @@
        IF ( humidity )  THEN
           q_init  = q_surface
-!
-!--       It is not allowed to choose initial profiles of rain water content 
-!--       and rain drop concentration. They are set to 0.0.
-          IF ( cloud_physics .AND. icloud_scheme == 0 ) THEN
-             qr_init = 0.0
-             nr_init = 0.0
-          ENDIF
        ENDIF
        IF ( ocean )           sa_init = sa_surface
@@ -1228,86 +1225,8 @@
              q_vertical_gradient_level(1) = 0.0
           ENDIF
-
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
-
-             i = 1
-             gradient = 0.0
-             qr_vertical_gradient_level_ind(1) = 0
-             DO  k = 1, nzt+1
-                IF ( i < 11 ) THEN
-                   IF ( qr_vertical_gradient_level(i) < zu(k)  .AND. &
-                        qr_vertical_gradient_level(i) >= 0.0 )  THEN
-                      gradient = qr_vertical_gradient(i) / 100.0
-                      qr_vertical_gradient_level_ind(i) = k - 1
-                      i = i + 1
-                   ENDIF
-                ENDIF
-                IF ( gradient /= 0.0 )  THEN
-                   IF ( k /= 1 )  THEN
-                      qr_init(k) = qr_init(k-1) + dzu(k) * gradient
-                   ELSE
-                      qr_init(k) = qr_init(k-1) + 0.5 * dzu(k) * gradient
-                   ENDIF
-                ELSE
-                   qr_init(k) = qr_init(k-1)
-                ENDIF
-!
-!--             Avoid negative rain water content
-                IF ( qr_init(k) < 0.0 )  THEN
-                   qr_init(k) = 0.0
-                ENDIF
-             ENDDO
-!
-!--          In case of no given rain water content gradients, choose zero gradient
-!--          conditions
-             IF ( qr_vertical_gradient_level(1) == -1.0 )  THEN
-                qr_vertical_gradient_level(1) = 0.0
-             ENDIF
-
-             i = 1
-             gradient = 0.0
-             nr_vertical_gradient_level_ind(1) = 0
-             DO  k = 1, nzt+1
-                IF ( i < 11 ) THEN
-                   IF ( nr_vertical_gradient_level(i) < zu(k)  .AND. &
-                        nr_vertical_gradient_level(i) >= 0.0 )  THEN
-                      gradient = nr_vertical_gradient(i) / 100.0
-                      nr_vertical_gradient_level_ind(i) = k - 1
-                      i = i + 1
-                   ENDIF
-                ENDIF
-                IF ( gradient /= 0.0 )  THEN
-                   IF ( k /= 1 )  THEN
-                      nr_init(k) = nr_init(k-1) + dzu(k) * gradient
-                   ELSE
-                      nr_init(k) = nr_init(k-1) + 0.5 * dzu(k) * gradient
-                   ENDIF
-                ELSE
-                   nr_init(k) = nr_init(k-1)
-                ENDIF
-!
-!--             Avoid negative rain water content
-                IF ( nr_init(k) < 0.0 )  THEN
-                   nr_init(k) = 0.0
-                ENDIF
-             ENDDO
-!
-!--          In case of no given rain water content gradients, choose zero gradient
-!--          conditions
-             IF ( nr_vertical_gradient_level(1) == -1.0 )  THEN
-                nr_vertical_gradient_level(1) = 0.0
-             ENDIF
-
-          ENDIF
 !
 !--       Store humidity, rain water content and rain drop concentration 
 !--       gradient at the top boundary for possile Neumann boundary condition
           bc_q_t_val  = ( q_init(nzt+1) - q_init(nzt) ) / dzu(nzt+1)
-            
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
-             bc_qr_t_val = ( qr_init(nzt+1) - qr_init(nzt) ) / dzu(nzt+1)
-             bc_nr_t_val = ( nr_init(nzt+1) - nr_init(nzt) ) / dzu(nzt+1)
-          ENDIF
-
        ENDIF
 
@@ -1787,45 +1706,5 @@
        ENDIF
 
-       IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
-          IF ( bc_qr_b == 'dirichlet' )  THEN
-             ibc_qr_b = 0
-          ELSEIF ( bc_qr_b == 'neumann' )  THEN
-             ibc_qr_b = 1
-          ELSE
-             message_string = 'unknown boundary condition: bc_qr_b ="' // &
-                              TRIM( bc_qr_b ) // '"'
-             CALL message( 'check_parameters', 'PA0352', 1, 2, 0, 6, 0 )
-          ENDIF
-          IF ( bc_qr_t == 'dirichlet' )  THEN
-             ibc_qr_t = 0
-          ELSEIF ( bc_qr_t == 'neumann' )  THEN
-             ibc_qr_t = 1
-          ELSE
-             message_string = 'unknown boundary condition: bc_qr_t ="' // &
-                              TRIM( bc_qr_t ) // '"'
-             CALL message( 'check_parameters', 'PA0353', 1, 2, 0, 6, 0 )
-          ENDIF
-          IF ( bc_nr_b == 'dirichlet' )  THEN
-             ibc_nr_b = 0
-          ELSEIF ( bc_nr_b == 'neumann' )  THEN
-             ibc_nr_b = 1
-          ELSE
-             message_string = 'unknown boundary condition: bc_nr_b ="' // &
-                              TRIM( bc_nr_b ) // '"'
-             CALL message( 'check_parameters', 'PA0355', 1, 2, 0, 6, 0 )
-          ENDIF
-          IF ( bc_nr_t == 'dirichlet' )  THEN
-             ibc_nr_t = 0
-          ELSEIF ( bc_nr_t == 'neumann' )  THEN
-             ibc_nr_t = 1
-          ELSE
-             message_string = 'unknown boundary condition: bc_nr_t ="' // &
-                              TRIM( bc_nr_t ) // '"'
-             CALL message( 'check_parameters', 'PA0356', 1, 2, 0, 6, 0 )
-          ENDIF
-       ENDIF       
-
-    ENDIF
-
+    ENDIF
 !
 !-- Boundary conditions for horizontal components of wind speed
@@ -2669,4 +2548,9 @@
                                  'lemented for cloud_scheme /= seifert_beheng'
                 CALL message( 'check_parameters', 'PA0358', 1, 2, 0, 6, 0 )
+             ELSEIF ( .NOT. precipitation )  THEN
+                message_string = 'data_output_pr = ' // &
+                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
+                                 'lemented for precipitation = .FALSE.'
+                CALL message( 'check_parameters', 'PA0361', 1, 2, 0, 6, 0 )
              ELSE
                 dopr_index(i) = 73
@@ -2686,4 +2570,9 @@
                                  'lemented for cloud_scheme /= seifert_beheng'
                 CALL message( 'check_parameters', 'PA0358', 1, 2, 0, 6, 0 )
+             ELSEIF ( .NOT. precipitation )  THEN
+                message_string = 'data_output_pr = ' // &
+                                 TRIM( data_output_pr(i) ) // ' is not imp' // &
+                                 'lemented for precipitation = .FALSE.'
+                CALL message( 'check_parameters', 'PA0361', 1, 2, 0, 6, 0 )
              ELSE
                 dopr_index(i) = 74
@@ -2906,4 +2795,8 @@
                          'res cloud_scheme = seifert_beheng'
                 CALL message( 'check_parameters', 'PA0359', 1, 2, 0, 6, 0 )
+             ELSEIF ( .NOT. precipitation )  THEN
+                message_string = 'output of "' // TRIM( var ) // '" requi' // &
+                                 'res precipitation = .TRUE.'
+                CALL message( 'check_parameters', 'PA0112', 1, 2, 0, 6, 0 )
              ENDIF
              unit = 'kg/kg'
Index: /palm/trunk/SOURCE/data_output_2d.f90
===================================================================
--- /palm/trunk/SOURCE/data_output_2d.f90	(revision 1114)
+++ /palm/trunk/SOURCE/data_output_2d.f90	(revision 1115)
@@ -20,4 +20,5 @@
 ! Current revisions:
 ! -----------------
+! ql is calculated by calc_liquid_water_content
 !
 ! Former revisions:
@@ -483,37 +484,15 @@
              CASE ( 'qc_xy', 'qc_xz', 'qc_yz' )
                 IF ( av == 0 )  THEN
+                   to_be_resorted => qc
+                ELSE
+                   to_be_resorted => qc_av
+                ENDIF
+                IF ( mode == 'xy' )  level_z = zu
+
+             CASE ( 'ql_xy', 'ql_xz', 'ql_yz' )
+                IF ( av == 0 )  THEN
                    to_be_resorted => ql
                 ELSE
                    to_be_resorted => ql_av
-                ENDIF
-                IF ( mode == 'xy' )  level_z = zu
-
-             CASE ( 'ql_xy', 'ql_xz', 'ql_yz' )
-                IF ( av == 0 )  THEN
-                   IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
-                      DO  i = nxlg, nxrg
-                         DO  j = nysg, nyng
-                            DO  k = nzb, nzt+1
-                               local_pf(i,j,k) = ql(k,j,i) + qr(k,j,i)
-                            ENDDO
-                         ENDDO
-                      ENDDO
-                      resorted = .TRUE.
-                   ELSE
-                      to_be_resorted => ql
-                   ENDIF
-                ELSE
-                   IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
-                      DO  i = nxlg, nxrg
-                         DO  j = nysg, nyng
-                            DO  k = nzb, nzt+1
-                               local_pf(i,j,k) = ql_av(k,j,i) + qr_av(k,j,i)
-                            ENDDO
-                         ENDDO
-                      ENDDO
-                      resorted = .TRUE.
-                   ELSE
-                      to_be_resorted => ql_av
-                   ENDIF
                 ENDIF
                 IF ( mode == 'xy' )  level_z = zu
Index: /palm/trunk/SOURCE/data_output_3d.f90
===================================================================
--- /palm/trunk/SOURCE/data_output_3d.f90	(revision 1114)
+++ /palm/trunk/SOURCE/data_output_3d.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! ql is calculated by calc_liquid_water_content
 !
 ! Former revisions:
@@ -337,36 +337,14 @@
           CASE ( 'qc' )
              IF ( av == 0 )  THEN
+                to_be_resorted => qc
+             ELSE
+                to_be_resorted => qc_av
+             ENDIF
+
+          CASE ( 'ql' )
+             IF ( av == 0 )  THEN
                 to_be_resorted => ql
              ELSE
                 to_be_resorted => ql_av
-             ENDIF
-
-          CASE ( 'ql' )
-             IF ( av == 0 )  THEN
-                IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
-                   DO  i = nxlg, nxrg
-                      DO  j = nysg, nyng
-                         DO  k = nzb, nz_do3d
-                            local_pf(i,j,k) = ql(k,j,i) + qr(k,j,i)
-                         ENDDO
-                      ENDDO
-                   ENDDO
-                   resorted = .TRUE.
-                ELSE
-                   to_be_resorted => ql
-                ENDIF
-             ELSE
-                IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
-                   DO  i = nxlg, nxrg
-                      DO  j = nysg, nyng
-                         DO  k = nzb, nz_do3d
-                            local_pf(i,j,k) = ql_av(k,j,i) + qr_av(k,j,i)
-                         ENDDO
-                      ENDDO
-                   ENDDO
-                   resorted = .TRUE.
-                ELSE
-                   to_be_resorted => ql_av
-                ENDIF
              ENDIF
 
Index: /palm/trunk/SOURCE/flow_statistics.f90
===================================================================
--- /palm/trunk/SOURCE/flow_statistics.f90	(revision 1114)
+++ /palm/trunk/SOURCE/flow_statistics.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! ql is calculated by calc_liquid_water_content
 !
 ! Former revisions:
@@ -766,29 +766,30 @@
                                                           rmask(j,i,sr)
                       IF ( .NOT. cloud_droplets )  THEN
-                         pts = 0.5 *                                          &
-                              ( ( q(k,j,i) - ql(k,j,i) ) -                    &
-                              hom(k,1,42,sr) +                                &
-                              ( q(k+1,j,i) - ql(k+1,j,i) ) -                  &
+                         pts = 0.5 *                                           &
+                              ( ( q(k,j,i) - ql(k,j,i) ) -                     &
+                              hom(k,1,42,sr) +                                 &
+                              ( q(k+1,j,i) - ql(k+1,j,i) ) -                   &
                               hom(k+1,1,42,sr) )
-                         sums_l(k,52,tn) = sums_l(k,52,tn) + pts * w(k,j,i) * &
+                         sums_l(k,52,tn) = sums_l(k,52,tn) + pts * w(k,j,i) *  &
                                                              rmask(j,i,sr)
                          IF ( icloud_scheme == 0  )  THEN
-                            sums_l(k,54,tn) = sums_l(k,54,tn) + ( ql(k,j,i) + &
-                                                                qr(k,j,i) ) * &
+                            sums_l(k,54,tn) = sums_l(k,54,tn) + ql(k,j,i) *    &
                                                                 rmask(j,i,sr)
-                            sums_l(k,73,tn) = sums_l(k,73,tn) + nr(k,j,i) *   &
+                            sums_l(k,75,tn) = sums_l(k,75,tn) + qc(k,j,i) *    &
                                                                 rmask(j,i,sr)
-                            sums_l(k,74,tn) = sums_l(k,74,tn) + qr(k,j,i) *   &
-                                                                rmask(j,i,sr)
-                            sums_l(k,75,tn) = sums_l(k,75,tn) + ql(k,j,i) *   &
-                                                                rmask(j,i,sr)
-                            sums_l(k,76,tn) = sums_l(k,76,tn) + prr(k,j,i) *  &
-                                                                rmask(j,i,sr)
+                            IF ( precipitation )  THEN
+                               sums_l(k,73,tn) = sums_l(k,73,tn) + nr(k,j,i) * &
+                                                                   rmask(j,i,sr)
+                               sums_l(k,74,tn) = sums_l(k,74,tn) + qr(k,j,i) * &
+                                                                   rmask(j,i,sr)
+                               sums_l(k,76,tn) = sums_l(k,76,tn) + prr(k,j,i) *&
+                                                                   rmask(j,i,sr)
+                            ENDIF
                          ELSE
-                            sums_l(k,54,tn) = sums_l(k,54,tn) + ql(k,j,i) *   &
+                            sums_l(k,54,tn) = sums_l(k,54,tn) + ql(k,j,i) *    &
                                                                 rmask(j,i,sr)
                          ENDIF
                       ELSE
-                         sums_l(k,54,tn) = sums_l(k,54,tn) + ql(k,j,i) *      &
+                         sums_l(k,54,tn) = sums_l(k,54,tn) + ql(k,j,i) *       &
                                                              rmask(j,i,sr)
                       ENDIF
Index: /palm/trunk/SOURCE/header.f90
===================================================================
--- /palm/trunk/SOURCE/header.f90	(revision 1114)
+++ /palm/trunk/SOURCE/header.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! descriptions for Seifert-Beheng-cloud-physics-scheme added
 !
 ! Former revisions:
@@ -1237,4 +1237,8 @@
        WRITE ( io, 415 )
        WRITE ( io, 416 ) surface_pressure, r_d, rho_surface, cp, l_v
+       IF ( icloud_scheme == 0 )  THEN
+          WRITE ( io, 510 ) 1.0E-6 * nc_const
+          IF ( precipitation )  WRITE ( io, 511 ) c_sedimentation
+       ENDIF
     ENDIF
 
@@ -1458,6 +1462,16 @@
     ELSEIF ( humidity  .AND.  cloud_physics )  THEN
        WRITE ( io, 432 )
-       IF ( radiation )      WRITE ( io, 132 )
-       IF ( precipitation )  WRITE ( io, 133 )
+       IF ( radiation )  WRITE ( io, 132 )
+       IF ( icloud_scheme == 1 )  THEN
+          IF ( precipitation )  WRITE ( io, 133 )
+       ELSEIF ( icloud_scheme == 0 )  THEN
+          IF ( drizzle )  WRITE ( io, 506 )
+          IF ( precipitation )  THEN
+             WRITE ( io, 505 )
+             IF ( turbulence )  WRITE ( io, 507 )
+             IF ( ventilation_effect )  WRITE ( io, 508 )
+             IF ( limiter_sedimentation )  WRITE ( io, 509 )
+          ENDIF
+       ENDIF
     ELSEIF ( humidity  .AND.  cloud_droplets )  THEN
        WRITE ( io, 433 )
@@ -2020,5 +2034,12 @@
 503 FORMAT (' --> Momentum advection via Wicker-Skamarock-Scheme 5th order')
 504 FORMAT (' --> Scalar advection via Wicker-Skamarock-Scheme 5th order')
-
+505 FORMAT ('    Precipitation parameterization via Seifert-Beheng-Scheme')
+506 FORMAT ('    Drizzle parameterization via Stokes law')
+507 FORMAT ('    Turbulence effects on precipitation process')
+508 FORMAT ('    Ventilation effects on evaporation of rain drops')
+509 FORMAT ('    Slope limiter used for sedimentation process')
+510 FORMAT ('        Droplet density    :   N_c   = ',F6.1,' 1/cm**3')
+511 FORMAT ('        Sedimentation Courant number:                  '/&
+            '                               C_s   = ',F3.1,'        ')
 
  END SUBROUTINE header
Index: /palm/trunk/SOURCE/init_3d_model.f90
===================================================================
--- /palm/trunk/SOURCE/init_3d_model.f90	(revision 1114)
+++ /palm/trunk/SOURCE/init_3d_model.f90	(revision 1115)
@@ -23,5 +23,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! unused variables removed
 !
 ! Former revisions:
@@ -392,44 +392,55 @@
              IF ( icloud_scheme == 0 )  THEN
 !
-!--             1D-rain sedimentation fluxes and rain drop size distribution
-!--             properties 
-                ALLOCATE ( dr(nzb:nzt+1), lambda_r(nzb:nzt+1),  &
-                           mu_r(nzb:nzt+1), sed_q(nzb:nzt+1),   &
-                           sed_qr(nzb:nzt+1), sed_nr(nzb:nzt+1),& 
-                           xr(nzb:nzt+1) ) 
-!
-!--             2D-rain water content and rain drop concentration arrays
-                ALLOCATE ( qrs(nysg:nyng,nxlg:nxrg),   &
-                   qrsws(nysg:nyng,nxlg:nxrg),         &
-                   qrswst(nysg:nyng,nxlg:nxrg),        &
-                   nrs(nysg:nyng,nxlg:nxrg),           &
-                   nrsws(nysg:nyng,nxlg:nxrg),         &
-                   nrswst(nysg:nyng,nxlg:nxrg) )
-!
-!--             3D-rain water content, rain drop concentration arrays
+!--             1D-arrays 
+                ALLOCATE ( nc_1d(nzb:nzt+1), pt_1d(nzb:nzt+1), &
+                           q_1d(nzb:nzt+1), qc_1d(nzb:nzt+1) ) 
+!
+!--             3D-cloud water content
 #if defined( __nopointer )
-                ALLOCATE( nr(nzb:nzt+1,nysg:nyng,nxlg:nxrg),         &
-                          nr_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
-                          qr(nzb:nzt+1,nysg:nyng,nxlg:nxrg),         &
-                          qr_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
-                          tnr_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg),      &
-                          tqr_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+                ALLOCATE( qc(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
 #else
-                ALLOCATE( nr_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
-                          nr_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
-                          nr_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
-                          qr_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
-                          qr_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
-                          qr_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+                ALLOCATE( qc_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
 #endif
 !
 !--             3D-tendency arrays
-                ALLOCATE( tend_nr(nzb:nzt+1,nysg:nyng,nxlg:nxrg),    &
-                          tend_pt(nzb:nzt+1,nysg:nyng,nxlg:nxrg),    &
-                          tend_q(nzb:nzt+1,nysg:nyng,nxlg:nxrg),     &
-                          tend_qr(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
-!
-!--             3d-precipitation rate
+                ALLOCATE( tend_pt(nzb:nzt+1,nysg:nyng,nxlg:nxrg),    &
+                          tend_q(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+
                 IF ( precipitation )  THEN
+!
+!--                1D-arrays 
+                   ALLOCATE ( nr_1d(nzb:nzt+1), qr_1d(nzb:nzt+1) ) 
+!
+!
+!--                3D-tendency arrays
+                   ALLOCATE( tend_nr(nzb:nzt+1,nysg:nyng,nxlg:nxrg),    &
+                             tend_qr(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+!
+!--                2D-rain water content and rain drop concentration arrays
+                   ALLOCATE ( qrs(nysg:nyng,nxlg:nxrg),                 &
+                              qrsws(nysg:nyng,nxlg:nxrg),               &
+                              qrswst(nysg:nyng,nxlg:nxrg),              &
+                              nrs(nysg:nyng,nxlg:nxrg),                 &
+                              nrsws(nysg:nyng,nxlg:nxrg),               &
+                              nrswst(nysg:nyng,nxlg:nxrg) )
+!
+!--                3D-rain water content, rain drop concentration arrays
+#if defined( __nopointer )
+                   ALLOCATE( nr(nzb:nzt+1,nysg:nyng,nxlg:nxrg),         &
+                             nr_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
+                             qr(nzb:nzt+1,nysg:nyng,nxlg:nxrg),         &
+                             qr_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
+                             tnr_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg),      &
+                             tqr_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+#else
+                   ALLOCATE( nr_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
+                             nr_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
+                             nr_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
+                             qr_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
+                             qr_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg),       &
+                             qr_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+#endif
+!
+!--                3d-precipitation rate
                    ALLOCATE( prr(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
                 ENDIF
@@ -582,6 +593,9 @@
              ql => ql_1
              IF ( icloud_scheme == 0 )  THEN
-                qr => qr_1;  qr_p  => qr_2;  tqr_m  => qr_3
-                nr => nr_1;  nr_p  => nr_2;  tnr_m  => nr_3
+                qc => qc_1
+                IF ( precipitation )  THEN
+                   qr => qr_1;  qr_p  => qr_2;  tqr_m  => qr_3
+                   nr => nr_1;  nr_p  => nr_2;  tnr_m  => nr_3
+                ENDIF
              ENDIF
           ENDIF
@@ -641,11 +655,16 @@
                 ENDDO
              ENDDO
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation )  THEN
                 DO  i = nxlg, nxrg
                    DO  j = nysg, nyng
-                      qr(:,j,i) = qr_init
-                      nr(:,j,i) = nr_init
+                      qr(:,j,i) = 0.0
+                      nr(:,j,i) = 0.0
                    ENDDO
                 ENDDO
+!
+!--             Initialze nc_1d with default value
+                nc_1d(:) = nc_const
+
              ENDIF
           ENDIF
@@ -693,5 +712,6 @@
           IF ( humidity  .OR.  passive_scalar )  THEN
              qs = 0.0
-             IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation )  THEN
                 qrs = 0.0
                 nrs = 0.0
@@ -762,11 +782,18 @@
                 ENDDO
              ENDDO
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
-                DO  i = nxlg, nxrg
-                   DO  j = nysg, nyng
-                      qr(:,j,i) = qr_init
-                      nr(:,j,i) = nr_init
+             IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
+!
+!--             Initialze nc_1d with default value
+                nc_1d(:) = nc_const
+
+                IF ( precipitation )  THEN
+                   DO  i = nxlg, nxrg
+                      DO  j = nysg, nyng
+                         qr(:,j,i) = 0.0
+                         nr(:,j,i) = 0.0
+                      ENDDO
                    ENDDO
-                ENDDO
+                ENDIF
+
              ENDIF
           ENDIF
@@ -910,11 +937,8 @@
 !--       Determine the near-surface water flux
           IF ( humidity  .OR.  passive_scalar )  THEN
-             IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
-                IF ( constant_waterflux_qr )  THEN
-                   qrsws = surface_waterflux_qr
-                ENDIF
-                IF (constant_waterflux_nr )  THEN
-                   nrsws = surface_waterflux_nr
-                ENDIF
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation )  THEN
+                qrsws = 0.0
+                nrsws = 0.0
              ENDIF
              IF ( constant_waterflux )  THEN
@@ -951,5 +975,6 @@
              IF ( humidity  .OR.  passive_scalar )  THEN
                 qswst = 0.0
-                IF ( cloud_physics  .AND.  icloud_scheme == 0 ) THEN
+                IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                     precipitation ) THEN
                    nrswst = 0.0
                    qrswst = 0.0
@@ -990,11 +1015,8 @@
           IF ( humidity  .OR.  passive_scalar )  THEN
              IF ( .NOT. constant_waterflux )  qsws   = 0.0
-             IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
-                IF ( .NOT. constant_waterflux_qr )  THEN
-                   qrsws = 0.0
-                ENDIF
-                IF ( .NOT. constant_waterflux_nr )  THEN
-                   nrsws = 0.0
-                ENDIF
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation )  THEN
+                qrsws = 0.0
+                nrsws = 0.0
              ENDIF
           ENDIF
@@ -1021,4 +1043,8 @@
           ql = 0.0
           IF ( precipitation )  precipitation_amount = 0.0
+          IF ( icloud_scheme == 0 )  THEN
+             qc = 0.0
+             nc_1d = nc_const
+          ENDIF 
        ENDIF
 !
@@ -1046,13 +1072,5 @@
             q_surface_initial_change /= 0.0 )  THEN
           q(nzb,:,:) = q(nzb,:,:) + q_surface_initial_change
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
-             IF ( qr_surface_initial_change /= 0.0 )  THEN
-                qr(nzb,:,:) = qr(nzb,:,:) + qr_surface_initial_change
-             ELSEIF ( nr_surface_initial_change /= 0.0 ) THEN
-                nr(nzb,:,:) = nr(nzb,:,:) + nr_surface_initial_change
-             ENDIF
-          ENDIF
-       ENDIF
-
+       ENDIF
 !
 !--    Initialize the random number generator (from numerical recipes)
@@ -1067,5 +1085,6 @@
           tq_m = 0.0
           q_p = q
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+          IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+               precipitation )  THEN
              tqr_m = 0.0
              qr_p = qr
@@ -1245,5 +1264,6 @@
        IF ( humidity  .OR.  passive_scalar )  THEN
           q_p = q
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+          IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+               precipitation )  THEN
              qr_p = qr
              nr_p = nr
@@ -1259,5 +1279,6 @@
        IF ( humidity  .OR.  passive_scalar )  THEN
           tq_m = 0.0
-          IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+          IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+               precipitation )  THEN
              tqr_m = 0.0
              tnr_m = 0.0
Index: /palm/trunk/SOURCE/microphysics.f90
===================================================================
--- /palm/trunk/SOURCE/microphysics.f90	(revision 1114)
+++ /palm/trunk/SOURCE/microphysics.f90	(revision 1115)
@@ -20,5 +20,6 @@
 ! Current revisions:
 ! ------------------
-! 
+! microphyical tendencies are calculated in microphysics_control in an optimized
+! way; unrealistic values are prevented; bugfix in evaporation; some reformatting
 !
 ! Former revisions:
@@ -35,5 +36,5 @@
 ! 1065 2012-11-22 17:42:36Z hoffmann
 ! Sedimentation process implemented according to Stevens and Seifert (2008).
-! Turbulence effects on autoconversion and accretion added (Seifert, Nuijens 
+! Turbulence effects on autoconversion and accretion added (Seifert, Nuijens
 ! and Stevens, 2010).
 !
@@ -48,11 +49,15 @@
 
     PRIVATE
-    PUBLIC dsd_properties, autoconversion, accretion, selfcollection_breakup, &
-           evaporation_rain, sedimentation_cloud, sedimentation_rain
-
-    INTERFACE dsd_properties
-       MODULE PROCEDURE dsd_properties
-       MODULE PROCEDURE dsd_properties_ij
-    END INTERFACE dsd_properties
+    PUBLIC microphysics_control
+
+    INTERFACE microphysics_control
+       MODULE PROCEDURE microphysics_control
+       MODULE PROCEDURE microphysics_control_ij
+    END INTERFACE microphysics_control
+
+    INTERFACE adjust_cloud
+       MODULE PROCEDURE adjust_cloud
+       MODULE PROCEDURE adjust_cloud_ij
+    END INTERFACE adjust_cloud
 
     INTERFACE autoconversion
@@ -92,10 +97,10 @@
 ! Call for all grid points
 !------------------------------------------------------------------------------!
-    SUBROUTINE dsd_properties
-
-       USE arrays_3d
-       USE cloud_parameters
-       USE constants
-       USE indices
+    SUBROUTINE microphysics_control
+
+       USE arrays_3d
+       USE control_parameters
+       USE indices
+       USE statistics
 
        IMPLICIT NONE
@@ -106,5 +111,5 @@
        DO  i = nxl, nxr
           DO  j = nys, nyn
-             DO  k = nzb_2d(j,i)+1, nzt
+             DO  k = nzb_s_inner(j,i)+1, nzt
 
              ENDDO
@@ -112,12 +117,10 @@
        ENDDO
 
-    END SUBROUTINE dsd_properties
-
-
-    SUBROUTINE autoconversion
-
-       USE arrays_3d
-       USE cloud_parameters
-       USE constants
+    END SUBROUTINE microphysics_control
+
+    SUBROUTINE adjust_cloud
+
+       USE arrays_3d
+       USE cloud_parameters
        USE indices
 
@@ -129,5 +132,5 @@
        DO  i = nxl, nxr
           DO  j = nys, nyn
-             DO  k = nzb_2d(j,i)+1, nzt
+             DO  k = nzb_s_inner(j,i)+1, nzt
 
              ENDDO
@@ -135,12 +138,13 @@
        ENDDO
 
-    END SUBROUTINE autoconversion
-
-
-    SUBROUTINE accretion
-
-       USE arrays_3d
-       USE cloud_parameters
-       USE constants
+    END SUBROUTINE adjust_cloud
+
+
+    SUBROUTINE autoconversion
+
+       USE arrays_3d
+       USE cloud_parameters
+       USE control_parameters
+       USE grid_variables
        USE indices
 
@@ -152,5 +156,5 @@
        DO  i = nxl, nxr
           DO  j = nys, nyn
-             DO  k = nzb_2d(j,i)+1, nzt
+             DO  k = nzb_s_inner(j,i)+1, nzt
 
              ENDDO
@@ -158,12 +162,12 @@
        ENDDO
 
-    END SUBROUTINE accretion
-
-
-    SUBROUTINE selfcollection_breakup
-
-       USE arrays_3d
-       USE cloud_parameters
-       USE constants
+    END SUBROUTINE autoconversion
+
+
+    SUBROUTINE accretion
+
+       USE arrays_3d
+       USE cloud_parameters
+       USE control_parameters
        USE indices
 
@@ -175,5 +179,5 @@
        DO  i = nxl, nxr
           DO  j = nys, nyn
-             DO  k = nzb_2d(j,i)+1, nzt
+             DO  k = nzb_s_inner(j,i)+1, nzt
 
              ENDDO
@@ -181,12 +185,12 @@
        ENDDO
 
-    END SUBROUTINE selfcollection_breakup
-
-
-    SUBROUTINE evaporation_rain
-
-       USE arrays_3d
-       USE cloud_parameters
-       USE constants
+    END SUBROUTINE accretion
+
+
+    SUBROUTINE selfcollection_breakup
+
+       USE arrays_3d
+       USE cloud_parameters
+       USE control_parameters
        USE indices
 
@@ -198,5 +202,5 @@
        DO  i = nxl, nxr
           DO  j = nys, nyn
-             DO  k = nzb_2d(j,i)+1, nzt
+             DO  k = nzb_s_inner(j,i)+1, nzt
 
              ENDDO
@@ -204,12 +208,13 @@
        ENDDO
 
-    END SUBROUTINE evaporation_rain
-
-
-    SUBROUTINE sedimentation_cloud
+    END SUBROUTINE selfcollection_breakup
+
+
+    SUBROUTINE evaporation_rain
 
        USE arrays_3d
        USE cloud_parameters
        USE constants
+       USE control_parameters
        USE indices
 
@@ -221,5 +226,5 @@
        DO  i = nxl, nxr
           DO  j = nys, nyn
-             DO  k = nzb_2d(j,i)+1, nzt
+             DO  k = nzb_s_inner(j,i)+1, nzt
 
              ENDDO
@@ -227,12 +232,13 @@
        ENDDO
 
-    END SUBROUTINE sedimentation_cloud
-
-
-    SUBROUTINE sedimentation_rain
+    END SUBROUTINE evaporation_rain
+
+
+    SUBROUTINE sedimentation_cloud
 
        USE arrays_3d
        USE cloud_parameters
        USE constants
+       USE control_parameters
        USE indices
 
@@ -244,5 +250,30 @@
        DO  i = nxl, nxr
           DO  j = nys, nyn
-             DO  k = nzb_2d(j,i)+1, nzt
+             DO  k = nzb_s_inner(j,i)+1, nzt
+
+             ENDDO
+          ENDDO
+       ENDDO
+
+    END SUBROUTINE sedimentation_cloud
+
+
+    SUBROUTINE sedimentation_rain
+
+       USE arrays_3d
+       USE cloud_parameters
+       USE constants
+       USE control_parameters
+       USE indices
+       USE statistics
+
+       IMPLICIT NONE
+
+       INTEGER ::  i, j, k
+
+ 
+       DO  i = nxl, nxr
+          DO  j = nys, nyn
+             DO  k = nzb_s_inner(j,i)+1, nzt
 
              ENDDO
@@ -256,21 +287,72 @@
 ! Call for grid point i,j
 !------------------------------------------------------------------------------!
-    SUBROUTINE dsd_properties_ij( i, j )
-
-       USE arrays_3d
-       USE cloud_parameters
-       USE constants
-       USE indices
-       USE control_parameters
-       USE user
-   
-       IMPLICIT NONE
-
-       INTEGER ::  i, j, k
-
-       DO  k = nzb_2d(j,i)+1, nzt
+
+    SUBROUTINE microphysics_control_ij( i, j )
+
+       USE arrays_3d
+       USE cloud_parameters
+       USE control_parameters
+       USE statistics
+
+       IMPLICIT NONE
+
+       INTEGER ::  i, j
+
+       dt_micro = dt_3d * weight_pres(intermediate_timestep_count)
+!
+!--    Adjust unrealistic values
+       IF ( precipitation )  CALL adjust_cloud( i,j ) 
+!
+!--    Use 1-d arrays
+       q_1d(:)  = q(:,j,i)
+       pt_1d(:) = pt(:,j,i)
+       qc_1d(:) = qc(:,j,i)
+       nc_1d(:) = nc_const
+       IF ( precipitation )  THEN
+          qr_1d(:) = qr(:,j,i)
+          nr_1d(:) = nr(:,j,i)
+       ENDIF
+!
+!--    Compute cloud physics
+       IF ( precipitation )  THEN
+          CALL autoconversion( i,j )
+          CALL accretion( i,j )
+          CALL selfcollection_breakup( i,j )
+          CALL evaporation_rain( i,j )
+          CALL sedimentation_rain( i,j )
+       ENDIF
+
+       IF ( drizzle )  CALL sedimentation_cloud( i,j )
+!
+!--    Derive tendencies
+       tend_q(:,j,i)  = ( q_1d(:) - q(:,j,i) ) / dt_micro
+       tend_pt(:,j,i) = ( pt_1d(:) - pt(:,j,i) ) / dt_micro
+       IF ( precipitation )  THEN
+          tend_qr(:,j,i) = ( qr_1d(:) - qr(:,j,i) ) / dt_micro
+          tend_nr(:,j,i) = ( nr_1d(:) - nr(:,j,i) ) / dt_micro
+       ENDIF
+
+    END SUBROUTINE microphysics_control_ij
+
+    SUBROUTINE adjust_cloud_ij( i, j )
+
+       USE arrays_3d
+       USE cloud_parameters
+       USE indices
+
+       IMPLICIT NONE
+
+       INTEGER ::  i, j, k
+!
+!--    Adjust number of raindrops to avoid nonlinear effects in 
+!--    sedimentation and evaporation of rain drops due to too small or
+!--    too big weights of rain drops (Stevens and Seifert, 2008).
+!--    The same procedure is applied to cloud droplets if they are determined
+!--    prognostically.  
+       DO  k = nzb_s_inner(j,i)+1, nzt
 
           IF ( qr(k,j,i) <= eps_sb )  THEN
              qr(k,j,i) = 0.0
+             nr(k,j,i) = 0.0
           ELSE
 !
@@ -283,21 +365,10 @@
                 nr(k,j,i) = qr(k,j,i) * hyrho(k) / xrmax
              ENDIF
-             xr(k) = hyrho(k) * qr(k,j,i) / nr(k,j,i)
-!
-!--          Weight averaged diameter of rain drops:
-             dr(k) = ( hyrho(k) * qr(k,j,i) / nr(k,j,i) *                     &
-                       dpirho_l )**( 1.0 / 3.0 )
-!
-!--          Shape parameter of gamma distribution (Milbrandt and Yau, 2005;
-!--          Stevens and Seifert, 2008):
-             mu_r(k) = 10.0 * ( 1.0 + TANH( 1.2E3 * ( dr(k) - 1.4E-3 ) ) )
-!
-!--          Slope parameter of gamma distribution (Seifert, 2008):
-             lambda_r(k) = ( ( mu_r(k) + 3.0 ) * ( mu_r(k) + 2.0 ) *          &
-                             ( mu_r(k) + 1.0 ) )**( 1.0 / 3.0 ) / dr(k)
+
           ENDIF
-       ENDDO
-
-    END SUBROUTINE dsd_properties_ij
+
+       ENDDO
+
+    END SUBROUTINE adjust_cloud_ij
 
 
@@ -306,26 +377,24 @@
        USE arrays_3d
        USE cloud_parameters
-       USE constants
-       USE indices
-       USE control_parameters
-       USE statistics
+       USE control_parameters
        USE grid_variables
-   
-       IMPLICIT NONE
-
-       INTEGER ::  i, j, k
-       REAL    ::  k_au, autocon, phi_au, tau_cloud, xc, nu_c, rc,   &
-                   l_mix, re_lambda, alpha_cc, r_cc, sigma_cc, epsilon
+       USE indices
+
+       IMPLICIT NONE
+
+       INTEGER ::  i, j, k
+       REAL    ::  alpha_cc, autocon, epsilon, k_au, l_mix, nu_c, phi_au,      &
+                   r_cc, rc, re_lambda, selfcoll, sigma_cc, tau_cloud, xc                      
 
 
        k_au = k_cc / ( 20.0 * x0 )
 
-       DO  k = nzb_2d(j,i)+1, nzt
-
-          IF ( ql(k,j,i) > 0.0 )  THEN
+       DO  k = nzb_s_inner(j,i)+1, nzt
+
+          IF ( qc_1d(k) > eps_sb )  THEN
 !
 !--          Intern time scale of coagulation (Seifert and Beheng, 2006):
-!--          (1.0 - ql(k,j,i) / ( ql(k,j,i) + qr(k,j,i) ))
-             tau_cloud = 1.0 - ql(k,j,i) / ( ql(k,j,i) + qr(k,j,i) + 1.0E-20 )
+!--          (1.0 - qc(k,j,i) / ( qc(k,j,i) + qr_1d(k) ))
+             tau_cloud = 1.0 - qc_1d(k) / ( qr_1d(k) + qc_1d(k) )
 !
 !--          Universal function for autoconversion process 
@@ -335,8 +404,8 @@
 !--          Shape parameter of gamma distribution (Geoffroy et al., 2010):
 !--          (Use constant nu_c = 1.0 instead?)
-             nu_c      = 1.0 !MAX( 0.0, 1580.0 * hyrho(k) * ql(k,j,i) - 0.28 )
+             nu_c      = 1.0 !MAX( 0.0, 1580.0 * hyrho(k) * qc(k,j,i) - 0.28 )
 !
 !--          Mean weight of cloud droplets:
-             xc        = hyrho(k) * ql(k,j,i) / nc
+             xc = hyrho(k) * qc_1d(k) / nc_1d(k)
 !
 !--          Parameterized turbulence effects on autoconversion (Seifert, 
@@ -371,16 +440,13 @@
 !
 !--          Autoconversion rate (Seifert and Beheng, 2006):
-             autocon   = k_au * ( nu_c + 2.0 ) * ( nu_c + 4.0 ) /              &
-                         ( nu_c + 1.0 )**2 * ql(k,j,i)**2 * xc**2 *            &
-                         ( 1.0 + phi_au / ( 1.0 - tau_cloud )**2 ) *           &
-                         rho_surface
-             autocon   = MIN( autocon, ql(k,j,i) / ( dt_3d *                   &
-                              weight_substep(intermediate_timestep_count) ) )
-!
-!--          Tendencies for q, qr, nr, pt:
-             tend_qr(k,j,i) = tend_qr(k,j,i) + autocon
-             tend_q(k,j,i)  = tend_q(k,j,i)  - autocon
-             tend_nr(k,j,i) = tend_nr(k,j,i) + autocon / x0 * hyrho(k)
-             tend_pt(k,j,i) = tend_pt(k,j,i) + autocon * l_d_cp * pt_d_t(k)
+             autocon = k_au * ( nu_c + 2.0 ) * ( nu_c + 4.0 ) /    &
+                       ( nu_c + 1.0 )**2 * qc_1d(k)**2 * xc**2 *   &
+                       ( 1.0 + phi_au / ( 1.0 - tau_cloud )**2 ) * &
+                       rho_surface
+             autocon = MIN( autocon, qc_1d(k) / dt_micro )
+
+             qr_1d(k) = qr_1d(k) + autocon * dt_micro
+             qc_1d(k) = qc_1d(k) - autocon * dt_micro 
+             nr_1d(k) = nr_1d(k) + autocon / x0 * hyrho(k) * dt_micro
 
           ENDIF
@@ -395,20 +461,17 @@
        USE arrays_3d
        USE cloud_parameters
-       USE constants
-       USE indices
-       USE control_parameters
-       USE statistics
-       
-       IMPLICIT NONE
-
-       INTEGER ::  i, j, k
-       REAL    ::  accr, phi_ac, tau_cloud, k_cr
-
-       DO  k = nzb_2d(j,i)+1, nzt
-
-          IF ( ( ql(k,j,i) > 0.0 )  .AND.  ( qr(k,j,i) > eps_sb ) )  THEN
+       USE control_parameters
+       USE indices
+
+       IMPLICIT NONE
+
+       INTEGER ::  i, j, k
+       REAL    ::  accr, k_cr, phi_ac, tau_cloud, xc
+
+       DO  k = nzb_s_inner(j,i)+1, nzt
+          IF ( ( qc_1d(k) > eps_sb )  .AND.  ( qr_1d(k) > eps_sb ) )  THEN
 !
 !--          Intern time scale of coagulation (Seifert and Beheng, 2006):
-             tau_cloud = 1.0 - ql(k,j,i) / ( ql(k,j,i) + qr(k,j,i) + 1.0E-20) 
+             tau_cloud = 1.0 - qc_1d(k) / ( qc_1d(k) + qr_1d(k) ) 
 !
 !--          Universal function for accretion process 
@@ -421,5 +484,5 @@
 !--          convert the dissipation (diss) from m2 s-3 to cm2 s-3.
              IF ( turbulence )  THEN
-                k_cr = k_cr0 * ( 1.0 + 0.05 *                            &
+                k_cr = k_cr0 * ( 1.0 + 0.05 *                             &
                                  MIN( 600.0, diss(k,j,i) * 1.0E4 )**0.25 )
              ELSE
@@ -428,13 +491,10 @@
 !
 !--          Accretion rate (Seifert and Beheng, 2006):
-             accr = k_cr * ql(k,j,i) * qr(k,j,i) * phi_ac *              &
+             accr = k_cr * qc_1d(k) * qr_1d(k) * phi_ac *                 &
                     SQRT( rho_surface * hyrho(k) )
-             accr = MIN( accr, ql(k,j,i) / ( dt_3d *                     &
-                    weight_substep(intermediate_timestep_count) ) )
-!
-!--          Tendencies for q, qr, pt:
-             tend_qr(k,j,i) = tend_qr(k,j,i) + accr
-             tend_q(k,j,i)  = tend_q(k,j,i)  - accr
-             tend_pt(k,j,i) = tend_pt(k,j,i) + accr * l_d_cp * pt_d_t(k) 
+             accr = MIN( accr, qc_1d(k) / dt_micro )
+
+             qr_1d(k) = qr_1d(k) + accr * dt_micro 
+             qc_1d(k) = qc_1d(k) - accr * dt_micro
 
           ENDIF
@@ -449,31 +509,25 @@
        USE arrays_3d
        USE cloud_parameters
-       USE constants
-       USE indices
-       USE control_parameters
-       USE statistics
+       USE control_parameters
+       USE indices
    
        IMPLICIT NONE
 
        INTEGER ::  i, j, k
-       REAL    ::  selfcoll, breakup, phi_br, phi_sc
-
-
-       DO  k = nzb_2d(j,i)+1, nzt
-
-          IF ( qr(k,j,i) > eps_sb )  THEN
-!
-!--          Selfcollection rate (Seifert and Beheng, 2006):
-!--          pirho_l**( 1.0 / 3.0 ) is necessary to convert [lambda_r] = m-1 to 
-!--          kg**( 1.0 / 3.0 ).
-             phi_sc = 1.0 !( 1.0 + kappa_rr / lambda_r(k) *                        &
-                      !pirho_l**( 1.0 / 3.0 ) )**( -9 )
-
-             selfcoll = k_rr * nr(k,j,i) * qr(k,j,i) * phi_sc *               &
+       REAL    ::  breakup, dr, phi_br, selfcoll
+
+       DO  k = nzb_s_inner(j,i)+1, nzt
+          IF ( qr_1d(k) > eps_sb )  THEN
+!
+!--          Selfcollection rate (Seifert and Beheng, 2001):
+             selfcoll = k_rr * nr_1d(k) * qr_1d(k) *         &
                         SQRT( hyrho(k) * rho_surface )
 !
+!--          Weight averaged diameter of rain drops:
+             dr = ( hyrho(k) * qr_1d(k) / nr_1d(k) * dpirho_l )**( 1.0 / 3.0 )
+!
 !--          Collisional breakup rate (Seifert, 2008):
-             IF ( dr(k) >= 0.3E-3 )  THEN
-                phi_br  = k_br * ( dr(k) - 1.1E-3 )
+             IF ( dr >= 0.3E-3 )  THEN
+                phi_br  = k_br * ( dr - 1.1E-3 )
                 breakup = selfcoll * ( phi_br + 1.0 )
              ELSE
@@ -481,12 +535,8 @@
              ENDIF
 
-             selfcoll =  MAX( breakup - selfcoll, -nr(k,j,i) / ( dt_3d *      &
-                              weight_substep(intermediate_timestep_count) ) )
-!
-!--          Tendency for nr:
-             tend_nr(k,j,i) = tend_nr(k,j,i) + selfcoll
+             selfcoll = MAX( breakup - selfcoll, -nr_1d(k) / dt_micro )
+             nr_1d(k) = nr_1d(k) + selfcoll * dt_micro
 
           ENDIF          
-
        ENDDO
 
@@ -502,21 +552,18 @@
        USE cloud_parameters
        USE constants
-       USE indices
-       USE control_parameters
-       USE statistics
-
-       IMPLICIT NONE
-
-       INTEGER ::  i, j, k
-       REAL    ::  evap, alpha, e_s, q_s, t_l, sat, temp, g_evap, f_vent, &
-                   mu_r_2, mu_r_5d2, nr_0
-
-
-       DO  k = nzb_2d(j,i)+1, nzt
-
-          IF ( qr(k,j,i) > eps_sb )  THEN
+       USE control_parameters
+       USE indices
+
+       IMPLICIT NONE
+
+       INTEGER ::  i, j, k
+       REAL    ::  alpha, dr, e_s, evap, evap_nr, f_vent, g_evap, lambda_r, &
+                   mu_r, mu_r_2, mu_r_5d2, nr_0, q_s, sat, t_l, temp, xr
+
+       DO  k = nzb_s_inner(j,i)+1, nzt
+          IF ( qr_1d(k) > eps_sb )  THEN
 !
 !--          Actual liquid water temperature:
-             t_l = t_d_pt(k) * pt(k,j,i)
+             t_l = t_d_pt(k) * pt_1d(k)
 !
 !--          Saturation vapor pressure at t_l:
@@ -526,45 +573,58 @@
              q_s = 0.622 * e_s / ( hyp(k) - 0.378 * e_s )
              alpha = 0.622 * l_d_r * l_d_cp / ( t_l * t_l )
-             q_s = q_s * ( 1.0 + alpha * q(k,j,i) ) / ( 1.0 + alpha * q_s )
+             q_s = q_s * ( 1.0 + alpha * q_1d(k) ) / ( 1.0 + alpha * q_s )
 !
 !--          Supersaturation:
-             sat = MIN( 0.0, ( q(k,j,i) - ql(k,j,i) ) / q_s - 1.0 )
+             sat = MIN( 0.0, ( q_1d(k) - qr_1d(k) - qc_1d(k) ) / q_s - 1.0 )
 !
 !--          Actual temperature:
-             temp = t_l + l_d_cp * ql(k,j,i)
-!
-!--          
-             g_evap = ( l_v / ( r_v * temp ) - 1.0 ) * l_v /                  &
-                      ( thermal_conductivity_l * temp ) + rho_l * r_v * temp /&
-                      ( diff_coeff_l * e_s )
-             g_evap = 1.0 / g_evap
+             temp = t_l + l_d_cp * ( qc_1d(k) + qr_1d(k) )
+    
+             g_evap = 1.0 / ( ( l_v / ( r_v * temp ) - 1.0 ) * l_v /   &
+                      ( thermal_conductivity_l * temp ) + r_v * temp / &
+                      ( diff_coeff_l * e_s ) )
+!
+!--          Mean weight of rain drops
+             xr = hyrho(k) * qr_1d(k) / nr_1d(k)
+!
+!--          Weight averaged diameter of rain drops:
+             dr = ( xr * dpirho_l )**( 1.0 / 3.0 )
 !
 !--          Compute ventilation factor and intercept parameter 
 !--          (Seifert and Beheng, 2006; Seifert, 2008):
              IF ( ventilation_effect )  THEN
-                mu_r_2   = mu_r(k) + 2.0
-                mu_r_5d2 = mu_r(k) + 2.5 
+!
+!--             Shape parameter of gamma distribution (Milbrandt and Yau, 2005;
+!--             Stevens and Seifert, 2008):
+                mu_r = 10.0 * ( 1.0 + TANH( 1.2E3 * ( dr - 1.4E-3 ) ) )
+!
+!--             Slope parameter of gamma distribution (Seifert, 2008):
+                lambda_r = ( ( mu_r + 3.0 ) * ( mu_r + 2.0 ) *       &
+                             ( mu_r + 1.0 ) )**( 1.0 / 3.0 ) / dr
+
+                mu_r_2   = mu_r + 2.0
+                mu_r_5d2 = mu_r + 2.5 
                 f_vent = a_vent * gamm( mu_r_2 ) *                            &
-                         lambda_r(k)**( -mu_r_2 ) +                           &
+                         lambda_r**( -mu_r_2 ) +                              &
                          b_vent * schmidt_p_1d3 *                             &
                          SQRT( a_term / kin_vis_air ) * gamm( mu_r_5d2 ) *    &
-                         lambda_r(k)**( -mu_r_5d2 ) *                         &
+                         lambda_r**( -mu_r_5d2 ) *                            &
                          ( 1.0 - 0.5 * ( b_term / a_term ) *                  &
-                         ( lambda_r(k) /                                      &
-                         (       c_term + lambda_r(k) ) )**mu_r_5d2 -         &
+                         ( lambda_r /                                         &
+                         (       c_term + lambda_r ) )**mu_r_5d2 -            &
                                  0.125 * ( b_term / a_term )**2 *             &
-                         ( lambda_r(k) /                                      &
-                         ( 2.0 * c_term + lambda_r(k) ) )**mu_r_5d2 -         &
+                         ( lambda_r /                                         &
+                         ( 2.0 * c_term + lambda_r ) )**mu_r_5d2 -            &
                                  0.0625 * ( b_term / a_term )**3 *            &
-                         ( lambda_r(k) /                                      &
-                         ( 3.0 * c_term + lambda_r(k) ) )**mu_r_5d2 -         &
+                         ( lambda_r /                                         &
+                         ( 3.0 * c_term + lambda_r ) )**mu_r_5d2 -            &
                                  0.0390625 * ( b_term / a_term )**4 *         &
-                         ( lambda_r(k) /                                      &
-                         ( 4.0 * c_term + lambda_r(k) ) )**mu_r_5d2 )
-                nr_0   = nr(k,j,i) * lambda_r(k)**( mu_r(k) + 1.0 ) /         &
-                         gamm( mu_r(k) + 1.0 ) 
+                         ( lambda_r /                                         &
+                         ( 4.0 * c_term + lambda_r ) )**mu_r_5d2 )
+                nr_0   = nr_1d(k) * lambda_r**( mu_r + 1.0 ) /                &
+                         gamm( mu_r + 1.0 ) 
              ELSE
                 f_vent = 1.0
-                nr_0   = nr(k,j,i) * dr(k)
+                nr_0   = nr_1d(k) * dr
              ENDIF
 !
@@ -572,13 +632,11 @@
              evap = 2.0 * pi * nr_0 * g_evap * f_vent * sat /    &
                     hyrho(k)
-             evap = MAX( evap, -qr(k,j,i) / ( dt_3d *            &
-                         weight_substep(intermediate_timestep_count) ) )
-!
-!--          Tendencies for q, qr, nr, pt:
-             tend_qr(k,j,i) = tend_qr(k,j,i) + evap
-             tend_q(k,j,i)  = tend_q(k,j,i)  - evap
-             tend_nr(k,j,i) = tend_nr(k,j,i) + c_evap * evap / xr(k) * hyrho(k)
-             tend_pt(k,j,i) = tend_pt(k,j,i) + evap * l_d_cp * pt_d_t(k) 
-
+
+             evap    = MAX( evap, -qr_1d(k) / dt_micro )
+             evap_nr = MAX( c_evap * evap / xr * hyrho(k), &
+                            -nr_1d(k) / dt_micro )
+
+             qr_1d(k) = qr_1d(k) + evap * dt_micro
+             nr_1d(k) = nr_1d(k) + evap_nr * dt_micro
           ENDIF          
 
@@ -593,32 +651,39 @@
        USE cloud_parameters
        USE constants
-       USE indices
-       USE control_parameters
+       USE control_parameters
+       USE indices
        
        IMPLICIT NONE
 
        INTEGER ::  i, j, k
-       REAL    ::  sed_q_const, sigma_gc = 1.3, k_st = 1.2E8
+       REAL    ::  sed_qc_const
+
+       REAL, DIMENSION(nzb:nzt+1) :: sed_qc
 
 !
 !--    Sedimentation of cloud droplets (Heus et al., 2010):
-       sed_q_const = k_st * ( 3.0 / ( 4.0 * pi * rho_l ))**( 2.0 / 3.0 ) *    &
+       sed_qc_const = k_st * ( 3.0 / ( 4.0 * pi * rho_l ))**( 2.0 / 3.0 ) *   &
                      EXP( 5.0 * LOG( sigma_gc )**2 )
 
-       sed_q = 0.0
-
-       DO  k = nzb_2d(j,i)+1, nzt
-          IF ( ql(k,j,i) > 0.0 )  THEN
-             sed_q(k) = sed_q_const * nc**( -2.0 / 3.0 ) *                    &
-                        ( ql(k,j,i) * hyrho(k) )**( 5.0 / 3.0 )
+       sed_qc(nzt+1) = 0.0
+
+       DO  k = nzt, nzb_s_inner(j,i)+1, -1
+          IF ( qc_1d(k) > eps_sb )  THEN
+             sed_qc(k) = sed_qc_const * nc_1d(k)**( -2.0 / 3.0 ) * &
+                        ( qc_1d(k) * hyrho(k) )**( 5.0 / 3.0 )
+          ELSE
+             sed_qc(k) = 0.0
           ENDIF
-       ENDDO
-!
-!--    Tendency for q, pt:
-       DO  k = nzb_2d(j,i)+1, nzt
-          tend_q(k,j,i)  = tend_q(k,j,i) + ( sed_q(k+1) - sed_q(k) ) *        &
-                           ddzu(k+1) / hyrho(k) 
-          tend_pt(k,j,i) = tend_pt(k,j,i) - ( sed_q(k+1) - sed_q(k) ) *       &
-                           ddzu(k+1) / hyrho(k) * l_d_cp * pt_d_t(k) 
+
+          sed_qc(k) = MIN( sed_qc(k), hyrho(k) * dzu(k+1) * q_1d(k) /     &
+                                      dt_micro + sed_qc(k+1) )
+
+          q_1d(k)  = q_1d(k)  + ( sed_qc(k+1) - sed_qc(k) ) * ddzu(k+1) /  &
+                                hyrho(k) * dt_micro
+          qc_1d(k) = qc_1d(k) + ( sed_qc(k+1) - sed_qc(k) ) * ddzu(k+1) / & 
+                                hyrho(k) * dt_micro
+          pt_1d(k) = pt_1d(k) - ( sed_qc(k+1) - sed_qc(k) ) * ddzu(k+1) / &
+                                hyrho(k) * l_d_cp * pt_d_t(k) * dt_micro
+
        ENDDO
 
@@ -631,6 +696,6 @@
        USE cloud_parameters
        USE constants
-       USE indices
-       USE control_parameters
+       USE control_parameters
+       USE indices
        USE statistics
        
@@ -638,25 +703,33 @@
 
        INTEGER ::  i, j, k, k_run
-       REAL    ::  c_run, d_max, d_mean, d_min, dt_sedi, flux, z_run
-
-       REAL, DIMENSION(nzb:nzt) :: c_nr, c_qr, nr_slope, qr_slope, w_nr, w_qr
-
+       REAL    ::  c_run, d_max, d_mean, d_min, dr, dt_sedi, flux, lambda_r,  &
+                   mu_r, z_run
+
+       REAL, DIMENSION(nzb:nzt+1) :: c_nr, c_qr, d_nr, d_qr, nr_slope,        &
+                                     qr_slope, sed_nr, sed_qr, w_nr, w_qr 
 !
 !--    Computation of sedimentation flux. Implementation according to Stevens
 !--    and Seifert (2008).
        IF ( intermediate_timestep_count == 1 )  prr(:,j,i) = 0.0
-
-       dt_sedi = dt_3d * weight_substep(intermediate_timestep_count)
-
-       w_nr = 0.0
-       w_qr = 0.0
 !
 !--    Compute velocities 
        DO  k = nzb_s_inner(j,i)+1, nzt
-          IF ( qr(k,j,i) > eps_sb )  THEN
+          IF ( qr_1d(k) > eps_sb )  THEN
+!
+!--          Weight averaged diameter of rain drops:
+             dr = ( hyrho(k) * qr_1d(k) / nr_1d(k) * dpirho_l )**( 1.0 / 3.0 )
+!
+!--          Shape parameter of gamma distribution (Milbrandt and Yau, 2005;
+!--          Stevens and Seifert, 2008):
+             mu_r = 10.0 * ( 1.0 + TANH( 1.2E3 * ( dr - 1.4E-3 ) ) )
+!
+!--          Slope parameter of gamma distribution (Seifert, 2008):
+             lambda_r = ( ( mu_r + 3.0 ) * ( mu_r + 2.0 ) *          &
+                        ( mu_r + 1.0 ) )**( 1.0 / 3.0 ) / dr
+
              w_nr(k) = MAX( 0.1, MIN( 20.0, a_term - b_term * ( 1.0 +          &
-                       c_term / lambda_r(k) )**( -1.0 * ( mu_r(k) + 1.0 ) ) ) )
+                       c_term / lambda_r )**( -1.0 * ( mu_r + 1.0 ) ) ) )
              w_qr(k) = MAX( 0.1, MIN( 20.0, a_term - b_term * ( 1.0 +          &
-                       c_term / lambda_r(k) )**( -1.0 * ( mu_r(k) + 4.0 ) ) ) )
+                       c_term / lambda_r )**( -1.0 * ( mu_r + 4.0 ) ) ) )
           ELSE
              w_nr(k) = 0.0
@@ -666,36 +739,31 @@
 !
 !--    Adjust boundary values
-       w_nr(nzb_2d(j,i)) = w_nr(nzb_2d(j,i)+1)
-       w_qr(nzb_2d(j,i)) = w_qr(nzb_2d(j,i)+1)
-       w_nr(nzt) = w_nr(nzt-1)
-       w_qr(nzt) = w_qr(nzt-1)
+       w_nr(nzb_s_inner(j,i)) = w_nr(nzb_s_inner(j,i)+1)
+       w_qr(nzb_s_inner(j,i)) = w_qr(nzb_s_inner(j,i)+1)
+       w_nr(nzt+1) = 0.0
+       w_qr(nzt+1) = 0.0
 !
 !--    Compute Courant number
-       DO  k = nzb_s_inner(j,i)+1, nzt-1
+       DO  k = nzb_s_inner(j,i)+1, nzt
           c_nr(k) = 0.25 * ( w_nr(k-1) + 2.0 * w_nr(k) + w_nr(k+1) ) * &
-                    dt_sedi * ddzu(k)
+                    dt_micro * ddzu(k)
           c_qr(k) = 0.25 * ( w_qr(k-1) + 2.0 * w_qr(k) + w_qr(k+1) ) * &
-                    dt_sedi * ddzu(k)
-       ENDDO
+                    dt_micro * ddzu(k)
+       ENDDO     
 !
 !--    Limit slopes with monotonized centered (MC) limiter (van Leer, 1977):
        IF ( limiter_sedimentation )  THEN
 
-          qr(nzb_s_inner(j,i),j,i) = 0.0
-          nr(nzb_s_inner(j,i),j,i) = 0.0
-          qr(nzt,j,i) = 0.0
-          nr(nzt,j,i) = 0.0
-
-          DO k = nzb_s_inner(j,i)+1, nzt-1
-             d_mean = 0.5 * ( qr(k+1,j,i) + qr(k-1,j,i) )
-             d_min  = qr(k,j,i) - MIN( qr(k+1,j,i), qr(k,j,i), qr(k-1,j,i) )
-             d_max  = MAX( qr(k+1,j,i), qr(k,j,i), qr(k-1,j,i) ) - qr(k,j,i)
+          DO k = nzb_s_inner(j,i)+1, nzt
+             d_mean = 0.5 * ( qr_1d(k+1) + qr_1d(k-1) )
+             d_min  = qr_1d(k) - MIN( qr_1d(k+1), qr_1d(k), qr_1d(k-1) )
+             d_max  = MAX( qr_1d(k+1), qr_1d(k), qr_1d(k-1) ) - qr_1d(k)
 
              qr_slope(k) = SIGN(1.0, d_mean) * MIN ( 2.0 * d_min, 2.0 * d_max, &
                                                      ABS( d_mean ) )
 
-             d_mean = 0.5 * ( nr(k+1,j,i) + nr(k-1,j,i) )
-             d_min  = nr(k,j,i) - MIN( nr(k+1,j,i), nr(k,j,i), nr(k-1,j,i) )
-             d_max  = MAX( nr(k+1,j,i), nr(k,j,i), nr(k-1,j,i) ) - nr(k,j,i)
+             d_mean = 0.5 * ( nr_1d(k+1) + nr_1d(k-1) )
+             d_min  = nr_1d(k) - MIN( nr_1d(k+1), nr_1d(k), nr_1d(k-1) )
+             d_max  = MAX( nr_1d(k+1), nr_1d(k), nr_1d(k-1) ) - nr_1d(k)
 
              nr_slope(k) = SIGN(1.0, d_mean) * MIN ( 2.0 * d_min, 2.0 * d_max, &
@@ -709,7 +777,10 @@
 
        ENDIF
+
+       sed_nr(nzt+1) = 0.0
+       sed_qr(nzt+1) = 0.0
 !
 !--    Compute sedimentation flux
-       DO  k = nzt-2, nzb_s_inner(j,i)+1, -1
+       DO  k = nzt, nzb_s_inner(j,i)+1, -1
 !
 !--       Sum up all rain drop number densities which contribute to the flux 
@@ -719,7 +790,7 @@
           k_run = k
           c_run = MIN( 1.0, c_nr(k) )
-          DO WHILE ( c_run > 0.0  .AND.  k_run <= nzt-1 )
+          DO WHILE ( c_run > 0.0  .AND.  k_run <= nzt )
              flux  = flux + hyrho(k_run) *                                    &
-                     ( nr(k_run,j,i) + nr_slope(k_run) * ( 1.0 - c_run ) *    &
+                     ( nr_1d(k_run) + nr_slope(k_run) * ( 1.0 - c_run ) *     &
                      0.5 ) * c_run * dzu(k_run)
              z_run = z_run + dzu(k_run)
@@ -731,9 +802,9 @@
 !--       available
           flux = MIN( flux,                                                   &
-                      hyrho(k) * dzu(k) * nr(k,j,i) + sed_nr(k+1) * dt_sedi )
-
-          sed_nr(k)      = flux / dt_sedi
-          tend_nr(k,j,i) = tend_nr(k,j,i) + ( sed_nr(k+1) - sed_nr(k) ) *     &
-                                            ddzu(k+1) / hyrho(k)
+                      hyrho(k) * dzu(k+1) * nr_1d(k) + sed_nr(k+1) * dt_micro )
+
+          sed_nr(k) = flux / dt_micro
+          nr_1d(k)  = nr_1d(k) + ( sed_nr(k+1) - sed_nr(k) ) * ddzu(k+1) /    &
+                                 hyrho(k) * dt_micro
 !
 !--       Sum up all rain water content which contributes to the flux 
@@ -747,5 +818,5 @@
 
              flux  = flux + hyrho(k_run) *                                    &
-                     ( qr(k_run,j,i) + qr_slope(k_run) * ( 1.0 - c_run ) *    &
+                     ( qr_1d(k_run) + qr_slope(k_run) * ( 1.0 - c_run ) *    &
                      0.5 ) * c_run * dzu(k_run)
              z_run = z_run + dzu(k_run)
@@ -757,14 +828,20 @@
 !--       It is not allowed to sediment more rain water content than available
           flux = MIN( flux,                                                   &
-                      hyrho(k) * dzu(k) * qr(k,j,i) + sed_qr(k+1) * dt_sedi )
-
-          sed_qr(k)      = flux / dt_sedi
-          tend_qr(k,j,i) = tend_qr(k,j,i) + ( sed_qr(k+1) - sed_qr(k) ) *     &
-                                            ddzu(k+1) / hyrho(k)
+                      hyrho(k) * dzu(k) * qr_1d(k) + sed_qr(k+1) * dt_micro )
+
+          sed_qr(k) = flux / dt_micro
+
+          qr_1d(k) = qr_1d(k) + ( sed_qr(k+1) - sed_qr(k) ) * ddzu(k+1) / &
+                                hyrho(k) * dt_micro
+          q_1d(k)  = q_1d(k)  + ( sed_qr(k+1) - sed_qr(k) ) * ddzu(k+1) / &
+                                hyrho(k) * dt_micro 
+          pt_1d(k) = pt_1d(k) - ( sed_qr(k+1) - sed_qr(k) ) * ddzu(k+1) / &
+                                hyrho(k) * l_d_cp * pt_d_t(k) * dt_micro
 !
 !--       Compute the rain rate
           prr(k,j,i) = prr(k,j,i) + sed_qr(k) / hyrho(k) *                    &
-                       weight_substep(intermediate_timestep_count) 
-       ENDDO
+                       weight_substep(intermediate_timestep_count)
+       ENDDO
+
 !
 !--    Precipitation amount
@@ -774,6 +851,6 @@
 
           precipitation_amount(j,i) = precipitation_amount(j,i) +   &
-                                      prr(nzb_2d(j,i)+1,j,i) *      &
-                                      hyrho(nzb_2d(j,i)+1) * dt_3d
+                                      prr(nzb_s_inner(j,i)+1,j,i) *      &
+                                      hyrho(nzb_s_inner(j,i)+1) * dt_3d
        ENDIF
 
Index: /palm/trunk/SOURCE/modules.f90
===================================================================
--- /palm/trunk/SOURCE/modules.f90	(revision 1114)
+++ /palm/trunk/SOURCE/modules.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! unused variables removed
 !
 ! Former revisions:
@@ -408,7 +408,7 @@
           c_u_m, c_u_m_l, c_v_m, c_v_m_l, c_w_m, c_w_m_l, ddzu, ddzu_pres,     &
           dd2zu, dzu, ddzw, dzw, hyp, inflow_damping_factor, lad, l_grid,      &
-          nr_init, ptdf_x, ptdf_y, pt_init, q_init, qr_init, rdf, rdf_sc,      &
-          sa_init, ug, u_init, u_nzb_p1_for_vfc, vg, v_init, v_nzb_p1_for_vfc, &
-          w_subs, zu, zw
+          nc_1d, nr_1d, ptdf_x, ptdf_y, pt_1d, pt_init, q_1d, q_init, qc_1d,   &
+          qr_1d, rdf, rdf_sc, sa_init, ug, u_init, u_nzb_p1_for_vfc, vg,       &
+          v_init, v_nzb_p1_for_vfc, w_subs, zu, zw
 
     REAL, DIMENSION(:,:), ALLOCATABLE ::                                       &
@@ -430,21 +430,20 @@
           tend_nr, tend_q, tend_qr, tric, 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
-            
 
 #if defined( __nopointer )
     REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET ::                             &
-          e, e_p, nr, nr_p, p, prho, pt, pt_p, q, q_p, ql, ql_c, ql_v, ql_vp,  &
-          qr, qr_p, rho, sa, sa_p, te_m, tnr_m, tpt_m, tq_m, tqr_m, tsa_m,     &
-          tu_m, tv_m, tw_m, u, u_p, v, v_p, vpt, w, w_p
+          e, e_p, nr, nr_p, p, prho, pt, pt_p, q, q_p, qc, ql, ql_c, ql_v,     &
+          ql_vp, qr, qr_p, rho, sa, sa_p, te_m, tnr_m, tpt_m, tq_m, tqr_m,     &
+          tsa_m, tu_m, tv_m, tw_m, u, u_p, v, v_p, vpt, w, w_p
 #else
     REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET ::                             &
           e_1, e_2, e_3, p, prho_1, nr_1, nr_2, nr_3, pt_1, pt_2, pt_3, q_1,   &
-          q_2, q_3, ql_v, ql_vp, ql_1, ql_2, qr_1, qr_2, qr_3, rho_1, sa_1,    &
-          sa_2, sa_3, u_1, u_2, u_3, v_1, v_2, v_3, vpt_1, w_1, w_2, w_3
+          q_2, q_3, qc_1, ql_v, ql_vp, ql_1, ql_2, qr_1, qr_2, qr_3, rho_1,    &
+          sa_1, sa_2, sa_3, u_1, u_2, u_3, v_1, v_2, v_3, vpt_1, w_1, w_2, w_3
 
     REAL, DIMENSION(:,:,:), POINTER ::                                         &
-          e, e_p, nr, nr_p, prho, pt, pt_p, q, q_p, ql, ql_c, qr, qr_p, rho,   &
-          sa, sa_p, te_m, tnr_m, tpt_m, tq_m, tqr_m, tsa_m, tu_m, tv_m, tw_m,  &
-          u, u_p, v, v_p, vpt, w, w_p
+          e, e_p, nr, nr_p, prho, pt, pt_p, q, q_p, qc, ql, ql_c, qr, qr_p,    &
+          rho, sa, sa_p, te_m, tnr_m, tpt_m, tq_m, tqr_m, tsa_m, tu_m, tv_m,   &
+          tw_m, u, u_p, v, v_p, vpt, w, w_p
 #endif
 
@@ -474,7 +473,7 @@
 
     REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: &
-          e_av, lpt_av, nr_av, p_av, pc_av, pr_av, prr_av, pt_av, q_av, ql_av, &
-          ql_c_av, ql_v_av, ql_vp_av, qr_av, qv_av, rho_av, s_av, sa_av, u_av, &
-          v_av, vpt_av, w_av
+          e_av, lpt_av, nr_av, p_av, pc_av, pr_av, prr_av, pt_av, q_av, qc_av, &
+          ql_av, ql_c_av, ql_v_av, ql_vp_av, qr_av, qv_av, rho_av, s_av, sa_av,&
+          u_av, v_av, vpt_av, w_av
   
  END MODULE averaging
@@ -525,8 +524,9 @@
              eps_ros = 1.0E-4,  & ! accuracy of Rosenbrock method
              eps_sb = 1.0E-20,  & ! threshold in two-moments scheme
-             k_cc = 9.44E09,    & ! const. rain-rain kernel (m3 kg-2 s-1) 
+             k_cc = 9.44E09,    & ! const. cloud-cloud kernel (m3 kg-2 s-1) 
              k_cr0 = 4.33,      & ! const. cloud-rain kernel (m3 kg-1 s-1) 
              k_rr = 7.12,       & ! const. rain-rain kernel (m3 kg-1 s-1) 
              k_br = 1000.,      & ! const. in breakup parametrization (m-1)
+             k_st = 1.2E8,      & ! const. in drizzle parametrization (m-1 s-1)
              kappa_rr = 60.7,   & ! const. in collision kernel (kg-1/3)
              kin_vis_air = 1.4086E-5, & ! kin. viscosity of air (m2 s-1)
@@ -536,5 +536,5 @@
              molecular_weight_of_solute = 0.05844, & ! mol. m. NaCl (kg mol-1)
              molecular_weight_of_water = 0.01801528, & ! mol. m. H2O (kg mol-1)
-             nc = 70.0E6,       & ! cloud droplet concentration 
+             nc_const = 70.0E6, & ! cloud droplet concentration
              prec_time_const = 0.001, & !coef. in Kessler scheme
              pirho_l, dpirho_l, & ! pi * rho_l / 6.0; 6.0 / ( pi * rho_l )
@@ -545,4 +545,5 @@
              schmidt = 0.71,    & ! Schmidt number
              schmidt_p_1d3,     & ! schmidt**( 1.0 / 3.0 )
+             sigma_gc = 1.3,    & ! log-normal geometric standard deviation
              stp = 2.5066282746310005, & ! parameter in gamma function
              thermal_conductivity_l = 2.43E-2, & ! therm. cond. air (J m-1 s-1 K-1)
@@ -554,6 +555,5 @@
              w_precipitation = 9.65      ! maximum terminal velocity (m s-1)
 
-    REAL, DIMENSION(:), ALLOCATABLE     ::  dr, hyrho, lambda_r, mu_r, pt_d_t, &
-                                            sed_nr, sed_q, sed_qr, t_d_pt, xr  
+    REAL, DIMENSION(:), ALLOCATABLE     ::  hyrho, pt_d_t, t_d_pt  
 
     REAL, DIMENSION(:,:), ALLOCATABLE   ::  precipitation_amount, &
@@ -630,10 +630,8 @@
                              scalar_advec = 'ws-scheme'
     CHARACTER (LEN=20)   ::  bc_e_b = 'neumann', bc_lr = 'cyclic', &
-                             bc_nr_b = 'neumann', bc_nr_t = 'neumann', &
                              bc_ns = 'cyclic', bc_p_b = 'neumann', &
                              bc_p_t = 'dirichlet', bc_pt_b = 'dirichlet', &
                              bc_pt_t = 'initial_gradient', &
                              bc_q_b = 'dirichlet', bc_q_t = 'neumann', &
-                             bc_qr_b = 'neumann', bc_qr_t = 'neumann',&
                              bc_s_b = 'dirichlet', bc_s_t = 'neumann', &
                              bc_sa_t = 'neumann', &
@@ -681,6 +679,6 @@
                 dp_level_ind_b = 0, dvrp_filecount = 0, &
                 dz_stretch_level_index, gamma_mg, gathered_size, &
-                grid_level, ibc_e_b, ibc_nr_b, ibc_nr_t, ibc_p_b, ibc_p_t, &
-                ibc_pt_b, ibc_pt_t, ibc_q_b, ibc_q_t, ibc_qr_b, ibc_qr_t, &
+                grid_level, ibc_e_b, ibc_p_b, ibc_p_t, &
+                ibc_pt_b, ibc_pt_t, ibc_q_b, ibc_q_t, &
                 ibc_sa_t, ibc_uv_b, ibc_uv_t, icloud_scheme, &
                 inflow_disturbance_begin = -1, inflow_disturbance_end = -1, &
@@ -705,8 +703,6 @@
                 mask_size(max_masks,3) = -1, mask_size_l(max_masks,3) = -1, &
                 mask_start_l(max_masks,3) = -1, &
-                nr_vertical_gradient_level_ind(10) = -9999, &
                 pt_vertical_gradient_level_ind(10) = -9999, &
                 q_vertical_gradient_level_ind(10) = -9999, &
-                qr_vertical_gradient_level_ind(10) = -9999, &
                 sa_vertical_gradient_level_ind(10) = -9999, &
                 section(100,3), section_xy(100) = -9999, &
@@ -738,6 +734,5 @@
                 constant_top_momentumflux = .FALSE., &
                 constant_top_salinityflux = .TRUE., &
-                constant_waterflux = .TRUE., constant_waterflux_nr = .TRUE., &
-                constant_waterflux_qr = .TRUE., create_disturbances = .TRUE., &
+                constant_waterflux = .TRUE., create_disturbances = .TRUE., &
                 data_output_2d_on_each_pe = .TRUE., & 
                 dissipation_control = .FALSE., disturbance_created = .FALSE., &
@@ -779,5 +774,5 @@
              averaging_interval = 0.0, averaging_interval_pr = 9999999.9, &
              averaging_interval_sp = 9999999.9, bc_pt_t_val, bc_q_t_val, &
-             bc_qr_t_val, bc_nr_t_val, bottom_salinityflux = 0.0, &
+             bottom_salinityflux = 0.0, &
              building_height = 50.0, building_length_x = 50.0, &
              building_length_y = 50.0, building_wall_left = 9999999.9, &
@@ -799,5 +794,5 @@
              dt_do2d_xy = 9999999.9, dt_do2d_xz = 9999999.9, &
              dt_do2d_yz = 9999999.9, dt_do3d = 9999999.9, dt_dvrp = 9999999.9, &
-             dt_max = 20.0, dt_restart = 9999999.9, &
+             dt_max = 20.0, dt_micro = -1.0, dt_restart = 9999999.9, &
              dt_run_control = 60.0, dt_3d = -1.0, dz = -1.0, &
              dz_max = 9999999.9, dz_stretch_factor = 1.08, &
@@ -810,5 +805,4 @@
              maximum_cpu_time_allowed = 0.0,  &
              molecular_viscosity = 1.461E-5, &
-             nr_surface = 0.0, nr_surface_initial_change = 0.0, &
              old_dt = 1.0E-10, omega = 7.29212E-5, omega_sor = 1.8, &
              particle_maximum_age = 9999999.9, &
@@ -819,5 +813,4 @@
              pt_surface = 300.0, pt_surface_initial_change = 0.0, &
              q_surface = 0.0, q_surface_initial_change = 0.0, &
-             qr_surface = 0.0, qr_surface_initial_change = 0.0, &
              rayleigh_damping_factor = -1.0, rayleigh_damping_height = -1.0, &
              recycling_width = 9999999.9, residual_limit = 1.0E-4, &
@@ -832,5 +825,4 @@
              surface_heatflux = 9999999.9, surface_pressure = 1013.25, &
              surface_scalarflux = 9999999.9, surface_waterflux = 9999999.9, &
-             surface_waterflux_nr = 0.0, surface_waterflux_qr = 0.0, &
              s_surface = 0.0, s_surface_initial_change = 0.0, &
              termination_time_needed = -1.0, time_coupling = 0.0, &
@@ -854,12 +846,8 @@
              lad_vertical_gradient_level(10) = -9999999.9, &
              mask_scale(3), &
-             nr_vertical_gradient(10) = 0.0, &
-             nr_vertical_gradient_level(10) = -1.0, &
              pt_vertical_gradient(10) = 0.0, &
              pt_vertical_gradient_level(10) = -9999999.9, &
              q_vertical_gradient(10) = 0.0, &
              q_vertical_gradient_level(10) = -1.0, &
-             qr_vertical_gradient(10) = 0.0, &
-             qr_vertical_gradient_level(10) = -1.0, &
              s_vertical_gradient(10) = 0.0, &
              s_vertical_gradient_level(10) = -1.0, &
@@ -877,5 +865,5 @@
              volume_flow(1:2) = 0.0, volume_flow_area(1:2) = 0.0, &
              volume_flow_initial(1:2) = 0.0, wall_heatflux(0:4) = 0.0, &
-             wall_humidityflux(0:4) = 0.0, wall_nrflux(0:4) = 0.0, & 
+             wall_humidityflux(0:4) = 0.0, wall_nrflux(0:4) = 0.0, &
              wall_qflux(0:4) = 0.0, wall_qrflux(0:4) = 0.0, &
              wall_salinityflux(0:4) = 0.0, wall_scalarflux(0:4) = 0.0, &
Index: /palm/trunk/SOURCE/parin.f90
===================================================================
--- /palm/trunk/SOURCE/parin.f90	(revision 1114)
+++ /palm/trunk/SOURCE/parin.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! unused variables removed
 !
 ! Former revisions:
@@ -224,7 +224,6 @@
              loop_optimization, masking_method, mg_cycles, &
              mg_switch_to_pe0_level, mixing_length_1d, momentum_advec, &
-             nc, netcdf_precision, neutral, ngsrb, nr_surface, &
-             nr_surface_initial_change, nr_vertical_gradient, &
-             nr_vertical_gradient_level, nsor, &
+             nc_const, netcdf_precision, neutral, ngsrb, &
+             nsor, &
              nsor_ini, nx, ny, nz, ocean, omega, omega_sor, &
              passive_scalar, pch_index, phi, plant_canopy, prandtl_layer, &
@@ -234,6 +233,5 @@
              pt_vertical_gradient_level, q_surface, q_surface_initial_change, &
              q_vertical_gradient, q_vertical_gradient_level, &
-             qr_surface, qr_surface_initial_change, &
-             qr_vertical_gradient, qr_vertical_gradient_level, radiation, &
+             radiation, &
              random_generator, random_heatflux, rayleigh_damping_factor, &
              rayleigh_damping_height, recycling_width, residual_limit, &
@@ -243,6 +241,6 @@
              statistic_regions, subs_vertical_gradient, &
              subs_vertical_gradient_level, surface_heatflux, surface_pressure, &
-             surface_scalarflux, surface_waterflux, surface_waterflux_nr, &
-             surface_waterflux_qr, s_surface, &
+             surface_scalarflux, surface_waterflux, &
+             s_surface, &
              s_surface_initial_change, s_vertical_gradient, &
              s_vertical_gradient_level, timestep_scheme, &
@@ -412,5 +410,4 @@
 !--                  also check the allocate statement in routine read_var_list
              ALLOCATE( lad(0:nz+1),pt_init(0:nz+1), q_init(0:nz+1),           &
-                       qr_init(0:nz+1), nr_init(0:nz+1),                      &
                        sa_init(0:nz+1), ug(0:nz+1), u_init(0:nz+1),           &
                        v_init(0:nz+1), vg(0:nz+1),                            &
Index: /palm/trunk/SOURCE/prognostic_equations.f90
===================================================================
--- /palm/trunk/SOURCE/prognostic_equations.f90	(revision 1114)
+++ /palm/trunk/SOURCE/prognostic_equations.f90	(revision 1115)
@@ -20,5 +20,6 @@
 ! Current revisions:
 ! ------------------
-! 
+! optimized cloud physics: calculation of microphysical tendencies transfered 
+! to microphysics.f90; qr and nr are only calculated if precipitation is required
 !
 ! Former revisions:
@@ -229,5 +230,4 @@
     IF ( ( ws_scheme_mom .OR. ws_scheme_sca )  .AND.  &
          intermediate_timestep_count == 1 )  CALL ws_statistics
-
 !
 !-- Loop over all prognostic equations
@@ -285,5 +285,4 @@
 
              CALL user_actions( i, j, 'u-tendency' )
-
 !
 !--          Prognostic equation for u-velocity component
@@ -341,5 +340,4 @@
 
              CALL user_actions( i, j, 'v-tendency' )
-
 !
 !--          Prognostic equation for v-velocity component
@@ -422,31 +420,8 @@
              ENDIF
           ENDIF
-
-!
-!--       If required, calculate tendencies for total water content, rain water 
-!--       content, rain drop concentration and liquid temperature
-          IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
-
-             tend_q(:,j,i)  = 0.0
-             tend_qr(:,j,i) = 0.0
-             tend_nr(:,j,i) = 0.0
-             tend_pt(:,j,i) = 0.0
-!
-!--          Droplet size distribution (dsd) properties are needed for the
-!--          computation of selfcollection, breakup, evaporation and
-!--          sedimentation of rain
-             IF ( precipitation )  THEN
-                CALL dsd_properties( i,j )
-                CALL autoconversion( i,j )
-                CALL accretion( i,j )
-                CALL selfcollection_breakup( i,j )
-                CALL evaporation_rain( i,j )
-                CALL sedimentation_rain( i,j )
-             ENDIF
-
-             IF ( drizzle )  CALL sedimentation_cloud( i,j )
-
-          ENDIF
-
+!
+!--       If required, calculate tendencies for total water content, liquid water 
+!--       potential temperature, rain water content and rain drop concentration
+          IF ( cloud_physics  .AND.  icloud_scheme == 0 )  CALL microphysics_control( i, j )
 !
 !--       If required, compute prognostic equation for potential temperature
@@ -650,5 +625,6 @@
 !--          If required, calculate prognostic equations for rain water content 
 !--          and rain drop concentration
-             IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.              &
+                  precipitation )  THEN
 !
 !--             Calculate prognostic equation for rain water content
@@ -667,5 +643,4 @@
                 ENDIF
                 CALL diffusion_s( i, j, qr, qrsws, qrswst, wall_qrflux )
-
 !
 !--             Using microphysical tendencies (autoconversion, accretion, 
@@ -673,20 +648,12 @@
                 tend(:,j,i) = tend(:,j,i) + tend_qr(:,j,i)
 
-!
-!--             If required, compute influence of large-scale subsidence/ascent
-                IF ( large_scale_subsidence )  THEN
-                   CALL subsidence( i, j, tend, qr, qr_init )
-                ENDIF
-
-!              CALL user_actions( i, j, 'qr-tendency' )
-
+                CALL user_actions( i, j, 'qr-tendency' )
 !
 !--             Prognostic equation for rain water content
                 DO  k = nzb_s_inner(j,i)+1, nzt
-                   qr_p(k,j,i) = qr(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +    &
-                                                     tsc(3) * tqr_m(k,j,i) )     &
-                                         - tsc(5) * rdf_sc(k) *                  &
-                                        ( qr(k,j,i) - qr_init(k) )
-                   IF ( qr_p(k,j,i) < 0.0 )  qr_p(k,j,i) = 0.1 * qr(k,j,i)
+                   qr_p(k,j,i) = qr(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +  &
+                                                       tsc(3) * tqr_m(k,j,i) ) &
+                                           - tsc(5) * rdf_sc(k) * qr(k,j,i)
+                   IF ( qr_p(k,j,i) < 0.0 )  qr_p(k,j,i) = 0.0
                 ENDDO
 !
@@ -711,7 +678,7 @@
                 IF ( timestep_scheme(1:5) == 'runge' )  THEN
                    IF ( ws_scheme_sca )  THEN
-                      CALL advec_s_ws( i, j, nr, 'nr', flux_s_nr,       & 
-                                       diss_s_nr, flux_l_nr, diss_l_nr, &
-                                       i_omp_start, tn )
+                      CALL advec_s_ws( i, j, nr, 'nr', flux_s_nr,    & 
+                                    diss_s_nr, flux_l_nr, diss_l_nr, &
+                                    i_omp_start, tn )
                    ELSE 
                       CALL advec_s_pw( i, j, nr )
@@ -721,5 +688,5 @@
                 ENDIF
                 CALL diffusion_s( i, j, nr, nrsws, nrswst, wall_nrflux )
-
+!
 !--             Using microphysical tendencies (autoconversion, accretion, 
 !--             selfcollection, breakup, evaporation; 
@@ -727,20 +694,12 @@
                 tend(:,j,i) = tend(:,j,i) + tend_nr(:,j,i)
 
-!
-!--             If required, compute influence of large-scale subsidence/ascent
-                IF ( large_scale_subsidence )  THEN
-                   CALL subsidence( i, j, tend, nr, nr_init )
-                ENDIF
-
-!                CALL user_actions( i, j, 'nr-tendency' )
-
+                CALL user_actions( i, j, 'nr-tendency' )
 !
 !--             Prognostic equation for rain drop concentration
                 DO  k = nzb_s_inner(j,i)+1, nzt
-                   nr_p(k,j,i) = nr(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) + &
-                                                  tsc(3) * tnr_m(k,j,i) )     &
-                                      - tsc(5) * rdf_sc(k) *                  &
-                                        ( nr(k,j,i) - nr_init(k) )
-                   IF ( nr_p(k,j,i) < 0.0 )  nr_p(k,j,i) = 0.1 * nr(k,j,i)
+                   nr_p(k,j,i) = nr(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +  &
+                                                       tsc(3) * tnr_m(k,j,i) ) &
+                                           - tsc(5) * rdf_sc(k) * nr(k,j,i)
+                   IF ( nr_p(k,j,i) < 0.0 )  nr_p(k,j,i) = 0.0
                 ENDDO
 !
@@ -762,6 +721,5 @@
              ENDIF
 
-          ENDIF
-
+         ENDIF
 !
 !--       If required, compute prognostic equation for turbulent kinetic 
Index: /palm/trunk/SOURCE/read_3d_binary.f90
===================================================================
--- /palm/trunk/SOURCE/read_3d_binary.f90	(revision 1114)
+++ /palm/trunk/SOURCE/read_3d_binary.f90	(revision 1115)
@@ -20,4 +20,5 @@
 ! Current revisions:
 ! -----------------
+! unused variables removed
 !
 ! Former revisions:
@@ -283,5 +284,5 @@
 !--    First compare the version numbers
        READ ( 13 )  version_on_file
-       binary_version = '3.3'
+       binary_version = '3.4'
        IF ( TRIM( version_on_file ) /= TRIM( binary_version ) )  THEN
           WRITE( message_string, * ) 'version mismatch concerning data ',      &
@@ -432,4 +433,32 @@
                                   tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
 
+                CASE ( 'nr' )
+                   IF ( k == 1 )  READ ( 13 )  tmp_3d
+                   nr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
+                                   tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
+
+                CASE ( 'nr_av' )
+                   IF ( .NOT. ALLOCATED( nr_av ) )  THEN
+                      ALLOCATE( nr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+                   ENDIF
+                   IF ( k == 1 )  READ ( 13 )  tmp_3d
+                   nr_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
+                                    tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
+
+                CASE ( 'nrs' )
+                   IF ( k == 1 )  READ ( 13 )  tmp_2d
+                   nrs(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  = &
+                                          tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
+
+                CASE ( 'nrsws' )
+                   IF ( k == 1 )  READ ( 13 )  tmp_2d
+                   nrsws(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  = &
+                                          tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
+
+                CASE ( 'nrswst' )
+                   IF ( k == 1 )  READ ( 13 )  tmp_2d
+                   nrswst(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  = &
+                                          tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
+
                 CASE ( 'p' )
                    IF ( k == 1 )  READ ( 13 )  tmp_3d
@@ -516,4 +545,17 @@
                                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
 
+                CASE ( 'qc' )
+                   IF ( k == 1 )  READ ( 13 )  tmp_3d
+                   qc(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
+                                       tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
+
+                CASE ( 'qc_av' )
+                   IF ( .NOT. ALLOCATED( qc_av ) )  THEN
+                      ALLOCATE( qc_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+                   ENDIF
+                   IF ( k == 1 )  READ ( 13 )  tmp_3d
+                   qc_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
+                                       tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
+
                 CASE ( 'ql' )
                    IF ( k == 1 )  READ ( 13 )  tmp_3d
@@ -534,9 +576,4 @@
                                    tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
 
-                CASE ( 'nr' )
-                   IF ( k == 1 )  READ ( 13 )  tmp_3d
-                   nr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
-                                   tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
-
                 CASE ( 'qr_av' )
                    IF ( .NOT. ALLOCATED( qr_av ) )  THEN
@@ -547,12 +584,4 @@
                                     tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
 
-                CASE ( 'nr_av' )
-                   IF ( .NOT. ALLOCATED( nr_av ) )  THEN
-                      ALLOCATE( nr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
-                   ENDIF
-                   IF ( k == 1 )  READ ( 13 )  tmp_3d
-                   nr_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &
-                                    tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
-
                 CASE ( 'qrs' )
                    IF ( k == 1 )  READ ( 13 )  tmp_2d
@@ -568,19 +597,4 @@
                    IF ( k == 1 )  READ ( 13 )  tmp_2d
                    qrswst(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  = &
-                                          tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
-
-                CASE ( 'nrs' )
-                   IF ( k == 1 )  READ ( 13 )  tmp_2d
-                   nrs(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  = &
-                                          tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
-
-                CASE ( 'nrsws' )
-                   IF ( k == 1 )  READ ( 13 )  tmp_2d
-                   nrsws(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  = &
-                                          tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
-
-                CASE ( 'nrswst' )
-                   IF ( k == 1 )  READ ( 13 )  tmp_2d
-                   nrswst(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  = &
                                           tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
 
Index: /palm/trunk/SOURCE/read_var_list.f90
===================================================================
--- /palm/trunk/SOURCE/read_var_list.f90	(revision 1114)
+++ /palm/trunk/SOURCE/read_var_list.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! unused variables removed
 !
 ! Former revisions:
@@ -185,5 +185,5 @@
 !-- Make version number check first
     READ ( 13 )  version_on_file
-    binary_version = '3.8'
+    binary_version = '3.9'
     IF ( TRIM( version_on_file ) /= TRIM( binary_version ) )  THEN
        WRITE( message_string, * ) 'version mismatch concerning control ', &
@@ -251,5 +251,5 @@
        ALLOCATE( lad(0:nz+1), ug(0:nz+1), u_init(0:nz+1), vg(0:nz+1),    &
                  v_init(0:nz+1), pt_init(0:nz+1), q_init(0:nz+1),        &
-                 qr_init(0:nz+1), nr_init(0:nz+1), sa_init(0:nz+1),      &
+                 sa_init(0:nz+1),                                        &
                  hom(0:nz+1,2,pr_palm+max_pr_user,0:statistic_regions),  &
                  hom_sum(0:nz+1,pr_palm+max_pr_user,0:statistic_regions) )
@@ -284,10 +284,4 @@
           CASE ( 'bc_ns' )
              READ ( 13 )  bc_ns
-          CASE ( 'bc_nr_b' )
-             READ ( 13 )  bc_nr_b
-          CASE ( 'bc_nr_t' )
-             READ ( 13 )  bc_nr_t
-          CASE ( 'bc_nr_t_val' )
-             READ ( 13 )  bc_nr_t_val
           CASE ( 'bc_p_b' )
              READ ( 13 )  bc_p_b
@@ -306,10 +300,4 @@
           CASE ( 'bc_q_t_val' )
              READ ( 13 )  bc_q_t_val
-          CASE ( 'bc_qr_b' )
-             READ ( 13 )  bc_qr_b
-          CASE ( 'bc_qr_t' )
-             READ ( 13 )  bc_qr_t
-          CASE ( 'bc_qr_t_val' )
-             READ ( 13 )  bc_qr_t_val
           CASE ( 'bc_s_b' )
              READ ( 13 )  bc_s_b
@@ -474,6 +462,6 @@
           CASE ( 'momentum_advec' )
              READ ( 13 )  momentum_advec
-          CASE ( 'nc' )
-             READ ( 13 )  nc
+          CASE ( 'nc_const' )
+             READ ( 13 )  nc_const
           CASE ( 'netcdf_precision' )
              READ ( 13 )  netcdf_precision
@@ -482,16 +470,4 @@
           CASE ( 'ngsrb' )
              READ ( 13 )  ngsrb
-          CASE ( 'nr_init' )
-             READ ( 13 )  nr_init
-          CASE ( 'nr_surface' )
-             READ ( 13 )  nr_surface
-          CASE ( 'nr_surface_initial_change' )
-             READ ( 13 )  nr_surface_initial_change
-          CASE ( 'nr_vertical_gradient' )
-             READ ( 13 )  nr_vertical_gradient
-          CASE ( 'nr_vertical_gradient_level' )
-             READ ( 13 )  nr_vertical_gradient_level
-          CASE ( 'nr_vertical_gradient_level_ind' )
-             READ ( 13 )  nr_vertical_gradient_level_ind
           CASE ( 'nsor' )
              READ ( 13 )  nsor
@@ -560,16 +536,4 @@
           CASE ( 'q_vertical_gradient_level_ind' )
              READ ( 13 )  q_vertical_gradient_level_ind
-          CASE ( 'qr_init' )
-             READ ( 13 )  qr_init
-          CASE ( 'qr_surface' )
-             READ ( 13 )  qr_surface
-          CASE ( 'qr_surface_initial_change' )
-             READ ( 13 )  qr_surface_initial_change
-          CASE ( 'qr_vertical_gradient' )
-             READ ( 13 )  qr_vertical_gradient
-          CASE ( 'qr_vertical_gradient_level' )
-             READ ( 13 )  qr_vertical_gradient_level
-          CASE ( 'qr_vertical_gradient_level_ind' )
-             READ ( 13 )  qr_vertical_gradient_level_ind
           CASE ( 'radiation' )
              READ ( 13 )  radiation
@@ -618,8 +582,4 @@
           CASE ( 'surface_waterflux' )
              READ ( 13 )  surface_waterflux      
-          CASE ( 'surface_waterflux_nr' )
-             READ ( 13 )  surface_waterflux_nr    
-          CASE ( 'surface_waterflux_qr' )
-             READ ( 13 )  surface_waterflux_qr    
           CASE ( 's_surface' )
              READ ( 13 )  s_surface
Index: /palm/trunk/SOURCE/sum_up_3d_data.f90
===================================================================
--- /palm/trunk/SOURCE/sum_up_3d_data.f90	(revision 1114)
+++ /palm/trunk/SOURCE/sum_up_3d_data.f90	(revision 1115)
@@ -20,4 +20,5 @@
 ! Current revisions:
 ! -----------------
+! ql is calculated by calc_liquid_water_content
 !
 ! Former revisions:
@@ -164,4 +165,10 @@
                 ENDIF
                 q_av = 0.0
+
+             CASE ( 'qc' )
+                IF ( .NOT. ALLOCATED( qc_av ) )  THEN
+                   ALLOCATE( qc_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
+                ENDIF
+                qc_av = 0.0
 
              CASE ( 'ql' )
@@ -411,4 +418,13 @@
              ENDDO
 
+          CASE ( 'qc' )
+             DO  i = nxlg, nxrg
+                DO  j = nysg, nyng
+                   DO  k = nzb, nzt+1
+                      qc_av(k,j,i) = qc_av(k,j,i) + qc(k,j,i)
+                   ENDDO
+                ENDDO
+             ENDDO
+
           CASE ( 'ql' )
              DO  i = nxlg, nxrg
Index: /palm/trunk/SOURCE/swap_timelevel.f90
===================================================================
--- /palm/trunk/SOURCE/swap_timelevel.f90	(revision 1114)
+++ /palm/trunk/SOURCE/swap_timelevel.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! calculation of qr and nr is restricted to precipitation
 ! 
 ! Former revisions:
@@ -130,5 +130,6 @@
           IF ( humidity  .OR.  passive_scalar )  THEN
              q => q_1;    q_p => q_2
-             IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation )  THEN
                 qr => qr_1;    qr_p => qr_2
                 nr => nr_1;    nr_p => nr_2
@@ -153,5 +154,6 @@
           IF ( humidity  .OR.  passive_scalar )  THEN
              q => q_2;    q_p => q_1
-             IF ( cloud_physics  .AND.  icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.  &
+                  precipitation)  THEN
                 qr => qr_2;    qr_p => qr_1
                 nr => nr_2;    nr_p => nr_1
Index: /palm/trunk/SOURCE/time_integration.f90
===================================================================
--- /palm/trunk/SOURCE/time_integration.f90	(revision 1114)
+++ /palm/trunk/SOURCE/time_integration.f90	(revision 1115)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! calculation of qr and nr is restricted to precipitation
 !
 ! Former revisions:
@@ -266,5 +266,6 @@
           IF (humidity  .OR.  passive_scalar)  THEN
              CALL exchange_horiz( q_p, nbgp )
-             IF ( cloud_physics .AND. icloud_scheme == 0 )  THEN
+             IF ( cloud_physics  .AND.  icloud_scheme == 0  .AND.              &
+                  precipitation)  THEN
                 CALL exchange_horiz( qr_p, nbgp )
                 CALL exchange_horiz( nr_p, nbgp )
@@ -349,11 +350,4 @@
 
 !
-!--       If required, compute virtuell potential temperature
-          IF ( humidity )  THEN
-             CALL compute_vpt
-             !$acc update device( vpt )
-          ENDIF
-
-!
 !--       If required, compute liquid water content
           IF ( cloud_physics )  THEN
@@ -361,5 +355,10 @@
              !$acc update device( ql )
           ENDIF
-
+! 
+!--       If required, compute virtual potential temperature 
+          IF ( humidity )  THEN 
+             CALL compute_vpt 
+             !$acc update device( vpt ) 
+          ENDIF 
 !
 !--       Compute the diffusion quantities
@@ -373,5 +372,4 @@
                 CALL cpu_log( log_point(19), 'prandtl_fluxes', 'stop' )
              ENDIF
-
 !
 !--          Compute the diffusion coefficients
Index: /palm/trunk/SOURCE/write_3d_binary.f90
===================================================================
--- /palm/trunk/SOURCE/write_3d_binary.f90	(revision 1114)
+++ /palm/trunk/SOURCE/write_3d_binary.f90	(revision 1115)
@@ -20,11 +20,9 @@
 ! Current revisions:
 ! -----------------
+! qr and nr are restricted to precipitation
 !
 ! Former revisions:
 ! -----------------
 ! $Id$
-!
-! 1074 2012-12-03 11:00:05Z hoffmann
-! Bugfix in writing qr, nr
 !
 ! 1053 2012-11-13 17:11:03Z hoffmann
@@ -117,5 +115,5 @@
 !
 !-- Write arrays.
-    binary_version = '3.3'
+    binary_version = '3.4'
 
     WRITE ( 14 )  binary_version
@@ -179,19 +177,25 @@
              WRITE ( 14 )  'ql_av               ';  WRITE ( 14 )  ql_av
           ENDIF
-          IF ( icloud_scheme == 0  .AND.  cloud_physics )  THEN
-             WRITE ( 14 )  'qr                  ';  WRITE ( 14 ) qr
-             WRITE ( 14 )  'nr                  ';  WRITE ( 14 ) nr
-             IF ( ALLOCATED( qr_av ) )  THEN
-                WRITE ( 14 )  'qr_av               ';  WRITE ( 14 )  qr_av
+          IF ( icloud_scheme == 0 )  THEN
+             WRITE ( 14 )  'qc                  ';  WRITE ( 14 ) qc
+             IF ( ALLOCATED( qc_av ) )  THEN
+                WRITE ( 14 )  'qc_av               ';  WRITE ( 14 )  qc_av
              ENDIF
-             IF ( ALLOCATED( nr_av ) )  THEN
-                WRITE ( 14 )  'nr_av               ';  WRITE ( 14 )  nr_av
+             IF ( precipitation )  THEN
+                WRITE ( 14 )  'nr                  ';  WRITE ( 14 ) nr
+                IF ( ALLOCATED( nr_av ) )  THEN
+                   WRITE ( 14 )  'nr_av               ';  WRITE ( 14 )  nr_av
+                ENDIF
+                WRITE ( 14 )  'nrs                 ';  WRITE ( 14 ) nrs
+                WRITE ( 14 )  'nrsws               ';  WRITE ( 14 ) nrsws
+                WRITE ( 14 )  'nrswst              ';  WRITE ( 14 ) nrswst
+                WRITE ( 14 )  'qr                  ';  WRITE ( 14 ) qr
+                IF ( ALLOCATED( qr_av ) )  THEN
+                   WRITE ( 14 )  'qr_av               ';  WRITE ( 14 )  qr_av
+                ENDIF
+                WRITE ( 14 )  'qrs                 ';  WRITE ( 14 ) qrs
+                WRITE ( 14 )  'qrsws               ';  WRITE ( 14 ) qrsws
+                WRITE ( 14 )  'qrswst              ';  WRITE ( 14 ) qrswst
              ENDIF
-             WRITE ( 14 )  'qrs                 ';  WRITE ( 14 ) qrs
-             WRITE ( 14 )  'qrsws               ';  WRITE ( 14 ) qrsws
-             WRITE ( 14 )  'qrswst              ';  WRITE ( 14 ) qrswst
-             WRITE ( 14 )  'nrs                 ';  WRITE ( 14 ) nrs
-             WRITE ( 14 )  'nrsws               ';  WRITE ( 14 ) nrsws
-             WRITE ( 14 )  'nrswst              ';  WRITE ( 14 ) nrswst
           ENDIF
        ENDIF
Index: /palm/trunk/SOURCE/write_var_list.f90
===================================================================
--- /palm/trunk/SOURCE/write_var_list.f90	(revision 1114)
+++ /palm/trunk/SOURCE/write_var_list.f90	(revision 1115)
@@ -20,4 +20,5 @@
 ! Current revisions:
 ! -----------------
+! unused variables removed
 !
 ! Former revisions:
@@ -165,5 +166,5 @@
 
 
-    binary_version = '3.8'
+    binary_version = '3.9'
 
     WRITE ( 14 )  binary_version
@@ -204,10 +205,4 @@
     WRITE ( 14 )  'bc_ns                         '
     WRITE ( 14 )  bc_ns
-    WRITE ( 14 )  'bc_nr_b                       '
-    WRITE ( 14 )  bc_nr_b
-    WRITE ( 14 )  'bc_nr_t                       '
-    WRITE ( 14 )  bc_nr_t
-    WRITE ( 14 )  'bc_nr_t_val                   '
-    WRITE ( 14 )  bc_nr_t_val
     WRITE ( 14 )  'bc_p_b                        '
     WRITE ( 14 )  bc_p_b
@@ -226,10 +221,4 @@
     WRITE ( 14 )  'bc_q_t_val                    '
     WRITE ( 14 )  bc_q_t_val
-    WRITE ( 14 )  'bc_qr_b                       '
-    WRITE ( 14 )  bc_qr_b
-    WRITE ( 14 )  'bc_qr_t                       '
-    WRITE ( 14 )  bc_qr_t
-    WRITE ( 14 )  'bc_qr_t_val                   '
-    WRITE ( 14 )  bc_qr_t_val
     WRITE ( 14 )  'bc_s_b                        '
     WRITE ( 14 )  bc_s_b
@@ -392,6 +381,6 @@
     WRITE ( 14 )  'momentum_advec                '
     WRITE ( 14 )  momentum_advec
-    WRITE ( 14 )  'nc                            '
-    WRITE ( 14 )  nc
+    WRITE ( 14 )  'nc_const                      '
+    WRITE ( 14 )  nc_const
     WRITE ( 14 )  'netcdf_precision              '
     WRITE ( 14 )  netcdf_precision
@@ -400,16 +389,4 @@
     WRITE ( 14 )  'ngsrb                         '
     WRITE ( 14 )  ngsrb
-    WRITE ( 14 )  'nr_init                       '
-    WRITE ( 14 )  nr_init
-    WRITE ( 14 )  'nr_surface                    '
-    WRITE ( 14 )  nr_surface
-    WRITE ( 14 )  'nr_surface_initial_change     '
-    WRITE ( 14 )  nr_surface_initial_change
-    WRITE ( 14 )  'nr_vertical_gradient          '
-    WRITE ( 14 )  nr_vertical_gradient
-    WRITE ( 14 )  'nr_vertical_gradient_level    '
-    WRITE ( 14 )  nr_vertical_gradient_level
-    WRITE ( 14 )  'nr_vertical_gradient_level_ind'
-    WRITE ( 14 )  nr_vertical_gradient_level_ind
     WRITE ( 14 )  'nsor                          '
     WRITE ( 14 )  nsor
@@ -476,16 +453,4 @@
     WRITE ( 14 )  'q_vertical_gradient_level_ind '
     WRITE ( 14 )  q_vertical_gradient_level_ind
-    WRITE ( 14 )  'qr_init                       '
-    WRITE ( 14 )  qr_init
-    WRITE ( 14 )  'qr_surface                    '
-    WRITE ( 14 )  qr_surface
-    WRITE ( 14 )  'qr_surface_initial_change     '
-    WRITE ( 14 )  qr_surface_initial_change
-    WRITE ( 14 )  'qr_vertical_gradient          '
-    WRITE ( 14 )  qr_vertical_gradient
-    WRITE ( 14 )  'qr_vertical_gradient_level    '
-    WRITE ( 14 )  qr_vertical_gradient_level
-    WRITE ( 14 )  'qr_vertical_gradient_level_ind'
-    WRITE ( 14 )  qr_vertical_gradient_level_ind
     WRITE ( 14 )  'radiation                     '
     WRITE ( 14 )  radiation
@@ -534,8 +499,4 @@
     WRITE ( 14 )  'surface_waterflux             '
     WRITE ( 14 )  surface_waterflux    
-    WRITE ( 14 )  'surface_waterflux_nr          '
-    WRITE ( 14 )  surface_waterflux_nr    
-    WRITE ( 14 )  'surface_waterflux_qr          '
-    WRITE ( 14 )  surface_waterflux_qr    
     WRITE ( 14 )  's_surface                     '
     WRITE ( 14 )  s_surface
