Dear Tangelo devs,
Hope all is good. This is me agian :)

I'm trying to create a circuit for Hamiltonian simulation and then export the resulting gateset to QASM or json.
The problem is that when the circuit is exported to QASM its gates are decomposed into 1-qubit and 2-qubit rotation gates. What I need is a QASM or json file which consists of Pauli rotation gates. For example, this is the description I'm looking for:

[
    {"gate": "PauliRot", "param": 0.045321883918106265, "pauli_word": "YXXY", "wires": [0, 1, 2, 3]},
    {"gate": "PauliRot", "param": -0.045321883918106265, "pauli_word": "YYXX", "wires": [0, 1, 2, 3]}
]

Is there a way of getting a QASM or jason like the one above from Tangelo?

Here is the code:

from tangelo import SecondQuantizedMolecule
from tangelo.toolboxes.qubit_mappings.mapping_transform import get_qubit_number, fermion_to_qubit_mapping
from tangelo.linq import Circuit, Gate
from tangelo.toolboxes.ansatz_generator.ansatz_utils import get_exponentiated_qubit_operator_circuit
from tangelo.linq.translator.translate_openqasm import translate_openqasm

# Molecule types: HH, C2HF3
molecule_name = "C2HF3"

if molecule_name == "C2HF3":
    molecule = [
        ("C", (0.0, 0.430631, 0.000000)),
        ("C", (-0.704353, -0.686847, 0.0)),
        ("F", (1.314791, 0.500587, 0.0)),
        ("F",(-0.554822, 1.630861, 0.0)),
        ("F", (-0.092517, -1.882156, 0.0)),
        ("H",(-1.780948, -0.706341, 0.0))   
    ]
elif molecule_name == "HH":
    molecule = [('H', (0, 0, 0)),('H', (0, 0, 0.74137727))]

second_quant_mol = SecondQuantizedMolecule(molecule, q=0, spin=0, basis="sto-3g")

hamiltonian = fermion_to_qubit_mapping(second_quant_mol.fermionic_hamiltonian, "JW")
circuit = get_exponentiated_qubit_operator_circuit(hamiltonian, time=0.5)

qasm_circuit_string = translate_openqasm(circuit)

# Write QASM file
qasm_file_name = f"./tangelo_qasm_circuits/circuit_molecule_{molecule_name}.qasm"
qasm_file_out = open(qasm_file_name, "w")
qasm_file_out.write(qasm_circuit_string)
qasm_file_out.close()

The code above works. It exports the QASM file but it breaks down the circuit into 1- and 2- qubit gates.

Issue: Bug ### Report

Expected Behavior
The stacked circuit feature should stack copies of a list of circuits width-wise. It seems that when the function tries to reindex_qubits a second circuit before stacking, it reindexes the first copy as well. The result is a circuit with two concatenated, and re-indexed copies. The reindex_qubits also fails when you put more than two circuits into stack due to wrong length, but that's okay for now.

Reproduce Issue

from tangelo.linq  import Gate, Circuit, stack
import copy

def circuit(theta=1.57):
    circuit=Circuit()
    circuit.add_gate(Gate("RX",0,parameter=3.14159165359))
    circuit.add_gate(Gate("RX",1,parameter=3.14159165359))    
    circuit.add_gate(Gate("CNOT",1, 0))
    circuit.add_gate(Gate("RZ",1,parameter=theta,is_variational=True))
    circuit.add_gate(Gate("CNOT",1, 0))
    circuit.add_gate(Gate("RX",0,parameter=10.995574287564))
    circuit.add_gate(Gate("RX",1,parameter=10.995574287564)) 
    return circuit

circuit = circuit()
#produces two concatenated, reindexed circuits as opposed to two stacked circuits
circ_stack=stack(circuit,circuit)
print(circ_stack)

# cannot put more than two circuits in
circ_stack3=stack(circuit,circuit,circuit)

Possible Solution
Make a deepcopy of the first circuit popped from list before reindexing

def stack_edit(*circuits):
    if not circuits:
        return Circuit()
    # Trim qubits of input circuit for maximum compactness
    circuits = [c.trim_qubits() for c in copy.deepcopy(list(circuits))]

    #stacked_circuit = circuits.pop(0)
    stacked_circuit=copy.deepcopy(circuits.pop(0))  #make a copy of initial circuit instead?!
    for c in circuits:
        c.reindex_qubits(list(range(stacked_circuit.width, stacked_circuit.width + c.width)))
        stacked_circuit += c
    return stacked_circuit

#produces stacked copies of circuits
circ_stack_edit=stack_edit(circuit,circuit)
print(circ_stack_edit)

Issue: Bug Report

Expected Behavior
When generating an ADAPT ansatz with scbk, there should be 2 qubits removed

Current Behavior
ADAPT returns a circuit with superfluous quits which are not actually used for anything. These can be removed using the trim_qubits method of the Circuit class.

Steps to Reproduce (minimal example)
The ansatz in the project Neon should be a 2 qubit circuit, but it clearly indicates a width of 4 qubits, even though there are only gates on the first 2.

Environment
MacOS 11.6.1

Issue: Bug Report

We wanted to provide a default draw method on the Circuit class in Tangelo for convenience, and decided that cirq 's would do: when we call draw() we implicitly ask cirq to draw it for us (Thanks cirq ❤️ )

Tangelo circuits can take almost anything as a value for gate.parameter, including a string. cirq is however not happy with that:

circuit = Circuit([Gate('H', 0), Gate('CNOT', 1, 0), Gate('RY', 1, parameter='theta')])
circuit.draw()

Can you help us ensure draw() works for strings as parameters ?

How to win

It turns out cirq is happy if we define a sympy symbol. Then we observe that both draw and print in Tangelo work just fine:

Which means one could ask draw to make a copy of the circuit on-the-fly (we don't want the method to modify the object itself) that replaces string values for parameters by equivalent sympy symbols.

• Try with an arbitrary circuit of 100 gates with parameters: if this takes less than 2-3 seconds, then your solution is performant enough.
• If not, compare the time of your new draw function and and show us the difference in performance with the old one (which you can keep around and rename: we'll only push the new one into the PR).

Validation

• Verify that draw and print give satisfying results on the circuit above
• give [email protected] access to a Google Colab notebook whose first cell installs Tangelo from your branch and then runs those draw / print cells with a few clicks and if it looks good, we're good !

Issue: Bug Report

I found a bug in the implementation of qubit tapering [2] where the optimal eigenvalues (aka eigenvalue sector) is determined incorrectly. I assume that the implementation follows the procedure as outlined in [3]. This means that after tapering, the smalles eigenvalue that can be obtained is higher than the smallest eigenvalue of the original hamiltonian! This is regarding the implementation of tangelo.toolboxes.operators.z2_tapering.get_eigenvalues.

This bug became apparent when using the HeH+ example as in the pennylane tutorial/demo on qubit tapering [1].
Most interestingly and importantly, this exact issue is present in the pennylane implementation (as one can see in the tutorial that the smallest eigenvalue becomes higher through tapering). Therefore, I have high doubts about the correctness of this entire approach of determining the optimal sector according to Setia et al., 2020 [3].
I am interested in investigating this concern in more details, and feel free to reach out to me to join.

Expected Behavior
When qubit tapering tapers qubits on which the hamiltonian has a non-trivial impact (i.e., Pauli-X instead of Identity), the eigenvalues (+1 or -1) of these qubits still have to be determined as they induce the tapered hamiltonian. An efficient implementation for this subroutine (to the best of my knowledge following the approach outlined in [3]) is part of Tangelo. It is expected that these eigenvalues are determined in a way that the tapered hamiltonian still has the same smallest eigenvalue as the original hamiltonian.
In the HeH+ example, this means that the two tapered qubits should either have eigenvalues -1, +1 or +1, -1, determined by tangelo.toolboxes.operators.z2_tapering.get_eigenvalues. The expected lowest eigenvalue of the (tapered) hamiltonian would be -3.2666.

Current Behavior
The current implementation of tangelo.toolboxes.operators.z2_tapering.get_eigenvalues does not seem to guarantee the computation of the correct/optimal eigenvalues. In the HeH+ example, the eigenvalues of -1, -1 are returned, which gives a lowest eigenvalue of the tapered hamiltonian of -2.8150.

Steps to Reproduce (minimal example)
The minimal example including the scenario as in the Pennylane tutorial [1] can be reproduced with the following code snippet. Note that it relies on the Psi4 backend as I could not make the SecondQuantizedMolecule initialization work using PySCF.

import numpy as np
from tangelo import SecondQuantizedMolecule
from tangelo.toolboxes.molecular_computation import IntegralSolverPsi4
from tangelo.toolboxes.qubit_mappings.mapping_transform import fermion_to_qubit_mapping
from tangelo.toolboxes.operators.operators import qubitop_to_qubitham
from tangelo.toolboxes.operators.operators import count_qubits
from openfermion.linalg import get_sparse_operator
from tangelo.toolboxes.operators.taper_qubits import QubitTapering
from tangelo.toolboxes.operators import MultiformOperator

symbols = ["He", "H"]
geometry = np.array([[0.00000000, 0.00000000, -0.87818361],
                     [0.00000000, 0.00000000,  0.87818362]])
xyz = list(zip(symbols, map(lambda x: tuple(x.tolist()), geometry)))
print(xyz)

mol = SecondQuantizedMolecule(xyz, q=+1, spin=0, solver=IntegralSolverPsi4(), frozen_orbitals=0)
qubit_h = qubitop_to_qubitham(fermion_to_qubit_mapping(mol.fermionic_hamiltonian, "jw", up_then_down=False), mapping="jw", up_then_down=False)
qubit_h_dense = get_sparse_operator(qubit_h).todense()
print("Qubits", count_qubits(qubit_h), "Terms", qubit_h.n_terms, "Constant", qubit_h.constant)
print("Eigenvalues:", np.linalg.eigvalsh(qubit_h_dense).tolist())

# ###### 🐛 BUG: Chooses incorrect sector/eigenvalues 🐛 ######
# 
# *** This leads to a higher smallest eigenvalue than for the untapered Hamiltonain ***

print("\n\nIncorrect current implementation to select correct eigenvalues for sector:")

taper = QubitTapering(qubit_h, n_qubits=4, n_electrons=2, spin=0, up_then_down=False)

taper_h = taper.z2_tapered_op
print("Tapered Hamiltonian:", taper.z2_tapered_op)
print(count_qubits(taper_h), taper_h.n_terms, taper_h.constant)
print("Sector:", taper.z2_properties["eigenvalues"])
taper_h_dense = get_sparse_operator(taper_h).todense()
print("Eigenvalues (tapered):", np.linalg.eigvalsh(taper_h_dense).tolist())

# ###### 🐛 Verify with correct optimal sector/eigenvalues 🐛 ######
# 
# *** (Manually) Choosing the eigenvalues to be -1, 1 (also applies to 1, -1) yields the matching correct smallest eigenvalue ***

print("\n\nVerify with manually chosen correct eigenvalues for sector:")

hybrid_op = MultiformOperator.from_qubitop(qubit_h, count_qubits(qubit_h))
sector = [-1, 1]
print("Sector:", sector)
taper_h = taper.z2_taper(hybrid_op, eigenvalues=sector)
taper_h_dense = get_sparse_operator(taper_h).todense()
print("Eigenvalues (tapered):", np.linalg.eigvalsh(taper_h_dense).tolist())

Environment
Please provide at least information about your OS, as well as the branch and version of Tangelo you are using. List other packages used and their version if relevant.

Possible Solution

Summary
The computation of the optimal eigenvalues/sector in qubit tapering leads to sub-optimal lowest eigenvalues, i.e., the tapered hamiltonian might has a higher lowest eigenvalue as the lowest one in the un-tapered hamiltonian. This is not the expected behavior.
However, since the same problem occurs in the Pennylane implementation, I highly suspect that there might be a deeper issue with the approach of determining the optimal eigenvalues/sector for the non-trivially tapered qubits as outlined in Setia et al. 2020 [3].

References

1. https://pennylane.ai/qml/demos/tutorial_qubit_tapering/
2. Bravyi, Sergey, Jay M. Gambetta, Antonio Mezzacapo, and Kristan Temme. “Tapering off Qubits to Simulate Fermionic Hamiltonians.” arXiv, January 27, 2017. https://doi.org/10.48550/arXiv.1701.08213.
3. Setia, Kanav, Richard Chen, Julia E. Rice, Antonio Mezzacapo, Marco Pistoia, and James D. Whitfield. “Reducing Qubit Requirements for Quantum Simulations Using Molecular Point Group Symmetries.” Journal of Chemical Theory and Computation 16, no. 10 (October 13, 2020): 6091–97. https://doi.org/10.1021/acs.jctc.0c00113.

Circuit depth code

The tangelo.linq.circuit.depth method:

One way to circumvent this is to use the cirq package (it uses moment during the construction, and it is somehow faster). For e.g.,

from tangelo.molecule_library import mol_H2O_sto3g, mol_H2O_321g
from tangelo.toolboxes.ansatz_generator import UCCSD

mol = mol_H2O_321g

ansatz = UCCSD(mol, mapping="JW")
ansatz.build_circuit()
circuit = ansatz.circuit

circuit.depth()

does not output anything in 20 min on my machine. However,

from tangelo.linq import translate_circuit

def depth_from_cirq(circuit):
    cirq_circuit = translate_circuit(circuit, source="tangelo", target="cirq")
    return len(cirq_circuit.moments)

depth_from_cirq(circuit)

outputs

CPU times: user 8.39 s, sys: 18.7 ms, total: 8.4 s
Wall time: 8.42 s
116869

Moreover, printing the depth in resource estimation should be optional.

I have identified some bugs in the current "main" state of the qSDK.

VQE

1. When using {"qubit_mapping": "bk", "verbose": True} or {"qubit_mapping": "scbk", "verbose": True}, the object returned does not have a count_qubits method. Simulation runs fine if verbose is set to False.

  File "QEMIST_qSDK/qsdk/electronic_structure_solvers/vqe_solver.py", line 128, in build
    n_qubits = self.qubit_hamiltonian.count_qubits()
AttributeError: 'QubitOperator' object has no attribute 'count_qubits'

2. When using {"initial_var_params": "random"} with UCCSD ansatsz, there is an unexpected dtype.

  File "QEMIST_qSDK/qsdk/toolboxes/ansatz_generator/uccsd.py", line 64, in set_var_params
    initial_var_params = np.random.random((self.n_var_params,), dtype=float)
  File "mtrand.pyx", line 426, in numpy.random.mtrand.RandomState.random
TypeError: random() got an unexpected keyword argument 'dtype'

Other

3. Even when freezing orbitals, bigger molecule cannot be simulated with the current code state. It is always doing an FCI computation. I would turn off those flags in qsdk/toolboxes/molecular_computation/molecular_data.py.

run_pyscf(self, run_scf=True, run_mp2=True, run_cisd=True, run_ccsd=True, run_fci=True)

Those are problems that can be fixed easily. I created an issue for documentation purposes.

Application of initial parameters seem problematic. When debugging DMET, MP2 initial parameters returns an error (maybe it is not correct to do them for fragments). When I tried to change to "ones", it was still trying to get MP2 parameters.

Here is a test file that I used:

from pyscf import gto
from qsdk.electronic_structure_solvers.vqe_solver import VQESolver

H2 = [
    ('H', (0., 0., 0.)),
    ('H', (0., 0., 0.75)),
    ]
mol = gto.Mole()
mol.atom = H2
mol.basis = "sto-3g"
mol.spin = 0
mol.build()

solver = VQESolver({"molecule": mol, "initial_var_params": [1., 2.]})

solver.build()

print(solver.initial_var_params)
print(solver.ansatz.var_params)

The last two lines confirm to me that the MP2 parameters (default) have been applied instead of [1., 2.] as stated. This is related to those two lines
https://github.com/1QB-Information-Technologies/QEMIST_qSDK/blob/4357c9fce35c260f5e07e8a2cb63a8cb0a59a58d/qsdk/electronic_structure_solvers/vqe_solver.py#L119-L120

I propose to change that to
self.ansatz.build_circuit(var_params=self.initial_var_params)

I still have problem with the initial parameters after that with DMET but it is related to the DMET code, not our VQESolver.

Incompatibility with PySCF >= 2.5.0

This error currently occurs in a single test: ansatz_generator.tests.test_uccsd.UCCSDTest.test_uccsd_H4_open, where unrestricted MP2 amplitudes are computed as initial parameters for a UCCSD ansatz in VQE.

With PySCF<=2.4.0, the test passes. However, with the latest versions of PySCF (2.5.0 and on) the test currently returns an error: it seems there was a change in a the eri data-structures. It may be that a simple fix could be done inside our MP2Solver class.

Expected Behavior
Find a way to make this test pass for both PySCF 2.5.0 while still supporting the older versions. Test locally in both environment (you can easily switch between versions with pip install and uninstall, by specifying the desired version).
The continuous integration environment yml file should then be bumped to pyscf latest version so that all automated tests run with it on.

Current Behavior
An IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed is raised when computing MP2 initial parameters for open-shell systems.

Steps to Reproduce (minimal example)

Minimal script:

from tangelo.molecule_library import mol_H4_cation_sto3g
from tangelo.toolboxes.ansatz_generator.uccsd import UCCSD

uccsd_ansatz = UCCSD(mol_H4_cation_sto3g)
uccsd_ansatz.build_circuit()

outputs

Traceback (most recent call last):
  File "/Users/alexandre.fleury/Work/Scratch/sandbox.py", line 6, in <module>
    uccsd_ansatz.build_circuit()
  File "/Users/alexandre.fleury/Repos/Tangelo/tangelo/toolboxes/ansatz_generator/uccsd.py", line 164, in build_circuit
    self.set_var_params()
  File "/Users/alexandre.fleury/Repos/Tangelo/tangelo/toolboxes/ansatz_generator/uccsd.py", line 127, in set_var_params
    initial_var_params = self._compute_mp2_params()
                         ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/alexandre.fleury/Repos/Tangelo/tangelo/toolboxes/ansatz_generator/uccsd.py", line 283, in _compute_mp2_params
    mp2.simulate()
  File "/Users/alexandre.fleury/Repos/Tangelo/tangelo/algorithms/classical/mp2_solver.py", line 316, in simulate
    return self.solver.simulate()
           ^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/alexandre.fleury/Repos/Tangelo/tangelo/algorithms/classical/mp2_solver.py", line 89, in simulate
    _, self.mp2_t2 = self.mp2_fragment.kernel()
                     ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/alexandre.fleury/.pyenv/versions/qsdk/lib/python3.11/site-packages/pyscf/mp/mp2.py", line 598, in kernel
    eris = self.ao2mo(mo_coeff)
           ^^^^^^^^^^^^^^^^^^^^
  File "/Users/alexandre.fleury/.pyenv/versions/qsdk/lib/python3.11/site-packages/pyscf/mp/ump2.py", line 437, in ao2mo
    return _make_eris(self, mo_coeff, verbose=self.verbose)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/alexandre.fleury/.pyenv/versions/qsdk/lib/python3.11/site-packages/pyscf/mp/ump2.py", line 499, in _make_eris
    eris._common_init_(mp, mo_coeff)
  File "/Users/alexandre.fleury/.pyenv/versions/qsdk/lib/python3.11/site-packages/pyscf/mp/ump2.py", line 475, in _common_init_
    mo_a = mo_coeff[0][:,mo_idx[0]]
IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed

The IndexError points out the fact that the mo_coeff data structures have changed from PySCF 2.4.0 to 2.5.0. The solution as of now is to restrict the PySCF support to 2.4.0 and older versions, but it would be great if this bug could be fixed for future releases.

Issue: Feature Request

What problem does this feature request help you overcome? Please describe.
In the DMET problem decomposition, after the bath orbitals are created and the dummy molecule is constructed and passed to the solver, I would like to have an option of frozen core / selection of active spaces, as this will reduce the computational complexity of the solver, especially VQE when you run it on small hardware devices.

Describe the solution you'd like
In the DMET options, provide a new option that set the frozen orbitals for each fragment.

DMETProblemDecomposition refactoring could be a long-term project, as it is an important module in our repo. Here is a task list I am creating and will be updating in the future for better transparency / if someone wants to tackle this endeavour.

• Improve overal code efficiency and data flow.
• Decouple PySCF from DMETProblemDecomposition.
• Implement MPI parallelization.
• Refactor "solvers" for better modularity.
• Remove RDMs purification and post-processing from the DMET code. This is a task related to the modularity of solvers.

Issue: Feature Request

What problem does this feature request help you overcome? Please describe.
Our current implementation to get circuits from the qasm format only supported a subset of OpenQASM, and was written early 2020.

Describe the solution you'd like
Implement this using the reverse qiskit translation combined with https://qiskit.org/documentation/stubs/qiskit.circuit.QuantumCircuit.from_qasm_str.html#qiskit.circuit.QuantumCircuit.from_qasm_str ?

The file is located here https://github.com/goodchemistryco/Tangelo/blob/develop/tangelo/linq/translator/translate_openqasm.py.

Issue: Feature Request

We recently went to a hackathon, and participants wanted to solve a challenge using a circuit with multi-controlled Ry gates. They wanted to translate the Tangelo circuit to Qiskit to simulate a fake IBM device or transpile to IBM native gate set, but got a NotImplemented error from Tangelo.

Indeed, we had written that function before Qiskit implemented multi-controlled gates, but Qiskit later supported them: here's the mcry gate https://docs.quantum.ibm.com/api/qiskit/0.43/qiskit.circuit.QuantumCircuit#mcry.

Can you modify https://github.com/ValentinS4t1qbit/Tangelo/blob/main/tangelo/linq/translator/translate_qiskit.py so that multi-controlled gates (CRX, CRY, CRZ) are supported in the translation from Tangelo to Qiskit and remove the error message ? (e.g, a list on integers can be passed for control qubits).

It may be that the crx, cry and crz gates in Qiskit can simply be seen as mcrx, mcry and mcrz with only one control qubit, which may help us resolve things elegantly with very little changes.

Test

Making this snippet work for CRY is a good starting point:

from tangelo.linq import Gate, Circuit, translate_circuit, get_backend

# Defines a multi-controlled RY gate in tangelo
c = Circuit([Gate('CRY', 2, control=[0,1], parameter=1.)])
c_qiskit = translate_circuit(c, target="qiskit")  # should not return an error

# Qiskit seems to parse this as a sequence of gate
print(c_qiskit)

You should see that Qiskit parses that gate into this:

To solve this issue:

• Create a test that looks at 3 separate circuits of one multi-controlled gate each (one for CRX, CRY, CRZ, like the one above. We could use a loop to do this elegantly)
• Define the circuit in Tangelo, translate it to qiskit.
• Write the equivalent native qiskit circuit using their mcr* gate.
• Test that the translated qiskit object is "equal" to the native one.

Issue: Bug Report

Expected Behavior
RZ gates when constructing the QMF circuit should be appended to the RX gates on the same qubits.

Current Behavior
RZ gates are appended to extra qubits.
https://github.com/quantumsimulation/QEMIST_qSDK/blob/f6187c853052a7a2f33dbeccf6b20a60bd72718e/tangelo/toolboxes/ansatz_generator/_qubit_mf.py#L159-L161

Possible Solution
target=idx%n_qubits

Here are some bugs that I found:

• self.qubit_hamiltonian.n_qubits should be removed.
  https://github.com/quantumsimulation/QEMIST_qSDK/blob/b6bac55a942c56c11af1ecb0dc6b96f7f6e7bda5/qsdk/qsdk/electronic_structure_solvers/vqe_solver.py#L125
• Adding mol.verbose = 0
  https://github.com/quantumsimulation/QEMIST_qSDK/blob/b6bac55a942c56c11af1ecb0dc6b96f7f6e7bda5/qsdk/qsdk/toolboxes/molecular_computation/molecule.py#L95-L112
• Adding copy_self.frozen_orbitals = frozen_orbitals
  https://github.com/quantumsimulation/QEMIST_qSDK/blob/b6bac55a942c56c11af1ecb0dc6b96f7f6e7bda5/qsdk/qsdk/toolboxes/molecular_computation/molecule.py#L319-L343
• In DMET code
  if purify and solver_fragment.qemist_molecule.n_electrons == 2:
  AttributeError: 'VQESolver' object has no attribute 'qemist_molecule'

Issue: Bug Report

Current Behavior
Computation of RDMS with CCSD raises an error:

Exception ignored in: <function H5TmpFile.__del__ at 0x7fd655bf93a0>
Traceback (most recent call last):
  File "/home/alex/VirtEnvs/qsdk/lib/python3.8/site-packages/pyscf/lib/misc.py", line 1004, in __del__
    self.close()
  File "/home/alex/VirtEnvs/qsdk/lib/python3.8/site-packages/h5py/_hl/files.py", line 552, in close
    self.id._close_open_objects(h5f.OBJ_LOCAL | h5f.OBJ_FILE)
  File "h5py/_objects.pyx", line 54, in h5py._objects.with_phil.wrapper
  File "h5py/_objects.pyx", line 55, in h5py._objects.with_phil.wrapper
  File "h5py/h5f.pyx", line 360, in h5py.h5f.FileID._close_open_objects
RuntimeError: Can't decrement id ref count (unable to extend file properly, errno = 2, error message = 'No such file or directory')
Segmentation fault (core dumped)

It may be related to pyscf config and h5py temporary file close method.

Steps to Reproduce (minimal example)

1. Create a new python environment.
2. Install tangelo (pip install tangelo-gc).
3. Compute the RDMs (with CCSDSolver) of any molecule (call the f = lib.H5TmpFile() line).

Environment
Virtual machine, Ubuntu 22.04.1 LTS (GNU/Linux 5.15.0-56-generic x86_64)

Summary of the h5py configuration
---------------------------------

h5py    3.7.0
HDF5    1.12.2
Python  3.8.15 (default, Dec  6 2022, 10:28:26)
[GCC 11.3.0]
sys.platform    linux
sys.maxsize     9223372036854775807
numpy   1.23.5
cython (built with) 0.29.30
numpy (built against) 1.17.3
HDF5 (built against) 1.12.2

I'm also getting the same error with pythons>=3.7, h5py>=2.8.

Possible Solution
(Temporary workaround) Setting the TMPDIR or the PYSCF_TMPDIR environment variables to /tmp help preventing the error. By default, PySCF set the scratch path as . when TMPDIR is not set. However, the H5 libraries does not seem to play well with relative paths.

Issue: Bug Report

Expected Behavior
The do_jkmn_transform(vector) function is supposed to return a numpy array of int.

Current Behavior
Instead, the do_jkmn_transform(vector) function returns a list of int. See the prep_jkmn_vector(vector) function in jkmn.py -- the output is a list.

QCC jobs fail if the JKMN encoding is used. The failure arises when initializing the qubit mean field variational parameters via the init_qmf_from_hf function: this function requires a numpy array of int from the statevector_mapping module.

Steps to Reproduce (minimal example)
Rename the attached test_qcc_jkmn_py.txt file to test_qcc_jkmn.py.
Then run python test_qcc_jkmn.py.
The resulting error message is attached as test_qcc_jkmn_error.txt

Environment
CentOS Linux release 7.8.2003 (Core)
Tangelo-gc v0.3.2 is installed.

Possible Solution
Change the output of the do_jkmn_transform(vector) function in statevector_mapping.py to np.array(prep_jkmn_vector(vector)). Alternatively, the output of prep_jkmn_vector(vector) can be changed to x_vector instead of list(x_vector).

test_qcc_jkmn_error.txt
test_qcc_jkmn_py.txt

Method angostic_simulator.simulator.Simulator._get_expectation_value_from_frequencies output weird results when playing with NoiseModel. Here is a minimal example I wrote:

from agnostic_simulator import Circuit, Gate
from agnostic_simulator import Simulator
from agnostic_simulator.noisy_simulation import NoiseModel

from qsdk.toolboxes.operators import QubitOperator

# Be Hamiltonian (scBK).
H = QubitOperator()
H.terms = {(): -14.41806525945003, ((0, 'Z'),): 0.0809953994342687, ((1, 'Z'),): 0.0809953994342687, ((0, 'Z'), (1, 'Z')): 0.0077184273651725865, ((0, 'X'), (1, 'X')): 0.0758664717894615}

# Hard coding of IQCC-ILC optimized circuit.
circuit=Circuit()
circuit.add_gate(Gate("RX",0,parameter=3.141595416808))
circuit.add_gate(Gate("RX",1,parameter=3.141588753134))

circuit.add_gate(Gate("H",0))
circuit.add_gate(Gate("RX",1,parameter=1.570796326795))
circuit.add_gate(Gate("CNOT",1, 0))
circuit.add_gate(Gate("RZ",1,parameter=0.43912793,is_variational=True))
circuit.add_gate(Gate("CNOT",1, 0))
circuit.add_gate(Gate("RX",1,parameter=10.995574287564))
circuit.add_gate(Gate("H",0))

# Noise model.
nmp = NoiseModel()
noise=0.01
nmp.add_quantum_error("CNOT", "pauli", [noise, noise, noise])
sim_noise = Simulator("qulacs", n_shots=10**5, noise_model=nmp)

# Noisy results
E_depol=sim_noise.get_expectation_value(H, circuit)

print("Ideal: -14.58922316")
print(f"Noisy original energy: {E_depol}")

In short:

• Changing the variable noise has little to no impact on the energy.
• Same results with qulacs and qiskit backend (cirq returns an error, ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all())
• HINT: Simulating with Q# instead gives very different frequencies (line 445 in Simulator.py).
  https://github.com/quantumsimulation/QEMIST_qSDK/blob/d0738278379f89fa3f4be7c2da4421ebba526408/agnostic_simulator/agnostic_simulator/simulator.py#L446

BAD EXPECTATION VALUE
qulacs with minimal noise (0.01) 10^5 shots
term -> frequencies
((0, 'Z'),) -> {'00': 1.0}
((1, 'Z'),) -> {'00': 1.0}
((0, 'Z'), (1, 'Z')) -> {'00': 1.0}
((0, 'X'), (1, 'X')) -> {'10': 0.24891, '11': 0.25068, '00': 0.2518, '01': 0.24861}

GOOD EXPECTATION VALUE
qdk noiseless 10^5 shots
term -> frequencies
((0, 'Z'),) -> {'00': 0.04832, '11': 0.95168}
((1, 'Z'),) -> {'00': 0.04731, '11': 0.95269}
((0, 'Z'), (1, 'Z')) -> {'00': 0.04792, '11': 0.95208}
((0, 'X'), (1, 'X')) -> {'00': 0.14334, '10': 0.35398, '01': 0.35948, '11': 0.1432}

Issue: Feature Request

What problem does this feature request help you overcome? Please describe.
Thousands of lines of Job info ID:: JOB_ID status :: ready are printed after submitting a job to a QPU on IonQ Cloud while it waits to be executed. The information is redundant and adds an extra step to processing and analyzing the data.

Describe the solution you'd like
The job_status = ready really only needs to be printed one time. This can be accomplished by modifying the job_results and job_status functions in tangelo/linq/qpu_connection/ionq_connection.py.

Describe alternatives you've considered
Instead of only printing once, the duration of time.sleep() in job_results could be increased to several minutes, hours, days, etc.

Issue: Bug Report

Expected Behavior
Running the attached Python script via "python minimum_example.py" is expected to

1. Create a SecondQuantizedMolecule for ground state H2 / STO-3G basis set.
2. Perform classical VQE to obtain the optimal QCC circuit for this system.
3. Connect with IonQ Cloud's API and submit the job to the simulator backend.
4. Return a histogram from the simulation.

Current Behavior
The job fails before connecting to IonQ Cloud and never submits the simulation to the simulator backend. Several errors related to the remote end disconnecting without response, a protocol error, and a connection error are reported. The full error message is attached to this report along with an example.

Steps to Reproduce (minimal example)
The job fails at Step 3 below.

1. Set the environment variable IONQ_APIKEY to its current value.
2. Activate the appropriate environment with the installation of tangelo-gc v0.3.2.
3. Run the script with "python minimum_example.py"

Environment
CentOS Linux release 7.8.2003 (Core)
The output of conda list is attached to this report. Briefly, tangelo-gc v0.3.2 and qemist-client v0.5.0 are installed.

Possible Solution
Commenting the line from .qemist_cloud_connection import job_submit, job_status, job_cancel, job_result, job_estimate in tangelo/linq/qpu_connection/__init__.py fixed the issue.

tangelo-gc_conda_env.txt
minimum_example_error.txt
minimum_example_py.txt

Expected Behavior

When running CCSD and DMET-CCSD on a Zn²⁺ atom using the LANL2DZ basis set and ECP, the energies obtained should be equal as there is only one fragment / one atom.

Current Behavior

The CCSD and DMET-CCSD energies are not the same. There is a significant difference (21.5307 Eh) between them, suggesting a possible bug when using ECPs with DMET.

Steps to Reproduce

import tangelo

from tangelo import SecondQuantizedMolecule
from tangelo.problem_decomposition import DMETProblemDecomposition
from tangelo.algorithms import CCSDSolver

zn = """
Zn 0.0 0.0 0.0
"""

mol_zn = SecondQuantizedMolecule(zn, q=2, spin=0, basis="lanl2dz", ecp="lanl2dz")

options_zn_dmet = {
    "molecule": mol_zn,
    "fragment_atoms": [1],
    "fragment_solvers": "ccsd",
    "verbose": True,
}

dmet_zn = DMETProblemDecomposition(options_zn_dmet)
dmet_zn.build()

energy_zn_dmet = dmet_zn.simulate()
energy_zn_hf = dmet_zn.mean_field.e_tot
print(f"SCF energy (hartree): \t {energy_zn_hf}")
print(f"DMET energy (hartree): \t {energy_zn_dmet}")

ccsd_zn = CCSDSolver(mol_zn)
energy_zn_ccsd = ccsd_zn.simulate()
print(f"CCSD energy (hartree): \t {energy_zn_ccsd}")

delta_zn_ccsd_hf = abs(energy_zn_ccsd - energy_zn_hf)
delta_zn_ccsd_dmet = abs(energy_zn_ccsd - energy_zn_dmet)
print(f"Difference CCSD vs DMET-CCSD energies (hartree): \t\t {delta_zn_ccsd_dmet}")

gives

 	Iteration =  10
 	----------------

		Fragment Number : #  1
		------------------------
		Fragment Energy = -41.24106640579692
		Number of Electrons in Fragment = 10.0
 	Iteration =  11
 	----------------

		Fragment Number : #  1
		------------------------
		Fragment Energy = -41.236798699871244
		Number of Electrons in Fragment = 10.0
 	*** DMET Cycle Done ***
 	DMET Energy ( a.u. ) =    -41.2367986999
 	Chemical Potential   =      0.0002653333
SCF energy (hartree): 	 -62.63753650865124
DMET energy (hartree): 	 -41.236798699871244
CCSD energy (hartree): 	 -62.77176583635452
Difference CCSD vs DMET-CCSD energies (hartree): 		 21.534967136483274

Environment

OS: Ventura 13.4 with Apple M1 Max
Python Version: Python 3.10.10 | packaged by conda-forge | (main, Mar 24 2023, 20:12:31) [Clang 14.0.6 ] on darwin
Package Version: Tangelo 0.3.4

Possible Solution

My guess is the fragment energy evaluation does not incorporate the core energy component correctly.

Summary

There seems to be a bug when using ECPs in the DMET calculation, leading to a difference between the CCSD and DMET-CCSD energies. The expected behavior is that these two energies are the same in cases where the fragment equals the full space.

Adding support for accessing and modifying variational gates in ansatzes

What problem does this feature request help you overcome? Please describe.

In PR #392 support was added for the Rotoselect optimizer, which modifies both the angle and axis of rotation for variational gates in a given ansatz. Currently, the Rotoselect algorithm can only be used by passing in a function of the form

def exp_rotoselect(var_params, var_rot_axes, ansatz, qubit_hamiltonian):
            ansatz.update_var_params(var_params)
            for i, axis in enumerate(var_rot_axes):
                ansatz.circuit._variational_gates[i].name = axis
            energy = sim.get_expectation_value(qubit_hamiltonian, ansatz.circuit)
            return energy

This function requires the user to manually modify the protected member ansatz.circuit._variational_gates, which is not ideal.

Describe the solution you'd like
It would be great if we could add a better interface for modifying variational gates in ansatzes that does not require directly modifying the ansatz circuit.

I can think of two ways of implementing this. The most invasive approach would involve extending the interface of the Ansatz abstract class (in tangelo.toolboxes.ansatz_generator.ansatz) to include the functions set_var_gates() and update_var_gates(). Then, each ansatz object which extends the abstract Ansatz class will need to be modified to accommodate it in a manner that "makes sense" in the context of each ansatz.

Describe alternatives you've considered
Another (less invasive) approach could be to create some methods that directly modify the variational gates in a Circuit, while performing some basic validation. A method to retrieve all of the variational gates in a circuit would also be necessary. This would be much easier to implement, but would place the burden of updating gates in a meaningful way on the user.

Additional context
I am currently working on a solution along the lines of the first approach, but I am encountering some difficulties with handling certain ansatzes (such as UCCGD) where there appears to be a one-to-many mapping of variational parameters to variational gates. If we want to modify the axis of rotation of gates in these ansatzes, I suppose it would make sense to respect this one-to-many mapping so that certain symmetries of the ansatz wave function are preserved. I am also concerned that adding this feature to some ansatzes may be "non-physical" in the sense that certain canonical commutation relations will be violated. For example, some ansatzes contain parameterized exponentials of fermionic creation/annihilation operators mapped to pauli strings with the same CCRs. If we were to allow for arbitrary axes of rotation that violate these CCRs, would this be problematic? Should we add this feature for only some ansatzes and not others?

I would like to see if anyone has any ideas or recommendations for whether or not this feature should be added and how it can be done in a manner that is useful and physically meaningful. ⚛️

When qemist_client is not installed in the python environment, the error raised when trying to submit a problem through qemist is not helpful:
Name Error: name util is not defined

This is because the ModuleNotFoundError is caught at the beginning of qemist_cloud_connection.py file. I think it was made like this to prevent unnecessary errors when testing. I suggest to modify the code so if the submission function is called without qemist_client, a ModuleNotFoundError is raised instead.

Issue: Feature Request

The Rotosolve optimizer is a straightforward method to obtain the optimal value for a single qubit rotation in variational quantum algorithms. The current implementation in Tangelo supports its use in ongoing research involving the Qubit Coupled Cluster algorithm where typically only a single rotation axis is considered. This feature can be easily augmented to the Rotoselect optimizer, which chooses both the optimal angle and axis of single qubit rotation.

• Paper describing Rotoselect: https://arxiv.org/pdf/1905.09692
• Implementation of Rotosolve in Tangelo: https://github.com/goodchemistryco/Tangelo/blob/main/tangelo/toolboxes/optimizers/rotosolve.py

Describe the solution you'd like

The solution would include a separate solver called Rotoselect which modifies the rotation axis (RX, RY, RZ) of the gates in the variational circuit. For each rotation axis, it could call rotosolve_step to find the optimal angle.

Finish line

Testing for this feature should be analogous to the rotosolve tests. The test would require an input variational circuit, and a Hamiltonian to minimize. These should be constructed such that the output optimized circuit has a different rotation axis for at least one rotation gate.

Automatic documentation deployment

Currently, updating the Tangelo documentation involves significant manual work. We would like to have an automatic deployment workflow, similar to what we have set up in the Tangelo-Examples repo.

What problem does this feature request help you overcome?
As of now, regenerating or updating the documentation requires performing the steps outlined in dev_tools/build_sphinx_docs.sh locally, and then uploading the files to an AWS bucket. It's cumbersome and manual, and sometimes it may not be updated right away.

Describe the solution you'd like
We would like an automatic workflow (similar to the one in the Tangelo-Examples repo, see files below) to compile and upload the documentation using GitHub Pages every time the main branch is updated, to guarantee that our documentation is always reliable !

• The documentation left behind by past versions of ourselves should be helpful here: https://github.com/goodchemistryco/Tangelo-Examples/blob/main/dev_tools/generate_documentation.rst
• An "automated workflow" file we should add to the repo is available here
• There are also other online resources that may be relevant.

The source files would reside in a gh-pages branch, and the website could be accessed via a GitHub Pages address https://goodchemistryco.github.io/Tangelo. After that, we would only need to point http://tangelo-docs.goodchemistry.com/ to the corresponding GitHub page. (this address will be changed to a SandboxAQ address later, don't worry about it 😄 )

Circuit as reference state in the ansatz definition

What problem does this feature request help you overcome?

In our stack of Ansätze, located in toolboxes/ansatz_generator, most of them take a string reference_state as an initialization parameter. This determines the reference state of the ansatze, which is commonly HF for an Hartree-Fock state or zero for an empty circuit. As an example, here is the code that handles the reference_state for the toolboxes.ansatz_generator.uccsd.UCCSD ansatz.

https://github.com/alexfleury/Tangelo/blob/9079cb2e8c4b893027a7ff70101e968548000714/tangelo/toolboxes/ansatz_generator/uccsd.py?plain=1#L136-L153

This is restrictive in situations where non-conventional mapping are used, i.e. where HF has no meaning.

Describe the solution you'd like

A workaround would be to add support for a tangelo.linq.circuit.Circuit objects as a reference state. Logics have to be implemented to handle the right number of qubits in some edge cases, and to consider optimization (or not) if the reference circuit contains variational gates.

For some ansatze, it would remove the need of specifying the mapping and up_then_down parameters, as they are only used for reference state preparation.

As a minimal use-case, exploring combinatorial mapping optimization could be easier if this option is available. For e.g., if the code below could work without error, it would be awesome!

from tangelo.molecule_library import mol_H2_sto3g
from tangelo.toolboxes.qubit_mappings import combinatorial
from tangelo.toolboxes.operators import count_qubits
from tangelo.linq import Circuit, Gate
from tangelo.toolboxes.ansatz_generator import HEA

H = combinatorial(mol_H2_sto3g.fermionic_hamiltonian, mol_H2_sto3g.n_active_mos, mol_H2_sto3g.n_active_electrons)
n_qubits = count_qubits(H)

ref_circ = Circuit([Gate("RY", 0, parameter=0.5), Gate("RY", 1, parameter=-0.5)])
ansatz = HEA(n_qubits=n_qubits, n_layers=2, reference_state=ref_circ)
ansatz.build_circuit()

It would be great to submit batch jobs to AWS Braket through QEMIST when running hardware experiments. Their API supports this. The feature would also need to keep track of the job ids for each circuit so we can call the results from QEMIST.

The input could be dictionary input like: {"job_title": {num_shots : circuit}...}
Which returns : {"job_title": job_id}

Issue: Bug Report

Expected Behavior
ADAPTSolver should update the optimal_* attributes after simulate has been called.

Current Behavior
When max_cycles>1 and the energy converges in less than max_cycles steps, the optimal attributes are not updated (var_params, circuit, etc.). They stay as None.

Steps to Reproduce (minimal example)

1. Perform ADAPT-VQE

from tangelo.algorithms import ADAPTSolver
from tangelo.molecule_library import mol_H2_sto3g

solver = ADAPTSolver({
    "molecule": mol_H2_sto3g,
    "max_cycles": 2,
    "verbose": True
})
solver.build()
solver.simulate()

circuit = solver.optimal_circuit
print(circuit)

2. Call optimal_circuit. In the previous example, circuit is None.

Environment
Ubuntu - python 3.10, latest commit of the develop branch

Small issue. An error should be raised when the number of target qubits does not match how many are involved in the gate. It looks like this check only exists for gates with a control qubit.

from tangelo.linq import Circuit, Gate

circ=Circuit()

swap_1=Gate('SWAP',[0])
swap_3=Gate('SWAP',[0,1,2])
circ.add_gate(swap_1)
circ.add_gate(swap_3)

XX_1=Gate('XX', [0])
XX_3=Gate('XX', [0,1,3])

circ.add_gate(XX_1)
circ.add_gate(XX_3)

RX_3=Gate('RX', [0,1,3])
circ.add_gate(RX_3)

print(circ)

returns:

SWAP      target : [0]   
SWAP      target : [0, 1, 2]   
XX        target : [0]   
XX        target : [0, 1, 3]   
RX        target : [0, 1, 3] 

Example code:

from numpy.random import randint
from qsdk.backendbuddy import Gate

# Randomly generates 2 random numbers in {0, 1}
i, j = randint(0, 2, 2)
g = Gate("CNOT",  i, j)

Currently triggers:

    if not (isinstance(target, int) and target >= 0):
        raise ValueError("Qubit index must be a positive integer.")

Suggestions:

• The error is triggered because the numpy int types are not captured by isinstance, here. To capture numpy int types, it should also include np.integer(https://stackoverflow.com/questions/37726830/how-to-determine-if-a-number-is-any-type-of-int-core-or-numpy-signed-or-not)
• Because the circuits are converted in different formats and type differences may lead to sneaky errors, I think for simplicity we should test for types (control, target, ...) and if everything is alright cast them as built-in ints in the init method.
• The error messages for target and control should be more useful and specific: Target qubit index must be a non-negative number (value: {target}) (likewise for control)
• Additional statement returning a warning if control == target. There's no reason to forbid it, but it's an odd thing to do. In this situation, the author of the code did not realize this could happen.

The following should be revised in UCCSD:

1. explicit declaration of molecule datatype requirements in documentation
2. error handling for initialization to zero-array (initial variational parameters)
3. alternative implementation of singlet_generator to track opt. params and gate params.

Feature Request: Deepcopy of SecondQuantizedMolecule

At this time, the SecondQuantizedMolecule can't be deepcopied:

from copy import deepcopy
from tangelo import SecondQuantizedMolecule

xyz_H2 = """
    H 0. 0. 0.00
    H 0. 0. 0.75
"""
mol = SecondQuantizedMolecule(xyz_H2)
deepcopy(mol)

outputs

Traceback (most recent call last):
  File "/home/alex/Scratch/mo_coeff_bug.py", line 14, in <module>
    deepcopy(mol)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 172, in deepcopy
    y = _reconstruct(x, memo, *rv)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 271, in _reconstruct
    state = deepcopy(state, memo)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 146, in deepcopy
    y = copier(x, memo)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 231, in _deepcopy_dict
    y[deepcopy(key, memo)] = deepcopy(value, memo)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 172, in deepcopy
    y = _reconstruct(x, memo, *rv)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 271, in _reconstruct
    state = deepcopy(state, memo)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 146, in deepcopy
    y = copier(x, memo)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 231, in _deepcopy_dict
    y[deepcopy(key, memo)] = deepcopy(value, memo)
  File "/home/alex/.pyenv/versions/3.10.13/lib/python3.10/copy.py", line 161, in deepcopy
    rv = reductor(4)
TypeError: cannot pickle 'module' object

This is not a critical bug at this moment, but it would be useful to be able to create a copy and update its attributes. For e.g.,

from copy import copy
from openfermion.linalg import eigenspectrum

def doing_something_with_molecule(sqmol):
    sqmol_copy = copy(sqmol)
    # Doing awesome stuff with mo_coeff and creating new_mo_coeff.
    sqmol_copy.mo_coeff = new_coeff
    return sqmol_copy

eigenvalues_original = eigenspectrum(mol.fermionic_hamiltonian)
updated_mol = doing_something_with_molecule(mol)
eigenvalues_verif = eigenspectrum(mol.fermionic_hamiltonian) # Should be the same as eigenvalues_original, but it is not.

# Doing other stuff with updated_mol...

The function doing_something_with_molecule is changing the MO coefficients of the sqmol_copy, but as the copy is shallow, it is also changing the MO coefficients of the mol object. This can introduce unexpected behaviour, as I experienced when working with FNO -> the MO coefficients were bound to the active space selection, so the fermionic Hamiltonian computed after the doing_something_with_molecule call wasn't trowing an error, but was physically wrong.

Here is a list of improvements that would upgrade the ONIOM offer:

• Geometry optimisation;
• Implementation of 3 or more layers (ex: FCI (high), CCSD (medium), RHF (low));
• Molecular mechanics solver;
• Semi-empirical solver(s);
• Capping broken bonds with functional groups instead of a single atom (ex: CH3).

Issue: Bug Report

Expected Behavior
Calling the inverse method on a S or T gates should return, at least in our implementation, a Phase gate with -pi/2 or -pi/4 as a parameter.

Current Behavior
It returns a S or T gate.

Steps to Reproduce (minimal example)

python
from tangelo.linq import translate_circuit, Gate
g = Gate("S", 1).inverse()
print(g.__dict__)

Possible Solution
This is do to the ordering of the conditionals in

EDIT: Parameters for S and T are also reversed.

The following code currently returns an error since UCC1 and 3 always return a circuit operating on 4 qubits (0 to 3), while the number of qubits in the Hamiltonian used for the expectation value may go beyond that.

Is this ansatz usable for any problem ? Should we indicate it in the documentation, or attempt to fix things so that it can be used on a broader variety of problem, if that makes sense at all ?

In particular, agnostic_simulator.Circuit can be instantiated with parameter n_qubits to "pad" a circuit until it reaches the desired number of qubits (we could do that in VQESolver.build, once we know the number of qubits of our Hamiltonian. This way fix some things, but I'm still not sure what makes sense for these ansatze (also, it means we're simulating qubits that don't need to be).

The code that crashes:

vqe_options = {"molecule": mol_H4, "ansatz": Ansatze.UCC1, "up_then_down": True}
vqe_solver = VQESolver(vqe_options)
vqe_solver.build()
vqe_solver.simulate()

Issue: Bug Report

Current Behavior
The H2 molecule in sto-3g basis with scBK qubit mapping should be solved with a single entangler with the UCCSD ansatz. Currently the algorithm seems to insert two copies of the same entangler.

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

from tangelo.algorithms.variational import VQESolver

vqe_options = {"molecule": mol_H2, "qubit_mapping": "scbk", "up_then_down":True}
vqe_solver = VQESolver(vqe_options)
vqe_solver.build()
vqe_solver.simulate()
print(vqe_solver.optimal_circuit)

returns

X         target : [0]   
X         target : [1]   
RX        target : [0]   parameter : 1.5707963267948966
RZ        target : [0]   parameter : 12.566261789411744	 (variational)
RX        target : [0]   parameter : 10.995574287564276
RX        target : [1]   parameter : 1.5707963267948966
RZ        target : [1]   parameter : 12.566261789411744	 (variational)
RX        target : [1]   parameter : 10.995574287564276
RX        target : [0]   parameter : 1.5707963267948966
H         target : [1]   
CNOT      target : [1]   control : [0]   
RZ        target : [1]   parameter : 0.11304437992941474	 (variational)
CNOT      target : [1]   control : [0]   
H         target : [1]   
RX        target : [0]   parameter : 10.995574287564276
H         target : [0]   
RX        target : [1]   parameter : 1.5707963267948966
CNOT      target : [1]   control : [0]   
RZ        target : [1]   parameter : 0.11304437992941474	 (variational)
CNOT      target : [1]   control : [0]   
RX        target : [1]   parameter : 10.995574287564276
H         target : [0]  

This circuit does return the right energy.

from tangelo.linq import Simulator 
sim=Simulator()
print(sim.get_expectation_value(vqe_solver.qubit_hamiltonian, vqe_solver.optimal_circuit))

returns -1.1372704157729252

Was this happening before? I can look into the BuiltInAnsatze.UCCSD or potentially the changes in the Circuit class

Here is a list where I am keeping documentation inconsistencies / small fixes in Tangelo. Before the next release, we could fix them all at once.

• The Pauli strings are reversed in (code is correct)

  Tangelo/tangelo/toolboxes/qubit_mappings/hcb.py

  Lines 67 to 69 in f4d46ad

  	In short, every creation operator b^{\dagger} (resp. annihilation b) are
  	mapped to X+iY strings (resp. X-iY), where X and Y are referring to the
  	Pauli matrices.

• lauren -> laurent
  

  Tangelo/tangelo/toolboxes/circuits/qsp.py

  Line 49 in 5d3e9da

  raise ModuleNotFoundError("pyqsp package is required to calculate QSP time-evolution phases using 'lauren' or 'tf' method.")

• It is not a pyscf molecule.
  

  Tangelo/tangelo/algorithms/variational/adapt_vqe_solver.py

  Line 138 in 3559208

  # Building molecule data with a pyscf molecule.

          **Update on the performance**

The code seems to be bottlenecked by the multiplication of QubitOperator. I tried several thing, like using MultiformOperator, the https://github.com/IntelLabs/mat2qubit package and several other options.

The very last thing I tried is splitting the multiplication with a divide-and-conquer strategy.

def element_to_qubitop(n_qubits, i, j, coeff=1.):
 
    # Must add 2 to the padding because of the "0b" prefix.
    bin_i = format(i, f"#0{n_qubits+2}b")
    bin_j = format(j, f"#0{n_qubits+2}b")

    qu_ops = [QubitOperator("", coeff)]
    for qubit, (bi, bj) in enumerate(zip(bin_i[2:][::-1], bin_j[2:][::-1])):
        if bi == "0" and bj == "0":
            qu_ops += [0.5 + QubitOperator(f"Z{qubit}", 0.5)]
        elif bi == "0" and bj == "1":
            qu_ops += [QubitOperator(f"X{qubit}", 0.5) + QubitOperator(f"Y{qubit}", 0.5j)]
        elif bi == "1" and bj == "0":
            qu_ops += [QubitOperator(f"X{qubit}", 0.5) + QubitOperator(f"Y{qubit}", -0.5j)]
        # The remaining case is 11.
        else:
            qu_ops += [0.5 + QubitOperator(f"Z{qubit}", -0.5)]

    qu_op = multiply_ops(qu_ops)

    return qu_op

def multiply_ops(qu_ops):

    if len(qu_ops) == 2:
        return qu_ops[0] * qu_ops[1]
    elif len(qu_ops) == 1:
        return qu_ops[0]
    else:
        return multiply_ops(qu_ops[:len(qu_ops)//2]) * multiply_ops(qu_ops[len(qu_ops)//2:])

However, this code is not speeding up things, in fact it is a little bit slower according to my manual tests. This suggest me that it is faster to do multiplication of big QubitOperator with a smaller one than doing the same thing with medium-size ones.

The next step is trying to leverage a faster language (like Julia or C). I have already begun working on an implementation using Julia.

Originally posted by @AlexandreF-1qbit in #286 (comment)

A picture (or a gif !) is worth a thousand words. Wouldn't it be nice if people could just see in practice how Tangelo works on simple tasks without having to dig into documentation, tutorials, etc ? They'd just "get it" and be able to decide if that's relevant to them or not.

Well, a gif in our README would be pretty neat. To win this:

• upload a gif in the docs folder on the repo
• alongside a simple text or rst file describing how you did it (don't worry about the README file, leave it alone). Mention the tools / website / techniques and the text used.

Here is what we want to see in that gif:

• a sequence of code snippets (provided below) gradually typed on screen, accompanied with the output
• try to keep the size of the text consistent across snippets.

Snippets and their results:

Snippet 1

# Run circuits on your favorite simulators

from tangelo.linq import Gate, Circuit, get_backend

c = Circuit([Gate('H', 0), Gate('X', 1)])
sim = get_backend("qiskit")  # or "qulacs", "cirq" ...
bitstring_frequencies, statevector = sim.simulate(c, return_statevector=True)

print(f'{bitstring_frequencies=} \n{statevector=}')

Snippet 2

# Use all the projects in the quantum ecosystem !

from tangelo.linq import translate_circuit

c_qiskit = translate_circuit(c, target="qiskit")  # to or from "qiskit", "pennylane", "cirq", ...
c_qiskit.draw(output='mpl', style={'backgroundcolor': '#EEEEEE'}) # it's a qiskit object !

(this is not the actual output, just run the above code with qiskit and pylatexenc in your environment)

Snippet 3

# Compute ground state of a molecule with VQE
# Check out the computational resources needed

from tangelo import SecondQuantizedMolecule
from tangelo.algorithms import VQESolver

xyz_H2 = [("H", (0., 0., 0.)), ("H", (0., 0., 0.7414))]
mol = SecondQuantizedMolecule(xyz_H2, q=0, spin=0, basis="sto-3g")

vqe_options = {"molecule": mol}  # Plenty options for qubit mappings etc
vqe_solver = VQESolver(vqe_options)
vqe_solver.build()
opt_energy = vqe_solver.simulate()

print(opt_energy)
print(vqe_solver.get_resources())

Snippet 4

from tangelo.linq.qpu_connection import IBMConnection

conn = IBMConnection(ibm_quantum_token='your_ibm_token')
job_id = conn.job_submit('sampler', 'ibmq_qasm_simulator', n_shots=100, circuit)
job_res = conn.job_results(job_id)  # Your job on IBM Quantum eventually completes

(no output needed)