using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using Azimuth;
using Azimuth.DenseLinearAlgebra;
using Azimuth.Approximation;

namespace Azimuth.UnitTests.Approximation
{
    public class ChebyshevTests : UnitTestSharp.TestFixture
    {
        private static Func<Scalar, Scalar> cheby_4 = x => ((8 * x * x) - 8) * x * x + 1;
        private static Func<Scalar, Scalar> cheby_3 = x => (4 * x * x - 3) * x;
        private static Func<Scalar, Scalar> cheby_2 = x => 2 * x * x - 1;
        private static Func<Scalar, Scalar> cheby_1 = x => x;
        private static Func<Scalar, Scalar> cheby_0 = x => 1;

        public class CtorTests : UnitTestSharp.TestFixture
        {
            public void Chebyshev_0()
            {            
                var window = new Interval(-1, 1);

                var cheby = new Chebyshev(cheby_0, window);

                CheckEqual(new DenseVector(new Scalar[] { 2, }), cheby.Coefficients);
            }

            public void Chebyshev_1()
            {            
                var window = new Interval(-1, 1);

                var cheby = new Chebyshev(cheby_1, window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 1, }), cheby.Coefficients);
            }

            public void Chebyshev_2()
            {
                var window = new Interval(-1, 1);

                var cheby = new Chebyshev(cheby_2, window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 0, 1, }), cheby.Coefficients);
            }

            public void Chebyshev_3()
            {
                var window = new Interval(-1, 1);

                var cheby = new Chebyshev(cheby_3, window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 0, 0, 1, }), cheby.Coefficients);
            }
        
            public void Chebyshev_4()
            {            
                var window = new Interval(-1, 1);

                var cheby = new Chebyshev(cheby_4, window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 0, 0, 0, 1, }), cheby.Coefficients);
            }

            public void x54321()
            {            
                var window = new Interval(-1, 1);

                var cheby = new Chebyshev(x => cheby_4(x) + 2 * cheby_3(x) + 3 * cheby_2(x) + 4 * cheby_1(x) + 2.5 * cheby_0(x), window);

                CheckEqual(new DenseVector(new Scalar[] { 5, 4, 3, 2, 1, }), cheby.Coefficients);
            }

            public void Chebyshev_0_Window()
            {            
                var window = new Interval(3, 7);

                var cheby = new Chebyshev(x => cheby_0(x*2+5), window);

                CheckEqual(new DenseVector(new Scalar[] { 2, }), cheby.Coefficients);
            }

            public void Chebyshev_1_Window()
            {
                var window = new Interval(-1, 4);

                var cheby = new Chebyshev(x => cheby_1((x-1.5)/2.5), window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 1, }), cheby.Coefficients);
            }

            public void Chebyshev_2_Window()
            {
                var window = new Interval(3, 6);

                var cheby = new Chebyshev(x => cheby_2((x-4.5)/1.5), window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 0, 1, }), cheby.Coefficients);
            }

            public void Chebyshev_3_Window()
            {
                var window = new Interval(-2, 11);

                var cheby = new Chebyshev(x => cheby_3((x-4.5)/6.5), window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 0, 0, 1, }), cheby.Coefficients);
            }

            public void Chebyshev_4_Window()
            {            
                var window = new Interval(2.5, 3.5);

                var cheby = new Chebyshev(x => cheby_4((x-3)/0.5), window);

                CheckEqual(new DenseVector(new Scalar[] { 0, 0, 0, 0, 1, }), cheby.Coefficients);
            }

            public void x54321_Window()
            {            
                var window = new Interval(-1, 0);

                var cheby = new Chebyshev(x => 
                {
                    x = ((x + 0.5) * 2);
                    return cheby_4(x) + 2 * cheby_3(x) + 3 * cheby_2(x) + 4 * cheby_1(x) + 2.5 * cheby_0(x);
                }, window);

                CheckEqual(new DenseVector(new Scalar[] { 5, 4, 3, 2, 1, }), cheby.Coefficients);
            }
        }

        public class EvaluateTests : UnitTestSharp.TestFixture
        {
            public void Basic()
            {
                Func<Scalar, Scalar> func = x => Math.Cos(10 * Math.PI * x) + 2 * Math.Sin(3 * Math.PI * x + 1) * x + 1;
                var window = new Interval(0, 2);

                var cheby = new Chebyshev(func, window);

                Scalar t = 0.25;
                var expected = func(t);
                CheckEqual(expected, cheby.Evaluate(t));
            }

            public void Chebyshev_0()
            {
                var window = new Interval(-1, 1);

                var func = cheby_0;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for(int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_1()
            {
                var window = new Interval(-1, 1);

                var func = cheby_1;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_2()
            {
                var window = new Interval(-1, 1);

                var func = cheby_2;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_3()
            {
                var window = new Interval(-1, 1);

                var func = cheby_3;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_4()
            {
                var window = new Interval(-1, 1);

                var func = cheby_4;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void x54321()
            {
                var window = new Interval(-1, 1);

                Func<Scalar, Scalar> func = x => cheby_4(x) + 2 * cheby_3(x) + 3 * cheby_2(x) + 4 * cheby_1(x) + 2.5 * cheby_0(x);
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_0_Window()
            {            
                var window = new Interval(3, 7);

                Func<Scalar, Scalar> func = x => cheby_0(x * 2 + 5);
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_1_Window()
            {
                var window = new Interval(-1, 4);

                Func<Scalar, Scalar> func = x => cheby_1((x - 1.5) / 2.5);
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_2_Window()
            {
                var window = new Interval(3, 6);

                Func<Scalar, Scalar> func = x => cheby_2((x - 4.5) / 1.5);
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_3_Window()
            {
                var window = new Interval(-2, 11);

                Func<Scalar, Scalar> func = x => cheby_3((x - 4.5) / 6.5);
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void Chebyshev_4_Window()
            {            
                var window = new Interval(2.5, 3.5);

                Func<Scalar, Scalar> func = x => cheby_4((x - 3) / 0.5);
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }

            public void x54321_Window()
            {            
                var window = new Interval(-1, 0);

                Func<Scalar, Scalar> func = x =>
                {
                    x = ((x + 0.5) * 2);
                    return cheby_4(x) + 2 * cheby_3(x) + 3 * cheby_2(x) + 4 * cheby_1(x) + 2.5 * cheby_0(x);
                };

                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(func(t), cheby.Evaluate(t));
                }
            }
        }

        public class FirstDerivativeTests : UnitTestSharp.TestFixture
        {
            public void Basic()
            {
                Func<Scalar, Scalar> func = x => Math.Cos(10 * Math.PI * x) + 2 * Math.Sin(3 * Math.PI * x + 1) * x + 1;
                Func<Scalar, Scalar> deri = x => -10 * Math.PI * Math.Sin(10 * Math.PI * x) + 6 * Math.Cos(3 * Math.PI * x + 1) * Math.PI * x + 2 * Math.Sin(3 * Math.PI * x + 1);
                var window = new Interval(0.25, 1.66);

                var cheby = new Chebyshev(func, window);

                for (int i = 0; i < 23; ++i)
                {
                    Scalar t = (i / 22.0) * (window.Max - window.Min) + window.Min;
                    var expected = deri(t);
                    CheckEqual(expected, cheby.FirstDerivative(t));
                }
            }

            public void Chebyshev_0()
            {
                var window = new Interval(-1, 1);

                var func = cheby_0;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(0, cheby.FirstDerivative(t));
                }
            }

            public void Chebyshev_1()
            {
                var window = new Interval(-1, 1);

                var func = cheby_1;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(1, cheby.FirstDerivative(t));
                }
            }

            public void Chebyshev_2()
            {
                var window = new Interval(-1, 1);

                var func = cheby_2;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(4 * t, cheby.FirstDerivative(t));
                }
            }

            public void Chebyshev_3()
            {
                var window = new Interval(-1, 1);

                var func = cheby_3;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(12 * t * t - 3, cheby.FirstDerivative(t));
                }
            }

            public void Chebyshev_4()
            {
                var window = new Interval(-1, 1);

                var func = cheby_4;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual((32 * t * t - 16) * t, cheby.FirstDerivative(t));
                }
            }
        }

        public class SecondDerivativeTests : UnitTestSharp.TestFixture
        {
            public void Basic()
            {
                Func<Scalar, Scalar> func = x => Math.Cos(10 * Math.PI * x) + 2 * Math.Sin(3 * Math.PI * x + 1) * x + 1;
                Func<Scalar, Scalar> deri2 = x =>
                {
                    Scalar w2 = 3 * Math.PI * x + 1;
                    Scalar w1 = 10 * Math.PI * x;
                    return Math.PI * (
                        -Math.PI * (100 * Math.Cos(w1) + 18 * x * Math.Sin(w2))
                        + 12 * Math.Cos(w2)
                    );
                };
                var window = new Interval(0.25, 1.66);

                var cheby = new Chebyshev(func, window);

                for (int i = 1; i < 22; ++i)
                {
                    Scalar t = (i / 22.0) * (window.Max - window.Min) + window.Min;
                    var expected = deri2(t);
                    CheckEqual(expected, cheby.SecondDerivative(t));
                }

                // i = 0, 22 have bad loss of precision, since they're at the edges of the interval, 
                // but check that they're approximately close at least

                {
                    Scalar t = window.Min;
                    var expected = deri2(t);
                    Check(Scalar.FuzzyEquality(expected, cheby.SecondDerivative(t), Math.Sqrt(Scalar.FuzzyEqualityEpsilon)));
                }

                {
                    Scalar t = window.Max;
                    var expected = deri2(t);
                    Check(Scalar.FuzzyEquality(expected, cheby.SecondDerivative(t), Math.Sqrt(Scalar.FuzzyEqualityEpsilon)));
                }
            }

            public void Quadratic_Unit()
            {
                Func<Scalar, Scalar> func = x => x * (4 * x + 3) + 1;
                Func<Scalar, Scalar> deri2 = x =>
                {
                    return 8;
                };
                var window = new Interval(-1, 1);

                var cheby = new Chebyshev(func, window);

                for (int i = 0; i < 23; ++i)
                {
                    Scalar t = (i / 22.0) * (window.Max - window.Min) + window.Min;
                    var expected = deri2(t);
                    CheckEqual(expected, cheby.SecondDerivative(t));
                }
            }

            public void Quadratic_OffUnit()
            {
                Func<Scalar, Scalar> func = x => x * (4 * x + 3) + 1;
                Func<Scalar, Scalar> deri2 = x =>
                {
                    return 8;
                };
                var window = new Interval(0.25, 1.75);

                var cheby = new Chebyshev(func, window);

                for (int i = 0; i < 23; ++i)
                {
                    Scalar t = (i / 22.0) * (window.Max - window.Min) + window.Min;
                    var expected = deri2(t);
                    CheckEqual(expected, cheby.SecondDerivative(t));
                }
            }

            public void Chebyshev_0()
            {
                var window = new Interval(-1, 1);

                var func = cheby_0;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(0, cheby.SecondDerivative(t));
                }
            }

            public void Chebyshev_1()
            {
                var window = new Interval(-1, 1);

                var func = cheby_1;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(0, cheby.SecondDerivative(t));
                }
            }

            public void Chebyshev_2()
            {
                var window = new Interval(-1, 1);

                var func = cheby_2;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(4, cheby.SecondDerivative(t));
                }
            }

            public void Chebyshev_3()
            {
                var window = new Interval(-1, 1);

                var func = cheby_3;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(24 * t, cheby.SecondDerivative(t));
                }
            }

            public void Chebyshev_4()
            {
                var window = new Interval(-1, 1);

                var func = cheby_4;
                var cheby = new Chebyshev(func, window);

                int max = 13;
                for (int i = 0; i < max; ++i)
                {
                    Scalar t = (Scalar)i / (Scalar)max * (window.Max - window.Min) + window.Min;
                    CheckEqual(96 * t * t - 16, cheby.SecondDerivative(t));
                }
            }
        }

        public class FindRealRootsInIntervalOfCreation : UnitTestSharp.TestFixture
        {
            public void OddDoubleRoot()
            {
                Func<Scalar, Scalar> function = (x) => -3 + 5 * Math.Sin(2.9562880943 + -1.4410157268 * x) + -7 * Math.Sin(-2.9562880943 + 2.4410157268 * x);

                var cheby = new Chebyshev(function, new Interval(0, 15));

                var roots = cheby.FindRealRootsInIntervalOfCreation();

                var expected = new Scalar[]
                {
                    0.0332997141421392,
                    1.2954858927,
                    4.9876994145,
                    6.2498855931,
                    8.0993275495,
                    8.8723870107,
                    10.173281710,
                    11.136292104,
                    12.968129194,
                    14.192681105,
                };

                CheckEqual(expected, roots);
            }
            
            public void FlatFunction_0()
            {
                Func<Scalar, Scalar> function = (x) => 0;

                var cheby = new Chebyshev(function, new Interval(10, 20));

                var roots = cheby.FindRealRootsInIntervalOfCreation();

                var expected = new Scalar[]
                {
                    Scalar.PositiveInfinity,
                };

                CheckEqual(expected, roots);
            }

            public void Empty()
            {
                var cheby = new Chebyshev(new DenseVector(0), new Interval(-1, 1));

                CheckEqual(new Scalar[] { Scalar.PositiveInfinity, }, cheby.FindRealRootsInIntervalOfCreation());
            }

            public void ExtraneousTerms()
            {
                var cheby = new Chebyshev(t =>
                {
                    return -10 + -10 - (-3.9355296627 + (-998.0302953003 + 250 * t) * t);
                }, new Interval(0.0161616033, 0.0166666667));

                var roots = cheby.FindRealRootsInIntervalOfCreation();

                CheckEqual(new Scalar[] { 0.0161616033 }, roots);
            }
        }

        public class Convolve : UnitTestSharp.TestFixture
        {
            public void SinCos()
            {
                var chebCos = new Chebyshev(x =>
                {
                    return Math.Cos(x);
                }, new Interval(-1, 1));

                var chebSin = new Chebyshev(x =>
                {
                    return Math.Sin(x+1);
                }, new Interval(-1, 1));

                var combined = new Chebyshev(x =>
                {
                    return Math.Cos(x) * Math.Sin(x+1);
                }, new Interval(-1, 1));

                var convolved1 = Chebyshev.Convolve(chebCos, chebSin);
                var convolved2 = Chebyshev.Convolve(chebSin, chebCos);

                CheckEqual(combined.Coefficients, convolved1.Coefficients);
                CheckEqual(combined.Coefficients, convolved2.Coefficients);
            }

            public void MismatchInterval()
            {
                var chebCos = new Chebyshev(x =>
                {
                    return Math.Cos(x);
                }, new Interval(-1, 1));

                var chebSin = new Chebyshev(x =>
                {
                    return Math.Sin(x + 1);
                }, new Interval(0, 1));

                CheckThrow(typeof(Exception));
                var convolved1 = Chebyshev.Convolve(chebCos, chebSin);
            }
        }

        public class Add : UnitTestSharp.TestFixture
        {
            public void SinCos()
            {
                var chebCos = new Chebyshev(x =>
                {
                    return Math.Cos(x);
                }, new Interval(-1, 1));

                var chebSin = new Chebyshev(x =>
                {
                    return Math.Sin(x + 1);
                }, new Interval(-1, 1));

                var combined = new Chebyshev(x =>
                {
                    return Math.Cos(x) + Math.Sin(x + 1);
                }, new Interval(-1, 1));

                var sum1 = Chebyshev.Add(chebCos, chebSin);
                var sum2 = Chebyshev.Add(chebSin, chebCos);

                CheckEqual(combined.Coefficients, sum1.Coefficients);
                CheckEqual(combined.Coefficients, sum2.Coefficients);
            }

            public void MismatchInterval()
            {
                var chebCos = new Chebyshev(x =>
                {
                    return Math.Cos(x);
                }, new Interval(-1, 1));

                var chebSin = new Chebyshev(x =>
                {
                    return Math.Sin(x + 1);
                }, new Interval(0, 1));

                CheckThrow(typeof(Exception));
                var convolved1 = Chebyshev.Add(chebCos, chebSin);
            }
        }

        public class FindHarmonic : UnitTestSharp.TestFixture
        {
            public void ReasonableTest_UnitInterval()
            {
                /// A * cos(w * t + theta) + B * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return 3 * Math.Cos(2 * x + 1) + 4 * Math.Sin(2 * x + 1);
                }, new Interval(-1, 1));

                var harmonic = Chebyshev.FindHarmonic(2, 1, 3, 4, new Interval(-1, 1));

                CheckEqual(cheby.Coefficients, harmonic.Coefficients);
            }

            public void ReasonableTest_CrazyInterval()
            {
                /// A * cos(w * t + theta) + B * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return 3 * Math.Cos(2 * x + 1) + 4 * Math.Sin(2 * x + 1);
                }, new Interval(17, 30));

                var harmonic = Chebyshev.FindHarmonic(2, 1, 3, 4, new Interval(17, 30));

                CheckEqual(cheby.Coefficients, harmonic.Coefficients);
            }

            public void SimpleCosine()
            {
                /// A * cos(w * t + theta) + B * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return Math.Cos(x);
                }, new Interval(-1, 1));

                var harmonic = Chebyshev.FindHarmonic(1, 0, 1, 0, new Interval(-1, 1));

                for (int i = 0; i < 100; ++i)
                {
                    var x = -1 + i/50.0;
                    CheckEqual(cheby.Evaluate(x), harmonic.Evaluate(x));
                }

                CheckEqual(cheby, harmonic);
            }

            public void SimpleSine()
            {
                /// A * cos(w * t + theta) + B * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return Math.Sin(x);
                }, new Interval(-1, 1));

                var harmonic = Chebyshev.FindHarmonic(1, 0, 0, 1, new Interval(-1, 1));

                for (int i = 0; i < 100; ++i)
                {
                    var x = -1 + i / 50.0;
                    CheckEqual(cheby.Evaluate(x), harmonic.Evaluate(x));
                }

                CheckEqual(cheby, harmonic);
            }

            public void SimpleCosine_Phase()
            {
                /// A * cos(w * t + theta) + B * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return Math.Cos(x+1);
                }, new Interval(-1, 1));

                var harmonic = Chebyshev.FindHarmonic(1, 1, 1, 0, new Interval(-1, 1));

                CheckEqual(cheby, harmonic);
            }

            public void SimpleSine_Phase()
            {
                /// A * cos(w * t + theta) + B * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return Math.Sin(x + 1);
                }, new Interval(-1, 1));

                var harmonic = Chebyshev.FindHarmonic(1, 1, 0, 1, new Interval(-1, 1));

                CheckEqual(cheby, harmonic);
            }

            public void SimpleCosSin_Phase()
            {
                /// A * cos(w * t + theta) + B * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return Math.Sin(x + 1) + Math.Cos(x + 1);
                }, new Interval(-1, 1));

                var harmonic = Chebyshev.FindHarmonic(1, 1, 1, 1, new Interval(-1, 1));

                CheckEqual(cheby, harmonic);
            }

            public void SimpleCosSin_PI_Phase()
            {
                var cheby = new Chebyshev(x =>
                {
                    return Math.Sin(Math.PI * x + 1) + Math.Cos(Math.PI * x + 1);
                }, new Interval(-1, 1));

                var harmonic = Chebyshev.FindHarmonic(Math.PI, 1, 1, 1, new Interval(-1, 1));
                
                CheckEqual(cheby, harmonic);
            }
        }

        public class ConvolveHarmonic : UnitTestSharp.TestFixture
        {
            public void Flat()
            {
                var interval = new Interval(-1, 1);
                var cheby = new Chebyshev(x =>
                {
                    return (5 * x + 3);
                }, interval);

                var X = new Chebyshev(x =>
                {
                    return 5 * x + 3;
                }, interval);

                var Y = new Chebyshev(x =>
                {
                    return (2 * x * x + 3 * x + 1);
                }, interval);

                var convolved = Chebyshev.HarmonicConvolution(0, 0, X, Y);

                CheckEqual(cheby, convolved);
            }

            public void Simple_Cos()
            {
                var interval = new Interval(-1, 1);
                var cheby = new Chebyshev(x =>
                {
                    return 5 * Math.Cos(Math.PI * x);
                }, interval);

                var X = new Chebyshev(x =>
                {
                    return 5;
                }, interval);

                var Y = new Chebyshev(x =>
                {
                    return 0;
                }, interval);

                var convolved = Chebyshev.HarmonicConvolution(Math.PI, 0, X, Y);

                for (int i = 0; i < 100; ++i)
                {
                    var x = interval.Min + i / 100.0 * (interval.Max - interval.Min);

                    CheckEqual(cheby.Evaluate(x), convolved.Evaluate(x));
                }

                CheckEqual(cheby, convolved);
            }

            public void Simple_Sin()
            {
                var interval = new Interval(-1, 1);
                var cheby = new Chebyshev(x =>
                {
                    return 5 * Math.Sin(Math.PI * x);
                }, interval);

                var X = new Chebyshev(x =>
                {
                    return 0;
                }, interval);

                var Y = new Chebyshev(x =>
                {
                    return 5;
                }, interval);

                var convolved = Chebyshev.HarmonicConvolution(Math.PI, 0, X, Y);

                for (int i = 0; i < 100; ++i)
                {
                    var x = interval.Min + i / 100.0 * (interval.Max - interval.Min);

                    CheckEqual(cheby.Evaluate(x), convolved.Evaluate(x));
                }

                CheckEqual(cheby, convolved);
            }

            public void Simple_SinCos()
            {
                var interval = new Interval(-1, 1);
                var cheby = new Chebyshev(x =>
                {
                    return 3 * Math.Sin(Math.PI * x) + 4 * Math.Cos(Math.PI * x);
                }, interval);

                var X = new Chebyshev(x =>
                {
                    return 4;
                }, interval);

                var Y = new Chebyshev(x =>
                {
                    return 3;
                }, interval);

                var convolved = Chebyshev.HarmonicConvolution(Math.PI, 0, X, Y);

                for (int i = 0; i < 100; ++i)
                {
                    var x = interval.Min + i / 100.0 * (interval.Max - interval.Min);

                    CheckEqual(cheby.Evaluate(x), convolved.Evaluate(x));
                }

                CheckEqual(cheby, convolved);
            }

            public void Simple_Cos_Phase()
            {
                var interval = new Interval(-1, 1);
                var cheby = new Chebyshev(x =>
                {
                    return Math.Cos(Math.PI * x + 1);
                }, interval);

                var X = new Chebyshev(x =>
                {
                    return 1;
                }, interval);

                var Y = new Chebyshev(x =>
                {
                    return 0;
                }, interval);

                var convolved = Chebyshev.HarmonicConvolution(Math.PI, 1, X, Y);

                for (int i = 0; i < 100; ++i)
                {
                    var x = interval.Min + i / 100.0 * (interval.Max - interval.Min);

                    //CheckEqual(cheby.Evaluate(x), convolved.Evaluate(x));
                }

                CheckEqual(cheby, convolved);
            }

            public void Simple_Sin_Phase()
            {
                var interval = new Interval(-1, 1);
                var cheby = new Chebyshev(x =>
                {
                    return 4 * Math.Sin(Math.PI * x + 1);
                }, interval);

                var X = new Chebyshev(x =>
                {
                    return 0;
                }, interval);

                var Y = new Chebyshev(x =>
                {
                    return 4;
                }, interval);

                var convolved = Chebyshev.HarmonicConvolution(Math.PI, 1, X, Y);

                for (int i = 0; i < 100; ++i)
                {
                    var x = interval.Min + i / 100.0 * (interval.Max - interval.Min);

                    CheckEqual(cheby.Evaluate(x), convolved.Evaluate(x));
                }

                // Slight differences end up compounding
                // CheckEqual(cheby, convolved);
            }

            public void ReasonableTest_UnitInterval()
            {
                var interval = new Interval(-1, 1);
                /// (5 * t + 3) * cos(w * t + theta) + (2t^2 + 3*t + 1) * sin(w * t + theta)
                var cheby = new Chebyshev(x =>
                {
                    return (5 * x + 3) * Math.Cos(2 * x + 1) + (2*x*x + 3*x + 1) * Math.Sin(2 * x + 1);
                }, interval);

                var X = new Chebyshev(x =>
                {
                    return 5 * x + 3;
                }, interval);

                var Y = new Chebyshev(x =>
                {
                    return (2*x*x + 3*x + 1);
                }, interval);

                var convolved = Chebyshev.HarmonicConvolution(2, 1, X, Y);

                for (int i = 0; i < 100; ++i)
                {
                    var x = interval.Min + i / 100.0 * (interval.Max - interval.Min);

                    // CheckEqual(cheby.Evaluate(x), convolved.Evaluate(x));
                }

                CheckEqual(cheby, convolved);
            }
        }

        public class CalculateFromMononomialBasis : UnitTestSharp.TestFixture
        {
            public void x3()
            {
                var expected = new DenseVector(new Scalar[]
                {
                    2, 1, 1,
                });

                var monomial = new DenseVector(new Scalar[] { 0, 1, 2, });
                
                var actual = Chebyshev.FromPolynomial(
                    monomial
                );

                CheckEqual(expected, actual);
            }

            public void x10()
            {
                var coef = new DenseVector(new Scalar[]
                {
                    1, 2, 3, 4, 5, 6, 7, 8, 9,
                });

                var expected = new DenseVector(new Scalar[]
                {
                    18.046875,
                    13.125,
                    11.21875,
                    5.5,
                    3.90625,
                    1.25,
                    0.78125,
                    0.125,
                    0.0703125,
                });

                var actual = Chebyshev.FromPolynomial(coef);

                CheckEqual(expected, actual);
            }
        }
    }
}