/*
Copyright 2015 Amazon.com, Inc.or its affiliates. All Rights Reserved.

Licensed under the Apache License, Version 2.0 (the "License").
You may not use this file except in compliance with the License.
A copy of the License is located at http ://aws.amazon.com/apache2.0/
or in the "license" file accompanying this file.
This file is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR
CONDITIONS OF ANY KIND, either express or implied.See the License
for the specific language governing permissions and limitations under
the License.
*/
using System;
using System.Collections.Generic;
using System.Linq;
using System.Web;
using System.Web.Mvc;
using System.Threading;
using CSharpMatrixMultiply.Models;
using System.Threading.Tasks;
using System.Runtime.InteropServices;
using System.Text;

namespace CSharpMatrixMultiply.Controllers
{
public class HomeController : Controller
{
enum Algorithms      // this must exactly duplicate enum in C++
{
    CSharp_Basic = 0,
    CSharp_ImprovedSerial = 1,
    CSharp_TPL = 2,
    CPP_Basic = 3,
    CPP_PPL = 4,
    CPP_AMP = 5,
    CPP_AMPTiling = 6
};

delegate float[][] CSharpMatrixMultiply(float[][] A, float[][] B, int N);

const int TESTLOOPS = 3;
const int N = 1024;          // matrix size must be multiple of C++ tilesize

public unsafe ActionResult Index()
{
    int NumAlgorithms = Enum.GetNames(typeof(Algorithms)).Length;

    var rand = new Random();
    double[] durations = new double[NumAlgorithms];

    var TaskResults = new TaskResults(NumAlgorithms);
    TaskResults.Description[0] = "C# Basic Serial (CPU)";
    TaskResults.Description[1] = "C# Improved Serial (CPU)";
    TaskResults.Description[2] = "C# Parallel with TPL (CPU)";
    TaskResults.Description[3] = "C++ Basic Serial (CPU)";
    TaskResults.Description[4] = "C++ Parallel with PPL (CPU)";
    TaskResults.Description[5] = "C++ Parallel with AMP (GPU)";
    TaskResults.Description[6] = "C++ Parallel with AMP Tiling (GPU)";
    TaskResults.Dimension = N;
    TaskResults.NumAlgorithms = NumAlgorithms;
    ViewData["TaskResults"] = TaskResults;

    // According to
    // http://www.heatonresearch.com/content/choosing-best-c-array-type-matrix-multiplication
    // ragged arrays perform better in C# than 2D arrays for matrix multiplication
    float[][] A = CreateRaggedMatrix(N);
    float[][] B = CreateRaggedMatrix(N);
    FillRaggedMatrix(A, N, rand);
    FillRaggedMatrix(B, N, rand);

    // C++ doesn't need ragged arrays for performance, and it's easier to marshall 
    // and process the data as 2D arrays
    float[,] A2 = new float[N, N];
    float[,] B2 = new float[N, N];
    // for comparing results, use the same random data in C++ as in C#
    CopyRaggedMatrixTo2D(A, A2, N);
    CopyRaggedMatrixTo2D(B, B2, N);

    // warm-up AMP and get GPU name before timing
    var sb = new StringBuilder(256);
    CPPWrapper.WarmUpAMP(sb, sb.Capacity);
    TaskResults.AMPDeviceName = sb.ToString();

    //*** Basic C#. Save this original result for future comparisons.
    long start = DateTime.Now.Ticks;
    float[][] original = MatrixMultiplyBasic(A, B, N);
    long stop = DateTime.Now.Ticks;
    durations[0] = (stop - start) / 10000000.0;

    if (!RunCSharpAlgorithm(
        original,
        A,
        B,
        N,
        MatrixMultiplySerial,
        "C# Improved Serial",
        (int)Algorithms.CSharp_ImprovedSerial,
        TaskResults,
        ref durations))
    {
        return PartialView();
    }

    if (!RunCSharpAlgorithm(
        original,
        A,
        B,
        N,
        MatrixMultiplyTPL,
        "C# TPL",
        (int)Algorithms.CSharp_TPL,
        TaskResults,
        ref durations))
    {
        return PartialView();
    }

    if (!RunCPPAlgorithm(original, A2, B2, N, "C++ Basic",
        (int)Algorithms.CPP_Basic, TaskResults, ref durations))
    {
        return PartialView();
    }

    if (!RunCPPAlgorithm(original, A2, B2, N, "C++ PPL",
        (int)Algorithms.CPP_PPL, TaskResults, ref durations))
    {
        return PartialView();
    }

    if (!RunCPPAlgorithm(original, A2, B2, N, "C++ AMP",
        (int)Algorithms.CPP_AMP, TaskResults, ref durations))
    {
        return PartialView();
    }

    if (!RunCPPAlgorithm(original, A2, B2, N, "C++ AMP Tiling",
        (int)Algorithms.CPP_AMPTiling, TaskResults, ref durations))
    {
        return PartialView();
    }

    var slowest = durations.Max();
    var fastest = durations.Min();

    // populate the Model for the HTML table in the View
    for (int k = 0; k < NumAlgorithms; k++)
    {
        TaskResults.Time[k] = string.Format("{0:0.000}", durations[k]);
        TaskResults.RelativeSpeedLabel[k] = string.Format("{0:0.0}X", slowest / durations[k]);
        TaskResults.PercentOfMax[k] = (int)(fastest / durations[k] * 100.0);
    }

    return PartialView();
}

bool RunCSharpAlgorithm(
        float[][] original,
        float[][] A,
        float[][] B,
        int N,
        CSharpMatrixMultiply function,
        string FunctionName,
        int AlgorithmIndex,
        TaskResults results,
        ref double[] durations)
{
    double[] test_durations = new double[TESTLOOPS];

    for (int k = 0; k < TESTLOOPS; k++)
    {
        long start = DateTime.Now.Ticks;
        float[][] C = function(A, B, N);
        test_durations[k] = (DateTime.Now.Ticks - start) / 10000000.0;

        if (!CompareMatrixes(original, C, N))
        {
            results.StatusMessage = "Error verifying " + FunctionName;
            return false;
        }
    }

    durations[AlgorithmIndex] = test_durations.Average();

    return true;
}

unsafe bool RunCPPAlgorithm(
    float[][] original, 
    float[,] A2, 
    float[,] B2, 
    int N,
    string FunctionName,
    int AlgorithmIndex,
    TaskResults results,
    ref double[] durations)
{
    double[] test_durations = new double[TESTLOOPS];

    for (int k = 0; k < TESTLOOPS; k++)
    {
        // allocate memory in C# to simplify marshalling/deallocation
        float[,] C2 = new float[N, N];

        long start = DateTime.Now.Ticks;
        fixed (float* pA2 = &A2[0, 0])
        fixed (float* pB2 = &B2[0, 0])
        fixed (float* pC2 = &C2[0, 0])
        {
            var error = new StringBuilder(1024);   // allocate string memory

            if (!CPPWrapper.CallCPPMatrixMultiply(AlgorithmIndex,
                pA2, pB2, pC2, N, error, error.Capacity))
            {
                results.StatusMessage = error.ToString();
                return false;
            }
        }

        if (!CompareRaggedMatrixTo2D(original, C2, N))
        {
            results.StatusMessage = "Error verifying " + FunctionName;
            return false;
        }

        test_durations[k] = (DateTime.Now.Ticks - start) / 10000000.0;
    }

    durations[AlgorithmIndex] = test_durations.Average();
    return true;
}

// Standard algorithm
float[][] MatrixMultiplyBasic(float[][] A, float[][] B, int N)
{
    float[][] C = CreateRaggedMatrix(N);   // C is the result matrix

    for (int i = 0; i < N; i++)
        for (int j = 0; j < N; j++)
            for (int k = 0; k < N; k++)
                C[i][j] += A[i][k] * B[k][j];

    return C;
}

// This function was developed by Heaton Research and is licensed under the Apache License, Version 2.0,
// available here: https://www.apache.org/licenses/LICENSE-2.0.html
// Improve the basic serial algorithm with optimized index order
float[][] MatrixMultiplySerial(float[][] A, float[][] B, int N)
{
    float[][] C = CreateRaggedMatrix(N);

    // according to http://www.heatonresearch.com/content/choosing-best-c-array-type-matrix-multiplication
    // this ikj index order performs the best for C# matrix multiplication
    for (int i = 0; i < N; i++)
    {
        float[] arowi = A[i];
        float[] crowi = C[i];
        for (int k = 0; k < N; k++)
        {
            float[] browk = B[k];
            float aik = arowi[k];
            for (int j = 0; j < N; j++)
            {
                crowi[j] += aik * browk[j];
            }
        }
    }

    return C;
}

// Parallel algorithm using TPL
float[][] MatrixMultiplyTPL(float[][] A, float[][] B, int N)
{
    float[][] C = CreateRaggedMatrix(N);

    Parallel.For(0, N, i =>
    {
        float[] arowi = A[i];
        float[] crowi = C[i];
        for (int k = 0; k < N; k++)
        {
            float[] browk = B[k];
            float aik = arowi[k];
            for (int j = 0; j < N; j++)
            {
                crowi[j] += aik * browk[j];
            }
        }
    });

    return C;
}

        float[][] CreateRaggedMatrix(int N)
        {
            float[][] result = new float[N][];
            for (int i = 0; i < N; i++)
            {
                result[i] = new float[N];
            }
            return result;
        }

        void FillRaggedMatrix(float[][] ary, int N, Random rand)
        {
            for (int i = 0; i < N; i++)
                for (int j = 0; j < N; j++)
                    ary[i][j] = (float)rand.NextDouble();
        }

        void CopyRaggedMatrixTo2D(float[][] src, float[,] dest, int N)
        {
            for (int i = 0; i < N; i++)
                for (int j = 0; j < N; j++)
                    dest[i, j] = src[i][j];
        }

        // return true if the two input matrixes are identical
        bool CompareMatrixes(float[][] A, float[][] B, int N)
        {
            for (int i = 0; i < N; i++)
                for (int j = 0; j < N; j++)
                    if (A[i][j] != B[i][j])
                        return false;

            return true;
        }

        bool CompareRaggedMatrixTo2D(float[][] A, float[,] B, int N)
        {
            for (int i = 0; i < N; i++)
                for (int j = 0; j < N; j++)
                    if (A[i][j] != B[i,j])
                        return false;

            return true;
        }
    }
}