Solving Ax=B using LAPACK, where x >= 0

Question

I'm currently working on an iOS app which handles chemical additions to water. In order to find the smallest possible additions, I'm solving Ax=B where A is a 6x6 matrix and B is one column. As far as I've been reading, LAPACK in the Accelerate Framework is the best way to handle this. I've been able to get it working for small additions using DGESV_, but larger ones have some negative values for x. This is an issue as you can't add a negative amount of a chemical.

So what I need to know is whether there is a function in LAPACK that will allow me to find the minimum solution to Ax=B, where x is always greater than or equal to 0? If not, is there an alternative solution (other than doing the math myself)?

Thank you

Answer 1

Unless A is rank deficient, the solution of Ax = B with A square is unique. So there is no way to get rid of negative components in the answer.

If you assume A and B might contain errors that cause a negative component but would like to find a "nearby" solution with all non-negative components, then you could cast this as:

minimize |Ax - b|  subject to x >= 0

This is a quadratic program. There are libraries to solve such problems, but LAPACK isn't one of them.

Edit

Your matrix is full rank. Wolfram Alpha is very nice for playing with little problems like this. Your 6x6 has a determinant of 5x10^11, so it's very nicely conditioned.

I'm no expert on QP, but this one has special structure. For standard form (see Wikipedia for notation), rewrite as:

minimize ( x'(A'A)x + (-2b'A)x )  subject to x >= 0

Here the quadratic coefficient matrix Q = A'A is symmetric positive definite, which (if I am recalling my matrix algebra correctly), makes the system convex and easier to solve: a single global solution is guarenteed. Note c = -2b'A.

Here is a page of QP libraries, some in C, and some optimized for the convex case. Perhaps one of them will work for you.