1 // cl_cossin_ratseries().
12 #include "cl_lfloat.h"
14 #include "cl_integer.h"
16 inline const cl_LF_cos_sin_t operator* (const cl_LF_cos_sin_t& a, const cl_LF_cos_sin_t& b)
18 return cl_LF_cos_sin_t(a.cos*b.cos-a.sin*b.sin,a.sin*b.cos+a.cos*b.sin);
21 const cl_LF_cos_sin_t cl_cossin_ratseries (const cl_LF& x)
23 // Similar to expx_ratseries.
24 var uintC len = TheLfloat(x)->len;
25 var cl_idecoded_float x_ = integer_decode_float(x);
26 // x = (-1)^sign * 2^exponent * mantissa
27 var uintL lq = cl_I_to_UL(- x_.exponent);
28 var const cl_I& p = x_.mantissa;
29 // Compute sin(p/2^lq) and cos(p/2^lq) by splitting into pieces.
30 var cl_boolean first_factor = cl_true;
31 var cl_LF_cos_sin_t product;
34 for (b1 = 0, b2 = 1; b1 < lq; b1 = b2, b2 = 2*b2) {
35 // Piece containing bits b1+1..b2 after "decimal point"
36 // in the binary representation of (p/2^lq).
37 var uintL lqk = (lq >= b2 ? b2 : lq);
38 var cl_I pk = ldb(p,cl_byte(lqk-b1,lq-lqk));
39 // Compute sin(pk/2^lqk) and cos(pk/2^lqk).
41 if (minusp(x_.sign)) { pk = -pk; }
42 var cl_LF_cos_sin_t factor = cl_cossin_aux(pk,lqk,len);
45 first_factor = cl_false;
47 product = product * factor;
51 return cl_LF_cos_sin_t(cl_I_to_LF(1,len),cl_I_to_LF(0,len));
55 // Bit complexity (N = length(x)): O(log(N)^2*M(N)).