diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 671eff440..db6f97ad8 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -5509,6 +5509,25 @@ the crash harness found. Both now agree with the solver | `"{ x : 1/0 = x }".ToEntity().Evaled` | `NaN` | `{ }` | | `"x subset 2".ToEntity().Solve("x")` | `x * subset ^ 2` had no solution: `{ }` — `subset` is new since | `{ }`, where the unreleased keyword threw | +### An inverse hyperbolic function beside a factor that is not a polynomial is integrated by parts against its logarithm + +**Answers where there were none.** `x (a + b acosh(c x))/(d + e x^2)^2`, `x (a + b acsch(c x))/(d + e x^2)^3`, +`asech(a + b x)/x^3` and `acsch(a + b x)/x^3`, Rubi's 7.2 to 7.6, were declined. The parser writes each inverse +hyperbolic function as a logarithm, `ln(A + sqrt(A^2 + k))` or `ln((1 + u)/(1 - u))`, and differentiated as written that +is a quotient with the root in both parts, which the search never brought down to the root alone. Its derivative is now +written as the algebraic function it is, `A'/sqrt(A^2 + k)` or `2 u'/(1 - u^2)`, and `P + Q L` is integrated as +`integral P + G L - integral G L'` with `G` the antiderivative of `Q`. Only where `Q` is rational and not a polynomial, +beside which the function is answered already; not where a factor below `Q`'s bar divides the function's own quadratic, +as in `atanh(a x)/(1 - a^2 x^2)^3`, which the substitution under the function answers; and not where `G` holds `L` +itself ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"x*(a + b*acosh(c*x))/(d + k*x^2)^2".ToEntity().Integrate("x")` | `integral(...)` | 1,624 characters | +| `"x*(a + b*acsch(c*x))/(d + k*x^2)^3".ToEntity().Integrate("x")` | `integral(...)` | 2,596 characters | +| `"asech(a + b*x)/x^3".ToEntity().Integrate("x")` | `integral(...)` | 2,831 characters | +| `"acsch(a + b*x)/x^3".ToEntity().Integrate("x")` | `integral(...)` | 2,746 characters | + ### `asech(a + b x)` and `acsch(a + b x)` beside a power of `x` are integrated By parts they leave `1/(x (a + b x)^2 sqrt(1/(a + b x)^2 - 1))`, and the chain never read it: diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index c91066d1e..3fd08cfbb 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -20738,6 +20738,159 @@ Entity TakeOut(Entity node) } } + /// + /// By parts against an inverse hyperbolic function written as the logarithm the parser makes of + /// it, the logarithm's derivative written as the algebraic function it is: + /// ln(A + sqrt(A^2 + k))' = A'/sqrt(A^2 + k) for a constant k -- asinh, + /// acosh, asech and acsch -- and ln((1 + u)/(1 - u))' = 2 u'/(1 - u^2), + /// the same for (u + 1)/(u - 1) -- atanh and acoth. With the integrand + /// P + Q L for the logarithm L and G = integral Q, it is + /// integral P + G L - integral G L'. Rubi's 7.1-7.6, x (a + b acosh(c x))/(d + e x^2)^2 + /// and asech(a + b x)/x^3, which timed out. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + /// Differentiated as written, the logarithm is a quotient with the root in both parts, which the + /// search never brought down to the root alone. The derivatives are identities of the algebra + /// with no branch in them. Only for a Q that is rational and not a polynomial: beside a + /// polynomial the inverse function is answered already, and asked here would come back respelled; + /// beside a root the rules for the inverse function answer, and asking the root's integral here + /// first took the time they answer in. Not + /// re-entered, and not where G holds L itself, whose question is the one asked: + /// asinh(x)/sqrt(x^2 + 1) is L L'. Nor where the answer holds on part of the + /// line only, an inequality in x among its conditions: Rubi's 7.5.1 (a + b asech(c x))/x^n + /// is answered for every x without this rule. Nor where a factor below Q's bar divides the + /// function's own quadratic, 1 - u^2 or the radicand, as in atanh(a x)/(1 - a^2 x^2)^3: + /// the substitution under the function answers those, and by parts here ran out of time. And an + /// artanh or arcoth only of a linear: beside atanh(x sqrt(e)/sqrt(d + e x^2))/x^n the + /// remainder was slow to decline. + /// + internal static Entity? SolveByPartsAgainstAnInverseHyperbolicFunctionAsItsLogarithm(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (againstAnInverseHyperbolicFunction) + return null; + // The one logarithm of x of either shape, its derivative, and the function's own quadratic. + Entity? logarithm = null, derivative = null; + List? quadratic = null; + foreach (var node in expr.Nodes) + { + if (node is not Logf(var @base, var argument) || @base != MathS.e || !argument.ContainsNode(x)) + continue; + if (ADerivativeOfTheLogarithm(argument, x) is not { } read) + continue; + if (logarithm is not null && logarithm != node) + return null; + (logarithm, derivative, quadratic) = (node, read.Derivative, read.Quadratic); + } + if (logarithm is null || derivative is null || quadratic is null) + return null; + // expr = P + Q L. + var t = Variable.CreateUnique(expr, "l_inv"); + var inT = expr.Replace(node => node == logarithm ? t : node); + if (inT.ContainsNode(logarithm) || !TreeAnalyzer.TryGetPolyLinear(inT, t, out var q, out var p) || TreeAnalyzer.IsZero(q) + || q.ContainsNode(t) || p.ContainsNode(t)) + return null; + if (TreeAnalyzer.TryGetPolynomial(q, x, out var asPolynomial) && asPolynomial.Keys.All(k => k.Sign >= 0)) + return null; + // And only for a rational one: beside a root, `(a + b asinh(c x)) sqrt(d + c^2 d x^2)` and the like, + // the rules for the inverse function answer at once, and the root's integral asked here took their time. + if (!IsARationalFunction(q, x)) + return null; + // Nor where a factor below Q's bar divides the function's own quadratic: the substitution + // under the function answers `atanh(a x)/(1 - a^2 x^2)^3` at once, and + // `atanh(a x)/(x^2 (1 - a^2 x^2)^2)`, where by parts here ran out of time. + foreach (var factor in Mulf.LinearChildren(Functions.SingleQuotient.Of(q).Denominator)) + { + var below = factor is Powf(var raised, Number.Integer { EInteger.Sign: > 0 }) ? raised : factor; + if (!below.ContainsNode(x)) + continue; + if (quadratic.Any(own => AreProportionalAtSampledPoints(below, own, x))) + return null; + // A linear one where the quadratic is zero at its root. + if (TreeAnalyzer.TryGetPolyLinear(below, x, out var slope, out var offset) && !VanishesIdentically(slope) + && VanishesIdentically(quadratic[0].Substitute(x, -offset / slope))) + return null; + } + againstAnInverseHyperbolicFunction = true; + try + { + if (Integration.ComputeAsAQuestionOfItsOwn(q, x, integrateByParts) is not { } g + || g.Nodes.Any(node => node is Integralf || node == MathS.NaN || node == logarithm)) + return null; + if (Integration.ComputeAsAQuestionOfItsOwn(g * derivative, x, integrateByParts) is not { } rest + || rest.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + Entity first = Number.Integer.Zero; + if (!TreeAnalyzer.IsZero(p)) + { + if (Integration.ComputeAsAQuestionOfItsOwn(p, x, integrateByParts) is not { } integralOfP + || integralOfP.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + first = integralOfP; + } + // Not an answer for part of the line only, where the function has one for all of it: + // beside 1/x^2 the remainder for `asech(c x)` was answered where `1 - c x >= 0`, and + // without this rule `(a + b asech(c x))/x^2` is answered for every x. An even power + // above zero leaves out a point, not a part: `(a + b x)^2 > 0`. + var answer = first + g * logarithm - rest; + if (answer.Nodes.Any(node => node is Greaterf or GreaterOrEqualf or Lessf or LessOrEqualf && node.ContainsNode(x) + && node is not Greaterf(Powf(_, Number.Integer { EInteger.IsEven: true }), Number.Integer { IsZero: true }))) + return null; + return answer; + } + finally + { + againstAnInverseHyperbolicFunction = false; + } + + // The derivative, and the function's own quadratic, 1 - u^2 or the radicand, as numerators + // in x: the quadratic first, then for 1 - u^2 its factors 1 - u and 1 + u. + static (Entity Derivative, List Quadratic)? ADerivativeOfTheLogarithm(Entity argument, Entity.Variable x) + { + // (1 + u)/(1 - u) or (u + 1)/(u - 1): the two sides add to 2 or differ by 2. + if (argument is Divf(var top, var bottom)) + { + Entity u; + if (Functions.PartialFractions.IsZeroAsAValue((top + bottom - 2).InnerSimplified)) + u = ((top - bottom) / 2).InnerSimplified; + else if (Functions.PartialFractions.IsZeroAsAValue((top - bottom - 2).InnerSimplified)) + u = ((top + bottom) / 2).InnerSimplified; + else + return null; + // Of a linear only: beside atanh(x sqrt(e)/sqrt(d + e x^2)) the remainder is slow to decline. + if (!TreeAnalyzer.TryGetPolyLinear(u, x, out _, out _)) + return null; + return (2 * u.Differentiate(x) / (1 - MathS.Sqr(u)), + new List { Above(1 - MathS.Sqr(u)), Above(1 - u), Above(1 + u) }); + } + // A + sqrt(R) with R - A^2 a nonzero constant. + if (argument is not (Sumf or Minusf)) + return null; + Entity? radicand = null; + Entity a = Number.Integer.Zero; + foreach (var term in Sumf.LinearChildren(argument)) + { + if (radicand is null && term is Powf(var r, Number.Rational half) && half == Number.Rational.Create(1, 2) && r.ContainsNode(x)) + radicand = r; + else + a += term; + } + if (radicand is null || !a.ContainsNode(x)) + return null; + var difference = radicand - a * a; + var atAPoint = Functions.PartialFractions.Bare(difference.Substitute(x, Number.Rational.Create(37, 23)).Simplify()); + if (atAPoint.ContainsNode(x) || VanishesIdentically(atAPoint) || !VanishesIdentically(difference - atAPoint)) + return null; + return (a.Differentiate(x) / MathS.Sqrt(radicand), new List { Above(radicand) }); + + // `1/(c x)^2 - 1` is `1 - (c x)^2` over `(c x)^2`, and the factors are above. + static Entity Above(Entity expr) => Functions.SingleQuotient.Of(expr).Numerator; + } + } + + /// Set while asks its questions. + [System.ThreadStatic] private static bool againstAnInverseHyperbolicFunction; + /// Set while asks its question. [System.ThreadStatic] private static bool takingTheWholePowersOut; diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index d55feb5a3..9223eb21e 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -795,6 +795,8 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) if ((answer = IndefiniteIntegralSolver.SolveAPowerOfAQuotientOfLinearsOverAHalfOddPowerOfTheirProduct(expr, x)) is { }) return answer; // A quotient by a radical sum whose conjugate product is a constant, over the conjugate. if ((answer = IndefiniteIntegralSolver.SolveByTheConjugateOfARadicalSumBelowTheBar(expr, x, integrateByParts)) is { }) return answer; + // By parts against an inverse hyperbolic function beside a non-polynomial, its logarithm's derivative written algebraically. + if ((answer = IndefiniteIntegralSolver.SolveByPartsAgainstAnInverseHyperbolicFunctionAsItsLogarithm(expr, x, integrateByParts)) is { }) return answer; // A power of an exponential with a positive base is the exponential of the product, // exactly, and only that spelling is one the exponential rules read. if ((answer = IndefiniteIntegralSolver.SolveByFlatteningAPowerOfAnExponential(expr, x, integrateByParts)) is { }) return answer; diff --git a/Sources/Tests/UnitTests/Calculus/ByPartsAgainstAnInverseHyperbolicFunctionAsItsLogarithmIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ByPartsAgainstAnInverseHyperbolicFunctionAsItsLogarithmIntegralTest.cs new file mode 100644 index 000000000..0349f1fe0 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ByPartsAgainstAnInverseHyperbolicFunctionAsItsLogarithmIntegralTest.cs @@ -0,0 +1,60 @@ +// +// 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 +{ + /// + /// By parts against an inverse hyperbolic function beside something that is not a polynomial, + /// the logarithm the parser writes it as differentiated as the algebraic function it is: + /// ln(A + sqrt(A^2 + k))' = A'/sqrt(A^2 + k). These timed out. Rubi's 7.2, 7.5 and 7.6. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ByPartsAgainstAnInverseHyperbolicFunctionAsItsLogarithmIntegralTest + { + [Theory] + [InlineData("x*(a + b*acosh(c*x))/(d + k*x^2)^2", 1.2, 2.9)] + [InlineData("x*(a + b*acsch(c*x))/(d + k*x^2)^3", 0.2, 2.9)] + [InlineData("x^3*(a + b*asech(c*x))/(d + k*x^2)^3", 0.05, 0.7)] + [InlineData("acsch(a + b*x)/x^3", 0.2, 2.9)] + // Rubi's 7.6.1:130: the remainder's coefficients reach c^294, and read at the pinned symbols + // they were taken for zero, which dropped two of its three terms. + [InlineData("(d + k*x^2)^2*(a + b*acsch(c*x))/x^8", 0.2, 2.9)] + // Rubi's 7.3.4:380 and 7.4.1:82: a power of the function's own quadratic below the bar, which + // the substitution under the function answers. The rule declines it; asked first, it timed out. + [InlineData("atanh(a*x)/(1 - a^2*x^2)^3", 0.05, 0.7)] + [InlineData("acoth(x)/(1 - x^2)^3", 1.2, 2.9)] + public void ByParts(string integrand, double from, double to) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Entity Pinned(Entity e) => e.Substitute("a", 1.3).Substitute("b", 0.7).Substitute("c", 1.1).Substitute("d", 0.4).Substitute("k", 0.3); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + // Where the function is real. + for (var i = 0; i < 6; i++) + { + var at = from + (to - from) * (i + 0.5) / 6; + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + if (want.IsNaN || got.IsNaN) + continue; + compared++; + Assert.True(Math.Abs((double)(got - want).RealPart) + Math.Abs((double)(got - want).ImaginaryPart) + < 1e-9 * Math.Max(1, Math.Abs((double)want.RealPart) + Math.Abs((double)want.ImaginaryPart)), + $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}"); + } + Assert.True(compared >= 5, $"only {compared} points could be compared for {integrand}"); + } + } +}