From f6c5ca55de8ffec4a77ef1cb0ed94b4a887bf22d Mon Sep 17 00:00:00 2001 From: Rafael Vuijk Date: Wed, 30 Sep 2026 07:48:08 +0000 Subject: [PATCH] The factorial, gamma and the totient compile under Compile Compile had no form for x!, gamma or phi, and Compile's own had defects. Both compilers now call one routine in SpecialFunctions. The factorial is exact at a whole number up to 170, rounded once to a double, and infinite past it; NaN at a negative whole number, one of the gamma function's poles, as in the interpreter; and elsewhere the Lanczos approximation Compile always used (g = 7, nine terms), taken in logarithms, the reflection formula included. Compile's overflowed its own power before the value did, so from just past 141!, about 2e243, it compiled to NaN, up to 170!. At a whole number it was a rounding off the integer, 23.999999999999996 at 4, which a compilation over integers truncates to 23. Off the real line the reflection's sin(pi z) overflowed where the value is tiny; past |Im(pi z)| = 20 it is a logarithm written out now. gamma(x) is (x - 1)!. The totient is the totient at a whole number, 0 at one that is not positive, and NaN at any other number, as the interpreter's is. Compile truncated a fraction to the whole number below it, so phi(2.5) was phi(2). Under Compile, long, int and BigInteger take it exactly. Against mpmath the compiled factorial is within 1.1e-13 of the value on a grid of [-14, 14]^2, 2.4e-13 up the imaginary axis to 10^4.5 and 3.2e-13 on circles out to radius 300. It is further out only where the conditioning forces it: next to a pole, and in the subnormal range past -177. SpecialFunctionCompilationTest gains the factorial and gamma against the interpreter, through both compilers and over doubles; whole numbers exact, including under Compile; 171 infinite and the poles NaN; and the totient at whole numbers over every integer type, NaN elsewhere. UserInvalidExceptions keeps binomial, gcd, lcm and prime as the nodes with no compiled form. Closes #1607. Co-Authored-By: Claude Opus 5.5 --- BREAKING-CHANGES.md | 21 ++++ .../Compilation/IntoFE/FastExpression.cs | 28 +---- .../IntoLinq/CompilationProtocol.cs | 2 + .../Compilation/IntoLinq/MathAllMethods.cs | 23 +++- .../Compilation/IntoLinq/MathAllMethods.tt | 17 ++- .../AngouriMath/Numerics/SpecialFunctions.cs | 114 ++++++++++++++++++ .../SpecialFunctionCompilationTest.cs | 93 +++++++++++++- .../UnitTests/Core/UserInvalidExceptions.cs | 11 +- 8 files changed, 271 insertions(+), 38 deletions(-) diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index b71dcb1c9..366f5f456 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -2354,6 +2354,27 @@ two of them, the points are samples of a continuum and not the set | `"floor(x)".Compile("x")`, and `ceil`, `round`, `max`, `min` | `UncompilableNodeException` | compiled | | `"max(x, 1)".ToEntity().Differentiate("x")` | `derivative(max(x, 1), x)` | `(1 + sgn(x - 1)) / 2 provided not x - 1 = 0` | +### The compiled factorial is exact at whole numbers and finite to `170!`, and the compiled totient has no value at a fraction + +Both compilers now evaluate the factorial, the gamma function and Euler's totient the same way, +and `Compile` has a form for them where it had none +([#1607](https://github.com/asc-community/AngouriMath/issues/1607)). The compiled factorial was +Lanczos's approximation everywhere, and its power overflowed before the value did: from just past +`141!`, about `2e243`, it compiled to NaN, up to `170!`. At a whole number it was off from the +integer by a rounding, which a compilation over integers would truncate. It is now the exact +factorial rounded once, up to `170!`; infinite past it; NaN at the negative whole numbers, the +poles; and the same approximation elsewhere, taken in logarithms so that it overflows only where +the value does. The compiled totient truncated a fraction to the whole number below it. There it +now has no value, as the interpreter's has none. + +| Input | Was (2.5.0) | Now | +|---|---|---| +| `"x!".Compile("x").Call(4)` | `23.999999999999996` | `24` | +| `"x!".Compile("x").Call(0)` | `0.9999999999999998` | `1` | +| `"x!".Compile("x").Call(142)`, and on to `170` | NaN | `2.695364137888163E+245` | +| `"phi(x)".Compile("x").Call(2.5)` | `1` | NaN | +| `"x!".ToEntity().Compile("x")`, and `gamma`, `phi` | `UncompilableNodeException` | compiled | + ### An equation the solver cannot invert is left unsolved, not answered with no roots `x! = 6` was answered `{ }`, a claim that it has no roots, and it has 3. The solver isolates `x` diff --git a/Sources/AngouriMath/Functions/Compilation/IntoFE/FastExpression.cs b/Sources/AngouriMath/Functions/Compilation/IntoFE/FastExpression.cs index 3d08fa302..15c4b1ec8 100644 --- a/Sources/AngouriMath/Functions/Compilation/IntoFE/FastExpression.cs +++ b/Sources/AngouriMath/Functions/Compilation/IntoFE/FastExpression.cs @@ -258,24 +258,7 @@ public System.Numerics.Complex Substitute(params System.Numerics.Complex[] value stack.Push(AngouriMath.Numerics.SpecialFunctions.Chi(stack.Pop())); break; case InstructionType.CALL_FACTORIAL: - // https://stackoverflow.com/a/15454784/5429648 - const int g = 7; - static System.Numerics.Complex Gamma(System.Numerics.Complex z) - { - if (z.Real < 0.5) return System.Math.PI / (System.Numerics.Complex.Sin(System.Math.PI * z) * Gamma(1 - z)); - else - { - z -= 1; - - System.Numerics.Complex x = gammaCoeffs[0]; - for (var i = 1; i < g + 2; i++) - x += gammaCoeffs[i] / (z + i); - - var t = z + g + 0.5; - return System.Math.Sqrt(2 * System.Math.PI) * System.Numerics.Complex.Pow(t, z + 0.5) * System.Numerics.Complex.Exp(-t) * x; - } - } - stack.Push(Gamma(stack.Pop() + 1)); + stack.Push(AngouriMath.Numerics.SpecialFunctions.Factorial(stack.Pop())); break; case InstructionType.CALL_SIGNUM: stack.Push(stack.Pop().Signum()); @@ -284,8 +267,7 @@ static System.Numerics.Complex Gamma(System.Numerics.Complex z) stack.Push(stack.Pop().Abs()); break; case InstructionType.CALL_PHI: - var n = (long) stack.Pop().Real; - stack.Push(n.Phi()); + stack.Push(AngouriMath.Numerics.SpecialFunctions.Phi(stack.Pop())); break; } if (stack.Count != 1) @@ -293,12 +275,6 @@ static System.Numerics.Complex Gamma(System.Numerics.Complex z) return stack.Pop(); } - [ConstantField] static readonly double[] gammaCoeffs = { - 0.99999999999980993, 676.5203681218851, -1259.1392167224028, - 771.32342877765313, -176.61502916214059, 12.507343278686905, - -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7 - }; - /// Might be useful for debug if a function works too slowly public override string ToString() => string.Join(" \n| ", instructions); } diff --git a/Sources/AngouriMath/Functions/Compilation/IntoLinq/CompilationProtocol.cs b/Sources/AngouriMath/Functions/Compilation/IntoLinq/CompilationProtocol.cs index 6db3fa3a4..adcded759 100644 --- a/Sources/AngouriMath/Functions/Compilation/IntoLinq/CompilationProtocol.cs +++ b/Sources/AngouriMath/Functions/Compilation/IntoLinq/CompilationProtocol.cs @@ -149,6 +149,8 @@ Sif when ShouldBeAtLeastDouble(e) is var newE => Expression.Call(GetDef("Si", Cif when ShouldBeAtLeastDouble(e) is var newE => Expression.Call(GetDef("Ci", 1, newE.Type), newE), Shif when ShouldBeAtLeastDouble(e) is var newE => Expression.Call(GetDef("Shi", 1, newE.Type), newE), Chif when ShouldBeAtLeastDouble(e) is var newE => Expression.Call(GetDef("Chi", 1, newE.Type), newE), + Factorialf when ShouldBeAtLeastDouble(e) is var newE => Expression.Call(GetDef("Factorial", 1, newE.Type), newE), + Phif => Expression.Call(GetDef("Phi", 1, e.Type), e), Notf => Expression.Not(e), diff --git a/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.cs b/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.cs index 46b76d12a..dea2e704e 100644 --- a/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.cs +++ b/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.cs @@ -288,8 +288,17 @@ public static System.Numerics.Complex Round(System.Numerics.Complex a) public static int Min(int a, int b) => Math.Min(a, b); public static BigInteger Min(BigInteger a, BigInteger b) => BigInteger.Min(a, b); - // The special functions, in double precision as the kernels define them. - // https://github.com/asc-community/AngouriMath/issues/1607 + // Euler's totient: exact over the integers, and over Complex and double a value at a whole + // number and NaN at any other. https://github.com/asc-community/AngouriMath/issues/1607 + public static System.Numerics.Complex Phi(System.Numerics.Complex a) => AngouriMath.Numerics.SpecialFunctions.Phi(a); + public static double Phi(double a) => AngouriMath.Numerics.SpecialFunctions.Phi(a); + public static float Phi(float a) => (float)AngouriMath.Numerics.SpecialFunctions.Phi(a); + public static long Phi(long a) => a.Phi(); + public static int Phi(int a) => (int)((long)a).Phi(); + public static BigInteger Phi(BigInteger a) => PeterO.Numbers.EInteger.FromBytes(a.ToByteArray(), littleEndian: true).Phi().ToBigInteger(); + + // The special functions, in double precision as the kernels define them, and the factorial, + // the gamma function one along. https://github.com/asc-community/AngouriMath/issues/1607 public static System.Numerics.Complex Erf(System.Numerics.Complex a) => AngouriMath.Numerics.SpecialFunctions.Erf(a); public static double Erf(double a) => AngouriMath.Numerics.SpecialFunctions.Erf(a); public static System.Numerics.Complex Erfc(System.Numerics.Complex a) => AngouriMath.Numerics.SpecialFunctions.Erfc(a); @@ -308,6 +317,8 @@ public static System.Numerics.Complex Round(System.Numerics.Complex a) public static double Shi(double a) => AngouriMath.Numerics.SpecialFunctions.Shi(a); public static System.Numerics.Complex Chi(System.Numerics.Complex a) => AngouriMath.Numerics.SpecialFunctions.Chi(a); public static double Chi(double a) => AngouriMath.Numerics.SpecialFunctions.Chi(a); + public static System.Numerics.Complex Factorial(System.Numerics.Complex a) => AngouriMath.Numerics.SpecialFunctions.Factorial(a); + public static double Factorial(double a) => AngouriMath.Numerics.SpecialFunctions.Factorial(a); /// /// True where System.Numerics.Complex.IsNaN is, spelled out because that overload is @@ -459,6 +470,12 @@ public static bool IsNaN(System.Numerics.Complex a) { ("Round", 1, typeof(long)), Def(a => Round(a)) }, { ("Round", 1, typeof(int)), Def(a => Round(a)) }, { ("Round", 1, typeof(BigInteger)), Def(a => Round(a)) }, + { ("Phi", 1, typeof(System.Numerics.Complex)), Def(a => Phi(a)) }, + { ("Phi", 1, typeof(double)), Def(a => Phi(a)) }, + { ("Phi", 1, typeof(float)), Def(a => Phi(a)) }, + { ("Phi", 1, typeof(long)), Def(a => Phi(a)) }, + { ("Phi", 1, typeof(int)), Def(a => Phi(a)) }, + { ("Phi", 1, typeof(BigInteger)), Def(a => Phi(a)) }, { ("Erf", 1, typeof(System.Numerics.Complex)), Def(a => Erf(a)) }, { ("Erf", 1, typeof(double)), Def(a => Erf(a)) }, { ("Erfc", 1, typeof(System.Numerics.Complex)), Def(a => Erfc(a)) }, @@ -477,6 +494,8 @@ public static bool IsNaN(System.Numerics.Complex a) { ("Shi", 1, typeof(double)), Def(a => Shi(a)) }, { ("Chi", 1, typeof(System.Numerics.Complex)), Def(a => Chi(a)) }, { ("Chi", 1, typeof(double)), Def(a => Chi(a)) }, + { ("Factorial", 1, typeof(System.Numerics.Complex)), Def(a => Factorial(a)) }, + { ("Factorial", 1, typeof(double)), Def(a => Factorial(a)) }, }; // An expression tree is how C# spells "the MethodInfo of this method" without a string: diff --git a/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.tt b/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.tt index 9fda9de21..809048155 100644 --- a/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.tt +++ b/Sources/AngouriMath/Functions/Compilation/IntoLinq/MathAllMethods.tt @@ -7,11 +7,11 @@ // to one and forgotten in the other. var unary = new[]{ "Sin", "Cos", "Tan", "Asin", "Acos", "Atan", "Cot", "Sec", "Csc", "Acot", "Asec", "Acsc", - "Abs", "Sgn", "Floor", "Ceil", "Round" }; + "Abs", "Sgn", "Floor", "Ceil", "Round", "Phi" }; var binary = new[]{ "Log", "Pow", "Max", "Min" }; // Defined for Complex and double alone: an integer or a float is widened to double first, as // for the trigonometry. - var special = new[]{ "Erf", "Erfc", "Erfi", "Ei", "Li", "Si", "Ci", "Shi", "Chi" }; + var special = new[]{ "Erf", "Erfc", "Erfi", "Ei", "Li", "Si", "Ci", "Shi", "Chi", "Factorial" }; var types = new[]{ "System.Numerics.Complex", "double", "float", "long", "int", "BigInteger" }; #> // @@ -191,8 +191,17 @@ namespace AngouriMath.Core.Compilation.IntoLinq public static int Min(int a, int b) => Math.Min(a, b); public static BigInteger Min(BigInteger a, BigInteger b) => BigInteger.Min(a, b); - // The special functions, in double precision as the kernels define them. - // https://github.com/asc-community/AngouriMath/issues/1607 + // Euler's totient: exact over the integers, and over Complex and double a value at a whole + // number and NaN at any other. https://github.com/asc-community/AngouriMath/issues/1607 + public static System.Numerics.Complex Phi(System.Numerics.Complex a) => AngouriMath.Numerics.SpecialFunctions.Phi(a); + public static double Phi(double a) => AngouriMath.Numerics.SpecialFunctions.Phi(a); + public static float Phi(float a) => (float)AngouriMath.Numerics.SpecialFunctions.Phi(a); + public static long Phi(long a) => a.Phi(); + public static int Phi(int a) => (int)((long)a).Phi(); + public static BigInteger Phi(BigInteger a) => PeterO.Numbers.EInteger.FromBytes(a.ToByteArray(), littleEndian: true).Phi().ToBigInteger(); + + // The special functions, in double precision as the kernels define them, and the factorial, + // the gamma function one along. https://github.com/asc-community/AngouriMath/issues/1607 <# foreach (var name in special) { #> public static System.Numerics.Complex <#= name #>(System.Numerics.Complex a) => AngouriMath.Numerics.SpecialFunctions.<#= name #>(a); public static double <#= name #>(double a) => AngouriMath.Numerics.SpecialFunctions.<#= name #>(a); diff --git a/Sources/AngouriMath/Numerics/SpecialFunctions.cs b/Sources/AngouriMath/Numerics/SpecialFunctions.cs index d19a454d1..966800902 100644 --- a/Sources/AngouriMath/Numerics/SpecialFunctions.cs +++ b/Sources/AngouriMath/Numerics/SpecialFunctions.cs @@ -380,6 +380,120 @@ private static NumericsComplex HalfOfTheExponentialSeries(NumericsComplex z, boo return sum; } + /// + /// x!, the gamma function one along. At a whole number up to 170 it is the exact + /// factorial rounded once, and past 170, where a double cannot hold it, it is infinite. At a + /// negative whole number, one of the gamma function's poles, it is NaN, as the interpreter's + /// is. Elsewhere it is Lanczos's approximation with g = 7 and nine terms, the one the + /// compiler has always used. Its t^(z + 1/2) e^(-t) is taken in logarithms, so that + /// it overflows only where the value does: the power alone passed the largest double just past + /// 141, where the factorial is about 2e243, and the compiled factorial was NaN from + /// there to 170. + /// https://github.com/asc-community/AngouriMath/issues/1607 + /// + public static NumericsComplex Factorial(NumericsComplex z) + { + if (z.Imaginary == 0) + return Factorial(z.Real); + if (IsNaN(z) || IsInfinite(z)) + return NaN; + return Gamma(new NumericsComplex(z.Real + 1, z.Imaginary)); + } + + public static double Factorial(double x) + { + if (double.IsNaN(x) || double.IsNegativeInfinity(x)) + return double.NaN; + if (double.IsPositiveInfinity(x)) + return double.PositiveInfinity; + if (x == Math.Floor(x)) + return x < 0 ? double.NaN : x <= 170 ? factorials[(int)x] : double.PositiveInfinity; + return Gamma(x + 1); + } + + /// 0! to 170!, each the exact integer rounded once to a double. + [AngouriMath.Core.ConstantField] private static readonly double[] factorials = ExactFactorials(); + + private static double[] ExactFactorials() + { + var table = new double[171]; + var factorial = PeterO.Numbers.EInteger.One; + for (var n = 0; n <= 170; n++) + { + if (n > 0) + factorial *= n; + table[n] = PeterO.Numbers.EFloat.FromEInteger(factorial).ToDouble(); + } + return table; + } + + [AngouriMath.Core.ConstantField] private static readonly double[] lanczos = { + 0.99999999999980993, 676.5203681218851, -1259.1392167224028, + 771.32342877765313, -176.61502916214059, 12.507343278686905, + -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7 }; + + private static double Gamma(double x) + { + // The reflection formula below 1/2. + if (x < 0.5) + return Math.PI / (Math.Sin(Math.PI * x) * Gamma(1 - x)); + x -= 1; + var sum = lanczos[0]; + for (var i = 1; i < 9; i++) + sum += lanczos[i] / (x + i); + var t = x + 7.5; + return Math.Exp((x + 0.5) * Math.Log(t) - t + Math.Log(Math.Sqrt(2 * Math.PI) * sum)); + } + + /// + /// The gamma function off the real line, as the exponential of its logarithm, so that it + /// overflows and underflows only where the value does. The reflection formula + /// pi/(sin(pi z) gamma(1 - z)) is taken in logarithms too: far from the real line + /// sin(pi z) overflows where the value it divides is tiny, and it made a NaN. + /// + private static NumericsComplex Gamma(NumericsComplex z) + => z.Real < 0.5 + ? ExpTimes(Math.Log(Math.PI) - LogSine(Math.PI * z) - LogGamma(1 - z), 1) + : ExpTimes(LogGamma(z), 1); + + /// A logarithm of the gamma function for Re z >= 1/2, by Lanczos's approximation. + private static NumericsComplex LogGamma(NumericsComplex z) + { + z -= 1; + NumericsComplex sum = lanczos[0]; + for (var i = 1; i < 9; i++) + sum += lanczos[i] / (z + i); + var t = z + 7.5; + return (z + 0.5) * NumericsComplex.Log(t) - t + NumericsComplex.Log(Math.Sqrt(2 * Math.PI) * sum); + } + + /// + /// A logarithm of sin w. Past |Im w| = 20 one of e^(i w) and + /// e^(-i w) in sin w = (e^(i w) - e^(-i w))/(2 i) is below 1e-17 of the + /// other. Its logarithm is then written out, where sin w itself would overflow. + /// + private static NumericsComplex LogSine(NumericsComplex w) + { + if (Math.Abs(w.Imaginary) < 20) + return NumericsComplex.Log(NumericsComplex.Sin(w)); + // (i/2) e^(-i w) above the real line, and -(i/2) e^(i w) below it. + return w.Imaginary > 0 + ? new NumericsComplex(w.Imaginary - Math.Log(2), Math.PI / 2 - w.Real) + : new NumericsComplex(-w.Imaginary - Math.Log(2), w.Real - Math.PI / 2); + } + + /// + /// Euler's totient at a whole number, 0 at one that is not positive, as the interpreter has + /// it, and NaN at any other number, where the interpreter's is NaN too. The compiler + /// truncated a fraction to the whole number below it, so phi(2.5) was phi(2). + /// Past long's range it is NaN: the interpreter's totient cannot go there either. + /// https://github.com/asc-community/AngouriMath/issues/1607 + /// + public static NumericsComplex Phi(NumericsComplex z) => z.Imaginary == 0 ? Phi(z.Real) : NaN; + + public static double Phi(double x) + => x == Math.Floor(x) && x >= long.MinValue && x < long.MaxValue ? ((long)x).Phi() : double.NaN; + /// /// Half of , part by part: dividing a complex number by 2 multiplies it /// by the complex number 2, and an infinite part times the other's zero makes a NaN. diff --git a/Sources/Tests/UnitTests/Convenience/SpecialFunctionCompilationTest.cs b/Sources/Tests/UnitTests/Convenience/SpecialFunctionCompilationTest.cs index 0732104c1..eeb9c608b 100644 --- a/Sources/Tests/UnitTests/Convenience/SpecialFunctionCompilationTest.cs +++ b/Sources/Tests/UnitTests/Convenience/SpecialFunctionCompilationTest.cs @@ -14,7 +14,7 @@ namespace AngouriMath.Tests.Convenience /// /// The special functions compile, through Compile and Compile<TIn, TOut> /// alike, in double precision and in agreement with the interpreter: its branch cuts, and its - /// values on the axes. + /// values on the axes. So do the factorial, the gamma function and Euler's totient. /// #1607 /// [Trait("Area", "Convenience")] @@ -114,5 +114,96 @@ public void AnInfiniteValueKeepsItsExactPart(string expression, double re, doubl [InlineData("Chi(x)", -1.7)] public void OverDoublesAValueOffTheRealLineIsNaN(string expression, double at) => Assert.True(double.IsNaN(expression.ToEntity().Compile("x")(at))); + + /// + /// The factorial and the gamma function, through both compilers and over doubles, against + /// the interpreter. The compiled factorial was NaN from just past 141 to 170, where its power + /// overflowed before the value did. + /// + [Theory] + [InlineData("x!", 0.5, 0)] + [InlineData("x!", -0.5, 0)] + [InlineData("x!", -2.5, 0)] + [InlineData("x!", 100.5, 0)] + [InlineData("x!", 142, 0)] + [InlineData("x!", 170, 0)] + [InlineData("x!", 2, 1)] + [InlineData("x!", -3.5, 2)] + [InlineData("gamma(x)", 5, 0)] + [InlineData("gamma(x)", 0.5, 0)] + [InlineData("gamma(x)", 1.5, -2)] + public void TheFactorialAgreesWithTheInterpreter(string expression, double re, double im) + { + var expected = Interpreted(expression, re, im); + var at = new NumericsComplex(re, im); + Close(expected, expression.ToEntity().Compile("x").Call(at)); + Close(expected, expression.ToEntity().Compile("x")(at)); + if (im == 0) + Close(expected, expression.ToEntity().Compile("x")(re)); + } + + /// + /// A whole number's factorial is the exact integer, rounded once: the compiled one was + /// 0.9999999999999998 at 0 and 23.999999999999996 at 4, which an integer-typed + /// compilation would have truncated to 23. + /// + [Theory] + [InlineData(0, 1)] + [InlineData(4, 24)] + [InlineData(20, 2432902008176640000)] + public void AWholeNumbersFactorialIsExact(long n, long factorial) + { + Assert.Equal(factorial, "x!".ToEntity().Compile("x").Call(n).Real); + Assert.Equal(factorial, "x!".ToEntity().Compile("x")(n)); + Assert.Equal(factorial, "gamma(x + 1)".ToEntity().Compile("x")(n)); + } + + /// Past 170 the factorial overflows a double; at a negative whole number, a pole, it has no value. + [Theory] + [InlineData(171, double.PositiveInfinity)] + [InlineData(-1, double.NaN)] + [InlineData(-3, double.NaN)] + public void PastTheLargestDoubleAndAtThePoles(double at, double expected) + { + Assert.Equal(expected, "x!".ToEntity().Compile("x").Call(at).Real); + Assert.Equal(expected, "x!".ToEntity().Compile("x")(at)); + } + + /// + /// Euler's totient at a whole number, and 0 at one that is not positive, as the interpreter + /// has it, through both compilers and over every integer type. + /// + [Theory] + [InlineData(12, 4)] + [InlineData(1, 1)] + [InlineData(97, 96)] + [InlineData(0, 0)] + [InlineData(-3, 0)] + public void ThePhiOfAWholeNumber(long n, long phi) + { + Assert.Equal(phi, "phi(x)".ToEntity().Compile("x").Call(n).Real); + Assert.Equal(phi, "phi(x)".ToEntity().Compile("x")(n)); + Assert.Equal(phi, "phi(x)".ToEntity().Compile("x")(n)); + Assert.Equal(phi, "phi(x)".ToEntity().Compile("x")((int)n)); + Assert.Equal(new System.Numerics.BigInteger(phi), "phi(x)".ToEntity().Compile("x")(n)); + // With the downcasting on, which reads the double as the integer it is: with it off, the + // interpreter's totient sees a decimal and has no value for it. + Assert.Equal(phi.ToString(), $"phi({n})".ToEntity().EvalNumerical().ToString()); + } + + /// + /// Anywhere else the totient has no value, and the interpreter says NaN. The compiler + /// truncated a fraction to the whole number below it, so phi(2.5) was phi(2), 1. + /// + [Theory] + [InlineData(2.5, 0)] + [InlineData(3, 1)] + public void ThePhiOfAnyOtherNumberIsNaN(double re, double im) + { + Assert.True(double.IsNaN("phi(x)".ToEntity().Compile("x").Call(new NumericsComplex(re, im)).Real)); + Assert.True(double.IsNaN("phi(x)".ToEntity().Compile("x")(new NumericsComplex(re, im)).Real)); + if (im == 0) + Assert.True(double.IsNaN("phi(x)".ToEntity().Compile("x")(re))); + } } } diff --git a/Sources/Tests/UnitTests/Core/UserInvalidExceptions.cs b/Sources/Tests/UnitTests/Core/UserInvalidExceptions.cs index 33f140962..bba878ab0 100644 --- a/Sources/Tests/UnitTests/Core/UserInvalidExceptions.cs +++ b/Sources/Tests/UnitTests/Core/UserInvalidExceptions.cs @@ -55,12 +55,13 @@ [Fact] public void UncompilableNode4() => // A node the Linq converter has no case for threw AngouriBugException, which asks the // caller to report a bug -- for a gap in coverage the library already knows about. Four // of these were nodes 2.0 added and never taught the compiler; floor, ceil and round have - // been taught since (https://github.com/asc-community/AngouriMath/issues/1603). + // been taught since (https://github.com/asc-community/AngouriMath/issues/1603), and phi, + // gamma and the factorial (https://github.com/asc-community/AngouriMath/issues/1607). [Theory] - [InlineData("phi(x)")] - [InlineData("gamma(x)")] - [InlineData("x!")] - [InlineData("(x + 1)! / x!")] + [InlineData("binomial(x, 2)")] + [InlineData("gcd(x, 6)")] + [InlineData("lcm(x, 6)")] + [InlineData("prime(x)")] public void ANodeWithNoCompiledFormSaysSoRatherThanAskingForABugReport(string expression) => Assert.Throws(() => expression.ToEntity().Compile("x"));