Pregunta

He implementado una solución para el siguiente sistema de ecuaciones

dy/dt = -t*y(t) - x(t)
dx/dt = 2*x(t) - y(t)^3

y(0) = x(0) = 1.
0 <= t <= 20

en primer lugar, en Mathematica y después en Python.

Mi código de 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}]

Desde que tengo la siguiente trama: Plot1 (si te da un mensaje 403 Prohibido por favor presione entrar dentro del campo de la url)

Más tarde me codificados de la misma en python:

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()

Y esta es la trama que obtengo:Plot2 (si te da un mensaje 403 Prohibido por favor presione entrar dentro del campo de la url)

Si uno traza de nuevo la solución de Mathematica en un pequeño campo:

ParametricPlot[Evaluate[{x[t], y[t]} /. s], {t, 0, 6}]

se obtendrá un resultado similar a la solución python.Sólo el eje " va a estar fuera de lugar.

¿Por qué hay una diferencia tan grande en las parcelas?¿Qué estoy haciendo mal?

Sospecho que mi python implementación del modelo es incorrecto, especialmente donde f1 es calculado.O tal vez la trama() la función no está al alcance de la mano para trazar ecuaciones paramétricas como en este caso.

Gracias.

ps:lo siento por hacer la vida difícil por no bofetadas de las imágenes dentro del texto;No tengo la reputación suficiente todavía.

¿Fue útil?

Solución

Estás usando t como su tercer parámetro en el vector de entrada, no como un parámetro independiente.El t en f(z,t) nunca se utiliza;en su lugar, utilice z[2], que no es igual el rango de t como definir antes (t=np.linspace(0,20.,1000)).El lambda función para g no ayuda aquí:usted lo utiliza solamente una vez para un t0, pero nunca después.

Simplificar el código, y quitar ese tercer parámetro a partir de su vector de entrada (así como la función lambda).Por ejemplo:

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()

Me han comentado que la trama real que usted desea, y en su lugar estoy tramando x frente a t, lo que tienen en los comentarios, ya que siento que hace que sea más fácil ver las cosas es correcta ahora.La figura de recibir:

enter image description here

Licenciado bajo: CC-BY-SA con atribución
No afiliado a StackOverflow
scroll top