diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 0909ba42e..b59eb4d7d 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -2963,6 +2963,46 @@ without a common factor | `"(c + d*tan(x))^(3/2)/(a + i*a*tan(x))^2".ToEntity().Integrate("x")` | `integral(...)`; 470,119 characters on the unreleased master | 77,377 characters | | `"(c + d*tan(x))^(1/2)/(a + i*a*tan(x))^2".ToEntity().Integrate("x")` | `integral(...)`; 258,943 characters on the unreleased master | 71,266 characters | +### A rational function of `x/sqrt(P)` for an even quartic `P`, high powers of two linears, and nine more algebraic shapes are integrated + +**Answers where there were none, and shorter ones.** Rubi's algebraic integrands of these shapes were declined or ran +out of time ([#233](https://github.com/asc-community/AngouriMath/issues/233)): +- A rational function beside a half-odd power of `P = a + b x^2 + c x^4`, where it is a rational function of + `t = x/sqrt(P)`: `dt/dx = (a - c x^4)/P^(3/2)`, so `(a - c x^4)/((a d + a e x^2 + c d x^4) sqrt(P))` is + `dt/(d + (a e - b d) t^2)`. Which integrands are is read in `y = x^2` modulo `y^2 - w y + a/c`, `w = y + a/(c y)`; + an elliptic integral does not pass. +- `(c x)^t`, `t` not whole, written `x^t` beside `(c x)^t x^(-t)`, whose derivative is zero wherever it is defined: + `(c x)^(-1 - 2 n - n p) (a + b x^n)^p`. +- Two square roots whose radicands have a constant times an even power of x for their quotient are one root, up to a + factor whose square is constant: `(g + h x) sqrt(a + b x + c x^2)/(a d + b d x + c d x^2)^(3/2)`. +- A polynomial over high whole powers of two linears, under `t = L1/L2`: `x^2/((1 - a x)^16 (1 + a x)^11)`. +- An algebraic function of `1/x` alone, under `x = 1/u`: `1/((a + b/x)^(5/2) (c + d/x)^2)`. +- `x^m (a + b x^n)^p` for whole `m` and `n` and `p ≤ -2`, by the reduction that raises `p` to `-1`: + `1/(x^2 (a + b x^3)^5)`. +- A power of a sum of powers of x beside a symbolic power of a binomial, multiplied out into binomial differentials: + `(a + b x^n)^3 (c + d x^n)^(-4 - 1/n)`. +- A difference of radicands that is one power of x, written as one: `x^4/(sqrt(a + b x) + sqrt(a + c x))^3`. +- A rational function beside an odd root of a linear past the square, split first where a linear below meets the root + at an irrational value: `1/(x^3 (1 - x^2)^(1/3) (3 + x^2))`. +- An odd integrand beside a high power of a quadratic is left to `u = x^2`, and a quartic in one shifted linear to + `u = d + k x`, where the partial fractions spent the budget or wrote pages. + +Answers the unreleased master gave change length: `x^2/((1 - a x)^16 (1 + a x)^11)` was 16,292 characters there, +after twenty seconds, and is 1,380 in a tenth of a second; `1/((a + b/x)^(5/2) (c + d/x)^2)` was 16,675 and is 4,667. +One is longer: `1/(x (1 - x^2)^(1/3) (3 + x^2))`, 342 characters there and 495 now. + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"(a - c*x^4)/((a*d + a*k*x^2 + c*d*x^4)*sqrt(a + b*x^2 + c*x^4))".ToEntity().Integrate("x")` | `integral(...)` | 1,949 characters | +| `"(g + h*x)*sqrt(a + b*x + c*x^2)/(a*d + b*d*x + c*d*x^2)^(3/2)".ToEntity().Integrate("x")` | `integral(...)` | 4,376 characters | +| `"x^2/((1 - a*x)^16*(1 + a*x)^11)".ToEntity().Integrate("x")` | `integral(...)` | 1,380 characters | +| `"1/((a + b/x)^(5/2)*(c + d/x)^2)".ToEntity().Integrate("x")` | `integral(...)` | 4,667 characters | +| `"(c*x)^(-1 - 2*n - n*p)*(a + b*x^n)^p".ToEntity().Integrate("x")` | `integral(...)` | 255 characters | +| `"1/(x^2*(a + b*x^3)^5)".ToEntity().Integrate("x")` | `integral(...)` | 741 characters | +| `"(d + k*x)^4/(a + b*(d + k*x)^2 + c*(d + k*x)^4)".ToEntity().Integrate("x")` | no answer in twenty seconds | 452 characters | +| `"1/(x^5*(a + b*x^2)^10)".ToEntity().Integrate("x")` | `integral(...)` | 652 characters, where the unreleased master gave 67,828 | +| `"1/(x*(1 - x^2)^(1/3)*(3 + x^2))".ToEntity().Integrate("x")` | `integral(...)` | 503 characters, where the unreleased master gave 350 | + ### A polynomial with symbols in it that is a binomial or an even quartic in a shifted variable is integrated in it **Answers where there were none.** `1/(c^2 x^3 + 3 b c x^2 + 3 b^2 x + 3 a b)` was left unevaluated, diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 31dc3156e..6f991d533 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -7909,6 +7909,315 @@ Entity Terms(bool positive, Entity s) ? (constant, cubic) : null; } + /// + /// A rational function beside a half-odd power of an even quartic, R(x) P^(m/2) with + /// P = a + b x^2 + c x^4, as a rational function of t = x/sqrt(P) where it is + /// one: Rubi's (a - c x^4)/((a d + a e x^2 + c d x^4) sqrt(a + b x^2 + c x^4)) is + /// int dt/(d + (a e - b d) t^2), and sqrt(a + b x^2 + c x^4)/(a d - c d x^4) is + /// int dt/(d ((1 - b t^2)^2 - 4 a c t^4)). + /// + /// + /// + /// dt/dx = (a - c x^4)/P^(3/2), so the integrand is H dt for + /// H = R P^((m + 3)/2)/(a - c x^4), a rational function of x. With y = x^2, + /// t^2 = y/(a + b y + c y^2) is 1/(c w + b) for w = y + a/(c y); so + /// H is a rational function of t where it is even, H = h(x^2), and h is a + /// rational function of w -- the same at y and at a/(c y), the two roots of + /// Y^2 - w Y + a/c. That is read by writing the numerator and the denominator of + /// h modulo that quadratic, each as u_1 y + u_0 with the u_i polynomials + /// in w: the quotient is free of y exactly where n_1 d_0 = n_0 d_1, decided as a value, + /// and is then n_1/d_1. An elliptic integral does not pass: 1/sqrt(P) has + /// h(a/(c y)) = -h(y). + /// + /// + /// Each step is an identity for the principal root wherever P is not zero, of either sign: + /// P^(m/2) P^(3/2) = P^((m + 3)/2) for an odd m, t^2 = x^2/P, and the derivative + /// of P^(-1/2) on an interval where P is negative. These were declined, as no rule + /// writes the variable over the root; with numbers for coefficients too, + /// x^2/((1 + x^4) sqrt(x^4 - 1)) is -int t^2/(1 + 4 t^4) dt. Asked late and at + /// the question asked, so that every rule reading the integrand as written goes first, and + /// the answer is differentiated back at sampled points before it is returned. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + internal static Entity? SolveByTheVariableOverTheRootOfAnEvenQuartic(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + // The one radicand, its total power, and the rest, which is to be rational in x. + Entity? quartic = null; + var power = ERational.Zero; + Entity rest = Number.Integer.One; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + if (factor is Powf(var @base, Number.Rational exponent) && exponent is not Number.Integer && @base.ContainsNode(x)) + { + if (!exponent.ERational.Denominator.Equals(EInteger.FromInt32(2)) || quartic is not null && quartic != @base) + return null; + quartic = @base; + power = (underneath ? power.Subtract(exponent.ERational) : power.Add(exponent.ERational)).ToLowestTerms(); + continue; + } + rest = underneath ? rest / factor : rest * factor; + } + if (quartic is null || !power.Denominator.Equals(EInteger.FromInt32(2)) || !power.Numerator.CanFitInInt32() + || rest.Nodes.Any(node => node.ContainsNode(x) && node is not (Sumf or Minusf or Mulf or Divf or Variable) + && !(node is Powf(_, Number.Integer))) + || !TreeAnalyzer.TryGetPolynomial(quartic, x, out var read) + || read.Keys.Any(degree => !(degree.IsZero || degree.Equals(EInteger.FromInt32(2)) || degree.Equals(EInteger.FromInt32(4)))) + || read.Values.Any(coefficient => coefficient.ContainsNode(x)) + || !read.TryGetValue(EInteger.Zero, out var a) || !read.TryGetValue(EInteger.FromInt32(4), out var c)) + return null; + var b = read.TryGetValue(EInteger.FromInt32(2), out var middle) ? middle : Number.Integer.Zero; + static bool Zero(Entity value) => Functions.PartialFractions.IsZeroAsAValue(value); + if (Zero(a) || Zero(c)) + return null; + // H, over one bar, its two sides even in x once a common power of x is taken off. + var raisedBy = Number.Integer.Create(power.Numerator.Add(EInteger.FromInt32(3)).ToInt32Checked() / 2); + var h = Functions.SingleQuotient.Combine(rest * MathS.Pow(quartic, raisedBy) / (a - c * MathS.Pow(x, 4))); + var (top, bottom) = Functions.SingleQuotient.Of(h); + if (!TreeAnalyzer.TryGetPolynomial(top.Expand(), x, out var above) || !TreeAnalyzer.TryGetPolynomial(bottom.Expand(), x, out var below) + || above.Count == 0 || below.Count == 0 + || above.Concat(below).Any(pair => pair.Key.Sign < 0 || !pair.Key.CanFitInInt32() || pair.Value.ContainsNode(x))) + return null; + var shift = EInteger.Min(above.Keys.Min()!, below.Keys.Min()!); + if (above.Keys.Concat(below.Keys).Any(degree => !(degree - shift).IsEven)) + return null; + // The two sides in y = x^2 modulo `y^2 - w y + a/c`, each as `u_1 y + u_0`, the u_i + // polynomials in w by their coefficients. + var ratio = Functions.PartialFractions.Bare((a / c).InnerSimplified); + var (n1, n0) = Reduced(above); + var (d1, d0) = Reduced(below); + if (!PolynomialDifference(PolynomialProduct(n1, d0), PolynomialProduct(n0, d1)).Values.All(Zero)) + return null; + var (overW, underW) = d1.Values.All(Zero) ? (n0, d0) : (n1, d1); + if (underW.Values.All(Zero)) + return null; + // At `w = (1 - b t^2)/(c t^2)`, each side times `(c t^2)` to its degree in w. + var t = Variable.CreateUnique(expr, "t_quartic"); + var (overT, overDegree) = InT(overW); + var (underT, underDegree) = InT(underW); + var inT = Functions.PartialFractions.Bare(Functions.SingleQuotient.Combine( + overT / underT * MathS.Pow(c * MathS.Sqr(t), Number.Integer.Create(underDegree - overDegree))).InnerSimplified); + if (inT.ContainsNode(x)) + return null; + if (Integration.ComputeAsAQuestionOfItsOwn(inT, t, integrateByParts) is not { } inTermsOfT) + return null; + var answer = inTermsOfT.Substitute(t, x / MathS.Sqrt(quartic)).InnerSimplified; + return answer.Nodes.Any(node => node == MathS.NaN) || !Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x) + ? null : answer; + + // A polynomial in w at `w = (1 - b t^2)/(c t^2)` times `(c t^2)^K`, K its degree. + (Entity, int) InT(Dictionary inW) + { + var degree = inW.Where(pair => !Zero(pair.Value)).Select(pair => pair.Key.ToInt32Checked()).DefaultIfEmpty(0).Max(); + Entity sum = Number.Integer.Zero; + foreach (var pair in inW) + { + if (Zero(pair.Value)) + continue; + var k = pair.Key.ToInt32Checked(); + sum += pair.Value * MathS.Pow(1 - b * MathS.Sqr(t), Number.Integer.Create(k)) * MathS.Pow(c * MathS.Sqr(t), Number.Integer.Create(degree - k)); + } + return (sum.Expand().InnerSimplified, degree); + } + + // The polynomial in y = x^2, read off the even powers past the shift, modulo + // `y^2 - w y + a/c`, by Horner's rule from the top: `(p y + q) y + k` is + // `(p w + q) y - p a/c + k` there. As `(p, q)` for `p y + q`. + (Dictionary, Dictionary) Reduced(Dictionary polynomial) + { + var p = new Dictionary(); + var q = new Dictionary(); + var byW = new Dictionary { [EInteger.One] = Number.Integer.One }; + var byRatio = new Dictionary { [EInteger.Zero] = (-ratio).InnerSimplified }; + var highest = (polynomial.Keys.Max()! - shift).ToInt32Checked() / 2; + for (var k = highest; k >= 0; k--) + { + var term = new Dictionary(); + if (polynomial.TryGetValue(shift + EInteger.FromInt32(2 * k), out var coefficient)) + term[EInteger.Zero] = coefficient; + (p, q) = (Bared(PolynomialSum(PolynomialProduct(p, byW), q)), Bared(PolynomialSum(term, PolynomialProduct(p, byRatio)))); + } + return (p, q); + } + + // Each coefficient without the condition the simplification of a quotient attaches. + static Dictionary Bared(Dictionary polynomial) + => polynomial.ToDictionary(pair => pair.Key, pair => Functions.PartialFractions.Bare(pair.Value)); + } + + /// + /// A factor (c x)^t of the integrand, t not whole, written as x^t with + /// (c x)^t x^(-t) in front: Rubi's (c x)^(-1 - 2 n - n p) (a + b x^n)^p and + /// (d x)^(-1 - 2 n (1 + p)) (a^2 + 2 a b x^n + b^2 x^(2 n))^p, whose powers of x alone are + /// answered. + /// + /// + /// The derivative of (c x)^t x^(-t) is t (c x)^t x^(-t) (c/(c x) - 1/x), zero + /// wherever it is defined, since (c x)^(t - 1) is (c x)^t/(c x) for the principal + /// power: it is a constant on each interval of the domain, for either sign of c and of x, and + /// goes in front of the integral of the rest. + /// writes c^t x^t for a positive c only, and reads (c x^n)^t for an n that is + /// not whole. Late and at the question asked, so that a rule reading (c x)^t as + /// written goes first; the answer is differentiated back at sampled points. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveByWritingAPowerOfAMultipleOfXAsAPowerOfX(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + Entity inFront = Number.Integer.One; + Entity written = Number.Integer.One; + var changed = false; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + var rewritten = factor; + if (factor is Powf(var multiple, var t) && multiple != x && !t.ContainsNode(x) && t.Evaled is not Number.Integer + && TreeAnalyzer.TryGetPolyLinear(multiple, x, out var slope, out var intercept) + && !slope.ContainsNode(x) && !intercept.ContainsNode(x) + && Functions.PartialFractions.IsZeroAsAValue(intercept) && !Functions.PartialFractions.IsZeroAsAValue(slope)) + { + rewritten = MathS.Pow(x, t); + var ratio = MathS.Pow(multiple, t) * MathS.Pow(x, -t); + inFront = underneath ? inFront / ratio : inFront * ratio; + changed = true; + } + written = underneath ? written / rewritten : written * rewritten; + } + if (!changed || Integration.ComputeAsAQuestionOfItsOwn(written, x, integrateByParts) is not { } integral) + return null; + var answer = (inFront * integral).InnerSimplified; + return answer.Nodes.Any(node => node == MathS.NaN) || !Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x) + ? null : answer; + } + + /// + /// Two square roots written as one up to a factor whose square is a constant: where the + /// quotient of the radicands is a constant times an even power of x, sqrt(A)/sqrt(B) + /// is K x^k, and where their product is free of a root of x that each holds, + /// sqrt(A) sqrt(B) is K sqrt(A B). Rubi's + /// (g + h x) sqrt(a + b x + c x^2)/(a d + b d x + c d x^2)^(3/2), + /// sqrt(b - a/x^2)/sqrt(a - b x^2) and + /// 1/(sqrt(1 + x) sqrt(-1 - sqrt(x)) sqrt(-1 + sqrt(x))). + /// + /// + /// + /// K is the quotient of the two sides, sqrt(A)/(x^k sqrt(B)) or + /// sqrt(A) sqrt(B)/sqrt(A B), so the rewriting is exact; and its square is + /// A/(x^(2k) B) or 1, free of x, so K is one of two values and constant on each + /// interval where it is continuous: it goes in front of the integral of the rest, as Rubi + /// writes it. For powers A^(p/2) B^(q/2) with odd p and q, the quotient is + /// B^((p + q)/2) (sqrt(A)/sqrt(B))^p, and the product needs p = q. Each was + /// declined: combines two roots of polynomials it + /// can sign, and reads no root inside a radicand, and + /// wants a logarithm to name + /// the radicand. Late and at the question asked, so that every rule reading the roots as + /// written goes first. The rewriting is checked at sampled points rather than the answer + /// differentiated back: the decimal evaluation of the derivative of K, zero, disagreed where + /// a root of x is imaginary, and the answer there was right. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + internal static Entity? SolveByWritingTwoRootsAsOneUpToAConstantFactor(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + var factors = FactorsOfTheIntegrand(expr); + var roots = new List<(int At, Entity Base, int Odd)>(); + for (var i = 0; i < factors.Count; i++) + if (factors[i].Factor is Powf(var @base, Number.Rational exponent) && exponent is not Number.Integer && @base.ContainsNode(x) + && exponent.ERational.Denominator.Equals(EInteger.FromInt32(2)) && exponent.ERational.Numerator.CanFitInInt32()) + { + var odd = exponent.ERational.Numerator.ToInt32Checked(); + roots.Add((i, @base, factors[i].Underneath ? -odd : odd)); + } + if (roots.Count < 2 || roots.Count > 4) + return null; + for (var i = 0; i < roots.Count; i++) + for (var j = 0; j < roots.Count; j++) + { + if (i == j || roots[i].Base == roots[j].Base) + continue; + var (a, p) = (roots[i].Base, roots[i].Odd); + var (b, q) = (roots[j].Base, roots[j].Odd); + Entity written, factor; + if (QuotientAsAConstantTimesAnEvenPowerOfX(a, b) is { } k) + { + // A^(p/2) B^(q/2) = B^((p + q)/2) (K x^k)^p. + written = Raised(b, (p + q) / 2) * Raised(x, k * p); + factor = Raised(MathS.Sqrt(a) / (Raised(x, k) * MathS.Sqrt(b)), p); + } + else if (i < j && p == q && HoldsARootOfX(a) && HoldsARootOfX(b) + && Functions.PartialFractions.Bare((a * b).Expand().InnerSimplified) is var product && !HoldsARootOfX(product) + && product.ContainsNode(x)) + { + // A^(p/2) B^(p/2) = (K sqrt(A B))^p. + written = MathS.Pow(product, Number.Rational.Create(p, 2)); + factor = Raised(MathS.Sqrt(a) * MathS.Sqrt(b) / MathS.Sqrt(product), p); + } + else + continue; + Entity rest = Number.Integer.One; + for (var m = 0; m < factors.Count; m++) + if (m != roots[i].At && m != roots[j].At) + rest = factors[m].Underneath ? rest / factors[m].Factor : rest * factors[m].Factor; + var rewritten = Functions.PartialFractions.Bare((rest * written).InnerSimplified); + // The rewriting checked, rather than the answer differentiated back: the decimal + // evaluation of the derivative of K, zero, disagreed at points where a root of + // x is imaginary, and the answer there is right. + bool rewrittenHolds; + using (MathS.Settings.DowncastingEnabled.Set(false)) + rewrittenHolds = Functions.PartialFractions.HoldsAtSampledPoints(factor * rewritten, expr, x); + if (!rewrittenHolds || Integration.ComputeAsAQuestionOfItsOwn(rewritten, x, integrateByParts) is not { } integral) + return null; + var answer = (factor * integral).InnerSimplified; + return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; + } + return null; + + static Entity Raised(Entity e, int power) => power switch + { + 0 => Number.Integer.One, + 1 => e, + _ => MathS.Pow(e, Number.Integer.Create(power)), + }; + + bool HoldsARootOfX(Entity e) => e.Nodes.Any(node => node is Powf(var @base, Number.Rational r) && r is not Number.Integer && @base.ContainsNode(x)); + + // k where A/B is a constant times x^(2k): each side a polynomial over a power of x, the + // polynomials a constant multiple of one another. + int? QuotientAsAConstantTimesAnEvenPowerOfX(Entity first, Entity second) + { + if (HoldsARootOfX(first) || HoldsARootOfX(second) + || OverAPowerOfX(first) is not var (firstAbove, firstPower) || OverAPowerOfX(second) is not var (secondAbove, secondPower) + || (secondPower - firstPower) % 2 != 0 + || TryReadAsAConstantMultiple(firstAbove, secondAbove, x, allowOne: true) is null) + return null; + return (secondPower - firstPower) / 2; + } + + // A polynomial over a whole power of x, as the polynomial and the power. + (Entity, int)? OverAPowerOfX(Entity e) + { + var (above, below) = Functions.SingleQuotient.Of(Functions.SingleQuotient.Combine(e)); + var power = 0; + foreach (var part in Mulf.LinearChildren(below)) + { + if (!part.ContainsNode(x)) + above /= part; + else if (part == x) + power++; + else if (part is Powf(var @base, Number.Integer whole) && @base == x && whole.EInteger.CanFitInInt32()) + power += whole.EInteger.ToInt32Checked(); + else + return null; + } + return TreeAnalyzer.TryGetPolynomial(above, x, out var read) && read.Keys.All(degree => degree.Sign >= 0) && read.Values.All(value => !value.ContainsNode(x)) + ? (above, power) : null; + } + } + /// /// A root of a quadratic binomial beside another, 1/((A + B x^2)^(1/3) (C + D x^2)) /// with B C + 3 A D = 0 or B C - 9 A D = 0, and @@ -18345,6 +18654,82 @@ private static bool IsAFunctionToRaise(Entity factor, Entity.Variable x) return Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x) ? answer : null; } + /// + /// A power of x beside a whole power below -1 of a binomial with a symbol in it, + /// x^m (a + b x^n)^p for whole m, n and p <= -2, by Rubi's reduction that + /// raises p by one, down to int x^m/(a + b x^n), which is asked on its own: + /// + /// J(p) = -x^(m + 1) (a + b x^n)^(p + 1)/(a n (p + 1)) + (m + 1 + n (p + 1))/(a n (p + 1)) J(p + 1) + /// + /// + /// + /// The derivative of x^(m + 1) B^(p + 1) is x^m B^p ((m + 1 + n (p + 1)) B - a n (p + 1)), + /// which is the step, exact for whole powers. Rubi's 1/(x^2 (a + b x^3)^5), which the + /// root of its perfect square leaves from 1/(x^2 (a^2 + 2 a b x^3 + b^2 x^6)^(5/2)), was + /// declined: the partial fractions over the cube roots of a symbolic binomial gave out past the + /// first power, which is answered. Late, so that every rule answering the power as written goes + /// first; the answer is differentiated back at sampled points. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveABinomialToAWholePowerBelowMinusOneByTheRecurrence(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (!expr.Vars.Any(symbol => symbol != x)) + return null; + Entity constant = Number.Integer.One; + var m = 0; + Entity? binomial = null; + var p = 0; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + if (!factor.ContainsNode(x)) + { + constant = underneath ? constant / factor : constant * factor; + continue; + } + var (@base, whole) = factor is Powf(var raised, Number.Integer { EInteger: var w }) && w.CanFitInInt32() + ? (raised, w.ToInt32Checked()) : (factor, 1); + if (underneath) + whole = -whole; + if (@base == x) + { + m += whole; + continue; + } + if (binomial is not null && binomial != @base || @base is not (Sumf or Minusf)) + return null; + binomial = @base; + p += whole; + } + if (binomial is null || p > -2 || p < -MostBinomialPowersReduced + || !TreeAnalyzer.TryGetPolynomial(binomial, x, out var read) || read.Count != 2 + || !read.TryGetValue(EInteger.Zero, out var a) || read.Keys.Max() is not { } degree || !degree.CanFitInInt32() + || degree.Sign <= 0 || read.Values.Any(value => value.ContainsNode(x)) || Functions.PartialFractions.IsZeroAsAValue(a)) + return null; + var n = degree.ToInt32Checked(); + Entity sum = Number.Integer.Zero; + Entity coefficient = Number.Integer.One; + var raisedX = m + 1 == 0 ? Number.Integer.One : m + 1 == 1 ? (Entity)x : MathS.Pow(x, Number.Integer.Create(m + 1)); + for (var q = p; q <= -2; q++) + { + var below = a * n * (q + 1); + sum += LowestOverTheSymbols(-coefficient / below) * raisedX * MathS.Pow(binomial, Number.Integer.Create(q + 1)); + coefficient = LowestOverTheSymbols(coefficient * (m + 1 + n * (q + 1)) / below); + } + if (!Functions.PartialFractions.IsZeroAsAValue(coefficient)) + { + var last = (m == 0 ? Number.Integer.One : MathS.Pow(x, Number.Integer.Create(m))) / binomial; + if (Integration.ComputeAsAQuestionOfItsOwn(last, x, integrateByParts) is not { } integral) + return null; + sum += coefficient * integral; + } + var answer = (constant * sum).InnerSimplified; + return answer.Nodes.Any(node => node == MathS.NaN) || !Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x) + ? null : answer; + } + + /// The most powers raises. + private const int MostBinomialPowersReduced = 16; + /// /// A binomial differential x^m (a + b x^n)^p with a symbol in its exponents, in /// Chebyshev's third case at a whole (m + 1)/n + p + 1 = -k of at most 0: Rubi's @@ -18575,6 +18960,127 @@ static Entity Simplified(Entity exponent) return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; } + /// + /// A polynomial over high whole powers of two linears, P(x)/(L1^m L2^n) with + /// m + n = k past and P of degree + /// at most k - 2, under t = L1/L2: a polynomial in t over a power of it. + /// Rubi's x^2/((1 - a x)^16 (1 + a x)^11). + /// + /// + /// + /// With D = a1 b2 - a2 b1, x = (a1 - a2 t)/(b2 t - b1), L2 = D/(b2 t - b1) + /// and dx = -D dt/(b2 t - b1)^2, the integrand is -D^(1 - k) t^(-m) S(t) dt for + /// S(t) = sum_j p_j (a1 - a2 t)^j (b2 t - b1)^(k - 2 - j), as + /// has it for exponents + /// that are not whole. Here they are, so each t^(i - m) is L1^(i - m) L2^(m - i) + /// exactly, and the one at i = m - 1 integrates to ln(L1) - ln(L2). + /// + /// + /// The partial fractions take the two linears apart with a symbol in them over k terms, and + /// ran past the budget at k = 27; below the bound they answer in their own form, and + /// this leaves them to it, as it leaves numbers to Hermite's reduction. The answer is + /// differentiated back at sampled points. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + internal static Entity? SolveAPolynomialOverHighWholePowersOfTwoLinears(Entity expr, Entity.Variable x) + { + // With numbers only, Hermite's reduction above has it. + if (!expr.Vars.Any(symbol => symbol != x)) + return null; + Entity constant = Number.Integer.One; + Entity above = Number.Integer.One; + var linears = new List<(Entity Linear, int Power)>(); + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + if (!factor.ContainsNode(x)) + { + constant = underneath ? constant / factor : constant * factor; + continue; + } + var (@base, whole) = factor is Powf(var raised, Number.Integer { EInteger: var w }) && w.CanFitInInt32() + ? (raised, w.ToInt32Checked()) : (factor, 1); + if (underneath) + whole = -whole; + if (whole < 0 && TreeAnalyzer.TryGetPolyLinear(@base, x, out var slope, out var intercept) + && !slope.ContainsNode(x) && !intercept.ContainsNode(x)) + { + var index = linears.FindIndex(pair => pair.Linear == @base); + if (index >= 0) + linears[index] = (@base, linears[index].Power - whole); + else if (linears.Count == 2) + return null; + else + linears.Add((@base, -whole)); + continue; + } + if (whole <= 0) + return null; + above *= factor; + } + if (linears.Count != 2) + return null; + var ((first, m), (second, n)) = (linears[0], linears[1]); + var k = m + n; + if (k < FewestLinearPowersWrittenInTheirQuotient || k > MostLinearPowersWrittenInTheirQuotient + || !TreeAnalyzer.TryGetPolynomial(above.Expand(), x, out var coefficients) + || coefficients.Any(pair => pair.Key.Sign < 0 || pair.Key.CompareTo(EInteger.FromInt32(k - 2)) > 0 || pair.Value.ContainsNode(x)) + || !TreeAnalyzer.TryGetPolyLinear(first, x, out var b1, out var a1) || !TreeAnalyzer.TryGetPolyLinear(second, x, out var b2, out var a2)) + return null; + // A power of x beside one linear is a binomial differential, which the rules for one + // answer at once. + var d = (a1 * b2 - a2 * b1).InnerSimplified; + if (Functions.PartialFractions.IsZeroAsAValue(d) || Functions.PartialFractions.IsZeroAsAValue(b1) || Functions.PartialFractions.IsZeroAsAValue(b2) + || Functions.PartialFractions.IsZeroAsAValue(a1) || Functions.PartialFractions.IsZeroAsAValue(a2)) + return null; + // S by its coefficients, from the binomial theorem: `(a1 - a2 t)^j (b2 t - b1)^(k - 2 - j)` + // has at t^i the sum over r + s = i of C(j, r) C(k - 2 - j, s) a1^(j - r) (-a2)^r b2^s (-b1)^(k - 2 - j - s). + // Expanded as a tree, a nineteenth power of `d t - b` came back unread. + var terms = new Dictionary>(); + foreach (var pair in coefficients) + { + var j = pair.Key.ToInt32Checked(); + var rest = k - 2 - j; + for (var r = 0; r <= j; r++) + for (var q = 0; q <= rest; q++) + { + var term = Number.Integer.Create(Binomial(j, r) * Binomial(rest, q)) * pair.Value + * Raised(a1, j - r) * Raised(-a2, r) * Raised(b2, q) * Raised(-b1, rest - q); + if (!terms.TryGetValue(r + q, out var list)) + terms[r + q] = list = new List(); + list.Add(term); + } + } + var ofT = terms.ToDictionary(pair => EInteger.FromInt32(pair.Key), pair => pair.Value.Skip(1).Aggregate(pair.Value[0], (sum, term) => sum + term)); + Entity sum = Number.Integer.Zero; + foreach (var pair in ofT.OrderBy(pair => pair.Key)) + { + var i = pair.Key.ToInt32Checked(); + var s = LowestOverTheSymbols(pair.Value); + if (s.Evaled is Number.Complex { IsZero: true }) + continue; + sum += i == m - 1 + ? s * (MathS.Ln(first) - MathS.Ln(second)) + : LowestOverTheSymbols(s / (i + 1 - m)) * MathS.Pow(first, Number.Integer.Create(i + 1 - m)) * MathS.Pow(second, Number.Integer.Create(m - 1 - i)); + } + var answer = (LowestOverTheSymbols(-constant * MathS.Pow(d, 1 - k)) * sum).InnerSimplified; + return answer.Nodes.Any(node => node == MathS.NaN) || !Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x) + ? null : answer; + + static Entity Raised(Entity e, int power) => power == 0 ? Number.Integer.One : power == 1 ? e : MathS.Pow(e, Number.Integer.Create(power)); + + static EInteger Binomial(int n, int r) + { + var value = EInteger.One; + for (var i = 1; i <= r; i++) + value = value * EInteger.FromInt32(n - r + i) / EInteger.FromInt32(i); + return value; + } + } + + /// The fewest and the most powers takes. + private const int FewestLinearPowersWrittenInTheirQuotient = 20, MostLinearPowersWrittenInTheirQuotient = 60; + /// /// The most negative sum of the two exponents /// takes: the polynomial @@ -23809,6 +24315,15 @@ Entity TakeOut(Entity node) // a second and ten levels of this ran past the budget. if (Functions.PartialFractions.HoldsTheImaginaryUnit(quadratic) || linear is null && VanishesIdentically(a1)) return null; + // Nor beside a power of x where the integrand is odd and the power of the quadratic high: + // it is then x times a function of x^2, the substitution's under u = x^2, and peeling + // fifteen powers of x off `(b + 2 c x^2)/(x^15 (b + c x^2)^8)` here ran past the budget, + // where in u it is two linears. An even one, and a low power, are left here, where they + // were answered. + if (linear is not null && power >= 6 && VanishesIdentically(a1) && TreeAnalyzer.TryGetPolyLinear(linear, x, out _, out var interceptOfTheLinear) + && VanishesIdentically(interceptOfTheLinear) + && top.Keys.All(degree => degree.IsEven != (linearPower % 2 == 0))) + return null; // A quadratic that factors over linears is the partial fractions' to take apart, and // they do it cheaply: `(d + e x)^10/((d + e x) (a e + c d x))^4`, written out, went from // two seconds past the budget here. @@ -26283,6 +26798,65 @@ private static (Entity Constant, Entity Base, Entity ExponentInU, Variable U, En return (constant, exponential.Base, exponentInU, u, linear, slope, power, beside ? (above.InnerSimplified, below.InnerSimplified, slopePower) : null); } + /// + /// An algebraic function of 1/x alone with a root in it, under x = 1/u: + /// f(1/u) (-1/u^2) du. Rubi's 1/((a + b/x)^(5/2) (c + d/x)^2) is + /// -int du/(u^2 (a + b u)^(5/2) (c + d u)^2), two linears in u beside a root of one. + /// + /// + /// x = 1/u is one to one on each side of zero, and d/dx F(1/x) = -F'(1/x)/x^2, + /// so the answer in u written back holds wherever the integrand is defined. Only where x + /// stands nowhere but below a bar, as b/x^k in two sums or more, each below the bar + /// itself, and some power of such a sum is not whole: a rational function of x is the + /// partial fractions', one sum beside a power of x is read by the binomial rules as it + /// stands, and with a sum above the bar the routes this pre-empts answer, where in u + /// sqrt(a + b u)/(u^2 (c + d u)^3) ran past the budget. wrote the + /// root over its linears, and the partial fractions it reached ran past the budget. At the + /// question asked, before that rule. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveAnAlgebraicFunctionOfTheReciprocalOfX(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + // Algebraic in x, with a power that is not whole of something holding it. + var aRoot = false; + foreach (var node in expr.Nodes) + { + if (node.ContainsNode(x) && node is not (Variable or Sumf or Minusf or Mulf or Divf or Powf)) + return null; + if (node is Powf(var @base, var exponent) && @base.ContainsNode(x) && exponent.ContainsNode(x)) + return null; + } + // A root to a power past the first below the bar: over a square root alone the routes in + // x answered, two of them faster than in u. + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + if (underneath && factor is Powf(var @base, Number.Rational { ERational: var power } exponent) && exponent is not Number.Integer + && @base.ContainsNode(x) && power.CompareTo(ERational.Create(3, 2)) >= 0) + aRoot = true; + // Two sums of 1/x or more, every factor with x below the bar, the shape the linears + // written once took apart: one beside powers of x is answered as it stands, and with a sum + // above the bar the routes it pre-empts answer, where in u two ran past the budget. + if (!aRoot || expr.Nodes.Where(node => node is Sumf or Minusf && node.ContainsNode(x)).Distinct().Count() < 2 + || FactorsOfTheIntegrand(expr).Any(pair => !pair.Underneath && pair.Factor.ContainsNode(x))) + return null; + var u = Variable.CreateUnique(expr, "u_reciprocal"); + var inU = expr.Replace(node => node switch + { + Divf(var above, Powf(var below, Number.Integer k)) when below == x && k.EInteger.Sign > 0 && !above.ContainsNode(x) => above * MathS.Pow(u, k), + Divf(var above, var below) when below == x && !above.ContainsNode(x) => above * u, + Powf(var below, Number.Integer { EInteger.Sign: < 0 } k) when below == x => MathS.Pow(u, (-k).InnerSimplified), + _ => node + }); + if (inU.ContainsNode(x)) + return null; + var integrand = Functions.PartialFractions.Bare((-inU / MathS.Sqr(u)).InnerSimplified); + if (Integration.ComputeAsAQuestionOfItsOwn(integrand, u, integrateByParts) is not { } inTermsOfU) + return null; + var answer = inTermsOfU.Substitute(u, 1 / x).InnerSimplified; + return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; + } + /// /// Sines, cosines and the rest of a function of the reciprocal of a linear, beside a whole /// power of the linear: L^m f(a + b/L^k) with L = c + d x. Under u = 1/L, @@ -30770,9 +31344,30 @@ bool IsAConstantTimesAPowerOfOneMinusUSquared(Entity polynomial, out Entity cons /// in it: with numbers the substitution's rational function is split over the rationals. /// internal static Entity? SolveARationalFunctionBesideARootOfALinearSplitFirst(Entity expr, Entity.Variable x, bool integrateByParts) + => expr.Vars.Any(symbol => symbol != x) ? SplitFirstBesideARootOfALinear(expr, x, integrateByParts, withNumbers: false) : null; + + /// + /// The same with numbers for coefficients, beside a cube root or a higher one, where a linear + /// below the bar has the radicand at its root no rational power of that order: Rubi's + /// 1/(u^2 (1 - u)^(1/3) (3 + u)), which x^3 (1 - x^2)^(1/3) (3 + x^2) is under + /// u = x^2. In v = (1 - u)^(1/3) the linear 3 + u is 4 - v^3, and + /// the partial fractions over the cube root of 4 beside the rest ran past the budget, where + /// split over the rationals first each term beside the root is answered in a tenth of a + /// second. Before the cancellation with functions as indeterminates, which spent five seconds + /// declining it, and the substitution, which this pre-empts only there. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveARationalFunctionBesideAnIrrationalOddRootOfALinearSplitFirst(Entity expr, Entity.Variable x, bool integrateByParts) { - if (!expr.Vars.Any(symbol => symbol != x)) + if (expr.Vars.Any(symbol => symbol != x)) return null; + var (numerator, denominator) = Functions.SingleQuotient.Of(expr); + return SplitFirstBesideARootOfALinear(numerator / denominator, x, integrateByParts, withNumbers: true); + } + + private static Entity? SplitFirstBesideARootOfALinear(Entity expr, Entity.Variable x, bool integrateByParts, bool withNumbers) + { + var linearsBelow = new List(); Entity? radicand = null; var exponent = ERational.Zero; Entity constant = Number.Integer.One; @@ -30803,12 +31398,21 @@ bool IsAConstantTimesAPowerOfOneMinusUSquared(Entity polynomial, out Entity cons { below = below == Number.Integer.One ? written : below * written; belowFactors++; + if (read.Keys.Max()!.Equals(EInteger.One)) + linearsBelow.Add(@base); } else above = above == Number.Integer.One ? written : above * written; } if (radicand is null) return null; + // With numbers, only an odd root past the square of a radicand that is no rational power + // of its order at the root of some linear below: there the substitution's rational function + // is over an irrational root, and elsewhere it is answered whole. + if (withNumbers && (exponent.Denominator.CompareTo(EInteger.FromInt32(3)) < 0 || exponent.Denominator.IsEven + || !exponent.Denominator.CanFitInInt32() + || !linearsBelow.Any(linear => RadicandIsNoRationalPowerAtTheRootOf(linear, exponent.Denominator.ToInt32Checked())))) + return null; // The whole part of the power goes into the rational function, and what is left of it // is a root below the bar, `L^q = L^k / L^(k - q)` with `k` the least whole number past // `q`: `1/(x (1 + x^2) (a + b x)^(3/2))` is split with `a + b x` among the factors and @@ -30844,6 +31448,25 @@ bool IsAConstantTimesAPowerOfOneMinusUSquared(Entity polynomial, out Entity cons sum += integral; } return constant * sum; + + // Whether the radicand, at the root of the linear, is a rational number with no rational + // root of the order: 4 at the root of `3 + u` beside `(1 - u)^(1/3)`. + bool RadicandIsNoRationalPowerAtTheRootOf(Entity linear, int order) + { + if (!TreeAnalyzer.TryGetPolyLinear(linear, x, out var slope, out var intercept) + || radicand!.Substitute(x, -intercept / slope).Evaled is not Number.Rational { ERational: var value } || value.IsZero) + return false; + return !IsAPowerOf(value.Numerator.Abs(), order) || !IsAPowerOf(value.Denominator, order); + } + + static bool IsAPowerOf(EInteger value, int order) + { + var root = EInteger.FromInt64((long)System.Math.Round(System.Math.Pow(EDecimal.FromEInteger(value).ToDouble(), 1.0 / order))); + for (var step = -1; step <= 1; step++) + if (root + EInteger.FromInt32(step) is var candidate && candidate.Sign >= 0 && candidate.Pow(order).Equals(value)) + return true; + return false; + } } /// @@ -30883,16 +31506,31 @@ bool IsAConstantTimesAPowerOfOneMinusUSquared(Entity polynomial, out Entity cons /// x^2 Shi(a + b x) sinh(a + b x) leaves a run of them beside powers of x, and /// written out each term is a search of its own. /// https://github.com/asc-community/AngouriMath/issues/718 + /// + /// And beside one power of a sum of powers of x with a symbol in its exponent and in the + /// powers, kept beside each term of eight at most: Rubi's + /// (a + b x^n)^3 (c + d x^n)^(-4 - 1/n) is four binomial differentials in Chebyshev's + /// third case, each answered, and as a product it was declined. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// /// internal static Entity? SolveByExpandingAProductOfSums(Entity expr, Entity.Variable x, bool integrateByParts) { if (expr is not Mulf || !Integration.AnsweringTheQuestionAsked) return null; var sums = 0; + // And one power of a sum of powers of x with a symbol in its exponent, kept beside each + // term: `(a + b x^n)^3 (c + d x^n)^(-4 - 1/n)` is four terms `x^(k n) (c + d x^n)^(-4 - 1/n)`, + // each a binomial differential in Chebyshev's third case, where the product is not. + Entity? kept = null; + Entity expanded = Number.Integer.One; foreach (var factor in Mulf.LinearChildren(expr)) { if (!factor.ContainsNode(x)) + { + expanded *= factor; continue; + } switch (factor) { case Sumf or Minusf when IsASumOfPowersOfX(factor, x): @@ -30902,15 +31540,27 @@ bool IsAConstantTimesAPowerOfOneMinusUSquared(Entity polynomial, out Entity cons case Entity.Variable when factor == x: case Powf(var @base, var exponent) when @base == x && !exponent.ContainsNode(x): break; + case Powf((Sumf or Minusf) and var sum, var exponent) when kept is null && exponent.Vars.Any() && !exponent.ContainsNode(x) + && sum.Nodes.Any(node => node is Powf(var @base, var power) && @base == x && power.Vars.Any()) && IsASumOfPowersOfX(sum, x): + kept = factor; + continue; default: return null; } + expanded *= factor; } if (sums == 0) return null; - var written = expr.Expand(); + var written = (kept is null ? expr : expanded).Expand(); if (written is not Sumf and not Minusf) return null; + if (kept is not null) + { + var terms = Sumf.LinearChildren(written).ToList(); + if (terms.Count > 8) + return null; + written = terms.Skip(1).Aggregate(terms[0] * kept, (sum, term) => sum + term * kept); + } return Integration.ComputeIndefiniteIntegral(written, x, integrateByParts); // Each term a constant times a power of x, the power's exponent free of x. @@ -33035,6 +33685,12 @@ private static (Entity Numerator, Entity Denominator)? WithASharedPowerCancelled // or a multiple of one root, and neither is this rule's. if (product.Evaled is Number.Complex { IsZero: true } || product.Simplify().Evaled is Number.Complex { IsZero: true }) return null; + // A difference that is one power of x, written as one: `a + b x - (a + c x)` is + // `(b - c) x`, and left as the inner simplification writes it, the cube below + // `x^4 (sqrt(a + b x) - sqrt(a + c x))^3` ran a substitution past the budget. + if (TreeAnalyzer.TryGetPolynomial(product.Expand(), x, out var difference) && difference.Count == 1 + && difference.Single() is { Key: var degree, Value: var coefficient } && degree.Sign > 0 && !coefficient.ContainsNode(x)) + product = Functions.PartialFractions.Bare(coefficient.InnerSimplified) * (degree.Equals(EInteger.One) ? x : MathS.Pow(x, Number.Integer.Create(degree))); // An even power of the conjugate has the root of the product of the radicands // in it, `2 sqrt(A B)` in the square, which over a constant difference leaves a // root of a quadratic where the power as written was answered: the squares over @@ -36331,6 +36987,15 @@ bool HoldsItsValue(Entity factor) if (Mulf.LinearChildren(denominator).Count(factor => factor.ContainsNode(x)) != 1 || Mulf.LinearChildren(denominator).First(factor => factor.ContainsNode(x)) is not (Sumf or Minusf)) return null; + // Nor one written in powers of a single shifted linear, `(d + e x)^4/(a + b (d + e x)^2 + c (d + e x)^4)`: + // in u = d + e x it is the even quartic the rules for a trinomial answer, and the shifted + // linear rule writes it so; split here first, the shift it asked of the split form ran on + // past the budget. + if (expr.Nodes.OfType().Where(power => power.Base.ContainsNode(x)).Select(power => power.Base).Distinct().ToList() is { Count: 1 } bases + && TreeAnalyzer.TryGetPolyLinear(bases[0], x, out _, out var offsetOfTheLinear) && !offsetOfTheLinear.ContainsNode(x) + && !VanishesIdentically(offsetOfTheLinear) + && !expr.Replace(node => node == bases[0] ? Variable.CreateUnique(expr, "u_shifted") : node).ContainsNode(x)) + return null; if (!TreeAnalyzer.TryGetPolynomial(denominator, x, out var read) || read.Keys.Any(key => key.Sign < 0 || key.CompareTo(EInteger.FromInt32(4)) > 0) || !read.ContainsKey(EInteger.FromInt32(4)) || read.Values.Any(value => value.ContainsNode(x)) || !TreeAnalyzer.TryGetPolynomial(numerator, x, out var readAbove) || readAbove.Values.Any(value => value.ContainsNode(x))) diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index b771c2129..5d00475e6 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -726,6 +726,9 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // A whole power of a linear in the inverse sine, cosine or hyperbolic cosine of // `c + d x^2` with `c^2 = 1`, by its recurrence, which closes: by parts it went round. if ((answer = IndefiniteIntegralSolver.SolveAPowerOfAnInverseFunctionOfAQuadraticByTheRecurrence(expr, x)) is { }) return answer; + // An algebraic function of 1/x alone with a root in it, under x = 1/u, before the linears + // are written once: `1/((a + b/x)^(5/2) (c + d/x)^2)` ran past the budget that way. + if ((answer = IndefiniteIntegralSolver.SolveAnAlgebraicFunctionOfTheReciprocalOfX(expr, x, integrateByParts)) is { }) return answer; // A product of powers of quotients of linears with one linear in two of them, written // over each linear once with the constant that takes in front. First, since it answers // the question asked only, and the polynomial term below writes `1/u` as `u^(-1)` and @@ -757,6 +760,9 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // A rational function with numbers for coefficients over a denominator past the twelfth // degree, by Hermite's reduction before anything searches: nothing below answered one. if ((answer = IndefiniteIntegralSolver.SolveAHighDegreeRationalFunctionByHermite(expr, x)) is { }) return answer; + // And a polynomial over high whole powers of two linears with a symbol in them, under + // t = L1/L2, where the partial fractions ran past the budget: `x^2/((1 - a x)^16 (1 + a x)^11)`. + if ((answer = IndefiniteIntegralSolver.SolveAPolynomialOverHighWholePowersOfTwoLinears(expr, x)) is { }) return answer; // A whole power of the square of a polynomial written out, written over its root: the rules // for a power of a binomial read `(a + b x^2)^4` where `(a^2 + 2 a b x^2 + b^2 x^4)^2` ran // past the budget. @@ -1022,6 +1028,11 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // split, since Hearn's `e^(x^2)/x + 2x e^(x^2) ln(x) + ...` is a sum whose terms are // not elementary apart; at the top only, and only for a tower no closed rule reads. if ((answer = IndefiniteIntegralSolver.SolveByRischNormanAnsatz(expr, x)) is { }) return answer; + // A rational function with numbers beside an odd root of a linear past the square, split over + // the rationals first where a linear below meets the root at an irrational value: + // `1/(u^2 (1 - u)^(1/3) (3 + u))`, whose substitution is over the cube root of 4, and which the + // cancellation below spent five seconds declining. + if ((answer = IndefiniteIntegralSolver.SolveARationalFunctionBesideAnIrrationalOddRootOfALinearSplitFirst(expr, x, integrateByParts)) is { }) return answer; // A quotient of polynomials in x and the functions of x in it, cancelled with those // functions as indeterminates and its denominator factored the same way: // `(x^2 + 2x sin(x) + sin(x)^2)/(x + sin(x))^2` is `1`, and nothing else reads it. @@ -1537,6 +1548,20 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // A half-odd power of a constant over a polynomial of either sign, written apart with // its sign in front. After the substitutions, which answer some of these without it. if ((answer = IndefiniteIntegralSolver.SolveByWritingAHalfOddPowerOfAReciprocalWithItsSign(expr, x, integrateByParts)) is { }) return answer; + // A power of x beside a whole power below -1 of a binomial with a symbol in it, by the + // reduction that raises the power to -1: `1/(x^2 (a + b x^3)^5)`, which the partial + // fractions declined past the first power. + if ((answer = IndefiniteIntegralSolver.SolveABinomialToAWholePowerBelowMinusOneByTheRecurrence(expr, x, integrateByParts)) is { }) return answer; + // A rational function beside a half-odd power of an even quartic P, as a rational + // function of x/sqrt(P) where it is one: Rubi's pseudo-elliptic + // `(a - c x^4)/((a d + a e x^2 + c d x^4) sqrt(a + b x^2 + c x^4))`. Late, since it fires + // on nothing a rule reading the root as written answers. + if ((answer = IndefiniteIntegralSolver.SolveByTheVariableOverTheRootOfAnEvenQuartic(expr, x, integrateByParts)) is { }) return answer; + // A factor `(c x)^t` as `x^t`, `(c x)^t x^(-t)` in front, constant on each interval. + if ((answer = IndefiniteIntegralSolver.SolveByWritingAPowerOfAMultipleOfXAsAPowerOfX(expr, x, integrateByParts)) is { }) return answer; + // Two square roots as one, up to a factor whose square is a constant, in front: + // `sqrt(b - a/x^2)/sqrt(a - b x^2)` is that factor over x. + if ((answer = IndefiniteIntegralSolver.SolveByWritingTwoRootsAsOneUpToAConstantFactor(expr, x, integrateByParts)) is { }) return answer; // The sign of a real-valued factor is constant between its zeros, and goes in // front of the antiderivative of the rest. After every rule that reads the sign // where it stands: `cos(x) sgn(sin(x))` is `|sin(x)|` by the substitution, and diff --git a/Sources/Tests/UnitTests/Calculus/ABinomialToAWholePowerBelowMinusOneIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ABinomialToAWholePowerBelowMinusOneIntegralTest.cs new file mode 100644 index 000000000..aa88687c8 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ABinomialToAWholePowerBelowMinusOneIntegralTest.cs @@ -0,0 +1,51 @@ +// +// 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 power of x beside a whole power below -1 of a binomial with a symbol in it, by the reduction + /// that raises the power to -1, and the root of a perfect square that leaves one: Rubi's 1.2.3.2. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ABinomialToAWholePowerBelowMinusOneIntegralTest + { + [Theory] + [InlineData("1/(x^2*(a + b*x^3)^5)")] + [InlineData("1/(x^3*(a + b*x^3)^5)")] + [InlineData("1/(x^2*(a^2 + 2*a*b*x^3 + b^2*x^6)^(5/2))")] + [InlineData("1/(x^3*(a^2 + 2*a*b*x^3 + b^2*x^6)^(5/2))")] + public void RaisedToTheFirstPower(string integrand) + { + 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); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/APowerOfAMultipleOfXAsAPowerOfXIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/APowerOfAMultipleOfXAsAPowerOfXIntegralTest.cs new file mode 100644 index 000000000..7fc06d3e8 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/APowerOfAMultipleOfXAsAPowerOfXIntegralTest.cs @@ -0,0 +1,53 @@ +// +// 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 factor (c x)^t written as x^t, with (c x)^t x^(-t) in front, whose + /// derivative is zero: Rubi's 1.1.3.2 and 1.2.3.2 with a symbol for the multiple of x. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class APowerOfAMultipleOfXAsAPowerOfXIntegralTest + { + [Theory] + [InlineData("(c*x)^(-1 - 2*n - n*p)*(a + b*x^n)^p")] + [InlineData("(c*x)^(-1 - 3*n - n*p)*(a + b*x^n)^p")] + [InlineData("(c*x)^(-1 - 4*n - n*p)*(a + b*x^n)^p")] + [InlineData("(d*x)^(-1 - 2*n*(1 + p))*(a^2 + 2*a*b*x^n + b^2*x^(2*n))^p")] + public void WrittenAsAPowerOfX(string integrand) + { + 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}"); + // n past the wrap, so that x^n at a negative x is not confined to a half-plane. + Entity Pinned(Entity e) => e.Substitute("a", 1.3).Substitute("b", 0.7).Substitute("c", -0.4).Substitute("d", 0.3) + .Substitute("n", 2.3).Substitute("p", 0.6); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/AProductOfSumsBesideASymbolicPowerIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AProductOfSumsBesideASymbolicPowerIntegralTest.cs new file mode 100644 index 000000000..6e30ebdb6 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AProductOfSumsBesideASymbolicPowerIntegralTest.cs @@ -0,0 +1,52 @@ +// +// 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 whole power of a binomial in x^n beside a power of another with a symbol in its + /// exponent, written out with the power kept beside each term: each term is a binomial + /// differential in Chebyshev's third case. Rubi's 1.1.3.3. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AProductOfSumsBesideASymbolicPowerIntegralTest + { + [Theory] + [InlineData("(a + b*x^n)^3*(c + d*x^n)^(-4 - 1/n)")] + [InlineData("(a + b*x^n)^2*(c + d*x^n)^(-4 - 1/n)")] + public void WrittenOutBesideThePower(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 10000, $"{text.Length} characters of answer for {integrand}"); + // n past the wrap, so that x^n at a negative x is not confined to a half-plane. + Entity Pinned(Entity e) => e.Substitute("a", 1.3).Substitute("b", 0.7).Substitute("c", 0.4).Substitute("d", 0.3) + .Substitute("n", 2.3); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/AQuarticInOneShiftedLinearIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AQuarticInOneShiftedLinearIntegralTest.cs new file mode 100644 index 000000000..de298ee17 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AQuarticInOneShiftedLinearIntegralTest.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 power of a linear over a trinomial in its square, written in the linear rather than split + /// over the quartic's quadratics first: Rubi's 1.2.3.2, which ran past the budget split. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AQuarticInOneShiftedLinearIntegralTest + { + [Theory] + [InlineData("(d + k*x)^4/(a + b*(d + k*x)^2 + c*(d + k*x)^4)")] + [InlineData("(d + k*x)^2/(a + b*(d + k*x)^2 + c*(d + k*x)^4)")] + public void WrittenInTheShiftedLinear(string integrand) + { + 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.9).Substitute("k", 1.1); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/ARationalFunctionBesideAnIrrationalOddRootIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ARationalFunctionBesideAnIrrationalOddRootIntegralTest.cs new file mode 100644 index 000000000..afd0d1a68 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ARationalFunctionBesideAnIrrationalOddRootIntegralTest.cs @@ -0,0 +1,53 @@ +// +// 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 rational function with numbers beside a cube root of a linear, split over the rationals + /// first where a linear below meets the root at an irrational value: Rubi's 1.1.2.4 and 1.1.3.4 + /// under the substitution of a power of x, which ran past the budget over the cube root of 4. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ARationalFunctionBesideAnIrrationalOddRootIntegralTest + { + [Theory] + [InlineData("1/(x^3*(1 - x^2)^(1/3)*(3 + x^2))")] + [InlineData("1/(x^5*(1 - x^2)^(1/3)*(3 + x^2))")] + [InlineData("1/(x^3*(1 - x^2)^(1/3)*(3 + x^2)^2)")] + [InlineData("1/(x^5*(1 - x^2)^(1/3)*(3 + x^2)^2)")] + [InlineData("1/(x^4*(1 - x^3)^(4/3)*(1 + x^3))")] + public void SplitOverTheRationalsFirst(string integrand) + { + 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; + 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.6, 0.8, 1.7, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/ASumOfRootsWithAMonomialDifferenceIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ASumOfRootsWithAMonomialDifferenceIntegralTest.cs new file mode 100644 index 000000000..ebf0a3a0a --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ASumOfRootsWithAMonomialDifferenceIntegralTest.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 power of x over a power of a sum of two roots whose radicands differ by a multiple of x, + /// rationalised with the difference written as that multiple: Rubi's 1.3.2. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ASumOfRootsWithAMonomialDifferenceIntegralTest + { + [Theory] + [InlineData("x^4/(sqrt(a + b*x) + sqrt(a + c*x))^3")] + [InlineData("x^2/(sqrt(a + b*x) + sqrt(a + c*x))^2")] + public void RationalisedOverAPowerOfX(string integrand) + { + 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); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/AnAlgebraicFunctionOfTheReciprocalOfXIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnAlgebraicFunctionOfTheReciprocalOfXIntegralTest.cs new file mode 100644 index 000000000..0bb84eed2 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnAlgebraicFunctionOfTheReciprocalOfXIntegralTest.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 +{ + /// + /// An algebraic function of 1/x alone with a root in it, under x = 1/u: Rubi's 1.1.3.3, which ran + /// past the budget with its linears written once. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnAlgebraicFunctionOfTheReciprocalOfXIntegralTest + { + [Theory] + [InlineData("1/((a + b/x)^(5/2)*(c + d/x)^2)")] + [InlineData("1/((a + b/x)^(5/2)*(c + d/x)^3)")] + public void UnderTheReciprocal(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 10000, $"{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.9); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -2.9, -1.2, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/AnOddRationalFunctionBesideAPowerOfXIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnOddRationalFunctionBesideAPowerOfXIntegralTest.cs new file mode 100644 index 000000000..899ad4fb3 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnOddRationalFunctionBesideAPowerOfXIntegralTest.cs @@ -0,0 +1,51 @@ +// +// 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 over a power of x and a whole power of a binomial in x^2, left to + /// the substitution under u = x^2 rather than lowered a power of x at a time: Rubi's 1.1.2.2 and + /// 1.1.3.4, which ran past the budget. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnOddRationalFunctionBesideAPowerOfXIntegralTest + { + [Theory] + [InlineData("1/(x^7*(a + b*x^2)^10)")] + [InlineData("(b + 2*c*x^2)/(x^15*(b + c*x^2)^8)")] + [InlineData("1/(x^5*(a + b*x^2)^10)")] + public void AnsweredUnderTheSquare(string integrand) + { + 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); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/HighWholePowersOfTwoLinearsIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/HighWholePowersOfTwoLinearsIntegralTest.cs new file mode 100644 index 000000000..5aa106fca --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/HighWholePowersOfTwoLinearsIntegralTest.cs @@ -0,0 +1,51 @@ +// +// 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 polynomial over high whole powers of two linears with a symbol in them, under the + /// quotient of the linears: Rubi's 1.1.1.3, where the partial fractions ran past the budget. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class HighWholePowersOfTwoLinearsIntegralTest + { + [Theory] + [InlineData("x^2/((1 - a*x)^16*(1 + a*x)^11)")] + [InlineData("x^2/((1 - a*x)^14*(1 + a*x)^10)")] + [InlineData("(c + x)^3/((a + b*x)^10*(c + d*x)^11)")] + [InlineData("1/((a + b*x)^10*(c + d*x)^11)")] + public void WrittenInTheQuotientOfTheLinears(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 10000, $"{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); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 0.3, 0.6, 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}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/TheVariableOverTheRootOfAnEvenQuarticIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/TheVariableOverTheRootOfAnEvenQuarticIntegralTest.cs new file mode 100644 index 000000000..f3aa02c86 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/TheVariableOverTheRootOfAnEvenQuarticIntegralTest.cs @@ -0,0 +1,62 @@ +// +// 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 rational function beside a half-odd power of an even quartic P, as a rational function + /// of t = x/sqrt(P): the pseudo-elliptic integrals of Rubi's 1.2.2.7 and 1.3.2, which + /// no rule wrote in t. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class TheVariableOverTheRootOfAnEvenQuarticIntegralTest + { + [Theory] + [InlineData("(a - c*x^4)/((a*d + a*k*x^2 + c*d*x^4)*sqrt(a + b*x^2 + c*x^4))")] + [InlineData("(a - c*x^4)/((a*d + a*k*x^2 + c*d*x^4)*sqrt(a - b*x^2 + c*x^4))")] + [InlineData("(a - c*x^4)/((a*k + c*d*x^2)*(d + k*x^2)*sqrt(a + b*x^2 + c*x^4))")] + [InlineData("sqrt(a + b*x^2 + c*x^4)/(a*d - c*d*x^4)")] + [InlineData("sqrt(a + b*x^2 - c*x^4)/(a*d + c*d*x^4)")] + [InlineData("sqrt(a + b*x^4)/(a*c - b*c*x^4)")] + [InlineData("sqrt(a - b*x^4)/(a*c + b*c*x^4)")] + [InlineData("x^2/((1 + x^4)*sqrt(-1 + x^4))")] + public void WrittenInTheVariableOverTheRoot(string integrand) + { + 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("k", 0.9); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -2.9, -1.2, -0.7, 0.3, 0.8, 1.7, 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}"); + } + + [Theory] + [InlineData("1/sqrt(a + b*x^2 + c*x^4)")] + [InlineData("x^2/sqrt(1 + x^4)")] + public void AnEllipticIntegralIsLeftAlone(string integrand) + => Assert.Contains("integral(", integrand.ToEntity().Integrate("x").Stringize()); + } +} diff --git a/Sources/Tests/UnitTests/Calculus/TwoRootsAsOneUpToAConstantFactorIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/TwoRootsAsOneUpToAConstantFactorIntegralTest.cs new file mode 100644 index 000000000..478c886c9 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/TwoRootsAsOneUpToAConstantFactorIntegralTest.cs @@ -0,0 +1,53 @@ +// +// 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 +{ + /// + /// Two square roots written as one up to a factor whose square is a constant, the factor in + /// front: where their radicands' quotient is a constant times an even power of x, or their + /// product is free of a root of x each holds. Rubi's 1.1.3.3, 1.2.1.6 and 1.3.2. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class TwoRootsAsOneUpToAConstantFactorIntegralTest + { + [Theory] + [InlineData("(g + h*x)*sqrt(a + b*x + c*x^2)/(a*d + b*d*x + c*d*x^2)^(3/2)")] + [InlineData("sqrt(b - a/x^2)/sqrt(a - b*x^2)")] + [InlineData("1/(sqrt(1 + x)*sqrt(-1 - sqrt(x))*sqrt(-1 + sqrt(x)))")] + [InlineData("1/(sqrt(a^2 + b^2*x)*sqrt(a - b*sqrt(x))*sqrt(a + b*sqrt(x)))")] + public void WrittenAsOneRoot(string integrand) + { + 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("g", 0.9).Substitute("h", 1.1); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 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}"); + } + } +}