srcsink_sandbox_wipp_gas.F90 coverage: 100.00 %func 98.53 %block
1) module SrcSink_Sandbox_WIPP_Gas_class
2)
3) ! Sandbox srcsink for WIPP gas generation source terms
4)
5) use PFLOTRAN_Constants_module
6) use SrcSink_Sandbox_Base_class
7)
8) implicit none
9)
10) private
11)
12) #include "petsc/finclude/petscsys.h"
13)
14) PetscInt, parameter, public :: WIPP_GAS_WATER_SATURATION_INDEX = 1
15) PetscInt, parameter, public :: WIPP_GAS_TEMPERATURE_INDEX = 2
16)
17) type, public, &
18) extends(srcsink_sandbox_base_type) :: srcsink_sandbox_wipp_gas_type
19) PetscReal :: inundated_corrosion_rate
20) PetscReal :: inundated_degradation_rate
21) PetscReal :: humid_corrosion_factor
22) PetscReal :: humid_degradation_factor
23) PetscReal :: h2_fe_ratio
24) PetscReal :: h2_ch2o_ratio
25) contains
26) procedure, public :: ReadInput => WIPPGasGenerationRead
27) procedure, public :: Setup => WIPPGasGenerationSetup
28) procedure, public :: Evaluate => WIPPGasGenerationSrcSink
29) procedure, public :: Destroy => WIPPGasGenerationDestroy
30) end type srcsink_sandbox_wipp_gas_type
31)
32) public :: WIPPGasGenerationCreate
33)
34) contains
35)
36) ! ************************************************************************** !
37)
38) function WIPPGasGenerationCreate()
39) !
40) ! Allocates WIPP gas generation src/sink object.
41) !
42) ! Author: Glenn Hammond, Edit: Heeho Park
43) ! Date: 04/11/14, 05/15/14
44) !
45)
46) implicit none
47)
48) class(srcsink_sandbox_wipp_gas_type), pointer :: WIPPGasGenerationCreate
49)
50) allocate(WIPPGasGenerationCreate)
51) call SSSandboxBaseInit(WIPPGasGenerationCreate)
52) WIPPGasGenerationCreate%inundated_corrosion_rate = 3.0d-8
53) WIPPGasGenerationCreate%inundated_degradation_rate = 1.5d-7
54) WIPPGasGenerationCreate%humid_corrosion_factor = 1.0d-3
55) WIPPGasGenerationCreate%humid_degradation_factor = 0.2d0
56) WIPPGasGenerationCreate%h2_fe_ratio = 1.3081d0
57) WIPPGasGenerationCreate%h2_ch2o_ratio = 1.1100d0
58) nullify(WIPPGasGenerationCreate%next)
59)
60) end function WIPPGasGenerationCreate
61)
62) ! ************************************************************************** !
63)
64) subroutine WIPPGasGenerationRead(this,input,option)
65) !
66) ! Reads input deck for WIPP gas generation src/sink parameters
67) !
68) ! Author: Glenn Hammond, Edit: Heeho Park
69) ! Date: 04/11/14, 05/15/14
70) !
71)
72) use Option_module
73) use String_module
74) use Input_Aux_module
75) use Units_module, only : UnitsConvertToInternal
76)
77) implicit none
78)
79) class(srcsink_sandbox_wipp_gas_type) :: this
80) type(input_type), pointer :: input
81) type(option_type) :: option
82)
83) PetscInt :: i
84) character(len=MAXWORDLENGTH) :: word, internal_units
85) PetscBool :: found
86)
87) do
88) call InputReadPflotranString(input,option)
89) if (InputError(input)) exit
90) if (InputCheckExit(input,option)) exit
91)
92) call InputReadWord(input,option,word,PETSC_TRUE)
93) call InputErrorMsg(input,option,'keyword', &
94) 'SRCSINK_SANDBOX,WIPP')
95) call StringToUpper(word)
96)
97) call SSSandboxBaseSelectCase(this,input,option,word,found)
98) if (found) cycle
99)
100) select case(trim(word))
101) case('INUNDATED_CORROSION_RATE')
102) internal_units = 'unitless/sec'
103) call InputReadDouble(input,option,this%inundated_corrosion_rate)
104) call InputDefaultMsg(input,option,'inundated_corrosion_rate')
105) case('INUNDATED_DEGRADATION_RATE')
106) internal_units = 'unitless/sec'
107) call InputReadDouble(input,option,this%inundated_degradation_rate)
108) call InputDefaultMsg(input,option,'inundated_degradation_rate')
109) case('HUMID_CORROSION_FACTOR')
110) call InputReadDouble(input,option,this%humid_corrosion_factor)
111) call InputDefaultMsg(input,option,'humid_corrosion_factor')
112) case('HUMID_DEGREDATION_FACTOR')
113) call InputReadDouble(input,option,this%humid_degradation_factor)
114) call InputDefaultMsg(input,option,'humid_degradation_factor')
115) case('H2_FE_RATIO')
116) call InputReadDouble(input,option,this%h2_fe_ratio)
117) call InputDefaultMsg(input,option,'h2_fe_ratio')
118) case('H2_CH2O_RATIO')
119) call InputReadDouble(input,option,this%h2_ch2o_ratio)
120) call InputDefaultMsg(input,option,'h2_ch2o_ratio')
121) case default
122) call InputKeywordUnrecognized(word, &
123) 'SRCSINK_SANDBOX,WIPP-GAS_GENERATION',option)
124) end select
125) enddo
126)
127) end subroutine WIPPGasGenerationRead
128)
129) ! ************************************************************************** !
130)
131) subroutine WIPPGasGenerationSetup(this,grid,option)
132) !
133) ! Sets up the WIPP gas generation src/sink
134) !
135) ! Author: Glenn Hammond
136) ! Date: 04/11/14
137) use Option_module
138) use Grid_module
139)
140) implicit none
141)
142) class(srcsink_sandbox_wipp_gas_type) :: this
143) type(grid_type) :: grid
144) type(option_type) :: option
145)
146) call SSSandboxBaseSetup(this,grid,option)
147)
148) end subroutine WIPPGasGenerationSetup
149)
150) ! ************************************************************************** !
151)
152) subroutine WIPPGasGenerationSrcSink(this,Residual,Jacobian, &
153) compute_derivative, &
154) material_auxvar,aux_real,option)
155) !
156) ! Evaluates src/sink storing residual and/or Jacobian
157) !
158) ! Author: Glenn Hammond, Edited: Heeho Park
159) ! Date: 04/11/14, 05/15/14
160) !
161)
162) use Option_module
163) use Reaction_Aux_module
164) use Material_Aux_class
165) use EOS_Gas_module
166)
167) implicit none
168)
169) class(srcsink_sandbox_wipp_gas_type) :: this
170) type(option_type) :: option
171) PetscBool :: compute_derivative
172) PetscReal :: Residual(option%nflowdof)
173) PetscReal :: Jacobian(option%nflowdof,option%nflowdof)
174) class(material_auxvar_type) :: material_auxvar
175) PetscReal :: aux_real(:)
176)
177) !gas generation calculation variables
178) PetscReal :: water_saturation
179) PetscReal :: humid_corrosion_rate
180) PetscReal :: humid_degradation_rate
181) PetscReal :: corrosion_gas_rate
182) PetscReal :: degradation_gas_rate
183) PetscReal :: gas_generation_rate
184)
185) !Enthalpy calculation variables
186) PetscReal :: T ! temperature [C]
187) PetscReal :: dummy_P ! pressure [Pa]
188) PetscReal :: H ! enthalpy [J/kmol]
189) PetscReal :: dH_dT ! derivative enthalpy wrt temperature
190) PetscReal :: dH_dP ! derivative enthalpy wrt pressure
191) PetscReal :: U ! internal energy [J/kmol]
192) PetscReal :: dU_dT ! deriv. internal energy wrt temperature
193) PetscReal :: dU_dP ! deriv. internal energy wrt pressure
194) PetscErrorCode :: ierr
195)
196) ! PetscReal :: test_value
197)
198) ! call water saturation on each cell
199) ! equations are given in V&V doc.
200) water_saturation = aux_real(WIPP_GAS_WATER_SATURATION_INDEX)
201) humid_corrosion_rate = this%humid_corrosion_factor * &
202) this%inundated_corrosion_rate
203) humid_degradation_rate = this%humid_degradation_factor * &
204) this%inundated_degradation_rate
205) corrosion_gas_rate = this%inundated_corrosion_rate * &
206) water_saturation + &
207) humid_corrosion_rate * (1.d0-water_saturation)
208) degradation_gas_rate = this%inundated_degradation_rate * &
209) water_saturation + humid_degradation_rate * &
210) (1.d0-water_saturation)
211) ! gas generation unit: mol of H2 / (m^3 * s)
212) gas_generation_rate = this%h2_fe_ratio * corrosion_gas_rate + &
213) this%h2_ch2o_ratio * degradation_gas_rate
214) ! convert to kmol/s -> mol/(m^3 * s) * (volume of the cell) *
215) ! (1 kmol/1000 mol)
216) gas_generation_rate = gas_generation_rate * material_auxvar%volume * 1.d-3
217) ! test_value = gas_generation_rate * 2.01588d0 / material_auxvar%volume
218)
219) ! positive is inflow
220) ! kmol/s
221) Residual(TWO_INTEGER) = gas_generation_rate
222)
223) T = aux_real(WIPP_GAS_TEMPERATURE_INDEX)
224) call EOSGasEnergy(T,dummy_P,H,dH_dT,dH_dP,U,dU_dT,dU_dP,ierr)
225) ! energy equation
226) ! units = MJ/s -> enthalpy(J/kmol) * (MJ/1E+6J) * gas_generation_rate (kmol/s)
227) Residual(THREE_INTEGER) = H * 1.d-6 * gas_generation_rate
228)
229) if (compute_derivative) then
230)
231) ! jacobian something
232)
233) endif
234)
235) end subroutine WIPPGasGenerationSrcSink
236)
237) ! ************************************************************************** !
238)
239) subroutine WIPPGasGenerationDestroy(this)
240) !
241) ! Destroys allocatable or pointer objects created in this
242) ! module
243) !
244) ! Author: Glenn Hammond
245) ! Date: 04/11/14
246) !
247)
248) implicit none
249)
250) class(srcsink_sandbox_wipp_gas_type) :: this
251)
252) call SSSandboxBaseDestroy(this)
253)
254) end subroutine WIPPGasGenerationDestroy
255)
256) end module SrcSink_Sandbox_WIPP_Gas_class