output_hdf5.F90 coverage: 71.43 %func 48.72 %block
1) module Output_HDF5_module
2)
3) use Logging_module
4) use Output_Aux_module
5) use Output_Common_module
6)
7) use PFLOTRAN_Constants_module
8)
9) implicit none
10)
11) private
12)
13) #include "petsc/finclude/petscsys.h"
14)
15) PetscMPIInt, private, parameter :: ON=1, OFF=0
16)
17) #if defined(SCORPIO_WRITE)
18) include "scorpiof.h"
19) #endif
20)
21) ! flags signifying the first time a routine is called during a given
22) ! simulation
23) PetscBool :: hdf5_first
24)
25) public :: OutputHDF5Init, &
26) OutputHDF5, &
27) OutputHDF5UGridXDMF, &
28) OutputHDF5FilenameID, &
29) OutputHDF5UGridXDMFExplicit, &
30) #if defined(PETSC_HAVE_HDF5)
31) OutputHDF5DatasetStringArray, &
32) OutputHDF5AttributeStringArray, &
33) #endif
34) OutputHDF5OpenFile, &
35) OutputHDF5CloseFile
36)
37) contains
38)
39) ! ************************************************************************** !
40)
41) subroutine OutputHDF5Init(num_steps)
42) !
43) ! Initializes module variables for HDF5 output
44) !
45) ! Author: Glenn Hammond
46) ! Date: 01/16/13
47) !
48)
49) use Option_module
50)
51) implicit none
52)
53) PetscInt :: num_steps
54)
55) if (num_steps == 0) then
56) hdf5_first = PETSC_TRUE
57) else
58) hdf5_first = PETSC_FALSE
59) endif
60)
61) end subroutine OutputHDF5Init
62)
63) ! ************************************************************************** !
64)
65) subroutine OutputHDF5(realization_base,var_list_type)
66) !
67) ! Print to HDF5 file
68) !
69) ! Author: Glenn Hammond
70) ! Date: 10/25/07
71) !
72)
73) use Realization_Base_class, only : realization_base_type
74) use Discretization_module
75) use Option_module
76) use Grid_module
77) use Field_module
78) use Patch_module
79) use Reaction_Aux_module
80) use String_module
81)
82) #if !defined(PETSC_HAVE_HDF5)
83) implicit none
84)
85) class(realization_base_type) :: realization_base
86) PetscInt :: var_list_type
87)
88) call printMsg(realization_base%option,'')
89) write(realization_base%option%io_buffer, &
90) '("PFLOTRAN must be compiled with HDF5 to &
91) &write HDF5 formatted structured grids Darn.")')
92) call printErrMsg(realization_base%option)
93) #else
94)
95) ! 64-bit stuff
96) #ifdef PETSC_USE_64BIT_INDICES
97) !#define HDF_NATIVE_INTEGER H5T_STD_I64LE
98) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
99) #else
100) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
101) #endif
102)
103) use hdf5
104) use HDF5_module
105) use HDF5_Aux_module
106)
107) implicit none
108)
109) #include "petsc/finclude/petscvec.h"
110) #include "petsc/finclude/petscvec.h90"
111) #include "petsc/finclude/petsclog.h"
112)
113) class(realization_base_type) :: realization_base
114) PetscInt :: var_list_type
115)
116) #if defined(SCORPIO_WRITE)
117) integer :: file_id
118) integer :: grp_id
119) integer :: file_space_id
120) integer :: realization_set_id
121) integer :: prop_id
122) PetscMPIInt :: rank
123) integer :: dims(3)
124) integer :: pio_dataset_groupid
125) #else
126) integer(HID_T) :: file_id
127) integer(HID_T) :: grp_id
128) integer(HID_T) :: file_space_id
129) integer(HID_T) :: realization_set_id
130) integer(HID_T) :: prop_id
131) PetscMPIInt :: rank
132) integer(HSIZE_T) :: dims(3)
133) #endif
134)
135) type(grid_type), pointer :: grid
136) type(discretization_type), pointer :: discretization
137) type(option_type), pointer :: option
138) type(field_type), pointer :: field
139) type(patch_type), pointer :: patch
140) type(output_option_type), pointer :: output_option
141) type(output_variable_type), pointer :: cur_variable
142)
143) Vec :: global_vec
144) Vec :: global_vec_vx
145) Vec :: global_vec_vy
146) Vec :: global_vec_vz
147) Vec :: natural_vec
148) PetscReal, pointer :: v_ptr
149)
150) character(len=MAXSTRINGLENGTH) :: string
151) character(len=MAXWORDLENGTH) :: word
152) character(len=2) :: free_mol_char, tot_mol_char, sec_mol_char
153) PetscReal, pointer :: array(:)
154) PetscInt :: i
155) PetscInt :: nviz_flow, nviz_tran, nviz_dof
156) PetscInt :: current_component
157) PetscMPIInt, parameter :: ON=1, OFF=0
158) PetscFortranAddr :: app_ptr
159) PetscMPIInt :: hdf5_err
160) PetscBool :: first
161) PetscInt :: ivar, isubvar, var_type
162) PetscErrorCode :: ierr
163)
164) discretization => realization_base%discretization
165) patch => realization_base%patch
166) option => realization_base%option
167) field => realization_base%field
168) output_option => realization_base%output_option
169)
170) call OutputHDF5OpenFile(option, output_option, var_list_type, file_id, first)
171)
172) grid => patch%grid
173) if (first) then
174) call OutputHDF5Provenance(option, output_option, file_id)
175)
176) ! create a group for the coordinates data set
177) #if defined(SCORPIO_WRITE)
178) string = "Coordinates" // CHAR(0)
179) call fscorpio_create_dataset_group(pio_dataset_groupid, string, file_id, &
180) option%iowrite_group_id, ierr)
181) ! set grp_id here
182) ! As we already created the group, we will use file_id as group_id
183) grp_id = file_id
184) #else
185) string = "Coordinates"
186) call h5gcreate_f(file_id,string,grp_id,hdf5_err,OBJECT_NAMELEN_DEFAULT_F)
187) #endif
188)
189) !GEH - Structured Grid Dependence - Begin
190) ! write out coordinates in x, y, and z directions
191) string = "X [m]"
192) allocate(array(grid%structured_grid%nx+1))
193) array(1) = discretization%origin_global(X_DIRECTION)
194) do i=2,grid%structured_grid%nx+1
195) array(i) = array(i-1) + grid%structured_grid%dx_global(i-1)
196) enddo
197) call WriteHDF5Coordinates(string,option,grid%structured_grid%nx+1, &
198) array,grp_id)
199) deallocate(array)
200)
201) string = "Y [m]"
202) allocate(array(grid%structured_grid%ny+1))
203) array(1) = discretization%origin_global(Y_DIRECTION)
204) do i=2,grid%structured_grid%ny+1
205) array(i) = array(i-1) + grid%structured_grid%dy_global(i-1)
206) enddo
207) call WriteHDF5Coordinates(string,option,grid%structured_grid%ny+1, &
208) array,grp_id)
209) deallocate(array)
210)
211) string = "Z [m]"
212) allocate(array(grid%structured_grid%nz+1))
213) array(1) = discretization%origin_global(Z_DIRECTION)
214) do i=2,grid%structured_grid%nz+1
215) array(i) = array(i-1) + grid%structured_grid%dz_global(i-1)
216) enddo
217) call WriteHDF5Coordinates(string,option,grid%structured_grid%nz+1, &
218) array,grp_id)
219) deallocate(array)
220) !GEH - Structured Grid Dependence - End
221)
222) #if defined(SCORPIO_WRITE)
223) call fscorpio_close_dataset_group(pio_dataset_groupid, file_id, &
224) option%iowrite_group_id, ierr)
225) #else
226) call h5gclose_f(grp_id,hdf5_err)
227) #endif
228)
229) endif
230)
231) ! create a group for the data set
232) write(string,'(''Time:'',es13.5,x,a1)') &
233) option%time/output_option%tconv,output_option%tunit
234) if (len_trim(output_option%plot_name) > 2) then
235) string = trim(string) // ' ' // output_option%plot_name
236) endif
237) !string = trim(string3) // ' ' // trim(string)
238) #if defined(SCORPIO_WRITE)
239) string = trim(string) //CHAR(0)
240) ! This opens existing dataset and creates it if needed
241) call fscorpio_create_dataset_group(pio_dataset_groupid, string, file_id, &
242) option%iowrite_group_id, ierr)
243) grp_id = file_id
244) #else
245) call h5eset_auto_f(OFF,hdf5_err)
246) call h5gopen_f(file_id,string,grp_id,hdf5_err)
247) if (hdf5_err /= 0) then
248) call h5gcreate_f(file_id,string,grp_id,hdf5_err,OBJECT_NAMELEN_DEFAULT_F)
249) endif
250) call h5eset_auto_f(ON,hdf5_err)
251) #endif
252) ! SCORPIO_WRITE
253)
254) ! write out data sets
255) call DiscretizationCreateVector(discretization,ONEDOF,global_vec,GLOBAL, &
256) option)
257) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vx)
258) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vy)
259) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vz)
260)
261) select case (var_list_type)
262)
263) case (INSTANTANEOUS_VARS)
264) ! loop over snapshot variables and write to file
265) cur_variable => output_option%output_snap_variable_list%first
266) do
267) if (.not.associated(cur_variable)) exit
268) call OutputGetVarFromArray(realization_base,global_vec, &
269) cur_variable%ivar, &
270) cur_variable%isubvar)
271) string = cur_variable%name
272) call StringSwapChar(string," ","_")
273) if (len_trim(cur_variable%units) > 0) then
274) word = cur_variable%units
275) call HDF5MakeStringCompatible(word)
276) string = trim(string) // ' [' // trim(word) // ']'
277) endif
278) if (cur_variable%iformat == 0) then
279) call HDF5WriteStructDataSetFromVec(string,realization_base, &
280) global_vec,grp_id, &
281) H5T_NATIVE_DOUBLE)
282) else
283) call HDF5WriteStructDataSetFromVec(string,realization_base, &
284) global_vec,grp_id, &
285) H5T_NATIVE_INTEGER)
286) endif
287) cur_variable => cur_variable%next
288) enddo
289)
290) case (AVERAGED_VARS)
291) cur_variable => output_option%aveg_output_variable_list%first
292) do ivar = 1,output_option%aveg_output_variable_list%nvars
293) string = 'Aveg. ' // cur_variable%name
294) if (len_trim(cur_variable%units) > 0) then
295) word = cur_variable%units
296) call HDF5MakeStringCompatible(word)
297) string = trim(string) // ' [' // trim(word) // ']'
298) endif
299)
300) call HDF5WriteStructDataSetFromVec(string,realization_base, &
301) field%avg_vars_vec(ivar),grp_id, &
302) H5T_NATIVE_DOUBLE)
303)
304) cur_variable => cur_variable%next
305) enddo
306)
307) end select
308)
309) if (output_option%print_hdf5_vel_cent .and. &
310) (var_list_type==INSTANTANEOUS_VARS)) then
311)
312) ! velocities
313) call OutputGetCellCenteredVelocities(realization_base, global_vec_vx, &
314) global_vec_vy,global_vec_vz, &
315) LIQUID_PHASE)
316)
317) string = "Liquid X-Velocity [m_per_" // trim(output_option%tunit) // "]"
318) call HDF5WriteStructDataSetFromVec(string,realization_base, &
319) global_vec_vx,grp_id,H5T_NATIVE_DOUBLE)
320)
321) string = "Liquid Y-Velocity [m_per_" // trim(output_option%tunit) // "]"
322) call HDF5WriteStructDataSetFromVec(string,realization_base, &
323) global_vec_vy,grp_id,H5T_NATIVE_DOUBLE)
324)
325) string = "Liquid Z-Velocity [m_per_" // trim(output_option%tunit) // "]"
326) call HDF5WriteStructDataSetFromVec(string,realization_base, &
327) global_vec_vz,grp_id,H5T_NATIVE_DOUBLE)
328)
329) if (option%nphase > 1) then
330) call OutputGetCellCenteredVelocities(realization_base,global_vec_vx, &
331) global_vec_vy,global_vec_vz, &
332) GAS_PHASE)
333) string = "Gas X-Velocity"
334) call HDF5WriteStructDataSetFromVec(string,realization_base, &
335) global_vec_vx,grp_id, &
336) H5T_NATIVE_DOUBLE)
337)
338) string = "Gas Y-Velocity"
339) call HDF5WriteStructDataSetFromVec(string,realization_base, &
340) global_vec_vy,grp_id, &
341) H5T_NATIVE_DOUBLE)
342)
343) string = "Gas Z-Velocity"
344) call HDF5WriteStructDataSetFromVec(string,realization_base, &
345) global_vec_vz,grp_id, &
346) H5T_NATIVE_DOUBLE)
347) endif
348) endif
349)
350) if (output_option%print_hdf5_vel_face .and. &
351) (var_list_type==INSTANTANEOUS_VARS)) then
352)
353) ! internal flux velocities
354) if (grid%structured_grid%nx > 1) then
355) string = "Liquid X-Flux Velocities"
356) call WriteHDF5FluxVelocities(string,realization_base,LIQUID_PHASE, &
357) X_DIRECTION,grp_id)
358) if (option%nphase > 1) then
359) string = "Gas X-Flux Velocities"
360) call WriteHDF5FluxVelocities(string,realization_base,GAS_PHASE, &
361) X_DIRECTION,grp_id)
362) endif
363) endif
364)
365) if (grid%structured_grid%ny > 1) then
366) string = "Liquid Y-Flux Velocities"
367) call WriteHDF5FluxVelocities(string,realization_base,LIQUID_PHASE, &
368) Y_DIRECTION,grp_id)
369) if (option%nphase > 1) then
370) string = "Gas Y-Flux Velocities"
371) call WriteHDF5FluxVelocities(string,realization_base,GAS_PHASE, &
372) Y_DIRECTION,grp_id)
373) endif
374) endif
375)
376) if (grid%structured_grid%nz > 1) then
377) string = "Liquid Z-Flux Velocities"
378) call WriteHDF5FluxVelocities(string,realization_base,LIQUID_PHASE, &
379) Z_DIRECTION,grp_id)
380) if (option%nphase > 1) then
381) string = "Gas Z-Flux Velocities"
382) call WriteHDF5FluxVelocities(string,realization_base,GAS_PHASE, &
383) Z_DIRECTION,grp_id)
384) endif
385) endif
386)
387) endif
388)
389) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
390) call VecDestroy(global_vec_vx,ierr);CHKERRQ(ierr)
391) call VecDestroy(global_vec_vy,ierr);CHKERRQ(ierr)
392) call VecDestroy(global_vec_vz,ierr);CHKERRQ(ierr)
393)
394) #if defined(SCORPIO_WRITE)
395) call fscorpio_close_dataset_group(pio_dataset_groupid, file_id, &
396) option%iowrite_group_id, ierr)
397) #else
398) call h5gclose_f(grp_id,hdf5_err)
399) #endif
400)
401) call OutputHDF5CloseFile(option, file_id)
402)
403) #endif
404) !PETSC_HAVE_HDF5
405)
406) hdf5_first = PETSC_FALSE
407)
408) end subroutine OutputHDF5
409)
410) ! ************************************************************************** !
411)
412) subroutine OutputHDF5OpenFile(option, output_option, var_list_type, file_id, &
413) first)
414) !
415) ! Determine the propper hdf5 output file name and open it.
416) !
417) ! Return the file handle and 'first' flag indicating if this is the
418) ! first time the file has been opened.
419) !
420) use Option_module, only : option_type, printMsg, printErrMsg
421)
422) #include "petsc/finclude/petscsysdef.h"
423)
424) #if !defined(PETSC_HAVE_HDF5)
425) implicit none
426)
427) type(option_type), intent(inout) :: option
428) type(output_option_type), intent(in) :: output_option
429) PetscInt, intent(in) :: var_list_type
430) character(len=MAXSTRINGLENGTH) :: filename
431) integer, intent(out) :: file_id
432) integer :: prop_id
433) PetscBool, intent(in) :: first
434)
435) call printMsg(option,'')
436) write(option%io_buffer, &
437) '("PFLOTRAN must be compiled with HDF5 to &
438) &write HDF5 formatted structured grids Darn.")')
439) call printErrMsg(option)
440) #else
441)
442) use hdf5
443)
444) implicit none
445)
446) type(option_type), intent(inout) :: option
447) type(output_option_type), intent(in) :: output_option
448) PetscInt, intent(in) :: var_list_type
449) character(len=MAXSTRINGLENGTH) :: filename
450) PetscBool, intent(out) :: first
451) PetscErrorCode :: ierr
452)
453) #if defined(SCORPIO_WRITE)
454) integer, intent(out) :: file_id
455) integer :: prop_id
456) #else
457) integer(HID_T), intent(out) :: file_id
458) integer(HID_T) :: prop_id
459) #endif
460)
461) character(len=MAXSTRINGLENGTH) :: string,string2,string3
462) PetscMPIInt :: hdf5_err
463)
464) select case (var_list_type)
465) case (INSTANTANEOUS_VARS)
466) string2=''
467) write(string3,'(i4)') output_option%plot_number
468) case (AVERAGED_VARS)
469) string2='-aveg'
470) write(string3,'(i4)') &
471) int(option%time/output_option%periodic_snap_output_time_incr)
472) end select
473)
474) if (output_option%print_single_h5_file) then
475) first = hdf5_first
476) filename = trim(option%global_prefix) // trim(option%group_prefix) // &
477) trim(string2) // '.h5'
478) else
479) string = OutputHDF5FilenameID(output_option,option,var_list_type)
480) select case (var_list_type)
481) case (INSTANTANEOUS_VARS)
482) if (mod(output_option%plot_number, &
483) output_option%times_per_h5_file) == 0) then
484) first = PETSC_TRUE
485) else
486) first = PETSC_FALSE
487) endif
488) case (AVERAGED_VARS)
489) if (mod((option%time-output_option%periodic_snap_output_time_incr)/ &
490) output_option%periodic_snap_output_time_incr, &
491) dble(output_option%times_per_h5_file))==0) then
492) first = PETSC_TRUE
493) else
494) first = PETSC_FALSE
495) endif
496) end select
497)
498) filename = trim(option%global_prefix) // trim(option%group_prefix) // &
499) '-' // trim(string) // trim(string2) // '.h5'
500) endif
501)
502) #if defined(SCORPIO_WRITE)
503) if (.not.first) then
504) filename = trim(filename) // CHAR(0)
505) call fscorpio_open_file(filename, option%iowrite_group_id, &
506) SCORPIO_FILE_READWRITE, file_id, ierr)
507) if (file_id == -1) first = PETSC_TRUE
508) endif
509) if (first) then
510) filename = trim(filename) // CHAR(0)
511) call fscorpio_open_file(filename, option%iowrite_group_id, &
512) SCORPIO_FILE_CREATE, file_id, ierr)
513) endif
514)
515) #else
516) ! SCORPIO_WRITE is not defined
517)
518) ! initialize fortran interface
519) call h5open_f(hdf5_err)
520)
521) call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
522) #ifndef SERIAL_HDF5
523) call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
524) #endif
525) if (.not.first) then
526) call h5eset_auto_f(OFF,hdf5_err)
527) call h5fopen_f(filename,H5F_ACC_RDWR_F,file_id,hdf5_err,prop_id)
528) if (hdf5_err /= 0) first = PETSC_TRUE
529) call h5eset_auto_f(ON,hdf5_err)
530) endif
531) if (first) then
532) call h5fcreate_f(filename,H5F_ACC_TRUNC_F,file_id,hdf5_err, &
533) H5P_DEFAULT_F,prop_id)
534) endif
535) call h5pclose_f(prop_id,hdf5_err)
536) #endif
537) ! SCORPIO_WRITE
538)
539) if (first) then
540) option%io_buffer = '--> creating hdf5 output file: ' // trim(filename)
541) else
542) option%io_buffer = '--> appending to hdf5 output file: ' // trim(filename)
543) endif
544) call printMsg(option)
545)
546) #endif
547) !PETSC_HAVE_HDF5
548)
549) end subroutine OutputHDF5OpenFile
550)
551) ! ************************************************************************** !
552)
553) subroutine OutputHDF5CloseFile(option, file_id)
554)
555) use Option_module, only : option_type, printMsg, printErrMsg
556)
557) #if !defined(PETSC_HAVE_HDF5)
558) implicit none
559)
560) type(option_type), intent(inout) :: option
561) integer, intent(in) :: file_id
562)
563) call printMsg(option,'')
564) write(option%io_buffer, &
565) '("PFLOTRAN must be compiled with HDF5 to &
566) &write HDF5 formatted structured grids Darn.")')
567) call printErrMsg(option)
568)
569) #else
570)
571) use hdf5
572)
573) implicit none
574)
575) type(option_type), intent(in) :: option
576) integer(HID_T), intent(in) :: file_id
577) integer :: hdf5_err
578) PetscErrorCode :: ierr
579)
580) #if defined(SCORPIO_WRITE)
581) call fscorpio_close_file(file_id, option%iowrite_group_id, ierr)
582) #else
583) call h5fclose_f(file_id, hdf5_err)
584) call h5close_f(hdf5_err)
585) #endif
586)
587) #endif
588) !PETSC_HAVE_HDF5
589)
590) end subroutine OutputHDF5CloseFile
591)
592) ! ************************************************************************** !
593)
594) subroutine OutputHDF5UGridXDMF(realization_base,var_list_type)
595) !
596) ! This routine writes unstructured grid data in HDF5 XDMF format.
597) !
598) ! Author: Gautam Bisht, LBNL
599) ! Date: 10/29/2012
600) !
601)
602) use Realization_Base_class, only : realization_base_type
603) use Discretization_module
604) use Option_module
605) use Grid_module
606) use Field_module
607) use Patch_module
608) use Reaction_Aux_module
609)
610) #if !defined(PETSC_HAVE_HDF5)
611)
612) implicit none
613)
614) class(realization_base_type) :: realization_base
615) PetscInt :: var_list_type
616)
617) call printMsg(realization_base%option,'')
618) write(realization_base%option%io_buffer, &
619) '("PFLOTRAN must be compiled with HDF5 to &
620) &write HDF5 formatted structured grids Darn.")')
621) call printErrMsg(realization_base%option)
622)
623) #else
624)
625) ! 64-bit stuff
626) #ifdef PETSC_USE_64BIT_INDICES
627) !#define HDF_NATIVE_INTEGER H5T_STD_I64LE
628) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
629) #else
630) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
631) #endif
632)
633) use hdf5
634) use HDF5_module, only : HDF5WriteDataSetFromVec
635) use HDF5_Aux_module
636)
637) implicit none
638)
639) #include "petsc/finclude/petscvec.h"
640) #include "petsc/finclude/petscvec.h90"
641) #include "petsc/finclude/petsclog.h"
642)
643) class(realization_base_type) :: realization_base
644) PetscInt :: var_list_type
645)
646) #if defined(SCORPIO_WRITE)
647) integer :: file_id
648) integer :: data_type
649) integer :: grp_id
650) integer :: file_space_id
651) integer :: memory_space_id
652) integer :: data_set_id
653) integer :: realization_set_id
654) integer :: prop_id
655) PetscMPIInt :: rank
656) integer :: rank_mpi,file_space_rank_mpi
657) integer :: dims(3)
658) integer :: start(3), length(3), stride(3)
659) #else
660) integer(HID_T) :: file_id
661) integer(HID_T) :: data_type
662) integer(HID_T) :: grp_id
663) integer(HID_T) :: file_space_id
664) integer(HID_T) :: realization_set_id
665) integer(HID_T) :: memory_space_id
666) integer(HID_T) :: data_set_id
667) integer(HID_T) :: prop_id
668) PetscMPIInt :: rank
669) PetscMPIInt :: rank_mpi,file_space_rank_mpi
670) integer(HSIZE_T) :: dims(3)
671) integer(HSIZE_T) :: start(3), length(3), stride(3)
672) #endif
673)
674) type(grid_type), pointer :: grid
675) type(discretization_type), pointer :: discretization
676) type(option_type), pointer :: option
677) type(field_type), pointer :: field
678) type(patch_type), pointer :: patch
679) type(output_option_type), pointer :: output_option
680) type(output_variable_type), pointer :: cur_variable
681)
682) Vec :: global_vec
683) Vec :: global_vec_vx,global_vec_vy,global_vec_vz
684) Vec :: natural_vec
685) PetscReal, pointer :: v_ptr
686)
687) character(len=MAXSTRINGLENGTH) :: filename
688) character(len=MAXSTRINGLENGTH) :: xmf_filename, att_datasetname, group_name
689) character(len=MAXSTRINGLENGTH) :: string, string2,string3
690) character(len=MAXWORDLENGTH) :: word
691) character(len=2) :: free_mol_char, tot_mol_char, sec_mol_char
692) PetscReal, pointer :: array(:)
693) PetscInt :: i
694) PetscInt :: nviz_flow, nviz_tran, nviz_dof
695) PetscInt :: current_component
696) PetscMPIInt, parameter :: ON=1, OFF=0
697) PetscFortranAddr :: app_ptr
698) PetscMPIInt :: hdf5_err
699) PetscBool :: first
700) PetscInt :: ivar, isubvar, var_type
701) PetscInt :: vert_count
702) PetscErrorCode :: ierr
703)
704) discretization => realization_base%discretization
705) patch => realization_base%patch
706) option => realization_base%option
707) field => realization_base%field
708) output_option => realization_base%output_option
709)
710) select case (var_list_type)
711) case (INSTANTANEOUS_VARS)
712) string2=''
713) write(string3,'(i4)') output_option%plot_number
714) xmf_filename = OutputFilename(output_option,option,'xmf','')
715) case (AVERAGED_VARS)
716) string2='-aveg'
717) write(string3,'(i4)') &
718) int(option%time/output_option%periodic_snap_output_time_incr)
719) xmf_filename = OutputFilename(output_option,option,'xmf','aveg')
720) end select
721)
722) if (output_option%print_single_h5_file) then
723) first = hdf5_first
724) filename = trim(option%global_prefix) // trim(string2) // &
725) trim(option%group_prefix) // '.h5'
726) else
727) string = OutputHDF5FilenameID(output_option,option,var_list_type)
728) select case (var_list_type)
729) case (INSTANTANEOUS_VARS)
730) if (mod(output_option%plot_number, &
731) output_option%times_per_h5_file) == 0) then
732) first = PETSC_TRUE
733) else
734) first = PETSC_FALSE
735) endif
736) case (AVERAGED_VARS)
737) if (mod((option%time-output_option%periodic_snap_output_time_incr)/ &
738) output_option%periodic_snap_output_time_incr, &
739) dble(output_option%times_per_h5_file))==0) then
740) first = PETSC_TRUE
741) else
742) first = PETSC_FALSE
743) endif
744) end select
745)
746) filename = trim(option%global_prefix) // trim(option%group_prefix) // &
747) trim(string2) // '-' // trim(string) // '.h5'
748) endif
749)
750) grid => patch%grid
751)
752) #ifdef SCORPIO_WRITE
753) option%io_buffer='OutputHDF5UGridXDMF not supported with SCORPIO_WRITE'
754) call printErrMsg(option)
755) #endif
756)
757) ! initialize fortran interface
758) call h5open_f(hdf5_err)
759)
760) call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
761) #ifndef SERIAL_HDF5
762) call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
763) #endif
764)
765) if (.not.first) then
766) call h5eset_auto_f(OFF,hdf5_err)
767) call h5fopen_f(filename,H5F_ACC_RDWR_F,file_id,hdf5_err,prop_id)
768) if (hdf5_err /= 0) first = PETSC_TRUE
769) call h5eset_auto_f(ON,hdf5_err)
770) endif
771) if (first) then
772) call h5fcreate_f(filename,H5F_ACC_TRUNC_F,file_id,hdf5_err, &
773) H5P_DEFAULT_F,prop_id)
774) else if (Uninitialized(realization_base%output_option%xmf_vert_len)) then
775) call DetermineNumVertices(realization_base,option)
776) endif
777) call h5pclose_f(prop_id,hdf5_err)
778)
779) if (first) then
780) option%io_buffer = '--> creating hdf5 output file: ' // trim(filename)
781) else
782) option%io_buffer = '--> appending to hdf5 output file: ' // trim(filename)
783) endif
784) call printMsg(option)
785)
786) if (first) then
787) ! create a group for the coordinates data set
788) string = "Domain"
789) call h5gcreate_f(file_id,string,grp_id,hdf5_err,OBJECT_NAMELEN_DEFAULT_F)
790) call WriteHDF5CoordinatesUGridXDMF(realization_base,option,grp_id)
791) call h5gclose_f(grp_id,hdf5_err)
792) endif
793)
794) if (option%myrank == option%io_rank) then
795) option%io_buffer = '--> write xmf output file: ' // trim(filename)
796) call printMsg(option)
797) open(unit=OUTPUT_UNIT,file=xmf_filename,action="write")
798) !call OutputXMFHeader(OUTPUT_UNIT,realization_base,filename)
799) call OutputXMFHeader(OUTPUT_UNIT, &
800) option%time/output_option%tconv, &
801) grid%nmax, &
802) realization_base%output_option%xmf_vert_len, &
803) grid%unstructured_grid%num_vertices_global,filename)
804) endif
805)
806) ! create a group for the data set
807) write(string,'(''Time'',es13.5,x,a1)') &
808) option%time/output_option%tconv,output_option%tunit
809) if (len_trim(output_option%plot_name) > 2) then
810) string = trim(string) // ' ' // output_option%plot_name
811) endif
812) string = trim(string3) // ' ' // trim(string)
813)
814) call h5eset_auto_f(OFF,hdf5_err)
815) call h5gopen_f(file_id,string,grp_id,hdf5_err)
816) group_name=string
817) if (hdf5_err /= 0) then
818) call h5gcreate_f(file_id,string,grp_id,hdf5_err,OBJECT_NAMELEN_DEFAULT_F)
819) endif
820) call h5eset_auto_f(ON,hdf5_err)
821)
822) ! write out data sets
823) call DiscretizationCreateVector(discretization,ONEDOF,global_vec,GLOBAL, &
824) option)
825) call DiscretizationCreateVector(discretization,ONEDOF,natural_vec,NATURAL, &
826) option)
827) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vx)
828) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vy)
829) call DiscretizationDuplicateVector(discretization,global_vec,global_vec_vz)
830)
831) select case (var_list_type)
832)
833) case (INSTANTANEOUS_VARS)
834) ! loop over snapshot variables and write to file
835) cur_variable => output_option%output_snap_variable_list%first
836) do
837) if (.not.associated(cur_variable)) exit
838) call OutputGetVarFromArray(realization_base,global_vec, &
839) cur_variable%ivar, &
840) cur_variable%isubvar)
841) call DiscretizationGlobalToNatural(discretization,global_vec, &
842) natural_vec,ONEDOF)
843) string = cur_variable%name
844) if (len_trim(cur_variable%units) > 0) then
845) word = cur_variable%units
846) call HDF5MakeStringCompatible(word)
847) string = trim(string) // ' [' // trim(word) // ']'
848) endif
849) if (cur_variable%iformat == 0) then
850) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
851) H5T_NATIVE_DOUBLE)
852) else
853) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
854) H5T_NATIVE_INTEGER)
855) endif
856) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
857) "/" // trim(string)
858) if (option%myrank == option%io_rank) then
859) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string,att_datasetname)
860) endif
861) cur_variable => cur_variable%next
862) enddo
863)
864) case (AVERAGED_VARS)
865) if (associated(output_option%aveg_output_variable_list%first)) then
866) cur_variable => output_option%aveg_output_variable_list%first
867) do ivar = 1,output_option%aveg_output_variable_list%nvars
868) string = 'Aveg. ' // cur_variable%name
869) if (len_trim(cur_variable%units) > 0) then
870) word = cur_variable%units
871) call HDF5MakeStringCompatible(word)
872) string = trim(string) // ' [' // trim(word) // ']'
873) endif
874)
875) call DiscretizationGlobalToNatural(discretization, &
876) field%avg_vars_vec(ivar), &
877) natural_vec,ONEDOF)
878) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
879) H5T_NATIVE_DOUBLE)
880) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
881) "/" // trim(string)
882) if (option%myrank == option%io_rank) then
883) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string, &
884) att_datasetname)
885) endif
886) cur_variable => cur_variable%next
887) enddo
888) endif
889)
890) end select
891)
892) !Output flowrates
893) if (output_option%print_hdf5_mass_flowrate.or. &
894) output_option%print_hdf5_energy_flowrate.or. &
895) output_option%print_hdf5_aveg_mass_flowrate.or. &
896) output_option%print_hdf5_aveg_energy_flowrate) then
897)
898) select case (var_list_type)
899) case (INSTANTANEOUS_VARS)
900) call OutputGetFaceFlowrateUGrid(realization_base)
901) if (output_option%print_hdf5_mass_flowrate.or.&
902) output_option%print_hdf5_energy_flowrate) then
903) call WriteHDF5FlowratesUGrid(realization_base,option,grp_id, &
904) var_list_type)
905) endif
906) case (AVERAGED_VARS)
907) if (output_option%print_hdf5_aveg_mass_flowrate.or.&
908) output_option%print_hdf5_aveg_energy_flowrate) then
909) call WriteHDF5FlowratesUGrid(realization_base,option,grp_id, &
910) var_list_type)
911) endif
912) end select
913) endif
914)
915) if (output_option%print_hdf5_vel_cent .and. &
916) (var_list_type==INSTANTANEOUS_VARS)) then
917)
918) ! velocities
919) call OutputGetCellCenteredVelocities(realization_base,global_vec_vx, &
920) global_vec_vy,global_vec_vz, &
921) LIQUID_PHASE)
922)
923) string = "Liquid X-Velocity [m_per_" // trim(output_option%tunit) // "]"
924) call DiscretizationGlobalToNatural(discretization,global_vec_vx, &
925) natural_vec,ONEDOF)
926) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
927) H5T_NATIVE_DOUBLE)
928) att_datasetname = trim(filename) // ":/" // trim(group_name) // "/" &
929) // trim(string)
930) if (option%myrank == option%io_rank) then
931) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string,att_datasetname)
932) endif
933)
934) string = "Liquid Y-Velocity [m_per_" // trim(output_option%tunit) // "]"
935) call DiscretizationGlobalToNatural(discretization,global_vec_vy, &
936) natural_vec,ONEDOF)
937) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
938) H5T_NATIVE_DOUBLE)
939) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
940) "/" // trim(string)
941) if (option%myrank == option%io_rank) then
942) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string,att_datasetname)
943) endif
944)
945) string = "Liquid Z-Velocity [m_per_" // trim(output_option%tunit) // "]"
946) call DiscretizationGlobalToNatural(discretization,global_vec_vz, &
947) natural_vec,ONEDOF)
948) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
949) H5T_NATIVE_DOUBLE)
950) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
951) "/" // trim(string)
952) if (option%myrank == option%io_rank) then
953) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string,att_datasetname)
954) endif
955)
956) if (option%nphase > 1) then
957) call OutputGetCellCenteredVelocities(realization_base,global_vec_vx, &
958) global_vec_vy,global_vec_vz, &
959) GAS_PHASE)
960)
961) string = "Gas X-Velocity [m_per_" // trim(output_option%tunit) // "]"
962) call DiscretizationGlobalToNatural(discretization,global_vec_vx, &
963) natural_vec,ONEDOF)
964) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
965) H5T_NATIVE_DOUBLE)
966) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
967) "/" // trim(string)
968) if (option%myrank == option%io_rank) then
969) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string,att_datasetname)
970) endif
971)
972) string = "Gas Y-Velocity [m_per_" // trim(output_option%tunit) // "]"
973) call DiscretizationGlobalToNatural(discretization,global_vec_vy, &
974) natural_vec,ONEDOF)
975) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
976) H5T_NATIVE_DOUBLE)
977) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
978) "/" // trim(string)
979) if (option%myrank == option%io_rank) then
980) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string,att_datasetname)
981) endif
982)
983) string = "Gas Z-Velocity [m_per_" // trim(output_option%tunit) // "]"
984) call DiscretizationGlobalToNatural(discretization,global_vec_vz, &
985) natural_vec,ONEDOF)
986) call HDF5WriteDataSetFromVec(string,option,natural_vec,grp_id, &
987) H5T_NATIVE_DOUBLE)
988) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
989) "/" // trim(string)
990) if (option%myrank == option%io_rank) then
991) call OutputXMFAttribute(OUTPUT_UNIT,grid%nmax,string,att_datasetname)
992) endif
993) endif
994) endif
995)
996) ! Output velocity at cell-face
997) if (output_option%print_hdf5_vel_face) then
998)
999) select case (var_list_type)
1000) case (INSTANTANEOUS_VARS)
1001) call OutputGetFaceVelUGrid(realization_base)
1002) if (output_option%print_hdf5_vel_face) then
1003) call WriteHDF5FaceVelUGrid(realization_base,option,grp_id, &
1004) var_list_type)
1005) endif
1006) case (AVERAGED_VARS)
1007) end select
1008) endif
1009)
1010) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
1011) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
1012) call VecDestroy(global_vec_vx,ierr);CHKERRQ(ierr)
1013) call VecDestroy(global_vec_vy,ierr);CHKERRQ(ierr)
1014) call VecDestroy(global_vec_vz,ierr);CHKERRQ(ierr)
1015)
1016) call h5gclose_f(grp_id,hdf5_err)
1017)
1018) call h5fclose_f(file_id,hdf5_err)
1019) call h5close_f(hdf5_err)
1020)
1021) if (option%myrank == option%io_rank) then
1022) call OutputXMFFooter(OUTPUT_UNIT)
1023) close(OUTPUT_UNIT)
1024) endif
1025)
1026) hdf5_first = PETSC_FALSE
1027)
1028) #endif
1029) ! !defined(PETSC_HAVE_HDF5)
1030)
1031) end subroutine OutputHDF5UGridXDMF
1032)
1033) ! ************************************************************************** !
1034)
1035) subroutine OutputHDF5UGridXDMFExplicit(realization_base,var_list_type)
1036) !
1037) ! This subroutine prints the explicit
1038) ! unstructured grid information in xdmf format
1039) !
1040) ! Author: Satish Karra, LANL
1041) ! Date: 07/17/2013
1042) !
1043)
1044) use Realization_Base_class, only : realization_base_type
1045) use Discretization_module
1046) use Option_module
1047) use Grid_module
1048) use Field_module
1049) use Patch_module
1050) use Reaction_Aux_module
1051)
1052) #if !defined(PETSC_HAVE_HDF5)
1053)
1054) implicit none
1055)
1056) class(realization_base_type) :: realization_base
1057) PetscInt :: var_list_type
1058)
1059) call printMsg(realization_base%option,'')
1060) write(realization_base%option%io_buffer, &
1061) '("PFLOTRAN must be compiled with HDF5 to &
1062) &write HDF5 formatted structured grids Darn.")')
1063) call printErrMsg(realization_base%option)
1064)
1065) #else
1066)
1067) ! 64-bit stuff
1068) #ifdef PETSC_USE_64BIT_INDICES
1069) !#define HDF_NATIVE_INTEGER H5T_STD_I64LE
1070) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
1071) #else
1072) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
1073) #endif
1074)
1075) use hdf5
1076) use HDF5_module, only : HDF5WriteDataSetFromVec
1077) use HDF5_Aux_module
1078)
1079) implicit none
1080)
1081) #include "petsc/finclude/petscvec.h"
1082) #include "petsc/finclude/petscvec.h90"
1083) #include "petsc/finclude/petsclog.h"
1084)
1085) class(realization_base_type) :: realization_base
1086) PetscInt :: var_list_type
1087)
1088) #if defined(SCORPIO_WRITE)
1089) integer :: file_id, new_file_id
1090) integer :: data_type
1091) integer :: grp_id, new_grp_id
1092) integer :: file_space_id
1093) integer :: memory_space_id
1094) integer :: data_set_id
1095) integer :: realization_set_id
1096) integer :: prop_id, new_prop_id
1097) PetscMPIInt :: rank
1098) integer :: rank_mpi,file_space_rank_mpi
1099) integer :: dims(3)
1100) integer :: start(3), length(3), stride(3)
1101) #else
1102) integer(HID_T) :: file_id, new_file_id
1103) integer(HID_T) :: data_type
1104) integer(HID_T) :: grp_id, new_grp_id
1105) integer(HID_T) :: file_space_id
1106) integer(HID_T) :: realization_set_id
1107) integer(HID_T) :: memory_space_id
1108) integer(HID_T) :: data_set_id
1109) integer(HID_T) :: prop_id, new_prop_id
1110) PetscMPIInt :: rank
1111) PetscMPIInt :: rank_mpi,file_space_rank_mpi
1112) integer(HSIZE_T) :: dims(3)
1113) integer(HSIZE_T) :: start(3), length(3), stride(3)
1114) #endif
1115)
1116) type(grid_type), pointer :: grid
1117) type(discretization_type), pointer :: discretization
1118) type(option_type), pointer :: option
1119) type(field_type), pointer :: field
1120) type(patch_type), pointer :: patch
1121) type(output_option_type), pointer :: output_option
1122) type(output_variable_type), pointer :: cur_variable
1123)
1124) Vec :: global_vec
1125) Vec :: natural_vec
1126) PetscReal, pointer :: v_ptr
1127)
1128) character(len=MAXSTRINGLENGTH) :: filename
1129) character(len=MAXSTRINGLENGTH) :: xmf_filename, att_datasetname, group_name
1130) character(len=MAXSTRINGLENGTH) :: new_filename
1131) character(len=MAXSTRINGLENGTH) :: string, string2,string3
1132) character(len=MAXWORDLENGTH) :: word
1133) character(len=2) :: free_mol_char, tot_mol_char, sec_mol_char
1134) PetscReal, pointer :: array(:)
1135) PetscInt :: i
1136) PetscInt :: nviz_flow, nviz_tran, nviz_dof
1137) PetscInt :: current_component
1138) PetscMPIInt, parameter :: ON=1, OFF=0
1139) PetscFortranAddr :: app_ptr
1140) PetscMPIInt :: hdf5_err
1141) PetscBool :: first
1142) PetscInt :: ivar, isubvar, var_type
1143) PetscInt :: istart
1144) PetscInt :: vert_count
1145) PetscErrorCode :: ierr
1146)
1147) discretization => realization_base%discretization
1148) patch => realization_base%patch
1149) option => realization_base%option
1150) field => realization_base%field
1151) output_option => realization_base%output_option
1152)
1153) select case (var_list_type)
1154) case (INSTANTANEOUS_VARS)
1155) string2=''
1156) write(string3,'(i4)') output_option%plot_number
1157) xmf_filename = OutputFilename(output_option,option,'xmf','')
1158) case (AVERAGED_VARS)
1159) string2='-aveg'
1160) write(string3,'(i4)') &
1161) int(option%time/output_option%periodic_snap_output_time_incr)
1162) xmf_filename = OutputFilename(output_option,option,'xmf','aveg')
1163) end select
1164)
1165) if (output_option%print_single_h5_file) then
1166) first = hdf5_first
1167) filename = trim(option%global_prefix) // trim(string2) // &
1168) trim(option%group_prefix) // '.h5'
1169) else
1170) string = OutputHDF5FilenameID(output_option,option,var_list_type)
1171) select case (var_list_type)
1172) case (INSTANTANEOUS_VARS)
1173) if (mod(output_option%plot_number, &
1174) output_option%times_per_h5_file) == 0) then
1175) first = PETSC_TRUE
1176) else
1177) first = PETSC_FALSE
1178) endif
1179) case (AVERAGED_VARS)
1180) if (mod((option%time-output_option%periodic_snap_output_time_incr)/ &
1181) output_option%periodic_snap_output_time_incr, &
1182) dble(output_option%times_per_h5_file))==0) then
1183) first = PETSC_TRUE
1184) else
1185) first = PETSC_FALSE
1186) endif
1187) end select
1188)
1189) filename = trim(option%global_prefix) // trim(option%group_prefix) // &
1190) trim(string2) // '-' // trim(string) // '.h5'
1191) endif
1192)
1193) grid => patch%grid
1194)
1195) #ifdef SCORPIO_WRITE
1196) option%io_buffer='OutputHDF5UGridXDMF not supported with SCORPIO_WRITE'
1197) call printErrMsg(option)
1198) #endif
1199)
1200) ! initialize fortran interface
1201) call h5open_f(hdf5_err)
1202)
1203) call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
1204) #ifndef SERIAL_HDF5
1205) call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
1206) #endif
1207)
1208) if (.not.first) then
1209) call h5eset_auto_f(OFF,hdf5_err)
1210) call h5fopen_f(filename,H5F_ACC_RDWR_F,file_id,hdf5_err,prop_id)
1211) if (hdf5_err /= 0) first = PETSC_TRUE
1212) call h5eset_auto_f(ON,hdf5_err)
1213) endif
1214) if (first) then
1215) call h5fcreate_f(filename,H5F_ACC_TRUNC_F,file_id,hdf5_err, &
1216) H5P_DEFAULT_F,prop_id)
1217) endif
1218) call h5pclose_f(prop_id,hdf5_err)
1219)
1220) if (first) then
1221) option%io_buffer = '--> creating hdf5 output file: ' // trim(filename)
1222) else
1223) option%io_buffer = '--> appending to hdf5 output file: ' // trim(filename)
1224) endif
1225) call printMsg(option)
1226)
1227) new_filename = trim(option%global_prefix) // '-domain.h5'
1228)
1229)
1230) if (first .and. option%print_explicit_primal_grid) then
1231) if (option%myrank == option%io_rank) then
1232) call h5pcreate_f(H5P_FILE_ACCESS_F,new_prop_id,hdf5_err)
1233) call h5fcreate_f(new_filename,H5F_ACC_TRUNC_F,new_file_id,hdf5_err, &
1234) H5P_DEFAULT_F,new_prop_id)
1235) call h5pclose_f(new_prop_id,hdf5_err)
1236) ! create a group for the coordinates data set
1237) string = "Domain"
1238) call h5gcreate_f(new_file_id,string,new_grp_id,hdf5_err, &
1239) OBJECT_NAMELEN_DEFAULT_F)
1240) call WriteHDF5CoordinatesUGridXDMFExplicit(realization_base,option, &
1241) new_grp_id)
1242) call h5gclose_f(new_grp_id,hdf5_err)
1243) call h5fclose_f(new_file_id,hdf5_err)
1244) endif
1245) endif
1246)
1247) if (option%myrank == option%io_rank .and. &
1248) option%print_explicit_primal_grid) then
1249) option%io_buffer = '--> write xmf output file: ' // trim(filename)
1250) call printMsg(option)
1251) open(unit=OUTPUT_UNIT,file=xmf_filename,action="write")
1252) call OutputXMFHeaderExplicit(OUTPUT_UNIT, &
1253) option%time/output_option%tconv, &
1254) grid%unstructured_grid%explicit_grid%num_elems, &
1255) realization_base%output_option%xmf_vert_len, &
1256) grid%unstructured_grid%explicit_grid%num_cells_global, &
1257) new_filename)
1258) endif
1259)
1260) ! create a group for the data set
1261) write(string,'(''Time'',es13.5,x,a1)') &
1262) option%time/output_option%tconv,output_option%tunit
1263) if (len_trim(output_option%plot_name) > 2) then
1264) string = trim(string) // ' ' // output_option%plot_name
1265) endif
1266) string = trim(string3) // ' ' // trim(string)
1267)
1268) call h5eset_auto_f(OFF,hdf5_err)
1269) call h5gopen_f(file_id,string,grp_id,hdf5_err)
1270) group_name=string
1271) if (hdf5_err /= 0) then
1272) call h5gcreate_f(file_id,string,grp_id,hdf5_err,OBJECT_NAMELEN_DEFAULT_F)
1273) endif
1274) call h5eset_auto_f(ON,hdf5_err)
1275)
1276) ! write out data sets
1277) call DiscretizationCreateVector(discretization,ONEDOF,global_vec,GLOBAL, &
1278) option)
1279) call DiscretizationCreateVector(discretization,ONEDOF,natural_vec,NATURAL, &
1280) option)
1281)
1282) select case (var_list_type)
1283)
1284) case (INSTANTANEOUS_VARS)
1285) ! loop over snapshot variables and write to file
1286) cur_variable => output_option%output_snap_variable_list%first
1287) do
1288) if (.not.associated(cur_variable)) exit
1289) call OutputGetVarFromArray(realization_base,global_vec, &
1290) cur_variable%ivar, &
1291) cur_variable%isubvar)
1292) call DiscretizationGlobalToNatural(discretization,global_vec, &
1293) natural_vec,ONEDOF)
1294) string = cur_variable%name
1295) if (len_trim(cur_variable%units) > 0) then
1296) word = cur_variable%units
1297) call HDF5MakeStringCompatible(word)
1298) string = trim(string) // ' [' // trim(word) // ']'
1299) endif
1300) if (cur_variable%iformat == 0) then
1301) call HDF5WriteDataSetFromVec(string,option,natural_vec, &
1302) grp_id,H5T_NATIVE_DOUBLE)
1303) else
1304) call HDF5WriteDataSetFromVec(string,option,natural_vec, &
1305) grp_id,H5T_NATIVE_INTEGER)
1306) endif
1307) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
1308) "/" // trim(string)
1309) if (option%myrank == option%io_rank .and. &
1310) option%print_explicit_primal_grid) then
1311) call OutputXMFAttributeExplicit(OUTPUT_UNIT,grid%nmax,string, &
1312) att_datasetname)
1313) endif
1314) cur_variable => cur_variable%next
1315) enddo
1316)
1317) case (AVERAGED_VARS)
1318) if (associated(output_option%aveg_output_variable_list%first)) then
1319) cur_variable => output_option%aveg_output_variable_list%first
1320) do ivar = 1,output_option%aveg_output_variable_list%nvars
1321) string = 'Aveg. ' // cur_variable%name
1322) if (len_trim(cur_variable%units) > 0) then
1323) word = cur_variable%units
1324) call HDF5MakeStringCompatible(word)
1325) string = trim(string) // ' [' // trim(word) // ']'
1326) endif
1327)
1328) call DiscretizationGlobalToNatural(discretization, &
1329) field%avg_vars_vec(ivar), &
1330) natural_vec,ONEDOF)
1331) call HDF5WriteDataSetFromVec(string,option,natural_vec, &
1332) grp_id,H5T_NATIVE_DOUBLE)
1333) att_datasetname = trim(filename) // ":/" // trim(group_name) // &
1334) "/" // trim(string)
1335) if (option%myrank == option%io_rank .and. &
1336) option%print_explicit_primal_grid) then
1337) call OutputXMFAttributeExplicit(OUTPUT_UNIT,grid%nmax,string, &
1338) att_datasetname)
1339) endif
1340) cur_variable => cur_variable%next
1341) enddo
1342) endif
1343)
1344) end select
1345)
1346) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
1347) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
1348) call h5gclose_f(grp_id,hdf5_err)
1349)
1350) call h5fclose_f(file_id,hdf5_err)
1351) call h5close_f(hdf5_err)
1352)
1353)
1354) if (option%myrank == option%io_rank) then
1355) call OutputXMFFooter(OUTPUT_UNIT)
1356) close(OUTPUT_UNIT)
1357) endif
1358)
1359) hdf5_first = PETSC_FALSE
1360)
1361) #endif
1362) ! !defined(PETSC_HAVE_HDF5)
1363)
1364) end subroutine OutputHDF5UGridXDMFExplicit
1365)
1366) ! ************************************************************************** !
1367)
1368) function OutputHDF5FilenameID(output_option,option,var_list_type)
1369) !
1370) ! This subroutine creates an ID for HDF5 filename for:
1371) ! - Instantaneous, or
1372) ! - Temporally averaged variables.
1373) !
1374) ! Author: Gautam Bisht, LBNL
1375) ! Date: 01/10/13
1376) !
1377)
1378) use Option_module
1379)
1380) implicit none
1381)
1382) type(option_type) :: option
1383) type(output_option_type) :: output_option
1384) PetscInt :: var_list_type
1385)
1386) character(len=MAXWORDLENGTH) :: OutputHDF5FilenameID
1387) PetscInt :: file_number
1388)
1389) select case(var_list_type)
1390) case (INSTANTANEOUS_VARS)
1391) file_number = floor(real(output_option%plot_number)/ &
1392) output_option%times_per_h5_file)
1393) case (AVERAGED_VARS)
1394) file_number = floor((option%time - &
1395) output_option%periodic_snap_output_time_incr)/ &
1396) output_option%periodic_snap_output_time_incr/ &
1397) output_option%times_per_h5_file)
1398) end select
1399)
1400) if (file_number < 10) then
1401) write(OutputHDF5FilenameID,'("00",i1)') file_number
1402) else if (output_option%plot_number < 100) then
1403) write(OutputHDF5FilenameID,'("0",i2)') file_number
1404) else if (output_option%plot_number < 1000) then
1405) write(OutputHDF5FilenameID,'(i3)') file_number
1406) else if (output_option%plot_number < 10000) then
1407) write(OutputHDF5FilenameID,'(i4)') file_number
1408) endif
1409)
1410) OutputHDF5FilenameID = adjustl(OutputHDF5FilenameID)
1411)
1412) end function OutputHDF5FilenameID
1413)
1414) #if defined(PETSC_HAVE_HDF5)
1415)
1416) ! ************************************************************************** !
1417)
1418) subroutine WriteHDF5FluxVelocities(name,realization_base,iphase,direction, &
1419) file_id)
1420) !
1421) ! Print flux velocities to HDF5 file
1422) !
1423) ! Author: Glenn Hammond
1424) ! Date: 10/25/07
1425) !
1426)
1427) use Realization_Base_class, only : realization_base_type
1428) use Discretization_module
1429) use Grid_module
1430) use Option_module
1431) use Field_module
1432) use Connection_module
1433) use Patch_module
1434) use hdf5
1435) use HDF5_module, only : HDF5WriteStructuredDataSet, trick_hdf5
1436)
1437) implicit none
1438)
1439) #include "petsc/finclude/petscvec.h"
1440) #include "petsc/finclude/petscvec.h90"
1441) #include "petsc/finclude/petsclog.h"
1442)
1443) character(len=32) :: name
1444) class(realization_base_type) :: realization_base
1445) PetscInt :: iphase
1446) PetscInt :: direction
1447) integer(HID_T) :: file_id
1448)
1449) PetscInt :: i, j, k
1450) PetscInt :: count, iconn
1451) PetscInt :: local_id, ghosted_id
1452) PetscInt :: nx_local, ny_local, nz_local
1453) PetscInt :: nx_global, ny_global, nz_global
1454) PetscErrorCode :: ierr
1455)
1456) type(grid_type), pointer :: grid
1457) type(discretization_type), pointer :: discretization
1458) type(option_type), pointer :: option
1459) type(field_type), pointer :: field
1460) type(patch_type), pointer :: patch
1461) type(output_option_type), pointer :: output_option
1462)
1463) PetscReal, allocatable :: array(:)
1464) PetscReal, pointer :: vec_ptr(:)
1465)
1466) PetscBool, save :: trick_flux_vel_x = PETSC_FALSE
1467) PetscBool, save :: trick_flux_vel_y = PETSC_FALSE
1468) PetscBool, save :: trick_flux_vel_z = PETSC_FALSE
1469)
1470) Vec :: global_vec
1471)
1472) type(connection_set_list_type), pointer :: connection_set_list
1473) type(connection_set_type), pointer :: cur_connection_set
1474)
1475) discretization => realization_base%discretization
1476) patch => realization_base%patch
1477) grid => patch%grid
1478) option => realization_base%option
1479) field => realization_base%field
1480) output_option => realization_base%output_option
1481)
1482) ! in a few cases (i.e. for small test problems), some processors may
1483) ! have no velocities to print. This results in zero-length arrays
1484) ! in collective H5Dwrite(). To avoid, we switch to independent
1485) ! H5Dwrite() and don't write from the zero-length procs.
1486) !GEH - Structured Grid Dependence - Begin
1487) if (hdf5_first) then
1488) trick_flux_vel_x = PETSC_FALSE
1489) trick_flux_vel_y = PETSC_FALSE
1490) trick_flux_vel_z = PETSC_FALSE
1491)
1492) nx_local = grid%structured_grid%nlx
1493) ny_local = grid%structured_grid%nly
1494) nz_local = grid%structured_grid%nlz
1495) if (grid%structured_grid%gxe-grid%structured_grid%lxe == 0) then
1496) nx_local = grid%structured_grid%nlx-1
1497) endif
1498) call MPI_Allreduce(nx_local,i,ONE_INTEGER_MPI,MPIU_INTEGER,MPI_MIN, &
1499) option%mycomm,ierr)
1500) if (i == 0) trick_flux_vel_x = PETSC_TRUE
1501) if (grid%structured_grid%gye-grid%structured_grid%lye == 0) then
1502) ny_local = grid%structured_grid%nly-1
1503) endif
1504) call MPI_Allreduce(ny_local,j,ONE_INTEGER_MPI,MPIU_INTEGER,MPI_MIN, &
1505) option%mycomm,ierr)
1506) if (j == 0) trick_flux_vel_y = PETSC_TRUE
1507) if (grid%structured_grid%gze-grid%structured_grid%lze == 0) then
1508) nz_local = grid%structured_grid%nlz-1
1509) endif
1510) call MPI_Allreduce(nz_local,k,ONE_INTEGER_MPI,MPIU_INTEGER,MPI_MIN, &
1511) option%mycomm,ierr)
1512) if (k == 0) trick_flux_vel_z = PETSC_TRUE
1513) endif
1514)
1515) nx_local = grid%structured_grid%nlx
1516) ny_local = grid%structured_grid%nly
1517) nz_local = grid%structured_grid%nlz
1518) nx_global = grid%structured_grid%nx
1519) ny_global = grid%structured_grid%ny
1520) nz_global = grid%structured_grid%nz
1521)
1522) select case(direction)
1523) case(X_DIRECTION)
1524) nx_global = grid%structured_grid%nx-1
1525) if (grid%structured_grid%gxe-grid%structured_grid%lxe == 0) then
1526) nx_local = grid%structured_grid%nlx-1
1527) endif
1528) if (trick_flux_vel_x) trick_hdf5 = PETSC_TRUE
1529) case(Y_DIRECTION)
1530) ny_global = grid%structured_grid%ny-1
1531) if (grid%structured_grid%gye-grid%structured_grid%lye == 0) then
1532) ny_local = grid%structured_grid%nly-1
1533) endif
1534) if (trick_flux_vel_y) trick_hdf5 = PETSC_TRUE
1535) case(Z_DIRECTION)
1536) nz_global = grid%structured_grid%nz-1
1537) if (grid%structured_grid%gze-grid%structured_grid%lze == 0) then
1538) nz_local = grid%structured_grid%nlz-1
1539) endif
1540) if (trick_flux_vel_z) trick_hdf5 = PETSC_TRUE
1541) end select
1542) allocate(array(nx_local*ny_local*nz_local))
1543)
1544)
1545) call DiscretizationCreateVector(discretization,ONEDOF,global_vec,GLOBAL, &
1546) option)
1547) call VecZeroEntries(global_vec,ierr);CHKERRQ(ierr)
1548) call VecGetArrayF90(global_vec,vec_ptr,ierr);CHKERRQ(ierr)
1549)
1550) ! place interior velocities in a vector
1551) connection_set_list => grid%internal_connection_set_list
1552) cur_connection_set => connection_set_list%first
1553) do
1554) if (.not.associated(cur_connection_set)) exit
1555) do iconn = 1, cur_connection_set%num_connections
1556) ghosted_id = cur_connection_set%id_up(iconn)
1557) local_id = grid%nG2L(ghosted_id) ! = zero for ghost nodes
1558) ! velocities are stored as the downwind face of the upwind cell
1559) if (local_id <= 0 .or. &
1560) dabs(cur_connection_set%dist(direction,iconn)) < 0.99d0) cycle
1561) vec_ptr(local_id) = patch%internal_velocities(iphase,iconn)
1562) enddo
1563) cur_connection_set => cur_connection_set%next
1564) enddo
1565)
1566) count = 0
1567) do k=1,nz_local
1568) do j=1,ny_local
1569) do i=1,nx_local
1570) count = count + 1
1571) local_id = i + (j-1)*grid%structured_grid%nlx + &
1572) (k-1)*grid%structured_grid%nlxy
1573) array(count) = vec_ptr(local_id)
1574) enddo
1575) enddo
1576) enddo
1577) call VecRestoreArrayF90(global_vec,vec_ptr,ierr);CHKERRQ(ierr)
1578)
1579) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
1580)
1581) array(1:nx_local*ny_local*nz_local) = & ! convert time units
1582) array(1:nx_local*ny_local*nz_local) * output_option%tconv
1583)
1584) call HDF5WriteStructuredDataSet(name,array,file_id,H5T_NATIVE_DOUBLE, &
1585) option,nx_global,ny_global,nz_global, &
1586) nx_local,ny_local,nz_local, &
1587) grid%structured_grid%lxs, &
1588) grid%structured_grid%lys, &
1589) grid%structured_grid%lzs)
1590) !GEH - Structured Grid Dependence - End
1591)
1592) deallocate(array)
1593) trick_hdf5 = PETSC_FALSE
1594)
1595) end subroutine WriteHDF5FluxVelocities
1596)
1597) ! ************************************************************************** !
1598)
1599) subroutine WriteHDF5Coordinates(name,option,length,array,file_id)
1600) !
1601) ! Writes structured coordinates to HDF5 file
1602) !
1603) ! Author: Glenn Hammond
1604) ! Date: 10/25/07
1605) !
1606)
1607) use hdf5
1608) use Option_module
1609)
1610) implicit none
1611)
1612) #if defined(SCORPIO_WRITE)
1613) character(len=32) :: name
1614) type(option_type) :: option
1615) PetscInt :: length
1616) PetscReal :: array(:)
1617) integer :: file_id
1618)
1619) integer :: file_space_id
1620) integer :: data_set_id
1621) integer :: prop_id
1622) integer :: dims(3)
1623) PetscMPIInt :: rank
1624) integer :: globaldims(3)
1625) #else
1626) character(len=32) :: name
1627) type(option_type) :: option
1628) PetscInt :: length
1629) PetscReal :: array(:)
1630) integer(HID_T) :: file_id
1631)
1632) integer(HID_T) :: file_space_id
1633) integer(HID_T) :: data_set_id
1634) integer(HID_T) :: prop_id
1635) integer(HSIZE_T) :: dims(3)
1636) PetscMPIInt :: rank
1637) #endif
1638)
1639) PetscMPIInt :: hdf5_flag
1640) PetscMPIInt :: hdf5_err
1641) PetscErrorCode :: ierr
1642)
1643) call PetscLogEventBegin(logging%event_output_coordinates_hdf5, &
1644) ierr);CHKERRQ(ierr)
1645) #if defined(SCORPIO_WRITE)
1646)
1647) name = trim(name) // CHAR(0)
1648) ! write out grid structure
1649) rank = 1
1650) dims = 0
1651) globaldims = 0
1652) ! x-direction
1653)
1654) ! Only process 0 writes coordinates
1655) if (option%myrank == 0 ) then
1656) dims(1) = length
1657) globaldims(1) = length
1658) else
1659) dims(1) = 0
1660) globaldims(1) = length
1661) endif
1662)
1663) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1664) call fscorpio_write_dataset(array, SCORPIO_DOUBLE, rank, globaldims, dims, &
1665) file_id, name, option%iowrite_group_id, SCORPIO_NONUNIFORM_CONTIGUOUS_WRITE, &
1666) ierr)
1667) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1668)
1669) #else
1670) !SCORPIO_WRITE is not defined
1671)
1672) ! write out grid structure
1673) rank = 1
1674) dims = 0
1675) ! x-direction
1676) dims(1) = length
1677) call h5screate_simple_f(rank,dims,file_space_id,hdf5_err,dims)
1678) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
1679) call h5dcreate_f(file_id,name,H5T_NATIVE_DOUBLE,file_space_id, &
1680) data_set_id,hdf5_err,prop_id)
1681) call h5pclose_f(prop_id,hdf5_err)
1682)
1683) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
1684) #ifndef SERIAL_HDF5
1685) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F,hdf5_err) ! must be independent and only from p0
1686) #endif
1687) if (option%myrank == option%io_rank) then
1688) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1689) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,array,dims, &
1690) hdf5_err,H5S_ALL_F,H5S_ALL_F,prop_id)
1691) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1692) endif
1693) call h5pclose_f(prop_id,hdf5_err)
1694) call h5dclose_f(data_set_id,hdf5_err)
1695) call h5sclose_f(file_space_id,hdf5_err)
1696)
1697) #endif
1698) ! SCORPIO_WRITE
1699)
1700) call PetscLogEventEnd(logging%event_output_coordinates_hdf5, &
1701) ierr);CHKERRQ(ierr)
1702)
1703) end subroutine WriteHDF5Coordinates
1704)
1705) ! ************************************************************************** !
1706)
1707) subroutine WriteHDF5CoordinatesUGrid(grid,option,file_id)
1708) !
1709) ! This subroutine writes unstructured coordinates to HDF5 file
1710) !
1711) ! Author: Gautam Bisht, ORNL
1712) ! Date: 05/31/12
1713) !
1714)
1715) use hdf5
1716) use Grid_module
1717) use Option_module
1718) use Grid_Unstructured_Aux_module
1719) use HDF5_module, only : trick_hdf5
1720) use Variables_module
1721)
1722) implicit none
1723)
1724) #include "petsc/finclude/petscvec.h"
1725) #include "petsc/finclude/petscvec.h90"
1726) #include "petsc/finclude/petsclog.h"
1727)
1728) type(grid_type), pointer :: grid
1729) type(option_type), pointer :: option
1730)
1731) #if defined(SCORPIO_WRITE)
1732) integer :: file_id
1733) integer :: data_type
1734) integer :: grp_id
1735) integer :: file_space_id
1736) integer :: memory_space_id
1737) integer :: data_set_id
1738) integer :: realization_set_id
1739) integer :: prop_id
1740) integer :: dims(3)
1741) integer :: start(3), length(3), stride(3)
1742) integer :: rank_mpi,file_space_rank_mpi
1743) integer :: hdf5_flag
1744) integer, parameter :: ON=1, OFF=0
1745) #else
1746) integer(HID_T) :: file_id
1747) integer(HID_T) :: data_type
1748) integer(HID_T) :: grp_id
1749) integer(HID_T) :: file_space_id
1750) integer(HID_T) :: realization_set_id
1751) integer(HID_T) :: memory_space_id
1752) integer(HID_T) :: data_set_id
1753) integer(HID_T) :: prop_id
1754) integer(HSIZE_T) :: dims(3)
1755) integer(HSIZE_T) :: start(3), length(3), stride(3)
1756) PetscMPIInt :: rank_mpi,file_space_rank_mpi
1757) PetscMPIInt :: hdf5_flag
1758) PetscMPIInt, parameter :: ON=1, OFF=0
1759) #endif
1760)
1761) character(len=MAXSTRINGLENGTH) :: string
1762) PetscMPIInt :: hdf5_err
1763)
1764) PetscInt :: local_size
1765) PetscInt :: istart
1766) PetscInt :: i,j
1767) PetscReal, pointer :: vec_x_ptr(:),vec_y_ptr(:),vec_z_ptr(:)
1768) PetscReal, pointer :: double_array(:)
1769) Vec :: global_x_vertex_vec,global_y_vertex_vec,global_z_vertex_vec
1770)
1771) PetscReal, pointer :: vec_ptr(:)
1772) Vec :: global_vec, natural_vec
1773) PetscInt, pointer :: int_array(:)
1774) type(ugdm_type),pointer :: ugdm_element
1775) PetscErrorCode :: ierr
1776)
1777) call VecCreateMPI(option%mycomm,PETSC_DECIDE, &
1778) grid%unstructured_grid%num_vertices_global, &
1779) global_x_vertex_vec,ierr);CHKERRQ(ierr)
1780) call VecCreateMPI(option%mycomm,PETSC_DECIDE, &
1781) grid%unstructured_grid%num_vertices_global, &
1782) global_y_vertex_vec,ierr);CHKERRQ(ierr)
1783) call VecCreateMPI(option%mycomm,PETSC_DECIDE, &
1784) grid%unstructured_grid%num_vertices_global, &
1785) global_z_vertex_vec,ierr);CHKERRQ(ierr)
1786)
1787) call VecGetLocalSize(global_x_vertex_vec,local_size,ierr);CHKERRQ(ierr)
1788) call VecGetLocalSize(global_y_vertex_vec,local_size,ierr);CHKERRQ(ierr)
1789) call VecGetLocalSize(global_z_vertex_vec,local_size,ierr);CHKERRQ(ierr)
1790)
1791) call GetVertexCoordinates(grid, global_x_vertex_vec,X_COORDINATE,option)
1792) call GetVertexCoordinates(grid, global_y_vertex_vec,Y_COORDINATE,option)
1793) call GetVertexCoordinates(grid, global_z_vertex_vec,Z_COORDINATE,option)
1794)
1795) call VecGetArrayF90(global_x_vertex_vec,vec_x_ptr,ierr);CHKERRQ(ierr)
1796) call VecGetArrayF90(global_y_vertex_vec,vec_y_ptr,ierr);CHKERRQ(ierr)
1797) call VecGetArrayF90(global_z_vertex_vec,vec_z_ptr,ierr);CHKERRQ(ierr)
1798)
1799) #if defined(SCORPIO_WRITE)
1800) write(*,*) 'SCORPIO_WRITE'
1801) option%io_buffer = 'WriteHDF5CoordinatesUGrid not supported for SCORPIO_WRITE'
1802) call printErrMsg(option)
1803) #else
1804)
1805) !
1806) ! not(SCORPIO_WRITE)
1807) !
1808)
1809) ! memory space which is a 1D vector
1810) rank_mpi = 1
1811) dims = 0
1812) dims(1) = local_size * 3
1813) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
1814)
1815) ! file space which is a 3D block
1816) rank_mpi = 2
1817) dims = 0
1818) dims(2) = grid%unstructured_grid%num_vertices_global
1819) dims(1) = 3
1820) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
1821)
1822) string = "Vertices" // CHAR(0)
1823)
1824) call h5eset_auto_f(OFF,hdf5_err)
1825) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
1826) hdf5_flag = hdf5_err
1827) call h5eset_auto_f(ON,hdf5_err)
1828) if (hdf5_flag < 0) then
1829) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
1830) call h5dcreate_f(file_id,string,H5T_NATIVE_DOUBLE,file_space_id, &
1831) data_set_id,hdf5_err,prop_id)
1832) else
1833) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
1834) endif
1835)
1836) call h5pclose_f(prop_id,hdf5_err)
1837)
1838) istart = 0
1839) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
1840) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
1841)
1842) start(2) = istart
1843) start(1) = 0
1844)
1845) length(2) = local_size
1846) length(1) = 3
1847)
1848) stride = 1
1849) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
1850) hdf5_err,stride,stride)
1851)
1852) ! write the data
1853) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
1854) #ifndef SERIAL_HDF5
1855) if (trick_hdf5) then
1856) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
1857) hdf5_err)
1858) else
1859) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_COLLECTIVE_F, &
1860) hdf5_err)
1861) endif
1862) #endif
1863)
1864) allocate(double_array(local_size*3))
1865)
1866) do i=1,local_size
1867) double_array((i-1)*3+1) = vec_x_ptr(i)
1868) double_array((i-1)*3+2) = vec_y_ptr(i)
1869) double_array((i-1)*3+3) = vec_z_ptr(i)
1870) enddo
1871)
1872) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1873) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,double_array,dims, &
1874) hdf5_err,memory_space_id,file_space_id,prop_id)
1875) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1876)
1877) deallocate(double_array)
1878) call h5pclose_f(prop_id,hdf5_err)
1879)
1880)
1881) call h5dclose_f(data_set_id,hdf5_err)
1882) call h5sclose_f(file_space_id,hdf5_err)
1883)
1884)
1885) #endif
1886)
1887) call VecRestoreArrayF90(global_x_vertex_vec,vec_x_ptr,ierr);CHKERRQ(ierr)
1888) call VecRestoreArrayF90(global_y_vertex_vec,vec_y_ptr,ierr);CHKERRQ(ierr)
1889) call VecRestoreArrayF90(global_z_vertex_vec,vec_z_ptr,ierr);CHKERRQ(ierr)
1890)
1891)
1892) call VecDestroy(global_x_vertex_vec,ierr);CHKERRQ(ierr)
1893) call VecDestroy(global_y_vertex_vec,ierr);CHKERRQ(ierr)
1894) call VecDestroy(global_z_vertex_vec,ierr);CHKERRQ(ierr)
1895)
1896)
1897) !
1898) ! Write elements
1899) !
1900) call UGridCreateUGDM(grid%unstructured_grid,ugdm_element,EIGHT_INTEGER,option)
1901) call UGridDMCreateVector(grid%unstructured_grid,ugdm_element,global_vec, &
1902) GLOBAL,option)
1903) call UGridDMCreateVector(grid%unstructured_grid,ugdm_element,natural_vec, &
1904) NATURAL,option)
1905) call GetCellConnections(grid,global_vec)
1906) call VecScatterBegin(ugdm_element%scatter_gton,global_vec,natural_vec, &
1907) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
1908) call VecScatterEnd(ugdm_element%scatter_gton,global_vec,natural_vec, &
1909) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
1910) call VecGetArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
1911)
1912) local_size = grid%unstructured_grid%nlmax
1913) #if defined(SCORPIO_WRITE)
1914) write(*,*) 'SCORPIO_WRITE'
1915) option%io_buffer = 'WriteHDF5CoordinatesUGrid not supported for SCORPIO_WRITE'
1916) call printErrMsg(option)
1917) #else
1918)
1919) !
1920) ! not(SCORPIO_WRITE)
1921) !
1922)
1923) ! memory space which is a 1D vector
1924) rank_mpi = 1
1925) dims = 0
1926) dims(1) = local_size*NINE_INTEGER
1927) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
1928)
1929) ! file space which is a 3D block
1930) rank_mpi = 2
1931) dims = 0
1932) dims(2) = grid%unstructured_grid%nmax
1933) dims(1) = NINE_INTEGER
1934) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
1935)
1936) string = "Cells" // CHAR(0)
1937)
1938) call h5eset_auto_f(OFF,hdf5_err)
1939) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
1940) hdf5_flag = hdf5_err
1941) call h5eset_auto_f(ON,hdf5_err)
1942) if (hdf5_flag < 0) then
1943) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
1944) call h5dcreate_f(file_id,string,H5T_NATIVE_INTEGER,file_space_id, &
1945) data_set_id,hdf5_err,prop_id)
1946) else
1947) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
1948) endif
1949)
1950) call h5pclose_f(prop_id,hdf5_err)
1951)
1952) istart = 0
1953) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
1954) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
1955)
1956) start(2) = istart
1957) start(1) = 0
1958)
1959) length(2) = local_size
1960) length(1) = NINE_INTEGER
1961)
1962) stride = 1
1963) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
1964) hdf5_err,stride,stride)
1965)
1966) ! write the data
1967) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
1968) #ifndef SERIAL_HDF5
1969) if (trick_hdf5) then
1970) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
1971) hdf5_err)
1972) else
1973) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_COLLECTIVE_F, &
1974) hdf5_err)
1975) endif
1976) #endif
1977)
1978) allocate(int_array(local_size*NINE_INTEGER))
1979)
1980) do i=1,local_size
1981) int_array((i-1)*9 + 1) = 0
1982) int_array((i-1)*9 + 2) = INT(vec_ptr((i-1)*8+1))
1983) int_array((i-1)*9 + 3) = INT(vec_ptr((i-1)*8+2))
1984) int_array((i-1)*9 + 4) = INT(vec_ptr((i-1)*8+3))
1985) int_array((i-1)*9 + 5) = INT(vec_ptr((i-1)*8+4))
1986) int_array((i-1)*9 + 6) = INT(vec_ptr((i-1)*8+5))
1987) int_array((i-1)*9 + 7) = INT(vec_ptr((i-1)*8+6))
1988) int_array((i-1)*9 + 8) = INT(vec_ptr((i-1)*8+7))
1989) int_array((i-1)*9 + 9) = INT(vec_ptr((i-1)*8+8))
1990) do j=2,9
1991) if (int_array((i-1)*9 + j)>0) int_array((i-1)*9 + 1)= int_array((i-1)*9 + 1) +1
1992) enddo
1993) enddo
1994)
1995) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1996) call h5dwrite_f(data_set_id,H5T_NATIVE_INTEGER,int_array,dims, &
1997) hdf5_err,memory_space_id,file_space_id,prop_id)
1998) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
1999)
2000) deallocate(int_array)
2001) call h5pclose_f(prop_id,hdf5_err)
2002)
2003)
2004) call h5dclose_f(data_set_id,hdf5_err)
2005) call h5sclose_f(file_space_id,hdf5_err)
2006)
2007) #endif
2008)
2009) call VecRestoreArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
2010) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
2011) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
2012) call UGridDMDestroy(ugdm_element)
2013)
2014) end subroutine WriteHDF5CoordinatesUGrid
2015)
2016) ! ************************************************************************** !
2017)
2018) subroutine WriteHDF5CoordinatesUGridXDMF(realization_base,option,file_id)
2019) !
2020) ! This routine writes unstructured coordinates to HDF5 file in XDMF format
2021) !
2022) ! Author: Gautam Bisht, LBNL
2023) ! Date: 10/29/2012
2024) !
2025)
2026) use hdf5
2027) use Realization_Base_class, only : realization_base_type
2028) use Grid_module
2029) use Option_module
2030) use Grid_Unstructured_Aux_module
2031) use Variables_module
2032)
2033) implicit none
2034)
2035) #include "petsc/finclude/petscvec.h"
2036) #include "petsc/finclude/petscvec.h90"
2037) #include "petsc/finclude/petsclog.h"
2038)
2039) class(realization_base_type) :: realization_base
2040) type(option_type), pointer :: option
2041)
2042) #if defined(SCORPIO_WRITE)
2043) integer :: file_id
2044) integer :: data_type
2045) integer :: grp_id
2046) integer :: file_space_id
2047) integer :: memory_space_id
2048) integer :: data_set_id
2049) integer :: realization_set_id
2050) integer :: prop_id
2051) integer :: dims(3)
2052) integer :: start(3), length(3), stride(3)
2053) integer :: rank_mpi,file_space_rank_mpi
2054) integer :: hdf5_flag
2055) integer, parameter :: ON=1, OFF=0
2056) #else
2057) integer(HID_T) :: file_id
2058) integer(HID_T) :: data_type
2059) integer(HID_T) :: grp_id
2060) integer(HID_T) :: file_space_id
2061) integer(HID_T) :: realization_set_id
2062) integer(HID_T) :: memory_space_id
2063) integer(HID_T) :: data_set_id
2064) integer(HID_T) :: prop_id
2065) integer(HSIZE_T) :: dims(3)
2066) integer(HSIZE_T) :: start(3), length(3), stride(3)
2067) PetscMPIInt :: rank_mpi,file_space_rank_mpi
2068) PetscMPIInt :: hdf5_flag
2069) PetscMPIInt, parameter :: ON=1, OFF=0
2070) #endif
2071)
2072) type(grid_type), pointer :: grid
2073) character(len=MAXSTRINGLENGTH) :: string
2074) PetscMPIInt :: hdf5_err
2075)
2076) PetscInt :: local_size,vert_count,nverts
2077) PetscInt :: i,j
2078) PetscInt :: temp_int, istart
2079) PetscReal, pointer :: vec_x_ptr(:),vec_y_ptr(:),vec_z_ptr(:)
2080) PetscReal, pointer :: double_array(:)
2081) Vec :: global_x_vertex_vec,global_y_vertex_vec,global_z_vertex_vec
2082) Vec :: global_x_cell_vec,global_y_cell_vec,global_z_cell_vec
2083) Vec :: natural_x_cell_vec,natural_y_cell_vec,natural_z_cell_vec
2084)
2085) PetscReal, pointer :: vec_ptr(:)
2086) Vec :: global_vec, natural_vec
2087) PetscInt, pointer :: int_array(:)
2088) type(ugdm_type),pointer :: ugdm_element, ugdm_cell
2089) PetscErrorCode :: ierr
2090)
2091) PetscInt :: TET_ID_XDMF = 6
2092) PetscInt :: PYR_ID_XDMF = 7
2093) PetscInt :: WED_ID_XDMF = 8
2094) PetscInt :: HEX_ID_XDMF = 9
2095)
2096) grid => realization_base%patch%grid
2097)
2098) call VecCreateMPI(option%mycomm,PETSC_DECIDE, &
2099) grid%unstructured_grid%num_vertices_global, &
2100) global_x_vertex_vec,ierr);CHKERRQ(ierr)
2101) call VecCreateMPI(option%mycomm,PETSC_DECIDE, &
2102) grid%unstructured_grid%num_vertices_global, &
2103) global_y_vertex_vec,ierr);CHKERRQ(ierr)
2104) call VecCreateMPI(option%mycomm,PETSC_DECIDE, &
2105) grid%unstructured_grid%num_vertices_global, &
2106) global_z_vertex_vec,ierr);CHKERRQ(ierr)
2107)
2108) call VecGetLocalSize(global_x_vertex_vec,local_size,ierr);CHKERRQ(ierr)
2109) call VecGetLocalSize(global_y_vertex_vec,local_size,ierr);CHKERRQ(ierr)
2110) call VecGetLocalSize(global_z_vertex_vec,local_size,ierr);CHKERRQ(ierr)
2111)
2112) call GetVertexCoordinates(grid, global_x_vertex_vec,X_COORDINATE,option)
2113) call GetVertexCoordinates(grid, global_y_vertex_vec,Y_COORDINATE,option)
2114) call GetVertexCoordinates(grid, global_z_vertex_vec,Z_COORDINATE,option)
2115)
2116) call VecGetArrayF90(global_x_vertex_vec,vec_x_ptr,ierr);CHKERRQ(ierr)
2117) call VecGetArrayF90(global_y_vertex_vec,vec_y_ptr,ierr);CHKERRQ(ierr)
2118) call VecGetArrayF90(global_z_vertex_vec,vec_z_ptr,ierr);CHKERRQ(ierr)
2119)
2120) #if defined(SCORPIO_WRITE)
2121) write(*,*) 'SCORPIO_WRITE'
2122) option%io_buffer = 'WriteHDF5CoordinatesUGrid not supported for SCORPIO_WRITE'
2123) call printErrMsg(option)
2124) #else
2125)
2126) !
2127) ! not(SCORPIO_WRITE)
2128) !
2129)
2130) ! memory space which is a 1D vector
2131) rank_mpi = 1
2132) dims = 0
2133) dims(1) = local_size * 3
2134) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
2135)
2136) ! file space which is a 2D block
2137) rank_mpi = 2
2138) dims = 0
2139) dims(2) = grid%unstructured_grid%num_vertices_global
2140) dims(1) = 3
2141) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
2142)
2143) string = "Vertices" // CHAR(0)
2144)
2145) call h5eset_auto_f(OFF,hdf5_err)
2146) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
2147) hdf5_flag = hdf5_err
2148) call h5eset_auto_f(ON,hdf5_err)
2149) if (hdf5_flag < 0) then
2150) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
2151) call h5dcreate_f(file_id,string,H5T_NATIVE_DOUBLE,file_space_id, &
2152) data_set_id,hdf5_err,prop_id)
2153) else
2154) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
2155) endif
2156)
2157) call h5pclose_f(prop_id,hdf5_err)
2158)
2159) istart = 0
2160) !geh: cannot use dims(1) in MPI_Allreduce as it causes errors on
2161) ! Juqueen
2162) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
2163) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
2164)
2165) start(2) = istart
2166) start(1) = 0
2167)
2168) length(2) = local_size
2169) length(1) = 3
2170)
2171) stride = 1
2172) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
2173) hdf5_err,stride,stride)
2174) ! write the data
2175) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
2176) #ifndef SERIAL_HDF5
2177) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
2178) hdf5_err)
2179) #endif
2180)
2181) allocate(double_array(local_size*3))
2182) do i=1,local_size
2183) double_array((i-1)*3+1) = vec_x_ptr(i)
2184) double_array((i-1)*3+2) = vec_y_ptr(i)
2185) double_array((i-1)*3+3) = vec_z_ptr(i)
2186) enddo
2187)
2188) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2189) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,double_array,dims, &
2190) hdf5_err,memory_space_id,file_space_id,prop_id)
2191) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2192)
2193) deallocate(double_array)
2194) call h5pclose_f(prop_id,hdf5_err)
2195)
2196) call h5dclose_f(data_set_id,hdf5_err)
2197) call h5sclose_f(file_space_id,hdf5_err)
2198)
2199) call VecRestoreArrayF90(global_x_vertex_vec,vec_x_ptr,ierr);CHKERRQ(ierr)
2200) call VecRestoreArrayF90(global_y_vertex_vec,vec_y_ptr,ierr);CHKERRQ(ierr)
2201) call VecRestoreArrayF90(global_z_vertex_vec,vec_z_ptr,ierr);CHKERRQ(ierr)
2202)
2203)
2204) call VecDestroy(global_x_vertex_vec,ierr);CHKERRQ(ierr)
2205) call VecDestroy(global_y_vertex_vec,ierr);CHKERRQ(ierr)
2206) call VecDestroy(global_z_vertex_vec,ierr);CHKERRQ(ierr)
2207)
2208) !
2209) ! Write elements
2210) !
2211) call UGridCreateUGDM(grid%unstructured_grid,ugdm_element,EIGHT_INTEGER,option)
2212) call UGridDMCreateVector(grid%unstructured_grid,ugdm_element,global_vec, &
2213) GLOBAL,option)
2214) call UGridDMCreateVector(grid%unstructured_grid,ugdm_element,natural_vec, &
2215) NATURAL,option)
2216) call GetCellConnections(grid,global_vec)
2217) call VecScatterBegin(ugdm_element%scatter_gton,global_vec,natural_vec, &
2218) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
2219) call VecScatterEnd(ugdm_element%scatter_gton,global_vec,natural_vec, &
2220) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
2221) call VecGetArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
2222)
2223) local_size = grid%unstructured_grid%nlmax
2224)
2225) vert_count=0
2226) do i=1,local_size*EIGHT_INTEGER
2227) if (int(vec_ptr(i)) >0 ) vert_count=vert_count+1
2228) enddo
2229) vert_count=vert_count+grid%nlmax
2230)
2231) ! memory space which is a 1D vector
2232) rank_mpi = 1
2233) dims(1) = vert_count
2234) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
2235)
2236) !geh: cannot use dims(1) in MPI_Allreduce as it causes errors on
2237) ! Juqueen
2238) call MPI_Allreduce(vert_count,temp_int,ONE_INTEGER_MPI,MPIU_INTEGER, &
2239) MPI_SUM,option%mycomm,ierr)
2240) dims(1) = temp_int
2241) realization_base%output_option%xmf_vert_len=int(dims(1))
2242)
2243) ! file space which is a 2D block
2244) rank_mpi = 1
2245) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
2246)
2247) string = "Cells" // CHAR(0)
2248)
2249) call h5eset_auto_f(OFF,hdf5_err)
2250) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
2251) hdf5_flag = hdf5_err
2252) call h5eset_auto_f(ON,hdf5_err)
2253) if (hdf5_flag < 0) then
2254) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
2255) call h5dcreate_f(file_id,string,H5T_NATIVE_INTEGER,file_space_id, &
2256) data_set_id,hdf5_err,prop_id)
2257) else
2258) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
2259) endif
2260)
2261) call h5pclose_f(prop_id,hdf5_err)
2262)
2263) istart = 0
2264) call MPI_Exscan(vert_count, istart, ONE_INTEGER_MPI, &
2265) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
2266)
2267) start(1) = istart
2268) length(1) = vert_count
2269) stride = 1
2270) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
2271) hdf5_err,stride,stride)
2272)
2273) ! write the data
2274) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
2275) #ifndef SERIAL_HDF5
2276) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
2277) hdf5_err)
2278) #endif
2279)
2280) allocate(int_array(vert_count))
2281)
2282) vert_count=0
2283) do i=1,local_size
2284) nverts=0
2285) do j=1,8
2286) if (vec_ptr((i-1)*8+j)>0) nverts=nverts+1
2287) enddo
2288) vert_count=vert_count+1
2289) select case (nverts)
2290) case (4) ! Tetrahedron
2291) int_array(vert_count) = TET_ID_XDMF
2292) case (5) ! Pyramid
2293) int_array(vert_count) = PYR_ID_XDMF
2294) case (6) ! Wedge
2295) int_array(vert_count) = WED_ID_XDMF
2296) case (8) ! Hexahedron
2297) int_array(vert_count) = HEX_ID_XDMF
2298) end select
2299)
2300) do j=1,8
2301) if (vec_ptr((i-1)*8+j)>0) then
2302) vert_count=vert_count+1
2303) int_array(vert_count) = INT(vec_ptr((i-1)*8+j))-1
2304) endif
2305) enddo
2306) enddo
2307)
2308) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2309) call h5dwrite_f(data_set_id,H5T_NATIVE_INTEGER,int_array,dims, &
2310) hdf5_err,memory_space_id,file_space_id,prop_id)
2311) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2312)
2313) deallocate(int_array)
2314) call h5pclose_f(prop_id,hdf5_err)
2315)
2316) call h5dclose_f(data_set_id,hdf5_err)
2317) call h5sclose_f(file_space_id,hdf5_err)
2318)
2319) call VecRestoreArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
2320) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
2321) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
2322) call UGridDMDestroy(ugdm_element)
2323)
2324) ! Cell center X/Y/Z
2325) call VecCreateMPI(option%mycomm,grid%nlmax, &
2326) PETSC_DETERMINE, &
2327) global_x_cell_vec,ierr);CHKERRQ(ierr)
2328) call VecCreateMPI(option%mycomm,grid%nlmax, &
2329) PETSC_DETERMINE, &
2330) global_y_cell_vec,ierr);CHKERRQ(ierr)
2331) call VecCreateMPI(option%mycomm,grid%nlmax, &
2332) PETSC_DETERMINE, &
2333) global_z_cell_vec,ierr);CHKERRQ(ierr)
2334)
2335) call GetCellCoordinates(grid, global_x_cell_vec,X_COORDINATE)
2336) call GetCellCoordinates(grid, global_y_cell_vec,Y_COORDINATE)
2337) call GetCellCoordinates(grid, global_z_cell_vec,Z_COORDINATE)
2338)
2339)
2340) call UGridCreateUGDM(grid%unstructured_grid,ugdm_cell,ONE_INTEGER,option)
2341) call UGridDMCreateVector(grid%unstructured_grid,ugdm_cell, &
2342) natural_x_cell_vec,NATURAL,option)
2343) call UGridDMCreateVector(grid%unstructured_grid,ugdm_cell, &
2344) natural_y_cell_vec,NATURAL,option)
2345) call UGridDMCreateVector(grid%unstructured_grid,ugdm_cell, &
2346) natural_z_cell_vec,NATURAL,option)
2347)
2348) call VecScatterBegin(ugdm_cell%scatter_gton,global_x_cell_vec, &
2349) natural_x_cell_vec,INSERT_VALUES,SCATTER_FORWARD, &
2350) ierr);CHKERRQ(ierr)
2351) call VecScatterEnd(ugdm_cell%scatter_gton,global_x_cell_vec, &
2352) natural_x_cell_vec,INSERT_VALUES,SCATTER_FORWARD, &
2353) ierr);CHKERRQ(ierr)
2354)
2355) call VecScatterBegin(ugdm_cell%scatter_gton,global_y_cell_vec, &
2356) natural_y_cell_vec,INSERT_VALUES,SCATTER_FORWARD, &
2357) ierr);CHKERRQ(ierr)
2358) call VecScatterEnd(ugdm_cell%scatter_gton,global_y_cell_vec, &
2359) natural_y_cell_vec,INSERT_VALUES,SCATTER_FORWARD, &
2360) ierr);CHKERRQ(ierr)
2361)
2362) call VecScatterBegin(ugdm_cell%scatter_gton,global_z_cell_vec, &
2363) natural_z_cell_vec,INSERT_VALUES,SCATTER_FORWARD, &
2364) ierr);CHKERRQ(ierr)
2365) call VecScatterEnd(ugdm_cell%scatter_gton,global_z_cell_vec, &
2366) natural_z_cell_vec,INSERT_VALUES,SCATTER_FORWARD, &
2367) ierr);CHKERRQ(ierr)
2368)
2369) call VecGetArrayF90(natural_x_cell_vec,vec_x_ptr,ierr);CHKERRQ(ierr)
2370) call VecGetArrayF90(natural_y_cell_vec,vec_y_ptr,ierr);CHKERRQ(ierr)
2371) call VecGetArrayF90(natural_z_cell_vec,vec_z_ptr,ierr);CHKERRQ(ierr)
2372) local_size = grid%unstructured_grid%nlmax
2373)
2374) ! XC
2375) ! memory space which is a 1D vector
2376) rank_mpi = 1
2377) dims = 0
2378) dims(1) = local_size
2379) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
2380)
2381) ! file space which is a 2D block
2382) rank_mpi = 1
2383) dims = 0
2384) dims(1) = grid%nmax
2385) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
2386)
2387) string = "XC" // CHAR(0)
2388)
2389) call h5eset_auto_f(OFF,hdf5_err)
2390) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
2391) hdf5_flag = hdf5_err
2392) call h5eset_auto_f(ON,hdf5_err)
2393) if (hdf5_flag < 0) then
2394) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
2395) call h5dcreate_f(file_id,string,H5T_NATIVE_DOUBLE,file_space_id, &
2396) data_set_id,hdf5_err,prop_id)
2397) else
2398) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
2399) endif
2400)
2401) call h5pclose_f(prop_id,hdf5_err)
2402)
2403) istart = 0
2404) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
2405) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
2406) start(1) = istart
2407) length(1) = local_size
2408)
2409) stride = 1
2410) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
2411) hdf5_err,stride,stride)
2412) ! write the data
2413) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
2414) #ifndef SERIAL_HDF5
2415) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
2416) hdf5_err)
2417) #endif
2418)
2419) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2420) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,vec_x_ptr,dims, &
2421) hdf5_err,memory_space_id,file_space_id,prop_id)
2422) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2423)
2424) call h5pclose_f(prop_id,hdf5_err)
2425)
2426) call h5dclose_f(data_set_id,hdf5_err)
2427) call h5sclose_f(file_space_id,hdf5_err)
2428)
2429) ! YC
2430) ! memory space which is a 1D vector
2431) rank_mpi = 1
2432) dims = 0
2433) dims(1) = local_size
2434) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
2435)
2436) ! file space which is a 2D block
2437) rank_mpi = 1
2438) dims = 0
2439) dims(1) = grid%nmax
2440) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
2441)
2442) string = "YC" // CHAR(0)
2443)
2444) call h5eset_auto_f(OFF,hdf5_err)
2445) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
2446) hdf5_flag = hdf5_err
2447) call h5eset_auto_f(ON,hdf5_err)
2448) if (hdf5_flag < 0) then
2449) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
2450) call h5dcreate_f(file_id,string,H5T_NATIVE_DOUBLE,file_space_id, &
2451) data_set_id,hdf5_err,prop_id)
2452) else
2453) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
2454) endif
2455)
2456) call h5pclose_f(prop_id,hdf5_err)
2457)
2458) istart = 0
2459) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
2460) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
2461) start(1) = istart
2462) length(1) = local_size
2463)
2464) stride = 1
2465) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
2466) hdf5_err,stride,stride)
2467) ! write the data
2468) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
2469) #ifndef SERIAL_HDF5
2470) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
2471) hdf5_err)
2472) #endif
2473)
2474) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2475) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,vec_y_ptr,dims, &
2476) hdf5_err,memory_space_id,file_space_id,prop_id)
2477) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2478)
2479) call h5pclose_f(prop_id,hdf5_err)
2480)
2481) call h5dclose_f(data_set_id,hdf5_err)
2482) call h5sclose_f(file_space_id,hdf5_err)
2483)
2484) ! ZC
2485) ! memory space which is a 1D vector
2486) rank_mpi = 1
2487) dims = 0
2488) dims(1) = local_size
2489) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
2490)
2491) ! file space which is a 2D block
2492) rank_mpi = 1
2493) dims = 0
2494) dims(1) = grid%nmax
2495) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
2496)
2497) string = "ZC" // CHAR(0)
2498)
2499) call h5eset_auto_f(OFF,hdf5_err)
2500) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
2501) hdf5_flag = hdf5_err
2502) call h5eset_auto_f(ON,hdf5_err)
2503) if (hdf5_flag < 0) then
2504) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
2505) call h5dcreate_f(file_id,string,H5T_NATIVE_DOUBLE,file_space_id, &
2506) data_set_id,hdf5_err,prop_id)
2507) else
2508) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
2509) endif
2510)
2511) call h5pclose_f(prop_id,hdf5_err)
2512)
2513) istart = 0
2514) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
2515) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
2516) start(1) = istart
2517) length(1) = local_size
2518)
2519) stride = 1
2520) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
2521) hdf5_err,stride,stride)
2522) ! write the data
2523) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
2524) #ifndef SERIAL_HDF5
2525) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
2526) hdf5_err)
2527) #endif
2528)
2529) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2530) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,vec_z_ptr,dims, &
2531) hdf5_err,memory_space_id,file_space_id,prop_id)
2532) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2533)
2534) call h5pclose_f(prop_id,hdf5_err)
2535)
2536) call h5dclose_f(data_set_id,hdf5_err)
2537) call h5sclose_f(file_space_id,hdf5_err)
2538)
2539)
2540) call VecRestoreArrayF90(natural_x_cell_vec,vec_x_ptr,ierr);CHKERRQ(ierr)
2541) call VecRestoreArrayF90(natural_y_cell_vec,vec_y_ptr,ierr);CHKERRQ(ierr)
2542) call VecRestoreArrayF90(natural_z_cell_vec,vec_z_ptr,ierr);CHKERRQ(ierr)
2543)
2544) call VecDestroy(global_x_cell_vec,ierr);CHKERRQ(ierr)
2545) call VecDestroy(global_y_cell_vec,ierr);CHKERRQ(ierr)
2546) call VecDestroy(global_z_cell_vec,ierr);CHKERRQ(ierr)
2547)
2548) call VecDestroy(natural_x_cell_vec,ierr);CHKERRQ(ierr)
2549) call VecDestroy(natural_y_cell_vec,ierr);CHKERRQ(ierr)
2550) call VecDestroy(natural_z_cell_vec,ierr);CHKERRQ(ierr)
2551)
2552) call UGridDMDestroy(ugdm_cell)
2553)
2554) #endif
2555) !if defined(SCORPIO_WRITE)
2556)
2557) end subroutine WriteHDF5CoordinatesUGridXDMF
2558)
2559) ! ************************************************************************** !
2560)
2561) subroutine DetermineNumVertices(realization_base,option)
2562) !
2563) ! Determine the number of vertices written out in the output HDF5 file
2564) !
2565) ! Author: Gautam Bisht, LBNL
2566) ! Date: 03/13/2015
2567) !
2568) use Realization_Base_class, only : realization_base_type
2569) use Grid_module
2570) use Option_module
2571) use Grid_Unstructured_Aux_module
2572) use Variables_module
2573)
2574) implicit none
2575)
2576) #include "petsc/finclude/petscvec.h"
2577) #include "petsc/finclude/petscvec.h90"
2578) #include "petsc/finclude/petsclog.h"
2579)
2580) class(realization_base_type) :: realization_base
2581) type(option_type), pointer :: option
2582)
2583) type(grid_type), pointer :: grid
2584) PetscInt :: local_size,vert_count
2585) PetscInt :: i
2586) PetscInt :: temp_int
2587)
2588) PetscReal, pointer :: vec_ptr(:)
2589) Vec :: global_vec, natural_vec
2590) type(ugdm_type),pointer :: ugdm_element
2591) PetscErrorCode :: ierr
2592)
2593) grid => realization_base%patch%grid
2594)
2595) call UGridCreateUGDM(grid%unstructured_grid,ugdm_element,EIGHT_INTEGER, &
2596) option)
2597) call UGridDMCreateVector(grid%unstructured_grid,ugdm_element,global_vec, &
2598) GLOBAL,option)
2599) call UGridDMCreateVector(grid%unstructured_grid,ugdm_element,natural_vec, &
2600) NATURAL,option)
2601) call GetCellConnections(grid,global_vec)
2602) call VecScatterBegin(ugdm_element%scatter_gton,global_vec,natural_vec, &
2603) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
2604) call VecScatterEnd(ugdm_element%scatter_gton,global_vec,natural_vec, &
2605) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
2606)
2607) local_size = grid%unstructured_grid%nlmax
2608)
2609) call VecGetArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
2610) vert_count=0
2611) do i=1,local_size*EIGHT_INTEGER
2612) if (int(vec_ptr(i)) >0 ) vert_count=vert_count+1
2613) enddo
2614) vert_count=vert_count+grid%nlmax
2615) call VecRestoreArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
2616)
2617) call MPI_Allreduce(vert_count,temp_int,ONE_INTEGER_MPI,MPIU_INTEGER, &
2618) MPI_SUM,option%mycomm,ierr)
2619) realization_base%output_option%xmf_vert_len=temp_int
2620)
2621) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
2622) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
2623) call UGridDMDestroy(ugdm_element)
2624)
2625) end subroutine DetermineNumVertices
2626)
2627) ! ************************************************************************** !
2628)
2629) subroutine WriteHDF5CoordinatesUGridXDMFExplicit(realization_base,option, &
2630) file_id)
2631) !
2632) ! Writes the coordinates of
2633) ! explicit grid to HDF5 file
2634) !
2635) ! Author: Satish Karra, LANL
2636) ! Date: 07/17/2013
2637) !
2638)
2639) use hdf5
2640) use Realization_Base_class, only : realization_base_type
2641) use Grid_module
2642) use Option_module
2643) use Grid_Unstructured_Aux_module
2644) use Variables_module
2645)
2646) implicit none
2647)
2648) #include "petsc/finclude/petscvec.h"
2649) #include "petsc/finclude/petscvec.h90"
2650) #include "petsc/finclude/petsclog.h"
2651)
2652) class(realization_base_type) :: realization_base
2653) type(option_type), pointer :: option
2654)
2655) #if defined(SCORPIO_WRITE)
2656) integer :: file_id
2657) integer :: data_type
2658) integer :: grp_id
2659) integer :: file_space_id
2660) integer :: memory_space_id
2661) integer :: data_set_id
2662) integer :: realization_set_id
2663) integer :: prop_id
2664) integer :: dims(3)
2665) integer :: start(3), length(3), stride(3)
2666) integer :: rank_mpi,file_space_rank_mpi
2667) integer :: hdf5_flag
2668) integer, parameter :: ON=1, OFF=0
2669) #else
2670) integer(HID_T) :: file_id
2671) integer(HID_T) :: data_type
2672) integer(HID_T) :: grp_id
2673) integer(HID_T) :: file_space_id
2674) integer(HID_T) :: realization_set_id
2675) integer(HID_T) :: memory_space_id
2676) integer(HID_T) :: data_set_id
2677) integer(HID_T) :: prop_id
2678) integer(HSIZE_T) :: dims(3)
2679) integer(HSIZE_T) :: start(3), length(3), stride(3)
2680) PetscMPIInt :: rank_mpi,file_space_rank_mpi
2681) PetscMPIInt :: hdf5_flag
2682) PetscMPIInt, parameter :: ON=1, OFF=0
2683) #endif
2684)
2685) type(grid_type), pointer :: grid
2686) character(len=MAXSTRINGLENGTH) :: string
2687) PetscMPIInt :: hdf5_err
2688)
2689) PetscInt :: local_size,vert_count,nverts
2690) PetscInt :: i,j
2691) PetscInt :: istart
2692) PetscReal, pointer :: vec_x_ptr(:),vec_y_ptr(:),vec_z_ptr(:)
2693) PetscReal, pointer :: double_array(:)
2694)
2695) PetscReal, pointer :: vec_ptr(:)
2696) Vec :: natural_vec
2697) PetscInt, pointer :: int_array(:)
2698) PetscErrorCode :: ierr
2699)
2700) PetscInt :: TRI_ID_XDMF = 4
2701) PetscInt :: QUAD_ID_XDMF = 5
2702) PetscInt :: TET_ID_XDMF = 6
2703) PetscInt :: PYR_ID_XDMF = 7
2704) PetscInt :: WED_ID_XDMF = 8
2705) PetscInt :: HEX_ID_XDMF = 9
2706)
2707) grid => realization_base%patch%grid
2708)
2709) allocate(vec_x_ptr(grid%unstructured_grid%explicit_grid%num_cells_global))
2710) allocate(vec_y_ptr(grid%unstructured_grid%explicit_grid%num_cells_global))
2711) allocate(vec_z_ptr(grid%unstructured_grid%explicit_grid%num_cells_global))
2712)
2713) do i = 1, grid%unstructured_grid%explicit_grid%num_cells_global
2714) vec_x_ptr(i) = grid%unstructured_grid%explicit_grid%vertex_coordinates(i)%x
2715) vec_y_ptr(i) = grid%unstructured_grid%explicit_grid%vertex_coordinates(i)%y
2716) vec_z_ptr(i) = grid%unstructured_grid%explicit_grid%vertex_coordinates(i)%z
2717) enddo
2718)
2719) #if defined(SCORPIO_WRITE)
2720) write(*,*) 'SCORPIO_WRITE'
2721) option%io_buffer = 'WriteHDF5CoordinatesUGrid not supported for &
2722) &SCORPIO_WRITE'
2723) call printErrMsg(option)
2724) #else
2725)
2726) !
2727) ! not(SCORPIO_WRITE)
2728) !
2729)
2730) local_size = grid%unstructured_grid%explicit_grid%num_cells_global
2731) ! memory space which is a 1D vector
2732) rank_mpi = 1
2733) dims = 0
2734) dims(1) = local_size * 3
2735) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
2736)
2737) ! file space which is a 2D block
2738) rank_mpi = 2
2739) dims = 0
2740) dims(2) = grid%unstructured_grid%explicit_grid%num_cells_global
2741) dims(1) = 3
2742) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
2743)
2744) string = "Vertices" // CHAR(0)
2745)
2746) call h5eset_auto_f(OFF,hdf5_err)
2747) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
2748) hdf5_flag = hdf5_err
2749) call h5eset_auto_f(ON,hdf5_err)
2750) if (hdf5_flag < 0) then
2751) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
2752) call h5dcreate_f(file_id,string,H5T_NATIVE_DOUBLE,file_space_id, &
2753) data_set_id,hdf5_err,prop_id)
2754) else
2755) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
2756) endif
2757)
2758) call h5pclose_f(prop_id,hdf5_err)
2759)
2760) istart = 0
2761) start(2) = istart
2762) start(1) = 0
2763)
2764) length(2) = local_size
2765) length(1) = 3
2766)
2767) stride = 1
2768) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
2769) hdf5_err,stride,stride)
2770) ! write the data
2771) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
2772)
2773)
2774) allocate(double_array(local_size*3))
2775) do i=1,local_size
2776) double_array((i-1)*3+1) = vec_x_ptr(i)
2777) double_array((i-1)*3+2) = vec_y_ptr(i)
2778) double_array((i-1)*3+3) = vec_z_ptr(i)
2779) enddo
2780)
2781) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2782) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,double_array,dims, &
2783) hdf5_err,memory_space_id,file_space_id,prop_id)
2784) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2785) deallocate(double_array)
2786)
2787) call h5pclose_f(prop_id,hdf5_err)
2788)
2789) call h5dclose_f(data_set_id,hdf5_err)
2790) call h5sclose_f(file_space_id,hdf5_err)
2791)
2792) deallocate(vec_x_ptr)
2793) deallocate(vec_y_ptr)
2794) deallocate(vec_z_ptr)
2795)
2796) !
2797) ! Write elements
2798) !
2799) local_size = grid%unstructured_grid%explicit_grid%num_elems
2800)
2801) call VecCreate(PETSC_COMM_SELF,natural_vec,ierr);CHKERRQ(ierr)
2802) call VecSetSizes(natural_vec,local_size*EIGHT_INTEGER, &
2803) PETSC_DECIDE,ierr);CHKERRQ(ierr)
2804) call VecSetBlockSize(natural_vec,EIGHT_INTEGER,ierr);CHKERRQ(ierr)
2805) call VecSetFromOptions(natural_vec,ierr);CHKERRQ(ierr)
2806)
2807) call GetCellConnectionsExplicit(grid,natural_vec)
2808) call VecGetArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
2809)
2810) vert_count=0
2811)
2812) do i=1,local_size*EIGHT_INTEGER
2813) if (int(vec_ptr(i)) >0 ) vert_count=vert_count+1
2814) enddo
2815) vert_count=vert_count+grid%unstructured_grid%explicit_grid%num_elems
2816)
2817) ! memory space which is a 1D vector
2818) rank_mpi = 1
2819) dims(1) = vert_count
2820) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
2821)
2822) realization_base%output_option%xmf_vert_len=int(dims(1))
2823)
2824) ! file space which is a 2D block
2825) rank_mpi = 1
2826) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
2827)
2828) string = "Cells" // CHAR(0)
2829)
2830) call h5eset_auto_f(OFF,hdf5_err)
2831) call h5dopen_f(file_id,string,data_set_id,hdf5_err)
2832) hdf5_flag = hdf5_err
2833) call h5eset_auto_f(ON,hdf5_err)
2834) if (hdf5_flag < 0) then
2835) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
2836) call h5dcreate_f(file_id,string,H5T_NATIVE_INTEGER,file_space_id, &
2837) data_set_id,hdf5_err,prop_id)
2838) else
2839) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
2840) endif
2841)
2842) call h5pclose_f(prop_id,hdf5_err)
2843)
2844) istart = 0
2845) start(1) = istart
2846) length(1) = vert_count
2847) stride = 1
2848) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
2849) hdf5_err,stride,stride)
2850)
2851) ! write the data
2852) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
2853)
2854) allocate(int_array(vert_count))
2855)
2856) vert_count=0
2857) do i=1,local_size
2858) nverts=0
2859) do j=1,8
2860) if (vec_ptr((i-1)*8+j)>0) nverts=nverts+1
2861) enddo
2862) vert_count=vert_count+1
2863) select case (nverts)
2864) case (3)
2865) int_array(vert_count) = TRI_ID_XDMF
2866) case (4)
2867) if (grid%unstructured_grid%grid_type /= TWO_DIM_GRID) then
2868) ! Tetrahedron
2869) int_array(vert_count) = TET_ID_XDMF
2870) else
2871) int_array(vert_count) = QUAD_ID_XDMF
2872) endif
2873) case (5) ! Pyramid
2874) int_array(vert_count) = PYR_ID_XDMF
2875) case (6) ! Wedge
2876) int_array(vert_count) = WED_ID_XDMF
2877) case (8) ! Hexahedron
2878) int_array(vert_count) = HEX_ID_XDMF
2879) end select
2880)
2881) do j=1,8
2882) if (vec_ptr((i-1)*8+j)>0) then
2883) vert_count=vert_count+1
2884) int_array(vert_count) = INT(vec_ptr((i-1)*8+j))-1
2885) endif
2886) enddo
2887) enddo
2888)
2889) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2890) call h5dwrite_f(data_set_id,H5T_NATIVE_INTEGER,int_array,dims, &
2891) hdf5_err,memory_space_id,file_space_id,prop_id)
2892) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
2893)
2894) deallocate(int_array)
2895)
2896) call h5pclose_f(prop_id,hdf5_err)
2897)
2898) call h5dclose_f(data_set_id,hdf5_err)
2899) call h5sclose_f(file_space_id,hdf5_err)
2900)
2901) call VecRestoreArrayF90(natural_vec,vec_ptr,ierr);CHKERRQ(ierr)
2902) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
2903)
2904) #endif
2905) !if defined(SCORPIO_WRITE)
2906)
2907) end subroutine WriteHDF5CoordinatesUGridXDMFExplicit
2908)
2909) ! ************************************************************************** !
2910)
2911) subroutine WriteHDF5FlowratesUGrid(realization_base,option,file_id, &
2912) var_list_type)
2913) !
2914) ! This routine writes (mass/energy) flowrate for unstructured grid.
2915) !
2916) ! Author: Gautam Bisht, LBNL
2917) ! Date: 03/19/2013
2918) !
2919)
2920) use hdf5
2921) use Realization_Base_class, only : realization_base_type
2922) use Patch_module
2923) use Grid_module
2924) use Option_module
2925) use Grid_Unstructured_Aux_module
2926) use Grid_Unstructured_Cell_module
2927) use Variables_module
2928) use Connection_module
2929) use Coupler_module
2930) use HDF5_Aux_module
2931) use Output_Aux_module
2932) use Field_module
2933)
2934) implicit none
2935)
2936) #include "petsc/finclude/petscvec.h"
2937) #include "petsc/finclude/petscvec.h90"
2938) #include "petsc/finclude/petsclog.h"
2939) #include "petsc/finclude/petscsys.h"
2940)
2941) class(realization_base_type) :: realization_base
2942) type(option_type), pointer :: option
2943) PetscInt :: var_list_type
2944)
2945) #if defined(SCORPIO_WRITE)
2946) integer :: file_id
2947) integer :: data_type
2948) integer :: grp_id
2949) integer :: file_space_id
2950) integer :: memory_space_id
2951) integer :: data_set_id
2952) integer :: realization_set_id
2953) integer :: prop_id
2954) integer :: dims(3)
2955) integer :: start(3), length(3), stride(3)
2956) integer :: rank_mpi,file_space_rank_mpi
2957) integer :: hdf5_flag
2958) integer, parameter :: ON=1, OFF=0
2959) #else
2960) integer(HID_T) :: file_id
2961) integer(HID_T) :: data_type
2962) integer(HID_T) :: grp_id
2963) integer(HID_T) :: file_space_id
2964) integer(HID_T) :: realization_set_id
2965) integer(HID_T) :: memory_space_id
2966) integer(HID_T) :: data_set_id
2967) integer(HID_T) :: prop_id
2968) integer(HSIZE_T) :: dims(3)
2969) integer(HSIZE_T) :: start(3), length(3), stride(3)
2970) PetscMPIInt :: rank_mpi,file_space_rank_mpi
2971) PetscMPIInt :: hdf5_flag
2972) PetscMPIInt, parameter :: ON=1, OFF=0
2973) #endif
2974)
2975) type(patch_type), pointer :: patch
2976) type(grid_type), pointer :: grid
2977) type(grid_unstructured_type),pointer :: ugrid
2978) type(connection_set_list_type), pointer :: connection_set_list
2979) type(connection_set_type), pointer :: cur_connection_set
2980) type(coupler_type), pointer :: boundary_condition
2981) type(ugdm_type),pointer :: ugdm
2982) type(output_option_type), pointer :: output_option
2983) type(field_type), pointer :: field
2984)
2985) PetscInt :: local_id
2986) PetscInt :: ghosted_id
2987) PetscInt :: idual
2988) PetscInt :: iconn
2989) PetscInt :: face_id
2990) PetscInt :: local_id_up,local_id_dn
2991) PetscInt :: ghosted_id_up,ghosted_id_dn
2992) PetscInt :: iface_up,iface_dn
2993) PetscInt :: dof
2994) PetscInt :: sum_connection
2995) PetscInt :: offset
2996) PetscInt :: cell_type
2997) PetscInt :: local_size
2998) PetscInt :: i
2999) PetscInt :: iface
3000) PetscInt :: ndof
3001)
3002) PetscReal, pointer :: flowrates(:,:,:)
3003) PetscReal, pointer :: vec_ptr1(:)
3004) PetscReal, pointer :: vec_ptr2(:)
3005) PetscReal, pointer :: double_array(:)
3006) PetscReal :: dtime
3007) PetscInt :: istart
3008)
3009) PetscBool :: mass_flowrate
3010) PetscBool :: energy_flowrate
3011)
3012) Vec :: global_flowrates_vec
3013) Vec :: natural_flowrates_vec
3014)
3015) PetscMPIInt :: hdf5_err
3016) PetscErrorCode :: ierr
3017)
3018) character(len=MAXSTRINGLENGTH) :: string
3019) character(len=MAXWORDLENGTH) :: unit_string
3020)
3021) patch => realization_base%patch
3022) grid => patch%grid
3023) ugrid => grid%unstructured_grid
3024) output_option =>realization_base%output_option
3025) field => realization_base%field
3026)
3027) #if defined(SCORPIO_WRITE)
3028) write(*,*) 'SCORPIO_WRITE'
3029) option%io_buffer = 'WriteHDF5FlowratesUGrid not supported for SCORPIO_WRITE'
3030) call printErrMsg(option)
3031) #else
3032) select case(option%iflowmode)
3033) case (RICHARDS_MODE)
3034) ndof=1
3035) case (TH_MODE)
3036) ndof=1
3037) if (output_option%print_hdf5_mass_flowrate .and. &
3038) output_option%print_hdf5_energy_flowrate) ndof = 2
3039) case default
3040) option%io_buffer='FLOWRATE output not supported in this mode'
3041) call printErrMsg(option)
3042) end select
3043)
3044) call VecGetLocalSize(field%flowrate_inst,local_size,ierr);CHKERRQ(ierr)
3045) local_size = local_size/(option%nflowdof*MAX_FACE_PER_CELL + 1)
3046)
3047) allocate(double_array(local_size*(MAX_FACE_PER_CELL+1)))
3048) double_array = 0.d0
3049)
3050) offset = option%nflowdof*MAX_FACE_PER_CELL+1
3051)
3052) select case (var_list_type)
3053) case (INSTANTANEOUS_VARS)
3054) call VecGetArrayF90(field%flowrate_inst,vec_ptr1,ierr);CHKERRQ(ierr)
3055) mass_flowrate = output_option%print_hdf5_mass_flowrate
3056) energy_flowrate = output_option%print_hdf5_energy_flowrate
3057) case (AVERAGED_VARS)
3058) call VecGetArrayF90(field%flowrate_inst,vec_ptr1,ierr);CHKERRQ(ierr)
3059) call VecGetArrayF90(field%flowrate_aveg,vec_ptr2,ierr);CHKERRQ(ierr)
3060) mass_flowrate = output_option%print_hdf5_aveg_mass_flowrate
3061) energy_flowrate = output_option%print_hdf5_aveg_energy_flowrate
3062) end select
3063)
3064)
3065) do dof = 1,option%nflowdof
3066)
3067) if (dof==1 .and. (.not.mass_flowrate)) exit
3068) if (dof==2 .and. (.not.energy_flowrate)) exit
3069)
3070) select case(option%iflowmode)
3071) case(RICHARDS_MODE)
3072) string = "Mass_Flowrate [kg_per_s]" // CHAR(0)
3073) case(TH_MODE)
3074) if (dof==1) then
3075) string = "Mass_Flowrate [kg_per_s]" // CHAR(0)
3076) else
3077) string = "Energy_Flowrate [MJ_per_s]" // CHAR(0)
3078) endif
3079) case default
3080) option%io_buffer='FLOWRATE output not implemented in this mode.'
3081) call printErrMsg(option)
3082) end select
3083)
3084) if (var_list_type==AVERAGED_VARS) string = 'Aveg_' // trim(string) // &
3085) char(0)
3086)
3087) ! memory space which is a 1D vector
3088) rank_mpi = 1
3089) dims = 0
3090) dims(1) = local_size*(MAX_FACE_PER_CELL+1)
3091) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
3092)
3093) ! file space which is a 2D block
3094) rank_mpi = 2
3095) dims = 0
3096) dims(2) = ugrid%nmax
3097) dims(1) = MAX_FACE_PER_CELL + 1
3098) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
3099)
3100) call h5eset_auto_f(OFF,hdf5_err)
3101) call h5dopen_f(file_id,trim(string),data_set_id,hdf5_err)
3102) hdf5_flag = hdf5_err
3103) call h5eset_auto_f(ON,hdf5_err)
3104) if (hdf5_flag < 0) then
3105) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
3106) call h5dcreate_f(file_id,trim(string),H5T_NATIVE_DOUBLE,file_space_id, &
3107) data_set_id,hdf5_err,prop_id)
3108) else
3109) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
3110) endif
3111)
3112) call h5pclose_f(prop_id,hdf5_err)
3113)
3114) istart = 0
3115) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
3116) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
3117)
3118) start(2) = istart
3119) start(1) = 0
3120)
3121) length(2) = local_size
3122) length(1) = MAX_FACE_PER_CELL + 1
3123)
3124) stride = 1
3125) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
3126) hdf5_err,stride,stride)
3127) ! write the data
3128) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
3129) #ifndef SERIAL_HDF5
3130) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
3131) hdf5_err)
3132) #endif
3133)
3134) select case (var_list_type)
3135) case (INSTANTANEOUS_VARS)
3136)
3137) do i=1,local_size
3138) ! Num. of faces for each cell (Note: Use vec_ptr1 not vec_ptr2)
3139) double_array((i-1)*(MAX_FACE_PER_CELL+1)+1) = &
3140) vec_ptr1((i-1)*offset+1)
3141) ! Flowrate values for each face
3142) do iface = 1,MAX_FACE_PER_CELL
3143) double_array((i-1)*(MAX_FACE_PER_CELL+1)+iface+1) = &
3144) vec_ptr1((i-1)*offset + (dof-1)*MAX_FACE_PER_CELL + iface + 1)
3145) enddo
3146) enddo
3147)
3148) case (AVERAGED_VARS)
3149)
3150) do i=1,local_size
3151) ! Num. of faces for each cell (Note: Use vec_ptr1 not vec_ptr2)
3152) double_array((i-1)*(MAX_FACE_PER_CELL+1)+1) = &
3153) vec_ptr1((i-1)*offset+1)
3154) ! Divide the flowrate values by integration 'time'
3155) do iface = 1,MAX_FACE_PER_CELL
3156) double_array((i-1)*(MAX_FACE_PER_CELL+1)+iface+1) = &
3157) vec_ptr2((i-1)*offset + (dof-1)*MAX_FACE_PER_CELL + iface + 1)/ &
3158) output_option%periodic_snap_output_time_incr
3159) enddo
3160) enddo
3161) end select
3162)
3163) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
3164) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,double_array,dims, &
3165) hdf5_err,memory_space_id,file_space_id,prop_id)
3166) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
3167)
3168) call h5pclose_f(prop_id,hdf5_err)
3169)
3170) call h5dclose_f(data_set_id,hdf5_err)
3171) call h5sclose_f(file_space_id,hdf5_err)
3172)
3173) enddo
3174)
3175) ! Free up memory
3176) deallocate(double_array)
3177) #endif
3178) ! #ifdef SCORPIO_WRITE
3179)
3180) end subroutine WriteHDF5FlowratesUGrid
3181)
3182) ! ************************************************************************** !
3183)
3184) subroutine WriteHDF5FaceVelUGrid(realization_base,option,file_id, &
3185) var_list_type)
3186) !
3187) ! This routine writes velocity at cell faces for unstructured grid.
3188) !
3189) ! Author: Gautam Bisht, LBNL
3190) ! Date: 05/25/2014
3191) !
3192)
3193) use hdf5
3194) use Realization_Base_class, only : realization_base_type
3195) use Patch_module
3196) use Grid_module
3197) use Option_module
3198) use Grid_Unstructured_Aux_module
3199) use Grid_Unstructured_Cell_module
3200) use Variables_module
3201) use Connection_module
3202) use Coupler_module
3203) use HDF5_Aux_module
3204) use Output_Aux_module
3205) use Field_module
3206)
3207) implicit none
3208)
3209) #include "petsc/finclude/petscvec.h"
3210) #include "petsc/finclude/petscvec.h90"
3211) #include "petsc/finclude/petsclog.h"
3212) #include "petsc/finclude/petscsys.h"
3213)
3214) class(realization_base_type) :: realization_base
3215) type(option_type), pointer :: option
3216) PetscInt :: var_list_type
3217)
3218) #if defined(SCORPIO_WRITE)
3219) integer :: file_id
3220) integer :: data_type
3221) integer :: grp_id
3222) integer :: file_space_id
3223) integer :: memory_space_id
3224) integer :: data_set_id
3225) integer :: realization_set_id
3226) integer :: prop_id
3227) integer :: dims(3)
3228) integer :: start(3), length(3), stride(3)
3229) integer :: rank_mpi,file_space_rank_mpi
3230) integer :: hdf5_flag
3231) integer, parameter :: ON=1, OFF=0
3232) #else
3233) integer(HID_T) :: file_id
3234) integer(HID_T) :: data_type
3235) integer(HID_T) :: grp_id
3236) integer(HID_T) :: file_space_id
3237) integer(HID_T) :: realization_set_id
3238) integer(HID_T) :: memory_space_id
3239) integer(HID_T) :: data_set_id
3240) integer(HID_T) :: prop_id
3241) integer(HSIZE_T) :: dims(3)
3242) integer(HSIZE_T) :: start(3), length(3), stride(3)
3243) PetscMPIInt :: rank_mpi,file_space_rank_mpi
3244) PetscMPIInt :: hdf5_flag
3245) PetscMPIInt, parameter :: ON=1, OFF=0
3246) #endif
3247)
3248) type(patch_type), pointer :: patch
3249) type(grid_type), pointer :: grid
3250) type(grid_unstructured_type),pointer :: ugrid
3251) type(connection_set_list_type), pointer :: connection_set_list
3252) type(connection_set_type), pointer :: cur_connection_set
3253) type(coupler_type), pointer :: boundary_condition
3254) type(ugdm_type),pointer :: ugdm
3255) type(output_option_type), pointer :: output_option
3256) type(field_type), pointer :: field
3257)
3258) PetscInt :: local_id
3259) PetscInt :: ghosted_id
3260) PetscInt :: idual
3261) PetscInt :: iconn
3262) PetscInt :: face_id
3263) PetscInt :: local_id_up,local_id_dn
3264) PetscInt :: ghosted_id_up,ghosted_id_dn
3265) PetscInt :: iface_up,iface_dn
3266) PetscInt :: iphase
3267) PetscInt :: sum_connection
3268) PetscInt :: offset
3269) PetscInt :: cell_type
3270) PetscInt :: local_size
3271) PetscInt :: i
3272) PetscInt :: idir
3273) PetscInt :: istart
3274) PetscInt :: iface
3275)
3276) PetscReal, pointer :: flowrates(:,:,:)
3277) PetscReal, pointer :: vx_ptr(:)
3278) PetscReal, pointer :: vy_ptr(:)
3279) PetscReal, pointer :: vz_ptr(:)
3280) PetscReal, pointer :: vec_ptr1(:)
3281) PetscReal, pointer :: double_array(:)
3282) PetscReal :: dtime
3283)
3284) Vec :: global_flowrates_vec
3285) Vec :: natural_flowrates_vec
3286)
3287) PetscMPIInt :: hdf5_err
3288) PetscErrorCode :: ierr
3289)
3290) character(len=MAXSTRINGLENGTH) :: string
3291) character(len=MAXWORDLENGTH) :: unit_string
3292) character(len=1) :: string_dir
3293)
3294) patch => realization_base%patch
3295) grid => patch%grid
3296) ugrid => grid%unstructured_grid
3297) output_option =>realization_base%output_option
3298) field => realization_base%field
3299)
3300) #if defined(SCORPIO_WRITE)
3301) write(*,*) 'SCORPIO_WRITE'
3302) option%io_buffer = 'WriteHDF5FaceVelUGrid not supported for SCORPIO_WRITE'
3303) call printErrMsg(option)
3304) #else
3305)
3306) call VecGetLocalSize(field%vx_face_inst,local_size,ierr);CHKERRQ(ierr)
3307) local_size = local_size/(option%nphase*MAX_FACE_PER_CELL + 1)
3308)
3309) allocate(double_array(local_size*(MAX_FACE_PER_CELL+1)))
3310) double_array = 0.d0
3311)
3312) offset = option%nphase*MAX_FACE_PER_CELL+1
3313)
3314) do idir = 1,3
3315)
3316) select case (idir)
3317) case (1)
3318) string_dir = 'X'
3319) call VecGetArrayF90(field%vx_face_inst,vec_ptr1,ierr);CHKERRQ(ierr)
3320) case (2)
3321) string_dir = 'Y'
3322) call VecGetArrayF90(field%vy_face_inst,vec_ptr1,ierr);CHKERRQ(ierr)
3323) case (3)
3324) string_dir = 'Z'
3325) call VecGetArrayF90(field%vz_face_inst,vec_ptr1,ierr);CHKERRQ(ierr)
3326) end select
3327)
3328) do iphase = 1,option%nphase
3329)
3330) select case (iphase)
3331) case (LIQUID_PHASE)
3332) string = "Liquid " // string_dir // "-Velocity at cell face [m_per_" &
3333) // trim(output_option%tunit) // "]"
3334) case (GAS_PHASE)
3335) string = "Gas " // string_dir // "-Velocity at cell face [m_per_" // &
3336) trim(output_option%tunit) // "]"
3337) end select
3338)
3339) ! memory space which is a 1D vector
3340) rank_mpi = 1
3341) dims = 0
3342) dims(1) = local_size*(MAX_FACE_PER_CELL+1)
3343) call h5screate_simple_f(rank_mpi,dims,memory_space_id,hdf5_err,dims)
3344)
3345) ! file space which is a 2D block
3346) rank_mpi = 2
3347) dims = 0
3348) dims(2) = ugrid%nmax
3349) dims(1) = MAX_FACE_PER_CELL + 1
3350) call h5pcreate_f(H5P_DATASET_CREATE_F,prop_id,hdf5_err)
3351)
3352) call h5eset_auto_f(OFF,hdf5_err)
3353) call h5dopen_f(file_id,trim(string),data_set_id,hdf5_err)
3354) hdf5_flag = hdf5_err
3355) call h5eset_auto_f(ON,hdf5_err)
3356) if (hdf5_flag < 0) then
3357) call h5screate_simple_f(rank_mpi,dims,file_space_id,hdf5_err,dims)
3358) call h5dcreate_f(file_id,trim(string),H5T_NATIVE_DOUBLE,file_space_id, &
3359) data_set_id,hdf5_err,prop_id)
3360) else
3361) call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
3362) endif
3363)
3364) call h5pclose_f(prop_id,hdf5_err)
3365)
3366) istart = 0
3367) call MPI_Exscan(local_size, istart, ONE_INTEGER_MPI, &
3368) MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
3369)
3370) start(2) = istart
3371) start(1) = 0
3372)
3373) length(2) = local_size
3374) length(1) = MAX_FACE_PER_CELL + 1
3375)
3376) stride = 1
3377) call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,start,length, &
3378) hdf5_err,stride,stride)
3379) ! write the data
3380) call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
3381) #ifndef SERIAL_HDF5
3382) call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F, &
3383) hdf5_err)
3384) #endif
3385)
3386) select case (var_list_type)
3387) case (INSTANTANEOUS_VARS)
3388)
3389) do i=1,local_size
3390) ! Num. of faces for each cell
3391) double_array((i-1)*(MAX_FACE_PER_CELL+1)+1) = &
3392) vec_ptr1((i-1)*offset+1)
3393) ! Flowrate values for each face
3394) do iface = 1,MAX_FACE_PER_CELL
3395) double_array((i-1)*(MAX_FACE_PER_CELL+1)+iface+1) = &
3396) vec_ptr1((i-1)*offset + (iphase-1)*MAX_FACE_PER_CELL + iface + 1)* &
3397) realization_base%output_option%tconv
3398) enddo
3399) enddo
3400)
3401) case (AVERAGED_VARS)
3402)
3403) end select
3404)
3405) call PetscLogEventBegin(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
3406) call h5dwrite_f(data_set_id,H5T_NATIVE_DOUBLE,double_array,dims, &
3407) hdf5_err,memory_space_id,file_space_id,prop_id)
3408) call PetscLogEventEnd(logging%event_h5dwrite_f,ierr);CHKERRQ(ierr)
3409)
3410) call h5pclose_f(prop_id,hdf5_err)
3411)
3412) call h5dclose_f(data_set_id,hdf5_err)
3413) call h5sclose_f(file_space_id,hdf5_err)
3414)
3415) enddo
3416)
3417) select case (idir)
3418) case (1)
3419) call VecRestoreArrayF90(field%vx_face_inst,vec_ptr1, &
3420) ierr);CHKERRQ(ierr)
3421) case (2)
3422) call VecRestoreArrayF90(field%vy_face_inst,vec_ptr1, &
3423) ierr);CHKERRQ(ierr)
3424) case (3)
3425) call VecRestoreArrayF90(field%vz_face_inst,vec_ptr1, &
3426) ierr);CHKERRQ(ierr)
3427) end select
3428)
3429) enddo
3430)
3431) ! Free up memory
3432) deallocate(double_array)
3433) #endif
3434) ! #ifdef SCORPIO_WRITE
3435)
3436) end subroutine WriteHDF5FaceVelUGrid
3437)
3438) ! ************************************************************************** !
3439)
3440) subroutine OutputHDF5Provenance(option, output_option, file_id)
3441) !
3442) ! write pflotran and petsc provenance information including a copy
3443) ! of the inputfile
3444) !
3445)
3446) use Option_module, only : option_type
3447) use Output_Aux_module, only : output_option_type
3448) use PFLOTRAN_Provenance_module, only : provenance_max_str_len
3449)
3450) #include "petsc/finclude/petscsysdef.h"
3451)
3452) use hdf5
3453)
3454) implicit none
3455)
3456) type(option_type), intent(in) :: option
3457) type(output_option_type), intent(in) :: output_option
3458) integer(HID_T), intent(in) :: file_id
3459)
3460) character(len=32) :: filename, name
3461) integer(HID_T) :: prop_id, provenance_id, string_type
3462) PetscMPIInt :: hdf5_err
3463) PetscBool :: first
3464) integer(SIZE_T) :: size_t_int
3465)
3466) ! create the provenance group
3467) name = "Provenance"
3468) call h5gcreate_f(file_id, name, provenance_id, hdf5_err, &
3469) OBJECT_NAMELEN_DEFAULT_F)
3470)
3471) ! create fixed length string datatype
3472) call h5tcopy_f(H5T_FORTRAN_S1, string_type, hdf5_err)
3473) size_t_int = provenance_max_str_len
3474) call h5tset_size_f(string_type, size_t_int, hdf5_err)
3475)
3476) call OutputHDF5Provenance_PFLOTRAN(option, provenance_id, string_type)
3477) call OutputHDF5Provenance_PETSc(provenance_id, string_type)
3478)
3479) ! close the provenance group
3480) call h5tclose_f(string_type, hdf5_err)
3481) call h5gclose_f(provenance_id, hdf5_err)
3482)
3483) end subroutine OutputHDF5Provenance
3484)
3485) ! ************************************************************************** !
3486)
3487) subroutine OutputHDF5Provenance_PFLOTRAN(option, provenance_id, string_type)
3488) !
3489) ! write the pflotran provenance data as attributes (small) or
3490) ! datasets (big details)
3491) !
3492)
3493) use Option_module, only : option_type
3494) use PFLOTRAN_Provenance_module
3495)
3496) use hdf5
3497)
3498) implicit none
3499)
3500) type(option_type), intent(in) :: option
3501) integer(HID_T), intent(in) :: provenance_id
3502) integer(HID_T), intent(in) :: string_type
3503)
3504) character(len=32) :: name
3505) integer(HID_T) :: pflotran_id
3506) PetscMPIInt :: hdf5_err
3507)
3508) ! Create the pflotran group under provenance
3509) name = "PFLOTRAN"
3510) call h5gcreate_f(provenance_id, name, pflotran_id, hdf5_err, &
3511) OBJECT_NAMELEN_DEFAULT_F)
3512)
3513) call OutputHDF5DatasetStringArray(pflotran_id, string_type, &
3514) "pflotran_compile_date_time", &
3515) ONE_INTEGER, pflotran_compile_date_time)
3516)
3517) call OutputHDF5DatasetStringArray(pflotran_id, string_type, &
3518) "pflotran_compile_user", &
3519) ONE_INTEGER, pflotran_compile_user)
3520)
3521) call OutputHDF5DatasetStringArray(pflotran_id, string_type, &
3522) "pflotran_compile_hostname", &
3523) ONE_INTEGER, pflotran_compile_hostname)
3524)
3525) call OutputHDF5AttributeStringArray(pflotran_id, string_type, &
3526) "pflotran_status", &
3527) ONE_INTEGER, pflotran_status)
3528)
3529) call OutputHDF5AttributeStringArray(pflotran_id, string_type, &
3530) "pflotran_changeset", &
3531) ONE_INTEGER, pflotran_changeset)
3532)
3533) call OutputHDF5DatasetStringArray(pflotran_id, string_type, &
3534) "detail_pflotran_fflags", &
3535) detail_pflotran_fflags_len, &
3536) detail_pflotran_fflags)
3537)
3538) call OutputHDF5DatasetStringArray(pflotran_id, string_type, &
3539) "detail_pflotran_status", &
3540) detail_pflotran_status_len, &
3541) detail_pflotran_status)
3542)
3543) call OutputHDF5DatasetStringArray(pflotran_id, string_type, &
3544) "detail_pflotran_parent", &
3545) detail_pflotran_parent_len, &
3546) detail_pflotran_parent)
3547)
3548) ! FIXME(bja, 2013-11-25): break gcc when diffs are present
3549) call OutputHDF5DatasetStringArray(pflotran_id, string_type, &
3550) "detail_pflotran_diff", &
3551) detail_pflotran_diff_len, &
3552) detail_pflotran_diff)
3553)
3554) call OutputHDF5Provenance_input(option, pflotran_id)
3555)
3556) ! close pflotran group
3557) call h5gclose_f(pflotran_id, hdf5_err)
3558)
3559) end subroutine OutputHDF5Provenance_PFLOTRAN
3560)
3561) ! ************************************************************************** !
3562)
3563) subroutine OutputHDF5Provenance_input(option, pflotran_id)
3564) !
3565) ! open the pflotran input file, figure out how long it is, read it
3566) ! into a buffer, then write the buffer as a pflotran provenance
3567) ! group dataset.
3568) !
3569) use hdf5
3570) use Input_Aux_module, only : input_type, InputCreate, InputDestroy, &
3571) InputGetLineCount, InputReadToBuffer
3572) use Option_module, only : option_type
3573) use PFLOTRAN_Constants_module, only : IN_UNIT, MAXSTRINGLENGTH
3574)
3575) implicit none
3576)
3577) type(option_type), intent(in) :: option
3578) integer(HID_T), intent(in) :: pflotran_id
3579)
3580) integer(HID_T) :: input_string_type
3581) type(input_type), pointer :: input
3582) PetscInt :: i, input_line_count
3583) character(len=MAXSTRINGLENGTH), allocatable :: input_buffer(:)
3584) PetscMPIInt :: hdf5_err
3585) integer(SIZE_T) :: size_t_int
3586)
3587) input => InputCreate(IN_UNIT, option%input_filename, option)
3588) input_line_count = InputGetLineCount(input)
3589) allocate(input_buffer(input_line_count))
3590) call InputReadToBuffer(input, input_buffer)
3591) call h5tcopy_f(H5T_FORTRAN_S1, input_string_type, hdf5_err)
3592) size_t_int = MAXWORDLENGTH
3593) call h5tset_size_f(input_string_type, size_t_int, hdf5_err)
3594) call OutputHDF5DatasetStringArray(pflotran_id, input_string_type, &
3595) "pflotran_input_file", &
3596) input_line_count, input_buffer)
3597) call h5tclose_f(input_string_type, hdf5_err)
3598) deallocate(input_buffer)
3599) call InputDestroy(input)
3600)
3601) end subroutine OutputHDF5Provenance_input
3602)
3603) ! ************************************************************************** !
3604)
3605) subroutine OutputHDF5Provenance_PETSc(provenance_id, string_type)
3606) !
3607) ! write the petsc provenance data as attributes (small) or datasets
3608) ! (big details)
3609) !
3610)
3611) use PFLOTRAN_Provenance_module
3612) use hdf5
3613)
3614) implicit none
3615)
3616) integer(HID_T), intent(in) :: provenance_id
3617) integer(HID_T), intent(in) :: string_type
3618)
3619) character(len=32) :: name
3620) integer(HID_T) :: petsc_id
3621) PetscMPIInt :: hdf5_err
3622)
3623) ! create the petsc group under provenance
3624) name = "PETSc"
3625) call h5gcreate_f(provenance_id, name, petsc_id, hdf5_err, &
3626) OBJECT_NAMELEN_DEFAULT_F)
3627)
3628) call OutputHDF5AttributeStringArray(petsc_id, string_type, "petsc_status", &
3629) ONE_INTEGER, petsc_status)
3630)
3631) call OutputHDF5AttributeStringArray(petsc_id, string_type, &
3632) "petsc_changeset", &
3633) ONE_INTEGER, petsc_changeset)
3634)
3635) call OutputHDF5DatasetStringArray(petsc_id, string_type, &
3636) "detail_petsc_status", &
3637) detail_petsc_status_len, &
3638) detail_petsc_status)
3639)
3640) call OutputHDF5DatasetStringArray(petsc_id, string_type, &
3641) "detail_petsc_parent", &
3642) detail_petsc_parent_len, &
3643) detail_petsc_parent)
3644)
3645) call OutputHDF5DatasetStringArray(petsc_id, string_type, &
3646) "detail_petsc_config", &
3647) detail_petsc_config_len, &
3648) detail_petsc_config)
3649)
3650) ! close the petsc group
3651) call h5gclose_f(petsc_id, hdf5_err)
3652)
3653) end subroutine OutputHDF5Provenance_PETSc
3654)
3655) ! ************************************************************************** !
3656)
3657) subroutine OutputHDF5AttributeStringArray(parent_id, type, name, length, data)
3658) ! create the dataspaces and attributes consisting of an array of
3659) ! strings, then write the data and cleanup
3660)
3661) use hdf5
3662)
3663) implicit none
3664)
3665) integer(HID_T), intent(in) :: parent_id, type
3666) character(len=*), intent(in) :: name
3667) PetscInt, intent(in) :: length
3668) character(len=*), intent(in) :: data(length)
3669)
3670) integer(HID_T) :: dataspace_id, attribute_id
3671) integer(HSIZE_T), dimension(1:1) :: dims
3672) PetscMPIInt :: hdf5_err
3673)
3674) dims = length
3675) call h5screate_simple_f(1, dims, dataspace_id, hdf5_err)
3676) call h5acreate_f(parent_id, name, type, dataspace_id, attribute_id, hdf5_err)
3677) call h5awrite_f(attribute_id, type, data, dims, hdf5_err)
3678) call h5aclose_f(attribute_id, hdf5_err)
3679) call h5sclose_f(dataspace_id, hdf5_err)
3680)
3681) end subroutine OutputHDF5AttributeStringArray
3682)
3683) ! ************************************************************************** !
3684)
3685) subroutine OutputHDF5DatasetStringArray(parent_id, type, name, length, data)
3686) ! create the dataspaces and dataset consisting of an array of
3687) ! strings, then write the data and cleanup
3688)
3689) use hdf5
3690)
3691) implicit none
3692)
3693) integer(HID_T), intent(in) :: parent_id, type
3694) character(len=*), intent(in) :: name
3695) PetscInt, intent(in) :: length
3696) character(len=*), intent(in) :: data(length)
3697)
3698) integer(HID_T) :: dataspace_id, attribute_id
3699) integer(HSIZE_T), dimension(1:1) :: dims
3700) PetscMPIInt :: hdf5_err
3701)
3702) dims = length
3703) call h5screate_simple_f(1, dims, dataspace_id, hdf5_err)
3704) call h5dcreate_f(parent_id, name, type, dataspace_id, attribute_id, hdf5_err)
3705) call h5dwrite_f(attribute_id, type, data, dims, hdf5_err)
3706) call h5dclose_f(attribute_id, hdf5_err)
3707) call h5sclose_f(dataspace_id, hdf5_err)
3708)
3709) end subroutine OutputHDF5DatasetStringArray
3710)
3711) ! PETSC_HAVE_HDF5
3712) #endif
3713)
3714) end module Output_HDF5_module