From 11601004583af4215df9cab40f0b098c06ab80a4 Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Wed, 21 Oct 2026 15:33:19 +0000 Subject: [PATCH] Decline four slow routes ahead of the answering one by value, for the hyperbolic integral functions and the tangent MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Four routes spent gigabytes declining before the one that answers, on Rubi's 8.5 and 4.3.3.1: - SolveBySplittingSum expanded Shi(u) or Chi(u) times (e^u ± e^-u)/x^n into exponential terms, each non-elementary on its own, and declined after 10-13 GB; SolveIntegratingByParts then answered the whole in 0.2-0.9 GB. IsAHyperbolicIntegralBesideItsExponentialsOverAPolynomial sends such an integrand to it first, the same call at the same depth. - SolveByScalingTheVariable ran three Simplify calls before finding that x = c t leaves c inside e^(b x) or Shi(u). TheScaleStaysInsideAFunction evaluates the scaled integrand at two scales and two points first; it can only decline. - SolveBySubstitution simplified each quotient by a special function's derivative, a gigabyte each for x^2 Chi(a + b x)^2. A special function is transcendental over the elementary ones, so with it held a quotient taking two values at two points keeps x through any rewriting: VariesWhereASpecialFunctionIsHeld refuses it before the simplification. - SolveByTangentSubstitution's cot^k leaves 1^k above the bar and u^k below, outside the lowest-terms class that skips Simplify; it is in it now. Answers are identical. Chain against master f7f2a7a8: pocket of 3,919 rows (4.5.x, 4.3.3.x and 8.x every 2nd row, the rest of 4.3 and 4.4 every 6th, every 30th elsewhere), verdicts identical, allocation over the rows solved on both 460.9 -> 435.0 GB, 8.5 alone 31.0 -> 11.8 GB; samples unchanged; no row's allocation grew past 1.3x. By the probe against d3a6bce9: 8.5:62 12.59 -> 0.17 GB, 8.5:154 12.68 -> 0.19 GB, 4.3.3.1:313 5.01 -> 0.20 GB. Part of #233. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- .../Integration/IndefiniteIntegralSolver.cs | 184 +++++++++++++++++- ...ralsBesideTheirExponentialsIntegralTest.cs | 84 ++++++++ 2 files changed, 266 insertions(+), 2 deletions(-) create mode 100644 Sources/Tests/UnitTests/Calculus/HyperbolicIntegralsBesideTheirExponentialsIntegralTest.cs diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 31dc3156e..d5e191da4 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -72,6 +72,14 @@ internal static class IndefiniteIntegralSolver // parts on what it hands down was a two-minute search. if (integrateByParts && IsAnInverseFunctionTimesARationalFunction(expr, x) && SolveIntegratingByParts(expr, x) is { } byParts) return byParts; + // And a hyperbolic sine or cosine integral of a linear u beside a sum of e^u and e^(-u), + // which is how a hyperbolic function of u beside it is written, over a polynomial: parts + // against the whole pair the integral's derivative, sinh(u)/u, with the whole hyperbolic + // function, where expanded each exponential leaves `Ei(u) e^(-u)/u` or the like, which + // has no closed form, and `cosh(b x) Shi(b x)/x^2` spent thirteen gigabytes declining its + // two terms before parts answered it whole in a quarter of one. + if (integrateByParts && IsAHyperbolicIntegralBesideItsExponentialsOverAPolynomial(expr, x) && SolveIntegratingByParts(expr, x) is { } byPartsWhole) + return byPartsWhole; var splitted = TreeAnalyzer.GatherLinearChildrenOverSumAndExpand(expr, e => e.ContainsNode(x)); if (splitted is null || splitted.Count < 2) return null; // nothing to do, let other solvers do the work // Each expanded term with a factor written on both sides of its bar cancelled: @@ -428,6 +436,44 @@ private static Entity WithTheFactorsEveryTermSharesOut(Entity sum, Entity.Variab /// The largest term a sum is split over as written, before the expanding gather. private const int LargestTermTakenAsWritten = 60; + /// + /// Whether is Shi(u) or Chi(u) for a linear u + /// times sums of constant multiples of e^u and e^(-u), over a polynomial in + /// that is not a constant, with nothing else in . + /// + private static bool IsAHyperbolicIntegralBesideItsExponentialsOverAPolynomial(Entity expr, Entity.Variable x) + { + if (expr is Sumf or Minusf) + return false; + var (above, below) = Functions.SingleQuotient.Of(expr); + if (!below.ContainsNode(x) || !TreeAnalyzer.TryGetPolynomial(below, x, out _)) + return false; + var factors = Mulf.LinearChildren(above).Where(factor => factor.ContainsNode(x)).ToList(); + if (factors.Count(factor => factor is Shif or Chif) != 1) + return false; + var argument = factors.First(factor => factor is Shif or Chif).DirectChildren.First(); + if (!TreeAnalyzer.TryGetPolyLinear(argument, x, out var slope, out var offset)) + return false; + bool AnExponentialOfTheArgument(Entity node) + => node is Powf(var @base, var exponent) && @base == MathS.e + && TreeAnalyzer.TryGetPolyLinear(exponent, x, out var rate, out var shift) + && (VanishesIdentically(rate - slope) && VanishesIdentically(shift - offset) + || VanishesIdentically(rate + slope) && VanishesIdentically(shift + offset)); + bool ATermInTheExponentials(Entity term) + => Mulf.LinearChildren(term).Count(factor => factor.ContainsNode(x)) == 1 + && Mulf.LinearChildren(term).Where(factor => factor.ContainsNode(x)).All(AnExponentialOfTheArgument); + var sums = 0; + foreach (var factor in factors) + { + if (factor is Shif or Chif) + continue; + if (factor is not (Sumf or Minusf) || !Sumf.LinearChildren(factor).All(term => !term.ContainsNode(x) || ATermInTheExponentials(term))) + return false; + sums++; + } + return sums > 0; + } + /// /// Whether the integrand is one logarithm or inverse trigonometric function of /// , or a whole power of one, times a rational function of @@ -13735,7 +13781,11 @@ static bool IsTwice(Entity argument, Entity.Variable x) : Functions.SingleQuotient.Combine(inU / (1 + MathS.Sqr(uSub))) is var overOneBar && IsAQuotientOfDistinctSumsOfSymbolsInLowestTerms(overOneBar, uSub) ? WithNeighboursToOnePowerGathered(overOneBar) - : overOneBar.Simplify(); + // And so beside a power of u below the bar, what a power of the cotangent leaves + // with its ones above: see AQuotientOfDistinctSumsOverAPowerOfTheVariable. + : AQuotientOfDistinctSumsOverAPowerOfTheVariable(overOneBar, uSub) is { } overThePower + ? overThePower + : overOneBar.Simplify(); if (integrand is Providedf(var withoutCondition, _)) integrand = withoutCondition; // Asked as the same question: the rewrite renames the variable and adds no step of @@ -13827,6 +13877,57 @@ private static bool IsAQuotientOfDistinctSumsOfSymbolsInLowestTerms(Entity expr, return holdsASymbol && ReferenceEquals(WithNeighboursToOnePowerGathered(top), top); } + /// + /// , a quotient of sums as + /// reads them beside one positive whole power of below the bar, with the + /// ones above it -- 1 or a power of it -- left out, as the simplifier hands it back; + /// null where it is not such a quotient. + /// + /// + /// A power of the cotangent is (1/u)^k under the tangent, and over one bar it is + /// 1^k above and u^k below: Rubi's 4.3.3.1 + /// cot(c + d x)^7 (a + b tan(c + d x))^4 (A + B tan(c + d x)) is + /// 1^7 (a + b u)^4 (A + B u)/(u^7 (1 + u^2)), which the simplifier spent 4.8 of the + /// integral's 5 GB on to hand back without the 1^7. The sums have constant terms, so + /// the power of u divides none of them and the quotient is in lowest terms as before. Only + /// with an exponent no other power below the bar has, since the gathering of neighbours to + /// one power was measured for sums alone. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + private static Entity? AQuotientOfDistinctSumsOverAPowerOfTheVariable(Entity expr, Entity.Variable u) + { + if (expr is not Divf(var top, var bottom)) + return null; + static bool IsOne(Entity factor) => factor == Number.Integer.One || factor is Powf(var one, Number.Integer) && one == Number.Integer.One; + var above = Mulf.LinearChildren(top).Where(factor => !IsOne(factor)).ToList(); + var below = Mulf.LinearChildren(bottom); + if (above.Count == 0) + return null; + Entity? power = null; + var exponent = EInteger.One; + var others = new List(); + foreach (var factor in below) + { + if (factor == u || factor is Powf(var @base, Number.Integer { EInteger.Sign: > 0 }) && @base == u) + { + if (power is not null) + return null; + power = factor; + exponent = factor is Powf(_, Number.Integer { EInteger: var whole }) ? whole : EInteger.One; + continue; + } + others.Add(factor); + } + if (power is null || others.Count == 0 + || others.Any(factor => factor is Powf(_, Number.Integer { EInteger: var whole }) && whole.Equals(exponent))) + return null; + static Entity Product(IEnumerable factors) => factors.Aggregate((left, right) => left * right); + var withoutThePower = Product(above) / Product(others); + return IsAQuotientOfDistinctSumsOfSymbolsInLowestTerms(withoutThePower, u) + ? WithNeighboursToOnePowerGathered(Product(above) / bottom) + : null; + } + /// /// with each product of two powers to one whole exponent written /// as that power of their product, from the leaves up: (P^2 Q^2) R^2 is @@ -30615,6 +30716,45 @@ bool IsAConstantTimesAPowerOfOneMinusUSquared(Entity polynomial, out Entity cons return scaled; } + /// + /// Whether under x = c t, for c the + /// , is shown by its values not to be a function of the scale times + /// one of t: where a function that is not algebraic has an argument holding the + /// variable, the quotient of the scaled integrand's values at two points of t is taken + /// at two values of the scale, and the two quotients differ. False wherever a value cannot be + /// had, the values agree, or no such function is written. + /// + private static bool TheScaleStaysInsideAFunction(Entity expr, Entity.Variable x, Entity.Variable scale) + { + if (!expr.Nodes.Any(node => node switch + { + Sumf or Minusf or Mulf or Divf or Number or Variable => false, + Powf(_, var exponent) => exponent.ContainsNode(x), + _ => node.DirectChildren.Any(child => child.ContainsNode(x)), + }) + || PowerOfXAsWritten(expr, x) > LargestPowerEvaluatedForARepeat) + return false; + // In decimals: downcast, every value is folded exactly. + using var _ = MathS.Settings.DowncastingEnabled.Set(false); + Number.Complex? At(string scaleValue, string point) + { + var c = Number.Real.Create(EDecimal.FromString(scaleValue)); + var value = (expr.Substitute(x, c * Number.Real.Create(EDecimal.FromString(point))).Substitute(scale, c) * c).Evaled; + return value is Number.Complex { IsNaN: false } number && number.IsFinite && !number.IsZero ? number : null; + } + if (At("1", "0.83") is not { } oneFirst || At("1", "1.29") is not { } oneSecond + || At("1.37", "0.83") is not { } scaledFirst || At("1.37", "1.29") is not { } scaledSecond) + return false; + var first = ToComplex(scaledFirst) / ToComplex(oneFirst); + var second = ToComplex(scaledSecond) / ToComplex(oneSecond); + var magnitude = System.Math.Max(1.0, System.Numerics.Complex.Abs(first)); + var difference = System.Numerics.Complex.Abs(first - second); + return !double.IsNaN(difference) && !double.IsInfinity(difference) && !double.IsNaN(magnitude) && !double.IsInfinity(magnitude) + && difference > 1e-6 * magnitude; + static System.Numerics.Complex ToComplex(Number.Complex number) + => new(number.RealPart.EDecimal.ToDouble(), number.ImaginaryPart.EDecimal.ToDouble()); + } + /// /// An integrand carrying one symbolic parameter, scaled by it — x = c t — so that /// what is left to integrate has the variable alone in it. @@ -30666,6 +30806,14 @@ bool IsAConstantTimesAPowerOfOneMinusUSquared(Entity polynomial, out Entity cons } if (scale is null) return null; + // A function that is not algebraic, of an argument holding the variable, takes the scale + // inside it under x = c t unless the argument divides it out -- `e^(b x)` is `e^(b^2 t)` + // and `Shi(x)` is `Shi(b t)` -- and then no power of it comes out, which two values of the + // scale say before the three simplifications below: they spent a gigabyte and more on + // `(b x e^(b x) + e^(b x)) (b x Shi(b x) - cosh(b x))`, what parts leaves of `x Shi(b x)^2`, to + // find it declined. Only a decline is decided so; agreeing values leave it to the check below. + if (TheScaleStaysInsideAFunction(expr, x, scale)) + return null; // dx = c dt. var t = Variable.CreateUnique(expr, "u_scale"); @@ -37845,7 +37993,7 @@ private static Dictionary PinsFor(IEnumerable symbol /// or a numeric power, where a symbol is not pinned, and where a value is not real and /// finite. A root of a negative number is not real here, odd or even. /// - private static (double Value, double Scale)? AtRealPoint(Entity expr, Dictionary pins, Entity.Variable x, double at) + private static (double Value, double Scale)? AtRealPoint(Entity expr, Dictionary pins, Entity.Variable x, double at, bool exponentials = false) { var scale = 0.0; double Of(Entity node) @@ -37866,6 +38014,7 @@ double Of(Entity node) => System.Math.Pow(Of(@base), whole.EInteger.ToInt32Checked()), Powf(var @base, Number.Rational fraction) => Of(@base) is var inner && inner >= 0 ? System.Math.Pow(inner, fraction.ERational.ToDouble()) : double.NaN, + Powf(var @base, var exponent) when exponentials && @base == MathS.e => System.Math.Exp(Of(exponent)), _ => double.NaN, }; if (System.Math.Abs(value) > scale) @@ -37934,6 +38083,29 @@ bool Constant(System.Func of) return false; } + /// + /// for a quotient by the derivative + /// of a special function, with that function written as : in doubles, + /// where it is a rational function of x, the symbols and exponentials, of any size. + /// + /// + /// A special function is transcendental over the elementary functions -- it is not one of + /// them, and what is algebraic over them is -- so a quotient elementary in x and rational in + /// the function is a function of the function alone only where, with it held, no x is left: + /// two values at two points then say that no rewriting that keeps the value takes the x out. + /// Real values at two positive points only, as there. + /// + private static bool VariesWhereASpecialFunctionIsHeld(Entity expr, Variable uSub, Entity.Variable x) + { + if (HoldsTheImaginaryUnit(expr) || !expr.ContainsNode(uSub) + || !expr.Nodes.All(node => node is Sumf or Minusf or Mulf or Divf or Number or Variable || node is Powf(_, Number) || node is Powf(var @base, _) && @base == MathS.e)) + return false; + var pins = PinsFor(expr.Vars, x); + return AtRealPoint(expr, pins, x, 0.29, exponentials: true) is (var first, var scaleFirst) + && AtRealPoint(expr, pins, x, 0.53, exponentials: true) is (var second, var scaleSecond) + && System.Math.Abs(first - second) > 1e-6 * System.Math.Max(1.0, System.Math.Max(scaleFirst, scaleSecond)); + } + /// /// Whether , written in and a candidate's variable /// , takes two real values at two points of x, @@ -41671,6 +41843,14 @@ private static bool IsAFunctionOfTheTangentAlone(Entity expr, Entity.Variable u) // simplification each time it was offered. if (KeepsTheVariableThroughEveryRewriting(inTermsOfU, u, uSub, x, source, duDx)) continue; + // And under a special function, where with it held the quotient takes two values at two + // points of x, refused before it is simplified, for the reason in + // VariesWhereASpecialFunctionIsHeld: `x^2 Chi(a + b x)^2` offered `Chi(a + b x)` to three + // remainders of its parts, and simplifying each quotient to refuse it was a gigabyte. + if (IsASpecialFunction(u) && inTermsOfU.ContainsNode(x) + && !inTermsOfU.Nodes.Any(node => IsASpecialFunction(node) && node.ContainsNode(x)) + && VariesWhereASpecialFunctionIsHeld(inTermsOfU, uSub, x)) + continue; integrandInU = !inTermsOfU.ContainsNode(x) && TreeAnalyzer.TryGetPolyLinear(u, x, out var slope, out var offset) && !slope.ContainsNode(x) && !offset.ContainsNode(x) && IsAFunctionOfTheTangentAlone(inTermsOfU, uSub) && !HoldsTheImaginaryUnit(inTermsOfU) diff --git a/Sources/Tests/UnitTests/Calculus/HyperbolicIntegralsBesideTheirExponentialsIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/HyperbolicIntegralsBesideTheirExponentialsIntegralTest.cs new file mode 100644 index 000000000..31f29db06 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/HyperbolicIntegralsBesideTheirExponentialsIntegralTest.cs @@ -0,0 +1,84 @@ +// +// 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 +{ + /// + /// Hyperbolic sine and cosine integrals beside a hyperbolic function of their argument, by + /// parts beside a power of x, and a power of the cotangent beside sums in the tangent. Each was + /// answered, after declines of five to sixteen gigabytes: the split wrote the hyperbolic + /// function in its two exponentials, neither elementary beside the integral apart, and asked + /// the second after the first had settled it; the scaling of the variable simplified three + /// times to find that `x = b t` leaves `b` inside `e^(b x)`; the substitution search simplified + /// quotients by the derivative of `Chi(a + b x)` that keep an x whatever is done to them; and + /// the tangent's quotient was simplified to drop the `1^7` the cotangent leaves. Rubi's 8.5 and + /// 4.3.3.1. #233 + /// + [Trait("Area", "Calculus")] + public sealed class HyperbolicIntegralsBesideTheirExponentialsIntegralTest + { + [Theory] + [InlineData("Shi(b*x)*sinh(b*x)/x^3")] + [InlineData("cosh(b*x)*Shi(b*x)/x^2")] + [InlineData("Chi(b*x)*cosh(b*x)/x^3")] + [InlineData("Chi(b*x)*sinh(b*x)/x^2")] + [InlineData("x*Chi(b*x)^2")] + [InlineData("x^2*Shi(a + b*x)^2")] + [InlineData("x^2*Chi(a + b*x)^2")] + [InlineData("cot(x)^7*(a + b*tan(x))^4*(A + B*tan(x))")] + public void DifferentiatesBack(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Entity Pinned(Entity e) => e.Substitute("a", 0.4).Substitute("b", 0.7).Substitute("A", 1.3).Substitute("B", 0.6); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 1.7, 2.9, 3.4 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + if (want.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}"); + } + + /// + /// The same shapes, bounded by what they allocate on this thread rather than by a clock, + /// in a symbol no other test here uses, since an integral already answered on the thread + /// is read from the cache and allocates nothing. Each took between three and a half and + /// twenty-three gigabytes before. + /// + [Theory] + [InlineData("Shi(w*x)*sinh(w*x)/x^3")] + [InlineData("cosh(w*x)*Shi(w*x)/x^2")] + [InlineData("Chi(w*x)*cosh(w*x)/x^3")] + [InlineData("Chi(w*x)*sinh(w*x)/x^2")] + [InlineData("x*Shi(w*x)^2")] + [InlineData("x*Chi(w*x)^2")] + [InlineData("cot(x)^7*(v + w*tan(x))^4*(P + Q*tan(x))")] + public void IsAnsweredWithoutTheLongDeclines(string integrand) + { + var expression = integrand.ToEntity(); + var before = GC.GetAllocatedBytesForCurrentThread(); + var text = expression.Integrate("x").Stringize(); + var allocated = GC.GetAllocatedBytesForCurrentThread() - before; + Assert.DoesNotContain("integral(", text); + Assert.True(allocated < 2_000_000_000L, $"integrating {integrand} allocated {allocated} bytes"); + } + } +}