i was trying to make a program that uses the chudnovsky algorithm to calculate pi with this code:

``````#include <cmath>
#include <iostream>
#include <cstdlib>

using namespace std;
float fac (float num)
{
float result=1;
for (float i=1; i<=num; ++i)
result=result*=i;
return result;
}

int main()
{
double pi=0;
float k=0;
float t[2];
cout << "Enter the number of terms:" <<endl;
cin >> k;

while (k<=10)
{
t[1]=fac(k);
t[2]=3*k+(3/2);
pi += fac(6*k)*(13591409+545140134*k)/fac(3*k)*pow(t[1],3)*pow((640320),t[2]);
k++;

}
cout << 1/(pi*12)<<endl;
system("pause");
return 0;
}``````

at first it outputted -1.#IND. Then after changing "cout << pi << endl;" to " cout << 1/(pi*12 << endl;" (which i had forgot about), got zero every time. any help?

6
Contributors
14
Replies
17
Views
8 Years
Discussion Span
Last Post by vmanes
Featured Replies
• 1
vmanes 1,165   8 Years Ago

It really helps if you implement the [URL="http://en.wikipedia.org/wiki/Chudnovsky_algorithm"]formula [/URL]correctly. [code] cin >> num; //the number of terms, or iterations to use for (k = 0; k < num; k++ ) //start k at 0, not the user's input { t[0]=fac(k); t[1]=3*k+3/2.0; //divide by a floating point, else 3/2 is just …

In C and C++ languages array subscripts begin with 0, not 1. So line 24 should be `t[0] = ...` . Other lines in your program need to be adjusted accordingly, such as use t[0] and t[1].

Edited by Ancient Dragon: n/a

you are getting a overflow and that is what the return value becomes. your call to pow when k=10 for instance returns 1.5556504132421497887188538979512e+174 but the max value for a float is 3.40282e+038. i would try using type double which has a max of 1.79769e+308.

Big, and I mean BIG, problem is first of all that you are using type float - which is not very precise in dealing with large numbers. And your numbers get really, really large.

What happens in this: `fac( 6 * k )` . At the limit of your loop, you're looking at 60! - wow! Any idea how many digits that is?

You should start by using doubles instead of floats - all the time. If your system and compiler actually provide a long double, use that (Visual C++ long doubles are only 8 bytes - no benefit.)

This really begs a math library supporting bigger numbers.

#IND is for indefinite.
You may also see #NAN and #INF from time to time as well.

All basically mean the same thing, and that is your calculations were way off (possibly the result of garbage data).

Here's a nice short intro to all the different things floats can represent.
http://www.randelshofer.ch/fhw/gri/float.html

Ok ,now the stuff with the infinite numbers is fixed, but now it gives values that aren't even close to pi! What have i done wrong now?

> What have i done wrong now?
You've changed your code, and we don't know what you've done.

Post the latest code.

``````#include <cmath>
#include <iostream>
#include <cstdlib>

using namespace std;
double fac (double num)
{
double result=1;
for (double i=1; i<=num; ++i)
result=result*=i;
return result;
}

int main()
{
long double pi=0;
long double k=0;
long double t[2];
cout << "Enter the number of terms:" <<endl;
cin >> k;

while (k<=10)
{
t[0]=fac(k);
t[1]=3*k+(3/2);
pi = fac(6*k)*(13591409+545140134*k)/fac(3*k)*pow(t[0],3)*pow((640320),t[1]) + pi;
k++;

}
cout << 1/(pi*12)<<endl;
system("pause");
return 0;
}``````

Edited by mebob: n/a

It really helps if you implement the formula correctly.

``````cin >> num;   //the number of terms, or iterations to use
for (k = 0; k < num; k++ ) //start k at 0, not the user's input
{
t[0]=fac(k);
t[1]=3*k+3/2.0;  //divide by a floating point, else 3/2 is just 1
//you forgot the -1^k term, which makes it alternating add/subtract
pi += pow( -1,k) * fac(6*k)*(13591409+545140134*k)/
(fac(3*k)*pow(t[0],3)*pow((640320),t[1])) ;
//and you forgot ( ) around the entire divisor
}``````

Be sure to set output precision to a large number of digits, else youi'll only see the default 6 digits C++ outputs.

for the most part, it works now, except when you enter 30 as the number of terms it outputs -0 and any more terms it outputs -1.#IND. how would i get it to display many digits?

Like I said earlier, you're pushing the limits of what a double can store. 170! is the upper limit, and then consider you're multiplying that by other large numbers.

So your term fac(6*k) is going to limit you to something less than 170/6 as the number of terms (iterations) you can successfully run, till you get on a system that truly supports long doubles.

Got a Sun system handy, they can do 16byte long doubles!

Just took a quick look at Fedora - in a VM running on top of 32 bit XP, g++ shows its long double as size 12 - I could run your code out to 50 iterations with no problem once all variables and function return updated to long double. Of course, it doesn't really get you any better pi due to the limit on the precision of the answer.

when i say i want long answers for pi, i mean like a program like pifast would output. is there any way to get it like that?

when i say i want long answers for pi, i mean like a program like pifast would output. is there any way to get it like that?

See jonsca's last post - you need special libraries/utilities that handle large precision numbers. Plain ole doubles/long doubles won't do it.