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

[Håkon Strandenes]

Did you come to any conclusion on this issue?

Regards,
Håkon

On 20. mars 2015 22:02, Håkon Strandenes wrote:
> On 20. mars 2015 20:48, Barry Smith wrote:
>> Why is 1 dimension a special case that is not worthy of its own
>> format? The same thing would hold for 2d and 3d. One could then argue
>> that we should have a single six dimensional format for the files for
>> all vectors that PETSc produces. Then a 1d problem has five of the
>> dimensions being 1.
>
> This is a very good point, and support my view.
>
> Let me come with two very simple example cases:
>
>
> Case 1:
> Create a list of grid points in an external preprocessor for the purpose
> of loading this into a Vec later:
>
>    x = np.linspace(0.0, 1.0, num=21)
>    f.create_dataset('/MESH/nodes/x', data=x)
>
> vs.
>
>    x = np.linspace(0.0, 1.0, num=21)
>    x = x.reshape((21,1))
>    f.create_dataset('/MESH/nodes/x', data=x)
>
>
> Case 2:
> Read three Vecs written to disk by PETSc, and calculate total "bounding
> box volume" of the grid:
>
>    g = h5py.File('grid.h5', 'r')
>    x = g['/MESH/nodes/x']
>    y = g['/MESH/nodes/y']
>    z = g['/MESH/nodes/z']
>    Vol = (xp[-1] - xp[0])*(yp[-1] - yp[0])*(zp[-1] - zp[0])
>
> vs.
>
>    g = h5py.File('grid.h5', 'r')
>    x = g['/MESH/nodes/x'][:,0]
>    y = g['/MESH/nodes/y'][:,0]
>    z = g['/MESH/nodes/z'][:,0]
>    Vol = (x[-1] - x[0])*(y[-1] - y[0])*(z[-1] - z[0])
>
>
> In both cases I think handling this extra, unnecessary dimension makes
> the code less attractive. It's not that either way is difficult,
> problematic or impossible, but it's just that 1D Vecs should intuitively
> be 1D datasets, and not 2D, 3D or 6D. This seriously confused me for
> quite a while until I figured this out, even after having written an
> entire Navier-Stokes DNS solver using the PETSc library for everything
> except time integration and filling these simple 1D coordinate arrays!
>
> Regards,
> Håkon