Quantum++  v1.0-rc2
A modern 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 {
36 class Gates final : public internal::Singleton<const Gates> // const Singleton
37 {
38  friend class internal::Singleton<const Gates>;
39 
40 public:
41  // One qubit gates
42  cmat Id2{cmat::Identity(2, 2)};
43  cmat H{cmat::Zero(2, 2)};
44  cmat X{cmat::Zero(2, 2)};
45  cmat Y{cmat::Zero(2, 2)};
46  cmat Z{cmat::Zero(2, 2)};
47  cmat S{cmat::Zero(2, 2)};
48  cmat T{cmat::Zero(2, 2)};
49 
50  // two qubit gates
51  cmat CNOT{cmat::Identity(4, 4)};
52  cmat CZ{cmat::Identity(4, 4)};
53  cmat CNOTba{cmat::Zero(4, 4)};
54  cmat SWAP{cmat::Identity(4, 4)};
55 
56  // three qubit gates
57  cmat TOF{cmat::Identity(8, 8)};
58  cmat FRED{cmat::Identity(8, 8)};
59 private:
64  {
65  H << 1 / std::sqrt(2.), 1 / std::sqrt(2.),
66  1 / std::sqrt(2.), -1 / std::sqrt(2.);
67  X << 0, 1, 1, 0;
68  Z << 1, 0, 0, -1;
69  Y << 0, -1_i, 1_i, 0;
70  S << 1, 0, 0, 1_i;
71  T << 1, 0, 0, std::exp(1_i * pi / 4.0);
72  CNOT.block(2, 2, 2, 2) = X;
73  CNOTba(0, 0) = 1;
74  CNOTba(1, 3) = 1;
75  CNOTba(2, 2) = 1;
76  CNOTba(3, 1) = 1;
77  CZ(3, 3) = -1;
78 
79  SWAP.block(1, 1, 2, 2) = X;
80  TOF.block(6, 6, 2, 2) = X;
81  FRED.block(4, 4, 4, 4) = SWAP;
82  }
83 
87  ~Gates() = default;
88 
89 public:
90  // variable gates
91 
92  // one qubit gates
93 
102  cmat Rn(double theta, const std::vector<double>& n) const
103  {
104  // EXCEPTION CHECKS
105 
106  // check 3-dimensional vector
107  if (n.size() != 3)
108  throw exception::CustomException("qpp::Gates::Rn()",
109  "n is not a 3-dimensional vector!");
110  // END EXCEPTION CHECKS
111 
112  cmat result(2, 2);
113  result = std::cos(theta / 2) * Id2
114  - 1_i * std::sin(theta / 2) * (n[0] * X + n[1] * Y + n[2] * Z);
115 
116  return result;
117  }
118 
119  // one quDit gates
120 
129  cmat Zd(idx D = 2) const
130  {
131  // EXCEPTION CHECKS
132 
133  if (D == 0)
134  throw exception::DimsInvalid("qpp::Gates::Zd()");
135  // END EXCEPTION CHECKS
136 
137  cmat result = cmat::Zero(D, D);
138  for (idx i = 0; i < D; ++i)
139  result(i, i) = std::pow(omega(D), i);
140 
141  return result;
142  }
143 
153  cmat Fd(idx D = 2) const
154  {
155  // EXCEPTION CHECKS
156 
157  if (D == 0)
158  throw exception::DimsInvalid("qpp::Gates::Fd()");
159  // END EXCEPTION CHECKS
160 
161  cmat result(D, D);
162 
163 #ifdef WITH_OPENMP_
164 #pragma omp parallel for collapse(2)
165 #endif // WITH_OPENMP_
166  for (idx j = 0; j < D; ++j) // column major order for speed
167  for (idx i = 0; i < D; ++i)
168  result(i, j) = 1 / std::sqrt(D) * std::pow(omega(D), i * j);
169 
170  return result;
171  }
172 
182  cmat Xd(idx D = 2) const
183  {
184  // EXCEPTION CHECKS
185 
186  if (D == 0)
187  throw exception::DimsInvalid("qpp::Gates::Xd()");
188  // END EXCEPTION CHECKS
189 
190  return Fd(D).inverse() * Zd(D) * Fd(D);
191  }
192 
202  template<typename Derived = Eigen::MatrixXcd>
203  Derived Id(idx D = 2) const
204  {
205  // EXCEPTION CHECKS
206 
207  if (D == 0)
208  throw exception::DimsInvalid("qpp::Gates::Id()");
209  // END EXCEPTION CHECKS
210 
211  return Derived::Identity(D, D);
212  }
213 
229  template<typename Derived>
230  dyn_mat<typename Derived::Scalar> CTRL(const Eigen::MatrixBase<Derived>& A,
231  const std::vector<idx>& ctrl,
232  const std::vector<idx>& subsys,
233  idx N, idx d = 2) const
234  {
235  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
236 
237  // EXCEPTION CHECKS
238 
239  // check matrix zero size
241  throw exception::ZeroSize("qpp::Gates::CTRL()");
242 
243  // check square matrix
245  throw exception::MatrixNotSquare("qpp::Gates::CTRL()");
246 
247  // check lists zero size
248  if (ctrl.size() == 0)
249  throw exception::ZeroSize("qpp::Gates::CTRL()");
250  if (subsys.size() == 0)
251  throw exception::ZeroSize("qpp::Gates::CTRL()");
252 
253  // check out of range
254  if (N == 0)
255  throw exception::OutOfRange("qpp::Gates::CTRL()");
256 
257  // check valid local dimension
258  if (d == 0)
259  throw exception::DimsInvalid("qpp::Gates::CTRL()");
260 
261  // ctrl + gate subsystem vector
262  std::vector<idx> ctrlgate = ctrl;
263  ctrlgate.insert(std::end(ctrlgate), std::begin(subsys),
264  std::end(subsys));
265  std::sort(std::begin(ctrlgate), std::end(ctrlgate));
266 
267  std::vector<idx> dims(N, d); // local dimensions vector
268 
269  // check that ctrl + gate subsystem is valid
270  // with respect to local dimensions
271  if (!internal::check_subsys_match_dims(ctrlgate, dims))
272  throw exception::SubsysMismatchDims("qpp::Gates::CTRL()");
273 
274  // check that subsys list match the dimension of the matrix
275  using Index = typename dyn_mat<typename Derived::Scalar>::Index;
276  if (rA.rows() !=
277  static_cast<Index>(std::llround(std::pow(d, subsys.size()))))
278  throw exception::DimsMismatchMatrix("qpp::Gates::CTRL()");
279  // END EXCEPTION CHECKS
280 
281  // Use static allocation for speed!
282  idx Cdims[maxn];
283  idx midx_row[maxn];
284  idx midx_col[maxn];
285 
286  idx CdimsA[maxn];
287  idx midxA_row[maxn];
288  idx midxA_col[maxn];
289 
290  idx Cdims_bar[maxn];
291  idx Csubsys_bar[maxn];
292  idx midx_bar[maxn];
293 
294  idx Ngate = subsys.size();
295  idx Nctrl = ctrl.size();
296  idx Nsubsys_bar = N - ctrlgate.size();
297  idx D = static_cast<idx>(std::llround(std::pow(d, N)));
298  idx DA = static_cast<idx>(rA.rows());
299  idx Dsubsys_bar = static_cast<idx>(
300  std::llround(std::pow(d, Nsubsys_bar)));
301 
302  // compute the complementary subsystem of ctrlgate w.r.t. dims
303  std::vector<idx> subsys_bar = complement(ctrlgate, N);
304  std::copy(std::begin(subsys_bar), std::end(subsys_bar),
305  std::begin(Csubsys_bar));
306 
307  for (idx k = 0; k < N; ++k)
308  {
309  midx_row[k] = midx_col[k] = 0;
310  Cdims[k] = d;
311  }
312 
313  for (idx k = 0; k < Nsubsys_bar; ++k)
314  {
315  Cdims_bar[k] = d;
316  midx_bar[k] = 0;
317  }
318 
319  for (idx k = 0; k < Ngate; ++k)
320  {
321  midxA_row[k] = midxA_col[k] = 0;
322  CdimsA[k] = d;
323  }
324 
326  typename Derived::Scalar>
327  ::Identity(D, D);
329 
330  // run over the complement indexes
331  for (idx i = 0; i < Dsubsys_bar; ++i)
332  {
333  // get the complement row multi-index
334  internal::n2multiidx(i, Nsubsys_bar, Cdims_bar, midx_bar);
335  for (idx k = 0; k < d; ++k)
336  {
337  Ak = powm(rA, k); // compute rA^k
338  // run over the subsys row multi-index
339  for (idx a = 0; a < DA; ++a)
340  {
341  // get the subsys row multi-index
342  internal::n2multiidx(a, Ngate, CdimsA, midxA_row);
343 
344  // construct the result row multi-index
345 
346  // first the ctrl part (equal for both row and column)
347  for (idx c = 0; c < Nctrl; ++c)
348  midx_row[ctrl[c]] = midx_col[ctrl[c]] = k;
349 
350  // then the complement part (equal for column)
351  for (idx c = 0; c < Nsubsys_bar; ++c)
352  midx_row[Csubsys_bar[c]] = midx_col[Csubsys_bar[c]] =
353  midx_bar[c];
354 
355  // then the subsys part
356  for (idx c = 0; c < Ngate; ++c)
357  midx_row[subsys[c]] = midxA_row[c];
358 
359  // run over the subsys column multi-index
360  for (idx b = 0; b < DA; ++b)
361  {
362  // get the subsys column multi-index
363  internal::n2multiidx(b, Ngate, CdimsA, midxA_col);
364 
365  // construct the result column multi-index
366  for (idx c = 0; c < Ngate; ++c)
367  midx_col[subsys[c]] = midxA_col[c];
368 
369  // finally write the values
370  result(internal::multiidx2n(midx_row, N, Cdims),
371  internal::multiidx2n(midx_col, N, Cdims))
372  = Ak(a, b);
373  }
374  }
375  }
376  }
377 
378  return result;
379  }
380 
396  template<typename Derived>
398  const Eigen::MatrixBase<Derived>& A, idx pos,
399  const std::vector<idx>& dims) const
400  {
401  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
402 
403  // EXCEPTION CHECKS
404 
405  // check zero-size
407  throw exception::ZeroSize("qpp::Gates::expandout()");
408 
409  // check that dims is a valid dimension vector
410  if (!internal::check_dims(dims))
411  throw exception::DimsInvalid("qpp::Gates::expandout()");
412 
413  // check square matrix
415  throw exception::MatrixNotSquare("qpp::Gates::expandout()");
416 
417  // check that position is valid
418  if (pos > dims.size() - 1)
419  throw exception::OutOfRange("qpp::Gates::expandout()");
420 
421  // check that dims[pos] match the dimension of A
422  if (static_cast<idx>(rA.rows()) != dims[pos])
423  throw exception::DimsMismatchMatrix("qpp::Gates::expandout()");
424  // END EXCEPTION CHECKS
425 
426  idx D = std::accumulate(std::begin(dims), std::end(dims),
427  static_cast<idx>(1), std::multiplies<idx>());
429  typename Derived::Scalar>
430  ::Identity(D, D);
431 
432  idx Cdims[maxn];
433  idx midx_row[maxn];
434  idx midx_col[maxn];
435 
436  for (idx k = 0; k < dims.size(); ++k)
437  {
438  midx_row[k] = midx_col[k] = 0;
439  Cdims[k] = dims[k];
440  }
441 
442  // run over the main diagonal multi-indexes
443  for (idx i = 0; i < D; ++i)
444  {
445  // get row multi_index
446  internal::n2multiidx(i, dims.size(), Cdims, midx_row);
447  // get column multi_index (same as row)
448  internal::n2multiidx(i, dims.size(), Cdims, midx_col);
449  // run over the gate row multi-index
450  for (idx a = 0; a < static_cast<idx>(rA.rows());
451  ++a)
452  {
453  // construct the total row multi-index
454  midx_row[pos] = a;
455 
456  // run over the gate column multi-index
457  for (idx b = 0;
458  b < static_cast<idx>(rA.cols()); ++b)
459  {
460  // construct the total column multi-index
461  midx_col[pos] = b;
462 
463  // finally write the values
464  result(internal::multiidx2n(midx_row, dims.size(), Cdims),
465  internal::multiidx2n(midx_col, dims.size(), Cdims))
466  = rA(a, b);
467  }
468  }
469  }
470 
471  return result;
472  }
473 
495  template<typename Derived>
497  const Eigen::MatrixBase<Derived>& A, idx pos,
498  const std::initializer_list<idx>& dims) const
499  {
500  return this->expandout(A, pos, std::vector<idx>(dims));
501  }
502 
519  template<typename Derived>
521  const Eigen::MatrixBase<Derived>& A, idx pos, idx N, idx d = 2)
522  const
523  {
524  // EXCEPTION CHECKS
525 
526  // check zero size
528  throw exception::ZeroSize("qpp::Gates::expandout()");
529 
530  // check valid dims
531  if (d == 0)
532  throw exception::DimsInvalid("qpp::Gates::expandout()");
533  // END EXCEPTION CHECKS
534 
535  std::vector<idx> dims(N, d); // local dimensions vector
536 
537  return this->expandout(A, pos, dims);
538  }
539 }; /* class Gates */
540 
541 } /* namespace qpp */
542 
543 #endif /* CLASSES_GATES_H_ */
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:129
Dimension(s) mismatch matrix size exception.
Definition: exception.h:322
cmat SWAP
SWAP gate.
Definition: gates.h:54
constexpr idx maxn
Maximum number of allowed qubits/qudits (subsystems)
Definition: constants.h:71
bool check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:221
Derived Id(idx D=2) const
Identity gate.
Definition: gates.h:203
~Gates()=default
Default destructor.
cplx omega(idx D)
D-th root of unity.
Definition: constants.h:93
Custom exception.
Definition: exception.h:636
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:77
idx multiidx2n(const idx *const midx, idx numdims, const idx *const dims) noexcept
Definition: util.h:70
dyn_mat< typename Derived::Scalar > expandout(const Eigen::MatrixBase< Derived > &A, idx pos, const std::initializer_list< idx > &dims) const
Expands out.
Definition: gates.h:496
void n2multiidx(idx n, idx numdims, const idx *const dims, idx *result) noexcept
Definition: util.h:40
Singleton policy class, used internally to implement the singleton pattern via CRTP (Curiously recurr...
Definition: singleton.h:76
Subsystems mismatch dimensions exception.
Definition: exception.h:395
const Singleton class that implements most commonly used gates
Definition: gates.h:36
Quantum++ main namespace.
Definition: codes.h:30
cmat Xd(idx D=2) const
Generalized X gate for qudits.
Definition: gates.h:182
Invalid dimension(s) exception.
Definition: exception.h:287
cmat S
S gate.
Definition: gates.h:47
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1811
cmat T
T gate.
Definition: gates.h:48
cmat FRED
Fredkin gate.
Definition: gates.h:58
cmat CNOTba
Controlled-NOT target control gate.
Definition: gates.h:53
dyn_mat< typename Derived::Scalar > expandout(const Eigen::MatrixBase< Derived > &A, idx pos, idx N, idx d=2) const
Expands out.
Definition: gates.h:520
cmat Zd(idx D=2) const
Generalized Z gate for qudits.
Definition: gates.h:129
cmat H
Hadamard gate.
Definition: gates.h:43
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:102
cmat Id2
Identity gate.
Definition: gates.h:42
cmat X
Pauli Sigma-X gate.
Definition: gates.h:44
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:777
cmat Fd(idx D=2) const
Fourier transform gate for qudits.
Definition: gates.h:153
dyn_mat< typename Derived::Scalar > expandout(const Eigen::MatrixBase< Derived > &A, idx pos, const std::vector< idx > &dims) const
Expands out.
Definition: gates.h:397
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:101
constexpr double pi
Definition: constants.h:76
Parameter out of range exception.
Definition: exception.h:567
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:230
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:60
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:142
cmat CNOT
Controlled-NOT control target gate.
Definition: gates.h:51
std::size_t idx
Non-negative integer index.
Definition: types.h:35
cmat Z
Pauli Sigma-Z gate.
Definition: gates.h:46
cmat CZ
Controlled-Phase gate.
Definition: gates.h:52
Matrix is not square exception.
Definition: exception.h:151
cmat Y
Pauli Sigma-Y gate.
Definition: gates.h:45
Gates()
Initializes the gates.
Definition: gates.h:63
cmat TOF
Toffoli gate.
Definition: gates.h:57
Object has zero size exception.
Definition: exception.h:134