# Why am I getting an error 'index greater than 1' when trying to simulate a circuit based on ode15s

3 views (last 30 days)
Meme Young on 27 Mar 2021
Commented: Walter Roberson on 3 Apr 2021
I am trying to simulate a circuit in a simulink model but I get an error saying that index > 1; if I change to other solvers like ode23tb or ode23s, it will have sigularity warning and the results are NaN. I have checked my code but I dont know why. Can someone help me? files are as attached.
Meme Young on 28 Mar 2021
@Walter Roberson Thank you Mr Roberson, I have also found this answer, but I still dont know how to solve this problem and fix my code

Meme Young on 28 Mar 2021
@Walter Roberson I have checked my matrix. The answer you posted says that the (M + lamda * dF/dy) has to be non-singular for all nonzero lamda, so I tried to construct the dF/dy and M:
M = diag([ones(1,5), zeros(1,5)]);
A = [-R12 / L12 0 0 0 0 1 / L12 0 0 0 0;
0 -R23_1 / L23_1 0 0 0 0 1 / L23_1 0 0 0;
0 0 -R23_2 / L23_2 0 0 0 0 1 / L23_2 0 0;
0 0 0 0 0 0 0 0 1 / C23_2 0;
0 0 0 0 -R34 / L34 0 0 0 0 1 / L34;
0 0 0 1 0 0 -1 1 0 0;
0 0 0 0 0 1 1 0 0 1;
0 0 -1 0 0 0 0 0 1 0;
1 -1 -1 0 0 0 0 0 0 0;
1 0 0 0 -1 0 0 0 0 0]; % dF/dy
But I dont know how to tell if (M + lamda * A) is non-singular for all lamda ≠ 0.
Walter Roberson on 3 Apr 2021
baseMVA = 100; baseKV = 500;
Zbase = baseKV ^ 2 / baseMVA;
fn = 60; wn = 2 * pi * fn;
%% line data p.u.
r12 = 0.0014; x12 = 0.03;
r23_1 = 0.0067; l23_1 = 0.0739;
r23_2 = 0.0074; l23_2 = 0.08; c23_2 = 0.55 * l23_2;
r34 = 2e-4; x34 = 0.02;
%% line data nominal
R12 = Zbase * r12; L12 = Zbase * x12 / wn;
R23_1 = r23_1 * Zbase; L23_1 = l23_1 * Zbase / wn;
R23_2 = r23_2 * Zbase; L23_2 = l23_2 * Zbase / wn; C23_2 = 1 / (c23_2 * Zbase * wn);
R34 = r34 * Zbase; L34 = x34 * Zbase / wn;
M = diag([ones(1,5), zeros(1,5)]);
A = [-R12 / L12 0 0 0 0 1 / L12 0 0 0 0;
0 -R23_1 / L23_1 0 0 0 0 1 / L23_1 0 0 0;
0 0 -R23_2 / L23_2 0 0 0 0 1 / L23_2 0 0;
0 0 0 0 0 0 0 0 1 / C23_2 0;
0 0 0 0 -R34 / L34 0 0 0 0 1 / L34;
0 0 0 1 0 0 -1 1 0 0;
0 0 0 0 0 1 1 0 0 1;
0 0 -1 0 0 0 0 0 1 0;
1 -1 -1 0 0 0 0 0 0 0;
1 0 0 0 -1 0 0 0 0 0]; % dF/dy
syms lambda
D = det(M + lambda * A);
R = solve(D./lambda^7, lambda);
vpa(R)
ans =
subs(M, lambda, R(3))
ans =
That is obviously singular.
Therefore, there are three lambda values (one of them real) for which M + lambda*A is singular.
If you are correct that "the (M + lamda * dF/dy) has to be non-singular for all nonzero lamda" then you do not meet the conditions, and that would imply that you cannot reduce the DAE order.
Since the slx file can run, it must be a set of solvable DAE
Your code might possibly not match the slx model.

R2020b

### Community Treasure Hunt

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

Start Hunting!