diff --git a/HartreeFock.ipynb b/HartreeFock.ipynb
index 08a4c6a..aabaa85 100644
--- a/HartreeFock.ipynb
+++ b/HartreeFock.ipynb
@@ -749,6 +749,97 @@
     "let () = iteration m_C 20\n"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "\n",
+    "let ee_ints = AOBasis.ee_ints ao_basis;;\n",
+    "let m_C = MOBasis.mo_coef mo_basis;;\n",
+    "let n_ao = Mat.dim1 m_C ;;\n",
+    "\n",
+    "let sum a = \n",
+    "  Array.fold_left (fun accu x -> accu +. x) 0. a\n",
+    ";;\n",
+    "\n",
+    "let i = 1;;\n",
+    "let j = 2;;\n",
+    "\n",
+    "Array.map (fun a ->\n",
+    "   let v = \n",
+    "    Array.map (fun b ->\n",
+    "        let u = \n",
+    "          Array.map (fun e ->\n",
+    "            let t = Array.map (fun f ->\n",
+    "            m_C.{f,j} *. ERI.get_phys ee_ints a b e f\n",
+    "           ) (Util.array_range 1 n_ao)\n",
+    "           in m_C.{e,i} *. sum t\n",
+    "        ) (Util.array_range 1 n_ao)\n",
+    "        in m_C.{b,j} *. sum u\n",
+    "    ) (Util.array_range 1 n_ao)\n",
+    "   in m_C.{a,i} *. sum v\n",
+    ") (Util.array_range 1 n_ao)\n",
+    "|> sum\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "$\\langle \\phi_i \\phi_j | \\phi_i \\phi_j \\rangle =\n",
+    "\\sum_a c_{ai} \\langle \\chi_a \\phi_j | \\phi_i \\phi_j \\rangle =\n",
+    "\\sum_a \\sum_b c_{ai} c_{bj} \\langle \\chi_a \\chi_b | \\phi_i \\phi_j \\rangle$\n",
+    "$\\sum_a \\sum_b \\sum_e \\sum_f \\left(c_{ai} \\left(c_{bj} \\left(c_{ei} \\left(c_{fj} \\langle \\chi_a \\chi_b | \\chi_e \\chi_f \\rangle \\right)\\right)\\right)\\right) $"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "$\\phi_i ^2 = \\left( \\sum_a c_{ai} \\chi_a \\right)^2$\n",
+    "$= \\sum_a c_{ai} \\sum_b c_{bi} \\chi_a \\chi_b$"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "$\\phi_i ^2 - \\phi_j ^2 = \\left( \\sum_a c_{ai} \\chi_a \\right)^2 - \\left( \\sum_a c_{aj} \\chi_a \\right)^2$\n",
+    "$= \\sum_a c_{ai} \\sum_b c_{bi} \\chi_a \\chi_b - \\sum_a c_{aj} \\sum_b c_{bj} \\chi_a \\chi_b$"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "$[\\phi_i ^2 - \\phi_j ^2 | \\phi_i \\phi_j ] = [\\left( \\sum_a c_{ai} \\chi_i \\right)^2 - \\left( \\sum_a c_{aj} \\chi_j \\right)^2| \\phi_i \\phi_j ] $\n",
+    "$= \\sum_a c_{ai} \\sum_b c_{bi} [\\chi_a \\chi_b| \\phi_i \\phi_j ]  - \\sum_a c_{aj} \\sum_b c_{bj} [ \\chi_a \\chi_b| \\phi_i \\phi_j ] $\n",
+    "$= \\left(\\sum_a c_{ai} \\sum_b c_{bi}  - \\sum_a c_{aj} \\sum_b c_{bj} \\right) [ \\chi_a \\chi_b| \\phi_i \\phi_j ] $\n",
+    "$= \\left(\\sum_a c_{ai} \\sum_b c_{bi}  - \\sum_a c_{aj} \\sum_b c_{bj} \\right) \\sum_e \\sum_f c_{ei} c_{fj} [ \\chi_a \\chi_b| \\chi_e \\chi_f ] $"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "$\\langle \\phi_i | x | \\phi_j \\rangle$\n",
+    "$=\\langle \\sum_a c_{ia} \\chi_a | x | \\sum_b c_{jb} \\chi_b\\rangle$\n",
+    "$=\\sum_a \\sum_b c_{jb} c_{ia} \\langle  \\chi_a | x | \\chi_b\\rangle$"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "let ints = MOBasis.ee_ints mo_basis in\n",
+    "ERI.get_phys ints 1 2 1 2;;"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,