美文网首页
追赶法求解三对角矩阵(C/C++/python/matlab/F

追赶法求解三对角矩阵(C/C++/python/matlab/F

作者: 平平又无奇 | 来源:发表于2017-11-02 14:48 被阅读197次

    原文地址
    Python地址

    追赶法求解过程
    import numpy as np
    
    ## Tri Diagonal Matrix Algorithm(a.k.a Thomas algorithm) solver
    def TDMAsolver(a, b, c, d):
        '''
        TDMA solver, a b c d can be NumPy array type or Python list type.
        refer to http://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm
        and to http://www.cfd-online.com/Wiki/Tridiagonal_matrix_algorithm_-_TDMA_(Thomas_algorithm)
        '''
        nf = len(d) # number of equations
        ac, bc, cc, dc = map(np.array, (a, b, c, d)) # copy arrays
        for it in range(1, nf):
            mc = ac[it-1]/bc[it-1]
            bc[it] = bc[it] - mc*cc[it-1] 
            dc[it] = dc[it] - mc*dc[it-1]
                    
        xc = bc
        xc[-1] = dc[-1]/bc[-1]
    
        for il in range(nf-2, -1, -1):
            xc[il] = (dc[il]-cc[il]*xc[il+1])/bc[il]
    
        return xc
    

    例如:

    A = np.array([[10,2,0,0],[3,10,4,0],[0,1,7,5],[0,0,3,4]],dtype=float)   
    
    a = np.array([3.,1,3]) 
    b = np.array([10.,10.,7.,4.])
    c = np.array([2.,4.,5.])
    d = np.array([3,4,5,6.])
    
    print(TDMAsolver(a, b, c, d))
    >> [ 0.14877589  0.75612053 -1.00188324  2.25141243]
    #compare against numpy linear algebra library
    print(np.linalg.solve(A, d))
    >> [ 0.14877589  0.75612053 -1.00188324  2.25141243]
    

    相关文章

      网友评论

          本文标题:追赶法求解三对角矩阵(C/C++/python/matlab/F

          本文链接:https://www.haomeiwen.com/subject/mygopxtx.html