srcsink_sandbox_base.F90 coverage: 71.43 %func 66.67 %block
1) module SrcSink_Sandbox_Base_class
2)
3) use PFLOTRAN_Constants_module
4) use Geometry_module
5)
6) implicit none
7)
8) private
9)
10) #include "petsc/finclude/petscsys.h"
11)
12) type, abstract, public :: srcsink_sandbox_base_type
13) PetscInt :: local_cell_id
14) type(point3d_type) :: coordinate
15) PetscReal, pointer :: instantaneous_mass_rate(:)
16) PetscReal, pointer :: cumulative_mass(:)
17) class(srcsink_sandbox_base_type), pointer :: next
18) contains
19) procedure, public :: ReadInput => SSSandboxBaseRead
20) procedure, public :: Setup => SSSandboxBaseSetup
21) procedure, public :: Update => SSSandboxBaseUpdate
22) procedure, public :: Evaluate => SSSandboxBaseEvaluate
23) procedure, public :: Destroy => SSSandboxBaseDestroy
24) end type srcsink_sandbox_base_type
25)
26) public :: SSSandboxBaseInit, &
27) SSSandboxBaseSetup, &
28) SSSandboxBaseRead, &
29) SSSandboxBaseSelectCase, &
30) SSSandboxBaseDestroy
31)
32) contains
33)
34) ! ************************************************************************** !
35)
36) subroutine SSSandboxBaseInit(this)
37)
38) implicit none
39)
40) class(srcsink_sandbox_base_type) :: this
41)
42) this%coordinate%x = UNINITIALIZED_DOUBLE
43) this%coordinate%y = UNINITIALIZED_DOUBLE
44) this%coordinate%z = UNINITIALIZED_DOUBLE
45) this%local_cell_id = UNINITIALIZED_INTEGER
46) nullify(this%instantaneous_mass_rate)
47) nullify(this%cumulative_mass)
48) nullify(this%next)
49)
50) end subroutine SSSandboxBaseInit
51)
52) ! ************************************************************************** !
53)
54) subroutine SSSandboxBaseSetup(this,grid,option)
55)
56) use Option_module
57) use Grid_module
58)
59) implicit none
60)
61) class(srcsink_sandbox_base_type) :: this
62) type(grid_type) :: grid
63) type(option_type) :: option
64)
65) PetscInt :: local_id
66) PetscInt :: i, iflag
67) PetscErrorCode :: ierr
68)
69) if (Initialized(this%coordinate%x)) then
70) call GridGetLocalIDFromCoordinate(grid,this%coordinate,option,local_id)
71) iflag = 0
72) if (local_id > 0) then
73) this%local_cell_id = local_id
74) iflag = 1
75) endif
76) call MPI_Allreduce(iflag,i,ONE_INTEGER_MPI,MPIU_INTEGER,MPI_MAX, &
77) option%mycomm,ierr)
78) iflag = i
79) if (iflag > 1) then
80) option%io_buffer = 'More than one grid cell mapped in SSSandboxBaseSetup.'
81) call printErrMsg(option)
82) else if (iflag == 0) then
83) option%io_buffer = 'No grid cells mapped in SSSandboxBaseSetup.'
84) call printErrMsg(option)
85) endif
86) else
87) option%io_buffer = 'Source/sink in SSSandbox not associate with the &
88) &domain thorugh either a REGION or COORDINATE.'
89) call printErrMsg(option)
90) endif
91)
92) end subroutine SSSandboxBaseSetup
93)
94) ! ************************************************************************** !
95)
96) subroutine SSSandboxBaseRead(this,input,option)
97)
98) use Option_module
99) use Input_Aux_module
100)
101) implicit none
102)
103) class(srcsink_sandbox_base_type) :: this
104) type(input_type), pointer :: input
105) type(option_type) :: option
106)
107) end subroutine SSSandboxBaseRead
108)
109) ! ************************************************************************** !
110)
111) subroutine SSSandboxBaseSelectCase(this,input,option,keyword,found)
112)
113) use Option_module
114) use Input_Aux_module
115) use Geometry_module
116)
117) implicit none
118)
119) class(srcsink_sandbox_base_type) :: this
120) type(input_type), pointer :: input
121) type(option_type) :: option
122) character(len=MAXWORDLENGTH) :: keyword
123) PetscBool :: found
124)
125) character(len=MAXSTRINGLENGTH) :: error_string
126)
127) error_string = 'SOURCE_SINK_SANDBOX'
128)
129) found = PETSC_TRUE
130) select case(trim(keyword))
131) case('REGION')
132) option%io_buffer = 'The REGION card has been deprecated in &
133) &Source/Sink Sandbox. Please switch to using a COORDINATE and &
134) &defining one Src/Sink block for each coordinate.'
135) call printErrMsg(option)
136) case('COORDINATE')
137) call GeometryReadCoordinate(input,option,this%coordinate,error_string)
138) case default
139) found = PETSC_FALSE
140) end select
141)
142) end subroutine SSSandboxBaseSelectCase
143)
144) ! ************************************************************************** !
145)
146) subroutine SSSandboxBaseUpdate(this,time,option)
147)
148) use Option_module
149)
150) implicit none
151)
152) class(srcsink_sandbox_base_type) :: this
153) PetscReal :: time
154) type(option_type) :: option
155)
156) if (associated(this%cumulative_mass)) then
157) this%cumulative_mass(:) = this%cumulative_mass(:) + &
158) option%flow_dt*this%instantaneous_mass_rate(:)
159) endif
160)
161) end subroutine SSSandboxBaseUpdate
162)
163) ! ************************************************************************** !
164)
165) subroutine SSSandboxBaseEvaluate(this,Residual,Jacobian,compute_derivative, &
166) material_auxvar,aux_real,option)
167)
168) use Option_module
169) use Material_Aux_class
170)
171) implicit none
172)
173) class(srcsink_sandbox_base_type) :: this
174) type(option_type) :: option
175) PetscBool :: compute_derivative
176) PetscReal :: Residual(option%nflowdof)
177) PetscReal :: Jacobian(option%nflowdof,option%nflowdof)
178) class(material_auxvar_type) :: material_auxvar
179) PetscReal :: aux_real(:)
180)
181) end subroutine SSSandboxBaseEvaluate
182)
183) ! ************************************************************************** !
184)
185) subroutine SSSandboxBaseDestroy(this)
186)
187) use Utility_module
188)
189) implicit none
190)
191) class(srcsink_sandbox_base_type) :: this
192)
193) call DeallocateArray(this%instantaneous_mass_rate)
194) call DeallocateArray(this%cumulative_mass)
195)
196) end subroutine SSSandboxBaseDestroy
197)
198) end module SrcSink_Sandbox_Base_class