Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

The module Matrices in SymPy doesn't support type mpq from gmpy2. #16883

Open
adukova opened this issue May 24, 2019 · 12 comments
Open

The module Matrices in SymPy doesn't support type mpq from gmpy2. #16883

adukova opened this issue May 24, 2019 · 12 comments
Labels

Comments

@adukova
Copy link

adukova commented May 24, 2019

Example

In[1]:   from sympy import Matrix,S,symbols,QQ,I,init_printing,Poly
            from sympy.polys.domains import QQ
            import gmpy2  as gm
            init_printing()
In[2]:    k2=QQ(1,4)
            k1=QQ(-3,1)
            k0=QQ(1345789308645,512855373753735374)
            type(k2),type(k1),type(k0)
Out[2]: (mpq, mpq, mpq)
In[3]:    A=Matrix([[k0,k1],[k2,k0]])
In[4]:    type(A[0,0]);type(A[0,1]),type(A[1,0]),type(A[1,1])
Out[4]: (sympy.core.numbers.Integer,
             sympy.core.numbers.Rational,
             sympy.core.numbers.Rational)

It does not allow to work with long arithmetic in linear algebra.

@jksuom
Copy link
Member

jksuom commented May 24, 2019

Matrix elements are converted to SymPy types by default. There is a variant class PolyMatrix whose instances do not convert the input data. That can be used, for instance, by from sympy.polys.polymatrix import PolyMatrix as Matrix.

@adukova
Copy link
Author

adukova commented May 24, 2019

@jksuom Thank you very much. I'll try to do it.

@adukova
Copy link
Author

adukova commented May 25, 2019

@jksuom
Thank you for your advice. Now matrix elements have type mpq. Unfortunately this does not solve my problems. From the module Matrices I use only rank and nullspace. The function rank works as long as for type Rational. May be it is needed to change the function rref in order to it can work with gmpy2. Do you have any thoughts on this? Thank you in advance.

@sylee957
Copy link
Member

sylee957 commented May 25, 2019

@adukova
Is the rank or nullspace giving the wrong answer?
I have got the result from the matrix above

>>> A.nullspace()
[]
>>> A.rank()
2

And I think Mathematica gives the same answer.

@jksuom
Copy link
Member

jksuom commented May 25, 2019

Can you describe the problems you meet?

@adukova
Copy link
Author

adukova commented May 25, 2019

I have matrices of sizes about 50 x 50. The entries are rational numbers with long numerators and denominators (about 50 digits and more). I need ranks for a finite sequence from such matrices (about 50 terms). Maple solves this problem in reasonable time (from 1 sec to 1 hour). The analogous routine in SymPy for some example can not solve the problem in five hours. I think the problem is that the module Matrices can not work with the library gmpy2. Moreover if entries were complex numbers the function rank in SymPy returned a wrong result. But this problem is solved now (see issue #16823).

@adukova
Copy link
Author

adukova commented May 25, 2019

Thank you all. It seems in the case of real matrices the problem is solved:

  1. with the help of PolyMatrix we create a matrix with mpq entries (thank you jksuom),
  2. the function rank() must be used with parameters r = M.rank(iszerofunc=lambda v: v==QQ.zero).
    Now the problem is solved with very high speed. It remains to check the case of matrices with complex elements. It is not clear what is considered as zero.

@sylee957
Copy link
Member

I have found the behavior that the matrix gets gradually filled with SymPy number types in each step of computation. So the booster may not last long.
I think that sympify automatically sanitizes the custom number types in the poly module to sympy number types. Is it an intended behavior?

@jksuom
Copy link
Member

jksuom commented May 28, 2019

sympify automatically sanitizes the custom number types in the poly module to sympy number types. Is it an intended behavior?

That should not happen. The existing variant matrix subclasses (RawMatrix, PolyMatrix) replace _sympify with identity mapping to avoid that. It seems to me that the best choice would be conversion to the ring data type. (I am working on that.)

oscarbenjamin pushed a commit that referenced this issue May 28, 2019
Matrix codebase currently makes use of the hardcoded constants
S.Zero and S.One. That is inconvenient when a non-Basic data
type is used to speed up computation (cf. #16883). This PR
replaces these constants with class attributes zero and one that
can be modified by subclassing.
@adukova
Copy link
Author

adukova commented Jun 9, 2019

@jksuom

Dear jksuom,

Unfortunately, I know you only as jksuom in the forum GitHub .

Thank you very much for your advice to my daughter Adukova Nataliya (GitHub issue # 16883).
Your advice help me to improve the work of my SymPy module for the Wiener-Hopf factorization
of matrix polynomials.

Can I thank you in my article for Journal of Symbolic Computation?
What your name can I note in my "Acknowledgments"?

Best wishes,

V. Adukov

@jksuom
Copy link
Member

jksuom commented Jun 9, 2019

Can I thank you in my article for Journal of Symbolic Computation?

Yes, of course, you can do that.

What your name can I note in my "Acknowledgments"?

My name is Kalevi Suominen.

@adukova
Copy link
Author

adukova commented Jun 9, 2019

@jksuom
Thank you

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

No branches or pull requests

3 participants