Lab 03: Gauss Elimination Function in C (Version 2)


The following code implements the gauss_eliminate() function to solve systems of equations. It does not have the code for back substitution which is left as an exercise! The main purpose of this code is to illustrate how 2D arrays are passed as parameters to functions in the C language.

One data file with data for 3 equations in 3 unknowns and another with data for 10 equations in 10 unknowns is given for testing the code. At the end of gauss_eliminate(), the resulting matrix must be strictly in upper diagonal form. Once the back substitution method is implemented, the solutions can be verified.

In the case of three unknowns, the values are $x_1=19, x_2=6, x_3=4$.


Data Files for gauss_eliminate() function

  • data.txt
    3
    1, 1, -1
    1, 2, 2
    -1, 5, -1
    21, 39, 7
  • data10.txt
    10
    12, -1, -1, -1, -1, -1, -1, -3, -1, -1
    -1, 10, -1, -1, -1, -3, -1, -1, -1, -1
    -1, -1, 7, 0, 0, 0, -2, 1, 1, 1
    -2, 1, -1, 15, -5, -1, 0, -2, 3, -1
    1, 1, 1, 1, -22, 1, 1, 1, -4, 1
    -3, -3, 0, -1, 0, 9, 0, 0, 0, 0
    5, -3, -2, -2, 0, 3, 15, -1, -3, -1
    1, -2, 1, -2, 1, -2, 1, 25, -1, 1
    1, 1, 1, 1, 1, -1, -1, -1, 18, -1
    -1, -1, -1, -1, 1, 1, 1, 1, 1, 10
    70, 0, 32, 5, 22, 0, 106, 30, 5, 69
In [ ]:
#include <stdio.h>
#include <stdlib.h>

void print_array(float **arr, int nrows, int ncols) {
    int i, j;
    
    for (i=0; i<nrows; i++) {
        for (j=0; j<ncols; j++)
            printf("%10.2f  ", arr[i][j]);
        printf("\n");
    }
}
/**************
*  This is a simpler version of the Gauss Elimination method 
*  where the pivot element is the diagonal element
*  The only sanity check done is to see if the pivot = 0
***************/
void gauss_eliminate(float **ip, int nrows, int ncols)
{
    int r_piv=0, c_piv=0;
    int crow, ccol;
    float mult_fact;
    
    while ((r_piv < nrows-1) && (c_piv < ncols-1)) {
        if (ip[r_piv][c_piv] == 0.0)
            c_piv += 1;
        else {
            for (crow=r_piv+1; crow<nrows; crow++) {
                mult_fact = ip[crow][c_piv] / ip[r_piv][c_piv];
                ip[crow][c_piv] = 0.0;
                for (ccol=c_piv+1;ccol<ncols; ccol++)
                    ip[crow][ccol] -= (ip[r_piv][ccol] * mult_fact);
            }
            r_piv += 1;
            c_piv += 1;
        }
    }
}
/****************
*   This function creates a new matrix with the given number
*   of rows and columns and returns a pointer to it. The matrix
*   is assumed to contain float values.
*****************/
float **create_matrix(int nrows, int ncols)
{
    float **A;
    int i;
    
    if ((A = (float **)malloc(nrows * sizeof(float *))) == NULL)
        return(NULL);
    
    for (i=0; i<nrows; i++) {
        A[i] = (float *)malloc(ncols * sizeof(float));
        if (A[i] == NULL)
            return(NULL);
    }
    return(A);
}

/********************
*  Read the number of unknowns, the coefficient matrix A and the
*  RHS matrix B from a file; do Gaussian Elimination and output
*  the upper triangular matrix.
*********************/
int main(void)
{
    float **A, *B, **C;
    int i, j, SIZE;
    FILE *ip;
    
    if ((ip = fopen("data.txt", "r")) == NULL) {
        printf("Error in reading input file\n");
        exit(1);
    }
    
    fscanf(ip, "%d", &SIZE);       /* Read number of unknowns */
    
    A = create_matrix(SIZE, SIZE);
    B = (float *)malloc(SIZE * sizeof(float));
    C = create_matrix(SIZE, SIZE+1);
    
    printf("Reading values into A matrix\n");
    for (i=0; i<SIZE; i++)
        for (j=0; j<SIZE; j++)
        fscanf(ip, "%f,", &A[i][j]);
    
    printf("Reading values of B matrix\n");
    for (i=0; i<SIZE; i++)
        fscanf(ip, "%f,", &B[i]);
    fclose(ip);
    
    for (i=0; i<SIZE; i++) {
        for (j=0; j<SIZE; j++) 
            C[i][j] = A[i][j];
        C[i][SIZE] = B[i];
    }
    /* Print augmented matrix to check if values are read
     * correctly from the data file                        */
    print_array(C, SIZE, SIZE+1);
    
    gauss_eliminate(C, SIZE, SIZE+1);
    
    printf("Reduced Upper Triangular Matrix\n");
    print_array(C, SIZE, SIZE+1);
    
    exit(0);
}


THE END