using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using Azimuth;
using Azimuth.DenseLinearAlgebra;

namespace Azimuth.UnitTests.DenseLinearAlgebra
{
    public class SymmetricToTridiagonalTests : UnitTestSharp.TestFixture
    {
        public void x00()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
                
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                
            });

            var expectedE = new DenseVector(new Scalar[]
            {
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void Basic_x11()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
                { 7, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                7,
            });

            var expectedE = new DenseVector(new Scalar[]
            {
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                { 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void Identity_x22()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 0, },
                { 0, 1, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                1, 1,
            });

            var expectedE = new DenseVector(new Scalar[]
            {
                0,
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 0, },
                { 0, 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void Identity_x33()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 0, 0, },
                { 0, 1, 0, },
                { 0, 0, 1, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                1, 1, 1,
            });

            var expectedE = new DenseVector(new Scalar[]
            {
                0, 0,
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 0, 0, },
                { 0, 1, 0, },
                { 0, 0, 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void x22()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
                { 3, 1, },
                { 1, 2, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                3, 2,
            });

            var expectedE = new DenseVector(new Scalar[]
            {
                1,
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 0, },
                { 0, 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void AlreadyTridiagonal()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 4, 0, },
                { 4, 2, 5, },
                { 0, 5, 3, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                1, 2, 3,
            });

            var expectedE = new DenseVector(new Scalar[]
            {
                -4, -5,
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 0, 0, },
                { 0, -1, 0, },
                { 0, 0, 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void x33()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
                { 1, 4, 6, },
                { 4, 2, 5, },
                { 6, 5, 3, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                -2.34426229508197, 5.34426229508197, 3
            });

            var expectedE = new DenseVector(new Scalar[]
            {
                1.21311475409836, -7.81024967590665
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                { 0.6401843997, -0.7682212796, 0, },
                { -0.7682212796, -0.6401843997, 0, },
                { 0, 0, 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void x44()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
	            { 60, 52, 43, 27, },
	            { 52, 76, 63, 53, },
	            { 43, 63, 99, 74, },
	            { 27, 53, 74, 74, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                14.913916894154, 41.0121980381735, 179.073885067672, 74
            });

            var expectedE = new DenseVector(new Scalar[]
            {
                6.01705801179844, 32.9902137571677, -94.9420876113434
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                { 0.4041311149, -0.8693697984, -0.2843838879, 0, },
                { -0.8070314555, -0.1925457629, -0.5582350392, 0, },
                { 0.4305557708, 0.4551068917, -0.7794225076, 0, },
                { 0, 0, 0, 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }

        public void x88()
        {
            var matrix = new SquareMatrix(new Scalar[,]
            {
	            { 60, 52, 43, 27, 53, 51, 24, 38, },
	            { 52, 76, 63, 53, 51, 65, 25, 46, },
	            { 43, 63, 99, 74, 41, 75, 29, 42, },
	            { 27, 53, 74, 74, 33, 53, 19, 35, },
	            { 53, 51, 41, 33, 55, 51, 25, 45, },
	            { 51, 65, 75, 53, 51, 83, 37, 49, },
	            { 24, 25, 29, 19, 25, 37, 33,  7, },
	            { 38, 46, 42, 35, 45, 49,  7, 61, },
            });

            DenseVector d = new DenseVector(matrix.Size);
            DenseVector e = new DenseVector(matrix.Size);
            SymmetricToTridiagonal.Decompose(ref matrix, ref d, ref e);

            var expectedD = new DenseVector(new Scalar[]
            {
                13.3088855839913,
                8.32949830932554,
                7.55354746370569,
                25.750437284575,
                43.1257691254079,
                37.8123956147073,
                344.119466618287, 
                61
            });

            var expectedE = new DenseVector(new Scalar[]
            {
                -0.450351437931699,
                -5.31171924146175,
                4.18153365016747,
                -16.8297440180236,
                21.4462411825045,
                -60.6274174794168,
                -104.995237987253
            });

            var expectedMatrix = new SquareMatrix(new Scalar[,]
            {
                 { 0.3875899626, 0.0466733281, 0.5232415045, -0.2169406808, -0.6219596097, -0.0944986591, -0.3619211759, 0, },
                 { -0.7095563267, -0.1243091066, -0.2558196642, -0.3851820496, -0.2726194037, -0.0316538546, -0.4381151077, 0, }, 
                 { 0.2756616821, -0.6808401498, 0.0176646259, -0.1266092815, 0.3516243632, 0.4005885317, -0.4000181418, 0, }, 
                 { 0.2087772149, 0.6766121838, -0.1298677779, -0.4251880277, 0.3820954323, 0.2093786502, -0.3333484515, 0, }, 
                 { 0.3395957296, -0.015973091, -0.5810648259, 0.346773705, -0.0759874327, -0.4868953131, -0.4285908662, 0, }, 
                 { -0.3335339754, 0.1699456959, 0.4839082315, 0.5649984785, 0.2938480163, -0.0483621877, -0.4666878321, 0, }, 
                 { 0.0125051764, 0.1785631835, -0.2679429795, 0.4102212543, -0.4208177108, 0.7391559368, -0.0666696903, 0, }, 
                 { 0, 0, 0, 0, 0, 0, 0, 1, },
            });

            CheckEqual(expectedD, d);
            CheckEqual(expectedE, e);
            CheckEqual(expectedMatrix, matrix);
        }
    }
}