<jupyter_start><jupyter_text>Tutoral for Quantum Chemistry 2020-06-11  ![jupyter](linux.jpeg)<jupyter_code>import projectq
from projectq.ops import All, Measure, QubitOperator, TimeEvolution
from projectq.ops import CNOT, H, Rz, Rx, Ry, X, Z, Y
import matplotlib.pyplot as plt
from scipy.optimize import minimize_scalar, minimize
import numpy as np<jupyter_output><empty_output><jupyter_text>small VQE demo $$H = 0.5 X_0 + 0.7 Y_0$$<jupyter_code>hamil = 0.5 * QubitOperator("X0") + 0.7 * QubitOperator("Y0")

hamil_matrix = 0.5* X.matrix + 0.7 * Y.matrix
hamil_matrix
np.linalg.eigh(hamil_matrix)
def small_vqe(theta, hamiltonian):

    energy = 0
    eng = projectq.MainEngine()
    qubits = eng.allocate_qureg(1)

    Z | qubits[0]
    
    Rx(theta[0]) | qubits[0]

    eng.flush()
    energy = eng.backend.get_expectation_value(hamiltonian, qubits)
    All(Measure) | qubits
    return energy
minimum = minimize(lambda theta: small_vqe(theta,hamil), [0])
minimum
def small_vqe(theta, hamiltonian):

    energy = 0
    eng = projectq.MainEngine()
    qubits = eng.allocate_qureg(1)

    H | qubits[0]
    
    Rz(theta[0]) | qubits[0]
#     Ry(theta[1]) | qubits[0]

    eng.flush()
    energy = eng.backend.get_expectation_value(hamiltonian, qubits)
    All(Measure) | qubits
    return energy
minimum = minimize(lambda theta: small_vqe(theta,hamil), [0.,0.])
minimum<jupyter_output><empty_output><jupyter_text>Hydrogen (H2)<jupyter_code>raw_data_table_1 = np.loadtxt("H2/H2_hamiltonians.txt")
target_energies = np.loadtxt("H2/H2_gs_energy.txt")
bond_lengths = [k[0] for k in raw_data_table_1];

plt.xlabel("length")
plt.ylabel("target energy")
plt.plot(bond_lengths, target_energies, "b.-",mew ='3')
plt.show()
def vqe(theta, hamiltonian):
    """ This is the function you need to test your parameterized quantum circuit.
    Args:
        theta (list): parameters
        hamiltonian (QubitOperator): Hamiltonian of the system
    Returns:
        energy of the wavefunction for parameters theta
    """
    
    try:
        len(theta)
    except:
        theta = [theta]
    
    # Initialize the energy
    energy = 0
    eng = projectq.MainEngine()
    qubits = eng.allocate_qureg(2)
    
    

    X | qubits[0]
    
    Rx(-np.pi/2) | qubits[0]
    Ry(np.pi/2) | qubits[1]
    CNOT |  (qubits[1],qubits[0])
    
    Rz(theta[0]) | qubits[0]

    CNOT |  (qubits[1],qubits[0])
    
    Rx(np.pi/2) | qubits[0]
    Ry(-np.pi/2) | qubits[1]
    
### ================================================
    
#     X | qubits[0]
    
#     ansatz_op = QubitOperator('X0 Y1')
#     # Apply the unitary e^{-i * ansatz_op * t}
#     TimeEvolution(theta[0], ansatz_op) | qubits
    
    
    eng.flush()
    energy = eng.backend.get_expectation_value(hamiltonian, qubits)
    All(Measure) | qubits
    return energy<jupyter_output><empty_output><jupyter_text>define some fuinctions<jupyter_code>def loss_H_H2(theta,config_H):
    
    hamiltonian = config_H[1] * QubitOperator(())  # == identity
    hamiltonian += config_H[2] * QubitOperator("Z0")
    hamiltonian += config_H[3] * QubitOperator("Z1")
    hamiltonian += config_H[4] * QubitOperator("Z0 Z1")
    hamiltonian += config_H[5] * QubitOperator("X0 X1")
    hamiltonian += config_H[6] * QubitOperator("Y0 Y1")
    
    return vqe(theta,hamiltonian)<jupyter_output><empty_output><jupyter_text>Solve VQE<jupyter_code>VQE_energies = []
VQE_thetas=[]
raw_data_hamiltonian = raw_data_table_1

for i in range(len(raw_data_hamiltonian)):
    
    # Use Scipy to perform the classical outerloop of the variational
    # eigensolver, i.e., the minimization of the parameter theta.
    
    # See documentation of Scipy for different optimizers.
    config_H = raw_data_hamiltonian[i,:]
    
    theta = 0
    minimum = minimize(loss_H_H2, theta, args = (config_H))

#     minimum = minimize_scalar(lambda theta: loss_H_H2(theta,config_H))

    VQE_energies.append(minimum.fun)
    VQE_thetas.append(minimum.x)
    
    if (i+1)%10 == 0:
        print("i=",i,"\t VQE_theta=",minimum.x)
        
VQE_thetas = np.array(VQE_thetas)
        
minimum
plt.xlabel("Bond length [Angstrom]")
plt.ylabel("Total Energy [Hartree]")
plt.title("Variational Quantum Eigensolver")
plt.plot(bond_lengths, VQE_energies, "r.-",mew ='3',
         label="VQE energy of molecular hydrogen")
plt.plot(bond_lengths, target_energies, "b-",
         label="target energy")
plt.legend()
plt.show()
plt.xlabel("Bond length [Angstrom]")
plt.ylabel("theta")
plt.title("variational parameter")
plt.plot(bond_lengths, VQE_thetas, "r.-",
         label="variational parameter of molecular hydrogen")
plt.legend()
plt.show()<jupyter_output><empty_output><jupyter_text>Let's try something else molecular LiH<jupyter_code>import projectq
from projectq.ops import All, Measure, QubitOperator, TimeEvolution, X
from projectq.ops import CNOT, H, Rz, Rx, Ry, X, Z, Y
import matplotlib.pyplot as plt
from scipy.optimize import minimize_scalar
import numpy as np
from numpy.linalg import norm
config_Hs= np.loadtxt("LiH/coeff_hams_more.txt")[:-2,:]
bond_lengths=np.loadtxt("LiH/spacings_more.txt")[:-2]
real_energies=np.loadtxt("LiH/gs_energies_more.txt")[:-2]
real_thetas=np.loadtxt("LiH/thetas_more.txt")[:-2]
n_samples,n_features=config_Hs.shape
operator_list=[(),"Z0","Z1","Z2", "Y0 Y1", "X0 X1","Y0 Y2", "X0 X2", "Z0 Z1","Z0 Z2", "Y1 Y2", "X1 X2", "Z1 Z2"]
plt.figure(0)
plt.plot(bond_lengths, real_energies, "r.-",
         label="Ground-state energy of molecular LiH")

plt.ylabel('Energy in Hartree')
plt.xlabel('Bond length in angstrom')
plt.legend()
plt.show()
def ham_LiH(config_H):
    #bond_distances.append(config_H[i,0])
    hamiltonian=config_H[0]*QubitOperator(operator_list[0])
    for j in range(1,n_features):
        hamiltonian+=config_H[j]*QubitOperator(operator_list[j])    
    return hamiltonian

hams = [ham_LiH(conf) for index, conf in enumerate(config_Hs)]
def variational_quantum_eigensolver(theta, hamiltonian):
    """
    Args:
        theta (float): variational parameter for ansatz wavefunction
        hamiltonian (QubitOperator): Hamiltonian of the system
    Returns:
        energy of the wavefunction for parameter theta
    """
    # Create a ProjectQ compiler with a simulator as a backend
    eng = projectq.MainEngine()
    # Allocate 3 qubits in state |00>
    wf = eng.allocate_qureg(3)
    # Initialize the Hartree Fock state |111>
    X | wf[0]
    X | wf[1]
    X | wf[2]
    
    Rx(-np.pi/2) | wf[0]
    Ry(np.pi/2) | wf[1]
    CNOT |  (wf[0],wf[1])
    
    Rz(theta[0]) | wf[1]

    CNOT |  (wf[0],wf[1])
    Rx(np.pi/2) | wf[0]
    Ry(-np.pi/2) | wf[1]
    
    ########
    Rx(-np.pi/2) | wf[1]
    Ry(np.pi/2) | wf[2]
    CNOT |  (wf[2],wf[1])
    
    Rz(theta[1]) | wf[1]

    CNOT |  (wf[2],wf[1])
    Rx(np.pi/2) | wf[1]
    Ry(-np.pi/2) | wf[2]

    # flush all gates
    eng.flush()
    # Calculate the energy.
    # The simulator can directly return expectation values, while on a
    # real quantum devices one would have to measure each term of the
    # Hamiltonian.
    energy = eng.backend.get_expectation_value(hamiltonian, wf)
    # Measure in order to return to return to a classical state
    # (as otherwise the simulator will give an error)
    All(Measure) | wf
    return energy

def loss_H(theta,config_H):
     
    hams = ham_LiH(config_H)
    
    return variational_quantum_eigensolver(theta,hams)
VQE_energies = []
thetas=[]
for i in range(len(bond_lengths)):
    
    # Use Scipy to perform the classical outerloop of the variational
    # eigensolver, i.e., the minimization of the parameter theta.
    
    # See documentation of Scipy for different optimizers.
    config_H = config_Hs[i,:]
#     minimum = minimize_scalar(lambda theta: ,[0.5,0.5])
    minimum = minimize(lambda theta: loss_H(theta,config_H), [0.5,0.5],method='nelder-mead',
               options={'xtol': 1e-5})
    VQE_energies.append(minimum.fun)
    thetas.append(minimum.x)
    if (i+1)%10 == 0:
        print("i=",i,"\t theta=",minimum.x)

# np.savetxt("VQE_energies.txt",np.array(VQE_energies))
%matplotlib inline
plt.figure()

plt.plot(bond_lengths, real_energies, "b-",
         label="Ground-state energy of molecular LiH")
plt.plot(bond_lengths, VQE_energies, "r.",
         label="VQE energy of molecular LiH")

plt.ylabel('Energy in Hartree')
plt.xlabel('Bond length in angstrom')
plt.legend()
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1409204423831609, 'fraction_numerical': 0.02081103579498157, 'mean_word_length': 3.711484593837535, 'pattern_counts': {'":': 0, '<': 25, '<?xml version=': 0, '>': 27, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '164859', 'license': 'Apache-2.0', 'n_tokens_mistral': 3386, 'n_tokens_neox': 3173, 'n_words': 666, 'path': '.ipynb_checkpoints/VQE_tutorial_1-checkpoint.ipynb', 'repo_name': 'yintaophy/Quantum_meetup_20200611', 'size': '8409'}