From 05a7e6d0d20f7e53688b4d250a317e002264a20f Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Thu, 1 Oct 2026 00:15:24 +0000 Subject: [PATCH] A special function beside an elementary factor is integrated by parts, and so is its square x Si(b x) sin(b x) was left unevaluated where Si(b x) sin(b x) was not. By parts it is Si(b x) against x sin(b x), and the integral of that factor was taken with parts off, which x sin(b x) needs. Beside a special function of a linear argument, a polynomial times something else is integrated by the polynomial's parts now, which end in as many steps as its degree, and what is left, the special function's derivative times polynomials and elementary functions, is given parts in turn. A square of a special function of b x is two rounds of parts, and the remainder of the first comes back as one product over a sum, (b x Ei(b x) - e^(b x)) e^(b x)/(b x) for Ei(b x)^2, that no rule reads whole, where its terms e^(b x) Ei(b x) and e^(2 b x)/(b x) are each answered when asked. So where the remainder after a special function is not answered, its terms are asked one at a time: the ordinary way first, and then as questions of their own. That is only at the top, where such asking reaches. It is done once, since a term asked so is at the top itself, and the same asking inside it lifted the search again with the depth reset each time: together with a rule that wrote a + b x as t, x Ei(a + b x)^2 grew a test host to twelve gigabytes. A decline that the once-only flag caused is not cached (Integration.DeclinedForItsScope), as a decline that ran out of depth is not. It covers a multiple of x only: of a + b x the remainder divides by the linear, its terms are the exponentials a hyperbolic function is written as, and the decline of x Shi(a + b x)^2 went from a second to past twenty asking them. Rubi's family 8, #1501's tranche of 420: 296 -> 367 solved, 0 wrong, 0 timeout, against master with #1638. The 71 gains are 0 of 71 on master and 71 of 71 here run alone, and the 53 still declined take 103 s against 117 s together, median ratio 1.00, the slowest change being x^2 Ci(a + b x) cos(a + b x), 0.25 s to 2 s. Nothing else moves: the independent suites are 1756 of 1814 and families 1 to 7 at five a file 813 of 912 on both, 0 wrong everywhere, and what is declined in both costs what it did (11.9 s against 11.5 s, 98.1 s against 98.9 s). The unit tests pass, 14,247, and the performance gate passes on ee5cb79e, whose tree is this commit's but for this entry in BREAKING-CHANGES.md: allocation is what the baseline says on all 19 gated benchmarks. SpecialFunctionsByPartsTest: six rows beside a power and the elementary factor of the derivative, eleven squares, each differentiated back with its parameters pinned. Part of #1501. Co-Authored-By: Claude Opus 5.5 --- BREAKING-CHANGES.md | 20 ++++ .../Integration/IndefiniteIntegralSolver.cs | 99 ++++++++++++++++++- .../Integration/Integration.Definition.cs | 7 ++ .../Calculus/SpecialFunctionsByPartsTest.cs | 38 +++++++ 4 files changed, 163 insertions(+), 1 deletion(-) diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 08fb6e2bb..5772f5ada 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -2295,6 +2295,26 @@ no integrand without one changes: Rubi's independent suites and a sample of its | `"e^(b*x)*Ei(b*x)/x".Integrate("x")` | `Ei * e ^ (b * x) + C`, with `Ei` a variable | `Ei(b * x) ^ 2 / 2 + C` | | `"Si(b*x)*sin(b*x)/x".Integrate("x")` | `Si * -cos(b * x) + C`, with `Si` a variable | `Si(b * x) ^ 2 / 2 + C` | +### A square of a special function, and one beside a power of `x` and its elementary derivative, are integrated by parts + +`x Si(b x) sin(b x)` was left unevaluated where `Si(b x) sin(b x)` was not: by parts it is +`Si(b x)` against `x sin(b x)`, whose integral needs parts of its own, and the integral of that +factor was taken with parts off. It is integrated by the polynomial's parts now, beside a special +function of a linear argument, and what is left, `sin(b x)/x` times sines and cosines, by parts in +turn. A square of one of `b x` is two rounds of parts, and what the first round leaves comes back +as one product over a sum -- `(b x Ei(b x) - e^(b x)) e^(b x)/(b x)` for `Ei(b x)^2` -- that no rule +reads whole: its terms are asked one at a time now +([#1501](https://github.com/asc-community/AngouriMath/issues/1501)). An integrand holding one of +these functions had no reading in 2.5.0, which the entries for the functions themselves record, and +no integrand without one changes: Rubi's independent suites and a sample of its families 1 to 7, +2726 problems, are answered alone as they were. + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"x*Si(b*x)*sin(b*x)".Integrate("x")` | a polynomial in a variable `Si` | an antiderivative in `Si(b x)`, `Ci(2 b x)` and `ln(x)` | +| `"Ei(b*x)^2".Integrate("x")` | `Ei * (b * x) ^ 3 / 3 / b + C`, with `Ei` a variable | an antiderivative in `Ei(b x)` and `Ei(2 b x)` | +| `"x*erf(b*x)^2".Integrate("x")` | `UnrecognizedFunctionParseException`: there is no function `erf` | an antiderivative in `erf(b x)` | + ### An inverse trigonometric function below the bar is integrated to the sine and cosine integrals `1/arcsin(x)` was left unintegrated. Under the substitution that undoes the inverse function, diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index b230c5e49..ff88f4070 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -1722,6 +1722,19 @@ static Entity WithTheConstantMatchedTo(Entity antiderivative, Entity derivativeO // a coefficient's sign whose conditions are decided, `0 = 0`, and a remainder // holding that piecewise is read by nothing. var integralOfU = Integration.ComputeIndefiniteIntegral(u, x, false)?.InnerSimplified; + // Beside a special function of a linear argument, a polynomial times something + // else is integrated by the polynomial's parts, which end in as many steps as its + // degree: `x sin(b x)` beside `Si(b x)`, which nothing answers with parts off. + // https://github.com/asc-community/AngouriMath/issues/1501 + var besideASpecialFunction = IsASpecialFunctionOfALinearOrAPowerOfOne(v, x); + var byThePolynomialsParts = false; + if (integralOfU is null && besideASpecialFunction + && ThePolynomialFactor(u, x) is var (polynomialBeside, restBeside) + && polynomialBeside is not null && restBeside is not null) + { + integralOfU = IntegrateByPartsPolynomial(polynomialBeside, restBeside, x)?.InnerSimplified; + byThePolynomialsParts = integralOfU is not null; + } if (integralOfU is null) return null; // Differentiate v @@ -1780,7 +1793,12 @@ static Entity WithTheConstantMatchedTo(Entity antiderivative, Entity derivativeO // https://github.com/asc-community/AngouriMath/issues/1501 var partsOnTheRemainder = remaining.Nodes.Count() < wholeSize || (remainingPower >= 1 && remainingPower < wholePower) - || IsASpecialFunctionOfALinear(v, x) && MathS.TryPolynomial(u, x, out _); + || IsASpecialFunctionOfALinear(v, x) && MathS.TryPolynomial(u, x, out _) + // And where what was integrated beside it was a polynomial times an + // elementary function, whose integral the polynomial's parts wrote: the + // remainder is then the special function's derivative times polynomials and + // elementary functions, and parts on it descend on their degrees. + || byThePolynomialsParts; // Spelled as the product its own next step reads, where there is one. `Simplify` // writes `2 arctan(x)/(1 + x^2) * x^2/2` as `arctan(x) x^2/(x^2 + 1)`, a // quotient, and this rule runs on a product. `x*arctan(x)^2` is answered by @@ -1804,6 +1822,46 @@ static Entity WithTheConstantMatchedTo(Entity antiderivative, Entity derivativeO // search above it carry on, thirty seconds where the gate leaves it at one. var remainingIntegral = (Integration.AnsweringTheQuestionAsked ? SolveByEulerSubstitution(remaining, x) : null) ?? Integration.ComputeIndefiniteIntegral(remaining, x, partsOnTheRemainder); + // After a special function, or a power of one, the remainder is its derivative + // times what was integrated, written as one product over a sum: for `Ei(b x)^2`, + // `(b x Ei(b x) - e^(b x)) e^(b x)/(b x)`, which no rule reads, where its two + // terms are `e^(b x) Ei(b x)` and `e^(2 b x)/(b x)`, each answered when asked. + // So its terms are asked, each as a question of its own. From the top only, + // which is where that asking reaches, and once: a term asked so is at the top + // itself, and its own remainder asked the same way would lift the search again, + // each level with the depth reset. And of a multiple of x only: of `a + b x`, the + // remainder is over the linear, and its terms are the exponentials a hyperbolic + // function is written as, each a harder question than the whole -- the decline of + // `x Shi(a + b x)^2` went from a second to past twenty asking them. + // https://github.com/asc-community/AngouriMath/issues/1501 + if (remainingIntegral is null && besideASpecialFunction && OfAMultipleOfTheVariable(v, x) + && Integration.AnsweringTheQuestionAsked && !askingTheTermsOfARemainder) + { + var terms = remaining is Divf(var remainingAbove, var remainingBelow) + ? DistributedOverTheSum(remainingAbove).Select(term => term / remainingBelow).ToList() + : DistributedOverTheSum(remaining); + askingTheTermsOfARemainder = true; + try + { + foreach (var term in terms) + { + if ((Integration.ComputeIndefiniteIntegral(term, x, integrateByParts: true) + ?? Integration.ComputeAsAQuestionOfItsOwn(term, x, integrateByParts: true)) is not { } termIntegral) + { + remainingIntegral = null; + break; + } + remainingIntegral = remainingIntegral is null ? termIntegral : remainingIntegral + termIntegral; + } + } + finally + { + askingTheTermsOfARemainder = false; + } + } + else if (remainingIntegral is null && besideASpecialFunction && OfAMultipleOfTheVariable(v, x) + && Integration.AnsweringTheQuestionAsked) + Integration.DeclinedForItsScope(); if (remainingIntegral is null) return null; return v * integralOfU - remainingIntegral; @@ -20155,10 +20213,49 @@ Entity Term(Entity? sum, Entity coefficient, int power) /// #1501: the error /// functions and the exponential, logarithmic, sine, cosine and hyperbolic integrals. /// + /// + /// Set while the terms of a remainder after a special function are asked as questions of + /// their own, so that the asking is not nested; see TryIntegrateByPartsOnce. + /// + [System.ThreadStatic] private static bool askingTheTermsOfARemainder; + private static bool IsASpecialFunction(Entity node) => node is Entity.Erff or Entity.Erfcf or Entity.Erfif or Entity.Eif or Entity.Lif or Entity.Sif or Entity.Cif or Entity.Shif or Entity.Chif; + /// + /// A special function of a linear argument, as + /// reads one, or a positive whole power of one: Si(b x)^2. + /// + private static bool IsASpecialFunctionOfALinearOrAPowerOfOne(Entity factor, Variable x) + => IsASpecialFunctionOfALinear(factor, x) + || factor is Powf(var @base, Number.Integer power) && power.EInteger.Sign > 0 && IsASpecialFunctionOfALinear(@base, x); + + /// + /// Whether the special function in , or in the base of its power, + /// is of a multiple of , b x, with no offset. + /// + private static bool OfAMultipleOfTheVariable(Entity factor, Variable x) + => (factor is Powf(var @base, _) ? @base : factor).DirectChildren.FirstOrDefault() is { } argument + && TreeAnalyzer.TryGetPolyLinear(argument, x, out _, out var offset) && TreeAnalyzer.IsZero(offset); + + /// + /// The factors of that are polynomials in of + /// positive degree, against the rest: x sin(b x) is x and sin(b x). + /// Both halves where either would be empty. + /// + private static (Entity? Polynomial, Entity? Others) ThePolynomialFactor(Entity expr, Variable x) + { + Entity? polynomial = null; + Entity? rest = null; + foreach (var factor in Mulf.LinearChildren(expr)) + if (factor.ContainsNode(x) && MathS.TryPolynomial(factor, x, out _)) + polynomial = polynomial is null ? factor : polynomial * factor; + else + rest = rest is null ? factor : rest * factor; + return polynomial is null || rest is null || !rest.ContainsNode(x) ? (null, null) : (polynomial, rest); + } + /// /// Whether has what the derivative of the special function /// is made of: an exponential of a quadratic in diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index a2770fba7..885414428 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -412,6 +412,13 @@ private static Entity Normalized(Entity expr, Entity.Variable x) => /// [System.ThreadStatic] private static bool descentTruncated; + /// + /// Marks what is being worked out as declined for the scope it was asked in rather than + /// for its mathematics, as running out of descent is, so that the decline is not + /// remembered: the same question asked where that scope does not bind may be answered. + /// + internal static void DeclinedForItsScope() => descentTruncated = true; + /// /// The integrals this thread is part-way through, so that asking for one again while it is /// still being worked out is recognised as a cycle rather than followed round again. diff --git a/Sources/Tests/UnitTests/Calculus/SpecialFunctionsByPartsTest.cs b/Sources/Tests/UnitTests/Calculus/SpecialFunctionsByPartsTest.cs index b9ead215c..0b774a01e 100644 --- a/Sources/Tests/UnitTests/Calculus/SpecialFunctionsByPartsTest.cs +++ b/Sources/Tests/UnitTests/Calculus/SpecialFunctionsByPartsTest.cs @@ -95,5 +95,43 @@ public void TheLogarithmicIntegralIsByParts(string integrand) [InlineData("x^3*Chi(b*x)")] public void APowerTimesASpecialFunctionIsByParts(string integrand) => DifferentiatesBack(integrand, Parameters); + + /// + /// Beside a power of x times the elementary function its derivative is made of, the + /// special function is still the factor differentiated: x sin(b x) is integrated by + /// its polynomial's parts, and what is left, sin(b x)/x times that, is products of + /// sines and cosines over powers of x. + /// + [Theory] + [InlineData("x*Si(b*x)*sin(b*x)")] + [InlineData("x^3*Si(b*x)*sin(b*x)")] + [InlineData("x^2*Ci(b*x)*cos(b*x)")] + [InlineData("x*Ci(b*x)*sin(b*x)")] + [InlineData("x*Si(a + b*x)*sin(a + b*x)")] + [InlineData("x*Si(c + d*x)*sin(a + b*x)")] + public void BesideAPowerAndTheElementaryFactorOfItsDerivative(string integrand) + => DifferentiatesBack(integrand, Parameters); + + /// + /// A square of a special function of b x is two rounds of parts: the first against + /// the power of x leaves the special function once, beside its derivative, and that + /// is the case above or a substitution. The remainder of a round is asked term by term, + /// since it comes back as one product over a sum -- (b x Ei(b x) - e^(b x)) e^(b x)/(b x) + /// for Ei(b x)^2 -- that no rule reads whole. + /// + [Theory] + [InlineData("Ei(b*x)^2")] + [InlineData("x*Ei(b*x)^2")] + [InlineData("x^2*Ei(b*x)^2")] + [InlineData("x*Si(b*x)^2")] + [InlineData("Ci(b*x)^2")] + [InlineData("x*Ci(b*x)^2")] + [InlineData("x^2*erf(b*x)^2")] + [InlineData("x*erfc(b*x)^2")] + [InlineData("erfi(b*x)^2/x^3")] + [InlineData("x*Shi(b*x)^2")] + [InlineData("Chi(b*x)^2")] + public void ASquareIsTwoRoundsOfParts(string integrand) + => DifferentiatesBack(integrand, Parameters); } }