[Fortran] 纯文本查看 复制代码
program mytp21
implicit none
integer i,j,n
double precision h
double precision,allocatable ::x(:,:),b(:),a(:),c(:),f(:),u(:)
write(6,*)'input un nombre n'
read(5,*)n
allocate(x(n,n),a(n),b(n),c(n),f(n),u(n))
h=1.d0/n
x(1,1)=1.d0/h+1.d0*h/3
f(1)=1.d0
do i=2,n
f(i)=0
!write(6,'("(",i3,")=",d15.8)')i,f(i)
enddo
do i=1,n
do j=1,n
if(i==j)then
x(i+1,j+1)=2.d0*x(1,1)
elseif(i==j+1)then
x(i,j)=-1.d0/h+1.d0*h/6.d0
elseif(j==i+1) then
x(i,j)=-1.d0/h+1.d0*h/6.d0
else
x(i,j)=0
endif
enddo
enddo
do i=1,n
b(i)=x(i,i)
!write(6,'("(",i3,")=",d15.8)')i,b(i)
enddo
do i=2,n
a(i)=x(i,i-1)
! write(6,'("(",i3,")=",d15.8)')i,a(i)
enddo
do i=1,n-1
c(i)=x(i,i+1)
! write(6,'("(",i3,")=",d15.8)')i,c(i)
end do
a(1)=0
c(n)=0
u(n)=0
open(30, file='le solution de u(x).txt')
call tridag(a,b,c,f,u,n)
do i=1,n
write(6,'("(",i3,")=",d15.8)')i,u(i)
write(30,*)'le solution de u(x)',u(i)
end do
close(30)
pause
contains
subroutine tridag(a,b,c,f,u,n)
implicit none
integer n
double precision gam(n),a(n),b(n),c(n),u(n),f(n),bet
integer j,n
if (b(1)==0.) pause 'b(1)=0 dans tridag'
bet=b(1)
u(1)=f(1)/bet
do j=2,n
gam(j)=c(j-1)/bet
bet=b(j)-a(j)*gam(j)
if (bet==0.) pause 'bet=0 dans tridag'
u(j)=(f(j)-a(j)*u(j-1))/bet
end do
do j=n-1,1,-1
u(j)=u(j)-gam(j+1)*u(j+1)
end do
end subroutine tridag
end program mytp21
! --------------------------------------------------
! Silverfrost FTN95 for Microsoft Visual Studio
! Free Format FTN95 Source File
! --------------------------------------------------
