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

Используйте в своих интересах структурированный Гессиан

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

В этом примере матрица Гессиана H имеет структуру H = B + A*A', где B является разреженным 512 512 симметрическая матрица, и A является 512 10 разреженной матрицей, состоявшей из многих плотных столбцов. Чтобы избежать чрезмерного использования памяти, которое могло произойти путем работы с H непосредственно, потому что H является плотным, пример обеспечивает, Гессиан умножают функцию, qpbox4mult. Эта функция, когда передано матричный Y, использует разреженные матрицы A и B, чтобы вычислить матричное произведение Гессиана       W = H*Y = (B + A*A')*Y.

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

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

Шаг 1: Решите что часть H передать quadprog в качестве первого аргумента.

Или A или B могут быть переданы в качестве первого аргумента quadprog. Пример принимает решение передать B в качестве первого аргумента, потому что это приводит к лучшему предварительному формирователю (см. Предварительное создание условий).

quadprog(B,f,[],[],[],[],l,u,xstart,options)

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

Теперь, задайте функциональный runqpbox4 это

  • Содержит вложенную функцию qpbox4mult, который использует A и B, чтобы вычислить матричное произведение Гессиана W, где       W = H*Y = (B + A*A')*Y. Вложенная функция должна иметь форму

    W = qpbox4mult(Hinfo,Y,...)

    Первые два аргумента Hinfo и Y требуются.

  • Загружает проблемные параметры от qpbox4.mat.

  • optimoptions использования, чтобы установить опцию HessianMultiplyFcn на указатель на функцию, который указывает на qpbox4mult.

  • Вызовы quadprog с B в качестве первого аргумента.

Первый аргумент к вложенной функции, qpbox4mult должен совпасть с первым аргументом, передал quadprog, который в этом случае является матрицей B.

Вторым аргументом к qpbox4mult является матричный Y (  W = H*Y). Поскольку quadprog ожидает, что Y будет использоваться, чтобы сформировать матричное произведение Гессиана, Y всегда является матрицей со строками n, где n является количеством размерностей в проблеме. Количество столбцов в Y может отличаться. Функциональный qpbox4mult вкладывается так, чтобы значение матричного A прибыло из внешней функции. Программное обеспечение Optimization Toolbox™ включает файл runqpbox4.m.

function [fval, exitflag, output, x] = runqpbox4
%RUNQPBOX4 demonstrates 'HessianMultiplyFcn' option for QUADPROG with bounds.

problem = load('qpbox4'); % Get xstart, u, l, B, A, f
xstart = problem.xstart; u = problem.u; l = problem.l;
B = problem.B; A = problem.A; f = problem.f;
mtxmpy = @qpbox4mult; % function handle to qpbox4mult nested function

% Choose algorithm and the HessianMultiplyFcn option
options = optimoptions(@quadprog,'Algorithm','trust-region-reflective','HessianMultiplyFcn',mtxmpy);

% Pass B to qpbox4mult via the H argument. Also, B will be used in
% computing a preconditioner for PCG.
[x, fval, exitflag, output] = quadprog(B,f,[],[],[],[],l,u,xstart,options);

    function W = qpbox4mult(B,Y)
        %QPBOX4MULT Hessian matrix product with dense structured Hessian.
        %   W = qpbox4mult(B,Y) computes W = (B + A*A')*Y where
        %   INPUT:
        %       B - sparse square matrix (512 by 512)
        %       Y - vector (or matrix) to be multiplied by B + A'*A.
        %   VARIABLES from outer function runqpbox4:
        %       A - sparse matrix with 512 rows and 10 columns.
        %
        %   OUTPUT:
        %       W - The product (B + A*A')*Y.
        %

        % Order multiplies to avoid forming A*A',
        %   which is large and dense
        W = B*Y + A*(A'*Y);
    end

end

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

Чтобы вызвать квадратичную стандартную программу минимизации, содержавшуюся в runqpbox4, войти

[fval,exitflag,output] = runqpbox4;

запускать предыдущий код. Затем отобразите значения для fval, exitflag, output.iterations и output.cgiterations.

fval,exitflag,output.iterations, output.cgiterations

fval =

  -1.0538e+03


exitflag =

     3


ans =

    18


ans =

    30

После 18 итераций с в общей сложности 30 итерациями PCG значение функции уменьшается до

fval
fval =
 -1.0538e+003

и оптимальность первого порядка

output.firstorderopt
ans =
    0.0043

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

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

Поскольку предварительный формирователь - больше аппроксимированное чем тогда, когда H доступен явным образом, настраивание параметра TolPCG к несколько меньшему значению может требоваться. Этот пример совпадает с предыдущим, но уменьшает TolPCG от значения по умолчанию 0.1 к 0,01.

function [fval, exitflag, output, x] = runqpbox4prec
%RUNQPBOX4PREC demonstrates 'HessianMultiplyFcn' option for QUADPROG with bounds.

problem = load('qpbox4'); % Get xstart, u, l, B, A, f
xstart = problem.xstart; u = problem.u; l = problem.l;
B = problem.B; A = problem.A; f = problem.f;
mtxmpy = @qpbox4mult; % function handle to qpbox4mult nested function

% Choose algorithm, the HessianMultiplyFcn option, and override the TolPCG option
options = optimoptions(@quadprog,'Algorithm','trust-region-reflective',...
    'HessianMultiplyFcn',mtxmpy,'TolPCG',0.01);

% Pass B to qpbox4mult via the H argument. Also, B will be used in
% computing a preconditioner for PCG.
% A is passed as an additional argument after 'options'
[x, fval, exitflag, output] = quadprog(B,f,[],[],[],[],l,u,xstart,options);

    function W = qpbox4mult(B,Y)
        %QPBOX4MULT Hessian matrix product with dense structured Hessian.
        %   W = qpbox4mult(B,Y) computes W = (B + A*A')*Y where
        %   INPUT:
        %       B - sparse square matrix (512 by 512)
        %       Y - vector (or matrix) to be multiplied by B + A'*A.
        %   VARIABLES from outer function runqpbox4prec:
        %       A - sparse matrix with 512 rows and 10 columns.
        %
        %   OUTPUT:
        %       W - The product (B + A*A')*Y.
        %

        % Order multiplies to avoid forming A*A',
        %   which is large and dense
        W = B*Y + A*(A'*Y);
    end

end

Теперь, войти

[fval,exitflag,output] = runqpbox4prec; 

запускать предыдущий код. После 18 итераций и 50 итераций PCG, значение функции имеет то же значение к пяти значительным цифрам

fval
fval = 
-1.0538e+003

и оптимальность первого порядка является по существу тем же самым.

output.firstorderopt
ans =
    0.0043

Примечание

Уменьшение TolPCG слишком много может существенно увеличить число итераций PCG.

Похожие темы