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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  Module name: CALC_RRATES_DES                                        !
4     !  Author: J. Musser                                  Date: 10-Oct-12  !
5     !                                                                      !
6     !  Purpose: Calculate reaction rates for various reactions for DES.    !
7     !                                                                      !
8     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
9           SUBROUTINE CALC_RRATES_DES(NP, pRgp, pRgc, pRPhase, pHoRg, pSUMRg)
10     
11           USE compar
12           USE constant
13           USE des_rxns
14           USE des_thermo
15           USE discretelement
16           USE energy
17           USE fldvar
18           USE funits
19           USE geometry
20           USE indices
21           USE parallel
22           USE param
23           USE param1
24           USE physprop
25           USE run
26           USE rxns
27           USE sendrecv
28           USE usr
29           Use parse
30           use functions
31           use toleranc, only: ZERO_X_gs, COMPARE
32           use physprop, only: NMAX
33           use param1, only: DIMENSION_LM
34           use param, only: DIMENSION_M
35     
36           IMPLICIT NONE
37     
38     ! Passed variables
39     !---------------------------------------------------------------------//
40           INTEGER, INTENT(IN) :: NP   ! particle index
41     
42     ! Local gas phase values.
43           DOUBLE PRECISION, INTENT(OUT) :: pRgp(NMAX(0)) ! Rate of production
44           DOUBLE PRECISION, INTENT(OUT) :: pRgc(NMAX(0)) ! Rate of consumption
45           DOUBLE PRECISION, INTENT(OUT) :: pHoRg   ! Heat of reaction
46           DOUBLE PRECISION, INTENT(OUT) :: pSUMRg
47           DOUBLE PRECISION, INTENT(OUT) :: pRPhase(DIMENSION_LM+DIMENSION_M-1)
48     
49     ! Local variables
50     !---------------------------------------------------------------------//
51           INTEGER :: pM   ! Global Solids Phase index
52           INTEGER :: IJK  ! fluid cell index
53     
54           INTEGER :: H    ! Reaction loop counter
55           INTEGER :: M    ! Global Phase index loop counter
56           INTEGER :: NN   ! Global species index
57           INTEGER :: lN   ! Local reaction speices index/loop counter
58           INTEGER :: LM   !
59     
60           INTEGER :: mXfr ! Global phase index for mass transfer
61     
62     ! User-defined reaction rates returned from USR_RATES
63           DOUBLE PRECISION :: DES_RATES(NO_OF_DES_RXNS)
64     
65           DOUBLE PRECISION :: lRate
66           DOUBLE PRECISION :: lTp
67           DOUBLE PRECISION :: lHoRs, llHoRg
68     
69           DOUBLE PRECISION :: RxH
70     
71     ! Reaction limiters. If a species mass fraction is less than this
72     ! value, then the reaction is suppressed.
73           DOUBLE PRECISION :: speciesLimiter
74     
75     ! External functions
76     !---------------------------------------------------------------------//
77     ! Enthalpy calculations (cal/gram)
78           DOUBLE PRECISION, EXTERNAL :: CALC_H
79     
80     ! Alias particle temperature.
81           lTp = DES_T_s(NP)
82           IJK = PIJK(NP,4)
83           pM = PIJK(NP,5)
84     
85     ! Initialize storage arrays
86           pRgp(:) = ZERO
87           pRgc(:) = ZERO
88           pHoRg = ZERO
89     
90           pRPhase = ZERO
91           pSUMRg = ZERO
92     
93     ! Set the species limiter:
94           speciesLimiter = ZERO_X_gs
95     
96     ! Calculate user defined reaction rates.
97           DES_RATES(:) = ZERO
98           CALL USR_RATES_DES(NP, pM, IJK, DES_RATES)
99     
100     ! Loop over reactions.
101           RXN_LP: DO H = 1, NO_OF_DES_RXNS
102     
103     ! Skip empty reactions
104              IF(DES_Reaction(H)%nSpecies == 0) CYCLE RXN_LP
105              IF(COMPARE(DES_RATES(H),ZERO)) CYCLE RXN_LP
106     
107     ! Initialize local loop arrays
108              llHoRg = ZERO
109              lHoRs = ZERO
110              RxH = ZERO
111     
112     ! Calculate the rate of formation/consumption for each species.
113     !---------------------------------------------------------------------//
114              DO lN = 1, DES_Reaction(H)%nSpecies
115     ! Global phase index.
116                 M = DES_Reaction(H)%Species(lN)%pMap
117     ! Global species index.
118                 NN = DES_Reaction(H)%Species(lN)%sMap
119     ! Index for interphase mass transfer. For a gas/solid reaction, the
120     ! index is stored with the gas phase.
121                 mXfr = DES_Reaction(H)%Species(lN)%mXfr
122                 lRate = DES_RATES(H) * DES_Reaction(H)%Species(lN)%MWxStoich
123     ! Gas Phase:
124                 IF(M == 0) THEN
125     ! Consumption of gas phase species.
126                    IF(lRate < ZERO) THEN
127                       IF(X_g(IJK,NN) > speciesLimiter) THEN
128                          pRgc(NN) = pRgc(NN) - lRate
129     ! Enthalpy transfer associated with mass transfer. (gas/solid)
130                          IF(M /= mXfr) RxH = RxH +                         &
131                             lRate*CALC_H(T_g(IJK),0,NN)
132                       ELSE
133     ! There is an insignificant amount of reactant. Skip this reaction.
134                          DES_RATES(H) = ZERO
135                          CYCLE RXN_LP
136                       ENDIF
137                    ELSE
138     ! Formation of gas phase species.
139                       pRgp(nn) = pRgp(nn) + lRate
140     ! Enthalpy transfer associated with mass transfer. (gas/solid)
141                       IF(M /= mXfr) RxH = RxH + lRate*CALC_H(lTp,0,NN)
142                    ENDIF
143     ! Discrete Solids Phase:
144                 ELSE
145     ! Formation/consumption of solids phase species.
146                    DES_R_s(NP,nn) = DES_R_s(NP,nn) + lRate
147                 ENDIF
148              ENDDO ! Loop of species
149     
150     
151     ! Calculate and store the heat of reaction.
152     !---------------------------------------------------------------------//
153              IF(ENERGY_EQ) THEN
154     ! Automated heat of reaction calculations
155                 IF(DES_Reaction(H)%Calc_DH) THEN
156     ! Loop over reaction species.
157                    DO lN = 1, DES_Reaction(H)%nSpecies
158     ! Global phase index.
159                       M = DES_Reaction(H)%Species(lN)%pMap
160     ! Global species index.
161                       NN = DES_Reaction(H)%Species(lN)%sMap
162     ! Rate of formation/consumption for speices N
163                       lRate = DES_RATES(H) * &
164                          DES_Reaction(H)%Species(lN)%MWxStoich
165     ! Gas phase enthalpy chnage from energy equation derivation.
166                       IF(M == 0) THEN
167                          llHORg = llHORg + CALC_H(T_g(IJK),0,nn) * lRate
168     ! Solid phase enthalpy change from energy equation derivation.
169                       ELSE
170                          lHORs = lHORs + CALC_H(lTp,M,nn) * lRate
171                       ENDIF
172                    ENDDO
173     
174     ! Apply enthalpy transfer associated with mass transfer to get the
175     ! complete heat of reaction for Reaction H.
176                    llHORg = llHORg - RxH
177                    lHORs = lHORs + RxH
178     
179     ! Convert the heat of reaction to the appropriate units (if SI), and
180     ! store in the global array.
181                    IF(UNITS == 'SI') THEN
182                       pHoRg = pHoRg + 4.183925d3*llHORg
183                       Q_Source(NP) = Q_Source(NP) - 4.183925d3*lHORs
184                    ELSE
185                       pHORg = pHORg + llHORg
186                       Q_Source(NP) = Q_Source(NP) - lHORs
187                    ENDIF
188                 ELSE
189     ! User-defined heat of reaction.
190                    HOR_g(IJK) = HOR_g(IJK) + &
191                       DES_Reaction(H)%HoR(0) * DES_RATES(H)
192                    Q_Source(NP) = Q_Source(NP) - &
193                       DES_Reaction(H)%HoR(pM) * DES_RATES(H)
194                 ENDIF
195              ENDIF
196     
197     ! Update rate of interphase mass transfer.
198     !---------------------------------------------------------------------//
199              LM = 1 + (pM - 1)*pM/2
200              pRPhase(LM) = pRPhase(LM) + &
201                 DES_RATES(H) * DES_Reaction(H)%rPHASE(LM)
202     
203           ENDDO RXN_LP ! Loop over reactions.
204     
205     ! Calculate the toal rate of formation and consumption for each species.
206     !---------------------------------------------------------------------//
207           IF(SPECIES_EQ(0)) THEN
208              pSUMRg = SUM(pRgp(:NMAX(0)) - pRgc(:NMAX(0)))
209           ELSE
210              DO H=1, NO_OF_DES_RXNS
211                 IF(DES_Reaction(H)%nPhases <= 0) CYCLE
212                 LM = 1 + ((pM-1)*pM)/2
213                 pSUMRg = pSUMRg + &
214                    DES_RATES(H) * DES_Reaction(H)%rPHASE(LM)
215              ENDDO
216           ENDIF
217     
218           RETURN
219           END SUBROUTINE CALC_RRATES_DES
220     
221     
222     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
223     !  Subroutine name: CALC_RRATE_DES                                     !
224     !                                                                      !
225     !  Purpose: This routine manages gas-solid reactions for the continuum !
226     !  phase.                                                              !
227     !                                                                      !
228     !  Author: J.Musser                                   Date: 16-May-11  !
229     !                                                                      !
230     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
231           SUBROUTINE ZERO_RRATE_DES
232     
233           USE des_rxns
234           USE param1, only: zero
235     
236           IMPLICIT NONE
237     
238           DES_R_gp(:,:) = ZERO
239           DES_R_gc(:,:) = ZERO
240           DES_R_PHASE(:,:) = ZERO
241           DES_HOR_G(:) = ZERO
242           DES_SUM_R_g(:) = ZERO
243     
244           RETURN
245           END SUBROUTINE ZERO_RRATE_DES
246