METODO DI JACOBI E CRITERI DI FERMATA....

[Fravilla1]

Ciao a tutti!!!!!! io ho fatto questo codice del metodo di jacobi con criterio di fermata a posteriori. solo che le soluzioni dei sistema non vengono... Dove ho sbagliato? grazie mille

/*METODO DI JACOBI*/

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

const int size=50;
typedef double matrice[size][size];
typedef double vettore[size];


void leggimatrice(matrice,int);
void leggivettore(vettore,int);
void jacobi(matrice,vettore,vettore,vettore,int);
double diagonaledominante(matrice,int);
double normamax(vettore, vettore, int);
void stampa(vettore,int);

main()
{
      matrice A; vettore b, vettk, vettkp1; double delta;
      
      int n, k=0; double toll, lambda, rho;
      
      printf("\n dammi la dimensione della matrice  n=");scanf("%d",&n);
      
      leggimatrice(A,n);
      
      lambda=diagonaledominante(A,n);//lamda=max-i(somm a[i][j])/a[i][i]
      
      if(lambda>=1)
                   {
                      printf("\n la matrice non e'a diagonale dominante!non puoi usare jacobi\n\n");
                      exit(1);
                   }
      
      leggivettore(b,n);
      
      printf("\n\n Inserire il limite di precisione delle soluzioni.        tolleranza=");
      scanf("%lf", &toll);
      
      for(int i=1; i<=n; i++)         vettk[i]=0;
      
      for(int i=1; i<=n; i++)    vettkp1[i]=0;
      
      delta=normamax(vettkp1, vettk, n);
      
      rho=(log(toll*((1-lambda)/delta))/log(lambda));
      
      do(k++);
      while(k<=rho);        
      
      while(delta>=toll)
             {
                        for(int j=1; j<n; j++)
                                {
                                        for(int i=1;i<=n;i++)   vettk[i]=vettkp1[i];
                   
                                        jacobi(A,b,vettk,vettkp1,n);
                                        
                                        delta=fabs(vettkp1[j]-vettk[j]);
                                 }
             }
     
     stampa(vettkp1,n);
     
     system("PAUSE");
     return 0;
}


void leggimatrice(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 leggivettore(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;
}


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 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;
}


double normamax(vettore vettkp1, vettore vettk, int n)
{
       
       double max=(fabs(vettkp1[1]-vettk[1]));
       
       for(int i=2;i<=n;i++)
               if(max<(fabs(vettkp1[i]-vettk[i])))     max=(fabs(vettkp1[i]-vettk[i]));
               
       return max;
}

        
void stampa(vettore y,int n)
{
     printf("\n Le soluzioni del sistema sono:");
     
     for(int i=1;i<=n;i++)   printf("\t%lf",y[i]);
     
     return;
}

[Fravilla1]

per favore, c'è qualcuno che conosce questo metodo? è importantissimo!!!!!!!!!!

[Riuzaki]

io farei cosi ma non so se è giusto ...ho scritto di getto senza controllare il funzionamento al compilatore:

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

#define EPS 0,00001   //tolleranza 
#define MAXINTER 100000  //numero max di interazioni 
using namespace std; 

const int MAX_DIM = 10; 

void leggi(double A[MAX_DIM][MAX_DIM],int &dim); 
void leggi_b(int dimensione, double b[]); 
void stampa(double A[MAX_DIM][MAX_DIM], double x[MAX_DIM], int dim); 

int jacobi(double A[MAX_DIM][MAX_DIM], double b[MAX_DIM],int dim) 
{ 
      int i,j; 
      int flag=0; 
      double sum; 
      double X[MAX_DIM]; 
      for( i=0; i<dim; i++){ 
          X[i]=b[i]/A[i][i]; 
          }    
          while(EPS && flag < MAXINTER) 
          { 
                
                for(i=0;i<dim;i++){ 
                    b[i]=X[i]; 
                    sum=b[i]; 
                    for(j=0;j<dim;++j) 
                    if(j!=i) 
                    sum-=A[i][j]*b[j]; 
                    X[i]=sum/A[i][i]; 
                    for (i=0;i<dim;i++) 
                    if(fabs(X[i]- b[i])>fabs(X[i])*EPS) 
                    ++flag; 
                    }          
           }                      
                
                
} 


void leggi( double A[MAX_DIM][MAX_DIM], int &dim) 
{ 
       int i,j; 
       cout << "Inserire la dimensione della matrice: "; 
       cin >> dim; 
       for( i=0; i<dim; i++) 
          for( j=0; j<dim; j++) { 
             cout << "A(" << i+1 << "," << j+1 << ")="; 
             cin >> A[i][j]; 
          } 
} 


void leggi_b( int dim, double b[MAX_DIM]) 
{ 
       int i; 
       for( i=0; i<dim; i++) { 
          cout << "b(" << i+1 << ")="; 
          cin >> b[i]; 
       } 
} 


void stampa(double A[MAX_DIM][MAX_DIM], double term_noti[MAX_DIM], int dim) 
{ 
       int i,j; 
       cout << "Matrice\n"; 
       for(i=0;i<dim;i++) { 
          for(j=0;j<dim;j++) 
              cout << "\tA(" 
              << i+1 << "," 
              << j+1 << ")=" 
              <<A[i][j]; 
              cout << "\n"; 
       } 
} 


void stampa_sol(double x[MAX_DIM], int dim) 
{ 
       int i; 
       cout << "Vettore soluzioni\n"; 
       for(i=0;i<dim;i++) {            
           cout << "\tx(" << i+1 << ")=" << x[i]; 
           cout << "\n"; 
       } 
} 


main() 
{ 
       double A[MAX_DIM][MAX_DIM],b[MAX_DIM], X[MAX_DIM]; 
       int dim; 
       int success; 

       cout << "\n\n\n"; 
       leggi(A,dim); 
       leggi_b(dim,b); 
       stampa(A,b,dim); 
       jacobi(A,b,dim); 
       cout << "\n\nIl vettore delle soluzioni e’:\n"; 
       stampa_sol(X,dim); 
        

       system ("PAUSE"); return 0; 
}