Skip to content

Commit 5dc25a6

Browse files
Rafael-SOWNetclaude
andcommitted
Write a biquadratic's terms by the sign of each root where both roots are real
Over a + b x^2 + c x^4, 1/(x^2 - r) was written atan(x/sqrt(-r))/sqrt(-r) for every root r. Where the discriminant is negative the roots are conjugate and the two arctangents sum to a real function, but where it is a positive number times even powers and the coefficients are real, both roots are real, and for a positive root the arctangent's argument is imaginary: on its cut wherever x^2 > r, where its value is decided by the sign of a rounding residue. There each term is now written by the sign of its root, as the quadratic rule writes 1/(x^2 - r): the arctangent for a negative root and ln((x - sqrt(r))/(x + sqrt(r)))/(2 sqrt(r)) for a positive one. The discriminant's square root is taken by halving its powers -- any root serves, the other swapping the two -- since sqrt(4 b^2) is 2 |b|, past which no root's sign can be read: cosh(x)^2/(a + b sech(x)) has -1 for a root with 2 b, and came to 6,909 characters with 2 |b|, 1,329 now. Where the roots' product a/c is a positive number times even powers too, the two share a sign, that of their sum, and one piecewise on it holds both terms. Over coefficients with the imaginary unit in them the roots are not real whatever the discriminant, and the arctangents are kept: Rubi's 1.1.1.2 1/((a - i a x)^(3/4) (a + i a x)^(1/4)) was answered NaN by a piecewise none of whose arms held. Rubi's 7.4.2 e^(3 acoth(a x)) sqrt(c - c/(a x))/x is u^4/((a - u^2)(2 a - u^2)) under u = sqrt((a x + 1)/x), discriminant a^2. The rule that writes each linear once checks its arms in interval arithmetic, and between 0 and 1/a both arctangents were on their cuts, so it declined and the answer that remained had no value there. Part of #718. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura
1 parent fc5d732 commit 5dc25a6

3 files changed

Lines changed: 175 additions & 5 deletions

File tree

‎BREAKING-CHANGES.md‎

Lines changed: 17 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -135,6 +135,23 @@ shares a power of `1 - w` with its numerator, and is left to the rules that answ
135135
| `"csch(c+d*x)/(a-b*sinh(c+d*x)^4)^2".ToEntity().Integrate("x")` | `integral(...)` | 2,511 characters |
136136
| `"1/((d+k*x^2)*(a+b*x^2+c*x^4)^2)".ToEntity().Integrate("x")` | `integral(...)` | 9,968 characters, in 0.9 s, where before this change master gave 18,154 in 11 s and 10.8 GB |
137137

138+
### A biquadratic with two real roots in `x^2` is integrated by the sign of each root
139+
140+
**Answers where there were none.** Over `a + b x^2 + c x^4`, each root `r` in `x^2` was written `atan(x/sqrt(-r))/sqrt(-r)`.
141+
Where the discriminant is a positive number times even powers, as `a^2` is, both roots are real, and for a positive
142+
root that arctangent has an imaginary argument, on its branch cut wherever `x^2 > r`. Each term is now written by the
143+
sign of its root, as `1/(x^2 - r)` is: the arctangent for a negative root, and `ln((x - sqrt(r))/(x + sqrt(r)))/(2 sqrt(r))`
144+
for a positive one, with the discriminant's square root taken by halving its powers, so that `4 b^2` gives `2 b` and
145+
not `2 |b|`. Where the product of the roots is a positive number times even powers too, the roots share a sign, and
146+
one piecewise holds both terms. Rubi's 7.4.2 `e^(3 acoth(a x)) sqrt(c - c/(a x))/x` had no value between 0 and
147+
`1/a`, where it is real ([#718](https://github.com/asc-community/AngouriMath/issues/718)).
148+
149+
| Input | Was (2.5.0) | Now |
150+
|---|---|---|
151+
| `"e^(3*acoth(a*x))*sqrt(c-c/(a*x))/x".ToEntity().Integrate("x")` | `integral(...)` | 3,677 characters, where before this change master gave 1,869 with no value between 0 and `1/a` |
152+
| `"1/((a-x^2)*(2*a-x^2))".ToEntity().Integrate("x")` | `integral(...)` | 289 characters, by the sign of `a`, where before this change master gave 218 with both arctangents on their cuts past `x^2 = 2 a` |
153+
| `"x^4/((a-x^2)*(2*a-x^2))".ToEntity().Integrate("x")` | `integral(...)` | 413 characters, where before this change master gave 332 |
154+
138155
### A polynomial over a symbolic quartic that is two quadratics in a hidden variable is split over them
139156

140157
**Answers where there were none.** `1/(a + b x + c x^2 + b x^3 + a x^4)` was declined, and so were Rubi's 1.3.1

‎Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs‎

Lines changed: 63 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -4820,9 +4820,20 @@ private static (Entity[] Polynomial, Entity[] Pole, Entity D, Entity E) SplitOve
48204820
/// </summary>
48214821
private static Entity ArctangentsOverABiquadratic(Entity total, Entity d, Entity e, Entity a, Entity b, Entity c, Entity discriminant, Entity.Variable x)
48224822
{
4823-
var q = MathS.Sqrt(discriminant);
4824-
var firstRoot = (-b + q) / (2 * c);
4825-
var secondRoot = (-b - q) / (2 * c);
4823+
// Real only over real coefficients: Rubi's 1.1.1.2 `1/((a - i a x)^(3/4) (a + i a x)^(1/4))` comes to a
4824+
// biquadratic with the imaginary unit in it and a positive discriminant, whose roots are complex,
4825+
// and the piecewise on their signs had no arm that held -- NaN.
4826+
var realRoots = SignOfANumberTimesEvenPowers(discriminant) is > 0
4827+
&& !new[] { a, b, c }.Any(coefficient => HoldsTheImaginaryUnit(coefficient.InnerSimplified));
4828+
// Any square root of the discriminant serves, the other swapping the two roots: of a number times
4829+
// even powers, the powers halved, since `sqrt(4 b^2)` is `2 |b|`, past which neither root's sign
4830+
// can be read, where with `2 b` one of `cosh(x)^2/(a + b sech(x))`'s roots is -1.
4831+
var rootOfTheDiscriminant = realRoots ? ARootOfANumberTimesEvenPowers(discriminant) : null;
4832+
var q = rootOfTheDiscriminant ?? MathS.Sqrt(discriminant);
4833+
Entity Root(Entity root) => rootOfTheDiscriminant is not null && !rootOfTheDiscriminant.Nodes.Any(node => node is Powf(_, Number.Rational and not Number.Integer))
4834+
? LowestOverTheSymbols(root) : root;
4835+
var firstRoot = Root((-b + q) / (2 * c));
4836+
var secondRoot = Root((-b - q) / (2 * c));
48264837
// `1/(x^2 - r)` is `atan(x/s)/s` with `s = sqrt(-r)` for every complex `r` but zero:
48274838
// `d/dx atan(x/s)/s` is `1/(s^2 + x^2)` whatever `s` is, and for a positive `r` the
48284839
// arctangent of an imaginary argument is the hyperbolic one, `-atanh(x/sqrt(r))/sqrt(r)`,
@@ -4831,16 +4842,43 @@ private static Entity ArctangentsOverABiquadratic(Entity total, Entity d, Entity
48314842
// case with real coefficients -- and the two terms then sum to a real function; a
48324843
// piecewise on the sign of a root, as the quadratic rule would write, has no
48334844
// value there at all.
4845+
// Where the discriminant is a positive number times even powers, both roots are real,
4846+
// and each term is written by the sign of its root, as the quadratic rule writes
4847+
// `1/(x^2 - r)`: for a positive `r` and `x^2 > r` the arctangent's argument is on its
4848+
// cut, the imaginary axis past `i`, where only the sign of a rounding residue decides
4849+
// its value. Rubi's 7.4.2 `e^(3 acoth(a x)) sqrt(c - c/(a x))/x` comes to
4850+
// `u^4/((a - u^2)(2 a - u^2))` under `u = sqrt((a x + 1)/x)`, whose discriminant is
4851+
// `a^2`, and between 0 and `1/a` both arctangents were on their cuts throughout.
4852+
// And where their product `a/c` is a positive number times even powers too, the two
4853+
// roots have one sign, that of their sum `-b/c`, and one piecewise on it holds both
4854+
// terms, where two would combine into four arms, two of which never hold.
4855+
// https://github.com/asc-community/AngouriMath/issues/718
4856+
var oneSign = realRoots && SignOfANumberTimesEvenPowers(a / c) is > 0;
4857+
Entity? arctangents = null, logarithms = null;
48344858
foreach (var (root, sign) in new[] { (firstRoot, 1), (secondRoot, -1) })
48354859
{
48364860
var coefficient = ((d + e * root) / q).InnerSimplified;
48374861
if (coefficient == Number.Integer.Zero)
48384862
continue;
48394863
var s = MathS.Sqrt(-root);
48404864
var part = MathS.Arctan(x / s) / s;
4841-
total += sign == 1 ? coefficient * part : -coefficient * part;
4865+
if (!realRoots)
4866+
{
4867+
total += sign == 1 ? coefficient * part : -coefficient * part;
4868+
continue;
4869+
}
4870+
var t = MathS.Sqrt(root);
4871+
var logarithm = MathS.Ln((x - t) / (x + t)) / (2 * t);
4872+
var signed = sign == 1 ? coefficient : -coefficient;
4873+
if (oneSign)
4874+
{
4875+
arctangents = arctangents is null ? signed * part : arctangents + signed * part;
4876+
logarithms = logarithms is null ? signed * logarithm : logarithms + signed * logarithm;
4877+
}
4878+
else
4879+
total += signed * BySign(-root, part, logarithm);
48424880
}
4843-
return total;
4881+
return arctangents is not null && logarithms is not null ? total + BySign((b / c).InnerSimplified, arctangents, logarithms) : total;
48444882
}
48454883

48464884
/// <summary>
@@ -17563,6 +17601,26 @@ private static Entity BySigns(Entity first, Entity second, Entity wherePositiveA
1756317601
(whereNegative, new Lessf(first, Number.Integer.Zero)));
1756417602
}
1756517603

17604+
// A square root of a positive number times even powers of symbols, `sqrt(k)` times each symbol to half
17605+
// its power; null where the quantity is anything else.
17606+
private static Entity? ARootOfANumberTimesEvenPowers(Entity quantity)
17607+
{
17608+
Entity root = Number.Integer.One;
17609+
foreach (var factor in Mulf.LinearChildren(quantity.InnerSimplified))
17610+
switch (factor)
17611+
{
17612+
case Number.Real { IsZero: false, IsNegative: false } number:
17613+
root *= MathS.Sqrt(number);
17614+
break;
17615+
case Powf(Variable symbol, Number.Integer power) when power.EInteger.IsEven && power.EInteger.Sign > 0 && power.EInteger.CanFitInInt32():
17616+
root *= WholePower(symbol, power.EInteger.ToInt32Unchecked() / 2);
17617+
break;
17618+
default:
17619+
return null;
17620+
}
17621+
return root.InnerSimplified;
17622+
}
17623+
1756617624
private static int? SignOfANumberTimesEvenPowers(Entity quantity)
1756717625
{
1756817626
var sign = 1;
Lines changed: 95 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,95 @@
1+
//
2+
// Copyright (c) 2019-2026 Angouri.
3+
// AngouriMath is licensed under MIT.
4+
// Details: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.
5+
// Website: https://am.angouri.org.
6+
//
7+
8+
using System;
9+
using AngouriMath.Extensions;
10+
using Xunit;
11+
12+
namespace AngouriMath.Tests.Calculus
13+
{
14+
/// <summary>
15+
/// Over a biquadratic whose discriminant is a positive number times even powers, both roots in
16+
/// <c>x^2</c> are real, and each term is written by the sign of its root: for a positive root
17+
/// <c>r</c> the arctangent <c>atan(x/sqrt(-r))</c> has an imaginary argument, on its cut wherever
18+
/// <c>x^2 > r</c>. Rubi's 7.4.2 <c>e^(3 acoth(a x)) sqrt(c - c/(a x))/x</c> comes to
19+
/// <c>u^4/((a - u^2)(2 a - u^2))</c> under <c>u = sqrt((a x + 1)/x)</c>, and its answer had no value
20+
/// between 0 and <c>1/a</c>, where the integrand is real.
21+
/// <a href="https://github.com/asc-community/AngouriMath/issues/718">#718</a>
22+
/// </summary>
23+
[Trait("Area", "Calculus")]
24+
public sealed class ABiquadraticWithRealRootsInTheSquareIntegralTest
25+
{
26+
[Theory]
27+
// Rubi 7.4.2:592.
28+
[InlineData("e^(3*acoth(a*x))*sqrt(c-c/(a*x))/x", new[] { 0.4, 0.9, 1.6, 2.5, -1.8, -3.1 })]
29+
[InlineData("e^(3*acoth(a*x))*sqrt(1-1/(a*x))/x", new[] { 0.4, 0.9, 1.6, 2.5, -1.8, -3.1 })]
30+
// Between and beyond the roots of x^2 = a and x^2 = 2 a, 1.14 and 1.61 at a = 1.3.
31+
[InlineData("1/((a - x^2)*(2*a - x^2))", new[] { 0.4, 1.4, 2.5, -0.7, -2.2 })]
32+
[InlineData("x^4/((a - x^2)*(2*a - x^2))", new[] { 0.4, 1.4, 2.5, -0.7, -2.2 })]
33+
public void HasAValueWhereTheIntegrandIsReal(string integrand, double[] points)
34+
{
35+
var integral = integrand.ToEntity().Integrate("x");
36+
var text = integral.Stringize();
37+
Assert.DoesNotContain("integral(", text);
38+
Entity Pinned(Entity e) => e.Substitute("a", 0.7).Substitute("c", 1.07);
39+
Entity PinnedForTheQuartic(Entity e) => e.Substitute("a", 1.3);
40+
Func<Entity, Entity> pin = integrand.Contains("acoth") ? Pinned : PinnedForTheQuartic;
41+
var antiderivative = pin(integral.Substitute("C", 0));
42+
var derivative = antiderivative.Differentiate("x");
43+
var original = pin(integrand.ToEntity());
44+
foreach (var at in points)
45+
{
46+
var want = original.Substitute("x", at).EvalNumerical();
47+
Assert.True(Math.Abs((double)want.ImaginaryPart) < 1e-9, $"{integrand} is {want} at x = {at}, which is not a point where it is real");
48+
Assert.False(antiderivative.Substitute("x", at).EvalNumerical().IsNaN, $"the antiderivative of {integrand} has no value at x = {at}, where the integrand is {want}");
49+
var got = derivative.Substitute("x", at).EvalNumerical();
50+
Assert.False(got.IsNaN, $"d/dx of the antiderivative of {integrand} has no value at x = {at}, where the integrand is {want}");
51+
Assert.True(Math.Abs((double)(got - want).RealPart) + Math.Abs((double)(got - want).ImaginaryPart)
52+
< 1e-9 * Math.Max(1, Math.Abs((double)want.RealPart)),
53+
$"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}");
54+
}
55+
}
56+
57+
[Theory]
58+
// Rubi 1.1.1.2:1343 and 1359: a biquadratic with the imaginary unit in it and a positive discriminant,
59+
// whose roots are not real, keeps its arctangents; split by their signs it had no arm that held.
60+
[InlineData("1/((a - i*a*x)^(3/4)*(a + i*a*x)^(1/4))")]
61+
[InlineData("(a - i*a*x)^(3/4)/(a + i*a*x)^(7/4)")]
62+
public void KeepsItsArctangentsOverComplexCoefficients(string integrand)
63+
{
64+
var integral = integrand.ToEntity().Integrate("x");
65+
Assert.DoesNotContain("integral(", integral.Stringize());
66+
var derivative = integral.Substitute("C", 0).Substitute("a", 1.3).Differentiate("x");
67+
var original = integrand.ToEntity().Substitute("a", 1.3);
68+
foreach (var at in new[] { -2.2, -0.9, 0.4, 1.4 })
69+
{
70+
var want = original.Substitute("x", at).EvalNumerical();
71+
var got = derivative.Substitute("x", at).EvalNumerical();
72+
Assert.False(got.IsNaN, $"d/dx of the antiderivative of {integrand} has no value at x = {at}, where the integrand is {want}");
73+
Assert.True(Math.Abs((double)(got - want).RealPart) + Math.Abs((double)(got - want).ImaginaryPart)
74+
< 1e-9 * Math.Max(1, Math.Abs((double)want.RealPart) + Math.Abs((double)want.ImaginaryPart)),
75+
$"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}");
76+
}
77+
}
78+
79+
[Theory]
80+
[InlineData("1/((a - x^2)*(2*a - x^2))")]
81+
[InlineData("x^4/((a - x^2)*(2*a - x^2))")]
82+
public void IsRealPastBothRoots(string integrand)
83+
{
84+
// Past both roots each logarithm's argument is positive, and the antiderivative is real;
85+
// as arctangents of imaginary arguments it was on their cuts there.
86+
var antiderivative = integrand.ToEntity().Integrate("x").Substitute("C", 0).Substitute("a", 1.3);
87+
foreach (var at in new[] { 2.5, -2.2 })
88+
{
89+
var value = antiderivative.Substitute("x", at).EvalNumerical();
90+
Assert.False(value.IsNaN, $"the antiderivative of {integrand} has no value at x = {at}");
91+
Assert.True(Math.Abs((double)value.ImaginaryPart) < 1e-9, $"the antiderivative of {integrand} is {value} at x = {at}");
92+
}
93+
}
94+
}
95+
}

0 commit comments

Comments
 (0)