QSVT Application to non-Hermitian matrices(with hermitian embedding)

Hello! If applicable, put your complete code example down below. Make sure that your code:

  • is 100% self-contained — someone can copy-paste exactly what is here and run it to
    reproduce the behaviour you are observing
  • includes comments
import pennylane as qml
from pennylane import numpy as np
import pyqsp
from pyqsp import angle_sequence

def generate_phases_pennylane_QSVT(normalised_matrix):
    kappa = np.linalg.cond(normalised_matrix) + 2
    pcoefs, s = pyqsp.poly.PolyOneOverX().generate(kappa, return_coef=True, ensure_bounded=True, return_scale=True)
    phi_pyqsp = pyqsp.angle_sequence.QuantumSignalProcessingPhases(pcoefs, signal_operator="Wx", tolerance=0.00001)
    phi_qsvt = qml.transform_angles(phi_pyqsp, "QSP", "QSVT")  # convert pyqsp angles to be compatible with QSVT
    return phi_qsvt, s, kappa
A = np.array(
    [
        [-1, 1, 0, 3],
        [ 0, -1, 4, 0],
        [ 0, 0, -1, 6],
        [ 1, 0, -1, 0]
    ]
)
H = np.block([[np.zeros((A.shape[0], A.shape[0]), dtype=complex), A],
              [A.conj().T, np.zeros((A.shape[0], A.shape[0]), dtype=complex)]])

# Normalize for block-encoding
alpha_H = np.linalg.norm(H, 2)
H_norm = H / alpha_H if alpha_H > 0 else H
A_norm = A / alpha_H

phases, s, kappa = generate_phases_pennylane_QSVT(A_norm)
def pad_b(b, H):
    ket1 = np.array([0, 1])
    if H.shape[0] / b.shape[0] == 2:
        return np.kron(ket1, b)
    else:
        print('b shape mismatch')
bb = np.array([1, 2, 3, 4])
b = pad_b(bb, H)
norm_b = np.linalg.norm(b)
normalised_b = b / norm_b
nH = int(np.log2(H.shape[0]))
dimH = 2 ** nH
sysH_wires = [f"sH{i}" for i in range(0, nH)]
ancH_q = nH + 1
ancH_wires = [f"aH{i}" for i in range(0, ancH_q)]
beH_wires = ancH_wires + sysH_wires
lcuH_wires = [f"lcuH{i}" for i in range(2)]
deviceH_wires = lcuH_wires + beH_wires
dev = qml.device("default.qubit", deviceH_wires)

be_HT = qml.FABLE(H_norm.conj().T, beH_wires)
be_H = qml.FABLE(H_norm, beH_wires)

def sum_even_odd(blockOp, phi_all, anc_parity):
    half = len(phi_all)//2
    phi_even, phi_odd = phi_all[:half], phi_all[half:]
    projectors_even = [qml.PCPhase(a, dim=2**nH, wires=beH_wires) for a in phi_even]
    projectors_odd  = [qml.PCPhase(a, dim=2**nH, wires=beH_wires) for a in phi_odd]
    qml.Hadamard(anc_parity)
    qml.ctrl(qml.QSVT, control=(anc_parity,), control_values=(0,))(blockOp, projectors_even)
    qml.ctrl(qml.QSVT, control=(anc_parity,), control_values=(0,))(blockOp, projectors_odd)
    qml.Hadamard(anc_parity)

def real_part_inverse(phi_all):
    # First branch: U on A^T
    qml.Hadamard(wires=lcuH_wires[0])
    qml.ctrl(sum_even_odd, control=(lcuH_wires[0],), control_values=(0,))(be_HT, phi_all, lcuH_wires[1])

    # Second branch: U* on A (adjoint of the same construction on A)
    qml.ctrl(qml.adjoint(sum_even_odd), control=(lcuH_wires[0],), control_values=(1,))(be_H, phi_all, lcuH_wires[1])
    qml.Hadamard(wires=lcuH_wires[0])

@qml.qnode(dev)
def qsvt_apply():
    # prepare |b> on the system; all ancillas start in |0>
    qml.StatePrep(normalised_b, sysH_wires)
    real_part_inverse(phases)
    return qml.state()

state = qsvt_apply()
psi = state.reshape(2, 2, 2**ancH_q, 2**nH)
sys_amp = psi[0, 0, 0, :]
x_est = (sys_amp * np.linalg.norm(b)) / s
x_est_n = x_est / np.linalg.norm(x_est)

x_true = np.linalg.solve(A, bb)
x_true_n = x_true / np.linalg.norm(x_true)

print("target x  :", np.round(x_true_n, 3))
print("computed x:", np.round(x_est_n, 3))

If you want help with diagnosing an error, please put the full error message below:

So what is happening is that for the system that I have specified, the algorithm is giving wrong result. But for the QSVT in Practice tutorial system, the code is working fine and giving correct result. At first I thought that it was because my matrix wasn’t hermitian. So I created a hermitian embedding using the standard trick: H = [[0, A], [A_dagger, 0]]. Then again for the tutorial system, this is giving correct result, but for the system I am considering, it is not giving correct result. Are there some more criteria for a matrix to satisfy for this to work?
I would be very grateful for your response and help.
Thank you!

And, finally, make sure to include the versions of your packages. Specifically, show us the output of qml.about().

Hi @Sxal ,

Welcome to the Forum!

It looks like you have modified the code quite a bit with respect to the demo. I would recommend going back to the implementation from the demo and changing just one thing at a time. If you make a change and it results in an error please let us know exactly what code you changed with respect to the demo. This is important so that we can find the source of the error!

Hi. So if I rephrase my question then it is: The code shown in the tutorial is only working for hermitian matrices. I tried to run for non-hermitian matrices, but got very poor results. So I tried doing the hermitian embedding: H = [[0, A], [A_dg, 0]], but this is also failing.

So is the algorithm only designed to work with hermitian matices? And if yes, then should I change how the hermitian embedding applied so something else?

Hi @Sxal , thanks for clarifying your question.

The QSVT algorithm itself works with non-Hermitian matrices. However, the method described in the section on solving a linear system with QSVT (within the QSVT in Practice demo), does assume that the matrix we will invert is hermitian.

I don’t know whether or not block encoding the matrix is enough to love the restriction for the part about solving the linear system. Is your goal specifically about solving linear systems? Or were you interested in using QSVT more generally for somehting else?