Skip to content

Commit

Permalink
externalize expansion
Browse files Browse the repository at this point in the history
to avoid local variables
  • Loading branch information
Expander committed Sep 17, 2023
1 parent 87177ba commit ce43981
Showing 1 changed file with 19 additions and 13 deletions.
32 changes: 19 additions & 13 deletions src/Li2.jl
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,15 @@ function reli2_approx(x::Float32)::Float32
end


# series expansion of Li2(z) for |z| <= 1 and Re(z) <= 0.5
# in terms of u = -log(1-z)
function li2_approx(u::ComplexF32)
B = (-1.0f0/4, 1.0f0/36, -1.0f0/3600, 1.0f0/211680)
u2 = u*u
u + u2*(B[1] + u*(B[2] + u2*(B[3] + u2*B[4])))
end


"""
reli2(x::Real)
Expand Down Expand Up @@ -148,6 +157,7 @@ li2(z::Real) = li2(Complex(z))

_li2(z::ComplexF16) = oftype(z, _li2(ComplexF32(z)))


function _li2(z::ComplexF32)::ComplexF32
clog(z) = log(abs(z)) + angle(z)*1.0f0im

Expand All @@ -157,7 +167,7 @@ function _li2(z::ComplexF32)::ComplexF32
if iz == 0.0f0
if rz <= 1.0f0
reli2(rz)
else # rz > 1.
else # Re(z) > 1
reli2(rz) - pi*clog(rz)*1.0f0im
end
else
Expand All @@ -166,24 +176,20 @@ function _li2(z::ComplexF32)::ComplexF32
if nz < eps(Float32)
z*(1.0f0 + 0.25f0*z)
else
(u, rest, sgn) = if rz <= 0.5f0
if rz <= 0.5f0
if nz > 1.0f0
(-clog(1.0f0 - inv(z)), -0.5f0*clog(-z)^2 - zeta2F32, -1.0f0)
else # nz <= 1.
(-clog(1.0f0 - z), 0.0f0 + 0.0f0im, 1.0f0)
-li2_approx(-clog(1.0f0 - inv(z))) - 0.5f0*clog(-z)^2 - zeta2F32
else # |z|^2 <= 1
li2_approx(-clog(1.0f0 - z))
end
else # rz > 0.5
else # Re(z) > 1/2
if nz <= 2.0f0*rz
l = -clog(z)
(l, l*clog(1.0f0 - z) + zeta2F32, -1.0f0)
else # nz > 2.0f0*rz
(-clog(1.0f0 - inv(z)), -0.5f0*clog(-z)^2 - zeta2F32, -1.0f0)
-li2_approx(l) + l*clog(1.0f0 - z) + zeta2F32
else # |z|^2 > 2*Re(z)
-li2_approx(-clog(1.0f0 - inv(z))) - 0.5f0*clog(-z)^2 - zeta2F32
end
end

B = (-1.0f0/4, 1.0f0/36, -1.0f0/3600, 1.0f0/211680)
u2 = u*u
rest + sgn*(u + u2*(B[1] + u*(B[2] + u2*(B[3] + u2*B[4]))))
end
end
end
Expand Down

0 comments on commit ce43981

Please sign in to comment.