10-07:Python线性方程组J、GS和超松弛迭代

为什么今天会这么冷。

今天上了数学的两节理论课，第一节课先是推导迭代法的通式，然后从迭代误差推导出收敛条件。第二节课介绍雅克比迭代和高斯-赛德尔迭代。
下面是我的雅克比迭代、高斯赛德尔迭代和超松弛迭代的代码，主要运行numpy的矩阵运算。

def Jocobi(A,b,initial,delta):
    '''
    输入：A是系数矩阵，N阶方阵
          b是N*1列向量
          initial是解的初始值，N*1大小

    输出：迭代后的解析解
    '''
    D = np.diag(np.diag(A)) # 获得D矩阵
    L = -np.tril(A,-1)      # 获得L矩阵
    U = -np.triu(A,1)       # 获得U矩阵
    d = np.linalg.inv(D)    # 对D矩阵求逆
    
    BJ = np.dot(d,L+U)      # 迭代矩阵BJ
    lamda,_ = np.linalg.eig(BJ)
    if np.max(lamda)<1:     # 谱半径小于1
        f = np.dot(d,b)
        X = np.dot(BJ,initial)+f    # 初次的解
        times = 0
        while np.linalg.norm(X - initial) > delta:
            initial = X
            X = np.dot(BJ,initial) + f
            times = times +1
        return X,times
    else:
        print('Sorry,不可收敛。')

def Seidel(A,b,initial,delta):
    '''
    输入：A是系数矩阵，N阶方阵
          b是N*1列向量
          initial是解的初始值，N*1大小
    输出：迭代后的解析解
    '''
    D = np.diag(np.diag(A)) 
    L = -np.tril(A,-1)
    U = -np.triu(A,1)
    d = np.linalg.inv( D - L )   # 这里是和雅克比不同的地方
    
    BG = np.dot(d,U)             # 迭代矩阵BG
    lamda,_ = np.linalg.eig(BG)
    if np.max(lamda)<1:          # 谱半径小于1
        f = np.dot(d,b)
        X = np.dot( BG ,initial ) + f
        times = 0
        while np.linalg.norm( X - initial ) > delta:
            initial = X
            X = np.dot( BG,initial )+f  
            times = times +1
        return X,times
    else:
        print('Sorry,不可收敛。')

def SOR(A,b,initial,delta,w,maxTimes):
    '''
    升级版高斯赛德尔迭代
    输入：A是系数矩阵，N阶方阵
          b是N*1列向量
          initial是解的初始值，N*1大小
    输出：迭代后的解析解
    '''
    D = np.diag(np.diag(A)) 
    L = -np.tril(A,-1)
    U = -np.triu(A,1)
    d = np.linalg.inv(D - w*L) 
    
    Bw = np.dot(d,(1-w)*D+w*U)         # 迭代矩阵Bw
    f = np.dot(d,w*b)
    X = np.dot( Bw ,initial ) + f
    for i in np.arange(maxTimes):
        if np.linalg.norm( X - initial ) > delta:
            initial = X
            X = np.dot( Bw,initial )+f  
    return X

今天课余时间很少，写的很粗糙，之后再自我迭代改进。