From 95f536a25235326dc0f180dbefa5dc2b26acad37 Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Thu, 15 Oct 2026 19:21:30 +0000 Subject: [PATCH 1/2] Split a power of x over a sum of two cubes over the cube roots of unity Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- .../Integration/IndefiniteIntegralSolver.cs | 89 +++++++++++++++++++ .../Integration/Integration.Definition.cs | 3 + .../Calculus/ASumOfTwoCubesIntegralTest.cs | 59 ++++++++++++ 3 files changed, 151 insertions(+) create mode 100644 Sources/Tests/UnitTests/Calculus/ASumOfTwoCubesIntegralTest.cs diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index c005eee53..f78fa02c3 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -19588,6 +19588,95 @@ void Add(int j, int k, Entity[] poly) return Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x) ? answer : null; } + /// + /// A power of x over a sextic that is a sum of two cubes, + /// d6 ((x^2 + s)^3 + q x^3), split over the cube roots of unity: + /// 1/(U^3 + V^3) = sum_j w_j/(3 V^2 (U + w_j V)) for U = x^2 + s, + /// V = m x, m^3 = q, so that x^k/(U^3 + V^3) is + /// sum_j w_j x^(k - 2)/(3 m^2 (x^2 + w_j m x + s)), three quadratics. + /// + /// + /// The sextic is read as d6 x^6 + d4 x^4 + d3 x^3 + d2 x^2 + d0 with + /// d4^2 = 3 d6 d2 and d2^2 = 3 d4 d0, which is what makes + /// d6 x^6 + d4 x^4 + d2 x^2 + d0 the cube d6 (x^2 + s)^3 with + /// s = d4/(3 d6). The identity is the partial fractions of 1/(U^3 + V^3) in + /// U, whose poles -w_j V have residues 1/(3 w_j^2 V^2) = w_j/(3 V^2); any + /// cube root m of q gives it. The three quadratics are one integral + /// x^n/(x^2 + B x + S) asked with B and S named. Rubi's 1.3.1 + /// x^k/(27 a^3 + 27 a^2 b x^2 + 27 a^2 c x^3 + 9 a b^2 x^4 + b^3 x^6) was declined, + /// the integrator's resultant needing numbers for coefficients. + /// https://github.com/asc-community/AngouriMath/issues/718 + /// + internal static Entity? SolveAPowerOfXOverASumOfTwoCubes(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (!Integration.AnsweringTheQuestionAsked || !expr.Vars.Any(symbol => symbol != x)) + return null; + var (numerator, denominator) = Functions.SingleQuotient.Of(expr); + Entity? sextic = null; + var power = 0; + Entity constant = Number.Integer.One; + foreach (var (side, below) in new[] { (numerator, false), (denominator, true) }) + foreach (var factor in Mulf.LinearChildren(side)) + { + if (!factor.ContainsNode(x)) + { + constant = below ? constant / factor : constant * factor; + continue; + } + if (factor == x || factor is Powf(var bx, Number.Integer { EInteger: var wx }) && bx == x && wx.CanFitInInt32()) + { + var whole = factor == x ? 1 : ((Number.Integer)((Powf)factor).Exponent).EInteger.ToInt32Unchecked(); + power += below ? -whole : whole; + continue; + } + if (below && sextic is null) + { + sextic = factor; + continue; + } + return null; + } + if (sextic is null || Functions.PartialFractions.HoldsTheImaginaryUnit(sextic) + || !TreeAnalyzer.TryGetPolynomial(sextic, x, out var read) || read.Keys.Any(key => key.Sign < 0 || key.CompareTo(EInteger.FromInt32(6)) > 0)) + return null; + Entity At(int degree) => read.TryGetValue(EInteger.FromInt32(degree), out var coefficient) ? coefficient.InnerSimplified : Number.Integer.Zero; + var (d6, d5, d4, d3, d2, d1, d0) = (At(6), At(5), At(4), At(3), At(2), At(1), At(0)); + if (!VanishesIdentically(d5) || !VanishesIdentically(d1) + || VanishesIdentically(d6) || VanishesIdentically(d4) || VanishesIdentically(d3) || VanishesIdentically(d0) + || !VanishesIdentically((d4 * d4 - 3 * d6 * d2).Expand()) || !VanishesIdentically((d2 * d2 - 3 * d4 * d0).Expand())) + return null; + var shift = Functions.PartialFractions.InLowestTermsOverTheSymbols(d4 / (3 * d6)); + var cube = Functions.PartialFractions.InLowestTermsOverTheSymbols(d3 / d6); + var root = MathS.Pow(cube, Number.Rational.Create(1, 3)); + var namedB = Variable.CreateUnique(expr, "B_c"); + var namedS = Variable.CreateUnique(expr + namedB, "S_c"); + var term = MathS.Pow(x, Number.Integer.Create(power - 2)) / (MathS.Pow(x, Number.Integer.Create(2)) + namedB * x + namedS); + if (Integration.ComputeAsAQuestionOfItsOwn(term, x, integrateByParts) is not { } each) + return null; + // The integral over a real quadratic splits on the sign of 4 S - B^2, which means nothing + // for the complex B of two of the three: its first arm, the arctangent over + // sqrt(4 S - B^2), differentiates back to the integrand for any B with 4 S != B^2, and + // the answer is checked below. + each = each.Replace(node => node switch + { + Piecewise piecewise when piecewise.Cases.Any() => piecewise.Cases.First().Expression, + Providedf provided => provided.Expression, + _ => node, + }); + var half = Number.Rational.Create(1, 2); + var unities = new Entity[] + { + Number.Integer.One, + -half + half * MathS.Sqrt(3) * MathS.i, + -half - half * MathS.Sqrt(3) * MathS.i, + }; + Entity answer = Number.Integer.Zero; + foreach (var unity in unities) + answer = answer + unity * each.Substitute(namedB, unity * root).Substitute(namedS, shift); + answer = (constant / (3 * d6 * MathS.Pow(root, Number.Integer.Create(2))) * answer).InnerSimplified; + return Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x) ? answer : null; + } + /// /// A power of a linear that is not whole, times a linear, over a whole power of a quadratic, /// with the power of the quadratic lowered one at a time: diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index ba13c868d..27552674c 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -678,6 +678,9 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // A rational function over a whole power of a symbolic quadratic, the power lowered one // at a time, before the partial fractions take the quadratic apart over its roots. if ((answer = IndefiniteIntegralSolver.SolveARationalFunctionOverAPowerOfASymbolicQuadratic(expr, x, integrateByParts)) is { }) return answer; + // A power of x over a sextic that is a sum of two cubes, split over the cube roots of + // unity into three quadratics before the partial fractions try the symbolic sextic. + if ((answer = IndefiniteIntegralSolver.SolveAPowerOfXOverASumOfTwoCubes(expr, x, integrateByParts)) is { }) return answer; // A root or a power of a product beside a power of one of its factors, the product // split: `(a + b x^2)^(5/2)/sqrt(a^2 - b^2 x^4)` is `K (a + b x^2)^2/sqrt(a - b x^2)`. // Early, since what reads the product whole is slower and longer: Rubi's 1.2.1.2:2218, diff --git a/Sources/Tests/UnitTests/Calculus/ASumOfTwoCubesIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ASumOfTwoCubesIntegralTest.cs new file mode 100644 index 000000000..4ef667494 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ASumOfTwoCubesIntegralTest.cs @@ -0,0 +1,59 @@ +// +// 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; +using AngouriMath.Extensions; +using Xunit; + +namespace AngouriMath.Tests.Calculus +{ + /// + /// A power of x over a sextic that is a sum of two cubes, + /// (3 a + b x^2)^3 + 27 a^2 c x^3, split over the cube roots of unity into three + /// quadratics. Rubi's 1.3.1 rows; the cube root of 27 a^2 c/b^3 is the principal one, + /// complex for a negative c, and the three quadratics are still the factors. + /// #718 + /// + [Trait("Area", "Calculus")] + public sealed class ASumOfTwoCubesIntegralTest + { + private static void DifferentiatesBack(string integrand, params (string, double)[] pins) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + + var pinned = integral.Substitute("C", 0); + Entity original = integrand.ToEntity(); + foreach (var (name, value) in pins) + { + pinned = pinned.Substitute(name, value); + original = original.Substitute(name, value); + } + var derivative = pinned.Differentiate("x"); + foreach (var at in new[] { -2.3, -0.7, 0.4, 1.9 }) + { + var got = derivative.Substitute("x", at).EvalNumerical(); + var want = original.Substitute("x", at).EvalNumerical(); + var difference = Math.Abs((double)(got - want).RealPart) + Math.Abs((double)(got - want).ImaginaryPart); + var scale = Math.Max(1.0, Math.Abs((double)want.RealPart) + Math.Abs((double)want.ImaginaryPart)); + Assert.True(difference / scale < 1e-8, + $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}"); + } + } + + [Theory] + [InlineData("x^4/(27*a^3+27*a^2*b*x^2+27*a^2*c*x^3+9*a*b^2*x^4+b^3*x^6)")] + [InlineData("x^2/(27*a^3+27*a^2*b*x^2+27*a^2*c*x^3+9*a*b^2*x^4+b^3*x^6)")] + [InlineData("1/(27*a^3+27*a^2*b*x^2+27*a^2*c*x^3+9*a*b^2*x^4+b^3*x^6)")] + [InlineData("1/(x^2*(27*a^3+27*a^2*b*x^2+27*a^2*c*x^3+9*a*b^2*x^4+b^3*x^6))")] + public void IsAnsweredAndDifferentiatesBack(string integrand) + { + DifferentiatesBack(integrand, ("a", 2), ("b", 3), ("c", 5)); + DifferentiatesBack(integrand, ("a", -1), ("b", 2), ("c", -3)); + } + } +} From 8f35b09956b81b00228a71b5d214f7bddee2570a Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Thu, 15 Oct 2026 20:29:52 +0000 Subject: [PATCH 2/2] Record the answers the split over the cube roots of unity gives Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- BREAKING-CHANGES.md | 13 +++++++++++++ 1 file changed, 13 insertions(+) diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 8e9ed4435..a2e79d456 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -447,6 +447,19 @@ case is most of its length ([#718](https://github.com/asc-community/AngouriMath/ | `"x^4/((a + b*x^2 + c*x^4)*sqrt(d + k*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 28,551 characters | | `"1/(x^2*(a + b*x^2 + c*x^4)*sqrt(d + k*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 10,194 characters | +### A power of x over a sum of two cubes is split over the cube roots of unity + +**Answers where there were none.** `x^k/((3 a + b x^2)^3 + 27 a^2 c x^3)` was declined: the integrator's resultant needs +numbers for coefficients. A sextic `d6 (x^2 + s)^3 + d3 x^3` is now split over the cube roots of unity into three +quadratics, `x^2 + w m x + s` with `m^3 = d3/d6`. For two of them `w` is complex, so the answer carries `i` the way +Rubi's does ([#718](https://github.com/asc-community/AngouriMath/issues/718)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"x^2/(27*a^3+27*a^2*b*x^2+27*a^2*c*x^3+9*a*b^2*x^4+b^3*x^6)".ToEntity().Integrate("x")` | `integral(...)` | 849 characters | +| `"1/(27*a^3+27*a^2*b*x^2+27*a^2*c*x^3+9*a*b^2*x^4+b^3*x^6)".ToEntity().Integrate("x")` | `integral(...)` | 2,091 characters | +| `"1/(x^2*(27*a^3+27*a^2*b*x^2+27*a^2*c*x^3+9*a*b^2*x^4+b^3*x^6))".ToEntity().Integrate("x")` | `integral(...)` | 3,820 characters | + ### A polynomial of degree six with a large leading coefficient is factored **Answers where there were none.** The factoring clears a leading coefficient by rescaling the variable, and declined