Ignore:
Timestamp:
Apr 13, 2016 1:12:11 PM (6 years ago)
Author:
hoffmann
Message:

Interpolation of collision kernels adjusted to more reasonable values

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/lpm_collision_kernels_mod.f90

    r1851 r1858  
    1919! Current revisions:
    2020! -----------------
    21 !
     21! Interpolation of collision kernels adjusted to more reasonable values.
     22! Reformatting of the code.
    2223!
    2324! Former revisions:
     
    2728! 1850 2016-04-08 13:29:27Z maronga
    2829! Module renamed
    29 !
    3030!
    3131! 1822 2016-04-07 07:49:42Z hoffmann
     
    195195!
    196196!--       Calculate the radius class bounds with logarithmic distances
    197 !--       in the interval [1.0E-6, 2.0E-4] m
     197!--       in the interval [1.0E-6, 1000.0E-6] m
    198198          rclass_lbound = LOG( 1.0E-6_wp )
    199           rclass_ubound = LOG( 2.0E-4_wp )
     199          rclass_ubound = LOG( 1000.0E-6_wp )
    200200          radclass(1)   = EXP( rclass_lbound )
    201201          DO  i = 2, radius_classes
     
    206206
    207207!
    208 !--       Set the class bounds for dissipation in interval [0.0, 0.1] m**2/s**3
     208!--       Set the class bounds for dissipation in interval [0.0, 600.0] cm**2/s**3
    209209          DO  i = 1, dissipation_classes
    210              epsclass(i) = 0.1_wp * REAL( i, KIND=wp ) / dissipation_classes
     210             epsclass(i) = 0.06_wp * REAL( i, KIND=wp ) / dissipation_classes
    211211          ENDDO
    212212!
     
    286286          DEALLOCATE( ec, ecf, epsclass, gck, hkernel, winf )
    287287
    288        ELSEIF( collision_kernel == 'hall'  .OR.  collision_kernel == 'wang' ) &
    289        THEN
    290 !
    291 !--       Initial settings for Hall- and Wang-Kernel
    292 !--       To be done: move here parts from turbsd, fallg, ecoll, etc.
    293288       ENDIF
    294289
     
    307302
    308303       USE particle_attributes,                                                &
    309            ONLY:  prt_count, radius_classes, wang_kernel
     304           ONLY:  number_of_particles, prt_count, radius_classes, wang_kernel
    310305
    311306       IMPLICIT NONE
     
    316311       INTEGER(iwp) ::  j1     !<
    317312       INTEGER(iwp) ::  k1     !<
    318        INTEGER(iwp) ::  pend   !<
    319        INTEGER(iwp) ::  pstart !<
    320 
    321 
    322        pstart = 1
    323        pend   = prt_count(k1,j1,i1)
    324        radius_classes = prt_count(k1,j1,i1)
    325 
    326        ALLOCATE( ec(1:radius_classes,1:radius_classes), &
    327                  radclass(1:radius_classes), winf(1:radius_classes) )
     313
     314
     315       number_of_particles = prt_count(k1,j1,i1)
     316       radius_classes      = number_of_particles   ! necessary to use the same
     317                                                   ! subroutines as for
     318                                                   ! precalculated kernels
     319
     320       ALLOCATE( ec(1:number_of_particles,1:number_of_particles), &
     321                 radclass(1:number_of_particles), winf(1:number_of_particles) )
    328322
    329323!
    330324!--    Store particle radii on the radclass array
    331        radclass(1:radius_classes) = particles(pstart:pend)%radius
     325       radclass(1:number_of_particles) = particles(1:number_of_particles)%radius
    332326
    333327       IF ( wang_kernel )  THEN
     
    341335!
    342336!--       Call routines to calculate efficiencies for the Wang kernel
    343           ALLOCATE( gck(1:radius_classes,1:radius_classes), &
    344                     ecf(1:radius_classes,1:radius_classes) )
     337          ALLOCATE( gck(1:number_of_particles,1:number_of_particles), &
     338                    ecf(1:number_of_particles,1:number_of_particles) )
    345339
    346340          CALL turbsd
     
    348342          CALL effic
    349343
    350           DO  j = 1, radius_classes
    351              DO  i =  1, radius_classes
    352                 ckernel(pstart+i-1,pstart+j-1,1) = ec(i,j) * gck(i,j) * ecf(i,j)
     344          DO  j = 1, number_of_particles
     345             DO  i =  1, number_of_particles
     346                ckernel(1+i-1,1+j-1,1) = ec(i,j) * gck(i,j) * ecf(i,j)
    353347             ENDDO
    354348          ENDDO
     
    362356          CALL effic
    363357
    364           DO  j = 1, radius_classes
    365              DO  i =  1, radius_classes
    366                 ckernel(pstart+i-1,pstart+j-1,1) = pi *                       &
    367                                           ( radclass(j) + radclass(i) )**2    &
    368                                           * ec(i,j) * ABS( winf(j) - winf(i) )
     358          DO  j = 1, number_of_particles
     359             DO  i =  1, number_of_particles
     360                ckernel(i,j,1) = pi * ( radclass(j) + radclass(i) )**2         &
     361                                    * ec(i,j) * ABS( winf(j) - winf(i) )
    369362             ENDDO
    370363          ENDDO
Note: See TracChangeset for help on using the changeset viewer.