From fe0887ef7ffb298a7caf06aacbb8625305e4616a Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Fri, 16 Oct 2026 06:22:56 +0000 Subject: [PATCH 1/3] Integrate a power of an inverse function of c + d x^2 by its recurrence For u the inverse sine or cosine of w = c + d x^2 with c^2 = 1, its hyperbolic cosine with c^2 = 1, or its hyperbolic sine with c^2 = -1, integral (a + b u)^n is x (a + b u)^n + 2 e n b H (a + b u)^(n - 1) - 4 e b^2 n (n - 1) integral (a + b u)^(n - 2), with H = 1/(d x u'(w)), and it closes at the first two powers. Rubi's 5.1.5, 5.2.5 and 7.2.5 (a + b acos(1 + d x^2))^4 and the like were declined. Part of https://github.com/asc-community/AngouriMath/issues/718 Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- .../Integration/IndefiniteIntegralSolver.cs | 79 +++++++++++++++++++ .../Integration/Integration.Definition.cs | 3 + ...InverseFunctionOfAQuadraticIntegralTest.cs | 68 ++++++++++++++++ 3 files changed, 150 insertions(+) create mode 100644 Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 5a0a8ccab..c6a85df42 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -28638,6 +28638,85 @@ static bool Through(Entity node, Entity.Variable x) } } + /// + /// A whole power of a + b u, where u is the inverse sine or cosine of + /// w = c + d x^2 with c^2 = 1, the inverse hyperbolic cosine of one, or the + /// inverse hyperbolic sine of one with c^2 = -1, by the recurrence + /// integral (a + b u)^n = x (a + b u)^n + 2 e n b H (a + b u)^(n - 1) - 4 e b^2 n (n - 1) integral (a + b u)^(n - 2), + /// with H = 1/(d x u'(w)), e = 1 for the circular functions and -1 for the + /// hyperbolic ones. It closes at n = 0 and n = 1, so nothing is asked below. + /// + /// + /// Each u has u'' = e w u'^3, and with c^2 = 1 (or -1) the radicand of + /// u' is d x^2 (2 c + d x^2) up to sign, so H u' = 2 and + /// (x (a + b u)^n + 2 e n b H (a + b u)^(n - 1))' is (a + b u)^n plus a multiple of + /// (a + b u)^(n - 2). Rubi's 5.1.5, 5.2.5 and 7.2.5 (a + b acos(1 + d x^2))^4 + /// and the like were declined. + /// https://github.com/asc-community/AngouriMath/issues/718 + /// + internal static Entity? SolveAPowerOfAnInverseFunctionOfAQuadraticByTheRecurrence(Entity expr, Entity.Variable x) + { + var (@base, power) = expr is Powf(var b0, Number.Integer p0) ? (b0, p0) : (expr, Number.Integer.One); + if (power.EInteger.Sign <= 0 || power.EInteger.CompareTo(EInteger.FromInt32(12)) > 0) + return null; + var n = power.EInteger.ToInt32Unchecked(); + // The inverse function, its argument w, the sign e, the radical H is written with, and + // whether w^2 is 1 or -1 at x = 0. + (Entity Node, Entity W, int Sign, Entity Root, int Square)? Read(Entity node) + { + switch (node) + { + case Arcsinf(var w): + return (node, w, 1, MathS.Sqrt(1 - MathS.Sqr(w)), 1); + case Arccosf(var w): + return (node, w, 1, -MathS.Sqrt(1 - MathS.Sqr(w)), 1); + // The parser writes acosh(w) and asinh(w) as ln(w + sqrt(w^2 -+ 1)). + case Logf(var logBase, Sumf(var w, Powf(var radicand, Number.Rational half))) + when logBase == MathS.e && half == Number.Rational.Create(1, 2) && w.ContainsNode(x): + foreach (var square in new[] { 1, -1 }) + if (VanishesIdentically((radicand - (MathS.Sqr(w) - square)).Expand())) + return (node, w, -1, MathS.Sqrt(radicand), square); + return null; + default: + return null; + } + } + var candidates = @base.Nodes.Select(Read).Where(read => read is not null).Select(read => read!.Value).ToList(); + if (candidates.Count == 0 || candidates.Any(read => read.Node != candidates[0].Node)) + return null; + var (inverse, argument, sign, root, squareAtZero) = candidates[0]; + var t = Variable.CreateUnique(expr, "u_inv"); + var inT = @base.Substitute(inverse, t); + if (inT.ContainsNode(x) || !TreeAnalyzer.TryGetPolyLinear(inT, t, out var b, out var a) + || Functions.PartialFractions.IsZeroAsAValue(b)) + return null; + if (!TreeAnalyzer.TryGetPolynomial(argument, x, out var read) || read.Keys.Any(key => !key.IsZero && !key.Equals(EInteger.FromInt32(2))) + || !read.TryGetValue(EInteger.FromInt32(2), out var d) || d.ContainsNode(x) || Functions.PartialFractions.IsZeroAsAValue(d)) + return null; + var c = read.TryGetValue(EInteger.Zero, out var c0) ? c0 : Number.Integer.Zero; + if (c.ContainsNode(x) || !Functions.PartialFractions.IsZeroAsAValue((MathS.Sqr(c) - squareAtZero).Expand().InnerSimplified)) + return null; + var h = root / (d * x); + var e = (Entity)Number.Integer.Create(sign); + var a_bu = a + b * inverse; + // I(k) by the recurrence, from I(0) = x and I(1) = x (a + b u) + 2 e b H. + var integrals = new Entity[n + 1]; + for (var k = 0; k <= n; k++) + { + Entity term = k == 0 ? x : x * MathS.Pow(a_bu, k) + 2 * e * k * b * h * (k == 1 ? Number.Integer.One : MathS.Pow(a_bu, k - 1)); + if (k >= 2) + term -= 4 * e * MathS.Sqr(b) * k * (k - 1) * integrals[k - 2]; + integrals[k] = term; + } + // The recurrence is exact where H u' = 2 with the derivative the library takes, the one + // step that depends on a branch; checked here rather than at sampled points, which for + // `c = 1` and a positive `d` put `w` past 1, on the cut of the inverse sine and cosine, + // and declined `asin(1 + d x^2)^2`, right for the `d < 0` where it is real. + if (!Functions.PartialFractions.IsZeroAsAValue((h * inverse.Differentiate(x) - 2).InnerSimplified)) + return null; + return integrals[n].InnerSimplified; + } + /// /// A power of a constant multiple of a radicand that a logarithm holds under a root, /// written over that radicand: acosh(c x) is ln(c x + sqrt(c^2 x^2 - 1)), diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index 606c7cfd4..5ee4e9a82 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -661,6 +661,9 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // dropped it re-enabled integration by parts one level below the call that // switched it off -- which is a cycle, since by parts calls back into here. // `x * ln(x)` went round it until the stack ran out. + // A whole power of a linear in the inverse sine, cosine or hyperbolic cosine of + // `c + d x^2` with `c^2 = 1`, by its recurrence, which closes: by parts it went round. + if ((answer = IndefiniteIntegralSolver.SolveAPowerOfAnInverseFunctionOfAQuadraticByTheRecurrence(expr, x)) is { }) return answer; // A product of powers of quotients of linears with one linear in two of them, written // over each linear once with the constant that takes in front. First, since it answers // the question asked only, and the polynomial term below writes `1/u` as `u^(-1)` and diff --git a/Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs new file mode 100644 index 000000000..45b735322 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs @@ -0,0 +1,68 @@ +// +// 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 whole power of a + b u for u the inverse sine, cosine or hyperbolic cosine of + /// c + d x^2 with c^2 = 1, or the inverse hyperbolic sine of one with c^2 = -1, + /// by the recurrence that closes at the first two powers: Rubi's 5.1.5, 5.2.5 and 7.2.5. Pinned + /// with d of either sign, on the real domain of each and off it. + /// #718 + /// + [Trait("Area", "Calculus")] + public sealed class APowerOfAnInverseFunctionOfAQuadraticIntegralTest + { + 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}"); + } + } + + // d with the sign that keeps the argument real at every sampled x, and |d| below 2/2.3^2 + // so that the inverse sine and cosine stay inside [-1, 1]; off it the exact checks hold too. + [Theory] + [InlineData("(a+b*acos(1+d*x^2))^4", -1)] + [InlineData("(a+b*acos(-1+d*x^2))^2", 1)] + [InlineData("(a+b*asin(1+d*x^2))^3", -1)] + [InlineData("(a+b*asin(-1+d*x^2))", 1)] + [InlineData("(a+b*acosh(-1+d*x^2))^2", -1)] + [InlineData("(a+b*acosh(1+d*x^2))^3", 1)] + [InlineData("(a+b*asinh(i+d*x^2))^2", 1)] + public void IsAnsweredAndDifferentiatesBack(string integrand, int sign) + { + DifferentiatesBack(integrand, ("a", 0.6), ("b", 1.25), ("d", 0.3 * sign)); + DifferentiatesBack(integrand, ("a", -0.75), ("b", -0.7), ("d", 0.15 * sign)); + } + + [Fact] + public void WithNumbersIsAnswered() + => DifferentiatesBack("asin(1-2*x^2)^2"); + } +} From 0da65f1e88236f9a68cd36346bdd800e1b097447 Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Fri, 16 Oct 2026 10:39:45 +0000 Subject: [PATCH 2/3] Integrate a negative power of an inverse function of c + d x^2 in sine and cosine integrals w - c is a constant times g(u/2 + beta)^2 for a sine or a cosine g, so x/g is constant on each side of 0 and dx = (x/g)(g'/2) du. The powers -1 and -2 are then sine and cosine integrals of (a + b u)/(2 b), hyperbolic ones for the inverse hyperbolic cosine, and the recurrence run backwards gives the rest. Rubi's 5.1.5, 5.2.5 and 7.2.5 1/(a + b acos(1 + d x^2)) and the like, to the first three negative powers, were declined. Part of https://github.com/asc-community/AngouriMath/issues/718 Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- .../Integration/IndefiniteIntegralSolver.cs | 60 ++++++++++++++++++- ...InverseFunctionOfAQuadraticIntegralTest.cs | 15 +++++ 2 files changed, 73 insertions(+), 2 deletions(-) diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index c6a85df42..c802a36f8 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -28644,7 +28644,12 @@ static bool Through(Entity node, Entity.Variable x) /// inverse hyperbolic sine of one with c^2 = -1, by the recurrence /// integral (a + b u)^n = x (a + b u)^n + 2 e n b H (a + b u)^(n - 1) - 4 e b^2 n (n - 1) integral (a + b u)^(n - 2), /// with H = 1/(d x u'(w)), e = 1 for the circular functions and -1 for the - /// hyperbolic ones. It closes at n = 0 and n = 1, so nothing is asked below. + /// hyperbolic ones. It closes at n = 0 and n = 1, so nothing is asked below. A + /// negative power goes through u: w - c is a constant times g(u/2 + beta)^2 + /// for a sine or a cosine g, so x/g is constant on each side of 0, and the + /// powers -1 and -2 are sine and cosine integrals of (a + b u)/(2 b) + /// (hyperbolic ones for the hyperbolic cosine), the rest the recurrence run backwards. Exact on + /// the real domain; off it x/g can change sign across a cut. /// /// /// Each u has u'' = e w u'^3, and with c^2 = 1 (or -1) the radicand of @@ -28657,7 +28662,7 @@ static bool Through(Entity node, Entity.Variable x) internal static Entity? SolveAPowerOfAnInverseFunctionOfAQuadraticByTheRecurrence(Entity expr, Entity.Variable x) { var (@base, power) = expr is Powf(var b0, Number.Integer p0) ? (b0, p0) : (expr, Number.Integer.One); - if (power.EInteger.Sign <= 0 || power.EInteger.CompareTo(EInteger.FromInt32(12)) > 0) + if (power.EInteger.IsZero || power.EInteger.Abs().CompareTo(EInteger.FromInt32(12)) > 0) return null; var n = power.EInteger.ToInt32Unchecked(); // The inverse function, its argument w, the sign e, the radical H is written with, and @@ -28699,6 +28704,8 @@ static bool Through(Entity node, Entity.Variable x) var h = root / (d * x); var e = (Entity)Number.Integer.Create(sign); var a_bu = a + b * inverse; + if (n < 0) + return NegativePower(); // I(k) by the recurrence, from I(0) = x and I(1) = x (a + b u) + 2 e b H. var integrals = new Entity[n + 1]; for (var k = 0; k <= n; k++) @@ -28715,6 +28722,55 @@ static bool Through(Entity node, Entity.Variable x) if (!Functions.PartialFractions.IsZeroAsAValue((h * inverse.Differentiate(x) - 2).InnerSimplified)) return null; return integrals[n].InnerSimplified; + + // A negative power, through u: w - c is a constant times g(theta)^2 with theta = u/2 + beta, + // so x = kappa g(theta) with kappa = x/g(theta) constant on each side of 0, and + // dx = (kappa/2) g'(theta) du. Then I(-1) and I(-2) are sine and cosine integrals of + // v = (a + b u)/(2 b), and the recurrence run backwards gives the rest. + Entity? NegativePower() + { + var plus = Functions.PartialFractions.IsZeroAsAValue((c - 1).InnerSimplified); + // Whether g is a sine, beta, and h = g' as (hyperbolic, cosine, sign). + (bool GIsSine, Entity Beta, bool Hyperbolic, bool Cosine, int Sign)? table = (inverse, sign, squareAtZero, plus) switch + { + (Arccosf, _, _, true) => (true, Number.Integer.Zero, false, true, 1), + (Arccosf, _, _, false) => (false, Number.Integer.Zero, false, false, -1), + (Arcsinf, _, _, true) => (true, -MathS.pi / 4, false, true, 1), + (Arcsinf, _, _, false) => (false, -MathS.pi / 4, false, false, -1), + (Logf, -1, 1, true) => (true, Number.Integer.Zero, true, true, 1), + (Logf, -1, 1, false) => (false, Number.Integer.Zero, true, false, 1), + _ => null, + }; + if (table is not var (gIsSine, beta, hyperbolic, cosine, hSign)) + return null; + var theta = inverse / 2 + beta; + var kappa = x / (gIsSine ? Sin(theta) : Cos(theta)); + var v = a_bu / (2 * b); + var gamma = a / (2 * b) - beta; + Entity Sin(Entity z) => hyperbolic ? MathS.Hyperbolic.Sinh(z) : MathS.Sin(z); + Entity Cos(Entity z) => hyperbolic ? MathS.Hyperbolic.Cosh(z) : MathS.Cos(z); + // The integral over v of h(v - gamma)/v for h a cosine or a sine, circular or hyperbolic. + Entity OverV(bool ofCosine) => (hyperbolic, ofCosine) switch + { + (false, true) => Cos(gamma) * MathS.Ci(v) + Sin(gamma) * MathS.Si(v), + (false, false) => Cos(gamma) * MathS.Si(v) - Sin(gamma) * MathS.Ci(v), + (true, true) => Cos(gamma) * MathS.Chi(v) - Sin(gamma) * MathS.Shi(v), + _ => Cos(gamma) * MathS.Shi(v) - Sin(gamma) * MathS.Chi(v), + }; + // h' is the other function: -sin for cos, cos for sin; sinh for cosh, cosh for sinh. + var hPrimeSign = hyperbolic || !cosine ? hSign : -hSign; + var hOfTheta = hSign * (cosine ? Cos(theta) : Sin(theta)); + var integrals = new Dictionary + { + [-1] = kappa / (2 * b) * hSign * OverV(cosine), + [-2] = kappa / 2 * (-hOfTheta / (b * a_bu) + hPrimeSign * OverV(!cosine) / (2 * MathS.Sqr(b))), + }; + for (var k = -1; k - 2 >= n; k--) + integrals[k - 2] = (x * MathS.Pow(a_bu, k) + 2 * e * k * b * h * MathS.Pow(a_bu, k - 1) - integrals[k]) + / (4 * e * MathS.Sqr(b) * k * (k - 1)); + var answer = integrals[n].InnerSimplified; + return Functions.PartialFractions.DerivativeDiffersWhereReal(answer, expr, x) ? null : answer; + } } /// diff --git a/Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs index 45b735322..0eb90e141 100644 --- a/Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs +++ b/Sources/Tests/UnitTests/Calculus/APowerOfAnInverseFunctionOfAQuadraticIntegralTest.cs @@ -61,6 +61,21 @@ public void IsAnsweredAndDifferentiatesBack(string integrand, int sign) DifferentiatesBack(integrand, ("a", -0.75), ("b", -0.7), ("d", 0.15 * sign)); } + // A negative power, in sine and cosine integrals (hyperbolic ones for acosh) of (a + b u)/(2 b): + // exact on the real domain, where d keeps the argument inside [-1, 1]. + [Theory] + [InlineData("1/(a+b*acos(1+d*x^2))", -1)] + [InlineData("1/(a+b*acos(-1+d*x^2))^2", 1)] + [InlineData("1/(a+b*asin(1+d*x^2))^3", -1)] + [InlineData("1/(a+b*asin(-1+d*x^2))", 1)] + [InlineData("1/(a+b*acosh(1+d*x^2))^2", 1)] + [InlineData("1/(a+b*acosh(-1+d*x^2))^3", -1)] + public void ANegativePowerIsAnsweredAndDifferentiatesBack(string integrand, int sign) + { + DifferentiatesBack(integrand, ("a", 0.6), ("b", 1.25), ("d", 0.3 * sign)); + DifferentiatesBack(integrand, ("a", -0.75), ("b", -0.7), ("d", 0.15 * sign)); + } + [Fact] public void WithNumbersIsAnswered() => DifferentiatesBack("asin(1-2*x^2)^2"); From d9c0f71a4afb55fcd3ddad14eca54a9aeba6eedf Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Fri, 16 Oct 2026 12:46:40 +0000 Subject: [PATCH 3/3] Record the negative powers in the breaking changes Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- BREAKING-CHANGES.md | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 60c75dea8..a34469655 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -396,6 +396,24 @@ in Rubi's 1.1.3.4 come out shorter as well, and those that took two seconds take | `"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 | +### A negative power of an inverse function of `c + d x^2` is integrated in sine and cosine integrals + +**Answers where there were none.** `1/(a + b acos(1 + d x^2))` was declined, and so were seventeen more of Rubi's +5.1.5, 5.2.5 and 7.2.5, to the first three negative powers. `w - c` is a constant times `g(u/2 + beta)^2` for a sine or +a cosine `g`, so `x/g` is constant on each side of 0; the powers `-1` and `-2` are then sine and cosine integrals of +`(a + b u)/(2 b)`, hyperbolic ones for the inverse hyperbolic cosine, and the recurrence run backwards gives the +rest. The answers hold on the real domain of the integrand +([#718](https://github.com/asc-community/AngouriMath/issues/718)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"1/(a+b*acos(1+d*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 184 characters | +| `"1/(a+b*acos(-1+d*x^2))^2".ToEntity().Integrate("x")` | `integral(...)` | 273 characters | +| `"1/(a+b*asin(1+d*x^2))^3".ToEntity().Integrate("x")` | `integral(...)` | 380 characters | +| `"1/(a+b*asin(-1+d*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 224 characters | +| `"1/(a+b*acosh(1+d*x^2))^2".ToEntity().Integrate("x")` | `integral(...)` | 608 characters | +| `"1/(a+b*acosh(-1+d*x^2))^3".ToEntity().Integrate("x")` | `integral(...)` | 633 characters | + ### The reciprocal of a binomial in a high power of a sine, cosine or hyperbolic one is split over its roots **Answers where there were none.** `1/(a + b sin(x)^5)` was declined, and so were twenty more of Rubi's 4.1.7, 4.2.7,