主题：关于三对角矩阵追赶法的程序！
LELEMA
[专家分：0] 发布于 2010-10-08 16:21:00
! 三对角矩阵追赶法
subroutine sdj(n,a,b,c,f)
doubleprecision a(0:n),b(0:n),c(0:n),f(0:n)
doubleprecision W
f(0)=f(0)/b(0)
w=b(0)
do k=0,n-1
b(k)=c(k)/w
w=b(k+1)-a(k)*b(k)
f(k+1)=(f(k+1)-a(k)*f(k))/w
end do
do j=1,n
k=n-j
f(k)=f(k)-B(K)*F(K+1)
END DO
RETURN
END
是不是应该把程序中的a(k)改成a(k+1)？谢谢
回复列表 （共10个回复）
沙发
jstzhurj [专家分：4680] 发布于 2010-10-08 16:59:00
为什么要用0下标？
板凳
jstzhurj [专家分：4680] 发布于 2010-10-08 17:08:00
program main
implicit none
integer,parameter::n=3 !变量x的个数
integer:: i,s
real(8):: a(n-1),b(n),c(n-1),f(n)
a=(/6,3/)
b=(/5,2,3/)
c=(/6,3/)
f=(/4,5,3/)
call sdj(n,a,b,c,f)
contains
!追赶法
subroutine sdj(n,a,b,c,f)
implicit none
integer:: n,i
real(8):: a(n-1),b(n),c(n-1),d(n-1),e(n-1),f(n)
real(8):: x(n)
!追赶法的核心
d(1)=f(1)/b(1)
e(1)=c(1)/b(1)
do i=2,n-1
e(i)=c(i)/(b(i)-a(i-1)*e(i-1))
d(i)=(f(i)-a(i-1)*d(i-1))/(b(i)-a(i-1)*e(i-1))
enddo
x(n)=(f(n)-a(n-1)*d(n-1))/(b(n)-a(n-1)*e(n-1))
do i=n-1,1,-1
x(i)=d(i)-e(i)*x(i+1)
end do
write(*,*) 'The roots is:'
write(*,*) x
return
end subroutine
end program
3 楼
LELEMA [专家分：0] 发布于 2010-10-08 17:15:00
我正在研究FORTRAN程序，我贴上去的这个程序是别人的，我现在需要首先给它看懂，然后改进，看的过程中我觉得应该把a(k)改成a(k+1)，不知道我理解的对，还是我朋友编的对，请各位大侠给个意见！谢谢
4 楼
LELEMA [专家分：0] 发布于 2010-10-09 09:14:00
e(i)=c(i)/(b(i)-a(i-1)*e(i-1))
为什么是a(i-1)？公式不是e(i)=c(i)/(b(i)-a(i)*e(i-1))么？
是我的公式导错了？
5 楼
jstzhurj [专家分：4680] 发布于 2010-10-09 09:39:00
[quote]e(i)=c(i)/(b(i)-a(i-1)*e(i-1))
为什么是a(i-1)？公式不是e(i)=c(i)/(b(i)-a(i)*e(i-1))么？
是我的公式导错了？[/quote]
注意a的下标起始点!若是从2开始不必减1。
6 楼
LELEMA [专家分：0] 发布于 2010-10-09 09:45:00
你的程序中i不是从2开始的么？怎么还减1了呢？a的下标不应该比e的下标高一个么？[em18]
7 楼
jstzhurj [专家分：4680] 发布于 2010-10-09 09:53:00
a(n-1) =a(1:n-1) ！
8 楼
LELEMA [专家分：0] 发布于 2010-10-09 10:45:00
[quote]
a(n-1) =a(1:n-1) ！[/quote]
这一点我知道啊，我再想想哦！[em8]
9 楼
jstzhurj [专家分：4680] 发布于 2010-10-09 10:59:00
[quote][quote]
a(n-1) =a(1:n-1) ！[/quote]
这一点我知道啊，我再想想哦！[em8][/quote]
为了便于你理解，我做了一点小的改动，效果一样!
program main
implicit none
integer,parameter::n=3 !变量x的个数
integer:: i,s
real(8):: a(2:n),b(n),c(n-1),f(n)
a=(/6,3/)
b=(/5,2,3/)
c=(/6,3/)
f=(/4,5,3/)
call sdj(n,a,b,c,f)
contains
!追赶法
subroutine sdj(n,a,b,c,f)
implicit none
integer:: n,i
real(8):: a(2:n),b(n),c(n-1),d(n-1),e(n-1),f(n)
real(8):: x(n)
!追赶法的核心
d(1)=f(1)/b(1)
e(1)=c(1)/b(1)
do i=2,n-1
e(i)=c(i)/(b(i)-a(i)*e(i-1))
d(i)=(f(i)-a(i)*d(i-1))/(b(i)-a(i)*e(i-1))
enddo
x(n)=(f(n)-a(n)*d(n-1))/(b(n)-a(n)*e(n-1))
do i=n-1,1,-1
x(i)=d(i)-e(i)*x(i+1)
end do
write(*,*) 'The roots is:'
write(*,*) x
return
end subroutine
end program
10 楼
LELEMA [专家分：0] 发布于 2010-10-09 15:48:00
现在的程序我是懂了
我也好像明白我的问题了，太谢谢了[em2]
我来回复
您所在位置：