diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 54b1efb85..d8d706432 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -4025,6 +4025,28 @@ in both quadratics, the cube's coefficients grow past what the reduction takes, | `"(7 + 13*x)/((5 + x + 2*x^2)^3*sqrt(2 + x + 3*x^2))".ToEntity().Integrate("x")` | `integral(...)` | an antiderivative | | `"sqrt(a + a*sec(x))/(c + d*sec(x))^2".ToEntity().Integrate("x")` | `integral(...)` | an antiderivative in `tan(x/2)` | +### Cube roots of two linears beside a power of a third, and an odd half power of a linear over a power of `x (b + c x)`, are integrated + +**Answers where there were none.** `(a + b x)^(1/3) (c + d x)^(2/3)/(e + f x)^2` and its kin, Rubi's 1.1.1.3, and +`(A + B x) sqrt(d + e x)/(b x + c x^2)^3` and its kin, 1.2.1.3, were declined. The first are rational in +`t = (a + b x)^(1/3)/(c + d x)^(1/3)` wherever the two powers add to a whole number, and exactly so for every complex +`x`; the substitution that reaches it left `k (b - d t^3) + f (c t^3 - a)` where a binomial in `t^3` stands, and the +rational rules ran out of time on it. It is now written as that binomial from the coefficients, each compound +coefficient named by a symbol of its own and written back into the answer, where the denominator in `t` is not +squarefree, as for a third linear below the bar at a power of `-2` or less. The second is rational in +`t = sqrt(d + e x)` with the quadratic kept as `x (b + c x)`, which the substitution that reaches it, late in the +chain, multiplied out; it is asked so now, named the same way, for the root above the bar and a third power or +higher. Each declines where its linears share a root, decided as a value +([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"(a + b*x)^(1/3)*(c + d*x)^(2/3)/(k + f*x)^2".ToEntity().Integrate("x")` | `integral(...)` | 2,243 characters | +| `"1/((a + b*x)^(1/3)*(c + d*x)^(2/3)*(k + f*x)^3)".ToEntity().Integrate("x")` | `integral(...)` | 2,242 characters | +| `"(a + b*x)^(4/3)/((c + d*x)^(4/3)*(k + f*x)^2)".ToEntity().Integrate("x")` | `integral(...)` | 1,327 characters | +| `"(d + k*x)^(1/2)/(b*x + c*x^2)^3".ToEntity().Integrate("x")` | `integral(...)` | 10,875 characters | +| `"(A + B*x)*(d + k*x)^(1/2)/(b*x + c*x^2)^3".ToEntity().Integrate("x")` | `integral(...)` | 14,763 characters | + ### Powers of two linears whose exponents sum to a whole number below `-2` are integrated **Answers where there were none.** `(a + b x)^m (c + d x)^(-3 - m)`, `x^(n - 4)/(a + b x)^n` and diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 8dc55b281..dba1f0072 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -13559,6 +13559,279 @@ static Entity WithWholePowersOfTheSecond(List<(EInteger Power, Entity Coefficien return Integration.ComputeAsTheSameQuestion(rewritten.InnerSimplified, x, integrateByParts); } + /// + /// Cube roots of two linears whose powers add up to a whole number, beside a whole power + /// of a third linear, (a + b x)^(p/3) (c + d x)^(r/3) (g + h x)^n with p + r + /// a multiple of 3, as the rational function the quotient of the two roots makes of it, + /// t = (a + b x)^(1/3)/(c + d x)^(1/3), with its compound coefficients named. + /// + /// + /// + /// With D = a d - b c, W = b - d t^3 and s = (p + r)/3: x is + /// (c t^3 - a)/W, c + d x = -D/W, g + h x = N/W for + /// N = (b g - a h) + (c h - d g) t^3, and dx = -3 D t^2 dt/W^2; and the two + /// roots are t^p (c + d x)^s. So the integrand is + /// 3 (-D)^(s + 1) t^(p + 2) N^n / W^(s + n + 2), rational in t, and exactly so: + /// a whole power of the quotient of the principal roots is the quotient of their powers + /// for every complex x, as says. + /// + /// + /// That rule makes the same substitution and declined these, past the budget: its + /// rational function is written as the substitution leaves it, with + /// k (b - d t^3) + f (c t^3 - a) where N is, and the rational rules ran out + /// of time on t^3/((b k - a f + (c f - d k) t^3) (b - d t^3)^2) where + /// t^3/((p + q t^3) (b - d t^3)^2) is answered in a quarter of a second. Here + /// N and W are written as binomials in t^3 from the coefficients, and a + /// coefficient that is a compound of the symbols is named by a symbol of its own and + /// written back into the answer. Naming loses no relation the decomposition needs: + /// W and N have the resultant h D, zero only where the third linear + /// is constant or the two under the roots are one, and those, and a third linear that is + /// a multiple of one of the two, are declined, each decided as a value. + /// + /// + /// Only cube roots, and only a third linear below the bar where the denominator in + /// t is not squarefree, n <= -2 or W^2 beside N: Rubi's + /// (a + b x)^m (c + d x)^n (e + f x)^p rows with m, n thirds. Of the + /// rest, the third linear above the bar or below it once beside W at most once, the + /// substitution rule answers most, and they are left to it as they are. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + internal static Entity? SolveCubeRootsOfTwoLinearsBesideAWholePowerOfAThird(Entity expr, Entity.Variable x, bool integrateByParts) + { + Entity constant = Number.Integer.One; + var powers = new List<(Entity Base, ERational Exponent)>(); + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + if (!factor.ContainsNode(x)) + { + constant = underneath ? constant / factor : constant * factor; + continue; + } + var (@base, exponent) = factor switch + { + Powf(var powerBase, Number.Rational power) => (powerBase, power.ERational), + Powf => (factor, (ERational?)null), + _ => (factor, ERational.One), + }; + if (exponent is null) + return null; + if (underneath) + exponent = exponent.Negate(); + var index = powers.FindIndex(pair => pair.Base == @base); + if (index >= 0) + powers[index] = (@base, powers[index].Exponent.Add(exponent)); + else + powers.Add((@base, exponent)); + } + powers.RemoveAll(pair => pair.Exponent.IsZero); + var roots = powers.Where(pair => !pair.Exponent.IsInteger()).ToList(); + var whole = powers.Where(pair => pair.Exponent.IsInteger()).ToList(); + if (roots.Count != 2 || whole.Count != 1) + return null; + var three = ERational.FromInt32(3); + var thirds = roots.Select(pair => pair.Exponent.Multiply(three)).ToList(); + if (thirds.Any(third => !third.IsInteger())) + return null; + var read = new[] { thirds[0].ToEInteger(), thirds[1].ToEInteger(), whole[0].Exponent.ToEInteger() }; + if (read.Any(power => !power.CanFitInInt32())) + return null; + var (p, r, n) = (read[0].ToInt32Checked(), read[1].ToInt32Checked(), read[2].ToInt32Checked()); + if ((p + r) % 3 != 0 || System.Math.Abs(p) > MaximumCubeRootNumerator || System.Math.Abs(r) > MaximumCubeRootNumerator + || System.Math.Abs(n) > MaximumThirdLinearPower) + return null; + var s = (p + r) / 3; + // Only the third linear below the bar with a repeated factor below it, N^2 or W^2 beside + // N: there the substitution rule ran out of time. Above, or once below beside W at most + // once, it answers most of them, and the question is left to it. + if (n >= 0 || n == -1 && s + n + 2 < 2) + return null; + var (first, second, third) = (roots[0].Base, roots[1].Base, whole[0].Base); + if (!TreeAnalyzer.TryGetPolyLinear(first, x, out var b, out var a) + || !TreeAnalyzer.TryGetPolyLinear(second, x, out var d, out var c) + || !TreeAnalyzer.TryGetPolyLinear(third, x, out var h, out var g) + || VanishesIdentically(b) || VanishesIdentically(d) || VanishesIdentically(h)) + return null; + var determinant = a * d - b * c; + var (free, slope) = (b * g - a * h, c * h - d * g); + if (VanishesIdentically(determinant) || VanishesIdentically(free) || VanishesIdentically(slope)) + return null; + + var t = Variable.CreateUnique(expr, "t_cube"); + var named = new List<(Variable Name, Entity Value)>(); + Entity Named(Entity value) + { + if (value is Number or Variable || value is Mulf(Number, Variable) || value is Powf(Variable, Number)) + return value; + var name = Variable.CreateUnique(expr + t + named.Aggregate((Entity)Number.Integer.Zero, (sum, pair) => sum + pair.Name), "k_named"); + named.Add((name, value)); + return name; + } + var cube = MathS.Pow(t, 3); + var w = Named(Reduced(b)) + Named(Reduced(-d)) * cube; + var nOfT = Named(Reduced(free)) + Named(Reduced(slope)) * cube; + // With nothing to name the rational function is the substitution rule's, which answers it. + if (named.Count == 0) + return null; + // One quotient, the powers whole on both sides, as the rational rules read it. + Entity above = Number.Integer.One, below = Number.Integer.One; + foreach (var (factor, power) in new[] { ((Entity)t, p + 2), (w, -(s + n + 2)), (nOfT, n) }) + { + if (power > 0) + above *= power == 1 ? factor : MathS.Pow(factor, power); + else if (power < 0) + below *= power == -1 ? factor : MathS.Pow(factor, -power); + } + var inT = (above / below).InnerSimplified; + if (Integration.ComputeAsAQuestionOfItsOwn(inT, t, integrateByParts) is not { } answer + || answer.Nodes.Any(node => node == MathS.NaN)) + return null; + foreach (var (name, value) in named) + answer = answer.Substitute(name, value); + var root = Number.Rational.Create(1, 3); + return Reduced(3 * constant * MathS.Pow(-determinant, s + 1)) * answer.Substitute(t, MathS.Pow(first, root) / MathS.Pow(second, root)); + + static Entity Reduced(Entity constant) + => constant.Vars.Any() ? Functions.PartialFractions.InLowestTermsOverTheSymbols(constant) : constant.InnerSimplified; + } + + /// + /// The largest numerator of a third takes. + /// + private const int MaximumCubeRootNumerator = 12; + + /// + /// The largest power of the third linear takes. + /// + private const int MaximumThirdLinearPower = 8; + + /// + /// A polynomial times an odd half power of a linear over a whole power of a quadratic with + /// a root at zero, P(x) (d + e x)^(j/2)/(b x + c x^2)^m with j odd and + /// positive and m >= 3, under t = sqrt(d + e x) with the quadratic written + /// as its two linears, x (b + c x). + /// + /// + /// + /// With x = (t^2 - d)/e, b + c x = (c t^2 + r)/e for r = b e - c d, and + /// dx = 2 t dt/e, the integrand is + /// 2 e^(2m - 1) P((t^2 - d)/e) t^(j + 1)/((t^2 - d)^m (c t^2 + r)^m), and exactly so: + /// t^j is (d + e x)^(j/2) for the principal root and every complex x. + /// The two factors below the bar have the resultant b e, and a compound coefficient + /// is named by a symbol of its own and written back into the answer, as in + /// . + /// + /// + /// Rubi's (A + B x) (d + e x)^(j/2)/(b x + c x^2)^3 ran past the budget: the rational + /// function in t is answered in a second, and the substitution that reaches it is + /// late in the chain and multiplies the quadratic out. Only the root above the bar and a + /// third power or higher, with a symbol to name: the root below the bar, the second power, + /// and numbers for coefficients are answered as they are, and left to that. Declined where + /// d or r is zero as a value, the linear then sharing a root with the quadratic. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + internal static Entity? SolveAnOddHalfPowerOfALinearOverAPowerOfAQuadraticWithARootAtZero(Entity expr, Entity.Variable x, bool integrateByParts) + { + Entity constant = Number.Integer.One; + Entity polynomial = Number.Integer.One; + (Entity Base, ERational Exponent)? root = null; + (Entity Base, int Power)? quadratic = null; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + if (!factor.ContainsNode(x)) + constant = underneath ? constant / factor : constant * factor; + else if (factor is Powf(var radicand, Number.Rational power) && power is not Number.Integer) + { + if (root is not null || underneath) + return null; + root = (radicand, power.ERational); + } + else if (underneath) + { + if (quadratic is not null || factor is not Powf(var @base, Number.Integer { EInteger: var m }) || !m.CanFitInInt32()) + return null; + quadratic = (@base, m.ToInt32Checked()); + } + else + polynomial *= factor; + } + if (root is not var (linear, half) || quadratic is not var (q, n) || half.Sign <= 0 + || half.Multiply(ERational.FromInt32(2)) is not { } twice || !twice.IsInteger() || !twice.ToEInteger().CanFitInInt32() + || n < 3 || n > MaximumQuadraticPowerUnderARoot) + return null; + var j = twice.ToEInteger().ToInt32Checked(); + if (j > MaximumHalfPowerOfTheLinear + || !TreeAnalyzer.TryGetPolyLinear(linear, x, out var e, out var d) + || !TreeAnalyzer.TryGetPolyQuadratic(q, x, out var c, out var b, out var free) + || !TreeAnalyzer.TryGetPolynomial(polynomial, x, out var read) + || read.Any(pair => pair.Key.Sign < 0 || pair.Value.ContainsNode(x))) + return null; + var degree = read.Keys.Max()!.ToInt32Checked(); + var r = b * e - c * d; + if (degree > MaximumQuadraticPowerUnderARoot || !VanishesIdentically(free) || VanishesIdentically(b) || VanishesIdentically(c) + || VanishesIdentically(e) || VanishesIdentically(d) || VanishesIdentically(r)) + return null; + + var t = Variable.CreateUnique(expr, "t_root"); + var named = new List<(Variable Name, Entity Value)>(); + Entity Named(Entity value) + { + if (value is Number or Variable || value is Mulf(Number, Variable) || value is Powf(Variable, Number)) + return value; + var name = Variable.CreateUnique(expr + t + named.Aggregate((Entity)Number.Integer.Zero, (sum, pair) => sum + pair.Name), "k_named"); + named.Add((name, value)); + return name; + } + // P((t^2 - d)/e) e^degree as a polynomial in t^2: the coefficient of t^(2l) is the sum + // over i of p_i e^(degree - i) binomial(i, l) (-d)^(i - l). + Entity above = Number.Integer.Zero; + for (var l = 0; l <= degree; l++) + { + Entity coefficient = Number.Integer.Zero; + foreach (var (i, p) in read.Select(pair => (pair.Key.ToInt32Checked(), pair.Value)).Where(pair => pair.Item1 >= l)) + coefficient += p * MathS.Pow(e, degree - i) * Binomial(i, l) * MathS.Pow(-d, i - l); + // t^0 is written as nothing: InnerSimplified makes it 1 provided t is not 0. + if (!VanishesIdentically(coefficient)) + above += l == 0 ? Named(Reduced(coefficient)) : Named(Reduced(coefficient)) * MathS.Pow(t, 2 * l); + } + var square = MathS.Sqr(t); + var below = MathS.Pow(square + Named(Reduced(-d)), n) * MathS.Pow(Named(Reduced(c)) * square + Named(Reduced(r)), n); + // With nothing to name the coefficients are numbers, and the rules after this answer them. + if (named.Count == 0) + return null; + var inT = (MathS.Pow(t, j + 1) * above / below).InnerSimplified; + if (Integration.ComputeAsAQuestionOfItsOwn(inT, t, integrateByParts) is not { } answer + || answer.Nodes.Any(node => node == MathS.NaN)) + return null; + foreach (var (name, value) in named) + answer = answer.Substitute(name, value); + return Reduced(2 * constant * MathS.Pow(e, 2 * n - 1 - degree)) * answer.Substitute(t, MathS.Sqrt(linear)); + + static Entity Reduced(Entity constant) + => constant.Vars.Any() ? Functions.PartialFractions.InLowestTermsOverTheSymbols(constant) : constant.InnerSimplified; + + static Entity Binomial(int i, int l) + { + var value = EInteger.One; + for (var k = 0; k < l; k++) + value = value.Multiply(i - k).Divide(k + 1); + return Number.Integer.Create(value); + } + } + + /// + /// The largest power of the quadratic, and degree of the polynomial, + /// takes. + /// + private const int MaximumQuadraticPowerUnderARoot = 8; + + /// + /// The largest odd j of (d + e x)^(j/2) + /// takes. + /// + private const int MaximumHalfPowerOfTheLinear = 15; + /// /// A polynomial over a power of a linear beside the square roots of two others, /// P(x) (a + b x)^(m/2) (c + d x)^(n/2) / (g + h x)^k with m and n diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index e978ecdcb..aeebf1123 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -732,6 +732,13 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // A linear and a quadratic to one power as one binomial, `(1 + x) (1 - x + x^2)` as // `1 + x^3`, which the rules for a power of a binomial answer. if ((answer = IndefiniteIntegralSolver.SolveALinearAndAQuadraticToOnePowerAsOneBinomial(expr, x, integrateByParts)) is { }) return answer; + // Cube roots of two linears beside a whole power of a third, under the quotient of the + // roots with the coefficients named, before the polynomial term takes it a level down to + // the substitution there, whose rational function ran past the budget. + if ((answer = IndefiniteIntegralSolver.SolveCubeRootsOfTwoLinearsBesideAWholePowerOfAThird(expr, x, integrateByParts)) is { }) return answer; + // And a root of a linear over a third power or higher of `b x + c x^2`, under the root, + // the quadratic kept as its two linears: the substitution there multiplies it out. + if ((answer = IndefiniteIntegralSolver.SolveAnOddHalfPowerOfALinearOverAPowerOfAQuadraticWithARootAtZero(expr, x, integrateByParts)) is { }) return answer; if ((answer = IndefiniteIntegralSolver.SolveByWritingEachLinearOnce(expr, x, integrateByParts)) is { }) return answer; // The same with a symbol for an exponent, where the product is one power that is not // whole beside whole ones: `e^atanh(a x) x^m sqrt(1 - a^2 x^2)` is `K x^m (a x + 1)`, and ran diff --git a/Sources/Tests/UnitTests/Calculus/AnOddHalfPowerOfALinearOverAPowerOfAQuadraticWithARootAtZeroIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnOddHalfPowerOfALinearOverAPowerOfAQuadraticWithARootAtZeroIntegralTest.cs new file mode 100644 index 000000000..13268e5cc --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnOddHalfPowerOfALinearOverAPowerOfAQuadraticWithARootAtZeroIntegralTest.cs @@ -0,0 +1,56 @@ +// +// 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 times an odd half power of a linear over a third power or higher of + /// b x + c x^2, under t = sqrt(d + k x) with the quadratic kept as + /// x (b + c x): the substitution reached the other way multiplies it out, and these + /// ran past the budget. Rubi's 1.2.1.2 and 1.2.1.3. Points on both sides of the zeros, where + /// the root is imaginary. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnOddHalfPowerOfALinearOverAPowerOfAQuadraticWithARootAtZeroIntegralTest + { + [Theory] + [InlineData("(d + k*x)^(1/2)/(b*x + c*x^2)^3")] + [InlineData("(d + k*x)^(5/2)/(b*x + c*x^2)^3")] + [InlineData("(A + B*x)*(d + k*x)^(1/2)/(b*x + c*x^2)^3")] + [InlineData("(A + B*x)*(d + k*x)^(7/2)/(b*x + c*x^2)^3")] + [InlineData("(2 + 3*x)^(3/2)/(b*x + c*x^2)^3")] + public void UnderTheRootWithTheQuadraticAsItsLinears(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 30000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) => e.Substitute("A", 1.1).Substitute("B", 0.6).Substitute("b", 1.3) + .Substitute("c", 0.7).Substitute("d", 0.4).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[] { -3.1, -1.7, -0.3, 0.4, 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/CubeRootsOfTwoLinearsBesideAWholePowerOfAThirdIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/CubeRootsOfTwoLinearsBesideAWholePowerOfAThirdIntegralTest.cs new file mode 100644 index 000000000..bea04824b --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/CubeRootsOfTwoLinearsBesideAWholePowerOfAThirdIntegralTest.cs @@ -0,0 +1,59 @@ +// +// 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 +{ + /// + /// Cube roots of two linears whose powers add up to a whole number beside a whole power of a + /// third linear, under t = (a + b x)^(1/3)/(c + d x)^(1/3) with the coefficients of the + /// rational function named: written as the substitution leaves them they ran past the budget. + /// Rubi's 1.1.1.3. Points on both sides of the roots' zeros, where the cube roots are complex. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class CubeRootsOfTwoLinearsBesideAWholePowerOfAThirdIntegralTest + { + [Theory] + [InlineData("(a + b*x)^(1/3)*(c + d*x)^(2/3)/(k + f*x)")] + [InlineData("(a + b*x)^(1/3)*(c + d*x)^(2/3)/(k + f*x)^2")] + [InlineData("(a + b*x)^(1/3)*(c + d*x)^(2/3)/(k + f*x)^4")] + [InlineData("(a + b*x)^(1/3)/((c + d*x)^(1/3)*(k + f*x)^3)")] + [InlineData("1/((a + b*x)^(1/3)*(c + d*x)^(2/3)*(k + f*x)^2)")] + [InlineData("1/((a + b*x)^(1/3)*(c + d*x)^(2/3)*(k + f*x)^3)")] + [InlineData("(a + b*x)^(4/3)/((c + d*x)^(4/3)*(k + f*x)^2)")] + [InlineData("(1 + x)^(1/3)*(2 - x)^(2/3)/(k + f*x)^2")] + [InlineData("(a + b*x)^(1/3)*(c + d*x)^(5/3)/(k + f*x)^3")] + public void UnderTheQuotientOfTheRoots(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 8000, $"{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).Substitute("f", 1.7); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -3.1, -1.7, -0.3, 0.4, 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}"); + } + } +}