Regarding pseudoinverse in Quantum Natural Gradient

PennyLane Help
1

Hello,

I am currently using QNGOptimizer to run optimization in the context of quantum chemistry. I recently came across with the error that notifies me of a singular matrix during the optimization, which occurs at some point during the optimization if I choose a finite number of shots in the device instead of using the analytic option. I overcame this issue by replacing the np.linalg.solve() function in the source code of the QNGOptimizer class by a matrix-vector multiplication where the matrix is a pseudoinverse of the geometric tensor (np.linalg.pinv()). I was wondering, since even the QNG tutorial uses a pseudoinverse instead of the usual matrix inverse, why the code still uses np.linalg.solve(). Am I just using the optimizer incorrectly?

Thank you very much.

2

Hi Tenzen! The source code originally also used np.linalg.pinvh(), and was written around the same time as the tutorial :slightly_smiling_face:

Later, we modified the source code to instead use np.linalg.solve(), as np.linalg.pinvh can be unstable for certain linear system. You can see the PR here: Replace pinv with np.linalg.solve in the QNG optimizer by josh146 · Pull Request #428 · PennyLaneAI/pennylane · GitHub

Your question is a very good one — we haven’t internally come across a case where pinv() provided a result, but solve() produced an error. Do you mind sharing minimal (non)-working PennyLane code that replicates this error? We can then explore what is exactly happening, and how best we can implement the QNG optimization.

4

Hi Josh, thank you for the response.
Here is a code that produces the error. It encounters a singular matrix most of the time. The circuit might not look simple, but it is just the alternating operator ansatz written out explicitly. I did this for a reason that is not related to the current issue, which is that the Multi-RZ gate used in the QAOAEmbedding does not seem to have a generator defined, so the metric tensor cannot be calculated without explicitly writing it out.

import pennylane as qml
from pennylane import numpy as np
from pennylane import expval
from pennylane.init import qaoa_embedding_normal
import autograd

n_q = 4 # number of qubits
n_l = 2 # number of layers
l_r = 0.01 # learninig rate / step size
n_it = 100 #number of iterations

dev = qml.device('default.qubit', wires=n_q, shots = 1024, analytic = False)

def qng_test(n_q, n_l, dev, l_r=0.01, max_iter=200, seed=0, diag_approx=False):
    @qml.qnode(dev)
    def circuit(params):
        # QAOA Ansatz
        for l in range(n_l):
            for q in range(n_q):
                qml.RX(features[q], wires=q)
            for q in range(n_q):
                qn = (q+1) % n_q
                p_idx = 2*l*n_q+q
                qml.CNOT(wires=[q,qn])
                qml.RZ(params[p_idx], wires=qn)
                qml.CNOT(wires=[q,qn])
            for q in range(n_q):
                p_idx = (2*l+1)*n_q+q
                qml.RY(params[p_idx], wires=q)
        for q in range(n_q):
            qml.RX(features[q], wires=q) 
        return qml.expval(qml.PauliZ(0))

    features = [i for i in range(n_q)]
    params = np.reshape(qaoa_embedding_normal(n_layers=n_l, n_wires=n_q, mean=0, std=0.2, seed=seed), 2*n_q*n_l)

    opt = qml.QNGOptimizer(stepsize=l_r, diag_approx=diag_approx)
    for i in range(max_iter):
        params = opt.step(circuit, params)

qng_test(n_q=n_q, n_l=n_l, dev=dev, l_r=l_r, max_iter=n_it, diag_approx=True)
5

Thanks @Tenzan_Araki!

In particular, thank you for pointing out that the MultiRZ operation currently doesn’t have a generator defined — the generator in this case is quite trivial, so we can work on adding support for that.

In the meantime, what you can do is ‘monkeypatch’ default.qubit so that it does not support the MultiRZ gate. This will cause PennyLane to automatically decompose the MultiRZ gate down to CNOTs and RZ gates:

dev = qml.device('default.qubit', wires=n_q, shots=1024, analytic=False)
dev.operations.remove("MultiRZ")

@qml.qnode(dev)
def circuit(params):
    QAOAEmbedding(features, params, wires=range(n_q), local_field="Y")
    return qml.expval(qml.PauliZ(0))

features = [i for i in range(n_q)]
params = qaoa_embedding_normal(n_layers=n_l, n_wires=n_q, mean=0, std=0.2)

opt = qml.QNGOptimizer(stepsize=l_r, diag_approx=True)
for i in range(n_it):
    params = opt.step(circuit, params)

This will create the circuit

 0: ──RX(0)──╭X──RZ(0.207)──╭X──────────────────────────────────────────────────────────────╭C─────────────╭C──RY(-1.222)──RX(0)──╭X──RZ(-0.032)──╭X─────────────────────────────────────────────────────────────╭C──────────────╭C──RY(-1.904)──RX(0)──┤ ⟨Z⟩
 1: ──RX(1)──╰C─────────────╰C──╭X──RZ(-0.607)──╭X───RY(-0.097)──RX(1)──────────────────────│──────────────│──────────────────────╰C──────────────╰C──╭X──RZ(-0.159)──╭X───RY(-0.025)──RX(1)─────────────────────│───────────────│──────────────────────┤
 2: ──RX(2)─────────────────────╰C──────────────╰C──╭X───────────RZ(0.082)──╭X──RY(-0.069)──│───RX(2)──────│──────────────────────────────────────────╰C──────────────╰C──╭X───────────RZ(0.154)──╭X──RY(0.189)──│───RX(2)───────│──────────────────────┤
 3: ──RX(3)─────────────────────────────────────────╰C──────────────────────╰C──────────────╰X──RZ(0.377)──╰X──RY(-0.315)──RX(3)──────────────────────────────────────────╰C──────────────────────╰C─────────────╰X──RZ(-0.017)──╰X──RY(-0.385)──RX(3)──┤

With respect to the singular matrix error, I suggest adding in a regularizer to the QNG optimizer, by passing the lam argument,

opt = QNGOptimizer(stepsize=0.01, diag_approx=True, lam=0.01)

This will add a small identity matrix to the metric tensor to avoid singularities (g\rightarrow g+\lambda I). You may find that you need to increase the regularizer further, sometimes even values of 0.1 are needed.

6

No problem, and thank you very much for the advice. Now it works fine for me.

1 Like