Минимизация с градиентом и гессианом

В этом примере показано, как решить нелинейную задачу минимизации с явной трехдиагональной матрицей Гессиана H(x). Проблема состоит в том, чтобы найти x минимизировать

f(x)=i=1n-1((xi2)(xi+12+1)+(xi+12)(xi2+1)),

где n = 1000.

Функция помощника brownfgh в конце этого примера вычисляет f(x), его градиент g(x), и его Гессиан H(x). Указывать что fminunc использование решателя производная информация, набор SpecifyObjectiveGradient и HessianFcn опции с помощью optimoptions. Использовать Гессиан с fminunc, необходимо использовать 'trust-region' алгоритм.

options = optimoptions(@fminunc,'Algorithm','trust-region',...
    'SpecifyObjectiveGradient',true,'HessianFcn','objective');

Установите параметр n к 1 000, и набор начальная точка xstart к –1 для нечетных компонентов и +1 для даже компонентов.

n = 1000;
xstart = -ones(n,1);
xstart(2:2:n) = 1;

Найдите минимальное значение f.

[x,fval,exitflag,output] = fminunc(@brownfgh,xstart,options);
Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.

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

disp(fval)
   2.8709e-17
disp(exitflag)
     1
disp(output)
         iterations: 7
          funcCount: 8
           stepsize: 0.0039
       cgiterations: 7
      firstorderopt: 4.7948e-10
          algorithm: 'trust-region'
            message: '...'
    constrviolation: []

Функция f(x) сумма степеней квадратов, и, поэтому, является неотрицательной. Решение fval почти нуль, таким образом, это - ясно минимум. Выходной флаг 1 также указывает на тот fminunc находит решение. output структура показывает тот fminunc берет только семь итераций, чтобы достигнуть решения.

Отобразите самые большие и самые маленькие элементы решения.

disp(max(x))
   1.1987e-10
disp(min(x))
  -1.1987e-10

Решение очень около точки где все элементы x = 0.

Функция помощника

Этот код создает brownfgh функция помощника.

function [f,g,H] = brownfgh(x)
%BROWNFGH  Nonlinear minimization problem (function, its gradients
% and Hessian)
% Documentation example        

%   Copyright 1990-2008 The MathWorks, Inc.

% Evaluate the function.
  n=length(x); y=zeros(n,1);
  i=1:(n-1);
  y(i)=(x(i).^2).^(x(i+1).^2+1)+(x(i+1).^2).^(x(i).^2+1);
  f=sum(y);
%
% Evaluate the gradient.
  if nargout > 1
     i=1:(n-1); g = zeros(n,1);
     g(i)= 2*(x(i+1).^2+1).*x(i).*((x(i).^2).^(x(i+1).^2))+...
             2*x(i).*((x(i+1).^2).^(x(i).^2+1)).*log(x(i+1).^2);
     g(i+1)=g(i+1)+...
              2*x(i+1).*((x(i).^2).^(x(i+1).^2+1)).*log(x(i).^2)+...
              2*(x(i).^2+1).*x(i+1).*((x(i+1).^2).^(x(i).^2));
  end
%
% Evaluate the (sparse, symmetric) Hessian matrix
  if nargout > 2
     v=zeros(n,1);
     i=1:(n-1);
     v(i)=2*(x(i+1).^2+1).*((x(i).^2).^(x(i+1).^2))+...
            4*(x(i+1).^2+1).*(x(i+1).^2).*(x(i).^2).*((x(i).^2).^((x(i+1).^2)-1))+...
            2*((x(i+1).^2).^(x(i).^2+1)).*(log(x(i+1).^2));
     v(i)=v(i)+4*(x(i).^2).*((x(i+1).^2).^(x(i).^2+1)).*((log(x(i+1).^2)).^2);
     v(i+1)=v(i+1)+...
              2*(x(i).^2).^(x(i+1).^2+1).*(log(x(i).^2))+...
              4*(x(i+1).^2).*((x(i).^2).^(x(i+1).^2+1)).*((log(x(i).^2)).^2)+...
              2*(x(i).^2+1).*((x(i+1).^2).^(x(i).^2));
     v(i+1)=v(i+1)+4*(x(i).^2+1).*(x(i+1).^2).*(x(i).^2).*((x(i+1).^2).^(x(i).^2-1));
     v0=v;
     v=zeros(n-1,1);
     v(i)=4*x(i+1).*x(i).*((x(i).^2).^(x(i+1).^2))+...
            4*x(i+1).*(x(i+1).^2+1).*x(i).*((x(i).^2).^(x(i+1).^2)).*log(x(i).^2);
     v(i)=v(i)+ 4*x(i+1).*x(i).*((x(i+1).^2).^(x(i).^2)).*log(x(i+1).^2);
     v(i)=v(i)+4*x(i).*((x(i+1).^2).^(x(i).^2)).*x(i+1);
     v1=v;
     i=[(1:n)';(1:(n-1))'];
     j=[(1:n)';(2:n)'];
     s=[v0;2*v1];
     H=sparse(i,j,s,n,n);
     H=(H+H')/2;
  end
end

Похожие темы