Ask Your Question

How to speed up the code?

asked 2022-10-10 16:59:47 +0200

mak3521 gravatar image

Given a matrix M, to find all 1,-1 matrices which satisfy the condition given in the code. For n=3, I am getting the outputs easily. For n=6, it ran for around 12 hours and got "flint exception error". I understand that len(T) is huge so there might be some memory issues. Is there anyway to get over it?

 A=[matrix(n,n,t) for t in T]
 for a in A:
     if 3*a^2==M:
edit retag flag offensive close merge delete


Can you give an example of matrix M resulting in an error?

Max Alekseyev gravatar imageMax Alekseyev ( 2022-10-10 20:27:40 +0200 )edit

As for speeding up the code it may be worth to compute the square root of matrix M/3and check if it composed of $\pm1$.

Max Alekseyev gravatar imageMax Alekseyev ( 2022-10-10 20:34:06 +0200 )edit

An example of M is:

  M = -n *ones_matrix(n) + (4*n)*(identity_matrix(n))
mak3521 gravatar imagemak3521 ( 2022-10-11 04:47:15 +0200 )edit

2 Answers

Sort by » oldest newest most voted

answered 2022-10-10 20:34:35 +0200

updated 2022-10-11 20:15:35 +0200

A simple modification would be to not create the list of matrices:

 n = M.nrows()
 for t in Tuples((1,-1), n^2):
     a = matrix(n, n, t)
     if 3*a^2 == M:

With n=5, this at least enters the loop immediately, whereas the original code spends a while creating the list of matrices. I gave up waiting for the original code after about 5 minutes. (You can add some print statements to see its progress.)

By the way, you can see that Tuples is fast, but creating a list from its entries can be slow:

sage: %time Tuples((1, -1), 16)
CPU times: user 24 µs, sys: 8 µs, total: 32 µs
Wall time: 46 µs
Tuples of (1, -1) of length 16

sage: %time L = [t for t in Tuples((1, -1), 16)] # slow
CPU times: user 453 ms, sys: 23.6 ms, total: 477 ms
Wall time: 479 ms

sage: %time L = [matrix(4, 4, t) for t in Tuples((1, -1), 16)] # slower
CPU times: user 988 ms, sys: 15.6 ms, total: 1 s
Wall time: 1.01 s

Going through all of these matrices is going to take a long time, in any case:

sage: len(Tuples((1, -1), 25))
sage: len(Tuples((1, -1), 36))

Edit: here is a modified version of your code:

def f(M):
    results = []
    n = M.nrows()
    j = 0
    for t in Tuples((1,-1), n^2):
        a = matrix(n, n, t)
        if j % 10000 == 0:
        j += 1
        if 3*a^2 == M:
            # if you want to stop after finding one solution, you could add "break" here to stop the loop.

You can then call f(M) for any given matrix M. It prints something after every 10,000 matrices, and it saves the results and returns them at the end.

edit flag offensive delete link more


Ok, thanks for the information

mak3521 gravatar imagemak3521 ( 2022-10-11 04:51:42 +0200 )edit

In fact I would guess that your original error occurred in the creation of the list A: I don't know if Python will successfully create a list of 2**36 6-by-6 matrices.

John Palmieri gravatar imageJohn Palmieri ( 2022-10-11 20:09:37 +0200 )edit

I started my computer on the 6x6 case with the matrix you mentioned above (-n *ones_matrix(n) + (4*n)*(identity_matrix(n)), and at the rate it's going, it looks like it would take a few months.

John Palmieri gravatar imageJohn Palmieri ( 2022-10-11 23:34:23 +0200 )edit

Your idea behind "modified version" is really very good . Thanks . I also ran it on my laptop a few hours back. Yes, it looks it is going to take a lot of time.

mak3521 gravatar imagemak3521 ( 2022-10-12 06:12:18 +0200 )edit

answered 2022-10-11 22:13:23 +0200

Max Alekseyev gravatar image

updated 2022-10-12 15:05:29 +0200

One can employ matrix square root to significantly speed up things. The following implementation illustrates this idea and supports only matrices with distinct eigenvalues:

def sqrt_pm1(M):
    if not M.det().is_square():    # quick test for insolubility
        return []

    n = M.nrows()

    d, L = M.eigenmatrix_right()
    d = d.diagonal()
    assert len(set(d)) == n          # only this case is supported
    L = L.change_ring(QQbar)  

    pm1 = {QQbar(-1),QQbar(1)}
    return [R for signs in Tuples((-1,1),n) if set((R := L * diagonal_matrix([sqrt(e)*s for e,s in zip(d,signs)]) * L.inverse()).list()).issubset(pm1)]

For example,

sqrt_pm1(Matrix([[0,-6],[6,0]]) / 3)

gives the following list of two $\pm 1$-matrices that are square roots of Matrix([[0,-6],[6,0]]) / 3:

[-1  1]  [ 1 -1]
[-1 -1], [ 1  1]
edit flag offensive delete link more


Okay...Thanks a lot for sparing your time. For M = -n ones_matrix(n) + (4n)*(identity_matrix(n)), mentioned above , eigen values are not distinct. Can anything be done in that case?

mak3521 gravatar imagemak3521 ( 2022-10-12 06:14:57 +0200 )edit

At very least you can quickly filter out unsolvable M by checking (M/3).det().is_square()

Max Alekseyev gravatar imageMax Alekseyev ( 2022-10-12 12:08:53 +0200 )edit

Okay. Thanks

mak3521 gravatar imagemak3521 ( 2022-10-12 19:00:01 +0200 )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

1 follower


Asked: 2022-10-10 16:59:47 +0200

Seen: 353 times

Last updated: Oct 12 '22