diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index de44812b0..ee194734d 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -55,6 +55,19 @@ numbers for every coefficient the solve is now done in exact rationals, to degre | `"1/((3-2*x+x^2)^(11/2)*(1+x+2*x^2)^5)".ToEntity().Integrate("x")` | `integral(...)` | 1,853 characters | | `"1/((3-2*x+x^2)^(21/2)*(1+x+2*x^2)^10)".ToEntity().Integrate("x")` | `integral(...)` | 3,475 characters | +### An exponential of something of the reciprocal of a linear is integrated beside another linear + +**Answers where there were none.** `f^(c/(a + b x)) x^2` was declined, and so were 24 more of Rubi's 2.3: +`F^(a + b/(c + d x))/(e + f x)^m`, `e^(e/(c + d x)^2) (a + b x)^m` and `F^(e + f (a + b x)/(c + d x))/(g + h x)^m`. +Under `u = 1/(c + d x)` another linear `p + q x` is `((d p - q c) u + q)/d` times `c + d x`, so the question in `u` +is an exponential beside a rational function, which the exponential integral and the moments answer. Only +powers of `c + d x` itself were read before ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"e^(k/(c+d*x))/(a+b*x)".ToEntity().Integrate("x")` | `integral(...)` | two `Ei`, of `k/(c + d x)` and of `k/(c + d x) - k b/(a d - b c)` | +| `"f^(c/(a+b*x)^2)*x^4".ToEntity().Integrate("x")` | `integral(...)` | 1,262 characters, the error function and powers of `1/(a + b x)` | + ### A polynomial over a symbolic quartic that is two quadratics in a hidden variable is split over them **Answers where there were none.** `1/(a + b x + c x^2 + b x^3 + a x^4)` was declined, and so were Rubi's 1.3.1 diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index bd16fd71b..b0b245da0 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -21983,7 +21983,8 @@ private static bool TryReadAPowerTimesARadicalQuadratic(Entity expr, Entity.Vari /// internal static Entity? SolveAGaussianInAReciprocal(Entity expr, Entity.Variable x) { - if (ReadAnExponentialInTheReciprocalOfALinear(expr, x) is not var (constant, @base, exponentInU, u, linear, slope, power) + // Beside another linear, the moments alone do not answer it: that is the next rule's. + if (ReadAnExponentialInTheReciprocalOfALinear(expr, x) is not (var constant, var @base, var exponentInU, var u, var linear, var slope, var power, null) || !TreeAnalyzer.TryGetPolyQuadratic(exponentInU, u, out var square, out _, out _) || TreeAnalyzer.IsZero(square)) return null; var powerInU = -power - EInteger.FromInt32(2); @@ -22013,7 +22014,7 @@ private static bool TryReadAPowerTimesARadicalQuadratic(Entity expr, Entity.Vari /// internal static Entity? SolveAnExponentialInTheReciprocalOfALinear(Entity expr, Entity.Variable x, bool integrateByParts) { - if (ReadAnExponentialInTheReciprocalOfALinear(expr, x) is not var (constant, @base, exponentInU, u, linear, slope, power) + if (ReadAnExponentialInTheReciprocalOfALinear(expr, x) is not var (constant, @base, exponentInU, u, linear, slope, power, beside) || !TreeAnalyzer.TryGetPolynomial(exponentInU, u, out var monomials) || monomials.Keys.Any(k => k.Sign < 0 || k.CompareTo(EInteger.FromInt32(8)) > 0) || !monomials.Keys.Any(k => k.Sign > 0)) @@ -22023,10 +22024,23 @@ private static bool TryReadAPowerTimesARadicalQuadratic(Entity expr, Entity.Vari return null; var exponential = MathS.Pow(@base, exponentInU); var inU = powerInU.IsZero ? exponential : MathS.Pow(u, Number.Integer.Create(powerInU)) * exponential; + // Another linear is a linear in u, and the rest one quotient, the way the rules read it. + var factor = -constant / slope; + if (beside is var (above, below, slopePower)) + { + // -constant d^(slopePower - 1), one power of d. + var raised = slopePower - EInteger.One; + factor = raised.IsZero ? -constant + : raised.Sign > 0 ? -constant * (raised.Equals(EInteger.One) ? slope : MathS.Pow(slope, Number.Integer.Create(raised))) + : -constant / (raised.Equals(EInteger.FromInt32(-1)) ? slope : MathS.Pow(slope, Number.Integer.Create(-raised))); + var top = powerInU.Sign > 0 ? MathS.Pow(u, Number.Integer.Create(powerInU)) * above : above; + var bottom = powerInU.Sign < 0 ? MathS.Pow(u, Number.Integer.Create(-powerInU)) * below : below; + inU = bottom == Number.Integer.One ? top * exponential : top * exponential / bottom; + } if (Integration.ComputeAsAQuestionOfItsOwn(inU, u, integrateByParts) is not { } answerInU || answerInU.Nodes.Any(node => node == MathS.NaN)) return null; - return (-constant / slope * answerInU.Substitute(u, 1 / linear)).InnerSimplified; + return (factor * answerInU.Substitute(u, 1 / linear)).InnerSimplified; } /// @@ -22037,12 +22051,23 @@ private static bool TryReadAPowerTimesARadicalQuadratic(Entity expr, Entity.Vari /// b/L^k is b u^k -- rather than by substituting x = (1/u - c)/d and /// asking the simplifier to see that 1/(1/u)^2 is u^2. /// - private static (Entity Constant, Entity Base, Entity ExponentInU, Variable U, Entity Linear, Entity Slope, EInteger Power)? + /// + /// A whole power of another linear M = p + q x is read too, as M/L written in + /// u times L: with L = c + d x, M/L = ((d p - q c) u + q)/d, a + /// linear in u, so M^k adds k to the power of L and + /// ((d p - q c) u + q)^k to Beside, and -k to its power of d, kept + /// apart from the constant so that the caller writes one power of d, not d/d. + /// Beside is when there is no such factor. A quotient N/L^k in the exponent with N linear is + /// (N/L) u^(k - 1) the same way. Rubi's 2.3 f^(c/(a + b x)) x^m, + /// F^(a + b/(c + d x))/(e + f x)^m and F^(e + f (a + b x)/(c + d x))/(g + h x)^m + /// were declined, for want of the other linear. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + private static (Entity Constant, Entity Base, Entity ExponentInU, Variable U, Entity Linear, Entity Slope, EInteger Power, (Entity Above, Entity Below, EInteger SlopePower)? Beside)? ReadAnExponentialInTheReciprocalOfALinear(Entity expr, Entity.Variable x) { Powf? exponential = null; - Entity? linear = null; - var power = EInteger.Zero; + var linears = new List<(Entity Linear, EInteger Power)>(); Entity constant = Number.Integer.One; foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) { @@ -22056,30 +22081,59 @@ private static (Entity Constant, Entity Base, Entity ExponentInU, Variable U, En exponential = (Powf)factor; continue; } - // The linear, or a whole power of it, above or below the bar. + // A linear, or a whole power of one, above or below the bar. var (candidate, k) = factor is Powf(var b, Number.Integer n) && n.EInteger.CanFitInInt32() ? (b, n.EInteger) : (factor, EInteger.One); - if (linear is not null && linear != candidate - || !TreeAnalyzer.TryGetPolyLinear(candidate, x, out var slopeOfIt, out _) || slopeOfIt is null || TreeAnalyzer.IsZero(slopeOfIt)) + if (!TreeAnalyzer.TryGetPolyLinear(candidate, x, out var slopeOfIt, out _) || slopeOfIt is null || TreeAnalyzer.IsZero(slopeOfIt)) return null; - linear = candidate; - power = underneath ? power - k : power + k; + linears.Add((candidate, underneath ? -k : k)); } if (exponential is null) return null; - linear ??= ALinearBelowABar(exponential.Exponent, x); - if (linear is null || !TreeAnalyzer.TryGetPolyLinear(linear, x, out var slope, out _) || slope is null || TreeAnalyzer.IsZero(slope)) + // L is the linear the exponent divides by. Wherever the exponent could be written in u + // before other linears were read, it is the linear beside it as well. + var linear = ALinearBelowABar(exponential.Exponent, x); + if (linear is null || !TreeAnalyzer.TryGetPolyLinear(linear, x, out var slope, out var intercept) || slope is null || TreeAnalyzer.IsZero(slope)) return null; var u = Variable.CreateUnique(expr, "u"); + // M/L in u, for a linear M: ((d p - q c) u + q)/d's numerator, the d left to the caller. + Entity? OverTheLinearInU(Entity other) + { + if (!TreeAnalyzer.TryGetPolyLinear(other, x, out var q, out var p) || q is null || TreeAnalyzer.IsZero(q)) + return null; + return ((slope * p - q * intercept).InnerSimplified * u + q).InnerSimplified; + } + var power = EInteger.Zero; + var slopePower = EInteger.Zero; + Entity above = Number.Integer.One, below = Number.Integer.One; + var beside = false; + foreach (var (other, k) in linears) + { + power += k; + if (other == linear) + continue; + if (k.Abs().CompareTo(EInteger.FromInt32(32)) > 0 || OverTheLinearInU(other) is not { } inU) + return null; + beside = true; + var raised = k.Abs().Equals(EInteger.One) ? inU : MathS.Pow(inU, Number.Integer.Create(k.Abs())); + if (k.Sign > 0) + above *= raised; + else + below *= raised; + slopePower -= k; + } var exponentInU = exponential.Exponent.Replace(node => node switch { - Divf(var numerator, Powf(var below, Number.Integer k)) when below == linear => numerator * MathS.Pow(u, k), - Divf(var numerator, var below) when below == linear => numerator * u, - Powf(var below, Number.Integer { EInteger.Sign: < 0 } k) when below == linear => MathS.Pow(u, -k), + Divf(var numerator, Powf(var under, Number.Integer k)) when under == linear && !numerator.ContainsNode(x) => numerator * MathS.Pow(u, k), + Divf(var numerator, var under) when under == linear && !numerator.ContainsNode(x) => numerator * u, + Divf(var numerator, Powf(var under, Number.Integer k)) when under == linear && OverTheLinearInU(numerator) is { } inU + => inU / slope * MathS.Pow(u, (k - 1).InnerSimplified), + Divf(var numerator, var under) when under == linear && OverTheLinearInU(numerator) is { } inU => inU / slope, + Powf(var under, Number.Integer { EInteger.Sign: < 0 } k) when under == linear => MathS.Pow(u, -k), _ => node }); if (exponentInU.ContainsNode(x)) return null; - return (constant, exponential.Base, exponentInU, u, linear, slope, power); + return (constant, exponential.Base, exponentInU, u, linear, slope, power, beside ? (above.InnerSimplified, below.InnerSimplified, slopePower) : null); } /// diff --git a/Sources/Tests/UnitTests/Calculus/AnExponentialInTheReciprocalOfALinearBesideAnotherLinearIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnExponentialInTheReciprocalOfALinearBesideAnotherLinearIntegralTest.cs new file mode 100644 index 000000000..05ba5c042 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnExponentialInTheReciprocalOfALinearBesideAnotherLinearIntegralTest.cs @@ -0,0 +1,61 @@ +// +// 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 something of 1/L, L = c + d x, beside a whole power of another + /// linear M. Under u = 1/L, M is ((d p - q c) u + q)/d times L, + /// so the question in u is an exponential beside a rational function, which the + /// exponential integral and the moments answer. They were declined, since only powers of + /// L itself were read. Rubi's 2.3. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnExponentialInTheReciprocalOfALinearBesideAnotherLinearIntegralTest + { + [Theory] + [InlineData("f^(c/(a + b*x))*x^2", 1.1)] + [InlineData("f^(c/(a + b*x))/x", 1.1)] + [InlineData("f^(c/(a + b*x)^2)*x", 1.1)] + [InlineData("f^(m + k/(c + d*x))/(g + h*x)^2", 1.1)] + [InlineData("e^(k/(c + d*x))*(a + b*x)^2", 1.1)] + [InlineData("e^(k/(c + d*x))/(a + b*x)", 1.1)] + [InlineData("e^(k/(c + d*x))/(a + b*x)", -0.9)] + [InlineData("e^(k/(c + d*x)^2)*(a + b*x)", -0.9)] + [InlineData("f^(m + k*(a + b*x)/(c + d*x))/(g + h*x)", 1.1)] + [InlineData("e^(1/x)/(x + 1)", 1.1)] + public void WrittenInTheReciprocal(string integrand, double k) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 5000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) => e.Substitute("a", 1.3).Substitute("b", 0.7).Substitute("c", 0.4).Substitute("d", 0.3) + .Substitute("f", 2.3).Substitute("g", 0.6).Substitute("h", 1.9).Substitute("m", 0.5).Substitute("k", k); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -0.9, -0.4, 0.3, 0.8, 1.3, 2.9 }) + { + 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}"); + } + } +}