gaussian elimination

gaussianelimination.c

/*****************************************************
 * 2017/10/27
 * input format
 * $./a.out
 * 4 8 16 24 32 2 7 12 17 6 17 32 59 7 22 46 105 160 70 198 291
 * ***************************************************/

#include <stdio.h>
#include <string.h>

int calculating_times = 0;

void inputMatrix(int row, int column, double data[row][column]);
void inputVector(int row, double data[row]);
void defaultPrint(int row, int column, double matrix[row][column], double data[row]);
void calculationPrint();
void unitMatrix(int row, int column, double data[row][column]);
void gaussianElimination(int num, int component, int row, int column, double data[row][column], double matrix[row][column]);
void multiplicationMA(int num, int size, int row, int column, double data[row][column], double matrix[row][column], double result[row][column]);
void multiplicationMb(int num, int size, int row, int column, double data[row], double matrix[row][column], double result[row]);
void answerPrint(int size, int row, int column, double matrix[row][column], double data[row]);

int main(void)
{
    int row, column, MATRIX_SIZE;
    int output_times = 1;
    int i, j, k;
    int diagonal_component;

    scanf("%d", &MATRIX_SIZE);

    double A[MATRIX_SIZE][MATRIX_SIZE];
    double b[MATRIX_SIZE];
    double M[MATRIX_SIZE][MATRIX_SIZE];
    double MA[MATRIX_SIZE][MATRIX_SIZE];
    double Mb[MATRIX_SIZE];

    /* input matrix */
    inputMatrix(MATRIX_SIZE, MATRIX_SIZE, A);
    inputVector(MATRIX_SIZE, b);
    /* output default value */
    defaultPrint(MATRIX_SIZE, MATRIX_SIZE, A, b);

    for (diagonal_component = 0; diagonal_component < MATRIX_SIZE - 1; diagonal_component++)
    {
        /* unit matrix */
        unitMatrix(MATRIX_SIZE, MATRIX_SIZE, M);
        /* Gaussian Elimination */
        gaussianElimination(output_times, diagonal_component, MATRIX_SIZE, MATRIX_SIZE, A, M);
        /* multiplication MA */
        multiplicationMA(output_times, MATRIX_SIZE, MATRIX_SIZE, MATRIX_SIZE, A, M, A);
        /* multiplication MB */
        multiplicationMb(output_times, MATRIX_SIZE, MATRIX_SIZE, MATRIX_SIZE, b, M, b);
        output_times++;
    }
    answerPrint(MATRIX_SIZE, MATRIX_SIZE, MATRIX_SIZE, A, b);
    return 0;
}

/* input matrix */
void inputMatrix(int row, int column, double data[row][column])
{
    int i, j;
    for (i = 0; i < row; i++)
    {
        for (j = 0; j < column; j++)
        {
            scanf("%lf", &data[i][j]);
        }
    }
}

/* input vector */
void inputVector(int row, double data[row])
{
    int i;
    for (i = 0; i < row; i++)
    {
        scanf("%lf", &data[i]);
    }
}

/* output default value */
void defaultPrint(int row, int column, double matrix[row][column], double data[row])
{
    int i, j;
    printf("Matrix A is\n");
    for (i = 0; i < row; i++)
    {
        for (j = 0; j < column; j++)
        {
            printf("%.3f ", matrix[i][j]);
        }
        printf("\n");
    }

    printf("Vector b is\n");
    for (i = 0; i < row; i++)
    {
        printf("%.3f ", data[i]);
    }
    printf("\n\n");
    printf("--------------------------------------------\n");
}

/* output calculating times */
void calculationPrint()
{
    printf("Num. of Mul./Div. = %d", calculating_times);
    printf("\n\n");
}

/* unit matrix */
void unitMatrix(int row, int column, double data[row][column])
{
    int i, j;
    for (i = 0; i < row; i++)
    {
        for (j = 0; j < column; j++)
        {
            if (i == j)
                data[i][j] = 1;
            else
            {
                data[i][j] = 0;
            }
        }
    }
}

/* Gaussian Elimination */
void gaussianElimination(int num, int component, int row, int column, double data[row][column], double matrix[row][column])
{
    int i, j;
    for (i = component; i < row; i++)
    {
        if (component == i)
        {
            matrix[i][component] = data[i][component] / data[component][component];
        }
        else
        {
            matrix[i][component] = -data[i][component] / data[component][component];
        }
        calculating_times++;
    }

    printf("M%d is\n", num);
    for (i = 0; i < row; i++)
    {
        for (j = 0; j < column; j++)
        {
            printf("%.3f ", matrix[i][j]);
        }
        printf("\n");
    }
    calculationPrint();
}

/* multiplication MA */
void multiplicationMA(int num, int size, int row, int column, double data[row][column], double matrix[row][column], double result[row][column])
{
    int i, j, k, sum = 0;
    for (i = 0; i < size; i++)
    {
        for (j = 0; j < size; j++)
        {
            for (k = 0; k < size; k++)
            {
                sum += matrix[i][k] * data[k][j];
                calculating_times++;
            }
            result[i][j] = sum;
            sum = 0;
        }
    }

    printf("M%dA is\n", num);
    for (i = 0; i < row; i++)
    {
        for (j = 0; j < column; j++)
        {
            printf("%.3f ", result[i][j]);
        }
        printf("\n");
    }
    calculationPrint();
}

/* multiplication MB */
void multiplicationMb(int num, int size, int row, int column, double data[row], double matrix[row][column], double result[row])
{
    int i, j, sum = 0;
    for (i = 0; i < size; i++)
    {
        for (j = 0; j < size; j++)
        {
            sum += matrix[i][j] * data[j];
            calculating_times++;
        }
        result[i] = sum;
        sum = 0;
    }

    printf("B%d is\n", num);
    for (i = 0; i < row; i++)
    {
        printf("%.3f ", result[i]);
    }
    printf("\n");

    calculationPrint();
    printf("--------------------------------------------\n");
}

/* output answer */
void answerPrint(int size, int row, int column, double matrix[row][column], double data[row])
{
    int i, j;
    double x[size];
    printf("Answers after  = %d mul/div is\n", calculating_times);

    for (i = 0; i < size; i++)
    {
        x[i] = 0;
    }

    for (i = size - 1; i >= 0; i--)
    {
        for (j = 0; j < size; j++)
        {
            data[i] = data[i] - matrix[i][j] * x[j];
        }
        x[i] = data[i]/matrix[i][i];
    }

    for (i = 0; i < row; i++)
    {
        printf("x[%d] = %.3f \n",i+1 ,x[i]);
    }
    printf("\n");
}