Error message: Invalid index.

Hi everyone, I am trying to solve two second order differential equations with two unknowns. However, matlab gives the warning that:

" Warning: Number of equations greater than number of indeterminates. Trying heuristics to reduce to square system. 
> In symengine
  In mupadengine/evalin (line 123)
  In mupadengine/feval (line 182)
  In dsolve>mupadDsolve (line 340)
  In dsolve (line 194)
  In InvalidIndexError (line 34)"

Then it later gives an error that:

" Error using mupadengine/feval (line 187)
Invalid index.

Error in dsolve>mupadDsolve (line 340)
T = feval(symengine,'symobj::dsolve',sys,x,options);

Error in dsolve (line 194)
sol = mupadDsolve(args, options);

Error in InvalidIndexError (line 34)
[aSol, bSol] = dsolve(odes,bc)".

What can I do? Please see the code below and on the attachment:

syms A B C E F G H I J K N O P a(x) b(x)

%%Parameters
A = 3.22e-9;   
B = 2.13e-9;   
C = 1.28e-9;   
E = 1.66e-8;   
F = 8.06e-9;  
G = 8.14e-5;   
H = 149;       
I = A/G;       
J = 6.24;      
K = 5.68e-8;   
N = 84.86;     
O = 0.63;      
P = 5.3e-8;    

%%Equations

ode1 = A * (1/J)* diff((a * b),x,2) + B * diff(a,x,2) + C * O * J * K * diff(((1)/(b^2 + J * b + J * K)),x,2)  == 0;
ode2 = E * diff(b,x,2) - B * diff(a,x,2) - 2 * C * O * J * K * diff(((1)/(b^2 + J * b + J * K)),x,2)- F * P* diff((1/b),x,2) ==0;

odes = [ode1; ode2];

%%Boundary conditions
bc1 = a(0) == (J * H * N)/3.47;
bc2 = E * diff(b,x) - B * diff(a,x) - 2 * C * O * J * K * diff(((1)/(b^2 + J * b + J * K)),x)- F * P* diff((1/b),x) ==0;
bc3 = a(I) == ((O * J * b)/(b.^2 + J * b + J * K));
bc4 = b(I)== 2.24e-4;
bc = [bc1; bc2; bc3; bc4];

[aSol, bSol] = dsolve(odes,bc)

fplot(aSol)
hold on
fplot(bSol)

legend('HSO_3^-','H^+','Location','best')

Answer 1

[Walter Roberson]

The invalid index is probably a MuPAD bug.

However, you also have a problem with your system.

Your bc1 and bc4 are good boundary conditions.

Your bc2 boundary condition is a differential equation, not a boundary condition. You have included it in the bc rather than in the odes, which is a logic error on your part. However, it does not matter, because MATLAB merges the odes and boundary conditions before passing on to MuPAD, so this just hints at an error in your thinking about the system. But this puts the number of ode equations up to 3 in two functions, which is part of the reason you get the message about the mismatch in count.

Your bc3 is a(I) == ((O * J * b)/(b.^2 + J * b + J * K)) . This locks in a() at one particular location I (which is A/G, approximately 3.9E-5) on the left side. But the right side involves b that is not evaluated at any particular location, so b(x) . That forces a(I) to be a function of x rather than function evaluated at a particular location. MATLAB is not smart enough to figure that out. If you want b to be evaluated at the same location then you need to specify that:

a(I) == ((O * J * b(I))/(b(I).^2 + J * b(I) + J * K))

To make it easier to compute, I suggest

bI = sym(2.24e-4);
bc3 = a(I) == ((O * J * bI)/(bI.^2 + J * bI + J * K));
bc4 = b(I) == bI;

which reminds me that I recommend putting sym() around all of the numeric values in your definitions of A through P: otherwise because of round-off error in calculating doubles, you would end up being slightly inconsistent in your definitions of the functions.

After this fix for a(I), because bc2 is still a differential equation, you would still have a mismatch between the number of equations and the number of functions, but at least you will get rid of the index error: instead you will quickly get the accurate

 Warning: Unable to find explicit solution.

Your system is somewhat complex, and I would not expect there to be any useful closed-form solution to it, even without the bc2 challenge.

[Dursman Mchabe]

Thanks a lot Walter, I now have a better understanding. I will try all your recommendations.

[Dursman Mchabe]

Wow! Your suggestion worked wonders. I don't get any error, however, the plot comes out blank. What can I try?

[Walter Roberson]

The solution is probably coming out empty, because MATLAB cannot resolve the system to a closed form because it is too complex.

Did you end up changing bc2 to not be a differential equation? That would help . But I think even if you did that the system would be too complex to resolve symbolically.

[Dursman Mchabe]

I set:

bc2 = b(0)==sym(3.47);

But still, the graph, comes out blank. I don't know what else to try for now.