hdf5_aux.F90       coverage:  66.67 %func     58.70 %block


     1) module HDF5_Aux_module
     2) 
     3) #if defined(PETSC_HAVE_HDF5)
     4)   use hdf5
     5) #endif
     6)   use Logging_module
     7) 
     8)   use PFLOTRAN_Constants_module
     9) 
    10)   implicit none
    11) 
    12) #include "petsc/finclude/petscsys.h"
    13) 
    14)   private
    15) 
    16)   PetscInt, parameter, public :: HDF5_READ_BUFFER_SIZE = 1000000
    17) !#define HDF5_BROADCAST
    18) 
    19)   PetscErrorCode :: ierr
    20) 
    21) #if defined(PETSC_HAVE_HDF5)
    22)   PetscMPIInt :: hdf5_err
    23)   PetscMPIInt :: io_rank_mpi
    24) ! 64-bit stuff
    25) #ifdef PETSC_USE_64BIT_INDICES
    26) !#define HDF_NATIVE_INTEGER H5T_STD_I64LE  
    27) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
    28) #else
    29) #define HDF_NATIVE_INTEGER H5T_NATIVE_INTEGER
    30) #endif
    31) 
    32)   public :: HDF5ReadNDimRealArray, &
    33) #ifdef SCORPIO
    34)             HDF5ReadDatasetInteger2D, &
    35)             HDF5ReadDatasetReal2D, &
    36)             HDF5GroupExists, &
    37)             HDF5DatasetExists, &
    38) #else
    39)             HDF5GroupExists, &
    40)             HDF5DatasetExists, &
    41) #endif
    42) ! SCORPIO
    43)             HDF5MakeStringCompatible, &
    44)             HDF5ReadDbase, &
    45)             HDF5OpenFileReadOnly
    46) 
    47) contains
    48) 
    49) #endif
    50)   
    51) 
    52) #if defined(PETSC_HAVE_HDF5)
    53) 
    54) ! ************************************************************************** !
    55) 
    56) subroutine HDF5ReadNDimRealArray(option,file_id,dataset_name,ndims,dims, &
    57)                                  real_array)
    58)   ! 
    59)   ! Read in an n-dimensional array from an hdf5 file
    60)   ! 
    61)   ! Author: Glenn Hammond
    62)   ! Date: 01/13/10
    63)   ! 
    64) 
    65)   use hdf5
    66)   
    67)   use Option_module
    68)   
    69)   implicit none
    70)   
    71)   type(option_type) :: option
    72)   character(len=MAXWORDLENGTH) :: dataset_name
    73)   integer(HID_T) :: file_id
    74)   PetscInt :: ndims
    75)   PetscInt, pointer :: dims(:)
    76)   PetscReal, pointer :: real_array(:)
    77)   
    78)   integer(HID_T) :: file_space_id
    79)   integer(HID_T) :: memory_space_id
    80)   integer(HID_T) :: data_set_id
    81)   integer(HID_T) :: prop_id
    82)   integer(HID_T) :: ndims_hdf5
    83)   integer(HSIZE_T), allocatable :: dims_h5(:), max_dims_h5(:)
    84)   integer(HSIZE_T) :: offset(1), length(1), stride(1)
    85)   PetscMPIInt :: rank_mpi
    86)   PetscInt :: index_count
    87)   integer(HSIZE_T) :: num_reals_in_dataset
    88)   PetscInt :: temp_int, i, index
    89)   PetscMPIInt :: int_mpi
    90)   
    91)   call PetscLogEventBegin(logging%event_read_ndim_real_array_hdf5, &
    92)                           ierr);CHKERRQ(ierr)
    93)                           
    94)   call h5dopen_f(file_id,dataset_name,data_set_id,hdf5_err)
    95)   call h5dget_space_f(data_set_id,file_space_id,hdf5_err)
    96)   ! should be a rank=1 data space
    97)   
    98)   call h5sget_simple_extent_ndims_f(file_space_id,ndims_hdf5,hdf5_err)
    99)   ndims = ndims_hdf5
   100)   allocate(dims_h5(ndims))
   101)   allocate(max_dims_h5(ndims))
   102)   allocate(dims(ndims))
   103)   call h5sget_simple_extent_dims_f(file_space_id,dims_h5,max_dims_h5,hdf5_err)
   104)   dims = int(dims_h5)
   105)   call h5sget_simple_extent_npoints_f(file_space_id,num_reals_in_dataset,hdf5_err)
   106)   temp_int = dims(1)
   107)   do i = 2, ndims
   108)     temp_int = temp_int * dims(i)
   109)   enddo
   110) 
   111)   rank_mpi = 1
   112)   offset = 0
   113)   length = num_reals_in_dataset
   114)   stride = 1
   115)   
   116)   call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
   117) #ifndef SERIAL_HDF5
   118)   call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F,hdf5_err)
   119) #endif
   120)   call h5screate_simple_f(rank_mpi,length,memory_space_id,hdf5_err,length)
   121) 
   122)   allocate(real_array(num_reals_in_dataset))
   123)   real_array = 0.d0
   124) #ifdef HDF5_BROADCAST
   125)   if (option%myrank == option%io_rank) then                           
   126) #endif
   127)     call PetscLogEventBegin(logging%event_h5dread_f,ierr);CHKERRQ(ierr)
   128)     call h5dread_f(data_set_id,H5T_NATIVE_DOUBLE,real_array,length, &
   129)                    hdf5_err,memory_space_id,file_space_id,prop_id)
   130)     call PetscLogEventEnd(logging%event_h5dread_f,ierr);CHKERRQ(ierr)
   131) #ifdef HDF5_BROADCAST
   132)   endif
   133)   if (option%mycommsize > 1) then
   134)     int_mpi = num_reals_in_dataset
   135)     call MPI_Bcast(real_array,int_mpi,MPI_DOUBLE_PRECISION, &
   136)                    option%io_rank,option%mycomm,ierr)
   137)   endif
   138) #endif
   139) 
   140)   call h5pclose_f(prop_id,hdf5_err)
   141)   if (memory_space_id > -1) call h5sclose_f(memory_space_id,hdf5_err)
   142)   call h5sclose_f(file_space_id,hdf5_err)
   143)   call h5dclose_f(data_set_id,hdf5_err)
   144)   
   145)   deallocate(dims_h5)
   146)   deallocate(max_dims_h5) 
   147) 
   148)   call PetscLogEventEnd(logging%event_read_ndim_real_array_hdf5, &
   149)                         ierr);CHKERRQ(ierr)
   150)                           
   151) end subroutine HDF5ReadNDimRealArray
   152) 
   153) #if defined(SCORPIO)
   154) 
   155) ! ************************************************************************** !
   156) 
   157) subroutine HDF5ReadDatasetInteger2D(filename,dataset_name,read_option,option, &
   158)            data,data_dims,dataset_dims)
   159)   ! 
   160)   ! Author: Gautam Bisht
   161)   ! Date: 05/13/2010
   162)   ! 
   163) 
   164)   use hdf5
   165)   use Option_module
   166)   
   167)   implicit none
   168)   
   169) #if defined(SCORPIO)
   170)   include "scorpiof.h"  
   171) #endif
   172) 
   173)   ! in
   174)   character(len=MAXSTRINGLENGTH) :: filename
   175)   character(len=MAXSTRINGLENGTH) :: dataset_name
   176)   integer :: read_option
   177)   type(option_type) :: option
   178)   
   179)   ! out
   180)   PetscInt,pointer :: data(:,:)
   181)   PetscInt :: data_dims(2)
   182)   PetscInt :: dataset_dims(2)
   183)   
   184)   ! local
   185)   PetscInt :: file_id
   186)   PetscInt :: ndims
   187)   PetscInt :: ii, remainder
   188) 
   189)   PetscErrorCode :: ierr
   190)   
   191)   ! Open file collectively
   192)   filename = trim(filename) // CHAR(0)
   193)   call fscorpio_open_file(filename, option%ioread_group_id, SCORPIO_FILE_READONLY, file_id, ierr)
   194) 
   195)   ! Get dataset dimnesions
   196)   call fscorpio_get_dataset_ndims(ndims, file_id, dataset_name, option%ioread_group_id, ierr)
   197)   if (ndims > 2) then
   198)     option%io_buffer='Dimension of ' // dataset_name // ' dataset in ' // filename // &
   199)     ' is greater than to 2.'
   200)     call printErrMsg(option)
   201)   endif
   202)   
   203)   ! Get size of each dimension
   204)   call fscorpio_get_dataset_dims(dataset_dims, file_id, dataset_name, option%ioread_group_id, ierr)
   205)   
   206)   data_dims(1) = dataset_dims(1)/option%mycommsize
   207)   data_dims(2) = dataset_dims(2)
   208) 
   209)   remainder = dataset_dims(1) - data_dims(1)*option%mycommsize
   210)   if (option%myrank < remainder) data_dims(1) = data_dims(1) + 1
   211) 
   212)   
   213)   allocate(data(data_dims(2),dataset_dims(1)))
   214)   
   215)   call fscorpio_get_dataset_dims(dataset_dims, file_id, dataset_name, option%ioread_group_id, ierr)
   216) 
   217)   ! Read the dataset collectively
   218)   call fscorpio_read_dataset( data, SCORPIO_INTEGER, ndims, dataset_dims, data_dims, & 
   219)             file_id, dataset_name, option%ioread_group_id, SCORPIO_NONUNIFORM_CONTIGUOUS_READ, ierr)
   220)   
   221)   data_dims(1) = data_dims(1) + data_dims(2)
   222)   data_dims(2) = data_dims(1) - data_dims(2)
   223)   data_dims(1) = data_dims(1) - data_dims(2)
   224)   
   225)   dataset_dims(1) = dataset_dims(1) + dataset_dims(2)
   226)   dataset_dims(2) = dataset_dims(1) - dataset_dims(2)
   227)   dataset_dims(1) = dataset_dims(1) - dataset_dims(2)
   228) 
   229)   ! Close file
   230)   call fscorpio_close_file( file_id, option%ioread_group_id, ierr)  
   231) 
   232) end subroutine HDF5ReadDatasetInteger2D
   233) #endif
   234) ! SCORPIO
   235) 
   236) #if defined(SCORPIO)
   237) 
   238) ! ************************************************************************** !
   239) 
   240) subroutine HDF5ReadDatasetReal2D(filename,dataset_name,read_option,option, &
   241)            data,data_dims,dataset_dims)
   242)   use hdf5
   243)   use Option_module
   244)   
   245)   implicit none
   246)   
   247) #if defined(SCORPIO)
   248)   include "scorpiof.h"  
   249) #endif
   250) 
   251)   ! in
   252)   character(len=MAXSTRINGLENGTH) :: filename
   253)   character(len=MAXSTRINGLENGTH) :: dataset_name
   254)   integer :: read_option
   255)   type(option_type) :: option
   256)   
   257)   ! out
   258)   PetscReal,pointer :: data(:,:)
   259)   PetscInt :: data_dims(2)
   260)   PetscInt :: dataset_dims(2)
   261)   
   262)   ! local
   263)   integer :: file_id
   264)   integer :: ndims
   265)   PetscInt :: ii, remainder
   266) 
   267)   PetscErrorCode :: ierr
   268)   
   269)   ! Open file collectively
   270)   filename = trim(filename) // CHAR(0)
   271)   call fscorpio_open_file(filename, option%ioread_group_id, SCORPIO_FILE_READONLY, file_id, ierr)
   272) 
   273)   ! Get dataset dimnesions
   274)   call fscorpio_get_dataset_ndims(ndims, file_id, dataset_name, option%ioread_group_id, ierr)
   275)   if (ndims > 2) then
   276)     option%io_buffer='Dimension of ' // dataset_name // ' dataset in ' // filename // &
   277)     ' is greater than to 2.'
   278)     call printErrMsg(option)
   279)   endif
   280)   
   281)   ! Get size of each dimension
   282)   call fscorpio_get_dataset_dims(dataset_dims, file_id, dataset_name, option%ioread_group_id, ierr)
   283)   
   284)   data_dims(1) = dataset_dims(1)/option%mycommsize
   285)   data_dims(2) = dataset_dims(2)
   286) 
   287)   remainder = dataset_dims(1) - data_dims(1)*option%mycommsize
   288)   if (option%myrank < remainder) data_dims(1) = data_dims(1) + 1
   289) 
   290)   
   291)   allocate(data(data_dims(2),dataset_dims(1)))
   292)   
   293)   call fscorpio_get_dataset_dims(dataset_dims, file_id, dataset_name, option%ioread_group_id, ierr)
   294) 
   295)   ! Read the dataset collectively
   296)   call fscorpio_read_dataset( data, SCORPIO_DOUBLE, ndims, dataset_dims, data_dims, & 
   297)             file_id, dataset_name, option%ioread_group_id, SCORPIO_NONUNIFORM_CONTIGUOUS_READ, ierr)
   298)   
   299)   data_dims(1) = data_dims(1) + data_dims(2)
   300)   data_dims(2) = data_dims(1) - data_dims(2)
   301)   data_dims(1) = data_dims(1) - data_dims(2)
   302) 
   303)   dataset_dims(1) = dataset_dims(1) + dataset_dims(2)
   304)   dataset_dims(2) = dataset_dims(1) - dataset_dims(2)
   305)   dataset_dims(1) = dataset_dims(1) - dataset_dims(2)
   306) 
   307)   ! Close file
   308)   call fscorpio_close_file( file_id, option%ioread_group_id, ierr)  
   309) 
   310) end subroutine HDF5ReadDatasetReal2D
   311) #endif
   312) 
   313) ! ************************************************************************** !
   314) 
   315) function HDF5GroupExists(filename,group_name,option)
   316)   ! 
   317)   ! SCORPIO
   318)   ! Returns true if a group exists
   319)   ! 
   320)   ! Author: Glenn Hammond
   321)   ! Date: 03/26/2012
   322)   ! 
   323) 
   324)   use hdf5
   325)   use Option_module
   326)   
   327)   implicit none
   328) 
   329) #if defined(SCORPIO)
   330)   include "scorpiof.h"  
   331) #endif
   332)   
   333)   character(len=MAXSTRINGLENGTH) :: filename
   334)   character(len=MAXWORDLENGTH) :: group_name
   335)   type(option_type) :: option
   336) 
   337)   integer(HID_T) :: file_id  
   338)   integer(HID_T) :: grp_id  
   339)   integer(HID_T) :: prop_id
   340)   PetscMPIInt, parameter :: ON=1, OFF=0 
   341)   PetscBool :: group_exists
   342)   
   343)   PetscBool :: HDF5GroupExists
   344) 
   345) #if defined(SCORPIO)
   346) 
   347)   ! Open file collectively
   348)   filename = trim(filename) // CHAR(0)
   349)   group_name = trim(group_name) // CHAR(0)
   350)   call fscorpio_open_file(filename, option%ioread_group_id, SCORPIO_FILE_READONLY, file_id, ierr)
   351)   call fscorpio_group_exists(group_name, file_id, option%ioread_group_id, ierr)
   352)   group_exists = (ierr == 1);
   353) 
   354)   if (group_exists) then
   355)     HDF5GroupExists = PETSC_TRUE
   356)     option%io_buffer = 'Group "' // trim(group_name) // '" in HDF5 file "' // &
   357)       trim(filename) // '" found in file.'
   358)   else
   359)     HDF5GroupExists = PETSC_FALSE
   360)     option%io_buffer = 'Group "' // trim(group_name) // '" in HDF5 file "' // &
   361)       trim(filename) // '" not found in file.  Therefore, assuming a ' // &
   362)       'cell-indexed dataset.'
   363)   endif
   364)   call printMsg(option)
   365) 
   366)   call fscorpio_close_file( file_id, option%ioread_group_id, ierr)  
   367) #else
   368)   ! open the file
   369)   call h5open_f(hdf5_err)
   370)   ! set read file access property
   371)   call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
   372) #ifndef SERIAL_HDF5
   373)   call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
   374) #endif
   375)   call HDF5OpenFileReadOnly(filename,file_id,prop_id,option)
   376)   call h5pclose_f(prop_id,hdf5_err)
   377) 
   378)   option%io_buffer = 'Testing group: ' // trim(group_name)
   379)   call printMsg(option)
   380)   ! I turn off error messaging since if the group does not exist, an error
   381)   ! will be printed, but the user does not need to see this.
   382)   call h5eset_auto_f(OFF,hdf5_err)
   383)   call h5gopen_f(file_id,group_name,grp_id,hdf5_err)
   384)   group_exists = .not.(hdf5_err < 0)
   385)   call h5eset_auto_f(ON,hdf5_err)  
   386) 
   387)   if (group_exists) then
   388)     HDF5GroupExists = PETSC_TRUE
   389)     call h5gclose_f(grp_id,hdf5_err)  
   390)     option%io_buffer = 'Group "' // trim(group_name) // '" in HDF5 file "' // &
   391)       trim(filename) // '" found in file.'
   392)   else
   393)     HDF5GroupExists = PETSC_FALSE
   394)     option%io_buffer = 'Group "' // trim(group_name) // '" in HDF5 file "' // &
   395)       trim(filename) // '" not found in file.  Therefore, assuming a ' // &
   396)       'cell-indexed dataset.'
   397)   endif
   398)   call printMsg(option)
   399) 
   400)   call h5fclose_f(file_id,hdf5_err)
   401)   call h5close_f(hdf5_err)  
   402) #endif 
   403) !SCORPIO
   404) 
   405) end function HDF5GroupExists
   406) 
   407) ! ************************************************************************** !
   408) 
   409) function HDF5DatasetExists(filename,group_name,dataset_name,option)
   410)   !
   411)   ! SCORPIO
   412)   ! Returns true if a dataset exists
   413)   !
   414)   ! Author: Gautam Bisht
   415)   ! Date: 04/30/2015
   416)   !
   417) 
   418)   use hdf5
   419)   use Option_module
   420) 
   421)   implicit none
   422) 
   423) #if defined(SCORPIO)
   424)   include "scorpiof.h"
   425) #endif
   426) 
   427)   character(len=MAXSTRINGLENGTH) :: filename
   428)   character(len=MAXWORDLENGTH) :: group_name
   429)   character(len=MAXWORDLENGTH) :: dataset_name
   430)   type(option_type) :: option
   431) 
   432)   character(len=MAXWORDLENGTH) :: group_name_local
   433)   integer(HID_T) :: file_id
   434)   integer(HID_T) :: grp_id
   435)   integer(HID_T) :: prop_id
   436)   integer(HID_T) :: dataset_id
   437)   PetscMPIInt, parameter :: ON=1, OFF=0
   438)   PetscBool :: group_exists
   439)   PetscBool :: dataset_exists
   440) 
   441)   PetscBool :: HDF5DatasetExists
   442) 
   443) #if defined(SCORPIO)
   444) 
   445)   option%io_buffer = 'Need to extend HDF5DatasetExists() for SCORPIO.'
   446)   call printErrMsg(option)
   447) 
   448) #else
   449) 
   450)   if (len_trim(group_name) == 0) then
   451)     group_name_local = "/" // CHAR(0)
   452)   else
   453)     group_name_local = trim(group_name) // "/" // CHAR(0)
   454)   endif
   455) 
   456)   ! open the file
   457)   call h5open_f(hdf5_err)
   458)   ! set read file access property
   459)   call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
   460) #ifndef SERIAL_HDF5
   461)   call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
   462) #endif
   463)   call HDF5OpenFileReadOnly(filename,file_id,prop_id,option)
   464)   call h5pclose_f(prop_id,hdf5_err)
   465) 
   466)   ! I turn off error messaging since if the group does not exist, an error
   467)   ! will be printed, but the user does not need to see this.
   468)   call h5eset_auto_f(OFF,hdf5_err)
   469)   call h5gopen_f(file_id,group_name_local,grp_id,hdf5_err)
   470)   group_exists = .not.(hdf5_err < 0)
   471)   call h5eset_auto_f(ON,hdf5_err)
   472) 
   473)   if (.not.group_exists) then
   474)     HDF5DatasetExists = PETSC_FALSE
   475)   endif
   476) 
   477)   call h5dopen_f(grp_id,dataset_name,dataset_id,hdf5_err)
   478)   dataset_exists = .not.(hdf5_err < 0)
   479) 
   480)   if (.not.dataset_exists) then
   481)     HDF5DatasetExists = PETSC_FALSE
   482)   else
   483)     HDF5DatasetExists = PETSC_TRUE
   484)     call h5dclose_f(dataset_id,hdf5_err)
   485)   endif
   486) 
   487)   if (group_exists) call h5gclose_f(grp_id,hdf5_err)
   488) 
   489)   call h5fclose_f(file_id,hdf5_err)
   490)   call h5close_f(hdf5_err)
   491) #endif
   492) !SCORPIO
   493) 
   494) end function HDF5DatasetExists
   495) 
   496) ! ************************************************************************** !
   497) 
   498) subroutine HDF5MakeStringCompatible(name)
   499)   ! 
   500)   ! Replaces '/' in string with '_' for hdf5 names
   501)   ! 
   502)   ! Author: Glenn Hammond
   503)   ! Date: 10/25/07
   504)   ! 
   505) 
   506)   implicit none
   507)   
   508)   character(len=*) :: name
   509)   
   510)   PetscInt :: len, ichar
   511)   
   512)   len = len_trim(name)
   513)   do ichar = 1, len
   514)     if (name(ichar:ichar) == '/') then
   515)       name(ichar:ichar) = '_'
   516)     endif
   517)   enddo
   518)   
   519)   name = trim(name)
   520) 
   521) end subroutine HDF5MakeStringCompatible
   522) 
   523) ! ************************************************************************** !
   524) 
   525) subroutine HDF5ReadDbase(filename,option)
   526)   ! 
   527)   ! Read in an ASCII database
   528)   ! 
   529)   ! Author: Glenn Hammond
   530)   ! Date: 08/19/14
   531)   ! 
   532)   use Option_module
   533)   use String_module
   534)   use Input_Aux_module, only : dbase
   535)   use h5lt
   536)   
   537)   implicit none
   538)   
   539)   character(len=MAXWORDLENGTH) :: filename
   540)   type(option_type) :: option
   541)   
   542)   character(len=MAXWORDLENGTH), allocatable :: wbuffer(:)
   543)   character(len=MAXWORDLENGTH) :: wbuffer_word
   544)   PetscReal, allocatable :: rbuffer(:)
   545)   PetscInt, allocatable :: ibuffer(:)
   546)   PetscInt :: dummy_int
   547)   PetscInt :: value_index
   548)   character(len=MAXWORDLENGTH) :: object_name, word
   549) #if defined(PETSC_HAVE_HDF5)  
   550)   integer(HID_T) :: file_id
   551)   integer(HID_T) :: num_objects
   552)   integer(HID_T) :: i_object
   553)   integer(HID_T) :: object_type
   554)   integer(HID_T) :: prop_id
   555)   integer(HID_T) :: dataset_id
   556)   integer(HID_T) :: class_id
   557)   integer(HID_T) :: datatype_id
   558)   integer(HID_T) :: datatype_id2
   559)   integer(HID_T) :: file_space_id
   560)   integer(HID_T) :: memory_space_id
   561)   integer(HSIZE_T) :: num_values_in_dataset
   562)   integer(SIZE_T) size_t_int
   563) !  integer(HSIZE_T) :: dims(1)
   564)   integer(SIZE_T) :: type_size
   565)   integer(HSIZE_T) :: offset(1), length(1), stride(1)
   566)   PetscMPIInt :: rank_mpi
   567)   PetscMPIInt :: int_mpi
   568)   PetscMPIInt :: hdf5_err
   569) #endif
   570)   PetscInt :: num_ints
   571)   PetscInt :: num_reals
   572)   PetscInt :: num_words
   573) 
   574)   call h5open_f(hdf5_err)
   575)   option%io_buffer = 'Opening hdf5 file: ' // trim(filename)
   576)   call printMsg(option)
   577)   call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
   578) #ifndef SERIAL_HDF5
   579)   call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
   580) #endif
   581)   call HDF5OpenFileReadOnly(filename,file_id,prop_id,option)
   582)   call h5pclose_f(prop_id,hdf5_err)
   583)   call h5gn_members_f(file_id, '.',num_objects, hdf5_err)
   584)   num_ints = 0
   585)   num_reals = 0
   586)   num_words = 0
   587)   ! index is zero-based
   588)   do i_object = 0, num_objects-1
   589)     call h5gget_obj_info_idx_f(file_id,'.',i_object,object_name, &
   590)                                object_type,hdf5_err)
   591)     if (object_type == H5G_DATASET_F) then
   592)       call h5dopen_f(file_id,object_name,dataset_id,hdf5_err)
   593)       call h5dget_type_f(dataset_id, datatype_id, hdf5_err)
   594)       call h5tget_class_f(datatype_id, class_id, hdf5_err)
   595)       ! cannot use a select case statement since the H5T definitions are not
   596)       ! guaranteed to be constant.  the preprocessor throws an error
   597)       if (class_id == H5T_INTEGER_F) then
   598)         num_ints = num_ints + 1
   599)       else if (class_id == H5T_FLOAT_F) then
   600)         num_reals = num_reals + 1
   601)       else if (class_id == H5T_STRING_F) then
   602)         num_words = num_words + 1
   603)       else
   604)         option%io_buffer = 'Unrecognized HDF5 datatype in Dbase: ' // &
   605)           trim(object_name)
   606)         call printErrMsg(option)
   607)       endif
   608)       call h5tclose_f(datatype_id, hdf5_err)
   609)       call h5dclose_f(dataset_id,hdf5_err)
   610)     endif
   611)   enddo
   612)   allocate(dbase)
   613)   if (num_ints > 0) then
   614)     allocate(dbase%icard(num_ints))
   615)     dbase%icard = ''
   616)     allocate(dbase%ivalue(num_ints))
   617)     dbase%ivalue = UNINITIALIZED_INTEGER
   618)   endif
   619)   if (num_reals > 0) then
   620)     allocate(dbase%rcard(num_reals))
   621)     dbase%rcard = ''
   622)     allocate(dbase%rvalue(num_reals))
   623)     dbase%rvalue = UNINITIALIZED_DOUBLE
   624)   endif
   625)   if (num_words > 0) then
   626)     allocate(dbase%ccard(num_words))
   627)     dbase%ccard = ''
   628)     allocate(dbase%cvalue(num_words))
   629)     dbase%cvalue = '-999'
   630)   endif
   631)   value_index = 1
   632)   if (option%id > 0) then
   633)     value_index = option%id
   634)   endif
   635)   num_ints = 0
   636)   num_reals = 0
   637)   num_words = 0
   638)   do i_object = 0, num_objects-1
   639)     call h5gget_obj_info_idx_f(file_id,'.',i_object,object_name, &
   640)                                object_type,hdf5_err)
   641)     if (object_type == H5G_DATASET_F) then
   642) ! use once HDF5 lite is linked in PETSc      
   643) !      call h5ltget_dataset_info_f(file_id,object_name,dims,dummy_int, &
   644) !                                  type_size,hdf5_err)
   645) !      allocate(buffer(dims(1)))
   646) !      buffer = 0.d0
   647) !      call h5ltread_dataset_double_f(file_id,object_name,buffer, &
   648) !                                     dims,hdf5_err)
   649) !      dbase%card(icount) = trim(object_name)
   650) !      if (option%id > 0) then
   651) !        if (option%id > dims(1)) then
   652) !          write(word,*) dims(1)
   653) !          option%io_buffer = 'DBASE dataset "' // trim(object_name) // &
   654) !            '" is too small (' // trim(adjustl(word)) // &
   655) !            ') for number of realizations.'
   656) !          call printErrMsg(option)
   657) !        endif
   658) !        dbase%value(icount) = buffer(option%id)
   659) !      else
   660) !        dbase%value(icount) = buffer(1)
   661) !      endif
   662) !      deallocate(buffer)
   663) 
   664)       call h5dopen_f(file_id,object_name,dataset_id,hdf5_err)
   665)       call h5dget_space_f(dataset_id,file_space_id,hdf5_err)
   666)       ! should be a rank=1 data space
   667)       call h5sget_simple_extent_npoints_f(file_space_id, &
   668)                                           num_values_in_dataset,hdf5_err)
   669)       if (option%id > 0) then
   670)         if (option%id > num_values_in_dataset) then
   671)           write(word,*) num_values_in_dataset
   672)           option%io_buffer = 'Data in DBASE_FILENAME "' // &
   673)             trim(object_name) // &
   674)             '" is too small (' // trim(adjustl(word)) // &
   675)             ') for number of realizations.'
   676)           call printErrMsg(option)
   677)         endif
   678)       endif
   679)       rank_mpi = 1
   680)       offset = 0
   681)       length = num_values_in_dataset
   682)       stride = 1
   683)       call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
   684) #ifndef SERIAL_HDF5
   685)       call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F,hdf5_err)
   686) #endif
   687)       call h5screate_simple_f(rank_mpi,length,memory_space_id,hdf5_err, &
   688)                               length)
   689) 
   690)       call h5dget_type_f(dataset_id, datatype_id, hdf5_err)
   691)       call h5tget_class_f(datatype_id, class_id, hdf5_err)
   692)       call h5tclose_f(datatype_id, hdf5_err)
   693) #ifdef HDF5_BROADCAST
   694)       if (option%myrank == option%io_rank) then                           
   695) #endif
   696)       call PetscLogEventBegin(logging%event_h5dread_f,ierr);CHKERRQ(ierr)
   697)       if (class_id == H5T_INTEGER_F) then
   698)         allocate(ibuffer(num_values_in_dataset))
   699)         ibuffer = 0
   700)         call h5dread_f(dataset_id,H5T_NATIVE_INTEGER,ibuffer,length, &
   701)                        hdf5_err,memory_space_id,file_space_id,prop_id)
   702)       else if (class_id == H5T_FLOAT_F) then
   703)         allocate(rbuffer(num_values_in_dataset))
   704)         rbuffer = 0.d0
   705)         call h5dread_f(dataset_id,H5T_NATIVE_DOUBLE,rbuffer,length, &
   706)                        hdf5_err,memory_space_id,file_space_id,prop_id)
   707)       else if (class_id == H5T_STRING_F) then
   708)         call h5tcopy_f(H5T_NATIVE_CHARACTER,datatype_id2,hdf5_err)
   709)         size_t_int = MAXWORDLENGTH
   710)         call h5tset_size_f(datatype_id2,size_t_int,hdf5_err)
   711)         allocate(wbuffer(num_values_in_dataset))
   712)         wbuffer = ''
   713)         call h5dread_f(dataset_id,datatype_id2,wbuffer,length, &
   714)                        hdf5_err,memory_space_id,file_space_id,prop_id)
   715)         wbuffer_word = wbuffer(value_index)
   716)         deallocate(wbuffer)
   717)         call h5tclose_f(datatype_id2,hdf5_err)
   718)       endif
   719)       call PetscLogEventEnd(logging%event_h5dread_f,ierr);CHKERRQ(ierr)
   720) #ifdef HDF5_BROADCAST
   721)       endif
   722)       if (option%mycommsize > 1) then
   723)         int_mpi = num_values_in_dataset
   724)         if (class_id == H5T_INTEGER_F) then
   725)           call MPI_Bcast(ibuffer,int_mpi,MPI_INTEGER, &
   726)                          option%io_rank,option%mycomm,ierr)
   727)         else if (class_id == H5T_FLOAT_F) then
   728)           call MPI_Bcast(rbuffer,int_mpi,MPI_DOUBLE_PRECISION, &
   729)                          option%io_rank,option%mycomm,ierr)
   730)         else if (class_id == H5T_STRING_F) then
   731)           int_mpi = MAXWORDLENGTH
   732)           call MPI_Bcast(wbuffer_word,int_mpi,MPI_CHARACTER, &
   733)                          option%io_rank,option%mycomm,ierr)
   734)         endif
   735)       endif
   736) #endif
   737)       call h5pclose_f(prop_id,hdf5_err)
   738)       if (memory_space_id > -1) call h5sclose_f(memory_space_id,hdf5_err)
   739)       call h5sclose_f(file_space_id,hdf5_err)
   740)       call h5dclose_f(dataset_id,hdf5_err)
   741)       call StringToUpper(object_name)
   742)       ! these conditionals must come after the bcasts above!!!
   743)       if (class_id == H5T_INTEGER_F) then
   744)         num_ints = num_ints + 1
   745)         dbase%icard(num_ints) = trim(object_name)
   746)         dbase%ivalue(num_ints) = ibuffer(value_index)
   747)         deallocate(ibuffer)
   748)       else if (class_id == H5T_FLOAT_F) then
   749)         num_reals = num_reals + 1
   750)         dbase%rcard(num_reals) = trim(object_name)
   751)         dbase%rvalue(num_reals) = rbuffer(value_index)
   752)         deallocate(rbuffer)
   753)       else if (class_id == H5T_STRING_F) then
   754)         num_words = num_words + 1
   755)         dbase%ccard(num_words) = trim(object_name)
   756)         dbase%cvalue(num_words) = wbuffer_word
   757)       endif
   758)     endif
   759)   enddo
   760)   call h5fclose_f(file_id,hdf5_err)
   761)   call h5close_f(hdf5_err)
   762)       
   763) end subroutine HDF5ReadDbase
   764) 
   765) ! ************************************************************************** !
   766) 
   767) subroutine HDF5OpenFileReadOnly(filename,file_id,prop_id,option)
   768)   ! 
   769)   ! Opens an HDF5 file.  This wrapper provides error messaging if the file
   770)   ! does not exist.
   771)   ! 
   772)   ! Author: Glenn Hammond
   773)   ! Date: 06/22/15
   774)   ! 
   775)   use hdf5
   776)   use Option_module
   777)   
   778)   character(len=*) :: filename  ! must be of variable length
   779)   integer(HID_T) :: file_id
   780)   integer(HID_T) :: prop_id
   781)   type(option_type) :: option
   782)   
   783)   PetscMPIInt :: hdf5_err
   784) 
   785)   call h5fopen_f(filename,H5F_ACC_RDONLY_F,file_id,hdf5_err,prop_id)
   786)   if (hdf5_err /= 0) then
   787)     option%io_buffer = 'HDF5 file "' // trim(filename) // '" not found.'
   788)     call printErrMsg(option)
   789)   endif
   790)   
   791) end subroutine HDF5OpenFileReadOnly
   792) 
   793) #endif
   794) ! defined(PETSC_HAVE_HDF5)
   795) 
   796) end module HDF5_Aux_module

generated by
Intel(R) C++/Fortran Compiler code-coverage tool
Web-Page Owner: Nobody