Clear Filters
Clear Filters

problem with solve y sym

1 view (last 30 days)
javiera de la carrera
javiera de la carrera on 20 Oct 2017
Commented: Walter Roberson on 23 Oct 2017
Hi, I am trying to solve a system of 11 equations with symbolic variables:
S=[q1 - q2.*(0.11.*x1aux) - q4.*(sumxxaux+0.17.*x2aux) - fc1aux + vpasaj1aux - lf1 - v1.*Delayxaux - q2*0.11.*x1aux - x1aux.*v1.*exp(0.3487*sumq0aux/6).*0.3487./6==0;
q6 - q7.*(0.17.*x2aux) - q4.*(sumxxaux+0.11.*x1aux) - fc2aux + vpasaj2aux - lf2 - v2.*Delayxaux - q7*0.17.*x2aux - x2aux.*v2.*exp(0.3487*sumq0aux/6).*0.3487./6==0;
q4/q2-dif2==0];
P= [q51 - q31.*xaux(1,1) - q4.*(0.11.*x1aux+0.17.*x2aux+sumxxaux-xaux(1,1)) - fcaux(1,1) + vpasajaux(1,1) - lfaux(1,1) - vaux(1,1).*Delayxaux - q31.*xaux(1,1) - xaux(1,1).*vaux(1,1).*exp(0.3487*sumq0aux/6).*0.3487./6==0;
q52 - q32.*xaux(1,2) - q4.*(0.11.*x1aux+0.17.*x2aux+sumxxaux-xaux(1,2)) - fcaux(1,2) + vpasajaux(1,2) - lfaux(1,2) - vaux(1,2).*Delayxaux - q32.*xaux(1,2) - xaux(1,2).*vaux(1,2).*exp(0.3487*sumq0aux/6).*0.3487./6==0;
q53 - q33.*xaux(1,3) - q4.*(0.11.*x1aux+0.17.*x2aux+sumxxaux-xaux(1,3)) - fcaux(1,3) + vpasajaux(1,3) - lfaux(1,3) - vaux(1,3).*Delayxaux - q33.*xaux(1,3) - xaux(1,3).*vaux(1,3).*exp(0.3487*sumq0aux/6).*0.3487./6==0];
R= [(-(3*q4^4 - 2*q4^3*q31 - 2*q4^3*q32 - 2*q4^3*q33 - 2*q4^3*q7 + q4^2*q7*q31 + q4^2*q7*q32 + q4^2*q7*q33 + q4^2*q31*q32 + q4^2*q31*q33 + q4^2*q32*q33 - q7*q31*q32*q33)/(3*q2*q4^4 + 3*q4^4*q7 + 3*q4^4*q31 + 3*q4^4*q32 + 3*q4^4*q33 - 4*q4^5 - 2*q2*q4^3*q7 - 2*q2*q4^3*q31 - 2*q2*q4^3*q32 - 2*q2*q4^3*q33 - 2*q4^3*q7*q31 - 2*q4^3*q7*q32 - 2*q4^3*q7*q33 - 2*q4^3*q31*q32 - 2*q4^3*q31*q33 - 2*q4^3*q32*q33 + q2*q4^2*q7*q31 + q2*q4^2*q7*q32 + q2*q4^2*q7*q33 + q2*q4^2*q31*q32 + q2*q4^2*q31*q33 + q2*q4^2*q32*q33 + q4^2*q7*q31*q32 + q4^2*q7*q31*q33 + q4^2*q7*q32*q33 + q4^2*q31*q32*q33 - q2*q7*q31*q32*q33))-(eta.*(0.11.*x1aux))./(q1 - q2.*(0.11.*x1aux) - q4.*(sumxxaux+0.17.*x2aux))==0;
(-(3*q4^4 - 2*q4^3*q31 - 2*q4^3*q32 - 2*q4^3*q33 - 2*q2*q4^3 + q2*q4^2*q31 + q2*q4^2*q32 + q2*q4^2*q33 + q4^2*q31*q32 + q4^2*q31*q33 + q4^2*q32*q33 - q2*q31*q32*q33)/(3*q2*q4^4 + 3*q4^4*q7 + 3*q4^4*q31 + 3*q4^4*q32 + 3*q4^4*q33 - 4*q4^5 - 2*q2*q4^3*q7 - 2*q2*q4^3*q31 - 2*q2*q4^3*q32 - 2*q2*q4^3*q33 - 2*q4^3*q7*q31 - 2*q4^3*q7*q32 - 2*q4^3*q7*q33 - 2*q4^3*q31*q32 - 2*q4^3*q31*q33 - 2*q4^3*q32*q33 + q2*q4^2*q7*q31 + q2*q4^2*q7*q32 + q2*q4^2*q7*q33 + q2*q4^2*q31*q32 + q2*q4^2*q31*q33 + q2*q4^2*q32*q33 + q4^2*q7*q31*q32 + q4^2*q7*q31*q33 + q4^2*q7*q32*q33 + q4^2*q31*q32*q33 - q2*q7*q31*q32*q33))-(eta.*(0.17.*x2aux))./(q6 - q7.*(0.17.*x2aux) - q4.*(0.11.*x1aux+sumxxaux))==0];
G= [(-(3*q4^4 - 2*q4^3*q7 - 2*q4^3*q32 - 2*q4^3*q33 - 2*q2*q4^3 + q2*q4^2*q7 + q2*q4^2*q32 + q2*q4^2*q33 + q4^2*q7*q32 + q4^2*q7*q33 + q4^2*q32*q33 - q2*q7*q32*q33)/(3*q2*q4^4 + 3*q4^4*q7 + 3*q4^4*q31 + 3*q4^4*q32 + 3*q4^4*q33 - 4*q4^5 - 2*q2*q4^3*q7 - 2*q2*q4^3*q31 - 2*q2*q4^3*q32 - 2*q2*q4^3*q33 - 2*q4^3*q7*q31 - 2*q4^3*q7*q32 - 2*q4^3*q7*q33 - 2*q4^3*q31*q32 - 2*q4^3*q31*q33 - 2*q4^3*q32*q33 + q2*q4^2*q7*q31 + q2*q4^2*q7*q32 + q2*q4^2*q7*q33 + q2*q4^2*q31*q32 + q2*q4^2*q31*q33 + q2*q4^2*q32*q33 + q4^2*q7*q31*q32 + q4^2*q7*q31*q33 + q4^2*q7*q32*q33 + q4^2*q31*q32*q33 - q2*q7*q31*q32*q33))-(eta.*(xaux(1,1)))./(q51 - q31.*xaux(1,1) - q4.*(0.11.*x1aux+0.17.*x2aux+sumxxaux-xaux(1,1)))==0;
(-(3*q4^4 - 2*q4^3*q7 - 2*q4^3*q31 - 2*q4^3*q33 - 2*q2*q4^3 + q2*q4^2*q7 + q2*q4^2*q31 + q2*q4^2*q33 + q4^2*q7*q31 + q4^2*q7*q33 + q4^2*q31*q33 - q2*q7*q31*q33)/(3*q2*q4^4 + 3*q4^4*q7 + 3*q4^4*q31 + 3*q4^4*q32 + 3*q4^4*q33 - 4*q4^5 - 2*q2*q4^3*q7 - 2*q2*q4^3*q31 - 2*q2*q4^3*q32 - 2*q2*q4^3*q33 - 2*q4^3*q7*q31 - 2*q4^3*q7*q32 - 2*q4^3*q7*q33 - 2*q4^3*q31*q32 - 2*q4^3*q31*q33 - 2*q4^3*q32*q33 + q2*q4^2*q7*q31 + q2*q4^2*q7*q32 + q2*q4^2*q7*q33 + q2*q4^2*q31*q32 + q2*q4^2*q31*q33 + q2*q4^2*q32*q33 + q4^2*q7*q31*q32 + q4^2*q7*q31*q33 + q4^2*q7*q32*q33 + q4^2*q31*q32*q33 - q2*q7*q31*q32*q33))-(eta.*(xaux(1,2)))./(q52 - q32.*xaux(1,2) - q4.*(0.11.*x1aux+0.17.*x2aux+sumxxaux-xaux(1,2)))==0;
(-(3*q4^4 - 2*q4^3*q7 - 2*q4^3*q31 - 2*q4^3*q32 - 2*q2*q4^3 + q2*q4^2*q7 + q2*q4^2*q31 + q2*q4^2*q32 + q4^2*q7*q31 + q4^2*q7*q32 + q4^2*q31*q32 - q2*q7*q31*q32)/(3*q2*q4^4 + 3*q4^4*q7 + 3*q4^4*q31 + 3*q4^4*q32 + 3*q4^4*q33 - 4*q4^5 - 2*q2*q4^3*q7 - 2*q2*q4^3*q31 - 2*q2*q4^3*q32 - 2*q2*q4^3*q33 - 2*q4^3*q7*q31 - 2*q4^3*q7*q32 - 2*q4^3*q7*q33 - 2*q4^3*q31*q32 - 2*q4^3*q31*q33 - 2*q4^3*q32*q33 + q2*q4^2*q7*q31 + q2*q4^2*q7*q32 + q2*q4^2*q7*q33 + q2*q4^2*q31*q32 + q2*q4^2*q31*q33 + q2*q4^2*q32*q33 + q4^2*q7*q31*q32 + q4^2*q7*q31*q33 + q4^2*q7*q32*q33 + q4^2*q31*q32*q33 - q2*q7*q31*q32*q33))-(eta.*(xaux(1,3)))./(q53 - q33.*xaux(1,3) - q4.*(0.11.*x1aux+0.17.*x2aux+sumxxaux-xaux(1,3)))==0];
D=[S;R;P;G];
[t11,t22,t331,t332,t333,t44,t551,t552,t553,t66,t77,param,cond]=solve(D,q1,q2,q31,q32,q33,q4,q51,q52,q53,q6,q7,'ReturnConditions',true,'Real',true);
But I obtain this answer:
Error using mupadengine/feval (line 163)
The row index is out of range.
Error in solve (line 369)
varargout{i} = transpose(eng.feval('map', solutions, '_index', i));
I don't know why! please help me!
Thanks
  6 Comments
Walter Roberson
Walter Roberson on 20 Oct 2017
Edited: Walter Roberson on 20 Oct 2017
It looks to me as if the code could be tested if it were preceded with
[Delayxaux, eta, fc1aux, fc2aux, lf1, lf2, sumxxaux, sumq0aux, v1, v2, vpasaj1aux, vpasaj2aux, x1aux, x2aux] = deal(rand, rand, rand, rand, rand, rand, rand, rand, rand, rand, rand, rand, rand, rand);
fcaux = rand(1,3);
lfaux = rand(1,3);
vaux = rand(1,3);
vpasajaux = rand(1,3);
xaux = rand(1,3);
syms dif2 q1 q2 q4 q6 q7 q31 q32 q33 q51 q52 q53
You appear to be solving 11 equations for 11 out of the 12 variables, so you would expect most of the outputs to be in terms of dif2
Can any constraints be put on the symbols? Such as real valued, or non-negative?
javiera de la carrera
javiera de la carrera on 20 Oct 2017
Yes, this is the code: syms q1 positive
syms q2 positive
syms q31 positive
syms q32 positive
syms q33 positive
syms q4 positive
syms q51 positive
syms q52 positive
syms q53 positive
syms q6 positive %a united
syms q7 positive %b united
syms dif2

Sign in to comment.

Answers (1)

Walter Roberson
Walter Roberson on 20 Oct 2017
With the above initialization, in R2017b, it works hard on it for a while and then says "Warning: Unable to find explicit solution." and returns empty symbols.
This does not surprise me. If you
collect(D, q4)
then you will see that you have q4^5 in some of the equations. There is no general closed form solution to equations of degree 5 and it is not possible for MATLAB to look for numeric approximations because of the free variable dif2 .
  5 Comments
javiera de la carrera
javiera de la carrera on 23 Oct 2017
Thanks! I think I have to change the problem. Or, Do you think another solution? Like a command to approximate equations or simplify to the solve comand? I tried with simplify D but nothing change.
Thanks!
Walter Roberson
Walter Roberson on 23 Oct 2017
Your R and your G both contain the ratio of quartics (polynomials of degree 4.) Just knowing that, we can predict that the solutions are going to be at least degree 7, or degree 8, or degree 16. You have 5 like that, and MATLAB needs to explore all of the possibilities for all of them. That is going to be slow and memory intensive.
Perhaps you can taylor with respect to q4 -- but to what order?
I extracted the division subexpression from the first R entry, and tossed in some random values for some of the variables
t2 = (3*q4^4 - 2*q4^3*q31 - 2*q4^3*q32 - 2*q4^3*q33 - 2*q4^3*q7 + q4^2*q7*q31 + q4^2*q7*q32 + q4^2*q7*q33 + q4^2*q31*q32 + q4^2*q31*q33 + q4^2*q32*q33 - q7*q31*q32*q33)/(3*q2*q4^4 + 3*q4^4*q7 + 3*q4^4*q31 + 3*q4^4*q32 + 3*q4^4*q33 - 4*q4^5 - 2*q2*q4^3*q7 - 2*q2*q4^3*q31 - 2*q2*q4^3*q32 - 2*q2*q4^3*q33 - 2*q4^3*q7*q31 - 2*q4^3*q7*q32 - 2*q4^3*q7*q33 - 2*q4^3*q31*q32 - 2*q4^3*q31*q33 - 2*q4^3*q32*q33 + q2*q4^2*q7*q31 + q2*q4^2*q7*q32 + q2*q4^2*q7*q33 + q2*q4^2*q31*q32 + q2*q4^2*q31*q33 + q2*q4^2*q32*q33 + q4^2*q7*q31*q32 + q4^2*q7*q31*q33 + q4^2*q7*q32*q33 + q4^2*q31*q32*q33 - q2*q7*q31*q32*q33)
t3 = subs(t2,[q2,q32,q33,q7],[2,32,33,7])
fsurf(t3,[-100 100 -80 40])
The result has at least 4 curved lines of discontinuity, indicating that if you cannot promise that the variables will stay in the "safe" ranges, that there cannot be any finite polynomial approximation.
I don't know... perhaps you could do something with continued fractions.

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!