Info
Questa domanda è chiusa. Riaprila per modificarla o per rispondere.
What is wrong with this code? Need help.
1 visualizzazione (ultimi 30 giorni)
Mostra commenti meno recenti
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
3 Commenti
KSSV
il 13 Dic 2018
YOur size of x is 5001*3.....but at that line you are trying to extract 50001 columns...you need to think on this line:
n=length(x(:,1));
Walter Roberson
il 13 Dic 2018
Your x is a 50001 x 3 array. m is 25000. n is 50001. What is your expectation as to what
y=diff(x(m,n,1))/dt;
will do for you? Even if x had that many columns, x(scalar,scalar,1) would be a scalar and diff() of a scalar is going to be empty.
Risposte (0)
Questa domanda è chiusa.
Vedere anche
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!