A polynomial times a fractional power of a base holding a function of x is answered as a polynomial times the next power of the base, by a linear system with the functions of x as indeterminates

Sources/.editorconfig

@@ -435,3 +435,6 @@ file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed
 
 [Tests/UnitTests/Calculus/RationalFunctionBesideTheRootOfAQuadraticTest.cs]
 file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed under MIT.\nDetails: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.\nWebsite: https://am.angouri.org.\n
+
+[Tests/UnitTests/Calculus/PolynomialTimesAPowerOfTheBaseTest.cs]
+file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed under MIT.\nDetails: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.\nWebsite: https://am.angouri.org.\n

Sources/AngouriMath/Functions/Algebra/Polynomials/MultivariatePolynomial.cs

@@ -251,6 +251,9 @@ internal MultivariatePolynomial ScaleBy(ERational factor)
             return new(VariableCount, result);
         }
 
+        /// <summary>The terms, each as its packed monomial and its coefficient.</summary>
+        internal IEnumerable<KeyValuePair<ulong, ERational>> Terms => terms;
+
         /// <summary>
         /// Read as a polynomial in <paramref name="variable"/> alone: the power of that
         /// variable mapped to the coefficient, itself a polynomial in the others.

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

@@ -3965,6 +3965,153 @@ private static (Entity Sine, Entity Cosine) SineAndCosineOfAMultiple(int n, Enti
             return ((sine * uPrevious).InnerSimplified, tCurrent);
         }
 
+        /// <summary>
+        /// <c>N B^r</c>, for a fractional <c>r</c> and an <c>N</c> and a <c>B</c> that are
+        /// polynomials in <paramref name="x"/> and the functions of it in them, answered as
+        /// <c>P(x) B^(r + 1)</c> for a polynomial <c>P</c> where there is one: the derivative of
+        /// that is <c>B^r (P' B + (r + 1) P B')</c>, and <c>P' B + (r + 1) P B' = N</c> is a linear
+        /// system in the coefficients of <c>P</c> once the functions of <c>x</c> are taken for
+        /// indeterminates.
+        /// </summary>
+        /// <remarks>
+        /// <para>
+        /// Bronstein's <c>(5x^2 + 3 (e^x + x)^(1/3) + e^x (3x + 2x^2))/(x (e^x + x)^(1/3))</c> is
+        /// <c>3/x + (5x + e^x (2x + 3))/(e^x + x)^(1/3)</c>, and the second is
+        /// <c>3x (e^x + x)^(2/3)</c>: with <c>P = c_0 + c_1 x</c> and <c>r = -1/3</c>,
+        /// <c>c_1 (e^x + x) + (2/3)(c_0 + c_1 x)(e^x + 1) = 5x + e^x (2x + 3)</c> is
+        /// <c>c_1 = 3</c>, <c>c_0 = 0</c>. Nothing else read it: the substitution
+        /// <c>u = e^x + x</c> wants <c>1 + e^x</c> beside the root and finds <c>5x + e^x (2x + 3)</c>.
+        /// </para>
+        /// <para>
+        /// The degree of <c>P</c> is the degree of <c>N</c> in <c>x</c>, the functions of <c>x</c>
+        /// held constant, and one more; the system is exact, over the rationals, and its
+        /// solution is checked by differentiating back at sampled points, since an identity
+        /// between two functions taken for independent indeterminates -- <c>e^x</c> and
+        /// <c>e^(2x)</c>, say -- can fail to be found and cannot be found wrongly, but a
+        /// derivative the reader did not expect can. Closed, and volunteered at any depth.
+        /// </para>
+        /// https://github.com/asc-community/AngouriMath/issues/718
+        /// </remarks>
+        internal static Entity? SolveByAPolynomialTimesAPowerOfTheBase(Entity expr, Entity.Variable x)
+        {
+            var (numerator, denominator) = Functions.SingleQuotient.Of(expr);
+            Entity? radicalBase = null;
+            Number.Rational? exponent = null;
+            Entity rest = Number.Integer.One;
+            Entity constantBelow = Number.Integer.One;
+            foreach (var (side, isBelow) in new[] { (numerator, false), (denominator, true) })
+                foreach (var factor in Mulf.LinearChildren(side))
+                {
+                    if (factor is Powf(var @base, Number.Rational power) && power is not Number.Integer && @base.ContainsNode(x))
+                    {
+                        if (radicalBase is not null)
+                            return null;
+                        radicalBase = @base;
+                        exponent = isBelow ? Number.Rational.Create(power.ERational.Negate()) : power;
+                        continue;
+                    }
+                    if (!factor.ContainsNode(x))
+                    {
+                        if (isBelow) constantBelow = constantBelow * factor;
+                        else rest = rest * factor;
+                        continue;
+                    }
+                    if (isBelow || HasARadicalOf(factor, x))
+                        return null;
+                    rest = rest * factor;
+                }
+            if (radicalBase is null || exponent is null || HasARadicalOf(radicalBase, x)
+                || !radicalBase.Nodes.Any(node => node.ContainsNode(x) && node is not (Variable or Sumf or Minusf or Mulf or Divf or Powf(_, Number.Integer))))
+                return null;
+
+            // The degree of N in x with the functions of x held constant.
+            var atoms = new Dictionary<Entity, Entity.Variable>();
+            var named = expr;
+            Entity WithAtoms(Entity e) => e.Replace(node =>
+            {
+                if (!node.ContainsNode(x) || node is Variable or Sumf or Minusf or Mulf or Divf || node is Powf(_, Number.Integer { IsNegative: false }))
+                    return node;
+                if (!atoms.TryGetValue(node, out var atom))
+                {
+                    atom = Variable.CreateUnique(named, "atom");
+                    named = named + atom;
+                    atoms[node] = atom;
+                }
+                return atom;
+            });
+            var n = WithAtoms(rest);
+            if (!TreeAnalyzer.TryGetPolynomial(n, x, out var nRead) || nRead.Keys.Any(k => !k.CanFitInInt32()))
+                return null;
+            var degree = (nRead.Count == 0 ? 0 : nRead.Keys.Max()!.ToInt32Unchecked()) + 1;
+            if (degree > 4)
+                return null;
+
+            // P' B + (r + 1) P B' - N, linear in the unknowns; written with the functions of x
+            // as indeterminates, as one polynomial, whose every coefficient must vanish.
+            var unknowns = new List<Entity.Variable>();
+            for (var k = 0; k <= degree; k++)
+            {
+                var unknown = Variable.CreateUnique(named, "c" + k);
+                named = named + unknown;
+                unknowns.Add(unknown);
+            }
+            Entity polynomial = Number.Integer.Zero;
+            for (var k = 0; k <= degree; k++)
+                polynomial = polynomial + unknowns[k] * (k == 0 ? Number.Integer.One : k == 1 ? x : MathS.Pow(x, k));
+            var rPlusOne = Number.Rational.Create(exponent.ERational.Add(ERational.One));
+            var identity = polynomial.Differentiate(x) * radicalBase + rPlusOne * polynomial * radicalBase.Differentiate(x) - rest / constantBelow;
+            var (top, bottom) = Functions.SingleQuotient.Of(Functions.PartialFractions.Bare(WithAtoms(Functions.PartialFractions.Bare(identity))));
+            if (bottom.ContainsNode(x) && !TreeAnalyzer.TryGetPolynomial(bottom, x, out _))
+                return null;
+            var variables = top.Vars.Where(v => !unknowns.Contains(v)).OrderBy(v => v.Name, System.StringComparer.Ordinal).ToList();
+            if (variables.Count == 0 || variables.Count > MultivariatePolynomial.MaxVariables)
+                return null;
+            var indices = new Dictionary<Variable, int>();
+            for (var i = 0; i < variables.Count; i++)
+                indices[variables[i]] = i;
+            // top = E_0 + sum_k c_k E_k.
+            Entity AtUnit(int which)
+            {
+                var e = top;
+                for (var k = 0; k <= degree; k++)
+                    e = e.Substitute(unknowns[k], k == which ? Number.Integer.One : Number.Integer.Zero);
+                return Functions.PartialFractions.Bare(e);
+            }
+            if (MultivariatePolynomial.TryParse(AtUnit(-1), indices) is not { } constantPart)
+                return null;
+            var columns = new List<MultivariatePolynomial>();
+            for (var k = 0; k <= degree; k++)
+            {
+                if (MultivariatePolynomial.TryParse(AtUnit(k), indices) is not { } column)
+                    return null;
+                columns.Add(column.Subtract(constantPart));
+            }
+            var monomials = columns.SelectMany(column => column.Terms.Select(term => term.Key)).Concat(constantPart.Terms.Select(term => term.Key)).Distinct().ToList();
+            var matrix = new Entity[monomials.Count][];
+            var rhs = new Entity[monomials.Count];
+            for (var row = 0; row < monomials.Count; row++)
+            {
+                matrix[row] = new Entity[degree + 1];
+                for (var k = 0; k <= degree; k++)
+                    matrix[row][k] = columns[k].Terms.FirstOrDefault(term => term.Key == monomials[row]).Value is { } c ? Number.Rational.Create(c) : Number.Integer.Zero;
+                rhs[row] = constantPart.Terms.FirstOrDefault(term => term.Key == monomials[row]).Value is { } d ? Number.Rational.Create(d.Negate()) : Number.Integer.Zero;
+            }
+            if (!Functions.PartialFractions.TrySolveLinear(matrix, rhs, out var values) || values is null)
+                return null;
+            Entity p = Number.Integer.Zero;
+            for (var k = 0; k <= degree; k++)
+                if (values[k].Evaled is Number.Complex { IsZero: false })
+                    p = p + values[k] * (k == 0 ? Number.Integer.One : k == 1 ? x : MathS.Pow(x, k));
+            if (p == Number.Integer.Zero)
+                return null;
+            var answer = (p * MathS.Pow(radicalBase, rPlusOne)).InnerSimplified;
+            // Checked as a fact about the function, since two functions taken for
+            // independent indeterminates may not be.
+            if (!Functions.PartialFractions.HoldsAtSampledPoints(answer.Differentiate(x), expr, x))
+                return null;
+            return answer;
+        }
+
         /// <summary>
         /// An integrand <c>e^h R</c>, with <c>R</c> rational in <c>x</c> and <c>h</c> rational
         /// in <c>x</c> or absent, integrated by the ansatz <c>F = e^h N/D</c>: <c>D</c> read

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

@@ -604,6 +604,10 @@ private static Entity Normalized(Entity expr) =>
             // an ansatz finds it or nothing does. After the substitutions, which answer the
             // linear-exponent cases in their own terms.
             if ((answer = IndefiniteIntegralSolver.SolveByExponentialAnsatz(expr, x)) is { }) return answer;
+            // A polynomial times a fractional power of a base that holds a function of x,
+            // answered as a polynomial times the next power of the base: a linear system in
+            // the polynomial's coefficients, exact, and volunteered at any depth.
+            if ((answer = IndefiniteIntegralSolver.SolveByAPolynomialTimesAPowerOfTheBase(expr, x)) is { }) return answer;
             // And an exponential times a rational function of the sine and cosine, closed by
             // the same kind of ansatz in the half-angle tangent: `e^x (1 - sin(x))/(1 - cos(x))`
             // is `-e^x cot(x/2)`.

Sources/Tests/UnitTests/Calculus/PolynomialTimesAPowerOfTheBaseTest.cs

@@ -0,0 +1,103 @@
+//
+// 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>
+    /// <c>N B^r</c>, for a fractional <c>r</c> and an <c>N</c> and a <c>B</c> that are polynomials in
+    /// <c>x</c> and the functions of it in them, answered as <c>P(x) B^(r + 1)</c> where there is a
+    /// polynomial <c>P</c>: <c>P' B + (r + 1) P B' = N</c> is a linear system in the coefficients of
+    /// <c>P</c> once the functions of <c>x</c> are taken for indeterminates, and it is solved exactly.
+    /// <a href="https://github.com/asc-community/AngouriMath/issues/718">#718</a>
+    /// </summary>
+    /// <remarks>
+    /// Bronstein's <c>(5x^2 + 3 (e^x + x)^(1/3) + e^x (3x + 2x^2))/(x (e^x + x)^(1/3))</c> is
+    /// <c>3/x + (5x + e^x (2x + 3))/(e^x + x)^(1/3)</c>, and the second is <c>3x (e^x + x)^(2/3)</c>;
+    /// nothing else read it, since the substitution <c>u = e^x + x</c> wants <c>1 + e^x</c> beside
+    /// the root and finds <c>5x + e^x (2x + 3)</c>.
+    /// </remarks>
+    [Trait("Area", "Calculus")]
+    public sealed class PolynomialTimesAPowerOfTheBaseTest
+    {
+        private static Entity DifferentiatesBack(string integrand, double[] points, params (string, double)[] pins)
+        {
+            var integral = integrand.ToEntity().Integrate("x");
+            Assert.DoesNotContain("integral(", integral.Stringize());
+            Assert.DoesNotContain("NaN", integral.Stringize());
+
+            Entity answer = integral.Substitute("C", 0);
+            Entity original = integrand.ToEntity();
+            foreach (var (name, value) in pins)
+            {
+                answer = answer.Substitute(name, value);
+                original = original.Substitute(name, value);
+            }
+            // With the symbols pinned the piecewise's conditions are decidable, and the
+            // simplification picks the arm before the derivative is taken.
+            var derivative = answer.Simplify().Differentiate("x");
+
+            var compared = 0;
+            foreach (var at in points)
+            {
+                var got = derivative.Substitute("x", at).EvalNumerical();
+                var want = original.Substitute("x", at).EvalNumerical();
+                if (got.IsNaN || want.IsNaN)
+                    continue;
+                compared++;
+                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));
+                Assert.True(difference / scale < 1e-8,
+                    $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, "
+                    + $"where the integrand is {want}");
+            }
+            Assert.True(compared >= 3,
+                $"only {compared} of {points.Length} points were comparable for {integrand}");
+            return integral;
+        }
+
+        /// <summary>
+        /// The exponential, a logarithm and a sine in the base, the root above the bar and
+        /// below it, and Bronstein's whole integrand, whose <c>3/x</c> the split takes off first.
+        /// </summary>
+        [Theory]
+        [InlineData("(5*x + e^x*(2*x + 3))/(e^x + x)^(1/3)", new[] { -1.5, -0.4, 0.3, 0.9, 2.6 })]
+        [InlineData("(5*x^2 + 3*(e^x + x)^(1/3) + e^x*(3*x + 2*x^2))/(x*(e^x + x)^(1/3))", new[] { 0.3, 0.9, 1.7, 2.6 })]
+        [InlineData("(5*x/2 + ln(x) + 3/2)*sqrt(x + ln(x))", new[] { 0.7, 0.9, 1.7, 2.6 })]
+        [InlineData("(3*x^3 + 2*x*sin(x) + x^2*cos(x)/2)/sqrt(x^2 + sin(x))", new[] { 0.5, 0.9, 1.7, 2.6 })]
+        public void APolynomialTimesTheNextPowerOfTheBase(string integrand, double[] points)
+            => DifferentiatesBack(integrand, points);
+
+        /// <summary>
+        /// The answer is the one it says: <c>3x (e^x + x)^(2/3)</c>, exactly, not a longer
+        /// spelling of it.
+        /// </summary>
+        [Fact]
+        public void TheAnswerIsTheProductItself()
+            => Assert.Equal("3 * x * (e ^ x + x) ^ (2/3) + C", "(5*x + e^x*(2*x + 3))/(e^x + x)^(1/3)".ToEntity().Integrate("x").Stringize());
+
+        /// <summary>
+        /// What is not of the form is declined by the system, not answered: a numerator one
+        /// coefficient off is no polynomial's derivative. Either verdict but a wrong answer.
+        /// </summary>
+        [Theory]
+        [InlineData("(2*x*ln(x) + x + ln(x))*sqrt(x + ln(x))")]
+        [InlineData("(3*x^2 + 2*x*sin(x) + x^2*cos(x))/sqrt(x^2 + sin(x))")]
+        public void WhatIsNotOfTheFormIsNotAnswered(string integrand)
+        {
+            var integral = integrand.ToEntity().Integrate("x");
+            Assert.DoesNotContain("NaN", integral.Stringize());
+            if (integral.Stringize().Contains("integral("))
+                return;   // unanswered is a legitimate verdict; a wrong answer is not
+            DifferentiatesBack(integrand, new[] { 0.7, 0.9, 1.7, 2.6 });
+        }
+    }
+}