diff --git a/Sources/.editorconfig b/Sources/.editorconfig index 01632d728..7febbe20b 100644 --- a/Sources/.editorconfig +++ b/Sources/.editorconfig @@ -432,3 +432,6 @@ file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed [Tests/UnitTests/Calculus/LinearBesideTheRootOfAQuadraticTest.cs] file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed under MIT.\nDetails: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.\nWebsite: https://am.angouri.org.\n + +[Tests/UnitTests/Calculus/RationalFunctionBesideTheRootOfAQuadraticTest.cs] +file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed under MIT.\nDetails: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.\nWebsite: https://am.angouri.org.\n diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs index ec3f98581..4c9c75dd8 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs @@ -4541,6 +4541,151 @@ private static Entity SquarefreePartAsWritten(Entity denominator, Entity.Variabl return answer.Nodes.Any(node => node is Number.Complex { IsNaN: true }) ? null : answer; } + /// + /// A polynomial over a product of distinct linears, beside the square root of a + /// quadratic above or below the bar, taken apart as the rational function it is over + /// that root: N sqrt(Q)/D is N Q/(D sqrt(Q)), and P/D is a polynomial + /// plus a constant over each linear, so the integrand is a polynomial over the root + /// plus one K/((x - p) sqrt(Q)) per linear -- the rule before this one's shape, + /// each. + /// + /// + /// + /// Moses' sqrt(A^2 + B^2 (1 - y^2))/(1 - y^2): the root is of A^2 + B^2 - B^2 y^2, + /// so this is (A^2 + B^2 - B^2 y^2)/((1 - y^2) sqrt(Q)), whose rational part is + /// B^2 + (A^2)/(1 - y^2), and 1/(1 - y^2) is a half over 1 - y and a + /// half over 1 + y. Each piece is a line of the table, and the whole was declined: + /// the Euler substitution takes the generic first substitution for a symbolic leading + /// coefficient, with sqrt(-B^2) in it. + /// + /// + /// Only where there is something to take apart -- at least two pieces, each strictly + /// smaller: a polynomial over the root alone, or one linear with nothing divided out, + /// is the rule before this one's or the chain's already, and asking again would be + /// asking the same question. The polynomial part goes to the chain, for + /// x^n/sqrt(Q) is the trigonometric substitution's or Euler's by its shape; the + /// linears are closed. + /// + /// https://github.com/asc-community/AngouriMath/issues/718 + /// + internal static Entity? SolveARationalFunctionBesideTheRootOfAQuadratic(Entity expr, Entity.Variable x) + { + var (numerator, denominator) = Functions.SingleQuotient.Of(expr); + Entity? radicand = null; + var rootBelow = false; + var wholePower = 0; + Entity above = Number.Integer.One; + Entity below = Number.Integer.One; + foreach (var (side, isBelow) in new[] { (numerator, false), (denominator, true) }) + foreach (var factor in Mulf.LinearChildren(side)) + { + if (factor is Powf(var @base, Number.Rational half) && half is not Number.Integer && @base.ContainsNode(x)) + { + if (radicand is not null || half.ERational.Denominator.CompareTo(EInteger.FromInt32(2)) != 0 + || !half.ERational.Numerator.CanFitInInt32()) + return null; + radicand = @base; + var n = half.ERational.Numerator.ToInt32Unchecked(); + rootBelow = isBelow != (n < 0); + // Q^(n/2) is Q^((|n| - 1)/2) sqrt(Q) on the side it is on. + wholePower = (System.Math.Abs(n) - 1) / 2; + if (wholePower > 0) + { + if (rootBelow) below = below * MathS.Pow(@base, wholePower); + else above = above * MathS.Pow(@base, wholePower); + } + continue; + } + if (factor.ContainsNode(x) && !TreeAnalyzer.TryGetPolynomial(factor, x, out _)) + return null; + if (isBelow) below = below * factor; + else above = above * factor; + } + if (radicand is null || !TreeAnalyzer.TryGetPolyQuadratic(radicand, x, out var a, out _, out _) + || a.Evaled is Number.Complex { IsZero: true }) + return null; + // Written over the root: N sqrt(Q)/D is N Q/(D sqrt(Q)). + if (!rootBelow) + above = above * radicand; + if (!below.ContainsNode(x)) + return null; + + // The denominator as distinct linears: factored over the integers where it is + // numeric, and as written otherwise. + var factors = new List(); + Entity constant = Number.Integer.One; + var written = Functions.PolynomialFactoring.TryFactor(below, x, out var factored) && factored is not null ? factored : below; + foreach (var factor in Mulf.LinearChildren(written)) + { + if (!factor.ContainsNode(x)) + { + constant = constant * factor; + continue; + } + if (!TreeAnalyzer.TryGetPolyLinear(factor, x, out var slope, out _) || slope.Evaled is Number.Complex { IsZero: true }) + return null; + factors.Add(factor); + } + // Distinct: two factors with the same root are one repeated, which is not this. + for (var i = 0; i < factors.Count; i++) + for (var j = i + 1; j < factors.Count; j++) + if (!TryReadAsQuotient(factors[i] / factors[j], out var top, out var bottom) + || Functions.PolynomialGcd.TryCancel(top, bottom, out _)) + return null; + var linears = factors.Aggregate(Number.Integer.One as Entity, (product, factor) => product * factor); + + // P/D as a polynomial plus a proper part, the proper part over each linear. + Entity polynomialPart = Number.Integer.Zero; + Entity properNumerator = above; + if (TreeAnalyzer.PolynomialLongDivision(above, linears, genericCase: true, inTermsOf: x) is var (quotient, proper) + && quotient.Evaled is not Number.Complex { IsZero: true }) + { + polynomialPart = quotient; + var (properTop, properBottom) = Functions.SingleQuotient.Of(proper); + if (properBottom != linears && !TryReadAsQuotient(proper, out properTop, out properBottom)) + return null; + properNumerator = properTop; + if (properNumerator.ContainsNode(x) && !TreeAnalyzer.TryGetPolynomial(properNumerator, x, out _)) + return null; + } + var pieces = new List(); + if (polynomialPart.Evaled is not Number.Complex { IsZero: true }) + pieces.Add(polynomialPart); + if (properNumerator.Evaled is not Number.Complex { IsZero: true }) + { + if (factors.Count == 1) + pieces.Add(properNumerator / factors[0]); + else if (Functions.PartialFractions.TrySplitOverWrittenFactors(properNumerator, linears, x, out var decomposition) && decomposition is not null) + { + // The split comes back as a sum over a constant; the constant goes to each term. + var (terms, over) = decomposition is Divf(var splitTop, var splitBottom) && !splitBottom.ContainsNode(x) + ? (splitTop, splitBottom) : (decomposition, Number.Integer.One as Entity); + foreach (var term in Sumf.LinearChildren(terms)) + if (term.Evaled is not Number.Complex { IsZero: true }) + pieces.Add(term / over); + } + else + return null; + } + if (pieces.Count < 2) + return null; // nothing was taken apart + + var root = MathS.Pow(radicand, Number.Rational.Create(1, 2)); + Entity answer = Number.Integer.Zero; + foreach (var piece in pieces) + { + // Written as `coefficient / sqrt(Q)`, the shape the table reads for a constant. + var overTheRoot = (piece / constant).InnerSimplified / root; + var integrated = SolveALinearBesideTheRootOfAQuadratic(overTheRoot, x) + ?? (piece.ContainsNode(x) && !TreeAnalyzer.TryGetPolynomial(piece, x, out _) ? null : Integration.ComputeIndefiniteIntegral(overTheRoot, x, integrateByParts: false)); + if (integrated is null) + return null; + answer = answer + integrated; + } + answer = answer.InnerSimplified; + return answer.Nodes.Any(node => node is Number.Complex { IsNaN: true }) ? null : answer; + } + /// /// A rational function of x and one square root of a quadratic in x, /// rationalised by an Euler substitution and handed to the rational integrator. diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs index 669aa2ac7..34a4ea15f 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/Integration.Definition.cs @@ -557,6 +557,7 @@ private static Entity Normalized(Entity expr) => // level down -- the nested radicals under `u = sqrt(1 + x)`, the remainders of by // parts -- each of which is closed here. if ((answer = IndefiniteIntegralSolver.SolveALinearBesideTheRootOfAQuadratic(expr, x)) is { }) return answer; + if ((answer = IndefiniteIntegralSolver.SolveARationalFunctionBesideTheRootOfAQuadratic(expr, x)) is { }) return answer; if ((answer = IndefiniteIntegralSolver.SolveByEulerSubstitution(expr, x)) is { }) return answer; // And a root of a palindromic quartic, which is a root of a quadratic in x -+ 1/x: // Charlwood's `(1 + x^2)/((1 - x^2) sqrt(1 + x^4))` is `-du/(u sqrt(u^2 + 2))`. diff --git a/Sources/AngouriMath/Functions/TreeAnalyzer/RealValued.cs b/Sources/AngouriMath/Functions/TreeAnalyzer/RealValued.cs index 63abd01f9..943a6292d 100644 --- a/Sources/AngouriMath/Functions/TreeAnalyzer/RealValued.cs +++ b/Sources/AngouriMath/Functions/TreeAnalyzer/RealValued.cs @@ -47,6 +47,12 @@ internal static bool IsRealValued(Entity expr, Variable x) // Only an integer exponent keeps a real base real: x ^ (1/2) is not real below 0. case Powf(var @base, Integer): return IsRealValued(@base, x); + // A positive constant base is real to any real exponent: e^x, and so + // `e^x - 1/e^x - 2`, whose sign the antiderivative of a hyperbolic integrand + // carries under `u = e^x`. + case Powf(var constantBase, var exponent) when !constantBase.ContainsNode(x) + && constantBase.Evaled is Real { IsNaN: false, IsPositive: true }: + return IsRealValued(exponent, x); case Sinf or Cosf or Tanf or Cotanf or Secantf or Cosecantf or Arctanf or Arccotanf or Signumf: return IsRealValued(expr.DirectChildren[0], x); diff --git a/Sources/Tests/UnitTests/Calculus/ModulusPhaseTest.cs b/Sources/Tests/UnitTests/Calculus/ModulusPhaseTest.cs index 4f44e5bb2..e5472ef12 100644 --- a/Sources/Tests/UnitTests/Calculus/ModulusPhaseTest.cs +++ b/Sources/Tests/UnitTests/Calculus/ModulusPhaseTest.cs @@ -87,6 +87,11 @@ public void ADecidedNaNAtInfinityIsKept() [InlineData("abs(x / i)", "sgn(x) provided not x = 0")] [InlineData("abs(x + 2)", "sgn(x + 2) provided not x + 2 = 0")] [InlineData("sgn(x + 2)", "0 provided not x + 2 = 0")] + // A positive constant base is real to a real exponent: the sign of e^x - 1/e^x - 2 is + // what the antiderivative of a hyperbolic integrand carries under u = e^x, and its + // derivative was left unevaluated, which no numeric check of that antiderivative got past. + [InlineData("sgn(e^x - 2)", "0 provided not e^x - 2 = 0")] + [InlineData("sgn(2^x + 1/2^x - 2)", "0 provided not 2^x + 1/2^x - 2 = 0")] public void TheDerivativeOfARealValuedModulusIsTheSignRule(string expr, string expected) => Assert.Equal(expected.ToEntity(), expr.ToEntity().Differentiate(x)); @@ -107,6 +112,8 @@ public void TheDerivativeOfTheModulusOfAnImaginaryQuotientHasTheRightSign(string [InlineData("abs(a * x)")] [InlineData("sgn(2 + i * x)")] [InlineData("sgn(a * x)")] + [InlineData("sgn((-2)^x)")] + [InlineData("sgn(x^x)")] public void TheDerivativeOfAModulusNotShownRealIsLeftUnevaluated(string expr) => Assert.IsType(expr.ToEntity().Differentiate(x)); diff --git a/Sources/Tests/UnitTests/Calculus/RationalFunctionBesideTheRootOfAQuadraticTest.cs b/Sources/Tests/UnitTests/Calculus/RationalFunctionBesideTheRootOfAQuadraticTest.cs new file mode 100644 index 000000000..db73b64f6 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/RationalFunctionBesideTheRootOfAQuadraticTest.cs @@ -0,0 +1,93 @@ +// +// 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 over a product of distinct linears beside the square root of a quadratic, + /// above or below the bar, taken apart as the rational function it is over that root: + /// N sqrt(Q)/D is N Q/(D sqrt(Q)), P/D is a polynomial plus a constant + /// over each linear, and each piece is closed -- the polynomial over the root by the rules + /// for it, each K/((x - p) sqrt(Q)) by the reciprocal of its linear. + /// #718 + /// + /// + /// Moses' sqrt(A^2 + B^2 (1 - y^2))/(1 - y^2) was declined: the Euler substitution + /// takes its generic first substitution for a symbolic leading coefficient, with + /// sqrt(-B^2) in it. Taken apart it is B^2 over the root, which is the table's + /// arcsine, and A^2/2 over each of 1 - y and 1 + y beside it. + /// + [Trait("Area", "Calculus")] + public sealed class RationalFunctionBesideTheRootOfAQuadraticTest + { + private static Entity DifferentiatesBack(string integrand, double[] points, params (string, double)[] pins) + { + var integral = integrand.ToEntity().Integrate("x"); + Assert.DoesNotContain("integral(", integral.Stringize()); + Assert.DoesNotContain("NaN", integral.Stringize()); + + Entity answer = integral.Substitute("C", 0); + Entity original = integrand.ToEntity(); + foreach (var (name, value) in pins) + { + answer = answer.Substitute(name, value); + original = original.Substitute(name, value); + } + // With the symbols pinned the piecewise's conditions are decidable, and the + // simplification picks the arm before the derivative is taken. + var derivative = answer.Simplify().Differentiate("x"); + + var compared = 0; + foreach (var at in points) + { + var got = derivative.Substitute("x", at).EvalNumerical(); + var want = original.Substitute("x", at).EvalNumerical(); + if (got.IsNaN || want.IsNaN) + continue; + compared++; + var difference = Math.Abs((double)(got - want).RealPart) + + Math.Abs((double)(got - want).ImaginaryPart); + var scale = Math.Max(1.0, Math.Abs((double)want.RealPart)); + Assert.True(difference / scale < 1e-8, + $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, " + + $"where the integrand is {want}"); + } + Assert.True(compared >= 3, + $"only {compared} of {points.Length} points were comparable for {integrand}"); + return integral; + } + + /// + /// The root above the bar and below it, with the linears written or found by + /// factoring; checked on both sides of every pole, where the sign of each + /// x - p is what makes one formula serve. + /// + [Theory] + [InlineData("sqrt(4 + 9*(1 - x^2))/(1 - x^2)", new[] { -1.1, -0.7, -0.3, 0.4, 0.8, 1.1 })] + [InlineData("sqrt(x^2 + x + 3)/(x*(x + 1))", new[] { -3.0, -1.5, -0.5, 0.4, 2.0 })] + [InlineData("1/((x - 1)*(x + 2)*sqrt(x^2 + 1))", new[] { -3.0, -1.0, 0.3, 2.0, 3.0 })] + [InlineData("x^2/((x - 1)*(x + 2)*sqrt(x^2 + 1))", new[] { -3.0, -1.0, 0.3, 2.0, 3.0 })] + [InlineData("sqrt(1 + x^2)/(1 + x)", new[] { -3.0, -1.7, -0.4, 0.9, 2.6 })] + public void TakenApartOverTheRoot(string integrand, double[] points) => DifferentiatesBack(integrand, points); + + /// + /// Moses', with symbols: B^2 over the root and A^2/2 over each linear. And + /// a symbolic quadratic under the root beside a written linear, on both signs of its + /// constant. + /// + [Theory] + [InlineData("-sqrt(A^2 + B^2*(1 - x^2))/(1 - x^2)", new[] { -1.7, -1.3, -0.5, 0.4, 1.3, 1.7 }, 1.3, 0.8, 3.0)] + [InlineData("(x + 1)/(x*sqrt(a*x^2 + b*x + c))", new[] { -2.5, -0.7, 0.4, 1.7 }, 1.0, 1.0, 3.0)] + [InlineData("(x + 1)/(x*sqrt(a*x^2 + b*x + c))", new[] { -4.0, -3.0, 1.6, 2.5 }, 1.0, 1.0, -3.0)] + public void TheCoefficientsMayBeSymbols(string integrand, double[] points, double first, double second, double third) + => DifferentiatesBack(integrand, points, ("A", first), ("B", second), ("a", first), ("b", second), ("c", third)); + } +}