yalmip.rar_1223.246domain com_by.1239 com_sdpvar\eig_yalmip lmi_

function output = bmibnb(p) %BMIBNB Branch-and-bound scheme for bilinear programs % % BMIBNB is never called by the user directly, but is called by % YALMIP from SOLVESDP, by choosing the solver tag 'bmibnb' in sdpsettings % % The behaviour of BMIBNB can be altered using the fields % in the field 'bmibnb' in SDPSETTINGS %s % WARNING: THIS IS EXPERIMENTAL CODE % % bmibnb.lowersolver - Solver for lower bound [standard solver tag ('')] % bmibnb.uppersolver - Solver for upper bound [standard solver tag ('')] % bmibnb.lpsolver - Solver for LP bound tightening [standard solver tag ('')] % bmibnb.branchmethod - Branch strategy ['maxvol' | 'best' ('best')] % bmibnb.branchrule - Branch position ['omega' | 'bisect' ('omega')] % bmibnb.lpreduce - Improve variable bounds using LP [ real [0,1] (0 means no reduction, 1 means all variables) % bmibnb.lowrank - partition variables into two disjoint sets and branch on smallest [ 0|1 (0)] % bmibnb.target - Exit if upper found<target [double (-inf)] % bmibnb.roottight - Improve variable bounds in root using full problem [ 0|1 (1)] % bmibnb.vartol - Cut tree when x_U-x_L < vartol on all branching variables % bmibnb.relgaptol - Tolerance on relative objective error (UPPER-LOWER)/(1+|UPPER|) [real (0.01)] % bmibnb.absgaptol - Tolerance on objective error (UPPER-LOWER) [real (0.01)] % bmibnb.pdtol - A number is declared non-negative if larger than...[ double (-1e-6)] % bmibnb.eqtol - A number is declared zero if abs(x) smaller than...[ double (1e-6)] % bmibnb.maxiter - Maximum number of solved nodes [int (100)] % bmibnb.maxtime - Maximum CPU time (sec.) [int (3600)] % Author Johan L�fberg % $Id: bmibnb.m,v 1.64 2006/03/28 21:47:09 joloef Exp $ % ************************************************************************* % INITIALIZE DIAGNOSTICS ETC % ************************************************************************* bnbsolvertime = clock; showprogress('Branch and bound started',p.options.showprogress); switch max(min(p.options.verbose,3),0) case 0 p.options.bmibnb.verbose = 0; case 1 p.options.bmibnb.verbose = 1; p.options.verbose = 0; case 2 p.options.bmibnb.verbose = 2; p.options.verbose = 0; case 3 p.options.bmibnb.verbose = 2; p.options.verbose = 1; otherwise p.options.bmibnb.verbose = 0; p.options.verbose = 0; end % ******************************* % Assume feasible % ******************************* p.feasible = 1; % ******************************* % Pre-calc linear/bilinear/nonlinear variables % (assumes polynomial program) % ******************************* [p.linears,p.bilinears,p.nonlinears] = compile_nonlinear_table(p); % ******************************* % Select branch variables % ******************************* p.branch_variables = decide_branch_variables(p); % ******************************* % Now reduce the branch variables by checking which quadratic % which only are used in the convex cost % This is experimental code.... % ******************************* if p.solver.lowersolver.objective.quadratic.convex % Setup quadratic Q_ = p.Q; for i = 1:size(p.bilinears,1) if p.c(p.bilinears(i,1)) Q_(p.bilinears(i,2),p.bilinears(i,3)) = p.c(p.bilinears(i,1))/2; Q_(p.bilinears(i,2),p.bilinears(i,3)) = Q_(p.bilinears(i,3),p.bilinears(i,2))+p.c(p.bilinears(i,1))/2; end end if nnz(Q_)>0 & all(eig(full(Q_))>-1e-12) Used_in_F = find(any(p.F_struc(:,2:end),1)); Used_in_F = intersect(Used_in_F,p.bilinears(:,1)); p.branch_variables = []; for i = 1:size(p.bilinears,1) j = p.bilinears(i,1); if ismember(j,Used_in_F) p.branch_variables = [p.branch_variables p.bilinears(i,2:3)]; end end p.branch_variables = unique( p.branch_variables); end end % ******************************** % Extract bounds from model % ******************************** p = preprocess_bounds(p); %p.branch_variables = setdiff(p.branch_variables ,find(p.lb == -10)); % ******************************** % Is the initial value feasible % ******************************** [p,x_min,upper] = initializesolution(p); % ******************************** % Simple pre-solve (empty rows etc) % Bound strenghtening % The loop is needed for variables % defined as w = x*y, x = t*u,y=.. % ******************************** start = [p.lb;p.ub]+1;i = 0;;goon = 1; while goon start = [p.lb;p.ub]; i = i+1; p = simple_presolve(p); p = updatenonlinearbounds(p); p = propagateequalities(p); p = updatenonlinearbounds(p); p = updatenonlinearbounds(p); i = i+1; goon = (norm(start-[p.lb;p.ub],'inf') > 1e-1) & i < 25; start = [p.lb;p.ub]; end % ******************************** % Default root bound improvement % ******************************** p = root_node_tighten(p,upper); p = updatenonlinearbounds(p); p = propagateequalities(p); p = updatenonlinearbounds(p); if p.feasible % ******************************* % Bounded domain? % ******************************* if ~isempty(p.bilinears) involved = unique(p.bilinears(:,2:3)); if any(isinf(p.lb(involved))) | any(isinf(p.ub(involved))) output.Primal = x_min; output.Dual = []; output.Slack = []; output.problem = -6; output.infostr = yalmiperror(-6); output.solved_nodes = 0; output.solverinput = 0; output.solveroutput = []; output.solvertime = 0; return end end % ******************************* % Save time & memory % ******************************* p.options.savesolverinput = 0; p.options.savesolveroutput = 0; p.options.saveduals = 0; p.options.dimacs = 0; p.options.dimacs = 0; % ******************************* % RUN BRANCH & BOUND % ******************************* [x_min,solved_nodes,lower,upper] = branch_and_bound(p,x_min,upper); % ********************************** % CREATE SOLUTION % ********************************** output.problem = 0; if isinf(upper) output.problem = 1; end if isinf(lower) & (lower<0) output.problem = 2; end if isinf(lower) & (lower>0) output.problem = 1; end if solved_nodes == p.options.bnb.maxiter output.problem = 3; end else output.problem = 1; x_min = repmat(nan,length(p.c),1); solved_nodes = 0; end output.solved_nodes = solved_nodes; output.Primal = x_min; output.Dual = []; output.Slack = []; output.infostr = yalmiperror(output.problem,'BNB'); output.solverinput = 0; output.solveroutput = []; output.solvertime = etime(clock,bnbsolvertime); function [x_min,solved_nodes,lower,upper] = branch_and_bound(p,x_min,upper) % *************************************** % Create function handles to solvers % *************************************** try lowersolver = eval(['@' p.solver.lowersolver.call]); % Solver for relaxed problem uppersolver = eval(['@' p.solver.uppersolver.call]); % Local nonlinear solver lpsolver = eval(['@' p.solver.lpsolver.call]); % LP solver catch disp(' '); disp('The internal branch & bound solver requires MATLAB 6') disp('I am too lazy too do the changes to make it compatible') disp('with MATLAB 5. If you really need it, contact me...'); disp(' '); error(lasterr); end % ************************************************ % GLOBAL PROBLEM DATA % ************************************************ c = p.c; Q