From e26f8f5a5beff5d3b64c516c37331a88cdb0286a Mon Sep 17 00:00:00 2001
From: Zeno Rogue <zeno@attnam.com>
Date: Wed, 16 Sep 2020 05:57:05 +0200
Subject: [PATCH] optimized inverses

---
 archimedean.cpp   |  4 +--
 binary-tiling.cpp |  4 +--
 control.cpp       |  4 +--
 drawing.cpp       |  4 +--
 geometry.cpp      |  4 +--
 geometry2.cpp     |  2 +-
 graph.cpp         | 32 ++++++++++----------
 hyperpoint.cpp    | 75 +++++++++++++++++++++++++++++++++++++++--------
 hypgraph.cpp      | 32 ++++++++++----------
 mapeditor.cpp     |  6 ++--
 models.cpp        |  2 +-
 nonisotropic.cpp  | 21 ++++++-------
 racing.cpp        | 11 +++----
 rug.cpp           | 14 ++++-----
 screenshot.cpp    |  2 +-
 shmup.cpp         |  6 ++--
 16 files changed, 138 insertions(+), 85 deletions(-)

diff --git a/archimedean.cpp b/archimedean.cpp
index f4b6fef1..2730bd33 100644
--- a/archimedean.cpp
+++ b/archimedean.cpp
@@ -625,7 +625,7 @@ struct hrmap_archimedean : hrmap {
       dynamicval<hrmap*> cm(currentmap, current_altmap);
       U = T;
       current_altmap->virtualRebase(alt, T);
-      U = U * inverse(T);
+      U = U * iso_inverse(T);
       }
     
     if(euclid) {
@@ -725,7 +725,7 @@ struct hrmap_archimedean : hrmap {
         h1 = h1->move(0);
         }
       else {
-        where = inverse(cof.adjcell_matrix(h2, 0)) * where;
+        where = iso_inverse(cof.adjcell_matrix(h2, 0)) * where;
         h2 = h2->move(0);
         }
       }
diff --git a/binary-tiling.cpp b/binary-tiling.cpp
index adf9e537..0c3f0a01 100644
--- a/binary-tiling.cpp
+++ b/binary-tiling.cpp
@@ -762,7 +762,7 @@ EX namespace bt {
       t[0] = parabolic3(horohex_scale, 0) * xpush(-l) * cspin(1, 2, M_PI/2);
       t[1] = cspin(1, 2, 2*M_PI/3) * t[0];
       t[2] = cspin(1, 2, 4*M_PI/3) * t[0];
-      auto it = inverse(t[0]);
+      auto it = iso_inverse(t[0]);
 
       t[5] = it * t[1] * t[1];
       t[6] = it * t[5];
@@ -782,7 +782,7 @@ EX namespace bt {
       use_direct >>= 1;
       }
     for(int i=0; i<S7; i++) if(use_direct_for(i))
-      inverse_tmatrix[i] = inverse(direct_tmatrix[i]);
+      inverse_tmatrix[i] = iso_inverse(direct_tmatrix[i]);
     }
   
   #if MAXMDIM == 4
diff --git a/control.cpp b/control.cpp
index 7769fb35..dd1d610a 100644
--- a/control.cpp
+++ b/control.cpp
@@ -331,9 +331,9 @@ EX void full_rotate_camera(int dir, ld val) {
   if(rug::rug_control() && lshiftclick) {
     hyperpoint h;
     if(nonisotropic) {
-      transmatrix T2 = eupush( tC0(inverse(View)) );
+      transmatrix T2 = eupush( tC0(view_inverse(View)) );
       transmatrix nlp = View * T2;  
-      auto rV = inverse(nlp) * View;
+      auto rV = ortho_inverse(nlp) * View;
       h = nlp * inverse_exp(shiftless(tC0(rV)));
       }
     else h = inverse_exp(shiftless(tC0(View)));
diff --git a/drawing.cpp b/drawing.cpp
index cc22a961..8c03e6aa 100644
--- a/drawing.cpp
+++ b/drawing.cpp
@@ -1278,10 +1278,10 @@ void draw_stretch(dqi_poly *p) {
   npoly.V = shiftless(Id);
   npoly.flags &= ~(POLY_INVERSE | POLY_FORCE_INVERTED);
   
-  transmatrix T2 = stretch::translate( tC0(inverse(View)) );
+  transmatrix T2 = stretch::translate( tC0(iso_inverse(View)) );
   transmatrix U = View * T2;
   
-  transmatrix iUV = inverse(U) * p->V.T;
+  transmatrix iUV = iso_inverse(U) * p->V.T;
   
   vector<hyperpoint> hs;
   vector<hyperpoint> ths;
diff --git a/geometry.cpp b/geometry.cpp
index 6f71ed2d..aa9d752d 100644
--- a/geometry.cpp
+++ b/geometry.cpp
@@ -577,7 +577,7 @@ void geometry_information::prepare_basics() {
   for(int d=0; d<S7; d++) 
     hexmove[d] = spin(hexshift-d * ALPHA) * xpush(-crossf)* spin(M_PI);  
 
-  for(int d=0; d<S7; d++) invhexmove[d] = inverse(hexmove[d]);
+  for(int d=0; d<S7; d++) invhexmove[d] = iso_inverse(hexmove[d]);
 
   hexvdist = hdist(xpush0(hexf), xspinpush0(ALPHA/2, hcrossf));
 
@@ -951,7 +951,7 @@ EX void switch_always3() {
     #endif
     if(dual::split(switch_fpp)) return;
     check_cgi(); cgi.require_basics();
-    View = inverse(models::rotmatrix()) * View;
+    View = iso_inverse(models::rotmatrix()) * View;
     if(!vid.always3) {
       vid.always3 = true;
       apply_always3();
diff --git a/geometry2.cpp b/geometry2.cpp
index 8ba5f9fd..7858264b 100644
--- a/geometry2.cpp
+++ b/geometry2.cpp
@@ -69,7 +69,7 @@ EX transmatrix master_relative(cell *c, bool get_inverse IS(false)) {
     else {
       auto li = gp::get_local_info(c);
       transmatrix T = spin(master_to_c7_angle()) * cgi.gpdata->Tf[li.last_dir][li.relative.first&31][li.relative.second&31][gp::fixg6(li.total_dir)];
-      if(get_inverse) T = inverse(T);
+      if(get_inverse) T = iso_inverse(T);
       return T;
       }
     }
diff --git a/graph.cpp b/graph.cpp
index 745493f7..375f5cb5 100644
--- a/graph.cpp
+++ b/graph.cpp
@@ -2315,7 +2315,7 @@ EX bool applyAnimation(cell *c, shiftmatrix& V, double& footphase, int layer) {
   
   if(sl2) {
     a.wherenow = slr::translate(tC0(a.wherenow));
-    hyperpoint h = tC0(inverse(a.wherenow));
+    hyperpoint h = tC0(iso_inverse(a.wherenow));
     hyperpoint ie = slr::get_inverse_exp(shiftless(h));
     auto R = hypot_d(3, ie);
     aspd *= (1+R+(shmup::on?1:0));
@@ -2346,9 +2346,9 @@ EX bool applyAnimation(cell *c, shiftmatrix& V, double& footphase, int layer) {
   else {
     hyperpoint wnow;
     if(a.attacking == 1)
-      wnow = tC0(inverse(a.wherenow) * a.attackat);
+      wnow = tC0(z_inverse(a.wherenow) * a.attackat);
     else
-      wnow = tC0(inverse(a.wherenow));
+      wnow = tC0(z_inverse(a.wherenow));
     
     if(prod) {
       auto d = product_decompose(wnow);
@@ -2778,7 +2778,7 @@ EX bool drawMonster(const shiftmatrix& Vparam, int ct, cell *c, color_t col, col
       }
     else if(NONSTDVAR) {
       transmatrix T = currentmap->adj(c, c->mondir);
-      Vb = Vb * T * rspintox(tC0(inverse(T))) * xpush(cgi.tentacle_length);
+      Vb = Vb * T * rspintox(tC0(iso_inverse(T))) * xpush(cgi.tentacle_length);
       }
     else {
       Vb = Vb * ddspin(c, c->mondir, M_PI);
@@ -3174,7 +3174,7 @@ void draw_movement_arrows(cell *c, const transmatrix& V, int df) {
         transmatrix Centered = rgpushxto0(unshift(tC0(cwtV)));
         int sd = md.subdir;
         
-        transmatrix T = inverse(Centered) * rgpushxto0(Centered * tC0(V)) * rspintox(Centered*tC0(V)) * spin(-sd * M_PI/S7) * xpush(0.2);
+        transmatrix T = iso_inverse(Centered) * rgpushxto0(Centered * tC0(V)) * rspintox(Centered*tC0(V)) * spin(-sd * M_PI/S7) * xpush(0.2);
         
         if(vid.axes >= 5)
           queuestr(shiftless(T), keysize, s0 + key, col >> 8, 1);
@@ -3477,7 +3477,7 @@ EX void pushdown(cell *c, int& q, const shiftmatrix &V, double down, bool rezoom
     // xyscale = xyscale + (zscale-xyscale) * (1+sin(ticks / 1000.0)) / 2;
     
     ptd.V.T = xyzscale( V.T, xyscale*zscale, zscale)
-      * inverse(V.T) * unshift(ptd.V, V.shift);
+      * z_inverse(V.T) * unshift(ptd.V, V.shift);
       
     if(!repriority) ;
     else if(nlev < -vid.lake_bottom-1e-3) {
@@ -3692,7 +3692,7 @@ void make_clipping_planes() {
     hyperpoint sx = point3(y1 * z2 - y2 * z1, z1 * x2 - z2 * x1, x1 * y2 - x2 * y1);
     sx /= hypot_d(3, sx);
     sx[3] = 0;
-    if(nisot::local_perspective_used()) sx = inverse(NLP) * sx;
+    if(nisot::local_perspective_used()) sx = ortho_inverse(NLP) * sx;
     clipping_planes.push_back(sx);
     };
   ld tx = current_display->tanfov;
@@ -4461,7 +4461,7 @@ EX transmatrix actual_view_transform;
 EX ld wall_radar(cell *c, transmatrix T, transmatrix LPe, ld max) {
   if(!in_perspective() || !vid.use_wall_radar) return max;
   transmatrix ori;
-  if(prod) ori = inverse(LPe);
+  if(prod) ori = ortho_inverse(LPe);
   ld step = max / 20;
   ld fixed_yshift = 0;
   for(int i=0; i<20; i++) {
@@ -4487,30 +4487,30 @@ EX void make_actual_view() {
   if(GDIM == 3) {
     ld max = WDIM == 2 ? vid.camera : vid.yshift;
     if(max) {
-      transmatrix Start = inverse(actual_view_transform * View);
+      transmatrix Start = view_inverse(actual_view_transform * View);
       ld d = wall_radar(centerover, Start, NLP, max);
-      actual_view_transform = get_shift_view_of(ztangent(d), actual_view_transform * View) * inverse(View); 
+      actual_view_transform = get_shift_view_of(ztangent(d), actual_view_transform * View) * view_inverse(View); 
       }
-    camera_level = asin_auto(tC0(inverse(actual_view_transform * View))[2]);
+    camera_level = asin_auto(tC0(view_inverse(actual_view_transform * View))[2]);
     }
   if(nonisotropic) {
     transmatrix T = actual_view_transform * View;
-    transmatrix T2 = eupush( tC0(inverse(T)) );
+    transmatrix T2 = eupush( tC0(view_inverse(T)) );
     NLP = T * T2;
-    actual_view_transform = inverse(NLP) * actual_view_transform;
+    actual_view_transform = ortho_inverse(NLP) * actual_view_transform;
     }
   #endif
   #if MAXMDIM >= 4
   if(GDIM == 3 && WDIM == 2) {
     transmatrix T = actual_view_transform * View;
-    transmatrix U = inverse(T);
+    transmatrix U = view_inverse(T);
     
     if(T[0][2]) 
       T = spin(-atan2(T[0][2], T[1][2])) * T;
     if(T[1][2] && T[2][2]) 
       T = cspin(1, 2, -atan2(T[1][2], T[2][2])) * T;
 
-    ld z = -asin_auto(tC0(inverse(T)) [2]);
+    ld z = -asin_auto(tC0(view_inverse(T)) [2]);
     T = zpush(-z) * T;
 
     radar_transform = T * U;
@@ -5297,7 +5297,7 @@ EX void animateAttack(const movei& m, int layer) {
   bool newanim = !animations[layer].count(m.s);
   animation& a = animations[layer][m.s];
   a.attacking = 1;
-  a.attackat = rspintox(tC0(inverse(T))) * xpush(hdist0(T*C0) / 3);
+  a.attackat = rspintox(tC0(iso_inverse(T))) * xpush(hdist0(T*C0) / 3);
   if(newanim) a.wherenow = Id, a.ltick = ticks, a.footphase = 0;
   }
 
diff --git a/hyperpoint.cpp b/hyperpoint.cpp
index 33bba721..40d03b15 100644
--- a/hyperpoint.cpp
+++ b/hyperpoint.cpp
@@ -567,10 +567,10 @@ EX transmatrix eupush(ld x, ld y) {
   return T;
   }
 
-EX transmatrix eupush(hyperpoint h) {
-  if(nonisotropic) return nisot::translate(h);
+EX transmatrix eupush(hyperpoint h, ld co IS(1)) {
+  if(nonisotropic) return nisot::translate(h, co);
   transmatrix T = Id;
-  for(int i=0; i<GDIM; i++) T[i][LDIM] = h[i];
+  for(int i=0; i<GDIM; i++) T[i][LDIM] = h[i] * co;
   return T;
   }
 
@@ -795,11 +795,8 @@ EX transmatrix rpushxto0(const hyperpoint& H) {
   }
 
 EX transmatrix ggpushxto0(const hyperpoint& H, ld co) {
-  if(translatable) {
-    if(nonisotropic)
-      return co > 0 ? eupush(H) : inverse(eupush(H));
-    return eupush(co * H);
-    }
+  if(translatable) 
+    return eupush(H, co);
   if(prod) {
     auto d = product_decompose(H);
     return mscale(PIU(ggpushxto0(d.second, co)), d.first * co);
@@ -924,7 +921,7 @@ void inverse_error(const transmatrix& T) {
   println(hlog, "Warning: inverting a singular matrix: ", T);
   }
 
-/** inverse */
+/** inverse of a 3x3 matrix */
 EX transmatrix inverse3(const transmatrix& T) {
   ld d = det(T);
   transmatrix T2;
@@ -939,6 +936,7 @@ EX transmatrix inverse3(const transmatrix& T) {
   return T2;
   }
 
+/** \brief inverse of a general matrix */
 EX transmatrix inverse(const transmatrix& T) {
   if(MDIM == 3) 
     return inverse3(T);
@@ -978,6 +976,59 @@ EX transmatrix inverse(const transmatrix& T) {
     }
   }
 
+/** \brief inverse of an orthogonal matrix, i.e., transposition */
+EX transmatrix ortho_inverse(transmatrix T) {
+  for(int i=1; i<MDIM; i++)
+    for(int j=0; j<i; j++)
+      swap(T[i][j], T[j][i]);
+  return T;
+  }
+
+/** \brief inverse of an orthogonal matrix in Minkowski space */
+EX transmatrix pseudo_ortho_inverse(transmatrix T) {
+  for(int i=1; i<MDIM; i++)
+    for(int j=0; j<i; j++)
+      swap(T[i][j], T[j][i]);
+  for(int i=0; i<MDIM-1; i++)
+    T[i][MDIM-1] = -T[i][MDIM-1],
+    T[MDIM-1][i] = -T[MDIM-1][i];
+  return T;
+  }
+
+/** \brief inverse of an isometry -- in most geometries this can be done more efficiently than using inverse */
+EX transmatrix iso_inverse(const transmatrix& T) {
+  if(hyperbolic)
+    return pseudo_ortho_inverse(T);
+  if(sphere) 
+    return ortho_inverse(T);
+  if(euclid && !(cgflags & qAFFINE)) {
+    transmatrix U = Id;
+    for(int i=0; i<MDIM-1; i++)
+      for(int j=0; j<MDIM-1; j++)
+        U[i][j] = T[j][i];
+    hyperpoint h = U * tC0(T);
+    for(int i=0; i<MDIM-1; i++)
+      U[i][MDIM-1] = -h[i];
+    return U;
+    }
+  return inverse(T);
+  }
+
+/** \brief T inverse a matrix T = O*S, where O is isometry and S is a scaling matrix (todo optimize) */
+EX transmatrix z_inverse(const transmatrix& T) {
+  return inverse(T);
+  }
+
+/** \brief T inverse a matrix T = O*P, where O is orthogonal and P is an isometry (todo optimize) */
+EX transmatrix view_inverse(transmatrix T) {
+  return inverse(T);
+  }
+
+/** \brief T inverse a matrix T = P*O, where O is orthogonal and P is an isometry (todo optimize) */
+EX transmatrix iview_inverse(transmatrix T) {
+  return inverse(T);
+  }
+
 EX pair<ld, hyperpoint> product_decompose(hyperpoint h) {
   ld z = zlevel(h);
   return make_pair(z, mscale(h, -z));
@@ -1050,7 +1101,7 @@ EX ld hdist(const hyperpoint& h1, const hyperpoint& h2) {
       return hypot(PIU(hdist(d1.second, d2.second)), d1.first - d2.first);
       }
     case gcSL2: 
-      return hdist0(inverse(slr::translate(h1)) * h2);
+      return hdist0(stretch::itranslate(h1) * h2);
     default:
       if(iv < 0) return 0;
       return sqrt(iv);
@@ -1349,12 +1400,12 @@ EX hyperpoint inverse_exp(const shiftpoint h, flagtype prec IS(pNORMAL)) {
 
 EX ld geo_dist(const hyperpoint h1, const hyperpoint h2, flagtype prec IS(pNORMAL)) {
   if(!nonisotropic) return hdist(h1, h2);
-  return hypot_d(3, inverse_exp(shiftless(inverse(nisot::translate(h1)) * h2, prec)));
+  return hypot_d(3, inverse_exp(shiftless(nisot::translate(h1, -1) * h2, prec)));
   }
 
 EX ld geo_dist(const shiftpoint h1, const shiftpoint h2, flagtype prec IS(pNORMAL)) {
   if(!nonisotropic) return hdist(h1, h2);
-  return hypot_d(3, inverse_exp(shiftless(inverse(nisot::translate(h1.h)) * h2.h, h2.shift - h1.shift), prec));
+  return hypot_d(3, inverse_exp(shiftless(nisot::translate(h1.h, -1) * h2.h, h2.shift - h1.shift), prec));
   }
 
 EX ld geo_dist_q(const hyperpoint h1, const hyperpoint h2, flagtype prec IS(pNORMAL)) {
diff --git a/hypgraph.cpp b/hypgraph.cpp
index 2f77b01f..30c88ae2 100644
--- a/hypgraph.cpp
+++ b/hypgraph.cpp
@@ -1523,11 +1523,11 @@ EX void spinEdge(ld aspd) {
     fixmatrix(our);
     fixmatrix(their);
     if(GDIM == 2) {
-      transmatrix T = their * inverse(our);
+      transmatrix T = their * iso_inverse(our);
       hyperpoint H = T * xpush0(1);
       downspin = -atan2(H[1], H[0]);
       }
-    else rotate_view(their * inverse(our));
+    else rotate_view(their * iso_inverse(our));
     }
   else if(playerfound && vid.fixed_facing) {
     hyperpoint H = gpushxto0(unshift(playerV) * C0) * unshift(playerV) * xpush0(5);
@@ -1576,8 +1576,8 @@ EX void centerpc(ld aspd) {
     transmatrix T = pc->at;
     int sl = snakelevel(cwt.at);
     if((sl || vid.eye) && WDIM == 2) T = T * zpush(cgi.SLEV[sl] - cgi.FLOOR + vid.eye);
-    View = inverse(T);
-    if(prod) NLP = inverse(pc->ori);
+    View = iso_inverse(T);
+    if(prod) NLP = ortho_inverse(pc->ori);
     if(WDIM == 2) rotate_view( cspin(0, 1, M_PI) * cspin(2, 1, M_PI/2 + shmup::playerturny[id]) * spin(-M_PI/2) );
     return;
     }
@@ -1610,7 +1610,7 @@ EX void centerpc(ld aspd) {
     if(sl || vid.eye) T = T * zpush(cgi.SLEV[sl] - cgi.FLOOR + vid.eye);
     }
   #endif
-  hyperpoint H = inverse(actual_view_transform) * tC0(T);
+  hyperpoint H = ortho_inverse(actual_view_transform) * tC0(T);
   ld R = (zero_d(GDIM, H) && !prod) ? 0 : hdist0(H);
   if(R < 1e-9) {
     // either already centered or direction unknown
@@ -1649,7 +1649,7 @@ EX void optimizeview() {
   if(dual::split(optimizeview)) return;
   
   cell *c = centerover;
-  transmatrix iView = inverse(View);
+  transmatrix iView = view_inverse(View);
   virtualRebase(centerover, iView);
   if(c != centerover && (sphere || sl2)) {
     transmatrix T = currentmap->relative_matrix(centerover, c, C0);
@@ -1657,7 +1657,7 @@ EX void optimizeview() {
     stretch::mstretch_matrix = T * stretch::mstretch_matrix;
     }
     
-  View = inverse(iView);
+  View = iview_inverse(iView);
   fixmatrix(View);
 
   #if CAP_ANIMATIONS
@@ -1732,7 +1732,7 @@ EX void resetview() {
     }
   cwtV = shiftless(View);
   current_display->which_copy = 
-    nonisotropic ? gpushxto0(tC0(inverse(View))) :
+    nonisotropic ? gpushxto0(tC0(view_inverse(View))) :
     View;
   // SDL_LockSurface(s);
   // SDL_UnlockSurface(s);
@@ -2266,11 +2266,11 @@ EX hyperpoint unshift(shiftpoint T, ld to IS(0)) {
   }
 
 EX transmatrix inverse_shift(const shiftmatrix& T1, const shiftmatrix& T2) {
-  return inverse(T1.T) * unshift(T2, T1.shift);
+  return iso_inverse(T1.T) * unshift(T2, T1.shift);
   }
 
 EX hyperpoint inverse_shift(const shiftmatrix& T1, const shiftpoint& T2) {
-  return inverse(T1.T) * unshift(T2, T1.shift);
+  return iso_inverse(T1.T) * unshift(T2, T1.shift);
   }
 
 EX void optimize_shift(shiftmatrix& T) {
@@ -2497,13 +2497,13 @@ EX transmatrix get_shift_view_of(const hyperpoint H, const transmatrix V) {
     return rgpushxto0(direct_exp(lp_iapply(H))) * V;
     }
   else if(!nisot::geodesic_movement) {
-    transmatrix IV = inverse(V);
+    transmatrix IV = view_inverse(V);
     nisot::fixmatrix(IV);
-    const transmatrix V1 = inverse(IV);
+    const transmatrix V1 = iview_inverse(IV);
     return V1 * eupush(IV * eupush(H) * V1 * C0);
     }
   else {
-    return inverse(nisot::parallel_transport(inverse(V), -H));
+    return iview_inverse(nisot::parallel_transport(view_inverse(V), -H));
     }
   }
 
@@ -2514,8 +2514,8 @@ EX void shift_view(hyperpoint H) {
   View = get_shift_view_of(H, View);
   auto& wc = current_display->which_copy;
   if(nonisotropic || stretch::in()) {
-    transmatrix ioldv = eupush(tC0(inverse(oView)));
-    transmatrix newv = inverse(eupush(tC0(inverse(View))));
+    transmatrix ioldv = eupush(tC0(view_inverse(oView)));
+    transmatrix newv = inverse(eupush(tC0(view_inverse(View))));
     wc = newv * ioldv * wc;
     }
   else
@@ -2586,7 +2586,7 @@ EX void set_view(hyperpoint camera, hyperpoint forward, hyperpoint upward) {
   
   if(det(rotator) < 0) rotator[0] = -rotator[0];
   
-  View = rotator * inverse(rgpushxto0(camera));
+  View = rotator * iso_inverse(rgpushxto0(camera));
   }
 
 }
diff --git a/mapeditor.cpp b/mapeditor.cpp
index d7d4b3c1..fad0ab0c 100644
--- a/mapeditor.cpp
+++ b/mapeditor.cpp
@@ -1565,7 +1565,7 @@ EX namespace mapeditor {
     ld step = 0.01;
     ld cdist = 0;
     
-    auto idt = inverse(unshift(drawtrans));
+    auto idt = z_inverse(unshift(drawtrans));
 
     auto qu = [&] (ld d) {
       ld d1 = front_edit;
@@ -1889,7 +1889,7 @@ EX namespace mapeditor {
     if(!mouseout()) {
       hyperpoint mh;
       if(GDIM == 2) {
-        transmatrix T = inverse(unshift(drawtrans)) * rgpushxto0(ccenter); /* todo? */
+        transmatrix T = z_inverse(unshift(drawtrans)) * rgpushxto0(ccenter); /* todo? */
         mh = spintox(gpushxto0(ccenter) * coldcenter) * T * unshift(mouseh);
         }
       else
@@ -2651,7 +2651,7 @@ EX namespace mapeditor {
     
     if(sio.id == sih.id) {
       if(c == ls)
-        textrans = inverse(V.T * applyPatterndir(ls, sio));
+        textrans = z_inverse(V.T * applyPatterndir(ls, sio));
       
       transmatrix mh = textrans * rgpushxto0(unshift(mouseh, V.shift));
       transmatrix ml = textrans * rgpushxto0(unshift(lstart, V.shift));
diff --git a/models.cpp b/models.cpp
index 98396961..42959fc1 100644
--- a/models.cpp
+++ b/models.cpp
@@ -173,7 +173,7 @@ EX namespace models {
       spiral_multiplier = cld(cos_spiral, sin_spiral) * cld(spiral_cone_rad * mul / 2., 0);
       }
     if(euclid) {
-      euclidean_spin = pispin * inverse(cview().T * master_relative(centerover, true));
+      euclidean_spin = pispin * iso_inverse(cview().T * master_relative(centerover, true));
       euclidean_spin = gpushxto0(euclidean_spin * C0) * euclidean_spin;
       hyperpoint h = inverse(euclidean_spin) * (C0 + (euc::eumove(gp::loc{1,0})*C0 - C0) * vpconf.spiral_x + (euc::eumove(gp::loc{0,1})*C0 - C0) * vpconf.spiral_y);
       spiral_multiplier = cld(0, 2 * M_PI) / cld(h[0], h[1]);
diff --git a/nonisotropic.cpp b/nonisotropic.cpp
index 8d611547..642a6fe3 100644
--- a/nonisotropic.cpp
+++ b/nonisotropic.cpp
@@ -16,9 +16,9 @@ EX namespace nisot {
   
   EX bool geodesic_movement = true;
 
-  EX transmatrix translate(hyperpoint h) {
-    if(sl2)
-      return slr::translate(h);
+  EX transmatrix translate(hyperpoint h, ld co IS(1)) {
+    if(sl2 || sphere)
+      return co > 0 ? stretch::translate(h) : stretch::itranslate(h);
     transmatrix T = Id;
     for(int i=0; i<GDIM; i++) T[i][LDIM] = h[i];
     if(sol && nih) {
@@ -35,6 +35,7 @@ EX namespace nisot {
       }
     if(nil) 
       T[2][1] = h[0];
+    if(co < 0) return iso_inverse(T);
     return T;
     }
 
@@ -927,7 +928,7 @@ EX namespace nilv {
     }
 
   EX hyperpoint on_geodesic(hyperpoint s0, hyperpoint s1, ld x) {
-    hyperpoint local = inverse(nisot::translate(s0)) * s1;
+    hyperpoint local = nisot::translate(s0, -1) * s1;
     hyperpoint h = get_inverse_exp(local);
     return nisot::translate(s0) * formula_exp(h * x);
     }
@@ -1372,7 +1373,7 @@ EX namespace hybrid {
           T = spintox(T * h2) * T;
           w = T * C0;
           w[1] = -w[1];
-          w = inverse(T) * w;
+          w = iso_inverse(T) * w;
           };
         if(prod)
           ((hrmap_hybrid*)currentmap)->in_underlying(f);
@@ -2046,7 +2047,7 @@ EX namespace rots {
     transmatrix Spin;
     hybrid::in_underlying_geometry([&] {
       hyperpoint h = tC0(T);
-      Spin = inverse(gpushxto0(h) * T);
+      Spin = iso_inverse(gpushxto0(h) * T);
       d = hr::inverse_exp(shiftless(h));
       alpha = atan2(Spin[0][1], Spin[0][0]);
       distance = hdist0(h);
@@ -2173,7 +2174,7 @@ EX namespace rots {
     auto g = std::move(gmatrix);
     auto g0 = std::move(gmatrix0);
     
-    ld alpha = atan2(inverse(NLP) * point3(1, 0, 0));
+    ld alpha = atan2(ortho_inverse(NLP) * point3(1, 0, 0));
     
     bool inprod = prod;
     transmatrix pView = View;
@@ -2762,7 +2763,7 @@ EX namespace nisot {
   EX void fixmatrix(transmatrix& T) {
     if(sphere) return hr::fixmatrix(T);
     transmatrix push = eupush( tC0(T) );
-    transmatrix push_back = inverse(push);
+    transmatrix push_back = eupush(tC0(T), -1);
     transmatrix gtl = push_back * T;
     { dynamicval<eGeometry> g(geometry, gSphere); hr::fixmatrix(gtl); }
     T = push * gtl;
@@ -2771,14 +2772,14 @@ EX namespace nisot {
   EX transmatrix parallel_transport(const transmatrix Position, const hyperpoint direction) {
     auto P = Position;
     nisot::fixmatrix(P);  
-    if(!geodesic_movement) return inverse(eupush(Position * translate(-direction) * inverse(Position) * C0)) * Position;
+    if(!geodesic_movement) return eupush(Position * translate(-direction) * inverse(Position) * C0, -1) * Position;
     return parallel_transport_bare(P, direction);
     }
   
   EX transmatrix spin_towards(const transmatrix Position, const hyperpoint goal, flagtype prec IS(pNORMAL)) {
 
     hyperpoint at = tC0(Position);
-    transmatrix push_back = inverse(translate(at));
+    transmatrix push_back = translate(at, -1);
     hyperpoint back_goal = push_back * goal;
     back_goal = inverse_exp(shiftless(back_goal), prec);
     
diff --git a/racing.cpp b/racing.cpp
index 31715501..58e0d316 100644
--- a/racing.cpp
+++ b/racing.cpp
@@ -778,7 +778,7 @@ EX bool set_view() {
   transmatrix at = ypush(-vid.yshift) * unshift(ggmatrix(who->base)) * who->at;
   
   if(racing::player_relative || quotient || (kite::in() && GDIM == 3)) {
-    View = inverse(at) * View;
+    View = iso_inverse(at) * View;
     }
   else {
     int z = get_info(who->base).completion;
@@ -787,14 +787,15 @@ EX bool set_view() {
     cell *c2 = racing::track[min(z+steps, isize(racing::track)-1)];
     transmatrix T1 = ypush(-vid.yshift) * unshift(ggmatrix(c1));
     transmatrix T2 = ypush(-vid.yshift) * unshift(ggmatrix(c2));
-    transmatrix T = spintox(inverse(T1) * T2 * C0);
-    hyperpoint h = T * inverse(T1) * at * C0;
+    transmatrix iT1 = iso_inverse(T1);
+    transmatrix T = spintox(iT1 * T2 * C0);
+    hyperpoint h = T * iT1 * at * C0;
     ld y = GDIM == 2 ? asin_auto(h[1]) : asin_auto(hypot(h[1], h[2]));
     ld x = asin_auto(h[0] / cos_auto(y));
     x += race_advance;
     if(GDIM == 3 && race_advance == 0 && pmodel == mdPerspective) race_advance = -1;
-    transmatrix Z = T1 * inverse(T) * xpush(x);
-    View = inverse(Z) * View;
+    transmatrix Z = T1 * iso_inverse(T) * xpush(x);
+    View = iso_inverse(Z) * View;
     if(GDIM == 3) View = cspin(2, 0, M_PI/2) * View;
     fixmatrix(View);
     }
diff --git a/rug.cpp b/rug.cpp
index ed52fc2f..0f994779 100644
--- a/rug.cpp
+++ b/rug.cpp
@@ -183,7 +183,7 @@ EX rugpoint *addRugpoint(shiftpoint h, double dist) {
   m->valid = false;
 
   if(euclid && quotient && !bounded) {
-    hyperpoint h1 = inverse(models::euclidean_spin) * eumove(euc::eu.user_axes[1]) * C0;
+    hyperpoint h1 = iso_inverse(models::euclidean_spin) * eumove(euc::eu.user_axes[1]) * C0;
     h1 /= sqhypot_d(2, h1);
     if(nonorientable) h1 /= 2;
     m->valid = good_shape = true;
@@ -626,7 +626,7 @@ bool force(rugpoint& m1, rugpoint& m2, double rd, bool is_anticusp=false, double
   bool nonzero = abs(t-rd) > err_zero_current;
   double forcev = (t - rd) / 2; // 20.0;
   
-  transmatrix T = inverse(rgpushxto0(m1.native));
+  transmatrix T = iso_inverse(rgpushxto0(m1.native));
   hyperpoint ie = inverse_exp(shiftless(T * m2.native));
 
   transmatrix iT = rgpushxto0(m1.native);
@@ -1086,9 +1086,9 @@ EX void drawRugScene() {
   auto& rug = queuecurve(shiftless(Id), 0, 0xFFFFFFFF, PPR::LINE);
 
   if(nonisotropic) {
-    transmatrix T2 = eupush( tC0(inverse(rugView)) );
+    transmatrix T2 = eupush( tC0(view_inverse(rugView)) );
     NLP = rugView * T2;  
-    rug.V = shiftless(inverse(NLP) * rugView);
+    rug.V = shiftless(ortho_inverse(NLP) * rugView);
     }
   else {
     rug.V = shiftless(rugView);
@@ -1361,9 +1361,9 @@ EX shiftpoint gethyper(ld x, ld y) {
         if(sp == 1 || sp == 2) continue;
         }
     
-      applymodel(shiftless(inverse(NLP) * rugView * r0->native), p0);
-      applymodel(shiftless(inverse(NLP) * rugView * r1->native), p1);
-      applymodel(shiftless(inverse(NLP) * rugView * r2->native), p2);
+      applymodel(shiftless(ortho_inverse(NLP) * rugView * r0->native), p0);
+      applymodel(shiftless(ortho_inverse(NLP) * rugView * r1->native), p1);
+      applymodel(shiftless(ortho_inverse(NLP) * rugView * r2->native), p2);
       }
 
     if(error || spherepoints == 1 || spherepoints == 2) continue;
diff --git a/screenshot.cpp b/screenshot.cpp
index 6071d2b6..c19324b4 100644
--- a/screenshot.cpp
+++ b/screenshot.cpp
@@ -1172,7 +1172,7 @@ EX void reflect_view() {
     cell *mbase = centerover;
     cell *c = centerover;
     if(shmup::reflect(c, mbase, T))
-      View = inverse(T.T) * View;
+      View = iso_inverse(T.T) * View;
     }
   }
 
diff --git a/shmup.cpp b/shmup.cpp
index 3c7ce2fa..f22c368b 100644
--- a/shmup.cpp
+++ b/shmup.cpp
@@ -197,11 +197,11 @@ void monster::rebasePat(const shiftmatrix& new_pat, cell *c2) {
       else {
         auto T = calc_relative_matrix(c2, base, tC0(at));
         base = c2;
-        at = inverse(T) * at;
+        at = iso_inverse(T) * at;
         }
       }
     if(multi::players == 1 && this == shmup::pc[0] && !eqmatrix(old_at, at))
-      current_display->which_copy = current_display->which_copy * old_at * inverse(at);
+      current_display->which_copy = current_display->which_copy * old_at * iso_inverse(at);
     return;
     }
   if(multi::players == 1 && this == shmup::pc[0])
@@ -1354,7 +1354,7 @@ bool reflectmatrix(shiftmatrix& M, cell *c1, cell *c2, bool onlypos) {
   ld d = hdist0(tC0(H));
   transmatrix T = xpush(-d/2) * S * inverse_shift(gmatrix[c1], M);
   if(onlypos && tC0(T)[0] < 0) return false;
-  M = gmatrix[c1] * inverse(S) * xpush(d/2) * MirrorX * T;
+  M = gmatrix[c1] * iso_inverse(S) * xpush(d/2) * MirrorX * T;
   return true;
   }