Deprecated: The each() function is deprecated. This message will be suppressed on further calls in /home/zhenxiangba/zhenxiangba.com/public_html/phproxy-improved-master/index.php on line 456
! 文字コードはutf-8
! 自由形式のfortranで書いています。
! ガウスの消去法で連立1次方程式a(i,j)*x(j)=b(i)を解くサブルーチン
! gausu の使用例。
implicit none
integer:: i,j,n
real:: a(100,100),x(100),b(100) !配列の次元を100まで切るけど使うのはnまで
!例1
n=3
a(1,1)=1.; a(1,2)=2.; a(1,3)=3.
a(2,1)=1.; a(2,2)=3.; a(2,3)=3.
a(3,1)=2.; a(3,2)=5.; a(3,3)=7.
b(1)=14.
b(2)=16.
b(3)=33.
!答え:x(1)=1.,x(2)=2.,x(3)=3.
!
!例2
!n=4
!a(1,1)= 4.;a(1,2)=-1.;a(1,3)=-1.;a(1,4)= 0.
!a(2,1)=-1.;a(2,2)= 4.;a(2,3)= 0.;a(2,4)=-1.
!a(3,1)=-1.;a(3,2)= 0.;a(3,3)= 4.;a(3,4)=-1.
!a(4,1)= 0.;a(4,2)=-1.;a(4,3)=-1.;a(4,4)= 4.
!b(1)=8.;b(2)=8.;b(3)=-16.;b(4)=8.
! 例2答え:x(1)=2.,x(2)=3.,x(3)=-3.,x(4)=2.
!
!例3
!n=5
!a(1,1)=1.;a(1,2)=1.;a(1,3)=1.;a(1,4)=0.;a(1,5)=0.
!a(2,1)=0.;a(2,2)=1.;a(2,3)=1.;a(2,4)=1.;a(2,5)=0.
!a(3,1)=0.;a(3,2)=0.;a(3,3)=1.;a(3,4)=1.;a(3,5)=1.
!a(4,1)=1.;a(4,2)=0.;a(4,3)=0.;a(4,4)=1.;a(4,5)=1.
!a(5,1)=1.;a(5,2)=1.;a(5,3)=0.;a(5,4)=0.;a(5,5)=1.
!b(1)=2.;b(2)=1.;b(3)=3.;b(4)=3.;b(5)=3.
!答え:x(1)=1.,x(2)=0.,x(3)=1.,x(4)=0.,x(5)=2.
!
!
print*,'a='
do i=1,n
print*, (a(i,j),j=1,n)
end do
!
print*,'b='
print*, (b(i),i=1,n)
! プログラムの中でa(n,n)*x(n)=b(n)を解きたければ、
! 上のようにa(n,n),b(n)を定義してから、
! 以下のように書けばサブルーチンで解いた答えがx(n)に入る
call gausu(n,a,x,b)
print*,'x='
print*, (x(i),i=1,n)
end
!
!以下はサブルーチン、プログラムの末尾の方に入れておく。
subroutine gausu(n,a,x,b)
implicit none
real:: a(100,100),x(100),b(100),akjxj
integer::i,j,k,n
!
!ガウスの消去法の参考としたのは、
!名取亮「すうがくぶっくす12 線形計算」(朝倉書店)p.10-15
!
do k=1,n-1 !a(k,k)を消去
do i=k+1,n !k+1行からn行まで
do j=k+1,n !k+1列からn列まで
a(i,j)=a(i,j)-a(k,j)*a(i,k)/a(k,k)
end do
b(i)=b(i)-b(k)*a(i,k)/a(k,k)
end do
end do
!
! 後退代入
x(n)=b(n)/a(n,n)
do k=n,1,-1
akjxj=0.d0
do j=k+1,n
akjxj=akjxj+a(k,j)*x(j)
end do
x(k)=(b(k)-akjxj)/a(k,k)
end do
!
return
end
!
!