[petsc-users] HDF5 vec output with DMPlex

[Adrian Croucher]

hi,

I am wanting to output time-stepped simulation results to an HDF5 file,
from a PETSc Vec created on a DMPlex using DMCreateGlobalVector().

I tried doing it the way it's done in
src/vec/vec/examples/tutorials/ex19.c, using PetscViewerHDF5Open(),
PetscViewerHDF5PushGroup() and PetscViewerHDF5SetTimestep(), but the
resulting HDF5 file only ever seems to have results for the last time
step in it.

It works OK for a Vec created simply using VecCreate(), as in the ex19
example. Is there any reason it shouldn't work with DMPlex?

Cheers, Adrian

[Matthew Knepley]

Post by Adrian Croucher
hi,
I am wanting to output time-stepped simulation results to an HDF5 file,
from a PETSc Vec created on a DMPlex using DMCreateGlobalVector().
I tried doing it the way it's done in
src/vec/vec/examples/tutorials/ex19.c, using PetscViewerHDF5Open(),
PetscViewerHDF5PushGroup() and PetscViewerHDF5SetTimestep(), but the
resulting HDF5 file only ever seems to have results for the last time step
in it.
It works OK for a Vec created simply using VecCreate(), as in the ex19
example. Is there any reason it shouldn't work with DMPlex?

No, that should work. This is how we do it in PyLith. There must be some
setup problem. When you look at the HDF5 file. Does
it have a dimension for timestep?

Thanks,

Matt

Post by Adrian Croucher
Cheers, Adrian
--
Dr Adrian Croucher
Senior Research Fellow
Department of Engineering Science
University of Auckland, New Zealand
tel: [phone]

[Adrian Croucher]

hi

Post by Matthew Knepley
No, that should work. This is how we do it in PyLith. There must be
some setup problem. When you look at the HDF5 file. Does
it have a dimension for timestep?

No, it doesn't.

Is there any extra setup I might need to do to get it to work? I've
attached a minimal example code, mesh and output.

Cheers, Adrian

[Matthew Knepley]

Post by Matthew Knepley
hi
No, that should work. This is how we do it in PyLith. There must be some
setup problem. When you look at the HDF5 file. Does
it have a dimension for timestep?
No, it doesn't.
Is there any extra setup I might need to do to get it to work? I've
attached a minimal example code, mesh and output.

Now I remember. I did not want the output to depend on the viewer.

Does your example work if you replace PetscViewerHDF5SetTimestep() with
DMSetOutputSequenceNumber()?

Thanks,

Matt

Post by Matthew Knepley
Cheers, Adrian
--
Dr Adrian Croucher
Senior Research Fellow
Department of Engineering Science
University of Auckland, New Zealand
tel: [phone]

[Adrian Croucher]

hi,

Post by Matthew Knepley
Now I remember. I did not want the output to depend on the viewer.
Does your example work if you replace PetscViewerHDF5SetTimestep()
with DMSetOutputSequenceNumber()?

Aha, yes, that appears to do the trick. It gives me the time values in
the output as well which was the next thing I was going to try and
figure out how to do.

Thanks...

Cheers, Adrian

[Adrian Croucher]

hi

When I read the HDF5 file back in again using PetscViewerHDF5Open() and
VecLoad(), is there any easy way to navigate to the last time step in
the results?

Or, equivalently, to find how many timesteps there are in there, so I
could then use PetscViewerHDF5SetTimestep() before doing VecLoad() ?

Cheers, Adrian

Post by Matthew Knepley
Now I remember. I did not want the output to depend on the viewer.
Does your example work if you replace PetscViewerHDF5SetTimestep()
with DMSetOutputSequenceNumber()?
Thanks,
Matt

[Matthew Knepley]

Post by Adrian Croucher
hi
When I read the HDF5 file back in again using PetscViewerHDF5Open() and
VecLoad(), is there any easy way to navigate to the last time step in the
results?
Or, equivalently, to find how many timesteps there are in there, so I
could then use PetscViewerHDF5SetTimestep() before doing VecLoad() ?

I am using

ierr = PetscViewerHDF5PushGroup(user.checkpoint,
"/fields");CHKERRQ(ierr);
ierr = DMSetOutputSequenceNumber(dm, user.restart, 0.0);CHKERRQ(ierr);
ierr = VecLoad(u, user.checkpoint);CHKERRQ(ierr);
ierr = DMOutputSequenceLoad(dm, user.checkpoint, "time",
user.restart, &time);CHKERRQ(ierr);
ierr = TSSetTime(ts, (PetscReal) time);CHKERRQ(ierr);
ierr = DMSetOutputSequenceNumber(dm, 0, 0.0);CHKERRQ(ierr);
ierr = PetscViewerHDF5PopGroup(user.checkpoint);CHKERRQ(ierr);
ierr = PetscViewerDestroy(&user.checkpoint);CHKERRQ(ierr);

but right now I do not have a query for the last timestep. We will have to
add an HDF5 query for this. I
am putting it on the list.

Thanks,

Matt

Post by Adrian Croucher
Cheers, Adrian
Now I remember. I did not want the output to depend on the viewer.
Does your example work if you replace PetscViewerHDF5SetTimestep() with
DMSetOutputSequenceNumber()?
Thanks,
Matt
--
Dr Adrian Croucher
Senior Research Fellow
Department of Engineering Science
University of Auckland, New Zealand
tel: [phone]