diff --git a/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver.cs b/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver.cs
index babf49b36..6416c9ec6 100644
--- a/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver.cs
+++ b/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver.cs
@@ -31,8 +31,16 @@ internal static Set Solve(Entity equation, Variable x)
if (solutions is FiniteSet finiteSet)
{
- return finiteSet.Select(simplifier)
+ var finite = finiteSet.Select(simplifier)
.Where(elem => elem.IsFinite && factorizer(equation.Substitute(x, elem)).IsFinite).ToSet();
+ // An equation written out as expr = 0 has had its answers checked against it
+ // since the extraneous root of ln(x) + ln(x+1) = 0 was fixed, but the same
+ // equation written as the bare expr came through here instead and was only
+ // checked for being finite. The two disagreed:
+ // (2^x + 2^(2x) - 6 = 0).Solve(x) gave { 1 } while
+ // "2^x + 2^(2x) - 6".SolveEquation(x) gave { 1, ln((-3) ^ (1 / ln(2))) },
+ // and the second of those is not a root of anything.
+ return StatementSolver.WithoutSpuriousRoots(finite, equation, x);
}
else
return solutions;
diff --git a/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/SolveStatement.cs b/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/SolveStatement.cs
index f0155a6b2..1bd603227 100644
--- a/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/SolveStatement.cs
+++ b/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/SolveStatement.cs
@@ -6,6 +6,7 @@
//
using AngouriMath.Extensions;
+using PeterO.Numbers;
using AngouriMath.Functions.Continuous.Solvers.SetSolver;
using static AngouriMath.Entity;
using static AngouriMath.Entity.Set;
@@ -36,11 +37,35 @@ private static Entity Minus(Entity left, Entity right)
/// substitutions that follows, so the answers are checked once, here, against the
/// equation as the caller wrote it.
///
- private static Set WithoutSpuriousRoots(Set roots, Entity equation, Variable x)
+ internal static Set WithoutSpuriousRoots(Set roots, Entity equation, Variable x)
=> roots is FiniteSet finite && finite.Any(root => IsSpurious(equation, x, root))
? finite.Where(root => !IsSpurious(equation, x, root)).ToSet()
: roots;
+ ///
+ /// How large the residual has to be next to the terms it is the sum of before it
+ /// counts as evidence against a root.
+ ///
+ ///
+ /// A residual that is merely non-zero is not evidence. A root found numerically is
+ /// only as accurate as the search that found it, and the terms of the equation at
+ /// such a root cancel to within that accuracy rather than exactly. Judging the
+ /// residual on its own size threw away every root of
+ /// 1/210 - 17x/210 + 101x^2/210 - 247x^3/210 + x^4, whose four roots come
+ /// back as decimals and leave 5.2e-6 against terms of about 1.2e-2, and answered
+ /// that quartic with the empty set.
+ ///
+ /// The two cases are far apart, which is what makes a threshold possible at all.
+ /// A spurious root is not a near miss: it comes from a rewrite that widened the
+ /// domain, so the residual is a whole quantity in its own right -- 2*pi*i for the
+ /// logarithm sum, or the difference between -3 and 9 for an exponential. Over the
+ /// equations measured, every genuine root sits at 4.5e-4 or below and most at
+ /// 1e-15, while every spurious one sits between 1.5 and 2. This is a hundredth:
+ /// twenty times above the worst genuine root and a hundred and fifty times below
+ /// the mildest spurious one.
+ ///
+ [ConstantField] private static readonly EDecimal RelativeResidualTolerance = EDecimal.FromString("0.01");
+
///
/// Whether a root demonstrably fails the equation. Anything that cannot be
/// evaluated to a number -- a root still carrying a parameter, say
@@ -53,14 +78,33 @@ private static bool IsSpurious(Entity equation, Variable x, Entity root)
return false;
try
{
- return equation.Substitute(x, root).Evaled is Number.Complex residual
- && residual.IsFinite
- && !residual.Abs().EDecimal.LessThan(MathS.Settings.PrecisionErrorCommon);
+ var substituted = equation.Substitute(x, root);
+ if (substituted.Evaled is not Number.Complex residual || !residual.IsFinite)
+ return false;
+ var size = residual.Abs().EDecimal;
+ if (size.LessThan(MathS.Settings.PrecisionErrorCommon))
+ return false;
+ return size.GreaterThan(LargestTerm(substituted).Multiply(RelativeResidualTolerance));
}
catch (Core.Exceptions.AngouriBugException) { throw; }
catch (System.Exception) { return false; }
}
+ ///
+ /// The largest of the terms the equation adds up to at the root, which is the scale
+ /// the residual has to be read against. A single term is its own largest, so an
+ /// equation that is not a sum is judged on the size of its residual alone -- there
+ /// is nothing there for a root to have cancelled inexactly.
+ ///
+ private static EDecimal LargestTerm(Entity substituted)
+ {
+ var largest = EDecimal.Zero;
+ foreach (var term in Sumf.LinearChildren(substituted))
+ if (term.Evaled is Number.Complex value && value.IsFinite)
+ largest = EDecimal.Max(largest, value.Abs().EDecimal);
+ return largest;
+ }
+
internal static Set Solve(Entity expr, Variable x)
=> expr switch
{
diff --git a/Sources/Tests/UnitTests/Algebra/SolveTest/SpuriousRootTest.cs b/Sources/Tests/UnitTests/Algebra/SolveTest/SpuriousRootTest.cs
new file mode 100644
index 000000000..f020f1d99
--- /dev/null
+++ b/Sources/Tests/UnitTests/Algebra/SolveTest/SpuriousRootTest.cs
@@ -0,0 +1,126 @@
+//
+// Copyright (c) 2019-2022 Angouri.
+// AngouriMath is licensed under MIT.
+// Details: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.
+// Website: https://am.angouri.org.
+//
+
+using System.Collections.Generic;
+using System.Linq;
+using AngouriMath;
+using AngouriMath.Extensions;
+using Xunit;
+
+namespace AngouriMath.Tests.Algebra.SolveTest
+{
+ ///
+ /// An equation and the same equation written out as = 0 have to be answered
+ /// alike, and neither answer may contain something that is not a root. No issue is
+ /// filed for this; it was found while checking whether
+ /// https://github.com/asc-community/AngouriMath/issues/214 could be closed.
+ ///
+ public sealed class SpuriousRootTest
+ {
+ ///
+ /// Whether every root returned satisfies the equation it came from. Checked
+ /// numerically, because the point is what the answers are worth as numbers and not
+ /// what form they are written in.
+ ///
+ private static void AssertEveryRootSatisfies(string expression, params string[] variables)
+ {
+ var expr = expression.ToEntity();
+ foreach (var variable in variables)
+ {
+ var roots = expr.SolveEquation(variable);
+ var finite = Assert.IsType(roots);
+ foreach (var root in finite)
+ {
+ if (root.Vars.Any())
+ continue;
+ var residual = expr.Substitute(variable, root).EvalNumerical();
+ Assert.True(
+ residual.Abs().EDecimal.ToDouble() < 1e-6,
+ $"{root.Stringize()} is not a root of {expression}; it leaves {residual.Stringize()}");
+ }
+ }
+ }
+
+ ///
+ /// The exponential solver substitutes t = e^x and inverts what it finds with
+ /// t = v^(1/p). Composing that with x = ln(t) is only the same number while the
+ /// imaginary part stays inside the principal branch, and for a negative v it does
+ /// not: ln((-3) ^ (1 / ln(2))) is 1.58 - 1.75i, a full turn away from the root
+ /// ln(-3) / ln(2) = 1.58 + 4.53i, and satisfies nothing.
+ ///
+ [Theory]
+ [InlineData("2 ^ x + 2 ^ (2 * x) - 6")]
+ [InlineData("2 ^ (x ^ 2 + 1) + 2 ^ (2 * x ^ 2 - 1) - 3")]
+ [InlineData("3 ^ x + 3 ^ (2 * x) - 12")]
+ [InlineData("e ^ x + e ^ (2 * x) - 2")]
+ [InlineData("2 ^ (x + 1) + 2 ^ (2 * x ^ 2 - 1) + 1")]
+ public void AnExponentialEquationReturnsOnlyRoots(string expression) =>
+ AssertEveryRootSatisfies(expression, "x");
+
+ // The answer must not depend on whether the caller wrote the equation out.
+ [Theory]
+ [InlineData("2 ^ x + 2 ^ (2 * x) - 6")]
+ [InlineData("ln(x) + ln(x + 1)")]
+ [InlineData("x ^ 2 - 4")]
+ [InlineData("sin(x)")]
+ [InlineData("x ^ 3 - 1")]
+ public void WritingOutTheEqualsSignChangesNothing(string expression)
+ {
+ var bare = Assert.IsType(expression.ToEntity().SolveEquation("x"));
+ var written = Assert.IsType($"{expression} = 0".ToEntity().Solve("x"));
+ // Compared as numbers rather than as trees: the two arrive at the same roots by
+ // different routes and write them differently, sqrt(3/4) * (-1) against
+ // (-1/2) * sqrt(3) and -1/2 against -0.5, while the claim here is only about
+ // which numbers are answered.
+ Assert.Equal(AsNumbers(written), AsNumbers(bare));
+ }
+
+ private static List AsNumbers(Entity.Set.FiniteSet roots)
+ {
+ var numbers = roots
+ .Where(root => !root.Vars.Any())
+ .Select(root => root.EvalNumerical())
+ .Select(value => $"{value.RealPart.EDecimal.ToDouble():G8} {value.ImaginaryPart.EDecimal.ToDouble():G8}")
+ .ToList();
+ numbers.Sort(System.StringComparer.Ordinal);
+ return numbers;
+ }
+
+ ///
+ /// A root found numerically is only as accurate as the search that found it, so the
+ /// terms of the equation cancel at it to within that accuracy rather than exactly.
+ /// This quartic's four roots come back as decimals and leave 5.2e-6 against terms of
+ /// about 1.2e-2, and judging that residual on its own size answered the whole
+ /// equation with the empty set.
+ ///
+ [Theory]
+ [InlineData("1/210 - (17*x)/210 + (101*x^2)/210 - (247*x^3)/210 + x^4", 4)]
+ [InlineData("x^4 - 10*x^3 + 35*x^2 - 50*x + 24", 4)]
+ public void RootsThatOnlyCancelToTheAccuracyTheyWereFoundWithAreKept(string expression, int count)
+ {
+ var bare = Assert.IsType(expression.ToEntity().SolveEquation("x"));
+ var written = Assert.IsType($"{expression} = 0".ToEntity().Solve("x"));
+ Assert.Equal(count, bare.Count);
+ Assert.Equal(count, written.Count);
+ }
+
+ // Nothing is dropped from the equations that were already answered correctly, and
+ // a root carrying a parameter is kept rather than judged.
+ [Theory]
+ [InlineData("x ^ 2 - 4", 2)]
+ [InlineData("x ^ 2 + 1", 2)]
+ [InlineData("x ^ 3 - 1", 3)]
+ [InlineData("2 ^ x - 8", 1)]
+ [InlineData("ln(x) ^ 2 - 3 * ln(x) + 2", 2)]
+ [InlineData("sin(x)", 2)]
+ public void TheRootsThatWereRightAreStillThere(string expression, int count)
+ {
+ var roots = Assert.IsType(expression.ToEntity().SolveEquation("x"));
+ Assert.Equal(count, roots.Count);
+ }
+ }
+}