using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Data;
using System.Drawing;
using System.Linq;
using System.Text;
using System.Windows.Forms;
using System.Threading;
using Azimuth.DenseLinearAlgebra;

namespace Azimuth.PerformanceTests
{
    public partial class MainForm : Form
    {
        bool QRDecompostion { get; set; }
        int progress;
        Random rand = new Random();
        
        public MainForm()
        {
            InitializeComponent();
        }

        private delegate void ProgressIncrement();

        private void StartTest_Click(object sender, EventArgs e)
        {
            Thread testThread = null;
            int workToBeDone = 0;
            progress = 0;

            if (QRDecompositionRadio.Checked)
            {
                testThread = RunQRTest(out workToBeDone);
            }
            else if (LinearLeastSquaresRadio.Checked)
            {
                testThread = RunLLSTest(out workToBeDone);
            }
            else if (BisectionRadio.Checked)
            {
                testThread = RunBisectionTest(out workToBeDone);
            }
            else if (radioButton3.Checked)
            {
                testThread = RunLLSInequalitiesTest(out workToBeDone);
            }
            else if (matMul.Checked)
            {
                testThread = RunMatMul(out workToBeDone);
            }
            else if (sparseMMCheckbox.Checked)
            {
                testThread = RunSparseMultiply(out workToBeDone);
            }
            else if (choleskyCheck.Checked)
            {
                testThread = RunCholesky(out workToBeDone);
            }
            else if (DantzigButton.Checked)
            {
                testThread = DantzigTest(out workToBeDone);
            }
            else if (radixsortButton.Checked)
            {
                testThread = RunRadixsort(out workToBeDone);
            }
            else if (quicksortButton.Checked)
            {
                testThread = RunQuicksort(out workToBeDone);
            }
            else
            {
                throw new Exception("Not yet implemented.");
            }

            progressBar.Minimum = 0;
            progressBar.Maximum = workToBeDone;

            var uiThread = new Thread(delegate()
            {
                var startTime = System.DateTime.Now;

                while (testThread.IsAlive)
                {
                    progressBar.Invoke((MethodInvoker) delegate()
                    {
                        progressBar.Value = Math.Min(progressBar.Maximum, Math.Max(progressBar.Minimum, Thread.VolatileRead(ref progress)));
                        progressLabel.Text = Thread.VolatileRead(ref progress).ToString() + "/" + workToBeDone.ToString();
                    });
                    Thread.Sleep(100);
                }

                progressBar.Invoke((MethodInvoker) delegate()
                {
                    progressBar.Value = progressBar.Maximum;
                    progressLabel.Text = Thread.VolatileRead(ref progress).ToString() + "/" + workToBeDone.ToString();
                });

                System.TimeSpan Duration = System.DateTime.Now - startTime;
                string operationCount = (!Scalar.IsCountingFlops) ? "" : "\nOperation count: " + Scalar.FLOPs;
                System.Windows.Forms.MessageBox.Show("Done.  Took " + Duration.TotalSeconds + " seconds." + operationCount);
            });

            testThread.IsBackground = true;
            testThread.Start();

            uiThread.IsBackground = true;
            uiThread.Start();
        }

        private Azimuth.DenseLinearAlgebra.RectangularMatrix RandomMatrix(int size)
        {
            Random rand = new Random();

            var matrix = new Azimuth.DenseLinearAlgebra.RectangularMatrix(size, size);
            for (int i = 0; i < matrix.RowCount; ++i)
            {
                for (int j = 0; j < matrix.ColumnCount; ++j)
                {
                    matrix[i, j] = (double)(rand.Next(-100, 100));
                }
            }

            return matrix;
        }

        private Azimuth.SparseLinearAlgebra.SparseRectangularMatrix RandomSparseMatrix(
            int rows, int cols, int nonZeroEntries)
        {
            Random rand = new Random();

            var matrix = new Azimuth.SparseLinearAlgebra.SparseRectangularMatrix(rows, cols);

            for(int i = 0; i < nonZeroEntries; ++i)
            {
                int col = rand.Next(cols);
                int row = rand.Next(rows);

                matrix[row, col] = rand.NextDouble() * 100000;
            }

            return matrix;
        }

        private Thread RunMatMul(out int doneValue)
        {
            /*
             * Here's Matlab/Octave code for comparison:
                rand1 = rand(1000,1000);
                rand2 = rand(1000,1000);
                t = cputime; mat*mat2; e = cputime - t
             * 
             * Octave w/out SSE took about 1.2 secs on a i7 @ 2.67 Ghz
             * (that is, each CPU tick can produce ~1 Flop,
             *  and we need 2 * N^3 Flops total).
             * */
            doneValue = 1000;

            var matrix1 = RandomMatrix(doneValue);
            var matrix2 = RandomMatrix(doneValue);

            var workerThread = new Thread(delegate()
            {
                var mat = matrix1.MatrixMultiply(matrix2, ref progress);
            });

            return workerThread;
        }

        private Thread RunQuicksort(out int doneValue)
        {
            doneValue = 4000;
            //doneValue *= doneValue;

            var list = new List<uint>(doneValue);

            for (int i = 0; i < doneValue; ++i)
            {
                list.Add((uint)rand.Next());
            }

            var workerThread = new Thread(delegate()
            {
                list.Sort();
            });

            return workerThread;
        }

        private Thread RunRadixsort(out int doneValue)
        {
            doneValue = 4000;
            //doneValue *= doneValue;

            var list = new uint[doneValue];

            for (int i = 0; i < doneValue; ++i)
            {
                list[i] = ((uint)rand.Next());
            }

            var workerThread = new Thread(delegate()
            {
                var newList = Azimuth.DataStructures.RadixSort.Sort(list.ToArray());
            });

            return workerThread;
        }

        private Thread DantzigTest(out int doneValue)
        {
            int size = 3;
            int iterCount = 1000000;
            doneValue = size;

            int unsolvable = 0;

            var mat1 = RandomPSDSymmetricMatrix(size);

            //var rand = new Random();

            var q = new DenseLinearAlgebra.DenseVector(size);
            for (int j = 0; j < size; ++j)
            {
                q[j] = -1;// (-2 * rand.NextDouble()) * 100000;
            }

            var workerThread = new Thread(delegate()
            {
                for (int i = 0; i < iterCount; ++i)
                {
                    if (iterCount > 1)
                    {
                        mat1 = RandomPSDSymmetricMatrix(size);
                    }

                    DenseVector z;

                    bool solved = DantzigsAlgorithm.SolveLCP(out z, mat1, q, ref progress);

                    if (solved)
                    {
                        var w = z.PreMatrixMultiplyBy(mat1);
                        w.AssignAdd(q);

                        for (int j = 0; j < w.Length; ++j)
                        {
                            if (w[j] == 0)
                            {
                                w[j] = 0;
                            }
                        }

                        for (int j = 0; j < w.Length; ++j)
                        {
                            if (w[j].Value < -Math.Sqrt(Scalar.FuzzyEqualityEpsilon))
                            {
                                throw new Exception("Not complementary!");
                            }
                            if (z[j].Value < 0)
                            {
                                throw new Exception("Not complementary!");
                            }
                            else if (Math.Abs(w[j] * z[j]) < Math.Sqrt(Scalar.FuzzyEqualityEpsilon))
                            {

                            }
                            else
                            {
                                throw new Exception("Not complementary!");
                            }
                        }
                    }
                    else
                    {
                        unsolvable++;
                    }

                    System.Threading.Interlocked.Increment(ref progress);
                }

                MessageBox.Show(String.Format("{0} unsolvable", (Scalar)unsolvable / (Scalar)iterCount));
            });

            return workerThread;
        }

        private Thread RunLLSInequalitiesTest(out int doneValue)
        {
            const int problemSize = 512;
            doneValue = problemSize * 3;
            progress = 0;

            Random rand = new Random();

            var matrix = RandomMatrix(problemSize);

            var constraints = RandomMatrix(problemSize);

            var B = RandomMatrix(problemSize);

            var b = new DenseLinearAlgebra.DenseVector(problemSize);
            for(int i = 0; i < problemSize; ++i)
            {
                b[i] = rand.NextDouble() * 100000;
            }

            var h = new DenseLinearAlgebra.DenseVector(problemSize);
            for (int i = 0; i < problemSize; ++i)
            {
                h[i] = rand.NextDouble() * 100000;
            }

            var workerThread = new Thread(delegate()
            {
                var lls = new Azimuth.DenseLinearAlgebra.LinearLeastSquares(matrix, ref progress);

                DenseLinearAlgebra.DenseVector vector;
                lls.SolveSystemWithConstraints(out vector, b, constraints, h, ref progress);
            });

            return workerThread;
        }

        private Thread RunQRTest(out int doneValue)
        {            
            doneValue = 512;
            progress = 0;

            var matrix = RandomMatrix(doneValue);
            
            var workerThread = new Thread(delegate()
            {
                Azimuth.DenseLinearAlgebra.QRDecomposition qr = 
                    new Azimuth.DenseLinearAlgebra.QRDecomposition(matrix, ref progress);
            });

            return workerThread;
        }

        private Thread RunLLSTest(out int doneValue)
        {
            int size = 512;
            doneValue = size * 2;
            progress = 0;

            var matrix = RandomMatrix(size);
            
            var workerThread = new Thread(delegate()
            {
                var lls = new Azimuth.DenseLinearAlgebra.LinearLeastSquares(matrix, ref progress);
            });

            return workerThread;
        }

        private Thread RunSparseMultiply(out int doneValue)
        {
            doneValue = 1000;

            var sparseMatrix1 = RandomSparseMatrix(doneValue, doneValue, doneValue);
            var sparseMatrix2 = RandomSparseMatrix(doneValue, doneValue, doneValue);

            var workerThread = new Thread(delegate()
            {
                sparseMatrix1.MatrixMultiplyIntoSparse(sparseMatrix2, ref progress);
            });

            return workerThread;
        }

        private Thread RunCholesky(out int doneValue)
        {
            doneValue = 1000;
            var mat1 = RandomPSDSymmetricMatrix(doneValue);

            var workerThread = new Thread(delegate()
            {
                var chol = new Azimuth.DenseLinearAlgebra.CholeskyDecomposition(
                    mat1, ref progress);
            });

            return workerThread;
        }

        private SquareMatrix RandomPSDSymmetricMatrix(int size)
        {
            var mat1 = RandomMatrix(size);

            return new SquareMatrix(mat1.MatrixMultiply(mat1.Transpose()));
        }

        private Thread RunBisectionTest(out int doneValue)
        {
            doneValue = 10000;
            progress = 0;

            var workerThread = new Thread(delegate()
            {
                for (progress = 0; progress < 10000; ++progress)
                {
                    Azimuth.RootFinding.Bisection.FindRoot(x => (x - 4) * (x + 4) * (x - 3) * (x + 2) * (x + 5),
                        new Interval(-10100, 45678));
                }
            });

            return workerThread;
        }
    }
}