diff --git a/simple_berggren.ipynb b/simple_berggren.ipynb
new file mode 100644
index 0000000..4e4e1d6
--- /dev/null
+++ b/simple_berggren.ipynb
@@ -0,0 +1,72 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "using Plots, LinearAlgebra\n",
+    "include(\"p_space.jl\")\n",
+    "\n",
+    "g0(R, p, q) = (exp(-(1/4)*(p + q)^2*R^2)*(-1 + exp(p*q*R^2))*R)/(2*sqrt(π))\n",
+    "g1(R, p, q) = (exp(-(1/4)*(p + q)^2*R^2)*(2 + p*q*R^2 + exp(p*q*R^2)*(-2 + p*q*R^2)))/(2*p*sqrt(π)*q*R)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "vertices = (0, 0.5 - 0.3im, 1, 6)\n",
+    "subdivisions = 64\n",
+    "p, w = get_mesh(vertices, subdivisions)\n",
+    "\n",
+    "# resonance example from my thesis\n",
+    "V_basis(p, q) = 2*g0(4, p, q) - 3*g0(2, p, q)\n",
+    "\n",
+    "basis_eig = eigen(get_H_matrix(V_basis, p, w))\n",
+    "basis = basis_eig.vectors .* w\n",
+    "\n",
+    "evals = basis_eig.values\n",
+    "E_target = 0.7\n",
+    "E = evals[argmin(norm.(evals .- E_target))]\n",
+    "\n",
+    "print(\"E = $E\")\n",
+    "scatter(real.(evals), imag.(evals), xlim = (0,2))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# ResonanceEC: Eq. (20)\n",
+    "V_system(c) = (p, q) -> c*(-5*g0(sqrt(3), p, q) + 2*g0(sqrt(10), p, q))\n",
+    "\n",
+    "H = get_H_matrix(V_system(0.45), p, w)\n",
+    "H_berggren = transpose(basis) * H * basis\n",
+    "\n",
+    "evals = eigvals(H)\n",
+    "scatter(real.(evals), imag.(evals), xlim = (0, 0.5))"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Julia 1.9.0",
+   "language": "julia",
+   "name": "julia-1.9"
+  },
+  "language_info": {
+   "file_extension": ".jl",
+   "mimetype": "application/julia",
+   "name": "julia",
+   "version": "1.9.0"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}