Help to find higher energy states of molecules

In this tutorial energy value of LiH , H2 is calculated . For this calculation molecule geometry is given there. If I want to run this for other materials like Be2, Be3 or B2 or others . Where can i get these geometry values of other molecules. Is there any data set for this?

Hey @Kaushik!

You’re in luck :sunglasses:. PennyLane has a datasets module: Quantum Data - PennyLane. There, you can find all sorts of molecular datasets. If there’s a molecule that PennyLane doesn’t have a dataset for, odds are you can get quick numbers from NIST (e.g., here: CCCBDB Experimental Geometry Data). Here’s the data for Be2: https://cccbdb.nist.gov/expgeom2x.asp

Let me know if that helps!

1 Like

Here I am facing some problems. When I am running this program for Be2 it is showing different energy each time .

Hey @Kaushik,

This could be due to a few things, but it’s hard to say without seeing what code you’re running (if you could attach a minimal example that would be lovely). My guess is that you’re randomly initializing your ansatz parameters. Each time you run your code, your parameters are different, which will have a huge impact on the final result :slight_smile:. You can get around this by setting a random seed.

Let me know if this helps!

from pennylane import numpy as np
import networkx as nx
import matplotlib.pyplot as plt
from math import pi

import pennylane as qml
from pennylane.templates import RandomLayers
#from pennylane.init import random_layers_uniform # modified - commented
from pennylane import broadcast
from pennylane import expval, var, device

symbols = ["Be", "Be"]
coordinates = np.array([0.0, 0.0, -1.230, 0.0, 0.0, 1.230])

h, qubits =  qml.qchem.molecular_hamiltonian(
        symbols,
        coordinates,
        charge=0,
        mult=1,
        basis='sto-3g',
        active_electrons=2,
        active_orbitals=2
)

print("Number of qubits = ", qubits)
print("The Hamiltonian is ", h)



energies = np.zeros(2)


# Get the number of qubits in the Hamiltonian
num_qubits = len(h.wires)

# Initialize the device
dev = qml.device('default.qubit', wires=num_qubits)


# Ansatz
def ansatz(params, wires, state_idx=0):
    # Need to prepare a different orthogonal state each time
    qml.PauliX(wires=wires[state_idx])
    qml.templates.StronglyEntanglingLayers(params, wires=wires)
    
#single_cost = qml.ExpvalCost(ansatz, h, dev) # modified - commented

# modified - added
@qml.qnode(dev)
def single_cost(params, state_idx=0):
  ansatz(params,wires=range(num_qubits),state_idx=state_idx)
  return qml.expval(h)

w = np.arange(num_qubits, 0, -1)


# The full cost - computes single_cost for each starting state
def total_cost(params):
    cost = 0
    for state_idx in range(num_qubits):
        cost += w[state_idx] * single_cost(params, state_idx=state_idx)
    return cost


# Set up a cost function and optimizer, and run the SSVQE
opt = qml.AdamOptimizer(stepsize=0.05)
max_iterations = 20 # modified from 200
costs = []

# Initial parameters for the ansatz
num_layers = 8
params = np.random.uniform(low=0, high=2*np.pi, size=(num_layers, num_qubits, 3))


# Optimize!
H2_callback_energies = []
H2_callback_energies_single_0 = []
H2_callback_energies_single_1 = []

for _ in range(max_iterations):
    params, prev_energy = opt.step_and_cost(total_cost, params)
    energy = total_cost(params)
    energy_ground_state = single_cost(params, state_idx=0)
    energy_first_excited_state = single_cost(params, state_idx=1)
    H2_callback_energies_single_0.append(energy_ground_state)
    H2_callback_energies_single_1.append(energy_first_excited_state)
    H2_callback_energies.append(energy)
    
    
    
# After optimization, get the energies for the original starting states
for state_idx in range(2):
    energies[state_idx] = single_cost(params, state_idx=state_idx)
    
    
energies = ",".join([str(E) for E in energies])


energies

This is the code I am using and I am getting various energy each time i run it

Hey @Kaushik,

Indeed it’s the lack of a random seed. If you’re looking for the same result every time, you can do the following:

np.random.seed(1234)
params = np.random.uniform(low=0, high=2*np.pi, size=(num_layers, num_qubits, 3))

You can change the integer 1234 to something different and get different results :slight_smile:. If you want more information on how random seeds work, the wikipedia entry isn’t all that bad (quick read as well): Random seed - Wikipedia.

But, I’d say it’s not necessary to dig deep into this if you’re just after a deterministic result from your code :slight_smile: