#pragma once

#ifndef TLIN_MATRIX_H
#define TLIN_MATRIX_H

#include "assert.h"

#include "tlin_basicops.h"

namespace tlin {

//----------------------------------------------------------------------------

/*!
  The tlin::matrix class represents matrices in tlin-compatible algorithms.
  This implementation both serves as a reference to other would-be matrix types,
  specifying the core methods they must provide to work in place of a
  tcg::matrix;
  plus, it provides the naive, universal dense-format of a matrix type.
  Please, observe that a tlin::matrix is stored in column-major order, unlike
  common
  C matrices, since it has to comply with the usual Fortran-style matrices
  supported
  by every BLAS implementation. In practice, this means that the values() array
  stores \b columns in consecutive data blocks.
*/

template <typename T>
class matrix {
  int m_rows;
  int m_cols;

  T *m_entries;

public:
  matrix() : m_rows(0), m_cols(0), m_entries(0) {}
  matrix(int rows, int cols)
      : m_rows(rows), m_cols(cols), m_entries(new T[rows * cols]) {
    memset(m_entries, 0, m_rows * m_cols * sizeof(T));
  }
  ~matrix() { delete[] m_entries; }

  matrix(const matrix &mat)
      : m_rows(mat.m_rows)
      , m_cols(mat.m_cols)
      , m_entries(new T[m_rows * m_cols]) {
    memcpy(m_entries, mat.m_entries, m_rows * m_cols * sizeof(T));
  }

  matrix &operator=(const matrix &mat) {
    if (m_rows != mat.m_rows || m_cols != mat.m_cols) {
      delete[] m_entries;
      m_rows = mat.m_rows, m_cols = mat.m_cols,
      m_entries = new T[m_rows * m_cols];
    }

    memcpy(m_entries, mat.m_entries, m_rows * m_cols * sizeof(T));
    return *this;
  }

  int rows() const { return m_rows; }
  int cols() const { return m_cols; }

  T &at(int row, int col) { return m_entries[m_rows * col + row]; }
  const T &get(int row, int col) const { return m_entries[m_rows * col + row]; }
  const T &operator()(int row, int col) const { return get(row, col); }

  //-----------------------------------------------------------------------------

  //  Dense-specific methods
  matrix(int rows, int cols, T val)
      : m_rows(rows), m_cols(cols), m_entries(new T[rows * cols]) {
    fill(val);
  }

  T *values() { return m_entries; }
  const T *values() const { return m_entries; }

  void fill(const T &val) {
    memset(m_entries, val, m_rows * m_cols * sizeof(T));
  }
};

//-------------------------------------------------------------------------

//  The Standard matrix data type in tlin is double

typedef tlin::matrix<double> mat;

}  // namespace tlin

#endif  // TLIN_MATRIX_H