CO2 active electrons & orbitals

I’m a high school student working on a science fair project and trying to simulate CO2 using PennyLane. There is a total of 22 electrons (C = 6, 2xO = 2x8) of which 16 are valence electrons (C = 4, 2xO = 2x6) and I think 15 molecular orbitals (C = 5, 2xO = 2x5). I have two questions: first - do I have the molecular orbital number correct, and second - is there a way to select an appropriate number of active electrons and active orbitals and do these need to be the same number? Is it just a case of determining how many qubits you have available (i.e. if I had 16 qubits then I’d choose 8 active orbitals etc)?

Hi @JoshLipman. Thanks for using PennyLane and sorry for the late reply.

The total number of molecular orbitals depends on the type of the basis set you use. For a minimal basis set such as STO-3G, where you have one basis function per each atomic orbital, you will have 15 molecular orbitals (2 S and 3 P orbitals for each element) and 15 * 2 = 30 spin-orbitals which means that you will need 30 qubits. Now, you can select some of the electrons and orbitals by defining an active space. Regarding your questions:

  1. You are correct, the total number of molecular orbitals is 15 if you are using a minimal basis such as STO-3G. You can find more info here.

  2. Defining a chemically meaningful active space depends on the problem and usually requires some experience. However, in practice, the only restriction is that the number of active spin-orbitals should be larger than the number of active electrons. For instance, you can have 2 active electrons and 2 active orbitals (i.e. 4 spin-orbitals) and you will need only 4 qubits for simulating CO2.

  3. Yes, you can define the number of active orbitals based on the number of available qubits.

I hope these help. Please feel free to let me know if you have any other questions.

Hi @JoshLipman,

I just wanted to add to @sjahangiri reply that the tutorial Building molecular Hamiltonians may be also useful to get familiar with this concepts and terminology.

Thank you!

Thanks very much for the great feedback. I’m simulating H2O, and I’ve created an adaptive circuit and optimized the geometry following the approach in the tutorial. It works well for 4 active electrons and 4 active orbitals, but when I step up to 6 active electrons and 5 active orbitals the bond length explodes, as you can see in cell 31 of my code. I’ve also tried this in google colab and got the same result, so I know it’s not an issue with my laptop environment. I’d welcome any suggestions, thanks.

Hi @JoshLipman. Nice work on building the adaptive circuit!

The issue is very likely due to an inconsistency in the finite difference Hamiltonian derivatives caused by pyscf. You can use psi4 to build the molecular Hamiltonian by specifying package='psi4' when you build molecular_hamiltonian. This will call for openfermionpsi4 to compute the Hamiltonian.

You might also consider using a smaller threshold for optimising the geometry, or use a smaller molecule such as BeH2 for the debuging process.

You can also have a look at the differentiable Hartree-Fock solver of PennyLane which allows computing the forces and circuit parameter gradients with the methods of automatic differentiation. It also allows optimising the basis set parameters. You can find an example here.

Please let us know if you face any other questions or face any issues.

Hi sjahangiri,

Thanks so much. I tried using psi4 but got the following error. I attempted install from source, through conda and through pip, but the same error each time. I guess I’m missing something simple, but would be grateful for any pointers.

Thanks!

In [8]:

active_electrons = 2
active_orbitals = 3
​
H, qubits = qchem.molecular_hamiltonian(
    symbols,
    geometry,
    active_electrons=active_electrons, # figure out how to get this
    active_orbitals=active_orbitals, # figure out how to get this
    package='psi4'
)
​
singles, doubles = qchem.excitations(active_electrons, qubits)
Psi4 calculation for H2-O1_sto-3g_singlet has failed.
---------------------------------------------------------------------------
FileNotFoundError                         Traceback (most recent call last)
/opt/anaconda3/lib/python3.8/site-packages/openfermionpsi4/_run_psi4.py in run_psi4(molecule, run_scf, run_mp2, run_cisd, run_ccsd, run_fci, verbose, tolerate_error, delete_input, delete_output, memory, template_file)
    208     try:
--> 209         process = subprocess.Popen(['psi4', input_file, output_file])
    210         process.wait()

/opt/anaconda3/lib/python3.8/subprocess.py in __init__(self, args, bufsize, executable, stdin, stdout, stderr, preexec_fn, close_fds, shell, cwd, env, universal_newlines, startupinfo, creationflags, restore_signals, start_new_session, pass_fds, encoding, errors, text)
    857 
--> 858             self._execute_child(args, executable, preexec_fn, close_fds,
    859                                 pass_fds, cwd, env,

/opt/anaconda3/lib/python3.8/subprocess.py in _execute_child(self, args, executable, preexec_fn, close_fds, pass_fds, cwd, env, startupinfo, creationflags, shell, p2cread, p2cwrite, c2pread, c2pwrite, errread, errwrite, restore_signals, start_new_session)
   1705                         err_msg = os.strerror(errno_num)
-> 1706                     raise child_exception_type(errno_num, err_msg, err_filename)
   1707                 raise child_exception_type(err_msg)

FileNotFoundError: [Errno 2] No such file or directory: 'psi4'

During handling of the above exception, another exception occurred:

UnboundLocalError                         Traceback (most recent call last)
<ipython-input-8-39f4437d2525> in <module>
      2 active_orbitals = 3
      3 
----> 4 H, qubits = qchem.molecular_hamiltonian(
      5     symbols,
      6     geometry,

/opt/anaconda3/lib/python3.8/site-packages/pennylane_qchem/qchem/structure.py in molecular_hamiltonian(symbols, coordinates, name, charge, mult, basis, package, active_electrons, active_orbitals, mapping, outpath, wires)
    836     """
    837 
--> 838     hf_file = meanfield(symbols, coordinates, name, charge, mult, basis, package, outpath)
    839 
    840     molecule = openfermion.MolecularData(filename=hf_file)

/opt/anaconda3/lib/python3.8/site-packages/pennylane_qchem/qchem/structure.py in meanfield(symbols, coordinates, name, charge, mult, basis, package, outpath)
    329         from openfermionpsi4 import run_psi4
    330 
--> 331         run_psi4(molecule, run_scf=1, verbose=0, tolerate_error=1)
    332 
    333     if package == "pyscf":

/opt/anaconda3/lib/python3.8/site-packages/openfermionpsi4/_run_psi4.py in run_psi4(molecule, run_scf, run_mp2, run_cisd, run_ccsd, run_fci, verbose, tolerate_error, delete_input, delete_output, memory, template_file)
    211     except:
    212         print('Psi4 calculation for {} has failed.'.format(molecule.name))
--> 213         process.kill()
    214         clean_up(molecule, delete_input, delete_output)
    215         if not tolerate_error:

UnboundLocalError: local variable 'process' referenced before assignment

Hi @JoshLipman and thanks again for bringing these issues to our attention.

Could you please verify that psi4 is installed properly by using import psi4? You can also run a small example to make sure psi4 works properly:

import psi4
psi4.set_memory('500 MB')
h2o = psi4.geometry("""
O
H 1 0.96
H 1 0.96 2 104.5
""")
psi4.energy('scf/sto-3g')

I recommend creating a fresh conda environment (python3.9 worked fine for me) and install psi4 with:

conda config --add channels http://conda.anaconda.org/psi4
python -m conda install psi4

When psi4 is properly installed and can be imported, you should also do a small modification in your code: the atomic coordinates are created with pennylane numpy and have the pennylane.numpy.tensor type. When this object is imported to psi4, with the H(x) function in your code, you will probably get an error when the psi4 input file is created (something like name 'tensor' is not defined). You can fix this by modifying the H(x) function such that the coordinates are converted to numpy.ndarray

import numpy
def H(x):
    return qml.qchem.molecular_hamiltonian(symbols,
                                           numpy.array(x),
                                           active_electrons = 6,
                                           active_orbitals = 5,
                                           package='psi4')[0]

I was able to run the modified version of your code and optimise the geometry of water. If you still face any issues please do not hesitate to let us know. Thanks!

Hi @sjahangiri,

I’m having trouble installing psi4. It’s not getting past the solving environment step. I tried pip also but no luck. Perhaps I need to build from source.

Thanks.

Hi @JoshLipman,
Just summarising what @sjahangiri nicely explained, assuming you are using conda.
You can create a new conda environment with:
conda create --name name_of_your_environment python=3.9
Activate the environment:
conda activate name_of_your_environment
and after this you can install the needed packages:
pip install numpy pennylane pennylane-qchem
Once you have performed these steps install psi4:
conda config --add channels http://conda.anaconda.org/psi4
conda install psi4
Installing psi4 can take some time, in some cases.

Once you have followed these steps you should be OK to run the example provided by @sjahangiri.
Let us know if this worked for you.

Hi @JoshLipman. We have seen this before, creating a new conda environment and installing psi4 with conda typically worked fine.

Hi @sjahangiri and @adusko, that worked. Thank you so much for the help on this! I used tensorflow to change the geometry tensor into an array and everything is working well.

Josh

Hey @JoshLipman!
Thank you for letting us know!