[math] Re: solving a linear optimization problem with commons math ?

[Luc Maisonobe <lu...@apache.org>]

Andreas Andreakis a écrit :
> hi,
> 
> sorry to send you an email regarding commons math directly. I tried to 
> subscribe to the mailing list, but for some reason, I was considered as 
> spam.
> Thus I chose two members randomly.
> 
> Well, anyway :) I want to find a solution to the following problem:
> 
> (explained by example)
> 1 *c1 + 4 * c2 + 4* c3 + 4 * c4 = 2.5
> 1 *c1 + 4 * c2 + 1* c3 + 4 * c4 = 3.0
> 
> here the constraints:
> 1) c1+c2+c3+c4 = 1 (of course this could be expressed as an additional 
> equation: c1 + c2 + c3 + c4 = 1.0)
> 2) c1,c2,c3 and c4 are all positive
> 
> Important is that I need only one and feasible solution, thus the 
> computation should be as fast as possible. So there is no need to search 
> for every possible result.
> And it needs to work with doubles (not only integers, because some 
> libraries work only with integers)
> 
> Actually the solution needs to scale for several dozens c´s ( there will 
> not be only 4 like c1,c2,c3,c4), but more like 50 c1, c2,...,c49, c50
> And there will not be only two equations (like in the example) but 
> several hundred.
> 
> 
> So the questions are:
> 1) can I find a solution using commons math. There is an optimization 
> package, but I could not figure out how to use it?

In this case, the estimation package is more suited to your needs. There 
are much better methods for such problems, but they are not available in 
commons-math.

The Levenberg-Marquardt solver implemented in thie estimation package 
should support both under- and over-determined problems (i.e. either 
more equations than variables or more variables than equations). The 
case with more variable than equations is still considered experimental, 
though. I did use it for such non homogeneous linear problems without 
problems myself.

This solver does not completely fit your needs as is because it 
currently do not handle constraints. You will have to set up an 
intermediate unconstrained problem that can be solved and from which the 
solution of your constrained problem can be derived. Here are a few 
tricks to do that, take care, they may induce numerical problems:

For equality constraints, either add the constraint to the equation set 
or remove one parameter (say c1) from the parameter set and compute it 
as needed by solving the constraint each time from the values of the 
other ones, i.e. compute c1 = 1 - f(c2, c3, ..., cn).

For simple bounds constraints on parameters (such as ck > ak), replace 
the real parameter ck by an unconstrained parameter pk using
   ck = ak + pk² (or ak - pk² if the constraint was ck < ak)

For interval constraints on parameters (such as lowk < ck < highk), 
replace the real parameter ck by an unconstrained parameter pk using
   ck = [ lowkk + highk ] / 2 + atan(pk) * [ highk - lowk ] / pk

Once you have changed your problem using these pk, it is an 
unconstrained problem that can be solved by our implementation of 
Levenberg-Marquardt. Once you have the pk values, you compute the ck values.


> 2) how will the performance be ?

Since the problem is linear, the performance should be fair.

> 3) could you give me a code snipped which is able to compute the given 
> example ? (I will generalize it then myself)

Try something along these lines. Take care, I wrote this exemple 
directly in the mail, I did not attempt to even compile it ;-)

Hope this helps,
Luc

import org.apache.commons.math.estimation.EstimatedParameter;
import org.apache.commons.math.estimation.EstimationException;
import org.apache.commons.math.estimation.LevenbergMarquardtEstimator;
import org.apache.commons.math.estimation.SimpleEstimationProblem;
import org.apache.commons.math.estimation.WeightedMeasurement;


public class LinearConstrainedProblem extends SimpleEstimationProblem {

     public static void main(String[] args) {

         try {
             // create the uninitialized problem
             LinearConstrainedProblem pb =
                new LinearConstrainedProblem(0.25, 0.25, 0.25, 0.25)

             // 1 *c1 + 4 * c2 + 4* c3 + 4 * c4 = 2.5
             pb.addEquation(1.0, 4.0, 4.0, 4.0, 2.5);

             // 1 *c1 + 4 * c2 + 1* c3 + 4 * c4 = 3.0
             pb.addEquation(1.0, 4.0, 1.0, 4.0, 3.0);

             // solve the problem, using a Levenberg-Marquardt algorithm
             // with default settings
             LevenbergMarquardtEstimator estimator =
               new LevenbergMarquardtEstimator();
             estimator.estimate(pb);
             System.out.println("RMS after solve = "
                              + estimator.getRMS(pb));

             // display the results
             System.out.println("c1 = " + pb.getC1()
                     + ", c2 = " + pb.getC2()
                     + ", c3 = " + pb.getC3()
                     + ", c4 = " + pb.getC4());

         } catch (EstimationException ee) {
             System.err.println(ee.getMessage());
         }

     }

     /** Fitter for a linear constrained model.
      * @param c1Guess start value for c1 (in fact, ignored)
      * @param c2Guess start value for c2
      * @param c3Guess start value for c3
      * @param c4Guess start value for c4
      */
     public LinearConstrainedProblem() {

         // one parameter (c1) is computed from the other ones
         // and the other ck are derived from pk using ck = pk*pk
         p2 = new EstimatedParameter("p2", Math.sqrt(c2Guess));
         p3 = new EstimatedParameter("p3", Math.sqrt(c3Guess));
         p4 = new EstimatedParameter("p4", Math.sqrt(c4Guess));

         // provide the parameters to the base class which
         // implements the getAllParameters and
         // getUnboundParameters methods
         addParameter(c2);
         addParameter(c3);
         addParameter(c4);

     }

     /** Add an equation.
      * @param a1 coeff for c1
      * @param a2 coeff for c2
      * @param a3 coeff for c3
      * @param a4 coeff for c4
      * @param rhs right hand side of the equation
      */
     public void addEquation(double a1, double a2, double a3,
                             double a4, double rhs) {
         addMeasurement(new Equation(a1, a2, a3, a4, rhs));
     }

     public double getC1() {
         // this solves the equality constraint
         return 1.0 - getC2() - getC3() - getC4();
     }

     public double getC2() {
       double p = p2.getEstimate();
       return p * p;
     }

     public double getC3() {
       double p = p3.getEstimate();
       return p * p;
     }

     public double getC4() {
       double p = p4.getEstimate();
       return p * p;
     }


     /** Internal measurements class, representing one equation. */
     private class Equation extends WeightedMeasurement {

         public Equation(double a1, double a2, double a3,
                         double a4, double rhs) {
             super(1.0, rhs);
             this.a1 = a1;
             this.a2 = a2;
             this.a3 = a3;
             this.a4 = a4;
         }

         public double getTheoreticalValue() {
             return a1 * getC1() + a2 * getC2()
                  + a3 * getC3() + a4 * getC4();
         }

         public double getPartial(EstimatedParameter parameter) {
           if (parameter == p2) {
             return 2.0 * p2.getEstimate() * (a2 - a1);
           } else if (parameter == p3) {
             return 2.0 * p3.getEstimate() * (a3 - a1);
           } else {
             return 2.0 * p4.getEstimate() * (a4 - a1);
           }
         }

         private final double a1;
         private final double a2;
         private final double a3;
         private final double a4;

     }

     private EstimatedParameter p2;
     private EstimatedParameter p3;
     private EstimatedParameter p4;

}

> 
> kind regards,
> Andreas



---------------------------------------------------------------------
To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
For additional commands, e-mail: user-help@commons.apache.org

Re: [math] Re: solving a linear optimization problem with commons math ?

[chengas123 <be...@gmail.com>]

Since this is a linear problem you can use SimplexSolver (only available in
Commons Math 2.0).  
You can reformulate the problem to the equivalent:
Maximize C1 (you can put anything here since in your case you don't care
about the maximization aspect of the problem but only the constraint
satisfaction)
Subject To
1 * c1 + 4 * c2 + 4 * c3 + 4 * c4 = 2.5
1 * c1 + 4 * c2 + 1 * c3 + 4 * c4 = 3.0
1 * c1 + 1 * c2 + 1 * c3 + 1 * c4 = 1.0 



> 1 *c1 + 4 * c2 + 4* c3 + 4 * c4 = 2.5
> 1 *c1 + 4 * c2 + 1* c3 + 4 * c4 = 3.0
> 
> here the constraints:
> 1) c1+c2+c3+c4 = 1 (of course this could be expressed as an additional
> equation: c1 + c2 + c3 + c4 = 1.0)
> 2) c1,c2,c3 and c4 are all positive 
> 

-- 
View this message in context: http://www.nabble.com/-math--Re%3A-solving-a-linear-optimization-problem-with-commons-math---tp14624075p23548298.html
Sent from the Commons - User mailing list archive at Nabble.com.


---------------------------------------------------------------------
To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
For additional commands, e-mail: user-help@commons.apache.org

Re: [math] Re: solving a linear optimization problem with commons math ?

[Luc Maisonobe <Lu...@free.fr>]

Andreas Andreakis a écrit :
> Hi Luc,
> 
> thank you very much for your quick and very detailed explanation. I´ll 
> try out your suggestions very soon.
> 
> You mentioned, that there are better method in order to deal with such 
> problems. Which are those methods ?

Constraint Programming. This is a mathematics fields by itself.

> Can you recommend a library which implements them ?

The major player in this field is ILOG. This company provides various 
proprietary libraries in various languages. The libraries are 
efficients, the guys behind them are *very* good, and the libraries are 
... really expensive. There are other commercial solutions apart from ILOG.

On the open-source side, you may try Choco (http://choco-solver.net/). I 
never used it myself but would be very interested if you could provide 
me feedback about it.

> 
> P.S.: don´t worry, I´m using commons math for other kind of calculations 
> anyway ;))

I am aware the solution I provided to you was a simple dirty trick, it 
is neither elegant nor efficient. By the way, rather than ck = pk², I 
now think ck = exp(pk) would be a better variable change (but it still 
is a dirty trick).

Luc


---------------------------------------------------------------------
To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
For additional commands, e-mail: user-help@commons.apache.org