Transformée inverse Laplace à pôles complexes

Version imprimable

1 pièce(s) jointe(s)

Bonjour,
Je travail sur de l'analyse système avec des équation différentielles (le bonheur !). Pour simplifier et automatiser toutes les calculs j'utilise Jupyter avec numpy et sympy.
J'ai une équation dans l'espace laplacien (variable s) : c**2/((c+m)(c*s**2+b*s−k)) − 1/(s**2*(c+m)) et impossible d'appliquer la transformée de Laplace !
Ce que je ne comprend pas c'est que quand j'utilise

Code:

sp.inverse_laplace_transform(-1/(s**2*(c+m)), s, t)

ça donne bien -t/(c+m)

Je pense que c'est un problème en rapport avec les du polynôme. Pourtant le résultat dans les tables de transformées de Laplace est connu. Je pourrais le faire à la main mais par curiosité j'aimerai savoir comment faire sur python.

Merci d'avance de m'aider, même si c'est peut être un peu pointu en math.

L'intégrale du code :

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
 
%matplotlib nbagg
import numpy as np
import sympy as sp
import matplotlib.pyplot as plt
sp.init_printing() #Affichage des outputs en LaTex
 
#Définition des variables avec sympy
M_m, M_c, g, b, k, s = sp.symbols('M_m M_c g b k s')
t = sp.Symbol('t', positive=True)
 
A = sp.Matrix([[0,0,1,0],
              [0,0,0,1],
              [0, k/(M_m+M_c), 0, -b/(M_m+M_c)],
              [0, k/M_c, 0, -b/M_c]])
B = sp.Matrix([[0,0],
              [0,0],
              [-1, 1/(M_m+M_c)],
              [-1, 0]])
phi=(s*sp.eye(4)-A)**-1
phi.simplify()
 
G = np.array(B)
for i in range(0, 4):
    for j in range(0, 2):
        G[i][j]=sp.apart((sp.eye(4).row(i)*phi*B)[j],s)
 
G_mat = sp.Matrix([[G[0][0], G[0][1]],
                   [G[1][0], G[1][1]],
                   [G[2][0], G[2][1]],
                   [G[3][0], G[3][1]]])
 
sp.inverse_laplace_transform(G[0][0], s, t)

Le code erreur :

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-22-c5aa2f068acb> in <module>()
----> 1 sp.inverse_laplace_transform(G[0][0],s,t)
 
~\Anaconda3\lib\site-packages\sympy\integrals\transforms.py in inverse_laplace_transform(F, s, t, plane, **hints)
   1269     if isinstance(F, MatrixBase) and hasattr(F, 'applyfunc'):
   1270         return F.applyfunc(lambda Fij: inverse_laplace_transform(Fij, s, t, plane, **hints))
-> 1271     return InverseLaplaceTransform(F, s, t, plane).doit(**hints)
   1272 
   1273 
 
~\Anaconda3\lib\site-packages\sympy\integrals\transforms.py in doit(self, **hints)
    115             try:
    116                 return self._compute_transform(self.function,
--> 117                     self.function_variable, self.transform_variable, **hints)
    118             except IntegralTransformError:
    119                 pass
 
~\Anaconda3\lib\site-packages\sympy\integrals\transforms.py in _compute_transform(self, F, s, t, **hints)
   1221 
   1222     def _compute_transform(self, F, s, t, **hints):
-> 1223         return _inverse_laplace_transform(F, s, t, self.fundamental_plane, **hints)
   1224 
   1225     def _as_integral(self, F, s, t):
 
~\Anaconda3\lib\site-packages\sympy\integrals\transforms.py in wrapper(*args, **kwargs)
    193         def wrapper(*args, **kwargs):
    194             noconds = kwargs.pop('noconds', default)
--> 195             res = func(*args, **kwargs)
    196             if noconds:
    197                 return res[0]
 
~\Anaconda3\lib\site-packages\sympy\integrals\transforms.py in _inverse_laplace_transform(F, s, t_, plane, simplify)
   1148     try:
   1149         f, cond = inverse_mellin_transform(F, s, exp(-t), (None, oo),
-> 1150                                            needeval=True, noconds=False)
   1151     except IntegralTransformError:
   1152         f = None
 
~\Anaconda3\lib\site-packages\sympy\integrals\transforms.py in inverse_mellin_transform(F, s, x, strip, **hints)
    861     hankel_transform, inverse_hankel_transform
    862     """
--> 863     return InverseMellinTransform(F, s, x, strip[0], strip[1]).doit(**hints)
    864 
    865 
 
~\Anaconda3\lib\site-packages\sympy\integrals\transforms.py in doit(self, **hints)
    140             try:
    141                 extra = self._collapse_extra(extra)
--> 142                 return tuple([res]) + tuple(extra)
    143             except IntegralTransformError:
    144                 pass
 
TypeError: 'NoneType' object is not iterable

25/12/2019, 13h45
robinechuca

réponsse qui n'aide pas

Bonjours,

Votre code tourne toujours sur mon ordinateur depuis maintenant plus d'une heure, il ne renvoi rien, donc pas d'erreur.
J'utilise la version '1.4' de sympy. Peut-être avez-vous une version d'avant?

Avez vous vraiment besoin d'une transformation inverse 'formelle'? Car si ce n'est pas le cas, vous pouvez repasser
par la définition avec l'intégrale affin de "prémâcher le travail de sympy". Bon ça résout pas tout car il faut encore réussir à trouver l’abscisse de sommabilité.

Si vous définissez M_m, M_c, g, b, k comme des réels positifs, votre code renvoi:
(M_c**2*exp(t*(b - sqrt(4*M_c*k + b**2))/(2*M_c)) - M_c**2*exp(t*(b + sqrt(4*M_c*k + b**2))/(2*M_c)) - M_m*t*sqrt(4*M_c*k + b**2)*exp(b*t/M_c))*exp(-b*t/M_c)/((M_c + M_m)*sqrt(4*M_c*k + b**2))

Mais je n'ai pas compris à quoi correspondent ces variables (M_m, M_c, g, b, k) donc je ne sait pas si cette hypothèse est censé.
29/12/2019, 21h17
flbos74

Remerciement et Solution peu élégantes

Merci beaucoup d'avoir essayé aussi et d'obtenir le même résultat que moi (infinite loop). Désormais je suis sur que le problème vient de la manière dont fonctionne la fonction de l'inverse de la transformée de Laplace dans Sympy (par identification car cela ne marche que pour des fractions simples).

Pour ceux qui aurait le même problème, j'ai réussi à contourner le problème en utilisant
- Transformation en éléments simples avec sp.apart(exp, s)
- Conversion en chaine de caractères avec str(exp)
- Découpage en plusieurs chaines de caractères en ciblant les " + " et les " - "
- Stockages des chaines dans un tableau
- Conversion en expression avec eval(exp)
- Application de la transformée de Laplace inverse sur chaque expression "simple" du tableau
- Somme de toutes les transformées inverses

C'est vraiment la galère mais efficace. Peut-être qu'une mise à jour de Sympy corrigera ce problème un jour même si je me doute que mathématiquement, cela doit être une horreur de définir une fonction universelle pour la transformée de Laplace.