TopologyKernel.cc 61.7 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42
/*===========================================================================*\
 *                                                                           *
 *                            OpenVolumeMesh                                 *
 *        Copyright (C) 2011 by Computer Graphics Group, RWTH Aachen         *
 *                        www.openvolumemesh.org                             *
 *                                                                           *
 *---------------------------------------------------------------------------*
 *  This file is part of OpenVolumeMesh.                                     *
 *                                                                           *
 *  OpenVolumeMesh is free software: you can redistribute it and/or modify   *
 *  it under the terms of the GNU Lesser General Public License as           *
 *  published by the Free Software Foundation, either version 3 of           *
 *  the License, or (at your option) any later version with the              *
 *  following exceptions:                                                    *
 *                                                                           *
 *  If other files instantiate templates or use macros                       *
 *  or inline functions from this file, or you compile this file and         *
 *  link it with other files to produce an executable, this file does        *
 *  not by itself cause the resulting executable to be covered by the        *
 *  GNU Lesser General Public License. This exception does not however       *
 *  invalidate any other reasons why the executable file might be            *
 *  covered by the GNU Lesser General Public License.                        *
 *                                                                           *
 *  OpenVolumeMesh is distributed in the hope that it will be useful,        *
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of           *
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the            *
 *  GNU Lesser General Public License for more details.                      *
 *                                                                           *
 *  You should have received a copy of the GNU LesserGeneral Public          *
 *  License along with OpenVolumeMesh.  If not,                              *
 *  see <http://www.gnu.org/licenses/>.                                      *
 *                                                                           *
\*===========================================================================*/

/*===========================================================================*\
 *                                                                           *
 *   $Revision$                                                         *
 *   $Date$                    *
 *   $LastChangedBy$                                                *
 *                                                                           *
\*===========================================================================*/

43 44 45 46
#ifndef NDEBUG
#include <iostream>
#endif

47 48
#include <queue>

49
#include "TopologyKernel.hh"
50
#include "../System/FunctionalInclude.hh"
51 52 53 54 55 56 57 58 59 60 61 62

namespace OpenVolumeMesh {

// Initialize constants
const VertexHandle      TopologyKernel::InvalidVertexHandle   = VertexHandle(-1);
const EdgeHandle        TopologyKernel::InvalidEdgeHandle     = EdgeHandle(-1);
const HalfEdgeHandle    TopologyKernel::InvalidHalfEdgeHandle = HalfEdgeHandle(-1);
const FaceHandle        TopologyKernel::InvalidFaceHandle     = FaceHandle(-1);
const HalfFaceHandle    TopologyKernel::InvalidHalfFaceHandle = HalfFaceHandle(-1);
const CellHandle        TopologyKernel::InvalidCellHandle     = CellHandle(-1);

TopologyKernel::TopologyKernel() :
63
    n_vertices_(0u),
64 65
    v_bottom_up_(true),
    e_bottom_up_(true),
66
    f_bottom_up_(true) {
67 68 69 70 71 72 73 74

}

TopologyKernel::~TopologyKernel() {
}

//========================================================================================

75 76 77 78
VertexHandle TopologyKernel::add_vertex() {

    ++n_vertices_;

79
    // Create item for vertex bottom-up incidences
80 81 82 83 84 85 86 87 88 89 90 91 92
    if(v_bottom_up_) {
        outgoing_hes_per_vertex_.resize(n_vertices_);
    }

    // Resize vertex props
    resize_vprops(n_vertices_);

    // Return 0-indexed handle
    return VertexHandle((int)(n_vertices_ - 1));
}

//========================================================================================

93 94
/// Add edge
EdgeHandle TopologyKernel::add_edge(const VertexHandle& _fromVertex,
95 96
                                    const VertexHandle& _toVertex,
                                    bool _allowDuplicates) {
97

98 99 100 101
    // If the conditions are not fulfilled, assert will fail (instead
	// of returning an invalid handle)
    assert(_fromVertex.is_valid() && (size_t)_fromVertex.idx() < n_vertices());
    assert(_toVertex.is_valid() && (size_t)_toVertex.idx() < n_vertices());
102 103

    // Test if edge does not exist, yet
104
    if(!_allowDuplicates) {
Mike Kremer's avatar
Mike Kremer committed
105 106
        if(v_bottom_up_) {

107
            assert((size_t)_fromVertex.idx() < outgoing_hes_per_vertex_.size());
Mike Kremer's avatar
Mike Kremer committed
108 109 110 111 112 113 114 115
            std::vector<HalfEdgeHandle>& ohes = outgoing_hes_per_vertex_[_fromVertex.idx()];
            for(std::vector<HalfEdgeHandle>::const_iterator he_it = ohes.begin(),
                    he_end = ohes.end(); he_it != he_end; ++he_it) {
                if(halfedge(*he_it).to_vertex() == _toVertex) {
                    return edge_handle(*he_it);
                }
            }
        } else {
116
            for(size_t i = 0; i < edges_.size(); ++i) {
Mike Kremer's avatar
Mike Kremer committed
117 118 119 120 121
                if(edge(EdgeHandle(i)).from_vertex() == _fromVertex && edge(EdgeHandle(i)).to_vertex() == _toVertex) {
                    return EdgeHandle(i);
                } else if(edge(EdgeHandle(i)).from_vertex() == _toVertex && edge(EdgeHandle(i)).to_vertex() == _fromVertex) {
                    return EdgeHandle(i);
                }
122
            }
123 124 125 126 127 128 129 130 131 132 133 134
        }
    }

    // Create edge object
    OpenVolumeMeshEdge e(_fromVertex, _toVertex);

    // Store edge locally
    edges_.push_back(e);

    // Resize props
    resize_eprops(n_edges());

135
    EdgeHandle eh((int)edges_.size()-1);
136

137
    // Update vertex bottom-up incidences
138
    if(v_bottom_up_) {
139 140 141
        assert((size_t)_fromVertex.idx() < outgoing_hes_per_vertex_.size());
        assert((size_t)_toVertex.idx() < outgoing_hes_per_vertex_.size());

142 143
        outgoing_hes_per_vertex_[_fromVertex.idx()].push_back(halfedge_handle(eh, 0));
        outgoing_hes_per_vertex_[_toVertex.idx()].push_back(halfedge_handle(eh, 1));
144 145
    }

146
    // Create item for edge bottom-up incidences
147 148 149 150
    if(e_bottom_up_) {
        incident_hfs_per_he_.resize(n_halfedges());
    }

151
    // Get handle of recently created edge
152
    return eh;
153 154 155 156 157 158 159
}

//========================================================================================

/// Add face via incident edges
FaceHandle TopologyKernel::add_face(const std::vector<HalfEdgeHandle>& _halfedges, bool _topologyCheck) {

160
#ifndef NDEBUG
161
    // Assert that halfedges are valid
162
    for(std::vector<HalfEdgeHandle>::const_iterator it = _halfedges.begin(),
163 164
            end = _halfedges.end(); it != end; ++it)
        assert(it->is_valid() && (size_t)it->idx() < edges_.size() * 2u);
165
#endif
166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184

    // Perform topology check
    if(_topologyCheck) {

        /*
         * Test if halfedges are connected
         *
         * The test works as follows:
         * For every edge in the parameter vector
         * put all incident vertices into a
         * set of either "from"-vertices or "to"-vertices,
         * respectively.
         * If and only if all edges are connected,
         * then both sets are identical.
         */

        std::set<VertexHandle> fromVertices;
        std::set<VertexHandle> toVertices;

185 186
        for(std::vector<HalfEdgeHandle>::const_iterator it = _halfedges.begin(),
            end = _halfedges.end(); it != end; ++it) {
187 188 189 190 191

            fromVertices.insert(halfedge(*it).from_vertex());
            toVertices.insert(halfedge(*it).to_vertex());
        }

192 193
        for(std::set<VertexHandle>::const_iterator v_it = fromVertices.begin(),
                v_end = fromVertices.end(); v_it != v_end; ++v_it) {
194
            if(toVertices.count(*v_it) != 1) {
195 196 197 198 199 200
                // The situation here is different, the caller has requested a
                // topology check and expects an invalid handle if the half-edges
                // are not connected. Give him a message in debug mode.
#ifndef NDEBUG
                std::cerr << "add_face(): The specified halfedges are not connected!" << std::endl;
#endif
201 202 203 204 205 206 207 208 209 210 211 212 213
                return InvalidFaceHandle;
            }
        }

        // The halfedges are now guaranteed to be connected
    }

    // Create face
    OpenVolumeMeshFace face(_halfedges);

    faces_.push_back(face);

    // Get added face's handle
214
    FaceHandle fh(faces_.size() - 1);
215 216 217 218

    // Resize props
    resize_fprops(n_faces());

219
    // Update edge bottom-up incidences
220 221 222 223
    if(e_bottom_up_) {

        for(std::vector<HalfEdgeHandle>::const_iterator it = _halfedges.begin(),
            end = _halfedges.end(); it != end; ++it) {
224 225 226 227

            assert((size_t)it->idx() < incident_hfs_per_he_.size());
            assert((size_t)opposite_halfedge_handle(*it).idx() < incident_hfs_per_he_.size());

228 229
            incident_hfs_per_he_[it->idx()].push_back(halfface_handle(fh, 0));
            incident_hfs_per_he_[opposite_halfedge_handle(*it).idx()].push_back(halfface_handle(fh, 1));
230 231 232
        }
    }

233
    // Create item for face bottom-up incidences
234 235 236 237
    if(f_bottom_up_) {
        incident_cell_per_hf_.resize(n_halffaces(), InvalidCellHandle);
    }

238 239 240 241 242 243 244 245 246 247
    // Return handle of recently created face
    return fh;
}

//========================================================================================

/// Add face via incident vertices
/// Define the _vertices in counter-clockwise order (from the "outside")
FaceHandle TopologyKernel::add_face(const std::vector<VertexHandle>& _vertices) {

248
#ifndef NDEBUG
249
    // Assert that all vertices have valid indices
250
    for(std::vector<VertexHandle>::const_iterator it = _vertices.begin(),
251 252
            end = _vertices.end(); it != end; ++it)
        assert(it->is_valid() && (size_t)it->idx() < n_vertices());
253
#endif
254 255 256 257

    // Add edge for each pair of vertices
    std::vector<HalfEdgeHandle> halfedges;
    std::vector<VertexHandle>::const_iterator it = _vertices.begin();
258 259
    std::vector<VertexHandle>::const_iterator end = _vertices.end();
    for(; (it+1) != end; ++it) {
260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283
        EdgeHandle e_idx = add_edge(*it, *(it+1));

        // Swap halfedge if edge already existed and
        // has been initially defined in reverse orientation
        int swap = 0;
        if(edge(e_idx).to_vertex() == *it) swap = 1;

        halfedges.push_back(halfedge_handle(e_idx, swap));
    }
    EdgeHandle e_idx = add_edge(*it, *_vertices.begin());
    int swap = 0;
    if(edge(e_idx).to_vertex() == *it) swap = 1;
    halfedges.push_back(halfedge_handle(e_idx, swap));

    // Add face
#ifndef NDEBUG
    return add_face(halfedges, true);
#else
    return add_face(halfedges, false);
#endif
}

//========================================================================================

284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309
void TopologyKernel::reorder_incident_halffaces(const EdgeHandle& _eh) {

    /* Put halffaces in clockwise order via the
     * same cell property which now exists.
     * Note, this only works for manifold configurations though.
     * Proceed as follows: Pick one starting halfface. Assuming
     * that all halfface normals point into the incident cell,
     * we find the adjacent halfface within the incident cell
     * along the considered halfedge. We set the found halfface
     * to be the one to be processed next. If we reach an outside
     * region, we try to go back from the starting halfface in reverse
     * order. If the complex is properly connected (the pairwise
     * intersection of two adjacent 3-dimensional cells is always
     * a 2-dimensional entity, namely a facet), such an ordering
     * always exists and will be found. If not, a correct order
     * can not be given and, as a result, the related iterators
     * will address the related entities in an arbitrary fashion.
     */

    for(unsigned char s = 0; s <= 1; s++) {

        HalfEdgeHandle cur_he = halfedge_handle(_eh, s);
        std::vector<HalfFaceHandle> new_halffaces;
        HalfFaceHandle start_hf = InvalidHalfFaceHandle;
        HalfFaceHandle cur_hf = InvalidHalfFaceHandle;

310 311
        // Start with one incident halfface and go into the first direction
        assert((size_t)cur_he.idx() < incident_hfs_per_he_.size());
312

313
        if(incident_hfs_per_he_[cur_he.idx()].size() != 0) {
314 315

            // Get start halfface
316
            cur_hf = *incident_hfs_per_he_[cur_he.idx()].begin();
317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337
            start_hf = cur_hf;

            while(cur_hf != InvalidHalfFaceHandle) {

                // Add halfface
                new_halffaces.push_back(cur_hf);

                // Go to next halfface
                cur_hf = adjacent_halfface_in_cell(cur_hf, cur_he);

                if(cur_hf != InvalidHalfFaceHandle)
                    cur_hf = opposite_halfface_handle(cur_hf);

                // End when we're through
                if(cur_hf == start_hf) break;
            }

            // First direction has terminated
            // If new_halffaces has the same size as old (unordered)
            // vector of incident halffaces, we are done here
            // If not, try the other way round
338
            if(new_halffaces.size() != incident_hfs_per_he_[cur_he.idx()].size()) {
339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356

                // Get opposite of start halfface
                cur_hf = start_hf;

                 while(cur_hf != InvalidHalfFaceHandle) {

                     cur_hf = opposite_halfface_handle(cur_hf);
                     cur_hf = adjacent_halfface_in_cell(cur_hf, cur_he);

                     if(cur_hf == start_hf) break;

                     if(cur_hf != InvalidHalfFaceHandle)
                         new_halffaces.insert(new_halffaces.begin(), cur_hf);
                     else break;
                }
            }

            // Everything worked just fine, set the new ordered vector
357 358
            if(new_halffaces.size() == incident_hfs_per_he_[cur_he.idx()].size()) {
                incident_hfs_per_he_[cur_he.idx()] = new_halffaces;
359 360 361 362 363 364 365
            }
        }
    }
}

//========================================================================================

366 367 368
/// Add cell via incident halffaces
CellHandle TopologyKernel::add_cell(const std::vector<HalfFaceHandle>& _halffaces, bool _topologyCheck) {

369
#ifndef NDEBUG
370
    // Assert that halffaces have valid indices
371
    for(std::vector<HalfFaceHandle>::const_iterator it = _halffaces.begin(),
372 373
            end = _halffaces.end(); it != end; ++it)
        assert(it->is_valid() && ((size_t)it->idx() < faces_.size() * 2u));
374
#endif
375 376 377 378 379 380 381 382 383 384 385 386 387 388 389

    // Perform topology check
    if(_topologyCheck) {

        /*
         * Test if all halffaces are connected and form a two-manifold
         * => Cell is closed
         *
         * This test is simple: The number of involved half-edges has to be
         * exactly twice the number of involved edges.
         */

        std::set<HalfEdgeHandle> incidentHalfedges;
        std::set<EdgeHandle>     incidentEdges;

390 391
        for(std::vector<HalfFaceHandle>::const_iterator it = _halffaces.begin(),
                end = _halffaces.end(); it != end; ++it) {
392 393

            OpenVolumeMeshFace hface = halfface(*it);
394 395
            for(std::vector<HalfEdgeHandle>::const_iterator he_it = hface.halfedges().begin(),
                    he_end = hface.halfedges().end(); he_it != he_end; ++he_it) {
396 397 398 399 400 401
                incidentHalfedges.insert(*he_it);
                incidentEdges.insert(edge_handle(*he_it));
            }
        }

        if(incidentHalfedges.size() != (incidentEdges.size() * 2u)) {
402 403 404
#ifndef NDEBUG
            std::cerr << "add_cell(): The specified half-faces are not connected!" << std::endl;
#endif
405 406 407 408 409 410 411 412 413 414 415 416 417 418
            return InvalidCellHandle;
        }

        // The halffaces are now guaranteed to form a two-manifold
    }

    // Create new cell
    OpenVolumeMeshCell cell(_halffaces);

    cells_.push_back(cell);

    // Resize props
    resize_cprops(n_cells());

419
    CellHandle ch((int)cells_.size()-1);
420

421
    // Update face bottom-up incidences
422 423 424 425 426
    if(f_bottom_up_) {

        std::set<EdgeHandle> cell_edges;
        for(std::vector<HalfFaceHandle>::const_iterator it = _halffaces.begin(),
                end = _halffaces.end(); it != end; ++it) {
427
            assert((size_t)it->idx() < incident_cell_per_hf_.size());
428

429
#ifndef NDEBUG
430 431
            if(_topologyCheck) {
                if(incident_cell_per_hf_[it->idx()] != InvalidCellHandle) {
432 433 434 435 436 437
                    // Shouldn't this situation be dealt with before adding the
                    // cell and return InvalidCellHandle in this case?
                	// Mike: Not if the user intends to add non-manifold
                	// configurations. Although, in this case, he should be
                	// warned about it.
                    std::cerr << "add_cell(): One of the specified half-faces is already incident to another cell!" << std::endl;
438 439
                }
            }
440
#endif
441 442

            // Overwrite incident cell for current half-face
443
            incident_cell_per_hf_[it->idx()] = ch;
444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470

            // Collect all edges of cell
            const std::vector<HalfEdgeHandle> hes = halfface(*it).halfedges();
            for(std::vector<HalfEdgeHandle>::const_iterator he_it = hes.begin(),
                    he_end = hes.end(); he_it != he_end; ++he_it) {
                cell_edges.insert(edge_handle(*he_it));
            }
        }

        if(e_bottom_up_) {

            // Try to reorder all half-faces w.r.t.
            // their incident half-edges such that all
            // half-faces are in cyclic order around
            // a half-edge
            for(std::set<EdgeHandle>::const_iterator e_it = cell_edges.begin(),
                    e_end = cell_edges.end(); e_it != e_end; ++e_it) {
                reorder_incident_halffaces(*e_it);
            }
        }
    }

    return ch;
}

//========================================================================================

471 472 473 474 475 476 477 478 479 480 481 482 483 484
/// Set the vertices of an edge
void TopologyKernel::set_edge(const EdgeHandle& _eh, const VertexHandle& _fromVertex, const VertexHandle& _toVertex) {

    Edge& e = edge(_eh);

    // Update bottom-up entries
    if(has_vertex_bottom_up_incidences()) {

        const VertexHandle& fv = e.from_vertex();
        const VertexHandle& tv = e.to_vertex();

        const HalfEdgeHandle heh0 = halfedge_handle(_eh, 0);
        const HalfEdgeHandle heh1 = halfedge_handle(_eh, 1);

Mike Kremer's avatar
Mike Kremer committed
485 486
        std::vector<HalfEdgeHandle>::iterator h_end =
        		std::remove(outgoing_hes_per_vertex_[fv.idx()].begin(), outgoing_hes_per_vertex_[fv.idx()].end(), heh0);
487 488 489 490
        outgoing_hes_per_vertex_[fv.idx()].resize(h_end - outgoing_hes_per_vertex_[fv.idx()].begin());

        h_end = std::remove(outgoing_hes_per_vertex_[tv.idx()].begin(), outgoing_hes_per_vertex_[tv.idx()].end(), heh1);
        outgoing_hes_per_vertex_[tv.idx()].resize(h_end - outgoing_hes_per_vertex_[tv.idx()].begin());
491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516

        outgoing_hes_per_vertex_[_fromVertex.idx()].push_back(heh0);
        outgoing_hes_per_vertex_[_toVertex.idx()].push_back(heh1);
    }

    e.set_from_vertex(_fromVertex);
    e.set_to_vertex(_toVertex);
}

//========================================================================================

/// Set the half-edges of a face
void TopologyKernel::set_face(const FaceHandle& _fh, const std::vector<HalfEdgeHandle>& _hes) {

    Face& f = face(_fh);

    if(has_edge_bottom_up_incidences()) {

        const HalfFaceHandle hf0 = halfface_handle(_fh, 0);
        const HalfFaceHandle hf1 = halfface_handle(_fh, 1);

        const std::vector<HalfEdgeHandle>& hes = f.halfedges();

        for(std::vector<HalfEdgeHandle>::const_iterator he_it = hes.begin(),
                he_end = hes.end(); he_it != he_end; ++he_it) {

Mike Kremer's avatar
Mike Kremer committed
517 518
        	std::vector<HalfFaceHandle>::iterator h_end =
        			std::remove(incident_hfs_per_he_[he_it->idx()].begin(),
519 520 521 522 523 524
                        		incident_hfs_per_he_[he_it->idx()].end(), hf0);
            incident_hfs_per_he_[he_it->idx()].resize(h_end - incident_hfs_per_he_[he_it->idx()].begin());

            h_end =  std::remove(incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].begin(),
                        		 incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].end(), hf1);
            incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].resize(h_end - incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].begin());
525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567
        }

        for(std::vector<HalfEdgeHandle>::const_iterator he_it = _hes.begin(),
                he_end = _hes.end(); he_it != he_end; ++he_it) {

            incident_hfs_per_he_[he_it->idx()].push_back(hf0);
            incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].push_back(hf1);
        }

        // TODO: Reorder incident half-faces
    }

    f.set_halfedges(_hes);
}

//========================================================================================

/// Set the half-faces of a cell
void TopologyKernel::set_cell(const CellHandle& _ch, const std::vector<HalfFaceHandle>& _hfs) {

    Cell& c = cell(_ch);

    if(has_face_bottom_up_incidences()) {

        const std::vector<HalfFaceHandle>& hfs = c.halffaces();
        for(std::vector<HalfFaceHandle>::const_iterator hf_it = hfs.begin(),
                hf_end = hfs.end(); hf_it != hf_end; ++hf_it) {

            incident_cell_per_hf_[*hf_it] = InvalidCellHandle;
        }

        for(std::vector<HalfFaceHandle>::const_iterator hf_it = _hfs.begin(),
                hf_end = _hfs.end(); hf_it != hf_end; ++hf_it) {

            incident_cell_per_hf_[*hf_it] = _ch;
        }
    }

    c.set_halffaces(_hfs);
}

//========================================================================================

568 569 570
/**
 * \brief Delete vertex from mesh
 *
Mike Kremer's avatar
Mike Kremer committed
571
 * Get all incident higher-dimensional entities and delete the complete
572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620
 * subtree of the mesh incident to vertex _h.
 * In this function all incident entities are gathered
 * and deleted using the delete_*_core functions
 * that do the actual deletion including the update
 * of the bottom-up incidences, etc.
 *
 * @param _h The handle to the vertex to be deleted
 */
VertexIter TopologyKernel::delete_vertex(const VertexHandle& _h) {

    std::vector<VertexHandle> vs;
    vs.push_back(_h);

    std::set<EdgeHandle> incidentEdges_s;
    get_incident_edges(vs, incidentEdges_s);

    std::set<FaceHandle> incidentFaces_s;
    get_incident_faces(incidentEdges_s, incidentFaces_s);

    std::set<CellHandle> incidentCells_s;
    get_incident_cells(incidentFaces_s, incidentCells_s);

    // Delete cells
    for(std::set<CellHandle>::const_reverse_iterator c_it = incidentCells_s.rbegin(),
            c_end = incidentCells_s.rend(); c_it != c_end; ++c_it) {
        delete_cell_core(*c_it);
    }

    // Delete faces
    for(std::set<FaceHandle>::const_reverse_iterator f_it = incidentFaces_s.rbegin(),
            f_end = incidentFaces_s.rend(); f_it != f_end; ++f_it) {
        delete_face_core(*f_it);
    }

    // Delete edges
    for(std::set<EdgeHandle>::const_reverse_iterator e_it = incidentEdges_s.rbegin(),
            e_end = incidentEdges_s.rend(); e_it != e_end; ++e_it) {
        delete_edge_core(*e_it);
    }

    // Delete vertex
    return delete_vertex_core(_h);
}

//========================================================================================

/**
 * \brief Delete edge from mesh
 *
Mike Kremer's avatar
Mike Kremer committed
621
 * Get all incident higher-dimensional entities and delete the complete
622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661
 * subtree of the mesh incident to edge _h.
 * In this function all incident entities are gathered
 * and deleted using the delete_*_core functions
 * that do the actual deletion including the update
 * of the bottom-up incidences, etc.
 *
 * @param _h The handle to the edge to be deleted
 */
EdgeIter TopologyKernel::delete_edge(const EdgeHandle& _h) {

    std::vector<EdgeHandle> es;
    es.push_back(_h);

    std::set<FaceHandle> incidentFaces_s;
    get_incident_faces(es, incidentFaces_s);

    std::set<CellHandle> incidentCells_s;
    get_incident_cells(incidentFaces_s, incidentCells_s);

    // Delete cells
    for(std::set<CellHandle>::const_reverse_iterator c_it = incidentCells_s.rbegin(),
            c_end = incidentCells_s.rend(); c_it != c_end; ++c_it) {
        delete_cell_core(*c_it);
    }

    // Delete faces
    for(std::set<FaceHandle>::const_reverse_iterator f_it = incidentFaces_s.rbegin(),
            f_end = incidentFaces_s.rend(); f_it != f_end; ++f_it) {
        delete_face_core(*f_it);
    }

    // Delete edge
    return delete_edge_core(_h);
}

//========================================================================================

/**
 * \brief Delete face from mesh
 *
Mike Kremer's avatar
Mike Kremer committed
662
 * Get all incident higher-dimensional entities and delete the complete
663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 807 808 809 810 811 812 813 814 815 816 817 818 819 820 821 822 823 824 825 826 827 828 829 830 831 832 833 834 835
 * subtree of the mesh incident to face _h.
 * In this function all incident entities are gathered
 * and deleted using the delete_*_core functions
 * that do the actual deletion including the update
 * of the bottom-up incidences, etc.
 *
 * @param _h The handle to the face to be deleted
 */
FaceIter TopologyKernel::delete_face(const FaceHandle& _h) {

    std::vector<FaceHandle> fs;
    fs.push_back(_h);

    std::set<CellHandle> incidentCells_s;
    get_incident_cells(fs, incidentCells_s);

    // Delete cells
    for(std::set<CellHandle>::const_reverse_iterator c_it = incidentCells_s.rbegin(),
            c_end = incidentCells_s.rend(); c_it != c_end; ++c_it) {
        delete_cell_core(*c_it);
    }

    // Delete face
    return delete_face_core(_h);
}

//========================================================================================

/**
 * \brief Delete cell from mesh
 *
 * Since there's no higher dimensional incident
 * entity to a cell, we can safely delete it from the
 * mesh.
 *
 * @param _h The handle to the cell to be deleted
 */
CellIter TopologyKernel::delete_cell(const CellHandle& _h) {

    return delete_cell_core(_h);
}

//========================================================================================

template <class ContainerT>
void TopologyKernel::get_incident_edges(const ContainerT& _vs,
                                        std::set<EdgeHandle>& _es) const {

    _es.clear();

    if(v_bottom_up_) {

        for(typename ContainerT::const_iterator v_it = _vs.begin(),
                v_end = _vs.end(); v_it != v_end; ++v_it) {

            const std::vector<HalfEdgeHandle>& inc_hes = outgoing_hes_per_vertex_[v_it->idx()];

            for(std::vector<HalfEdgeHandle>::const_iterator he_it = inc_hes.begin(),
                    he_end = inc_hes.end(); he_it != he_end; ++he_it) {

                _es.insert(edge_handle(*he_it));
            }
        }
    } else {

        for(typename ContainerT::const_iterator v_it = _vs.begin(),
                v_end = _vs.end(); v_it != v_end; ++v_it) {

            for(EdgeIter e_it = edges_begin(), e_end = edges_end(); e_it != e_end; ++e_it) {

                const Edge& e = edge(*e_it);

                if(e.from_vertex() == *v_it || e.to_vertex() == *v_it) {
                    _es.insert(*e_it);
                }
            }
        }
    }
}

//========================================================================================

template <class ContainerT>
void TopologyKernel::get_incident_faces(const ContainerT& _es,
                                        std::set<FaceHandle>& _fs) const {

    _fs.clear();

    if(e_bottom_up_) {

        for(typename ContainerT::const_iterator e_it = _es.begin(),
                e_end = _es.end(); e_it != e_end; ++e_it) {

            for(HalfEdgeHalfFaceIter hehf_it = hehf_iter(halfedge_handle(*e_it, 0));
                    hehf_it.valid(); ++hehf_it) {

                const FaceHandle fh = face_handle(*hehf_it);

                if(_fs.count(fh) == 0) {
                    _fs.insert(fh);
                }
            }
        }
    } else {

        for(typename ContainerT::const_iterator e_it = _es.begin(),
                e_end = _es.end(); e_it != e_end; ++e_it) {

            for(FaceIter f_it = faces_begin(),
                    f_end = faces_end(); f_it != f_end; ++f_it) {

                const std::vector<HalfEdgeHandle>& hes = face(*f_it).halfedges();

                for(std::vector<HalfEdgeHandle>::const_iterator he_it = hes.begin(),
                        he_end = hes.end(); he_it != he_end; ++he_it) {

                    if(edge_handle(*he_it) == *e_it) {
                        _fs.insert(*f_it);
                        break;
                    }
                }
            }
        }
    }
}

//========================================================================================

template <class ContainerT>
void TopologyKernel::get_incident_cells(const ContainerT& _fs,
                                        std::set<CellHandle>& _cs) const {

    _cs.clear();

    if(f_bottom_up_) {

        for(typename ContainerT::const_iterator f_it = _fs.begin(),
            f_end = _fs.end(); f_it != f_end; ++f_it) {

            const HalfFaceHandle hfh0 = halfface_handle(*f_it, 0);
            const HalfFaceHandle hfh1 = halfface_handle(*f_it, 1);

            const CellHandle c0 = incident_cell(hfh0);
            const CellHandle c1 = incident_cell(hfh1);

            if(c0.is_valid()) _cs.insert(c0);
            if(c1.is_valid()) _cs.insert(c1);
        }
    } else {

        for(typename ContainerT::const_iterator f_it = _fs.begin(),
            f_end = _fs.end(); f_it != f_end; ++f_it) {

            for(CellIter c_it = cells_begin(), c_end = cells_end();
                c_it != c_end; ++c_it) {

                const std::vector<HalfFaceHandle>& hfs = cell(*c_it).halffaces();

                for(std::vector<HalfFaceHandle>::const_iterator hf_it = hfs.begin(),
                        hf_end = hfs.end(); hf_it != hf_end; ++hf_it) {

                    if(face_handle(*hf_it) == *f_it) {
                        _cs.insert(*c_it);
                        break;
                    }
                }
            }
        }
    }
}

//========================================================================================

836 837 838 839 840 841 842 843 844
/**
 * \brief Delete vertex from mesh
 *
 * After performing this operation, all vertices
 * following vertex _h in the array will be accessible
 * through their old handle decreased by one.
 * This function directly fixes the vertex links
 * in all edges. These steps are performed:
 *
845 846
 * 1) Decrease all vertex handles > _h in incident edges
 * 2) Delete entry in bottom-up list: V -> HE
847 848 849 850
 * 3) Delete vertex itself (not necessary here since
 *    a vertex is only represented by a number)
 * 4) Delete property entry
 *
851
 * @param _h A vertex's handle
852
 */
853
VertexIter TopologyKernel::delete_vertex_core(const VertexHandle& _h) {
854

855
    assert(_h.is_valid() && (size_t)_h.idx() < n_vertices());
856 857 858 859

    // 1)
    if(v_bottom_up_) {

860 861
        // Decrease all vertex handles >= _h in all edge definitions
        for(int i = _h.idx(), end = n_vertices(); i < end; ++i) {
862 863 864
            const std::vector<HalfEdgeHandle>& hes = outgoing_hes_per_vertex_[i];
            for(std::vector<HalfEdgeHandle>::const_iterator he_it = hes.begin(),
                    he_end = hes.end(); he_it != he_end; ++he_it) {
865

866 867 868 869 870 871
                Edge& e = edge(edge_handle(*he_it));
                if(e.from_vertex().idx() == i) {
                    e.set_from_vertex(VertexHandle(i-1));
                }
                if(e.to_vertex().idx() == i) {
                    e.set_to_vertex(VertexHandle(i-1));
872 873 874
                }
            }
        }
875

876 877 878 879 880 881 882 883 884 885 886 887 888 889 890 891 892 893
    } else {

        // Iterate over all edges
        for(EdgeIter e_it = edges_begin(), e_end = edges_end();
                e_it != e_end; ++e_it) {

            // Decrease all vertex handles in edge definitions that are greater than _h
            if(edge(*e_it).from_vertex() > _h) {
                edge(*e_it).set_from_vertex(VertexHandle(edge(*e_it).from_vertex().idx() - 1));
            }
            if(edge(*e_it).to_vertex() > _h) {
                edge(*e_it).set_to_vertex(VertexHandle(edge(*e_it).to_vertex().idx() - 1));
            }
        }
    }

    // 2)
    if(v_bottom_up_) {
894
        assert((size_t)_h.idx() < outgoing_hes_per_vertex_.size());
895 896 897 898 899 900 901 902 903 904 905 906 907 908 909 910 911 912 913 914 915 916 917 918
        outgoing_hes_per_vertex_.erase(outgoing_hes_per_vertex_.begin() + _h.idx());
    }

    // 3)
    --n_vertices_;

    // 4)
    vertex_deleted(_h);

    // Iterator to next element in vertex list
    return (vertices_begin() + _h.idx());
}

//========================================================================================

/**
 * \brief Delete edge from mesh
 *
 * After performing this operation, all edges
 * following edge _h in the array will be accessible
 * through their old handle decreased by one.
 * This function directly fixes the edge links
 * in all faces. These steps are performed:
 *
919 920 921 922 923 924 925
 * 1) Delete bottom-up links from incident vertices
 * 2) Decrease all half-edge handles > _h in incident faces
 * 3) Delete entry in bottom-up list: HE -> HF
 * 4) Decrease all half-edge handles > 2*_h.idx() in
 *    vertex bottom-up list
 * 5) Delete edge itself
 * 6) Delete property entry
926
 *
927
 * @param _h An edge's handle
928
 */
929
EdgeIter TopologyKernel::delete_edge_core(const EdgeHandle& _h) {
930

931
    assert(_h.is_valid() && (size_t)_h.idx() < edges_.size());
932 933 934 935 936 937

    // 1)
    if(v_bottom_up_) {

        VertexHandle v0 = edge(_h).from_vertex();
        VertexHandle v1 = edge(_h).to_vertex();
938 939
        assert(v0.is_valid() && (size_t)v0.idx() < outgoing_hes_per_vertex_.size());
        assert(v1.is_valid() && (size_t)v1.idx() < outgoing_hes_per_vertex_.size());
940

941 942 943
        outgoing_hes_per_vertex_[v0.idx()].erase(
                std::remove(outgoing_hes_per_vertex_[v0.idx()].begin(),
                            outgoing_hes_per_vertex_[v0.idx()].end(),
944
                            halfedge_handle(_h, 0)),
945
                            outgoing_hes_per_vertex_[v0.idx()].end());
946

947 948 949
        outgoing_hes_per_vertex_[v1.idx()].erase(
                std::remove(outgoing_hes_per_vertex_[v1.idx()].begin(),
                            outgoing_hes_per_vertex_[v1.idx()].end(),
950
                            halfedge_handle(_h, 1)),
951
                            outgoing_hes_per_vertex_[v1.idx()].end());
952 953 954 955 956
    }

    // 2)
    if(e_bottom_up_) {

957
        assert((size_t)halfedge_handle(_h, 0).idx() < incident_hfs_per_he_.size());
958

959 960
        // Decrease all half-edge handles > he and
        // delete all half-edge handles == he in face definitions
961 962 963
        // Get all faces that need updates
        std::set<FaceHandle> update_faces;
        for(std::vector<std::vector<HalfFaceHandle> >::const_iterator iit =
964
                (incident_hfs_per_he_.begin() + halfedge_handle(_h, 0).idx()),
965 966 967 968 969 970 971 972 973 974 975 976 977 978 979 980 981
                iit_end = incident_hfs_per_he_.end(); iit != iit_end; ++iit) {
            for(std::vector<HalfFaceHandle>::const_iterator it = iit->begin(),
                    end = iit->end(); it != end; ++it) {
                update_faces.insert(face_handle(*it));
            }
        }
        // Update respective handles
        HEHandleCorrection cor(halfedge_handle(_h, 1));
        for(std::set<FaceHandle>::iterator f_it = update_faces.begin(),
                f_end = update_faces.end(); f_it != f_end; ++f_it) {

            std::vector<HalfEdgeHandle> hes = face(*f_it).halfedges();

            // Delete current half-edge from face's half-edge list
            hes.erase(std::remove(hes.begin(), hes.end(), halfedge_handle(_h, 0)), hes.end());
            hes.erase(std::remove(hes.begin(), hes.end(), halfedge_handle(_h, 1)), hes.end());

Jan Möbius's avatar
Jan Möbius committed
982
#if defined(__clang_major__) && (__clang_major__ >= 5)
983 984 985 986 987
            for(std::vector<HalfEdgeHandle>::iterator it = hes.begin(), end = hes.end();
                it != end; ++it) {
                cor.correctValue(*it);
            }
#else
988
            std::for_each(hes.begin(), hes.end(),
989
                          fun::bind(&HEHandleCorrection::correctValue, &cor, fun::placeholders::_1));
990
#endif
991 992 993 994 995 996 997 998
            face(*f_it).set_halfedges(hes);
        }
    } else {

        // Iterate over all faces
        for(FaceIter f_it = faces_begin(), f_end = faces_end();
                f_it != f_end; ++f_it) {

999
            // Get face's half-edges
1000 1001 1002 1003 1004 1005 1006 1007
            std::vector<HalfEdgeHandle> hes = face(*f_it).halfedges();

            // Delete current half-edge from face's half-edge list
            hes.erase(std::remove(hes.begin(), hes.end(), halfedge_handle(_h, 0)), hes.end());
            hes.erase(std::remove(hes.begin(), hes.end(), halfedge_handle(_h, 1)), hes.end());

            // Decrease all half-edge handles greater than _h in face
            HEHandleCorrection cor(halfedge_handle(_h, 1));
Jan Möbius's avatar
Jan Möbius committed
1008
#if defined(__clang_major__) && (__clang_major__ >= 5)
1009 1010 1011 1012 1013
            for(std::vector<HalfEdgeHandle>::iterator it = hes.begin(), end = hes.end();
                    it != end; ++it) {
                cor.correctValue(*it);
            }
#else
1014
            std::for_each(hes.begin(), hes.end(),
1015
                          fun::bind(&HEHandleCorrection::correctValue, &cor, fun::placeholders::_1));
1016
#endif
1017 1018 1019 1020 1021 1022
            face(*f_it).set_halfedges(hes);
        }
    }

    // 3)
    if(e_bottom_up_) {
1023 1024
        assert((size_t)halfedge_handle(_h, 1).idx() < incident_hfs_per_he_.size());

1025 1026 1027 1028 1029 1030 1031
        incident_hfs_per_he_.erase(incident_hfs_per_he_.begin() + halfedge_handle(_h, 1).idx());
        incident_hfs_per_he_.erase(incident_hfs_per_he_.begin() + halfedge_handle(_h, 0).idx());
    }

    // 4)
    if(v_bottom_up_) {
        HEHandleCorrection cor(halfedge_handle(_h, 1));
Jan Möbius's avatar
Jan Möbius committed
1032
#if defined(__clang_major__) && (__clang_major__ >= 5)
1033 1034 1035 1036 1037
        for(std::vector<std::vector<HalfEdgeHandle> >::iterator it = outgoing_hes_per_vertex_.begin(),
                end = outgoing_hes_per_vertex_.end(); it != end; ++it) {
            cor.correctVecValue(*it);
        }
#else
1038 1039
        std::for_each(outgoing_hes_per_vertex_.begin(),
                      outgoing_hes_per_vertex_.end(),
1040
                      fun::bind(&HEHandleCorrection::correctVecValue, &cor, fun::placeholders::_1));
1041
#endif
1042 1043 1044 1045 1046 1047 1048 1049 1050 1051 1052 1053 1054 1055 1056 1057 1058 1059 1060 1061 1062 1063 1064
    }

    // 5)
    edges_.erase(edges_.begin() + _h.idx());

    // 6)
    edge_deleted(_h);

    // Return iterator to next element in list
    return (edges_begin() + _h.idx());
}

//========================================================================================

/**
 * \brief Delete face from mesh
 *
 * After performing this operation, all faces
 * following face _h in the array will be accessible
 * through their old handle decreased by one.
 * This function directly fixes the face links
 * in all cells. These steps are performed:
 *
1065 1066 1067 1068 1069 1070 1071
 * 1) Delete bottom-up links from incident edges
 * 2) Decrease all half-face handles > _h in incident cells
 * 3) Delete entry in bottom-up list: HF -> C
 * 4) Decrease all half-face handles > 2*_h.idx() in
 *    half-edge bottom-up list
 * 5) Delete face itself
 * 6) Delete property entry
1072
 *
1073
 * @param _h An face's handle
1074
 */
1075
FaceIter TopologyKernel::delete_face_core(const FaceHandle& _h) {
1076

1077
    assert(_h.is_valid() && (size_t)_h.idx() < faces_.size());
1078 1079 1080 1081 1082 1083 1084 1085

    // 1)
    if(e_bottom_up_) {

        const std::vector<HalfEdgeHandle>& hes = face(_h).halfedges();
        for(std::vector<HalfEdgeHandle>::const_iterator he_it = hes.begin(),
                he_end = hes.end(); he_it != he_end; ++he_it) {

1086
            assert((size_t)std::max(he_it->idx(), opposite_halfedge_handle(*he_it).idx()) < incident_hfs_per_he_.size());
1087

1088 1089 1090 1091
            incident_hfs_per_he_[he_it->idx()].erase(
                    std::remove(incident_hfs_per_he_[he_it->idx()].begin(),
                                incident_hfs_per_he_[he_it->idx()].end(),
                                halfface_handle(_h, 0)), incident_hfs_per_he_[he_it->idx()].end());
1092 1093


1094 1095 1096 1097
            incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].erase(
                    std::remove(incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].begin(),
                                incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].end(),
                                halfface_handle(_h, 1)), incident_hfs_per_he_[opposite_halfedge_handle(*he_it).idx()].end());
1098 1099 1100 1101 1102 1103 1104
        }
    }

    // 2)
    if(f_bottom_up_) {

        // Decrease all half-face handles > _h in all cells
1105
        // and delete all half-face handles == _h
1106
        std::set<CellHandle> update_cells;
1107
        for(std::vector<CellHandle>::const_iterator c_it = (incident_cell_per_hf_.begin() + halfface_handle(_h, 0).idx()),
1108 1109 1110 1111 1112 1113 1114 1115 1116 1117 1118 1119 1120 1121
                c_end = incident_cell_per_hf_.end(); c_it != c_end; ++c_it) {
            if(!c_it->is_valid()) continue;
            update_cells.insert(*c_it);
        }
        for(std::set<CellHandle>::const_iterator c_it = update_cells.begin(),
                c_end = update_cells.end(); c_it != c_end; ++c_it) {

            std::vector<HalfFaceHandle> hfs = cell(*c_it).halffaces();

            // Delete current half-faces from cell's half-face list
            hfs.erase(std::remove(hfs.begin(), hfs.end(), halfface_handle(_h, 0)), hfs.end());
            hfs.erase(std::remove(hfs.begin(), hfs.end(), halfface_handle(_h, 1)), hfs.end());

            HFHandleCorrection cor(halfface_handle(_h, 1));
Jan Möbius's avatar
Jan Möbius committed
1122
#if defined(__clang_major__) && (__clang_major__ >= 5)
1123 1124 1125 1126 1127
            for(std::vector<HalfFaceHandle>::iterator it = hfs.begin(),
                    end = hfs.end(); it != end; ++it) {
                cor.correctValue(*it);
            }
#else
1128
            std::for_each(hfs.begin(), hfs.end(),
1129
                          fun::bind(&HFHandleCorrection::correctValue, &cor, fun::placeholders::_1));
1130
#endif
1131 1132 1133 1134 1135 1136 1137 1138 1139
            cell(*c_it).set_halffaces(hfs);
        }

    } else {

        // Iterate over all cells
        for(CellIter c_it = cells_begin(), c_end = cells_end(); c_it != c_end; ++c_it) {

            std::vector<HalfFaceHandle> hfs = cell(*c_it).halffaces();
1140

1141 1142 1143 1144 1145
            // Delete current half-faces from cell's half-face list
            hfs.erase(std::remove(hfs.begin(), hfs.end(), halfface_handle(_h, 0)), hfs.end());
            hfs.erase(std::remove(hfs.begin(), hfs.end(), halfface_handle(_h, 1)), hfs.end());

            HFHandleCorrection cor(halfface_handle(_h, 1));
Jan Möbius's avatar
Jan Möbius committed
1146
#if defined(__clang_major__) && (__clang_major__ >= 5)
1147 1148 1149 1150 1151
            for(std::vector<HalfFaceHandle>::iterator it = hfs.begin(),
                    end = hfs.end(); it != end; ++it) {
                cor.correctValue(*it);
            }
#else
1152
            std::for_each(hfs.begin(), hfs.end(),
1153
                          fun::bind(&HFHandleCorrection::correctValue, &cor, fun::placeholders::_1));
1154
#endif
1155 1156 1157 1158 1159 1160
            cell(*c_it).set_halffaces(hfs);
        }
    }

    // 3)
    if(f_bottom_up_) {
1161 1162
        assert((size_t)halfface_handle(_h, 1).idx() < incident_cell_per_hf_.size());

1163 1164 1165 1166 1167 1168 1169
        incident_cell_per_hf_.erase(incident_cell_per_hf_.begin() + halfface_handle(_h, 1).idx());
        incident_cell_per_hf_.erase(incident_cell_per_hf_.begin() + halfface_handle(_h, 0).idx());
    }

    // 4)
    if(e_bottom_up_) {
        HFHandleCorrection cor(halfface_handle(_h, 1));
Jan Möbius's avatar
Jan Möbius committed
1170
#if defined(__clang_major__) && (__clang_major__ >= 5)
1171 1172 1173 1174
        for(std::vector<std::vector<HalfFaceHandle> >::iterator it = incident_hfs_per_he_.begin(), end = incident_hfs_per_he_.end(); it != end; ++it) {
            cor.correctVecValue(*it);
        }
#else
1175 1176
        std::for_each(incident_hfs_per_he_.begin(),
                      incident_hfs_per_he_.end(),
1177
                      fun::bind(&HFHandleCorrection::correctVecValue, &cor, fun::placeholders::_1));
1178
#endif
1179 1180 1181 1182 1183 1184 1185 1186 1187 1188 1189 1190 1191 1192 1193 1194 1195 1196 1197 1198 1199 1200 1201 1202 1203 1204 1205 1206 1207
    }

    // 5)
    faces_.erase(faces_.begin() + _h.idx());

    // 6)
    face_deleted(_h);

    // Return iterator to next element in list
    return (faces_begin() + _h.idx());
}

//========================================================================================

/**
 * \brief Delete cell from mesh
 *
 * After performing this operation, all cells
 * following cell _h in the array will be accessible
 * through their old handle decreased by one.
 * These steps are performed:
 *
 * 1) Delete links in BU: HF -> C
 * 2) Decrease all entries > c in BU: HF -> C
 * 3) Delete cell from storage array
 * 4) Delete property item
 *
 * @param _h A cell handle
 */
1208
CellIter TopologyKernel::delete_cell_core(const CellHandle& _h) {
1209

1210
    assert(_h.is_valid() && (size_t)_h.idx() < cells_.size());
1211 1212 1213 1214 1215 1216

    // 1)
    if(f_bottom_up_) {
        const std::vector<HalfFaceHandle>& hfs = cell(_h).halffaces();
        for(std::vector<HalfFaceHandle>::const_iterator hf_it = hfs.begin(),
                hf_end = hfs.end(); hf_it != hf_end; ++hf_it) {
1217
            assert((size_t)hf_it->idx() < incident_cell_per_hf_.size());
1218

1219
            incident_cell_per_hf_[hf_it->idx()] = InvalidCellHandle;
1220 1221 1222 1223 1224 1225
        }
    }

    // 2)
    if(f_bottom_up_) {
        CHandleCorrection cor(_h);
Jan Möbius's avatar
Jan Möbius committed
1226
#if defined(__clang_major__) && (__clang_major__ >= 5)
1227 1228 1229 1230 1231
        for(std::vector<CellHandle>::iterator it = incident_cell_per_hf_.begin(),
                end = incident_cell_per_hf_.end(); it != end; ++it) {
            cor.correctValue(*it);
        }
#else
1232 1233
        std::for_each(incident_cell_per_hf_.begin(),
                      incident_cell_per_hf_.end(),
1234
                      fun::bind(&CHandleCorrection::correctValue, &cor, fun::placeholders::_1));
1235
#endif
1236 1237 1238 1239 1240 1241 1242 1243 1244
    }

    // 3)
    cells_.erase(cells_.begin() + _h.idx());

    // 4)
    cell_deleted(_h);

    return (cells_begin() + _h.idx());
1245 1246 1247 1248
}

//========================================================================================

1249 1250
void TopologyKernel::delete_multiple_vertices(const std::vector<bool>& _tag) {

1251
    assert(_tag.size() == n_vertices());
1252 1253 1254 1255 1256 1257 1258 1259 1260 1261 1262 1263 1264 1265 1266 1267 1268 1269 1270 1271 1272 1273 1274 1275 1276 1277 1278 1279 1280 1281 1282 1283 1284 1285 1286 1287 1288 1289 1290 1291 1292 1293 1294 1295 1296 1297 1298 1299 1300 1301 1302 1303 1304 1305 1306 1307 1308 1309 1310 1311 1312 1313 1314 1315 1316 1317 1318 1319 1320 1321 1322 1323 1324 1325 1326 1327 1328 1329 1330 1331 1332 1333 1334 1335 1336 1337 1338 1339 1340 1341 1342 1343 1344 1345 1346 1347 1348 1349 1350 1351 1352 1353 1354 1355 1356 1357 1358 1359 1360 1361 1362 1363 1364 1365 1366 1367 1368 1369 1370 1371 1372 1373 1374 1375 1376 1377 1378

    std::vector<int> newIndices(n_vertices(), -1);
    int curIdx = 0;

    std::vector<int>::iterator idx_it = newIndices.begin();
    for(std::vector<bool>::const_iterator t_it = _tag.begin(),
            t_end = _tag.end(); t_it != t_end; ++t_it, ++idx_it) {

        if(!(*t_it)) {
            // Not marked as deleted
            *idx_it = curIdx;
            ++curIdx;
        } else {
            --n_vertices_;
        }
    }

    // Delete properties accordingly
    delete_multiple_vertex_props(_tag);

    EdgeCorrector corrector(newIndices);
    std::for_each(edges_.begin(), edges_.end(), corrector);
}

//========================================================================================

void TopologyKernel::delete_multiple_edges(const std::vector<bool>& _tag) {

    assert(_tag.size() == n_edges());

    std::vector<int> newIndices(n_edges(), -1);
    int curIdx = 0;

    std::vector<Edge> newEdges;

    std::vector<int>::iterator idx_it = newIndices.begin();
    std::vector<Edge>::const_iterator e_it = edges_.begin();

    for(std::vector<bool>::const_iterator t_it = _tag.begin(),
            t_end = _tag.end(); t_it != t_end; ++t_it, ++idx_it, ++e_it) {

        if(!(*t_it)) {
            // Not marked as deleted

            newEdges.push_back(*e_it);

            *idx_it = curIdx;
            ++curIdx;
        }
    }

    // Swap edges
    edges_.swap(newEdges);

    // Delete properties accordingly
    delete_multiple_edge_props(_tag);

    FaceCorrector corrector(newIndices);
    std::for_each(faces_.begin(), faces_.end(), corrector);
}

//========================================================================================

void TopologyKernel::delete_multiple_faces(const std::vector<bool>& _tag) {

    assert(_tag.size() == n_faces());

    std::vector<int> newIndices(n_faces(), -1);
    int curIdx = 0;

    std::vector<Face> newFaces;

    std::vector<int>::iterator idx_it = newIndices.begin();
    std::vector<Face>::const_iterator f_it = faces_.begin();

    for(std::vector<bool>::const_iterator t_it = _tag.begin(),
            t_end = _tag.end(); t_it != t_end; ++t_it, ++idx_it, ++f_it) {

        if(!(*t_it)) {
            // Not marked as deleted

            newFaces.push_back(*f_it);

            *idx_it = curIdx;
            ++curIdx;
        }
    }

    // Swap faces
    faces_.swap(newFaces);

    // Delete properties accordingly
    delete_multiple_face_props(_tag);

    CellCorrector corrector(newIndices);
    std::for_each(cells_.begin(), cells_.end(), corrector);
}

//========================================================================================

void TopologyKernel::delete_multiple_cells(const std::vector<bool>& _tag) {

    assert(_tag.size() == n_cells());

    std::vector<Cell> newCells;

    std::vector<Cell>::const_iterator c_it = cells_.begin();

    for(std::vector<bool>::const_iterator t_it = _tag.begin(),
            t_end = _tag.end(); t_it != t_end; ++t_it, ++c_it) {

        if(!(*t_it)) {
            // Not marked as deleted