Add two-qubit lambda in qubitvector

[hhorii]

Summary

This PR optimizes two-qubit gate operations in qubitvector

Details and comments

Currently two-qubit operations use index. An index is generated for each iteration in a loop.

This PR introduces a new apply_lambda for two-qubits operations. cnot is implemented as follows:

  const int_t GAP_ctrl = BITS[qubit_ctrl];
  const int_t GAP_trgt = BITS[qubit_trgt];
  auto lambda = [&](const int_t i00)->void {
    auto i10 = i00 | GAP_ctrl;
    std::swap(data_[i10], data_[i10 | GAP_trgt]);
  };
  apply_lambda(lambda, qubit_ctrl, qubit_trgt);

In the above, i01 = i00 | GAP_trgt and i11 = i 10 | GAP_trgt. This PR may increases arithmetic operations, but they are cheap comparing with index generation.

My local environment (MacBook Pro (Retina, 13-inch, Early 2015)), cnot x 25 for a 25-qubits circuit took 3.19 sec with this optimization and 9.55 sec without this optimization.

[chriseclectic]

@hhorii PR #124 removes these two qubit gates and replaces them with general multi-controlled versions.

If performance is a lot better can you start from that PR and look at adding this back in as special case for N=2 qubit gates.

[hhorii]

Static index with std::array can show the similar performance. This latest version add static index on top of #124.

[chriseclectic]

EDIT: I just realised i was running the stabilizer simulator, not statevector doh, now i see a big difference. I think PR where we removed static indexing earlier was flawed.

Test case

Using a single-shot test circuit of N CNOTS for N qubits (N=10,15,20,25)

from qiskit import *
sim = Aer.get_backend('qasm_simulator')
opts = {'method': 'statevector'}

def make_qobj(nq):
    q = QuantumRegister(nq)
    qc = QuantumCircuit(q)
    for j in range(nq):
        qc.cx(q[j], q[(j + 1) % nq])
    return compile(qc, Aer.get_backend('qasm_simulator'), shots=1)

nqs = [10, 15, 20, 25]
qobjs = [make_qobj(nq) for nq in nqs]

and simple timing using Jupyter

%%timeit
result = sim.run(qobj, validate=False, backend_options=opts).result()

EDIT: UPDATED
I get results on my laptop of

Using dynamic indexes:

• N=10: 16.2 ms ± 2.65 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
• N=15: 38.6 ms ± 6.95 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
• N=20: 435 ms ± 43.3 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
• N=25: 11.5 s ± 404 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

Using static indexes:

• N=10: 14 ms ± 554 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
• N=15: 30.2 ms ± 6.03 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
• N=20: 99.4 ms ± 11.3 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
• N=25: 1.51 s ± 155 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

[chriseclectic]

Some further testing and I think we can remove the single-qubit apply matrix functions and use the old static indexing. Changing the above example to apply and x gate instead of cx and comparing the current 1-qubit apply_lambda to an N=1 static indexed case I get that they are about the same (or that static indexing is slightly faster).