From ccf89f8186cdd114f54bac51d8a1e848769ed51f Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Sat, 4 May 2013 13:47:31 +0200 Subject: [PATCH] RootFinding: protect newton-raphson better against singularities, more tests --- src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs | 4 ++-- src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs | 7 +++++++ 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs b/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs index 3fd9bd82..11093b1a 100644 --- a/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs +++ b/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs @@ -129,7 +129,7 @@ namespace MathNet.Numerics.RootFinding.Algorithms continue; } - if (Math.Abs(step) < accuracy) + if (Math.Abs(step) < accuracy && Math.Abs(fx) < accuracy) { return true; } @@ -149,7 +149,7 @@ namespace MathNet.Numerics.RootFinding.Algorithms lowerBound = root; fmin = fx; } - else if (Math.Sign(fmin) != Math.Sign(fmax)) + else if (Math.Sign(fmin) != Math.Sign(fmax) && Math.Abs(fx) < accuracy) { return true; } diff --git a/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs b/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs index 808cd76d..a045842e 100644 --- a/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs +++ b/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs @@ -91,6 +91,13 @@ namespace MathNet.Numerics.UnitTests.RootFindingTests Assert.AreEqual(Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f3, df3, 1, 1.99, 1e-14, 100, 20)); Assert.AreEqual(Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f3, df3, -1.5, 1.99, 1e-14, 100, 20)); Assert.AreEqual(Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f3, df3, 1, 6, 1e-14, 100, 20)); + + Func f4 = x => 1/(2 - x) - x + 6; + Func df4 = x => 1/(x*x - 4*x + 4) - 1; + Assert.AreEqual(4 + Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f4, df4, 5, 6, 1e-14, 100, 20), 1e-14); + Assert.AreEqual(4 - Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f4, df4, 2.01, 3, 1e-14, 100, 20)); + Assert.AreEqual(4 - Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f4, df4, 2.01, 5, 1e-14, 100, 20)); + Assert.AreEqual(4 - Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f4, df4, -2, 4, 1e-14, 100, 20)); } [Test]