matlab解常微分方程组数值解法(二元常微分方程组的解法)

上篇博客介绍了Matlab求解常微分方程组解析解的方法:博客地址微分方程组复杂时,无法求出解析解时,就需要求其数值解,这里来介绍。以下内容按照Matlab官方文档提供的方程来展开(提议多看官方文档)介绍一下核心函数ode45()一般形式:[t,y]=ode45(odefun,tspan,y0) 其中tspan=[t0tf]功能介绍:求微分方程组y′=f(t,y)从t0…

大家好,又见面了,我是你们的朋友全栈君。

上篇博客介绍了Matlab求解常微分方程组解析解的方法:博客地址
微分方程组复杂时,无法求出解析解时,就需要求其数值解,这里来介绍。
以下内容按照Matlab官方文档提供的方程来展开(提议多看官方文档)

介绍一下核心函数ode45()

一般形式:[t,y] = ode45(odefun,tspan,y0) 其中 tspan = [t0 tf]
功能介绍:求微分方程组 y′=f(t,y) 从 t0 到 tf 的积分,初始条件为 y0。解数组 y 中的每一行都与列向量 t 中返回的值相对应。

1. 一阶微分方程求解(简单调用即可)

方程:y’=2*t
代码

tspan=[1 6]; %定义自变量x的取值空间为1-6
y0=0;%定义因变量的初值,当x=1(x取值空间的第一个数)时,y0=0
[t,y]=ode45(@(t,y) 2*t,tspan,y0); %定义函数y'=2*t,使用ode45求解
plot(t,y,'-o'); %绘制求得的数值曲线

说明:简单的odefun参数就是这个形式,@(x,y) fun 中间是空格,不能用逗号
结果
在这里插入图片描述

2. 二阶微分方程求解(引入函数文件)

方程:范德波尔方程 y1’’-u(1-y1²)*y1’+y1=0;(这里设u=1)
代码
定义输入的方程,以函数形式定义

function dydt=odefun(t,y)
%二阶方程为y1''-(1-y1²)*y1'+y1=0;
%降阶为两个方程:y1'=y2; 
%               y2'=(1-y1²)*y2-y1;
%t虽然没有使用,但必须要作为参数写入
dydt=[y(2);(1-y(1)^2)*y(2)-y(1)];

求解作图

tspan=[0 20]; %定义自变量x的取值空间为0-20
y0=[2;0];%定义因变量的初值,当x=0时,y1=2,y2=y1'=0;
[t,y]=ode45(@odefun,tspan,y0); %使用ode45求解
%%下面为作图过程,不解释
plot(t,y(:,1),'-o',t,y(:,2),'-o')
title('Solution of van der Pol Equation (\mu = 1) with ODE45');
xlabel('Time t');
ylabel('Solution y');
legend('y_1','y_2')

说明:先进行降阶,使y2=y1’,这时候y2’=y1’’,然后把两个方程写成列向量形式就行了
结果
在这里插入图片描述

3. 求解微分方程组(和2类似)

这里就和求解二阶方程类似的,只不过不需要降阶,仍旧需要一个函数来定义方程组。我们这里不用官方文档的例子,用同学的循坏摆问题来进行演示。
方程在这里插入图片描述
给定的初值(w接近0,但实际上不能设置为0):
在这里插入图片描述
代码
定义输入的方程

function dRvw=func(t,Rvw)
%% 函数功能:为ode45提供微分方程
%输入:t:时间序列,就是θ;Rvw:因变量,Rvw(1)代表R,Rvw(2)代表v,Rvw(3)代表w
%输出:dRvw:因变量的一阶微分,dRvw(1)代表dR,dRvw(2)代表dv,dRvw(3)代表dw
%% 初始化因变量的一阶微分,3×1的向量
dRvw=zeros(3,1);
%% 参数初始化
r=0.01;u=0.1;g=9.8;M=10;m=1;
%% 输入微分方程式
dRvw(1)=-Rvw(2)/Rvw(3)-r;
dRvw(2)=(M*g-(m*Rvw(3)^2*Rvw(1)+m*g*sin(t)) *exp(u*(t+pi/2)))/(Rvw(3)*(M+m));
dRvw(3)=Rvw(3)*r/Rvw(1)+2*Rvw(2)/Rvw(1)+g*cos(t)/(Rvw(3)*Rvw(1));

求解作图

%% 参数初始化 start_Theta是θ的初始值 end_Theta是θ的结束值
%R是半径初值;v是线速度初值;w是角速度初值
start_Theta=0;end_Theta=2*pi;R=1;v=0;w=1e-5;
%% 使用ode45方法计算微分方程组func的数值解
%func是带有方程组的函数 
%[start_Theta end_Theta]是自变量范围
%[R;v;w]是方程初值
%T是自变量的数组,Rvw是对应的因变量的数值。Rvw(:,1)=R;Rvw(:,2)=v;Rvw(:,3)=w;
[T,Rvw]=ode45(@func,[start_Theta end_Theta],[R;v;w]);
%% 组合自变量和因变量。TRvw(:,1)=θ;TRvw(:,2)=R;TRvw(:,3)=v;Trvw(:,4);
TRvw=[T,Rvw];
%% 绘制自变量-因变量图像.figure1是R-θ,figure2是v-θ,figure3是w-θ
plot(T,Rvw(:,1),'-o',T,Rvw(:,2),'-o',T,Rvw(:,3),'-o')
title('自变量-因变量图像');
xlabel('Angle θ');
ylabel('Solution');
legend('R','v','w')

说明:注释的应该是比较清楚的,把三个方程写成列向量的形式就行

PS:有些人和我说不能运行,然后我看了他们出错的原因,有点儿哭笑不得。出错的基本上都是运行上面的dRvw=func(t,Rvw)这个函数的。说明一下,这是有参数的函数,不给参数不能直接运行的。下面的求解作图脚本才是需要运行的哈,它调用了函数,才得到的结果。如果大家还发现什么问题,欢迎私戳或评论。

PS+ 有了PS之后,还是很多人问我参数的问题,我在这里直接把文件给大家:cupt.zip
提取码:6k8n。
解压后,有两个文件,大家直接运行cupt.m即可。

结果
在这里插入图片描述

4. 更多形式

讲到这里,大部分我们用到的微分方程形式都可以求解了,Matlab还支持带有时变项和额外参数的微分方程求解,这里不再赘述,大家可以自行参阅官方文档。

版权声明:本文内容由互联网用户自发贡献,该文观点仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 举报,一经查实,本站将立刻删除。

发布者:全栈程序员-用户IM,转载请注明出处:https://javaforall.cn/126077.html原文链接:https://javaforall.cn

【正版授权,激活自己账号】: Jetbrains全家桶Ide使用,1年售后保障,每天仅需1毛

【官方授权 正版激活】: 官方授权 正版激活 支持Jetbrains家族下所有IDE 使用个人JB账号...

(2)
blank

相关推荐

  • from_unixtime函数类型_localtime_r函数

    from_unixtime函数类型_localtime_r函数Unix时间戳(Unixtimestamp),是一种时间表示方式,定义为从格林威治时间1970年01月01日00时00分00秒起至现在的总秒数。在MySQL中如何格式化时间戳?在mysql中因为

  • 利用Nginx构建负载均衡server[通俗易懂]

    利用Nginx构建负载均衡server

  • OpenSSL密码库算法笔记——第5.1.2章 椭圆曲线算法集

    OpenSSL密码库算法笔记——第5.1.2章 椭圆曲线算法集在定义椭圆曲线点群时出现了描述曲线所用算法的参数constEC_METHOD*meth,这一节就来看看这个参数有什么用处。椭圆曲线算法集的定义如下。typedefstructec_method_stEC_METHOD;structec_method_st{(具体定义略,详情可参见代码文件ec_lcl.h)};//EC_METHOD在结构体ec…

  • 操作系统虚拟存储技术_虚拟存储

    操作系统虚拟存储技术_虚拟存储虚拟存储管理   在前面总结了集中存储管理的刚上,要求作业的逻辑地址空间连续的存放主存储器的某个区域中。当主存储器中没有足够大的区域是,则作业是无法装入的,或必须移动某些作业后才能装入。是否有可能吧作业的连续逻辑地址空间分散到几个不连续的主存区域,且仍能使作业正确执行呢?若可行的话,则可充分利用主存空间有可减少移动所花费的开销。不仅如此,还可采用虚拟存储管理技

  • for循环break和continue[通俗易懂]

    for循环break和continue[通俗易懂]for循环像while循环一样,for可以完成循环的功能。在Python中for循环可以遍历任何序列的项目,如一个列表或者一个字符串等。for循环的格式for临时变量in列表或者字符串等:循环满足条件时执行的代码demo1name=‘itheima’forxinname:print(x)运行结果如下:itheimademo2name=‘h…

  • 微信开放平台网站应用扫码登陆[通俗易懂]

    微信开放平台网站应用扫码登陆[通俗易懂]cookie问题现象:接收微信回调数据后,后端存储cookie,未超出cookie有效期,但cookie在浏览器关闭页面重新打开后消失原因:在微信开放平台填写授权回调域名为xx.com,但在浏览器重新打开页面时,域名默认为www.xx.com,而www.xx.cn无法获取到xx.com的cookie信息,看起来就像cookie消失了一样解决办法:在微信开放平台填写授权回调域名为www….

发表回复

您的电子邮箱地址不会被公开。

关注全栈程序员社区公众号