数理最適化セミナーのご案内

5.7 対称行列クラスSymmetricMatrix

　対称行列はSymmetricMatrixというクラスで表現されます．対称行列自体の定義と，対称行列の構造の定義は別々に行う必要があります．例えば，次のような二次正方対称行列を定義したいとします．

\[X = \left( \begin{array}{cc} x+3 & 4y+1.5z \\ 4y+1.5z & 2x+10y \end{array} \right)\]

　この場合，以下のように記述します．

Set S = "1 2";
Element i(set = S), j(set = S);
Variable x, y, z;
SymmetricMatrix X((i, j));
X["1, 1"] = x + 3;
X["1, 2"] = 4 * y + 1.5 * z;
X["2, 2"] = 2 * x + 10 * y;

　対角要素以外の成分の定義は上下三角部分いずれか一方に関してのみ定義してください．上記の例では[1,2]成分が定義されているので，[2,1]成分を定義する必要はありません．

　対称成分が重複して定義された場合は，先に定義された方は無視されます．例えば，

X["1, 1"] = x + 3;
X["1, 2"] = 1.5 * y + 4 * z; // 次の[2,1]要素の定義で打ち消される
X["2, 1"] = 4 * y + 1.5 * z;
X["2, 2"] = 2 * x + y;

は，次の行列を定義している事になります．

\[X = \left( \begin{array}{cc} x+3 & 4y+1.5z \\ 4y+1.5z & 2x+y \end{array} \right)\]

　Variable x,y,zに添字を付けて，係数に対してもParameterを導入すれば，行列の定義を一括して行う事もできます．以下の例をご覧ください． 

Set S = "1 2";
Set T = "1 2 3";
Element i(set = S), j(set = S);
Element k(set = T);
Variable x(index = k);
Parameter a(index = (k, i, j));
Parameter b(index = (i, j));
SymmetricMatrix X((i, j));
// 上三角部分を定義
a["1, 1, 1"] = 1;
a["1, 2, 2"] = 2;
a["2, 1, 2"] = 4;
a["2, 2, 2"] = 1.5;
a["3, 1, 2"] = 10;
b["1, 1"] = 3;
X[i, j] = sum(a[k, i, j] * x[k], k) + b[i, j], i <= j;

　この例では，行列$X = \left( \begin{array}{cc} x_{1} + 3 & 4x_{2} + 10x_{3} \\ 4x_{2} + 10x_{3} & 2x_{1} + 1.5x_{2} \end{array} \right)$を定義している事になります．

　対称行列に対しては，制約として半正定値制約を課す事ができます．半正定値制約を課す場合，右辺に>= 0を記述する必要があります．次の例では対称行列Xの半正定値制約を記述しています． 

SymmetricMatrix X((i, j));
X >= 0;

　右辺には0以外のスカラー値を用いることもできます．この場合，右辺値は左辺行列の最小固有値を意味します．例えば，次の例は行列Xの最小固有値が2つまり，$X - 2E \succeq 0$であることを意味します． 

X >= 2;

　次の記述は誤りです．

SymmetricMatrix X((i, j));
SymmetricMatrix Y((i, j));
X >= Y;

　複数の対称行列を一括して定義する事もできます．例えば，次の例では対称行列$X_{1}, .., X_{10}$を一括して定義しています． 

Set S = "1 2";
Element i(set = S), j(set = S);
Set N = "1 .. 10";
Element n(set = N);
SymmetricMatrix X(index = n, (i, j));

　対称行列自体の添字nにはindex =を付ける必要があります．個別の対称行列内部の添字(i,j)にはindex =を付与してはいけません．

　複数の対称行列に対して一気に半正定値制約を設定するには，次のように記述します．

SymmetricMatrix X(index = n, (i, j));
X[n] >= 0;

　大規模な行列を考える場合，定数並びに行列成分を疎形式で定義することで高速な求解が可能です．次の例は少し難解ですが，Parameter a, bが出現する添字に絞って行列成分を定義しています． 

Set S, T;
Element i(set = S), j(set = S);
Element k(set = T);
Variable x(index = k);
Set A(dim = 3, superSet = (T, S, S));
Element m(set = A);
Set B(dim = 2, superSet = (S, S));
Element n(set = B);
Parameter a(name = "a", index = m);
Parameter b(name = "b", index = n);
A = A | "1" * B;
B = A.slice(2, 3);
S = A.slice(1);
T = A.slice(2) | A.slice(3);
SymmetricMatrix X((i, j));
X[i, j] = sum(a[k, i, j] * x[k], (k, (k, i, j) < A)) + b[i, j], (i, j) < B;

　上記の例ではParameter a, bの値を外部ファイルから与えるケースを想定しています．外部ファイルから自動代入によって定まるa, bの添字の範囲がそれぞれ集合A, Bであると定められています．行列内部での疎形式定義の為に，集合A, Bをそれぞれ加工します．また，変数の範囲を示す集合T, 行列内の添字の範囲を示す集合Sは，いずれも集合Aから作成されています．

　具体的に以下のようなデータファイルが与えられたとします．

a =
[1, 1, 1] 1
[1, 2, 3] 3
[2, 2, 2] 2
[3, 2, 2] 0.5
[3, 3, 3] 5
;
b =
[1, 1] 10
[1, 2] 4
;

　これは，次の対称行列を定義している事に相当します．

\[X = \left( \begin{array}{ccc} x_{1} + 10 & 4 & \\ 4 & 2x_{2} + 0.5x_{3} & 3x_{1} \\ & 3x_{1} & 5x_{3} \end{array} \right)\]

　自動代入の段階では，集合Aの要素は

\[\{(1,1,1), (1,2,3), (2,2,2), (3,2,2), (3,3,3)\}\]

　集合Bの要素は

\[\{(1,1), (1,2)\}\]

です．これに対して以下の加工を施す事により，

A = A | "1" * B;
B = A.slice(2, 3);

　集合Aの要素は

\[\{(1,1,1), (1,1,2), (1,2,3), (2,2,2), (3,2,2), (3,3,3)\}\]

　集合Bの要素は

\[\{(1,1), (1,2), (2,2), (2,3), (3,3)\}\]

となります．なお，"1" * Bでは，1個の要素（1）からなる集合と集合Bの直積を求めています．また，A | "1" * Bにより集合Aと（直積）集合"1" * Bの和集合を求めています．

　関数slice()は，集合の一部を射影して切り出す機能を有しています．なお，関数slice()については引数を最大5個とることができます．

　上記の例では，集合Aの第二第三成分を切り出して，集合Bに渡しています．この結果|A| = 6, |B| = 5となりましたが，特に何の工夫もしない場合|A| = 27, |B| = 9となり，内部で余分な領域を必要とします．特に行列の次元が大きい場合には，パフォーマンスにかなりの違いが生じます．

5.6 範囲演算関数sum, prod
5.8 式クラスExpression

 

 

上に戻る