Une très courte introduction à SymPy (Python)
Vincent Jalby
Université de Limoges
25 août 2023
1 Introduction mais l’affichage restera basique, et nécessitera d’utiliser
la fonction print() pour afficher le moindre résultat.
Le langage de programmation Python est particulière-
ment adapté aux mathématiques. De nombreuses biblio- La méthode moderne adaptée aux sciences est l’utilisa-
thèques (libraries) permettent d’étendre ses fonctions à de tion d’un notebook de type Jupyter. On retrouve l’inter-
nombreux domaines d’application. Nous nous intéres- activité du mode console, avec un affichage amélioré à
sons ici plus particulièrement à Sympy permettant d’ef- travers un navigateur web. L’enregistrement est facilité,
fectuer des calculs formels (symboliques) « comme en toujours à travers une interface web. C’est cette méthode
cours ». que nous allons privilégier par la suite.
D’autres bibliothèques seront utilisées explicitement ou
implicitement : matplotlib, jupyterlab.
4 Utilisation de Jupyter
2 Installation Si vous avez tenté une installation standard de Python
et de ses bibliothèques, vous lancez Jupyter (et donc Py-
Les plus aguerris pourront installer Python 3.x en le télé-
thon) à l’aide de l’instruction jupyter-lab. Mais peu
chargeant directement depuis le site [Link] puis en
de chance que cela fonctionne au premier coup !
installant les bibliothèques citées précédemment (il fau-
dra pour cela utiliser la commande pip ou pip3). Si vous avez installé la distribution de Anaconda, lan-
Mais il est sans aucun doute beaucoup plus facile d’ins- cez l’application Anaconda-Navigator puis cliquez sur
taller la distribution Anaconda à partir du site ana- le bouton Launch de la case JupyterLab.
[Link]. Elle permet d’installer Python ainsi que
Si vous utilisez la version cloud de anaconda,
toutes les bibliothèques nécessaires (et beaucoup plus
connectez-vous sur [Link] puis cliquez
encore) en un simple clic. Inconvénient, elle va utiliser
sur le menu Notebooks.
environ 7 gigaoctets d’espace sur votre disque dur.
La dernière méthode consiste à ne rien installer, mais à Quelque soit la méthode utilisée, vous obtiendrez dans
utiliser Anaconda dans le cloud (en étant donc connecté votre navigateur un écran proche du suivant :
à Internet). L’utilisation d’un compte gratuit sera large-
ment suffisante à notre niveau. Pour cela, il suffit de rem-
plir le formulaire Sign Up sur le site [Link].
3 Utilisation
Il existe (au moins) trois façons d’utiliser python.
Le mode console, par exemple avec l’application IDLE
inclus dans l’installation standard de Python, consiste à
taper des instructions et à obtenir immédiatement le ré-
sultat correspondant. L’affichage des résultats est basique
La partie de gauche permet de naviguer et de gérer les
(ASCII) et il y a peu de possibilité d’enregistrer son travail
fichiers (notebooks) Python créés. Celle de droite affiche
si ce n’est le copier-coller.
les notebooks.
On peut aussi créer un programme complet (fichier .py) Pour créer votre premier notebook, il suffit de cliquer sur
puis l’executer avec Python, par exemple avec IDLE, le premier icône de la liste dans la zone Notebooks.
5 Notions de base 7 Calcul dans ℝ
Un notebook est composé d’une suite de lignes de com- On peut utiliser les opérations standards de ℝ : somme
mandes (dans la suite en gris) et de résultats (dans la suite (+), différence (-), produit (*), puissance (**) :
en bleu) : 1+2*3+2**4
[1]: 1+1 23
[1]: 2 La division (/) donne la valeur numérique du résultat :
Chaque ligne de commande peut contenir une ou plu- 3/9
sieurs instructions. Dans ce cas, il faudra les séparer par 0.3333333333333333
des points-virgules (;). Pour exécuter une ligne de com- Pour obtenir la version rationnelle, on utilise la fonction
mande, il suffit de placer le curseur dans la ligne (inutile Rational() :
de se positionner à la fin de la ligne !) contenant la ou Rational(3,9)
les commandes et d’appuyer sur les touches Majuscule
1
et Retour ou Contrôle et Entrée.
3
Pour rappeler le dernier résultat (c’est-à-dire, le dernier
calcul effectué par Python lors de la session courante),
on utilise l’instruction ( _ ) (tiret bas ou underscore) : 8 Fonctions usuelles
_ + 3 SymPy définit les fonctions usuelles avec les notations
5 mathématiques standards : sin(), cos(), tan(),
On peut enregistrer une valeur (nombre) ou une expres- ln (), exp(), etc.
sion (fonction, équation, etc) à l’aide de l’opérateur d’af- sin(pi/3)
fectation (=) : √3
x0 = 12 2
puis l’utiliser comme un symbole mathématique : La constante π = 3.14 … s’obtient avec l’instruction pi.
x0 + x0 + 3 De même, e = 2.71 … s’obtient avec l’instruction E.
27 Racine carrée et valeur absolue sont obtenues via les
Attention de n’utiliser que des caractères alphanumé- fonctions sqrt() (SQuare RooT) et abs() :
riques ASCII (A-Z, 0-9) dans les noms d’affectation. sqrt(8)
2√2
6 SymPy On peut forcer l’évaluation numérique (floating-point)
d’un résultat avec la propriété evalf() :
La bibliothèque Sympy permet d’effectuer de nombreux
sqrt(8).evalf()
calculs formels dans Python tels que calculs de dérivées,
de limites, d’intégrales ou encore résolutions d’équa- 2.82842712474619
tions. La factorielle (!) est obtenue avec la fonction
factorial() et les coefficients binomiaux (combi-
Il faut d’abord commencer par importer la bibliothèque
naisons) avec binomial(n,k).
à l’aide de l’instruction suivante :
from sympy import *
On déclare ensuite les variables (au sens mathématique) 9 Simplification
que l’on souhaite utiliser, par exemple x et y : Les instructions expand(), factor(),
x, y = symbols('x y') simplify() permettent de développer, factoriser,
Il est alors possible d’effectuer des calculs avec ces va- simplifier des expressions mathématiques :
riables, par exemple : expand((x-1)**2)
3*x - x + 1 + x*y/x x2 − 2x + 1
2x + y + 1 factor(x**2-1)
Lors de la déclaration des variables, il est parfois utile de (x − 1)(x + 1)
préciser leur type (réel ou entier) :
simplify(x**2/(x+x**3))
x, y = symbols('x y', real=True)
x
n = symbols('n', integer=True) 2
x +1
2
10 Résolution d’équations limit(f(x),x,oo)
∞
On utilise l’instruction solve() pour trouver les solu-
tions exactes d’une équation, en indiquant l’expression Noter l’utilisation de « oo » (deux « o » minuscules) pour
devant s’annuler sans le « = 0 » (ici x2 + x − 2) et la représenter le symbole ∞.
variable (ici x) : Les limites à droite et à gauche s’obtiennent en rajoutant
solve(x**2+x-2, x) un argument '+' ou '-' à la limites :
[−2, 1] limit(1/x,x,0,'-')
Cela fonctionne de même pour les systèmes d’équations, −∞
en regroupant les équations (et variables) entre crochets : L’instruction diff() permet de calculer la dérivée
solve([x+2*y-7,x-y-1],[x,y]) d’une fonction :
{x ∶ 3, y ∶ 2} diff(f(x),x)
Pour obtenir les racines d’un polynôme avec leur ordre 4x
de multiplicité, on utilise roots() : Les dérivées successives s’obtiennent en précisant l’ordre
roots(x**3-3*x**2-9*x+27,x) la dérivée :
{−3 ∶ 1, 3 ∶ 2} diff(f(x),x,2)
4
Lorsqu’il n’est pas possible de résoudre formellement
l’équation (ou le système), on peut utiliser la fonction On peut tester la convexité d’une fonction sur ℝ ou un
nsolve() pour trouver une solution numérique dans intervalle avec is_convex() :
un intervalle (ici, [0, 1]) : is_convex(x**2,x)
nsolve(cos(x)-x,x,(0,1)) True
0.739085133215161
is_convex(x**3,x)
False
11 Inégalités
is_convex(x**3,x,
On procède de même pour résoudre les inéquations : domain = Interval(0,oo))
reduce_inequalities(x**2-1 >= 0, x) True
1 ⩽ x ∨ x ⩽ −1 L’instruction integrate() permet de calculer des pri-
La lecture du résultat est moins aisée. Le symbole ∨ si- mitives ou des intégrales :
gnifie « ou » (et ∧ signifie « et »). Le résultat est donc integrate(f(x),x)
1 ⩽ x ou x ⩽ −1. Soit x ∈ ]−∞, −1] ∪ [1, +∞[.
2x3
− 3x
3
12 Etude d’une fonction
integrate(f(x),(x,0,1))
On peut facilement définir une fonction (mathématique)
7
comme une fonction Python standard : −
3
def f(x):
return 2*x**2-3 Pour obtenir le développement limité d’une fonction
f(x) à l’ordre n au voisinage de x0 , on utilise l’instruc-
On peut ensuite utiliser les notations mathématiques
tion series(f(x),x,x0,n+1) :
standards :
series(exp(x),x,0,4)
f(3)
15 x2 x3
1+x+ + + O(x4 )
2 6
f(x+1)
2(x + 1)2 − 3
Les limites sont obtenues avec l’instruction limit : 13 Représentations graphiques
limit(f(x),x,0) L’instruction plot() permet d’obtenir une représenta-
−3 tion graphique d’une fonction sur un intervalle :
3
plot(f(x),(x,-3,5)) 14 Fonctions de plusieurs variables
Comme pour les fonctions d’une variable, on utilise une
fonction Python pour définir une fonction mathématique
40
de plusieurs variables :
30 def f(x,y):
return x*y-x**2
f(x)
20
On utilise alors la fonction avec les notations naturelles :
10 f(2,3)
2
0
3 2 1 0 1 2 3 4 5 f(x,2*x)
x
On peut aussi préciser les valeurs extrêmes de l’axe y et x2
même la couleur de la courbe : Les dérivées partielles premières et secondes s’ob-
plot(exp(x),(x,-10,5),ylim=(0,20), tiennent avec :
line_color='red',ylabel='Exp(x)',
title="Graphique") diff(f(x,y),x)
Graphique −2x + y
20.0
17.5 diff(f(x,y),y)
x
15.0
12.5 diff(f(x,y),x,2)
Exp(x)
10.0 −2
7.5
diff(f(x,y),x,y)
5.0 1
2.5 Le gradient s’obtient de manière indirecte :
0.0 Matrix(derive_by_array(f(x,y),[x,y]))
10 8 6 4 2 0 2 4
x
−2x + y
Il est possible de représenter plusieurs fonctions sur le [ ]
x
même graphique :
plot(exp(x),x**2,2*x,(x,-1,3), C’est plus simple pour la matrice hessienne et le hessien :
legend=True) hessian(f(x,y),(x,y))
ex 20 −2 1
x2 [ ]
2x 1 0
15
det( hessian(f(x,y),(x,y)) )
10 −1
f(x)
Il est possible de représenter f graphiquement en 3 di-
5 mensions. Pour cela, il est nécessaire d’importer la fonc-
tion plot3d() :
0 from [Link] import plot3d
1.0 0.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0
x
plot3d(f(x,y),(x,-1,1),(y,-1,1))
4
15 Courbes paramétrées
L’instruction plot_parametric permet de représen-
0.0 ter facilement des courbes paramétrées :
0.5
f(x, y)
t = symbols('t')
1.0 plot_parametric((t*cos(t),t*sin(t)),
1.5 (t,0,6*pi), aspect_ratio=(1,1))
2.0
15
1.00
0.75
0.50
1.000.75 0.25 10
0.00
0.500.25 0.25 y
0.000.25 0.50
x 0.500.75 0.75
1.00 1.00 5
On peut également tracer les courbes de niveau :
0
plot_implicit(f(x,y)-1, (x,-9,9), 15 10 5 0 5 10 15 20
(y,-9,9) )
5
8
6 10
4
2 15
0
y
8 6 4 2 0 2 4 6 8
2x
4
6 16 Suites et séries
8
Une suite (un )n se définit de la même façon qu’une
Pour représenter plusieurs lignes de niveau, il faut tra- fonction :
vailler un peu plus : def u(n):
plot1 = plot_implicit(f(x,y)-1, return 1/n**2
(x,-9,9), (y,-9,9), show=False)
On peut alors calculer les termes successifs de la suite et
plot2 = plot_implicit(f(x,y)-2,
déterminer sa limite :
(x,-9,9),(y,-9,9),
line_color='red', show=False) u(2)
plot3 = plot_implicit(f(x,y)-3, 0.25
(x,-9,9),(y,-9,9),
line_color='green', show=False) limit(u(n),n=oo)
[Link](plot2[0]) 0
[Link](plot3[0])
[Link]()
Pour représenter graphiquement une suite (les couples
8 (n, un )), il est nécessaire de charger la bibliothèque
6 matplotlib :
4 import [Link] as plt
2
On utilise ensuite une des deux méthodes suivantes :
0
y
8 6 4 2 0 2 4 6 8 for k in range(1,10):
2x [Link](n,u(k),'r.')
4
ou
6
[Link](range(1,10),
8 [u(k) for k in range(1,10)],'r.')
5
1.0
La condition initiale x(0) = 1 est précisée via l’option
ics :
0.8
dsolve(diff(x(t),t)-t*x(t)-t,x(t),
ics={x(0):1})
t2
0.6 x(t) = 2e 2 − 1
On résout de même les équations différentielles
0.4 d’ordre 2, par exemple x″ (t) + 2x′ (t) + x(t) = 1 +
t, x(0) = 1, x′ (0) = 2 :
0.2
dsolve(diff(x(t),t,2)+2*diff(x(t),t)
+x(t)-(1+t),x(t), ics={x(0):1,
0.0 diff(x(t),t).subs(t,0):2})
1 2 3 4 5 6 7 8 9
x(t) = t + (3t + 2)e−t − 1
La somme de la série (Σun ) s’obtient avec
summation() :
summation(u(n),(n,1,oo)) 19 Vecteurs et matrices
2
π
Un vecteur est défini comme une matrice à une dimen-
6 sion de la façon suivante :
X = Matrix([1,2,3]); X
17 Equations récurrentes linéaires 1
2
Pour résoudre l’équation un+1 = 3un + n , on com- ⎡2⎤
mence par définir la suite inconnue (un ) avec l’instruc- ⎣3⎦
tion Function() (attention au « F » majuscule) :
On peut alors calculer sa norme avec norm() :
u = Function('u')
[Link]()
On peut alors résoudre l’équation en l’écrivant sous la
forme un+1 − 3un − n2 (sans le = 0) : √14
rsolve(u(n+1)-3*u(n)-n**2, u(n)) Le calculs basiques se font assez naturellement :
2
n n 1 Y = Matrix([4,5,6]); X + 2*Y
3n C0 − − −
2 2 2 9
Des conditions initiales
peuvent être indi- ⎡12⎤
quées, en utilisant la
syntaxe particulière ⎣15⎦
{u(n0):u0,u(n1):u1,…} :
Le produit scalaire utilise une syntaxe un peu particu-
rsolve(u(n+1)-3*u(n)-n**2, u(n), lière :
{u(0):2})
[Link](Y)
5 ⋅ 3n n2 n 1
− − − 32
2 2 2 2
Plus généralement, une matrice se définit par :
18 Equations différentielles linéaires A = Matrix([[1,2,3],[4,5,6]]); A
La syntaxe pour résoudre une équation différentielle est 1 2 3
[ ]
très proche de celle des équations récurrentes. 4 5 6
Pour résoudre l’équation x′ (t) − tx(t) = t, on définit Les instructions eye(n) et diag() permettent d’obte-
d’abord la variable et la fonction inconnue : nir la matrice identité d’ordre n et une matrice diagonale
t = symbols('t') quelconque :
x = Function('x')
eye(3)
puis on utilise l’instruction dsolve() :
1 0 0
dsolve(diff(x(t),t)-t*x(t)-t,x(t))
⎡0 1 0⎤
t2
x(t) = C1 e 2 − 1 ⎣0 0 1⎦
6
diag(3,4,5) ou directement en utilisant l’instruction
3 0 0 diagonalize() :
⎡0 4 0⎤ P, D = [Link]()
⎣0 0 5⎦ D
La transposée est obtenue par tranpose() et le déter- 2 0 0
minant (d’une matrice carrée) par det() : ⎡0 3 0⎤
transpose(A) ⎣0 0 3⎦
1 4 P
⎡2 5⎤ −1 0 1
⎣3 6⎦ ⎡ 1 1 0⎤
⎣ 0 0 1⎦
B = Matrix([[2,5],[1,3]]);B
La puissance n-ième de la matrice s’obtient alors facile-
2 5
[ ] ment :
1 3
P*D**n*P**(-1)
det(B) 2n 0 −2n + 3n
1 ⎡−2 + 3 3n 2n − 3n ⎤
n n
⎣ 0 0 3n ⎦
Evidemment, le produit de deux matrices est obtenu avec
la multiplication standard * : Mais on aurait pu juste utiliser
B*A A**n
22 29 36
[ ]
13 17 21 20 Nombres complexes
De même, l’inverse d’une matrice (inversible !) est ob- Un nombre complexe est noté avec la notation naturelle
tenu en la mettant à la puissance −1 : (mais en utilisant un « I » majuscule).
B**(-1) z = 3+4*I; z
3 −5 3 + 4i
[ ]
−1 2 Les instructions re(), im(), abs(), arg() ren-
voient, respectivement, la partie réelle et imaginaire, le
Pour diagonaliser une matrice module et l’argument du nombre complexe.
A = Matrix([[2,0,1],[1,3,- abs(z)
1],[0,0,3]]);A 5
2 0 1 Le conjugué d’un nombre est obtenu avec l’instruction
⎡1 3 −1⎤ conjugate() :
⎣0 0 3 ⎦ conjugate(z)
on peut calculer le polynôme caractéristique : 3 − 4i
lamda = symbols('lamda', real=True)
[Link](lamda).as_expr() 21 Fonctions (encore)
λ3 − 8λ2 + 21λ − 18 Pour définir une fonction mathématique f, outre l’utilisa-
puis les valeurs propres : tion d’une fonction Python vue plus haut, il est possible
d’utiliser l’instruction lambdify() de la manière sui-
[Link]()
vante :
{3 ∶ 2, 2 ∶ 1} f = lambdify(x,1+x**3)
et les vecteurs propres :
f(2)
[Link]()
9
[(2,1,[Matrix([[-1],[1],[0]])]),
(3,2,[Matrix([[0],[1],[0]]), diff(f(x),x)
Matrix([[1],[0],[1]])])] 3x2
7
Cela est aussi possible avec les fonctions de plusieurs va- La fonction f admet donc un maximum local en x = −1
riables : et un minimum local en x = +1. Ces extremums sont
f = lambdify((x,y),1+x*y) f(-1)
2
f(3,4)
13 f(1)
Parfois, il est juste nécessaire de nommer une expression, −2
sans en faire véritablement une fonction : Alternative pour les conditions suffisantes :
f=3*x**2 candidats = solve(diff(f(x),x),x)
On utilise alors la propriété subs() pour donner des for x0 in candidats:
valeurs à la variable : print("f''(", x0, ") = ",
[Link](x,2) diff(f(x),x,2).subs(x,x0))
12 f''( -1 ) = -6
f''( 1 ) = 6
Cela, y compris pour des expressions de plusieurs va-
riables : Finalement, on peut faire une représentation graphique
pour illustrer le résultat :
f = x*y
plot(f(x),(x,-2.5,2.5))
[Link]([(x,3),(y,4)]) 8
12 6
4
22 Exemple : Optimisation 1 variable 2
0
f(x)
On souhaite optimiser la fonction f(x) = x3 − 3x sur ℝ. 2 1 0 1 2
On commence par définir la fonction : 2x
4
def f(x):
return x**3 - 3*x 6
Les limites en ±∞ indiquent que f n’admet pas d’extre- 8
mum global :
limit(f(x),x,+oo)
∞
limit(f(x),x,-oo)
−∞
En utilisant la dérivée de f :
diff(f(x),x)
3x2 − 2
on détermine les points critiques (candidats) :
solve(diff(f(x),x),x)
[−1, 1]
A l’aide de la dérivée seconde de f :
diff(f(x),x,2)
6x
on détermine la nature des points critiques :
diff(f(x),x,2).subs(x,-1)
−6
diff(f(x),x,2).subs(x,1)
6