2

Im having a problem with a function i wrote. The idea was to calculate sin and cosin values (operating on radians) using taylor expansion instead of js math objects. These are the equations:

sin(x) = (x^1)/1! - (x^3)/3! + (x^5)/5! - (x^7)/7! + (x^9)/9! - (x^11)/11! + ...

cos(x) = (x^0)/0! - (x^2)/2! + (x^4)/4! - (x^6)/6! + (x^8)/8! - (x^10)/10! + ...

I understand that when i type something like myCos(10,2) the result is going to be inaccurate because of low amount of iterations, however i dont understand why for (for example) x = 10 results start to be real at specifically iterNum = 6, and becomes NaN at iterNum = 80. The point is that for ranges like myCos/Sin(1-40, 5-50) (more or less) function works, but for higher numbers result becomes NaN. Not sure if my explanation is understandable, but i hope it is, just go ahead and play with the function in the console and youll see whats the problem
Here is my code:

function power(a,n) {
    var result = 1;
    for (var i = 0; i < n; i++) {
        result = result * a;
    }
    return result;
}
function factorial(z) {
    var result = 1;
    for (var i = 1; i <= z; i++) {
        result = result * i;
    }
    return result;
}
function mySin(x, iterNum) {
    var sin = 0;
    var n = 1;
    for (var i = 0; i <= iterNum; i++) {
        sin = sin + (power(x,n)/factorial(n) - power(x,n+2)/factorial(n+2));
        n = n + 4;
    }
    console.log(sin + " = my function.");
    console.log(Math.sin(x)  + " math.sin");
}
function myCos(x, iterNum) {
    var cos = 0;
    var n = 0;
    for (var i = 0; i <= iterNum; i++) {
        cos = cos + (power(x,n)/factorial(n) - power(x,n+2)/factorial(n+2));
        n = n + 4;
    }
    console.log(cos + " = my function.");
    console.log(Math.cos(x) + " math.cos");
}
ishmaelMakitla
  • 3,784
  • 3
  • 26
  • 32
mkrds
  • 35
  • 1
  • 7
  • See http://stackoverflow.com/a/22791396/3088138, http://stackoverflow.com/a/28227419/3088138 for factorial and power free implementations, which avoids the NaN, but not the catastrophic cancellation of very large terms for larger `x`. – Lutz Lehmann Sep 09 '16 at 07:56

2 Answers2

1

The problem you're running into is

> Infinity/Infinity
NaN

So when both power and factorial produce a value larger than JavaScript's floating point numbers can represent, you are adding NaN to your sum, and it propagates to the result.

You should be able to fix it by using

for (var i = 0; i <= iterNum; i++) {
    var member = power(x,n)/factorial(n) - power(x,n+2)/factorial(n+2);
    if (isNaN(member)) // maybe || member === 0
        break;
    sin += member;
    n += 4;
}
Bergi
  • 630,263
  • 148
  • 957
  • 1,375
  • Correct rationale and constructive solution recommendation. Additionally the reason for starting with whole numbers and then suddenly getting real numbers is due to the range of your integer data type. You are using JavaScript, which is like a chameleon about data types. At some point one of your function results of power or factorial becomes too big to fit into an integer, so it becomes a floating point number instead. From now on you work with floating point numbers. – Chris Tophski Sep 08 '16 at 10:46
  • @ChrisTophski: Wouldn't you need to cast the integers to floats for the division in any language? – Bergi Sep 08 '16 at 10:49
  • Not necessarily, or at least only implicitly. In strongly typed languages you would define a function accepting a float. Then you work with floats inside the function anyway and giving an integer as argument it would be converted implicitly to a float. One solution for the specific problem here would be to multiply at least one operand with `1.0`, just to make sure, a given integer is interpreted as float. If this is what you mean, I completely agree with you. – Chris Tophski Sep 08 '16 at 10:55
  • @ChrisTophski: You seemed to imply that `power` and `factorial` would return integers, which one would then need to cast (as we don't want integer division). If they return floats already, it works just the same as in JavaScript. – Bergi Sep 08 '16 at 10:57
  • Both the `power` and `factorial` functions here could just `return 1.0 * result;` to make them become a float. This would guarantee the same output data type for every input. The division would then be done on floats anyway and should work as you wrote in your answer. – Chris Tophski Sep 08 '16 at 11:04
  • @Bergi Thanks for reply, that solution does indeed solve the NaN problem, however for large numbers results are still drastically different mySin(50,50); 56931.192037791116 = my function. -0.26237485370392877 math.sin Now im not sure if im just using wrong tool (JS due to its limitations with numbers) or my code still has some flaws? – mkrds Sep 08 '16 at 11:05
  • 1
    @mkrds The taylor series just doesn't work that well for large values, it eventually diverges. You'd need a very large `iterNum`, larger than what you still can get useful power and factorial values for in JavaScript. That's why every real implementation uses `x = x % Math.PI` first. – Bergi Sep 08 '16 at 11:20
  • @ChrisTophski and Bergi: thanks guys for guiding me in a proper direction. It seems ill have to dig deeper into integers and floating point numbers in JS. – mkrds Sep 08 '16 at 11:35
  • @bergi better make that (2* Math.PI) – Jeremy Kahan Sep 08 '16 at 11:41
  • @JeremyKahan Oops, of course. I guess I wanted to write `x = x % Math.TAU` :-D – Bergi Sep 08 '16 at 11:43
  • @mkrds: mySin(50,50) fails as the largest term involved is around power 50, i.e., 50^49/49! or 50^51/51!. As these terms with values around e^50 cancel, they leave rather large cancellation errors of size exp(50)*10^(-16) ~ 5e5 that again may cancel or add to the total error of the evaluation. – Lutz Lehmann Sep 09 '16 at 08:04
  • @Bergi: The Taylor series for exp, sin, cos converge absolutely for all values of `x`. The evaluation error for floating point evaluation of sin and cos will grow proportional to `exp(x)`, as `x^n/n!~(ex/n)^n` with `n=floor(x)` or `n=ceil(x)` will be the largest term in the alternating series and thus `exp(x)*mu` will dominate the cancellation errors. – Lutz Lehmann Sep 09 '16 at 08:12
  • @LutzL I'm sorry, I must have used the wrong wording. Yes, the infinite series converges for all `x`. What I meant was that a finite expansion will have a large and growing error for `x` far away from `0`, [this animation](https://de.wikipedia.org/wiki/Datei:Taylor-polynomial-sinus-n-0-40.gif) displays it very well. – Bergi Sep 09 '16 at 10:21
0

The quotient between term k-1 and k of the sine series has the nice simple form

-x*x / ( (2*k)*(2*k+1) )

This allows for a very simple implementation using a minimum of floating point operations:

function mySin(x, iterNum) {
    var mxx = -x*x;
    var sin = 1;
    var n = 0;
    var term = 1;
    for (var i = 1; i <= 2*iterNum; i++) {
        n = n + 2;
        term = term * mxx / ( n*(n+1) );
        sin = sin + term
    }
    sin = x*sin;
    console.log(sin + " = my function.");
    console.log(Math.sin(x)  + " math.sin");
}

Since neither the power nor the factorial are actually computed, the largest value in this computation is the term where n-0.5 <= x <= n+1.5 with size about exp(x). This can still lead to values NaN, but only above x>705.

Lutz Lehmann
  • 25,219
  • 2
  • 22
  • 51