MATLAB: Using solve with inequalities

linprogMATLABsolve

I would like to use Solve to check if there are solutions to various systems of equalities/inequalities. Here's a toy example:
syms x y
S = solve(x+y==1,x<0.4,y<0.4,[x,y]);
Obviously there's no solution to that one, so Matlab gives me empty S.x and S.y. Interestingly, when I do the following, even though there are obviously a continuum of solutions, it spits out one answer:
S = solve(x+y==1,x<0.6,y<0.6,[x,y]);
It gives S.x = S.y = 1/2, which is, in some sense, the most intuitive possible solution. I like that it gives me that. Unfortunately, when I do the following, it's more concerned about that fact that there is a continuum of solutions:
S = solve(x+y==1,x<0.6,y<0.7,[x,y]);
In this case, it gives me:
Warning: Cannot find explicit solution.
What I would like to do is to have it spit out just one possible solution — I don't really care which it chooses. I know I can get a parametrized solution doing the following:
S = solve(x+y==1,x<0.6,y<0.7,[x,y],'ReturnConditions', true)
But then the mathworks page for Solve suggests how I can then propose a solution and have Matlab check whether it meets those conditions. But it doesn't offer advice on how to have Matlab propose a solution. Is there a way to do so?
Incidentally, I've also considered doing this numerically with linprog, which seems a little faster. The issue I have is with the tolerance on the constraints. Since I have strict inequalities, I've been making them weak inequalities by saying a <= b+1e-4 is roughly equivalent to a<b. However, this causes some issues for me (where it fails to find solutions where some exist). I'd really need to be able to set the constraint tolerances to zero to do it that way I think, and I don't believe that's possible.
Any advice (on how to do this numerically using linprog or symbolically) would be greatly appreciated.

Best Answer

My suggestion would be to use LCON2VERT ( Download ) if your system of in/equalities isn't too large (it couldn't be, I guess, if you were hoping to do the analysis symbolically) and if they define a bounded region.
LCON2VERT will give you the vertices of the polyhedral region corresponding to the in/equalities. If you just want an arbitrary point that strictly satisfies the inequalities, just take the mean of all the vertices.
As an example, instead of
S = solve(x+y==1,x<0.6,y<0.7,[x,y]);
you can do this,
>> Aeq=[1 1]; beq=1; A=eye(2); b=[0.6;0.7];
>> vertices=lcon2vert(A,b,Aeq,beq)
vertices =
0.6000 0.4000
0.3000 0.7000
>> mean(vertices)
ans =
0.4500 0.5500