Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 12 additions & 0 deletions BREAKING-CHANGES.md
Original file line number Diff line number Diff line change
Expand Up @@ -447,6 +447,18 @@ case is most of its length ([#718](https://github.com/asc-community/AngouriMath/
| `"x^4/((a + b*x^2 + c*x^4)*sqrt(d + k*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 28,551 characters |
| `"1/(x^2*(a + b*x^2 + c*x^4)*sqrt(d + k*x^2))".ToEntity().Integrate("x")` | `integral(...)` | 10,194 characters |

### A rational function past the twelfth degree is reduced by Hermite's method

**Answers where there were none.** The integrator read a denominator of degree twelve at most, so a derivative over
powers of `x^4 + x + 3` up to the fourth was declined, after seven to fifteen seconds. A rational function with numbers
for coefficients over a denominator past the twelfth degree now takes Hermite's reduction, as long as what is left for the
logarithms stays within twelve ([#718](https://github.com/asc-community/AngouriMath/issues/718)).

| Input | Was (2.5.0) | Now |
|---|---|---|
| `"3*(-47+228*x+120*x^2+19*x^3)/(3+x+x^4)^4+(42-320*x-75*x^2-8*x^3)/(3+x+x^4)^3+30*x/(3+x+x^4)^2".ToEntity().Integrate("x")` | `integral(...)` | 142 characters |
| `"(-3+10*x+4*x^3-30*x^5)/(3+x+x^4)^3-3*(1+4*x^3)*(2-3*x+5*x^2+x^4-5*x^6)/(3+x+x^4)^4".ToEntity().Integrate("x")` | `integral(...)` | 142 characters |

### A power of a linear over a power of a quadratic is reduced a power at a time

**Answers where there were none.** `sqrt(d + k x)/(a + b x + c x^2)^3` took thirty seconds for 26,437 characters on the
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2621,6 +2621,50 @@ node is Powf(var radicand, Number.Rational r) && r is not Number.Integer && radi
/// <see cref="RothsteinTrager"/>. After <see cref="SolveByPartialFractions"/>, so that
/// everything that answers keeps the form it gives.
/// </summary>
internal static Entity? SolveAHighDegreeRationalFunctionByHermite(Entity expr, Entity.Variable x)
{
// A sum of quotients over powers of one denominator, Rubi's spelling of a derivative,
// is one quotient for this. At the question asked only, which is where nothing else
// answered one.
if (!Integration.AnsweringTheQuestionAsked || expr.Vars.Any(symbol => symbol != x))
return null;
// Over the least common multiple of the terms' denominators, each base to its highest
// power, rather than their product: three terms over the fourth, third and second
// powers of one quartic are one quotient over its fourth power.
var terms = Sumf.LinearChildren(expr).Select(term => Functions.SingleQuotient.Of(term)).ToList();
var highest = new Dictionary<Entity, int>();
var powersOf = new List<Dictionary<Entity, int>>();
foreach (var (_, below) in terms)
{
var powers = new Dictionary<Entity, int>();
foreach (var factor in Mulf.LinearChildren(below))
{
var (@base, whole) = factor is Powf(var b0, Number.Integer { EInteger: var w }) && w.Sign > 0 && w.CanFitInInt32()
? (b0, w.ToInt32Unchecked()) : (factor, 1);
powers[@base] = powers.TryGetValue(@base, out var already) ? already + whole : whole;
}
foreach (var pair in powers)
highest[pair.Key] = System.Math.Max(highest.TryGetValue(pair.Key, out var already) ? already : 0, pair.Value);
powersOf.Add(powers);
}
Entity denominator = Number.Integer.One;
foreach (var pair in highest)
denominator = denominator * MathS.Pow(pair.Key, Number.Integer.Create(pair.Value));
Entity numerator = Number.Integer.Zero;
for (var i = 0; i < terms.Count; i++)
{
Entity cofactor = Number.Integer.One;
foreach (var pair in highest)
{
var missing = pair.Value - (powersOf[i].TryGetValue(pair.Key, out var have) ? have : 0);
if (missing > 0)
cofactor = cofactor * MathS.Pow(pair.Key, Number.Integer.Create(missing));
}
numerator = numerator + terms[i].Numerator * cofactor;
}
return RothsteinTrager.IntegrateAHighDegree(numerator, denominator, x);
}

internal static Entity? SolveByRothsteinTrager(Entity expr, Entity.Variable x)
{
if (!TryReadAsQuotient(expr, out var numerator, out var denominator))
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -666,6 +666,9 @@ private static Entity WithItsPiecewiseSumsCombined(Entity expr)
// the question asked only, and the polynomial term below writes `1/u` as `u^(-1)` and
// asks that a level down: `1/(e^atanh(a x) (c - c/(a x))^(3/2))` reached it nowhere.
if ((answer = IndefiniteIntegralSolver.SolveByWritingEachLinearOnce(expr, x, integrateByParts)) is { }) return answer;
// A rational function with numbers for coefficients over a denominator past the twelfth
// degree, by Hermite's reduction before anything searches: nothing below answered one.
if ((answer = IndefiniteIntegralSolver.SolveAHighDegreeRationalFunctionByHermite(expr, x)) is { }) return answer;
// A power of a linear beside the derivative of a quadratic over a whole power of it, by
// parts once, before anything takes the symbolic quadratic apart.
if ((answer = IndefiniteIntegralSolver.SolveAPowerOfALinearBesideTheDerivativeOfAQuadraticOverItsPower(expr, x, integrateByParts)) is { }) return answer;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -63,16 +63,44 @@ internal static class RothsteinTrager
/// </summary>
private const int MaxDegree = 12;

/// <summary>
/// The largest denominator degree read at all: Hermite's reduction is a few Diophantine
/// solves over the rationals, and what it leaves for the resultant is the squarefree part,
/// held to <see cref="MaxDegree"/>. Rubi's 1.3.1 <c>3 (19 x^3 + 120 x^2 + 228 x - 47)/(x^4 + x + 3)^4 + ...</c>
/// is a derivative over a sixteenth-degree denominator with nothing left for the logarithms,
/// and was declined after a minute and a half of searching elsewhere.
/// </summary>
private const int MaxReducedDegree = 24;

/// <summary>
/// The integral of <paramref name="numerator"/> over <paramref name="denominator"/>,
/// both polynomials in <paramref name="x"/> with rational coefficients and the fraction
/// proper; <see langword="null"/> where they are not, or where the answer does not
/// differentiate back to the integrand at the sampled points.
/// </summary>
/// <summary>
/// <see cref="Integrate"/>, for a denominator past <see cref="MaxDegree"/> only: what the
/// rest of the integrator never answered, asked before it searches.
/// </summary>
internal static Entity? IntegrateAHighDegree(Entity numerator, Entity denominator, Variable x)
=> PolynomialFactoring.TryGetRationalCoefficients(denominator, x, 1, 1, MaxReducedDegree, out var belowCoefficients)
&& RationalPolynomial.Create(belowCoefficients).Degree > MaxDegree
? IntegrateUpToDegree(numerator, denominator, x, MaxReducedDegree)
: null;

/// <summary>
/// The integral of <paramref name="numerator"/> over <paramref name="denominator"/>,
/// held to <see cref="MaxDegree"/> unless asked otherwise: below a substitution a
/// denominator past it -- <c>1/(cos(x) + cos(3x))^5</c> under the tangent's half angle --
/// was reduced past the budget where it used to be declined at once.
/// </summary>
internal static Entity? Integrate(Entity numerator, Entity denominator, Variable x)
=> IntegrateUpToDegree(numerator, denominator, x, MaxDegree);

private static Entity? IntegrateUpToDegree(Entity numerator, Entity denominator, Variable x, int degree)
{
if (!PolynomialFactoring.TryGetRationalCoefficients(numerator, x, 1, 0, MaxDegree, out var aboveCoefficients)
|| !PolynomialFactoring.TryGetRationalCoefficients(denominator, x, 1, 1, MaxDegree, out var belowCoefficients))
if (!PolynomialFactoring.TryGetRationalCoefficients(numerator, x, 1, 0, degree, out var aboveCoefficients)
|| !PolynomialFactoring.TryGetRationalCoefficients(denominator, x, 1, 1, degree, out var belowCoefficients))
return null;
var above = RationalPolynomial.Create(aboveCoefficients);
var below = RationalPolynomial.Create(belowCoefficients);
Expand All @@ -87,7 +115,7 @@ internal static class RothsteinTrager
var answer = rationalPart;
if (!logAbove.IsZero)
{
if (LogarithmicPart(logAbove, logBelow, x) is not { } logarithmic)
if (logBelow.Degree > MaxDegree || LogarithmicPart(logAbove, logBelow, x) is not { } logarithmic)
return null;
answer = answer + logarithmic;
}
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,41 @@
//
// 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
{
/// <summary>
/// A rational function with numbers for coefficients over a denominator past the twelfth
/// degree, by Hermite's reduction: Rubi's 1.3.1 derivatives over powers of
/// <c>x^4 + x + 3</c>, written as sums of quotients over different powers of it.
/// <a href="https://github.com/asc-community/AngouriMath/issues/718">#718</a>
/// </summary>
[Trait("Area", "Calculus")]
public sealed class AHighDegreeRationalFunctionIntegralTest
{
[Theory]
[InlineData("3*(-47+228*x+120*x^2+19*x^3)/(3+x+x^4)^4+(42-320*x-75*x^2-8*x^3)/(3+x+x^4)^3+30*x/(3+x+x^4)^2")]
[InlineData("(-3+10*x+4*x^3-30*x^5)/(3+x+x^4)^3-3*(1+4*x^3)*(2-3*x+5*x^2+x^4-5*x^6)/(3+x+x^4)^4")]
public void IsAnsweredAndDifferentiatesBack(string integrand)
{
var integral = integrand.ToEntity().Integrate("x");
Assert.DoesNotContain("integral(", integral.Stringize());
var derivative = integral.Substitute("C", 0).Differentiate("x");
foreach (var at in new[] { -2.3, -0.7, 0.4, 1.9 })
{
var got = derivative.Substitute("x", at).EvalNumerical();
var want = integrand.ToEntity().Substitute("x", at).EvalNumerical();
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}");
}
}
}
}
Loading