Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

FacetBasis outward normal vectors are false #1158

Closed
TNBaran opened this issue Aug 8, 2024 · 9 comments · Fixed by #1159
Closed

FacetBasis outward normal vectors are false #1158

TNBaran opened this issue Aug 8, 2024 · 9 comments · Fixed by #1159

Comments

@TNBaran
Copy link

TNBaran commented Aug 8, 2024

Hello
I am new user, and I began with 8.1 version, I encountered the following issue.
FacetBasis didn't work well since upgrading to 10.0.0-3, the outward normal vectors at boundary are sometimes false.
Here an example in order to show this issue, using an annular domain, there are 2 boundaries 'Int' and 'Ext'.
In G_Vh normal vectors are OK on two boundary 'Int' and 'Ext'
In Gm_Vh normal vectors are OK on boundary 'Ext'
In Gu_Vh normal vectors are not OK on boundary 'Int', this could be shown on the figure only two vectors are given on the circular boundary (red points) and in Gu_Vh.normals.value there are only two normal vectors.
array([[[ 0.89926879, 0.89926879],
[-0.89926879, -0.89926879],
........
[ 0.89926879, 0.89926879],
[-0.89926879, -0.89926879]],

   [[-0.43739644, -0.43739644],
    [ 0.43739644,  0.43739644],
    .......
    [-0.43739644, -0.43739644],
    [ 0.43739644,  0.43739644]]])

Here attached the mesh file created with gmsh
Anneau01.zip

import os
print(os.environ['PATH'])
from pathlib import Path
import numpy as np
from skfem import *
from skfem.helpers import dot, grad
from skfem.io import from_meshio
import matplotlib.pyplot as plt

mesh = Mesh.load('Anneau01.msh')
element = ElementTriP1()
Vh = Basis(mesh, element)
G_Vh = FacetBasis(mesh, element)
Gm_Vh = FacetBasis(mesh, element, facets=mesh.boundaries['Ext'])
Gu_Vh = FacetBasis(mesh, element, facets=mesh.boundaries['Int'])

plt.figure(1)
plt.plot(G_Vh.normals.value[0],G_Vh.normals.value[1],'',color='blue')
plt.plot(Gm_Vh.normals.value[0],Gm_Vh.normals.value[1],'
',color='green')
plt.plot(Gu_Vh.normals.value[0],Gu_Vh.normals.value[1],'*',color='red')
plt.show()

@kinnala
Copy link
Owner

kinnala commented Aug 8, 2024

Is it fixed by 10.0.1 bugfix or a different issue?

@kinnala
Copy link
Owner

kinnala commented Aug 8, 2024

This is the release from few days ago:

[10.0.1] - 2024-08-06

Fixed: Mesh.load returned incorrect orientation for some Gmsh meshes with tagged interfaces

@kinnala
Copy link
Owner

kinnala commented Aug 9, 2024

It seems to be a different issue, if you do mesh.draw(boundaries=True).show() you can see the following:

image

It seems that the element associated with the internal boundary is not correct. I wonder if there is something that has been misunderstood concerning the Gmsh orientation.

@kinnala
Copy link
Owner

kinnala commented Aug 9, 2024

Before this is solved, in this case the workaround is the following:

mesh = MeshTri.load('Anneau01.msh', ignore_orientation=True)

which avoids the orientation of the boundaries according to the Gmsh logic. There is another logic in native scikit-fem meshes which assures that on exterior boundaries normal vector is always correct.

@TNBaran
Copy link
Author

TNBaran commented Aug 9, 2024

Thanks,
it is OK when using this workaround.
However, there is something abnormal, it worked fine with version 8.1 and if you examine normal vectors obtained by G_Vh = FacetBasis(mesh, element) for all boundaries, they are OK for all boundaries, so Mesh.load returned correct orientation.
I think that may be it's facets option that's the issue.
Best regards

@kinnala
Copy link
Owner

kinnala commented Aug 9, 2024

Yes, it seems that mesh.boundaries['Int'].ori is incorrect. This is a variable which is supposed to be inferred from the Gmsh file in skfem.io.meshio.from_meshio.

The value of the variable is inferred incorrectly while reading your example mesh, and I will need to study in detail why this happens. There was a clear bug in 10.0.0 but as we already fixed one closely related bug in 10.0.1 there must be another one remaining.

@kinnala
Copy link
Owner

kinnala commented Aug 9, 2024

I have not been able to find Gmsh documentation on the topic of normal vectors so I am just finding many of these so-called "conventions" related to the orientation of normal vectors by trial-and-error.

@kinnala
Copy link
Owner

kinnala commented Aug 9, 2024

Alright, I think there is actually no remaining bug in skfem.io.meshio.from_meshio. It is trying to do an outward normal for a closed loop as it is supposed to do but now the issue is that the mesh has no elements inside of the closed loop, and orientation for such boundary does not make sense.

@kinnala
Copy link
Owner

kinnala commented Aug 9, 2024

Thanks for the report, this has been fixed in #1159. The fix will be included in the next release.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging a pull request may close this issue.

2 participants