www.digitalmars.com         C & C++   DMDScript  

digitalmars.D.bugs - Corrected pi.d

reply Berin Loritsch <bloritsch d-haven.org> writes:
Content-Type: text/plain; charset=us-ascii; format=flowed
Content-Transfer-Encoding: 7bit

There were a few deficiencies in the pi.d sample, including the display 
error pointed out by Sampsa.

I took the opportunity to modify the sample to show off dynamic array
size language feature of D.  I also noticed that there is a calculation
error on the last digit, so I modified the source to calculate one more
digit than it would display.  This seems to solve the most egregious and
obvious errors.  (i.e. PI does not equal 3.12, it should be 3.14).

I added a little checking so that we don't try to compute a negative
number of digits--and print out a nice message.  Also, instead of
forcing someone to never compute more than 4000 digits, I had it able
to compute any number of digits the user entered.  Although, it does
print out a warning "Be prepared to wait a while..." if the number is
larger than 4000.

While I can't vouch for the perfectness of the math (I did nothing to
change it), it no longer outputs obviously wrong info when low precision
is used.

Please include it in the next version release samples.
Jul 23 2004
parent Berin Loritsch <bloritsch d-haven.org> writes:
Is there no interest in correcting the sample?  It looks kind of bad
if it is broken, and the file is included here...

Berin Loritsch wrote:

 There were a few deficiencies in the pi.d sample, including the display 
 error pointed out by Sampsa.
 
 I took the opportunity to modify the sample to show off dynamic array
 size language feature of D.  I also noticed that there is a calculation
 error on the last digit, so I modified the source to calculate one more
 digit than it would display.  This seems to solve the most egregious and
 obvious errors.  (i.e. PI does not equal 3.12, it should be 3.14).
 
 I added a little checking so that we don't try to compute a negative
 number of digits--and print out a nice message.  Also, instead of
 forcing someone to never compute more than 4000 digits, I had it able
 to compute any number of digits the user entered.  Although, it does
 print out a warning "Be prepared to wait a while..." if the number is
 larger than 4000.
 
 While I can't vouch for the perfectness of the math (I did nothing to
 change it), it no longer outputs obviously wrong info when low precision
 is used.
 
 Please include it in the next version release samples.
 
 
 ------------------------------------------------------------------------
 
 import std.c.stdio;
 import std.c.stdlib;
 import std.c.time;
 
 const int LONG_TIME=4000;
 
 byte[] p;
 byte[] t;
 int q;
 
 int main(char[][] args)
 {
 	int startime, endtime;
 	int i;
 
 	if (args.length == 2) {
 		sscanf(&args[1][0],"%d",&q);
 	} else {
 		printf("Usage: pi [precision]\n");
 		exit(55);
 	}
 
         if (q < 0)
         {
                 printf("Precision was too low, running with precision of
0.\n");
                 q = 0;
         }
 
         if (q > LONG_TIME)
         {
                 printf("Be prepared to wait a while...\n");
         }
 
         // Compute one more digit than we display to compensate for rounding
         q++;
 
         p.length = q + 1;
         t.length = q + 1;
 
 	/* compute pi */
 
 	std.c.time.time(&startime);
 	arctan(2);
 	arctan(3);
 	mul4();
 	std.c.time.time(&endtime);
 
         // Return to the number of digits we want to display
         q--;
 
 	/* print pi */
 
 	printf("pi = %d.",cast(int)(p[0]));
 	for (i = 1; i <= q; i++)
 		printf("%d",cast(int)(p[i]));
 	printf("\n");
 	printf("%ld seconds to compute pi with a precision of %d
digits.\n",endtime-startime,q);
 
 	return 0;
 }
 
 void arctan(int s)
 {
 	int n;
 
 	t[0] = 1;
 	div(s);			/* t[] = 1/s */
 	add();
 	n = 1;
 	do {
 		mul(n);
 		div(s * s);
 		div(n += 2);
 		if (((n-1) / 2) % 2 == 0)
 			add();
 		else
 			sub();
 	} while (!tiszero());
 }
 
 void add()
 {
 	int j;
 
 	for (j = q; j >= 0; j--)
 	{
 		if (t[j] + p[j] > 9) {
 			p[j] += t[j] - 10;
 			p[j-1] += 1;
 		} else
 			p[j] += t[j];
 	}
 }
 
 void sub()
 {
 	int j;
 
 	for (j = q; j >= 0; j--)
 		if (p[j] < t[j]) {
 			p[j] -= t[j] - 10;
 			p[j-1] -= 1;
 		} else
 			p[j] -= t[j];
 }
 
 void mul(int multiplier)
 {
 	int b;
 	int i;
 	int carry = 0, digit = 0;
 
 	for (i = q; i >= 0; i--) {
 		b = (t[i] * multiplier + carry);
 		digit = b % 10;
 		carry = b / 10;
 		t[i] = digit;
 	}
 }
 
 /* t[] /= l */
 
 void div(int divisor)
 {
 	int i, b;
 	int quotient, remainder = 0;
 
 	for (i = 0; i <= q; i++) {
 		b = (10 * remainder + t[i]);
 		quotient = b / divisor;
 		remainder = b % divisor;
 		t[i] = quotient;
 	}
 }
 
 void div4()
 {
 	int i, c, d = 0;
 
 	for (i = 0; i <= q; i++) {
 		c = (10 * d + p[i]) / 4;
 		d = (10 * d + p[i]) % 4;
 		p[i] = c;
 	}
 }
 
 void mul4()
 {
 	int i, c, d;
 
 	d = c = 0;
 
 	for (i = q; i >= 0; i--) {
 		d = (p[i] * 4 + c) % 10;
 		c = (p[i] * 4 + c) / 10;
 		p[i] = d;
 	}
 }
 
 int tiszero()
 {
 	int k;
 
 	for (k = 0; k <= q; k++)
 		if (t[k] != 0)
 			return false;
 	return true;
 }
 
 

Jul 28 2004