Odd bugs with Constructive Solid Geometry 2D

Main Forum Old-Forum - do not post here
1

I’ve noticed some odd behavior when generating 2D meshes using Constructive Solid Geometry. Consider the following code, which should generate the intersection of two overlapping circles:

import netgen.gui
from netgen.geom2d import CSG2d, Circle
from ngsolve import Mesh

geo = CSG2d()
c0 = Circle(center=(0,0),radius=1, bc='a')
c1 = Circle(center=(1,0),radius=1, bc='b')
domain = c0 * c1
geo.Add(domain)
mesh = Mesh(geo.GenerateMesh())

We instead get a single circle centered at at the origin. The same thing happens with a union. An intersection yields an empty mesh; it’s like the two circles are the same. If we shift the center a small amount, the code works as expected.

I’ve done this using a conda install on my Macbook Pro as well as with an apt-get install on Ubuntu 20.04 LTS, both with the same result.

2

Hello,

I cannot reproduce the issue with Netgen 6.2.2009-69-gea7f6c1e.
Which version are you running? If it’s older than 6.2.2009, I suggest an update.

Best,
Matthias

3

Thanks for the response. I am running 6.2.2009 on both my Mac and on my Ubuntu server. I’ve converted a Jupyter notebook to a webpage so you can see my version and full code here:

https://marksmath.org/ProblematicMeshIllustration.html

4

Here’s a similar problem. To produce this, I set up a virtual machine running Ubuntu 20.04 LTS on Digital Ocean. Python version is 3.8. I then logged in as root and followed almost exactly the Ubuntu instructions on this page. - thus:

apt-add-repository universe
add-apt-repository ppa:ngsolve/ngsolve
apt-get update
apt-get install ngsolve

(Since I was root, I didn’t need sudo.)

I then used my favorite, console based text editor to create a file named “ng_bad_script.py” with the following content:

#!/usr/bin/python3
from ngsolve import Mesh
from netgen.geom2d import Solid2d, Circle, CSG2d

dom = Solid2d([(0,0),(1,0),(1,1),(0.5,0.5)])-Circle((0.5,0.5),0.2)
geo = CSG2d()
geo.Add(dom)
mesh = Mesh(geo.GenerateMesh())

If I run this script, meshing then fails:

root@netgen:~# python3 ng_bad_script.py 
importing NGSolve-6.2.2009
 Generate Mesh from spline geometry
 Boundary mesh done, np = 30
 CalcLocalH: 30 Points 0 Elements 0 Surface Elements
 Meshing domain 1 / 1
 load internal triangle rules
 give up with qualclass 201
 number of frontlines = 27
 Surface meshing done
Traceback (most recent call last):
  File "ng_script.py", line 9, in <module>
    mesh = Mesh(geo.GenerateMesh())
netgen.libngpy._meshing.NgException: meshing failed

The situation is worse if I delete the point (0.5,0.5) from the Solid2d object and store the program in “ng_really_bad_script.py”, for then I enter an infinite loop:

root@netgen:~# python3 ng_really_bad_script.py 
importing NGSolve-6.2.2009
 Generate Mesh from spline geometry
 Boundary mesh done, np = 30
 CalcLocalH: 30 Points 0 Elements 0 Surface Elements
 Meshing domain 1 / 1
 load internal triangle rules
loclines.Size = 1
loclines.Size = 1
loclines.Size = 1
...
...

Note that loclines.Size=1 printout repeats until I terminate.

Again, this is a totally fresh install and the exact same thing happens on my Mac with an Anaconda based installation.

To be clear, though, I’m generally quite impressed with the software and I’m looking forward to using it to generate examples for the PDE class that I’ll be teaching this semester. Examples where the points are in generic position tend to work; thus, problems like the above tend to go away if the points are perturbed a bit.

5

Thanks for the detailed information. I discovered the reason for this bug in the source code. I will let you know, when it’s fixed.

Best,
Matthias

6

The issues are fixed on the master branch. Nightly builds will be available tomorrow.

Best,
Matthias

7

I look forward to trying it - thanks!!

8

I wasn’t able to get to this right away but I can confirm that this fixed all the problems I was having.

Thank you!!

9

Hello there, just to let you know I hit the same bug with the very simple MWE below. This is with NGSolve 6.2.210x with x=1,2,3,4.

import netgen.gui
from netgen.geom2d import SplineGeometry
geo = SplineGeometry()

# Control points for a (near-)circle drawn with quadratic rational splines
pts = [[(297.79, 347.97), (280.5658347837522, 348.1496523099468)],
       [(280.47875, 366.16625000000005),
        (280.39167008489085, 384.1818404370362)],
       [(297.48, 384.52000000000004), (315.803620573473, 384.49016911587523)],
       [(315.75124999999997, 366.35), (315.6988839573919, 348.211400289232)]]

# Create NGSolve point objects
ngpts = []
for pt in pts:
    ngpts.append([geo.AppendPoint(*lopt) for lopt in pt])

# Connect the end point of every spline with the starting point of the next
for i, pt in enumerate(ngpts[:-1]):
    ngpt.append(ngpts[i + 1][0])
ngpts[-1].append(ngpts[0][0])

# Create splines
for pt in ngpts:
    geo.Append(["spline3", *pt])

# This crashes and prints "loclines.Size = 1" over and over.
geo.GenerateMesh()
10

You are having a different issue here (not using 2d CSG, but SplineGeometry).

Since you are passing the points in clockwise order, the mesher is meshing the outside (infinite) domain.

Either pass the points in counterclockwise order or define left/rightdomain (0 == dont mesh)

# Create splines
for pt in ngpts:
    geo.Append(["spline3", *pt], leftdomain=0, rightdomain=1)

Best,
Matthias

11

Gurh you are right. Thanks a bunch for the swift reply!!!