From d512e445653825617fda768404ece6df427c32b5 Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Mon, 2 Nov 2026 13:29:57 +0000 Subject: [PATCH] Halve a sign change down to the root before Newton's method polishes it Newton's method, which Solve falls back to, also starts from each bracket where the expression changes sign along the real line, and it started from the bracket's midpoint. (x^2 - 3x - 3)^(1/3) + sqrt(3x - 4) = 3 has its one real root at 3.79617, just past the cube root's zero at 3.79129, where the slope is unbounded. From the midpoint 3.7 the iteration left the real line, the root was never answered, and over the reals the answer was { }, which claims there is none. The bracket is now halved down to the root in double precision, and Newton starts from there. A pole between the two ends is halved down to as well, and its value rejects it as before: 1/x = 0 is still { }. radcheck, 400 random radical equations at seed 51: 1 missed root on master, 0 here; seeds 22 and 53 (pow) are the same on both. Co-Authored-By: Claude Opus 5.5 --- BREAKING-CHANGES.md | 13 +++++ .../Solvers/NumericalSolving/NewtonSolver.cs | 34 ++++++++++++-- .../SolveTest/ARootPastARadicandsZeroTest.cs | 47 +++++++++++++++++++ 3 files changed, 91 insertions(+), 3 deletions(-) create mode 100644 Sources/Tests/UnitTests/Algebra/SolveTest/ARootPastARadicandsZeroTest.cs diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index ed8e67ccc..c7451f5b4 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -105,6 +105,19 @@ root. An equation in a square root whose squared form is such a quartic was answ | `"31x^4 - 62x^3 + 33x^2 - 2x - 1 = 0".ToEntity().Solve("x")` | `{ }` | `{ 1/2 + -sqrt((27/62 + -8/31 * sqrt(2)) / 2), … }`, four real roots | | `"x^2 - x - 1 = 4x*sqrt(2x^2 - 4x + 2)".ToEntity().Solve("x")` | `{ }` | `{ 1/2 + -sqrt((27/62 + 8/31 * sqrt(2)) / 2) }` | +### A real root just past where a radical's argument changes sign is found + +**Silent: an equation Newton's method answers gains a real root it missed, and over the reals is no longer answered +`{ }`.** Newton's method, which `Solve` falls back to, also starts from each bracket where the expression changes sign +along the real line, and started from the bracket's midpoint. Where a radical's argument crosses zero inside the +bracket the slope there is unbounded, and from the midpoint the iteration left the real line and never came back, so a +real root sitting just past that zero was not answered. The bracket is halved down to the root first. + +| Input | Was (2.6.0) | Now | +|---|---|---| +| `"(x^2 - 3x - 3)^(1/3) + sqrt(3x - 4) = 3".ToEntity().Solve("x")` | `{ 1.3099 - 2.3063i, 1.3099 + 2.3063i }` | also `3.7961670056822383` | +| `"(x^2 - 3x - 3)^(1/3) + sqrt(3x - 4) = 3 provided x in RR".ToEntity().Solve("x")` | `{ }` | `{ 3.7961670056822383 }` | + ### A quotient of logarithms is differentiated only where both diverge **Silent: a limit of a quotient of logarithms whose divisor tends to 0 is no longer answered by l'Hopital's rule, which diff --git a/Sources/AngouriMath/Functions/Continuous/Solvers/NumericalSolving/NewtonSolver.cs b/Sources/AngouriMath/Functions/Continuous/Solvers/NumericalSolving/NewtonSolver.cs index 72c218789..3719a1555 100644 --- a/Sources/AngouriMath/Functions/Continuous/Solvers/NumericalSolving/NewtonSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Solvers/NumericalSolving/NewtonSolver.cs @@ -207,14 +207,42 @@ private static IEnumerable RealSignChanges( yield return new NumericsComplex(at, 0); else if (value is { } here && previous is { } there && (here < 0) != (there < 0)) - // The midpoint, rather than either end: the bracket holds a root and - // Newton converges fastest from as near it as we can say it is. - yield return new NumericsComplex((at + previousAt) / 2, 0); + // The bracket halved down to the root, rather than its midpoint: the + // bracket holds a root and Newton converges fastest from as near it as we + // can say it is. From the midpoint it could not converge at all where the + // slope is unbounded inside the bracket -- (x^2 - 3x - 3)^(1/3) + // + sqrt(3x - 4) = 3 has its root at 3.7962, past the cube root's zero at + // 3.7913, and from 3.7 the iteration went off the real line and its one + // real root was never answered. A pole between the two is halved down to as + // well, and its value there rejects it as it did before. + yield return new NumericsComplex(Halved(previousAt, there, at, RealValueAt), 0); previousAt = at; previous = value; } } + /// + /// A point of [low, high], or [high, low], as near a sign change of + /// as halving the interval gets in double precision, given + /// the value at . Stops where a midpoint has no real value. + /// + private static double Halved(double low, double atLow, double high, Func valueAt) + { + for (var step = 0; step < 64; step++) + { + var middle = (low + high) / 2; + if (middle == low || middle == high || valueAt(middle) is not { } here) + return middle; + if (here == 0) + return middle; + if ((here < 0) == (atLow < 0)) + (low, atLow) = (middle, here); + else + high = middle; + } + return (low + high) / 2; + } + /// /// How far apart two candidates may be, relative to their own size, and still be /// the same root. diff --git a/Sources/Tests/UnitTests/Algebra/SolveTest/ARootPastARadicandsZeroTest.cs b/Sources/Tests/UnitTests/Algebra/SolveTest/ARootPastARadicandsZeroTest.cs new file mode 100644 index 000000000..a28e16c5b --- /dev/null +++ b/Sources/Tests/UnitTests/Algebra/SolveTest/ARootPastARadicandsZeroTest.cs @@ -0,0 +1,47 @@ +// +// Copyright (c) 2019-2026 Angouri. +// AngouriMath is licensed under MIT. +// Details: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md. +// Website: https://am.angouri.org. +// + +using System.Linq; +using AngouriMath; +using AngouriMath.Extensions; +using Xunit; +using static AngouriMath.Entity.Number; + +namespace AngouriMath.Tests.Algebra +{ + /// + /// Newton's method starts from each bracket where the expression changes sign along the real + /// line, and started from the bracket's midpoint. (x^2 - 3x - 3)^(1/3) + sqrt(3x - 4) = 3 + /// has its one real root at 3.79617, just past the cube root's zero at 3.79129, where the + /// slope is unbounded; from the midpoint 3.7 the iteration left the real line, and the root + /// was not answered, so over the reals the answer was { }. No issue was filed; found + /// by a cross-check of random radical equations against mpmath. + /// + [Trait("Area", "Algebra")] + public sealed class ARootPastARadicandsZeroTest + { + // mpmath's findroot at 40 digits. + private static readonly Real Root = Real.Create(PeterO.Numbers.EDecimal.FromString("3.7961670056822382211489100492838228303")); + + private static bool IsTheRoot(Entity candidate) + => candidate.Evaled is Complex value && ((Complex)(value - Root)).Abs().EDecimal.CompareTo(PeterO.Numbers.EDecimal.FromString("1e-12")) < 0; + + [Fact] + public void IsAnswered() + { + var roots = Assert.IsType("(x^2 - 3x - 3)^(1/3) + sqrt(3x - 4) = 3".ToEntity().Solve("x")); + Assert.Contains(roots, IsTheRoot); + } + + [Fact] + public void IsTheOneRealRoot() + { + var roots = Assert.IsType("(x^2 - 3x - 3)^(1/3) + sqrt(3x - 4) = 3 provided x in RR".ToEntity().Solve("x")); + Assert.True(IsTheRoot(Assert.Single(roots))); + } + } +}