diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index be5134275..18bcbfe71 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -4505,6 +4505,90 @@ Entity Integral(ERational power) return expr * inTermsOfW.Substitute(w, wOfX) / inX; } + /// + /// sec(y) (a ± a sec(y))^m (c ∓ c sec(y))^n, with n half-odd and m anything, + /// integrated in u = a ± a sec(y): the two sums multiply to -a c tan(y)^2, so the + /// root of the second is tan(y) over the root of u times a constant, + /// sec(y) tan(y) dy is du over one, and what is left is + /// u^(m - 1/2) (2 - u/a)^(n - 1/2). The cosecant's are the same with the cotangent. + /// + /// + /// Rubi's sec(e + f x) (a + a sec(e + f x))^m sqrt(c - c sec(e + f x)) and two more of + /// 4.5.2.3 were declined. The constants are not written: the answer is the integrand times + /// the antiderivative in u over what that differentiates back to, a quotient whose + /// square is one; and it is kept only where that square is one at the sampled points. + /// https://github.com/asc-community/AngouriMath/issues/718 + /// + internal static Entity? SolveASecantBesidePowersOfItsConjugateSums(Entity expr, Entity.Variable x, bool integrateByParts) + { + if (!Integration.AnsweringTheQuestionAskedOrOneBelow) + return null; + Entity? argument = null; + foreach (var node in expr.Nodes) + { + if (TrigonometricArgument(node) is not { } thisArgument || !thisArgument.ContainsNode(x)) + continue; + if (argument is null) + argument = thisArgument; + else if (argument != thisArgument) + return null; + } + if (argument is null || !TreeAnalyzer.TryGetPolyLinear(argument, x, out var rate, out _) + || rate.ContainsNode(x) || TreeAnalyzer.IsZero(rate)) + return null; + var secant = MathS.Sec(argument); + var cosecant = new Cosecantf(argument); + Entity? function = null; + var sums = new List<(Entity Radicand, Entity Exponent, Entity A, bool Plus)>(); + foreach (var (factor, underneath) in FactorsOfTheIntegrand(expr)) + { + if (!factor.ContainsNode(x)) + continue; + if ((factor == secant || factor == cosecant) && function is null && !underneath) + { + function = factor; + continue; + } + var (radicand, exponent) = factor is Powf(var @base, var power) ? (@base, power) : (factor, (Entity)Number.Integer.One); + if (exponent.ContainsNode(x) || ReadAsOnePlusMinusAFunction(radicand, secant, cosecant) is not var (sumConstant, plus, _)) + return null; + sums.Add((radicand, underneath ? (-exponent).InnerSimplified : exponent, sumConstant, plus)); + } + if (function is null || sums.Count != 2 || sums[0].Plus == sums[1].Plus + || ReadAsOnePlusMinusAFunction(sums[0].Radicand, secant, cosecant) is not var (_, _, firstIsSecant) + || ReadAsOnePlusMinusAFunction(sums[1].Radicand, secant, cosecant) is not var (_, _, secondIsSecant) + || firstIsSecant != secondIsSecant || firstIsSecant != (function == secant)) + return null; + // The half-odd one, which is written through the other; a power of it at least a half + // first, so that what is left of it is a whole power to expand. + static bool IsHalfOdd(Entity exponent) => exponent is Number.Rational rational and not Number.Integer + && rational.ERational.Denominator.Equals(EInteger.FromInt32(2)) && rational.ERational.Numerator.CanFitInInt32(); + var order = sums[0].Exponent is Number.Rational { ERational.Sign: > 0 } && IsHalfOdd(sums[0].Exponent) ? new[] { 0, 1 } : new[] { 1, 0 }; + var (half, other) = (sums[order[0]], sums[order[1]]); + if (!IsHalfOdd(half.Exponent)) + return null; + var n = ((Number.Rational)half.Exponent).ERational; + // With u the other sum, the half-odd one is c (2 - u/a), whatever the signs. + var u = Variable.CreateUnique(expr, "u_conjugate"); + var c = half.A; + var a = other.A; + var leftOver = n.Subtract(ERational.FromInt32(1).Divide(ERational.FromInt32(2))); + var inU = (MathS.Pow(u, other.Exponent - Number.Rational.Create(1, 2)) + * MathS.Pow(2 - u / a, Number.Integer.Create(leftOver.ToEInteger()))).InnerSimplified; + inU = Functions.PartialFractions.Bare(inU); + if (inU.ContainsNode(x) || inU.Nodes.Any(node => node == MathS.NaN)) + return null; + var inX = inU.Substitute(u, other.Radicand) * other.Radicand.Differentiate(x); + // The square of the constant: -c^(2n)/(a rate^2). + var constantSquared = -MathS.Pow(c, Number.Integer.Create(n.Multiply(ERational.FromInt32(2)).ToEInteger())) / (a * MathS.Sqr(rate)); + if (!Functions.PartialFractions.HoldsAtSampledPoints(constantSquared * MathS.Sqr(inX), MathS.Sqr(expr), x)) + return null; + if (Integration.ComputeAsAQuestionOfItsOwn(inU, u, integrateByParts) is not { } inTermsOfU + || inTermsOfU.Nodes.Any(node => node == MathS.NaN)) + return null; + return expr * inTermsOfU.Substitute(u, other.Radicand) / inX; + } + /// /// tan(y) tan(2y) written sec(2y) - 1, and the integrand asked again in the /// one argument 2y. diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index c408d4161..105c4a753 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -908,6 +908,9 @@ private static Entity Normalized(Entity expr, Entity.Variable x) => // And beside a half-odd power of `c + d sec(y)`, which the half-angle tangent leaves as a // second root: in sin(y/2), where the two powers of cos(y) they bring make a whole one. if ((answer = IndefiniteIntegralSolver.SolveTwoHalfOddPowersOfSecantSumsByTheHalfAngleSine(expr, x, integrateByParts)) is { }) return answer; + // And the secant beside powers of a + a sec(y) and c - c sec(y), whose product is a + // square of the tangent: in u = a + a sec(y), a symbolic power of it too. + if ((answer = IndefiniteIntegralSolver.SolveASecantBesidePowersOfItsConjugateSums(expr, x, integrateByParts)) is { }) return answer; // And `a ± a cosh(y)` under a fractional power: `2a cosh(y/2)^2`, `-2a sinh(y/2)^2`. if ((answer = IndefiniteIntegralSolver.SolveByTheHalfAngleWhereOnePlusAHyperbolicCosineIsASquare(expr, x, integrateByParts)) is { }) return answer; // The hyperbolic sine's, off the real line: 1 + i sinh(y) is (cosh(y/2) + i sinh(y/2))^2. diff --git a/Sources/Tests/UnitTests/Calculus/ASecantBesideConjugateSumsIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ASecantBesideConjugateSumsIntegralTest.cs new file mode 100644 index 000000000..22c1bb9ca --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ASecantBesideConjugateSumsIntegralTest.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 +{ + /// + /// sec(y) (a ± a sec(y))^m (c ∓ c sec(y))^n, with n half-odd and m a + /// symbol, in u = a ± a sec(y), and the cosecant's. Rubi's 4.5.2.3. + /// #718 + /// + [Trait("Area", "Calculus")] + public sealed class ASecantBesideConjugateSumsIntegralTest + { + [Theory] + [InlineData("sec(g + f*x)*(a + a*sec(g + f*x))^m*(c - c*sec(g + f*x))^(5/2)")] + [InlineData("sec(g + f*x)*(a + a*sec(g + f*x))^m*sqrt(c - c*sec(g + f*x))")] + [InlineData("sec(x)*(a - a*sec(x))^m*sqrt(c + c*sec(x))")] + [InlineData("csc(x)*(a + a*csc(x))^m*sqrt(c - c*csc(x))")] + public void InTheSum(string integrand) + { + 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", 1.3).Substitute("c", 0.9).Substitute("m", 2.3).Substitute("g", 0.2).Substitute("f", 0.7); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.2, -0.7, 0.3, 0.8, 1.3, 2.9 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + if (want.IsNaN) + continue; + compared++; + Assert.True(Math.Abs((double)(got - want).RealPart) + Math.Abs((double)(got - want).ImaginaryPart) + < 1e-9 * Math.Max(1, Math.Abs((double)want.RealPart) + Math.Abs((double)want.ImaginaryPart)), + $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}"); + } + Assert.True(compared >= 5, $"only {compared} points could be compared for {integrand}"); + } + } +}