用Python解線性方程組——Scipy包和自己寫
## 一、基於SymPy庫
用Python解決方程組、微積分等問題,主要是用到Python的一個庫——SymPy庫。可以說這個專案也主要是學習SymPy庫的用法。
解二元一次方程功能實現
解方程的功能主要是使用Sympy中solve函式實現。
示例題目是:
方程表示
程式碼表示與手寫還是有區別的,下面列出常用的:
加號 +
減號 -
除號 /
乘號 *
指數 **
對數 log()
e的指數次冪 exp()
題目中表達式可表示為:
2 * x - y - 3 = 0
3 * x + y - 7 = 0
完整程式碼為:
from sympy import *
x = Symbol('x')
y = Symbol('y')
print solve([2 * x - y - 3, 3 * x + y - 7],[x, y])
麻麻,我跟正確答案一樣哦~
二、自己試圖寫個簡單的函式
求解方程組需要先消元,我隨手寫了一個這樣的方程組:
x + 2y +3z = 6
3x + y = 3
-x + y = 3
2.1 先消元
按照最原始的方法消元(初中學的方法),程式碼如下:
# 寫一個n*n非奇異線性方程組求解函式
# Elimination:消元
def Elimination1(A = A):
for i in range(1, A.shape[0]):
A[i] = A[i] - A[0]*(A[i,0]/A[0,0])
return(A)
def Elimination(A = A):
j = A.shape[0] - 1
i = 0
while j > 0:
A[i:,i:] = Elimination1(A = A[i:,i:])
i = i + 1
j = j - 1
return(A)
第一個函式是消一次,第二個函式通過呼叫第一個函式,迴圈消元。
把方程組寫進來,執行函式輸出。
import numpy as np
A = np.matrix([[1,2,3,6],[3,1,0,3],[-1,1,0,3]])
Elimination(A=A)
輸出結果
matrix([[ 1, 2, 3, 6],
[ 0, -5, -9, -15],
[ 0, 0, -2, 0]])
消元后,下一步就可以求解了。
2.2 求解方程組
程式碼有問題,求解有錯誤
但是還沒查出來是哪裡出了問題
# 求解方程組
def SolveX(A = None):
A = Elimination(A=A)
X = {}
m = A.shape[0]
n = A.shape[1]
for i in range(m):
j = m-1-i
Y = A[j,-1]
k = i
while(k>0):
count = 0 # 問題在這!這樣內部迴圈count永遠不變。
Y = Y - (A[2-k-count,-2-count]*X[2-count])
count = count + 1
k = k-1
X[j] = Y/A[j,-2-i]
return X
第三個迴圈錯了,明天一步一步的核對,看哪個地方錯了。
分開一步一步核對是對的,然而成了函式又錯了,很奇怪哪裡出了問題。
修改如下:
# 寫一個n*n非奇異線性方程組求解函式
# Elimination:消元
def Elimination1(A = None):
for i in range(1, A.shape[0]):
A[i] = A[i] - A[0]*(A[i,0]/A[0,0])
return(A)
def Elimination(A = None):
j = A.shape[0] - 1
i = 0
while j > 0:
A[i:,i:] = Elimination1(A = A[i:,i:])
i = i + 1
j = j - 1
return(A)
# 求解方程組
def SolveX(A = None):
A = Elimination(A=A)
X = {}
m = A.shape[0]
for i in range(m):
j = m-1-i
Y = A[j,-1]
k = i
count = 0
while(k>0):
Y = Y - A[2-k-count,-2-count]*X[2-count]
count = count + 1
k = k-1
X[j] = Y/A[j,-2-i]
return X
import numpy as np
A = np.matrix([[1,2,3,6],[3,1,0,3],[-1,1,0,3]])
SolveX(A = A)
結果:
{0: 0.0, 1: 3.0, 2: -0.0}
正確了終於。
---------------------