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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
4           SUBROUTINE COMP_MEAN_FIELDS1
5     
6     !-----------------------------------------------
7     ! Modules
8     !-----------------------------------------------
9           USE param
10           USE param1
11           USE fldvar
12           USE geometry
13           USE indices
14           USE compar
15           USE parallel
16           USE sendrecv
17           USE discretelement
18           use desgrid
19           use desmpi
20           USE mfix_pic
21           USE functions
22           use particle_filter, only: FILTER_WEIGHT
23           use particle_filter, only: FILTER_CELL
24           IMPLICIT NONE
25     !-----------------------------------------------
26     ! Local variables
27     !-----------------------------------------------
28     ! Loop counters: partciles, filter cells, phases
29           INTEGER NP, LC, M
30     ! Fluid cell index
31           INTEGER IJK
32     ! Total Mth solids phase volume in IJK
33           DOUBLE PRECISION :: SOLVOLINC(DIMENSION_3,DES_MMAX)
34     ! One divided by the total solids volume.
35           DOUBLE PRECISION :: OoSOLVOL
36     ! PVOL times statistical weight, and times filter weight
37           DOUBLE PRECISION :: VOL_WT, VOLxWEIGHT
38     ! Loop bound for filter
39           INTEGER :: LP_BND
40     
41     
42     !-----------------------------------------------
43     
44     
45           SOLVOLINC(:,:) = ZERO
46     
47           IF(MPPIC) THEN
48              DES_U_s(:,:) = ZERO
49              DES_V_s(:,:) = ZERO
50              DES_W_s(:,:) = ZERO
51           ENDIF
52     
53     ! Loop bounds for interpolation.
54           LP_BND = merge(27,9,DO_K)
55     
56     ! Calculate the gas phase forces acting on each particle.
57     !$omp parallel default(none)                                           &
58     !$omp private(NP, VOL_WT, M, LC, IJK, VOLXWEIGHT)                      &
59     !$omp shared(MAX_PIP, PEA, PVOL, DES_STAT_WT, PIJK, LP_BND, MPPIC,     &
60     !$omp    FILTER_WEIGHT, SOLVOLINC,DES_U_S, DES_V_S, DES_W_S, DO_K,     &
61     !$omp    FILTER_CELL,DES_VEL_NEW)
62     !$omp do
63           do NP=1,MAX_PIP
64              IF(.NOT.PEA(NP,1)) CYCLE
65              IF(any(PEA(NP,2:3))) CYCLE
66     
67              VOL_WT = PVOL(NP)
68              IF(MPPIC) VOL_WT = VOL_WT*DES_STAT_WT(NP)
69     ! Particle phase for data binning.
70              M = PIJK(NP,5)
71     
72              DO LC=1,LP_BND
73                 IJK = FILTER_CELL(LC,NP)
74     ! Particle volume times the weight for this cell.
75                 VOLxWEIGHT = VOL_WT*FILTER_WEIGHT(LC,NP)
76     ! Accumulate total solids volume (by phase)
77                 !$omp atomic
78                 SOLVOLINC(IJK,M) = SOLVOLINC(IJK,M) + VOLxWEIGHT
79                 IF(MPPIC) THEN
80     ! Accumulate total solids momentum (by phase)
81                    !$omp atomic
82                    DES_U_S(IJK,M) = DES_U_S(IJK,M) + &
83                       DES_VEL_NEW(1,NP)*VOLxWEIGHT
84                    !$omp atomic
85                    DES_V_S(IJK,M) = DES_V_S(IJK,M) + &
86                       DES_VEL_NEW(2,NP)*VOLxWEIGHT
87                    IF(DO_K) THEN
88                       !$omp atomic
89                       DES_W_S(IJK,M) = DES_W_S(IJK,M) + &
90                          DES_VEL_NEW(3,NP)*VOLxWEIGHT
91                    ENDIF
92                 ENDIF
93              ENDDO
94           ENDDO
95     !$omp end do
96     !$omp end parallel
97     
98     ! Calculate the cell average solids velocity, the bulk density,
99     ! and the void fraction.
100     !---------------------------------------------------------------------//
101     !$omp parallel do if(ijkend3 .ge. 2000) default(shared)                &
102     !$omp private(IJK,M,OoSOLVOL)
103           DO IJK = IJKSTART3, IJKEND3
104              IF(.NOT.FLUID_AT(IJK)) CYCLE
105     
106     ! calculating the cell average solids velocity for each solids phase
107              DO M = 1, DES_MMAX
108                 IF(SOLVOLINC(IJK,M).GT.ZERO) THEN
109                    OoSOLVOL = ONE/SOLVOLINC(IJK,M)
110                    DES_U_s(IJK,M) = DES_U_s(IJK,M)*OoSOLVOL
111                    DES_V_s(IJK,M) = DES_V_s(IJK,M)*OoSOLVOL
112                    IF(DO_K) DES_W_s(IJK,M) = DES_W_s(IJK,M)*OoSOLVOL
113                 ENDIF
114     
115     ! calculating the bulk density of solids phase m based on the total
116     ! number of particles having their center in the cell
117                 DES_ROP_S(IJK,M) = DES_RO_S(M)*SOLVOLINC(IJK,M)/VOL(IJK)
118     
119              ENDDO   ! end loop over M=1,DES_MMAX
120     
121           ENDDO     ! end loop over IJK=ijkstart3,ijkend3
122     !$omp end parallel do
123     
124     
125     ! Halo exchange of solids volume fraction data.
126           calL SEND_RECV(DES_ROP_S,2)
127     
128           end SUBROUTINE COMP_MEAN_FIELDS1
129