CLASE 5a: SymPy¶

_ SymPy es una biblioteca de Python para matemática simbólica. Apunta a convertirse en un sistema de algebra computacional (CAS) con todas sus prestaciones manteniendo el código tan simple como sea posible para manterlo comprensible y fácilmente extensible. SymPy está escrito totalmente en Python y no requiere bibliotecas adicionales. Este proyecto comenzó en 2005, fue lanzado al público en 2007 y a él han contribuido durante estos años cientos de personas.

_ Otros CAS conocidos son Mathematica y Maple, sin embargo ambos son software privativo y de pago. Aquí puedes encontrar una comparativa de SymPy con Maple. _

Hoy veremos cómo:

Crear símbolos y expresiones.
Manipular expresiones (simplificación, expansión...)
Calcular derivadas e integrales.
Límites y desarrollos en serie.
Resolución de ecuaciones.
Resolción de EDOs.
Matrices

Sin embargo, SymPy no acaba aquí ni mucho menos...

Documentación & SymPy Live Shell¶

In [1]:

from IPython.display import HTML
HTML('<iframe src="http://docs.sympy.org/latest/index.html" width="700" height="400"></iframe>')

Out[1]:

SymPy Gamma¶

In [2]:

HTML('<iframe src="http://www.sympygamma.com/input/?i=integrate%281+%2F+%281+%2B+x^2%29%29" width="700" height="400"></iframe>')

Out[2]:

Creación de símbolos¶

Lo primero, como siempre, es importar aquello que vayamos a necesitar:

In [3]:

# Importación
from sympy import init_session

In [4]:

init_session(use_latex=True)

IPython console for SymPy 0.7.5 (Python 3.4.2-64-bit) (ground types: python)

These commands were executed:
>>> from __future__ import division
>>> from sympy import *
>>> x, y, z, t = symbols('x y z t')
>>> k, m, n = symbols('k m n', integer=True)
>>> f, g, h = symbols('f g h', cls=Function)

Documentation can be found at http://www.sympy.org

WARNING: Hook shutdown_hook is deprecated. Use the atexit module instead.

Lo primero que vemos es que el comando init_session ha llevado a cabo algunas acciones por nostros:

Gracias a use_latex=True obtenemos la salida en $\LaTeX$.
Ha creado una serie de variables para que podamos ponernos a trabajar en el momento.

Nota: En Python, no se declaran las variables, sin embargo, no puedes usar una hasta que no le hayas asignado un valor. Si ahora intentamos crear una variable `a` que sea `a = 2 * b`, veamos qué ocurre:

In [5]:

# Intentamos usar un símbolo que no hemos creado
a = 2 * b

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-5-9024d2cb83d8> in <module>()
      1 # Intentamos usar un símbolo que no hemos creado
----> 2 a = 2 * b

NameError: name 'b' is not defined

Como en b no había sido creada, Python no sabe qué es b.

Esto mismo nos ocurre con los símbolos de SymPy. Antes de usar una variable, debo decir que es un símbolo y asignárselo:

In [6]:

# Creamos el símbolo a
a = symbols('a')
a

Out[6]:

$$a$$

In [7]:

# Número pi
(a + pi) ** 2

Out[7]:

$$\left(a + \pi\right)^{2}$$

In [8]:

# Unidad imaginaria
a + 2 * I

Out[8]:

$$a + 2 i$$

In [9]:

# Número e
E

Out[9]:

$$e$$

In [10]:

# Vemos qué tipo de variable es a
type(a)

Out[10]:

sympy.core.symbol.Symbol

Ahora ya podría crear b = 2 * a:

In [11]:

b = 2 * a
b

Out[11]:

$$2 a$$

In [12]:

type(b)

Out[12]:

sympy.core.mul.Mul

¿Qué está ocurriendo? Python detecta que a es una variable de tipo Symbol y al multiplicarla por 2 devuelve una variable de Sympy.

Como Python permite que el tipo de una variable cambie, si ahora le asigno a a un valor float deja de ser un símbolo.

In [13]:

a = 2.26492
a

Out[13]:

$$2.26492$$

In [14]:

type(a)

Out[14]:

float

Las conclusiones son:

Si quiero usar una variable como símbolo debo crearla previamente.
Las operaciones con símbolos devuelven símbolos.
Si una varibale que almacenaba un símbolo recibe otra asignación, cambia de tipo.

Las variables de tipo Symbol actúan como contenedores en los que no sabemos qué hay (un real, un complejo, una lista...). Hay que tener en cuenta que: una cosa es el nombre de la variable y otra el símbolo con el que se representa.

In [15]:

#creación de símbolos
coef_traccion = symbols('c_T')
coef_traccion

Out[15]:

$$c_{T}$$

Incluso puedo hacer cosas raras como:

In [16]:

# Diferencia entre variable y símbolo
a = symbols('b')
a

Out[16]:

$$b$$

Además, se pueden crear varos símbolos a la vez:

In [17]:

x, y, z, t = symbols('x y z t')

y símbolos griegos:

In [18]:

w = symbols('omega')
W = symbols('Omega')
w, W

Out[18]:

$$\begin{pmatrix}\omega, & \Omega\end{pmatrix}$$

Fuente: Documentación oficial de SymPy

Por defecto, SymPy entiende que los símbolos son números complejos. Esto puede producir resultados inesperados ante determinadas operaciones como, por ejemplo, lo logaritmos. Podemos indicar que la variable es real, entera... en el momento de la creación:

In [19]:

# Creamos símbolos reales
x, y, z, t = symbols('x y z t', real=True)

In [20]:

# Podemos ver las asunciones de un símbolo
x.assumptions0

Out[20]:

{'complex': True,
 'real': True,
 'hermitian': True,
 'imaginary': False,
 'commutative': True}

Expresiones¶

Comencemos por crear una expresión como: $\cos(x)^2+\sin(x)^2$

In [21]:

expr = cos(x)**2 + sin(x)**2
expr

Out[21]:

$$\sin^{2}{\left (x \right )} + \cos^{2}{\left (x \right )}$$

`simplify()`¶

Podemos pedirle que simplifique la expresión anterior:

In [22]:

simplify(expr)

Out[22]:

$$1$$

En este caso parece estar claro lo que quiere decir más simple, pero como en cualquier CAS el comando simplify puede no devolvernos la expresión que nosotros queremos. Cuando esto ocurra necesitaremos usar otras instrucciones.

`.subs()`¶

En algunas ocasiones necesitaremos sustituir una variable por otra, por otra expresión o por un valor.

In [23]:

expr

Out[23]:

$$\sin^{2}{\left (x \right )} + \cos^{2}{\left (x \right )}$$

In [24]:

# Sustituimos x por y ** 2
expr.subs(x, y**2)

Out[24]:

$$\sin^{2}{\left (y^{2} \right )} + \cos^{2}{\left (y^{2} \right )}$$

In [25]:

# ¡Pero la expresión no cambia!
expr

Out[25]:

$$\sin^{2}{\left (x \right )} + \cos^{2}{\left (x \right )}$$

In [26]:

# Para que cambie
expr = expr.subs(x, y**2)
expr

Out[26]:

$$\sin^{2}{\left (y^{2} \right )} + \cos^{2}{\left (y^{2} \right )}$$

Cambia el sin(x) por exp(x)

In [27]:

expr.subs(sin(x), exp(x))

Out[27]:

$$\sin^{2}{\left (y^{2} \right )} + \cos^{2}{\left (y^{2} \right )}$$

Particulariza la expresión $sin(x) + 3 x $ en $x = \pi$

In [28]:

(sin(x) + 3 * x).subs(x, pi)

Out[28]:

$$3 \pi$$

Aunque si lo que queremos es obtener el valor numérico lo mejor es .evalf()

In [29]:

(sin(x) + 3 * x).subs(x, pi).evalf(25)

Out[29]:

$$9.42477796076937971538793$$

In [30]:

#ver pi con 25 decimales
pi.evalf(25)

Out[30]:

$$3.141592653589793238462643$$

In [31]:

#el mismo resultado se obtiene ocn la función N()
N(pi,25)

Out[31]:

$$3.141592653589793238462643$$

Simplificación¶

SymPy ofrece numerosas funciones para simplificar y manipular expresiones. Entre otras, destacan:

expand()
factor()
collect()
apart()
cancel()

Puedes consultar en la documentación de SymPy lo que hace cada una y algunos ejemplos. Existen también funciones específicas de simplificación para funciones trigonométricas, potencias y logaritmos. Abre esta documentación si lo necesitas.

¡Te toca!¶

Pasaremos rápidamente por esta parte, para hacer cosas "más interesantes". Te proponemos algunos ejemplos para que te familiarices con el manejor de expresiones:

Crea las expresiones de la izquierda y averigua qué función te hace obtener la de la derecha:

expresión 1	expresión 2
$\left(x^{3} + 3 y + 2\right)^{2}$	$x^{6} + 6 x^{3} y + 4 x^{3} + 9 y^{2} + 12 y + 4$
$\frac{\left(3 x^{2} - 2 x + 1\right)}{\left(x - 1\right)^{2}} $	$3 + \frac{4}{x - 1} + \frac{2}{\left(x - 1\right)^{2}}$
$x^{3} + 9 x^{2} + 27 x + 27$	$\left(x + 3\right)^{3}$
$\sin(x+2y)$	$\left(2 \cos^{2}{\left (y \right )} - 1\right) \sin{\left (x \right )} + 2 \sin{\left (y \right )} \cos{\left (x \right )} \cos{\left (y \right )}$

In [32]:

#1
expr1 = (x ** 3 + 3 * y + 2) ** 2
expr1

Out[32]:

$$\left(x^{3} + 3 y + 2\right)^{2}$$

In [33]:

expr1_exp = expr1.expand()
expr1_exp

Out[33]:

$$x^{6} + 6 x^{3} y + 4 x^{3} + 9 y^{2} + 12 y + 4$$

In [34]:

#2
expr2 = (3 * x ** 2 - 2 * x + 1) / (x - 1) ** 2
expr2

Out[34]:

$$\frac{1}{\left(x - 1\right)^{2}} \left(3 x^{2} - 2 x + 1\right)$$

In [35]:

expr2.apart()

Out[35]:

$$3 + \frac{4}{x - 1} + \frac{2}{\left(x - 1\right)^{2}}$$

In [36]:

#3
expr3 = x ** 3 + 9 * x ** 2 + 27 * x + 27
expr3

Out[36]:

$$x^{3} + 9 x^{2} + 27 x + 27$$

In [37]:

expr3.factor()

Out[37]:

$$\left(x + 3\right)^{3}$$

In [38]:

#4
expr4 = sin(x + 2 * y)
expr4

Out[38]:

$$\sin{\left (x + 2 y \right )}$$

In [39]:

expand(expr4)

Out[39]:

$$\sin{\left (x + 2 y \right )}$$

In [40]:

expand_trig(expr4)

Out[40]:

$$\left(2 \cos^{2}{\left (y \right )} - 1\right) \sin{\left (x \right )} + 2 \sin{\left (y \right )} \cos{\left (x \right )} \cos{\left (y \right )}$$

In [41]:

expand(expr4, trig=True)

Out[41]:

$$2 \sin{\left (x \right )} \cos^{2}{\left (y \right )} - \sin{\left (x \right )} + 2 \sin{\left (y \right )} \cos{\left (x \right )} \cos{\left (y \right )}$$

Derivadas e integrales¶

Puedes derivar una expresion usando el método .diff() y la función dif()

In [42]:

#creamos una expresión
expr = cos(x)

#obtenemos la derivada primera con funcion
diff(expr, x)

Out[42]:

$$- \sin{\left (x \right )}$$

In [43]:

#utilizando método
expr.diff(x)

Out[43]:

$$- \sin{\left (x \right )}$$

¿derivada tercera?

In [44]:

expr.diff(x, x, x)

Out[44]:

$$\sin{\left (x \right )}$$

In [45]:

expr.diff(x, 3)

Out[45]:

$$\sin{\left (x \right )}$$

¿varias variables?

In [46]:

expr_xy = y ** 3 * sin(x) ** 2 + x ** 2 * cos(y)
expr_xy

Out[46]:

$$x^{2} \cos{\left (y \right )} + y^{3} \sin^{2}{\left (x \right )}$$

In [47]:

diff(expr_xy, x, 2, y, 2)

Out[47]:

$$2 \left(- 6 y \sin^{2}{\left (x \right )} + 6 y \cos^{2}{\left (x \right )} - \cos{\left (y \right )}\right)$$

Queremos que la deje indicada, usamos Derivative()

In [48]:

Derivative(expr_xy, x, 2, y)

Out[48]:

$$\frac{\partial^{3}}{\partial x^{2}\partial y} \left(x^{2} \cos{\left (y \right )} + y^{3} \sin^{2}{\left (x \right )}\right)$$

¿Será capaz SymPy de aplicar la regla de la cadena?

In [49]:

# Creamos una función F
F = Function('F')
F(x)

Out[49]:

$$F{\left (x \right )}$$

In [50]:

# Creamos una función G
G = Function('G')
G(x)

Out[50]:

$$G{\left (x \right )}$$

$$\frac{d}{d x} F{\left (G(x) \right )} $$

In [51]:

# Derivamos la función compuesta F(G(x))
F(G(x)).diff(x)

Out[51]:

$$\frac{d}{d x} G{\left (x \right )} \left. \frac{d}{d \xi_{1}} F{\left (\xi_{1} \right )} \right|_{\substack{ \xi_{1}=G{\left (x \right )} }}$$

En un caso en el que conocemos las funciones:

In [52]:

# definimos una f
f = 2 * y * exp(x)
f

Out[52]:

$$2 y e^{x}$$

In [53]:

# definimos una g(f)
g = f **2 * cos(x) + f
g

Out[53]:

$$4 y^{2} e^{2 x} \cos{\left (x \right )} + 2 y e^{x}$$

In [54]:

#la derivamos
diff(g,x)

Out[54]:

$$- 4 y^{2} e^{2 x} \sin{\left (x \right )} + 8 y^{2} e^{2 x} \cos{\left (x \right )} + 2 y e^{x}$$

Te toca integrar¶

Si te digo que se integra usando el método .integrate() o la función integrate(). ¿Te atreves a integrar estas casi inmediatas...?:

$$\int{\cos(x)^2}dx$$$$\int{\frac{dx}{\sin(x)}}$$$$\int{\frac{dx}{(x^2+a^2)^2}}$$

In [55]:

int1 = cos(x) ** 2
integrate(int1)

Out[55]:

$$\frac{x}{2} + \frac{1}{2} \sin{\left (x \right )} \cos{\left (x \right )}$$

In [56]:

int2 =  1 / sin(x)
integrate(int2)

Out[56]:

$$\frac{1}{2} \log{\left (\cos{\left (x \right )} - 1 \right )} - \frac{1}{2} \log{\left (\cos{\left (x \right )} + 1 \right )}$$

In [57]:

x, a = symbols('x a', real=True)

int3 = 1 / (x**2 + a**2)**2
integrate(int3, x)

Out[57]:

$$\frac{x}{2 a^{4} + 2 a^{2} x^{2}}$$

Límites¶

Calculemos este límite sacado del libro Cálculo: definiciones, teoremas y resultados, de Juan de Burgos:

$$\lim_{x \to 0} \left(\frac{x}{\tan{\left (x \right )}}\right)^{\frac{1}{x^{2}}}$$

Primero creamos la expresión:

In [58]:

x = symbols('x', real=True)
expr = (x / tan(x)) ** (1 / x**2)
expr

Out[58]:

$$\left(\frac{x}{\tan{\left (x \right )}}\right)^{\frac{1}{x^{2}}}$$

Obtenemos el límite con la función limit() y si queremos dejarlo indicado, podemos usar Limit():

In [59]:

limit(expr, x, 0)

Out[59]:

$$e^{- \frac{1}{3}}$$

Series¶

Los desarrollos en serie se pueden llevar a cabo con el método .series() o la función series()

In [60]:

#creamos la expresión
expr = exp(x)
expr

Out[60]:

$$e^{x}$$

In [61]:

#la desarrollamos en serie
series(expr)

Out[61]:

$$1 + x + \frac{x^{2}}{2} + \frac{x^{3}}{6} + \frac{x^{4}}{24} + \frac{x^{5}}{120} + \mathcal{O}\left(x^{6}\right)$$

Se puede especificar el número de términos pasándole un argumento n=.... El número que le pasemos será el primer término que desprecie.

In [62]:

# Indicando el número de términos
series(expr, n=10)

Out[62]:

$$1 + x + \frac{x^{2}}{2} + \frac{x^{3}}{6} + \frac{x^{4}}{24} + \frac{x^{5}}{120} + \frac{x^{6}}{720} + \frac{x^{7}}{5040} + \frac{x^{8}}{40320} + \frac{x^{9}}{362880} + \mathcal{O}\left(x^{10}\right)$$

Si nos molesta el $\mathcal{O}(x^{10})$ lo podemos quitar con removeO():

In [63]:

series(expr, n=10).removeO()

Out[63]:

$$\frac{x^{9}}{362880} + \frac{x^{8}}{40320} + \frac{x^{7}}{5040} + \frac{x^{6}}{720} + \frac{x^{5}}{120} + \frac{x^{4}}{24} + \frac{x^{3}}{6} + \frac{x^{2}}{2} + x + 1$$

In [64]:

series(sin(x), n=8, x0=pi/3).removeO().subs(x, x-pi/3)

Out[64]:

$$\frac{x}{2} - \frac{1}{10080} \left(x - \frac{\pi}{3}\right)^{7} - \frac{\sqrt{3}}{1440} \left(x - \frac{\pi}{3}\right)^{6} + \frac{1}{240} \left(x - \frac{\pi}{3}\right)^{5} + \frac{\sqrt{3}}{48} \left(x - \frac{\pi}{3}\right)^{4} - \frac{1}{12} \left(x - \frac{\pi}{3}\right)^{3} - \frac{\sqrt{3}}{4} \left(x - \frac{\pi}{3}\right)^{2} - \frac{\pi}{6} + \frac{\sqrt{3}}{2}$$

Resolución de ecuaciones¶

Como se ha mencionado anteriormente las ecuaciones no se pueden crear con el =

In [65]:

#creamos la ecuación
ecuacion = Eq(x ** 2 - x, 3)
ecuacion

Out[65]:

$$x^{2} - x = 3$$

In [66]:

# También la podemos crear como
Eq(x ** 2 - x -3)

Out[66]:

$$x^{2} - x - 3 = 0$$

In [67]:

#la resolvemos
solve(ecuacion)

Out[67]:

$$\begin{bmatrix}\frac{1}{2} + \frac{\sqrt{13}}{2}, & - \frac{\sqrt{13}}{2} + \frac{1}{2}\end{bmatrix}$$

Pero la gracia es resolver con símbolos, ¿no? $$a e^{\frac{x}{t}} = C$$

In [68]:

# Creamos los símbolos y la ecuación
a, x, t, C = symbols('a, x, t, C', real=True)
ecuacion = Eq(a * exp(x/t), C)
ecuacion

Out[68]:

$$a e^{\frac{x}{t}} = C$$

In [69]:

# La resolvemos
solve(ecuacion ,x)

Out[69]:

$$\begin{bmatrix}t \log{\left (\frac{C}{a} \right )}\end{bmatrix}$$

Si consultamos la ayuda, vemos que las posibilidades y el número de parámetros son muchos, no vamos a entrar ahora en ellos, pero ¿se ve la potencia?

Ecuaciones diferenciales¶

Tratemos de resolver, por ejemplo:

$$y{\left (x \right )} + \frac{d}{d x} y{\left (x \right )} + \frac{d^{2}}{d x^{2}} y{\left (x \right )} = \cos{\left (x \right )}$$

In [70]:

x = symbols('x')
f = Function('y')
ecuacion_dif = Eq(y(x).diff(x,2) + y(x).diff(x) + y(x), cos(x))
ecuacion_dif

Out[70]:

$$y{\left (x \right )} + \frac{d}{d x} y{\left (x \right )} + \frac{d^{2}}{d x^{2}} y{\left (x \right )} = \cos{\left (x \right )}$$

In [71]:

#resolvemos
dsolve(ecuacion_dif, f(x))

Out[71]:

$$y{\left (x \right )} = \frac{1}{\sqrt{e^{x}}} \left(C_{1} \sin{\left (\frac{\sqrt{3} x}{2} \right )} + C_{2} \cos{\left (\frac{\sqrt{3} x}{2} \right )}\right) + \sin{\left (x \right )}$$

Matrices¶

In [72]:

#creamos una matriz llena de símbolos
a, b, c, d = symbols('a b c d')
A = Matrix([
            [a, b],
            [c, d]
    ])
A

Out[72]:

$$\left[\begin{matrix}a & b\\c & d\end{matrix}\right]$$

In [73]:

#sacamos autovalores
A.eigenvals()

Out[73]:

$$\begin{Bmatrix}\frac{a}{2} + \frac{d}{2} - \frac{1}{2} \sqrt{a^{2} - 2 a d + 4 b c + d^{2}} : 1, & \frac{a}{2} + \frac{d}{2} + \frac{1}{2} \sqrt{a^{2} - 2 a d + 4 b c + d^{2}} : 1\end{Bmatrix}$$

In [74]:

#inversa
A.inv()

Out[74]:

$$\left[\begin{matrix}\frac{1}{a} + \frac{b c}{a^{2} \left(d - \frac{b c}{a}\right)} & - \frac{b}{a \left(d - \frac{b c}{a}\right)}\\- \frac{c}{a \left(d - \frac{b c}{a}\right)} & \frac{1}{d - \frac{b c}{a}}\end{matrix}\right]$$

In [75]:

#elevamos al cuadrado la matriz
A ** 2

Out[75]:

$$\left[\begin{matrix}a^{2} + b c & a b + b d\\a c + c d & b c + d^{2}\end{matrix}\right]$$

_ Esto ha sido un rápido recorrido por algunas de las posibilidades que ofrece SymPy . El cálculo simbólico es un terreno díficil y este joven paquete avanza a pasos agigantados gracias a un grupo de desarrolladores siempre dispuestos a mejorar y escuchar sugerencias. Sus posibilidades no acaban aquí. En la siguiente clase presentaremos el paquete mechanics, pero además cuenta con herramientas para geometría, mecánica cuántica, teoría de números, combinatoria... Puedes echar un ojo aquí. _

Si te ha gustado esta clase:

¡Síguenos en Twitter!

Follow @Pybonacci <a href="https://twitter.com/AlexS12" class="twitter-follow-button" data-show-count="false" align="right";>Follow @AlexS12 <a href="https://twitter.com/newlawrence" class="twitter-follow-button" data-show-count="false" align="right";>Follow @newlawrence

CLASE 5a: SymPy¶

Documentación & SymPy Live Shell¶

SymPy Gamma¶

Creación de símbolos¶

Expresiones¶

simplify()¶

.subs()¶

Simplificación¶

¡Te toca!¶

Derivadas e integrales¶

Te toca integrar¶

Límites¶

Series¶

Resolución de ecuaciones¶

Ecuaciones diferenciales¶

Matrices¶

¡Síguenos en Twitter!

Curso AeroPython por Juan Luis Cano Rodriguez y Alejandro Sáez Mollejo se distribuye bajo una Licencia Creative Commons Atribución 4.0 Internacional.¶

`simplify()`¶

`.subs()`¶