<div dir="ltr"><div dir="ltr">On Fri, Dec 30, 2022 at 5:59 AM Berend van Wachem <<a href="mailto:berend.vanwachem@ovgu.de">berend.vanwachem@ovgu.de</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Dear Mark,<br>
<br>
Yes, I have tried that. That will only work if I convert the DMForest to <br>
a DMPlex first, and then write the DMPlex to file. But then, I cannot <br>
read in a DMForest when I want to continue the calculations later on.<br>
<br>
Below is the code I use to write a DMPlex to file. When I call this <br>
routine with a DMForest, I get the error:<br></blockquote><div><br></div><div>This is true. I am going to have to look at this with Toby. I will get back to you.</div><div><br></div><div>  THanks,</div><div><br></div><div>     Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
[0]PETSC ERROR: Invalid argument<br>
[0]PETSC ERROR: Wrong type of object: Parameter # 1<br>
[0]PETSC ERROR: See <a href="https://petsc.org/release/faq/" rel="noreferrer" target="_blank">https://petsc.org/release/faq/</a> for trouble shooting.<br>
[0]PETSC ERROR: Petsc Release Version 3.18.2, Nov 28, 2022<br>
[0]PETSC ERROR: ./bin/write_periodic on a linux-gcc-dev named <br>
<a href="http://multiflow.multiflow.org" rel="noreferrer" target="_blank">multiflow.multiflow.org</a> by berend Fri Dec 30 11:54:51 2022<br>
[0]PETSC ERROR: Configure options --with-clean --download-metis=yes <br>
--download-parmetis=yes --download-hdf5 --download-p4est <br>
--download-triangle --download-tetgen --with-zlib-lib=/usr/lib64/libz.a <br>
--with-zlib-include=/usr/include --with-mpi=yes --with-mpi-dir=/usr <br>
--with-mpiexec=/usr/bin/mpiexec --download-slepc=yes --download-fftw=yes<br>
[0]PETSC ERROR: #1 PetscSectionGetChart() at <br>
/usr/local/petsc-3.18.2/src/vec/is/section/interface/section.c:592<br>
[0]PETSC ERROR: #2 DMPlexGetChart() at <br>
/usr/local/petsc-3.18.2/src/dm/impls/plex/plex.c:2702<br>
[0]PETSC ERROR: #3 DMPlexGetDepthStratum() at <br>
/usr/local/petsc-3.18.2/src/dm/impls/plex/plex.c:4882<br>
[0]PETSC ERROR: #4 DMPlexCreatePointNumbering() at <br>
/usr/local/petsc-3.18.2/src/dm/impls/plex/plex.c:8303<br>
[0]PETSC ERROR: #5 DMPlexTopologyView() at <br>
/usr/local/petsc-3.18.2/src/dm/impls/plex/plex.c:1854<br>
[0]PETSC ERROR: #6 WriteDMAndVectortoHDF5File() at <br>
/home/berend/git/Code/petsc-dmplex-restart-test/src/createandwrite.c:176<br>
[0]PETSC ERROR: #7 main() at <br>
/home/berend/git/Code/petsc-dmplex-restart-test/src/createandwrite.c:555<br>
<br>
The code I use to write is:<br>
<br>
PetscErrorCode WriteDMAndVectortoHDF5File(DM dm, Vec *DataVector, const <br>
char *HDF5file)<br>
{<br>
   PetscViewer H5Viewer;<br>
   PetscErrorCode ierr;<br>
   DM dmCopy;<br>
   PetscSection sectionCopy;<br>
   PetscInt i;<br>
   PetscScalar *xVecArray;<br>
   PetscInt numPoints, numPointsCopy;<br>
   Vec vectorCopy;<br>
   PetscScalar *array;<br>
<br>
   PetscFunctionBegin;<br>
   ierr = PetscViewerHDF5Open(PETSC_COMM_WORLD, HDF5file, <br>
FILE_MODE_WRITE, &H5Viewer);<br>
   CHKERRQ(ierr);<br>
   ierr = PetscObjectSetName((PetscObject) dm, "DM");<br>
   CHKERRQ(ierr);<br>
   ierr = PetscViewerPushFormat(H5Viewer, PETSC_VIEWER_HDF5_PETSC);<br>
   CHKERRQ(ierr);<br>
   ierr = DMPlexTopologyView(dm, H5Viewer);<br>
   CHKERRQ(ierr);<br>
   ierr = DMPlexLabelsView(dm, H5Viewer);<br>
   CHKERRQ(ierr);<br>
   ierr = DMPlexCoordinatesView(dm, H5Viewer);<br>
   CHKERRQ(ierr);<br>
   ierr = PetscViewerPopFormat(H5Viewer);<br>
   CHKERRQ(ierr);<br>
<br>
   ierr = DMClone(dm, &dmCopy);<br>
   CHKERRQ(ierr);<br>
   ierr = PetscObjectSetName((PetscObject) dmCopy, "DM");<br>
   CHKERRQ(ierr);<br>
   ierr = DMGetLocalSection(dm, &sectionCopy);<br>
   CHKERRQ(ierr);<br>
   ierr = DMSetLocalSection(dmCopy, sectionCopy);<br>
   CHKERRQ(ierr);<br>
<br>
   /* Write the section to the file */<br>
   ierr = DMPlexSectionView(dm, H5Viewer, dmCopy);<br>
   CHKERRQ(ierr);<br>
<br>
   ierr = DMGetGlobalVector(dmCopy, &vectorCopy);<br>
   CHKERRQ(ierr);<br>
<br>
   /*** We have to copy the vector into the new vector ... ***/<br>
   ierr = VecGetArray(vectorCopy, &array);<br>
   CHKERRQ(ierr);<br>
   ierr = VecGetLocalSize(*DataVector, &numPoints);<br>
   CHKERRQ(ierr);<br>
   ierr = VecGetLocalSize(vectorCopy, &numPointsCopy);<br>
   CHKERRQ(ierr);<br>
   assert(numPoints == numPointsCopy);<br>
   ierr = VecGetArray(*DataVector, &xVecArray);<br>
   CHKERRQ(ierr);<br>
<br>
   for (i = 0; i < numPoints; i++) /* Loop over all internal cells */<br>
   {<br>
     array[i] = xVecArray[i];<br>
   }<br>
<br>
   ierr = VecRestoreArray(vectorCopy, &array);<br>
   CHKERRQ(ierr);<br>
   ierr = VecRestoreArray(*DataVector, &xVecArray);<br>
   CHKERRQ(ierr);<br>
<br>
   ierr = PetscObjectSetName((PetscObject) vectorCopy, "DataVector");<br>
   CHKERRQ(ierr);<br>
<br>
   /* Write the vector to the file */<br>
   ierr = DMPlexGlobalVectorView(dm, H5Viewer, dmCopy, vectorCopy);<br>
   CHKERRQ(ierr);<br>
<br>
   /* Close the file */<br>
   ierr = PetscViewerDestroy(&H5Viewer);<br>
   CHKERRQ(ierr);<br>
<br>
   ierr = DMDestroy(&dmCopy);<br>
   CHKERRQ(ierr);<br>
   /*** End of writing ****/<br>
   PetscFunctionReturn(0);<br>
}<br>
<br>
<br>
<br>
<br>
<br>
On 29/12/2022 21:54, Mark Adams wrote:<br>
> Have you tried using the DMForest as you would use a DMPLex?<br>
> That should work. Forst has extra stuff but it just makes a Plex and the <br>
> end of the day.<br>
> <br>
> On Tue, Dec 27, 2022 at 5:21 AM Berend van Wachem <br>
> <<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a> <mailto:<a href="mailto:berend.vanwachem@ovgu.de" target="_blank">berend.vanwachem@ovgu.de</a>>> wrote:<br>
> <br>
>     Dear Petsc-team/users,<br>
> <br>
>     I am trying to save a DMForest which has been used for a calculation<br>
>     (so<br>
>     refined/coarsened in places) to a file and later read it from the file<br>
>     to continue a calculation with.<br>
>     My question is: how do I do this? I've tried a few things, such as<br>
>     using<br>
>     DMView and DMLoad, and saving the BaseDM, but that doesn't seem to save<br>
>     the refining/coarsening which has taken place in the initial<br>
>     calculation.<br>
>     I haven't been able to find any example or documentation for this. Any<br>
>     pointers or examples would be very much appreciated!<br>
> <br>
>     Best regards, Berend.<br>
> <br>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>