Matrix Manipulation Using C

            Matrix  is a great invention of Mathematics. It is very useful in solving system of equations through computer. Of course famous numerical computing platform, MATLAB is the well known application of matrix.
            In MATLAB everything is matrix. The images are also represented as a matrix in ease. It is easy to perform operations over matrices. This is why the matrix manipulation is an important.

Download the code here

            Here is my matrix library. It may be inefficient, since it is indented for educational purpose. I want to keep things simple so that the beginner also be able to understand. This article is also helpful for learning memory allocation in C, Pointer and its applications, recursion,  modularity.
                                                                
CodeListing

1. /*----------------------------------------*/
2. /* Program: matrix.h                      */
3. /* Author: Vanangamudi                    */
4. /* Date: 15/05/2011                       */
5. /*----------------------------------------*/
6. #include<stdio.h>
7. #ifndef _MATRIX_H_
8. #define _MATRIX_H_
9. 
   
10. typedef struct matrix {
11. float **mat;
12. int m,n;
13. }Matrix;
14. 
    
15. typedef struct matrixEquation {
16. Matrix *lhs;
17. Matrix *rhs;
18. Matrix *Soln;
19.  }MatrixEqn;
20. 
    
21. Matrix * matrix_new_Matrix(int m,int n);
22. MatrixEqn * matrix_new_MatrixEqn(Matrix* lhs,Matrix *rhs);
23. 
    
24. float matrix_Det(Matrix *matA);
25. Matrix* matrix_Cofactor(Matrix *matA);
26. Matrix* matrix_Minor(Matrix *matA, int m, int n);
27. Matrix* matrix_Adjoint(Matrix *matA);
28. Matrix* matrix_Inverse(Matrix *matA);
29. 
    
30. int matrix_Solve_Det(MatrixEqn *matEqn);
31. int matrix_Solve_Inv(MatrixEqn *matEqn);
32. 
    
33. Matrix* matrix_Mult(Matrix *matA , Matrix *matB);
34. Matrix* matrix_Add(Matrix *matA , Matrix *matB);
35. Matrix* matrix_Sub(Matrix *matA , Matrix *matB);
36. Matrix* matrix_tran(Matrix *matA);
37. Matrix* matrix_ScalarMultiply(float X , Matrix *matA);
38. Matrix* matrix_Duplicate(Matrix *matA);
39. Matrix* matrix_ColumnReplace(Matrix *matA, Matrix* matB ,int IndA,int IndB );
40. 
    
41. int matrix_Print(Matrix *matA);
42. int matrix_Read(Matrix *matA);
43. 
    
44. #endif       //_MATRIX_H_

1. /*-----------------------------------------------*/
2. /* Program: matrix.c                             */
3. /* Author: Vanangamudi                           */
4. /* Date: 15/05/2011                              */
5. /*-----------------------------------------------*/
6. 
   
7. #include "matrix.h"
8. #include<malloc.h>
9. #include<stdio.h>
10. 
    
11. //Constructors
12. Matrix* matrix_new_Matrix(int m, int n)
13. {
14.        Matrix *newMat = (Matrix*)malloc(sizeof(float));
15.        newMat->mat = (float*)malloc(sizeof(float)*m);
16.     int i = 0;
17.     for(i =0 ; i < m ; i++){
18.           newMat->mat[i] = (float*)malloc(sizeof(float)*n);
19.     }
20.       newMat->m = m;
21.       newMat->n = n;
22.     return newMat;
23. }
24. 
    
25. MatrixEqn * matrix_new_MatrixEqn(Matrix *lhs,Matrix *rhs)
26. {
27.       MatrixEqn* newEqn = (MatrixEqn*)malloc(sizeof(float));
28.       newEqn->lhs = matrix_Duplicate(lhs); 
29.       newEqn->rhs = matrix_Duplicate(rhs);
30.       newEqn->Soln = matrix_new_Matrix(rhs->m,1);
31.     return newEqn;
32. }
33. 
    
34. // Matrix Operations
35. Matrix* matrix_Minor(Matrix *matA, int m, int n)
36. {
37.     int i,j;
38.     int p,q;
39.     Matrix *minor = matrix_new_Matrix(matA->m-1,matA->n-1);
40.     for(i = 0,p = 0; i < matA->m , p < minor->m; i++){
41.       if( i != m-1){
42.          for(j=0, q=0; j <matA->n , q < minor->n; j++){
43.             if(j != n-1){
44.               minor->mat[p][q] = matA->mat[i][j];
45.               q++;
46.             }
47.          }
48.          p++;
49.       }
50.     }
51.     return minor;
52. }
53. 
    
54. float matrix_Det(Matrix *matA)
55. {
56.     float det=0;
57.     if((matA->m==1)&&(matA->n==1)){
58.           det = matA->mat[0][0];
59.     }else {
60.       int i;
61.       int sign=1;
62.       for(i = 0; i<matA->m ; i++){
63.          det += sign*matA->mat[0][i]*matrix_Det(matrix_Minor(matA,0,i));
64.           sign = -1*sign;
65.       }
66.     }
67.     return det;
68. }
69. 
    
70. Matrix* matrix_Cofactor(Matrix *matA)
71. {
72.     Matrix *Cofact = matrix_new_Matrix(matA->m,matA->n);
73.     int i,j;
74.     int sign =1;
75.     for(i =0 ; i<matA->m;i++){
76.        for(j=0 ; j<matA->n ; j++){
77.           Cofact->mat[i][j] = sign*matrix_Det(matrix_Minor(matA,i,j));
78.           sign = sign * (-1);
79.        }
80.        if(matA->m%2 == 0){
81.               sign = (-1)*sign;
82.        }
83.     }
84.     return Cofact;
85. }
86. 
    
87. Matrix* matrix_Adjoint(Matrix *matA)
88. {
89.      Matrix *Cofact = matrix_Cofactor(matA);
90.      Matrix *Adj = matrix_tran(Cofact);
91.      return Adj;
92. }
93. 
    
94. Matrix* matrix_Inverse(Matrix *matA)
95. {
96.     Matrix *Inv = matrix_Adjoint(matA);
97.     float Det = matrix_Det(matA);
98.     Inv = matrix_ScalarMultiply(1/Det, Inv);
99.     return Inv;
100. }
101. 
     
102. //Arithmetics
103. Matrix* matrix_Mult(Matrix *matA , Matrix *matB)
104. {
105.      if(matA->n != matB->m){
106.         printf("\nMatrices cannot be multiplied");
107.         return NULL;
108.      }else{
109.              Matrix *matC = matrix_new_Matrix(matA->m,matB->n);
110.         int i,j,k;
111.         for(i = 0 ;i< matB->m; i++){
112.            for(j=0; j < matB->n ; j++){
113.               matC->mat[i][j] = 0;
114.               for( k =0;k < matA->n;k++){
115.                  matC->mat[i][j]+=matA->mat[i][k]*matB->mat[k][j];
116.               }
117.             }
118.          }
119.         return matC;
120.      }
121. }
122. 
     
123. Matrix* matrix_ScalarMultiply(float X , Matrix *matA)
124. {
125.         Matrix *matB = matrix_new_Matrix(matA->m,matA->n);
126.      int i,j;
127.      for(i = 0 ;i< matA->m; i++){
128.         for(j=0; j < matA->n ; j++){
129.            matB->mat[i][j] = X * matA->mat[i][j];
130.         }
131.      }
132.      return matB;
133. }
134. 
     
135. Matrix* matrix_Add(Matrix *matA , Matrix *matB)
136. {
137.         Matrix *matC = matrix_new_Matrix(matA->m,matA->n);
138.      int i,j;
139.      for(i = 0 ;i< matA->m; i++){
140.         for(j=0; j < matA->n ; j++){
141.           matC->mat[i][j] = matA->mat[i][j]+matB->mat[i][j];
142.         }
143.      }
144.       return matC;
145. }
146. 
     
147. Matrix* matrix_Sub(Matrix *matA , Matrix *matB)
148. {
149.         Matrix *matC = matrix_new_Matrix(matA->m,matA->n);
150.      int i,j;
151.      for(i = 0 ;i< matA->m; i++){
152.         for(j=0; j < matA->n ; j++){
153.            matC->mat[i][j] = matA->mat[i][j]-matB->mat[i][j];
154.         }
155.      }
156.      return matC;
157. }
158. 
     
159. Matrix* matrix_tran(Matrix *matA)
160. {
161.     Matrix *matB = matrix_new_Matrix(matA->n,matA->m);
162.     int i,j;
163.     for(i = 0 ;i< matB->m; i++){
164.        for(j=0; j < matB->n ; j++){
165.           matB->mat[i][j] = matA->mat[j][i];
166.        }
167.     }
168.     return matB;
169. }
170. 
     
171. //Matrix Eqn Solvers
172. int matrix_Solve_Inv(MatrixEqn *matEqn)
173. {
174.     Matrix* Inv = matrix_Inverse(matEqn->lhs);
175.     matEqn->Soln = matrix_Mult(Inv,matEqn->rhs);
176. }
177. 
     
178. int matrix_Solve_Det(MatrixEqn *matEqn)
179. {
180.     float Det = matrix_Det(matEqn->lhs);
181.     Matrix* matA;
182.     float DetXYZ;
183.     int i;
184.     for(i = 0;i <matEqn->Soln->m;i++){
185.       matA = matrix_ColumnReplace(matEqn->lhs,matEqn->rhs,i+1,1);
186.       DetXYZ = matrix_Det(matA);
187.       matEqn->Soln->mat[i][0] = DetXYZ/Det;
188.     }
189. }
190. 
     
191. //Misc
192. Matrix* matrix_ColumnReplace(Matrix *matA, Matrix* matB ,int IndA,int IndB )
193. {
194.        Matrix* matC = matrix_Duplicate(matA);
195.     int i;
196.     for(i = 0;i<matA->n; i++){
197.             matC->mat[i][IndA-1] = matB->mat[i][IndB-1];
198.     }
199.     return matC;
200. }
201. 
     
202. Matrix* matrix_Duplicate(Matrix *matA)
203. {
204.        Matrix *matB = matrix_new_Matrix(matA->m,matA->n);
205.     int i,j;
206.     for(i = 0 ;i< matB->m; i++){
207.        for(j=0; j < matB->n ; j++){
208.                 matB->mat[i][j] = matA->mat[i][j];
209.        }
210.     }
211.     return matB;
212. }
213. 
     
214. //IO Operations
215. int matrix_Read(Matrix *matA)
216. {
217.      int i=0 , j=0;
218.         printf("\nEnter the Elements:\n");
219.      for(i = 0; i < matA->m; i++){
220.         for(j = 0; j < matA->n; j++ ){
221.                  scanf("%f",&matA->mat[i][j]);
222.         }
223.      }
224. }
225. 
     
226. int matrix_Print(Matrix *matA)
227. {
228.      int i=0 , j=0;
229.         printf("\n");
230.      for(i = 0; i < matA->m; i++){
231.         for(j = 0; j < matA->n; j++ ){
232.                  printf("%0.3f\t",matA->mat[i][j]);
233.         }
234.              printf("\n\n");
235.      }
236.      return 0;
237. }

Sample Program to Test

1. /*-----------------------------------------------------*/
2. /* Program: matrix_test.c                              */
3. /* Author  : Vanangamudi                               */
4. /* Date     :15/05/2011                                */
5. /*-----------------------------------------------------*/
6. 
   
7. #include "matrix.h"
8. #include<stdio.h>
9. 
   
10. int main()
11. {
12.     Matrix *matA  = matrix_new_Matrix(3,3);
13.     Matrix *matB  = matrix_new_Matrix(3,3);
14. 
    
15.     printf("\t\t|--------------------------|");
16.     printf("\n\t\t|**************************|");
17.     printf("\n\t\t|*** MATRIX ILLUSTRATION **|\n");
18.     printf("\t\t|**************************|\n");
19.     printf("\t\t|--------------------------|\n\n");
20. 
    
21.     matrix_Read(matA);
22.     matrix_Read(matB);
23. 
    
24.     printf("\nmatA + matB:\n");
25.     matrix_Print(matrix_Add(matA,matB));
26. 
    
27.     printf("\nmatA - matB:\n");
28.     matrix_Print(matrix_Sub(matA,matB));
29. 
    
30.     printf("\nmatA * matB:\n");
31.     matrix_Print(matrix_Mult(matA,matB));
32. 
    
33.     printf("\n5 * matA:\n");
34.     matrix_Print(matrix_ScalarMultiply(5,matA));
35. 
    
36.     printf("\nTranspose ofmatA :\n");
37.     matrix_Print(matrix_tran(matA));
38. 
    
39.     printf("\nDeterminant of matA:\t");
40.     printf("%f",matrix_Det(matA));
41. 
    
42.     printf("\nCofactor matrix of matA :\n");
43.     matrix_Print(matrix_Cofactor(matA));
44.     printf("\nInverse matrix of matA:\n");
45.     matrix_Print(matrix_Inverse(matA));
46.     printf("\nMinor matrix of ( i , j) = (1,1) of matA:\n");
47.     matrix_Print(matrix_Minor(matA,1,1));
48. 
    
49.     return 0;
50. }

        Once we wrote a library with generic functions, it easy to develop some projects easily based on our own library. Even though, there are many efficient scientific and mathematical libraries available free of charge.  We learn Programming easily through writing programs like this. 

Happy coding!!!