Merge pull request #28 from yangguohao/task81

【PaddlePaddle Hackathon】81 时间演化电路的性能优化

documents/单测文件.py

 from paddle_quantum.circuit import UAnsatz
+from paddle_quantum.utils import Hamiltonian,NKron, gate_fidelity,SpinOps,dagger
+from paddle_quantum.trotter import construct_trotter_circuit, get_1d_heisenberg_hamiltonian
+from paddle import matmul, transpose, trace
+import paddle
+import numpy as np
+import scipy
+from scipy import linalg
+import matplotlib.pyplot as plt
+
+def get_evolve_op(t,h): return scipy.linalg.expm(-1j * t * h.construct_h_matrix())
+
+def test(h,n):
+    t = 2
+    r = 1
+    cir = UAnsatz(n)
+    construct_trotter_circuit(cir, h, tau=t/r, steps=r) 
+    print('系统的哈密顿量为：')
+    print(h)
+    print('电路的酉矩阵与正确的演化算符之间的保真度为：%.2f' % gate_fidelity(cir.U.numpy(), get_evolve_op(t,h)))
+    print(cir)
+
+
+print('--------------test1------------')
+h1 = get_1d_heisenberg_hamiltonian(length=2, j_x=1, j_y=1, j_z=2,h_z=2 * np.random.rand(2) - 1,periodic_boundary_condition=False)# 
+test(h1,2)
+print('--------------test2------------')
+h2 = Hamiltonian([[1., 'X0, X1'], [1., 'Z2, Z3'], [1., 'Y0, Y1'], [1., 'X1, X2'], [1., 'Y2, Y3'], [1., 'Z0, Z1']])
+test(h2,4)
+print('--------------test3------------')
+h3 = Hamiltonian([ [1, 'Y0, Y1'], [1, 'X1, X0'],[1, 'X0, Y1'],[1, 'Z0, Z1']])
+test(h3,2)
+
+"""
+
+运行耗时: 130毫秒
+--------------test1------------
+([1.0, 1.0, 2.0, 0.8812020131972615, 0.7453483155128535], ['XX', 'YY', 'ZZ', 'Z', 'Z'], [[0, 1], [0, 1], [0, 1], [0], [1]])
+系统的哈密顿量为：
+1.0 X0, X1
+1.0 Y0, Y1
+2.0 Z0, Z1
+0.8812020131972615 Z0
+0.7453483155128535 Z1
+电路的酉矩阵与正确的演化算符之间的保真度为：0.99
+---------------x----Rz(6.429)----*-----------------x----Rz(-1.57)----Rz(3.525)--
+               |                 |                 |                            
+--Rz(1.571)----*----Ry(-3.85)----x----Ry(3.854)----*----Rz(2.981)---------------
+                                                                                
+--------------test2------------
+([1.0, 1.0, 1.0, 1.0, 1.0, 1.0], ['XX', 'YY', 'ZZ', 'ZZ', 'XX', 'YY'], [[0, 1], [0, 1], [0, 1], [2, 3], [1, 2], [2, 3]])
+系统的哈密顿量为：
+1.0 X0, X1
+1.0 Z2, Z3
+1.0 Y0, Y1
+1.0 X1, X2
+1.0 Y2, Y3
+1.0 Z0, Z1
+电路的酉矩阵与正确的演化算符之间的保真度为：0.51
+-------------------x--------Rz(2.429)--------*---------------------x----Rz(-1.57)-------------------------------------------------------------------------------
+                   |                         |                     |                                                                                            
+--Rz(1.571)--------*--------Ry(-3.85)--------x--------Ry(3.854)----*--------H--------*-----------------*----H---------------------------------------------------
+                                                                                     |                 |                                                        
+------*-------------------------*------------H---------------------------------------x----Rz(4.000)----x----H----Rx(1.571)----*-----------------*----Rx(-1.57)--
+      |                         |                                                                                             |                 |               
+------x--------Rz(4.000)--------x--------Rx(1.571)----------------------------------------------------------------------------x----Rz(4.000)----x----Rx(-1.57)--
+                                                                                                                                                                
+--------------test3------------
+([1.0, 1.0, 1.0, 1.0], ['YY', 'XX', 'ZZ', 'XY'], [[0, 1], [1, 0], [0, 1], [0, 1]])
+系统的哈密顿量为：
+1.0 Y0, Y1
+1.0 X1, X0
+1.0 X0, Y1
+1.0 Z0, Z1
+电路的酉矩阵与正确的演化算符之间的保真度为：0.46
+---------------x----Rz(2.429)----*-----------------x----Rz(-1.57)----H----*-----------------*--------H------
+               |                 |                 |                      |                 |               
+--Rz(1.571)----*----Ry(-3.85)----x----Ry(3.854)----*----Rx(1.571)---------x----Rz(4.000)----x----Rx(-1.57)--
+"""
+
 from paddle_quantum.utils import partial_trace,plot_state_in_bloch_sphere,partial_trace_discontiguous,NKron,plot_n_qubit_state_in_bloch_sphere
 from mpl_toolkits.mplot3d import Axes3D
 import numpy as np
@@ -19,3 +98,4 @@ state = rho
 plot_n_qubit_state_in_bloch_sphere(state,show_arrow=True)
 plot_n_qubit_state_in_bloch_sphere(cir2.run_density_matrix(),show_arrow=True)
 plot_n_qubit_state_in_bloch_sphere(cir1.run_state_vector(),show_arrow=True)
+

paddle_quantum/trotter.py

@@ -18,6 +18,7 @@ Trotter Hamiltonian time evolution circuit module
 
 from paddle_quantum.utils import Hamiltonian
 from paddle_quantum.circuit import UAnsatz
+from collections import defaultdict
 import warnings
 import numpy as np
 import re
@@ -257,13 +258,54 @@ def __add_first_order_trotter_block(circuit, tau, grouped_hamiltonian, reverse=F
     if not reverse:
         for hamiltonian in grouped_hamiltonian:
             assert isinstance(hamiltonian, Hamiltonian)
+            
+            #将原哈密顿量中相同site的XX，YY，ZZ组合到一起
+            grouped_hamiltonian = []
+            coeffs, pauli_words, sites = hamiltonian.decompose_with_sites()
+            grouped_terms_indices = []
+            left_over_terms_indices = []
+            d = defaultdict(list)
+            #合并相同site的XX,YY,ZZ
+            for term_index in range(len(coeffs)):
+                site = sites[term_index]
+                pauli_word = pauli_words[term_index]
+                for pauli in ['XX', 'YY', 'ZZ']:
+                    assert isinstance(pauli_word, str), "Each pauli word should be a string type"
+                    if (pauli_word==pauli or pauli_word==pauli.lower()):
+                        key = tuple(sorted(site))
+                        d[key].append((pauli,term_index))
+                        if len(d[key])==3:
+                            terms_indices_to_be_grouped = [x[1] for x in d[key]]
+                            grouped_terms_indices.extend(terms_indices_to_be_grouped)
+                            grouped_hamiltonian.append(hamiltonian[terms_indices_to_be_grouped])
+            #其他的剩余项
+            for term_index in range(len(coeffs)):
+                if term_index not in grouped_terms_indices:
+                    left_over_terms_indices.append(term_index)
+            if len(left_over_terms_indices):
+                for term_index in left_over_terms_indices:
+                    grouped_hamiltonian.append(hamiltonian[term_index])
+            #得到新的哈密顿量
+            res = grouped_hamiltonian[0]
+            for i in range(1,len(grouped_hamiltonian)):
+                res+=grouped_hamiltonian[i]
+            hamiltonian = res
+            
             # decompose the Hamiltonian into 3 lists
             coeffs, pauli_words, sites = hamiltonian.decompose_with_sites()
             # apply rotational gate of each term
-            for term_index in range(len(coeffs)):
-                # get the sorted pauli_word and site (an array of qubit indices) according to their qubit indices
-                pauli_word, site = __sort_pauli_word(pauli_words[term_index], sites[term_index])
-                add_n_pauli_gate(circuit, 2 * tau * coeffs[term_index], pauli_word, site)
+            term_index = 0
+            while term_index <len(coeffs):
+                if term_index+3<=len(coeffs) and \
+                len(set(y for x in sites[term_index:term_index+3] for y in x ))==2 and\
+                set(pauli_words[term_index:term_index+3])=={'XX','YY','ZZ'}:
+                    optimal_circuit(circuit,[tau*i for i in coeffs[term_index:term_index+3]],sites[term_index])
+                    term_index+=3
+                else:
+                    # get the sorted pauli_word and site (an array of qubit indices) according to their qubit indices
+                    pauli_word, site = __sort_pauli_word(pauli_words[term_index], sites[term_index])
+                    add_n_pauli_gate(circuit, 2 * tau * coeffs[term_index], pauli_word, site)
+                    term_index+=1
     # in the reverse mode, if the Hamiltonian is a single element list, reverse the order its each term
     else:
         if len(grouped_hamiltonian) == 1:
@@ -361,7 +403,30 @@ def add_n_pauli_gate(circuit, theta, pauli_word, which_qubits):
                 circuit.h(which_qubits[qubit_index])
             elif re.match(r'Y', pauli_word[qubit_index], flags=re.I):
                 circuit.rx(- PI / 2, which_qubits[qubit_index])
+                
+def optimal_circuit(circuit,theta,which_qubits):
+    r""" 添加一个优化电路，哈密顿量为'XXYYZZ'`
 
+    Args:
+        circuit (UAnsatz): 需要添加门的电路
+        theta list(tensor or float): 旋转角度需要传入三个参数
+        which_qubits (list or np.ndarray): ``pauli_word`` 中的每个算符所作用的量子比特编号
+    """
+    p = np.pi/2
+    x,y,z = theta
+    alpha = paddle.to_tensor(3*p-4*x*p+2*x,dtype='float64')
+    beta = paddle.to_tensor(-3*p+4*y*p-2*y,dtype='float64')
+    gamma = paddle.to_tensor(2*z-p,dtype='float64')
+    which_qubits.sort()
+    a,b = which_qubits
+    circuit.rz(paddle.to_tensor(p,dtype='float64'),b)
+    circuit.cnot([b,a])
+    circuit.rz(gamma,a)
+    circuit.ry(alpha,b)
+    circuit.cnot([a,b])
+    circuit.ry(beta,b)
+    circuit.cnot([b,a])
+    circuit.rz(paddle.to_tensor(-p,dtype='float64'),a)
 
 def __group_hamiltonian_xyz(hamiltonian):
     r""" 将哈密顿量拆分成 X、Y、Z 以及剩余项四个部分，并返回由他们组成的列表