Integrate a symbolic secant power, a half-odd power of a ± a sin, and an even cosine power (#1975)

Three seams of Rubi's family 4:

- (b sec(y))^m tan(y)^3 was declined only because m was not a number; in
  u = sec(y) it is (b u)^m (u^2 - 1)/u, and the rule for an odd power of the
  tangent beside a function of the secant reads a power free of x whatever its
  exponent now. A negative power of the tangent beside it is hypergeometric and
  still declined.
- 1/(sqrt(a + a sin(y)) (c + d sin(y))^2) and its kin are answered by the half
  angle at which a + a sin(y) is a square. They never reached it: the polynomial
  term asked the substitution search first, which spent the budget.
- cos(y)^4 sin(y)^3/(a + b sin(y))^2: an even whole power of the cosine is a
  power of 1 - sin(y)^2, so beside a power of the sine above the bar the
  integrand is a polynomial in the sine over the power. Without a sine above
  the bar the rules after it answer as before: cos(y)^4/(a + b sin(y))^8 is 40K
  characters there and was 106K here.

Rubi pocket (2,061 rows of family 4): 1967 -> 1994 solved; the moves re-run
alone gain 28 rows (4.1.2.1 x16, 4.1.2.2 x7, 4.1.3.1:364, 4.3.0 x4), none lost.

Part of #233.


Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura

Co-authored-by: Claude Opus 5.5 (1M context) <noreply@anthropic.com>

Author: Rafael-SOWNet

BREAKING-CHANGES.md

@@ -350,6 +350,30 @@ cache is now keyed on the downcasting setting too
 |---|---|---|
 | `MathS.FromString("1/3").Evaled` under downcasting off, after the same parse with it on | `1/3` | `0.3333…` (one hundred digits) |
 
+### A symbolic power of the secant, a half-odd power of `a ± a sin` beside other linears, and an even power of the cosine over a sine's binomial are integrated
+
+**Answers where there were none.** `(b sec(g + f x))^m tan(g + f x)^3` and three more of Rubi's 4.3.0 were declined only
+because `m` was not a number: in `u = sec(y)` it is `(b u)^m (u^2 - 1)/u`, and the rule for an odd power of the tangent
+beside a function of the secant now reads a power free of x whatever its exponent.
+`1/(sqrt(a + a sin(y)) (c + d sin(y))^2)`, `1/((a + a sin(y))^2 sqrt(c - c sin(y)))` and fourteen more of 4.1.2.1 and
+4.1.3.1 are answered by the half angle at which `a + a sin(y)` is a square, which they never reached: the polynomial
+term wrote the reciprocal as a power and asked the substitution search, which spent the budget. A half-odd power of
+`a ± a sin(y)` beside at least one other linear in the same sine or cosine is asked it first now. And
+`cos(y)^4 sin(y)^3/(a + b sin(y))^2` and six more of 4.1.2.2 are read by the rule below, an even whole power of the
+cosine being a power of `1 - sin(y)^2` exactly, where a power of the sine stands above the bar beside it. That file's
+neighbours the unreleased master answered are respelled, some shorter and some longer ([#233](https://github.com/asc-community/AngouriMath/issues/233)).
+
+| Input | Was (2.5.0) | Now |
+|---|---|---|
+| `"(b*sec(g+f*x))^m*tan(g+f*x)^3".ToEntity().Integrate("x")` | `integral(...)` | 144 characters |
+| `"cot(g+f*x)^5*(b*csc(g+f*x))^m".ToEntity().Integrate("x")` | `integral(...)` | 221 characters |
+| `"1/(sqrt(a+a*sin(g+f*x))*(c+d*sin(g+f*x))^2)".ToEntity().Integrate("x")` | `integral(...)` | 2,679 characters |
+| `"1/((a+a*sin(g+f*x))^2*sqrt(c-c*sin(g+f*x)))".ToEntity().Integrate("x")` | `integral(...)` | 329 characters |
+| `"cos(c+d*x)^4*sin(c+d*x)^3/(a+b*sin(c+d*x))^2".ToEntity().Integrate("x")` | `integral(...)` | 4,174 characters |
+| `"cos(c+d*x)^6*sin(c+d*x)^3/(a+b*sin(c+d*x))^3".ToEntity().Integrate("x")` | `integral(...)` | 9,029 characters |
+| `"cos(c+d*x)^4*sin(c+d*x)^2/(a+b*sin(c+d*x))^2".ToEntity().Integrate("x")` | `integral(...)`; 7,068 characters, in 1.3 s, on the unreleased master | 3,677 characters, in 154 ms |
+| `"cos(c+d*x)^2*sin(c+d*x)^2/(a+b*sin(c+d*x))^3".ToEntity().Integrate("x")` | `integral(...)`; 4,621 characters, in 431 ms, on the unreleased master | 6,584 characters, in 570 ms |
+
 ### A polynomial in a sine over a power of a binomial in it is written in powers of the binomial
 
 **Answers where there were none.** `(A + B cos(c + d x) + C cos(c + d x)^2)/(a + b cos(c + d x))^4` ran out of time, and

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

@@ -5364,7 +5364,8 @@ Entity Integral(ERational power)
         /// </para>
         /// https://github.com/asc-community/AngouriMath/issues/718
         /// </remarks>
-        internal static Entity? SolveByTheHalfAngleWhereOnePlusASineIsASquare(Entity expr, Entity.Variable x, bool integrateByParts)
+        internal static Entity? SolveByTheHalfAngleWhereOnePlusASineIsASquare(Entity expr, Entity.Variable x, bool integrateByParts,
+            bool collectTheLinears = false)
         {
             if (!Integration.AnsweringTheQuestionAsked)
                 return null;
@@ -5445,6 +5446,41 @@ Entity Integral(ERational power)
                     constant = MathS.Pow(a, exponent) * constant;
                 return constant * MathS.Pow(square, Number.Integer.Create(twoP));
             });
+            // `c + d sin(y)` as `(c - d) + 2d sin(u)^2`, collected, where asked: written
+            // `c + d (2 sin(u)^2 - 1)`, `sin(u)^5/(c + d (2 sin(u)^2 - 1))^2` ran past the budget.
+            if (collectTheLinears)
+            {
+                var (function, square) = sineKind ? (sine, MathS.Sqr(sinU)) : (cosine, MathS.Sqr(cosU));
+                rewritten = rewritten.Replace(node =>
+                {
+                    if (node is not Sumf || !node.ContainsNode(function))
+                        return node;
+                    Entity constant = Number.Integer.Zero, coefficient = Number.Integer.Zero;
+                    foreach (var term in Sumf.LinearChildren(node))
+                    {
+                        if (!term.ContainsNode(x))
+                        {
+                            constant += term;
+                            continue;
+                        }
+                        Entity factor = Number.Integer.One;
+                        var seen = false;
+                        foreach (var part in Mulf.LinearChildren(term))
+                            if (part == function && !seen)
+                                seen = true;
+                            else if (part.ContainsNode(x))
+                                return node;
+                            else
+                                factor *= part;
+                        if (!seen)
+                            return node;
+                        coefficient += factor;
+                    }
+                    var shifted = (constant - coefficient).InnerSimplified;
+                    var collected = (2 * coefficient).InnerSimplified * square;
+                    return VanishesIdentically(shifted) ? collected : shifted + collected;
+                });
+            }
             var inU = rewritten.Substitute(sine, sineInU).Substitute(cosine, cosineInU);
             // What is left of the variable is the variable itself: the angle is `rate x + offset`
             // and `u` is half of it (a quarter turn on for the sine), so `x` is
@@ -5468,6 +5504,82 @@ Entity Integral(ERational power)
             return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer;
         }
 
+        /// <summary>
+        /// A half-odd power of <c>a ± a sin(y)</c> times whole or half-odd powers of other linears
+        /// in the same sine, <c>c + d sin(y)</c>, and nothing else of <c>x</c>: by
+        /// <see cref="SolveByTheHalfAngleWhereOnePlusASineIsASquare"/>, asked first. The cosine's
+        /// the same.
+        /// </summary>
+        /// <remarks>
+        /// <para>
+        /// Rubi's <c>1/(sqrt(a + a sin(e + f x)) (c + d sin(e + f x))^2)</c>,
+        /// <c>1/((a + a sin(e + f x))^2 sqrt(c - c sin(e + f x)))</c> and their kin in 4.1.2.1
+        /// ran past the budget. The half angle answers each, its question in <c>u</c> taking
+        /// tens of milliseconds, but it was not reached: the polynomial term above it wrote the
+        /// reciprocal as a power, took the constant out of the root, and asked the substitution
+        /// search below that, whose simplification of what it found took the whole budget.
+        /// </para>
+        /// <para>
+        /// The answer is that rule's, so exact on the same terms. Read narrowly: a product of
+        /// powers of linears in one sine or one cosine, one of them half-odd and of
+        /// <c>a ± a sin(y)</c>, and at least two of them, since one alone is
+        /// <see cref="SolveAHalfPowerOfOnePlusASine"/>'s.
+        /// https://github.com/asc-community/AngouriMath/issues/233
+        /// </para>
+        /// </remarks>
+        internal static Entity? SolveAHalfOddPowerOfOnePlusMinusASineBesidePowersOfLinearsInIt(Entity expr, Entity.Variable x, bool integrateByParts)
+        {
+            if (!Integration.AnsweringTheQuestionAsked)
+                return null;
+            if (!TryReadAsQuotient(expr, out var above, out var below))
+                (above, below) = (expr, Number.Integer.One);
+            Entity? function = null;
+            var (factors, halfOdd) = (0, false);
+            foreach (var side in new[] { above, below })
+                foreach (var factor in Mulf.LinearChildren(side))
+                {
+                    if (!factor.ContainsNode(x))
+                        continue;
+                    var (@base, exponent) = factor is Powf(var inner, Number.Rational power) ? (inner, power.ERational) : (factor, ERational.One);
+                    var isHalfOdd = exponent.Denominator.Equals(EInteger.FromInt32(2));
+                    if (!isHalfOdd && !exponent.Denominator.Equals(EInteger.One))
+                        return null;
+                    // c + d f(y), the one f throughout, and c there: a power of f itself is not this shape.
+                    Entity? found = null;
+                    var constant = false;
+                    foreach (var term in Sumf.LinearChildren(@base))
+                    {
+                        if (!term.ContainsNode(x))
+                        {
+                            constant = true;
+                            continue;
+                        }
+                        foreach (var part in Mulf.LinearChildren(term))
+                            if (part.ContainsNode(x))
+                            {
+                                if (found is not null || part is not (Sinf or Cosf))
+                                    return null;
+                                found = part;
+                            }
+                    }
+                    if (found is null || !constant || function is not null && found != function)
+                        return null;
+                    function = found;
+                    if (isHalfOdd)
+                    {
+                        if (ReadAsOnePlusMinusAFunction(@base, function is Sinf ? function : MathS.Sin(function.DirectChildren.First()),
+                                function is Cosf ? function : MathS.Cos(function.DirectChildren.First())) is null)
+                            return null;
+                        halfOdd = true;
+                    }
+                    factors++;
+                }
+            if (!halfOdd || factors < 2 || function is null
+                || !TreeAnalyzer.TryGetPolyLinear(function.DirectChildren.First(), x, out var rate, out _) || rate.ContainsNode(x))
+                return null;
+            return SolveByTheHalfAngleWhereOnePlusASineIsASquare(expr, x, integrateByParts, collectTheLinears: true);
+        }
+
         /// <summary>
         /// Powers of <c>a ± a sin(y)</c> with a symbolic exponent, or beside a fractional power of
         /// <c>g cos(y)</c>, and anything else in the sine alone, under <c>u = sin(y)</c>:
@@ -6154,8 +6266,12 @@ static bool IsHalfOdd(Entity exponent) => exponent is Number.Rational rational a
         /// <remarks>
         /// Rubi's <c>cot(c + d x)/(a + b sec(c + d x))^(3/2)</c> was declined and
         /// <c>cot(c + d x)^3 sqrt(a + b sec(c + d x))</c> ran past thirty seconds, where in the
-        /// secant each is a rational function beside the root of a linear.
+        /// secant each is a rational function beside the root of a linear. A power free of
+        /// <c>x</c> is read whatever its exponent: Rubi's <c>(b sec(e + f x))^m tan(e + f x)^3</c>,
+        /// 4.3.0, was declined only because <c>m</c> was not a number, where in <c>u</c> it is
+        /// <c>(b u)^m (u^2 - 1)/u</c>.
         /// https://github.com/asc-community/AngouriMath/issues/718
+        /// https://github.com/asc-community/AngouriMath/issues/233
         /// </remarks>
         internal static Entity? SolveAnOddPowerOfTheTangentBesideAFunctionOfTheSecant(Entity expr, Entity.Variable x, bool integrateByParts)
         {
@@ -6174,10 +6290,24 @@ static bool IsHalfOdd(Entity exponent) => exponent is Number.Rational rational a
             if (argument is null || !TreeAnalyzer.TryGetPolyLinear(argument, x, out var rate, out _)
                 || rate.ContainsNode(x) || TreeAnalyzer.IsZero(rate))
                 return null;
-            if (!expr.Nodes.All(node => !node.ContainsNode(x)
-                    || node is Variable or Sumf or Minusf or Mulf or Divf or Sinf or Cosf or Secantf or Cosecantf or Tanf or Cotanf
-                    || node is Powf(_, Number.Rational)))
-                return null;
+            // A symbolic power too, of a product of constants and secants or cosecants:
+            // `(b sec(y))^m` is `(b u)^m`. Of a sum it brings `(a + b u)^m/u`, not elementary.
+            var z = Variable.CreateUnique(expr, "z_secant");
+            bool IsAProductOfSecants(Entity @base)
+            {
+                var inZ = @base.Replace(node => node is Secantf(var a) && a == argument || node is Cosecantf(var b) && b == argument ? z : node);
+                return !inZ.ContainsNode(x) && inZ.Nodes.All(node => !node.ContainsNode(z) || node == z || node is Mulf or Divf or Powf(_, Number.Rational));
+            }
+            var symbolic = false;
+            foreach (var node in expr.Nodes)
+            {
+                if (!node.ContainsNode(x) || node is Variable or Sumf or Minusf or Mulf or Divf or Sinf or Cosf or Secantf or Cosecantf or Tanf or Cotanf
+                    || node is Powf(_, Number.Rational))
+                    continue;
+                if (node is not Powf(var @base, var exponent) || exponent.ContainsNode(x) || !IsAProductOfSecants(@base))
+                    return null;
+                symbolic = true;
+            }
             var t = Variable.CreateUnique(expr, "t_odd_tangent");
             var u = Variable.CreateUnique(expr, "u_secant");
             foreach (var ofTheSecant in new[] { true, false })
@@ -6193,6 +6323,9 @@ static bool IsHalfOdd(Entity exponent) => exponent is Number.Rational rational a
                 var n = exponents[0];
                 if (rest.ContainsNode(t) || !n.IsInteger() || n.ToEInteger().IsEven || !n.ToEInteger().CanFitInInt32())
                     continue;
+                // Beside a symbolic power a negative one leaves `u^m/(u (u^2 - 1)^k)`, hypergeometric.
+                if (symbolic && n.Sign < 0)
+                    continue;
                 var inU = rest.Replace(node => node switch
                 {
                     Secantf(var a) when ofTheSecant && a == argument => u,
@@ -27949,22 +28082,57 @@ when two.Equals(EInteger.FromInt32(2)) && e1 == MathS.e && e2 == MathS.e
         /// The rules after it took the whole quotient through the half-angle and wrote 30 million
         /// characters for <c>(A + B cos(x) + C cos(x)^2)/(a + b cos(x))^4</c> in eighteen seconds, where
         /// <c>1/(a + b cos(x))^4</c> alone is 44 thousand. Rubi's 4.1.4 and 4.2.4
-        /// <c>(a + b cos)^m (A + B cos + C cos^2)</c> at <c>m = -4</c> timed out.
+        /// <c>(a + b cos)^m (A + B cos + C cos^2)</c> at <c>m = -4</c> timed out. An even whole power
+        /// of the other function above the bar is a power of <c>1 - f^2</c>, exactly, and read as
+        /// one: Rubi's 4.1.2.2 <c>cos(c + d x)^6 sin(c + d x)^3/(a + b sin(c + d x))^2</c> and its
+        /// kin timed out too, where in the sine each is the polynomial over the power.
         /// https://github.com/asc-community/AngouriMath/issues/718
+        /// https://github.com/asc-community/AngouriMath/issues/233
         /// </remarks>
         internal static Entity? SolveAPolynomialInASineOverAPowerOfABinomialInItInPowersOfTheBinomial(Entity expr, Entity.Variable x)
         {
             if (!Integration.AnsweringTheQuestionAsked)
                 return null;
             var functions = expr.Nodes.Where(node => node is Sinf or Cosf && node.ContainsNode(x)).Distinct().ToList();
-            if (functions.Count != 1)
+            if (functions.Count is not (1 or 2))
+                return null;
+            var s = Variable.CreateUnique(expr, "s_pow");
+            Entity? function = null, inS = null;
+            foreach (var candidate in functions)
+            {
+                // The other of the two only to even whole powers above the bar, each a power of
+                // 1 - f^2: `cos(y)^4 sin(y)^3/(a + b sin(y))^2` is a polynomial in the sine over the power.
+                var written = candidate;
+                if (functions.Count == 2)
+                {
+                    var other = functions.Single(f => f != candidate);
+                    if (candidate is Sinf && other is not Cosf || candidate is Cosf && other is not Sinf
+                        || other.DirectChildren.First() != candidate.DirectChildren.First())
+                        return null;
+                    // And only beside the function itself above the bar: `cos(y)^4/(a + b sin(y))^8`, with
+                    // no sine above it, is answered by the rules after it in 40 thousand characters, and
+                    // was 106 thousand here, in 60% more time.
+                    var (top, _) = Functions.SingleQuotient.Of(expr);
+                    if (!top.ContainsNode(candidate))
+                        continue;
+                    written = expr.Replace(node => node is Powf(var @base, Number.Integer { EInteger: var power }) && @base == other
+                        && power.Sign > 0 && power.IsEven && power.CompareTo(EInteger.FromInt32(12)) <= 0
+                        ? MathS.Pow(1 - MathS.Sqr(candidate), Number.Integer.Create(power / 2))
+                        : node);
+                    if (written.ContainsNode(other))
+                        continue;
+                }
+                else
+                    written = expr;
+                function = candidate;
+                inS = written.Substitute(candidate, s);
+                break;
+            }
+            if (function is null || inS is null)
                 return null;
-            var function = functions[0];
             var argument = function switch { Sinf(var y) => y, Cosf(var y) => y, _ => null! };
             if (!TreeAnalyzer.TryGetPolyLinear(argument, x, out var slope, out _) || Functions.PartialFractions.IsZeroAsAValue(slope))
                 return null;
-            var s = Variable.CreateUnique(expr, "s_pow");
-            var inS = expr.Substitute(function, s);
             if (inS.ContainsNode(x))
                 return null;
             var (above, below) = Functions.SingleQuotient.Of(inS);

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

@@ -776,6 +776,9 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr)
             // A polynomial beside error functions of one shifted linear, written in it: by parts in
             // x leaves a Gaussian of the shifted linear beside a polynomial, which nothing reads.
             if ((answer = IndefiniteIntegralSolver.SolveAPolynomialBesideErrorFunctionsOfAShiftedLinear(expr, x, integrateByParts)) is { }) return answer;
+            // A half-odd power of a ± a sin(y) beside powers of other linears in the sine, by the half
+            // angle, before the polynomial term: beneath it the substitution search ran out of time.
+            if ((answer = IndefiniteIntegralSolver.SolveAHalfOddPowerOfOnePlusMinusASineBesidePowersOfLinearsInIt(expr, x, integrateByParts)) is { }) return answer;
             if ((answer = IndefiniteIntegralSolver.SolveAsPolynomialTerm(expr, x, integrateByParts)) is { }) return answer;
             // A product of powers of the variable with a power of a constant multiple of it
             // among them, `(c x)^m x^n`, by the power rule with the written power kept as