Changeset 790 for palm/trunk

Timestamp:

Nov 29, 2011 3:11:20 AM (13 years ago)

Author:

raasch

Message:

Bugfix for output of mean particle radius + preliminary works for implementing the Wang collision kernel

Location:

palm/trunk/SOURCE

Files:

• Makefile (modified) (5 diffs)
• data_output_2d.f90 (modified) (2 diffs)
• data_output_3d.f90 (modified) (3 diffs)
• diffusion_e.f90 (modified) (5 diffs)
• init_3d_model.f90 (modified) (2 diffs)
• modules.f90 (modified) (2 diffs)
• package_parin.f90 (modified) (2 diffs)
• sum_up_3d_data.f90 (modified) (2 diffs)
• time_integration.f90 (modified) (3 diffs)
• wang_kernel.f90 (added)
• write_3d_binary.f90 (modified) (2 diffs)

palm/trunk/SOURCE/Makefile

@@ -4 +4 @@
 # Current revisions:
 # ------------------
+# +wang_kernel
 #
 # Former revisions:
@@ -104 +105 @@
         user_particle_attributes.f90 user_read_restart_data.f90 \
         user_spectra.f90 user_statistics.f90 wall_fluxes.f90 \
-        write_3d_binary.f90 write_compressed.f90 write_var_list.f90
+        wang_kernel.f90 write_3d_binary.f90 write_compressed.f90 \
+        write_var_list.f90
 
 OBJS =  advec_particles.o advec_s_bc.o advec_s_pw.o advec_s_up.o \
@@ -146 +148 @@
         user_module.o user_parin.o user_particle_attributes.o \
         user_read_restart_data.o user_spectra.o user_statistics.o \
-        wall_fluxes.o write_3d_binary.o write_compressed.o write_var_list.o
+        wall_fluxes.o wang_kernel.o write_3d_binary.o write_compressed.o \
+        write_var_list.o
 
 CC = cc
@@ -175 +178 @@
 
 
-advec_particles.o: modules.o random_function.o
+advec_particles.o: modules.o random_function.o wang_kernel.o
 advec_s_bc.o: modules.o
 advec_s_pw.o: modules.o
@@ -322 +325 @@
 user_statistics.o: modules.o user_module.o
 wall_fluxes.o: modules.o
+wang_kernel.o: modules.o user_module.o
 write_3d_binary.o: modules.o random_function.o
 write_compressed.o: modules.o

palm/trunk/SOURCE/data_output_2d.f90

@@ -4 +4 @@
 ! Current revisions:
 ! -----------------
+! bugfix: calculation of 'pr' must depend on the particle weighting factor
 !
 ! Former revisions:
@@ -303 +304 @@
                                s_r4 = 0.0
                                DO  n = psi, psi+prt_count(k,j,i)-1
-                                  s_r3 = s_r3 + particles(n)%radius**3
-                                  s_r4 = s_r4 + particles(n)%radius**4
+                                  s_r3 = s_r3 + particles(n)%radius**3 * &
+                                                particles(n)%weight_factor
+                                  s_r4 = s_r4 + particles(n)%radius**4 * &
+                                                particles(n)%weight_factor
                                ENDDO
                                IF ( s_r3 /= 0.0 )  THEN

palm/trunk/SOURCE/data_output_3d.f90

@@ -4 +4 @@
 ! Current revisions:
 ! -----------------
+! bugfix: calculation of 'pr' must depend on the particle weighting factor,
+! nzt+1 replaced by nz_do3d for 'pr'
 !
 ! Former revisions:
@@ -207 +209 @@
                 DO  i = nxl, nxr
                    DO  j = nys, nyn
-                      DO  k = nzb, nzt+1
+                      DO  k = nzb, nz_do3d
                          psi = prt_start_index(k,j,i)
                          s_r3 = 0.0
                          s_r4 = 0.0
                          DO  n = psi, psi+prt_count(k,j,i)-1
-                            s_r3 = s_r3 + particles(n)%radius**3
-                            s_r4 = s_r4 + particles(n)%radius**4
+                         s_r3 = s_r3 + particles(n)%radius**3 * &
+                                       particles(n)%weight_factor
+                         s_r4 = s_r4 + particles(n)%radius**4 * &
+                                       particles(n)%weight_factor
                          ENDDO
                          IF ( s_r3 /= 0.0 )  THEN
@@ -227 +231 @@
                 DO  i = nxlg, nxrg
                    DO  j = nysg, nyng
-                      DO  k = nzb, nzt+1
+                      DO  k = nzb, nz_do3d
                          local_pf(i,j,k) = tend(k,j,i)
                       ENDDO

palm/trunk/SOURCE/diffusion_e.f90

@@ -4 +4 @@
 ! Current revisions:
 ! -----------------
+! diss is also calculated in case that the Wang kernel is used
 !
 ! Former revisions:
@@ -160 +161 @@
 !--          Store dissipation if needed for calculating the sgs particle
 !--          velocities
-             IF ( use_sgs_for_particles )  THEN
+             IF ( use_sgs_for_particles  .OR.  wang_collision_kernel )  THEN
                 DO  j = nys, nyn
                    DO  k = nzb_s_inner(j,i)+1, nzt
@@ -250 +251 @@
 !--          Store dissipation if needed for calculating the sgs particle
 !--          velocities
-             IF ( use_sgs_for_particles )  THEN
+             IF ( use_sgs_for_particles  .OR.  wang_collision_kernel )  THEN
                 DO  j = nys, nyn
                    DO  k = nzb_s_inner(j,i)+1, nzt
@@ -264 +265 @@
 !
 !--    Boundary condition for dissipation
-       IF ( use_sgs_for_particles )  THEN
+       IF ( use_sgs_for_particles  .OR.  wang_collision_kernel )  THEN
           DO  i = nxl, nxr
              DO  j = nys, nyn
@@ -368 +369 @@
 !
 !--    Store dissipation if needed for calculating the sgs particle velocities
-       IF ( use_sgs_for_particles )  THEN
+       IF ( use_sgs_for_particles  .OR.  wang_collision_kernel )  THEN
           DO  k = nzb_s_inner(j,i)+1, nzt
              diss(k,j,i) = dissipation(k)

palm/trunk/SOURCE/init_3d_model.f90

@@ -7 +7 @@
 ! Current revisions:
 ! ------------------
-! 
+! diss is also allocated in case that the Wang kernel is used
 !
 ! Former revisions:
@@ -331 +331 @@
 !-- 3D-array for storing the dissipation, needed for calculating the sgs
 !-- particle velocities
-    IF ( use_sgs_for_particles )  THEN
+    IF ( use_sgs_for_particles  .OR.  wang_collision_kernel )  THEN
        ALLOCATE ( diss(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
     ELSE

palm/trunk/SOURCE/modules.f90

@@ -5 +5 @@
 ! Current revisions:
 ! -----------------
-! 
+! +turbulence_effects_on_collision, wang_collision_kernel
 !
 ! Former revisions:
@@ -1170 +1170 @@
     INTEGER, DIMENSION(:,:,:), ALLOCATABLE ::  prt_count, prt_start_index
 
-    LOGICAL ::  particle_advection = .FALSE., random_start_position = .FALSE., &
+    LOGICAL ::  turbulence_effects_on_collision = .FALSE.,                     &
+                particle_advection = .FALSE., random_start_position = .FALSE., &
                 read_particles_from_restartfile = .TRUE.,                      &
                 uniform_particles = .TRUE., use_particle_tails = .FALSE.,      &
                 use_sgs_for_particles = .FALSE.,                               &
+                wang_collision_kernel = .FALSE.,                               &
                 write_particle_statistics = .FALSE.
 

palm/trunk/SOURCE/package_parin.f90

@@ -4 +4 @@
 ! Current revisions:
 ! -----------------
-! 
+! +turbulence_effects_on_collision, wang_collision_kernel in particles_par
 !
 ! Former revisions:
@@ -83 +83 @@
                                   random_start_position,                       &
                                   read_particles_from_restartfile,             &
-                                  skip_particles_for_tail, use_particle_tails, &
-                                  use_sgs_for_particles,                       &
+                                  skip_particles_for_tail,                     &
+                                  turbulence_effects_on_collision,             &
+                                  use_particle_tails, use_sgs_for_particles,   &
                                   vertical_particle_advection,                 &
+                                  wang_collision_kernel,                       &
                                   write_particle_statistics
     NAMELIST /spectra_par/        averaging_interval_sp, comp_spectra_level,   &

palm/trunk/SOURCE/sum_up_3d_data.f90

@@ -4 +4 @@
 ! Current revisions:
 ! -----------------
+! bugfix: calculation of 'pr' must depend on the particle weighting factor
 !
 ! Former revisions:
@@ -295 +296 @@
                       s_r4 = 0.0
                       DO  n = psi, psi+prt_count(k,j,i)-1
-                         s_r3 = s_r3 + particles(n)%radius**3
-                         s_r4 = s_r4 + particles(n)%radius**4
+                         s_r3 = s_r3 + particles(n)%radius**3 * &
+                                       particles(n)%weight_factor
+                         s_r4 = s_r4 + particles(n)%radius**4 * &
+                                       particles(n)%weight_factor
                       ENDDO
                       IF ( s_r3 /= 0.0 )  THEN

palm/trunk/SOURCE/time_integration.f90

@@ -4 +4 @@
 ! Current revisions:
 ! -----------------
-! 
+! exchange of ghostpoints for array diss
 !
 ! Former revisions:
@@ -208 +208 @@
              CALL exchange_horiz( prho, nbgp )
           ENDIF
-          IF (humidity  .OR.  passive_scalar) CALL exchange_horiz( q_p, nbgp )
+          IF (humidity  .OR.  passive_scalar)  CALL exchange_horiz( q_p, nbgp )
           IF ( cloud_droplets )  THEN
              CALL exchange_horiz( ql, nbgp )
@@ -215 +215 @@
              CALL exchange_horiz( ql_vp, nbgp )
           ENDIF
+          IF ( wang_collision_kernel )  CALL exchange_horiz( diss, nbgp )
 
           CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'stop' )

palm/trunk/SOURCE/write_3d_binary.f90

@@ -4 +4 @@
 ! Current revisions:
 ! -----------------
+! bugfix: output of 'ql' also if cloud droplets are used
 !
 ! Former revisions:
@@ -141 +142 @@
        ENDIF
        WRITE ( 14 )  'q_m                 ';  WRITE ( 14 ) q_m
-       IF ( cloud_physics ) THEN
+       IF ( cloud_physics  .OR.  cloud_droplets )  THEN
           WRITE ( 14 )  'ql                  ';  WRITE ( 14 ) ql
           IF ( ALLOCATED( ql_av ) )  THEN