various fixes related to product spaces, especially product+inverse

cell.cpp

@@ -292,10 +292,10 @@ EX void initcells() {
   
   hrmap* res = callhandlers((hrmap*)nullptr, hooks_newmap);
   if(res) currentmap = res;
+  else if(nonisotropic || hybri) currentmap = nisot::new_map();
   else if(INVERSE) currentmap = gp::new_inverse();
   else if(fake::in()) currentmap = fake::new_map();
   else if(asonov::in()) currentmap = asonov::new_map();
-  else if(nonisotropic || hybri) currentmap = nisot::new_map();
   #if CAP_CRYSTAL
   else if(cryst) currentmap = crystal::new_map();
   #endif
@@ -938,11 +938,11 @@ EX cdata *arcmCdata(cell *c) {
 
 EX int getCdata(cell *c, int j) {
   if(fake::in()) return FPIU(getCdata(c, j));
-  if(INVERSE) {
+  if(hybri) { c = hybrid::get_where(c).first; return PIU(getBits(c)); }
+  else if(INVERSE) {
     cell *c1 = gp::get_mapped(c);
     return UIU(getCdata(c1, j));
     }
-  if(hybri) { c = hybrid::get_where(c).first; return PIU(getBits(c)); }
   else if(euc::in()) return getEuclidCdata(euc2_coordinates(c))->val[j];
   else if(arcm::in() && euclid)
     return getEuclidCdata(pseudocoords(c))->val[j];
@@ -961,11 +961,11 @@ EX int getCdata(cell *c, int j) {
 
 EX int getBits(cell *c) {
   if(fake::in()) return FPIU(getBits(c));
-  if(INVERSE) {
+  if(hybri) { c = hybrid::get_where(c).first; return PIU(getBits(c)); }
+  else if(INVERSE) {
     cell *c1 = gp::get_mapped(c);
     return UIU(getBits(c1));
     }
-  if(hybri) { c = hybrid::get_where(c).first; return PIU(getBits(c)); }
   else if(euc::in()) return getEuclidCdata(euc2_coordinates(c))->bits;
   else if(arcm::in() && euclid)
     return getEuclidCdata(pseudocoords(c))->bits;

geometry.cpp

@@ -679,6 +679,7 @@ void geometry_information::prepare_basics() {
       single_step = S3 * S7 - 2 * S7 - 2 * S3;
       psl_steps = 2 * S7;    
       if(BITRUNCATED) psl_steps *= S3;
+      if(inv) psl_steps = 2 * S3;
       if(single_step < 0) single_step = -single_step;
       }
     DEBB(DF_GEOM | DF_POLY, ("steps = ", psl_steps, " / ", single_step));

geometry2.cpp

@@ -594,6 +594,7 @@ EX bool approx_nearcorner = false;
 
 EX hyperpoint nearcorner(cell *c, int i) {
   if(GOLDBERG_INV) {
+    i = gmod(i, c->type);
     cellwalker cw(c, i);
     cw += wstep;
     transmatrix cwm = currentmap->adj(c, i);

goldberg.cpp

@@ -1126,6 +1126,8 @@ EX namespace gp {
       }
 
     transmatrix relative_matrix(cell *c2, cell *c1, const hyperpoint& hint) override {
+      c1 = mapping[c1];
+      c2 = mapping[c2];
       return in_underlying([&] { return currentmap->relative_matrix(c2, c1, hint); });
       }
 

nonisotropic.cpp

@@ -1231,6 +1231,8 @@ EX namespace hybrid {
       dynamicval<hrmap*> gpm(pmap, this);
       dynamicval<eGeometry> gag(actual_geometry, geometry);
       dynamicval<eGeometry> g(geometry, underlying);
+      dynamicval<int> gss(underlying_cgip->single_step, cgi.single_step);
+      dynamicval<int> gsp(underlying_cgip->psl_steps, cgi.psl_steps);
       dynamicval<geometry_information*> gc(cgip, underlying_cgip);
       dynamicval<hrmap*> gu(currentmap, underlying_map);
       return t();
@@ -1308,6 +1310,8 @@ EX namespace hybrid {
     if(!hybri) return f();
     dynamicval<eGeometry> g(geometry, underlying);
     dynamicval<eGeometry> gag(actual_geometry, geometry);
+    dynamicval<int> gss(underlying_cgip->single_step, cgi.single_step);
+    dynamicval<int> gsp(underlying_cgip->psl_steps, cgi.psl_steps);
     dynamicval<geometry_information*> gc(cgip, underlying_cgip);
     dynamicval<hrmap*> gpm(pmap, currentmap);
     dynamicval<hrmap*> gm(currentmap, get_umap());
@@ -1345,7 +1349,7 @@ EX namespace hybrid {
     }
   
   EX int wall_offset(cell *c) {
-    if(GOLDBERG) {
+    if(GOLDBERG || INVERSE) {
       /* a bit slow... */
       cell *c1 = WDIM == 2 ? c : get_where(c).first;
       gp::draw_li = WDIM == 2 ? gp::get_local_info(c1) : PIU(gp::get_local_info(c1));

system.cpp

@@ -1318,7 +1318,7 @@ EX void set_geometry(eGeometry target) {
     if(bt::in() || WDIM == 3 || kite::in() || arb::in()) if(!hybri) variation = eVariation::pure;
     #endif
     if(S3 >= OINF) variation = eVariation::pure;
-    if(INVERSE) variation = gp::variation_for(gp::param);
+    if(INVERSE && !hybri) variation = gp::variation_for(gp::param);
     if(ginf[target].default_variation == eVariation::pure && geometry != gArchimedean)
       variation = eVariation::pure;
     if(was_default) pmodel = default_model();