[petsc-users] realCoords for DOFs
Yann Jobic
yann.jobic at univ-amu.fr
Tue Dec 13 07:03:28 CST 2022
That a really smart trick !
Thanks for sharing it.
I previously looked closeley to the DMProjectFunction without success
yet as it has the solution, but yours is just too good.
Yann
Le 12/13/2022 à 1:52 PM, Mark Adams a écrit :
> You should be able to use PetscFECreateDefault instead of
> PetscFECreateLagrange here and set the "order" on the command
> line, which is recommended, but start with this and low order to debug.
>
> Mark
>
> static PetscErrorCode crd_func(PetscInt dim, PetscReal time, const
> PetscReal x[], PetscInt Nf_dummy, PetscScalar *u, void *actx)
> {
> int i;
> PetscFunctionBeginUser;
> for (i = 0; i < dim; ++i) u[i] = x[i];
> PetscFunctionReturn(0);
> }
>
> PetscErrorCode (*initu[1])(PetscInt, PetscReal, const PetscReal [],
> PetscInt, PetscScalar [], void *);
> PetscFE fe;
>
> /* create coordinate DM */
> ierr = DMClone(dm, &crddm);CHKERRV(ierr);
> ierr = PetscFECreateLagrange(PETSC_COMM_SELF, dim, dim, PETSC_FALSE,
> order, PETSC_DECIDE, &fe);CHKERRV(ierr);
> ierr = PetscFESetFromOptions(fe);CHKERRV(ierr);
> ierr = DMSetField(crddm, 0, NULL, (PetscObject)fe);CHKERRV(ierr);
> ierr = DMCreateDS(crddm);CHKERRV(ierr);
> ierr = PetscFEDestroy(&fe);CHKERRV(ierr);
> /* project coordinates to vertices */
> ierr = DMCreateGlobalVector(crddm, &crd_vec);CHKERRV(ierr);
> initu[0] = crd_func;
> ierr = DMProjectFunction(crddm, 0.0, initu, NULL, INSERT_ALL_VALUES,
> crd_vec);CHKERRV(ierr);
> ierr = VecViewFromOptions(crd_vec, NULL, "-coord_view");CHKERRV(ierr);
>
> On Tue, Dec 13, 2022 at 3:38 AM Yann Jobic <yann.jobic at univ-amu.fr
> <mailto:yann.jobic at univ-amu.fr>> wrote:
>
>
>
> Le 12/13/2022 à 2:04 AM, Mark Adams a écrit :
> > PETSc does not store the coordinates for high order elements (ie,
> the
> > "midside nodes").
> > I get them by projecting a f(x) = x, function.
> > Not sure if that is what you want but I can give you a code
> snippet if
> > there are no better ideas.
>
> It could really help me !
> If i have the node coordinates in the reference element, then it's easy
> to project them to the real space. But i couldn't find a way to have
> them, so if you could give me some guidance, it could be really helpful.
> Thanks,
> Yann
>
> >
> > Mark
> >
> >
> > On Mon, Dec 12, 2022 at 6:06 PM Yann Jobic
> <yann.jobic at univ-amu.fr <mailto:yann.jobic at univ-amu.fr>
> > <mailto:yann.jobic at univ-amu.fr <mailto:yann.jobic at univ-amu.fr>>>
> wrote:
> >
> > Hi all,
> >
> > I'm trying to get the coords of the dofs of a DMPlex for a
> PetscFE
> > discretization, for orders greater than 1.
> >
> > I'm struggling to run dm/impls/plex/tutorials/ex8.c
> > I've got the following error (with option -view_coord) :
> >
> > [0]PETSC ERROR: --------------------- Error Message
> > --------------------------------------------------------------
> > [0]PETSC ERROR: Object is in wrong state
> > [0]PETSC ERROR: DMGetCoordinatesLocalSetUp() has not been called
> > [0]PETSC ERROR: See https://petsc.org/release/faq/
> <https://petsc.org/release/faq/>
> > <https://petsc.org/release/faq/
> <https://petsc.org/release/faq/>> for trouble shooting.
> > [0]PETSC ERROR: Petsc Release Version 3.18.2, unknown
> > [0]PETSC ERROR:
> >
> /home/jobic/projet/fe-utils/petsc/fetools/cmake-build-debug/ex_closure_petsc
> > on a named luke by jobic Mon Dec 12 23:34:37 2022
> > [0]PETSC ERROR: Configure options
> > --prefix=/local/lib/petsc/3.18/p3/gcc/nompi_hdf5 --with-mpi=0
> > --with-debugging=1 --with-blacs=1
> --download-zlib,--download-p4est
> > --download-hdf5=1 --download-triangle=1 --with-single-library=0
> > --with-large-file-io=1 --with-shared-libraries=0 -CFLAGS=" -g
> -O0"
> > -CXXFLAGS=" -g -O0" -FFLAGS=" -g -O0" PETSC_ARCH=nompi_gcc_hdf5
> > [0]PETSC ERROR: #1 DMGetCoordinatesLocalNoncollective() at
> >
> /home/devel/src_linux/petsc-3.18.0/src/dm/interface/dmcoordinates.c:621
> > [0]PETSC ERROR: #2 DMPlexGetCellCoordinates() at
> >
> /home/devel/src_linux/petsc-3.18.0/src/dm/impls/plex/plexgeometry.c:1291
> > [0]PETSC ERROR: #3 main() at
> >
> /home/jobic/projet/fe-utils/petsc/fetools/src/ex_closure_petsc.c:86
> > [0]PETSC ERROR: PETSc Option Table entries:
> > [0]PETSC ERROR: -dm_plex_box_faces 2,2
> > [0]PETSC ERROR: -dm_plex_dim 2
> > [0]PETSC ERROR: -dm_plex_simplex 0
> > [0]PETSC ERROR: -petscspace_degree 1
> > [0]PETSC ERROR: -view_coord
> > [0]PETSC ERROR: ----------------End of Error Message
> -------send entire
> > error message to petsc-maint at mcs.anl.gov----------
> >
> > Maybe i've done something wrong ?
> >
> >
> > Moreover, i don't quite understand the function
> DMPlexGetLocalOffsets,
> > and how to use it with DMGetCoordinatesLocal. It seems that
> > DMGetCoordinatesLocal do not have the coords of the dofs, but
> only the
> > nodes defining the geometry.
> >
> > I've made some small modifications of ex8.c, but i still have
> an error :
> > [0]PETSC ERROR: --------------------- Error Message
> > --------------------------------------------------------------
> > [0]PETSC ERROR: Invalid argument
> > [0]PETSC ERROR: Wrong type of object: Parameter # 1
> > [0]PETSC ERROR: WARNING! There are option(s) set that were
> not used!
> > Could be the program crashed before they were used or a spelling
> > mistake, etc!
> > [0]PETSC ERROR: Option left: name:-sol value: vtk:sol.vtu
> > [0]PETSC ERROR: See https://petsc.org/release/faq/
> <https://petsc.org/release/faq/>
> > <https://petsc.org/release/faq/
> <https://petsc.org/release/faq/>> for trouble shooting.
> > [0]PETSC ERROR: Petsc Release Version 3.18.2, unknown
> > [0]PETSC ERROR:
> >
> /home/jobic/projet/fe-utils/petsc/fetools/cmake-build-debug/view_coords
> > on a named luke by jobic Mon Dec 12 23:51:05 2022
> > [0]PETSC ERROR: Configure options
> > --prefix=/local/lib/petsc/3.18/p3/gcc/nompi_hdf5 --with-mpi=0
> > --with-debugging=1 --with-blacs=1
> --download-zlib,--download-p4est
> > --download-hdf5=1 --download-triangle=1 --with-single-library=0
> > --with-large-file-io=1 --with-shared-libraries=0 -CFLAGS=" -g
> -O0"
> > -CXXFLAGS=" -g -O0" -FFLAGS=" -g -O0" PETSC_ARCH=nompi_gcc_hdf5
> > [0]PETSC ERROR: #1 PetscFEGetHeightSubspace() at
> >
> /home/devel/src_linux/petsc-3.18.0/src/dm/dt/fe/interface/fe.c:1692
> > [0]PETSC ERROR: #2 DMPlexGetLocalOffsets() at
> >
> /home/devel/src_linux/petsc-3.18.0/src/dm/impls/plex/plexceed.c:98
> > [0]PETSC ERROR: #3 ViewOffsets() at
> > /home/jobic/projet/fe-utils/petsc/fetools/src/view_coords.c:28
> > [0]PETSC ERROR: #4 main() at
> > /home/jobic/projet/fe-utils/petsc/fetools/src/view_coords.c:99
> > [0]PETSC ERROR: PETSc Option Table entries:
> > [0]PETSC ERROR: -heat_petscspace_degree 2
> > [0]PETSC ERROR: -sol vtk:sol.vtu
> > [0]PETSC ERROR: ----------------End of Error Message
> -------send entire
> > error message to petsc-maint at mcs.anl.gov----------
> >
> >
> > Is dm/impls/plex/tutorials/ex8.c a good example for viewing
> the coords
> > of the dofs of a DMPlex ?
> >
> >
> > Thanks,
> >
> > Yann
> >
>
More information about the petsc-users
mailing list