Ask Your Question
2

Issue with inversion of complex symbolic array initialized with numpy

asked 2021-04-28 13:40:29 +0100

yorgos_sot gravatar image

Hello!

I am trying to invert a symbolic matrix that is initialized with a combination of numpy arrays and I get an error. Below i present a simple code that gives the error.

import numpy as np

a = np.zeros((5,5) , dtype = 'complex')
np.fill_diagonal(a,1)

b = var('x')*a

c = matrix(b)
c.inverse()

Gives the error:

ECL says: THROW: The catch MACSYMA-QUIT is undefined.

What I have noticed up until now is that the issue stems from the data type. When I try casting the numpy array to float before turning it into a matrix, it works. However, my actual code makes use of complex coefficients. I believe it might have to do with how the imaginary part is represented in SageMath in comparison to numpy. ( j vs I )

When I manually create the symbolic array with I for the imaginary part, the .inverse() has no issue.

edit retag flag offensive close merge delete

1 Answer

Sort by » oldest newest most voted
2

answered 2021-04-28 14:45:50 +0100

slelievre gravatar image

Thanks for reporting.

One alternative approach would be, after creating the array a, instead of multiplying by x and then converting to a matrix, to first convert to a matrix and then multiply by x.

Create the array a:

sage: import numpy as np
sage: a = np.zeros((5, 5), dtype='complex')
sage: np.fill_diagonal(a, 1)
sage: x = SR.var('x')
sage: b = x*a
sage: a
array([[1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j],
       [0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j],
       [0.+0.j, 0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j],
       [0.+0.j, 0.+0.j, 0.+0.j, 1.+0.j, 0.+0.j],
       [0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 1.+0.j]])

From the question: multiply array by x, then convert to matrix:

sage: b
array([[((1+0j))*x, 0j, 0j, 0j, 0j],
       [0j, ((1+0j))*x, 0j, 0j, 0j],
       [0j, 0j, ((1+0j))*x, 0j, 0j],
       [0j, 0j, 0j, ((1+0j))*x, 0j],
       [0j, 0j, 0j, 0j, ((1+0j))*x]], dtype=object)
sage: c = matrix(b)
sage: c
[((1+0j))*x         0j         0j         0j         0j]
[        0j ((1+0j))*x         0j         0j         0j]
[        0j         0j ((1+0j))*x         0j         0j]
[        0j         0j         0j ((1+0j))*x         0j]
[        0j         0j         0j         0j ((1+0j))*x]
sage: c.parent()
Full MatrixSpace of 5 by 5 dense matrices over Symbolic Ring

sage: c.inverse()
Traceback (most recent call last)
...
TypeError: ECL says: THROW: The catch MACSYMA-QUIT is undefined.

Alternative approach: convert to matrix, then multiply by x:

sage: b = matrix(a)
sage: b
[1.0 0.0 0.0 0.0 0.0]
[0.0 1.0 0.0 0.0 0.0]
[0.0 0.0 1.0 0.0 0.0]
[0.0 0.0 0.0 1.0 0.0]
[0.0 0.0 0.0 0.0 1.0]
sage: b.parent()
Full MatrixSpace of 5 by 5 dense matrices over Complex Double Field

sage: c = x * b
sage: c
[1.0*x   0.0   0.0   0.0   0.0]
[  0.0 1.0*x   0.0   0.0   0.0]
[  0.0   0.0 1.0*x   0.0   0.0]
[  0.0   0.0   0.0 1.0*x   0.0]
[  0.0   0.0   0.0   0.0 1.0*x]
sage: c.parent()
Full MatrixSpace of 5 by 5 dense matrices over Symbolic Ring

sage: c.inverse()
[1.0/x     0     0     0     0]
[    0 1.0/x     0     0     0]
[    0     0 1.0/x     0     0]
[    0     0     0 1.0/x     0]
[    0     0     0     0 1.0/x]
edit flag offensive delete link more

Comments

Thank you so much. I would like to also add something that I had to modify in my code.

Because I use python, I had to also define exponential variables as:

import sage.all as sg
e_Iphi = sg.exp(sg.I * sg.var('phi', domain='real') )
yorgos_sot gravatar imageyorgos_sot ( 2021-04-28 17:37:06 +0100 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

Stats

Asked: 2021-04-28 13:31:54 +0100

Seen: 395 times

Last updated: Apr 28 '21