CLS
  REM ==CYCLIC JACOBI METHOD==
  
  DIM A(10, 10)
  DIM B(10)
  DIM U(10, 10)
  
  REM ==INPUT DATA==
  GOSUB 70
  REM ==INITIALIZE==
  GOSUB 660
  GOSUB 980
  REM ==OUTPUT RESULTS==
  GOSUB 750
  GOTO 3000
  REM ==INPUT DATA==
70 :
  PRINT
  PRINT "CYCLIC JACOBI METHOD"
  PRINT
  PRINT
  PRINT "ENTER SIZE OF PROBLEM"
  INPUT "N=  "; N
  PRINT
  INPUT "NUMBER OF SIGNIFICANT FIGURES "; S1
  PRINT
  PRINT "ENTER UPPER TRIANGLE "
  PRINT
  PRINT "OF MATRIX A , BY COLUMNS"
  PRINT
  FOR J = 1 TO N
  PRINT
  PRINT "ENTER UPPER PART OF COLUMN "; J
  PRINT
  FOR I = 1 TO J
  PRINT "A("; I; ","; J; ") =";
  INPUT A(I, J)
  NEXT I
  PRINT
  NEXT J
  FOR I = 1 TO N
  I1 = I - 1
  FOR J = 1 TO I1
  A(I, J) = A(J, I)
  NEXT J
  NEXT I
  PRINT
  PRINT "ENTER ELEMENTS  OF "
  PRINT "DIAGONAL MATRIX B"
  PRINT
  FOR I = 1 TO N
  PRINT "B("; I; ") = ";
  INPUT B(I)
  NEXT I
  PRINT
  RETURN
  REM ==INITIALIZE TOLERANCE AND MAX NO. OF ROTATIONS==
660 :
  Z = 2 * S1
  T1 = 1 / (10 ^ Z)
  PRINT
  PRINT " TOLERANCE= "; T1
  R = 5 * N * N
  R1 = 0
  T2 = .1
  N1 = N - 1
  RETURN
  REM ==OUTPUT EIGENSOLUTION==
750 :
  PRINT
  PRINT
  PRINT "SOLUTION"
  PRINT "********"
  PRINT
  PRINT "NO. OF ROTATIONS REQUIRED = "; R1
  PRINT
  PRINT
  FOR J = 1 TO N
  PRINT
  PRINT "EIGENVALUE "; J; " IS "; B(J)
  PRINT
  PRINT "ITS EIGENVECTOR IS "
  FOR I = 1 TO N
  PRINT U(I, J)
  NEXT I
  PRINT
  PRINT "PRESS THE ESC KEY TO CONTINUE"
  DO
  LOOP UNTIL INKEY$ = CHR$(27)
  PRINT
  NEXT J
  RETURN
  REM ==EIGENPROBLEM SOLUTION==
980 :
  GOSUB 1130
  REM ==PERFORM ONE CYCLE OF ROTATIONS==
990 :
  GOSUB 1290
  REM ==CHECK TOLERANCE==
  IF X1 < T1 THEN 1110
  REM CHECK NO. OF ROTATIONS
  IF R1 > R THEN 1040
  T2 = .1 * X1
  GOTO 990

1040 :
  PRINT
  PRINT " *** ERROR ***"
  PRINT
  PRINT "NO CONVERGENCE ATTAINED "
  PRINT
  PRINT "WITH "; R1; " ROTATIONS"
  END
1110 :
  GOSUB 1830
  RETURN
 
1130 :
  FOR I = 1 TO N
  FOR J = 1 TO N
  U(I, J) = 0
  NEXT J
  U(I, I) = 1
  NEXT I
  FOR I = 1 TO N
  B1 = SQR(B(I))
  B(I) = 1 / B1
  NEXT I
  FOR I = 1 TO N
  FOR J = 1 TO N
  A(I, J) = B(I) * A(I, J) * B(J)
  NEXT J
  NEXT I
  RETURN

1290 :
  X1 = 0
  FOR K = 1 TO N1
  K1 = K + 1
  FOR L = K1 TO N
  A1 = A(K, K)
  A2 = A(K, L)
  A3 = A(L, L)
  X = A2 * A2 / (A1 * A3)
  REM ==CHECK IF ROTATION IS NEEDED==
  IF X > X1 THEN 1390
  GOTO 1400
1390 :
  X1 = X
1400 :
  IF X < T2 THEN 1800
  R1 = R1 + 1
  REM ==COMPUTE ANGLE=
  IF A1 = A3 THEN 1470
  Z = .5 * (A1 - A3) / A2
  Z1 = 1 + 1 / (Z * Z)
  T = -Z * (1 + SQR(Z1))
  GOTO 1480
1470 :
  T = 1
1480 :
  C = 1 / SQR(1 + T * T)
  S = C * T
  S2 = S * S
  C2 = C * C
  A(K, L) = 0
  REM ==TRANSFORM DIAGONAL ELEMENTS==
  A0 = 2 * A2 * C * S
  A(K, K) = A1 * C2 + A0 + A3 * S2
  A(L, L) = A1 * S2 - A0 + A3 * C2
  REM ==TRANSFORM OFF DIAGONAL ELEMENTS==
  FOR I = 1 TO N
  IF I < K THEN 1600
  IF I > K THEN 1640
  GOTO 1740
1600 :
  A0 = A(I, K)
  A(I, K) = C * A0 + S * A(I, L)
  A(I, L) = -S * A0 + C * A(I, L)
  GOTO 1740
1640 :
  IF I < L THEN 1670
  IF I > L THEN 1710
  GOTO 1740
1670 :
  A0 = A(K, I)
  A(K, I) = C * A0 + S * A(I, L)
  A(I, L) = -S * A0 + C * A(I, L)
  GOTO 1740
1710 :
  A0 = A(K, I)
  A(K, I) = C * A0 + S * A(L, I)
  A(L, I) = -S * A0 + C * A(L, I)
1740 :
  NEXT I
  REM ==TRANSFORM MATRIX U TO GENERATE EIGENVECTORS==
  FOR I = 1 TO N
  U0 = U(I, K)
  U(I, K) = C * U0 + S * U(I, L)
  U(I, L) = -S * U0 + C * U(I, L)
  NEXT I
1800 :
  NEXT L
  NEXT K
  RETURN
  
  REM ==NORMALIZE EIGENVECTORS==
1830 :
  FOR I = 1 TO N
  FOR J = 1 TO N
  U(I, J) = U(I, J) * B(I)
  NEXT J
  NEXT I
  FOR I = 1 TO N
  B(I) = A(I, I)
  NEXT I
  REM ==ORDER EIGENSOLUTION==
  FOR I = 1 TO N1
  I1 = I + 1
  Z = B(I)
  M = I
  FOR J = I1 TO N
  IF Z < B(J) THEN 1990
  Z = B(J)
  M = J
  NEXT J
1990 :
  B(M) = B(I)
  B(I) = Z
  FOR J = 1 TO N
  Z = U(J, I)
  U(J, I) = U(J, M)
  U(J, M) = Z
  NEXT J
  NEXT I
  RETURN
3000 :
  PRINT "END OF PROGRAM"
  END