Quantum++  v1.0
A modern C++11 quantum computing library
gates.h
Go to the documentation of this file.
1 /*
2  * This file is part of Quantum++.
3  *
4  * MIT License
5  *
6  * Copyright (c) 2013 - 2018 Vlad Gheorghiu (vgheorgh@gmail.com)
7  *
8  * Permission is hereby granted, free of charge, to any person obtaining a copy
9  * of this software and associated documentation files (the "Software"), to deal
10  * in the Software without restriction, including without limitation the rights
11  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
12  * copies of the Software, and to permit persons to whom the Software is
13  * furnished to do so, subject to the following conditions:
14  *
15  * The above copyright notice and this permission notice shall be included in
16  * all copies or substantial portions of the Software.
17  *
18  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
19  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
20  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
21  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
22  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
23  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
24  * SOFTWARE.
25  */
26 
32 #ifndef CLASSES_GATES_H_
33 #define CLASSES_GATES_H_
34 
35 namespace qpp {
40 class Gates final : public internal::Singleton<const Gates> // const Singleton
41 {
42  friend class internal::Singleton<const Gates>;
43 
44  public:
45  // One qubit gates
46  cmat Id2{cmat::Identity(2, 2)};
47  cmat H{cmat::Zero(2, 2)};
48  cmat X{cmat::Zero(2, 2)};
49  cmat Y{cmat::Zero(2, 2)};
50  cmat Z{cmat::Zero(2, 2)};
51  cmat S{cmat::Zero(2, 2)};
52  cmat T{cmat::Zero(2, 2)};
53 
54  // two qubit gates
55  cmat CNOT{cmat::Identity(4, 4)};
56  cmat CZ{cmat::Identity(4, 4)};
57  cmat CNOTba{cmat::Zero(4, 4)};
58  cmat SWAP{cmat::Identity(4, 4)};
59 
60  // three qubit gates
61  cmat TOF{cmat::Identity(8, 8)};
62  cmat FRED{cmat::Identity(8, 8)};
63  private:
67  Gates() {
68  H << 1 / std::sqrt(2.), 1 / std::sqrt(2.), 1 / std::sqrt(2.),
69  -1 / std::sqrt(2.);
70  X << 0, 1, 1, 0;
71  Z << 1, 0, 0, -1;
72  Y << 0, -1_i, 1_i, 0;
73  S << 1, 0, 0, 1_i;
74  T << 1, 0, 0, std::exp(1_i * pi / 4.0);
75  CNOT.block(2, 2, 2, 2) = X;
76  CNOTba(0, 0) = 1;
77  CNOTba(1, 3) = 1;
78  CNOTba(2, 2) = 1;
79  CNOTba(3, 1) = 1;
80  CZ(3, 3) = -1;
81 
82  SWAP.block(1, 1, 2, 2) = X;
83  TOF.block(6, 6, 2, 2) = X;
84  FRED.block(4, 4, 4, 4) = SWAP;
85  }
86 
90  ~Gates() = default;
91 
92  public:
93  // variable gates
94 
95  // one qubit gates
96 
105  cmat Rn(double theta, const std::vector<double>& n) const {
106  // EXCEPTION CHECKS
107 
108  // check 3-dimensional vector
109  if (n.size() != 3)
111  "qpp::Gates::Rn()", "n is not a 3-dimensional vector!");
112  // END EXCEPTION CHECKS
113 
114  cmat result(2, 2);
115  result = std::cos(theta / 2) * Id2 -
116  1_i * std::sin(theta / 2) * (n[0] * X + n[1] * Y + n[2] * Z);
117 
118  return result;
119  }
120 
121  // one quDit gates
122 
132  cmat Zd(idx D = 2) const {
133  // EXCEPTION CHECKS
134 
135  if (D == 0)
136  throw exception::DimsInvalid("qpp::Gates::Zd()");
137  // END EXCEPTION CHECKS
138 
139  cmat result = cmat::Zero(D, D);
140  for (idx i = 0; i < D; ++i)
141  result(i, i) = std::pow(omega(D), static_cast<double>(i));
142 
143  return result;
144  }
145 
156  cmat Fd(idx D = 2) const {
157  // EXCEPTION CHECKS
158 
159  if (D == 0)
160  throw exception::DimsInvalid("qpp::Gates::Fd()");
161  // END EXCEPTION CHECKS
162 
163  cmat result(D, D);
164 
165 #ifdef WITH_OPENMP_
166 #pragma omp parallel for collapse(2)
167 #endif // WITH_OPENMP_
168  // column major order for speed
169  for (idx j = 0; j < D; ++j)
170  for (idx i = 0; i < D; ++i)
171  result(i, j) = 1 / std::sqrt(D) *
172  std::pow(omega(D), static_cast<double>(i * j));
173 
174  return result;
175  }
176 
186  cmat Xd(idx D = 2) const {
187  // EXCEPTION CHECKS
188 
189  if (D == 0)
190  throw exception::DimsInvalid("qpp::Gates::Xd()");
191  // END EXCEPTION CHECKS
192 
193  return Fd(D).inverse() * Zd(D) * Fd(D);
194  }
195 
205  template <typename Derived = Eigen::MatrixXcd>
206  Derived Id(idx D = 2) const {
207  // EXCEPTION CHECKS
208 
209  if (D == 0)
210  throw exception::DimsInvalid("qpp::Gates::Id()");
211  // END EXCEPTION CHECKS
212 
213  return Derived::Identity(D, D);
214  }
215 
231  template <typename Derived>
233  CTRL(const Eigen::MatrixBase<Derived>& A, const std::vector<idx>& ctrl,
234  const std::vector<idx>& subsys, idx N, idx d = 2) const {
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 =
300  static_cast<idx>(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  midx_row[k] = midx_col[k] = 0;
309  Cdims[k] = d;
310  }
311 
312  for (idx k = 0; k < Nsubsys_bar; ++k) {
313  Cdims_bar[k] = d;
314  midx_bar[k] = 0;
315  }
316 
317  for (idx k = 0; k < Ngate; ++k) {
318  midxA_row[k] = midxA_col[k] = 0;
319  CdimsA[k] = d;
320  }
321 
325 
326  // run over the complement indexes
327  for (idx i = 0; i < Dsubsys_bar; ++i) {
328  // get the complement row multi-index
329  internal::n2multiidx(i, Nsubsys_bar, Cdims_bar, midx_bar);
330  for (idx k = 0; k < d; ++k) {
331  Ak = powm(rA, k); // compute rA^k
332  // run over the subsys row multi-index
333  for (idx a = 0; a < DA; ++a) {
334  // get the subsys row multi-index
335  internal::n2multiidx(a, Ngate, CdimsA, midxA_row);
336 
337  // construct the result row multi-index
338 
339  // first the ctrl part (equal for both row and column)
340  for (idx c = 0; c < Nctrl; ++c)
341  midx_row[ctrl[c]] = midx_col[ctrl[c]] = k;
342 
343  // then the complement part (equal for column)
344  for (idx c = 0; c < Nsubsys_bar; ++c)
345  midx_row[Csubsys_bar[c]] = midx_col[Csubsys_bar[c]] =
346  midx_bar[c];
347 
348  // then the subsys part
349  for (idx c = 0; c < Ngate; ++c)
350  midx_row[subsys[c]] = midxA_row[c];
351 
352  // run over the subsys column multi-index
353  for (idx b = 0; b < DA; ++b) {
354  // get the subsys column multi-index
355  internal::n2multiidx(b, Ngate, CdimsA, midxA_col);
356 
357  // construct the result column multi-index
358  for (idx c = 0; c < Ngate; ++c)
359  midx_col[subsys[c]] = midxA_col[c];
360 
361  // finally write the values
362  result(internal::multiidx2n(midx_row, N, Cdims),
363  internal::multiidx2n(midx_col, N, Cdims)) =
364  Ak(a, b);
365  }
366  }
367  }
368  }
369 
370  return result;
371  }
372 
388  template <typename Derived>
390  expandout(const Eigen::MatrixBase<Derived>& A, idx pos,
391  const std::vector<idx>& dims) const {
392  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
393 
394  // EXCEPTION CHECKS
395 
396  // check zero-size
398  throw exception::ZeroSize("qpp::Gates::expandout()");
399 
400  // check that dims is a valid dimension vector
401  if (!internal::check_dims(dims))
402  throw exception::DimsInvalid("qpp::Gates::expandout()");
403 
404  // check square matrix
406  throw exception::MatrixNotSquare("qpp::Gates::expandout()");
407 
408  // check that position is valid
409  if (pos > dims.size() - 1)
410  throw exception::OutOfRange("qpp::Gates::expandout()");
411 
412  // check that dims[pos] match the dimension of A
413  if (static_cast<idx>(rA.rows()) != dims[pos])
414  throw exception::DimsMismatchMatrix("qpp::Gates::expandout()");
415  // END EXCEPTION CHECKS
416 
417  idx D = std::accumulate(std::begin(dims), std::end(dims),
418  static_cast<idx>(1), std::multiplies<idx>());
421 
422  idx Cdims[maxn];
423  idx midx_row[maxn];
424  idx midx_col[maxn];
425 
426  for (idx k = 0; k < dims.size(); ++k) {
427  midx_row[k] = midx_col[k] = 0;
428  Cdims[k] = dims[k];
429  }
430 
431  // run over the main diagonal multi-indexes
432  for (idx i = 0; i < D; ++i) {
433  // get row multi_index
434  internal::n2multiidx(i, dims.size(), Cdims, midx_row);
435  // get column multi_index (same as row)
436  internal::n2multiidx(i, dims.size(), Cdims, midx_col);
437  // run over the gate row multi-index
438  for (idx a = 0; a < static_cast<idx>(rA.rows()); ++a) {
439  // construct the total row multi-index
440  midx_row[pos] = a;
441 
442  // run over the gate column multi-index
443  for (idx b = 0; b < static_cast<idx>(rA.cols()); ++b) {
444  // construct the total column multi-index
445  midx_col[pos] = b;
446 
447  // finally write the values
448  result(internal::multiidx2n(midx_row, dims.size(), Cdims),
449  internal::multiidx2n(midx_col, dims.size(), Cdims)) =
450  rA(a, b);
451  }
452  }
453  }
454 
455  return result;
456  }
457 
479  template <typename Derived>
481  expandout(const Eigen::MatrixBase<Derived>& A, idx pos,
482  const std::initializer_list<idx>& dims) const {
483  return this->expandout(A, pos, std::vector<idx>(dims));
484  }
485 
502  template <typename Derived>
504  expandout(const Eigen::MatrixBase<Derived>& A, idx pos, idx N,
505  idx d = 2) const {
506  // EXCEPTION CHECKS
507 
508  // check zero size
510  throw exception::ZeroSize("qpp::Gates::expandout()");
511 
512  // check valid dims
513  if (d == 0)
514  throw exception::DimsInvalid("qpp::Gates::expandout()");
515  // END EXCEPTION CHECKS
516 
517  std::vector<idx> dims(N, d); // local dimensions vector
518 
519  return this->expandout(A, pos, dims);
520  }
521 }; /* class Gates */
522 
523 } /* namespace qpp */
524 
525 #endif /* CLASSES_GATES_H_ */
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:123
Dimension(s) mismatch matrix size exception.
Definition: exception.h:300
cmat SWAP
SWAP gate.
Definition: gates.h:58
constexpr idx maxn
Maximum number of allowed qubits/qudits (subsystems)
Definition: constants.h:75
bool check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:210
Derived Id(idx D=2) const
Identity gate.
Definition: gates.h:206
~Gates()=default
Default destructor.
cplx omega(idx D)
D-th root of unity.
Definition: constants.h:97
Custom exception.
Definition: exception.h:571
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:81
idx multiidx2n(const idx *const midx, idx numdims, const idx *const dims) noexcept
Definition: util.h:71
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:481
void n2multiidx(idx n, idx numdims, const idx *const dims, idx *result) noexcept
Definition: util.h:43
Singleton policy class, used internally to implement the singleton pattern via CRTP (Curiously recurr...
Definition: singleton.h:80
Subsystems mismatch dimensions exception.
Definition: exception.h:365
const Singleton class that implements most commonly used gates
Definition: gates.h:40
Quantum++ main namespace.
Definition: codes.h:35
cmat Xd(idx D=2) const
Generalized X gate for qudits.
Definition: gates.h:186
Invalid dimension(s) exception.
Definition: exception.h:269
cmat S
S gate.
Definition: gates.h:51
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1733
cmat T
T gate.
Definition: gates.h:52
cmat FRED
Fredkin gate.
Definition: gates.h:62
cmat CNOTba
Controlled-NOT target control gate.
Definition: gates.h:57
dyn_mat< typename Derived::Scalar > expandout(const Eigen::MatrixBase< Derived > &A, idx pos, idx N, idx d=2) const
Expands out.
Definition: gates.h:504
cmat Zd(idx D=2) const
Generalized Z gate for qudits.
Definition: gates.h:132
cmat H
Hadamard gate.
Definition: gates.h:47
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:105
cmat Id2
Identity gate.
Definition: gates.h:46
cmat X
Pauli Sigma-X gate.
Definition: gates.h:48
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:751
cmat Fd(idx D=2) const
Fourier transform gate for qudits.
Definition: gates.h:156
dyn_mat< typename Derived::Scalar > expandout(const Eigen::MatrixBase< Derived > &A, idx pos, const std::vector< idx > &dims) const
Expands out.
Definition: gates.h:390
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:99
constexpr double pi
Definition: constants.h:80
Parameter out of range exception.
Definition: exception.h:515
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:233
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:64
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:134
cmat CNOT
Controlled-NOT control target gate.
Definition: gates.h:55
std::size_t idx
Non-negative integer index.
Definition: types.h:39
cmat Z
Pauli Sigma-Z gate.
Definition: gates.h:50
cmat CZ
Controlled-Phase gate.
Definition: gates.h:56
Matrix is not square exception.
Definition: exception.h:149
cmat Y
Pauli Sigma-Y gate.
Definition: gates.h:49
Gates()
Initializes the gates.
Definition: gates.h:67
cmat TOF
Toffoli gate.
Definition: gates.h:61
Object has zero size exception.
Definition: exception.h:134