matrix problem

[eduard77]

How can I define this algorithm to calculate the matrix that has uneven nr of rows and columns?

//==============================================================================
//
// Linear System Solution by Gauss method
//
// 
//
//==============================================================================
#include <iostream>
#include <stdio.h>
#include <windows.h>
#include <cmath>

//==============================================================================
void VectorPrint(int nDim, double* pfVect)
{
  int i;

  printf("------------------------------------------------------------------\n");
  for(i=0; i<nDim; i++)
  {
    printf("%9.2lf ", pfVect[i]);
  }
  printf("\n");
}
//==============================================================================
void MatrixPrint(int nDim, double* pfMatr)
{
  int i,j;

  printf("------------------------------------------------------------------\n");
  for(i=0; i<nDim; i++)
  {
    for(j=0; j<nDim; j++)
    {
      printf("%9.2lf ", pfMatr[nDim*i + j]);
    }
    printf("\n");
  }
}

//==============================================================================
// return 1 if system not solving
// nDim - system dimension
// pfMatr - matrix with coefficients
// pfVect - vector with free members
// pfSolution - vector with system solution
// pfMatr becames trianglular after function call
// pfVect changes after function call
//
//
//
//==============================================================================
int LinearEquationsSolving(int nDim, double* pfMatr, double* pfVect, double* pfSolution)
{
  double fMaxElem;
  double fAcc;

  int i , j, k, m;




  for(k=0; k<(nDim-1); k++) // base row of matrix
  {
    // search of line with max element
 fMaxElem= fabs(pfMatr[k*nDim + k]);
    m = k;
    for(i=k+1; i<nDim; i++)
    {
      if(fMaxElem < fabs(pfMatr[i*nDim + k]) )
      {
        fMaxElem = pfMatr[i*nDim + k];
        m = i;
      }
    }

    // permutation of base line (index k) and max element line(index m)
    if(m != k)
    {
      for(i=k; i<nDim; i++)
      {
        fAcc               = pfMatr[k*nDim + i];
        pfMatr[k*nDim + i] = pfMatr[m*nDim + i];
        pfMatr[m*nDim + i] = fAcc;
      }
      fAcc = pfVect[k];
      pfVect[k] = pfVect[m];
      pfVect[m] = fAcc;
    }

    if( pfMatr[k*nDim + k] == 0.) return 1; // needs improvement !!!

    // triangulation of matrix with coefficients
    for(j=(k+1); j<nDim; j++) // current row of matrix
    {
      fAcc = - pfMatr[j*nDim + k] / pfMatr[k*nDim + k];
      for(i=k; i<nDim; i++)
      {
        pfMatr[j*nDim + i] = pfMatr[j*nDim + i] + fAcc*pfMatr[k*nDim + i];
      }
      pfVect[j] = pfVect[j] + fAcc*pfVect[k]; // free member recalculation
    }
  }

  for(k=(nDim-1); k>=0; k--)
  {
    pfSolution[k] = pfVect[k];
    for(i=(k+1); i<nDim; i++)
    {
      pfSolution[k] -= (pfMatr[k*nDim + i]*pfSolution[i]);
    }
    pfSolution[k] = pfSolution[k] / pfMatr[k*nDim + k];
  }

  return 0;
}
//==============================================================================
// testing of function
//==============================================================================

#define MATRIX_DIMENSION 11

int main(int nArgs, char** pArgs)
{
  int nDim = MATRIX_DIMENSION;
  double fMatr[MATRIX_DIMENSION*MATRIX_DIMENSION] =
  {
   9.00,  12.00, 11.00, 15.00, 27.00, 46.00, 34.00, 50.00, 64.00,  2.71, 35.00,
   8.00,  10.45,  9.11,  12.1, 21.76, 39.42, 30.19, 39.00, 36.83,  0.00,  0.00,
   4.00,   2.00,  2.00,  4.00,  8.00,  1.50,  1.50, 10.00, 39.00,  0.00,  0.00,  
   1.50,   1.50,  2.50,  5.50,  4.50,  6.00,  0.00, 30.00, 15.50,  0.00,  0.00,
  33.80,  30.80, 27.10, 14.50, 23.48, 22.30, 13.00, 23.50, 29.40, 88.00, 18.14,
   2.00,   2.00,  5.00, 10.50,  8.50, 10.50, 23.00,  0.00,  0.00,  0.00,  4.00,
   0.26,   0.34,  0.41,  0.63,  0.60,  2.89,  1.21,  2.42,  5.05,  0.00,  0.95,
   0.20,   0.28,  0.33,  0.46,  0.36,  2.56,  1.05,  1.90,  4.32,  0.00,  0.00,

  };
  double fVec[MATRIX_DIMENSION] = {};

  double fSolution[MATRIX_DIMENSION];
  int res;
  int i;

  res = LinearEquationsSolving(nDim, fMatr, fVec, fSolution); // !!!

  if(res)
  {
    printf("No solution!\n");
    return 1;
  }
  else
  {
    printf("Solution:\n");
    VectorPrint(nDim, fSolution);

  }


  return 0;
}