dataset_base.F90 coverage: 84.62 %func 73.96 %block
1) module Dataset_Base_class
2)
3) use Time_Storage_module
4)
5) use PFLOTRAN_Constants_module
6)
7) implicit none
8)
9) private
10)
11) #include "petsc/finclude/petscsys.h"
12)
13) type, public :: dataset_base_type
14) character(len=MAXWORDLENGTH) :: name
15) character(len=MAXSTRINGLENGTH) :: filename
16) type(time_storage_type), pointer :: time_storage ! stores transient times
17) character(len=MAXSTRINGLENGTH) :: header
18) PetscInt :: rank ! size of dims(:)
19) PetscInt, pointer :: dims(:) ! dimensions of arrays (excludes time)
20) PetscInt :: data_type
21) PetscInt, pointer :: iarray(:)
22) PetscInt, pointer :: ibuffer(:)
23) PetscReal, pointer :: rarray(:)
24) PetscReal, pointer :: rbuffer(:)
25) PetscInt :: buffer_slice_offset ! index of the first time slice in the buffer
26) PetscInt :: buffer_nslice ! # of time slices stored in buffer
27) class(dataset_base_type), pointer :: next
28) end type dataset_base_type
29)
30) ! dataset types
31) PetscInt, public, parameter :: DATASET_INTEGER = 1
32) PetscInt, public, parameter :: DATASET_REAL = 2
33)
34) public :: DatasetBaseCreate, &
35) DatasetBaseInit, &
36) DatasetBaseCopy, &
37) DatasetBaseVerify, &
38) DatasetBaseInterpolateTime, &
39) DatasetBaseReorder, &
40) DatasetBaseGetPointer, &
41) DatasetBaseAddToList, &
42) DatasetBasePrint, &
43) DatasetBaseStrip, &
44) DatasetBaseDestroy
45) contains
46)
47) ! ************************************************************************** !
48)
49) function DatasetBaseCreate()
50) !
51) ! Creates members of base database class
52) !
53) ! Author: Glenn Hammond
54) ! Date: 05/03/13
55) !
56)
57) implicit none
58)
59) class(dataset_base_type), pointer :: dataset
60)
61) class(dataset_base_type), pointer :: DatasetBaseCreate
62)
63) allocate(dataset)
64) call DatasetBaseInit(dataset)
65)
66) DatasetBaseCreate => dataset
67)
68) end function DatasetBaseCreate
69)
70) ! ************************************************************************** !
71)
72) subroutine DatasetBaseInit(this)
73) !
74) ! Initializes members of base database class
75) !
76) ! Author: Glenn Hammond
77) ! Date: 05/03/13
78) !
79)
80) implicit none
81)
82) class(dataset_base_type) :: this
83)
84) this%name = ''
85) this%filename = ''
86) this%header = ''
87) this%rank = 0
88) this%data_type = 0
89) nullify(this%time_storage)
90) nullify(this%dims)
91) nullify(this%iarray)
92) nullify(this%ibuffer)
93) nullify(this%rarray)
94) nullify(this%rbuffer)
95) this%buffer_slice_offset = 0
96) this%buffer_nslice = 0
97) nullify(this%next)
98)
99) end subroutine DatasetBaseInit
100)
101) ! ************************************************************************** !
102)
103) subroutine DatasetBaseCopy(this, that)
104) !
105) ! Copies members of base database class
106) !
107) ! Author: Glenn Hammond
108) ! Date: 05/03/13
109) !
110)
111) implicit none
112)
113) class(dataset_base_type) :: this
114) class(dataset_base_type) :: that
115)
116) that%name = this%name
117) that%filename = this%filename
118) that%rank = this%rank
119) that%data_type = this%data_type
120) if (associated(this%time_storage)) then
121) that%time_storage => this%time_storage
122) nullify(this%time_storage)
123) endif
124) if (associated(this%dims)) then
125) that%dims => this%dims
126) nullify(this%dims)
127) endif
128) if (associated(this%iarray)) then
129) that%iarray => this%iarray
130) nullify(this%iarray)
131) endif
132) if (associated(this%ibuffer)) then
133) that%ibuffer => this%ibuffer
134) nullify(this%ibuffer)
135) endif
136) if (associated(this%rarray)) then
137) that%rarray => this%rarray
138) nullify(this%rarray)
139) endif
140) if (associated(this%rbuffer)) then
141) that%rbuffer => this%rbuffer
142) nullify(this%rbuffer)
143) endif
144) that%buffer_slice_offset = this%buffer_slice_offset
145) that%buffer_nslice = this%buffer_nslice
146) that%next => this%next
147) nullify(this%next)
148)
149) end subroutine DatasetBaseCopy
150)
151) ! ************************************************************************** !
152)
153) subroutine DatasetBaseVerify(this,option)
154) !
155) ! Verifies that data structure is properly set up.
156) !
157) ! Author: Glenn Hammond
158) ! Date: 10/08/13
159) !
160)
161) use Option_module
162)
163) implicit none
164)
165) class(dataset_base_type) :: this
166) type(option_type) :: option
167)
168) if (associated(this%time_storage) .or. associated(this%iarray) .or. &
169) associated(this%iarray) .or. associated(this%ibuffer) .or. &
170) associated(this%rarray) .or. associated(this%rbuffer)) then
171) if (len_trim(this%name) < 1) then
172) this%name = 'Unknown Dataset'
173) endif
174) if (this%data_type == 0) then
175) option%io_buffer = '"data_type" not defined in dataset: ' // &
176) trim(this%name)
177) call printErrMsg(option)
178) endif
179) if (associated(this%ibuffer) .or. associated(this%rbuffer)) then
180) if (.not.associated(this%dims)) then
181) option%io_buffer = '"dims" not allocated in dataset: ' // &
182) trim(this%name)
183) call printErrMsg(option)
184) endif
185) if (this%dims(1) == 0) then
186) option%io_buffer = '"dims" not defined in dataset: ' // &
187) trim(this%name)
188) call printErrMsg(option)
189) endif
190) if (size(this%dims) /= this%rank) then
191) option%io_buffer = 'Size of "dims" not match "rank" in dataset: ' // &
192) trim(this%name)
193) call printErrMsg(option)
194) endif
195) endif
196) if (associated(this%ibuffer) .and. .not.associated(this%iarray)) then
197) allocate(this%iarray(this%dims(1)))
198) this%iarray = 0
199) endif
200) if (associated(this%rbuffer) .and. .not.associated(this%rarray)) then
201) allocate(this%rarray(this%dims(1)))
202) this%rarray = 0.d0
203) endif
204) else if (len_trim(this%name) < 1) then
205) option%io_buffer = 'NULL dataset in input deck.'
206) call printErrMsg(option)
207) endif
208)
209) end subroutine DatasetBaseVerify
210)
211) ! ************************************************************************** !
212)
213) subroutine DatasetBaseInterpolateTime(this)
214) !
215) ! Interpolates dataset between two buffer times
216) !
217) ! Author: Glenn Hammond
218) ! Date: 10/26/11
219) !
220)
221) use Option_module
222)
223) implicit none
224)
225) class(dataset_base_type) :: this
226) type(option_type) :: option
227)
228) PetscInt :: array_size,i
229) PetscInt :: time_interpolation_method
230) PetscReal :: weight2
231) PetscInt :: time1_start, time1_end, time2_start, time2_end
232)
233) if (.not.associated(this%rbuffer)) return
234)
235) time_interpolation_method = this%time_storage%time_interpolation_method
236)
237) if (this%time_storage%cur_time_index >= &
238) this%time_storage%max_time_index) then
239) ! dataset has reached the end of the time array and is not cyclic.
240) ! no interpolation needed
241) if (this%time_storage%cur_time_index_changed) then
242) array_size = size(this%rarray)
243) time1_end = array_size*(this%time_storage%max_time_index - &
244) this%buffer_slice_offset)
245) time1_start = time1_end - array_size + 1
246) this%rarray = this%rbuffer(time1_start:time1_end)
247) endif
248) return
249) endif
250)
251) array_size = size(this%rarray)
252) select case(time_interpolation_method)
253) case(INTERPOLATION_NULL)
254) ! it possible that the time_interpolation_method is null during
255) ! the load process. in that case, we need to set the array to an
256) ! uninitialized value (UNINITIALIZED_DOUBLE).
257) this%rarray = UNINITIALIZED_DOUBLE
258) case(INTERPOLATION_STEP)
259) ! if time index has not changed skip
260) if (.not.this%time_storage%cur_time_index_changed) return
261) time1_end = array_size*(this%time_storage%cur_time_index - &
262) this%buffer_slice_offset)
263) time1_start = time1_end - array_size + 1
264) this%rarray = this%rbuffer(time1_start:time1_end)
265) case(INTERPOLATION_LINEAR)
266) ! if fraction has not changed skip
267) if (.not.this%time_storage%cur_time_fraction_changed) return
268) time1_end = array_size*(this%time_storage%cur_time_index - &
269) this%buffer_slice_offset)
270) time1_start = time1_end - array_size + 1
271) time2_end = time1_end + array_size
272) time2_start = time1_start + array_size
273) this%rarray = (1.d0-this%time_storage%cur_time_fraction) * &
274) this%rbuffer(time1_start:time1_end) + &
275) this%time_storage%cur_time_fraction * &
276) this%rbuffer(time2_start:time2_end)
277) end select
278)
279) end subroutine DatasetBaseInterpolateTime
280)
281) ! ************************************************************************** !
282)
283) subroutine DatasetBaseInterpolateSpace(this,xx,yy,zz,time,real_value,option)
284) !
285) ! Interpolates data from the dataset
286) !
287) ! Author: Glenn Hammond
288) ! Date: 10/26/11
289) !
290)
291) use Utility_module, only : InterpolateBilinear
292) use Option_module
293)
294) implicit none
295)
296) class(dataset_base_type) :: this
297) PetscReal, intent(in) :: xx, yy, zz
298) PetscReal :: time
299) PetscReal :: real_value
300) type(option_type) :: option
301)
302) end subroutine DatasetBaseInterpolateSpace
303)
304) ! ************************************************************************** !
305)
306) subroutine DatasetBaseReorder(this,option)
307) !
308) ! If a dataset is loaded from an HDF5 file, and it was
309) ! multidimensional in the HDF5 file, the array needs to be
310) ! reordered fro Fortran -> C indexing. This subroutine
311) ! takes care of the reordering.
312) !
313) ! Author: Glenn Hammond
314) ! Date: 10/26/11
315) !
316)
317) use Option_module
318)
319) implicit none
320)
321) class(dataset_base_type) :: this
322) type(option_type) :: option
323)
324) PetscReal, allocatable :: temp_real(:)
325) PetscInt :: i, j, k, l
326) PetscInt :: dims(4), n1, n1Xn2, n1Xn2Xn3
327) PetscInt :: count, index
328) PetscReal, pointer :: rarray(:)
329)
330) if (this%data_type == DATASET_INTEGER) then
331) option%io_buffer = 'Reordering of integer data sets not yet supported.'
332) call printErrMsg(option)
333) endif
334)
335) ! set each dim to 1 by default for loop below
336) dims(:) = 1
337) dims(1:this%rank) = this%dims(1:this%rank)
338) if (associated(this%rbuffer)) then
339) rarray => this%rbuffer
340) dims(this%rank+1) = this%buffer_nslice
341) else
342) rarray => this%rarray
343) endif
344)
345) ! Not necessary for 1D arrays
346) if (maxval(dims(2:)) == 1) return
347)
348) l = 1
349) do i = 1, size(dims)
350) l = l*dims(i)
351) enddo
352) allocate(temp_real(l))
353)
354) n1 = dims(1)
355) n1Xn2 = n1*dims(2)
356) n1Xn2Xn3 = n1Xn2*dims(3)
357) count = 0
358) do i = 1, dims(1)
359) do j = 0, dims(2)-1
360) do k = 0, dims(3)-1
361) do l = 0, dims(4)-1
362) index = l*n1Xn2Xn3+k*n1Xn2+j*n1+i
363) count = count+1
364) temp_real(index) = rarray(count)
365) enddo
366) enddo
367) enddo
368) enddo
369)
370) rarray = temp_real
371) deallocate(temp_real)
372)
373) end subroutine DatasetBaseReorder
374)
375) ! ************************************************************************** !
376)
377) subroutine DatasetBaseGetTimes(this, option, max_sim_time, time_array)
378) !
379) ! Fills an array of times based on a dataset
380) !
381) ! Author: Glenn Hammond
382) ! Date: 10/26/11
383) !
384)
385) use Option_module
386)
387) implicit none
388)
389) class(dataset_base_type) :: this
390) type(option_type) :: option
391) PetscReal :: max_sim_time
392) PetscReal, pointer :: time_array(:)
393)
394)
395) if (associated(this%time_storage)) then
396) call TimeStorageGetTimes(this%time_storage, option, max_sim_time, &
397) time_array)
398) endif
399)
400) end subroutine DatasetBaseGetTimes
401)
402) ! ************************************************************************** !
403)
404) subroutine DatasetBaseAddToList(dataset,list)
405) !
406) ! Adds a dataset to linked list
407) !
408) ! Author: Glenn Hammond
409) ! Date: 01/12/11
410) !
411)
412) implicit none
413)
414) class(dataset_base_type), pointer :: dataset
415) class(dataset_base_type), pointer :: list
416)
417) class(dataset_base_type), pointer :: cur_dataset
418)
419) if (associated(list)) then
420) cur_dataset => list
421) ! loop to end of list
422) do
423) if (.not.associated(cur_dataset%next)) exit
424) cur_dataset => cur_dataset%next
425) enddo
426) cur_dataset%next => dataset
427) else
428) list => dataset
429) endif
430)
431) end subroutine DatasetBaseAddToList
432)
433) ! ************************************************************************** !
434)
435) function DatasetBaseGetPointer(dataset_list, dataset_name, debug_string, &
436) option)
437) !
438) ! Returns the pointer to the dataset named "name"
439) !
440) ! Author: Glenn Hammond
441) ! Date: 01/12/11
442) !
443)
444) use Option_module
445) use String_module
446)
447) class(dataset_base_type), pointer :: dataset_list
448) character(len=MAXWORDLENGTH) :: dataset_name
449) character(len=MAXSTRINGLENGTH) :: debug_string
450) type(option_type) :: option
451)
452) class(dataset_base_type), pointer :: DatasetBaseGetPointer
453) PetscBool :: found
454) class(dataset_base_type), pointer :: cur_dataset
455)
456) found = PETSC_FALSE
457) cur_dataset => dataset_list
458) do
459) if (.not.associated(cur_dataset)) exit
460) if (StringCompare(dataset_name, &
461) cur_dataset%name,MAXWORDLENGTH)) then
462) found = PETSC_TRUE
463) DatasetBaseGetPointer => cur_dataset
464) return
465) endif
466) cur_dataset => cur_dataset%next
467) enddo
468) if (.not.found) then
469) option%io_buffer = 'Dataset "' // trim(dataset_name) // '" in "' // &
470) trim(debug_string) // '" not found among available datasets.'
471) call printErrMsgByRank(option)
472) endif
473)
474) end function DatasetBaseGetPointer
475)
476) ! ************************************************************************** !
477)
478) subroutine DatasetBasePrint(this,option)
479) !
480) ! Prints dataset info
481) !
482) ! Author: Glenn Hammond
483) ! Date: 10/22/13
484) !
485)
486) use Option_module
487)
488) implicit none
489)
490) class(dataset_base_type) :: this
491) type(option_type) :: option
492)
493) if (len_trim(this%filename) > 0) then
494) write(option%fid_out,'(10x,''Filename: '',a)') trim(this%filename)
495) endif
496) if (associated(this%time_storage)) then
497) write(option%fid_out,'(10x,''Is transient?: yes'')')
498) write(option%fid_out,'(10x,''Number of times: '',i6)') &
499) this%time_storage%max_time_index
500) if (this%time_storage%is_cyclic) then
501) write(option%fid_out,'(10x,''Is cyclic?: yes'')')
502) else
503) write(option%fid_out,'(10x,''Is cyclic?: no'')')
504) endif
505) else
506) write(option%fid_out,'(10x,''Transient: no'')')
507) endif
508) if (associated(this%rbuffer)) then
509) write(option%fid_out,'(10x,''Buffer:'')')
510) write(option%fid_out,'(12x,''Rank: '',i2)') this%rank
511) if (associated(this%dims)) then
512) write(option%fid_out,'(12x,''Dims: '',10i4)') this%dims
513) endif
514) write(option%fid_out,'(12x,''Buffer Slice Size: '',i3)') this%buffer_nslice
515) endif
516)
517) end subroutine DatasetBasePrint
518)
519) ! ************************************************************************** !
520)
521) subroutine DatasetBaseStrip(this)
522) !
523) ! Strips allocated objects within base dataset object
524) !
525) ! Author: Glenn Hammond
526) ! Date: 05/03/13
527) !
528)
529) use Utility_module, only : DeallocateArray
530)
531) implicit none
532)
533) class(dataset_base_type) :: this
534)
535) call DeallocateArray(this%iarray)
536) call DeallocateArray(this%rarray)
537) call DeallocateArray(this%ibuffer)
538) call DeallocateArray(this%rbuffer)
539) call DeallocateArray(this%dims)
540)
541) call TimeStorageDestroy(this%time_storage)
542)
543) end subroutine DatasetBaseStrip
544)
545) ! ************************************************************************** !
546)
547) subroutine DatasetBaseDestroy(dataset)
548) !
549) ! Destroys a dataset
550) !
551) ! Author: Glenn Hammond
552) ! Date: 01/12/11, 05/03/13
553) !
554)
555) implicit none
556)
557) class(dataset_base_type), pointer :: dataset
558)
559) if (.not.associated(dataset)) return
560)
561) call DatasetBaseStrip(dataset)
562)
563) deallocate(dataset)
564) nullify(dataset)
565)
566) end subroutine DatasetBaseDestroy
567)
568) end module Dataset_Base_class