transcendental equation solver for solveset: transolve

sympy/solvers/solveset.py

@@ -14,10 +14,12 @@
 from sympy.core.containers import Tuple
 from sympy.core.facts import InconsistentAssumptions
 from sympy.core.numbers import I, Number, Rational, oo
-from sympy.core.function import (Lambda, expand_complex, AppliedUndef, Function)
+from sympy.core.function import (Lambda, expand_complex, AppliedUndef, Function,
+                                expand_log)
 from sympy.core.relational import Eq
 from sympy.core.symbol import Symbol
 from sympy.simplify.simplify import simplify, fraction, trigsimp
+from sympy.simplify import powdenest
 from sympy.functions import (log, Abs, tan, cot, sin, cos, sec, csc, exp,
                              acos, asin, acsc, asec, arg,
                              piecewise_fold, Piecewise)
@@ -30,9 +32,9 @@
 from sympy.sets.sets import Set
 from sympy.matrices import Matrix
 from sympy.polys import (roots, Poly, degree, together, PolynomialError,
-                         RootOf)
+                         RootOf, factor)
 from sympy.solvers.solvers import (checksol, denoms, unrad,
-    _simple_dens, recast_to_symbols)
+    _simple_dens, recast_to_symbols, _ispow)
 from sympy.solvers.polysys import solve_poly_system
 from sympy.solvers.inequalities import solve_univariate_inequality
 from sympy.utilities import filldedent
@@ -943,7 +945,12 @@ def _solveset(f, symbol, domain, _check=False):
                     elif equation.has(Abs):
                         result += _solve_abs(f, symbol, domain)
                     else:
-                        result += _solve_as_rational(equation, symbol, domain)
+                        new_result = _solve_as_rational(equation, symbol, domain)
+                        if isinstance(new_result, ConditionSet):
+                            # may be a transcendental type equation
+                            result += transolve(equation, symbol, domain)
+                        else:
+                            result += new_result
                 else:
                     result += solver(equation, symbol)
 
@@ -978,6 +985,180 @@ def _solveset(f, symbol, domain, _check=False):
     return result
 
 
+def _expo_solver(f, symbol):
+    """
+    Helper function for solving exponential equations.
+
+    This function solves exponential equation by using logarithms.
+    Exponents are converted to a more general form of logs which can further be
+    solved by passing to _solveset.
+    It deals with equations of type `a*f(x) + b*g(x)`, where f(x) and g(x)
+    are power terms.
+
+    eg: 3**(2*x) - 2**(x + 3) can be transformed to a better log form,
+    2*x*log(3) - (x+3)*log(2), which is easily solvable.
+    """
+
+    a, b = ordered(f.args)
+
+    lhs = a
+    rhs = -b
+
+    lhs = expand_log(log(lhs), force=True)
+    rhs = expand_log(log(rhs), force=True)
+
+    return (lhs - rhs)
+
+
+def _check_expo(f, symbol):
+    """
+    Helper to check whether an equation is exponential or not.
+    Returns True if it is of exponential type otherwise False.
+    """
+    try:
+        a, b = ordered(f.args)
+    except ValueError:
+        return False
+
+    ad = a.as_independent(symbol)[1]
+    bd = b.as_independent(symbol)[1]
+
+    return any(_ispow(i) for i in (ad, bd))
+
+
+def transolve(f, symbol, domain, **flags):
+    """
+    Function to solve transcendental equations. It is a
+    helper to solveset and should be used internally as of now.
+    Handles the following class of transcendental equations:
+
+        - Exponential equations
+        - Logarithmic equations
+        - LambertW type equations.
+        - Trigonometric equations
+
+
+    Parameters
+    ==========
+
+    f: Expr
+       The target equation.
+    symbol: Symbol
+        The variable for which the eqquation is solved.
+    doamin: Set
+        The domain over which the equation is solved.
+    flags: Dictionary
+        Takes care of the recursive calls.
+
+
+    How is transolve better than _tsolve.
+    ====================================
+
+    1) Improved Output
+
+       The output of many types of equations is much better and easy to
+       understand than the ones computed by _tsolve.
+
+       eg: for 3**(2*x) - 2**(x + 3) transolve would return
+       FiniteSet(-3*log(2)/(-2*log(3) + log(2))), whereas _tsolve
+       would return [-log(2**(3/log(2/9)))]. Both are same but the former
+       increases readability and simplicity.
+
+
+    2) Less Complex API
+
+       transolve's API and its flow is simple to understand. Unlike _tsolve which
+       computes by solving with lots of recursive calls.
+       transolve eases the task by reducing it to a two step procedure as
+       dicussed below.
+
+
+    3) Extensible
+
+       transolve is easily extensible, unlike _tsolve which requires in depth knowledge
+       of the method and adding new class of equation is difficult due to its complex
+       structure.
+       To add a new class of equation in transolve one needs to figure out a way to
+       identify the equation and a generalised way to solve that particular
+       class of eqaution.
+
+
+    How transolve works
+    ===================
+
+    The way transolve solves any transcendental equation is very much
+    different from the old solve way of solving as pointed out above.
+
+    transolve solves equations in two step procedure.
+
+    i)  Identification of the type of equation.
+
+        Helpers are used with heuristics implemented to determine
+        if the equation is of a certain type.
+
+    ii) Invoking the respective helper to solve the equation.
+
+        Once identified what family the equation belongs, respective
+        helpers are called and the equation is solved by either returning
+        the solution or by simplifying the original one to the one that can be
+        handled by _solveset.
+
+
+    Examples
+    ========
+
+    >>> from sympy.solvers.solveset import transolve
+
+    >>> x = symbols('x')
+    >>> transolve(5**(x-3) - 3**(2*x + 1), x, S.Reals)
+    FiniteSet(-log(375)/(-log(5) + 2*log(3)))
+    """
+
+    if 'tsolve_saw' not in flags:
+        flags['tsolve_saw'] = []
+    if f in flags['tsolve_saw']:
+        return ConditionSet(symbol, Eq(f, 0), domain)
+    else:
+        flags['tsolve_saw'].append(f)
+
+    inverter = invert_real if domain.is_subset(S.Reals) else invert_complex
+    lhs, rhs_s = inverter(f, 0, symbol, domain)
+
+    result = ConditionSet(symbol, Eq(f, 0), domain)
+
+    if rhs_s is S.EmptySet:
+        result = S.EmptySet
+    else:
+        rhs = rhs_s.args[0]
+        # do we need a loop for rhs_s?
+        # Can't determine a case as of now.
+        if lhs == symbol:
+            result = rhs_s
+
+        elif lhs.is_Add:
+            equation = factor(powdenest(lhs-rhs))
+
+            if equation.is_Mul:
+                result = _solveset(equation, symbol, domain)
+
+            # check if it is exponential type equation
+            elif _check_expo(equation, symbol):
+                new_f = _expo_solver(equation, symbol)
+                result = _solveset(new_f, symbol, domain)
+
+                if isinstance(result, ConditionSet):
+                    result = ConditionSet(symbol, Eq(f, 0), domain)
+
+            else:
+                result = transolve(f, symbol, domain, **flags)
+
+        elif lhs.is_Pow:
+            new_f = _expo_solver(lhs, symbol)
+            result = _solveset(new_f, symbol, domain)
+
+    return result
+
+
 def solveset(f, symbol=None, domain=S.Complexes):
     r"""Solves a given inequality or equation with set as output
 

sympy/solvers/tests/test_solveset.py

@@ -38,7 +38,7 @@
     solveset_real, domain_check, solveset_complex, linear_eq_to_matrix,
     linsolve, _is_function_class_equation, invert_real, invert_complex,
     solveset, solve_decomposition, substitution, nonlinsolve, solvify,
-    _is_finite_with_finite_vars)
+    _is_finite_with_finite_vars, transolve)
 
 
 a = Symbol('a', real=True)
@@ -1613,9 +1613,9 @@ def test_issue_10555():
     f = Function('f')
     g = Function('g')
     assert solveset(f(x) - pi/2, x, S.Reals) == \
-        ConditionSet(x, Eq(2*f(x) - pi, 0), S.Reals)
+        ConditionSet(x, Eq(f(x) - pi/2, 0), S.Reals)
     assert solveset(f(g(x)) - pi/2, g(x), S.Reals) == \
-        ConditionSet(g(x), Eq(2*f(g(x)) - pi, 0), S.Reals)
+        ConditionSet(g(x), Eq(f(g(x)) - pi/2, 0), S.Reals)
 
 
 def test_issue_8715():
@@ -1749,3 +1749,64 @@ def test_issue_14454():
     number = CRootOf(x**4 + x - 1, 2)
     raises(ValueError, lambda: invert_real(number, 0, x, S.Reals))
     assert invert_real(x**2, number, x, S.Reals)  # no error
+
+
+def test_transolve():
+    from sympy.abc import x
+    from sympy.simplify import simplify
+
+    assert transolve(2**x - 32, x, S.Reals) == FiniteSet(log(32)/log(2))
+    assert transolve(3**x, x, S.Reals) == S.EmptySet
+    assert simplify(transolve(3**x - 9**(x+5), x, S.Reals)) == FiniteSet(-10)
+
+
+def test_expo_solver():
+    from sympy.abc import x, y, z
+    from sympy.simplify import simplify
+
+    e1 = 3**(2*x) - 2**(x + 3)
+    e2 = 4**(5 - 9*x) - 8**(2 - x)
+    e3 = 2**x + 4**x
+    e4 = exp(log(5)*x) - 2**x
+    e5 = x**(y*z) - x     # issue 10864
+    e6 = 5**(x/2) - 2**(x/3)
+    e7 = 2**(x) + 4**(x) + 8**(x) - 84
+    e8 = 4**(x+1) + 4**(x+2) + 4**(x-1) - 3**(x+2) -3**(x+3)
+
+    f1 = exp(x/y)*exp(-z/y) - 2
+    f2 = x**x
+
+    assert solveset(e1, x, S.Reals) == FiniteSet(-3*log(2)/(-2*log(3) + log(2)))
+    assert simplify(solveset(e2, x, S.Reals)) == FiniteSet(4/S(15))
+    assert solveset(e3, x, S.Reals) == S.EmptySet
+    assert solveset(e4, x, S.Reals) == FiniteSet(0)
+    assert solveset(e5, x, S.Reals) == FiniteSet(1)
+    assert solveset(e6, x, S.Reals) == FiniteSet(0)
+    assert simplify(solveset(e7, x, S.Reals)) == FiniteSet(2)
+    assert solveset(e8, x, S.Reals) == FiniteSet(2)
+    assert solveset(f1, x, S.Reals) == Intersection(S.Reals, FiniteSet(y*log(2*exp(z/y))))
+    assert solveset(f2, x, S.Reals) == S.EmptySet
+
+
+def test_expo_conditionset():
+    from sympy.abc import x, y
+
+    f1 = (exp(x) + 1)**x - 2
+    f2 = (x + 2)**y*x - 3
+
+    assert solveset(f1, x, S.Reals) == ConditionSet(x, Eq(log(x) - log(-exp(x) - 1), 0), S.Reals)
+    assert solveset(f2, x ,S.Reals) == ConditionSet(x, Eq(x*(x + 2)**y - 3, 0), S.Reals)
+
+@XFAIL
+def test_expo_fail():
+    from sympy.abc import x, y, z, a, b
+
+    assert solveset(z**x - y, x, S.Reals) == Intersection(S.Reals, FiniteSet(log(y)/log(z)))
+    assert solveset(y - a*x**b, x) == FiniteSet((y/a)**(1/b))
+
+    w = symbols('w', integer=True)
+    f1 = 2*x**w - 4*y**w
+    f2 = (x/y)**w - 2
+    ans1 = solveset(f1, w, S.Reals)
+    ans2 = solveset(f2, w, S.Reals)
+    assert ans1 == ans2