Adding Fmpz mpoly

[deinst]

This is a more work in progress than a request to change the code. I'd like to get some comments before building this out completely. A few questions:

Do we want some way of specifying a polynomial with a list of tuples of the form (fmpz, (exp1, ..., expn)) ? What would be the best format for this?

I'm not certain when I should use NotImplemented vs ValueError when for example adding polynomials from different context and similar situations.

If this looks ok, I'll forge on and get the rest of fmpz_mpoly working and do the same for fmpq_mpoly and fmpz_mpoly_q

[oscarbenjamin]

Do we want some way of specifying a polynomial with a list of tuples of the form (fmpz, (exp1, ..., expn)) ? What would be the best format for this?

For SymPy the corresponding Python format looks like:

In [4]: p = QQ[x,y].convert(x**2*y + 1)

In [5]: p
Out[5]: x**2*y + 1

In [6]: dict(p)
Out[6]: {(0, 0): 1, (2, 1): 1}

In [7]: type(dict(p)[(0,0)])
Out[7]: flint._flint.fmpq

Here p is a dict (actually a subclass of dict). The keys are tuples of python int and the values are fmpq if python-flint is being used.

I would definitely like to be able to just pass that p straight to fmpz_mpoly somehow like maybe

fmpz_mpoly(p, ctx)

or perhaps a method on the context object:

ctx.from_dict(p)

I discussed the idea of a higher-level interface for all of this in gh-53 but we do not need to do that here e.g. this is like level 1 referred to at #55 (comment)

I'm not certain when I should use NotImplemented vs ValueError when for example adding polynomials from different context and similar situations.

NotImplemented is used for any types that you do not recognise. It allows the other type to potentially overload the method instead and it is good practice to do that for any type that you have not explicitly contemplated. Note that if both operands give NotImplemented then in the end a TypeError will be raised so returning NotImplemented usually results in TypeError.

ValueError would be when you recognise the type but the operation is not valid e.g. fmpz_mpoly with different contexts.

So generally just check for the types that you want to allow and otherwise return NotImplemented. If you find a type that is expected and allowed but the value is not compatible (different context) raise ValueError.

get the rest of fmpz_mpoly working and do the same for fmpq_mpoly and fmpz_mpoly_q

I think it probably makes sense to have a common base class for the mpoly types. The situation would be like:

cdef class mpoly_base:
    # any generic code that just calls other methods
    # no C functions are called here

cdef class fmpz_mpoly(mpoly_base):
     cdef fmpz_mpoly_t val

    # methods that actually call C functions.

It would be easier to ensure a common interface like this and the idea would be that some of the code can be shared between the three classes.

[oscarbenjamin]

There seem to be various unrelated changes bundled in here. This one is possibly what breaks building the wheels in CI.

[oscarbenjamin]

It is better to keep different changes like this in separate pull requests so this one can be for adding fmpz_mpoly but changes for flint 3 would go somewhere else (actually we already have a gh-43 for that).

Any changes to the build system should definitely be kept separate from anything else so that we can see clearly the impact of each change independently.

[deinst]

Thanks for the info on a coefficient/exponent dict.
Sorry about the build system confusion (I forgot that I did that I'll see about backing it out)
I agree about a common base class (one exists, but it is currently empty)

[deinst]

There are actually three unrelated changes, include changes for flint 3.0, build system changes, and fmpz_changes. Sorry, I'll get to it.

[deinst]

Ok, I've backed out the flint-3 and setuptools and will see if that makes the CI happy. I am building a copy of flint-2.9 so I can test it. What is the version of arb that corresponds to flint-2.9?

[oscarbenjamin]

Looks like CI is working again.

ValueError would be when you recognise the type but the operation is not valid e.g. fmpz_mpoly with different contexts.

Not that it needs to be fixed here but actually it would be better to have more of an exception hierarchy like

class FlintError(Exception):
    pass

class FlintDomainError(FlintError):
    """Raised when the domains of operands are not compatible."""
    pass

class FlintIncompatibleGeneratorsError(FlintDomainError):
    """Raised when polynomial rings have incompatible generators."""
    pass

or something like that.

Generally it is better to have something more specific than ValueError and especially so if we think that anyone might ever want to catch the exception.

[oscarbenjamin]

@GiacomoPope will the changes in #61 impact on this PR?

I guess maybe ctx/thectx might need to be imported or used differently but it doesn't look like it's being used here.

It is probably worth being aware of the work @deinst is doing here while refactoring things to try to reduce potential conflicts.

[GiacomoPope]

Yeah the change should be minimal.

I think the only issues that will arise are when this code is using something which is previously global and will need to be imported to get everything working.

[oscarbenjamin]

What is the status of this now. Should we merge it and make further improvements later?

There are many more functions defined in fmpz_mpoly.h but we don't need to add them all in one PR. The functionality that is here should be tested though.

Here we have ctx.fmpz_mpoly_from_dict but there is no corresponding to_dict function. Probably we need a primitive for converting mpoly into other things.

[oscarbenjamin]

Probably the number one missing feature here that I would want for mpoly is factor.

[oscarbenjamin]

At some point when some version of this is merged I will try to integrate it with SymPy which and will probably discover a large number of methods that need to be added. For now the most important thing missing to be able to do that is just a to_dict method/function: with that I can at least implement temporary fallback methods for anything SymPy needs that is not yet provided by fmpz_mpoly.

[deinst]

This is not yet done, I need to touch up the call code, and I will add the factorization code. I got distracted by cython wierdness when attempting fmpq_mpoly, which I think is caused by a combination of the giant single module, and the trick of cythoning pyflint to make a _flint submodule, thereby dragging the _.flint pxd into view.

I have a mostly modularized codebase in my fork in the mod_experiment branch (I still need to extract the dirichlet module). It is not ready for a PR, as it includes a couple of poor design decisions, the most glaring of which is that I renamed all the submodules to having a prefixed underscore to avoid conflicts between submodule names and class names, when I should have them in a subdirectory/module like flint.flint. The result of this experiment has convinced me that splitting the code into modules is worth doing, and is worth doing before we expand the code further.

[oscarbenjamin]

The result of this experiment has convinced me that splitting the code into modules is worth doing, and is worth doing before we expand the code further.

Fair enough. A start was made in #61 and #63 and there was more discussion in #15 about where to put the files.

[oscarbenjamin]

We need a more general decision on how to handle a / b in non-field rings.

Currently this is a TypeError in all other cases:

In [2]: fmpz(2) / fmpz(2)
...
TypeError: unsupported operand type(s) for /: 'flint.types.fmpz.fmpz' and 'flint.types.fmpz.fmpz'

Likewise all polynomials and matrices. Division by a scalar is allowed where the base ring is a field e.g.:

In [4]: fmpz_poly([1, 0]) / fmpz(1)
...
TypeError: unsupported operand type(s) for /: 'flint.types.fmpz_poly.fmpz_poly' and 'flint.types.fmpz.fmpz'

In [5]: fmpq_poly([1, 0]) / fmpq(1)
Out[5]: 1

There is one exception:

In [8]: fmpz_mat([[1], [0]]) / fmpz(2)
Out[8]: 
[1/2]
[  0]

In other words scalar division of an fmpz_mat by an fmpz returns an fmpq_mat.

There is one other case of implicit conversion to a larger domain:

In [15]: fmpz_mat([[2, 0], [0, 2]]).inv()
Out[15]: 
[1/2,   0]
[  0, 1/2]

[deinst]

Ok, I think you're probably right. The gains in convenience are probably not worth the loss of sanity.

[oscarbenjamin]

I'm not sure what is best here in general. At least in the places where SymPy would use this the code is very careful not to use a / b if it is not known that we are working in a field. It would look like:

In [1]: from sympy import *

In [2]: x, y = symbols('x, y')

In [3]: R = QQ[x,y]

In [4]: R
Out[4]: QQ[x,y]

In [5]: p1 = R(x - y)

In [6]: p2 = R(x + y)

In [7]: K = R.get_field()

In [8]: K
Out[8]: QQ(x,y)

In [9]: K.convert_from(p1, R) / K.convert_from(p2, R)
Out[9]: (x - y)/(x + y)

Currently lacking in python-flint to make this work is:

• A systematic way to get the field of fractions from a ring.
• A systematic way to convert from the ring to that field.
• A systematic way to get the numerator and denominator back from the field as elements of the ring of integers.

It is not an immediate problem for SymPy that these things are lacking but I think that for someone using python-flint directly there should be some way to do this.

Probably what is needed is to have ring objects for every domain (including those that don't have a context in Flint like fmpz, fmpq, nmod). Probably this is all something that would be easier with generic rings so I am unsure right now what is the best approach to take in python-flint for conversions.

It might be good to have a meeting at some point to talk about some of these basic design questions.

[fredrik-johansson]

IMO the / operator should work as div in the new generics system: preserve the type, allow divisions when the quotient exists in the ring, raise an exception otherwise.

[oscarbenjamin]

preserve the type, allow divisions when the quotient exists in the ring, raise an exception otherwise.

That seems reasonable to me. Also if it is something that is used consistently in generic rings then it would be good to match that here and to use it consistently also across python-flint.

In SymPy this kind of division is called exquo:

In [10]: R.exquo(p1, p2)
...
ExactQuotientFailed: x + y does not divide x - y in QQ[x,y]

In [11]: R.exquo(p1, p1)
Out[11]: 1

This is a useful operation that is used in e.g. fraction-free algorithms. Exposing it as a / b would give SymPy a direct way to access the operation from python-flint rather than needing to use divmod.

It is unfortunate that SymPy's domain types cannot be trusted to use a / b for exquo but it all starts from using int or gmpy2.mpz for ZZ. Having this for fmpz leaves it still a little incompatible with int and gmpy2.mpz but that is unavoidable unless we want to match what they do:

In [12]: 2/3
Out[12]: 0.6666666666666666

I think it is fine for fmpz to be incompatible with this because gmpy2.mpz already exists for someone who wants a dropin faster replacement for int.

[deinst]

A meeting would be great. This week is not great as I have people fixing my roof, but after that anytime after 1200 GMT after that would be great.

I also have type related questions in reference to the call method of the various mpoly types: should I insist that all the parameters have identical types, or that at least one of the parameters is of a type that all the others can be promoted to (for example if one was an arb_poly then the other parameters could be arb, fmpz or fmpz_poly), or that there was some type, not necessarily occuring in the list, that all the parameters could be converted to.

I do think that we should probably find a way of inserting the flint generics under python-flint. Issue #86 addresses this as well.

[oscarbenjamin]

@GiacomoPope or @fredrik-johansson would you like to join a meeting with me and @deinst (possibly next week some time)?

[fredrik-johansson]

Sure.

[GiacomoPope]

Yeah that sounds good. Happy to make time next week. I'm in the UK re: timezones and planning but happy to be flexible :)

[oscarbenjamin]

Great, I'll follow up by email.

I will just note here some of the questions that I would like to resolve or at least make progress on:

• How exactly python-flint should use Flint's generics and also how to handle things that are not really generic.
• What a higher level python-flint interface might look like.
• How to handle explicit conversions between different types/domains.
• What to do with operations that do not necessarily make sense in a given domain like division and inverse.

[oscarbenjamin]

@deinst do you have any update on this?

I could take over soon if you are not going to work on this right now.

[deinst]

@deinst do you have any update on this?

I could take over soon if you are not going to work on this right now.

Sorry, I got distracted by life. I'll get it done this weekend.

[Jake-Moss]

Hi @deinst , I'm looking at picking up this PR as part of some contributions I'm planning on making to SymPy (mailing list thread here). Just wanted to let you know I'm planning on continuing your work here in the coming weeks. If you have any comments or objections please let me know.