# coding:utf-8

import numpy as np
import sys

# 設定矩陣
def set_matrix():
    # 設定係數矩陣A
    matrix_a =np.mat([
        [2.0, 1.0, 2.0],
        [5.0, -1.0, 1.0],
        [1.0, -3.0, -4.0]],dtype=float)
    matrix_b = np.mat([5, 8, -4],dtype=float).T
    return matrix_a, matrix_b


# 高斯順序消去法
def gauss_shunxu(mat):
    for i in range(0,mat.shape[0]-1):
        # 判斷順序順序主子式的首個元素不為0
        if mat[i,i] == 0:
            print mat
            print 'break:(', i, '，', i, ')元素為0'
            break
        else:
            # i行下面的每一行分別跟i行做計算，消掉第i個元素
            mat[i+1:,:] = mat[i+1:,:] - (mat[i+1:,i]/mat[i,i])*mat[i,:]
    return mat

def huidai(mat):
    x = np.mat(np.zeros(mat.shape[0],dtype=float))
    # 先算x(n)   用 b(n)/a(nn)
    n = x.shape[1]-1
    x[0,n] = mat[n,n+1]/mat[n,n]
    # 再分別帶入上一部式子算出n-1
    for i in range(n):
        n -= 1
        x[0,n] = (mat[n,mat.shape[1]-1] - np.sum(np.multiply(x[0, n+1:], mat[n,n+1:mat.shape[1]-1])))/mat[n,n]
    return x


if __name__ == "__main__":
    # 增廣矩陣m為係數矩陣A加上列矩陣b
    m = np.hstack(set_matrix())            #按列合併：vstack()   按行合併：hstack()
    print '原矩陣：'
    print m

    # 順序消去過程
    m1 = gauss_shunxu(m)
    print '\n上三角矩陣：'
    print m1

    # 迴帶過程
    x = huidai(m1)
    print '\n\nX的值為:', x