# Gaussian Elimination C Code

The internet is a wild place full of college homework. It appears that there is a lot of incorrect C code online for Gaussian Elimination with off by one indexing errors. This is an attempt to remedy that. Remember folks, C indexing starts at zero. The pseudocode you will find has indexing that starts at one. Feel free to verify the indexing.

``````#include <stdio.h>

#define n 3
#define my_fabs(a) (((a) < 0) ? -(a) : (a))

void gaussian(float A[n][n], float v[n]) {
float *b = v;
float *x = v;
for (int k = 0; k < n-1; k++) {
// Partial pivot
float cur_max = my_fabs(A[k][k]);
int m = k;
for (int i = k+1; i < n; i++){
float potential_max = my_fabs(A[i][k]);
if (potential_max > cur_max) {
cur_max = potential_max;
m = i;
}
}
if (m != k) {
float b_temp = b[k];
b[k]  = b[m];
b[m]  = b_temp;
for (int j = k; j < n; j++) {
float A_temp = A[k][j];
A[k][j] = A[m][j];
A[m][j] = A_temp;
}
}
// Forward elimination
for (int i = k+1; i < n; i++) {
float factor = A[i][k]/A[k][k];
for (int j = k+1; j < n; j++) {
A[i][j] -= factor*A[k][j];
}
b[i] -= factor*b[k];
}
}
// Back substitution
for (int i = n-1; i >= 0; i--) {
x[i] = b[i];
for (int j = i+1; j < n; j++) {
x[i] -= A[i][j]*x[j];
}
x[i] /= A[i][i];
}
}

int main() {
float A[n][n] = {{ 0,  1,  1},
{ 1, -2,  2},
{ 1,  1, -1}};
float b[n] = {11, 4, 1};
float *v = b;
gaussian(A, v);
float *x = v;
printf("The solution is:\n");
for (int i = 0; i < n; i++) {
printf("x[%d] = %f\n", i, x[i]);
}
return 0;
}
``````