mirror_circuit_rb.cpp 15.2 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
200
201
202
203
204
205
206
207
208
209
210
211
std::pair<bool, xacc::HeterogeneousMap> MirrorCircuitValidator::validate(
    std::shared_ptr<xacc::Accelerator> qpu,
    std::shared_ptr<qcor::CompositeInstruction> program,
    xacc::HeterogeneousMap options) {
  // Some default values
  int n_trials = 1000;
  // 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");
  }
212

213
214
215
216
  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
217
218
219
    auto mirror_data = qcor::MirrorCircuitValidator::createMirrorCircuit(program);
    auto mirror_circuit = mirror_data.first;
    auto expected_result = mirror_data.second;
220
221
    const std::string expectedBitString = [&]() {
      std::string bitStr;
222
223
224
      if (qpu->getBitOrder() == xacc::Accelerator::BitOrder::MSB) {
        std::reverse(expected_result.begin(), expected_result.end());
      }
225
226
227
228
229
230
231
232
233
234
235
236
237
      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());
238
239
240
241
242
243
244
245
    {
      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());
    }
246
247
248
249
250
251
252
253
254
255
256
257
    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);
}

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

264
265
266
267
268
269
270
271
272
273
274
275
  // 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();
276
277
278
279
  // Tracking the Pauli layer as it is commuted through
  std::vector<qcor::utils::PauliLabel> net_paulis(n,
                                                  qcor::utils::PauliLabel::I);

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

  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:
290
            result.emplace_back(std::make_pair("I", std::vector<int>{i}));
291
292
            break;
          case qcor::utils::PauliLabel::X:
293
            result.emplace_back(std::make_pair("X", std::vector<int>{i}));
294
295
            break;
          case qcor::utils::PauliLabel::Y:
296
            result.emplace_back(std::make_pair("Y", std::vector<int>{i}));
297
298
            break;
          case qcor::utils::PauliLabel::Z:
299
            result.emplace_back(std::make_pair("Z", std::vector<int>{i}));
300
301
302
303
304
305
306
            break;
          default:
            __builtin_unreachable();
          }
        }
        return result;
      };
307
308

  // program->as_xacc()->toGraph()->write(std::cout);
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
  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});
  };
326
327

  const auto d = program->depth();
328
  for (int layer = d - 1; layer >= 0; --layer) {
329
    auto current_layers = getLayer(program, layer);
330
    for (const auto &gate : current_layers) {
331
332
333
      if (gate->bits().size() == 1) {
        assert(gate->name() == "U");
        const auto u3_angles = decomposeU3Angle(gate);
334
335
336
        const auto [theta1_inv, theta2_inv, theta3_inv] =
            qcor::utils::invU3Gate(u3_angles);
        const size_t qubit = gate->bits()[0];
337
        program->addInstruction(gateProvider->createInstruction(
338
339
            "U", {qubit},
            {theta2_inv - M_PI, theta1_inv - 3.0 * M_PI, theta3_inv}));
340
      } else {
341
342
        assert(gate->name() == "CNOT");
        program->addInstruction(gate->clone());
343
344
      }
    }
345
346
  }
  const int newDepth = program->depth();
347
  for (int layer = 0; layer < newDepth; ++layer) {
348
    auto current_layers = getLayer(program, layer);
349
350
351
352
353
354
355
356
357
358
    // 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)]);
      }
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
359
360
361
362
363
364
365
366
      {
        std::stringstream ss;
        ss << "Random Pauli: ";
        for (const auto &p : random_paulis) {
          ss << p << " ";
        }
        xacc::info(ss.str());
      }
367
368
369
370

      return random_paulis;
    }(n);

Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
371
372
373
374
375
376
377
378
379
380
381
    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
382
383
      }

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

    const auto current_net_paulis_as_layer = pauliListToLayer(net_paulis);
389
    for (const auto &gate : current_layers) {
390
391
      if (gate->bits().size() == 1) {
        assert(gate->name() == "U");
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
        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);
        {
          std::stringstream ss;
          ss << "Net Pauli: ";
          for (const auto &p : net_paulis) {
            ss << p << " ";
          }
          xacc::info(ss.str());
        }

409
        const size_t qubit = gate->bits()[0];
410
        const auto [theta1, theta2, theta3] = decomposeU3Angle(gate);
411
412
413
414
415
        // 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]);
416
417
        mirrorCircuit.emplace_back(
            createU3GateFromAngle(qubit, theta1_new, theta2_new, theta3_new));
418
      } else {
Nguyen, Thien Minh's avatar
Nguyen, Thien Minh committed
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
        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);
        {
          std::stringstream ss;
          ss << "Net Pauli: ";
          for (const auto &p : net_paulis) {
            ss << p << " ";
          }
          xacc::info(ss.str());
        }
438
439
440
441
442
443
444
445
446
447
448
449
      }
    }
  }

  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);
  }

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