Algorithme de la Méthode de Pivot de Gauss en Langage C :
Le code C ci-dessous illustre une implémentation de cette méthode de résolution ; il est
constitué des blocs suivants :
Le code C ci-dessous illustre une implémentation de cette méthode de résolution ; il est
constitué des blocs suivants :
- les quatre premières fonctions alloc_vecteur, free_vecteur, alloc_matrice et free_matrice sont des fonctions générales permettant de gérer dynamiquement le stockage des matrices et vecteurs ; ces fonctions sont utilisées ici, mais ne font pas partie de l’algorithme du pivot de Gauss,
- les fonctions trouve_pivot, permute_lignes et gauss correspondent à la résolution par le pivot de Gauss,
- la fonction controle n’est là que pour vérifier que la solution trouvée satisfait bien le système initial,
- enfin, les fonctions random et main permettent de tester avec un exemple particulier.
| #include <math.h> |
#include <stdio.h>
#include <stdlib.h>
#include <malloc.h>
#define EPS 1.0e-12
//--------------------------------------------------------------
// Fonction d’allocation d’un vecteur (n)
//--------------------------------------------------------------
double * alloc_vecteur (int n)
{
return (double *)malloc(n*sizeof(double));
}
//--------------------------------------------------------------
// Fonction de désallocation d'un vecteur (n)
//--------------------------------------------------------------
void free_vecteur (double *v)
{
if (v!=NULL) free((void *)v);
}
//--------------------------------------------------------------
// Fonction d'allocation d'une matrice (n,n)
// Remarque : on désalloue en cas d’échec en cours !
//--------------------------------------------------------------
double ** alloc_matrice (int n)
{
double **a;
a=(double **)malloc(n*sizeof(double *));
if (a!=NULL)
{
for (int i=0; i<n; i++)
{
a[i]=(double *)malloc(n*sizeof(double));
if (a[i]==NULL)
{
for (int j=0; j<i; j++)
free((void *)a[j]);
free((void *)a);
return NULL;
}
}
}
return a;
}
/*--------------------------------------------------------------
Fonction de désallocation d'une matrice (n,n)
--------------------------------------------------------------*/
void free_matrice (double **a, int n)
{
if (a!=NULL)
{
for (int i=0; i<n; i++)
if (a[i]!=NULL)
free((void *)a[i]);
free((void *)a);
}
}
//--------------------------------------------------------------
// Fonction de recherche du pivot maximum sur une colonne
// et à partir d'une ligne spécifiée
//--------------------------------------------------------------
int trouve_pivot (double **A, int n, int ligne, int colonne)
{
double v,max;
inti,pivot;
for (i=ligne, max=0.0; i<n; i++)
{
v=fabs(A[i][colonne]);
if (v>max)
{
pivot=i; // pivot identifie la ligne contenant le pivot max.
max=v;
}
}
if (max<EPS) pivot=-1; // pivot trop petit !
return pivot;
}
//--------------------------------------------------------------
// Fonction de permutation de 2 lignes de la matrice
//--------------------------------------------------------------
void permute_lignes (double **A, double *b, intn,int ligne1, int ligne2)
{
double x;
for (int colonne=0; colonne<n; colonne++)
{
x=A[ligne1][colonne];
A[ligne1][colonne]=A[ligne2][colonne];
A[ligne2][colonne]=x;
}
x=b[ligne1];
b[ligne1]=b[ligne2];
b[ligne2]=x;
}
//--------------------------------------------------------------
// Cette fonction cherche la solution du système ax=b
// a est la matrice (n,n), b le second membre (n) et x la solution
// trouvée (n), n est la dimension du système.
// La méthode de Gauss modifie a et b. Pour éviter ce problème, on
// duplique a dans A, b dans B, et on utilise A et B pour la méthode
// du pivot de Gauss. On ne modifie donc ni a ni b.
// Valeur retournée : 0 en cas d'erreur
// 1 en cas de succès
//--------------------------------------------------------------
int gauss (double **a, double *b, double *x, int n)
{
interr,pivot,indpivot,ligne,colonne;
double **A,*B,coef;
A=alloc_matrice(n);
if (A==NULL) return 0; // allocation matrice A
B=alloc_vecteur(n); // allocation vecteur B
if (B==NULL)
{
free_matrice(A,n);
return 0;
}
for (ligne=0; ligne<n; ligne++) // copie de a dans A
{
for (colonne=0; colonne<n; colonne++) // et de b dans B
A[ligne][colonne]=a[ligne][colonne];
B[ligne]=b[ligne];
}
err=1; // code d’erreur
for (pivot=0; pivot<n-1; pivot++)
{
indpivot=trouve_pivot(A,n,pivot,pivot); // recherche du pivot max.
if (indpivot==-1)
{ // problème : pas de pivot satisfaisant
err=0;
break;
}
if (pivot!=indpivot) // permutation lignes si nécessaire
permute_lignes(A,B,n,pivot,indpivot);
for (ligne=1+pivot; ligne<n; ligne++) // calcul de la nouvelle matrice
{
coef=A[ligne][pivot]/A[pivot][pivot];
A[ligne][pivot]=0.0;
for (colonne=1+pivot; colonne<n; colonne++)
A[ligne][colonne]-=A[pivot][colonne]*coef;
B[ligne]-=B[pivot]*coef; // et du nouveau second membre
}
if (fabs(A[pivot][pivot])<EPS) // pivot trop petit : problème dans
{
err=0; // la résolution finale du système
break;
}
}
if (err==1) // si on n’a pas rencontré d’erreur
{
for (ligne=n-1; ligne>=0; ligne--) // calcul des solutions, en remontant
de la
{ // dernière jusqu’à la première
coef=B[ligne];
for (colonne=1+ligne; colonne<n; colonne++)
coef-=A[ligne][colonne]*x[colonne];
x[ligne]=coef/A[ligne][ligne];
}
}
free_matrice(A,n); // désallocation de A
free_vecteur(B); // désallocation de B
return err;
}
//--------------------------------------------------------------
// Cette fonction calcule (et affiche) la différence Ax-b (contrôle)
//--------------------------------------------------------------
void controle (double **a, double *b, double *x, int n)
{
double d;
printf("Vecteur Ax-b :\n");
for (int i=0; i<n; i++)
{
d=-b[i];
for (int j=0; j<n; j++) d+=a[i][j]*x[j];
printf("%9.2le \n ",d);
}
}
//--------------------------------------------------------------
// Cette fonction renvoie un nombre aléatoire entre -range et +range
//--------------------------------------------------------------
double random (double range)
{
return range*(1.0-2.0*(double)rand()/RAND_MAX);
}
//--------------------------------------------------------------
// Exemple d’appel de la fonction gauss
// 1. on alloue dynamiquement a et b (x=b+n)
// 2. la matrice a est aléatoire entre -1 et +1, idem pour b
// 3. on affiche a et b
// 4. on calcule la solution x par la fonction gauss
// 5. on affiche x, puis la différence (ax-b)
// 6. on désalloue a et b
main ()
{
double **a,*b,*x;
int n=5;
a=alloc_matrice(n); if (a==NULL) return 0;
b=alloc_vecteur(2*n);
if (b==NULL)
{
free_matrice(a,n);
return 0;
}
x=b+n;
for (int i=0; i<n; i++)
{
for (int j=0; j<n; j++)
a[i][j]=random(1.0); // Génértion aléatoire de la matrice a
b[i]=random(1.0); // Génértion aléatoire du vecteur b
}
printf("\n Matrice a (GENEREE DE FACON ALEATOIRE):\n");
for (int k=0; k<n; k++)
{
for (int j=0; j<n; j++)
printf("%9.6lf ",a[k][j]);
printf("\n");
}
printf("\n Vecteur b (GENERE DE FACON ALEATOIRE):\n");
for (int m=0; m<n; m++)
printf("%9.6lf ",b[m]);
printf("\n\n");
interr=gauss(a,b,x,n);
printf("Retour de la fonction gauss : %d\n",err);
printf("Solution x :\n");
for (int p=0; p<n; p++)
printf("%9.6lf \n",x[p]);
controle(a,b,x,n); // Permet de vérifier que A.x - b =0
free_matrice(a,n);
free_vecteur(b);
return 0;
}
Aucun commentaire:
Enregistrer un commentaire