技術標籤：數值分析非線性方程組的求解

import numpy as np


def fgauss(A, b):
    n = A.shape[0]
    zengguan = np.hstack((A, b.T))
    ra = np.linalg.matrix_rank(A)
    rz = np.linalg.matrix_rank(zengguan)
    temp1 = rz - ra
    if temp1 > 0:
        print("無一般意義下的解，係數矩陣與增廣矩陣的秩不同")
        return
    if
 
 ra == rz:
        if ra == n:
            x = np.mat(np.zeros(zengguan.shape[0], dtype=float))
            for p in range(0, n):
                for k in range(p+1, n):
                    m = zengguan[k, p]/zengguan[p, p]
                    zengguan[k, p:n+1] = zengguan[k, p:n+1] - m * zengguan[p, p:
 
n+1]
            b1 = zengguan[0:n, n]
            a1 = zengguan[0:n, 0:n]
            x[0, n-1] = b1[n-1]/a1[n-1, n-1]
            for i in range(n - 2, -1, -1):
                try:
                    x[0, i] = (b1[i] - np.sum(np.multiply(a1[i, i+1:n], x[0, i+1:n]))) / (a1[i, i])
                except:
                    print
 
("錯誤")
    return x


if __name__ == "__main__":
    x = [0, 0]
    err = 1
    k = 0
    while err > 0.0000001:
        k = k+1
        A1 = np.mat([[4+0.1*np.exp(x[0]), -1],
                    [-1+0.25*x[0], 4]], dtype=float)
        b2 = np.mat([-1*(4*x[0]-x[1]+0.1*np.exp(x[0])-1), -1*(-x[0]+4*x[1]+x[0]*x[0]/8)], dtype=float)
        dx = fgauss(A1, b2)
        x[0] = x[0] + dx[0, 0]
        x[1] = x[1] + dx[0, 1]
        err = max(abs(4*x[0]-x[1]+0.1*np.exp(x[0])-1), abs(-x[0]+4*x[1]+x[0]*x[0]/8))
        print('第{}次迭代，誤差為{}'.format(k, err))
        print('迭代結果為', x)