#include <math.h>
#include <iomanip>
#include "Polinomio.h"
#include "Matriz.h"
using namespace std;
#ifndef _MATRIZ_H_
#define _MATRIZ_H_
class Matriz
{
double **a; //array de punteros a array double, contiene los terminos de la matriz
protected:
int length; //orden de la matriz (n filas x ncolumnas)
double det;
Polinomio ct; //Polinomio caracteristico de la matriz
public:
Matriz (int n); // reserva memoria, lenght=n: MatrizTridiag A(n);
~Matriz (); // destructor
bool redim(int n);
void operator=(const Matriz& B); // A=B;
double operator () (int i,int j) const ; // a=m(i,j);
double& operator () (int i,int j); // function: m(i,j)=a;
Matriz operator + (const Matriz& B); // C=A+B;
void read(); // A.read();
void fill(); // A.fill();
void print () const; // A.print();
double menores(int i);
void factdet();
void ccaracteristico();
void mcaracteristico();
void cambiarfilas(int i,int j);
double determinante() {return det;}
int size () const {return length;} // n = A.size();
};
#endif
// Matriz A(n);
Matriz::Matriz(int n)
{
int i,j;
length = n;
det=1;
ct.q= new double [length];
for(i=0;i<n;i++) ct.q[i]=0;
// reserva memoria
a = new double *[n]; // reserva memoria para el vector de punteros
for (i=0;i<n;i++) a[i] = new double[n]; // reserva memoria para cada fila
for(i=0;i<n;i++) for(j=0;j<n;j++) a[i][j]=0; // inicializa a cero
}
// destructor
Matriz::~Matriz()
{
if( a != NULL)
{ // libera memoria
for(int i=0;i<length;i++) delete [] a[i];
delete [] a;
}
if( ct.q != NULL)
{
delete [] ct.q;
}
}
// A = B
void Matriz::operator = (const Matriz& B)
{
redim(B.size());
for(int i=0;i<length;i++) for(int j=0;j<length;j++) a[i][j]=B.a[i][j];
}
// function: a=m[i,j];
inline double Matriz::operator() (int i,int j) const
{ return a[i][j];}
// function: m[i,j]=a;
double& Matriz::operator() (int i,int j)
{ return a[i][j];}
// operador suma: C=A+B
Matriz Matriz::operator+(const Matriz& B)
{
if(length != B.length) cout<< " error: matrices de distinta dimension "<<endl;
Matriz temp(length);
for(int i=0;i<length;i++) for(int j=0;j<length;i++) temp.a[i][j]=a[i][j]+B.a[i][j];
return temp;
}
// producto matriz por vector: A.prod(v,w);
// redimensionamiento de una matriz
bool Matriz::redim(int n)
{
if (length == n)
return false; // no es necesario redimensionar
else
{
if( a!= NULL)
{
// primero libera memoria
for(int i=0;i<length;i++) delete [] a[i];
delete [] a;
}
length = n;
// reserva memoria
a = new double *[n]; // reserva memoria para el vector de punteros
for (int i=0;i<n;i++) a[i] = new double [n]; // reserva memoria para cada fila
return true;
}
}
// A.read()
void Matriz::read ()
{
for(int i=0;i<length;i++)
{
cout << "fila "<<i<<"-esima" << "n";
for(int j=0;j<length;j++)
{
cout << "termino a["<<i<<"]["<<j<<"] = ?" << "n";
cin >> a[i][j];
}
}
}
// A.fill()
void Matriz::fill()
{
for(int i=0;i<length;i++)
{
a[i][i]=4;
if(i<length-1)
{
a[i+1][i]=-1;
a[i][i+1]=-1;
}
if(i<length-3)
{
a[i+3][i]=-1;
a[i][i+3]=-1;
}
}
}
// A.print()
void Matriz::print () const
{
for(int i=0;i<length;i++)
{
for(int j=0;j<length;j++)
cout << a[i][j]<<" ";
cout<<endl;
}
}
// A.factLU()
void Matriz::factdet()
{
int permutaciones=0;
int i,j;
for(int cp=0;cp<length-1;cp++)
{
for(i=1;i<length-cp;i++){
if(a[cp][cp]==0) //si el pivote es cero cambiar filas
{
cambiarfilas(cp,cp+i);
permutaciones++;
}
}
if(a[cp][cp]==0)
det= 0;//si el pivote sigue siendo cero, el det es cero.
//ahora ya aplicar formula de ferragut
for(i=cp+1;i<length;i++)
{
double pivote=a[i][cp]/a[cp][cp];
for(j=0;j<length;j++)
{
a[i][j]=a[i][j]-pivote*a[cp][j];
}
}
}
//se calcula el determinante
for(i=0;i<length;i++)
{det=det*a[i][i];}
det=det*pow(-1,permutaciones);
}
void Matriz::cambiarfilas(int i, int j){
double tem;
for(int l=0;l<length;l++){
tem=a[i][l];
a[i][l]=a[j][l];
a[j][l]=tem;
}
}
double Matriz::menores(int i){
double coef=0;
int j,k,l;
for (j=0;j<(length-i+1);j++);{ //recorre la diagonal en busca de menores
Matriz menor(i); //crea menores de orden deseado
for (k=0;k<i;k++);
{
for (l=0;l<i;l++); //los iguala a la matriz menor para no perturbar la principal
menor.a[k+j][l+j]=a[k+j][l+j];
}
menor.factdet(); //calcula el determinante
coef+=menor.det; //los va sumando
menor.~Matriz();
}
return coef;
}//en caso de pasar i=length devuelve el determinante de la matriz
void Matriz::ccaracteristico(){
ct.q[length-1]=pow(-1,length);
for(int i=0;i<length-1;i++){
ct.q[i]=menores(length-i)*pow(-1,i);
}
}
void Matriz::mcaracteristico(){
ct.print1();
}
int main() {
Matriz prueba(2);
prueba.fill();
Matriz p(2);
p.fill();
Matriz c(2);
c=prueba+p; //al llegar aqui
c.print();
p.ccaracteristico(); //o aqui, sale del programa sin terminar de ejecutarlo
p.mcaracteristico();
return 0;
}
