Skip to content

Commit

Permalink
Minor refactoring mistake
Browse files Browse the repository at this point in the history
  • Loading branch information
auroranockert committed Mar 17, 2012
1 parent af9ae86 commit ec630d8
Show file tree
Hide file tree
Showing 2 changed files with 33 additions and 20 deletions.
51 changes: 32 additions & 19 deletions lib/csmath.js
Expand Up @@ -237,7 +237,7 @@ function(global) {
CSMath.min = Math.min
CSMath.max = Math.max

CSMath.hypot = norm2
CSMath.hypot = CSMath.norm2
/* TODO: The trigonometric functions have low precision */

CSMath.round = Math.round
Expand All @@ -249,34 +249,47 @@ function(global) {
if (CSMath.accurate['sqrt']) {
CSMath.sqrt = Math.sqrt
} else {
/* Algorithm from 'Divide, Square Root and Remainder Algorithms for the IA64 Architecture */
/*
* Algorithm from 'Divide, Square Root and Remainder Algorithms for the IA64
* Architecture' by Intel
*
* - It is proven correct under a few assumptions
* 1. That Math.sqrt is at least as accurate as the IA64 approximation
* instruction.
* 2. That we're using the IA64 extended floating point format. (Which
* we are not)
* - I haven't checked if this algorithm actually works in double-precision,
* but my intuition says that it should be fine.
*/

Math.sqrt = function(a) {
var y0 = 1.0 / Math.sqrta(a)
var sqrtApproximation = Math.sqrt

var H0 = 0.5 * y0,
S0 = (a * y0)
CSMath.sqrt = function(a) {
var y0 = 1.0 / sqrtApproximation(a)

var d0 = -Math.fms(S0, H0, 0.5)
var H0 = 0.5 * y0,
S0 = (a * y0)

var H1 = Math.fma(d0, H0, H0),
S1 = Math.fma(d0, S0, S0)
var d0 = -CSMath.fms(S0, H0, 0.5)

var d1 = -Math.fms(S1, H1, 0.5)
var H1 = CSMath.fma(d0, H0, H0),
S1 = CSMath.fma(d0, S0, S0)

var H2 = Math.fma(d1, H1, H1),
S2 = Math.fma(d1, S1, S1)
var d1 = -CSMath.fms(S1, H1, 0.5)

var d2 = -Math.fms(S2, H2, 0.5),
e2 = -Math.fms(S2, S2, a)
var H2 = CSMath.fma(d1, H1, H1),
S2 = CSMath.fma(d1, S1, S1)

var H3 = Math.fma(d2, H2, H2),
S3 = Math.fma(e2, H2, S2)
var d2 = -CSMath.fms(S2, H2, 0.5),
e2 = -CSMath.fms(S2, S2, a)

var e3 = -Math.fms(S3, S3, a)
var H3 = CSMath.fma(d2, H2, H2),
S3 = CSMath.fma(e2, H2, S2)

return Math.fma(e3, H3, S3)
}
var e3 = -CSMath.fms(S3, S3, a)

return CSMath.fma(e3, H3, S3)
}
}

/* TODO: The root functions have low precision */
Expand Down
2 changes: 1 addition & 1 deletion src/norm.erb.js
Expand Up @@ -51,4 +51,4 @@ CSMath.hypot2 = function() {
CSMath.min = Math.min
CSMath.max = Math.max

CSMath.hypot = norm2
CSMath.hypot = CSMath.norm2

0 comments on commit ec630d8

Please sign in to comment.