数值分析—列主元消去法解线性方程组—FORTRAN程序

数值分析—列主元消去法解线性方程组—FORTRAN程序

program main
implicit none
integer::i,n
!用矩阵（实型动态数组）将线性方程组表示出来
real8,dimension( :,: ),allocatable::a,b
real8,dimension( : ),allocatable::x
print*,‘请输入线性方程组的元数n:’
read*,n
!给动态数组分配内存
allocate(a(n,n+1),b(n,n+1),x(n))
open(111,file=‘A.txt’)
open(333,file=‘ans.txt’)
do i=1,n
read(111,)A(i,1:n+1) !请依次按行输入线性方程组的增广矩阵A
end do
print,“A=”
do i=1,n
print’(1X,5F6.1)’,A(i,1:n+1)
end do
call Gauss_lzy(n,a,b,x)
write(333,*)‘x=’,x
end

subroutine Gauss_lzy(n,a,b,x)
implicit none
integer::i,j,k,imax,n
real8:?,p
real8,dimension(n,n+1)::a,b
real*8,dimension(n)::x
do k=1,n-1
!下面选取每列的最大列主元素
z=abs(a(k,k))
imax=k
do i=k+1,n
if(abs(a(i,k))>z) then
z=abs(a(i,k))
imax=i
endif
enddo
!将最大元素所在的行与第k行进行交换
do j=k,n+1
b(k,j)=a(k,j)
a(k,j)=a(imax,j)
a(imax,j)=b(k,j)
enddo
!对方程组按x1,x2,…,xn的顺序进行依次消元，将增广矩阵a化为上三角矩阵
do i=k+1,n
b(i,k)=a(i,k)/a(k,k)
do j=k+1,n+1
a(i,j)=a(i,j)-b(i,k)a(k,j)
enddo
enddo
enddo
!列主元消去法的第一步消去至此结束
!先计算出x(n)
x(n)=a(n,n+1)/a(n,n)
!下面进行回代计算
do k=n-1,1,-1
p=0
do j=n,k+1,-1
p=p+a(k,j)x(j)
enddo
x(k)=(a(k,n+1)-p)/a(k,k)
enddo
print,‘线性方程组的解为：’
do i=1,n
print,x(i)
enddo
end subroutine
!输入：n=4
!输入：a=1 1 1 2 10
! 3 2 -1 -1 4
! 3 4 3 1 12
! 2 5 4 3 20
!输出：解为 1.0000000000000002
! 1.9999999999999996
! -0.99999999999999967
! 3.9999999999999996*