ガウスの消去
ガウスの消去法のプログラムを作ったのですがうまく動きません
どこが間違っているのでしょう
ちなみに連立方程式を解くプログラムです
* the numerical solution of linear equation
* by gauss method
parameter(ll=10)
dimension a(ll,ll+1), x(ll)
read(5,100) n
100 format(i3)
do 150 i = 1, n
150 read(5,200) (a(i,j), j=1, n+1)
200 format(11f4.1)
write(6,300)
300 format(' ', 10x, 'coefficient')
do 10 i = 1, n
write(6,400) (a(i,j), j=1, n+1)
400 format(' ', 11 (5x, f4.1))
10 continue
esp=10.0e-19
call rgaule(a, x, ll, n, esp, ipivot)
if (ipivot. eq. 1) go to 20
write(6,500)
500 format(' ', 10x, 'the pivot is little '/
1 'so the solution is singular')
go to 110
20 write(6,600)
600 format(/' ',2x, 8hsolution)
write (6,700) (i, x(i), i = 1, n)
700 format (' ',5x, 'x (', i2, ' ) =', 2x, e14.7)
110 stop
end
subroutine rgaule(a, x, ll, n, esp, ipivot)
dimension a(ll, ll+1), x(ll), ln(100)
ipivot = 1
* the order of x
do 5 i = 1, n
ln(i) = i
5 continue
do 100 m = 1, n-1
* the selection of pivot
amax = 0.0
do 10 i = m, n
do 20 j = m, n
if (amax. ge. abs(a(i,j))) go to 20
amax = abs(a(i, j))
irow = i
icolum = j
20 continue
10 continue
if (amax. le. eps) go to 70
if (m. eq. irow) go to 22
* the exchange of row
do 27 l = m, n+1
swap = a(irow, l)
a(irow, l) = a(m, l)
a(m, l) = swap
27 continue
22 if (m. eq. icolum) go to 30
* the exchange of colum
do 25 l =1, n
swap = a(l, icolum)
a(l, icolum) = a(l,m)
a(l, m) = swap
25 continue
* the exchange of x
iswap = ln(m)
ln(m) = ln(icolum)
ln(icolum) = iswap
* gaussian elimination
30 do 35 i = m+1, n
do 37 j = m+1, n+1
a(i,j) = a(i,j) - a(i,m) * a(m,j) / a(m,m)
37 continue
35 continue
100 continue
if (abs(a(n,n)). le. eps) go to 70
* back subsititution
x(n) = a(n, n+1) / a(n,n)
kk = ln(n)
a(n, kk) = x(n)
do 50 i = n-1, 1, -1
k = n-i
x(i) = 0.0
do 52 j = 1, k
ll = i + j
x(i) = a(i, ll) * x(ll) + x(i)
52 continue
x(i) = (a(i, n+1) - x(i)) / a(i,i)
kk = ln(i)
a(n, kk) = x(i)
50 continue
do 60 i = 1, n
x(i) = a(n, i)
60 continue
return
70 ipivot = 0
return
end
