// Matrice.h
#pragma once
#include <iostream>
#include <vector>
#include <iomanip>
using namespace std;

template<class T>
class Matrice
{
private:
	int ordin;
	vector<vector<T>> elemente;

	Matrice<T> getMinor(int, int);

	void interschimbareLinii(int, int, int);
	void multiplicareLinie(int, int, T);
	void operatieLinii(int, int, int, T);
	int liberColoana(int, int, int);

public:
	Matrice() : Matrice(0) {}
	Matrice(int);
	Matrice(const Matrice<T>&);
	~Matrice();

	friend istream& operator >>(istream& in, Matrice<T>& M)
	{
		for (int i = 0; i < M.ordin; ++i)
			for (int j = 0; j < M.ordin; ++j)
				in >> M[i][j];
		return in;
	}
	friend ostream& operator <<(ostream& out, Matrice<T> M)
	{
		for (int i = 0; i < M.ordin; ++i, out << endl)
			for (int j = 0; j < M.ordin; ++j)
				out << setw(10) << M[i][j] << ' ';
		return out;
	}

	vector<T>& operator[](int indexLinie)
	{
		return elemente[indexLinie];
	}

	Matrice<T>& operator =(const Matrice<T>&);

	bool operator ==(const Matrice<T>& A)
	{
		if (ordin != A.ordin)
			return false;
		for (int i = 0; i < ordin; ++i)
			for (int j = 0; j < ordin; ++j)
				if (elemente[i][j] != A.elemente[i][j])
					return false;
		return true;
	}
	bool operator !=(const Matrice<T>& A)
	{
		if (ordin != A.ordin)
			return true;
		int nrElemEgale = 0;
		for (int i = 0; i < ordin; ++i)
			for (int j = 0; j < ordin; ++j)
				if (elemente[i][j] == A.elemente[i][j])
					nrElemEgale++;
		if (nrElemEgale == ordin * ordin)
			return false;

		return true;
	}

	Matrice<T> operator +(Matrice<T>);
	Matrice<T> operator -(Matrice<T>);
	Matrice<T> operator *(Matrice<T>);
	Matrice<T> operator *(T);

	Matrice<T>& operator += (const Matrice<T>&);
	Matrice<T>& operator -= (const Matrice<T>&);
	Matrice<T>& operator *= (const Matrice<T>&);

	Matrice<T> operator^(int);

	Matrice<T> getTranspusa();

	T operator !();
	Matrice<T> getInversa();
};

template<class T>
Matrice<T>::Matrice(int ordin) : ordin(ordin)
{
	elemente.resize(ordin);
	for (int i = 0; i < ordin; ++i)
		for (int j = 0; j < ordin; ++j)
			elemente[i].push_back(0);
}

template<class T>
Matrice<T>::Matrice(const Matrice<T>& M)
{
	ordin = M.ordin;
	elemente.resize(ordin);

	for (int i = 0; i < ordin; ++i)
		for (int j = 0; j < ordin; ++j)
			elemente[i].push_back(M.elemente[i][j]);
}

template<class T>
Matrice<T>::~Matrice()
{
	elemente.clear();
	elemente.shrink_to_fit();
}

template<class T>
Matrice<T>& Matrice<T>::operator=(const Matrice<T>& M)
{
	if (this != &M)
	{
		this->~Matrice();

		ordin = M.ordin;
		elemente.resize(ordin);

		for (int i = 0; i < M.ordin; ++i)
			for (int j = 0; j < M.ordin; ++j)
				elemente[i].push_back(M.elemente[i][j]);
	}

	return *this;
}

template<class T>
Matrice<T> Matrice<T>::operator+(Matrice<T> M)
{
	Matrice<T> Suma(M.ordin);

	for (int i = 0; i < M.ordin; ++i)
		for (int j = 0; j < M.ordin; ++j)
			Suma[i][j] = elemente[i][j] + M[i][j];

	return Suma;
}

template<class T>
Matrice<T> Matrice<T>::operator-(Matrice<T> M)
{
	Matrice Diferenta(M.ordin);

	for (int i = 0; i < M.ordin; ++i)
		for (int j = 0; j < M.ordin; ++j)
			Diferenta[i][j] = elemente[i][j] - M[i][j];

	return Diferenta;
}

template<class T>
Matrice<T> Matrice<T>::operator*(Matrice<T> M)
{
	Matrice Produs(ordin);

	for (int i = 0; i < ordin; ++i)
		for (int j = 0; j < M.ordin; ++j)
			for (int k = 0; k < ordin; ++k)
				Produs[i][j] += elemente[i][k] * M[k][j];

	return Produs;
}

template<class T>
Matrice<T> Matrice<T>::operator*(T scalar)
{
	Matrice Produs(ordin);
	for (int i = 0; i < ordin; ++i)
		for (int j = 0; j < ordin; ++j)
			Produs[i][j] = elemente[i][j] * scalar;

	return Produs;
}

template<class T>
Matrice<T>& Matrice<T>::operator+=(const Matrice<T>& M)
{
	*this = *this + M;
	return *this;
}

template<class T>
Matrice<T>& Matrice<T>::operator-=(const Matrice<T>& M)
{
	*this = *this - M;
	return *this;
}

template<class T>
Matrice<T>& Matrice<T>::operator*=(const Matrice<T>& M)
{
	*this = *this * M;
	return *this;
}

template<class T>
Matrice<T> Matrice<T>::operator^(int putere)
{
	Matrice<T> matr(*this);
	if (putere == 0)
		return matr * matr.getInversa();

	for (int i = 1; i < putere; ++i)
		matr *= *this;

	return matr;
}

template<class T>
Matrice<T> Matrice<T>::getTranspusa()
{
	Matrice<T> transp(*this);

	for (int i = 0; i < ordin; ++i)
		for (int j = 0; j < ordin; ++j)
			transp.elemente[i][j] = elemente[j][i];

	return transp;
}

template<class T>
Matrice<T> Matrice<T>::getMinor(int linie, int coloana)
{
	Matrice<T> minor(ordin - 1);

	int nr = 0;
	for (int i = 0; i < ordin; ++i)
		for (int j = 0; j < ordin; ++j)
			if (i != linie && j != coloana)
			{
				minor[nr / (ordin - 1)][nr % (ordin - 1)] = elemente[i][j];
				nr++;
			}
	return minor;
}

template<class T>
T Matrice<T>::operator!()
{
	if (ordin == 1)
		return elemente[0][0];

	if (ordin == 2)
		return elemente[0][0] * elemente[1][1] - elemente[0][1] * elemente[1][0];

	float suma = 0;
	int semn = 1;
	for (int j = 0; j < ordin; ++j)
	{
		Matrice<T> minor = this->getMinor(0, j);
		suma += elemente[0][j] * (!minor) * semn;
		semn = -semn;
	}

	return suma;
}

template<class T>
void Matrice<T>::multiplicareLinie(int nrC, int linie, T a) // M_i(a) sau L_i -> a*L_i
{
	for (int i = 0; i < nrC; ++i)
		elemente[linie][i] = elemente[linie][i] * a;
}

template<class T>
void Matrice<T>::interschimbareLinii(int nrC, int linie1, int linie2) //S_ij sau L_i <-> L_j
{
	for (int i = 0; i < nrC; ++i)
		std::swap(elemente[linie1][i], elemente[linie2][i]);
}

template<class T>
void Matrice<T>::operatieLinii(int nrC, int linie1, int linie2, T a) //T_i(a*L_j) sau L_i -> a*L_j + L_i
{
	for (int j = 0; j < nrC; ++j)
		if ((T)(elemente[linie1][j] - a * elemente[linie2][j]) < (T)0.000001 &&
			(T)(elemente[linie1][j] - a * elemente[linie2][j]) > (T)-0.000001)
			elemente[linie1][j] = 0;
		else
			elemente[linie1][j] = elemente[linie1][j] - a * elemente[linie2][j];
}

template<class T>
int Matrice<T>::liberColoana(int nrL, int linie, int coloana)
{
	for (int i = linie; i < nrL; ++i)
		if (elemente[i][coloana] != 0)
			return i;
	return -1;
}

template<class T>
Matrice<T> Matrice<T>::getInversa() // Metoda Jordan-Gauss
{
	Matrice<T> mat(*this);
	for (int i = 0; i < ordin; i++)
		for (int j = 0; j < ordin; j++)
			if (i == j)
				mat[i].push_back(1);
			else
				mat[i].push_back(0);

	int contorLinie = 0, contorColoana = 0;
	while (contorLinie < ordin && contorColoana < ordin)
	{
		if (mat[contorLinie][contorColoana] != 1 && mat[contorLinie][contorColoana] != 0)
			mat.multiplicareLinie(2 * ordin, contorLinie, 1 / mat[contorLinie][contorColoana]);
		else
			if (mat[contorLinie][contorColoana] == 0)
			{
				int indiceLinie = mat.liberColoana(ordin, contorLinie, contorColoana);
				if (indiceLinie != -1)
					mat.interschimbareLinii(2 * ordin, contorLinie, indiceLinie);
				else
					contorColoana++;
			}
			else
				if (mat[contorLinie][contorColoana] == 1)
				{
					for (int i = 0; i < ordin; i++)
						if (i != contorLinie && mat[i][contorColoana] != 0)
							mat.operatieLinii(2 * ordin, i, contorLinie, mat[i][contorColoana]);

					contorLinie++;
					contorColoana++;
				}
	}

	Matrice<T> inversa(ordin);

	for (int i = 0; i < ordin; ++i)
		for (int j = ordin; j < 2 * ordin; ++j)
		{
			if (mat[i][j] == -0) mat[i][j] = 0;
			inversa[i][j - ordin] = mat[i][j];
		}

	return inversa;
}

// Matrice.cpp
#include "Matrice.h"

template<class T>
Matrice<T> suma(Matrice<T> A, int ordin, int putere)
{
	Matrice<T> tA = A.getTranspusa(), Suma(ordin);
	for (int i = 1; i <= putere; ++i)
	{
		Suma += (tA ^ i);

		cout << "(tA)^" << i;
		if (i < putere)
			cout << " + ";
		else
			cout << " = " << endl;
	}

	return Suma;
}

int main()
{
	int putere, ordin;
	cout << "Introduceti ordinul matricei : ";
	cin >> ordin;

	Matrice<float> A(ordin);
	cout << "Introduceti matricea : \n";
	cin >> A;
	cout << "Introduceti puterea : ";
	cin >> putere;

	cout << suma(A, ordin, putere) << endl;
	cout << "Determinantul matricei A : " << (!A) << endl;
	if (!A != 0)
		cout << "Inversa matricei A : \n" << A.getInversa();
	else
		cout << "Matricea nu este inversabila. \n";

	return 0;
}