2
/**************************************************************************
4
* Regina - A Normal Surface Theory Calculator *
5
* Computational Engine *
7
* Copyright (c) 1999-2014, Ben Burton *
8
* For further details contact Ben Burton (bab@debian.org). *
10
* This program is free software; you can redistribute it and/or *
11
* modify it under the terms of the GNU General Public License as *
12
* published by the Free Software Foundation; either version 2 of the *
13
* License, or (at your option) any later version. *
15
* As an exception, when this program is distributed through (i) the *
16
* App Store by Apple Inc.; (ii) the Mac App Store by Apple Inc.; or *
17
* (iii) Google Play by Google Inc., then that store may impose any *
18
* digital rights management, device limits and/or redistribution *
19
* restrictions that are required by its terms of service. *
21
* This program is distributed in the hope that it will be useful, but *
22
* WITHOUT ANY WARRANTY; without even the implied warranty of *
23
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
24
* General Public License for more details. *
26
* You should have received a copy of the GNU General Public *
27
* License along with this program; if not, write to the Free *
28
* Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, *
29
* MA 02110-1301, USA. *
31
**************************************************************************/
36
#include "census/ngluingpermsearcher.h"
37
#include "triangulation/nedge.h"
38
#include "triangulation/nfacepair.h"
39
#include "triangulation/ntriangulation.h"
40
#include "utilities/boostutils.h"
41
#include "utilities/memutils.h"
45
const char NHyperbolicMinSearcher::ECLASS_TWISTED = 1;
46
const char NHyperbolicMinSearcher::ECLASS_LOWDEG = 2;
48
const char NHyperbolicMinSearcher::dataTag_ = 'h';
50
NHyperbolicMinSearcher::NHyperbolicMinSearcher(const NFacePairing* pairing,
51
const NFacePairing::IsoList* autos, bool orientableOnly,
52
UseGluingPerms use, void* useArgs) :
53
NEulerSearcher(0, pairing, autos, orientableOnly,
54
PURGE_NON_MINIMAL_HYP, use, useArgs) {
57
void NHyperbolicMinSearcher::runSearch(long maxDepth) {
58
unsigned nTets = getNumberOfTetrahedra();
60
// Larger than we will ever see (and in fact grossly so).
61
maxDepth = nTets * 4 + 1;
65
// Search initialisation.
68
// Do we in fact have no permutation at all to choose?
69
if (maxDepth == 0 || pairing_->dest(0, 0).isBoundary(nTets)) {
79
// Is it a partial search that has already finished?
80
if (orderElt == orderSize) {
87
// ---------- Selecting the individual gluing permutations ----------
89
int minOrder = orderElt;
90
int maxOrder = orderElt + maxDepth;
94
while (orderElt >= minOrder) {
95
face = order[orderElt];
96
adj = (*pairing_)[face];
98
// TODO (long-term): Check for cancellation.
100
// Move to the next permutation.
102
// Be sure to preserve the orientation of the permutation if necessary.
103
if ((! orientableOnly_) || adj.facet == 0)
106
permIndex(face) += 2;
108
// Are we out of ideas for this face?
109
if (permIndex(face) >= 6) {
110
// Yep. Head back down to the previous face.
111
permIndex(face) = -1;
115
// Pull apart vertex and edge links at the previous level.
116
if (orderElt >= minOrder) {
117
splitVertexClasses();
124
// We are sitting on a new permutation to try.
125
permIndex(adj) = NPerm4::invS3[permIndex(face)];
127
// Merge edge links and run corresponding tests.
128
if (mergeEdgeClasses()) {
129
// We created an invalid edge.
134
// Merge vertex links and run corresponding tests.
135
mergeResult = mergeVertexClasses();
136
if (mergeResult & VLINK_BAD_EULER) {
137
// Our vertex link can never obtain the correct
138
// Euler characteristic. Stop now.
139
splitVertexClasses();
144
// Fix the orientation if appropriate.
145
if (adj.facet == 0 && orientableOnly_) {
146
// It's the first time we've hit this tetrahedron.
147
if ((permIndex(face) + (face.facet == 3 ? 0 : 1) +
148
(adj.facet == 3 ? 0 : 1)) % 2 == 0)
149
orientation[adj.simp] = -orientation[face.simp];
151
orientation[adj.simp] = orientation[face.simp];
154
// Move on to the next face.
157
// If we're at the end, try the solution and step back.
158
if (orderElt == orderSize) {
159
// We in fact have an entire triangulation.
160
// Run through the automorphisms and check whether our
161
// permutations are in canonical form.
163
use_(this, useArgs_);
165
// Back to the previous face.
168
// Pull apart vertex and edge links at the previous level.
169
if (orderElt >= minOrder) {
170
splitVertexClasses();
174
// Not a full triangulation; just one level deeper.
176
// We've moved onto a new face.
177
// Be sure to get the orientation right.
178
face = order[orderElt];
179
if (orientableOnly_ && pairing_->dest(face).facet > 0) {
180
// permIndex(face) will be set to -1 or -2 as appropriate.
181
adj = (*pairing_)[face];
182
if (orientation[face.simp] == orientation[adj.simp])
187
if ((face.facet == 3 ? 0 : 1) + (adj.facet == 3 ? 0 : 1) == 1)
188
permIndex(face) = (permIndex(face) + 1) % 2;
190
permIndex(face) -= 2;
193
if (orderElt == maxOrder) {
194
// We haven't found an entire triangulation, but we've
195
// gone as far as we need to.
196
// Process it, then step back.
197
use_(this, useArgs_);
199
// Back to the previous face.
200
permIndex(face) = -1;
203
// Pull apart vertex and edge links at the previous level.
204
if (orderElt >= minOrder) {
205
splitVertexClasses();
212
// And the search is over.
214
// Some extra sanity checking.
216
// Our vertex classes had better be 4n standalone vertices.
217
if (nVertexClasses != 4 * nTets)
218
std::cerr << "ERROR: nVertexClasses == "
219
<< nVertexClasses << " at end of search!" << std::endl;
220
for (int i = 0; i < static_cast<int>(nTets) * 4; i++) {
221
if (vertexState[i].parent != -1)
222
std::cerr << "ERROR: vertexState[" << i << "].parent == "
223
<< vertexState[i].parent << " at end of search!"
225
if (vertexState[i].rank != 0)
226
std::cerr << "ERROR: vertexState[" << i << "].rank == "
227
<< vertexState[i].rank << " at end of search!" << std::endl;
228
if (vertexState[i].bdry != 3)
229
std::cerr << "ERROR: vertexState[" << i << "].bdry == "
230
<< vertexState[i].bdry << " at end of search!" << std::endl;
231
if (vertexState[i].euler != 2)
232
std::cerr << "ERROR: vertexState[" << i << "].euler == "
233
<< vertexState[i].euler << " at end of search!"
235
if (vertexState[i].hadEqualRank)
236
std::cerr << "ERROR: vertexState[" << i << "].hadEqualRank == "
237
"true at end of search!" << std::endl;
238
if (vertexState[i].bdryEdges != 3)
239
std::cerr << "ERROR: vertexState[" << i << "].bdryEdges == "
240
<< static_cast<int>(vertexState[i].bdryEdges)
241
<< " at end of search!" << std::endl;
242
if (vertexState[i].bdryNext[0] != i)
243
std::cerr << "ERROR: vertexState[" << i << "].bdryNext[0] == "
244
<< vertexState[i].bdryNext[0] << " at end of search!"
246
if (vertexState[i].bdryNext[1] != i)
247
std::cerr << "ERROR: vertexState[" << i << "].bdryNext[1] == "
248
<< vertexState[i].bdryNext[1] << " at end of search!"
250
if (vertexState[i].bdryTwist[0])
251
std::cerr << "ERROR: vertexState[" << i << "].bdryTwist == "
252
"true at end of search!" << std::endl;
253
if (vertexState[i].bdryTwist[1])
254
std::cerr << "ERROR: vertexState[" << i << "].bdryTwist == "
255
"true at end of search!" << std::endl;
257
for (unsigned i = 0; i < nTets * 8; i++)
258
if (vertexStateChanged[i] != VLINK_JOIN_INIT)
259
std::cerr << "ERROR: vertexStateChanged[" << i << "] == "
260
<< vertexStateChanged[i] << " at end of search!"
263
// And our edge classes had better be 6n standalone edges.
264
if (nEdgeClasses != 6 * nTets)
265
std::cerr << "ERROR: nEdgeClasses == "
266
<< nEdgeClasses << " at end of search!" << std::endl;
267
for (unsigned i = 0; i < nTets * 6; i++) {
268
if (edgeState[i].parent != -1)
269
std::cerr << "ERROR: edgeState[" << i << "].parent == "
270
<< edgeState[i].parent << " at end of search!"
272
if (edgeState[i].rank != 0)
273
std::cerr << "ERROR: edgeState[" << i << "].rank == "
274
<< edgeState[i].rank << " at end of search!" << std::endl;
275
if (edgeState[i].size != 1)
276
std::cerr << "ERROR: edgeState[" << i << "].size == "
277
<< edgeState[i].size << " at end of search!" << std::endl;
278
if (! edgeState[i].bounded)
279
std::cerr << "ERROR: edgeState[" << i << "].bounded == "
280
"false at end of search!" << std::endl;
281
if (edgeState[i].hadEqualRank)
282
std::cerr << "ERROR: edgeState[" << i << "].hadEqualRank == "
283
"true at end of search!" << std::endl;
285
for (unsigned i = 0; i < nTets * 8; i++)
286
if (edgeStateChanged[i] != -1)
287
std::cerr << "ERROR: edgeStateChanged[" << i << "] == "
288
<< edgeStateChanged[i] << " at end of search!"
295
void NHyperbolicMinSearcher::dumpData(std::ostream& out) const {
296
NEulerSearcher::dumpData(out);
299
NHyperbolicMinSearcher::NHyperbolicMinSearcher(std::istream& in,
300
UseGluingPerms use, void* useArgs) :
301
NEulerSearcher(in, use, useArgs) {
305
// Did we hit an unexpected EOF?
310
int NHyperbolicMinSearcher::mergeEdgeClasses() {
312
* As well as detecting edges that are self-identified in reverse,
313
* we strip out low-degree edges here. Although we are also interested
314
* in non-geometric triangulations, we can still ignore triangulations
315
* with low-degree edges because (with a little work) they can be
316
* proven to be non-minimal. For details see:
317
* "The cusped hyperbolic census is complete", B.B.
319
NTetFace face = order[orderElt];
320
NTetFace adj = (*pairing_)[face];
324
NPerm4 p = gluingPerm(face);
334
char parentTwists, hasTwist;
335
for (v2 = 0; v2 < 4; v2++) {
341
// Look at the edge opposite v1-v2.
342
e = 5 - NEdge::edgeNumber[v1][v2];
343
f = 5 - NEdge::edgeNumber[w1][w2];
345
orderIdx = v2 + 4 * orderElt;
347
// We declare the natural orientation of an edge to be smaller
348
// vertex to larger vertex.
349
hasTwist = (p[NEdge::edgeVertex[e][0]] > p[NEdge::edgeVertex[e][1]] ?
353
eRep = findEdgeClass(e + 6 * face.simp, parentTwists);
354
fRep = findEdgeClass(f + 6 * adj.simp, parentTwists);
357
edgeState[eRep].bounded = false;
359
if (edgeState[eRep].size <= 2)
360
retVal |= ECLASS_LOWDEG;
361
else if (edgeState[eRep].size == 3) {
362
// Flag as LOWDEG only if three distinct tetrahedra are used.
363
middleTet = pairing_->dest(face.simp, v2).simp;
364
if (face.simp != adj.simp && adj.simp != middleTet &&
365
middleTet != face.simp)
366
retVal |= ECLASS_LOWDEG;
368
if (hasTwist ^ parentTwists)
369
retVal |= ECLASS_TWISTED;
371
edgeStateChanged[orderIdx] = -1;
373
if (edgeState[eRep].rank < edgeState[fRep].rank) {
374
// Join eRep beneath fRep.
375
edgeState[eRep].parent = fRep;
376
edgeState[eRep].twistUp = hasTwist ^ parentTwists;
377
edgeState[fRep].size += edgeState[eRep].size;
379
edgeStateChanged[orderIdx] = eRep;
381
// Join fRep beneath eRep.
382
edgeState[fRep].parent = eRep;
383
edgeState[fRep].twistUp = hasTwist ^ parentTwists;
384
if (edgeState[eRep].rank == edgeState[fRep].rank) {
385
edgeState[eRep].rank++;
386
edgeState[fRep].hadEqualRank = true;
388
edgeState[eRep].size += edgeState[fRep].size;
390
edgeStateChanged[orderIdx] = fRep;
400
void NHyperbolicMinSearcher::splitEdgeClasses() {
401
NTetFace face = order[orderElt];
410
for (v2 = 3; v2 >= 0; v2--) {
414
// Look at the edge opposite v1-v2.
415
e = 5 - NEdge::edgeNumber[v1][v2];
417
eIdx = e + 6 * face.simp;
418
orderIdx = v2 + 4 * orderElt;
420
if (edgeStateChanged[orderIdx] < 0)
421
edgeState[findEdgeClass(eIdx)].bounded = true;
423
subRep = edgeStateChanged[orderIdx];
424
rep = edgeState[subRep].parent;
426
edgeState[subRep].parent = -1;
427
if (edgeState[subRep].hadEqualRank) {
428
edgeState[subRep].hadEqualRank = false;
429
edgeState[rep].rank--;
432
edgeState[rep].size -= edgeState[subRep].size;
434
edgeStateChanged[orderIdx] = -1;
440
} // namespace regina
added
removed
engine/census/hyperbolic.cpp
