File: RELATIVE:/../../../mfix.git/model/chem/stiff_chem_dbg_mod.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  Module name: STIFF_CHEM_DEBUG                                       !
4     !                                                                      !
5     !  Purpose:                                                            !
6     !                                                                      !
7     !  Author: J.Musser                                   Date:            !
8     !                                                                      !
9     !  Comments:                                                           !
10     !                                                                      !
11     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
12           MODULE STIFF_CHEM_DBG
13     
14           PRIVATE
15     
16     ! Variable Access:
17     !---------------------------------------------------------------------//
18           PUBLIC :: ODE_DEBUG_LEVEL
19     
20     ! Subroutine Access:
21     !---------------------------------------------------------------------//
22           PUBLIC :: ALLOCATE_STIFF_CHEM_DBG
23     
24           PUBLIC :: UPDATE_ODE_OLD
25           PUBLIC :: RESET_ODE
26     
27           PUBLIC :: CHECK_ODE_DATA
28           PUBLIC :: WRITE_ODE_LOG
29     
30     ! Routine used to compare to values.
31           LOGICAL, external :: COMPARE
32     
33     ! Static variables/parameters.
34     !---------------------------------------------------------------------//
35     ! Debug Level: (Messages)
36     ! 0 - None
37     ! 1 - Limited
38     ! 2 - Aggressive
39           INTEGER :: ODE_DEBUG_LEVEL = 1
40     
41     ! File unit for ODE Error Log.
42           INTEGER, parameter :: OEL_Unit = 6589
43     
44     ! Original field variables. Used to reset after failed integration.
45           DOUBLE PRECISION, allocatable :: ODE_VARS_o(:)
46     
47           LOGICAL :: PRINT_ERR_HEADER
48     
49           contains
50     
51     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
52     !                                                                      !
53     !  Module name: ALLOCATE_STIFF_CHEM_STATS                              !
54     !                                                                      !
55     !  Purpose:                                                            !
56     !                                                                      !
57     !  Author: J.Musser                                   Date:            !
58     !                                                                      !
59     !  Comments:                                                           !
60     !                                                                      !
61     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
62           SUBROUTINE ALLOCATE_STIFF_CHEM_DBG(lODE_DIMN)
63     
64           implicit none
65     
66     ! The number ODEs (maximum).
67           INTEGER, intent(in) :: lODE_DIMN
68     
69           allocate( ODE_VARS_o(lODE_DIMN) )
70     
71           RETURN
72           END SUBROUTINE ALLOCATE_STIFF_CHEM_DBG
73     
74     
75     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
76     !  Module name: ODE_UPDATE_OLD                                         !
77     !                                                                      !
78     !  Purpose:                                                            !
79     !                                                                      !
80     !  Author: J.Musser                                   Date:            !
81     !                                                                      !
82     !  Comments:                                                           !
83     !                                                                      !
84     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
85           SUBROUTINE UPDATE_ODE_OLD(lODE_DIMN, lODE_VARS)
86     
87           implicit none
88     
89     ! The number ODEs (maximum).
90           INTEGER, intent(in) :: lODE_DIMN
91     
92     ! MFIX field variables mapped into the ODE format.
93           DOUBLE PRECISION, intent(in) :: lODE_VARS(lODE_DIMN)
94     
95           ODE_VARS_o = lODE_VARS
96     
97           RETURN
98           END SUBROUTINE UPDATE_ODE_OLD
99     
100     
101     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
102     !  Module name: ODE_RESET                                              !
103     !                                                                      !
104     !  Purpose:                                                            !
105     !                                                                      !
106     !  Author: J.Musser                                   Date:            !
107     !                                                                      !
108     !  Comments:                                                           !
109     !                                                                      !
110     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
111           SUBROUTINE RESET_ODE(lDIMN, lVARS, lAtps)
112     
113           use compar, only: myPE
114     
115           use stiff_chem_stats, only: failedCount
116     
117           implicit none
118     
119     ! The number ODEs (maximum).
120           INTEGER, intent(in) :: lDIMN
121     
122     ! MFIX field variables mapped into the ODE format.
123           DOUBLE PRECISION, intent(inout) :: lVARS(lDIMN)
124     
125           INTEGER, intent(in) :: lAtps
126     
127           lVARS = ODE_VARS_o
128     
129           IF(lAtps == 3) failedCount(myPE) = failedCount(myPE) + 1
130     
131           RETURN
132           END SUBROUTINE RESET_ODE
133     
134     
135     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
136     !                                                                      !
137     !  Module name: CHECK_ODE_DATA                                         !
138     !                                                                      !
139     !  Purpose:                                                            !
140     !                                                                      !
141     !  Author: J.Musser                                   Date:            !
142     !                                                                      !
143     !  Comments:                                                           !
144     !                                                                      !
145     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
146           SUBROUTINE CHECK_ODE_DATA(lnD, lNEQ, loD, VARS, lUNLIMITED,     &
147              lState, lErr)
148     
149           use param1,   only : ZERO
150           use physprop, only : NMAX, MMAX
151           use toleranc, only : TMax, TMin
152     
153           implicit none
154     
155           INTEGER, intent(in) :: lnD
156     
157           INTEGER, intent(in) :: lNEQ(lnD)
158     
159     ! The number ODEs (maximum).
160           INTEGER, intent(in) :: loD
161     
162     ! MFIX field variables mapped into the ODE format.
163           DOUBLE PRECISION, intent(in) :: VARS(loD)
164     ! Flag for unlimited steps.
165           LOGICAL, intent(in) :: lUNLIMITED
166     ! Sate value returned from ODEPACK
167           INTEGER, intent(in)  :: lState
168     ! Error Flag
169           INTEGER, intent(out) :: lErr
170     
171     
172     ! Loop indicies:
173           INTEGER :: M  ! phase
174           INTEGER :: N  ! species
175           INTEGER :: Node  ! other
176     
177     
178           lErr = 0
179     
180     ! Check ODE_VAR for NaNs.
181     !---------------------------------------------------------------------//
182           DO Node=1,loD
183              IF(VARS(node).NE.VARS(node)) THEN
184                 lErr = -8
185                 return
186              ENDIF
187           ENDDO
188     
189     
190     ! Check ODEPACK State
191     !---------------------------------------------------------------------//
192           IF(lState == -1) THEN
193              lErr = lState
194              if(lUNLIMITED) return
195           ELSEIF(lState == 2) THEN
196              lErr = 0
197           ELSE
198              lErr = lState
199           ENDIF
200     
201     
202     ! Initialize.
203           Node = 1
204     
205     ! Check variables for physical values.
206     !---------------------------------------------------------------------//
207     
208     ! Gas phase density.
209           IF(VARS(Node) <= ZERO) THEN
210              lErr = -(100 + Node)
211              return
212           ENDIF
213           Node = Node + 1
214     
215     
216     ! Gas phase temperature.
217           IF(VARS(Node) > Tmax .OR. VARS(Node) < Tmin) THEN
218              lErr = -(100 + Node)
219              return
220           ENDIF
221           Node = Node + 1
222     
223     ! Gas phase species mass fractions.
224           DO N=1,NMAX(0)
225              IF(VARS(Node) > 1.009d0) THEN
226                 lErr = -(100 + Node)
227                 return
228              ENDIF
229              IF(VARS(Node) < 0.0d0) THEN
230                 IF(abs(VARS(Node)) > 1.0d-5) THEN
231                    lErr = -(100 + Node)
232                    return
233                 ENDIF
234              ENDIF
235              Node = Node + 1
236           ENDDO
237     
238     ! Solids temperature.
239           DO M = 1, MMAX
240              IF(VARS(Node) > Tmax .OR. VARS(Node) < Tmin) THEN
241                 lErr = -(100 + Node)
242                 return
243              ENDIF
244              Node = Node + 1
245           ENDDO
246     
247     
248           DO M = 1, MMAX
249              IF(lNEQ(2+M) == 1) THEN
250     
251     ! Solids volume fraction.
252                 IF(VARS(Node) <= ZERO) THEN
253                    lErr = -(100 + Node)
254                    return
255                 ENDIF
256                 Node = Node + 1
257     
258     ! Solids phase species mass fractions.
259                 DO N=1,NMAX(M)
260     
261                    IF(VARS(Node) > 1.009d0) THEN
262                       lErr = -(100 + Node)
263                       return
264                    ENDIF
265                    IF(VARS(Node) < 0.0d0) THEN
266                       IF(abs(VARS(Node)) > 1.0d-5) THEN
267                          lErr = -(100 + Node)
268                          return
269                       ENDIF
270                    ENDIF
271                    Node = Node + 1
272     
273                 ENDDO
274              ENDIF
275           ENDDO
276     
277           RETURN
278           END SUBROUTINE CHECK_ODE_DATA
279     
280     
281     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
282     !                                                                      !
283     !  Module name: ODE_ErrorLog                                           !
284     !                                                                      !
285     !  Purpose:                                                            !
286     !                                                                      !
287     !  Author: J.Musser                                   Date:            !
288     !                                                                      !
289     !  Comments:                                                           !
290     !                                                                      !
291     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
292           SUBROUTINE WRITE_ODE_LOG(lErr, lnD, lNEQ, loD, lVARS)
293     
294           use compar,   only : myPE
295           use run,      only : TIME
296           use physprop, only : MMAX
297     
298           use indices
299     
300           implicit none
301     
302           INTEGER, intent(in) :: lErr
303     
304           INTEGER, intent(in) :: lnD
305           INTEGER, intent(in) :: lNEQ(lnD)
306     
307     ! The number ODEs (maximum).
308           INTEGER, intent(in) :: loD
309           DOUBLE PRECISION, intent(in) :: lVARS(loD)
310     
311           INTEGER :: lc1
312           INTEGER :: IJK
313     
314           DOUBLE PRECISION :: ddt_lVARS(loD)
315     
316           CHARACTER(LEN=255) :: lFile
317     
318           LOGICAL :: lExist
319     
320           IJK = lNEQ(2)
321     
322           lFile = ''; write(lFile,"('StiffChem_',I2.2,'.log')") myPE
323     
324           inquire(file=lFile,exist=lExist)
325           if(lExist) then
326              open(OEL_Unit,file=trim(adjustl(lFile)), position='append',convert='big_endian')
327           else
328              open(OEL_Unit,file=trim(adjustl(lFile)), status='new',convert='big_endian')
329           endif
330     
331     
332           if(PRINT_ERR_HEADER) then
333              WRITE(OEL_Unit,9000) Time
334              PRINT_ERR_HEADER = .FALSE.
335           endif
336     
337     
338           WRITE(OEL_Unit,9001) IJK, I_OF(IJK), J_OF(IJK), myPE
339           WRITE(OEL_Unit,9002) lErr
340           WRITE(OEL_Unit,9003) lNEQ(1)
341     
342           DO lc1=1,MMAX
343              IF(lNEQ(lc1+2) == 1) THEN
344                 WRITE(OEL_Unit,9004) lc1
345              ELSEIF(lNEQ(lc1+2) == 0) THEN
346                 WRITE(OEL_Unit,9005) lc1
347              ELSE
348                 WRITE(OEL_Unit,9006) lc1, lNEQ(lc1+2)
349              ENDIF
350           ENDDO
351     
352           WRITE(OEL_Unit,"(/5x,'Field Variables:')")
353           WRITE(OEL_Unit,9007)
354           DO lc1=1, loD
355              WRITE(OEL_Unit,9008) lc1, ODE_VARS_o(lc1), lVARS(lc1), &
356                 (lVARS(lc1) - ODE_VARS_o(lc1))
357           ENDDO
358     
359           CALL STIFF_CHEM_RRATES(lNEQ, Time, lVARS, ddt_lVARS)
360     
361           WRITE(OEL_Unit,"(/5x,'Derivatives:')")
362           WRITE(OEL_Unit,9009)
363           DO lc1=1, loD
364              WRITE(OEL_Unit,9010) lc1, ddt_lVARS(lc1)
365           ENDDO
366     
367     
368           RETURN
369     
370      9000 FORMAT(//3x,'Time: ',g15.8)
371      9001 FORMAT(//5x,'Fluid Cell: ',I4,4x,'(',I3,',',I3,')',10x,&
372              'Process: ',I4)
373      9002 FORMAT(  5x,'Error Flag/ODEPACK Status: ',I4,/)
374     
375      9003 FORMAT(  5x,'Number of ODEs Solved: ', I4,/)
376      9004 FORMAT(  5x,'Solving solids phase ',I2,&
377              ' bulk density and species equations.')
378     
379      9005 FORMAT(  5x,'NOT solving solids phase ',I2,&
380              ' bulk density and species equations.')
381     
382      9006 FORMAT(  5x,'Unknown flag for solids phase ',I2,'  Flag = ',I4)
383     
384     
385      9007 FORMAT(  5x,'Var',6x,'Incoming',13x,'Returned',13x,'Difference')
386      9008 FORMAT(  5x,I3,3(3x,G18.8))
387     
388     
389      9009 FORMAT(  5x,'Var',6x,'  YDOT  ')
390      9010 FORMAT(  5x,I3,3x,G18.8)
391     
392           END SUBROUTINE WRITE_ODE_LOG
393     
394     
395           END MODULE STIFF_CHEM_DBG
396