[petsc-users] Load Vec from 1D HDF5 dataset MATT READ THIS EMAIL!

[Barry Smith]

  MATT READ THIS EMAIL!

  Ok here is the scoop.  In  VecView_MPI_HDF5() is the code (use meta x vc-annotate in emacs)

df09a399 src/vec/vec/impls/mpi/pdvec.c (Matthew Knepley      2011-01-20 16:27:23 -0600  682)   if (bs >= 1) {
272345b0 src/vec/vec/impls/mpi/pdvec.c (Karl Rupp            2013-01-26 20:27:50 -0600  683)     dims[dim]      = bs;
272345b0 src/vec/vec/impls/mpi/pdvec.c (Karl Rupp            2013-01-26 20:27:50 -0600  684)     maxDims[dim]   = dims[dim];
61ea93a1 src/vec/vec/impls/mpi/pdvec.c (Brad Aagaard         2010-11-06 15:31:18 -0700  685)     chunkDims[dim] = dims[dim];
43c8710a src/vec/vec/impls/mpi/pdvec.c (Matthew Knepley      2010-11-05 19:21:58 -0500  686)     ++dim;
5464dd29 src/vec/vec/impls/mpi/pdvec.c (Matthew Knepley      2010-10-31 19:04:52 -0500  687)   }  if (bs >= 1) {

  then in PetscErrorCode VecLoad_HDF5(Vec xin, PetscViewer viewer)

  f6e5521d src/vec/vec/utils/vecio.c     (Karl Rupp            2013-01-28 13:32:07 -0600 277) if (bs >= 1) ++dim;

I check the commit message for df09a399 at bitbucket.com/petsc/petsc/ and it is 

Small change to HDF5 output (now Vec always has blocksize dim)

So what is happening is even the simplest vector (with a block size of 1) is treated as a 2d HDF5 object 

Matt, what was the rational for making vectors always have a dimension of block size even when it is one? Seems bad to me. Maybe makes HDF5 readers simpler since they don't need to handle bs 1 differently?

 Barry




> On Mar 19, 2015, at 12:31 PM, Håkon Strandenes <haakon at hakostra.net> wrote:
> 
> I'm sorry I didn't include this in the first email. Here is an example HDF5-file, together with a sample program that demonstrate the problem.
> 
> I have been trying to dig into VecLoad_HDF5() in vecio.c in order to suggest a solution myself, however, I do not think I fully understand the concept of the Vec block size and how that affects how the Vec is stored in memory and on disk (HDF5 file).
> 
> Regards,
> Håkon Strandenes
> 
> 
> On 19. mars 2015 18:17, Barry Smith wrote:
>> 
>>   Please email the file in.h5 so we can try to reproduce the problem and debug it. Are you using petsc version 3.5.?
>> 
>>   Barry
>> 
>>> On Mar 19, 2015, at 4:19 AM, Håkon Strandenes <haakon at hakostra.net> wrote:
>>> 
>>> Hi,
>>> 
>>> I am again in trouble with PETSc and its HDF5-related functions.
>>> 
>>> I have a code in which I want to read a 1D dataset from a HDF5 file into a PETSc Vec. The array are "small", and each MPI process should read the entire dataset and store it in its own sequential Vec.
>>> 
>>> Currently my code is as follows (somewhat simplified):
>>> 
>>> VecCreate(PETSC_COMM_SELF, &v);
>>> VecSetType(v, VECSEQ);
>>> PetscObjectSetName((PetscObject)v, "x");
>>> PetscViewerHDF5Open(PETSC_COMM_WORLD, "in.h5", FILE_MODE_READ, &viewer);
>>> VecLoad(v, viewer);
>>> 
>>> This works, but only if the dataset in the HDF5 file is a 2D dataset, of size Nx1, where N is the number of elements I read. If the dataset in the HDF5 file is of size N (i.e. a 1D dataset) this does not work.
>>> 
>>> The error message I get is the following:
>>> [0]PETSC ERROR: Unexpected data in file
>>> [0]PETSC ERROR: Dimension of array in file 1 not 2 as expected
>>> etc. and it points back to the line 296 in file vecio.c, function VecLoad_HDF5().
>>> 
>>> Technically this is not a problem. I could of course just go ahead and "reshape" all my datasets into 2D datasets. I do however find this behaviour very un-logical. Why can't I read a 1D dataset into a Vec?
>>> 
>>> So my question are:
>>> 1) Do I do anything wrong when creating my Vec and reading my datasets?
>>> 2) Is there anything I can do to read my 1D datasets into Vecs?
>>> 3) Could PETSc perhaps be a little less restrictive in these cases?
>>> 
>>> Thanks in advance.
>>> 
>>> Regards,
>>> Håkon Strandenes
>> 
>> 
> <VecLoadTest.tar.gz>