tsai
/
mathnet-numerics
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

style cop 

Signed-off-by: Christoph Ruegg <git@cdrnet.ch>
pull/2/head
Christoph Ruegg 17 years ago
parent
795b406854
commit
ddd3441649
6 changed files with 163 additions and 149 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Unified View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 29
      COPYRIGHT.markdown
  2. 2
      src/Managed/Complex.cs
  3. 4
      src/Managed/Control.cs
  4. 78
      src/Managed/Distributions/Continuous/ContinuousUniform.cs
  5. 111
      src/Managed/Distributions/Continuous/Normal.cs
  6. 88
      src/Managed/SpecialFunctions/Erf.cs

29
COPYRIGHT.markdown
View File

@ -61,4 +61,31 @@ LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
BOOST
-----
Boost Software License - Version 1.0 - August 17th, 2003
Permission is hereby granted, free of charge, to any person or organization
obtaining a copy of the software and accompanying documentation covered by
this license (the "Software") to use, reproduce, display, distribute,
execute, and transmit the Software, and to prepare derivative works of the
Software, and to permit third-parties to whom the Software is furnished to
do so, all subject to the following:
The copyright notices in the Software and this entire statement, including
the above license grant, this restriction and the following disclaimer,
must be included in all copies of the Software, in whole or in part, and
all derivative works of the Software, unless such copies or derivative
works are solely in the form of machine-executable object code generated by
a source language processor.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE, TITLE AND NON-INFRINGEMENT. IN NO EVENT
SHALL THE COPYRIGHT HOLDERS OR ANYONE DISTRIBUTING THE SOFTWARE BE LIABLE
FOR ANY DAMAGES OR OTHER LIABILITY, WHETHER IN CONTRACT, TORT OR OTHERWISE,
ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
DEALINGS IN THE SOFTWARE.

2
src/Managed/Complex.cs
View File

@ -286,7 +286,7 @@ namespace MathNet.Numerics
/// </value> /// </value>
public bool IsImaginary public bool IsImaginary
{ {
get { return _real.AlmostZero(); } get { return _real.AlmostZero(); }
} }
#region Static Initializers #region Static Initializers

4
src/Managed/Control.cs
View File

@ -1,4 +1,4 @@
// <copyright file="ContinuousUniformTests.cs" company="Math.NET"> // <copyright file="Control.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project // Math.NET Numerics, part of the Math.NET Project
// http://mathnet.opensourcedotnet.info // http://mathnet.opensourcedotnet.info
// //
@ -56,7 +56,7 @@ namespace MathNet.Numerics
/// Thread safe RNG about two and half time slower than non-thread safe RNG. /// Thread safe RNG about two and half time slower than non-thread safe RNG.
/// </summary> /// </summary>
/// <value> /// <value>
/// <c>true</c> to use thread safe random number generators ; otherwise, <c>false</c>. /// <c>true</c> to use thread safe random number generators ; otherwise, <c>false</c>.
/// </value> /// </value>
public static bool ThreadSafeRandomNumberGenerators { get; set; } public static bool ThreadSafeRandomNumberGenerators { get; set; }
} }

78
src/Managed/Distributions/Continuous/ContinuousUniform.cs
View File

@ -30,7 +30,7 @@ namespace MathNet.Numerics.Distributions
{ {
using System; using System;
using System.Collections.Generic; using System.Collections.Generic;
using MathNet.Numerics.Properties; using Properties;
/// <summary> /// <summary>
/// The continuous uniform distribution is a distribution over real numbers. For details about this distribution, see /// The continuous uniform distribution is a distribution over real numbers. For details about this distribution, see
@ -46,12 +46,12 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The distribution's lower bound. /// The distribution's lower bound.
/// </summary> /// </summary>
private double mLower; private double _lower;
/// <summary> /// <summary>
/// The distribution's upper bound. /// The distribution's upper bound.
/// </summary> /// </summary>
private double mUpper; private double _upper;
/// <summary> /// <summary>
/// Initializes a new instance of the ContinuousUniform class with lower bound 0 and upper bound 1. /// Initializes a new instance of the ContinuousUniform class with lower bound 0 and upper bound 1.
@ -69,7 +69,7 @@ namespace MathNet.Numerics.Distributions
public ContinuousUniform(double lower, double upper) public ContinuousUniform(double lower, double upper)
{ {
SetParameters(lower, upper); SetParameters(lower, upper);
RandomSource = new System.Random(); RandomSource = new Random();
} }
/// <summary> /// <summary>
@ -78,7 +78,7 @@ namespace MathNet.Numerics.Distributions
/// <returns>a string representation of the distribution.</returns> /// <returns>a string representation of the distribution.</returns>
public override string ToString() public override string ToString()
{ {
return "ContinuousUniform(Lower = " + mLower + ", Upper = " + mUpper + ")"; return "ContinuousUniform(Lower = " + _lower + ", Upper = " + _upper + ")";
} }
/// <summary> /// <summary>
@ -93,7 +93,8 @@ namespace MathNet.Numerics.Distributions
{ {
return false; return false;
} }
else if(Double.IsNaN(upper) || Double.IsNaN(lower))
if (Double.IsNaN(upper) || Double.IsNaN(lower))
{ {
return false; return false;
} }
@ -109,15 +110,13 @@ namespace MathNet.Numerics.Distributions
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double lower, double upper) private void SetParameters(double lower, double upper)
{ {
if (!Control.CheckDistributionParameters || IsValidParameterSet(lower, upper)) if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper))
{
mLower = lower;
mUpper = upper;
}
else
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
_lower = lower;
_upper = upper;
} }
/// <summary> /// <summary>
@ -127,12 +126,12 @@ namespace MathNet.Numerics.Distributions
{ {
get get
{ {
return mLower; return _lower;
} }
set set
{ {
SetParameters(value, mUpper); SetParameters(value, _upper);
} }
} }
@ -143,12 +142,12 @@ namespace MathNet.Numerics.Distributions
{ {
get get
{ {
return mUpper; return _upper;
} }
set set
{ {
SetParameters(mLower, value); SetParameters(_lower, value);
} }
} }
@ -164,7 +163,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get { return (mLower + mUpper)/2.0; } get { return (_lower + _upper) / 2.0; }
} }
/// <summary> /// <summary>
@ -172,7 +171,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get { return (mUpper - mLower)*(mUpper - mLower)/12.0; } get { return (_upper - _lower) * (_upper - _lower) / 12.0; }
} }
/// <summary> /// <summary>
@ -180,7 +179,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get { return (mUpper - mLower) / System.Math.Sqrt(12.0); } get { return (_upper - _lower) / Math.Sqrt(12.0); }
} }
/// <summary> /// <summary>
@ -189,7 +188,7 @@ namespace MathNet.Numerics.Distributions
/// <value></value> /// <value></value>
public double Entropy public double Entropy
{ {
get { return System.Math.Log(mUpper - mLower); } get { return Math.Log(_upper - _lower); }
} }
/// <summary> /// <summary>
@ -209,7 +208,7 @@ namespace MathNet.Numerics.Distributions
/// <value></value> /// <value></value>
public double Mode public double Mode
{ {
get { return (mLower + mUpper)/2.0; } get { return (_lower + _upper) / 2.0; }
} }
/// <summary> /// <summary>
@ -218,7 +217,7 @@ namespace MathNet.Numerics.Distributions
/// <value></value> /// <value></value>
public double Median public double Median
{ {
get { return (mLower + mUpper)/2.0; } get { return (_lower + _upper) / 2.0; }
} }
/// <summary> /// <summary>
@ -226,7 +225,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get { return mLower; } get { return _lower; }
} }
/// <summary> /// <summary>
@ -234,7 +233,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get { return mUpper; } get { return _upper; }
} }
/// <summary> /// <summary>
@ -244,9 +243,9 @@ namespace MathNet.Numerics.Distributions
/// <returns>the density at <paramref name="x"/>.</returns> /// <returns>the density at <paramref name="x"/>.</returns>
public double Density(double x) public double Density(double x)
{ {
if (x >= mLower && x <= mUpper) if (x >= _lower && x <= _upper)
{ {
return 1.0/(mUpper - mLower); return 1.0 / (_upper - _lower);
} }
return 0.0; return 0.0;
@ -259,9 +258,9 @@ namespace MathNet.Numerics.Distributions
/// <returns>the log density at <paramref name="x"/>.</returns> /// <returns>the log density at <paramref name="x"/>.</returns>
public double DensityLn(double x) public double DensityLn(double x)
{ {
if (x >= mLower && x <= mUpper) if (x >= _lower && x <= _upper)
{ {
return -Math.Log(mUpper - mLower); return -Math.Log(_upper - _lower);
} }
return Double.NegativeInfinity; return Double.NegativeInfinity;
@ -274,16 +273,17 @@ namespace MathNet.Numerics.Distributions
/// <returns>the cumulative density at <paramref name="x"/>.</returns> /// <returns>the cumulative density at <paramref name="x"/>.</returns>
public double CumulativeDistribution(double x) public double CumulativeDistribution(double x)
{ {
if(x <= mLower) if (x <= _lower)
{ {
return 0.0; return 0.0;
} }
else if(x >= mUpper)
if (x >= _upper)
{ {
return 1.0; return 1.0;
} }
return (x - mLower)/(mUpper - mLower); return (x - _lower) / (_upper - _lower);
} }
/// <summary> /// <summary>
@ -292,7 +292,7 @@ namespace MathNet.Numerics.Distributions
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, mLower, mUpper); return DoSample(RandomSource, _lower, _upper);
} }
/// <summary> /// <summary>
@ -301,9 +301,9 @@ namespace MathNet.Numerics.Distributions
/// <returns>a sequence of samples from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples() public IEnumerable<double> Samples()
{ {
while(true) while (true)
{ {
yield return DoSample(RandomSource, mLower, mUpper); yield return DoSample(RandomSource, _lower, _upper);
} }
} }
@ -316,7 +316,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="lower">The lower bound of the uniform random variable.</param> /// <param name="lower">The lower bound of the uniform random variable.</param>
/// <param name="upper">The upper bound of the uniform random variable.</param> /// <param name="upper">The upper bound of the uniform random variable.</param>
/// <returns>a uniformly distributed sample.</returns> /// <returns>a uniformly distributed sample.</returns>
public static double Sample(System.Random rnd, double lower, double upper) public static double Sample(Random rnd, double lower, double upper)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper)) if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper))
{ {
@ -333,14 +333,14 @@ namespace MathNet.Numerics.Distributions
/// <param name="lower">The lower bound of the uniform random variable.</param> /// <param name="lower">The lower bound of the uniform random variable.</param>
/// <param name="upper">The upper bound of the uniform random variable.</param> /// <param name="upper">The upper bound of the uniform random variable.</param>
/// <returns>a sequence of uniformly distributed samples.</returns> /// <returns>a sequence of uniformly distributed samples.</returns>
public static IEnumerable<double> Samples(System.Random rnd, double lower, double upper) public static IEnumerable<double> Samples(Random rnd, double lower, double upper)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper)) if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper))
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
while(true) while (true)
{ {
yield return DoSample(rnd, lower, upper); yield return DoSample(rnd, lower, upper);
} }
@ -353,9 +353,9 @@ namespace MathNet.Numerics.Distributions
/// <param name="lower">The lower bound of the uniform random variable.</param> /// <param name="lower">The lower bound of the uniform random variable.</param>
/// <param name="upper">The upper bound of the uniform random variable.</param> /// <param name="upper">The upper bound of the uniform random variable.</param>
/// <returns>a uniformly distributed random number.</returns> /// <returns>a uniformly distributed random number.</returns>
private static double DoSample(System.Random rnd, double lower, double upper) private static double DoSample(Random rnd, double lower, double upper)
{ {
return lower + (rnd.NextDouble()*(upper - lower)); return lower + (rnd.NextDouble() * (upper - lower));
} }
} }
} }

111
src/Managed/Distributions/Continuous/Normal.cs
View File

@ -30,8 +30,7 @@ namespace MathNet.Numerics.Distributions
{ {
using System; using System;
using System.Collections.Generic; using System.Collections.Generic;
using MathNet.Numerics; using Properties;
using MathNet.Numerics.Properties;
/// <summary> /// <summary>
/// Implements the univariate Normal (or Gaussian) distribution. For details about this distribution, see /// Implements the univariate Normal (or Gaussian) distribution. For details about this distribution, see
@ -47,12 +46,12 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Keeps track of the mean of the normal distribution. /// Keeps track of the mean of the normal distribution.
/// </summary> /// </summary>
private double mMean; private double _mean;
/// <summary> /// <summary>
/// Keeps track of the standard deviation of the normal distribution. /// Keeps track of the standard deviation of the normal distribution.
/// </summary> /// </summary>
private double mStdDev; private double _stdDev;
/// <summary> /// <summary>
/// Initializes a new instance of the Normal class. This is a normal distribution with mean 0.0 /// Initializes a new instance of the Normal class. This is a normal distribution with mean 0.0
@ -96,7 +95,7 @@ namespace MathNet.Numerics.Distributions
/// <returns>a normal distribution.</returns> /// <returns>a normal distribution.</returns>
public static Normal WithMeanVariance(double mean, double var) public static Normal WithMeanVariance(double mean, double var)
{ {
return new Normal(mean, System.Math.Sqrt(var)); return new Normal(mean, Math.Sqrt(var));
} }
/// <summary> /// <summary>
@ -108,7 +107,7 @@ namespace MathNet.Numerics.Distributions
/// <returns>a normal distribution.</returns> /// <returns>a normal distribution.</returns>
public static Normal WithMeanAndPrecision(double mean, double prec) public static Normal WithMeanAndPrecision(double mean, double prec)
{ {
return new Normal(mean, 1.0 / System.Math.Sqrt(prec)); return new Normal(mean, 1.0 / Math.Sqrt(prec));
} }
/// <summary> /// <summary>
@ -117,7 +116,7 @@ namespace MathNet.Numerics.Distributions
/// <returns>a string representation of the distribution.</returns> /// <returns>a string representation of the distribution.</returns>
public override string ToString() public override string ToString()
{ {
return "Normal(Mean = " + mMean + ", StdDev = " + mStdDev + ")"; return "Normal(Mean = " + _mean + ", StdDev = " + _stdDev + ")";
} }
/// <summary> /// <summary>
@ -128,15 +127,7 @@ namespace MathNet.Numerics.Distributions
/// <returns>True when the parameters are valid, false otherwise.</returns> /// <returns>True when the parameters are valid, false otherwise.</returns>
private static bool IsValidParameterSet(double mean, double stddev) private static bool IsValidParameterSet(double mean, double stddev)
{ {
if (stddev < 0.0) if (stddev < 0.0 || Double.IsNaN(mean) || Double.IsNaN(stddev))
{
return false;
}
else if (System.Double.IsNaN(mean))
{
return false;
}
else if (System.Double.IsNaN(stddev))
{ {
return false; return false;
} }
@ -152,17 +143,15 @@ namespace MathNet.Numerics.Distributions
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double mean, double stddev) private void SetParameters(double mean, double stddev)
{ {
if (!Control.CheckDistributionParameters || IsValidParameterSet(mean, stddev)) if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev))
{
mMean = mean;
mStdDev = stddev;
}
else
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
_mean = mean;
_stdDev = stddev;
} }
/// <summary> /// <summary>
/// Gets or sets the precision of the normal distribution. /// Gets or sets the precision of the normal distribution.
/// </summary> /// </summary>
@ -170,20 +159,20 @@ namespace MathNet.Numerics.Distributions
{ {
get get
{ {
return 1.0 / (mStdDev * mStdDev); return 1.0 / (_stdDev * _stdDev);
} }
set set
{ {
double sdev = 1.0/Math.Sqrt(value); double sdev = 1.0 / Math.Sqrt(value);
// Handle the case when the precision is -0. // Handle the case when the precision is -0.
if(Double.IsInfinity(sdev)) if (Double.IsInfinity(sdev))
{ {
sdev = Double.PositiveInfinity; sdev = Double.PositiveInfinity;
} }
SetParameters(mMean, sdev); SetParameters(_mean, sdev);
} }
} }
@ -199,8 +188,8 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get { return mMean; } get { return _mean; }
set { SetParameters(value, mStdDev); } set { SetParameters(value, _stdDev); }
} }
/// <summary> /// <summary>
@ -208,8 +197,8 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get { return mStdDev * mStdDev; } get { return _stdDev * _stdDev; }
set { SetParameters(mMean, value); } set { SetParameters(_mean, value); }
} }
/// <summary> /// <summary>
@ -217,8 +206,8 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get { return mStdDev; } get { return _stdDev; }
set { SetParameters(mMean, value); } set { SetParameters(_mean, value); }
} }
/// <summary> /// <summary>
@ -226,7 +215,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get { return Math.Log(mStdDev) + Constants.LogSqrt2PiE; } get { return Math.Log(_stdDev) + Constants.LogSqrt2PiE; }
} }
/// <summary> /// <summary>
@ -245,7 +234,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get { return mMean; } get { return _mean; }
} }
/// <summary> /// <summary>
@ -253,7 +242,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get { return mMean; } get { return _mean; }
} }
/// <summary> /// <summary>
@ -261,7 +250,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get { return System.Double.NegativeInfinity; } get { return Double.NegativeInfinity; }
} }
/// <summary> /// <summary>
@ -269,7 +258,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get { return System.Double.PositiveInfinity; } get { return Double.PositiveInfinity; }
} }
/// <summary> /// <summary>
@ -279,8 +268,8 @@ namespace MathNet.Numerics.Distributions
/// <returns>the density at <paramref name="x"/>.</returns> /// <returns>the density at <paramref name="x"/>.</returns>
public double Density(double x) public double Density(double x)
{ {
double d = (x - mMean) / mStdDev; double d = (x - _mean) / _stdDev;
return Math.Exp(-0.5*d*d) / (Constants.Sqrt2Pi*mStdDev); return Math.Exp(-0.5 * d * d) / (Constants.Sqrt2Pi * _stdDev);
} }
/// <summary> /// <summary>
@ -290,8 +279,8 @@ namespace MathNet.Numerics.Distributions
/// <returns>the log density at <paramref name="x"/>.</returns> /// <returns>the log density at <paramref name="x"/>.</returns>
public double DensityLn(double x) public double DensityLn(double x)
{ {
double d = (x - mMean) / mStdDev; double d = (x - _mean) / _stdDev;
return (-0.5 * d * d) - Math.Log(mStdDev) - Constants.LogSqrt2Pi; return (-0.5 * d * d) - Math.Log(_stdDev) - Constants.LogSqrt2Pi;
} }
/// <summary> /// <summary>
@ -301,7 +290,7 @@ namespace MathNet.Numerics.Distributions
/// <returns>the cumulative density at <paramref name="x"/>.</returns> /// <returns>the cumulative density at <paramref name="x"/>.</returns>
public double CumulativeDistribution(double x) public double CumulativeDistribution(double x)
{ {
return 0.5 * (1.0 + SpecialFunctions.Erf((x - mMean) / (mStdDev * System.Math.Sqrt(2.0)))); return 0.5 * (1.0 + SpecialFunctions.Erf((x - _mean) / (_stdDev * Math.Sqrt(2.0))));
} }
/// <summary> /// <summary>
@ -311,7 +300,7 @@ namespace MathNet.Numerics.Distributions
public double Sample() public double Sample()
{ {
double r2; double r2;
return mMean + (mStdDev * SampleBoxMuller(RandomSource, out r2)); return _mean + (_stdDev * SampleBoxMuller(RandomSource, out r2));
} }
/// <summary> /// <summary>
@ -325,8 +314,8 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
double r1 = SampleBoxMuller(RandomSource, out r2); double r1 = SampleBoxMuller(RandomSource, out r2);
yield return mMean + (mStdDev * r1); yield return _mean + (_stdDev * r1);
yield return mMean + (mStdDev * r2); yield return _mean + (_stdDev * r2);
} }
} }
#endregion #endregion
@ -335,10 +324,10 @@ namespace MathNet.Numerics.Distributions
/// Computes the inverse cumulative distribution function of the normal distribution. /// Computes the inverse cumulative distribution function of the normal distribution.
/// </summary> /// </summary>
/// <param name="p">The location at which to compute the inverse cumulative density.</param> /// <param name="p">The location at which to compute the inverse cumulative density.</param>
/// <returns>the inverse cumulative density at <paramref name="x"/>.</returns> /// <returns>the inverse cumulative density at <paramref name="p"/>.</returns>
public double InverseCumulativeDistribution(double p) public double InverseCumulativeDistribution(double p)
{ {
return mMean - (mStdDev * System.Math.Sqrt(2.0) * SpecialFunctions.ErfcInv(2.0 * p)); return _mean - (_stdDev * Math.Sqrt(2.0) * SpecialFunctions.ErfcInv(2.0 * p));
} }
/// <summary> /// <summary>
@ -348,9 +337,9 @@ namespace MathNet.Numerics.Distributions
/// <param name="mean">The mean of the normal distribution from which to generate samples.</param> /// <param name="mean">The mean of the normal distribution from which to generate samples.</param>
/// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param> /// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public static double Sample(System.Random rng, double mean, double stddev) public static double Sample(Random rng, double mean, double stddev)
{ {
if(Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev)) if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev))
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
@ -366,22 +355,20 @@ namespace MathNet.Numerics.Distributions
/// <param name="mean">The mean of the normal distribution from which to generate samples.</param> /// <param name="mean">The mean of the normal distribution from which to generate samples.</param>
/// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param> /// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param>
/// <returns>a sequence of samples from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(System.Random rng, double mean, double stddev) public static IEnumerable<double> Samples(Random rng, double mean, double stddev)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev)) if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev))
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
else
{ double r2;
double r2;
while (true) while (true)
{ {
double r1 = SampleBoxMuller(rng, out r2); double r1 = SampleBoxMuller(rng, out r2);
yield return mean + (stddev * r1); yield return mean + (stddev * r1);
yield return mean + (stddev * r2); yield return mean + (stddev * r2);
}
} }
} }
@ -391,7 +378,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="r2">A second random number from the standard normal distribution computed as a side product.</param> /// <param name="r2">A second random number from the standard normal distribution computed as a side product.</param>
/// <returns>a random number from the standard normal distribution.</returns> /// <returns>a random number from the standard normal distribution.</returns>
internal static double SampleBoxMuller(System.Random rnd, out double r2) internal static double SampleBoxMuller(Random rnd, out double r2)
{ {
double v1 = (2.0 * rnd.NextDouble()) - 1.0; double v1 = (2.0 * rnd.NextDouble()) - 1.0;
double v2 = (2.0 * rnd.NextDouble()) - 1.0; double v2 = (2.0 * rnd.NextDouble()) - 1.0;
@ -403,7 +390,7 @@ namespace MathNet.Numerics.Distributions
r = (v1 * v1) + (v2 * v2); r = (v1 * v1) + (v2 * v2);
} }
double fac = System.Math.Sqrt(-2.0 * System.Math.Log(r) / r); double fac = Math.Sqrt(-2.0 * Math.Log(r) / r);
r2 = v2 * fac; r2 = v2 * fac;
return v1 * fac; return v1 * fac;
} }

88
src/Managed/SpecialFunctions/Erf.cs
View File

@ -49,10 +49,10 @@ namespace MathNet.Numerics
/// <param name="x">The value to evaluate.</param> /// <param name="x">The value to evaluate.</param>
/// <returns>the error function evaluated at given value.</returns> /// <returns>the error function evaluated at given value.</returns>
/// <remarks> /// <remarks>
/// <list type="bullet"> /// <list type="bullet">
/// <item>returns 1 if <c>x == Double.PositiveInfinity</c>.</item> /// <item>returns 1 if <c>x == Double.PositiveInfinity</c>.</item>
/// <item>returns -1 if <c>x == Double.NegativeInfinity</c>.</item> /// <item>returns -1 if <c>x == Double.NegativeInfinity</c>.</item>
/// </list> /// </list>
/// </remarks> /// </remarks>
public static double Erf(double x) public static double Erf(double x)
{ {
@ -83,10 +83,10 @@ namespace MathNet.Numerics
/// <param name="x">The value to evaluate.</param> /// <param name="x">The value to evaluate.</param>
/// <returns>the complementary error function evaluated at given value.</returns> /// <returns>the complementary error function evaluated at given value.</returns>
/// <remarks> /// <remarks>
/// <list type="bullet"> /// <list type="bullet">
/// <item>returns 0 if <c>x == Double.PositiveInfinity</c>.</item> /// <item>returns 0 if <c>x == Double.PositiveInfinity</c>.</item>
/// <item>returns 2 if <c>x == Double.NegativeInfinity</c>.</item> /// <item>returns 2 if <c>x == Double.NegativeInfinity</c>.</item>
/// </list> /// </list>
/// </remarks> /// </remarks>
public static double Erfc(double x) public static double Erfc(double x)
{ {
@ -110,7 +110,7 @@ namespace MathNet.Numerics
return Double.NaN; return Double.NaN;
} }
return ErfImp(x,true); return ErfImp(x, true);
} }
/// <summary> /// <summary>
@ -158,7 +158,7 @@ namespace MathNet.Numerics
double[] n = new double[] { 0.00337916709551257388990745, -0.00073695653048167948530905, -0.374732337392919607868241, 0.0817442448733587196071743, -0.0421089319936548595203468, 0.0070165709512095756344528, -0.00495091255982435110337458, 0.000871646599037922480317225 }; double[] n = new double[] { 0.00337916709551257388990745, -0.00073695653048167948530905, -0.374732337392919607868241, 0.0817442448733587196071743, -0.0421089319936548595203468, 0.0070165709512095756344528, -0.00495091255982435110337458, 0.000871646599037922480317225 };
double[] d = new double[] { 1, -0.218088218087924645390535, 0.412542972725442099083918, -0.0841891147873106755410271, 0.0655338856400241519690695, -0.0120019604454941768171266, 0.00408165558926174048329689, -0.000615900721557769691924509 }; double[] d = new double[] { 1, -0.218088218087924645390535, 0.412542972725442099083918, -0.0841891147873106755410271, 0.0655338856400241519690695, -0.0120019604454941768171266, 0.00408165558926174048329689, -0.000615900721557769691924509 };
result = (z * 1.125) + (z * evaluate_polynomial(n, z) / evaluate_polynomial(d, z)); result = (z * 1.125) + (z * EvaluatePolynomial(n, z) / EvaluatePolynomial(d, z));
} }
} }
else if ((z < 110) || ((z < 110) && invert)) else if ((z < 110) || ((z < 110) && invert))
@ -173,15 +173,15 @@ namespace MathNet.Numerics
// Worst case absolute error found: 5.582813374e-21 // Worst case absolute error found: 5.582813374e-21
double[] n = new double[] { -0.0361790390718262471360258, 0.292251883444882683221149, 0.281447041797604512774415, 0.125610208862766947294894, 0.0274135028268930549240776, 0.00250839672168065762786937 }; double[] n = new double[] { -0.0361790390718262471360258, 0.292251883444882683221149, 0.281447041797604512774415, 0.125610208862766947294894, 0.0274135028268930549240776, 0.00250839672168065762786937 };
double[] d = new double[] { 1, 1.8545005897903486499845, 1.43575803037831418074962, 0.582827658753036572454135, 0.124810476932949746447682, 0.0113724176546353285778481 }; double[] d = new double[] { 1, 1.8545005897903486499845, 1.43575803037831418074962, 0.582827658753036572454135, 0.124810476932949746447682, 0.0113724176546353285778481 };
r = evaluate_polynomial(n, z - 0.5) / evaluate_polynomial(d, z - 0.5); r = EvaluatePolynomial(n, z - 0.5) / EvaluatePolynomial(d, z - 0.5);
b = 0.3440242112F; b = 0.3440242112F;
} }
else if (z < 1.25) else if (z < 1.25)
{ {
// Worst case absolute error found: 4.01854729e-21 // Worst case absolute error found: 4.01854729e-21
double[] n = new double[] { -0.0397876892611136856954425, 0.153165212467878293257683, 0.191260295600936245503129, 0.10276327061989304213645, 0.029637090615738836726027, 0.0046093486780275489468812, 0.000307607820348680180548455 }; double[] n = new double[] { -0.0397876892611136856954425, 0.153165212467878293257683, 0.191260295600936245503129, 0.10276327061989304213645, 0.029637090615738836726027, 0.0046093486780275489468812, 0.000307607820348680180548455 };
double[] d = new double[] { 1, 1.95520072987627704987886, 1.64762317199384860109595, 0.768238607022126250082483, 0.209793185936509782784315, 0.0319569316899913392596356, 0.00213363160895785378615014 }; double[] d = new double[] { 1, 1.95520072987627704987886, 1.64762317199384860109595, 0.768238607022126250082483, 0.209793185936509782784315, 0.0319569316899913392596356, 0.00213363160895785378615014 };
r = evaluate_polynomial(n, z - 0.75) / evaluate_polynomial(d, z - 0.75); r = EvaluatePolynomial(n, z - 0.75) / EvaluatePolynomial(d, z - 0.75);
b = 0.419990927F; b = 0.419990927F;
} }
else if (z < 2.25) else if (z < 2.25)
@ -189,15 +189,15 @@ namespace MathNet.Numerics
// Worst case absolute error found: 2.866005373e-21 // Worst case absolute error found: 2.866005373e-21
double[] n = new double[] { -0.0300838560557949717328341, 0.0538578829844454508530552, 0.0726211541651914182692959, 0.0367628469888049348429018, 0.00964629015572527529605267, 0.00133453480075291076745275, 0.778087599782504251917881e-4 }; double[] n = new double[] { -0.0300838560557949717328341, 0.0538578829844454508530552, 0.0726211541651914182692959, 0.0367628469888049348429018, 0.00964629015572527529605267, 0.00133453480075291076745275, 0.778087599782504251917881e-4 };
double[] d = new double[] { 1, 1.75967098147167528287343, 1.32883571437961120556307, 0.552528596508757581287907, 0.133793056941332861912279, 0.0179509645176280768640766, 0.00104712440019937356634038, -0.106640381820357337177643e-7 }; double[] d = new double[] { 1, 1.75967098147167528287343, 1.32883571437961120556307, 0.552528596508757581287907, 0.133793056941332861912279, 0.0179509645176280768640766, 0.00104712440019937356634038, -0.106640381820357337177643e-7 };
r = evaluate_polynomial(n, z - 1.25) / evaluate_polynomial(d, z - 1.25); r = EvaluatePolynomial(n, z - 1.25) / EvaluatePolynomial(d, z - 1.25);
b = 0.4898625016F; ; b = 0.4898625016F;
} }
else if (z < 3.5) else if (z < 3.5)
{ {
// Worst case absolute error found: 1.045355789e-21 // Worst case absolute error found: 1.045355789e-21
double[] n = new double[] { -0.0117907570137227847827732, 0.014262132090538809896674, 0.0202234435902960820020765, 0.00930668299990432009042239, 0.00213357802422065994322516, 0.00025022987386460102395382, 0.120534912219588189822126e-4 }; double[] n = new double[] { -0.0117907570137227847827732, 0.014262132090538809896674, 0.0202234435902960820020765, 0.00930668299990432009042239, 0.00213357802422065994322516, 0.00025022987386460102395382, 0.120534912219588189822126e-4 };
double[] d = new double[] { 1, 1.50376225203620482047419, 0.965397786204462896346934, 0.339265230476796681555511, 0.0689740649541569716897427, 0.00771060262491768307365526, 0.000371421101531069302990367 }; double[] d = new double[] { 1, 1.50376225203620482047419, 0.965397786204462896346934, 0.339265230476796681555511, 0.0689740649541569716897427, 0.00771060262491768307365526, 0.000371421101531069302990367 };
r = evaluate_polynomial(n, z - 2.25) / evaluate_polynomial(d, z - 2.25); r = EvaluatePolynomial(n, z - 2.25) / EvaluatePolynomial(d, z - 2.25);
b = 0.5317370892F; b = 0.5317370892F;
} }
else if (z < 5.25) else if (z < 5.25)
@ -205,7 +205,7 @@ namespace MathNet.Numerics
// Worst case absolute error found: 8.300028706e-22 // Worst case absolute error found: 8.300028706e-22
double[] n = new double[] { -0.00546954795538729307482955, 0.00404190278731707110245394, 0.0054963369553161170521356, 0.00212616472603945399437862, 0.000394984014495083900689956, 0.365565477064442377259271e-4, 0.135485897109932323253786e-5 }; double[] n = new double[] { -0.00546954795538729307482955, 0.00404190278731707110245394, 0.0054963369553161170521356, 0.00212616472603945399437862, 0.000394984014495083900689956, 0.365565477064442377259271e-4, 0.135485897109932323253786e-5 };
double[] d = new double[] { 1, 1.21019697773630784832251, 0.620914668221143886601045, 0.173038430661142762569515, 0.0276550813773432047594539, 0.00240625974424309709745382, 0.891811817251336577241006e-4, -0.465528836283382684461025e-11 }; double[] d = new double[] { 1, 1.21019697773630784832251, 0.620914668221143886601045, 0.173038430661142762569515, 0.0276550813773432047594539, 0.00240625974424309709745382, 0.891811817251336577241006e-4, -0.465528836283382684461025e-11 };
r = evaluate_polynomial(n, z - 3.5) / evaluate_polynomial(d, z - 3.5); r = EvaluatePolynomial(n, z - 3.5) / EvaluatePolynomial(d, z - 3.5);
b = 0.5489973426F; b = 0.5489973426F;
} }
else if (z < 8) else if (z < 8)
@ -213,23 +213,23 @@ namespace MathNet.Numerics
// Worst case absolute error found: 1.700157534e-21 // Worst case absolute error found: 1.700157534e-21
double[] n = new double[] { -0.00270722535905778347999196, 0.0013187563425029400461378, 0.00119925933261002333923989, 0.00027849619811344664248235, 0.267822988218331849989363e-4, 0.923043672315028197865066e-6 }; double[] n = new double[] { -0.00270722535905778347999196, 0.0013187563425029400461378, 0.00119925933261002333923989, 0.00027849619811344664248235, 0.267822988218331849989363e-4, 0.923043672315028197865066e-6 };
double[] d = new double[] { 1, 0.814632808543141591118279, 0.268901665856299542168425, 0.0449877216103041118694989, 0.00381759663320248459168994, 0.000131571897888596914350697, 0.404815359675764138445257e-11 }; double[] d = new double[] { 1, 0.814632808543141591118279, 0.268901665856299542168425, 0.0449877216103041118694989, 0.00381759663320248459168994, 0.000131571897888596914350697, 0.404815359675764138445257e-11 };
r = evaluate_polynomial(n, z - 5.25) / evaluate_polynomial(d, z - 5.25); r = EvaluatePolynomial(n, z - 5.25) / EvaluatePolynomial(d, z - 5.25);
b = 0.5571740866F; b = 0.5571740866F;
} }
else if (z < 11.5) else if (z < 11.5)
{ {
//Worst case absolute error found: 3.002278011e-22 // Worst case absolute error found: 3.002278011e-22
double[] n = new double[] { -0.00109946720691742196814323, 0.000406425442750422675169153, 0.000274499489416900707787024, 0.465293770646659383436343e-4, 0.320955425395767463401993e-5, 0.778286018145020892261936e-7 }; double[] n = new double[] { -0.00109946720691742196814323, 0.000406425442750422675169153, 0.000274499489416900707787024, 0.465293770646659383436343e-4, 0.320955425395767463401993e-5, 0.778286018145020892261936e-7 };
double[] d = new double[] { 1, 0.588173710611846046373373, 0.139363331289409746077541, 0.0166329340417083678763028, 0.00100023921310234908642639, 0.24254837521587225125068e-4 }; double[] d = new double[] { 1, 0.588173710611846046373373, 0.139363331289409746077541, 0.0166329340417083678763028, 0.00100023921310234908642639, 0.24254837521587225125068e-4 };
r = evaluate_polynomial(n, z - 8) / evaluate_polynomial(d, z - 8); r = EvaluatePolynomial(n, z - 8) / EvaluatePolynomial(d, z - 8);
b = 0.5609807968F; b = 0.5609807968F;
} }
else if (z < 17) else if (z < 17)
{ {
//Worst case absolute error found: 6.741114695e-21 // Worst case absolute error found: 6.741114695e-21
double[] n = new double[] { -0.00056907993601094962855594, 0.000169498540373762264416984, 0.518472354581100890120501e-4, 0.382819312231928859704678e-5, 0.824989931281894431781794e-7 }; double[] n = new double[] { -0.00056907993601094962855594, 0.000169498540373762264416984, 0.518472354581100890120501e-4, 0.382819312231928859704678e-5, 0.824989931281894431781794e-7 };
double[] d = new double[] { 1, 0.339637250051139347430323, 0.043472647870310663055044, 0.00248549335224637114641629, 0.535633305337152900549536e-4, -0.117490944405459578783846e-12 }; double[] d = new double[] { 1, 0.339637250051139347430323, 0.043472647870310663055044, 0.00248549335224637114641629, 0.535633305337152900549536e-4, -0.117490944405459578783846e-12 };
r = evaluate_polynomial(n, z - 11.5) / evaluate_polynomial(d, z - 11.5); r = EvaluatePolynomial(n, z - 11.5) / EvaluatePolynomial(d, z - 11.5);
b = 0.5626493692F; b = 0.5626493692F;
} }
else if (z < 24) else if (z < 24)
@ -237,7 +237,7 @@ namespace MathNet.Numerics
// Worst case absolute error found: 7.802346984e-22 // Worst case absolute error found: 7.802346984e-22
double[] n = new double[] { -0.000241313599483991337479091, 0.574224975202501512365975e-4, 0.115998962927383778460557e-4, 0.581762134402593739370875e-6, 0.853971555085673614607418e-8 }; double[] n = new double[] { -0.000241313599483991337479091, 0.574224975202501512365975e-4, 0.115998962927383778460557e-4, 0.581762134402593739370875e-6, 0.853971555085673614607418e-8 };
double[] d = new double[] { 1, 0.233044138299687841018015, 0.0204186940546440312625597, 0.000797185647564398289151125, 0.117019281670172327758019e-4 }; double[] d = new double[] { 1, 0.233044138299687841018015, 0.0204186940546440312625597, 0.000797185647564398289151125, 0.117019281670172327758019e-4 };
r = evaluate_polynomial(n, z - 17) / evaluate_polynomial(d, z - 17); r = EvaluatePolynomial(n, z - 17) / EvaluatePolynomial(d, z - 17);
b = 0.5634598136F; b = 0.5634598136F;
} }
else if (z < 38) else if (z < 38)
@ -245,7 +245,7 @@ namespace MathNet.Numerics
// Worst case absolute error found: 2.414228989e-22 // Worst case absolute error found: 2.414228989e-22
double[] n = new double[] { -0.000146674699277760365803642, 0.162666552112280519955647e-4, 0.269116248509165239294897e-5, 0.979584479468091935086972e-7, 0.101994647625723465722285e-8 }; double[] n = new double[] { -0.000146674699277760365803642, 0.162666552112280519955647e-4, 0.269116248509165239294897e-5, 0.979584479468091935086972e-7, 0.101994647625723465722285e-8 };
double[] d = new double[] { 1, 0.165907812944847226546036, 0.0103361716191505884359634, 0.000286593026373868366935721, 0.298401570840900340874568e-5 }; double[] d = new double[] { 1, 0.165907812944847226546036, 0.0103361716191505884359634, 0.000286593026373868366935721, 0.298401570840900340874568e-5 };
r = evaluate_polynomial(n, z - 24) / evaluate_polynomial(d, z - 24); r = EvaluatePolynomial(n, z - 24) / EvaluatePolynomial(d, z - 24);
b = 0.5638477802F; b = 0.5638477802F;
} }
else if (z < 60) else if (z < 60)
@ -253,7 +253,7 @@ namespace MathNet.Numerics
// Worst case absolute error found: 5.896543869e-24 // Worst case absolute error found: 5.896543869e-24
double[] n = new double[] { -0.583905797629771786720406e-4, 0.412510325105496173512992e-5, 0.431790922420250949096906e-6, 0.993365155590013193345569e-8, 0.653480510020104699270084e-10 }; double[] n = new double[] { -0.583905797629771786720406e-4, 0.412510325105496173512992e-5, 0.431790922420250949096906e-6, 0.993365155590013193345569e-8, 0.653480510020104699270084e-10 };
double[] d = new double[] { 1, 0.105077086072039915406159, 0.00414278428675475620830226, 0.726338754644523769144108e-4, 0.477818471047398785369849e-6 }; double[] d = new double[] { 1, 0.105077086072039915406159, 0.00414278428675475620830226, 0.726338754644523769144108e-4, 0.477818471047398785369849e-6 };
r = evaluate_polynomial(n, z - 38) / evaluate_polynomial(d, z - 38); r = EvaluatePolynomial(n, z - 38) / EvaluatePolynomial(d, z - 38);
b = 0.5640528202F; b = 0.5640528202F;
} }
else if (z < 85) else if (z < 85)
@ -261,7 +261,7 @@ namespace MathNet.Numerics
// Worst case absolute error found: 3.080612264e-21 // Worst case absolute error found: 3.080612264e-21
double[] n = new double[] { -0.196457797609229579459841e-4, 0.157243887666800692441195e-5, 0.543902511192700878690335e-7, 0.317472492369117710852685e-9 }; double[] n = new double[] { -0.196457797609229579459841e-4, 0.157243887666800692441195e-5, 0.543902511192700878690335e-7, 0.317472492369117710852685e-9 };
double[] d = new double[] { 1, 0.052803989240957632204885, 0.000926876069151753290378112, 0.541011723226630257077328e-5, 0.535093845803642394908747e-15 }; double[] d = new double[] { 1, 0.052803989240957632204885, 0.000926876069151753290378112, 0.541011723226630257077328e-5, 0.535093845803642394908747e-15 };
r = evaluate_polynomial(n, z - 60) / evaluate_polynomial(d, z - 60); r = EvaluatePolynomial(n, z - 60) / EvaluatePolynomial(d, z - 60);
b = 0.5641309023F; b = 0.5641309023F;
} }
else else
@ -269,11 +269,11 @@ namespace MathNet.Numerics
// Worst case absolute error found: 8.094633491e-22 // Worst case absolute error found: 8.094633491e-22
double[] n = new double[] { -0.789224703978722689089794e-5, 0.622088451660986955124162e-6, 0.145728445676882396797184e-7, 0.603715505542715364529243e-10 }; double[] n = new double[] { -0.789224703978722689089794e-5, 0.622088451660986955124162e-6, 0.145728445676882396797184e-7, 0.603715505542715364529243e-10 };
double[] d = new double[] { 1, 0.0375328846356293715248719, 0.000467919535974625308126054, 0.193847039275845656900547e-5 }; double[] d = new double[] { 1, 0.0375328846356293715248719, 0.000467919535974625308126054, 0.193847039275845656900547e-5 };
r = evaluate_polynomial(n, z - 85) / evaluate_polynomial(d, z - 85); r = EvaluatePolynomial(n, z - 85) / EvaluatePolynomial(d, z - 85);
b = 0.5641584396F; b = 0.5641584396F;
} }
double g = System.Math.Exp(-z * z) / z; double g = Math.Exp(-z * z) / z;
result = (g * b) + (g * r); result = (g * b) + (g * r);
} }
else else
@ -293,18 +293,18 @@ namespace MathNet.Numerics
return result; return result;
} }
///<summary>Calculates the complementary inverse error function evaluated at z.</summary> /// <summary>Calculates the complementary inverse error function evaluated at z.</summary>
/// <returns>The complementary inverse error function evaluated at given value.</returns> /// <returns>The complementary inverse error function evaluated at given value.</returns>
/// <remarks> We have tested this implementation against the arbitrary precision mpmath library /// <remarks> We have tested this implementation against the arbitrary precision mpmath library
/// and found cases where we can only guarantee 9 significant figures correct. /// and found cases where we can only guarantee 9 significant figures correct.
/// <list type="bullet"> /// <list type="bullet">
/// <item>returns Double.PositiveInfinity if <c>z &lt;= 0.0</c>.</item> /// <item>returns Double.PositiveInfinity if <c>z &lt;= 0.0</c>.</item>
/// <item>returns Double.NegativeInfinity if <c>z &gt;= 2.0</c>.</item> /// <item>returns Double.NegativeInfinity if <c>z &gt;= 2.0</c>.</item>
/// </list> /// </list>
/// </remarks> /// </remarks>
///<summary>Calculates the complementary inverse error function evaluated at z.</summary> /// <summary>calculates the complementary inverse error function evaluated at z.</summary>
///<param name="z">value to evaluate.</param> /// <param name="z">value to evaluate.</param>
///<returns>the complementary inverse error function evaluated at Z.</returns> /// <returns>the complementary inverse error function evaluated at Z.</returns>
public static double ErfcInv(double z) public static double ErfcInv(double z)
{ {
if (z <= 0.0) if (z <= 0.0)
@ -363,7 +363,7 @@ namespace MathNet.Numerics
double[] P = new double[] { -0.000508781949658280665617, -0.00836874819741736770379, 0.0334806625409744615033, -0.0126926147662974029034, -0.0365637971411762664006, 0.0219878681111168899165, 0.00822687874676915743155, -0.00538772965071242932965 }; double[] P = new double[] { -0.000508781949658280665617, -0.00836874819741736770379, 0.0334806625409744615033, -0.0126926147662974029034, -0.0365637971411762664006, 0.0219878681111168899165, 0.00822687874676915743155, -0.00538772965071242932965 };
double[] Q = new double[] { 1, -0.970005043303290640362, -1.56574558234175846809, 1.56221558398423026363, 0.662328840472002992063, -0.71228902341542847553, -0.0527396382340099713954, 0.0795283687341571680018, -0.00233393759374190016776, 0.000886216390456424707504 }; double[] Q = new double[] { 1, -0.970005043303290640362, -1.56574558234175846809, 1.56221558398423026363, 0.662328840472002992063, -0.71228902341542847553, -0.0527396382340099713954, 0.0795283687341571680018, -0.00233393759374190016776, 0.000886216390456424707504 };
double g = p * (p + 10); double g = p * (p + 10);
double r = evaluate_polynomial(P, p) / evaluate_polynomial(Q, p); double r = EvaluatePolynomial(P, p) / EvaluatePolynomial(Q, p);
result = (g * Y) + (g * r); result = (g * Y) + (g * r);
} }
else if (q >= 0.25) else if (q >= 0.25)
@ -385,7 +385,7 @@ namespace MathNet.Numerics
double[] Q = new double[] { 1, 6.24264124854247537712, 3.9713437953343869095, -28.6608180499800029974, -20.1432634680485188801, 48.5609213108739935468, 10.8268667355460159008, -22.6436933413139721736, 1.72114765761200282724 }; double[] Q = new double[] { 1, 6.24264124854247537712, 3.9713437953343869095, -28.6608180499800029974, -20.1432634680485188801, 48.5609213108739935468, 10.8268667355460159008, -22.6436933413139721736, 1.72114765761200282724 };
double g = System.Math.Sqrt(-2 * System.Math.Log(q)); double g = System.Math.Sqrt(-2 * System.Math.Log(q));
double xs = q - 0.25; double xs = q - 0.25;
double r = evaluate_polynomial(P, xs) / evaluate_polynomial(Q, xs); double r = EvaluatePolynomial(P, xs) / EvaluatePolynomial(Q, xs);
result = g / (Y + r); result = g / (Y + r);
} }
else else
@ -417,7 +417,7 @@ namespace MathNet.Numerics
double[] P = new double[] { -0.131102781679951906451, -0.163794047193317060787, 0.117030156341995252019, 0.387079738972604337464, 0.337785538912035898924, 0.142869534408157156766, 0.0290157910005329060432, 0.00214558995388805277169, -0.679465575181126350155e-6, 0.285225331782217055858e-7, -0.681149956853776992068e-9 }; double[] P = new double[] { -0.131102781679951906451, -0.163794047193317060787, 0.117030156341995252019, 0.387079738972604337464, 0.337785538912035898924, 0.142869534408157156766, 0.0290157910005329060432, 0.00214558995388805277169, -0.679465575181126350155e-6, 0.285225331782217055858e-7, -0.681149956853776992068e-9 };
double[] Q = new double[] { 1, 3.46625407242567245975, 5.38168345707006855425, 4.77846592945843778382, 2.59301921623620271374, 0.848854343457902036425, 0.152264338295331783612, 0.01105924229346489121 }; double[] Q = new double[] { 1, 3.46625407242567245975, 5.38168345707006855425, 4.77846592945843778382, 2.59301921623620271374, 0.848854343457902036425, 0.152264338295331783612, 0.01105924229346489121 };
double xs = x - 1.125; double xs = x - 1.125;
double R = evaluate_polynomial(P, xs) / evaluate_polynomial(Q, xs); double R = EvaluatePolynomial(P, xs) / EvaluatePolynomial(Q, xs);
result = (Y * x) + (R * x); result = (Y * x) + (R * x);
} }
else if (x < 6) else if (x < 6)
@ -427,7 +427,7 @@ namespace MathNet.Numerics
double[] P = new double[] { -0.0350353787183177984712, -0.00222426529213447927281, 0.0185573306514231072324, 0.00950804701325919603619, 0.00187123492819559223345, 0.000157544617424960554631, 0.460469890584317994083e-5, -0.230404776911882601748e-9, 0.266339227425782031962e-11 }; double[] P = new double[] { -0.0350353787183177984712, -0.00222426529213447927281, 0.0185573306514231072324, 0.00950804701325919603619, 0.00187123492819559223345, 0.000157544617424960554631, 0.460469890584317994083e-5, -0.230404776911882601748e-9, 0.266339227425782031962e-11 };
double[] Q = new double[] { 1, 1.3653349817554063097, 0.762059164553623404043, 0.220091105764131249824, 0.0341589143670947727934, 0.00263861676657015992959, 0.764675292302794483503e-4 }; double[] Q = new double[] { 1, 1.3653349817554063097, 0.762059164553623404043, 0.220091105764131249824, 0.0341589143670947727934, 0.00263861676657015992959, 0.764675292302794483503e-4 };
double xs = x - 3; double xs = x - 3;
double R = evaluate_polynomial(P, xs) / evaluate_polynomial(Q, xs); double R = EvaluatePolynomial(P, xs) / EvaluatePolynomial(Q, xs);
result = (Y * x) + (R * x); result = (Y * x) + (R * x);
} }
else if (x < 18) else if (x < 18)
@ -437,7 +437,7 @@ namespace MathNet.Numerics
double[] P = new double[] { -0.0167431005076633737133, -0.00112951438745580278863, 0.00105628862152492910091, 0.000209386317487588078668, 0.149624783758342370182e-4, 0.449696789927706453732e-6, 0.462596163522878599135e-8, -0.281128735628831791805e-13, 0.99055709973310326855e-16 }; double[] P = new double[] { -0.0167431005076633737133, -0.00112951438745580278863, 0.00105628862152492910091, 0.000209386317487588078668, 0.149624783758342370182e-4, 0.449696789927706453732e-6, 0.462596163522878599135e-8, -0.281128735628831791805e-13, 0.99055709973310326855e-16 };
double[] Q = new double[] { 1, 0.591429344886417493481, 0.138151865749083321638, 0.0160746087093676504695, 0.000964011807005165528527, 0.275335474764726041141e-4, 0.282243172016108031869e-6 }; double[] Q = new double[] { 1, 0.591429344886417493481, 0.138151865749083321638, 0.0160746087093676504695, 0.000964011807005165528527, 0.275335474764726041141e-4, 0.282243172016108031869e-6 };
double xs = x - 6; double xs = x - 6;
double R = evaluate_polynomial(P, xs) / evaluate_polynomial(Q, xs); double R = EvaluatePolynomial(P, xs) / EvaluatePolynomial(Q, xs);
result = (Y * x) + (R * x); result = (Y * x) + (R * x);
} }
else if (x < 44) else if (x < 44)
@ -447,7 +447,7 @@ namespace MathNet.Numerics
double[] P = new double[] { -0.0024978212791898131227, -0.779190719229053954292e-5, 0.254723037413027451751e-4, 0.162397777342510920873e-5, 0.396341011304801168516e-7, 0.411632831190944208473e-9, 0.145596286718675035587e-11, -0.116765012397184275695e-17 }; double[] P = new double[] { -0.0024978212791898131227, -0.779190719229053954292e-5, 0.254723037413027451751e-4, 0.162397777342510920873e-5, 0.396341011304801168516e-7, 0.411632831190944208473e-9, 0.145596286718675035587e-11, -0.116765012397184275695e-17 };
double[] Q = new double[] { 1, 0.207123112214422517181, 0.0169410838120975906478, 0.000690538265622684595676, 0.145007359818232637924e-4, 0.144437756628144157666e-6, 0.509761276599778486139e-9 }; double[] Q = new double[] { 1, 0.207123112214422517181, 0.0169410838120975906478, 0.000690538265622684595676, 0.145007359818232637924e-4, 0.144437756628144157666e-6, 0.509761276599778486139e-9 };
double xs = x - 18; double xs = x - 18;
double R = evaluate_polynomial(P, xs) / evaluate_polynomial(Q, xs); double R = EvaluatePolynomial(P, xs) / EvaluatePolynomial(Q, xs);
result = (Y * x) + (R * x); result = (Y * x) + (R * x);
} }
else else
@ -457,7 +457,7 @@ namespace MathNet.Numerics
double[] P = new double[] { -0.000539042911019078575891, -0.28398759004727721098e-6, 0.899465114892291446442e-6, 0.229345859265920864296e-7, 0.225561444863500149219e-9, 0.947846627503022684216e-12, 0.135880130108924861008e-14, -0.348890393399948882918e-21 }; double[] P = new double[] { -0.000539042911019078575891, -0.28398759004727721098e-6, 0.899465114892291446442e-6, 0.229345859265920864296e-7, 0.225561444863500149219e-9, 0.947846627503022684216e-12, 0.135880130108924861008e-14, -0.348890393399948882918e-21 };
double[] Q = new double[] { 1, 0.0845746234001899436914, 0.00282092984726264681981, 0.468292921940894236786e-4, 0.399968812193862100054e-6, 0.161809290887904476097e-8, 0.231558608310259605225e-11 }; double[] Q = new double[] { 1, 0.0845746234001899436914, 0.00282092984726264681981, 0.468292921940894236786e-4, 0.399968812193862100054e-6, 0.161809290887904476097e-8, 0.231558608310259605225e-11 };
double xs = x - 44; double xs = x - 44;
double R = evaluate_polynomial(P, xs) / evaluate_polynomial(Q, xs); double R = EvaluatePolynomial(P, xs) / EvaluatePolynomial(Q, xs);
result = (Y * x) + (R * x); result = (Y * x) + (R * x);
} }
} }
@ -471,7 +471,7 @@ namespace MathNet.Numerics
/// <param name="poly">The coefficients of the polynomial.</param> /// <param name="poly">The coefficients of the polynomial.</param>
/// <param name="z">The location where to evaluate the polynomial at.</param> /// <param name="z">The location where to evaluate the polynomial at.</param>
/// <returns>the evaluation of the polynomial.</returns> /// <returns>the evaluation of the polynomial.</returns>
private static double evaluate_polynomial(double[] poly, double z) private static double EvaluatePolynomial(double[] poly, double z)
{ {
int count = poly.Length; int count = poly.Length;
double sum = poly[count - 1]; double sum = poly[count - 1];

Write Preview
Loading…
Cancel
Save