CGAL Mesh_3 ERROR: assertion violation! Expr: minimal_size_ > 0 || sq_d > 0

Question

My application uses CGAL to create a tetrahedral mesh. The input are six patches (each of them in the form of an OFF-file) forming the boundary of the domain as well as a file with the feature edges (which are the interfaces between the patches).

Thanks to the help I got here and here again, it has worked perfectly many times. However, recently I came across this seemingly inocuous example:

Note that the bottom patch is a mesh of a cylinder, i.e., the domain is not convex.

This fails with the following error message

Error: CGAL ERROR: assertion violation!
Expr: minimal_size_ > 0 || sq_d > 0
File: /path/to/CGAL-5.2/include/CGAL/Mesh_3/Protect_edges_sizing_field.h
Line: 750

I have activated the debugging output above the assertion. Just before crashing, it prints the following

Nearest power vertex of (16.1 2.74455 159.976) is 0x44fafe0 (16.1 2.74455 159.976 5.47731e-28) at distance: 8.08582e-28
Nearest power vertex of (16.1 2.74455 159.976) is 0x44fb3d0 (16.1 2.74455 159.976 4.84338e-28) at distance: 9.66355e-30
Nearest power vertex of (16.1 2.74455 159.976) is 0x44fb050 (16.1 2.74455 159.976 2.73616e-27) at distance: 3.23196e-27
Nearest power vertex of (16.1 2.74455 159.976) is 0x44fb590 (16.1 2.74455 159.976 7.82747e-30) at distance: 0

To see the point (16.1 2.74455 159.976), one has to zoom very close to the corner where the three feature polylines meet.

Question 1

What is wrong?

My guess is that while constructing the protecting balls, the required squared radius somehow falls down to something like 10^-60, which is considered zero. But that should not happen, right? I tried helping it by changing several input parameters (target edge size, surface tolerance, etc.) but to no avail.

Question 2 (if Question 1 cannot be answered)

Can you maybe point me to the place in the papers where exactly this assertion comes into play?

Thinking that I have decent enough knowledge of Delaunay triangulations, I have spent a few hours reading [7] and [8]. However, it would take much more time for me to understand the details of the algorithm and the implementation enough to find out what's wrong. So any help (even if it just restricts the area I should search) is appreciated.

Related questions

I have googled around but haven't found much.

• This user got the same error because his meshes were not co-refined (mine are).
• Another user got the same error after using detect_features() (whereas I use add_features()).
• There was also a discussion in the mailing list (I cannot find it right now), where someone has got the same error because his triangle had zero height (which is not my case).

Answer 1

As @sloriot correctly pointed out in the comments, the problem was that my feature polylines do not meet exactly.

To help future readers, here is the relevant part of my polyline file:

(...)
    1.6100000e+1    2.2933305e-1    1.5999984e+2
    1.6100000e+1    1.4022496e+0    1.5999386e+2
    1.6100000e+1    2.7445457e+0    1.5997646e+2
34
    1.6100000e+1    2.7445456e+0    1.5997646e+2
    1.4946531e+1    3.3592107e+0    1.5996473e+2
    1.4337449e+1    3.6838449e+0    1.5995759e+2
(...)

Notice the difference of the y-coordinates of the last point of the previous polyline and the first point of the second polyline.