exponenta event banner

Минимизация с помощью шаблона градиента и гессенского разрежения

В этом примере показано, как решить нелинейную задачу минимизации с помощью тридиагональной гессенской матрицы, аппроксимированной разреженными конечными разностями, вместо явного вычисления.

Проблема заключается в поиске x для минимизации

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

где n = 1000.

n = 1000;

Для использования trust-region метод в fminunc, необходимо вычислить градиент в целевой функции; он не является необязательным, как в quasi-newton способ.

Вспомогательная функция brownfg в конце этого примера вычисляет целевую функцию и градиент.

Чтобы обеспечить эффективное вычисление разреженной конечной разности аппроксимации гессенской матрицы H (x), структура разреженности H должна быть предопределена. В этом случае структураHstr, разреженная матрица, доступна в файле brownhstr.mat. Использование spy команда, вы можете видеть, что Hstr действительно, является разреженным (всего 2998 ненулевых).

load brownhstr
spy(Hstr)

Figure contains an axes. The axes contains an object of type line.

Установите HessPattern опция для Hstr использование optimoptions. Когда такая большая проблема имеет очевидную структуру разреженности, не устанавливая HessPattern опция использует большой объем памяти и вычисления без необходимости, потому что fminunc пытается использовать конечную разностность для полной гессенской матрицы из одного миллиона ненулевых записей.

Чтобы использовать шаблон разреженности Гессена, необходимо использовать trust-region алгоритм fminunc. Этот алгоритм также требует установки SpecifyObjectiveGradient опция для true использование optimoptions.

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

Задайте для целевой функции значение @brownfg. Установите начальную точку в -1 для нечетных x компонентов и + 1 для четных x компонентов.

xstart = -ones(n,1); 
xstart(2:2:n,1) = 1;
fun = @brownfg;

Решить проблему, позвонив fminunc используя начальную точку xstart и варианты options.

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

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

Проверьте решение и процесс решения.

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

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

Отображение наибольших и наименьших элементов решения.

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

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

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

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

function [f,g] = brownfg(x)
% BROWNFG Nonlinear minimization test problem
% 
% 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
  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
end

Связанные темы