File: /nfs/home/0/users/jenkins/mfix.git/model/chem/stiff_chem_maps_mod.f

1           MODULE STIFF_CHEM_MAPS
2     
3           PRIVATE
4     
5     ! Variable Access:
6     !---------------------------------------------------------------------//
7     
8     ! Subroutine Access:
9     !---------------------------------------------------------------------//
10           PUBLIC :: mapMFIXtoODE,   &
11                     mapODEtoMFIX
12     
13           contains
14     
15     
16     
17     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
18     !                                                                      !
19     !  Module name: mapMFIXtoODE                                           !
20     !                                                                      !
21     !  Purpose: This routine maps MFIX variables into the ODE array.       !
22     !                                                                      !
23     !  Author: J.Musser                                   Date: 07-Feb-13  !
24     !                                                                      !
25     !  Comments:                                                           !
26     !                                                                      !
27     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
28           SUBROUTINE mapMFIXtoODE(lnD, lNEQ, loD, lVars)
29     
30     ! Global Variables:
31     !---------------------------------------------------------------------//
32           use fldvar, only: ROP_g
33           use fldvar, only: ROP_s
34           use fldvar, only: T_g
35           use fldvar, only: T_s
36           use fldvar, only: X_g
37           use fldvar, only: X_s
38     
39           use physprop, only: NMAX
40           use physprop, only: MMAX
41     
42     
43           implicit none
44     
45     
46     ! Passed Variables:
47     !---------------------------------------------------------------------//
48     ! Passed array dimensions
49           INTEGER, intent(in) :: lnD  ! lNEQ
50           INTEGER, intent(in) :: loD  ! lVars
51     
52     ! (1) Number of ODEs to be solve
53     ! (2) Fluid cell index
54           INTEGER, intent(in) :: lNEQ(lnD)
55     ! Array of dependent variable initial values.
56           DOUBLE PRECISION, intent(out)  :: lVars(loD)
57     
58     
59     ! Local Variables:
60     !---------------------------------------------------------------------//
61     ! Fluid cell index
62           INTEGER :: IJK
63     ! Loop indicies:
64           INTEGER :: M  ! phase
65           INTEGER :: N  ! species
66     ! ODE Equation Counter
67           INTEGER :: Node
68     
69     
70     ! Initialize.
71           IJK = lNEQ(2)
72           lVars = 0.0d0
73           Node = 1
74     
75     ! Gas phase density.
76           lVars(Node) = ROP_G(IJK);             Node = Node + 1
77     ! Gas phase temperature.
78           lVars(Node) = T_G(IJK);               Node = Node + 1
79     ! Gas phase species mass fractions.
80           DO N=1,NMAX(0)
81              lVars(Node) = X_G(IJK,N);          Node = Node + 1
82           ENDDO
83     
84     ! Solids temperature.
85           DO M = 1, MMAX
86              lVars(Node) = T_S(IJK,M);          Node = Node + 1
87           ENDDO
88     
89           DO M = 1, MMAX
90              IF(lNEQ(2+M) == 1) THEN
91     ! Solids volume fraction.
92                 lVars(Node) = ROP_S(IJK,M);     Node = Node + 1
93     ! Solids phase species mass fractions.
94                 DO N=1,NMAX(M)
95                    lVars(Node) = X_S(IJK,M,N);  Node = Node + 1
96                 ENDDO
97              ENDIF
98           ENDDO
99     
100           RETURN
101           END SUBROUTINE mapMFIXtoODE
102     
103     
104     
105     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
106     !                                                                      !
107     !  Module name: mapODEtoMFIX                                           !
108     !                                                                      !
109     !  Purpose: This is a driver routine for mapping variables stored in   !
110     !  the ODE array back to MFIX field variables.                         !
111     !                                                                      !
112     !  Author: J.Musser                                   Date: 07-Feb-13  !
113     !                                                                      !
114     !  Comments:                                                           !
115     !                                                                      !
116     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
117           SUBROUTINE mapODEtoMFIX(lnD, lNEQ, loD, lVars)
118     
119     ! Global Variables:
120     !---------------------------------------------------------------------//
121     ! Material density.
122           use fldvar, only : RO_g, RO_s
123     ! Bulk density.
124           use fldvar, only : ROP_g, ROP_s
125     ! Gas temperature. (K)
126           use fldvar, only : T_g, T_s
127     ! Species Mass frcations.
128           use fldvar, only : X_g, X_s
129     ! Number of solids phases.
130           use physprop, only : MMAX
131     ! Number of species comprising each phase.
132           use physprop, only : NMAX
133     ! Gas constant (cal/g.K)
134           use constant, only : GAS_CONST
135     ! Gas pressure
136           use fldvar, only : P_g
137     ! Gas phase volume fraction.
138           use fldvar,   only : EP_g
139     ! Species molecular weights
140           use physprop, only : MW_g, MW_s
141     ! Gas phase mixture molecular weight.
142           use physprop, only : MW_MIX_g
143     ! Baseline/Unreaced solids density.
144           use physprop, only: BASE_ROs
145     ! Initial mass fraction of inert solids species.
146           use physprop, only: X_S0
147     ! Index of inert solids phase species.
148           use physprop, only: INERT_SPECIES
149     ! Rank ID and Rank of IO process.
150           use compar,   only : myPE, PE_IO
151     ! Runtime flag for variable solids density.
152           use run, only: SOLVE_ROs
153     
154     ! Variable solids density calculation.
155           USE eos, ONLY: EOSS
156     
157     ! Global Parameters:
158     !---------------------------------------------------------------------//
159           use param1,   only : ONE
160           use param1,   only : LARGE_NUMBER
161           use param1,   only : SMALL_NUMBER
162     
163     
164           implicit none
165     
166     
167     ! Passed Variables:
168     !---------------------------------------------------------------------//
169     ! (1) Number of ODEs to be solve
170     ! (2) Fluid cell index
171           INTEGER, intent(in) :: lnD
172           INTEGER, intent(in) :: lNEQ(lnD)
173     
174     ! Array of dependent variable initial values.
175           INTEGER, intent(in) :: loD
176           DOUBLE PRECISION, intent(in)  :: lVars(loD)
177     
178     
179     ! Local Variables:
180     !---------------------------------------------------------------------//
181     ! Fluid Cell index.
182           INTEGER :: IJK
183           INTEGER :: L
184     
185     ! Loop indicies:
186           INTEGER :: M    ! phase
187           INTEGER :: N    ! species
188           INTEGER :: Node ! ODE Equation Counter
189     
190     ! Error flags/counters.
191           INTEGER :: countNaN
192           LOGICAL :: writeMsg
193     
194           IJK = lNEQ(2)
195     
196     ! Check if any NaNs were returned from ODEPACK.  This shouldn't happend
197     ! but errors in the usr_rates file can cause this to occur.
198     !-----------------------------------------------------------------------
199           countNaN = 0
200           writeMsg = .FALSE.
201           NaN_lp: do l=1, loD
202              if(lVars(l).NE.lVars(l)) then
203                 countNaN = countNan + 1
204                 writeMsg = .TRUE.
205              endif
206           enddo NaN_lp
207     
208           if(writeMsg) then
209              write(*,"(3x,'From MapODEtoMFIX: NaNs Found! :: ',3(3x,I4))") &
210                 myPE, IJK, countNaN
211     
212              if(countNaN < loD) then
213                 do l=1, loD
214                    if(lVars(l).NE.lVars(l))                                &
215                       write(*,"(5x,' NaN in Var ',I2)") l
216                 enddo
217              endif
218           endif
219     
220     ! Directly map the ODE values into the field variable names.
221     !-----------------------------------------------------------------------
222     ! Initialize the loop counter for ODEs.
223           Node = 1
224     
225     ! Gas phase density.
226           ROP_G(IJK) = lVars(Node);                    Node = Node + 1
227     ! Gas phase temperature.
228           T_G(IJK) = lVars(Node);                      Node = Node + 1
229     ! Gas phase species mass fractions.
230           DO N=1,NMAX(0)
231              X_G(IJK,N) = lVars(Node);                 Node = Node + 1
232           ENDDO
233     
234     ! Solids temperature.
235           DO M = 1, MMAX
236              IF(ROP_s(IJK,M) > 1.0d-8) &
237                 T_S(IJK,M) = lVars(Node);              Node = Node + 1
238           ENDDO
239     
240     ! Only map back what was calculated.
241           DO M = 1, MMAX
242              IF(lNEQ(2+M) == 1) THEN
243     ! Solids bulk density.
244                 ROP_S(IJK,M) = lVars(Node);            Node = Node + 1
245     ! Solids phase species mass fractions.
246                 DO N=1,NMAX(M)
247                    X_S(IJK,M,N) = lVars(Node);         Node = Node + 1
248                 ENDDO
249     ! Update the solids density from speices composition. This update must
250     ! come after the species data is mapped back into the MFIX variables.
251                 IF(SOLVE_ROs(M)) RO_S(IJK,M) = EOSS(BASE_ROs(M), &
252                    X_s0(M,INERT_SPECIES(M)), X_s(IJK,M,INERT_SPECIES(M)))
253              ENDIF
254           ENDDO
255     
256     ! Calculate the gas volume fraction from solids volume fractions. Only
257     ! update it's value if the solids equations are being solved.
258           IF(sum(lNEQ(3:)) > 0) EP_G(IJK) = &
259              ONE - sum(ROP_S(IJK,1:MMAX)/RO_S(IJK,1:MMAX))
260     
261     ! Gas phase bulk density is updated within the stiff solver (lVar(1)).
262     ! Now that the gas phase volume fraction is updated, the gas phase
263     ! density can be backed out. RO_g * EP_g = ROP_g
264           IF(EP_g(IJK) > small_number) THEN
265              RO_g(IJK) = ROP_g(IJK) / EP_g(IJK)
266           ELSE
267     ! This case shouldn't happen, however 'LARGE_NUMBER' is used to aid
268     ! in tracking errors should this somehow become and issue.
269              RO_g(IJK) = LARGE_NUMBER
270           ENDIF
271     
272     ! Calculate the mixture molecular weight.
273           MW_MIX_G(IJK) = sum(X_G(IJK,1:NMAX(0))/MW_g(1:NMAX(0)))
274           MW_MIX_G(IJK) = ONE/MW_MIX_G(IJK)
275     
276     ! Calculate the gas phase pressure.
277           P_G(IJK) = (RO_G(IJK)*GAS_CONST*T_G(IJK))/MW_MIX_G(IJK)
278     
279           RETURN
280           END SUBROUTINE mapODEtoMFIX
281     
282           END MODULE STIFF_CHEM_MAPS
283