Why do solutions get lost when I add a assumption where I know, that the correct solution would satisfie this assumption. Same if I add the inequality instead of the assume

Hey,

I have a question refering to the the symbolic solution of a system of equations and unequations.

If I run my system of equations with solve, it finds a solution which satisfies x+y>1, as expected. But if I add x+y>1 as an inequation in solve or I demand this by assume(x+y>1), than the solution either gets lost or I get a warning "not able to find an explicit solution". Why is that so?

Thank to everybody who trys to help.

[Torsten]

Maybe you didn't specify x and y as real-valued:

syms x y real

?

We cannot tell because we don't know your problem.

[Felix]

Hey, thank you for you response,

x and y are defined as real. It seems like it only accepts inequalitys like x>0 or y>0 but he struggles with x+y>1

I'm new here, thats why I dont knwo how everything works but I try to show you my code.

clear all
close all
clc
syms A11 A12 A13 A21 A22 A23 B11 B12 B13 B21 B22 B23 C11 C12 C13 C21 C22 C23 w z lambda my x y real
% assume(x>0)
% assume(y>0)
% assumeAlso(x+y>1)
% assume(w>0)
% assume(z>0)
% assume(lambda>0 & lambda<1)
% assume(my>0 & my<1)
% assume((w-1)*y+z*(1-x)<0)
% assume(w*(1-y)+(z-1)*x<0)
% assumptions
F = [A11 A12; A21 A22];
f = [A13; A23];
G = [B11 B12; B21 B22];
g = [B13; B23];
H = [C11 C12; C21 C22];
h = [C13; C23];
P1 = [0;0];
P2 = [1;0];
P3 = [x;y];
P4 = [0;1];
% A1 = [0;0]
% A2 = [lambda;0]
% A3 = [w;z]
% A4 = [0;my]
A = [0 lambda w 0;...
    0 0 z my];
% B1 = [1;0]
% B2 = [x;y]
% B3 = [w;z]
% B4 = [lambda;0]
B = [1 x w lambda;...
    0 y z 0];
% C1 = [0;1]
% C2 = [0;my]
% C3 = [w,z]
% C4 = [x;y]
C = [0 0 w x;...
    1 my z y];
D = [1 2 3 4; 2 3 4 1; 3 4 1 2; 4 1 2 3;...
    1 4 3 2; 2 1 4 3; 3 2 1 4; 4 3 2 1];
i = 4
A1b=A(:,D(i,1));
A2b=A(:,D(i,2));
A3b=A(:,D(i,3));
A4b=A(:,D(i,4));
for j = 8
    B1b=B(:,D(j,1));
    B2b=B(:,D(j,2));
    B3b=B(:,D(j,3));
    B4b=B(:,D(j,4));
    for k = 5
        C1b=C(:,D(k,1));
        C2b=C(:,D(k,2));
        C3b=C(:,D(k,3));
        C4b=C(:,D(k,4));
        p=D(i,:); 
        q=D(j,:);
        r=D(k,:);
        eqns = [cat(1,F*P1+f==A1b,  F*P2+f==A2b, F*P3+f==A3b, F*P4+f==A4b, G*P1+g==B1b, G*P2+g==B2b, G*P3+g==B3b, G*P4+g==B4b, H*P1+h==C1b, H*P2+h==C2b, H*P3+h==C3b, H*P4+h==C4b, lambda<1, lambda>0, my>0, my<1, x>0, y>0, w>0, z>0)]
        txt = ['Perm.F ist   ' num2str(p) '   Perm.G ist   ' num2str(q) '   Perm.H ist   ' num2str(r)]
        sol = solve(eqns, [A11 A12 A13 A21 A22 A23 B11 B12 B13 B21 B22 B23 C11 C12 C13 C21 C22 C23 lambda my w z x y],'ReturnConditions',true);
        A11Sol = sol.A11;
        A12Sol = sol.A12;
        A13Sol = sol.A13;
        A21Sol = sol.A21;
        A22Sol = sol.A22;
        A23Sol = sol.A23;
        B11Sol = sol.B11;
        B12Sol = sol.B12;
        B13Sol = sol.B13;
        B21Sol = sol.B21;
        B22Sol = sol.B22;
        B23Sol = sol.B23;
        C11Sol = sol.C11;
        C12Sol = sol.C12;
        C13Sol = sol.C13;
        C21Sol = sol.C21;
        C22Sol = sol.C22;
        C23Sol = sol.C23;
        lambdaSol = sol.lambda
        mySol = sol.my
        wSol = sol.w
        zSol = sol.z
        xSol = sol.x
        ySol = sol.y
        conditionsSol = sol.conditions
        % lambdaSol = eval(lambdaSol)
        % mySol=eval(mySol)
        % wSol=eval(wSol)
        % zSol=eval(zSol)
        xSol=eval(xSol)
        ySol=eval(ySol)
    end
end

There reason for the while construction is, that I want to test multiple combinations how points get mapped on each other. Thats also why I added conditions, because sometime the solution will be a family.

The solution I get are

lambdaSol =

23/(144*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) - ((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)/4 + 5/24

mySol =

23/(36*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) - ((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3) + 5/6

wSol =

(((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3) - 23/(36*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) + 1/6)^2/2 + 23/(144*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) - ((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)/4 + 17/24

zSol =

23/(72*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) - (((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3) - 23/(36*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) + 1/6)^2 - ((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)/2 + 5/12

xSol =

(((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3) - 23/(36*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) + 1/6)^2/2 - 23/(144*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) + ((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)/4 + 19/24

ySol =

((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3) - 23/(36*((29^(1/2)*108^(1/2))/108 + 19/216)^(1/3)) + 1/6

conditionsSol =

symtrue

xSol =

0.8478

ySol =

0.2581

So it is clear, that x+y>1 is true. However, when I now add x+y>1 at eqns than the output is:

Warning: Unable to find explicit solution. For options, see help.

> In sym/solve (line 317)

Answer 1

[Walter Roberson]

First of all, the Symbolic Engine is not very good about reasoning with inequalities.

Second of all, when the Symbolic Engine returns an indication of an solution inequality (range) and you do not have ReturnConditions set, then when you go to display, the MATLAB level examines the range and picks one representative value from the range to present, without indication that an inequality exists.

For example,

syms X
ineq = (X-3)^2+2<5

ineq = 

sol1 = solve(ineq)

sol1 = 

sol2 = solve(ineq, 'returnconditions', true)

sol2 = struct with fields:
             X: [2×1 sym]
    parameters: x
    conditions: [2×1 sym]

sol2.X

ans = 

sol2.parameters

ans = 

x

sol2.conditions

ans = 

In sol2, the first solution for X involves values of a parameter x such that particular conditions are met. In sol1, the MATLAB level's rules about picking one particular solution first checked 0, found that did not work, then checked π and found that did work so it arbitrarily picked π . In sol2, the second solution for X involves all real-valued x; the value 0 is first checked and found to be an acceptable real so 0 is substitued into the 3+x*i giving 3 as the solution for sol1 .

Both cases are families of solutions, but the MATLAB level has presented π and 3 as if they are the only solutions. It is rather misleading.

Thus you have two situations: the internal computation engine can only handle relatively simple inequalities; and the MATLAB level gives weak or misleading answers even when the computation engine is able to reason about the inequalities.

Often, in order to deal with inequalities, you have to convert the expression into a equality involving an extra variable that you constrain to be > 0 or >= 0 . For example,

syms delta_xy positive
eqn = x + y == 1 + delta_xy

with delta_xy constrained to be non-negative, x + y > 1 has to be true.

[Felix]

Hey Walter Roberson, thank you very much for your answer. Its a really nice trick to maka inequations to equations.

That is what it looks right now.

clear all
close all
clc


syms A11 A12 A13 A21 A22 A23 B11 B12 B13 B21 B22 B23 C11 C12 C13 C21 C22 C23 
syms x y lambda my w z delta_xy delta_lambda1 delta_my1 delta_xywz1 delta_xywz2 positive 

% assume(x>0)
% assume(y>0)
% assumeAlso(x+y>1)
% assume(w>0)
% assume(z>0)
% assume(lambda>0 & lambda<1)
% assume(my>0 & my<1)
% assume((w-1)*y+z*(1-x)<0)
% assume(w*(1-y)+(z-1)*x<0)
% assumptions


F = [A11 A12; A21 A22];
f = [A13; A23];
G = [B11 B12; B21 B22];
g = [B13; B23];
H = [C11 C12; C21 C22];
h = [C13; C23];

P1 = [0;0];
P2 = [1;0];
P3 = [x;y];
P4 = [0;1];

% A1 = [0;0]
% A2 = [lambda;0]
% A3 = [w;z]
% A4 = [0;my]

A = [0 lambda w 0;...
    0 0 z my];

% B1 = [1;0]
% B2 = [x;y]
% B3 = [w;z]
% B4 = [lambda;0]

B = [1 x w lambda;...
    0 y z 0];

% C1 = [0;1]
% C2 = [0;my]
% C3 = [w,z]
% C4 = [x;y]

C = [0 0 w x;...
    1 my z y];

D = [1 2 3 4; 2 3 4 1; 3 4 1 2; 4 1 2 3;...
    1 4 3 2; 2 1 4 3; 3 2 1 4; 4 3 2 1];


for i = 5
 A1b=A(:,D(i,1));
 A2b=A(:,D(i,2));
 A3b=A(:,D(i,3));
 A4b=A(:,D(i,4));
 AnzahlLoesungen=0
    for j = 1:8
        B1b=B(:,D(j,1));
        B2b=B(:,D(j,2));
        B3b=B(:,D(j,3));
        B4b=B(:,D(j,4));
        for k = 1:8
            C1b=C(:,D(k,1));
            C2b=C(:,D(k,2));
            C3b=C(:,D(k,3));
            C4b=C(:,D(k,4));


            p=D(i,:); 
            q=D(j,:);
            r=D(k,:);

            eqns = [cat(1,F*P1+f==A1b,  F*P2+f==A2b, F*P3+f==A3b, F*P4+f==A4b, G*P1+g==B1b, G*P2+g==B2b, G*P3+g==B3b, G*P4+g==B4b, H*P1+h==C1b, H*P2+h==C2b, H*P3+h==C3b, H*P4+h==C4b, x+y==1+delta_xy, lambda==1-delta_lambda1, my==1-delta_my1, (w-1)*y+z*(1-x)== - delta_xywz1, w*(1-y)+(z-1)*x== - delta_xywz2)];

            txt = ['Perm.F ist   ' num2str(p) '   Perm.G ist   ' num2str(q) '   Perm.H ist   ' num2str(r)]

            sol = solve(eqns,[A11 A12 A13 A21 A22 A23 B11 B12 B13 B21 B22 B23 C11 C12 C13 C21 C22 C23 lambda my w z delta_xy delta_lambda1 delta_my1 delta_xywz1 delta_xywz2 x y],'ReturnConditions',true);

What I still dont understand: If I say "syms A11 ... C23 real" , what they should be (per defintions of a affine map but also because I already know the solution), than the problem remains. Than the programm will tell me the same warning: cant find explicite solution. When I cancel real from syms than it works and tells me the true answer, that there is no solution. Why is that?

Furthermore I would like to know why the following happens: in solve(eqns,var,conditions) I tell the programm for whch variables it should solve. Why is there a diffrence if I say [A11 ... x y] or if I say [A11 ... y x]? In one case it can determine that there is no solution in a few seconds, in the other case the programm is working for houres and never stops. Is it possible to find out why the programm cant answer the question if the system of equations has a solution or not?

[Walter Roberson]

The order of variables passed to solve() determines the order that solve does eliminations.

When you are working with multiple variables and non-linear equations, you often get to a point in the calculations where the Symbolic Engine has to decide which variable it is going to solve for next.

The Symbolic Engine does not create parallel processes to attempt to solve for all of the possibilities simultaneously (to the limit of the number of available cores): it choses one possibility and tries to chase it down to a conclusion, not backtracking until it has proven that the entire tree from that point cannot work. So the order it looks at variables can be important.

syms A11 ... C23 real

There are algebraic rules and identies that exist over the complex numbers that might not always be true over reals. The Symbolic Engine sometimes just skips trying to apply those rules to systems involving reals, instead of doing something like checking to see whether the rule applies anyhow in the circumstances, or instead of injecting an assumption that would allow the rule to be true and then reasoning about the augmented system.

[Felix]

Thanks again @Walter Roberson.

But what can I do? If I say A11 ... C23 real than it doenst catch all solutions and gives me a warning that I maybe miss some solutions. When I let that open, than the programm is computing to death.

How long can it take until the programm realises, that chasing down one variable leeds to nothing? It worked over night while in all the other cases it needs a maximum of maybe 15 sec. The systems also dont differ much. It is just the allocation of mapped points. How can I help the programm to not get lost in the wrong tress? Does it help the programm if i kill some trivial equations? Some are really easy. when ever I put (0,0) into the affin map than the entries of the translations vektors are simply there. so I could reduce maybe 5 -10 equations.

[Walter Roberson]

If I say A11 ... C23 real than it doenst catch all solutions

Yes, that does indeed happen with the Symbolic Engine.

How long can it take until the programm realises, that chasing down one variable leeds to nothing

Let me put it this way: I ran one symbolic calculation for 13 days. When I put in a breakpoint and checked to see how far it had reached, I could see that it was only a small fraction of the way through.

As far as I know, the algorithms in the Symbolic Engine do not have any checks along the lines of "It's been 3 days already, this is too hard, I'm going to give up!"

[Walter Roberson]

I have sometimes solved for a few variables "manually", substituted the values in, and asked to solve() the result.

• Sometimes the result is something that the Symbolic Engine gives up on, or takes much longer on than if I had not done the pre-calculation -- that is, sometimes my steps make the situation worse.
• Sometimes the result is something that the Symbolic Engine is able to calculate with, when otherwise it would have given up. This is especially the case when there are fractional powers involved: the Symbolic Engine has a tendancy to give up on fractional powers earlier than necessary.
• Sometimes the result is just to be able to track down what, more precisely, is preventing the calculation from working. When solve() cannot find a solution, it never explains what is wrong; sometimes you need to guide it through steps to figure out the mathematical issue being encountered.

If you can eliminate variables appearing linearly then it might well be worth doing so. Beyond that, elimination tends to help more often than not, but elimination can cause problems as well.