It also has serious applications obviously. www.sympy.org/scipy-2017-codegen-tutorial/ mentions code generation capabilities, which sounds super cool!

The code in this section was tested on sympy==1.8 and Python 3.9.5.

Let's start with some basics. fractions:

from sympy import *
sympify(2)/3 + sympify(1)/2

outputs:

7/6

Note that this is an exact value, it does not get converted to floating-point numbers where precision could be lost!

We can also do everything with symbols:

from sympy import *
x, y = symbols('x y')
expr = x/3 + y/2
print(expr)

outputs:

x/3 + y/2

We can now evaluate that expression object at any time:

expr.subs({x: 1, y: 2})

outputs:

4/3

How about a square root?

x = sqrt(2)
print(x)

outputs:

sqrt(2)

so we understand that the value was kept without simplification. And of course:

sqrt(2)**2

outputs 2. Also:

sqrt(-1)

outputs:

I is the imaginary unit. We can use that symbol directly as well, e.g.:

I*I

gives:

-1

Let's do some trigonometry:

cos(pi)

gives:

-1

and:

cos(pi/4)

gives:

sqrt(2)/2

The exponential also works:

exp(I*pi)

gives;

-1

Now for some calculus. To find the derivative of the natural logarithm:

from sympy import *
x = symbols('x')
diff(ln(x), x)

outputs:

1/x

Just read that. One over x. Beauty.

Let's do some more. Let's solve a simple differential equation:

y''(t) - 2y'(t) + y(t) = sin(t)

Doing:

from sympy import *
x = symbols('x')
f, g = symbols('f g', cls=Function)
diffeq = Eq(f(x).diff(x, x) - 2*f(x).diff(x) + f(x), sin(x)**4)
print(dsolve(diffeq, f(x)))

outputs:

Eq(f(x), (C1 + C2*x)*exp(x) + cos(x)/2)

which means:

f (x) = C_{1} + C_{2} x e^{x} + c o s (x) / 2

(1)

To be fair though, it can't do anything crazy, it likely just goes over known patterns that it has solvers for, e.g. if we change it to:

diffeq = Eq(f(x).diff(x, x)**2 + f(x), 0)

it just blows up:

NotImplementedError: solve: Cannot solve f(x) + Derivative(f(x), (x, 2))**2

Sad.

Let's try some polynomial equations:

from sympy import *
x, a, b, c = symbols('x a b c d e f')
eq = Eq(a*x**2 + b*x + c, 0)
sol = solveset(eq, x)
print(sol)

which outputs:

FiniteSet(-b/(2*a) - sqrt(-4*a*c + b**2)/(2*a), -b/(2*a) + sqrt(-4*a*c + b**2)/(2*a))

which is a not amazingly nice version of the quadratic formula. Let's evaluate with some specific constants after the fact:

sol.subs({a: 1, b: 2, c: 3})

which outputs

FiniteSet(-1 + sqrt(2)*I, -1 - sqrt(2)*I)

Let's see if it handles the quartic equation:

x, a, b, c, d, e, f = symbols('x a b c d e f')
eq = Eq(e*x**4 + d*x**3 + c*x**2 + b*x + a, 0)
solveset(eq, x)

Something comes out. It takes up the entire terminal. Naughty. And now let's try to mess with it:

x, a, b, c, d, e, f = symbols('x a b c d e f')
eq = Eq(f*x**5 + e*x**4 + d*x**3 + c*x**2 + b*x + a, 0)
solveset(eq, x)

and this time it spits out something more magic:

ConditionSet(x, Eq(a + b*x + c*x**2 + d*x**3 + e*x**4 + f*x**5, 0), Complexes)

Oh well.

Let's try some linear algebra.

m = Matrix([[1, 2], [3, 4]])

Let's invert it:

m**-1

outputs:

Matrix([
[ -2,    1],
[3/2, -1/2]])

SymPy special function

0 0

@cirosantilli/_file/python/sympy_cheat/python/sympy_cheat/logarithm_integral.py (`li`)

0 0

 Tagged (1)

Partial differential equation solver

 Ancestors (8)

 Incoming links (1)

SageMath

View article source

Discussion (0)

New discussion

There are no discussions about this article yet.

Articles by others on the same topic (1)

Numerical software by

Wikipedia Bot 0 1970-01-01

Read the full article

See all articles in the same topic Create my own version

Discussion (0) Subscribe (1)

Discussion (0)