Метод вращения Якоби (ошибка в цикле)

[crawter]

Всем доброго времени суток. Цель программы - поиск собственных чисел матрицы 4х4 методом вращения Якоби.
Суть алгоритма :
находим наибольший не диагональный элемент
находим fi
заполняем матрицу вращения U и U транспонированную
ищем собственные числа.
Первые итерации у меня считает правильно, но потом программка зацикливается.

 class Systems_of_the_linear_equations
    {
        private int n { get; set; }

        private double[,] A;

   public Systems_of_the_linear_equations()
        {
            My_Func2();
        }
 public void My_Func2()
        {
            n = 4;

            A = new double[n, n];

            A[0, 0] = 7;   A[1, 0] = 0;   A[2, 0] = 7;  A[3, 0] = 10;
            A[0, 1] = 0;   A[1, 1] = 7;   A[2, 1] = 1;   A[3, 1] = 1;
            A[0, 2] = 7;  A[1, 2] = 1;   A[2, 2] = 1;   A[3, 2] = 2;
            A[0, 3] = 10;   A[1, 3] = 1;   A[2, 3] = 2;   A[3, 3] = 8;


        }
 public void Jacobi_Own_Num()
        {
            double fi = 0, o = 1;
            int i1 = 0, j1 = 0, iter = 0;
            double[,] U = new double[n, n];
            double[,] Ut = new double[n, n];
            double[,] temp = new double[n, n];
            double[,] Ak = new double[n, n ];

            while (Math.Abs(o) > 0.0001) 
            {

                double max = -1;

                for (int i = 0; i < n; i++)
                {
                    for (int j = 0; j < n; j++)
                    {
                        if ((i != j) && (A[j, j] >= max))
                        {
                            max = Math.Abs(A[i, j]);
                            i1 = i; j1 = j;
                        }
                        max = A[i1, j1];
                    }
                }
              

                if (A[i1, i1] != A[j1, j1])
                {
                    fi = Math.Atan(2 * max / (A[i1, i1] - A[j1, j1])) / 2;
                }
                else
                {
                    fi = Math.PI / 4;
                }


                for (int i = 0; i < n; i++)
                {
                    for (int j = 0; j < n; j++)
                    {
                        if (i != j)
                        {
                            U[i, j] = 0;
                            Ut[i, j] = 0;
                        }
                        else
                        {
                            U[i, j] = 1;
                            Ut[i, j] = 1;
                        }
                    }
                }

                U[i1, i1] = Math.Cos(fi); U[i1, j1] = -(Math.Sin(fi)); U[j1, i1] = Math.Sin(fi); U[j1, j1] = Math.Cos(fi);
                Ut[i1, j1] = Math.Sin(fi); Ut[j1, i1] = -(Math.Sin(fi)); Ut[i1, i1] = Math.Cos(fi); Ut[j1, j1] = Math.Cos(fi);
                


                double a = 0;

                for (int i = 0; i < n; i++)
                {
                    for (int j = 0; j < n; j++)
                    {
                        for (int k = 0; k < n; k++)
                        {
                            temp[i, j] = a + Ut[i, k] * A[k, j];
                            a = temp[i, j];
                        }
                        a = 0;
                    }
                }

                a = 0;
                for (int i = 0; i < n; i++)
                {
                    for (int j = 0; j < n; j++)
                    {
                        for (int k = 0; k < n; k++)
                        {
                            Ak[i, j] = a + temp[i, k] * U[k, j];
                            a = Ak[i, j];
                        }
                        a = 0;
                    }
                }

                o = 0;
                for (int i = 0; i < n; i++)
                {
                    for (int j = 0; j < n; j++)
                    {
                        if (i != j)
                        {
                            o += Ak[i, j] * Ak[i, j];
                            A[i, j] = Ak[i, j];
                        }
                    }
                }

                iter++;
                Console.WriteLine("Iteration {0} \n\n", iter);
                for (int i = 0; i < n; i++)
                {
                    
                    for (int j = 0; j < n; j++)
                    {
                        Console.Write("{0} ", Ak[i, j]);
                    }
                    Console.WriteLine("\n");
                }
                Console.Write("\n\nOwn number {0} ", o);

            }

            Console.WriteLine("\n");
            for (int i = 0; i < n; i++)
            {
                Console.WriteLine("L{0} = {1}\n", i+1, Ak[i, i]);
            }

        }
}