Минимизация с плотным структурированным гессианом, линейными равенствами

Гессиан умножает функцию для более низкой памяти

interior-point fmincon и алгоритмы trust-region-reflective и алгоритм trust-region fminunc могут решить проблемы, где Гессиан является плотным, но структурирован. Для этих проблем fmincon и fminunc не вычисляют H*Y с Гессианом H непосредственно, потому что формирование H было бы интенсивно использующим память. Вместо этого необходимо предоставить fmincon или fminunc с функцией, которая, учитывая матричный Y и информацию о H, вычисляет W = H*Y.

В этом примере целевая функция нелинейна, и линейные равенства существуют так, fmincon используется. Описание применяется к доверительной области отражающий алгоритм; алгоритм trust-region fminunc подобен. Для алгоритма внутренней точки см., что опция HessianMultiplyFcn в Гессиане Умножает Функцию. Целевая функция имеет структуру

f(x)=f^(x)12xTVVTx,

где V 1000 2 матрица. Гессиан f является плотным, но Гессиан f^ issparse . Если Гессиан H^, затем H, Гессиан f,

H=H^VVT.

Чтобы избежать чрезмерного использования памяти, которое могло произойти путем работы с H непосредственно, пример обеспечивает, Гессиан умножают функцию, hmfleq1. Эта функция, когда передано матричный Y, использует разреженные матрицы Hinfo, который соответствует, и V, чтобы вычислить матричное произведение Гессиана

W = H*Y = (Hinfo - V*V')*Y

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

Однако не константа и должна быть вычислена в текущем x. Можно сделать это путем вычисления в целевой функции и возврата как Hinfo в третьем выходном аргументе. При помощи optimoptions, чтобы установить опции 'Hessian' на 'on', fmincon знает, чтобы получить значение Hinfo от целевой функции, и передать его Гессиану умножают функциональный hmfleq1.

Шаг 1: Запишите файл brownvv.m, который вычисляет целевую функцию, градиент и разреженную часть Гессиана.

Пример передает brownvv fmincon как целевая функция. Файл brownvv.m длинен и не включен здесь. Можно просмотреть код с командой

type brownvv

Поскольку brownvv вычисляет градиент, а также целевую функцию, пример (Шаг 3) использует optimoptions, чтобы установить опцию SpecifyObjectiveGradient на true.

Шаг 2: Запишите функцию, чтобы вычислить матричные произведения Гессиана для H, учитывая матрицу Y.

Теперь, задайте функциональный hmfleq1, который использует Hinfo, который вычисляется в brownvv и V, который можно получить в указателе на функцию к анонимной функции, чтобы вычислить матричное произведение Гессиана W где       W = H*Y = (Hinfo - V*V')*Y. Эта функция должна иметь форму

W = hmfleq1(Hinfo,Y)

Первый аргумент должен совпасть с третьим аргументом, возвращенным целевой функцией brownvv. Второй аргумент к Гессиану умножается, функцией является матричный Y (W = H*Y).

Поскольку fmincon ожидает, что второй аргумент Y будет использоваться, чтобы сформировать матричное произведение Гессиана, Y всегда является матрицей со строками n, где n является количеством размерностей в проблеме. Количество столбцов в Y может отличаться. Наконец, можно использовать указатель на функцию для анонимной функции, чтобы получить V, таким образом, V может быть третий аргумент к 'hmfleqq'.

function W = hmfleq1(Hinfo,Y,V);
%HMFLEQ1 Hessian-matrix product function for BROWNVV objective.
%   W = hmfleq1(Hinfo,Y,V) computes W = (Hinfo-V*V')*Y
%   where Hinfo is a sparse matrix computed by BROWNVV 
%   and V is a 2 column matrix.
W = Hinfo*Y - V*(V'*Y);

Примечание

Функциональный hmfleq1 доступен в папке optimdemos как hmfleq1.m.

Шаг 3: Вызовите нелинейную стандартную программу минимизации с отправной точкой и линейными ограничениями равенства.

Загрузите проблемный параметр, V, и разреженные ограничительные матрицы равенства, Aeq и beq, от fleq1.mat, который доступен в папке optimdemos. Используйте optimoptions, чтобы установить опцию SpecifyObjectiveGradient на true и установить опцию HessianMultiplyFcn на указатель на функцию, который указывает на hmfleq1. Вызовите fmincon с целевой функцией brownvv и с V как дополнительный параметр:

function [fval,exitflag,output,x] = runfleq1
% RUNFLEQ1 demonstrates 'HessMult' option for FMINCON with linear
% equalities.

problem = load('fleq1');   % Get V, Aeq, beq
V = problem.V; Aeq = problem.Aeq; beq = problem.beq;
n = 1000;             % problem dimension
xstart = -ones(n,1); xstart(2:2:n,1) = ones(length(2:2:n),1); % starting point
options = optimoptions(@fmincon,...
    'Algorithm','trust-region-reflective',...
    'SpecifyObjectiveGradient',true, ...
    'HessianMultiplyFcn',@(Hinfo,Y)hmfleq1(Hinfo,Y,V),...
    'Display','iter',...
    'OptimalityTolerance',1e-9,...
    'FunctionTolerance',1e-9);
[x,fval,exitflag,output] = fmincon(@(x)brownvv(x,V),xstart,[],[],Aeq,beq,[],[], ...
                                    [],options);

Чтобы запустить предыдущий код, войти

[fval,exitflag,output,x] = runfleq1;

Поскольку итеративное отображение было установлено с помощью optimoptions, эта команда генерирует следующее итеративное отображение:

                                Norm of      First-order 
 Iteration        f(x)          step          optimality   CG-iterations
     0            2297.63                      1.41e+03                
     1            1084.59         6.3903            578           1
     2            1084.59            100            578           3
     3            1084.59             25            578           0
     4            1084.59           6.25            578           0
     5            1047.61         1.5625            240           0
     6            761.592          3.125           62.4           2
     7            761.592           6.25           62.4           4
     8            746.478         1.5625            163           0
     9            546.578          3.125           84.1           2
    10            274.311           6.25           26.9           2
    11            55.6193        11.6597             40           2
    12            55.6193             25             40           3
    13            22.2964           6.25           26.3           0
    14            -49.516           6.25             78           1
    15           -93.2772         1.5625             68           1
    16           -207.204          3.125           86.5           1
    17           -434.162           6.25           70.7           1
    18           -681.359           6.25           43.7           2
    19           -681.359           6.25           43.7           4
    20           -698.041         1.5625            191           0
    21           -723.959          3.125            256           7
    22            -751.33        0.78125            154           3
    23           -793.974         1.5625           24.4           3
    24           -820.831        2.51937           6.11           3
    25           -823.069       0.562132           2.87           3
    26           -823.237       0.196753          0.486           3
    27           -823.245      0.0621202          0.386           3
    28           -823.246      0.0199951           0.11           6
    29           -823.246     0.00731333         0.0404           7
    30           -823.246     0.00505883         0.0185           8
    31           -823.246     0.00126471        0.00268           9
    32           -823.246     0.00149326        0.00521           9
    33           -823.246    0.000373314        0.00091           9

Local minimum possible.

fmincon stopped because the final change in function value relative to 
its initial value is less than the value of the function tolerance.

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

problem = load('fleq1');   % Get V, Aeq, beq
V = problem.V; Aeq = problem.Aeq; beq = problem.beq;
norm(Aeq*x-beq,inf)

ans =
    1.8874e-14

Предварительное создание условий

В этом примере fmincon не может использовать H, чтобы вычислить предварительный формирователь, потому что H только существует неявно. Вместо H fmincon использует Hinfo, третий аргумент, возвращенный brownvv, чтобы вычислить предварительный формирователь. Hinfo является хорошим выбором, потому что он одного размера как H и аппроксимирует H до некоторой степени. Если бы Hinfo не был тем же размером как H, fmincon вычислил бы предварительный формирователь на основе некоторых диагональных матриц масштабирования, определенных из алгоритма. Как правило, это не выполнило бы также.