diff --git a/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs b/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs index 0e05b4fe..1a3a1325 100644 --- a/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs +++ b/src/Numerics/RootFinding/Algorithms/HybridNewtonRaphson.cs @@ -93,6 +93,11 @@ namespace MathNet.Numerics.RootFinding.Algorithms double step = fx/dfx; root -= step; + if (Math.Abs(step) < accuracy && Math.Abs(fx) < accuracy) + { + return true; + } + bool overshoot = root > upperBound, undershoot = root < lowerBound; if (overshoot || undershoot || Math.Abs(2*fx) > Math.Abs(lastStep*dfx)) { @@ -131,11 +136,6 @@ namespace MathNet.Numerics.RootFinding.Algorithms continue; } - if (Math.Abs(step) < accuracy && Math.Abs(fx) < accuracy) - { - return true; - } - // Evaluation fx = f(root); lastStep = step; diff --git a/src/UnitTests/RootFindingTests/BrentTest.cs b/src/UnitTests/RootFindingTests/BrentTest.cs index 6dea542f..94067a1a 100644 --- a/src/UnitTests/RootFindingTests/BrentTest.cs +++ b/src/UnitTests/RootFindingTests/BrentTest.cs @@ -66,6 +66,20 @@ namespace MathNet.Numerics.UnitTests.RootFindingTests Assert.AreEqual(0, f1(Brent.FindRoot(f1, -2, 4, 1e-14, 100)), 1e-14); } + [Test] + public void Cubic() + { + // with complex roots (looking for the real root only) + Func f1 = x => 3 * x * x * x + 4 * x * x + 5 * x + 6; + Assert.AreEqual(-1.265328088928, Brent.FindRoot(f1, -2, -1, 1e-8, 100), 1e-6); + + // real roots only + Func f2 = x => 2 * x * x * x + 4 * x * x - 50 * x + 6; + Assert.AreEqual(-6.1466562197069, Brent.FindRoot(f2, -6.5, -5.5, 1e-8, 100), 1e-6); + Assert.AreEqual(0.12124737195841, Brent.FindRoot(f2, -0.5, 0.5, 1e-8, 100), 1e-6); + Assert.AreEqual(4.0254088477485, Brent.FindRoot(f2, 3.5, 4.5, 1e-8, 100), 1e-6); + } + [Test] public void NoRoot() { diff --git a/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs b/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs index 2d32f75d..9fa69da7 100644 --- a/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs +++ b/src/UnitTests/RootFindingTests/NewtonRaphsonTest.cs @@ -100,6 +100,23 @@ namespace MathNet.Numerics.UnitTests.RootFindingTests Assert.AreEqual(4 - Math.Sqrt(3), HybridNewtonRaphson.FindRoot(f4, df4, -2, 4, 1e-14, 100, 20)); } + [Test] + public void Cubic() + { + // with complex roots (looking for the real root only) + Func f1 = x => 3*x*x*x + 4*x*x + 5*x + 6; + Func df1 = x => 9*x*x + 8*x + 5; + Assert.AreEqual(-1.265328088928, HybridNewtonRaphson.FindRoot(f1, df1, -2, -1, 1e-10, 100, 20), 1e-6); + Assert.AreEqual(-1.265328088928, HybridNewtonRaphson.FindRoot(f1, df1, -5, 5, 1e-10, 100, 20), 1e-6); + + // real roots only + Func f2 = x => 2*x*x*x + 4*x*x - 50*x + 6; + Func df2 = x => 6*x*x + 8*x - 50; + Assert.AreEqual(-6.1466562197069, HybridNewtonRaphson.FindRoot(f2, df2, -8, -5, 1e-10, 100, 20), 1e-6); + Assert.AreEqual(0.12124737195841, HybridNewtonRaphson.FindRoot(f2, df2, -1, 1, 1e-10, 100, 20), 1e-6); + Assert.AreEqual(4.0254088477485, HybridNewtonRaphson.FindRoot(f2, df2, 3, 5, 1e-10, 100, 20), 1e-6); + } + [Test] public void NoRoot() {