いぬおさんのおもしろ数学実験室

おいしい紅茶でも飲みながら数学、物理、工学、プログラミング、そして読書を楽しみましょう

写真から立体を再現(2) 一般逆行列を試す

 前回、一般逆行列について書きました。

変数の数より式の本数の方が多いとき(行列を使った表現で言えば行が多いとき)、連立方程式を満たすは存在しない可能性があります。そんなときでも

f:id:Inuosann:20200212225126p:plain:w100

の両辺の左から一般逆行列(疑似逆行列

f:id:Inuosann:20200212225151p:plain:w150

かけると一応、が求まります。この x は|A|が最小になる、言わば「ベストの」解です。
 実験してみます。Pythonのライブラリに一般逆行列を求める関数があるので、今回は上の式を使うのではなくこれを使ってみます。
連立方程式

x + y = 3 ……①

x - y = 1 ……②

     y = 2 ……③

を解きます。

行列で書けば

f:id:Inuosann:20200212230538p:plain:w200

 です。両辺の左から
f:id:Inuosann:20200212230627p:plain:w100
の一般逆行列をかけてみます。

#=========================
#一般逆行列を試す
#
import numpy as np
import math
import sys
A = np.array([[1, 1], [1, -1], [0, 1]]) #行列をセットする
ginv = np.linalg.pinv(A) #Aの一般逆行列を求める
print(ginv)
arr = np.dot(ginv, [3, 1, 2]) #一般逆行列を左からかける
print(arr)  #かけた結果の表示
sys.exit()
#=========================

 以下が結果です。
 f:id:Inuosann:20200212231459p:plain:w400

 土台無理な方程式を解かせているのでもちろんこの「解」が3本の式を満たすわけではありません。こうして計算された「解」をAにかけてみます。

#=========================
#一般逆行列を試す
#
import numpy as np
import math
import sys
#行列をセットする
A = np.array([[1, 1], [1, -1], [0, 1]])
ginv = np.linalg.pinv(A) #Aの一般逆行列
print(ginv)
arr = np.dot(ginv, [3, 1, 2]) #一般逆行列を左からかける
print(arr) 

print("\n計算された「解」をAにかけてみると……")
arr1 = np.dot(A, arr)
print(arr1)
sys.exit()

f:id:Inuosann:20200212232318p:plain:w400
となって、[3, 1, 2]には戻りません。もちろんこれは仕方ないのですね。成立させるのは土台不可能な3本の方程式をまあまあ満たす、「ベストの」解が
f:id:Inuosann:20200212233702p:plain:w300
なのです。