connection.F90 coverage: 87.50 %func 89.20 %block
1) module Connection_module
2)
3) use PFLOTRAN_Constants_module
4)
5) implicit none
6)
7) #include "petsc/finclude/petscsys.h"
8)
9) private
10)
11) type, public :: connection_set_type
12) PetscInt :: id
13) PetscInt :: itype ! connection type (boundary, internal, source sink
14) PetscInt :: num_connections
15) PetscInt :: offset
16) PetscInt, pointer :: local(:) ! 1 if connection is local, 0 if connection is ghosted
17) PetscInt, pointer :: id_up(:) ! list of ids of upwind cells
18) PetscInt, pointer :: id_dn(:) ! list of ids of downwind cells
19) PetscInt, pointer :: id_up2(:) ! list of ids of 2nd upwind cells
20) PetscInt, pointer :: id_dn2(:) ! list of ids of 2nd downwind cells
21) PetscReal, pointer :: dist(:,:) ! list of distance vectors, size(-1:3,num_connections) where
22) ! -1 = fraction upwind
23) ! 0 = magnitude of distance
24) ! 1-3 = components of unit vector
25) PetscReal, pointer :: intercp(:,:) ! x,y,z location of intercept between the line connecting
26) ! upwind and downwind cells with the face shared by the cells
27) PetscReal, pointer :: area(:) ! list of areas of faces normal to distance vectors
28) PetscReal, pointer :: cntr(:,:) ! coordinates (1:3, num_connections) of the mass center of the face
29) PetscInt, pointer :: face_id(:) ! list of ids of faces (in local order)
30) type(connection_set_type), pointer :: next
31) end type connection_set_type
32)
33)
34) ! pointer data structure required for making an array of region pointers in F90
35) type, public :: connection_set_ptr_type
36) type(connection_set_type), pointer :: ptr ! pointer to the connection_set_type
37) end type connection_set_ptr_type
38)
39) type, public :: connection_set_list_type
40) PetscInt :: num_connection_objects
41) type(connection_set_type), pointer :: first
42) type(connection_set_type), pointer :: last
43) type(connection_set_ptr_type), pointer :: array(:)
44) end type connection_set_list_type
45)
46) public :: ConnectionCreate, &
47) ConnectionAddToList, &
48) ConnectionGetNumberInList, &
49) ConnectionInitList, &
50) ConnectionCalculateDistances, &
51) ConnectionDestroyList, &
52) ConnectionDestroy
53)
54) contains
55)
56) ! ************************************************************************** !
57)
58) function ConnectionCreate(num_connections,connection_itype)
59) !
60) ! Allocates and initializes a new connection
61) !
62) ! Author: Glenn Hammond
63) ! Date: 10/15/07
64) !
65)
66) implicit none
67)
68) PetscInt :: num_connections
69) PetscInt :: num_dof
70) PetscInt :: connection_itype
71)
72) type(connection_set_type), pointer :: ConnectionCreate
73)
74) type(connection_set_type), pointer :: connection
75)
76) allocate(connection)
77) connection%id = 0
78) connection%itype = connection_itype
79) connection%offset = 0
80) connection%num_connections = num_connections
81) nullify(connection%local)
82) nullify(connection%id_up)
83) nullify(connection%id_dn)
84) nullify(connection%id_up2)
85) nullify(connection%id_dn2)
86) nullify(connection%face_id)
87) nullify(connection%dist)
88) nullify(connection%intercp)
89) nullify(connection%area)
90) nullify(connection%cntr)
91) select case(connection_itype)
92) case(INTERNAL_CONNECTION_TYPE)
93) allocate(connection%id_up(num_connections))
94) allocate(connection%id_dn(num_connections))
95) allocate(connection%dist(-1:3,num_connections))
96) allocate(connection%intercp(1:3,num_connections))
97) allocate(connection%area(num_connections))
98) allocate(connection%face_id(num_connections))
99) connection%id_up = 0
100) connection%id_dn = 0
101) connection%face_id = 0
102) connection%dist = 0.d0
103) connection%intercp = 0.d0
104) connection%area = 0.d0
105) case(BOUNDARY_CONNECTION_TYPE)
106) allocate(connection%id_dn(num_connections))
107) allocate(connection%dist(-1:3,num_connections))
108) allocate(connection%intercp(1:3,num_connections))
109) allocate(connection%area(num_connections))
110) allocate(connection%face_id(num_connections))
111) connection%id_dn = 0
112) connection%dist = 0.d0
113) connection%intercp = 0.d0
114) connection%area = 0.d0
115) case(SRC_SINK_CONNECTION_TYPE,INITIAL_CONNECTION_TYPE)
116) allocate(connection%id_dn(num_connections))
117) connection%id_dn = 0
118) end select
119) nullify(connection%next)
120)
121) ConnectionCreate => connection
122)
123) end function ConnectionCreate
124)
125) ! ************************************************************************** !
126)
127) function ConnectionGetNumberInList(list)
128) !
129) ! Returns the number of connections in a list
130) !
131) ! Author: Glenn Hammond
132) ! Date: 11/19/07
133) !
134)
135) implicit none
136)
137) type(connection_set_list_type) :: list
138)
139) PetscInt :: ConnectionGetNumberInList
140) type(connection_set_type), pointer :: cur_connection_set
141)
142) ConnectionGetNumberInList = 0
143) cur_connection_set => list%first
144) do
145) if (.not.associated(cur_connection_set)) exit
146) ConnectionGetNumberInList = ConnectionGetNumberInList + &
147) cur_connection_set%num_connections
148) cur_connection_set => cur_connection_set%next
149) enddo
150)
151) end function ConnectionGetNumberInList
152)
153) ! ************************************************************************** !
154)
155) subroutine ConnectionInitList(list)
156) !
157) ! InitConnectionModule: Initializes module variables, lists, arrays.
158) !
159) ! Author: Glenn Hammond
160) ! Date: 10/15/07
161) !
162)
163) implicit none
164)
165) type(connection_set_list_type) :: list
166)
167) nullify(list%first)
168) nullify(list%last)
169) nullify(list%array)
170) list%num_connection_objects = 0
171)
172) end subroutine ConnectionInitList
173)
174) ! ************************************************************************** !
175)
176) subroutine ConnectionAddToList(new_connection_set,list)
177) !
178) ! Adds a new connection of the module global list of
179) ! connections
180) !
181) ! Author: Glenn Hammond
182) ! Date: 10/15/07
183) !
184)
185) implicit none
186)
187) type(connection_set_type), pointer :: new_connection_set
188) type(connection_set_list_type) :: list
189)
190) list%num_connection_objects = list%num_connection_objects + 1
191) new_connection_set%id = list%num_connection_objects
192) if (.not.associated(list%first)) list%first => new_connection_set
193) if (associated(list%last)) list%last%next => new_connection_set
194) list%last => new_connection_set
195)
196) end subroutine ConnectionAddToList
197)
198) ! ************************************************************************** !
199)
200) subroutine ConnectionConvertListToArray(list)
201) !
202) ! Creates an array of pointers to the
203) ! connections in the connection list
204) !
205) ! Author: Glenn Hammond
206) ! Date: 10/15/07
207) !
208)
209) implicit none
210)
211) type(connection_set_list_type) :: list
212)
213) PetscInt :: count
214) type(connection_set_type), pointer :: cur_connection_set
215)
216)
217) allocate(list%array(list%num_connection_objects))
218)
219) cur_connection_set => list%first
220) do
221) if (.not.associated(cur_connection_set)) exit
222) list%array(cur_connection_set%id)%ptr => cur_connection_set
223) cur_connection_set => cur_connection_set%next
224) enddo
225)
226) end subroutine ConnectionConvertListToArray
227)
228) ! ************************************************************************** !
229)
230) subroutine ConnectionCalculateDistances(dist,gravity,distance_upwind, &
231) distance_downwind,distance_gravity, &
232) upwind_weight)
233) !
234) ! Calculates the various distances and weights used in a flux calculation.
235) !
236) ! Author: Glenn Hammond
237) ! Date: 01/09/14
238) !
239)
240) implicit none
241)
242) PetscReal, intent(in) :: dist(-1:3)
243) PetscReal, intent(in) :: gravity(3)
244)
245) PetscReal, intent(out) :: distance_upwind
246) PetscReal, intent(out) :: distance_downwind
247) PetscReal, intent(out) :: distance_gravity
248) PetscReal, intent(out) :: upwind_weight
249)
250) ! dist(-1) = scalar - fraction upwind
251) ! dist(0) = scalar - magnitude of distance
252) ! gravity = vector(3)
253) ! dist(1:3) = vector(3) - unit vector
254) distance_gravity = dist(0) * & ! distance_gravity = dx*g*n
255) dot_product(gravity,dist(1:3))
256) distance_upwind = dist(0)*dist(-1)
257) distance_downwind = dist(0)-distance_upwind ! should avoid truncation error
258) ! upweight could be calculated as 1.d0-fraction_upwind
259) ! however, this introduces ever so slight error causing pflow-overhaul not
260) ! to match pflow-orig. This can be changed to 1.d0-fraction_upwind
261) upwind_weight = distance_downwind/(distance_upwind+distance_downwind)
262)
263) end subroutine ConnectionCalculateDistances
264)
265) ! ************************************************************************** !
266)
267) subroutine ConnectionDestroy(connection)
268) !
269) ! Deallocates a connection
270) !
271) ! Author: Glenn Hammond
272) ! Date: 10/23/07
273) !
274) use Utility_module, only : DeallocateArray
275)
276) implicit none
277)
278) type(connection_set_type), pointer :: connection
279)
280) if (.not.associated(connection)) return
281)
282) call DeallocateArray(connection%local)
283) call DeallocateArray(connection%id_up)
284) call DeallocateArray(connection%id_dn)
285) call DeallocateArray(connection%id_up2)
286) call DeallocateArray(connection%id_dn2)
287) call DeallocateArray(connection%face_id)
288) call DeallocateArray(connection%dist)
289) call DeallocateArray(connection%intercp)
290) call DeallocateArray(connection%area)
291) call DeallocateArray(connection%cntr)
292)
293) nullify(connection%next)
294)
295) deallocate(connection)
296) nullify(connection)
297)
298) end subroutine ConnectionDestroy
299)
300) ! ************************************************************************** !
301)
302) subroutine ConnectionDestroyList(list)
303) !
304) ! Deallocates the module global list and array of regions
305) !
306) ! Author: Glenn Hammond
307) ! Date: 10/15/07
308) !
309)
310) implicit none
311)
312) type(connection_set_list_type), pointer :: list
313)
314) type(connection_set_type), pointer :: cur_connection_set, prev_connection_set
315)
316) if (.not.associated(list)) return
317)
318) if (associated(list%array)) deallocate(list%array)
319) nullify(list%array)
320)
321) cur_connection_set => list%first
322) do
323) if (.not.associated(cur_connection_set)) exit
324) prev_connection_set => cur_connection_set
325) cur_connection_set => cur_connection_set%next
326) call ConnectionDestroy(prev_connection_set)
327) enddo
328)
329) nullify(list%first)
330) nullify(list%last)
331) list%num_connection_objects = 0
332)
333) deallocate(list)
334) nullify(list)
335)
336) end subroutine ConnectionDestroyList
337)
338) end module Connection_module