Merge pull request #288 from tylerjarvis/main

Pulling from Main

Author: wlgns0330

CombinedNotebook.ipynb

@@ -34,12 +34,12 @@ def test_univariate():
     assert len(roots) == 128
     assert np.max(np.abs(f(roots))) < tol1
 
-def test_univariate_power():
+def test_univariate_power():    
     coeff = np.zeros(5)
     coeff[0], coeff[1], coeff[2], coeff[3], coeff[4] = -2, 2, 3, 4, 5
     f = yr.MultiPower(coeff)
 
-    roots = yr.solve(f, -2, 1)
+    roots = yr.solve(f, -1, 1)
     assert len(roots) == 2
     assert np.max(np.abs(f(roots))) < tol1
 
@@ -48,8 +48,8 @@ def test_univariate_cheb():
     coeff[0], coeff[1], coeff[2], coeff[3] = 0, 1, 2, 3
     f = yr.MultiCheb(coeff)
 
-    roots = yr.solve(f, -0.5, 1)
-    assert len(roots) == 2
+    roots = yr.solve(f, -1, 1)
+    assert len(roots) == 3
     assert np.max(np.abs(f(roots))) < tol1
 
 # Test Multidimensional Examples
@@ -103,6 +103,28 @@ def test_multiCheb_multiPower():
     assert np.max(np.abs(f(roots))) < tol2
     assert np.max(np.abs(g(roots))) < tol2
 
+# Test MultiCheb and MultiPower
+def test_multiCheb_multiPower_non_unit_box():
+    """
+    f(x,y) = 5x^3 + 4 xy^2 + 3x^2 + 2y^2 + 1
+    g(x,y) = 5 T_2(x) + 3T_1(x)T_2(y) + 2
+
+    """
+
+    coeff = np.zeros((4,4))
+    coeff[3,0], coeff[1,2], coeff[2,0], coeff[0,2], coeff[0,0] = 5, 4, 3, 2, 1
+    f = yr.MultiPower(coeff)
+
+    coeff = np.zeros((3,3))
+    coeff[2,0], coeff[1,2], coeff[0,0] = 5, 3, 2
+    g = yr.MultiCheb(coeff)
+
+    roots = yr.solve([f,g],[-2,-2],[2,2])
+
+    assert len(roots) == 2
+    assert np.max(np.abs(f(roots))) < tol2
+    assert np.max(np.abs(g(roots))) < tol2
+
 def test_multiPower():
     """
     f(x,y) = 5x^3 + 4 xy^2 + 3x^2 + 2y^2 - 5
@@ -117,9 +139,9 @@ def test_multiPower():
     coeff[2,1], coeff[1,2], coeff[2,0], coeff[0,2], coeff[0,0] = 3, -4, 3, 2, -1
     g = yr.MultiPower(coeff)
 
-    roots = yr.solve([f,g],[-2,-2],[2,2])
+    roots = yr.solve([f,g],[-1,-1],[1,1])
 
-    assert len(roots) == 2
+    assert len(roots) == 1
     assert np.max(np.abs(f(roots))) < tol2
     assert np.max(np.abs(g(roots))) < tol2
 
@@ -208,8 +230,8 @@ def test_exact_option():
     yroots_non_exact = yr.solve(funcs,a,b,exact=False)
     yroots_exact = yr.solve(funcs,a,b,exact=True)
 
-    actual_roots = np.load('../../Polished_results/polished_2.3.npy')
-    chebfun_roots = np.loadtxt('../../Chebfun_results/test_roots_2.3.csv', delimiter=',')
+    actual_roots = np.load('../Polished_results/polished_2.3.npy')
+    chebfun_roots = np.loadtxt('../Chebfun_results/test_roots_2.3.csv', delimiter=',')
 
     assert len(yroots_non_exact) == len(actual_roots)
     assert len(yroots_exact) == len(actual_roots)

tests/test_Combined_Solver.py

@@ -441,9 +441,9 @@ def chebApproximate(f, a, b, relApproxTol=1e-10):
     except TypeError as e:
         raise ValueError("Invalid input: length of the upper/lower bound lists must match the dimension of the function")
 
-    # If the function is a MultiCheb object on [-1,1]^n, then return its matrix as the approximation
-    if isinstance(f,MultiCheb) and np.allclose(a,-np.ones_like(a)) and np.allclose(b,np.ones_like(b)):
-        return f.coeff.astype(float), 0
+    # # If the function is a MultiCheb object on [-1,1]^n, then return its matrix as the approximation
+    # if isinstance(f,MultiCheb) and np.allclose(a,-np.ones_like(a)) and np.allclose(b,np.ones_like(b)):
+    #     return f.coeff.astype(float), 0
 
     # Generate and return the approximation
     degs, epsilons, rhos = getChebyshevDegrees(f, a, b, relApproxTol)

yroots/ChebyshevApproximator.py

@@ -5,6 +5,7 @@
 from scipy.spatial import HalfspaceIntersection, QhullError
 from scipy.optimize import linprog
 from yroots.QuadraticCheck import quadratic_check
+from time import time
 import copy
 import warnings
 
@@ -1238,6 +1239,7 @@ def solvePolyRecursive(Ms, trackedInterval, errors, solverOptions):
     originalIntervalSize = trackedInterval.size()
     #Zoom in while we can
     lastSizes = trackedInterval.dimSize()
+    start_time = time()
     while changed and zoomCount <= solverOptions.maxZoomCount:
         #Zoom in until we stop changing or we hit machine epsilon
         Ms, errors, trackedInterval, changed, should_stop = zoomInOnIntervalIter(Ms, errors, trackedInterval, solverOptions.exact)
@@ -1248,6 +1250,7 @@ def solvePolyRecursive(Ms, trackedInterval, errors, solverOptions):
         if np.all(newSizes >= lastSizes / 2): #Check all dims and use >= to account for a dimension being 0.
             zoomCount += 1
         lastSizes = newSizes
+    finish_time = time()
     if should_stop:
         #Start the final step if the is in the options and we aren't already in it.
         if trackedInterval.finalStep or not solverOptions.useFinalStep:

yroots/ChebyshevSubdivisionSolver.py

@@ -4,7 +4,8 @@
 import functools
 import yroots.ChebyshevSubdivisionSolver as ChebyshevSubdivisionSolver
 import yroots.ChebyshevApproximator as ChebyshevApproximator
-from yroots.polynomial import MultiCheb, MultiPower
+from yroots.polynomial import MultiCheb,MultiPower
+from time import time
 
 def solve(funcs,a=-1,b=1, verbose = False, returnBoundingBoxes = False, exact=False, minBoundingIntervalSize=1e-5):
     """Finds and returns the roots of a system of functions on the search interval [a,b].
@@ -106,18 +107,24 @@ def solve(funcs,a=-1,b=1, verbose = False, returnBoundingBoxes = False, exact=Fa
     polys = np.array(funcs)
     errs = np.array([0.]*dim)
     macheps = 2**-52
+    unit_box = True
+    # Check if original region is in the unit box
+    if not np.allclose(a,-np.ones_like(a)) or not np.allclose(b,np.ones_like(b)):
+        unit_box = False
     # Get an approximation for each function.
     if verbose:
         print("Approximation shapes:", end=" ")
     for i in range(dim):
-        if isinstance(funcs[i], MultiPower):
+        # t = time()
+        if unit_box and isinstance(funcs[i], MultiPower):
             polys[i] = funcs[i].to_cheb()
             errs[i] = macheps
-        elif isinstance(funcs[i], MultiCheb):
+        elif unit_box and isinstance(funcs[i], MultiCheb):
             polys[i] = funcs[i].coeff
             errs[i] = macheps
         else:
             polys[i], errs[i] = ChebyshevApproximator.chebApproximate(funcs[i],a,b)
+        # return time() - t
         if verbose:
             print(f"{i}: {polys[i].shape}", end = " " if i != dim-1 else '\n')
     if verbose: