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 :

scipy et resolution F(x)=0


Sujet :

Calcul scientifique Python

  1. #1
    Futur Membre du Club
    Profil pro
    Inscrit en
    Septembre 2007
    Messages
    10
    Détails du profil
    Informations personnelles :
    Localisation : France

    Informations forums :
    Inscription : Septembre 2007
    Messages : 10
    Points : 5
    Points
    5
    Par défaut scipy et resolution F(x)=0
    Bonjour
    j'ai réalisé un programme en python qui résout 2 fois un problème du type F(x)=0. j'ai utilisé plusieurs méthode (brentq, ec) de root_scalar de scipy.optimize. Je vois que la méthode ne trouve pas toujours la solution car les fonctions ne sont pas toujours monotones croissante. Je dois donc souvent ajuster l'intervalle d'étude pour que le système converge. J'ai testé de résoudre le même problème avec le solveur excel avec des contraintes et là aucun problème de convergence. Je ne pense pas que le solder d'excel soi plus puissant que python. Sauriez-vous quelle est la meilleure bibliothèque de python résoudre de type de problèmes. Je précise que mes connaissance sont très faibles en python. je suis prêt à vous montrer des portions du code mais je ne sais pas comment on peut les insérer dans le message.
    Merci pour votre aide.

  2. #2
    Expert éminent
    Avatar de fred1599
    Homme Profil pro
    Lead Dev Python
    Inscrit en
    Juillet 2006
    Messages
    3 978
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : France, Meurthe et Moselle (Lorraine)

    Informations professionnelles :
    Activité : Lead Dev Python
    Secteur : Arts - Culture

    Informations forums :
    Inscription : Juillet 2006
    Messages : 3 978
    Points : 7 407
    Points
    7 407
    Par défaut
    Bonjour,

    j'ai réalisé un programme en python qui résout 2 fois un problème du type F(x)=0
    Hein ? Pourquoi 2 fois ?

    j'ai utilisé plusieurs méthode (brentq, ec) de root_scalar de scipy.optimize.
    Dommage que vous n'ayez pas donné toutes les méthodes testées

    Je ne pense pas que le solder d'excel soi plus puissant que python.
    Vous visez juste

    Sauriez-vous quelle est la meilleure bibliothèque de python résoudre de type de problèmes.
    scipy semble très bien, avez vous testé avec fsolve ?

    je suis prêt à vous montrer des portions du code mais je ne sais pas comment on peut les insérer dans le message
    En regardant les règles du forum sur comment présenter un code. Je vous aide, utilisez le bouton #
    Celui qui trouve sans chercher est celui qui a longtemps cherché sans trouver.(Bachelard)
    La connaissance s'acquiert par l'expérience, tout le reste n'est que de l'information.(Einstein)

  3. #3
    Expert éminent sénior
    Homme Profil pro
    Architecte technique retraité
    Inscrit en
    Juin 2008
    Messages
    21 548
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : France, Manche (Basse Normandie)

    Informations professionnelles :
    Activité : Architecte technique retraité
    Secteur : Industrie

    Informations forums :
    Inscription : Juin 2008
    Messages : 21 548
    Points : 37 198
    Points
    37 198
    Par défaut
    Salut,

    Citation Envoyé par fulgar Voir le message
    Je dois donc souvent ajuster l'intervalle d'étude pour que le système converge.
    Chaque méthode de calcul numérique aura avantages et inconvénients.
    A vous de (sa)voir laquelle sera la plus pertinente dans chaque cas de figure et quelles données fournir pour que ça "converge".
    Par exemple avec:
    • dichotomie: ce sera un intervalle [a, b] tel que f(a)*f(b) < 0... et ça trouvera une des racines de f dans cet intervalle.
    • newton: une valeur proche d'une racine (à calculer)
    • brent: des intervalles contenant une seule racine (sinon ca peut ne pas converger).


    Savoir quel outil de calcul numérique utiliser (et comment) n'a rien à voir avec la programmation Python.
    (Quel algorithme est utilisé par le solveur Excel peut être une question intéressante mais à poser dans un forum Excel...)

    Citation Envoyé par fulgar Voir le message
    Je précise que mes connaissance sont très faibles en python.
    Ce qu'il faudrait booster, dans le cas présent, c'est côté connaissances en calcul numérique.

    - W
    Architectures post-modernes.
    Python sur DVP c'est aussi des FAQs, des cours et tutoriels

  4. #4
    Futur Membre du Club
    Profil pro
    Inscrit en
    Septembre 2007
    Messages
    10
    Détails du profil
    Informations personnelles :
    Localisation : France

    Informations forums :
    Inscription : Septembre 2007
    Messages : 10
    Points : 5
    Points
    5
    Par défaut
    Bonjour,
    merci pour vos retour. pour répondre à vos questionnement, on programme est structure en deux boucle imbriquées l'une dans l'autre la première itère sur la recherche d'un débit d'air QA qui est utilisée dans la deuxième qui itère pour la recherche d'une température TE2. La recherche de TE2 qui implique d'avoir une valeur de QA se fait en cherchant la racine de la fonction DeltaMerkel :
    Code : Sélectionner tout - Visualiser dans une fenêtre à part
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    def DeltaMerkel(teau_froide):
        global DonE,Dontemp
        TE1 = DonE.dte + teau_froide # Calcul de TE1 en fonction de TE2
        Dontemp.QE=calculer_QE(DonE,TE1)
        Dontemp.RapQ=Dontemp.QE/Dontemp.QA
        Merkelconstruteur = Me_etoile(DonR)
        Merkel_calcul,_= Merkel(teau_froide, TE1)
        delta = Merkel_calcul-Merkelconstruteur
        Dontemp.Me= Merkelconstruteur
        return delta
    pour la recherche de cette racine j'utilise cette fonction :
    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
    def Calcul_TE2():
        global Dontemp, air_ent, DonE
        # Calculer la valeur de th
        th, _, _, _, _, _, _ = th_de_t_et_P_et_r(air_ent.tair, DonE.P, air_ent.r1)
        # Définir les bornes de l'intervalle
        intervalle_temperature_minimum = 19 #th[0] + 9
        intervalle_temperature_maximum = 30#th[0] + 20 + DonE.dte
        print(f"min={intervalle_temperature_minimum:.2f}, max={intervalle_temperature_maximum:.2f}")
        try:
            # Utiliser root_scalar avec la méthode spécifiée
            # techniques possibles :bisect,brentq,brenth,ridder,toms748 
            sol = root_scalar(DeltaMerkel, method='brentq', bracket=[intervalle_temperature_minimum, intervalle_temperature_maximum], xtol=1e-4)
            if sol.converged:
                Dontemp.TE2 = sol.root
                #print(f"Racine trouvée: {sol.root}")
            else:
                print(f"La méthode {method} n'a pas convergé.")
        except ValueError as e:
            print(f"Erreur lors de la recherche de la racine: {e}")
    le débit QA est utilisée dans la variable Dontemp.RapQ
    ensuit pour une valeur de TE2 solution de DeltaMerkel=0, j'itère ensuite sur le débit QA via ces deux fonctions selon la même philosophie. Je cherche le débit QA pour lequel Delta_ro2 s'annule :
    Code : Sélectionner tout - Visualiser dans une fenêtre à part
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    def Delta_ro2(Q_air):
        global Dontemp,air_ent, DonR
        Dontemp.QA=Q_air
        Calcul_vit_debitante()
        Calcul_TE2()
        print(f"La température pour laquelle DeltaMerkel(teau_froide) s'annule est {Dontemp.TE2:.4f} °C")
        calcul_air_chaud()
        Dontemp.CF=calcul_Cf()
        ro2_constructeur = air_ent.ro1-air_ent.ro1*((Dontemp.CF*(Dontemp.VD)**2)/(2*DonR.g*DonR.H))
        print((f"ro2 constructeur = {ro2_constructeur:.4f}"))
        deltaro2=ro2_constructeur-Dontemp.ro2
        return deltaro2
    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
    def Calcul_QA():
        global Dontemp
     
        # Définir les bornes de l'intervalle pour Q_air
        Q_air_min = 24000 # Valeur minimum arbitraire, mais positive
        Q_air_max = 30000   # Valeur maximum estimée pour QA, à ajuster si nécessaire
     
        try:
            # Utiliser root_scalar pour rechercher la valeur de Q_air pour laquelle Delta_ro2(Q_air) = 0
            solution = root_scalar(Delta_ro2, method='brentq', bracket=[Q_air_min, Q_air_max], xtol=1e-4)
     
            if solution.converged:
                Dontemp.QA = solution.root
                print(f"Solution trouvée : Dontemp.QA = {Dontemp.QA:.4f}")
            else:
                print("La méthode brentq n'a pas convergé.")
     
        except ValueError as e:
            print(f"Erreur lors de la recherche de la racine pour Q_air : {e}")
    je vois par exemple que lorsque je modifie Qair_min=2400 convergence OK avec les bons résultats : par contre en mettant Qair_min=2400 j'ai le bon résultat avec ce message
    Y:\FCT\04_PYTHON\Calcul aéro_1.py:126: IntegrationWarning: The integral is probably divergent, or slowly convergent.
    result, error = integrate.quad(integrandMe, TE2, TE1,epsabs=precision)

    je vois bien que ce sont mes bornes de la recherche de fonction qui pose problème notamment pour le calcul du zéro de fonction. peut-être que le produit des bornes de l'intervalle lors de la recherche du zero ne sont pas toujours négatifs. Je peux peut-être essayé une méthode de balayage en amont de la recherche du zero de fonction pour trouver 2 bornes dont le produit serait négatif. j'ai aussi essayé bisect et brenth mais j'ai le même problème. J'ai pense à utiliser Newton-Raphson mais cela me semblait plus complexe car je devais écrire une nouvelle fonction décrivant la dérivée.
    Quelqu'un aurait-il une idée pour rendre ce code moins sensible au borne d"initialisation de QA et TE2. Merci pour votre aide.

  5. #5
    Expert éminent sénior
    Homme Profil pro
    Architecte technique retraité
    Inscrit en
    Juin 2008
    Messages
    21 548
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : France, Manche (Basse Normandie)

    Informations professionnelles :
    Activité : Architecte technique retraité
    Secteur : Industrie

    Informations forums :
    Inscription : Juin 2008
    Messages : 21 548
    Points : 37 198
    Points
    37 198
    Par défaut
    Citation Envoyé par fulgar Voir le message
    Quelqu'un aurait-il une idée pour rendre ce code moins sensible au borne d"initialisation de QA et TE2.
    Il faut étudier la fonction pour savoir quelle sera la moins mauvaise méthode et les données initiales à fournir. C'est un problème de calcul numérique à résoudre en demandant de l'aide du côté des forums algorithmique & C°. En tout cas tant que vous restez scotché à vouloir que la solution soit côté code, vous n'avancerez pas.

    - W
    Architectures post-modernes.
    Python sur DVP c'est aussi des FAQs, des cours et tutoriels

  6. #6
    Futur Membre du Club
    Profil pro
    Inscrit en
    Septembre 2007
    Messages
    10
    Détails du profil
    Informations personnelles :
    Localisation : France

    Informations forums :
    Inscription : Septembre 2007
    Messages : 10
    Points : 5
    Points
    5
    Par défaut
    Citation Envoyé par wiztricks Voir le message
    Il faut étudier la fonction pour savoir quelle sera la moins mauvaise méthode et les données initiales à fournir. C'est un problème de calcul numérique à résoudre en demandant de l'aide du côté des forums algorithmique & C°. En tout cas tant que vous restez scotché à vouloir que la solution soit côté code, vous n'avancerez pas.

    - W
    Bonjour Wiztrick
    j'ai effectivement retraillé le code en mettant un test sur f(a)*f(b) pour calcule le zero de fonction comme ci-dessus :
    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
    def Calcul_TE2():
        global Dontemp, air_ent, DonE
        # Calculer la valeur de th
        th, _, _, _, _, _, _ = th_de_t_et_P_et_r(air_ent.tair, DonE.P, air_ent.r1)
        # Définir les bornes de l'intervalle
        intervalle_temperature_minimum = 18 #th[0] + 9
        intervalle_temperature_maximum = 30 #th[0] + 20 + DonE.dte
        print(f"min={intervalle_temperature_minimum:.2f}, max={intervalle_temperature_maximum:.2f}")
         # Balayage de l'intervalle pour vérifier  e signe de DeltaMerkel aux bornes
        delta_min = DeltaMerkel(intervalle_temperature_minimum)
        delta_max = DeltaMerkel(intervalle_temperature_maximum)
     
        # On vérifie que les signes sont opposés, sinon on réduit l'intervalle
        while delta_min * delta_max > 0 and intervalle_temperature_maximum - intervalle_temperature_minimum > 0.1:
            intervalle_temperature_minimum += 0.5
            intervalle_temperature_maximum -= 0.5
            delta_min = DeltaMerkel(intervalle_temperature_minimum)
            delta_max = DeltaMerkel(intervalle_temperature_maximum)
            print(f"Nouveau min={intervalle_temperature_minimum:.2f}, max={intervalle_temperature_maximum:.2f}")
     
        # Vérification finale avant d'utiliser root_scalar
        if delta_min * delta_max > 0:
            print("Impossible de trouver une racine avec cet intervalle, DeltaMerkel a le même signe aux bornes.")
            return None
     
        try:
            # Utiliser root_scalar avec la méthode spécifiée
            # techniques possibles :bisect,brentq,brenth,ridder,toms748 
            sol = root_scalar(DeltaMerkel, method='brentq', bracket=[intervalle_temperature_minimum, intervalle_temperature_maximum], xtol=1e-4)
            if sol.converged:
                Dontemp.TE2 = sol.root
                #print(f"Racine trouvée: {sol.root}")
            else:
                print(f"La méthode {method} n'a pas convergé.")
        except ValueError as e:
            print(f"Erreur lors de la recherche de la racine: {e}")
    mais j'ai toujours un problème de convergence de l'integrand avec ce message d'erreur même si le code me fournit malgré tout les bonnes valeurs (j'ai testé le même calcul sous excel")
    je vais suivre votre conseil et je vais essayer de poser mon problème chez les mathématiciens. merci pour vos conseils et votre réactivité dans les réponses.

Discussions similaires

  1. [SciPy] Résolution d'un système d'équas diff
    Par captainspaulding dans le forum Calcul scientifique
    Réponses: 5
    Dernier message: 20/07/2014, 17h33
  2. Erreur de code (résolution d'équadiff métode scipy)
    Par frozzen dans le forum Calcul scientifique
    Réponses: 1
    Dernier message: 01/06/2014, 16h23
  3. résolution d'equation f(x) = 0
    Par magicien dans le forum C
    Réponses: 8
    Dernier message: 06/05/2003, 17h06
  4. [Kylix] Resolution e Kylix
    Par ulisse dans le forum EDI
    Réponses: 1
    Dernier message: 02/03/2003, 16h57
  5. recuperer la résolution de l'écran
    Par florent dans le forum C++Builder
    Réponses: 11
    Dernier message: 07/06/2002, 16h01

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