diff --git a/.github/workflows/ci.yaml b/.github/workflows/ci.yaml
index 464b63e98..c752cd075 100644
--- a/.github/workflows/ci.yaml
+++ b/.github/workflows/ci.yaml
@@ -102,6 +102,8 @@ jobs:
 
       shell: ${{steps.filter.outputs.shell || steps.filter.outputs.ci}}
       shell_files: ${{steps.filter.outputs.shell_files}}
+
+      notebooks: ${{steps.filter.outputs.notebooks || steps.filter.outputs.ci}}
     steps:
       - name: Check out a copy of the OpenFermion git repository
         uses: actions/checkout@de0fac2e4500dabe0009e67214ff5f5447ce83dd # v6.0.2
@@ -135,6 +137,8 @@ jobs:
             shell:
               - '**/*.sh'
               - 'check/*'
+            notebooks:
+              - '**/*.ipynb'
 
   python-checks:
     if: needs.changes.outputs.python == 'true'
@@ -418,6 +422,32 @@ jobs:
           # shellcheck disable=SC2086
           shellcheck ${CHANGED_FILES}
 
+  notebook-checks:
+    if: needs.changes.outputs.notebooks == 'true'
+    name: Notebook format checks
+    needs: changes
+    runs-on: ubuntu-24.04
+    timeout-minutes: 10
+    steps:
+      - name: Check out a copy of the git repository
+        uses: actions/checkout@de0fac2e4500dabe0009e67214ff5f5447ce83dd # v6.0.2
+
+      - name: Set up Python and restore cache
+        uses: actions/setup-python@a309ff8b426b58ec0e2a45f0f869d46889d02405 # v6.2.0
+        with:
+          python-version: ${{inputs.python_ver || env.python_ver}}
+          cache: pip
+          cache-dependency-path: ${{env.python_dep_files}}
+
+      - name: Upgrade pip
+        run: python -m pip install --upgrade pip
+
+      - name: Install requirements
+        run: pip install -r dev_tools/requirements/envs/pytest.env.txt
+
+      - name: Check notebook format
+        run: check/nbformat
+
   report-results:
     name: CI
     if: always()
@@ -425,6 +455,7 @@ jobs:
       - coverage
       - docker-lint
       - json-lint
+      - notebook-checks
       - pytest
       - pytest-compat
       - pytest-extra
diff --git a/check/all b/check/all
index 76a03c436..1a09b79cc 100755
--- a/check/all
+++ b/check/all
@@ -99,6 +99,7 @@ fi
 run check/mypy
 run check/pytest-and-incremental-coverage "${rev[@]}"
 run check/shellcheck
+run check/nbformat
 
 # ~~~~ Summarize the results and exit with a status code ~~~~
 
diff --git a/check/nbformat b/check/nbformat
new file mode 100755
index 000000000..10eabd6e5
--- /dev/null
+++ b/check/nbformat
@@ -0,0 +1,78 @@
+#!/usr/bin/env bash
+# Copyright 2025 Google LLC
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#     https://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+set -euo pipefail
+
+declare -r usage="\
+Usage: ${0##*/} [--help] [--apply]
+
+Checks Jupyter notebook (.ipynb) format using tensorflow_docs.tools.nbfmt.
+
+Without --apply, exits 0 if all notebooks are formatted and 1 otherwise.
+With --apply, reformats notebooks in place.
+"
+
+# Get the working directory to the repo root.
+thisdir=$(dirname "${BASH_SOURCE[0]:?}")
+repo_dir=$(git -C "${thisdir}" rev-parse --show-toplevel)
+cd "${repo_dir}"
+
+opt_apply=0
+for arg in "$@"; do
+    case "${arg}" in
+        -h | --help)
+            echo "${usage}"
+            exit 0
+            ;;
+        --apply)
+            opt_apply=1
+            ;;
+        *)
+            echo "Unknown argument: '${arg}'" >&2
+            echo "See '${0} --help' for usage." >&2
+            exit 1
+            ;;
+    esac
+done
+
+if ! python -c "import tensorflow_docs.tools.nbfmt"; then
+    echo "tensorflow-docs is not installed." >&2
+    exit 2
+fi
+
+declare -a notebook_paths=(
+    docs
+    src/openfermion/resource_estimates/pbc/notebooks
+)
+
+my_nbfmt() {
+    python -m tensorflow_docs.tools.nbfmt --indent=1 "$@"
+}
+
+exit_code=0
+if output=$(my_nbfmt --test "${notebook_paths[@]}" 2>&1); then
+    echo "Notebooks are formatted."
+elif (( opt_apply )); then
+    my_nbfmt "${notebook_paths[@]}" &&
+    echo "Reformatted notebooks." ||
+    exit_code=$?
+else
+    echo "The following notebooks require formatting:" >&2
+    grep '^-' <<<"${output}" >&2 || echo "${output}" >&2
+    echo "Run 'check/nbformat --apply' to fix notebooks." >&2
+    exit_code=1
+fi
+
+exit "${exit_code}"
diff --git a/check/shellcheck b/check/shellcheck
index f5d36bfdf..1c796506c 100755
--- a/check/shellcheck
+++ b/check/shellcheck
@@ -57,6 +57,7 @@ required_shell_scripts=(
     check/all
     check/format-incremental
     check/mypy
+    check/nbformat
     check/pylint
     check/pylint-changed-files
     check/pytest
diff --git a/dev_tools/requirements/deps/pytest.txt b/dev_tools/requirements/deps/pytest.txt
index e6d43b3d5..0dca16a9f 100644
--- a/dev_tools/requirements/deps/pytest.txt
+++ b/dev_tools/requirements/deps/pytest.txt
@@ -6,3 +6,4 @@ pytest-retry
 pytest-xdist
 
 nbformat
+tensorflow-docs
diff --git a/dev_tools/requirements/envs/dev.env.txt b/dev_tools/requirements/envs/dev.env.txt
index 183e714a0..1de8ea73f 100644
--- a/dev_tools/requirements/envs/dev.env.txt
+++ b/dev_tools/requirements/envs/dev.env.txt
@@ -1,13 +1,17 @@
 #
-# This file is autogenerated by pip-compile with Python 3.10
+# This file is autogenerated by pip-compile with Python 3.13
 # by the following command:
 #
 #    pip-compile --output-file=envs/dev.env.txt deps/format.txt deps/mypy.txt deps/pip-tools.txt deps/pylint.txt deps/pytest.txt deps/resource_estimates_runtime.txt deps/runtime.txt deps/shellcheck.txt
 #
+absl-py==2.4.0
+    # via tensorflow-docs
 ase==3.28.0
     # via -r deps/resource_estimates_runtime.txt
 ast-serialize==0.5.0
     # via mypy
+astor==0.8.1
+    # via tensorflow-docs
 astroid==3.3.11
     # via pylint
 attrs==26.1.0
@@ -15,10 +19,6 @@ attrs==26.1.0
     #   cirq-core
     #   jsonschema
     #   referencing
-backports-asyncio-runner==1.2.0 ; python_full_version < "3.11"
-    # via
-    #   -r deps/runtime.txt
-    #   pytest-asyncio
 black==26.5.1
     # via -r deps/format.txt
 build==1.5.0
@@ -33,6 +33,13 @@ click==8.4.1
     # via
     #   black
     #   pip-tools
+colorama==0.4.6
+    # via
+    #   build
+    #   click
+    #   pylint
+    #   pytest
+    #   tqdm
 contourpy==1.3.2
     # via
     #   -r deps/runtime.txt
@@ -47,8 +54,6 @@ dill==0.4.1
     # via pylint
 duet==0.2.9
     # via cirq-core
-exceptiongroup==1.3.1
-    # via pytest
 execnet==2.1.2
     # via pytest-xdist
 fastjsonschema==2.21.2
@@ -71,6 +76,8 @@ jaxlib==0.4.38
     # via
     #   -r deps/resource_estimates_runtime.txt
     #   jax
+jinja2==3.1.6
+    # via tensorflow-docs
 jsonschema==4.26.0
     # via nbformat
 jsonschema-specifications==2025.9.1
@@ -81,6 +88,8 @@ kiwisolver==1.5.0
     # via matplotlib
 librt==0.11.0
     # via mypy
+markupsafe==3.0.3
+    # via jinja2
 matplotlib==3.10.9
     # via
     #   ase
@@ -100,7 +109,9 @@ mypy-extensions==1.1.0
     #   black
     #   mypy
 nbformat==5.10.4
-    # via -r deps/pytest.txt
+    # via
+    #   -r deps/pytest.txt
+    #   tensorflow-docs
 networkx==3.4.2
     # via
     #   -r deps/runtime.txt
@@ -152,6 +163,8 @@ pluggy==1.6.0
     # via
     #   pytest
     #   pytest-cov
+protobuf==7.35.0
+    # via tensorflow-docs
 pubchempy==1.0.5
     # via -r deps/runtime.txt
 pygments==2.20.0
@@ -192,6 +205,8 @@ pytokens==0.4.1
     # via black
 pytz==2026.2
     # via pandas
+pyyaml==6.0.3
+    # via tensorflow-docs
 referencing==0.37.0
     # via
     #   jsonschema
@@ -220,15 +235,8 @@ sympy==1.14.0
     # via
     #   -r deps/runtime.txt
     #   cirq-core
-tomli==2.4.1
-    # via
-    #   black
-    #   build
-    #   coverage
-    #   mypy
-    #   pip-tools
-    #   pylint
-    #   pytest
+tensorflow-docs==2025.2.19.33219
+    # via -r deps/pytest.txt
 tomlkit==0.15.0
     # via pylint
 tqdm==4.67.3
@@ -247,13 +255,8 @@ types-setuptools==82.0.0.20260518
     # via -r deps/mypy.txt
 typing-extensions==4.15.0
     # via
-    #   astroid
-    #   black
     #   cirq-core
-    #   exceptiongroup
     #   mypy
-    #   pytest-asyncio
-    #   referencing
 tzdata==2026.2
     # via pandas
 urllib3==2.7.0
diff --git a/dev_tools/requirements/envs/pytest.env.txt b/dev_tools/requirements/envs/pytest.env.txt
index c2f399485..8f6275d99 100644
--- a/dev_tools/requirements/envs/pytest.env.txt
+++ b/dev_tools/requirements/envs/pytest.env.txt
@@ -1,20 +1,23 @@
 #
-# This file is autogenerated by pip-compile with Python 3.10
+# This file is autogenerated by pip-compile with Python 3.13
 # by the following command:
 #
 #    pip-compile --constraint=envs/dev.env.txt --output-file=envs/pytest.env.txt deps/pytest.txt deps/runtime.txt
 #
+absl-py==2.4.0
+    # via
+    #   -c envs/dev.env.txt
+    #   tensorflow-docs
+astor==0.8.1
+    # via
+    #   -c envs/dev.env.txt
+    #   tensorflow-docs
 attrs==26.1.0
     # via
     #   -c envs/dev.env.txt
     #   cirq-core
     #   jsonschema
     #   referencing
-backports-asyncio-runner==1.2.0 ; python_full_version < "3.11"
-    # via
-    #   -c envs/dev.env.txt
-    #   -r deps/runtime.txt
-    #   pytest-asyncio
 certifi==2026.5.20
     # via
     #   -c envs/dev.env.txt
@@ -27,6 +30,11 @@ cirq-core==1.5.0
     # via
     #   -c envs/dev.env.txt
     #   -r deps/runtime.txt
+colorama==0.4.6
+    # via
+    #   -c envs/dev.env.txt
+    #   pytest
+    #   tqdm
 contourpy==1.3.2
     # via
     #   -c envs/dev.env.txt
@@ -48,10 +56,6 @@ duet==0.2.9
     # via
     #   -c envs/dev.env.txt
     #   cirq-core
-exceptiongroup==1.3.1
-    # via
-    #   -c envs/dev.env.txt
-    #   pytest
 execnet==2.1.2
     # via
     #   -c envs/dev.env.txt
@@ -76,6 +80,10 @@ iniconfig==2.3.0
     # via
     #   -c envs/dev.env.txt
     #   pytest
+jinja2==3.1.6
+    # via
+    #   -c envs/dev.env.txt
+    #   tensorflow-docs
 jsonschema==4.26.0
     # via
     #   -c envs/dev.env.txt
@@ -92,6 +100,10 @@ kiwisolver==1.5.0
     # via
     #   -c envs/dev.env.txt
     #   matplotlib
+markupsafe==3.0.3
+    # via
+    #   -c envs/dev.env.txt
+    #   jinja2
 matplotlib==3.10.9
     # via
     #   -c envs/dev.env.txt
@@ -104,6 +116,7 @@ nbformat==5.10.4
     # via
     #   -c envs/dev.env.txt
     #   -r deps/pytest.txt
+    #   tensorflow-docs
 networkx==3.4.2
     # via
     #   -c envs/dev.env.txt
@@ -142,6 +155,10 @@ pluggy==1.6.0
     #   -c envs/dev.env.txt
     #   pytest
     #   pytest-cov
+protobuf==7.35.0
+    # via
+    #   -c envs/dev.env.txt
+    #   tensorflow-docs
 pubchempy==1.0.5
     # via
     #   -c envs/dev.env.txt
@@ -192,6 +209,10 @@ pytz==2026.2
     # via
     #   -c envs/dev.env.txt
     #   pandas
+pyyaml==6.0.3
+    # via
+    #   -c envs/dev.env.txt
+    #   tensorflow-docs
 referencing==0.37.0
     # via
     #   -c envs/dev.env.txt
@@ -224,11 +245,10 @@ sympy==1.14.0
     #   -c envs/dev.env.txt
     #   -r deps/runtime.txt
     #   cirq-core
-tomli==2.4.1
+tensorflow-docs==2025.2.19.33219
     # via
     #   -c envs/dev.env.txt
-    #   coverage
-    #   pytest
+    #   -r deps/pytest.txt
 tqdm==4.67.3
     # via
     #   -c envs/dev.env.txt
@@ -242,9 +262,6 @@ typing-extensions==4.15.0
     # via
     #   -c envs/dev.env.txt
     #   cirq-core
-    #   exceptiongroup
-    #   pytest-asyncio
-    #   referencing
 tzdata==2026.2
     # via
     #   -c envs/dev.env.txt
diff --git a/docs/fqe/guide/introduction.ipynb b/docs/fqe/guide/introduction.ipynb
index 8c6191843..266c4280e 100644
--- a/docs/fqe/guide/introduction.ipynb
+++ b/docs/fqe/guide/introduction.ipynb
@@ -2,7 +2,9 @@
  "cells": [
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "b000658d0d2b"
+   },
    "source": [
     "##### Copyright 2020 The OpenFermion Developers"
    ]
@@ -10,7 +12,10 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "cellView": "form",
+    "id": "906e07f6e562"
+   },
    "outputs": [],
    "source": [
     "#@title Licensed under the Apache License, Version 2.0 (the \"License\");\n",
@@ -28,14 +33,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "e6e2b37a3cea"
+   },
    "source": [
     "# Introduction to FQE"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "761ea25c35e5"
+   },
    "source": [
     "<table class=\"tfo-notebook-buttons\" align=\"left\">\n",
     "  <td>\n",
@@ -55,7 +64,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "7240f48e0397"
+   },
    "source": [
     "OpenFermion-FQE is an emulator for quantum computing specialized for simulations of Fermionic many-body problems, where FQE stands for 'Fermionic quantum emulator.' By focusing on Fermionic physics, OpenFermion-FQE realizes far more efficient emulation of quantum computing than generic quantum emulators such as Cirq, both in computation and memory costs; the speed-up and improved memory footprint originate from the use of the spin and number symmetries as well as highly optimized special algorithms.\n",
     "\n",
@@ -69,7 +80,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "1b44f1bf44c3"
+   },
    "outputs": [],
    "source": [
     "try:\n",
@@ -81,7 +94,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "a9308de9e5d7"
+   },
    "outputs": [],
    "source": [
     "import fqe\n",
@@ -90,14 +105,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "4c902bfe31c5"
+   },
    "source": [
     "## The FQE Wavefunction"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "ce6f7941f2a4"
+   },
    "source": [
     "The `Wavefunction` is an interface to the objects that hold the actual wavefunction data.  As mentioned, the wavefunction is partitioned into sectors with fixed particle and $Sz$ quantum number.  This partitioning information is the necessary information for initializing a `Wavefunction` object.  \n",
     "\n",
@@ -107,7 +126,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "3db9fbaacd17"
+   },
    "outputs": [],
    "source": [
     "wfn = fqe.Wavefunction([[4, 4, 4], [4, 2, 4], [4, 0, 4], [4, -2, 4], [4, -4, 4]])"
@@ -115,21 +136,27 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "b8554465dd55"
+   },
    "source": [
     "This command initializes a wavefunction with the following block structure:"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "c142731388ad"
+   },
    "source": [
     "![Image of wf sectors](./wf_init_sectors.png)"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "daf8dd3e1448"
+   },
    "source": [
     "Each sector corresponds to a set of bit strings \n",
     "\n",
@@ -142,14 +169,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "ac7422ec233a"
+   },
    "source": [
     "## Methods to initialize wavefunctions"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "785013468950"
+   },
    "source": [
     "FQE wavefunctions can be initialized by calling the constructor directly."
    ]
@@ -157,7 +188,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "5194a4b1575a"
+   },
    "outputs": [],
    "source": [
     "wfn_fqe = fqe.Wavefunction([[2, -2, 4]], broken=None)"
@@ -165,7 +198,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "74e94543a067"
+   },
    "source": [
     "When wavefunctions are first created, they are initialized to empty values. We can see this by printing out the wavefunction. "
    ]
@@ -173,7 +208,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "88e63b8991c3"
+   },
    "outputs": [],
    "source": [
     "wfn_fqe.print_wfn()"
@@ -181,7 +218,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "a6e99cce0b75"
+   },
    "source": [
     "To set the values of a wavefunction, we can use the `set_wfn` method with a provided `strategy`."
    ]
@@ -189,7 +228,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "d61c657172b1"
+   },
    "outputs": [],
    "source": [
     "wfn_fqe.set_wfn(strategy=\"hartree-fock\")\n",
@@ -198,7 +239,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "a4953b85f7f4"
+   },
    "source": [
     "Users can access the wavefunction through the `get_sector` method.  This returns the entire matrix of data representing the specified sector of the wavefunction.  For example, we can grab the sector corresponding to $n = 2$ and $sz = -2$ by doing the following."
    ]
@@ -206,7 +249,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "3a099aa1414a"
+   },
    "outputs": [],
    "source": [
     "interesting_states = wfn_fqe.get_coeff((2, -2))\n",
@@ -215,7 +260,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "f53313dce74f"
+   },
    "source": [
     "Other than the `Wavefunction` constructor, several utility methods are available to initialize wavefunctions. The function `fqe.get_wavefunction` builds a wavefunction with definite particle number and spin."
    ]
@@ -223,7 +270,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "d01012da95d5"
+   },
    "outputs": [],
    "source": [
     "wfn_fqe = fqe.get_wavefunction(4, -2, 10)"
@@ -231,7 +280,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "4c0c5fcbd16e"
+   },
    "source": [
     "The function `fqe.get_wavefunction_multiple` constructs multiple wavefunctions with different particle number, spin, and orbital number."
    ]
@@ -239,7 +290,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "0de4e8538d37"
+   },
    "outputs": [],
    "source": [
     "wfn_fqe1, wfn_fqe2 = fqe.get_wavefunction_multiple([[4, 0, 10], [5, -5, 10]])"
@@ -247,7 +300,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "d37154398c2b"
+   },
    "source": [
     "There are also functions like `fqe.get_number_conserving_wavefunction` and  `fqe.get_spin_conserving_wavefunction` to get number or spin conserving wavefunctions, respectively."
    ]
@@ -255,7 +310,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "d6ce0faddcc8"
+   },
    "outputs": [],
    "source": [
     "# Get a spin conserving wavefunction.\n",
@@ -267,14 +324,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "49e30a484711"
+   },
    "source": [
     "### Conversions between FQE and Cirq wavefunction representations"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "d386dba3b341"
+   },
    "source": [
     "Wavefunctions on $n$ qubits in Cirq are represented by Numpy arrays with $2^n$ amplitudes."
    ]
@@ -282,7 +343,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "f6bebeaea688"
+   },
    "outputs": [],
    "source": [
     "nqubits = 4\n",
@@ -296,7 +359,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "e55bd63e10d8"
+   },
    "source": [
     "To convert from this representation to the FQE representation, the function `fqe.from_cirq` can be used."
    ]
@@ -304,7 +369,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "64241a2c5480"
+   },
    "outputs": [],
    "source": [
     "fqe_wfn = fqe.from_cirq(cirq_wfn, thresh=0.0001)\n",
@@ -313,14 +380,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "df2d4caef151"
+   },
    "source": [
     "> *Note*: The `thresh` argument is the value below which amplitudes are considered zero."
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "bcda03e151d0"
+   },
    "source": [
     "We can convert back to the Cirq representation using `fqe._to_cirq`."
    ]
@@ -328,7 +399,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "052cdc627835"
+   },
    "outputs": [],
    "source": [
     "cirq_wfn_from_fqe = fqe.to_cirq(fqe_wfn)\n",
@@ -340,7 +413,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "06f2c59885f6"
+   },
    "outputs": [],
    "source": [
     "assert np.allclose(cirq_wfn_from_fqe, cirq_wfn)"
@@ -348,21 +423,27 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "aa7e202a9975"
+   },
    "source": [
     "An important thing to note in these conversions is the ordering of the $\\alpha$ and $\\beta$ strings in the converted wavefunctions.  The FQE uses the OpenFermion convention of interleaved $\\alpha$ and $\\beta$ orbitals.  Thus when converting to Cirq we first convert each bitstring into an OpenFermion operator and then call normal ordering."
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "593b7a9e2b25"
+   },
    "source": [
     "### Printing and saving wavefunctions"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "2488c8d21178"
+   },
    "source": [
     "Printing is currently available as alpha beta strings followed by the coefficient as well as orbital occupation representation."
    ]
@@ -371,7 +452,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "5d67ec2a25ef"
    },
    "outputs": [],
    "source": [
@@ -383,7 +464,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "e8e97254580d"
    },
    "outputs": [],
    "source": [
@@ -393,21 +474,27 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "96d1f944ac19"
+   },
    "source": [
     "Wavefunctions can also be saved to disk using the `save` method which takes a filename and optional path."
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "a6096cef2e03"
+   },
    "source": [
     "# Action on Wavefunctions: Fermionic algebra operations and their unitaries on the state"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "86e96919e84f"
+   },
    "source": [
     "FermionOperators can be directly passed in to create a new wavefunction based on application of the operators.  The FermionOperators are parsed according to the interleaved $\\alpha$ $\\beta$ indexing of the spin-orbitals.  This means that odd index FermionOperators correspond to $\\beta$-spin orbitals and even are $\\alpha$-spin orbitals.  \n",
     "\n",
@@ -419,7 +506,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "4f217ab9ad12"
+   },
    "outputs": [],
    "source": [
     "from openfermion import FermionOperator, hermitian_conjugated\n",
@@ -431,14 +520,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "9767dc4eea74"
+   },
    "source": [
     "## Unitary operations"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "01ae90adda3f"
+   },
    "source": [
     "Any simulator backend must be able to perform unitary evolution on a state.  The FQE accomplishes this by implementing code for evolving a state via the action of a unitary generated by fermionic generators.  Given a fermion operator $g$, the unitary \n",
     "\n",
@@ -460,7 +553,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "1164bb0810b6"
+   },
    "outputs": [],
    "source": [
     "i, j, theta = 0, 1, np.pi / 3\n",
@@ -474,31 +569,24 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "f6607620f34c"
+   },
    "source": [
     "In other tutorials we will do a deeper dive into supported time-evolution operations.  To serve a full functioning time-evolution platform the FQE also implements arbitrary time-evolution of full Hamiltonian operators consisting of sums of non-commuting terms.  The Taylor and Chebyshev expansion methods are used to do the exact time evolution.  "
    ]
   }
  ],
  "metadata": {
+  "colab": {
+   "name": "introduction.ipynb",
+   "toc_visible": true
+  },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.8"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 4
+ "nbformat_minor": 0
 }
diff --git a/docs/fqe/tutorials/diagonal_coulomb_evolution.ipynb b/docs/fqe/tutorials/diagonal_coulomb_evolution.ipynb
index d97c1ba0f..11a4acb37 100644
--- a/docs/fqe/tutorials/diagonal_coulomb_evolution.ipynb
+++ b/docs/fqe/tutorials/diagonal_coulomb_evolution.ipynb
@@ -2,7 +2,9 @@
  "cells": [
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "b000658d0d2b"
+   },
    "source": [
     "##### Copyright 2020 The OpenFermion Developers"
    ]
@@ -10,7 +12,10 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "cellView": "form",
+    "id": "906e07f6e562"
+   },
    "outputs": [],
    "source": [
     "#@title Licensed under the Apache License, Version 2.0 (the \"License\");\n",
@@ -28,14 +33,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "7cac3ba4e0b6"
+   },
    "source": [
     "# FQE vs OpenFermion vs Cirq: Diagonal Coulomb Operators"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "194e5ebc2e40"
+   },
    "source": [
     "<table class=\"tfo-notebook-buttons\" align=\"left\">\n",
     "  <td>\n",
@@ -55,7 +64,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "027a78b74a41"
+   },
    "source": [
     "Special routines are available for evolving under a diagonal Coulomb operator.  This notebook describes how to use these built in routines and how they work."
    ]
@@ -63,7 +74,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "1b44f1bf44c3"
+   },
    "outputs": [],
    "source": [
     "try:\n",
@@ -75,7 +88,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "a623fa4c377a"
+   },
    "outputs": [],
    "source": [
     "from itertools import product\n",
@@ -92,7 +107,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "07eb37e66cea"
+   },
    "outputs": [],
    "source": [
     "#Utility function\n",
@@ -151,7 +168,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "4f4ed62f910c"
+   },
    "source": [
     "The first example we will perform is diagonal Coulomb evolution on the Hartree-Fock state.  The diagonal Coulomb operator is defined as\n",
     "\n",
@@ -170,7 +189,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "2d22924fcb4e"
+   },
    "outputs": [],
    "source": [
     "norbs = 4\n",
@@ -192,7 +213,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "f4b8c93472be"
+   },
    "source": [
     "Now we can define a random 2-electron operator $V$.  To define $V$ we need a $4 \\times 4$ matrix.  We will generate this matrix by making a full random two-electron integral matrix and then just take the diagonal elements"
    ]
@@ -200,7 +223,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "56aa8604a016"
+   },
    "outputs": [],
    "source": [
     "tei_compressed = np.random.randn(tedim**2).reshape((tedim, tedim))\n",
@@ -221,7 +246,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "5270f0939273"
+   },
    "source": [
     "Evolution under $V$ can be computed by looking at each bitstring, seeing if $n_{p\\alpha}n_{q\\beta}$ is non-zero and then phasing that string by $V_{pq}$.  For the Hartree-Fock state we can easily calculate this phase accumulation.  The alpha and beta bitstrings are \"0001\" and \"0001\". "
    ]
@@ -229,7 +256,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "d7ee0994bdc2"
+   },
    "outputs": [],
    "source": [
     "alpha_occs = [list(range(fci_graph.nalpha()))]\n",
@@ -251,7 +280,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "9797ff1de370"
+   },
    "source": [
     "We can now try this out for more than 2 electrons.  Let's reinitialize a wavefunction on 6-orbitals with 4-electrons $S_{z} = 0$ to a random state."
    ]
@@ -259,7 +290,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "7a26d169c115"
+   },
    "outputs": [],
    "source": [
     "norbs = 6\n",
@@ -278,7 +311,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "89354765d52a"
+   },
    "source": [
     "We need to build our Diagoanl Coulomb operator For this bigger system."
    ]
@@ -286,7 +321,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "eec8fa7be891"
+   },
    "outputs": [],
    "source": [
     "tei_compressed = np.random.randn(tedim**2).reshape((tedim, tedim))\n",
@@ -307,7 +344,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "0bab89d13edf"
+   },
    "source": [
     "Now we can convert our wavefunction to a cirq wavefunction, evolve under the diagonal_coulomb operator we constructed and then compare the outputs."
    ]
@@ -315,7 +354,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "e28ac36062b7"
+   },
    "outputs": [],
    "source": [
     "cirq_wfn = fqe.to_cirq(fqe_wfn).reshape((-1, 1))\n",
@@ -327,7 +368,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "5a1a81b3c86a"
+   },
    "outputs": [],
    "source": [
     "fqe_wfn = fqe_wfn.time_evolve(1, dc_ham)\n",
@@ -338,7 +381,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "866925b171a5"
+   },
    "outputs": [],
    "source": [
     "print(\"From Cirq Evolution\")\n",
@@ -350,7 +395,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "7108f16f2401"
+   },
    "source": [
     "Finally, we can compare against evolving each term of $V$ individually."
    ]
@@ -358,7 +405,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "8a5588d11373"
+   },
    "outputs": [],
    "source": [
     "fqe_wfn = fqe.Wavefunction([[n_elec, sz, norbs]])\n",
@@ -372,34 +421,18 @@
     "               fqe_wfn.get_coeff((n_elec, sz)))\n",
     "print(\"Individual term evolution is equivalent\")"
    ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": []
   }
  ],
  "metadata": {
+  "colab": {
+   "name": "diagonal_coulomb_evolution.ipynb",
+   "toc_visible": true
+  },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.8"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 4
+ "nbformat_minor": 0
 }
diff --git a/docs/fqe/tutorials/fqe_vs_openfermion_quadratic_hamiltonians.ipynb b/docs/fqe/tutorials/fqe_vs_openfermion_quadratic_hamiltonians.ipynb
index 97eaa33de..7a0f317b8 100644
--- a/docs/fqe/tutorials/fqe_vs_openfermion_quadratic_hamiltonians.ipynb
+++ b/docs/fqe/tutorials/fqe_vs_openfermion_quadratic_hamiltonians.ipynb
@@ -2,7 +2,9 @@
  "cells": [
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "b000658d0d2b"
+   },
    "source": [
     "##### Copyright 2020 The OpenFermion Developers"
    ]
@@ -10,7 +12,10 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "cellView": "form",
+    "id": "906e07f6e562"
+   },
    "outputs": [],
    "source": [
     "#@title Licensed under the Apache License, Version 2.0 (the \"License\");\n",
@@ -28,14 +33,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "7cb64bac2714"
+   },
    "source": [
     "# FQE vs OpenFermion vs Cirq: Quadratic Hamiltonian Evolution"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "8cea787bddc0"
+   },
    "source": [
     "<table class=\"tfo-notebook-buttons\" align=\"left\">\n",
     "  <td>\n",
@@ -55,7 +64,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "37235f2b18ee"
+   },
    "source": [
     "This notebook demonstrates how evolve a state under a quadratic Hamiltonian with FQE, Cirq, and OpenFermion.\n",
     "\n",
@@ -65,7 +76,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "1b44f1bf44c3"
+   },
    "outputs": [],
    "source": [
     "try:\n",
@@ -77,7 +90,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "f84b86b1d410"
+   },
    "outputs": [],
    "source": [
     "from typing import Union\n",
@@ -103,7 +118,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "32b50215aa8b"
+   },
    "outputs": [],
    "source": [
     "# set up the number of orbs and a random quadratic Hamiltonian\n",
@@ -115,7 +132,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "1cbcb8dfed3f"
+   },
    "outputs": [],
    "source": [
     "# set up initial full Hilbert space wavefunction\n",
@@ -126,7 +145,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "dff29dcd0134"
+   },
    "outputs": [],
    "source": [
     "# full space unitary\n",
@@ -137,7 +158,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "880d190f075e"
+   },
    "outputs": [],
    "source": [
     "# build a FQE operator corresponding to the Hamiltonian evolution\n",
@@ -148,7 +171,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "7790282deccb"
+   },
    "outputs": [],
    "source": [
     "# initialize the FQE wavefunction [n-elec, sz, norbs]\n",
@@ -162,7 +187,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "96c80609b6f5"
+   },
    "source": [
     "# Aside\n",
     "\n",
@@ -172,7 +199,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "543a2aa79f4a"
+   },
    "outputs": [],
    "source": [
     "# initialize wavefunction [n-elec, sz, norbs]\n",
@@ -192,7 +221,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "9fc2809420e8"
+   },
    "outputs": [],
    "source": [
     "fqe_wfn.print_wfn()  # before\n",
@@ -203,7 +234,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "0e67f43cb4af"
+   },
    "outputs": [],
    "source": [
     "final_wfn = hamiltonian_evolution @ wfn  # Full space time-evolution"
@@ -211,7 +244,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "a77ff82a140b"
+   },
    "source": [
     "## Compare 1-RDMs\n",
     "\n",
@@ -221,10 +256,11 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "90c2c9d56205"
+   },
    "outputs": [],
    "source": [
-    "\n",
     "def get_opdm(wfn: Union[sp.sparse.coo_matrix, sp.sparse.csr_matrix,\n",
     "                        sp.sparse.csc_matrix], nso: int) -> np.ndarray:\n",
     "    \"\"\"\n",
@@ -257,7 +293,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "370f800e4494"
+   },
    "outputs": [],
    "source": [
     "fqe_opdm = fqe_wfn.rdm('i^ j') # get the FQE-opdm\n",
@@ -300,7 +338,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "d210103b1b78"
+   },
    "source": [
     "# Givens rotations with FQE\n",
     "\n",
@@ -315,7 +355,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "716314c7285b"
+   },
    "outputs": [],
    "source": [
     "# Get the basis rotation circuit\n",
@@ -333,7 +375,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "7286f75b718d"
+   },
    "outputs": [],
    "source": [
     "# Iterate through each layer and time evolve by the appropriate fermion operators\n",
@@ -350,14 +394,15 @@
     "for idx, final_phase in enumerate(diagonal):\n",
     "    if not np.isclose(final_phase, 1.0):\n",
     "        op = of.FermionOperator(((2 * idx, 1), (2 * idx, 0)), -np.angle(final_phase))\n",
-    "        fqe_wfn = fqe_wfn.time_evolve(1.0, op)\n",
-    "\n"
+    "        fqe_wfn = fqe_wfn.time_evolve(1.0, op)\n"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "6e3ee3bbc751"
+   },
    "outputs": [],
    "source": [
     "# map the wavefunction to a cirq wavefunction and get the opdm\n",
@@ -371,24 +416,15 @@
   }
  ],
  "metadata": {
+  "colab": {
+   "name": "fqe_vs_openfermion_quadratic_hamiltonians.ipynb",
+   "toc_visible": true
+  },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.8"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 4
+ "nbformat_minor": 0
 }
diff --git a/docs/fqe/tutorials/hamiltonian_time_evolution_and_expectation_estimation.ipynb b/docs/fqe/tutorials/hamiltonian_time_evolution_and_expectation_estimation.ipynb
index 92697d4c2..2de518cbe 100644
--- a/docs/fqe/tutorials/hamiltonian_time_evolution_and_expectation_estimation.ipynb
+++ b/docs/fqe/tutorials/hamiltonian_time_evolution_and_expectation_estimation.ipynb
@@ -2,7 +2,9 @@
  "cells": [
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "b000658d0d2b"
+   },
    "source": [
     "##### Copyright 2020 The OpenFermion Developers"
    ]
@@ -10,7 +12,10 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "cellView": "form",
+    "id": "906e07f6e562"
+   },
    "outputs": [],
    "source": [
     "#@title Licensed under the Apache License, Version 2.0 (the \"License\");\n",
@@ -28,14 +33,18 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "12cd0703cdd5"
+   },
    "source": [
     "# Hamiltonian Time Evolution and Expectation Value Computation"
    ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "e8b4b6017009"
+   },
    "source": [
     "<table class=\"tfo-notebook-buttons\" align=\"left\">\n",
     "  <td>\n",
@@ -55,7 +64,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "704fa8bc7617"
+   },
    "source": [
     "This tutorial describes the FQE's capabilities for Hamiltonian time-evolution and expectation value estimation\n",
     "\n",
@@ -65,7 +76,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "1b44f1bf44c3"
+   },
    "outputs": [],
    "source": [
     "try:\n",
@@ -75,65 +88,66 @@
    ]
   },
   {
-    "cell_type": "code",
-    "metadata": {
-      "scrolled": true,
-      "id": "3i55Vs4SYsm5"
-    },
-    "source": [
-      "Print = True\n",
-      "from openfermion import FermionOperator, MolecularData\n",
-      "from openfermion.utils import hermitian_conjugated\n",
-      "import numpy\n",
-      "import fqe\n",
-      "\n",
-      "numpy.set_printoptions(floatmode='fixed', precision=6, linewidth=80, suppress=True)\n",
-      "numpy.random.seed(seed=409)"
-    ],
-    "execution_count": null,
-    "outputs": []
-  },
-  {
-    "cell_type": "code",
-    "metadata": {
-      "id": "obExleOZZSAq"
-    },
-    "source": [
-      "!curl -O https://raw.githubusercontent.com/quantumlib/OpenFermion-FQE/main/tests/unittest_data/build_lih_data.py"
-    ],
-    "execution_count": null,
-    "outputs": []
-  },
-  {
-    "cell_type": "code",
-    "metadata": {
-      "id": "E3RaAyBkZuov"
-    },
-    "source": [
-      "import build_lih_data"
-    ],
-    "execution_count": null,
-    "outputs": []
-  },
-  {
-    "cell_type": "code",
-    "metadata": {
-	       "id": "g8XMRUw6Zb8W"
-    },
-    "source": [
-      "h1e, h2e, wfn = build_lih_data.build_lih_data('energy')\n",
-      "lih_hamiltonian = fqe.get_restricted_hamiltonian(([h1e, h2e]))\n",
-      "lihwfn = fqe.Wavefunction([[4, 0, 6]])\n",
-      "lihwfn.set_wfn(strategy='from_data', raw_data={(4, 0): wfn})\n",
-      "if Print:\n",
-      "    lihwfn.print_wfn()"
-    ],
-    "execution_count": null,
-    "outputs": []
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "id": "3i55Vs4SYsm5"
+   },
+   "outputs": [],
+   "source": [
+    "Print = True\n",
+    "from openfermion import FermionOperator, MolecularData\n",
+    "from openfermion.utils import hermitian_conjugated\n",
+    "import numpy\n",
+    "import fqe\n",
+    "\n",
+    "numpy.set_printoptions(floatmode='fixed', precision=6, linewidth=80, suppress=True)\n",
+    "numpy.random.seed(seed=409)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "id": "obExleOZZSAq"
+   },
+   "outputs": [],
+   "source": [
+    "!curl -O https://raw.githubusercontent.com/quantumlib/OpenFermion-FQE/main/tests/unittest_data/build_lih_data.py"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "id": "E3RaAyBkZuov"
+   },
+   "outputs": [],
+   "source": [
+    "import build_lih_data"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "id": "g8XMRUw6Zb8W"
+   },
+   "outputs": [],
+   "source": [
+    "h1e, h2e, wfn = build_lih_data.build_lih_data('energy')\n",
+    "lih_hamiltonian = fqe.get_restricted_hamiltonian(([h1e, h2e]))\n",
+    "lihwfn = fqe.Wavefunction([[4, 0, 6]])\n",
+    "lihwfn.set_wfn(strategy='from_data', raw_data={(4, 0): wfn})\n",
+    "if Print:\n",
+    "    lihwfn.print_wfn()"
+   ]
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "0afb126bbbfb"
+   },
    "source": [
     "## Application of one- and two-body fermionic gates\n",
     "\n",
@@ -143,7 +157,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "cba870e1720c"
+   },
    "outputs": [],
    "source": [
     "# dummy geometry\n",
@@ -160,7 +176,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "a7cf9f35b47f"
    },
    "outputs": [],
    "source": [
@@ -181,7 +197,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "3817bccb304d"
+   },
    "source": [
     "## Exact evolution implementation of quadratic Hamiltonians\n",
     "\n",
@@ -194,7 +212,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": false
+    "id": "c4b6a6165223"
    },
    "outputs": [],
    "source": [
@@ -221,7 +239,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "eb25c93086c9"
+   },
    "source": [
     "Exact evolution of dense quadratic hamiltonians is supported.  Here is an evolution example using a spin restricted Hamiltonian on a number and spin conserving wavefunction"
    ]
@@ -229,7 +249,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "01e066ed856a"
+   },
    "outputs": [],
    "source": [
     "norb = 4 \n",
@@ -251,7 +273,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "c305659829ee"
+   },
    "source": [
     "The GSO Hamiltonian is for evolution of quadratic hamiltonians that are spin broken and number conserving."
    ]
@@ -260,7 +284,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "a7466488a4d4"
    },
    "outputs": [],
    "source": [
@@ -283,7 +307,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "9ff1417ce5a5"
+   },
    "source": [
     "The BCS hamiltonian evolves spin conserved and number broken wavefunctions."
    ]
@@ -292,7 +318,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "ef9ac0c1ef2f"
    },
    "outputs": [],
    "source": [
@@ -321,7 +347,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "2fd96475486a"
+   },
    "source": [
     "Exact Evolution Implementation of Diagonal Coulomb terms"
    ]
@@ -330,7 +358,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "d6d63f7dc1db"
    },
    "outputs": [],
    "source": [
@@ -355,7 +383,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "69adcf1e8f58"
+   },
    "source": [
     "Exact evolution of individual n-body anti-Hermitian generators"
    ]
@@ -364,7 +394,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "95f3bf5a0dda"
    },
    "outputs": [],
    "source": [
@@ -384,7 +414,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "7ac8b740a309"
+   },
    "source": [
     " Approximate evolution of sums of n-body generators\n",
     "\n",
@@ -395,7 +427,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": false
+    "id": "4304d19db236"
    },
    "outputs": [],
    "source": [
@@ -408,7 +440,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "fab3d2d5c632"
    },
    "outputs": [],
    "source": [
@@ -435,7 +467,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "517df9eb2ea8"
+   },
    "source": [
     "API for determining desired expectation values"
    ]
@@ -444,7 +478,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "932fda9dea0c"
    },
    "outputs": [],
    "source": [
@@ -464,7 +498,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "869e982eed37"
+   },
    "source": [
     "2.B.1 RDMs \n",
     "In addition to the above API higher order density matrices in addition to hole densities can be calculated."
@@ -474,7 +510,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "e46fb551ae48"
    },
    "outputs": [],
    "source": [
@@ -488,7 +524,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "f6ef4e9fa2d8"
+   },
    "source": [
     "2.B.2 Hamiltonian expectations (or any expectation values)"
    ]
@@ -497,7 +535,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "c651e1e1e4cd"
    },
    "outputs": [],
    "source": [
@@ -511,7 +549,9 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "f3ddbe06e5e6"
+   },
    "source": [
     "2.B.3 Symmetry operations"
    ]
@@ -520,7 +560,7 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "scrolled": true
+    "id": "a61ed4b79b6e"
    },
    "outputs": [],
    "source": [
@@ -536,24 +576,15 @@
   }
  ],
  "metadata": {
+  "colab": {
+   "name": "hamiltonian_time_evolution_and_expectation_estimation.ipynb",
+   "toc_visible": true
+  },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.8"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 2
+ "nbformat_minor": 0
 }
diff --git a/docs/tutorials/binary_code_transforms.ipynb b/docs/tutorials/binary_code_transforms.ipynb
index dad891066..eb2ef5a5c 100644
--- a/docs/tutorials/binary_code_transforms.ipynb
+++ b/docs/tutorials/binary_code_transforms.ipynb
@@ -293,22 +293,9 @@
   },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.10"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 1
+ "nbformat_minor": 0
 }
diff --git a/docs/tutorials/bosonic_operators.ipynb b/docs/tutorials/bosonic_operators.ipynb
index d20696732..babb63208 100644
--- a/docs/tutorials/bosonic_operators.ipynb
+++ b/docs/tutorials/bosonic_operators.ipynb
@@ -986,22 +986,9 @@
   },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.10"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 1
+ "nbformat_minor": 0
 }
diff --git a/docs/tutorials/circuits_1_basis_change.ipynb b/docs/tutorials/circuits_1_basis_change.ipynb
index 651c6dc37..b48268c98 100644
--- a/docs/tutorials/circuits_1_basis_change.ipynb
+++ b/docs/tutorials/circuits_1_basis_change.ipynb
@@ -294,22 +294,9 @@
   },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.10"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 1
+ "nbformat_minor": 0
 }
diff --git a/docs/tutorials/circuits_2_diagonal_coulomb_trotter.ipynb b/docs/tutorials/circuits_2_diagonal_coulomb_trotter.ipynb
index 5ad5624e3..0bd2c99a8 100644
--- a/docs/tutorials/circuits_2_diagonal_coulomb_trotter.ipynb
+++ b/docs/tutorials/circuits_2_diagonal_coulomb_trotter.ipynb
@@ -481,22 +481,9 @@
   },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.10"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 1
+ "nbformat_minor": 0
 }
diff --git a/docs/tutorials/circuits_3_arbitrary_basis_trotter.ipynb b/docs/tutorials/circuits_3_arbitrary_basis_trotter.ipynb
index 78b97bbf4..a1e9f2611 100644
--- a/docs/tutorials/circuits_3_arbitrary_basis_trotter.ipynb
+++ b/docs/tutorials/circuits_3_arbitrary_basis_trotter.ipynb
@@ -281,22 +281,9 @@
   },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.10"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 1
+ "nbformat_minor": 0
 }
diff --git a/docs/tutorials/intro_to_openfermion.ipynb b/docs/tutorials/intro_to_openfermion.ipynb
index 0414e2567..2084b175b 100644
--- a/docs/tutorials/intro_to_openfermion.ipynb
+++ b/docs/tutorials/intro_to_openfermion.ipynb
@@ -894,22 +894,9 @@
   },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.10"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 1
+ "nbformat_minor": 0
 }
diff --git a/docs/tutorials/intro_workshop_exercises.ipynb b/docs/tutorials/intro_workshop_exercises.ipynb
index f23779c06..c4be74a09 100644
--- a/docs/tutorials/intro_workshop_exercises.ipynb
+++ b/docs/tutorials/intro_workshop_exercises.ipynb
@@ -2236,7 +2236,7 @@
    "collapsed_sections": [
     "Wp1a_Zf8EUen"
    ],
-   "name": "chemistry.ipynb",
+   "name": "intro_workshop_exercises.ipynb",
    "toc_visible": true
   },
   "kernelspec": {
diff --git a/docs/tutorials/jordan_wigner_and_bravyi_kitaev_transforms.ipynb b/docs/tutorials/jordan_wigner_and_bravyi_kitaev_transforms.ipynb
index 3eef5b4c2..740bcfd88 100644
--- a/docs/tutorials/jordan_wigner_and_bravyi_kitaev_transforms.ipynb
+++ b/docs/tutorials/jordan_wigner_and_bravyi_kitaev_transforms.ipynb
@@ -554,22 +554,9 @@
   },
   "kernelspec": {
    "display_name": "Python 3",
-   "language": "python",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.6.10"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 1
+ "nbformat_minor": 0
 }
diff --git a/src/openfermion/resource_estimates/pbc/notebooks/isdf.ipynb b/src/openfermion/resource_estimates/pbc/notebooks/isdf.ipynb
index 014c021ab..564c3ffee 100644
--- a/src/openfermion/resource_estimates/pbc/notebooks/isdf.ipynb
+++ b/src/openfermion/resource_estimates/pbc/notebooks/isdf.ipynb
@@ -3,7 +3,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "2c80abc258c1"
+   },
    "source": [
     "# Interpolative Seperable Density Fitting (ISDF) "
    ]
@@ -11,7 +13,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "46d00234458a"
+   },
    "source": [
     "The ISDF implementation currently provides a THC-like factorization of the two electron integrals which should converge to the FFTDF representation of the ERIs in the limit of large THC rank. This differs from the assumption of using RSGDF throughout the rest of the resource estimation scripts. However, we typically are only interested in ISDF as an initial guess for the THC factors which are then subsequently reoptimized to regularize $\\lambda$. The assumption here is that FFTDF / ISDF is a good enough approximation to the RSGDF ERIs and thus serves as a good initial guess.\n",
     "\n",
@@ -33,7 +37,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "706a3437e034"
+   },
    "outputs": [],
    "source": [
     "import os\n",
@@ -70,7 +76,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "47b66c6397f2"
+   },
    "source": [
     "Now let's find the ISDF THC factors using the KMeans-CVT algorithm to find the interpolating points. It's easiest to use the helper function `solve_kmeans_kpisdf` which will perform the necessary steps. "
    ]
@@ -78,7 +86,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "943ce7196a5d"
+   },
    "outputs": [],
    "source": [
     "import openfermion.resource_estimates.pbc.thc as kthc\n",
@@ -91,7 +101,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "48eaa5eb78da"
+   },
    "source": [
     "We see that the kpt_thc class has 4 attributes, `chi`, `xi`, `zeta` and `G_mapping`. `chi` corresponds to the cell periodic part of the Bloch orbital (i.e. $u_{p\\mathbf{k}}(\\mathbf{r}_\\mu))$. `xi` corresponds to $\\xi_{\\mu}(\\mathbf{r})$ in Eq. (1) above. To understand `zeta` and `G_mapping` it is helpful to recall we want to build\n",
     "\n",
@@ -109,7 +121,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "d603b7c193fd"
+   },
    "outputs": [],
    "source": [
     "from openfermion.resource_estimates.pbc.hamiltonian import build_momentum_transfer_mapping \n",
@@ -143,7 +157,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "06facefd6fab"
+   },
    "source": [
     "Now let's check convergence of the integral error, and corresponding MP2 error, with the THC dimension or equivalently the THC rank parameter."
    ]
@@ -151,7 +167,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "e71049388472"
+   },
    "outputs": [],
    "source": [
     "from openfermion.resource_estimates.pbc.hamiltonian.cc_extensions import compute_emp2_approx, build_cc_inst\n",
@@ -177,7 +195,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "80b58f25abcc"
+   },
    "source": [
     "Let's look at the convergence of the integral error"
    ]
@@ -185,7 +205,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "657a292b9342"
+   },
    "outputs": [],
    "source": [
     "import matplotlib.pyplot as plt\n",
@@ -199,7 +221,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "9d304a910c99"
+   },
    "source": [
     "Let's see how this corresponds to the MP2 error."
    ]
@@ -207,7 +231,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "c0895c6e0ed1"
+   },
    "outputs": [],
    "source": [
     "plt.cla()\n",
@@ -220,7 +246,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "be726ade0225"
+   },
    "source": [
     "We see that apart from some non-monotonic behaviour (which is expected due to the non-linear nature of the ISDF procedure), that a relatively large rank parameter is required to obtain say $< 0.1$ mHa error per cell. Note this could likely be reduced by carefully selecting the orbital sets we perform ISDF on as oo, ov, and vv blocks exhibit different low-rank behaviour, but for quantum algorithms this is not relevant. "
    ]
@@ -228,7 +256,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "dc31c542f689"
+   },
    "source": [
     "### Optional: Effect of Mesh Density\n",
     "\n",
@@ -240,7 +270,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "75e24f8cf796"
+   },
    "outputs": [],
    "source": [
     "plt.plot(thc_ranks, 100*np.array(thc_ranks)*num_spatial_orbs/(11**3), marker=\"o\")\n",
@@ -251,7 +283,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "27a120c10f35"
+   },
    "source": [
     "Now let us crank up the FFTDF accuracy and see if the results change significantly. This cell will take around 10 minutes to run."
    ]
@@ -259,10 +293,11 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "c461ea0c2ad7"
+   },
    "outputs": [],
    "source": [
-    "\n",
     "results = {11: [], 15: [], 19: [], 21: [], 28: []}\n",
     "#results = {11: [], 15: []}\n",
     "for mesh in list(results.keys()):\n",
@@ -308,9 +343,10 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "ac4d3f67419e"
+   },
    "source": [
-    "\n",
     "### Optional: Effect of Basis Set \n",
     "\n",
     "Another concern is that the basis set size is tiny so maybe things get worse as the basis set increases. Let's look into it by increasing the basis set size but still use a fairly coarse FFT grid, as we've seen it's not super important. This will also take several minutes to run.\n"
@@ -319,10 +355,11 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "b4f30404e56a"
+   },
    "outputs": [],
    "source": [
-    "\n",
     "basis_results = {\"gth-szv\": [], \"gth-dzvp\": [], \"gth-tzvp\": []}\n",
     "for basis in list(basis_results.keys()):\n",
     "    ase_atom = bulk(\"C\", \"diamond\", a=3.5)\n",
@@ -369,10 +406,11 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "d72642d90ae8"
+   },
    "outputs": [],
    "source": [
-    "\n",
     "for basis in list(basis_results.keys()):\n",
     "    plt.plot(thc_ranks, np.abs(basis_results[basis]), marker=\"o\", label=f\"basis = {basis}\")\n",
     "plt.yscale(\"log\")\n",
@@ -384,7 +422,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "09b912a14890"
+   },
    "source": [
     "## Effect on $\\lambda$ \n",
     "\n",
@@ -394,10 +434,11 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "b0844c27eb01"
+   },
    "outputs": [],
    "source": [
-    "\n",
     "ase_atom = bulk(\"C\", \"diamond\", a=3.5)\n",
     "cell = gto.Cell()\n",
     "cell.atom = pyscf_ase.ase_atoms_to_pyscf(ase_atom)\n",
@@ -443,7 +484,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "046319b60929"
+   },
    "source": [
     "Now we can try to improve $\\lambda$ by reoptimizing the THC factors using our ISDF factors as an initial guess. Practically this should mean we can use a smaller THC rank parameter for comparable MP2 accuracy. Note reoptimizing the THC factors is quite expensive. This cell may take 20 minutes to run. You should see that for a $c_\\mathrm{THC}=6$ the MP2 error is reduced by an order of magnitude. The optimization can be sped up by running on a GPU."
    ]
@@ -451,10 +494,11 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "72c4ab6084cb"
+   },
    "outputs": [],
    "source": [
-    "\n",
     "from openfermion.resource_estimates.pbc import utils\n",
     "# Recall we need RSGDF integrals to fit to.\n",
     "rsmf = scf.KRHF(mf.cell, mf.kpts).rs_density_fit()\n",
@@ -479,30 +523,15 @@
   }
  ],
  "metadata": {
+  "colab": {
+   "name": "isdf.ipynb",
+   "toc_visible": true
+  },
   "kernelspec": {
-   "display_name": "pyscf_pip",
-   "language": "python",
+   "display_name": "Python 3",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.11.4"
-  },
-  "orig_nbformat": 4,
-  "vscode": {
-   "interpreter": {
-    "hash": "61e04ed291edad8fb55208ca4954976506ff083febf358c737688b49027371c1"
-   }
   }
  },
  "nbformat": 4,
- "nbformat_minor": 2
+ "nbformat_minor": 0
 }
diff --git a/src/openfermion/resource_estimates/pbc/notebooks/resource_estimates.ipynb b/src/openfermion/resource_estimates/pbc/notebooks/resource_estimates.ipynb
index 20b0fce01..6bf2a3fd9 100644
--- a/src/openfermion/resource_estimates/pbc/notebooks/resource_estimates.ipynb
+++ b/src/openfermion/resource_estimates/pbc/notebooks/resource_estimates.ipynb
@@ -3,7 +3,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "f575d27f7440"
+   },
    "source": [
     "# Estimating Fault Tolerant Resources for Periodic Systems\n",
     "\n",
@@ -15,7 +17,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "9570df2a4539"
+   },
    "outputs": [],
    "source": [
     "import matplotlib.pyplot as plt\n",
@@ -54,7 +58,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "f59b6f8ef7b1"
+   },
    "source": [
     "Armed with our SCF solution we can now generate the one and two-electron integrals required to compute $\\lambda$."
    ]
@@ -62,7 +68,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "d36eb9ca42e5"
+   },
    "outputs": [],
    "source": [
     "from openfermion.resource_estimates.pbc.hamiltonian import build_hamiltonian\n",
@@ -79,7 +87,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "caae0e41d067"
+   },
    "source": [
     "# Hamiltonian Representation\n",
     "The resource estimation module provides four representations for the Hamiltonian: sparse, single-factorization (SF), double-factorization (SF), and tensor hypercontraction (THC). Each of these approaches introduces a different parameter which controls the accuracy of the factorization, and results need to be monitored with respect to these truncation parameters.\n",
@@ -112,7 +122,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "2a6ca742b64e"
+   },
    "outputs": [],
    "source": [
     "from openfermion.resource_estimates.pbc import sparse\n",
@@ -131,7 +143,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "832b23c73209"
+   },
    "source": [
     "With the Hamiltonian at hand we can compute the 1-norm (called $\\lambda$) of the Hamiltonian which is essential for computing resource estimates:"
    ]
@@ -139,7 +153,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "cc244a55f502"
+   },
    "outputs": [],
    "source": [
     "sparse_lambda = sparse.compute_lambda(hcore_mo, sparse_ham)\n",
@@ -149,7 +165,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "7a432f9e2a17"
+   },
    "source": [
     "Finally we can compute the total resource costs. In particular the following code will compute the number of Toffoli gates required for a single step of phase estimation (`toffolis_per_step`), the total Toffoli count (`total_toffolis`), and the number of logical qubits (`logical_qubits`).\n",
     "\n",
@@ -167,7 +185,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "037d545d5073"
+   },
    "outputs": [],
    "source": [
     "num_spin_orbs = 2*num_mo\n",
@@ -178,7 +198,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "293cc2fa0196"
+   },
    "source": [
     "We see there are the following steps required:\n",
     "\n",
@@ -195,7 +217,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "fb6341c69e6c"
+   },
    "outputs": [],
    "source": [
     "thresholds = np.logspace(-1, -5, 5)\n",
@@ -205,7 +229,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "a8a02370f36d"
+   },
    "source": [
     "The costing table is a `pandas.DataFrame` which can be convenient for saving to `csv` `(to_csv)`, string `(to_sting())` or outputting results to $\\LaTeX$ (`to_latex()`)."
    ]
@@ -213,7 +239,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "530218f1822a"
+   },
    "outputs": [],
    "source": [
     "print(sparse_costing_table.to_string(index=False))"
@@ -222,7 +250,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "a47d01eb89ac"
+   },
    "source": [
     "A crucial step in resource estimation is determining when the Hamiltonian is sufficiently accurate given a certain truncation. One way to check this is to use a correlated wavefunction method and monitor the convergence of the correlation energy with the truncation parameter. Let's look at the convergence of the MP2 error with the sparsity threshold. Note that MP2 is not expected to be a particularly faithful model to monitor convergence as it only requires a subset of the integral blocks. CCSD is a better option but comes with a considerable overhead."
    ]
@@ -230,7 +260,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "a86b462b7188"
+   },
    "outputs": [],
    "source": [
     "plt.plot(sparse_costing_table.cutoff, np.abs(sparse_costing_table.approx_energy-sparse_costing_table.exact_energy), marker=\"o\")\n",
@@ -243,7 +275,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "629cea0a9b01"
+   },
    "outputs": [],
    "source": [
     "plt.plot(sparse_costing_table.cutoff, np.abs(sparse_costing_table.lambda_total), marker=\"o\")\n",
@@ -255,7 +289,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "236da33eec3d"
+   },
    "outputs": [],
    "source": [
     "plt.plot(sparse_costing_table.cutoff, np.abs(sparse_costing_table.toffolis_per_step), marker=\"o\")\n",
@@ -267,7 +303,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "750982e4558a"
+   },
    "source": [
     "## Single Factorization \n",
     "\n",
@@ -310,7 +348,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "a342eb6b44a8"
+   },
    "outputs": [],
    "source": [
     "from openfermion.resource_estimates.pbc import sf \n",
@@ -321,7 +361,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "df7cd6d4807c"
+   },
    "outputs": [],
    "source": [
     "print(sf_costing_table.to_string(index=None))"
@@ -330,7 +372,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "fd8a8732f7fa"
+   },
    "outputs": [],
    "source": [
     "plt.plot(sf_costing_table.cutoff, np.abs(sf_costing_table.approx_energy-sf_costing_table.exact_energy), marker=\"o\")\n",
@@ -342,9 +386,10 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "1e72c4cb7572"
+   },
    "source": [
-    "\n",
     "## Double Factorization \n",
     "\n",
     "The double factorization follows on from the single-factorized Hamiltonian, and one performs a further factorization on the $A$ and $B$ operators defined above. \n",
@@ -360,7 +405,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "9ae37c3a0557"
+   },
    "outputs": [],
    "source": [
     "from openfermion.resource_estimates.pbc import df\n",
@@ -371,7 +418,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "64e0e9aeed5b"
+   },
    "outputs": [],
    "source": [
     "print(df_costing_table.to_string(index=False))"
@@ -380,7 +429,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "137fadd333c5"
+   },
    "outputs": [],
    "source": [
     "plt.rcParams['figure.dpi'] = 130 \n",
@@ -394,7 +445,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "9f74dbbba3d6"
+   },
    "source": [
     "## Tensor hypercontraction factorization\n",
     "\n",
@@ -404,7 +457,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "7407f17e3c21"
+   },
    "outputs": [],
    "source": [
     "from openfermion.resource_estimates.pbc import thc \n",
@@ -415,7 +470,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "4a4f279a3db8"
+   },
    "outputs": [],
    "source": [
     "print(thc_costing_table.to_string(index=False))"
@@ -424,7 +481,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "ca5f8bc86401"
+   },
    "outputs": [],
    "source": [
     "plt.rcParams['figure.dpi'] = 130 \n",
@@ -437,7 +496,9 @@
   {
    "attachments": {},
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "74229e6a1cdc"
+   },
    "source": [
     "## Other Considerations\n",
     "\n",
@@ -449,7 +510,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "495f12ac2752"
+   },
    "outputs": [],
    "source": [
     "from ase.build import bulk\n",
@@ -483,7 +546,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "610ff06485c0"
+   },
    "outputs": [],
    "source": [
     "cutoffs = np.logspace(-1, -5, 5)\n",
@@ -493,9 +558,10 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {},
+   "metadata": {
+    "id": "12367e81abd2"
+   },
    "source": [
-    "\n",
     "### Threshold convergence\n",
     "\n",
     "So far we have used MP2 to monitor the convergence of the factorization, however we can instead use CCSD as a better model chemistry."
@@ -504,7 +570,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "cdfe1f595615"
+   },
    "outputs": [],
    "source": [
     "thresholds = np.logspace(-1, -5, 5)\n",
@@ -514,7 +582,9 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {},
+   "metadata": {
+    "id": "28c3c63d0d17"
+   },
    "outputs": [],
    "source": [
     "print(sparse_costing_table.to_string(index=False))"
@@ -522,25 +592,15 @@
   }
  ],
  "metadata": {
+  "colab": {
+   "name": "resource_estimates.ipynb",
+   "toc_visible": true
+  },
   "kernelspec": {
-   "display_name": "kpoint_eri",
-   "language": "python",
+   "display_name": "Python 3",
    "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.11.4"
-  },
-  "orig_nbformat": 4
+  }
  },
  "nbformat": 4,
- "nbformat_minor": 2
+ "nbformat_minor": 0
 }