我已經想出了如何撰寫一個函式,以便您可以行歸約和解決線性代數問題,我遇到的唯一問題是設定 if 條件,當用于行歸約的常數時進行部分旋轉等于0。我在下面嘗試過,邏輯是有道理的,但我很難理解為什么我的解決方案不起作用。
import numpy as np
def gaussElim(A,B):
M = np.concatenate((A,B), axis=1)# Combines the two matrices (assuming they have the same ammount of rows and matrix B has been )
nr, nc = M.shape
for r in range (nr):
const = M[r][r]
if const == 0: # **This is the condition that is tripping me up**
for i in range (nr-1):
M[r][i]=M[r 1][i]
M[r 1][i] = M[r][i]
const = M[r][r]
for c in range (r,nc):
M[r][c] = M[r][c]/const
for rr in range(nr):
if rr!= r :
const = M[rr][r]
for c in range(r,nc):
M[rr][c] = M[rr][c] - const * M[r][c]
return M[:, nc-1]
Mrx = np.array([ [1.0,3,2,4,3,1], [-4,0,3,2,3,4], [3,-1,3,2,2,5], [3,3,12,2,-
6,-4], [-1,-2,-3,7,6,4], [7,5,0,0,4,2] ])
Rhs = np.array([[ 4, 5, 6, 10, 6, -8 ]]) # This is a row vectorv
RhsT = Rhs.T # Rhs.T is transpose of Rhs and a column vector
S = gaussElim(Mrx,RhsT)
print(S)
A1 = np.array([[2.0,1,-1],[2,1,-2],[1,-1,1]])
b1 = np.array([[1.0],[-2],[2]])
S1 = gaussElim(A1,b1)
print (S1)
x = np.linalg.solve(A1,b1)
print(x)
是的,我已經查看了其他人對高斯消除的解決方案,但我想更具體地了解為什么我的部分旋轉解決方案不起作用。謝謝!
將 A1,b1 輸入我的函式給了我
[1, 0.8, 1.8]
正確答案是
[1, 2, 3]
列印陳述句是為了讓我可以看到我的功能是否按預期作業。
uj5u.com熱心網友回復:
一個問題是這段代碼:
for i in range (nr-1):
M[r][i]=M[r 1][i] # line 1
M[r 1][i] = M[r][i] # line 2
const = M[r][r] # line 3
我評論的那一行line 2只是撤銷了line 1. 如果您嘗試交換值,請嘗試將第 1 行和第 2 行替換為:
M[r][i], M[r 1][i] = M[r 1][i], M[r][i]
...或者,如果您更喜歡明確說明交換:
temp = M[r][i]
M[r][i]=M[r 1][i] # line 1
M[r 1][i] = temp # line 2
您的代碼中的第二個(可能是良性的)問題是line 3上面的 ( const = M[r][r]) 不需要像當前那樣位于回圈內,我相信您可以在最終結果沒有變化的情況下減少一個級別。
第三個(良性)問題是M[r][c] = M[r][c]/const可以簡化為M[r][c] /= const,假設原始陣列 A 和 B(因此 M)是浮點數。
同樣,M[rr][c] = M[rr][c] - const * M[r][c]可以簡化為M[rr][c] -= const * M[r][c]。
把它們放在一起(最重要的是,糾正上面的第一個問題):
import numpy as np
def gaussElim(A,B):
M = np.concatenate((A,B), axis=1)# Combines the two matrices (assuming they have the same ammount of rows and matrix B has been )
nr, nc = M.shape
for r in range (nr):
const = M[r][r]
if const == 0: # **This is the condition that is tripping me up**
for i in range (nr-1):
M[r][i], M[r 1][i] = M[r 1][i], M[r][i]
const = M[r][r]
for c in range (r,nc):
M[r][c] /= const
for rr in range(nr):
if rr!= r :
const = M[rr][r]
for c in range(r,nc):
M[rr][c] -= const * M[r][c]
return M[:, nc-1]
Mrx = np.array([ [1.0,3,2,4,3,1], [-4,0,3,2,3,4], [3,-1,3,2,2,5], [3,3,12,2,-
6,-4], [-1,-2,-3,7,6,4], [7,5,0,0,4,2] ])
Rhs = np.array([[ 4, 5, 6, 10, 6, -8 ]]) # This is a row vectorv
RhsT = Rhs.T # Rhs.T is transpose of Rhs and a column vector
S = gaussElim(Mrx,RhsT)
print(S)
A1 = np.array([[2.0,1,-1],[2,1,-2],[1,-1,1]])
b1 = np.array([[1.0],[-2],[2]])
S1 = gaussElim(A1,b1)
print (S1)
x = np.linalg.solve(A1,b1)
print(x)
輸出:
[-0.97124946 2.88607869 -1.80198876 3.47492434 -6.16688284 4.51794207]
[0.33333333 1.33333333 1. ]
[[1.]
[2.]
[3.]]
轉載請註明出處,本文鏈接:https://www.uj5u.com/qiye/473664.html
