Chromium Code Reviews
chromiumcodereview-hr@appspot.gserviceaccount.com (chromiumcodereview-hr) | Please choose your nickname with Settings | Help | Chromium Project | Gerrit Changes | Sign out
(625)

Side by Side Diff: src/math.js

Issue 78873006: Embed trigonometric lookup table. (Closed) Base URL: https://v8.googlecode.com/svn/branches/bleeding_edge
Patch Set: tiny fixes Created 7 years, 1 month ago
Use n/p to move between diff chunks; N/P to move between comments. Draft comments are only viewable by you.
Jump to:
View unified diff | Download patch | Annotate | Revision Log
OLDNEW
1 // Copyright 2012 the V8 project authors. All rights reserved. 1 // Copyright 2012 the V8 project authors. All rights reserved.
2 // Redistribution and use in source and binary forms, with or without 2 // Redistribution and use in source and binary forms, with or without
3 // modification, are permitted provided that the following conditions are 3 // modification, are permitted provided that the following conditions are
4 // met: 4 // met:
5 // 5 //
6 // * Redistributions of source code must retain the above copyright 6 // * Redistributions of source code must retain the above copyright
7 // notice, this list of conditions and the following disclaimer. 7 // notice, this list of conditions and the following disclaimer.
8 // * Redistributions in binary form must reproduce the above 8 // * Redistributions in binary form must reproduce the above
9 // copyright notice, this list of conditions and the following 9 // copyright notice, this list of conditions and the following
10 // disclaimer in the documentation and/or other materials provided 10 // disclaimer in the documentation and/or other materials provided
(...skipping 61 matching lines...) Expand 10 before | Expand all | Expand 10 after
72 return %Math_atan2(TO_NUMBER_INLINE(y), TO_NUMBER_INLINE(x)); 72 return %Math_atan2(TO_NUMBER_INLINE(y), TO_NUMBER_INLINE(x));
73 } 73 }
74 74
75 // ECMA 262 - 15.8.2.6 75 // ECMA 262 - 15.8.2.6
76 function MathCeil(x) { 76 function MathCeil(x) {
77 return %Math_ceil(TO_NUMBER_INLINE(x)); 77 return %Math_ceil(TO_NUMBER_INLINE(x));
78 } 78 }
79 79
80 // ECMA 262 - 15.8.2.7 80 // ECMA 262 - 15.8.2.7
81 function MathCos(x) { 81 function MathCos(x) {
82 return MathCosImpl(x); 82 x = MathAbs(x); // Convert to number and get rid of -0.
83 return TrigonometricInterpolation(x, 1);
83 } 84 }
84 85
85 // ECMA 262 - 15.8.2.8 86 // ECMA 262 - 15.8.2.8
86 function MathExp(x) { 87 function MathExp(x) {
87 return %Math_exp(TO_NUMBER_INLINE(x)); 88 return %Math_exp(TO_NUMBER_INLINE(x));
88 } 89 }
89 90
90 // ECMA 262 - 15.8.2.9 91 // ECMA 262 - 15.8.2.9
91 function MathFloor(x) { 92 function MathFloor(x) {
92 x = TO_NUMBER_INLINE(x); 93 x = TO_NUMBER_INLINE(x);
(...skipping 86 matching lines...) Expand 10 before | Expand all | Expand 10 after
179 return (x < 0 ? (x + 0x100000000) : x) * 2.3283064365386962890625e-10; 180 return (x < 0 ? (x + 0x100000000) : x) * 2.3283064365386962890625e-10;
180 } 181 }
181 182
182 // ECMA 262 - 15.8.2.15 183 // ECMA 262 - 15.8.2.15
183 function MathRound(x) { 184 function MathRound(x) {
184 return %RoundNumber(TO_NUMBER_INLINE(x)); 185 return %RoundNumber(TO_NUMBER_INLINE(x));
185 } 186 }
186 187
187 // ECMA 262 - 15.8.2.16 188 // ECMA 262 - 15.8.2.16
188 function MathSin(x) { 189 function MathSin(x) {
189 return MathSinImpl(x); 190 x = x * 1; // Convert to number and deal with -0.
191 if (%_IsMinusZero(x)) return x;
192 return TrigonometricInterpolation(x, 0);
190 } 193 }
191 194
192 // ECMA 262 - 15.8.2.17 195 // ECMA 262 - 15.8.2.17
193 function MathSqrt(x) { 196 function MathSqrt(x) {
194 return %_MathSqrt(TO_NUMBER_INLINE(x)); 197 return %_MathSqrt(TO_NUMBER_INLINE(x));
195 } 198 }
196 199
197 // ECMA 262 - 15.8.2.18 200 // ECMA 262 - 15.8.2.18
198 function MathTan(x) { 201 function MathTan(x) {
199 return MathSinImpl(x) / MathCosImpl(x); 202 return MathSin(x) / MathCos(x);
200 } 203 }
201 204
202 // Non-standard extension. 205 // Non-standard extension.
203 function MathImul(x, y) { 206 function MathImul(x, y) {
204 return %NumberImul(TO_NUMBER_INLINE(x), TO_NUMBER_INLINE(y)); 207 return %NumberImul(TO_NUMBER_INLINE(x), TO_NUMBER_INLINE(y));
205 } 208 }
206 209
207 210
208 var MathSinImpl = function(x) { 211 var kInversePiHalf = 0.636619772367581343; // 2 / pi
209 InitTrigonometricFunctions(); 212 var kInversePiHalfS26 = 9.48637384723993156e-9; // 2 / pi / (2^26)
210 return MathSinImpl(x); 213 var kS26 = 1 << 26;
214 var kTwoStepThreshold = 1 << 27;
215 // pi / 2 rounded up
216 var kPiHalf = 1.570796326794896780; // 0x192d4454fb21f93f
217 // We use two parts for pi/2 to emulate a higher precision.
218 // pi_half_1 only has 26 significant bits for mantissa.
219 // Note that pi_half > pi_half_1 + pi_half_2
220 var kPiHalf1 = 1.570796325802803040; // 0x00000054fb21f93f
221 var kPiHalf2 = 9.920935796805404252e-10; // 0x3326a611460b113e
222
223 var kSamples; // Initialized to a number during genesis.
224 var kIndexConvert; // Initialized to a kSamples / (pi/2) during genesis.
225 var kSinTable; // Initialized to a Float64Array during genesis.
226 var kCosXIntervalTable; // Initialized to a Float64Array during genesis.
227
228 // This implements sine using the following algorithm.
229 // 1) Multiplication takes care of to-number conversion.
230 // 2) Reduce x to the first quadrant [0, pi/2].
231 // Conveniently enough, in case of +/-Infinity, we get NaN.
232 // Note that we try to use only 26 instead of 52 significant bits for
233 // mantissa to avoid rounding errors when multiplying. For very large
234 // input we therefore have additional steps.
235 // 3) Replace x by (pi/2-x) if x was in the 2nd or 4th quadrant.
236 // 4) Do a table lookup for the closest samples to the left and right of x.
237 // 5) Find the derivatives at those sampling points by table lookup:
238 // dsin(x)/dx = cos(x) = sin(pi/2-x) for x in [0, pi/2].
239 // 6) Use cubic spline interpolation to approximate sin(x).
240 // 7) Negate the result if x was in the 3rd or 4th quadrant.
241 // 8) Get rid of -0 by adding 0.
242 function TrigonometricInterpolation(x, phase) {
243 if (x < 0 || x > kPiHalf) {
244 var multiple;
245 while (x < -kTwoStepThreshold || x > kTwoStepThreshold) {
246 // Let's assume this loop does not terminate.
247 // All numbers x in each loop forms a set S.
248 // (1) abs(x) > 2^27 for all x in S.
249 // (2) abs(multiple) != 0 since (2^27 * inverse_pi_half_s26) > 1
250 // (3) multiple is rounded down in 2^26 steps, so the rounding error is
251 // at most max(ulp, 2^26).
252 // (4) so for x > 2^27, we subtract at most (1+pi/4)x and at least
253 // (1-pi/4)x
254 // (5) The subtraction results in x' so that abs(x') <= abs(x)*pi/4.
255 // Note that this difference cannot be simply rounded off.
256 // Set S cannot exist since (5) violates (1). Loop must terminate.
257 multiple = MathFloor(x * kInversePiHalfS26) * kS26;
258 x = x - multiple * kPiHalf1 - multiple * kPiHalf2;
259 }
260 multiple = MathFloor(x * kInversePiHalf);
261 x = x - multiple * kPiHalf1 - multiple * kPiHalf2;
262 phase += multiple;
263 }
264 var double_index = x * kIndexConvert;
265 if (phase & 1) double_index = kSamples - double_index;
266 var index = double_index | 0;
267 var t1 = double_index - index;
268 var t2 = 1 - t1;
269 var y1 = kSinTable[index];
270 var y2 = kSinTable[index + 1];
271 var dy = y2 - y1;
272 return (t2 * y1 + t1 * y2 +
273 t1 * t2 * ((kCosXIntervalTable[index] - dy) * t2 +
274 (dy - kCosXIntervalTable[index + 1]) * t1))
275 * (1 - (phase & 2)) + 0;
211 } 276 }
212 277
213
214 var MathCosImpl = function(x) {
215 InitTrigonometricFunctions();
216 return MathCosImpl(x);
217 }
218
219
220 var InitTrigonometricFunctions;
221
222
223 // Define constants and interpolation functions.
224 // Also define the initialization function that populates the lookup table
225 // and then wires up the function definitions.
226 function SetupTrigonometricFunctions() {
227 var samples = 1800; // Table size. Do not change arbitrarily.
228 var inverse_pi_half = 0.636619772367581343; // 2 / pi
229 var inverse_pi_half_s_26 = 9.48637384723993156e-9; // 2 / pi / (2^26)
230 var s_26 = 1 << 26;
231 var two_step_threshold = 1 << 27;
232 var index_convert = 1145.915590261646418; // samples / (pi / 2)
233 // pi / 2 rounded up
234 var pi_half = 1.570796326794896780; // 0x192d4454fb21f93f
235 // We use two parts for pi/2 to emulate a higher precision.
236 // pi_half_1 only has 26 significant bits for mantissa.
237 // Note that pi_half > pi_half_1 + pi_half_2
238 var pi_half_1 = 1.570796325802803040; // 0x00000054fb21f93f
239 var pi_half_2 = 9.920935796805404252e-10; // 0x3326a611460b113e
240 var table_sin;
241 var table_cos_interval;
242
243 // This implements sine using the following algorithm.
244 // 1) Multiplication takes care of to-number conversion.
245 // 2) Reduce x to the first quadrant [0, pi/2].
246 // Conveniently enough, in case of +/-Infinity, we get NaN.
247 // Note that we try to use only 26 instead of 52 significant bits for
248 // mantissa to avoid rounding errors when multiplying. For very large
249 // input we therefore have additional steps.
250 // 3) Replace x by (pi/2-x) if x was in the 2nd or 4th quadrant.
251 // 4) Do a table lookup for the closest samples to the left and right of x.
252 // 5) Find the derivatives at those sampling points by table lookup:
253 // dsin(x)/dx = cos(x) = sin(pi/2-x) for x in [0, pi/2].
254 // 6) Use cubic spline interpolation to approximate sin(x).
255 // 7) Negate the result if x was in the 3rd or 4th quadrant.
256 // 8) Get rid of -0 by adding 0.
257 var Interpolation = function(x, phase) {
258 if (x < 0 || x > pi_half) {
259 var multiple;
260 while (x < -two_step_threshold || x > two_step_threshold) {
261 // Let's assume this loop does not terminate.
262 // All numbers x in each loop forms a set S.
263 // (1) abs(x) > 2^27 for all x in S.
264 // (2) abs(multiple) != 0 since (2^27 * inverse_pi_half_s26) > 1
265 // (3) multiple is rounded down in 2^26 steps, so the rounding error is
266 // at most max(ulp, 2^26).
267 // (4) so for x > 2^27, we subtract at most (1+pi/4)x and at least
268 // (1-pi/4)x
269 // (5) The subtraction results in x' so that abs(x') <= abs(x)*pi/4.
270 // Note that this difference cannot be simply rounded off.
271 // Set S cannot exist since (5) violates (1). Loop must terminate.
272 multiple = MathFloor(x * inverse_pi_half_s_26) * s_26;
273 x = x - multiple * pi_half_1 - multiple * pi_half_2;
274 }
275 multiple = MathFloor(x * inverse_pi_half);
276 x = x - multiple * pi_half_1 - multiple * pi_half_2;
277 phase += multiple;
278 }
279 var double_index = x * index_convert;
280 if (phase & 1) double_index = samples - double_index;
281 var index = double_index | 0;
282 var t1 = double_index - index;
283 var t2 = 1 - t1;
284 var y1 = table_sin[index];
285 var y2 = table_sin[index + 1];
286 var dy = y2 - y1;
287 return (t2 * y1 + t1 * y2 +
288 t1 * t2 * ((table_cos_interval[index] - dy) * t2 +
289 (dy - table_cos_interval[index + 1]) * t1))
290 * (1 - (phase & 2)) + 0;
291 }
292
293 var MathSinInterpolation = function(x) {
294 x = x * 1; // Convert to number and deal with -0.
295 if (%_IsMinusZero(x)) return x;
296 return Interpolation(x, 0);
297 }
298
299 // Cosine is sine with a phase offset.
300 var MathCosInterpolation = function(x) {
301 x = MathAbs(x); // Convert to number and get rid of -0.
302 return Interpolation(x, 1);
303 };
304
305 %SetInlineBuiltinFlag(Interpolation);
306 %SetInlineBuiltinFlag(MathSinInterpolation);
307 %SetInlineBuiltinFlag(MathCosInterpolation);
308
309 InitTrigonometricFunctions = function() {
310 table_sin = new global.Float64Array(samples + 2);
311 table_cos_interval = new global.Float64Array(samples + 2);
312 %PopulateTrigonometricTable(table_sin, table_cos_interval, samples);
313 MathSinImpl = MathSinInterpolation;
314 MathCosImpl = MathCosInterpolation;
315 }
316 }
317
318 SetupTrigonometricFunctions();
319
320
321 // ------------------------------------------------------------------- 278 // -------------------------------------------------------------------
322 279
323 function SetUpMath() { 280 function SetUpMath() {
324 %CheckIsBootstrapping(); 281 %CheckIsBootstrapping();
325 282
326 %SetPrototype($Math, $Object.prototype); 283 %SetPrototype($Math, $Object.prototype);
327 %SetProperty(global, "Math", $Math, DONT_ENUM); 284 %SetProperty(global, "Math", $Math, DONT_ENUM);
328 %FunctionSetInstanceClassName(MathConstructor, 'Math'); 285 %FunctionSetInstanceClassName(MathConstructor, 'Math');
329 286
330 // Set up math constants. 287 // Set up math constants.
(...skipping 56 matching lines...) Expand 10 before | Expand all | Expand 10 after
387 "atan2", MathAtan2, 344 "atan2", MathAtan2,
388 "pow", MathPow, 345 "pow", MathPow,
389 "max", MathMax, 346 "max", MathMax,
390 "min", MathMin, 347 "min", MathMin,
391 "imul", MathImul 348 "imul", MathImul
392 )); 349 ));
393 350
394 %SetInlineBuiltinFlag(MathSin); 351 %SetInlineBuiltinFlag(MathSin);
395 %SetInlineBuiltinFlag(MathCos); 352 %SetInlineBuiltinFlag(MathCos);
396 %SetInlineBuiltinFlag(MathTan); 353 %SetInlineBuiltinFlag(MathTan);
354 %SetInlineBuiltinFlag(TrigonometricInterpolation);
397 } 355 }
398 356
399 SetUpMath(); 357 SetUpMath();
OLDNEW
« no previous file with comments | « src/bootstrapper.cc ('k') | src/runtime.h » ('j') | src/trig-table.h » ('J')

Powered by Google App Engine
This is Rietveld 408576698