In python, when computing fixed points and conjugating matrix by unit upper triangular matrix to avoid division by zero, consider two conjugations.

Author: unhyperbolic

python/drilling/test_cases.py

@@ -38,12 +38,12 @@
 Tests drilling geodesics that are entirely in the 2-skeleton::
 
     >>> M.drill_words(['a','acAADa']).canonical_retriangulation().triangulation_isosig(ignore_orientation=False)
-    'ivvPQQcfhghgfghfaaaaaaaaa_BabBBbBaBBbabbab'
+    'ivvPQQcfhghgfghfaaaaaaaaa_BabBBbBaBBaBbbab'
 
 Same test as verified computation::
 
-    sage: M.drill_words(['a','acAADa'], verified = True).canonical_retriangulation().triangulation_isosig(ignore_orientation=False)
-    'ivvPQQcfhghgfghfaaaaaaaaa_BabBBbBaBBbabbab'
+    sage: M.drill_words(['a','acAADa'], verified = True, bits_prec = 60).canonical_retriangulation().triangulation_isosig(ignore_orientation=False)
+    'ivvPQQcfhghgfghfaaaaaaaaa_BabBBbBaBBaBbbab'
 
 Test error when drilling something close to core curve::
 
@@ -70,7 +70,7 @@
 handled by 2-3 moves (in particular, commit f9879d04 introduced a bug):
 
     >>> Manifold("m004").drill_words(['CAC','CCbC']).canonical_retriangulation().triangulation_isosig(ignore_orientation=False)
-    'qLvvLvAMQQQkcgimopkllmpkonnnpixcaelchapewetvrn_bcaaBbBBbaBaBbB'
+    'qLvvLvAMQQQkcgimopkllmpkonnnpixcaelchapewetvrn_bcaaBbBBbBaaBbB'
 
 
 An interesting case where geodesic intersects triangulation in only one tetrahedron:
@@ -91,9 +91,9 @@
     ...         except RuntimeError:
     ...             pass
     >>> drilled_isosig(Manifold('K11n34(0,1)'), ['iFcdbEiFJ', 'iFJ'])
-    'zLLvLLwzAwPQMQzzQkcdgijkjplssrnrotqruvwyxyxyhsgnnighueqdniblsipklpxgcr_BcbDbBba'
+    'zLLvLLwzAwPQMQzzQkcdgijkjplssrnrotqruvwyxyxyhsgnnighueqdniblsipklpxgcr_BcaBbBcB'
     >>> drilled_isosig(Manifold('K11n34(0,1)'), ['iFJ', 'iFcdbEiFJ'])
-    'zLLvLLwzAwPQMQzzQkcdgijkjplssrnrotqruvwyxyxyhsgnnighueqdniblsipklpxgcr_babBbaBcaB'
+    'zLLvLLwzAwPQMQzzQkcdgijkjplssrnrotqruvwyxyxyhsgnnighueqdniblsipklpxgcr_babBcBBcaB'
     >>> sys.setrecursionlimit(original_limit)
 
 Stress test by using large perturbation. In particular, this is testing the
@@ -104,18 +104,18 @@
     >>> from snappy.drilling import perturb
     >>> original_radius = perturb._tube_developing_radius
     >>> perturb._tube_developing_radius = 1
-    >>> Manifold("m307").drill_word('dadadabCdada', verbose=True).isometry_signature(of_link=True) # doctest: +NUMERIC9
-    Tubes lower bound injectivity radius: 0.380575727320247
-    Number of geodesic pieces: [9]
-    Number of tets after subdividing: 45
+    >>> Manifold("m209").drill_word('aaaEDBd', verbose=True).isometry_signature(of_link=True) # doctest: +NUMERIC9
+    Tubes lower bound injectivity radius: 0.505170045966405
+    Number of geodesic pieces: [12]
+    Number of tets after subdividing: 69
     Shortening geodesic by sweeping across triangle.
-    'oLLwQvvPQQcbeefgemnllnmnmlhhaaaaaahaaaaah_bBbabaab'
-    >>> Manifold("m320").drill_word('daaacDA', verbose=True).isometry_signature(of_link=True) # doctest: +NUMERIC9
-    Tubes lower bound injectivity radius: 0.397319067589326
+    'mLALzLQMPbcbefjhhijkllxxnnxagbxnhhn_cbBabBcB'
+    >>> Manifold("m307").drill_word('ADADcBADADAD', verbose=True).isometry_signature(of_link=True) # doctest: +NUMERIC9
+    Tubes lower bound injectivity radius: 0.380575727319372
     Number of geodesic pieces: [9]
     Number of tets after subdividing: 49
     Shortening geodesic by sweeping across triangle.
-    'rLLPwAPvvPQQcccdfehgjiqpooqppqoqffaaaaaaaqaaaqaaa_bBbabaab'
+    'oLLwQvvPQQcbeefgemnllnmnmlhhaaaaaahaaaaah_bBbabaab'
     >>> perturb._tube_developing_radius = original_radius
 
 """

python/geometric_structure/geodesic/fixed_points.py

@@ -29,14 +29,18 @@ def r13_fixed_points_of_psl2c_matrix(m):
     bc = m[1,0] - m[0,1]
     ad = m[1,1] - m[0,0]
     abs_c_pm_m_pp = _lower_bound_abs(bc + ad)
+    abs_c_pp_m_pm = _lower_bound_abs(bc - ad)
 
-    if abs_c > abs_c_pm_m_pp:
+    if abs_c > abs_c_pm_m_pp and abs_c > abs_c_pp_m_pm:
         return _r13_fixed_points_of_psl2c_matrix(m)
 
     pp = make_matrix([[ 1, 0],[ 1, 1]], ring=m.base_ring())
     pm = make_matrix([[ 1, 0],[-1, 1]], ring=m.base_ring())
 
-    tinv, t = pm, pp
+    if abs_c_pm_m_pp > abs_c_pp_m_pm:
+        tinv, t = pm, pp
+    else:
+        tinv, t = pp, pm
 
     pts = _r13_fixed_points_of_psl2c_matrix(tinv * m * t)
     o13_t = psl2c_to_o13(t)