srcsink_sandbox_mass_rate.F90 coverage: 0.00 %func 0.00 %block
1) module SrcSink_Sandbox_Mass_Rate_class
2)
3) ! Sandbox srcsink for mass rate source terms. This source/sink is identical
4) ! to the standard mass rate source/sink in PFLOTRAN, but is more for
5) ! illustrative purposes
6)
7) use PFLOTRAN_Constants_module
8) use SrcSink_Sandbox_Base_class
9)
10) implicit none
11)
12) private
13)
14) #include "petsc/finclude/petscsys.h"
15)
16) type, public, &
17) extends(srcsink_sandbox_base_type) :: srcsink_sandbox_mass_rate_type
18) PetscReal, pointer :: rate(:)
19) contains
20) procedure, public :: ReadInput => MassRateRead
21) procedure, public :: Setup => MassRateSetup
22) procedure, public :: Evaluate => MassRateSrcSink
23) procedure, public :: Destroy => MassRateDestroy
24) end type srcsink_sandbox_mass_rate_type
25)
26) public :: MassRateCreate
27)
28) contains
29)
30) ! ************************************************************************** !
31)
32) function MassRateCreate()
33) !
34) ! Allocates mass rate src/sink object.
35) !
36) ! Author: Glenn Hammond
37) ! Date: 05/06/14
38)
39) implicit none
40)
41) class(srcsink_sandbox_mass_rate_type), pointer :: MassRateCreate
42)
43) allocate(MassRateCreate)
44) call SSSandboxBaseInit(MassRateCreate)
45) nullify(MassRateCreate%rate)
46)
47)
48) end function MassRateCreate
49)
50) ! ************************************************************************** !
51)
52) subroutine MassRateRead(this,input,option)
53) !
54) ! Reads input deck for mass rate src/sink parameters
55) !
56) ! Author: Glenn Hammond
57) ! Date: 05/06/14
58)
59) use Option_module
60) use String_module
61) use Input_Aux_module
62) use Units_module, only : UnitsConvertToInternal
63)
64) implicit none
65)
66) class(srcsink_sandbox_mass_rate_type) :: this
67) type(input_type), pointer :: input
68) type(option_type) :: option
69)
70) PetscInt :: i
71) character(len=MAXWORDLENGTH) :: word, internal_units
72) PetscBool :: found
73)
74) do
75) call InputReadPflotranString(input,option)
76) if (InputError(input)) exit
77) if (InputCheckExit(input,option)) exit
78)
79) call InputReadWord(input,option,word,PETSC_TRUE)
80) call InputErrorMsg(input,option,'keyword', &
81) 'SOURCE_SINK_SANDBOX,MASS_RATE')
82) call StringToUpper(word)
83)
84) call SSSandboxBaseSelectCase(this,input,option,word,found)
85) if (found) cycle
86)
87) select case(trim(word))
88)
89) case('RATE')
90) allocate(this%rate(option%nflowdof))
91) do i = 1, option%nflowdof
92) word = ''
93) select case(i)
94) case(ONE_INTEGER)
95) word = 'Liquid Component Rate'
96) internal_units = 'kg/sec'
97) case(TWO_INTEGER)
98) word = 'Gas Component Rate'
99) internal_units = 'kg/sec'
100) case(THREE_INTEGER)
101) word = 'Energy Rate'
102) internal_units = 'MW|MJ/sec'
103) case default
104) write(word,*) i
105) option%io_buffer = 'Unknown dof #' // trim(adjustl(word)) // &
106) ' in MassRateRead.'
107) call printErrMsg(option)
108) end select
109) call InputReadDouble(input,option,this%rate(i))
110) call InputErrorMsg(input,option,word,'SOURCE_SINK_SANDBOX,MASS_RATE')
111) call InputReadWord(input,option,word,PETSC_TRUE)
112) if (input%ierr == 0) then
113) this%rate(i) = this%rate(i) * &
114) UnitsConvertToInternal(word,internal_units,option)
115) endif
116) enddo
117) case default
118) call InputKeywordUnrecognized(word,'SRCSINK_SANDBOX,MASS_RATE',option)
119) end select
120) enddo
121)
122) end subroutine MassRateRead
123)
124) ! ************************************************************************** !
125)
126) subroutine MassRateSetup(this,grid,option)
127) !
128) ! Sets up the mass rate src/sink
129) !
130) ! Author: Glenn Hammond
131) ! Date: 05/06/14
132)
133) use Option_module
134) use Grid_module
135) use General_Aux_module, only : general_fmw_com => fmw_comp
136)
137) implicit none
138)
139) class(srcsink_sandbox_mass_rate_type) :: this
140) type(grid_type) :: grid
141) type(option_type) :: option
142)
143) call SSSandboxBaseSetup(this,grid,option)
144) ! convert rate from kg/s to mol/s
145) select case(option%iflowmode)
146) case(RICHARDS_MODE)
147) this%rate(1) = this%rate(1) / FMWH2O
148) case(G_MODE)
149) this%rate(1) = this%rate(1) / general_fmw_com(1)
150) this%rate(2) = this%rate(2) / general_fmw_com(2)
151) case default
152) option%io_buffer = 'Rate conversion not set up for flow mode in ' // &
153) 'MassRateSetup'
154) call printErrMsg(option)
155) end select
156)
157) end subroutine MassRateSetup
158)
159) ! ************************************************************************** !
160)
161) subroutine MassRateSrcSink(this,Residual,Jacobian,compute_derivative, &
162) material_auxvar,aux_real,option)
163) !
164) ! Evaluates src/sink storing residual and/or Jacobian
165) !
166) ! Author: Glenn Hammond
167) ! Date: 05/06/14
168)
169) use Option_module
170) use Reaction_Aux_module
171) use Material_Aux_class
172)
173) implicit none
174)
175) class(srcsink_sandbox_mass_rate_type) :: this
176) type(option_type) :: option
177) PetscBool :: compute_derivative
178) PetscReal :: Residual(option%nflowdof)
179) PetscReal :: Jacobian(option%nflowdof,option%nflowdof)
180) class(material_auxvar_type) :: material_auxvar
181) PetscReal :: aux_real(:)
182)
183) PetscInt :: idof
184)
185) do idof = 1, option%nflowdof
186) Residual(idof) = this%rate(idof)
187) enddo
188)
189) if (compute_derivative) then
190)
191) ! since the rates are constant, there is no derivative
192)
193) endif
194)
195) end subroutine MassRateSrcSink
196)
197) ! ************************************************************************** !
198)
199) subroutine MassRateDestroy(this)
200) !
201) ! Destroys allocatable or pointer objects created in this module
202) !
203) ! Author: Glenn Hammond
204) ! Date: 05/06/14
205)
206) implicit none
207)
208) class(srcsink_sandbox_mass_rate_type) :: this
209)
210) call SSSandboxBaseDestroy(this)
211) deallocate(this%rate)
212) nullify(this%rate)
213)
214) end subroutine MassRateDestroy
215)
216) end module SrcSink_Sandbox_Mass_Rate_class