Branch data Line data Source code
1 : : /*
2 : : * primitives.cpp
3 : : *
4 : : * Created on: Mar 22, 2010
5 : : * Author: iulian
6 : : */
7 : :
8 : : #include "primitives.h"
9 : :
10 : 65628 : void filterValid(moab::Interface * mb, std::vector<moab::EntityHandle> & io) {
11 : 65628 : int next = 0, size = io.size();
12 [ - + ]: 65628 : if (size==0)
13 : 0 : return;
14 [ + - ]: 65628 : std::vector<unsigned char> tags(size);
15 [ + - ][ + - ]: 65628 : moab::ErrorCode rval = mb->tag_get_data(validTag, &io[0], size, &tags[0] );
[ + - ]
16 [ - + ]: 65628 : if (moab::MB_SUCCESS != rval)
17 : 0 : return;
18 [ + + ]: 800422 : for (int i = 0; i < size; i++) {
19 : : //moab::EntityHandle eh = io[i];
20 [ + - ][ + + ]: 734794 : if (tags[i]) {
21 [ + - ][ + - ]: 629382 : io[next++] = io[i];
22 : : }
23 : : }
24 [ + + ]: 65628 : if (next<size)
25 [ + - ]: 28488 : io.resize(next);
26 : 65628 : return;
27 : : }
28 : :
29 : 61160 : moab::ErrorCode contractionRegion(moab::Interface * mb, moab::EntityHandle v1,
30 : : moab::EntityHandle v2, std::vector<moab::EntityHandle>& changed) {
31 : : moab::EntityHandle vlist[2];
32 : 61160 : vlist[0] = v1;
33 : 61160 : vlist[1] = v2;
34 : : // it makes more sense to use vector, than range
35 : : // the entities could be very disjoint at some point
36 : : // moab::Range adj_ents;
37 : : moab::ErrorCode rval = mb->get_adjacencies(vlist, 2, 2, false, changed,
38 [ + - ]: 61160 : moab::Interface::UNION);
39 [ + - ]: 61160 : if (opts.useDelayedDeletion)
40 [ + - ]: 61160 : filterValid(mb, changed);
41 : 61160 : return rval;
42 : : }
43 : :
44 : : //
45 : 0 : int classifyVertex(moab::Interface * mb, moab::EntityHandle v1) {
46 : : // return 0, 1, or 2 if the vertex is interior, on the border, or
47 : : // border only (corner)
48 : : // to do
49 [ # # ]: 0 : std::vector<moab::EntityHandle> adjEdges;
50 : : moab::ErrorCode rval = mb->get_adjacencies(&v1, 1, 1, false, adjEdges,
51 [ # # ]: 0 : moab::Interface::UNION);
52 [ # # ]: 0 : if (moab::MB_SUCCESS != rval)
53 : 0 : return 0; // interior??
54 [ # # ]: 0 : if (opts.useDelayedDeletion)
55 [ # # ]: 0 : filterValid(mb, adjEdges);
56 : 0 : unsigned int nBorder = 0;
57 [ # # ]: 0 : for (unsigned int i = 0; i < adjEdges.size(); i++) {
58 [ # # ]: 0 : moab::EntityHandle edg = adjEdges[i];
59 [ # # ]: 0 : std::vector<moab::EntityHandle> adjFaces;
60 : : rval = mb->get_adjacencies(&edg, 1, 2, false, adjFaces,
61 [ # # ]: 0 : moab::Interface::UNION);
62 [ # # ]: 0 : if (opts.useDelayedDeletion)
63 [ # # ]: 0 : filterValid(mb, adjFaces);
64 [ # # ]: 0 : if (adjFaces.size() == 1)
65 : 0 : nBorder++;
66 : 0 : }
67 [ # # ]: 0 : if (nBorder == 0)
68 : 0 : return 0;// everything interior
69 [ # # ]: 0 : else if (nBorder == adjEdges.size())
70 : 0 : return 2;
71 : : else
72 : 0 : return 1; // some edges are interior
73 : : }
74 : :
75 : 1154548 : Vec3 getVec3FromMBVertex(moab::Interface * mbi, moab::EntityHandle v) {
76 : : double c[3];
77 [ + - ]: 1154548 : mbi->get_coords(&v, 1, c);
78 [ + - ]: 1154548 : return Vec3(c[0], c[1], c[2]);
79 : : }
80 : : // every time we are getting the normal, we compute a new plane
81 : : // maybe we should store it??
82 : : // No debate needed
83 : : // it will be much cheaper to store it, for "-m" option
84 : : // there, we will need it a lot
85 : 10400 : Plane trianglePlane(moab::Interface * mb, moab::EntityHandle tri) {
86 : : // retrieve it from tag
87 : 10400 : Plane pl;
88 : 10400 : mb->tag_get_data(planeDataTag, &tri, 1, &pl);
89 : 10400 : return pl;
90 : : }
91 : :
92 : 10000 : void computeTrianglePlane (moab::Interface * mb, moab::EntityHandle tri)
93 : : {
94 : : // get connectivity of triangle
95 : : const moab::EntityHandle * conn;
96 : : int num_nodes;
97 [ + - ]: 10000 : moab::ErrorCode rval = mb->get_connectivity(tri, conn, num_nodes);
98 [ + - ][ - + ]: 10000 : assert(3==num_nodes && rval == moab::MB_SUCCESS);
[ # # ]
99 [ + - ]: 10000 : Vec3 ve1 = getVec3FromMBVertex(mb, conn[0]);
100 [ + - ]: 10000 : Vec3 ve2 = getVec3FromMBVertex(mb, conn[1]);
101 [ + - ]: 10000 : Vec3 ve3 = getVec3FromMBVertex(mb, conn[2]);
102 [ + - ]: 10000 : Plane pl = Plane(ve1, ve2, ve3);
103 [ + - ]: 10000 : mb->tag_set_data(planeDataTag, &tri, 1, &pl);
104 : 10000 : return;
105 : : }
106 : :
107 : 544 : moab::ErrorCode contract(moab::Interface * mb, moab::EntityHandle v0, moab::EntityHandle v1,
108 : : Vec3 & vnew, std::vector<moab::EntityHandle> & changed) {
109 : :
110 : : //
111 : : //// Collect all the faces that are going to be changed
112 : : ////
113 [ + - ]: 544 : std::vector<moab::EntityHandle> adj_entities;
114 [ + - ]: 544 : contractionRegion(mb, v0, v1, adj_entities);
115 : : // those are all triangles that are affected
116 : : // find also all edges that are affect
117 : : moab::EntityHandle vlist[2];
118 : 544 : vlist[0] = v0;
119 : 544 : vlist[1] = v1;
120 : : // it makes more sense to use vector, than range
121 : : // the entities could be very disjoint at some point
122 : : // moab::Range adj_ents;
123 [ + - ]: 1088 : std::vector<moab::EntityHandle> edges;
124 : : moab::ErrorCode rval = mb->get_adjacencies(vlist, 2, 1, false, edges,
125 [ + - ]: 544 : moab::Interface::UNION);
126 [ - + ]: 544 : if (moab::MB_SUCCESS != rval) return rval;
127 [ + - ]: 544 : if (opts.useDelayedDeletion)
128 [ + - ]: 544 : filterValid(mb, edges);
129 [ - + ][ # # ]: 544 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
130 [ # # ][ # # ]: 0 : *opts.logfile << "Edges Adjacent:" << edges.size();
131 [ # # ]: 0 : for (unsigned int i = 0; i < edges.size(); i++)
132 [ # # ][ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(edges[i]);
[ # # ][ # # ]
133 [ # # ]: 0 : *opts.logfile << std::endl;
134 : : }
135 : : // we have the edges and the triangles that are affected
136 : : // 2 situations
137 : : // 1) edge v0 v1 is existing
138 : : // we will delete edge (v0, v1), and triangles formed
139 : : // with edge (v0, v1)
140 : : // 2) edge v0 v1 is not existing, but due to proximity
141 : : // only edges v2 v1 , v2, v0 will need to merge
142 : : // more important is case 1)
143 : :
144 : : // first, find edge v0, v1
145 : : moab::EntityHandle ev0v1;
146 : 544 : int foundEdge = 0;
147 [ + - ]: 3596 : for (unsigned int i = 0; i < edges.size(); i++) {
148 [ + - ]: 3596 : moab::EntityHandle e = edges[i];
149 : : int nnodes;
150 : : const moab::EntityHandle * conn2;
151 [ + - ]: 3596 : mb->get_connectivity(e, conn2, nnodes);
152 [ - + ][ # # ]: 3596 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS)
153 [ # # ][ # # ]: 0 : *opts.logfile << "Edge: " << mb->id_from_handle(e) << " nnodes:"
[ # # ][ # # ]
154 [ # # ][ # # ]: 0 : << nnodes << " vertices:" << mb->id_from_handle(conn2[0])
[ # # ][ # # ]
155 [ # # ][ # # ]: 0 : << " " << mb->id_from_handle(conn2[1]) << std::endl;
[ # # ][ # # ]
156 [ + + ][ + + ]: 3596 : if ((conn2[0] == v0 && conn2[1] == v1) || (conn2[0] == v1 && conn2[1]
[ + + ][ + + ]
157 : 728 : == v0)) {
158 : 544 : foundEdge = 1;
159 : 544 : ev0v1 = e; // could be ev1v0, but who cares?
160 [ - + ][ # # ]: 544 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS)
161 [ # # ][ # # ]: 0 : *opts.logfile << "Edge found " << mb->id_from_handle(e)
[ # # ]
162 [ # # ]: 0 : << std::endl;
163 : 544 : break;
164 : : }
165 : :
166 : : }
167 : : // set the position of new vertices in vnew
168 : : double newCoords[3];
169 [ + - ]: 544 : newCoords[0] = vnew[0];
170 [ + - ]: 544 : newCoords[1] = vnew[1];
171 [ + - ]: 544 : newCoords[2] = vnew[2];
172 [ + - ]: 544 : mb->set_coords(&v0, 1, newCoords);
173 [ + - ]: 544 : mb->set_coords(&v1, 1, newCoords);
174 : : // although, vertex v1 will be deleted in the end; do we really need to set it?
175 : : // yes, for merging purposes
176 : : //
177 : :
178 [ + - ]: 544 : if (opts.useDelayedDeletion) {
179 : : // big copy from version 3512
180 : : // although vertex v1 will be merged!!
181 [ + - ]: 544 : std::vector<moab::EntityHandle> edgePairs; // the one that has v0 will
182 : : // be kept
183 [ + - ][ + - ]: 1088 : std::vector<moab::EntityHandle> edgesWithV1;
184 [ + - ]: 544 : if (foundEdge) {
185 : : // this is case 1, the most complicated
186 : : // get triangles connected to edge ev0v1
187 [ + - ]: 544 : std::vector<moab::EntityHandle> tris;
188 : : rval = mb->get_adjacencies(&ev0v1, 1, 2, false, tris,
189 [ + - ]: 544 : moab::Interface::UNION);
190 [ - + ]: 544 : if (moab::MB_SUCCESS != rval) return rval;
191 [ + - ]: 544 : filterValid(mb, tris);
192 [ - + ][ # # ]: 544 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
193 [ # # ]: 0 : *opts.logfile << "Triangles adjacent to found edge:"
194 [ # # ][ # # ]: 0 : << tris.size() << ":";
195 [ # # ]: 0 : for (unsigned int i = 0; i < tris.size(); i++)
196 [ # # ][ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(tris[i]);
[ # # ][ # # ]
197 [ # # ]: 0 : *opts.logfile << std::endl;
198 : : }
199 [ + + ]: 1546 : for (unsigned int i = 0; i < tris.size(); i++) {
200 [ + - ]: 1002 : moab::EntityHandle triangleThatCollapses = tris[i];
201 [ + - ]: 1002 : std::vector<moab::EntityHandle> localEdges;
202 : : rval = mb->get_adjacencies(&triangleThatCollapses, 1, 1, false,
203 [ + - ]: 1002 : localEdges, moab::Interface::UNION);
204 [ - + ]: 1002 : if (moab::MB_SUCCESS != rval) return rval;
205 [ + - ]: 1002 : filterValid(mb, localEdges);
206 [ - + ][ # # ]: 1002 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
207 [ # # ]: 0 : *opts.logfile << "Triangle " << mb->id_from_handle(
208 [ # # ][ # # ]: 0 : triangleThatCollapses) << " Edges: "
[ # # ]
209 [ # # ]: 0 : << localEdges.size();
210 : :
211 [ # # ]: 0 : for (unsigned int i = 0; i < localEdges.size(); i++)
212 [ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(
213 [ # # ][ # # ]: 0 : localEdges[i]);
[ # # ]
214 [ # # ]: 0 : *opts.logfile << std::endl;
215 : : }
216 : : // find the edges that contains v0
217 : : moab::EntityHandle e[2];// the 2 edges e0, e1, that are not ev0v1;
218 [ - + ]: 1002 : if (localEdges.size() != 3)
219 : 0 : return moab::MB_FAILURE; // failure
220 : 1002 : int index = 0;
221 [ + + ]: 4008 : for (int k = 0; k < 3; k++)
222 [ + - ][ + + ]: 3006 : if (localEdges[k] != ev0v1)
223 [ + - ]: 2004 : e[index++] = localEdges[k];
224 : : // among those 2 edges, find out which one has v0, and which one v1
225 [ - + ]: 1002 : if (index != 2)
226 : 0 : return moab::MB_FAILURE; // failure
227 [ + - ][ + - ]: 2554 : for (int j = 0; j < 2; j++) {
228 : : int nn;
229 : : const moab::EntityHandle * conn2;
230 [ + - ]: 1552 : mb->get_connectivity(e[j], conn2, nn);
231 [ + + ][ + + ]: 1552 : if (conn2[0] == v0 || conn2[1] == v0) {
232 : : // this is the edge that will be kept, the other one collapsed
233 [ + - ]: 1002 : edgePairs.push_back(e[j]);
234 : : //j = (j + 1) % 2;// the other one
235 : 1002 : int j1 = (j + 1) % 2; // do not modify j, as someone might do something crazy
236 : : // with that index in a for_loop (optimizations?)
237 [ + - ]: 1002 : edgePairs.push_back(e[j1]);
238 [ + - ]: 1002 : edgesWithV1.push_back(e[j1]);
239 : 1002 : break; // no need to check the other one. it
240 : : // will contain v1
241 : : }
242 : : }
243 : 1002 : }
244 : : // look at all triangles that are adjacent
245 : : // invalidate first tris
246 : 544 : unsigned char invalid = 0;
247 [ + + ]: 544 : if (tris.size() > 0)
248 : : {
249 : : // set the invalid tag one by one, we do not want to create an array of invalids
250 [ + + ]: 1544 : for (unsigned int k = 0; k < tris.size(); k++)
251 : : {
252 [ + - ]: 1002 : moab::EntityHandle ct = tris[k];
253 [ + - ]: 1002 : mb->tag_set_data(validTag, &ct, 1, &invalid);
254 : : }
255 : : }
256 [ - + ][ # # ]: 544 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
257 [ # # ][ # # ]: 0 : *opts.logfile << "Triangles invalidated: " << tris.size()
258 [ # # ]: 0 : << ":";
259 [ # # ]: 0 : for (unsigned int i = 0; i < tris.size(); i++)
260 [ # # ][ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(tris[i]);
[ # # ][ # # ]
261 [ # # ]: 0 : *opts.logfile << std::endl;
262 : : }
263 : : // then look at all triangles that are not in tris (adj_entities), and then
264 : : // replace connectivity of v1 with v0
265 [ + - ][ + - ]: 1088 : std::vector<moab::EntityHandle> trisChanged;
266 [ + + ]: 6034 : for (unsigned int i = 0; i < adj_entities.size(); i++) {
267 [ + - ]: 5490 : moab::EntityHandle tr = adj_entities[i];
268 [ + - ][ + + ]: 5490 : if (!ehIsValid(tr))
269 : 1002 : continue;
270 : : // see the connectivity of tr; if there is a v1, replace it with v0
271 : : // will that propagate to edges? or do we have to set it separately?
272 : : int nn;
273 : : const moab::EntityHandle * conn3;
274 [ + - ]: 4488 : mb->get_connectivity(tr, conn3, nn);
275 : :
276 [ - + ][ # # ]: 4488 : assert(3==nn);
277 [ + + ]: 17952 : for (int j = 0; j < 3; j++) {
278 [ + + ]: 13464 : if (conn3[j] == v1) {
279 : : // replace it with v0, and reset it
280 : : moab::EntityHandle connNew[3];
281 : 1834 : connNew[0] = conn3[0];
282 : 1834 : connNew[1] = conn3[1];
283 : 1834 : connNew[2] = conn3[2];
284 : 1834 : connNew[j] = v0;
285 [ - + ][ # # ]: 1834 : if (opts.logfile && opts.selected_output
286 : 0 : & OUTPUT_CONTRACTIONS) {
287 [ # # ]: 0 : std::vector<moab::EntityHandle> localEdges;
288 : : rval = mb->get_adjacencies(&tr, 1, 1, false,
289 [ # # ]: 0 : localEdges, moab::Interface::UNION);
290 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
291 [ # # ]: 0 : filterValid(mb, localEdges);
292 [ # # ]: 0 : *opts.logfile << "Triangle t"
293 [ # # ][ # # ]: 0 : << mb->id_from_handle(tr)
294 [ # # ][ # # ]: 0 : << " filtered : " << localEdges.size();
295 [ # # ]: 0 : for (unsigned int j = 0; j < localEdges.size(); j++)
296 [ # # ]: 0 : *opts.logfile << " e" << mb->id_from_handle(
297 [ # # ][ # # ]: 0 : localEdges[j]);
[ # # ]
298 [ # # ][ # # ]: 0 : *opts.logfile << std::endl;
299 : : }
300 [ - + ][ # # ]: 1834 : if (opts.logfile && opts.selected_output
301 : 0 : & OUTPUT_CONTRACTIONS) {
302 [ # # ]: 0 : *opts.logfile << "replace connectivity t"
303 [ # # ][ # # ]: 0 : << mb->id_from_handle(tr) << " v"
[ # # ]
304 [ # # ][ # # ]: 0 : << mb->id_from_handle(conn3[0]) << " v"
[ # # ]
305 [ # # ][ # # ]: 0 : << mb->id_from_handle(conn3[1]) << " v"
[ # # ]
306 [ # # ][ # # ]: 0 : << mb->id_from_handle(conn3[2])
307 [ # # ]: 0 : << " to: v";
308 : : }
309 [ + - ]: 1834 : rval = mb->set_connectivity(tr, connNew, 3);
310 [ - + ]: 1834 : if (moab::MB_SUCCESS != rval) return rval;
311 [ - + ][ # # ]: 1834 : if (opts.logfile && opts.selected_output
312 : 0 : & OUTPUT_CONTRACTIONS) {
313 [ # # ][ # # ]: 0 : *opts.logfile << mb->id_from_handle(connNew[0])
314 [ # # ][ # # ]: 0 : << " v" << mb->id_from_handle(connNew[1])
[ # # ]
315 [ # # ][ # # ]: 0 : << " v" << mb->id_from_handle(connNew[2])
[ # # ]
316 [ # # ]: 0 : << std::endl;
317 : : }
318 [ + - ]: 1834 : trisChanged.push_back(tr);
319 : : }
320 : : }
321 : :
322 : : }
323 : 544 : validFaceCount -= tris.size();
324 [ + - ]: 544 : rval = mb->tag_set_data(validTag, &ev0v1, 1, &invalid);
325 [ - + ]: 544 : if (moab::MB_SUCCESS != rval) return rval;
326 : : // invalidate the edges connected for sure to v1
327 [ + + ]: 544 : if (edgesWithV1.size() > 0)
328 : : {
329 [ + + ]: 1544 : for(unsigned int k=0; k<edgesWithV1.size(); k++)
330 : : {
331 [ + - ]: 1002 : moab::EntityHandle eV1=edgesWithV1[k];
332 [ + - ]: 1002 : mb->tag_set_data(validTag, &eV1, 1, &invalid);
333 : : }
334 : : }
335 : : // reset the connectivity of some edges (from v1 to v0)
336 [ + + ]: 6900 : for (unsigned int i = 0; i < edges.size(); i++) {
337 [ + - ]: 6356 : moab::EntityHandle e1 = edges[i];
338 [ + - ][ + + ]: 6356 : if (!ehIsValid(e1)) // it could be the ones invalidated
339 : 1546 : continue;
340 : :
341 : : int nn;
342 : : const moab::EntityHandle * conn;
343 [ + - ]: 4810 : mb->get_connectivity(e1, conn, nn);
344 : :
345 [ - + ][ # # ]: 4810 : assert(2==nn);
346 [ + + ]: 14430 : for (int j = 0; j < 2; j++) {
347 [ + + ]: 9620 : if (conn[j] == v1) {
348 : : // replace it with v0, and reset it
349 : : moab::EntityHandle connNew[2];
350 : 1402 : connNew[0] = conn[0];
351 : 1402 : connNew[1] = conn[1];
352 : 1402 : connNew[j] = v0;
353 [ - + ][ # # ]: 1402 : if (opts.logfile && opts.selected_output
354 : 0 : & OUTPUT_CONTRACTIONS) {
355 [ # # ]: 0 : *opts.logfile << "replace connectivity edge: "
356 [ # # ][ # # ]: 0 : << mb->id_from_handle(e1) << " "
[ # # ]
357 [ # # ][ # # ]: 0 : << mb->id_from_handle(conn[0]) << " "
[ # # ]
358 [ # # ][ # # ]: 0 : << mb->id_from_handle(conn[1]) << " to: ";
[ # # ]
359 : : }
360 [ + - ]: 1402 : rval = mb->set_connectivity(e1, connNew, 2);
361 [ - + ]: 1402 : if (moab::MB_SUCCESS != rval) return rval;
362 [ - + ][ # # ]: 1402 : if (opts.logfile && opts.selected_output
363 : 0 : & OUTPUT_CONTRACTIONS) {
364 [ # # ][ # # ]: 0 : *opts.logfile << mb->id_from_handle(connNew[0])
365 [ # # ][ # # ]: 0 : << " " << mb->id_from_handle(connNew[1])
[ # # ]
366 [ # # ]: 1402 : << std::endl;
367 : : }
368 : : }
369 : : }
370 : : }
371 : : // the question: is the adjacency between triangles and edges restored?
372 : : // yes, it is; check set_connectivity logic
373 : : // we need to remove adjacencies between triangles and edges that are not valid
374 : : // and add some more adjacencies to the edges that are now part of the triangle
375 [ + + ]: 2378 : for (unsigned int i = 0; i < trisChanged.size(); i++) {
376 : : // get adjacencies now, and bail out early; maybe we should create if missing
377 [ + - ]: 1834 : std::vector<moab::EntityHandle> localEdges;
378 [ + - ]: 1834 : moab::EntityHandle tr = trisChanged[i];
379 : : rval = mb->get_adjacencies(&tr, 1, 1, false, localEdges,
380 [ + - ]: 1834 : moab::Interface::UNION);
381 [ - + ]: 1834 : if (moab::MB_SUCCESS != rval) return rval;
382 [ + - ]: 1834 : filterValid(mb, localEdges);
383 [ - + ][ # # ]: 1834 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
384 [ # # ][ # # ]: 0 : *opts.logfile << "Triangle t" << mb->id_from_handle(tr)
[ # # ]
385 [ # # ][ # # ]: 0 : << " filtered : " << localEdges.size();
386 [ # # ]: 0 : for (unsigned int j = 0; j < localEdges.size(); j++)
387 [ # # ]: 0 : *opts.logfile << " e" << mb->id_from_handle(
388 [ # # ][ # # ]: 0 : localEdges[j]);
[ # # ]
389 [ # # ]: 0 : *opts.logfile << std::endl;
390 : : }
391 [ - + ][ # # ]: 1834 : assert(localEdges.size()==3);
[ + - ]
392 : 1834 : }
393 : :
394 [ + - ]: 544 : filterValid(mb, adj_entities);
395 [ + - ]: 544 : changed = adj_entities; // deep copy
396 [ - + ][ # # ]: 544 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
397 [ # # ][ # # ]: 0 : *opts.logfile << "Triangles changed:" << changed.size();
398 [ # # ]: 0 : for (unsigned int i = 0; i < changed.size(); i++)
399 [ # # ][ # # ]: 0 : *opts.logfile << " t" << mb->id_from_handle(changed[i]);
[ # # ][ # # ]
400 [ # # ][ + - ]: 544 : *opts.logfile << std::endl;
401 : 544 : }
402 : : } else // this should appear only when proximity limit > 0
403 : : {
404 : : // in this case, we only need to worry if vertex 0 and 1 are connected to the same
405 : : // vertex 2; then we need to merge edges v0-v2 and v1 - v2
406 : : // no triangles get deleted, only those edges;
407 : : // the crack v0-v2-v1 gets seamed
408 : : // get edges connected to vertex v0
409 [ # # ]: 0 : std::vector<moab::EntityHandle> edges0;
410 : : rval = mb->get_adjacencies(&v0, 1, 1, false, edges0,
411 [ # # ]: 0 : moab::Interface::UNION);
412 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
413 [ # # ]: 0 : filterValid(mb, edges0);
414 : :
415 : : // get edges connected to vertex v1
416 [ # # ][ # # ]: 0 : std::vector<moab::EntityHandle> edges1;
417 : : rval = mb->get_adjacencies(&v1, 1, 1, false, edges1,
418 [ # # ]: 0 : moab::Interface::UNION);
419 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
420 [ # # ]: 0 : filterValid(mb, edges1);
421 : : // find all edges that will be merged, of type v0-v2 v1-v2 (so that they have a
422 : : // common vertex v2
423 : : // in that case, we will have to merge them as before, and delete the
424 : : // one that contains v1 before
425 : : // for sure, there is no edge from v0 to v1 !!!
426 : : // keep all the vertices from edges1, different from v1
427 [ # # ][ # # ]: 0 : std::vector<moab::EntityHandle> otherVertices1;
428 [ # # ]: 0 : for (unsigned int i = 0; i < edges1.size(); i++) {
429 [ # # ]: 0 : moab::EntityHandle edgeThatGoes = edges1[1];
430 : : // find the other vertex, not v1
431 : : int nn;
432 : : const moab::EntityHandle * conn2;
433 [ # # ]: 0 : mb->get_connectivity(edgeThatGoes, conn2, nn);
434 : 0 : int other = 0;
435 [ # # ]: 0 : if (conn2[0] == v1)
436 : 0 : other = 1;
437 [ # # ]: 0 : else if (conn2[1] == v1)// is this really necessary in a correct model?
438 : 0 : other = 0;
439 [ # # ]: 0 : otherVertices1.push_back(conn2[other]);
440 : : }
441 [ # # ]: 0 : for (unsigned int i = 0; i < edges0.size(); i++) {
442 [ # # ]: 0 : moab::EntityHandle edgeThatStays = edges0[i];
443 : : // find the other vertex, not v0
444 : : int nn;
445 : : const moab::EntityHandle * conn2;
446 [ # # ]: 0 : mb->get_connectivity(edgeThatStays, conn2, nn);
447 : 0 : int other = 0;
448 [ # # ]: 0 : if (conn2[0] == v1)
449 : 0 : other = 1;
450 [ # # ]: 0 : else if (conn2[1] == v1)// is this really necessary in a correct model?
451 : 0 : other = 0;
452 : 0 : moab::EntityHandle v2 = conn2[other];
453 : : // let's see now if v2 is among vertices from otherVertices1
454 [ # # ]: 0 : for (unsigned int i1 = 0; i1 < otherVertices1.size(); i1++) {
455 [ # # ][ # # ]: 0 : if (v2 == otherVertices1[i1]) {
456 : : // we have a match, some work to do
457 : : // invalidate the edge edges1[i1]
458 : 0 : unsigned char invalid = 0;
459 [ # # ][ # # ]: 0 : mb->tag_set_data(validTag, &(edges1[i1]), 1, &invalid);
460 : 0 : break; // we stop looking for a match, only one possible
461 : : }
462 : : }
463 : : }
464 : : // triangles that need reconnected
465 [ # # ][ # # ]: 0 : std::vector<moab::EntityHandle> tri1;
466 : : rval = mb->get_adjacencies(&v1, 1, 2, false, tri1,
467 [ # # ]: 0 : moab::Interface::UNION);
468 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
469 : : // start copy tri reconnect
470 : 0 : unsigned int i = 0;
471 [ # # ]: 0 : for (i=0; i < tri1.size(); i++) {
472 [ # # ]: 0 : moab::EntityHandle tr = tri1[i];
473 [ # # ][ # # ]: 0 : if (!ehIsValid(tr))
474 : 0 : continue;
475 : : // see the connectivity of tr; if there is a v1, replace it with v0
476 : : // will that propagate to edges? or do we have to set it separately?
477 : : int nn;
478 : : const moab::EntityHandle * conn3;
479 [ # # ]: 0 : mb->get_connectivity(tr, conn3, nn);
480 : :
481 [ # # ][ # # ]: 0 : assert(3==nn);
482 [ # # ]: 0 : for (int j = 0; j < 3; j++) {
483 [ # # ]: 0 : if (conn3[j] == v1) {
484 : : // replace it with v0, and reset it
485 : : moab::EntityHandle connNew[3];
486 : 0 : connNew[0] = conn3[0];
487 : 0 : connNew[1] = conn3[1];
488 : 0 : connNew[2] = conn3[2];
489 : 0 : connNew[j] = v0;
490 [ # # ]: 0 : rval = mb->set_connectivity(tr, connNew, 3);
491 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
492 [ # # ][ # # ]: 0 : if (opts.logfile && opts.selected_output
493 : 0 : & OUTPUT_CONTRACTIONS) {
494 [ # # ][ # # ]: 0 : *opts.logfile << mb->id_from_handle(connNew[0])
495 [ # # ][ # # ]: 0 : << " v" << mb->id_from_handle(connNew[1])
[ # # ]
496 [ # # ][ # # ]: 0 : << " v" << mb->id_from_handle(connNew[2])
[ # # ]
497 [ # # ]: 0 : << std::endl;
498 : : }
499 : : }
500 : : }
501 : :
502 : : }
503 : : // now reconnect edges1 that are still valid
504 [ # # ]: 0 : for (i = 0; i < edges1.size(); i++) {
505 [ # # ]: 0 : moab::EntityHandle e1 = edges1[i];
506 [ # # ][ # # ]: 0 : if (!ehIsValid(e1))
507 : 0 : continue;
508 : : // see the connectivity of e1; if there is a v1, replace it with v0
509 : : // will that propagate to edges? or do we have to set it separately?
510 : : int nn;
511 : : const moab::EntityHandle * conn3;
512 [ # # ]: 0 : mb->get_connectivity(e1, conn3, nn);
513 : :
514 [ # # ][ # # ]: 0 : assert(2==nn);
515 [ # # ]: 0 : for (int j = 0; j < 2; j++) {
516 [ # # ]: 0 : if (conn3[j] == v1) {
517 : : // replace it with v0, and reset it
518 : : moab::EntityHandle connNew[2];
519 : 0 : connNew[0] = conn3[0];
520 : 0 : connNew[1] = conn3[1];
521 : 0 : connNew[j] = v0;
522 [ # # ]: 0 : rval = mb->set_connectivity(e1, connNew, 2);
523 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
524 : :
525 : : }
526 : : }
527 : :
528 : : }
529 : : // all the triangles connected to v0 are now changed, need recomputation
530 : : //
531 : : rval = mb->get_adjacencies(&v0, 1, 2, false, changed,
532 [ # # ]: 0 : moab::Interface::UNION);
533 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
534 [ # # ][ # # ]: 544 : filterValid(mb, changed);
[ + - ]
535 : :
536 : 544 : }
537 : : // end big copy from version 3512
538 : : } else {
539 : : // vertex v1 will be merged!!
540 [ # # ]: 0 : std::vector<moab::EntityHandle> edgePairsToMerge; // the one that has v0 will
541 : : // be kept
542 [ # # ]: 0 : if (foundEdge) {
543 : : // this is case 1, the most complicated
544 : : // get triangles connected to edge ev0v1
545 [ # # ]: 0 : std::vector<moab::EntityHandle> tris;
546 : : rval = mb->get_adjacencies(&ev0v1, 1, 2, false, tris,
547 [ # # ]: 0 : moab::Interface::UNION);
548 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
549 : : // find all edges that will be merged ( xv0, xv1, etc)
550 [ # # ][ # # ]: 0 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
551 [ # # ]: 0 : *opts.logfile << "Triangles adjacent to found edge:"
552 [ # # ]: 0 : << tris.size();
553 [ # # ]: 0 : for (unsigned int i = 0; i < tris.size(); i++)
554 [ # # ][ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(tris[i]);
[ # # ][ # # ]
555 [ # # ]: 0 : *opts.logfile << std::endl;
556 : : }
557 : 0 : unsigned int i=0;
558 [ # # ]: 0 : for (i=0; i < tris.size(); i++) {
559 [ # # ]: 0 : moab::EntityHandle triangleThatCollapses = tris[i];
560 [ # # ]: 0 : std::vector<moab::EntityHandle> localEdges;
561 : : rval = mb->get_adjacencies(&triangleThatCollapses, 1, 1, false,
562 [ # # ]: 0 : localEdges, moab::Interface::UNION);
563 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
564 [ # # ][ # # ]: 0 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
565 [ # # ]: 0 : *opts.logfile << "Triangle " << mb->id_from_handle(
566 [ # # ][ # # ]: 0 : triangleThatCollapses) << " Edges: "
[ # # ]
567 [ # # ]: 0 : << localEdges.size();
568 : :
569 [ # # ]: 0 : for (unsigned int i = 0; i < localEdges.size(); i++)
570 [ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(
571 [ # # ][ # # ]: 0 : localEdges[i]);
[ # # ]
572 [ # # ]: 0 : *opts.logfile << std::endl;
573 : : }
574 : : // find the edges that contains v0
575 : : moab::EntityHandle e[2];// the 2 edges e0, e1, that are not ev0v1;
576 [ # # ]: 0 : if (localEdges.size() != 3)
577 : 0 : return moab::MB_FAILURE; // failure
578 : 0 : int index = 0;
579 [ # # ]: 0 : for (int k = 0; k < 3; k++)
580 [ # # ][ # # ]: 0 : if (localEdges[k] != ev0v1)
581 [ # # ]: 0 : e[index++] = localEdges[k];
582 : : // among those 2 edges, find out which one has v0, and which one v1
583 [ # # ]: 0 : if (index != 2)
584 : 0 : return moab::MB_FAILURE; // failure
585 [ # # ][ # # ]: 0 : for (int j = 0; j < 2; j++) {
586 : : int nn;
587 : : const moab::EntityHandle * conn2;
588 [ # # ]: 0 : mb->get_connectivity(e[j], conn2, nn);
589 [ # # ][ # # ]: 0 : if (conn2[0] == v0 || conn2[1] == v0) {
590 : : // this is the edge that will be kept, the other one collapsed
591 [ # # ]: 0 : edgePairsToMerge.push_back(e[j]);
592 : : //j = (j + 1) % 2;
593 : 0 : int j1 = (j + 1)%2;// do not modify original j
594 [ # # ]: 0 : edgePairsToMerge.push_back(e[j1]);
595 : 0 : break; // no need to check the other one. it
596 : : // will contain v1
597 : : }
598 : : }
599 : 0 : }
600 : : // first merge vertices v0 and v1 : will also NOT delete v1 (yet)
601 : : // the tag on v1 will be deleted too, and we do not want that, at least until
602 : : // after the merging of edges, and deleting the pair
603 [ # # ]: 0 : rval = mb->merge_entities(v0, v1, false, false);
604 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
605 : : // merge edgePairsToMerge // now, v0 and v1 should be collapsed!
606 [ # # ]: 0 : for (unsigned int j = 0; j < edgePairsToMerge.size(); j += 2) {
607 : : // will also delete edges that contained v1 before
608 [ # # ][ # # ]: 0 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS)
609 [ # # ]: 0 : *opts.logfile << "Edges merged:" << mb->id_from_handle(
610 [ # # ][ # # ]: 0 : edgePairsToMerge[j]) << " " << mb->id_from_handle(
[ # # ][ # # ]
611 [ # # ][ # # ]: 0 : edgePairsToMerge[j + 1]) << std::endl;
[ # # ][ # # ]
612 [ # # ]: 0 : mb->merge_entities(edgePairsToMerge[j],
613 [ # # ][ # # ]: 0 : edgePairsToMerge[j + 1], false, true);
614 : :
615 : : }
616 : : // the only things that need deleted are triangles
617 [ # # ][ # # ]: 0 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
618 [ # # ][ # # ]: 0 : *opts.logfile << "Triangles invalidated:" << tris.size();
619 [ # # ]: 0 : for (unsigned int i = 0; i < tris.size(); i++)
620 [ # # ][ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(tris[i]);
[ # # ][ # # ]
621 [ # # ]: 0 : *opts.logfile << std::endl;
622 : : }
623 [ # # ][ # # ]: 0 : rval = mb->delete_entities(&(tris[0]), tris.size());
624 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
625 [ # # ]: 0 : if (iniSet)
626 [ # # ][ # # ]: 0 : mb->remove_entities(iniSet, &(tris[0]), tris.size());
627 : 0 : validFaceCount -= tris.size();
628 : : // hopefully, all adjacencies are preserved
629 : : // delete now the edge ev0v1
630 [ # # ][ # # ]: 0 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS)
631 [ # # ]: 0 : *opts.logfile << "Edge invalidated"
632 [ # # ][ # # ]: 0 : << mb->id_from_handle(ev0v1) << std::endl;
[ # # ]
633 [ # # ]: 0 : rval = mb->delete_entities(&ev0v1, 1);
634 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
635 [ # # ]: 0 : if (iniSet)
636 [ # # ]: 0 : mb->remove_entities(iniSet, &ev0v1, 1);// it may not be part of it, but
637 : : // delete it anyway
638 : : // among adj_entities, remove the tris triangles, and return them
639 [ # # ]: 0 : for (unsigned int j = 0; j < adj_entities.size(); j++) {
640 [ # # ]: 0 : moab::EntityHandle F = adj_entities[j];
641 : 0 : int inTris = 0;
642 [ # # ]: 0 : for (unsigned int k = 0; k < tris.size(); k++)
643 [ # # ][ # # ]: 0 : if (F == tris[k]) {
644 : 0 : inTris = 1;
645 : 0 : break;
646 : : }
647 [ # # ]: 0 : if (!inTris)
648 [ # # ]: 0 : changed.push_back(F);
649 : : }
650 [ # # ][ # # ]: 0 : if (opts.logfile && opts.selected_output & OUTPUT_CONTRACTIONS) {
651 [ # # ][ # # ]: 0 : *opts.logfile << "Triangles changed:" << changed.size();
652 [ # # ]: 0 : for (unsigned int i = 0; i < changed.size(); i++)
653 [ # # ][ # # ]: 0 : *opts.logfile << " " << mb->id_from_handle(changed[i]);
[ # # ][ # # ]
654 [ # # ][ # # ]: 0 : *opts.logfile << std::endl;
655 : 0 : }
656 : :
657 : : } else // this should appear only when proximity limit > 0
658 : : {
659 : : // in this case, we only need to worry if vertex 0 and 1 are connected to the same
660 : : // vertex 2; then we need to merge edges v0-v2 and v1 - v2
661 : : // no triangles get deleted, only those edges
662 : : // get edges connected to vertex v0
663 [ # # ]: 0 : std::vector<moab::EntityHandle> edges0;
664 : : rval = mb->get_adjacencies(&v0, 1, 1, false, edges0,
665 [ # # ]: 0 : moab::Interface::UNION);
666 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
667 : :
668 : : // get edges connected to vertex v1
669 [ # # ][ # # ]: 0 : std::vector<moab::EntityHandle> edges1;
670 : : rval = mb->get_adjacencies(&v1, 1, 1, false, edges1,
671 [ # # ]: 0 : moab::Interface::UNION);
672 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
673 : : // find all edges that will be merged, of type v0-v2 v1-v2 (so that they have a
674 : : // common vertex v2
675 : : // in that case, we will have to merge them as before, and delete the
676 : : // one that contains v1 before
677 : : // for sure, there is no edge from v0 to v1 !!!
678 : : // keep all the vertices from edges1, different from v1
679 [ # # ][ # # ]: 0 : std::vector<moab::EntityHandle> otherVertices1;
680 [ # # ]: 0 : for (unsigned int i1 = 0; i1 < edges1.size(); i1++) {
681 [ # # ]: 0 : moab::EntityHandle edgeThatGoes = edges1[i1];
682 : : // find the other vertex, not v1
683 : : int nn;
684 : : const moab::EntityHandle * conn2;
685 [ # # ]: 0 : mb->get_connectivity(edgeThatGoes, conn2, nn);
686 : 0 : int other = 0;
687 [ # # ]: 0 : if (conn2[0] == v1)
688 : 0 : other = 1;
689 [ # # ]: 0 : else if (conn2[1] == v1)// is this really necessary in a correct model?
690 : 0 : other = 0;
691 [ # # ]: 0 : otherVertices1.push_back(conn2[other]);
692 : : }
693 [ # # ]: 0 : for (unsigned int i = 0; i < edges0.size(); i++) {
694 [ # # ]: 0 : moab::EntityHandle edgeThatStays = edges0[i];
695 : : // find the other vertex, not v0
696 : : int nn;
697 : : const moab::EntityHandle * conn2;
698 [ # # ]: 0 : mb->get_connectivity(edgeThatStays, conn2, nn);
699 : 0 : int other = 0;
700 [ # # ]: 0 : if (conn2[0] == v1)
701 : 0 : other = 1;
702 [ # # ]: 0 : else if (conn2[1] == v1)// is this really necessary in a correct model?
703 : 0 : other = 0;
704 : 0 : moab::EntityHandle v2 = conn2[other];
705 : : // let's see now if v2 is among vertices from otherVertices1
706 : : // if yes, then we have a match, edges that need collapsed
707 [ # # ]: 0 : for (unsigned int i1 = 0; i1 < otherVertices1.size(); i1++) {
708 [ # # ][ # # ]: 0 : if (v2 == otherVertices1[i1]) {
709 : : // we have a match, some work to do
710 [ # # ]: 0 : edgePairsToMerge.push_back(edgeThatStays);
711 [ # # ][ # # ]: 0 : edgePairsToMerge.push_back(edges1[i1]);
712 : 0 : break; // we stopp looking for a match, only one possible
713 : : }
714 : : }
715 : : }
716 : :
717 : : // first merge vertices v0 and v1 : will also NOT delete v1 (yet)
718 : : // the tag on v1 will be deleted too, and we do not want that, at least until
719 : : // after the merging of edges, and deleting the pair
720 [ # # ]: 0 : rval = mb->merge_entities(v0, v1, false, false);
721 [ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
722 : : // merge edgePairsToMerge // now, v0 and v1 should be collapsed!
723 [ # # ]: 0 : for (unsigned int j = 0; j < edgePairsToMerge.size(); j += 2) {
724 : : // will also delete edges that contained v1 before
725 [ # # ]: 0 : mb->merge_entities(edgePairsToMerge[j],
726 [ # # ][ # # ]: 0 : edgePairsToMerge[j + 1], false, true);
727 : : }
728 : : // all the triangles connected to v0 are now changed, need recomputation
729 : : //
730 : : rval = mb->get_adjacencies(&v0, 1, 2, false, changed,
731 [ # # ]: 0 : moab::Interface::UNION);
732 [ # # ][ # # ]: 0 : if (moab::MB_SUCCESS != rval) return rval;
[ # # ]
733 : :
734 : 0 : }
735 : : }
736 : 1088 : return moab::MB_SUCCESS;
737 : :
738 [ + - ][ + - ]: 312 : }
739 : :
|
