我试图用欧拉方法求解一个微分方程组,首先,我将系统编码成一个向量,然后将初始条件传递给函数ode_Euler。
但是,我的尝试出现了一些错误。我收到以下错误:
>> nm06p03a
Unable to perform assignment because the size of the left side is 1-by-2 and the size of the right side is 2-by-2.
Error in ode_Euler (line 11)
y(k+1,:)= y(k,:) +h*feval(f,t(k),y(k,:));
Error in nm06p03a (line 12)
tic, [tE,xE]=ode_Euler(f,tspan,x0,N); t_Euler=toc;
这是我的代码:
clear, clf
f=@(t,x)[-x(2)+x(1)*x(2); x(1)-(0.5.*x(1).^2)+(0.5.*x(2).^2)]; %Encoding system of differential equations into a vector
t0=0; tf=10;
tspan=[t0 tf];
N=100;
x0s=[0.2 0]; % A matrix consisting of initial points
for iter=1:size(x0s,1)
x0=x0s(iter,:);
tic, [tE,xE]=ode_Euler(f,tspan,x0,N); t_Euler=toc;
subplot(220+iter),
plot(tE,xE,'r:')
legend('ode_ Euler')
end
下面是欧拉方法:
function [t,y]=ode_Euler(f,tspan,y0,N)
if nargin<4|N<=0, N=100; end
if nargin<3, y0=0; end
h=(tspan(2)-tspan(1))/N;
t=tspan(1)+[0:N]'*h;
y(1,:)=y0(:)'; %make it a row vector
for k=1:N
y(k+1,:)= y(k,:) +h*feval(f,t(k),y(k,:));
end
当我使用另一个方法ode_Heun
时,我得到了同样的错误:
function [t,y]=ode_Heun(f,tspan,y0,N)
if nargin<4|N<=0, N=100; end
if nargin<3, y0=0; end
h=(tspan(2)-tspan(1))/N; % Step-size
t=tspan(1)+[0:N]'*h; % Time vector
y(1,:)=y0(:)'; % make the initial value a row vector
for k=1:N
fk= feval(f,t(k),y(k,:)); y(k+1,:)= y(k,:)+h*fk; % Eq.(6.2.3)
y(k+1,:)= y(k,:) +h/2*(fk +feval(f,t(k+1),y(k+1,:))); % Eq.(6.2.4)
end
我可以得到一些帮助来了解我的代码问题吗?
1条答案
按热度按时间t1rydlwq1#
y(k,:)
是一个行向量,而f
的返回值是一个列向量。根据广播规则,行向量和列向量之和是一个矩阵,即重复的行向量和列向量的矩阵之和。这在向量和矩阵运算的上下文中不是很合乎逻辑,但是在处理向量的(有限)序列时是有意义的,不幸的是,这种区别在类型系统中没有实现和实施。