Changeset 266 for palm/trunk/SOURCE/set_particle_attributes.f90
 Timestamp:
 Mar 25, 2009 7:43:59 AM (14 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

palm/trunk/SOURCE/set_particle_attributes.f90
r264 r266 32 32 33 33 INTEGER :: i, j, k, n 34 REAL :: aa, bb, cc, dd, gg, pt_int, pt_int_l, pt_int_u, w_int, & 35 w_int_l, w_int_u, x, y 34 REAL :: aa, absuv, bb, cc, dd, gg, pt_int, pt_int_l, pt_int_u, u_int, & 35 u_int_l, u_int_u, v_int, v_int_l, v_int_u, w_int, w_int_l, & 36 w_int_u, x, y 36 37 37 38 … … 41 42 ! Set particle color 42 43 IF ( particle_color == 'absuv' ) THEN 44 45 ! 46 ! Set particle color depending on the absolute value of the horizontal 47 ! velocity 48 DO n = 1, number_of_particles 49 ! 50 ! Interpolate u velocitycomponent, determine left, front, bottom 51 ! index of uarray 52 i = ( particles(n)%x + 0.5 * dx ) * ddx 53 j = particles(n)%y * ddy 54 k = ( particles(n)%z + 0.5 * dz * atmos_ocean_sign ) / dz & 55 + offset_ocean_nzt ! only exact if equidistant 56 57 ! 58 ! Interpolation of the velocity components in the xyplane 59 x = particles(n)%x + ( 0.5  i ) * dx 60 y = particles(n)%y  j * dy 61 aa = x**2 + y**2 62 bb = ( dx  x )**2 + y**2 63 cc = x**2 + ( dy  y )**2 64 dd = ( dx  x )**2 + ( dy  y )**2 65 gg = aa + bb + cc + dd 66 67 u_int_l = ( ( gg  aa ) * u(k,j,i) + ( gg  bb ) * u(k,j,i+1) & 68 + ( gg  cc ) * u(k,j+1,i) + ( gg  dd ) * u(k,j+1,i+1) & 69 ) / ( 3.0 * gg )  u_gtrans 70 IF ( k+1 == nzt+1 ) THEN 71 u_int = u_int_l 72 ELSE 73 u_int_u = ( ( ggaa ) * u(k+1,j,i) + ( ggbb ) * u(k+1,j,i+1) & 74 + ( ggcc ) * u(k+1,j+1,i) + ( ggdd ) * u(k+1,j+1,i+1) & 75 ) / ( 3.0 * gg )  u_gtrans 76 u_int = u_int_l + ( particles(n)%z  zu(k) ) / dz * & 77 ( u_int_u  u_int_l ) 78 ENDIF 79 80 ! 81 ! Same procedure for interpolation of the v velocitycomponent (adopt 82 ! index k from u velocitycomponent) 83 i = particles(n)%x * ddx 84 j = ( particles(n)%y + 0.5 * dy ) * ddy 85 86 x = particles(n)%x  i * dx 87 y = particles(n)%y + ( 0.5  j ) * dy 88 aa = x**2 + y**2 89 bb = ( dx  x )**2 + y**2 90 cc = x**2 + ( dy  y )**2 91 dd = ( dx  x )**2 + ( dy  y )**2 92 gg = aa + bb + cc + dd 93 94 v_int_l = ( ( gg  aa ) * v(k,j,i) + ( gg  bb ) * v(k,j,i+1) & 95 + ( gg  cc ) * v(k,j+1,i) + ( gg  dd ) * v(k,j+1,i+1) & 96 ) / ( 3.0 * gg )  v_gtrans 97 IF ( k+1 == nzt+1 ) THEN 98 v_int = v_int_l 99 ELSE 100 v_int_u = ( ( ggaa ) * v(k+1,j,i) + ( ggbb ) * v(k+1,j,i+1) & 101 + ( ggcc ) * v(k+1,j+1,i) + ( ggdd ) * v(k+1,j+1,i+1) & 102 ) / ( 3.0 * gg )  v_gtrans 103 v_int = v_int_l + ( particles(n)%z  zu(k) ) / dz * & 104 ( v_int_u  v_int_l ) 105 ENDIF 106 107 absuv = SQRT( u_int**2 + v_int**2 ) 108 109 ! 110 ! Limit values by the given interval and normalize to interval [0,1] 111 absuv = MIN( absuv, color_interval(2) ) 112 absuv = MAX( absuv, color_interval(1) ) 113 114 absuv = ( absuv  color_interval(1) ) / & 115 ( color_interval(2)  color_interval(1) ) 116 117 ! 118 ! Number of available colors is defined in init_dvrp 119 particles(n)%color = 1 + absuv * ( dvrp_colortable_entries_prt  1 ) 120 121 ENDDO 43 122 44 123 ELSEIF ( particle_color == 'pt*' ) THEN … … 110 189 particles(n)%color = 1 + pt_int * ( dvrp_colortable_entries_prt  1 ) 111 190 112 ENDDO191 ENDDO 113 192 114 193 ELSEIF ( particle_color == 'z' ) THEN 194 ! 195 ! Set particle color depending on the height above the bottom 196 ! boundary (z=0) 197 DO n = 1, number_of_particles 198 199 height = particles(n)%z 200 ! 201 ! Limit values by the given interval and normalize to interval [0,1] 202 height = MIN( height, color_interval(2) ) 203 height = MAX( height, color_interval(1) ) 204 205 height = ( height  color_interval(1) ) / & 206 ( color_interval(2)  color_interval(1) ) 207 208 ! 209 ! Number of available colors is defined in init_dvrp 210 particles(n)%color = 1 + height * ( dvrp_colortable_entries_prt  1 ) 211 212 ENDDO 115 213 116 214 ENDIF … … 149 247 ENDIF 150 248 151 particles(n)%dvrp_psize = ( 0.2 + MIN( 0.6, ABS( w_int ) / 6.0 ) ) * & 152 dx * 1.4 249 ! 250 ! Limit values by the given interval and normalize to interval [0,1] 251 w_int = ABS( w_int ) 252 w_int = MIN( w_int, dvrpsize_interval(2) ) 253 w_int = MAX( w_int, dvrpsize_interval(1) ) 254 255 w_int = ( w_int  dvrpsize_interval(1) ) / & 256 ( dvrpsize_interval(2)  dvrpsize_interval(1) ) 257 258 particles(n)%dvrp_psize = ( 0.25 + w_int * 0.6 ) * dx 153 259 154 260 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.