#include "SparseMatrix.h"
#include <iostream>
#include <cassert>
#include <vector>
using namespace std;

SparseMatrix::SparseMatrix(void)
{
}


SparseMatrix::~SparseMatrix(void)
{
}


void SparseMatrix::unitTest(void)
{

    SparseMatrix A, At, B, Bt, C;
    A.setValue(0,0,0,0);
    At.setValue(0,0,0,0);
    B.setValue(0,0,0,0);
    Bt.setValue(0,0,0,0);
    C.setValue(0,0,0,0);

    SparseMatrix C0, C1, C2;
    struct Entry C0A[] = {{0,0,0},{0,0,3},{0,1,-1},{0,2,5},{1,0,1},{1,2,2},{3,1,2},{4,2,-1}};
    C0.setValue(5,3,7,C0A);
    struct Entry C1A[] = {{0,0,0},{0,0,3},{0,1,2},{0,2,6},{1,0,3},{1,1,7},{1,2,5},{2,2,-2}};
    C1.setValue(3,3,7,C1A);
    struct Entry C2A[] = {{0,0,0},{0,0,-1},{0,3,2},{1,0,9},{1,1,3},{2,0,-13},
                                  {2,1,-4},{2,3,6},{2,4,2},{3,0,-1},{3,1,1},
                                  {3,3,8},{4,0,10},{4,1,4},{4,4,-2}};
    C2.setValue(5,5,14,C2A);

    A.read();
    A.m_A[A.m_nElems+1].row = -1;

    B.read();
    B.m_A[B.m_nElems+1].row = -1;

    if ((A.m_nRows != B.m_nRows)||(A.m_nCols != B.m_nCols))
    {
        cout << "Matrix A and B are not of the same dimension!\n";
        return;
    }

    cout << "\n";
    cout << "Matrix A:\n";
    A.print();

    cout << "Matrix B:\n";
    B.print();

    Bt.transpose(B);
    assert(Bt.equal(C0));
    cout << "Transpose Matrix B:\n";
    Bt.print();

    At.transpose(A);

    A.multiply(Bt);
    assert(A.equal(C1));

    cout << "Matrix A * Matrix B':\n";
    A.print();

    cout << "Transpose Matrix A:\n";
    At.print();

    At.multiply(B);
    assert(At.equal(C2));

    cout << "Matrix A' * Matrix B:\n";
    At.print();

    cin.get();cin.get();
}

void SparseMatrix::setValue(int nRows, int nCols, int nElems, struct Entry *A)
{
    m_nRows = nRows;
    m_nCols = nCols;
    m_nElems = nElems;
    m_A = A;
}


void SparseMatrix::read(void)
{
    int i;
    cin >> m_nRows >> m_nCols >> m_nElems;
    delete[] m_A;
    m_A = new Entry[m_nElems+2];
    for (i=1; i<=m_nElems; i++)
        cin >> m_A[i].row >> m_A[i].col >> m_A[i].value;
}


void SparseMatrix::print(void) const
{
    int i;
    cout << m_nRows << " * " << m_nCols << " matrix with "
         << m_nElems << " non-zero item:\n";
    for (i=1; i<=m_nElems; i++)
        cout << '<' << m_A[i].row << ',' << m_A[i].col << ':' << m_A[i].value << ">\n";
    cout << "====================\n";
}


void SparseMatrix::transpose(const SparseMatrix &A)
{
    int i, iB=0;
    vector<int> start(A.m_nCols+1);

    m_nRows = A.m_nCols; m_nCols = A.m_nRows; m_nElems = A.m_nElems;
    delete[] m_A;
    m_A = new Entry[m_nElems+2];

    for (start[0]=1, i=1; i<=A.m_nCols; i++)
        start[i] = 0;
    for (i=1; i<=A.m_nElems; i++)
        start[A.m_A[i].col+1]++; // counting element's in a row
    for (i=1; i<=A.m_nCols; i++)
        start[i] += start[i-1]; // accumulating

    for (i=1; i<=A.m_nElems; i++)
    {
        iB = start[A.m_A[i].col]++;
        m_A[iB].row = A.m_A[i].col;
        m_A[iB].col = A.m_A[i].row;
        m_A[iB].value = A.m_A[i].value;
    }
    m_A[m_nElems+1].row = -1; // extra entry for multiply to test end of B1's row
}

bool SparseMatrix::equal(const SparseMatrix &B) const
{
    int i;
    if ((m_nRows!=B.m_nRows)||(m_nCols!=B.m_nCols)||(m_nElems!=B.m_nElems))
        return false;
    for (i=1; i<=m_nElems; i++)
        if ((m_A[i].row!=B.m_A[i].row)||(m_A[i].col!=B.m_A[i].col)||
            (m_A[i].value!=B.m_A[i].value))
            return false;
    return true;
}

void SparseMatrix::multiply(const SparseMatrix &B)
{
    int i, iA, iC, iA1, iB1, iA1row, iB1row, sum;

    SparseMatrix B1;
    B1.setValue(0,0,0,0);
    B1.transpose(B);

    vector<Entry> result;
    iA = iC = 1;
    while (iA <= m_nElems)
    {
        sum = 0;
        iA1 = iA; iB1 = 1;
        iA1row = m_A[iA1].row; iB1row = B1.m_A[iB1].row;
        while (iB1 <= B1.m_nElems)
        {
            if (m_A[iA1].col == B1.m_A[iB1].col)
            {
                sum += m_A[iA1].value * B1.m_A[iB1].value;
                iA1++; iB1++;
            }
            else if (m_A[iA1].col > B1.m_A[iB1].col)
                iB1++;
            else // if (m_A[iA1].col < B1.m_A[iB1].col)
                iA1++;
            if (m_A[iA1].row != iA1row)
                while (B1.m_A[iB1].row == iB1row) iB1++;
            if (B1.m_A[iB1].row != iB1row)
            {
                if (iB1 <= B1.m_nElems)
                    iA1 = iA;
                else
                    while (m_A[iA1].row == iA1row) iA1++;
                if (sum!=0)
                {
                    Entry e = {iA1row, iB1row, sum};
                    result.push_back(e);
                    sum = 0;
                }
                iB1row = B1.m_A[iB1].row;
            }
        }
        iA = iA1;
    }
    delete[] B1.m_A;

    m_nElems = result.size();
    m_A = new Entry[m_nElems+2];
    for (i=0; i<m_nElems; i++)
        m_A[i+1] = result[i];
    m_A[0].row = m_A[0].col = m_A[0].value = 0;
    m_nCols = B.m_nCols;
}