[petsc-users] output DMDA to hdf5 file?

Matteo Semplice matteo.semplice at uninsubria.it
Mon Jul 12 10:40:36 CDT 2021


Dear all,

     I am experimenting with hdf5+xdmf output. At 
https://www.xdmf.org/index.php/XDMF_Model_and_Format I read that "XDMF 
uses XML to store Light data and to describe the data Model. Either 
HDF5[3] <https://www.hdfgroup.org/HDF5> or binary files can be used to 
store Heavy data. The data Format is stored redundantly in both XML and 
HDF5."

However, if I call DMView(dmda,hdf5viewer) and then I run h5ls or h5stat 
on the resulting h5 file, I see no "geometry" section in the file. How 
should I write the geometry to the HDF5 file?

Here below is what I have tried.

Best

     Matteo

//Setup

ierr = DMDACreate2d(PETSC_COMM_WORLD,
                               DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,
                               DMDA_STENCIL_STAR,
                               ctx.Nx,ctx.Ny, //global dim
                               PETSC_DECIDE,PETSC_DECIDE, //n proc on 
each dim
                               2,stWidth, //dof, stencil width
                               NULL, NULL, //n nodes per direction on 
each cpu
                               &(ctx.daAll));

ierr = DMDASetUniformCoordinates(ctx.daAll, 0.0, 1.0, 0.0, 1.0, 0.0, 
1.0); CHKERRQ(ierr);
ierr = DMDASetFieldName(ctx.daAll,0,"first"); CHKERRQ(ierr);
ierr = DMDASetFieldName(ctx.daAll,1,"second"); CHKERRQ(ierr);
ierr = DMDAGetLocalInfo(ctx.daAll,&ctx.daInfo); CHKERRQ(ierr);
ierr = DMCreateFieldDecomposition(ctx.daAll,NULL, NULL, &ctx.is, 
&ctx.daField); CHKERRQ(ierr);

//Initial data
ierr = DMCreateGlobalVector(ctx.daAll,&ctx.U0); CHKERRQ(ierr);
ierr = VecISSet(ctx.U0,ctx.is[0],1.); CHKERRQ(ierr);
ierr = VecISSet(ctx.U0,ctx.is[1],100.0); CHKERRQ(ierr);


Write mesh and Vec to hdf5:

PetscViewer viewer;
   ierr = 
PetscViewerHDF5Open(PETSC_COMM_WORLD,"solution.h5",FILE_MODE_WRITE,&viewer); 
CHKERRQ(ierr);
   ierr = DMView(ctx.daAll , viewer); CHKERRQ(ierr);  //does not output 
anything to solution.h5??
   ierr = VecView(ctx.U0,viewer); CHKERRQ(ierr);

ierr = PetscViewerDestroy(&viewer); CHKERRQ(ierr);

Attempt to save the two fields separately::

PetscViewer viewer;
ierr = 
PetscViewerHDF5Open(PETSC_COMM_WORLD,"solution.h5",FILE_MODE_WRITE,&viewer); 
CHKERRQ(ierr);
ierr = DMView(ctx.daField[0] , viewer); CHKERRQ(ierr); //does not output 
anything to solution.h5??

Vec uF;
ierr = VecGetSubVector(ctx.U,ctx.is[0],&uF); CHKERRQ(ierr);
PetscObjectSetName((PetscObject) uF, "first");
ierr = VecView(uF,viewer); CHKERRQ(ierr);
ierr = VecRestoreSubVector(ctx.U,ctx.is[0],&uF); CHKERRQ(ierr);

ierr = VecGetSubVector(ctx.U,ctx.is[1],&uF); CHKERRQ(ierr);
PetscObjectSetName((PetscObject) uF, "second");
ierr = VecView(uF,viewer); CHKERRQ(ierr);
ierr = VecRestoreSubVector(ctx.U,ctx.is[1],&uF); CHKERRQ(ierr);

ierr = PetscViewerDestroy(&viewer); CHKERRQ(ierr);

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210712/1cdd85c8/attachment.html>


More information about the petsc-users mailing list