c c Numerical Analysis: c The Mathematics of Scientific Computing c D.R. Kincaid & E.W. Cheney c Brooks/Cole Publ., 1990 c c Section 4.3 c c Example of gaussian elimination with scaled row pivoting c Solving U(trans)z=c and L(trans)Py=z c c c file: yaec.f c parameter (n = 3) dimension a(n,n),s(n),c(n),y(n),z(n) integer p(n) data (a(1,j),j=1,n) /2.0,3.0,-6.0/ data (a(2,j),j=1,n) /1.0,-6.0,8.0/ data (a(3,j),j=1,n) /3.0,-2.0,1.0/ data (c(i),i=1,n) /6.0,-5.0,3.0/ c print * print *,' Gaussian elimination with scaled row pivoting' print *,' Section 4.3, Kincaid-Cheney' print * c call gauss(n,a,s,p) c do 3 j=1,n sum = c(j) do 2 i=1,j-1 sum = sum - a(p(i),j)*z(i) 2 continue z(j) = sum/a(p(j),j) 3 continue do 5 j=n,1,-1 sum = z(j) do 4 i=j+1,n sum = sum - a(p(i),j)*y(p(i)) 4 continue y(p(j)) = sum 5 continue c do 6 i=1,n print 7,i,y(i) 6 continue c 7 format (3x,'y(',i2,') =',e13.6) stop end c subroutine gauss(n,a,s,p) c c gaussian elimination c dimension a(n,n),s(n) integer p(n) real max c do 9 i=1,n p(i) = i smax = 0.0 do 8 j=1,n smax = max(smax,abs(a(i,j))) 8 continue s(i) = smax 9 continue c do 13 k=1,n-1 rmax = 0.0 do 10 i=k,n r = abs(a(p(i),k))/s(p(i)) if (r .gt. rmax) then j = i rmax = r endif 10 continue c pk = p(j) p(j) = p(k) p(k) = pk c do 12 i=k+1,n z = a(p(i),k)/a(p(k),k) a(p(i),k) = z do 11 j=k+1,n a(p(i),j) = a(p(i),j) - z*a(p(k),j) 11 continue 12 continue 13 continue c return end