On Mon, Dec 15, 2025 at 3:43 AM Aldo Bonfiglioli <[email protected]> wrote:
> On 12/2/25 15:02, Matthew Knepley wrote: > > On Mon, Dec 1, 2025 at 8:22 AM Aldo Bonfiglioli < > [email protected]> wrote: > >> Dear developers, >> >> I wrote a code that extracts subDMs corresponding to the various strata >> in the Face Sets. >> >> I run into troubles when I view a subDM or a Vec attached to the subDM >> using the VTK format. >> >> More precisely, the problem only occurs on more than one processor when >> the rank=0 processor has no points on a given subDM. >> >> For instance, when the attached reproducer is run on 2 procs, the >> u_01.vtu file (global u Vec mapped to the subDM corresponding to >> stratum=1) >> >> only includes the header, but no data. All other u_0?.vtu files can >> successfully be loaded and viewed in paraview. >> >> The problem does NOT arise when I view the same objects in HDF5 format. >> >> However, my problem in using the HDF5 lies in the fact that: >> >> while the hdf5 file obtained with DMView can be post-processed with >> "petsc/lib/petsc/bin/petsc_gen_xdmf.py" to create a xmf file readable by >> paraview >> > > Hi Aldo, > > Sorry about this, I would like to make it more intuitive. First, the > solution (I think) > > -dm_plex_view_hdf5_storage_version 1.1.0 > > will write the Viz field by default, so that PAraview will see it. Can you > try this? > > Why do we need this? I have now made version-controlled output formats. > There is something about > this in the manual, but not enough. Paraview only supports vertex-based > fields and cell-based fields > (at least that I understand), so we need to write a separate copy of the > field (since Plex supports any > layout). Lots of people do not want a separate copy, since they are > checkpointing, so we control this > with a format (PETSC_VIEWER_HDF5_VIZ). You can pass this for specific > output, or use the format > version that does it automatically. > > Let me know if this works. > > Thanks, > > Matt > > >> I do not know how to view the field(s) associated with the DM when the >> hdf5 file is obtained from VecView. >> >> The reproducer compiles with the latest petsc release. >> >> Thanks, >> >> Aldo >> >> -- >> Dr. Aldo Bonfiglioli >> Associate professor of Fluid Mechanics >> Dipartimento di Ingegneria >> Universita' della Basilicata >> V.le dell'Ateneo Lucano, 10 85100 Potenza ITALY >> tel:+39.0971.205203 fax:+39.0971.205215 >> web: >> https://urldefense.us/v3/__http://docenti.unibas.it/site/home/docente.html?m=002423__;!!G_uCfscf7eWS!eXE2csxUb1J5bnRosf9LIGxLs17TdstMxQcGs0mbXFroRjzLN81jVmeAWfMr41X8JGEuVm286lVTmDgmi5s1zfsfegJNuWVVWTM$ >> > > > -- > 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 > > https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!e3ZEg0fMK55wMqFxwA8I9trdBUMCc27ap7Q8V8beAiyV258mz8N43KSQa8MUBZMDUTVJw33Yay2EgEb4psDr$ > > <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!e3ZEg0fMK55wMqFxwA8I9trdBUMCc27ap7Q8V8beAiyV258mz8N43KSQa8MUBZMDUTVJw33Yay2EgGGjuvYF$ > > > > Matt, > > use of the option "-dm_plex_view_hdf5_storage_version 1.1.0" is indeed > necessary, > > but I also realized, thanks to a suggestion from > [email protected], that I have to View BOTH the dm and vec in > the same hdf5 file, i.e. > > ! > ! dump the dm+u to the same HDF5 file > ! > > filename = "test.h5" > PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, trim(filename), > FILE_MODE_WRITE, viewer, ierr)) > PetscCall(DMView(dm, viewer, ierr)) > PetscCall(PetscViewerDestroy(viewer, ierr)) > PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, trim(filename), > FILE_MODE_APPEND, viewer, ierr)) > PetscCall(VecView(u, viewer, ierr)) > PetscCall(PetscViewerDestroy(viewer, ierr)) > > > Once test.h5 has been processed with "petsc_gen_xdmf.py", I can load the > xmf file in paraview e see the solution (there is NO solution unless > -dm_plex_view_hdf5_storage_version 1.1.0 is in the options db). > > I was probably misled by the fact that a single VecView in VTK format > gives both the mesh and solution in the same file. > > Does this make sense? > Ah, yes. VTK does not have a way to construct the file without the DM, so we force it. Our HDF5 format can actually handle multiple DMs (so the adaptive refinement can be visualized), so you need to specify what to put in. > Final question: > > is it possible to specify, using command line options or the options db, > that vecview should be appended to an existing file? > > Yes. You use the mode -vec_view hdf5:sol.h5::append Thanks, Matt > Thanks, > > Aldo > > -- > Dr. Aldo Bonfiglioli > Associate professor of Fluid Mechanics > Dipartimento di Ingegneria > Universita' della Basilicata > V.le dell'Ateneo Lucano, 10 85100 Potenza ITALY > tel:+39.0971.205203 fax:+39.0971.205215 > web: > https://urldefense.us/v3/__http://docenti.unibas.it/site/home/docente.html?m=002423__;!!G_uCfscf7eWS!e3ZEg0fMK55wMqFxwA8I9trdBUMCc27ap7Q8V8beAiyV258mz8N43KSQa8MUBZMDUTVJw33Yay2EgLcMeEyC$ > > > -- 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 https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!e3ZEg0fMK55wMqFxwA8I9trdBUMCc27ap7Q8V8beAiyV258mz8N43KSQa8MUBZMDUTVJw33Yay2EgEb4psDr$ <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!e3ZEg0fMK55wMqFxwA8I9trdBUMCc27ap7Q8V8beAiyV258mz8N43KSQa8MUBZMDUTVJw33Yay2EgGGjuvYF$ >
