Tangelo VQE: Custom Ansatz and qubit Hamiltonian Tutorial

Tangelo comes packaged with an implementation of several standard ansatz circuits for the user to take advantage of. In this tutorial, we’ll explore how you can incorporate the built in VQESolver into your own workflow, by introducing a user-defined custom ansatz circuit and/or a qubit Hamiltonian. We’ll base our work here on the VQESolver class, and take advantage of tools readily available through Tangelo.

# Installation of tangelo if not already installed.
try:
    import tangelo
except ModuleNotFoundError:
    !pip install git+https://github.com/goodchemistryco/Tangelo.git@develop --quiet
import numpy as np

from tangelo.algorithms import VQESolver, FCISolver
from tangelo import SecondQuantizedMolecule
from tangelo.toolboxes.ansatz_generator.ansatz import Ansatz
from tangelo.toolboxes.qubit_mappings.statevector_mapping import get_reference_circuit
from tangelo.toolboxes.qubit_mappings.mapping_transform import get_qubit_number, fermion_to_qubit_mapping
from tangelo.linq import Circuit, Gate

Hardware Efficient Ansatz

For our example, we’re going to implement the so-called Hardware Efficient Ansatz (HEA), developed by Kandala et al at IBM. In this ansatz, a circuit is constructed with repeated layers of a simple structure. Each layer consists of entangling gates (e.g. CNOT or CZ) which couple neighbouring qubits, followed by a series of Euler rotations carried out as single-qubit rotations \(e^{i\theta_i^1 Z_i}e^{i\theta_i^2 X_i}e^{i\theta_i^3 Z_i}\). We’ll start by initializing our Ansatz class, and then fill in the functionality required to implement this ansatz in VQESolver.

To construct our HEA ansatz, we’re going to make use of the three helper functions defined here. The first will go through a register of qubits and add a layer of Euler-rotations as prescribed above. The second adds two columns of alternating CNOT gates, establishing long-range entanglement. The third brings these together into a sequence of alternating entanglers and Euler rotations.

def EulerCircuit(n_qubits):
    """Construct a circuit applying an Euler Z-X-Z rotation to each qubit."""
    circuit = Circuit()
    for target in range(n_qubits):
        circuit.add_gate(Gate("RZ" , target, parameter=0.0, is_variational=True))
        circuit.add_gate(Gate("RX", target, parameter=0.0, is_variational=True))
        circuit.add_gate(Gate("RZ", target, parameter=0.0, is_variational=True))
    return circuit

def EntanglerCircuit(n_qubits):
    """Construct a circuit applying two columns of staggered CNOT gates to all qubits
     and their neighbours"""
    circuit = Circuit()
    for ii in range(n_qubits//2):
        circuit.add_gate(Gate("CNOT", control=2*ii, target=2*ii + 1))
    for ii in range(n_qubits//2 - 1):
        circuit.add_gate(Gate("CNOT", control=2*ii + 1, target=2*(ii+1)))
    return circuit

def HEACircuit(n_qubits, n_layers):
    """Construct a circuit consisting of alternating sequence of Euler rotations and entanglers"""
    circuit = EulerCircuit(n_qubits)
    for ii in range(n_layers):
        circuit += EntanglerCircuit(n_qubits)
        circuit += EulerCircuit(n_qubits)
    return circuit

Ansatz Class

In the VQESolver, we are expecting an instance of an abstract Ansatz class, which will be responsible for constructing the variational circuit we use to minimize the energy of our problem. To build up our own Ansatz class, we’ll require the following.

  1. init: an initialization function to instantiate the class.
  2. set_var_params: initialize the variational circuit parameters
  3. update_var_params: update the parametric gates in the circuit
  4. prepare_reference_state: get fixed circuit for initializing the reference, e.g. HF state.
  5. build_circuit: instantiate the variational circuit object

Below, we’re going to type out the entire class as we will use it. This is a lot of code in one place. So afterwards, we’ll break it down into each of the relevant member methods.

class HEA(Ansatz):

    def __init__(self, n_spinorbitals, n_electrons, n_layers, mapping='jw'):

        self.n_spinorbitals = n_spinorbitals
        self.n_qubits = get_qubit_number(mapping, n_spinorbitals)
        self.n_electrons = n_electrons
        #number of layers of repeated entangler + Euler rotations
        self.n_layers = n_layers
        
        #specify fermion-to-qubit mapping (required for the initial reference state)
        self.mapping = mapping
        
        #Each layer has 3 variational parameters per qubit, and one non-variational entangler
        #There is an additional layer with no entangler.
        self.n_var_params = self.n_qubits * 3 * (self.n_layers + 1)

        self.var_params = None
        self.circuit = None

    def set_var_params(self, var_params=None):
        """Set initial variational parameter values"""
        if var_params is None:
            var_params = np.random.random(self.n_var_params)
        elif isinstance(var_params, str) and var_params == "ones":
            var_params = np.ones(self.n_var_params, dtype=float)
        elif len(var_params) != self.n_var_params:
            raise ValueError('Invalid number of parameters.')
        self.var_params = var_params
        return var_params

    def update_var_params(self, var_params):
        """Update variational parameters (done repeatedly during VQE)"""
        for param_index in range(self.n_var_params):
            self.circuit._variational_gates[param_index].parameter = var_params[param_index]
    
    def prepare_reference_state(self):
        """Prepare a circuit generating the HF reference state."""
        return get_reference_circuit(n_spinorbitals=self.n_spinorbitals, n_electrons=self.n_electrons,mapping=self.mapping)

    def build_circuit(self, var_params=None):
        """Construct the variational circuit to be used as our ansatz."""
        self.var_params = self.set_var_params(var_params)

        reference_state_circuit = self.prepare_reference_state()
        hea_circuit = HEACircuit(self.n_qubits, self.n_layers)

        if reference_state_circuit.size != 0:
            self.circuit = reference_state_circuit + hea_circuit
        else:
            self.circuit = hea_circuit
        return self.circuit

Very briefly, we’ll go through the member methods required to construct an Ansatz class. These code blocks duplicate the code above. We emphasize here that these member functions can be really as simple or as elaborate as you like.

Let’s start with set_var_params. We’re going to do something very basic and just force this to be a random or all ones numpy array. We add some error handling in case the number of parameters is incompatible with the number of variational gates in the ansatz circuit. Have a look at the implementation of UCCSD to see how you can make this more fancy and interesting.

def set_var_params(self, var_params=None):
    if var_params is None:
        var_params = np.random.random(self.n_var_params)
    elif isinstance(var_params, str) and var_params == "ones":
        var_params = np.ones(self.n_var_params, dtype=float)
    elif len(var_params) != self.n_var_params:
        raise ValueError('Invalid number of parameters.')
    self.var_params = var_params
    return var_params

Next, we’ll implement update_var_params, where the circuit is updated with a new batch of variational parameters. The tangelo.linq Circuit class keeps a record of the variational gates in the circuit, making this update very straightforward, and avoids having to rebuild the circuit from scratch. All variational gates in the circuit are updated as per the var_params argument.

def update_var_params(self, var_params):
    for param_index in range(self.n_var_params):
        self.circuit._variational_gates[param_index].parameter = var_params[param_index]

We’ll use the methods from the qubit_mappings toolbox to construct a Hartree-Fock reference state. This will just generate a circuit with an X-gate applied to each qubit which we want to begin in the \(|1\rangle\) state.

def prepare_reference_state(self):
    circuit = get_reference_circuit(n_spinorbitals=self.n_spinorbitals, n_electrons=self.n_electrons, mapping=self.mapping)
    return circuit

Finally, we’ll implement the build_circuit method. As compared to the three others here, this is really the only method in the present case that requires much effort–everything else above has followed pretty boilerplate code. For this, we’re just going to alternate between entanglers and Euler rotations, using the HEACircuit helper function we defined earlier. We then combine this with the Hartree Fock reference circuit. In the event that no qubits are instantiated as \(|1\rangle\), we skip this empty reference circuit.

def build_circuit(self, var_params=None):
    """Construct the variational circuit to be used as our ansatz."""
    self.var_params = self.set_var_params(var_params)

    reference_state_circuit = self.prepare_reference_state()
    hea_circuit = HEACircuit(self.n_qubits, self.n_layers)

    if reference_state_circuit.size != 0:
        self.circuit = reference_state_circuit + hea_circuit
    else:
        self.circuit = hea_circuit
    return self.circuit

HEA-VQE on an H\(_2\)-dimer

With the Ansatz so defined, we’re ready to go ahead and build our VQE solver class, and run a calculation on a molecule of interest. I’m going to use the SecondQuantizedMolecule class to build a hydrogen dimer and run VQE on it using my custom ansatz. In the next section, I will also show how to run VQE with a custom qubit Hamiltonian that did not require passing a molecule.

Using a molecule object as input

H2 = [('H', (0, 0, 0)),('H', (0, 0, 0.74137727))]
mol_H2 = SecondQuantizedMolecule(H2, q=0, spin=0, basis="sto-3g")

With this molecule prepared, we’re ready to instantiate our ansatz, and feed it into VQE. I’ll access details of the molecule required to build our ansatz circuit (i.e. number of spin-orbitals and number of electrons) from the molecule object. Feel free to change the number of layers in the circuit, and explore how this changes VQE results, and timing.

n_spinorbitals = mol_H2.n_active_sos
n_electrons = mol_H2.n_active_electrons
hea_layers = 4
HEA_ansatz = HEA(n_spinorbitals=n_spinorbitals, n_electrons=n_electrons, n_layers=hea_layers)

Finally, we can instantiate the VQESolver, and run. We have decided to initialize all the variational parameters to ones through an option supported in the set_var_params of our HEA ansatz class.

vqe_options = {"molecule": mol_H2, "qubit_mapping": 'JW', 'ansatz': HEA_ansatz, "initial_var_params": "ones"}

HEA_VQE = VQESolver(vqe_options)
HEA_VQE.build()
HEA_VQE.simulate()
-1.1372335844997628

Using a custom qubit operator as input

VQESolver can also directly take as input a custom qubit operator instead of computing a Hamiltonian starting from a molecule. This situation can be relevant in the case where one wants to perform VQE using a qubit Hamiltonian: * That does not correspond to a molecular system; * That has been tailored by the user (to reduce complexity or to study something specific); * That is too expensive to recompute.

As written in the previous section, it is possible to get the same result by providing VQESolver with a qubit Hamiltonian instead of a molecule. First, let’s store the optimal variational parameters to use them as a starting point for comparison.

Next step is to generate the qubit Hamiltonian. With the help of MolecularData and fermion_to_qubit_mapping, we are able to generate a fermionic Hamiltonian and transform it into operations doable on a quantum computer.

fermionic_hamiltonian_H2 = mol_H2.fermionic_hamiltonian
qubit_hamiltonian_H2 = fermion_to_qubit_mapping(fermionic_hamiltonian_H2, mapping="jw", n_spinorbitals=n_spinorbitals, n_electrons=n_electrons)
print(qubit_hamiltonian_H2)
(-0.09883484730799605+0j) [] +
(-0.0453218839181063+0j) [X0 X1 Y2 Y3] +
(0.0453218839181063+0j) [X0 Y1 Y2 X3] +
(0.0453218839181063+0j) [Y0 X1 X2 Y3] +
(-0.0453218839181063+0j) [Y0 Y1 X2 X3] +
(0.17120123806595933+0j) [Z0] +
(0.16862327595071597+0j) [Z0 Z1] +
(0.1205461274055685+0j) [Z0 Z2] +
(0.1658680113236748+0j) [Z0 Z3] +
(0.1712012380659593+0j) [Z1] +
(0.1658680113236748+0j) [Z1 Z2] +
(0.1205461274055685+0j) [Z1 Z3] +
(-0.22279639651093203+0j) [Z2] +
(0.17434948757007063+0j) [Z2 Z3] +
(-0.22279639651093205+0j) [Z3]

Notes: Users can directly retrieve the qubit Hamiltonian object used in VQESolver by accessing its qubit_hamiltonian attribute.

The dihydrogen molecule ground state energy is computed again from this qubit Hamiltonian (which is identical to the one used in the first part computed from the molecule).

vqe_alternative_options = {"qubit_hamiltonian": qubit_hamiltonian_H2, 'ansatz': HEA_ansatz, "initial_var_params": "ones"}

HEA_VQE_HAMILTONIAN = VQESolver(vqe_alternative_options)
HEA_VQE_HAMILTONIAN.build()
HEA_VQE_HAMILTONIAN.simulate()
-1.1372335844997628

The behaviour of VQE is identical to the previous one. Regarding chemical accuracy, how well did we do here? Let’s compare against Hartree Fock and FCI.

energy_fci = FCISolver(mol_H2).simulate()
energy_hf = mol_H2.mf_energy
energy_vqe = HEA_VQE.optimal_energy
energy_vqe_hamiltonian = HEA_VQE_HAMILTONIAN.optimal_energy
print(f'FCI ENERGY: {energy_fci :.7f} Ha')
print(f'HF ENERGY: {energy_hf :.7f} Ha')
print(f'HEA-VQE ENERGY (from molecule): {energy_vqe :.7f} Ha')
print(f'HEA-VQE ENERGY (from qubit Hamiltonian): {energy_vqe_hamiltonian :.7f} Ha')
FCI ENERGY: -1.1372704 Ha
HF ENERGY: -1.1166856 Ha
HEA-VQE ENERGY (from molecule): -1.1372336 Ha
HEA-VQE ENERGY (from qubit Hamiltonian): -1.1372336 Ha

Using a custom circuit as input

Lastly, it is possible to pass a pre-built arbitrary circuit with variational gates as input to VQESolver. In that case, VQESolver will optimize the variational parameters of these gates in order to minimize the cost function, which still requires a qubit Hamiltonian. This situation can be relevant to users who wish to apply variational approaches to arbitrary circuits or a circuit they got from other collaborators.

circuit = HEA_ansatz.build_circuit()
print(circuit)
Circuit object. Size 74 

X         target : [0]   
X         target : [1]   
RZ        target : [0]   parameter : 0.0     (variational)
RX        target : [0]   parameter : 0.0     (variational)
RZ        target : [0]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RX        target : [1]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RX        target : [2]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
RX        target : [3]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
CNOT      target : [1]   control : [0]   
CNOT      target : [3]   control : [2]   
CNOT      target : [2]   control : [1]   
RZ        target : [0]   parameter : 0.0     (variational)
RX        target : [0]   parameter : 0.0     (variational)
RZ        target : [0]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RX        target : [1]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RX        target : [2]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
RX        target : [3]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
CNOT      target : [1]   control : [0]   
CNOT      target : [3]   control : [2]   
CNOT      target : [2]   control : [1]   
RZ        target : [0]   parameter : 0.0     (variational)
RX        target : [0]   parameter : 0.0     (variational)
RZ        target : [0]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RX        target : [1]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RX        target : [2]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
RX        target : [3]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
CNOT      target : [1]   control : [0]   
CNOT      target : [3]   control : [2]   
CNOT      target : [2]   control : [1]   
RZ        target : [0]   parameter : 0.0     (variational)
RX        target : [0]   parameter : 0.0     (variational)
RZ        target : [0]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RX        target : [1]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RX        target : [2]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
RX        target : [3]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
CNOT      target : [1]   control : [0]   
CNOT      target : [3]   control : [2]   
CNOT      target : [2]   control : [1]   
RZ        target : [0]   parameter : 0.0     (variational)
RX        target : [0]   parameter : 0.0     (variational)
RZ        target : [0]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RX        target : [1]   parameter : 0.0     (variational)
RZ        target : [1]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RX        target : [2]   parameter : 0.0     (variational)
RZ        target : [2]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)
RX        target : [3]   parameter : 0.0     (variational)
RZ        target : [3]   parameter : 0.0     (variational)

The dihydrogen energy is recomputed for the third time by using its qubit Hamiltonian and the abstract circuit provided. We are here optimizing the variational parameters of the pre-built circuit, in an attempt to return the minimal energy for this ansatz.

vqe_second_alternative_options = {"qubit_hamiltonian": qubit_hamiltonian_H2, "ansatz": circuit, "initial_var_params": "random"}

HEA_VQE_CIRCUIT = VQESolver(vqe_second_alternative_options)
HEA_VQE_CIRCUIT.build()
HEA_VQE_CIRCUIT.simulate()
-1.1372416894541575

Note: Users may also give VQESolver a SecondQuantizedMolecule with a custom circuit. The energy at the end is the very same as the previous calculations with a molecule + custom ansatz and a qubit Hamiltonian.

energy_vqe_circuit = HEA_VQE_CIRCUIT.optimal_energy
print(f'HEA-VQE ENERGY (from a circuit): {energy_vqe_circuit :.7f} Ha')
HEA-VQE ENERGY (from a circuit): -1.1372417 Ha

Conclusion

In this tutorial, we’ve seen how to implement a custom ansatz circuit for VQE using the tools from Tangelo and how to combine it with a built-in or a custom qubit Hamiltonian. Hopefully, this gives some impression of how this platform is designed to help users construct their own workflows easily, focusing on the specific issues they are interested in studying without the distraction of building the supporting framework from scratch.