diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 0954fd435..dc47247c6 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -35,6 +35,8 @@ read first. | **Silent** | `limit(i, i, 0)` | unevaluated | `0` | | **Silent** | `integral(i, i)` | `-1/2 + C` | `i_1 ^ 2 / 2 + C` | | | `lambda(i, i + 1)` | `InvalidArgumentParseException` | the lambda | +| **Silent** | the symbolic determinant of a matrix, substituted where a pivot vanishes — `[[0,1,2],[3,4,5],[6,7,8]]` through `[[a,b,c],[d,e,f],[g,h,i]]` | `NaN` | `0` | +| **Silent** | `((Entity.Matrix)"[[x, 1], [2, y]]".ToEntity()).Determinant.Simplify()` | `x * y - 2 provided not x = 0` | `x * y - 2` | | **Silent** | `[1, 2] in RR`, and a matrix or a finite set against any of `BB`, `ZZ`, `QQ`, `RR`, `CC` | `True` | `False` | | **Silent** | `"sin(pi)".ToEntity().Differentiate(MathS.pi)`, and every derivative over `pi` or `e` | `-1`, the chain rule run over a symbol that cannot change | `0` | | **Silent** | `"x ^ 3".ToEntity().Differentiate(x, 2)`, and every `Differentiate(x, n)` with `n >= 1` | `(0 * x ^ 2 + 2 * x ^ 1 * 1 * 3) * 1 + 0 * 3 * x ^ 2` | `2 * x * 3` | @@ -55,6 +57,64 @@ read first. | **Silent** | `integral(x, [x, y]T)` | `[[C + x ^ 2, C + x * y]]` | `[[x ^ 2 / 2 + C, x * y + C]]` | | **Silent** | `derivative(e ^ 2, e)`, over a named constant | `0` | `2 * e` | +### The symbolic determinant is a polynomial, not a quotient by its pivots + +`Matrix.Determinant` was computed by Gaussian elimination, which leaves the pivots as literal +divisions. The expression it returned was therefore undefined wherever a pivot vanishes — at points +where the determinant itself is perfectly well defined. + +```csharp +var m = (Entity.Matrix)"[[x, 1], [2, y]]".ToEntity(); +m.Determinant // was: x * (y * x + -2) / x now: x * y + -2 +m.Determinant.Simplify() // was: x * y - 2 provided not x = 0 now: x * y - 2 +m.Determinant.Substitute("x", 0).Substitute("y", 5).Evaled // was: NaN now: -2 +``` + +`x = 0` is an ordinary point of that matrix — its determinant there is `-2` — and `NaN` asserts the +value **does not exist**. At 3×3 it stops being an edge case, because the denominator is +`a ^ 4 * (a * e - b * d)` and so two conditions have to miss: + +```csharp +// [[a, b, c], [d, e, f], [g, h, i]].Determinant, substituted: +// [[0, 1, 2], [3, 4, 5], [6, 7, 8]] was: NaN now: 0 (the pivot a is 0) +// [[1, 2, 3], [2, 4, 6], [1, 1, 1]] was: NaN now: 0 (a * e = b * d) +// [[1, 2, 3], [4, 5, 6], [7, 8, 10]] was: -3 now: -3 +// [[2, 1, 0], [1, 2, 1], [0, 1, 2]] was: 4 now: 4 +``` + +Two of those four are ordinary matrices, and the first is the singular example every linear-algebra +course opens with. Both wrong answers were **silent**: the call succeeded and returned `NaN`, which +is exactly the answer a caller checking for singularity was looking for. + +The determinant of a matrix over a commutative ring is a polynomial in its entries, so Laplace +expansion — which never divides — needs no condition at all. It is what the property's own +documentation and the comment above it already claimed was in use. + +**This is not a performance trade.** Measured on the same machine, property only, both arms built +from source: Laplace returns a *smaller* expression at every size, and the elimination it replaces +was the slower of the two on numeric matrices by a wide margin. + +| entries | Gaussian complexity | Laplace complexity | Gaussian, numeric | Laplace, numeric | +|---|---|---|---|---| +| 2×2 | 13 | 9 | | | +| 3×3 | 79 | 37 | | | +| 4×4 | 443 | 163 | | | +| 5×5 | 2461 | 833 | | | +| 6×6 | 13673 | 5021 | 164 ms | 126 ms | +| 7×7 | | 35173 | | | +| 8×8 | | | over 120 s | 358 ms | +| 10×10 | | | over 120 s | 4478 ms | + +Laplace expansion is `O(n!)`, and a fully symbolic determinant has `n!` terms however it is +computed, so that is the size of the answer rather than an overhead. The practical ceiling on a +numeric matrix moves from 7×7 to 10×10; 11×11 does not return, where under the elimination 8×8 +already did not. A fraction-free elimination (Bareiss) would raise it further and is worth having, +but it is a separate change and this one is not blocked on it. + +Measured: the whole suite passes with the fix in, and the corpus is unchanged — 116/119 with 0 +wrong, 0 error, 0 timeout, and no case's verdict or answer altered. +[#992](https://github.com/asc-community/AngouriMath/issues/992). + ### A matrix is no longer a member of every special set at once `SpecialSet.TryContains` decided membership by asking `MayContain`, which is deliberately permissive diff --git a/Sources/AngouriMath/Core/Entity/Omni/Entity.Matrix.cs b/Sources/AngouriMath/Core/Entity/Omni/Entity.Matrix.cs index 4c917bb9e..40be77e9b 100644 --- a/Sources/AngouriMath/Core/Entity/Omni/Entity.Matrix.cs +++ b/Sources/AngouriMath/Core/Entity/Omni/Entity.Matrix.cs @@ -264,7 +264,14 @@ public Entity AsScalar() private LazyPropertyA t; // We do not need to use Gaussian elimination here - // since we anyway get N! memory use + // since we anyway get N! memory use. + // Gaussian elimination is also not merely no cheaper, it is wrong here: + // it leaves the pivots as literal divisions, so the expression it returns + // is undefined wherever a pivot vanishes -- at points where the + // determinant itself is perfectly well defined. The determinant of a + // matrix over a commutative ring is a polynomial in its entries, and + // Laplace expansion never divides, so there is nothing to exclude. + // https://github.com/asc-community/AngouriMath/issues/992 /// /// Finds the symbolical determinant via Laplace's method /// @@ -273,7 +280,7 @@ public Entity AsScalar() { if (!@this.IsSquare) return null; - return @this.InnerMatrix.DeterminantGaussianSafeDivision().InnerSimplified; + return @this.InnerMatrix.DeterminantLaplace().InnerSimplified; }, this ); diff --git a/Sources/Tests/UnitTests/Algebra/MatrixTest.cs b/Sources/Tests/UnitTests/Algebra/MatrixTest.cs index 2f1ed0a8d..fc4715fc7 100644 --- a/Sources/Tests/UnitTests/Algebra/MatrixTest.cs +++ b/Sources/Tests/UnitTests/Algebra/MatrixTest.cs @@ -5,6 +5,7 @@ // Website: https://am.angouri.org. // +using System.Linq; using AngouriMath; using AngouriMath.Core.Exceptions; using Xunit; @@ -871,5 +872,56 @@ [Fact] public void Concat13() ) ) ); + + // The determinant of a matrix over a commutative ring is a polynomial in + // its entries, so it is defined wherever the entries are. Computing it by + // Gaussian elimination left the pivots as literal divisions, which made the + // returned expression undefined wherever a pivot vanishes -- so substituting + // an ordinary matrix into a symbolic determinant gave NaN for a determinant + // that exists. Laplace expansion never divides, so there is nothing to + // exclude. https://github.com/asc-community/AngouriMath/issues/992 + [Theory] + // the 2x2 of the issue: the pivot is x, and x = 0 is an ordinary point + [InlineData("[[x, 1], [2, y]]", "x", 0, "y", 5, -2)] + [InlineData("[[x, 1], [2, y]]", "x", 3, "y", 4, 10)] + // a vanishing pivot in the other position + [InlineData("[[1, x], [y, 0]]", "x", 2, "y", 0, 0)] + public void SymbolicDeterminantIsDefinedWhereAPivotVanishes( + string matrixRaw, string a, int aVal, string b, int bVal, int expected) + { + Matrix m = matrixRaw; + var det = TestExtensions.AsNotNull(m.Determinant); + Assert.Equal(expected, det.Substitute(a, aVal).Substitute(b, bVal).EvalNumerical()); + } + + [Theory] + // both of these are singular, and both defeated the divided form: the first + // has a zero in the pivot position, the second a vanishing 2x2 leading minor + [InlineData("[[0, 1, 2], [3, 4, 5], [6, 7, 8]]", 0)] + [InlineData("[[1, 2, 3], [2, 4, 6], [1, 1, 1]]", 0)] + [InlineData("[[1, 2, 3], [4, 5, 6], [7, 8, 10]]", -3)] + [InlineData("[[2, 1, 0], [1, 2, 1], [0, 1, 2]]", 4)] + public void GeneralSymbolicDeterminantAgreesWithTheNumericOne(string matrixRaw, int expected) + { + Matrix numeric = matrixRaw; + Matrix symbolic = "[[a, b, c], [d, e, f], [g, h, i]]"; + var det = TestExtensions.AsNotNull(symbolic.Determinant); + foreach (var (name, value) in new[] { "a", "b", "c", "d", "e", "f", "g", "h", "i" } + .Select((name, index) => (name, numeric[index / 3, index % 3]))) + det = det.Substitute(name, value); + Assert.Equal(expected, det.EvalNumerical()); + Assert.Equal(expected, TestExtensions.AsNotNull(numeric.Determinant).EvalNumerical()); + } + + [Theory] + [InlineData("[[x, 1], [2, y]]", "x * y - 2")] + [InlineData("[[a, b], [c, d]]", "a * d - b * c")] + public void SymbolicDeterminantCarriesNoCondition(string matrixRaw, string expected) + { + Matrix m = matrixRaw; + var det = TestExtensions.AsNotNull(m.Determinant).Simplify(); + Assert.DoesNotContain(det.Nodes, node => node is Providedf); + Assert.Equal(0, ((Entity)expected - det).Simplify().EvalNumerical()); + } } }