dataset.F90 coverage: 100.00 %func 82.72 %block
1) module Dataset_module
2)
3) use Dataset_Base_class
4) use Dataset_Ascii_class
5) use Dataset_Common_HDF5_class
6) use Dataset_Gridded_HDF5_class
7) use Dataset_Map_HDF5_class
8) use Dataset_Global_HDF5_class
9)
10) use PFLOTRAN_Constants_module
11)
12) implicit none
13)
14) private
15)
16) #include "petsc/finclude/petscsys.h"
17)
18) public :: DatasetRead, &
19) DatasetScreenForNonCellIndexed, &
20) DatasetLoad, &
21) DatasetVerify, &
22) DatasetUpdate, &
23) DatasetFindInList, &
24) DatasetIsTransient, &
25) DatasetGetClass, &
26) DatasetPrint, &
27) DatasetReadDoubleOrDataset, &
28) DatasetDestroy
29)
30) contains
31)
32) ! ************************************************************************** !
33)
34) subroutine DatasetRead(input,dataset,option)
35) !
36) ! Reads a dataset from the input file
37) !
38) ! Author: Glenn Hammond
39) ! Date: 03/26/12
40) !
41)
42) use Input_Aux_module
43) use Option_module
44) use String_module
45)
46) implicit none
47)
48) type(input_type), pointer :: input
49) class(dataset_base_type), pointer :: dataset
50) class(dataset_map_hdf5_type), pointer :: dataset_map_hdf5
51) class(dataset_global_hdf5_type), pointer :: dataset_global_hdf5
52) type(option_type) :: option
53)
54) character(len=MAXWORDLENGTH) :: word, word2
55)
56) ! read from the buffer on the first line the type of dataset. If it does
57) ! not match any of type in the select case, assume default and use the
58) ! word as the name.
59)
60) call InputReadWord(input,option,word,PETSC_TRUE)
61) word2 = word
62) call StringToUpper(word2)
63)
64) select case(word2)
65) case('MAPPED')
66) dataset_map_hdf5 => DatasetMapHDF5Create()
67) call InputReadWord(input,option,dataset_map_hdf5%name,PETSC_TRUE)
68) call InputDefaultMsg(input,option,'DATASET name')
69) call DatasetMapHDF5Read(dataset_map_hdf5,input,option)
70) dataset => dataset_map_hdf5
71) case('GLOBAL')
72) dataset_global_hdf5 => DatasetGlobalHDF5Create()
73) call InputReadWord(input,option,dataset_global_hdf5%name,PETSC_TRUE)
74) call InputDefaultMsg(input,option,'DATASET name')
75) call DatasetCommonHDF5Read(dataset_global_hdf5,input,option)
76) dataset => dataset_global_hdf5
77) case default ! CELL_INDEXED, GLOBAL, GRIDDED
78) dataset => DatasetCommonHDF5Create()
79) select case(word2)
80) case('CELL_INDEXED', 'GRIDDED')
81) call InputReadWord(input,option,dataset%name,PETSC_TRUE)
82) case default
83) if (.not.InputError(input)) then
84) dataset%name = trim(word)
85) endif
86) end select
87) call InputDefaultMsg(input,option,'DATASET name')
88) call DatasetCommonHDF5Read(DatasetCommonHDF5Cast(dataset),input, &
89) option)
90) end select
91)
92) end subroutine DatasetRead
93)
94) ! ************************************************************************** !
95)
96) subroutine DatasetScreenForNonCellIndexed(datasets,option)
97) !
98) ! Recasts datasets of dataset_common_hdf5_type
99) ! to dataset_gridded_hdf5_type
100) !
101) ! Author: Glenn Hammond
102) ! Date: 03/26/12
103) !
104)
105) use Option_module
106)
107) implicit none
108)
109) class(dataset_base_type), pointer :: datasets
110) type(option_type) :: option
111)
112) class(dataset_base_type), pointer :: cur_dataset
113) class(dataset_base_type), pointer :: prev_dataset
114) class(dataset_gridded_hdf5_type), pointer :: dataset_gridded_hdf5
115) PetscBool :: swapped
116)
117) cur_dataset => datasets
118) nullify(prev_dataset)
119) do
120) if (.not.associated(cur_dataset)) exit
121) nullify(dataset_gridded_hdf5)
122) swapped = PETSC_FALSE
123) select type(cur_dataset)
124) class is(dataset_map_hdf5_type)
125) ! do nothing
126) class is(dataset_global_hdf5_type)
127) ! do nothing
128) class is(dataset_common_hdf5_type)
129) ! if class is dataset_common_hdf5_type, the dataset should really
130) ! be dataset_gridded_hdf5 unless cell indexed
131) cur_dataset%is_cell_indexed = &
132) DatasetCommonHDF5IsCellIndexed(cur_dataset,option)
133) ! if not cell indexed, we need to change the type to the extended
134) ! dataset_gridded_hdf5 type
135) if (.not.cur_dataset%is_cell_indexed) then
136) swapped = PETSC_TRUE
137) dataset_gridded_hdf5 => DatasetGriddedHDF5Create()
138) call DatasetCommonHDF5Copy(cur_dataset,dataset_gridded_hdf5)
139) endif
140) class default
141) option%io_buffer = &
142) 'Unknown dataset type in DatasetScreenForNonCellIndexed.'
143) call printErrMsg(option)
144) end select
145) ! if we changed the derived type, we need to replace the old dataset
146) ! in the linked list of datasets and destroy the old dataset.
147) if (swapped) then
148) ! dataset_gridded_hdf5%next is already set to cur_dataset%next
149) if (associated(prev_dataset)) then
150) prev_dataset%next => dataset_gridded_hdf5
151) else
152) datasets => dataset_gridded_hdf5
153) endif
154) ! just to be sure
155) nullify(cur_dataset%next)
156) call DatasetDestroy(cur_dataset)
157) cur_dataset => dataset_gridded_hdf5
158) endif
159) prev_dataset => cur_dataset
160) cur_dataset => cur_dataset%next
161) enddo
162)
163) end subroutine DatasetScreenForNonCellIndexed
164)
165) ! ************************************************************************** !
166)
167) subroutine DatasetVerify(dataset,default_time_storage,option)
168) !
169) ! Verifies that a dataset is intact and useable.
170) !
171) ! Author: Glenn Hammond
172) ! Date: 10/08/13
173) !
174)
175) use Time_Storage_module
176) use Option_module
177)
178) implicit none
179)
180) class(dataset_base_type), pointer :: dataset
181) type(time_storage_type), pointer :: default_time_storage
182) type(option_type) :: option
183)
184) type(time_storage_type), pointer :: default
185)
186) if (.not.associated(dataset)) return
187)
188) call TimeStorageVerify(0.d0,dataset%time_storage,default_time_storage,option)
189) select type(dataset_ptr => dataset)
190) class is(dataset_ascii_type)
191) call DatasetAsciiVerify(dataset_ptr,option)
192) class is(dataset_base_type)
193) call DatasetBaseVerify(dataset_ptr,option)
194) class default
195) option%io_buffer = 'DatasetXXXVerify needed for unknown dataset type'
196) call printErrMsg(option)
197) end select
198)
199) end subroutine DatasetVerify
200)
201) ! ************************************************************************** !
202)
203) recursive subroutine DatasetUpdate(dataset,time,option)
204) !
205) ! Updates a dataset based on type
206) !
207) ! Author: Glenn Hammond
208) ! Date: 10/09/13
209) !
210)
211) use Option_module
212)
213) implicit none
214)
215) PetscReal :: time
216) class(dataset_base_type), pointer :: dataset
217) type(option_type) :: option
218)
219) class(dataset_ascii_type), pointer :: dataset_ascii
220)
221) if (.not.associated(dataset)) return
222)
223) if (associated(dataset%time_storage)) then
224) dataset%time_storage%cur_time = time
225) endif
226) select type (selector => dataset)
227) class is (dataset_ascii_type)
228) dataset_ascii => selector
229) call DatasetAsciiUpdate(dataset_ascii,option)
230) class is (dataset_common_hdf5_type)
231) if (associated(selector%time_storage) .and. &
232) .not.selector%is_cell_indexed) then
233) call DatasetLoad(dataset,option)
234) endif
235) end select
236)
237) end subroutine DatasetUpdate
238)
239) ! ************************************************************************** !
240)
241) recursive subroutine DatasetLoad(dataset,option)
242) !
243) ! Loads a dataset based on type
244) !
245) ! Author: Glenn Hammond
246) ! Date: 06/03/13
247) !
248)
249) use DM_Kludge_module
250) use Option_module
251)
252) implicit none
253)
254) class(dataset_base_type), pointer :: dataset
255) type(option_type) :: option
256)
257) class(dataset_global_hdf5_type), pointer :: dataset_global_hdf5
258) class(dataset_gridded_hdf5_type), pointer :: dataset_gridded_hdf5
259) class(dataset_map_hdf5_type), pointer :: dataset_map_hdf5
260) class(dataset_common_hdf5_type), pointer :: dataset_common_hdf5
261) class(dataset_base_type), pointer :: dataset_base
262)
263) select type (selector => dataset)
264) class is (dataset_global_hdf5_type)
265) dataset_global_hdf5 => selector
266) call DatasetGlobalHDF5Load(dataset_global_hdf5,option)
267) class is (dataset_gridded_hdf5_type)
268) dataset_gridded_hdf5 => selector
269) call DatasetGriddedHDF5Load(dataset_gridded_hdf5,option)
270) class is (dataset_map_hdf5_type)
271) dataset_map_hdf5 => selector
272) call DatasetMapHDF5Load(dataset_map_hdf5,option)
273) class is (dataset_common_hdf5_type)
274) dataset_common_hdf5 => selector
275) if (dataset_common_hdf5%is_cell_indexed) then
276) option%io_buffer = 'Cell Indexed datasets opened later.'
277) call printMsg(option)
278) else
279) option%io_buffer = 'Unrecognized dataset that extends ' // &
280) 'dataset_common_hdf5 in DatasetLoad.'
281) call printErrMsg(option)
282) endif
283) class is (dataset_ascii_type)
284) ! do nothing. dataset should already be in memory at this point
285) class is (dataset_base_type)
286) dataset_base => selector
287) option%io_buffer = 'DatasetLoad not yet supported for base dataset class.'
288) call printErrMsg(option)
289) end select
290)
291) end subroutine DatasetLoad
292)
293) ! ************************************************************************** !
294)
295) subroutine DatasetFindInList(list,dataset_base,default_time_storage, &
296) error_string,option)
297) !
298) ! Uses a dummy dataset with name to find the actual
299) ! dataset in a list of datasets
300) !
301) ! Author: Glenn Hammond
302) ! Date: 10/07/13
303) !
304)
305) use Dataset_Base_class
306) use Dataset_Ascii_class
307) use Time_Storage_module
308) use Option_module
309)
310) implicit none
311)
312) class(dataset_base_type), pointer :: list
313) class(dataset_base_type), pointer :: dataset_base
314) type(time_storage_type), pointer :: default_time_storage
315) character(len=MAXSTRINGLENGTH) :: error_string
316) type(option_type) :: option
317)
318) character(len=MAXWORDLENGTH) :: dataset_name
319) PetscReal, parameter :: time = 0.d0
320)
321) ! check for dataset in flow_dataset
322) if (associated(dataset_base)) then
323) select type(dataset => dataset_base)
324) class is(dataset_ascii_type)
325) ! do nothing as the correct dataset already exists
326) class default
327) dataset_name = dataset%name
328) ! delete the dataset since it is solely a placeholder
329) call DatasetDestroy(dataset_base)
330) ! get dataset from list
331) dataset_base => &
332) DatasetBaseGetPointer(list,dataset_name,error_string,option)
333) ! once a dataset is linked to the dataset list, it needs to be loaded
334) ! immediately
335) call DatasetLoad(dataset_base,option)
336) call DatasetVerify(dataset_base,default_time_storage,option)
337) ! must update after DatasetVerify since the time interpolation method
338) ! may not have been properly set during the load! Force the update.
339) if (associated(dataset_base%time_storage)) then
340) dataset_base%time_storage%force_update = PETSC_TRUE
341) endif
342) call DatasetUpdate(dataset_base,time,option)
343) end select
344) endif
345)
346) end subroutine DatasetFindInList
347)
348) ! ************************************************************************** !
349)
350) function DatasetIsTransient(dataset)
351) !
352) ! Determines whether a dataset is transient
353) !
354) ! Author: Glenn Hammond
355) ! Date: 10/10/13
356) !
357)
358) use Time_Storage_module
359)
360) implicit none
361)
362) class(dataset_base_type), pointer :: dataset
363)
364) PetscBool :: DatasetIsTransient
365)
366) DatasetIsTransient = PETSC_FALSE
367) if (associated(dataset)) then
368) if (associated(dataset%time_storage)) then
369) DatasetIsTransient = PETSC_TRUE
370) endif
371) endif
372)
373) end function DatasetIsTransient
374)
375) ! ************************************************************************** !
376)
377) subroutine DatasetPrint(this,option)
378) !
379) ! Prints dataset info
380) !
381) ! Author: Glenn Hammond
382) ! Date: 10/22/13
383) !
384)
385) use Option_module
386)
387) implicit none
388)
389) class(dataset_base_type) :: this
390) type(option_type) :: option
391)
392) write(option%fid_out,'(8x,''Dataset: '',a)') trim(this%name)
393) write(option%fid_out,'(10x,''Type: '',a)') trim(DatasetGetClass(this))
394)
395) call DatasetBasePrint(this,option)
396)
397) select type (d=>this)
398) class is (dataset_ascii_type)
399) call DatasetAsciiPrint(d,option)
400) class is (dataset_global_hdf5_type)
401) call DatasetGlobalHDF5Print(d,option)
402) class is (dataset_gridded_hdf5_type)
403) call DatasetGriddedHDF5Print(d,option)
404) class is (dataset_map_hdf5_type)
405) call DatasetMapHDF5Print(d,option)
406) class is (dataset_common_hdf5_type)
407) call DatasetCommonHDF5Print(d,option)
408) class default
409) option%io_buffer = 'Unknown dataset type for dataset "' // &
410) trim(this%name) // '" in DatasetPrint()'
411) end select
412)
413) end subroutine DatasetPrint
414)
415) ! ************************************************************************** !
416)
417) function DatasetGetClass(dataset)
418) !
419) ! Returns a string defining the class of dataset
420) !
421) ! Author: Glenn Hammond
422) ! Date: 10/10/13
423) !
424)
425) use Option_module
426)
427) implicit none
428)
429) class(dataset_base_type) :: dataset
430)
431) character(len=MAXSTRINGLENGTH) :: DatasetGetClass
432)
433) select type (dataset)
434) class is (dataset_ascii_type)
435) DatasetGetClass = 'dataset_ascii_type'
436) class is (dataset_global_hdf5_type)
437) DatasetGetClass = 'dataset_global_hdf5_type'
438) class is (dataset_gridded_hdf5_type)
439) DatasetGetClass = 'dataset_gridded_hdf5_type'
440) class is (dataset_map_hdf5_type)
441) DatasetGetClass = 'dataset_map_hdf5_type'
442) class is (dataset_common_hdf5_type)
443) DatasetGetClass = 'dataset_hdf5_type'
444) class is (dataset_base_type)
445) DatasetGetClass = 'dataset_base_type'
446) class default
447) DatasetGetClass = 'dataset class unknown'
448) end select
449)
450) end function DatasetGetClass
451)
452) ! ************************************************************************** !
453)
454) subroutine DatasetReadDoubleOrDataset(input,double_value,dataset, &
455) error_variable,error_card,option)
456) !
457) ! Prints dataset info
458) !
459) ! Author: Glenn Hammond
460) ! Date: 10/22/13
461) !
462)
463) use Option_module
464) use Input_Aux_module
465) use String_module
466)
467) implicit none
468)
469) type(input_type), pointer :: input
470) PetscReal :: double_value
471) class(dataset_base_type), pointer :: dataset
472) character(len=*) :: error_variable
473) character(len=*) :: error_card
474) type(option_type) :: option
475)
476) character(len=MAXSTRINGLENGTH) :: string
477) character(len=MAXSTRINGLENGTH) :: buffer_save
478)
479) buffer_save = input%buf
480) call InputReadNChars(input,option,string,MAXSTRINGLENGTH,PETSC_TRUE)
481) call InputErrorMsg(input,option,error_variable,error_card)
482) call StringToUpper(string)
483) if (StringCompare(string,'DATASET',SEVEN_INTEGER)) then
484) dataset => DatasetBaseCreate()
485) call InputReadWord(input,option,dataset%name,PETSC_TRUE)
486) string = trim(error_card) // ',' // trim(error_variable)
487) call InputErrorMsg(input,option,'DATASET,NAME',string)
488) else
489) input%buf = buffer_save
490) call InputReadDouble(input,option,double_value)
491) call InputErrorMsg(input,option,error_variable,error_card)
492) endif
493)
494) end subroutine DatasetReadDoubleOrDataset
495)
496) ! ************************************************************************** !
497)
498) recursive subroutine DatasetDestroy(dataset)
499) !
500) ! Destroys a dataset
501) !
502) ! Author: Glenn Hammond
503) ! Date: 01/12/11
504) !
505)
506) implicit none
507)
508) class(dataset_base_type), pointer :: dataset
509)
510) class(dataset_global_hdf5_type), pointer :: dataset_global_hdf5
511) class(dataset_gridded_hdf5_type), pointer :: dataset_gridded_hdf5
512) class(dataset_map_hdf5_type), pointer :: dataset_map_hdf5
513) class(dataset_common_hdf5_type), pointer :: dataset_common_hdf5
514) class(dataset_ascii_type), pointer :: dataset_ascii
515) class(dataset_base_type), pointer :: dataset_base
516)
517) if (.not.associated(dataset)) return
518)
519) if (associated(dataset%next)) then
520) call DatasetDestroy(dataset%next)
521) endif
522)
523) select type (selector => dataset)
524) class is (dataset_ascii_type)
525) dataset_ascii => selector
526) call DatasetAsciiDestroy(dataset_ascii)
527) class is (dataset_global_hdf5_type)
528) dataset_global_hdf5 => selector
529) call DatasetGlobalHDF5Destroy(dataset_global_hdf5)
530) class is (dataset_gridded_hdf5_type)
531) dataset_gridded_hdf5 => selector
532) call DatasetGriddedHDF5Destroy(dataset_gridded_hdf5)
533) class is (dataset_map_hdf5_type)
534) dataset_map_hdf5 => selector
535) call DatasetMapHDF5Destroy(dataset_map_hdf5)
536) class is (dataset_common_hdf5_type)
537) dataset_common_hdf5 => selector
538) call DatasetCommonHDF5Destroy(dataset_common_hdf5)
539) class is (dataset_base_type)
540) dataset_base => selector
541) call DatasetBaseDestroy(dataset_base)
542) end select
543)
544) end subroutine DatasetDestroy
545)
546) end module Dataset_module