Quantum++  v0.8.8.2
C++11 quantum computing library
gates.h
Go to the documentation of this file.
1 /*
2  * Quantum++
3  *
4  * Copyright (c) 2013 - 2016 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 
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("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) const
130  {
131  // EXCEPTION CHECKS
132 
133  if ( D == 0 )
134  throw Exception("qpp::Gates::Zd()", Exception::Type::DIMS_INVALID);
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) const
154  {
155  // EXCEPTION CHECKS
156 
157  if ( D == 0 )
158  throw Exception("qpp::Gates::Fd()", Exception::Type::DIMS_INVALID);
159  // END EXCEPTION CHECKS
160 
161  cmat result(D, D);
162 
163 #ifdef _WITH_OPENMP_
164 #pragma omp parallel for collapse(2)
165 #endif
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(static_cast<double>(D))
169  * std::pow(omega(D), i * j);
170 
171  return result;
172  }
173 
182  cmat Xd(idx D) const
183  {
184  // EXCEPTION CHECKS
185 
186  if ( D == 0 )
187  throw Exception("qpp::Gates::Xd()", Exception::Type::DIMS_INVALID);
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) const
204  {
205  // EXCEPTION CHECKS
206 
207  if ( D == 0 )
208  throw Exception("qpp::Gates::Id()", Exception::Type::DIMS_INVALID);
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("qpp::Gates::CTRL()", Exception::Type::ZERO_SIZE);
242 
243  // check square matrix
244  if ( !internal::_check_square_mat(rA))
245  throw Exception("qpp::Gates::CTRL()",
247 
248  // check lists zero size
249  if ( ctrl.size() == 0 )
250  throw Exception("qpp::Gates::CTRL()", Exception::Type::ZERO_SIZE);
251  if ( subsys.size() == 0 )
252  throw Exception("qpp::Gates::CTRL()", Exception::Type::ZERO_SIZE);
253 
254  // check out of range
255  if ( n == 0 )
256  throw Exception("qpp::Gates::CTRL()",
258 
259  // check valid local dimension
260  if ( d == 0 )
261  throw Exception("qpp::Gates::CTRL()",
263 
264  // ctrl + gate subsystem vector
265  std::vector<idx> ctrlgate = ctrl;
266  ctrlgate.insert(std::end(ctrlgate), std::begin(subsys),
267  std::end(subsys));
268  std::sort(std::begin(ctrlgate), std::end(ctrlgate));
269 
270  std::vector<idx> dims(n, d); // local dimensions vector
271 
272  // check that ctrl + gate subsystem is valid
273  // with respect to local dimensions
274  if ( !internal::_check_subsys_match_dims(ctrlgate, dims))
275  throw Exception("qpp::Gates::CTRL()",
277 
278  // check that subsys list match the dimension of the matrix
279  if ( rA.rows() != std::llround(std::pow(d, subsys.size())))
280  throw Exception("qpp::Gates::CTRL()",
282  // END EXCEPTION CHECKS
283 
284  // Use static allocation for speed!
285  idx Cdims[maxn];
286  idx midx_row[maxn];
287  idx midx_col[maxn];
288 
289  idx CdimsA[maxn];
290  idx midxA_row[maxn];
291  idx midxA_col[maxn];
292 
293  idx Cdims_bar[maxn];
294  idx Csubsys_bar[maxn];
295  idx midx_bar[maxn];
296 
297  idx ngate = subsys.size();
298  idx nctrl = ctrl.size();
299  idx nsubsys_bar = n - ctrlgate.size();
300  idx D = static_cast<idx>(std::llround(std::pow(d, n)));
301  idx DA = static_cast<idx>(rA.rows());
302  idx Dsubsys_bar = static_cast<idx>(
303  std::llround(std::pow(d, nsubsys_bar)));
304 
305  // compute the complementary subsystem of ctrlgate w.r.t. dims
306  std::vector<idx> subsys_bar = complement(ctrlgate, n);
307  std::copy(std::begin(subsys_bar), std::end(subsys_bar),
308  std::begin(Csubsys_bar));
309 
310  for ( idx k = 0; k < n; ++k )
311  {
312  midx_row[k] = midx_col[k] = 0;
313  Cdims[k] = d;
314  }
315 
316  for ( idx k = 0; k < nsubsys_bar; ++k )
317  {
318  Cdims_bar[k] = d;
319  midx_bar[k] = 0;
320  }
321 
322  for ( idx k = 0; k < ngate; ++k )
323  {
324  midxA_row[k] = midxA_col[k] = 0;
325  CdimsA[k] = d;
326  }
327 
329  typename Derived::Scalar>
330  ::Identity(D, D);
332 
333  // run over the complement indexes
334  for ( idx i = 0; i < Dsubsys_bar; ++i )
335  {
336  // get the complement row multi-index
337  internal::_n2multiidx(i, nsubsys_bar, Cdims_bar, midx_bar);
338  for ( idx k = 0; k < d; ++k )
339  {
340  Ak = powm(rA, k); // compute rA^k
341  // run over the subsys row multi-index
342  for ( idx a = 0; a < DA; ++a )
343  {
344  // get the subsys row multi-index
345  internal::_n2multiidx(a, ngate, CdimsA, midxA_row);
346 
347  // construct the result row multi-index
348 
349  // first the ctrl part (equal for both row and column)
350  for ( idx c = 0; c < nctrl; ++c )
351  midx_row[ctrl[c]] = midx_col[ctrl[c]] = k;
352 
353  // then the complement part (equal for column)
354  for ( idx c = 0; c < nsubsys_bar; ++c )
355  midx_row[Csubsys_bar[c]] = midx_col[Csubsys_bar[c]] =
356  midx_bar[c];
357 
358  // then the subsys part
359  for ( idx c = 0; c < ngate; ++c )
360  midx_row[subsys[c]] = midxA_row[c];
361 
362  // run over the subsys column multi-index
363  for ( idx b = 0; b < DA; ++b )
364  {
365  // get the subsys column multi-index
366  internal::_n2multiidx(b, ngate, CdimsA, midxA_col);
367 
368  // construct the result column multi-index
369  for ( idx c = 0; c < ngate; ++c )
370  midx_col[subsys[c]] = midxA_col[c];
371 
372  // finally write the values
373  result(internal::_multiidx2n(midx_row, n, Cdims),
374  internal::_multiidx2n(midx_col, n, Cdims))
375  = Ak(a, b);
376  }
377  }
378  }
379  }
380 
381  return result;
382  }
383 
399  template<typename Derived>
401  const Eigen::MatrixBase<Derived>& A, idx pos,
402  const std::vector<idx>& dims) const
403  {
404  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
405 
406  // EXCEPTION CHECKS
407 
408  // check zero-size
410  throw Exception("qpp::Gates::expandout()",
412 
413  // check that dims is a valid dimension vector
414  if ( !internal::_check_dims(dims))
415  throw Exception("qpp::Gates::expandout()",
417 
418  // check square matrix
419  if ( !internal::_check_square_mat(rA))
420  throw Exception("qpp::Gates::expandout()",
422 
423  // check that position is valid
424  if ( pos > dims.size() - 1 )
425  throw Exception("qpp::Gates::expandout()",
427 
428  // check that dims[pos] match the dimension of A
429  if ( static_cast<idx>(rA.rows()) != dims[pos] )
430  throw Exception("qpp::Gates::expandout()",
432  // END EXCEPTION CHECKS
433 
434  idx D = std::accumulate(std::begin(dims), std::end(dims),
435  static_cast<idx>(1), std::multiplies<idx>());
437  typename Derived::Scalar>
438  ::Identity(D, D);
439 
440  idx Cdims[maxn];
441  idx midx_row[maxn];
442  idx midx_col[maxn];
443 
444  for ( idx k = 0; k < dims.size(); ++k )
445  {
446  midx_row[k] = midx_col[k] = 0;
447  Cdims[k] = dims[k];
448  }
449 
450  // run over the main diagonal multi-indexes
451  for ( idx i = 0; i < D; ++i )
452  {
453  // get row multi_index
454  internal::_n2multiidx(i, dims.size(), Cdims, midx_row);
455  // get column multi_index (same as row)
456  internal::_n2multiidx(i, dims.size(), Cdims, midx_col);
457  // run over the gate row multi-index
458  for ( idx a = 0; a < static_cast<idx>(rA.rows());
459  ++a )
460  {
461  // construct the total row multi-index
462  midx_row[pos] = a;
463 
464  // run over the gate column multi-index
465  for ( idx b = 0;
466  b < static_cast<idx>(rA.cols()); ++b )
467  {
468  // construct the total column multi-index
469  midx_col[pos] = b;
470 
471  // finally write the values
472  result(internal::_multiidx2n(midx_row, dims.size(), Cdims),
473  internal::_multiidx2n(midx_col, dims.size(), Cdims))
474  = rA(a, b);
475  }
476  }
477  }
478 
479  return result;
480  }
481 
482 }; /* class Gates */
483 
484 } /* namespace qpp */
485 
486 #endif /* CLASSES_GATES_H_ */
dyn_mat< typename Derived::Scalar > expandout(const Eigen::MatrixBase< Derived > &A, idx pos, const std::vector< idx > &dims) const
Expands out.
Definition: gates.h:400
cmat SWAP
SWAP gate.
Definition: gates.h:55
constexpr idx maxn
Maximum number of allowed qu(d)its (subsystems)
Definition: constants.h:74
bool _check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:183
~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:83
cmat Zd(idx D) const
Generalized Z gate for qudits.
Definition: gates.h:129
Singleton policy class, used internally to implement the singleton pattern via CRTP (Curiously recurr...
Definition: singleton.h:77
const Singleton class that implements most commonly used gates
Definition: gates.h:37
Quantum++ main namespace.
Definition: codes.h:30
bool _check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:84
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:1818
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:153
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
idx _multiidx2n(const idx *midx, idx numdims, const idx *dims) noexcept
Definition: util.h:61
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)
Matrix power.
Definition: functions.h:783
bool _check_nonzero_size(const T &x) noexcept
Definition: util.h:112
Derived Id(idx D) const
Identity gate.
Definition: gates.h:203
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
constexpr double pi
Definition: constants.h:79
cmat Rn(double theta, const std::vector< double > &n) const
Rotation of theta about the 3-dimensional real unit vector n.
Definition: gates.h:102
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:66
void _n2multiidx(idx n, idx numdims, const idx *dims, idx *result) noexcept
Definition: util.h:48
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
bool _check_dims(const std::vector< idx > &dims)
Definition: util.h:125
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:182