A power of x beside two blocks in one power of x is split over the blocks first

BREAKING-CHANGES.md

@@ -294,6 +294,23 @@ takes 6 s for 265,499 now ([#718](https://github.com/asc-community/AngouriMath/i
 | `"(d + f*x + g*x^2 + h*x^3 + k*x^4 + m*x^5 + n*x^8 + p*x^11)/(a + b*x^2 + c*x^4)^3".ToEntity().Integrate("x")` | `integral(...)` | 77,389 characters |
 | `"csc(c + d*x)^2/(a - b*sin(c + d*x)^4)^3".ToEntity().Integrate("x")` | `integral(...)` | 157,787 characters |
 
+### A power of x beside two blocks in one power of x is split over the blocks first
+
+**Shorter answers, and sooner.** `1/(x^7 (a + b x^3)(c + d x^3))` came to 265,499 characters on the unreleased
+master: the Hermite reduction read `x^7` as a repeated factor, and what it left went to the split at the power of x,
+whose rest over both blocks had coefficients pages long. The blocks are split in `u = x^3` first now, with the power
+of x kept outside, and each piece is a power of x over one block: 351 characters. The other `1/(x^m (a + b x^n)(c + d x^n))`
+in Rubi's 1.1.3.4 come out shorter as well, and those that took two seconds take milliseconds
+([#718](https://github.com/asc-community/AngouriMath/issues/718)).
+
+| Input | Was (2.5.0) | Now |
+|---|---|---|
+| `"1/(x*(a + b*x^3)*(c + d*x^3))".ToEntity().Integrate("x")` | `integral(...)` | 185 characters |
+| `"1/(x^4*(a + b*x^3)*(c + d*x^3))".ToEntity().Integrate("x")` | `integral(...)` | 259 characters |
+| `"1/(x^7*(a + b*x^3)*(c + d*x^3))".ToEntity().Integrate("x")` | `integral(...)` | 351 characters |
+| `"1/(x^8*(a + b*x^3)*(c + d*x^3))".ToEntity().Integrate("x")` | `integral(...)` | 1,217 characters |
+| `"1/(x^6*(a + b*x^4)*(c + d*x^4))".ToEntity().Integrate("x")` | `integral(...)` | 795 characters |
+
 ### One linear twice is not two linear powers
 
 **Answers that had no value.** `(a c + b c x)^(-3 - 2 p) (f + g x) (a^2 + 2 a b x + b^2 x^2)^p` was answered

Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs

@@ -1057,6 +1057,14 @@ is var (multiple, leftover)
                     ?? Integration.ComputeIndefiniteIntegral(numerator / overTheConjugates, x, integrateByParts)) is { } overConjugates)
                 return overConjugates;
 
+            // A power of x beside two blocks or more in one power of x, `1/(x^7 (a + b x^3)(c + d x^3))`:
+            // split in `u = x^3` over the blocks, each fraction then a power of x over one block.
+            // Before the Hermite reduction, which reads the power of x as a repeated factor: what it
+            // leaves went to the split at the power of x, whose rest over both blocks had coefficients
+            // pages long, and the answer was 265,499 characters.
+            if (IntegrateOverBlocksInAPowerOfX(numerator, denominator, x, integrateByParts, besideAPowerOfX: true) is { } besideBlocksInAPower)
+                return besideBlocksInAPower;
+
             // A denominator with a written repeated factor takes the Hermite reduction first:
             // the rational part of the answer in one linear solve, and what is left is a proper
             // fraction over a squarefree denominator for the splits below. `(1 + x^2)/(x (1 + x^3)^2)`
@@ -1622,19 +1630,27 @@ static Entity Binomial(int n, int k)
         /// <c>n >= 2</c>, linear or quadratic in it, with a symbol among them: split in
         /// <c>u = x^n</c> by <see cref="Functions.PartialFractions.TrySplitOverWrittenFactors"/>,
         /// one residue of the numerator's powers modulo <c>n</c> at a time, each fraction then
-        /// over one block in x.
+        /// over one block in x. With <paramref name="besideAPowerOfX"/>, a power of x beside two
+        /// blocks or more as well, and only then.
         /// </summary>
         /// <remarks>
         /// <c>1/((a + b x^3)(c + d x^3))</c> was declined: over the rationals it does not factor,
         /// and the split over written factors reads a block of the third degree as nothing. In
         /// <c>u = x^3</c> it is <c>1/((a + b u)(c + d u))</c>, two linear factors. A numerator
         /// term <c>x^r M(x^n)</c> keeps its <c>x^r</c> outside the split, which is exact since
         /// the blocks are polynomials in <c>x^n</c> alone.
+        /// <para>
+        /// A power of x beside the blocks stays outside the split: <c>1/(x^7 (a + b x^3)(c + d x^3))</c>
+        /// is <c>(b/(x^7 (a + b x^3)) - d/(x^7 (c + d x^3)))/(b c - a d)</c>, a power of x over one
+        /// block twice. Split at the power of x first, the rest over both blocks had coefficients
+        /// that were pages of symbols, and the answer came to 265,499 characters.
+        /// </para>
         /// https://github.com/asc-community/AngouriMath/issues/718
         /// </remarks>
-        private static Entity? IntegrateOverBlocksInAPowerOfX(Entity numerator, Entity denominator, Entity.Variable x, bool integrateByParts, bool oneBlockAtItsRoots = false)
+        private static Entity? IntegrateOverBlocksInAPowerOfX(Entity numerator, Entity denominator, Entity.Variable x, bool integrateByParts, bool oneBlockAtItsRoots = false, bool besideAPowerOfX = false)
         {
             var blocks = new List<(Dictionary<EInteger, Entity> Polynomial, int Power)>();
+            var powerOfX = 0;
             Entity constant = Number.Integer.One;
             foreach (var factor in Mulf.LinearChildren(denominator))
             {
@@ -1646,6 +1662,13 @@ static Entity Binomial(int n, int k)
                 var (@base, power) = factor is Powf(var b, Number.Integer e) && e.EInteger.Sign > 0 && e.EInteger.CanFitInInt32()
                     ? (b, e.EInteger.ToInt32Unchecked())
                     : (factor, 1);
+                if (@base == x)
+                {
+                    if (!besideAPowerOfX)
+                        return null;
+                    powerOfX += power;
+                    continue;
+                }
                 if (!TreeAnalyzer.TryGetPolynomial(@base, x, out var polynomial) || polynomial.Count < 2)
                     return null;
                 foreach (var term in polynomial)
@@ -1660,6 +1683,8 @@ static Entity Binomial(int n, int k)
             }
             if (blocks.Count == 0 || !(numerator + denominator).Vars.Any(symbol => symbol != x))
                 return null;
+            if (besideAPowerOfX && (powerOfX == 0 || blocks.Count < 2 || powerOfX > MaximumPowerOfXOverABlock))
+                return null;
             // n, the greatest common divisor of every power the blocks hold.
             var n = EInteger.Zero;
             foreach (var (polynomial, _) in blocks)
@@ -1678,6 +1703,8 @@ static Entity Binomial(int n, int k)
             foreach (var term in above)
                 if (term.Key.Sign < 0 || !term.Key.CanFitInInt32() || term.Value.ContainsNode(x))
                     return null;
+            if (besideAPowerOfX && above.Values.Any(coefficient => coefficient is not (Number or Variable)))
+                return null;
 
             var u = Entity.Variable.CreateUnique(numerator + denominator, "u");
             Entity InU(Dictionary<EInteger, Entity> polynomial)
@@ -1735,11 +1762,33 @@ Entity InU(Dictionary<EInteger, Entity> polynomial)
                 }
                 if (!Functions.PartialFractions.TrySplitOverWrittenFactors(inU, below, u, out var decomposition))
                     return null;
-                var outside = residue.Key == 0 ? (Entity)Number.Integer.One : residue.Key == 1 ? x : MathS.Pow(x, residue.Key);
-                foreach (var fraction in Sumf.LinearChildren(decomposition))
+                // The power of x below stays outside the split, so each fraction is a power of x
+                // over one block, which the split at a power of x beside a block reads.
+                var outside = (residue.Key - powerOfX) switch
                 {
+                    0 => (Entity)Number.Integer.One,
+                    1 => x,
+                    var power => MathS.Pow(x, power),
+                };
+                // The split comes as a sum over a constant. Beside a power of x it is taken term by
+                // term: whole, it is the power of x over both blocks again, for the split at that
+                // power, whose coefficients run to `a^54 b^6 c^46 d^14`.
+                var (terms, scale) = besideAPowerOfX && decomposition is Divf(var sum, var divisor) && !divisor.ContainsNode(u)
+                    ? (sum, divisor)
+                    : (decomposition, (Entity)Number.Integer.One);
+                foreach (var term in Sumf.LinearChildren(terms))
+                {
+                    if (term == Number.Integer.Zero)
+                        continue;
+                    var fraction = scale == Number.Integer.One ? term : term / scale;
                     var inX = outside * fraction.Substitute(u, MathS.Pow(x, step));
-                    if (Integration.ComputeIndefiniteIntegral(inX, x, integrateByParts) is not { } integrated)
+                    // Beside a power of x, each piece is a part of the question asked, over one block
+                    // fewer, and answered at its level: the split at the power of x that reads it
+                    // answers only the question asked or one below.
+                    var integrated = besideAPowerOfX
+                        ? Integration.ComputeAsTheSameQuestion(inX, x, integrateByParts)
+                        : Integration.ComputeIndefiniteIntegral(inX, x, integrateByParts);
+                    if (integrated is null)
                         return null;
                     total += integrated;
                 }

Sources/Tests/UnitTests/Calculus/APowerOfXBesideBlocksInAPowerOfXIntegralTest.cs

@@ -0,0 +1,53 @@
+//
+// 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
+{
+    /// <summary>
+    /// A power of x beside two blocks in one power of x is split over the blocks in
+    /// <c>u = x^n</c> with the power of x kept outside: <c>1/(x^7 (a + b x^3)(c + d x^3))</c> is a
+    /// power of x over one block twice. Rubi's 1.1.3.4.
+    /// <a href="https://github.com/asc-community/AngouriMath/issues/718">#718</a>
+    /// </summary>
+    [Trait("Area", "Calculus")]
+    public sealed class APowerOfXBesideBlocksInAPowerOfXIntegralTest
+    {
+        [Theory]
+        [InlineData("1/(x*(a + b*x^3)*(c + d*x^3))")]
+        [InlineData("1/(x^4*(a + b*x^3)*(c + d*x^3))")]
+        [InlineData("1/(x^7*(a + b*x^3)*(c + d*x^3))")]
+        [InlineData("1/(x^6*(a + b*x^4)*(c + d*x^4))")]
+        public void DifferentiatesBack(string integrand)
+        {
+            var integral = integrand.ToEntity().Integrate("x");
+            Assert.DoesNotContain("integral(", integral.Stringize());
+            Entity Pinned(Entity e) => e.Substitute("a", 2).Substitute("b", 3).Substitute("c", 5).Substitute("d", 7);
+            var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x");
+            var original = Pinned(integrand.ToEntity());
+            foreach (var at in new[] { 0.35, 0.6, 0.85, 1.4 })
+            {
+                var want = original.Substitute("x", at).EvalNumerical();
+                var got = derivative.Substitute("x", at).EvalNumerical();
+                Assert.True(Math.Abs((double)(got - want).RealPart) + Math.Abs((double)(got - want).ImaginaryPart)
+                            < 1e-9 * Math.Max(1, Math.Abs((double)want.RealPart)),
+                    $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}");
+            }
+        }
+
+        /// <summary>
+        /// Split at the power of x first, the rest over both blocks had coefficients pages long,
+        /// and the answer was 265,499 characters; over the blocks first it is 351.
+        /// </summary>
+        [Fact]
+        public void TheAnswerIsTheSizeOfItsPieces()
+            => Assert.True("1/(x^7*(a + b*x^3)*(c + d*x^3))".ToEntity().Integrate("x").Stringize().Length < 2000);
+    }
+}