compileinterfacedata.m

matlab波形优化算法经常要用到的matlab toolbox工具箱:yalmip

        for i = 1:length(logdetStruct)
            Ptemp = blkdiag(Ptemp,logdetStruct.P{i});
        end
        P = {Ptemp};
        if length(F)>0
            if isequal(P,sdpvar(F(end)))
                F = F(1:end-1);
            end
        end
        F = F + detset(t,P{1});
        if isempty(h)
            h = -t;
        else
            h = h-t;
            % Warn about logdet -> det^1/m
            if options.verbose>0 & options.warning>0
                disp(' ')
                disp('Objective c''x-logdet(P) has been changed to c''x-det(P)^(1/(2^ceil(log2(length(X))))).')
                disp('See the MAXDET section in the manual for details.')
                disp(' ')
            end
        end
        P = [];
    end
end

% *************************************************************************
%% Change binary variables to integer?
% *************************************************************************
old_binary_variables = binary_variables;
if ~isempty(binary_variables) & (solver.constraint.binary==0)
    x_bin = recover(binary_variables(ismember(binary_variables,unique([getvariables(h) getvariables(F)]))));
    F = F + set(x_bin<1)+set(x_bin>0);
    integer_variables = union(binary_variables,integer_variables);
    binary_variables = [];
end

% *************************************************************************
%% Model quadratics using SOCP?
% Should not be done when using PENNLP or BMIBNB or FMINCON, or if we have relaxed the
% monmial terms or...Ouch, need to clean up all special cases, this sucks.
% *************************************************************************
convertQuadraticObjective = ~strcmpi(solver.tag,'pennlp-standard');
convertQuadraticObjective = convertQuadraticObjective & ~strcmpi(solver.tag,'bmibnb');
convertQuadraticObjective = convertQuadraticObjective & ( ~(options.relax==1 | options.relax==3)  &(~isempty(quad_info) & solver.objective.quadratic.convex==0 | (~isempty(quad_info) & strcmp(solver.tag,'bnb') & localsolver.objective.quadratic.convex==0)));
convertQuadraticObjective = convertQuadraticObjective; % | strcmpi(solver.tag,'cutsdp');
if strcmpi(solver.tag,'bnb') & ~isempty(quad_info)
    if solver.lower.objective.quadratic.convex==0
        convertQuadraticObjective = 1;
    end
end
if convertQuadraticObjective
    t = sdpvar(1,1);
    x = quad_info.x;
    R = quad_info.R;
    if ~isempty(R)
        c = quad_info.c;
        f = quad_info.f;
        F = F + lmi(cone([2*R*x;1-(t-c'*x-f)],1+t-c'*x-f));
        h = t;
    end
    quad_info = [];
end
if solver.constraint.inequalities.rotatedsecondordercone == 0
    [F,changed] = convertlorentz(F);
    if changed
        options.saveduals = 0; % We cannot calculate duals since we change the problem
    end
end
% Whoa, horrible tests to find out when to convert SOCP to SDP
% This should not be done if :
%   1. Problem is actually a QCQP and solver supports this
%   2. Problem is integer, local solver supports SOCC
%   3. Solver supports SOCC
if ~((solver.constraint.inequalities.elementwise.quadratic.convex == 1) & socp_are_really_qc)
    if ~(strcmp(solver.tag,'bnb') & socp_are_really_qc & localsolver.constraint.inequalities.elementwise.quadratic.convex==1 )
        if ((solver.constraint.inequalities.secondordercone == 0) | (strcmpi(solver.tag,'bnb') & localsolver.constraint.inequalities.secondordercone==0))
            [F,changed] = convertsocp(F);
            if changed
                options.saveduals = 0; % We cannot calculate duals since we change the problem
            end
        end
    end
end

% *************************************************************************
%% Add logaritmic barrier cost/constraint for MAXDET and SDPT3-4. Note we
% have to add it her in order for a complex valued matrix to be converted.
% *************************************************************************
if ~isempty(logdetStruct) & solver.objective.maxdet==1
    for i = 1:length(logdetStruct.P)
        F = F + set(logdetStruct.P{i} > 0);
        if ~isreal(logdetStruct.P{i})
            logdetStruct.gain(i) = logdetStruct.gain(i)/2;
            ProblemClass.complex = 1;
        end
    end
end

if ((solver.complex==0) & ProblemClass.complex) | ((strcmp(solver.tag,'bnb') & localsolver.complex==0) & ProblemClass.complex)
    showprogress('Converting to real constraints',options.showprogress)
    F = imag2reallmi(F);
    if ~isempty(logdetStruct) 
        for i = 1:length(logdetStruct.P)
            P{i} = sdpvar(F(end-length(logdetStruct.P)+i));
        end
    end
    options.saveduals = 0; % We cannot calculate duals since we change the problem
%else
%    complex_logdet = zeros(length(P),1);
end

% *************************************************************************
%% CREATE OBJECTIVE FUNCTION c'*x+x'*Q*x
% *************************************************************************
showprogress('Processing objective h(x)',options.showprogress);
try
    % If these solvers, the Q term is placed in c, hence quadratic terms
    % are treated as any other nonlinear term
    geometric = strcmpi(solver.tag,'fmincon-geometric')| strcmpi(solver.tag,'gpposy') | strcmpi(solver.tag,'mosek-geometric') | strcmpi(solver.tag,'snopt-geometric') | strcmpi(solver.tag,'ipopt-geometric') | strcmpi(solver.tag,'pennon-geometric');
    if strcmpi(solver.tag,'bnb')
        lowersolver = lower([solver.lower.tag '-' solver.lower.version]);
        if strcmpi(lowersolver,'fmincon-geometric')| strcmpi(lowersolver,'gpposy-') |  strcmpi(lowersolver,'mosek-geometric')
            geometric = 1;
        end
    end
    if strcmpi(solver.tag,'bmibnb') | strcmpi(solver.tag,'pennlp-standard') | geometric | evaluation_based ;
        tempoptions = options;
        tempoptions.relax = 1;
        [c,Q,f]=createobjective(h,logdetStruct,tempoptions,quad_info);
    else
        [c,Q,f]=createobjective(h,logdetStruct,options,quad_info);
    end
catch
    error(lasterr)
end

% *************************************************************************
%% Convert {F(x),G(x)} to a numerical SeDuMi-like format
% *************************************************************************
showprogress('Processing F(x)',options.showprogress);
F = lmi(F);
[F_struc,K,KCut] = lmi2sedumistruct(F);
% We add a field to remember the dimmesion of the logarithmic cost.
% Actually, the actually value is not interesting, we know that the
% logarithmic cost corresponds to the last SDP constraint anyway
if isempty(logdetStruct)
    K.m = 0;
else
    for i = 1:length(logdetStruct.P)
        K.m(i) = length(logdetStruct.P{i});        
    end
    K.maxdetgain = logdetStruct.gain;
end

% *************************************************************************
%% SOME HORRIBLE CODE TO DETERMINE USED VARIABLES
% *************************************************************************
% Which sdpvar variables are actually in the problem
used_variables_LMI = find(any(F_struc(:,2:end),1));
used_variables_obj = find(any(c',1) | any(Q));
used_variables = uniquestripped([used_variables_LMI used_variables_obj]);
% The problem is that linear terms might be missing in problems with only
% nonlinear expressions
[monomtable,variabletype] = yalmip('monomtable');
if (options.relax==1)|(options.relax==3)
    monomtable = [];
    nonlinearvariables = [];
    linearvariables = used_variables;
else
    nonlinearvariables = find(variabletype);
    linearvariables = used_variables(find(variabletype(used_variables)==0));
end
needednonlinear = nonlinearvariables(ismembc(nonlinearvariables,used_variables));
linearinnonlinear = find(sum(abs(monomtable(needednonlinear,:)),1));
missinglinear = setdiff(linearinnonlinear(:),linearvariables);
used_variables = uniquestripped([used_variables(:);missinglinear(:)]);

% *************************************************************************
%% So are we done now? No... What about variables hiding inside so called
% evaluation variables. We detect these, and at the same time set up the
% structures needed to support general functions such as exp, log, etc
% NOTE : This is experimental code
% FIX  : Clean up...
% *************************************************************************
[evalMap,evalVariables,used_variables,nonlinearvariables,linearvariables] = detectHiddenNonlinear(used_variables,options,nonlinearvariables,linearvariables,evalVariables);

% Attach information on the evaluation based variables that was generated
% when the model was expanded
if ~isempty(evalMap)
    for i = 1:length(operators)
        index = find(operators{i}.properties.models(1) == used_variables(evalVariables));
        if ~isempty(index)
            evalMap{index}.properties = operators{i}.properties;
        end
    end
    for i = 1:length(evalMap)
        for j = 1:length(evalMap{i}.computes)
            evalMap{i}.computes(j) = find(evalMap{i}.computes(j) == used_variables);
        end
    end
end

% *************************************************************************
%% REMOVE UNNECESSARY VARIABLES FROM PROBLEM
% *************************************************************************
if length(used_variables)<yalmip('nvars')
    c = c(used_variables);
    Q = Q(:,used_variables);Q = Q(used_variables,:);
    if ~isempty(F_struc)
        F_struc = sparse(F_struc(:,[1 1+used_variables]));
    end
end

% *************************************************************************
%% Map variables and constraints in low-rank definition to local stuff
% *************************************************************************
if ~isempty(lowrankdetails)
    % Identifiers of the SDP constraints
    lmiid = getlmiid(F);
    for i = 1:length(lowrankdetails)
        lowrankdetails{i}.id = find(ismember(lmiid,lowrankdetails{i}.id));
        if ~isempty(lowrankdetails{i}.variables)
            index = ismember(used_variables,lowrankdetails{i}.variables);
            lowrankdetails{i}.variables = find(index);
        end
    end
end

% *************************************************************************
%% SPECIAL VARIABLES
% Relax = 1 : relax both integers and nonlinear stuff
% Relax = 2 : relax integers
% Relax = 3 : relax nonlinear stuff
% *************************************************************************
if (options.relax==1) | (options.relax==3)
    nonlins = [];
end
if (options.relax) | (options.relax==2)
    integer_variables = [];
    binary_variables  = [];
    old_binary_variables  = find(ismember(used_variables,old_binary_variables));
else
    integer_variables = find(ismember(used_variables,integer_variables));
    binary_variables  = find(ismember(used_variables,binary_variables));
    old_binary_variables  = find(ismember(used_variables,old_binary_variables));
end
parametric_variables  = find(ismember(used_variables,parametric_variables));
extended_variables =  find(ismember(used_variables,yalmip('extvariables')));


% *************************************************************************
%% Equality constraints not supported or supposed to be removed
% *************************************************************************
% We may save some data in order to reconstruct
% dual variables related to equality constraints that
% have been removed.
oldF_struc = [];
oldQ = [];
oldc = [];
oldK = K;
Fremoved = [];
if (K.f>0)
    % reduce if user explicitely says remove, or user says nothing but
    % solverdefinitions does, and there are no nonlinear variables
    if ((options.removeequalities==1 | options.removeequalities==2) & isempty(intersect(used_variables,nonlinearvariables))) | ((options.removeequalities==0) & (solver.constraint.equalities.linear==-1))
        showprogress('Solving equalities',options.showprogress);
        [x_equ,H,A_equ,b_equ,factors] = solveequalities(F_struc,K,options.removeequalities==1);
        % Exit if no consistent solution exist
        if (norm(A_equ*x_equ-b_equ,'inf')>1e-5)%sqrt(eps)*size(A_equ,2))
            diagnostic.solvertime = 0;
            diagnostic.info = yalmiperror(1,'YALMIP');
            diagnostic.problem = 1;