求助,用 MATLAB 的 ode45 求解微分方程组

 我来答
lhmhz
高粉答主

推荐于2017-09-28 · 专注matlab等在各领域中的应用。
lhmhz
采纳数:7263 获赞数:16999

向TA提问 私信TA
展开全部

例如:求解下列微分方程组

求解步骤:

①自定义函数  rigid。m

function dy = rigid(t,y)
dy = zeros(3,1);    % a column vector
dy(1) = y(2) * y(3);
dy(2) = -y(1) * y(3);
dy(3) = -0.51 * y(1) * y(2);

②在执行窗口下,执行下列命令

options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-5]);
[t,Y] = ode45(@rigid,[0 12],[0 1 1],options);

plot(T,Y(:,1),'-',T,Y(:,2),'-.',T,Y(:,3),'.')

③求解结果

Sievers分析仪
2024-10-13 广告
是的。传统上,对于符合要求的内毒素检测,最终用户必须从标准内毒素库存瓶中构建至少一式两份三点标准曲线;必须有重复的阴性控制;每个样品和PPC必须一式两份。有了Sievers Eclipse内毒素检测仪,这些步骤可以通过使用预嵌入的内毒素标准... 点击进入详情页
本回答由Sievers分析仪提供
tianxiawulang
2015-11-24 · TA获得超过2.7万个赞
知道大有可为答主
回答量:4732
采纳率:89%
帮助的人:2593万
展开全部

1、你贴出来的报错信息和代码对不上号:前面显示错误的那行代码和你贴出来的完全不一样;而后面的错误(Input argument 'u1' is undefined)也不可能是目前的代码所导致的——的确是有错,但错误应该是iL未定义才对。

 

2、使用ode*系列函数解常微分方程,用于描述微分方程的函数(例如你这里的current)输入输出参数是有固定格式要求的,应该是

    dy = current( t, y )

其中t是时间,y是t时刻的状态变量。这两二个参数即使你在函数中用不上,也必须列在参数表中。当然,变量的具体名字可以自定,但含义就是上面说的。

 

例如,按照你现在的写法,传递到current函数的u1就是时间,而u2则是微分方程的状态变量,相当于y,是一个向量,有3个元素。如果按照上面的固定格式,current函数不允许有三个输入参数,那么,既然你写了三个输入参数,第三个参数iL自然就没有定义,所以会出错。

 

我猜测,你的方程中u1、u2和iL其实就是状态变量y——如果是这样,那么在current函数中将其分别以y(1)、y(2)、y(3)代替就可以了,后面的大部分内容也就没必要看了。

 

之所以对描述微分方程的函数(例如current)有这种固定的格式要求,是因为 ode*系列函数属于所谓“函数的函数”(Function Functions),也就是说,它的输入参数当中有其它函数(可以是函数文件名,或者函数句柄,也可以是inline函数、匿名函数)。ode* 函数在求解过程中,会反复调用作为参数传递给它的函数,而调用的过程并非由我们决定的,而是固定地写在ode*函数代码中的,所以,对描述微分方程的函数有固定格式要求也就不足为奇了。

 

3、假如描述微分方程的函数除了t和y之外,的确还需要其它数据,应该怎么办?

 

首先,请认真考虑一下问题本身,这些数据到底是什么性质?

  • 如果是常数,可以考虑直接写在current函数里面;

  • 如果变化的,但仅与时间t和状态变量y及其导数有关,也可以直接在current函数里面计算出来,而无需由外部传递。

 

如果上面两种情况都不适合,那么就需要用到传递附加参数了。传递附加参数的常用手段包括使用匿名函数、嵌套函数以及通过函数参数传递三种方式,下面介绍一下第三种方式。

 

ode系列函数较为一般的调用格式为

    [T,Y] =solver(odefun,tspan,y0,options,p1,p2...)

这些参数中,p1、p2等就是要额外传递的参数,而相应的微分方程函数应该定义成

    dy = current( t, y,p1, p2... )

options为求解器选项,如果不知道(同时也不想知道)它有什么用,可以不用管它,调用的时候用空数组([])代替即可。

 

4、还有一点小问题:

    dy=zeros(1,3);
应改成

    dy=zeros(3,1);
描述微分方程的函数要求返回列向量。

 

参考改动

function dy = current( t, y )
u1 = y(1);
u2 = y(2);
iL = y(3);
dy=zeros(3,1);
dy(1)=u1;
dy(2)=u1/(1.96*10^-4)+u2/(1.96*10^-4)+iL/(10^-7);
dy(3)=-u2/(17.2*10^-3)-0.5*iL/(17.2*10^-3);

 

这样修改后,程序可以运行,但求出来的结果是发散的,请你再仔细检查一下方程是否正确。

 

由于不确定你的方程中u1、u2和iL到底是什么,所以只能帮你说到这里了。

已赞过 已踩过<
你对这个回答的评价是?
评论 收起
推荐律师服务: 若未解决您的问题,请您详细描述您的问题,通过百度律临进行免费专业咨询

为你推荐:

下载百度知道APP,抢鲜体验
使用百度知道APP,立即抢鲜体验。你的手机镜头里或许有别人想知道的答案。
扫描二维码下载
×

类别

我们会通过消息、邮箱等方式尽快将举报结果通知您。

说明

0/200

提交
取消

辅 助

模 式