Линейная или квадратичная цель с квадратичными ограничениями

В этом примере показано, как решить задачу оптимизации, которая имеет линейные или квадратичные объективные и квадратичные ограничения неравенства. Это показывает, как сгенерировать и использовать градиент и Гессиан ограничительных функций и цели.

Квадратичная ограниченная проблема

Предположим, что можно поместить проблему в форму

minx(12xTQx+fTx+c)

при ограничениях

12xTHix+kiTx+di0,

где 1 ≤ i ≤ m. Примите в наименьшем количестве, один Hi является ненулевым; в противном случае можно использовать quadprog или linprog решать эту задачу. С ненулевым Hi ограничения нелинейны, и Таблица решений Оптимизации утверждает тот fmincon соответствующий решатель.

Пример принимает, что квадратичные матрицы симметричны. Это без потери общности; можно заменить несимметричный H (или Q) матрица с эквивалентной симметрированной версией (H + HT)/2.

Если x имеет компоненты N, то Q и Hi является N-by-N матрицы, f и ki является N-by-1 векторы, и c и di являются скалярами.

Целевая функция

Сформулируйте проблему с помощью fmincon синтаксис. Примите тот x и f вектор-столбцы. X вектор-столбец когда начальный векторный x0 .)

function [y,grady] = quadobj(x,Q,f,c)
y = 1/2*x'*Q*x + f'*x + c;
if nargout > 1
    grady = Q*x + f;
end

Ограничительная функция

Для непротиворечивости и легкой индексации, поместите каждую квадратичную матрицу ограничений в один массив ячеек. Точно так же поместите линейные и постоянные условия в массивы ячеек.

function [y,yeq,grady,gradyeq] = quadconstr(x,H,k,d)
jj = length(H); % jj is the number of inequality constraints
y = zeros(1,jj);
for i = 1:jj
    y(i) = 1/2*x'*H{i}*x + k{i}'*x + d{i};
end
yeq = [];
    
if nargout > 2
    grady = zeros(length(x),jj);
    for i = 1:jj
        grady(:,i) = H{i}*x + k{i};
    end
end
gradyeq = [];

Числовой пример

Например, предположите, что у вас есть следующая проблема.

Q = [3,2,1;
     2,4,0;
     1,0,5];
f = [-24;-48;-130];
c = -2;

rng default % for reproducibility
% Two sets of random quadratic constraints:
H{1} = gallery('randcorr',3); % random positive definite matrix
H{2} = gallery('randcorr',3);
k{1} = randn(3,1);
k{2} = randn(3,1);
d{1} = randn;
d{2} = randn;

Гессиан

Создайте функцию Гессиана. Гессиан функции Лагранжа дан уравнением

xx2L(x,λ)=2f(x)+λi2ci(x)+λi2ceqi(x).

fmincon вычисляет аппроксимированный набор множителей Лагранжа λi и группирует их в структуре. Чтобы включать Гессиан, используйте следующую функцию.

function hess = quadhess(x,lambda,Q,H)
hess = Q;
jj = length(H); % jj is the number of inequality constraints
for i = 1:jj
    hess = hess + lambda.ineqnonlin(i)*H{i};
end

Решение

Используйте fmincon interior-point алгоритм, чтобы решить задачу наиболее эффективно. Этот алгоритм принимает функцию Гессиана, которую вы предоставляете. Установите эти опции.

options = optimoptions(@fmincon,'Algorithm','interior-point',...
    'SpecifyObjectiveGradient',true,'SpecifyConstraintGradient',true,...
    'HessianFcn',@(x,lambda)quadhess(x,lambda,Q,H));

Вызовите fmincon решать задачу.

fun = @(x)quadobj(x,Q,f,c);
nonlconstr = @(x)quadconstr(x,H,k,d);
x0 = [0;0;0]; % column vector
[x,fval,eflag,output,lambda] = fmincon(fun,x0,...
    [],[],[],[],[],[],nonlconstr,options);

Исследуйте множители Лагранжа.

lambda.ineqnonlin
ans =

   12.8412
   39.2337

Оба нелинейных множителя неравенства являются ненулевыми, таким образом, оба квадратичных ограничения активны в решении.

КПД, когда обеспечение гессиана

Алгоритм внутренней точки с градиентами и Гессианом эффективен. Исследуйте количество вычислений функции.

output
output = 

         iterations: 9
          funcCount: 10
    constrviolation: 0
           stepsize: 5.3547e-04
          algorithm: 'interior-point'
      firstorderopt: 1.5851e-05
       cgiterations: 0
            message: 'Local minimum found that satisfies the constraints.

Optimization compl...'

fmincon используемый всего 10 вычислений функции, чтобы решить задачу.

Сравните это с решением без Гессиана.

options.HessianFcn = [];
[x2,fval2,eflag2,output2,lambda2] = fmincon(fun,[0;0;0],...
    [],[],[],[],[],[],nonlconstr,options);
output2
output2 = 

         iterations: 17
          funcCount: 22
    constrviolation: 0
           stepsize: 2.8475e-04
          algorithm: 'interior-point'
      firstorderopt: 1.7680e-05
       cgiterations: 0
            message: 'Local minimum found that satisfies the constraints.

Optimization compl...'

На этот раз fmincon используемый о вдвое большем количестве итераций и вычислений функции. Решения являются тем же самым к в допусках.

Расширение квадратичных ограничений равенства

Если у вас также есть квадратичные ограничения равенства, можно использовать по существу тот же метод. Проблемой является то же самое с дополнительными ограничениями

12xTJix+piTx+qi=0.

Переформулируйте свои ограничения, чтобы использовать Ji, pi и переменные qi. lambda.eqnonlin(i) структура имеет множители Лагранжа для ограничений равенства.

Связанные примеры

Больше о