2

I am attempting to write a function that takes a value x and calculates the cos x by the series expansion, I'm always getting -inf, indifferent which value is read in.

#include <iostream>
#include <math.h>
#include <iomanip>

using namespace std;

int fac(int n){
        return n == 0 || n == 1 ? 1 : n * fac(n-1);
}

int main(){
        double eps = 1e-15;
        double x;
        cin >> x;
        long double ak = 1, sn = 0;
        for(int k=1; abs(ak) > eps * abs(sn); k++){
            double sgn = k % 2 == 0 ? 1 : -1;
            sn += ak;
            ak = sgn * pow(x, 2 * k) / fac(2*k);
        }
        cout << setprecision(4) << setw(5) << "x" << setprecision(15) << setw(20) << "cos x" << endl;
        cout << setprecision(4) << setw(5)<< x << " " << setprecision(15) << setw(20) << sn << endl;
        cout << setw(26) << cos(x) << endl;
        return 0;
}

I debugged the code and at a certain point ak gets -inf. But why?

ATW
  • 233
  • 2
  • 14
  • 1
    [What is a debugger and how can it help me diagnose problems?](https://stackoverflow.com/q/25385173/5910058) – Jesper Juhl May 23 '20 at 15:37
  • 7
    `fac(2*k)` would overflow an `int` real fast. `pow(x, 2 * k)` will lose precision real fast, and possibly also overflow in the end. You are computing two very large numbers and then dividing one by another - that ends up accumulating a huge error, one that dwarfs any useful result. You need a computationally stable formula for the series - I imagine one exists, but I don't remember it off the top. – Igor Tandetnik May 23 '20 at 15:38
  • 1
    With values smaller than 1 (say .5) it works better. – Marc Glisse May 23 '20 at 15:38
  • 1
    Print fac(2*k) at each step instead of pow, sorry. Or use -fsanitize=undefined. – Marc Glisse May 23 '20 at 15:40
  • @IgorTandetnik the overflow of fac is indeed a problem (making it return double helps a lot), but dividing huge numbers (as long as they don't overflow) doesn't necessarily give a bad precision. – Marc Glisse May 23 '20 at 15:43
  • 6
    At the very least, instead of computing each member from scratch, derive it from the previous member: `ak = -ak * x * x / (2*k - 1) / (2*k)` – Igor Tandetnik May 23 '20 at 15:43
  • 1
    You should include `SafeInt.hpp` it will help you detect integer overflow https://learn.microsoft.com/en-us/cpp/safeint/safeint-library?view=vs-2019 – metablaster May 23 '20 at 16:05

0 Answers0