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

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

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

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

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

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

H=H^VVT.

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

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

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

Однако H^ не константа и должна быть вычислена в текущем x. Можно сделать это путем вычисления H^ в целевой функции и возврате H^ как 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 вычислил бы предварительный формирователь на основе некоторых диагональных матриц масштабирования, определенных из алгоритма. Как правило, это не выполнило бы также.