Why does replacing zeros in the mass matrix by small numbers solve: "DAE appears to be of index greater than 1"?

1 次查看（过去 30 天）

Danny Helwegen 2022-4-16

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1698115-why-does-replacing-zeros-in-the-mass-matrix-by-small-numbers-solve-dae-appears-to-be-of-index-gre

评论： Torsten 2022-4-19

Hello,

I've a set of stiff differential equations that I wanted to solve using ode15s. At small initial concentrations the system worked fine, no problems. However, when I increased the initial values by a factor 10 the error following error popped up:

This DAE appears to be of index greater than 1.

In order to speed up the simulations I used a Mass Matrix M with zeros on the diagonal line except for the last cell, a simple example:

As stated previously, this worked fine for small inital values but not for larger initial values. From some literature and other questions I saw that by replacing the zeros on the Mass Matrix diagonal line by small values, e.g. 1e-12, the system would work. Which it did, so this is the mass matrix that I currently use:

1e-12 0     0     0     0 0
   1e-12 0     0     0 0
   0     1e-12 0     0 0
   0     0     1e-12 0 0
   0     0     0     1 0 
   0     0     0     0 1

So the problem is solved, only I don't know why. What's the reason that the DAE index is greater than 1 when using zeros, but when using really small values it is fine?

Does somebody have an explanation for this?

As an example, the simplified version of my system (where parameter values are stored in .par):

options = odeset('Mass',M,'RelTol',1e-6,'AbsTol',1e-6);
[t,c]=ode15s(@Equations,par.tspan,par.Initial,options,par); 
function dcdt = Equations(c,par)
% Create a storage
dcdt = zeros(size(c));
% Determine Rate
Rate = par.k .* c(1,1:end-1); 
% Determine Reaction Profile
    dcdt(1,1:end-1) = - Rate;
    dcdt(2,1:end-1) = + Rate;
% Determine Liquid Profile
dcdt(1,end) = dcdt(1,par.Nr+1) + par.A * (par.c0 - c(1,end) ) + par.B * ( par.C - c(1,end) );
dcdt(2,end) = dcdt(2,par.Nr+1) + par.A * (par.c0 - c(2,end) ) + par.B * ( par.C - c(2,end) );
end

Kind regards,

Danny

15 个评论
显示 13更早的评论隐藏 13更早的评论

Danny Helwegen 2022-4-19

在 MATLAB Online 中打开

Even if I say that no conditions are fixed, fixed_y0 = fixed_yp0 = [ ], I recieve the following error:

Error using decic>sls (line 170)
Index may be greater than one.

Torsten 2022-4-19

Fix the differential variables to their initial values and give initial guesses for the rest (algebraic variables and derivatives). This should fix your system.

If decic is not able to compute consistent initial values, I could look over your equations if you post them in a comprehensible way.

请先登录，再进行评论。

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