reaction_sandbox.F90 coverage: 100.00 %func 81.74 %block
1) module Reaction_Sandbox_module
2)
3) use Reaction_Sandbox_Base_class
4) use Reaction_Sandbox_CLM_CN_class
5) use Reaction_Sandbox_UFD_WP_class
6) use Reaction_Sandbox_Example_class
7) use Reaction_Sandbox_Simple_class
8) use Reaction_Sandbox_Cyber_class
9)
10) ! Add new reacton sandbox classes here.
11)
12) use PFLOTRAN_Constants_module
13)
14) implicit none
15)
16) private
17)
18) #include "petsc/finclude/petscsys.h"
19)
20) class(reaction_sandbox_base_type), pointer, public :: rxn_sandbox_list
21)
22) interface RSandboxRead
23) module procedure RSandboxRead1
24) module procedure RSandboxRead2
25) end interface
26)
27) interface RSandboxDestroy
28) module procedure RSandboxDestroy1
29) module procedure RSandboxDestroy2
30) end interface
31)
32) public :: RSandboxInit, &
33) RSandboxRead, &
34) RSandboxSkipInput, &
35) RSandboxSetup, &
36) RSandbox, &
37) RSandboxDestroy
38)
39) contains
40)
41) ! ************************************************************************** !
42)
43) subroutine RSandboxInit(option)
44) !
45) ! Initializes the sandbox list
46) !
47) ! Author: Glenn Hammond
48) ! Date: 01/28/13
49) !
50) use Option_module
51) implicit none
52) type(option_type) :: option
53)
54) if (associated(rxn_sandbox_list)) then
55) call RSandboxDestroy()
56) endif
57) nullify(rxn_sandbox_list)
58)
59) end subroutine RSandboxInit
60)
61) ! ************************************************************************** !
62)
63) subroutine RSandboxSetup(reaction,option)
64) !
65) ! Calls all the initialization routines for all reactions in
66) ! the sandbox list
67) !
68) ! Author: Glenn Hammond
69) ! Date: 01/28/13
70) !
71)
72) use Option_module
73) use Reaction_Aux_module, only : reaction_type
74)
75) implicit none
76)
77) type(reaction_type) :: reaction
78) type(option_type) :: option
79)
80) class(reaction_sandbox_base_type), pointer :: cur_sandbox
81)
82) ! sandbox reactions
83) cur_sandbox => rxn_sandbox_list
84) do
85) if (.not.associated(cur_sandbox)) exit
86) call cur_sandbox%Setup(reaction,option)
87) cur_sandbox => cur_sandbox%next
88) enddo
89)
90) end subroutine RSandboxSetup
91)
92) ! ************************************************************************** !
93)
94) subroutine RSandboxRead1(input,option)
95) !
96) ! Reads input deck for reaction sandbox parameters
97) !
98) ! Author: Glenn Hammond
99) ! Date: 05/16/13
100) !
101)
102) use Option_module
103) use String_module
104) use Input_Aux_module
105) use Utility_module
106)
107) implicit none
108)
109) type(input_type), pointer :: input
110) type(option_type) :: option
111)
112) call RSandboxRead(rxn_sandbox_list,input,option)
113)
114) end subroutine RSandboxRead1
115)
116) ! ************************************************************************** !
117)
118) subroutine RSandboxRead2(local_sandbox_list,input,option)
119) !
120) ! RSandboxRead: Reads input deck for reaction sandbox parameters
121) !
122) ! Author: Glenn Hammond
123) ! Date: 11/08/12
124) !
125)
126) use Option_module
127) use String_module
128) use Input_Aux_module
129) use Utility_module
130)
131) implicit none
132)
133) class(reaction_sandbox_base_type), pointer :: local_sandbox_list
134) type(input_type), pointer :: input
135) type(option_type) :: option
136)
137) character(len=MAXSTRINGLENGTH) :: string
138) character(len=MAXWORDLENGTH) :: word
139) class(reaction_sandbox_base_type), pointer :: new_sandbox, cur_sandbox
140)
141) nullify(new_sandbox)
142) do
143) call InputReadPflotranString(input,option)
144) if (InputError(input)) exit
145) if (InputCheckExit(input,option)) exit
146)
147) call InputReadWord(input,option,word,PETSC_TRUE)
148) call InputErrorMsg(input,option,'keyword','CHEMISTRY,REACTION_SANDBOX')
149) call StringToUpper(word)
150)
151) select case(trim(word))
152) case('CLM-CN')
153) new_sandbox => CLM_CN_Create()
154) ! Add new cases statements for new reacton sandbox classes here.
155) case('UFD-WP')
156) new_sandbox => WastePackageCreate()
157) case('EXAMPLE')
158) new_sandbox => EXAMPLECreate()
159) case('SIMPLE')
160) new_sandbox => SimpleCreate()
161) case('CYBERNETIC')
162) new_sandbox => CyberCreate()
163) case default
164) call InputKeywordUnrecognized(word,'CHEMISTRY,REACTION_SANDBOX',option)
165) end select
166)
167) call new_sandbox%ReadInput(input,option)
168)
169) if (.not.associated(local_sandbox_list)) then
170) local_sandbox_list => new_sandbox
171) else
172) cur_sandbox => local_sandbox_list
173) do
174) if (.not.associated(cur_sandbox%next)) exit
175) cur_sandbox => cur_sandbox%next
176) enddo
177) cur_sandbox%next => new_sandbox
178) endif
179) enddo
180)
181) end subroutine RSandboxRead2
182)
183) ! ************************************************************************** !
184)
185) subroutine RSandboxSkipInput(input,option)
186) !
187) ! Intelligently skips over REACTION_SANDBOX block
188) !
189) ! Author: Glenn Hammond
190) ! Date: 02/04/13
191) !
192)
193) use Option_module
194) use String_module
195) use Input_Aux_module
196) use Utility_module
197)
198) implicit none
199)
200) type(input_type), pointer :: input
201) type(option_type) :: option
202)
203) class(reaction_sandbox_base_type), pointer :: dummy_list
204)
205) nullify(dummy_list)
206) call RSandboxRead(dummy_list,input,option)
207) call RSandboxDestroy(dummy_list)
208)
209) end subroutine RSandboxSkipInput
210)
211) ! ************************************************************************** !
212)
213) subroutine RSandbox(Residual,Jacobian,compute_derivative,rt_auxvar, &
214) global_auxvar,material_auxvar,reaction,option)
215) !
216) ! Evaluates reaction storing residual and/or Jacobian
217) !
218) ! Author: Glenn Hammond
219) ! Date: 11/08/12
220) !
221)
222) use Option_module
223) use Reaction_Aux_module
224) use Reactive_Transport_Aux_module
225) use Global_Aux_module
226) use Material_Aux_class, only: material_auxvar_type
227)
228) implicit none
229)
230) type(option_type) :: option
231) type(reaction_type) :: reaction
232) PetscBool :: compute_derivative
233) PetscReal :: Residual(reaction%ncomp)
234) PetscReal :: Jacobian(reaction%ncomp,reaction%ncomp)
235) type(reactive_transport_auxvar_type) :: rt_auxvar
236) type(global_auxvar_type) :: global_auxvar
237) class(material_auxvar_type) :: material_auxvar
238)
239) class(reaction_sandbox_base_type), pointer :: cur_reaction
240)
241) cur_reaction => rxn_sandbox_list
242) do
243) if (.not.associated(cur_reaction)) exit
244) ! select type(cur_reaction)
245) ! class is(reaction_sandbox_clm_cn_type)
246) call cur_reaction%Evaluate(Residual,Jacobian,compute_derivative, &
247) rt_auxvar,global_auxvar,material_auxvar, &
248) reaction,option)
249) ! end select
250) cur_reaction => cur_reaction%next
251) enddo
252)
253) end subroutine RSandbox
254)
255) ! ************************************************************************** !
256)
257) subroutine RSandboxDestroy1()
258) !
259) ! Destroys master sandbox list
260) !
261) ! Author: Glenn Hammond
262) ! Date: 05/16/13
263) !
264)
265) implicit none
266)
267) call RSandboxDestroy(rxn_sandbox_list)
268)
269) end subroutine RSandboxDestroy1
270)
271) ! ************************************************************************** !
272)
273) subroutine RSandboxDestroy2(local_sandbox_list)
274) !
275) ! Destroys arbitrary sandbox list
276) !
277) ! Author: Glenn Hammond
278) ! Date: 11/08/12
279) !
280)
281) implicit none
282)
283) class(reaction_sandbox_base_type), pointer :: local_sandbox_list
284)
285) class(reaction_sandbox_base_type), pointer :: cur_sandbox, prev_sandbox
286)
287) ! sandbox reactions
288) cur_sandbox => local_sandbox_list
289) do
290) if (.not.associated(cur_sandbox)) exit
291) prev_sandbox => cur_sandbox%next
292) call cur_sandbox%Destroy()
293) deallocate(cur_sandbox)
294) cur_sandbox => prev_sandbox
295) enddo
296)
297) end subroutine RSandboxDestroy2
298)
299) end module Reaction_Sandbox_module