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))); + } + } +}