File: N:\mfix\model\des\des_thermo_conv.f

1     #include "version.inc"
2     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
3     !                                                                      !
4     !  Module name: CONV_GS_DES1                                           !
5     !  Author: J.Musser: 16-Jun-10                                         !
6     !                                                                      !
7     !  Purpose: This routine is called from the DISCRETE side to calculate !
8     !  the gas-particle convective heat transfer.                          !
9     !                                                                      !
10     !  Comments: Explicitly coupled simulations use a stored convective    !
11     !  heat transfer coefficient. Otherwise, the convective heat transfer  !
12     !  coeff is calculated every time step and the total interphase energy !
13     !  transfered is 'stored' and used explictly in the gas phase. The     !
14     !  latter conserves all energy
15     !                                                                      !
16     !  REF: Zhou, Yu, and Zulli, "Particle scale study of heat transfer in !
17     !       packed and bubbling fluidized beds," AIChE Journal, Vol. 55,   !
18     !       no 4, pp 868-884, 2009.                                        !
19     !                                                                      !
20     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
21           SUBROUTINE CONV_GS_DES1
22     
23           use constant, only: Pi
24           Use des_thermo
25           Use discretelement
26           Use fldvar
27           Use interpolation
28           Use param1
29           use des_thermo, only: GAMMAxSA
30           use geometry, only: NO_K
31           use particle_filter, only: DES_INTERP_ON
32           use particle_filter, only: FILTER_CELL
33           use particle_filter, only: FILTER_WEIGHT
34           use particle_filter, only: FILTER_SIZE
35           use functions, only: FLUID_AT
36           use functions, only: IS_NORMAL
37           IMPLICIT NONE
38     
39           DOUBLE PRECISION :: lTg, GAMMA
40           DOUBLE PRECISION :: Qcv_DT, Qcv
41           DOUBLE PRECISION :: l4Pi
42           INTEGER :: IJK, LC, NP
43     
44           l4Pi = 4.0d0*Pi
45     
46           DO NP=1,MAX_PIP
47              IF(.NOT.IS_NORMAL(NP)) CYCLE
48     
49     ! Calculate the gas temperature.
50              IF(DES_INTERP_ON) THEN
51                 lTg = ZERO
52                 DO LC=1,FILTER_SIZE
53                    IJK = FILTER_CELL(LC,NP)
54                    lTg = lTg + T_G(IJK)*FILTER_WEIGHT(LC,NP)
55                 ENDDO
56              ELSE
57                 IJK = PIJK(NP,4)
58                 lTg = T_G(IJK)
59              ENDIF
60     
61              IJK = PIJK(NP,4)
62     
63     ! Avoid convection calculations in cells without fluid (cut-cell)
64              IF(.NOT.FLUID_AT(IJK)) THEN
65                 GAMMAxSA(NP) = ZERO
66                 CONV_Qs(NP) = ZERO
67     
68     ! For explicit coupling, use the heat transfer coefficient calculated
69     ! for the gas phase heat transfer calculations.
70              ELSEIF(DES_EXPLICITLY_COUPLED) THEN
71                 CONV_Qs(NP) = GAMMAxSA(NP)*(lTg - DES_T_s(NP))
72     
73              ELSE
74     
75     ! Calculate the heat transfer coefficient.
76                 CALL CALC_GAMMA_DES(NP, GAMMA)
77                 GAMMAxSA(NP) = GAMMA* l4Pi*DES_RADIUS(NP)*DES_RADIUS(NP)
78     
79     ! Calculate the rate of heat transfer to the particle
80                 Qcv = GAMMAxSA(NP)*(lTg - DES_T_s(NP))
81     ! Store convection source in global energy source array.
82                 Q_Source(NP) = Q_Source(NP) + Qcv
83     
84     ! Calculate the gas phase source term components.
85                 Qcv_DT = Qcv*DTSOLID
86                 IF(DES_INTERP_ON) THEN
87                    DO LC=1,FILTER_SIZE
88                       CONV_SC(IJK)=CONV_Sc(IJK)-Qcv_DT*FILTER_WEIGHT(LC,NP)
89                    ENDDO
90                 ELSE
91                    CONV_SC(IJK) = CONV_Sc(IJK) - Qcv_DT
92                 ENDIF
93              ENDIF
94           ENDDO
95     
96     ! Note that MPI sync is managed at the end of des_time_march for
97     ! non-explicitly coupled cases that use interpolation.
98     
99           RETURN
100           END SUBROUTINE CONV_GS_DES1
101     
102     
103     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
104     !                                                                      !
105     !  Subroutine: CONV_GS_GAS1                                            !
106     !  Author: J.Musser                                   Date: 21-NOV-14  !
107     !                                                                      !
108     !                                                                      !
109     !  Purpose: This routine is called from the CONTINUUM. It calculates   !
110     !  the scalar cell center drag force acting on the fluid using         !
111     !  interpolated values for the gas velocity and volume fraction. The   !
112     !  The resulting sources are interpolated back to the fluid grid.      !
113     !                                                                      !
114     !  NOTE: The loop over particles includes ghost particles so that MPI  !
115     !  communications are needed to distribute overlapping force between   !
116     !  neighboring grid cells. This is possible because only cells "owned" !
117     !  by the current process will have non-zero weights.                  !
118     !                                                                      !
119     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
120           SUBROUTINE CONV_GS_GAS1
121     
122     ! Flag: The fluid and discrete solids are explicitly coupled.
123           use discretelement, only: DES_EXPLICITLY_COUPLED
124     ! Size of particle array on this process.
125           use discretelement, only: MAX_PIP
126     ! Flag to use interpolation
127           use particle_filter, only: DES_INTERP_ON
128     ! Interpolation cells and weights
129           use particle_filter, only: FILTER_CELL, FILTER_WEIGHT, FILTER_SIZE
130     ! IJK of fluid cell containing particles center
131           use discretelement, only: PIJK
132     ! Particle temperature
133           use des_thermo, only: DES_T_s
134     ! Gas phase energy equation sources
135           use des_thermo, only: CONV_Sp, CONV_Sc
136           Use discretelement, only: DES_RADIUS
137     ! Heat transfer coefficint (GAMMA) multiplied by sufrace area
138           use des_thermo, only: GAMMAxSA
139     ! Funtion for identifying fluid cells and normal particles.
140           use functions, only: FLUID_AT
141           use functions, only: IS_NORMAL
142     ! MPI function for collecting interpolated data from ghost cells.
143           use sendrecvnode, only: DES_COLLECT_gDATA
144     ! MPI wrapper for halo exchange.
145           use sendrecv, only: SEND_RECV
146     
147     ! Global Parameters:
148     !---------------------------------------------------------------------//
149     ! Double precision values.
150           use param1, only: ZERO, ONE
151           use constant, only: Pi
152     
153           IMPLICIT NONE
154     
155     ! Loop counters: Particle, fluid cell, neighbor cells
156           INTEGER :: NP, IJK, LC
157     ! Interpolation weight
158           DOUBLE PRECISION :: WEIGHT
159           DOUBLE PRECISION :: GAMMAxSAxTp, GAMMA
160           DOUBLE PRECISION :: l4Pi
161     
162           l4Pi = 4.0d0*Pi
163     
164     ! Initialize fluid cell values.
165           CONV_Sc = ZERO
166           CONV_Sp = ZERO
167     
168     ! Calculate the gas phase forces acting on each particle.
169           DO NP=1,MAX_PIP
170     
171              IF(.NOT.IS_NORMAL(NP)) CYCLE
172              IF(.NOT.FLUID_AT(PIJK(NP,4))) CYCLE
173     
174     ! Calculate the heat transfer coefficient.
175              CALL CALC_GAMMA_DES(NP, GAMMA)
176     
177     ! Calculate the surface area of the particle
178     
179              GAMMAxSA(NP) = GAMMA*l4Pi*DES_RADIUS(NP)*DES_RADIUS(NP)
180              GAMMAxSAxTp = GAMMAxSA(NP)*DES_T_s(NP)
181     
182              IF(DES_INTERP_ON) THEN
183                 DO LC=1,FILTER_SIZE
184                    IJK = FILTER_CELL(LC,NP)
185                    WEIGHT = FILTER_WEIGHT(LC,NP)
186     
187                    CONV_Sc(IJK) = CONV_Sc(IJK) + WEIGHT*GAMMAxSAxTp
188                    CONV_Sp(IJK) = CONV_Sp(IJK) + WEIGHT*GAMMAxSA(NP)
189                 ENDDO
190              ELSE
191                 IJK = PIJK(NP,4)
192     
193                 CONV_Sc(IJK) = CONV_Sc(IJK) + GAMMAxSAxTp
194                 CONV_Sp(IJK) = CONV_Sp(IJK) + GAMMAxSA(NP)
195              ENDIF
196     
197           ENDDO
198     
199     ! Add in data stored in ghost cells from interpolation. This call must
200     ! preceed the SEND_RECV to avoid overwriting ghost cell data.
201           IF(DES_INTERP_ON) THEN
202              CALL DES_COLLECT_gDATA(CONV_SC)
203              CALL DES_COLLECT_gDATA(CONV_SP)
204           ENDIF
205     
206     ! Update the drag force and sources in ghost layers.
207           CALL SEND_RECV(CONV_SC, 2)
208           CALL SEND_RECV(CONV_SP, 2)
209     
210           RETURN
211           END SUBROUTINE CONV_GS_GAS1
212     
213     
214     
215     
216     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
217     !  Subroutine: ZERO_ENERGY_SOURCE                                      !
218     !                                                                      !
219     !  Purpose: ZERO out the array that passes energy source terms back to !
220     !  the continuum model. Additional entries may be needed to include    !
221     !  heat transfer to the hybrid mode.                                   !
222     !                                                                      !
223     !  Author: J.Musser                                   Date: 15-Jan-11  !
224     !                                                                      !
225     !  Comments:                                                           !
226     !                                                                      !
227     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
228           SUBROUTINE ZERO_ENERGY_SOURCE
229     
230           Use des_thermo
231           Use param1
232     
233           IMPLICIT NONE
234     
235           CONV_Sc = ZERO
236           CONV_Sp = ZERO
237     
238           RETURN
239           END SUBROUTINE ZERO_ENERGY_SOURCE
240     
241