Skip to content

Instantly share code, notes, and snippets.

/CheckKahanAngle.jl Secret
Created Aug 23, 2017

Embed
What would you like to do?
module CheckKahanAngle
using PyPlot
function kahanAngle(a, b, c)
b >= c >= 0 || c > b >= 0 || error("not a triangle")
μ = ifelse(b >= c >= 0, c - (a - b), b - (a - c))
2*atan(sqrt(((a-b)+c)*μ/((a+(b+c))*((a-c)+b))))
end
function f(v1, v2)
a, b, c = norm(v1), norm(v2), norm(v1 - v2)
kahanAngle(a, b, c)
end
function f2(v1, v2)
v1, v2 = big.(v1), big.(v2)
acos(dot(v1./norm(v1), v2./norm(v2)))
end
function fe{T<:AbstractFloat}::T)
v1, v2 = [1.0, 0.0], T.([cos(big(θ)), sin(big(θ))])
T(abs(big(f(v1, v2)) / f2(big.(v1), big.(v2)) - 1))
end
function main()
x = 10.^linspace(-16, -4, 1024)
figure()
loglog(x, fe.(pi/2-x))
loglog(x, fe.(x))
end
end
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.