Calculation of Dynamical Lie Algebra

Hello! With great interest I read the demo about the dynamical lie algebra. I would like to calculate the dimension of the DLA for arbitrary generator sets. However, I struggle finding out if the elements of my current set of operators are linearly independent. In the demo example, this was straight forward and visible by eye. However, I was wondering if there is a smooth “pennylane” way of finding out if two operators are linearly independent, such that this algorithm below actually calculates the DLA:

import numpy as np
import pennylane as qml
from pennylane import X, Y, Z, I
from copy import deepcopy
qml.operation.enable_new_opmath()


# Generator Set for TFIM circuit with 3 qubits (also see https://arxiv.org/pdf/2105.14377.pdf Eq. 26  )
generators = [(1j + (Z(0) @ Z(1)) + (Z(1) @ Z(2))), 1j * X(0) + 1j * X(1) + 1j * X(2)]

# Following the terminology of the demo and the algorithm introduced in section E of this paper https://arxiv.org/abs/2105.14377
S = generators.copy()
"""Basis of the Lie algebra"""
S_new = (None,)
"""Set of matrices that was added to the Basis during the current iteration"""
S_prev = generators.copy()
"""Set of matrices that was added to the Basis during the last iteration"""

# algorithm:
while len(S_new) > 0:
    # reset the array of matrices we have added to the basis
    S_new = ()

    for op1 in generators:
        for op2 in S_prev:
            # calculating the nested commutator
            C = qml.commutator(op1, op2)/2
            C = C.simplify() 

            # Pseudo Code: 
            # Question would be how to check for linear independence. 
            if LinearIndependance(C,*S):
                S += (C,)
                S_new += (C,)

    S_prev = S_new

print('DLA:' S, 'and dim DLA:' len(S))

Thank you!

Hi @georg !

Glad you’ve liked my demo. This is an excellent question and something I have been actually working on in the meantime since releasing the demo. The problem comes down to asking generally for a vector v to be linearly independent of the columns of a matrix A. Take for example this collection of 4-dimensional vectors and two candidate vectors v1 and v2

A = np.array([
    [1, 1, 1,],
    [1, 1, 0,],
    [1, 0, 0,],
    [1, 0, 1,]
])

v1 = np.array([0, 0, 1, 0,]) # is linearly independent
v2 = np.array([1, 1, 0, 0,]) # is linearly dependent

If A is orthonormal, we can check if it spans v by checking if vres = v - A@(np.conj(A).T @ v) has non-zero norm. If A is not orthonormal (as is the case above) we can instead orthonormalize it and check:

>>> v1res = v1 - A @ np.linalg.inv(np.conj(A.T) @ A) @ np.conj(A).T @ v1
>>> np.linalg.norm(v1res) > 0.0
True
>>> v2res = v2 - A @ np.linalg.inv(np.conj(A.T) @ A) @ np.conj(A).T @ v2
>>> np.linalg.norm(v2res) > 0.0
False

This is rather general. For Pauli operators we just need to translate their pauli representation in a matrix form A by mapping Pauli words (tensor products of paulis) to rows that build up the matrix (this is kind of keyword based sparse representation of the Pauli sentence (a linear combination of pauli words).

I am currently implementing this in PennyLane. It already works, just hasnt been merged to the master branch yet. It will be part of the next release in ~2 months. You can already use it if you install PennyLane from source (see pennylane.ai/install >> GitHub) and then git checkout lieclosure2. On that branch, you can then do (I modified your example to the correct generators from the linked paper)

# Generator Set for TFIM circuit with 3 qubits (also see https://arxiv.org/pdf/2105.14377.pdf Eq. 26  )
n = 3
generators = [Z(i) @ Z(i+1) for i in range(n-1)]
generators += [X(i) for i in range(n)]

dla = qml.dla.lie_closure(generators)
dla = [ps.operation() for ps in dla] # currently the function returns not Operators but qml.pauli.PauliSentence instances, this will likely change
dla # note also that we dont care about minus signs at the moment, also something that may change
[Z(1) @ Z(0),
 Z(2) @ Z(1),
 X(0),
 X(1),
 X(2),
 Z(1) @ Y(0),
 Y(1) @ Z(0),
 Z(2) @ Y(1),
 Y(2) @ Z(1),
 -1.0 * (Z(2) @ X(1) @ Z(0)),
 -1.0 * (Y(1) @ Y(0)),
 -1.0 * (Y(2) @ Y(1)),
 -1.0 * (Z(2) @ X(1) @ Y(0)),
 Y(2) @ X(1) @ Z(0),
 -1.0 * (Y(2) @ X(1) @ Y(0))]

Two words of caution:

  • this is an experimental feature, the API and functionality may change and bugs may be undetected yet. Development is happening in #5161 and 5169
  • since you are looking into barren plateaus and all that jazz, you’ll find that almost all DLAs are exponentially sized. So dont be surprised if things for n>8 qubits wont work / take forever (the scaling here is 4**n). Only few exceptions, like the Ising model, have a polynomial size. In 1D there is a full characterization of all possible DLAs in this paper.

Hope that helps! Let me know if you have any further questions :slight_smile: