Browse Source

Polynomial: more integration into FindRoots

ridge-regression
Christoph Ruegg 8 years ago
parent
commit
bd8742aec3
  1. 1
      MathNet.Numerics.sln.DotSettings
  2. 31
      src/Numerics.Tests/RootFindingTests/CubicTest.cs
  3. 22
      src/Numerics/FindRoots.cs
  4. 43
      src/Numerics/Polynomial.cs

1
MathNet.Numerics.sln.DotSettings

@ -81,6 +81,7 @@ OTHER DEALINGS IN THE SOFTWARE.
<s:Boolean x:Key="/Default/Environment/SettingsMigration/IsMigratorApplied/=JetBrains_002EReSharper_002EPsi_002ECSharp_002ECodeStyle_002ESettingsUpgrade_002ECSharpPlaceAttributeOnSameLineMigration/@EntryIndexedValue">True</s:Boolean>
<s:Boolean x:Key="/Default/Environment/SettingsMigration/IsMigratorApplied/=JetBrains_002EReSharper_002EPsi_002ECSharp_002ECodeStyle_002ESettingsUpgrade_002EMigrateBlankLinesAroundFieldToBlankLinesAroundProperty/@EntryIndexedValue">True</s:Boolean>
<s:Boolean x:Key="/Default/Environment/SettingsMigration/IsMigratorApplied/=JetBrains_002EReSharper_002EPsi_002ECSharp_002ECodeStyle_002ESettingsUpgrade_002EMigrateThisQualifierSettings/@EntryIndexedValue">True</s:Boolean>
<s:Boolean x:Key="/Default/UserDictionary/Words/=Chebychev/@EntryIndexedValue">True</s:Boolean>
<s:Boolean x:Key="/Default/UserDictionary/Words/=Numerics/@EntryIndexedValue">True</s:Boolean>
</wpf:ResourceDictionary>

31
src/Numerics.Tests/RootFindingTests/CubicTest.cs

@ -97,6 +97,19 @@ namespace MathNet.Numerics.UnitTests.RootFindingTests
Assert.That(roots.Item3.Imaginary, Is.EqualTo(0).Within(1e-14));
}
[TestCase(6.0, -5.0, -2.0, 1.0, 3.0, -2.0, 1.0)]
public void ComplexRoots_TripleReal_AsPolynomial(double d, double c, double b, double a, double x1, double x2, double x3)
{
var roots = FindRoots.Polynomial(new[] { d, c, b, a });
Assert.That(roots.Length, Is.EqualTo(3));
Assert.That(roots[0].Real, Is.EqualTo(x1).Within(1e-14));
Assert.That(roots[0].Imaginary, Is.EqualTo(0).Within(1e-14));
Assert.That(roots[1].Real, Is.EqualTo(x2).Within(1e-14));
Assert.That(roots[1].Imaginary, Is.EqualTo(0).Within(1e-14));
Assert.That(roots[2].Real, Is.EqualTo(x3).Within(1e-14));
Assert.That(roots[2].Imaginary, Is.EqualTo(0).Within(1e-14));
}
[TestCase(-350.0, 162.0, -30.0, 2.0)]
[TestCase(6.0, -5.0, -2.0, 1.0)]
[TestCase(1.0, 5.0, 2.0, 1.0)]
@ -113,5 +126,23 @@ namespace MathNet.Numerics.UnitTests.RootFindingTests
Assert.That(Polynomial.Evaluate(roots.Item3, d, c, b, a).Real, Is.EqualTo(0).Within(1e-12));
Assert.That(Polynomial.Evaluate(roots.Item3, d, c, b, a).Imaginary, Is.EqualTo(0).Within(1e-12));
}
[TestCase(-350.0, 162.0, -30.0, 2.0)]
[TestCase(6.0, -5.0, -2.0, 1.0)]
[TestCase(1.0, 5.0, 2.0, 1.0)]
[TestCase(1.0, 5.0, 0.0, 1.0)]
[TestCase(1.0, 0.0, 2.0, 1.0)]
[TestCase(0.0, 0.0, 0.0, 2.0)]
public void ComplexRootsAreRoots_AsPolynomial(double d, double c, double b, double a)
{
var roots = FindRoots.Polynomial(new[] { d, c, b, a });
Assert.That(roots.Length, Is.EqualTo(3));
Assert.That(Polynomial.Evaluate(roots[0], d, c, b, a).Real, Is.EqualTo(0).Within(1e-12));
Assert.That(Polynomial.Evaluate(roots[0], d, c, b, a).Imaginary, Is.EqualTo(0).Within(1e-12));
Assert.That(Polynomial.Evaluate(roots[1], d, c, b, a).Real, Is.EqualTo(0).Within(1e-12));
Assert.That(Polynomial.Evaluate(roots[1], d, c, b, a).Imaginary, Is.EqualTo(0).Within(1e-12));
Assert.That(Polynomial.Evaluate(roots[2], d, c, b, a).Real, Is.EqualTo(0).Within(1e-12));
Assert.That(Polynomial.Evaluate(roots[2], d, c, b, a).Imaginary, Is.EqualTo(0).Within(1e-12));
}
}
}

22
src/Numerics/FindRoots.cs

@ -114,11 +114,21 @@ namespace MathNet.Numerics
/// <summary>
/// Find all roots of a polynomial by calculating the characteristic polynomial of the companion matrix
/// </summary>
/// <param name="poly">the values for the polynomial in ascending order e.G new double[] {5, 0, 2} = "5 + 0 x^1 + 2 x^2"</param>
/// <returns>the roots of the polynomial</returns>
public static Complex[] Polynomial(double[] poly)
/// <param name="coefficients">The coefficients of the polynomial in ascending order, e.g. new double[] {5, 0, 2} = "5 + 0 x^1 + 2 x^2"</param>
/// <returns>The roots of the polynomial</returns>
public static Complex[] Polynomial(double[] coefficients)
{
return new Polynomial(poly).Roots();
return new Polynomial(coefficients).Roots();
}
/// <summary>
/// Find all roots of a polynomial by calculating the characteristic polynomial of the companion matrix
/// </summary>
/// <param name="polynomial">The polynomial.</param>
/// <returns>The roots of the polynomial</returns>
public static Complex[] Polynomial(Polynomial polynomial)
{
return polynomial.Roots();
}
/// <summary>
@ -139,7 +149,7 @@ namespace MathNet.Numerics
double location = 0.5*(intervalBegin + intervalEnd);
double scale = 0.5*(intervalEnd - intervalBegin);
// evaluate first kind chebyshev nodes
// evaluate first kind chebychev nodes
double angleFactor = Constants.Pi/(2*degree);
var samples = new double[degree];
@ -168,7 +178,7 @@ namespace MathNet.Numerics
double location = 0.5*(intervalBegin + intervalEnd);
double scale = 0.5*(intervalEnd - intervalBegin);
// evaluate second kind chebyshev nodes
// evaluate second kind chebychev nodes
double angleFactor = Constants.Pi/(degree + 1);
var samples = new double[degree];

43
src/Numerics/Polynomial.cs

@ -76,10 +76,21 @@ namespace MathNet.Numerics
/// Create a constant polynomial.
/// Example: 3.0 -> "p : x -> 3.0"
/// </summary>
/// <param name="coefficient">just the "x^0" part</param>
/// <param name="coefficient">The coefficient of the "x^0" monomial.</param>
public Polynomial(double coefficient)
{
Coefficients = coefficient == 0.0 ? new double[0] : new[] { coefficient };
if (coefficient == 0.0)
{
#if NET40
Coefficients = new double[0];
#else
Coefficients = Array.Empty<double>();
#endif
}
else
{
Coefficients = new[] { coefficient };
}
}
/// <summary>
@ -101,6 +112,17 @@ namespace MathNet.Numerics
{
}
public static Polynomial Zero => new Polynomial();
/// <summary>
/// Least-Squares fitting the points (x,y) to a k-order polynomial y : x -> p0 + p1*x + p2*x^2 + ... + pk*x^k
/// </summary>
public static Polynomial Fit(double[] x, double[] y, int order, DirectRegressionMethod method = DirectRegressionMethod.QR)
{
var coefficients = Numerics.Fit.Polynomial(x, y, order, method);
return new Polynomial(coefficients);
}
static int EvaluateDegree(double[] coefficients)
{
for (int i = coefficients.Length - 1; i >= 0; i--)
@ -114,22 +136,11 @@ namespace MathNet.Numerics
return -1;
}
public static Polynomial Zero => new Polynomial();
/// <summary>
/// Least-Squares fitting the points (x,y) to a k-order polynomial y : x -> p0 + p1*x + p2*x^2 + ... + pk*x^k
/// </summary>
public static Polynomial Fit(double[] x, double[] y, int order, DirectRegressionMethod method = DirectRegressionMethod.QR)
{
var coefficients = Numerics.Fit.Polynomial(x, y, order, method);
return new Polynomial(coefficients);
}
#region Evaluation
/// <summary>
/// Evaluate a polynomial at point x.
/// Coefficients are ordered by power with power k at index k.
/// Coefficients are ordered ascending by power with power k at index k.
/// Example: coefficients [3,-1,2] represent y=2x^2-x+3.
/// </summary>
/// <param name="z">The location where to evaluate the polynomial at.</param>
@ -148,7 +159,7 @@ namespace MathNet.Numerics
/// <summary>
/// Evaluate a polynomial at point x.
/// Coefficients are ordered by power with power k at index k.
/// Coefficients are ordered ascending by power with power k at index k.
/// Example: coefficients [3,-1,2] represent y=2x^2-x+3.
/// </summary>
/// <param name="z">The location where to evaluate the polynomial at.</param>
@ -167,7 +178,7 @@ namespace MathNet.Numerics
/// <summary>
/// Evaluate a polynomial at point x.
/// Coefficients are ordered by power with power k at index k.
/// Coefficients are ordered ascending by power with power k at index k.
/// Example: coefficients [3,-1,2] represent y=2x^2-x+3.
/// </summary>
/// <param name="z">The location where to evaluate the polynomial at.</param>

Loading…
Cancel
Save