HexahedralMeshTopologyKernel.cc 16.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
43
44
45
46
/*===========================================================================*\
 *                                                                           *
 *                            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$                                                *
 *                                                                           *
\*===========================================================================*/

#include "HexahedralMeshTopologyKernel.hh"

namespace OpenVolumeMesh {

47
48

FaceHandle HexahedralMeshTopologyKernel::add_face(const std::vector<HalfEdgeHandle>& _halfedges, bool _topologyCheck) {
49
50

    if(_halfedges.size() != 4) {
51
52
53
54
#ifndef NDEBUG
        std::cerr << "HexahedralMeshTopologyKernel::add_face(): Face valence is not four! Returning" << std::endl;
        std::cerr << "invalid handle." << std::endl;
#endif
55
        return TopologyKernel::InvalidFaceHandle;
56
57
    }

58
    return TopologyKernel::add_face(_halfedges, _topologyCheck);
59
60
61
62
}

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

63

64
FaceHandle
65
HexahedralMeshTopologyKernel::add_face(const std::vector<VertexHandle>& _vertices) {
66
67

    if(_vertices.size() != 4) {
68
69
70
71
#ifndef NDEBUG
        std::cerr << "HexahedralMeshTopologyKernel::add_face(): Face valence is not four! Returning" << std::endl;
        std::cerr << "invalid handle." << std::endl;
#endif
72
        return TopologyKernel::InvalidFaceHandle;
73
74
    }

75
    return TopologyKernel::add_face(_vertices);
76
77
78
79
}

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

80

81
CellHandle
82
HexahedralMeshTopologyKernel::add_cell(const std::vector<HalfFaceHandle>& _halffaces, bool _topologyCheck) {
83
84

    if(_halffaces.size() != 6) {
85
86
// To make this consistent with add_face
#ifndef NDEBUG
87
        std::cerr << "Cell valence is not six! Aborting." << std::endl;
88
#endif
89
        return TopologyKernel::InvalidCellHandle;
90
    }
91
    for(std::vector<HalfFaceHandle>::const_iterator it = _halffaces.begin();
92
            it != _halffaces.end(); ++it) {
93
        if(TopologyKernel::halfface(*it).halfedges().size() != 4) {
94
#ifndef NDEBUG
95
            std::cerr << "Incident face does not have valence four! Aborting." << std::endl;
96
#endif
97
            return TopologyKernel::InvalidCellHandle;
98
99
100
101
102
103
104
        }
    }

    // Create new halffaces vector
    std::vector<HalfFaceHandle> ordered_halffaces;

    // The user wants the faces to be reordered
105
106
107
108
109
110
111
    if(_topologyCheck) {

        // Are faces in correct ordering?
        bool ordered = check_halfface_ordering(_halffaces);

        if(!ordered) {

112
#ifndef NDEBUG
113
            std::cerr << "The specified half-faces are not in correct order. Trying automatic re-ordering." << std::endl;
114
#endif
115
116
117

            // Ordering array (see below for details)
            const int orderTop[] = {2, 4, 3, 5};
118
            //const int orderBot[] = {3, 4, 2, 5};
119

120
            ordered_halffaces.resize(6, TopologyKernel::InvalidHalfFaceHandle);
121

122
123
            // Create top side
            ordered_halffaces[0] = _halffaces[0];
124

125
126
127
128
129
            // Go over all incident halfedges
            std::vector<HalfEdgeHandle> hes = TopologyKernel::halfface(ordered_halffaces[0]).halfedges();
            unsigned int idx = 0;
            for(std::vector<HalfEdgeHandle>::const_iterator he_it = hes.begin();
                    he_it != hes.end(); ++he_it) {
130

131
132
                HalfFaceHandle ahfh = get_adjacent_halfface(ordered_halffaces[0], *he_it, _halffaces);
                if(ahfh == TopologyKernel::InvalidHalfFaceHandle) {
133
#ifndef NDEBUG
134
                    std::cerr << "The current halfface is invalid!" << std::endl;
135
#endif
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
                    continue;
                }
                ordered_halffaces[orderTop[idx]] = ahfh;
                ++idx;
            }

            // Now set bottom-halfface
            HalfFaceHandle cur_hf = ordered_halffaces[0];
            HalfEdgeHandle cur_he = *(TopologyKernel::halfface(cur_hf).halfedges().begin());
            cur_hf = get_adjacent_halfface(cur_hf, cur_he, _halffaces);
            cur_he = TopologyKernel::opposite_halfedge_handle(cur_he);
            cur_he = TopologyKernel::next_halfedge_in_halfface(cur_he, cur_hf);
            cur_he = TopologyKernel::next_halfedge_in_halfface(cur_he, cur_hf);
            cur_hf = get_adjacent_halfface(cur_hf, cur_he, _halffaces);

            if(cur_hf != TopologyKernel::InvalidHalfFaceHandle) {
                ordered_halffaces[1] = cur_hf;
            } else {
154
#ifndef NDEBUG
155
                std::cerr << "The current halfface is invalid!" << std::endl;
156
#endif
157
                return TopologyKernel::InvalidCellHandle;
158
159
160
            }

        } else {
161
162
            // Assume right ordering at the user's risk
            ordered_halffaces = _halffaces;
163
164
165
166
167
168
169
        }

    } else {
        // Just copy the original ones
        ordered_halffaces = _halffaces;
    }

170
171
172
173
174
175
176
    return TopologyKernel::add_cell(ordered_halffaces, _topologyCheck);
}

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

bool HexahedralMeshTopologyKernel::check_halfface_ordering(const std::vector<HalfFaceHandle>& _hfs) const {

177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
    /*
     * The test works as follows: Test for both the first and second face in the list,
     * whether the following order holds (clockwise):
     *
     * View from above, outside.
     *           ____
     *          | 4  |
     *      ____|____|____
     *     | 5  | 1  | 6  |
     *     |____|____|____|
     *          | 3  |
     *          |____|
     *
     * View from below, outside.
     *           ____
     *          | 3  |
     *      ____|____|____
     *     | 5  | 2  | 6  |
     *     |____|____|____|
     *          | 4  |
     *          |____|
     */

200
201
202
203
204
    const int orderTop[] = {2, 4, 3, 5};
    const int orderBot[] = {3, 4, 2, 5};

    HalfFaceHandle hfhTop = _hfs[0];
    HalfFaceHandle hfhBot = _hfs[1];
205

206
207
    std::vector<HalfEdgeHandle> halfedgesTop = TopologyKernel::halfface(_hfs[0]).halfedges();
    std::vector<HalfEdgeHandle> halfedgesBot = TopologyKernel::halfface(_hfs[1]).halfedges();
208
209
210
211
212

    int offsetTop = -1;
    int offsetBot = -1;

    // Traverse halfedges top
213
    for(std::vector<HalfEdgeHandle>::const_iterator it = halfedgesTop.begin();
214
215
            it != halfedgesTop.end(); ++it) {

216
        HalfFaceHandle ahfh = get_adjacent_halfface(hfhTop, *it, _hfs);
217
218

        if(offsetTop == -1) {
219
220
221
222
            if(ahfh == _hfs[2])       offsetTop = 0;
            else if(ahfh == _hfs[4])  offsetTop = 1;
            else if(ahfh == _hfs[3])  offsetTop = 2;
            else if(ahfh == _hfs[5])  offsetTop = 3;
223
224
        } else {
            offsetTop = (offsetTop + 1) % 4;
225
            if(ahfh != _hfs[orderTop[offsetTop]]) {
226
#ifndef NDEBUG
227
                std::cerr << "Faces not in right order!" << std::endl;
228
#endif
229
                return false;
230
231
232
233
234
            }
        }
    }

    if(offsetTop == -1) {
235
#ifndef NDEBUG
236
        std::cerr << "Faces not in right order!" << std::endl;
237
#endif
238
        return false;
239
240
241
    }

    // Traverse halfedges bottom
242
    for(std::vector<HalfEdgeHandle>::const_iterator it = halfedgesBot.begin();
243
244
            it != halfedgesBot.end(); ++it) {

245
        HalfFaceHandle ahfh = get_adjacent_halfface(hfhBot, *it, _hfs);
246
247

        if(offsetBot == -1) {
248
249
250
251
            if(ahfh == _hfs[3])       offsetBot = 0;
            else if(ahfh == _hfs[4])  offsetBot = 1;
            else if(ahfh == _hfs[2])  offsetBot = 2;
            else if(ahfh == _hfs[5])  offsetBot = 3;
252
253
        } else {
            offsetBot = (offsetBot + 1) % 4;
254
            if(ahfh != _hfs[orderBot[offsetBot]]) {
255
#ifndef NDEBUG
256
                std::cerr << "Faces not in right order!" << std::endl;
257
#endif
258
                return false;
259
260
261
262
263
            }
        }
    }

    if(offsetBot == -1) {
264
#ifndef NDEBUG
265
        std::cerr << "Faces not in right order!" << std::endl;
266
#endif
267
        return false;
268
269
    }

270
    return true;
271
272
273
274
}

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

275
CellHandle
276
HexahedralMeshTopologyKernel::add_cell(const std::vector<VertexHandle>& _vertices, bool _topologyCheck) {
277

278
279
    // debug mode checks
    assert(TopologyKernel::has_full_bottom_up_incidences());
280
281
    assert(_vertices.size() == 8);

282
    // release mode checks
283
    if(!TopologyKernel::has_full_bottom_up_incidences()) {
284
285
286
287
288
289
290
        return CellHandle(-1);
    }

    if(_vertices.size() != 8) {
        return CellHandle(-1);
    }

291
    HalfFaceHandle hf0, hf1, hf2, hf3, hf4, hf5;
292
293
294
295

    std::vector<VertexHandle> vs;

    // Half-face XF
296
    vs.push_back(_vertices[3]);
297
298
299
    vs.push_back(_vertices[2]);
    vs.push_back(_vertices[1]);
    vs.push_back(_vertices[0]);
300
    hf0 = TopologyKernel::halfface_extensive(vs); vs.clear();
301
302

    // Half-face XB
303
    vs.push_back(_vertices[7]);
304
305
306
    vs.push_back(_vertices[6]);
    vs.push_back(_vertices[5]);
    vs.push_back(_vertices[4]);
307
    hf1 = TopologyKernel::halfface_extensive(vs); vs.clear();
308
309
310
311
312

    // Half-face YF
    vs.push_back(_vertices[1]);
    vs.push_back(_vertices[2]);
    vs.push_back(_vertices[6]);
313
314
    vs.push_back(_vertices[7]);
    hf2 = TopologyKernel::halfface_extensive(vs); vs.clear();
315
316

    // Half-face YB
317
    vs.push_back(_vertices[4]);
318
319
320
    vs.push_back(_vertices[5]);
    vs.push_back(_vertices[3]);
    vs.push_back(_vertices[0]);
321
    hf3 = TopologyKernel::halfface_extensive(vs); vs.clear();
322
323

    // Half-face ZF
324
    vs.push_back(_vertices[1]);
325
326
327
    vs.push_back(_vertices[7]);
    vs.push_back(_vertices[4]);
    vs.push_back(_vertices[0]);
328
    hf4 = TopologyKernel::halfface_extensive(vs); vs.clear();
329
330

    // Half-face ZB
331
    vs.push_back(_vertices[2]);
332
333
334
    vs.push_back(_vertices[3]);
    vs.push_back(_vertices[5]);
    vs.push_back(_vertices[6]);
335
    hf5 = TopologyKernel::halfface_extensive(vs);
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390

    if(!hf0.is_valid()) {

        vs.clear();
        vs.push_back(_vertices[3]); vs.push_back(_vertices[2]);
        vs.push_back(_vertices[1]); vs.push_back(_vertices[0]);
        FaceHandle fh = TopologyKernel::add_face(vs);
        hf0 = halfface_handle(fh, 0);
    }

    if(!hf1.is_valid()) {

        vs.clear();
        vs.push_back(_vertices[7]); vs.push_back(_vertices[6]);
        vs.push_back(_vertices[5]); vs.push_back(_vertices[4]);
        FaceHandle fh = TopologyKernel::add_face(vs);
        hf1 = halfface_handle(fh, 0);
    }

    if(!hf2.is_valid()) {

        vs.clear();
        vs.push_back(_vertices[1]); vs.push_back(_vertices[2]);
        vs.push_back(_vertices[6]); vs.push_back(_vertices[7]);
        FaceHandle fh = TopologyKernel::add_face(vs);
        hf2 = halfface_handle(fh, 0);
    }

    if(!hf3.is_valid()) {

        vs.clear();
        vs.push_back(_vertices[4]); vs.push_back(_vertices[5]);
        vs.push_back(_vertices[3]); vs.push_back(_vertices[0]);
        FaceHandle fh = TopologyKernel::add_face(vs);
        hf3 = halfface_handle(fh, 0);
    }

    if(!hf4.is_valid()) {

        vs.clear();
        vs.push_back(_vertices[1]); vs.push_back(_vertices[7]);
        vs.push_back(_vertices[4]); vs.push_back(_vertices[0]);
        FaceHandle fh = TopologyKernel::add_face(vs);
        hf4 = halfface_handle(fh, 0);
    }

    if(!hf5.is_valid()) {

        vs.clear();
        vs.push_back(_vertices[2]); vs.push_back(_vertices[3]);
        vs.push_back(_vertices[5]); vs.push_back(_vertices[6]);
        FaceHandle fh = TopologyKernel::add_face(vs);
        hf5 = halfface_handle(fh, 0);
    }

Mike Kremer's avatar
Mike Kremer committed
391
392
    assert(hf0.is_valid()); assert(hf1.is_valid()); assert(hf2.is_valid());
    assert(hf3.is_valid()); assert(hf4.is_valid()); assert(hf5.is_valid());
393

394

395
396
397
398
    std::vector<HalfFaceHandle> hfs;
    hfs.push_back(hf0); hfs.push_back(hf1); hfs.push_back(hf2);
    hfs.push_back(hf3); hfs.push_back(hf4); hfs.push_back(hf5);

399
    if (_topologyCheck) {
400
        /*
401
402
403
404
405
406
407
408
409
410
411
        * 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;

        for(std::vector<HalfFaceHandle>::const_iterator it = hfs.begin(),
412
413
414
415
416
417
418
419
                end = hfs.end(); it != end; ++it) {

            OpenVolumeMeshFace hface = halfface(*it);
            for(std::vector<HalfEdgeHandle>::const_iterator he_it = hface.halfedges().begin(),
                    he_end = hface.halfedges().end(); he_it != he_end; ++he_it) {
                incidentHalfedges.insert(*he_it);
                incidentEdges.insert(edge_handle(*he_it));
            }
420
421
422
        }

        if(incidentHalfedges.size() != (incidentEdges.size() * 2u)) {
423
#ifndef NDEBUG
424
            std::cerr << "The specified halffaces are not connected!" << std::endl;
425
#endif
426
            return InvalidCellHandle;
427
428
        }
        // The halffaces are now guaranteed to form a two-manifold
429
430
431
432
433
434

        if(has_face_bottom_up_incidences()) {

            for(std::vector<HalfFaceHandle>::const_iterator it = hfs.begin(),
                    end = hfs.end(); it != end; ++it) {
                if(incident_cell(*it) != InvalidCellHandle) {
435
#ifndef NDEBUG
436
                    std::cerr << "Warning: One of the specified half-faces is already incident to another cell!" << std::endl;
437
#endif
438
439
440
441
442
443
                    return InvalidCellHandle;
                }
            }

        }

444
445
    }

446
    return TopologyKernel::add_cell(hfs, false);
447
448
449
}

//========================================================================================
450

451
const HalfFaceHandle&
452
HexahedralMeshTopologyKernel::get_adjacent_halfface(const HalfFaceHandle& _hfh, const HalfEdgeHandle& _heh,
453
454
455
456
        const std::vector<HalfFaceHandle>& _halffaces) const {

    // Search for halfface that is incident to the opposite
    // halfedge of _heh
457
    HalfEdgeHandle o_he = TopologyKernel::opposite_halfedge_handle(_heh);
458

459
    for(std::vector<HalfFaceHandle>::const_iterator it = _halffaces.begin();
460
461
            it != _halffaces.end(); ++it) {
        if(*it == _hfh) continue;
462
        std::vector<HalfEdgeHandle> halfedges = TopologyKernel::halfface(*it).halfedges();
463
        for(std::vector<HalfEdgeHandle>::const_iterator h_it = halfedges.begin();
464
465
466
467
468
                h_it != halfedges.end(); ++h_it) {
            if(*h_it == o_he) return *it;
        }
    }

469
    return TopologyKernel::InvalidHalfFaceHandle;
470
471
472
}

} // Namespace OpenVolumeMesh