Error en índices sistema EDOs
Publicado por Geroge (1 intervención) el 14/07/2014 01:44:40
Hola, pretendo resolver un sistema de EDO´s pero no he podido. Les dejo mi código:
El error que arroja matlab es:
Attempted to access y(3); index out of bounds because numel(y)=1.
Error in concentraciones (line 40)
v2=k2*(y(r1+1)*y(r1)-(Keq2*y(r1-1)*(0.055-2*y(r1+1))^2)/y(r1));
Error in odearguments (line 88)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in Solve_concentraciones (line 30)
[t,y] = ode45(@concentraciones,[0 50],y0);
No entiendo el error. Alguna sugerencia por favor?
Muchas gracias
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
function dy = concentraciones(t,y)
% Constantes
VI=2*10^(-5);
VII=11*10^(-5);
VIII=8*10^(-4);
Def=1.7*10^(-10);
S1=2.83*10^(-11);
N2=11944.5;
S2=2.37*10^(-5);
RfI=1.5*10^(-6);
RfII=1.375*10^(-3);
RfIII=RfI+RfII;
Keq1=106;k1=2.5*10^(-5);
Keq2=10;k2=1.1*10^(-6);
kA=3.0*10^(-5);
V=(VI+VII)/VII;
% Dominio de resolución
dr=0.5*RfI;
dt=10^(-3);
RI=dr:dr:RfI;
RII=RfI+dr:dr:RfII;
RIII=RfII+dr:dr:RfIII;
R=[RI RII RIII];
r1=length(RI);
r2=length(RI)+length(RII);
r3=length(R);
dy=ones(length(R),1);%matriz de concentraciones
for r=1:1:r3
if r<r1v2=k2*(y(r1+1)*y(r1)-(Keq2*y(r1-1)*(0.055-2*y(r1+1))^2)/y(r1));
dy(r)=(S1/VI)*v2;
else
if r==r1
y(r)=2*7.65*10^(-4)+2.5-2*y(r-1);
else
if r==r1+1
v2=k2*(y(r)*y(r1)-(Keq2*y(r1-1)*(0.055-2*y(r))^2)/y(r1));
dy(r)=(V*Def/dr^2)*(2*y(r+1)-2*y(r))-(S1/VII)*v2;
else
if (r1+1<r)&&(r<r2)
v2=k2*(y(r1+1)*y(r1)-(Keq2*y(r1-1)*(0.055-2*y(r1+1))^2)/y(r1));
dy(r)=(V*Def/dr^2)*(y(r+1)-2*y(r)+y(r-1))-(S1/VII)*v2;
else
if r==r2
v2=k2*(y(r1+1)*y(r1)-(Keq2*y(r1-1)*(0.055-2*y(r1+1))^2)/y(r1));
dy(r)=(2*V*Def/dr^2)*((dr*kA/Def)*(y(r2+2)-y(r2+1))-y(r)+y(r-1))-(S1/VII)*v2;
else
if r==r2+1
v1=k1*((y(r)*(0.055-2*y(r2))/(10^(-3)))-((y(r2)*10^(-3))/(Keq1*(0.055-2*y(r2)))));
y(r)=y(r+1)+((N2*S2)/(VIII*kA))*v1;
else
if r2+1<r
dy(r)=-(N2*S2/VIII)*kA*(y(r)-y(r2+1));
end
end
end
end
end
end
end
end
end
%Para resolver utilizo:
RfI=1.5*10^(-6);
RfII=1.375*10^(-3);
RfIII=RfI+RfII;
dr=0.5*RfI;
RI=dr:dr:RfI;
RII=RfI+dr:dr:RfII;
RIII=RfII+dr:dr:RfIII;
R=[RI RII RIII];
y0=zeros(1,length(R));
for i=1:length(R)
if i==1
y0=0;else
if i==2
y0=2.5;
else
if i==length(R)-1
y0=0;else
if i==length(R)
y0=7.65*10^(-4);
end
end
end
end
end
[t,y] = ode45(@concentraciones,[0 50],y0);
El error que arroja matlab es:
Attempted to access y(3); index out of bounds because numel(y)=1.
Error in concentraciones (line 40)
v2=k2*(y(r1+1)*y(r1)-(Keq2*y(r1-1)*(0.055-2*y(r1+1))^2)/y(r1));
Error in odearguments (line 88)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 114)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in Solve_concentraciones (line 30)
[t,y] = ode45(@concentraciones,[0 50],y0);
No entiendo el error. Alguna sugerencia por favor?
Muchas gracias
Valora esta pregunta
0