diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 1d9a1ab9e..16c411cb1 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -247,6 +247,25 @@ after it |---|---|---| | `"(a*c + b*c*x)^(-3-2*p)*(f + g*x)*(a^2 + 2*a*b*x + b^2*x^2)^p".ToEntity().Integrate("x")` | `integral(...)`; an answer with no value on the unreleased master | the antiderivative | +### A symbolic cubic's linear factors are found before the partial fractions + +**Answers where there were none.** `1/(a c h + (b c h + a d h + a c f) x + (b d h + b c f + a d f) x^2 + b d f x^3)`, +`(a + b x)(c + d x)(h + f x)` multiplied out, was declined, with its square and cube, Rubi's 1.3.1. A polynomial of +the third degree or more with symbols in it, whose leading and constant coefficients are single terms, is tried +for linear factors `u ± v x`, with `u` a divisor of the constant term and `v` one of the leading coefficient, as for +the rational roots of a polynomial over the integers. The factors that divide it exactly are written apart, and +those the numerator also holds are cancelled first. The same split answers three rows of Rubi's 4.7.2 and 4.7.7, +whose denominator in the tangent is a product of linears +([#718](https://github.com/asc-community/AngouriMath/issues/718)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"1/(a*c*h + (b*c*h + a*d*h + a*c*f)*x + (b*d*h + b*c*f + a*d*f)*x^2 + b*d*f*x^3)".ToEntity().Integrate("x")` | `integral(...)` | the antiderivative, in logarithms of the three linears, 288 characters | +| the same cubic squared, and cubed | `integral(...)` | answered, in 3,047 and 10,662 characters | +| `"sec(c + d*x)^6/(a*cos(c + d*x) + b*sin(c + d*x))".ToEntity().Integrate("x")` | `integral(...)` | answered, in 5,573 characters | +| `"sec(c + d*x)^3/(a*cos(c + d*x) + b*sin(c + d*x))^4".ToEntity().Integrate("x")` | `integral(...)` | answered, in 9,322 characters | +| `"1/(a*sec(x) + b*tan(x))^4".ToEntity().Integrate("x")` | `integral(...)` | answered, in 6,548 characters | + ### An odd power of a quadratic's derivative over a power of the quadratic is integrated in the quadratic **Shorter answers, and answers where there were none.** `(b + 2c x)^9/(a + b x + c x^2)^3` was answered through the diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 92edc7dc5..878771429 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -809,6 +809,18 @@ Entity Back(Entity mapped) if (SolveByReducingThePolynomialOverALinearFactor(expr, x, integrateByParts, alone: true) is { } overThePowersOfTheLinear) return overThePowersOfTheLinear; + // A symbolic polynomial of the third degree or more whose leading and constant coefficients + // are monomials, with the linear factors that divide it written apart: Rubi's 1.3.1 has + // `(a + b x)(c + d x)(e + f x)` multiplied out below the bar, which nothing read. + // The linears the numerator holds as well are divided out of both: the half-angle form of + // 4.7.7's `1/(a sec(x) + b tan(x))^4` has `1 - t` and `1 + t` four times each below the bar + // and `2 (1 - t^2)^8` above it, and uncancelled the partial fractions ran past two minutes + // where it was declined in three seconds. Through the partial fractions only: handed to the + // whole integrator where they decline, it went to the substitution search. + if (WithSymbolicLinearFactorsFound(numerator, denominator, x) is var (aboveTheFound, belowTheFound) + && SolveByPartialFractions(aboveTheFound / belowTheFound, x, integrateByParts) is { } overTheFoundLinears) + return overTheFoundLinears; + // A symbolic quadratic whose discriminant is a square in the symbols, written as its two // linears before the division: the numerator may share one of them, and over the // quadratic as written the division of `(d + e x)^8` ran past a minute. @@ -15177,6 +15189,115 @@ private static bool SharesAFactor(Entity above, Entity below, Entity.Variable x) return false; } + /// + /// with each written factor of the third degree or more in + /// , with symbols in it, whose leading and constant coefficients are + /// single terms, written as the linear factors u ± v x that divide it exactly and what + /// is left: u a divisor of the constant term and v one of the leading + /// coefficient, as for the rational roots of a polynomial over the integers. + /// where no factor has one. + /// + private static (Entity Numerator, Entity Denominator)? WithSymbolicLinearFactorsFound(Entity numerator, Entity denominator, Entity.Variable x) + { + var changed = false; + var cancelled = false; + Entity product = Number.Integer.One; + foreach (var factor in Mulf.LinearChildren(denominator)) + { + var (@base, power) = factor is Powf(var b, Number.Integer { EInteger.Sign: > 0 } p) ? (b, (Entity)p) : (factor, (Entity)Number.Integer.One); + if (@base is not (Sumf or Minusf) || !@base.ContainsNode(x) || !@base.Vars.Any(v => v != x)) + { + product *= factor; + continue; + } + var variables = @base.Vars.Concat(numerator.Vars).Distinct().OrderBy(v => v.Name, System.StringComparer.Ordinal).ToList(); + if (variables.Count > MultivariatePolynomial.MaxVariables) + { + product *= factor; + continue; + } + var indices = new Dictionary(); + for (var i = 0; i < variables.Count; i++) + indices[variables[i]] = i; + var at = indices[x]; + if (MultivariatePolynomial.TryParse(@base, indices) is not { } polynomial || polynomial.DegreeIn(at) < 3 + || !polynomial.CoefficientsIn(at).TryGetValue(polynomial.DegreeIn(at), out var leading) || leading.TermCount != 1 + || !polynomial.CoefficientsIn(at).TryGetValue(0, out var constant) || constant.TermCount != 1) + { + product *= factor; + continue; + } + var remaining = polynomial; + var found = new List(); + var xMonomial = MultivariatePolynomial.PackMonomial(at, 1); + foreach (var u in MonomialDivisors(constant, variables.Count)) + foreach (var v in MonomialDivisors(leading, variables.Count)) + foreach (var sign in new[] { ERational.One, ERational.One.Negate() }) + { + if (!v.Monomials.Any()) + continue; + var linear = u.Add(MultivariatePolynomial.Term(variables.Count, v.Monomials.First() | xMonomial, v.CoefficientOf(v.Monomials.First()).Multiply(sign))); + while (remaining.DegreeIn(at) >= 2 && remaining.DivideExact(linear) is { } quotient) + { + found.Add(linear); + remaining = quotient; + } + } + if (found.Count == 0) + { + product *= factor; + continue; + } + // Each found linear the numerator holds is divided out of it, once per power it stands to. + var multiplicity = ((Number.Integer)power).EInteger.ToInt32Unchecked(); + var above = MultivariatePolynomial.TryParse(numerator, indices); + foreach (var linear in found) + { + var left = multiplicity; + while (left > 0 && above?.DivideExact(linear) is { } fewer) + { + above = fewer; + left--; + cancelled = true; + } + if (left > 0) + product *= left == 1 ? linear.ToEntity(variables) : MathS.Pow(linear.ToEntity(variables), left); + } + if (cancelled && above is { }) + numerator = above.ToEntity(variables); + product *= power == Number.Integer.One ? remaining.ToEntity(variables) : MathS.Pow(remaining.ToEntity(variables), power); + changed = true; + } + return changed ? (numerator, product) : null; + + // The divisors of a single term: every product of powers of its variables up to theirs, + // with the coefficient one or the term's own. + static IEnumerable MonomialDivisors(MultivariatePolynomial term, int count) + { + var monomial = term.Monomials.First(); + var coefficient = term.CoefficientOf(monomial); + var powers = Enumerable.Range(0, count).Select(v => MultivariatePolynomial.PowerOfMonomial(monomial, v)).ToArray(); + var divisors = new List { 0 }; + for (var v = 0; v < count; v++) + { + var next = new List(); + foreach (var partial in divisors) + for (var k = 0; k <= powers[v]; k++) + next.Add(partial | (k == 0 ? 0 : MultivariatePolynomial.PackMonomial(v, k))); + divisors = next; + if (divisors.Count > 256) + yield break; + } + var magnitude = coefficient.Sign < 0 ? coefficient.Negate() : coefficient; + foreach (var divisor in divisors) + { + yield return MultivariatePolynomial.Term(count, divisor, ERational.One); + if (!magnitude.Equals(ERational.One)) + yield return MultivariatePolynomial.Term(count, divisor, magnitude); + } + } + } + /// /// with each sum among its factors, or the base of a power of one, /// written as writes it: the spelling diff --git a/Sources/Tests/UnitTests/Calculus/SymbolicLinearFactorsOfACubicIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/SymbolicLinearFactorsOfACubicIntegralTest.cs new file mode 100644 index 000000000..443c2c75e --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/SymbolicLinearFactorsOfACubicIntegralTest.cs @@ -0,0 +1,49 @@ +// +// 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 +{ + /// + /// A symbolic cubic whose leading and constant coefficients are single terms, written as the + /// linear factors that divide it before the partial fractions. Rubi's 1.3.1. + /// #718 + /// + [Trait("Area", "Calculus")] + public sealed class SymbolicLinearFactorsOfACubicIntegralTest + { + [Theory] + [InlineData("1/(a*c*h + (b*c*h + a*d*h + a*c*f)*x + (b*d*h + b*c*f + a*d*f)*x^2 + b*d*f*x^3)")] + [InlineData("1/(a*c*h + (b*c*h + a*d*h + a*c*f)*x + (b*d*h + b*c*f + a*d*f)*x^2 + b*d*f*x^3)^2")] + public void OverTheFoundLinears(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 40000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) => e.Substitute("a", 1.3).Substitute("b", 0.4).Substitute("c", 1.9).Substitute("d", 1.6).Substitute("f", 0.7).Substitute("h", 0.45); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -2.1, -0.4, 0.3, 0.6, 1.6, 2.5 }) + { + 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}"); + } + } +}