AX=B 的天际线求解器,其中 A 是对称天际线矩阵

计算科学 线性求解器 正则 数据存储
2021-12-15 00:17:47

我正在 Fortran 90(GNU 编译器)中寻找一个简单的子程序来求解类型的线性方程,其中是以对称天际线矩阵的形式存储的对称矩阵。我想要这种类型的解决方案:AX=BAn×n

X = skysolve(A,p,X,B,n)
A = {11,22,13,0,33,24,34,44,55,16,0,0,46,56,66}

p对角线元素的索引在哪里p = {1,2,5,8,9,15}

如果我能对这种情况有一个简单的子程序类型的解决方案而不是一个库,那就更好了。

1个回答

我建议看一下 Robert Taylor 为本书编写的 FEAPpv 有限元代码:

http://www.amazon.com/The-Finite-Element-Method-Fundamentals/dp/1856176339

您可以从此页面下载 fortran 源代码:

http://www.ce.berkeley.edu/projects/feap/feappv/

该代码包含两个子程序,datri分别dasol执行以天际线格式存储的矩阵的分解和求解;您首先调用datri矩阵的因子,然后调用dasol执行前向和后向替换。

两个子例程中都有足够详细的注释以便能够使用它们。这是一个包含五个方程和一个对称矩阵的简单示例(不包括对角线下方的条目):

      program test
C
C  Solve A*U = b      
C        [2   -2    0    0   -1]
C        [0    3   -2    0    0]
C  A =   [0    0    5   -3    0]
C        [0    0    0   10    4]
C        [0    0    0    0   10] 
C
C The solution, U, is: 636.,619.,292.,74.,34.
C
        parameter (neqs = 5)
        real*8 au(7), ad(neqs), b(neqs), energy
        integer jp(neqs)
        data au/-2.,-2.,-3.,-1.,0.,0.,4./
        data ad/2.,3.,5.,10.,10./
        data b/0.,1.,0.,0.,0./
        data jp/0,1,2,3,7/
        call datri(0, au, ad, jp, neqs, neqs)
        call dasol(0, au, ad, b, jp, neqs, neqs, energy)
        write(6,*) 'solution', (b(i), i=1,neqs)
      end