I think I found a hint in src/dm/impls/plex/examples/tests/ex5.c to get my code
corrected.
Apparently the contiguous numbering is
NOT: cells, vertices, edges, faces (as you said),
BUT: cells, vertices, faces, edges.
Then in numPoints, you don't put depth 0..3 or depth 3..0,
BUT: #vertices, #faces, #edges, #cells.
This seems completely random, but works in the attached example.
Maybe you should finally take the time to document this mess if even the
developers no longer understand it.
regards
Chris
Mar 12, 2019, 1:38 PM by [email protected]:
> On Tue, Mar 12, 2019 at 9:25 AM <> [email protected]
> <mailto:[email protected]>> > wrote:
>
>> Thank you,
>> I added DMPlexInterpolate(). Now I get a different error in
>> DMPlexInterpolate():
>> ...
>>
>
> Okay, you did not need interpolate. You already specified all the levels.
> However, your orientations are wrong.
> Get rid of that code.
>
>
>>
>> [0]PETSC ERROR: Dimension 0 not supported
>> ...
>>
>> The source code and full output are attached. Anybody able to fix this?
>>
>
> The problem is a numbering convention in the library. Plex will accept any
> consistent DAG. However, if you
> want to use other things in the library, like geometry routines, then there
> is a convention on numbering (which
> makes many things simpler). We require that you number contiguously:
>
> Cells: [0, Nc)
> Vertices: [Nc, Nc+Nv)
> Edges: [Nc+Nv, Nc+Nv+Ne)
> Faces: [Nc+Nv+Ne, Nc+Nv+Ne+Nf)
>
> You numbered the faces in the vertices slot, so the geometry routines got
> confused.
>
> Thanks,
>
> Matt
>
>
>>
>> regards
>> Chris
>>
>> --
>> Securely sent with Tutanota. Get your own encrypted, ad-free mailbox:
>> https://tutanota.com <https://tutanota.com>
>>
>>
>> Mar 12, 2019, 1:08 PM by >> [email protected] <mailto:[email protected]>>> :
>>
>>> On Tue, Mar 12, 2019 at 9:03 AM Chris Finn via petsc-users <>>>
>>> [email protected] <mailto:[email protected]>>>> > wrote:
>>>
>>>> Hello,
>>>> with the code below, I create a tetrahedron using DMPlexCreateFromDAG,
>>>> then I try to run DMPlexComputeCellGeometryFVM on this cell. The latter
>>>> call fails with output:
>>>>
>>>
>>> All the geometry stuff requires that you interpolate the mesh I think, Just
>>> use DMPlexInterpolate().
>>>
>>> Thanks,
>>>
>>> Matt
>>>
>>>
>>>>
>>>> ...
>>>> [0]PETSC ERROR: --------------------- Error Message -----------------
>>>> [0]PETSC ERROR: Argument out of range
>>>> [0]PETSC ERROR: Cannot handle faces with 1 vertices
>>>> ...
>>>> (full output is attached).
>>>>
>>>> What am I doing wrong?
>>>> regards
>>>> Chris
>>>>
>>>> Here is the code:
>>>> static char help[] = "No help \n";
>>>>
>>>> #include <petscdmplex.h>
>>>>
>>>> #undef __FUNCT__
>>>> #define __FUNCT__ "main"
>>>> int main(int argc,char **args)
>>>> {
>>>> PetscErrorCode ierr;
>>>> PetscInitialize(&argc,&args,(char*)0,help);
>>>>
>>>> DM dm;
>>>> int cStart,cEnd;
>>>> int fStart,fEnd;
>>>>
>>>> int depth = 3;
>>>> int dim = 3;
>>>>
>>>> PetscInt numPoints[4] = {1,4,6,4};
>>>> PetscInt coneSize[15] = {4,3,3,3,3,2,2,2,2,2,2,0,0,0,0};
>>>> PetscInt cones[28] = {1,2,3,4, 5,9,8, 9,6,10, 10,8,7,
>>>> 5,6,7, 11,12, 12,13, 13,11, 11,14, 12,14, 13,14};
>>>> PetscInt coneOrientations[28] = {0 };
>>>> PetscScalar vertexCoords[12] = {0,0,0, 1,0,0, 0,1,0, 0,0,1};
>>>>
>>>> DMCreate(PETSC_COMM_WORLD, &dm);
>>>> DMSetType(dm, DMPLEX);
>>>> DMSetDimension(dm,dim);
>>>> DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones,
>>>> coneOrientations, vertexCoords);
>>>>
>>>> DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);
>>>> DMPlexGetHeightStratum(dm, 1, &fStart, &fEnd);
>>>>
>>>> for (int k =cStart;k<cEnd;k++){
>>>> double vol;
>>>> double centroid[3];
>>>> double normal[3];
>>>> ierr = DMPlexComputeCellGeometryFVM(dm, k, &vol,
>>>> centroid,NULL);CHKERRQ(ierr);
>>>> printf("FVM: V=%f c=(%f %f %f) n=(%f %f
>>>> %f)\n",vol,centroid[0],centroid[1],centroid[2],
>>>> normal[0],normal[1],normal[2]);
>>>> }
>>>> ierr = PetscFinalize();
>>>> return 0;
>>>> }
>>>>
>>>>
>>>>
>>>
>>>
>>> --
>>> 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://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
>>>
>>
>>
>
>
> --
> 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://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
>
static char help[] = "No help \n";
#include <petscdmplex.h>
#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc,char **args)
{
PetscErrorCode ierr;
PetscInitialize(&argc,&args,(char*)0,help);
DM dm;
int cStart,cEnd;
int fStart,fEnd;
int depth = 3;
int dim = 3;
PetscInt numPoints[4] = {4,4,6,1};
//v f e c
PetscInt coneSize[15] = {4,0,0,0,0,3,3,3,3,2,2,2,2,2,2};
PetscInt cones[28] = {5,6,7,8,
9,13,12, 13,10,14, 14,12,11, 9,10,11,
1,2, 2,3, 3,1, 1,4, 2,4, 3,4};
PetscInt coneOrientations[28] = {0 };
PetscScalar vertexCoords[12] = {0,0,0, 1,0,0, 0,1,0, 0,0,1};
DMCreate(PETSC_COMM_WORLD, &dm);
DMSetType(dm, DMPLEX);
DMSetDimension(dm,dim);
ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords); CHKERRQ(ierr);
ierr = DMPlexGetDepth(dm, &depth);CHKERRQ(ierr);
ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr);
fprintf(stderr,"reconstructed depth,dim: %d %d\n",depth,dim);
DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);
DMPlexGetHeightStratum(dm, 1, &fStart, &fEnd);
for (int k =cStart;k<cEnd;k++){
double vol;
double centroid[3];
double normal[3];
ierr = DMPlexComputeCellGeometryFVM(dm, k, &vol, centroid,NULL);CHKERRQ(ierr);
printf("FVM: V=%f c=(%f %f %f) n=(%f %f %f)\n",vol,centroid[0],centroid[1],centroid[2],
normal[0],normal[1],normal[2]);
}
for (int k =fStart;k<fEnd;k++){
double vol;
double centroid[3];
double normal[3];
ierr = DMPlexComputeCellGeometryFVM(dm, k, &vol, centroid,normal);CHKERRQ(ierr);
printf("FVM: V=%f c=(%f %f %f) n=(%f %f %f)\n",vol,centroid[0],centroid[1],centroid[2],
normal[0],normal[1],normal[2]);
}
ierr = PetscFinalize();
return 0;
}
