Gradient in differentiable Hartree-Fock demo

Thanks @Raul_Guerrero, I created a simple working example that you can use to debug the original code. Please also have a look at this discussion which is a bout a similar topic. Please let me know if you have any further questions. Thanks!

from autograd import grad
import pennylane as qml
from pennylane import numpy as np

symbols = ["O", "H", "H"]
geometry = np.array([[0.028,  0.054,  0.10],
                     [0.986,  1.610, -0.10],
                     [1.855,  0.002,  0.20]], requires_grad=True)

mol = qml.qchem.Molecule(symbols, geometry)

dev = qml.device("default.qubit")

def energy(mol):
    @qml.qnode(dev, interface="autograd")
    def circuit(*args):
        
        # note that active_electrons=2, active_orbitals=2 in this example
        qml.BasisState(np.array([1, 1, 0, 0]), wires=range(4))
        
        # note that this simple circuit has only one gate
        qml.DoubleExcitation(*args[0][0], wires=[0, 1, 2, 3])
        
        # note that active_electrons=2, active_orbitals=2 in this example
        H = qml.qchem.molecular_hamiltonian(mol.symbols, mol.coordinates,
            active_electrons=2, active_orbitals=2, args=args[1:])[0]
        
        return qml.expval(H)
    return circuit

# number of zeros should match the number of circuit gates
circuit_param = [np.array([0.0], requires_grad=True)]

for n in range(36):

    args = [circuit_param, geometry]
    mol = qml.qchem.Molecule(symbols, geometry)

    g_param = grad(energy(mol), argnum = 0)(*args)
    circuit_param = circuit_param - 0.25 * g_param[0]

    forces = -grad(energy(mol), argnum = 1)(*args)
    geometry = geometry + 0.5 * forces

    print(f'n: {n}, E: {energy(mol)(*args):.8f}, Force-max: {abs(forces).max():.8f}')
1 Like