Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
15 changes: 10 additions & 5 deletions base/complex.jl
Original file line number Diff line number Diff line change
Expand Up @@ -472,9 +472,13 @@ function inv(z::Complex{T}) where T<:Union{Float16,Float32}
end
function inv(w::ComplexF64)
c, d = reim(w)
(isinf(c) | isinf(d)) && return complex(copysign(0.0, c), flipsign(-0.0, d))
absc, absd = abs(c), abs(d)
cd = ifelse(absc>absd, absc, absd) # cheap `max`: don't need sign- and nan-checks here
cd, dc = ifelse(absc>absd, (absc, absd), (absd, absc))
# no overflow from abs2
if sqrt(floatmin(Float64)/2) <= cd <= sqrt(floatmax(Float64)/2)
return conj(w) / muladd(cd, cd, dc*dc)
end
(isinf(c) | isinf(d)) && return complex(copysign(0.0, c), flipsign(-0.0, d))

ϵ = eps(Float64)
bs = 2/(ϵ*ϵ)
Expand All @@ -493,12 +497,13 @@ function inv(w::ComplexF64)
else
q, p = robust_cinv(-d, -c)
end
return ComplexF64(p*s, q*s) # undo scaling
return ComplexF64(p*s, q*s)
end
function robust_cinv(c::Float64, d::Float64)
r = d/c
p = inv(muladd(d, r, c))
q = -r*p
z = muladd(d, r, c)
p = 1.0/z
q = -r/z
return p, q
end

Expand Down
2 changes: 1 addition & 1 deletion doc/src/manual/complex-and-rational-numbers.md
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,7 @@ julia> 3(2 - 5im)^2
-63 - 60im

julia> 3(2 - 5im)^-1.0
0.20689655172413796 + 0.5172413793103449im
0.20689655172413793 + 0.5172413793103449im
```

The promotion mechanism ensures that combinations of operands of different types just work:
Expand Down