diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 7b8cb4a18..0cb4c08b6 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -83,6 +83,24 @@ powers of `c + d x` itself were read before ([#233](https://github.com/asc-commu | `"e^(k/(c+d*x))/(a+b*x)".ToEntity().Integrate("x")` | `integral(...)` | two `Ei`, of `k/(c + d x)` and of `k/(c + d x) - k b/(a d - b c)` | | `"f^(c/(a+b*x)^2)*x^4".ToEntity().Integrate("x")` | `integral(...)` | 1,262 characters, the error function and powers of `1/(a + b x)` | +### A transcendental function of one shifted linear is integrated in the linear + +**Answers where there were none.** Rubi's 2.3 `F^(a + b/(c + d x)^3)/(c + d x)^19` and 3.4 `ln(c (d + e/(f + g x)^3)^q)` +ran past the budget, and 3.3's `(e g + f g x)^(p - 1)/ln(d (e + f x)^p)` was declined. Where every x stands in one linear +`c + d x`, up to constant multiples, and either a power of it stands below the bar inside a sum or it is the base of powers +written two ways, the integrand is now written in `u = c + d x` before the polynomial term takes the slope out of a power +of it: that wrote `(c/d + x)^19` beside the `(c + d x)^3` in the exponent, and the exponential ansatz on the two spellings +took longer the higher the power, 24 s at the 16th. In u it is 20 ms. Rows of these shapes that were answered are answered +with the slope in another place, `F^(a + b/(c + d x)^2)/(c + d x)^2` in 15 ms where it took 1.6 s +([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"F^(a+b/(c+d*x)^3)/(c+d*x)^19".ToEntity().Integrate("x")` | `integral(...)` | 386 characters, `F^(a + b/(c + d x)^3)` times a polynomial in `(c + d x)^3` | +| `"ln(c*(d+k/(f+g*x)^3)^q)".ToEntity().Integrate("x")` | `integral(...)` | 626 characters, by parts, the cube roots of `-k/d` | +| `"(k*g+f*g*x)^(p-1)/ln(d*(k+f*x)^p)".ToEntity().Integrate("x")` | `integral(...)` | 150 characters, `Ei(ln(d (k + f x)^p))` | +| `"F^(a+b/(c+d*x)^2)/(c+d*x)^2".ToEntity().Integrate("x")` | `integral(...)` | `-F^a sqrt(pi) erf(sqrt(-b ln F)/(c + d x))/(2 d sqrt(-b ln F))` | + ### A polynomial over a symbolic quartic that is two quadratics in a hidden variable is split over them **Answers where there were none.** `1/(a + b x + c x^2 + b x^3 + a x^4)` was declined, and so were Rubi's 1.3.1 @@ -216,6 +234,22 @@ it is still declined ([#233](https://github.com/asc-community/AngouriMath/issues | `"e^(n*acoth(a*x))*x/(c-a^2*c*x^2)^(3/2)".ToEntity().Integrate("x")` | `integral(...)` | 155 characters | | `"e^(n*acoth(a*x))*x^3/(c-a^2*c*x^2)^(5/2)".ToEntity().Integrate("x")` | `integral(...)` | 397 characters | +### A product of powers whose logarithmic derivative is a constant is its own derivative's multiple + +**Answers where there were none.** Rubi's 2.3 `(a F^(c + d x))^m (b F^(e + f x))^n` and `x^(-2 - 1/x) (1 - ln(x))`, and 6.7.1's +`cosh(x)^x (ln(cosh(x)) + x tanh(x))` and `F^(a + b x) (cosh(c + d x) - sinh(c + d x))^n`, were declined. The derivative of +a product of powers was read with some of them raised by one, and never with none: where the product's own logarithmic +derivative is a constant, the product over that constant is the answer. A power with x in its base and its exponent is an +exponential, `f^g = e^(g ln f)`, its terms free of x apart, and where the rest of the integrand is a constant times the +exponentials' logarithmic derivative the exponentials are the answer: `x^(-2 - 1/x) (1 - ln(x))` is `-(x^(-1/x))'` +([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"(a*F^(c+d*x))^m*(b*F^(k+f*x))^n".ToEntity().Integrate("x")` | `integral(...)` | `(a F^(c + d x))^m (b F^(k + f x))^n/((d m + f n) ln(F))` | +| `"x^(-2-1/x)*(1-ln(x))".ToEntity().Integrate("x")` | `integral(...)` | `-x^(-1/x)` | +| `"cosh(x)^x*(ln(cosh(x))+x*tanh(x))".ToEntity().Integrate("x")` | `integral(...)` | `((e^x + e^(-x))/2)^x`, which is `cosh(x)^x` | + ### A power of a quotient of linears beside whole powers of them is integrated under the quotient **Answers where there were none.** `e^(n acoth(a x))/(c - a^2 c x^2)^2` was declined, and so were its `^3` and `^4`, Rubi's @@ -245,6 +279,26 @@ found by one linear solve whose last equation is checked; where it fails the rul | `"e^(1/2*atanh(a*x))*x^2/(c-a^2*c*x^2)^(9/8)".ToEntity().Integrate("x")` | `integral(...)` | 121 characters | | `"e^(2*atanh(a*x))*x^2/(c-a^2*c*x^2)^3".ToEntity().Integrate("x")` | `integral(...)` | 109 characters | +### Powers of two linears beside a negative power of a logarithm of their quotient are integrated onto the exponential integral + +**Answers where there were none.** Rubi's 3.2.1 `1/((a g + b g x)^2 (A + B ln(e (a + b x)/(c + d x)))^2)` and three more +ran past the budget, and the rest of its rows and 3.2.2's with the powers of the linears summing to `-2` or `-3` took from +0.8 s to past 8 s: twenty of them past 8 s on a loaded machine, now under a quarter of a second. Under `t = (a + b x)/(c + d x)` each is a power of `t` and a polynomial beside a negative power of +`A + B ln(e t^n)`, which the rule for a power over a power of a logarithm answers onto `Ei`; that substitution was the +rule's already, far down the chain, and the polynomial term reached first took `b^2 g^2` out of `(a g + b g x)^2` and +asked again, past the budget. It is now tried first for these shapes. And `1/(t^2 (A + B ln(k t))^2)`, the question it +asks, had an answer with no value at any point: under `u = ln(k t)` the square `B^2 u^2 + 2 A B u + A^2` was split over +its two roots, which are one, and the partial fractions divided by their difference. A quadratic below the bar whose +discriminant is zero is now one linear twice ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"1/((a*g+b*g*x)^2*(A+B*ln(k*(a+b*x)/(c+d*x)))^2)".ToEntity().Integrate("x")` | no answer in twenty seconds | 340 characters, `Ei(-(A + B ln(k (a + b x)/(c + d x)))/B)` and the elementary term | +| `"1/((a*g+b*g*x)^2*(A+B*ln(k*((a+b*x)/(c+d*x))^n))^2)".ToEntity().Integrate("x")` | `integral(...)` | 394 characters | +| `"1/((a*g+b*g*x)^2*(A+B*ln(k*(c+d*x)^2/(a+b*x)^2)))".ToEntity().Integrate("x")` | no answer in twenty seconds | 231 characters, one `Ei` | +| `"1/((a*g+b*g*x)^3*(A+B*ln(k*(a+b*x)/(c+d*x)))^2)".ToEntity().Integrate("x")` | no answer in twenty seconds | 644 characters, two powers of t | +| `"1/(x^2*(A+B*ln(k*x))^2)".ToEntity().Integrate("x")` | `integral(...)` | 133 characters, `Ei(-(A/B + ln(k x)))`, where before this change master gave an answer with no value at any point | + ### The last digits of an evaluated sine, cosine, tangent, arcsine or arctangent **Silent, and in the hundredth digit.** `EvalNumerical` at the default hundred digits gives the @@ -323,12 +377,29 @@ no order. Five integrals whose working squares such a sum were past the budget a | Input | Was (2.5.0) | Now | |---|---|---| -| `"1/((a*g+b*g*x)^2*(A+B*ln(k*((a+b*x)/(c+d*x))^n))^2)".ToEntity().Integrate("x")` | `integral(...)` | 412 characters | -| `"1/((a*g+b*g*x)^2*(A+B*ln(k*(a+b*x)/(c+d*x)))^2)".ToEntity().Integrate("x")` | past 60 s | 358 characters | +| `"1/((a*g+b*g*x)^2*(A+B*ln(k*((a+b*x)/(c+d*x))^n))^2)".ToEntity().Integrate("x")` | `integral(...)` | 394 characters | +| `"1/((a*g+b*g*x)^2*(A+B*ln(k*(a+b*x)/(c+d*x)))^2)".ToEntity().Integrate("x")` | past 60 s | 340 characters | | `"(a+a*sin(k+f*x))^(3/2)*(A+B*sin(k+f*x))/(c+d*sin(k+f*x))^2".ToEntity().Integrate("x")` | `integral(...)` | 5,373 characters | | `"sech(x)^2*tanh(x)^2/(2+tanh(x)^3)^2".ToEntity().Integrate("x")` | past 60 s | 117 characters | | `"coth(x)^3*csch(x)*sqrt(1+csch(x))".ToEntity().Integrate("x")` | `integral(...)` | 325 characters | +### The derivative of a function over a power of a logarithm, an inverse tangent or a sum with a logarithm in it is read as one + +**Answers where there were none.** Rubi's 3.2.3 `1/((c + d x)(c - a + (d - b) x) ln(Q)) + ln(1 - Q)/((a + b x)(c + d x) ln(Q)^2)` +for `Q = (a + b x)/(c + d x)`, its twin, 5.3.4's `x^3/((1 + a^2 x^2) atan(a x)^3) - 3 x^2/(2 a atan(a x)^2)` and 3.5's +`(a d n x^m - a d m x^m ln(c x^n) - b d n (q - 1) ln(c x^n)^q)/(x (a x^m + b ln(c x^n)^q)^2)` ran past the budget. Each is +`(U/V^k)' = U'/V^k - k U V'/V^(k + 1)`, whose terms are not elementary apart, so every split lost it and by parts on either +term asked the other's question again. For `V` a logarithm or an inverse trigonometric function, the terms over +`V^(k + 1)` give `U = -B/(k V')`, its rational part in lowest terms; for `V` a sum with a logarithm in it, `U` is a +monomial of the numerator times a constant read monomial by monomial. The answer is given only where its derivative is +the whole integrand at sampled points ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"1/((c+d*x)*(-a+c+(-b+d)*x)*ln((a+b*x)/(c+d*x)))+ln(1+(-a-b*x)/(c+d*x))/((a+b*x)*(c+d*x)*ln((a+b*x)/(c+d*x))^2)".ToEntity().Integrate("x")` | no answer in twenty seconds | `ln(1 + (-a - b x)/(c + d x))/((a d - b c) ln((a + b x)/(c + d x)))` | +| `"x^3/((1+a^2*x^2)*atan(a*x)^3)-3/2*x^2/(a*atan(a*x)^2)".ToEntity().Integrate("x")` | `integral(...)` | `-x^3/(2 a atan(a x)^2)` | +| `"(a*d*n*x^m-a*d*m*x^m*ln(c*x^n)-b*d*n*(q-1)*ln(c*x^n)^q)/(x*(a*x^m+b*ln(c*x^n)^q)^2)".ToEntity().Integrate("x")` | `integral(...)` | `d ln(x^n c)/(a x^m + b ln(c x^n)^q)` | + ### A rational evaluated at a higher precision than before has all of its digits **More digits where there were a hundred.** A rational kept its decimal value in the precision current when it was @@ -373,6 +444,23 @@ parts what is left is a rational function and one over the root |---|---|---| | `"x^3*ln(-1+4*x+4*sqrt((x-1)*x))".ToEntity().Integrate("x")` | `integral(...)` | `x^4 ln(...)/4` minus a polynomial, `ln(8 x + 1)` and terms over `sqrt(x^2 - x)` | +### The Risch-Norman ansatz reads `ln(x)` beside x, with the logarithmic integral and the error function among its terms + +**Answers where there were none.** Rubi's 3.5 `(1 - ln(x))/(x (x + ln(x)))` and `(1 + x)/(ln(x) (x + ln(x)))`, and 2.3's +`e^((a + x)^2)/x^2 - 2 a e^((a + x)^2)/x`, were declined. The ansatz left a plain `ln(x)` to the logarithm tower's, a +polynomial in it, which has no logarithm of `x + ln(x)`; it now reads one where a sum below the bar holds x beside +`ln(x)` to two or more powers. Its answer may now hold `li(x)`, whose derivative is `1/ln(x)`, where `ln(x)` divides the +denominator, and for an exponential of a quadratic `erfi(sqrt(A) (x + B/(2 A)))`, whose derivative is a constant times +the exponential: `(1 + x)/(ln(x) (x + ln(x)))` is `(li(x) + ln(ln(x)) - ln(x + ln(x)))'`, and by parts the first term +of the third leaves the second's negative and `2 e^((a + x)^2)`. With either among the terms the elimination pins the +symbols first where it declines with them in ([#233](https://github.com/asc-community/AngouriMath/issues/233)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"(1-ln(x))/(x*(x+ln(x)))".ToEntity().Integrate("x")` | `integral(...)` | `ln(x + ln(x)) - ln(x)` | +| `"(1+x)/(ln(x)*(x+ln(x)))".ToEntity().Integrate("x")` | `integral(...)` | `li(x) + ln(ln(x)) - ln(x + ln(x))` | +| `"e^((a+x)^2)/x^2-2*e^((a+x)^2)*a/x".ToEntity().Integrate("x")` | `integral(...)` | `-e^((a + x)^2)/x + sqrt(pi) erfi(x + a)` | + ### The last digits of an evaluated logarithm, exponential, real power or complex power **Silent, and in the hundredth digit.** `ln`, `log`, `e^x`, a real base to a real power and the diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 5a953b156..2aabecbf1 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -17575,6 +17575,108 @@ node is Powf(var @base, var power) && @base.ContainsNode(u) _ => false, }; + /// The largest integrand reads; the shapes it answers are 17 to 25. + private const int MostComplexityWrittenInOneShiftedLinear = 80; + + /// + /// A function with something other than sums, products and powers in it, of one shifted + /// linear c + d x alone, written in u = c + d x: Rubi's + /// F^(a + b/(c + d x)^3)/(c + d x)^19 is F^(a + b/u^3)/(d u^19), + /// ln(c (d + e/(f + g x)^3)^q) is ln(c (d + e/u^3)^q)/g, and + /// (e g + f g x)^(p - 1)/ln(d (e + f x)^p) is (g u)^(p - 1)/(f ln(d u^p)). + /// + /// + /// In x, the polynomial term takes the slope out of a power of the linear and writes + /// (c/d + x)^19 beside the (c + d x)^3 in the exponent, two spellings of one + /// linear, and the exponential ansatz on that took longer the higher the power -- 24 s at + /// (c + d x)^16, past the budget at 19 -- for what is 20 ms in u. By parts on the + /// logarithm leaves a rational function in x of degree six with the symbols in it, where in u + /// it is 1/(d u^3 + e). Every x stands in a sum linear in it, each sum a constant + /// multiple of one of them, and either a power of the linear stands below the bar inside a + /// sum or the linear is written two ways: a function of the linear otherwise -- an + /// exponential of it, ln(c (d + e x)^n)^p -- is read as written by the rules for it, + /// and was left so. The algebraic integrands are the rule above's. At the question asked only: + /// a level down, the polynomial term's own monic spelling `(c/d + x)` beside `c + d x` is two + /// spellings, and the substitution rewrote answers the chain had given as written. In u the + /// linear is u itself, and this declines. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveByWritingATranscendentalFunctionOfOneShiftedLinear(Entity expr, Entity.Variable x, bool integrateByParts) + { + // A written integrand, not what a substitution expands one into: the inverse sine's of + // `(d + e x^2)^3 (a + b asinh(c x))^2` is a sum of exponentials of linears in u at the same + // level, and reading every sum in it for one was 1.4 s and 370 MB of a decline. + if (!Integration.AnsweringTheQuestionAsked || expr.Complexity > MostComplexityWrittenInOneShiftedLinear) + return null; + if (!expr.Nodes.Any(node => node.ContainsNode(x) && node is not (Variable or Sumf or Minusf or Mulf or Divf) + && !(node is Powf(_, var exponent) && !exponent.ContainsNode(x)))) + return null; + // The sums linear in x, outermost: every x must stand in one. + var linears = new List<(Entity Linear, Entity Slope, Entity Offset)>(); + var asBases = new List(); + var belowTheBarInASum = false; + bool IsLinear(Entity node) => node is Sumf or Minusf && TreeAnalyzer.TryGetPolyLinear(node, x, out var slope, out var offset) + && !slope.ContainsNode(x) && !offset.ContainsNode(x) && !TreeAnalyzer.IsZero(slope); + bool Read(Entity node, bool insideASum) + { + if (!node.ContainsNode(x)) + return true; + if (IsLinear(node)) + { + if (!linears.Any(entry => entry.Linear == node)) + { + TreeAnalyzer.TryGetPolyLinear(node, x, out var slope, out var offset); + linears.Add((node, slope!, offset!)); + } + return true; + } + if (node == x) + return false; + if (insideASum && (node is Divf(_, var under) && (IsLinear(under) || under is Powf(var raised, var p) && !p.ContainsNode(x) && IsLinear(raised)) + || node is Powf(var @base, var power) && power.Evaled is Number.Real { IsNegative: true } && IsLinear(@base))) + belowTheBarInASum = true; + if (node is Powf(var powered, var exponentOfIt) && !exponentOfIt.ContainsNode(x) && IsLinear(powered) && !asBases.Contains(powered)) + asBases.Add(powered); + var inside = insideASum || node is Sumf or Minusf; + return node.DirectChildren.All(child => Read(child, inside)); + } + // A power of the linear below the bar inside a sum -- `b/(c + d x)^3` in the exponent, or + // in the logarithm -- which in x is a rational function of a degree the power multiplies, + // or two spellings of one linear. A function of the linear alone otherwise, `ln(c (d + e x)^n)^p`, + // is read as written by the rules for the logarithm of a linear. + // Two spellings: two linears, neither monic, each the base of a power. A monic one beside + // its multiple is the polynomial term's own spelling, asked as the same question, and a + // linear in an exponent, `f^(c + d x)` beside `f^(2 c + 2 d x)`, is the exponential rules'. + if (!Read(expr, insideASum: false) || linears.Count == 0) + return null; + var twoSpellings = asBases.Count >= 2 && asBases.All(linear => + TreeAnalyzer.TryGetPolyLinear(linear, x, out var slope, out _) && !VanishesIdentically(slope - 1)); + if (!belowTheBarInASum && !twoSpellings) + return null; + // The one in the fewest symbols is the variable, so that the others are its constant + // multiples as written: `e g + f g x` is `g u` for `u = e + f x`. + var (linear, slope0, offset0) = linears.OrderBy(entry => entry.Linear.Complexity).First(); + if (VanishesIdentically(offset0)) + return null; + var u = Variable.CreateUnique(expr, "u_linear"); + var inU = new Dictionary(); + foreach (var (other, slope, offset) in linears) + { + if (!VanishesIdentically(slope * offset0 - slope0 * offset)) + return null; + var ratio = LowestOverTheSymbols((slope / slope0).InnerSimplified); + inU[other] = ratio == Number.Integer.One ? u : ratio * u; + } + var rewritten = expr.Replace(node => inU.TryGetValue(node, out var written) ? written : node); + if (rewritten.ContainsNode(x)) + return null; + // dx = du/d. + if (Integration.ComputeAsAQuestionOfItsOwn((rewritten / slope0).InnerSimplified, u, integrateByParts) is not { } answer) + return null; + var back = answer.Substitute(u, linear); + return back.Nodes.Any(node => node == MathS.NaN) ? null : back; + } + /// /// A polynomial times powers of error functions of one shifted linear, written in /// u = a + b x: x erf(a + b x)^2 is (u - a) erf(u)^2/b^2, which by parts @@ -19255,6 +19357,7 @@ static void Join(List<(Number.Rational Power, Entity Cofactor)> groups, Number.R var raisable = new List<(Entity Base, Entity Exponent)>(); var logarithmicDerivativesKept = new List(); Entity kept = Number.Integer.One; + var varyingInBoth = false; foreach (var (written, writtenUnderneath) in FactorsOfTheIntegrand(expr)) { // A secant is a power of the cosine, and a cosecant of the sine, below the bar, @@ -19294,18 +19397,63 @@ static void Join(List<(Number.Rational Power, Entity Cofactor)> groups, Number.R kept = underneath ? kept / factor : kept * factor; continue; } - // A power with x in its base and its exponent is the exponential of the exponent - // times the base's logarithm, kept as the exponentials are, with that derivative: - // `cosh(x)^x (ln(cosh(x)) + x tanh(x))` is `(cosh(x)^x)'`, and raised it was squared. - if (factor is Powf(var variableBase, var variableExponent) && variableBase.ContainsNode(x) && variableExponent.ContainsNode(x)) + // And a power with x in its base and its exponent, `f^g = e^(g ln f)`, its terms + // free of x apart as a power to raise: `x^(-2 - 1/x)` is `x^(-2) x^(-1/x)`. + if (factor is Powf(var varyingBase, var varyingExponent) && varyingBase.ContainsNode(x)) { - var derivative = Functions.PartialFractions.Bare((variableExponent * MathS.Ln(variableBase)).Differentiate(x).InnerSimplified); + varyingInBoth = true; + Entity free = Number.Integer.Zero, varying = Number.Integer.Zero; + foreach (var term in Sumf.LinearChildren(varyingExponent)) + if (term.ContainsNode(x)) + varying += term; + else + free += term; + varying = varying.InnerSimplified; + free = free.InnerSimplified; + if (free.Evaled is not Number.Complex { IsZero: true }) + raisable.Add((varyingBase, underneath ? (-free).InnerSimplified : free)); + var exponential = MathS.Pow(varyingBase, varying); + var derivative = Functions.PartialFractions.Bare((varying * MathS.Ln(varyingBase)).Differentiate(x).InnerSimplified); logarithmicDerivativesKept.Add(underneath ? -derivative : derivative); - kept = underneath ? kept / factor : kept * factor; + kept = underneath ? kept / exponential : kept * exponential; continue; } raisable.Add((factor, underneath ? Number.Integer.MinusOne : Number.Integer.One)); } + // The exponentials alone, the rest of the integrand a constant multiple of their + // logarithmic derivative, where one of them is a power with x in its base and its exponent: + // `x^(-2 - 1/x) (1 - ln(x))` is `-(x^(-1/x))'`, and `cosh(x)^x (ln(cosh(x)) + x tanh(x))` is + // `(cosh(x)^x)'`. No power is raised, the rest is the derivative's, and raising x by one + // could not reach the `x^2` that `(-ln(x)/x)'` lowers. With constant bases only it is the + // substitution's, which answers it as written. + if (varyingInBoth && kept != Number.Integer.One && logarithmicDerivativesKept.Count > 0 && raisable.Count <= 4) + { + Entity rest = bracket ?? Number.Integer.One; + foreach (var (@base, exponent) in raisable) + rest *= MathS.Pow(@base, exponent); + Entity sum = Number.Integer.Zero; + foreach (var derivative in logarithmicDerivativesKept) + sum += derivative; + sum = Functions.PartialFractions.Bare(sum.InnerSimplified); + if (AreProportionalAtSampledPoints(rest, sum, x)) + { + var ratio = Functions.PartialFractions.Bare(Functions.SingleQuotient.Combine(rest / sum).Simplify()); + if (ratio is Providedf(var inner, _)) + ratio = inner; + // With no symbol in it a number, read at a point where the simplification left + // `tanh(x)` beside `sinh(x)/cosh(x)`: a small fraction, and checked by the derivative. + if (ratio.ContainsNode(x) && !(rest / sum).Vars.Any(v => v != x) + && (rest / sum).Substitute(x, Number.Rational.Create(7, 5)).Evaled is Number.Complex { ImaginaryPart.IsZero: true } at + && SmallFractionNear(at.RealPart.EDecimal) is { } near) + ratio = near; + if (!ratio.ContainsNode(x) && !ratio.Nodes.Any(node => node == MathS.NaN)) + { + var answer = (constant * ratio * kept).InnerSimplified; + if (Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x)) + return answer; + } + } + } // A function inside the bracket, raised from nothing: the bracket of `(K sin(u))'` is // `K'/K sin(u) + u' cos(u)`, so Rubi's // `F^(c (a + b x)) (f x)^m (e x cos(d + e x) + (m + b c x ln F) sin(d + e x))/x` is @@ -19325,9 +19473,20 @@ static void Join(List<(Number.Rational Power, Entity Cofactor)> groups, Number.R // `(a + b x^n + c x^(2 n))^p (a + b (1 + n + n p) x^n + c (1 + 2 n (1 + p)) x^(2 n))` // is `(x (a + b x^n + c x^(2 n))^(p + 1))'`. Whole powers only are a rational // function's, or a polynomial's beside a function, and the rules for those answer them. + // With no sum and no power raised, the product itself, where its logarithmic derivative is + // a constant: `(a F^(c + d x))^m (b F^(e + f x))^n` is `((a F^(c + d x))^m (b F^(e + f x))^n)'` + // over `(d m + f n) ln(F)`, and `F^(a + b x) (cosh(c + d x) + sinh(c + d x))^n` over + // `b ln(F) + d n`. Raising a power of either went nowhere, and the rules for one power of + // an exponential read one. + var itself = false; + var onlyItself = false; if (raisable.Any(pair => pair.Exponent is not Number.Integer)) { - if (raisable.Count >= 2) + // Beside an exponential, one power is the product itself only: raised, it was the + // substitution's question, answered as written. + if (bracket is null && (raisable.Count >= 2 || kept != Number.Integer.One)) + (itself, onlyItself) = (true, raisable.Count < 2); + if (raisable.Count >= 2 || kept != Number.Integer.One) bracket ??= Number.Integer.One; if (!raisable.Any(pair => pair.Base == x)) raisable.Insert(0, (x, Number.Integer.Zero)); @@ -19340,7 +19499,7 @@ static void Join(List<(Number.Rational Power, Entity Cofactor)> groups, Number.R var besideAFunction = new List<(Entity BracketOfTheDerivative, Entity Candidate)>(); // With a power of a variable base kept, the candidate may be the integrand's powers as // they are, nothing raised. - for (var subset = variablePowerKept ? 0 : 1; subset < 1 << raisable.Count; subset++) + for (var subset = variablePowerKept || itself ? 0 : 1; subset < (onlyItself ? 1 : 1 << raisable.Count); subset++) { var kernelsRaised = 0; var powersRaised = 0; @@ -19385,6 +19544,16 @@ static void Join(List<(Number.Rational Power, Entity Cofactor)> groups, Number.R candidate *= newExponent == Number.Integer.One ? @base : MathS.Pow(@base, newExponent); if (!raised) { + // In the product itself, a sum that is a constant multiple of its own derivative, + // `cosh(c + d x) + sinh(c + d x)`, has that constant for its logarithmic derivative, + // read monomial by monomial: as a quotient it was left as written, and nothing + // read a constant off it. + if (subset == 0 && @base is Sumf or Minusf && MonomialsOf(@base.Differentiate(x), x) is { } derived + && MonomialsOf(@base, x) is { } plain && TheRatioOfTheMonomials(derived, plain) is { } rate) + { + bracketOfTheDerivative += newExponent * rate * productOfRaisedBases; + continue; + } bracketOfTheDerivative += newExponent * @base.Differentiate(x) * productOfRaisedBases / @base; continue; } @@ -19604,6 +19773,178 @@ private static bool IsAFunctionToRaise(Entity factor, Entity.Variable x) return null; } + /// + /// The derivative of U/V^k for V a sum with a logarithm of x in it, standing to the + /// power k + 1 below the bar of a single quotient: Rubi's + /// (a d n x^m - a d m x^m ln(c x^n) - b d n (q - 1) ln(c x^n)^q)/(x (a x^m + b ln(c x^n)^q)^2) + /// is (d ln(c x^n)/(a x^m + b ln(c x^n)^q))'. U is a constant times one of the + /// monomials above the bar, tried in turn: the numerator W over the power must be + /// lambda (U' V - k U V'), decided at sampled points and then read monomial by monomial + /// with each power of x and of each function merged, so that ln(c x^n) ln(c x^n)^(q - 1) + /// is ln(c x^n)^q as the integrand writes it; or null. + /// + private static Entity? AsTheDerivativeOfAFunctionOverAPowerOfASum(Entity expr, Entity.Variable x) + { + if (expr is Sumf or Minusf) + return null; + var factors = FactorsOfTheIntegrand(expr); + foreach (var (factor, underneath) in factors) + { + // A sum with a logarithm of x in it: over a sum of exponentials or of trigonometric + // functions the rules for those answer as they did, the substitution among them. + if (!underneath || factor is not Powf(var sum, Number.Integer n) || sum is not (Sumf or Minusf) || !sum.ContainsNode(x) + || !sum.Nodes.Any(node => node is Logf(_, var argument) && argument.ContainsNode(x)) + || n.EInteger.CompareTo(EInteger.FromInt32(2)) < 0 || n.EInteger.CompareTo(EInteger.FromInt32(8)) > 0) + continue; + var k = n.EInteger.ToInt32Checked() - 1; + Entity above = Number.Integer.One, below = Number.Integer.One; + var skipped = false; + foreach (var (other, otherUnderneath) in factors) + { + if (!skipped && other == factor && otherUnderneath) + { + skipped = true; + continue; + } + if (otherUnderneath) below *= other; + else above *= other; + } + var w = above / below; + if (Monomials(w) is not { } read || read.Count == 0 || read.Count > 6) + continue; + // The candidates: each monomial above the bar with its constant dropped, and each of its + // functions of x on its own. Not 1: `V'/V^(k + 1)` is a power of a sum beside its + // derivative, which the rules for that answer as they always have. + var candidates = new List(); + foreach (var term in Sumf.LinearChildren(above).Concat(Sumf.LinearChildren(above.Expand()))) + { + Entity product = Number.Integer.One; + foreach (var part in Mulf.LinearChildren(term)) + if (part.ContainsNode(x)) + { + product = product == Number.Integer.One ? part : product * part; + if (!candidates.Contains(part)) + candidates.Add(part); + } + if (!candidates.Contains(product)) + candidates.Add(product); + } + var derivativeOfTheSum = sum.Differentiate(x); + foreach (var candidate in candidates.Take(12)) + { + if (candidate.ContainsNode(sum)) + continue; + var ofTheCandidate = candidate.Differentiate(x) * sum - k * candidate * derivativeOfTheSum; + if (!AreProportionalAtSampledPoints(w, ofTheCandidate, x) || Monomials(ofTheCandidate) is not { } theirs + || TheRatioOfTheMonomials(read, theirs) is not { } lambda) + continue; + var answer = lambda * candidate / (k == 1 ? sum : MathS.Pow(sum, k)); + if (Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x)) + return answer; + } + } + return null; + + List<(Entity Coefficient, Entity Power, List<(Entity Base, Entity Power)> Atoms)>? Monomials(Entity e) => MonomialsOf(e, x); + } + + /// + /// Each term of the expansion of as a coefficient, a power of + /// and powers of functions of it, each base's powers merged and like terms + /// summed; null where a factor is a sum holding x. + /// + private static List<(Entity Coefficient, Entity Power, List<(Entity Base, Entity Power)> Atoms)>? MonomialsOf(Entity e, Entity.Variable x) + { + var terms = new List<(Entity Coefficient, Entity Power, List<(Entity Base, Entity Power)> Atoms)>(); + // The conditions a derivative carries, `provided not x = 0`, are the generic case's. + e = e.Replace(node => node is Providedf(var inner, _) ? inner : node); + foreach (var term in Sumf.LinearChildren(e.Expand())) + { + var (numerator, denominator) = Functions.SingleQuotient.Of(term); + Entity coefficient = Number.Integer.One, power = Number.Integer.Zero; + var atoms = new List<(Entity Base, Entity Power)>(); + foreach (var (side, sign) in new[] { (numerator, 1), (denominator, -1) }) + foreach (var part in Mulf.LinearChildren(side)) + { + if (!part.ContainsNode(x)) + { + coefficient = sign > 0 ? coefficient * part : coefficient / part; + continue; + } + var (@base, exponent) = part is Powf(var b, var p) && !p.ContainsNode(x) ? (b, p) : (part, (Entity)Number.Integer.One); + if (@base is Sumf or Minusf) + return null; + if (sign < 0) + exponent = -exponent; + if (@base == x) + { + power += exponent; + continue; + } + var at = atoms.FindIndex(atom => atom.Base == @base); + if (at >= 0) + atoms[at] = (@base, atoms[at].Power + exponent); + else + atoms.Add((@base, exponent)); + } + power = Functions.PartialFractions.Bare(power.Simplify()); + atoms = atoms.Select(atom => (atom.Base, Functions.PartialFractions.Bare(atom.Power.Simplify()))) + .Where(atom => atom.Item2.Evaled is not Number.Complex { IsZero: true }).ToList(); + var like = terms.FindIndex(t => TheSameValue(t.Power, power) && t.Atoms.Count == atoms.Count + && t.Atoms.All(mine => atoms.Any(other => other.Base == mine.Base && TheSameValue(other.Power, mine.Power)))); + if (like >= 0) + terms[like] = (terms[like].Coefficient + coefficient, terms[like].Power, terms[like].Atoms); + else + terms.Add((coefficient, power, atoms)); + } + return terms.Where(t => !VanishesIdentically(t.Coefficient)).ToList(); + } + + /// + /// The one constant every monomial of is of the same monomial of + /// , or null: the two are then proportional, read exactly. + /// + private static Entity? TheRatioOfTheMonomials(List<(Entity Coefficient, Entity Power, List<(Entity Base, Entity Power)> Atoms)> mine, + List<(Entity Coefficient, Entity Power, List<(Entity Base, Entity Power)> Atoms)> theirs) + { + if (mine.Count != theirs.Count) + return null; + Entity? ratio = null; + foreach (var term in mine) + { + var matching = theirs.Where(t => TheSameValue(t.Power, term.Power) && t.Atoms.Count == term.Atoms.Count + && t.Atoms.All(atom => term.Atoms.Any(other => other.Base == atom.Base && TheSameValue(other.Power, atom.Power)))).ToList(); + if (matching.Count != 1) + return null; + var here = Functions.PartialFractions.Bare((term.Coefficient / matching[0].Coefficient).Simplify()); + if (here.Nodes.Any(node => node == MathS.NaN)) + return null; + if (ratio is null) + ratio = here; + else if (!VanishesIdentically(ratio - here)) + return null; + } + return ratio is null ? null : LowestOverTheSymbols(ratio); + } + + private static bool TheSameValue(Entity left, Entity right) => left == right || VanishesIdentically(left - right); + + /// + /// The fraction with a denominator of at most 24 within 1e-20 of , + /// or null. + /// + private static Number.Rational? SmallFractionNear(EDecimal value) + { + for (var denominator = 1; denominator <= 24; denominator++) + { + var scaled = value.Multiply(EDecimal.FromInt32(denominator)); + var nearest = scaled.RoundToExponent(0, ERounding.HalfEven); + if (scaled.Subtract(nearest).Abs().CompareTo(EDecimal.FromString("1e-20")) < 0) + return Number.Rational.Create(nearest.ToEInteger(), EInteger.FromInt32(denominator)); + } + return null; + } + /// Set while asks its question. [System.ThreadStatic] private static bool derivativeBesideAFunction; @@ -19816,6 +20157,128 @@ private static bool IsAFunctionToRaise(Entity factor, Entity.Variable x) /// The most powers raises. private const int MostBinomialPowersReduced = 16; + /// + /// The derivative of a quotient U/V^k for V a logarithm or an inverse + /// trigonometric function: (U/V^k)' = U'/V^k - k U V'/V^(k + 1), so the terms of the + /// integrand with V^(k + 1) below the bar are -k U V' over it, and + /// U = -B/(k V') for their sum B. Rubi's + /// x^3/((1 + a^2 x^2) atan(a x)^3) - 3 x^2/(2 a atan(a x)^2) is + /// (-x^3/(2 a atan(a x)^2))', and + /// 1/((c + d x)(c - a + (d - b) x) ln(Q)) + ln(1 - Q)/((a + b x)(c + d x) ln(Q)^2), for + /// Q = (a + b x)/(c + d x), is (-ln(1 - Q)/((b c - a d) ln(Q)))'. + /// + /// + /// Neither term is elementary on its own, so every split loses it, and by parts on either + /// leaves the other's question again; both were searches past the budget. The rational part + /// of B/(k V') is brought to lowest terms, where it is a constant or a polynomial for + /// these, and the answer is U/V^k only where its derivative is the whole integrand at + /// sampled points, the terms with lower powers of V included: where it is not, the + /// integrand is no such derivative, and this declines. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveAsTheDerivativeOfAFunctionOverAPowerOfALogarithmOrAnInverse(Entity expr, Entity.Variable x) + { + if (!Integration.AnsweringTheQuestionAskedOrOneBelow || expr is not (Sumf or Minusf or Mulf or Divf)) + return null; + static bool AnAtom(Entity node, Entity.Variable x) + => node.ContainsNode(x) && node is Logf or Arctanf or Arcsinf or Arccosf or Arccotanf; + var terms = Sumf.LinearChildren(expr).Select(FactorsOfTheIntegrand).ToList(); + // The atom and the highest power of it below the bar in any one term. + var candidates = new List<(Entity Atom, int Power)>(); + foreach (var factors in terms) + foreach (var (factor, underneath) in factors) + if (underneath && factor is Powf(var atom, Number.Integer n) && AnAtom(atom, x) + && n.EInteger.CompareTo(EInteger.FromInt32(2)) >= 0 && n.EInteger.CompareTo(EInteger.FromInt32(12)) <= 0) + { + var power = n.EInteger.ToInt32Checked(); + var at = candidates.FindIndex(pair => pair.Atom == atom); + if (at < 0) + candidates.Add((atom, power)); + else if (candidates[at].Power < power) + candidates[at] = (atom, power); + } + foreach (var (atom, highest) in candidates) + { + var derivative = atom.Differentiate(x); + if (derivative.ContainsNode(atom)) + continue; + // B and its factors that are no rational function of x, kept as they are. + Entity rational = Number.Integer.Zero; + Entity? transcendental = null; + var read = true; + foreach (var factors in terms) + { + var power = 0; + Entity above = Number.Integer.One, below = Number.Integer.One, others = Number.Integer.One; + foreach (var (factor, underneath) in factors) + { + var (@base, times) = factor is Powf(var b, Number.Integer n) && n.EInteger.CanFitInInt32() ? (b, n.EInteger.ToInt32Checked()) : (factor, 1); + if (@base == atom) + power += underneath ? times : -times; + else if (!factor.ContainsNode(x) || IsARationalFunction(factor, x)) + { + if (underneath) below *= factor; + else above *= factor; + } + else if (underneath || factor.ContainsNode(atom)) + read = false; + else + others *= factor; + } + if (!read) + break; + if (power != highest) + continue; + // One transcendental factor among the terms of B, the same in each: with two + // different ones the quotient is no product of one and a rational function. + if (others != Number.Integer.One) + { + if (transcendental is not null && transcendental != others) + { + read = false; + break; + } + transcendental = others; + } + rational += above / below; + } + if (!read) + continue; + if (rational == Number.Integer.Zero) + continue; + var k = highest - 1; + var quotient = InLowestTerms((-rational / (k * derivative)).InnerSimplified); + // A constant U is a power of the function beside its derivative, which the rules for + // that answer as they did: `1/(x ln(c x)^2)` is `-1/ln(c x)` there. + if (transcendental is null && !quotient.ContainsNode(x)) + continue; + var answer = (transcendental is null ? quotient : quotient * transcendental) / MathS.Pow(atom, k); + if (answer.Nodes.Any(node => node == MathS.NaN)) + continue; + if (Functions.PartialFractions.DerivativeHoldsAtSampledPoints(answer, expr, x)) + return answer; + } + return AsTheDerivativeOfAFunctionOverAPowerOfASum(expr, x); + + // The rational part in lowest terms: a constant -- `-1/(b c - a d)` for the logarithm of a + // quotient of linears, whose derivative over one bar is a quotient of quadratics -- as its + // value at a point, simplified, and otherwise cancelled over one bar. + Entity InLowestTerms(Entity r) + { + if (AreProportionalAtSampledPoints(r, Number.Integer.One, x)) + foreach (var at in new[] { 0, 1, 2 }) + { + var value = Functions.PartialFractions.Bare(r.Substitute(x, at).Simplify()); + if (!value.ContainsNode(x) && !value.Nodes.Any(node => node == MathS.NaN) && value.Evaled is not Number.Complex { IsZero: true }) + return LowestOverTheSymbols(value); + } + var combined = Functions.SingleQuotient.Combine(r).InnerSimplified; + var (numerator, denominator) = Functions.SingleQuotient.Of(combined); + return Functions.PolynomialGcd.TryCancel(numerator.Expand().InnerSimplified, denominator.Expand().InnerSimplified, out var cancelled) && cancelled is not null + ? Functions.PartialFractions.Bare(cancelled) : combined; + } + } + /// /// A binomial differential x^m (a + b x^n)^p with a symbol in its exponents, in /// Chebyshev's third case at a whole (m + 1)/n + p + 1 = -k of at most 0: Rubi's @@ -20387,7 +20850,24 @@ private static bool AreProportionalAtSampledPoints(Entity left, Entity right, En static bool IsPlain(Entity tower, Entity.Variable x) => tower is Powf(_, var g) && TreeAnalyzer.TryGetPolyLinear(g, x, out _, out _) || tower is Logf(_, var argument) && argument == x; - if (towers.All(tower => IsPlain(tower, x))) + // Unless a sum below the bar holds x beside the logarithm, to two or more powers: + // `(1 - ln(x))/(x (x + ln(x)))` is `(ln(x + ln(x)) - ln(x))'`, a logarithm of that sum, which + // the tower's ansatz, a polynomial in `ln(x)`, does not write. Where one power of x is a + // factor of every term, `x + 4 x ln(x)^2`, the sum is that power times one in the logarithm, + // which the rules for a function of `ln(x)` answer as written. + static bool MixesTheVariableWithItsLogarithm(Entity expr, Entity.Variable x) + { + var logarithm = MathS.Ln(x); + var placeholder = Variable.CreateUnique(expr, "t_mixed"); + foreach (var term in Sumf.LinearChildren(expr)) + foreach (var factor in Mulf.LinearChildren(Functions.SingleQuotient.Of(term).Denominator)) + if ((factor is Powf(var b, Number.Integer) ? b : factor) is var @base && @base is Sumf or Minusf && @base.ContainsNode(logarithm) + && TreeAnalyzer.TryGetPolynomial(@base.Replace(node => node == logarithm ? placeholder : node), x, out var inX) + && inX.Count > 1) + return true; + return false; + } + if (towers.All(tower => IsPlain(tower, x)) && !MixesTheVariableWithItsLogarithm(expr, x)) return null; var monomials = new List<(Entity.Variable Indeterminate, Entity Written, Entity Derivative)>(); var map = new Dictionary(); @@ -20475,6 +20955,40 @@ Entity Differentiated(Entity generator) return result; } var exponentialIndeterminates = monomials.Where(m => m.Written is Powf).Select(m => m.Indeterminate).ToList(); + Entity.Variable? logarithmOfX = monomials.Where(m => m.Written == MathS.Ln(x)).Select(m => (Entity.Variable?)m.Indeterminate).FirstOrDefault(); + // Each exponential of a quadratic in x alone, with the error function its integral is. + var gaussians = new List<(Entity Generator, Entity Kappa, Entity.Variable Indeterminate)>(); + foreach (var (indeterminate, written, _) in monomials) + if (written is Powf(_, var exponent) && !monomials.Any(m => exponent.ContainsNode(m.Indeterminate)) + && TreeAnalyzer.TryGetPolyQuadratic(exponent, x, out var leading, out var middle, out var last) + && !leading.ContainsNode(x) && !middle.ContainsNode(x) && !last.ContainsNode(x) && !VanishesIdentically(leading)) + { + var rootOfTheLeading = MathS.Sqrt(leading); + var shift = Functions.PartialFractions.Bare((middle / (2 * leading)).Simplify()); + var offset = Functions.PartialFractions.Bare((MathS.Sqr(middle) / (4 * leading) - last).Simplify()); + var kappa = 2 * rootOfTheLeading / MathS.Sqrt(MathS.pi) * (VanishesIdentically(offset) ? Number.Integer.One : MathS.Pow(MathS.e, offset)); + gaussians.Add((MathS.Erfi(rootOfTheLeading * (x + shift)), kappa, indeterminate)); + } + // The denominator with one factor `t` taken out, where it has one as a factor. + static Entity? TheDenominatorWithoutOne(Entity below, Entity.Variable t) + { + var factors = Mulf.LinearChildren(below).ToList(); + for (var i = 0; i < factors.Count; i++) + { + Entity? rest = factors[i] == t ? Number.Integer.One + : factors[i] is Powf(var b, Number.Integer n) && b == t && n.EInteger.CompareTo(EInteger.One) > 0 + ? (n.EInteger.Equals(EInteger.FromInt32(2)) ? t : MathS.Pow(t, Number.Integer.Create(n.EInteger.Subtract(EInteger.One)))) + : null; + if (rest is null) + continue; + Entity product = rest; + for (var j = 0; j < factors.Count; j++) + if (j != i) + product = product * factors[j]; + return product; + } + return null; + } // Q: the factors with their powers lowered by one, and then as they are; each // exponential indeterminate to every power from zero up to the numerator's degree // in it plus one, since a power of an exponential above the bar is one below in @@ -20499,6 +21013,18 @@ Entity Differentiated(Entity generator) return System.Math.Min(MaximumTowerDegree, aboveDegree + 2); }).ToList(); var logarithmGenerators = factors.Select(f => f.Factor).Distinct().ToList(); + // The special functions in a second pass only, where the first found nothing: beside them + // every system has a column more, and Hearn's `e^(x^2)` integrands, answered without one, + // took five to ten times as long. The Gaussian's only beside a factor of x alone below the + // bar, as in `e^((a + x)^2)/x^2`: a polynomial times exponentials of a quadratic is the + // rules' for those, which answer it after this, and the pass declined `x^2 sinh(a + b x + c x^2)` + // and `x/(a + b csch(c + d x^2))` first at three times their allocation. And one Gaussian: a + // hyperbolic function of a quadratic is two, `sinh(a + b x^2)^2/x^2` among them, which the pass + // declined at three times the allocation of the rules after it. + var aGaussianBesideAPowerOfX = gaussians.Count == 1 + && Mulf.LinearChildren(below).Any(factor => factor.ContainsNode(x) && !monomials.Any(m => factor.ContainsNode(m.Indeterminate))); + var specialOnes = aGaussianBesideAPowerOfX || logarithmOfX is { } lnX && TheDenominatorWithoutOne(below, lnX) is not null; + foreach (var withTheSpecialFunctions in specialOnes ? new[] { false, true } : new[] { false }) foreach (var q in qCandidates) foreach (var exponentPowers in Combinations(exponentialBounds)) { @@ -20568,6 +21094,33 @@ Entity Differentiated(Entity generator) } columns.Add(column); } + // And the logarithmic integral, `li(x)' = 1/ln(x)`, where `ln(x)` is a factor of the + // integrand's denominator: `(1 + x)/(ln(x) (x + ln(x)))` is + // `(li(x) + ln(ln(x)) - ln(x + ln(x)))'`, and without it no rational part and no + // logarithms answer it. + if (withTheSpecialFunctions && !failed && logarithmOfX is { } t && TheDenominatorWithoutOne(below, t) is { } belowWithoutIt) + { + var numerator = MathS.Sqr(denominator) * belowWithoutIt * logArguments * towersBelow; + if (TryReadInAll(Functions.PartialFractions.Bare(numerator.Expand().InnerSimplified), variables, out var column)) + { + generators.Add(MathS.Li(x)); + columns.Add(column); + } + } + // And the Gaussian's, for an exponential `t = e^Q` of a quadratic: `erfi(sqrt(A) (x + B/(2 A)))'` + // is `kappa t`. `e^((a + x)^2)/x^2 - 2 a e^((a + x)^2)/x` is + // `(-e^((a + x)^2)/x + sqrt(pi) erfi(a + x))'`: by parts the first term leaves the second's + // negative and `2 e^((a + x)^2)`, and neither term is elementary on its own. + if (withTheSpecialFunctions && !failed && aGaussianBesideAPowerOfX) + foreach (var (generator, kappa, indeterminate) in gaussians) + { + var numerator = kappa * indeterminate * MathS.Sqr(denominator) * below * logArguments * towersBelow; + if (TryReadInAll(Functions.PartialFractions.Bare(numerator.Expand().InnerSimplified), variables, out var column)) + { + generators.Add(generator); + columns.Add(column); + } + } if (failed) continue; var targetNumerator = above * MathS.Sqr(denominator) * logArguments * towersBelow; @@ -20583,7 +21136,12 @@ Entity Differentiated(Entity generator) matrix[row][c] = columns[c].TryGetValue(keys[row], out var entry) ? entry : Number.Integer.Zero; rhs[row] = target.TryGetValue(keys[row], out var wanted) ? wanted : Number.Integer.Zero; } - if (!Functions.PartialFractions.TrySolveLinear(matrix, rhs, out var values) || values is null) + // With a special function among the generators, whose constant carries the symbols + // into the columns, an elimination that leaves `2 a - 2 a` as written declines; then + // with the symbols pinned first, and the columns the pinned solutions need. + var special = generators.Count > monomialGenerators.Count + logarithmGenerators.Count; + if ((!Functions.PartialFractions.TrySolveLinear(matrix, rhs, out var values) || values is null) + && (!special || !Functions.PartialFractions.TrySolveLinearWithSymbols(matrix, rhs, out values) || values is null)) continue; Entity answer = Number.Integer.Zero; for (var c = 0; c < values.Length; c++) @@ -20591,6 +21149,10 @@ Entity Differentiated(Entity generator) var value = values[c].InnerSimplified; if (value.Evaled is Number.Complex { IsZero: true }) continue; + // A special function's coefficient is a quotient of its constant, `1/(2/sqrt(pi))` + // as the elimination leaves it: small, and simplified. + if (c >= monomialGenerators.Count + logarithmGenerators.Count && value.Complexity < 80) + value = Functions.PartialFractions.Bare(value.Simplify()); answer += value * generators[c]; } if (answer == Number.Integer.Zero) @@ -31693,6 +32255,38 @@ private static Entity WithLogarithmArgumentsOverOneBar(Entity expr, Entity.Varia /// against the integrand at sampled points before it is given. /// internal static Entity? SolveByTheQuotientOfTheLogarithmsLinears(Entity expr, Entity.Variable x, bool integrateByParts) + => SolveByTheQuotientOfTheLogarithmsLinears(expr, x, integrateByParts, ontoTheExponentialIntegral: false); + + /// + /// where + /// the logarithm's factor is a negative whole power and the integrand in t is a power + /// of t and a polynomial beside it, t^q (b - d t)^j (A + B ln(K t^n))^(-k): the + /// powers of the linears sum to -2 - j, and + /// answers it onto the exponential integral, a power of t at a time. Rubi's + /// 1/((a g + b g x)^2 (A + B ln(e (a + b x)/(c + d x)))^2) is + /// t^(-2) (A + B ln(e t))^(-2)/(g^2 D). + /// + /// + /// Before the polynomial term, which takes b^2 g^2 out of (a g + b g x)^2 and asks + /// 1/((a/b + x)^2 (A + B ln(...))^2) again, and a level below that the substitution + /// search ran past the budget: the rule above, far down the chain, was not reached. Of Rubi's + /// 3.2.1 those rows were from 1.3 s to 23 s, and four past the budget. + /// https://github.com/asc-community/AngouriMath/issues/233 + /// + internal static Entity? SolveByTheQuotientOfTheLogarithmsLinearsOntoTheExponentialIntegral(Entity expr, Entity.Variable x) + { + // Read first, as cheaply as it can be: a factor with a logarithm in it below the bar, the + // logarithm of a quotient. Early in the chain, the rule's reading of the linears simplified + // a derivative for every logarithm it was shown. + static bool OfAQuotient(Entity node) => node is Logf(_, var argument) + && argument.Nodes.Any(inner => inner is Divf || inner is Powf(_, Number.Integer { IsNegative: true })); + if (!FactorsOfTheIntegrand(expr).Any(pair => pair.Factor.Nodes.Any(OfAQuotient) + && (pair.Underneath ? pair.Factor is not Powf(_, Number.Integer { IsNegative: true }) : pair.Factor is Powf(_, Number.Integer { IsNegative: true })))) + return null; + return SolveByTheQuotientOfTheLogarithmsLinears(expr, x, integrateByParts: false, ontoTheExponentialIntegral: true); + } + + private static Entity? SolveByTheQuotientOfTheLogarithmsLinears(Entity expr, Entity.Variable x, bool integrateByParts, bool ontoTheExponentialIntegral) { if (!Integration.AnsweringTheQuestionAskedOrOneBelow) return null; @@ -31785,6 +32379,7 @@ private static Entity WithLogarithmArgumentsOverOneBar(Entity expr, Entity.Varia Entity powerOfT = Number.Integer.Zero; Entity powerOfDenominator = Number.Integer.Create(-2); // dx = D dt/(b - d t)^2 Entity? logarithmsFactor = null; + Entity logarithmsPower = Number.Integer.One; // A polynomial beside the logarithm is by parts in x as it stands, one closed // step; in t it is a power of b - d t below the bar, which is the slower way. var aPolynomialBesideTheLogarithm = true; @@ -31808,6 +32403,7 @@ private static Entity WithLogarithmArgumentsOverOneBar(Entity expr, Entity.Varia if (inT.ContainsNode(x)) return null; logarithmsFactor = MathS.Pow(inT, power); + logarithmsPower = power; continue; } if (!TryReadLinear(@base, out var slope, out var offset)) @@ -31836,6 +32432,12 @@ private static Entity WithLogarithmArgumentsOverOneBar(Entity expr, Entity.Varia if (readFromTheReciprocal && powerOfDenominator.Evaled is Number.Complex { IsZero: true } && powerOfT.Evaled is not Number.Complex { IsZero: true }) return null; + // In t a power of t beside a polynomial and the logarithm's power: the linears' powers sum + // to -2 or less, and b - d t stands above the bar to a whole power. + if (ontoTheExponentialIntegral && (logarithmsPower.Evaled is not Number.Integer { IsNegative: true } + || Functions.PartialFractions.Bare(powerOfDenominator.Simplify()).Evaled is not Number.Integer { IsNegative: false } whole + || whole.EInteger.CompareTo(EInteger.FromInt32(12)) > 0)) + return null; Entity above = constantInFront * logarithmsFactor, below = Number.Integer.One; foreach (var (@base, power) in new[] { (t, Functions.PartialFractions.Bare(powerOfT.Simplify())), (denominatorInT, Functions.PartialFractions.Bare(powerOfDenominator.Simplify())) }) { @@ -31849,7 +32451,8 @@ private static Entity WithLogarithmArgumentsOverOneBar(Entity expr, Entity.Varia var rewritten = (below == Number.Integer.One ? above : above / below).InnerSimplified; if (rewritten.ContainsNode(x)) return null; - if (Integration.ComputeAsAQuestionOfItsOwn(rewritten, t, integrateByParts) is not { } integrated) + if ((ontoTheExponentialIntegral ? SolveAPowerOverAPowerOfALogarithm(rewritten, t) + : Integration.ComputeAsAQuestionOfItsOwn(rewritten, t, integrateByParts)) is not { } integrated) return null; var answer = integrated.Substitute(t, l1.Linear / l2.Linear); if (answer.Nodes.Any(node => node == MathS.NaN)) @@ -31951,7 +32554,7 @@ static bool AreEqualAsPolynomials(Entity left, Entity right) /// with m <= k - 2: P written in powers of L at its root, /// P = p_0 + p_1 L + ... + p_m L^m, makes the integrand sum p_j G/L^(k - j) /// exactly, and each term is a power of one of the logarithm's linears beside a function of - /// it, which answers under the + /// it, which answers under the /// quotient. Rubi's 7.3.4 x (a + b atanh(c x))^2/(d + c d x)^3 is /// ((a + b atanh(c x))^2/(d + c d x)^2 - d (a + b atanh(c x))^2/(d + c d x)^3)/(c d). /// https://github.com/asc-community/AngouriMath/issues/233 diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index 2e33e1ec9..6f03cb88e 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -766,6 +766,11 @@ internal static bool ASumOverRootsWouldAnswer // since it fires only where that factor holds x, and the rules before it spent 28 seconds on // `atan(e x/sqrt(-a e^2/b - e^2 x^2))/sqrt(a + b x^2)` declining. if ((answer = IndefiniteIntegralSolver.SolveAFunctionOfAnInverseTimesItsDerivativeUpToAConstant(expr, x, integrateByParts)) is { }) return answer; + // The derivative of a function over a power of a logarithm, an inverse trigonometric + // function, or a sum with a logarithm in it, whose terms are not elementary apart: early, + // since it reads only a power of one below the bar, and by parts on either term asked the + // other's question again. + if ((answer = IndefiniteIntegralSolver.SolveAsTheDerivativeOfAFunctionOverAPowerOfALogarithmOrAnInverse(expr, x)) is { }) return answer; // A whole power of a linear in the inverse sine, cosine or hyperbolic cosine of // `c + d x^2` with `c^2 = 1`, by its recurrence, which closes: by parts it went round. if ((answer = IndefiniteIntegralSolver.SolveAPowerOfAnInverseFunctionOfAQuadraticByTheRecurrence(expr, x)) is { }) return answer; @@ -849,6 +854,13 @@ internal static bool ASumOverRootsWouldAnswer // the sum and writes what is left out in powers of x, where the linear is gone; and // before the substitution search, which reads one function at a time. if ((answer = IndefiniteIntegralSolver.SolveByWritingAFunctionOfOneShiftedLinear(expr, x, integrateByParts)) is { }) return answer; + // And a transcendental function of one shifted linear alone, for the same reason: the + // polynomial term writes `(c/d + x)^19` beside `F^(a + b/(c + d x)^3)`, two spellings. + if ((answer = IndefiniteIntegralSolver.SolveByWritingATranscendentalFunctionOfOneShiftedLinear(expr, x, integrateByParts)) is { }) return answer; + // Powers of two linears beside a negative power of a logarithm of their quotient, which in + // the quotient are a power of it beside the logarithm alone, onto the exponential integral: + // the polynomial term takes the slope out of the linear and asks it again, past the budget. + if ((answer = IndefiniteIntegralSolver.SolveByTheQuotientOfTheLogarithmsLinearsOntoTheExponentialIntegral(expr, x)) is { }) return answer; // A polynomial beside error functions of one shifted linear, written in it: by parts in // x leaves a Gaussian of the shifted linear beside a polynomial, which nothing reads. if ((answer = IndefiniteIntegralSolver.SolveAPolynomialBesideErrorFunctionsOfAShiftedLinear(expr, x, integrateByParts)) is { }) return answer; diff --git a/Sources/Tests/UnitTests/Calculus/AFunctionOverAPowerOfALogarithmOrAnInverseIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/AFunctionOverAPowerOfALogarithmOrAnInverseIntegralTest.cs new file mode 100644 index 000000000..46cb5ab02 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/AFunctionOverAPowerOfALogarithmOrAnInverseIntegralTest.cs @@ -0,0 +1,67 @@ +// +// 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 derivative of U/V^k for V a logarithm or an inverse trigonometric function, whose + /// terms are not elementary apart: U = -B/(k V') for the terms B over V^(k + 1); and + /// for V a sum with a logarithm in it, U a monomial of the numerator, read monomial by + /// monomial. Rubi's 3.2.3, 3.5 and 5.3.4 rows ran past the budget. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class AFunctionOverAPowerOfALogarithmOrAnInverseIntegralTest + { + private static readonly (string Symbol, double Value)[] Pins = + { + ("a", 1.3), ("b", 0.7), ("c", 0.4), ("d", 0.3), ("f", 0.9), ("g", 1.1), ("k", 0.6), ("h", 0.5), + ("m", 1.7), ("n", 2.3), ("p", 1.6), ("q", 2.2), ("A", 0.8), ("B", 1.2), ("F", 2.4), + }; + + [Theory] + [InlineData("1/((c + d*x)*(-a + c + (-b + d)*x)*ln((a + b*x)/(c + d*x))) + ln(1 + (-a - b*x)/(c + d*x))/((a + b*x)*(c + d*x)*ln((a + b*x)/(c + d*x))^2)")] + [InlineData("(-1)/((a + b*x)*(a - c + (b - d)*x)*ln((a + b*x)/(c + d*x))) + ln(1 + (-c - d*x)/(a + b*x))/((a + b*x)*(c + d*x)*ln((a + b*x)/(c + d*x))^2)")] + [InlineData("x^3/((1 + a^2*x^2)*atan(a*x)^3) - 3/2*x^2/(a*atan(a*x)^2)")] + [InlineData("(a*d*n*x^m - a*d*m*x^m*ln(c*x^n) - b*d*n*(q - 1)*ln(c*x^n)^q)/(x*(a*x^m + b*ln(c*x^n)^q)^2)")] + public void AsTheDerivativeOfAQuotient(string integrand) => AssertAnswered(integrand); + + private static void AssertAnswered(string integrand, string variable = "x") + { + var integral = integrand.ToEntity().Integrate(variable); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 5000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) + { + foreach (var (symbol, value) in Pins) + if (symbol != variable) + e = e.Substitute(symbol, value); + return e; + } + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate(variable); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 1.9, 2.9, 3.7 }) + { + var want = original.Substitute(variable, at).EvalNumerical(); + var got = derivative.Substitute(variable, 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/d{variable} of the antiderivative of {integrand} is {got} at {variable} = {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/ARischNormanAnsatzBesideTheLogarithmOfXIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ARischNormanAnsatzBesideTheLogarithmOfXIntegralTest.cs new file mode 100644 index 000000000..c9d76d22d --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ARischNormanAnsatzBesideTheLogarithmOfXIntegralTest.cs @@ -0,0 +1,68 @@ +// +// 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 Risch-Norman ansatz on ln(x) alone where a sum below the bar holds x beside it to two + /// or more powers, with the logarithmic integral among its generators: (1 - ln(x))/(x (x + ln(x))) + /// is (ln(x + ln(x)) - ln(x))', and (1 + x)/(ln(x) (x + ln(x))) is + /// (li(x) + ln(ln(x)) - ln(x + ln(x)))'. Rubi's 3.5. And the error function among them for an + /// exponential of a quadratic: e^((a + x)^2)/x^2 - 2 a e^((a + x)^2)/x is + /// (-e^((a + x)^2)/x + sqrt(pi) erfi(a + x))', Rubi's 2.3, whose terms are not elementary apart. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ARischNormanAnsatzBesideTheLogarithmOfXIntegralTest + { + private static readonly (string Symbol, double Value)[] Pins = + { + ("a", 1.3), ("b", 0.7), ("c", 0.4), ("d", 0.3), ("f", 0.9), ("g", 1.1), ("k", 0.6), ("h", 0.5), + ("m", 1.7), ("n", 2.3), ("p", 1.6), ("q", 2.2), ("A", 0.8), ("B", 1.2), ("F", 2.4), + }; + + [Theory] + [InlineData("(1 - ln(x))/(x*(x + ln(x)))")] + [InlineData("(1 + x)/(ln(x)*(x + ln(x)))")] + [InlineData("e^((a + x)^2)/x^2 - 2*e^((a + x)^2)*a/x")] + public void BesideTheLogarithmOfX(string integrand) => AssertAnswered(integrand); + + private static void AssertAnswered(string integrand, string variable = "x") + { + var integral = integrand.ToEntity().Integrate(variable); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 5000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) + { + foreach (var (symbol, value) in Pins) + if (symbol != variable) + e = e.Substitute(symbol, value); + return e; + } + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate(variable); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 1.9, 2.9, 3.7 }) + { + var want = original.Substitute(variable, at).EvalNumerical(); + var got = derivative.Substitute(variable, 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/d{variable} of the antiderivative of {integrand} is {got} at {variable} = {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/ATranscendentalFunctionOfOneShiftedLinearIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/ATranscendentalFunctionOfOneShiftedLinearIntegralTest.cs new file mode 100644 index 000000000..e507d9e32 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/ATranscendentalFunctionOfOneShiftedLinearIntegralTest.cs @@ -0,0 +1,67 @@ +// +// 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 transcendental function of one shifted linear alone, written in u = c + d x, where a + /// power of the linear stands below the bar inside a sum or the linear is written two ways. Rubi's + /// F^(a + b/(c + d x)^3)/(c + d x)^19 ran past the budget in x, and is 20 ms in u. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class ATranscendentalFunctionOfOneShiftedLinearIntegralTest + { + private static readonly (string Symbol, double Value)[] Pins = + { + ("a", 1.3), ("b", 0.7), ("c", 0.4), ("d", 0.3), ("f", 0.9), ("g", 1.1), ("k", 0.6), ("h", 0.5), + ("m", 1.7), ("n", 2.3), ("p", 1.6), ("q", 2.2), ("A", 0.8), ("B", 1.2), ("F", 2.4), + }; + + [Theory] + [InlineData("F^(a + b/(c + d*x)^3)/(c + d*x)^19")] + [InlineData("F^(a + b/(c + d*x)^3)/(c + d*x)^16")] + [InlineData("F^(a + b/(c + d*x)^2)/(c + d*x)^10")] + [InlineData("ln(c*(d + k/(f + g*x)^3)^q)")] + [InlineData("(k*g + f*g*x)^(p - 1)/ln(d*(k + f*x)^p)")] + public void WrittenInTheShiftedLinear(string integrand) => AssertAnswered(integrand); + + private static void AssertAnswered(string integrand, string variable = "x") + { + var integral = integrand.ToEntity().Integrate(variable); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 5000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) + { + foreach (var (symbol, value) in Pins) + if (symbol != variable) + e = e.Substitute(symbol, value); + return e; + } + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate(variable); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 1.9, 2.9, 3.7 }) + { + var want = original.Substitute(variable, at).EvalNumerical(); + var got = derivative.Substitute(variable, 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/d{variable} of the antiderivative of {integrand} is {got} at {variable} = {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/TheExponentialsAloneAsTheirOwnDerivativeIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/TheExponentialsAloneAsTheirOwnDerivativeIntegralTest.cs new file mode 100644 index 000000000..60b31e5e3 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/TheExponentialsAloneAsTheirOwnDerivativeIntegralTest.cs @@ -0,0 +1,68 @@ +// +// 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 product whose logarithmic derivative is a constant, or whose exponentials' logarithmic derivative + /// is the rest of it up to a constant -- with f^g = e^(g ln f) for x in the base and the exponent: + /// (a F^(c + d x))^m (b F^(k + f x))^n, x^(-2 - 1/x) (1 - ln(x)), cosh(x)^x (ln(cosh(x)) + x tanh(x)). + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class TheExponentialsAloneAsTheirOwnDerivativeIntegralTest + { + private static readonly (string Symbol, double Value)[] Pins = + { + ("a", 1.3), ("b", 0.7), ("c", 0.4), ("d", 0.3), ("f", 0.9), ("g", 1.1), ("k", 0.6), ("h", 0.5), + ("m", 1.7), ("n", 2.3), ("p", 1.6), ("q", 2.2), ("A", 0.8), ("B", 1.2), ("F", 2.4), + }; + + [Theory] + [InlineData("(a*F^(c + d*x))^m*(b*F^(k + f*x))^n")] + [InlineData("F^(a + b*x)*(cosh(c + d*x) + sinh(c + d*x))^n")] + [InlineData("F^(a + b*x)*(cosh(c + d*x) - sinh(c + d*x))^n")] + [InlineData("x^(-2 - 1/x)*(1 - ln(x))")] + [InlineData("cosh(x)^x*(ln(cosh(x)) + x*tanh(x))")] + [InlineData("x^x*(1 + ln(x))")] + public void TheirOwnDerivative(string integrand) => AssertAnswered(integrand); + + private static void AssertAnswered(string integrand, string variable = "x") + { + var integral = integrand.ToEntity().Integrate(variable); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 5000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) + { + foreach (var (symbol, value) in Pins) + if (symbol != variable) + e = e.Substitute(symbol, value); + return e; + } + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate(variable); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 1.9, 2.9, 3.7 }) + { + var want = original.Substitute(variable, at).EvalNumerical(); + var got = derivative.Substitute(variable, 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/d{variable} of the antiderivative of {integrand} is {got} at {variable} = {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/TheQuotientOfTheLogarithmsLinearsOntoTheExponentialIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/TheQuotientOfTheLogarithmsLinearsOntoTheExponentialIntegralTest.cs new file mode 100644 index 000000000..279d2feb7 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/TheQuotientOfTheLogarithmsLinearsOntoTheExponentialIntegralTest.cs @@ -0,0 +1,75 @@ +// +// 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 +{ + /// + /// Powers of two linears beside a negative power of a logarithm of their quotient, which under + /// t = (a + b x)/(c + d x) are a power of t and a polynomial beside the logarithm alone, + /// onto the exponential integral: Rubi's 3.2.1 rows with the powers of the linears summing to + /// -2 or -3, which ran past the budget. And the square below the bar the substitution + /// u = ln(k t) writes, whose discriminant is zero: its two roots were one, and the answer had + /// no value at any point. + /// #233 + /// + [Trait("Area", "Calculus")] + public sealed class TheQuotientOfTheLogarithmsLinearsOntoTheExponentialIntegralTest + { + private static readonly (string Symbol, double Value)[] Pins = + { + ("a", 1.3), ("b", 0.7), ("c", 0.4), ("d", 0.3), ("f", 0.9), ("g", 1.1), ("k", 0.6), ("h", 0.5), + ("m", 1.7), ("n", 2.3), ("p", 1.6), ("q", 2.2), ("A", 0.8), ("B", 1.2), ("F", 2.4), + }; + + [Theory] + [InlineData("1/((a*g + b*g*x)^2*(A + B*ln(k*((a + b*x)/(c + d*x))^n))^2)")] + [InlineData("1/((a*g + b*g*x)^2*(A + B*ln(k*(a + b*x)/(c + d*x)))^2)")] + [InlineData("1/((a*g + b*g*x)^2*(A + B*ln(k*(a + b*x)^2/(c + d*x)^2))^2)")] + [InlineData("1/((a*g + b*g*x)^2*(A + B*ln(k*(c + d*x)^2/(a + b*x)^2)))")] + [InlineData("1/((a*g + b*g*x)^3*(A + B*ln(k*(a + b*x)/(c + d*x)))^2)")] + public void OntoTheExponentialIntegral(string integrand) => AssertAnswered(integrand); + + [Theory] + [InlineData("1/(x^2*(A + B*ln(k*x))^2)")] + [InlineData("1/(x^2*(A + B*ln(k*x^2))^2)")] + public void OverASquareWithAZeroDiscriminant(string integrand) => AssertAnswered(integrand); + + private static void AssertAnswered(string integrand, string variable = "x") + { + var integral = integrand.ToEntity().Integrate(variable); + var text = integral.Stringize(); + Assert.DoesNotContain("integral(", text); + Assert.True(text.Length < 5000, $"{text.Length} characters of answer for {integrand}"); + Entity Pinned(Entity e) + { + foreach (var (symbol, value) in Pins) + if (symbol != variable) + e = e.Substitute(symbol, value); + return e; + } + var derivative = Pinned(integral.Substitute("C", 0)).Differentiate(variable); + var original = Pinned(integrand.ToEntity()); + var compared = 0; + foreach (var at in new[] { 0.3, 0.8, 1.3, 1.9, 2.9, 3.7 }) + { + var want = original.Substitute(variable, at).EvalNumerical(); + var got = derivative.Substitute(variable, 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/d{variable} of the antiderivative of {integrand} is {got} at {variable} = {at}, where the integrand is {want}"); + } + Assert.True(compared >= 5, $"only {compared} points could be compared for {integrand}"); + } + } +}