Fortran
program main
c main program to illustrate eigenproblem solvers
c ndim array dimension, ndim = 6 in this example
c n number of equations, n
c a coefficient matrix, A(i,j)
c x eigenvector, x(i)
c y intermediate vector, y(i)
c norm specifies unity component of eigenvecor
c iter number of iterations allowed
c tol convergence tolerance
c shift amount by which eigenvalue is shifted
c iw intermediate results output flag: 0 no, 1 yes
dimension a(6,6) ,x(6) ,y(6)
data ndim = 6, n = 3,norm =3 , iter = 50, tol = 1.e-6, shift = 0.0, iw=1 / 6, 3, 3, 50,1.e-6,0.0,1/
c data ndim, n,norm, iter, tol, shift, iw/6, 3, 3, 50, 1.e-6, 13.870584, 2/
data (a(i,1),i=l,3) / 8.0, -2.0, -2.0 /
data (a(i,2),i=i,3) /-2.0, 4.0, -2.0 /
data (a(i,3),i=i,3) /-2.0, -2.0, 13.0 /
data (x(i),i=l,3) / 1.0, 1.0, 1.0 /
write (6,1000)
do i=1, n
write (6,1010) i, (a(i,j),j=l,n),x(i)
end do
if (shift.gt.0.0) then
write (6,1005) shift
do i=i,n
a(i, i)=a(i,i) -shift
write (6,1010) i, (a(i,j),j=i,n)
end do
end if
call power ( ndim, n, a, x, y, norm, i ter , t ol , shift, i w , k, ev2 )
write (6,1020)
write (6,1010) k, ev2, (x(i),i=l,n)
stop
1000 format (' The power method'/' '/' A and x(O) '/' ')
1005 format (' '/' A shifted by shift = ',flO. 6/' ')
1010 format (ix, i$, 6f12.6)
1020 format (' '/' k lambda and eigenvector components'/' ')
end
subroutine power (ndim, n, a, x, y, norm, i ter, tol , shift, iw, k, ev2 )
c the direct power method
dimensiOn a(ndim, ndim),x(ndim),y(ndim)
evl = 0.0
if (iw. eq.1) write (6,1000)
if (iw. eq.1) write (6,1010) k, evl, (x(i),i=l,n)
do k=1,iter
c calculate y(i)
do i=l,n
y(i) =0.0
do j=l,n
y(i) =y(i) +a (i, j) *x(j)
end do
end do
c calculate lambda and x(i)
ev2 =y(norm)
if (abs (ev2) .le.1. Oe-3) then
write (6,1020) ev2
return
else
do i=l,n
x (i) =y(i )/ev2
end do
end if
if (iw. eq. 1) write (6,1010) k, ev2, (x(i),i=i,n)
c check for convergence
if (abs(ev2-evl) .Ie. tol) then
if (shift.ne.O.O) then
evl =ev2
ev2=ev2+shi ft
write (6, 1040) evl, ev2
end if
return
else
evl =ev2
end if
end do
write (6, 1030)
return
1000 format (' '/' k lambda and eigenvector components'/' ')
I010 format (ix, i3,6f12.6)
1020 forma (' '/' lambda = ',elO.2, ' approaching zero, stop')
1030 format (' '/' Iterations did not converge, stop']
1040 format (' '/' lambda shifted =',f12.6,' and lambda =',fi2.6)
end
C PROGRAM
#include <stdio.h>
#include <conio.h>
#include <math.h>
int main()
{
double A[3][3] = {8,-2,-2,-2,4,-2,-2,-2,13};
double X[6], Y[6];
int x,y,ndim = 6,n = 3,norm = 3, iter = 50, tol = 1, iw=1;
double shift = 0.0,tol = 1e-6;
for (x=0;x<=2;x++)
{
for (y=0;y<=2;y++)
{
printf("%7.4f\t",A[x][y]);
}
printf("\n");
}
getch();
}

New Topic/Question
This topic is locked



MultiQuote








|