Roots of polynomial , eigenvalue method

using System;
using System.IO;
using System.Globalization;
namespace NamespaceName
{
    public class ClassName
    {
        const int maxiter = 1000;
        const double eps = 1e-12;
        public static double hypot(double a,double b)
        {
            double absa,absb;
            absa = Math.Abs(a);
            absb = Math.Abs(b);
            if(absa > absb) return absa * Math.Sqrt(1+(double)(absb/absa)*(double)(absb/absa));
            else return (absb == 0?0: absb * Math.Sqrt(1 + (double)(absa/absb)*(double)(absa/absb)));

        }
        public static void printMatrix(int m,int n, double[,] A)
        {
            NumberFormatInfo nfi = new NumberFormatInfo();
            nfi.NumberDecimalDigits = 12;
            for(int i = 0;i < m;i++)
            {
                for(int j = 0; j < n; j++)
                {
                    Console.Write("{0} , ",A[i,j].ToString("N",nfi));
                }
                Console.WriteLine();
            }
            Console.WriteLine();
        }
        public static void QR_Givens(int m,int n,double [,] A,double[,] Q)
        {
            for(int i = 0; i < m; i++)
                for(int j = 0; j < m; j++)
                    Q[i,j] = (i == j ? 1 : 0);
            int min = (m < n ? m : n);
            for(int i = 0; i < min; i++)
            {
                for(int j = i + 1; j < m; j++)
                {
                    if(A[j,i] != 0)
                    {
                        double r = hypot(A[i,i],A[j,i]);
                        double c = (double) (A[i,i]/r);
                        double s = (double) (A[j,i]/r);
                        for(int k = 0;k < n; k++)
                        {
                            double temp = A[i,k];
                            A[i,k] = c * A[i,k] + s * A[j,k];
                            A[j,k] = -s * temp + c * A[j,k];
                        }
                        for(int k = 0;k < m;k++)
                        {
                            double temp = Q[k,i];
                            Q[k,i] = c * Q[k,i] + s * Q[k,j];
                            Q[k,j] = -s * temp + c * Q[k,j];
                        }
                    }
                }
            }

        }
        public static void copyMatrix(double[,] A,double[,] B,int m,int n)
        {
            for(int i = 0;i < m;i++)
                for(int j = 0;j < n;j++)
                    B[i,j] = A[i,j];
        }
        public static void multiplyMatrix(double[,] A,double[,] B,double[,] C,int m,int n,int p)
        {
            for(int i = 0;i < m;i++)
                for(int j = 0;j < p;j++)
                {
                    double sum = 0;
                    for(int k = 0;k < n;k++)
                        sum += A[i,k] * B[k,j];
                    C[i,j] = sum;
                }
        }
        public static double rayleigh(int n, double[,] A,double[] q)
        {
            double norm = 0;
            double sum;
            double[] v = new double[n];
            for(int i = 0;i < n;i++)
                norm += q[i] * q[i];
            for(int i = 0;i < n;i++)
            {
                sum = 0;
                for(int j = 0;j < n;j++)
                    sum += A[i,j]*q[j];
                v[i] = sum;
            }
            double r = 0;
            for(int i = 0;i < n;i++)
                r += q[i]*v[i];
            r /= (double)norm;
            return r;
        }
        public static void Main(string[] args)
        {
            char esc;
            int n;
            double[,] A,Q,R;
            double[] v;
            double[] a;
            double r;
            Random rnd = new Random();
            NumberFormatInfo nfi = new NumberFormatInfo();
            nfi.NumberDecimalDigits = 12;


            using(StreamWriter sw = new StreamWriter("polyroots.txt",true))
            {
                do
                {
                Console.WriteLine("Podaj stopien wielomianu");
                int.TryParse(Console.ReadLine(),out n);
                A = new double[n,n];
                Q = new double[n,n];
                R = new double[n,n];
                v = new double[n];
                a = new double[n + 1];
                for(int i = n;i >= 0;i--)
                {
                    Console.Write("Podaj a[{0}]= ", i);
                    double.TryParse(Console.ReadLine(),out a[i]);
                }
                for(int i = n;i >= 0;i--)
                    if(a[i] < 0)
                        sw.WriteLine("-{0}x^{1} ",(-a[i]).ToString("N",nfi),i);
                    else
                        sw.WriteLine("+{0}x^{1} ",a[i].ToString("N",nfi),i);
                    sw.WriteLine();
                for(int i=0;i<n;i++)
                    A[0,i] = (double)(-a[n-i-1]/a[n]);
                for(int i = 1;i < n;i++)
                    for(int j = 0;j<n;j++)
                        A[i,j] = (i == j+1)?1:0;
                printMatrix(n,n,A);
                for(int i = 0;i < n;i++)
                {
                    for(int j = 0; j < n; j++)
                        sw.Write("{0} ",A[i,j].ToString("N",nfi));
                    sw.WriteLine();
                }
                sw.WriteLine();
                for(int i = 0;i < n;i++)
                    v[i] = i + rnd.NextDouble();
                for(int i = 0;i < maxiter;i++)
                {
                    r = rayleigh(n,A,v);
                    for(int j = 0;j < n;j++)
                        A[j,j] -= r;
                    QR_Givens(n,n,A,Q);
                    copyMatrix(A,R,n,n);
                    multiplyMatrix(R,Q,A,n,n,n);
                    for(int j = 0;j < n;j++)
                        A[j,j] += r;
                    for(int j = 0;j < n;j++)
                        v[j] = Q[n-1,j];
                }
                printMatrix(n,n,A);
                for(int i = 0;i < n;i++)
                {
                    for(int j = 0; j < n; j++)
                        sw.Write("{0} ",A[i,j].ToString("N",nfi));
                    sw.WriteLine();
                }
                sw.WriteLine();
                Console.WriteLine("Pierwiastki danego rownania wielomianowego to: ");
                sw.WriteLine("Pierwiastki danego rownania wielomianowego to: ");
                int k = 0;
                while(k<n)
                {
                    if(k + 1 < n && Math.Abs(A[k+1,k])>eps)
                    {
                        double p = 0.5*(A[k,k]+A[k+1,k+1]);
                        double q = A[k,k] * A[k+1,k+1] - A[k,k+1] * A[k + 1,k];
                        double d = q - p * p;
                        Console.WriteLine("x[{0}]={1}-{2}i",k,p.ToString("N",nfi),Math.Sqrt(d).ToString("N",nfi));
                        Console.WriteLine("x[{0}]={1}+{2}i",k+1,p.ToString("N",nfi),Math.Sqrt(d).ToString("N",nfi));
                        sw.WriteLine("x[{0}]={1}-{2}i",k,p.ToString("N",nfi),Math.Sqrt(d).ToString("N",nfi));
                        sw.WriteLine("x[{0}]={1}+{2}i",k+1,p.ToString("N",nfi),Math.Sqrt(d).ToString("N",nfi));
                        k += 2;
                    }
                    else
                    {
                        Console.WriteLine("x[{0}]={1}",k,A[k,k].ToString("N",nfi));
                        sw.WriteLine("x[{0}]={1}",k,A[k,k].ToString("N",nfi));
                        k++;
                    }
                }
                esc = (char) Console.ReadKey().Key;
            }
            while(esc != (char)ConsoleKey.Escape);
        }
        }
    }
}