From ff232fe993a9370dd37eb480dd18091b39b72f66 Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Wed, 21 Oct 2026 15:28:12 +0000 Subject: [PATCH] Read hyperbolic integrands through cosh^2 = (1 + i sinh)(1 - i sinh), a symbolic power's own variable, and a phase MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Five shapes of Rubi's 6.1.1 and 6.7.1, each through an identity the substitutions did not use: - (g + f x)^m cosh(c + d x)^k/(a + i a sinh(c + d x)) and its sinh^k rows: cosh^2 = (1 + i sinh)(1 - i sinh), so the numerator divides by the denominator, leaving a polynomial and a constant over it. - 1/(sech(x) ± i tanh(x))^4 is cosh(x)^4/(1 ± i sinh(x))^4, which the same identity cancels to (1 ∓ i sinh)^2/(1 ± i sinh)^2; cosh^(2m) over (1 ± i sinh)^(2m) and the other way up only. - A symbolic power of sinh, cosh, sech or csch beside an odd power of the other, under u = sinh, cosh, sech or csch; a power of csch stays a power of u, since it is not sinh^(-n) where the sine is negative. - (cosh(x) coth(x))^(3/2) and (csch(x) + sinh(x))^(3/2) are (cosh(x)^2/sinh(x))^(3/2), cosh^3 sinh^(-3/2) times a constant on each side of zero. - Half-odd powers of S = b cosh(x) + c sinh(x) ± sqrt(b^2 - c^2): with T = b sinh + c cosh, d/dx (T S^(n-1)) = n S^n - a (2n - 1) S^(n-1) steps the power to ±1/2. The whole powers are left to the exponential substitution, which writes 1/S^2, 1/S^3 and 1/S^4 in 182, 272 and 416 characters, where the recurrence gave 350, 539 and 728. Chain against master 4877ed2b: pocket of 1,756 rows, 1,686 -> 1,723 solved, 0 wrong; samples unchanged or better; the 43 rows that moved, alone, 4 -> 42 solved. Re-measured on master f7f2a7a8 after the rebase: 14 of the 43 are new there (12 of (e + f x)^k cosh^m/(a + i a sinh) past 20 s -> 0.1-0.9 s, (csch(x) + sinh(x))^(3/2) and ^(5/2) declined -> answered), and the other 29 answer on both, most much cheaper here: sech(a + b x)^n tanh^3 9.12 -> 0.02 GB, 1/(sech(x) + i tanh(x))^4 10.5 -> 0.07 GB and 551 -> 241 characters. Part of #233. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura --- BREAKING-CHANGES.md | 35 ++ .../Integration/IndefiniteIntegralSolver.cs | 568 ++++++++++++++++++ .../Integration/Integration.Definition.cs | 16 + ...neAndSinePlusTheirAmplitudeIntegralTest.cs | 94 +++ ...nBesideAnOddPowerOfTheOtherIntegralTest.cs | 80 +++ ...SineAgainstTheCosineSquaredIntegralTest.cs | 50 ++ ...usAnImaginaryHyperbolicSineIntegralTest.cs | 54 ++ ...licCosineSquaredOutOfAPowerIntegralTest.cs | 50 ++ 8 files changed, 947 insertions(+) create mode 100644 Sources/Tests/UnitTests/Calculus/APowerOfAHyperbolicCosineAndSinePlusTheirAmplitudeIntegralTest.cs create mode 100644 Sources/Tests/UnitTests/Calculus/ASymbolicPowerOfAHyperbolicFunctionBesideAnOddPowerOfTheOtherIntegralTest.cs create mode 100644 Sources/Tests/UnitTests/Calculus/CancellingOnePlusAnImaginaryHyperbolicSineAgainstTheCosineSquaredIntegralTest.cs create mode 100644 Sources/Tests/UnitTests/Calculus/DividingOutOnePlusAnImaginaryHyperbolicSineIntegralTest.cs create mode 100644 Sources/Tests/UnitTests/Calculus/TakingTheHyperbolicCosineSquaredOutOfAPowerIntegralTest.cs diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 0909ba42e..b4b2bab18 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -741,6 +741,41 @@ power, and the integrand is answered in a fraction of a second | `"csch(c+d*x)*(a+b*sinh(c+d*x)^4)^2".ToEntity().Integrate("x")` | past fifteen minutes | 426 characters | | `"csch(c+d*x)*(a+b*sinh(c+d*x)^4)^3".ToEntity().Integrate("x")` | past fifteen minutes | 646 characters | +### Hyperbolic integrands are read through `cosh^2 = (1 + i sinh)(1 - i sinh)`, a symbolic power's own variable, and a phase + +**Answers where there were none, and shorter ones.** Rows of these shapes in Rubi's 6.1.1 and 6.7.1 were declined or ran +out of time. Each is answered now through an identity of the hyperbolic functions that the substitutions did not use +([#233](https://github.com/asc-community/AngouriMath/issues/233)): +- `(g + f x)^m cosh(c + d x)^k/(a + i a sinh(c + d x))` and its `sinh^k` rows. `cosh^2 = (1 + i sinh)(1 - i sinh)`, so + the numerator divides by the denominator, leaving a polynomial and a constant over it. +- `1/(sech(x) ± i tanh(x))^4`. It is `cosh(x)^4/(1 ± i sinh(x))^4`, which the same identity cancels to + `(1 ∓ i sinh(x))^2/(1 ± i sinh(x))^2`. Only `cosh^(2m)` over `(1 ± i sinh)^(2m)` and the other way up: cancelling + where a polynomial or another power of the binomial is left made other answers longer. +- A symbolic power of `sinh`, `cosh`, `sech` or `csch` beside an odd power of the other function, with `tanh` and + `coth` read as quotients. `cosh(a + b x)^3 sinh(a + b x)^n` is integrated under `u = sinh`, and + `sech(a + b x)^n tanh(a + b x)^3` under `u = sech`. The power of `csch` stays a power of `u`, because it is not + `sinh^(-n)` where the sine is negative. +- `(cosh(x) coth(x))^(3/2)` and `(csch(x) + sinh(x))^(3/2)`. Both are `(cosh(x)^2/sinh(x))^(3/2)`, which is + `cosh(x)^3 sinh(x)^(-3/2)` times a constant on each side of zero, odd in the cosine. +- Half-odd powers of `S = b cosh(x) + c sinh(x) ± sqrt(b^2 - c^2)`. With `T = b sinh(x) + c cosh(x)` and `a` the + constant term, `d/dx (T S^(n - 1)) = n S^n - a (2n - 1) S^(n - 1)` steps the power to `±1/2`. The whole powers stay + the exponential substitution's, which writes them shorter. + +2.5.0 declined or ran out of time on all of the inputs below. The unreleased master answered some of them, and +those answers are now two to twelve times shorter: `(sech(x) + i tanh(x))^2` and `^4` were 437 and 2,916 +characters, and `1/(sech(x) + i tanh(x))^2` and `^4` were 372 and 551. + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"(g+f*x)*cosh(c+d*x)^2/(a+i*a*sinh(c+d*x))".ToEntity().Integrate("x")` | `integral(...)` | 366 characters | +| `"cosh(a+b*x)^3*sinh(a+b*x)^n".ToEntity().Integrate("x")` | `integral(...)` | 149 characters | +| `"sech(a+b*x)^n*tanh(a+b*x)^3".ToEntity().Integrate("x")` | past two minutes | 183 characters | +| `"coth(a+b*x)^3*csch(a+b*x)^n".ToEntity().Integrate("x")` | past a minute | 242 characters | +| `"(cosh(x)*coth(x))^(3/2)".ToEntity().Integrate("x")` | `integral(...)` | 247 characters | +| `"1/sqrt(b*cosh(x)+c*sinh(x)+sqrt(b^2-c^2))".ToEntity().Integrate("x")` | `integral(...)` | 409 characters | +| `"1/(sech(x)+i*tanh(x))^4".ToEntity().Integrate("x")` | `integral(...)` | 241 characters | +| `"(sech(x)+i*tanh(x))^4".ToEntity().Integrate("x")` | `integral(...)`, after 54 seconds | 251 characters, where the unreleased master gave 2,916 | + ### A power of a polynomial whose leading coefficient is zero spelled otherwise is not divided by it **Answers where there were none.** `1/((a (1 + x^2))^2 - (a (1 - x^2))^2)^3` was left unevaluated. That is diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 31dc3156e..f6c7a8fcb 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -5649,6 +5649,101 @@ Entity Integral(ERational power) /// private const int MaximumAmplitudePower = 8; + /// + /// (a + b cosh(y) + c sinh(y))^n for a^2 = b^2 - c^2, y linear in + /// x and n half-odd, in closed form: + /// 's recurrence for the + /// hyperbolic functions. With S the base and T = b sinh(y) + c cosh(y), + /// T' = S - a, S' = T and T^2 = S (S - 2a), the last of which is + /// a^2 = b^2 - c^2, and from them + /// + /// d/dy (T S^(n-1)) = n S^n - a (2n - 1) S^(n-1) + /// + /// which steps n down to 1/2, where the integral is 2T/sqrt(S), or up + /// to -1/2, where it is (ln(T - i q sqrt(S)) - ln(T + i q sqrt(S)))/(i q) for + /// q = sqrt(2a). + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + /// + /// Rubi's sqrt(b cosh(x) + c sinh(x) + sqrt(b^2 - c^2)) and the rest of the twelve in + /// its 6.7.1 were declined or past the budget: b cosh(y) + c sinh(y) is + /// sqrt(b^2 - c^2) cosh(y + phi), so the base is a square of the half angle turned by a + /// phase, and reads + /// a ± a cosh(y) without one. The recurrence needs no phase. + /// + /// + /// The base case is two logarithms rather than the arctangent of T/(q sqrt(S)): for + /// real b > |c|, S is not negative whichever sign a has, and for a + /// negative a that argument is imaginary and past i, on the arctangent's cut. + /// The two logarithms are of complex conjugates for a positive a and of two numbers of + /// one sign for a negative one, real and continuous between the zeros of S either way. + /// The derivations use only S^(3/2) = S sqrt(S) and (q sqrt(S))^2 = 2a S, which + /// hold for the principal powers. a^2 = b^2 - c^2 is decided, not assumed, and all three + /// are there: b = 0 is a ± i a sinh(y) and c = 0 is a ± a cosh(y), + /// which the half-angle rules answer. + /// + /// + /// Half-odd powers only. The whole ones are the exponential substitution's, which answers + /// 1/S^2, 1/S^3 and 1/S^4 in 182, 272 and 416 characters, where the + /// same recurrence stepped up to T/(a S) gives 350, 539 and 728. + /// + /// + internal static Entity? SolveAPowerOfAHyperbolicCosineAndSinePlusTheirAmplitude(Entity expr, Entity.Variable x) + { + var sh = Variable.CreateUnique(expr, "s_hyp"); + var ch = Variable.CreateUnique(expr, "c_hyp"); + if (ReadTheHyperbolicFunctions(expr, x, sh, ch) is not var (read, argument) || read.ContainsNode(x) + || !TreeAnalyzer.TryGetPolyLinear(argument, x, out var rate, out _) || rate.ContainsNode(x) || TreeAnalyzer.IsZero(rate)) + return null; + // A constant times one power of the base, above the bar or below it. + Entity constant = Number.Integer.One; + Entity? @base = null; + ERational? n = null; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(read)) + { + if (!factor.ContainsNode(sh) && !factor.ContainsNode(ch)) + { + constant = underneath ? constant / factor : constant * factor; + continue; + } + if (@base is not null || factor is not Powf(var inner, var power) || power.Evaled is not Number.Rational rational) + return null; + (@base, n) = (inner, underneath ? rational.ERational.Negate() : rational.ERational); + } + // Half-odd: the whole powers are answered shorter by the exponential substitution. + if (@base is null || n is null || !n.Denominator.Equals(EInteger.FromInt32(2)) + || n.Abs().CompareTo(ERational.FromInt32(MaximumAmplitudePower)) > 0) + return null; + // a + b cosh(y) + c sinh(y), each coefficient free of both. + if (!TreeAnalyzer.TryGetPolyLinear(@base, ch, out var b, out var rest) || !TreeAnalyzer.TryGetPolyLinear(rest, sh, out var c, out var a) + || b.ContainsNode(sh) || b.ContainsNode(ch) || c.ContainsNode(sh) || c.ContainsNode(ch) || a.ContainsNode(sh) || a.ContainsNode(ch)) + return null; + if (VanishesIdentically(a) || VanishesIdentically(b) || VanishesIdentically(c) + || !VanishesIdentically(MathS.Sqr(a) - MathS.Sqr(b) + MathS.Sqr(c))) + return null; + + var s = a + b * MathS.Hyperbolic.Cosh(argument) + c * MathS.Hyperbolic.Sinh(argument); + var t = b * MathS.Hyperbolic.Sinh(argument) + c * MathS.Hyperbolic.Cosh(argument); + var half = ERational.Create(1, 2); + Entity Integral(ERational power) + { + if (power.CompareTo(half) == 0) + return 2 * t / MathS.Sqrt(s); + if (power.CompareTo(half.Negate()) == 0) + { + var q = MathS.Sqrt(2 * a); + return (MathS.Ln(t - MathS.i * q * MathS.Sqrt(s)) - MathS.Ln(t + MathS.i * q * MathS.Sqrt(s))) / (MathS.i * q); + } + var current = Number.Rational.Create(power); + if (power.Sign > 0) + return t * MathS.Pow(s, Number.Rational.Create(power.Subtract(ERational.One))) / current + + a * (2 * current - 1) / current * Integral(power.Subtract(ERational.One)); + return ((current + 1) * Integral(power.Add(ERational.One)) - t * MathS.Pow(s, current)) / (a * (2 * current + 1)); + } + return (constant * Integral(n) / rate).InnerSimplified; + } + /// /// Fractional powers of a ± a sin(y), or of a ± a cos(y), beside anything /// rational in the sine and cosine of y, by the half angle at which they are @@ -35038,6 +35133,479 @@ static bool IsZero(Entity value) => value.InnerSimplified.Evaled is Number.Compl return null; } + /// + /// A polynomial in sinh(y) and cosh(y) over a ± i a sinh(y), beside + /// anything else in x, divided out: cosh(y)^2 = (1 + i sinh(y))(1 - i sinh(y)), + /// so with cosh(y)^2 written as 1 + sinh(y)^2 the numerator is + /// N0(sinh) + cosh N1(sinh), and each part divides by the linear + /// a ± i a sinh(y) to a polynomial and a constant. (e + f x) cosh(c + d x)^2/(a + i a sinh(c + d x)) + /// is (e + f x)(1 - i sinh(c + d x))/a, and (e + f x) sinh(y)^2/(a + i a sinh(y)) + /// is (e + f x)(1 - i sinh(y))/a - (e + f x)/(a + i a sinh(y)). + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + /// Rubi's (e + f x)^m cosh(c + d x)^k/(a + i a sinh(c + d x)) and its sinh^k + /// rows were declined or past the budget, where the pieces the division leaves are answered + /// in a fraction of a second. Exact wherever the denominator is not zero. Only where the remainder beside + /// the cosine vanishes, since cosh(y)/(a ± i a sinh(y)) beside a power of x is + /// a polylogarithm by parts; where the quotient keeps a hyperbolic function, a constant + /// there being one step from the quotient the rules already answer; and where something + /// else in x multiplies it, a rational function of the exponential alone being the + /// exponential substitution's. + /// + internal static Entity? SolveByDividingOutOnePlusAnImaginaryHyperbolicSine(Entity expr, Entity.Variable x, bool integrateByParts) + { + var sh = Variable.CreateUnique(expr, "s_hyp"); + var ch = Variable.CreateUnique(expr, "c_hyp"); + if (ReadTheHyperbolicFunctions(expr, x, sh, ch) is not var (read, argument)) + return null; + // One linear in the sine below the bar, a polynomial in the two above, the rest beside. + Entity? denominator = null; + Entity numerator = Number.Integer.One; + Entity beside = Number.Integer.One; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(read)) + { + if (!factor.ContainsNode(sh) && !factor.ContainsNode(ch)) + beside = underneath ? beside / factor : beside * factor; + else if (!underneath) + numerator = numerator == Number.Integer.One ? factor : numerator * factor; + else if (denominator is null) + denominator = factor; + else + return null; + } + if (denominator is null || numerator == Number.Integer.One || !beside.ContainsNode(x) + || numerator.ContainsNode(x) || denominator.ContainsNode(x) || denominator.ContainsNode(ch) + || ReadAsOnePlusMinusAnImaginarySine(denominator, sh) is not var (a, plus)) + return null; + // a (1 + i s) = i a (s - i), a (1 - i s) = -i a (s + i). + var root = plus ? MathS.i : -MathS.i; + var leading = plus ? MathS.i * a : -MathS.i * a; + // N0 + c N1 with c^2 = 1 + s^2. + if (!TreeAnalyzer.TryGetPolynomial(numerator, ch, out var byCosine)) + return null; + Entity even = Number.Integer.Zero, odd = Number.Integer.Zero; + foreach (var pair in byCosine) + { + var (power, coefficient) = (pair.Key, pair.Value); + if (power.Sign < 0 || power.CompareTo(EInteger.FromInt32(MaximumDividedPower)) > 0 || coefficient.ContainsNode(ch)) + return null; + var half = power.ToInt32Checked() / 2; + var term = half == 0 ? coefficient : coefficient * MathS.Pow(1 + MathS.Sqr(sh), half); + if (power.IsEven) + even += term; + else + odd += term; + } + if (Divided(even) is not var (evenQuotient, evenRemainder) || Divided(odd) is not var (oddQuotient, oddRemainder) + || !VanishesIdentically(oddRemainder) || evenQuotient.Count + oddQuotient.Count == 0 + || evenQuotient.Count <= 1 && oddQuotient.Count == 0) + return null; + var sinh = MathS.Hyperbolic.Sinh(argument); + var cosh = MathS.Hyperbolic.Cosh(argument); + Entity? rewritten = null; + void Add(Entity term) => rewritten = rewritten is null ? term : rewritten + term; + foreach (var (part, function) in new[] { (evenQuotient, (Entity)Number.Integer.One), (oddQuotient, cosh) }) + for (var k = 0; k < part.Count; k++) + if (!VanishesIdentically(part[k])) + { + Entity power = k == 0 ? Number.Integer.One : k == 1 ? sinh : MathS.Pow(sinh, k); + var hyperbolic = function == Number.Integer.One ? power : power == Number.Integer.One ? function : function * power; + // Free of x and small, so simplified: `i/(i a)` is otherwise left as written. + Add(Functions.PartialFractions.Bare((part[k] / leading).Simplify()) * hyperbolic * beside); + } + if (!VanishesIdentically(evenRemainder)) + Add(evenRemainder * beside / denominator.Substitute(sh, sinh)); + if (rewritten is null || Integration.ComputeAsTheSameQuestion(rewritten, x, integrateByParts) is not { } answer + || answer.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + return answer; + + // p(s) = (s - root) q(s) + r, by synthetic division: q's coefficients from the + // constant up, and r. + (List Quotient, Entity Remainder)? Divided(Entity polynomial) + { + if (!TreeAnalyzer.TryGetPolynomial(polynomial, sh, out var bySine) || bySine.Keys.Any(k => k.Sign < 0) || bySine.Values.Any(v => v.ContainsNode(sh))) + return null; + var highest = bySine.Keys.Aggregate(EInteger.Zero, (most, k) => k.CompareTo(most) > 0 ? k : most); + if (highest.CompareTo(EInteger.FromInt32(2 * MaximumDividedPower)) > 0) + return null; + var degree = highest.ToInt32Checked(); + Entity Coefficient(int k) => bySine.TryGetValue(EInteger.FromInt32(k), out var value) ? value : Number.Integer.Zero; + var quotient = new Entity[degree]; + Entity carried = Number.Integer.Zero; + for (var k = degree; k >= 1; k--) + { + carried = (Coefficient(k) + root * carried).InnerSimplified; + quotient[k - 1] = carried; + } + return (quotient.ToList(), (Coefficient(0) + root * carried).InnerSimplified); + } + } + + /// + /// The largest power of the hyperbolic cosine + /// writes in the sine: each is a binomial expanded. + /// + private const int MaximumDividedPower = 12; + + /// + /// A symbolic power of sinh(y) beside an odd whole power of cosh(y), under + /// u = sinh(y), where cosh(y)^(2j + 1) dy is (1 + u^2)^j du; and a + /// symbolic power of cosh(y) beside an odd power of sinh(y), under + /// u = cosh(y), with (u^2 - 1)^j. cosh(a + b x)^3 sinh(a + b x)^n is + /// (1 + u^2) u^n/b. A symbolic power of sech(y) or csch(y) the same way + /// under u = sech(y) or u = csch(y), with tanh(y) and coth(y) + /// beside it read as quotients: sech(a + b x)^n tanh(a + b x)^3 is + /// sinh^3 cosh^(-3) sech^n, which is -(1 - u^2) u^(n - 1)/b. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + /// Rubi's four in 6.7.1 were past the budget: + /// reads whole and rational powers, and + /// a symbolic one is neither. Exact: d/dx F(sinh(y)) = y' cosh(y) F'(sinh(y)) wherever + /// F is differentiable, whatever the power, and (1 + sinh(y)^2)^j = cosh(y)^(2j) + /// for a whole j. Odd powers from three: the first is the derivative of a power beside + /// it, which the rules for that answer. So were its sech(a + b x)^n tanh(a + b x)^3 and + /// coth(a + b x)^3 csch(a + b x)^n. The reciprocal is the variable, not its reciprocal + /// to the power negated: csch(y)^n is e^(2 pi i n) sinh(y)^(-n) where the sine is + /// negative, and u^n under u = csch(y) is csch(y)^n everywhere. + /// + internal static Entity? SolveASymbolicPowerOfAHyperbolicFunctionBesideAnOddPowerOfTheOther(Entity expr, Entity.Variable x, bool integrateByParts) + { + var sh = Variable.CreateUnique(expr, "s_hyp"); + var ch = Variable.CreateUnique(expr, "c_hyp"); + if (ReadTheHyperbolicFunctions(expr, x, sh, ch) is not var (read, argument) || read.ContainsNode(x) + || !TreeAnalyzer.TryGetPolyLinear(argument, x, out var rate, out _) || rate.ContainsNode(x) || TreeAnalyzer.IsZero(rate)) + return null; + Entity constant = Number.Integer.One; + Entity? symbolic = null, exponent = null, odd = null; + var oddPower = 0; + foreach (var (factor, underneath) in FactorsOfTheIntegrand(read)) + { + if (!factor.ContainsNode(sh) && !factor.ContainsNode(ch)) + { + constant = underneath ? constant / factor : constant * factor; + continue; + } + if (factor is Powf(var function, var power) && (function == sh || function == ch) && power.Evaled is not Number) + { + if (symbolic is not null) + return null; + (symbolic, exponent) = (function, underneath ? -power : power); + continue; + } + var (other, whole) = factor is Powf(var b, Number.Integer k) ? (b, k.EInteger) : (factor, EInteger.One); + if (underneath || odd is not null || other != sh && other != ch + || whole.IsEven || whole.Sign <= 0 || whole.CompareTo(EInteger.FromInt32(MaximumDividedPower)) > 0) + return BesideTheReciprocalsAndTheQuotients(); + (odd, oddPower) = (other, whole.ToInt32Checked()); + } + if (symbolic is null || exponent is null || odd is null || symbolic == odd || oddPower < 3) + return null; + var u = Variable.CreateUnique(expr, "u_hyp"); + var ofTheCosine = odd == ch; + var square = ofTheCosine ? 1 + MathS.Sqr(u) : MathS.Sqr(u) - 1; + var half = (oddPower - 1) / 2; + var question = constant * (half == 1 ? square : MathS.Pow(square, half)) * MathS.Pow(u, exponent) / rate; + if (Integration.ComputeAsAQuestionOfItsOwn(question, u, integrateByParts) is not { } answer + || answer.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + return answer.Substitute(u, ofTheCosine ? MathS.Hyperbolic.Sinh(argument) : MathS.Hyperbolic.Cosh(argument)); + + // The symbolic power of sinh(y), cosh(y), csch(y) or sech(y), beside whole powers of the + // two -- tanh(y)^3 is sinh(y)^3 cosh(y)^(-3) -- under u = that function. Of the + // reciprocals: du = -cosh(y) sinh(y)^(-2) dy under u = csch(y), where cosh(y)^2 is + // (1 + u^2)/u^2, and du = -sinh(y) cosh(y)^(-2) dy under u = sech(y), where sinh(y)^2 is + // (1 - u^2)/u^2. The power stays the power of u: csch(y)^n is not sinh(y)^(-n) where + // the sine is negative. + Entity? BesideTheReciprocalsAndTheQuotients() + { + Entity coefficient = Number.Integer.One; + Entity? function = null, power = null; + var (sines, cosines, functionSines, functionCosines) = (0, 0, 0, 0); + foreach (var (factor, underneath) in FactorsOfTheIntegrand(read)) + { + if (factor is Powf(var @base, var p) && p.Evaled is not Number && (@base.ContainsNode(sh) || @base.ContainsNode(ch))) + { + if (function is not null || PowersOfTheTwo(@base) is not var (one, s, c) || one != Number.Integer.One || System.Math.Abs(s) + System.Math.Abs(c) != 1) + return null; + (function, power, functionSines, functionCosines) = (@base, underneath ? -p : p, s, c); + continue; + } + if (PowersOfTheTwo(factor) is not var (k, ks, kc)) + return null; + var sign = underneath ? -1 : 1; + (sines, cosines) = (sines + sign * ks, cosines + sign * kc); + if (k != Number.Integer.One) + coefficient = underneath ? coefficient / k : coefficient * k; + } + // The other's power after du takes one of it, odd, so that what is left is even. + var otherPower = functionSines != 0 ? cosines : sines; + if (function is null || power is null || otherPower < 3 || otherPower % 2 == 0 || otherPower > MaximumDividedPower) + return null; + var j = (otherPower - 1) / 2; + var w = Variable.CreateUnique(expr, "u_hyp"); + Entity Power(Entity of, int n) => n == 0 ? Number.Integer.One : n == 1 ? of : MathS.Pow(of, n); + var (inW, back) = (functionSines, functionCosines) switch + { + (1, _) => (coefficient * Power(1 + MathS.Sqr(w), j) * Power(w, sines), MathS.Hyperbolic.Sinh(argument)), + (_, 1) => (coefficient * Power(MathS.Sqr(w) - 1, j) * Power(w, cosines), MathS.Hyperbolic.Cosh(argument)), + (-1, _) => (-coefficient * Power(1 + MathS.Sqr(w), j) * Power(w, -sines - 2 - 2 * j), MathS.Hyperbolic.Cosech(argument)), + _ => (-coefficient * Power(1 - MathS.Sqr(w), j) * Power(w, -cosines - 2 - 2 * j), MathS.Hyperbolic.Sech(argument)), + }; + if (Integration.ComputeAsAQuestionOfItsOwn(inW * MathS.Pow(w, power) / rate, w, integrateByParts) is not { } result + || result.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + return result.Substitute(w, back); + } + + // A product of whole powers of the two and a constant: the constant, and the powers. + (Entity Constant, int Sines, int Cosines)? PowersOfTheTwo(Entity node) + { + if (node == sh) + return (Number.Integer.One, 1, 0); + if (node == ch) + return (Number.Integer.One, 0, 1); + if (!node.ContainsNode(sh) && !node.ContainsNode(ch)) + return (node, 0, 0); + switch (node) + { + case Mulf(var left, var right) when PowersOfTheTwo(left) is var (k, s, c) && PowersOfTheTwo(right) is var (l, t, d): + return (k == Number.Integer.One ? l : l == Number.Integer.One ? k : k * l, s + t, c + d); + case Divf(var left, var right) when PowersOfTheTwo(left) is var (k, s, c) && PowersOfTheTwo(right) is var (l, t, d): + return (l == Number.Integer.One ? k : k / l, s - t, c - d); + case Powf(var @base, Number.Integer { EInteger: var n }) when n.Abs().CompareTo(EInteger.FromInt32(2 * MaximumDividedPower)) <= 0 + && PowersOfTheTwo(@base) is var (k, s, c): + var whole = n.ToInt32Checked(); + return (k == Number.Integer.One ? k : MathS.Pow(k, whole), whole * s, whole * c); + } + return null; + } + } + + /// + /// A fractional power of a rational function of sinh(y) and cosh(y) whose + /// numerator has cosh(y)^2 as a factor -- with cosh(y)^2 written + /// 1 + sinh(y)^2 -- taken out of the power: (cosh(y)^2 A/B)^p is + /// cosh(y)^(2p) A^p B^(-p) times a constant on every interval where both are + /// continuous. (cosh(x) coth(x))^(3/2) and (csch(x) + sinh(x))^(3/2) are both + /// (cosh(x)^2/sinh(x))^(3/2), which is cosh(x)^3 sinh(x)^(-3/2) on + /// x > 0, odd in the hyperbolic cosine. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + /// Rubi's six in 6.7.1 were declined: the radicand is a quotient, and + /// reads the radicands that are + /// polynomials. The constant is not written: each power as written over its form apart is a + /// quotient whose logarithmic derivative is zero, so it is constant wherever it is + /// continuous, and the answer is the antiderivative of the rewritten integrand times their + /// product. Only where an odd power of the sine is left beside it. + /// + internal static Entity? SolveByTakingTheHyperbolicCosineSquaredOutOfAPower(Entity expr, Entity.Variable x, bool integrateByParts) + { + var sh = Variable.CreateUnique(expr, "s_hyp"); + var ch = Variable.CreateUnique(expr, "c_hyp"); + if (ReadTheHyperbolicFunctions(expr, x, sh, ch) is not var (read, argument)) + return null; + var sinh = MathS.Hyperbolic.Sinh(argument); + var cosh = MathS.Hyperbolic.Cosh(argument); + Entity constant = Number.Integer.One; + var found = false; + var rewritten = read.Replace(node => + { + if (node is not Powf(var radicand, Number.Rational exponent) || exponent is Number.Integer + || !radicand.ContainsNode(sh) && !radicand.ContainsNode(ch) || radicand.ContainsNode(x)) + return node; + var (top, bottom) = Functions.SingleQuotient.Of(Functions.SingleQuotient.Combine(radicand)); + if (InTheSine(top) is not { } above || InTheSine(bottom) is not { } below) + return node; + // cosh(y)^2 above the bar, and an odd degree in the sine left beside it: what is left + // is then no square, whose root -- `sqrt(1 + csch(x)^2)` -- is a sign the rules for + // the power of a square answer. + var squares = TakeTheSquaresOut(above) - TakeTheSquaresOut(below); + if (squares <= 0 || (above.Count + below.Count) % 2 != 1) + return node; + found = true; + Entity written = MathS.Pow(cosh, (2 * squares * exponent).InnerSimplified); + if (above.Count > 1 || above[0] != Number.Integer.One) + written *= MathS.Pow(InX(above), exponent); + if (below.Count > 1 || below[0] != Number.Integer.One) + written *= MathS.Pow(InX(below), -exponent); + constant *= MathS.Pow(radicand.Substitute(ch, cosh).Substitute(sh, sinh), exponent) / written; + return written; + }); + if (!found) + return null; + var inX = rewritten.Substitute(ch, cosh).Substitute(sh, sinh); + if (Integration.ComputeAsTheSameQuestion(inX, x, integrateByParts) is not { } result + || result.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + return constant * result; + + // The coefficients of a polynomial in the two, the cosine to even powers only, as one in + // the sine, from the constant up. + List? InTheSine(Entity polynomial) + { + if (!TreeAnalyzer.TryGetPolynomial(polynomial, ch, out var byCosine)) + return null; + Entity inSine = Number.Integer.Zero; + foreach (var pair in byCosine) + { + var (power, coefficient) = (pair.Key, pair.Value); + if (power.Sign < 0 || !power.IsEven || power.CompareTo(EInteger.FromInt32(MaximumDividedPower)) > 0 || coefficient.ContainsNode(ch)) + return null; + var half = power.ToInt32Checked() / 2; + inSine += half == 0 ? coefficient : coefficient * MathS.Pow(1 + MathS.Sqr(sh), half); + } + if (!TreeAnalyzer.TryGetPolynomial(inSine, sh, out var bySine) || bySine.Keys.Any(k => k.Sign < 0) || bySine.Values.Any(v => v.ContainsNode(sh))) + return null; + var highest = bySine.Keys.Aggregate(EInteger.Zero, (most, k) => k.CompareTo(most) > 0 ? k : most); + if (highest.CompareTo(EInteger.FromInt32(2 * MaximumDividedPower)) > 0) + return null; + var coefficients = new List(); + for (var k = 0; k <= highest.ToInt32Checked(); k++) + coefficients.Add(bySine.TryGetValue(EInteger.FromInt32(k), out var value) ? value.InnerSimplified : Number.Integer.Zero); + while (coefficients.Count > 1 && VanishesIdentically(coefficients[^1])) + coefficients.RemoveAt(coefficients.Count - 1); + return coefficients; + } + + // Divides 1 + s^2 out as often as it goes, in place, and says how often. + static int TakeTheSquaresOut(List coefficients) + { + var times = 0; + while (coefficients.Count > 2) + { + var quotient = new List(coefficients); + for (var k = quotient.Count - 1; k >= 2; k--) + quotient[k - 2] = (quotient[k - 2] - quotient[k]).InnerSimplified; + if (!VanishesIdentically(quotient[0]) || !VanishesIdentically(quotient[1])) + break; + coefficients.Clear(); + coefficients.AddRange(quotient.Skip(2)); + times++; + } + return times; + } + + Entity InX(List coefficients) + { + Entity? sum = null; + for (var k = 0; k < coefficients.Count; k++) + { + if (VanishesIdentically(coefficients[k])) + continue; + var term = k == 0 ? coefficients[k] : coefficients[k] * (k == 1 ? sinh : MathS.Pow(sinh, k)); + sum = sum is null ? term : sum + term; + } + return sum ?? Number.Integer.Zero; + } + } + + /// + /// cosh(y)^(2m) over (1 ± i sinh(y))^(2m), or the other way up, times a + /// constant, cancelled to ((1 ∓ i sinh(y))/(1 ± i sinh(y)))^m: + /// cosh(y)^2 = (1 + i sinh(y))(1 - i sinh(y)). 1/(sech(y) + i tanh(y))^4 is + /// cosh(y)^4/(1 + i sinh(y))^4, which is (1 - i sinh(y))^2/(1 + i sinh(y))^2. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + /// + /// Rubi's 1/(sech(x) ± i tanh(x))^4 in 6.7.1 were declined: the quotient is even in + /// the hyperbolic cosine, so does not + /// read it, and the substitutions spent the budget on a denominator of degree eight in the + /// exponential, where the quotient the cancellation leaves is answered in a fifth of a + /// second. Exact on the real line, where 1 ± i sinh(y) is not zero, and the points + /// off it where it is are removable. Only this shape, read by the multiplicities of + /// sinh(y) ∓ i on either side once the squares of the cosine are written in the sine: + /// cancelling where a polynomial or another power of 1 ± i sinh(y) is left made the + /// answers to cosh(x)^8/(i + sinh(x)) and cosh(x)^2/(1 + i sinh(x))^3 longer, + /// and writing the cosh(y)^2 a quotient carries on both sides as the two binomials + /// made (a sech(x) + b tanh(x))^4 a search past the budget. + /// + internal static Entity? SolveByCancellingOnePlusAnImaginaryHyperbolicSineAgainstTheCosineSquared(Entity expr, Entity.Variable x, bool integrateByParts) + { + var sh = Variable.CreateUnique(expr, "s_hyp"); + var ch = Variable.CreateUnique(expr, "c_hyp"); + if (ReadTheHyperbolicFunctions(expr, x, sh, ch) is not var (read, argument) || read.ContainsNode(x)) + return null; + var (top, bottom) = Functions.SingleQuotient.Of(Functions.SingleQuotient.Combine(read)); + if (bottom == Number.Integer.One || InTheSine(top) is not { } above || InTheSine(bottom) is not { } below) + return null; + // s - i = -i (1 + i s) and s + i = i (1 - i s): what is left of either side must be a + // constant, and the two powers opposite, so that the quotient is ((s - i)/(s + i))^e. + var (aboveAtI, belowAtI) = (TakenOut(above, MathS.i), TakenOut(below, MathS.i)); + var (aboveAtMinusI, belowAtMinusI) = (TakenOut(above, -MathS.i), TakenOut(below, -MathS.i)); + var e = aboveAtI - belowAtI; + if (above.Count != 1 || below.Count != 1 || e == 0 || e != belowAtMinusI - aboveAtMinusI || System.Math.Abs(e) > MaximumDividedPower) + return null; + // ((s - i)/(s + i))^e = (-1)^e ((1 + i s)/(1 - i s))^e. + var constant = ((e % 2 == 0 ? above[0] : -above[0]) / below[0]).InnerSimplified; + if (constant.Vars.Any()) + constant = Functions.PartialFractions.Bare(constant.Simplify()); + var sinh = MathS.Hyperbolic.Sinh(argument); + var (up, down) = e > 0 ? (1 + MathS.i * sinh, 1 - MathS.i * sinh) : (1 - MathS.i * sinh, 1 + MathS.i * sinh); + var m = System.Math.Abs(e); + Entity rewritten = m == 1 ? up / down : MathS.Pow(up, m) / MathS.Pow(down, m); + if (constant != Number.Integer.One) + rewritten = constant * rewritten; + if (Integration.ComputeAsTheSameQuestion(rewritten, x, integrateByParts) is not { } answer + || answer.Nodes.Any(node => node is Integralf || node == MathS.NaN)) + return null; + return answer; + + // A polynomial in the two even in the cosine, as one in the sine from the constant up, + // without its zero top coefficients; null where it is odd in the cosine anywhere. + List? InTheSine(Entity polynomial) + { + if (!TreeAnalyzer.TryGetPolynomial(polynomial, ch, out var byCosine)) + return null; + Entity inSine = Number.Integer.Zero; + foreach (var pair in byCosine) + { + var (power, coefficient) = (pair.Key, pair.Value); + if (power.Sign < 0 || !power.IsEven || power.CompareTo(EInteger.FromInt32(2 * MaximumDividedPower)) > 0 || coefficient.ContainsNode(ch)) + return null; + var half = power.ToInt32Checked() / 2; + inSine += half == 0 ? coefficient : coefficient * MathS.Pow(1 + MathS.Sqr(sh), half); + } + if (!TreeAnalyzer.TryGetPolynomial(inSine, sh, out var bySine) || bySine.Keys.Any(k => k.Sign < 0) || bySine.Values.Any(v => v.ContainsNode(sh))) + return null; + var highest = bySine.Keys.Aggregate(EInteger.Zero, (most, k) => k.CompareTo(most) > 0 ? k : most); + if (highest.CompareTo(EInteger.FromInt32(4 * MaximumDividedPower)) > 0) + return null; + var coefficients = new List(); + for (var k = 0; k <= highest.ToInt32Checked(); k++) + coefficients.Add(bySine.TryGetValue(EInteger.FromInt32(k), out var value) ? value.InnerSimplified : Number.Integer.Zero); + while (coefficients.Count > 0 && VanishesIdentically(coefficients[^1])) + coefficients.RemoveAt(coefficients.Count - 1); + return coefficients.Count == 0 ? null : coefficients; + } + + // Divides s - root out of the polynomial, in place, as often as it goes, by synthetic + // division, and says how often. + static int TakenOut(List polynomial, Entity root) + { + var times = 0; + while (polynomial.Count > 1) + { + var quotient = new Entity[polynomial.Count - 1]; + var carried = polynomial[^1]; + for (var k = polynomial.Count - 2; k >= 0; k--) + { + quotient[k] = carried; + carried = (polynomial[k] + root * carried).InnerSimplified; + } + if (!VanishesIdentically(carried)) + break; + polynomial.Clear(); + polynomial.AddRange(quotient); + times++; + } + return times; + } + } + /// /// A sum every term of which carries the same power of one linear form in /// , written as that power times the sum of the rest: diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index b771c2129..8a8e906d8 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -1134,6 +1134,8 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) // And a power of `a + b cos(y) + c sin(y)` with `a^2 = b^2 + c^2`, the same square // turned by a phase, in closed form: the substitution search spent the budget on it. if ((answer = IndefiniteIntegralSolver.SolveAPowerOfACosineAndASinePlusTheirAmplitude(expr, x)) is { }) return answer; + // And `a + b cosh(y) + c sinh(y)` with `a^2 = b^2 - c^2`, the hyperbolic one. + if ((answer = IndefiniteIntegralSolver.SolveAPowerOfAHyperbolicCosineAndSinePlusTheirAmplitude(expr, x)) is { }) return answer; // And symbolic powers of `a ± a sin(y)` beside a power of `g cos(y)`, by u = sin(y), where // `(1 + u)(1 - u)` is the cosine's square: the half angle wants numeric powers, and the // substitution search spent the budget on these. Numeric ones too beside a fractional @@ -1158,6 +1160,20 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr) if ((answer = IndefiniteIntegralSolver.SolveAHalfOddPowerOfOnePlusAnImaginaryHyperbolicSine(expr, x, integrateByParts)) is { }) return answer; // And the cosecant's: a ± i a csch(y) is a (w ± i)^2/(w^2 - 1) in w = e^y. if ((answer = IndefiniteIntegralSolver.SolveAPowerOfAnImaginaryHyperbolicCosecantSumInTheExponential(expr, x, integrateByParts)) is { }) return answer; + // And a polynomial in sinh(y) and cosh(y) over a ± i a sinh(y) beside a power of x, + // divided out, cosh(y)^2 being (1 + i sinh(y))(1 - i sinh(y)). + if ((answer = IndefiniteIntegralSolver.SolveByDividingOutOnePlusAnImaginaryHyperbolicSine(expr, x, integrateByParts)) is { }) return answer; + // A symbolic power of one hyperbolic function beside an odd power of the other, under + // u = the first: Bioche's substitution reads whole and rational powers. The first may be + // sech(y) or csch(y), and tanh(y) and coth(y) beside it are read as quotients. + if ((answer = IndefiniteIntegralSolver.SolveASymbolicPowerOfAHyperbolicFunctionBesideAnOddPowerOfTheOther(expr, x, integrateByParts)) is { }) return answer; + // And cosh(y)^2 taken out of a fractional power of a quotient, `cosh(x) coth(x)` being + // `cosh(x)^2/sinh(x)`: Bioche's substitution reads the radicands that are polynomials. + if ((answer = IndefiniteIntegralSolver.SolveByTakingTheHyperbolicCosineSquaredOutOfAPower(expr, x, integrateByParts)) is { }) return answer; + // And cosh(y)^(2m) over (1 ± i sinh(y))^(2m), or the other way up, cancelled through + // cosh(y)^2 = (1 + i sinh(y))(1 - i sinh(y)): `1/(sech(x) + i tanh(x))^4` is + // `(1 - i sinh(x))^2/(1 + i sinh(x))^2`, where the substitutions spent the budget. + if ((answer = IndefiniteIntegralSolver.SolveByCancellingOnePlusAnImaginaryHyperbolicSineAgainstTheCosineSquared(expr, x, integrateByParts)) is { }) return answer; // A quotient of linears as a function's argument, written over its denominator: a // constant plus a multiple of the reciprocal of a linear, which the rules read. if ((answer = IndefiniteIntegralSolver.SolveByWritingAQuotientOfLinearsOverItsDenominator(expr, x, integrateByParts)) is { }) return answer; diff --git a/Sources/Tests/UnitTests/Calculus/APowerOfAHyperbolicCosineAndSinePlusTheirAmplitudeIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/APowerOfAHyperbolicCosineAndSinePlusTheirAmplitudeIntegralTest.cs new file mode 100644 index 000000000..f041b93d6 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/APowerOfAHyperbolicCosineAndSinePlusTheirAmplitudeIntegralTest.cs @@ -0,0 +1,94 @@ +// +// 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 half-odd power of a + b cosh(y) + c sinh(y) with a^2 = b^2 - c^2, a square of + /// the half angle turned by a phase, in closed form through T = b sinh(y) + c cosh(y). + /// Rubi's 6.7.1. + /// #233 + /// + /// + /// Compared as complex numbers: with b < |c| the constant is imaginary and so is the + /// base everywhere. Away from the zero of the base where a is negative. + /// + [Trait("Area", "Calculus")] + public sealed class APowerOfAHyperbolicCosineAndSinePlusTheirAmplitudeIntegralTest + { + [Theory] + [InlineData("(b*cosh(x) + c*sinh(x) + sqrt(b^2 - c^2))^(5/2)", 2.0, 1.0)] + [InlineData("sqrt(b*cosh(x) + c*sinh(x) - sqrt(b^2 - c^2))", 2.0, 1.0)] + [InlineData("1/sqrt(b*cosh(x) + c*sinh(x) + sqrt(b^2 - c^2))", 2.0, 1.0)] + [InlineData("1/(b*cosh(x) + c*sinh(x) - sqrt(b^2 - c^2))^(5/2)", 2.0, -1.0)] + [InlineData("1/(b*cosh(g + f*x) + c*sinh(g + f*x) + sqrt(b^2 - c^2))^(3/2)", 0.7, 1.1)] + [InlineData("(b*cosh(g + f*x) + c*sinh(g + f*x) - sqrt(b^2 - c^2))^(3/2)", 0.7, 1.1)] + [InlineData("1/sqrt(5*cosh(2*x + 1) + 3*sinh(2*x + 1) + 4)", 0.0, 0.0)] + [InlineData("3*(2*cosh(x) - sinh(x) - sqrt(3))^(-3/2)", 0.0, 0.0)] + public void IsWrittenThroughTheDerivativeOfItsBase(string integrand, double b, double c) + => DifferentiatesBack(integrand, b, c); + + // A whole power is left to the exponential substitution, which writes `1/S^2`, `1/S^3` and + // `1/S^4` shorter than the same recurrence stepped up to `T/(a S)`: answered all the same. + [Theory] + [InlineData("1/(b*cosh(x) + c*sinh(x) + sqrt(b^2 - c^2))^2", 2.0, 1.0)] + [InlineData("1/(b*cosh(x) + c*sinh(x) + sqrt(b^2 - c^2))^3", 2.0, -1.0)] + [InlineData("1/(b*cosh(x) + c*sinh(x) + sqrt(b^2 - c^2))^4", 0.7, 1.1)] + [InlineData("1/(b*cosh(g + f*x) + c*sinh(g + f*x) - sqrt(b^2 - c^2))^3", 2.0, 1.0)] + [InlineData("3*(2*cosh(x) - sinh(x) - sqrt(3))^(-2)", 0.0, 0.0)] + [InlineData("1/(5*cosh(2*x + 1) + 3*sinh(2*x + 1) + 4)^3", 0.0, 0.0)] + public void AWholePowerIsAnsweredToo(string integrand, double b, double c) + => DifferentiatesBack(integrand, b, c); + + private static void DifferentiatesBack(string integrand, double b, double c) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Entity Pinned(Entity e) => e.Substitute("b", b).Substitute("c", c).Substitute("g", 0.4).Substitute("f", 1.3); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + foreach (var at in new[] { -1.9, -1.1, -0.2, 0.4, 1.1, 1.9 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + 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}"); + } + } + + /// + /// Real and continuous where the integrand is: across the point where T vanishes and + /// the base does not, and on either side of the base's zero where a is negative. + /// The difference of the answer against Simpson's rule. + /// + [Theory] + [InlineData(1.0, -2.0, 1.0)] + [InlineData(-1.0, 0.1, 2.0)] + [InlineData(-1.0, -3.0, -1.0)] + public void IsAContinuousAntiderivativeOnTheRealLine(double sign, double from, double to) + { + var root = sign > 0 ? "sqrt(3)" : "-sqrt(3)"; + var integral = $"1/sqrt(2*cosh(x) + sinh(x) + {root})".ToEntity().Integrate("x").Substitute("C", 0); + Assert.DoesNotContain("integral(", integral.Stringize()); + double Integrand(double t) => 1 / Math.Sqrt(2 * Math.Cosh(t) + Math.Sinh(t) + sign * Math.Sqrt(3)); + const int steps = 400; + var h = (to - from) / steps; + var simpson = Integrand(from) + Integrand(to); + for (var k = 1; k < steps; k++) + simpson += (k % 2 == 1 ? 4 : 2) * Integrand(from + k * h); + simpson *= h / 3; + var difference = integral.Substitute("x", to).EvalNumerical() - integral.Substitute("x", from).EvalNumerical(); + Assert.True(Math.Abs((double)difference.ImaginaryPart) < 1e-9, $"F({to}) - F({from}) is {difference}"); + Assert.True(Math.Abs((double)difference.RealPart - simpson) < 1e-7, $"F({to}) - F({from}) is {difference}, Simpson's rule {simpson}"); + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/ASymbolicPowerOfAHyperbolicFunctionBesideAnOddPowerOfTheOtherIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ASymbolicPowerOfAHyperbolicFunctionBesideAnOddPowerOfTheOtherIntegralTest.cs new file mode 100644 index 000000000..f22f0bb39 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ASymbolicPowerOfAHyperbolicFunctionBesideAnOddPowerOfTheOtherIntegralTest.cs @@ -0,0 +1,80 @@ +// +// Copyright (c) 2019-2026 Angouri. +// AngouriMath is licensed under MIT. +// Details: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md. +// Website: https://am.angouri.org. +// + +using System; +using AngouriMath.Extensions; +using Xunit; + +namespace AngouriMath.Tests.Calculus +{ + /// + /// A symbolic power of sinh(y) beside an odd power of cosh(y), under + /// u = sinh(y), and the other way about under u = cosh(y). Rubi's 6.7.1. + /// #233 + /// + /// Compared as complex numbers: a symbolic power of a negative sine is complex. + [Trait("Area", "Calculus")] + public sealed class ASymbolicPowerOfAHyperbolicFunctionBesideAnOddPowerOfTheOtherIntegralTest + { + [Theory] + [InlineData("cosh(a + b*x)^3*sinh(a + b*x)^n")] + [InlineData("cosh(a + b*x)^5*sinh(a + b*x)^n")] + [InlineData("cosh(a + b*x)^m*sinh(a + b*x)^3")] + [InlineData("cosh(a + b*x)^m*sinh(a + b*x)^5")] + [InlineData("2*cosh(x)^3/sinh(x)^n")] + public void IsAPolynomialTimesAPowerInTheOther(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Entity Pinned(Entity e) => e.Substitute("a", 0.4).Substitute("b", 1.3).Substitute("n", 0.37).Substitute("m", -1.6); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + foreach (var at in new[] { -1.9, -1.1, -0.7, 0.4, 1.1, 1.9 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + 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}"); + } + } + + /// + /// The reciprocals to the symbolic power and the tangent and cotangent beside them, read as + /// powers of the two: sech(y)^n tanh(y)^3 under u = sech(y), and + /// coth(y)^3 csch(y)^n under u = csch(y), whose power is complex where the sine + /// is negative and not that of sinh(y)^(-n). Rubi's 6.7.1. Both signs of the + /// argument at each slope. + /// + [Theory] + [InlineData("sech(a + b*x)^n*tanh(a + b*x)^3", 0.4, 1.3)] + [InlineData("coth(a + b*x)^3*csch(a + b*x)^n", 0.4, 1.3)] + [InlineData("coth(a + b*x)^3*csch(a + b*x)^n", 0.9, -0.8)] + [InlineData("sech(a + b*x)^m*tanh(a + b*x)^5", -0.4, 1.3)] + [InlineData("coth(a + b*x)^5*csch(a + b*x)^m", 0.9, -0.8)] + [InlineData("cosh(a + b*x)^n*tanh(a + b*x)^3", 0.4, 1.3)] + [InlineData("sinh(a + b*x)^m*coth(a + b*x)^3", -0.4, 1.3)] + [InlineData("csch(a + b*x)^n*cosh(a + b*x)^3", 0.4, 1.3)] + [InlineData("2*sech(x)^n*sinh(x)^3", 0.0, 1.0)] + public void ReadsTheReciprocalsAndTheQuotients(string integrand, double a, double b) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Entity Pinned(Entity e) => e.Substitute("a", a).Substitute("b", b).Substitute("n", 0.37).Substitute("m", -1.6); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + foreach (var at in new[] { -1.9, -1.1, -0.7, 0.4, 1.1, 1.9 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + 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}"); + } + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/CancellingOnePlusAnImaginaryHyperbolicSineAgainstTheCosineSquaredIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/CancellingOnePlusAnImaginaryHyperbolicSineAgainstTheCosineSquaredIntegralTest.cs new file mode 100644 index 000000000..4fc93b2e9 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/CancellingOnePlusAnImaginaryHyperbolicSineAgainstTheCosineSquaredIntegralTest.cs @@ -0,0 +1,50 @@ +// +// 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 +{ + /// + /// cosh(y)^(2m) over (1 ± i sinh(y))^(2m), or the other way up, cancelled through + /// cosh(y)^2 = (1 + i sinh(y))(1 - i sinh(y)): 1/(sech(x) + i tanh(x))^4 is + /// (1 - i sinh(x))^2/(1 + i sinh(x))^2. Rubi's 6.7.1. + /// #233 + /// + /// Compared as complex numbers: the integrand is complex wherever it is defined. + [Trait("Area", "Calculus")] + public sealed class CancellingOnePlusAnImaginaryHyperbolicSineAgainstTheCosineSquaredIntegralTest + { + [Theory] + [InlineData("1/(sech(x) + i*tanh(x))^4")] + [InlineData("1/(sech(x) - i*tanh(x))^4")] + [InlineData("1/(sech(x) + i*tanh(x))^6")] + [InlineData("1/(sech(a + b*x) - i*tanh(a + b*x))^4")] + [InlineData("cosh(x)^4/(1 + i*sinh(x))^4")] + [InlineData("3*cosh(a + b*x)^6/(a - i*a*sinh(a + b*x))^6")] + [InlineData("(sech(x) + i*tanh(x))^6")] + [InlineData("cosh(x)^2/(i + sinh(x))^2")] + public void IsTheQuotientTheCancellationLeaves(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Entity Pinned(Entity e) => e.Substitute("a", -0.4).Substitute("b", 1.3); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + foreach (var at in new[] { -1.9, -1.1, -0.2, 0.4, 1.1, 1.9 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + 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}"); + } + } + } +} diff --git a/Sources/Tests/UnitTests/Calculus/DividingOutOnePlusAnImaginaryHyperbolicSineIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/DividingOutOnePlusAnImaginaryHyperbolicSineIntegralTest.cs new file mode 100644 index 000000000..f8570286c --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/DividingOutOnePlusAnImaginaryHyperbolicSineIntegralTest.cs @@ -0,0 +1,54 @@ +// +// 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 in sinh(y) and cosh(y) over a ± i a sinh(y) beside a power + /// of a linear, divided out: cosh(y)^2 = (1 + i sinh(y))(1 - i sinh(y)). Rubi's 6.1.1. + /// #233 + /// + /// Compared as complex numbers: the integrand is complex wherever it is defined. + [Trait("Area", "Calculus")] + public sealed class DividingOutOnePlusAnImaginaryHyperbolicSineIntegralTest + { + [Theory] + [InlineData("(g + f*x)*cosh(c + d*x)^2/(a + i*a*sinh(c + d*x))")] + [InlineData("(g + f*x)^3*cosh(c + d*x)^3/(a + i*a*sinh(c + d*x))")] + [InlineData("cosh(c + d*x)^2/((g + f*x)^2*(a + i*a*sinh(c + d*x)))")] + [InlineData("cosh(c + d*x)^3/((g + f*x)*(a + i*a*sinh(c + d*x)))")] + [InlineData("(g + f*x)*sinh(c + d*x)^3/(a + i*a*sinh(c + d*x))")] + [InlineData("(g + f*x)^2*cosh(c + d*x)^3/(a - i*a*sinh(c + d*x))")] + [InlineData("x^2*cosh(2*x + 1)^4/(3 - 3*i*sinh(2*x + 1))")] + public void IsAPolynomialInTheHyperbolicFunctionsAndAQuotient(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Entity Pinned(Entity e) => e.Substitute("a", 1.3).Substitute("c", 0.4).Substitute("d", 0.3).Substitute("g", 0.6).Substitute("f", 1.1); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { -1.9, -1.1, -0.2, 0.4, 1.1, 1.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/TakingTheHyperbolicCosineSquaredOutOfAPowerIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/TakingTheHyperbolicCosineSquaredOutOfAPowerIntegralTest.cs new file mode 100644 index 000000000..91e81bd79 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/TakingTheHyperbolicCosineSquaredOutOfAPowerIntegralTest.cs @@ -0,0 +1,50 @@ +// +// 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 +{ + /// + /// cosh(y)^2 taken out of a fractional power of a quotient: cosh(x) coth(x) and + /// csch(x) + sinh(x) are both cosh(x)^2/sinh(x), whose power is then odd in the + /// hyperbolic cosine. Rubi's 6.7.1. + /// #233 + /// + /// + /// Compared as complex numbers on both sides of zero: for x < 0 the radicand is + /// negative, and the constant the power was written apart with is -1 there. + /// + [Trait("Area", "Calculus")] + public sealed class TakingTheHyperbolicCosineSquaredOutOfAPowerIntegralTest + { + [Theory] + [InlineData("(cosh(x)*coth(x))^(1/2)")] + [InlineData("(cosh(x)*coth(x))^(5/2)")] + [InlineData("(csch(x) + sinh(x))^(3/2)")] + [InlineData("1/(csch(x) + sinh(x))^(1/2)")] + [InlineData("(cosh(a + b*x)^2/(c + d*sinh(a + b*x)))^(3/2)")] + public void IsOddInTheCosineOnceItIsOut(string integrand) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Entity Pinned(Entity e) => e.Substitute("a", 0.4).Substitute("b", 1.3).Substitute("c", 0.7).Substitute("d", 1.1); + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x"); + var original = Pinned(integrand.ToEntity()); + foreach (var at in new[] { -1.9, -1.1, -0.2, 0.4, 1.1, 1.9 }) + { + var want = original.Substitute("x", at).EvalNumerical(); + var got = derivative.Substitute("x", at).EvalNumerical(); + 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}"); + } + } + } +}