用python解决微分方程

一.用python求解一阶常微分方程

标准形式:

准备工作:

1.首先需要导入scipy库,只需在终端输入pip install scipy然后回车即可

同理导入numpy,matplotlib 包

2. scipy.integrate.odeint()是求解的具体方法,这个函数完整版是这样的

scipy.integrate.odeint(func, y0, t, args=(), Dfun=None, col_deriv=0, full_output=0, ml=None, mu=None, rtol=None, atol=None, tcrit=None, h0=0.0, hmax=0.0, hmin=0.0, ixpr=0, mxstep=0, mxhnil=0, mxordn=12, mxords=5, printmessg=0, tfirst=False)

但是我们只需要用到前三个参数:

func : 导数函数f(y,t),即y在t处的导数,用函数形式表示

y0 :初始条件y0,用数组的形式表示

t : 求解函数值对应的时间点的序列。序列的第一个元素是与初始条件 y0 对应的初始时间 t0;时间序列必须是单调递增或单调递减的,允许重复值。

好啦,现在我们开始正式编程求解:

首先导入库函数:

from scipy.integrate import odeint

import numpy as np

import matplotlib.pyplot as plt 

接着定义函数f(y,t):

def dy_dt(y,t):

    return np.sin(t**2)

使用odeint()方法:

y0=[1]

t = np.arange(-10,10,0.01)   #(start,stop,step)

y=odeint(dy_dt,y0,t)

最后就是绘图啦:

plt.plot(t, y)

plt.title("picture")

plt.show()

完整代码贴一贴:

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

def dy_dt(y,t):
    return np.sin(t**2)

y0=[1]
t = np.arange(-10,10,0.01)
y=odeint(dy_dt,y0,t)

plt.plot(t, y)
plt.title("picture")
plt.show()

得到结果:

二.求解高阶微分方程:

高阶常微分方程,必须做变量替换,化为一阶微分方程组,再用 odeint 求数值解

我们来看这个例子:

直接上代码:

import matplotlib.pyplot as plt
from scipy import linspace,exp
from scipy.integrate import odeint
import numpy as np

def fvdp(t,y):
    '''
    要把y看出一个向量,y = [dy0,dy1,dy2,...]分别表示y的n阶导,那么
    y[0]就是需要求解的函数,y[1]表示一阶导,y[2]表示二阶导,以此类推
    '''
    dy1 = y[1]      # y[1]=dy/dt,一阶导
    dy2 = -4 * y[1] - 5* y[0]  
    # y[0]是最初始,也就是需要求解的函数
    # 注意返回的顺序是[一阶导, 二阶导],这就形成了一阶微分方程组
    return [dy1,dy2] 


def solve_second_order_ode():
    '''
    求解二阶ODE
    '''
    x = linspace(0,20,1000)#给x规定范围
    y0 = [3.0, -5.0] # 初值条件
    # 初值[3.0, -5.0]表示y(0)=3,y'(0)=-5
    # 返回y,其中y[:,0]是y[0]的值,就是最终解,y[:,1]是y'(x)的值
    y = odeint(fvdp, y0, x, tfirst=True)
    
    y1, = plt.plot(x,y[:,0],label='y')
    y1_1, = plt.plot(x,y[:,1],label='y‘')             
    plt.legend(handles=[y1,y1_1])   #创建图例
    
    plt.show()
    
solve_second_order_ode()

 得到结果:

  • 19
    点赞
  • 223
    收藏
    觉得还不错? 一键收藏
  • 5
    评论
评论 5
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值