[tombouctou]

Stacked on #2073: the base is its branch. Once #2073 is merged, this one moves to master. It removes the 40 regressions that #2073 lists, and with them the xfail that #2073 added. Part of #2071.

In WMA, every factor of a product is a variable of the monomial, not only the symbols. Mathics3 skipped factors such as f[x] or Erf[1/x] when it built the key of a product. This PR counts them, as WMA 15 does:

• An expression that is not a symbol is a variable (GeneralVariable in keycomparable.py). These variables come after the symbols and are ordered among themselves as general expressions (head, length, elements). A general expression on its own is a monomial in itself.

  Sort[{x f[x], Pi g[x], y f[x], a f[x], f[x], z, g[x], f[x]^2, 1 + f[x]}]
     {z, f[x], a f[x], x f[x], y f[x], f[x]^2, 1 + f[x], g[x], Pi g[x]}
  Sort[{a Sin[x], b Cos[x]}]                 {b Cos[x], a Sin[x]}
  Integrate[Exp[-1/x^2], x]                  E^(-1/x^2) x + Sqrt[Pi] Erf[1/x]
  

• A power with an exponent that is not a number comes after the numeric powers of its base: Sort[{E^x, E, E^2, F, x}] is {E, E^2, E^x, F, x}, and x E^x is E^x x.
• A sum, alone or as a factor, sorts right after the largest variable of its largest term. Its key extends that variable's key with the terms. This replaces the infinite exponent from Canonical order: sums, symbols from other contexts and coefficients, as in WMA #2073, so that a sum as a factor goes there too:

  Collect[3 x y + 2 Sin[x z] (x + 2 y^2 + x) + (x + y)^3, y]
     x^3 + (3 x + 3 x^2) y + y^3 + 4 x Sin[x z] + y^2 (3 x + 4 Sin[x z])
  Expand[(x + a)^2 + (y + a)^2 + (x + y) (x + a), y]
     a^2 + x (a + x) + (a + x)^2 + 2 a y + (a + x) y + y^2
  

  Both are now exactly what WMA gives.

How much closer to WMA

Order on the same 3003 pairs as in #2073, compared with WMA 15:

	pairs that differ from WMA
master	1152
#2073	425
this PR	20

No pair that master or #2073 gets right goes wrong. The 20 that remain are cases the WMA rule for which is not clear yet (see #2071):

• a sum whose largest term has a negative coefficient, or with a negative constant (1 - x, -1 + x, x - y);
• "s" against Sqrt[2].

Also still as before: a power with a numeric base and a symbolic exponent (2^x), and numeric expressions such as -Sqrt[3] against Sqrt[3].

Doctests and tests

• 47 doctests had the old order as expected output. I ran every new output back through WMA 15 (parsed held, evaluated by WMA, compared as InputForm). In 35 of them WMA keeps exactly the Mathics3 order; only 5 of the old expectations had it. In the other 12, WMA simplifies differently, e.g. Sqrt[2/Pi], 1/E^x or a third argument of Root, and the order of the terms is the same. Two differ in order and belong to the known cases above (2^Re[z], -Sqrt[3]).
• pytest: expectations updated in test_basic.py, test_algebra.py, test_calculus.py, test_diffeqns.py and test_linalg.py. I checked each new output against WMA 15; the test_linalg eigenvalues match WMA up to parentheses. test_canonical_order gets 7 more cases; 6 of them fail on Canonical order: sums, symbols from other contexts and coefficients, as in WMA #2073.
• Full run on kimmy: doctests 4622 passed, 0 failed. pytest with -n 4 --dist loadfile: 20 failed, in test_tensor.py, test_series.py, test_list.py, test_testing_expressions.py and test_latex.py. Each of these files passes on its own; the Max[x] failure, for example, is an x = 3 left by another file in the same worker.

Speed

With general expressions now monomials, comparisons got deeper. To keep this from slowing things down:

• element_order is kept on the expression, together with the head and elements it was built from, so a changed expression rebuilds it;
• Sort and Orderless sort with element_order as a key, which builds each key once instead of at every comparison;
• Monomial has a fast path for a single variable and an == that compares the exponents. A tuple comparison calls == before <, and before this both ran the full comparison.

On kimmy, best of 7, master vs this PR:

	master	this PR
Sort of 3000 f[n, g[k]]	125 ms	216 ms
Table[f[Mod[7 k, 101], g[k]], {k, 3000}]	2017 ms	2022 ms
Plus @@ Table[Sin[k x], {k, 1500}]	1528 ms	1471 ms
Total[Table[a[k] x^k, {k, 1500}]]	382 ms	378 ms
Expand[(1 + a + b + c + x + y)^6]	229 ms	171 ms
Expand[(a + b + c + d)^8]	99 ms	77 ms
Times @@ Table[f[k]^k, {k, 1000}]	749 ms	762 ms
Sort of 3000 x^i y^j	4156 ms	3866 ms

The one slowdown left is a bare Sort of general expressions, 1.7×. Sums and products of such expressions are as fast as before or faster.

🤖 Generated with Claude Code

[rocky]

Right now, the slow behavior I do not feel is a big problem. First, get it right, and then get it fast.

My current belief is that if we get rule matching fixed and we spend $O(1)$ time instead of $O(n)$ time for $n$ fixed rule values of a definition, this will more than compensate for this.

The other big wins I believe will be in handling NumericArray as NumPy data (needed in plotting) and expression compilation.

So LGTM. In fact, this is a huge help!

@mmatera - your thoughts?

[rocky]

One small problem:

The Expression.elements setter invalidates _cache, but not _element_order_cache.

set_head() invalidates head changes, but not _element_order_cache; element_order is used as the comparison key for Sort, Orderless, and ordered expressions. If an expression is rebuilt/restructured and then reused without a fresh cache, it has an old key and sorts incorrectly.

So it would be good to clear _element_order_cache in Expression.elements, and clear it in set_head().

[tombouctou]

Done in 7575db2: both the elements setter and set_head() now clear _element_order_cache. The cached key was already checked against the head and elements it was built from, so a new head or elements tuple was rebuilt anyway. That check stays, for code that sets _elements or _head directly (box layouts, ListExpression). Added a test that changes elements, head and a single element.