File: RELATIVE:/../../../mfix.git/model/check_data/check_odepack_stiff_chem.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  Module name: CHECK_DATA_CHEM                                        !
4     !  Author: J.Musser                                   Date: 02-Aug-13  !
5     !                                                                      !
6     !  Purpose: Check the chemical rxns namelist variables.                !
7     !                                                                      !
8     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
9           SUBROUTINE CHECK_ODEPACK_STIFF_CHEM
10     
11     ! Global Variables:
12     !---------------------------------------------------------------------//
13     ! Dimension of IJK arrays.
14           use param, only: DIMENSION_3
15     ! Double precision zero.
16           use param1, only: ZERO
17     ! Double precision value for undefined variables.
18           use param1, only: UNDEFINED
19     ! Constant gas phase density
20           use physprop, only : RO_G0
21     ! Runtime logical for solving energy equations.
22           use run, only: ENERGY_EQ
23     ! Run time logical for solving species equations.
24           use run, only: SPECIES_EQ
25     ! Net rate of gas phase production/consumption
26           use rxns, only: SUM_R_g
27     ! Net rate of solids phase production/consumption
28           use rxns, only: SUM_R_s
29     ! Run time logical for using stiff chemistry solver
30           use stiff_chem, only: STIFF_CHEMISTRY
31     ! Run time logicals for identifying cells owned by myPE
32           use stiff_chem, only: notOwner
33     
34     ! Full access to the following modules:
35     !---------------------------------------------------------------------//
36           use compar
37           use geometry
38           use indices
39     
40           use error_manager
41           use functions
42     
43           implicit none
44     
45     ! Local Variables:
46     !---------------------------------------------------------------------//
47           INTEGER :: I, J, K, IJK
48     
49           CALL INIT_ERR_MSG('CHECK_ODEPACK_STIFF_CHEM')
50     
51     ! Verify that there is sufficient run complexity to use the stiff solver
52           IF(STIFF_CHEMISTRY) THEN
53     
54     ! Energy equations must be solved.
55              IF(.NOT.ENERGY_EQ) THEN
56                 WRITE(ERR_MSG,1004)'ENERGY_EQ = .FALSE.'
57                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
58              ENDIF
59     
60     ! Must be compressible
61              IF(RO_G0 /= UNDEFINED) THEN
62                 WRITE(ERR_MSG,1004)'RO_G0 /= UNDEFINED'
63                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
64              ENDIF
65     
66              IF(.NOT.SPECIES_EQ(0)) THEN
67                 WRITE(ERR_MSG,1004)'SPECIES_EQ(0) = .FALSE.'
68                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
69              ENDIF
70     
71     ! The stiff chemistry solver only needs to loop over physical cells
72     ! owned by a process (e.g., not ghost cells). To avoid having a
73     ! triple do loop, this array is populated to identify the cells that
74     ! are not owned.
75              ALLOCATE( notOwner(DIMENSION_3) ); notOwner = .TRUE.
76              do k=kstart, kend
77              do j=jstart, jend
78              do i=istart, iend
79                 ijk = funijk(i,j,k)
80                 notOwner(IJK) = .FALSE.
81              enddo
82              enddo
83              enddo
84     
85     ! Initialize ODEPACK operating parameters.
86              CALL ODEPACK_INIT
87     
88     ! Clear the interphase mass transfer terms as the stiff solver
89     ! does no use them.
90              IF(allocated(SUM_R_g)) SUM_R_g = ZERO
91              IF(allocated(SUM_R_s)) SUM_R_s = ZERO
92     
93           ENDIF
94     
95           CALL FINL_ERR_MSG
96     
97     
98     
99      1004 FORMAT('Error 1004: ',                                           &
100              'Invalid parameters for stiff chemistry solver!',//           &
101              ' The following criteria must be satisfied:',/                &
102              '   > Solving the energy equations.',/                        &
103              '   > Compressible gas phase.',/                              &
104              '   > Solving gas phase species equations.',//                &
105              ' >>> Invalid Parameter: ',A,//                               &
106              'Check the user documentation and correct the mfix.dat file.')
107     
108           RETURN
109           END SUBROUTINE CHECK_ODEPACK_STIFF_CHEM
110     
111     
112     
113     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
114     !                                                                         C
115     !     Module name: MCHEM_ODEPACK_INIT                                     C
116     !     Purpose: controlling values for ODEAPCK(reference to ODEPACK manual)C
117     !                                                                         C
118     !     Author: Nan Xie                                   Date: 02-Aug-04   C
119     !     Reviewer:                                         Date:             C
120     !                                                                         C
121     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
122           SUBROUTINE ODEPACK_INIT
123     
124     ! Global Variables:
125     !---------------------------------------------------------------------//
126           use physprop, only: MMAX
127           use physprop, only: NMAX
128     
129           use run, only: SPECIES_EQ
130     
131           use stiff_chem, only: NEQ_DIMN
132           use stiff_chem, only: ODE_DIMN_g
133           use stiff_chem, only: ODE_DIMN_all
134     
135           use stiff_chem, only: ODE_ITOL
136           use stiff_chem, only: ODE_RTOL
137           use stiff_chem, only: ODE_ATOL
138     
139           use stiff_chem, only: ODE_LRW
140           use stiff_chem, only: ODE_LIW
141           use stiff_chem, only: ODE_JT
142     
143           use stiff_chem, only: STIFF_CHEM_MAX_STEPS
144     
145           use stiff_chem_dbg, only: ALLOCATE_STIFF_CHEM_DBG
146           use stiff_chem_stats, only: ALLOCATE_STIFF_CHEM_STATS
147     
148     ! Double precision value for undefined variables.
149           use param1, only: UNDEFINED_I
150     
151           implicit none
152     
153     
154     ! Local Variables:
155     !---------------------------------------------------------------------//
156           INTEGER :: LRN, LRS
157           INTEGER :: M
158     
159           LOGICAL :: LIMIT_MAX_STEPS
160     
161     ! Number of ODEs (maximum)
162           NEQ_DIMN = 2 + MMAX
163     
164     ! Dimension of ODEs solved if only the gas phase is present:
165     ! Gas density, temperature, and species
166           ODE_DIMN_g = 2 + NMAX(0)
167     ! Solids temperature excluding for all phases.
168           ODE_DIMN_g = ODE_DIMN_g + MMAX
169     
170     ! Calculate the total number of ODEs that are solve.
171           ODE_DIMN_all = ODE_DIMN_g
172           DO M=1, MMAX
173     ! Solids bulk density and species.
174              IF(SPECIES_EQ(M)) ODE_DIMN_all = ODE_DIMN_all + (1 + NMAX(M))
175           ENDDO
176     
177     ! Indicates type of Error control.
178           ODE_ITOL = 2 ! :: EWT(i) = RTOL * ABS(Y(i)) * ATOL(i)
179     
180     ! Relative error tolerance parameter.
181           ODE_RTOL(1) = 1.0D-5
182     
183     ! Absolute error tolerance parameter.
184           IF(.NOT.(allocated(ODE_ATOL))) allocate(ODE_ATOL(ODE_DIMN_all))
185           ODE_ATOL(:) = 1.0D-6  ! All Equations
186     
187     ! Declared length of RWORK.
188           LRN = 20 + 16*ODE_DIMN_all
189           LRS = 22 + 9*ODE_DIMN_all + (ODE_DIMN_all**2)
190           ODE_LRW = max(LRN, LRS)
191     
192     ! Declared length of IWORK.
193           ODE_LIW = 20 + ODE_DIMN_all
194     
195     ! Jacobian type indicator.
196           ODE_JT = 2 ! Internally generated.
197     
198     
199           IF(STIFF_CHEM_MAX_STEPS == UNDEFINED_I) THEN
200              STIFF_CHEM_MAX_STEPS = 500000
201              LIMIT_MAX_STEPS = .FALSE.
202           ELSE
203              LIMIT_MAX_STEPS = .TRUE.
204           ENDIF
205     
206     
207           CALL ALLOCATE_STIFF_CHEM_DBG(ODE_DIMN_all)
208           CALL ALLOCATE_STIFF_CHEM_STATS
209     
210           return
211           END SUBROUTINE ODEPACK_INIT
212