どうも、木村（@kimu3_slime）です。
今回は、Julia（SymPy）で行列の和・積、逆行列、線形方程式を解く方法を紹介します。
予備知識：Julia（SymPy）でベクトルの計算をする方法（内積・ノルム・距離・角度、外積）
準備
SymPyを使うので、持っていなければインストールしておきましょう。
1 2 | using Pkg Pkg.add("SymPy") |
準備として、以下のコードを実行しておきます。
1 | using SymPy |
Julia（SymPy）でベクトルの計算をする方法
行列の用意
「Sym[1行目; 2行目]」で、2行の行列を表すことができます。セミコロン;が改行の印です。
1 | B = Sym[1 2; 3 4] |
B = Sym[1 2; 3 4]
より一般に、文字を使って行列を表しましょう。
1 2 | @vars a11 a12 a21 a22 x y z n α β A = Sym[a11 a12; a21 a22] |
\[ \begin{aligned}\left[ \begin{array}{rr}a_{11}&a_{12}\\a_{21}&a_{22}\end{array}\right]\end{aligned} \]
行列の和、積、転置
通常の数と同じ演算記号「+,*,^」が行列に対しても使えます。
1 2 3 4 | A+B A*B B*A A*[x,y] |
\[ \begin{aligned}\left[ \begin{array}{rr}a_{11} + 1&a_{12} + 2\\a_{21} + 3&a_{22} + 4\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{rr}a_{11} + 3 a_{12}&2 a_{11} + 4 a_{12}\\a_{21} + 3 a_{22}&2 a_{21} + 4 a_{22}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{rr}a_{11} + 2 a_{21}&a_{12} + 2 a_{22}\\3 a_{11} + 4 a_{21}&3 a_{12} + 4 a_{22}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{r}a_{11} x + a_{12} y\\a_{21} x + a_{22} y\end{array} \right]\end{aligned} \]
行列の積では、一般に\(AB\)と\(BA\)が異なることが、見てわかりますね。
1 2 3 4 5 | A^2 A^3 B^2 B^3 B^n |
\[ \begin{aligned}\left[ \begin{array}{rr}a_{11}^{2} + a_{12} a_{21}&a_{11} a_{12} + a_{12} a_{22}\\a_{11} a_{21} + a_{21} a_{22}&a_{12} a_{21} + a_{22}^{2}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{rr}a_{11} \left(a_{11}^{2} + a_{12} a_{21}\right) + a_{12} \left(a_{11} a_{21} + a_{21} a_{22}\right)&a_{11} \left(a_{11} a_{12} + a_{12} a_{22}\right) + a_{12} \left(a_{12} a_{21} + a_{22}^{2}\right)\\a_{21} \left(a_{11}^{2} + a_{12} a_{21}\right) + a_{22} \left(a_{11} a_{21} + a_{21} a_{22}\right)&a_{21} \left(a_{11} a_{12} + a_{12} a_{22}\right) + a_{22} \left(a_{12} a_{21} + a_{22}^{2}\right)\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{rr}7&10\\15&22\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{rr}37&54\\81&118\end{array}\right]\end{aligned} \]
MethodError: no method matching schur!(::Matrix{Sym})
具体的なべき乗も、素早く計算してくれます。しかし、\(B^n\)のような一般的なべき乗は計算できませんでした。
「B^-1」のように書くことで、逆行列も計算できます。これについて詳しくは別記事で。
1 | B^-1 |
\[ \begin{aligned}\left[ \begin{array}{rr}-2&1\\\frac{3}{2}&- \frac{1}{2}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}R_\theta = \begin{pmatrix} \cos \theta& -\sin \theta\\ \sin\theta & \cos \theta \end{pmatrix}\end{aligned} \]
の性質を、計算で確かめてみましょう。
回転行列の積はまた、回転行列になります。「sympy.simplify(数式)」で数式を単純化できます。
1 2 3 4 | R_α= Sym[cos(α) -sin(α); sin(α) cos(α)] R_β= Sym[cos(β) -sin(β); sin(β) cos(β)] R_α*R_β sympy.simplify(R_α*R_β) |
\[ \begin{aligned}\left[ \begin{array}{rr}\cos{\left(β \right)}&- \sin{\left(β \right)}\\\sin{\left(β \right)}&\cos{\left(β \right)}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{rr}\cos{\left(β \right)}&- \sin{\left(β \right)}\\\sin{\left(β \right)}&\cos{\left(β \right)}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[ \begin{array}{rr}- \sin{\left(α \right)} \sin{\left(β \right)} + \cos{\left(α \right)} \cos{\left(β \right)}&- \sin{\left(α \right)} \cos{\left(β \right)} – \sin{\left(β \right)} \cos{\left(α \right)}\\\sin{\left(α \right)} \cos{\left(β \right)} + \sin{\left(β \right)} \cos{\left(α \right)}&- \sin{\left(α \right)} \sin{\left(β \right)} + \cos{\left(α \right)} \cos{\left(β \right)}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[\begin{smallmatrix}\cos{\left(α + β \right)} & – \sin{\left(α + β \right)}\\\sin{\left(α + β \right)} & \cos{\left(α + β \right)}\end{smallmatrix}\right]\end{aligned} \]
転置行列\(A^{\top}\)は、「行列.T」で求めることができます。
1 | R_α.T |
\[ \begin{aligned}\left[ \begin{array}{rr}\cos{\left(α \right)}&\sin{\left(α \right)}\\- \sin{\left(α \right)}&\cos{\left(α \right)}\end{array}\right]\end{aligned} \]
回転行列\(R_{\alpha}\)の転置行列は、\(R_{-\alpha}\)に等しいことを確かめてみましょう。「R_α.subs(α,-α)」で、変数\(\alpha\)に\(-\alpha\)を代入して置き換えることができます。
1 2 | R_α.subs(α,-α) R_α.T == R_α.subs(α,-α) |
\[ \begin{aligned}\left[ \begin{array}{rr}\cos{\left(α \right)}&\sin{\left(α \right)}\\- \sin{\left(α \right)}&\cos{\left(α \right)}\end{array}\right]\end{aligned} \]
回転行列の逆行列は、その転置行列であることを計算で確かめてみましょう。
1 2 | R_α*R_α.T sympy.simplify(R_α*R_α.T) |
\[ \begin{aligned}\left[ \begin{array}{rr}\sin^{2}{\left(α \right)} + \cos^{2}{\left(α \right)}&0\\0&\sin^{2}{\left(α \right)} + \cos^{2}{\left(α \right)}\end{array}\right]\end{aligned} \]
\[ \begin{aligned}\left[\begin{smallmatrix}1 & 0\\0 & 1\end{smallmatrix}\right]\end{aligned} \]
確かにそれらの積が、単位行列になっていますね。
線形方程式を解く
連立方程式を行列で表した、\(Ax=b\)という形の線形方程式を解きましょう。「linsolve((A, b),変数)」で解けます。
例えば、
\[ \begin{aligned}\left\{ \begin{array} {l} x-y=1\\ -2x+y=4 \end{array} \right.\end{aligned} \]
ならば、左辺を表す行列は「Sym[1 -1; -2 1]」で、右辺は「Sym[1,4]」です。
1 | linsolve((Sym[1 -1; -2 1], Sym[1,4]),x,y) |
\[ \begin{aligned}\left\{\left( -5, \ -6\right)\right\}\end{aligned} \]
解に不定性（自由度）がある場合は、変数として文字を使って解を表してくれます。
\[ \begin{aligned} \left\{ \begin{array} {l} x-y=1\\ 0x+0y=0 \end{array} \right.\end{aligned} \]
\[ \begin{aligned} \left\{ \begin{array} {l} x+2y+3z=4\\ 5x+6y+7z=8\\9x+10y+11z=12 \end{array} \right.\end{aligned} \]
1 2 | linsolve((Sym[1 -1; 0 0], Sym[1,0]),x,y) linsolve((Sym[1 2 3; 5 6 7; 9 10 11], Sym[4,8,12]),x,y,z) |
\[ \begin{aligned}\left\{\left( y + 1, \ y\right)\right\}\end{aligned} \]
\[ \begin{aligned}\left\{\left( z – 2, \ 3 – 2 z, \ z\right)\right\}\end{aligned} \]
参考：Julia（SymPy）で連立方程式を解く方法、1次方程式を行列で解くメリット・方法・条件について、幾何学的に見る
以上、Julia（SymPy）で行列の和・積、逆行列、線形方程式を解く方法を紹介してきました。
SymPyでは記号を扱えるので、回転行列の例のように文字を含んだ行列の計算ができるのが嬉しいですね。
木村すらいむ（@kimu3_slime）でした。ではでは。
コロナ社 (2020-03-26T00:00:01Z)
￥7,353 (コレクター商品)
こちらもおすすめ
Julia（SymPy）でベクトルの計算をする方法（内積・ノルム・距離・角度、外積）
1次方程式を行列で解くメリット・方法・条件について、幾何学的に見る
