From cdb2a001ca262686ed9c8f9871893b6f41d11974 Mon Sep 17 00:00:00 2001
From: Zeno Rogue <zeno@attnam.com>
Date: Mon, 11 Jul 2022 17:32:27 +0200
Subject: [PATCH] rulegen3

---
 hyper.cpp    |   1 +
 reg3.cpp     |  14 +++++++
 rulegen.cpp  |  47 +++++++++++++++++----
 rulegen3.cpp | 113 +++++++++++++++++++++++++++++++++++++++++++++++++++
 4 files changed, 168 insertions(+), 7 deletions(-)
 create mode 100644 rulegen3.cpp

diff --git a/hyper.cpp b/hyper.cpp
index d0c988ec..e46b5a3f 100644
--- a/hyper.cpp
+++ b/hyper.cpp
@@ -40,6 +40,7 @@
 #include "archimedean.cpp"
 #include "arbitrile.cpp"
 #include "rulegen.cpp"
+#include "rulegen3.cpp"
 #include "euclid.cpp"
 #include "sphere.cpp"
 #include "fake.cpp"
diff --git a/reg3.cpp b/reg3.cpp
index 49d8bb07..c21565ce 100644
--- a/reg3.cpp
+++ b/reg3.cpp
@@ -1732,6 +1732,16 @@ EX namespace reg3 {
       }
     };
 
+  EX int get_aid(cell *c) {
+    auto m = dynamic_cast<hrmap_h3*> (currentmap);
+    if(!m) throw hr_exception("get_aid incorrect");
+    return m->cell_id[c];
+    }
+
+  EX int get_size_of_aid(int aid) {
+    auto m = dynamic_cast<hrmap_h3*> (currentmap);
+    if(!m) throw hr_exception("get_size_of_fv incorrect");
+    return m->quotient_map->acells[aid]->type;
     };
 
   struct hrmap_sphere3 : hrmap_closed3 {
@@ -2206,6 +2216,10 @@ EX int quotient_count() {
   return isize(hypmap()->quotient_map->allh);
   }
 
+EX int quotient_count_sub() {
+  return isize(hypmap()->quotient_map->acells);
+  }
+
 /** This is a generalization of hyperbolic_celldistance in expansion.cpp to three dimensions.
     It still assumes that there are at most 4 cells around every edge, and that distances from
     the origin are known, so it works only in {5,3,4}.
diff --git a/rulegen.cpp b/rulegen.cpp
index ab964d8c..d464b0c1 100644
--- a/rulegen.cpp
+++ b/rulegen.cpp
@@ -184,21 +184,25 @@ twalker addstep(twalker x) {
 
 int number_of_types() {
   if(arb::in()) return isize(arb::current.shapes);
+  if(WDIM == 3) return reg3::quotient_count_sub();
   throw hr_exception("unknown number_of_types");
   }
 
 int get_id(cell *c) {
   if(arb::in()) return shvid(c);
+  if(GDIM == 3) return reg3::get_aid(c);
   throw hr_exception("unknown get_id");
   }
 
 int shape_size(int id) {
   if(arb::in()) return isize(arb::current.shapes[id].connections);
+  if(GDIM == 3) return reg3::get_size_of_aid(id);
   throw hr_exception("unknown shape_size");
   }
 
 int cycle_size(int id) {
   if(arb::in()) return arb::current.shapes[id].cycle_length;
+  if(GDIM == 3) return reg3::get_size_of_aid(id);
   throw hr_exception("unknown shape size");
   }
 
@@ -222,8 +226,8 @@ tcell *gen_tcell(int id) {
   return c;
   }
 
-map<cell*, tcell*> cell_to_tcell;
-map<tcell*, cell*> tcell_to_cell;
+EX map<cell*, tcell*> cell_to_tcell;
+EX map<tcell*, cell*> tcell_to_cell;
 
 void numerical_fix(twalker pw) {
   auto& shs = arb::current.shapes;
@@ -335,6 +339,7 @@ tcell* tmove(tcell *c, int d) {
 
 /** check whether we have completed the vertex to the right of edge d of c */
 void check_loops(twalker pw) {
+  if(GDIM == 3) throw hr_exception("check_loops called");
   ufind(pw);
   auto& shs = arb::current.shapes;
   int id = pw.at->id;
@@ -391,6 +396,7 @@ EX void unify(twalker pw1, twalker pw2) {
   ufind(pw1);
   ufind(pw2);
   if(pw1 == pw2) return;
+  if(GDIM == 3) throw hr_exception("check_loops called");
   callhooks(hooks_gen_tcell, 3, pw1);
   callhooks(hooks_gen_tcell, 4, pw2);
   if(pw1.at->unified_to.at != pw1.at)
@@ -1016,8 +1022,8 @@ using aid_t = pair<int, int>;
 struct analyzer_state {
   int analyzer_id;
   int id, dir;
-  array<analyzer_state*, 10> substates;
-  analyzer_state() { id = MYSTERY; dir = MYSTERY; for(int i=0; i<10; i++) substates[i] = nullptr; }
+  map<int, analyzer_state*> substates;
+  analyzer_state() { id = MYSTERY; dir = MYSTERY; } // for(int i=0; i<10; i++) substates[i] = nullptr; }
   vector<twalker> inhabitants;
   };
 
@@ -1169,6 +1175,7 @@ int get_sidecache(twalker what) {
   }
 
 int get_side(twalker what) {
+  if(WDIM == 3) throw hr_exception("called get_side");
 
   bool side = !(flags & w_no_sidecache);
   bool fast = (flags & w_slow_side);
@@ -1287,6 +1294,11 @@ EX int move_code(twalker cs) {
 
      int x;
 
+     if(WDIM == 3) {
+       if(cs2.at->parent_dir == cs2.spin) return C_PARENT;
+       else return get_roadsign(cs);
+       }
+
      if(!(flags & w_no_relative_distance)) x = C_EQUAL;
      else if(y == 1) x = C_NEPHEW;
      else if(y == 0) x = C_EQUAL;
@@ -1313,7 +1325,15 @@ EX void id_at_spin(twalker cw, vector<twalker>& sprawl, vector<analyzer_state*>&
       a = alloc_analyzer();
       }
     states.push_back(a);
-    if(isize(sprawl) <= cw.at->type) {
+    if(WDIM == 3) {
+      auto& ae = check_all_edges(cw, a, isize(sprawl));
+      int id = isize(sprawl);
+      if(id < isize(ae)) {
+        a->id = ae[id].first;
+        a->dir = ae[id].second;
+        }
+      }
+    else if(isize(sprawl) <= cw.at->type) {
       a->id = 0, a->dir = isize(sprawl)-1;
       // println(hlog, "need to go in direction ", a->dir);
       }
@@ -1411,9 +1431,12 @@ vector<int> gen_rule(twalker cwmain, int id) {
     cids.push_back(id1);
     }
 
-  for(int i=0; i<isize(cids); i++) if(cids[i] == DIR_UNKNOWN)
+  if(WDIM != 3) for(int i=0; i<isize(cids); i++) if(cids[i] == DIR_UNKNOWN)
     cids[i] = get_side(cwmain+i) < 0 ? DIR_RIGHT : DIR_LEFT;
 
+  if(WDIM == 3) for(int i=0; i<isize(cids); i++) if(cids[i] == DIR_UNKNOWN)
+    cids[i] = get_roadsign(cwmain+i);
+
   return cids;
   }
 
@@ -1662,6 +1685,7 @@ void verified_treewalk(twalker& tw, int id, int dir) {
 EX bool view_examine_branch = false;
 
 bool examine_branch(int id, int left, int right) {
+  if(WDIM == 3) return true;
   auto rg = treestates[id].giver;
 
   if(debugflags & DF_GEOM)
@@ -2012,6 +2036,7 @@ EX void cleanup() {
   important.clear();
   shortcuts.clear();
   single_live_branch_close_to_root.clear();
+  cleanup3();
   }
 
 EX void clear_all() {  
@@ -2033,7 +2058,15 @@ EX void generate_rules() {
   start_time = SDL_GetTicks();
   delete_tmap();
 
-  if(!arb::in()) try {
+  if(WDIM == 3) {
+    stop_game();
+    reg3::reg3_rule_available = false;
+    fieldpattern::use_rule_fp = true;
+    fieldpattern::use_quotient_fp = true;
+    flags |= w_numerical;
+    start_game();
+    }
+  else if(!arb::in()) try {
     arb::convert::convert();
     if(flags & w_numerical) arb::convert::activate();
     }
diff --git a/rulegen3.cpp b/rulegen3.cpp
new file mode 100644
index 00000000..b164c038
--- /dev/null
+++ b/rulegen3.cpp
@@ -0,0 +1,113 @@
+// Hyperbolic Rogue -- rule generator
+// Copyright (C) 2011-2021 Zeno Rogue, see 'hyper.cpp' for details
+
+/** \file rulegen3.cpp 
+ *  \brief An algorithm to create strict tree rules for arb tessellations -- 3D parts
+ */
+
+#include "hyper.h"
+
+namespace hr {
+
+EX namespace rulegen {
+
+/** next roadsign ID -- they start at -100 and go downwards */
+int next_roadsign_id = -100;
+
+/** get the ID of a roadsign path */
+EX map<vector<int>, int> roadsign_id;
+
+EX int get_roadsign(twalker what) {
+  int dlimit = what.at->dist - 1;
+  tcell *s = what.at, *t = what.peek();
+  vector<int> result;
+  while(s->dist > dlimit) { 
+    twalker s0 = s;
+    get_parent_dir(s0);
+    if(s->parent_dir == MYSTERY) throw hr_exception("parent_dir unknown");
+    result.push_back(s->parent_dir); s = s->move(s->parent_dir);
+    result.push_back(s->dist - dlimit);
+    }
+  vector<int> tail;
+  while(t->dist > dlimit) {
+    twalker t0 = t;
+    get_parent_dir(t0);
+    if(t->parent_dir == MYSTERY) throw hr_exception("parent_dir unknown");
+    tail.push_back(t->dist - dlimit);
+    tail.push_back(t->c.spin(t->parent_dir));
+    t = t->move(t->parent_dir);
+    }
+  map<tcell*, int> visited;
+  queue<tcell*> vqueue;
+  auto visit = [&] (tcell *c, int dir) {
+    if(visited.count(c)) return;
+    visited[c] = dir;
+    vqueue.push(c);
+    };
+  visit(s, MYSTERY);
+  while(true) {
+    if(vqueue.empty()) throw hr_exception("vqueue empty");
+    tcell *c = vqueue.front();
+    if(c == t) break;
+    vqueue.pop();
+    for(int i=0; i<c->type; i++)
+      if(c->move(i) && c->move(i)->dist <= dlimit)
+        visit(c->move(i), c->c.spin(i));
+    }
+  while(t != s) {
+    int d = visited.at(t);
+    tail.push_back(t->dist - dlimit);
+    tail.push_back(t->c.spin(d));
+    t = t->move(d);
+    }
+  reverse(tail.begin(), tail.end());
+  for(auto t: tail) result.push_back(t);
+  if(roadsign_id.count(result)) return roadsign_id[result];
+  return roadsign_id[result] = next_roadsign_id--;
+  }
+
+map<pair<int, int>, vector<pair<int, int>> > all_edges;
+
+EX vector<pair<int, int>>& check_all_edges(twalker cw, analyzer_state* a, int id) {
+  auto& ae = all_edges[{cw.at->id, cw.spin}];
+  if(ae.empty()) {
+    set<tcell*> seen;
+    vector<pair<twalker, transmatrix> > visited;
+    auto visit = [&] (twalker tw, const transmatrix& T, int id, int dir) {
+      if(seen.count(tw.at)) return;
+      seen.insert(tw.at);
+      auto& sh0 = currentmap->get_cellshape(tcell_to_cell[cw.at]);
+      auto& sh1 = currentmap->get_cellshape(tcell_to_cell[tw.at]);
+      int common = 0;
+      vector<hyperpoint> rotated;
+      for(auto w: sh1.vertices_only_local) rotated.push_back(T*w);
+      
+      for(auto v: sh0.vertices_only_local)      
+      for(auto w: rotated)
+        if(sqhypot_d(MDIM, v-w) < 1e-6)
+          common++;
+      if(common < 2) return;
+      visited.emplace_back(tw, T);
+      ae.emplace_back(id, dir);
+      };
+    visit(cw, Id, -1, -1);
+    for(int i=0; i<isize(visited); i++) {
+      auto tw = visited[i].first;      
+      for(int j=0; j<tw.at->type; j++) {
+        visit(tw + j + wstep, visited[i].second * currentmap->adj(tcell_to_cell[tw.at], (tw+j).spin), i, j);
+        }
+      }
+    println(hlog, "for ", tie(cw.at->id, cw.spin), " generated all_edges structure: ", ae, " of size ", isize(ae));
+    }
+  return ae;
+  }
+
+EX void cleanup3() {
+  all_edges.clear();
+  roadsign_id.clear();
+  next_roadsign_id = -100;
+  }
+
+}
+
+}
\ No newline at end of file