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

[Matthew Knepley]

On Fri, Mar 20, 2015 at 12:32 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:

>
> > On Mar 20, 2015, at 1:17 PM, Matthew Knepley <knepley at gmail.com> wrote:
> >
> > On Fri, Mar 20, 2015 at 12:08 PM, Barry Smith <bsmith at mcs.anl.gov>
> wrote:
> >
> > > On Mar 20, 2015, at 11:31 AM, Håkon Strandenes <haakon at hakostra.net>
> wrote:
> > >
> > > Yes, in the HDF5 file. Currently a dof=1 DMDA Vec is saved as a NZ x
> NY x NX dataset, wile a dof > 1 is saved as NZ x NY x NX x dof in the HDF5
> file.
> >
> >   Actually a 1d DMDA is stored as a NX dataset, a 2d DMDA is stored as a
> NY x NX data set and a 3d DMDA is stored as a NZ x NY x NX data set and an
> extra dimension is added if dof > 1.   And I think this is great and
> intuitively obvious for users.
> >
> >    I agree with Hakon and advocate NOT sticking the extra dimension in
> when bs is 1. Here is Matt's entire commit
> >
> >
> >    maxDims[dim] = dims[dim];
> >    chunkDims[dim] = dims[dim];
> >    ++dim;
> > -  if (bs > 1) {
> > +  if (bs >= 1) {
> >      dims[dim]    = bs;
> >      maxDims[dim] = dims[dim];
> >      chunkDims[dim] = dims[dim];
> >    }
> >    count[dim] = PetscHDF5IntCast(xin->map->n)/bs;
> >    ++dim;
> > -  if (bs > 1) {
> > +  if (bs >= 1) {
> >      count[dim] = bs;
> >      ++dim;
> >    }
> >    }
> >    offset[dim] = PetscHDF5IntCast(low/bs);
> >    ++dim;
> > -  if (bs > 1) {
> > +  if (bs >= 1) {
> >      offset[dim] = 0;
> >      ++dim;
> >    }
> >
> >   (Then there is Karl's one line fix for reading HDF5 files.)
> >
> >   Arguing that handling bs = 1 as a special case greatly (or even
> modestly) increases the complexity of the PETSc source code is nonsense.
> >
> >   So I want a stronger case presented then
> >
> > >> My rationale was to simplify our code path, and have one way to
> handle fields with components. I did not
> > >> want to special case a field with one component.
> >
> > If we special case here, we must special case in everything we do, such
> as automatic generation of Xdmf.
> > I would rather say that the PETSc format always has the bs in it, which
> makes it easier to process and
> > consistent.
>
>    From what I see the "easier" is just trivially easier; having a tiny
> bit more complicated code that is intuitively clear for users seems a
> better choice.
>
> > However, we allowing reading of of a set without block size.
>
>   Hmm, according to the bug report which I confirmed, VecLoad() refuses to
> load in a 1d dataset
>

Yes, I was saying this is the part we change.


> >
> > We would have the ability to write without blocksize by setting the
> blocksize to -1.
>
>   Ok, so now you are still supporting that "special case" but with some
> totally non-intuitive value for a block size that will crash if used with
> the vector in any other usage of the vector. How is this either simpler
> code or more intuitive code.
>

I am supporting a special case in this function for output, but not for our
accepted format. I want
to avoid a special case in the format, not in the function itself. Then I
can use the format from
other tools without making them complicated too.


>  As Hakon notes also, VecView from a DMDA is handled differently with
> regard to the bs dimension then VecView from a simple "vector". Which is
> terrible.


Yes, I think we should change it so that all vectors come with the
blocksize.

  Matt


>
> Barry
>
> >
> >    Matt
> >
> >   otherwise I think we should undo Matt's change.
> >
> >   Barry
> >
> >
> >
> > >
> > > In my humble opinion, it is best to leave the DMDA's as it is, and fix
> up the code related to the plain Vecs. That code is already full of "if (bs
> >= 1) dim++" etc. which is useless (as Barry points out).
> > >
> > > My final point is that one should handle this in a CONSISTENT manner
> for both "plain Vecs" and "DMDA Vecs", i.e. decide wither one want to
> always add one dimension for bs and dof or leave that dimension out for the
> special case with bs=1/dof=1.
> > >
> > > Håkon
> > >
> > >
> > > On 20. mars 2015 16:57, Jed Brown wrote:
> > >> Håkon Strandenes <haakon at hakostra.net> writes:
> > >>> I mean that if you decide to stick with the current way of
> > >>> reading/writing a "plain 1D" Vec with N elements as a N x 1
> > >>> array/dataset, you should treat a dof=1 DMDA Vec as an NZ x NY x NX
> x 1
> > >>> array/dataset.
> > >>
> > >> In the HDF5 file?  I agree.
> > >>
> >
> >
> >
> >
> > --
> > What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which their
> experiments lead.
> > -- Norbert Wiener
>
>


-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150320/ae46aacd/attachment-0001.html>