Qulacs Python Advanced Guide

This chapter is for those who are familiar with the terminology of quantum information and want fine tuning for numerical calculations. For more information on terminology, please see the textbook Quantum Computation and Quantum Information by M.A. Nielsen et al..

Quantum States

This class allocates and manages \(2^n\) complex arrays on the CPU/GPU with the precision of complex128.

Create and Destroy

Necessary memory is secured at the time of instance creation and released when the instance is destroyed, but you can explicitly destroy it to release the memory with del. A formatted display of the state vector is provided through __repr__ function overrides.

[1]:
from qulacs import QuantumState
n = 2
state = QuantumState(n)
print(state)
del state
 *** Quantum State ***
 * Qubit Count : 2
 * Dimension   : 4
 * State vector :
(1,0)
(0,0)
(0,0)
(0,0)

Transform between quantum state and numpy array

Transformation between quantum state and numpy array can be realized by get_vector or load function. In principle, it is not checked whether the norm is saved.

[2]:
from qulacs import QuantumState

state = QuantumState(2)
vec = state.get_vector()
print(vec)
state.load([0,1,2,3])
print(state.get_vector())
[1.+0.j 0.+0.j 0.+0.j 0.+0.j]
[0.+0.j 1.+0.j 2.+0.j 3.+0.j]

Copy of quantum state

Quantum states can generate new instances of the same state by copy. By giving a quantum state to the load function, it is possible to copy a quantum vector of another quantum state without securing a new area in the existing quantum state. This allows you to reuse the already allocated space. You can use the allocate_buffer function if you want to secure a state vector of the same size as the quantum state you already have and without copying the state.

[3]:
from qulacs import QuantumState

initial_state = QuantumState(3)
buffer = initial_state.allocate_buffer()
for ind in range(10):
    buffer.load(initial_state)
    # some computation and get results

Initialization of quantum state

The following is an function initializing a quantum state to a specific state.

[4]:
from qulacs import QuantumState

n = 3
state = QuantumState(n)
# Initialize as |0> state
state.set_zero_state()
print(state.get_vector())

# Initialize the specified value to the calculation base in binary notation
state.set_computational_basis(0b101)
print(state.get_vector())

# Initialize to random pure state with Haar measure using argument value as seed
# If no value is specified, the time function is used as a seed. Pseudo random number uses xorshift.
state.set_Haar_random_state(0)
print(state.get_vector())
[1.+0.j 0.+0.j 0.+0.j 0.+0.j 0.+0.j 0.+0.j 0.+0.j 0.+0.j]
[0.+0.j 0.+0.j 0.+0.j 0.+0.j 0.+0.j 1.+0.j 0.+0.j 0.+0.j]
[ 0.17875932+0.12507972j  0.28161322-0.20162504j  0.21880416-0.04754893j
  0.01570731+0.46979145j -0.0216192 +0.46023702j -0.20700612+0.09443371j
 -0.03203312+0.3305667j   0.43080876+0.03745933j]

Check quantum state

The following is a list of functions for checking quantum state information without changing the quantum state.

[5]:
from qulacs import QuantumState

n = 5
state = QuantumState(n)
state.set_Haar_random_state(0)

# Get quantum bit numbers
qubit_count = state.get_qubit_count()
print("qubit_count", qubit_count)

# Get the probability that the specified qubit will be measured as 0
prob = state.get_zero_probability(1)
print("zero_prob_1", prob)

# Get arbitrary marginal probabilities
# Argument is an array of the same length as the number of qubits
# Specify 0,1,2. 0,1 is the probability of the subscript measured at that value
# 2 means that bit is peripheralized.
# For example, calculation of the probability that the third is measured as 0 and the 0th is measured as 1:
prob = state.get_marginal_probability([1,2,2,0,2])
print("marginal_prob", prob)

# Get the entropy of the probability distribution when measured on the Z basis
ent = state.get_entropy()
print("entropy", ent)

# Get squared norm (<a|a>)
# Because the operation maybe not Trace preserving, the norm of state does not necessarily to be 1.
sq_norm = state.get_squared_norm()
print("sqaured_norm", sq_norm)

# The number of measurements and sampling of all qubits Z-basis is given by the argument.
# Get a list of integers converted from the resulting binaries.
samples = state.sampling(10)
print("sampling", samples)

# Get a character string indicating whether the state vector is on CPU or GPU
dev_type = state.get_device_name()
print("device", dev_type)
qubit_count 5
zero_prob_1 0.6177236324405834
marginal_prob 0.3266606041128599
entropy 3.1462789915919283
sqaured_norm 0.9999999999999998
sampling [10, 4, 4, 28, 21, 9, 22, 1, 21, 21]
device cpu

Deformation of quantum state

The following function rewrites the quantum state.

[6]:
from qulacs import QuantumState
state = QuantumState(2)
state.set_computational_basis(0)
buffer = QuantumState(2)
buffer.set_computational_basis(2)
print("state" , state.get_vector())
print("buffer", buffer.get_vector())

# Sum of quantum state (state <- state+buffer)
# Add the buffer state to the state to create a superimposed state.
# The norm after operation generally is not 1.
state.add_state(buffer)
print("added", state.get_vector())

# Product of quantum state and complex number
# Multiplies all elements by the complex number of the argument.
# The norm after operation generally is not 1.
coef = 0.5 + 0.1j
state.multiply_coef(coef)
print("mul_coef", state.get_vector())

# Normalization of quantum states
# The current squared norm needs to be provided as an argument.
squared_norm = state.get_squared_norm()
print("sq_norm", squared_norm)
state.normalize(squared_norm)
print("normalized", state.get_vector())
print("sq_norm", state.get_squared_norm())
state [1.+0.j 0.+0.j 0.+0.j 0.+0.j]
buffer [0.+0.j 0.+0.j 1.+0.j 0.+0.j]
added [1.+0.j 0.+0.j 1.+0.j 0.+0.j]
mul_coef [0.5+0.1j 0. +0.j  0.5+0.1j 0. +0.j ]
sq_norm 0.52
normalized [0.69337525+0.13867505j 0.        +0.j         0.69337525+0.13867505j
 0.        +0.j        ]
sq_norm 0.9999999999999998

Opetation on classic register

Quantum states have classical registers as arrays of integers with variable length. The classical register is used to write the result of the Instrument operation or to describe a gate that executes conditions as the result of the classical register. The value of a classic register that has not yet been written is 0. The classical register is copied at the same time when the quantum state is copied by the copy and load functions.

[7]:
from qulacs import QuantumState
state = QuantumState(3)
position = 0
# Write the value to position number
value = 20
state.set_classical_value(position, value)
# Get the value of the register at position position number
obtained = state.get_classical_value(position)
print(obtained)
20

Calculation between quantum states

The inner product between quantum states can be obtained by inner_product.

[8]:
from qulacs import QuantumState
from qulacs.state import inner_product

n = 5
state_bra = QuantumState(n)
state_ket = QuantumState(n)
state_bra.set_Haar_random_state()
state_ket.set_computational_basis(0)

# Calculation of inner product
value = inner_product(state_bra, state_ket)
print(value)
(-0.06763387457138917-0.02885458593002302j)

Though the usage is the same as QuantumState, there are two aspects to keep in mind: * The get_vector function takes a long time because it requires copying between the GPU and CPU. This function should be avoided whenever possible. * inner_product between CPU / GPU states cannot be calculated. It is possible to load a state vector between the GPU and CPU state vectors, but it is time consuming and should be avoided.

Calculation using GPU

When Qulacs is installed from qulacs-gpu package, QuantumStateGpu is available. Except the different class name, the usage is the same as QuantumState.

from qulacs import QuantumStateGpu
state = QuantumStateGpu(2)
print(state)
# print(state.get_device_name())
# gpu

Quantum gates

Types of quantum gate

Quantum gates are divided into two types: special gates and general gates. In Qulacs, quantum gates are not limited to unitary operators, and the operation of updating an arbitrary quantum state, such as Instrument and CPTP-map, is also called a gate.

Special gates are those that have a pre-specified gate matrix and can only perform limited deformations on quantum gates. For example, Pauli gate, rotation Pauli gate, projection measurement, etc. are supported. The advantage of a special gate is that the update function of the quantum state is more efficient than a general gate with limited properties. Also, at the time of definition, it holds information on whether or not it is diagonalized by each qubit at the basis of some Pauli, and this information is used for circuit optimization. The disadvantage of special gates is that the possible operations on the gates are limited for the reasons mentioned above.

A gate that has an explicit operation gate matrix is called a general gate. The advantage of general gates is that you can specify the gate matrix as you like, but the disadvantage is that the updates are slower than special gates.

Common operations for quantum gates

The gate matrix of the generated quantum gate can be obtained with the get_matrix function. Control qubits are not included in the gate matrix. In particular, be careful when obtain gates that do not have a gate matrix (for example, n-qubit Pauli rotating gates), which require a very large amount of memory and time. print function displays the gate information.

[9]:
import numpy as np
from qulacs.gate import X
gate = X(2)
mat = gate.get_matrix()
print(mat)
print(gate)
[[0.+0.j 1.+0.j]
 [1.+0.j 0.+0.j]]
 *** gate info ***
 * gate name : X
 * target    :
 2 : commute X
 * control   :
 * Pauli     : yes
 * Clifford  : yes
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no

Special gate

The special gates are listed below. #### 1 qubit gate Takes the index of the target bit as the first argument.

[10]:
from qulacs.gate import Identity # Identity matrix
from qulacs.gate import X, Y, Z     # Pauli
from qulacs.gate import H, S, Sdag, sqrtX, sqrtXdag, sqrtY, sqrtYdag # クリフォード
from qulacs.gate import T, Tdag # T gate
from qulacs.gate import P0, P1 # Projection to 0,1 (not normalized)
target = 3
gate = T(target)
print(gate)
 *** gate info ***
 * gate name : T
 * target    :
 3 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : yes
 * Parametric: no
 * Diagonal  : no

Identity does not update the quantum state, but when it is included into the quantum circuit, it is counted as a gate that consumes 1 step. #### 1 qubit rotating gate Take the index of the target bit as the first argument and the rotation angle as the second argument.

[11]:
import numpy as np
from qulacs.gate import RX, RY, RZ
target = 0
angle = 0.1
gate = RX(target, angle)
print(gate)
print(gate.get_matrix())
 *** gate info ***
 * gate name : X-rotation
 * target    :
 0 : commute X
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no

[[0.99875026+0.j         0.        +0.04997917j]
 [0.        +0.04997917j 0.99875026+0.j        ]]

The definition of rotating operation is \(R_X(\theta) = \exp(i\frac{\theta}{2} X)\)です。 #### IBMQ basis gate IBMQ basis gate is a gate based on the virtual-Z decomposition defined by IBMQ’s OpenQASM.

[12]:
from qulacs.gate import U1,U2,U3
print(U3(0, 0.1, 0.2, 0.3))
 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 0 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
            (0.99875,0) (-0.0477469,-0.0147699)
 (0.0489829,0.00992933)     (0.876486,0.478826)

Definitions are: * \(U_1(\lambda) = R_Z(\lambda)\) * \(U_2(\phi, \lambda) = R_Z(\phi+\frac{\pi}{2}) R_X(\frac{\pi}{2}) R_Z(\lambda-\frac{\pi}{2})\) * \(U_3(\theta, \phi, \lambda) = R_Z(\phi+3\pi) R_X(\pi/2) R_Z(\theta+\pi) R_X(\pi/2) R_Z(\lambda)\) U3 matches the degree of freedom of any 1-qubit unitary operation.

2 qubit gate

Take the indexes of the target bit in the first and second arguments. The first argument of the CNOT gate is a control qubit. The remaining gates are symmetric operations.

[13]:
from qulacs.gate import CNOT, CZ, SWAP
control = 5
target = 2
target2 = 3
gate = CNOT(control, target)
print(gate)
gate = CZ(control, target)
gate = SWAP(target, target2)
 *** gate info ***
 * gate name : CNOT
 * target    :
 2 : commute X
 * control   :
 5 : value 1
 * Pauli     : no
 * Clifford  : yes
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no

Multi-bit Pauli operation

Multi-bit Pauli operations define a gate with arguments as a list of target qubits and a list of Pauli operators. Since the update speed of an n-qubit Pauli operation has the same order as the update cost of a 1-qubit Pauli operation, Pauli’s tensor product is often defined the gate in this form. In the Pauli operator specification, 1, 2, and 3 correspond to X, Y, and Z, respectively.

[14]:
from qulacs.gate import Pauli
target_list = [0,3,5]
pauli_index = [1,3,1] # 1:X , 2:Y, 3:Z
gate = Pauli(target_list, pauli_index) # = X_0 Z_3 X_5
print(gate)
print(gate.get_matrix())
 *** gate info ***
 * gate name : Pauli
 * target    :
 0 : commute X
 3 : commute     Z
 5 : commute X
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no

[[ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  1.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  1.+0.j  0.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j -1.-0.j]
 [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j -1.-0.j  0.+0.j]
 [ 0.+0.j  1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]
 [ 1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j  0.+0.j -1.-0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]
 [ 0.+0.j  0.+0.j -1.-0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]]

Multi-bit Pauli rotating operation

Multi-bit Pauli rotating operation rotates multi-bit Pauli operator. The multi-bit Pauli rotation becomes heavy calculation when the gate matrix is calculated intuitively, but it can be updated efficiently if it is defined in this form.

[15]:
from qulacs.gate import PauliRotation
target_list = [0,3,5]
pauli_index = [1,3,1] # 1:X , 2:Y, 3:Z
angle = 0.5
gate = PauliRotation(target_list, pauli_index, angle) # = exp(i angle/2 X_0 Z_3 X_5)
print(gate)
print(gate.get_matrix().shape)
 *** gate info ***
 * gate name : Pauli-rotation
 * target    :
 0 : commute X
 3 : commute     Z
 5 : commute X
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no

(8, 8)

Reversible circuit

Reversible circuit performs permutation operation between basis by giving a bijective function to a total number of \(2^n\) index. This is equivalent to the gate matrix being a permutation matrix. Please note that it will not work properly unless it is bijective.

[16]:
from qulacs.gate import ReversibleBoolean
def upper(val, dim):
    return (val+1)%dim
target_list = [0,1]
gate = ReversibleBoolean(target_list, upper)
print(gate)
state = QuantumState(3)
state.load(np.arange(2**3))
print(state.get_vector())
gate.update_quantum_state(state)
print(state.get_vector())
 *** gate info ***
 * gate name : ReversibleBoolean
 * target    :
 0 : commute
 1 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no

[0.+0.j 1.+0.j 2.+0.j 3.+0.j 4.+0.j 5.+0.j 6.+0.j 7.+0.j]
[3.+0.j 0.+0.j 1.+0.j 2.+0.j 7.+0.j 4.+0.j 5.+0.j 6.+0.j]

The above code moves the elements of the vector down one by one in the subspace of the qubit of interest (the bottom element moves to the top). ### State reflection This gate is \((I-2\left|a\right> \left<a \right|)\), which is defined with the quantum state \(\left|a\right>\) as an argument. This corresponds to the operation of reflecting based on the quantum state \(\left|a\right>\). This gate appears in Grover search. The number of qubits on which this gate operates must match the number of qubits in the quantum state given as an argument.

[17]:
from qulacs.gate import StateReflection
from qulacs import QuantumState
axis = QuantumState(2)
axis.set_Haar_random_state(0)
state = QuantumState(2)
gate = StateReflection(axis)
gate.update_quantum_state(state)
print("axis", axis.get_vector())
print("reflected", state.get_vector())
gate.update_quantum_state(state)
print("two reflection", state.get_vector())
axis [0.26990561+0.18885571j 0.42520294-0.30443016j 0.33036863-0.07179331j
 0.02371619+0.70933j   ]
reflected [-0.78296897+0.j          0.11454257-0.32493882j  0.15121954-0.16353884j
  0.28072431+0.37394642j]
two reflection [ 1.00000000e+00-2.77555756e-17j -5.55111512e-17+0.00000000e+00j
 -2.77555756e-17+2.77555756e-17j  0.00000000e+00-1.11022302e-16j]

General gate

This is a quantum gate with a gate matrix.

Dense matrix gate

Gate defined based on a dense matrix.

[18]:
from qulacs.gate import DenseMatrix

# 1-qubit gate
gate = DenseMatrix(0, [[0,1],[1,0]])
print(gate)

# 2-qubit gate
gate = DenseMatrix([0,1], [[1,0,0,0],[0,1,0,0],[0,0,0,1],[0,0,1,0]])
print(gate)
 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 0 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
(0,0) (1,0)
(1,0) (0,0)

 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 0 : commute
 1 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
(1,0) (0,0) (0,0) (0,0)
(0,0) (1,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (1,0)
(0,0) (0,0) (1,0) (0,0)

Sparse matrix gate

A gate defined based on a sparse matrix. If the elements are sparse enough, they can be updated faster than a dense matrix. Sparse matrices should be defined using scipy’s csc_matrix.

[20]:
from qulacs import QuantumState
from qulacs.gate import SparseMatrix
from scipy.sparse import csc_matrix
mat = csc_matrix((2,2))
mat[1,1] = 1
print("sparse matrix", mat)

gate = SparseMatrix([0], mat)
print(gate)

qs = QuantumState(2)
qs.load([1,2,3,4])
gate.update_quantum_state(qs)
print(qs.get_vector())
sparse matrix   (1, 1)  1.0
 *** gate info ***
 * gate name : SparseMatrix
 * target    :
 0 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
0 0
0 (1,0)


[0.+0.j 2.+0.j 0.+0.j 4.+0.j]

Add control bit

General gates can add a control qubit using the add_control_qubit function. You can also specify whether the target qubit has an effect when control qubit is 0 or 1.

[21]:
import numpy as np
from qulacs.gate import to_matrix_gate, X

index = 0
x_gate = X(index)
x_mat_gate = to_matrix_gate(x_gate)

# Operate only when 1st-qubit is 0
control_index = 1
control_with_value = 0
x_mat_gate.add_control_qubit(control_index, control_with_value)
print(x_mat_gate)

from qulacs import QuantumState
state = QuantumState(3)
state.load(np.arange(2**3))
print(state.get_vector())

x_mat_gate.update_quantum_state(state)
print(state.get_vector())
 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 0 : commute X
 * control   :
 1 : value 0
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
(0,0) (1,0)
(1,0) (0,0)

[0.+0.j 1.+0.j 2.+0.j 3.+0.j 4.+0.j 5.+0.j 6.+0.j 7.+0.j]
[1.+0.j 0.+0.j 2.+0.j 3.+0.j 5.+0.j 4.+0.j 6.+0.j 7.+0.j]

Operation to create a new gate from multiple gates

Gate product

Combine successively the operating quantum gates to create a new single quantum gate. This reduces access to quantum states.

[22]:
import numpy as np
from qulacs import QuantumState
from qulacs.gate import X, RY, merge

n = 3
state = QuantumState(n)
state.set_zero_state()

index = 1
x_gate = X(index)
angle = np.pi / 4.0
ry_gate = RY(index, angle)

# Create the new gate by combining gates
# First argument operates first
x_and_ry_gate = merge(x_gate, ry_gate)
print(x_and_ry_gate)
 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 1 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
 (0.382683,0)   (0.92388,0)
  (0.92388,0) (-0.382683,0)

Gate sum

You can add multiple gates to create a new gate. This is useful for making the projection of the Pauli operator \(P\) onto the +1 eigenvalue space such as $ (I + P) / 2$.

[23]:
import numpy as np
from qulacs.gate import P0,P1,add, merge, Identity, X, Z

gate00 = merge(P0(0),P0(1))
gate11 = merge(P1(0),P1(1))
# |00><00| + |11><11|
proj_00_or_11 = add(gate00, gate11)
print(proj_00_or_11)

gate_ii_zz = add(Identity(0), merge(Z(0),Z(1)))
gate_ii_xx = add(Identity(0), merge(X(0),X(1)))
proj_00_plus_11 = 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)
print(proj_00_plus_11)
 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 0 : commute
 1 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
(1,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (1,0)

 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 0 : commute
 1 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
(0.5,0)   (0,0)   (0,0) (0.5,0)
  (0,0)   (0,0)   (0,0)   (0,0)
  (0,0)   (0,0)   (0,0)   (0,0)
(0.5,0)   (0,0)   (0,0) (0.5,0)

Random unitary

Use the RandomUnitary function to sample a random unitary matrix with the Haar measure and generate a dense matrix gate.

[24]:
from qulacs.gate import RandomUnitary
target_list = [2,3]
gate = RandomUnitary(target_list)
print(gate)
 *** gate info ***
 * gate name : DenseMatrix
 * target    :
 2 : commute
 3 : commute
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : no
 * Matrix
 (-0.549679,-0.199309)    (-0.41282,0.191987)  (-0.430215,-0.393607)  (-0.152234,-0.296078)
   (0.255527,0.290869)   (0.136685,-0.154401)  (0.0998426,-0.553577) (-0.700763,-0.0097331)
   (0.056699,0.690511) (-0.488394,-0.0743973)  (0.0480177,-0.260356)    (0.439362,0.113074)
 (0.156442,-0.0611172)  (-0.695389,-0.150231)    (0.177095,0.492055) (-0.433419,-0.0657552)

Stochastic operation

Using the Probabilistic function, create the operation by giving multiple gate operations and probability distribution. If the sum of the given probability distributions is less than 1, Identity will operate with a probability less than 1.

[25]:
from qulacs.gate import Probabilistic, H, Z
distribution = [0.2, 0.2, 0.2]
gate_list = [H(0), Z(0), X(1)]
gate = Probabilistic(distribution, gate_list)
print(gate)

from qulacs import QuantumState
state = QuantumState(2)
for _ in range(10):
    gate.update_quantum_state(state)
    print(state.get_vector())
 *** gate info ***
 * gate name : Generic gate
 * target    :
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : yes

[0.+0.j 0.+0.j 1.+0.j 0.+0.j]
[0.+0.j 0.+0.j 1.+0.j 0.+0.j]
[0.+0.j 0.+0.j 1.+0.j 0.+0.j]
[ 0.+0.j -0.-0.j  1.+0.j -0.-0.j]
[0.        +0.j 0.        +0.j 0.70710678+0.j 0.70710678+0.j]
[0.        +0.j 0.        +0.j 0.70710678+0.j 0.70710678+0.j]
[ 0.        +0.j -0.        -0.j  0.70710678+0.j -0.70710678-0.j]
[ 0.70710678+0.j -0.70710678-0.j  0.        +0.j -0.        -0.j]
[0.+0.j 1.+0.j 0.+0.j 0.+0.j]
[ 0.+0.j -1.-0.j  0.+0.j -0.-0.j]

BitFlipNoise, DephasingNoise, IndependentXZNoise, DepolarizingNoise, and TwoQubitDepolarizingNoise gates are defined as stochastic gates. Probabilistic instances are generated by entering the error probabilities.

[26]:
from qulacs.gate import BitFlipNoise, DephasingNoise, IndependentXZNoise, DepolarizingNoise, TwoQubitDepolarizingNoise
target = 0
second_target = 1
error_prob = 0.8
gate = BitFlipNoise(target, error_prob) # X: prob
gate = DephasingNoise(target, error_prob) # Z: prob
gate = IndependentXZNoise(target, error_prob) # X,Z : prob*(1-prob), Y: prob*prob
gate = DepolarizingNoise(target, error_prob) # X,Y,Z : prob/3
gate = TwoQubitDepolarizingNoise(target, second_target, error_prob) # {I,X,Y,Z} \times {I,X,Y,Z} \setminus {II} : prob/15

from qulacs import QuantumState
state = QuantumState(2)
for _ in range(10):
    gate.update_quantum_state(state)
    print(state.get_vector())
[1.+0.j 0.+0.j 0.+0.j 0.+0.j]
[0.-0.j 0.+1.j 0.-0.j 0.+0.j]
[0.-0.j 0.+1.j 0.-0.j 0.+0.j]
[ 0.-0.j -0.-1.j  0.-0.j -0.-0.j]
[0.+0.j 0.-0.j 0.+0.j 0.-1.j]
[0.+0.j 0.-0.j 0.+0.j 0.-1.j]
[-1.+0.j  0.+0.j -0.+0.j  0.+0.j]
[ 0.-0.j  0.+0.j  0.+0.j -0.-1.j]
[ 0.-0.j  0.+0.j  0.+0.j -0.-1.j]
[ 0.-0.j  0.+0.j  0.+0.j -0.-1.j]

CPTP mapping

CPTP operates by giving a list of Claus operators that satisfy the completeness.

[27]:
from qulacs.gate import merge,CPTP, P0,P1

gate00 = merge(P0(0),P0(1))
gate01 = merge(P0(0),P1(1))
gate10 = merge(P1(0),P0(1))
gate11 = merge(P1(0),P1(1))

gate_list = [gate00, gate01, gate10, gate11]
gate = CPTP(gate_list)

from qulacs import QuantumState
from qulacs.gate import H,merge
state = QuantumState(2)
for _ in range(10):
    state.set_zero_state()
    merge(H(0),H(1)).update_quantum_state(state)
    gate.update_quantum_state(state)
    print(state.get_vector())
[0.+0.j 0.+0.j 1.+0.j 0.+0.j]
[0.+0.j 1.+0.j 0.+0.j 0.+0.j]
[0.+0.j 0.+0.j 1.+0.j 0.+0.j]
[1.+0.j 0.+0.j 0.+0.j 0.+0.j]
[0.+0.j 0.+0.j 0.+0.j 1.+0.j]
[1.+0.j 0.+0.j 0.+0.j 0.+0.j]
[0.+0.j 0.+0.j 1.+0.j 0.+0.j]
[0.+0.j 1.+0.j 0.+0.j 0.+0.j]
[0.+0.j 0.+0.j 0.+0.j 1.+0.j]
[0.+0.j 0.+0.j 1.+0.j 0.+0.j]

The AmplitudeDampingNoise gate is defined as CPTP-map.

[28]:
from qulacs.gate import AmplitudeDampingNoise
target = 0
damping_rate = 0.1
AmplitudeDampingNoise(target, damping_rate)
#K_0: [[1,0],[0,sqrt(1-p)]], K_1: [[0,sqrt(p)], [0,0]]
[28]:
 *** gate info ***
 * gate name : Generic gate
 * target    :
 * control   :
 * Pauli     : no
 * Clifford  : no
 * Gaussian  : no
 * Parametric: no
 * Diagonal  : yes

Instrument

Instrument is an operation to get the subscript of the Claus operator that acts randomly in addition to the general CPTP-map operation. For example, measurement on the Z basis is equivalent to applying a CPTP-map consisting of P0 and P1 and knowing which one has acted. In cppsim, this is achieved in the Instrument function, by specifying the information of the CPTP-map and the address of the classic register in which the subscript of the operated Claus operator is written.

[29]:
from qulacs import QuantumState
from qulacs.gate import merge,Instrument, P0,P1

gate00 = merge(P0(0),P0(1))
gate01 = merge(P1(0),P0(1))
gate10 = merge(P0(0),P1(1))
gate11 = merge(P1(0),P1(1))

gate_list = [gate00, gate01, gate10, gate11]
classical_pos = 0
gate = Instrument(gate_list, classical_pos)

from qulacs import QuantumState
from qulacs.gate import H,merge
state = QuantumState(2)
for index in range(10):
    state.set_zero_state()
    merge(H(0),H(1)).update_quantum_state(state)
    gate.update_quantum_state(state)
    result = state.get_classical_value(classical_pos)
    print(index, format(result,"b").zfill(2), state.get_vector())
0 01 [0.+0.j 1.+0.j 0.+0.j 0.+0.j]
1 01 [0.+0.j 1.+0.j 0.+0.j 0.+0.j]
2 01 [0.+0.j 1.+0.j 0.+0.j 0.+0.j]
3 10 [0.+0.j 0.+0.j 1.+0.j 0.+0.j]
4 00 [1.+0.j 0.+0.j 0.+0.j 0.+0.j]
5 01 [0.+0.j 1.+0.j 0.+0.j 0.+0.j]
6 11 [0.+0.j 0.+0.j 0.+0.j 1.+0.j]
7 11 [0.+0.j 0.+0.j 0.+0.j 1.+0.j]
8 00 [1.+0.j 0.+0.j 0.+0.j 0.+0.j]
9 10 [0.+0.j 0.+0.j 1.+0.j 0.+0.j]

Note that Measurement gate is defined as Instrument.

[30]:
from qulacs.gate import Measurement
target = 0
classical_pos = 0
gate = Measurement(target, classical_pos)

Adaptive operation

This is a gate that uses a function that returns a Boolean value with variable-length list of classical register values as an argument, and determines whether to perform an operation according to the conditions obtained from the classical register. Conditions can be written as python functions. A python function must be a function that takes a list of type unsigned int as an argument and returns a bool type value.

[31]:
from qulacs.gate import Adaptive, X

def func(list):
    print("func is called! content is ",list)
    return list[0]==1


gate = Adaptive(X(0), func)

state = QuantumState(1)
state.set_zero_state()

# func returns False, so X does not operate
state.set_classical_value(0,0)
gate.update_quantum_state(state)
print(state.get_vector())

# func returns True, so X operates
state.set_classical_value(0,1)
gate.update_quantum_state(state)
print(state.get_vector())
func is called! content is  [0]
[1.+0.j 0.+0.j]
func is called! content is  [1]
[0.+0.j 1.+0.j]

Operators

Pauli operator

Observables are represented as linear combinations of Pauli operators with real coefficients. PauliOperator class is a class that expresses each term with the coefficient added to the $ n $ -qubit Pauli operator. Unlike gates, quantum states cannot be updated.

Create Pauli operator and obtain state

[32]:
from qulacs import PauliOperator
coef = 0.1
s = "X 0 Y 1 Z 3"
pauli = PauliOperator(s, coef)

# Added pauli symbol later
pauli.add_single_Pauli(3, 2)

# Get the subscript of each pauli symbol
index_list = pauli.get_index_list()

# Get pauli symbols (I,X,Y,Z -> 0,1,2,3)
pauli_id_list = pauli.get_pauli_id_list()

# Get pauli coefficient
coef = pauli.get_coef()

# Create a copy of pauli operator
another_pauli = pauli.copy()

s = ["I","X","Y","Z"]
pauli_str = [s[i] for i in pauli_id_list]
terms_str = [item[0]+str(item[1]) for item in zip(pauli_str,index_list)]
full_str = str(coef) + " " + " ".join(terms_str)
print(full_str)
(0.1+0j) X0 Y1 Z3 Y3

Expected value of Pauli operator

You can evaluate the expected value and transition moment of the Pauli operator for one state.

[33]:
from qulacs import PauliOperator, QuantumState

n = 5
coef = 2.0
Pauli_string = "X 0 X 1 Y 2 Z 4"
pauli = PauliOperator(Pauli_string,coef)

# Calculate expectation value <a|H|a>
state = QuantumState(n)
state.set_Haar_random_state()
value = pauli.get_expectation_value(state)
print("expect", value)

# Calculate transition moment <a|H|b>
# The first arguments comes to the bra side
bra = QuantumState(n)
bra.set_Haar_random_state()
value = pauli.get_transition_amplitude(bra, state)
print("transition", value)
expect (0.6132028798856664+0j)
transition (0.5279713648817308-0.6040063466314933j)

General linear operators

The linear operator GeneralQuantumOperator is represented by a linear combination of the complex numbers of the Pauli operators. PauliOperator with coefficient can be added as a term with add_operator.

[34]:
from qulacs import GeneralQuantumOperator, PauliOperator, QuantumState

n = 5
operator = GeneralQuantumOperator(n)

# Pauli operator can be added
coef = 2.0+0.5j
Pauli_string = "X 0 X 1 Y 2 Z 4"
pauli = PauliOperator(Pauli_string,coef)
operator.add_operator(pauli)
# Pauli operator can also be added directly from coefficients and strings
operator.add_operator(0.5j, "Y 1 Z 4")

# Get number of terms
term_count = operator.get_term_count()

# Get number of quantum bit
qubit_count = operator.get_qubit_count()

# Get specific terms as PauliOperator
index = 1
pauli = operator.get_term(index)


# Expected value calculation <a|H|a>
## Generally not self-adjoint, can return complex value
state = QuantumState(n)
state.set_Haar_random_state()
value = operator.get_expectation_value(state)
print("expect", value)

# Transition moment calculation <a|H|b>
# The first arguments comes to the bra side
bra = QuantumState(n)
bra.set_Haar_random_state()
value = operator.get_transition_amplitude(bra, state)
print("transition", value)
expect (-0.12123883067235244+0.004874745096201068j)
transition (-0.19076380005803-0.03521147847936834j)

Generation of observables using OpenFermion

OpenFermion is a tool that gives Hamiltonian to be solved by chemical calculation in the form of Pauli operator. The output of this tool can be read in the form of a file or a string and can be used in the form of an operator.

[35]:
from qulacs.quantum_operator import create_quantum_operator_from_openfermion_file
from qulacs.quantum_operator import create_quantum_operator_from_openfermion_text

open_fermion_text = """
(-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]
"""

operator = create_quantum_operator_from_openfermion_text(open_fermion_text)
print(operator.get_term_count())
print(operator.get_qubit_count())
# In the case of create_quantum_operator_from_openfermion_file, specify the path of the file where the above is written in the argument.
15
4

Hermite operator / observable

The Hermite operator is represented by a real linear combination of the Pauli operators. Equivalent to the GeneralQuatnumOperator class, except that eigenvalues or expected values are guaranteed to be real.

The function to read and process from an external file is possible by replacing the quantum_operator with the observable, and using functions such as create_observable_from_openfermion_file create_observable_from_openfermion_text create_split_observable.

Separates operators into diagonal and off-diagonal terms

When reading an operator from a file, you can separate it into diagonal and off-diagonal components with the create_split_observable function.

[ ]:
from qulacs.observable import create_split_observable, create_observable_from_openfermion_file

# H2.txt must be placed in openfermon format beforehand.
operator = create_observable_from_openfermion_file("./H2.txt")
diag, nondiag = create_split_observable("./H2.txt")
print(operator.get_term_count(), diag.get_term_count(), nondiag.get_term_count())
print(operator.get_qubit_count(), diag.get_qubit_count(), nondiag.get_qubit_count())

Quantum Circuits

Structure of a quantum circuit

A quantum circuit is represented as a set of quantum gates. For example, a quantum circuit can be configured as follows.

[37]:
from qulacs import QuantumState, QuantumCircuit
from qulacs.gate import Z
n = 3
state = QuantumState(n)
state.set_zero_state()

circuit = QuantumCircuit(n)

# Add hadamard gate to quantum circuit
for i in range(n):
    circuit.add_H_gate(i)

# Create gate, which can also be added
for i in range(n):
    circuit.add_gate(Z(i))

# Operate quantum circuit to state
circuit.update_quantum_state(state)

print(state.get_vector())
[ 0.35355339+0.j -0.35355339-0.j -0.35355339-0.j  0.35355339+0.j
 -0.35355339-0.j  0.35355339+0.j  0.35355339+0.j -0.35355339-0.j]

Calculation and optimization of depth of quantum circuits

By combining quantum gates into a single quantum gate, the number of quantum gates can be reduced and the time required for numerical calculations can be reduced. (Of course, if the number of target qubits increases, or if a quantum gate with a dedicated function is synthesized into a quantum gate without a dedicated function, the total calculation time may not decrease. It depends.)

The code below uses the optimize function to repeat the greedy synthesis of the quantum gate of the quantum circuit until the target qubit becomes three.

[38]:
from qulacs import QuantumCircuit
from qulacs.circuit import QuantumCircuitOptimizer
n = 5
depth = 10
circuit = QuantumCircuit(n)
for d in range(depth):
    for i in range(n):
        circuit.add_H_gate(i)

# Calculate the depth (depth=10)
print(circuit.calculate_depth())

# Optimization
opt = QuantumCircuitOptimizer()
# The maximum quantum gate size allowed to be created
max_block_size = 1
opt.optimize(circuit, max_block_size)

# Calculate the depth (depth=1へ)
print(circuit.calculate_depth())
10
1

Debugging quantum circuits

When print a quantum circuit, statistical information about the gates included in the quantum circuit will be displayed.

[39]:
from qulacs import QuantumCircuit
from qulacs.circuit import QuantumCircuitOptimizer
n = 5
depth = 10
circuit = QuantumCircuit(n)
for d in range(depth):
    for i in range(n):
        circuit.add_H_gate(i)

print(circuit)
*** Quantum Circuit Info ***
# of qubit: 5
# of step : 10
# of gate : 50
# of 1 qubit gate: 50
Clifford  : yes
Gaussian  : no


Parametric Quantum Circuits

Defining a quantum circuit as a ParametricQuantumCircuit class enables some functions that are useful for optimizing quantum circuits using variational methods, in addition to the usual functions of the QuantumCircuit class.

Application 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 a quantum gate added as a parametric gate, the number of parametric gates can be extracted after the quantum circuit is configured, and the rotation angle can be changed later.

[40]:
from qulacs import ParametricQuantumCircuit
from qulacs import QuantumState
import numpy as np

n = 5
depth = 10

# construct parametric quantum circuit with random rotation
circuit = ParametricQuantumCircuit(n)
for d in range(depth):
    for i in range(n):
        angle = np.random.rand()
        circuit.add_parametric_RX_gate(i,angle)
        angle = np.random.rand()
        circuit.add_parametric_RY_gate(i,angle)
        angle = np.random.rand()
        circuit.add_parametric_RZ_gate(i,angle)
    for i in range(d%2, n-1, 2):
        circuit.add_CNOT_gate(i,i+1)

# add multi-qubit Pauli rotation gate as parametric gate (X_0 Y_3 Y_1 X_4)
target = [0,3,1,4]
pauli_ids = [1,2,2,1]
angle = np.random.rand()
circuit.add_parametric_multi_Pauli_rotation_gate(target, pauli_ids, angle)

# get variable parameter count, and get current parameter
parameter_count = circuit.get_parameter_count()
param = [circuit.get_parameter(ind) for ind in range(parameter_count)]

# set 3rd parameter to 0
circuit.set_parameter(3, 0.)

# update quantum state
state = QuantumState(n)
circuit.update_quantum_state(state)

# output state and circuit info
print(state)
print(circuit)
 *** Quantum State ***
 * Qubit Count : 5
 * Dimension   : 32
 * State vector :
   (0.0292587,0.0842756)
   (-0.0199201,0.242442)
    (0.181468,-0.104986)
  (0.0995482,-0.0233467)
     (0.16811,0.0345288)
   (0.0562338,-0.225101)
    (-0.131558,0.100174)
    (0.144257,-0.120167)
  (-0.173724,-0.0457571)
  (-0.142581,-0.0917222)
     (0.210296,0.176269)
   (0.0014233,0.0511538)
   (-0.00440611,0.26444)
   (0.0175043,0.0854573)
    (-0.197366,0.133705)
     (0.191306,0.316103)
   (-0.132828,-0.128943)
   (-0.00599114,0.03745)
  (-0.203638,-0.0249769)
   (-0.157577,0.0372862)
   (-0.132844,-0.105019)
  (-0.136556,-0.0918098)
   (-0.0584003,0.043396)
    (-0.02902,0.0680596)
   (0.0952735,-0.052631)
   (0.0738972,-0.127688)
   (0.165229,-0.0591462)
(0.0608202,-0.000749349)
  (-0.0873545,0.0887971)
    (-0.12346,0.0469901)
 (0.0344535,-0.00645451)
   (-0.216478,0.0550732)

*** Quantum Circuit Info ***
# of qubit: 5
# of step : 41
# of gate : 171
# of 1 qubit gate: 150
# of 2 qubit gate: 20
# of 3 qubit gate: 0
# of 4 qubit gate: 1
Clifford  : no
Gaussian  : no

*** Parameter Info ***
# of parameter: 151