- Thank you received: 0
Bug? report: specialcf.normal orientation
4 years 1 month ago - 4 years 1 month ago #3219
by cpfeiler
Bug? report: specialcf.normal orientation was created by cpfeiler
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
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
Last edit: 4 years 1 month ago by cpfeiler. Reason: somehow attachment got lost first try
4 years 1 month ago #3220
by cpfeiler
[attachment=undefined]normal_bug.jpg[/attachment] [attachment=undefined]normal_bug.jpg[/attachment]
Replied by cpfeiler on topic Bug? report: specialcf.normal orientation
Attachment not found
4 years 1 month ago #3221
by cpfeiler
Replied by cpfeiler on topic Bug? report: specialcf.normal orientation
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()
geo.Add(big-small)
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")
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()
geo.Add(big-small)
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")
4 years 1 month ago #3223
by joachim
Replied by joachim on topic Bug? report: specialcf.normal orientation
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).
Joachim
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).
Code:
for i in range(ngmesh.GetNFaceDescriptors()):
fd = ngmesh.FaceDescriptor(i+1)
print (fd)
print (fd.domin, fd.domout)
Joachim
4 years 1 month ago #3224
by cpfeiler
Replied by cpfeiler on topic Bug? report: specialcf.normal orientation
Hey Joachim,
thanks for the quick reply.
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
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:
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
for some boundary faces. Anyhow, I still wanted to report this, since I could imagine other problems implied by the incosistent orientation.)
Best, Carl
thanks for the quick reply.
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
Code:
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
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()
geo.Add(big-small)
ngmesh = geo.GenerateMesh(maxh=0.5)
mesh = ngsMesh(ngmesh)
### new from here ###
from netgen.meshing import FaceDescriptor, Element2D, Mesh
new_ngmesh = Mesh()
fd_outside = new_ngmesh.Add(FaceDescriptor(bc=1,domin=1,surfnr=1))
pmap = {}
for e in mesh.ngmesh.Elements2D():
for v in e.vertices:
if (v not in pmap):
pmap[v] = new_ngmesh.Add(mesh.ngmesh[v])
for e in mesh.ngmesh.Elements2D():
new_ngmesh.Add(Element2D(fd_outside, [pmap[v] for v in e.vertices]))
new_ngmesh.GenerateVolumeMesh()
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
Code:
[pmap[v] for v in e.vertices]
Best, Carl
4 years 1 month ago #3225
by joachim
Replied by joachim on topic Bug? report: specialcf.normal orientation
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.
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.
Time to create page: 0.102 seconds