codice:
/* ALGORITMO DI JACOBI CON CRITERIO DI ARRESTO A PRIORI*/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
const int size=50;
typedef double matrice [size][size];
typedef double vettore [size];
void contenuto();
int leggidim();
void leggimat(matrice, int);
void leggivet(vettore, int);
double diagonaledominante(matrice,int);
void jacobi(matrice, vettore, vettore, vettore, int);
double norma_max (vettore, int);
main()
{
int n,i,cont=0; matrice A; vettore b,x0,xp1, diff;
double lambda,delta,eps, rho;
contenuto();
n=leggidim();
leggimat(A,n);
lambda=diagonaledominante(A,n);
if (lambda>=1)//Non rispetta la condizione necessaria per la convergenza
{
printf("\n ATTENZIONE!!La matrice A non e'a DIAGONALE STRETTAMENTE DOMINANTE\n");
printf("\n quindi il metodo di Jacobi potrebbe non convergere!Ricompila\n");
printf("\n il programma inserendo una matrice a DIAGONALE STRETTAMENTE DOMINANTE!\n\n");
exit(1);
}
leggivet(b,n);
printf("\n\n Inserire il limite di precisione delle soluzioni. tolleranza=");
scanf("%lf", &eps);
for(i=1; i<=n; i++) x0[i]=0;//Pongo la prima approssimazione del vettore soluzione come il vettore nullo
do
{
jacobi(A,b,x0,xp1, n);
for (i=1;i<=n;i++) diff[i]=x0[i]-xp1[i];
/*Il vettore errore e' il vettore che ha come elementi la differenza tra le componenti del vettore dell'approssimazione
precedente e quelle dell'approssimazione successiva*/
delta=norma_max(diff,n);//L'errore e'il modulo della componente massima del vettore errore
for (i=1;i<=n;i++) x0[i]=xp1[i];
/*Aggiorno il vettore dell'approssimazione precedente che assume il valore del vettore dell'approssimazione
successiva e riinizio Jacobi*/
rho=(log(eps*((1-lambda)/delta))/log(lambda));
cont++;
}
while(cont<rho);
printf("\n\n Il vettore soluzione e'");
for(i=1; i<=n; i++) printf(" %lf ",x0[i]);
printf("\n Sono state eseguite %d iterazioni\n\n",cont);
printf("\n L'errore e': %.15f\n\n",delta);
system("PAUSE");
return 0;
}
void contenuto()
{
printf("\n ALGORITMO DI JACOBI CON IL CRITERIO DI ARRESTO A PRIORI\n\n");
return;
}
int leggidim()
{
int n;
do
{
printf("\n\n Qual è la dimensione della matrice A e del vettore b? n=");
scanf("%d", &n);
}
while(n>size);
return n;
}
void leggimat(matrice A,int n)
{
printf("\n La matrice A e' costituita dai seguenti elementi:");
for(int i=1;i<=n;i++)
for(int j=1;j<=n;j++)
{
printf("\n Inserire l'elemento A[%d][%d]=",i,j);
scanf("%lf",&A[i][j]);
}
return;
}
void leggivet(vettore b,int n)
{
for(int i=1;i<=n;i++)
{
printf("\n Inserire il termine noto b[%d]=",i);
scanf("%lf",&b[i]);
}
return;
}
double diagonaledominante(matrice A,int n)
{
double sum, lambda=0;
for(int i=1;i<=n;i++)
{
sum=0;
for(int j=1;j<=n;j++)
if(i!=j) sum=sum+(fabs(A[i][j]));
sum=sum/fabs(A[i][i]);
if(lambda<sum) lambda=sum;
}
return lambda;
}
void jacobi(matrice A, vettore b,vettore x,vettore y,int n)
{
for(int i=1;i<=n;i++)
{
y[i]=b[i];
for(int j=1;j<=n;j++)
if(j!=i) y[i]=y[i]-A[i][j]*x[j];
y[i]=y[i]/A[i][i];
}
return;
}
double norma_max (vettore diff,int n)
{
double max=fabs(diff[1]);
for(int i=2; i<=n; i++)
{
if (fabs(diff[i])>max)
max=fabs(diff[i]);
}
return max;
}