Probleme mit odeint

mit matplotlib, NumPy, pandas, SciPy, SymPy und weiteren mathematischen Programmbibliotheken.
Antworten
Pumpkinpy
User
Beiträge: 10
Registriert: Freitag 26. Mai 2017, 11:45

Hallo alle zusammen.

Ich beschäschäftige mich gerade mit der numerischen Lösung von DGL. Für Python hab ich dazu die odeint-Funktion gefunden.
Ich bekomme meine DGL aus der vorangegangen Rechnung als Vektor übergeben. Sprich jede Zeile ist quasi eine der DGL, welche zum Teil gekoppelt sind.

Der Code zur berechnung sieht dann quasi so aus:

Code: Alles auswählen

#geht
import matplotlib.pyplot as plt
import scipy.integrate as sci
import numpy as np


def con(c, t):
        
    a=0
    cbdt = np.matrix([[-b_1 * c[1]], [b_1 * c[1] - b_2 * c[2]], [b_2 * c[2]]])

    return(a, cbdt[0], cbdt[1], cbdt[2])

b_1 = .5
b_2 = .1
c0 = [0, 1, 0, 0]


t = np.linspace(0 , 10, 1000)

sol = sci.odeint(con, c0 , t)


plt.plot(t, sol[:,1], 'r-' )
plt.plot(t, sol[:,2], 'g-' )
plt.plot(t, sol[:,3], 'b-')
plt.show()
Das Problem auf das ich gestoßen bin ist, dass sobald ich den a=0 Parameter entferne (natürlich auch die Ausgabe in der FKT und den Parameter beim Imput), den brauche ich eigenlich gar nicht, das Programm mit den Fehler "RuntimeError: The array return by func must be one-dimensional, but got ndim=3." wirft. Ich verstehe aber nicht wieso der Fehler mit der Nullfunktion nicht auftaucht.

Code: Alles auswählen

#geht nicht
import matplotlib.pyplot as plt
import scipy.integrate as sci
import numpy as np


def con(c, t):
        
    cbdt = np.matrix([[-b_1 * c[1]], [b_1 * c[1] - b_2 * c[2]], [b_2 * c[2]]])

    return(cbdt[0], cbdt[1], cbdt[2])

b_1 = .5
b_2 = .1
c0 = [1, 0, 0]


t = np.linspace(0 , 10, 1000)

sol = sci.odeint(con, c0 , t)


plt.plot(t, sol[:,1], 'r-' )
plt.plot(t, sol[:,2], 'g-' )
plt.plot(t, sol[:,3], 'b-')
plt.show()
Wenn ich die Gleichungen einzeln aufschtreibe kommt das gleiche Schaubild raus wie wenn ich es als Vektor aufschreibe, daher Verstehe ich das Problem nicht das mit Python hier zurückgibt. Da das ganze aber einen Fehler produziert den ich durch den "Pfusch a=0" unterdrücke, würde ich gerne wissen woran es liegt, nicht dass mir da später noch Fehler auftauchen die darauf zurückzuführen sind.
Wenn mir jemand das Phänomen erklären könnte und vll auch eine korrekte Lösung dafür anbietet wäre ich sehr dankbar.

Mit freundlichen Grüßen
Pumpkinpy
Zuletzt geändert von Anonymous am Samstag 27. Mai 2017, 11:30, insgesamt 1-mal geändert.
Grund: Quelltext in Python-Codebox-Tags gesetzt.
Benutzeravatar
BigZ
User
Beiträge: 30
Registriert: Mittwoch 1. Juli 2015, 21:18
Wohnort: Hamburg
Kontaktdaten:

Das Problem wird ersichtlich wenn du dir den Return Value von con(c0, t) anschaust.
con(c0, t) gibt ein Tupel zurück.
Wenn du die 0 weg lässt hast du im Tupel 3 matrizen (oder wie auch immer die Mehrzahl von Matrix ist).
Das ergibt dann so was wie:
(
( [[-0.5]] ),
( [[ 0.5]] ),
( [[ 0.]] )
)
was faktisch gesehen ein drei dimensionales Array ist, odeint verlangt aber offenbar ein eindimensionales array.
Warum odeint darin ein eindimensionales Array mit der 0 vorweg erkennt, ist mir auch noch ein bisschen unklar,
aber du könntest eventuell mit der np.reshape aus dem drei dimensionalem Array ein ein dimensionales Array machen ohne dabei die Daten zu korrumpieren.

edit sagt:
np.matrix hat auch die funktion flatten(), aber dabei kommt dann ein zwei dimensionales Array raus, aus dem man noch ein eindimensionales basteln müsste, für odeint zum weiterverarbeiten. :)

edit² sagt:
cbdt.flat gibt ein ein dimensionales Array zurück, aber offenbar fehlt dann immer noch ein weiteres Element im Array :/
"Ist noch λ?"
"Ja, aber das ϕ ist noch ϱ"
Pumpkinpy
User
Beiträge: 10
Registriert: Freitag 26. Mai 2017, 11:45

Ok habe mal nochmal etwas rumprobiert. Da odeint ja anscheinend keine Vektren oder Matrix verarbeiten kann (zumindest als erster Übergabewert :shock: ) habe ich mittels Transponierung und umschreiben in ein Tupel das ganze Umgangen. Wenn jemand eine schönere Lösung hat bin ich über weitere Vorschläge glücklich. Meine Lösung sieht zur Zeit so aus:

Code: Alles auswählen

import matplotlib.pyplot as plt
import scipy.integrate as sci
import numpy as np


def con(c, t):
        
    adt = np.matrix([[-b_1 * c[0]], [b_1 * c[0] - b_2 * c[1]], [b_2 * c[1]]])

    bdt = np.array(adt.T)

    cdt = tuple(bdt[0])

    return(cdt)

b_1 = 1
b_2 = 1
c0 = [1, 0, 0]



t = np.linspace(0 , 10, 1000)

sol = sci.odeint(con, c0 , t)

plt.plot(t, sol[:,0], 'r-' )
plt.plot(t, sol[:,1], 'g-' )
plt.plot(t, sol[:,2], 'b-')
plt.show()
Danke euch und schönes Wochenende

Pumpkinpy
Zuletzt geändert von Anonymous am Samstag 27. Mai 2017, 11:44, insgesamt 1-mal geändert.
Grund: Quelltext in Python-Codebox-Tags gesetzt.
Antworten