diff --git a/src/db/db/dbTriangle.cc b/src/db/db/dbTriangle.cc index 490f27aa1..bb7c73910 100644 --- a/src/db/db/dbTriangle.cc +++ b/src/db/db/dbTriangle.cc @@ -494,4 +494,45 @@ Triangle::contains (const db::DPoint &point) const return res; } +double +Triangle::min_edge_length () const +{ + double lmin = edge (0)->d ().length (); + for (int i = 1; i < 3; ++i) { + lmin = std::min (lmin, edge (i)->d ().length ()); + } + return lmin; +} + +double +Triangle::b () const +{ + double lmin = min_edge_length (); + auto cr = circumcircle (); + return lmin / cr.second; +} + +bool +Triangle::has_segment () const +{ + for (int i = 0; i < 3; ++i) { + if (edge (i)->is_segment ()) { + return true; + } + } + return false; +} + +unsigned int +Triangle::num_segments () const +{ + unsigned int n = 0; + for (int i = 0; i < 3; ++i) { + if (edge (i)->is_segment ()) { + ++n; + } + } + return n; +} + } diff --git a/src/db/db/dbTriangle.h b/src/db/db/dbTriangle.h index 89a4f7bc6..ff39ab837 100644 --- a/src/db/db/dbTriangle.h +++ b/src/db/db/dbTriangle.h @@ -497,6 +497,26 @@ public: return mp_e1.get () == e || mp_e2.get () == e || mp_e3.get () == e; } + /** + * @brief Returns the minimum edge length + */ + double min_edge_length () const; + + /** + * @brief Returns the min edge length to circumcircle radius ratio + */ + double b () const; + + /** + * @brief Returns a value indicating whether the triangle borders to a segment + */ + bool has_segment () const; + + /** + * @brief Returns the number of segments the triangle borders to + */ + unsigned int num_segments () const; + private: bool m_is_outside; tl::weak_ptr mp_e1, mp_e2, mp_e3; diff --git a/src/db/db/dbTriangles.cc b/src/db/db/dbTriangles.cc index b3e6a057d..dd04c56cf 100644 --- a/src/db/db/dbTriangles.cc +++ b/src/db/db/dbTriangles.cc @@ -26,6 +26,7 @@ #include "dbWriter.h" #include "tlStream.h" #include "tlLog.h" +#include "tlTimer.h" #include @@ -1389,16 +1390,12 @@ Triangles::create_constrained_delaunay (const db::Region ®ion, double dbu) static bool is_skinny (db::Triangle *tri, const Triangles::TriangulateParameters ¶m) { - if (param.b < db::epsilon) { + if (param.min_b < db::epsilon) { return false; } else { - auto cr = tri->circumcircle (); - double lmin = tri->edge (0)->d ().length (); - for (int i = 1; i < 3; ++i) { - lmin = std::min (lmin, tri->edge (i)->d ().length ()); - } - double delta = (fabs (lmin / cr.second) + fabs (param.b)) * db::epsilon; - return lmin / cr.second < param.b - delta; + double b = tri->b (); + double delta = (b + param.min_b) * db::epsilon; + return b < param.min_b - delta; } } @@ -1410,43 +1407,28 @@ static bool is_invalid (db::Triangle *tri, const Triangles::TriangulateParameter double amax = param.max_area; if (param.max_area_border > db::epsilon) { - bool at_border = false; - for (int i = 0; i < 3; ++i) { - if (tri->edge (i)->is_segment ()) { - at_border = true; - } - } - if (at_border) { + if (tri->has_segment ()) { amax = param.max_area_border; } } if (amax > db::epsilon) { double a = tri->area (); - double delta = (fabs (a) + fabs (amax)) * db::epsilon; + double delta = (a + amax) * db::epsilon; return tri->area () > amax + delta; } else { return false; } } -static unsigned int num_segments (db::Triangle *tri) -{ - unsigned int n = 0; - for (int i = 0; i < 3; ++i) { - if (tri->edge (i)->is_segment ()) { - ++n; - } - } - return n; -} - void Triangles::triangulate (const db::Region ®ion, const TriangulateParameters ¶meters, double dbu) { + tl::SelfTimer timer (tl::verbosity () > parameters.base_verbosity, "Triangles::triangulate"); + create_constrained_delaunay (region, dbu); - if (parameters.b < db::epsilon && parameters.max_area < db::epsilon && parameters.max_area_border < db::epsilon) { + if (parameters.min_b < db::epsilon && parameters.max_area < db::epsilon && parameters.max_area_border < db::epsilon) { // no refinement requested - we're done. remove_outside_triangles (); @@ -1464,7 +1446,9 @@ Triangles::triangulate (const db::Region ®ion, const TriangulateParameters &p while (nloop < parameters.max_iterations) { // @@@ ++nloop; - tl::info << "Iteration " << nloop << " .."; + if (tl::verbosity () >= parameters.base_verbosity + 10) { + tl::info << "Iteration " << nloop << " .."; + } std::list > to_consider; for (auto t = new_triangles.begin (); t != new_triangles.end (); ++t) { @@ -1472,12 +1456,15 @@ Triangles::triangulate (const db::Region ®ion, const TriangulateParameters &p to_consider.push_back (*t); } } - tl::info << "@@@ to_consider " << to_consider.size(); if (to_consider.empty ()) { break; } + if (tl::verbosity () >= parameters.base_verbosity + 10) { + tl::info << to_consider.size() << " triangles to consider"; + } + new_triangles.clear (); for (auto t = to_consider.begin (); t != to_consider.end (); ++t) { @@ -1493,8 +1480,10 @@ Triangles::triangulate (const db::Region ®ion, const TriangulateParameters &p if ((*t)->contains (center) >= 0) { // heuristics #1: never insert a point into a triangle with more than two segments - if (num_segments (t->get ()) <= 1) { - // @@@ print(f"Inserting in-triangle center {repr(center)} of {repr(tri)}") + if (t->get ()->num_segments () <= 1) { + if (tl::verbosity () >= parameters.base_verbosity + 20) { + tl::info << "Inserting in-triangle center " << center.to_string () << " of " << (*t)->to_string (true); + } insert_point (center, &new_triangles); } @@ -1514,7 +1503,9 @@ Triangles::triangulate (const db::Region ®ion, const TriangulateParameters &p if (! edge->is_segment () || edge->side_of (*vstart) * edge->side_of (center) >= 0) { - // @@@ print(f"Inserting out-of-triangle center {repr(center)} of {repr(tri)}") + if (tl::verbosity () >= parameters.base_verbosity + 20) { + tl::info << "Inserting out-of-triangle center " << center << " of " << (*t)->to_string (true); + } insert_point (center, &new_triangles); } else { @@ -1522,7 +1513,9 @@ Triangles::triangulate (const db::Region ®ion, const TriangulateParameters &p double sr = edge->d ().length () * 0.5; db::DPoint pnew = *edge->v1 () + edge->d () * 0.5; - // @@@ print(f"split edge {repr(edge)} at {repr(pnew)}") + if (tl::verbosity () >= parameters.base_verbosity + 20) { + tl::info << "Split edge " << edge->to_string (true) << " at " << pnew.to_string (); + } db::Vertex *vnew = insert_point (pnew, &new_triangles); auto vertexes_in_diametral_circle = find_points_around (vnew, sr); @@ -1537,7 +1530,9 @@ Triangles::triangulate (const db::Region ®ion, const TriangulateParameters &p } } - // @@@ print(f" -> found {len(to_delete)} vertexes to remove") + if (tl::verbosity () >= parameters.base_verbosity + 20) { + tl::info << " -> found " << to_delete.size () << " vertexes to remove"; + } for (auto v = to_delete.begin (); v != to_delete.end (); ++v) { remove (*v, &new_triangles); } @@ -1548,10 +1543,11 @@ Triangles::triangulate (const db::Region ®ion, const TriangulateParameters &p } - // @@@ dump ("debug2.gds"); // @@@ - } + if (tl::verbosity () >= parameters.base_verbosity + 20) { + tl::info << "Finishing .."; + } remove_outside_triangles (); } diff --git a/src/db/db/dbTriangles.h b/src/db/db/dbTriangles.h index 2cd5d8d62..c246ea544 100644 --- a/src/db/db/dbTriangles.h +++ b/src/db/db/dbTriangles.h @@ -43,16 +43,17 @@ public: struct TriangulateParameters { TriangulateParameters () - : b (1.0), + : min_b (1.0), max_area (0.0), max_area_border (0.0), - max_iterations (std::numeric_limits::max ()) + max_iterations (std::numeric_limits::max ()), + base_verbosity (30) { } /** * @brief Min. readius-to-shortest edge ratio */ - double b; + double min_b; /** * @brief Max area or zero for "no constraint" @@ -68,6 +69,11 @@ public: * @brief Max number of iterations */ size_t max_iterations; + + /** + * @brief The verbosity level above which triangulation reports details + */ + int base_verbosity; }; typedef tl::shared_collection triangles_type; @@ -114,6 +120,9 @@ public: /** * @brief Creates a refined Delaunay triangulation for the given region + * + * The database unit should be chosen in a way that target area values are "in the order of 1". + * The algorithm becomes numerically unstable area constraints below 1e-4. */ void triangulate (const db::Region ®ion, const TriangulateParameters ¶meters, double dbu = 1.0); diff --git a/src/db/unit_tests/dbTrianglesTests.cc b/src/db/unit_tests/dbTrianglesTests.cc index 6bfaa1825..591e842a8 100644 --- a/src/db/unit_tests/dbTrianglesTests.cc +++ b/src/db/unit_tests/dbTrianglesTests.cc @@ -222,7 +222,7 @@ TEST(Triangles_illegal_edge2) db::Triangle t2 (&ee1, &ee2, &e1); - EXPECT_EQ (db::Triangles::is_illegal_edge (&e1), true); + EXPECT_EQ (db::Triangles::is_illegal_edge (&e1), false); } } @@ -631,35 +631,33 @@ TEST(Triangles_test_triangulate) r -= r2; db::Triangles::TriangulateParameters param; - param.b = 1.21; + param.min_b = 1.2; param.max_area = 1.0; - param.max_area_border = 0.0; db::Triangles tri; tri.triangulate (r, param, 0.001); - tri.dump ("debug.gds"); + EXPECT_EQ (tri.check (), true); + + for (auto t = tri.begin (); t != tri.end (); ++t) { + EXPECT_EQ (t->area () <= param.max_area, true); + EXPECT_EQ (t->b () >= param.min_b, true); + } + + EXPECT_EQ (tri.num_triangles () > 100 && tri.num_triangles () < 150, true); + tl::info << tri.num_triangles (); + + param.min_b = 1.0; + param.max_area = 0.1; + + tri.triangulate (r, param, 0.001); EXPECT_EQ (tri.check (), true); + for (auto t = tri.begin (); t != tri.end (); ++t) { + EXPECT_EQ (t->area () <= param.max_area, true); + EXPECT_EQ (t->b () >= param.min_b, true); + } + + EXPECT_EQ (tri.num_triangles () > 900 && tri.num_triangles () < 1000, true); } - -#if 0 - -assert tris.check(check_delaunay = False) - -amax = 0.0 -l2rmin = 2.0 -for tri in tris.triangles: -if not tri.is_outside: - _, radius = tri.circumcircle() - lmin = min([math.sqrt(t.square(s.d())) for s in tri.edges()]) - l2rmin = min(l2rmin, lmin / radius) - amax = max(amax, tri.area()) -print(f"max. area = {'%.5g'%amax}") -print(f"l/R min = {'%.5g'%l2rmin}") - -tris.dump_as_gdstxt("out.txt") - - -#endif