From 6cdafe62be2b3fc80714be5e7aaaa668fa0a5fb7 Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Thu, 25 Dec 2014 20:01:22 +0100 Subject: [PATCH] Minor special functions cleanup --- CONTRIBUTORS.md | 2 + src/Numerics/Numerics.csproj | 2 +- src/Numerics/SpecialFunctions/Beta.cs | 5 +- src/Numerics/SpecialFunctions/Erf.cs | 4 +- src/Numerics/SpecialFunctions/Evaluate.cs | 38 ++++++------- .../SpecialFunctions/ExponentialIntegral.cs | 56 ++++++++++--------- src/Numerics/SpecialFunctions/Factorial.cs | 10 ++-- src/Numerics/SpecialFunctions/Gamma.cs | 4 +- src/Numerics/SpecialFunctions/Harmonic.cs | 4 +- src/Numerics/SpecialFunctions/Logistic.cs | 10 ++-- .../SpecialFunctions/ModifiedBessel.cs | 4 +- .../SpecialFunctions/ModifiedStruve.cs | 38 ++++++------- src/Numerics/SpecialFunctions/Stability.cs | 44 +++++++-------- .../{Test.cs => TestFunctions.cs} | 0 14 files changed, 114 insertions(+), 107 deletions(-) rename src/Numerics/SpecialFunctions/{Test.cs => TestFunctions.cs} (100%) diff --git a/CONTRIBUTORS.md b/CONTRIBUTORS.md index 0adeace9..06953a54 100644 --- a/CONTRIBUTORS.md +++ b/CONTRIBUTORS.md @@ -33,6 +33,7 @@ Feel free to add a link to your personal site/blog and/or twitter handle.* - Scott Stephens - Superbest - Anders Gustafsson +- Ashley Messer - Gauthier Segay - Patrick van der Velde - Robin Neatherway @@ -49,6 +50,7 @@ Feel free to add a link to your personal site/blog and/or twitter handle.* - Kosei ABE - Martin Posch - Paul Varkey +- Sunny Ahuwanya - Till Hoffmann - Tomas Petricek - VicPara diff --git a/src/Numerics/Numerics.csproj b/src/Numerics/Numerics.csproj index c741b621..d88a75ae 100644 --- a/src/Numerics/Numerics.csproj +++ b/src/Numerics/Numerics.csproj @@ -203,7 +203,7 @@ - + diff --git a/src/Numerics/SpecialFunctions/Beta.cs b/src/Numerics/SpecialFunctions/Beta.cs index 9e3064a0..71ed278a 100644 --- a/src/Numerics/SpecialFunctions/Beta.cs +++ b/src/Numerics/SpecialFunctions/Beta.cs @@ -33,12 +33,13 @@ // ALGLIB 2.0.1, Sergey Bochkanov // +using System; +using MathNet.Numerics.Properties; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - using Properties; public static partial class SpecialFunctions { diff --git a/src/Numerics/SpecialFunctions/Erf.cs b/src/Numerics/SpecialFunctions/Erf.cs index b67d3705..7428d20d 100644 --- a/src/Numerics/SpecialFunctions/Erf.cs +++ b/src/Numerics/SpecialFunctions/Erf.cs @@ -37,12 +37,12 @@ // * double ErfInvImpl(double p, double q, double s) // +using System; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - /// /// This partial implementation of the SpecialFunctions class contains all methods related to the error function. /// diff --git a/src/Numerics/SpecialFunctions/Evaluate.cs b/src/Numerics/SpecialFunctions/Evaluate.cs index 52d6dd86..d36be692 100644 --- a/src/Numerics/SpecialFunctions/Evaluate.cs +++ b/src/Numerics/SpecialFunctions/Evaluate.cs @@ -32,10 +32,10 @@ // CERN - European Laboratory for Particle Physics // http://www.docjar.com/html/api/cern/jet/math/Bessel.java.html // Copyright 1999 CERN - European Laboratory for Particle Physics. -// Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose -// is hereby granted without fee, provided that the above copyright notice appear in all copies and -// that both that copyright notice and this permission notice appear in supporting documentation. -// CERN makes no representations about the suitability of this software for any purpose. +// Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose +// is hereby granted without fee, provided that the above copyright notice appear in all copies and +// that both that copyright notice and this permission notice appear in supporting documentation. +// CERN makes no representations about the suitability of this software for any purpose. // It is provided "as is" without expressed or implied warranty. // TOMS757 - Uncommon Special Functions (Fortran77) by Allan McLeod // http://people.sc.fsu.edu/~jburkardt/f77_src/toms757/toms757.html @@ -44,16 +44,16 @@ // ALGLIB 2.0.1, Sergey Bochkanov // -// ReSharper disable CheckNamespace -namespace MathNet.Numerics -// ReSharper restore CheckNamespace -{ - using System; +using System; #if !NOSYSNUMERICS - using Complex = System.Numerics.Complex; +using System.Numerics; #endif +// ReSharper disable CheckNamespace +namespace MathNet.Numerics +// ReSharper restore CheckNamespace +{ /// /// Evaluation functions, useful for function approximation. /// @@ -140,7 +140,7 @@ namespace MathNet.Numerics compensation -= y; sum = t; } - while (Math.Abs(sum) < Math.Abs(factor * current)); + while (Math.Abs(sum) < Math.Abs(factor*current)); return sum; } @@ -198,11 +198,11 @@ namespace MathNet.Numerics { b2 = b1; b1 = b0; - b0 = x * b1 - b2 + coefficients[p++]; + b0 = x*b1 - b2 + coefficients[p++]; } while (--i > 0); - return 0.5 * (b0 - b2); + return 0.5*(b0 - b2); } /// @@ -234,10 +234,10 @@ namespace MathNet.Numerics { u2 = u1; u1 = u0; - u0 = xx * u1 + coefficients[i] - u2; + u0 = xx*u1 + coefficients[i] - u2; } - return (u0 - u2) / 2.0; + return (u0 - u2)/2.0; } // If ABS ( T ) > = 0.6 use the Reinsch modification @@ -254,11 +254,11 @@ namespace MathNet.Numerics { d2 = d1; double u2 = u1; - d1 = xx * u2 + coefficients[i] + d2; + d1 = xx*u2 + coefficients[i] + d2; u1 = d1 + u2; } - return (d1 + d2) / 2.0; + return (d1 + d2)/2.0; } else { @@ -273,11 +273,11 @@ namespace MathNet.Numerics { d2 = d1; double u2 = u1; - d1 = xx * u2 + coefficients[i] - d2; + d1 = xx*u2 + coefficients[i] - d2; u1 = d1 - u2; } - return (d1 - d2) / 2.0; + return (d1 - d2)/2.0; } } } diff --git a/src/Numerics/SpecialFunctions/ExponentialIntegral.cs b/src/Numerics/SpecialFunctions/ExponentialIntegral.cs index 228ccf7f..66a159e0 100644 --- a/src/Numerics/SpecialFunctions/ExponentialIntegral.cs +++ b/src/Numerics/SpecialFunctions/ExponentialIntegral.cs @@ -32,12 +32,12 @@ // Ashley Messer // +using System; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - public static partial class SpecialFunctions { /// @@ -53,16 +53,17 @@ namespace MathNet.Numerics /// "Advanced mathematical methods for scientists and engineers", Bender, Carl M.; Steven A. Orszag (1978). page 253 /// /// - /// for x > 1 uses continued fraction approach that is often used to compute incomplete gamma. - /// for 0 < x <= 1 uses taylor series expansion + /// for x > 1 uses continued fraction approach that is often used to compute incomplete gamma. + /// for 0 < x <= 1 uses Taylor series expansion /// /// Our unit tests suggest that the accuracy of the Exponential Integral function is correct up to 13 floating point digits. /// public static double ExponentialIntegral(double x, int n) { //parameter validation - if (n < 0 || x < 0.0 ) { - throw new ArgumentOutOfRangeException(string.Format("x and n must be positive: x={0}, n={1}", x, n)); + if (n < 0 || x < 0.0) + { + throw new ArgumentOutOfRangeException(string.Format("x and n must be positive: x={0}, n={1}", x, n)); } const double epsilon = 0.00000000000000001; @@ -79,57 +80,60 @@ namespace MathNet.Numerics //special cases if (n == 0) { - result = Math.Exp( -1.0d * x ) / x; + result = Math.Exp(-1.0d*x)/x; return result; } else if (x == 0.0d) { - result = 1.0d / (ndbl - 1.0d); + result = 1.0d/(ndbl - 1.0d); return result; } //general cases //continued fraction for large x - if (x > 1.0d) - { + if (x > 1.0d) + { b = x + ((double)n); - c = 1.0d / nearDoubleMin; - d = 1.0d / b; + c = 1.0d/nearDoubleMin; + d = 1.0d/b; h = d; for (i = 1; i <= maxIterations; i++) { - a = -1.0d * ((double)i) * ((ndbl - 1.0d) + (double)i); + a = -1.0d*((double)i)*((ndbl - 1.0d) + (double)i); b += 2.0d; - d = 1.0d / (a * d + b); - c = b + a / c; - del = c * d; - h = h * del; + d = 1.0d/(a*d + b); + c = b + a/c; + del = c*d; + h = h*del; if (Math.Abs(del - 1.0d) < epsilon) { - result = h * Math.Exp( -x ); + result = h*Math.Exp(-x); return result; } } throw new ArithmeticException(string.Format("continued fraction failed to converge for x={0}, n={1})", x, n)); } - //series computation for small x + //series computation for small x else { - result = ((ndbl - 1.0d) != 0 ? 1.0 / (ndbl - 1.0d) : (-1.0d * Math.Log(x) - Constants.EulerMascheroni)); //Set first term. + result = ((ndbl - 1.0d) != 0 ? 1.0/(ndbl - 1.0d) : (-1.0d*Math.Log(x) - Constants.EulerMascheroni)); //Set first term. for (i = 1; i <= maxIterations; i++) { - factorial *= (-1.0d * x / ((double)i)); - if (i != (ndbl - 1.0d)) { del = -factorial / (i - (ndbl - 1.0d)); } + factorial *= (-1.0d*x/((double)i)); + if (i != (ndbl - 1.0d)) + { + del = -factorial/(i - (ndbl - 1.0d)); + } else { - psi = -1.0d * Constants.EulerMascheroni; + psi = -1.0d*Constants.EulerMascheroni; for (ii = 1; ii <= (ndbl - 1.0d); ii++) { - psi += (1.0d / ((double)ii)); + psi += (1.0d/((double)ii)); } - del = factorial * (-1.0d * Math.Log(x) + psi); + del = factorial*(-1.0d*Math.Log(x) + psi); } result += del; - if (Math.Abs(del) < Math.Abs(result) * epsilon) + if (Math.Abs(del) < Math.Abs(result)*epsilon) { return result; } diff --git a/src/Numerics/SpecialFunctions/Factorial.cs b/src/Numerics/SpecialFunctions/Factorial.cs index 20e75b65..9676d0e0 100644 --- a/src/Numerics/SpecialFunctions/Factorial.cs +++ b/src/Numerics/SpecialFunctions/Factorial.cs @@ -41,8 +41,8 @@ namespace MathNet.Numerics { public partial class SpecialFunctions { - private const int FactorialMaxArgument = 170; - private static double[] _factorialCache; + const int FactorialMaxArgument = 170; + static double[] _factorialCache; /// /// Initializes static members of the SpecialFunctions class. @@ -52,13 +52,13 @@ namespace MathNet.Numerics InitializeFactorial(); } - private static void InitializeFactorial() + static void InitializeFactorial() { _factorialCache = new double[FactorialMaxArgument + 1]; _factorialCache[0] = 1.0; for (int i = 1; i < _factorialCache.Length; i++) { - _factorialCache[i] = _factorialCache[i - 1] * i; + _factorialCache[i] = _factorialCache[i - 1]*i; } } @@ -214,4 +214,4 @@ namespace MathNet.Numerics return Math.Floor(0.5 + Math.Exp(ret)); } } -} \ No newline at end of file +} diff --git a/src/Numerics/SpecialFunctions/Gamma.cs b/src/Numerics/SpecialFunctions/Gamma.cs index f2d4e4b7..043c8d5f 100644 --- a/src/Numerics/SpecialFunctions/Gamma.cs +++ b/src/Numerics/SpecialFunctions/Gamma.cs @@ -33,12 +33,12 @@ // ALGLIB 2.0.1, Sergey Bochkanov // +using System; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - public static partial class SpecialFunctions { /// diff --git a/src/Numerics/SpecialFunctions/Harmonic.cs b/src/Numerics/SpecialFunctions/Harmonic.cs index 0f723f73..233c1e2d 100644 --- a/src/Numerics/SpecialFunctions/Harmonic.cs +++ b/src/Numerics/SpecialFunctions/Harmonic.cs @@ -33,12 +33,12 @@ // ALGLIB 2.0.1, Sergey Bochkanov // +using System; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - /// /// This partial implementation of the SpecialFunctions class contains all methods related to the harmonic function. /// diff --git a/src/Numerics/SpecialFunctions/Logistic.cs b/src/Numerics/SpecialFunctions/Logistic.cs index 0149c818..2bf764b9 100644 --- a/src/Numerics/SpecialFunctions/Logistic.cs +++ b/src/Numerics/SpecialFunctions/Logistic.cs @@ -33,13 +33,13 @@ // ALGLIB 2.0.1, Sergey Bochkanov // +using System; +using MathNet.Numerics.Properties; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - using Properties; - /// /// This partial implementation of the SpecialFunctions class contains all methods related to the logistic function. /// @@ -52,7 +52,7 @@ namespace MathNet.Numerics /// The logistic function of . public static double Logistic(double p) { - return 1.0 / (Math.Exp(-p) + 1.0); + return 1.0/(Math.Exp(-p) + 1.0); } /// @@ -68,7 +68,7 @@ namespace MathNet.Numerics throw new ArgumentOutOfRangeException("p", Resources.ArgumentBetween0And1); } - return Math.Log(p / (1.0 - p)); + return Math.Log(p/(1.0 - p)); } } } diff --git a/src/Numerics/SpecialFunctions/ModifiedBessel.cs b/src/Numerics/SpecialFunctions/ModifiedBessel.cs index 9d13f3c1..b28a5a23 100644 --- a/src/Numerics/SpecialFunctions/ModifiedBessel.cs +++ b/src/Numerics/SpecialFunctions/ModifiedBessel.cs @@ -44,12 +44,12 @@ // ALGLIB 2.0.1, Sergey Bochkanov // +using System; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - /// /// This partial implementation of the SpecialFunctions class contains all methods related to the modified bessel function. /// diff --git a/src/Numerics/SpecialFunctions/ModifiedStruve.cs b/src/Numerics/SpecialFunctions/ModifiedStruve.cs index 35fab64a..d3e5455a 100644 --- a/src/Numerics/SpecialFunctions/ModifiedStruve.cs +++ b/src/Numerics/SpecialFunctions/ModifiedStruve.cs @@ -44,12 +44,12 @@ // ALGLIB 2.0.1, Sergey Bochkanov // +using System; + // ReSharper disable CheckNamespace namespace MathNet.Numerics // ReSharper restore CheckNamespace { - using System; - /// /// This partial implementation of the SpecialFunctions class contains all methods related to the modified bessel function. /// @@ -246,11 +246,11 @@ namespace MathNet.Numerics { if (x < xlow) { - return Constants.TwoInvPi * x; + return Constants.TwoInvPi*x; } - double T = (4.0 * x - 24.0) / (x + 24.0); - return Constants.TwoInvPi * x * Evaluate.ChebyshevSum(nterm1, ARL0, T) * Math.Exp(x); + double T = (4.0*x - 24.0)/(x + 24.0); + return Constants.TwoInvPi*x*Evaluate.ChebyshevSum(nterm1, ARL0, T)*Math.Exp(x); } // Code for |xvalue| > 16 @@ -261,7 +261,7 @@ namespace MathNet.Numerics } else { - double T = (x - 28.0) / (4.0 - x); + double T = (x - 28.0)/(4.0 - x); ch1 = Evaluate.ChebyshevSum(nterm2, ARL0AS, T); } @@ -272,18 +272,18 @@ namespace MathNet.Numerics } else { - double xsq = x * x; - double T = (800.0 - xsq) / (288.0 + xsq); + double xsq = x*x; + double T = (800.0 - xsq)/(288.0 + xsq); ch2 = Evaluate.ChebyshevSum(nterm3, AI0ML0, T); } - double test = Math.Log(ch1) - Constants.LogSqrt2Pi - Math.Log(x) / 2.0 + x; + double test = Math.Log(ch1) - Constants.LogSqrt2Pi - Math.Log(x)/2.0 + x; if (test > Math.Log(xmax)) { throw new ArithmeticException("ERROR IN MISCFUN FUNCTION STRVL0: ARGUMENT CAUSES OVERFLOW"); } - return Math.Exp(test) - Constants.TwoInvPi * ch2 / x; + return Math.Exp(test) - Constants.TwoInvPi*ch2/x; } /// @@ -488,14 +488,14 @@ namespace MathNet.Numerics return 0.0; } - double xsq = x * x; + double xsq = x*x; if (x < xlow1) { - return xsq / Constants.Pi3Over2; + return xsq/Constants.Pi3Over2; } - double t = (4.0 * x - 24.0) / (x + 24.0); - return xsq * Evaluate.ChebyshevSum(nterm1, ARL1, t) * Math.Exp(x) / Constants.Pi3Over2; + double t = (4.0*x - 24.0)/(x + 24.0); + return xsq*Evaluate.ChebyshevSum(nterm1, ARL1, t)*Math.Exp(x)/Constants.Pi3Over2; } // CODE FOR |x| > 16 @@ -506,7 +506,7 @@ namespace MathNet.Numerics } else { - double t = (x - 30.0) / (2.0 - x); + double t = (x - 30.0)/(2.0 - x); ch1 = Evaluate.ChebyshevSum(nterm2, ARL1AS, t); } @@ -517,18 +517,18 @@ namespace MathNet.Numerics } else { - double xsq = x * x; - double t = (800.0 - xsq) / (288.0 + xsq); + double xsq = x*x; + double t = (800.0 - xsq)/(288.0 + xsq); ch2 = Evaluate.ChebyshevSum(nterm3, AI1ML1, t); } - double test = Math.Log(ch1) - Constants.LogSqrt2Pi - Math.Log(x) / 2.0 + x; + double test = Math.Log(ch1) - Constants.LogSqrt2Pi - Math.Log(x)/2.0 + x; if (test > Math.Log(xmax)) { throw new ArithmeticException("ERROR IN MISCFUN FUNCTION STRVL1: ARGUMENT CAUSES OVERFLOW"); } - return Math.Exp(test) - Constants.TwoInvPi * ch2; + return Math.Exp(test) - Constants.TwoInvPi*ch2; } /// diff --git a/src/Numerics/SpecialFunctions/Stability.cs b/src/Numerics/SpecialFunctions/Stability.cs index 471959fb..4940748c 100644 --- a/src/Numerics/SpecialFunctions/Stability.cs +++ b/src/Numerics/SpecialFunctions/Stability.cs @@ -28,16 +28,16 @@ // OTHER DEALINGS IN THE SOFTWARE. // -// ReSharper disable CheckNamespace -namespace MathNet.Numerics -// ReSharper restore CheckNamespace -{ - using System; +using System; #if !NOSYSNUMERICS - using Complex = System.Numerics.Complex; +using System.Numerics; #endif +// ReSharper disable CheckNamespace +namespace MathNet.Numerics +// ReSharper restore CheckNamespace +{ public partial class SpecialFunctions { /// @@ -81,15 +81,15 @@ namespace MathNet.Numerics { if (a.Magnitude > b.Magnitude) { - var r = b.Magnitude / a.Magnitude; - return a.Magnitude * Math.Sqrt(1 + (r * r)); + var r = b.Magnitude/a.Magnitude; + return a.Magnitude*Math.Sqrt(1 + (r*r)); } if (b != 0.0) { // NOTE (ruegg): not "!b.AlmostZero()" to avoid convergence issues (e.g. in SVD algorithm) - var r = a.Magnitude / b.Magnitude; - return b.Magnitude * Math.Sqrt(1 + (r * r)); + var r = a.Magnitude/b.Magnitude; + return b.Magnitude*Math.Sqrt(1 + (r*r)); } return 0d; @@ -105,15 +105,15 @@ namespace MathNet.Numerics { if (a.Magnitude > b.Magnitude) { - var r = b.Magnitude / a.Magnitude; - return a.Magnitude * (float)Math.Sqrt(1 + (r * r)); + var r = b.Magnitude/a.Magnitude; + return a.Magnitude*(float)Math.Sqrt(1 + (r*r)); } if (b != 0.0f) { // NOTE (ruegg): not "!b.AlmostZero()" to avoid convergence issues (e.g. in SVD algorithm) - var r = a.Magnitude / b.Magnitude; - return b.Magnitude * (float)Math.Sqrt(1 + (r * r)); + var r = a.Magnitude/b.Magnitude; + return b.Magnitude*(float)Math.Sqrt(1 + (r*r)); } return 0f; @@ -129,15 +129,15 @@ namespace MathNet.Numerics { if (Math.Abs(a) > Math.Abs(b)) { - double r = b / a; - return Math.Abs(a) * Math.Sqrt(1 + (r * r)); + double r = b/a; + return Math.Abs(a)*Math.Sqrt(1 + (r*r)); } if (b != 0.0) { // NOTE (ruegg): not "!b.AlmostZero()" to avoid convergence issues (e.g. in SVD algorithm) - double r = a / b; - return Math.Abs(b) * Math.Sqrt(1 + (r * r)); + double r = a/b; + return Math.Abs(b)*Math.Sqrt(1 + (r*r)); } return 0d; @@ -153,15 +153,15 @@ namespace MathNet.Numerics { if (Math.Abs(a) > Math.Abs(b)) { - float r = b / a; - return Math.Abs(a) * (float)Math.Sqrt(1 + (r * r)); + float r = b/a; + return Math.Abs(a)*(float)Math.Sqrt(1 + (r*r)); } if (b != 0.0) { // NOTE (ruegg): not "!b.AlmostZero()" to avoid convergence issues (e.g. in SVD algorithm) - float r = a / b; - return Math.Abs(b) * (float)Math.Sqrt(1 + (r * r)); + float r = a/b; + return Math.Abs(b)*(float)Math.Sqrt(1 + (r*r)); } return 0f; diff --git a/src/Numerics/SpecialFunctions/Test.cs b/src/Numerics/SpecialFunctions/TestFunctions.cs similarity index 100% rename from src/Numerics/SpecialFunctions/Test.cs rename to src/Numerics/SpecialFunctions/TestFunctions.cs