Generates tables of indices and matrix elements to build electronic Hamiltonians #809

agran2018 · 2020-09-17T20:15:34Z

Context:
This PR extends PennyLane-QChem functionalities to construct observables of many-body quantum systems.

Description of the Change:
The functions one_particle and two_particle have been implemented to construct the tables of indices and
matrix elements of arbitrary one- and two-particle operators, respectively. These functions
can be used in conjunction with the observable function to construct electronic hamiltonians.
For example, in the case of a molecular Hamiltonian:

from pennylane import qchem
from openfermion.hamiltonians import MolecularData

# load the HF elect. structure
hf_data = MolecularData(filename="./lih")

# create lists with the indices and matrix elements of the one- and two-particle operators
t, t_core = qchem.one_particle(hf_data.one_body_integrals, core=[0], active=[1, 2])
v, v_core = qchem.two_particle(hf_data.two_body_integrals, core=[0], active=[1, 2])

# build the second-quantized Hamiltonian and map it to the qubit basis
H = qchem.observable([t, v], init_term=t_core+v_core)

Benefits:
These functions allow to build any electronic Hamiltonian based on the representations of the one- and two-particle operators in the basis of single-particle states.

Possible Drawbacks:
None

Related GitHub Issues:
None

…of one-particle operators

codecov · 2020-09-17T20:25:40Z

Codecov Report

Merging #809 into master will increase coverage by 0.08%.
The diff coverage is 100.00%.

@@            Coverage Diff             @@
##           master     #809      +/-   ##
==========================================
+ Coverage   90.64%   90.72%   +0.08%     
==========================================
  Files         127      127              
  Lines        8275     8347      +72     
==========================================
+ Hits         7501     7573      +72     
  Misses        774      774

Impacted Files	Coverage Δ
qchem/pennylane_qchem/qchem/obs.py	`100.00% <100.00%> (ø)`

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 99bd629...eedceda. Read the comment docs.

ixfoduap · 2020-09-21T19:53:55Z

qchem/pennylane_qchem/qchem/obs.py

@@ -444,10 +444,308 @@ def particle_number(orbitals, mapping="jordan_wigner", wires=None):
    return observable([table], mapping=mapping, wires=wires)


+def one_particle(t_matrix_elements, core=None, active=None, cutoff=1.0e-12):


Let's avoid using more than one underscore ( _ ) when naming!
For example, here, we can rename the argument to just matrix_elements. It's clear from context (the function is called one-particle) that these are elements for one-particle operators

The reason I added the t_ prefix was to differentiate from the matrix_elements argument of the observable function. But I agree that it must be clear from the docstring.

ixfoduap · 2020-09-21T19:56:46Z

qchem/pennylane_qchem/qchem/obs.py

+    .. math::
+
+        \hat{T} = \sum_{\alpha, \beta} \langle \alpha \vert \hat{t} \vert \beta \rangle
+        [\hat{c}_{\alpha\uparrow}^\dagger \hat{c}_{\beta\uparrow} +


What is the meaning of the arrows in alpha and beta? I found this confusing while reading the docstring

This is explained just below in the docstring.
"the spin quantum numbers are indicated explicitly with the up/down arrows."
This notation is very convenient to stress that these operators do not act on the spin degrees of freedom and also to read the code.
See here https://pennylane--809.org.readthedocs.build/en/809/code/api/pennylane_qchem.qchem.obs.one_particle.html

qchem/pennylane_qchem/qchem/obs.py

ixfoduap · 2020-09-21T20:00:17Z

qchem/pennylane_qchem/qchem/obs.py

+        \hat{c}_{\alpha\downarrow}^\dagger \hat{c}_{\beta\downarrow}].
+
+    In the equation above the indices :math:`\alpha, \beta` run over the basis of spatial
+    orbitals :math:`\vert \alpha \rangle = \phi_\alpha(r)`. Since the operator :math:`t`


This equation |alpha> = phi(r) doesn't make sense, since we are equating two different mathematical objects, namely a vector and a function. We could either use <r|alpha> or, what I think is better, say "spatial orbitals |alpha> with associated wavefunctions phi(r)"

ixfoduap · 2020-09-21T20:12:04Z

qchem/pennylane_qchem/qchem/obs.py

+
+    If an active space is defined (see :func:`~.active_space`), the summation indices
+    run over the active orbitals and the contribution due to core orbitals is computed as
+    :math:`T_\mathrm{core} = 2 \sum_{\alpha\in \mathrm{core}}


Should this be a sum over alpha and beta? Currently only sum over alpha

No. This is correct except that there is a typo!, thanks!. See the notes below.

ixfoduap · 2020-09-21T20:50:16Z

qchem/pennylane_qchem/qchem/obs.py

+        )
+
+    # Indices of the matrix elements with absolute values >= cutoff
+    indices = np.nonzero(np.abs(t_matrix_elements) >= cutoff)


ixfoduap · 2020-09-21T22:06:21Z

qchem/pennylane_qchem/qchem/obs.py

+
+    In the equation above the indices :math:`\alpha, \beta, \gamma, \delta` run over the basis
+    of spatial orbitals :math:`\vert \alpha \rangle = \phi_\alpha(r)`. Since the operator
+    :math:`v` acts only on the spatial the spin quantum numbers are indicated explicitly


This sentence is missing a part, it now reads "acts only on the spatial the spin quantum...". I guess it should be "acts only on the spatial orbitals, the spin quantum..."

Good catch!, thanks @ixfoduap

ixfoduap · 2020-09-21T22:07:08Z

qchem/pennylane_qchem/qchem/obs.py

+        \hat{c}_{\beta\downarrow}^\dagger \hat{c}_{\gamma\downarrow} \hat{c}_{\delta\downarrow}~].
+
+    In the equation above the indices :math:`\alpha, \beta, \gamma, \delta` run over the basis
+    of spatial orbitals :math:`\vert \alpha \rangle = \phi_\alpha(r)`. Since the operator


Same comment as above regarding |alpha> = phi

ixfoduap · 2020-09-21T22:08:34Z

qchem/pennylane_qchem/qchem/obs.py

+    \langle \alpha, \beta \vert \hat{v} \vert \alpha, \beta \rangle`.
+
+    Args:
+        v_matrix_elements (array[float]): 4D NumPy array with the matrix elements


One question for me is: where is the user supposed to get access to these matrix elements? Is this a standard output from quantum chemistry packages?

Yes. It is not straightforward to get them printed out but there is always a possibility. In the particular cases of PySCF and Psi4 this is straightforward thanks to the interfaces with OpenFermion.

ixfoduap · 2020-09-21T22:09:02Z

qchem/pennylane_qchem/qchem/obs.py

+
+    **Example**
+
+    >>> v_matrix_elements = np.array([[[[ 6.82389533e-01, -1.45716772e-16],


remove underscores from this example usage :-)

soranjh · 2020-09-22T13:12:09Z

qchem/pennylane_qchem/qchem/obs.py

+            not correlated in the many-body wave function
+        active (list): indices of active orbitals, i.e., the orbitals used to
+            build the correlated many-body wave function
+        cutoff (float): Cutoff value for including matrix elements. The


Cutoff might be changed to cutoff for consistency.

Here I am following guidelines for PennyLane docstrings.
See here https://pennylane.readthedocs.io/en/stable/development/guide/documentation.html#punctuation-and-spelling

soranjh · 2020-09-22T13:14:23Z

qchem/pennylane_qchem/qchem/obs.py

+                )
+            )
+
+        # Compute contribution due to core orbitals


Is this comment necessary?

I think it helps reading the code. If it is not a violation of the guidelines I would vote for keeping it.

It's fine to keep, but usually the need for comments is a sign that the code could be made clearer, so as to not need comments

It's fine to keep, but usually the need for comments is a sign that the code could be made clearer, so as to not need comments

Completely tangential, but while this might be the case in a lot of 'commodity' code, I find it's usually the opposite with scientific coding. Scientific code always has that extra layer of needing to interpret (potentially highly complex) equations, and implement them in a highly optimized manner.

So I always welcome commenting like this one here 🙂

Yes. Scientific codes without comments are very difficult to read/maintain. This is the case of GAMESS, for example. More recent packages are much easier to develop.

Co-authored-by: ixfoduap <40441298+ixfoduap@users.noreply.github.com>

agran2018

Thanks @ixfoduap for all suggestions!. Left few notes answering some of your questions.

agran2018 · 2020-09-22T15:49:38Z

qchem/pennylane_qchem/qchem/obs.py

+    .. math::
+
+        \hat{T} = \sum_{\alpha, \beta} \langle \alpha \vert \hat{t} \vert \beta \rangle
+        [\hat{c}_{\alpha\uparrow}^\dagger \hat{c}_{\beta\uparrow} +


This is explained just below in the docstring.
"the spin quantum numbers are indicated explicitly with the up/down arrows."
This notation is very convenient to stress that these operators do not act on the spin degrees of freedom and also to read the code.
See here https://pennylane--809.org.readthedocs.build/en/809/code/api/pennylane_qchem.qchem.obs.one_particle.html

agran2018 · 2020-09-22T16:28:46Z

qchem/pennylane_qchem/qchem/obs.py

+
+    If an active space is defined (see :func:`~.active_space`), the summation indices
+    run over the active orbitals and the contribution due to core orbitals is computed as
+    :math:`T_\mathrm{core} = 2 \sum_{\alpha\in \mathrm{core}}


No. This is correct except that there is a typo!, thanks!. See the notes below.

agran2018 · 2020-09-22T16:30:41Z

qchem/pennylane_qchem/qchem/obs.py

+            not correlated in the many-body wave function
+        active (list): indices of active orbitals, i.e., the orbitals used to
+            build the correlated many-body wave function
+        cutoff (float): Cutoff value for including matrix elements. The


Here I am following guidelines for PennyLane docstrings.
See here https://pennylane.readthedocs.io/en/stable/development/guide/documentation.html#punctuation-and-spelling

agran2018 · 2020-09-22T16:33:33Z

qchem/pennylane_qchem/qchem/obs.py

+                )
+            )
+
+        # Compute contribution due to core orbitals


I think it helps reading the code. If it is not a violation of the guidelines I would vote for keeping it.

agran2018 · 2020-09-22T16:34:46Z

qchem/pennylane_qchem/qchem/obs.py

+
+    In the equation above the indices :math:`\alpha, \beta, \gamma, \delta` run over the basis
+    of spatial orbitals :math:`\vert \alpha \rangle = \phi_\alpha(r)`. Since the operator
+    :math:`v` acts only on the spatial the spin quantum numbers are indicated explicitly


Good catch!, thanks @ixfoduap

agran2018 · 2020-09-22T16:37:07Z

qchem/pennylane_qchem/qchem/obs.py

+    \langle \alpha, \beta \vert \hat{v} \vert \alpha, \beta \rangle`.
+
+    Args:
+        v_matrix_elements (array[float]): 4D NumPy array with the matrix elements


Yes. It is not straightforward to get them printed out but there is always a possibility. In the particular cases of PySCF and Psi4 this is straightforward thanks to the interfaces with OpenFermion.

ixfoduap · 2020-09-22T21:41:56Z

qchem/pennylane_qchem/qchem/obs.py

+                )
+            )
+
+        # Compute contribution due to core orbitals


It's fine to keep, but usually the need for comments is a sign that the code could be made clearer, so as to not need comments

josh146

Hi @agran2018 --- looks great! I didn't do a thorough review however, just left comments on the parts that stood out to me. Happy for this to be merged once the other reviewers have approved.

josh146 · 2020-09-23T11:16:11Z

qchem/pennylane_qchem/qchem/obs.py

+                )
+            )
+
+        # Compute contribution due to core orbitals


It's fine to keep, but usually the need for comments is a sign that the code could be made clearer, so as to not need comments

Completely tangential, but while this might be the case in a lot of 'commodity' code, I find it's usually the opposite with scientific coding. Scientific code always has that extra layer of needing to interpret (potentially highly complex) equations, and implement them in a highly optimized manner.

So I always welcome commenting like this one here 🙂

josh146 · 2020-09-23T11:16:27Z

qchem/pennylane_qchem/qchem/obs.py

+     [ 2.          2.         -1.12554473]
+     [ 3.          3.         -1.12554473]]
+    >>> print(t_core)
+    -9.45478626


Great docstring @agran2018!

josh146 · 2020-09-23T11:18:16Z

qchem/pennylane_qchem/qchem/obs.py

+            )
+
+        # Compute contribution due to core orbitals
+        t_core = 2 * sum([matrix_elements[alpha, alpha] for alpha in core])


This looks like something that can be vectorized :)

not entirely sure what it would look like (you might have to play around with it), but something like the following:

2 * np.sum(matrix_elements[np.array(core), np.array(core)])

(some other permutation of the indexing might be required)

Thanks @josh146, this works fine.

josh146 · 2020-09-23T11:18:35Z

qchem/pennylane_qchem/qchem/obs.py

+        )
+
+    # Indices of the matrix elements with absolute values >= cutoff
+    indices = np.nonzero(np.abs(matrix_elements) >= cutoff)


nice use of numpy functions!

josh146 · 2020-09-23T11:20:15Z

qchem/pennylane_qchem/qchem/obs.py

+    pairs = [
+        [indices[0][i], indices[1][i]]
+        for i in range(num_indices)
+        if all(indices[j][i] in active for j in range(len(indices)))
+    ]
+
+    # Building the table of indices and matrix elements
+    table = np.zeros((2 * len(pairs), 3))
+
+    for i, pair in enumerate(pairs):
+        alpha, beta = pair
+        element = matrix_elements[alpha, beta]
+
+        # spin-up term
+        table[2 * i, 0] = 2 * active.index(alpha)
+        table[2 * i, 1] = 2 * active.index(beta)
+        table[2 * i, 2] = element
+
+        # spin-down term
+        table[2 * i + 1, 0] = 2 * active.index(alpha) + 1
+        table[2 * i + 1, 1] = 2 * active.index(beta) + 1
+        table[2 * i + 1, 2] = element


another potential vectorization 😆 It will probably be very ugly though, so it depends on the speed up needed

At this point we have already dropped the negligible matrix elements and selected only a few active orbitals. I don't think it is critical here and keeping it explicitly will be helpful to remember the spin index mappings.

Sounds good, thanks @agran2018!

trbromley

Looks awesome! I left some very minor/optional comments.

qchem/pennylane_qchem/qchem/obs.py

trbromley · 2020-09-23T21:07:45Z

qchem/pennylane_qchem/qchem/obs.py

+    if core is None:
+        t_core = 0
+    else:
+        if True in [i > orbitals - 1 or i < 0 for i in core]:


Minor thing: could it be any([i > orbitals - 1 or i < 0 for i in core]) which might be easier to read

qchem/pennylane_qchem/qchem/obs.py

Co-authored-by: Tom Bromley <49409390+trbromley@users.noreply.github.com>

agran2018 added 16 commits September 11, 2020 14:52

Add function 'one_particle' to generate the table of matrix elements …

43e198a

…of one-particle operators

Black 'one_particle'

c0db0da

Polish function 'one_particle'

264138d

Add unit tests for 'one_particle' function

f07ad99

Merge branch 'master' into qchem-table-matrix-elements

d379262

Add 'two_particle' function

2ed234d

Polish functions 'one_particle' and 'two_particle'

291c671

Add unit test for the 'two_particle' function

6ba5c3e

Update Example section of 'two_particle' function

fb739b8

Update docstring of unit tests for the 'one_particle' function

1a72c58

Add unit tests for the 'two_particle' function

22f5d44

Update 'test_one_particle'

a8c9ebe

Update error message in 'one_particle' function

4e65ab3

Merge branch 'master' into qchem-table-matrix-elements

60f74c2

Polish docstring

d32f940

Update 'CHANGELOG.md'

a6cbd10

agran2018 added review-ready 👌 PRs which are ready for review by someone from the core team. qchem ⚛️ Related to the QChem package labels Sep 17, 2020

agran2018 requested review from josh146, soranjh, ixfoduap and trbromley September 17, 2020 20:15

agran2018 self-assigned this Sep 17, 2020

Remove trailing spaces and update 'CHANGELOG.md'

f603981

agran2018 and others added 4 commits September 17, 2020 16:58

Polish the docstrings

13f979d

Merge branch 'master' into qchem-table-matrix-elements

f5d16e5

Merge branch 'master' into qchem-table-matrix-elements

b218557

Merge branch 'master' into qchem-table-matrix-elements

ed9e4d9

ixfoduap reviewed Sep 21, 2020

View reviewed changes

Merge branch 'master' into qchem-table-matrix-elements

9d74cda

soranjh reviewed Sep 22, 2020

View reviewed changes

Apply suggestions from code review

bf13233

Co-authored-by: ixfoduap <40441298+ixfoduap@users.noreply.github.com>

agran2018 commented Sep 22, 2020

View reviewed changes

agran2018 added 2 commits September 22, 2020 14:08

Apply suggestions from code review

91b37b7

Polish docstrings

76d5993

agran2018 requested review from ixfoduap and soranjh September 22, 2020 21:33

ixfoduap approved these changes Sep 22, 2020

View reviewed changes

josh146 reviewed Sep 23, 2020

View reviewed changes

Apply suggestions from code review

5bf5641

soranjh approved these changes Sep 23, 2020

View reviewed changes

Merge branch 'master' into qchem-table-matrix-elements

e0221f6

trbromley approved these changes Sep 23, 2020

View reviewed changes

agran2018 and others added 2 commits September 23, 2020 17:35

Apply suggestions from code review

32813c5

Co-authored-by: Tom Bromley <49409390+trbromley@users.noreply.github.com>

Apply suggestions from code review

eedceda

agran2018 merged commit 79d50e7 into master Sep 23, 2020

agran2018 deleted the qchem-table-matrix-elements branch September 23, 2020 21:50

		@@ -444,10 +444,308 @@ def particle_number(orbitals, mapping="jordan_wigner", wires=None):
		return observable([table], mapping=mapping, wires=wires)


		def one_particle(t_matrix_elements, core=None, active=None, cutoff=1.0e-12):


		Example

		>>> v_matrix_elements = np.array([[[[ 6.82389533e-01, -1.45716772e-16],

Generates tables of indices and matrix elements to build electronic Hamiltonians #809

Generates tables of indices and matrix elements to build electronic Hamiltonians #809

Conversation

agran2018 commented Sep 17, 2020 • edited

codecov bot commented Sep 17, 2020 • edited

Codecov Report

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

agran2018 Sep 22, 2020 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

soranjh Sep 22, 2020 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

agran2018 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

agran2018 Sep 22, 2020 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

josh146 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

trbromley left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

agran2018 commented Sep 17, 2020 •

edited

codecov bot commented Sep 17, 2020 •

edited

agran2018 Sep 22, 2020 •

edited

soranjh Sep 22, 2020 •

edited

agran2018 Sep 22, 2020 •

edited