2015年12月29日 星期二

Jacobi Iterative (雅可比迭代) C語言實現



#define _CRT_SECURE_NO_WARNINGS
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

void main(){
int i, j, k;
int n, N;
double **arrA, *arrB, *arrXO, *arrX, TOL, tem;

printf("***********  input  **********\n");

printf("n=");
scanf("%d", &n);

arrA = (double **)malloc(sizeof(double) * n * (n + 1));
for (i = 0; i < n; i++)
arrA[i] = (double *)arrA + n * (i + 1);
arrB = (double *)malloc(sizeof(double) * n);
arrXO = (double *)malloc(sizeof(double) * n);
arrX = (double *)malloc(sizeof(double) * n);

printf("\n");
for (i = 0; i < n; i++){
for (j = 0; j < n; j++){
printf("A[%d][%d]=", i, j);
scanf("%lf", &arrA[i][j]);
}
printf("\n");
}

printf("\n");
for (i = 0; i < n; i++){
printf("B[%d]=", i);
scanf("%lf", &arrB[i]);
}

printf("\n");
for (i = 0; i < n; i++){
printf("XO[%d]=", i);
scanf("%lf", &arrXO[i]);
}

printf("\nTOL=");
scanf("%lf", &TOL);

printf("\nN=");
scanf("%d", &N);

/*printf("***********  information  **********\n");
printf("\nA=\n");
for (i = 0; i < n; i++){
for (j = 0; j < n; j++)
printf(" %lf", arrA[i][j]);
printf("\n");
}
printf("\nB=\n");
for (i = 0; i < n; i++)
printf(" %lf", arrB[i]);
printf("\nXO=\n");
for (i = 0; i < n; i++)
printf(" %lf", arrXO[i]);*/

printf("\n\n***********  Jacobi Iterative (雅可比迭代?)  **********\n");

k = 0;
while (k <= N){
for (i = 0; i < n; i++){
arrX[i] = 0;
for (j = 0; j < n; j++){
if (i == j)
continue;
arrX[i] -= arrA[i][j] * arrXO[j];
}
arrX[i] = (arrX[i] + arrB[i]) / arrA[i][i];
}
tem = 0;
for (i = 0; i < n; i++)
tem += abs(arrX[i] - arrXO[i]);

if (tem < TOL){
printf("\n\nX =");
for (i = 0; i < n; i++)
printf(" %lf", arrX[i]);
break;
}
k++;
for (i = 0; i < n; i++)
arrXO[i] = arrX[i];

if (k > N)
printf("\n\nMaximum number of iterations exceeded\n\n");
}

free(arrA);
free(arrB);
free(arrXO);
free(arrX);

printf("\n\n***********  end  **********\n");
system("pause");
}


沒有留言:

張貼留言