[petsc-users] dmplex face normals orientation

Fri Apr 14 04:28:12 CDT 2017

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],
fgeom[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],
coneOrientation[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],
coneOrientation[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 <knepley at gmail.com>:

> On Wed, Apr 12, 2017 at 10:52 AM, Ingo Gaertner <
> ingogaertner.tus at gmail.com> 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_-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
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20170414/b0bdd93f/attachment-0001.html>