Calcul de pi avec la bibliothèque gmp

Commenter

Contenu du snippet

Un calcul du nombre PI grâce à la librairie GMP
Vous pouvez trouver un très grand nombre de décimales

Source / Exemple : 


#include <iostream>
#include <cstdlib>
#include <ctime>
#include <gmpxx.h>
#include <iomanip>
#include <vector>
#include <math.h>

using namespace std;
mpz_class MoinsUnPuissanceN(mpz_class n);

inline mpz_class MoinsUnPuissanceN(mpz_class n)
{
     mpz_class retour;
     mpz_class un(1);

     mpz_and(retour.get_mpz_t(),n.get_mpz_t(),un.get_mpz_t());
     if (retour==1)
	  retour=-1;
     else
	  retour=1;
     
     return retour;
}

mpf_class Process_MadhavaLeibniz(mpz_class iterations)
{
     mpf_class retour(0);
     mpf_class deux(2), un(1);

     for (mpz_class i=0; i<=iterations;i++)
     {
	  retour+=(MoinsUnPuissanceN(i)/(deux*i+un));
     }

     return retour*deux*deux;
}

void MadhavaLeibniz()
{
     mpz_class iterations;
     time_t temps;
     unsigned long int precision;

     cout<<endl<<"Nombre d'iterations : ";
     cin>>iterations;
     cout<<endl<<"Nombre de bits pour chaque variable : ";
     cin>>precision;
     mpf_set_default_prec(precision);
     temps=time(NULL);
     cout<<endl<<"Valeur approchée de PI : "<< setprecision(precision)<<Process_MadhavaLeibniz(iterations)<<endl;
     cout<<"en "<<time(NULL)-temps<<" secondes"<<endl;
}

mpf_class Process_Wallis(mpz_class iterations)
{
     mpf_class retour(1);
     mpf_class deux(2), un(1);
     mpf_class puissance,donnee;

     for (mpz_class i=1; i<=iterations; i++)
     {
	  donnee=deux*i+1;
	  mpf_pow_ui(puissance.get_mpf_t(),donnee.get_mpf_t(),2);
	  retour*=(un-un/(puissance));
     }
     return retour*deux*deux;
     
}

void Wallis()
{
     mpz_class iterations;
     time_t temps;
     unsigned long int precision;

     cout<<endl<<"Nombre d'iterations : ";
     cin>>iterations;
     cout<<endl<<"Nombre de bits pour chaque variable : ";
     cin>>precision;
     mpf_set_default_prec(precision);
     temps=time(NULL);
     cout<<endl<<"Valeur approchée de PI : "<< setprecision(precision)<<Process_Wallis(iterations)<<endl;
     cout<<"en "<<time(NULL)-temps<<" secondes"<<endl;

}

mpf_class Process_Viete(mpz_class iterations)
{
     mpf_class retour, encours, addition;
     mpf_class deux(2);
     mpf_sqrt_ui(retour.get_mpf_t(),2);
     encours=retour;
     retour=(deux/retour)*deux;

     for (mpz_class i=1; i<=iterations; i++)
     {
	  addition=2+encours;
	  mpf_sqrt(encours.get_mpf_t(),addition.get_mpf_t());
	  retour*=(deux/encours);
     }
     return retour;
     
}

void Viete()
{
     mpz_class iterations;
     time_t temps;
     unsigned long int precision;

     cout<<endl<<"Nombre d'iterations : ";
     cin>>iterations;
     cout<<endl<<"Nombre de bits pour chaque variable : ";
     cin>>precision;
     mpf_set_default_prec(precision);
     temps=time(NULL);
     cout<<endl<<"Valeur approchée de PI : "<< setprecision(precision)<<Process_Viete(iterations)<<endl;
     cout<<"en "<<time(NULL)-temps<<" secondes"<<endl;

}

mpf_class Process_SuiteLeibniz(mpz_class iterations)
{
     mpf_class retour;
    
     vector<mpf_class> Vecteur1(iterations.get_ui()*2+1);
     vector< vector<mpf_class> > Vecteur2(iterations.get_ui()*2+1,Vecteur1);

     for (mpz_class j=1; j<Vecteur2.size();j++)
     {
	  Vecteur2[0][j.get_ui()]=Process_MadhavaLeibniz(j);
     }

     for (mpz_class i=1; i<Vecteur1.size();i++)
     {
	  for (mpz_class j=1; j<Vecteur2.size();j++)
	  {
	       Vecteur2[i.get_ui()][j.get_ui()]=(Vecteur2[i.get_ui()-1][j.get_ui()]+Vecteur2[i.get_ui()-1][j.get_ui()+1])/2;
	  }
     }
     
     return Vecteur2[iterations.get_ui()][iterations.get_ui()];
     
}

void SuiteLeibniz()
{
     mpz_class iterations;
     time_t temps;
     unsigned long int precision;

     cout<<endl<<"Nombre d'iterations : ";
     cin>>iterations;
     cout<<endl<<"Nombre de bits pour chaque variable : ";
     cin>>precision;
     mpf_set_default_prec(precision);
     temps=time(NULL);
     cout<<endl<<"Valeur approchée de PI : "<< setprecision(precision)<<Process_SuiteLeibniz(iterations)<<endl;
     cout<<"en "<<time(NULL)-temps<<" secondes"<<endl;
}

mpf_class Process_Salamin_Brent(mpz_class iterations)
{
     mpf_class retour, an, bn, tn, pn, b0,temp;
     mpf_class a0(1), p0(1),t0(0.25);
     mpf_sqrt_ui(b0.get_mpf_t(),2);
     b0=1/b0;
    
     for (mpz_class i=1; i<=iterations;i++)
     {
	  an=(a0+b0)/2;
	  temp=a0*b0;
	  mpf_sqrt(bn.get_mpf_t(),temp.get_mpf_t());
	  temp=a0-an;
	  mpf_pow_ui(temp.get_mpf_t(),temp.get_mpf_t(),2);
	  tn=t0-temp*p0;
	  pn=2*p0;

	  a0=an;
	  b0=bn;
	  t0=tn;
	  p0=pn;

     }
     temp=an+bn;
     mpf_pow_ui(retour.get_mpf_t(),temp.get_mpf_t(),2);
     retour/=(4*tn);

     return retour;
     
}

void Algo_Salamin_Brent()
{
     mpz_class iterations;
     time_t temps;
     unsigned long int precision;

     cout<<endl<<"Nombre d'iterations : ";
     cin>>iterations;
     cout<<endl<<"Nombre de bits pour chaque variable : ";
     cin>>precision;
     mpf_set_default_prec(precision);
     temps=time(NULL);
     cout<<endl<<"Valeur approchée de PI : "<< setprecision(precision)<<Process_Salamin_Brent(iterations)<<endl;
     cout<<"en "<<time(NULL)-temps<<" secondes"<<endl;
}

mpf_class Process_BBP(mpz_class iterations)
{
     mpf_class retour(0);

     if (iterations<pow(2,sizeof(int)*8))
     {
#define avecentier
     }
     else
     {
#define avecGMP
     }

#if defined avecentier
     for (unsigned int i=0; i<=iterations.get_ui(); i++)
#elif defined avecGMP
     for( mpf_class i=0; i<iterations;i++)
#endif
     {
	  mpf_class ajout, temp8i, temp/*, huit(8), quatre(4), deux(2), cinq(5), six(6), un(1)*/, seize(16);

	  temp8i=8*i;
	  temp=temp8i+1;
	  ajout=4/temp;

	  temp=temp8i+4;
	  ajout-=(2/temp);

	  temp=temp8i+5;
	  ajout-=(1/temp);

	  temp=temp8i+6;
	  ajout-=(1/temp);

	 
#if defined avecentier

	  mpf_pow_ui(temp.get_mpf_t(),seize.get_mpf_t(),i);
#elif defined avecGMP

	  temp=1;

	  for (mpf_class j=0; j<i; j++)
	  {
	       mpf_mul_ui(temp.get_mpf_t(),temp.get_mpf_t(),16);
	  } 

#endif

	  ajout*=(1/temp);

	  retour+=ajout;

     }

     return retour;

}

void BBP()
{
mpz_class iterations;
     time_t temps;
     unsigned long int precision;

     cout<<endl<<"Nombre d'iterations : ";
     cin>>iterations;
     cout<<endl<<"Nombre de bits pour chaque variable : ";
     cin>>precision;
     mpf_set_default_prec(precision);
     temps=time(NULL);
     cout<<endl<<"Valeur approchée de PI : "<< setprecision(precision)<<Process_BBP(iterations)<<endl;
     cout<<"en "<<time(NULL)-temps<<" secondes"<<endl;

}

int main(void)
{
     int choix;

     do
     {
	  cout<<endl<<"Calcul du nombre PI par différentes méthodes"<<endl<<endl;
	  cout<<"1 : Par la formule de Madhava-Leibniz"<<endl;
	  cout<<"2 : Par la formule de Wallis"<<endl;
	  cout<<"3 : Par la formule de Viete"<<endl;
	  cout<<"4 : Par la suite de Leibniz"<<endl;
	  cout<<"5 : Par l'algorithme de Salamin et Brent"<<endl;
	  cout<<"6 : Par la formule BBP"<<endl;
	  cout<<"0 : Arret du programme"<<endl;
	  cout<<endl<<"Votre choix ?"<<endl;
	  cin>>choix;

	  switch (choix)
	  {
	  case 1:
	       MadhavaLeibniz();
	       break;
	 
	  case 2:
	       Wallis();
	       break;
	  case 3:
	       Viete();
	       break;
	  case 4:
	       SuiteLeibniz();
	       break;
	  case 5:
	       Algo_Salamin_Brent();
	       break;
	  case 6:
	       BBP();
	       break;
	  }

     }
     while (choix!=0);

     return EXIT_SUCCESS;
}

Conclusion :


Il faut récupérer la librairie sur ce site : http://gmplib.org/
et compiler le projet en mettant les liens de la librairie.
De cette façon par exemple avec le compilateur gcc : g++ -o PI_GMP -lgmpxx -lgmp -g *.cpp

Pour calculer les 128000 premiers chiffres, j'ai mis environ 2 heures sur un core i7 avec l'algorithme de Salamin et Brent.
Les chiffres sont disponibles sur ce site : http://jlsigrist.com/pi.html

A voir également

S'inscrire maintenant

Vous n'êtes pas encore membre ?

inscrivez-vous, c'est gratuit et ça prend moins d'une minute !

Les membres obtiennent plus de réponses que les utilisateurs anonymes.

Le fait d'être membre vous permet d'avoir un suivi détaillé de vos demandes et codes sources.

Le fait d'être membre vous permet d'avoir des options supplémentaires.

Du même auteur (cs_lann)