qr
シンボリック行列の QR 分解
構文
説明
___ = qr(___,'real')
は、入力引数と中間結果が実数であることを仮定するため、abs
および conj
への呼び出しを無効にします。このフラグを使用すると、qr
はすべてのシンボリック変数が実数を表すことを仮定します。このフラグを使用する際は、すべての数値引数が実数であることを確認してください。
'real'
を使用して、結果に複素共役が含まれないようにします。
例
QR 分解の R 部
4
行 4
列のウィルキンソンの固有値テスト行列 QR 分解の R 部を計算します。
4
行 4
列のウィルキンソンの固有値テスト行列を作成します。
A = sym(wilkinson(4))
A = [ 3/2, 1, 0, 0] [ 1, 1/2, 1, 0] [ 0, 1, 1/2, 1] [ 0, 0, 1, 3/2]
単一の出力引数をもつ構文を使用して、QR 分解の Q 部を返さずに R 部を返します。
R = qr(A)
R = [ 13^(1/2)/2, (4*13^(1/2))/13, (2*13^(1/2))/13, 0] [ 0, (13^(1/2)*53^(1/2))/26, (10*13^(1/2)*53^(1/2))/689, (2*13^(1/2)*53^(1/2))/53] [ 0, 0, (53^(1/2)*381^(1/2))/106, (172*53^(1/2)*381^(1/2))/20193] [ 0, 0, 0, (35*381^(1/2))/762]
パスカル行列の QR 分解
3
行 3
列のパスカル行列の QR 分解を計算します。
3
行 3
列のパスカル行列を作成します。
A = sym(pascal(3))
A = [ 1, 1, 1] [ 1, 2, 3] [ 1, 3, 6]
A
の QR 分解を表現する Q
行列と R
行列を計算します。
[Q,R] = qr(A)
Q = [ 3^(1/2)/3, -2^(1/2)/2, 6^(1/2)/6] [ 3^(1/2)/3, 0, -6^(1/2)/3] [ 3^(1/2)/3, 2^(1/2)/2, 6^(1/2)/6] R = [ 3^(1/2), 2*3^(1/2), (10*3^(1/2))/3] [ 0, 2^(1/2), (5*2^(1/2))/2] [ 0, 0, 6^(1/2)/6]
isAlways
を使用して A = Q*R
であることを検証します。
isAlways(A == Q*R)
ans = 3×3 logical array 1 1 1 1 1 1 1 1 1
置換情報
置換を使用すると、浮動小数点行列の QR 分解の数値安定性を上げることができます。関数 qr
は、置換情報を行列またはベクトルとして返します。
可変精度の演算を行うための有効小数桁数を 10 に設定します。3
行 3
列のシンボリック ヒルベルト行列を浮動小数点数で近似します。
previoussetting = digits(10); A = vpa(hilb(3))
A = [ 1.0, 0.5, 0.3333333333] [ 0.5, 0.3333333333, 0.25] [ 0.3333333333, 0.25, 0.2]
まず、A
の QR 分解を置換なしで計算します。
[Q,R] = qr(A)
Q = [ 0.8571428571, -0.5016049166, 0.1170411472] [ 0.4285714286, 0.5684855721, -0.7022468832] [ 0.2857142857, 0.6520863915, 0.7022468832] R = [ 1.166666667, 0.6428571429, 0.45] [ 0, 0.1017143303, 0.1053370325] [ 0, 0, 0.003901371573]
A
と Q*R
の差を計算します。丸め誤差のため、計算された行列 Q
と R
は、厳密に等式 A*P = Q*R
を満たしません。
A - Q*R
ans = [ -1.387778781e-16, -3.989863995e-16, -2.064320936e-16] [ -3.469446952e-18, -8.847089727e-17, -1.084202172e-16] [ -2.602085214e-18, -6.591949209e-17, -6.678685383e-17]
QR 分解の数値安定性を上げるには、3 つの出力引数を使用する構文を指定することで置換を使用します。シンボリック変数、式または関数を含まない行列では、この構文は、R
に返された行列の abs(diag(R))
が減少するようにピボットを起動します。
[Q,R,P] = qr(A)
Q = [ 0.8571428571, -0.4969293466, -0.1355261854] [ 0.4285714286, 0.5421047417, 0.7228063223] [ 0.2857142857, 0.6776309272, -0.6776309272] R = [ 1.166666667, 0.45, 0.6428571429] [ 0, 0.1054092553, 0.1016446391] [ 0, 0, 0.003764616262] P = 1 0 0 0 0 1 0 1 0
再び、等式 A*P = Q*R
を確認します。置換を使用すると QR 分解の丸め誤差が低くなります。
A*P - Q*R
ans = [ -3.469446952e-18, -4.33680869e-18, -6.938893904e-18] [ 0, -8.67361738e-19, -1.734723476e-18] [ 0, -4.33680869e-19, -1.734723476e-18]
次に、'vector'
引数を使用して置換情報をベクトルとして返します。
[Q,R,p] = qr(A,'vector')
Q = [ 0.8571428571, -0.4969293466, -0.1355261854] [ 0.4285714286, 0.5421047417, 0.7228063223] [ 0.2857142857, 0.6776309272, -0.6776309272] R = [ 1.166666667, 0.45, 0.6428571429] [ 0, 0.1054092553, 0.1016446391] [ 0, 0, 0.003764616262] p = 1 3 2
A(:,p) = Q*R
であることを検証します。
A(:,p) - Q*R
ans = [ -3.469446952e-18, -4.33680869e-18, -6.938893904e-18] [ 0, -8.67361738e-19, -1.734723476e-18] [ 0, -4.33680869e-19, -1.734723476e-18]
正確なシンボリック計算により丸め誤差を回避することができます。
A = sym(hilb(3)); [Q,R] = qr(A); A - Q*R
ans = [ 0, 0, 0] [ 0, 0, 0] [ 0, 0, 0]
小数点以下の有効桁数を既定の設定に戻します。
digits(previoussetting)
QR 分解を使用した行列方程式の求解
qr
を使用して、方程式系を行列形式で解くことができます。
たとえば、方程式系 A*X = b
を解く必要があるとします。A
と b
は、次の行列とベクトルです。
A = sym(invhilb(5)) b = sym([1:5]')
A = [ 25, -300, 1050, -1400, 630] [ -300, 4800, -18900, 26880, -12600] [ 1050, -18900, 79380, -117600, 56700] [ -1400, 26880, -117600, 179200, -88200] [ 630, -12600, 56700, -88200, 44100] b = 1 2 3 4 5
qr
を使用して、C = Q'*B
と A = Q*R
になるように、行列 C
と R
を求めます。
[C,R] = qr(A,b);
解 X
を計算します。
X = R\C
X = 5 71/20 197/70 657/280 1271/630
isAlways
を使用して、X
が系 A*X = b
の解であることを検証します。
isAlways(A*X == b)
ans = 5×1 logical array 1 1 1 1 1
置換情報をもつ QR 分解を使用した行列方程式の求解
浮動小数点数をもつ方程式系を解く際に、置換行列または置換ベクトルをもつ QR 分解を使用します。
たとえば、方程式系 A*X = b
を解く必要があるとします。A
と b
は、次の行列とベクトルです。
previoussetting = digits(10); A = vpa([2 -3 -1; 1 1 -1; 0 1 -1]); b = vpa([2; 0; -1]);
qr
を使用して、C = Q'*B
と A = Q*R
になるように、行列 C
と R
を求めます。
[C,R,P] = qr(A,b)
C = -2.110579412 -0.2132007164 0.7071067812 R = [ 3.31662479, 0.3015113446, -1.507556723] [ 0, 1.705605731, -1.492405014] [ 0, 0, 0.7071067812] P = 0 0 1 1 0 0 0 1 0
解 X
を計算します。
X = P*(R\C)
X = 1.0 -0.25 0.75
または、置換情報をベクトルとして返します。
[C,R,p] = qr(A,b,'vector')
C = -2.110579412 -0.2132007164 0.7071067812 R = [ 3.31662479, 0.3015113446, -1.507556723] [ 0, 1.705605731, -1.492405014] [ 0, 0, 0.7071067812] p = 2 3 1
この場合、次のように解 X
を計算します。
X(p,:) = R\C
X = 1.0 -0.25 0.75
小数点以下の有効桁数を既定の設定に戻します。
digits(previoussetting)
"エコノミー サイズ" 分解
'econ'
を使用して、"エコノミー サイズ" QR 分解を行います。
4
行 4
列のパスカル行列の最初の 2 列の行列を作成します。
A = sym(pascal(4)); A = A(:,1:2)
A = [ 1, 1] [ 1, 2] [ 1, 3] [ 1, 4]
この行列の QR 分解を計算します。
[Q,R] = qr(A)
Q = [ 1/2, -(3*5^(1/2))/10, (3^(1/2)*10^(1/2))/10, 0] [ 1/2, -5^(1/2)/10, -(2*3^(1/2)*10^(1/2))/15, 6^(1/2)/6] [ 1/2, 5^(1/2)/10, -(3^(1/2)*10^(1/2))/30, -6^(1/2)/3] [ 1/2, (3*5^(1/2))/10, (3^(1/2)*10^(1/2))/15, 6^(1/2)/6] R = [ 2, 5] [ 0, 5^(1/2)] [ 0, 0] [ 0, 0]
次に、この行列の "エコノミー サイズ" QR 分解を計算します。行数は列数を上回るので、qr
は Q
の最初の 2
列と R
の最初の2
行のみを計算します。
[Q,R] = qr(A,'econ')
Q = [ 1/2, -(3*5^(1/2))/10] [ 1/2, -5^(1/2)/10] [ 1/2, 5^(1/2)/10] [ 1/2, (3*5^(1/2))/10] R = [ 2, 5] [ 0, 5^(1/2)]
複素共役の回避
'real'
フラグを使用して、結果に複素共役が含まれないようにします。
行列を作成します。その行列の 1 つの要素は変数です。
syms x A = [1 2; 3 x]
A = [ 1, 2] [ 3, x]
この行列の QR 分解を計算します。既定の設定では、qr
は x
が複素数を表現すると仮定するので、結果には関数 abs
を使用する式が含まれます。
[Q,R] = qr(A)
Q = [ 10^(1/2)/10, -((3*x)/10 - 9/5)/(abs(x/10 - 3/5)^2... + abs((3*x)/10 - 9/5)^2)^(1/2)] [ (3*10^(1/2))/10, (x/10 - 3/5)/(abs(x/10 - 3/5)^2... + abs((3*x)/10 - 9/5)^2)^(1/2)] R = [ 10^(1/2), (10^(1/2)*(3*x + 2))/10] [ 0, (abs(x/10 - 3/5)^2 + abs((3*x)/10 - 9/5)^2)^(1/2)]
'real'
を使用すると、qr
はすべてのシンボリック変数が実数を表現すると仮定するので、短い結果を返すことができます。
[Q,R] = qr(A,'real')
Q = [ 10^(1/2)/10, -((3*x)/10 - 9/5)/(x^2/10 - (6*x)/5... + 18/5)^(1/2)] [ (3*10^(1/2))/10, (x/10 - 3/5)/(x^2/10 - (6*x)/5... + 18/5)^(1/2)] R = [ 10^(1/2), (10^(1/2)*(3*x + 2))/10] [ 0, (x^2/10 - (6*x)/5 + 18/5)^(1/2)]
入力引数
出力引数
詳細
ヒント
上三角行列
R
は次の条件を満たします。R = chol(A'*A)
。引数
'econ'
と0
は返される行列の形状のみに影響を与えます。シンボリック オブジェクトではない数値行列 (
sym
、syms
またはvpa
により作成されていない行列) でqr
を呼び出すと、MATLAB® 関数qr
が起動します。'vector'
の代わりに'matrix'
を使用した場合、qr
は既定の場合と同様に置換行列を返します。'matrix'
と'econ'
を使用すると、qr
はエラーをスローします。多くのシンボリック変数が含まれる行列計算は低速になる可能性があります。計算速度を向上させるには、特定の値を変数に代入することでシンボリック変数の数を減らします。