dataset_map_hdf5.F90       coverage:  0.00 %func     0.00 %block


     1) module Dataset_Map_HDF5_class
     2)  
     3)   use Dataset_Common_HDF5_class
     4)   
     5)   use PFLOTRAN_Constants_module
     6) 
     7)   implicit none
     8) 
     9)   private
    10) 
    11) #include "petsc/finclude/petscsys.h"
    12) 
    13)   type, public, extends(dataset_common_hdf5_type) :: dataset_map_hdf5_type
    14)     character(len=MAXSTRINGLENGTH) :: h5_dataset_map_name
    15)     character(len=MAXSTRINGLENGTH) :: map_filename
    16)     PetscInt, pointer :: mapping(:,:)
    17)     PetscInt :: map_dims_global(2)
    18)     PetscInt :: map_dims_local(2)
    19)     PetscInt, pointer :: datatocell_ids(:)
    20)     PetscInt, pointer :: cell_ids_local(:)
    21)     PetscBool :: first_time
    22)   end type dataset_map_hdf5_type
    23)   
    24)   PetscInt, parameter :: MAX_NSLICE = 100
    25)   
    26)   public :: DatasetMapHDF5Create, &
    27)             DatasetMapHDF5Init, &
    28)             DatasetMapHDF5Cast, &
    29)             DatasetMapHDF5Read, &
    30)             DatasetMapHDF5Load, &
    31)             DatasetMapHDF5Print, &
    32)             DatasetMapHDF5Destroy
    33)   
    34) contains
    35) 
    36) ! ************************************************************************** !
    37) 
    38) function DatasetMapHDF5Create()
    39)   ! 
    40)   ! Creates global dataset class
    41)   ! 
    42)   ! Author: Glenn Hammond
    43)   ! Date: 05/29/13
    44)   ! 
    45)   
    46)   implicit none
    47)   
    48)   class(dataset_map_hdf5_type), pointer :: dataset
    49) 
    50)   class(dataset_map_hdf5_type), pointer :: DatasetMapHDF5Create
    51)   
    52)   allocate(dataset)
    53)   call DatasetMapHDF5Init(dataset)
    54) 
    55)   DatasetMapHDF5Create => dataset
    56)     
    57) end function DatasetMapHDF5Create
    58) 
    59) ! ************************************************************************** !
    60) 
    61) subroutine DatasetMapHDF5Init(this)
    62)   ! 
    63)   ! Initializes members of global dataset class
    64)   ! 
    65)   ! Author: Glenn Hammond
    66)   ! Date: 05/29/13
    67)   ! 
    68)   
    69)   implicit none
    70)   
    71)   class(dataset_map_hdf5_type) :: this
    72)   
    73)   call DatasetCommonHDF5Init(this)
    74)   this%h5_dataset_map_name = ''
    75)   this%map_filename = ''
    76)   nullify(this%mapping)
    77)   this%map_dims_global = 0
    78)   this%map_dims_local = 0
    79)   nullify(this%datatocell_ids)
    80)   nullify(this%cell_ids_local)
    81)   this%first_time = PETSC_TRUE
    82)     
    83) end subroutine DatasetMapHDF5Init
    84) 
    85) ! ************************************************************************** !
    86) 
    87) function DatasetMapHDF5Cast(this)
    88)   ! 
    89)   ! Casts a dataset_base_type to dataset_map_hdf5_type
    90)   ! 
    91)   ! Date: 05/03/13
    92)   ! 
    93) 
    94)   use Dataset_Base_class
    95)   
    96)   implicit none
    97) 
    98)   class(dataset_base_type), pointer :: this
    99) 
   100)   class(dataset_map_hdf5_type), pointer :: DatasetMapHDF5Cast
   101)   
   102)   nullify(DatasetMapHDF5Cast)
   103)   select type (this)
   104)     class is (dataset_map_hdf5_type)
   105)       DatasetMapHDF5Cast => this
   106)   end select
   107)     
   108) end function DatasetMapHDF5Cast
   109) 
   110) ! ************************************************************************** !
   111) 
   112) subroutine DatasetMapHDF5Read(this,input,option)
   113)   ! 
   114)   ! Reads in contents of a dataset card
   115)   ! 
   116)   ! Author: Glenn Hammond
   117)   ! Date: 01/12/11, 06/04/13
   118)   ! 
   119) 
   120)   use Option_module
   121)   use Input_Aux_module
   122)   use String_module
   123) 
   124)   implicit none
   125)   
   126)   class(dataset_map_hdf5_type) :: this
   127)   type(input_type), pointer :: input
   128)   type(option_type) :: option
   129)   
   130)   character(len=MAXWORDLENGTH) :: keyword
   131)   PetscBool :: found
   132) 
   133)   input%ierr = 0
   134)   do
   135)   
   136)     call InputReadPflotranString(input,option)
   137) 
   138)     if (InputCheckExit(input,option)) exit  
   139) 
   140)     call InputReadWord(input,option,keyword,PETSC_TRUE)
   141)     call InputErrorMsg(input,option,'keyword','DATASET')
   142)     call StringToUpper(keyword)   
   143)       
   144)     call DatasetCommonHDF5ReadSelectCase(this,input,keyword,found,option)
   145)     if (.not.found) then
   146)       select case(trim(keyword))
   147)         case('MAP_HDF5_DATASET_NAME') 
   148)           call InputReadWord(input,option,this%h5_dataset_map_name,PETSC_TRUE)
   149)           call InputErrorMsg(input,option,'map dataset name','DATASET')
   150)         case('MAP_FILENAME_NAME') 
   151)           call InputReadWord(input,option,this%map_filename,PETSC_TRUE)
   152)           call InputErrorMsg(input,option,'map filename','DATASET')
   153)         case default
   154)           call InputKeywordUnrecognized(keyword,'dataset',option)
   155)       end select
   156)     endif
   157)   
   158)   enddo
   159)   
   160)   if (len_trim(this%hdf5_dataset_name) < 1) then
   161)     this%hdf5_dataset_name = this%name
   162)   endif
   163) 
   164)   if (len_trim(this%map_filename) < 1) then
   165)     this%map_filename = this%filename
   166)   endif
   167)   
   168) end subroutine DatasetMapHDF5Read
   169) 
   170) ! ************************************************************************** !
   171) 
   172) subroutine DatasetMapHDF5Load(this,option)
   173)   ! 
   174)   ! Load new data into dataset buffer
   175)   ! 
   176)   ! Author: Glenn Hammond
   177)   ! Date: 05/29/13
   178)   ! 
   179)   
   180)   use Option_module
   181)   use Time_Storage_module
   182)   use Dataset_Base_class
   183) 
   184)   implicit none
   185)   
   186)   class(dataset_map_hdf5_type) :: this
   187)   type(option_type) :: option
   188)   
   189)   if (DatasetCommonHDF5Load(this,option)) then
   190) #if defined(PETSC_HAVE_HDF5)    
   191)     if (.not.associated(this%mapping)) then
   192)       call DatasetMapHDF5ReadMap(this,option)
   193)     endif
   194)     call DatasetMapHDF5ReadData(this,option)
   195) #endif    
   196) !    call this%Reorder(option)
   197)     call DatasetBaseReorder(this,option)
   198)   endif
   199)   call DatasetBaseInterpolateTime(this)
   200)     
   201) end subroutine DatasetMapHDF5Load
   202) 
   203) #if defined(PETSC_HAVE_HDF5)
   204) 
   205) ! ************************************************************************** !
   206) 
   207) subroutine DatasetMapHDF5ReadData(this,option)
   208)   ! 
   209)   ! Read an hdf5 data into a array
   210)   ! 
   211)   ! Author: Glenn Hammond
   212)   ! Date: 10/25/11, 05/29/13
   213)   ! 
   214) 
   215)   use hdf5
   216)   use Option_module
   217)   use Units_module
   218)   use Logging_module
   219)   use HDF5_Aux_module
   220)   
   221)   implicit none
   222)   
   223)   class(dataset_map_hdf5_type) :: this
   224)   type(option_type) :: option
   225)   
   226)   integer(HID_T) :: file_id
   227)   integer(HID_T) :: file_space_id
   228)   integer(HID_T) :: memory_space_id
   229)   integer(HID_T) :: dataset_id
   230)   integer(HID_T) :: prop_id
   231)   integer(HID_T) :: grp_id
   232)   integer(HID_T) :: attribute_id
   233)   integer(HID_T) :: ndims_h5
   234)   integer(HID_T) :: atype_id
   235)   integer(HSIZE_T), allocatable :: dims_h5(:), max_dims_h5(:)
   236)   integer(HSIZE_T) :: attribute_dim(3)
   237)   integer(HSIZE_T) :: offset(4), length(4), stride(4)
   238)   integer(HSIZE_T) :: num_data_values
   239)   PetscInt :: i, temp_array(4)
   240)   PetscInt :: time_dim, num_times
   241)   PetscInt :: num_dims_in_h5_file, num_times_in_h5_file
   242)   PetscMPIInt :: array_rank_mpi
   243)   PetscMPIInt :: hdf5_err
   244)   PetscErrorCode :: ierr
   245)   character(len=MAXWORDLENGTH) :: dataset_name, word
   246) 
   247)   call PetscLogEventBegin(logging%event_dataset_map_hdf5_read, &
   248)                           ierr);CHKERRQ(ierr)
   249) 
   250)   ! open the file
   251)   call h5open_f(hdf5_err)
   252)   option%io_buffer = 'Opening hdf5 file: ' // trim(this%filename)
   253)   call printMsg(option)
   254)   
   255)   ! set read file access property
   256)   call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
   257) #ifndef SERIAL_HDF5
   258)   call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
   259) #endif
   260)   call HDF5OpenFileReadOnly(this%filename,file_id,prop_id,option)
   261)   call h5pclose_f(prop_id,hdf5_err)
   262) 
   263)   ! the dataset is actually stored in a group.  the group contains
   264)   ! a "data" dataset and optionally a "time" dataset.
   265)   option%io_buffer = 'Opening group: ' // trim(this%hdf5_dataset_name)
   266)   call printMsg(option)  
   267)   call h5gopen_f(file_id,this%hdf5_dataset_name,grp_id,hdf5_err)
   268) 
   269)   time_dim = -1
   270)   num_times = 1
   271)   if (associated(this%time_storage)) then
   272)     num_times = this%time_storage%max_time_index
   273)     time_dim = 2
   274)   endif
   275)   
   276)   ! open the "data" dataset
   277)   dataset_name = 'Data'
   278)   if (this%realization_dependent) then
   279)     write(word,'(i9)') option%id
   280)     dataset_name = trim(dataset_name) // trim(adjustl(word))
   281)   endif
   282)   call h5dopen_f(grp_id,dataset_name,dataset_id,hdf5_err)
   283)   call h5dget_space_f(dataset_id,file_space_id,hdf5_err)
   284) 
   285)   ! get dataset dimensions
   286)   if (.not.associated(this%dims)) then
   287)     call h5sget_simple_extent_ndims_f(file_space_id,ndims_h5,hdf5_err)
   288)     allocate(dims_h5(ndims_h5))
   289)     allocate(max_dims_h5(ndims_h5))
   290)     call h5sget_simple_extent_dims_f(file_space_id,dims_h5,max_dims_h5,hdf5_err)
   291)     if (associated(this%time_storage)) then
   292)       num_times_in_h5_file = dims_h5(1)
   293)     else
   294)       num_times_in_h5_file = 0
   295)     endif
   296)     if ((ndims_h5 == 2 .and. .not.associated(this%time_storage)) .or. &
   297)         (ndims_h5 == 1 .and. associated(this%time_storage))) then
   298)       option%io_buffer = 'Inconsistent dimensions in dataset file.'
   299)       call printErrMsg(option)
   300)     endif
   301)     ! rank in space will always be 1
   302)     this%rank = 1
   303)     allocate(this%dims(this%rank))
   304)     ! have to invert dimensions, but do not count the time dimension
   305)     this%dims(1) = int(dims_h5(ndims_h5))
   306)     deallocate(dims_h5)
   307)     deallocate(max_dims_h5) 
   308) 
   309)     ! check to ensure that dataset is properly sized
   310)     call h5sget_simple_extent_npoints_f(file_space_id,num_data_values,hdf5_err)
   311)     if (num_data_values/num_times /= this%dims(1)) then
   312)       option%io_buffer = &
   313)         'Number of values in dataset does not match dimensions.'
   314)       call printErrMsg(option)
   315)     endif
   316)     if (associated(this%time_storage) .and. &
   317)         num_times_in_h5_file /= num_times) then
   318)       option%io_buffer = &
   319)         'Number of times does not match last dimension of data array.'
   320)       call printErrMsg(option)
   321)     endif
   322)     if (.not.associated(this%rarray)) then
   323)       allocate(this%rarray(this%dims(1)))
   324)     endif
   325)     this%rarray = 0.d0
   326)     if (associated(this%time_storage) .and. .not.associated(this%rbuffer)) then
   327)       ! buffered array
   328)       allocate(this%rbuffer(size(this%rarray)*MAX_NSLICE))
   329)       this%rbuffer = 0.d0
   330)     endif
   331)   endif
   332) 
   333)   call PetscLogEventBegin(logging%event_h5dread_f,ierr);CHKERRQ(ierr)
   334) 
   335)   if (associated(this%time_storage)) then
   336)     num_dims_in_h5_file = this%rank + 1
   337)   else
   338)     num_dims_in_h5_file = this%rank
   339)   endif
   340)   
   341)   array_rank_mpi = 1
   342)   length = 1
   343)   ! length (or size) must be adjusted according to the size of the 
   344)   ! remaining data in the file
   345)   this%buffer_nslice = min(MAX_NSLICE,(num_times-this%buffer_slice_offset))
   346)   if (time_dim > 0) then
   347)     length(1) = size(this%rarray) * this%buffer_nslice
   348)   else
   349)     length(1) = size(this%rarray)
   350)   endif
   351)   call h5screate_simple_f(array_rank_mpi,length,memory_space_id,hdf5_err,length)    
   352) 
   353)   length = 1
   354)   stride = 1
   355)   offset = 0
   356)   length(1) = this%dims(1)
   357)   ! cannot read beyond end of the buffer
   358)   if (time_dim > 0) then
   359)     length(time_dim) = this%buffer_nslice
   360)     offset(time_dim) = this%buffer_slice_offset
   361)   endif
   362)   
   363)   !geh: for some reason, we have to invert here.  Perhaps because the
   364)   !     dataset was generated in C???
   365)   temp_array(1:num_dims_in_h5_file) = int(length(1:num_dims_in_h5_file))
   366)   do i = 1, num_dims_in_h5_file
   367)     length(i) = temp_array(num_dims_in_h5_file-i+1)
   368)   enddo
   369)   temp_array(1:num_dims_in_h5_file) = int(offset(1:num_dims_in_h5_file))
   370)   do i = 1, num_dims_in_h5_file
   371)     offset(i) = temp_array(num_dims_in_h5_file-i+1)
   372)   enddo
   373)   ! stride is fine
   374)   
   375)   call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
   376) #ifndef SERIAL_HDF5
   377)   call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_INDEPENDENT_F,hdf5_err)
   378) #endif
   379)   call h5sselect_hyperslab_f(file_space_id, H5S_SELECT_SET_F,offset,length, &
   380)                              hdf5_err,stride,stride) 
   381)   if (associated(this%rbuffer)) then
   382)     call h5dread_f(dataset_id,H5T_NATIVE_DOUBLE,this%rbuffer,length, &
   383)                    hdf5_err,memory_space_id,file_space_id,prop_id)
   384)   else
   385)     call h5dread_f(dataset_id,H5T_NATIVE_DOUBLE,this%rarray,length, &
   386)                    hdf5_err,memory_space_id,file_space_id,prop_id)
   387) !    this%rmax = maxval(this%rarray)
   388) !    this%rmin = minval(this%rarray)
   389)   endif
   390)   
   391)   call h5pclose_f(prop_id,hdf5_err)
   392)   if (memory_space_id > -1) call h5sclose_f(memory_space_id,hdf5_err)
   393)   call h5sclose_f(file_space_id,hdf5_err)
   394)   call h5dclose_f(dataset_id,hdf5_err)  
   395) 
   396)   call PetscLogEventEnd(logging%event_h5dread_f,ierr);CHKERRQ(ierr)
   397) 
   398)   option%io_buffer = 'Closing group: ' // trim(this%hdf5_dataset_name)
   399)   call printMsg(option)  
   400)   call h5gclose_f(grp_id,hdf5_err)  
   401)   option%io_buffer = 'Closing hdf5 file: ' // trim(this%filename)
   402)   call printMsg(option)  
   403)   call h5fclose_f(file_id,hdf5_err)
   404)   call h5close_f(hdf5_err)
   405)   
   406)   call PetscLogEventEnd(logging%event_dataset_map_hdf5_read, &
   407)                         ierr);CHKERRQ(ierr)
   408)                           
   409) end subroutine DatasetMapHDF5ReadData
   410) 
   411) ! ************************************************************************** !
   412) 
   413) subroutine DatasetMapHDF5ReadMap(this,option)
   414)   ! 
   415)   ! Read an hdf5 array
   416)   ! 
   417)   ! Author: Glenn Hammond
   418)   ! Date: 10/25/11, 05/29/13
   419)   ! 
   420) 
   421)   use hdf5
   422)   use Option_module
   423)   use Units_module
   424)   use Logging_module
   425)   use HDF5_Aux_module
   426)   
   427)   implicit none
   428)   
   429)   class(dataset_map_hdf5_type) :: this
   430)   type(option_type) :: option
   431)   
   432)   integer(HID_T) :: file_id
   433)   integer(HID_T) :: file_space_id
   434)   integer(HID_T) :: memory_space_id
   435)   integer(HID_T) :: dataset_id
   436)   integer(HID_T) :: prop_id
   437)   integer(HID_T) :: grp_id
   438)   integer(HID_T) :: ndims_hdf5
   439)   integer(HSIZE_T), allocatable :: dims_h5(:), max_dims_h5(:)
   440)   integer(HSIZE_T) :: offset(2), length(2)
   441)   PetscInt :: i
   442)   PetscMPIInt :: array_rank_mpi
   443)   PetscMPIInt :: hdf5_err
   444)   PetscInt :: nids_local, remainder, istart, iend
   445)   PetscErrorCode :: ierr
   446)   character(len=MAXWORDLENGTH) :: dataset_name
   447) 
   448)   call PetscLogEventBegin(logging%event_dataset_map_hdf5_read, &
   449)                           ierr);CHKERRQ(ierr)
   450) 
   451)   ! open the file
   452)   call h5open_f(hdf5_err)
   453)   option%io_buffer = 'Opening hdf5 file: ' // trim(this%map_filename)
   454)   call printMsg(option)
   455)   
   456)   ! set read file access property
   457)   call h5pcreate_f(H5P_FILE_ACCESS_F,prop_id,hdf5_err)
   458) #ifndef SERIAL_HDF5
   459)   call h5pset_fapl_mpio_f(prop_id,option%mycomm,MPI_INFO_NULL,hdf5_err)
   460) #endif
   461)   call HDF5OpenFileReadOnly(this%map_filename,file_id,prop_id,option)
   462)   call h5pclose_f(prop_id,hdf5_err)
   463) 
   464)   ! the dataset is actually stored in a group.  the group contains
   465)   ! a "data" dataset and optionally a "time" dataset.
   466)   option%io_buffer = 'Opening group: ' // trim(this%h5_dataset_map_name)
   467)   call printMsg(option)  
   468)   call h5gopen_f(file_id,this%h5_dataset_map_name,grp_id,hdf5_err)
   469)   
   470)   ! Open the "data" dataset
   471)   dataset_name = 'Data'
   472)   call h5dopen_f(grp_id,dataset_name,dataset_id,hdf5_err)
   473)   call h5dget_space_f(dataset_id,file_space_id,hdf5_err)
   474) 
   475)   ! Get number of dimensions and check
   476)   call h5sget_simple_extent_ndims_f(file_space_id,ndims_hdf5,hdf5_err)
   477)   if (ndims_hdf5 /= 2) then
   478)     option%io_buffer='Dimension of '// trim(this%h5_dataset_map_name) // &
   479)       '/Data dataset in ' // trim(this%filename) // ' is not equal to 2.'
   480)     call printErrMsg(option)
   481)   endif
   482) 
   483)   ! Get dimensions of dataset
   484)   allocate(dims_h5(ndims_hdf5))
   485)   allocate(max_dims_h5(ndims_hdf5))
   486)   call h5sget_simple_extent_dims_f(file_space_id,dims_h5,max_dims_h5,hdf5_err)
   487)   
   488)   nids_local=int(dims_h5(2)/option%mycommsize)
   489)   remainder =int(dims_h5(2))-nids_local*option%mycommsize
   490)   if (option%myrank<remainder) nids_local=nids_local+1
   491) 
   492)   ! Find istart and iend
   493)   istart = 0
   494)   iend   = 0
   495)   call MPI_Exscan(nids_local, istart, ONE_INTEGER_MPI, &
   496)                   MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
   497)   call MPI_Scan(nids_local, iend, ONE_INTEGER_MPI, &
   498)                 MPIU_INTEGER, MPI_SUM, option%mycomm, ierr)
   499)   
   500)   ! Determine the length and offset of data to be read by each processor
   501)   length(1) = dims_h5(1)
   502)   length(2) = iend-istart
   503)   offset(1) = 0
   504)   offset(2) = istart
   505) 
   506)   ! TOdo: Only read part of data
   507) !  option%io_buffer='Gautam: Modify code for reading HDF5 to generate mapping '//&
   508) !   'of dataset'
   509) !  call printMsg(option)
   510) !  nids_local=dims_h5(2)
   511) !  length(:) = dims_h5(:)
   512) !  offset(:) = 0
   513)   
   514)   ! Save dimension size
   515)   this%map_dims_global(:) = int(dims_h5(:))
   516)   this%map_dims_local(:) = int(length(:))
   517)   
   518)   ! Create data space for dataset
   519)   array_rank_mpi=2
   520)   call h5screate_simple_f(array_rank_mpi,length,memory_space_id,hdf5_err)
   521) 
   522)   ! Create property list
   523)   call h5pcreate_f(H5P_DATASET_XFER_F,prop_id,hdf5_err)
   524) #ifndef SERIAL_HDF5
   525)   call h5pset_dxpl_mpio_f(prop_id,H5FD_MPIO_COLLECTIVE_F, hdf5_err)
   526) #endif
   527) 
   528)   ! Select hyperslab
   529)   call h5sselect_hyperslab_f(file_space_id,H5S_SELECT_SET_F,offset,length, &
   530)                              hdf5_err)
   531)   
   532)   ! Initialize data buffer
   533)   allocate(this%mapping(length(1), length(2)))
   534) 
   535)   ! Read the dataset collectively
   536)   call h5dread_f(dataset_id, H5T_NATIVE_INTEGER, this%mapping, &
   537)                  dims_h5, hdf5_err, memory_space_id, file_space_id,prop_id)
   538) 
   539)   call h5pclose_f(prop_id,hdf5_err)
   540) 
   541)   if (memory_space_id > -1) call h5sclose_f(memory_space_id,hdf5_err)
   542)   call h5sclose_f(file_space_id,hdf5_err)
   543)   call h5dclose_f(dataset_id,hdf5_err)  
   544)   
   545)   option%io_buffer = 'Closing group: ' // trim(this%hdf5_dataset_name)
   546)   call printMsg(option)  
   547)   call h5gclose_f(grp_id,hdf5_err)  
   548)   option%io_buffer = 'Closing hdf5 file: ' // trim(this%filename)
   549)   call printMsg(option)  
   550)   call h5fclose_f(file_id,hdf5_err)
   551)   call h5close_f(hdf5_err)  
   552)   
   553)   call PetscLogEventEnd(logging%event_dataset_map_hdf5_read, &
   554)                         ierr);CHKERRQ(ierr)
   555)   
   556) end subroutine DatasetMapHDF5ReadMap
   557) #endif
   558) 
   559) ! ************************************************************************** !
   560) 
   561) subroutine DatasetMapHDF5Print(this,option)
   562)   ! 
   563)   ! Prints dataset info
   564)   ! 
   565)   ! Author: Glenn Hammond
   566)   ! Date: 10/22/13
   567)   ! 
   568) 
   569)   use Option_module
   570) 
   571)   implicit none
   572)   
   573)   class(dataset_map_hdf5_type), target :: this
   574)   type(option_type) :: option
   575)   
   576)   class(dataset_common_hdf5_type), pointer :: dataset_hdf5
   577) 
   578)   dataset_hdf5 => this
   579)   call DatasetCommonHDF5Print(this,option)
   580) 
   581)   if (len_trim(this%h5_dataset_map_name) > 0) then
   582)     write(option%fid_out,'(10x,''HDF5 Dataset Map Name: '',a)') &
   583)       trim(this%h5_dataset_map_name)
   584)   endif
   585)   if (len_trim(this%map_filename) > 0) then
   586)     write(option%fid_out,'(10x,''Map Filename: '',a)') &
   587)       trim(this%map_filename)
   588)   endif
   589)   write(option%fid_out,'(10x,''Global Dimensions: '',2i8)') &
   590)     this%map_dims_global(:)
   591)   
   592) end subroutine DatasetMapHDF5Print
   593) 
   594) ! ************************************************************************** !
   595) 
   596) subroutine DatasetMapHDF5Strip(this)
   597)   ! 
   598)   ! Strips allocated objects within Map dataset object
   599)   ! 
   600)   ! Author: Glenn Hammond
   601)   ! Date: 05/03/13
   602)   ! 
   603) 
   604)   use Utility_module, only : DeallocateArray
   605) 
   606)   implicit none
   607)   
   608)   class(dataset_map_hdf5_type) :: this
   609)   
   610)   call DatasetCommonHDF5Strip(this)
   611)   
   612)   call DeallocateArray(this%mapping)
   613)   call DeallocateArray(this%datatocell_ids)
   614)   call DeallocateArray(this%cell_ids_local)
   615)   
   616) end subroutine DatasetMapHDF5Strip
   617) 
   618) ! ************************************************************************** !
   619) 
   620) subroutine DatasetMapHDF5Destroy(this)
   621)   ! 
   622)   ! Destroys a dataset
   623)   ! 
   624)   ! Author: Glenn Hammond
   625)   ! Date: 05/29/13
   626)   ! 
   627) 
   628)   implicit none
   629)   
   630)   class(dataset_map_hdf5_type), pointer :: this
   631)   
   632)   if (.not.associated(this)) return
   633)   
   634)   call DatasetMapHDF5Strip(this)
   635)   
   636)   deallocate(this)
   637)   nullify(this)
   638)   
   639) end subroutine DatasetMapHDF5Destroy
   640) 
   641) end module Dataset_Map_HDF5_class

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