Molecular Hamiltonian with PennyLane and Qiskit

PennyLane Help
1

Hi!

I am trying to compare VQE when used in PennyLane and Qiskit with the H2 molecule, but I’m finding trouble in obtaining the same qubit Hamiltonian with both libraries.

What I am doing in PennyLane is this:

from pennylane import numpy as np
import pennylane as qml

symbols = [“H”, “H”]
coordinates = np.array([0.0, 0.0, 0.0, 0.0, 0.0, 0.735])
H, qubits = qml.qchem.molecular_hamiltonian(symbols, coordinates)
print("Number of qubits = ", qubits)
print("The Hamiltonian is ", H)

which gives me this Hamiltonian:

Number of qubits = 4
The Hamiltonian is (-0.4636962339610332) [Z2]
+ (-0.4636962339610332) [Z3]
+ (0.2377881241080374) [Z1]
+ (0.23778812410803746) [Z0]
+ (0.7905416330599709) [I0]
+ (0.1407985649568116) [Z0 Z2]
+ (0.1407985649568116) [Z1 Z3]
+ (0.18181634392642693) [Z0 Z3]
+ (0.18181634392642693) [Z1 Z2]
+ (0.18466765822807416) [Z0 Z1]
+ (0.19192132833854414) [Z2 Z3]
+ (-0.041017778969615365) [Y0 Y1 X2 X3]
+ (-0.041017778969615365) [X0 X1 Y2 Y3]
+ (0.041017778969615365) [Y0 X1 X2 Y3]
+ (0.041017778969615365) [X0 Y1 Y2 X3]

In Qiskit, I am doing this:

from qiskit_nature.drivers import Molecule
from qiskit_nature.drivers.second_quantization import ElectronicStructureMoleculeDriver, ElectronicStructureDriverType
from qiskit_nature.problems.second_quantization import ElectronicStructureProblem
from qiskit_nature.converters.second_quantization import QubitConverter
from qiskit_nature.mappers.second_quantization import JordanWignerMapper

molecule = Molecule(geometry=[[‘H’, [0., 0., 0.]],
[‘H’, [0., 0., 0.735]]],
charge=0, multiplicity=1)

driver = ElectronicStructureMoleculeDriver(molecule, basis=‘sto3g’, driver_type=ElectronicStructureDriverType.PYSCF)
es_problem = ElectronicStructureProblem(driver)
qubit_converter = QubitConverter(JordanWignerMapper())
second_q_op = es_problem.second_q_ops()
qubit_op = qubit_converter.convert(second_q_op[0])
print(“Qubit Hamiltonian”)
print(qubit_op)

This gives me the following Hamiltonian:
Qubit Hamiltonian
-0.8105479805373279 * IIII
+ 0.1721839326191554 * IIIZ
- 0.22575349222402372 * IIZI
+ 0.17218393261915543 * IZII
- 0.2257534922240237 * ZIII
+ 0.12091263261776627 * IIZZ
+ 0.16892753870087907 * IZIZ
+ 0.045232799946057826 * YYYY
+ 0.045232799946057826 * XXYY
+ 0.045232799946057826 * YYXX
+ 0.045232799946057826 * XXXX
+ 0.1661454325638241 * ZIIZ
+ 0.1661454325638241 * IZZI
+ 0.17464343068300453 * ZIZI
+ 0.12091263261776627 * ZZII

As you can see, they are very different, and so are the ground state energies that I obtain when running VQE on them.

I would really appreciate it if you could shed some light here. What am I doing wrong?

Thanks in advance!

1 Like
2

Hi @combarro,

In the case of PennyLane I can see that the coordinates aren’t right. They should be the nuclear coordinates in atomic units so maybe the units you used are wrong. The right coordinates are:

coordinates = np.array([0.0, 0.0, -0.6614, 0.0, 0.0, 0.6614])

Our introductory demo to VQE uses the H2 molecule too. You can use it to compare whether you’re getting the same results as us.

Regarding the Qiskit case I would check with them what units you should use for the nuclear coordinates.

I hope this helps!

3

Thank you, Catalina. My main interest is in the comparison of PennyLane and Qiskit for this task. That is why I used the [0.0, 0.0, 0.0, 0.0, 0.0, 0.735] coordinates (those were the values in the Qiskit tutorial). I didn’t think there were some “right coordinates” in the sense that, in principle, you can compute the energy for different distances and positions of the H atoms (for instance, if you are computing the dissociation profile).

In any case, I will check what units Qiskit is using.

Thanks!

4

Oh I see. I think that Qiskit may be using Angstrom. We use atomic units (Bohr). So you writing coordinates = np.array([0.0, 0.0, 0.0, 0.0, 0.0, 1.3889487]) in PennyLane should be equivalent to [0.0, 0.0, 0.0, 0.0, 0.0, 0.735] in Qiskit.

5

Thank you, Catalina. I had already tried that, with no luck. Then I noticed that Qiskit’s Hamiltonian is only for the electronic structure, while PannyLane’s accounts for the energy of the whole molecule (if I am not mistaken).

6

Hi @combarro,

Oh that could be it. I’m not sure how it is in Qiskit but in PennyLane we do account for the total energy of the molecule. You can learn more about building molecular Hamiltonians with PennyLane in this demo :slight_smile:

7

Yes, I think that it is the difference. Thanks!

8

I’m glad I could help!

Enjoy using PennyLane!

9

Hi @CatalinaAlbornoz,
I am also trying to follow demo referenced above but can’t figure out from where we take the following coordinates?

symbols = ["H", "H"]
coordinates = np.array([0.0, 0.0, -0.6614, 0.0, 0.0, 0.6614])
  1. According to DatasetsSource on github we have 0.3710 Angstrom that equals to 0.3710 / 0.5292 = 0.7010 atomic units (Bohr).
2
H2 0
CCSD(T)/aug-cc-pVQZ
0.7420
H 0.3710 0.0 0.0
H -0.3710 0.0 0.0
  1. From another source Quantum Datasets bond length also equals 0.742 so 0.742 / 2 = 0.371 Angstrom.
qml.data.load("qchem", molname="H2", bondlength=0.742, basis="STO-3G")
  1. The same values were provide in documentation by reading values from h2.xyz but as we saw above they are different.

Let’s take as a source of true 0.742 Angstrom then by formula 1 bohr - 0.53 Angstrom 0.742 / 0.53 = 1.4 atomic units and it’s very close to number that was mentioned here

So you writing coordinates = np.array([0.0, 0.0, 0.0, 0.0, 0.0, 1.3889487]) in PennyLane

but not to 0.6614 * 2 = 1.3228. So interesting how we get 1.3889487 then?

Probably I have missed something obvious, is there some reference to the source from where 0.6614 was taken?

Thanks in advance!

10

Hey @TerraVenil! Welcome to the forum :slight_smile:

I wouldn’t read too much into the coordinates listed there. They’re just placeholder coordinates to get the demo rolling!

11

Hey hey @isaacdevlugt! Thanks for welcoming me! :raised_hands:

According to the comment above

They should be the nuclear coordinates in atomic units so maybe the units you used are wrong. The right coordinates are:
coordinates = np.array([0.0, 0.0, -0.6614, 0.0, 0.0, 0.6614])

so it’s not just a placeholder but the right coordinates of the atoms.

12

I think that conversation above is being taken a little bit out of context. Those coordinates (coordinates = np.array([0.0, 0.0, -0.6614, 0.0, 0.0, 0.6614])) are kinda close to the exact, true values of the coordinates of both hydrogen atoms. With less computationally-costly basis sets, if you try to find the equilibrium H-H bond length it might come out to be close to that.