# C++ Tutorial ## Quantum states ### Generate quantum states Generate $n$ qubit quantum states using `QuantumState` class and initialize it as $|0\rangle^{\otimes n}$. ``` cpp #include int main(){ unsigned int n = 5; QuantumState state(n); state.set_zero_state(); return 0; } ``` You can not generate the quantum states if the memory is not sufficient. With the memory of a typical laptop or desktop computer, the limit of qubits you can create is about 26, 27 qubits. ### Obtain data of quantum states The quantum state is expressed as an array of length $2^n$. Note that if the quantum state is formed by GPU, and if the $n$ is large, it can become an extremely heavy operation. ``` cpp #include int main(){ unsigned int n = 5; QuantumState state(n); state.set_zero_state(); // If use GNU C++, obtain double_Complex array // If use MSVC, obtain std::complex array const CTYPE* raw_data_c = state.data_c(); // obtain std::complex array const CPPCTYPE* raw_data_cpp = state.data_cpp(); } ``` If you want to set the quantum state directly to the specified array, it is recommended to create the corresponding quantum gate and perform it as an operation of the quantum gate. ### Initialize quantum states The generated quantum state can be initialized to a computational basis using the `set_computational_basis` or to a random state using the `set_Harr_random_state`. Note that in Qulacs, the subscripts of the qubits start from 0, and the rightmost bit is the 0-th qubit when written as $\ket{0000}$ (In other libraries and textbooks, the leftmost bit may be the 0-th qubit). ``` cpp #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); // Initialize as |00101> state.set_computational_basis(0b00101); // Generate random initial states state.set_Haar_random_state(); // Generate random initial state with specifying seed state.set_Haar_random_state(0); return 0; } ``` ### Copy and load quantum state data The quantum state can be copied and loaded from other quantum state data. ``` cpp #include int main(){ unsigned int n = 5; QuantumState state(n); state.set_computational_basis(0b00101); // Create new quantum state by copying auto second_state = state.copy(); // Create a new quantum state and copy the existing state vector QuantumState third_state(n); third_state.load(&state); return 0; } ``` ### Operate classic registers The quantum state can be read and written as a classic register. ``` cpp #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); // Set the 3rd classical register as 1 int register_position = 3; int register_value = 1; state.set_classical_bit(register_position, register_value); // Obtain the value of the 3rd classical register int obtained_value; obtained_value = state.get_classical_bit(register_position); return 0; } ``` ### Calculate quantum states The following operations do not change quantum states. Calculations that changes quantum states are performed by quantum gates or quantum circuits. ``` cpp #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); // Calculate norm double norm = state.get_squared_norm(); // Calculate the entropy when measured in z-basis double entropy = state.get_entropy(); // Calculate the probability that the index-th qubit is 0 when measured in z-basis unsigned int index = 3; double zero_prob = state.get_zero_probability(index); // Calculate marginal probabilities // (Here is an example of the probability that 0,3-th qubit is measured as 0 and 1,2-th qubit is measured as 1) std::vector value_list = { 0,1,1,0,2 }; double marginal_prob = state.get_marginal_probability(value_list); return 0; } ``` ### Calculate the inner product of quantum states The inner product of quantum states can be calculated by the `inner_product` function. ``` cpp #include int main(){ unsigned int n = 5; QuantumState state_ket(n); state_ket.set_zero_state(); QuantumState state_bra(n); state_bra.set_Haar_random_state(); // Calculate the inner product std::complex value = state::inner_product(&state_ket, &state_bra); return 0; } ``` ## Quantum Gate ### Generate and operate quantum gates The quantum gate implemented by default is generated through the function of `gate_factory`. This quantum gate can operate on the quantum state specified by the argument of `update_quantum_state`. The quantum gate generated by `gate_factory` is not released automatically, so you must release it. ``` cpp #define _USE_MATH_DEFINES #include #include #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); // operation of X gata unsigned int index = 3; auto x_gate = gate::X(index); x_gate->update_quantum_state(&state); // Rotation by PI/2 using Y gate double angle = M_PI / 2.0; auto ry_gate = gate::RY(index, angle); ry_gate->update_quantum_state(&state); delete x_gate; delete ry_gate; return 0; } ``` Quantum gates defined in the `gate` namespace are as following: - Single-qubit Pauli operation: `Identity`, `X`, `Y`, `Z` - Single-qubit Clifford operation: `H`, `S`, `Sdag`, `T`, `Tdag`, `sqrtX`, `sqrtXdag`, `sqrtY`, `sqrtYdag` - Two-qubit Clifford operation: `CNOT`, `CZ`, `SWAP` - Single-qubit Pauli rotation: `RX`, `RY`, `RZ`, `RotX`, `RotY`, `RotZ`, `RotInvX`, `RotInvY`, `RotInvZ` - General Pauli operation: `Pauli`, `PauliRotation` - IBMQ basis-gate: `U1`, `U2`, `U3` - General gate: `DenseMatrix` - Measurement : `Measurement` - Noise : `BitFlipNoise`, `DephasingNoise`, `IndependentXZNoise`, `DepolarizingNoise` Rotation gates `RX`, `RY`, and `RZ` operate as Pauli rotation $\exp(i\frac{\theta}{2}P)$ based on corresponding Pauli operator $P$ and argument $\theta$. Please check the API documents for details of each gate. Definitions are as follows: `RX`:$R_X(\theta) = \exp(i\frac{\theta}{2} X) = \begin{pmatrix} \cos(\frac{\theta}{2}) & i\sin(\frac{\theta}{2}) \\ i\sin(\frac{\theta}{2}) & \cos(\frac{\theta}{2}) \end{pmatrix} $ `RY`:$R_Y(\theta) = \exp(i\frac{\theta}{2} Y) = \begin{pmatrix} \cos(\frac{\theta}{2}) & \sin(\frac{\theta}{2}) \\ -\sin(\frac{\theta}{2}) & \cos(\frac{\theta}{2}) \end{pmatrix} $ `RZ`:$R_Z(\theta) = \exp(i\frac{\theta}{2} Z) = \begin{pmatrix} e^{i\frac{\theta}{2}} & 0 \\ 0 & e^{-i\frac{\theta}{2}} \end{pmatrix} $ The rotation gates used in Qulacs are basically `RX`, `RY`, and `RZ`, but there are also `RotX`, `RotY`, and `RotZ` with opposite rotation directions. The definitions of these operations are as follows: `RotX`:$RotX(\theta) = \exp(-i\frac{\theta}{2}X)$ `RotY`:$RotX(\theta) = \exp(-i\frac{\theta}{2}Y)$ `RotZ`:$RotX(\theta) = \exp(-i\frac{\theta}{2}Z)$ Gates with opposite rotation direction to `RotX`, `RotY` and `RotZ` are `RotInvX`, `RotInvY` and `RotInvZ`. ### Merge quantum gates You can generate a new single quantum gate by merging the quantum gates that continue to operate in sequence. You have to release the synthesized gate by yourself. ``` cpp #define _USE_MATH_DEFINES #include #include #include #include #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); unsigned int index = 3; auto x_gate = gate::X(index); double angle = M_PI / 2.0; auto ry_gate = gate::RY(index, angle); // combine X, RY in the successive operation order auto x_and_ry_gate = gate::merge(x_gate, ry_gate); x_and_ry_gate->update_quantum_state(&state); delete x_gate; delete ry_gate; delete x_and_ry_gate; return 0; } ``` ### Sum of quantum gate matrices A new gate can be generated by summing gates. (Not available for gates with control-qubit, because that operation is undefined yet.) ``` cpp #define _USE_MATH_DEFINES #include #include #include #include #include int main() { auto gate00 = gate::merge(gate::P0(0), gate::P0(1)); auto gate11 = gate::merge(gate::P1(0), gate::P1(1)); // |00><00| + |11><11| auto proj_00_or_11 = gate::add(gate00, gate11); std::cout << proj_00_or_11 << std::endl; auto gate_ii_zz = gate::add(gate::Identity(0), gate::merge(gate::Z(0), gate::Z(1))); auto gate_ii_xx = gate::add(gate::Identity(0), gate::merge(gate::X(0), gate::X(1))); auto proj_00_plus_11 = gate::merge(gate_ii_zz, gate_ii_xx); // ((|00>+|11>)(<00|+<11|))/2 = (II + ZZ)(II + XX)/4 proj_00_plus_11->multiply_scalar(0.25); std::cout << proj_00_plus_11 << std::endl; return 0; } ``` ### Special quantum gate and common quantum gate In cppsim, the basic quantum gates are divided into the following two ways: - Special gate: There are dedicated speed-up functions for utilizing the special gate. - Common gate: The gate holds the gate matrix and operates by multiplying the matrix. The special gate is faster than the common gate because of the dedicated functions. But in a special gate, operations that change the function of a quantum gate, such as increasing the number of control qubits, cannot be performed later. This kind of change can be made only when the special gate is transformed into a common gate, which can be realized by `gate::to_matrix_gate`. Here's an example: ``` cpp #include #include #include #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); unsigned int index = 3; auto x_gate = gate::X(index); // Add control qubit so that only operates when 1st-qubit is 0 auto x_mat_gate = gate::to_matrix_gate(x_gate); unsigned int control_index = 1; unsigned int control_with_value = 0; x_mat_gate->add_control_qubit(control_index, control_with_value); x_mat_gate->update_quantum_state(&state); delete x_gate; delete x_mat_gate; return 0; } ``` Please check the API documents for details of the special quantum gate. ### Obtain the gate matrix of the quantum gate The gate matrix of the generated quantum gate can be obtained, but gate matrices do not include the control qubit. Especially, note that gates without gate matrix (ex. $n$-qubit Pauli rotation gate) require a very long time and memory to obtain the matrix. ``` cpp #include #include #include #include int main(){ unsigned int index = 3; auto x_gate = gate::X(index); // obtain the matrix element // ComplexMatrix is a complex matrix type with RowMajor by Eigen::MatrixXcd ComplexMatrix matrix; x_gate->set_matrix(matrix); std::cout << matrix << std::endl; return 0; } ``` ### Obtain information about quantum gate Debug information of quantum gate can be shown by using ostream. Note that if the gate matrix of the quantum gate is very large, it takes a long time. Quantum gates with dedicated functions do not display their gate matrix. ``` cpp #include #include #include #include int main(){ unsigned int index = 3; auto x_gate = gate::X(index); std::cout << x_gate << std::endl; delete x_gate; return 0; } ``` ### Implement of common quantum gate Cppsim implements various maps of quantum information in the following forms. #### Unitary operation Implemented as the quantum gate. #### Projection operator and Kraus operator, etc. Implemented as the quantum gate. In general, the norm of the quantum state is not preserved after the operation. The gate can be generated by `DenseMatrix`. ``` cpp #define _USE_MATH_DEFINES #include #include #include #include #include #include int main() { ComplexMatrix one_qubit_matrix(2, 2); one_qubit_matrix << 0, 1, 1, 0; auto one_qubit_gate = gate::DenseMatrix(0, one_qubit_matrix); std::cout << one_qubit_gate << std::endl; ComplexMatrix two_qubit_matrix(4,4); two_qubit_matrix << 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0; auto two_qubit_gate = gate::DenseMatrix({0,1}, two_qubit_matrix); std::cout << two_qubit_gate << std::endl; return 0; } ``` #### Probabilistic unitary operations With given multiple unitary operations and probability distributions, stochastic unitary operations can be created by `Probabilistic` function. ``` cpp #define _USE_MATH_DEFINES #include #include #include #include #include #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); unsigned int index = 3; auto x_gate = gate::X(index); auto z_gate = gate::Z(index); auto probabilistic_xz_gate = gate::Probabilistic({ 0.1,0.2 } , { x_gate,z_gate }); auto depolarizing_gate = gate::DepolarizingNoise(index, 0.3); depolarizing_gate->update_quantum_state(&state); probabilistic_xz_gate->update_quantum_state(&state); return 0; } ``` #### CPTP-map CPTP-map can be created by giving the `CPTP` function a list of Kraus operators satisfying completeness. ``` cpp #define _USE_MATH_DEFINES #include #include #include #include #include #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); unsigned int index = 3; auto p0 = gate::P0(index); auto p1_fix = gate::merge(gate::P1(index), gate::X(index)); auto correction = gate::CPTP({p0,p1_fix}); auto noise = gate::BitFlipNoise(index,0.1); noise->update_quantum_state(&state); correction->update_quantum_state(&state); return 0; } ``` #### POVM Since it is the same as Instrument in numerical calculation, it is realized as Instrument. #### Instrument In addition to the general CPTP-map operation, `Instrument` is an operation that obtains the array subscripts of the randomly acting Kraus operator. For example, a measurement on the Z basis is to operate on the CPTP-map consisting of `P0` and `P1` and know which one was operated. In cppsim, this is achieved by specifying the information of the CPTP-map and the address of the classic register in which the subscripts of the operated Kraus operator are written in the `Instrument` function. ``` cpp #define _USE_MATH_DEFINES #include #include #include #include #include #include int main() { auto gate00 = gate::merge(gate::P0(0), gate::P0(1)); auto gate01 = gate::merge(gate::P0(0), gate::P1(1)); auto gate10 = gate::merge(gate::P1(0), gate::P0(1)); auto gate11 = gate::merge(gate::P1(0), gate::P1(1)); std::vector gate_list = { gate00, gate01, gate10, gate11 }; unsigned int classical_pos = 0; auto gate = gate::Instrument(gate_list, classical_pos); QuantumState state(2); state.set_Haar_random_state(); std::cout << state << std::endl; gate->update_quantum_state(&state); unsigned int result = state.get_classical_value(classical_pos); std::cout << state << std::endl; std::cout << result << std::endl; return 0; } ``` #### Adaptive The operation is performed or not performed depending on the value written to the classical register. In cppsim, this is achieved by specifying a function that takes a register of type `std::vector` as an argument and returns a `bool` type in the `Adaptive` function. ``` cpp #define _USE_MATH_DEFINES #include #include #include #include #include #include int main() { unsigned int n = 5; QuantumState state(n); state.set_zero_state(); unsigned int index = 3; auto h = gate::H(index); h->update_quantum_state(&state); auto meas = gate::Measurement(index,0); meas->update_quantum_state(&state); auto condition = [](const std::vector reg){ return reg[0]==1; }; auto correction = gate::Adaptive(gate::X(index), condition); correction->update_quantum_state(&state); return 0; } ``` #### CP-map If Kraus rank is 1, please treat it as a single Kraus operator as described above. In other cases, please adjust the Kraus operator so that it becomes TP, and then adjust it by applying the `Identity` operator multiplied by a constant with the `multiply_scalar()` function. ## Quantum Circuits ### Construct the quantum circuit A quantum circuit is represented as a set of quantum gates. For example, you can construct a quantum circuit as follows: ``` cpp #include #include #include int main(){ unsigned int n = 5; QuantumState state(n); state.set_zero_state(); // Define quantum circuit QuantumCircuit circuit(n); // Add gate to quantum circuit for(int i=0;i #include #include #include int main() { unsigned int n = 5; unsigned int depth = 10; QuantumCircuit circuit(n); for (int d = 0; d < depth; ++d) { for (int i = 0; i < n; ++i) { circuit.add_gate(gate::H(i)); } } // 量子回路の最適化 QuantumCircuitOptimizer opt; unsigned int max_block_size = 3; opt.optimize(&circuit, max_block_size); return 0; } ``` ### Obtain debug information of quantum circuits The same as the quantum gate, debug information of the quantum circuit can be shown by using `ostream`. ``` cpp #include #include #include int main() { unsigned int n = 5; unsigned int depth = 10; QuantumCircuit circuit(n); for (int d = 0; d < depth; ++d) { for (int i = 0; i < n; ++i) { circuit.add_gate(gate::H(i)); } } // Output information of quantum circuit std::cout << circuit << std::endl; return 0; } ``` ## Observables ### Generate observables Observables are represented as a set of Pauli operators. The Pauli operator can be defined as follows: ``` cpp #include #include int main() { unsigned int n = 5; double coef = 2.0; std::string Pauli_string = "X 0 X 1 Y 2 Z 4"; Observable observable(n); observable.add_operator(coef,Pauli_string.c_str()); return 0; } ``` ### Cooperate with OpenFermion Observables can be generated from files in the following formats generated using OpenFermion. At this time, the observable has the minimum size necessary to generate it. For example, it is possible to generate an observable by reading an observable obtained using OpenFermion such as the following. ``` python from openfermion.ops import FermionOperator from openfermion.transforms import bravyi_kitaev h_00 = h_11 = -1.252477 h_22 = h_33 = -0.475934 h_0110 = h_1001 = 0.674493 h_2332 = h_3323 = 0.697397 h_0220 = h_0330 = h_1221 = h_1331 = h_2002 = h_3003 = h_2112 = h_3113 = 0.663472 h_0202 = h_1313 = h_2130 = h_2310 = h_0312 = h_0132 = 0.181287 fermion_operator = FermionOperator('0^ 0', h_00) fermion_operator += FermionOperator('1^ 1', h_11) fermion_operator += FermionOperator('2^ 2', h_22) fermion_operator += FermionOperator('3^ 3', h_33) fermion_operator += FermionOperator('0^ 1^ 1 0', h_0110) fermion_operator += FermionOperator('2^ 3^ 3 2', h_2332) fermion_operator += FermionOperator('0^ 3^ 3 0', h_0330) fermion_operator += FermionOperator('1^ 2^ 2 1', h_1221) fermion_operator += FermionOperator('0^ 2^ 2 0', h_0220-h_0202) fermion_operator += FermionOperator('1^ 3^ 3 1', h_1331-h_1313) fermion_operator += FermionOperator('0^ 1^ 3 2', h_0132) fermion_operator += FermionOperator('2^ 3^ 1 0', h_0132) fermion_operator += FermionOperator('0^ 3^ 1 2', h_0312) fermion_operator += FermionOperator('2^ 1^ 3 0', h_0312) ## Bravyi-Kitaev transformation bk_operator = bravyi_kitaev(fermion_operator) ## output fp = open("H2.txt", 'w') fp.write(str(bk_operator)) fp.close() ``` The `H2.txt` file generated by the above Python code has the following format. ``` txt (-0.8126100000000005+0j) [] + (0.04532175+0j) [X0 Z1 X2] + (0.04532175+0j) [X0 Z1 X2 Z3] + (0.04532175+0j) [Y0 Z1 Y2] + (0.04532175+0j) [Y0 Z1 Y2 Z3] + (0.17120100000000002+0j) [Z0] + (0.17120100000000002+0j) [Z0 Z1] + (0.165868+0j) [Z0 Z1 Z2] + (0.165868+0j) [Z0 Z1 Z2 Z3] + (0.12054625+0j) [Z0 Z2] + (0.12054625+0j) [Z0 Z2 Z3] + (0.16862325+0j) [Z1] + (-0.22279649999999998+0j) [Z1 Z2 Z3] + (0.17434925+0j) [Z1 Z3] + (-0.22279649999999998+0j) [Z2] ``` You can create an observable from such a file through a function as follows: ``` cpp #include #include int main() { unsigned int n = 5; std::string filename = "H2.txt"; Observable* observable = observable::create_observable_from_openfermion_file(filename); delete observable; return 0; } ``` ### Evaluate observable An evaluation of the expected value of the observable for the state can be obtained. ``` cpp #include #include #include int main() { unsigned int n = 5; double coef = 2.0; std::string Pauli_string = "X 0 X 1 Y 2 Z 4"; Observable observable(n); observable.add_operator(coef, Pauli_string.c_str()); QuantumState state(n); observable.get_expectation_value(&state); return 0; } ``` ### Rotation of Observable The rotation of Observable $H$, $e^{i\theta H}$, is performed by Trotter expansion. `num_repeats` defaults as the following code, but can be optionally specified by the user. ``` cpp #include #include #include int main() { UINT n; UINT num_repeats; double theta = 0.1; Observable* observable = observable::create_observable_from_openfermion_file("../test/cppsim/H2.txt"); n = observable->get_qubit_count(); QuantumState state(n); state.set_computational_basis(0); QuantumCircuit circuit(n); num_repeats = (UINT)std::ceil(theta * (double)n* 100.); circuit.add_observable_rotation_gate(*observable, theta, num_repeats); circuit.update_quantum_state(&state); auto result = observable->get_expectation_value(&state); std::cout << result << std::endl; delete observable; return 0; } ``` ## Parametric Quantum Circuits Defining a quantum circuit as the `ParametricQuantumCircuit` class allows you to use some functions that are useful for optimizing quantum circuits using variational methods, in addition to the usual functions of the `QuantumCircuit` class. ### Examples of parametric quantum circuits Quantum gates with one rotation angle (X-rot, Y-rot, Z-rot, multi_qubit_pauli_rotation) can be added to quantum circuits as parametric quantum gates. For quantum gates added as parametric gates, the number of parametric gates can be extracted after the quantum circuit is constructed, and the rotation angle can be changed later. ``` cpp #include #include #include int main(){ const UINT n = 3; const UINT depth = 10; // create n-qubit parametric circuit ParametricQuantumCircuit* circuit = new ParametricQuantumCircuit(n); Random random; for (UINT d = 0; d < depth; ++d) { // add parametric X,Y,Z gate with random initial rotation angle for (UINT i = 0; i < n; ++i) { circuit->add_parametric_RX_gate(i, random.uniform()); circuit->add_parametric_RY_gate(i, random.uniform()); circuit->add_parametric_RZ_gate(i, random.uniform()); } // add neighboring two-qubit ZZ rotation for (UINT i = d % 2; i + 1 < n; i+=2) { circuit->add_parametric_multi_Pauli_rotation_gate({ i,i + 1 }, { 3,3 }, random.uniform()); } } // get parameter count UINT param_count = circuit->get_parameter_count(); // get current parameter, and set shifted parameter for (UINT p = 0; p < param_count; ++p) { double current_angle = circuit->get_parameter(p); circuit->set_parameter(p, current_angle + random.uniform()); } // create quantum state and update QuantumState state(n); circuit->update_quantum_state(&state); // output state and circuit info std::cout << state << std::endl; std::cout << circuit << std::endl; // release quantum circuit delete circuit; } ```