Segfault when creating 1D mesh

Hey team,

I’m trying to build a 1D mesh with two subdomains, somewhat like 1.5 Spaces and forms on subdomains — NGS-Py 6.2.2302-87-ga5a5eff3b documentation

I was unable to find helper methods for 1D meshes, so I tried something more manual, like

from netgen.meshing import Element0D, Element1D, Element2D, MeshPoint, Pnt
from netgen.meshing import Mesh as NetGenMesh
import ngsolve

N = 10
ngmesh = NetGenMesh(dim=1)
pids = []
# Add points to the mesh.
for i in range(N + 1):
    pids.append(ngmesh.Add(MeshPoint(Pnt(i / N, 0, 0))))
# Left half of the domain is material 1, right half material 2.
for i in range(N):
    ngmesh.Add(Element1D([pids[i], pids[i + 1]], index=1 if i <= 5 else 2))
ngmesh.SetMaterial(1, "domain_left")
ngmesh.SetMaterial(2, "domain_right")

# Add BC to the mesh.
ngmesh.Add(Element0D(pids[0], index=1))  # Added index is 1
ngmesh.Add(Element0D(pids[N], index=2))
ngmesh.SetBCName(0, "l")  # Index is retrieved as 0
ngmesh.SetBCName(1, "r")

# Create NGSolve mesh and refine it once uniformly.
mesh = ngsolve.Mesh(ngmesh)
for el in mesh.Elements():
    mesh.SetRefinementFlag(el, True)
mesh.Refine()

Now, this code segfaults:

$ netgen mwe.py
NETGEN-6.2-dev
Developed by Joachim Schoeberl at
2010-xxxx Vienna University of Technology
2006-2010 RWTH Aachen University
1996-2006 Johannes Kepler University Linz
optfile ./ng.opt does not exist - using default values
togl-version : 2
loading ngsolve library
NGSolve-6.2.1907
Using Lapack
Including sparse direct solver UMFPACK
Running parallel using 8 thread(s)
importing NGSolve-6.2.1907
(should) load python file 'mwe.py'
 Mesh bisection
Segmentation fault: 11

I was unable to run it through valgrind (no support yet on newer MacOS systems, :-() so I have little idea what is going wrong. In any case, I have two questions:

  • how would I make this code work, or at least achieve the same goal?
  • In setting the boundary conditions, there is a discrepancy between the indices (see comment) – what am I doing wrong?

Thank you!!

Hi Jan,

thank you for the report …

using Netgen/NGSolve in 1D is not the main focus, and some things are a bit crooked.
Yes, we have this 0/1-based discrepancy, and the way out is to use labels wherever possible.

I just added an exception telling that mesh refinement in 1D is not supported, instead of segfaulting.

Joachim

Hi,

we have now an updated method to get rid of the 0/1-based indexing for the element-index.

We first create a 0/1/2/3-dimensional region, and get the region-index back.
This index is used for the elements. We should not worry whether it is 0 or 1-based:

idx_doml = ngmesh.AddRegion("domain_left", dim=1)
idx_domr = ngmesh.AddRegion("domain_right", dim=1)
    
for i in range(N):
    ngmesh.Add(Element1D([pids[i], pids[i + 1]], index=idx_doml if i < 5 else idx_domr))

the whole example is attached.

the tutorials are also updated.
you need the latest nightly build.

Joachim

Attachment: test1d.py