Этот пример показывает, как решить нелинейную проблему минимизации с трехдиагональной матрицей Гессиана, аппроксимированной разреженными конечными разностями вместо явного вычисления.
Проблема состоит в том, чтобы найти, что x минимизирует
где n = 1000.
Чтобы использовать метод trust-region в fminunc, необходимо вычислить градиент в fun; это не является дополнительным как в методе quasi-newton.
Файл 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Чтобы позволить эффективное вычисление разреженного приближения конечной разности матрицы Гессиана H (x), структура разреженности H должна быть предопределена. В этом случае примите эту структуру, Hstr, разреженная матрица, доступен в файле brownhstr.mat. Используя команду spy вы видите, что Hstr действительно разрежен (только 2 998 ненулей). Используйте optimoptions, чтобы установить опцию HessPattern на Hstr. Когда проблема, столь большая, как это имеет очевидную структуру разреженности, не устанавливая опцию HessPattern, требует огромной суммы ненужной памяти и вычисления, потому что fminunc пытается использовать конечное дифференцирование на полной матрице Гессиана одного миллиона ненулевых записей.
Необходимо также установить опцию SpecifyObjectiveGradient на true с помощью optimoptions, поскольку градиент вычисляется в brownfg.m. Затем выполните fminunc как показано на Шаге 2.
fun = @brownfg; load brownhstr % Get Hstr, structure of the Hessian spy(Hstr) % View the sparsity structure of Hstr

n = 1000;
xstart = -ones(n,1);
xstart(2:2:n,1) = 1;
options = optimoptions(@fminunc,'Algorithm','trust-region',...
'SpecifyObjectiveGradient',true,'HessPattern',Hstr);
[x,fval,exitflag,output] = fminunc(fun,xstart,options); Эта проблема с 1000 переменными решена в семи итерациях и семи итерациях метода сопряженных градиентов с положительной сходимостью указания exitflag. Итоговое значение функции и мера оптимальности в решении, которое x оба близко к нулю (для fminunc, оптимальность первого порядка является нормой бесконечности градиента функции, которая является нулем в локальном минимуме):
exitflag,fval,output.firstorderopt
exitflag =
1
fval =
7.4738e-17
ans =
7.9822e-10