diff --git a/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/Gruntz.cs b/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/Gruntz.cs index 3944a9538..c5dfb0f43 100644 --- a/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/Gruntz.cs +++ b/Sources/AngouriMath/Functions/Continuous/Limits/Gruntz/Gruntz.cs @@ -78,6 +78,7 @@ private static void Trace(string what, object? value) { if (depth > 0) return null; // already inside; the caller is the entry + expr = AsExponentials(expr, x); try { return LimitInf(expr, x); } catch (Core.Exceptions.AngouriBugException) { throw; } catch (OperationCanceledException) { throw; } @@ -164,6 +165,31 @@ private static bool IsExponential(Entity e, out Entity exponent) private static Entity Exponential(Entity exponent) => MathS.Pow(MathS.e, exponent); + /// + /// Every power whose exponent moves rewritten as an exponential, which is how + /// reads one in any case. + /// + /// + /// The mrv set holds subexpressions of the expression and + /// substitutes them by name, so a member has to occur in the expression as it stands. + /// Reading b^p as exp(p * ln(b)) inside Mrv alone put a member in the set that was + /// nowhere in the expression: for x^x / e^(x * ln(x)) the set came back holding the + /// same exponential twice, once as the constructed e^(x * ln(x)) and once as the + /// denominator's own e^(ln(x) * x), and the substitution found only the second. The + /// numerator went into the series as x^x, whose leading exponent reads as +1, and the + /// limit came back 0 where the two sides are equal and the answer is 1. + /// https://github.com/asc-community/AngouriMath/issues/735 + /// This assumes b is positive, which is what Mrv's own reading of the same node + /// already assumed; the algorithm is scoped to the exp-log functions, and a moving + /// exponent over a base that changes sign is outside that class either way. + /// + private static Entity AsExponentials(Entity e, Variable x) => e.Replace(node => + node is Powf(var @base, var power) + && @base != MathS.e + && power.ContainsNode(x) + ? Exponential((power * MathS.Ln(@base)).InnerSimplified) + : node); + /// /// The expression without the domain conditions simplification leaves behind. A limit /// is taken of an expression read as continuous, so a condition that excludes a point @@ -379,6 +405,12 @@ int Contained(Entity member) } if (rewritten.ContainsNode(x) && !logarithmOfW.ContainsNode(x)) return null; + // A member the substitution did not find is a member left in the series, where it + // is read as part of a coefficient and the conclusion is drawn from a leading term + // that is not the leading term. There is nothing to salvage from that, and saying + // nothing is the only safe reading -- this is where x^x / e^(x * ln(x)) answered 0. + if (members.Any(member => rewritten.ContainsNode(member.Member))) + return null; return (rewritten, logarithmOfW); } } diff --git a/Sources/Tests/UnitTests/Calculus/GruntzMovingExponentTest.cs b/Sources/Tests/UnitTests/Calculus/GruntzMovingExponentTest.cs new file mode 100644 index 000000000..5a7694426 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/GruntzMovingExponentTest.cs @@ -0,0 +1,82 @@ +// +// Copyright (c) 2019-2022 Angouri. +// AngouriMath is licensed under MIT. +// Details: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md. +// Website: https://am.angouri.org. +// + +using AngouriMath; +using AngouriMath.Extensions; +using Xunit; + +namespace AngouriMath.Tests.Calculus +{ + /// + /// Gruntz's algorithm puts the subexpressions of the fastest comparability class in a set + /// and then rewrites the expression by substituting each of them by name, so a member of + /// that set has to occur in the expression as it stands. A power whose exponent moves was + /// read as exp(p * ln(b)) when the set was built but left as written in the expression, so + /// the set came back holding one exponential twice -- once as the constructed + /// e^(x * ln(x)) and once as the denominator's own e^(ln(x) * x), the same product with + /// its factors the other way round -- and the substitution found only the second of them. + /// The numerator then went into the series as x^x, whose leading exponent reads as +1, + /// which the algorithm concludes from as "tends to zero". + /// + /// No issue exists for this; it was found while measuring what a factorial's Stirling + /// expansion would need, where x^x is the term the factorial is compared against. + /// + public sealed class GruntzMovingExponentTest + { + private static void AssertLimit(string expression, string expected) => + Assert.Equal( + expected.ToEntity().Evaled, + expression.ToEntity().Limit("x", "+oo".ToEntity()).Evaled); + + /// + /// x^x and e^(x * ln(x)) are one function, so each of these is an expression whose + /// value is known exactly rather than only asymptotically -- which is what makes the + /// expected answers here beyond argument. Every one of them answered 0 before. + /// + [Theory] + [InlineData("x ^ x / e ^ (x * ln(x))", "1")] + [InlineData("e ^ (x * ln(x)) / x ^ x", "1")] + [InlineData("x ^ (2 * x) / e ^ (2 * x * ln(x))", "1")] + [InlineData("(x ^ 2) ^ x / e ^ (2 * x * ln(x))", "1")] + public void APowerAndItsExponentialAreOneFunction(string expression, string expected) => + AssertLimit(expression, expected); + + /// + /// The same cancellation with something left over, so that the answer is not 1 and a + /// rule which merely stopped saying 0 would not pass. The quotients are e^x, x and + /// e^(-x) exactly. + /// + [Theory] + [InlineData("x ^ x / e ^ (x * ln(x) - x)", "+oo")] + [InlineData("x ^ x / e ^ (x * ln(x) - ln(x))", "+oo")] + [InlineData("x ^ x / e ^ (x * ln(x) + x)", "0")] + public void WhatIsLeftOverDecidesIt(string expression, string expected) => + AssertLimit(expression, expected); + + /// + /// The claim the expected values above rest on, checked at a point rather than argued: + /// the ratio is not merely close to 1, it is 1. + /// + [Theory] + [InlineData("(50 ^ 50) / e ^ (50 * ln(50))")] + [InlineData("(20 ^ 20) / e ^ (20 * ln(20))")] + public void ThePowerAndTheExponentialAgreeAtAPoint(string written) => + Assert.Equal(Entity.Number.Integer.Create(1), written.ToEntity().EvalNumerical()); + + // Growths that were already right and have to stay so, including the ones where a + // moving exponent competes with a fixed one. + [Theory] + [InlineData("x ^ x / e ^ x", "+oo")] + [InlineData("x ^ x / 2 ^ x", "+oo")] + [InlineData("2 ^ x / x ^ 2", "+oo")] + [InlineData("x ^ 2 / e ^ x", "0")] + [InlineData("x ^ x / x ^ (2 * x)", "0")] + [InlineData("x ^ (2 * x) / x ^ x", "+oo")] + public void TheOrdinaryOnesAreUnchanged(string expression, string expected) => + AssertLimit(expression, expected); + } +}