mirror_circuit_rb.cpp 15.5 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
/*******************************************************************************
 * Copyright (c) 2018-, UT-Battelle, LLC.
 * All rights reserved. This program and the accompanying materials
 * are made available under the terms of the BSD 3-Clause License 
 * which accompanies this distribution. 
 *
 * Contributors:
 *   Alexander J. McCaskey - initial API and implementation
 *   Thien Nguyen - implementation
 *******************************************************************************/
11
#include "mirror_circuit_rb.hpp"
12
#include "AllGateVisitor.hpp"
13
14
15
16
#include "clifford_gate_utils.hpp"
#include "qcor_ir.hpp"
#include "qcor_pimpl_impl.hpp"
#include "xacc.hpp"
17
#include "xacc_plugin.hpp"
18
19
20
#include "xacc_service.hpp"
#include <cassert>
#include <random>
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
namespace {
std::vector<std::shared_ptr<xacc::Instruction>>
getLayer(std::shared_ptr<xacc::CompositeInstruction> circuit, int layerId) {
  std::vector<std::shared_ptr<xacc::Instruction>> result;
  assert(layerId < circuit->depth());
  auto graphView = circuit->toGraph();
  for (int i = 1; i < graphView->order() - 1; i++) {
    auto node = graphView->getVertexProperties(i);
    if (node.get<int>("layer") == layerId) {
      result.emplace_back(
          circuit->getInstruction(node.get<std::size_t>("id") - 1)->clone());
    }
  }
  assert(!result.empty());
  return result;
}
} // namespace
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
namespace xacc {
namespace quantum {
// Helper to convert a gate
class GateConverterVisitor : public AllGateVisitor {
public:
  GateConverterVisitor() {
    m_gateRegistry = xacc::getService<xacc::IRProvider>("quantum");
    m_program = m_gateRegistry->createComposite("temp_composite");
  }

  // Keep these 2 gates:
  void visit(CNOT &cnot) override { m_program->addInstruction(cnot.clone()); }
  void visit(U &u) override { m_program->addInstruction(u.clone()); }

  // Rotation gates:
  void visit(Ry &ry) override {
    const double theta = InstructionParameterToDouble(ry.getParameter(0));
    m_program->addInstruction(m_gateRegistry->createInstruction(
        "U", {ry.bits()[0]}, {theta, 0.0, 0.0}));
  }

  void visit(Rx &rx) override {
    const double theta = InstructionParameterToDouble(rx.getParameter(0));
    m_program->addInstruction(m_gateRegistry->createInstruction(
        "U", {rx.bits()[0]}, {theta, -1.0 * M_PI / 2.0, M_PI / 2.0}));
  }

  void visit(Rz &rz) override {
    const double theta = InstructionParameterToDouble(rz.getParameter(0));
    m_program->addInstruction(m_gateRegistry->createInstruction(
        "U", {rz.bits()[0]}, {0.0, theta, 0.0}));
  }

  void visit(X &x) override {
    Rx rx(x.bits()[0], M_PI);
    visit(rx);
  }
  void visit(Y &y) override {
    Ry ry(y.bits()[0], M_PI);
    visit(ry);
  }
  void visit(Z &z) override {
    Rz rz(z.bits()[0], M_PI);
    visit(rz);
  }
  void visit(S &s) override {
    Rz rz(s.bits()[0], M_PI / 2.0);
    visit(rz);
  }
  void visit(Sdg &sdg) override {
    Rz rz(sdg.bits()[0], -M_PI / 2.0);
    visit(rz);
  }
  void visit(T &t) override {
    Rz rz(t.bits()[0], M_PI / 4.0);
    visit(rz);
  }
  void visit(Tdg &tdg) override {
    Rz rz(tdg.bits()[0], -M_PI / 4.0);
    visit(rz);
  }

  void visit(Hadamard &h) override {
    m_program->addInstruction(m_gateRegistry->createInstruction(
        "U", {h.bits()[0]}, {M_PI / 2.0, 0.0, M_PI}));
  }

  void visit(Measure &measure) override {
106
    // Ignore measure
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
  }

  void visit(Identity &i) override {}

  void visit(CY &cy) override {
    // controlled-Y = Sdg(target) - CX - S(target)
    CNOT c1(cy.bits());
    Sdg sdg(cy.bits()[1]);
    S s(cy.bits()[1]);

    visit(sdg);
    visit(c1);
    visit(s);
  }

  void visit(CZ &cz) override {
    // CZ = H(target) - CX - H(target)
    CNOT c1(cz.bits());
    Hadamard h1(cz.bits()[1]);
    Hadamard h2(cz.bits()[1]);

    visit(h1);
    visit(c1);
    visit(h2);
  }

  void visit(CRZ &crz) override {
    const auto theta = InstructionParameterToDouble(crz.getParameter(0));
    // Decompose
    Rz rz1(crz.bits()[1], theta / 2);
    CNOT c1(crz.bits());
    Rz rz2(crz.bits()[1], -theta / 2);
    CNOT c2(crz.bits());

    // Revisit:
    visit(rz1);
    visit(c1);
    visit(rz2);
    visit(c2);
  }

  void visit(CH &ch) override {
    // controlled-H = Ry(pi/4, target) - CX - Ry(-pi/4, target)
    CNOT c1(ch.bits());
    Ry ry1(ch.bits()[1], M_PI_4);
    Ry ry2(ch.bits()[1], -M_PI_4);

    visit(ry1);
    visit(c1);
    visit(ry2);
  }

159
  std::shared_ptr<CompositeInstruction> getProgram() { return balanceLayer(); }
160
161

private:
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
  std::shared_ptr<CompositeInstruction> balanceLayer() {
    // Balance all layers
    // The mirroring protocol doesn't work well when there are layers that has
    // no gate on a qubit line, hence, just add an Identity there:
    auto program = m_gateRegistry->createComposite("temp_composite");
    const auto d = m_program->depth();
    const auto nbQubits = m_program->nPhysicalBits();
    for (int layer = 0; layer < d; ++layer) {
      std::set<int> qubits;
      auto current_layers = getLayer(m_program, layer);
      for (const auto &gate : current_layers) {
        program->addInstruction(gate);
        for (const auto &bit : gate->bits()) {
          qubits.emplace(bit);
        }
      }
      if (qubits.size() < nbQubits) {
        for (int i = 0; i < nbQubits; ++i) {
          if (!xacc::container::contains(qubits, i)) {
            program->addInstruction(
                m_gateRegistry->createInstruction("U", {(size_t)i}, {0.0, 0.0, 0.0}));
          }
        }
      }
    }
    return program;
  }
189
190
191
192
193
194
  std::shared_ptr<CompositeInstruction> m_program;
  std::shared_ptr<xacc::IRProvider> m_gateRegistry;
};
} // namespace quantum
} // namespace xacc

195
namespace qcor {
196
197
198
199
std::pair<bool, xacc::HeterogeneousMap> MirrorCircuitValidator::validate(
    std::shared_ptr<xacc::Accelerator> qpu,
    std::shared_ptr<qcor::CompositeInstruction> program,
    xacc::HeterogeneousMap options) {
200
201
202
203
  // Some default values: 4 Paulis for each qubit (double that to make sure we
  // cover more cases). Max is 1000 trials for many qubits.
  int n_trials = std::min(
      2 * static_cast<int>(std::pow(4, program->nPhysicalBits())), 1000);
204
205
206
207
208
209
210
211
212
213
  // 10% error allowed... (away from the expected bitstring)
  double eps = 0.1;
  int n_shots = 1024;

  if (options.keyExists<int>("trials")) {
    n_trials = options.get<int>("trials");
  }
  if (options.keyExists<double>("epsilon")) {
    eps = options.get<double>("epsilon");
  }
214

215
216
217
218
  qpu->updateConfiguration({{"shots", n_shots}});
  std::vector<double> trial_success_probs;
  auto provider = xacc::getIRProvider("quantum");
  for (int i = 0; i < n_trials; ++i) {
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
219
220
221
    auto mirror_data = qcor::MirrorCircuitValidator::createMirrorCircuit(program);
    auto mirror_circuit = mirror_data.first;
    auto expected_result = mirror_data.second;
222
223
    const std::string expectedBitString = [&]() {
      std::string bitStr;
224
225
226
      if (qpu->getBitOrder() == xacc::Accelerator::BitOrder::MSB) {
        std::reverse(expected_result.begin(), expected_result.end());
      }
227
228
229
230
231
232
233
234
235
236
237
238
239
      for (const auto &bit : expected_result) {
        bitStr += std::to_string(bit);
      }
      return bitStr;
    }();

    for (size_t qId = 0; qId < mirror_circuit->nPhysicalBits(); ++qId) {
      mirror_circuit->addInstruction(
          provider->createInstruction("Measure", {qId}));
    }

    auto mc_buffer = xacc::qalloc(mirror_circuit->nPhysicalBits());
    qpu->execute(mc_buffer, mirror_circuit->as_xacc());
240
241
242
243
244
245
246
247
    {
      std::stringstream ss;
      ss << "Trial " << i << "\n";
      ss << "Circuit:\n" << mirror_circuit->toString() << "\n";
      ss << "Result:\n" << mc_buffer->toString() << "\n";
      ss << "Expected bitstring:" << expectedBitString << "\n";
      xacc::info(ss.str());
    }
248
249
250
251
252
253
254
255
256
257
258
259
    const auto bitStrProb =
        mc_buffer->computeMeasurementProbability(expectedBitString);
    trial_success_probs.emplace_back(bitStrProb);
  }
  xacc::HeterogeneousMap data;
  data.insert("trial-success-probabilities", trial_success_probs);
  const bool pass_fail_result =
      std::all_of(trial_success_probs.begin(), trial_success_probs.end(),
                  [&eps](double val) { return val > 1.0 - eps; });
  return std::make_pair(pass_fail_result, data);
}

260
std::pair<std::shared_ptr<CompositeInstruction>, std::vector<bool>>
261
262
MirrorCircuitValidator::createMirrorCircuit(
    std::shared_ptr<CompositeInstruction> in_circuit) {
263
264
265
  std::vector<std::shared_ptr<xacc::Instruction>> mirrorCircuit;
  auto gateProvider = xacc::getService<xacc::IRProvider>("quantum");

266
267
268
269
270
271
272
273
274
275
276
277
  // Gate conversion:
  xacc::quantum::GateConverterVisitor visitor;
  xacc::InstructionIterator it(in_circuit->as_xacc());
  while (it.hasNext()) {
    auto nextInst = it.next();
    if (nextInst->isEnabled() && !nextInst->isComposite()) {
      nextInst->accept(&visitor);
    }
  }

  auto program = visitor.getProgram();
  const int n = program->nPhysicalBits();
278
279
280
281
  // Tracking the Pauli layer as it is commuted through
  std::vector<qcor::utils::PauliLabel> net_paulis(n,
                                                  qcor::utils::PauliLabel::I);

282
  // Sympletic group
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
283
  const auto srep_dict = qcor::utils::computeGateSymplecticRepresentations();
284
285
286
287
288
289
290
291

  const auto pauliListToLayer =
      [](const std::vector<qcor::utils::PauliLabel> &in_paulis) {
        qcor::utils::CliffordGateLayer_t result;
        for (int i = 0; i < in_paulis.size(); ++i) {
          const auto pauli = in_paulis[i];
          switch (pauli) {
          case qcor::utils::PauliLabel::I:
292
            result.emplace_back(std::make_pair("I", std::vector<int>{i}));
293
294
            break;
          case qcor::utils::PauliLabel::X:
295
            result.emplace_back(std::make_pair("X", std::vector<int>{i}));
296
297
            break;
          case qcor::utils::PauliLabel::Y:
298
            result.emplace_back(std::make_pair("Y", std::vector<int>{i}));
299
300
            break;
          case qcor::utils::PauliLabel::Z:
301
            result.emplace_back(std::make_pair("Z", std::vector<int>{i}));
302
303
304
305
306
307
308
            break;
          default:
            __builtin_unreachable();
          }
        }
        return result;
      };
309
310

  // program->as_xacc()->toGraph()->write(std::cout);
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
  const auto decomposeU3Angle = [](xacc::InstPtr u3_gate) {
    const double theta = InstructionParameterToDouble(u3_gate->getParameter(0));
    const double phi = InstructionParameterToDouble(u3_gate->getParameter(1));
    const double lam = InstructionParameterToDouble(u3_gate->getParameter(2));
    // Convert to 3 rz angles:
    const double theta1 = lam;
    const double theta2 = theta + M_PI;
    const double theta3 = phi + 3.0 * M_PI;
    return std::make_tuple(theta1, theta2, theta3);
  };

  const auto createU3GateFromAngle = [](size_t qubit, double theta1,
                                        double theta2, double theta3) {
    auto gateProvider = xacc::getService<xacc::IRProvider>("quantum");
    return gateProvider->createInstruction(
        "U", {qubit}, {theta2 - M_PI, theta3 - 3.0 * M_PI, theta1});
  };
328
329

  const auto d = program->depth();
330
  for (int layer = d - 1; layer >= 0; --layer) {
331
    auto current_layers = getLayer(program, layer);
332
    for (const auto &gate : current_layers) {
333
334
335
      if (gate->bits().size() == 1) {
        assert(gate->name() == "U");
        const auto u3_angles = decomposeU3Angle(gate);
336
337
338
        const auto [theta1_inv, theta2_inv, theta3_inv] =
            qcor::utils::invU3Gate(u3_angles);
        const size_t qubit = gate->bits()[0];
339
        program->addInstruction(gateProvider->createInstruction(
340
341
            "U", {qubit},
            {theta2_inv - M_PI, theta1_inv - 3.0 * M_PI, theta3_inv}));
342
      } else {
343
344
        assert(gate->name() == "CNOT");
        program->addInstruction(gate->clone());
345
346
      }
    }
347
348
  }
  const int newDepth = program->depth();
349
  for (int layer = 0; layer < newDepth; ++layer) {
350
    auto current_layers = getLayer(program, layer);
351
352
353
354
355
356
357
358
359
360
    // New random Pauli layer
    const std::vector<qcor::utils::PauliLabel> new_paulis = [](int nQubits) {
      static std::random_device rd;
      static std::mt19937 gen(rd());
      static std::uniform_int_distribution<size_t> dis(
          0, qcor::utils::ALL_PAULI_OPS.size() - 1);
      std::vector<qcor::utils::PauliLabel> random_paulis;
      for (int i = 0; i < nQubits; ++i) {
        random_paulis.emplace_back(qcor::utils::ALL_PAULI_OPS[dis(gen)]);
      }
361
362
363
364
365
366
367
368
      // {
      //   std::stringstream ss;
      //   ss << "Random Pauli: ";
      //   for (const auto &p : random_paulis) {
      //     ss << p << " ";
      //   }
      //   xacc::info(ss.str());
      // }
369
370
371
372

      return random_paulis;
    }(n);

Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
373
374
375
376
377
378
379
380
381
382
383
    const auto gateToLayerInfo = [](xacc::InstPtr gate, int nbQubits) {
      qcor::utils::CliffordGateLayer_t result;
      std::vector<int> operands;
      for (const auto &bit : gate->bits()) {
        operands.emplace_back(bit);
      }

      for (int i = 0; i < nbQubits; ++i) {
        if (!xacc::container::contains(operands, i)) {
          result.emplace_back(std::make_pair("I", std::vector<int>{i}));
        }
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
384
385
      }

Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
386
387
388
389
390
      result.emplace_back(std::make_pair(gate->name(), operands));
      return result;
    };

    const auto current_net_paulis_as_layer = pauliListToLayer(net_paulis);
391
    for (const auto &gate : current_layers) {
392
393
      if (gate->bits().size() == 1) {
        assert(gate->name() == "U");
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
394
395
396
397
398
399
400
401
        const auto new_paulis_as_layer = pauliListToLayer(new_paulis);
        const auto new_net_paulis_reps =
            qcor::utils::computeCircuitSymplecticRepresentations(
                {new_paulis_as_layer, current_net_paulis_as_layer}, n,
                srep_dict);

        // Update the tracking net
        net_paulis = qcor::utils::find_pauli_labels(new_net_paulis_reps.second);
402
403
404
405
406
407
408
409
        // {
        //   std::stringstream ss;
        //   ss << "Net Pauli: ";
        //   for (const auto &p : net_paulis) {
        //     ss << p << " ";
        //   }
        //   xacc::info(ss.str());
        // }
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
410

411
        const size_t qubit = gate->bits()[0];
412
        const auto [theta1, theta2, theta3] = decomposeU3Angle(gate);
413
414
415
416
417
        // Compute the pseudo_inverse gate:
        const auto [theta1_new, theta2_new, theta3_new] =
            qcor::utils::computeRotationInPauliFrame(
                std::make_tuple(theta1, theta2, theta3), new_paulis[qubit],
                net_paulis[qubit]);
418
419
        mirrorCircuit.emplace_back(
            createU3GateFromAngle(qubit, theta1_new, theta2_new, theta3_new));
420
      } else {
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
421
422
423
424
425
426
427
428
429
430
431
        mirrorCircuit.emplace_back(gate->clone());
        // we need to account for how the net pauli changes when it gets passed
        // through the clifford layers
        const auto new_net_paulis_reps =
            qcor::utils::computeCircuitSymplecticRepresentations(
                {gateToLayerInfo(gate, n), current_net_paulis_as_layer,
                 gateToLayerInfo(gate, n)},
                n, srep_dict);

        // Update the tracking net
        net_paulis = qcor::utils::find_pauli_labels(new_net_paulis_reps.second);
432
433
434
435
436
437
438
439
        // {
        //   std::stringstream ss;
        //   ss << "Net Pauli: ";
        //   for (const auto &p : net_paulis) {
        //     ss << p << " ";
        //   }
        //   xacc::info(ss.str());
        // }
440
441
442
443
444
445
446
447
448
449
450
451
      }
    }
  }

  const auto [telp_s, telp_p] =
      qcor::utils::computeLayerSymplecticRepresentations(
          pauliListToLayer(net_paulis), n, srep_dict);
  std::vector<bool> target_bitString;
  for (int i = n; i < telp_p.size(); ++i) {
    target_bitString.emplace_back(telp_p[i] == 2);
  }

452
  auto mirror_comp = gateProvider->createComposite(program->name() + "_MIRROR");
453
454
455
456
457
  mirror_comp->addInstructions(mirrorCircuit);
  return std::make_pair(std::make_shared<CompositeInstruction>(mirror_comp),
                        target_bitString);
}
} // namespace qcor
458
REGISTER_PLUGIN(qcor::MirrorCircuitValidator, qcor::BackendValidator)