Question

matlab 一直出错，求助

程序一直出错
为了解一个多元方程，编辑了两个M文件，其中方程6中的一个系数X为随着t变化的曲线上纵坐标值：
In order to solve the equations， I edited 2 M files. In equation (6), X is coefficient which varies with t(t is time)
【ODE451_main】
clc;clear
tspan=[0,180];
y0=[1e-5;0;0;0;0];
[t,y]=ode23('ODE45_fun',tspan,y0);
[m,n]=size(y);
for i=1:m
t=[0 25 50 70 80 100 180];
x=[0 5 7 9 10 11 14];
%plot(t,x,'r*'); hold on;
k = 3; %多项式的阶次
p=polyfit(t,x,k); %得到拟合的多项式系数
%dt = linspace(0,100,100);
%dx = polyval(p,dt);
%plot(dt,dx);
%xlabel('t');
%ylabel('x');
%title('x-t的拟合曲线');
syms t;%syms t Y;
X = [t^3 t^2 t 1]* p'; %得到X-t的函数
z =-178075801.6*y(i,4)*y(i,5)+X*y(i,3)*y(i,2)+1119396.815471592+7481.104004*sin(2.047*t(i)); %得到z(t,Y）的符号表达式

Z = matlabFunction(z); %将符号表达式转化为函数
t = linspace(0,180,180);
y(i,6) = Z(t); %Dz-t的方程;dz = Z(3,t);y(i,6) = Z(t);
end
data=[t,y];
save ODE45_data.txt data -ascii
subplot(2,3,1),plot(t,y(:,1)),title('y(1)')
xlabel('t');ylabel('y');
subplot(2,3,2),plot(t,y(:,2)),title('y(2)')
xlabel('t');ylabel('y');
subplot(2,3,3),plot(t,y(:,3)),title('y(3)')
xlabel('t');ylabel('y');
subplot(2,3,4),plot(t,y(:,4)),title('y(4)')
xlabel('t');ylabel('y');
subplot(2,3,5),plot(t,y(:,5)),title('y(5)')
xlabel('t');ylabel('y');
subplot(2,3,6),plot(t,y(:,6)),title('y(6)')
xlabel('t');ylabel('y');
% plot(t,y(:,1),'bo',t,y(:,2),'rx',t,y(:,3),'gv',t,y(:,4),'r-');
grid on

【ODE45_fun】
function dy=ODE45_fun(t,y)
dy(1)=-1.918298553*y(3)*y(4)-121.6697369*y(5)*y(2)+0.006472085*y(2)*y(2)+15.25250926*y(5)*y(5)-0.518363603*sin(2.047*t)+0.001124759;

dy(2)=0.007229182*y(5)*y(1)-0.013867729*y(3)-0.005151943*sin(2.047*t)+33.43424564*y(4)*y(5)-0.092169794*y(3)*y(2)-0.698266828;

dy(3)=72.10986245*y(2)+0.52129529*y(4)*y(1)+0.025471074*y(3)+921.886526*y(4)/y(1)-0.47870471*y(4)+0.025471074*y(1)*t-0.38220722*cos(2.047*t)-4.62279911;

dy(4)=-57.37263009*y(5)+0.001053501*y(3)*y(1)+20.28825016/y(1)-0.001053501*y(3)+0.064741605*y(4)+20.28825016*t-0.006201915*sin(2.047*t)-3651.885374303154;

dy(5)=0.017284293*y(4)+0.00278644*y(1)*y(2)-0.551218454*y(2)*y(5)*y(5)+0.010839281*y(2)*y(2)*y(5)+0.020353114*cos(2.047*t)+0.110594984;

dy=[dy(1);dy(2);dy(3);dy(4);dy(5)];

Answer 1

【ODE451_main】：

.....

t = linspace(0,180,180);

y(i,6) = Z(t);  %Dz-t的方程;dz = Z(3,t);y(i,6) = Z(t); 

......

----------------------------

1. Z(t)算出来是1行180列的矩阵，而y(i, 6)只是一个元素，左右两边不匹配，无法赋值

2. 从workspace看到，y是一个75行5列的矩阵，而这里却要引用y(i, 6)，下标超出范围

这两个问题你根据自己的需求先做修改吧，改完运行下看哪里还有问题再逐步解决

追问

关键是不会改，呜呜.............

追答

额，要先弄懂原理吧，比如为什么这里要放一句y(i,6) = Z(t); ？