Les paso unas rutinas que creé para el cálculo matricial en C++.

Se las puede utilizar como vector o como estructura (está sobrecargado)

CABECERA

Código:

struct stMatrices
{
 int rx;
 int ry;
 vector<float> Mat;

 public:
 //FUNCIONES
 stMatrices(void); //Constructor sobrecargado
 stMatrices(int x, int y); //Constructor
 ~stMatrices(void); //Destructor
 float get(int x,int y); //Retorna el valor de la posición (x,y)
 void set(int x, int y, float val); //Asigna el valor val a la posición (x,y)
 void redim (int rx, int ry); //Redimensiona la matriz
 int ubound(int col); //Retorna el máximo valor de x o y
};


//FUNCIONES CON VECTORES
void MMuestre(vector<float> Mat,int R1,int R2);
vector<float> MInversa(vector<float> Mat, int R1);
vector<float> MMult(vector<float> MatA, vector<float> MatB, int RAF, int RAC, int RBC);
vector<float> MTransp(vector<float> Mat, int RF, int RC);
vector<float> MkMult(vector<float> Mat, int RF, int RC,float factK);
vector<float> MSuma(vector<float> MatA, vector<float> MatB,int RF, int RC,float factK);

//FUNCIONES CON ESTRUCTURAS
void MMuestre(stMatrices Mat);
stMatrices MInversa(stMatrices Mat);
stMatrices MMult(stMatrices MatA, stMatrices MatB);
stMatrices MTransp(stMatrices Mat);
stMatrices MkMult(stMatrices Mat,float factK);
stMatrices MSuma(stMatrices MatA, stMatrices MatB,float factK);

MÓDULO

Código:

//[O]
//Herramientas matriciales

#include <iostream>
#include <string>
#include <cstdio>
#include <vector>
#include "mat.h"

using namespace std;

//*********************************
// FUNCIONES DEL OBJETO
//*********************************

stMatrices::stMatrices(void) //Constructor sobrecargado
{
 redim(0,0);
}

stMatrices::stMatrices(int x, int y) //Constructor
{
 redim(x,y);
}

stMatrices::~stMatrices(void) //Destructor
{
}

float stMatrices::get(int x,int y)
{
 if (x <= rx)
  if (y <= ry)
   return Mat[x*ry+y];
}

void stMatrices::set(int x, int y, float val)
{
 if (x <= rx)
  if (y <= ry)
   Mat[x*ry+y] = val;
}

void stMatrices::redim (int x, int y)
{
 rx = x;
 ry = y;
 Mat.resize(rx*ry);
}

int stMatrices::ubound(int col)
{
 if (col == 1)
  return rx;
 else if (col == 2)
  return ry;
 else
  return 0;
}


//************************************************
// FUNCIONES MATRICIALES CON VECTORES
//************************************************

void MMuestre (vector<float> Mat,int R1,int R2)
{
 for(int i = 0;i < R1; i++)
 {
  for (int j = 0; j < R2; j++)
   printf("%f ",Mat[i*R2+j]);
  printf("\n");
 }
}


vector<float> MInversa(vector<float> Mat, int R1)
{
//Inversión de la Matriz 'Mat', de rango R1, por el método de Gauss-Jordan
//int R1 = El rango de la Matriz

 vector<float> MatM;
 int R12 = R1*2;
 MatM.resize(R1*R12);

//Creo la Matriz complementaria de la forma [M|I], donde M = Mat e I = Identidad
//1º la parte de MatM = [Mat|...]

 for (int i = 0; i < R1; i++)
  for (int j = 0; j < R1; j++)
   MatM[i * R12 + j] = Mat[i * R1 + j];


//2º La parte de MatM = [...|I]

 for (int i = 0; i < R1; i++)
  for (int j = 0; j < R1; j++)
   if (i == j)
    MatM[i * R12 + R1 + j] = 1;
   else
    MatM[i * R12 + R1 + j] = 0;



//INVERSIÓN PROPIAMENTE DICHA
 float Pibot=0; //Pibot
 int iP = 0 ; //Índice del Pibot
 float factK = 0; //Factor K
 for (int i = 0; i < R1; i++)
 {
  iP = i; //Índice del Pibot
  Pibot = MatM[iP * R12 + iP]; //Obtiene el Pibot

  for (int k = iP; k < (R1 * 2 + 1); k++) //Reduce la FILA del pibot
   MatM[i * R12 + k] = MatM[i * R12 + k] / Pibot;

 //Fila n = Fila n -k.Fila n
 //Obtiene el 0 delante del pibot
  for (int h = 0; h < R1; h++)
   if (h != iP) //La fila del pibot no se toca
   {
    factK = MatM[h * R12 + iP]; //Obtiene el factor k
    for (int k = 0; k < R1 * 2 + 1; k++)
     MatM[h*R12+k] = MatM[h*R12+k] - factK * MatM[iP*R12+k];
   }

 }

//FINAL PRESENTACIÓN DE LOS DATOS
 vector<float> MatT;
 MatT.resize(R1*R1);
 for (int i = 0; i < R1; i++)
  for (int j = 0; j < R1; j++)
   MatT[i*R1+j] = MatM[i * R12 + j + R1];

 return MatT;
}


vector<float> MMult(vector<float> MatA, vector<float> MatB, int RAF, int RAC, 
int RBC){
//Multiplica MatA x MatB, siendo
//RAF = Cantidad de filas de MatA
//RAC = Cantidad de columnas de MatA
//RBC = Cantidad de columnas de MatB
//La cantidad de filas de MatB es la cantidad de columnas de MatA
//La respuesta tiene la cantidad de filas de MatA y de columnas de MatB
 int RBF = RAC;
 vector<float> MatR;
 float val=0;
 MatR.resize(RAF*RBC);

 for (int i = 0; i < RAF; i++)
  for (int j = 0; j < RBC; j++)
  {
   for (int k = 0; k < RAC; k++)
    val += MatA[i * RAC + k] * MatB[k * RBC + j];
   MatR[i * RBC + j] = val;
   val = 0;
  }
 return MatR;
}


vector<float> MTransp(vector<float> Mat, int RF, int RC)
{
//Transpone la matriz MatA
 vector<float> MatR;
 MatR.resize(RF*RC);

 for (int i = 0; i < RF; i++)
  for (int j = 0; j < RC; j++)
   MatR[j * RF + i] = Mat[i * RC + j];
 return MatR;
}

vector<float> MkMult(vector<float> Mat, int RF, int RC,float factK)
{
//Multiplica una matriz por una constante
 vector<float> MatR;
 MatR.resize(RF*RC);
 for (int i = 0; i < RF; i++)
  for (int j = 0; j < RC; j++)
   MatR[i*RC+j]= Mat[i*RC+j] * factK;
 return MatR;
}

vector<float> MSuma(vector<float> MatA, vector<float> MatB,int RF, int RC,float 
factK)
{
//Multiplica suma 2 matrices multiplicando a la segunda por un factor
//MatR = MatA + k.MatB
//Ej. para una resta de matrices, factK = -1
 vector<float> MatR;
 MatR.resize(RF*RC);
 for (int i = 0; i < RF; i++)
  for (int j = 0; j < RC; j++)
   MatR[i*RC+j]= MatA[i*RC+j] + factK * MatB[i*RC+j];
 return MatR;
}




//****************************************************
// FUNCIONES MATRICIALES CON ESTRUCTURAS
//****************************************************

void MMuestre (stMatrices Mat)
{
 for(int i = 0;i < Mat.ubound(1) ; i++)
 {
  for (int j = 0; j < Mat.ubound(2); j++)
   printf("%f ",Mat.get(i,j));
  printf("\n");
 }
}

stMatrices MInversa(stMatrices Mat)
{
//Inversión de la Matriz 'Mat', de rango R1, por el método de Gauss-Jordan


 int R1 = Mat.ubound(1);
 int R12 = R1*2;
 stMatrices MatM(R1,R12);
 //MatM.resize(R1,R12);

//Creo la Matriz complementaria de la forma [M|I], donde M = Mat e I = Identidad
//1º la parte de MatM = [Mat|...]

 for (int i = 0; i < R1; i++)
  for (int j = 0; j < R1; j++)
   MatM.set(i,j, Mat.get(i,j));

//2º La parte de MatM = [...|I]

 for (int i = 0; i < R1; i++)
  for (int j = 0; j < R1; j++)
   if (i == j)
    MatM.set(i,j+R1,1);
   else
    MatM.set(i,j+R1,0);


//INVERSIÓN PROPIAMENTE DICHA
 float Pibot=0; //Pibot
 int iP = 0 ; //Índice del Pibot
 float factK = 0; //Factor K
 for (int i = 0; i < R1; i++)
 {
  iP = i; //Índice del Pibot
  Pibot = MatM.get(iP,iP); //Obtiene el Pibot

  for (int k = iP; k < (R1 * 2 + 1); k++) //Reduce la FILA del pibot
   MatM.set(i,k, MatM.get(i,k) / Pibot);

 //Fila n = Fila n -k.Fila n
 //Obtiene el 0 delante del pibot
  for (int h = 0; h < R1; h++)
   if (h != iP) //La fila del pibot no se toca
   {
    factK = MatM.get(h,iP); //Obtiene el factor k
    for (int k = 0; k < R1 * 2 + 1; k++)
     MatM.set(h,k, MatM.get(h,k) - factK * MatM.get(iP,k));
   }

 }

//FINAL PRESENTACIÓN DE LOS DATOS
 stMatrices MatT(R1,R1);
 for (int i = 0; i < R1; i++)
  for (int j = 0; j < R1; j++)
   MatT.set(i,j, MatM.get(i,j+R1));

 return MatT;
}


stMatrices MMult(stMatrices MatA, stMatrices MatB)
{
//Multiplica MatA x MatB, siendo
//RAF = Cantidad de filas de MatA
//RAC = Cantidad de columnas de MatA
//RBC = Cantidad de columnas de MatB
//La cantidad de filas de MatB es la cantidad de columnas de MatA
//La respuesta tiene la cantidad de filas de MatA y de columnas de MatB
 int RAF = MatA.ubound(1);
 int RAC = MatA.ubound(2);
 int RBF = MatB.ubound(1);
 int RBC = MatB.ubound(2);
//  if (RBF == RAC)
//   return;

 stMatrices MatR(RAF,RBC);
 float val = 0;

 for (int i = 0; i < RAF; i++)
  for (int j = 0; j < RBC; j++)
  {
   for (int k = 0; k < RAC; k++)
    val += MatA.get(i,k) * MatB.get(k,j);
   MatR.set(i,j,val);
   val = 0;
  }
 return MatR;
}


stMatrices MTransp(stMatrices Mat)
{
//Transpone la matriz Mat
 int RF = Mat.ubound(1);
 int RC = Mat.ubound(2);
 stMatrices MatR(RC,RF);

 for (int i = 0; i < RF; i++)
  for (int j = 0; j < RC; j++)
   MatR.set(j,i, Mat.get(i,j));
 return MatR;
}

stMatrices MkMult(stMatrices Mat,float factK)
{
//Multiplica una matriz por una constante
 int RF = Mat.ubound(1);
 int RC = Mat.ubound(2);
 stMatrices MatR(RF,RC);

 for (int i = 0; i < RF; i++)
  for (int j = 0; j < RC; j++)
   MatR.set(i,j, Mat.get(i,j)*factK);
 return MatR;
}

stMatrices MSuma(stMatrices MatA, stMatrices MatB, float factK)
{
//Multiplica suma 2 matrices multiplicando a la segunda por un factor
//MatR = MatA + k.MatB
//Ej. para una resta de matrices, factK = -1
 int RF = MatA.ubound(1);
 int RC = MatA.ubound(2);
 stMatrices MatR(RF,RC);
 for (int i = 0; i < RF; i++)
  for (int j = 0; j < RC; j++)
   MatR.set(i,j, MatA.get(i,j) + factK * MatB.get(i,j));
 return MatR;
}

Saludos

ahhh esta bien loco tu codigo, pero la verdad sera de gran ayuda!!!

gracias Vanchi!!!

Att. player1