Hola, Necesito resolver un problema de optimización por el método de ramificación y poda(branch and bound). No se nada sobre matlab y encontré este programa en internet, pero no se ni si quiera como meter las variables para ver si funciona o no.
Espero me puedan ayudar y decirme como puedo meter las variables.
Este es el problema a optimizar (Minimizar la función objetivo):
Min F= x1-2*x2
subject to
-4x1+6x2 <= 9
x1+x2 <= 4
x1,x2 >= 0
Este es el código:
% [x,val,status]=MILP_bb(c,A,b,M,eps)
% this function solves a mixed-integer linear programming problem
% using the branch and bound algorithm
% The code uses Kartik's revised simplex and dual simplex codes
% to solve the LP relaxations at each node of the branch and bound tree.
% max c'*x
% subject to
% A*x = b
% x >= 0
% M is a vector of indices for the variables that are constrained to be integers
% eps is the user tolerance
% the output variables are :
% x : the optimal solution
% val: value of the objective function at the optimal solution
% status =1 if successful
% =2 if the original problem is infeasible
% Kartik's MATLAB code
% Last updated: 19th April, 2006
% The B&B driver routine based on a MATLAB code by Sherif Tawfik, Cairo University
% downloaded from the MATLAB File Exchange.
function [x,val,status]=MILP_bb(c,A,b,M,eps)
[m,n] = size(A);
% Use the phase 1 procedure to find an initial basic feasible solution
% for the LP relaxation at the root node.
Aphase1 = A;
bphase1 = b;
index = find(b < 0);
bphase1(index) = -b(index);
Aphase1(index,:) = -A(index,:);
Aphase1 = [Aphase1 eye(m)];
cphase1 = [zeros(n,1); -ones(m,1)];
x0phase1 = [zeros(n,1); bphase1];
B0 = [n+1:n+m]';
[val,x0,B] = revised_simplex_bb(cphase1,Aphase1,bphase1,eps,x0phase1,B0);
clear Aphase1 cphase1 bphase1 x0phase1 B0 index;
% If the root LP relaxation is infeasible, then just quit!
if val < 0,
x = [];
status = 2;
val = [];
return;
end
bound = -inf; % the initial lower bound is set to -ve infinity
root = 1;
[x,B,status,b]=branch(c,A,b,x0,B,M,eps,bound,root); % a recursive function that processes the BB tree
val = c'*x(1:n); % objective value of the optimal solution
function [xx,B,status,bb]=branch(c,A,b,x,B,M,eps,bound,root)
% x is an initial solution
% Solve the LP relaxation at the current node
% If it is the LP relaxation at the root node use the revised simplex method
% Else use the dual simplex method
if root == 1,
[val0,x0,B] = revised_simplex_bb(c,A,b,eps,x,B);
% The root LP relaxation is solved using Kartik's revised simplex code
% Initial bfs is the one found using the phase 1 procedure.
status0 = 1;
root = 0; % We are now dealing with descendents of the root node
else
[val0,x0,B,status0] = dual_simplex_bb(c,A,b,eps,x,B);
% All LP relaxations at descendent nodes are solved using Kartik's dual
% simplex code with warm-start information.
end
% If the LP relaxation is infeasible then PRUNE THE NODE BY INFEASIBILITY
% If the LP objective value is less than our current integer bound then PRUNE THE NODE BY BOUNDS
if status0 == 2 | val0 < bound
xx=x;
status=status0; bb=bound;
% We are pruning this node, so there is no need to branch further on
% this node
return;
end
% If the solution to the LP relaxation is feasible in the MILP problem, then check its objective value
% against that of the incumbent solution
% If the new feasible solution has a greater objective value then update the lower bound
% In any case, PRUNE THE NODE BY OPTIMALITY
ind=find(abs(x0(M)-round(x0(M))) > eps );
if isempty(ind)
status=1;
if val0 > bound % replace the incumbent solution
x0(M)=round(x0(M));
xx=x0;
bb=val0;
else
xx=x; % we have pruned the node and there is no need to branch further on this node
bb=bound;
end
return;
end
% The solution of the LP relaxation is not feasible in the MILP problem.
% However, the objective value of the LP relaxation is greater than the current lower bound.
% So we branch on this node to create two subproblems.
% We will solve the two subproblems recursively by calling the same branching function.
% The first LP problem with the added constraint that x_i <= floor(x_i) , i=ind(1)
br_var=M(ind(1));
br_value=x0(br_var);
[m,n] = size(A);
A1 = [A zeros(m,1); zeros(1,n+1)];
A1(m+1,br_var) = 1;
A1(m+1,n+1) = 1;
b1 = [b; floor(br_value)];
c1 = [c; 0];
x01 = [x0; b1(m+1)-x0(br_var)];
% The basis size grows by 1 with the new slack variable in the basis
B01 = [B; n+1];
% second LP problem with the added constraint that x_i >= ceil(x_i) , i=ind(1)
A2 = [A zeros(m,1); zeros(1,n+1)];
A2(m+1,br_var) = 1;
A2(m+1,n+1) = -1;
b2 = [b; ceil(br_value)];
c2 = [c; 0];
x02 = [x0; x0(br_var)-b2(m+1)];
% The basis size grows by 1 with the new slack variable in the basis
B02 = [B; n+1];
% solve the first child problem
[x1,B1,status1,bound1]=branch(c1,A1,b1,x01,B01,M,eps,bound,root);
status=status1;
if status1 == 1 & bound1 > bound % if the solution was successful and gives a better bound
xx=x1;
B=B1;
bound=bound1;
bb=bound1;
else
xx=x0;
bb=bound;
end
% solve the second child problem
[x2,B2,status2,bound2]=branch(c2,A2,b2,x02,B02,M,eps,bound,root);
if status2 == 1 & bound2 > bound % if the solution was successful and gives a better bound
status=status2;
xx=x2;
B=B2;
bb=bound2;
end
Espero puedan ayudarme.
Gracias.
