Parametrisches Diagramm der Lösung von 2x2 diff.System in Python, Mathematica
Frage
Ich habe eine Lösung für das folgende Gleichungssystem implementiert
dy/dt = -t*y(t) - x(t)
dx/dt = 2*x(t) - y(t)^3
y(0) = x(0) = 1.
0 <= t <= 20
zunächst in Mathematica und anschließend in Python.
Mein Code in Mathematica:
s = NDSolve[
{x'[t] == -t*y[t] - x[t], y'[t] == 2 x[t] - y[t]^3, x[0] == y[0] == 1},
{x, y}, {t, 20}]
ParametricPlot[Evaluate[{x[t], y[t]} /. s], {t, 0, 20}]
Daraus ergebe ich folgende Handlung: Handlung1 (Wenn die Meldung 403 Forbidden angezeigt wird, drücken Sie bitte die Eingabetaste im URL-Feld.)
Später habe ich dasselbe in Python codiert:
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
g = lambda t: t
def f(z,t):
xi = z[0]
yi = z[1]
gi = z[2]
f1 = -gi*yi-xi
f2 = 2*xi-yi**3
return [f1,f2]
# Initial Conditions
x0 = 1.
y0 = 1.
g0 = g(0)
z0 = [x0,y0,g0]
t= np.linspace(0,20.,1000)
# Solve the ODEs
soln = odeint(f,z0,t)
x = soln[:,0]
y = soln[:,1]
plt.plot(x,y)
plt.show()
Und das ist die Handlung, die ich bekomme:Handlung2 (Wenn die Meldung 403 Forbidden angezeigt wird, drücken Sie bitte die Eingabetaste im URL-Feld.)
Wenn man die Mathematica-Lösung noch einmal in einem kleineren Feld aufträgt:
ParametricPlot[Evaluate[{x[t], y[t]} /. s], {t, 0, 6}]
er wird ein ähnliches Ergebnis wie die Python-Lösung erhalten.Nur die Achse wird verlegt.
Warum gibt es so große Unterschiede in den Handlungssträngen?Was mache ich falsch?
Ich vermute, dass meine Python-Implementierung des Modells falsch ist, insbesondere dort, wo f1 berechnet wird.Oder vielleicht ist die Funktion plot() überhaupt nicht zum Zeichnen parametrischer Gleichungen geeignet, wie in diesem Fall.
Danke.
PS:Es tut mir leid, dass ich Ihnen das Leben schwer gemacht habe, indem ich die Bilder nicht in den Text eingefügt habe.Ich habe noch nicht genug Ruf.
Lösung
Du verwendest t
als Ihr dritter Parameter im Eingabevektor, nicht als separater Parameter.Der t
In f(z,t)
wird nie verwendet;Stattdessen verwenden Sie z[2]
, was nicht dem Bereich von entspricht t
wie Sie es zuvor definiert haben (t=np.linspace(0,20.,1000)
).Der lambda
Funktion für g
hilft hier nicht weiter:Sie verwenden es nur einmal, um eine festzulegen t0
, aber nie danach.
Vereinfachen Sie Ihren Code und entfernen Sie diesen dritten Parameter aus Ihrem Eingabevektor (sowie aus der Lambda-Funktion).Zum Beispiel:
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
def f(z,t):
xi = z[0]
yi = z[1]
f1 = -t*yi-xi
f2 = 2*xi-yi**3
return [f1,f2]
# Initial Conditions
x0 = 1.
y0 = 1.
#t= np.linspace(0,20.,1000)
t = np.linspace(0, 10., 100)
# Solve the ODEs
soln = odeint(f,[x0,y0],t)
x = soln[:,0]
y = soln[:,1]
ax = plt.axes()
#plt.plot(x,y)
plt.plot(t,x)
# Put those axes at their 0 value position
ax.spines['left'].set_position('zero')
ax.spines['bottom'].set_position('zero')
ax.spines['right'].set_color('none')
ax.spines['top'].set_color('none')
ax.xaxis.set_ticks_position('bottom')
ax.yaxis.set_ticks_position('left')
#plt.axis([-0.085, 0.085, -0.05, 0.07])
plt.show()
Ich habe die eigentliche Handlung, die Sie wollen, auskommentiert und schreibe stattdessen eine Handlung x
gegen t
, was Sie in den Kommentaren haben, da ich der Meinung bin, dass es jetzt einfacher ist zu erkennen, dass die Dinge richtig sind.Die Zahl, die ich bekomme: