Python solvers.solve函数代码示例

def test_idiff():
    # the use of idiff in ellipse also provides coverage
    circ = x**2 + y**2 - 4
    ans = -3*x*(x**2 + y**2)/y**5
    assert ans == idiff(circ, y, x, 3).simplify()
    assert ans == idiff(circ, [y], x, 3).simplify()
    assert idiff(circ, y, x, 3).simplify() == ans
    explicit  = 12*x/sqrt(-x**2 + 4)**5
    assert ans.subs(y, solve(circ, y)[0]).equals(explicit)
    assert True in [sol.diff(x, 3).equals(explicit) for sol in solve(circ, y)]
    assert idiff(x + t + y, [y, t], x) == -Derivative(t, x) - 1

    def parameter_value(self, other, t):
        """Return the parameter corresponding to the given point.
        Evaluating an arbitrary point of the entity at this parameter
        value will return the given point.

        Examples
        ========

        >>> from sympy import Line, Point
        >>> from sympy.abc import t
        >>> a = Point(0, 0)
        >>> b = Point(2, 2)
        >>> Line(a, b).parameter_value((1, 1), t)
        {t: 1/2}
        >>> Line(a, b).arbitrary_point(t).subs(_)
        Point2D(1, 1)
        """
        from sympy.geometry.point import Point
        from sympy.core.symbol import Dummy
        from sympy.solvers.solvers import solve
        if not isinstance(other, GeometryEntity):
            other = Point(other, dim=self.ambient_dimension)
        if not isinstance(other, Point):
            raise ValueError("other must be a point")
        T = Dummy('t', real=True)
        sol = solve(self.arbitrary_point(T) - other, T, dict=True)
        if not sol:
            raise ValueError("Given point is not on %s" % func_name(self))
        return {t: sol[0][T]}

def test_issue_8235():
    assert reduce_inequalities(x**2 - 1 < 0) == \
        And(S(-1) < x, x < S(1))
    assert reduce_inequalities(x**2 - 1 <= 0) == \
        And(S(-1) <= x, x <= 1)
    assert reduce_inequalities(x**2 - 1 > 0) == \
        Or(And(-oo < x, x < -1), And(x < oo, S(1) < x))
    assert reduce_inequalities(x**2 - 1 >= 0) == \
        Or(And(-oo < x, x <= S(-1)), And(S(1) <= x, x < oo))

    eq = x**8 + x - 9  # we want CRootOf solns here
    sol = solve(eq >= 0)
    tru = Or(And(rootof(eq, 1) <= x, x < oo), And(-oo < x, x <= rootof(eq, 0)))
    assert sol == tru

    # recast vanilla as real
    assert solve(sqrt((-x + 1)**2) < 1) == And(S(0) < x, x < 2)

def test_idiff():
    x = Symbol('x', real=True)
    y = Symbol('y', real=True)
    t = Symbol('t', real=True)
    f = Function('f')
    g = Function('g')
    # the use of idiff in ellipse also provides coverage
    circ = x**2 + y**2 - 4
    ans = -3*x*(x**2 + y**2)/y**5
    assert ans == idiff(circ, y, x, 3).simplify()
    assert ans == idiff(circ, [y], x, 3).simplify()
    assert idiff(circ, y, x, 3).simplify() == ans
    explicit  = 12*x/sqrt(-x**2 + 4)**5
    assert ans.subs(y, solve(circ, y)[0]).equals(explicit)
    assert True in [sol.diff(x, 3).equals(explicit) for sol in solve(circ, y)]
    assert idiff(x + t + y, [y, t], x) == -Derivative(t, x) - 1
    assert idiff(f(x) * exp(f(x)) - x * exp(x), f(x), x) == (x + 1) * exp(x - f(x))/(f(x) + 1)
    assert idiff(f(x) - y * exp(x), [f(x), y], x) == (y + Derivative(y, x)) * exp(x)
    assert idiff(f(x) - y * exp(x), [y, f(x)], x) == -y + exp(-x) * Derivative(f(x), x)
    assert idiff(f(x) - g(x), [f(x), g(x)], x) == Derivative(g(x), x)

    def analyticalButtonClicked(self):
        self.derTableWidget.clear()
        self.tableWidget.clear()
        self.tableWidget.setRowCount(0)
        self.derTableWidget.setColumnCount(0)
        expr = parse_expr(self.fEdit.text())
        f = lambdify(x, expr, 'numpy')
        diff = expr.diff(x)
        diffroots = solve(diff, x)
        diffroots.insert(0, float("-inf"))
        diffroots.append(float("inf"))
        colPosition = self.derTableWidget.columnCount()
        self.derTableWidget.setRowCount(2)
        self.derTableWidget.insertColumn(colPosition)
        self.derTableWidget.setItem(0, colPosition, QTableWidgetItem("DerInterval"))
        self.derTableWidget.setItem(1, colPosition, QTableWidgetItem("Sign f(x)"))
        for i in range(len(diffroots)):
            colPosition = self.derTableWidget.columnCount()
            self.derTableWidget.insertColumn(colPosition)
            self.derTableWidget.setItem(0, colPosition, QTableWidgetItem(("%f"%(diffroots[i]))))
            self.derTableWidget.setItem(1, colPosition, QTableWidgetItem(str(sign(f(diffroots[i])))))

        xspace = np.linspace(int(diffroots[1]-100), int(diffroots[len(diffroots)-2]+100), int(diffroots[len(diffroots)-2]-diffroots[1]+201))
        #TODO: to end work with long table and algorithm with step
        dct = {}
        for i in xspace:
            dct[i] = sign(f(i))
        self.tableWidget.setColumnCount(2)
        rowCount = self.tableWidget.rowCount()
        self.tableWidget.insertRow(rowCount)
        self.tableWidget.setItem(rowCount, 0, QTableWidgetItem("x"))
        self.tableWidget.setItem(rowCount, 1, QTableWidgetItem("Sign f(x)"))
        keylist = []
        keylist = dct.keys()
        for key in sorted(keylist):
            rowCount = self.tableWidget.rowCount()
            self.tableWidget.insertRow(rowCount)
            self.tableWidget.setItem(rowCount, 0, QTableWidgetItem("%f" %key))
            self.tableWidget.setItem(rowCount, 1, QTableWidgetItem(dct[key]))

        self.derEdit.setText(str(diff))

        intervals = ''
        keylist = sorted(keylist)
        for i in range(len(keylist)-1):
            if (dct[keylist[i]] == '-' and dct[keylist[i+1]] == '+') or (dct[keylist[i]] == '+' and dct[keylist[i+1]] == '-'):
                intervals =( "%s%s" %(intervals, makeinterval(keylist[i], keylist[i+1], f)))
            elif (dct[keylist[i]] == '0'):
                intervals = ('%s{%f} and' % (intervals, keylist[i]))

        self.rootIntervalEdit.setText(intervals)

def test_issue_8235():
    x = Symbol('x', real=True)
    assert reduce_inequalities(x**2 - 1 < 0) == \
        And(S(-1) < x, x < S(1))
    assert reduce_inequalities(x**2 - 1 <= 0) == \
        And(S(-1) <= x, x <= 1)
    assert reduce_inequalities(x**2 - 1 > 0) == \
        Or(And(-oo < x, x < -1), And(x < oo, S(1) < x))
    assert reduce_inequalities(x**2 - 1 >= 0) == \
        Or(And(-oo < x, x <= S(-1)), And(S(1) <= x, x < oo))

    eq = x**8 + x**2 - 9
    sol = solve(eq >= 0)
    known_sol = Or(And(-oo < x, RootOf(x**8 + x**2 - 9, 1) <= x, x < oo), \
            And(-oo < x, x < oo, x <= RootOf(x**8 + x**2 - 9, 0)))
    assert sol == known_sol

def _lambert(eq, x):
    """
    Given an expression assumed to be in the form
        ``F(X, a..f) = a*log(b*X + c) + d*X + f = 0``
    where X = g(x) and x = g^-1(X), return the Lambert solution if possible:
        ``x = g^-1(-c/b + (a/d)*W(d/(a*b)*exp(c*d/a/b)*exp(-f/a)))``.
    """
    eq = _mexpand(expand_log(eq))
    mainlog = _mostfunc(eq, log, x)
    if not mainlog:
        return []  # violated assumptions
    other = eq.subs(mainlog, 0)
    if (-other).func is log:
        eq = (eq - other).subs(mainlog, mainlog.args[0])
        mainlog = mainlog.args[0]
        if mainlog.func is not log:
            return []  # violated assumptions
        other = -(-other).args[0]
        eq += other
    if not x in other.free_symbols:
        return [] # violated assumptions
    d, f, X2 = _linab(other, x)
    logterm = collect(eq - other, mainlog)
    a = logterm.as_coefficient(mainlog)
    if a is None or x in a.free_symbols:
        return []  # violated assumptions
    logarg = mainlog.args[0]
    b, c, X1 = _linab(logarg, x)
    if X1 != X2:
        return []  # violated assumptions

    u = Dummy('rhs')
    sol = []
    # check only real solutions:
    for k in [-1, 0]:
        l = LambertW(d/(a*b)*exp(c*d/a/b)*exp(-f/a), k)
        # if W's arg is between -1/e and 0 there is
        # a -1 branch real solution, too.
        if k and not l.is_real:
            continue
        rhs = -c/b + (a/d)*l

        solns = solve(X1 - u, x)
        for i, tmp in enumerate(solns):
            solns[i] = tmp.subs(u, rhs)
            sol.append(solns[i])
    return sol

def solve_univariate_inequality(expr, gen, assume=True, relational=True):
    """Solves a real univariate inequality.

    Examples
    ========

    >>> from sympy.solvers.inequalities import solve_univariate_inequality
    >>> from sympy.core.symbol import Symbol
    >>> x = Symbol('x', real=True)

    >>> solve_univariate_inequality(x**2 >= 4, x)
    Or(And(-oo < x, x <= -2), And(2 <= x, x < oo))

    >>> solve_univariate_inequality(x**2 >= 4, x, relational=False)
    (-oo, -2] U [2, oo)

    """

    # Implementation for continous functions

    from sympy.solvers.solvers import solve

    solns = solve(expr.lhs - expr.rhs, gen, assume=assume)
    oo = S.Infinity

    start = -oo

    sol_sets = [S.EmptySet]

    for x in sorted(s for s in solns if s.is_real):
        end = x
        if simplify(expr.subs(gen, (start + end)/2 if start != -oo else end - 1)):
            sol_sets.append(Interval(start, end, True, True))

        if simplify(expr.subs(gen, x)):
            sol_sets.append(FiniteSet(x))

        start = end

    end = oo

    if simplify(expr.subs(gen, start + 1)):
        sol_sets.append(Interval(start, end, True, True))

    rv = Union(*sol_sets)
    return rv if not relational else rv.as_relational(gen)

    def _integrate(field=None):
        irreducibles = set()

        for poly in reducibles:
            for z in poly.atoms(Symbol):
                if z in V:
                    break
            else:
                continue

            irreducibles |= set(root_factors(poly, z, filter=field))

        log_coeffs, log_part = [], []
        B = _symbols('B', len(irreducibles))

        for i, poly in enumerate(irreducibles):
            if poly.has(*V):
                log_coeffs.append(B[i])
                log_part.append(log_coeffs[-1] * log(poly))

        coeffs = poly_coeffs + log_coeffs

        candidate = poly_part/poly_denom + Add(*log_part)

        h = F - _derivation(candidate) / denom

        numer = h.as_numer_denom()[0].expand(force=True)

        equations = defaultdict(lambda: S.Zero)

        for term in Add.make_args(numer):
            coeff, dependent = term.as_independent(*V)
            equations[dependent] += coeff

        solution = solve(equations.values(), *coeffs)

        return (solution, candidate, coeffs) if solution else None

def test_slow_general_univariate():
    r = rootof(x**5 - x**2 + 1, 0)
    assert solve(sqrt(x) + 1/root(x, 3) > 1) == \
        Or(And(S(0) < x, x < r**6), And(r**6 < x, x < oo))

def test_issue_10047():
    assert solve(sin(x) < 2) == And(-oo < x, x < oo)

def pde_1st_linear_variable_coeff(eq, func, order, match, solvefun):
    r"""
    Solves a first order linear partial differential equation
    with variable coefficients. The general form of this partial differential equation is

    .. math:: a(x, y) \frac{df(x, y)}{dx} + a(x, y) \frac{df(x, y)}{dy}
                + c(x, y) f(x, y) - G(x, y)

    where `a(x, y)`, `b(x, y)`, `c(x, y)` and `G(x, y)` are arbitrary functions
    in `x` and `y`. This PDE is converted into an ODE by making the following transformation.

    1] `\xi` as `x`

    2] `\eta` as the constant in the solution to the differential equation
    `\frac{dy}{dx} = -\frac{b}{a}`

    Making the following substitutions reduces it to the linear ODE

    .. math:: a(\xi, \eta)\frac{du}{d\xi} + c(\xi, \eta)u - d(\xi, \eta) = 0

    which can be solved using dsolve.

    The general form of this PDE is::

        >>> from sympy.solvers.pde import pdsolve
        >>> from sympy.abc import x, y
        >>> from sympy import Function, pprint
        >>> a, b, c, G, f= [Function(i) for i in ['a', 'b', 'c', 'G', 'f']]
        >>> u = f(x,y)
        >>> ux = u.diff(x)
        >>> uy = u.diff(y)
        >>> genform = a(x, y)*u + b(x, y)*ux + c(x, y)*uy - G(x,y)
        >>> pprint(genform)
                                             d                     d
        -G(x, y) + a(x, y)*f(x, y) + b(x, y)*--(f(x, y)) + c(x, y)*--(f(x, y))
                                             dx                    dy

    Examples
    ========

    >>> from sympy.solvers.pde import pdsolve
    >>> from sympy import Function, diff, pprint, exp
    >>> from sympy.abc import x,y
    >>> f = Function('f')
    >>> eq =  x*(u.diff(x)) - y*(u.diff(y)) + y**2*u - y**2
    >>> pdsolve(eq)
    f(x, y) == F(x*y)*exp(y**2/2) + 1

    References
    ==========

    - Viktor Grigoryan, "Partial Differential Equations"
      Math 124A - Fall 2010, pp.7

    """
    from sympy.integrals.integrals import integrate
    from sympy.solvers.ode import dsolve

    xi, eta = symbols("xi eta")
    f = func.func
    x = func.args[0]
    y = func.args[1]
    b = match[match['b']]
    c = match[match['c']]
    d = match[match['d']]
    e = -match[match['e']]


    if not d:
         # To deal with cases like b*ux = e or c*uy = e
         if not (b and c):
            if c:
                try:
                    tsol = integrate(e/c, y)
                except NotImplementedError:
                    raise NotImplementedError("Unable to find a solution"
                        " due to inability of integrate")
                else:
                    return Eq(f(x,y), solvefun(x) + tsol)
            if b:
                try:
                    tsol = integrate(e/b, x)
                except NotImplementedError:
                    raise NotImplementedError("Unable to find a solution"
                        " due to inability of integrate")
                else:
                    return Eq(f(x,y), solvefun(y) + tsol)

    if not c:
        # To deal with cases when c is 0, a simpler method is used.
        # The PDE reduces to b*(u.diff(x)) + d*u = e, which is a linear ODE in x
        plode = f(x).diff(x)*b + d*f(x) - e
        sol = dsolve(plode, f(x))
        syms = sol.free_symbols - plode.free_symbols - set([x, y])
        rhs = _simplify_variable_coeff(sol.rhs, syms, solvefun, y)
        return Eq(f(x, y), rhs)

    if not b:
        # To deal with cases when b is 0, a simpler method is used.
        # The PDE reduces to c*(u.diff(y)) + d*u = e, which is a linear ODE in y
#.........这里部分代码省略.........

#.........这里部分代码省略.........
                else:
                    p,op,e = subeq
                    #print(p,op,e)
                    p = objs[p][1]
                    e = objs[e][1]
                    op = op.strip()
                    pute = compute(p,op,e,target,problem)
                    #print("OPERATION SELECTED: ",op)
                    #p.details()
                    #e.details()
                    #print(substr,pute[1].num)
                    objs[substr]=pute
                if pute == -1:
                    exit()
                score,c = pute
                thisscore.append(score)
            if target[1][1].entity != c.entity:
                thisscore.append(-0.2)
            #print("WAT",thisscore,c.ent,c.num)
            
            scores.append(sum(thisscore))

            #print(compound)
        m = np.argmax(scores)
        #print(scores[m],ST.equations[m].toString())
        srt = sorted([(x,i) for i,x in enumerate(scores)],reverse=True)
        print('\n Top scoring 3 equations: ')
        for x,i in srt[:3]:
            print(x,ST.equations[i].toString())

        '''
        try:
            if target.ent=='dozen':
                guess = solve('('+numlist[0].num+'/12)'+"-"+target.num,'x')[0]
                print(numlist[0].num+"/12="+target.num)
            else:
                guess = solve(numlist[0].num+"-"+target.num,'x')[0]
                print(numlist[0].num+"="+target.num)
        '''
        eqidxs = [y[0] for y in sorted(enumerate(scores),key=lambda x:x[1],reverse=True)]
        eqnidsx = [x[1] for x in srt]
        seen = []
        tright = 0
        for i in eqidxs:
            if len(seen)>=1:break
            eq = ST.equations[i].toString()
            #eq = eq.replace("=",'-')
            splitEquation = eq.split('=')
            eq = splitEquation[0] + '- (' + splitEquation[1] + ')'
            #print(scores[i], eq)
            try:
                guess = solve(eq,'x')[0]
            except: guess = -1

            # This is the non-negative constraint
            # wrapped in a "check for complex number" try statement :/
            try:
                if guess < 0:
                    continue
            except:
                continue

            #this is a constraint agianst fractional answers when the problem is integers
            if not guess.is_integer:
                if integerproblem:
                    continue

def checkpdesol(pde, sol, func=None, solve_for_func=True):
    """
    Checks if the given solution satisfies the partial differential
    equation.

    pde is the partial differential equation which can be given in the
    form of an equation or an expression. sol is the solution for which
    the pde is to be checked. This can also be given in an equation or
    an expression form. If the function is not provided, the helper
    function _preprocess from deutils is used to identify the function.

    If a sequence of solutions is passed, the same sort of container will be
    used to return the result for each solution.

    The following methods are currently being implemented to check if the
    solution satisfies the PDE:

        1. Directly substitute the solution in the PDE and check. If the
           solution hasn't been solved for f, then it will solve for f
           provided solve_for_func hasn't been set to False.

    If the solution satisfies the PDE, then a tuple (True, 0) is returned.
    Otherwise a tuple (False, expr) where expr is the value obtained
    after substituting the solution in the PDE. However if a known solution
    returns False, it may be due to the inability of doit() to simplify it to zero.

    Examples
    ========
    >>> from sympy import Function, symbols, diff
    >>> from sympy.solvers.pde import checkpdesol, pdsolve
    >>> x, y = symbols('x y')
    >>> f = Function('f')
    >>> eq = 2*f(x,y) + 3*f(x,y).diff(x) + 4*f(x,y).diff(y)
    >>> sol = pdsolve(eq)
    >>> assert checkpdesol(eq, sol)[0]
    >>> eq = x*f(x,y) + f(x,y).diff(x)
    >>> checkpdesol(eq, sol)
    (False, (x*F(4*x - 3*y) - 6*F(4*x - 3*y)/25 + 4*Subs(Derivative(F(_xi_1), _xi_1), (_xi_1,), (4*x - 3*y,)))*exp(-6*x/25 - 8*y/25))
    """

    # Converting the pde into an equation
    if not isinstance(pde, Equality):
        pde = Eq(pde, 0)

    # If no function is given, try finding the function present.
    if func is None:
        try:
            _, func = _preprocess(pde.lhs)
        except ValueError:
            funcs = [s.atoms(AppliedUndef) for s in (
                sol if is_sequence(sol, set) else [sol])]
            funcs = reduce(set.union, funcs, set())
            if len(funcs) != 1:
                raise ValueError(
                    'must pass func arg to checkpdesol for this case.')
            func = funcs.pop()

    # If the given solution is in the form of a list or a set
    # then return a list or set of tuples.
    if is_sequence(sol, set):
        return type(sol)(map(lambda i: checkpdesol(pde, i,
            solve_for_func=solve_for_func), sol))

    # Convert solution into an equation
    if not isinstance(sol, Equality):
        sol = Eq(func, sol)

    # Try solving for the function
    if solve_for_func and not (sol.lhs == func and not sol.rhs.has(func)) and not \
            (sol.rhs == func and not sol.lhs.has(func)):
        try:
            solved = solve(sol, func)
            if not solved:
                raise NotImplementedError
        except NotImplementedError:
            pass
        else:
            if len(solved) == 1:
                result = checkpdesol(pde, Eq(func, solved[0]),
                    order=order, solve_for_func=False)
            else:
                result = checkpdesol(pde, [Eq(func, t) for t in solved],
                order=order, solve_for_func=False)

    # The first method includes direct substitution of the solution in
    # the PDE and simplifying.
    pde = pde.lhs - pde.rhs
    if sol.lhs == func:
        s = pde.subs(func, sol.rhs).doit()
    elif sol.rhs == func:
        s = pde.subs(func, sol.lhs).doit()
    if s:
        ss = simplify(s)
        if ss:
            return False, ss
        else:
            return True, 0
    else:
        return True, 0

def pde_1st_linear_constant_coeff(eq, func, order, match, solvefun):
    r"""
    Solves a first order linear partial differential equation
    with constant coefficients.

    The general form of this partial differential equation is

    .. math:: a \frac{df(x,y)}{dx} + b \frac{df(x,y)}{dy} + c f(x,y) = G(x,y)

    where `a`, `b` and `c` are constants and `G(x, y)` can be an arbitrary
    function in `x` and `y`.

    The general solution of the PDE is::

        >>> from sympy.solvers import pdsolve
        >>> from sympy.abc import x, y, a, b, c
        >>> from sympy import Function, pprint
        >>> f = Function('f')
        >>> G = Function('G')
        >>> u = f(x,y)
        >>> ux = u.diff(x)
        >>> uy = u.diff(y)
        >>> genform = a*u + b*ux + c*uy - G(x,y)
        >>> pprint(genform)
                  d               d
        a*f(x, y) + b*--(f(x, y)) + c*--(f(x, y)) - G(x, y)
                  dx              dy
        >>> pprint(pdsolve(genform, hint='1st_linear_constant_coeff_Integral'))
                  //          b*x + c*y                                             \
                  ||              /                                                 |
                  ||             |                                                  |
                  ||             |                                       a*xi       |
                  ||             |                                     -------      |
                  ||             |                                      2    2      |
                  ||             |      /b*xi + c*eta  -b*eta + c*xi\  b  + c       |
                  ||             |     G|------------, -------------|*e        d(xi)|
                  ||             |      |   2    2         2    2   |               |
                  ||             |      \  b  + c         b  + c    /               |
                  ||             |                                                  |
                  ||            /                                                   |
                  ||                                                                |
        f(x, y) = ||F(eta) + -------------------------------------------------------|*
                  ||                                  2    2                        |
                  \\                                 b  + c                         /
        <BLANKLINE>
                \|
                ||
                ||
                ||
                ||
                ||
                ||
                ||
                ||
          -a*xi ||
         -------||
          2    2||
         b  + c ||
        e       ||
                ||
                /|eta=-b*y + c*x, xi=b*x + c*y


    Examples
    ========

    >>> from sympy.solvers.pde import pdsolve
    >>> from sympy import Function, diff, pprint, exp
    >>> from sympy.abc import x,y
    >>> f = Function('f')
    >>> eq = -2*f(x,y).diff(x) + 4*f(x,y).diff(y) + 5*f(x,y) - exp(x + 3*y)
    >>> pdsolve(eq)
    f(x, y) == (F(4*x + 2*y) + exp(x/2 + 4*y)/15)*exp(x/2 - y)

    References
    ==========

    - Viktor Grigoryan, "Partial Differential Equations"
      Math 124A - Fall 2010, pp.7

    """

    # TODO : For now homogeneous first order linear PDE's having
    # two variables are implemented. Once there is support for
    # solving systems of ODE's, this can be extended to n variables.

    xi, eta = symbols("xi eta")
    f = func.func
    x = func.args[0]
    y = func.args[1]
    b = match[match['b']]
    c = match[match['c']]
    d = match[match['d']]
    e = -match[match['e']]
    expterm = exp(-S(d)/(b**2 + c**2)*xi)
    functerm = solvefun(eta)
    solvedict = solve((b*x + c*y - xi, c*x - b*y - eta), x, y)
    # Integral should remain as it is in terms of xi,
    # doit() should be done in _handle_Integral.
    genterm = (1/S(b**2 + c**2))*C.Integral(
#.........这里部分代码省略.........