integral_flux.F90 coverage: 30.00 %func 17.54 %block
1) module Integral_Flux_module
2)
3) use Geometry_module
4)
5) use PFLOTRAN_Constants_module
6)
7) implicit none
8)
9) private
10)
11) #include "petsc/finclude/petscsys.h"
12)
13) PetscInt, parameter, public :: INTEGRATE_FLOW = 1
14) PetscInt, parameter, public :: INTEGRATE_TRANSPORT = 2
15)
16) type, public :: integral_flux_type
17) PetscInt :: id
18) character(len=MAXWORDLENGTH) :: name
19) type(point3d_type), pointer :: coordinates(:)
20) PetscBool :: invert_direction
21) PetscInt, pointer :: connections(:)
22) PetscReal, pointer :: integral_value(:)
23) type(integral_flux_type), pointer :: next
24) end type integral_flux_type
25)
26) type, public :: integral_flux_list_type
27) PetscInt :: num_integral_fluxes
28) type(integral_flux_type), pointer :: first
29) type(integral_flux_type), pointer :: last
30) type(integral_flux_type), pointer :: array(:)
31) end type integral_flux_list_type
32)
33) public :: IntegralFluxCreate, &
34) IntegralFluxDestroy, &
35) IntegralFluxRead, &
36) IntegralFluxAddToList, &
37) IntegralFluxInitList, &
38) IntegralFluxDestroyList, &
39) IntegralFluxGetPtrFromList, &
40) IntegralFluxSizeStorage, &
41) IntegralFluxUpdate, &
42) IntegralFluxGetInstantaneous
43)
44) contains
45)
46) ! ************************************************************************** !
47)
48) function IntegralFluxCreate()
49) !
50) ! Create object that stores integral flux information
51) !
52) ! Author: Glenn Hammond
53) ! Date: 10/20/14
54) !
55)
56) implicit none
57)
58) type(integral_flux_type), pointer :: IntegralFluxCreate
59)
60) type(integral_flux_type), pointer :: integral_flux
61)
62) allocate(integral_flux)
63)
64) integral_flux%name = ''
65) integral_flux%id = 0
66) integral_flux%invert_direction = PETSC_FALSE
67) nullify(integral_flux%connections)
68) nullify(integral_flux%coordinates)
69) nullify(integral_flux%integral_value)
70) nullify(integral_flux%next)
71)
72) IntegralFluxCreate => integral_flux
73)
74) end function IntegralFluxCreate
75)
76) ! ************************************************************************** !
77)
78) subroutine IntegralFluxRead(integral_flux,input,option)
79) !
80) ! Reads integral flux data from the input file
81) !
82) ! Author: Glenn Hammond
83) ! Date: 10/20/14
84) !
85)
86) use Input_Aux_module
87) use String_module
88) use Option_module
89)
90) implicit none
91)
92) type(integral_flux_type) :: integral_flux
93) type(input_type), pointer :: input
94) type(option_type) :: option
95)
96) character(len=MAXWORDLENGTH) :: keyword
97)
98) input%ierr = 0
99) do
100)
101) call InputReadPflotranString(input,option)
102)
103) if (InputCheckExit(input,option)) exit
104)
105) call InputReadWord(input,option,keyword,PETSC_TRUE)
106) call InputErrorMsg(input,option,'keyword','integral_flux')
107)
108) select case(trim(keyword))
109) case('NAME')
110) call InputReadWord(input,option,integral_flux%name,PETSC_TRUE)
111) call InputErrorMsg(input,option,'name','INTEGRAL_FLUX')
112) case('INVERT_DIRECTION')
113) integral_flux%invert_direction = PETSC_TRUE
114) case('COORDINATES')
115) call GeometryReadCoordinates(input,option,integral_flux%name, &
116) integral_flux%coordinates)
117) case default
118) call InputKeywordUnrecognized(keyword,'INTEGRAL_FLUX',option)
119) end select
120)
121) enddo
122)
123) if (len_trim(integral_flux%name) < 1) then
124) option%io_buffer = 'All INTEGRAL_FLUXes must have a name.'
125) call printErrMsg(option)
126) endif
127)
128) end subroutine IntegralFluxRead
129)
130) ! ************************************************************************** !
131)
132) subroutine IntegralFluxSizeStorage(integral_flux,option)
133) !
134) ! Sizes the arrays that store the integrated flux
135) !
136) ! Author: Glenn Hammond
137) ! Date: 10/20/14
138) !
139)
140) use Option_module
141)
142) implicit none
143)
144) type(integral_flux_type) :: integral_flux
145) type(option_type) :: option
146)
147) allocate(integral_flux%integral_value(option%nflowdof+option%ntrandof))
148) integral_flux%integral_value = 0.d0
149)
150) end subroutine IntegralFluxSizeStorage
151)
152) ! ************************************************************************** !
153)
154) subroutine IntegralFluxUpdate(integral_flux_list,internal_fluxes, &
155) boundary_fluxes,iflag,option)
156) !
157) ! Updates the stored integrated value of each integral flux measurement
158) !
159) ! Author: Glenn Hammond
160) ! Date: 10/20/14
161) !
162)
163) use Option_module
164)
165) implicit none
166)
167) type(integral_flux_list_type) :: integral_flux_list
168) PetscReal :: internal_fluxes(:,:)
169) PetscReal :: boundary_fluxes(:,:)
170) PetscInt :: iflag
171) type(option_type) :: option
172)
173) type(integral_flux_type), pointer :: integral_flux
174) PetscReal, allocatable :: sum_array(:)
175) PetscInt :: offset
176) PetscInt :: num_values
177) PetscReal :: dt
178)
179) if (.not.associated(integral_flux_list%first)) return
180)
181) select case(iflag)
182) case(INTEGRATE_FLOW)
183) offset = 0
184) num_values = option%nflowdof
185) dt = option%flow_dt
186) case(INTEGRATE_TRANSPORT)
187) offset = option%nflowdof
188) num_values = option%ntrandof
189) dt = option%tran_dt
190) case default
191) offset = -1 ! to catch bugs
192) end select
193)
194) allocate(sum_array(num_values))
195) integral_flux => integral_flux_list%first
196) do
197) if (.not.associated(integral_flux)) exit
198) sum_array = 0.d0
199) call IntegralFluxGetInstantaneous(integral_flux, internal_fluxes, &
200) boundary_fluxes,num_values, &
201) sum_array,option)
202) integral_flux%integral_value(offset+1:offset+num_values) = &
203) integral_flux%integral_value(offset+1:offset+num_values) + &
204) sum_array(1:num_values)*dt
205) integral_flux => integral_flux%next
206) enddo
207) deallocate(sum_array)
208)
209) end subroutine IntegralFluxUpdate
210)
211)
212) ! ************************************************************************** !
213)
214) subroutine IntegralFluxGetInstantaneous(integral_flux, internal_fluxes, &
215) boundary_fluxes,num_values, &
216) sum_array,option)
217) !
218) ! Returns the instantaneous mole flux for an integral flux object
219) !
220) ! Author: Glenn Hammond
221) ! Date: 10/20/14
222) !
223)
224) use Option_module
225)
226) implicit none
227)
228) type(integral_flux_type) :: integral_flux
229) PetscReal :: internal_fluxes(:,:)
230) PetscReal :: boundary_fluxes(:,:)
231) PetscInt :: num_values
232) PetscReal :: sum_array(:)
233) type(option_type) :: option
234)
235) PetscInt :: i
236) PetscInt :: iconn
237)
238) sum_array = 0.d0
239) if (.not.associated(integral_flux%connections)) return
240)
241) do i = 1, size(integral_flux%connections)
242) iconn = integral_flux%connections(i)
243) if (iconn > 0) then ! internal
244) sum_array(1:num_values) = sum_array(1:num_values) + &
245) internal_fluxes(1:num_values,iconn)
246) else ! boundary
247) sum_array(1:num_values) = sum_array(1:num_values) + &
248) boundary_fluxes(1:num_values,-iconn)
249) endif
250) enddo
251) if (integral_flux%invert_direction) then
252) sum_array = -1.d0 * sum_array
253) endif
254)
255) end subroutine IntegralFluxGetInstantaneous
256)
257) ! ************************************************************************** !
258)
259) subroutine IntegralFluxInitList(list)
260) !
261) ! Initializes a integral flux list
262) !
263) ! Author: Glenn Hammond
264) ! Date: 10/20/14
265) !
266)
267) implicit none
268)
269) type(integral_flux_list_type) :: list
270)
271) nullify(list%first)
272) nullify(list%last)
273) nullify(list%array)
274) list%num_integral_fluxes = 0
275)
276) end subroutine IntegralFluxInitList
277)
278) ! ************************************************************************** !
279)
280) subroutine IntegralFluxAddToList(new_integral_flux,list)
281) !
282) ! Adds a new integral_flux to a integral flux list
283) !
284) ! Author: Glenn Hammond
285) ! Date: 10/20/14
286) !
287)
288) implicit none
289)
290) type(integral_flux_type), pointer :: new_integral_flux
291) type(integral_flux_list_type) :: list
292)
293) list%num_integral_fluxes = list%num_integral_fluxes + 1
294) new_integral_flux%id = list%num_integral_fluxes
295) if (.not.associated(list%first)) list%first => new_integral_flux
296) if (associated(list%last)) list%last%next => new_integral_flux
297) list%last => new_integral_flux
298)
299) end subroutine IntegralFluxAddToList
300)
301) ! ************************************************************************** !
302)
303) function IntegralFluxGetPtrFromList(integral_flux_name,integral_flux_list)
304) !
305) ! Returns a pointer to the integral flux matching integral_flux_name
306) !
307) ! Author: Glenn Hammond
308) ! Date: 10/20/14
309) !
310)
311) use String_module
312)
313) implicit none
314)
315) type(integral_flux_type), pointer :: IntegralFluxGetPtrFromList
316) character(len=MAXWORDLENGTH) :: integral_flux_name
317) type(integral_flux_list_type) :: integral_flux_list
318)
319) PetscInt :: length
320) type(integral_flux_type), pointer :: integral_flux
321)
322) nullify(IntegralFluxGetPtrFromList)
323) integral_flux => integral_flux_list%first
324)
325) do
326) if (.not.associated(integral_flux)) exit
327) length = len_trim(integral_flux_name)
328) if (length == len_trim(integral_flux%name) .and. &
329) StringCompare(integral_flux%name,integral_flux_name, &
330) length)) then
331) IntegralFluxGetPtrFromList => integral_flux
332) return
333) endif
334) integral_flux => integral_flux%next
335) enddo
336)
337) end function IntegralFluxGetPtrFromList
338)
339) ! ************************************************************************** !
340)
341) subroutine IntegralFluxDestroyList(integral_flux_list)
342) !
343) ! Deallocates a list of integral fluxes
344) !
345) ! Author: Glenn Hammond
346) ! Date: 10/20/14
347) !
348)
349) implicit none
350)
351) type(integral_flux_list_type), pointer :: integral_flux_list
352)
353) type(integral_flux_type), pointer :: integral_flux, prev_integral_flux
354)
355) if (.not.associated(integral_flux_list)) return
356)
357) integral_flux => integral_flux_list%first
358) do
359) if (.not.associated(integral_flux)) exit
360) prev_integral_flux => integral_flux
361) integral_flux => integral_flux%next
362) call IntegralFluxDestroy(prev_integral_flux)
363) enddo
364)
365) integral_flux_list%num_integral_fluxes = 0
366) nullify(integral_flux_list%first)
367) nullify(integral_flux_list%last)
368) if (associated(integral_flux_list%array)) deallocate(integral_flux_list%array)
369) nullify(integral_flux_list%array)
370)
371) deallocate(integral_flux_list)
372) nullify(integral_flux_list)
373)
374) end subroutine IntegralFluxDestroyList
375)
376) ! ************************************************************************** !
377)
378) subroutine IntegralFluxDestroy(integral_flux)
379) !
380) ! Deallocates a integral flux
381) !
382) ! Author: Glenn Hammond
383) ! Date: 10/20/14
384) !
385) use Utility_module
386)
387) implicit none
388)
389) type(integral_flux_type), pointer :: integral_flux
390)
391) PetscInt :: i
392)
393) if (.not.associated(integral_flux)) return
394)
395) if (associated(integral_flux%coordinates)) &
396) deallocate(integral_flux%coordinates)
397) nullify(integral_flux%coordinates)
398) call DeallocateArray(integral_flux%connections)
399) call DeallocateArray(integral_flux%integral_value)
400) deallocate(integral_flux)
401) nullify(integral_flux)
402)
403) end subroutine IntegralFluxDestroy
404)
405) end module Integral_Flux_module