はじめに
Aを行列としたN次元連立一次方程式 $A x = b$ を$x$について解く。
内容
(1) Aの逆行列$A^-1$をnumpyのlinalg.invメソッドにより計算し，
解ベクトル$x=A^{-1} b$を求める。
(2) np.linalg.solveメソッドから直接$x$を計算する。こちらの方が一般的に速い。
コード(1)
%%time
import numpy as np
"""
A x = bを満たすxベクトルを linalg.solveを使って求める
"""
# 行列Aの生成
a1_lis = [1, 0, 0]
a2_lis = [0, 2, 0]
a3_lis = [0, 0, 4]
A_matrix=np.array([a1_lis, a2_lis, a3_lis])
# b vector
b = np.array([4,5,6]) # 行ベクトルとしてｂを生成し，それの転置行列として"列ベクトル”を生成する
b_vec = b.T
A_inv = np.linalg.inv(A_matrix)
x_vec= np.dot(A_inv,b_vec) # xベクトルを計算する。行列の積に@演算子を使っている(参考)
print(x_vec)
または"matrix"を使う下の書き方もある。
可読性は高くなる。しかし，なにかと不便なことが多いのであまり推奨しない。
import numpy as np
"""
A x = bを満たすxベクトルを x = A^-1 bの計算により求める
"""
# 行列Aの生成
a1_lis = [1, 0, 0]
a2_lis = [0, 2, 0]
a3_lis = [0, 0, 4]
A_matrix=np.matrix([a1_lis, a2_lis, a3_lis])
# b vector
bb = np.matrix([4,5,6]) # 行ベクトルとしてｂを生成し，それの転置行列として"列ベクトル”を生成する
b_vec = bb.T
x_vec=A_matrix.I@b_vec # xベクトルを計算する。行列の積に@演算子を使っている(参考)
print(x_vec)
結果(1)
[[ 4. ]
[ 2.5]
[ 1.5]]
(2) numpy.linalg.solveを使う場合(推奨)
import numpy as np
"""
A x = bを満たすxベクトルを linalg.solveを使って求める
"""
# 行列Aの生成
a1_lis = [1, 0, 0]
a2_lis = [0, 2, 0]
a3_lis = [0, 0, 4]
A_matrix=np.array([a1_lis, a2_lis, a3_lis])
# b vector
b = np.array([4,5,6]) # 行ベクトルとしてｂを生成し，それの転置行列として"列ベクトル”を生成する
x_vec= np.linalg.solve(A_matrix,b) # xベクトルを計算する。行列の積に@演算子を使っている(参考)
print(x_vec)
結果(2)
[ 4. 2.5 1.5]
補遺
numpy.linalg.solveメソッドは内部的にはLapackのLAPACKのdgesv(実連立一次方程式の解)とzgesv(複素連立一次方程式の解)のラッパーである[1]。
● 計算物理学等では連立一次方程式を求めることが多い。行列$A$の性質によって逆行列を求めるためのさまざまなスキームが存在するが，コードを書くときはさしあたりsolveメソッドを使うだけで良いと思われる...。
連立一次方程式を解く初歩的な数値計算の手法としてガウスの消去法，ヤコビ法，ガウス-ザイデル法，共役勾配法などがある[2]。必要になったときに調べると良いだろう。
参考文献
[1] numpy.linalg.solveについての説明: https://docs.scipy.org/doc/numpy-1.7.0/reference/generated/numpy.linalg.solve.html
[2] 川上一郎，「数値計算」，岩波書店，1989.
@演算子について: [Pythonによる科学・技術計算] @演算子による行列積の計算, python3.5以降, numpy