diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index f13d7d318..3f68253aa 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -557,6 +557,24 @@ Every spelling but `lcm` was a parse error in 2.5.0, the congruence node being n | `forall n in ZZ : 6 divides n^3 + 5 n` | `UnhandledParseException` | `True` | | `exists x, y in ZZ : 3 x^2 - 5 y^2 = 1` | `UnhandledParseException` | `False` | +### An antiderivative built through an even root says where the root is real + +The linear-radical substitution, `u = (a x + b)^(1/q)`, and the general substitution for an even +root both take a power of `u` out of a radical on the strength of `u` being a non-negative real — +`sqrt(u^2 (1 + u))` as `u sqrt(1 + u)` — which it is only where the radical is. Beyond the +radicand's zero `u` is imaginary, and the integrand can still be real there: Rubi's +`x sqrt(c - a c x) / e^(3 atanh(a x))` above `a x = 1` is the product of two imaginary factors, +and the answer built for a real `u` is not its antiderivative — 14% off at `x = 0.59` with +`a = 3.1`, `c = 0.7`, and it shipped as an answer. An answer that used the step now carries the +condition, `provided a x + b >= 0`: `NaN` beyond it, which is no claim, and unchanged where it +holds. An answer the step did not touch — `x/sqrt(1 + x)`, `1/(sqrt(x) - x^(-1/3))` — is as it +was ([#718](https://github.com/asc-community/AngouriMath/issues/718)). + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"x * sqrt(c - a*c*x) / e^(3 * atanh(a*x))".Integrate("x")` | `integral(…)` — left unevaluated, with `e` evaluated to a hundred digits inside it | an antiderivative in `sqrt(c - a c x)`, `provided c - a * c * x >= 0` | +| `"1/sqrt(x + x^(3/2))".Integrate("x")` | `integral(1 / sqrt(x + x ^ (3/2)), x)` — left unevaluated | `2 * sqrt(1 + sqrt(x)) / (1/2) + C provided x >= 0` | + ### `binomial(n, k)` is a function **Addition, not silent.** The binomial coefficient is a node, `Entity.Binomialf`, spelled diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index 780c59f02..1b61a3c75 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -5169,10 +5169,16 @@ static bool IsTwice(Entity argument, Entity.Variable x) /// two substitutions at once and is declined rather than half-rewritten. /// /// - /// No condition is owed by the rewrite. a is non-zero, since a zero + /// No condition is owed by the rewrite itself. a is non-zero, since a zero /// coefficient is not a linear expression in x and the reader below rejects it, /// and u^q = a*x + b is invertible wherever the radical it came from is defined. - /// The answer inherits the radical's own domain and adds nothing. + /// The answer inherits the radical's own domain and adds nothing -- except where the + /// even-root step below has read u as a non-negative real, which it is only + /// where the radical is real. Beyond that, u is imaginary and the integrand can + /// still be real: x sqrt(c - a c x) / e^(3 atanh(a x)) above a x = 1 is the + /// product of two imaginary factors, and the answer built for a real u is not its + /// antiderivative there (Rubi 7.3.6). So the answer says what it used: it is given + /// provided a x + b >= 0 exactly when that step changed the integrand. /// /// https://github.com/asc-community/AngouriMath/issues/718 /// @@ -5304,13 +5310,24 @@ node is Powf(var radical, Number.Rational e) && e is not Number.Integer integrand = inner; // For an even q the principal root is not negative wherever it is real, and a root // holding a power of u gives that power up: `1/sqrt(t + t^(3/2))` under `u = sqrt(t)` - // is `2u/sqrt(u^2 + u^3)`, a root of a cubic, and is `2/sqrt(1 + u)`. + // is `2u/sqrt(u^2 + u^3)`, a root of a cubic, and is `2/sqrt(1 + u)`. Where the + // radical is not real, u is not a non-negative real and the identities do not hold, + // so an answer that used them is conditioned on the radicand -- unless the caller + // knows its variable to be non-negative, as an exponential is. + // https://github.com/asc-community/AngouriMath/issues/718 + Entity? condition = null; if (q % 2 == 0) - integrand = FactorANonnegativeVariableOutOfRadicals(integrand, u); + { + var factored = FactorANonnegativeVariableOutOfRadicals(integrand, u); + if (factored != integrand && !variableIsNonnegative) + condition = radicalBase >= Number.Integer.Zero; + integrand = factored; + } - return Integration.ComputeIndefiniteIntegral(integrand, u, integrateByParts) is { } result - ? result.Substitute(u, MathS.Pow(radicalBase, Number.Rational.Create(1, q))) - : null; + if (Integration.ComputeIndefiniteIntegral(integrand, u, integrateByParts) is not { } result) + return null; + var back = result.Substitute(u, MathS.Pow(radicalBase, Number.Rational.Create(1, q))); + return condition is null ? back : back.Provided(condition); } /// @@ -13787,11 +13804,16 @@ private static bool IsARationalFunctionOfExponentials(Entity expr, Entity.Variab && Integration.ComputeIndefiniteIntegral(factored, uSub, integrateByParts) is { } forPositiveU && !forPositiveU.Nodes.Any(node => node == MathS.NaN)) { - // An even root is not negative, and what holds for u > 0 is the answer as - // it stands; `sqrt(x^(1/3))` under `1/(sqrt(x) - x^(-1/3))` was extended by - // parity it did not need, to an answer with a sign function in it. - if (u is Powf(_, Number.Rational root) && root.ERational.Denominator.IsEven) - return Functions.PartialFractions.Bare(forPositiveU).Substitute(uSub, u); + // An even root is not negative where it is real, and what holds for u > 0 + // is the answer as it stands there; `sqrt(x^(1/3))` under + // `1/(sqrt(x) - x^(-1/3))` was extended by parity it did not need, to an + // answer with a sign function in it. Where the radicand is negative the root + // is imaginary and the identities the factoring used do not hold, so the + // answer is conditioned on the radicand, as the linear-radical substitution + // conditions its own. https://github.com/asc-community/AngouriMath/issues/718 + if (u is Powf(var radicand, Number.Rational root) && root.ERational.Denominator.IsEven) + return Functions.PartialFractions.Bare(forPositiveU).Substitute(uSub, u) + .Provided(radicand >= Number.Integer.Zero); if (ExtendedByParity(integrandInU, forPositiveU, uSub) is { } onBothSides) return onBothSides.Substitute(uSub, u); } diff --git a/Sources/Tests/UnitTests/Calculus/LinearRadicalSubstitutionIntegralTest.cs b/Sources/Tests/UnitTests/Calculus/LinearRadicalSubstitutionIntegralTest.cs index bbb8e8eb9..065bf99b3 100644 --- a/Sources/Tests/UnitTests/Calculus/LinearRadicalSubstitutionIntegralTest.cs +++ b/Sources/Tests/UnitTests/Calculus/LinearRadicalSubstitutionIntegralTest.cs @@ -77,6 +77,46 @@ private static void DifferentiatesBack(string integrand) [InlineData("x/sqrt(1 + x^2 + (1 + x^2)^(3/2))")] public void APowerOfTheRootLeavesTheRadical(string integrand) => DifferentiatesBack(integrand); + /// + /// That step reads u as a non-negative real, which it is only where the radical + /// is real, and the answer says so: 1/sqrt(x + x^(3/2)) comes back + /// provided x >= 0. Beyond the radicand's zero u is imaginary, and the + /// integrand can still be real there -- Rubi's x sqrt(c - a c x) / e^(3 atanh(a x)) + /// above a x = 1 is the product of two imaginary factors -- while the answer built + /// for a real u is not its antiderivative: with a = 3.1, c = 0.7 at + /// x = 0.59 its derivative was 14% off, and it shipped as an answer. Now it is + /// NaN there, which is no claim, and right where the condition holds. + /// https://github.com/asc-community/AngouriMath/issues/718 + /// + [Fact] + public void TheEvenRootStepSaysWhereItHolds() + { + var conditioned = "1/sqrt(x + x^(3/2))".ToEntity().Integrate("x"); + Assert.Contains("provided x >= 0", conditioned.Stringize()); + + var integrand = "x * sqrt(c - a*c*x) / e^(3 * atanh(a*x))".ToEntity(); + var integral = integrand.Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Assert.Contains("provided c - a * c * x >= 0", integral.Stringize()); + var pinned = integral.Substitute("C", 0).Substitute("a", 3.1).Substitute("c", 0.7); + var derivative = pinned.Differentiate("x"); + var original = integrand.Substitute("a", 3.1).Substitute("c", 0.7); + // Inside the condition, a x < 1: the derivative is the integrand. + foreach (var at in new[] { 0.05, 0.17, 0.31 }) + { + var got = derivative.Substitute("x", at).EvalNumerical(); + var want = original.Substitute("x", at).EvalNumerical(); + Assert.False(got.IsNaN, $"no answer at x = {at}, inside the condition"); + Assert.True(Math.Abs((double)(got - want).RealPart) < 1e-9 * Math.Max(1, Math.Abs((double)want.RealPart)), + $"d/dx of the antiderivative is {got} at x = {at}, where the integrand is {want}"); + } + // Beyond it the integrand is real, 3.1 * 0.59 being above 1, and the answer claims + // nothing rather than something wrong. + var outside = original.Substitute("x", 0.59).EvalNumerical(); + Assert.True(Math.Abs((double)outside.ImaginaryPart) < 1e-12, $"the integrand at x = 0.59 is {outside}, not real"); + Assert.True(derivative.Substitute("x", 0.59).Evaled.IsNaN, "the answer claims a derivative outside its condition"); + } + /// /// A polynomial over a square root of something linear. None of these had an /// antiderivative, and each is a first-year exercise.