From b204d304208c6f3c0d115906e4b335867c23b19b Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Wed, 30 Sep 2026 13:11:06 +0000 Subject: [PATCH 1/2] Gruntz answers the example it is introduced with lim e^x (exp(1/x + e^(-x)) - exp(1/x)), x -> oo, is the example Gonnet and Gruntz introduce the algorithm with, and it came back as written. The algorithm did its steps right -- the mrv set {e^(-x), e^x}, the rewrite (e^(w + 1/x) - e^(1/x))/w -- and the series engine then read a leading term where there is none, in three places: - a constant term was taken off by subtracting it, which left 1/x - 1/x at w^0 in the argument of an exponential, and Exponentiate declined an argument it took to run off to infinity; - normalising a series divided the leading coefficient by itself, and the inner simplification leaves x/x as written, so the rest carried x/x - 1 at w^0 and the geometric series multiplied it out term by term; - e^(1/x) - e^(1/x), the two exponentials' constant terms, was not known to be zero: collecting like terms gives 0 provided not x = 0, a zero wherever the coefficient is defined. The constant term is now removed exactly, the normalised leading coefficient is 1 exactly, and the zero test collects like terms and reads a conditional zero as zero. Part of #353. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- .../Limits/Gruntz/AsymptoticSeries.cs | 44 ++++++++++++++++--- .../Tests/UnitTests/Calculus/GruntzTest.cs | 16 +++++++ 2 files changed, 54 insertions(+), 6 deletions(-) diff --git a/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/AsymptoticSeries.cs b/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/AsymptoticSeries.cs index 6528f3d6e..52832404b 100644 --- a/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/AsymptoticSeries.cs +++ b/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/AsymptoticSeries.cs @@ -97,7 +97,21 @@ private static bool IsKnownZero(Entity coefficient) var simplified = coefficient.InnerSimplified; if (simplified == Integer.Zero) return true; - return simplified.Evaled is Complex { IsZero: true }; + if (simplified.Evaled is Complex { IsZero: true }) + return true; + // A sum can cancel term for term where the inner simplification leaves it written + // out: e^(1/x) - e^(1/x) is what the constant terms of the two exponentials leave in + // e^x (e^(1/x + e^(-x)) - e^(1/x)), the example Gruntz's algorithm is introduced with. + // Collecting like terms settles that without a search. It answers 0 provided not + // x = 0, a zero wherever the coefficient is defined, and the condition excludes a + // point that says nothing about x running off to infinity -- the same reading as + // Gruntz.Bare's. + if (simplified is not (Sumf or Minusf)) + return false; + var collected = Simplificator.SimplifyChildren(simplified); + while (collected is Providedf(var value, _)) + collected = value; + return collected == Integer.Zero; } private static SortedDictionary Fresh() => new(); @@ -173,10 +187,28 @@ internal AsymptoticSeries Multiply(AsymptoticSeries other) var terms = Fresh(); foreach (var term in Terms) if (term.Key.CompareTo(e) >= 0) - terms[term.Key.Subtract(e)] = (term.Value / c).InnerSimplified; + // The leading coefficient over itself is 1 exactly, and is written so: the + // inner simplification leaves x/x as it is, and the rest of the series is + // this minus 1, where x/x - 1 would read as a term at w^0 that nothing can + // tell from zero. + terms[term.Key.Subtract(e)] = term.Key.CompareTo(e) == 0 ? Integer.One : (term.Value / c).InnerSimplified; return new AsymptoticSeries(terms, Order.Subtract(e)); } + /// + /// The series without its constant term: the series minus that term, taken off + /// exactly rather than subtracted, so that nothing is left at w^0 that the zero + /// test would have to recognise. + /// + internal AsymptoticSeries WithoutConstant() + { + var terms = Fresh(); + foreach (var term in Terms) + if (term.Key.CompareTo(ERational.Zero) != 0) + terms[term.Key] = term.Value; + return new AsymptoticSeries(terms, Order); + } + /// /// The reciprocal, by the geometric series on what is left once the leading term is /// taken out. Declines if the leading term cannot be found, since dividing by a @@ -187,7 +219,7 @@ internal AsymptoticSeries Multiply(AsymptoticSeries other) if (Normalised(out var coefficient, out var power) is not { } unit) return null; // 1/(1 + d) = 1 - d + d^2 - ..., which terminates because d starts above w^0. - var rest = unit.Add(Constant(-1)); + var rest = unit.WithoutConstant(); var sum = Constant(1); var term = Constant(1); for (var i = 0; i < MaxExpansionTerms; i++) @@ -212,7 +244,7 @@ internal AsymptoticSeries Multiply(AsymptoticSeries other) internal AsymptoticSeries? Exponentiate(ERational requested) { var constant = Terms.TryGetValue(ERational.Zero, out var atZero) ? atZero : Integer.Zero; - var rest = Add(Constant(-constant)); + var rest = WithoutConstant(); if (rest.LeadingTerm() is var (_, least) && least.CompareTo(ERational.Zero) <= 0) return null; // exp(constant) is left unexpanded on purpose: it is a coefficient, and opening @@ -241,7 +273,7 @@ internal AsymptoticSeries Multiply(AsymptoticSeries other) return null; var head = (MathS.Ln(coefficient) + Rational.Create(power) * logarithmOfW).InnerSimplified; // log(1 + d) = d - d^2/2 + d^3/3 - ... - var rest = unit.Add(Constant(-1)); + var rest = unit.WithoutConstant(); var sum = Constant(head); var term = Constant(1); for (var i = 1; i < MaxExpansionTerms; i++) @@ -324,7 +356,7 @@ internal AsymptoticSeries Multiply(AsymptoticSeries other) // Otherwise take the leading term out and use the binomial series on the rest. if (expanded.Normalised(out var coefficient, out var leading) is not { } unit) return null; - var rest = unit.Add(Constant(-1)); + var rest = unit.WithoutConstant(); var sum = Constant(1); var term = Constant(1); for (var i = 1; i < MaxExpansionTerms; i++) diff --git a/Sources/Tests/UnitTests/Calculus/GruntzTest.cs b/Sources/Tests/UnitTests/Calculus/GruntzTest.cs index f11b2dcf2..dc8d3e5ae 100644 --- a/Sources/Tests/UnitTests/Calculus/GruntzTest.cs +++ b/Sources/Tests/UnitTests/Calculus/GruntzTest.cs @@ -27,6 +27,22 @@ private static void AssertLimit(string expression, string destination, string ex expected.ToEntity().Evaled, expression.ToEntity().Limit("x", destination.ToEntity()).Evaled); + /// + /// The example the algorithm is introduced with -- Gonnet and Gruntz, Limit + /// Computation in Computer Algebra, ETH technical report 187, 1992, in Salvy's summary + /// of Gruntz's 1993 seminar -- where expanding in powers of 1/x gives + /// O(x^-k) for every k and settles nothing. In w = e^(-x) it is + /// e^(1/x) (e^w - 1)/w, whose leading coefficient e^(1/x) tends to 1. The + /// series has to carry 1/x as a parameter and see the two constant terms + /// e^(1/x) cancel. https://github.com/asc-community/AngouriMath/issues/353 + /// + [Theory] + [InlineData("e ^ x * (e ^ (1/x + e ^ (-x)) - e ^ (1/x))", "+oo", "1")] + [InlineData("e ^ x * (e ^ (2/x + e ^ (-x)) - e ^ (2/x))", "+oo", "1")] + [InlineData("e ^ (2 * x) * (e ^ (1/x + e ^ (-2 * x)) - e ^ (1/x))", "+oo", "1")] + public void TheExampleTheAlgorithmIsIntroducedWith(string expression, string destination, string expected) => + AssertLimit(expression, destination, expected); + /// /// The example the algorithm is usually shown with. Expanding the two exponentials /// separately gives two divergent series whose difference cancels entirely; rewriting From ca71bee36bf95587db1b55c8198ed10f49ccb9fe Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Wed, 30 Sep 2026 13:17:57 +0000 Subject: [PATCH 2/2] BREAKING-CHANGES: the Gruntz example is answered (measured on v2.5.0) Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- BREAKING-CHANGES.md | 19 +++++++++++++++++++ 1 file changed, 19 insertions(+) diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 0cb6e1ad6..e4bc750a0 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -2663,6 +2663,25 @@ arctangent, arccotangent or arccosecant whose argument diverges, which is now an | `integral(1/(2 + cos(x)), x, 0, pi)` | `integral(1 / (2 + cos(x)), x, 0, pi)` | `sqrt(3) * pi / 3` | | `limitleft(arctan(tan(x / 2)), x, pi)` | `limitleft(arctan(tan(x / 2)), x, pi)` | `pi / 2` | +### A limit at infinity whose terms cancel around a slowly varying factor is answered + +`lim e^x (exp(1/x + e^(-x)) - exp(1/x))` as `x → ∞` is the example Gonnet and Gruntz introduce +their algorithm with, and the implementation left it as written. The algorithm's steps were right, +and the series engine then read a leading term where there is none. It took a constant term off by +subtracting it, and divided a leading coefficient by itself, which left `1/x - 1/x` and `x/x - 1` at +`w^0`; it also did not know `e^(1/x) - e^(1/x)` for zero. The constant term is now removed exactly, +the normalised leading coefficient is `1` exactly, and the zero test collects like terms and reads a +conditional zero -- `0 provided not x = 0` -- as zero. +[#353](https://github.com/asc-community/AngouriMath/issues/353). Both columns measured on a build, +`v2.5.0` against this change. + +| `"….".ToEntity().Limit("x", "+oo")` of | Was (2.5.0) | Is | +|---|---|---| +| `e ^ x * (e ^ (1/x + e ^ (-x)) - e ^ (1/x))` | left as written | `1` | +| `e ^ x * (e ^ (2/x + e ^ (-x)) - e ^ (2/x))` | left as written | `1` | +| `e ^ (2 * x) * (e ^ (1/x + e ^ (-2 * x)) - e ^ (1/x))` | left as written | `1` | +| `e ^ (x + e ^ (-x)) - e ^ x` | `1` | the same | + ### With the downcasting off, integration and limits compute on exact numbers With `DowncastingEnabled` off, `1/(1 + c^2*x^2)` integrated to `NaN + C` and `limit(sin(c*x)/x, x, 0)`