Give affine schemes unique representation (needed for elliptic curves and forking)

[Jan Keitel]

The following used to work in Sage, but does not anymore.

Define

@fork
def compute_E():
    E = EllipticCurve([2,3])
    p = E(3,6,1)
    return p

Then

sage: p = compute_E()
sage: 2*p
---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
<ipython-input-5-a9d49bd49c00> in <module>()
----> 1 Integer(2)*p

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/element.so in sage.structure.element.RingElement.__mul__ (build/cythonized/sage/structure/element.c:16558)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/coerce.so in sage.structure.coerce.CoercionModel_cache_maps.bin_op (build/cythonized/sage/structure/coerce.c:7839)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/coerce.so in sage.structure.coerce.CoercionModel_cache_maps.get_action (build/cythonized/sage/structure/coerce.c:13470)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/coerce.so in sage.structure.coerce.CoercionModel_cache_maps.discover_action (build/cythonized/sage/structure/coerce.c:14918)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/parent.so in sage.structure.parent.Parent.get_action (build/cythonized/sage/structure/parent.c:20486)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/parent.so in sage.structure.parent.Parent.discover_action (build/cythonized/sage/structure/parent.c:21783)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/coerce_actions.so in sage.structure.coerce_actions.IntegerMulAction.__init__ (build/cythonized/sage/structure/coerce_actions.c:8432)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/element.so in sage.structure.element.ModuleElement.__add__ (build/cythonized/sage/structure/element.c:11295)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/coerce.so in sage.structure.coerce.CoercionModel_cache_maps.bin_op (build/cythonized/sage/structure/coerce.c:7904)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/coerce.so in sage.structure.coerce.CoercionModel_cache_maps.canonical_coercion (build/cythonized/sage/structure/coerce.c:9490)()

/home/pcl337b/jkeitel/sage/sage/local/lib/python2.7/site-packages/sage/structure/coerce.so in sage.structure.coerce.CoercionModel_cache_maps._coercion_error (build/cythonized/sage/structure/coerce.c:16376)()

RuntimeError: There is a bug in the coercion code in Sage.
Both x (=(3 : 6 : 1)) and y (=(3 : -6 : 1)) are supposed to have identical parents but they don't.
In fact, x has parent 'Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field'
whereas y has parent 'Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field'
Original elements (3 : 6 : 1) (parent Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field) and (3 : -6 : 1) (parent Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field) and maps
<type 'NoneType'> None
<type 'sage.structure.coerce_maps.DefaultConvertMap_unique'> (map internal to coercion system -- copy before use)
Conversion map:
  From: Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field
  To:   Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field
sage: 

and

sage: (-p).parent() is p.parent()
sage: False

Note that this depends crucially on the @fork, without it things work just fine. Intuitively, I suppose that since E is created in the forked process, the main process where I'm negating p does not know about E and therefore creates a new parent. However, I'm not sure how to fix this.

Since this definitely worked a few months ago, my gut says that it might have something to do with the changes introduced in #11474 and I've taken the liberty to include some of you who worked on this. I hope that's okay for you.

[Nils Bruin]

It looks like the problem is in the point sets, not the elliptic curves:

sage: A=parent(p)
sage: B=parent(-p)
sage: id(A)==id(B)
False
sage: A._codomain
Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field
sage: id(A._codomain)==id(B._codomain)
True

In fact, even the "Spectrum of Ration Field" is multiply instantiated (it arises as domain of the point set. Do we really get any mileage out of carrying around the formal homset stuff?

sage: A._domain
Spectrum of Rational Field
sage: id(A._domain)==id(B._domain)
False

We can create nonidentical-but-equal pointsets without fork:

sage: E=EllipticCurve([2,3])
sage: A=E.point_set()
sage: B=E.point_set(QQ)
sage: parent(A(3,6,1)+B(3,6,1)) is A
True
sage: parent(B(3,6,1)+B(3,6,1)) is A
True
sage: A is B
False
sage: A==B
True

[John Cremona]

Replying to nbruin:

Do we really get any mileage out of carrying around the formal homset stuff?

I use the elliptic curve & point stuff in Sage all the time, and never have any use for that. I don't use point-sets at all.

[Nils Bruin]

Replying to cremona:

I use the elliptic curve & point stuff in Sage all the time, and never have any use for that. I don't use point-sets at all.

Not explicitly perhaps, but it does get created (unless "parent" triggers something lazy here):

sage: parent(EllipticCurve([2,3])(3,6,1))
Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field

which is why I'm slightly concerned about "carrying around" this stuff. It may be OK: The point set should obviously carry around the field (in the form of the ._domain attribute) and it may be that "Spectrum of Rational Field" is a sufficiently lightweight/permanent wrapper that it doesn't hurt too much.

[Nils Bruin]

It looks to me the problem with point sets indeed explains the problem: point_homset is a cached_method, but the "fork" construction somehow circumvents the cache:

sage: p=compute_E()
sage: p.parent()._codomain.point_homset.cache
{}
sage: 2*p
...
sage: p.parent()._codomain.point_homset.cache
{((None,),
  ()): Abelian group of points on Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field}

A new pointset was created! If we avoid this, everything works fine:

sage: p.parent()._codomain.point_homset.set_cache(parent(p))
sage: 2*p
(-23/144 : 2827/1728 : 1)

This is probably explained by:

sage: p.parent().__reduce_ex__(2)
(<class 'sage.schemes.generic.homset.SchemeHomsetFactory'>,
 (Spectrum of Rational Field,
  Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Rational Field,
  Category of schemes over Integer Ring,
  Integer Ring,
  False,
  True))
sage: p.parent()._codomain.__reduce_ex__(2)
(<function sage.structure.factory.generic_factory_unpickle>,
 (<class 'sage.schemes.elliptic_curves.constructor.EllipticCurveFactory'>,
  (6, 4, 'beta2'),
  (Rational Field, (0, 0, 0, 2, 3)),
  {}))

so, when p gets pickled, then p.parent() gets pickled and as a result the elliptic curve as well, but the cache doesn't get pickled (as you can see, the elliptic curve gets created "fresh").

The problem could be solved by making point sets unique, but probably you'd have to make spec unique in the process as well; it is not at the moment:

sage: a=Spec(QQ); b=Spec(QQ)
sage: a is b , a == b
(False, True)

Incidentally, this also explains:

sage: E=EllipticCurve([2,3])
sage: A=E.point_homset()
sage: B=E.point_homset(QQ)
sage: A is B
False
sage: [A is c for c in E.point_homset.cache.values()]
[True, False]

two entries in the cache: one for each input type. However, these point sets seem a little better coordinated than the ones resulting from unpickling:

sage: 2*A(3,6,1)
(-23/144 : 2827/1728 : 1)
sage: 2*B(3,6,1)
(-23/144 : 2827/1728 : 1)
sage: A(3,6,1)+B(3,6,1)
(-23/144 : 2827/1728 : 1)

[Peter Bruin]

Replying to nbruin:

The problem could be solved by making point sets unique, but probably you'd have to make spec unique in the process as well; it is not at the moment:

sage: a=Spec(QQ); b=Spec(QQ)
sage: a is b , a == b
(False, True)

Indeed. I am now experimenting with making AffineScheme inherit from UniqueRepresentation, and in fact this already gives the desired results without any special changes to point sets:

sage: @fork
....: def compute_E():
....:     E = EllipticCurve([2,3])
....:     p = E(3,6,1)
....:     return p
....:
sage: p = compute_E()
sage: 2*p
(-23/144 : 2827/1728 : 1)

sage: a=Spec(QQ); b=Spec(QQ)
sage: a is b , a == b
(True, True)

sage: E=EllipticCurve([2,3])
sage: A=E.point_homset()
sage: B=E.point_homset(QQ)
sage: A is B
True

[Simon King]

Replying to pbruin:

Indeed. I am now experimenting with making AffineScheme inherit from UniqueRepresentation, and in fact this already gives the desired results without any special changes to point sets:

Isn't UniqueRepresentation cool? :-D

But as you know, it can introduce difficult to debug memory leaks in some cases...

[Peter Bruin]

Replying to SimonKing:

Replying to pbruin:

Indeed. I am now experimenting with making AffineScheme inherit from UniqueRepresentation, and in fact this already gives the desired results without any special changes to point sets:

Isn't UniqueRepresentation cool? :-D

It seems to work quite well, except it took me some time to figure out the following caveat:

sage: class X(UniqueRepresentation):
....:     def __init__(self, x, y=None):
....:         self.x = x
....:         self.y = y
....:
sage: A = X(1)
sage: B = X(1)
sage: C = X(1, None)
sage: A is B
True
sage: B is C
False

Apparently it makes a difference whether the default None argument is explicitly given or not. I accept that this may be a technical limitation that isn't worth fixing, though.

[Simon King]

Replying to pbruin:

Apparently it makes a difference whether the default None argument is explicitly given or not. I accept that this may be a technical limitation that isn't worth fixing, though.

Certainly not here. Anyway, the problem is that caching happens in the __classcall__ method that is inherited from CachedRepresentation, and this is unaware of default arguments to the __init__ method.

The point is that @cached_function, that is used to decorate the __classcall__, only takes into account the argspec of __classcall__ (which is generic), but not the argspec of __init__.

Default arguments are taken into account by a so-called argument fixer. Normally, this is created using the argspec of the to-be-wrapped method (here: __classcall__). But it could be possible to find a mechanism to create an argument fixer using the argspec of a different method.

Problem: Suppose you have a hierarchy of classes. They all inherit the __classcall__ from CachedRepresentation, but have their own __init__ (or some just inherit it), and thus the classes should inherit different instances of the cached __classcall__ function, each with its own argument fixer.

One may think that this could be possible to do, for example, in sage.misc.classcall_metaclass.ClasscallMetaclass.__cinit__: This is where the classcall is assigned to a class. But I am not sure if __init__ is available at that point. After all, the class for which we want to inspect __init__ is just to-be-created.

Perhaps this is worth a different ticket, perhaps not (because of the bad ration of complication and gain).

[Travis Scrimshaw]

@pbruin

So what Simon is recommending (I believe) is to add a method:

@staticmethod
def __classcall__(cls, R, S=None):
    # Do any other parsing/standardization that's needed here
    # In particular, if S is None, then set it to ZZ following
    #   the __init__ of Scheme
    return super(AffineScheme, cls).__classcall__(R, S)

or perhaps call it __classcall_private__ (so it doesn't get accidentally used by subclasses). It almost seems like this should go one further to make Scheme a UniqueRepresentation, but IDK how much more work that would be.

However I disagree with this change:

src/sage/categories/homset.py

@@ -505 +505 @@
         sage: loads(H.dumps()) is H
         True
 
-    Homsets of non-unique parents are non-unique as well::
+    Homsets of unique parents are unique as well::
 
         sage: H = End(AffineSpace(2, names='x,y'))
         sage: loads(dumps(AffineSpace(2, names='x,y'))) is AffineSpace(2, names='x,y')
-        False
-        sage: loads(dumps(AffineSpace(2, names='x,y'))) == AffineSpace(2, names='x,y')
         True
         sage: loads(dumps(H)) is H
-        False
-        sage: loads(dumps(H)) == H
         True
 
     """

I think it would be much better to instead change it to a non-unique parent.

EDIT (clarification) - By 'it' above I meant the parent in the test.

[Simon King]

Replying to tscrim:

So what Simon is recommending (I believe) is to add a method:

@staticmethod
def __classcall__(cls, R, S=None):
    # Do any other parsing/standardization that's needed here
    # In particular, if S is None, then set it to ZZ following
    #   the __init__ of Scheme
    return super(AffineScheme, cls).__classcall__(R, S)

or perhaps call it __classcall_private__ (so it doesn't get accidentally used by subclasses).

That's the way to work around the limitation of UniqueRepresentation in special cases. Note that what we want here is to make __classcall__ aware of the default arguments of __init__ (or rather: Provide these default arguments). Hence, we want that this behaviour is inherited for all sub-classes that also inherit __init__. Therefore, __classcall__ should be better than __classcall_private__ in this case.

However, what I was mentioning was an approach to leverage the limitation generally, without the need to implement __classcall__ or __classcall_private__.

[Simon King]

Replying to tscrim:

However I disagree with this change:

src/sage/categories/homset.py

@@ -505 +505 @@
         sage: loads(H.dumps()) is H
         True
 
-    Homsets of non-unique parents are non-unique as well::
+    Homsets of unique parents are unique as well::
 
         sage: H = End(AffineSpace(2, names='x,y'))
         sage: loads(dumps(AffineSpace(2, names='x,y'))) is AffineSpace(2, names='x,y')
-        False
-        sage: loads(dumps(AffineSpace(2, names='x,y'))) == AffineSpace(2, names='x,y')
         True
         sage: loads(dumps(H)) is H
-        False
-        sage: loads(dumps(H)) == H
         True
 
     """

I think it would be much better to instead change it to a non-unique parent.

Homsets are cached by identity of domain and codomain. Hence, if domain or codomain are non-unique parents, then homsets automatically are non-unique parents as well. However, I think it is true (and already documented, that's why I disagree with this change) that homsets of unique parents are unique parents.

[Nils Bruin]

I'd expect that caching strategies need amending if schemes move over to UniqueRepresentations. As we know, if one has an object A from which another object b=B(A,...) is constructed then if b is UniqueRepresentation, it should NOT be cached on A. That's because the global UniqueRepresentation cache holds a strong reference to A for the lifetime of b (i.e., A will not die before b).

If A also holds a strong reference to b (in its local cache) then b will also not die before A: they're immortal! The garbage collector does not recognize this, because the reference to A is truly global, so it's not a cycle (it's key in a WeakValueDict, so the reference is artificially removed should b be deallocated).

I haven't checked, but the fact that point sets are cached on the scheme seems likely to instantiate a scenario as above.

UniqueRepresentation is nice for coercion, but it's a pain for memory management. Most memory management issues in Python are alleviated by the cyclic garbage collector (CGC) (reference counting simply doesn't cut it), but the mechanisms involved in UniqueRepresentation can defeat the CGC. That's not pleasant.

You'd better have good reasons for condemning data structures to UniqueRepresentation. The original problem on this ticket is more a pickling problem. It can also be solved by amending the pickling data structures (e.g., by making point sets pickle to point set constructors rather than a generic type __new__ and filling the __dict__, or by making points pickle differently)

[Nils Bruin]

It seems we're currently not leaking pointsets and elliptic curves themselves. We are leaking all kinds of other structures, though, so constructing many reductions of an elliptic curve will not be possible. It's also scary that we're creating multivariate polynomial rings (and leaking them!) just to construct a point on an elliptic curve (in fact all this junk already is created if just the elliptic curve is created).

sage: import gc
sage: from collections import Counter
sage: pre = {id(a) for a in gc.get_objects()}
sage: for p in prime_range(5,2000):
....:       k=GF(p)
....:       E=EllipticCurve(k,[0,1])
....:       V=E(0,1)
....:     
sage: gc.collect()
2297
sage: gc.collect()
0
sage: new = Counter(str(type(c)) for c in gc.get_objects() if id(c) not in pre)
sage: for a,b in new.iteritems(): print a,b
<type 'sage.misc.cachefunc.CachedFunction'> 1
<type 'sage.misc.cachefunc.CachedMethodCallerNoArgs'> 5
<type 'frame'> 2
<class 'weakref.KeyedRef'> 8751
<type 'weakref'> 3384
<type 'set'> 1
<type 'method_descriptor'> 19
<type 'sage.misc.constant_function.ConstantFunction'> 1204
<class '_ast.comprehension'> 1
<class '_ast.Compare'> 1
<type 'sage.structure.coerce_dict.TripleDict'> 907
<class '_ast.Attribute'> 1
<class '_ast.Interactive'> 1
<class 'sage.categories.schemes.Schemes_over_base_with_category'> 1
<type 'sage.structure.coerce_actions.LeftModuleAction'> 301
<class 'sage.categories.commutative_algebras.CommutativeAlgebras_with_category'> 1
<class 'sage.rings.ideal.Ideal_pid'> 301
<type 'tuple'> 4504
<type 'sage.rings.finite_rings.integer_mod.NativeIntStruct'> 301
<class 'sage.structure.dynamic_class.DynamicMetaclass'> 14
<type 'sage.misc.cachefunc.CachedMethodCaller'> 3
<type 'sage.rings.polynomial.multi_polynomial_libsingular.MPolynomialRing_libsingular'> 301
<class 'sage.rings.polynomial.term_order.TermOrder'> 602
<type 'builtin_function_or_method'> 608
<type 'module'> 3
<type 'wrapper_descriptor'> 18
<class 'sage.rings.finite_rings.finite_field_prime_modn.FiniteField_prime_modn_with_category'> 301
<class 'sage.rings.homset.RingHomset_quo_ring_with_category'> 301
<class 'sage.schemes.projective.projective_space.ProjectiveSpace_finite_field_with_category'> 1
<type 'dict'> 1576
<type 'sage.rings.polynomial.polynomial_element.PolynomialBaseringInjection'> 301
<type 'cell'> 647
<type 'sage.rings.polynomial.multi_polynomial_libsingular.MPolynomial_libsingular'> 1807
<type 'sage.structure.coerce_dict.MonoDictEraser'> 1814
<type 'list'> 6653
<type 'sage.misc.cachefunc.CachedMethod'> 1
<class '_ast.Call'> 5
<class 'sage.structure.dynamic_class.DynamicClasscallMetaclass'> 309
<type 'sage.structure.coerce_dict.TripleDictEraser'> 907
<type 'instancemethod'> 302
<type 'sage.rings.finite_rings.integer_mod.IntegerMod_int'> 22373
<type 'function'> 359
<class 'sage.categories.category.JoinCategory_with_category'> 5
<class '_ast.Module'> 1
<class 'sage.schemes.elliptic_curves.ell_finite_field.EllipticCurve_finite_field_with_category'> 1
<class '_ast.Name'> 9
<type 'listiterator'> 1
<type 'sage.structure.coerce_dict.MonoDict'> 1814
<type 'type'> 2
<type 'staticmethod'> 323
<class '_ast.Assign'> 1
<class 'sage.schemes.generic.scheme.AffineScheme_with_category'> 2
<type 'sage.rings.finite_rings.integer_mod.Int_to_IntegerMod'> 301
<type 'frozenset'> 2
<type 'sage.rings.finite_rings.integer_mod.Integer_to_IntegerMod'> 301
<class 'sage.categories.groupoid.Groupoid_with_category'> 301

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​daaaee8

Trac 17008: add doctest for homsets of non-unique parents

[Peter Bruin]

Replying to tscrim:

It almost seems like this should go one further to make Scheme a UniqueRepresentation, but IDK how much more work that would be.

I thought about this, but I think it does not make sense, because Scheme is just an abstract base class for general schemes, which don't usually have a set of defining data that determine them uniquely up to isomorphism. Affine schemes are uniquely determined by their coordinate ring, elliptic curves by their base ring and Weierstrass coefficients, but general schemes are just too general. We could do this for affine and projective spaces Aⁿ and Pⁿ, though.

However I disagree with this change: [...] I think it would be much better to instead change it to a non-unique parent.

Done in previous commit.

[Peter Bruin]

Replying to nbruin:

It seems we're currently not leaking pointsets and elliptic curves themselves. We are leaking all kinds of other structures, though, so constructing many reductions of an elliptic curve will not be possible.

You are right, and this is a reason why people are forced to use forking for long computations...

It's also scary that we're creating multivariate polynomial rings (and leaking them!) just to construct a point on an elliptic curve (in fact all this junk already is created if just the elliptic curve is created).

In this case it is probably because an elliptic curve is also a projective curve. From EllipticCurve_generic:

def __init__(self, K, ainvs):
    self.__base_ring = K
    self.__ainvs = tuple(K(a) for a in ainvs)
    if self.discriminant() == 0:
        raise ArithmeticError("invariants " + str(ainvs) + " define a singular curve")
    PP = projective_space.ProjectiveSpace(2, K, names='xyz');
    x, y, z = PP.coordinate_ring().gens()
    a1, a2, a3, a4, a6 = ainvs
    f = y**2*z + (a1*x + a3*z)*y*z \
        - (x**3 + a2*x**2*z + a4*x*z**2 + a6*z**3)
    plane_curve.ProjectiveCurve_generic.__init__(self, PP, f)

[Travis Scrimshaw]

Replying to pbruin:

I thought about this, but I think it does not make sense, because Scheme is just an abstract base class for general schemes, which don't usually have a set of defining data that determine them uniquely up to isomorphism. Affine schemes are uniquely determined by their coordinate ring, elliptic curves by their base ring and Weierstrass coefficients, but general schemes are just too general. We could do this for affine and projective spaces Aⁿ and Pⁿ, though.

Ah okay (I don't know that much about the math in this area), I was looking at the code and it seems like all of input standardization (the _base_* attributes) could be done at initialization for very low cost (really beforehand). I think == doesn't check for isomorphism (at least in the base class and there are only 3 classes that I can see that directly inherit from Scheme), so I feel like making it into a UniqueRepresentation would be okay. (In fact, it seems like no __eq__ is defined on any scheme, so it basically always returns false.) *shrugs*

However I disagree with this change: [...] I think it would be much better to instead change it to a non-unique parent.

Done in previous commit.

Thanks!

Replying to nbruin:

It seems we're currently not leaking pointsets and elliptic curves themselves. We are leaking all kinds of other structures, though, so constructing many reductions of an elliptic curve will not be possible.

I thought UniqueRepresentation's are weakly cached, so these aren't truly leaked, but instead waiting for Sage to use up all it's memory before collecting them? I know polynomial rings use their own (weak) cache and exhibit similar output to the above.

Related question: Do we have a (documented) mechanism for freeing the UniqueRepresentation cache? Should we if we don't?

[Nils Bruin]

Replying to tscrim:

I thought UniqueRepresentation's are weakly cached, so these aren't truly leaked, but instead waiting for Sage to use up all it's memory before collecting them? I know polynomial rings use their own (weak) cache and exhibit similar output to the above.

They are: as long as the object lives, we remember it in a WeakValueDict, keyed on construction parameters. When the object is deallocated, we remove the entry from the cache. The WeakValueDict reference by itself doesn't prevent deallocation of the object. However, the key IS strongly referenced, so if the keys somehow hold a reference to the original object, the cache is will keep the object alive. See the scenario I described in comment 13 or ​sage-devel for details.

Of course, normally one wouldn't expect that construction parameters refer to the object that is constructed from them. However, local caching strategies tend to do exactly that.

Another issue is that libsingular polynomial rings are, unfortunately, immortal. See #13447. That may well explain all the leaked objects for the particular example. The main issue is: there are plenty of things one can do with elliptic curves without explicitly referring to a plane projective model of it, so it seems a little worrying that we're always creating that.

Related question: Do we have a (documented) mechanism for freeing the UniqueRepresentation cache? Should we if we don't?

Do you mean: removing an entry from the cache even when the object still exists? That sounds like a bad idea to sanction with a formal API. However, you can do it if you want, since you can find the cache if you work at it:

sage: A=AbelianGroup([0])
sage: k=[k for k,v in sage.structure.unique_representation.CachedRepresentation.__classcall__.cache.iteritems() if v is A][0]; k
((sage.groups.abelian_gps.abelian_group.AbelianGroup_class, (0,), 'f'), ())
sage: del sage.structure.unique_representation.CachedRepresentation.__classcall__.cache[k]
sage: B=AbelianGroup([0])
sage: A is B
False

However, sage behaviour after this is likely to be questionable (so, basically no worse than before).

[Travis Scrimshaw]

Replying to nbruin:

They are: as long as the object lives, we remember it in a WeakValueDict, keyed on construction parameters. When the object is deallocated, we remove the entry from the cache. The WeakValueDict reference by itself doesn't prevent deallocation of the object. However, the key IS strongly referenced, so if the keys somehow hold a reference to the original object, the cache is will keep the object alive. See the scenario I described in comment 13 or ​sage-devel for details.

Ahh, I see what the problem is. One of these days, this will have been explained to me enough times I'll remember it. :P

Another issue is that libsingular polynomial rings are, unfortunately, immortal. See #13447. That may well explain all the leaked objects for the particular example. The main issue is: there are plenty of things one can do with elliptic curves without explicitly referring to a plane projective model of it, so it seems a little worrying that we're always creating that.

So then would you expect #13447 to help with the leak here? Is it worthwhile to try and finish #13447?

How bad do we expect the leak here to be in practice? Mainly, what's the lesser evil: the error this ticket's about or the memory leak, and should we just have a follow-up ticket if the latter is the lesser evil?

Related question: Do we have a (documented) mechanism for freeing the UniqueRepresentation cache? Should we if we don't?

Do you mean: removing an entry from the cache even when the object still exists? That sounds like a bad idea to sanction with a formal API. However, you can do it if you want, since you can find the cache if you work at it:

sage: A=AbelianGroup([0])
sage: k=[k for k,v in sage.structure.unique_representation.CachedRepresentation.__classcall__.cache.iteritems() if v is A][0]; k
((sage.groups.abelian_gps.abelian_group.AbelianGroup_class, (0,), 'f'), ())
sage: del sage.structure.unique_representation.CachedRepresentation.__classcall__.cache[k]
sage: B=AbelianGroup([0])
sage: A is B
False

However, sage behaviour after this is likely to be questionable (so, basically no worse than before).

I was misunderstanding the weak cache as I thought they didn't get garbage collected until space was needed (even if there were no strong references). Nevermind.

[Peter Bruin]

Replying to tscrim:

So then would you expect #13447 to help with the leak here? Is it worthwhile to try and finish #13447?

Just to clarify, #13447 is orthogonal to this ticket; the leak occurs both with and without the attached branch. I checked that the leak goes away after applying

src/sage/rings/polynomial/multi_polynomial_libsingular.pyx

@@ -245 +245 @@
 import sage.libs.pari.gen
 import polynomial_element
 
-permstore=[]
 cdef class MPolynomialRing_libsingular(MPolynomialRing_generic):
 
     def __cinit__(self):
@@ -367 +366 @@
         from sage.rings.polynomial.polynomial_element import PolynomialBaseringInjection
         base_inject = PolynomialBaseringInjection(base_ring, self)
         self.register_coercion(base_inject)
-        #permanently store a reference to this ring until deallocation works reliably
-        permstore.append(self)
 
     def __dealloc__(self):
         r"""

which should be done as part of #13447. My branch does not introduce any additional leak as far as I can see.

How bad do we expect the leak here to be in practice? Mainly, what's the lesser evil: the error this ticket's about or the memory leak, and should we just have a follow-up ticket if the latter is the lesser evil?

It seems to me that #13447 and this ticket can and should be fixed independently of each other, and there is currently no separate problem requiring another ticket.

[Nils Bruin]

Replying to pbruin:

Just to clarify, #13447 is orthogonal to this ticket; the leak occurs both with and without the attached branch. I checked that the leak goes away after (removing permstore). which should be done as part of #13447. My branch does not introduce any additional leak as far as I can see.

That doesn't prove that properly solving #13447 would leave the memory leak in place. permstore was only introduced because we had examples where the deletion of polynomial rings caused memory corruption in libsingular. We don't know that permstore is the only link that prevents them from being deleted (and since we're not testing that, it's quite probable it's not). Memory management of polynomial rings is just broken.

It seems to me that #13447 and this ticket can and should be fixed independently of each other, and there is currently no separate problem requiring another ticket.

I do agree with that.

I do find it funny that with the fork/pickle idiom, we're basically reintroducing an extremely heavy-handed version of PARI's grepile memory management (except that pickling is an even harder problem to solve and, as we can see here, is probably broken for pretty much any sufficiently complicated data structure).