#include <iostream>
#include <stdlib.h>
#include <math.h>
#define MaxD 10
#define Nmax 15

using namespace std;

float Matrix[MaxD][MaxD];
float X[2][MaxD];
float b[MaxD];
float e;
int N;

float Sumof(int CurLine,int MLine)
{
    float RetVal = 0;
    for(int i=0; i<N; i++)
    {
        if(MLine != i)
            RetVal += Matrix[MLine][i]*X[CurLine][i];
    }
    return RetVal;
}

float CompareX()
{
    float temp = 0;
    float RetVal = fabs(X[0][0]-X[1][0]);
    for(int i=1; i<N; i++)
    {
        temp = fabs(X[0][i]-X[1][i]);
        if(RetVal<temp)
            RetVal = temp;
    }
    return RetVal;
}
void Jacobi(int Xinit)
{
    int i=0,j=0;
    for(i=0; i<MaxD; i++)
    {
        X[0][i] = Xinit;
        X[1][i] = Xinit;
    }
    
    for(i=0; i<Nmax; i++)
    {
        for(j=0; j<N; j++)
        {
            if(Matrix[j][j]==0)
            {
                printf("JACOBI迭代失效");
                return;
            }
            X[i%2][j] = (b[j] - Sumof((i+1)%2,j))/Matrix[j][j];
        }
       
        printf("第%d步迭代結果:\n",i);
        for(int k=0; k<N; k++)
        {
            printf("X[%d]=%5f\n",k,X[i%2][k]);
        }
        printf("--------------------------------------------\n");
        if(CompareX()<e)
        {
            printf("最終迭代結果:\n");
            for(int k=0; k<N; k++)
            {
                printf("X[%d]=%5f\n",k,X[i%2][k]);
            }
            return;
        }
    }
    printf("JACOBI迭代失效");
}

int main(int argc, char *argv[])
{
    int i=0,j=0;
    e = 0.01;
    N = 3;
    do
    {
        printf("請輸入矩陣維數 不大于10\n");
        scanf("%d",&N);
    }while(N > 10); 
    
    printf("請輸入允許誤差 e=");
    scanf("%f",&e);  
    
    printf("請輸入矩陣元素:\n");
    for(i=0; i<N; i++)
    {
        for(j=0; j<N; j++)
        {
            printf("Matrix[%d][%d]=",i,j);
            scanf("%f",&Matrix[i][j]);
        }
    }
    
    printf("請輸入常數項向量\n");
    for(i=0; i<N; i++)
    {
        printf("b[%d]=",i);
        scanf("%f",&b[i]);
    }
    
    for(i=0; i<N; i++)
    {
        for(j=0; j<N; j++)
        {
            printf("Matrix[%d][%d]=%5f    ",i,j,Matrix[i][j]);
        }
        printf("\n");
    }
    
    Jacobi(1);
    system("PAUSE");	
    return 0;
}