diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index ee194734d..a14f230c9 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -84,6 +84,19 @@ coefficients are the roots of that quadratic, and the polynomial over it is spli | `"1/(4*a*c+4*c^2*x^2+4*c*d*x^3+d^2*x^4)".ToEntity().Integrate("x")` | `integral(...)` | 3,920 characters | | `"x/(4*a*c+4*c^2*x^2+4*c*d*x^3+d^2*x^4)".ToEntity().Integrate("x")` | `integral(...)` | 2,315 characters | +### An exponential of a quadratic over a power of the linear through its vertex is a moment about it + +**Answers where there were none.** `f^(a + b x + c x^2)/(b + 2 c x)` was declined, and so were five more of Rubi's 2.3, +`f^(a + b x + c x^2)/(b + 2 c x)^n` and `f^(b x + c x^2)/(b + 2 c x)^n` for `n` up to 3. Under `u = b + 2 c x` the exponent is +`u^2/(4 c) + a - b^2/(4 c)`, with no linear part, and the integrand is a moment about 0 of the Gaussian: the +exponential integral of the square for an odd `n`, the error function for an even one +([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"f^(a+b*x+c*x^2)/(b+2*c*x)".ToEntity().Integrate("x")` | `integral(...)` | `f^(a - b^2/(4 c)) Ei(ln(f) (b + 2 c x)^2/(4 c))/(4 c)` | +| `"e^(x^2+2*x)/(x+2)".ToEntity().Integrate("x")` | `integral(...)` | `integral(...)`, not through the vertex | + ### `|` is divisibility, and was disjunction **Read this one first.** It is the only entry here that changes what existing input *means* rather diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index b0b245da0..148c4e3f3 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -6995,6 +6995,65 @@ private static bool TryReadAMonomial(Entity term, Entity.Variable x, out ERation return Functions.PartialFractions.Bare((constant * shifted / d * (sum + ei * MathS.Ei(rate * u)).Substitute(u, linear)).InnerSimplified); } + /// + /// An exponential of a quadratic over a whole power of the linear through its vertex: + /// F^(A x^2 + B x + C)/(p + q x)^n with 2 A p = q B, which is how + /// f^(a + b x + c x^2)/(b + 2 c x)^n arrives. Under u = p + q x the exponent is + /// A/q^2 u^2 + C - B^2/(4 A), with no linear part, and the integrand + /// F^(A/q^2 u^2 + C - B^2/(4 A))/(q u^n) is a moment about 0 of the Gaussian: the + /// exponential integral of the square for an odd n, the error function for an even + /// one. Rubi's 2.3. https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + /// Only through the vertex: anywhere else the shift leaves a linear part in the exponent, + /// and the moments about 0 do not answer it. The exponent in u is written from the + /// coefficients, not by substituting x = (u - p)/q and asking the simplifier to see + /// that the linear part cancels. + /// + internal static Entity? SolveAnExponentialOfAQuadraticOverTheLinearThroughItsVertex(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (expr is not Divf) + return null; + Powf? exponential = null; + Entity? linear = null; + var power = EInteger.Zero; + Entity constant = Number.Integer.One; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + if (!factor.ContainsNode(x)) + { + constant = underneath ? constant / factor : constant * factor; + continue; + } + if (!underneath && exponential is null && factor is Powf(var @base, var exponent) && !@base.ContainsNode(x) && exponent.ContainsNode(x)) + { + exponential = (Powf)factor; + continue; + } + // The linear, or a whole power of it, below the bar; nothing else. + if (!underneath) + return null; + var (candidate, k) = factor is Powf(var below, Number.Integer n) && n.EInteger.Sign > 0 && n.EInteger.CanFitInInt32() + ? (below, n.EInteger) : (factor, EInteger.One); + if (linear is not null && linear != candidate) + return null; + linear = candidate; + power += k; + } + if (exponential is null || linear is null || power.CompareTo(EInteger.FromInt32(32)) > 0 + || !TreeAnalyzer.TryGetPolyLinear(linear, x, out var q, out var p) || TreeAnalyzer.IsZero(q) + || !TreeAnalyzer.TryGetPolyQuadratic(exponential.Exponent, x, out var a, out var b, out var c) || TreeAnalyzer.IsZero(a) + || !VanishesIdentically(2 * a * p - q * b)) + return null; + var u = Variable.CreateUnique(expr, "u"); + var exponentInU = (a / MathS.Sqr(q)).InnerSimplified * MathS.Sqr(u) + (c - MathS.Sqr(b) / (4 * a)).InnerSimplified; + var question = MathS.Pow(exponential.Base, exponentInU) / (power.Equals(EInteger.One) ? u : MathS.Pow(u, Number.Integer.Create(power))); + if (Integration.ComputeAsAQuestionOfItsOwn(question, u, integrateByParts) is not { } answer + || answer.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + return (constant / q * answer.Substitute(u, linear)).InnerSimplified; + } + /// /// int u^m e^(k u) du, where is e^(k u), as its /// elementary part and the coefficient of Ei(k u) in it: diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index 4f34b7cdc..f6c14bc2a 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -827,6 +827,8 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) if ((answer = IndefiniteIntegralSolver.SolveAProductOfExponentialsOfQuadratics(expr, x)) is { }) return answer; // An exponential of a linear over a power of a linear, onto the exponential integral. if ((answer = IndefiniteIntegralSolver.SolveAnExponentialOfALinearOverAPowerOfALinear(expr, x)) is { }) return answer; + // And of a quadratic over a power of the linear through its vertex, onto the moments about 0. + if ((answer = IndefiniteIntegralSolver.SolveAnExponentialOfAQuadraticOverTheLinearThroughItsVertex(expr, x, integrateByParts)) is { }) return answer; // And a sum of exponentials, which is how sinh and cosh arrive, onto Shi and Chi. if ((answer = IndefiniteIntegralSolver.SolveAHyperbolicOfALinearOverAPowerOfALinear(expr, x)) is { }) return answer; // And either over several linears, split into partial fractions over them first, as the diff --git a/Sources/Tests/UnitTests/Calculus/AnExponentialOfAQuadraticOverTheLinearThroughItsVertexIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnExponentialOfAQuadraticOverTheLinearThroughItsVertexIntegralTest.cs new file mode 100644 index 000000000..1b9aa4099 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnExponentialOfAQuadraticOverTheLinearThroughItsVertexIntegralTest.cs @@ -0,0 +1,57 @@ +// +// 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 +{ + /// + /// An exponential of a quadratic over a power of the linear through its vertex, + /// f^(a + b x + c x^2)/(b + 2 c x)^n. Under u = b + 2 c x the exponent has no linear + /// part, and the integrand is a moment about 0 of the Gaussian, which the integrator answers. + /// They were declined. Rubi's 2.3. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnExponentialOfAQuadraticOverTheLinearThroughItsVertexIntegralTest + { + [Theory] + [InlineData("f^(a + b*x + c*x^2)/(b + 2*c*x)")] + [InlineData("f^(a + b*x + c*x^2)/(b + 2*c*x)^2")] + [InlineData("f^(a + b*x + c*x^2)/(b + 2*c*x)^3")] + [InlineData("f^(b*x + c*x^2)/(b + 2*c*x)^3")] + [InlineData("e^(x^2 + 2*x)/(x + 1)")] + public void AMomentAboutTheVertex(string integrand) + { + 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", 0.4).Substitute("f", 2.3); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -2.9, -0.4, 0.3, 0.8, 1.3, 2.1 }) + { + 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}"); + } + + [Fact] + public void NotThroughTheVertexIsDeclined() + => Assert.Contains("integral(", "e^(x^2 + 2*x)/(x + 2)".ToEntity().Integrate("x").Stringize()); + } +}