diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 16c411cb1..2b46a0c8d 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -293,6 +293,23 @@ left is rational beside one root ([#718](https://github.com/asc-community/Angour | `"sec(e + f*x)*sqrt(a + a*sec(e + f*x))/sqrt(c + d*sec(e + f*x))".ToEntity().Integrate("x")` | `integral(...)` | 1,294 characters | | `"sec(e + f*x)^2/(sqrt(a + a*sec(e + f*x))*sqrt(c + d*sec(e + f*x)))".ToEntity().Integrate("x")` | `integral(...)` | 13,563 characters, a case for each sign of the root's leading coefficient | +### An odd rational function with symbols in it is integrated in the square of the variable + +**Answers where there were none, shorter ones, and sooner.** `(A + B x^2)/(x (a + b x^2 + c x^4)^3)` was declined, and so were +`x^3` in place of `x`, `1/(x^5 (d + g x^2)(a + b x^2 + c x^4))` and `1/(x (d + g x^2)(a + c x^4)^2)`, Rubi's 1.2.2.4, and +`(g x + f x^3)/(a + b x^2 + c x^4)^3`, 1.2.2.5, and `1/(x^5 (a + b x^2)^10)`, 1.1.2.2; on the unreleased master each ran past its budget, and several past a minute. A rational function +with symbols in it, past the third degree below the bar, whose one side holds only odd powers of x and the other only +even ones, is integrated in `u = x^2` before the partial fractions: `x N(x^2)/D(x^2) dx` is `N(u)/(2 D(u)) du`, and +`N(x^2)/(x D(x^2)) dx` is `N(u)/(2u D(u)) du`, where the written factors keep their shape. Each of these takes under a +second, in 650 to 4,200 characters. Where the unreleased master answered already, the answers are a third to a half +the size: `(A + B x^2)/(x^3 (a + b x^2))` was 656 characters and is 174 +([#718](https://github.com/asc-community/AngouriMath/issues/718)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"(A + B*x^2)/(x*(a + b*x^2 + c*x^4)^3)".ToEntity().Integrate("x")` | `integral(...)`; no answer within 60 s on the unreleased master | the antiderivative, in 0.7 s | +| `"(A + B*x^2)/(x^3*(a + b*x^2))".ToEntity().Integrate("x")` | `integral(...)` | `(-A/(a x^2) + (B a - A b)/a^2 ln(x^2) + (A b^2 - B a b)/(a^2 b) ln(b x^2 + a))/2` | + ### The secant beside powers of `a + a sec` and `c - c sec` is integrated in the first sum **Answers where there were none.** `sec(e + f x) (a + a sec(e + f x))^m sqrt(c - c sec(e + f x))` was declined, with diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 878771429..316a6759e 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -720,6 +720,75 @@ Entity Back(Entity mapped) return inTermsOfU.Substitute(u, quadratic); } + /// + /// An odd rational function N/D with symbols in it, past the third degree below the + /// bar, in u = x^2: one side holds only odd powers of x and the other only even + /// ones, so N/D dx is N1(u)/(2 D(u)) du for an odd numerator x N1(x^2), + /// and N(u)/(2u D1(u)) du for an odd denominator x D1(x^2). + /// + /// + /// In x the partial fractions split `(A + B x^2)/(x (a + b x^2 + c x^4)^3)` at x and carried + /// a polynomial of the eleventh degree over the cube of the quartic through the reduction, + /// past a minute; in u it is `(A + B u)/(2u (a + b u + c u^2)^3)`, answered in a tenth of a + /// second. Rubi's 1.2.2.4 and 1.2.2.5. + /// https://github.com/asc-community/AngouriMath/issues/718 + /// + private static Entity? InTheSquareAnOddRationalFunction(Entity numerator, Entity denominator, Entity.Variable x, bool integrateByParts) + { + if (!(numerator + denominator).Vars.Any(symbol => symbol != x) + || !TreeAnalyzer.TryGetPolynomial(numerator, x, out var above) || above.Count == 0 + || !TreeAnalyzer.TryGetPolynomial(denominator, x, out var below) || below.Count == 0 + || above.Values.Concat(below.Values).Any(coefficient => coefficient.ContainsNode(x)) + || above.Keys.Concat(below.Keys).Any(power => power.Sign < 0 || !power.CanFitInInt32()) + || below.Keys.Max()!.CompareTo(EInteger.FromInt32(4)) < 0) + return null; + static bool AllOfParity(Dictionary polynomial, bool odd) + => polynomial.Keys.All(power => power.IsEven != odd); + var oddAbove = AllOfParity(above, odd: true) && AllOfParity(below, odd: false); + var oddBelow = AllOfParity(above, odd: false) && AllOfParity(below, odd: true); + if (!oddAbove && !oddBelow) + return null; + var u = Variable.CreateUnique(numerator + denominator, "u_sq"); + // An even side as it is written, its powers of x halved; one with x otherwise, and an + // odd side, from its coefficients, a power of x taken off first. + Entity? AsWritten(Entity side) + { + var halved = side.Replace(node => node is Powf(var b, Number.Integer e) && b == x && e.EInteger.IsEven + ? MathS.Pow(u, Number.Integer.Create(e.EInteger.Divide(2))) : node); + return halved.ContainsNode(x) ? null : halved; + } + Entity FromCoefficients(Dictionary polynomial, int shift) + { + Entity sum = Number.Integer.Zero; + foreach (var pair in polynomial.OrderBy(pair => pair.Key)) + { + var k = (pair.Key.ToInt32Unchecked() - shift) / 2; + sum += k == 0 ? pair.Value : pair.Value * (k == 1 ? u : MathS.Pow(u, k)); + } + return sum; + } + // An odd side with an odd power of x written among its factors has one x taken off it + // there, and the rest is read as it is written: `x^5 (d + e x^2)(a + b x^2 + c x^4)` is + // `x` times `u^2 (d + e u)(a + b u + c u^2)`, where from its coefficients it is a sextic + // in u that nothing factors. + Entity? WithAnXTakenOff(Entity side) + { + var factors = Mulf.LinearChildren(side).ToList(); + var at = factors.FindIndex(factor => factor == x || factor is Powf(var b, Number.Integer e) && b == x && !e.EInteger.IsEven); + if (at < 0) + return null; + factors[at] = factors[at] == x ? Number.Integer.One + : MathS.Pow(x, Number.Integer.Create(((Number.Integer)((Powf)factors[at]).Exponent).EInteger.Subtract(EInteger.One))); + return AsWritten(factors.Aggregate((Entity)Number.Integer.One, (product, factor) => factor == Number.Integer.One ? product : product == Number.Integer.One ? factor : product * factor)); + } + var top = oddAbove ? WithAnXTakenOff(numerator) ?? FromCoefficients(above, 1) : AsWritten(numerator) ?? FromCoefficients(above, 0); + var bottom = oddBelow ? WithAnXTakenOff(denominator) ?? FromCoefficients(below, 1) : AsWritten(denominator) ?? FromCoefficients(below, 0); + var inU = oddAbove ? top / (2 * bottom) : top / (2 * u * bottom); + return Integration.ComputeAsAQuestionOfItsOwn(inU, u, integrateByParts) is { } answer + ? answer.Substitute(u, MathS.Pow(x, 2)) + : null; + } + /// /// A quotient of polynomials, split into two smaller quotients and integrated in two /// parts — at a rational root of the denominator where it has one, and otherwise at a @@ -771,6 +840,10 @@ Entity Back(Entity mapped) if (InTheQuadraticAnOddPowerOfItsDerivative(numerator, denominator, x, integrateByParts) is { } inTheQuadratic) return inTheQuadratic; + // An odd rational function with symbols in it, in u = x^2: see InTheSquareAnOddRationalFunction. + if (InTheSquareAnOddRationalFunction(numerator, denominator, x, integrateByParts) is { } inTheSquare) + return inTheSquare; + // A power of a sum in the numerator whose coefficients share a symbol, `(b d + 2c d x)^9`, // is that symbol's power times one of a sum without it: written so, it is a power of the // denominator's derivative, and the division of the improper fraction as written ran past diff --git a/Sources/Tests/UnitTests/Calculus/AnOddRationalFunctionInTheSquareIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnOddRationalFunctionInTheSquareIntegralTest.cs new file mode 100644 index 000000000..44c5a6854 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnOddRationalFunctionInTheSquareIntegralTest.cs @@ -0,0 +1,44 @@ +// +// 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 odd rational function with symbols in it is integrated in u = x^2 before the + /// partial fractions, which in x took a minute or more. Rubi's 1.2.2.4 and 1.2.2.5. + /// #718 + /// + [Trait("Area", "Calculus")] + public sealed class AnOddRationalFunctionInTheSquareIntegralTest + { + [Theory] + [InlineData("(A + B*x^2)/(x*(a + b*x^2 + c*x^4)^3)")] + [InlineData("1/(x^5*(d + g*x^2)*(a + b*x^2 + c*x^4))")] + [InlineData("(g*x + f*x^3)/(a + b*x^2 + c*x^4)^3")] + public void DifferentiatesBack(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Entity Pinned(Entity e) => e.Substitute("a", 2).Substitute("b", 3).Substitute("c", 5) + .Substitute("d", -1).Substitute("g", 7).Substitute("f", 1).Substitute("A", 4).Substitute("B", -2); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + foreach (var at in new[] { -1.7, -0.6, 0.45, 1.3 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + Assert.True(Math.Abs((double)(got - want).RealPart) + Math.Abs((double)(got - want).ImaginaryPart) + < 1e-9 * Math.Max(1, Math.Abs((double)want.RealPart)), + $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}"); + } + } + } +}