From 7589568a8a575b740cecd38325e2c5a546369993 Mon Sep 17 00:00:00 2001 From: diluculo Date: Tue, 16 Jan 2018 14:45:23 +0900 Subject: [PATCH 1/6] Trig: updated Tan(z), Cot(z), Tanh(z), Coth(z), Csch(z), and Sech(z) based on the Kahan algorithm. --- src/Numerics/Trigonometry.cs | 183 +++++++++++++++++++++++++----- src/UnitTests/TrigonometryTest.cs | 48 ++++++++ 2 files changed, 203 insertions(+), 28 deletions(-) diff --git a/src/Numerics/Trigonometry.cs b/src/Numerics/Trigonometry.cs index 70c9534e..7fe39db9 100644 --- a/src/Numerics/Trigonometry.cs +++ b/src/Numerics/Trigonometry.cs @@ -192,11 +192,39 @@ namespace MathNet.Numerics return new Complex(Tan(value.Real), 0.0); } - var cosr = Cos(value.Real); - var sinhi = Sinh(value.Imaginary); - var denom = (cosr * cosr) + (sinhi * sinhi); + // tan(x + j*y) = (tan(x) + j*cosh(y)*sinh(y)/cos^2(x))/(1 + sinh^2(y)/cos^2(x)) + // if |y| > asinh(sqrt(max)), tan(z) = 4*cos(x)*sin(x)*exp(-2*|y|) + j*sign(y) + // if exp(-|y|) = 0, tan(z) = j*sign(y) + // if tan(x) = +/- oo or 1/cos^2(x) = 1 + tan^2(x) = oo, tan(z) = j*cosh(y)/sinh(y) + // + // The algorithm is from: + // Kahan, W. Branch Cuts for Complex Elementary Functions, or Much Ado + // About Nothing's Sign Bit." In The State of the Art in Numerical Analysis: + // Proceedings of the Joint IMA / SIAM Conference on the State of the Art in + // Numerical Analysis Held at the UN (Ed.A.Iserles and M.J.D.Powell). + // New York: Clarendon Press, pp. 165 - 211, 1987. + + double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 + + if (Math.Abs(value.Imaginary) > upperLimit) + { + double e = Math.Exp(-Math.Abs(value.Imaginary)); + return e == 0.0 + ? new Complex(0.0, Math.Sign(value.Imaginary)) + : new Complex(4.0 * Math.Cos(value.Real) * Math.Sin(value.Real) * e * e, Math.Sign(value.Imaginary)); + } + + double tanr = Math.Tan(value.Real); + double beta = 1 + tanr * tanr; // beta = 1/cos(x)^2 = 1 + t^2 + double sinhi = Math.Sinh(value.Imaginary); + double coshi = Math.Cosh(value.Imaginary); + + if (double.IsInfinity(tanr)) + return new Complex(0.0, coshi / sinhi); - return new Complex(Sin(value.Real) * cosr / denom, sinhi * Cosh(value.Imaginary) / denom); + double denom = 1.0 + beta * sinhi * sinhi; + + return new Complex(tanr / denom, beta * coshi * sinhi / denom); } /// @@ -221,11 +249,34 @@ namespace MathNet.Numerics return new Complex(Cot(value.Real), 0d); } - var sinr = Sin(value.Real); - var sinhi = Sinh(value.Imaginary); - var denom = (sinr * sinr) + (sinhi * sinhi); + // cot(x + j*y) = (cot(x) - j*cosh(y)*sinh(y)/sin^2(x))/(1 + sinh^2(y)/sin^2(x)) + // if |y| > asinh(sqrt(max)), cot(z) = 4*cos(x)*sin(x)*exp(-2*|y|) - j*sign(y) + // if exp(-|y|) = 0, cot(z) = -j*sign(y) + // if cot(x) = +/- oo or 1/sin^2(x) = 1 + cot^2(x) = oo, cot(z) = -j*cosh(y)/sinh(y) + // + // The algorithm is based on Kahan. + + double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 + + if (Math.Abs(value.Imaginary) > upperLimit) + { + double e = Math.Exp(-Math.Abs(value.Imaginary)); + return e == 0.0 + ? new Complex(0.0, -Math.Sign(value.Imaginary)) + : new Complex(4.0 * Math.Cos(value.Real) * Math.Sin(value.Real) * e * e, -Math.Sign(value.Imaginary)); + } + + double cotr = Cot(value.Real); + double beta = 1 + cotr * cotr; // beta = 1/cos(x)^2 = 1 + t^2 + double sinhi = Sinh(value.Imaginary); + double coshi = Cosh(value.Imaginary); + + if (double.IsInfinity(cotr)) + return new Complex(0.0, - coshi / sinhi); - return new Complex(sinr * Cos(value.Real) / denom, -sinhi * Cosh(value.Imaginary) / denom); + double denom = 1.0 + beta * sinhi * sinhi; + + return new Complex(cotr / denom, - beta * coshi * sinhi / denom); } /// @@ -515,17 +566,33 @@ namespace MathNet.Numerics return new Complex(Tanh(value.Real), 0.0); } - var cosi = Cos(value.Imaginary); - var sinhr = Sinh(value.Real); + // tanh(x + j*y) = (cosh(x)*sinh(x)/cos^2(y) + j*tan(y))/(1 + sinh^2(x)/cos^2(y)) + // if |x| > asinh(sqrt(max)), tanh(z) = sign(x) + j*4*cos(y)*sin(y)*exp(-2*|x|) + // if exp(-|x|) = 0, tanh(z) = sign(x) + // if tan(y) = +/- oo or 1/cos^2(y) = 1 + tan^2(y) = oo, tanh(z) = cosh(x)/sinh(x) + // + // The algorithm is based on Kahan. + + double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 - if (double.IsInfinity(sinhr)) + if (Math.Abs(value.Real) > upperLimit) { - return new Complex(double.IsPositiveInfinity(sinhr) ? 1 : -1, 0.0); + double e = Math.Exp(-Math.Abs(value.Real)); + return e == 0.0 + ? new Complex(Math.Sign(value.Real), 0.0) + : new Complex(Math.Sign(value.Real), 4.0 * Math.Cos(value.Imaginary) * Math.Sin(value.Imaginary) * e * e); } - var denom = (cosi * cosi) + (sinhr * sinhr); + double tani = Tan(value.Imaginary); + double beta = 1 + tani * tani; // beta = 1/cos^2(y) = 1 + t^2 + double sinhr = Sinh(value.Real); + double coshr = Cosh(value.Real); + + if (double.IsInfinity(tani)) + return new Complex(coshr / sinhr, 0.0); - return new Complex(Cosh(value.Real) * sinhr / denom, cosi * Sin(value.Imaginary) / denom); + double denom = 1.0 + beta * sinhr * sinhr; + return new Complex(beta * coshr * sinhr / denom, tani / denom); } /// @@ -562,17 +629,33 @@ namespace MathNet.Numerics return new Complex(Coth(value.Real), 0.0); } - var sini = Sin(value.Imaginary); - var sinhr = Sinh(value.Real); + // coth(x + j*y) = (cosh(x)*sinh(x)/sin^2(y) - j*cot(y))/(1 + sinh^2(x)/sin^2(y)) + // if |x| > asinh(sqrt(max)), coth(z) = sign(x) - j*4*cos(y)*sin(y)*exp(-2*|x|) + // if exp(-|x|) = 0, coth(z) = sign(x) + // if cot(y) = +/- oo or 1/sin^2(y) = 1 + cot^2(y) = oo, coth(z) = cosh(x)/sinh(x) + // + // The algorithm is based on Kahan. - if (double.IsInfinity(sinhr)) + double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 + + if (Math.Abs(value.Real) > upperLimit) { - return new Complex(double.IsPositiveInfinity(sinhr) ? 1 : -1, 0.0); + double e = Math.Exp(-Math.Abs(value.Real)); + return e == 0.0 + ? new Complex(Math.Sign(value.Real), 0.0) + : new Complex(Math.Sign(value.Real), -4.0 * Math.Cos(value.Imaginary) * Math.Sin(value.Imaginary) * e * e); } - var denom = (sini * sini) + (sinhr * sinhr); + double coti = Cot(value.Imaginary); + double beta = 1 + coti * coti; + double sinhr = Math.Sinh(value.Real); + double coshr = Math.Cosh(value.Real); + + if (double.IsInfinity(coti)) + return new Complex(coshr / sinhr, 0.0); - return new Complex(sinhr * Cosh(value.Real) / denom, -sini * Cos(value.Imaginary) / denom); + double denom = 1.0 + beta * sinhr * sinhr; + return new Complex(beta * coshr * sinhr / denom, -coti / denom); } /// @@ -597,14 +680,36 @@ namespace MathNet.Numerics return new Complex(Sech(value.Real), 0.0); } - var exp = value.Exp(); + // sech(x + j*y) = (cosh(x)/cos(y) - j*sinh(x)*tan(y)/cos(y))/(1 + sinh^2(x)/cos^2(y)) + // if |x| > asinh(sqrt(max)), sech(z) = 4*cosh(x)*cos(y)*exp(-2*|x|) - j*4*sinh(x)*tan(y)*cos(y)*exp(-2*|x|) + // if exp(-|x|) = 0, sech(z) = 0 + // if tan(y) = +/- oo or 1/cos^2(y) = 1 + tan^2(y) = oo, sech(z) = -j*sign(tan(y))/sinh(x) + // + // The algorithm is based on Kahan. + + double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 + + double tani = Tan(value.Imaginary); + double cosi = Cos(value.Imaginary); + double beta = 1.0 + tani * tani; + double sinhr = Math.Sinh(value.Real); + double coshr = Math.Cosh(value.Real); + + if (Math.Abs(value.Real) > upperLimit) + { + double e = Math.Exp(-Math.Abs(value.Real)); + return e == 0.0 + ? new Complex(0, 0) + : new Complex(4.0 * coshr * cosi * e * e, -4.0 * sinhr * tani * cosi * e * e); + } - if (exp.IsInfinity()) + if (double.IsInfinity(tani)) { - return Complex.Zero; + return new Complex(0.0, -Math.Sign(tani) / sinhr); } - return 2 * exp / (exp.Square() + 1); + double denom = 1.0 + beta * sinhr * sinhr; + return new Complex(coshr / cosi / denom, -sinhr * tani / cosi / denom); } /// @@ -629,14 +734,36 @@ namespace MathNet.Numerics return new Complex(Csch(value.Real), 0.0); } - var exp = value.Exp(); + // csch(x + j*y) = (sinh(x)*cot(y)/sin(y) - j*cosh(x)/sin(y))/(1 + sinh^2(x)/sin^2(y)) + // if |x| > asinh(sqrt(max)), csch(z) = 4*sinh(x)*cot(y)*sin(y)*exp(-2*|x|) - j*4*cosh(x)*sin(y)*exp(-2*|x|) + // if exp(-|x|) = 0, csch(z) = 0 + // if cot(y) = +/- oo or 1/sin^2(x) = 1 + cot^2(x) = oo, csch(z) = sign(cot(y))/sinh(x) + // + // The algorithm is based on Kahan. + + double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); + + double coti = Cot(value.Imaginary); + double sini = Sin(value.Imaginary); + double beta = 1 + coti * coti; + double sinhr = Sinh(value.Real); + double coshr = Cosh(value.Real); + + if (Math.Abs(value.Real) > upperLimit) + { + double e = Math.Exp(-Math.Abs(value.Real)); + return e == 0.0 + ? new Complex(0, 0) + : new Complex(4.0 * sinhr * coti * sini * e * e, -4.0 * coshr * sini * e * e); + } - if (exp.IsInfinity()) + if (double.IsInfinity(coti)) { - return Complex.Zero; + return new Complex(Math.Sign(coti) / sinhr, 0.0); } - return 2 * exp / (exp.Square() - 1); + double denom = 1.0 + beta * sinhr * sinhr; + return new Complex(sinhr * coti / sini / denom, -coshr / sini / denom); } diff --git a/src/UnitTests/TrigonometryTest.cs b/src/UnitTests/TrigonometryTest.cs index a428cda6..2fb3194b 100644 --- a/src/UnitTests/TrigonometryTest.cs +++ b/src/UnitTests/TrigonometryTest.cs @@ -58,6 +58,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(-1.19209289550780998537e-7, 0.0, 0.99999999999999289, 0.0)] [TestCase(8.388608e6, 1.19209289550780998537e-7, -0.90175467375876572, -5.1528001100635277e-8)] [TestCase(-1.19209289550780998537e-7, -8.388608e6, double.PositiveInfinity, double.NegativeInfinity)] + [TestCase(512.0, 32.0, -39356457675156.6, -3139507838262.74)] + [TestCase(-512.0, -32.0, -39356457675156.6, -3139507838262.74)] + [TestCase(32.0, -512.0, +9.52855589474963e+221, +6.29843301304523e+221)] + [TestCase(-32.0, -512.0, +9.52855589474963e+221, -6.29843301304523e+221)] public void CanComputeComplexCosine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Cos(); @@ -79,6 +83,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(-1.19209289550780998537e-7, 0.0, -1.19209289550780998537e-7, 0.0)] [TestCase(8.388608e6, 1.19209289550780998537e-7, 0.43224820225680083, -1.0749753400787824e-7)] [TestCase(-1.19209289550780998537e-7, -8.388608e6, double.NegativeInfinity, double.NegativeInfinity)] + [TestCase(512.0, 32.0, +3139507838262.74, -39356457675156.6)] + [TestCase(-512.0, -32.0, -3139507838262.74, +39356457675156.6)] + [TestCase(32.0, -512.0, +6.29843301304523e+221, -9.52855589474963e+221)] + [TestCase(-32.0, -512.0, -6.29843301304523e+221, -9.52855589474963e+221)] public void CanComputeComplexSine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sin(); @@ -100,6 +108,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(-1.19209289550780998537e-7, 0.0, -1.1920928955078157e-7, 0.0)] [TestCase(8.388608e6, 1.19209289550780998537e-7, -0.47934123862653449, 1.4659977233982276e-7)] [TestCase(-8.388608e6, -1.19209289550780998537e-7, 0.47934123862653449, -1.4659977233982276e-7)] + [TestCase(512.0, 32.0, -5.08515122860094E-29, 1)] + [TestCase(-512.0, -32.0, +5.08515122860094E-29, -1)] + [TestCase(32.0, 512.0, +3.52598650243739e-445, 1)] + [TestCase(32.0, -512.0, +3.52598650243739e-445, -1)] public void CanComputeComplexTangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Tan(); @@ -549,6 +561,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(-1.19209289550780998537e-7, 0.0, -8388607.999999978, 0.0)] [TestCase(8.388608e6, 1.19209289550780998537e-7, -2.0861964701080704, -6.3803383253713457e-7)] [TestCase(-8.388608e6, -1.19209289550780998537e-7, 2.0861964701080704, 6.3803383253713457e-7)] + [TestCase(512.0, 32.0, -5.08515122860094E-29, -1.0)] + [TestCase(-512.0, -32.0, 5.08515122860094E-29, +1.0)] + [TestCase(32.0, 512.0, +3.52598650243739e-445, -1.0)] + [TestCase(32.0, -512.0, +3.52598650243739e-445, +1.0)] public void CanComputeComplexCotangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Cot(); @@ -570,6 +586,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(-1.19209289550780998537e-7, 0.0, 1.0000000000000071, 0.0)] [TestCase(8.388608e6, 1.19209289550780998537e-7, -1.1089490624226177, 6.3367488045143761e-8)] [TestCase(-8.388608e6, -1.19209289550780998537e-7, -1.1089490624226177, 6.3367488045143761e-8)] + [TestCase(512.0, 32.0, -2.52481261731330570134e-14, +2.01407074478744036413e-15)] + [TestCase(-512.0, -32.0, -2.52481261731330570134e-14, +2.01407074478744036413e-15)] + [TestCase(32.0, 512.0, +7.30361056703505e-223, +4.82773070945482e-223)] + [TestCase(32.0, -512.0, +7.30361056703505e-223, -4.82773070945482e-223)] public void CanComputeComplexSecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sec(); @@ -591,6 +611,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(-1.19209289550780998537e-7, 0.0, -8388608.0000000376, 0.0)] [TestCase(8.388608e6, 1.19209289550780998537e-7, 2.3134856195557596, 5.7534999050657057e-7)] [TestCase(-8.388608e6, -1.19209289550780998537e-7, -2.3134856195557596, -5.7534999050657057e-7)] + [TestCase(512.0, 32.0, 2.01407074478744e-15, 2.52481261731331e-14)] + [TestCase(-512.0, -32.0, -2.01407074478744e-15, -2.52481261731331e-14)] + [TestCase(32.0, 512.0, +4.82773070945482e-223, -7.30361056703505e-223)] + [TestCase(32.0, -512.0, +4.82773070945482e-223, +7.30361056703505e-223)] public void CanComputeComplexCosecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Csc(); @@ -616,6 +640,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(0.5, -0.5, 0.45730415318424922, -0.54061268571315335)] [TestCase(-0.5, 0.5, -0.45730415318424922, 0.54061268571315335)] [TestCase(-0.5, -0.5, -0.45730415318424922, -0.54061268571315335)] + [TestCase(512.0, 32.0, 9.52855589474962752845e221, 6.29843301304522728211e221)] + [TestCase(-512.0, -32.0, -9.52855589474962752845e221, -6.29843301304522728211e221)] + [TestCase(-32.0, 512.0, +3.93564576751566e13, +3.13950783826274e12)] + [TestCase(32.0, -512.0, -3.93564576751566e13, -3.13950783826274e12)] public void CanComputeComplexHyperbolicSine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sinh(); @@ -641,6 +669,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(0.5, -0.5, 0.9895848833999199, -0.24982639750046154)] [TestCase(-0.5, 0.5, 0.9895848833999199, -0.24982639750046154)] [TestCase(-0.5, -0.5, 0.9895848833999199, 0.24982639750046154)] + [TestCase(512.0, 32.0, 9.52855589474962752845e221, 6.29843301304522728211e221)] + [TestCase(-512.0, -32.0, 9.52855589474962752845e221, 6.29843301304522728211e221)] + [TestCase(-32.0, 512.0, -3.93564576751566e13, -3.13950783826274e12)] + [TestCase(32.0, -512.0, -3.93564576751566e13, -3.13950783826274e12)] public void CanComputeComplexHyperbolicCosine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Cosh(); @@ -666,6 +698,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(0.5, -0.5, 0.56408314126749848, -0.40389645531602575)] [TestCase(-0.5, 0.5, -0.56408314126749848, 0.40389645531602575)] [TestCase(-0.5, -0.5, -0.56408314126749848, -0.40389645531602575)] + [TestCase(512.0, 32.0, 1.0, 0.0)] + [TestCase(-512.0, -32.0, -1.0, 0.0)] + [TestCase(-32.0, 512.0, -1.0, -5.08515122860093626173e-29)] + [TestCase(32.0, -512.0, +1.0, +5.08515122860093626173e-29)] public void CanComputeComplexHyperbolicTangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Tanh(); @@ -691,6 +727,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(0.5, -0.5, 1.1719451445243514, 0.8391395790248311)] [TestCase(-0.5, 0.5, -1.1719451445243514, -0.8391395790248311)] [TestCase(-0.5, -0.5, -1.1719451445243514, 0.8391395790248311)] + [TestCase(512.0, 32.0, 1.0, 0.0)] + [TestCase(-512.0, -32.0, -1.0, 0.0)] + [TestCase(-32.0, 512.0, -1.0, 5.08515122860093626173e-29)] + [TestCase(32.0, -512.0, +1.0, -5.08515122860093626173e-29)] public void CanComputeComplexHyperbolicCotangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Coth(); @@ -716,6 +756,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(0.5, -0.5, 0.94997886761549463, 0.23982763093808804)] [TestCase(-0.5, 0.5, 0.94997886761549463, 0.23982763093808804)] [TestCase(-0.5, -0.5, 0.94997886761549463, -0.23982763093808804)] + [TestCase(512.0, 32.0, 7.30361056703505051822e-223, -4.82773070945482080706e-223)] + [TestCase(-512.0, -32.0, 7.30361056703505051822e-223, -4.82773070945482080706e-223)] + [TestCase(-32.0, 512.0, -2.52481261731330570134e-14, +2.01407074478744036413e-15)] + [TestCase(-32.0, -512.0, -2.52481261731330570134e-14, -2.01407074478744036413e-15)] public void CanComputeComplexHyperbolicSecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sech(); @@ -741,6 +785,10 @@ namespace MathNet.Numerics.UnitTests [TestCase(0.5, -0.5, 0.91207426403881078, 1.0782296946540223)] [TestCase(-0.5, 0.5, -0.91207426403881078, -1.0782296946540223)] [TestCase(-0.5, -0.5, -0.91207426403881078, 1.0782296946540223)] + [TestCase(512.0, 32.0, 7.30361056703505051822e-223, -4.82773070945482080706e-223)] + [TestCase(-512.0, -32.0, -7.30361056703505051822e-223, +4.82773070945482080706e-223)] + [TestCase(-32.0, 512.0, +2.52481261731330570134e-14, -2.01407074478744036413e-15)] + [TestCase(32.0, -512.0, -2.52481261731330570134e-14, +2.01407074478744036413e-15)] public void CanComputeComplexHyperbolicCosecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Csch(); From 9a67b5d9fd612731c0003a933552296261e397d4 Mon Sep 17 00:00:00 2001 From: diluculo Date: Tue, 16 Jan 2018 20:24:51 +0900 Subject: [PATCH 2/6] To avoid overflow, widen the range of x to apply the approximation of sinh^2(x). Simplify Coth(z), Tan(z), and Cot(z) --- src/Numerics/Trigonometry.cs | 105 ++++-------------------------- src/UnitTests/TrigonometryTest.cs | 28 ++++++++ 2 files changed, 42 insertions(+), 91 deletions(-) diff --git a/src/Numerics/Trigonometry.cs b/src/Numerics/Trigonometry.cs index 7fe39db9..45b4b8fa 100644 --- a/src/Numerics/Trigonometry.cs +++ b/src/Numerics/Trigonometry.cs @@ -192,39 +192,10 @@ namespace MathNet.Numerics return new Complex(Tan(value.Real), 0.0); } - // tan(x + j*y) = (tan(x) + j*cosh(y)*sinh(y)/cos^2(x))/(1 + sinh^2(y)/cos^2(x)) - // if |y| > asinh(sqrt(max)), tan(z) = 4*cos(x)*sin(x)*exp(-2*|y|) + j*sign(y) - // if exp(-|y|) = 0, tan(z) = j*sign(y) - // if tan(x) = +/- oo or 1/cos^2(x) = 1 + tan^2(x) = oo, tan(z) = j*cosh(y)/sinh(y) - // - // The algorithm is from: - // Kahan, W. Branch Cuts for Complex Elementary Functions, or Much Ado - // About Nothing's Sign Bit." In The State of the Art in Numerical Analysis: - // Proceedings of the Joint IMA / SIAM Conference on the State of the Art in - // Numerical Analysis Held at the UN (Ed.A.Iserles and M.J.D.Powell). - // New York: Clarendon Press, pp. 165 - 211, 1987. - - double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 - - if (Math.Abs(value.Imaginary) > upperLimit) - { - double e = Math.Exp(-Math.Abs(value.Imaginary)); - return e == 0.0 - ? new Complex(0.0, Math.Sign(value.Imaginary)) - : new Complex(4.0 * Math.Cos(value.Real) * Math.Sin(value.Real) * e * e, Math.Sign(value.Imaginary)); - } - - double tanr = Math.Tan(value.Real); - double beta = 1 + tanr * tanr; // beta = 1/cos(x)^2 = 1 + t^2 - double sinhi = Math.Sinh(value.Imaginary); - double coshi = Math.Cosh(value.Imaginary); - - if (double.IsInfinity(tanr)) - return new Complex(0.0, coshi / sinhi); + // tan(z) = - j*tanh(j*z) - double denom = 1.0 + beta * sinhi * sinhi; - - return new Complex(tanr / denom, beta * coshi * sinhi / denom); + Complex z = Tanh(new Complex(-value.Imaginary, value.Real)); + return new Complex(z.Imaginary, -z.Real); } /// @@ -249,34 +220,10 @@ namespace MathNet.Numerics return new Complex(Cot(value.Real), 0d); } - // cot(x + j*y) = (cot(x) - j*cosh(y)*sinh(y)/sin^2(x))/(1 + sinh^2(y)/sin^2(x)) - // if |y| > asinh(sqrt(max)), cot(z) = 4*cos(x)*sin(x)*exp(-2*|y|) - j*sign(y) - // if exp(-|y|) = 0, cot(z) = -j*sign(y) - // if cot(x) = +/- oo or 1/sin^2(x) = 1 + cot^2(x) = oo, cot(z) = -j*cosh(y)/sinh(y) - // - // The algorithm is based on Kahan. - - double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 - - if (Math.Abs(value.Imaginary) > upperLimit) - { - double e = Math.Exp(-Math.Abs(value.Imaginary)); - return e == 0.0 - ? new Complex(0.0, -Math.Sign(value.Imaginary)) - : new Complex(4.0 * Math.Cos(value.Real) * Math.Sin(value.Real) * e * e, -Math.Sign(value.Imaginary)); - } - - double cotr = Cot(value.Real); - double beta = 1 + cotr * cotr; // beta = 1/cos(x)^2 = 1 + t^2 - double sinhi = Sinh(value.Imaginary); - double coshi = Cosh(value.Imaginary); - - if (double.IsInfinity(cotr)) - return new Complex(0.0, - coshi / sinhi); - - double denom = 1.0 + beta * sinhi * sinhi; + // cot(z) = - j*coth(-j*z) - return new Complex(cotr / denom, - beta * coshi * sinhi / denom); + Complex z = Coth(new Complex(value.Imaginary, -value.Real)); + return new Complex(z.Imaginary, -z.Real); } /// @@ -567,13 +514,13 @@ namespace MathNet.Numerics } // tanh(x + j*y) = (cosh(x)*sinh(x)/cos^2(y) + j*tan(y))/(1 + sinh^2(x)/cos^2(y)) - // if |x| > asinh(sqrt(max)), tanh(z) = sign(x) + j*4*cos(y)*sin(y)*exp(-2*|x|) + // if |x| > huge, tanh(z) = sign(x) + j*4*cos(y)*sin(y)*exp(-2*|x|) // if exp(-|x|) = 0, tanh(z) = sign(x) // if tan(y) = +/- oo or 1/cos^2(y) = 1 + tan^2(y) = oo, tanh(z) = cosh(x)/sinh(x) // // The algorithm is based on Kahan. - double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 + double upperLimit = 177.6189650; // Asinh(double.MaxValue) / 4.0; if (Math.Abs(value.Real) > upperLimit) { @@ -629,33 +576,9 @@ namespace MathNet.Numerics return new Complex(Coth(value.Real), 0.0); } - // coth(x + j*y) = (cosh(x)*sinh(x)/sin^2(y) - j*cot(y))/(1 + sinh^2(x)/sin^2(y)) - // if |x| > asinh(sqrt(max)), coth(z) = sign(x) - j*4*cos(y)*sin(y)*exp(-2*|x|) - // if exp(-|x|) = 0, coth(z) = sign(x) - // if cot(y) = +/- oo or 1/sin^2(y) = 1 + cot^2(y) = oo, coth(z) = cosh(x)/sinh(x) - // - // The algorithm is based on Kahan. - - double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 - - if (Math.Abs(value.Real) > upperLimit) - { - double e = Math.Exp(-Math.Abs(value.Real)); - return e == 0.0 - ? new Complex(Math.Sign(value.Real), 0.0) - : new Complex(Math.Sign(value.Real), -4.0 * Math.Cos(value.Imaginary) * Math.Sin(value.Imaginary) * e * e); - } + // Coth(z) = 1/tanh(z) - double coti = Cot(value.Imaginary); - double beta = 1 + coti * coti; - double sinhr = Math.Sinh(value.Real); - double coshr = Math.Cosh(value.Real); - - if (double.IsInfinity(coti)) - return new Complex(coshr / sinhr, 0.0); - - double denom = 1.0 + beta * sinhr * sinhr; - return new Complex(beta * coshr * sinhr / denom, -coti / denom); + return Complex.One / Tanh(value); } /// @@ -681,13 +604,13 @@ namespace MathNet.Numerics } // sech(x + j*y) = (cosh(x)/cos(y) - j*sinh(x)*tan(y)/cos(y))/(1 + sinh^2(x)/cos^2(y)) - // if |x| > asinh(sqrt(max)), sech(z) = 4*cosh(x)*cos(y)*exp(-2*|x|) - j*4*sinh(x)*tan(y)*cos(y)*exp(-2*|x|) + // if |x| > huge, sech(z) = 4*cosh(x)*cos(y)*exp(-2*|x|) - j*4*sinh(x)*tan(y)*cos(y)*exp(-2*|x|) // if exp(-|x|) = 0, sech(z) = 0 // if tan(y) = +/- oo or 1/cos^2(y) = 1 + tan^2(y) = oo, sech(z) = -j*sign(tan(y))/sinh(x) // // The algorithm is based on Kahan. - double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); // 355.58450362725193 + double upperLimit = 177.6189650; // Asinh(double.MaxValue) / 4.0; double tani = Tan(value.Imaginary); double cosi = Cos(value.Imaginary); @@ -735,13 +658,13 @@ namespace MathNet.Numerics } // csch(x + j*y) = (sinh(x)*cot(y)/sin(y) - j*cosh(x)/sin(y))/(1 + sinh^2(x)/sin^2(y)) - // if |x| > asinh(sqrt(max)), csch(z) = 4*sinh(x)*cot(y)*sin(y)*exp(-2*|x|) - j*4*cosh(x)*sin(y)*exp(-2*|x|) + // if |x| > huge, csch(z) = 4*sinh(x)*cot(y)*sin(y)*exp(-2*|x|) - j*4*cosh(x)*sin(y)*exp(-2*|x|) // if exp(-|x|) = 0, csch(z) = 0 // if cot(y) = +/- oo or 1/sin^2(x) = 1 + cot^2(x) = oo, csch(z) = sign(cot(y))/sinh(x) // // The algorithm is based on Kahan. - double upperLimit = Asinh(Math.Sqrt(double.MaxValue)); + double upperLimit = 177.6189650; // Asinh(double.MaxValue) / 4.0; double coti = Cot(value.Imaginary); double sini = Sin(value.Imaginary); diff --git a/src/UnitTests/TrigonometryTest.cs b/src/UnitTests/TrigonometryTest.cs index 2fb3194b..a7c4bc82 100644 --- a/src/UnitTests/TrigonometryTest.cs +++ b/src/UnitTests/TrigonometryTest.cs @@ -62,6 +62,8 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -39356457675156.6, -3139507838262.74)] [TestCase(32.0, -512.0, +9.52855589474963e+221, +6.29843301304523e+221)] [TestCase(-32.0, -512.0, +9.52855589474963e+221, -6.29843301304523e+221)] + [TestCase(355.0, 355.0, -7.4732769989672537746e153, +2.25277102712558920722e149)] + [TestCase(355.58, 355.58, -1.11645148443766450218e154, +7.31511888973169437082e153)] public void CanComputeComplexCosine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Cos(); @@ -87,6 +89,8 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -3139507838262.74, +39356457675156.6)] [TestCase(32.0, -512.0, +6.29843301304523e+221, -9.52855589474963e+221)] [TestCase(-32.0, -512.0, -6.29843301304523e+221, -9.52855589474963e+221)] + [TestCase(355.0, 355.0, -2.25277102712558920722e149, -7.4732769989672537746e153)] + [TestCase(355.58, 355.58, -+7.31511888973169437082e153, -1.11645148443766450218e154)] public void CanComputeComplexSine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sin(); @@ -112,6 +116,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, +5.08515122860094E-29, -1)] [TestCase(32.0, 512.0, +3.52598650243739e-445, 1)] [TestCase(32.0, -512.0, +3.52598650243739e-445, -1)] + [TestCase(355.0, 355.0, +5.39739014654622322424e-313, +1.0)] + [TestCase(355.58, 355.58, +2.57308259096360322163e-309, +1.0)] + [TestCase(196.0, 16.0, +1.62991625102252204721e-14, +1.00000000000001938715)] public void CanComputeComplexTangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Tan(); @@ -565,6 +572,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, 5.08515122860094E-29, +1.0)] [TestCase(32.0, 512.0, +3.52598650243739e-445, -1.0)] [TestCase(32.0, -512.0, +3.52598650243739e-445, +1.0)] + [TestCase(355.0, 355.0, +5.39739014654622322424e-313, -1.0)] + [TestCase(355.58, 355.58, +2.57308259096360322163e-309, -1.0)] + [TestCase(196.0, 16.0, +1.62991625102252204721e-14, -1.00000000000001938715)] public void CanComputeComplexCotangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Cot(); @@ -590,6 +600,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -2.52481261731330570134e-14, +2.01407074478744036413e-15)] [TestCase(32.0, 512.0, +7.30361056703505e-223, +4.82773070945482e-223)] [TestCase(32.0, -512.0, +7.30361056703505e-223, -4.82773070945482e-223)] + [TestCase(355.0, 355.0, -1.33810107564527561878e-154, -4.03361916732891580715e-159)] + [TestCase(355.58, 355.58, -6.26665948010988663559e-155, -4.10598756663013051681e-155)] + [TestCase(196.0, 16.0, +7.70790453157574313253e-8, +2.11460357915116320189e-7)] public void CanComputeComplexSecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sec(); @@ -615,6 +628,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -2.01407074478744e-15, -2.52481261731331e-14)] [TestCase(32.0, 512.0, +4.82773070945482e-223, -7.30361056703505e-223)] [TestCase(32.0, -512.0, +4.82773070945482e-223, +7.30361056703505e-223)] + [TestCase(355.0, 355.0, -4.03361916732891580715e-159, +1.33810107564527561878e-154)] + [TestCase(355.58, 355.58, -4.10598756663013051681e-155, +6.26665948010988663559e-155)] + [TestCase(196.0, 16.0, +2.11460357915116320189e-7, -7.70790453157574313253e-8)] public void CanComputeComplexCosecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Csc(); @@ -702,6 +718,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -1.0, 0.0)] [TestCase(-32.0, 512.0, -1.0, -5.08515122860093626173e-29)] [TestCase(32.0, -512.0, +1.0, +5.08515122860093626173e-29)] + [TestCase(355.0, 355.0, +1.0, +5.39739014654622322424e-313)] + [TestCase(355.58, 355.58, +1.0, +2.57308259096360322163e-309)] + [TestCase(196.0, 16.0, +1.0, +6.29623407730470174659e-171)] public void CanComputeComplexHyperbolicTangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Tanh(); @@ -731,6 +750,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -1.0, 0.0)] [TestCase(-32.0, 512.0, -1.0, 5.08515122860093626173e-29)] [TestCase(32.0, -512.0, +1.0, -5.08515122860093626173e-29)] + [TestCase(355.0, 355.0, +1.0, -5.39739014654622322424e-313)] + [TestCase(355.58, 355.58, +1.0, -2.57308259096360322163e-309)] + [TestCase(196.0, 16.0, +1.0, -6.29623407730470174659e-171)] public void CanComputeComplexHyperbolicCotangent(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Coth(); @@ -760,6 +782,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, 7.30361056703505051822e-223, -4.82773070945482080706e-223)] [TestCase(-32.0, 512.0, -2.52481261731330570134e-14, +2.01407074478744036413e-15)] [TestCase(-32.0, -512.0, -2.52481261731330570134e-14, -2.01407074478744036413e-15)] + [TestCase(355.0, 355.0, -1.33810107564527561878e-154, +4.03361916732891580715e-159)] + [TestCase(355.58, 355.58, -6.26665948010988663559e-155, +4.10598756663013051681e-155)] + [TestCase(196.0, 16.0, -1.447180343166980306269e-85, +4.350690711792111509525e-86)] public void CanComputeComplexHyperbolicSecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sech(); @@ -789,6 +814,9 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -7.30361056703505051822e-223, +4.82773070945482080706e-223)] [TestCase(-32.0, 512.0, +2.52481261731330570134e-14, -2.01407074478744036413e-15)] [TestCase(32.0, -512.0, -2.52481261731330570134e-14, +2.01407074478744036413e-15)] + [TestCase(355.0, 355.0, -1.33810107564527561878e-154, +4.03361916732891580715e-159)] + [TestCase(355.58, 355.58, -6.26665948010988663559e-155, +4.10598756663013051681e-155)] + [TestCase(196.0, 16.0, -1.447180343166980306269e-85, +4.350690711792111509525e-86)] public void CanComputeComplexHyperbolicCosecant(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Csch(); From 7d28c170826cdced3e5507a319cb321e678bba29 Mon Sep 17 00:00:00 2001 From: diluculo Date: Wed, 17 Jan 2018 18:45:16 +0900 Subject: [PATCH 3/6] Set huge = 22.0 Apply approximation to Sinh(z) and Cosh(z). --- src/Numerics/Trigonometry.cs | 38 ++++++++++++++++++++++--------- src/UnitTests/TrigonometryTest.cs | 2 ++ 2 files changed, 29 insertions(+), 11 deletions(-) diff --git a/src/Numerics/Trigonometry.cs b/src/Numerics/Trigonometry.cs index 45b4b8fa..fc75847d 100644 --- a/src/Numerics/Trigonometry.cs +++ b/src/Numerics/Trigonometry.cs @@ -447,9 +447,20 @@ namespace MathNet.Numerics return new Complex(Sinh(value.Real), 0.0); } + // sinh(x + j y) = sinh(x)*cos(y) + j*cosh(x)*sin(y) + // if x > huge, sinh(x + jy) = sign(x)*exp(|x|)/2*cos(y) + j*exp(|x|)/2*sin(y) + + if (Math.Abs(value.Real) >= 22.0) // Taken from the msun library in FreeBSD + { + double h = Math.Exp(Math.Abs(value.Real)) * 0.5; + return new Complex( + Math.Sign(value.Real)*h*Cos(value.Imaginary), + h*Sin(value.Imaginary)); + } + return new Complex( Sinh(value.Real) * Cos(value.Imaginary), - Cosh(value.Real) * Sin(value.Imaginary)); + Cosh(value.Real) * Sin(value.Imaginary)); } /// @@ -474,6 +485,17 @@ namespace MathNet.Numerics return new Complex(Cosh(value.Real), 0.0); } + // cosh(x + j*y) = cosh(x)*cos(y) + j*sinh(x)*sin(y) + // if x > huge, cosh(x + j*y) = exp(|x|)/2*cos(y) + j*sign(x)*exp(|x|)/2*sin(y) + + if (Math.Abs(value.Real) >= 22.0) // Taken from the msun library in FreeBSD + { + double h = Math.Exp(Math.Abs(value.Real)) * 0.5; + return new Complex( + h * Cos(value.Imaginary), + Math.Sign(value.Real) * h * Sin(value.Imaginary)); + } + return new Complex( Cosh(value.Real) * Cos(value.Imaginary), Sinh(value.Real) * Sin(value.Imaginary)); @@ -519,10 +541,8 @@ namespace MathNet.Numerics // if tan(y) = +/- oo or 1/cos^2(y) = 1 + tan^2(y) = oo, tanh(z) = cosh(x)/sinh(x) // // The algorithm is based on Kahan. - - double upperLimit = 177.6189650; // Asinh(double.MaxValue) / 4.0; - - if (Math.Abs(value.Real) > upperLimit) + + if (Math.Abs(value.Real) >= 22.0) // Taken from the msun library in FreeBSD { double e = Math.Exp(-Math.Abs(value.Real)); return e == 0.0 @@ -610,15 +630,13 @@ namespace MathNet.Numerics // // The algorithm is based on Kahan. - double upperLimit = 177.6189650; // Asinh(double.MaxValue) / 4.0; - double tani = Tan(value.Imaginary); double cosi = Cos(value.Imaginary); double beta = 1.0 + tani * tani; double sinhr = Math.Sinh(value.Real); double coshr = Math.Cosh(value.Real); - if (Math.Abs(value.Real) > upperLimit) + if (Math.Abs(value.Real) >= 22.0) // Taken from the msun library in FreeBSD { double e = Math.Exp(-Math.Abs(value.Real)); return e == 0.0 @@ -664,15 +682,13 @@ namespace MathNet.Numerics // // The algorithm is based on Kahan. - double upperLimit = 177.6189650; // Asinh(double.MaxValue) / 4.0; - double coti = Cot(value.Imaginary); double sini = Sin(value.Imaginary); double beta = 1 + coti * coti; double sinhr = Sinh(value.Real); double coshr = Cosh(value.Real); - if (Math.Abs(value.Real) > upperLimit) + if (Math.Abs(value.Real) >= 22.0) // Taken from the msun library in FreeBSD { double e = Math.Exp(-Math.Abs(value.Real)); return e == 0.0 diff --git a/src/UnitTests/TrigonometryTest.cs b/src/UnitTests/TrigonometryTest.cs index a7c4bc82..c55e0d27 100644 --- a/src/UnitTests/TrigonometryTest.cs +++ b/src/UnitTests/TrigonometryTest.cs @@ -660,6 +660,7 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, -9.52855589474962752845e221, -6.29843301304522728211e221)] [TestCase(-32.0, 512.0, +3.93564576751566e13, +3.13950783826274e12)] [TestCase(32.0, -512.0, -3.93564576751566e13, -3.13950783826274e12)] + [TestCase(709.0, 709.0, +2.22042071181169647963e307, -3.45764185010438140812e307)] public void CanComputeComplexHyperbolicSine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Sinh(); @@ -689,6 +690,7 @@ namespace MathNet.Numerics.UnitTests [TestCase(-512.0, -32.0, 9.52855589474962752845e221, 6.29843301304522728211e221)] [TestCase(-32.0, 512.0, -3.93564576751566e13, -3.13950783826274e12)] [TestCase(32.0, -512.0, -3.93564576751566e13, -3.13950783826274e12)] + [TestCase(709.0, 709.0, +2.22042071181169647963e307, -3.45764185010438140812e307)] public void CanComputeComplexHyperbolicCosine(double real, double imag, double expectedReal, double expectedImag) { var actual = new Complex(real, imag).Cosh(); From 118a3c028e42709b8b593ad5c40085e5c4a9a332 Mon Sep 17 00:00:00 2001 From: diluculo Date: Wed, 17 Jan 2018 21:26:15 +0900 Subject: [PATCH 4/6] Apply approximation to Asinh(x) and Acosh(x) --- src/Numerics/Trigonometry.cs | 14 +++++++++++++- src/UnitTests/TrigonometryTest.cs | 3 +++ 2 files changed, 16 insertions(+), 1 deletion(-) diff --git a/src/Numerics/Trigonometry.cs b/src/Numerics/Trigonometry.cs index fc75847d..73c6c5c6 100644 --- a/src/Numerics/Trigonometry.cs +++ b/src/Numerics/Trigonometry.cs @@ -713,7 +713,13 @@ namespace MathNet.Numerics /// The hyperbolic angle, i.e. the area of its hyperbolic sector. public static double Asinh(double value) { - return Math.Log(value + Math.Sqrt((value * value) + 1), Math.E); + // asinh(x) = Sign(x) * ln(|x| + sqrt(x*x + 1)) + // if |x| > huge, asinh(x) ~= Sign(x) * ln(2|x|) + + if (Math.Abs(value) >= 268435456.0) // 2^28, taken from freeBSD + return Math.Sign(value) * (Math.Log(Math.Abs(value)) + Math.Log(2.0)); + + return Math.Sign(value) * Math.Log(Math.Abs(value) + Math.Sqrt((value * value) + 1)); } /// @@ -733,6 +739,12 @@ namespace MathNet.Numerics /// The hyperbolic angle, i.e. the area of its hyperbolic sector. public static double Acosh(double value) { + // acosh(x) = ln(x + sqrt(x*x - 1)) + // if |x| >= 2^28, acosh(x) ~ ln(x) + ln(2) + + if (Math.Abs(value) >= 268435456.0) // 2^28, taken from freeBSD + return Math.Log(value) + Math.Log(2.0); + return Math.Log(value + (Math.Sqrt(value - 1) * Math.Sqrt(value + 1)), Math.E); } diff --git a/src/UnitTests/TrigonometryTest.cs b/src/UnitTests/TrigonometryTest.cs index c55e0d27..f853d401 100644 --- a/src/UnitTests/TrigonometryTest.cs +++ b/src/UnitTests/TrigonometryTest.cs @@ -338,6 +338,7 @@ namespace MathNet.Numerics.UnitTests /// Expected value. [TestCase(1.0, 0.0)] [TestCase(8388608, 16.635532333438682)] + [TestCase(1.7976931348623157E+308, 710.47586007394394203711)] public void CanComputeInverseHyperbolicCosine(double value, double expected) { var actual = Trig.Acosh(value); @@ -383,6 +384,8 @@ namespace MathNet.Numerics.UnitTests [TestCase(-8388608, -16.63553233343869)] [TestCase(1.19209289550780998537e-7, 1.1920928955078072e-7)] [TestCase(-1.19209289550780998537e-7, -1.1920928955078072e-7)] + [TestCase(1.7976931348623157E+308, 710.47586007394394203711)] + [TestCase(-1.7976931348623157E+308, -710.47586007394394203711)] public void CanComputeInverseHyperbolicSine(double value, double expected) { var actual = Trig.Asinh(value); From 0263f20c3d86dd4daa11136bb7d46f1535ec05c6 Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Mon, 22 Jan 2018 19:34:12 +0100 Subject: [PATCH 5/6] BUG: Distributions: BetaScaled no longer ignores optional random source parameter #547 --- src/Numerics/Distributions/BetaScaled.cs | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Numerics/Distributions/BetaScaled.cs b/src/Numerics/Distributions/BetaScaled.cs index 9d234b4b..58eeaa3b 100644 --- a/src/Numerics/Distributions/BetaScaled.cs +++ b/src/Numerics/Distributions/BetaScaled.cs @@ -80,7 +80,7 @@ namespace MathNet.Numerics.Distributions throw new ArgumentException(Resources.InvalidDistributionParameters); } - _random = SystemRandomSource.Default; + _random = randomSource ?? SystemRandomSource.Default; _shapeA = a; _shapeB = b; _location = location; From a7cc9def57d624c5519f5108c5cc5ebe5f75da29 Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Mon, 22 Jan 2018 19:47:19 +0100 Subject: [PATCH 6/6] Release v3.20.2 --- .paket/Paket.Restore.targets | 6 +++++- RELEASENOTES.md | 5 +++++ src/FSharp/AssemblyInfo.fs | 6 +++--- src/FSharpUnitTests/AssemblyInfo.fs | 6 +++--- src/Numerics/Properties/AssemblyInfo.cs | 6 +++--- src/UnitTests/Properties/AssemblyInfo.cs | 6 +++--- 6 files changed, 22 insertions(+), 13 deletions(-) diff --git a/.paket/Paket.Restore.targets b/.paket/Paket.Restore.targets index a86be3a1..830e5699 100644 --- a/.paket/Paket.Restore.targets +++ b/.paket/Paket.Restore.targets @@ -23,6 +23,9 @@ <_PaketExeExtension>$([System.IO.Path]::GetExtension("$(PaketExePath)")) dotnet "$(PaketExePath)" + + "$(PaketExePath)" + $(PaketRootPath)paket.bootstrapper.exe $(PaketToolsPath)paket.bootstrapper.exe "$(PaketBootStrapperExePath)" @@ -145,9 +148,10 @@ + diff --git a/RELEASENOTES.md b/RELEASENOTES.md index f702528b..514aab09 100644 --- a/RELEASENOTES.md +++ b/RELEASENOTES.md @@ -1,3 +1,8 @@ +### 3.20.2 - 2018-01-22 +* Bugfixes backported from v4: +* BUG: Distributions: BetaScaled no longer ignores optional random source parameter +* BUG: Trigonometry: Complex trigonometric functions behave on large imput *~diluculo* + ### 3.20.1 - 2018-01-13 * Bugfixes backported from v4: * BUG: Trigonometry: Fix imaginary part sign of complex hyperbolic cotangent diff --git a/src/FSharp/AssemblyInfo.fs b/src/FSharp/AssemblyInfo.fs index 4535d34b..172b8250 100644 --- a/src/FSharp/AssemblyInfo.fs +++ b/src/FSharp/AssemblyInfo.fs @@ -44,9 +44,9 @@ open System.Runtime.InteropServices [] [] -[] -[] -[] +[] +[] +[] #if PORTABLE #else diff --git a/src/FSharpUnitTests/AssemblyInfo.fs b/src/FSharpUnitTests/AssemblyInfo.fs index b889ff7e..6e53935c 100644 --- a/src/FSharpUnitTests/AssemblyInfo.fs +++ b/src/FSharpUnitTests/AssemblyInfo.fs @@ -10,9 +10,9 @@ open System.Runtime.InteropServices [] [] -[] -[] -[] +[] +[] +[] #if PORTABLE #else diff --git a/src/Numerics/Properties/AssemblyInfo.cs b/src/Numerics/Properties/AssemblyInfo.cs index ecf98445..5dcc5f16 100644 --- a/src/Numerics/Properties/AssemblyInfo.cs +++ b/src/Numerics/Properties/AssemblyInfo.cs @@ -45,9 +45,9 @@ using System.Runtime.InteropServices; [assembly: CLSCompliant(true)] [assembly: NeutralResourcesLanguage("en")] -[assembly: AssemblyVersion("3.20.1.0")] -[assembly: AssemblyFileVersion("3.20.1.0")] -[assembly: AssemblyInformationalVersion("3.20.1")] +[assembly: AssemblyVersion("3.20.2.0")] +[assembly: AssemblyFileVersion("3.20.2.0")] +[assembly: AssemblyInformationalVersion("3.20.2")] #if PORTABLE diff --git a/src/UnitTests/Properties/AssemblyInfo.cs b/src/UnitTests/Properties/AssemblyInfo.cs index cd983ca4..4ede7c9b 100644 --- a/src/UnitTests/Properties/AssemblyInfo.cs +++ b/src/UnitTests/Properties/AssemblyInfo.cs @@ -9,8 +9,8 @@ using MathNet.Numerics.UnitTests; [assembly: ComVisible(false)] [assembly: Guid("04157581-63f3-447b-a277-83c6e69126a4")] -[assembly: AssemblyVersion("3.20.1.0")] -[assembly: AssemblyFileVersion("3.20.1.0")] -[assembly: AssemblyInformationalVersion("3.20.1")] +[assembly: AssemblyVersion("3.20.2.0")] +[assembly: AssemblyFileVersion("3.20.2.0")] +[assembly: AssemblyInformationalVersion("3.20.2")] [assembly: UseLinearAlgebraProvider]