IdentifiantMot de passe
Loading...
Mot de passe oublié ?Je m'inscris ! (gratuit)
Navigation

Inscrivez-vous gratuitement
pour pouvoir participer, suivre les réponses en temps réel, voter pour les messages, poser vos propres questions et recevoir la newsletter

Calcul scientifique Python Discussion :

Transformée inverse Laplace à pôles complexes


Sujet :

Calcul scientifique Python

  1. #1
    Futur Membre du Club
    Homme Profil pro
    Étudiant
    Inscrit en
    Novembre 2019
    Messages
    2
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Âge : 26
    Localisation : France, Haute Savoie (Rhône Alpes)

    Informations professionnelles :
    Activité : Étudiant

    Informations forums :
    Inscription : Novembre 2019
    Messages : 2
    Points : 5
    Points
    5
    Par défaut Transformée inverse Laplace à pôles complexes
    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 : Sélectionner tout - Visualiser dans une fenêtre à part
    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 : Sélectionner tout - Visualiser dans une fenêtre à part
    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 : Sélectionner tout - Visualiser dans une fenêtre à part
    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
    Images attachées Images attachées  

  2. #2
    Membre habitué
    Homme Profil pro
    Vagabong étudiant en annalyse du signal.
    Inscrit en
    Avril 2019
    Messages
    123
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Âge : 25
    Localisation : France, Isère (Rhône Alpes)

    Informations professionnelles :
    Activité : Vagabong étudiant en annalyse du signal.
    Secteur : High Tech - Multimédia et Internet

    Informations forums :
    Inscription : Avril 2019
    Messages : 123
    Points : 162
    Points
    162
    Par défaut 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é.

  3. #3
    Futur Membre du Club
    Homme Profil pro
    Étudiant
    Inscrit en
    Novembre 2019
    Messages
    2
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Âge : 26
    Localisation : France, Haute Savoie (Rhône Alpes)

    Informations professionnelles :
    Activité : Étudiant

    Informations forums :
    Inscription : Novembre 2019
    Messages : 2
    Points : 5
    Points
    5
    Par défaut 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.

+ Répondre à la discussion
Cette discussion est résolue.

Discussions similaires

  1. [DirectX9/HLSL] Transformation inverse.
    Par F-fisher dans le forum DirectX
    Réponses: 5
    Dernier message: 10/01/2010, 21h01
  2. transformée de laplace
    Par fabrice_francois dans le forum MATLAB
    Réponses: 5
    Dernier message: 07/08/2009, 16h08
  3. [Débutant] transformée de laplace inverse
    Par slyphe dans le forum Signal
    Réponses: 2
    Dernier message: 06/04/2009, 14h08
  4. [Débutant] Transformation de Laplace discrète sur un vecteur de type double ?
    Par Grandblabla26 dans le forum MATLAB
    Réponses: 3
    Dernier message: 21/12/2008, 20h02
  5. transformation inverse ondelette
    Par harafado dans le forum MATLAB
    Réponses: 3
    Dernier message: 03/12/2008, 17h34

Partager

Partager
  • Envoyer la discussion sur Viadeo
  • Envoyer la discussion sur Twitter
  • Envoyer la discussion sur Google
  • Envoyer la discussion sur Facebook
  • Envoyer la discussion sur Digg
  • Envoyer la discussion sur Delicious
  • Envoyer la discussion sur MySpace
  • Envoyer la discussion sur Yahoo