如何使用Python内置函数odeint求解微分方程? [英] How to solve differential equation using Python builtin function odeint?

问题描述

我想在给定的初始条件下求解该微分方程:

I want to solve this differential equations with the given initial conditions:

(3x-1)y''-(3x+2)y'+(6x-8)y=0, y(0)=2, y'(0)=3

ans应该是

y=2*exp(2*x)-x*exp(-x)

the ans should be

y=2*exp(2*x)-x*exp(-x)

这是我的代码:

def g(y,x):
    y0 = y[0]
    y1 = y[1]
    y2 = (6*x-8)*y0/(3*x-1)+(3*x+2)*y1/(3*x-1)
    return [y1,y2]

init = [2.0, 3.0]
x=np.linspace(-2,2,100)
sol=spi.odeint(g,init,x)
plt.plot(x,sol[:,0])
plt.show()

但是我得到的答案却有所不同. 我做错了什么?

but what I get is different from the answer. what have I done wrong?

推荐答案

这里有几处错误.首先，你的方程式显然是

There are several things wrong here. Firstly, your equation is apparently

(3x-1)y''-(3x + 2)y'-(6x-8)y = 0; y(0)= 2，y'(0)= 3

(3x-1)y''-(3x+2)y'-(6x-8)y=0; y(0)=2, y'(0)=3

(请注意y中项的符号).对于这个方程，您的解析解和y2的定义是正确的.

(note the sign of the term in y). For this equation, your analytical solution and definition of y2 are correct.

其次，就像@Warren Weckesser所说的那样，您必须将2个参数作为y传递给g:y[0](y)，y[1](y')并返回它们的导数y'和y' '.

Secondly, as the @Warren Weckesser says, you must pass 2 parameters as y to g: y[0] (y), y[1] (y') and return their derivatives, y' and y''.

第三，给定x = 0的初始条件，但是要集成的x网格从-2开始.在odeint的文档中，此参数在呼叫签名说明中为t:

Thirdly, your initial conditions are given for x=0, but your x-grid to integrate on starts at -2. From the docs for odeint, this parameter, t in their call signature description:

odeint(func, y0, t, args=(),...):

t:数组 求解y的时间点序列.最初的 价值点应该是该序列的第一个元素.

t : array A sequence of time points for which to solve for y. The initial value point should be the first element of this sequence.

因此，您必须从0开始积分或提供从-2开始的初始条件.

So you must integrate starting at 0 or provide initial conditions starting at -2.

最后，您的积分范围包括x = 1/3的奇点. odeint在这里可能会很糟糕(但显然不是).

Finally, your range of integration covers a singularity at x=1/3. odeint may have a bad time here (but apparently doesn't).

这是一种可行的方法:

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

def g(y, x):
    y0 = y[0]
    y1 = y[1]
    y2 = ((3*x+2)*y1 + (6*x-8)*y0)/(3*x-1)
    return y1, y2

# Initial conditions on y, y' at x=0
init = 2.0, 3.0
# First integrate from 0 to 2
x = np.linspace(0,2,100)
sol=odeint(g, init, x)
# Then integrate from 0 to -2
plt.plot(x, sol[:,0], color='b')
x = np.linspace(0,-2,100)
sol=odeint(g, init, x)
plt.plot(x, sol[:,0], color='b')

# The analytical answer in red dots
exact_x = np.linspace(-2,2,10)
exact_y = 2*np.exp(2*exact_x)-exact_x*np.exp(-exact_x)
plt.plot(exact_x,exact_y, 'o', color='r', label='exact')
plt.legend()

plt.show()

这篇关于如何使用Python内置函数odeint求解微分方程?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！