Создайте функцию Гессиана

Использование информации о производной второго порядка в fmincon решатель, вы должны создать Гессиана, который является Гессианом Лагранжа. Гессиан Лагранжа задается уравнением

Вот, $$f(x)$$ является bigtoleft функции, и $$c_i(x)$$ являются двумя функциями ограничения конуса. The hessinterior вспомогательная функция в конце этого примера вычисляет Гессиана Лагранжа в точке x со структурой множителя Лагранжа lambda. Функция сначала вычисляет $${\nabla }^{2}f(x)$$. Затем он вычисляет два ограничительных Гессиана $${\nabla }^2{c}_1(x)$$ и $${\nabla }^2{c}_2(x)$$, умножает их на соответствующие им множители Лагранжа lambda.ineqnonlin(1) и lambda.ineqnonlin(2), и добавляет их. Вы можете увидеть из определения twocone функция ограничения, которая $${\nabla }^2{c}_1(x)={\nabla }^2{c}_2(x)$$, что упрощает вычисление.

Создайте опции для использования производных

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

options = optimoptions('fmincon','Algorithm','interior-point',...
    "SpecifyConstraintGradient",true,"SpecifyObjectiveGradient",true,...
    'HessianFcn',@hessinterior);

Минимизируйте использование всей производной информации

Установите начальную точку x0 = [-1,-1,-1].

x0 = [-1,-1,-1];

Задача не имеет линейных ограничений или границ. Установите эти аргументы равными [].

A = [];
b = [];
Aeq = [];
beq = [];
lb = [];
ub = [];

Решите проблему.

[x,fval,eflag,output] = fmincon(@bigtoleft,x0,...
           A,b,Aeq,beq,lb,ub,@twocone,options);

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

Исследуйте процесс решения и решения

Исследуйте решение, значение целевой функции, выходной флаг и количество вычислений функции и итераций.

disp(x)

   -6.5000   -0.0000   -3.5000

disp(fval)

  -2.8941e+03

disp(eflag)

     1

disp([output.funcCount,output.iterations])

     7     6

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

options.HessianFcn = [];
[x2,fval2,eflag2,output2] = fmincon(@bigtoleft,x0,...
           A,b,Aeq,beq,lb,ub,@twocone,options);

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

disp([output2.funcCount,output2.iterations])

    13     9

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

options.SpecifyConstraintGradient = false;
options.SpecifyObjectiveGradient = false;
[x3,fval3,eflag3,output3] = fmincon(@bigtoleft,x0,...
           A,b,Aeq,beq,lb,ub,@twocone,options);

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

disp([output3.funcCount,output3.iterations])

    43     9

Вспомогательные функции

Этот код создает bigtoleft вспомогательная функция.

function [f gradf] = bigtoleft(x)
% This is a simple function that grows rapidly negative
% as x(1) becomes negative
%
f = 10*x(:,1).^3+x(:,1).*x(:,2).^2+x(:,3).*(x(:,1).^2+x(:,2).^2);

if nargout > 1

   gradf=[30*x(1)^2+x(2)^2+2*x(3)*x(1);
       2*x(1)*x(2)+2*x(3)*x(2);
       (x(1)^2+x(2)^2)];

end
end

Этот код создает twocone вспомогательная функция.

function [c ceq gradc gradceq] = twocone(x)
% This constraint is two cones, z > -10 + r
% and z < 3 - r

ceq = [];
r = sqrt(x(1)^2 + x(2)^2);
c = [-10+r-x(3);
    x(3)-3+r];

if nargout > 2

    gradceq = [];
    gradc = [x(1)/r,x(1)/r;
       x(2)/r,x(2)/r;
       -1,1];

end
end

Этот код создает hessinterior вспомогательная функция.

function h = hessinterior(x,lambda)

h = [60*x(1)+2*x(3),2*x(2),2*x(1);
    2*x(2),2*(x(1)+x(3)),2*x(2);
    2*x(1),2*x(2),0];% Hessian of f
r = sqrt(x(1)^2+x(2)^2);% radius
rinv3 = 1/r^3;
hessc = [(x(2))^2*rinv3,-x(1)*x(2)*rinv3,0;
    -x(1)*x(2)*rinv3,x(1)^2*rinv3,0;
    0,0,0];% Hessian of both c(1) and c(2)
h = h + lambda.ineqnonlin(1)*hessc + lambda.ineqnonlin(2)*hessc;
end