主要内容

このページの翻訳は最新ではありません。ここをクリックして,英語の最新版を参照してください。

線形方程式

計算上の考慮事項

技術計算の中で最も重要な問題の1つは,線形方程式を解くことです。

行列表記では一般的な問題は次のような形式になります。[2つの行列Aとがあるとき,Ax = bまたはxA = bの条件のいずれかを満たす一意の行列xは存在しますか。]

1行1列の簡単な例を考えます。たとえば次の方程式について考えます。

7 x = 21

上記に一意の解は存在しますか。

もちろん存在します。この方程式には一意の解x = 3が存在します。解は除算から簡単に得られます。

X = 21/7 = 3

解を求める場合7の逆数すなわち7-1= 0.142857……を計算し,この7-1に21を乗算する方法は通常“使いません”。この方法では余計な計算ステップが必要となり,7-1を有限小数値で表すと精度が低くなります。複数の未知数がある線形方程式にも同様の考え方で,MATLAB®は逆行列を使わずにこのような方程式を解きます。

標準の数学的な表記ではありませんが,MATLABではスカラーの場合と同じ除算記号を使用して一般的な連立方程式の解を記述します。2つの除算記号”スラッシュ”(/)および”バックスラッシュ”(\)が2つのMATLAB関数mrdivideおよびmldivideに対応します。これらの演算子は,係数行列の左右どちらかが未知の行列となる2つの状況で使用されます。

x = b / A

mrdivideを使用して得られる行列方程式xA = bの解です。

x = A \ b

mldivideを使用して得られる行列方程式Ax = bの解です。

Ax = bまたはxA = bの方程式の両辺を一で割ると考えます。係数行列一个は常に“分”母になります。

x = A \ bに対する次元の整合性の条件は,2つの行列一个bの行数が同じであることです。その結果,解xbと同じ列数になり,行数は一个の列数と同じになります。x = b / Aに対しては,行と列の処理が反転します。

実際にはAx = bの形式の線形方程式の方がxA = bの形式より頻繁に使用されます。そのためバックスラッシュはスラッシュよりも頻繁に使用されます。この節の以降ではバックスラッシュ演算子を中心に説明します。スラッシュ演算子の特性は次の等式から推定できます。

(b / A) =(“\ b”)。

係数行列一个は正方である必要はありません。一个のサイズがm行n列である場合,次の3つのケースがあります。

m = n

正方システム。厳密解が得られます。

m > n

未知数より方程式が多い過決定システム。最小二乗解を求めます。

m < n

未知数より方程式が少ない劣決定システム。最大でm個のゼロでない構成要素をもつ基本解を求めます。

mldivideアルゴリズム

mldivide演算子は,さまざまなソルバーを使って各種の係数行列を処理します。実行するアルゴリズムは係数行列を調べて自動的に決められます。詳細については,関数mldivideのリファレンスページの”アルゴリズム”の節を参照してください。

一般解

線形方程式系Ax = bの一般解は,すべての解を示します。次のようにして一般解を求めることができます。

  1. 対応する同次方程式Ax = 0を解きます。コマンドを使用して”零(A)“と入力します。これによって,Ax = 0の解空間の基底が返されます。どの解も基底ベクトルの線形結合になります。

  2. 非同次方程式Ax = bの特殊解を導出します。

Ax = bの解はすべて,手順2のAx = bの特殊解に手順1の基底ベクトルの線形結合を加えたものとして表現できます。

この節の以降では,MATLABを使用して手順2のAx = bの特殊解を導出する方法を説明します。

正方システム

一般的な例は,正方係数行列一个と右辺に単一の列ベクトルbがある場合です。

特異でない係数行列

行列一个が特異でない場合,x = A \ bの解はbと同じサイズになります。たとえば,

一个=帕斯卡(3);u = [3;1;4);x = A\u x = 10 -12 5

* xuと等価であることを確認してください。

一个bが正方で同じサイズの場合,x = A \ bも同じサイズになります。

b =魔法(3);X = A \ b X = 19 13 6 0 6 3 1 -17 4

* xbと等価であることを確認してください。

上記2つの例は,厳密に整数の解になります。これは係数行列がフルランク行列(特異でない)の帕斯卡(3)であるためです。

特異な係数行列

正方行列一个が線形に独立した列をもっていない場合を“特”異といいます。一个が特異な場合、Ax = b の解は存在しないかまたは一意ではありません。バックスラッシュ演算一个\ bは,一个が特異に近い場合か,厳密な特異性が検出された場合に警告を発します。

一个が特異でAx = bが解をもつ場合は,次のように入力することで一意ではない特殊解を導出できます。

P = pinv b (A) *

pinv (A)は一の疑似逆行列です。一个x=bが厳密解をもたない場合、pinv (A)は最小二乗解を返します。

たとえば,

A = [1 3 7 -1 4 4 1 10 18]

が特異であることは以下のように入力することでわかります。

rank(A) ans = 2

一个はフルランクではないため,いくつかの特異値がゼロに等しくなります。

厳密解。b = (5; 2; 12)に対して,方程式Ax = bは次のような厳密解をもちます。

pinv(A)*b ans = 0.3850 -0.1103 0.7066

以下のように入力すると,pinv b (A) *が厳密解であることを確かめることができます。

A*pinv(A)*b ans = 5.0000、2.0000、12.0000

最小二乗解。ただし,b =[3、6 0]の場合,Ax = bは厳密解をもちません。この例の場合,pinv b (A) *は最小二乗解を返します。以下のように入力すると,

A*pinv(A)*b ans = -1.0000 4.0000 2.0000

元のベクトルbには戻りません。

Ax = bが厳密解をもつかどうかは,拡大行列[b]を行の階段型に変換することによって判定できます。例では以下のようになります。

rref([A b]) ans = 1.0000 0 2.2857 0 1.0000 1.5714 000 1.0000

最終行は最後の要素以外はすべてゼロであるため,方程式は解をもちません。この例の場合,pinv (A)は最小二乗解を返します。

過決定システム

この例では,実験データでのさまざまな曲線近似において過決定システムが頻繁に発生することを説明します。

yは,時間tのいくつかの異なる値で測定され,次の観測値を生成します。次のステートメントで,データを入力して表内に表示できます。

T = [0.3 .8 1.1 . 1.6 . 2.3]';y =[。83.63.60.60.50]';表(t, B = y)
B =6×2表T y ___ ____ 0 0.82 0.3 0.72 0.8 0.63 1.1 0.6 1.6 0.55 2.3 0.5

減衰指数関数を使用してデータのモデル化を試します。

y t c 1 + c 2 e - t

上記の式は,ベクトルyが他の2つのベクトルの線形結合で近似されることを示しています。一方はすべて1の定数ベクトルであり,もう一方は成分exp (- t)を含むベクトルです。未知の係数 c 1 および c 2 は,モデルのデータの偏差の二乗和を最小化する最小二乗近似を実行することで計算可能です。2つの未知数に対し6つの方程式があり,6行2列の行列で表されます。

E = [ones(size(t)) exp(-t)]
E =6×21.0000 1.0000 0.7408 1.0000 0.4493 1.0000 0.3329 1.0000 0.2019 1.0000 0.1003

最小二乗解を求めるには,バックスラッシュ演算子を使用します。

c = E \ y
c =2×10.4760 - 0.3413

つまり,データの最小二乗近似は次のようになります。

y t 0 4 7 6 0 + 0 3. 4 1 3. e - t

次のステートメントでは,tを等間隔にインクリメントしてモデルを評価し,結果を元のデータとともにプロットします。

T =(0:0.1:2.5)”;Y = [ones(size(T)) exp(-T)]*c;情节(T Y“- - -”、t、y,“o”

图中包含一个轴。坐标轴包含两个类型线对象。

E * cyと完全に等しくはありませんが,その差は元のデータの測定誤差より小さくなっている可能性があります。

方形行列一个は,線形独立の列がない場合はランク落ちとなります。一个がランク落ちの場合,AX = Bの最小二乗解は一意ではありません。一个\ B一个がランク落ちであり,最小二乗解を生成する場合に警告を発します。lsqminnormを使用して,すべての解の中で最小のノルムをもつ解Xを求めることができます。

劣決定システム

この例では,劣決定システムの解が一意にならない理由について説明します。劣決定の線形システムは方程式より未知数の数が多くなります。MATLABの行列の左除算演算で基底最小二乗解を求めます。つまり、n列の係数行列には,最大個の非ゼロの成分があります。

ここでは小さな乱数の例をあげます。

R = [6 8 7 3;3 5 4 1] rng(0);b =兰迪(2,1)
R = 6 8 7 3 3 5 4 1 b = 7 8

線形システムRp = bには4つの未知数の中に2つの式があります。この係数行列に含まれる整数は小さいため,格式コマンドを使用して有理数形式の解を表示するのが適切です。特殊解は次のステートメントで得られます。

格式老鼠p = R \ b
P = 0 17/7 0 -29/7

R (: 2)Rの最大ノルムをもつ列であるため,ゼロ以外の要素の1つは(2页)になります。R (: 2)を取り除くとR (: 4)が最大ノルムをもつため,他のゼロ以外の要素は(4页)になります。

劣決定システムの完全な一般解は,ヌル空間ベクトルの任意の線形結合にpを付加することで特性化され,有理数で表示するオプションにより関数を使用して求められます。

Z = null (R,“r”
Z = -1/2 -7/6 -1/2 / 1/2

R * Zがゼロになります。残差R * x - bは,以下に示す任意のベクトルxの場合は小さくなります。

x = p + Z*q

Zの列はヌル空間ベクトルであるため,積Z *问はこれらのベクトルの線形結合です。

Z x 1 x 2 u w u x 1 + w x 2

説明のため,任意のを選択し,xを作成します。

q = [2;1);x = p + Z*q;

残差のノルムを計算します。

格式规范(R * x - b)
ans = 2.6645 e15汽油

解を無限に利用できる場合,最小ノルムをもつ解が特に重要となります。lsqminnormを使用して最小ノルムをもつ最小二乗解を計算することができます。この解は常模(p)に対する可能な最小値をもちます。

p = lsqminnorm (R, b)
P = -207/137 365/137 79/137 -424/137

複数の右辺の解法

同一の係数行列一个をもつが,異なる右辺bをもつ線形システムを解く問題もあります。bの異なる値が同時に使用できる場合は,bを複数の列をもつ行列として構成し,単一のバックスラッシュコマンドを使用してX = A\[b1 b2 b3…]のように方程式系全体を一度に解くことができます。

ただし,bの異なる値がすべて同時に使用できない場合は,いくつかの方程式系を連続して解かなくてはなりません。これらの方程式系のいずれかをスラッシュ(/)またはバックスラッシュ(\)を使用して解く場合,演算子は係数行列一个を因数分解し,この行列の分解を使用して解を計算します。ただし,異なるbをもつ類似の方程式系を連続して解く場合,演算子は毎回同じ一个の分解を演算するため冗長演算となります。

この問題を解決するには,一个の分解を事前計算してから,その因子を再利用して異なるbの値の解を求めます。しかし実際には,この方法で分解を事前に計算することは困難です。その理由は,問題を解決するためにどの分解(LU,低密度脂蛋白,コレスキーなど)を計算するかだけでなく,因子をどのように乗算するかについても知る必要があるからです。たとえば,陆分解では,元のシステムAx = bを解くために次の2つの線形システムを解かなければなりません。

陆[L U] =(一个);x = U \ (L \ b);

代わりに,連続する複数の右辺をもつ線形システムを解くために推奨される方法は,分解オブジェクトを使用する方法です。これらのオブジェクトによって,行列分解の事前計算によるパフォーマンス上の利点を活用できます。それらに,行列分解の使用方法に関する知識は必要的“ありません”。前述の陆分解を以下のように置き換えることができます。

dA =分解(,“陆”);x = dA \ b;

どの分解を使用すればよいかわからない場合は,分解(一)によりバックスラッシュと同様に一个のプロパティに基づいて正しいタイプを選択します。

この方法によるパフォーマンス上の利点を測定する簡単なテストを次に示します。このテストでは,バックスラッシュと(\)分解の両方を使用して,同じスパース線形システムの解決を100回繰り返します。

n = 1 e3;A = sprand(n,n,0.2) + speye(n);b = 1 (n, 1);%反斜杠的解决方案抽搐k = 1:100 x = A\b;结束toc
运行时间是9.006156秒。
%分解方案tic dA =分解(A);k = 1:100 x = dA\b;结束toc
经过时间是0.374347秒。

この問題では,分解を使用した方法がバックスラッシュのみを使用した方法に比べてはるかに高速ですが,構文はシンプルです。

反復法

係数行列一个が大規模でスパース性がある場合,分解する方法は一般に効果がありません。反復法は近似解を生成できます。MATLABは大規模でスパース性のある行列を扱う、いくつかの反復法を用意しています。

関数 説明
pcg

前処理付き共役勾配法。この方法はエルミート正定値係数行列一个に適しています。

bicg

双共役傾斜法

bicgstab

双共役傾斜安定化法

bicgstabl

BiCGStab (l)法

研究生院理事会

共役勾配二乗法

巨磁电阻

一般化最小残差法

lsqr

LSQR法

minres

最小残差法。この方法はエルミート係数行列一个に適しています。

qmr

準最小残差法

symmlq

対称LQ法

tfqmr

転置なし準最小残差法

マルチスレッド計算

MATLABは多数の線形代数関数と要素ごとに計算する数値関数のマルチスレッド計算をサポートしています。これらの関数は自動的にマルチスレッドで実行されます。関数や式に対し、複数の CPU を使用して計算を高速化するにはさまざまな条件があります。

  1. 関数は同時実行可能な部分へ簡単に分割できる処理を実行します。この分割部分は各処理間で通信をほとんど行わずに実行されるものでなければなりません。すなわちシーケンシャルな処理がほとんどないものを扱います。

  2. データを分割し,複数の実行スレッドを管理する時間を含めても同時実行する価値があるように,データサイズは十分に大きいものを扱います。たとえば配列が数千個以上の要素を含む場合は,ほとんどの関数に対して高速化の効果があります。

  3. 処理がメモリに拘束されないものを扱います。すなわち処理時間の大部分がメモリのアクセス時間にならないものを扱います。一般的には,シンプルな処理より複雑な処理をする関数の方が,高速化の効果があります。

関数发票lscovlinsolvemldivideはマルチスレッド計算を有効にすると,大規模な倍精度配列(10000個以上の要素)に対する計算速度が大幅に上昇します。

参考

||||

関連するトピック