GCC Code Coverage Report

Directory: ./
Coverage: low: ≥ 0% medium: ≥ 75.0% high: ≥ 90.0%
Coverage Exec / Excl / Total
Lines: 98.1% 105 / 0 / 107
Functions: 100.0% 9 / 0 / 9
Branches: 96.3% 104 / 0 / 108
src/geo3d/ds/cvh/lib.hpp
Line Branch Exec Source
1 #pragma once
2
3 #include "../../../math/kahan/lib.hpp"
4 #include "../../dis/pl/lib.hpp"
5 #include "../../dis/plp_s/lib.hpp"
6 #include "../../dis/pp/lib.hpp"
7 #include "../../rel/plp/lib.hpp"
8 #include "../pl/lib.hpp"
9
10 namespace tifa_libs {
11
12 template <class FP>
13 class cvh3d {
14 struct facet {
15 // NOLINTNEXTLINE(modernize-avoid-c-arrays)
16 u32 n[3], id, vistime = 0;
17 bool isdel = false;
18 planev<FP> p;
19 // NOLINTNEXTLINE(modernize-avoid-c-arrays)
20 u32 pid[3];
21
22 CEXP facet() = default;
23 490 CEXP facet(u32 id, vec<point3d<FP>> CR vp, u32 id0, u32 id1, u32 id2) NE : id(id), p(vp[id0], vp[id1], vp[id2]), pid{id0, id1, id2} {}
24 454 CEXP facet(vec<point3d<FP>> CR vp, u32 id0, u32 id1, u32 id2) NE : facet(0, vp, id0, id1, id2) {}
25 36 CEXP void in(u32 n1, u32 n2, u32 n3)NE { n[0] = n1, n[1] = n2, n[2] = n3; }
26 };
27 u32 tm = 0, idx = 0;
28 kahan<FP> suf_area = 0;
29
30 public:
31 //! 1-indexed
32 vec<facet> faces;
33 //! 1-indexed, no co-linear or co-plane
34 9 cvh3d(vec<point3d<FP>> CR vp) NE {
35 9 faces.emplace_back(vp, 0, 0, 0), faces[0].isdel = 1;
36 9 vvecu ptsid(5);
37 {
38 #define v(n, i) vp[n##id[i]]
39 9 arr<u32, 6> pid{1, 1, 1, 1, 1, 1};
40
2/2
✓ Branch 1 taken 441 times.
✓ Branch 2 taken 9 times.
 450 flt_ (u32, i, 2, (u32)vp.size()) {
41
2/2
✓ Branch 4 taken 27 times.
✓ Branch 5 taken 414 times.
 441 if (is_gt(vp[i].x, v(p, 0).x)) pid[0] = i;
42
2/2
✓ Branch 4 taken 9 times.
✓ Branch 5 taken 432 times.
 441 if (is_gt(v(p, 1).x, vp[i].x)) pid[1] = i;
43
2/2
✓ Branch 4 taken 25 times.
✓ Branch 5 taken 416 times.
 441 if (is_gt(vp[i].y, v(p, 2).y)) pid[2] = i;
44
2/2
✓ Branch 4 taken 13 times.
✓ Branch 5 taken 428 times.
 441 if (is_gt(v(p, 3).y, vp[i].y)) pid[3] = i;
45
2/2
✓ Branch 4 taken 26 times.
✓ Branch 5 taken 415 times.
 441 if (is_gt(vp[i].z, v(p, 4).z)) pid[4] = i;
46
2/2
✓ Branch 4 taken 16 times.
✓ Branch 5 taken 425 times.
 441 if (is_gt(v(p, 5).z, vp[i].z)) pid[5] = i;
47 }
48 9 arr<u32, 4> sid{pid[0], pid[0], pid[0], pid[0]};
49
2/2
✓ Branch 0 taken 54 times.
✓ Branch 1 taken 9 times.
 63 flt_ (u32, i, 0, 6)
50
2/2
✓ Branch 0 taken 135 times.
✓ Branch 1 taken 54 times.
 189 flt_ (u32, j, i + 1, 6)
51
2/2
✓ Branch 11 taken 20 times.
✓ Branch 12 taken 115 times.
 135 if (is_gt(dist3_PP(v(p, i), v(p, j)), dist3_PP(v(s, 0), v(s, 1)))) sid[0] = pid[i], sid[1] = pid[j];
52
2/2
✓ Branch 0 taken 54 times.
✓ Branch 1 taken 9 times.
 63 flt_ (u32, i, 0, 6)
53
2/2
✓ Branch 17 taken 10 times.
✓ Branch 18 taken 44 times.
 54 if (is_gt(dist3_PL(v(p, i), {v(s, 0), v(s, 1)}), dist3_PL(v(s, 2), {v(s, 0), v(s, 1)}))) sid[2] = pid[i];
54
2/2
✓ Branch 1 taken 450 times.
✓ Branch 2 taken 9 times.
 459 flt_ (u32, i, 1, (u32)vp.size())
55
2/2
✓ Branch 14 taken 23 times.
✓ Branch 15 taken 427 times.
 450 if (is_gt(abs(sdist3_PlP(planev<FP>{v(s, 0), v(s, 1), v(s, 2)}, vp[i])),
56 473 abs(sdist3_PlP(planev<FP>{v(s, 0), v(s, 1), v(s, 2)}, v(s, 3))))) sid[3] = i;
57
2/2
✓ Branch 11 taken 3 times.
✓ Branch 12 taken 6 times.
 9 if (is_neg(sdist3_PlP(planev<FP>{v(s, 0), v(s, 1), v(s, 2)}, v(s, 3)))) swap(sid[1], sid[2]);
58 9 faces.emplace_back(1, vp, sid[0], sid[2], sid[1]), faces.emplace_back(2, vp, sid[0], sid[1], sid[3]);
59 9 faces.emplace_back(3, vp, sid[1], sid[2], sid[3]), faces.emplace_back(4, vp, sid[2], sid[0], sid[3]);
60 9 faces[1].in(4, 3, 2), faces[2].in(1, 3, 4), faces[3].in(1, 4, 2), faces[4].in(1, 2, 3);
61
2/2
✓ Branch 1 taken 450 times.
✓ Branch 2 taken 9 times.
 459 flt_ (u32, i, 1, (u32)vp.size()) {
62
10/10
✓ Branch 4 taken 441 times.
✓ Branch 5 taken 9 times.
✓ Branch 10 taken 432 times.
✓ Branch 11 taken 9 times.
✓ Branch 16 taken 423 times.
✓ Branch 17 taken 9 times.
✓ Branch 22 taken 9 times.
✓ Branch 23 taken 414 times.
✓ Branch 24 taken 36 times.
✓ Branch 25 taken 414 times.
 450 if (vp[i] == v(s, 0) || vp[i] == v(s, 1) || vp[i] == v(s, 2) || vp[i] == v(s, 3)) continue;
63
2/2
✓ Branch 0 taken 959 times.
✓ Branch 1 taken 56 times.
 1015 flt_ (u32, j, 1, 5)
64
2/2
✓ Branch 3 taken 358 times.
✓ Branch 4 taken 601 times.
 959 if (relation_PlP(faces[j].p, vp[i]) == RELPLP::above) {
65 358 ptsid[j].push_back(i);
66 358 break;
67 }
68 }
69 9 idx = 1;
70 #undef v
71 }
72 struct edge {
73 u32 netid, facetid;
74 };
75 36 vec<edge> e1(vp.size()), e2(vp.size());
76 63 vecu vistime(vp.size()), resfdel(vp.size()), resfnew(vp.size()), resptid(vp.size());
77 807 auto horizon = [&](auto&& f, u32 id, point3d<FP> CR p) NE -> u32 {
78
2/2
✓ Branch 2 taken 445 times.
✓ Branch 3 taken 353 times.
 798 if (relation_PlP(faces[id].p, p) != RELPLP::above) return 0;
79
2/2
✓ Branch 1 taken 124 times.
✓ Branch 2 taken 229 times.
 353 if (faces[id].vistime == tm) return -1_u32;
80 229 faces[id].vistime = tm, faces[id].isdel = 1, resfdel.push_back(faces[id].id);
81 229 u32 ret = -2_u32;
82
2/2
✓ Branch 0 taken 687 times.
✓ Branch 1 taken 229 times.
 916 flt_ (u32, i, 0, 3)
83
2/2
✓ Branch 2 taken 445 times.
✓ Branch 3 taken 242 times.
 687 if (u32 res = f(f, faces[id].n[i], p); !res) {
84 445 arr<u32, 2> pt{faces[id].pid[i], faces[id].pid[(i + 1) % 3]};
85
2/2
✓ Branch 0 taken 890 times.
✓ Branch 1 taken 445 times.
 1335 flt_ (u32, j, 0, 2)
86
2/2
✓ Branch 2 taken 445 times.
✓ Branch 3 taken 445 times.
 890 if (vistime[pt[j]] != tm) vistime[pt[j]] = tm, e1[pt[j]].netid = pt[j ^ 1], e1[pt[j]].facetid = faces[id].n[i];
87 445 else e2[pt[j]].netid = pt[j ^ 1], e2[pt[j]].facetid = faces[id].n[i];
88 445 ret = pt[0];
89
2/2
✓ Branch 0 taken 118 times.
✓ Branch 1 taken 124 times.
 242 } else if ((int)res > 0) ret = res;
90 229 return ret;
91 };
92 18 std::queue<u32> que({idx});
93
2/2
✓ Branch 2 taken 27 times.
✓ Branch 3 taken 9 times.
 36 for (auto i : faces[idx].n) que.push(i);
94 9 u32 snew = 0;
95
2/2
✓ Branch 1 taken 481 times.
✓ Branch 2 taken 9 times.
 490 while (!que.empty()) {
96 481 u32 nf = que.front();
97 481 que.pop();
98
2/2
✓ Branch 1 taken 104 times.
✓ Branch 2 taken 377 times.
 747 if (faces[nf].isdel) continue;
99
2/2
✓ Branch 2 taken 266 times.
✓ Branch 3 taken 111 times.
 377 if (ptsid[nf].empty()) {
100 266 snew = nf;
101 266 continue;
102 }
103 111 u32 pid = ptsid[nf][0];
104 111 point3d<FP> p = vp[pid];
105
2/2
✓ Branch 2 taken 548 times.
✓ Branch 3 taken 111 times.
 659 flt_ (u32, i, 1, (u32)ptsid[nf].size())
106
2/2
✓ Branch 7 taken 98 times.
✓ Branch 8 taken 450 times.
 548 if (is_gt(sdist3_PlP(faces[nf].p, vp[ptsid[nf][i]]),
107 646 sdist3_PlP(faces[nf].p, p))) p = vp[pid = ptsid[nf][i]];
108 111 ++tm, resfdel.clear();
109 111 u32 s = horizon(horizon, nf, p);
110 111 ++tm, resfnew.clear();
111 111 u32 from = 0, lastf = 0, fstf = 0;
112
2/2
✓ Branch 1 taken 445 times.
✓ Branch 2 taken 111 times.
 556 while (vistime[s] != tm) {
113 445 vistime[s] = tm;
114 u32 net, f, fnew;
115
2/2
✓ Branch 1 taken 148 times.
✓ Branch 2 taken 297 times.
 445 if (e1[s].netid == from) net = e2[s].netid, f = e2[s].facetid;
116 297 else net = e1[s].netid, f = e1[s].facetid;
117 445 u32 pt1 = -1_u32, pt2 = -1_u32, pt1id = -1_u32, pt2id = -1_u32;
118
2/2
✓ Branch 0 taken 1335 times.
✓ Branch 1 taken 445 times.
 1780 flt_ (u32, i, 0, 3) {
119
2/2
✓ Branch 4 taken 445 times.
✓ Branch 5 taken 890 times.
 1335 if (vp[s] == faces[f].p.get(i)) pt1 = i, pt1id = s;
120
2/2
✓ Branch 4 taken 445 times.
✓ Branch 5 taken 890 times.
 1335 if (vp[net] == faces[f].p.get(i)) pt2 = i, pt2id = net;
121 }
122
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 445 times.
 445 if ((pt1 + 1) % 3 != pt2) swap(pt1, pt2);
123 445 faces.emplace_back(vp, pt2id, pt1id, pid), faces[fnew = u32(faces.size() - 1)].id = u32(faces.size() - 1);
124
2/2
✓ Branch 4 taken 334 times.
✓ Branch 5 taken 111 times.
 445 if (ptsid.emplace_back(), resfnew.push_back(fnew), faces[faces[fnew].n[0] = f].n[pt1] = fnew; lastf) {
125
1/2
✓ Branch 3 taken 334 times.
✗ Branch 4 not taken.
 334 if (*faces[fnew].p.v == *faces[lastf].p.u) faces[fnew].n[1] = lastf, faces[lastf].n[2] = fnew;
126 ✗ else faces[fnew].n[2] = lastf, faces[lastf].n[1] = fnew;
127 111 } else fstf = fnew;
128 445 lastf = fnew, from = s, s = net;
129 }
130
1/2
✓ Branch 3 taken 111 times.
✗ Branch 4 not taken.
 111 if (*faces[fstf].p.v == *faces[lastf].p.u) faces[fstf].n[1] = lastf, faces[lastf].n[2] = fstf;
131 ✗ else faces[fstf].n[2] = lastf, faces[lastf].n[1] = fstf;
132 111 resptid.clear();
133
2/2
✓ Branch 10 taken 229 times.
✓ Branch 11 taken 111 times.
 340 for (auto i : resfdel) move(ptsid[i], std::back_inserter(resptid)), ptsid[i].clear();
134
2/2
✓ Branch 5 taken 838 times.
✓ Branch 6 taken 111 times.
 949 for (auto i : resptid) {
135
2/2
✓ Branch 2 taken 111 times.
✓ Branch 3 taken 727 times.
 838 if (vp[i] == p) continue;
136
2/2
✓ Branch 5 taken 2072 times.
✓ Branch 6 taken 247 times.
 2319 for (auto j : resfnew)
137
2/2
✓ Branch 3 taken 480 times.
✓ Branch 4 taken 1592 times.
 2072 if (relation_PlP(faces[j].p, vp[i]) == RELPLP::above) {
138 480 ptsid[j].push_back(i);
139 480 break;
140 }
141 }
142
2/2
✓ Branch 6 taken 445 times.
✓ Branch 7 taken 111 times.
 556 for (auto i : resfnew) que.push(i);
143 }
144 9 idx = snew;
145 9 }
146 template <class F>
147 requires requires(F op, kahan<FP>& v, planev<FP> p) { op(v, p); }
148 9 CEXP void dfs(F&& op) NE {
149 774 auto f = [&](auto&& f, u32 nf) NE -> void {
150
2/2
✓ Branch 1 taken 513 times.
✓ Branch 2 taken 252 times.
 765 if (faces[nf].vistime == tm) return;
151
2/2
✓ Branch 5 taken 756 times.
✓ Branch 6 taken 252 times.
 1008 for (faces[nf].vistime = tm, op(suf_area, faces[nf].p); auto i : faces[nf].n) f(f, i);
152 };
153 9 ++tm, f(f, idx);
154 9 }
155 9 CEXP FP surface_area() NE {
156
1/2
✗ Branch 2 not taken.
✓ Branch 3 taken 9 times.
 9 if (!is_zero((FP)suf_area)) return suf_area;
157 261 dfs([](kahan<FP>& v, planev<FP> p) NE { v += p.area(); });
158 9 return suf_area;
159 }
160 };
161
162 } // namespace tifa_libs
163