CODEWITHSUNDEEP
c++intdoubleprecision
oobrbzn
oobrbzn

Reputation: 103

Gaussian integral in C++ is not working. Why?

Here is my code:

const double kEps(0.00000001);
    double gaussianIntegral(double x)
    {
        double res(x), numerator(x);
        for(unsigned int i(1), k(3); (abs(numerator/k) > kEps) || (res < 0); ++i, k+=2)
        {
            numerator*=(-x*x/i);
            res+=numerator/k;
        }
        return res;
    }

Here is what I am trying to compute:

enter image description here

When I try to pass 30 as an argument my computations go forever. What is wrong? I am very stuck, it seems for me like there is no error and everything should work just fine, but still it is not so.

Upvotes: 0

Views: 1018

Answers (1)

lomereiter
lomereiter

Reputation: 356

Although formally Taylor series converges, in practice you will run into machine precision limits even for the argument as small as 0.1 (and that is using kEps=0)

It's better to use (scaled appropriately) std::erf (C++11), or if it's a homework, look up an algorithm for computing erf function, e.g. here: https://math.stackexchange.com/questions/97/how-to-accurately-calculate-the-error-function-erfx-with-a-computer

Upvotes: 2

Related Questions