File: /nfs/home/0/users/jenkins/mfix.git/model/des/comp_mean_fields0.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
3           SUBROUTINE COMP_MEAN_FIELDS0
4     
5     !-----------------------------------------------
6     ! Modules
7     !-----------------------------------------------
8           USE param
9           USE param1
10           USE parallel
11           USE constant
12           USE physprop
13           USE fldvar
14           USE run
15           USE geometry
16           USE indices
17           USE bc
18           USE compar
19           USE sendrecv
20           USE discretelement
21           USE drag
22           USE interpolation
23           use desmpi
24           USE cutcell
25           USE mfix_pic
26           USE mpi_utility
27     
28     
29           use mpi_node_des, only: des_addnodevalues_mean_fields
30           use particle_filter, only: DES_REPORT_MASS_INTERP
31     
32     
33           use functions, only: FLUID_AT
34           use functions, only: FUNIJK
35           use functions, only: IS_ON_myPE_wobnd
36     
37           IMPLICIT NONE
38     !-----------------------------------------------
39     ! Local variables
40     !-----------------------------------------------
41     ! general i, j, k indices
42           INTEGER :: I, J, K, IJK, &
43                      II, JJ, KK
44           INTEGER :: I1, I2, J1, J2, K1, K2
45           INTEGER :: IDIM, IJK2
46           INTEGER :: ICUR, JCUR, KCUR, CUR_IJK
47     ! indices used for interpolation stencil (unclear why IE, JN, KTP are
48     ! needed)
49           INTEGER :: IW, IE, JS, JN, KB, KTP
50     ! i,j,k indices of the fluid cell the particle resides in minus 1
51     ! (e.g., shifted 1 in west, south, bottom direction)
52           INTEGER, DIMENSION(3):: PCELL
53     ! order of interpolation set in the call to set_interpolation_scheme
54     ! unless it is re/set later through the call to
55     ! set_interpolation_stencil
56           INTEGER :: ONEW
57     ! constant whose value depends on dimension of system
58     ! avg_factor=0.250 (in 3D) or =0.50 (in 2D)
59           DOUBLE PRECISION :: AVG_FACTOR
60     ! index of solid phase that particle NP belongs to
61           INTEGER :: M
62     ! particle number index, used for looping
63           INTEGER :: NP, NINDX
64     ! Statistical weight of the particle. Equal to one for DEM
65           DOUBLE PRECISION :: WTP
66     ! one over the solids volume fraction
67           DOUBLE PRECISION :: OEPS
68     
69           DOUBLE PRECISION :: MASS_SOL1, MASS_SOL2
70     ! sum of mass_sol1 and mass_sol2 across all processors
71           DOUBLE PRECISION :: MASS_SOL1_ALL, MASS_SOL2_ALL
72     
73           DOUBLE PRECISION :: TEMP1, TEMP2
74     
75           DOUBLE PRECISION, DIMENSION(3) :: DES_VEL_DENSITY
76           DOUBLE PRECISION :: DES_ROP_DENSITY
77     
78           INTEGER :: COUNT_NODES_OUTSIDE, COUNT_NODES_INSIDE, &
79                      COUNT_NODES_INSIDE_MAX
80           double precision :: RESID_ROPS(DES_MMAX), &
81                               RESID_VEL(3, DES_MMAX)
82           double precision :: NORM_FACTOR, VOL_RATIO
83     !Handan Liu added on Jan 17 2013
84               DOUBLE PRECISION, DIMENSION(2,2,2,3) :: gst_tmp
85               DOUBLE PRECISION, DIMENSION(2,2,2) :: weight_ft
86     !-----------------------------------------------
87     
88     ! initializing
89           MASS_SOL1 = ZERO
90           MASS_SOL2 = ZERO
91           MASS_SOL1_ALL = ZERO
92           MASS_SOL2_ALL = ZERO
93     ! avg_factor=0.25 (in 3D) or =0.5 (in 2D)
94           AVG_FACTOR = merge(0.5d0, 0.25D0, NO_K)
95     
96     ! cartesian_grid related quantities
97           COUNT_NODES_INSIDE_MAX = merge(4, 8, NO_K)
98     
99     
100     ! Initialize entire arrays to zero
101           DES_VEL_NODE = ZERO
102           DES_ROPS_NODE = ZERO
103           DES_ROP_S = zero
104           DES_U_S = ZERO
105           DES_V_S = ZERO
106           IF(DO_K) DES_W_S = ZERO
107     
108     
109     ! sets several quantities including interp_scheme, scheme, and
110     ! order and allocates arrays necessary for interpolation
111           CALL SET_INTERPOLATION_SCHEME(2)
112     
113     !$omp parallel default(shared)                                             &
114     !$omp private(IJK, I, J, K, PCELL, IW, IE, JS, JN, KB, KTP, ONEW, GST_TMP, &
115     !$omp    COUNT_NODES_INSIDE, II, JJ, KK, CUR_IJK, NINDX, NP, WTP, M, ICUR, &
116     !$omp    JCUR, KCUR, WEIGHT_FT, VOL_RATIO, I1, I2, J1, J2, K1, K2, IDIM,   &
117     !$omp    IJK2, NORM_FACTOR, RESID_ROPS, RESID_VEL,COUNT_NODES_OUTSIDE, TEMP1)
118     !$omp do reduction(+:MASS_SOL1) reduction(+:DES_ROPS_NODE,DES_VEL_NODE)
119           DO IJK = IJKSTART3,IJKEND3
120     
121     ! Cycle this cell if not in the fluid domain or if it contains no
122     ! particle/parcel
123              IF(.NOT.FLUID_AT(IJK)) CYCLE
124              IF( PINC(IJK) == 0) CYCLE
125     
126              PCELL(1) = I_OF(IJK)-1
127              PCELL(2) = J_OF(IJK)-1
128              PCELL(3) = merge(K_OF(IJK)-1, 1, DO_K)
129     
130     ! setup the stencil based on the order of interpolation and factoring in
131     ! whether the system has any periodic boundaries. sets onew to order.
132              CALL SET_INTERPOLATION_STENCIL(PCELL, IW, IE, JS, JN, KB, KTP,&
133                 INTERP_SCHEME, DIMN, ORDERNEW=ONEW)
134     
135              COUNT_NODES_OUTSIDE = 0
136     ! Computing/setting the geometric stencil
137              DO K=1, merge(1, ONEW, NO_K)
138              DO J=1, ONEW
139              DO I=1, ONEW
140     
141                 II = IW + I-1
142                 JJ = JS + J-1
143                 KK = KB + K-1
144                 CUR_IJK = funijk_map_c(II,JJ,KK)
145     
146                 GST_TMP(I,J,K,1) = XE(II)
147                 GST_TMP(I,J,K,2) = YN(JJ)
148                 GST_TMP(I,J,K,3) = merge(DZ(1), ZT(KK), NO_K)
149     
150                 IF(CARTESIAN_GRID) THEN
151                    IF(SCALAR_NODE_ATWALL(CUR_IJK))                   &
152                       COUNT_NODES_OUTSIDE = COUNT_NODES_OUTSIDE + 1
153                 ENDIF
154     
155              ENDDO
156              ENDDO
157              ENDDO
158     
159     
160     ! Calculate des_rops_node so des_rop_s, and in turn, ep_g can be updated
161     !----------------------------------------------------------------->>>
162     
163     ! looping through particles in the cell
164              DO NINDX=1, PINC(IJK)
165                 NP = PIC(IJK)%P(NINDX)
166     
167                 call DRAG_WEIGHTFACTOR(gst_tmp,des_pos_new(:,np),weight_ft)
168     
169                 M = PIJK(NP,5)
170                 WTP = ONE
171                 IF(MPPIC) WTP = DES_STAT_WT(NP)
172     
173                 MASS_SOL1 = MASS_SOL1 + PMASS(NP)*WTP
174     
175                 TEMP2 = DES_RO_S(M)*PVOL(NP)*WTP
176     
177                 DO K = 1, merge(1, ONEW, NO_K)
178                 DO J = 1, ONEW
179                 DO I = 1, ONEW
180     ! shift loop index to new variables for manipulation
181                    II = IW + I-1
182                    JJ = JS + J-1
183                    KK = KB + K-1
184     
185                    CUR_IJK = FUNIJK_MAP_C(II,JJ,KK)
186     
187                    TEMP1 = WEIGHT_FT(I,J,K)*TEMP2
188     
189                    DES_ROPS_NODE(CUR_IJK,M) = DES_ROPS_NODE(CUR_IJK,M) +   &
190                       TEMP1
191     
192                    DES_VEL_NODE(CUR_IJK,:,M) = DES_VEL_NODE(CUR_IJK,:,M) + &
193                       TEMP1*DES_VEL_NEW(:,NP)
194                 ENDDO
195                 ENDDO
196                 ENDDO
197              ENDDO   ! end do (nindx=1,pinc(ijk))
198     !-----------------------------------------------------------------<<<
199     
200     
201     ! Only for cutcell cases may count_nodes_inside become less than its
202     ! original set value. In such an event, the contribution of scalar nodes
203     ! that do not reside in the domain is added to a residual array. This
204     ! array is then redistribited equally to the nodes that are in the fluid
205     ! domain. These steps are done to conserve mass.
206     !----------------------------------------------------------------->>>
207              IF (CARTESIAN_GRID) THEN
208     
209     ! only for cartesian_grid will count_nodes_outside be modified from zero
210                 COUNT_NODES_INSIDE = &
211                    COUNT_NODES_INSIDE_MAX - COUNT_NODES_OUTSIDE
212     
213                 IF(COUNT_NODES_INSIDE.LT.COUNT_NODES_INSIDE_MAX) THEN
214     
215     ! initializing
216                    RESID_ROPS(1:DES_MMAX) = ZERO
217                    RESID_VEL(:, 1:DES_MMAX) = ZERO
218     
219     ! Convention used to number node numbers
220     ! i=1, j=2           i=2, j=2
221     !   _____________________
222     !   |                   |
223     !   |  I = 2, J = 2     |
224     !   |___________________|
225     ! i=1, j=1           i=2, j=1
226     ! setting indices based on convention
227                    I = I_OF(IJK)
228                    J = J_OF(IJK)
229                    K = K_OF(IJK)
230                    I1 = I-1
231                    I2 = I
232                    J1 = J-1
233                    J2 = J
234                    K1 = merge(K, K-1, NO_K)
235                    K2 = K
236     ! first calculate the residual des_rops_node and des_vel_node that was
237     ! computed on nodes that do not belong to the domain
238     
239                    DO KK = K1, K2
240                    DO JJ = J1, J2
241                    DO II = I1, I2
242     
243                       IJK2 = funijk(II, JJ, KK)
244     
245                       IF(SCALAR_NODE_ATWALL(IJK2)) THEN
246                          RESID_ROPS(1:DES_MMAX) = RESID_ROPS(1:DES_MMAX) + &
247                             DES_ROPS_NODE(IJK2,1:DES_MMAX)
248     
249                          DES_ROPS_NODE(IJK2,1:DES_MMAX) = ZERO
250                          DO IDIM = 1, merge(2,3,NO_K)
251                             RESID_VEL(IDIM, 1:DES_MMAX) =                  &
252                                 RESID_VEL(IDIM, 1:DES_MMAX) +              &
253                                DES_VEL_NODE(IJK2,IDIM, 1:DES_MMAX)
254                             DES_VEL_NODE(IJK2,IDIM, 1:DES_MMAX) = ZERO
255                          ENDDO
256                       ENDIF
257                    ENDDO
258                    ENDDO
259                    ENDDO
260     
261     ! now add this residual equally to the remaining nodes
262                    NORM_FACTOR = ONE/REAL(COUNT_NODES_INSIDE)
263                    DO KK = K1, K2
264                    DO JJ = J1, J2
265                    DO II = I1, I2
266                       IJK2 = funijk(II, JJ, KK)
267     
268                       IF(.NOT.SCALAR_NODE_ATWALL(IJK2)) THEN
269                          DES_ROPS_NODE(IJK2,1:DES_MMAX) =                  &
270                             DES_ROPS_NODE(IJK2,1:DES_MMAX) +               &
271                             RESID_ROPS(1:DES_MMAX)*NORM_FACTOR
272                          DO IDIM = 1, merge(2,3,NO_K)
273                             DES_VEL_NODE(IJK2,IDIM, 1:DES_MMAX) =          &
274                                DES_VEL_NODE(IJK2,IDIM, 1:DES_MMAX) +       &
275                                RESID_VEL(IDIM, 1:DES_MMAX)*NORM_FACTOR
276                          ENDDO
277                       ENDIF
278     
279                    ENDDO
280                    ENDDO
281                    ENDDO
282                 ENDIF
283              ENDIF   ! end if (cartesian_grid)
284           ENDDO
285     !$omp end parallel
286     
287     
288     ! At the interface des_rops_node has to be added since particles
289     ! across the processors will contribute to the same scalar node.
290     ! sendrecv will be called and the node values will be added
291     ! at the junction. des_rops_node is altered by the routine when
292     ! periodic boundaries are invoked
293           CALL DES_ADDNODEVALUES_MEAN_FIELDS
294     
295     
296     ! Now go from node to scalar center. Same convention as sketched
297     ! earlier
298     !----------------------------------------------------------------->>>
299     ! Explanation by RG: 08/17/2012
300     ! the approach used here is to make it general enough for cutcells to be
301     ! included as well. The new changes do not alter earlier calculations
302     ! but make the technique general as to include cartesian grid (cut-cell)
303     ! simulations.
304     ! Previously, the volume of the node (by array des_vol_node) was used to
305     ! first scale the nodal values. Subsequently, these nodal values were
306     ! equally added to compute the cell centered values for the scalar cell.
307     
308     ! Consider an internal node next to an edge node (a node adjacent to a
309     ! boundary). In 2D, the volume of an edge node will be half that of an
310     ! internal node. And an edge node will contribute double compared to
311     ! an internal node to the value of the scalar cell they share. These
312     ! calculations were previously accomplished via the variable volume of
313     ! node.  Now this is accomplished by the ratio vol(ijk2)/vol_sur, where
314     ! vol(ijk2) is the volume of the scalar cell in consideration and
315     ! vol_sur is the sum of all the scalar cell volumes that have this node
316     ! as the common node.
317     !---------------------------------------------------------------------//
318     !$omp parallel do default(none) collapse (3)                           &
319     !$omp shared(KSTART2, KEND1, JSTART2, JEND1, ISTART2, IEND1, DO_K, VOL,&
320     !$omp   DEAD_CELL_AT, FUNIJK_MAP_C, VOL_SURR, DES_MMAX, DES_ROPS_NODE, &
321     !$omp   DES_VEL_NODE)                                                  &
322     !$omp private(I, J, K, IJK, M, II, JJ, KK, IJK2, DES_ROP_DENSITY,      &
323     !$omp   DES_VEL_DENSITY)                                               &
324     !$omp reduction(+:DES_ROP_S, DES_U_S, DES_V_S, DES_W_S)
325           DO K = KSTART2, KEND1
326           DO J = JSTART2, JEND1
327           DO I = ISTART2, IEND1
328              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
329              IJK = funijk(I,J,K)
330              if (VOL_SURR(IJK).eq.ZERO) CYCLE ! no FLUID_AT any of the stencil points have
331     
332     ! looping over stencil points (NODE VALUES)
333              DO M = 1, DES_MMAX
334     
335                 DES_ROP_DENSITY = DES_ROPS_NODE(IJK, M)/VOL_SURR(IJK)
336                 DES_VEL_DENSITY(:) = DES_VEL_NODE(IJK, :, M)/VOL_SURR(IJK)
337     
338                 DO KK = K, merge(K+1, K, DO_K)
339                 DO JJ = J, J+1
340                 DO II = I, I+1
341                    IF (DEAD_CELL_AT(II,JJ,KK)) CYCLE  ! skip dead cells
342     
343                    IJK2 = funijk_map_c(II, JJ, KK)
344                    IF(FLUID_AT(IJK2).and.(IS_ON_myPE_wobnd(II, JJ, KK))) THEN
345     ! Since the data in the ghost cells is spurious anyway and overwritten during
346     ! subsequent send receives, do not compute any value here as this will
347     ! mess up the total mass value that is computed below to ensure mass conservation
348     ! between Lagrangian and continuum representations
349                       DES_ROP_S(IJK2, M) = DES_ROP_S(IJK2, M) + DES_ROP_DENSITY*VOL(IJK2)
350                       DES_U_S(IJK2, M) = DES_U_S(IJK2, M) + DES_VEL_DENSITY(1)*VOL(IJK2)
351                       DES_V_S(IJK2, M) = DES_V_S(IJK2, M) + DES_VEL_DENSITY(2)*VOL(IJK2)
352                       IF(DO_K) DES_W_S(IJK2, M) = DES_W_S(IJK2, M) + DES_VEL_DENSITY(3)*VOL(IJK2)
353                    ENDIF
354                 ENDDO  ! end do (ii=i1,i2)
355                 ENDDO  ! end do (jj=j1,j2)
356                 ENDDO  ! end do (kk=k1,k2)
357              ENDDO
358           ENDDO   ! end do (i=istart2,iend1)
359           ENDDO   ! end do (j=jstart2,jend1)
360           ENDDO   ! end do (k=kstart2,kend1)
361     !omp end parallel do
362     
363     !-----------------------------------------------------------------<<<
364     
365     
366     !$omp parallel do default(none) private(IJK, M)                        &
367     !$omp shared(IJKSTART3, IJKEND3, DO_K, DES_MMAX, DES_ROP_s, DES_U_S,   &
368     !$omp   DES_V_S, DES_W_S, VOL)
369           DO IJK = IJKSTART3, IJKEND3
370              IF(.NOT.FLUID_AT(IJK)) CYCLE
371     
372              DO M = 1, DES_MMAX
373                 IF(DES_ROP_S(IJK, M).GT.ZERO) THEN
374                    DES_U_S(IJK, M) = DES_U_S(IJK,M)/DES_ROP_S(IJK, M)
375                    DES_V_S(IJK, M) = DES_V_S(IJK,M)/DES_ROP_S(IJK, M)
376                    IF(DO_K) DES_W_S(IJK, M) = DES_W_S(IJK,M)/DES_ROP_S(IJK, M)
377     
378     ! Divide by scalar cell volume to obtain the bulk density
379                    DES_ROP_S(IJK, M) = DES_ROP_S(IJK, M)/VOL(IJK)
380     
381                 ENDIF
382              ENDDO   ! end loop over M=1,DES_MMAX
383           ENDDO  ! end loop over IJK=ijkstart3,ijkend3
384     !omp end parallel do
385     
386     
387           IF (MPPIC) CALL SEND_RECV(DES_ROP_S,2)
388     
389           CALL CALC_EPG_DES
390     
391     
392           IF(MPPIC) THEN
393     
394     ! Now calculate Eulerian mean velocity fields like U_S, V_S, and W_S.
395              CALL SEND_RECV(DES_U_S,2)
396              CALL SEND_RECV(DES_V_S,2)
397              IF(DO_K) CALL SEND_RECV(DES_W_S,2)
398     
399     ! The Eulerian velocity field is used to set up the stencil to interpolate
400     ! mean solid velocity at the parcel's location. DES_U_S could have also been
401     ! used, but that also would have require the communication at this stage.
402     ! The final interpolated value does not change if the stencil is formed by
403     ! first obtaining face centered Eulerian velocities (U_S, etc.)
404     ! and then computing the node velocities from them or directly computing
405     ! the node velocities from cell centered average velocity field (DES_U_S,
406     ! etc.). We are using the first approach as it is more natural to set
407     ! BC's on solid velocity field in the face centered represenation (U_S,
408     ! etc.)
409     
410              IF(.NOT.CARTESIAN_GRID) THEN
411                 CALL MPPIC_COMP_EULERIAN_VELS_NON_CG
412              ELSE
413                 CALL MPPIC_COMP_EULERIAN_VELS_CG
414              ENDIF
415           ENDIF   ! end if (.not.mppic)
416     
417     ! turn on the below statements to check if the mass is conserved
418     ! between discrete and continuum representations. Should be turned to
419     ! false for any production runs.
420           IF(DES_REPORT_MASS_INTERP) THEN
421     
422     
423              DO IJK = IJKSTART3, IJKEND3
424                 IF(.NOT.FLUID_AT(IJK)) CYCLE
425     
426                 I = I_OF(IJK)
427                 J = J_OF(IJK)
428                 K = K_OF(IJK)
429     
430     ! It is important to check both FLUID_AT and IS_ON_MYPE_WOBND.
431                 IF(IS_ON_myPE_wobnd(I,J,K)) MASS_SOL2 = MASS_SOL2 +        &
432                    sum(DES_ROP_S(IJK,1:DES_MMAX))*VOL(IJK)
433              ENDDO
434     
435     
436              CALL GLOBAL_SUM(MASS_SOL1, MASS_SOL1_ALL)
437              CALL GLOBAL_SUM(MASS_SOL2, MASS_SOL2_ALL)
438              if(myPE.eq.pe_IO) THEN
439                 WRITE(*,'(/,5x,A,4(2x,g17.8),/)') &
440                      'SOLIDS MASS DISCRETE AND CONTINUUM =  ', &
441                      MASS_SOL1_ALL, MASS_SOL2_ALL
442              ENDIF
443           ENDIF
444           END SUBROUTINE COMP_MEAN_FIELDS0
445     
446     
447     
448     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
449     !  Subroutine: DRAG_WEIGHTFACTOR                                        C
450     !  Purpose: DES - Calculate the fluid velocity interpolated at the      C
451     !           particle's location and weights. Replace 'interpolator'     C
452     !                       interface for OpenMP implementation.            C
453     !                                                                       C
454     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
455     
456           SUBROUTINE DRAG_WEIGHTFACTOR(GSTEN,DESPOS,WEIGHTFACTOR)
457     
458           use geometry, only: NO_K
459     
460             IMPLICIT NONE
461     
462     !-----------------------------------------------
463     ! Local Variables
464     !-----------------------------------------------
465             DOUBLE PRECISION, DIMENSION(2,2,2,3), INTENT(IN):: GSTEN
466             DOUBLE PRECISION, DIMENSION(3), INTENT(IN):: DESPOS
467             DOUBLE PRECISION, DIMENSION(2,2,2), INTENT(OUT) :: WEIGHTFACTOR
468             INTEGER :: II, JJ, KK
469     
470             DOUBLE PRECISION, DIMENSION(2) :: XXVAL, YYVAL, ZZVAL
471             DOUBLE PRECISION :: DXX, DYY, DZZ
472             DOUBLE PRECISION, DIMENSION(3) :: ZETAA
473     
474             DXX = GSTEN(2,1,1,1) - GSTEN(1,1,1,1)
475             DYY = GSTEN(1,2,1,2) - GSTEN(1,1,1,2)
476     
477             ZETAA(1:2) = DESPOS(1:2) - GSTEN(1,1,1,1:2)
478     
479             ZETAA(1) = ZETAA(1)/DXX
480             ZETAA(2) = ZETAA(2)/DYY
481     
482             XXVAL(1)=1-ZETAA(1)
483             YYVAL(1)=1-ZETAA(2)
484             XXVAL(2)=ZETAA(1)
485             YYVAL(2)=ZETAA(2)
486     
487             IF(NO_K) THEN
488                DO JJ=1,2
489                   DO II=1,2
490                      WEIGHTFACTOR(II,JJ,1) = XXVAL(II)*YYVAL(JJ)
491                   ENDDO
492                ENDDO
493             ELSE
494                DZZ = GSTEN(1,1,2,3) - GSTEN(1,1,1,3)
495                ZETAA(3) = DESPOS(3) - GSTEN(1,1,1,3)
496                ZETAA(3) = ZETAA(3)/DZZ
497                ZZVAL(1)=1-ZETAA(3)
498                ZZVAL(2)=ZETAA(3)
499                DO KK=1,2
500                   DO JJ=1,2
501                      DO II=1,2
502                         WEIGHTFACTOR(II,JJ,KK) = XXVAL(II)*YYVAL(JJ)*ZZVAL(KK)
503                      ENDDO
504                   ENDDO
505                ENDDO
506             ENDIF
507     
508           END SUBROUTINE DRAG_WEIGHTFACTOR
509