data_mediator_dataset.F90 coverage: 100.00 %func 80.36 %block
1) module Data_Mediator_Dataset_class
2)
3) use PFLOTRAN_Constants_module
4) use Data_Mediator_Base_class
5) use Dataset_Global_HDF5_class
6)
7) implicit none
8)
9) private
10)
11) #include "petsc/finclude/petscsys.h"
12) #include "petsc/finclude/petscvec.h"
13) #include "petsc/finclude/petscvec.h90"
14)
15) type, public, extends(data_mediator_base_type) :: data_mediator_dataset_type
16) PetscInt :: idof
17) class(dataset_global_hdf5_type), pointer :: dataset
18) contains
19) procedure, public :: Update => DataMediatorDatasetUpdate
20) procedure, public :: Strip => DataMediatorDatasetStrip
21) end type data_mediator_dataset_type
22)
23) public :: DataMediatorDatasetCreate, &
24) DataMediatorDatasetRead, &
25) DataMediatorDatasetInit
26)
27) contains
28)
29) ! ************************************************************************** !
30)
31) function DataMediatorDatasetCreate()
32) !
33) ! Creates a data mediator object
34) !
35) ! Author: Glenn Hammond
36) ! Date: 05/01/13
37) !
38)
39) implicit none
40)
41) class(data_mediator_dataset_type), pointer :: DataMediatorDatasetCreate
42)
43) class(data_mediator_dataset_type), pointer :: data_mediator
44)
45) allocate(data_mediator)
46) call DataMediatorBaseCreate(data_mediator)
47) data_mediator%idof = 0
48) nullify(data_mediator%dataset)
49) DataMediatorDatasetCreate => data_mediator
50)
51) end function DataMediatorDatasetCreate
52)
53) ! ************************************************************************** !
54)
55) subroutine DataMediatorDatasetRead(data_mediator,input,option)
56) !
57) ! Reads in contents of a data mediator card
58) !
59) ! Author: Glenn Hammond
60) ! Date: 05/01/13
61) !
62)
63) use Option_module
64) use Input_Aux_module
65) use String_module
66)
67) implicit none
68)
69) class(data_mediator_dataset_type) :: data_mediator
70) type(input_type), pointer :: input
71) type(option_type) :: option
72)
73) character(len=MAXWORDLENGTH) :: keyword, word
74)
75) input%ierr = 0
76) do
77)
78) call InputReadPflotranString(input,option)
79)
80) if (InputCheckExit(input,option)) exit
81)
82) call InputReadWord(input,option,keyword,PETSC_TRUE)
83) call InputErrorMsg(input,option,'keyword','MASS_TRANSFER')
84) call StringToUpper(keyword)
85)
86) select case(trim(keyword))
87) case('IDOF')
88) call InputReadInt(input,option,data_mediator%idof)
89) call InputErrorMsg(input,option,'idof','MASS_TRANSFER')
90) case('DATASET')
91) data_mediator%dataset => DatasetGlobalHDF5Create()
92) call InputReadNChars(input,option, &
93) data_mediator%dataset%name,&
94) MAXWORDLENGTH,PETSC_TRUE)
95) call InputErrorMsg(input,option,'DATASET,NAME','MASS_TRANSFER')
96) case default
97) call InputKeywordUnrecognized(keyword,'MASS_TRANSFER',option)
98) end select
99)
100) enddo
101)
102) end subroutine DataMediatorDatasetRead
103)
104) ! ************************************************************************** !
105)
106) subroutine DataMediatorDatasetInit(data_mediator, discretization, &
107) available_datasets, option)
108) !
109) ! Initializes data mediator object opening dataset to
110) ! set up times, vectors, etc.
111) !
112) ! Author: Glenn Hammond
113) ! Date: 05/09/13
114) !
115)
116) use Discretization_module
117) use Dataset_Base_class
118) use Dataset_Common_HDF5_class
119) use Option_module
120)
121) implicit none
122)
123) #include "petsc/finclude/petscvec.h"
124) #include "petsc/finclude/petscvec.h90"
125)
126) class(data_mediator_dataset_type) :: data_mediator
127) type(discretization_type) :: discretization
128) class(dataset_base_type), pointer :: available_datasets
129) type(option_type) :: option
130)
131) class(dataset_base_type), pointer :: dataset_base_ptr
132) character(len=MAXSTRINGLENGTH) :: string
133) PetscErrorCode :: ierr
134)
135) if (.not.associated(data_mediator%dataset)) then
136) option%io_buffer = 'A "global" DATASET does not exist for ' // &
137) 'MASS_TRANSFER object "' // trim(data_mediator%name) // '".'
138) call printErrMsg(option)
139) endif
140)
141) string = 'Data Mediator ' // trim(data_mediator%name)
142) dataset_base_ptr => &
143) DatasetBaseGetPointer(available_datasets,data_mediator%dataset%name, &
144) string,option)
145) call DatasetGlobalHDF5Destroy(data_mediator%dataset)
146) select type(dataset => dataset_base_ptr)
147) class is(dataset_global_hdf5_type)
148) data_mediator%dataset => dataset
149) class default
150) option%io_buffer = 'DATASET ' // trim(dataset%name) // 'is not of ' // &
151) 'GLOBAL type, which is necessary for all MASS_TRANSFER objects.'
152) call printErrMsg(option)
153) end select
154) ! dm_wrapper is solely a pointer; it should not be allocated
155) data_mediator%dataset%dm_wrapper => discretization%dm_1dof
156) data_mediator%dataset%local_size = discretization%grid%nlmax
157) data_mediator%dataset%global_size = discretization%grid%nmax
158)
159) if (.not.associated(data_mediator%dataset%time_storage)) then
160) #if defined(PETSC_HAVE_HDF5)
161) call DatasetCommonHDF5ReadTimes(data_mediator%dataset%filename, &
162) data_mediator%dataset%hdf5_dataset_name, &
163) data_mediator%dataset%time_storage,option)
164) #endif
165) ! if time interpolation methods not set in hdf5 file, set to default of STEP
166) if (data_mediator%dataset%time_storage%time_interpolation_method == &
167) INTERPOLATION_NULL) then
168) data_mediator%dataset%time_storage%time_interpolation_method = &
169) INTERPOLATION_STEP
170) endif
171) endif
172)
173) end subroutine DataMediatorDatasetInit
174)
175) ! ************************************************************************** !
176)
177) recursive subroutine DataMediatorDatasetUpdate(this,data_mediator_vec,option)
178) !
179) ! Updates a data mediator object transfering data from
180) ! the buffer into the PETSc Vec
181) !
182) ! Author: Glenn Hammond
183) ! Date: 05/01/13
184) !
185) use Option_module
186)
187) implicit none
188)
189) #include "petsc/finclude/petscvec.h"
190) #include "petsc/finclude/petscvec.h90"
191)
192) class(data_mediator_dataset_type) :: this
193) Vec :: data_mediator_vec
194) type(option_type) :: option
195)
196) PetscReal, pointer :: vec_ptr(:)
197) PetscInt :: ndof_per_cell
198) PetscInt :: mdof_local_size
199) PetscInt :: offset
200) PetscInt :: i
201) PetscErrorCode :: ierr
202)
203) call DatasetGlobalHDF5Load(this%dataset,option)
204)
205) call VecGetLocalSize(data_mediator_vec,mdof_local_size,ierr);CHKERRQ(ierr)
206) ndof_per_cell = mdof_local_size / this%dataset%local_size
207) if (mod(mdof_local_size,this%dataset%local_size) > 0) then
208) option%io_buffer = 'Mismatched vector size in MassTransferUpdate.'
209) call printErrMsg(option)
210) endif
211) call VecGetArrayF90(data_mediator_vec,vec_ptr,ierr);CHKERRQ(ierr)
212) offset = this%idof
213) do i = 1, this%dataset%local_size
214) vec_ptr(offset) = vec_ptr(offset) + this%dataset%rarray(i)
215) offset = offset + ndof_per_cell
216) enddo
217) call VecRestoreArrayF90(data_mediator_vec,vec_ptr,ierr);CHKERRQ(ierr)
218)
219) end subroutine DataMediatorDatasetUpdate
220)
221) ! ************************************************************************** !
222)
223) recursive subroutine DataMediatorDatasetStrip(this)
224) !
225) ! Destroys a data mediator object
226) !
227) ! Author: Glenn Hammond
228) ! Date: 05/01/13
229) !
230)
231) implicit none
232)
233) class(data_mediator_dataset_type) :: this
234)
235) PetscErrorCode :: ierr
236)
237) ! update the next one
238) if (associated(this%next)) then
239) call this%next%Strip()
240) deallocate(this%next)
241) nullify(this%next)
242) endif
243)
244) ! Simply nullify the pointer as the dataset resides in a list to be
245) ! destroyed separately.
246) nullify(this%dataset)
247)
248) end subroutine DataMediatorDatasetStrip
249)
250) end module Data_Mediator_Dataset_class