• 0

Matrix inverse


Question

Hello,

I've written the following program to invert a NxN matrix. It's based on this : http://math.uww.edu/~mcfarlat/inverse.htm

#include <stdio.h>
#include <stdlib.h>
#include <assert.h>

/*Function Prototypes*/
int** createMatrix(int size);
void readMatrix(int** array,int size,FILE *fin);
double** inverseMatrix(int** array,int size);
void deallocate(void** array,int size);

int main()
{
	int** inMatrix;
	double** outMatrix;
	FILE *fin;
	int size,i,j;
	fin = fopen("in.dat","r");
	assert(fin); /*Check if call to fopen failed*/
	fscanf(fin,"%d %d",&size,&size); /*O pinakas einai tetragonikos ara oi grammes kai oi stiles einai ises*/
	inMatrix = createMatrix(size);
	readMatrix(inMatrix,size,fin);
	outMatrix = inverseMatrix(inMatrix,size);
	for (i=0; i < size; i++)
	{
		for (j=0; j < size; j++)
		{
			printf("%8.2f",outMatrix[i][j]);
		}
		printf("\n");
	}
	deallocate(outMatrix,size);
	deallocate(inMatrix,size);
	printf("Press any key to continue...");
	getchar();
	return 0;
}

int** createMatrix(int size)
{
	int** matrix;
	int i;
	matrix = malloc(sizeof(int*)*size);
	assert(matrix); /*Malloc success check*/
	for (i=0; i < size; i++)
	{
		matrix[i] = malloc(sizeof(int)*size);
		assert(matrix[i]);
	}
	return matrix;
}

void readMatrix(int** array,int size,FILE *fin)
{
	int data;
	int i,j;
	for (i=0; i < size; i++)
	{
		for (j=0; j < size; j++)
		{
			fscanf(fin,"%d",&data);
			array[i][j] = data;
		}
	}
}

double** inverseMatrix(int** array,int size)
{
	double data;
	double** augMatrix;
	double** tempMatrix;
	int i,j,k;
	double x;
	tempMatrix = malloc(sizeof(double*)*size);
	assert(tempMatrix);
	augMatrix = malloc(sizeof(double*)*size);
	assert(augMatrix);
	for (i=0; i < size; i++)
	{
		tempMatrix[i] = malloc(sizeof(double)*size);
		assert(tempMatrix[i]);
		augMatrix[i] = malloc(sizeof(double)*size);
		assert(augMatrix[i]);
	}
	/*Dimiourgia tou epauksimenou pinaka*/
	for (i=0; i < size; i++)
	{
		for (j = 0; j< size; j++)
		{
			if (i==j)
			{
				augMatrix[i][j] = 1;
			}
			else
			{
				augMatrix[i][j] = 0;
			}
			 tempMatrix[i][j] = array[i][j];
		}

	}
	for (i = 0; i < size; i++)
	{
		data = tempMatrix[i][i];
		for (j=0; j < size; j++)
		{
			tempMatrix[i][j] =  tempMatrix[i][j]/data;
			augMatrix[i][j] = augMatrix[i][j]/data;
		}

		for (j = 0; j < size; j++)
		{
			if (j != i)
			{

			   x = tempMatrix[j][i]/tempMatrix[i][i];
			   for (k = 0; k < size; k++)
			   {
				   tempMatrix[j][k] = tempMatrix[j][k] - x*tempMatrix[i][k];
				   augMatrix[j][k] = augMatrix[j][k] - x*augMatrix[i][k];
			   }

			}
		}
	}
	deallocate(tempMatrix,size);
	return augMatrix;
}

void deallocate(void** array,int size)
{
	int i;
	for (i = 0; i < size; i++)
	{
		free(array[i]);
	}
	free(array);
}

Normally I would debug it on my own but I am a bit tight on time right now.

Can someone tell me if there is something wrong with my algorithm or have I messed something up with the language?

EDIT : nvm I found some logic problems, I am working on it.

EDIT 2 : fixed it.

Edited by gianpan
Link to comment
Share on other sites

4 answers to this question

Recommended Posts

  • 0

OK, I fixed this! However I still can't get it to work correctly (the results are all wrong:P)

Thanks for the help on this one it was a very silly mistake..

Link to comment
Share on other sites

This topic is now closed to further replies.
  • Recently Browsing   0 members

    • No registered users viewing this page.