# Bug? report: specialcf.normal orientation

Hello everybody,

I’m curious whether it is a bug and/or known, that specialcf.normal(3) doesn’t necessarily point outwards for non-trivial geometries.

Attached you find a minimal (non) working example.
Best, Carl

https://ngsolve.org/media/kunena/attachments/972/normal_bug.py [attachment=undefined]normal_bug.jpg[/attachment] [attachment=undefined]normal_bug.jpg[/attachment]

Attachment doesn’t work…
Here the example:

import netgen.gui
from ngsolve import Draw, Mesh
from netgen.csg import OrthoBrick, Pnt, CSGeometry

big = OrthoBrick(Pnt(0,0,0), Pnt(4,4,4))
small = OrthoBrick(Pnt(1, 1, 2), Pnt(3, 3, 4))

geo = CSGeometry()
ngmesh = geo.GenerateMesh(maxh=0.5)
mesh = Mesh(ngmesh)

from ngsolve import specialcf
Draw(specialcf.normal(3), mesh, “nu”)
input(“drawn; enable Visual >> Draw Surface Vectors”)

Hi Carl,

that’s not a bug. The normal vector is defined outward for you ‘small’ brick.
Later you subtract it to get the subdomain, this does not change the normal vector.

You can use this code snipped to obtain the sub-domain numbers in front / behind the individual faces (taking care of 1-based enumeration in Netgen).

``````for i in range(ngmesh.GetNFaceDescriptors()):
fd = ngmesh.FaceDescriptor(i+1)
print (fd)
print (fd.domin, fd.domout)``````

Joachim

Hey Joachim,

Naturally, I suppose one wants to consider (orientation of) normal vectors with respect to the meshed domain. Of course, there opinions might differ, BUT:

Normal vectors not pointing outwards everywhere, (I suppose) is a consequence of boundary faces not being oriented consistently. And therefore, this simple example of big_brick minus small_brick above, breaks your tutorial on Working with meshes:

Since boundary faces are not oriented consistently, creating a volume mesh from the boundary mesh fails with the correct error message

``` ERROR: Edge 307 - 308 multiple times in surface mesh ... same error for some other edges ... netgen.libngpy._meshing.NgException: Stop meshing since surface mesh not consistent ```
Since boundary face orientation is incosistent, alle edges belonging to the boundary of the top face of the small brick appear twice in the mesh – instead of once, plus once in reverse direction!

Using the mesh above, and a snippet from the tutorial, reproduce error via:

[code]from ngsolve import Mesh as ngsMesh
from netgen.csg import OrthoBrick, Pnt, CSGeometry

big = OrthoBrick(Pnt(0,0,0), Pnt(4,4,4))
small = OrthoBrick(Pnt(1, 1, 2), Pnt(3, 3, 4))

geo = CSGeometry()
ngmesh = geo.GenerateMesh(maxh=0.5)
mesh = ngsMesh(ngmesh)

### new from here

from netgen.meshing import FaceDescriptor, Element2D, Mesh
new_ngmesh = Mesh()

pmap = {}
for e in mesh.ngmesh.Elements2D():
for v in e.vertices:
if (v not in pmap):

for e in mesh.ngmesh.Elements2D():
new_ngmesh.Add(Element2D(fd_outside, [pmap[v] for v in e.vertices]))

new_ngmesh.GenerateVolumeMesh()
[/code]

If it’s not a bug, there’s definitely an inconsistency.
(Nothing I couldn’t work around with the badhack of changing the list-order of

``[pmap[v] for v in e.vertices]``

for some boundary faces. Anyhow, I still wanted to report this, since I could imagine other problems implied by the incosistent orientation.)

Best, Carl

yes, that’s a good example where you can use the code snipped from my previous answer:

The surface element has an index pointing into the array of FaceDescriptors.
Whenever the FaceDescriptor has domout = 0, you use the orientation as is, if fd.domin = 0 you flip the orientation.

If both are non-zero, you found an element belonging to an internal interface, which we don’t want to copy to the new mesh.