Offset method robustness for complex shapes about gdspy

[rosivda]

Furthermore to @okianus comment above (and related to issue #1), gdspy offset algorithm does not have an expected behavior. For example, suppose we have two structures that are initially separated but, for a sufficiently large offset number, would merge into a single structure. That is not what is happening, as the gdspy offset method treats each structure separately and enlarges them individually, but only if they are simple, preferably convex hull shapes. Knowing this, one of the solutions to @okianus problem is to fragment complex shapes into a set of simpler shapes and subsequently using 'OR' to unite the shapes. This is not a generic solution though. Because of this gdspy offset algorithm limitation, we are forced to deploy external algorithms for our production layouts that handle arbitrary shapes and unite them into a smaller set of resulting shapes. One note here is that we use offsets not only to make small adjustments to the shape but also to produce large exclusion zones for further boolean operations. It often happens that the amount of offset exceeds the size of the feature, which is exactly when gdspy offset algorithm fails, as demonstrated in @okianus example. Error messages:
Ubuntu 14.04 LTS: Process finished with exit code 139
Windows 10: Process finished with exit code -1073741819 (0xC0000005)

[heitzmann]

@okianus : This problem is definitely a bug in the clipping code. To be honest, it might be a while before I can fix it, because those are usually difficult to find. I have to check...

@rosivda : The issue you're describing comes from the fact that the offset function was initially thought of a function operating in individual polygons. Indeed the idea for whole sets of polygons was to merge them afterwards using a boolean OR, as you mention. I cannot see a different approach to simplify this process, but if you have any ideas they'll be very welcome! Otherwise, once the problem above has been fixed, this solution should at least work.

[rosivda]

@heitzmann : I am afraid that the problem that we are facing is considerably more complex. The fundamental flaw of the present offset() implementation is that it assumes that:

OFFSET( polygon1 OR polygon2) == OFFSET(polygon1) OR OFFSET(polygon2).

This is incorrect for adjacent or overlapping shapes. I have created an example to illustrate this. Please try zero and non-zero offset. You will observe that, even though the shape on the right is the same thing as the two shapes on the left, the outcome of the offset() algorithm leads to different results. This creates a failure opportunity for any shape that consists of numerous polygons.

import gdspy
import numpy

cell = gdspy.Cell('CELL')

pl1 = numpy.array([[ -1.,  0.], [0., 1.],  [0., 0.]])
pl2 = numpy.array([[ 0., 0.], [ 0., 1.], [ 1., 0.]])
pl3 = numpy.array([[ 2., 0.], [ 3., 1.], [4., 0.]])

cell.add(gdspy.Polygon(pl1))
cell.add(gdspy.Polygon(pl2))
cell.add(gdspy.Polygon(pl3))

offset_result = gdspy.offset(gdspy.CellReference(cell), distance=0.2)
offset_cell = gdspy.Cell('CELL_OFFSET')
offset_cell.add(offset_result)
offset_cell.add(cell)

gdspy.LayoutViewer(cells=offset_cell, depth=3)

[rosivda]

@heitzmann : I would like to share a bit more general insight about gdspy - from the industrial perspective. gdspy has now been used by our R&D team to create large-scale production layouts, start-to-end. We have adopted 'natively' gdspy object structure and architecture with its approach that uses a dictionary with unique naming convention and cell referencing. gdspy has proven to be robust and scalable to large projects. It also yielded nearly theoretically-optimal gdsii layout files. gdspy has allowed us to bring all of our design tools into a single IDE.
Besides @okianus fixes, there are two areas that continue limiting gdspy from wider adoption as a full-fledged and completely self-contained module:
(1) fundamental implementation issues of offset(), as described above
(2) lack of a more powerful viewer that is able to handle large layouts
While we use separate tools to visualize layout files, the lack of robust offset() algorithm has been, perhaps, the only issue that continues to disrupt our workflow.

[heitzmann]

I've changed how the offset computation is performed in commit bb43726. The previous implementation was trying to be smart about the shapes but it lacked generality. The algorithm now should be general and more mathematically rigorous with the definition of offset.
The price for that is the computational cost, which has increased significantly due to the new rounded corners, whereas before they were single points.
Also note that similar to boolean, offset also has an eps parameter that can help in case of segfaults.

@rosivda's example now gives the same output in both cases for positive offsets (expansion). Negative offsets (shrinking) take into account internal boundaries, and that is by design. If you don't want internal boundaries showing up in the result, you have to merge all polygons beforehand. Merging the geometry is a good idea for expansion as well, since it reduces the size of the problem.

Please let me know if the new offset works and I'll make an official release.

[okianus]

Thank you for the quick response! Will try the new offset with my test cases as soon as possible.

[okianus]

I tried the modified version of offset on a typical layout. The layout contains some simple but realistic device structures. Prior to the latest commit, the offset operation (enlargement of a collection of structures by a relatively small amount) took approximately 0.5 seconds. After installing the latest commit bb43726, the operation takes an unrealistically long time (about 460 seconds with eps = 1e-9, while the default eps = 1e-13 results in a segmentation fault after about 7 minutes).

[rosivda]

@heitzmann : Thanks for looking into this critical issue. Your proposed offset algorithm is based on a rounding corner approach. This may be able to cover some of the aspects of the design but is not generally suitable for our design flow. We rely on miter-style offset rather than rounding. I would not recommend putting this algo in a new release of gdspy. I have no illusions about the complexity of a generic offset algorithm but have encountered successful implementations from vendors that are able to handle all of our cases, even the most complex ones. The link below highlights the complexity of the algorithm but I would not know the exact algorithms that our external packages are based on:
http://doc.cgal.org/Manual/3.2/doc_html/cgal_manual/Straight_skeleton_2/Chapter_main.html

[heitzmann]

Yes, I thought the new algorithm might end up being too slow for any practical use. The better solution seems to be using an external algorithm dedicated to polygon clipping and offset. I'm currently considering clipper, but I'm open to suggestions as long as they are compatible with the gdspy's GPL licence.

[rosivda]

I have looked at the examples of the clipper algorithm - they look quite general with no obvious limitations. However, no one in our group has tried it in practice. Should you decide to proceed with its implementation, we would be glad to test our libraries with your new algorithm and provide you with feedback. Please note that our individual chips (not a wafer) typically contain tens of thousands polygons, and some of them hundreds of thousands polygons. This would provide a hard test for the robustness of the algorithm.

[heitzmann]

Sounds good! I'll let you know when I have something to test.

[heitzmann]

@rosivda, @okianus: The clipper library integration is coming slowly. Meanwhile, I've tried to fix the offset function in the master branch to provide the 3 joint types: round (still slow), miter and bevel. In particular, miter should be similar to what the previous offset function used to do (if you set the tolerance high enough), but without the errors you first reported.
If you have time to check it out, it worked fine with my (small) test cases.

[okianus]

Thank you for the recent update. I've tried the updated offset function using the default joint type (miter), and it no longer generates segmentation errors for the letter C case above (once enlarged by distance>=0.18 so that features are lost).

The algorithm still lacks generality though when looking at a case like presented above by @rosivda or when enlarging letter C by a larger amount (e.g. distance=2). I think the skeleton approach to construct an offset curve of a polygon, as described in the link provided by @rosivda, is a rigorous approach that would be able to handle all cases. We look forward to a more generic implementation of the offset algorithm. Thank you for your work on this.

[okianus]

@heitzmann, I have recently tested the new offset algorithm in release 0.8.1 and compared it to the older algorithm in release 0.7.1. With default parameters passed to the offset() function, it is substantially slower than the older algorithm (by at least a factor of 6). Is this expected? If so, is it possible to keep the option of the older (more limited, but faster) algorithm in the upcoming releases?

[heitzmann]

@okianus, unfortunately that version of offset does not implement a proper offset, as discussed above, so I don't think it is a good idea to leave it available in the official release. If it is only for your use, you could always keep it in a separate module for yourself and import it as needed.

[okianus]

@heitzmann, could you please advise what optional parameters (joint type, tolerances, etc.) in the latest release would resolve the problem described by @rosivda above?

If I understand correctly, you mentioned in your comment on Feb 3 that the new implementation of offset resolves this problem, and the clipper library follows a similar, but more efficient algorithm. However, I've ran some basic experiments with the different joint types, and I still see different outputs whether there is one isosceles triangle, or the same shape formed by two adjacent triangles for positive offsets.

In the case of this simple problem, I could obtain the desired output by running fast_boolean prior to offset. However, for large scale problems, fast_boolean would result in re-segmentation. The segmented structures, when enlarged through offset, would still produce an incorrect result. Would appreciate your input on this.

[heitzmann]

@okianus, the problem you have is that offset is performed independently in each polygon. So if your shape is composed of many polygons (as in the case of 2 adjacent triangles) each will be offset independent of the other, so their joints might misbehave. The only setting which will give similar results for a single polygon vs a shape made of adjacent polygons is round joins (that is, if you don't care about overlaps).

I can think of only one feasible solution for this problem, which would be to add a flag in offset to allow it to internally join all polygons before the actual offset is performed. That will have to be implemented in the c extension, though, because that is the only place where polygons with holes can be represented without segmentation.

That should not be hard to do, so as soon as I have the time, I'll try to implement it.

UPDATE

@okianus could you test out commit 81a37ea and see if it solves your issue? I've added a join_first flag to offset that should be all you'll need to set. Internally it runs offset twice, so it will be a little slower unfortunately.

[okianus]

@heitzmann I've tested commit 81a37ea with simple cases that had problems similar to what was described above by @rosivda. With join_first=True, all the cases produced the correct result for positive offset. Thank you for your quick action on this!

A related question -- I've noticed that the performance of the offset algorithm depends on the distance parameter (for example, enlarging the same set of features by distance=20 takes roughly twice the amount of time compared to distance=0.3). Do you know what is the reason for this? (In contrast, the older algorithm in release 0.7.1 did not depend on the value of the required offset.)

[rosivda]

@heitzmann Thanks a lot for providing improvements to the offset algorithm! The join_first=True tag in your 81a37ea commit appears to have been a great leap. Preliminary testing indicates that the offset algorithm is now functioning correctly for adjacent geometries; we will now do testing of this algorithm in production environment. The performance issue has been noted but it is not yet clear if it is critical for practical-size layouts.

[okianus]

@heitzmann, could you please offer some insight as to why the execution time of the new offset algorithm increases for larger values of distance (with join_first=True)? Intuitively, larger distance would cause more polygons to merge, reducing the feature set to be passed to the clipper algorithm and presumably make the enlargement job easier, but that's not what I am observing.

[heitzmann]

Does it also happen when join_first = False? This flag forces the offset algorithm to run first with distance = 0 to perform all merging, and then with the actual distance. So, it seems strange to me that larger distances would have this effect only when join_first = True.

[okianus]

@heitzmann, I stand corrected, my previous comment was incorrect. The observation of longer execution times for larger distance values is with join_first=False (not with join_first=True as I had originally written). Here are some rough numbers for a large test case:

distance=0.45, join_first=False: 76 sec
distance=4.5, join_first=False: 444 sec
distance=45, join_first=False: 22472 sec

Even though there is no merging involved, why does the algorithm depend on the value of distance? When join_first=True is used, the times for the three cases above are about the same.

[rosivda]

@heitzmann, We have conducted testing of the offset algorithm for polygon sets containing about million vertices each (in join_first=True configuration). The algorithm produces expected results. The introduction of the offset algorithm has dramatically simplified our workflow, which now relies only on gdspy. We took a note that the offset algorithm was the longest (time-wise) operation for us but it was acceptable for the assembly of even the most complex chips. We would like to thank you for implementing it.
This has been discussed with @okianus. Our joint recommendation is to close this issue.

[heitzmann]

Sounds great! Thank you for your input and taking the time to improve gdspy!