diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 2b46a0c8d..91cf2dca3 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -193,6 +193,27 @@ of a linear over a linear, which is answered in a second | `"1/(sqrt(a + i*a*tan(x))*(c + d*tan(x))^(3/2))".ToEntity().Integrate("x")` | `integral(...)`; past thirty seconds on the unreleased master | 2,217 characters | | `"1/((a + i*a*tan(x))^(3/2)*(c + d*tan(x))^(5/2))".ToEntity().Integrate("x")` | `integral(...)`; past thirty seconds on the unreleased master | 3,333 characters | +### Every linear beside a block is split off, a quadratic with rational roots among them + +**Answers where there were none.** `sec(c + d x)/(a + b sin(c + d x)^3)` was declined, and so were +`sec(c + d x)^3` over the same block, both over its square (Rubi's 4.1.7), `cot(f x + g)` and `cot(f x + g)^3` over +`a + b sec(f x + g)^3` (4.5.7), and `1/((1 - x^2)(a + b x^3))` and `1/(x (1 - x^2)(a + b x^3))`. The split that writes +`N/(L^k B)` as `P/L^k + R/B` at the root of a linear `L` took one linear only, and the rest over the block was +integrated a level down, where neither split answers. A second linear in the block is split off at the same level +now, where its answer differentiates back, and with no linear written a quadratic over the rationals with rational +roots, `1 - x^2`, is taken as its two ([#718](https://github.com/asc-community/AngouriMath/issues/718)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"sec(c + d*x)/(a + b*sin(c + d*x)^3)".ToEntity().Integrate("x")` | `integral(...)` | 1,788 characters | +| `"sec(c + d*x)^3/(a + b*sin(c + d*x)^3)".ToEntity().Integrate("x")` | `integral(...)` | 2,743 characters | +| `"sec(c + d*x)/(a + b*sin(c + d*x)^3)^2".ToEntity().Integrate("x")` | `integral(...)` | 2,851 characters | +| `"sec(c + d*x)^3/(a + b*sin(c + d*x)^3)^2".ToEntity().Integrate("x")` | `integral(...)` | 4,571 characters | +| `"cot(f*x + g)/(a + b*sec(f*x + g)^3)".ToEntity().Integrate("x")` | `integral(...)` | 1,838 characters | +| `"cot(f*x + g)^3/(a + b*sec(f*x + g)^3)".ToEntity().Integrate("x")` | `integral(...)` | 2,829 characters | +| `"1/((1 - x^2)*(a + b*x^3))".ToEntity().Integrate("x")` | `integral(...)` | 1,600 characters | +| `"1/(x*(1 - x^2)*(a + b*x^3))".ToEntity().Integrate("x")` | `integral(...)` | 1,637 characters | + ### A symbol shared by a numerator sum's coefficients is taken out before the partial fractions **Answers where there were none.** `(b d + 2c d x)^9/(a + b x + c x^2)^3` ran past a minute, where diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 316a6759e..bb53654c7 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -1381,11 +1381,28 @@ Entity Coefficient(Dictionary polynomial, int power) : series[j] * MathS.Pow(x, j - k + 1) / (j - k + 1); if (remainder == Number.Integer.Zero) return overX; - return Integration.ComputeIndefiniteIntegral(remainder / block, x, integrateByParts) is { } overTheBlock + return OverTheRestOfTheBlock(remainder, block, x, integrateByParts) is { } overTheBlock ? overX + overTheBlock : null; } + /// + /// The rest R/B of a split at a power of x or of a linear: split again at a linear + /// the block still holds, at the same level, where that answer differentiates back, and + /// otherwise integrated a level down. + /// + /// + /// The splits answer only the question asked or one substitution below it, so through the + /// integrator a second linear was out of reach under the sine: `sec(x)/(a + b sin(x)^3)` is + /// `1/((1 - u^2)(a + b u^3))`, and the second linear of `1 - u^2` was two levels down. + /// https://github.com/asc-community/AngouriMath/issues/718 + /// + private static Entity? OverTheRestOfTheBlock(Entity remainder, Entity block, Entity.Variable x, bool integrateByParts) + => IntegrateOverAPowerOfALinearBesideABlock(remainder, block, x, integrateByParts) is { } atAnotherRoot + && Functions.PartialFractions.DerivativeHoldsAtSampledPoints(atAnotherRoot, remainder / block, x) + ? atAnotherRoot + : Integration.ComputeIndefiniteIntegral(remainder / block, x, integrateByParts); + /// /// N/(L^k B) for a linear L = g + h x other than x, beside a block B /// with a symbol in it, past the second degree, that does not vanish at the root of L: @@ -1417,11 +1434,11 @@ Entity Coefficient(Dictionary polynomial, int power) var (@base, power) = factor is Powf(var raised, Number.Integer e) && e.EInteger.Sign > 0 && e.EInteger.CanFitInInt32() ? (raised, e.EInteger.ToInt32Unchecked()) : (factor, 1); + // The first linear is split off, and any other is left in the block, where the rest + // over the block meets it again: `1/((1 - x)(1 + x)(a + b x^3))` was declined. if (@base != x && @base.ContainsNode(x) && TreeAnalyzer.TryGetPolyLinear(@base, x, out var slope, out var offset) - && !slope.ContainsNode(x) && !offset.ContainsNode(x)) + && !slope.ContainsNode(x) && !offset.ContainsNode(x) && (linear is null || linear == @base)) { - if (linear is { } && linear != @base) - return null; linear = @base; k += power; } @@ -1431,7 +1448,13 @@ Entity Coefficient(Dictionary polynomial, int power) pastTheSecondDegree |= TreeAnalyzer.TryGetPolynomial(@base, x, out var terms) && terms.Keys.Any(degree => degree.CompareTo(EInteger.FromInt32(2)) > 0); } } - if (linear is null || k > MaximumPowerOfXBesideABlock || !pastTheSecondDegree || !block.ContainsNode(x) || !block.Vars.Any(symbol => symbol != x) + // With no linear written, a quadratic over the rationals with rational roots is two: + // `1/((1 - x^2)(a + b x^3))` is what the sine makes of `sec(x)/(a + b sin(x)^3)`. + if (linear is null) + return WithQuadraticsOverTheRationalsSplit(denominator, x) is { } split + ? IntegrateOverAPowerOfALinearBesideABlock(numerator, split, x, integrateByParts) + : null; + if (k > MaximumPowerOfXBesideABlock || !pastTheSecondDegree || !block.ContainsNode(x) || !block.Vars.Any(symbol => symbol != x) || !TreeAnalyzer.TryGetPolyLinear(linear, x, out var h, out var g) || VanishesIdentically(h) || !TreeAnalyzer.TryGetPolynomial(block, x, out var below) || !TreeAnalyzer.TryGetPolynomial(numerator, x, out var above)) return null; @@ -1522,7 +1545,7 @@ Entity[] InS(Dictionary polynomial, int degree) var scale = k == 1 ? h : MathS.Pow(h, k); if (remainder == Number.Integer.Zero) return overTheLinear / scale; - return Integration.ComputeIndefiniteIntegral(remainder / block, x, integrateByParts) is { } overTheBlock + return OverTheRestOfTheBlock(remainder, block, x, integrateByParts) is { } overTheBlock ? (overTheLinear + overTheBlock) / scale : null; @@ -22031,6 +22054,37 @@ private static bool IsAProductOfSymbolicLinearFactors(Entity denominator, Entity return (linears >= 1 || !aLinearAmongThem) && linears + quadratics >= 2 && symbolic; } + /// + /// with each written factor that is a quadratic over the + /// rationals with rational roots, or a power of one, written as its leading coefficient + /// times its two linears; where there is none. + /// + private static Entity? WithQuadraticsOverTheRationalsSplit(Entity denominator, Entity.Variable x) + { + var changed = false; + Entity product = Number.Integer.One; + foreach (var factor in Mulf.LinearChildren(denominator)) + { + var (@base, power) = factor is Powf(var raised, Number.Integer { EInteger.Sign: > 0 } p) ? (raised, p) : (factor, Number.Integer.One); + if (!TreeAnalyzer.TryGetPolynomial(@base, x, out var read) || read.Count == 0 + || !read.Keys.Max()!.Equals(EInteger.FromInt32(2)) || !HasRationalRoots(read)) + { + product *= factor; + continue; + } + Entity Coefficient(int degree) => read.TryGetValue(EInteger.FromInt32(degree), out var c) ? c : Number.Integer.Zero; + var (a, b, c) = (Coefficient(2), Coefficient(1), Coefficient(0)); + var root = MathS.Sqrt(b * b - 4 * a * c).Evaled; + var first = (x - ((-b - root) / (2 * a)).Evaled).InnerSimplified; + var second = (x - ((-b + root) / (2 * a)).Evaled).InnerSimplified; + product *= power == Number.Integer.One + ? a.Evaled * first * second + : MathS.Pow(a.Evaled, power) * MathS.Pow(first, power) * MathS.Pow(second, power); + changed = true; + } + return changed ? product : null; + } + /// /// Whether a quadratic read as its coefficients has rational coefficients and a /// discriminant that is the square of a rational. diff --git a/Sources/Tests/UnitTests/Calculus/EveryLinearBesideABlockIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/EveryLinearBesideABlockIntegralTest.cs new file mode 100644 index 000000000..5ee158f60 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/EveryLinearBesideABlockIntegralTest.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 +{ + /// + /// Every linear beside a block with symbols in it is split off, at the same level, and a + /// quadratic over the rationals with rational roots is taken as its two linears: under the + /// sine, sec(x)/(a + b sin(x)^3) is 1/((1 - u^2)(a + b u^3)). Rubi's 4.1.7. + /// #718 + /// + [Trait("Area", "Calculus")] + public sealed class EveryLinearBesideABlockIntegralTest + { + [Theory] + [InlineData("1/((1 - x)*(1 + x)*(a + b*x^3))")] + [InlineData("1/(x*(1 - x^2)*(a + b*x^3))")] + [InlineData("sec(x)/(a + b*sin(x)^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); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + foreach (var at in new[] { 0.15, 0.35, 0.6, 0.85 }) + { + 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}"); + } + } + } +}