#include <iostream>
#include <iomanip>
#include <fstream>
using namespace std;

class Matrix
{
public:
    Matrix(ifstream &is);
    Matrix();
    ~Matrix();
    Matrix *multiply(Matrix &rhs) const;
    void print(ostream &os) const;
    void clone(Matrix &src);
private:
    static double **allocate(int nr, int nc);
    static void deallocate(int nr, double **data);
    
    int m_nRows;
    int m_nColumns;
    double **m_data;
};

Matrix::Matrix()
    : m_nRows(0), m_nColumns(0), m_data(0)
{
}

double **Matrix::allocate(int nr, int nc)
{
    int ir;
    double **mat;
    mat = new double *[nr];
    for (ir=0; ir<nr; ir++)
        mat[ir] = new double[nc];
    return mat;
}

Matrix::Matrix(ifstream &is)
{
    int ir, ic;
    is >> m_nRows >> m_nColumns;
    m_data = allocate(m_nRows, m_nColumns);
    for (ir=0; ir<m_nRows; ir++)
        for (ic=0; ic<m_nColumns; ic++)
            is >> m_data[ir][ic];
}

void Matrix::deallocate(int nr, double **data)
{
    int ir;
    for (ir=0; ir<nr; ir++)
        delete[] data[ir];
    delete[] data;
}

Matrix::~Matrix()
{
    deallocate(m_nRows, m_data);
}

Matrix *Matrix::multiply(Matrix &rhs) const
{
    if (m_nColumns != rhs.m_nRows)
        return 0;

    Matrix *c = new Matrix;
    c->m_nRows = m_nRows;
    c->m_nColumns = rhs.m_nColumns;
    c->m_data = allocate(c->m_nRows, c->m_nColumns);

    int ir, ic, ip;
    for (ir=0; ir<c->m_nRows; ir++)
        for (ic=0; ic<c->m_nColumns; ic++)
        {
            c->m_data[ir][ic] = 0;
            for (ip=0; ip<m_nColumns; ip++)
                c->m_data[ir][ic] += m_data[ir][ip] * rhs.m_data[ip][ic];
        }
    return c;
}

void Matrix::print(ostream &os) const
{
    int ir, ic;
    os << "     # rows=" << m_nRows << "   # columns=" << m_nColumns << endl;
    for (ir=0; ir<m_nRows; ir++)
    {
        os << "              ";
        for (ic=0; ic<m_nColumns; ic++)
            os << setw(10) << m_data[ir][ic];
        os << endl;
    }
}

void Matrix::clone(Matrix &src)
{
    if (m_data)
        deallocate(m_nRows, m_data);
    m_nRows = src.m_nRows;
    m_nColumns = src.m_nColumns;
    m_data = allocate(m_nRows, m_nColumns);

    int ir, ic;
    for (ir=0; ir<m_nRows; ir++)
        for (ic=0; ic<m_nColumns; ic++)
            m_data[ir][ic] = src.m_data[ir][ic];
}

void main()
{
    ifstream infile("matrix.dat");
    Matrix matA(infile);
    Matrix matB(infile);
    
    Matrix *matC=0;
    matC = matA.multiply(matB);
    
    cout << "The first matrix A is " << endl;
    matA.print(cout);
    cout << "The second matrix B is " << endl;
    matB.print(cout);
    cout << "The product of A x B is " << endl;
    matC->print(cout);
    
    Matrix matD;
    matD.clone(*matC);
    cout << "The copy of matrix C is " << endl;
    matD.print(cout);
    delete matC;
}