#include <cmath>
#include <vector>
#include <string>
#include <ostream>
#include <iostream>
#include <cstdlib>

using namespace std;

void MatrixError(const string&);



class Matrix
{
    private:  // These are the internal data of the Matrix
    
    int rows;
    int columns;
    vector<double> elements;
    
    public:   // Public interface methods

    // Matrix A(2,2)
    Matrix(int r,int c)
       {
          SetSize(r,c);
       }

    // Matrix B(A);
    Matrix(const Matrix& A)   // Make a new matrix which is a copy of A.
       {
          SetSize(A.Rows(), A.Columns());   
          for (int r=0; r<rows; r++)
            for (int c=0; c<columns; c++)
               Element(r,c)=A(r,c);
       }

    // Common code to (re)define size of matrix, and allocate space.
    void SetSize(int r, int c)
    {
       rows=r;
       columns=c;
       // Allocate space for elements, setting them to zero.
       elements.resize(rows*columns, 0.0);
    }

    // Return number of rows in matrix
    int Rows() const
    {
       return rows;
    }
    
    // Return number of columns in matrix
    int Columns() const
    {
       return columns;
    }


    // Declare methods that are later defined in Matrix.cpp
    
    bool hasElement(int r, int c) const;                // True if r,c is an element in the matrix.
    int index(int r, int c) const;                      // Get linear index corresponding to r,c

    const double&  Element(int r,int c) const;          // Return the constant value of the element at position r,c
    double& Element(int r,int c);                       // Return a modifiable reference to the element at position r,c
    void SetElement(int r,int c, double e);             // Set the element at position r,c to value e
    
    void Identity();                                    // Turn matrix into identiy matrix (only for symmetric matrices!)
    
    // Row operations to help implementing Gauss elimination
    void AddScaledRow(int r1, int r2,double scale);     // Add row r2 *scale to row r1
    void ScaleRow(int r,double scale);                  // Multiply row r by c
    void SwapRows(int r1,int r2);                       // Swap rows r1 and r2


    // Stuff for you to implement    
    void GaussJordan(Matrix& Inverse);                  // Apply Gauss Jordan turning the matrix into the identity, and return its inverse.

/*
    void Transpose();                                   // Transposes the matrix
    double Trace();                                     // Returns the trace of a matrix

    // For 2x2 matrices
    void EigenValues(double& e1, double& e2 );          // Calculate the eigenvalues    

    // Somewhat complicated:
    double Determinant();                               // Determinant (hint create a Constructor that takes a matrix and removes a column and row)
*/



    // Defining a lot of operators
    const double& operator()(int r,int c) const;        // return a constant value of the element A(r,c).    E.g. cout << A(0,0);     The value is unchanged
    double& operator()(int r,int c);                    // return modifiable reference to the element,       E.g. A(r,c) = 0.0;       The reference is changed.

    // All these operators modify the internal state of A, but does not change B.
    Matrix& operator=(const Matrix &B);                 // A=B   => A.operator=(B)
    Matrix& operator+=(const Matrix &B);                // A+=B  => A.operator+=(B)
    Matrix& operator-=(const Matrix &B);                // A-=B  => A.operator-=(B)
    Matrix& operator*=(const Matrix &B);                // A*=B  => A.operator*=(B)
    Matrix& operator*=(const double &x);                // A*=x  => A.operator*=(x)

    // These are binary operators which return a new matrix, the original matrices are not modified.
    const Matrix operator+(const Matrix &B) const;      // A+B  => A.operator+(B)
    const Matrix operator-(const Matrix &B) const;      // A-B  => A.operator-(B)
    const Matrix operator*(const Matrix &B) const;      // A*B  => A.operator*(B)
    const Matrix operator*(double) const;               // A*x  => A.operator*(double x)
    
    // x*A is special, since the left hand value is a double, and we need to define it with two arguments as follows:
    friend  Matrix operator*(const double& x, const Matrix &M)
    {
       Matrix copy(M);
       return copy*x;  // Calls existing Matrix::operator*(double);
    }

    friend ostream& operator<<(ostream&, Matrix&);
};