Processing math: 11%
Ask Your Question
1

Unimodular matrices with additional restrictions

asked 5 years ago

Daniel L gravatar image

updated 5 years ago

I'd like to generate some unimodular matrices over ZZ with some restrictions:

  1. the entries are between -2 and 2
  2. the matrix is symmetric: A=At

Item (1) above can be addressed with the option upper_bound. Simply searching for matrices satisfying (1) and (2),

for j in [1..1000]:
A = random_matrix(ZZ,3,3, algorithm = 'unimodular', upper_bound = 3, max_tries = 1000)
    if A == A.transpose():
        A

becomes extremely inefficient, especially as the dimension grows. Is there a way to enforce items (1),(2) within the unimodular algorithm? Or another workaround?

Edit: I originally also wanted the diagonal entries fixed at 2, but this is a bit too restrictive.

Preview: (hide)

Comments

3

Are you sure that your restrictions are correct? For example, let A=(2aba2cbc2), with a,b,cZ. Since det A will be unimodular if and only if a^2+b^2+c^2-abc=\frac{7}{2}, which is impossible because a^2+b^2+c^2-abc\in\mathbb{Z}.

Juanjo gravatar imageJuanjo ( 5 years ago )

Good point. This certainly works if the dimension is 8 e.g. the E_8 root lattice gram matrix: E8 lattice wiki. I edited the post. Thank you

Daniel L gravatar imageDaniel L ( 5 years ago )

1 Answer

Sort by » oldest newest most voted
1

answered 5 years ago

Juanjo gravatar image

updated 5 years ago

I have adapted an idea I have found in Stackexchange. It is based on the following facts:

  • Let J_k be the identity matrix of order n with the null elements in the k-th row replaced by random integers. Then J_k is unimodular, i.e. \lvert J_k\rvert=1.

  • If A is an unimodular symmetric matrix, so it is J_k^T A J_k.

Thus, one can start with A=I_n, choose a row k at random, construct J_k, compute J_k^T A J_k and replace A by that matrix. This process can be repeated as many times as wanted. The result is an unimodular symmetric matrix. The following code implements this idea:

def unimod_sym_matrix(n, niter=1, min_int=-1, max_int=1):
    A = identity_matrix(n)
    rows = random_vector(niter, x=0, y=n)
    for k in rows:
        J = identity_matrix(n)
        J[k,:] = random_vector(n, x=min_int, y=max_int+1)
        J[k,k] = 1
        A = J.transpose()*A*J
    return A

Observe that, to get J_k , random integers are taken from min_int to max_int, both included. Likewise, niter is the number of matrices J_k that are constructed.

For example, the following code

set_random_seed(1000)
A = unimod_sym_matrix(6,3); show(A)
A = unimod_sym_matrix(6); show(A)

yields

\begin{pmatrix} 2 & -2 & -1 & 0 & 1 & 0 \\ -2 & 9 & 4 & 1 & -3 & 3 \\ -1 & 4 & 2 & 0 & -1 & 2 \\ 0 & 1 & 0 & 2 & -1 & -1 \\ 1 & -3 & -1 & -1 & 2 & 1 \\ 0 & 3 & 2 & -1 & 1 & 6 \end{pmatrix}

\begin{pmatrix} 2 & 0 & 1 & -1 & -1 & -1 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ 1 & 0 & 2 & -1 & -1 & -1 \\ -1 & 0 & -1 & 1 & 1 & 1 \\ -1 & 0 & -1 & 1 & 2 & 1 \\ -1 & 0 & -1 & 1 & 1 & 2 \end{pmatrix}

I have included set_random_seed just for the sake of reproducibility of the above output.

Since you want matrix entries between -2 and 2, you had better use the default values for niter, min_int and max_int.

Preview: (hide)
link

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

Stats

Asked: 5 years ago

Seen: 640 times

Last updated: Mar 11 '20