A simple idea regarding groundtypes for Matrix

6 views
Skip to first unread message

SherjilOzair

unread,
May 29, 2011, 3:04:26 PM5/29/11
to sympy
Why didn't anyone think of this before ?!

Q : What benefits does groundtypes give ?
A : Operations are faster, as one does not have to type-check. Each
groundtype has its own __add__, __mul__ installed, and global Add, Mul
are avoided. The Matrix doesn't have to worry about the zero-
equivalence. x belong to type MYTYPE for example, the question x == 0,
should be solved by the groundtype and not by the Matrix class.

Q : So, all the elements have to be of the same groundtype ? How do
you plan to do it ?
A : Simple. Instead of sympify-ing the vals before putting it in the
matrix, I typify them, typify being a variable function, which the
user/caller will provide, like this A = Matrix(..., typify = sympify)
or A = Matrix(..., typify = poly). The Matrix constructor would then
use A[i,j] = typify(val).

Q : So the matrix just needs to have a typify function ?
A : Well, It could have another variable which says "This Matrix
belongs to the Poly domain !". But Its just for the user. The matrix
doesn't need to know which domain it is in (there are exceptions, keep
reading). The Matrix shouldn't care what A[6,5] and A[7,8] types are.
When they need to be multiplied or added, A[6,5] * A[7,8] will
automatically the right __mul__ and similarly for addition.

Q : That easy ?! :O
A : No, this can only work for 'domain-preserving algorithms'. An
algorithm which for example employs the division or the square root,
shouldn't/wouldn't operate on Polys. It should return an Error,
complaining about the bad groundtype.
OR
We could have type-conversion. For division-employing algorithms, we
could convert the groundtype to RationalFunction, for squareroot-
employing algorithms, we could convert to the Expr.

Q : So the user has to be smart enough to know why typify function to
use for what ?
A : Maybe or maybe not. We can always have automatic identification
added in later.

Q : What do you need from Polys for this ?
A : Basically, I need an exhaustive set of groundtypes which in
essence have the fundamental operators coded along with __eq__ for the
zero-equivalence check. If someone (mattpap ?) could help me compile a
list of groundtypes which I can extract from Polys with description of
each and what operations it supports and, equally important, what it
does not suuport.

I may have used the terms 'ground type' to mean both 'ground domain'
and 'ground type'. But I think my meaning is clear.

This is just a simple idea that came in my mind which made me think
'why not ?'. Did I leave out/forget anything ? What does this idea
lack ? Am I completely wrong ?

I think, this idea is good enough to iterate over and build on, and
would love to see the communities' modifications and suggestions.

Cheers,
Sherjil Ozair

Tom Bachmann

unread,
May 29, 2011, 3:13:06 PM5/29/11
to sy...@googlegroups.com
How is this at all different from what Polys does? I'm not saying it's
bad, I'm just not seeing your point. Basically what you call ground
types are called domains in polys, and they are in polys/domains ...

As far as I understand, Poly is mostly a smart wrapper that figures out
what domains to use, and translates highlevel operations on polynomials
into low-level operations.

So the Matrix class should not /use/ Poly instances, but parallel Poly
in the sense that it is a container, figuring out what domains to use,
and dispatching operation to low-level implementations.

Again, I should have prefixed this by saying that Mateusz is the expert
and I'm just using polys a bit...

SherjilOzair

unread,
May 29, 2011, 3:22:26 PM5/29/11
to sympy
I see. I must confess I'm not deeply familiar with Polys, and also
haven't even used it much yet.
This idea is infact very similar to Polys. Its just that I've been
able to figure out a translation to the Matrix code, thus my
excitement.

On May 30, 12:13 am, Tom Bachmann <e_mc...@web.de> wrote:
> How is this at all different from what Polys does? I'm not saying it's
> bad, I'm just not seeing your point. Basically what you call ground
> types are called domains in polys, and they are in polys/domains ...

I need usable types. For example, you mentioned that one usable type
is DMF, and it is in poly/polyclasses.
How many other such types exist and where ?

>
> As far as I understand, Poly is mostly a smart wrapper that figures out
> what domains to use, and translates highlevel operations on polynomials
> into low-level operations.
>
> So the Matrix class should not /use/ Poly instances, but parallel Poly
> in the sense that it is a container, figuring out what domains to use,
> and dispatching operation to low-level implementations.

I see your point.

Tom Bachmann

unread,
May 29, 2011, 5:10:55 PM5/29/11
to sy...@googlegroups.com
> On May 30, 12:13 am, Tom Bachmann<e_mc...@web.de> wrote:
>> How is this at all different from what Polys does? I'm not saying it's
>> bad, I'm just not seeing your point. Basically what you call ground
>> types are called domains in polys, and they are in polys/domains ...
>
> I need usable types. For example, you mentioned that one usable type
> is DMF, and it is in poly/polyclasses.
> How many other such types exist and where ?
>

DMF is already higher level. The things you should probably be looking
at are in polys/domains. There's a lot there, and it's a bit of a shame
that there's not much documentation. I think construct_domain should
somehow construct a domain that can hold your data. If you need to
divide, or take square roots, etc, you presumably have to figure out
what larger domain you need. There are probably functions for this, but
I don't know about them. Ask Mateusz about anything specific.

As for ground types, they work roughly like this:

>>> from sympy.polys.domains import FF

This imports a constructor for finite fields. As I understand it, this
will automatically use gmpy or python types depending on what is available.


Construct a domain for arithmetic mod 5:
>>> F5 = FF(5)

Do some arithmetic:
>>> F5(3) + F5(8)
1 mod 5
>>> F5(3)/F5(2)
4 mod 5

Funny things can happen if you do division where you should not:
>>> from sympy.polys.domains import ZZ
>>> ZZ(2)/ZZ(3)
0.666666666667

But this works:
>>> Q = ZZ.get_field()
>>> Q(2)/Q(3)
2/3

I suppose to see how to do this sort of stuff automatically, you have to
read polytools.py.

Mateusz Paprocki

unread,
May 29, 2011, 5:44:16 PM5/29/11
to sy...@googlegroups.com
Hi,

On 29 May 2011 23:10, Tom Bachmann <e_m...@web.de> wrote:
On May 30, 12:13 am, Tom Bachmann<e_mc...@web.de>  wrote:
How is this at all different from what Polys does? I'm not saying it's
bad, I'm just not seeing your point. Basically what you call ground
types are called domains in polys, and they are in polys/domains ...

I need usable types. For example, you mentioned that one usable type
is DMF, and it is in poly/polyclasses.
How many other such types exist and where ?


DMF is already higher level. The things you should probably be looking at are in polys/domains.

Well, DMF is low-level. In domains you will find FractionField domain that uses DMF a ground type.
 
There's a lot there, and it's a bit of a shame that there's not much documentation.

There isn't that much, but main ideas are described in my thesis (ch. 2). If more information is needed, I can always provide it (as long as a question is specific).
 
I think construct_domain should somehow construct a domain that can hold your data. If you need to divide, or take square roots, etc, you presumably have to figure out what larger domain you need. There are probably functions for this, but I don't know about them. Ask Mateusz about anything specific.

As for ground types, they work roughly like this:

>>> from sympy.polys.domains import FF

This imports a constructor for finite fields. As I understand it, this will automatically use gmpy or python types depending on what is available.


Construct a domain for arithmetic mod 5:
>>> F5 = FF(5)

Do some arithmetic:
>>> F5(3) + F5(8)
1 mod 5
>>> F5(3)/F5(2)
4 mod 5

Funny things can happen if you do division where you should not:
>>> from sympy.polys.domains import ZZ
>>> ZZ(2)/ZZ(3)
0.666666666667

But this works:
>>> Q = ZZ.get_field()
>>> Q(2)/Q(3)
2/3

Division in Python is tricky because there are / and //. The meaning of / can differ (2/3 can give either 0 or a float). If you want to make your code independent of this use quo(), rem(), div() methods of an appropriate domain, e.g.:

In [1]: ZZ.quo(ZZ(2), ZZ(3))
Out[1]: 0

In [2]: ZZ.rem(ZZ(2), ZZ(3))
Out[2]: 2

In [3]: ZZ.div(ZZ(2), ZZ(3))
Out[3]: (0, 2)

Domains and ground types in polys work as follows: +, - (unary/binary) and * must be implemented in ground types and must solve zero equivalence problem. Division, gcd, lcm and other can be implemented in ground types but code using those ground types can't take advantage of them and must use domain.div(), domain.gcd(), domain.lcm() and so on.

Using +, -, * in code directly makes the code fast. Division isn't very common so we can use a thin layer above / and //, and as division is generally a mess (e.g. extra qdiv() in gmpy), the extra layer in necessary to use single code for multiple ground types (here I'm not only speaking about handling numbers and polynomials in one code, but also different implementations of numbers (Integer, int, mpz, ...), because they have very different APIs).


I suppose to see how to do this sort of stuff automatically, you have to read polytools.py.


--
You received this message because you are subscribed to the Google Groups "sympy" group.
To post to this group, send email to sy...@googlegroups.com.
To unsubscribe from this group, send email to sympy+un...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/sympy?hl=en.


Mateusz

Aaron Meurer

unread,
Jun 1, 2011, 9:31:26 PM6/1/11
to sy...@googlegroups.com
I think the key point that should be made here is that the matrices
should try to use the *exact* same ground types/domains as the polys,
so that there is no code duplication. This will probably involve
improving and even changing a little bit the design of the polys
classes (but this is good, because it will make them more modular,
which is the goal anyway).

On Sun, May 29, 2011 at 3:44 PM, Mateusz Paprocki <mat...@gmail.com> wrote:
> Hi,
>
> On 29 May 2011 23:10, Tom Bachmann <e_m...@web.de> wrote:
>>>
>>> On May 30, 12:13 am, Tom Bachmann<e_mc...@web.de>  wrote:
>>>>
>>>> How is this at all different from what Polys does? I'm not saying it's
>>>> bad, I'm just not seeing your point. Basically what you call ground
>>>> types are called domains in polys, and they are in polys/domains ...
>>>
>>> I need usable types. For example, you mentioned that one usable type
>>> is DMF, and it is in poly/polyclasses.
>>> How many other such types exist and where ?
>>>
>>
>> DMF is already higher level. The things you should probably be looking at
>> are in polys/domains.
>
> Well, DMF is low-level. In domains you will find FractionField domain that
> uses DMF a ground type.
>
>>
>> There's a lot there, and it's a bit of a shame that there's not much
>> documentation.
>
> There isn't that much, but main ideas are described in my thesis (ch. 2). If
> more information is needed, I can always provide it (as long as a question
> is specific).

There would be more, as I wrote doctests for all those things back in
https://github.com/asmeurer/sympy/tree/polydocs. But the domains code
was all moved to polys/domains and the API changed so much that I
never bothered to transfer the doctests.

If you need to know how a module works and there is no documentation,
writing doctests for all the functions/methods is a great way to fix
both problems. This is how I learned how the polys in general worked
at the beginning of last summer (see the above linked branch). In
this case, you could just work on transferring my doctests that were
never transferred.

Aaron Meurer

Brian Granger

unread,
Jun 1, 2011, 10:09:18 PM6/1/11
to sy...@googlegroups.com
Is there a link somewhere to documentation about the groundtypes in
polys? I am interested why such things are needed? It sort of sounds
like there are two semi-independent code bases living here...

--
Brian E. Granger
Cal Poly State University, San Luis Obispo
bgra...@calpoly.edu and elli...@gmail.com

SherjilOzair

unread,
Jun 2, 2011, 12:27:22 AM6/2/11
to sympy
scipy.sparse implements a dtype kwarg argument, but which currently
cannot take in arbitrary unknown types though.
One thing that can be done about this, is to define an interface for
the dtype. It would be taken for granted that it will have +, *, /
defined. Checks will be used in algorithms if it has pow, inverse
defined or not.
The caller will be provided with a dtype function argument in the
Matrix constructor.

I list some built-in dtypes that sympy has and can provide. With only
these 6 dtypes, Matrix should suffice for 90% of symbolic matrix
needs. Possibly, if the caller doesn't want any of these dtypes, then
he should specify 'other' in the dtype argument.

Int, numeric, can employ addition, multiplication, raising to positive
integral power.
Rational, numeric, can employ addition, multiplication, division,
inverse, raising to integral power.
Real, numeric, can employ addition, multiplication, division, inverse,
raising to any power.
Poly (Or one of its internals), symbolic, to support addition,
multiplication, division by scalar, *not* inverse, raising to positive
integral power.
Rational Function, symbolic, to support addition, multiplication,
division by scalar, inverse, raising to integral power.
Expr, symbolic, to support addition, multiplication, division by
scalar, inverse, raising to all powers.

SherjilOzair

unread,
Jun 2, 2011, 12:28:58 AM6/2/11
to sympy, mat...@gmail.com
@mattpap, are you writing RationalFunction ?
If not, I'll be glad to. Just tell me how should I go about it ? Which
interface to follow and all.

SherjilOzair

unread,
Jun 2, 2011, 12:34:20 AM6/2/11
to sympy
I need review on this style for matrix groundtypes, i.e. having a list
of builtin dtypes, and an other option for misc backends.
Do you think the list of builtins is good enough ? Should I add/remove
stuff ?

The Matrix will still support seamlessness. I mean, if a Matrix of
Ints comes up, and user wants to take its inverse, then the dtype will
automatically be changed to Rational, and the inverse calculated in
the Rational field. The original Matrix however would remain the same.
( In accordance to immutability).

SherjilOzair

unread,
Jun 2, 2011, 3:58:22 AM6/2/11
to sympy

Mateusz Paprocki

unread,
Jun 2, 2011, 4:28:52 AM6/2/11
to sy...@googlegroups.com
Hi,

On 2 June 2011 04:09, Brian Granger <elli...@gmail.com> wrote:
Is there a link somewhere to documentation about the groundtypes in
polys?  I am interested why such things are needed?  It sort of sounds
like there are two semi-independent code bases living here...

See my previous response. Unfortunately, currently there isn't much documentation on this topic. As to two code bases, there are actually two cores, one (should be) purely symbolic (sympy.core) and the other algebraic (sympy.polys and sympy.polys.domains, later should be extracted to a separate top-level module, to make it reusable in other parts of SymPy).

Mateusz

Mateusz Paprocki

unread,
Jun 2, 2011, 4:47:09 AM6/2/11
to sy...@googlegroups.com
Hi,

On 2 June 2011 06:27, SherjilOzair <sherji...@gmail.com> wrote:
scipy.sparse implements a dtype kwarg argument, but which currently
cannot take in arbitrary unknown types though.

Comparing two entities from two very different libraries may be misleading, because SciPy has different needs and SymPy has different.
I'm not sure where you are going with this. In sympy.polys (and specifically domains submodule) there is everything (well, almost) that you need to start working on matrices. We ask you not to reinvent the wheel, but to use the framework that is already in SymPy, report bugs and issues, and, if necessary, improve it.

In another thread I asked you to start with something simple, mainly with making solve_linear_system() aware of sympy.polys.domains. This way you should be able to quickly learn how to interact with domains and ground types, and what are the possible issues, and why this is needed at all. In some very old local branch I did this exercise my self. Until this is done (I mean until you use the framework for something practical), all of the discussions will be purely theoretical.


--
You received this message because you are subscribed to the Google Groups "sympy" group.
To post to this group, send email to sy...@googlegroups.com.
To unsubscribe from this group, send email to sympy+un...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/sympy?hl=en.


Mateusz

Brian Granger

unread,
Jun 2, 2011, 11:48:29 AM6/2/11
to sy...@googlegroups.com
On Thu, Jun 2, 2011 at 1:28 AM, Mateusz Paprocki <mat...@gmail.com> wrote:
> Hi,
>
> On 2 June 2011 04:09, Brian Granger <elli...@gmail.com> wrote:
>>
>> Is there a link somewhere to documentation about the groundtypes in
>> polys?  I am interested why such things are needed?  It sort of sounds
>> like there are two semi-independent code bases living here...
>
> See my previous response. Unfortunately, currently there isn't much
> documentation on this topic. As to two code bases, there are actually two
> cores, one (should be) purely symbolic (sympy.core) and the other algebraic
> (sympy.polys and sympy.polys.domains, later should be extracted to a
> separate top-level module, to make it reusable in other parts of SymPy).

Then I don't understand what polys has to do with Matrices. Surely
matrices can contain expressions other than polynomials. This is
really what I am wondering about. Are we talking about a special
Matrix subclass that *only* contains polynomials?

Brian Granger

unread,
Jun 2, 2011, 11:56:29 AM6/2/11
to sy...@googlegroups.com
I am a bit confused here as well.

Are you considering adding a dtype argument to sympy.Matrix?
Are you consider making sympy.Expr objects work inside numpy matrices?
Are you considering making sympy.Expr work inside scipy.sparse matrices?

Cheers,

Brian

sherji...@gmail.com

unread,
Jun 2, 2011, 12:13:42 PM6/2/11
to sy...@googlegroups.com
Yes and no and no.
Numpy/scipy matrices have their backend in C/fortran types. Arbitrary objects can't be put inside them.

That's where sympy comes in and finds a market for its use. I'm adding a dtype(the usgae of this name might be wrong) argument to sympy matrices. 6 basic types, ints, rationals, reals, polys, rational functions, and finally Exprs will be built in. It will have an 'other' option which will the user to put in any arbitrary type as elements of the matrix which support the fundamental operations required for matrix algorithms.

Possibly, a template could be provided to the user to have provide the sympy matrix with domain.sum, domain.typify, etc..

Still, 90% of symbolic matrix needs would be covered by the 6 builtins.
Sent on my BlackBerry® from Vodafone

Aaron S. Meurer

unread,
Jun 2, 2011, 9:41:14 PM6/2/11
to sy...@googlegroups.com

On Jun 2, 2011, at 9:48 AM, Brian Granger wrote:

> On Thu, Jun 2, 2011 at 1:28 AM, Mateusz Paprocki <mat...@gmail.com> wrote:
>> Hi,
>>
>> On 2 June 2011 04:09, Brian Granger <elli...@gmail.com> wrote:
>>>
>>> Is there a link somewhere to documentation about the groundtypes in
>>> polys? I am interested why such things are needed? It sort of sounds
>>> like there are two semi-independent code bases living here...
>>
>> See my previous response. Unfortunately, currently there isn't much
>> documentation on this topic. As to two code bases, there are actually two
>> cores, one (should be) purely symbolic (sympy.core) and the other algebraic
>> (sympy.polys and sympy.polys.domains, later should be extracted to a
>> separate top-level module, to make it reusable in other parts of SymPy).
>
> Then I don't understand what polys has to do with Matrices. Surely
> matrices can contain expressions other than polynomials. This is
> really what I am wondering about. Are we talking about a special
> Matrix subclass that *only* contains polynomials?

It isn't related to the polys in the sense the that matrix elements are polynomials. Rather, think of the coefficients of a polynomial. Poly(2*x**3 + x + 4) is stored internally as [2, 0, 1, 4]. The problems that arise when working with coefficients of polynomials are very similar to the ones that arise when working with matrices. Mateusz has modularized the coefficients of Poly pretty well so that there are different classes for integer coefficients or rational coefficients or even rational function coefficients (like Poly((1 + y)/y**2*x, x)). The integer and rational number coefficient domains are written so that they can use gmpy if it's installed, or else it uses Python, and all automatically. This should all be adapted to the matrices, because it makes the code cleaner and, more importantly, it makes things much faster.

Aaron Meurer

Aaron S. Meurer

unread,
Jun 2, 2011, 9:41:45 PM6/2/11
to sy...@googlegroups.com

On Jun 2, 2011, at 10:13 AM, sherji...@gmail.com wrote:

> Yes and no and no.
> Numpy/scipy matrices have their backend in C/fortran types. Arbitrary objects can't be put inside them.
>
> That's where sympy comes in and finds a market for its use. I'm adding a dtype(the usgae of this name might be wrong) argument to sympy matrices. 6 basic types, ints, rationals, reals, polys, rational functions, and finally Exprs will be built in. It will have an 'other' option which will the user to put in any arbitrary type as elements of the matrix which support the fundamental operations required for matrix algorithms.
>
> Possibly, a template could be provided to the user to have provide the sympy matrix with domain.sum, domain.typify, etc..

I'd say this already exists as the super class of the poly domains.

Aaron Meurer

Reply all
Reply to author
Forward
0 new messages