NGSolve / ngsolve

Netgen/NGSolve is a high performance multiphysics finite element software. It is widely used to analyze models from solid mechanics, fluid dynamics and electromagnetics. Due to its flexible Python interface new physical equations and solution algorithms can be implemented easily.
https://ngsolve.org/
GNU Lesser General Public License v2.1
436 stars 79 forks source link

Mesh.Refine() does not update mesh.faces properly #52

Closed Alex-Vasile closed 2 years ago

Alex-Vasile commented 2 years ago

The issue I'm having is that the information provided by mesh.faces (edges, points, and facets may also have this issue) is not update correctly after a call to mesh.Refine(). Simple test case included below.

from ngsolve import Mesh
from netgen.geom2d import unit_square

mesh = Mesh(unit_square.GenerateMesh())

assert len(mesh.faces[0].edges) > 0

mesh.Refine()

empty_faces = []
for face in mesh.faces:
    if len(face.edges) == 0:
        empty_faces.append(face)
assert len(empty_faces) == 0  # This will fail
sigmundholm commented 2 years ago

The empty faces does seem very strange. After refinement, the mesh in your example has eight elements/faces, which make sense. However, the iterable mesh.faces has 10 elements, where two of these are empty. The numbering of the edges make sense with respect to the position of the elements, if you disregard the two empty elements/faces.

schruste commented 2 years ago

Hi, this is not a bug. Old facets are kept after refinement although they are no longer geometrically directly involved in the new mesh. But for setting up the grid hierarchy in Multigrid methods this is needed. Indeed you can simply disregard the facets without neighbors, they are simply kept from the coarser mesh.

Alex-Vasile commented 2 years ago

Hi, this is not a bug. Old facets are kept after refinement although they are no longer geometrically directly involved in the new mesh. But for setting up the grid hierarchy in Multigrid methods this is needed. Indeed you can simply disregard the facets without neighbors, they are simply kept from the coarser mesh.

Okay. That makes sense. I make my workarounds permanent then. (I refine then immediately save to file and load from file)