Simplest way to solve mathematical equations in Python

Question

I want to solve a set of equations, linear, or sometimes quadratic. I don't have a specific problem, but often, I have been in this situation often.

It is simple to use wolframalpha.com, the web equivalent of Mathematica, to solve them. But that doesn't provide the comfort and convenience of an iPython shell.

Is there a simple library to work on linear and quadratic equations from a python shell?

Personally, I find it extremely convenient to use the Casio 991 MS scientific calculator. I know how to set variables, solve equations, and do a lot. I want such a tool preferably usable from within an ipython shell. I am surprised not to have found any. I'm not impressed enough by sage; perhaps I am missing something.

Answer 1

sympy is exactly what you're looking for.

Answer 2

You discount the best answer as unacceptable.

Your question is "I want a free Computer Algebra System that I can use in Python."

The answer is "SAGE does that."

Have you looked at maxima/macsyma? SAGE provides bindings for it, and that's one of the more powerful free ones.

http://maxima.sourceforge.net/

Answer 3

Here is how to solve your original question using Python (via Sage). This basically clarifies the remark Paul McMillan makes above.

sage: a,b,c = var('a,b,c')
sage: solve([a+b+c==1000, a^2+b^2==c^2], a,b,c)
[[a == 1000*(r1 + sqrt(r1^2 + 2000*r1 - 1000000))/(r1 + sqrt(r1^2 + 2000*r1 - 1000000) + 1000), b == -1/2*r1 - 1/2*sqrt(r1^2 + 2000*r1 - 1000000) + 500, c == r1], [a == 1000*(r2 - sqrt(r2^2 + 2000*r2 - 1000000))/(r2 - sqrt(r2^2 + 2000*r2 - 1000000) + 1000), b == -1/2*r2 + 1/2*sqrt(r2^2 + 2000*r2 - 1000000) + 500, c == r2]]

Answer 4

For inexact solutions, read up on linear programming and sequential quadratic optimization, then search for Python libraries that performs such optimizations for you.

If the equations require integer solutions, you should search for Diophantine equation solvers for Python.

Just note that using a simple solver for Project Euler is missing the point. The fun, and educational part, is learning how to solve it yourself using primitive methods!

Answer 5

A free web-service for solving large-scale systems of nonlinear equations (1 million+) is APMonitor.com. There is a browser interface and an API to Python / MATLAB. The API to Python is a single script (apm.py) that is available for download from the apmonitor.com homepage. Once the script is loaded into a Python code, it gives the ability to solve problems of:

• Nonlinear equations
• Differential and algebraic equations
• Least squares model fitting
• Moving horizon estimation
• Nonlinear model predictive control
• etc.

For the new user, the APM Python software has a Google Groups forum where a user can post questions. There are bi-weekly webinars that showcase optimization problems in operations research and engineering.

Below is an example of an optimization problem (hs71.apm).

Model

  Variables

    x[1] = 1, >=1, <=5

    x[2] = 5, >=1, <=5

    x[3] = 5, >=1, <=5

    x[4] = 1, >=1, <=5

  End Variables



  Equations

    x[1] * x[2] * x[3] * x[4] > 25

    x[1]^2 + x[2]^2 + x[3]^2 + x[4]^2 = 40



    minimize  x[1] * x[4] * (x[1]+x[2]+x[3]) + x[3]

  End Equations

End Model

The optimization problem is solved with the following Python script:

# Import

from apm import *

# Select server

server = 'http://xps.apmonitor.com'

# Application name

app = 'eqn'

# Clear previous application

apm(server,app,'clear all')

# Load model file

apm_load(server,app,'hs71.apm')

# Option to select solver (1=APOPT, 2=BPOPT, 3=IPOPT)

apm_option(server,app,'nlc.solver',3)

# Solve on APM server

solver_output = apm(server,app,'solve')


# Display solver output

print solver_output


# Retrieve results

results = apm_sol(server,app)

# Display results

print '--- Results of the Optimization Problem ---'

print results

# Display Results in Web Viewer 

url = apm_var(server,app)

print "Opened Web Viewer: " + url

Answer 6

Have you looked at SciPy?

It has an example in the tutorials on solving linear algebra:

http://docs.scipy.org/doc/scipy/reference/tutorial/linalg.html#solving-linear-system

Answer 7

For reference: Wolfram Alpha's solution:

a-1000!=0,   b = (1000 (a-500))/(a-1000),   c = (-a^2+1000 a-500000)/(a-1000)

In python, using sympy's solver module (note that it assumes all equations are set equal to zero):

>>> import sympy
>>> a, b, c = sympy.symbols('a, b, c')
>>> sympy.solve([a + b + c - 1000, a**2 + b**2 - c**2], b, c)
[(1000*(a - 500)/(a - 1000), (-a**2 + 1000*a - 500000)/(a - 1000))]

And of course, a != 1000, as a-1000 is the denominator of the two equations.

Answer 8

I have just started using GNU Scientific Library, which however is C library. Looks like there are Python bindings too. So, it might be worth looking at.

Answer 9

I'd use Octave for this but I agree, the syntax of Octave isn't what I'd call thrilling (and the docs always confuse me more than they help, too).

Answer 10

SymPy symbolic Python library demo: symbolically solving math and integrals & pretty-printing the output

From @Autoplectic:

sympy is exactly what you're looking for.

While I have a tendency to write some of the longest answers on Stack Overflow, that is the shortest answer I've seen on Stack Overflow.

Let's add a basic demo.

References and tutorials you'll need:

1. https://www.sympy.org/en/index.html --> click "Get started with the tutorial"
   1. Tutorial index: https://docs.sympy.org/latest/tutorial/index.html Example pages:
      1. Basic Operations
      2. Printing
      3. Simplification
      4. Calculus

Calculus demo I came up with (see here for the main tutorial I looked at to get started):

Short version:

#!/usr/bin/env python3

from sympy import *
x, y, z = symbols('x y z')
init_printing(use_unicode=True)

integral = Integral(exp(-x**2 - y**2), (x, -oo, oo), (y, -oo, oo))
print(pretty(integral))
result = integral.doit() # "do it": evaluate the integral
print(pretty(result))

Longer version:

eRCaGuy_hello_world/math/sympy_integral_and_printing.py:

#!/usr/bin/env python3

from sympy import *
x, y, z = symbols('x y z')
init_printing(use_unicode=True)

# Store an integral
integral = Integral(exp(-x**2 - y**2), (x, -oo, oo), (y, -oo, oo))
# print it in various ways
print("ORIGINAL INTEGRAL:")
print(pretty(integral))
pprint(integral)        # "pretty print": same as above
print(pretty(integral, use_unicode=False))  # pretty print with ASCII instead 
                                            # of unicode
print(integral)         # plain text
print(latex(integral))  # in LaTeX format

# Now solve the integral and output the result by telling it to
# "do it".
result = integral.doit()
# print it in various ways
print("\nRESULT:")
print(pretty(result))
pprint(result)
print(pretty(result, use_unicode=False))
print(result)
print(latex(result))

Note that using pprint(integral) is the same as print(pretty(integral)).

Output of running the above commands:

ORIGINAL INTEGRAL:
∞  ∞                  
⌠  ⌠                  
⎮  ⎮      2    2      
⎮  ⎮   - x  - y       
⎮  ⎮  ℯ          dx dy
⌡  ⌡                  
-∞ -∞                 
∞  ∞                  
⌠  ⌠                  
⎮  ⎮      2    2      
⎮  ⎮   - x  - y       
⎮  ⎮  ℯ          dx dy
⌡  ⌡                  
-∞ -∞                 
 oo  oo                 
  /   /                 
 |   |                  
 |   |      2    2      
 |   |   - x  - y       
 |   |  e          dx dy
 |   |                  
/   /                   
-oo -oo                 
Integral(exp(-x**2 - y**2), (x, -oo, oo), (y, -oo, oo))
\int\limits_{-\infty}^{\infty}\int\limits_{-\infty}^{\infty} e^{- x^{2} - y^{2}}\, dx\, dy

RESULT:
π
π
pi
pi
\pi

The SymPy symbolic math library in Python can do pretty much any kind of math, solving equations, simplifying, factoring, substituting values for variables, pretty printing, converting to LaTeX format, etc. etc. It seems to be a pretty robust solver in my very limited use so far. I recommend trying it out.

Installing it, for me (tested on Linux Ubuntu), was as simple as:

pip3 install sympy

Answer 11

I don't think there is a unified way of dealing with both linear and quadratic (or generally nonlinear) equations simultaneously. With linear systems, python has bindings to linear algebra and matrix packages. Nonlinear problems tend to be solved on a case by case basis.

Answer 12

It depends on your needs:

If you want an interactive graphical interface, then sage is probably the best solution.

If you want to avoid using a graphical interface, but you still want to do computer algebra, then sympy or maxima may cover your needs. (sympy looks very promising, but it still have a long way to go before they can replace mathematica).

If you don't really need symbolic algrebra, but you need a way to program with matrices, solve differential equations, and minimize functions, then scipy or octave are excelent starting points.

Answer 13

Try applying Bisection method in py to find the root given an interval:

def f(x, rhs): # f(x) = e^x
    return math.e ** x - rhs # e^x = rhs -> e^x - rhs = 0

def solve(rhs, a = 0, b = 100, tol = 1e-3):
    while True:
        c  = (a + b) / 2.0
        if(f(a, rhs) * f(c, rhs) > 0):
            a = c
        else:
            b = c
        if(abs(f(c, rhs)) < tol):
            break
    return c

y = math.e ** 3.75 # x = 3.75
print(solve(y)) # 3.7499..

Answer 14

Take a look at this:

http://openopt.org/FuncDesignerDoc#Solving_systems_of_nonlinear_equations

It is extremely easy to use and quite powerful

Answer 15

Well, I just googled into this page by accident. I see many suggestions regarding this and that software tool, but does any tool actually provide an answer? The actual answer is:

[a,b,c] = [200,375,425]

How did I get this? By writing a quick program in the Maxima programming language to find it via "brute force" searching. It only took about 10 minutes to write, seeing as how I'm familiar with the Maxima language. It took a few seconds for the program to run. Here is the program:

euler_solve():= block ( [ a, b, A, B, end:1000],

for a thru end do
    (
    for b thru end do
        (
        c: 1000 -a -b,
        if c < 0 then
            b:end
        else if a^2 + b^2 = c^2 then
            (
            A:a,
            B:b,
            a:end,
            b:end
            )
        )
    ),
return( [A,B,c])
);

You can just cut and paste the above code into the wxMaxima user interface, which I run under Ubuntu and not MS Windows. Then you just enter the function name: euler_solve(), hit return, wait a few seconds, and out pops the answer. This particular kind of problem is so simple that you could use any general-purpose programming language to do the search.

Answer 16

On second thoughts, I went through sage in detail and clearly it is the best math free software available.

Just some of the different python math related libraries, it integrates is absolutely awesome.

Mathematics packages contained in Sage:

Algebra GAP, Maxima, Singular 
Algebraic Geometry  Singular
Arbitrary    Precision
Arithmetic  GMP, MPFR, MPFI,    NTL
Arithmetic Geometry PARI, NTL,      
mwrank, ecm Calculus    Maxima, SymPy, 
GiNaC Combinatorics Symmetrica,     
Sage-Combinat Linear Algebra    Linbox,
IML Graph Theory    NetworkX Group     
Theory  GAP Numerical
computation GSL,    SciPy, NumPy,
ATLAS

Other packages contained in Sage:

Command line    IPython Database    ZODB,
Python Pickles, SQLite Graphical
Interface   Sage Notebook, jsmath
Graphics    Matplotlib, Tachyon3d, GD,
Jmol Interactive 
programming language  Python 
Networking  Twisted