Long double overflows but value smaller than maximum representable value

Question

I'm trying to compute a series using C++. The series is:
(for those wondering)

My code is the following:

#include <iostream>
#include <fstream>
#include <cmath> // exp 
#include <iomanip> //setprecision, setw 
#include <limits> //numeric_limits (http://en.cppreference.com/w/cpp/types/numeric_limits)

long double SminOneCenter(long double gamma)
{
    using std::endl; using std::cout;
    long double result=0.0l;
    for (long double k = 1; k < 1000 ; k++)
    {   
            if(isinf(pow(1.0l+pow(gamma,k),6.0l/4.0l)))
            {   
                    cout << "infinity for reached for gamma equals:   " << gamma <<  "value of k:  " << k ; 
                    cout << "maximum allowed:   " <<  std::numeric_limits<long double>::max()<< endl;
                    break;
            }   

                    // CAS PAIR: -1^n = 1
                    if ((int)k%2 == 0)
                    {   
                            result += pow(4.0l*pow(gamma,k),3.0l/4.0l) /(pow(1+pow(gamma,k)),6.0l/4.0l);
                    }   
                    // CAS IMPAIR:-1^n = -1
                    else if ((int)k%2!=0)
                    {   
                            result -= pow(4.0l*pow(gamma,k),3.0l/4.0l) /(pow(1+pow(gamma,k)),6.0l/4.0l);

                            //if (!isinf(pow(k,2.0l)*zeta/2.0l))
                    }   
                    //              cout << result << endl;
    }    


    return 1.0l + 2.0l*result;
}

Output will be, for instance with gamma = 1.7 : infinity reached for gamma equals: 1.7 value of k: 892

The maximum value a long double can represent, as provided by the STL numeric_limits, is: 1.18973e+4932.

However (1+1.7^892)= 2.19.... × 10^308 which is way lower than 10^4932, so it shouldn't be considered as infinity.

Provided my code is not wrong (but it very well might be), could anyone tell me why the discussed code evals to infinity when it should not?

To start with: You do know that floating point arithmetic on computers will leads to compounded rounding errors? (See e.g. [Is floating point math broken?](http://stackoverflow.com/questions/588004/is-floating-point-math-broken)) Secondly, have you tried stepping through your code, line by line, in a debugger? — Some programmer dude, Mar 20 '17 at 15:56
Yeah, i knew floating point math isn't safe from errors, but that wasn't the issue, see below. — Devolution, Mar 20 '17 at 16:33

Bathsheba · Accepted Answer · 2017-03-20T15:59:05.593

4

You need to use powl rather than pow if you want to supply long double arguments.

Currently you are hitting the numeric_limits<double>::max() in your pow calls.

As an alternative, consider using std::pow which has appropriate overloads.

Reference http://en.cppreference.com/w/c/numeric/math/pow

edited Mar 20 '17 at 15:59

answered Mar 20 '17 at 15:57

Bathsheba

231,907
34
361
483

Isn't pow supposed to be overloaded for long double? http://en.cppreference.com/w/cpp/numeric/math/pow – Devolution Mar 20 '17 at 15:58
No. See the link. – Bathsheba Mar 20 '17 at 15:58
2

Your link is for ``. The OP is using `` which is overloaded: http://en.cppreference.com/w/cpp/numeric/math/pow – NathanOliver Mar 20 '17 at 15:58
Well, it can't be much else. It's either that sizeof(long double) is sizeof(double) or the compiler is picking up `pow`. But let's keep an eye on it: please retain the downvote if you think I'm incorrect as it helps the peer review process. I will delete the answer if it's shown to be incorrect. – Bathsheba Mar 20 '17 at 15:59
2

@NathanOliver But what magic is allowing that `pow` to be called unqualified without `using namespace std;`? (Unless that was omitted in OP's code, that is.) Something here doesn't quite make sense... – cdhowie Mar 20 '17 at 16:00
@cdhowie Dang and blast. I thought I saw a `using namespace std;`. They must also be using the `isinf` from `` as well – NathanOliver Mar 20 '17 at 16:01
I have a special key to type `std::` on my DasKeyboard. I *never* use `using namespace std;`. – Bathsheba Mar 20 '17 at 16:02
@cdhowie `` must put its names into `::std`, but it *may* also put them into `::`. `` must put its names into `::`, but it *may* also put them into `::std`. – Angew is no longer proud of SO Mar 20 '17 at 16:02
Damnit. So isinf also exist outside namespace std? I'll try to add std:: everywhere i use math function and see if it was it. Strange thing is i have almost the same function for another series and it's working perfectly without using std::pow in lieu of pow. (And this second series uses fractional exponent as well) – Devolution Mar 20 '17 at 16:04
1

Okay, seems it was it. Calling std::pow and std::isinf doesn't trigger the if unless the long double is actually overflowing. Thank you very much for your help – Devolution Mar 20 '17 at 16:10
1

As a final point, always add the smaller floating point terms first in a series. – Bathsheba Mar 20 '17 at 16:12
Shouldn't you have got a compiler warning about loss of precision when you passed a long double to a ::pow that took a double? Unless, of course, your warning level is not set high enough. – Klitos Kyriacou Mar 20 '17 at 16:44

Long double overflows but value smaller than maximum representable value

1 Answers1

Linked