This was discovered while digging for the underlying reason behind #33156:

sage: E = EllipticCurve('11a3')
sage: R.<X> = QQ[]
sage: S.<Y> = R.quotient(X^2)
sage: E.division_polynomial(5, x=Y)
-5*Y
sage: E.division_polynomial(5, x=X)
-5*Y

The problem is that Y == X, which makes the method return an incorrect result from its cache.

In this patch, we streamline caching in three major ways:

• Only cache the generic division polynomials, rather than evaluations at all inputs. (Currently all input-output pairs are stored, and this is both a source of errors — as above — and a source of memory leaks.)
• Share the caches of .division_polynomial_0() and .division_polynomial(). This makes later calls to .division_polynomial() (with x==None) faster and faster as the cache is populated with more intermediate nodes appearing in the recursive formula for division polynomials.
• In a similar vein, we make opportunistic use of already computed division polynomials, for instance when evaluating at base-field elements (and the degree is small enough for this to be beneficial), or when computing another division polynomial with a different two-torsion multiplicity.

Depends on #33165

CC: @defeo @JohnCremona

Component: elliptic curves

Author: Lorenz Panny

Branch/Commit: 0aaec67

Reviewer: Travis Scrimshaw

Issue created by migration from https://trac.sagemath.org/ticket/33164