What is wrong with this code? Need help.

1 次查看（过去 30 天）
ZAINULABADIN ZAFAR
ZAINULABADIN ZAFAR 2019-8-31
评论： ZAINULABADIN ZAFAR 2019-9-1
Main code:
clear all
close all
clc
M = [1 0 0
0 1 0
0 0 1];
x0=[-0.5; 0.2; 0.5];
options = odeset('Mass',M,'RelTol',1e-12,'AbsTol',[1e-14 1e-14 1e-14], 'Vectorized','on');
global t x y z dt alpha
dt=0.01;
for alpha=0.8800:0.0005:1.6000
alpha
clear n
clear m
[t,x]=ode15s(@equations,0:dt:500,x0);
n=length(x(:,1));
m=floor(n/2);
y=diff(x(m,n,1))/dt; %% Error Message at this line: Index exceeds matrix dimensions
z=diff(y)/dt;
k=1;
clear aa
for i=m:n
t0=t(i); %Comp. local max. pts for m<t<n
option=optimset('display','off');
zer=fsolve(@differ,t0,option);
if interp1(t(m:n-2),z,zer)>0
aa(k)=interp1(t(m:n),x(m:n,1),zer);
k=k+1;
end
end
kmax=k-1;
h=plot(alpha.*ones(1,kmax),aa,'r.');
hold on
set(h,'MarkerSize',0.1);
end
Function equations.m:
function xdot=equations(t,x)
global alpha
a=0.0005; b=0.01; eps=0.01; beta=-1; R=0.3;
xdot(1) = (-x(2) + alpha*x(1)^2 + beta*x(1)^3)/eps;
xdot(2) = x(1) - x(3) - R*x(2);
xdot(3) = a - b*x(2);
xdot = xdot';
end
Function differ.m:
function f=differ(a)
global t x dt
s=diff(x(m:n,1))/dt;
f=interp1(t(m:n-1),s,aa);
end
The author of this program plot the attached picture. But when I tried cannot get that plot. Please help!!!!

请先登录，再回答此问题。

回答（1 个）

KALYAN ACHARJYA
KALYAN ACHARJYA 2019-8-31
y=diff(x(m,n,1))/dt; %% Error Message at this line: Index exceeds matrix dimensions
Yes, see the reason:
>> m
m =
25000
>> n
n =
50001
>> whos x
Name Size Bytes Class Attributes
x 50001x3 1200024 double global
The size of x is 5000x1 and you are trying to acess of x as follows (In first iterations)
x(25000,50001,1)
%..^.........^ no x values
Thats why it shows the error index exceeds matrix dimensions. That means you trying to acess the indices more than it have.
  3 个评论
显示 1更早的评论
Walter Roberson
Walter Roberson 2019-9-1
Considering your function differ() uses
s=diff(x(m:n,1))/dt;
then I speculate that your earlier line
y=diff(x(m,n,1))/dt;
should be
y=diff(x(m:n,1))/dt;
ZAINULABADIN ZAFAR
ZAINULABADIN ZAFAR 2019-9-1
Respected Sir,
I did tried for that also. It gives error. Can you help and plot. Plot is given in pdf file.
Thanks

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by