Experimento: cálculo de cifras del número Pi

Me se ha ocurrido medir el rendimiento del ordenador, simplemente con un programa que calcula las cifras del número pi que están en x posición decimal.

El programa no lo he creado yo, no tengo suficiente nivel para ello, pero sí lo he compilado para GNU/Linux y Windows, y poder comparar el rendimiendo de las dos plataformas en el mismo ordenador.

GNU/Linux (compilador GCC 4.6.1):
En primer lugar, lo he compilado con la optimización por defecto, luego se verá que con la máxima optimización el programa calcula mucho más rápido.

g++ -lm pi.c -o pisinopti
Luego con la máxima optimización:
g++ -lm -O3 pi.c -o piopti
He aquí la diferencia, siendo el código fuente el mismo:

juan@siemens:~/Compilar/Pi$ ./piopti 5000
Decimal digits of pi at position 5000: 156951623
12478.83699997328 milliseconds

juan@siemens:~/Compilar/Pi$ ./pisinopti 5000
Decimal digits of pi at position 5000: 156951623
37960.17299999949 milliseconds

El programa optimizado lo consigue 12s y el no optimizado tarda más del doble, casi 38s. Los dos ejecutables pesan 11,4 KiB. El optimizado utiliza 64K de memoria, y el otro, 60K.

Windows (compilador GCC 4.6.1):

Hacemos lo mismo que en GNU/Linux en MinGW32. Los resultados (12,5s y 39,2s respectivamente), ligeramente más altos que en GNU/Linux, pero despreciable la diferencia.

Windows (compilador Visual Studio 2010):

He probado también con el compilador de Microsoft porque pensaba que era más eficiente, pero por lo que vemos, resulta ser bastante inferior en rendimiendo (tarda el doble en el cálculo, incluso con la máxima optimización), aunque los ejecutables que genera son más pequeños (8KiB VS2010 vs. 56 KiB GCC) y utilizan un poco de menos memoria durante la ejecución (192KiB VS2010 vs. 232 KiB GCC).

piopti = 12,5 segundos (g++ -O3)
pisin = 39,2 segundos (g++)
pivsopti = 24,1 segundos (flag /Ox)
pivssin = 24,1 segundos (flag /O2)

El programa venía sin contador de tiempo de ejecucación, por lo que le tenido que agregar uno como mejor he podido básandome en el código de la siguiente página: http://cmasomenos.blogspot.com.es/2008/03/medir-el-tiempo-de-una-rutina.html
Utilizando gettimeofday(), que no es muy preciso en Windows, pero para esto ya nos vale.

Para compilar el código en Visual Studio 2010… además he tenido que meter código para que haga la función POSIX gettimeofday(), ya que la marca como desconocida. A parte de otro error en el código original de cálculo con las raíces cuadradas y los logaritmos, el error C2668, “ambiguous call to overloaded function”, que he solucionado tomando como referencia el ejemplo d.

El programa se podría compilar para Android e inclusive iOS, y comprobar el rendimiendo de los teléfonos inteligentes y ver si algún teléfono tiene más potencia que mi PC antiguo (Pentium 4 a 2,53 GHz).

He aquí el código con el contador:

/*
 * Computation of the n'th decimal digit of \pi with very little memory.
 * Written by Fabrice Bellard on January 8, 1997.
 * 
 * We use a slightly modified version of the method described by Simon
 * Plouffe in "On the Computation of the n'th decimal digit of various
 * transcendental numbers" (November 1996). We have modified the algorithm
 * to get a running time of O(n^2) instead of O(n^3log(n)^3).
 * 
 * This program uses mostly integer arithmetic. It may be slow on some
 * hardwares where integer multiplications and divisons must be done
 * by software. We have supposed that 'int' has a size of 32 bits. If
 * your compiler supports 'long long' integers of 64 bits, you may use
 * the integer version of 'mul_mod' (see HAS_LONG_LONG). 
 */

#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <time.h>
#include <sys/time.h>

/* uncomment the following line to use 'long long' integers */
/* #define HAS_LONG_LONG */

#ifdef HAS_LONG_LONG
#define mul_mod(a,b,m) (( (long long) (a) * (long long) (b) ) % (m))
#else
#define mul_mod(a,b,m) fmod( (double) a * (double) b, m)
#endif

/* retorna "a - b" en segundos */
double timeval_diff(struct timeval *a, struct timeval *b)
{
 return
 (double)(a->tv_sec + (double)a->tv_usec/1000000) -
 (double)(b->tv_sec + (double)b->tv_usec/1000000);
}

/* return the inverse of x mod y */
int inv_mod(int x,int y) {
 int q,u,v,a,c,t;

 u=x;
 v=y;
 c=1;
 a=0;
 do {
 q=v/u;

 t=c;
 c=a-q*c;
 a=t;

 t=u;
 u=v-q*u;
 v=t;
 } while (u!=0);
 a=a%y;
 if (a<0) a=y+a;
 return a;
}

/* return (a^b) mod m */
int pow_mod(int a,int b,int m)
{
 int r,aa;

 r=1;
 aa=a;
 while (1) {
 if (b&1) r=mul_mod(r,aa,m);
 b=b>>1;
 if (b == 0) break;
 aa=mul_mod(aa,aa,m);
 }
 return r;
}

/* return true if n is prime */
int is_prime(int n)
{
 int r,i;
 if ((n % 2) == 0) return 0;

 r=(int)(sqrt(n));
 for(i=3;i<=r;i+=2) if ((n % i) == 0) return 0;
 return 1;
}

/* return the prime number immediatly after n */
int next_prime(int n)
{
 do {
 n++;
 } while (!is_prime(n));
 return n;
}

int main(int argc,char *argv[])
{
 int av,a,vmax,N,n,num,den,k,kq,kq2,t,v,s,i;
 double sum;
 struct timeval t_ini, t_fin;
 double secs;

 if (argc<2 || (n=atoi(argv[1])) <= 0) {
 printf("This program computes the n'th decimal digit of \\pi\n"
 "usage: pi n , where n is the digit you want\n"
 );
 exit(1);
 }

 gettimeofday(&t_ini, NULL);

 N=(int)((n+20)*log(10)/log(2));

 sum=0;

 for(a=3;a<=(2*N);a=next_prime(a)) {

 vmax=(int)(log(2*N)/log(a));
 av=1;
 for(i=0;i<vmax;i++) av=av*a;

 s=0;
 num=1;
 den=1;
 v=0;
 kq=1;
 kq2=1;

 for(k=1;k<=N;k++) {

 t=k;
 if (kq >= a) {
 do {
 t=t/a;
 v--;
 } while ((t % a) == 0);
 kq=0;
 }
 kq++;
 num=mul_mod(num,t,av);

 t=(2*k-1);
 if (kq2 >= a) {
 if (kq2 == a) {
 do {
 t=t/a;
 v++;
 } while ((t % a) == 0);
 }
 kq2-=a;
 }
 den=mul_mod(den,t,av);
 kq2+=2;

 if (v > 0) {
 t=inv_mod(den,av);
 t=mul_mod(t,num,av);
 t=mul_mod(t,k,av);
 for(i=v;i<vmax;i++) t=mul_mod(t,a,av);
 s+=t;
 if (s>=av) s-=av;
 }

 }

 t=pow_mod(10,n-1,av);
 s=mul_mod(s,t,av);
 sum=fmod(sum+(double) s/ (double) av,1.0);
 }
 printf("Decimal digits of pi at position %d: %09d\n",n,(int)(sum*1e9));
 gettimeofday(&t_fin, NULL);
 secs = timeval_diff(&t_fin, &t_ini);
 printf("%.16g milliseconds\n", secs * 1000.0);
 return 0;
}
Anuncis
Aquesta entrada ha esta publicada en Informática e Internet. Afegeix a les adreces d'interès l'enllaç permanent.

Deixa un comentari

Fill in your details below or click an icon to log in:

WordPress.com Logo

Esteu comentant fent servir el compte WordPress.com. Log Out / Canvia )

Twitter picture

Esteu comentant fent servir el compte Twitter. Log Out / Canvia )

Facebook photo

Esteu comentant fent servir el compte Facebook. Log Out / Canvia )

Google+ photo

Esteu comentant fent servir el compte Google+. Log Out / Canvia )

Connecting to %s