diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 44d906a18..f2ef6dbd3 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -532,6 +532,27 @@ written in powers of `u = a + b f`, with coefficients `c_j = sum_i p_i C(i, j) ( | `"(a+b*sin(e+f*x))^2*(A+B*sin(e+f*x))/(c+d*sin(e+f*x))^2".ToEntity().Integrate("x")` | `integral(...)` | 3,337 characters | | `"cos(c+d*x)^4*(A+C*cos(c+d*x)^2)/(a+b*cos(c+d*x))^4".ToEntity().Integrate("x")` | `integral(...)` | 10,859 characters | +### An even power of a cosine over a power of `a + b sin` is reduced in the power of the binomial + +**Answers where there were none.** `cos(c + d x)^8/(a + b sin(c + d x))^8` ran out of time, and so did +`cos^6`, `sec^2` and `sec^4` over the same eighth power and `sec^4/(a + b sin)^3`, Rubi's 4.1.1.2: with symbols for `a` +and `b` the half-angle's partial fractions did not finish. With `u = a + b sin(y)`, `b^2 cos(y)^2 = (b^2 - a^2) + 2 a u - u^2` +lowers a power of the cosine above the bar, and read for `a^2 - b^2` it lowers one below it; what is left is the integral +`I_n` of `1/u^n`, and `(n - 1)(a^2 - b^2) I_n = b cos(y) u^(1 - n) + (2 n - 3) a I_(n - 1) - (n - 2) I_(n - 2)` takes it +down to `I_1`. A sine over a binomial in the cosine is the same with `-b sin(y)`. The coefficients are summed exactly, so +the answer holds `I_1` once. For `a^2 != b^2`, the generic case, as Rubi answers it: `a + a sin` declines. `1/(a + b sin)^n` +and a polynomial in the sine over it are answered this way too, more shortly than by the half-angle +([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"cos(c+d*x)^8/(a+b*sin(c+d*x))^8".ToEntity().Integrate("x")` | `integral(...)` | 4,788 characters | +| `"sec(c+d*x)^4/(a+b*sin(c+d*x))^8".ToEntity().Integrate("x")` | `integral(...)` | 6,910 characters | +| `"sec(c+d*x)^4/(a+b*sin(c+d*x))^3".ToEntity().Integrate("x")` | `integral(...)` | 3,030 characters | +| `"cos(c+d*x)^4/(a+b*sin(c+d*x))^8".ToEntity().Integrate("x")` | `integral(...)` | 4,092 characters | +| `"sin(x)^4/(a+b*cos(x))".ToEntity().Integrate("x")` | `integral(...)` | 1,309 characters | +| `"1/(a+b*cos(c+d*x))^4".ToEntity().Integrate("x")` | `integral(...)` | 1,851 characters | + ### A difference of a negated term is printed as a sum **A printed form.** A solution and a difference were printed with the sign of a negative factor left inside: @@ -658,6 +679,18 @@ Only where the product `T B` is in the integrand; without it the rules after thi | `"x^5*(a+b*acoth(c*x))*(d+k*ln(1-c^2*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 717 characters | | `"(a+b*acoth(c*x))*(d+k*ln(1-c^2*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 354 characters | +### An odd power of a tangent beside a root of `a + b sin^2` is integrated under `v = sin^2` + +**Answers where there were none.** `tan(e + f x)^5/(a + b sin(e + f x)^2)^(5/2)`, Rubi's 4.1.7, ran out of time: the +substitution in the secant wrote the root as a power of `a + b - b/u^2` and spent eight seconds writing it apart. With +`v = sin(y)^2`, `tan(y)^(2 k + 1) dy = v^k/(1 - v)^(k + 1) dv/2`, and in `v` the root is of a linear. Only where a root of a +sum with a constant in it holds the sine's square ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"tan(k+f*x)^5/(a+b*sin(k+f*x)^2)^(5/2)".ToEntity().Integrate("x")` | `integral(...)` | 4,109 characters | +| `"tan(k+f*x)^5/sqrt(a-a*sin(k+f*x)^2)".ToEntity().Integrate("x")` | `integral(...)` | 190 characters | + ### A power of x times an inverse tangent and the logarithm of its quadratic is integrated by parts **Answers where there were none.** `x^4 (a + b atan(c x)) (d + e ln(1 + c^2 x^2))` ran out of time, and so did three more @@ -834,6 +867,17 @@ over a power of a linear. That is answered with `Si` and `Ci` | `"cos(c+d*x)^3/((pe+f*x)*(a+a*sin(c+d*x)))".ToEntity().Integrate("x")` | `integral(...)` | 246 characters | | `"cos(c+d*x)^3/((pe+f*x)^2*(a+a*sin(c+d*x)))".ToEntity().Integrate("x")` | `integral(...)` | 790 characters | +### From the twelfth power, a power of a tangent over a power of `a + a sec` is cancelled + +**Answers where there were none.** `tan(c + d x)^12/(a + a sec(c + d x))^3`, Rubi's 4.5.1.4, ran out of time. `tan^2` is +`(sec - 1)(sec + 1)`, so `tan(y)^m/(k (1 + sec(y)))^n` is `tan(y)^(m - 2 n) (sec(y) - 1)^n/k^n`, and over `k (1 - sec(y))` +the same with `sec(y) + 1` over `(-k)^n`. From the twelfth power only, where the rules after it stop answering in time: +below it they answer more shortly ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"tan(c+d*x)^12/(a+a*sec(c+d*x))^3".ToEntity().Integrate("x")` | `integral(...)` | 1,140 characters | + ### Every linear beside a block is split off, a quadratic with rational roots among them **Answers where there were none.** `sec(c + d x)/(a + b sin(c + d x)^3)` was declined, and so were @@ -1006,6 +1050,20 @@ split. It is written `b (s + r)(s^2 - r s + r^2)` now, with `r = (a/b)^(1/3)`, w | `"sec(c + d*x)^2/(a + b*sin(c + d*x)^3)".ToEntity().Integrate("x")` | `integral(...)` | 7,010 characters | | `"sec(c + d*x)^4/(a + b*sin(c + d*x)^3)".ToEntity().Integrate("x")` | `integral(...)` | 11,646 characters | +### A root of `1 - cos` beside another root is integrated under the cosine + +**Answers where there were none.** `sqrt(1 - cos(x))/sqrt(a - cos(x))`, Rubi's 4.2.2.1, was declined and +`sqrt((1 - cos(x))/(a - cos(x)))` ran out of time: the half angle at which `1 - cos` is a square reads only roots of +that kind. With `g` the other function, `sqrt(k (1 + s f)/D)` is `K g/(sqrt(D) sqrt(1 - s f))` where `K` has the square +`k` and so is constant between the zeros of `g`, and under `u = f(y)` what is left is +`1/(sqrt(1 + u) sqrt(a - u))`. For a number `k` and another root with a symbol in it only: the rules after it answer +the rest, and more shortly ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"sqrt(1-cos(x))/sqrt(a-cos(x))".ToEntity().Integrate("x")` | `integral(...)` | 98 characters | +| `"sqrt((1-cos(x))/(a-cos(x)))".ToEntity().Integrate("x")` | `integral(...)` | 138 characters | + ### A linear factor common to a sum below the bar is taken out of it **Answers where there were none.** `1/(a x + b x/log(c x^n))` was left unevaluated. With x taken out of the sum below the @@ -1424,6 +1482,18 @@ integral of the function ([#718](https://github.com/asc-community/AngouriMath/is | `"atan(k*x/sqrt(-a*k^2/b-k^2*x^2))^m/sqrt(a+b*x^2)".ToEntity().Integrate("x")` | `integral(...)` | 320 characters | | `"1/(atan(k*x/sqrt(-a*k^2/b-k^2*x^2))^3*sqrt(a+b*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 314 characters | +### Roots that are a constant times a power of x are taken out as that constant + +**Answers where there were none.** `sin(x) sqrt(b - a/x^2)/sqrt(a - b x^2)`, Rubi's 4.7.7, was declined. The roots are +`K/x` with `K = x sqrt(b - a/x^2)/sqrt(a - b x^2)`, whose square is `-1`: constant between the zeros of the roots, where +it changes sign. With `r` the degree of the roots both at infinity and at zero, `K = A/x^r` is taken as constant where its +derivative is zero at sampled points, and the integral is `K` times that of `x^r` times the rest, here `K Si(x)`. Last +but two, where everything reading the roots as written has declined ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"sin(x)*sqrt(b-a/x^2)/sqrt(a-b*x^2)".ToEntity().Integrate("x")` | `integral(...)` | 107 characters | + ### Half an odd power of an inverse hyperbolic cosine of `c + d x^2` is integrated in error functions **Answers where there were none.** `(a + b acosh(1 + d x^2))^(1/2)` was declined, and so were thirteen more of Rubi's @@ -1504,6 +1574,21 @@ answer these integrands; the split now ends at the first term declined. More slo | `"1/(a*coth(x)+b*csch(x))^2".ToEntity().Integrate("x")` | `integral(...)` | 1,831 characters | | `"1/(1-cosh(x)^8)".ToEntity().Integrate("x")` | `integral(...)` | 1,232 characters | +### A power of a secant over a power of `a cos + b sin` is written in powers of `a + b tan` + +**Answers where there were none.** `sec(c + d x)^4/(a cos(c + d x) + b sin(c + d x))^3` and `sec^3/(a cos + b sin)^4`, Rubi's +4.7.2, ran out of time. Where the two powers together are even the quotient is a rational function of the tangent beside +`sec^2`; odd, it went to the half-angle. It is `sec(y) (1 + tan(y)^2)^k/(a + b tan(y))^n` with `2 k + 1 = m + n`, and in +powers of `w = a + b tan(y)` a term below the bar, `sec(y)/w^i`, is `cos(y)^(i - 1)/(a cos(y) + b sin(y))^i`, and one above +it a polynomial in the tangent beside the secant. From `m + n = 7`: below it the half-angle answers in time +([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"sec(c+d*x)^4/(a*cos(c+d*x)+b*sin(c+d*x))^3".ToEntity().Integrate("x")` | `integral(...)` | 2,175 characters | +| `"sec(c+d*x)^3/(a*cos(c+d*x)+b*sin(c+d*x))^4".ToEntity().Integrate("x")` | `integral(...)` | 3,401 characters | +| `"sec(c+d*x)^6/(a*cos(c+d*x)+b*sin(c+d*x))".ToEntity().Integrate("x")` | `integral(...)` | 1,151 characters | + ### A whole power of an inverse function of `c + d x^2` is integrated by its recurrence **Answers where there were none.** `(a + b acos(1 + d x^2))^4` was declined, and so were twenty-five more of Rubi's @@ -1629,6 +1714,18 @@ twice ([#1877](https://github.com/asc-community/AngouriMath/issues/1877), | `"sin(a + b*x)/((d*e - c*f)*x^3 + f)".ToEntity().Integrate("x")` | `integral(...)` | 1,416 characters | | `"sin(x)/(a*x^2 + 2*a*x + a)".ToEntity().Integrate("x")` | `integral(...)` | 246 characters | +### An odd power of a sine over a power of `a + b tan^2` is integrated under the cosine with `a - b` named + +**Answers where there were none.** `sin(e + f x)^3/(a + b tan(e + f x)^2)^3`, Rubi's 4.3.7, ran out of time in the +substitution search before Bioche's substitution was reached, and `sin^5` over the same cube took fifteen seconds. Under +`u = cos(y)`, `a + b tan(y)^2 = ((a - b) u^2 + b)/u^2`, and the rational function in `u` is asked with `a - b` named by +a symbol of its own: written in, it took seven seconds, named a tenth of one ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"sin(k+f*x)^3/(a+b*tan(k+f*x)^2)^3".ToEntity().Integrate("x")` | `integral(...)` | 2,667 characters | +| `"sin(k+f*x)^5/(a+b*tan(k+f*x)^2)^3".ToEntity().Integrate("x")` | `integral(...)` | 6,183 characters | + ### A trigonometric function of `a + b/x` beside a rational function is integrated in `1/x` **Answers where there were none.** `sin(c + d/x)/(e + f x)` was declined: the substitution search reads a candidate's @@ -3842,7 +3939,7 @@ of the six functions read inside a power as it is | `"(c*sin(a+b*x)^n)^(1/n)".ToEntity().Integrate("x")` | `integral(...)` | 77 characters | | `"sin(c+d*x)^6/(a+b*sin(c+d*x)^2)^2".ToEntity().Integrate("x")` | `integral(...)` | 2,169 characters | | `"cot(c+d*x)^2/(a+b*sec(c+d*x)^2)^2".ToEntity().Integrate("x")` | `integral(...)` | 1,963 characters | -| `"sqrt(1-cos(x))/sqrt(a-cos(x))".ToEntity().Integrate("x")` | `integral(...)` | 70 characters | +| `"sqrt(1-cos(x))/sqrt(a-cos(x))".ToEntity().Integrate("x")` | `integral(...)` | 98 characters | ### A product of powers of the trigonometric functions, some not whole, is integrated as sine and cosine powers diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index c309bc8d9..741f4eff6 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -33453,6 +33453,834 @@ when two.Equals(EInteger.FromInt32(2)) && e1 == MathS.e && e2 == MathS.e return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; } + /// + /// N/((a^2 - b^2)^d b^e) with N a polynomial over the rationals in a, + /// b and a third placeholder: the coefficients the reduction formulas in a power of + /// a + b sin(y) leave, kept exact, so that summing the many paths that reach one + /// term adds numerators instead of growing an expression. + /// + private readonly record struct OverTheDifferenceOfSquares(MultivariatePolynomial Numerator, int Difference, int B) + { + /// a^2 - b^2 in the placeholders, a the first and b the second. + internal static MultivariatePolynomial DifferenceOfSquares + => MultivariatePolynomial.Monomial(3, 0).Multiply(MultivariatePolynomial.Monomial(3, 0))! + .Subtract(MultivariatePolynomial.Monomial(3, 1).Multiply(MultivariatePolynomial.Monomial(3, 1))!); + + internal static OverTheDifferenceOfSquares Of(MultivariatePolynomial numerator) => new(numerator, 0, 0); + + /// Over one bar, or null where the numerator grows past what the polynomials hold. + internal OverTheDifferenceOfSquares? Plus(OverTheDifferenceOfSquares other) + { + var difference = System.Math.Max(Difference, other.Difference); + var b = System.Math.Max(B, other.B); + return Raised(difference - Difference, b - B) is { } left && other.Raised(difference - other.Difference, b - other.B) is { } right + ? new OverTheDifferenceOfSquares(left.Add(right), difference, b) : null; + } + + internal OverTheDifferenceOfSquares? Times(MultivariatePolynomial factor) + => Numerator.Multiply(factor) is { } product ? this with { Numerator = product } : null; + + private MultivariatePolynomial? Raised(int differences, int bs) + { + MultivariatePolynomial? raised = Numerator; + for (var i = 0; i < differences && raised is not null; i++) + raised = raised.Multiply(DifferenceOfSquares); + return raised is null || bs == 0 ? raised : raised.ShiftedBy(1, bs); + } + + /// In and , what divides the numerator divided out first. + internal Entity ToEntity(IReadOnlyList placeholders, Entity a, Entity b) + { + var (numerator, difference, power) = (Numerator, Difference, B); + while (difference > 0 && numerator.DivideExact(DifferenceOfSquares) is { } lower) + (numerator, difference) = (lower, difference - 1); + while (power > 0 && numerator.DivideExact(MultivariatePolynomial.Monomial(3, 1)) is { } lower) + (numerator, power) = (lower, power - 1); + Entity value = numerator.ToEntity(placeholders).Substitute(placeholders[0], a).Substitute(placeholders[1], b); + if (difference > 0) + value /= difference == 1 ? a * a - b * b : MathS.Pow(a * a - b * b, Number.Integer.Create(difference)); + if (power > 0) + value /= power == 1 ? b : MathS.Pow(b, Number.Integer.Create(power)); + return value; + } + } + + /// + /// An even power of a cosine over a whole power of a binomial in the sine of the same + /// linear, cos(y)^p/(a + b sin(y))^n, or of a sine over one in the cosine, with + /// p of either sign -- a power of the secant is one below the bar -- by the + /// reduction formulas in the power of the binomial. With u = a + b sin(y), + /// b^2 cos(y)^2 = (b^2 - a^2) + 2 a u - u^2 lowers a power of the cosine above the + /// bar, and read for a^2 - b^2, (a^2 - b^2) = 2 a u - u^2 - b^2 cos(y)^2 + /// lowers one below it; what is left is I_n, the integral of 1/u^n, and + /// (n - 1)(a^2 - b^2) I_n = b cos(y) u^(1 - n) + (2 n - 3) a I_(n - 1) - (n - 2) I_(n - 2) + /// takes it down to I_1, the one integral asked of the rest of the integrator + /// besides the powers of the sine and the secant the reductions end on. For the cosine, + /// b sin(y) is -b sin(y) in the last. + /// + /// + /// The rules after it took the whole quotient through the half-angle, where the + /// coefficients are symbols: Rubi's 4.1.1.2 cos(c + d x)^8/(a + b sin(c + d x))^8, + /// cos^6/(a + b sin)^8, sec^2/(a + b sin)^8, sec^4/(a + b sin)^8 + /// and sec^4/(a + b sin)^3 ran out of time, and cos^4/(a + b sin)^8 and + /// cos^6/(a + b sin), answered, were 40 thousand characters each. For + /// a^2 - b^2 != 0, the generic case Rubi answers: a + a sin(y) is + /// another shape, and declines. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveAnEvenPowerOfACosineOverAPowerOfABinomialInTheSineByReduction(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + static IEnumerable<(Entity Factor, bool Below)> Factors(Entity e, bool below) => e switch + { + Mulf(var l, var r) => Factors(l, below).Concat(Factors(r, below)), + Divf(var l, var r) => Factors(l, below).Concat(Factors(r, !below)), + _ => new[] { (e, below) } + }; + Entity constant = Number.Integer.One; + Entity? binomial = null; + var n = 0; + var trigonometric = new List<(Entity Function, int Power)>(); + foreach (var (factor, below) in Factors(expr, false)) + { + if (!factor.ContainsNode(x)) + { + constant = below ? constant / factor : constant * factor; + continue; + } + var (@base, power) = factor is Powf(var b0, Number.Integer { EInteger: var e0 }) && e0.CanFitInInt32() + ? (b0, e0.ToInt32Unchecked()) : (factor, 1); + if (below) + power = -power; + if (@base is Sinf or Cosf or Secantf or Cosecantf) + trigonometric.Add((@base, power)); + else if (@base is Sumf or Minusf && (binomial is null || binomial == @base)) + { + binomial = @base; + n -= power; + } + else + return null; + } + if (binomial is null || n < 1 || n > 12) + return null; + var functions = binomial.Nodes.Where(node => node is Sinf or Cosf && node.ContainsNode(x)).Distinct().ToList(); + if (functions.Count != 1) + return null; + var function = functions[0]; + var sine = function is Sinf; + var argument = function is Sinf(var ys) ? ys : ((Cosf)function).Argument; + // The power of the co-function, its reciprocal counted below the bar. + Entity cofunction = sine ? MathS.Cos(argument) : MathS.Sin(argument); + Entity reciprocal = sine ? MathS.Sec(argument) : MathS.Cosec(argument); + var p = 0; + foreach (var (written, power) in trigonometric) + { + if (written == cofunction) + p += power; + else if (written == reciprocal) + p -= power; + else + return null; + } + // An odd power is a rational function of the binomial's own function, and 1/u is the base case. + if (p % 2 != 0 || p == 0 && n < 2 || System.Math.Abs(p) > 12) + return null; + if (!TreeAnalyzer.TryGetPolyLinear(argument, x, out var slope, out _) || Functions.PartialFractions.IsZeroAsAValue(slope)) + return null; + var s = Variable.CreateUnique(expr, "s_red"); + if (!TreeAnalyzer.TryGetPolyLinear(binomial.Substitute(function, s), s, out var b, out var a) + || a.ContainsNode(x) || b.ContainsNode(x) || a.ContainsNode(s) || b.ContainsNode(s) + || Functions.PartialFractions.IsZeroAsAValue(a) || Functions.PartialFractions.IsZeroAsAValue(b) + || VanishesIdentically(a * a - b * b)) + return null; + + // The placeholders a, b and u, in that order. + var placeholders = new[] { Variable.CreateUnique(expr, "a_red"), Variable.CreateUnique(expr, "b_red"), Variable.CreateUnique(expr, "u_red") }; + var pa = MultivariatePolynomial.Monomial(3, 0); + var pb = MultivariatePolynomial.Monomial(3, 1); + var pu = MultivariatePolynomial.Monomial(3, 2); + MultivariatePolynomial Whole(int value) => MultivariatePolynomial.Constant(3, ERational.FromInt32(value)); + // The terms: (0, j) for cos(y) u^(-j), (1, m) for I_m, and (2, q, i) for the integral of + // cos(y)^q sin(y)^i, q zero or negative, asked of the integrator. + var combination = new Dictionary<(int Kind, int First, int Second), OverTheDifferenceOfSquares>(); + bool Add(Dictionary<(int, int, int), OverTheDifferenceOfSquares> into, (int, int, int) key, OverTheDifferenceOfSquares? value) + { + if (value is not { } added) + return false; + if (!into.TryGetValue(key, out var existing)) + { + into[key] = added; + return true; + } + if (existing.Plus(added) is not { } sum) + return false; + into[key] = sum; + return true; + } + // cos(y)^q u^k, k >= 0, in powers of the sine: u^k = sum C(k, i) a^(k - i) b^i sin(y)^i. + bool AddPowers(Dictionary<(int, int, int), OverTheDifferenceOfSquares> into, int q, int k, OverTheDifferenceOfSquares factor) + { + var choose = EInteger.One; + for (var i = 0; i <= k; i++) + { + if (i > 0) + choose = choose * (k - i + 1) / i; + var monomial = pa.Power(k - i)?.Multiply(pb.Power(i)!)?.ScaleBy(ERational.FromEInteger(choose)); + if (monomial is null || !Add(into, (2, q, i), factor.Times(monomial))) + return false; + } + return true; + } + var reduced = new Dictionary<(int, int), Dictionary<(int, int, int), OverTheDifferenceOfSquares>?>(); + // cos(y)^q/u^m, q even and at most zero, in I_m and the powers. + Dictionary<(int, int, int), OverTheDifferenceOfSquares>? Reduced(int q, int m) + { + if (reduced.TryGetValue((q, m), out var known)) + return known; + Dictionary<(int, int, int), OverTheDifferenceOfSquares>? terms = new(); + var one = OverTheDifferenceOfSquares.Of(Whole(1)); + if (m <= 0) + terms = AddPowers(terms, q, -m, one) ? terms : null; + else if (q == 0) + terms[(1, m, 0)] = one; + else + { + // (a^2 - b^2)/(cos^(-q) u^m) = 2 a/(cos^(-q) u^(m - 1)) - 1/(cos^(-q) u^(m - 2)) - b^2/(cos^(-q - 2) u^m). + foreach (var (term, factor) in new[] { ((q, m - 1), pa.ScaleBy(ERational.FromInt32(2))), ((q, m - 2), Whole(-1)), ((q + 2, m), pb.Multiply(pb)!.ScaleBy(ERational.FromInt32(-1))) }) + { + if (terms is null || Reduced(term.Item1, term.Item2) is not { } below) + { + terms = null; + break; + } + foreach (var entry in below) + if (!Add(terms, entry.Key, entry.Value.Times(factor) is { } scaled ? scaled with { Difference = scaled.Difference + 1 } : null)) + { + terms = null; + break; + } + } + } + reduced[(q, m)] = terms; + return terms; + } + if (p > 0) + { + // cos(y)^p = ((b^2 - a^2) + 2 a u - u^2)^(p/2)/b^p, in powers of u over u^n. + var square = pb.Multiply(pb)!.Subtract(pa.Multiply(pa)!).Add(pa.Multiply(pu)!.ScaleBy(ERational.FromInt32(2))).Subtract(pu.Multiply(pu)!); + if (square.Power(p / 2) is not { } expanded) + return null; + // .Key and .Value: KeyValuePair has no Deconstruct on netstandard2.0. + foreach (var power in expanded.CoefficientsIn(2)) + { + var j = power.Key; + var factor = OverTheDifferenceOfSquares.Of(power.Value) with { B = p }; + if (n - j >= 1 ? !Add(combination, (1, n - j, 0), factor) : !AddPowers(combination, 0, j - n, factor)) + return null; + } + } + else if (Reduced(p, n) is { } reducedTerms) + combination = reducedTerms; + else + return null; + // I_m, m >= 2, in cos(y) u^(-j) and I_1: (m - 1)(a^2 - b^2) I_m = sigma b cos(y) u^(1 - m) + (2 m - 3) a I_(m - 1) - (m - 2) I_(m - 2). + var sigma = sine ? 1 : -1; + var recurrences = new Dictionary?>(); + Dictionary<(int, int, int), OverTheDifferenceOfSquares>? Recurrence(int m) + { + if (recurrences.TryGetValue(m, out var known)) + return known; + Dictionary<(int, int, int), OverTheDifferenceOfSquares>? terms = new(); + if (m == 1) + terms[(1, 1, 0)] = OverTheDifferenceOfSquares.Of(Whole(1)); + else + { + var over = ERational.Create(EInteger.One, EInteger.FromInt32(m - 1)); + OverTheDifferenceOfSquares Lowered(MultivariatePolynomial numerator) => new(numerator.ScaleBy(over), 1, 0); + var ok = Add(terms, (0, m - 1, 0), Lowered(pb.ScaleBy(ERational.FromInt32(sigma)))); + foreach (var (lower, factor) in new[] { (m - 1, pa.ScaleBy(ERational.FromInt32(2 * m - 3))), (m - 2, Whole(-(m - 2))) }) + { + if (!ok || lower < 1 || factor.IsZero) + continue; + if (Recurrence(lower) is not { } below) + { + ok = false; + break; + } + foreach (var entry in below) + ok &= Add(terms, entry.Key, entry.Value.Times(factor.ScaleBy(over)) is { } scaled ? scaled with { Difference = scaled.Difference + 1 } : null); + } + if (!ok) + terms = null; + } + recurrences[m] = terms; + return terms; + } + var final = new Dictionary<(int, int, int), OverTheDifferenceOfSquares>(); + foreach (var entry in combination) + { + if (entry.Key.Item1 != 1) + { + if (!Add(final, entry.Key, entry.Value)) + return null; + continue; + } + if (Recurrence(entry.Key.Item2) is not { } expandedPower) + return null; + foreach (var term in expandedPower) + if (!Add(final, term.Key, term.Value.Numerator.Multiply(entry.Value.Numerator) is { } product + ? new OverTheDifferenceOfSquares(product, term.Value.Difference + entry.Value.Difference, term.Value.B + entry.Value.B) : null)) + return null; + } + // The algebraic terms simplified apart from the integrals, which a piecewise + // answer for I_1 would otherwise take into each of its branches. + Entity algebraic = Number.Integer.Zero; + Entity integrals = Number.Integer.Zero; + foreach (var entry in final.OrderBy(entry => entry.Key)) + { + if (entry.Value.Numerator.IsZero) + continue; + var coefficient = entry.Value.ToEntity(placeholders, a, b); + var (kind, first, second) = entry.Key; + if (kind == 0) + { + algebraic += coefficient * cofunction / (first == 1 ? binomial : MathS.Pow(binomial, Number.Integer.Create(first))); + continue; + } + Entity question = kind == 1 ? 1 / binomial + : (first == 0 ? Number.Integer.One : MathS.Pow(reciprocal, Number.Integer.Create(-first))) + * (second == 0 ? Number.Integer.One : second == 1 ? function : MathS.Pow(function, Number.Integer.Create(second))); + if (Integration.ComputeAsAQuestionOfItsOwn(question.InnerSimplified, x, integrateByParts: true) is not { } integrated) + return null; + integrals += coefficient * integrated; + } + var answer = constant * ((algebraic / slope).InnerSimplified + integrals); + return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; + } + + /// + /// A power of a secant over a whole power of a cos(y) + b sin(y), the two powers + /// together odd, sec(y)^m/(a cos(y) + b sin(y))^n, in powers of the binomial + /// w = a + b tan(y): the quotient is sec(y) (1 + tan(y)^2)^k/w^n with + /// 2 k + 1 = m + n, and with tan(y) = (w - a)/b the power + /// (1 + tan(y)^2)^k = ((w - a)^2 + b^2)^k/b^(2 k) = sum d_j w^j. A term below the + /// bar, sec(y)/w^i, is cos(y)^(i - 1)/(a cos(y) + b sin(y))^i, and one above + /// it a polynomial in the tangent beside the secant; each is its own question. + /// + /// + /// Where m + n is even the quotient is a rational function of the tangent beside + /// sec(y)^2, which the substitution answers; odd, it took the half-angle, and with + /// symbols for a and b Rubi's 4.7.2 sec(c + d x)^4/(a cos + b sin)^3 + /// and sec^3/(a cos + b sin)^4 ran out of time, and sec^6/(a cos + b sin), + /// answered, was 5.6 thousand characters. From m + n = 7, where the half-angle + /// stops answering in time. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveASecantPowerOverAPowerOfACosineAndASineInPowersOfABinomialInTheTangent(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + static IEnumerable<(Entity Factor, bool Below)> Factors(Entity e, bool below) => e switch + { + Mulf(var l, var r) => Factors(l, below).Concat(Factors(r, below)), + Divf(var l, var r) => Factors(l, below).Concat(Factors(r, !below)), + _ => new[] { (e, below) } + }; + Entity constant = Number.Integer.One; + Entity? combination = null; + Entity? argument = null; + var n = 0; + var m = 0; + foreach (var (factor, below) in Factors(expr, false)) + { + if (!factor.ContainsNode(x)) + { + constant = below ? constant / factor : constant * factor; + continue; + } + var (@base, power) = factor is Powf(var b0, Number.Integer { EInteger: var e0 }) && e0.CanFitInInt32() + ? (b0, e0.ToInt32Unchecked()) : (factor, 1); + if (below) + power = -power; + switch (@base) + { + case Secantf(var y) when argument is null || argument == y: + argument = y; + m += power; + break; + case Cosf(var y) when argument is null || argument == y: + argument = y; + m -= power; + break; + case Sumf or Minusf when combination is null || combination == @base: + combination = @base; + n -= power; + break; + default: + return null; + } + } + // From m + n = 7: below it the half-angle answers in time, and this, asking up to five + // questions, was no faster on `sec(y)/(a cos(y) + b sin(y))^4`. + if (combination is null || argument is null || m < 1 || n < 1 || (m + n) % 2 == 0 || m + n < 7 || m + n > 13) + return null; + if (!TreeAnalyzer.TryGetPolyLinear(argument, x, out var slope, out _) || Functions.PartialFractions.IsZeroAsAValue(slope)) + return null; + // a cos(y) + b sin(y), nothing else. + var cosine = Variable.CreateUnique(expr, "c_tan"); + var sine = Variable.CreateUnique(expr, "s_tan"); + var inBoth = combination.Substitute(MathS.Cos(argument), cosine).Substitute(MathS.Sin(argument), sine); + if (inBoth.ContainsNode(x) + || !TreeAnalyzer.TryGetPolyLinear(inBoth, cosine, out var a, out var rest) + || !TreeAnalyzer.TryGetPolyLinear(rest, sine, out var b, out var free) + || a.ContainsNode(sine) || !Functions.PartialFractions.IsZeroAsAValue(free) + || Functions.PartialFractions.IsZeroAsAValue(a) || Functions.PartialFractions.IsZeroAsAValue(b)) + return null; + var k = (m + n - 1) / 2; + // The placeholders a, b and w, in that order; every coefficient is over b^(2 k). + var placeholders = new[] { Variable.CreateUnique(expr, "a_tan"), Variable.CreateUnique(expr, "b_tan"), Variable.CreateUnique(expr, "w_tan") }; + var pa = MultivariatePolynomial.Monomial(3, 0); + var pb = MultivariatePolynomial.Monomial(3, 1); + var pw = MultivariatePolynomial.Monomial(3, 2); + var shifted = pw.Subtract(pa); + if (shifted.Multiply(shifted)?.Add(pb.Multiply(pb)!).Power(k) is not { } expanded) + return null; + // The terms: (true, i) for cos(y)^(i - 1)/(a cos(y) + b sin(y))^i, (false, i) for sec(y) tan(y)^i. + var terms = new Dictionary<(bool Below, int Power), MultivariatePolynomial>(); + void Add((bool, int) key, MultivariatePolynomial value) + => terms[key] = terms.TryGetValue(key, out var existing) ? existing.Add(value) : value; + // .Key and .Value: KeyValuePair has no Deconstruct on netstandard2.0. + foreach (var power in expanded.CoefficientsIn(2)) + { + var (j, d) = (power.Key, power.Value); + if (j < n) + { + Add((true, n - j), d); + continue; + } + // sec(y) w^(j - n), w = a + b tan(y), in powers of the tangent. + var choose = EInteger.One; + for (var i = 0; i <= j - n; i++) + { + if (i > 0) + choose = choose * (j - n - i + 1) / i; + if (pa.Power(j - n - i)?.Multiply(pb.Power(i)!)?.Multiply(d)?.ScaleBy(ERational.FromEInteger(choose)) is not { } term) + return null; + Add((false, i), term); + } + } + Entity answer = Number.Integer.Zero; + foreach (var entry in terms.OrderBy(entry => entry.Key)) + { + if (entry.Value.IsZero) + continue; + // Over b^(2 k), with the powers of b the numerator holds divided out first. + var (numerator, power) = (entry.Value, 2 * k); + while (power > 0 && numerator.DivideExact(pb) is { } lower) + (numerator, power) = (lower, power - 1); + Entity coefficient = numerator.ToEntity(placeholders).Substitute(placeholders[0], a).Substitute(placeholders[1], b); + if (power > 0) + coefficient /= power == 1 ? b : MathS.Pow(b, Number.Integer.Create(power)); + var (isBelow, i) = entry.Key; + Entity question = isBelow + ? (i == 1 ? Number.Integer.One : i == 2 ? MathS.Cos(argument) : MathS.Pow(MathS.Cos(argument), Number.Integer.Create(i - 1))) + / (i == 1 ? combination : MathS.Pow(combination, Number.Integer.Create(i))) + : MathS.Sec(argument) * (i == 0 ? Number.Integer.One : i == 1 ? MathS.Tan(argument) : MathS.Pow(MathS.Tan(argument), Number.Integer.Create(i))); + if (Integration.ComputeAsAQuestionOfItsOwn(question.InnerSimplified, x, integrateByParts: true) is not { } integrated) + return null; + answer += coefficient * integrated; + } + answer = constant * answer; + return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; + } + + /// + /// An odd power of a sine over a whole power of a + b tan(y)^2, under + /// u = cos(y): a + b tan(y)^2 = ((a - b) u^2 + b)/u^2 and + /// sin(y)^m dy = -(1 - u^2)^((m - 1)/2) du, so the integrand is the rational + /// function -(1 - u^2)^((m - 1)/2) u^(2 n)/((a - b) u^2 + b)^n, asked with + /// a - b named by a symbol of its own and written back into the answer. + /// + /// + /// Bioche's substitution, which comes after the substitution search: the search ran out + /// of time on Rubi's 4.3.7 sin(e + f x)^3/(a + b tan(e + f x)^2)^3, and took + /// fifteen seconds on sin^5/(a + b tan^2)^3. The rational function with + /// a - b written in it took seven seconds, where named it takes a tenth of one. + /// The quadratic has no linear term, so the name hides no relation between its + /// coefficients; a = b is another shape, and declines. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveAnOddPowerOfASineOverAPowerOfABinomialInTheSquaredTangentUnderTheCosine(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + static IEnumerable<(Entity Factor, bool Below)> Factors(Entity e, bool below) => e switch + { + Mulf(var l, var r) => Factors(l, below).Concat(Factors(r, below)), + Divf(var l, var r) => Factors(l, below).Concat(Factors(r, !below)), + _ => new[] { (e, below) } + }; + Entity constant = Number.Integer.One; + Entity? binomial = null; + Entity? argument = null; + var n = 0; + var m = 0; + foreach (var (factor, below) in Factors(expr, false)) + { + if (!factor.ContainsNode(x)) + { + constant = below ? constant / factor : constant * factor; + continue; + } + var (@base, power) = factor is Powf(var b0, Number.Integer { EInteger: var e0 }) && e0.CanFitInInt32() + ? (b0, e0.ToInt32Unchecked()) : (factor, 1); + if (below) + power = -power; + switch (@base) + { + case Sinf(var y) when argument is null || argument == y: + argument = y; + m += power; + break; + case Sumf or Minusf when binomial is null || binomial == @base: + binomial = @base; + n -= power; + break; + default: + return null; + } + } + if (binomial is null || argument is null || m < 1 || m % 2 == 0 || m > 15 || n < 1 || n > 8) + return null; + if (!TreeAnalyzer.TryGetPolyLinear(argument, x, out var slope, out _) || Functions.PartialFractions.IsZeroAsAValue(slope)) + return null; + var tangent = MathS.Tan(argument); + var t = Variable.CreateUnique(expr, "t_sq"); + var inT = binomial.Substitute(tangent, t); + if (inT.ContainsNode(x) || !TreeAnalyzer.TryGetPolynomial(inT, t, out var read) + || read.Count != 2 || !read.TryGetValue(EInteger.Zero, out var a) || !read.TryGetValue(EInteger.FromInt32(2), out var b) + || Functions.PartialFractions.IsZeroAsAValue(b)) + return null; + var difference = (a - b).InnerSimplified; + if (VanishesIdentically(difference)) + return null; + var named = difference is Number or Variable ? difference : Variable.CreateUnique(expr, "g_sq"); + var u = Variable.CreateUnique(expr + named, "u_sq"); + var k = (m - 1) / 2; + Entity inU = MathS.Pow(u, Number.Integer.Create(2 * n)) / (n == 1 ? named * MathS.Sqr(u) + b : MathS.Pow(named * MathS.Sqr(u) + b, Number.Integer.Create(n))); + inU = -(k == 0 ? inU : MathS.Pow(1 - MathS.Sqr(u), Number.Integer.Create(k)) * inU); + if (Integration.ComputeAsAQuestionOfItsOwn(inU, u, integrateByParts: true) is not { } integrated) + return null; + if (named is Variable name && name != difference) + integrated = integrated.Substitute(name, difference); + var answer = constant * integrated.Substitute(u, MathS.Cos(argument)) / slope; + return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; + } + + /// + /// An odd power of a tangent beside a function of the square of the sine of the same + /// linear, tan(y)^(2 k + 1) F(sin(y)^2), under v = sin(y)^2: with + /// dv = 2 sin(y) cos(y) dy, tan(y)^(2 k + 1) dy = v^k/(1 - v)^(k + 1) dv/2, + /// and v^k F(v)/(2 (1 - v)^(k + 1)) is asked of the rest of the integrator. + /// + /// + /// The substitution in the secant after it writes a + b sin(y)^2 as + /// a + b (1 - 1/u^2), and writing the power of that quotient apart took eight + /// seconds on Rubi's 4.1.7 tan(e + f x)^5/(a + b sin(e + f x)^2)^(5/2), which + /// ran out of time; in v it is a power of a linear beside a rational function, + /// answered in a tenth of a second. Only where a root holds the sine's square: a + /// rational function of it the secant's substitution answers as shortly. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveAnOddPowerOfATangentBesideAFunctionOfTheSquaredSineUnderIt(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + static IEnumerable<(Entity Factor, bool Below)> Factors(Entity e, bool below) => e switch + { + Mulf(var l, var r) => Factors(l, below).Concat(Factors(r, below)), + Divf(var l, var r) => Factors(l, below).Concat(Factors(r, !below)), + _ => new[] { (e, below) } + }; + Entity constant = Number.Integer.One; + Entity rest = Number.Integer.One; + Entity? argument = null; + var j = 0; + foreach (var (factor, below) in Factors(expr, false)) + { + if (!factor.ContainsNode(x)) + { + constant = below ? constant / factor : constant * factor; + continue; + } + var (@base, power) = factor is Powf(var b0, Number.Integer { EInteger: var e0 }) && e0.CanFitInInt32() + ? (b0, e0.ToInt32Unchecked()) : (factor, 1); + if (@base is Tanf(var y) && (argument is null || argument == y)) + { + argument = y; + j += below ? -power : power; + } + else + rest = below ? rest / factor : rest * factor; + } + if (argument is null || j < 1 || j % 2 == 0 || j > 15 || rest == Number.Integer.One) + return null; + if (!TreeAnalyzer.TryGetPolyLinear(argument, x, out var slope, out _) || Functions.PartialFractions.IsZeroAsAValue(slope)) + return null; + var v = Variable.CreateUnique(expr, "v_sq"); + var sine = MathS.Sin(argument); + var inV = rest.Replace(node => node is Powf(var s0, Number.Integer { EInteger: var e1 }) && s0 == sine && !e1.IsZero && e1.IsEven && e1.CanFitInInt32() + ? MathS.Pow(v, Number.Integer.Create(e1.ToInt32Unchecked() / 2)) : node); + // A root of a sum in v with a constant in it, which the secant's substitution writes + // apart: a rational function of the sine's square it answers as one in the secant, and as + // shortly, and `tan(y)/sqrt(a sin(y)^2)` in 110 characters, where in v it is 627. + if (inV.ContainsNode(x) || !inV.ContainsNode(v) + || !inV.Nodes.Any(node => node is Powf(var radicand, Number.Rational and not Number.Integer) && radicand.ContainsNode(v) + && TreeAnalyzer.TryGetPolynomial(radicand, v, out var inPowers) && inPowers.Count >= 2 && inPowers.ContainsKey(EInteger.Zero))) + return null; + var k = (j - 1) / 2; + Entity measure = (k == 0 ? Number.Integer.One : k == 1 ? v : MathS.Pow(v, Number.Integer.Create(k))) + / (2 * (k == 0 ? 1 - v : MathS.Pow(1 - v, Number.Integer.Create(k + 1)))); + if (Integration.ComputeAsAQuestionOfItsOwn(measure * inV, v, integrateByParts: true) is not { } integrated) + return null; + var answer = constant * integrated.Substitute(v, MathS.Sqr(sine)) / slope; + return answer.Nodes.Any(node => node == MathS.NaN) ? null : answer; + } + + /// + /// An even power of a tangent over a power of k (1 +- sec(y)), cancelled where the + /// power above is at least twice the one below: tan^2 = (sec - 1)(sec + 1), so + /// tan(y)^m/(k (1 + sec(y)))^n is tan(y)^(m - 2 n) (sec(y) - 1)^n/k^n, and + /// over k (1 - sec(y)) it is tan(y)^(m - 2 n) (sec(y) + 1)^n/(-k)^n. + /// + /// + /// Rubi's 4.5.1.4 tan(c + d x)^12/(a + a sec(c + d x))^3 ran out of time, and + /// tan(c + d x)^6 (sec(c + d x) - 1)^3/a^3 is answered in a third of a second. + /// From the twelfth power only, where the rules after it stop answering in time -- + /// tan^12/(a + a sec) and tan^14/(a + a sec)^3 ran past a minute -- since + /// below it they answer more shortly: tan^10/(a + a sec)^3 in 507 characters + /// where cancelled it is 851. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveAnEvenPowerOfATangentOverAPowerOfOnePlusASecantByCancelling(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + static IEnumerable<(Entity Factor, bool Below)> Factors(Entity e, bool below) => e switch + { + Mulf(var l, var r) => Factors(l, below).Concat(Factors(r, below)), + Divf(var l, var r) => Factors(l, below).Concat(Factors(r, !below)), + _ => new[] { (e, below) } + }; + Entity constant = Number.Integer.One; + Entity? binomial = null; + Entity? argument = null; + var n = 0; + var m = 0; + foreach (var (factor, below) in Factors(expr, false)) + { + if (!factor.ContainsNode(x)) + { + constant = below ? constant / factor : constant * factor; + continue; + } + var (@base, power) = factor is Powf(var b0, Number.Integer { EInteger: var e0 }) && e0.CanFitInInt32() + ? (b0, e0.ToInt32Unchecked()) : (factor, 1); + if (below) + power = -power; + switch (@base) + { + case Tanf(var y) when argument is null || argument == y: + argument = y; + m += power; + break; + case Sumf or Minusf when binomial is null || binomial == @base: + binomial = @base; + n -= power; + break; + default: + return null; + } + } + if (binomial is null || argument is null || n < 1 || m % 2 != 0 || m < 2 * n || m < 12 || m > 24) + return null; + // k (1 + sigma sec(y)). + var secant = MathS.Sec(argument); + var s = Variable.CreateUnique(expr, "s_tc"); + if (!TreeAnalyzer.TryGetPolyLinear(binomial.Substitute(secant, s), s, out var slope, out var k) + || slope.ContainsNode(x) || k.ContainsNode(x) || Functions.PartialFractions.IsZeroAsAValue(k)) + return null; + var sigma = (slope / k).InnerSimplified; + int sign; + if (Functions.PartialFractions.IsZeroAsAValue((sigma - 1).InnerSimplified)) sign = 1; + else if (Functions.PartialFractions.IsZeroAsAValue((sigma + 1).InnerSimplified)) sign = -1; + else return null; + var left = m - 2 * n; + var tangent = MathS.Tan(argument); + Entity cancelled = (left == 0 ? Number.Integer.One : MathS.Pow(tangent, Number.Integer.Create(left))) + * (n == 1 ? secant - sign : MathS.Pow(secant - sign, Number.Integer.Create(n))) + / (n == 1 ? sign * k : MathS.Pow(sign * k, Number.Integer.Create(n))); + return Integration.ComputeAsAQuestionOfItsOwn((constant * cancelled).InnerSimplified, x, integrateByParts: true); + } + + /// + /// A product of roots that is a constant times a whole power of the variable, though + /// written with it, beside something else: with A the roots and r the sum of + /// their exponents times the degrees of what is under them, K = A/x^r has the + /// derivative zero wherever it is real, and integral A T dx = K integral x^r T dx. + /// sqrt(b - a/x^2)/sqrt(a - b x^2) is K/x, with K^2 = -1. + /// + /// + /// r is the degree of A both at infinity and at zero, which a constant times + /// x^r needs, and K is taken as constant only where its derivative is zero at + /// the sampled points; the answer is checked. Rubi's 4.7.7 + /// sin(x) sqrt(b - a/x^2)/sqrt(a - b x^2), which is K Si(x), was declined. + /// Late, where everything reading the integrand as written has declined. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveByTakingOutRootsConstantUpToAPowerOfTheVariable(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + static IEnumerable<(Entity Factor, bool Below)> Factors(Entity e, bool below) => e switch + { + Mulf(var l, var r) => Factors(l, below).Concat(Factors(r, below)), + Divf(var l, var r) => Factors(l, below).Concat(Factors(r, !below)), + _ => new[] { (e, below) } + }; + Entity roots = Number.Integer.One; + Entity rest = Number.Integer.One; + var atInfinity = ERational.Zero; + var atZero = ERational.Zero; + foreach (var (factor, below) in Factors(expr, false)) + { + if (factor is Powf(var radicand, Number.Rational { ERational: var exponent }) && radicand.ContainsNode(x) + && exponent.ToLowestTerms() is var q && !q.Denominator.Equals(EInteger.One)) + { + // The degrees of what is under the root at infinity and at zero, negative powers read. + if (!TreeAnalyzer.TryGetPolynomial(radicand, x, out var read) || read.Count == 0 + || read.Values.Any(value => value.ContainsNode(x) || Functions.PartialFractions.IsZeroAsAValue(value))) + return null; + var signed = below ? q.Negate() : q; + atInfinity += signed * ERational.FromEInteger(read.Keys.Max()!); + atZero += signed * ERational.FromEInteger(read.Keys.Min()!); + roots = below ? roots / factor : roots * factor; + } + else + rest = below ? rest / factor : rest * factor; + } + atInfinity = atInfinity.ToLowestTerms(); + if (roots == Number.Integer.One || !rest.ContainsNode(x) || atInfinity.CompareTo(atZero) != 0 + || !atInfinity.Denominator.Equals(EInteger.One) || !atInfinity.Numerator.CanFitInInt32()) + return null; + var r = atInfinity.Numerator.ToInt32Unchecked(); + if (System.Math.Abs(r) > 6) + return null; + // As written, x^(-r) beside the roots: simplified, the roots were taken apart and put back. + var constant = r == 0 ? roots : (r == -1 ? x : MathS.Pow(x, Number.Integer.Create(-r))) * roots; + // Constant where its derivative is zero at the sampled points: its value is not one + // number across them, since `sqrt(b - a/x^2) x/sqrt(a - b x^2)` is i on one side of + // sqrt(a/b), where both roots vanish, and -i on the other. + if (!constant.ContainsNode(x) + || !Functions.PartialFractions.HoldsAtSampledPoints(constant.Differentiate(x), Number.Integer.Zero, x)) + return null; + var power = r == 0 ? rest : r == 1 ? x * rest : MathS.Pow(x, Number.Integer.Create(r)) * rest; + if (Integration.ComputeAsAQuestionOfItsOwn(power.InnerSimplified, x, integrateByParts: true) is not { } integrated) + return null; + var answer = constant * integrated; + return answer.Nodes.Any(node => node == MathS.NaN) || Functions.PartialFractions.DerivativeDiffersWhereReal(answer, expr, x) ? null : answer; + } + + /// + /// The root of k (1 +- f(y)), or of it over something in f(y), beside another + /// root of something in f(y), for f a sine or a cosine, under + /// u = f(y): with g the other function, R = sqrt(k (1 + s f)/D) is + /// K g/(sqrt(D) sqrt(1 - s f)) where K = R sqrt(D) sqrt(1 - s f)/g has the square + /// k, and so is constant between the zeros of g; and g dy is + /// du for the sine and -du for the cosine. What is asked in u is the + /// rest over sqrt(D(u)) sqrt(1 - s u). + /// + /// + /// The half angle at which 1 +- f is a square reads only roots of that kind, and + /// Rubi's 4.2.2.1 sqrt(1 - cos(x))/sqrt(a - cos(x)) was declined and + /// sqrt((1 - cos(x))/(a - cos(x))) ran out of time; under u = cos(x) they are + /// 1/(sqrt(1 + u) sqrt(a - u)), answered in a tenth of a second. Only for a number + /// k, and where another root of something in the function with a symbol in it is + /// there: the rules after it answer the rest, and more shortly. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveARootOfOnePlusASineBesideAnotherRootUnderTheSine(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAsked) + return null; + var functions = expr.Nodes.Where(node => node.ContainsNode(x) && node is Sinf or Cosf or Tanf or Cotanf or Secantf or Cosecantf).Distinct().ToList(); + if (functions.Count != 1 || functions[0] is not (Sinf or Cosf)) + return null; + var function = functions[0]; + var sine = function is Sinf; + var argument = function is Sinf(var ys) ? ys : ((Cosf)function).Argument; + if (!TreeAnalyzer.TryGetPolyLinear(argument, x, out var slope, out _) || Functions.PartialFractions.IsZeroAsAValue(slope)) + return null; + var u = Variable.CreateUnique(expr, "u_root"); + // k (1 + s f) as written, the sign s, or null. + int? Sign(Entity radicand) + { + if (!TreeAnalyzer.TryGetPolyLinear(radicand.Substitute(function, u), u, out var coefficient, out var k) + || coefficient.ContainsNode(x) || k.ContainsNode(x) || Functions.PartialFractions.IsZeroAsAValue(k)) + return null; + var ratio = (coefficient / k).InnerSimplified; + return Functions.PartialFractions.IsZeroAsAValue((ratio - 1).InnerSimplified) ? 1 + : Functions.PartialFractions.IsZeroAsAValue((ratio + 1).InnerSimplified) ? -1 : null; + } + static IEnumerable<(Entity Factor, bool Below)> Factors(Entity e, bool below) => e switch + { + Mulf(var l, var r) => Factors(l, below).Concat(Factors(r, below)), + Divf(var l, var r) => Factors(l, below).Concat(Factors(r, !below)), + _ => new[] { (e, below) } + }; + var factors = Factors(expr, false).ToList(); + for (var at = 0; at < factors.Count; at++) + { + var (root, rootBelow) = factors[at]; + if (rootBelow || root is not Powf(var radicand, Number.Rational { ERational: var half }) + || half.CompareTo(ERational.Create(1, 2)) != 0 || !radicand.ContainsNode(x)) + continue; + var (above, below) = Functions.SingleQuotient.Of(radicand); + // A number for k: `sqrt(a + a sin(y))` beside another root the rules after it answer, + // and more shortly than this, `sqrt(a + a sin(y)) sqrt(c + d sin(y))/sin(y)` in 1,071 + // characters where this wrote 1,879. + if (above.Vars.Any(v => v != x) || Sign(above) is not { } s) + continue; + // Everything else, and another root of something in the function among it or below the root. + Entity rest = Number.Integer.One; + for (var other = 0; other < factors.Count; other++) + if (other != at) + rest = factors[other].Below ? rest / factors[other].Factor : rest * factors[other].Factor; + // With a symbol in it: over numbers the half angle answers, `sqrt(1 - cos(x))/sqrt(2 - cos(x))` + // in 95 characters where this wrote over five hundred. + var anotherRoot = below.ContainsNode(x) && below.Vars.Any(v => v != x) + || rest.Nodes.Any(node => node is Powf(var r0, Number.Rational and not Number.Integer) + && r0.ContainsNode(x) && r0.Vars.Any(v => v != x) && Sign(r0) is null); + if (!anotherRoot) + continue; + Entity complement = sine ? MathS.Cos(argument) : MathS.Sin(argument); + Entity conjugate = s == 1 ? 1 - function : 1 + function; + var constant = root * (below.ContainsNode(x) ? MathS.Sqrt(below) : Number.Integer.One) * MathS.Sqrt(conjugate) / complement; + // g dy = du for the sine, -du for the cosine. + var inU = ((sine ? rest : -rest) / ((below.ContainsNode(x) ? MathS.Sqrt(below) : Number.Integer.One) * MathS.Sqrt(conjugate))) + .Substitute(function, u); + if (inU.ContainsNode(x)) + continue; + if (Integration.ComputeAsAQuestionOfItsOwn(inU.InnerSimplified, u, integrateByParts: true) is not { } integrated) + return null; + var answer = constant * integrated.Substitute(u, function) / slope; + return answer.Nodes.Any(node => node == MathS.NaN) || Functions.PartialFractions.DerivativeDiffersWhereReal(answer, expr, x) ? null : answer; + } + return null; + } + /// /// A power of x times a whole power, past the first, of a sum holding an inverse function or a /// logarithm of a power of x, x^m (a + b atanh(c x^n))^p, under t = x^n where diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index 466ed1a97..2112f2f44 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -719,6 +719,26 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // A power of x times a power of a sum holding an inverse function of x^n, under t = x^n, // before the rules that timed out on `x^3 (a + b atanh(c x^2))^2`. if ((answer = IndefiniteIntegralSolver.SolveAPowerOfXTimesAFunctionOfAnotherUnderThatPower(expr, x)) is { }) return answer; + // An even power of a cosine over a whole power of a binomial in the sine, or the other way + // round, by the reduction formulas in the power of the binomial, before the half-angle ran out + // of time on `cos(c + d x)^8/(a + b sin(c + d x))^8` with symbols for coefficients. + if ((answer = IndefiniteIntegralSolver.SolveAnEvenPowerOfACosineOverAPowerOfABinomialInTheSineByReduction(expr, x)) is { }) return answer; + // A power of a secant over one of a cos(y) + b sin(y), the two together odd, in powers of + // a + b tan(y), before the half-angle ran out of time on `sec(c + d x)^4/(a cos + b sin)^3`. + if ((answer = IndefiniteIntegralSolver.SolveASecantPowerOverAPowerOfACosineAndASineInPowersOfABinomialInTheTangent(expr, x)) is { }) return answer; + // An odd power of a sine over a power of a + b tan(y)^2 under u = cos(y), a - b named, before + // the substitution search ran out of time on `sin(e + f x)^3/(a + b tan(e + f x)^2)^3`. + if ((answer = IndefiniteIntegralSolver.SolveAnOddPowerOfASineOverAPowerOfABinomialInTheSquaredTangentUnderTheCosine(expr, x)) is { }) return answer; + // An odd power of a tangent beside a function of sin(y)^2 under v = sin(y)^2, before the + // secant's substitution took eight seconds on `tan(e + f x)^5/(a + b sin(e + f x)^2)^(5/2)`. + if ((answer = IndefiniteIntegralSolver.SolveAnOddPowerOfATangentBesideAFunctionOfTheSquaredSineUnderIt(expr, x)) is { }) return answer; + // From the twelfth power, an even power of a tangent over a power of k (1 +- sec), cancelled: + // the rules after it ran out of time on `tan(c + d x)^12/(a + a sec(c + d x))^3`. + if ((answer = IndefiniteIntegralSolver.SolveAnEvenPowerOfATangentOverAPowerOfOnePlusASecantByCancelling(expr, x)) is { }) return answer; + // A root of k (1 +- sin) or of k (1 +- cos), or of it over something in the function, beside + // another root, under u = sin or u = cos: the half angle reads only roots of the first kind, + // and `sqrt((1 - cos(x))/(a - cos(x)))` ran out of time. + if ((answer = IndefiniteIntegralSolver.SolveARootOfOnePlusASineBesideAnotherRootUnderTheSine(expr, x)) is { }) return answer; // A polynomial in a sine or a cosine over a whole power of a binomial in it, in powers of the // binomial, before the half-angle wrote 30 million characters for `(A + B cos + C cos^2)/(a + b cos)^4`. if ((answer = IndefiniteIntegralSolver.SolveAPolynomialInASineOverAPowerOfABinomialInItInPowersOfTheBinomial(expr, x)) is { }) return answer; @@ -1612,6 +1632,10 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // 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; + // Roots that are a constant times a power of the variable, though written with it, taken + // out as that constant: `sin(x) sqrt(b - a/x^2)/sqrt(a - b x^2)` is K sin(x)/x. As late, + // since it fires only where everything reading the roots as written has declined. + if ((answer = IndefiniteIntegralSolver.SolveByTakingOutRootsConstantUpToAPowerOfTheVariable(expr, x)) 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/ARootOfOnePlusASineBesideAnotherRootUnderTheSineIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ARootOfOnePlusASineBesideAnotherRootUnderTheSineIntegralTest.cs new file mode 100644 index 000000000..97c609572 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ARootOfOnePlusASineBesideAnotherRootUnderTheSineIntegralTest.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 +{ + /// + /// The root of k (1 +- f(y)), or of it over something in f(y), beside another root of + /// something in the function, under u = f(y), the factor constant between the zeros of the + /// other function in front. Rubi's 4.2.2.1 sqrt((1 - cos(x))/(a - cos(x))) ran out of time. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ARootOfOnePlusASineBesideAnotherRootUnderTheSineIntegralTest + { + [Theory] + [InlineData("sqrt(1 - cos(x))/sqrt(a - cos(x))", 1.3, 0.7)] + [InlineData("sqrt((1 - cos(x))/(a - cos(x)))", 2.9, 0.7)] + [InlineData("sqrt(1 - sin(c + d*x))/sqrt(a + b*sin(c + d*x))", 1.3, 0.7)] + [InlineData("sqrt(1 + cos(x))/sqrt(a + b*cos(x))", 1.3, 0.7)] + public void UnderTheFunction(string integrand, double a, double b) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 3000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) => e.Substitute("a", a).Substitute("b", b).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("f", 0.6).Substitute("k", 0.2); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 2.9, 4.1, 5.3 }) + { + 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-8 * 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/ASecantPowerOverAPowerOfACosineAndASineInPowersOfABinomialInTheTangentIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ASecantPowerOverAPowerOfACosineAndASineInPowersOfABinomialInTheTangentIntegralTest.cs new file mode 100644 index 000000000..eea2d68a1 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ASecantPowerOverAPowerOfACosineAndASineInPowersOfABinomialInTheTangentIntegralTest.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 power of a secant over a whole power of a cos(y) + b sin(y), the two powers together + /// odd and at least seven, in powers of a + b tan(y). The half-angle ran out of time on + /// Rubi's 4.7.2 sec(c + d x)^4/(a cos(c + d x) + b sin(c + d x))^3. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ASecantPowerOverAPowerOfACosineAndASineInPowersOfABinomialInTheTangentIntegralTest + { + [Theory] + [InlineData("sec(c + d*x)^4/(a*cos(c + d*x) + b*sin(c + d*x))^3", 1.3, 0.7)] + [InlineData("sec(c + d*x)^3/(a*cos(c + d*x) + b*sin(c + d*x))^4", 0.6, -1.9)] + [InlineData("sec(c + d*x)^6/(a*cos(c + d*x) + b*sin(c + d*x))", 1.3, 0.7)] + [InlineData("sec(c + d*x)^2/(a*cos(c + d*x) + b*sin(c + d*x))^5", -0.6, 1.9)] + public void InPowersOfTheBinomialInTheTangent(string integrand, double a, double b) + { + 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", a).Substitute("b", b).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("f", 0.6).Substitute("k", 0.2); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 2.9, 4.1, 5.3 }) + { + 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-8 * 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/AnEvenPowerOfACosineOverAPowerOfABinomialInTheSineByReductionIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnEvenPowerOfACosineOverAPowerOfABinomialInTheSineByReductionIntegralTest.cs new file mode 100644 index 000000000..78b3879cc --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnEvenPowerOfACosineOverAPowerOfABinomialInTheSineByReductionIntegralTest.cs @@ -0,0 +1,57 @@ +// +// 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 even power of a cosine over a whole power of a + b sin(y), or of a sine over one of + /// a + b cos(y), by the reduction formulas in the power of the binomial down to + /// 1/(a + b sin(y)). The half-angle ran out of time on Rubi's 4.1.1.2 + /// cos(c + d x)^8/(a + b sin(c + d x))^8 and sec^4/(a + b sin)^8. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnEvenPowerOfACosineOverAPowerOfABinomialInTheSineByReductionIntegralTest + { + [Theory] + [InlineData("cos(c + d*x)^8/(a + b*sin(c + d*x))^8", 1.3, 0.7)] + [InlineData("cos(c + d*x)^6/(a + b*sin(c + d*x))^8", 0.6, 1.9)] + [InlineData("sec(c + d*x)^2/(a + b*sin(c + d*x))^8", 1.3, 0.7)] + [InlineData("sec(c + d*x)^4/(a + b*sin(c + d*x))^8", 0.6, 1.9)] + [InlineData("sec(c + d*x)^4/(a + b*sin(c + d*x))^3", 1.3, -0.7)] + [InlineData("sin(c + d*x)^4/(a + b*cos(c + d*x))^5", 1.3, 0.7)] + [InlineData("csc(c + d*x)^2/(a + b*cos(c + d*x))^3", 0.6, 1.9)] + [InlineData("1/(a + b*cos(c + d*x))^6", 1.3, 0.7)] + public void ByTheReductionFormulas(string integrand, double a, double b) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 12000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) => e.Substitute("a", a).Substitute("b", b).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("f", 0.6).Substitute("k", 0.2); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 2.9, 4.1, 5.3 }) + { + 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-8 * 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/AnEvenPowerOfATangentOverAPowerOfOnePlusASecantByCancellingIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnEvenPowerOfATangentOverAPowerOfOnePlusASecantByCancellingIntegralTest.cs new file mode 100644 index 000000000..97d5f6ece --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnEvenPowerOfATangentOverAPowerOfOnePlusASecantByCancellingIntegralTest.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 +{ + /// + /// From the twelfth power, an even power of a tangent over a power of k (1 +- sec(y)), + /// cancelled by tan^2 = (sec - 1)(sec + 1). The rules after it ran out of time on Rubi's + /// 4.5.1.4 tan(c + d x)^12/(a + a sec(c + d x))^3. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnEvenPowerOfATangentOverAPowerOfOnePlusASecantByCancellingIntegralTest + { + [Theory] + [InlineData("tan(c + d*x)^12/(a + a*sec(c + d*x))^3", 1.3, 0.7)] + [InlineData("tan(c + d*x)^14/(a + a*sec(c + d*x))", 1.3, 0.7)] + [InlineData("tan(c + d*x)^12/(a - a*sec(c + d*x))^3", 1.3, 0.7)] + [InlineData("tan(c + d*x)^14/(a + a*sec(c + d*x))^3", 0.6, 0.7)] + public void Cancelled(string integrand, double a, double b) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 6000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) => e.Substitute("a", a).Substitute("b", b).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("f", 0.6).Substitute("k", 0.2); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 2.9, 4.1, 5.3 }) + { + 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-8 * 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/AnOddPowerOfASineOverAPowerOfABinomialInTheSquaredTangentUnderTheCosineIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnOddPowerOfASineOverAPowerOfABinomialInTheSquaredTangentUnderTheCosineIntegralTest.cs new file mode 100644 index 000000000..5e5599059 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnOddPowerOfASineOverAPowerOfABinomialInTheSquaredTangentUnderTheCosineIntegralTest.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 +{ + /// + /// An odd power of a sine over a whole power of a + b tan(y)^2, under u = cos(y) with + /// a - b named. The substitution search ran out of time on Rubi's 4.3.7 + /// sin(e + f x)^3/(a + b tan(e + f x)^2)^3. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnOddPowerOfASineOverAPowerOfABinomialInTheSquaredTangentUnderTheCosineIntegralTest + { + [Theory] + [InlineData("sin(k + f*x)^3/(a + b*tan(k + f*x)^2)^3", 1.3, 0.7)] + [InlineData("sin(k + f*x)^5/(a + b*tan(k + f*x)^2)^3", 0.6, 1.9)] + [InlineData("sin(k + f*x)/(a + b*tan(k + f*x)^2)^2", 1.3, 0.7)] + [InlineData("sin(k + f*x)^3/(a + b*tan(k + f*x)^2)", 0.6, 1.9)] + public void UnderTheCosine(string integrand, double a, double b) + { + 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", a).Substitute("b", b).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("f", 0.6).Substitute("k", 0.2); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 2.9, 4.1, 5.3 }) + { + 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-8 * 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/AnOddPowerOfATangentBesideAFunctionOfTheSquaredSineUnderItIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AnOddPowerOfATangentBesideAFunctionOfTheSquaredSineUnderItIntegralTest.cs new file mode 100644 index 000000000..08975651a --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AnOddPowerOfATangentBesideAFunctionOfTheSquaredSineUnderItIntegralTest.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 +{ + /// + /// An odd power of a tangent beside a root of a function of sin(y)^2, under + /// v = sin(y)^2. The secant's substitution ran out of time on Rubi's 4.1.7 + /// tan(e + f x)^5/(a + b sin(e + f x)^2)^(5/2). + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AnOddPowerOfATangentBesideAFunctionOfTheSquaredSineUnderItIntegralTest + { + [Theory] + [InlineData("tan(k + f*x)^5/(a + b*sin(k + f*x)^2)^(5/2)", 1.3, 0.7)] + [InlineData("tan(k + f*x)^3/(a + b*sin(k + f*x)^2)^(5/2)", 1.3, 0.7)] + [InlineData("tan(k + f*x)^5/sqrt(a + b*sin(k + f*x)^2)", 1.9, -0.6)] + [InlineData("tan(k + f*x)/(a - a*sin(k + f*x)^2)^(3/2)", 1.3, 0.7)] + public void UnderTheSquaredSine(string integrand, double a, double b) + { + 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", a).Substitute("b", b).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("f", 0.6).Substitute("k", 0.2); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 2.9, 4.1, 5.3 }) + { + 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-8 * 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/ByTakingOutRootsConstantUpToAPowerOfTheVariableIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ByTakingOutRootsConstantUpToAPowerOfTheVariableIntegralTest.cs new file mode 100644 index 000000000..333cd0be5 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ByTakingOutRootsConstantUpToAPowerOfTheVariableIntegralTest.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 +{ + /// + /// Roots that are a constant times a whole power of the variable, though written with it, taken + /// out: sqrt(b - a/x^2)/sqrt(a - b x^2) is K/x with K^2 = -1, and Rubi's 4.7.7 + /// sin(x) sqrt(b - a/x^2)/sqrt(a - b x^2) is K Si(x). + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ByTakingOutRootsConstantUpToAPowerOfTheVariableIntegralTest + { + [Theory] + [InlineData("sin(x)*sqrt(b - a/x^2)/sqrt(a - b*x^2)", 1.3, 0.7)] + [InlineData("sin(x)*sqrt(b - a/x^2)/sqrt(a - b*x^2)", 0.6, 1.9)] + [InlineData("cos(x)*sqrt(a*x^2 + b)/sqrt(a + b/x^2)", 1.3, 0.7)] + public void TakenOut(string integrand, double a, double b) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 2000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) => e.Substitute("a", a).Substitute("b", b).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("f", 0.6).Substitute("k", 0.2); + 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.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-8 * 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}"); + } + } +}