Finding linear combination of two vectors such as every element is positive
13 次查看(过去 30 天)
显示 更早的评论
Lets say I have two vectors:

How to determinate existence and eventually how to find
and
such as



where every element of
satisfies:


0 个评论
采纳的回答
John D'Errico
2023-3-9
编辑:John D'Errico
2023-3-9
Hmm. Is this homework? I wonder, but I have learned not to trust that students won't post homework, and this is a cute problem.
In general, this is not always possible to do. For the trivial counter-example:
V1 = [1 -1]; V2 = [0 0];
Clearly there is no solution such that a*V1 + b*V2 > 0, though there is a trivial solution such that the sum is identically [0 0]. That happens at a=b=0. So while it may be possible in some circumstances, but as I showed, it is trivially easy to choose a pair of vectors such that no solution exists where the result is always strictly positive.
Can you find some solution, when one exists? It would seem this is then just a linear programming problem, although if the goal is to find a strictly positive solution, that will not fall under linear programming, since it will only allow you to specify a constraint that is greater than or equal to zero. Equality will be allowed. And then you have the issue of finding a NON-trivial solution, since a==b==0 will be the solution a solver will always find.
And, whenever a solution exists, it will almost always be true there will be infinitely many solutions, and there will be no trivial way to describe them all.
Having said all that, I would first point out that if ANY solution exists, such that a*V1 + b*V2 > 0 for some pair [a,b], then [k*a,k*b] is also a solution for any positive value of k. But that also allows you to resolve the issue of the non-solution where the sum could be zero. How? Suppose you found a solution (a,b) where every element of the sum was at least epsilon, for some specific arbitrarily small value of epsilon? Then the solution (a/epsilon,b/epsilon) is at least 1 for every element of the sum. And this gives you way to find A solution using linear programming. Just find the solution that has the sum at least 1 for EVERY element. Now we can use the greater than or equal to inequality that linprog assumes.
As example, I'll pick two vectors V1 and V2. In this example, to insure a solution will always exist, I chose a vector V2 which was ALWAYS positive. This will insure a solution exists.
V1 = randn(8,1);
V2 = rand(8,1);
A = [V1,V2]
You now need to find a vector of length 2, I'll call it X such that A*X>1 for every element. I'll use linear programming for the problem. Next, I want to solve the problem such that the elements of X are as small as possible. I'll do that using slack variables. Essentially, I want to solve the problem that abs(a)+abs(b) is as small as possible. This will give me a unique solution, the solution with minimum 1-norm. There will now be 4 unknowns to solve for, because of the slack variables I had to introduce. I'll write it as:
a1*V1 - a2*V1 + b1*V2 - b2*V2 > 1
where each of a1,a2,b1,b2 are all positive numbers, AND such that the sum a1+a2+b1+b2 is as small as possible. I could formulate this using an optimproblem, or I could just set up the problem directly into linprog.
n = length(V1);
a = optimvar('a',[1 2]);
a.LowerBound = [0 0];
b = optimvar('b',[1 2]);
b.LowerBound = [0 0];
prob = optimproblem;
prob.Objective = sum(a) + sum(b);
prob.Constraints.C1 = V1*a(1) - V1*a(2) + V2*b(1) - V2*b(2) >= ones(n,1)
absol = solve(prob)
In this case, the solution was found. Again, the use of slack variables will ALWAYS have one of the variables a(1) or a(2) equal to 0, similarly for b(1) and b(2). One member of each pair will always be zero. This is a classical way to solve for a 1-norm solution. So the solution is
asol = dot(absol.a,[1 -1])
bsol = dot(absol.b,[1 -1])
Does this solution work?
V3 = asol*V1 + bsol*V2
Of course. So all STRICTLY positive elements. Note that the smallest such element is exactly 1 in V3. Again, you can scale that solution by any positive constant k.
How can you know if no solution exists? No solurion exists if linprog returns with no feasible solution found. For example...
V1 = randn(8,1);
V2 = randn(8,1); % note my use of randn for b here. That will usually generate a problem with no solution.
[V1,V2]
prob.Constraints.C1 = V1*a(1) - V1*a(2) + V2*b(1) - V2*b(2) >= ones(n,1);
solve(prob)
So no solution can ever exist for that pair of vectors.
Would other methods exist to solve the problem? Well, yes. We could probably write this where
a = r*cos(theta)
b = r*sin(theta)
Essentially, we want to solve the problem now in polar coordinates. So, does a solution exist of that form? If r>0, then we will be solving the problem for the case where the solution lies on a circle of radius r. That too would, I think, be solvable, under a different norm on the vector (a,b). In this case, it would be the 2-norm. As well, a nice thing about the polar coordinate solution is it allows us to find the set of ALL solutions, and to so so moderately easily. In this last case, we knew no solution exists. How can we learn that is the case using the polar coordinate approach? As a reminder, here are the vectors V1 and V2.
A = [V1,V2]
syms r theta real
ax = r*cos(theta)
by = r*sin(theta)
Now, form the desired sum, as
V3 = ax*V1 + by*V2;
vpa(V3,5)
We could trivially factor out r from there. Again, remember we can always scale any solution by any positive constant. But as well, we can recombine each element of those vectors into one term, using an old trig identity I hope you remember from high school trig.
The one you need to use is:
sin(u +/- v) = sin(u)*cos(v) +/- cos(u)*sin(v)
Essentially, this allows us to re-write an expression with a sum of sines plus cosines, into a single term. For example, we would have this result:
combine(0.8*cos(theta) + 0.3*sin(theta),'sincos')
Can we use that idea here? Of course. A little algebra would give us this:
V3 = r*cos(theta - atan2(V2,V1)./sqrt((V1.^2 + V2.^2)));
vpa(V3,5)
Does that make sense? Can we use it? OF COURSE! On what interval is cos(u) greater than zero? This happens when u is in the interval [-pi/2,pi/2]. What does that tell us? It tells us that We need to find theta such that both of these inequalities always hold:
theta - atan2(V2,V1) > -pi/2
AND
theta - atan2(V2,V1) < pi/2
If such a value of theta exists, then a solution will exist. AND even nicer, it gives us an entire interval on theta, such that ANY value of theta in that interval, a solution will exist. So all we need to do now is to find the intersection of all of those intervals. If that intersection is non-empty, then we have the complete solution set. We can see how this works now.
T = [max(-pi/2 + atan2(V2,V1)),min(pi/2 + atan2(V2,V1))]
If T(1) is less than T(2), then a solution exists, and then the solution set is given by all values of theta in that interval. As you can see here, T(1) > T(2), so no solution was possible for this problem. The test for a solution is nothing more than diff(T) > 0.
Going back to a case where I know a solution will always exist, try this variation again, to insure a solution must exist:
V1 = randn(8,1);
V2 = rand(8,1);
T = [max(-pi/2 + atan2(V2,V1)),min(pi/2 + atan2(V2,V1))]
Do you see that now the interval for theta allows us a valid solution? ANY value of theta in that interval (Remember, T is in radians) is a solution. So I'll just choose the mean of T, as the most logical choice.
Tmean = mean(T)
Asol = cos(Tmean)
Bsol = sin(Tmean)
Asol*V1 + Bsol*V2
And a solution exists, as predicted. Note that linprog was NEVER even necessary using this formulation. The nice thing is, identifying a solution took only one line to compute T. If diff(T)>0, then a solution exists.
ANY other choice in that interval is as equally valid a solution. For example, here the interval on theta was [1.4343,1.8854]. So theta=1.5 would have been a solution too. Anything in that interval should work. If you go near the endpoints of the interval, then one of the elements of V3 will probably start to wander closer to 0. With some additional effort, we might even be able to choose the best value from that interval, under some definition of best. The question at hand did not specify what the best possible solution would be in this context.
Hmm. I'm still wondering if I just did your homework assignment. I hope not. Too bad if it was, because this would have been a nice homework problem to assign. As you can see, the slack variable solution using linprog is nice, but the solution using a trig formulation is nicer yet. It gives you directly the entire set of values of theta for which a solution will exist, and there is no linear programming problem that needs be solved. Optimization is a good thing. But mathematics is sometimes better than resorting to optimization.
5 个评论
Bruno Luong
2023-3-10
Why don't you use Torsen's formulation. Coding in done in 4 lines, extensible in n vectors and simple enough no?
John D'Errico
2023-3-20
编辑:John D'Errico
2023-3-20
Here is an interestign way of looking at the problem for more than 2 variables.
V = randn(6,3)
V =
0.66973 -2.2619 0.51456
0.11795 0.65938 -0.41276
-1.1796 0.51772 -1.0951
-0.74004 -1.8528 0.66123
-0.86157 -0.42123 0.54267
0.50844 -0.42586 -0.11782
X = randn(3,1e7);
VX = any(V*X <= 0,1);
X(:,VX) = [];
size(X)
ans =
3 14467
So only 14467 solutions found, out of an initial sample size of 1e7. Normalize then to the unit sphere.
Xnorm = normalize(X,'norm');
One arbitrarily chosen member of the set does this:
V*Xnorm(:,1)
ans =
0.24612
0.045562
1.1226
0.43443
7.9712e-05
0.099349
plot3(Xnorm(1,:),Xnorm(2,:),Xnorm(3,:),'.')

Normalized to the surface of a unit 3-sphere. It is a small patch, where the boundaries of the patch are as I suggested, great circles on the sphere. From another point of view, we can see it as a polygonal region. Again, this must be, since if you have any two valid solutions to the problem, then any convex linear combination of the two solutions will also be a solution.

In higher dimensions though, is it simplest to use linprog to find a solution? Probably yes. I still like the slack variable formulation I wrote, but any solution that insures the minimum element of the product is 1 will do.
Vsol = Vpos(V)
Solving problem using linprog.
Optimal solution found.
Vsol =
-7.8283
-7.2159
-16.187
>> V*Vsol
ans =
2.7496
1
23.226
8.4597
1
1
For smaller problems with only 3 variables, there may be a method one could put together that would construct a polygonal region of all solutions, but I'm not sure I see the need for that. Only you know if there would be any value.
更多回答(1 个)
Torsten
2023-3-9
编辑:Torsten
2023-3-9
n = 5;
v1 = -10 + rand(1,n)*20
v2 = -20 + rand(1,n)*40
% If there is a solution x for which A*x > zeros(size(A,1),1), then there is also one with
% A*x >= ones(size(A,1),1)
% Thus solve feasibility problem A*x >= ones(size(A,1),1)
A = -[v1;v2].';
b = -ones(n,1);
f = [0;0];
x = linprog(f,A,b)
if ~isempty(x)
x(1)*v1 + x(2)*v2
end
8 个评论
Bruno Luong
2023-3-10
编辑:Bruno Luong
2023-3-10
No, I'm not sure. But it looks to me we need for 2 angles in 3D, 3 angles in 4D etc...
In 3D probably express formulation with quaternion can give some semi-analytic formula. I have impression it is more simpler than trigonometric. But I guess more complex than 2D... I'm not sure with ND.
LINPROG extension however is quite straightforward.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Linear Programming and Mixed-Integer Linear Programming 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!