diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IntegralPatterns.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IntegralPatterns.cs index 6211b028c..dcd697031 100644 --- a/Sources/AngouriMath/Functions/Continuous/Integration/IntegralPatterns.cs +++ b/Sources/AngouriMath/Functions/Continuous/Integration/IntegralPatterns.cs @@ -124,6 +124,14 @@ _ when TryReadSineCosinePowers(expr, x, out var trigArg, out var sinePower, out && qa.Evaled is Entity.Number.Complex { IsZero: false } => IntegrateRootOfQuadratic(qa, qb, qc, radicand, x), + // ∫ k / (x^2 * sqrt(ax^2 + c)) dx, the shape a trigonometric substitution is + // usually taught for. Differentiating sqrt(ax^2 + c)/x gives exactly + // -c/(x^2 * sqrt(ax^2 + c)), so every sign of a and c is the one formula and + // there is no case analysis to get wrong. Without it 1/(x^2 * sqrt(x^2 - 1)) + // had no antiderivative at all. + _ when TryReadOverSquareTimesRoot(expr, x, out var overFactor, out var overRadicand, out var overConstant) + => -overFactor * MathS.Sqrt(overRadicand) / (overConstant * x), + // ∫ k / sqrt(ax^2 + bx + c) dx -- the arcsine and logarithm forms. Without // these, 1/sqrt(1 - x^2) had no antiderivative at all. Entity.Divf(var numerator, @@ -273,6 +281,48 @@ private static bool IsExponentialRate(Entity expr, Entity.Variable x, out Entity /// is the square root of something linear, which the ordinary power rule already /// integrates, and dividing by a would not be allowed anyway. /// + /// + /// Reads k / (x^2 * sqrt(ax^2 + c)), giving back k, the radicand and its + /// constant term. The radicand has to be a quadratic in x with no linear term and + /// with neither of its two coefficients zero: a zero constant makes the formula + /// below divide by it, and with a zero a there is no root of x left to speak of. + /// + private static bool TryReadOverSquareTimesRoot( + Entity expr, Entity.Variable x, + out Entity factor, out Entity radicand, out Entity constantTerm) + { + factor = radicand = constantTerm = 0; + if (expr is not Entity.Divf(var numerator, var denominator) || numerator.ContainsNode(x)) + return false; + Entity coefficient = numerator; + var squares = 0; + Entity? root = null, constant = null; + foreach (var part in Entity.Mulf.LinearChildren(denominator)) + switch (part) + { + case Entity.Powf(var square, Entity.Number.Integer(2)) when square == x: + squares++; + break; + case Entity.Powf(var under, Entity.Number.Rational(Entity.Number.Integer(1), Entity.Number.Integer(2))) + when root is null + && TreeAnalyzer.TryGetPolyQuadratic(under, x, out var qa, out var qb, out var qc) + && qa.Evaled is Entity.Number.Complex { IsZero: false } + && qb.Evaled is Entity.Number.Complex { IsZero: true } + && qc.Evaled is Entity.Number.Complex { IsZero: false }: + (root, constant) = (under, qc); + break; + case var other when !other.ContainsNode(x): + coefficient /= other; + break; + default: + return false; + } + if (squares != 1 || root is null || constant is null) + return false; + (factor, radicand, constantTerm) = (coefficient, root, constant); + return true; + } + private static Entity IntegrateRootOfQuadratic( Entity a, Entity b, Entity c, Entity radicand, Entity.Variable x) => (2 * a * x + b) * MathS.Sqrt(radicand) / (4 * a) diff --git a/Sources/Tests/UnitTests/Calculus/RootOverSquareIntegralsTest.cs b/Sources/Tests/UnitTests/Calculus/RootOverSquareIntegralsTest.cs new file mode 100644 index 000000000..d3ca736a2 --- /dev/null +++ b/Sources/Tests/UnitTests/Calculus/RootOverSquareIntegralsTest.cs @@ -0,0 +1,80 @@ +// +// Copyright (c) 2019-2022 Angouri. +// AngouriMath is licensed under MIT. +// Details: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md. +// Website: https://am.angouri.org. +// + +using AngouriMath; +using AngouriMath.Extensions; +using Xunit; + +namespace AngouriMath.Tests.Calculus +{ + /// + /// k / (x^2 * sqrt(ax^2 + c)), the shape a trigonometric substitution is usually + /// taught for. It had no antiderivative at all. Named in the issue's own list of what + /// is missing: https://github.com/asc-community/AngouriMath/issues/233. + /// Each answer is checked by differentiating it back and comparing at points, since + /// what matters is that it is an antiderivative and not what form it is written in. + /// + public sealed class RootOverSquareIntegralsTest + { + private static void AssertIsAntiderivative(string integrand, params double[] points) + { + var f = integrand.ToEntity(); + var antiderivative = f.Integrate("x"); + Assert.DoesNotContain("integral(", antiderivative.Stringize()); + var derivative = antiderivative.Substitute("C", 0).Differentiate("x"); + foreach (var point in points) + { + var expected = f.Substitute("x", point).EvalNumerical().RealPart.EDecimal.ToDouble(); + var actual = derivative.Substitute("x", point).EvalNumerical().RealPart.EDecimal.ToDouble(); + Assert.Equal(expected, actual, 7); + } + } + + // Differentiating sqrt(ax^2 + c)/x gives -c/(x^2 * sqrt(ax^2 + c)), so every sign + // of a and c is the same formula: the root of a sum, of a difference, and the two + // with the sign of x^2 the other way round. + [Theory] + [InlineData("1 / (x ^ 2 * sqrt(x ^ 2 - 1))", new[] { 1.4, 2.6, 4.1 })] + [InlineData("1 / (x ^ 2 * sqrt(x ^ 2 + 1))", new[] { 0.4, 1.6, 3.1, -2.2 })] + [InlineData("1 / (x ^ 2 * sqrt(1 - x ^ 2))", new[] { 0.4, 0.8, -0.6 })] + [InlineData("1 / (x ^ 2 * sqrt(4 - x ^ 2))", new[] { 0.4, 1.6, -1.2 })] + [InlineData("1 / (x ^ 2 * sqrt(2 * x ^ 2 + 3))", new[] { 0.4, 1.6, -2.2 })] + public void EverySignOfTheQuadraticUnderTheRoot(string integrand, double[] points) => + AssertIsAntiderivative(integrand, points); + + // A constant anywhere in the quotient is carried through rather than refused. + [Theory] + [InlineData("3 / (x ^ 2 * sqrt(x ^ 2 - 4))", new[] { 2.4, 3.6 })] + [InlineData("1 / (2 * x ^ 2 * sqrt(x ^ 2 + 9))", new[] { 0.4, 1.6 })] + [InlineData("(-1) / (x ^ 2 * sqrt(x ^ 2 + 1))", new[] { 0.4, 1.6 })] + public void AConstantFactorIsCarriedThrough(string integrand, double[] points) => + AssertIsAntiderivative(integrand, points); + + // The shapes this sits next to have to keep working. + [Theory] + [InlineData("1 / sqrt(x ^ 2 - 1)", new[] { 1.4, 2.6 })] + [InlineData("1 / sqrt(1 - x ^ 2)", new[] { 0.4, 0.8 })] + [InlineData("sqrt(x ^ 2 + 1)", new[] { 0.4, 1.6 })] + [InlineData("1 / x ^ 2", new[] { 0.4, 1.6 })] + [InlineData("1 / (x ^ 2 + 1)", new[] { 0.4, 1.6 })] + [InlineData("x / sqrt(x ^ 2 + 1)", new[] { 0.4, 1.6 })] + [InlineData("sin(x)", new[] { 0.4, 1.6 })] + public void NeighbouringFormsAreUnaffected(string integrand, double[] points) => + AssertIsAntiderivative(integrand, points); + + /// + /// Outside the family the formula does not apply and nothing is claimed: a power of + /// x other than the square, and a radicand whose constant term is zero, which the + /// formula would divide by. + /// + [Theory] + [InlineData("1 / (x ^ 3 * sqrt(x ^ 2 - 1))")] + [InlineData("1 / (x ^ 2 * sqrt(x ^ 2 + x + 1))")] + public void OutsideTheFamilyNothingIsClaimed(string integrand) => + Assert.Contains("integral(", integrand.ToEntity().Integrate("x").Stringize()); + } +}