From 3228cf945a2d5a782d156034623860c86cce21fd Mon Sep 17 00:00:00 2001
From: Yichao Yu <yyc1992@gmail.com>
Date: Sat, 6 May 2017 12:46:38 -0400
Subject: [PATCH] Use sincos in cis and other complex math functions

---
 base/complex.jl  | 35 ++++++++++++++++++++++++-----------
 base/fastmath.jl |  5 ++++-
 2 files changed, 28 insertions(+), 12 deletions(-)

diff --git a/base/complex.jl b/base/complex.jl
index ae9a191b181d4..b2f50aa0955b9 100644
--- a/base/complex.jl
+++ b/base/complex.jl
@@ -421,7 +421,10 @@ sqrt(z::Complex) = sqrt(float(z))
 # end
 
 # compute exp(im*theta)
-cis(theta::Real) = Complex(cos(theta),sin(theta))
+function cis(theta::Real)
+    s, c = sincos(theta)
+    Complex(c, s)
+end
 
 """
     cis(z)
@@ -430,7 +433,8 @@ Return ``\\exp(iz)``.
 """
 function cis(z::Complex)
     v = exp(-imag(z))
-    Complex(v*cos(real(z)), v*sin(real(z)))
+    s, c = sincos(real(z))
+    Complex(v * c, v * s)
 end
 
 """
@@ -507,7 +511,8 @@ function exp(z::Complex)
         if iszero(zi)
             Complex(er, zi)
         else
-            Complex(er*cos(zi), er*sin(zi))
+            s, c = sincos(zi)
+            Complex(er * c, er * s)
         end
     end
 end
@@ -530,8 +535,13 @@ function expm1(z::Complex)
             Complex(erm1, zi)
         else
             er = erm1+one(erm1)
-            wr = isfinite(er) ? erm1 - 2.0*er*(sin(0.5*zi))^2 : er*cos(zi)
-            Complex(wr, er*sin(zi))
+            if isfinite(er)
+                wr = erm1 - 2.0*er*(sin(0.5*zi))^2
+                return Complex(wr, er * sin(zi))
+            else
+                s, c = sincos(zi)
+                return Complex(er * c, er * s)
+            end
         end
     end
 end
@@ -592,13 +602,15 @@ end
 function exp2(z::Complex{T}) where T
     er = exp2(real(z))
     theta = imag(z) * log(convert(T, 2))
-    Complex(er*cos(theta), er*sin(theta))
+    s, c = sincos(theta)
+    Complex(er * c, er * s)
 end
 
 function exp10(z::Complex{T}) where T
     er = exp10(real(z))
     theta = imag(z) * log(convert(T, 10))
-    Complex(er*cos(theta), er*sin(theta))
+    s, c = sincos(theta)
+    Complex(er * c, er * s)
 end
 
 function ^(z::T, p::T) where T<:Complex
@@ -620,8 +632,7 @@ function ^(z::T, p::T) where T<:Complex
         rp = rp*exp(-pim*theta)
         ntheta = ntheta + pim*log(r)
     end
-    cosntheta = cos(ntheta)
-    sinntheta = sin(ntheta)
+    sinntheta, cosntheta = sincos(ntheta)
     re, im = rp*cosntheta, rp*sinntheta
     if isinf(rp)
         if isnan(re)
@@ -681,7 +692,8 @@ function sin(z::Complex{T}) where T
             Complex(F(NaN), F(NaN))
         end
     else
-        Complex(sin(zr)*cosh(zi), cos(zr)*sinh(zi))
+        s, c = sincos(zr)
+        Complex(s * cosh(zi), c * sinh(zi))
     end
 end
 
@@ -700,7 +712,8 @@ function cos(z::Complex{T}) where T
             Complex(F(NaN), F(NaN))
         end
     else
-        Complex(cos(zr)*cosh(zi), -sin(zr)*sinh(zi))
+        s, c = sincos(zr)
+        Complex(c * cosh(zi), -s * sinh(zi))
     end
 end
 
diff --git a/base/fastmath.jl b/base/fastmath.jl
index d05a80131d971..6eca987d0703b 100644
--- a/base/fastmath.jl
+++ b/base/fastmath.jl
@@ -327,7 +327,10 @@ sincos_fast(v) = (sin_fast(v), cos_fast(v))
 
     # complex numbers
 
-    cis_fast(x::T) where {T<:FloatTypes} = Complex{T}(cos(x), sin(x))
+    function cis_fast(x::T) where {T<:FloatTypes}
+        s, c = sincos_fast(x)
+        Complex{T}(c, s)
+    end
 
     # See <http://en.cppreference.com/w/cpp/numeric/complex>
     pow_fast(x::T, y::T) where {T<:ComplexTypes} = exp(y*log(x))