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
17 changes: 17 additions & 0 deletions BREAKING-CHANGES.md
Original file line number Diff line number Diff line change
Expand Up @@ -552,6 +552,23 @@ power of the variable every term holds comes out the same way.
| `"x^2/((a*x + a)*(1 - x^2))".Integrate("x")` | unevaluated | the partial fractions over `(x + 1)^2 (1 - x)`, over `a` |
| `"cot(x)^3/(a + b*csc(x))".Integrate("x")` | unevaluated | `-1/(b sin(x)) - a ln(sin(x))/b^2 + (b^2 - a^2) ln(a sin(x) + b)/(a b^2)`, once `a u^4 + b u^3` is read as `u^3 (a u + b)` |

### Symbolic factors of a denominator that share a factor are written over it

**Improvement, not silent.** Partial fractions over factors with symbols in them read the written
factors as coprime and squarefree, and the refactoring after them reads rational coefficients only.
A quadratic that shares a root with a linear beside it, or that is a square, was read as an
irreducible quadratic and declined. The written factors are now taken apart over their greatest
common divisors, in every variable, with one another and each with its derivative, and a factor
shared with the numerator is cancelled. Rubi's 1.2.1.2
([#718](https://github.com/asc-community/AngouriMath/issues/718)).

| Input | Was (2.5.0) | Now |
|---|---|---|
| `"1/((d+e*x)*(a*d*e+(c*d^2+a*e^2)*x+c*d*e*x^2)^2)".Integrate("x")` | unevaluated | the partial fractions over `(d + e x)^3 (a e + c d x)^2` |
| `"(d+e*x)^6/(a*d*e+(c*d^2+a*e^2)*x+c*d*e*x^2)^4".Integrate("x")` | no answer in 90 s | the antiderivative of `(d + e x)^2/(a e + c d x)^4` |
| `"1/((d+e*x)^2*(c*d^2+2*c*d*e*x+c*e^2*x^2)^3)".Integrate("x")` | unevaluated | `-1/(7 c^3 e^8 (d/e + x)^7)` |
| `"x^2/((a+b*x)*(a^2+2*a*b*x+b^2*x^2))".Integrate("x")` | unevaluated | the partial fractions over `(a + b x)^3` |

### Fractional powers of `a ± a sin` are integrated by the half angle at which they are squares

**Improvement, not silent.** `(A + C sin(y)^2)/((c - c sin(y))^(3/2) sqrt(a + a sin(y)))` and
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -665,6 +665,19 @@ is var (multiple, leftover)
?? Integration.ComputeIndefiniteIntegral(numerator / primitive, x, integrateByParts)) is { } overPrimitives)
return overPrimitives;

// Written factors with a symbol in them that share a factor, or repeat one inside,
// written over what they share: `(d + e x)(a d e + (c d^2 + a e^2) x + c d e x^2)^2`
// is `(d + e x)^3 (a e + c d x)^2`, and `a^2 + 2 a b x + b^2 x^2` beside `a + b x`
// is a third power of it. The splits below read the written factors as coprime and
// squarefree, and the refactoring after them reads rational coefficients only:
// Rubi's 1.2.1.2 `1/((d + e x)(a d e + (c d^2 + a e^2) x + c d e x^2)^2)` was declined
// after four seconds, where written so it is answered in a tenth of one. Once: what
// this writes is coprime and squarefree already.
if (OverSymbolicFactorsWrittenApart(numerator, denominator, x) is var (above, below)
&& (SolveByPartialFractions(above / below, x, integrateByParts)
?? Integration.ComputeIndefiniteIntegral(above / below, x, integrateByParts)) is { } overWhatTheyShare)
return overWhatTheyShare;

// Written linear factors with symbols in their coefficients, two or more, one of
// them to a power: decomposed over the written factors, the coefficients read
// off derivatives at the roots and each a line. First, before the respellings
Expand Down Expand Up @@ -10275,6 +10288,221 @@ private static bool IsASumOfMonomials(Entity expr, Entity.Variable x)
return changed ? product : null;
}

/// <summary>
/// <paramref name="numerator"/> over <paramref name="denominator"/>, the denominator's
/// written factors taken apart over their greatest common divisors -- with one another,
/// and each with its own derivative in <paramref name="x"/> -- until they are coprime and
/// squarefree as written, and what the numerator shares with them cancelled;
/// <see langword="null"/> where no factor holds a symbol, they are coprime and squarefree
/// already, or one is not a polynomial.
/// </summary>
/// <remarks>
/// The divisors are taken over the rationals in every variable, so a factor of the symbols
/// alone can come with one: that goes in front, as a constant, and each factor kept is
/// primitive in <paramref name="x"/>. `c d^2 + 2 c d e x + c e^2 x^2` shares `c (d + e x)`
/// with its derivative, and is `c (d + e x)^2`.
/// </remarks>
private static (Entity Numerator, Entity Denominator)? OverSymbolicFactorsWrittenApart(Entity numerator, Entity denominator, Entity.Variable x)
{
var written = new List<(Entity Base, int Power)>();
Entity constant = Number.Integer.One;
foreach (var factor in Mulf.LinearChildren(denominator))
{
if (!factor.ContainsNode(x))
{
constant = constant == Number.Integer.One ? factor : constant * factor;
continue;
}
var (@base, power) = factor is Powf(var b, Number.Integer e) && e.EInteger.Sign > 0 && e.EInteger.CanFitInInt32()
? (b, e.EInteger.ToInt32Unchecked())
: (factor, 1);
// Over one bar, the constant below it in front: the integrator writes a quadratic
// monic, `a/c + (c d^2 + a e^2)/(c d e) x + x^2`, and it is a polynomial in the
// symbols over `c d e`.
var (overOneBar, belowTheBar) = Functions.SingleQuotient.Of(Functions.SingleQuotient.Combine(@base));
if (belowTheBar.ContainsNode(x))
return null;
if (belowTheBar != Number.Integer.One)
{
var front = power == 1 ? belowTheBar : MathS.Pow(belowTheBar, power);
constant = constant == Number.Integer.One ? 1 / front : constant / front;
}
written.Add((overOneBar, power));
}
if (!written.Any(pair => pair.Base.Vars.Any(v => v != x)))
return null;
var (numeratorOverOneBar, belowTheNumerator) = Functions.SingleQuotient.Of(Functions.SingleQuotient.Combine(numerator));
if (!belowTheNumerator.ContainsNode(x))
{
numerator = numeratorOverOneBar;
if (belowTheNumerator != Number.Integer.One)
constant = constant == Number.Integer.One ? belowTheNumerator : constant * belowTheNumerator;
}
var variables = numerator.Vars.Concat(denominator.Vars).Distinct().OrderBy(v => v.Name, System.StringComparer.Ordinal).ToList();
if (variables.Count > MultivariatePolynomial.MaxVariables)
return null;
var indices = new Dictionary<Variable, int>();
for (var i = 0; i < variables.Count; i++)
indices[variables[i]] = i;
var all = Enumerable.Range(0, variables.Count).ToList();
var at = indices[x];
var rest = all.Where(i => i != at).ToList();
var parts = new List<(MultivariatePolynomial Poly, int Power)>();
var degree = 0;
foreach (var (@base, power) in written)
{
if (MultivariatePolynomial.TryParse(@base, indices) is not { } poly || poly.DegreeIn(at) == 0)
return null;
degree += poly.DegreeIn(at);
parts.Add((poly, power));
}
if (degree > 12 || parts.Count < 2 && parts[0].Poly.DegreeIn(at) < 2)
return null;

// A power of x is what a factor shares with x^k beside it, or with its own derivative
// where its lowest term is x^k, and taking it out is the content's, before this. Split
// here, `b x^2 + c x^4` came apart as `x^2 (b/c + x^2) c`, and over that monic quadratic
// with a symbol below its bar Rubi's `x^12 (A + B x^2)/(b x^2 + c x^4)^3` had an answer
// of four million characters, where with the content taken out it has seven thousand.
bool IsAPowerOfX(MultivariatePolynomial poly)
=> poly.Monomials.Take(2).Count() == 1
&& rest.All(variable => MultivariatePolynomial.PowerOfMonomial(poly.Monomials.First(), variable) == 0);

// What is free of x in a factor goes in front, to the factor's power.
void Front(MultivariatePolynomial free, int power)
{
if (free.IsConstant && free.SameAs(MultivariatePolynomial.One(free.VariableCount)))
return;
var front = free.ToEntity(variables);
var raised = power == 1 ? front : MathS.Pow(front, power);
constant = constant == Number.Integer.One ? raised : constant * raised;
}
bool Keep(MultivariatePolynomial poly, int power)
{
// A part free of x goes in front whole, its sign with it: a divisor is found up to
// a unit, so taking -(b u^2 - 1)^2 apart leaves a part -1, whose content is 1.
if (poly.DegreeIn(at) == 0)
{
Front(poly, power);
return true;
}
if (Functions.PolynomialGcd.ContentIn(poly, at, rest, 0) is not { } content || poly.DivideExact(content) is not { } primitive)
return false;
Front(content, power);
parts.Add((primitive, power));
return true;
}

var changed = false;
for (var round = 0; round < 64; round++)
{
var split = false;
for (var i = 0; i < parts.Count && !split; i++)
{
var (poly, power) = parts[i];
// A factor repeated inside: what it shares with its derivative.
if (poly.DegreeIn(at) >= 2 && Functions.PolynomialGcd.Gcd(poly, poly.DerivativeIn(at), all, 0) is { } repeated
&& repeated.DegreeIn(at) >= 1 && !IsAPowerOfX(repeated) && poly.DivideExact(repeated) is { } remaining)
{
parts.RemoveAt(i);
if (!Keep(repeated, power) || !Keep(remaining, power))
return null;
split = true;
break;
}
// Two factors that share one.
for (var j = i + 1; j < parts.Count; j++)
{
if (Functions.PolynomialGcd.Gcd(poly, parts[j].Poly, all, 0) is not { } common || common.DegreeIn(at) == 0 || IsAPowerOfX(common))
continue;
if (poly.DivideExact(common) is not { } left || parts[j].Poly.DivideExact(common) is not { } right)
return null;
var other = parts[j].Power;
parts.RemoveAt(j);
parts.RemoveAt(i);
if (!Keep(common, power + other) || !Keep(left, power) || !Keep(right, other))
return null;
split = true;
break;
}
}
if (!split)
break;
changed = true;
}

// A factor that shares one with the numerator, taken apart the same way:
// `(d + e x)^6/(a d e + (c d^2 + a e^2) x + c d e x^2)^4` is `(d + e x)^2/(a e + c d x)^4`.
var top = MultivariatePolynomial.TryParse(numerator, indices);
for (var round = 0; top is not null && round < 16; round++)
{
var split = false;
for (var i = 0; i < parts.Count; i++)
{
var (poly, power) = parts[i];
if (poly.DegreeIn(at) < 2 || Functions.PolynomialGcd.Gcd(top, poly, all, 0) is not { } common
|| common.DegreeIn(at) == 0 || IsAPowerOfX(common) || common.DegreeIn(at) == poly.DegreeIn(at)
|| poly.DivideExact(common) is not { } remaining)
continue;
parts.RemoveAt(i);
if (!Keep(common, power) || !Keep(remaining, power))
return null;
split = true;
break;
}
if (!split)
break;
changed = true;
}
if (!changed)
return null;

// What the numerator shares with the factors, cancelled.
Entity above = numerator;
if (top is not null)
{
var cancelled = false;
for (var i = 0; i < parts.Count; i++)
while (parts[i].Power > 0 && top.DivideExact(parts[i].Poly) is { } fewer)
{
top = fewer;
parts[i] = (parts[i].Poly, parts[i].Power - 1);
cancelled = true;
}
if (cancelled)
above = top.ToEntity(variables);
}
// Each factor monic in x, the leading coefficient to its power in front, which is how
// the integrator writes a factor: `d + e x` is `e (d/e + x)`, and written the other
// way the next normalisation took `e` out again and the answer carried `e^(-8) e^8`.
var monicFactors = new List<Entity>();
foreach (var (poly, power) in parts)
{
if (power == 0)
continue;
var lead = poly.LeadingCoefficientIn(at).ToEntity(variables);
Entity monic = Number.Integer.Zero;
foreach (var pair in poly.CoefficientsIn(at).OrderBy(pair => pair.Key))
{
var coefficient = pair.Key == poly.DegreeIn(at) ? Number.Integer.One : Reduced(pair.Value.ToEntity(variables) / lead);
var term = pair.Key == 0 ? coefficient : coefficient * (pair.Key == 1 ? x : MathS.Pow(x, pair.Key));
monic = monic == Number.Integer.Zero ? term : monic + term;
}
monicFactors.Add(power == 1 ? monic : MathS.Pow(monic, power));
var raised = power == 1 ? lead : MathS.Pow(lead, power);
constant = constant == Number.Integer.One ? raised : constant * raised;
}
// The constant in lowest terms: the bars and the contents taken out of the factors
// and the leading coefficients put back leave `e^(-8) e^8` in it as written.
Entity below = Reduced(constant);
foreach (var factor in monicFactors)
below = below == Number.Integer.One ? factor : below * factor;
return (above, below);

static Entity Reduced(Entity constant)
=> constant.Vars.Any() ? Functions.PartialFractions.InLowestTermsOverTheSymbols(constant) : constant.InnerSimplified;
}

/// <summary>
/// <paramref name="denominator"/> with every written factor that has rational coefficients
/// and splits over the rationals written in its factors, where another written factor
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -278,6 +278,35 @@ public void TheCoefficientsAreInLowestTerms(string integrand, int atMost)
public void ProportionalLinearFactorsAreOneFactor(string integrand)
=> DifferentiatesBack(integrand, "x", ("a", 1.7), ("b", 2.3), ("c", 0.4), ("d", 1.1), ("g", 0.7));

/// <summary>
/// Written factors that share a factor, or repeat one inside, are written over what they
/// share, the greatest common divisors taken in every variable: Rubi's 1.2.1.2,
/// <c>(d + e x)^m (a d e + (c d^2 + a e^2) x + c d e x^2)^p</c>, whose quadratic is
/// <c>(d + e x)(a e + c d x)</c>, and the square <c>c d^2 + 2 c d e x + c e^2 x^2</c>
/// beside a power of its linear, or with the numerator sharing it. Each was declined or
/// ran past the corpus's budget.
/// </summary>
[Theory]
[InlineData("1/((d + h*x)*(a*d*h + (c*d^2 + a*h^2)*x + c*d*h*x^2)^2)")]
[InlineData("(d + h*x)^6/(a*d*h + (c*d^2 + a*h^2)*x + c*d*h*x^2)^4")]
[InlineData("1/((d + h*x)^2*(c*d^2 + 2*c*d*h*x + c*h^2*x^2)^3)")]
[InlineData("x^2/((a + b*x)*(a^2 + 2*a*b*x + b^2*x^2))")]
public void FactorsThatShareOneAreWrittenOverIt(string integrand)
=> DifferentiatesBack(integrand, "x", ("a", 1.7), ("b", 2.3), ("c", 0.4), ("d", 1.1), ("h", 0.7));

/// <summary>
/// The negative of a square, <c>-(b x^2 - 1)^2</c>, is the square with its sign in front:
/// a divisor is found up to a unit, and the part the split left was <c>-1</c>. Rubi's
/// 1.1.2.2 <c>x^4/(a + b x^2)^(3/2)</c> is this under <c>t = x/sqrt(a + b x^2)</c>, and was
/// answered with the wrong sign.
/// </summary>
[Theory]
[InlineData("(2*x^2/b - 1/b^2)/(-b^2*x^4 + 2*b*x^2 - 1)")]
[InlineData("x^4/(a + b*x^2)^(3/2)")]
[InlineData("(c + d*x^3)^2/(a + b*x^3)^(4/3)")]
public void TheNegativeOfASquareKeepsItsSign(string integrand)
=> DifferentiatesBack(integrand, "x", ("a", 1.7), ("b", 2.3), ("c", 0.4), ("d", 1.1));

/// <summary>
/// A written power of a linear is integrated as the power it is: expanded,
/// <c>(f/g + x)^2</c> is a quadratic whose discriminant is zero in a spelling the
Expand Down
Loading