快速验证矩阵求MP广义逆及最小范数解或最小二乘解是否正确
生活随笔
收集整理的這篇文章主要介紹了
快速验证矩阵求MP广义逆及最小范数解或最小二乘解是否正确
小編覺得挺不錯的,現在分享給大家,幫大家做個參考.
\quad設計了一個測試函數,傳入原矩陣AAA和其最大秩分解B,D,A=BDB,D,A=BDB,D,A=BD,列向量bbb。函數依次輸出以下內容
- 1、判斷當前給出的最大秩分解是否正確
- 2、給出AAA的M-P廣義逆A+A^+A+
- 3、判斷矩陣是否相容
- 4、若矩陣相容則輸出Ax=bAx=bAx=b的最小范數解和通解
- 5、若不相容則輸出Ax=bAx=bAx=b的最小二乘解、通解和最佳逼近解
\quad運行環境:Python3+numpy庫
import numpy as npdef same(A, B):m = len(A)n = len(A[0])A = A.tolist()B = B.tolist()for i in range(m):for j in range(n):if(abs(A[i][j]-B[i][j])>0.0001):return False;return True def run(A, B, D, b):"""A為原矩陣,A=BD為其最大秩分解,Ax=b"""A = np.mat(A)B = np.mat(B)D = np.mat(D)b = np.mat(b).reshape(len(b), 1)if(same(B*D, A)):print("最大秩分解正確")else:print("最大秩分解錯誤")Ap = D.T * np.linalg.inv(D*D.T) * np.linalg.inv(B.T*B) * B.Tprint("M-P廣義逆為:\n", Ap)if(same(A*Ap*b, b)):print("該方程相容")x = Ap*bprint("最小范數解為:\n", x)E = np.eye(len(x))print("通解為:\n{}+{}U".format(x, E-Ap*A))else:print("該方程不相容")x = Ap*bprint("最小二乘解為:\n", x)E = np.eye(len(x))print("最小二乘解通解為:\n{}+{}U".format(x, E-Ap*A))print("最佳逼近解為:\n", x)def test1():A = [[0, 1, -1],[-1, 0, 1],[1, -1, 0],[0, 1, -1]]B = [[0, 1], [-1, 0], [1, -1], [0, 1]]D = [[1, 0, -1], [0, 1, -1]]b = [1, -1, 0, 1]run(A, B, D, b)def test2():A = [[0, 0, 2], [1, 1, 0], [0, 0, 1], [1, 1, 1]]B = [[0, 2], [1, 0], [0, 1], [1, 1]]D = [[1, 1, 0], [0, 0, 1]]b = [1, 1, 1, 1]run(A, B, D, b)if __name__ == '__main__':test1()test2()\quad給出的測試用例test1來自于2014年第8題。輸出結果為
\quad第二個題為書上P233 21題第一小題,運行結果為:
總結
以上是生活随笔為你收集整理的快速验证矩阵求MP广义逆及最小范数解或最小二乘解是否正确的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 详解React的Transition工作
- 下一篇: 微信小程序 swiper 自定义组件