3d:: all the regular honeycombs

cell.cpp

@@ -232,7 +232,7 @@ void initcells() {
   else if(fulltorus) currentmap = new hrmap_torus;
   else if(euclid && DIM == 3) currentmap = euclid3::new_map();
   else if(euclid) currentmap = new hrmap_euclidean;
-  else if(sphere && DIM == 3) currentmap = new sphere3::hrmap_spherical3;
+  else if(DIM == 3 && !binarytiling) currentmap = reg3::new_map();
   else if(sphere) currentmap = new hrmap_spherical;
   else if(quotient) currentmap = new quotientspace::hrmap_quotient;
   else currentmap = new hrmap_hyperbolic;
@@ -442,6 +442,7 @@ int celldistAlt(cell *c) {
   #if MAXMDIM == 4
   if(euclid && DIM == 3) return euclid3::dist_alt(c);
   #endif
+  if(hyperbolic && DIM == 3) return reg3::dist_alt(c);
   if(!c->master->alt) return 0;
   #if CAP_IRR
   if(IRREGULAR) return irr::celldist(c, true);
@@ -875,6 +876,8 @@ int celldistance(cell *c1, cell *c2) {
   if(euclid && DIM == 3) 
     return euclid3::celldistance(c1, c2);
   
+  if(hyperbolic && DIM == 3) return reg3::celldistance(c1, c2);
+
   return hyperbolic_celldistance(c1, c2);
   }
 

classes.cpp

@@ -1771,11 +1771,22 @@ vector<geometryinfo> ginf = {
   {"{6,4}", "Crystal",  "dimensional crystal",                        "Crystal",  6, 4, qANYQ,     gcHyperbolic, 0x28000, {{5, 3}}, eVariation::pure},
   {"{3,4}", "none",     "{3,4} (octahedron)",                         "4x3",      3, 4, qsSMALLB,  gcSphere,     0x28200, {{SEE_ALL, SEE_ALL}}, eVariation::bitruncated},
   {"bin3",  "none",     "3D binary tiling",                           "binary3",  9, 4, 0,         gcHyperbolic, 0x30000, {{7, 3}}, eVariation::pure},
-  {"cube",  "none",     "3D cube tiling",                             "cube",     6, 4, 0,         gcEuclid,     0x30200, {{7, 5}}, eVariation::pure},
+  {"{4,3,4}","none",     "3D cube tiling",                            "cube",     6, 4, 0,         gcEuclid,     0x30200, {{7, 5}}, eVariation::pure},
-  {"120c",  "none",     "120-cell",                                   "120c",    12, 4, qsSMALLB,  gcSphere,     0x30400, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{5,3,3}","none",     "{5,3,3} 120-cell",                          "533",     12, 3, qsSMALLB,  gcSphere,     0x30400, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
-  {"e120c", "elliptic", "120-cell (elliptic space)",                  "e120c",   12, 4, qsSMALLBE, gcSphere,     0x30600, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{5,3,3}", "elliptic", "{5,3,3} 120-cell (elliptic space)",        "e533",    12, 3, qsSMALLBE, gcSphere,     0x30600, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
   {"rhombic","none",    "rhombic dodecahedral honeycomb",             "rhombic", 12, 4, 0,         gcEuclid,     0x31000, {{7, 5}}, eVariation::pure},
   {"bitrunc","none",    "bitruncated cubic honeycomb",                "bitrunc", 14, 3, 0,         gcEuclid,     0x31200, {{7, 5}}, eVariation::pure},
+  {"{5,3,4}","none",    "{5,3,4} hyperbolic honeycomb",               "534",     12, 4, 0,         gcHyperbolic, 0x31400, {{7, 1}}, eVariation::pure},
+  {"{4,3,5}","none",    "{4,3,5} hyperbolic honeycomb",               "435",      6, 4, 0,         gcHyperbolic, 0x31600, {{7, 1}}, eVariation::pure},
+  {"{3,3,3}","none",    "{3,3,3} 5-cell",                             "333",      4, 3, qsSMALLB,  gcSphere,     0x38000, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{3,3,4}","none",    "{3,3,4} 16-cell",                            "334",      4, 3, qsSMALLB,  gcSphere,     0x38200, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{3,3,4}","elliptic","{3,3,4} 16-cell (elliptic)",                 "e334",     4, 3, qsSMALLBE, gcSphere,     0x39200, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{4,3,3}","none",    "{3,3,4} 8-cell",                             "433",      6, 3, qsSMALLB,  gcSphere,     0x38400, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{4,3,3}","elliptic","{3,3,4} 8-cell (elliptic)",                  "e433",     6, 3, qsSMALLBE, gcSphere,     0x39400, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{3,4,3}","none",    "{3,4,3} 24-cell",                            "343",      8, 3, qsSMALLB,  gcSphere,     0x38600, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{3,4,3}","elliptic","{3,4,3} 24-cell (elliptic)",                 "e343",     8, 3, qsSMALLBE, gcSphere,     0x39600, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{3,3,5}","none",    "{3,3,5} 600-cell",                           "335",      4, 3, qsSMALLB,  gcSphere,     0x30800, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
+  {"{3,3,5}","elliptic","{3,3,5} 600-cell (elliptic)",                "e335",     4, 3, qsSMALLBE, gcSphere,     0x31800, {{SEE_ALL, SEE_ALL}}, eVariation::pure},
   };
 
   // bits: 9, 10, 15, 16, (reserved for later) 17, 18

classes.h

@@ -216,6 +216,12 @@ enum eGeometry {
   gKleinQuartic, gBolza, gBolza2, gMinimal, gBinaryTiling, gArchimedean, 
   gMacbeath, gBring, gSchmutzM2, gSchmutzM3, gCrystal, gOctahedron, 
   gBinary3, gCubeTiling, gCell120, gECell120, gRhombic3, gBitrunc3, 
+  gSpace534, gSpace435, 
+  gCell5, 
+  gCell8, gECell8,
+  gCell16, gECell16,
+  gCell24, gECell24,
+  gCell600, gECell600,
   gGUARD};
 
 enum eGeometryClass { gcHyperbolic, gcEuclid, gcSphere };

compileunits.h

@@ -35,6 +35,7 @@
 #include "sphere.cpp"
 #include "quotient.cpp"
 #include "crystal.cpp"
+#include "reg3.cpp"
 #include "language.cpp"
 #include "cell.cpp"
 #include "expansion.cpp"

game.cpp

@@ -134,7 +134,7 @@ vector<int> reachedfrom;
 // it remembers from where we have got to this location
 // the opposite cell will be added to the queue first,
 // which helps the AI
-vector<cell*> movesofgood[8];
+vector<cell*> movesofgood[MAX_EDGE+1];
 
 int first7;           // the position of the first monster at distance 7 in dcal
 
@@ -5119,7 +5119,7 @@ void moveshadow() {
 void moveghosts() {
 
   if(invismove) return;
-  for(int d=0; d<8; d++) movesofgood[d].clear();  
+  for(int d=0; d<=MAX_EDGE; d++) movesofgood[d].clear();  
 
   for(int i=0; i<isize(ghosts); i++) {
     cell *c = ghosts[i];
@@ -5137,7 +5137,7 @@ void moveghosts() {
       }
     }
   
-  for(int d=0; d<8; d++) for(int i=0; i<isize(movesofgood[d]); i++) {
+  for(int d=0; d<=MAX_EDGE; d++) for(int i=0; i<isize(movesofgood[d]); i++) {
     cell *c = movesofgood[d][i];
     if(c->stuntime) return;
     
@@ -6053,7 +6053,7 @@ void consMove(cell *c, eMonster param) {
 void moveNormals(eMonster param) {
   pathdata pd(param);
 
-  for(int d=0; d<8; d++) movesofgood[d].clear();
+  for(int d=0; d<=MAX_EDGE; d++) movesofgood[d].clear();
 
   for(int i=0; i<isize(pathqm); i++) 
     consMove(pathqm[i], param);
@@ -6064,7 +6064,7 @@ void moveNormals(eMonster param) {
     if(c->pathdist == PINFD) consMove(c, param);
     }
 
-  for(int d=0; d<8; d++) for(int i=0; i<isize(movesofgood[d]); i++) {
+  for(int d=0; d<=MAX_EDGE; d++) for(int i=0; i<isize(movesofgood[d]); i++) {
     cell *c = movesofgood[d][i];
     if(normalMover(c->monst)) {
       moveNormal(c, MF_PATHDIST);

geom-exp.cpp

@@ -383,7 +383,15 @@ vector<eGeometry> quotientlist = {
   };
 
 vector<eGeometry> list3d = {
-  gBinary3, gCubeTiling, gRhombic3, gBitrunc3, gCell120, gECell120
+  gBinary3,
+  gSpace534, gSpace435, 
+  gCubeTiling, gRhombic3, gBitrunc3, 
+  gCell120, gECell120,
+  gCell600, gECell600,
+  gCell24, gECell24,
+  gCell16, gECell16,
+  gCell8, gECell8,
+  gCell5
   };
 
 void ge_select_tiling(const vector<eGeometry>& lst) {

geometry.cpp

@@ -88,7 +88,7 @@ void precalc() {
     goto finish;
     }
   
-  if(sphere && DIM == 3) {
+  if((sphere || hyperbolic) && DIM == 3 && !binarytiling) {
     rhexf = hexf = 0.378077;
     crossf = hcrossf = 0.620672;
     tessf = 1.090550;

geometry2.cpp

@@ -96,7 +96,7 @@ transmatrix calc_relative_matrix(cell *c2, cell *c1, const hyperpoint& point_hin
   #endif
   #if MAXMDIM == 4
   if(euclid && DIM == 3) return euclid3::relative_matrix(c2->master, c1->master);
-  if(sphere && DIM == 3) return sphere3::relative_matrix(c2->master, c1->master);
+  if(DIM == 3) return reg3::relative_matrix(c2->master, c1->master);
   #endif
   #if CAP_ARCM
   if(archimedean) return arcm::relative_matrix(c2->master, c1->master);

graph.cpp

@@ -3842,16 +3842,19 @@ bool isWall3(cell *c, color_t& wcol) {
 
 // how much should be the d-th wall darkened in 3D
 int get_darkval(int d) {
-  const int darkval_h[9] = {0,2,2,0,6,6,8,8,0};
+  const int darkval_hbt[9] = {0,2,2,0,6,6,8,8,0};
-  const int darkval_s[12] = {0,1,2,3,4,5,0,1,2,3,4,5};
+  const int darkval_s12[12] = {0,1,2,3,4,5,0,1,2,3,4,5};
   const int darkval_e6[6] = {0,4,6,0,4,6};
   const int darkval_e12[12] = {0,4,6,0,4,6,0,4,6,0,4,6};
   const int darkval_e14[14] = {0,0,0,4,6,4,6,0,0,0,6,4,6,4};
-  if(sphere) return darkval_s[d];
+  if(sphere) return darkval_s12[d];
   if(euclid && S7 == 6) return darkval_e6[d];
   if(euclid && S7 == 12) return darkval_e12[d];
   if(euclid && S7 == 14) return darkval_e14[d];
-  return darkval_h[d];
+  if(binarytiling) return darkval_hbt[d];
+  if(hyperbolic && S7 == 6) return darkval_e6[d];
+  if(hyperbolic && S7 == 12) return darkval_s12[d];
+  return 0;
   }
 
 void drawcell(cell *c, transmatrix V, int spinv, bool mirrored) {
@@ -4689,7 +4692,7 @@ void drawcell(cell *c, transmatrix V, int spinv, bool mirrored) {
 
           for(int a=0; a<c->type; a++)
             if(c->move(a) && !isWall3(c->move(a), dummy)) {
-              if(a < 4 && hyperbolic) {
+              if(a < 4 && binarytiling) {
                 if(celldistAlt(c) >= celldistAlt(viewctr.at->c7)) continue;
                 queuepoly(V, shWall3D[a], darkena(wcol - d * get_darkval(a), 0, 0xFF));
                 }
@@ -5245,8 +5248,8 @@ void drawcell(cell *c, transmatrix V, int spinv, bool mirrored) {
       else if(DIM == 3) {
         for(int t=0; t<c->type; t++) {
           if(!c->move(t)) continue;
-          if(hyperbolic && !among(t, 5, 6, 8)) continue;
+          if(binarytiling && !among(t, 5, 6, 8)) continue;
-          if(!hyperbolic && c->move(t) < c) continue;
+          if(!binarytiling && c->move(t) < c) continue;
           queuepoly(V, shWall3D[t], 0);
           }
         }
@@ -5872,8 +5875,8 @@ void drawthemap() {
   #if MAXMDIM == 4
   else if(euclid && DIM == 3)
     euclid3::draw();
-  else if(sphere && DIM == 3)
+  else if(DIM == 3)
-    sphere3::draw();
+    reg3::draw();
   #endif
   else
     drawStandard();

heptagon.cpp

@@ -234,6 +234,8 @@ heptagon *createStep(heptagon *h, int d) {
   #if MAXMDIM == 4
   if(!h->move(d) && euclid && DIM == 3) 
     return euclid3::createStep(h, d);
+  if(!h->move(d) && DIM == 3) 
+    return reg3::createStep(h, d);
   #endif
   #if CAP_ARCM
   if(!h->move(d) && archimedean) {

hyper.h

@@ -96,7 +96,7 @@ void addMessage(string s, char spamtype = 0);
 
 #define binarytiling (geometry == gBinaryTiling || geometry == gBinary3)
 #define archimedean (geometry == gArchimedean)
-#define eubinary (euclid || binarytiling || geometry == gCrystal)
+#define eubinary (euclid || binarytiling || geometry == gCrystal || (DIM == 3 && hyperbolic))
 
 #define cgclass (ginf[geometry].cclass)
 #define euclid (cgclass == gcEuclid)
@@ -1518,7 +1518,7 @@ bool bearsCamelot(eLand l);
 extern bool safety;
 
 #define SAGEMELT .1
-#define TEMPLE_EACH ((DIM == 3 && hyperbolic) ? 2 : 6)
+#define TEMPLE_EACH ((DIM == 3 && binarytiling) ? 2 : (DIM == 3 && hyperbolic) ? 3 : 6)
 #define PT(x, y) ((tactic::on || quotient == 2 || daily::on) ? (y) : inv::on ? min(2*(y),x) : (x))
 #define ROCKSNAKELENGTH 50
 #define WORMLENGTH 15
@@ -4257,6 +4257,7 @@ namespace binary {
   transmatrix parabolic(ld u);
   transmatrix parabolic3(ld u, ld v);
   extern ld btrange, btrange_cosh;
+  transmatrix relative_matrix(heptagon *h2, heptagon *h1);
   }
 #endif
 
@@ -4266,6 +4267,18 @@ namespace euclid3 {
   hrmap* new_map();
   void draw();
   }
+
+namespace reg3 {
+  void generate();
+  hrmap* new_map();
+  heptagon *createStep(heptagon *parent, int d);
+  extern vector<hyperpoint> cellshape;
+  void draw();
+  transmatrix relative_matrix(heptagon *h2, heptagon *h1);  
+  int dist_alt(cell *c);
+  int celldistance(cell *c1, cell *c2);
+  bool pseudohept(cell *c);
+  }
 #endif
 
 namespace arcm {

hypgraph.cpp

@@ -1132,8 +1132,8 @@ void optimizeview() {
       if(binarytiling) T = binary::relative_matrix(h2, viewctr.at);
       #endif
       #if MAXMDIM == 4
-      if(euclid && DIM == 3) T = euclid3::relative_matrix(h2, viewctr.at);
+      else if(euclid && DIM == 3) T = euclid3::relative_matrix(h2, viewctr.at);
-      if(sphere && DIM == 3) T = sphere3::relative_matrix(h2, viewctr.at);
+      else if(DIM == 3) T = reg3::relative_matrix(h2, viewctr.at);
       #endif
       #if CAP_ARCM
       if(archimedean)  T = arcm::relative_matrix(h2, viewctr.at);

pattern2.cpp

@@ -450,8 +450,8 @@ int getHemisphere(heptagon *h, int which) {
 int getHemisphere(cell *c, int which) {
   if(euwrap) return 0;
   if(DIM == 3) {
-    ld z = sphere3::vertices120[c->master->zebraval][which];
+    hyperpoint p = tC0(calc_relative_matrix(c, currentmap->gamestart(), C0));
-    return int(z * 6  + 10.5) - 10;
+    return int(p[which] * 6  + 10.5) - 10;
     }
   if(which == 0 && GOLDBERG && has_nice_dual()) {
     set<cell*> visited;
@@ -765,7 +765,7 @@ namespace patterns {
     }
 
   void val_all(cell *c, patterninfo &si, int sub, int pat) {
-    if(IRREGULAR || archimedean || binarytiling) si.symmetries = 1;
+    if(IRREGULAR || archimedean || binarytiling || DIM == 3) si.symmetries = 1;
     else if(a46) val46(c, si, sub, pat);
     else if(a38) val38(c, si, sub, pat);
     else if(sphere && S3 == 3) valSibling(c, si, sub, pat);
@@ -1297,7 +1297,7 @@ bool pseudohept(cell *c) {
   #if MAXMDIM == 4
   if(DIM == 3) {
     if(euclid) return euclid3::pseudohept(c);
-    if(sphere) return sphere3::pseudohept(c);
+    else return reg3::pseudohept(c);
     }
   #endif
   #if CAP_ARCM

patterns.cpp

@@ -363,7 +363,7 @@ int emerald_heptagon(int parent, int dir) {
   
   #undef RULE
 
-  if(weirdhyperbolic) return 0;
+  if(weirdhyperbolic || DIM == 3) return 0;
   printf("HEPTAGONAL RULE MISSING for (%d,%d)\n", parent,dir);
   exit(1);
   }

polygons.cpp

@@ -2414,12 +2414,13 @@ void create_wall3d() {
       }
     }
   
-  if(DIM == 3 && sphere) {
+  if(DIM == 3 && !euclid && !binarytiling) {
-    sphere3::gen600();
+    reg3::generate();
-    for(int w=0; w<12; w++) {
+    int facesize = isize(reg3::cellshape) / S7;
+    for(int w=0; w<S7; w++) {
       bshape(shWall3D[w], PPR::WALL);
-      for(int a=0; a<=5; a++) 
+      for(int a=0; a<=facesize; a++) 
-        hpcpush(sphere3::dodefaces[w*5+a%5]);
+        hpcpush(reg3::cellshape[w*facesize+a%facesize]);
       }
     }
   

reg3.cpp

@@ -0,0 +1,462 @@
+// Hyperbolic Rogue -- regular honeycombs
+// works with spherical and hyperbolic ones -- Euclidean cubic tiling implemented in euclid.cpp
+// hyperbolic honeycombs rely on binary:: to deal with floating point errors (just like archimedean)
+// Copyright (C) 2011-2019 Zeno Rogue, see 'hyper.cpp' for details
+
+namespace hr {
+#if MAXMDIM >= 4
+transmatrix cpush(int cid, ld alpha);
+transmatrix cspin(int a, int b, ld alpha);
+extern
+vector<hpcshape> shWall3D, shMiniWall3D;
+namespace binary {   
+  void build_tmatrix(); 
+  void virtualRebaseSimple(heptagon*& base, transmatrix& at);
+  int celldistance3(heptagon *c1, heptagon *c2);
+  hyperpoint deparabolic3(hyperpoint h);
+  }
+
+namespace reg3 {
+
+  int loop=4, face=5;
+
+  vector<hyperpoint> cellshape;
+  
+  transmatrix spins[12], adjmoves[12];
+
+  template<class T> ld binsearch(ld dmin, ld dmax, const T& f) {
+    for(int i=0; i<200; i++) {
+      ld d = (dmin + dmax) / 2;
+      if(f(d)) dmax = d;
+      else dmin = d;
+      }
+    return dmin;
+    } 
+
+  void generate() {
+  
+    using namespace hyperpoint_vec;
+
+    if(S7 == 4) face = 3;
+    if(S7 == 6) face = 4;
+    if(S7 == 12) face = 5;
+    if(S7 == 8) face = 3;
+    /* icosahedron not implemented */
+    loop = ginf[geometry].tiling_name[5] - '0';
+    println(hlog, "face = ", face, " loop = ", loop, " S7 = ", S7);
+    
+    ld dual_angle = binsearch(0, M_PI, [&] (ld d) {
+      hyperpoint h0 = cpush(0, 1) * C0;
+      hyperpoint h1 = cspin(0, 1, d) * h0;
+      hyperpoint h2 = cspin(1, 2, 2*M_PI/loop) * h1;
+      return hdist(h0, h1) > hdist(h1, h2);
+      });
+
+    ld dodecahedron_angle = binsearch(0, M_PI, [&] (ld d) {
+      hyperpoint h0 = cpush(0, 1) * C0;
+      hyperpoint h1 = cspin(0, 1, d) * h0;
+      hyperpoint h2 = cspin(1, 2, 2*M_PI/face) * h1;
+      return hdist(h0, h1) > hdist(h1, h2);
+      });
+    
+    if(S7 == 8) {
+      /* 24-cell is a special case because it is the only one with '4' in the middle of the Schlaefli symbol. */
+      /* The computations above assume 3 */
+      hyperpoint h1 = hpxy3(.5,.5,.5);
+      hyperpoint h2 = hpxy3(.5,.5,-.5);
+      dual_angle = hdist(h1, h2);
+      }
+    
+    println(hlog, "dodecahedron angle = ", dodecahedron_angle);
+    println(hlog, "dual angle = ", dual_angle);
+    
+    ld inp_length = binsearch(0, 1.55, [&] (ld d) {
+      hyperpoint h = xpush(-d) * spin(2*M_PI/face) * xpush0(d);
+      ld alpha = M_PI - atan2(-h[1], h[0]);
+      return (alpha < dual_angle / 2) ? hyperbolic : sphere;
+      });
+    
+    println(hlog, "inp length = ", inp_length);
+    
+    ld edge_length = hdist(xpush0(inp_length), spin(2*M_PI/face) * xpush0(inp_length));
+    if(S7 == 8) edge_length = hdist(normalize(hpxyz3(1,1,0,0)), normalize(hpxyz3(1,0,1,0))); 
+    println(hlog, "edge length = ", edge_length);
+    
+    hyperpoint h0 = cpush(0, 1) * C0;
+    hyperpoint h1 = cspin(0, 1, dodecahedron_angle) * h0;
+    hyperpoint h2 = cspin(1, 2, 2*M_PI/face) * h1;
+    hyperpoint h3 = cspin(1, 2, -2*M_PI/face) * h1;
+    hyperpoint a2 = S7 == 8 ? normalize(h1 + h2) : normalize(h0 + h1 + h2);
+    hyperpoint a3 = S7 == 8 ? normalize(h1 + h3) : normalize(h0 + h1 + h3);
+    
+    println(hlog, "S7 = ", S7);
+
+    ld whereonline = binsearch(0, 5, [&] (ld d) {
+      // sometimes breaks in elliptic
+      dynamicval<eGeometry> g(geometry, elliptic ? gCell120 : geometry);
+      hyperpoint z2 = a2 * d + C0 * (1-d);
+      if(hyperbolic && intval(z2, Hypc) >= 0) return true;
+      hyperpoint b2 = normalize(z2);
+      hyperpoint z3 = a3 * d + C0 * (1-d);
+      hyperpoint b3 = normalize(z3);
+      return hdist(b2, b3) >= edge_length;
+      });
+    
+    println(hlog, "whereonline = ", whereonline);
+    a2 = normalize(a2 * whereonline + C0 * (1-whereonline));
+    a3 = normalize(a3 * whereonline + C0 * (1-whereonline));
+
+    hyperpoint mid = Hypc;
+    for(int i=0; i<face; i++) mid += cspin(1, 2, 2*i*M_PI/face) * a2;
+    mid = normalize(mid);
+    ld between_centers = 2 * hdist0(mid);
+    println(hlog, "between_centers = ", between_centers);
+
+    if(S7 == 12 || S7 == 8) {
+      spins[0] = Id;
+      spins[1] = cspin(0, 1, dodecahedron_angle) * cspin(1, 2, M_PI);
+      for(int a=2; a<face+1; a++) spins[a] = cspin(1, 2, 2*M_PI*(a-1)/face) * spins[1];      
+      for(int a=S7/2; a<S7; a++) spins[a] = cspin(0, 1, M_PI) * spins[a-S7/2];
+      }
+    
+    if(S7 == 6) {
+      spins[0] = Id;
+      spins[1] = cspin(0, 1, dodecahedron_angle) * cspin(1, 2, M_PI);
+      spins[2] = cspin(1, 2, M_PI/2) * spins[1];
+      for(int a=S7/2; a<S7; a++) spins[a] = spins[a-S7/2] * cspin(0, 1, M_PI);
+      }
+    
+    if(S7 == 4) {
+      spins[0] = Id;
+      spins[1] = cspin(0, 1, dodecahedron_angle) * cspin(1, 2, M_PI);
+      for(int a=2; a<face+1; a++) spins[a] = cspin(1, 2, 2*M_PI*(a-1)/face) * spins[1];      
+      }
+    
+    cellshape.clear();
+    for(int a=0; a<S7; a++)
+    for(int b=0; b<face; b++)
+      cellshape.push_back(spins[a] * cspin(1, 2, 2*M_PI*b/face) * a2);
+    
+    adjmoves[0] = cpush(0, between_centers) * cspin(0, 2, M_PI);
+    for(int i=1; i<S7; i++) adjmoves[i] = spins[i] * adjmoves[0];
+
+    for(int a=0; a<S7; a++)
+      println(hlog, "center of ", a, " is ", tC0(adjmoves[a]));    
+    
+    println(hlog, "doublemove = ", tC0(adjmoves[0] * adjmoves[0]));
+    // exit(1);
+    }
+
+  void binary_rebase(heptagon *h, const transmatrix& V) {
+    }
+  
+  void test();
+    
+  struct hrmap_reg3 : hrmap {
+  
+    heptagon *origin;
+    hrmap *binary_map;
+    
+    int mycellcount = 1;
+
+    unordered_map<heptagon*, pair<heptagon*, transmatrix>> reg_gmatrix;
+    unordered_map<heptagon*, vector<pair<heptagon*, transmatrix> > > altmap;
+
+    hrmap_reg3() {
+      generate();
+      origin = tailored_alloc<heptagon> (S7);
+      heptagon& h = *origin;
+      h.s = hsOrigin;
+      h.cdata = NULL;
+      h.alt = NULL;
+      h.distance = 0;
+      h.c7 = newCell(S7, origin);
+      sightranges[geometry] = 3;
+      
+      dynamicval<hrmap*> cr(currentmap, this);
+      
+      heptagon *alt = NULL;
+      transmatrix T = Id;
+      
+      if(hyperbolic) {
+        dynamicval<eGeometry> g(geometry, gBinary3); 
+        binary::build_tmatrix();
+        alt = tailored_alloc<heptagon> (S7);
+        alt->s = hsOrigin;
+        alt->emeraldval = 0;
+        alt->zebraval = 0;
+        alt->distance = 0;
+        alt->alt = alt;
+        alt->cdata = NULL;
+        binary_map = newAltMap(alt); 
+        T = xpush(.01241) * spin(1.4117) * xpush(0.1241) * cspin(0, 2, 1.1249) * xpush(0.07) * Id;
+        }
+      
+      reg_gmatrix[origin] = make_pair(alt, T);
+      altmap[alt].emplace_back(origin, T);
+      }
+    
+    ld worst_error1 = 0, worst_error2 = 0;
+
+    heptagon *getOrigin() {
+      return origin;
+      }
+    
+    #define DEB 0
+
+    heptagon *createStep(heptagon *parent, int d) {
+      auto& p1 = reg_gmatrix[parent];
+      if(DEB) println(hlog, "creating step ", parent, ":", d, ", at ", p1.first, tC0(p1.second));
+      heptagon *alt = p1.first;
+      transmatrix T = p1.second * adjmoves[d];
+      transmatrix T1 = T;
+      if(hyperbolic) {
+        dynamicval<eGeometry> g(geometry, gBinary3);
+        binary::virtualRebaseSimple(alt, T);
+        }
+      
+      fixmatrix(T);
+      auto hT = tC0(T);
+      
+      if(DEB) println(hlog, "searching at ", alt, ":", hT);
+
+      if(DEB) for(auto& p2: altmap[alt]) println(hlog, "for ", tC0(p2.second), " intval is ", intval(tC0(p2.second), hT));
+      
+      ld err;
+      
+      for(auto& p2: altmap[alt]) if((err = intval(tC0(p2.second), hT)) < 1e-3) {
+        if(err > worst_error1) println(hlog, format("worst_error1 = %lg", double(worst_error1 = err)));
+        // println(hlog, "YES found in ", isize(altmap[alt]));
+        if(DEB) println(hlog, "-> found ", p2.first);
+        int fb = 0;
+        hyperpoint old = T * (inverse(T1) * tC0(p1.second));
+        for(int d2=0; d2<S7; d2++) {
+          hyperpoint back = p2.second * tC0(adjmoves[d2]);
+          if((err = intval(back, old)) < 1e-3) {
+            if(err > worst_error2) println(hlog, format("worst_error2 = %lg", double(worst_error2 = err)));
+            if(p2.first->move(d2)) println(hlog, "error: repeated edge");
+            p2.first->c.connect(d2, parent, d, false);
+            fb++;
+            }
+          }
+        if(fb != 1) { 
+          println(hlog, "found fb = ", fb); 
+          println(hlog, old);
+          for(int d2=0; d2<S7; d2++) {
+            println(hlog, p2.second * tC0(adjmoves[d2]), " in distance ", intval(p2.second * tC0(adjmoves[d2]), old));
+            }
+          parent->c.connect(d, parent, d, false);
+          return parent;
+          }
+        return p2.first;
+        }
+      
+      // println(hlog, "NOT found in ", isize(altmap[alt]), " intval = ", intval(tC0(p1.second),Hypc), " to = ", intval(hT,Hypc), " mcc = ", mycellcount);
+      if(DEB) println(hlog, "-> not found");
+      mycellcount++;
+      heptagon *created = tailored_alloc<heptagon> (S7);
+      created->c7 = newCell(S7, created);
+      created->alt = NULL;
+      created->zebraval = hrand(10);
+      fixmatrix(T);
+      reg_gmatrix[created] = make_pair(alt, T);
+      altmap[alt].emplace_back(created, T);
+      created->c.connect(0, parent, d, false);
+      return created;
+      }
+
+    };
+  
+hrmap* new_map() {
+  return new hrmap_reg3;
+  }
+
+hrmap_reg3* regmap() {
+  return ((hrmap_reg3*) currentmap);
+  }
+
+heptagon *createStep(heptagon *parent, int d) {
+  return regmap()->createStep(parent, d);
+  }
+
+transmatrix relative_matrix(heptagon *h2, heptagon *h1) {
+  auto m = regmap();
+  auto p1 = m->reg_gmatrix[h1];
+  auto p2 = m->reg_gmatrix[h2];
+  transmatrix T = Id;
+  if(hyperbolic) { 
+    dynamicval<eGeometry> g(geometry, gBinary3);
+    T = binary::relative_matrix(p2.first, p1.first);
+    }
+  return inverse(p1.second) * T * p2.second;
+  }
+
+void draw() {
+  sphereflip = Id;
+  
+  // for(int i=0; i<S6; i++) queuepoly(ggmatrix(cwt.at), shWall3D[i], 0xFF0000FF);
+  
+  dq::visited.clear();
+  dq::enqueue(viewctr.at, cview());
+  
+  while(!dq::drawqueue.empty()) {
+    auto& p = dq::drawqueue.front();
+    heptagon *h = get<0>(p);
+    transmatrix V = get<1>(p);
+    dynamicval<ld> b(band_shift, get<2>(p));
+    bandfixer bf(V);
+    dq::drawqueue.pop();
+    
+    
+    cell *c = h->c7;
+    if(!do_draw(c, V)) continue;
+    drawcell(c, V, 0, false);
+
+    for(int i=0; i<S7; i++)
+      dq::enqueue(h->move(i), V * relative_matrix(h->move(i), h));
+    }
+  }
+
+int celldistance(cell *c1, cell *c2) {
+  if(c1 == c2) return 0;
+  auto r = regmap();
+  dynamicval<eGeometry> g(geometry, gBinary3);  
+  return 1 + binary::celldistance3(r->reg_gmatrix[c1->master].first, r->reg_gmatrix[c2->master].first);
+  }
+
+bool pseudohept(cell *c) {
+  if(sphere) {
+    hyperpoint h = tC0(relative_matrix(c->master, regmap()->origin));
+    if(S7 == 12) {
+      hyperpoint h1 = cspin(0, 1, atan2(16, 69) + M_PI/4) * h;
+      for(int i=0; i<4; i++) if(abs(abs(h1[i]) - .5) > .01) return false;
+      return true;
+      }
+    if(S7 == 8)
+      return h[3] >= .99 || h[3] <= -.99 || abs(h[3]) < .01; 
+    if(loop == 3 && face == 3 && S7 == 4)
+      return c == currentmap->gamestart();
+    if(loop == 4 && face == 3)
+      return abs(h[3]) > .9;
+    if(loop == 3 && face == 4)
+      return abs(h[3]) > .9;
+    if(loop == 5 && face == 3)
+      return abs(h[3]) > .99 || abs(h[0]) > .99 || abs(h[1]) > .99 || abs(h[2]) > .99;
+    }
+  if(hyperbolic) {
+    heptagon *h = regmap()->reg_gmatrix[c->master].first;
+    return (h->zebraval == 1) && (h->distance & 1);
+    }    
+  return false;
+  }
+
+int dist_alt(cell *c) {
+  return regmap()->reg_gmatrix[c->master].first->distance;
+  }
+#endif
+
+#if 0
+/* More precise, but very slow distance. Not used/optimized for now */
+
+ld adistance(cell *c) {  
+  hyperpoint h = tC0(regmap()->reg_gmatrix[c->master].second);
+  h = binary::deparabolic3(h);
+  return regmap()->reg_gmatrix[c->master].first->distance * log(2) - h[0];
+  }
+
+int bucketer(ld x) {
+  return int(x * 10 + 100000.5) - 100000;
+  }
+
+int bucketer(hyperpoint h) {
+  return bucketer(h[0]) + 1000 * bucketer(h[1]) + 1000000 * bucketer(h[2]);
+  }
+
+map<int, int> close_distances;
+
+unordered_map<pair<cell*, cell*>, int> memo;
+
+bool cdd;
+
+int celldistance(cell *c1, cell *c2) {
+  if(memo.count(make_pair(c1, c2))) return memo[make_pair(c1, c2)];
+  if(c1 == c2) return 0;
+  vector<cell*> v[2];
+  v[0].push_back(c1);
+  v[1].push_back(c2);
+  
+  int steps = 0;
+  
+  map<cell*, int> visited;
+  visited[c1] = 1;
+  visited[c2] = 2;
+  
+  while(true) {
+    if(cdd) {
+      println(hlog, "state ", steps, "/",isize(v[0]), "/", isize(v[1]));
+      println(hlog, "  A: ", v[0]);
+      println(hlog, "  B: ", v[1]); 
+      }
+    for(int i: {0,1}) {
+      vector<cell*> new_v;
+      for(cell *c: v[i]) forCellCM(cn, c) if(adistance(cn) < adistance(c)) {
+        auto &vi = visited[cn];
+        if((vi&3) == 0) {
+          vi = 4 * (steps+1);
+          vi |= (1<<i);
+          new_v.push_back(cn);
+          }
+        else if((vi&3) == 2-i) {
+          vector<pair<cell*, int>> ca1, ca2;
+          int b1 = 4*steps-4;
+          int b2 = ((vi>>2)<<2) - 4;
+          for(auto p: visited) {
+            if(cdd) println(hlog, p);
+            int ps = p.second & 3;
+            if(ps == 1+i && p.second >= b1)
+              ca1.emplace_back(p.first, p.second/4);
+            if(ps == 2-i && p.second >= b2 && p.second <= b2+8)
+              ca2.emplace_back(p.first, p.second/4);
+            }
+          int bound = 1<<16;
+          for(auto p1: ca1) for(auto p2: ca2) {
+            hyperpoint h = tC0(relative_matrix(p1.first->master, p2.first->master));
+            int b = bucketer(h);
+            if(close_distances.count(b)) {
+              int d = close_distances[b] + p1.second + p2.second;
+              if(cdd) println(hlog, "candidate: close=", close_distances[b], p1, p2, "; h = ", h);
+              if(d < bound) bound = d;
+              }
+            else if(cdd) println(hlog, "bucket missing");
+            }
+          return memo[make_pair(c1, c2)] = bound;
+          return bound;
+          }
+        }
+      v[i] = std::move(new_v);
+      }
+    steps++;
+    }
+  }
+
+cellwalker target;
+int tsteps;
+
+int dist_alt(cell *c) {
+  if(!target.at) {
+    target = cellwalker(currentmap->gamestart(), 0);
+    tsteps = 0;
+    for(int i=0; i<30; i++) target += wstep, target += rev, tsteps++;
+    }
+  if(specialland == laCamelot) return reg3::celldistance(c, target.at);
+  else {
+    int d = reg3::celldistance(c, target.at) - tsteps;
+    if(d < 10) target += wstep, target += rev, tsteps++;
+    return d;
+    }
+  }
+#endif
+
+}
+}
+

sphere.cpp

@@ -164,271 +164,4 @@ heptagon *getDodecahedron(int i) {
   return s->dodecahedron[i];
   }
 
-namespace sphere3 {
-
-vector<hyperpoint> vertices120;
-array<transmatrix, 120> vmatrix120;
-vector<int> adj0;
-array<array<int, 4>, 120> js;
-array<hyperpoint, 60> dodefaces;
-
-int opposite[120];
-
-hyperpoint zero4;
-
-int root;
-
-ld norm(hyperpoint a, hyperpoint b) {
-  ld res = 0;
-  for(int i=0; i<4; i++) res += pow(a[i]-b[i], 2);
-  return res;
-  }
-
-void gen600() {
-  dynamicval<eGeometry> gp(geometry, gCell120);
-  vertices120.clear();
-  root = 23;
-  
-  /// coordinates taken from Wikipedia
-  
-  for(int a=0; a<16; a++) {
-    hyperpoint v = zero4;
-    for(int i=0; i<4; i++) v[i] = ((a >> i) & 1) ? .5 : -.5;
-    vertices120.push_back(v);
-    }
-  
-  for(int i=0; i<4; i++) for(int q: {-1, 1}) {
-    hyperpoint v = zero4;
-    v[i]=q;
-    vertices120.push_back(v);
-    }
-  
-  ld phi = (1 + sqrt(5)) / 2;
-  
-  array<ld, 4> coo = make_array<ld>(1, phi, 1/phi, 0);
-  
-  // all permutations
-  array<int, 4> tab;
-  for(int i=0; i<4; i++) tab[i] = i;
-  
-  do {
-
-    // check the permutation's sign
-    auto tabs = tab;
-    int inv = 0;
-
-    for(int i=0; i<4; i++) while(tabs[i] != i) {
-      swap(tabs[i], tabs[tabs[i]]);
-      inv++;
-      }
-    
-    if(inv&1) goto again;
-    
-    // 8 vertices for each permutation
-    
-    for(int sg=0; sg<8; sg++) {
-      hyperpoint v;
-      for(int i=0; i<4; i++) 
-        v[i] = (((sg >> tab[i])&1) ? 1 : -1) * coo[tab[i]]/2;
-      vertices120.push_back(v);
-      }
-      
-    again: ;
-    }
-  while(std::next_permutation(tab.begin(), tab.end()));
-  
-  if(isize(vertices120) != 120) {
-    printf("error: wrong number of vertices\n");
-    exit(1);
-    }
-  // we add edges between vertices which are close to each other
-  // ((specifically in distance 1/phi/phi)
-  
-  bool inedge[120][120];
-    
-  for(int i=0; i<120; i++)
-  for(int j=0; j<120; j++) {
-    ld d = hdist(vertices120[i], vertices120[j]);
-    inedge[i][j] = (i != j) && d < sqrt(.4);
-    }
-  
-  vector<hyperpoint> cellvertices;
-  
-  for(int i=0; i<120; i++)
-    for(int j=0; j<120; j++) if(inedge[i][j])
-    for(int k=0; k<120; k++) if(inedge[i][k] && inedge[k][j])
-    for(int l=0; l<120; l++) if(inedge[i][l] && inedge[j][l] && inedge[k][l]) {
-      array<int, 4> ijkl = make_array(i, j, k, l);
-      transmatrix T;
-      for(int z=0; z<4; z++) set_column(T, z, vertices120[ijkl[z]]);
-      if(det(T) > 0) js[i] = ijkl;
-      }
-    
-  /* transmatrix src;
-  for(int z=0; z<4; z++) set_column(src, z, vertices120[js[0][z]]); */
-
-  for(int i=0; i<120; i++)
-    for(int z=0; z<4; z++) set_column(vmatrix120[i], z, vertices120[js[i][z]]);
-
-  for(int i=0; i<120; i++) {
-    println(hlog, i, ": ", js[i], " -> ", vmatrix120[i]);  
-    println(hlog, vmatrix120[i] * hyperpoint(1,0,0,0), " should be ", vertices120[i]);  
-    }
-  
-  adj0.clear();
-  for(int i=0; i<120; i++) if(inedge[root][i]) adj0.push_back(i);
-  
-  using namespace hyperpoint_vec;
-  for(int i=0; i<6; i++) for(int j=i+1; j<12; j++) if(zero_d(3, vertices120[adj0[i]] + vertices120[adj0[j]]))
-    swap(adj0[j], adj0[i+6]);
-
-  for(int i=0; i<120; i++) for(int j=0; j<120; j++)
-    if(hdist(vertices120[i], vertices120[j]) > 3)
-      opposite[i] = j;
-  
-  using namespace hyperpoint_vec;
-  
-  int id = 0;
-  for(int i=0; i<12; i++) {
-    int ot = adj0[i];
-    
-    vector<int> pentagon;
-    for(int j: adj0) if(inedge[ot][j]) pentagon.push_back(j);
-    println(hlog, i, ": ", pentagon);
-    
-    int illegal = -1;
-    int at = pentagon[0];
-    for(int d=0; d<5; d++) {
-      for(int s: pentagon) if(inedge[at][s] && s != illegal) {
-        hyperpoint m = vertices120[root] + vertices120[ot] + vertices120[at] + vertices120[s];
-        m = mid(m, m);
-        println(hlog, id, ": ", m);
-        dodefaces[id++] = m;
-        illegal = at;
-        at = s;
-        break;
-        }
-      }
-    }
-  
-  printf("id = %d\n", id);
-  }
-
-bool goodside(int i) {
-  if(!elliptic) return true;
-  hyperpoint& h = vertices120[i];
-  for(int k=3; k>=0; k--) {
-    if(h[k] > 1e-3) return true;
-    if(h[k] < -1e-3) return false;
-    }
-  return false;
-  }
-
-struct hrmap_spherical3 : hrmap {
-  heptagon* cells[120];
-
-  hrmap_spherical3() {
-    gen600();
-    for(int i=0; i<120; i++) {
-      if(!goodside(i)) {
-        cells[i] = NULL;
-        continue;
-        }
-      cells[i] = tailored_alloc<heptagon> (12);
-      heptagon& h = *(cells[i]);
-      h.s = hsOrigin;
-      h.emeraldval = i;
-      h.zebraval = i;
-      h.fiftyval = i;
-      h.rval0 = h.rval1 = 0;
-      h.alt = NULL;
-      h.cdata = NULL;
-      h.c.fullclear();
-      h.fieldval = i;
-      h.c7 = newCell(12, &h);      
-      }
-
-    for(int i=0; i<120; i++) if(cells[i]) {
-      for(int k=0; k<12; k++) {
-        hyperpoint which = vmatrix120[i] * inverse(vmatrix120[root]) * vertices120[adj0[k]];
-        for(int s=0; s<120; s++) if(hdist(which, vertices120[s]) < 1e-6) {
-          int s1 = s;
-          if(!cells[s1]) continue;
-          cells[i]->move(k) = cells[s1];
-          println(hlog, i,".",k, " -> ", s1, " ; ", js[i], " distance = ", hdist(vertices120[i], vertices120[s]));
-          }
-        }
-      }
-
-    for(int i=0; i<120; i++)
-      for(int k=0; k<12; k++) 
-        for(int l=0; l<12; l++)
-          if(cells[i] && cells[i]->move(k)->move(l) == cells[i])
-            cells[i]->c.setspin(k, l, false);
-    }
-
-  heptagon *getOrigin() { return cells[root]; }
-
-  ~hrmap_spherical3() {
-    for(int i=0; i<120; i++) if(cells[i]) tailored_delete(cells[i]);
-    }    
-  };
-
-transmatrix gmatr(heptagon *h) {
-  return vmatrix120[h->zebraval] * inverse(vmatrix120[root]);
-  }
-
-transmatrix relative_matrix(heptagon *h2, heptagon *h1) {
-  return inverse(gmatr(h1)) * gmatr(h2);
-  }
-
-bool pseudohept(cell *c) {
-  return c->master->zebraval < 16;
-  }
-
-void draw() {
-  sphereflip = Id;
-  auto m = (hrmap_spherical3*) currentmap;
-
-  for(int i=0; i<120; i++) if(m->cells[i])
-    drawcell(m->cells[i]->c7, View * relative_matrix(m->cells[i], viewctr.at), 0, false);
-  }
-
-void makewax(int x) {
-  int waxcenter = 63;
-  auto m = (hrmap_spherical3*) currentmap;
-  for(int i=0; i<120; i++) if(m->cells[i]) m->cells[i]->c7->wall = waNone;
-  m->cells[waxcenter]->c7->wall = waDune;
-  int cols[16] = {0x202020, 0x2020A0, 0x20A020, 0x20A0A0, 0xA02020, 0xA020A0, 0xA0A020, 0xA0A0A0,
-    0x606060, 0x6060FF, 0x60FF60, 0x60FFFF, 0xFF6060, 0xFF60FF, 0xFFFF60, 0xFFFFFF };
-  if(x) for(int i=0; i<12; i++) {
-    m->cells[waxcenter]->c7->move(i)->wall = waWaxWall;
-    m->cells[waxcenter]->c7->move(i)->landparam = cols[i];
-    }
-  }
-
-#if CAP_COMMANDLINE  
-int readArgs() {
-  using namespace arg;
-           
-  if(argis("-wax1")) {
-    PHASE(3);
-    start_game();
-    makewax(1);
-    }
-  else if(argis("-wax0")) {
-    PHASE(3);
-    start_game();
-    makewax(0);
-    }
-  else return 1;
-  return 0;
-  }
-
-auto hook = 
-  addHook(hooks_args, 100, readArgs);
-#endif
-
-}
-
 }