Dear Matt,
please explain your previous answer ("The normal should be outward, unless
the face has orientation < 0") with respect to my example.
As I think, the example shows that the face normals for faces 11 and 12 are
pointing outward, although they have orientation=-2. For all other faces
the normal direction and their orientation sign agree with what you said.Thanks Ingo 2017-04-14 11:28 GMT+02:00 Ingo Gaertner <[email protected]>: > Thank you, Matt, > as you say, the face orientations do change sign when switching between > the two adjacent cells. (I confused my program output. You are right.) But > it seems not to be always correct to keep the normal direction for > orientation>=0 and invert it for orientation<0: > > I include sample code below to make my question more clear. The program > creates a HexBoxMesh split horizontally in two cells (at x=0.5) It produces > this output: > > "Face centroids (c) and normals(n): > face #008 c=(0.250000 0.000000 0.000000) n=(-0.000000 -0.500000 0.000000) > face #009 c=(0.750000 0.000000 0.000000) n=(-0.000000 -0.500000 0.000000) > face #010 c=(0.250000 1.000000 0.000000) n=(0.000000 0.500000 0.000000) > face #011 c=(0.750000 1.000000 0.000000) n=(0.000000 0.500000 0.000000) > face #012 c=(0.000000 0.500000 0.000000) n=(-1.000000 -0.000000 0.000000) > face #013 c=(0.500000 0.500000 0.000000) n=(1.000000 0.000000 0.000000) > face #014 c=(1.000000 0.500000 0.000000) n=(1.000000 0.000000 0.000000) > Cell faces orientations: > cell #0, faces:[8 13 10 12] orientations:[0 0 -2 -2] > cell #1, faces:[9 14 11 13] orientations:[0 0 -2 -2]" > > Looking at the face normals, all boundary normals point outside (good). > The normal of face #013 points outside with respect to the left cell #0, > but inside w.r.t. the right cell #1. > > Face 13 is shared between both cells. It has orientation 0 for cell #0, > but orientation -2 for cell #1 (good). > What I don't understand is the orientation of face 12 (cell 0) and of face > 11 (cell 1). These are negative, which would make them point into the cell. > Have I done some other stupid mistake? > > Thanks > Ingo > > Here is the code: > > static char help[] = "Check face normals orientations.\n\n"; > #include <petscdmplex.h> > #undef __FUNCT__ > #define __FUNCT__ "main" > > int main(int argc, char **argv) > { > DM dm; > PetscErrorCode ierr; > PetscFVCellGeom *cgeom; > PetscFVFaceGeom *fgeom; > Vec cellgeom,facegeom; > int dim=2; > int cells[]={2,1}; > int cStart,cEnd,fStart,fEnd; > int coneSize,supportSize; > const int *cone,*coneOrientation; > > ierr = PetscInitialize(&argc, &argv, NULL, help);CHKERRQ(ierr); > ierr = PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL);CHKERRQ(ierr); > ierr = DMPlexCreateHexBoxMesh(PETSC_COMM_WORLD, dim, > cells,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE, > &dm);CHKERRQ(ierr); > > ierr = DMPlexComputeGeometryFVM(dm, &cellgeom,&facegeom);CHKERRQ(ierr); > ierr = VecGetArray(cellgeom, (PetscScalar**)&cgeom);CHKERRQ(ierr); > ierr = VecGetArray(facegeom, (PetscScalar**)&fgeom);CHKERRQ(ierr); > ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr); > ierr = DMPlexGetHeightStratum(dm, 1, &fStart, &fEnd);CHKERRQ(ierr); > > fprintf(stderr,"Face centroids (c) and normals(n):\n"); > for (int f=fStart;f<fEnd;f++){ > fprintf(stderr,"face #%03d c=(%03f %03f %03f) n=(%03f %03f %03f)\n",f, > fgeom[f-fStart].centroid[0],fgeom[f-fStart].centroid[1],fgeo > m[f-fStart].centroid[2], > fgeom[f-fStart].normal[0],fgeom[f-fStart].normal[1],fgeom[f- > fStart].normal[2]); > } > > fprintf(stderr,"Cell faces orientations:\n"); > for (int c=cStart;c<cEnd;c++){ > ierr = DMPlexGetConeSize(dm, c, &coneSize);CHKERRQ(ierr); > ierr = DMPlexGetCone(dm, c, &cone);CHKERRQ(ierr); > ierr = DMPlexGetConeOrientation(dm, c, &coneOrientation);CHKERRQ(ierr > ); > if (dim==2){ > if (coneSize!=4){ > fprintf(stderr,"Expected coneSize 4, got %d.\n",coneSize); > exit(1); > } > fprintf(stderr,"cell #%d, faces:[%d %d %d %d] orientations:[%d %d > %d %d]\n",c, > cone[0],cone[1],cone[2],cone[3], > coneOrientation[0],coneOrientation[1],coneOrientation[2],con > eOrientation[3] > ); > } else if (dim==3){ > if (coneSize!=6){ > fprintf(stderr,"Expected coneSize 6, got %d.\n",coneSize); > exit(1); > } > fprintf(stderr,"cell #%d, faces:[%d %d %d %d %d %d] > orientations:[%d %d %d %d %d %d]\n",c, > cone[0],cone[1],cone[2],cone[3],cone[4],cone[5], > coneOrientation[0],coneOrientation[1],coneOrientation[2],con > eOrientation[3],coneOrientation[4],coneOrientation[5] > ); > } else { > fprintf(stderr,"Dimension %d not implemented.\n",dim); > exit(1); > } > } > ierr = PetscFinalize(); > > } > > > 2017-04-14 11:00 GMT+02:00 Matthew Knepley <[email protected]>: > >> On Wed, Apr 12, 2017 at 10:52 AM, Ingo Gaertner < >> [email protected]> wrote: >> >>> Hello, >>> I have problems determining the orientation of the face normals of a >>> DMPlex. >>> >>> I create a DMPlex, for example with DMPlexCreateHexBoxMesh(). >>> Next, I get the face normals using DMPlexComputeGeometryFVM(DM dm, Vec >>> *cellgeom, Vec *facegeom). facegeom gives the correct normals, but I don't >>> know how the inside/outside is defined with respect to the adjacant cells? >>> >> >> The normal should be outward, unless the face has orientation < 0. >> >> >>> Finally, I iterate over all cells. For each cell I iterate over the >>> bounding faces (obtained from DMPlexGetCone) and try to obtain their >>> orientation with respect to the current cell using >>> DMPlexGetConeOrientation(). However, the six integers for the orientation >>> are the same for each cell. I expect them to flip between neighbour cells, >>> because if a face normal is pointing outside for any cell, the same normal >>> is pointing inside for its neighbour. Apparently I have a misunderstanding >>> here. >>> >> >> I see the orientations changing sign for adjacent cells. Want to send a >> simple code? You should see this >> for examples. You can run SNES ex12 with -dm_view ::ascii_info_detail to >> see the change in sign. >> >> Thanks, >> >> Matt >> >> >>> How can I make use of the face normals in facegeom and the orientation >>> values from DMPlexGetConeOrientation() to get the outside face normals for >>> each cell? >>> >>> Thank you >>> Ingo >>> >>> >>> <https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail> >>> Virenfrei. >>> www.avast.com >>> <https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail> >>> <#m_2360740100696784902_m_-1897645241821352774_m_-6432344443275881332_DAB4FAD8-2DD7-40BB-A1B8-4E2AA1F9FDF2> >>> >> >> >> >> -- >> 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 >> > >
