csycsy 发表于 2014-6-3 14:09:07

求助三对角矩阵的追赶法

哪位大神帮帮忙啊!!急啊

楚香饭 发表于 2014-6-3 14:37:44

三对角矩阵,用追赶法。(取自徐世良的算法集)

其中,数据文件内容为:
32 -8
12 -22 14
14 -22 16
16 -22 18
18 -22 20
20 -22 22
22 -22 24
24 -22 26
26 -22 28
28 -22
20 0 0 0 0 0 0 0 0 0

解为
解为:0.490707967461251      -0.537168130154994       -1.26472817663892
-1.26897912899290      -0.426771667312277       0.672632382050108
   1.06060662506127       0.355643056093564      -0.678092760285062
-0.863027149453715

以下是代码:

Program www_fcode_cn
Implicit None
Integer , parameter :: N = 10 !// 10*10 矩阵
Integer , parameter :: M = 3*N-2
Integer , parameter :: DP = Kind(0.0D0)
Real( Kind = DP ) :: B(M) , D(N)
integer :: L
Open( 12 , File = "B.txt" )
read( 12 , * ) B(:) !// 三对角系数
read( 12 , * ) D(:) !// 右端向量
call ATRDE( B , N , M , D , L )
write(*,*) '解为:',D(:)
End Program www_fcode_cn

SUBROUTINE ATRDE(B,N,M,D,L)
DIMENSION B(M),D(N)
DOUBLE PRECISION B,D
L=1
IF (M.NE.(3*N-2)) THEN
    L=-1
    WRITE(*,10)
    RETURN
END IF
10FORMAT(1X,'ERR')
DO 20 K=1,N-1
    J=3*K-2
    IF (ABS(B(J))+1.0.EQ.1.0) THEN
      L=0
      WRITE(*,10)
      RETURN
    END IF
    B(J+1)=B(J+1)/B(J)
    D(K)=D(K)/B(J)
    B(J+3)=B(J+3)-B(J+2)*B(J+1)
    D(K+1)=D(K+1)-B(J+2)*D(K)
20CONTINUE
IF (ABS(B(3*N-2))+1.0.EQ.1.0) THEN
    L=0
    WRITE(*,10)
    RETURN
END IF
D(N)=D(N)/B(3*N-2)
DO 30 K=N-1,1,-1
    D(K)=D(K)-B(3*K-1)*D(K+1)
30CONTINUE
RETURN
END

csycsy 发表于 2014-6-3 14:59:53

chuxf 发表于 2014-6-3 14:37
三对角矩阵,用追赶法。(取自徐世良的算法集)

其中,数据文件内容为:


谢谢啊,你这是先把数字文件存为B。txt吗?然后再用这个代码解吗?

楚香饭 发表于 2014-6-3 15:17:21

csycsy 发表于 2014-6-3 14:59
谢谢啊,你这是先把数字文件存为B。txt吗?然后再用这个代码解吗?

是的{:3_48:}

csycsy 发表于 2014-6-3 16:33:42

chuxf 发表于 2014-6-3 15:17
是的

大神大神,错了错了,,,不好意思啊,那个数据文本为32 -8
12 -22 6
14 -22 4
16 -22 2
18 -22 0
20 -22 -2
22 -22 -4
24 -22 -6
26 -22 -8
28 -22
20 0 0 0 0 0 0 0 0 0


能不能帮我运行下结果???????????万分感谢!!!!!!

楚香饭 发表于 2014-6-3 17:07:33

csycsy 发表于 2014-6-3 16:33
大神大神,错了错了,,,不好意思啊,那个数据文本为32 -8
12 -22 6
14 -22 4

你这要求过分了点。我不是你的保姆,自己运行吧。

csycsy 发表于 2014-6-3 17:11:35

chuxf 发表于 2014-6-3 17:07
你这要求过分了点。我不是你的保姆,自己运行吧。

大神啊,主要我不是这个专业的,隔行如隔山啊!!万分感谢啊@@!!

楚香饭 发表于 2014-6-3 17:18:28

csycsy 发表于 2014-6-3 17:11
大神啊,主要我不是这个专业的,隔行如隔山啊!!万分感谢啊@@!!

中国没有大学开设 Fortran 专业。

1.你一开始强调需要源代码。现在给你源代码了,你又不用,改而要结果。这是对回答者的不尊重。
2.这里是讨论技术问题的,不是免费劳动力市场。

楚香饭 发表于 2014-6-3 17:21:32

解为:0.753906250000000       0.515625000000000       0.382812500000000
0.300781250000000       0.246093750000000       0.207495629370629
0.178485576923077       0.159555288461538       0.128906250000000
0.164062500000000

csycsy 发表于 2014-6-3 17:22:22

chuxf 发表于 2014-6-3 17:18
中国没有大学开设 Fortran 专业。

1.你一开始强调需要源代码。现在给你源代码了,你又不用,改而要结果 ...

万分感谢!!不好意思!!!!来重庆请你吃饭!
页: [1] 2 3
查看完整版本: 求助三对角矩阵的追赶法