文档建议反馈控制台
首页
学习
活动
专区
圈层
工具
MCP广场
文章/答案/技术大牛
发布
首页
学习
活动
专区
圈层
工具
MCP广场
返回腾讯云官网
社区首页 >问答首页 >无法使用ODE45求解器解决ODE问题

问无法使用ODE45求解器解决ODE问题
EN

Stack Overflow用户
提问于 2022-08-01 17:14:19
回答 1查看 36关注 0票数 0

这是我的原始密码。

代码语言：javascript
运行
复制
function dydt = Biogearsmcdaniel(t,y)
SP= 0.60;
kPT=3.7*10^4;
thetaP =1.35*10^-4;
kPTB= 3.1;
kPTMT=6.3*10^-3;
kPTNT=6.1*10^-4;
SMT=2.6*10^-2;
Mv=0.3;
kMT=6.43*10^-5;
kMTB=36.0;
SNT=7.0*10^-7;
Nv =1*10^8;
kNTB=36.0;
kNTMT=0.16;
kNT=6.1*10^-2;
SB =4.6*10^-2;
kBPT=26.0;
kBR= 0.14;
kBNT=4.0*10^-8;
kPBNA=5.8;
xPN= 0.5;
kPS =6.9*10^3;
xPS=1.3*10^4;
kMP=1.01;
xMP =-37.5;
kMD=5.0*10^-2;
xMD = 0.75;
xMTNF=0.4;
kM6=0.1;
xM6 =1.0;
xM10 =0.297;
kMR=0.05;
SM= 1.0;
kMA=0.2;
kNP= 33.75;
xNP=56.25;
kND= 0.05;
xND =0.4;
kNTNF= 0.2;
xENOSP =1.015; 
k10MA=0.1;
xNTNF= 2.0; 
kENOS= 4.0;
k10TNF= 1.485;
kN6 =1.5;
kNO3= 0.46;
x10TNF=0.05;
xN6 =1.0;
kNOMN= 2.0; 
k106= 5.1*10^-2;
xN10= 0.2; 
kTNFN =2.97;
x106 =8.0*10^-2;
kNR=0.05;
kTNFM=0.1;
k10R= 0.1;
SN =1.0; 
xTNF10= 7.9*10^-2;
x1012=1.0*10^-2;
kNA=0.5;
xTNF6= 5.9*10^-2;
k10= 0.35;
kINOSN =1.5;
kTNF= 1.4;
S10= 1.0*10^-2;
kINOSM =0.1;
k6N= 0.2;
k12M= 0.303;
kINOSEC= 0.1;
k6M= 3.03;
x1210= 0.2525;
xINOSTNF =0.05;
k6TNF= 1.0;
k12= 5.0*10^-2;
kINOSd= 0.05;
x6TNF= 0.1; 
kD= 0.15;
kINOS6= 2.0;
k6NO= 2.97;
kD6= 0.125;
xINOS6 =0.1;
x6NO= 0.4; 
xD6= 0.85;
xINOS10= 0.1;
x610= 0.1782;
xDNO= 0.5;
xINOSNO= 0.3;
x66= 0.5;
kINOS= 0.101;
k6= 0.7;
kENOSEC= 0.05;
xENOSTNF= 0.4;
k10N=0.1;
S6= 1.0*10^-3;

dydt=zeros(19,1);%DAEs 

function f=HU1(x,n,h)

    f=x^h/1+(x/n)^h;

end

function g=HU2(x,n,h)

    g=x^h/(x^h+n^h);

end

function r=HUD(x,n,h)

    r=1/1+(x/n)^h;

end

R=1;

dydt(1)=(SP/kPT)*y(1)*(1-y(1))-(thetaP*y(1)/1+kPTB*y(1))-(kPTMT*y(2)*y(1))-kPTNT*y(3)*y(1);
dydt(2)=SMT*y(3)*Mv/1+kMTB*y(4)-kMT*y(2);
dydt(3)=SNT*R*Nv/((1+kNTB*y(4))*(1+kNTMT*y(2)))-kNT*y(3);
dydt(4)=(SB/1+kBPT)*y(4)*(1-y(4))-kBR*R*y(4)-kBNT*y(3)*y(4);
dydt(5)=SP*y(5)+thetaP*y(1)/1+kPTB*y(1)-kPS*y(5)/xPS+y(5)-kPBNA*y(9)*HU2(y(5),xPN,2);
dydt(6)=(-(kMP*HU2(y(5),xMP,2)+kMD*HU2(1-y(18),xMD,4)*(HU2(y(14),xMTNF,2)+kM6*HU2(y(15),xM6,2)))*HUD(y(16),xM10,2)*y(6)-kMR*(y(6)-SM));
dydt(7)=(kMP*HU2(y(5),xMP,2)+kMD*HU2(1-y(18),xMD,4))*(HU2(y(14),xMTNF,2)+kM6*HU2(y(15),xM6,2))*HUD(y(16),xM10,2)*y(6)-kMA*y(7);
dydt(8)=-((kNP*HU2(y(5),xNP,1))+kND*HU1(1-y(18),xND,2)+kNTNF*HU1(y(14),xNTNF,1)+kN6*HU1(y(15),xN6,2))*HUD(y(16),xN10,2)*y(8)-kNR*(y(8)-SN);
dydt(9)=(kNP*HU2(y(5),xNP,1))+kND*HU1(1-
y(18),xND,2)+kNTNF*HU1(y(14),xNTNF,1)+kN6*HU1(y(15),xN6,2)*HUD(y(16),xN10,2)*y(8)-kNA*y(9);
dydt(10=(kINOSN*y(9)+kINOSM*y(7)+kINOSEC*HU1(y(14),xINOSTNF,2)+kINOS6*HU1(y(15),xINOS6,2))*HUD(y(16),xINOS10,2)*HUD(y(19),xINOSNO,4)-kINOSd*y(10);
dydt(11)=kINOS*(y(10)-y(11));
dydt(12)=kENOSEC*HUD(y(14),xENOSTNF,1)*HUD(y(5),xENOSP,1)-kENOS*y(12);
dydt(13)=kNO3*(y(19)-y(13));
dydt(14)=(kTNFN*y(9)+kTNFM*y(7))*HUD(y(16),xTNF10,2)*HUD(y(15),xTNF6,3)-kTNF*y(14);
dydt(15)=(k6N*y(9)+y(7))*(k6M+k6TNF*HU2(y(14),x6TNF,2)+k6NO*HU2(y(19),x6NO,2))*HUD(y(16),x610,2)*HUD(y(16),x66,1)-k6*(y(15)-S6);
dydt(16)=(k10N*y(9)+y(7))*(k10MA+k10TNF*HU2(y(14),x10TNF,4)+k106*HU2(y(15),x106,4))*((1-k10R)*HUD(y(17),x1012,4)+k10R)-k10*(y(16)-S10);
dydt(17)=k12M*y(7)*HUD(y(16),x1210,2)-k12*y(17);
dydt(18)=kD*(1-y(18))*(y(18)-0.05)-(y(18)-0.05)*kD6*HU2(y(15),xD6,6)*(1/xDNO^2+y(19)^2);
dydt(19)=y(11)*(1+kNOMN*(y(7)+y(9)))+y(12);

end

下面的代码是我为代码的求解者编写的语法。

代码语言：javascript
运行
复制
tspan = 0:120;

yo=[1 1000 1000 0 0 1000 0 1000 0 0 0 0 0 1.1 1.02 1.05 1.2 0 0];

[t,y]=ode15s(@(t,y)Biogearsmcdaniel,tspan,yo);

for i=1:19

    figure(i) 

    plot(t,y(:,i))

    hold on 

end

当我运行上面的文件没有足够的输入参数时，我会得到以下错误。

Biogearsmcdaniel (第109行) dydt(1)=(SP/kPT)y(1)(1-y(1))-(thetaP_y(1)/1+kPTB_y(1))-(kPTMT*y(2)_y(1))-kPTNT_y(3)*y(1)；的

误差

untitled148>@(t，y)Biogearsmcdaniel (第3行) t，y=ode15s(@(t，y)Biogearsmcdaniel，tspan，yo)中的误差；

f0 = ode(t0，y0，args{:})；% ODE15I将args{1}设置为yp0。

错误在ode15s (第153行)的参数(odeIsFuncHandle，odeTreatAsMFile，solver_name，ode，tspan，y0，options，varargin)；

untitled148 (第3行) t，y=ode15s(@(t，y)Biogearsmcdaniel，tspan，yo)中的误差；

请帮帮我。

matlab
EN

回答 1

Stack Overflow用户

发布于 2022-08-06 17:41:33

我是John BG jgb2012@sky.com

我修补了您的启动脚本，使其运行时没有错误。

但是，不确定所获得的结果是否是您要寻找的结果。

​

​

启动并运行脚本的步骤：

1.-您没有正确地输入ode15s，在告诉函数的名称时，必须定义该函数的输入。

2.-我把每个函数放在一个不同的文件中

3.-我已使大量参数全局化，在主支持函数之外可能没有必要，但如果从这里开始，您可以采取

4.-以下是不崩溃的主要支持函数的代码

代码语言：javascript
运行
复制
function dydt = Biogearsmcdaniel(t,y)
global SP
global kPT
global thetaP
global kPTB
global kPTMT
global kPTNT
global SMT
global Mv
global kMT
global kMTB
global SNT
global Nv
global kNTB
global kNTMT
global kNT
global SB
global kBPT
global kBR
global kBNT
global kPBNA
global xPN
global kPS
global xPS
global kMP
global xMP
global kMD
global xMD
global xMTNF
global kM6
global xM6
global xM10
global kMR
global SM
global kMA
global kNP
global xNP
global kND
global xND
global kNTNF
global xENOSP
global k10MA
global xNTNF
global kENOS
global k10TNF
global kN6global
global kNO3
global x10TNF
global xN6global
global kNOMN
global k106
global xN10
global kTNFNglobal
global x106global 
global kNR
global kTNFM
global k10R
global SN
global xTNF10
global x1012
global kNA
global xTNF6
global k10
global kINOSN
global kTNF
global S10
global kINOSM
global k6N
global k12M
global kINOSEC
global k6M
global x1210
global xINOSTNF
global k6TNF
global k12
global kINOSd
global x6TNF
global kD
global kINOS6
global k6NO
global kD6
global xINOS6
global x6NO
global xD6
global xINOS10
global x610
global xDNO
global xINOSNO
global x66
global kINOS
global k6
global kENOSEC
global xENOSTNF
global k10N
global S6

SP= 0.60;
kPT=3.7*10^4;
thetaP =1.35*10^-4;
kPTB= 3.1;
kPTMT=6.3*10^-3;
kPTNT=6.1*10^-4;
SMT=2.6*10^-2;
Mv=0.3;
kMT=6.43*10^-5;
kMTB=36.0;
SNT=7.0*10^-7;
Nv =1*10^8;
kNTB=36.0;
kNTMT=0.16;
kNT=6.1*10^-2;
SB =4.6*10^-2;
kBPT=26.0;
kBR= 0.14;
kBNT=4.0*10^-8;
kPBNA=5.8;
xPN= 0.5;
kPS =6.9*10^3;
xPS=1.3*10^4;
kMP=1.01;
xMP =-37.5;
kMD=5.0*10^-2;
xMD = 0.75;
xMTNF=0.4;
kM6=0.1;
xM6 =1.0;
xM10 =0.297;
kMR=0.05;
SM= 1.0;
kMA=0.2;
kNP= 33.75;
xNP=56.25;
kND= 0.05;
xND =0.4;
kNTNF= 0.2;
xENOSP =1.015; 
k10MA=0.1;
xNTNF= 2.0; 
kENOS= 4.0;
k10TNF= 1.485;
kN6 =1.5;
kNO3= 0.46;
x10TNF=0.05;
xN6 =1.0;
kNOMN= 2.0; 
k106= 5.1*10^-2;
xN10= 0.2; 
kTNFN =2.97;
x106 =8.0*10^-2;
kNR=0.05;
kTNFM=0.1;
k10R= 0.1;
SN =1.0; 
xTNF10= 7.9*10^-2;
x1012=1.0*10^-2;
kNA=0.5;
xTNF6= 5.9*10^-2;
k10= 0.35;
kINOSN =1.5;
kTNF= 1.4;
S10= 1.0*10^-2;
kINOSM =0.1;
k6N= 0.2;
k12M= 0.303;
kINOSEC= 0.1;
k6M= 3.03;
x1210= 0.2525;
xINOSTNF =0.05;
k6TNF= 1.0;
k12= 5.0*10^-2;
kINOSd= 0.05;
x6TNF= 0.1; 
kD= 0.15;
kINOS6= 2.0;
k6NO= 2.97;
kD6= 0.125;
xINOS6 =0.1;
x6NO= 0.4; 
xD6= 0.85;
xINOS10= 0.1;
x610= 0.1782;
xDNO= 0.5;
xINOSNO= 0.3;
x66= 0.5;
kINOS= 0.101;
k6= 0.7;
kENOSEC= 0.05;
xENOSTNF= 0.4;
k10N=0.1;
S6= 1.0*10^-3;


dydt=zeros(19,1);%DAEs 

% function f=HU1(x,n,h)
% these functions in separate files
%     f=x^h/1+(x/n)^h;
% 
% end

% function g=HU2(x,n,h)
% 
%     g=x^h/(x^h+n^h);
% 
% end

% function r=HUD(x,n,h)
% 
%     r=1/1+(x/n)^h;
% 
% end

R=1;

dydt(1)=(SP/kPT)*y(1)*(1-y(1))-(thetaP*y(1)/1+kPTB*y(1))-(kPTMT*y(2)*y(1))-kPTNT*y(3)*y(1);
dydt(2)=SMT*y(3)*Mv/1+kMTB*y(4)-kMT*y(2);
dydt(3)=SNT*R*Nv/((1+kNTB*y(4))*(1+kNTMT*y(2)))-kNT*y(3);
dydt(4)=(SB/1+kBPT)*y(4)*(1-y(4))-kBR*R*y(4)-kBNT*y(3)*y(4);
dydt(5)=SP*y(5)+thetaP*y(1)/1+kPTB*y(1)-kPS*y(5)/xPS+y(5)-kPBNA*y(9)*HU2(y(5),xPN,2);
dydt(6)=(-(kMP*HU2(y(5),xMP,2)+kMD*HU2(1-y(18),xMD,4)*(HU2(y(14),xMTNF,2)+kM6*HU2(y(15),xM6,2)))*HUD(y(16),xM10,2)*y(6)-kMR*(y(6)-SM));
dydt(7)=(kMP*HU2(y(5),xMP,2)+kMD*HU2(1-y(18),xMD,4))*(HU2(y(14),xMTNF,2)+kM6*HU2(y(15),xM6,2))*HUD(y(16),xM10,2)*y(6)-kMA*y(7);
dydt(8)=-((kNP*HU2(y(5),xNP,1))+kND*HU1(1-y(18),xND,2)+kNTNF*HU1(y(14),xNTNF,1)+kN6*HU1(y(15),xN6,2))*HUD(y(16),xN10,2)*y(8)-kNR*(y(8)-SN);
dydt(9)=(kNP*HU2(y(5),xNP,1))+kND*HU1(1-y(18),xND,2)+kNTNF*HU1(y(14),xNTNF,1)+kN6*HU1(y(15),xN6,2)*HUD(y(16),xN10,2)*y(8)-kNA*y(9);
dydt(10)=(kINOSN*y(9)+kINOSM*y(7)+kINOSEC*HU1(y(14),xINOSTNF,2)+kINOS6*HU1(y(15),xINOS6,2))*HUD(y(16),xINOS10,2)*HUD(y(19),xINOSNO,4)-kINOSd*y(10);
dydt(11)=kINOS*(y(10)-y(11));
dydt(12)=kENOSEC*HUD(y(14),xENOSTNF,1)*HUD(y(5),xENOSP,1)-kENOS*y(12);
dydt(13)=kNO3*(y(19)-y(13));
dydt(14)=(kTNFN*y(9)+kTNFM*y(7))*HUD(y(16),xTNF10,2)*HUD(y(15),xTNF6,3)-kTNF*y(14);
dydt(15)=(k6N*y(9)+y(7))*(k6M+k6TNF*HU2(y(14),x6TNF,2)+k6NO*HU2(y(19),x6NO,2))*HUD(y(16),x610,2)*HUD(y(16),x66,1)-k6*(y(15)-S6);
dydt(16)=(k10N*y(9)+y(7))*(k10MA+k10TNF*HU2(y(14),x10TNF,4)+k106*HU2(y(15),x106,4))*((1-k10R)*HUD(y(17),x1012,4)+k10R)-k10*(y(16)-S10);
dydt(17)=k12M*y(7)*HUD(y(16),x1210,2)-k12*y(17);
dydt(18)=kD*(1-y(18))*(y(18)-0.05)-(y(18)-0.05)*kD6*HU2(y(15),xD6,6)*(1/xDNO^2+y(19)^2);
dydt(19)=y(11)*(1+kNOMN*(y(7)+y(9)))+y(12);

end

5.-这是打电话的人

代码语言：javascript
运行
复制
clear all;clc;close all

global SP
global kPT
global thetaP
global kPTB
global kPTMT
global kPTNT
global SMT
global Mv
global kMT
global kMTB
global SNT
global Nv
global kNTB
global kNTMT
global kNT
global SB
global kBPT
global kBR
global kBNT
global kPBNA
global xPN
global kPS
global xPS
global kMP
global xMP
global kMD
global xMD
global xMTNF
global kM6
global xM6
global xM10
global kMR
global SM
global kMA
global kNP
global xNP
global kND
global xND
global kNTNF
global xENOSP
global k10MA
global xNTNF
global kENOS
global k10TNF
global kN6global
global kNO3
global x10TNF
global xN6global
global kNOMN
global k106
global xN10
global kTNFNglobal
global x106global 
global kNR
global kTNFM
global k10R
global SN
global xTNF10
global x1012
global kNA
global xTNF6
global k10
global kINOSN
global kTNF
global S10
global kINOSM
global k6N
global k12M
global kINOSEC
global k6M
global x1210
global xINOSTNF
global k6TNF
global k12
global kINOSd
global x6TNF
global kD
global kINOS6
global k6NO
global kD6
global xINOS6
global x6NO
global xD6
global xINOS10
global x610
global xDNO
global xINOSNO
global x66
global kINOS
global k6
global kENOSEC
global xENOSTNF
global k10N
global S6

SP= 0.60;
kPT=3.7*10^4;
thetaP =1.35*10^-4;
kPTB= 3.1;
kPTMT=6.3*10^-3;
kPTNT=6.1*10^-4;
SMT=2.6*10^-2;
Mv=0.3;
kMT=6.43*10^-5;
kMTB=36.0;
SNT=7.0*10^-7;
Nv =1*10^8;
kNTB=36.0;
kNTMT=0.16;
kNT=6.1*10^-2;
SB =4.6*10^-2;
kBPT=26.0;
kBR= 0.14;
kBNT=4.0*10^-8;
kPBNA=5.8;
xPN= 0.5;
kPS =6.9*10^3;
xPS=1.3*10^4;
kMP=1.01;
xMP =-37.5;
kMD=5.0*10^-2;
xMD = 0.75;
xMTNF=0.4;
kM6=0.1;
xM6 =1.0;
xM10 =0.297;
kMR=0.05;
SM= 1.0;
kMA=0.2;
kNP= 33.75;
xNP=56.25;
kND= 0.05;
xND =0.4;
kNTNF= 0.2;
xENOSP =1.015; 
k10MA=0.1;
xNTNF= 2.0; 
kENOS= 4.0;
k10TNF= 1.485;
kN6 =1.5;
kNO3= 0.46;
x10TNF=0.05;
xN6 =1.0;
kNOMN= 2.0; 
k106= 5.1*10^-2;
xN10= 0.2; 
kTNFN =2.97;
x106 =8.0*10^-2;
kNR=0.05;
kTNFM=0.1;
k10R= 0.1;
SN =1.0; 
xTNF10= 7.9*10^-2;
x1012=1.0*10^-2;
kNA=0.5;
xTNF6= 5.9*10^-2;
k10= 0.35;
kINOSN =1.5;
kTNF= 1.4;
S10= 1.0*10^-2;
kINOSM =0.1;
k6N= 0.2;
k12M= 0.303;
kINOSEC= 0.1;
k6M= 3.03;
x1210= 0.2525;
xINOSTNF =0.05;
k6TNF= 1.0;
k12= 5.0*10^-2;
kINOSd= 0.05;
x6TNF= 0.1; 
kD= 0.15;
kINOS6= 2.0;
k6NO= 2.97;
kD6= 0.125;
xINOS6 =0.1;
x6NO= 0.4; 
xD6= 0.85;
xINOS10= 0.1;
x610= 0.1782;
xDNO= 0.5;
xINOSNO= 0.3;
x66= 0.5;
kINOS= 0.101;
k6= 0.7;
kENOSEC= 0.05;
xENOSTNF= 0.4;
k10N=0.1;
S6= 1.0*10^-3;

tspan = [0:120];

yo=[1 1000 1000 0 0 1000 0 1000 0 0 0 0 0 1.1 1.02 1.05 1.2 0 0];


figure(1)
hold all

for i=1:19

%     figure(i) 

[t,y]=ode15s(@(t,y) Biogearsmcdaniel(t,yo),tspan,yo);

    plot(t,y(:,i))

    hold on 

end

grid on

6.-没有必要在for循环中重复等待，只需重叠每个生成的图

是的，所有的重复，太多的行，但这是一个很好的做法，尽快发布的脚本工作；过度工作的脚本往往回到最初解决的错误。

希望它能帮上忙

约翰·BG jgb2012@sky.com

票数 0
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：

https://stackoverflow.com/questions/73197325

复制
相关文章

相似问题

领券

腾讯云开发者

扫码关注腾讯云开发者

扫码关注腾讯云开发者

领取腾讯云代金券

热门产品

热门推荐

更多推荐

Copyright © 2013 - 2025 Tencent Cloud. All Rights Reserved. 腾讯云 版权所有 

深圳市腾讯计算机系统有限公司 ICP备案/许可证号：粤B2-20090059 深公网安备号 44030502008569

腾讯云计算（北京）有限责任公司 京ICP证150476号 |  京ICP备11018762号 | 京公网安备号11010802020287

问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档