Reformulate quadratic program to linear program

I have a quadratic program:

minimize

want to have a linear program by substitution:

and after substitution I have the following constraints:

,

My are binary integer values (just 0 or 1) and also my should be binary as well.

I have no idea how to implement the new constraints.

Best regards !

Answer 1

[Matt J]

Organize x into an N-vector X=x(:) and organize q into an NxN matrix Q=reshape(q,[N,N]).

X=optimvar('X',[N,1],'Type','integer','LowerBound',0,'UpperBound',1);
Q=optimvar('Q',[N,N],'Type','integer','LowerBound',0,'UpperBound',1);

Now you can use the problem-based approach, expressing the constraints as follows:

    e=ones(N,1);
Constraint1= Q<=e*X.';
Constraint2= Q<=X*e.';
Constraint3= Q>=e*X.'+X*e.' - 1;

[Matt J]

Kernel7364's comment moved here:

Thank you very much for the quick answer !

But I didnt get it at all. How do I have to use intprog now ?

[x,fval] = intprog(p,.....)

I am sorry that I didnt gave you all the information, but I also have constraints to my vector x:

And now I have this constraints to my x, but I have to make a substitution and the constraints to x still hold but additionally there are these 4 constraints to my q ... I have no clue how to do this and soory that I dont understand your first answer.

[Matt J]

How do I have to use intprog now ?

You use solve(), as explained at the link I gave you.

I am sorry that I didnt gave you all the information, but I also have constraints to my vector x:

You can freely add linear constraints to either X or Q:

Constraint4= A*X<=b;
Constraint5= Aeq*X==beq;

[Kernel7364]

Hey Matt,

thanks again for your help. It seems like you will rescue me today, but I am an absolut beginner in matlab and I always get some errors. That is my code:

Y = [7,3,10;5,0,12;10,6,1];
T= [0,2,4,6;2,0,3,5;4,3,0,7;6,5,7,0];
prob = optimproblem('ObjectiveSense','minimize')
X=optimvar('X',[size(Y,1)*size(T,1),1],'Type','integer','LowerBound',0,'UpperBound',1);
Q=optimvar('Q',[size(Y,1)*size(T,1),size(Y,1)*size(T,1)],'Type','integer','LowerBound',0,'UpperBound',1);
prob.Objective = X*M*X'

    e=ones(size(Y,1)*size(T,1),1);
prob.Constraints1= Q <= e*X.';
prob.Constraints2= Q <= X*e.';
prob.Constraints3= Q >= e*X.'+X*e.' - 1;
prob.Constraints4= A*X <= b;
prob.Constraints5= Aeq*X==beq;

sol = solve(prob)

But I always get this error

prob = 

  OptimizationProblem with properties:

       Description: ''
    ObjectiveSense: 'minimize'
         Variables: [0×0 struct] containing 0 OptimizationVariables
         Objective: [0×0 OptimizationExpression]
       Constraints: [0×0 struct] containing 0 OptimizationConstraints

  No problem defined.

Error using optim.internal.problemdef.MatrixOperator
Incorrect dimensions for matrix multiplication. Check that the number of columns in the first matrix matches the
number of rows in the second matrix. To perform elementwise multiplication, use '.*'.

Error in optim.internal.problemdef.Mtimes

Error in *

Error in Untitled3 (line 10)
prob.Objective = X*M*X'
 

I tried to use your tipps and the tipps from the link you gave me, but I dont get it.

Best regards !

[Matt J]

I think it's just this line. The transpose was in the wrong place.

prob.Objective = X.'*M*X;

[Matt J]

Also, the constraints need to be input like this

prob.Constraints.con1= Q <= e*X.';
prob.Constraints.con2= Q <= X*e.';
prob.Constraints.con3= Q >= e*X.'+X*e.' - 1;
prob.Constraints.con4= A*X <= b;
prob.Constraints.con5= Aeq*X==beq;

The particular names con1,con2,..con5 are not important and could be any other names that you like.

[Kernel7364]

I think it has the problem with

prob.Constraints.con4= A*X <= b;
prob.Constraints.con5= Aeq*X==beq;

It says that

prob = 

  OptimizationProblem with properties:

       Description: ''
    ObjectiveSense: 'minimize'
         Variables: [0×0 struct] containing 0 OptimizationVariables
         Objective: [0×0 OptimizationExpression]
       Constraints: [0×0 struct] containing 0 OptimizationConstraints

  No problem defined.

Error using <=
Argument dimensions 4-by-1 and 1-by-4 must agree.

Error in Untitled3 (line 17)
prob.Constraints.con4= A*X <= b;

I am wondering becaus A is matrix of dimension (4x12) and X is a vector of dimension (12x1) and b is of dimension (4x1). That should actually work.

And I am also wondering that there is "No problem defined."

I hope if I fix the error with the dimension in

prob.Constraints.con4= A*X <= b; and I also guess in

prob.Constraints.con5= Aeq*X <= beq;

that the program will work.

But at the moment I cant see where the error with the dimension is.

[Kernel7364]

And I am also wondering because I never implemented that

I think I have to tell matlab that in any way, right ?!

I am sorry if there are some stupid questions, but that is the first time that I want to solve optimization problems and espacially that I use the solve() function.

I am very thankful about your advices and your help.

[Matt J]

The following code (with some fake data) does run for me without errors. See what happens when you plug in the actual M, A,b, Aeq,beq.

Y = [7,3,10;5,0,12;10,6,1];
T= [0,2,4,6;2,0,3,5;4,3,0,7;6,5,7,0];

A=rand(4,12);b=rand(4,1); %Fake additional data
[Aeq,beq]=deal(A,b);
M=rand(12);

N=size(Y,1)*size(T,1);


prob = optimproblem('ObjectiveSense','minimize');
X=optimvar('X',[N,1],'Type','integer','LowerBound',0,'UpperBound',1);
Q=optimvar('Q',[N,N],'Type','integer','LowerBound',0,'UpperBound',1);
prob.Objective = sum(sum(Q.*M));

    e=ones(N,1);
    
prob.Constraints.con1= Q <= e*X.';
prob.Constraints.con2= Q <= X*e.';
prob.Constraints.con3= Q >= e*X.'+X*e.' - 1;
prob.Constraints.con4= A*X <= b;
prob.Constraints.con5= Aeq*X==beq;


sol = solve(prob)

[Kernel7364]

I cant belive it, it is working !!!

Thank you very much !!

LP:                Optimal objective value is 0.000000.                                             

Cut Generation:    Applied 2 strong CG cuts,                                                        
                   and 2 zero-half cuts.                                                            
                   Lower bound is 0.000000.                                                         

Heuristics:        Found 1 solution using diving.                                                   
                   Upper bound is 150.000000.                                                       
                   Relative gap is 99.34%.                                                         

Heuristics:        Found 1 solution using rss.                                                      
                   Upper bound is 148.000000.                                                       
                   Relative gap is 99.33%.                                                         

Cut Generation:    Applied 42 implication cuts.                                                     
                   Lower bound is 126.000000.                                                       
                   Relative gap is 0.00%.                                                          


Optimal solution found.

Intlinprog stopped at the root node because the objective value is within a gap tolerance of the optimal value,
options.AbsoluteGapTolerance = 0 (the default value). The intcon variables are integer within tolerance,
options.IntegerTolerance = 1e-05 (the default value).


sol = 

  struct with fields:

    Q: [12×12 double]
    X: [12×1 double]

But why does it said that the Optimal objective value is 0.000000 ????

The solution is 126 and that what it also said. I also get the right vector X in sol. Everything is perfectly working.

And hopefully my last question: Why do I dont have to define ??????

[Matt J]

But why does it said that the Optimal objective value is 0.000000 ????

If you know you are getting the correct X, just plug the optimal X into the objective function to verify its value.

And hopefully my last question: Why do I dont have to define ????

Your constraints

,

and the fact that q and x are binary already imply this. You can verify this by checking if your solutions for X and Q satisfy Q=X*X.'

[Kernel7364]

Yes when I check it with the objective function everything is correct, I was just wondering about the 0,0000.

I am using this one for the first time thats why I dont understand everything.

What does CG cuts mean in this comment ?

Cut Generation: Applied 2 strong CG cuts,

and 2 zero-half cuts.

Lower bound is 0.000000.

And when I try to formulate this as a quadratic problem, then it writes:

Error using optim.problemdef.OptimizationProblem/solve
Quadratic problem with integer variables not supported.

Thats why I have to define my q. I start to understand everything.

Thank you so much for your help ! You are an absolut hero !

I will ask some more questions if I try something out and dont understand it. :)

[Matt J]

You're welcome, but please Accept-click the answer to signify that your question has been addressed.

[Kernel7364]

It takes so long when I use it for a bigger matrix.

Do you know how to solve this with cplex ?

I hope that this is faster when I do it with cplex.

Best regards !

[Matt J]

Kernel7364's comment moved here:

Hey,

i tried this one out with some bigger matrix and this one took almost 15 minutes to solve it and it is even a wrong solution.

Thats why I would like to solve it with Cplex. I have heard that this program is much better and faster.

But I have no idea how to use the Cplex solver to solve it. Where do I have to write Cplex in my program ?

How is it possible that I get a wrong solution and that it takes 15 minutes if my N is 140 ?

(T is a 14x14 matrix and Y 14x14)

Best regards !

[Matt J]

I'm afraid I have no experience with CPLEX, but this thread looks like it might have releveant information

https://www.mathworks.com/matlabcentral/answers/6291-call-cplex-from-matlab

As for why a wrong solution is reached, you should check

(1) The exitflag output of the solver - mabe it terminated prematurely.

(2) That there may be multiple solutions - see if the final objective values are nearly the same for both the expected solution and the one you obtained.

[Kernel7364]

Oh no, if you dont no a solution to my problem I am really lost.

I read this older chat and I tried the advices there but it is still not running. :(

[Matt J]

What do you mean, "don't know"? I gave you 3 pieces of advice...

[Kernel7364]

Sorry for the missunderstanding.

Your advice for the wrong solution makes sense and this helped me. After around 15 minutes the solver finished and the exitflag was -6 at first so I knew that there was a problem with it.

I am trying to find a way how to tell matlab that he please has to use CPLEX for my problem and not the solver in matlab itself. I just want a command for this but I am still looking. I know how to do it if I have a linear problem, then I have to call

cplexlp(f,A,b,Aeq,beq,lb,ub)

And Matlab uses cplex. But when i use the solve() function, I dont know how to tell matlab to solve it with Cplex.