如何从dy_dt(y,x,t)方程绘制y,x图?

问题描述 投票:0回答:1

我从dy_dt(y,t,x)方程中绘制了x线变化的y,t图。像这样。

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
import math

def dy_dt(y,t,z):
    c_1 = 5.0 / (1.38 * 1223.0 * pow(10.0, 28.0)*pow(z,3.0))
    c_2 = 6*pow(10.0, 5.0)
    c_3 = 9.5*pow(10,19)
    c_4 = 9.5*pow(10,10)
    c_5 = pow(600000,2.7)
    return c_1 * ( y/math.sqrt(1.0+pow(y, 2.0)) ) * ( ((1.0-pow(y, 3.0))/(z* pow(y,(1.0/3.0)))) - (y * c_2))-(((3*y*(c_3)*np.exp(-36*1570/1223)))/((2*pow(c_4,2.7)))*(c_5))

t = np.linspace(0, pow(10.0, 12.0), pow(10.0, 6.0))
y_0 = 1.0
z = 0.00001
y1 = odeint(dy_dt, y_0, t, args=(z,))
z = 0.000015
y2 = odeint(dy_dt, y_0, t, args=(z,))
y1 = np.array(y1).flatten()
y2 = np.array(y2).flatten()

plt.rcParams.update({'font.size': 10})
plt.ylim(0,1.0)
plt.xlabel("x")
plt.ylabel("y")
plt.plot(t, y1, 'r-')
plt.plot(t, y2, 'b-')
plt.grid(True)
plt.show()

但是,我不确切地知道如何从等式中绘制出t变化的y,x图。我就是这样,但是有ZeroDivisionError。

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
import math

def dy_dt(y,x,t):
    c_1 = 5.0 / (1.38 * 1223.0 * pow(10.0, 28.0)*pow(x,3.0))
    c_2 = 600000
    c_3 = 9.5*pow(10,19)
    c_4 = 9.5*pow(10,10)
    c_5 = pow(600000,2.7)
    return c_1 * (y/math.sqrt(1.0+pow(y, 2.0)) ) * ( ((1.0-pow(y, 3.0))/(x* pow(y,(1.0/3.0)))) - (y * c_2))-(((3*y*(c_3)*np.exp(-36*1570/1223)))/((2*pow(c_4,2.7)))*(c_5))

xs = np.linspace(0, 0.0006, 100)
y_0 = 1.0
t = pow(10, 11)
y1 = odeint(dy_dt, y_0, xs, args=(t,))
t = pow(10, 12)
y2 = odeint(dy_dt, y_0, xs, args=(t,))
y1 = np.array(y1).flatten()
y2 = np.array(y2).flatten()



plt.rcParams.update({'font.size': 10})
plt.ylim(0,1.0)
plt.xlabel("x")
plt.ylabel("y")
plt.plot(xs, y1, 'r-')
plt.plot(xs, y2, 'b-')
plt.grid(True)
plt.show()

在前一张图中和后一张图中,初始条件均为y(0)= 1。(t = 0-> y = 1,x = 0-> y = 1)

python plot odeint
1个回答
0
投票

您正在c1中的第二个函数和返回行中除以x。当您的xs从0到0.0006时,这将导致零除错误。

要么从linspace中的一个很小的数字开始,要么检查x == 0并引入避免零除的特殊情况。

© www.soinside.com 2019 - 2024. All rights reserved.