Adding upstream version 1.16.10.upstream/1.16.10upstream

Signed-off-by: Daniel Baumann <daniel.baumann@progress-linux.org>

17 files changed, 2769 insertions, 0 deletions

diff --git a/src/math/cmplx/abs.go b/src/math/cmplx/abs.go
new file mode 100644
index 0000000..2f89d1b
--- /dev/null
+++ b/src/math/cmplx/abs.go
@@ -0,0 +1,13 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+// Package cmplx provides basic constants and mathematical functions for
+// complex numbers. Special case handling conforms to the C99 standard
+// Annex G IEC 60559-compatible complex arithmetic.
+package cmplx
+
+import "math"
+
+// Abs returns the absolute value (also called the modulus) of x.
+func Abs(x complex128) float64 { return math.Hypot(real(x), imag(x)) }
diff --git a/src/math/cmplx/asin.go b/src/math/cmplx/asin.go
new file mode 100644
index 0000000..30d019e
--- /dev/null
+++ b/src/math/cmplx/asin.go
@@ -0,0 +1,221 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// The original C code, the long comment, and the constants
+// below are from http://netlib.sandia.gov/cephes/c9x-complex/clog.c.
+// The go code is a simplified version of the original C.
+//
+// Cephes Math Library Release 2.8:  June, 2000
+// Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
+//
+// The readme file at http://netlib.sandia.gov/cephes/ says:
+//    Some software in this archive may be from the book _Methods and
+// Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
+// International, 1989) or from the Cephes Mathematical Library, a
+// commercial product. In either event, it is copyrighted by the author.
+// What you see here may be used freely but it comes with no support or
+// guarantee.
+//
+//   The two known misprints in the book are repaired here in the
+// source listings for the gamma function and the incomplete beta
+// integral.
+//
+//   Stephen L. Moshier
+//   moshier@na-net.ornl.gov
+
+// Complex circular arc sine
+//
+// DESCRIPTION:
+//
+// Inverse complex sine:
+//                               2
+// w = -i clog( iz + csqrt( 1 - z ) ).
+//
+// casin(z) = -i casinh(iz)
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10     10100       2.1e-15     3.4e-16
+//    IEEE      -10,+10     30000       2.2e-14     2.7e-15
+// Larger relative error can be observed for z near zero.
+// Also tested by csin(casin(z)) = z.
+
+// Asin returns the inverse sine of x.
+func Asin(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case im == 0 && math.Abs(re) <= 1:
+		return complex(math.Asin(re), im)
+	case re == 0 && math.Abs(im) <= 1:
+		return complex(re, math.Asinh(im))
+	case math.IsNaN(im):
+		switch {
+		case re == 0:
+			return complex(re, math.NaN())
+		case math.IsInf(re, 0):
+			return complex(math.NaN(), re)
+		default:
+			return NaN()
+		}
+	case math.IsInf(im, 0):
+		switch {
+		case math.IsNaN(re):
+			return x
+		case math.IsInf(re, 0):
+			return complex(math.Copysign(math.Pi/4, re), im)
+		default:
+			return complex(math.Copysign(0, re), im)
+		}
+	case math.IsInf(re, 0):
+		return complex(math.Copysign(math.Pi/2, re), math.Copysign(re, im))
+	}
+	ct := complex(-imag(x), real(x)) // i * x
+	xx := x * x
+	x1 := complex(1-real(xx), -imag(xx)) // 1 - x*x
+	x2 := Sqrt(x1)                       // x2 = sqrt(1 - x*x)
+	w := Log(ct + x2)
+	return complex(imag(w), -real(w)) // -i * w
+}
+
+// Asinh returns the inverse hyperbolic sine of x.
+func Asinh(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case im == 0 && math.Abs(re) <= 1:
+		return complex(math.Asinh(re), im)
+	case re == 0 && math.Abs(im) <= 1:
+		return complex(re, math.Asin(im))
+	case math.IsInf(re, 0):
+		switch {
+		case math.IsInf(im, 0):
+			return complex(re, math.Copysign(math.Pi/4, im))
+		case math.IsNaN(im):
+			return x
+		default:
+			return complex(re, math.Copysign(0.0, im))
+		}
+	case math.IsNaN(re):
+		switch {
+		case im == 0:
+			return x
+		case math.IsInf(im, 0):
+			return complex(im, re)
+		default:
+			return NaN()
+		}
+	case math.IsInf(im, 0):
+		return complex(math.Copysign(im, re), math.Copysign(math.Pi/2, im))
+	}
+	xx := x * x
+	x1 := complex(1+real(xx), imag(xx)) // 1 + x*x
+	return Log(x + Sqrt(x1))            // log(x + sqrt(1 + x*x))
+}
+
+// Complex circular arc cosine
+//
+// DESCRIPTION:
+//
+// w = arccos z  =  PI/2 - arcsin z.
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      5200      1.6e-15      2.8e-16
+//    IEEE      -10,+10     30000      1.8e-14      2.2e-15
+
+// Acos returns the inverse cosine of x.
+func Acos(x complex128) complex128 {
+	w := Asin(x)
+	return complex(math.Pi/2-real(w), -imag(w))
+}
+
+// Acosh returns the inverse hyperbolic cosine of x.
+func Acosh(x complex128) complex128 {
+	if x == 0 {
+		return complex(0, math.Copysign(math.Pi/2, imag(x)))
+	}
+	w := Acos(x)
+	if imag(w) <= 0 {
+		return complex(-imag(w), real(w)) // i * w
+	}
+	return complex(imag(w), -real(w)) // -i * w
+}
+
+// Complex circular arc tangent
+//
+// DESCRIPTION:
+//
+// If
+//     z = x + iy,
+//
+// then
+//          1       (    2x     )
+// Re w  =  - arctan(-----------)  +  k PI
+//          2       (     2    2)
+//                  (1 - x  - y )
+//
+//               ( 2         2)
+//          1    (x  +  (y+1) )
+// Im w  =  - log(------------)
+//          4    ( 2         2)
+//               (x  +  (y-1) )
+//
+// Where k is an arbitrary integer.
+//
+// catan(z) = -i catanh(iz).
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      5900       1.3e-16     7.8e-18
+//    IEEE      -10,+10     30000       2.3e-15     8.5e-17
+// The check catan( ctan(z) )  =  z, with |x| and |y| < PI/2,
+// had peak relative error 1.5e-16, rms relative error
+// 2.9e-17.  See also clog().
+
+// Atan returns the inverse tangent of x.
+func Atan(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case im == 0:
+		return complex(math.Atan(re), im)
+	case re == 0 && math.Abs(im) <= 1:
+		return complex(re, math.Atanh(im))
+	case math.IsInf(im, 0) || math.IsInf(re, 0):
+		if math.IsNaN(re) {
+			return complex(math.NaN(), math.Copysign(0, im))
+		}
+		return complex(math.Copysign(math.Pi/2, re), math.Copysign(0, im))
+	case math.IsNaN(re) || math.IsNaN(im):
+		return NaN()
+	}
+	x2 := real(x) * real(x)
+	a := 1 - x2 - imag(x)*imag(x)
+	if a == 0 {
+		return NaN()
+	}
+	t := 0.5 * math.Atan2(2*real(x), a)
+	w := reducePi(t)
+
+	t = imag(x) - 1
+	b := x2 + t*t
+	if b == 0 {
+		return NaN()
+	}
+	t = imag(x) + 1
+	c := (x2 + t*t) / b
+	return complex(w, 0.25*math.Log(c))
+}
+
+// Atanh returns the inverse hyperbolic tangent of x.
+func Atanh(x complex128) complex128 {
+	z := complex(-imag(x), real(x)) // z = i * x
+	z = Atan(z)
+	return complex(imag(z), -real(z)) // z = -i * z
+}
diff --git a/src/math/cmplx/cmath_test.go b/src/math/cmplx/cmath_test.go
new file mode 100644
index 0000000..3011e83
--- /dev/null
+++ b/src/math/cmplx/cmath_test.go
@@ -0,0 +1,1589 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import (
+	"math"
+	"testing"
+)
+
+// The higher-precision values in vc26 were used to derive the
+// input arguments vc (see also comment below). For reference
+// only (do not delete).
+var vc26 = []complex128{
+	(4.97901192488367350108546816 + 7.73887247457810456552351752i),
+	(7.73887247457810456552351752 - 0.27688005719200159404635997i),
+	(-0.27688005719200159404635997 - 5.01060361827107492160848778i),
+	(-5.01060361827107492160848778 + 9.63629370719841737980004837i),
+	(9.63629370719841737980004837 + 2.92637723924396464525443662i),
+	(2.92637723924396464525443662 + 5.22908343145930665230025625i),
+	(5.22908343145930665230025625 + 2.72793991043601025126008608i),
+	(2.72793991043601025126008608 + 1.82530809168085506044576505i),
+	(1.82530809168085506044576505 - 8.68592476857560136238589621i),
+	(-8.68592476857560136238589621 + 4.97901192488367350108546816i),
+}
+
+var vc = []complex128{
+	(4.9790119248836735e+00 + 7.7388724745781045e+00i),
+	(7.7388724745781045e+00 - 2.7688005719200159e-01i),
+	(-2.7688005719200159e-01 - 5.0106036182710749e+00i),
+	(-5.0106036182710749e+00 + 9.6362937071984173e+00i),
+	(9.6362937071984173e+00 + 2.9263772392439646e+00i),
+	(2.9263772392439646e+00 + 5.2290834314593066e+00i),
+	(5.2290834314593066e+00 + 2.7279399104360102e+00i),
+	(2.7279399104360102e+00 + 1.8253080916808550e+00i),
+	(1.8253080916808550e+00 - 8.6859247685756013e+00i),
+	(-8.6859247685756013e+00 + 4.9790119248836735e+00i),
+}
+
+// The expected results below were computed by the high precision calculators
+// at https://keisan.casio.com/.  More exact input values (array vc[], above)
+// were obtained by printing them with "%.26f".  The answers were calculated
+// to 26 digits (by using the "Digit number" drop-down control of each
+// calculator).
+
+var abs = []float64{
+	9.2022120669932650313380972e+00,
+	7.7438239742296106616261394e+00,
+	5.0182478202557746902556648e+00,
+	1.0861137372799545160704002e+01,
+	1.0070841084922199607011905e+01,
+	5.9922447613166942183705192e+00,
+	5.8978784056736762299945176e+00,
+	3.2822866700678709020367184e+00,
+	8.8756430028990417290744307e+00,
+	1.0011785496777731986390856e+01,
+}
+
+var acos = []complex128{
+	(1.0017679804707456328694569 - 2.9138232718554953784519807i),
+	(0.03606427612041407369636057 + 2.7358584434576260925091256i),
+	(1.6249365462333796703711823 + 2.3159537454335901187730929i),
+	(2.0485650849650740120660391 - 3.0795576791204117911123886i),
+	(0.29621132089073067282488147 - 3.0007392508200622519398814i),
+	(1.0664555914934156601503632 - 2.4872865024796011364747111i),
+	(0.48681307452231387690013905 - 2.463655912283054555225301i),
+	(0.6116977071277574248407752 - 1.8734458851737055262693056i),
+	(1.3649311280370181331184214 + 2.8793528632328795424123832i),
+	(2.6189310485682988308904501 - 2.9956543302898767795858704i),
+}
+var acosh = []complex128{
+	(2.9138232718554953784519807 + 1.0017679804707456328694569i),
+	(2.7358584434576260925091256 - 0.03606427612041407369636057i),
+	(2.3159537454335901187730929 - 1.6249365462333796703711823i),
+	(3.0795576791204117911123886 + 2.0485650849650740120660391i),
+	(3.0007392508200622519398814 + 0.29621132089073067282488147i),
+	(2.4872865024796011364747111 + 1.0664555914934156601503632i),
+	(2.463655912283054555225301 + 0.48681307452231387690013905i),
+	(1.8734458851737055262693056 + 0.6116977071277574248407752i),
+	(2.8793528632328795424123832 - 1.3649311280370181331184214i),
+	(2.9956543302898767795858704 + 2.6189310485682988308904501i),
+}
+var asin = []complex128{
+	(0.56902834632415098636186476 + 2.9138232718554953784519807i),
+	(1.5347320506744825455349611 - 2.7358584434576260925091256i),
+	(-0.054140219438483051139860579 - 2.3159537454335901187730929i),
+	(-0.47776875817017739283471738 + 3.0795576791204117911123886i),
+	(1.2745850059041659464064402 + 3.0007392508200622519398814i),
+	(0.50434073530148095908095852 + 2.4872865024796011364747111i),
+	(1.0839832522725827423311826 + 2.463655912283054555225301i),
+	(0.9590986196671391943905465 + 1.8734458851737055262693056i),
+	(0.20586519875787848611290031 - 2.8793528632328795424123832i),
+	(-1.0481347217734022116591284 + 2.9956543302898767795858704i),
+}
+var asinh = []complex128{
+	(2.9113760469415295679342185 + 0.99639459545704326759805893i),
+	(2.7441755423994259061579029 - 0.035468308789000500601119392i),
+	(-2.2962136462520690506126678 - 1.5144663565690151885726707i),
+	(-3.0771233459295725965402455 + 1.0895577967194013849422294i),
+	(3.0048366100923647417557027 + 0.29346979169819220036454168i),
+	(2.4800059370795363157364643 + 1.0545868606049165710424232i),
+	(2.4718773838309585611141821 + 0.47502344364250803363708842i),
+	(1.8910743588080159144378396 + 0.56882925572563602341139174i),
+	(2.8735426423367341878069406 - 1.362376149648891420997548i),
+	(-2.9981750586172477217567878 + 0.5183571985225367505624207i),
+}
+var atan = []complex128{
+	(1.5115747079332741358607654 + 0.091324403603954494382276776i),
+	(1.4424504323482602560806727 - 0.0045416132642803911503770933i),
+	(-1.5593488703630532674484026 - 0.20163295409248362456446431i),
+	(-1.5280619472445889867794105 + 0.081721556230672003746956324i),
+	(1.4759909163240799678221039 + 0.028602969320691644358773586i),
+	(1.4877353772046548932715555 + 0.14566877153207281663773599i),
+	(1.4206983927779191889826 + 0.076830486127880702249439993i),
+	(1.3162236060498933364869556 + 0.16031313000467530644933363i),
+	(1.5473450684303703578810093 - 0.11064907507939082484935782i),
+	(-1.4841462340185253987375812 + 0.049341850305024399493142411i),
+}
+var atanh = []complex128{
+	(0.058375027938968509064640438 + 1.4793488495105334458167782i),
+	(0.12977343497790381229915667 - 1.5661009410463561327262499i),
+	(-0.010576456067347252072200088 - 1.3743698658402284549750563i),
+	(-0.042218595678688358882784918 + 1.4891433968166405606692604i),
+	(0.095218997991316722061828397 + 1.5416884098777110330499698i),
+	(0.079965459366890323857556487 + 1.4252510353873192700350435i),
+	(0.15051245471980726221708301 + 1.4907432533016303804884461i),
+	(0.25082072933993987714470373 + 1.392057665392187516442986i),
+	(0.022896108815797135846276662 - 1.4609224989282864208963021i),
+	(-0.08665624101841876130537396 + 1.5207902036935093480142159i),
+}
+var conj = []complex128{
+	(4.9790119248836735e+00 - 7.7388724745781045e+00i),
+	(7.7388724745781045e+00 + 2.7688005719200159e-01i),
+	(-2.7688005719200159e-01 + 5.0106036182710749e+00i),
+	(-5.0106036182710749e+00 - 9.6362937071984173e+00i),
+	(9.6362937071984173e+00 - 2.9263772392439646e+00i),
+	(2.9263772392439646e+00 - 5.2290834314593066e+00i),
+	(5.2290834314593066e+00 - 2.7279399104360102e+00i),
+	(2.7279399104360102e+00 - 1.8253080916808550e+00i),
+	(1.8253080916808550e+00 + 8.6859247685756013e+00i),
+	(-8.6859247685756013e+00 - 4.9790119248836735e+00i),
+}
+var cos = []complex128{
+	(3.024540920601483938336569e+02 + 1.1073797572517071650045357e+03i),
+	(1.192858682649064973252758e-01 + 2.7857554122333065540970207e-01i),
+	(7.2144394304528306603857962e+01 - 2.0500129667076044169954205e+01i),
+	(2.24921952538403984190541e+03 - 7.317363745602773587049329e+03i),
+	(-9.148222970032421760015498e+00 + 1.953124661113563541862227e+00i),
+	(-9.116081175857732248227078e+01 - 1.992669213569952232487371e+01i),
+	(3.795639179042704640002918e+00 + 6.623513350981458399309662e+00i),
+	(-2.9144840732498869560679084e+00 - 1.214620271628002917638748e+00i),
+	(-7.45123482501299743872481e+02 + 2.8641692314488080814066734e+03i),
+	(-5.371977967039319076416747e+01 + 4.893348341339375830564624e+01i),
+}
+var cosh = []complex128{
+	(8.34638383523018249366948e+00 + 7.2181057886425846415112064e+01i),
+	(1.10421967379919366952251e+03 - 3.1379638689277575379469861e+02i),
+	(3.051485206773701584738512e-01 - 2.6805384730105297848044485e-01i),
+	(-7.33294728684187933370938e+01 + 1.574445942284918251038144e+01i),
+	(-7.478643293945957535757355e+03 + 1.6348382209913353929473321e+03i),
+	(4.622316522966235701630926e+00 - 8.088695185566375256093098e+00i),
+	(-8.544333183278877406197712e+01 + 3.7505836120128166455231717e+01i),
+	(-1.934457815021493925115198e+00 + 7.3725859611767228178358673e+00i),
+	(-2.352958770061749348353548e+00 - 2.034982010440878358915409e+00i),
+	(7.79756457532134748165069e+02 + 2.8549350716819176560377717e+03i),
+}
+var exp = []complex128{
+	(1.669197736864670815125146e+01 + 1.4436895109507663689174096e+02i),
+	(2.2084389286252583447276212e+03 - 6.2759289284909211238261917e+02i),
+	(2.227538273122775173434327e-01 + 7.2468284028334191250470034e-01i),
+	(-6.5182985958153548997881627e-03 - 1.39965837915193860879044e-03i),
+	(-1.4957286524084015746110777e+04 + 3.269676455931135688988042e+03i),
+	(9.218158701983105935659273e+00 - 1.6223985291084956009304582e+01i),
+	(-1.7088175716853040841444505e+02 + 7.501382609870410713795546e+01i),
+	(-3.852461315830959613132505e+00 + 1.4808420423156073221970892e+01i),
+	(-4.586775503301407379786695e+00 - 4.178501081246873415144744e+00i),
+	(4.451337963005453491095747e-05 - 1.62977574205442915935263e-04i),
+}
+var log = []complex128{
+	(2.2194438972179194425697051e+00 + 9.9909115046919291062461269e-01i),
+	(2.0468956191154167256337289e+00 - 3.5762575021856971295156489e-02i),
+	(1.6130808329853860438751244e+00 - 1.6259990074019058442232221e+00i),
+	(2.3851910394823008710032651e+00 + 2.0502936359659111755031062e+00i),
+	(2.3096442270679923004800651e+00 + 2.9483213155446756211881774e-01i),
+	(1.7904660933974656106951860e+00 + 1.0605860367252556281902109e+00i),
+	(1.7745926939841751666177512e+00 + 4.8084556083358307819310911e-01i),
+	(1.1885403350045342425648780e+00 + 5.8969634164776659423195222e-01i),
+	(2.1833107837679082586772505e+00 - 1.3636647724582455028314573e+00i),
+	(2.3037629487273259170991671e+00 + 2.6210913895386013290915234e+00i),
+}
+var log10 = []complex128{
+	(9.6389223745559042474184943e-01 + 4.338997735671419492599631e-01i),
+	(8.8895547241376579493490892e-01 - 1.5531488990643548254864806e-02i),
+	(7.0055210462945412305244578e-01 - 7.0616239649481243222248404e-01i),
+	(1.0358753067322445311676952e+00 + 8.9043121238134980156490909e-01i),
+	(1.003065742975330237172029e+00 + 1.2804396782187887479857811e-01i),
+	(7.7758954439739162532085157e-01 + 4.6060666333341810869055108e-01i),
+	(7.7069581462315327037689152e-01 + 2.0882857371769952195512475e-01i),
+	(5.1617650901191156135137239e-01 + 2.5610186717615977620363299e-01i),
+	(9.4819982567026639742663212e-01 - 5.9223208584446952284914289e-01i),
+	(1.0005115362454417135973429e+00 + 1.1383255270407412817250921e+00i),
+}
+
+type ff struct {
+	r, theta float64
+}
+
+var polar = []ff{
+	{9.2022120669932650313380972e+00, 9.9909115046919291062461269e-01},
+	{7.7438239742296106616261394e+00, -3.5762575021856971295156489e-02},
+	{5.0182478202557746902556648e+00, -1.6259990074019058442232221e+00},
+	{1.0861137372799545160704002e+01, 2.0502936359659111755031062e+00},
+	{1.0070841084922199607011905e+01, 2.9483213155446756211881774e-01},
+	{5.9922447613166942183705192e+00, 1.0605860367252556281902109e+00},
+	{5.8978784056736762299945176e+00, 4.8084556083358307819310911e-01},
+	{3.2822866700678709020367184e+00, 5.8969634164776659423195222e-01},
+	{8.8756430028990417290744307e+00, -1.3636647724582455028314573e+00},
+	{1.0011785496777731986390856e+01, 2.6210913895386013290915234e+00},
+}
+var pow = []complex128{
+	(-2.499956739197529585028819e+00 + 1.759751724335650228957144e+00i),
+	(7.357094338218116311191939e+04 - 5.089973412479151648145882e+04i),
+	(1.320777296067768517259592e+01 - 3.165621914333901498921986e+01i),
+	(-3.123287828297300934072149e-07 - 1.9849567521490553032502223e-7i),
+	(8.0622651468477229614813e+04 - 7.80028727944573092944363e+04i),
+	(-1.0268824572103165858577141e+00 - 4.716844738244989776610672e-01i),
+	(-4.35953819012244175753187e+01 + 2.2036445974645306917648585e+02i),
+	(8.3556092283250594950239e-01 - 1.2261571947167240272593282e+01i),
+	(1.582292972120769306069625e+03 + 1.273564263524278244782512e+04i),
+	(6.592208301642122149025369e-08 + 2.584887236651661903526389e-08i),
+}
+var sin = []complex128{
+	(-1.1073801774240233539648544e+03 + 3.024539773002502192425231e+02i),
+	(1.0317037521400759359744682e+00 - 3.2208979799929570242818e-02i),
+	(-2.0501952097271429804261058e+01 - 7.2137981348240798841800967e+01i),
+	(7.3173638080346338642193078e+03 + 2.249219506193664342566248e+03i),
+	(-1.964375633631808177565226e+00 - 9.0958264713870404464159683e+00i),
+	(1.992783647158514838337674e+01 - 9.11555769410191350416942e+01i),
+	(-6.680335650741921444300349e+00 + 3.763353833142432513086117e+00i),
+	(1.2794028166657459148245993e+00 - 2.7669092099795781155109602e+00i),
+	(2.8641693949535259594188879e+03 + 7.451234399649871202841615e+02i),
+	(-4.893811726244659135553033e+01 - 5.371469305562194635957655e+01i),
+}
+var sinh = []complex128{
+	(8.34559353341652565758198e+00 + 7.2187893208650790476628899e+01i),
+	(1.1042192548260646752051112e+03 - 3.1379650595631635858792056e+02i),
+	(-8.239469336509264113041849e-02 + 9.9273668758439489098514519e-01i),
+	(7.332295456982297798219401e+01 - 1.574585908122833444899023e+01i),
+	(-7.4786432301380582103534216e+03 + 1.63483823493980029604071e+03i),
+	(4.595842179016870234028347e+00 - 8.135290105518580753211484e+00i),
+	(-8.543842533574163435246793e+01 + 3.750798997857594068272375e+01i),
+	(-1.918003500809465688017307e+00 + 7.4358344619793504041350251e+00i),
+	(-2.233816733239658031433147e+00 - 2.143519070805995056229335e+00i),
+	(-7.797564130187551181105341e+02 - 2.8549352346594918614806877e+03i),
+}
+var sqrt = []complex128{
+	(2.6628203086086130543813948e+00 + 1.4531345674282185229796902e+00i),
+	(2.7823278427251986247149295e+00 - 4.9756907317005224529115567e-02i),
+	(1.5397025302089642757361015e+00 - 1.6271336573016637535695727e+00i),
+	(1.7103411581506875260277898e+00 + 2.8170677122737589676157029e+00i),
+	(3.1390392472953103383607947e+00 + 4.6612625849858653248980849e-01i),
+	(2.1117080764822417640789287e+00 + 1.2381170223514273234967850e+00i),
+	(2.3587032281672256703926939e+00 + 5.7827111903257349935720172e-01i),
+	(1.7335262588873410476661577e+00 + 5.2647258220721269141550382e-01i),
+	(2.3131094974708716531499282e+00 - 1.8775429304303785570775490e+00i),
+	(8.1420535745048086240947359e-01 + 3.0575897587277248522656113e+00i),
+}
+var tan = []complex128{
+	(-1.928757919086441129134525e-07 + 1.0000003267499169073251826e+00i),
+	(1.242412685364183792138948e+00 - 3.17149693883133370106696e+00i),
+	(-4.6745126251587795225571826e-05 - 9.9992439225263959286114298e-01i),
+	(4.792363401193648192887116e-09 + 1.0000000070589333451557723e+00i),
+	(2.345740824080089140287315e-03 + 9.947733046570988661022763e-01i),
+	(-2.396030789494815566088809e-05 + 9.9994781345418591429826779e-01i),
+	(-7.370204836644931340905303e-03 + 1.0043553413417138987717748e+00i),
+	(-3.691803847992048527007457e-02 + 9.6475071993469548066328894e-01i),
+	(-2.781955256713729368401878e-08 - 1.000000049848910609006646e+00i),
+	(9.4281590064030478879791249e-05 + 9.9999119340863718183758545e-01i),
+}
+var tanh = []complex128{
+	(1.0000921981225144748819918e+00 + 2.160986245871518020231507e-05i),
+	(9.9999967727531993209562591e-01 - 1.9953763222959658873657676e-07i),
+	(-1.765485739548037260789686e+00 + 1.7024216325552852445168471e+00i),
+	(-9.999189442732736452807108e-01 + 3.64906070494473701938098e-05i),
+	(9.9999999224622333738729767e-01 - 3.560088949517914774813046e-09i),
+	(1.0029324933367326862499343e+00 - 4.948790309797102353137528e-03i),
+	(9.9996113064788012488693567e-01 - 4.226995742097032481451259e-05i),
+	(1.0074784189316340029873945e+00 - 4.194050814891697808029407e-03i),
+	(9.9385534229718327109131502e-01 + 5.144217985914355502713437e-02i),
+	(-1.0000000491604982429364892e+00 - 2.901873195374433112227349e-08i),
+}
+
+// huge values along the real axis for testing reducePi in Tan
+var hugeIn = []complex128{
+	1 << 28,
+	1 << 29,
+	1 << 30,
+	1 << 35,
+	-1 << 120,
+	1 << 240,
+	1 << 300,
+	-1 << 480,
+	1234567891234567 << 180,
+	-1234567891234567 << 300,
+}
+
+// Results for tanHuge[i] calculated with https://github.com/robpike/ivy
+// using 4096 bits of working precision.
+var tanHuge = []complex128{
+	5.95641897939639421,
+	-0.34551069233430392,
+	-0.78469661331920043,
+	0.84276385870875983,
+	0.40806638884180424,
+	-0.37603456702698076,
+	4.60901287677810962,
+	3.39135965054779932,
+	-6.76813854009065030,
+	-0.76417695016604922,
+}
+
+// special cases conform to C99 standard appendix G.6 Complex arithmetic
+var inf, nan = math.Inf(1), math.NaN()
+
+var vcAbsSC = []complex128{
+	NaN(),
+}
+var absSC = []float64{
+	math.NaN(),
+}
+var acosSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.1.1
+	{complex(zero, zero),
+		complex(math.Pi/2, -zero)},
+	{complex(-zero, zero),
+		complex(math.Pi/2, -zero)},
+	{complex(zero, nan),
+		complex(math.Pi/2, nan)},
+	{complex(-zero, nan),
+		complex(math.Pi/2, nan)},
+	{complex(1.0, inf),
+		complex(math.Pi/2, -inf)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(-inf, 1.0),
+		complex(math.Pi, -inf)},
+	{complex(inf, 1.0),
+		complex(0.0, -inf)},
+	{complex(-inf, inf),
+		complex(3*math.Pi/4, -inf)},
+	{complex(inf, inf),
+		complex(math.Pi/4, -inf)},
+	{complex(inf, nan),
+		complex(nan, -inf)}, // imaginary sign unspecified
+	{complex(-inf, nan),
+		complex(nan, inf)}, // imaginary sign unspecified
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(nan, -inf)},
+	{NaN(),
+		NaN()},
+}
+var acoshSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.2.1
+	{complex(zero, zero),
+		complex(zero, math.Pi/2)},
+	{complex(-zero, zero),
+		complex(zero, math.Pi/2)},
+	{complex(1.0, inf),
+		complex(inf, math.Pi/2)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(-inf, 1.0),
+		complex(inf, math.Pi)},
+	{complex(inf, 1.0),
+		complex(inf, zero)},
+	{complex(-inf, inf),
+		complex(inf, 3*math.Pi/4)},
+	{complex(inf, inf),
+		complex(inf, math.Pi/4)},
+	{complex(inf, nan),
+		complex(inf, nan)},
+	{complex(-inf, nan),
+		complex(inf, nan)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(inf, nan)},
+	{NaN(),
+		NaN()},
+}
+var asinSC = []struct {
+	in,
+	want complex128
+}{
+	// Derived from Asin(z) = -i * Asinh(i * z), G.6 #7
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(1.0, inf),
+		complex(0, inf)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, 1),
+		complex(math.Pi/2, inf)},
+	{complex(inf, inf),
+		complex(math.Pi/4, inf)},
+	{complex(inf, nan),
+		complex(nan, inf)}, // imaginary sign unspecified
+	{complex(nan, zero),
+		NaN()},
+	{complex(nan, 1),
+		NaN()},
+	{complex(nan, inf),
+		complex(nan, inf)},
+	{NaN(),
+		NaN()},
+}
+var asinhSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.2.2
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(1.0, inf),
+		complex(inf, math.Pi/2)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, 1.0),
+		complex(inf, zero)},
+	{complex(inf, inf),
+		complex(inf, math.Pi/4)},
+	{complex(inf, nan),
+		complex(inf, nan)},
+	{complex(nan, zero),
+		complex(nan, zero)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(inf, nan)}, // sign of real part unspecified
+	{NaN(),
+		NaN()},
+}
+var atanSC = []struct {
+	in,
+	want complex128
+}{
+	// Derived from Atan(z) = -i * Atanh(i * z), G.6 #7
+	{complex(0, zero),
+		complex(0, zero)},
+	{complex(0, nan),
+		NaN()},
+	{complex(1.0, zero),
+		complex(math.Pi/4, zero)},
+	{complex(1.0, inf),
+		complex(math.Pi/2, zero)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, 1),
+		complex(math.Pi/2, zero)},
+	{complex(inf, inf),
+		complex(math.Pi/2, zero)},
+	{complex(inf, nan),
+		complex(math.Pi/2, zero)},
+	{complex(nan, 1),
+		NaN()},
+	{complex(nan, inf),
+		complex(nan, zero)},
+	{NaN(),
+		NaN()},
+}
+var atanhSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.2.3
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(zero, nan),
+		complex(zero, nan)},
+	{complex(1.0, zero),
+		complex(inf, zero)},
+	{complex(1.0, inf),
+		complex(0, math.Pi/2)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, 1.0),
+		complex(zero, math.Pi/2)},
+	{complex(inf, inf),
+		complex(zero, math.Pi/2)},
+	{complex(inf, nan),
+		complex(0, nan)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(zero, math.Pi/2)}, // sign of real part not specified.
+	{NaN(),
+		NaN()},
+}
+var vcConjSC = []complex128{
+	NaN(),
+}
+var conjSC = []complex128{
+	NaN(),
+}
+var cosSC = []struct {
+	in,
+	want complex128
+}{
+	// Derived from Cos(z) = Cosh(i * z), G.6 #7
+	{complex(zero, zero),
+		complex(1.0, -zero)},
+	{complex(zero, inf),
+		complex(inf, -zero)},
+	{complex(zero, nan),
+		complex(nan, zero)}, // imaginary sign unspecified
+	{complex(1.0, inf),
+		complex(inf, -inf)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, zero),
+		complex(nan, -zero)},
+	{complex(inf, 1.0),
+		NaN()},
+	{complex(inf, inf),
+		complex(inf, nan)}, // real sign unspecified
+	{complex(inf, nan),
+		NaN()},
+	{complex(nan, zero),
+		complex(nan, -zero)}, // imaginary sign unspecified
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(inf, nan)},
+	{NaN(),
+		NaN()},
+}
+var coshSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.2.4
+	{complex(zero, zero),
+		complex(1.0, zero)},
+	{complex(zero, inf),
+		complex(nan, zero)}, // imaginary sign unspecified
+	{complex(zero, nan),
+		complex(nan, zero)}, // imaginary sign unspecified
+	{complex(1.0, inf),
+		NaN()},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, zero),
+		complex(inf, zero)},
+	{complex(inf, 1.0),
+		complex(inf*math.Cos(1.0), inf*math.Sin(1.0))}, // +inf  cis(y)
+	{complex(inf, inf),
+		complex(inf, nan)}, // real sign unspecified
+	{complex(inf, nan),
+		complex(inf, nan)},
+	{complex(nan, zero),
+		complex(nan, zero)}, // imaginary sign unspecified
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		NaN()},
+	{NaN(),
+		NaN()},
+}
+var expSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.3.1
+	{complex(zero, zero),
+		complex(1.0, zero)},
+	{complex(-zero, zero),
+		complex(1.0, zero)},
+	{complex(1.0, inf),
+		NaN()},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, zero),
+		complex(inf, zero)},
+	{complex(-inf, 1.0),
+		complex(math.Copysign(0.0, math.Cos(1.0)), math.Copysign(0.0, math.Sin(1.0)))}, // +0 cis(y)
+	{complex(inf, 1.0),
+		complex(inf*math.Cos(1.0), inf*math.Sin(1.0))}, // +inf  cis(y)
+	{complex(-inf, inf),
+		complex(zero, zero)}, // real and imaginary sign unspecified
+	{complex(inf, inf),
+		complex(inf, nan)}, // real sign unspecified
+	{complex(-inf, nan),
+		complex(zero, zero)}, // real and imaginary sign unspecified
+	{complex(inf, nan),
+		complex(inf, nan)}, // real sign unspecified
+	{complex(nan, zero),
+		complex(nan, zero)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		NaN()},
+	{NaN(),
+		NaN()},
+}
+var vcIsNaNSC = []complex128{
+	complex(math.Inf(-1), math.Inf(-1)),
+	complex(math.Inf(-1), math.NaN()),
+	complex(math.NaN(), math.Inf(-1)),
+	complex(0, math.NaN()),
+	complex(math.NaN(), 0),
+	complex(math.Inf(1), math.Inf(1)),
+	complex(math.Inf(1), math.NaN()),
+	complex(math.NaN(), math.Inf(1)),
+	complex(math.NaN(), math.NaN()),
+}
+var isNaNSC = []bool{
+	false,
+	false,
+	false,
+	true,
+	true,
+	false,
+	false,
+	false,
+	true,
+}
+
+var logSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.3.2
+	{complex(zero, zero),
+		complex(-inf, zero)},
+	{complex(-zero, zero),
+		complex(-inf, math.Pi)},
+	{complex(1.0, inf),
+		complex(inf, math.Pi/2)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(-inf, 1.0),
+		complex(inf, math.Pi)},
+	{complex(inf, 1.0),
+		complex(inf, 0.0)},
+	{complex(-inf, inf),
+		complex(inf, 3*math.Pi/4)},
+	{complex(inf, inf),
+		complex(inf, math.Pi/4)},
+	{complex(-inf, nan),
+		complex(inf, nan)},
+	{complex(inf, nan),
+		complex(inf, nan)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(inf, nan)},
+	{NaN(),
+		NaN()},
+}
+var log10SC = []struct {
+	in,
+	want complex128
+}{
+	// derived from Log special cases via Log10(x) = math.Log10E*Log(x)
+	{complex(zero, zero),
+		complex(-inf, zero)},
+	{complex(-zero, zero),
+		complex(-inf, float64(math.Log10E)*float64(math.Pi))},
+	{complex(1.0, inf),
+		complex(inf, float64(math.Log10E)*float64(math.Pi/2))},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(-inf, 1.0),
+		complex(inf, float64(math.Log10E)*float64(math.Pi))},
+	{complex(inf, 1.0),
+		complex(inf, 0.0)},
+	{complex(-inf, inf),
+		complex(inf, float64(math.Log10E)*float64(3*math.Pi/4))},
+	{complex(inf, inf),
+		complex(inf, float64(math.Log10E)*float64(math.Pi/4))},
+	{complex(-inf, nan),
+		complex(inf, nan)},
+	{complex(inf, nan),
+		complex(inf, nan)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(inf, nan)},
+	{NaN(),
+		NaN()},
+}
+var vcPolarSC = []complex128{
+	NaN(),
+}
+var polarSC = []ff{
+	{math.NaN(), math.NaN()},
+}
+var vcPowSC = [][2]complex128{
+	{NaN(), NaN()},
+	{0, NaN()},
+}
+var powSC = []complex128{
+	NaN(),
+	NaN(),
+}
+var sinSC = []struct {
+	in,
+	want complex128
+}{
+	// Derived from Sin(z) = -i * Sinh(i * z), G.6 #7
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(zero, inf),
+		complex(zero, inf)},
+	{complex(zero, nan),
+		complex(zero, nan)},
+	{complex(1.0, inf),
+		complex(inf, inf)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, zero),
+		complex(nan, zero)},
+	{complex(inf, 1.0),
+		NaN()},
+	{complex(inf, inf),
+		complex(nan, inf)},
+	{complex(inf, nan),
+		NaN()},
+	{complex(nan, zero),
+		complex(nan, zero)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(nan, inf)},
+	{NaN(),
+		NaN()},
+}
+
+var sinhSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.2.5
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(zero, inf),
+		complex(zero, nan)}, // real sign unspecified
+	{complex(zero, nan),
+		complex(zero, nan)}, // real sign unspecified
+	{complex(1.0, inf),
+		NaN()},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, zero),
+		complex(inf, zero)},
+	{complex(inf, 1.0),
+		complex(inf*math.Cos(1.0), inf*math.Sin(1.0))}, // +inf  cis(y)
+	{complex(inf, inf),
+		complex(inf, nan)}, // real sign unspecified
+	{complex(inf, nan),
+		complex(inf, nan)}, // real sign unspecified
+	{complex(nan, zero),
+		complex(nan, zero)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		NaN()},
+	{NaN(),
+		NaN()},
+}
+
+var sqrtSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.4.2
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(-zero, zero),
+		complex(zero, zero)},
+	{complex(1.0, inf),
+		complex(inf, inf)},
+	{complex(nan, inf),
+		complex(inf, inf)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(-inf, 1.0),
+		complex(zero, inf)},
+	{complex(inf, 1.0),
+		complex(inf, zero)},
+	{complex(-inf, nan),
+		complex(nan, inf)}, // imaginary sign unspecified
+	{complex(inf, nan),
+		complex(inf, nan)},
+	{complex(nan, 1.0),
+		NaN()},
+	{NaN(),
+		NaN()},
+}
+var tanSC = []struct {
+	in,
+	want complex128
+}{
+	// Derived from Tan(z) = -i * Tanh(i * z), G.6 #7
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(zero, nan),
+		complex(zero, nan)},
+	{complex(1.0, inf),
+		complex(zero, 1.0)},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, 1.0),
+		NaN()},
+	{complex(inf, inf),
+		complex(zero, 1.0)},
+	{complex(inf, nan),
+		NaN()},
+	{complex(nan, zero),
+		NaN()},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		complex(zero, 1.0)},
+	{NaN(),
+		NaN()},
+}
+var tanhSC = []struct {
+	in,
+	want complex128
+}{
+	// G.6.2.6
+	{complex(zero, zero),
+		complex(zero, zero)},
+	{complex(1.0, inf),
+		NaN()},
+	{complex(1.0, nan),
+		NaN()},
+	{complex(inf, 1.0),
+		complex(1.0, math.Copysign(0.0, math.Sin(2*1.0)))}, // 1 + i 0 sin(2y)
+	{complex(inf, inf),
+		complex(1.0, zero)}, // imaginary sign unspecified
+	{complex(inf, nan),
+		complex(1.0, zero)}, // imaginary sign unspecified
+	{complex(nan, zero),
+		complex(nan, zero)},
+	{complex(nan, 1.0),
+		NaN()},
+	{complex(nan, inf),
+		NaN()},
+	{NaN(),
+		NaN()},
+}
+
+// branch cut continuity checks
+// points on each axis at |z| > 1 are checked for one-sided continuity from both the positive and negative side
+// all possible branch cuts for the elementary functions are at one of these points
+
+var zero = 0.0
+var eps = 1.0 / (1 << 53)
+
+var branchPoints = [][2]complex128{
+	{complex(2.0, zero), complex(2.0, eps)},
+	{complex(2.0, -zero), complex(2.0, -eps)},
+	{complex(-2.0, zero), complex(-2.0, eps)},
+	{complex(-2.0, -zero), complex(-2.0, -eps)},
+	{complex(zero, 2.0), complex(eps, 2.0)},
+	{complex(-zero, 2.0), complex(-eps, 2.0)},
+	{complex(zero, -2.0), complex(eps, -2.0)},
+	{complex(-zero, -2.0), complex(-eps, -2.0)},
+}
+
+// functions borrowed from pkg/math/all_test.go
+func tolerance(a, b, e float64) bool {
+	d := a - b
+	if d < 0 {
+		d = -d
+	}
+
+	// note: b is correct (expected) value, a is actual value.
+	// make error tolerance a fraction of b, not a.
+	if b != 0 {
+		e = e * b
+		if e < 0 {
+			e = -e
+		}
+	}
+	return d < e
+}
+func veryclose(a, b float64) bool { return tolerance(a, b, 4e-16) }
+func alike(a, b float64) bool {
+	switch {
+	case a != a && b != b: // math.IsNaN(a) && math.IsNaN(b):
+		return true
+	case a == b:
+		return math.Signbit(a) == math.Signbit(b)
+	}
+	return false
+}
+
+func cTolerance(a, b complex128, e float64) bool {
+	d := Abs(a - b)
+	if b != 0 {
+		e = e * Abs(b)
+		if e < 0 {
+			e = -e
+		}
+	}
+	return d < e
+}
+func cSoclose(a, b complex128, e float64) bool { return cTolerance(a, b, e) }
+func cVeryclose(a, b complex128) bool          { return cTolerance(a, b, 4e-16) }
+func cAlike(a, b complex128) bool {
+	var realAlike, imagAlike bool
+	if isExact(real(b)) {
+		realAlike = alike(real(a), real(b))
+	} else {
+		// Allow non-exact special cases to have errors in ULP.
+		realAlike = veryclose(real(a), real(b))
+	}
+	if isExact(imag(b)) {
+		imagAlike = alike(imag(a), imag(b))
+	} else {
+		// Allow non-exact special cases to have errors in ULP.
+		imagAlike = veryclose(imag(a), imag(b))
+	}
+	return realAlike && imagAlike
+}
+func isExact(x float64) bool {
+	// Special cases that should match exactly.  Other cases are multiples
+	// of Pi that may not be last bit identical on all platforms.
+	return math.IsNaN(x) || math.IsInf(x, 0) || x == 0 || x == 1 || x == -1
+}
+
+func TestAbs(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Abs(vc[i]); !veryclose(abs[i], f) {
+			t.Errorf("Abs(%g) = %g, want %g", vc[i], f, abs[i])
+		}
+	}
+	for i := 0; i < len(vcAbsSC); i++ {
+		if f := Abs(vcAbsSC[i]); !alike(absSC[i], f) {
+			t.Errorf("Abs(%g) = %g, want %g", vcAbsSC[i], f, absSC[i])
+		}
+	}
+}
+func TestAcos(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Acos(vc[i]); !cSoclose(acos[i], f, 1e-14) {
+			t.Errorf("Acos(%g) = %g, want %g", vc[i], f, acos[i])
+		}
+	}
+	for _, v := range acosSC {
+		if f := Acos(v.in); !cAlike(v.want, f) {
+			t.Errorf("Acos(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Acos(Conj(z))  == Conj(Acos(z))
+		if f := Acos(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Acos(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Acos(pt[0]), Acos(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Acos(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestAcosh(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Acosh(vc[i]); !cSoclose(acosh[i], f, 1e-14) {
+			t.Errorf("Acosh(%g) = %g, want %g", vc[i], f, acosh[i])
+		}
+	}
+	for _, v := range acoshSC {
+		if f := Acosh(v.in); !cAlike(v.want, f) {
+			t.Errorf("Acosh(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Acosh(Conj(z))  == Conj(Acosh(z))
+		if f := Acosh(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Acosh(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Acosh(pt[0]), Acosh(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Acosh(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestAsin(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Asin(vc[i]); !cSoclose(asin[i], f, 1e-14) {
+			t.Errorf("Asin(%g) = %g, want %g", vc[i], f, asin[i])
+		}
+	}
+	for _, v := range asinSC {
+		if f := Asin(v.in); !cAlike(v.want, f) {
+			t.Errorf("Asin(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Asin(Conj(z))  == Asin(Sinh(z))
+		if f := Asin(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Asin(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Asin(-z)  == -Asin(z)
+		if f := Asin(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Asin(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Asin(pt[0]), Asin(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Asin(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestAsinh(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Asinh(vc[i]); !cSoclose(asinh[i], f, 4e-15) {
+			t.Errorf("Asinh(%g) = %g, want %g", vc[i], f, asinh[i])
+		}
+	}
+	for _, v := range asinhSC {
+		if f := Asinh(v.in); !cAlike(v.want, f) {
+			t.Errorf("Asinh(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Asinh(Conj(z))  == Asinh(Sinh(z))
+		if f := Asinh(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Asinh(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Asinh(-z)  == -Asinh(z)
+		if f := Asinh(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Asinh(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Asinh(pt[0]), Asinh(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Asinh(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestAtan(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Atan(vc[i]); !cVeryclose(atan[i], f) {
+			t.Errorf("Atan(%g) = %g, want %g", vc[i], f, atan[i])
+		}
+	}
+	for _, v := range atanSC {
+		if f := Atan(v.in); !cAlike(v.want, f) {
+			t.Errorf("Atan(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Atan(Conj(z))  == Conj(Atan(z))
+		if f := Atan(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Atan(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Atan(-z)  == -Atan(z)
+		if f := Atan(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Atan(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Atan(pt[0]), Atan(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Atan(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestAtanh(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Atanh(vc[i]); !cVeryclose(atanh[i], f) {
+			t.Errorf("Atanh(%g) = %g, want %g", vc[i], f, atanh[i])
+		}
+	}
+	for _, v := range atanhSC {
+		if f := Atanh(v.in); !cAlike(v.want, f) {
+			t.Errorf("Atanh(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Atanh(Conj(z))  == Conj(Atanh(z))
+		if f := Atanh(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Atanh(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Atanh(-z)  == -Atanh(z)
+		if f := Atanh(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Atanh(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Atanh(pt[0]), Atanh(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Atanh(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestConj(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Conj(vc[i]); !cVeryclose(conj[i], f) {
+			t.Errorf("Conj(%g) = %g, want %g", vc[i], f, conj[i])
+		}
+	}
+	for i := 0; i < len(vcConjSC); i++ {
+		if f := Conj(vcConjSC[i]); !cAlike(conjSC[i], f) {
+			t.Errorf("Conj(%g) = %g, want %g", vcConjSC[i], f, conjSC[i])
+		}
+	}
+}
+func TestCos(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Cos(vc[i]); !cSoclose(cos[i], f, 3e-15) {
+			t.Errorf("Cos(%g) = %g, want %g", vc[i], f, cos[i])
+		}
+	}
+	for _, v := range cosSC {
+		if f := Cos(v.in); !cAlike(v.want, f) {
+			t.Errorf("Cos(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Cos(Conj(z))  == Cos(Cosh(z))
+		if f := Cos(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Cos(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Cos(-z)  == Cos(z)
+		if f := Cos(-v.in); !cAlike(v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Cos(%g) = %g, want %g", -v.in, f, v.want)
+		}
+	}
+}
+func TestCosh(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Cosh(vc[i]); !cSoclose(cosh[i], f, 2e-15) {
+			t.Errorf("Cosh(%g) = %g, want %g", vc[i], f, cosh[i])
+		}
+	}
+	for _, v := range coshSC {
+		if f := Cosh(v.in); !cAlike(v.want, f) {
+			t.Errorf("Cosh(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Cosh(Conj(z))  == Conj(Cosh(z))
+		if f := Cosh(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Cosh(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Cosh(-z)  == Cosh(z)
+		if f := Cosh(-v.in); !cAlike(v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Cosh(%g) = %g, want %g", -v.in, f, v.want)
+		}
+	}
+}
+func TestExp(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Exp(vc[i]); !cSoclose(exp[i], f, 1e-15) {
+			t.Errorf("Exp(%g) = %g, want %g", vc[i], f, exp[i])
+		}
+	}
+	for _, v := range expSC {
+		if f := Exp(v.in); !cAlike(v.want, f) {
+			t.Errorf("Exp(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Exp(Conj(z))  == Exp(Cosh(z))
+		if f := Exp(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Exp(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+	}
+}
+func TestIsNaN(t *testing.T) {
+	for i := 0; i < len(vcIsNaNSC); i++ {
+		if f := IsNaN(vcIsNaNSC[i]); isNaNSC[i] != f {
+			t.Errorf("IsNaN(%v) = %v, want %v", vcIsNaNSC[i], f, isNaNSC[i])
+		}
+	}
+}
+func TestLog(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Log(vc[i]); !cVeryclose(log[i], f) {
+			t.Errorf("Log(%g) = %g, want %g", vc[i], f, log[i])
+		}
+	}
+	for _, v := range logSC {
+		if f := Log(v.in); !cAlike(v.want, f) {
+			t.Errorf("Log(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Log(Conj(z))  == Conj(Log(z))
+		if f := Log(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Log(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Log(pt[0]), Log(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Log(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestLog10(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Log10(vc[i]); !cVeryclose(log10[i], f) {
+			t.Errorf("Log10(%g) = %g, want %g", vc[i], f, log10[i])
+		}
+	}
+	for _, v := range log10SC {
+		if f := Log10(v.in); !cAlike(v.want, f) {
+			t.Errorf("Log10(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Log10(Conj(z))  == Conj(Log10(z))
+		if f := Log10(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Log10(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+	}
+}
+func TestPolar(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if r, theta := Polar(vc[i]); !veryclose(polar[i].r, r) && !veryclose(polar[i].theta, theta) {
+			t.Errorf("Polar(%g) = %g, %g want %g, %g", vc[i], r, theta, polar[i].r, polar[i].theta)
+		}
+	}
+	for i := 0; i < len(vcPolarSC); i++ {
+		if r, theta := Polar(vcPolarSC[i]); !alike(polarSC[i].r, r) && !alike(polarSC[i].theta, theta) {
+			t.Errorf("Polar(%g) = %g, %g, want %g, %g", vcPolarSC[i], r, theta, polarSC[i].r, polarSC[i].theta)
+		}
+	}
+}
+func TestPow(t *testing.T) {
+	// Special cases for Pow(0, c).
+	var zero = complex(0, 0)
+	zeroPowers := [][2]complex128{
+		{0, 1 + 0i},
+		{1.5, 0 + 0i},
+		{-1.5, complex(math.Inf(0), 0)},
+		{-1.5 + 1.5i, Inf()},
+	}
+	for _, zp := range zeroPowers {
+		if f := Pow(zero, zp[0]); f != zp[1] {
+			t.Errorf("Pow(%g, %g) = %g, want %g", zero, zp[0], f, zp[1])
+		}
+	}
+	var a = complex(3.0, 3.0)
+	for i := 0; i < len(vc); i++ {
+		if f := Pow(a, vc[i]); !cSoclose(pow[i], f, 4e-15) {
+			t.Errorf("Pow(%g, %g) = %g, want %g", a, vc[i], f, pow[i])
+		}
+	}
+	for i := 0; i < len(vcPowSC); i++ {
+		if f := Pow(vcPowSC[i][0], vcPowSC[i][1]); !cAlike(powSC[i], f) {
+			t.Errorf("Pow(%g, %g) = %g, want %g", vcPowSC[i][0], vcPowSC[i][1], f, powSC[i])
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Pow(pt[0], 0.1), Pow(pt[1], 0.1); !cVeryclose(f0, f1) {
+			t.Errorf("Pow(%g, 0.1) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestRect(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Rect(polar[i].r, polar[i].theta); !cVeryclose(vc[i], f) {
+			t.Errorf("Rect(%g, %g) = %g want %g", polar[i].r, polar[i].theta, f, vc[i])
+		}
+	}
+	for i := 0; i < len(vcPolarSC); i++ {
+		if f := Rect(polarSC[i].r, polarSC[i].theta); !cAlike(vcPolarSC[i], f) {
+			t.Errorf("Rect(%g, %g) = %g, want %g", polarSC[i].r, polarSC[i].theta, f, vcPolarSC[i])
+		}
+	}
+}
+func TestSin(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Sin(vc[i]); !cSoclose(sin[i], f, 2e-15) {
+			t.Errorf("Sin(%g) = %g, want %g", vc[i], f, sin[i])
+		}
+	}
+	for _, v := range sinSC {
+		if f := Sin(v.in); !cAlike(v.want, f) {
+			t.Errorf("Sin(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Sin(Conj(z))  == Conj(Sin(z))
+		if f := Sin(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Sinh(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Sin(-z)  == -Sin(z)
+		if f := Sin(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Sinh(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+}
+func TestSinh(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Sinh(vc[i]); !cSoclose(sinh[i], f, 2e-15) {
+			t.Errorf("Sinh(%g) = %g, want %g", vc[i], f, sinh[i])
+		}
+	}
+	for _, v := range sinhSC {
+		if f := Sinh(v.in); !cAlike(v.want, f) {
+			t.Errorf("Sinh(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Sinh(Conj(z))  == Conj(Sinh(z))
+		if f := Sinh(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Sinh(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Sinh(-z)  == -Sinh(z)
+		if f := Sinh(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Sinh(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+}
+func TestSqrt(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Sqrt(vc[i]); !cVeryclose(sqrt[i], f) {
+			t.Errorf("Sqrt(%g) = %g, want %g", vc[i], f, sqrt[i])
+		}
+	}
+	for _, v := range sqrtSC {
+		if f := Sqrt(v.in); !cAlike(v.want, f) {
+			t.Errorf("Sqrt(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Sqrt(Conj(z)) == Conj(Sqrt(z))
+		if f := Sqrt(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Sqrt(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+	}
+	for _, pt := range branchPoints {
+		if f0, f1 := Sqrt(pt[0]), Sqrt(pt[1]); !cVeryclose(f0, f1) {
+			t.Errorf("Sqrt(%g) not continuous, got %g want %g", pt[0], f0, f1)
+		}
+	}
+}
+func TestTan(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Tan(vc[i]); !cSoclose(tan[i], f, 3e-15) {
+			t.Errorf("Tan(%g) = %g, want %g", vc[i], f, tan[i])
+		}
+	}
+	for _, v := range tanSC {
+		if f := Tan(v.in); !cAlike(v.want, f) {
+			t.Errorf("Tan(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Tan(Conj(z))  == Conj(Tan(z))
+		if f := Tan(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Tan(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Tan(-z)  == -Tan(z)
+		if f := Tan(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Tan(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+}
+func TestTanh(t *testing.T) {
+	for i := 0; i < len(vc); i++ {
+		if f := Tanh(vc[i]); !cSoclose(tanh[i], f, 2e-15) {
+			t.Errorf("Tanh(%g) = %g, want %g", vc[i], f, tanh[i])
+		}
+	}
+	for _, v := range tanhSC {
+		if f := Tanh(v.in); !cAlike(v.want, f) {
+			t.Errorf("Tanh(%g) = %g, want %g", v.in, f, v.want)
+		}
+		if math.IsNaN(imag(v.in)) || math.IsNaN(imag(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Tanh(Conj(z))  == Conj(Tanh(z))
+		if f := Tanh(Conj(v.in)); !cAlike(Conj(v.want), f) && !cAlike(v.in, Conj(v.in)) {
+			t.Errorf("Tanh(%g) = %g, want %g", Conj(v.in), f, Conj(v.want))
+		}
+		if math.IsNaN(real(v.in)) || math.IsNaN(real(v.want)) {
+			// Negating NaN is undefined with regard to the sign bit produced.
+			continue
+		}
+		// Tanh(-z)  == -Tanh(z)
+		if f := Tanh(-v.in); !cAlike(-v.want, f) && !cAlike(v.in, -v.in) {
+			t.Errorf("Tanh(%g) = %g, want %g", -v.in, f, -v.want)
+		}
+	}
+}
+
+// See issue 17577
+func TestInfiniteLoopIntanSeries(t *testing.T) {
+	want := Inf()
+	if got := Cot(0); got != want {
+		t.Errorf("Cot(0): got %g, want %g", got, want)
+	}
+}
+
+func BenchmarkAbs(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Abs(complex(2.5, 3.5))
+	}
+}
+func BenchmarkAcos(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Acos(complex(2.5, 3.5))
+	}
+}
+func BenchmarkAcosh(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Acosh(complex(2.5, 3.5))
+	}
+}
+func BenchmarkAsin(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Asin(complex(2.5, 3.5))
+	}
+}
+func BenchmarkAsinh(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Asinh(complex(2.5, 3.5))
+	}
+}
+func BenchmarkAtan(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Atan(complex(2.5, 3.5))
+	}
+}
+func BenchmarkAtanh(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Atanh(complex(2.5, 3.5))
+	}
+}
+func BenchmarkConj(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Conj(complex(2.5, 3.5))
+	}
+}
+func BenchmarkCos(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Cos(complex(2.5, 3.5))
+	}
+}
+func BenchmarkCosh(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Cosh(complex(2.5, 3.5))
+	}
+}
+func BenchmarkExp(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Exp(complex(2.5, 3.5))
+	}
+}
+func BenchmarkLog(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Log(complex(2.5, 3.5))
+	}
+}
+func BenchmarkLog10(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Log10(complex(2.5, 3.5))
+	}
+}
+func BenchmarkPhase(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Phase(complex(2.5, 3.5))
+	}
+}
+func BenchmarkPolar(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Polar(complex(2.5, 3.5))
+	}
+}
+func BenchmarkPow(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Pow(complex(2.5, 3.5), complex(2.5, 3.5))
+	}
+}
+func BenchmarkRect(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Rect(2.5, 1.5)
+	}
+}
+func BenchmarkSin(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Sin(complex(2.5, 3.5))
+	}
+}
+func BenchmarkSinh(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Sinh(complex(2.5, 3.5))
+	}
+}
+func BenchmarkSqrt(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Sqrt(complex(2.5, 3.5))
+	}
+}
+func BenchmarkTan(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Tan(complex(2.5, 3.5))
+	}
+}
+func BenchmarkTanh(b *testing.B) {
+	for i := 0; i < b.N; i++ {
+		Tanh(complex(2.5, 3.5))
+	}
+}
diff --git a/src/math/cmplx/conj.go b/src/math/cmplx/conj.go
new file mode 100644
index 0000000..34a4277
--- /dev/null
+++ b/src/math/cmplx/conj.go
@@ -0,0 +1,8 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+// Conj returns the complex conjugate of x.
+func Conj(x complex128) complex128 { return complex(real(x), -imag(x)) }
diff --git a/src/math/cmplx/example_test.go b/src/math/cmplx/example_test.go
new file mode 100644
index 0000000..f0ed963
--- /dev/null
+++ b/src/math/cmplx/example_test.go
@@ -0,0 +1,28 @@
+// Copyright 2016 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx_test
+
+import (
+	"fmt"
+	"math"
+	"math/cmplx"
+)
+
+func ExampleAbs() {
+	fmt.Printf("%.1f", cmplx.Abs(3+4i))
+	// Output: 5.0
+}
+
+// ExampleExp computes Euler's identity.
+func ExampleExp() {
+	fmt.Printf("%.1f", cmplx.Exp(1i*math.Pi)+1)
+	// Output: (0.0+0.0i)
+}
+
+func ExamplePolar() {
+	r, theta := cmplx.Polar(2i)
+	fmt.Printf("r: %.1f, θ: %.1f*π", r, theta/math.Pi)
+	// Output: r: 2.0, θ: 0.5*π
+}
diff --git a/src/math/cmplx/exp.go b/src/math/cmplx/exp.go
new file mode 100644
index 0000000..d5d0a5d
--- /dev/null
+++ b/src/math/cmplx/exp.go
@@ -0,0 +1,72 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// The original C code, the long comment, and the constants
+// below are from http://netlib.sandia.gov/cephes/c9x-complex/clog.c.
+// The go code is a simplified version of the original C.
+//
+// Cephes Math Library Release 2.8:  June, 2000
+// Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
+//
+// The readme file at http://netlib.sandia.gov/cephes/ says:
+//    Some software in this archive may be from the book _Methods and
+// Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
+// International, 1989) or from the Cephes Mathematical Library, a
+// commercial product. In either event, it is copyrighted by the author.
+// What you see here may be used freely but it comes with no support or
+// guarantee.
+//
+//   The two known misprints in the book are repaired here in the
+// source listings for the gamma function and the incomplete beta
+// integral.
+//
+//   Stephen L. Moshier
+//   moshier@na-net.ornl.gov
+
+// Complex exponential function
+//
+// DESCRIPTION:
+//
+// Returns the complex exponential of the complex argument z.
+//
+// If
+//     z = x + iy,
+//     r = exp(x),
+// then
+//     w = r cos y + i r sin y.
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      8700       3.7e-17     1.1e-17
+//    IEEE      -10,+10     30000       3.0e-16     8.7e-17
+
+// Exp returns e**x, the base-e exponential of x.
+func Exp(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case math.IsInf(re, 0):
+		switch {
+		case re > 0 && im == 0:
+			return x
+		case math.IsInf(im, 0) || math.IsNaN(im):
+			if re < 0 {
+				return complex(0, math.Copysign(0, im))
+			} else {
+				return complex(math.Inf(1.0), math.NaN())
+			}
+		}
+	case math.IsNaN(re):
+		if im == 0 {
+			return complex(math.NaN(), im)
+		}
+	}
+	r := math.Exp(real(x))
+	s, c := math.Sincos(imag(x))
+	return complex(r*c, r*s)
+}
diff --git a/src/math/cmplx/huge_test.go b/src/math/cmplx/huge_test.go
new file mode 100644
index 0000000..f8e60c2
--- /dev/null
+++ b/src/math/cmplx/huge_test.go
@@ -0,0 +1,22 @@
+// Copyright 2020 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+// Disabled for s390x because it uses assembly routines that are not
+// accurate for huge arguments.
+
+// +build !s390x
+
+package cmplx
+
+import (
+	"testing"
+)
+
+func TestTanHuge(t *testing.T) {
+	for i, x := range hugeIn {
+		if f := Tan(x); !cSoclose(tanHuge[i], f, 3e-15) {
+			t.Errorf("Tan(%g) = %g, want %g", x, f, tanHuge[i])
+		}
+	}
+}
diff --git a/src/math/cmplx/isinf.go b/src/math/cmplx/isinf.go
new file mode 100644
index 0000000..6273cd3
--- /dev/null
+++ b/src/math/cmplx/isinf.go
@@ -0,0 +1,21 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// IsInf reports whether either real(x) or imag(x) is an infinity.
+func IsInf(x complex128) bool {
+	if math.IsInf(real(x), 0) || math.IsInf(imag(x), 0) {
+		return true
+	}
+	return false
+}
+
+// Inf returns a complex infinity, complex(+Inf, +Inf).
+func Inf() complex128 {
+	inf := math.Inf(1)
+	return complex(inf, inf)
+}
diff --git a/src/math/cmplx/isnan.go b/src/math/cmplx/isnan.go
new file mode 100644
index 0000000..d3382c0
--- /dev/null
+++ b/src/math/cmplx/isnan.go
@@ -0,0 +1,25 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// IsNaN reports whether either real(x) or imag(x) is NaN
+// and neither is an infinity.
+func IsNaN(x complex128) bool {
+	switch {
+	case math.IsInf(real(x), 0) || math.IsInf(imag(x), 0):
+		return false
+	case math.IsNaN(real(x)) || math.IsNaN(imag(x)):
+		return true
+	}
+	return false
+}
+
+// NaN returns a complex ``not-a-number'' value.
+func NaN() complex128 {
+	nan := math.NaN()
+	return complex(nan, nan)
+}
diff --git a/src/math/cmplx/log.go b/src/math/cmplx/log.go
new file mode 100644
index 0000000..fd39c76
--- /dev/null
+++ b/src/math/cmplx/log.go
@@ -0,0 +1,65 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// The original C code, the long comment, and the constants
+// below are from http://netlib.sandia.gov/cephes/c9x-complex/clog.c.
+// The go code is a simplified version of the original C.
+//
+// Cephes Math Library Release 2.8:  June, 2000
+// Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
+//
+// The readme file at http://netlib.sandia.gov/cephes/ says:
+//    Some software in this archive may be from the book _Methods and
+// Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
+// International, 1989) or from the Cephes Mathematical Library, a
+// commercial product. In either event, it is copyrighted by the author.
+// What you see here may be used freely but it comes with no support or
+// guarantee.
+//
+//   The two known misprints in the book are repaired here in the
+// source listings for the gamma function and the incomplete beta
+// integral.
+//
+//   Stephen L. Moshier
+//   moshier@na-net.ornl.gov
+
+// Complex natural logarithm
+//
+// DESCRIPTION:
+//
+// Returns complex logarithm to the base e (2.718...) of
+// the complex argument z.
+//
+// If
+//       z = x + iy, r = sqrt( x**2 + y**2 ),
+// then
+//       w = log(r) + i arctan(y/x).
+//
+// The arctangent ranges from -PI to +PI.
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      7000       8.5e-17     1.9e-17
+//    IEEE      -10,+10     30000       5.0e-15     1.1e-16
+//
+// Larger relative error can be observed for z near 1 +i0.
+// In IEEE arithmetic the peak absolute error is 5.2e-16, rms
+// absolute error 1.0e-16.
+
+// Log returns the natural logarithm of x.
+func Log(x complex128) complex128 {
+	return complex(math.Log(Abs(x)), Phase(x))
+}
+
+// Log10 returns the decimal logarithm of x.
+func Log10(x complex128) complex128 {
+	z := Log(x)
+	return complex(math.Log10E*real(z), math.Log10E*imag(z))
+}
diff --git a/src/math/cmplx/phase.go b/src/math/cmplx/phase.go
new file mode 100644
index 0000000..03cece8
--- /dev/null
+++ b/src/math/cmplx/phase.go
@@ -0,0 +1,11 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// Phase returns the phase (also called the argument) of x.
+// The returned value is in the range [-Pi, Pi].
+func Phase(x complex128) float64 { return math.Atan2(imag(x), real(x)) }
diff --git a/src/math/cmplx/polar.go b/src/math/cmplx/polar.go
new file mode 100644
index 0000000..9b192bc
--- /dev/null
+++ b/src/math/cmplx/polar.go
@@ -0,0 +1,12 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+// Polar returns the absolute value r and phase θ of x,
+// such that x = r * e**θi.
+// The phase is in the range [-Pi, Pi].
+func Polar(x complex128) (r, θ float64) {
+	return Abs(x), Phase(x)
+}
diff --git a/src/math/cmplx/pow.go b/src/math/cmplx/pow.go
new file mode 100644
index 0000000..5a405f8
--- /dev/null
+++ b/src/math/cmplx/pow.go
@@ -0,0 +1,81 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// The original C code, the long comment, and the constants
+// below are from http://netlib.sandia.gov/cephes/c9x-complex/clog.c.
+// The go code is a simplified version of the original C.
+//
+// Cephes Math Library Release 2.8:  June, 2000
+// Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
+//
+// The readme file at http://netlib.sandia.gov/cephes/ says:
+//    Some software in this archive may be from the book _Methods and
+// Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
+// International, 1989) or from the Cephes Mathematical Library, a
+// commercial product. In either event, it is copyrighted by the author.
+// What you see here may be used freely but it comes with no support or
+// guarantee.
+//
+//   The two known misprints in the book are repaired here in the
+// source listings for the gamma function and the incomplete beta
+// integral.
+//
+//   Stephen L. Moshier
+//   moshier@na-net.ornl.gov
+
+// Complex power function
+//
+// DESCRIPTION:
+//
+// Raises complex A to the complex Zth power.
+// Definition is per AMS55 # 4.2.8,
+// analytically equivalent to cpow(a,z) = cexp(z clog(a)).
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    IEEE      -10,+10     30000       9.4e-15     1.5e-15
+
+// Pow returns x**y, the base-x exponential of y.
+// For generalized compatibility with math.Pow:
+//	Pow(0, ±0) returns 1+0i
+//	Pow(0, c) for real(c)<0 returns Inf+0i if imag(c) is zero, otherwise Inf+Inf i.
+func Pow(x, y complex128) complex128 {
+	if x == 0 { // Guaranteed also true for x == -0.
+		if IsNaN(y) {
+			return NaN()
+		}
+		r, i := real(y), imag(y)
+		switch {
+		case r == 0:
+			return 1
+		case r < 0:
+			if i == 0 {
+				return complex(math.Inf(1), 0)
+			}
+			return Inf()
+		case r > 0:
+			return 0
+		}
+		panic("not reached")
+	}
+	modulus := Abs(x)
+	if modulus == 0 {
+		return complex(0, 0)
+	}
+	r := math.Pow(modulus, real(y))
+	arg := Phase(x)
+	theta := real(y) * arg
+	if imag(y) != 0 {
+		r *= math.Exp(-imag(y) * arg)
+		theta += imag(y) * math.Log(modulus)
+	}
+	s, c := math.Sincos(theta)
+	return complex(r*c, r*s)
+}
diff --git a/src/math/cmplx/rect.go b/src/math/cmplx/rect.go
new file mode 100644
index 0000000..bf94d78
--- /dev/null
+++ b/src/math/cmplx/rect.go
@@ -0,0 +1,13 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// Rect returns the complex number x with polar coordinates r, θ.
+func Rect(r, θ float64) complex128 {
+	s, c := math.Sincos(θ)
+	return complex(r*c, r*s)
+}
diff --git a/src/math/cmplx/sin.go b/src/math/cmplx/sin.go
new file mode 100644
index 0000000..febac0e
--- /dev/null
+++ b/src/math/cmplx/sin.go
@@ -0,0 +1,184 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// The original C code, the long comment, and the constants
+// below are from http://netlib.sandia.gov/cephes/c9x-complex/clog.c.
+// The go code is a simplified version of the original C.
+//
+// Cephes Math Library Release 2.8:  June, 2000
+// Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
+//
+// The readme file at http://netlib.sandia.gov/cephes/ says:
+//    Some software in this archive may be from the book _Methods and
+// Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
+// International, 1989) or from the Cephes Mathematical Library, a
+// commercial product. In either event, it is copyrighted by the author.
+// What you see here may be used freely but it comes with no support or
+// guarantee.
+//
+//   The two known misprints in the book are repaired here in the
+// source listings for the gamma function and the incomplete beta
+// integral.
+//
+//   Stephen L. Moshier
+//   moshier@na-net.ornl.gov
+
+// Complex circular sine
+//
+// DESCRIPTION:
+//
+// If
+//     z = x + iy,
+//
+// then
+//
+//     w = sin x  cosh y  +  i cos x sinh y.
+//
+// csin(z) = -i csinh(iz).
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      8400       5.3e-17     1.3e-17
+//    IEEE      -10,+10     30000       3.8e-16     1.0e-16
+// Also tested by csin(casin(z)) = z.
+
+// Sin returns the sine of x.
+func Sin(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case im == 0 && (math.IsInf(re, 0) || math.IsNaN(re)):
+		return complex(math.NaN(), im)
+	case math.IsInf(im, 0):
+		switch {
+		case re == 0:
+			return x
+		case math.IsInf(re, 0) || math.IsNaN(re):
+			return complex(math.NaN(), im)
+		}
+	case re == 0 && math.IsNaN(im):
+		return x
+	}
+	s, c := math.Sincos(real(x))
+	sh, ch := sinhcosh(imag(x))
+	return complex(s*ch, c*sh)
+}
+
+// Complex hyperbolic sine
+//
+// DESCRIPTION:
+//
+// csinh z = (cexp(z) - cexp(-z))/2
+//         = sinh x * cos y  +  i cosh x * sin y .
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    IEEE      -10,+10     30000       3.1e-16     8.2e-17
+
+// Sinh returns the hyperbolic sine of x.
+func Sinh(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case re == 0 && (math.IsInf(im, 0) || math.IsNaN(im)):
+		return complex(re, math.NaN())
+	case math.IsInf(re, 0):
+		switch {
+		case im == 0:
+			return complex(re, im)
+		case math.IsInf(im, 0) || math.IsNaN(im):
+			return complex(re, math.NaN())
+		}
+	case im == 0 && math.IsNaN(re):
+		return complex(math.NaN(), im)
+	}
+	s, c := math.Sincos(imag(x))
+	sh, ch := sinhcosh(real(x))
+	return complex(c*sh, s*ch)
+}
+
+// Complex circular cosine
+//
+// DESCRIPTION:
+//
+// If
+//     z = x + iy,
+//
+// then
+//
+//     w = cos x  cosh y  -  i sin x sinh y.
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      8400       4.5e-17     1.3e-17
+//    IEEE      -10,+10     30000       3.8e-16     1.0e-16
+
+// Cos returns the cosine of x.
+func Cos(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case im == 0 && (math.IsInf(re, 0) || math.IsNaN(re)):
+		return complex(math.NaN(), -im*math.Copysign(0, re))
+	case math.IsInf(im, 0):
+		switch {
+		case re == 0:
+			return complex(math.Inf(1), -re*math.Copysign(0, im))
+		case math.IsInf(re, 0) || math.IsNaN(re):
+			return complex(math.Inf(1), math.NaN())
+		}
+	case re == 0 && math.IsNaN(im):
+		return complex(math.NaN(), 0)
+	}
+	s, c := math.Sincos(real(x))
+	sh, ch := sinhcosh(imag(x))
+	return complex(c*ch, -s*sh)
+}
+
+// Complex hyperbolic cosine
+//
+// DESCRIPTION:
+//
+// ccosh(z) = cosh x  cos y + i sinh x sin y .
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    IEEE      -10,+10     30000       2.9e-16     8.1e-17
+
+// Cosh returns the hyperbolic cosine of x.
+func Cosh(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case re == 0 && (math.IsInf(im, 0) || math.IsNaN(im)):
+		return complex(math.NaN(), re*math.Copysign(0, im))
+	case math.IsInf(re, 0):
+		switch {
+		case im == 0:
+			return complex(math.Inf(1), im*math.Copysign(0, re))
+		case math.IsInf(im, 0) || math.IsNaN(im):
+			return complex(math.Inf(1), math.NaN())
+		}
+	case im == 0 && math.IsNaN(re):
+		return complex(math.NaN(), im)
+	}
+	s, c := math.Sincos(imag(x))
+	sh, ch := sinhcosh(real(x))
+	return complex(c*ch, s*sh)
+}
+
+// calculate sinh and cosh
+func sinhcosh(x float64) (sh, ch float64) {
+	if math.Abs(x) <= 0.5 {
+		return math.Sinh(x), math.Cosh(x)
+	}
+	e := math.Exp(x)
+	ei := 0.5 / e
+	e *= 0.5
+	return e - ei, e + ei
+}
diff --git a/src/math/cmplx/sqrt.go b/src/math/cmplx/sqrt.go
new file mode 100644
index 0000000..eddce2f
--- /dev/null
+++ b/src/math/cmplx/sqrt.go
@@ -0,0 +1,107 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import "math"
+
+// The original C code, the long comment, and the constants
+// below are from http://netlib.sandia.gov/cephes/c9x-complex/clog.c.
+// The go code is a simplified version of the original C.
+//
+// Cephes Math Library Release 2.8:  June, 2000
+// Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
+//
+// The readme file at http://netlib.sandia.gov/cephes/ says:
+//    Some software in this archive may be from the book _Methods and
+// Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
+// International, 1989) or from the Cephes Mathematical Library, a
+// commercial product. In either event, it is copyrighted by the author.
+// What you see here may be used freely but it comes with no support or
+// guarantee.
+//
+//   The two known misprints in the book are repaired here in the
+// source listings for the gamma function and the incomplete beta
+// integral.
+//
+//   Stephen L. Moshier
+//   moshier@na-net.ornl.gov
+
+// Complex square root
+//
+// DESCRIPTION:
+//
+// If z = x + iy,  r = |z|, then
+//
+//                       1/2
+// Re w  =  [ (r + x)/2 ]   ,
+//
+//                       1/2
+// Im w  =  [ (r - x)/2 ]   .
+//
+// Cancellation error in r-x or r+x is avoided by using the
+// identity  2 Re w Im w  =  y.
+//
+// Note that -w is also a square root of z. The root chosen
+// is always in the right half plane and Im w has the same sign as y.
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10     25000       3.2e-17     9.6e-18
+//    IEEE      -10,+10   1,000,000     2.9e-16     6.1e-17
+
+// Sqrt returns the square root of x.
+// The result r is chosen so that real(r) ≥ 0 and imag(r) has the same sign as imag(x).
+func Sqrt(x complex128) complex128 {
+	if imag(x) == 0 {
+		// Ensure that imag(r) has the same sign as imag(x) for imag(x) == signed zero.
+		if real(x) == 0 {
+			return complex(0, imag(x))
+		}
+		if real(x) < 0 {
+			return complex(0, math.Copysign(math.Sqrt(-real(x)), imag(x)))
+		}
+		return complex(math.Sqrt(real(x)), imag(x))
+	} else if math.IsInf(imag(x), 0) {
+		return complex(math.Inf(1.0), imag(x))
+	}
+	if real(x) == 0 {
+		if imag(x) < 0 {
+			r := math.Sqrt(-0.5 * imag(x))
+			return complex(r, -r)
+		}
+		r := math.Sqrt(0.5 * imag(x))
+		return complex(r, r)
+	}
+	a := real(x)
+	b := imag(x)
+	var scale float64
+	// Rescale to avoid internal overflow or underflow.
+	if math.Abs(a) > 4 || math.Abs(b) > 4 {
+		a *= 0.25
+		b *= 0.25
+		scale = 2
+	} else {
+		a *= 1.8014398509481984e16 // 2**54
+		b *= 1.8014398509481984e16
+		scale = 7.450580596923828125e-9 // 2**-27
+	}
+	r := math.Hypot(a, b)
+	var t float64
+	if a > 0 {
+		t = math.Sqrt(0.5*r + 0.5*a)
+		r = scale * math.Abs((0.5*b)/t)
+		t *= scale
+	} else {
+		r = math.Sqrt(0.5*r - 0.5*a)
+		t = scale * math.Abs((0.5*b)/r)
+		r *= scale
+	}
+	if b < 0 {
+		return complex(t, -r)
+	}
+	return complex(t, r)
+}
diff --git a/src/math/cmplx/tan.go b/src/math/cmplx/tan.go
new file mode 100644
index 0000000..67a1133
--- /dev/null
+++ b/src/math/cmplx/tan.go
@@ -0,0 +1,297 @@
+// Copyright 2010 The Go Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package cmplx
+
+import (
+	"math"
+	"math/bits"
+)
+
+// The original C code, the long comment, and the constants
+// below are from http://netlib.sandia.gov/cephes/c9x-complex/clog.c.
+// The go code is a simplified version of the original C.
+//
+// Cephes Math Library Release 2.8:  June, 2000
+// Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
+//
+// The readme file at http://netlib.sandia.gov/cephes/ says:
+//    Some software in this archive may be from the book _Methods and
+// Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
+// International, 1989) or from the Cephes Mathematical Library, a
+// commercial product. In either event, it is copyrighted by the author.
+// What you see here may be used freely but it comes with no support or
+// guarantee.
+//
+//   The two known misprints in the book are repaired here in the
+// source listings for the gamma function and the incomplete beta
+// integral.
+//
+//   Stephen L. Moshier
+//   moshier@na-net.ornl.gov
+
+// Complex circular tangent
+//
+// DESCRIPTION:
+//
+// If
+//     z = x + iy,
+//
+// then
+//
+//           sin 2x  +  i sinh 2y
+//     w  =  --------------------.
+//            cos 2x  +  cosh 2y
+//
+// On the real axis the denominator is zero at odd multiples
+// of PI/2. The denominator is evaluated by its Taylor
+// series near these points.
+//
+// ctan(z) = -i ctanh(iz).
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      5200       7.1e-17     1.6e-17
+//    IEEE      -10,+10     30000       7.2e-16     1.2e-16
+// Also tested by ctan * ccot = 1 and catan(ctan(z))  =  z.
+
+// Tan returns the tangent of x.
+func Tan(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case math.IsInf(im, 0):
+		switch {
+		case math.IsInf(re, 0) || math.IsNaN(re):
+			return complex(math.Copysign(0, re), math.Copysign(1, im))
+		}
+		return complex(math.Copysign(0, math.Sin(2*re)), math.Copysign(1, im))
+	case re == 0 && math.IsNaN(im):
+		return x
+	}
+	d := math.Cos(2*real(x)) + math.Cosh(2*imag(x))
+	if math.Abs(d) < 0.25 {
+		d = tanSeries(x)
+	}
+	if d == 0 {
+		return Inf()
+	}
+	return complex(math.Sin(2*real(x))/d, math.Sinh(2*imag(x))/d)
+}
+
+// Complex hyperbolic tangent
+//
+// DESCRIPTION:
+//
+// tanh z = (sinh 2x  +  i sin 2y) / (cosh 2x + cos 2y) .
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    IEEE      -10,+10     30000       1.7e-14     2.4e-16
+
+// Tanh returns the hyperbolic tangent of x.
+func Tanh(x complex128) complex128 {
+	switch re, im := real(x), imag(x); {
+	case math.IsInf(re, 0):
+		switch {
+		case math.IsInf(im, 0) || math.IsNaN(im):
+			return complex(math.Copysign(1, re), math.Copysign(0, im))
+		}
+		return complex(math.Copysign(1, re), math.Copysign(0, math.Sin(2*im)))
+	case im == 0 && math.IsNaN(re):
+		return x
+	}
+	d := math.Cosh(2*real(x)) + math.Cos(2*imag(x))
+	if d == 0 {
+		return Inf()
+	}
+	return complex(math.Sinh(2*real(x))/d, math.Sin(2*imag(x))/d)
+}
+
+// reducePi reduces the input argument x to the range (-Pi/2, Pi/2].
+// x must be greater than or equal to 0. For small arguments it
+// uses Cody-Waite reduction in 3 float64 parts based on:
+// "Elementary Function Evaluation:  Algorithms and Implementation"
+// Jean-Michel Muller, 1997.
+// For very large arguments it uses Payne-Hanek range reduction based on:
+// "ARGUMENT REDUCTION FOR HUGE ARGUMENTS: Good to the Last Bit"
+// K. C. Ng et al, March 24, 1992.
+func reducePi(x float64) float64 {
+	// reduceThreshold is the maximum value of x where the reduction using
+	// Cody-Waite reduction still gives accurate results. This threshold
+	// is set by t*PIn being representable as a float64 without error
+	// where t is given by t = floor(x * (1 / Pi)) and PIn are the leading partial
+	// terms of Pi. Since the leading terms, PI1 and PI2 below, have 30 and 32
+	// trailing zero bits respectively, t should have less than 30 significant bits.
+	//	t < 1<<30  -> floor(x*(1/Pi)+0.5) < 1<<30 -> x < (1<<30-1) * Pi - 0.5
+	// So, conservatively we can take x < 1<<30.
+	const reduceThreshold float64 = 1 << 30
+	if math.Abs(x) < reduceThreshold {
+		// Use Cody-Waite reduction in three parts.
+		const (
+			// PI1, PI2 and PI3 comprise an extended precision value of PI
+			// such that PI ~= PI1 + PI2 + PI3. The parts are chosen so
+			// that PI1 and PI2 have an approximately equal number of trailing
+			// zero bits. This ensures that t*PI1 and t*PI2 are exact for
+			// large integer values of t. The full precision PI3 ensures the
+			// approximation of PI is accurate to 102 bits to handle cancellation
+			// during subtraction.
+			PI1 = 3.141592502593994      // 0x400921fb40000000
+			PI2 = 1.5099578831723193e-07 // 0x3e84442d00000000
+			PI3 = 1.0780605716316238e-14 // 0x3d08469898cc5170
+		)
+		t := x / math.Pi
+		t += 0.5
+		t = float64(int64(t)) // int64(t) = the multiple
+		return ((x - t*PI1) - t*PI2) - t*PI3
+	}
+	// Must apply Payne-Hanek range reduction
+	const (
+		mask     = 0x7FF
+		shift    = 64 - 11 - 1
+		bias     = 1023
+		fracMask = 1<<shift - 1
+	)
+	// Extract out the integer and exponent such that,
+	// x = ix * 2 ** exp.
+	ix := math.Float64bits(x)
+	exp := int(ix>>shift&mask) - bias - shift
+	ix &= fracMask
+	ix |= 1 << shift
+
+	// mPi is the binary digits of 1/Pi as a uint64 array,
+	// that is, 1/Pi = Sum mPi[i]*2^(-64*i).
+	// 19 64-bit digits give 1216 bits of precision
+	// to handle the largest possible float64 exponent.
+	var mPi = [...]uint64{
+		0x0000000000000000,
+		0x517cc1b727220a94,
+		0xfe13abe8fa9a6ee0,
+		0x6db14acc9e21c820,
+		0xff28b1d5ef5de2b0,
+		0xdb92371d2126e970,
+		0x0324977504e8c90e,
+		0x7f0ef58e5894d39f,
+		0x74411afa975da242,
+		0x74ce38135a2fbf20,
+		0x9cc8eb1cc1a99cfa,
+		0x4e422fc5defc941d,
+		0x8ffc4bffef02cc07,
+		0xf79788c5ad05368f,
+		0xb69b3f6793e584db,
+		0xa7a31fb34f2ff516,
+		0xba93dd63f5f2f8bd,
+		0x9e839cfbc5294975,
+		0x35fdafd88fc6ae84,
+		0x2b0198237e3db5d5,
+	}
+	// Use the exponent to extract the 3 appropriate uint64 digits from mPi,
+	// B ~ (z0, z1, z2), such that the product leading digit has the exponent -64.
+	// Note, exp >= 50 since x >= reduceThreshold and exp < 971 for maximum float64.
+	digit, bitshift := uint(exp+64)/64, uint(exp+64)%64
+	z0 := (mPi[digit] << bitshift) | (mPi[digit+1] >> (64 - bitshift))
+	z1 := (mPi[digit+1] << bitshift) | (mPi[digit+2] >> (64 - bitshift))
+	z2 := (mPi[digit+2] << bitshift) | (mPi[digit+3] >> (64 - bitshift))
+	// Multiply mantissa by the digits and extract the upper two digits (hi, lo).
+	z2hi, _ := bits.Mul64(z2, ix)
+	z1hi, z1lo := bits.Mul64(z1, ix)
+	z0lo := z0 * ix
+	lo, c := bits.Add64(z1lo, z2hi, 0)
+	hi, _ := bits.Add64(z0lo, z1hi, c)
+	// Find the magnitude of the fraction.
+	lz := uint(bits.LeadingZeros64(hi))
+	e := uint64(bias - (lz + 1))
+	// Clear implicit mantissa bit and shift into place.
+	hi = (hi << (lz + 1)) | (lo >> (64 - (lz + 1)))
+	hi >>= 64 - shift
+	// Include the exponent and convert to a float.
+	hi |= e << shift
+	x = math.Float64frombits(hi)
+	// map to (-Pi/2, Pi/2]
+	if x > 0.5 {
+		x--
+	}
+	return math.Pi * x
+}
+
+// Taylor series expansion for cosh(2y) - cos(2x)
+func tanSeries(z complex128) float64 {
+	const MACHEP = 1.0 / (1 << 53)
+	x := math.Abs(2 * real(z))
+	y := math.Abs(2 * imag(z))
+	x = reducePi(x)
+	x = x * x
+	y = y * y
+	x2 := 1.0
+	y2 := 1.0
+	f := 1.0
+	rn := 0.0
+	d := 0.0
+	for {
+		rn++
+		f *= rn
+		rn++
+		f *= rn
+		x2 *= x
+		y2 *= y
+		t := y2 + x2
+		t /= f
+		d += t
+
+		rn++
+		f *= rn
+		rn++
+		f *= rn
+		x2 *= x
+		y2 *= y
+		t = y2 - x2
+		t /= f
+		d += t
+		if !(math.Abs(t/d) > MACHEP) {
+			// Caution: Use ! and > instead of <= for correct behavior if t/d is NaN.
+			// See issue 17577.
+			break
+		}
+	}
+	return d
+}
+
+// Complex circular cotangent
+//
+// DESCRIPTION:
+//
+// If
+//     z = x + iy,
+//
+// then
+//
+//           sin 2x  -  i sinh 2y
+//     w  =  --------------------.
+//            cosh 2y  -  cos 2x
+//
+// On the real axis, the denominator has zeros at even
+// multiples of PI/2.  Near these points it is evaluated
+// by a Taylor series.
+//
+// ACCURACY:
+//
+//                      Relative error:
+// arithmetic   domain     # trials      peak         rms
+//    DEC       -10,+10      3000       6.5e-17     1.6e-17
+//    IEEE      -10,+10     30000       9.2e-16     1.2e-16
+// Also tested by ctan * ccot = 1 + i0.
+
+// Cot returns the cotangent of x.
+func Cot(x complex128) complex128 {
+	d := math.Cosh(2*imag(x)) - math.Cos(2*real(x))
+	if math.Abs(d) < 0.25 {
+		d = tanSeries(x)
+	}
+	if d == 0 {
+		return Inf()
+	}
+	return complex(math.Sin(2*real(x))/d, -math.Sinh(2*imag(x))/d)
+}