Sequential Quantum Gate Decomposer  v1.9.3
Powerful decomposition of general unitarias into one- and two-qubit gates gates
Gate.cpp
Go to the documentation of this file.
1 /*
2 Created on Fri Jun 26 14:13:26 2020
3 Copyright 2020 Peter Rakyta, Ph.D.
4 
5 Licensed under the Apache License, Version 2.0 (the "License");
6 you may not use this file except in compliance with the License.
7 You may obtain a copy of the License at
8 
9  http://www.apache.org/licenses/LICENSE-2.0
10 
11 Unless required by applicable law or agreed to in writing, software
12 distributed under the License is distributed on an "AS IS" BASIS,
13 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 See the License for the specific language governing permissions and
15 limitations under the License.
16 
17 @author: Peter Rakyta, Ph.D.
18 */
24 #include "Gate.h"
25 #include "common.h"
26 
27 #ifdef USE_AVX
30 #endif
31 
32 #include "apply_kernel_to_input.h"
34 
40 
41  // A string labeling the gate operation
42  name = "Gate";
43  // number of qubits spanning the matrix of the operation
44  qbit_num = -1;
45  // The size N of the NxN matrix associated with the operations.
46  matrix_size = -1;
47  // The type of the operation (see enumeration gate_type)
49  // The index of the qubit on which the operation acts (target_qbit >= 0)
50  target_qbit = -1;
51  // The index of the qubit which acts as a control qubit (control_qbit >= 0) in controlled operations
52  control_qbit = -1;
53  // the number of free parameters of the operation
54  parameter_num = 0;
55  // the index in the parameter array (corrensponding to the encapsulated circuit) where the gate parameters begin (if gates are placed into a circuit a single parameter array is used to execute the whole circuit)
57 }
58 
59 
60 
66 Gate::Gate(int qbit_num_in) {
67 
68  if (qbit_num_in > 30) {
69  std::string err("Gate::Gate: Number of qubits supported up to 30");
70  throw err;
71  }
72 
73  // A string labeling the gate operation
74  name = "Gate";
75  // number of qubits spanning the matrix of the operation
76  qbit_num = qbit_num_in;
77  // the size of the matrix
79  // A string describing the type of the operation
81  // The index of the qubit on which the operation acts (target_qbit >= 0)
82  target_qbit = -1;
83  // The index of the qubit which acts as a control qubit (control_qbit >= 0) in controlled operations
84  control_qbit = -1;
85  // The number of parameters
86  parameter_num = 0;
87  // the index in the parameter array (corrensponding to the encapsulated circuit) where the gate parameters begin (if gates are placed into a circuit a single parameter array is used to execute the whole circuit)
89 }
90 
91 
96 }
97 
102 void Gate::set_qbit_num( int qbit_num_in ) {
103 
104  if (qbit_num_in > 30) {
105  std::string err("Gate::set_qbit_num: Number of qubits supported up to 30");
106  throw err;
107  }
108 
109  if ( qbit_num_in <= target_qbit || qbit_num_in <= control_qbit ) {
110  std::string err("Gate::set_qbit_num: The number of qbits is too small, conflicting with either target_qbit os control_qbit");
111  throw err;
112  }
113 
114 
115  // setting the number of qubits
116  qbit_num = qbit_num_in;
117 
118  // update the size of the matrix
120 
121 }
122 
123 
128 Matrix
130 
131  return matrix_alloc;
132 }
133 
134 
140 Matrix
141 Gate::get_matrix(int parallel) {
142 
143  std::string err("Gate::get_matrix: Unimplemented function");
144  throw err;
145 
146 }
147 
148 
154 Matrix
156 
157  std::string err("Gate::get_matrix: Unimplemented function");
158  throw err;
159 
160 }
161 
162 
163 
170 Matrix
172 
173  std::string err("Gate::get_matrix: Unimplemented function");
174  throw err;
175 
176 }
177 
178 
179 
180 
186 void
187 Gate::apply_to_list( std::vector<Matrix>& inputs, int parallel ) {
188 
189  int work_batch = 1;
190  if ( parallel == 0 ) {
191  work_batch = inputs.size();
192  }
193  else {
194  work_batch = 1;
195  }
196 
197 
198  tbb::parallel_for( tbb::blocked_range<int>(0,inputs.size(),work_batch), [&](tbb::blocked_range<int> r) {
199  for (int idx=r.begin(); idx<r.end(); ++idx) {
200 
201  Matrix* input = &inputs[idx];
202 
203  apply_to( *input, parallel );
204 
205  }
206 
207  });
208 
209 
210 
211 
212 }
213 
214 
215 
222 void
223 Gate::apply_to_list( Matrix_real& parameters_mtx, std::vector<Matrix>& inputs, int parallel ) {
224 
225  return;
226 
227 }
228 
229 
230 
236 void
237 Gate::apply_to( Matrix& input, int parallel ) {
238 
239  if (input.rows != matrix_size ) {
240  std::string err("Gate::apply_to: Wrong matrix size in Gate gate apply.");
241  throw err;
242  }
243 
244  Matrix ret = dot(matrix_alloc, input);
245  memcpy( input.get_data(), ret.get_data(), ret.size()*sizeof(QGD_Complex16) );
246  //input = ret;
247 }
248 
249 
256 void
257 Gate::apply_to( Matrix_real& parameter_mtx, Matrix& input, int parallel ) {
258 
259  return;
260 }
261 
262 
263 
270 std::vector<Matrix>
271 Gate::apply_derivate_to( Matrix_real& parameters_mtx_in, Matrix& input, int parallel ) {
272 
273  std::vector<Matrix> ret;
274  return ret;
275 
276 }
277 
278 
279 
284 void
286 
287  Matrix ret = dot(input, matrix_alloc);
288  memcpy( input.get_data(), ret.get_data(), ret.size()*sizeof(QGD_Complex16) );
289 
290 }
291 
292 
298 void
300  matrix_alloc = input;
301 }
302 
303 
308 void Gate::set_control_qbit(int control_qbit_in){
309 
310  if ( control_qbit_in >= qbit_num ) {
311  std::string err("Gate::set_target_qbit: Wrong value of the control qbit: of out of the range given by qbit_num");
312  throw err;
313  }
314 
315 
316  control_qbit = control_qbit_in;
317 }
318 
319 
324 void Gate::set_target_qbit(int target_qbit_in){
325 
326  if ( target_qbit_in >= qbit_num ) {
327  std::string err("Gate::set_target_qbit: Wrong value of the target qbit: out of the range given by qbit_num");
328  throw err;
329  }
330 
331 
332  target_qbit = target_qbit_in;
333 }
334 
339 void Gate::reorder_qubits( std::vector<int> qbit_list ) {
340 
341  // check the number of qubits
342  if ((int)qbit_list.size() != qbit_num ) {
343  std::string err("Gate::reorder_qubits: Wrong number of qubits.");
344  throw err;
345  }
346 
347 
348  int control_qbit_new = control_qbit;
349  int target_qbit_new = target_qbit;
350 
351  // setting the new value for the target qubit
352  for (int idx=0; idx<qbit_num; idx++) {
353  if (target_qbit == qbit_list[idx]) {
354  target_qbit_new = qbit_num-1-idx;
355  }
356  if (control_qbit == qbit_list[idx]) {
357  control_qbit_new = qbit_num-1-idx;
358  }
359  }
360 
361  control_qbit = control_qbit_new;
362  target_qbit = target_qbit_new;
363 }
364 
365 
371  return target_qbit;
372 }
373 
379  return control_qbit;
380 }
381 
386 std::vector<int> Gate::get_involved_qubits() {
387 
388  std::vector<int> involved_qbits;
389 
390  if( target_qbit != -1 ) {
391  involved_qbits.push_back( target_qbit );
392  }
393 
394  if( control_qbit != -1 ) {
395  involved_qbits.push_back( control_qbit );
396  }
397 
398 
399  return involved_qbits;
400 
401 
402 }
403 
404 
409 void Gate::add_parent( Gate* parent ) {
410 
411  // check if parent already present in th elist of parents
412  if ( std::count(parents.begin(), parents.end(), parent) > 0 ) {
413  return;
414  }
415 
416  parents.push_back( parent );
417 
418 }
419 
420 
421 
426 void Gate::add_child( Gate* child ) {
427 
428  // check if parent already present in th elist of parents
429  if ( std::count(children.begin(), children.end(), child) > 0 ) {
430  return;
431  }
432 
433  children.push_back( child );
434 
435 }
436 
437 
442 
443  children.clear();
444 
445 }
446 
447 
452 
453  parents.clear();
454 
455 }
456 
457 
458 
463 std::vector<Gate*> Gate::get_parents() {
464 
465  return parents;
466 
467 }
468 
469 
474 std::vector<Gate*> Gate::get_children() {
475 
476  return children;
477 
478 }
479 
480 
481 
487  return parameter_num;
488 }
489 
490 
496  return type;
497 }
498 
499 
505  return qbit_num;
506 }
507 
508 
514 
515  Gate* ret = new Gate( qbit_num );
516  ret->set_matrix( matrix_alloc );
517 
519  ret->set_parents( parents );
520  ret->set_children( children );
521 
522  return ret;
523 
524 }
525 
526 
527 
536 void
537 Gate::apply_kernel_to(Matrix& u3_1qbit, Matrix& input, bool deriv, int parallel) {
538 
539 #ifdef USE_AVX
540 
541  // apply kernel on state vector
542  if ( input.cols == 1 && (qbit_num<14 || !parallel) ) {
544  return;
545  }
546  else if ( input.cols == 1 ) {
547  if ( parallel == 1 ) {
549  }
550  else if ( parallel == 2 ) {
552  }
553  else {
554  std::string err("Gate::apply_kernel_to: the argument parallel should be either 0,1 or 2. Set 0 for sequential execution (default), 1 for parallel execution with OpenMP and 2 for parallel with TBB");
555  throw err;
556  }
557  return;
558  }
559 
560 
561 
562  // unitary transform kernels
563  if ( qbit_num < 4 ) {
565  return;
566  }
567  else if ( qbit_num < 10 || !parallel) {
569  return;
570  }
571  else {
573  return;
574  }
575 
576 
577 #else
578 
579  // apply kernel on state vector
580  if ( input.cols == 1 && (qbit_num<10 || !parallel) ) {
582  return;
583  }
584  else if ( input.cols == 1 ) {
586  return;
587  }
588 
589 
590  // apply kernel on unitary
591  apply_kernel_to_input(u3_1qbit, input, deriv, target_qbit, control_qbit, matrix_size);
592 
593 
594 
595 
596 
597 #endif // USE_AVX
598 
599 
600 }
601 
602 
603 
604 
605 
612 void
614 
615 
616  int index_step_target = 1 << target_qbit;
617  int current_idx = 0;
618  int current_idx_pair = current_idx+index_step_target;
619 
620 //std::cout << "target qbit: " << target_qbit << std::endl;
621 
622  while ( current_idx_pair < input.cols ) {
623 
624  for(int idx=0; idx<index_step_target; idx++) {
625  //tbb::parallel_for(0, index_step_target, 1, [&](int idx) {
626 
627  int current_idx_loc = current_idx + idx;
628  int current_idx_pair_loc = current_idx_pair + idx;
629 
630  // determine the action according to the state of the control qubit
631  if ( control_qbit<0 || ((current_idx_loc >> control_qbit) & 1) ) {
632 
633  for ( int row_idx=0; row_idx<matrix_size; row_idx++) {
634 
635  int row_offset = row_idx*input.stride;
636 
637 
638  int index = row_offset+current_idx_loc;
639  int index_pair = row_offset+current_idx_pair_loc;
640 
641  QGD_Complex16 element = input[index];
642  QGD_Complex16 element_pair = input[index_pair];
643 
644  QGD_Complex16 tmp1 = mult(u3_1qbit[0], element);
645  QGD_Complex16 tmp2 = mult(u3_1qbit[2], element_pair);
646  input[index].real = tmp1.real + tmp2.real;
647  input[index].imag = tmp1.imag + tmp2.imag;
648 
649  tmp1 = mult(u3_1qbit[1], element);
650  tmp2 = mult(u3_1qbit[3], element_pair);
651  input[index_pair].real = tmp1.real + tmp2.real;
652  input[index_pair].imag = tmp1.imag + tmp2.imag;
653 
654  }
655 
656  }
657  else {
658  // leave the state as it is
659  continue;
660  }
661 
662 
663 //std::cout << current_idx_target << " " << current_idx_target_pair << std::endl;
664 
665 
666  //});
667  }
668 
669 
670  current_idx = current_idx + (index_step_target << 1);
671  current_idx_pair = current_idx_pair + (index_step_target << 1);
672 
673 
674  }
675 
676 
677 }
678 
679 #ifdef _WIN32
680 void sincos(double x, double *s, double *c)
681 {
682  *s = sin(x), *c = cos(x);
683 }
684 #elif defined(__APPLE__)
685 #define sincos __sincos
686 #endif
687 
695 Matrix Gate::calc_one_qubit_u3(double ThetaOver2, double Phi, double Lambda ) {
696 
697  Matrix u3_1qbit = Matrix(2,2);
698 #ifdef DEBUG
699  if (isnan(ThetaOver2)) {
700  std::stringstream sstream;
701  sstream << "Matrix U3::calc_one_qubit_u3: ThetaOver2 is NaN." << std::endl;
702  print(sstream, 1);
703  }
704  if (isnan(Phi)) {
705  std::stringstream sstream;
706  sstream << "Matrix U3::calc_one_qubit_u3: Phi is NaN." << std::endl;
707  print(sstream, 1);
708  }
709  if (isnan(Lambda)) {
710  std::stringstream sstream;
711  sstream << "Matrix U3::calc_one_qubit_u3: Lambda is NaN." << std::endl;
712  print(sstream, 1);
713  }
714 #endif // DEBUG
715 
716  double cos_theta = 1.0, sin_theta = 0.0;
717  double cos_phi = 1.0, sin_phi = 0.0;
718  double cos_lambda = 1.0, sin_lambda = 0.0;
719 
720  if (ThetaOver2!=0.0) sincos(ThetaOver2, &sin_theta, &cos_theta);
721  if (Phi!=0.0) sincos(Phi, &sin_phi, &cos_phi);
722  if (Lambda!=0.0) sincos(Lambda, &sin_lambda, &cos_lambda);
723 
724  // the 1,1 element
725  u3_1qbit[0].real = cos_theta;
726  u3_1qbit[0].imag = 0;
727  // the 1,2 element
728  u3_1qbit[1].real = -cos_lambda*sin_theta;
729  u3_1qbit[1].imag = -sin_lambda*sin_theta;
730  // the 2,1 element
731  u3_1qbit[2].real = cos_phi*sin_theta;
732  u3_1qbit[2].imag = sin_phi*sin_theta;
733  // the 2,2 element
734  //cos(a+b)=cos(a)cos(b)-sin(a)sin(b)
735  //sin(a+b)=sin(a)cos(b)+cos(a)sin(b)
736  u3_1qbit[3].real = (cos_phi*cos_lambda-sin_phi*sin_lambda)*cos_theta;
737  u3_1qbit[3].imag = (sin_phi*cos_lambda+cos_phi*sin_lambda)*cos_theta;
738  //u3_1qbit[3].real = cos(Phi+Lambda)*cos_theta;
739  //u3_1qbit[3].imag = sin(Phi+Lambda)*cos_theta;
740 
741 
742  return u3_1qbit;
743 
744 }
745 
751 
752  std::string err("Gate::calc_one_qubit_u3: Unimplemented abstract function");
753  throw err;
754 
755  Matrix u3_1qbit = Matrix(2,2);
756  return u3_1qbit;
757 
758 }
759 
766 void
767 Gate::parameters_for_calc_one_qubit(double& ThetaOver2, double& Phi, double& Lambda ) {
768 
769  return;
770 
771 }
772 
773 
778 void
780 
781  parameter_start_idx = start_idx;
782 
783 }
784 
789 void
790 Gate::set_parents( std::vector<Gate*>& parents_ ) {
791 
792  parents = parents_;
793 
794 }
795 
796 
801 void
802 Gate::set_children( std::vector<Gate*>& children_ ) {
803 
804  children = children_;
805 
806 }
807 
808 
813 int
815 
816  return parameter_start_idx;
817 
818 }
819 
820 
821 
828 
829  return Matrix_real(0,0);
830 
831 }
832 
833 
838 std::string
840 
841  return name;
842 
843 }
844 
845 
846 
847 
std::vector< Gate * > parents
list of parent gates to be applied in the circuit prior to this current gate
Definition: Gate.h:95
Matrix dot(Matrix &A, Matrix &B)
Call to calculate the product of two complex matrices by calling method zgemm3m from the CBLAS librar...
Definition: dot.cpp:38
Gate()
Default constructor of the class.
Definition: Gate.cpp:39
virtual void parameters_for_calc_one_qubit(double &ThetaOver2, double &Phi, double &Lambda)
Calculate the matrix of a U3 gate gate corresponding to the given parameters acting on a single qbit ...
Definition: Gate.cpp:767
void print(const std::stringstream &sstream, int verbose_level=1) const
Call to print output messages in the function of the verbosity level.
Definition: logging.cpp:55
Header file for a class for the representation of general gate operations.
void add_child(Gate *child)
Call to add a child gate to the current gate.
Definition: Gate.cpp:426
virtual Matrix get_matrix()
Call to retrieve the operation matrix.
Definition: Gate.cpp:129
void clear_children()
Call to erase data on children.
Definition: Gate.cpp:441
int control_qbit
The index of the qubit which acts as a control qubit (control_qbit >= 0) in controlled operations...
Definition: Gate.h:87
virtual std::vector< Matrix > apply_derivate_to(Matrix_real &parameters_mtx_in, Matrix &input, int parallel)
Call to evaluate the derivate of the circuit on an inout with respect to all of the free parameters...
Definition: Gate.cpp:271
virtual Gate * clone()
Call to create a clone of the present class.
Definition: Gate.cpp:513
int stride
The column stride of the array. (The array elements in one row are a_0, a_1, ... a_{cols-1}, 0, 0, 0, 0. The number of zeros is stride-cols)
Definition: matrix_base.hpp:46
void set_children(std::vector< Gate *> &children_)
Call to set the children of the current gate.
Definition: Gate.cpp:802
virtual void set_qbit_num(int qbit_num_in)
Set the number of qubits spanning the matrix of the operation.
Definition: Gate.cpp:102
void apply_kernel_to_input_AVX_small(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
AVX kernel to apply single qubit gate kernel on an input matrix (efficient for small inputs) ...
int target_qbit
The index of the qubit on which the operation acts (target_qbit >= 0)
Definition: Gate.h:85
void apply_kernel_to_input_AVX(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
AVX kernel to apply single qubit gate kernel on an input matrix (single threaded) ...
void apply_kernel_to_state_vector_input_parallel_OpenMP_AVX(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
Parallel AVX kernel on a state vector (parallelized with OpenMP)
void apply_kernel_to_state_vector_input_AVX(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
AVX kernel on a state vector.
void apply_kernel_from_right(Matrix &u3_1qbit, Matrix &input)
Call to apply the gate kernel on the input state or unitary from right (no AVX support) ...
Definition: Gate.cpp:613
virtual ~Gate()
Destructor of the class.
Definition: Gate.cpp:95
QGD_Complex16 mult(QGD_Complex16 &a, QGD_Complex16 &b)
Call to calculate the product of two complex scalars.
Definition: common.cpp:259
void apply_kernel_to_state_vector_input_parallel_AVX(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
Parallel AVX kernel on a state vector (parallelized with Intel TBB)
int matrix_size
The size N of the NxN matrix associated with the operations.
Definition: Gate.h:89
scalar * get_data() const
Call to get the pointer to the stored data.
int parameter_start_idx
the index in the parameter array (corrensponding to the encapsulated circuit) where the gate paramete...
Definition: Gate.h:93
void apply_kernel_to_input_AVX_parallel(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
Parallel AVX kernel to apply single qubit gate kernel on an input matrix.
gate_type type
The type of the operation (see enumeration gate_type)
Definition: Gate.h:83
virtual void apply_from_right(Matrix &input)
Call to apply the gate on the input array/matrix by input*Gate.
Definition: Gate.cpp:285
int rows
The number of rows.
Definition: matrix_base.hpp:42
int cols
The number of columns.
Definition: matrix_base.hpp:44
void apply_kernel_to_input(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
Call to apply kernel to apply single qubit gate kernel on an input matrix.
gate_type get_type()
Call to get the type of the operation.
Definition: Gate.cpp:495
void set_control_qbit(int control_qbit_in)
Call to set the control qubit for the gate operation.
Definition: Gate.cpp:308
void set_parameter_start_idx(int start_idx)
Call to set the starting index of the parameters in the parameter array corresponding to the circuit ...
Definition: Gate.cpp:779
std::vector< Gate * > get_parents()
Call to get the parents of the current gate.
Definition: Gate.cpp:463
int get_parameter_start_idx()
Call to get the starting index of the parameters in the parameter array corresponding to the circuit ...
Definition: Gate.cpp:814
void add_parent(Gate *parent)
Call to add a parent gate to the current gate.
Definition: Gate.cpp:409
void set_target_qbit(int target_qbit_in)
Call to set the target qubit for the gate operation.
Definition: Gate.cpp:324
Structure type representing complex numbers in the SQUANDER package.
Definition: QGDTypes.h:38
virtual Matrix_real extract_parameters(Matrix_real &parameters)
Call to extract parameters from the parameter array corresponding to the circuit, in which the gate i...
Definition: Gate.cpp:827
void clear_parents()
Call to erase data on parents.
Definition: Gate.cpp:451
virtual Matrix calc_one_qubit_u3()
Calculate the matrix of the constans gates.
Definition: Gate.cpp:750
void set_matrix(Matrix input)
Call to set the stored matrix in the operation.
Definition: Gate.cpp:299
std::vector< Gate * > get_children()
Call to get the children of the current gate.
Definition: Gate.cpp:474
int Power_of_2(int n)
Calculates the n-th power of 2.
Definition: common.cpp:117
Class to store data of complex arrays and its properties.
Definition: matrix.h:38
int size() const
Call to get the number of the allocated elements.
int get_parameter_num()
Call to get the number of free parameters.
Definition: Gate.cpp:486
std::string name
A string labeling the gate operation.
Definition: Gate.h:79
std::vector< Gate * > children
list of child gates to be applied after this current gate
Definition: Gate.h:97
int get_target_qbit()
Call to get the index of the target qubit.
Definition: Gate.cpp:370
void apply_kernel_to(Matrix &u3_1qbit, Matrix &input, bool deriv=false, int parallel=0)
Call to apply the gate kernel on the input state or unitary with optional AVX support.
Definition: Gate.cpp:537
Base class for the representation of general gate operations.
Definition: Gate.h:73
std::string get_name()
Call to get the name label of the gate.
Definition: Gate.cpp:839
virtual std::vector< int > get_involved_qubits()
Call to get the qubits involved in the gate operation.
Definition: Gate.cpp:386
int parameter_num
the number of free parameters of the operation
Definition: Gate.h:91
Header file for commonly used functions and wrappers to CBLAS functions.
virtual void apply_to_list(std::vector< Matrix > &inputs, int parallel)
Call to apply the gate on a list of inputs.
Definition: Gate.cpp:187
void set_parents(std::vector< Gate *> &parents_)
Call to set the parents of the current gate.
Definition: Gate.cpp:790
double real
the real part of a complex number
Definition: QGDTypes.h:40
void apply_kernel_to_state_vector_input(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
Call to apply a gate kernel on a state vector.
int qbit_num
number of qubits spanning the matrix of the operation
Definition: Gate.h:81
int get_qbit_num()
Call to get the number of qubits composing the unitary.
Definition: Gate.cpp:504
gate_type
Type definition of operation types (also generalized for decomposition classes derived from the class...
Definition: Gate.h:34
virtual void apply_to(Matrix &input, int parallel)
Call to apply the gate on the input array/matrix.
Definition: Gate.cpp:237
void apply_kernel_to_state_vector_input_parallel(Matrix &u3_1qbit, Matrix &input, const bool &deriv, const int &target_qbit, const int &control_qbit, const int &matrix_size)
Call to apply a gate kernel on a state vector.
virtual void reorder_qubits(std::vector< int > qbit_list)
Call to reorder the qubits in the matrix of the operation.
Definition: Gate.cpp:339
Matrix matrix_alloc
Matrix of the operation.
Definition: Gate.h:102
int get_control_qbit()
Call to get the index of the control qubit.
Definition: Gate.cpp:378
Class to store data of complex arrays and its properties.
Definition: matrix_real.h:39
double imag
the imaginary part of a complex number
Definition: QGDTypes.h:42