From d8b2a3ce61c623f843c5d5ca17955bbe845a6eb5 Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Wed, 5 Aug 2026 18:41:38 +0000 Subject: [PATCH 1/2] Solve a polynomial through its factors rather than by its degree (#272) The solver read a polynomial equation only as a polynomial: it gathered the monomials, took the degree, and applied the formula for that degree -- Cardano at three, Ferrari at four, nothing at all above four. A factorization was never looked for, and one handed over outright was expanded away before the degree was taken. So (x - 1)(x^2 - 3) = 0, whose roots are written on its face, came back as two nested cube roots of 26 + 18i, and a quintic that splits into a cubic's worth of whole roots and a quadratic fell through to the numeric solver. Two things are added. A product is now solved factor by factor, since a product is zero exactly where one of its factors is; and where the equation is written out as a sum, its rational roots are divided out first, which leaves a factor of lower degree for the formulas to answer exactly. Zeroing one factor only makes the product zero where the others are defined, so each root is checked against the whole product before it is kept -- x * (1/x) is undefined at 0, not zero there, and the verification the solvers share cannot catch that, since it drops a root only on positive evidence and a NaN residual is not evidence. Two-term polynomials are left alone. a*x^n + b is answered whole by inverting x^n = -b/a, which gives the roots as a + bi; splitting it off at its rational roots leaves the rest to be dug out of a quotient, and x^3 - 8 would read (-2 -+ sqrt(-12))/2 rather than (-1/2 +- i*sqrt(3)/2)*2. A test pins that. lim of the measured effect, against 40ce1093: x^3 - x^2 - 3x + 3 = 0 -(-1 + (26 + 18i)^(1/3) + 10/(26 + 18i)^(1/3))/3, ... -> { 1, sqrt(3), -sqrt(3) } 160 ms -> 9 ms x^5 - x^4 - 7x^3 + 7x^2 + 12x - 12 = 0 { -2, -1.73205080756887741... - 2.3e-15i, 1, ... } -> { 1, 2, -2, sqrt(3), -sqrt(3) } 153 ms -> 12 ms (x - 1)(x^2 - 3) = 0 -> { 1, sqrt(3), -sqrt(3) } 138 ms -> 8 ms x^3 - 6x^2 + 11x - 6 = 0-> { 1, 2, 3 } 160 ms -> 1 ms 4929 unit tests and 130 F# tests pass; corpus 112/117 with 0 wrong, 0 error and 0 timeout, unchanged. Co-Authored-By: Claude Opus 5 --- .../AnalyticalEquationSolver.cs | 43 ++++++ .../Simplification/PolynomialFactoring.cs | 56 +++++++ .../SolveTest/FactoredPolynomialSolveTest.cs | 144 ++++++++++++++++++ 3 files changed, 243 insertions(+) create mode 100644 Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs diff --git a/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/AnalyticalEquationSolver.cs b/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/AnalyticalEquationSolver.cs index dbe75bbfd..37a9c24ed 100644 --- a/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/AnalyticalEquationSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Solvers/EquationSolver/AnalyticalEquationSolver.cs @@ -140,10 +140,53 @@ static Entity TryDowncast(Entity equation, Variable x, Entity root) return Solve(expression, x, compensateSolving).Filter(predicate, x); case Piecewise p: return EquationSolver.SolvePiecewise(p, x, (e, x) => Solve(e, x, compensateSolving)); + // A product is zero exactly where one of its factors is. Without this the + // factors are expanded back into a polynomial and the equation is answered + // by the general formula for whatever degree that turns out to be, so + // (x - 1)(x^2 - 3) = 0 -- handed over already factored -- came back as two + // nested cube roots of 26 + 18i rather than as 1 and +-sqrt(3). + // https://github.com/asc-community/AngouriMath/issues/272 + // + // A factor free of x contributes no roots. It could still be zero itself, + // which would make every x a solution, but that is a reading the solver + // does not take anywhere else either: a * (x^2 - 3) = 0 already answered + // +-sqrt(3) by dividing a out, so nothing changes here. + case Mulf(var left, var right) when left.ContainsNode(x) || right.ContainsNode(x): + { + var factors = new List(2); + if (left.ContainsNode(x)) factors.Add(left); + if (right.ContainsNode(x)) factors.Add(right); + var product = expr; + // Zeroing one factor makes the product zero only where the others are + // defined: x * (1/x) is undefined at 0, not zero there. This has to be + // checked here rather than left to the verification every solver shares, + // because that deliberately keeps a root whose residual does not + // evaluate to a number at all -- a root is only dropped on positive + // evidence, and NaN is not evidence (see SolveStatement.IsSpurious). + bool ProductIsDefinedAt(Entity root) + => root.Vars.Any() + || product.Substitute(x, root).Evaled is not Complex value + || value.IsFinite; + var fromFactors = factors.Select(factor => Solve(factor, x, compensateSolving)).Unite(); + // Collapsed first: uniting sets builds a Unionf, and the roots cannot be + // looked at one at a time until it is a finite set again. + return fromFactors.InnerSimplified is FiniteSet found + ? found.Where(ProductIsDefinedAt).ToSet() + : fromFactors; + } default: break; } + // The same equation written out as a sum. Where a rational root can be divided + // out, what is left is a polynomial of lower degree, and the factors are then + // solved one at a time by the case above -- which is how x^3 - x^2 - 3x + 3 is + // answered exactly, and how a quintic is answered at all. + // https://github.com/asc-community/AngouriMath/issues/272 + if (Functions.PolynomialFactoring.TrySplitOffRationalRoots(expr, x, out var split) + && split is Mulf) + return Solve(split, x, compensateSolving); + if (PolynomialSolver.SolveAsPolynomial(expr, x, out var isIdentity) is { } poly) return poly.Select(e => TryDowncast(expr, x, e.InnerSimplified)).ToSet(); if (isIdentity) diff --git a/Sources/AngouriMath/Functions/Simplification/PolynomialFactoring.cs b/Sources/AngouriMath/Functions/Simplification/PolynomialFactoring.cs index 159551cd7..08228cd40 100644 --- a/Sources/AngouriMath/Functions/Simplification/PolynomialFactoring.cs +++ b/Sources/AngouriMath/Functions/Simplification/PolynomialFactoring.cs @@ -105,6 +105,62 @@ internal static bool TryFactor(Entity expr, Variable x, [NotNullWhen(true)] out return degree > 1; } + /// + /// written as the linear factors of its rational roots + /// times whatever will not divide, or where it has no + /// rational root at all. + /// + /// + /// + /// answers a different question and answers it more + /// strictly: it offers a factorization only when the polynomial splits completely + /// into whole linear factors, because a partly factored sum is not obviously a + /// better way of *writing* the same thing. Solving is not a question of how to + /// write it. (x - 1)(x^2 - 3) is no nicer to read than + /// x^3 - x^2 - 3x + 3, but it is the difference between answering with + /// 1 and +-sqrt(3) and answering with two nested cube roots of + /// 26 + 18i -- and above degree four, between answering at all and handing + /// the equation to the numeric solver. + /// https://github.com/asc-community/AngouriMath/issues/272 + /// + /// + /// Fractional roots are kept here where declines them, for + /// the same reason: a root of 1/2 is as exact an answer as any, whatever + /// the factored form looks like on the page. + /// + /// + internal static bool TrySplitOffRationalRoots( + Entity expr, Variable x, [NotNullWhen(true)] out Entity? factored) + { + factored = null; + if (!TryGetRationalCoefficients(expr, x, out var coefficients)) + return false; + + // Two terms is a x^n + b, which is answered whole by inverting x^n = -b/a, and + // that gives the n roots in the form a + bi. Splitting it would divide out the + // rational ones and leave the rest to be dug out of a quotient: x^3 - 8 reads + // as 2 and (-1/2 +- i*sqrt(3)/2)*2 that way, and as 2 and (-2 -+ sqrt(-12))/2 + // this way. Nothing here improves on that, so it is left alone. + if (coefficients.Count(coefficient => !coefficient.IsZero) <= 2) + return false; + + var factors = new List(); + foreach (var candidate in RootCandidates(coefficients)) + while (coefficients.Length > 1 && Evaluate(coefficients, candidate).IsZero) + { + factors.Add(LinearFactor(candidate, x)); + coefficients = DivideByLinear(coefficients, candidate); + } + if (factors.Count == 0) + return false; + + var remainder = BuildPolynomial(coefficients, x); + if (remainder != Integer.One) + factors.Add(remainder); + factored = factors.Aggregate((left, right) => left * right); + return true; + } + /// /// One step of a partial fraction decomposition, at a rational root of the /// denominator: N/D becomes A/(x - r) + R/Q, where D = (x - r)Q. diff --git a/Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs b/Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs new file mode 100644 index 000000000..0e1fd9925 --- /dev/null +++ b/Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs @@ -0,0 +1,144 @@ +// +// 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.Linq; +using AngouriMath; +using AngouriMath.Extensions; +using Xunit; + +namespace AngouriMath.Tests.Algebra +{ + /// + /// A polynomial equation that has a factor the solver can answer exactly must be + /// answered through that factor, rather than by handing the whole polynomial to the + /// general formula for its degree -- or, above degree four, to the numeric solver. + /// https://github.com/asc-community/AngouriMath/issues/272 + /// + public sealed class FactoredPolynomialSolveTest + { + private static Entity.Set.FiniteSet Roots(string equation) => + (Entity.Set.FiniteSet)equation.ToEntity().Solve("x"); + + private static void AssertRootsAre(string equation, params string[] expected) + { + var roots = Roots(equation); + foreach (var root in expected) + Assert.Contains(root.ToEntity().InnerSimplified, roots); + Assert.Equal(expected.Length, roots.Count); + } + + /// + /// The clearest statement of the defect, and the one that does not depend on how a + /// radical is spelled: an equation whose roots are all real must not be answered + /// with expressions built out of complex numbers. Cardano's formula reaches the + /// real roots of x^3 - x^2 - 3x + 3 through (26 + 18i)^(1/3). + /// + private static void AssertNoImaginaryPartsInTheWorking(string equation) + { + foreach (var root in Roots(equation)) + foreach (var node in root.Nodes) + if (node is Entity.Number.Complex and not Entity.Number.Real) + Assert.True(false, + $"{equation} has only real roots, but {root.Stringize()} reaches one " + + $"through {node.Stringize()}"); + } + + /// A root that is a decimal is a root the solver gave up on. + private static void AssertEveryRootIsExact(string equation) + { + foreach (var root in Roots(equation)) + foreach (var node in root.Nodes) + if (node is Entity.Number.Real and not Entity.Number.Rational) + Assert.True(false, + $"{equation} was answered with the decimal {node.Stringize()} " + + $"inside {root.Stringize()}"); + } + + // Handed the factors outright, the solver still expanded them into a cubic and ran + // Cardano over it, so the answer to an equation whose roots are 1 and +-sqrt(3) came + // back as two nested cube roots of 26 + 18i. + [Fact] + public void AProductIsSolvedThroughItsFactors() => + AssertRootsAre("(x - 1) * (x ^ 2 - 3) = 0", "1", "sqrt(3)", "-sqrt(3)"); + + [Fact] + public void AProductOfSeveralFactorsIsSolvedThroughThemToo() => + AssertRootsAre("(x - 1) * (x - 2) * (x + 2) * (x ^ 2 - 3) = 0", + "1", "2", "-2", "sqrt(3)", "-sqrt(3)"); + + // The same equation written out. A rational root is there to be found, and dividing + // it out leaves a quadratic that is answered exactly. + [Fact] + public void AnExpandedCubicWithARationalRootIsSplitAtIt() => + AssertRootsAre("x ^ 3 - x ^ 2 - 3 * x + 3 = 0", "1", "sqrt(3)", "-sqrt(3)"); + + // Above degree four there is no general formula, so the whole equation used to fall + // through to the numeric solver -- which returned decimals, one of them carrying a + // spurious -2.3e-15i. + [Fact] + public void AQuinticThatSplitsIntoSolvableFactorsIsSolvedExactly() => + AssertRootsAre("x ^ 5 - x ^ 4 - 7 * x ^ 3 + 7 * x ^ 2 + 12 * x - 12 = 0", + "1", "2", "-2", "sqrt(3)", "-sqrt(3)"); + + [Theory] + [InlineData("(x - 1) * (x ^ 2 - 3) = 0")] + [InlineData("x ^ 3 - x ^ 2 - 3 * x + 3 = 0")] + [InlineData("x ^ 3 - 6 * x ^ 2 + 11 * x - 6 = 0")] + public void RealRootsAreNotReachedThroughComplexNumbers(string equation) => + AssertNoImaginaryPartsInTheWorking(equation); + + [Theory] + [InlineData("x ^ 3 - x ^ 2 - 3 * x + 3 = 0")] + [InlineData("x ^ 5 - x ^ 4 - 7 * x ^ 3 + 7 * x ^ 2 + 12 * x - 12 = 0")] + [InlineData("(x - 1) * (x - 2) * (x + 2) * (x ^ 2 - 3) = 0")] + public void RootsThatCanBeWrittenExactlyAre(string equation) => + AssertEveryRootIsExact(equation); + + // What must not change. A polynomial with no rational root has nothing to split off, + // and its roots are what they always were. + [Fact] + public void APolynomialWithoutARationalRootIsUnaffected() => + AssertRootsAre("x ^ 2 - 2 = 0", "sqrt(2)", "-sqrt(2)"); + + // A two-term polynomial is answered whole by inverting x^n = c, which gives its + // roots as a + bi. Splitting off the rational ones instead leaves the others to be + // dug out of the quotient, and x^3 - 8 comes back as (-2 -+ sqrt(-12))/2 rather + // than as (-1/2 +- i*sqrt(3)/2)*2. Nothing is gained, so those are left alone. + [Theory] + [InlineData("x ^ 3 - 8 = 0")] + [InlineData("x ^ 6 - 1 = 0")] + public void ABinomialKeepsTheAnswerInversionGivesIt(string equation) + { + foreach (var root in Roots(equation)) + foreach (var node in root.Nodes) + if (node is Entity.Powf(Entity.Number.Real { IsNegative: true } radicand, Entity.Number.Rational)) + Assert.True(false, + $"{equation} was answered through a root of the negative " + + $"{radicand.Stringize()}, in {root.Stringize()}"); + } + + // Repeated factors are one root, not several. + [Fact] + public void ARepeatedFactorContributesItsRootOnce() => + AssertRootsAre("(x - 1) ^ 2 * (x + 3) = 0", "1", "-3"); + + // A factor free of the variable has no roots to contribute, and must not make the + // equation look unsolvable either. + [Fact] + public void AConstantFactorContributesNothing() => + AssertRootsAre("5 * (x ^ 2 - 3) = 0", "sqrt(3)", "-sqrt(3)"); + + // Solving factor by factor must not smuggle in a value that is not a root of the + // whole equation. Here x = 0 zeroes the first factor but makes the second undefined. + [Fact] + public void AFactorsRootThatTheWholeEquationRejectsIsDropped() + { + var roots = Roots("x * (1 / x) = 0"); + Assert.DoesNotContain(Entity.Number.Integer.Create(0), roots); + } + } +} From 38fcae788cb68773c3e8ee462e34bca2aef50438 Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Wed, 5 Aug 2026 18:54:19 +0000 Subject: [PATCH 2/2] Pin the roots the split turns out to recover A sweep over 572 polynomials built from factors whose roots are known in advance -- every product and every expansion of two or three factors drawn from a pool of linear and irreducible quadratic ones -- says the fix does more than tidy the output. Before it, ten of those came back **short of roots**, and after it none do; no case gained a root that is not one, in either direction. Every one of the ten had the same shape, a rational root times a quadratic times x^2 + x + 1: x^5 - 2x^3 - x^2 + 2 = 0 { 1, sqrt(2), -sqrt(2) } -> all five, the pair (-1 +- i*sqrt(3))/2 included Above degree four the equation went to the numeric solver whole, and that searches the real axis, so the complex conjugate pair was simply absent -- and a set two roots short reads exactly like a complete one. Splitting off the rational root leaves a quartic, which is answered in one piece. This is the stronger half of #272 and it was not visible from the issue, which reports the answers as ugly rather than as incomplete. Co-Authored-By: Claude Opus 5 --- .../SolveTest/FactoredPolynomialSolveTest.cs | 31 +++++++++++++++++++ 1 file changed, 31 insertions(+) diff --git a/Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs b/Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs index 0e1fd9925..9e04cefda 100644 --- a/Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs +++ b/Sources/Tests/UnitTests/Algebra/SolveTest/FactoredPolynomialSolveTest.cs @@ -98,6 +98,37 @@ public void RealRootsAreNotReachedThroughComplexNumbers(string equation) => public void RootsThatCanBeWrittenExactlyAre(string equation) => AssertEveryRootIsExact(equation); + private static bool IsNear(Entity root, double re, double im) + { + if (root.Vars.Any()) + return false; + var value = root.EvalNumerical(); + return System.Math.Abs(value.RealPart.EDecimal.ToDouble() - re) < 1e-9 + && System.Math.Abs(value.ImaginaryPart.EDecimal.ToDouble() - im) < 1e-9; + } + + /// + /// Worse than an unreadable answer: an incomplete one. Above degree four the whole + /// equation went to the numeric solver, which searches for real roots, so + /// (x - 1)(x^2 - 2)(x^2 + x + 1) = 0 came back as { 1, sqrt(2), -sqrt(2) } -- the + /// two roots of its x^2 + x + 1 factor were missing, and a set that is short of two + /// roots does not look any different from one that is not. + /// Splitting off the rational root leaves a quartic, which is solved whole. + /// + [Theory] + [InlineData("x ^ 5 - 2 * x ^ 3 - x ^ 2 + 2 = 0")] + [InlineData("(x - 1) * (x ^ 2 - 2) * (x ^ 2 + x + 1) = 0")] + public void AQuinticKeepsTheRootsOfItsIrreducibleFactor(string equation) + { + var roots = Roots(equation); + var sqrt2 = System.Math.Sqrt(2); + var halfSqrt3 = System.Math.Sqrt(3) / 2; + foreach (var (re, im) in new[] { (1.0, 0.0), (sqrt2, 0.0), (-sqrt2, 0.0), + (-0.5, halfSqrt3), (-0.5, -halfSqrt3) }) + Assert.Contains(roots, root => IsNear(root, re, im)); + Assert.Equal(5, roots.Count); + } + // What must not change. A polynomial with no rational root has nothing to split off, // and its roots are what they always were. [Fact]