I need to calculate the following:

S= 1- x^2 / 2! x^4 / 4! - x^6 / 6! ... (-1)^n * x^2n / (2n)!

Where n is between 1 and 100, and x is a double.

I have the following code:

unsigned int factorial (unsigned int n)
{
    if (n == 0)
        return 1;
    return n * factorial(n - 1);
}

double exFive(int n, double x)
{
    double s = 1;

    for (int i = 1; i <= n;   i)
    {
        int j = 2 * i; 
        s = s   pow(-1, i) * pow(x, 2*i) / factorial(j); //problem is here I guess
    }

    return s;
}

void fiveIO()
{
    int n = 1;
    double x;
   
    cout << "Input n: ";
    cin >> n;
    while ((n < 1) || (n > 100))
    {
        cout << "Wrong number, input again: ";
        cin >> n;
    } 

    cout << "Input: ";
    cin >> x;

    cout << "Result is " << fixed << setprecision(2) << exFive(n, x);
}

It works however the result is nan where n is above ~15.. but I don't know why. I would presume it’s the FiveX function.

So, for instance, n = 3, x = 5 outputs -7.16 (which is correct), but n = 50, x = 5 outputs "nan".. is it because the output is too large of a number? But then how am I supposed to do it?

CodePudding user response：

Avoid int overflow (undefined behavior (UB)) in factorial(j) Typical 32-bit int can only hold result of factorial(12).

Improve loop computation by calculating the term based on its prior value.

//for (int i = 1; i <= n;   i) {
//    int j = 2 * i; 
//    s = s   pow(-1, i) * pow(x, 2*i) / factorial(j); //problem is here I guess
//}

// Something like
double xx = x*x;
double term = 1.0;
for (int i = 1; i <= n;   i) {
    int j = 2 * i; 
    term *= -xx / ((j-1)*j);
    s  = term;
}

A more advanced approach would calculate the sum in reverse to minimize computation errors, yet a lot depends on the allowable range of x, something not constrained by OP.

CodePudding user response：

For the formula:

S= 1- x^2 / 2! x^4 / 4! - x^6 / 6! ... (-1)^n * x^2n / (2n)!

I submit this implementation scales slightly better since it avoids multiplying everything at once.

double exFive(int n, double x)
{
    double s = 1;

    bool subtract = true;

    for (int i = 1; i <= n; i  )
    {
        double inner = 1;
        int n2 = i * 2;
        for (int j = 1; j <= n2; j  )
        {
            inner *= (x / j);
        }

        if (subtract)
        {
            s = s - inner;
        }
        else
        {
            subtract = s   inner;
        }
        subtract = !subtract;
    }
    return s;
}

Sample Run:

Input n: 100
Input: 94
Result is -4417.00

The implementation above, when it needs to compute say (x^n)/(n!). Instead of computing

T  = (x*x*x*x....*x*x*x) / (1*2*3*4*5*6*....*(N-1)*N)`

In the above expression, both the numerator and denominator will overflow when attempting to compute either using pow of factorial functions. So instead it will compute it as:

T = (x/1) * (x/2) * (x/3) * (x/4) * (x/5) * (x/6)

This avoids using factorial and pow functions. And since each iteration of the outer loop is alternating subtraction or addition operations, the chance of overflow remains smaller.