Quantum++  v1.0.0-beta2
C++11 quantum computing library
gates.h
Go to the documentation of this file.
1 /*
2  * Quantum++
3  *
4  * Copyright (c) 2013 - 2017 Vlad Gheorghiu (vgheorgh@gmail.com)
5  *
6  * This file is part of Quantum++.
7  *
8  * Quantum++ is free software: you can redistribute it and/or modify
9  * it under the terms of the GNU General Public License as published by
10  * the Free Software Foundation, either version 3 of the License, or
11  * (at your option) any later version.
12  *
13  * Quantum++ is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  * GNU General Public License for more details.
17  *
18  * You should have received a copy of the GNU General Public License
19  * along with Quantum++. If not, see <http://www.gnu.org/licenses/>.
20  */
21 
27 #ifndef CLASSES_GATES_H_
28 #define CLASSES_GATES_H_
29 
30 namespace qpp
31 {
32 
37 class Gates final : public internal::Singleton<const Gates> // const Singleton
38 {
39  friend class internal::Singleton<const Gates>;
40 
41 public:
42  // One qubit gates
43  cmat Id2{cmat::Identity(2, 2)};
44  cmat H{cmat::Zero(2, 2)};
45  cmat X{cmat::Zero(2, 2)};
46  cmat Y{cmat::Zero(2, 2)};
47  cmat Z{cmat::Zero(2, 2)};
48  cmat S{cmat::Zero(2, 2)};
49  cmat T{cmat::Zero(2, 2)};
50 
51  // two qubit gates
52  cmat CNOT{cmat::Identity(4, 4)};
53  cmat CZ{cmat::Identity(4, 4)};
54  cmat CNOTba{cmat::Zero(4, 4)};
55  cmat SWAP{cmat::Identity(4, 4)};
56 
57  // three qubit gates
58  cmat TOF{cmat::Identity(8, 8)};
59  cmat FRED{cmat::Identity(8, 8)};
60 private:
65  {
66  H << 1 / std::sqrt(2.), 1 / std::sqrt(2.),
67  1 / std::sqrt(2.), -1 / std::sqrt(2.);
68  X << 0, 1, 1, 0;
69  Z << 1, 0, 0, -1;
70  Y << 0, -1_i, 1_i, 0;
71  S << 1, 0, 0, 1_i;
72  T << 1, 0, 0, std::exp(1_i * pi / 4.0);
73  CNOT.block(2, 2, 2, 2) = X;
74  CNOTba(0, 0) = 1;
75  CNOTba(1, 3) = 1;
76  CNOTba(2, 2) = 1;
77  CNOTba(3, 1) = 1;
78  CZ(3, 3) = -1;
79 
80  SWAP.block(1, 1, 2, 2) = X;
81  TOF.block(6, 6, 2, 2) = X;
82  FRED.block(4, 4, 4, 4) = SWAP;
83  }
84 
88  ~Gates() = default;
89 
90 public:
91  // variable gates
92 
93  // one qubit gates
94 
103  cmat Rn(double theta, const std::vector<double>& n) const
104  {
105  // EXCEPTION CHECKS
106 
107  // check 3-dimensional vector
108  if (n.size() != 3)
109  throw Exception("qpp::Gates::Rn()",
110  "n is not a 3-dimensional vector!");
111  // END EXCEPTION CHECKS
112 
113  cmat result(2, 2);
114  result = std::cos(theta / 2) * Id2
115  - 1_i * std::sin(theta / 2) * (n[0] * X + n[1] * Y + n[2] * Z);
116 
117  return result;
118  }
119 
120  // one quDit gates
121 
130  cmat Zd(idx D) const
131  {
132  // EXCEPTION CHECKS
133 
134  if (D == 0)
135  throw Exception("qpp::Gates::Zd()", Exception::Type::DIMS_INVALID);
136  // END EXCEPTION CHECKS
137 
138  cmat result = cmat::Zero(D, D);
139  for (idx i = 0; i < D; ++i)
140  result(i, i) = std::pow(omega(D), i);
141 
142  return result;
143  }
144 
154  cmat Fd(idx D) const
155  {
156  // EXCEPTION CHECKS
157 
158  if (D == 0)
159  throw Exception("qpp::Gates::Fd()", Exception::Type::DIMS_INVALID);
160  // END EXCEPTION CHECKS
161 
162  cmat result(D, D);
163 
164 #ifdef WITH_OPENMP_
165 #pragma omp parallel for collapse(2)
166 #endif // WITH_OPENMP_
167  for (idx j = 0; j < D; ++j) // column major order for speed
168  for (idx i = 0; i < D; ++i)
169  result(i, j) = 1 / std::sqrt(D) * std::pow(omega(D), i * j);
170 
171  return result;
172  }
173 
183  cmat Xd(idx D) const
184  {
185  // EXCEPTION CHECKS
186 
187  if (D == 0)
188  throw Exception("qpp::Gates::Xd()", Exception::Type::DIMS_INVALID);
189  // END EXCEPTION CHECKS
190 
191  return Fd(D).inverse() * Zd(D) * Fd(D);
192  }
193 
203  template<typename Derived = Eigen::MatrixXcd>
204  Derived Id(idx D) const
205  {
206  // EXCEPTION CHECKS
207 
208  if (D == 0)
209  throw Exception("qpp::Gates::Id()", Exception::Type::DIMS_INVALID);
210  // END EXCEPTION CHECKS
211 
212  return Derived::Identity(D, D);
213  }
214 
230  template<typename Derived>
231  dyn_mat<typename Derived::Scalar> CTRL(const Eigen::MatrixBase<Derived>& A,
232  const std::vector<idx>& ctrl,
233  const std::vector<idx>& subsys,
234  idx N, idx d = 2) const
235  {
236  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
237 
238  // EXCEPTION CHECKS
239 
240  // check matrix zero size
242  throw Exception("qpp::Gates::CTRL()", Exception::Type::ZERO_SIZE);
243 
244  // check square matrix
246  throw Exception("qpp::Gates::CTRL()",
248 
249  // check lists zero size
250  if (ctrl.size() == 0)
251  throw Exception("qpp::Gates::CTRL()", Exception::Type::ZERO_SIZE);
252  if (subsys.size() == 0)
253  throw Exception("qpp::Gates::CTRL()", Exception::Type::ZERO_SIZE);
254 
255  // check out of range
256  if (N == 0)
257  throw Exception("qpp::Gates::CTRL()",
259 
260  // check valid local dimension
261  if (d == 0)
262  throw Exception("qpp::Gates::CTRL()",
264 
265  // ctrl + gate subsystem vector
266  std::vector<idx> ctrlgate = ctrl;
267  ctrlgate.insert(std::end(ctrlgate), std::begin(subsys),
268  std::end(subsys));
269  std::sort(std::begin(ctrlgate), std::end(ctrlgate));
270 
271  std::vector<idx> dims(N, d); // local dimensions vector
272 
273  // check that ctrl + gate subsystem is valid
274  // with respect to local dimensions
275  if (!internal::check_subsys_match_dims(ctrlgate, dims))
276  throw Exception("qpp::Gates::CTRL()",
278 
279  // check that subsys list match the dimension of the matrix
280  if (rA.rows() != std::llround(std::pow(d, subsys.size())))
281  throw Exception("qpp::Gates::CTRL()",
283  // END EXCEPTION CHECKS
284 
285  // Use static allocation for speed!
286  idx Cdims[maxn];
287  idx midx_row[maxn];
288  idx midx_col[maxn];
289 
290  idx CdimsA[maxn];
291  idx midxA_row[maxn];
292  idx midxA_col[maxn];
293 
294  idx Cdims_bar[maxn];
295  idx Csubsys_bar[maxn];
296  idx midx_bar[maxn];
297 
298  idx Ngate = subsys.size();
299  idx Nctrl = ctrl.size();
300  idx Nsubsys_bar = N - ctrlgate.size();
301  idx D = static_cast<idx>(std::llround(std::pow(d, N)));
302  idx DA = static_cast<idx>(rA.rows());
303  idx Dsubsys_bar = static_cast<idx>(
304  std::llround(std::pow(d, Nsubsys_bar)));
305 
306  // compute the complementary subsystem of ctrlgate w.r.t. dims
307  std::vector<idx> subsys_bar = complement(ctrlgate, N);
308  std::copy(std::begin(subsys_bar), std::end(subsys_bar),
309  std::begin(Csubsys_bar));
310 
311  for (idx k = 0; k < N; ++k)
312  {
313  midx_row[k] = midx_col[k] = 0;
314  Cdims[k] = d;
315  }
316 
317  for (idx k = 0; k < Nsubsys_bar; ++k)
318  {
319  Cdims_bar[k] = d;
320  midx_bar[k] = 0;
321  }
322 
323  for (idx k = 0; k < Ngate; ++k)
324  {
325  midxA_row[k] = midxA_col[k] = 0;
326  CdimsA[k] = d;
327  }
328 
330  typename Derived::Scalar>
331  ::Identity(D, D);
333 
334  // run over the complement indexes
335  for (idx i = 0; i < Dsubsys_bar; ++i)
336  {
337  // get the complement row multi-index
338  internal::n2multiidx(i, Nsubsys_bar, Cdims_bar, midx_bar);
339  for (idx k = 0; k < d; ++k)
340  {
341  Ak = powm(rA, k); // compute rA^k
342  // run over the subsys row multi-index
343  for (idx a = 0; a < DA; ++a)
344  {
345  // get the subsys row multi-index
346  internal::n2multiidx(a, Ngate, CdimsA, midxA_row);
347 
348  // construct the result row multi-index
349 
350  // first the ctrl part (equal for both row and column)
351  for (idx c = 0; c < Nctrl; ++c)
352  midx_row[ctrl[c]] = midx_col[ctrl[c]] = k;
353 
354  // then the complement part (equal for column)
355  for (idx c = 0; c < Nsubsys_bar; ++c)
356  midx_row[Csubsys_bar[c]] = midx_col[Csubsys_bar[c]] =
357  midx_bar[c];
358 
359  // then the subsys part
360  for (idx c = 0; c < Ngate; ++c)
361  midx_row[subsys[c]] = midxA_row[c];
362 
363  // run over the subsys column multi-index
364  for (idx b = 0; b < DA; ++b)
365  {
366  // get the subsys column multi-index
367  internal::n2multiidx(b, Ngate, CdimsA, midxA_col);
368 
369  // construct the result column multi-index
370  for (idx c = 0; c < Ngate; ++c)
371  midx_col[subsys[c]] = midxA_col[c];
372 
373  // finally write the values
374  result(internal::multiidx2n(midx_row, N, Cdims),
375  internal::multiidx2n(midx_col, N, Cdims))
376  = Ak(a, b);
377  }
378  }
379  }
380  }
381 
382  return result;
383  }
384 
400  template<typename Derived>
402  const Eigen::MatrixBase<Derived>& A, idx pos,
403  const std::vector<idx>& dims) const
404  {
405  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
406 
407  // EXCEPTION CHECKS
408 
409  // check zero-size
411  throw Exception("qpp::Gates::expandout()",
413 
414  // check that dims is a valid dimension vector
415  if (!internal::check_dims(dims))
416  throw Exception("qpp::Gates::expandout()",
418 
419  // check square matrix
421  throw Exception("qpp::Gates::expandout()",
423 
424  // check that position is valid
425  if (pos > dims.size() - 1)
426  throw Exception("qpp::Gates::expandout()",
428 
429  // check that dims[pos] match the dimension of A
430  if (static_cast<idx>(rA.rows()) != dims[pos])
431  throw Exception("qpp::Gates::expandout()",
433  // END EXCEPTION CHECKS
434 
435  idx D = std::accumulate(std::begin(dims), std::end(dims),
436  static_cast<idx>(1), std::multiplies<idx>());
438  typename Derived::Scalar>
439  ::Identity(D, D);
440 
441  idx Cdims[maxn];
442  idx midx_row[maxn];
443  idx midx_col[maxn];
444 
445  for (idx k = 0; k < dims.size(); ++k)
446  {
447  midx_row[k] = midx_col[k] = 0;
448  Cdims[k] = dims[k];
449  }
450 
451  // run over the main diagonal multi-indexes
452  for (idx i = 0; i < D; ++i)
453  {
454  // get row multi_index
455  internal::n2multiidx(i, dims.size(), Cdims, midx_row);
456  // get column multi_index (same as row)
457  internal::n2multiidx(i, dims.size(), Cdims, midx_col);
458  // run over the gate row multi-index
459  for (idx a = 0; a < static_cast<idx>(rA.rows());
460  ++a)
461  {
462  // construct the total row multi-index
463  midx_row[pos] = a;
464 
465  // run over the gate column multi-index
466  for (idx b = 0;
467  b < static_cast<idx>(rA.cols()); ++b)
468  {
469  // construct the total column multi-index
470  midx_col[pos] = b;
471 
472  // finally write the values
473  result(internal::multiidx2n(midx_row, dims.size(), Cdims),
474  internal::multiidx2n(midx_col, dims.size(), Cdims))
475  = rA(a, b);
476  }
477  }
478  }
479 
480  return result;
481  }
482 
483 }; /* class Gates */
484 
485 } /* namespace qpp */
486 
487 #endif /* CLASSES_GATES_H_ */
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:127
dyn_mat< typename Derived::Scalar > expandout(const Eigen::MatrixBase< Derived > &A, idx pos, const std::vector< idx > &dims) const
Expands out.
Definition: gates.h:401
cmat SWAP
SWAP gate.
Definition: gates.h:55
constexpr idx maxn
Maximum number of allowed qubits/qudits (subsystems)
Definition: constants.h:74
bool check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:198
~Gates()=default
Default destructor.
cplx omega(idx D)
D-th root of unity.
Definition: constants.h:96
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:78
idx multiidx2n(const idx *const midx, idx numdims, const idx *const dims) noexcept
Definition: util.h:71
void n2multiidx(idx n, idx numdims, const idx *const dims, idx *result) noexcept
Definition: util.h:48
cmat Zd(idx D) const
Generalized Z gate for qudits.
Definition: gates.h:130
Singleton policy class, used internally to implement the singleton pattern via CRTP (Curiously recurr...
Definition: singleton.h:77
dyn_mat< typename Derived::Scalar > CTRL(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &ctrl, const std::vector< idx > &subsys, idx N, idx d=2) const
Generates the multi-partite multiple-controlled-A gate in matrix form.
Definition: gates.h:231
const Singleton class that implements most commonly used gates
Definition: gates.h:37
Quantum++ main namespace.
Definition: codes.h:30
cmat S
S gate.
Definition: gates.h:48
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1814
cmat T
T gate.
Definition: gates.h:49
cmat FRED
Fredkin gate.
Definition: gates.h:59
cmat CNOTba
Controlled-NOT target control gate.
Definition: gates.h:54
cmat Fd(idx D) const
Fourier transform gate for qudits.
Definition: gates.h:154
cmat H
Hadamard gate.
Definition: gates.h:44
cmat Id2
Identity gate.
Definition: gates.h:43
cmat X
Pauli Sigma-X gate.
Definition: gates.h:45
Generates custom exceptions, used when validating function parameters.
Definition: exception.h:39
dyn_mat< typename Derived::Scalar > powm(const Eigen::MatrixBase< Derived > &A, idx n)
Fast matrix power based on the SQUARE-AND-MULTIPLY algorithm.
Definition: functions.h:778
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:99
Derived Id(idx D) const
Identity gate.
Definition: gates.h:204
constexpr double pi
Definition: constants.h:79
cmat Rn(double theta, const std::vector< double > &n) const
Qubit rotation of theta about the 3-dimensional real (unit) vector n.
Definition: gates.h:103
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:61
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:140
cmat CNOT
Controlled-NOT control target gate.
Definition: gates.h:52
std::size_t idx
Non-negative integer index.
Definition: types.h:36
cmat Z
Pauli Sigma-Z gate.
Definition: gates.h:47
cmat CZ
Controlled-Phase gate.
Definition: gates.h:53
cmat Y
Pauli Sigma-Y gate.
Definition: gates.h:46
Gates()
Initializes the gates.
Definition: gates.h:64
cmat TOF
Toffoli gate.
Definition: gates.h:58
cmat Xd(idx D) const
Generalized X gate for qudits.
Definition: gates.h:183