File: /nfs/home/0/users/jenkins/mfix.git/model/rxn_com_mod.f
1 MODULE RXN_COM
2
3 Use param
4 Use param1
5 USE compar
6 Use funits
7
8
9
10
11
12 TYPE SPECIES_
13
14
15 INTEGER pMap
16
17
18 INTEGER sMap
19
20 DOUBLE PRECISION Coeff
21
22 DOUBLE PRECISION MW
23
24 DOUBLE PRECISION xXfr
25
26 INTEGER mXfr
27
28 DOUBLE PRECISION MWxStoich
29 END TYPE SPECIES_
30
31
32 TYPE REACTION_BLOCK
33
34 CHARACTER(LEN=32) :: Name
35
36 CHARACTER(LEN=512) :: ChemEq
37
38 CHARACTER(LEN=16) :: Classification
39
40
41 LOGICAL Calc_DH
42
43 INTEGER nPhases
44
45 INTEGER nSpecies
46
47 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: HoR
48
49 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: rPhase
50
51 TYPE(SPECIES_), DIMENSION(:), ALLOCATABLE :: Species
52
53 END TYPE REACTION_BLOCK
54
55 CONTAINS
56
57
58
59
60
61
62
63
64
65
66 SUBROUTINE checkDuplicateAliases(lNg, SA_g, lMMx, lNs, SA_s)
67
68 use error_manager
69
70 IMPLICIT NONE
71
72
73 INTEGER, INTENT(IN) :: lNg
74
75 CHARACTER(len=32), DIMENSION(DIM_N_g), INTENT(IN) :: SA_g
76
77 INTEGER, INTENT(IN) :: lMMx
78
79 INTEGER, DIMENSION(DIM_M), INTENT(IN) :: lNs
80
81 CHARACTER(len=32), DIMENSION(DIM_M, DIM_N_s), INTENT(IN) :: SA_s
82
83
84 INTEGER M1, N1
85 INTEGER M2, N2
86
87 CHARACTER(len=32) SA1, SA2
88
89 CALL INIT_ERR_MSG("RXN_COM --> checkDuplicateAliases")
90
91
92 = 0
93 M2 = 0
94
95 DO N1 = 1, lNg
96 SA1 =SA_g(N1)
97 IF(len_trim(SA1) == 0) CYCLE
98 DO N2=N1+1,lNg
99 SA2 = SA_g(N2)
100 IF(len_trim(SA2) == 0) CYCLE
101 IF(compareAliases(SA1, SA2)) GoTo 100
102 ENDDO
103
104 DO M2 = 1, lMMx
105 DO N2 = 1, lNs(M2)
106 SA2 = SA_s(M2,N2)
107 IF(len_trim(SA2) == 0) CYCLE
108 IF(compareAliases(SA1, SA2)) GoTo 100
109 ENDDO
110 ENDDO
111 ENDDO
112
113 DO M1 = 1, lMMx
114 DO N1 = 1, lNs(M1)
115 SA1 = SA_s(M1,N1)
116 IF(len_trim(SA1) == 0) CYCLE
117
118 = M1
119 DO N2=N1+1, lNs(M2)
120 SA2 = SA_s(M2,N2)
121 IF(len_trim(SA2) == 0) CYCLE
122 IF(compareAliases(SA1, SA2)) GoTo 100
123 ENDDO
124
125 DO M2 = M1+1, lMMx
126 DO N2 = 1, lNs(M2)
127 SA2 = SA_s(M2,N2)
128 IF(len_trim(SA2) == 0) CYCLE
129 IF(compareAliases(SA1, SA2)) GoTo 100
130 ENDDO
131 ENDDO
132 ENDDO
133 ENDDO
134
135 CALL FINL_ERR_MSG
136 RETURN
137
138 100 WRITE(ERR_MSG, 1100) M1, N1, SA1, M2, N2, SA2
139 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
140
141 1100 FORMAT('Error 1100: Non-unique species aliases detected.',/ &
142 3x,'Phase: ',I2,', Species: ',I3,' - Alias: ',A,/ &
143 3x,'Phase: ',I2,', Species: ',I3,' - Alias: ',A,// &
144 'Please correct the mfix.dat file.')
145
146 END SUBROUTINE checkDuplicateAliases
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161 SUBROUTINE checkSpeciesInc(lNg, SA_g, lMMx, lNs, SA_s, &
162 lNRxn, lRNames, lNRxn_DES, lRNames_DES)
163
164 use error_manager
165 use toleranc
166
167 IMPLICIT NONE
168
169
170 INTEGER, INTENT(IN) :: lNg
171
172 CHARACTER(len=32), DIMENSION(DIM_N_g), INTENT(IN) :: SA_g
173
174 INTEGER, INTENT(IN) :: lMMx
175
176 INTEGER, DIMENSION(DIM_M), INTENT(IN) :: lNs
177
178 CHARACTER(len=32), DIMENSION(DIM_M, DIM_N_s), INTENT(IN) :: SA_s
179
180 INTEGER, INTENT(IN) :: lNRxn
181
182 CHARACTER(len=32), INTENT(IN) :: lRNames(DIMENSION_RXN)
183
184 INTEGER, INTENT(IN) :: lNRxn_DES
185
186 CHARACTER(len=32), INTENT(IN) :: lRNames_DES(DIMENSION_RXN)
187
188
189 INTEGER :: IOS
190
191 INTEGER, PARAMETER :: FUNIT = 167
192
193 CHARACTER(len=256) :: FILENAME
194 CHARACTER(len=128) :: INPUT
195
196 INTEGER :: SRC, M
197
198 INTEGER :: POS
199
200 INTEGER :: lIndex
201 CHARACTER(len=64) :: lName
202 CHARACTER(len=32) :: tName
203
204 INTEGER :: LINE_LEN
205
206 INTEGER, EXTERNAL :: SEEK_COMMENT
207
208 LOGICAL, EXTERNAL :: BLANK_LINE
209
210
211 INCLUDE 'mfix_directory_path.inc'
212
213 CALL INIT_ERR_MSG("RXN_COM --> checkDuplicateAliases")
214
215 SRC = 0
216
217
218 SRC_LP: DO
219 SRC = SRC + 1
220 SELECT CASE(SRC)
221
222
223 CASE(1); FILENAME = 'species.inc'
224 OPEN(UNIT=FUNIT,FILE=trim(FILENAME),STATUS='OLD',IOSTAT=IOS)
225 IF(IOS /= 0) CYCLE SRC_LP
226 WRITE(ERR_MSG, 1000)'species.inc'
227 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
228
229
230 CASE(2); FILENAME = trim(MFIX_PATH)//'/species.inc'
231 OPEN(UNIT=FUNIT,FILE=trim(FILENAME),STATUS='OLD',IOSTAT=IOS)
232 IF(IOS /= 0) CYCLE SRC_LP
233 WRITE(ERR_MSG, 1000)'mfix/model/species.inc'
234 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
235
236 1000 FORMAT(/2X,'Verifying reaction aliases in ',A)
237
238
239 CASE DEFAULT
240 WRITE(ERR_MSG, 1004)
241 CALL FLUSH_ERR_MSG
242 EXIT SRC_LP
243 END SELECT
244
245 1004 FORMAT('Warning 1004: Unable to locate the species.inc file. No ',&
246 'verification',/'of mfix.dat species aliases or reaction ', &
247 'names can be preformed.')
248
249 REWIND(FUNIT)
250 READ_LP: DO
251 READ(FUNIT,"(A)",IOSTAT=IOS) INPUT
252
253
254
255 IF(IOS > 0) THEN
256 WRITE(ERR_MSG,1200) trim(adjustl(FILENAME))
257 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
258 1200 FORMAT('Error 1200: There was a problem reading file: ',A)
259
260
261 ELSEIF(IOS<0)THEN
262 EXIT SRC_LP
263 ENDIF
264
265
266 = SEEK_COMMENT(INPUT,LEN(INPUT)) - 1
267 CALL REMOVE_COMMENT(INPUT, LINE_LEN + 1, LEN(INPUT))
268 CALL MAKE_UPPER_CASE(INPUT, LINE_LEN)
269 CALL REPLACE_TAB(INPUT, LINE_LEN)
270
271
272 IF(LINE_LEN <= 0) CYCLE READ_LP
273 IF(BLANK_LINE(INPUT)) CYCLE READ_LP
274
275 POS = INDEX(INPUT,"INTEGER, PARAMETER ::")
276 IF(POS /= 0) THEN
277 INPUT = INPUT((POS + 21):)
278 ELSE
279 CYCLE READ_LP
280 ENDIF
281
282
283 = INDEX(INPUT,"=")
284 IF(POS == 0) CYCLE READ_LP
285
286
287 WRITE(lName,"(A)") trim(adjustl(INPUT(:(POS-1))))
288
289
290 READ(INPUT((POS+1):),*,IOSTAT=IOS) lIndex
291 IF(IOS /= 0) THEN
292 WRITE(ERR_MSG,1205) 'index', trim(INPUT)
293 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
294 ENDIF
295
296 1205 FORMAT('Error 1205: Unable to obtain alias index from species.', &
297 'inc file.',//' INPUT: ',A)
298
299
300
301 IF(lIndex <= lNg) THEN
302 tName = SA_g(lIndex)
303 IF(compareAliases(tName, lName)) CYCLE READ_LP
304 ENDIF
305
306
307 DO M = 1, lMMx
308 IF(lIndex <= lNs(M))THEN
309 tName = SA_s(M, lIndex)
310 IF(compareAliases(tName, lName)) CYCLE READ_LP
311 ENDIF
312 ENDDO
313
314
315 IF(lIndex <= lNRxn)THEN
316 tName = lRNames(lIndex)
317 IF(compareAliases(tName, lName)) CYCLE READ_LP
318 ENDIF
319
320
321 IF(lIndex <= lNRxn_DES)THEN
322 tName = lRNames_DES(lIndex)
323 IF(compareAliases(tName, lName)) CYCLE READ_LP
324 ENDIF
325
326 WRITE(ERR_MSG,1300) trim(lName), lIndex
327 CALL FLUSH_ERR_MSG
328
329 1300 FORMAT('Error 1300: An entry in the species.inc file does not ', &
330 'match any inputs',/'in the mfix.dat file.'/3x,'Name: ',A,4x, &
331 'Index: ',I3,/'If the quantity or order of gas species, ', &
332 'solids species, or chemical',/'reactions has changed, then ',&
333 'the executable must be re-build. Please',/'see the document',&
334 'ation for specifying chemical reactions.')
335
336 ENDDO READ_LP
337 ENDDO SRC_LP
338
339 CLOSE(FUNIT)
340 CALL FINL_ERR_MSG
341 RETURN
342
343 END SUBROUTINE checkSpeciesInc
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358 LOGICAL FUNCTION compareAliases(lS1, lS2, N1, N2)
359
360 IMPLICIT NONE
361
362 CHARACTER(len=*), INTENT(IN) :: lS1, lS2
363
364 INTEGER, OPTIONAL, INTENT(IN) :: N1, N2
365
366 CALL MAKE_UPPER_CASE (lS1, len(lS1))
367 CALL MAKE_UPPER_CASE (lS2, len(lS2))
368
369 compareAliases = .FALSE.
370 IF(trim(lS1) == trim(lS2)) compareAliases = .TRUE.
371
372 IF(.NOT.compareAliases) RETURN
373
374 IF(PRESENT(N1) .AND. PRESENT(N2)) THEN
375 IF(N1 == N2) THEN
376 compareAliases = .TRUE.
377 ELSE
378 compareAliases = .FALSE.
379 ENDIF
380 ENDIF
381
382 RETURN
383 END FUNCTION compareAliases
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399 SUBROUTINE WRITE_RXN_SUMMARY(RxN, lSAg, lSAs, ABORT, fUNIT)
400
401 USE toleranc
402
403 IMPLICIT NONE
404
405
406 TYPE(REACTION_BLOCK), POINTER, INTENT(IN) :: RxN
407
408
409 CHARACTER(len=32), DIMENSION(DIM_N_g), INTENT(IN) :: lSAg
410
411 CHARACTER(len=32), DIMENSION(DIM_M, DIM_N_s), INTENT(IN) :: lSAs
412
413 LOGICAL, INTENT(IN) :: ABORT
414
415
416 INTEGER, OPTIONAL :: fUNIT
417
418 CHARACTER(LEN=72) :: OUTPUT, full, divided, empty
419
420 CHARACTER(LEN=32) :: lSP
421
422 INTEGER lN, M, N
423 INTEGER lS, lE
424
425 INTEGER UNIT_FLAG
426
427 IF(present(fUnit)) THEN
428 UNIT_FLAG = fUNIT
429 ELSE
430 UNIT_FLAG = -1
431 ENDIF
432
433 empty = ' '
434 CALL WRITE_RS0(empty, UNIT_FLAG)
435
436 full = ''
437 WRITE(full,2000)
438
439 divided = ''
440 WRITE(divided,2005)
441
442
443 CALL WRITE_RS0(full, UNIT_FLAG)
444
445 = ''
446 WRITE(OUTPUT, 2001)trim(RxN%Name)
447 OUTPUT(72:72) = '|'
448 CALL WRITE_RS0(OUTPUT, UNIT_FLAG)
449
450
451 CALL WRITE_RS0(full, UNIT_FLAG)
452
453 OUTPUT = ''
454 WRITE(OUTPUT, 2002)trim(RxN%ChemEq(1:54))
455 OUTPUT(72:72) = '|'
456 CALL WRITE_RS0(OUTPUT, UNIT_FLAG)
457
458 CALL WRITE_RS0(full, UNIT_FLAG)
459
460 IF(RxN%nSpecies > 0) THEN
461
462 OUTPUT = ''
463 WRITE(OUTPUT, 2007)trim(RxN%Classification)
464 OUTPUT(72:72) = '|'
465 CALL WRITE_RS0(OUTPUT, UNIT_FLAG)
466
467 CALL WRITE_RS0(full, UNIT_FLAG)
468
469 WRITE(OUTPUT,2003); CALL WRITE_RS0(OUTPUT, UNIT_FLAG)
470 WRITE(OUTPUT,2004); CALL WRITE_RS0(OUTPUT, UNIT_FLAG)
471 CALL WRITE_RS0(divided, UNIT_FLAG)
472 ENDIF
473
474
475 DO lN = 1, RxN%nSpecies
476
477 M = RxN%Species(lN)%pMap
478 N = RxN%Species(lN)%sMap
479
480 WRITE(OUTPUT,2006)
481
482 IF(M == 0) THEN
483 IF(len_trim(lSAg(N)) > 8) THEN
484 lSP = lSAg(N)
485 OUTPUT(5:13) = lSP(1:8)
486 ELSE
487 lS = (9-int(len_trim(lSAg(N))/2))
488 lE = lS + len_trim(lSAg(N))
489 OUTPUT(lS:lE) = trim(lSAg(N))
490 ENDIF
491 WRITE(OUTPUT(32:35),"(A)") 'Gas'
492 ELSE
493 IF(len_trim(lSAs(M,N)) > 8) THEN
494 lSP = lSAs(M,N)
495 OUTPUT(5:13) = lSP(1:8)
496 ELSE
497 lS = (9-int(len_trim(lSAs(M,N))/2))
498 lE = lS + len_trim(lSAs(M,N))
499 OUTPUT(lS:lE) = trim(lSAs(M,N))
500 ENDIF
501 WRITE(OUTPUT(30:36),"(A,I2)") 'Solid',M
502 ENDIF
503 WRITE(OUTPUT(43:44),"(I2)") N
504 WRITE(OUTPUT(51:60),"(F9.4)") RxN%Species(lN)%MW
505
506 IF(COMPARE(RxN%Species(lN)%Coeff, ZERO)) THEN
507 WRITE(OUTPUT(17:26),"(F9.4)") ZERO
508 WRITE(OUTPUT(63:71),"(A)") 'Catalyst'
509 ELSEIF(RxN%Species(lN)%Coeff < ZERO) THEN
510 WRITE(OUTPUT(17:26),"(F9.4)") -RxN%Species(lN)%Coeff
511 WRITE(OUTPUT(63:71),"(A)") 'Reactant'
512 ELSE
513 WRITE(OUTPUT(17:26),"(F9.4)") RxN%Species(lN)%Coeff
514 WRITE(OUTPUT(63:70),"(A)") 'Product'
515 ENDIF
516 CALL WRITE_RS0(OUTPUT, UNIT_FLAG)
517 CALL WRITE_RS0(divided, UNIT_FLAG)
518
519 ENDDO
520
521 CALL WRITE_RS0(empty, UNIT_FLAG)
522
523 IF(ABORT) CALL MFIX_EXIT(myPE)
524 RETURN
525
526
527 2000 FORMAT(2X,'|',68('-'),'|')
528
529 2001 FORMAT(2X,'| Name: ',A)
530 2002 FORMAT(2x,'| Chemical Eq: ',A)
531
532 2003 FORMAT(' | Species | Stoich | | Species |', &
533 ' Molecular | |')
534
535 2004 FORMAT(' | Alias | Coeff. | Phase | Index |', &
536 ' Weight | Type |')
537
538
539 2005 FORMAT(2X,'|',10('-'),'|',13('-'),'|',9('-'),'|',9('-'),'|', &
540 12('-'),'|',10('-'),'|')
541
542 2006 FORMAT(2X,'|',10(' '),'|',13(' '),'|',9(' '),'|',9(' '),'|', &
543 12(' '),'|',10(' '),'|')
544
545
546 2007 FORMAT(2X,'| Classification: ',A)
547
548 contains
549
550
551
552
553
554
555
556
557
558
559
560
561
562 SUBROUTINE WRITE_RS0(LINE, UFLAG)
563
564 use error_manager
565
566 IMPLICIT NONE
567
568 CHARACTER(len=*), INTENT(IN) :: LINE
569 INTEGER, INTENT(IN) :: UFLAG
570
571 CALL INIT_ERR_MSG("WRITE_RXN_SUMMARY --> WRITE_RS0")
572
573 IF(UFLAG == -1)THEN
574 WRITE(ERR_MSG,*) LINE
575 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
576 ELSE
577 WRITE(UFLAG,*) LINE
578 ENDIF
579 CALL FINL_ERR_MSG
580
581 RETURN
582 END SUBROUTINE WRITE_RS0
583 END SUBROUTINE WRITE_RXN_SUMMARY
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598 SUBROUTINE checkThermoReqs(RxN, S_g, S_s, rDB, MWg, MWs, Cpg0, Cps0)
599
600 use error_manager
601 use toleranc
602
603 IMPLICIT NONE
604
605
606 TYPE(REACTION_BLOCK), POINTER, INTENT(INOUT) :: RxN
607
608 CHARACTER(len=18), INTENT(IN) :: S_g(DIM_N_g)
609 CHARACTER(len=18), INTENT(in) :: S_s(DIM_M, DIM_N_s)
610 LOGICAL, INTENT(inout) :: rDB(0:DIM_M, DIM_N_g)
611 DOUBLE PRECISION, INTENT(in) :: Cpg0
612 DOUBLE PRECISION, INTENT(in) :: Cps0(DIM_M)
613 DOUBLE PRECISION, INTENT(inout) :: MWg(DIM_N_g)
614 DOUBLE PRECISION, INTENT(inout) :: MWs(DIM_M, DIM_N_s)
615
616 LOGICAL :: CP_FATAL
617 LOGICAL :: CHECK_DATABASE
618
619 INTEGER :: M, N, lN
620
621
622 CALL INIT_ERR_MSG("RXN_COM --> checkThermoReqs")
623
624 CHECK_DATABASE = .FALSE.
625 CP_FATAL = .FALSE.
626
627
628
629 DO lN = 1, RxN%nSpecies
630 M = RxN%Species(lN)%pMap
631 N = RxN%Species(lN)%sMap
632 IF(M == 0) THEN
633 IF(Cpg0 /= UNDEFINED) THEN
634 CP_FATAL = .TRUE.
635 ELSEIF((RxN%Calc_DH .AND. .NOT.rDB(M,N)) .OR. &
636 (MWg(N) == UNDEFINED)) THEN
637 CHECK_DATABASE = .TRUE.
638 ENDIF
639 ELSE
640 IF(Cps0(M) /= UNDEFINED) THEN
641 CP_FATAL = .TRUE.
642 ELSEIF((RxN%Calc_DH .AND. .NOT.rDB(M,N)) .OR. &
643 (MWs(M,N) == UNDEFINED)) THEN
644 CHECK_DATABASE = .TRUE.
645 ENDIF
646 ENDIF
647 ENDDO
648
649
650 IF(CP_FATAL) THEN
651
652 WRITE(ERR_MSG, 1100) trim(RxN%Name)
653 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
654
655 1100 FORMAT('Error 1100: One or more phases associated with ', &
656 'reaction ',A,/'has specified constant specific heat (C_PG0/',&
657 'Cps0). This is',/'not permitted for reacting phases. ', &
658 'Please correct the mfix.dat file.')
659
660 ELSEIF(CHECK_DATABASE) THEN
661
662 WRITE(ERR_MSG, 1101) trim(RxN%Name)
663 CALL FLUSH_ERR_MSG
664
665 1101 FORMAT('Message 1101: One or more molecular weights and/or ', &
666 'thermochemical data',/'is needed for reaction ',A,'. The ', &
667 'thermochemical database',/'will be used to gather the ', &
668 'necessary data.')
669
670 ENDIF
671
672 IF(CHECK_DATABASE) THEN
673 WRITE(ERR_MSG, 1200)
674 CALL FLUSH_ERR_MSG(FOOTER=.FALSE.)
675 ENDIF
676
677 1200 FORMAT('Message 1200: Searching thermochemical databases for ',&
678 'species data.',/' ')
679
680 DO lN = 1, RxN%nSpecies
681 M = RxN%Species(lN)%pMap
682 N = RxN%Species(lN)%sMap
683 IF(M == 0) THEN
684 IF((RxN%Calc_DH .AND. .NOT.rDB(M,N)) .OR. &
685 (MWg(N) == UNDEFINED)) THEN
686
687
688 IF(S_g(N) == UNDEFINED_C) THEN
689 WRITE(ERR_MSG,1000) trim(iVar('SPECIES_g',N))
690 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
691 ENDIF
692
693
694 WRITE(ERR_MSG, 3001) N, trim(S_g(N))
695 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
696
697 CALL READ_DATABASE(0, N, S_g(N), MWg(N))
698
699 (0,N) = .TRUE.
700 ENDIF
701 RxN%Species(lN)%MW = MWg(N)
702 ELSE
703 IF((RxN%Calc_DH .AND. .NOT.rDB(M,N)) .OR. &
704 (MWs(M,N) == UNDEFINED)) THEN
705
706
707 IF(S_s(M,N) == UNDEFINED_C) THEN
708 WRITE(ERR_MSG,1000) trim(iVar('SPECIES_s',M,N))
709 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
710 ENDIF
711
712 WRITE(ERR_MSG, 3001) N, trim(S_s(M,N))
713 CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
714 CALL READ_DATABASE(M,N,S_s(M,N),MWs(M,N))
715
716 (M,N) = .TRUE.
717 ENDIF
718 RxN%Species(lN)%MW = MWs(M,N)
719 ENDIF
720 ENDDO
721
722 IF(CHECK_DATABASE) CALL FLUSH_ERR_MSG(HEADER=.FALSE.)
723
724 3001 FORMAT(/2x,'>',I3,': Species: ',A)
725
726 CALL FINL_ERR_MSG
727
728 RETURN
729
730 1000 FORMAT('Error 1000: Required input not specified: ',A,/'Please ',&
731 'correct the mfix.dat file.')
732
733 END SUBROUTINE checkThermoReqs
734
735
736
737
738
739
740
741
742
743
744
745
746
747 SUBROUTINE checkMassBalance(CALLER, RxN, lnMT, IER)
748
749 USE toleranc
750
751 IMPLICIT NONE
752
753 CHARACTER(len=*), INTENT(IN) :: CALLER
754
755
756 TYPE(REACTION_BLOCK), POINTER, INTENT(INOUT) :: RxN
757
758 DOUBLE PRECISION, INTENT(OUT) :: lnMT(0:DIM_M)
759 INTEGER, INTENT(OUT) :: IER
760
761 INTEGER M, N, lN
762 DOUBLE PRECISION rSUM, pSUM
763 DOUBLE PRECISION MWxStoich
764
765 INTEGER sprCount, sprIndex
766
767 DOUBLE PRECISION, PARAMETER :: massBalanceTol = 1.0d-3
768
769
770
771
772 = 0
773 rSUM = ZERO
774 pSUM = ZERO
775 lnMT(:) = ZERO
776 sprCount = 0
777
778
779
780 DO lN = 1, RxN%nSpecies
781 M = RxN%Species(lN)%pMap
782 N = RxN%Species(lN)%sMap
783
784
785 = RxN%Species(lN)%MW * RxN%Species(lN)%Coeff
786 RxN%Species(lN)%MWxStoich = MWxStoich
787
788
789
790
791 (M) = lnMT(M) + MWxStoich
792
793 IF(MWxStoich < ZERO) THEN
794 rSUM = rSUM - MWxStoich
795 IF(M /= 0) THEN
796 sprCount = sprCount + 1
797 IF(sprCount == 1) THEN
798 sprIndex = M
799
800 ELSEIF( M /= sprIndex) THEN
801 IF(DMP_LOG) THEN
802 WRITE(*,1002) trim(CALLER), trim(RxN%Name)
803 WRITE(UNIT_LOG,1002) trim(CALLER), trim(RxN%Name)
804 IER = 1
805 ENDIF
806 ENDIF
807 ENDIF
808 ELSE
809 pSUM = pSUM + MWxStoich
810 ENDIF
811 ENDDO
812
813 IF (.NOT.COMPARE(rSUM,pSUM)) THEN
814 IF(DMP_LOG) THEN
815 WRITE(*,1001) trim(CALLER), trim(RxN%Name), rSUM, pSUM
816 WRITE(UNIT_LOG,1001) trim(CALLER), trim(RxN%Name), rSUM,pSUM
817 IER = 1
818 ENDIF
819 ENDIF
820
821 RETURN
822
823
824
825
826 FORMAT(/1X,70('*')/' From: ',A,' --> RXN_COM -->', &
827 ' checkMassBalance',/' Error 1001: Stoichiometry is not', &
828 ' consistent with molecular weights',/' for reaction ',A,'.',/&
829 ' Mass of reactants: ',F12.4,/' Mass of products: ',F12.4,/ &
830 1X,70('*')/)
831
832 1002 FORMAT(/1X,70('*')/' From: ',A,' --> RXN_COM -->', &
833 ' checkMassBalance',/' Error 1002: More than one solids', &
834 ' phase fules was detected. Unable to',/' determine solids/', &
835 'solids heat of reaction unambiguously for',/' reaction ',A, &
836 '.',/1X,70('*')/)
837
838 END SUBROUTINE checkMassBalance
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854 SUBROUTINE calcInterphaseTxfr(CALLER, RxN, lnMT, lEEq, lSEq, &
855 lSAg, lMMx, lSAs)
856
857 USE toleranc
858
859 IMPLICIT NONE
860
861 CHARACTER(len=*), INTENT(IN) :: CALLER
862
863
864 TYPE(REACTION_BLOCK), POINTER, INTENT(INOUT) :: RxN
865
866 DOUBLE PRECISION, INTENT(IN) :: lnMT(0:DIM_M)
867
868 LOGICAL, INTENT(IN) :: lEEq
869
870 LOGICAL, INTENT(IN) :: lSEq(0:DIM_M)
871
872 CHARACTER(len=32), DIMENSION(DIM_N_g), INTENT(IN) :: lSAg
873
874 INTEGER, INTENT(IN) :: lMMx
875
876 CHARACTER(len=32), DIMENSION(DIM_M, DIM_N_s), INTENT(IN) :: lSAs
877
878 INTEGER toPhase, toPhaseCount, mCount
879 INTEGER fromPhase, fromPhaseCount
880 INTEGER catPhase
881
882 INTEGER M, MM, LL
883 INTEGER lM, lN
884
885 DOUBLE PRECISION, PARAMETER :: massBalanceTol = 1.0d-3
886
887
888 IF(Allocated(RxN%rPhase)) RxN%rPhase(:) = ZERO
889
890
891
892 IF(RxN%nPhases == 1) THEN
893
894
895
896 %rPhase(:) = ZERO
897
898 DO lN = 1, RxN%nSpecies
899 M = RxN%Species(lN)%pMap
900 RxN%Species(lN)%mXfr = M
901 ENDDO
902 RxN%Classification = "Homogeneous"
903
904 ELSE
905
906 = 0
907 fromPhaseCount = 0
908 DO M = 0, lMMx
909
910
911 IF (lnMT(M) > massBalanceTol) THEN
912 toPhaseCount = toPhaseCount + 1
913 toPhase = M
914
915
916 ELSEIF(lnMT(M) < -massBalanceTol) THEN
917 fromPhaseCount = fromPhaseCount + 1
918 fromPhase = M
919 ENDIF
920 ENDDO
921
922
923 IF(toPhaseCount == 1) THEN
924
925 %Classification = "Heterogeneous"
926 DO M = 0, lMMx
927 IF(M /= toPhase) THEN
928 IF (toPhase < M) THEN
929 LM = 1 + toPhase + ((M-1)*M)/2
930 RxN%rPhase(LM) = -lnMT(M)
931 ELSE
932 LM = 1 + M + ((toPhase-1)*toPhase)/2
933 RxN%rPhase(LM) = lnMT(M)
934 ENDIF
935
936
937
938
939 IF(abs(RxN%rPhase(LM)) > SMALL_NUMBER) THEN
940 IF((lSEq(toPhase) .AND. .NOT.lSEq(M)) .OR. &
941 (.NOT.lSEq(toPhase) .AND. lSEq(M))) THEN
942 IF(DMP_LOG) THEN
943 WRITE(*,1001) trim(CALLER)
944 WRITE(UNIT_LOG,1001) trim(CALLER)
945 IF(lSEq(M)) THEN
946 WRITE(*,1101) M, 'Solving'
947 WRITE(*,1101) toPhase, 'Not Solving'
948 WRITE(UNIT_LOG,1101) M, 'Solving'
949 WRITE(UNIT_LOG,1101) toPhase, &
950 'Not Solving'
951 ELSE
952 WRITE(*,1101) toPhase, 'Solving'
953 WRITE(*,1101) M, 'Not Solving'
954 WRITE(UNIT_LOG,1101) M, 'Solving'
955 WRITE(UNIT_LOG,1101) toPhase, &
956 'Not Solving'
957 ENDIF
958 WRITE(*,1000)
959 WRITE(UNIT_LOG,1000)
960 ENDIF
961 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), &
962 lSAs(:,:), .TRUE.)
963 ENDIF
964 ENDIF
965 ENDIF
966 ENDDO
967
968
969 IF(lEEq .AND. RxN%Calc_DH) THEN
970 DO lN = 1, RxN%nSpecies
971 M = RxN%Species(lN)%pMap
972
973 IF(M == 0) THEN
974
975 IF(toPhase == 0) THEN
976
977
978 = 0
979
980 DO MM = 1, lMMx
981 LM = 1 + (MM-1)*MM/2
982
983 IF(RxN%rPhase(LM) > 0) THEN
984
985 %Species(lN)%mXfr = MM
986
987
988 %Species(lN)%xXfr = ZERO
989
990 = mCount + 1
991 ENDIF
992 ENDDO
993 IF(mCount /= 1) THEN
994 IF(DMP_LOG) THEN
995 WRITE(*,1002) trim(CALLER), &
996 trim(RxN%ChemEq)
997 WRITE(*,1000)
998 WRITE(UNIT_LOG,1002) trim(CALLER), &
999 trim(RxN%ChemEq)
1000 WRITE(UNIT_LOG,1000)
1001 ENDIF
1002 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), &
1003 lSAs(:,:), .TRUE.)
1004 ENDIF
1005
1006
1007 ELSE
1008
1009
1010
1011 %Species(lN)%mXfr = toPhase
1012
1013 %Species(lN)%xXfr = ZERO
1014 ENDIF
1015
1016
1017 ELSE
1018
1019 DO LL = 1, lMMx-1
1020 DO MM = LL + 1, lMMx
1021 IF(M /= LL .AND. M /= MM) CYCLE
1022 LM = LL + 1 + (MM-1)*MM/2
1023 IF(RxN%rPhase(LM) == ZERO) CYCLE
1024
1025
1026 IF( M == LL .AND. &
1027 RxN%rPhase(LM) < ZERO) THEN
1028
1029 %Species(lN)%mXfr = MM
1030
1031
1032 %Species(lN)%xXfr = &
1033 abs(lnMT(MM) / lnMT(LL))
1034
1035
1036 ELSEIF( M == MM .AND. &
1037 RxN%rPhase(LM) > ZERO) THEN
1038
1039 %Species(lN)%mXfr = LL
1040
1041
1042 %Species(lN)%xXfr = &
1043 abs(lnMT(LL) / lnMT(MM))
1044 ENDIF
1045 ENDDO
1046 ENDDO
1047 ENDIF
1048 ENDDO
1049 ENDIF
1050
1051
1052 ELSEIF(fromPhaseCount == 1) THEN
1053 RxN%Classification = "Heterogeneous"
1054 DO M = 0, lMMx
1055 IF (M /= fromPhase) THEN
1056 IF(M < fromPhase) THEN
1057 LM = 1 + M + ((fromPhase-1)*fromPhase)/2
1058 RxN%rPhase(LM) = lnMT(M)
1059 ELSE
1060 LM = 1 + fromPhase + ((M-1)*M)/2
1061 RxN%rPhase(LM) = -lnMT(M)
1062 ENDIF
1063
1064
1065
1066 IF(abs(RxN%rPhase(LM)) > SMALL_NUMBER) THEN
1067 IF((lSEq(fromPhase) .AND. .NOT.lSEq(M)) .OR. &
1068 (.NOT.lSEq(fromPhase) .AND. lSEq(M))) THEN
1069 IF(DMP_LOG) THEN
1070 WRITE(*,1001) trim(CALLER)
1071 WRITE(UNIT_LOG,1001) trim(CALLER)
1072 IF(lSEq(M)) THEN
1073 WRITE(*,1101) M, 'Solving'
1074 WRITE(*,1101) fromPhase, 'Not Solving'
1075 WRITE(UNIT_LOG,1101) M, 'Solving'
1076 WRITE(UNIT_LOG,1101) fromPhase, &
1077 'Not Solving'
1078 ELSE
1079 WRITE(*,1101) toPhase, 'Solving'
1080 WRITE(*,1101) M, 'Not Solving'
1081 WRITE(UNIT_LOG,1101) fromPhase, 'Solving'
1082 WRITE(UNIT_LOG,1101) M, 'Not Solving'
1083 ENDIF
1084 WRITE(*,1000)
1085 WRITE(UNIT_LOG,1000)
1086 ENDIF
1087 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), &
1088 lSAs(:,:), .TRUE.)
1089 ENDIF
1090 ENDIF
1091 ENDIF
1092 END DO
1093
1094
1095 IF(lEEq .AND. RxN%Calc_DH) THEN
1096 DO lN = 1, RxN%nSpecies
1097 M = RxN%Species(lN)%pMap
1098
1099
1100 IF(M == 0) THEN
1101
1102 IF(fromPhase == 0) THEN
1103
1104
1105 = 0
1106
1107 DO MM = 1, lMMx
1108 LM = 1 + (MM-1)*MM/2
1109
1110 IF(RxN%rPhase(LM) < 0) THEN
1111
1112 %Species(lN)%mXfr = MM
1113
1114
1115 %Species(lN)%xXfr = ZERO
1116
1117 = mCount + 1
1118 ENDIF
1119 ENDDO
1120 IF(mCount /=1 ) THEN
1121 IF(DMP_LOG) THEN
1122 WRITE(*,1002) trim(CALLER), &
1123 trim(RxN%ChemEq)
1124 WRITE(*,1000)
1125 WRITE(UNIT_LOG,1002) trim(CALLER), &
1126 trim(RxN%ChemEq)
1127 WRITE(UNIT_LOG,1000)
1128 ENDIF
1129 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), &
1130 lSAs(:,:), .TRUE.)
1131 ENDIF
1132 ELSE
1133
1134
1135 %Species(lN)%mXfr = fromPhase
1136
1137
1138 %Species(lN)%xXfr = ZERO
1139 ENDIF
1140
1141
1142 ELSE
1143
1144 DO LL = 1, lMMx-1
1145 DO MM = LL + 1, lMMx
1146 IF(M /= LL .AND. M /= MM) CYCLE
1147 LM = LL + 1 + (MM-1)*MM/2
1148 IF(RxN%rPhase(LM) == ZERO) CYCLE
1149
1150
1151 IF( M == LL .AND. &
1152 RxN%rPhase(LM) < ZERO) THEN
1153
1154 %Species(lN)%mXfr = MM
1155
1156
1157 %Species(lN)%xXfr = &
1158 abs(lnMT(MM) / lnMT(LL))
1159
1160
1161 ELSEIF( M == MM .AND. &
1162 RxN%rPhase(LM) > ZERO) THEN
1163
1164 %Species(lN)%mXfr = LL
1165
1166
1167 %Species(lN)%xXfr = &
1168 abs(lnMT(LL) / lnMT(MM))
1169 ENDIF
1170 ENDDO
1171 ENDDO
1172 ENDIF
1173 ENDDO
1174 ENDIF
1175
1176
1177
1178 ELSEIF(toPhaseCount == 0 .AND. fromPhaseCount == 0) THEN
1179
1180
1181 IF(RxN%nPhases > 0) RxN%Classification = "Catalytic"
1182 RxN%rPhase(:) = ZERO
1183
1184 IF(lEEq .AND. RxN%Calc_DH) THEN
1185
1186
1187 = -1
1188 DO lN= 1, RxN%nSpecies
1189 IF(COMPARE(RxN%Species(lN)%Coeff,ZERO)) THEN
1190 IF(catPhase /= -1) THEN
1191 IF(catPhase /= RxN%Species(lN)%pMap) THEN
1192 IF(DMP_LOG) THEN
1193 WRITE(*,1002) trim(CALLER), &
1194 trim(RxN%Name)
1195 WRITE(*,1000)
1196 WRITE(UNIT_LOG,1002) trim(CALLER), &
1197 trim(RxN%Name)
1198 WRITE(UNIT_LOG,1000)
1199 ENDIF
1200 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), &
1201 lSAs(:,:), .TRUE.)
1202 ENDIF
1203 ELSE
1204 catPhase = RxN%Species(lN)%pMap
1205 ENDIF
1206 ENDIF
1207 ENDDO
1208
1209 IF(catPhase == -1) THEN
1210 IF(DMP_LOG) THEN
1211 WRITE(*,1003) trim(CALLER), 'catalyst', &
1212 trim(RxN%Name)
1213 WRITE(*,1000)
1214 WRITE(UNIT_LOG,1003) trim(CALLER), &
1215 'catalyst', trim(RxN%Name)
1216 WRITE(UNIT_LOG,1000)
1217 ENDIF
1218 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), &
1219 lSAs(:,:), .TRUE.)
1220 ENDIF
1221
1222
1223 = -1
1224 DO lN = 1, RxN%nSpecies
1225 IF(.NOT.COMPARE(RxN%Species(lN)%Coeff,ZERO)) THEN
1226 IF(toPhase /= -1) THEN
1227 IF(toPhase /= RxN%Species(lN)%pMap) THEN
1228 IF(DMP_LOG) THEN
1229 WRITE(*,1002) trim(CALLER), &
1230 trim(RxN%Name)
1231 WRITE(*,1000)
1232 WRITE(UNIT_LOG,1002) trim(CALLER), &
1233 trim(RxN%Name)
1234 WRITE(UNIT_LOG,1000)
1235 ENDIF
1236 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), &
1237 lSAs(:,:), .TRUE.)
1238 ENDIF
1239 ELSE
1240 toPhase = RxN%Species(lN)%pMap
1241 ENDIF
1242 ENDIF
1243 ENDDO
1244
1245 IF(toPhase == -1) THEN
1246 IF(DMP_LOG) THEN
1247 WRITE(*,1003) trim(CALLER), 'reacting', &
1248 trim(RxN%Name)
1249 WRITE(*,1000)
1250 WRITE(UNIT_LOG,1003) trim(CALLER), 'reacting', &
1251 trim(RxN%Name)
1252 WRITE(UNIT_LOG,1000)
1253 ENDIF
1254 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), lSAs(:,:),.TRUE.)
1255 ENDIF
1256
1257
1258 IF(catPhase == toPhase) THEN
1259 IF(DMP_LOG) THEN
1260 WRITE(*,1004) trim(CALLER), trim(RxN%Name)
1261 WRITE(*,1000)
1262 WRITE(UNIT_LOG,1004) trim(CALLER),trim(RxN%Name)
1263 WRITE(UNIT_LOG,1000)
1264 ENDIF
1265 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), lSAs(:,:),.TRUE.)
1266
1267 ELSEIF(toPhase == 0) THEN
1268 DO lN = 1, RxN%nSpecies
1269 IF(RxN%Species(lN)%pMap == 0) THEN
1270
1271 %Species(lN)%mXfr = catPhase
1272
1273
1274 %Species(lN)%xXfr = ZERO
1275 ENDIF
1276 ENDDO
1277 ELSEIF(catPhase == 0) THEN
1278 DO lN = 1, RxN%nSpecies
1279 IF(RxN%Species(lN)%pMap == 0) THEN
1280
1281 %Species(lN)%mXfr = toPhase
1282
1283
1284 %Species(lN)%xXfr = ZERO
1285 ENDIF
1286 ENDDO
1287 ENDIF
1288 ENDIF
1289 ELSE
1290
1291
1292
1293 CALL WRITE_RXN_SUMMARY(RxN, lSAg(:), lSAs(:,:),.FALSE.)
1294 WRITE(*,1002) trim(CALLER), trim(RxN%ChemEq)
1295 WRITE(*,1000)
1296 WRITE(UNIT_LOG,1002) trim(CALLER), trim(RxN%ChemEq)
1297 WRITE(UNIT_LOG,1000)
1298 CALL MFiX_EXIT(myPE)
1299 ENDIF
1300 ENDIF
1301
1302 RETURN
1303
1304
1305
1306
1307 FORMAT(/' Please refer to the Readme file on the required input',&
1308 ' format and make',/' the necessary corrections to the data', &
1309 ' file.',/1X,70('*')//)
1310
1311 1001 FORMAT(//1X,70('*')/' From: ',A,' --> RXN_COM -->', &
1312 ' calcInterphaseTxfr',/' Error 1001: A chemical reaction or', &
1313 ' phase change was detected between',/' a phases solving', &
1314 ' species equations and another phase not solving',/ &
1315 ' species equations.',/)
1316
1317 1101 FORMAT(' Phase ',I2,': ',A,' species equations.')
1318
1319 1002 FORMAT(//1X,70('*')/' From: ',A,' --> RXN_COM -->', &
1320 ' calcInterphaseTxfr',/' Error 1002: Reaction complexity', &
1321 ' exceeds implementation capabilities.',/' Unable to', &
1322 ' determine unambiguously interphase heat or mass transfer.', &
1323 //' Reaction: ',A,//' Consider splitting the chemical', &
1324 ' reaction equation into two or more',/' separate equations.',&
1325 ' The same reaction rate calculated in usr_rates',/' can be', &
1326 ' used for the multiple reactions to ensure mass')
1327
1328 1003 FORMAT(//1X,70('*')/' From: ',A,' --> RXN_COM -->', &
1329 ' calcInterphaseTxfr',/' Error 1003: Unable to determine ',A, &
1330 ' phase for catalytic reaction'/1X,A,'.')
1331
1332 1004 FORMAT(//1X,70('*')/' From: ',A,' --> RXN_COM -->', &
1333 ' calcInterphaseTxfr',/' Error 1004: Unable to distinguish', &
1334 ' catalyst phase from reacting phase',/' for catalytic', &
1335 ' reaction ',A,'.')
1336
1337 END SUBROUTINE calcInterphaseTxfr
1338
1339
1340 END MODULE RXN_COM
1341