|
|
1.1 ! root 1: /* ! 2: C program for floating point sin/cos. ! 3: Calls modf. ! 4: There are no error exits. ! 5: Coefficients are #3370 from Hart & Cheney (18.80D). ! 6: */ ! 7: ! 8: static double twoopi = 0.63661977236758134308; ! 9: static double p0 = .1357884097877375669092680e8; ! 10: static double p1 = -.4942908100902844161158627e7; ! 11: static double p2 = .4401030535375266501944918e6; ! 12: static double p3 = -.1384727249982452873054457e5; ! 13: static double p4 = .1459688406665768722226959e3; ! 14: static double q0 = .8644558652922534429915149e7; ! 15: static double q1 = .4081792252343299749395779e6; ! 16: static double q2 = .9463096101538208180571257e4; ! 17: static double q3 = .1326534908786136358911494e3; ! 18: ! 19: double ! 20: cos(arg) ! 21: double arg; ! 22: { ! 23: double sinus(); ! 24: if(arg<0) ! 25: arg = -arg; ! 26: return(sinus(arg, 1)); ! 27: } ! 28: ! 29: double ! 30: sin(arg) ! 31: double arg; ! 32: { ! 33: double sinus(); ! 34: return(sinus(arg, 0)); ! 35: } ! 36: ! 37: static double ! 38: sinus(arg, quad) ! 39: double arg; ! 40: int quad; ! 41: { ! 42: double modf(); ! 43: double e, f; ! 44: double ysq; ! 45: double x,y; ! 46: int k; ! 47: double temp1, temp2; ! 48: ! 49: x = arg; ! 50: if(x<0) { ! 51: x = -x; ! 52: quad = quad + 2; ! 53: } ! 54: x = x*twoopi; /*underflow?*/ ! 55: if(x>32764){ ! 56: y = modf(x,&e); ! 57: e = e + quad; ! 58: modf(0.25*e,&f); ! 59: quad = e - 4*f; ! 60: }else{ ! 61: k = x; ! 62: y = x - k; ! 63: quad = (quad + k) & 03; ! 64: } ! 65: if (quad & 01) ! 66: y = 1-y; ! 67: if(quad > 1) ! 68: y = -y; ! 69: ! 70: ysq = y*y; ! 71: temp1 = ((((p4*ysq+p3)*ysq+p2)*ysq+p1)*ysq+p0)*y; ! 72: temp2 = ((((ysq+q3)*ysq+q2)*ysq+q1)*ysq+q0); ! 73: return(temp1/temp2); ! 74: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.