構造化されたヘッシアンの利用

quadprog信頼領域 反思的法では、ヘッシアンが密であっても構造化されていれば、大規模な問題を解くことができます。これらの問題に対しquadprogは、スパース行列 Hを使った信頼領域 反思的法の問題で行うように、ヘッシアン Hを使用して直接 H*Yを計算するという方法は取りません。Hの形成には多大なメモリを要するからです。その代わり、行列 Yと Hの情報をもとに W=H*Yを計算する関数をquadprogに与えなければなりません。

この例ではヘッセ行列HはH=B+A*A'という構造をしています。ここで、Bは 512行 512列の対称スパース行列でA.は密な列から構成される512行10列のスパース行列です。Hが密であるため、直接Hを使うことで起こる過度のメモリ使用を避けるために、この例ではヘッセ乗算関数qpbox4multを与えます。この関数は行列Yが渡されると、スパース行列A.とBを使用してヘッセ行列の積W=H*Y=（B+A*A'）*Yを計算します。

この例の最初の部分では,ヘッセ乗算関数qpbox4multに、行列A.とBを渡す必要があります。quadprogの第 1.引数に、ヘッセ乗算関数に渡す行列を入力できます。第 2.の行列の値を与えるために入れ子関数を使用できます。

この例の 2.番目の部分では、正確な行列Hではなく近似の前提条件子を使用することによる影響を補正するために托尔普克許容誤差を厳しくする方法を説明します。

手順1:1番目の引数としてquadprogに渡すHの一部を決定する

A.またはBをquadprogの 1.番目の引数として渡すことができます。例では良い前提条件子 (前処理を参照) の結果から 1.番目の引数としてBを渡すことにします。

quadprog（B、f、[]、[]、[]、[]、[]、l、u、xstart、选项）

手順2:Hに対するヘッセ行列の乗算を計算する関数を記述する

ここで、次のような関数runqpbox4を定義します。

入れ子関数qpbox4multの第 1.引数はquadprogに渡される第 1.引数と同じである必要があります。この場合は行列 Bです。

qpbox4multに対する第 2.引数は行列Y(W=H*Yのもの) です。quadprogではヘッセ行列の積の構成に使われるYを想定しているため、問題の次元数をNとした場合、Yは常にN行の行列になります。Yの列数は可変です。関数qpbox4multは、行列A.の値が外側の関数から得られる入れ子になります。优化工具箱™ ソフトウェアにはrunqpbox4.mファイルが付属しています。

function [fval, exitflag, output, x] = runqpbox4 % runqpbox4演示了QUADPROG的'HessianMultiplyFcn'选项。问题=负载(“qpbox4”);% Get xstart, u, l, B, A, f xstart = problem.xstart;u = problem.u;l = problem.l;B = problem.B;一个=实例计算;f = problem.f;mtxmpy = @qpbox4mult;% select algorithm and the HessianMultiplyFcn option options = optimoptions(@quadprog，' algorithm '，'信任区域-反射'，'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

手順 三：開始点を使った二次の最小化ルーチンを呼び出す

runqpbox4内の二次の最小化ルーチンを呼び出すには次のように入力してください。

[fval，exitflag，output]=runqpbox4；

上記のコードを実行します。次に,未来值、出口滞后、输出迭代次数およびoutput.cg迭代次数の値を表示します。

fval，exitflag，output.iterations，output.cgtations fval=-1.0538e+03 exitflag=3 ans=18 ans=30

30の PCG反復による 18回の反復の後、関数値が減ります。

fval fval=-1.0538e+003

そして、1.次の最適性は、以下のようになります。

output.firstorderopt ans=0.0043

前処理

Hが陰的に存在するためquadprogが前提条件子の計算にHを使用できないことがあります。その代わり、quadprogはHの代わりに渡す引数として前提条件子の計算にBを使用します。BはHとサイズが同じであり、またHとある程度近似しているため適切な選択となります。BがHと同じサイズでない場合,quadprogはアルゴリズムで決められた対角行列をベースに前提条件子を計算します。通常、これはうまく実行されません。

Hが陽的に使用可能な場合、前提条件子はより良い近似になるため、托尔普克パラメーターの調整はやや小さめの値にしなければなりません。この例は前の例と同じですが、托尔普克は既定の 0.1から 0.01に減らしています。

作用[fval，exitflag，output，x]=runqpbox4prec%RUNQPBOX4PREC演示了带边界的QUADPROG的“HessianMultiplyFcn”选项。问题=负载(“qpbox4”);u, l, B, A, fxstart=problem.xstart；u=problem.u；l=problem.l；B=problem.B；A=problem.A；f=problem.f；mtxmpy=@qpbox4mult；%qpbox4mult嵌套函数的函数句柄%选择algorithm、HessianMultiplyFcn选项，并覆盖TolPCG选项options=options（@quadprog，“算法”,“trust-region-reflective”,...“HessianMultiplyFcn”，mtxmpy，“托尔普克”,0.01);%通过H参数将B传递给qpbox4mult。此外，B将在%计算PCG的预处理器。%A作为“选项”后面的附加参数传递[x，fval，exitflag，output]=quadprog（B，f，[]，[]，[]，[]，[]，[]，[]，l，u，xstart，options）；作用W = qpbox4mult (B, Y)%QPBOX4MULT Hessian矩阵积，具有稠密结构Hessian。% W = qpbox4mult(B,Y)计算W = (B + A*A')*Y%输入：%B-稀疏平方矩阵（512 x 512）%Y-向量（或矩阵）乘以B+A'*A。%外部函数runqpbox4prec中的变量：% A -稀疏矩阵，512行10列。%%输出：%W-产品（B+A*A'）*Y。%%为了避免形成*A'，对倍数进行排序，%它又大又密W = b * y + a *(a '* y);终止终止

ここで以下を入力し,

[fval，exitflag，output]=runqpbox4prec；

上記のコードを実行します。18回の反復と 50回の PCG反復の後、関数値は 5.桁の有効数字となる同じ値になります。

fval fval=-1.0538e+003

さらに、1.次の最適性は基本的には同じです。

output.firstorderopt ans=0.0043