Solving DAE - Index may be greater than one.

This is my program to solve DAEs:
function btp_25Jan4
syms cm(x) Pe(x) g(x) delpie(x) Rad(x) Re Sc de Rr C0 D A B kp t Peo delp Rm Rg Mw kad n L;
eqs = [(1/(2*(cm-1)))*(g-(g^2)/5)*diff(cm,1)-(1-3*g/10)*diff(g,1)==(8/(Re*Sc*de/L))*(8-Pe*g)/(g^2),...
Pe*cm*Rr*g + 8*(cm-1)==(de/(C0*D))*g*(A*cm/(1+B*cm))*kp*(exp(-kp*t)),...
Pe==Peo*(1-delpie/delp)/( 1+Rad/Rm),...
delpie==2*Rg*C0*Rr*cm/Mw,...
Rad==kad*((A*cm/(1+B*cm))*(1-exp(-kp*t)))^n];
vars = [cm(x); Pe(x); g(x); delpie(x); Rad(x)];
Mu = 0.001;
Re = 600;
D = 1.4*(10^-9);
Sc = Mu/D;
de = (4*0.08/pi)^0.5;
Rr = 0.72;
C0 = 4*(10^-3);
A = 0.26;
B = 130;
kp = 0.58;
t = 2*3600;%2 hours
Peo = 24;
delp = 400000;
Rm = 1/(0.001*1.6*10^-11);
Rg = 8.314;
Mw = 20000;
kad = 4*(10^13);
n = 0.15;
L = 0.145;
go = (128/(3*Re*Sc*de/L))^(1/3);
[eqs, vars, ~] = reduceDifferentialOrder(eqs, vars);
[DAEs,DAEvars] = reduceDAEIndex(eqs,vars);
[DAEs,DAEvars] = reduceRedundancies(DAEs,DAEvars)
DAEs = subs(DAEs);
f = daeFunction(DAEs, DAEvars);
F = @(x, Y, YP) f(x, Y, YP);
y0est = [0; 1];% have used only 2 because after reduceRedundancies only 2 DAEvars and DAEs are left
yp0est = zeros(2,1);
opt = odeset('RelTol', 10.0^(-2), 'AbsTol' , 10.0^(-2));
[y0, yp0] = decic(F, 0.0005, y0est,[], yp0est,[], opt);
end
But I get this error
Error using decic>sls (line 170)
Index may be greater than one.
Error in decic (line 77)
[dy,dyp] = sls(res,dfdy,dfdyp,neq,free_y,free_yp);
Error in btp_25Jan4 (line 48)
[y0, yp0] = decic(F, 0.0005, y0est,[], yp0est,[], opt);
Please help me to know where have I done the mistake. Thanks.

3 commentaires

akshay
akshay le 12 Mar 2016
hey were you able to resolve the error
I am getting same error...don't know how to proceed
Ghaith Salah
Ghaith Salah le 6 Mai 2018
Im having the same problem, kindly let me know if you resolved the error
Ramanuja Srinivasan
Ramanuja Srinivasan le 31 Août 2020
Hi, have you sorted out this problem, because i am encountering the same error. I would really appreciate your help.

Connectez-vous pour commenter.

Réponses (1)

Szar
Szar le 13 Oct 2017

1 vote

Try fixing some of the initial guesses.

Catégories

En savoir plus sur Mathematics dans Centre d'aide et File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by