[수치계산] 연립선형방정식의 해 - Gauss Seidel 반복법

Gauss Seidel 반복법

● Gauss 소거법 & Gauss-Jordan 소거법

- 방정식 개수가 수십개인 작은 크기의 연립방정식에서 실근에 매우 근접한 해를 구할 수 있다.

- O(n^3)의 덧셈/곱셈 연산을 필요로 하기 때문에 매우 큰 규모의 연립방정식을 풀기 힘들다

● Gauss-Seidel 반복법

- 큰 규모의 연립방정식 근사해를 빠르게 구할 수 있다

- 과정

1. 초기 근사해를 설정

2. 1의 해를 이용해 더 나은 근사해를 구한다

3. 허용가능한 오차범위보다 작아질때까지 2를 반복한다

#include <stdio.h>
#include <math.h>

#define n 4
#define THRESHOLD 0.001

void main(){
    float a[5][5], b[5], x_now[5], x_previous[5];
    int i, j, k;
    float e;

    printf("Input A : \n");

    for(i = 1; i <= 4; i++)
        scanf("%f %f %f %f", &a[i][1], &a[i][2], &a[i][3], &a[i][4]);

    printf("Input b : \n");
    scanf("%f %f %f %f", &b[1], &b[2], &b[3], &b[4]);

    x_previous[1] = 0;
    x_previous[2] = 0;
    x_previous[3] = 0;
    x_previous[4] = 0;

    printf(" 0 : x1 = %8.3f \t x2 = %8.3f \t x3 = %8.3f \t x4 = %8.3f\n",
            x_previous[1], x_previous[2], x_previous[3], x_previous[4]);

    k = 0;
    while(1){
        k++;

        for(i = 1; i <= n; i++){
            x_now[i] = 0;
            for(j = 1; j <= i - 1; j++)
                x_now[i] = x_now[i] - a[i][j] * x_now[j];

            for(j = i + 1; j <= n; j++)
                x_now[i] = x_now[i] - a[i][j] * x_previous[j];

            x_now[i] = (x_now[i] + b[i]) / a[i][i];
        }

        printf("%2d : x1 = %8.3f \t x2 = %8.3f \t x3 = %8.3f \t x4 = %8.3f\n",
                k, x_now[1], x_now[2], x_now[3], x_now[4]);

        e = 0;

        for(i = 1; i <= n; i++){
            e = e + fabs(x_now[i] - x_previous[i]);
        }

        e = e / n;
        if(e <= THRESHOLD)  break;
        for(i = 1; i <= n; i++)
            x_previous[i] = x_now[i];
    }

    printf("\nResult : \n");
    for(i = 1; i <= n; i++)
        printf("x%d = %8.3f\n", i, x_now[i]);
}