59  Cálculo simbólico con SymPy

Fecha: 10 y 12 de noviembre, 2026

Invitado: Dr. José Crispin Ruiz Pantaleon (JCRP)

SymPy es una biblioteca de Python para hacer matemáticas simbólicas: manipula expresiones con variables como lo haríamos con lápiz y papel, y obtiene resultados exactos en lugar de aproximaciones decimales.

En este capítulo veremos:

• La diferencia entre el cálculo numérico (NumPy) y el simbólico (SymPy).
• Cómo definir símbolos y manipular expresiones: sustituir, simplificar, factorizar y expandir.
• Cómo resolver ecuaciones y sistemas de ecuaciones.
• Cálculo: límites, derivadas, integrales, series y sumas.
• Matrices y ecuaciones diferenciales.
• Cómo convertir resultados a LaTeX para tus documentos de Quarto.

59.1 Cálculo numérico vs. cálculo simbólico

import numpy as np
import sympy as sp

print("NumPy: ", np.sqrt(8))   # Aproximación decimal

NumPy:  2.8284271247461903

print("SymPy: ", sp.sqrt(8))   # Resultado exacto

SymPy:  2*sqrt(2)

	Cálculo numérico (NumPy)	Cálculo simbólico (SymPy)
Trabaja con	Números decimales.	Símbolos y expresiones exactas.
\(\sqrt{8}\)	2.8284271247461903	2*sqrt(2)
\(\frac{1}{3} + \frac{1}{6}\)	0.5 (con posibles errores de redondeo)	1/2 (fracción exacta)
Derivar \(x^2\)	Aproxima la pendiente en puntos.	Obtiene \(2x\).
Ideal para	Datos, simulaciones, grandes cantidades de números.	Álgebra, cálculo, demostrar y verificar fórmulas.

print(0.1 + 0.2)                                  # Error de redondeo de los decimales

0.30000000000000004

print(sp.Rational(1, 10) + sp.Rational(2, 10))    # Fracciones exactas

3/10

print(sp.pi.evalf(30))                            # Pi con 30 dígitos

3.14159265358979323846264338328

59.2 Instalar e importar SymPy

En tu computadora, instálalo en tu entorno virtual (consulta el apéndice Entornos para trabajar con Python):

python -m pip install sympy

En Google Colab y en los editores de esta página ya está instalado. Por convención se importa con el alias sp:

import sympy as sp
print(sp.__version__)

1.14.0

59.3 Símbolos y expresiones

En Python, una variable debe tener un valor antes de usarse. En SymPy primero declaramos qué nombres son símbolos matemáticos con sp.symbols():

x, y = sp.symbols("x y")

expresion = x**2 + 2*x*y + y**2
print(expresion)

x**2 + 2*x*y + y**2

print(type(x))

<class 'sympy.core.symbol.Symbol'>

Importante

Para elevar a una potencia se usa **, no ^. En Python, x^2 es otra operación (el “o exclusivo” entre bits) y produce un error con los símbolos.

Podemos indicar supuestos sobre los símbolos, que SymPy usa para simplificar:

t = sp.symbols("t")
r = sp.symbols("r", positive=True)

print(sp.sqrt(t**2))   # No sabe si t es negativo

sqrt(t**2)

print(sp.sqrt(r**2))   # Sabe que r es positivo

r

59.3.1 Sustituir y evaluar

Método	¿Qué hace?	Ejemplo
expr.subs(x, valor)	Sustituye un símbolo por un valor u otra expresión.	(x**2).subs(x, 3) → 9
expr.evalf() o sp.N(expr)	Calcula el valor decimal.	sp.sqrt(2).evalf() → 1.414...
sp.N(expr, n)	Valor decimal con n dígitos.	sp.N(sp.pi, 50)

f = x**2 + 3*x - 1

print(f.subs(x, 2))                  # Evaluar en x = 2

9

print(f.subs(x, y + 1))              # Sustituir por otra expresión

3*y + (y + 1)**2 + 2

print(f.subs(x, sp.sqrt(2)))         # Resultado exacto

1 + 3*sqrt(2)

print(f.subs(x, sp.sqrt(2)).evalf()) # Resultado decimal

5.24264068711928

59.3.2 Fracciones exactas

Si escribes 1/2 dentro de una expresión, Python calcula primero el decimal 0.5. Para obtener una fracción exacta usa sp.Rational():

print(x + 1/2)                 # 0.5: decimal

x + 0.5

print(x + sp.Rational(1, 2))   # 1/2: fracción exacta

x + 1/2

59.4 Simplificar, expandir y factorizar

Función	¿Qué hace?	Ejemplo	Resultado
sp.expand()	Desarrolla productos y potencias.	sp.expand((x + 1)**2)	x**2 + 2*x + 1
sp.factor()	Factoriza.	sp.factor(x**2 - 1)	(x - 1)*(x + 1)
sp.simplify()	Intenta encontrar la forma más simple.	sp.simplify(sp.sin(x)**2 + sp.cos(x)**2)	1
sp.cancel()	Simplifica fracciones de polinomios.	sp.cancel((x**2 - 1)/(x - 1))	x + 1
sp.apart()	Fracciones parciales.	sp.apart(1/(x**2 - 1))	-1/(2*(x + 1)) + 1/(2*(x - 1))
sp.collect()	Agrupa términos por un símbolo.	sp.collect(x*y + x, x)	x*(y + 1)
sp.trigsimp()	Simplifica expresiones trigonométricas.	sp.trigsimp(2*sp.sin(x)*sp.cos(x))	sin(2*x)

print(sp.expand((x + y)**3))

x**3 + 3*x**2*y + 3*x*y**2 + y**3

print(sp.factor(x**3 - x**2 - x + 1))

(x - 1)**2*(x + 1)

print(sp.cancel((x**2 - 1) / (x - 1)))

x + 1

print(sp.apart(1 / (x**2 - 1)))

-1/(2*(x + 1)) + 1/(2*(x - 1))

59.5 Ecuaciones

Para escribir una ecuación no se usa = (que asigna) ni == (que compara), sino sp.Eq(lado_izquierdo, lado_derecho). Para resolverla usamos sp.solve():

ecuacion = sp.Eq(x**2 - 5*x, -6)

print(ecuacion)

Eq(x**2 - 5*x, -6)

print(sp.solve(ecuacion, x))

[2, 3]

Si pasamos solo una expresión, sp.solve() supone que está igualada a cero:

print(sp.solve(x**2 - 4, x))

[-2, 2]

print(sp.solve(x**2 + 1, x))   # Soluciones complejas: I es la unidad imaginaria

[-I, I]

Con sp.solveset() podemos indicar en qué conjunto buscar soluciones:

print(sp.solveset(x**2 + 1, x, domain=sp.S.Reals))   # Sin soluciones reales

EmptySet

59.5.1 La fórmula general

SymPy puede resolver la ecuación cuadrática con coeficientes simbólicos y obtener la fórmula general que usamos en el capítulo de condicionales:

a, b, c = sp.symbols("a b c")
print(sp.solve(a*x**2 + b*x + c, x))

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

59.5.2 Sistemas de ecuaciones

Es el mismo sistema que resolvimos con np.linalg.solve() en el capítulo de NumPy, pero ahora con solución exacta:

x, y, z = sp.symbols("x y z")

sistema = [sp.Eq(2*x + y - z, 3),
           sp.Eq(x - y + 2*z, 5),
           sp.Eq(3*x + 2*y + z, 10)]

print(sp.solve(sistema, [x, y, z]))

{x: 2, y: 1, z: 2}

59.6 Cálculo

Operación	Notación	SymPy
Límite	\(\lim_{x \to a} f(x)\)	sp.limit(f, x, a)
Límite al infinito	\(\lim_{x \to \infty} f(x)\)	sp.limit(f, x, sp.oo)
Derivada	\(\frac{df}{dx}\)	sp.diff(f, x)
Derivada de orden \(n\)	\(\frac{d^n f}{dx^n}\)	sp.diff(f, x, n)
Derivada parcial	\(\frac{\partial f}{\partial y}\)	sp.diff(f, y)
Integral indefinida	\(\int f\,dx\)	sp.integrate(f, x)
Integral definida	\(\int_a^b f\,dx\)	sp.integrate(f, (x, a, b))
Serie de Taylor	\(f(x) \approx \sum \frac{f^{(n)}(a)}{n!}(x-a)^n\)	sp.series(f, x, a, n)
Suma	\(\sum_{i=1}^{n} f(i)\)	sp.summation(f, (i, 1, n))

El infinito se escribe sp.oo (dos letras “o”), y las constantes son sp.pi, sp.E y sp.I.

59.6.1 Límites

x, n = sp.symbols("x n")

print(sp.limit(sp.sin(x) / x, x, 0))

1

print(sp.limit((1 + 1/n)**n, n, sp.oo))   # El número e

E

print(sp.limit(1/x, x, 0, dir="+"))       # Límite por la derecha

oo

59.6.2 Derivadas

f = x**2 * sp.sin(x)

print(sp.diff(f, x))         # Primera derivada

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

print(sp.diff(f, x, 2))      # Segunda derivada

-x**2*sin(x) + 4*x*cos(x) + 2*sin(x)

x, y = sp.symbols("x y")
g = x**2 * y + sp.sin(y)
print(sp.diff(g, y))         # Derivada parcial respecto a y

x**2 + cos(y)

59.6.3 Integrales

print(sp.integrate(x * sp.exp(x), x))                    # Indefinida (sin la constante C)

(x - 1)*exp(x)

print(sp.integrate(x * sp.exp(x), (x, 0, 1)))            # Definida

1

print(sp.integrate(sp.exp(-x**2), (x, -sp.oo, sp.oo)))   # Integral de Gauss

sqrt(pi)

59.6.4 Series de Taylor

print(sp.series(sp.sin(x), x, 0, 8))

x - x**3/6 + x**5/120 - x**7/5040 + O(x**8)

print(sp.series(sp.exp(x), x, 0, 5).removeO())   # Sin el término O(x^n)

x**4/24 + x**3/6 + x**2/2 + x + 1

59.6.5 Sumas

Para sumar los números del 1 al 100 podríamos sumarlos uno por uno; en el capítulo de ciclos for aprenderemos a hacerlo con un ciclo. SymPy, en cambio, obtiene directamente la fórmula para cualquier \(n\):

i, n = sp.symbols("i n", integer=True, positive=True)

print(sp.factor(sp.summation(i, (i, 1, n))))      # 1 + 2 + ... + n

n*(n + 1)/2

print(sp.summation(i, (i, 1, 100)))               # Con n = 100

5050

59.7 Matrices

sp.Matrix trabaja con valores exactos, incluidas fracciones y símbolos:

M = sp.Matrix([[2, 1],
               [1, 3]])

print("det(M) =", M.det())

det(M) = 5

print("M^-1 =", M.inv())

M^-1 = Matrix([[3/5, -1/5], [-1/5, 2/5]])

print("Valores propios:", M.eigenvals())

Valores propios: {5/2 - sqrt(5)/2: 1, sqrt(5)/2 + 5/2: 1}

Operación	SymPy
Producto	A * B (en SymPy, * sí es el producto de matrices)
Transpuesta	A.T
Determinante	A.det()
Inversa	A.inv()
Forma escalonada reducida	A.rref()
Valores y vectores propios	A.eigenvals(), A.eigenvects()

59.8 Ecuaciones diferenciales

Para resolver una ecuación diferencial ordinaria declaramos una función desconocida con sp.Function() y usamos sp.dsolve(). Por ejemplo, \(y'' + y = 0\):

y = sp.Function("y")
x = sp.symbols("x")

edo = sp.Eq(y(x).diff(x, 2) + y(x), 0)
print(sp.dsolve(edo, y(x)))

Eq(y(x), C1*sin(x) + C2*cos(x))

C1 y C2 son las constantes de integración.

59.9 De SymPy a NumPy: lambdify()

Para evaluar una expresión en muchos puntos (por ejemplo, para graficarla), conviértela en una función de NumPy con sp.lambdify():

expr = sp.sin(x)**2
f_numpy = sp.lambdify(x, expr, "numpy")

valores = np.linspace(0, np.pi, 5)
print(f_numpy(valores).round(4))

[0.  0.5 1.  0.5 0. ]

Usaremos esta herramienta en el capítulo de visualización de datos.

59.10 De SymPy a LaTeX

sp.latex() convierte cualquier expresión al código LaTeX que usamos en Markdown y Quarto:

integral = sp.Integral(sp.exp(-x**2), (x, -sp.oo, sp.oo))
print(sp.latex(integral), "=", sp.latex(integral.doit()))

\int\limits_{-\infty}^{\infty} e^{- x^{2}}\, dx = \sqrt{\pi}

En un documento de Quarto puedes mostrar el resultado como fórmula con la opción #| output: asis:

print(f"$${sp.latex(integral)} = {sp.latex(integral.doit())}$$")

\[\int\limits_{-\infty}^{\infty} e^{- x^{2}}\, dx = \sqrt{\pi}\]

59.11 Errores comunes

Error	Ejemplo	Corrección
Usar un símbolo sin declararlo	sp.diff(x**2, x) sin sp.symbols (NameError)	x = sp.symbols("x")
Usar ^ para potencias	x^2 (TypeError)	x**2
Fracciones con / entre enteros	x + 1/3 da x + 0.333...	x + sp.Rational(1, 3)
Usar funciones de NumPy con símbolos	np.sin(x) (TypeError)	sp.sin(x)
Escribir ecuaciones con =	x**2 = 4 (SyntaxError)	sp.Eq(x**2, 4)
Olvidar el símbolo en solve o diff	sp.diff(x*y) (ambiguo)	sp.diff(x*y, x)

59.12 Ejercicios

Resuelve los ejercicios en los editores interactivos. Recuerda que lo que escribas no se guarda al recargar la página. La primera vez que ejecutes un editor, SymPy puede tardar unos segundos en cargar.

59.12.1 Ejercicio 1. Factorizar y expandir

1. Expande \((x + 1)^3\).
2. Factoriza \(x^2 - 5x + 6\) y deja el resultado en la última línea.

NotaPista

Usa sp.expand((x + 1)**3) y sp.factor(x**2 - 5*x + 6). Recuerda que la potencia se escribe con **.

TipSolución

import sympy as sp
x = sp.symbols("x")

print(sp.expand((x + 1)**3))   # x**3 + 3*x**2 + 3*x + 1

sp.factor(x**2 - 5*x + 6)      # (x - 3)*(x - 2)

59.12.2 Ejercicio 2. Raíces de un polinomio

Encuentra todas las raíces del polinomio \(p(x) = x^3 - 6x^2 + 11x - 6\). Después verifica que cada raíz cumple \(p(r) = 0\) sustituyéndola con subs(). Deja la lista de raíces en la última línea.

NotaPista

sp.solve(p, x) regresa una lista con las raíces. Para verificar: p.subs(x, 1).

TipSolución

import sympy as sp
x = sp.symbols("x")

p = x**3 - 6*x**2 + 11*x - 6
raices = sp.solve(p, x)

print("Raíces:", raices)     # [1, 2, 3]
print("p(1) =", p.subs(x, 1)) # 0

raices

59.12.3 Ejercicio 3. El número \(e\) como límite

Calcula \(\displaystyle\lim_{n \to \infty} \left(1 + \frac{1}{n}\right)^n\) y muestra su valor con 20 dígitos. Deja el límite en la última línea.

NotaPista

El infinito es sp.oo. Para ver el valor con 20 dígitos usa sp.N(limite, 20).

TipSolución

import sympy as sp
n = sp.symbols("n")

limite = sp.limit((1 + 1/n)**n, n, sp.oo)

print("Límite:", limite)                    # E
print("Valor decimal:", sp.N(limite, 20))   # 2.7182818284590452354
limite

59.12.4 Ejercicio 4. Derivadas

Para \(f(x) = x^2 \sin(x)\):

1. Calcula la primera y la segunda derivada.
2. Evalúa la primera derivada en \(x = \pi\) y deja el resultado en la última línea.

NotaPista

sp.diff(f, x) y sp.diff(f, x, 2). Para evaluar: primera.subs(x, sp.pi). Usa sp.pi, no 3.1416, para obtener el resultado exacto.

TipSolución

import sympy as sp
x = sp.symbols("x")

f = x**2 * sp.sin(x)

primera = sp.diff(f, x)       # x**2*cos(x) + 2*x*sin(x)
segunda = sp.diff(f, x, 2)    # -x**2*sin(x) + 4*x*cos(x) + 2*sin(x)

print("f'(x) =", primera)
print("f''(x) =", segunda)

primera.subs(x, sp.pi)        # -pi**2

59.12.5 Ejercicio 5. Integrales

1. Calcula la integral indefinida \(\int x e^x \, dx\).
2. Verifica tu resultado derivándolo: debes obtener de nuevo \(x e^x\).
3. Calcula la integral definida \(\int_0^1 x e^x \, dx\) y déjala en la última línea.

NotaPista

sp.integrate(x * sp.exp(x), x) para la indefinida y sp.integrate(x * sp.exp(x), (x, 0, 1)) para la definida. Para verificar, deriva: sp.diff(indefinida, x).

TipSolución

import sympy as sp
x = sp.symbols("x")

indefinida = sp.integrate(x * sp.exp(x), x)
print("Integral indefinida:", indefinida)                     # (x - 1)*exp(x)
print("Verificación:", sp.simplify(sp.diff(indefinida, x)))   # x*exp(x)

sp.integrate(x * sp.exp(x), (x, 0, 1))                        # 1

59.12.6 Ejercicio 6. Fórmula de la suma de cuadrados

Usa sp.summation() para obtener una fórmula de \(\displaystyle\sum_{i=1}^{n} i^2 = 1^2 + 2^2 + \cdots + n^2\). Factoriza el resultado y déjalo en la última línea. Después úsala para calcular la suma con \(n = 100\).

NotaPista

sp.factor(sp.summation(i**2, (i, 1, n))). Para evaluar la fórmula: formula.subs(n, 100).

TipSolución

import sympy as sp
i, n = sp.symbols("i n", integer=True, positive=True)

formula = sp.factor(sp.summation(i**2, (i, 1, n)))
print("Fórmula:", formula)                    # n*(n + 1)*(2*n + 1)/6
print("Con n = 100:", formula.subs(n, 100))   # 338350

formula

59.13 Reto final: dos series famosas

La serie de Leibniz permite aproximar \(\pi\); en el capítulo de ciclos for la usaremos para aproximarlo sumando 1000 de sus términos. Con SymPy podemos calcular la suma infinita de forma exacta:

1. Calcula \(\displaystyle\sum_{k=0}^{\infty} \frac{(-1)^k}{2k + 1}\) (serie de Leibniz).
2. Calcula \(\displaystyle\sum_{k=1}^{\infty} \frac{1}{k^2}\) (el problema de Basilea, resuelto por Euler en 1734) y déjala en la última línea.

NotaPista

sp.summation((-1)**k / (2*k + 1), (k, 0, sp.oo)) y sp.summation(1 / k**2, (k, 1, sp.oo)).

TipSolución

import sympy as sp
k = sp.symbols("k", integer=True)

leibniz = sp.summation((-1)**k / (2*k + 1), (k, 0, sp.oo))
print("Serie de Leibniz:", leibniz)                          # pi/4

basilea = sp.summation(1 / k**2, (k, 1, sp.oo))
print("Problema de Basilea:", basilea, "≈", sp.N(basilea))   # pi**2/6 ≈ 1.64493406684823
basilea

Por eso, para aproximar \(\pi\) con esta serie, hay que multiplicar la suma por 4: \(4 \cdot \frac{\pi}{4} = \pi\). Lo haremos en el capítulo de ciclos for.

59.14 Resumen

Tema	SymPy
Símbolos	x, y = sp.symbols("x y")
Sustituir / evaluar	expr.subs(x, 2), expr.evalf(), sp.N(expr, n)
Fracciones y constantes	sp.Rational(1, 3), sp.pi, sp.E, sp.I, sp.oo
Manipular expresiones	expand, factor, simplify, cancel, apart, trigsimp
Ecuaciones	sp.Eq(izq, der), sp.solve(), sp.solveset()
Cálculo	limit, diff, integrate, series, summation
Matrices	sp.Matrix, .det(), .inv(), .eigenvals()
Ecuaciones diferenciales	sp.Function, sp.dsolve()
Conversión	sp.lambdify() (a NumPy), sp.latex() (a LaTeX)

59.15 Referencias

• Documentación de SymPy
• Tutorial introductorio de SymPy
• Scientific Python Lectures: SymPy