diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md
index 7cf40b66e..08fb6e2bb 100644
--- a/BREAKING-CHANGES.md
+++ b/BREAKING-CHANGES.md
@@ -2278,6 +2278,23 @@ elementary integrand whose antiderivative is reached through one is answered whe
|---|---|---|
| `"(a+b*ln(c*x^n))/(x^2*(d+e*ln(f*x^m)))".Integrate("x")`, Rubi's 3.1.5 row 213 | `integral(...)` | an antiderivative in `Ei`, provided `f > 0` and `e^d f^e > 0` |
+### A special function beside its derivative is a substitution
+
+`e^(c - b^2 x^2) erf(b x)^n`, `Ei(b x) e^(b x)/x` and `Si(b x) sin(b x)/x` are each a power of a
+special function beside its derivative, and `u = erf(b x)`, `u = Ei(b x)` or `u = Si(b x)` writes
+them as a power of `u`. None of the nine special functions was a candidate for that substitution.
+Each is one now, as the logarithm is, wherever its derivative can be the differential
+([#1501](https://github.com/asc-community/AngouriMath/issues/1501)). An integrand holding one of
+these functions had no reading in 2.5.0, which the entries for the functions themselves record, and
+no integrand without one changes: Rubi's independent suites and a sample of its families 1 to 7,
+2726 problems, are answered alone as they were.
+
+| Input | Was (2.5.0) | Now |
+|---|---|---|
+| `"e^(c - b^2*x^2)*erf(b*x)".Integrate("x")` | `UnrecognizedFunctionParseException`: there is no function `erf` | `e ^ c * pi ^ (1/2) * erf(b * x) ^ 2 / 2 / (2 * b) + C` |
+| `"e^(b*x)*Ei(b*x)/x".Integrate("x")` | `Ei * e ^ (b * x) + C`, with `Ei` a variable | `Ei(b * x) ^ 2 / 2 + C` |
+| `"Si(b*x)*sin(b*x)/x".Integrate("x")` | `Si * -cos(b * x) + C`, with `Si` a variable | `Si(b * x) ^ 2 / 2 + C` |
+
### An inverse trigonometric function below the bar is integrated to the sine and cosine integrals
`1/arcsin(x)` was left unintegrated. Under the substitution that undoes the inverse function,
diff --git a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs
index 4b57434f7..b230c5e49 100644
--- a/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs
+++ b/Sources/AngouriMath/Functions/Continuous/Integration/IndefiniteIntegralSolver.cs
@@ -1610,8 +1610,7 @@ private static bool IsDifferentiatedBeforeAPolynomial(Entity factor)
Logf or Entity.Arcsinf or Entity.Arccosf
or Entity.Arctanf or Entity.Arccotanf
or Entity.Arcsecantf or Entity.Arccosecantf => true,
- Entity.Erff or Entity.Erfcf or Entity.Erfif or Entity.Eif or Entity.Lif
- or Entity.Sif or Entity.Cif or Entity.Shif or Entity.Chif => true,
+ _ when IsASpecialFunction(factor) => true,
// And a whole power of one, which is the same function for this purpose:
// differentiating ln(x)^2 gives 2ln(x)/x, whose x cancels against the integrated
// polynomial exactly as ln(x)'s does, leaving x*ln(x) -- one step simpler, and
@@ -19624,6 +19623,33 @@ private static bool HoldsAnEvenRootBesides(Entity expr, Entity radicand, Entity.
integrandInU = SimplifiedWithoutTheImaginaryUnit(cancelled, expr);
}
}
+ // Powers of one base on the two sides of the bar: `e^(c - b^2 x^2)` over the
+ // `e^(-b^2 x^2)` of du/dx is `e^c`, which the one-level simplification leaves as
+ // written, and `e^(c - b^2 x^2) erf(b x)` was refused under `u = erf(b x)` for
+ // the x in it, where `e^(-b^2 x^2) erf(b x)` was answered.
+ // https://github.com/asc-community/AngouriMath/issues/1501
+ // The divisor's factors are spread first, each to the power -1: gathered as it
+ // stands, the whole divisor `2 e^(-(b x)^2) b` is one factor below the bar, with
+ // a base of its own.
+ // For a special function only, which is where the constant in the exponent comes
+ // from: done for every candidate that left an x, it gathered the exponentials of
+ // `1/((c + d x)^3 (a + a tanh(e + f x)))` and simplified them for every one, and
+ // a decline in a third of a second became a timeout.
+ if (integrandInU.ContainsNode(x) && IsASpecialFunction(u))
+ {
+ var (above, below) = Functions.SingleQuotient.Of(Functions.SingleQuotient.Combine(integrandInU));
+ var spread = above;
+ foreach (var factor in Mulf.LinearChildren(below))
+ spread *= MathS.Pow(factor, -1);
+ // Simplified a level, not only inner-simplified, since the gathered exponent
+ // `c - (b x)^2 + (b x)^2` is left standing by the inner simplification -- and
+ // only where a base was gathered at all. Simplified for every candidate that
+ // left an x, it took `x cosh(a + b x) Shi(a + b x)` from a second to a timeout.
+ var gatheredRaw = Patterns.GatherPowersOfOneBase(spread);
+ if (!ReferenceEquals(gatheredRaw, spread)
+ && SimplifiedWithoutTheImaginaryUnit(gatheredRaw, expr) is var gathered && !gathered.ContainsNode(x))
+ integrandInU = gathered;
+ }
// A polynomial in x left over under a candidate that is itself a polynomial
// is written in the candidate where it is one in it: `(1 - x)^2` is
// `1 - 2x + x^2`, and that is `u` for Apostol's `(1 - 2x + x^2)^(1/5)/(1 - x)`,
@@ -20124,6 +20150,34 @@ Entity Term(Entity? sum, Entity coefficient, int power)
/// 2. f(x^n) * x^(n-1) -> u = x^n
/// 3. f(g(x)) * g'(x) -> u = g(x)
///
+ ///
+ /// The special functions of
+ /// #1501: the error
+ /// functions and the exponential, logarithmic, sine, cosine and hyperbolic integrals.
+ ///
+ private static bool IsASpecialFunction(Entity node)
+ => node is Entity.Erff or Entity.Erfcf or Entity.Erfif or Entity.Eif or Entity.Lif
+ or Entity.Sif or Entity.Cif or Entity.Shif or Entity.Chif;
+
+ ///
+ /// Whether has what the derivative of the special function
+ /// is made of: an exponential of a quadratic in
+ /// for the error functions, something in below the
+ /// bar for the exponential, trigonometric and hyperbolic integrals, whose derivatives
+ /// divide by their argument, and a logarithm below the bar for the logarithmic integral.
+ ///
+ private static bool TheDifferentialCanBeThere(Entity special, Entity expr, Entity.Variable x)
+ {
+ var below = Functions.SingleQuotient.Of(Functions.SingleQuotient.Combine(expr)).Denominator;
+ return special switch
+ {
+ Entity.Erff or Entity.Erfcf or Entity.Erfif => expr.Nodes.Any(node => node is Powf(var @base, var exponent)
+ && !@base.ContainsNode(x) && TreeAnalyzer.TryGetPolyQuadratic(exponent, x, out var square, out _, out _) && !TreeAnalyzer.IsZero(square)),
+ Entity.Lif => below.Nodes.Any(node => node is Logf(_, var antilogarithm) && antilogarithm.ContainsNode(x)),
+ _ => below.ContainsNode(x),
+ };
+ }
+
private static IEnumerable FindSubstitutionCandidates(Entity expr, Entity.Variable x)
{
var candidates = new List();
@@ -20234,6 +20288,18 @@ bool APowerCanBeExact(int k)
candidates.Add(node); // Logarithm itself (for cases like 1/(x*ln(x)))
if (antilog != x && antilog.ContainsNode(x)) candidates.Add(antilog); // Also add the argument if it's not just x
break;
+ // A special function itself, as the logarithm is: each has an elementary
+ // derivative, so beside a power of it that derivative is the differential.
+ // `e^(c - b^2 x^2) erf(b x)^n` is `sqrt(pi) e^c/(2b) u^n` under `u = erf(b x)`.
+ // https://github.com/asc-community/AngouriMath/issues/1501
+ // Only where that derivative can be there: an exponential of a quadratic for
+ // the error functions, a divisor in x for the u below the bar of e^u/u, sin(u)/u
+ // and the others, a logarithm below the bar for li. Offered beside anything,
+ // `x cosh(a + b x) Shi(a + b x)` went from a second to ten, simplifying quotients
+ // of exponentials that were never going to lose their x.
+ case var special when IsASpecialFunction(special) && TheDifferentialCanBeThere(special, expr, x):
+ candidates.Add(node);
+ break;
case Sumf(var aug, var add) when !rational && !large && node.Complexity <= LargestSumOffered:
if (aug.ContainsNode(x) || add.ContainsNode(x)) candidates.Add(node); // Linear expressions ax + b
break;
diff --git a/Sources/Tests/UnitTests/Calculus/SpecialFunctionSubstitutionTest.cs b/Sources/Tests/UnitTests/Calculus/SpecialFunctionSubstitutionTest.cs
new file mode 100644
index 000000000..3d7f80d67
--- /dev/null
+++ b/Sources/Tests/UnitTests/Calculus/SpecialFunctionSubstitutionTest.cs
@@ -0,0 +1,66 @@
+//
+// 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 System.Linq;
+using AngouriMath.Extensions;
+using Xunit;
+
+namespace AngouriMath.Tests.Calculus
+{
+ ///
+ /// A special function as the substitution: beside a power of it, its derivative is the
+ /// differential, so e^(c - b^2 x^2) erf(b x)^n is a power of u = erf(b x), and
+ /// Si(b x) sin(b x)/x is u du under u = Si(b x). The rows are Rubi's, from
+ /// its files 8.1, 8.3, 8.4 and 8.5.
+ /// https://github.com/asc-community/AngouriMath/issues/1501
+ ///
+ [Trait("Area", "Calculus")]
+ public sealed class SpecialFunctionSubstitutionTest
+ {
+ /// Off 0, where the reciprocals are undefined.
+ private static readonly double[] Points = { -1.7, -0.6, 0.35, 0.9, 1.45 };
+
+ private static readonly (string, string)[] Parameters = { ("b", "13/10"), ("c", "7/10") };
+
+ ///
+ /// Integrates, pins the parameters, and compares the derivative of the answer with the
+ /// integrand at . The parameters are pinned after integrating, so the
+ /// rule is asked the symbolic question.
+ ///
+ private static void DifferentiatesBack(string integrand)
+ {
+ var integral = integrand.ToEntity().Integrate("x");
+ Assert.DoesNotContain("integral(", integral.Stringize());
+ Entity Pinned(Entity e) => Parameters.Aggregate(e, (current, pin) => current.Substitute(pin.Item1, pin.Item2.ToEntity()));
+ var derivative = Pinned(integral.Substitute("C", 0)).Differentiate("x");
+ var original = Pinned(integrand.ToEntity());
+ foreach (var at in Points)
+ {
+ var got = derivative.Substitute("x", at).EvalNumerical();
+ var want = original.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) + Math.Abs((double)want.ImaginaryPart));
+ Assert.True(difference / scale < 1e-9,
+ $"d/dx of the antiderivative of {integrand} is {got} at x = {at}, where the integrand is {want}");
+ }
+ }
+
+ [Theory]
+ [InlineData("e^(c - b^2*x^2)*erf(b*x)")]
+ [InlineData("e^(c - b^2*x^2)*erf(b*x)^2")]
+ [InlineData("e^(c - b^2*x^2)/erf(b*x)")]
+ [InlineData("e^(c - b^2*x^2)/erf(b*x)^2")]
+ [InlineData("e^(c - b^2*x^2)*erfc(b*x)^3")]
+ [InlineData("e^(c + b^2*x^2)*erfi(b*x)")]
+ [InlineData("e^(b*x)*Ei(b*x)/x")]
+ [InlineData("Si(b*x)*sin(b*x)/x")]
+ [InlineData("cos(b*x)*Ci(b*x)/x")]
+ public void ThePowerBesideTheDerivativeIsASubstitution(string integrand)
+ => DifferentiatesBack(integrand);
+ }
+}