diff --git a/examples/example10 - growth_rates.ipynb b/examples/example10 - growth_rates.ipynb
index 3e6e78ce..c64c6046 100644
--- a/examples/example10 - growth_rates.ipynb	
+++ b/examples/example10 - growth_rates.ipynb	
@@ -246,6 +246,13 @@
     "Here we will be using **all** available metal modes, however, later on we will show how to exclude those (via additional options to `nl_solver`) with specific frequencies, timescales, coupling, etc to reduce computational runtime (without losing too much information)."
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -350,7 +357,7 @@
     "\n",
     "# multiply each metal current (from the eigenvector) with its corresponding Greens matrix and sum\n",
     "# (don't forget to omit the plasma current mode, i.e. the final element)\n",
-    "flux = np.sum(mode_currents[0:-1, np.newaxis, np.newaxis] * nonlinear_solver.vessel_modes_greens, axis=0)\n",
+    "flux = np.sum(mode_currents[0:-1, np.newaxis, np.newaxis] * nonlinear_solver.vessel_modes_greens[i], axis=0)\n",
     "\n",
     "# plot\n",
     "fig, ax = plt.subplots(1, 1, figsize=(5, 8), dpi=60)\n",
@@ -369,6 +376,13 @@
     "cbar = plt.colorbar(im, ax=ax, fraction=0.09)"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  },
   {
    "cell_type": "markdown",
    "metadata": {},
@@ -503,7 +517,7 @@
     "\n",
     "This parameter quantifies the ability of the active coil system to stabilise vertical displacements of the plasma. It is purely inductive as it does not contain any resistance terms and is given by\n",
     "\n",
-    "$$ m_s := \\lambda \\left[ -M_{m,m}^{-1} \\left( M_{m,m} + M_{m,y} \\frac{\\partial \\bold{I}_y}{\\partial \\bold{I}_m} \\Bigg|_{GS} \\right) \\right], $$\n",
+    "$$ m_s := \\lambda \\left[ -M_{m,m}^{-1} \\left( M_{m,m} + M_{m,y} \\frac{\\partial \\vec{I}_y}{\\partial \\vec{I}_m} \\Bigg|_{GS} \\right) \\right], $$\n",
     "\n",
     "where the operator $\\lambda[A]$ returns the largest non-negative eigenvalue of the matrix $A$. This is the counterpart to the \"resistive\" stability margin (i.e. the deformable timescale) that is described in the \"Removing the plasma equation\" - where in this case the plasma equation has been omitted (i.e. $\\dot{I_p} = 0$).\n",
     "\n",
@@ -530,7 +544,7 @@
     "The (generalised) Leuer stability parameter is defined as the ratio of stabilising force gradient to de-stabilising force gradient acting on the system:\n",
     "\n",
     "\\begin{equation}\n",
-    "    f := \\frac{f_{stab}}{f_{destab}} = \\frac{\\bold{I}_y^T M'_{y,m} M^{-1}_{m,m} M'_{m,y} \\bold{I}_y}{\\bold{I}_y^T M''_{y,m} \\bold{I}_m},\n",
+    "    f := \\frac{f_{stab}}{f_{destab}} = \\frac{\\vec{I}_y^T M'_{y,m} M^{-1}_{m,m} M'_{m,y} \\vec{I}_y}{\\vec{I}_y^T M''_{y,m} \\vec{I}_m},\n",
     "\\end{equation}\n",
     "\n",
     "where $M'$ and $M''$ denote the first and second derivatives of the mutual inductances with respect to the $Z$ coordinate, respectively. Other matrices are defined as they are above. We can find these derivatives using the relations:\n",
@@ -580,10 +594,10 @@
    "source": [
     "### Back to the deformable model: Removing the plasma equation\n",
     "\n",
-    "The deformable growth rates can also be computed in the case where we assume $\\bold{\\dot{I}}_p = 0$. For this, the final row/column of the state space matrices above are omitted and the system becomes:\n",
+    "The deformable growth rates can also be computed in the case where we assume $\\vec{\\dot{I}}_p = 0$. For this, the final row/column of the state space matrices above are omitted and the system becomes:\n",
     "\n",
     "\\begin{equation*}\n",
-    "    \\left( M_{m,m} + M_{m,y} \\frac{\\partial \\bold{I}_y}{\\partial \\bold{I}_m} \\Bigg|_{GS} \\right) \\dot{\\bold{I}}_m + R_{m,m} \\bold{I}_m = \\bold{V}_m - M_{m,y} \\frac{\\partial \\bold{I}_y}{\\partial \\bold{\\theta}}\\Bigg|_{GS} \\dot{\\bold{\\theta}},\n",
+    "    \\left( M_{m,m} + M_{m,y} \\frac{\\partial \\vec{I}_y}{\\partial \\vec{I}_m} \\Bigg|_{GS} \\right) \\dot{\\vec{I}}_m + R_{m,m} \\vec{I}_m = \\vec{V}_m - M_{m,y} \\frac{\\partial \\vec{I}_y}{\\partial \\vec{\\theta}}\\Bigg|_{GS} \\dot{\\vec{\\theta}},\n",
     "\\end{equation*}\n",
     "\n",
     "and we find the eigenvaluesas before. \n",
@@ -642,20 +656,20 @@
     "\\end{equation*}\n",
     "where $\\mathcal{I}_{act,act}$ is the identity matrix (with size equal to number of active coils). Supposing we retain $k$ columns of $P$ (including all of the active coils and a subset of the \"most strongly coupled\" passive structures), we can write the full vector of metal currents as \n",
     "\\begin{equation*}\n",
-    "    \\bold{I}_m = P \\bold{I}_d,\n",
+    "    \\vec{I}_m = P \\vec{I}_d,\n",
     "\\end{equation*}\n",
-    "where $\\bold{I}_d$ is the \"reduced\" set of $k$ mode currents.\n",
+    "where $\\vec{I}_d$ is the \"reduced\" set of $k$ mode currents.\n",
     "\n",
-    "Now, recall the state space matrix system from above. If we apply our expression for $\\bold{I}_m$, we get:\n",
+    "Now, recall the state space matrix system from above. If we apply our expression for $\\vec{I}_m$, we get:\n",
     "\n",
     "\\begin{equation*}\n",
     "\\begin{aligned}\n",
     "\\left[\\begin{array}{cc}\n",
-    " M_{m,m}P + M_{m,y} \\frac{\\partial \\bold{I}_y}{\\partial \\bold{I}_m} \\Bigg|_{GS}P & M_{m,y} \\frac{\\partial \\bold{I}_y}{\\partial I_p}\\Bigg|_{GS} \\\\ \n",
-    "\\frac{\\hat{\\bold{I}}_y}{R_p} \\left( M_{y,y} \\frac{\\partial \\bold{I}_y}{\\partial \\bold{I}_m} \\Bigg|_{GS}P + M_{y,m}P \\right) & \\frac{\\hat{\\bold{I}}_y}{R_p} M_{y,y} \\frac{\\partial \\bold{I}_y}{\\partial I_p}\\Bigg|_{GS}\n",
+    " M_{m,m}P + M_{m,y} \\frac{\\partial \\vec{I}_y}{\\partial \\vec{I}_d} \\Bigg|_{GS} & M_{m,y} \\frac{\\partial \\vec{I}_y}{\\partial I_p}\\Bigg|_{GS} \\\\ \n",
+    "\\frac{\\hat{\\vec{I}}_y}{R_p} \\left( M_{y,y} \\frac{\\partial \\vec{I}_y}{\\partial \\vec{I}_m} \\Bigg|_{GS}P + M_{y,m}P \\right) & \\frac{\\hat{\\vec{I}}_y}{R_p} M_{y,y} \\frac{\\partial \\vec{I}_y}{\\partial I_p}\\Bigg|_{GS}\n",
     "\\end{array}\\right]\n",
     "\\left[\\begin{array}{c}\n",
-    "\\dot{\\bold{I}}_d \\\\ \n",
+    "\\dot{\\vec{I}}_d \\\\ \n",
     "\\dot{I}_p\n",
     "\\end{array}\\right] +\n",
     "\n",
@@ -664,19 +678,19 @@
     "0 & 1\n",
     "\\end{array}\\right]\n",
     "\\left[\\begin{array}{c}\n",
-    "\\bold{I}_d \\\\ \n",
+    "\\vec{I}_d \\\\ \n",
     "I_p\n",
     "\\end{array}\\right]\n",
     "\n",
     "&=\n",
     "\\left[\\begin{array}{c}\n",
-    "\\bold{V}_m \\\\ \n",
+    "\\vec{V}_m \\\\ \n",
     "0 \n",
     "\\end{array}\\right]\n",
     "-\n",
     "\\left[\\begin{array}{c}\n",
-    "M_{m,y} \\frac{\\partial \\bold{I}_y}{\\partial \\bold{\\theta}}\\Bigg|_{GS} \\dot{\\bold{\\theta}} \\\\ \n",
-    "\\frac{\\hat{\\bold{I}}_y}{R_p} M_{y,y} \\frac{\\partial \\bold{I}_y}{\\partial \\bold{\\theta}}\\Bigg|_{GS} \\dot{\\bold{\\theta}}\n",
+    "M_{m,y} \\frac{\\partial \\vec{I}_y}{\\partial \\vec{\\theta}}\\Bigg|_{GS} \\dot{\\vec{\\theta}} \\\\ \n",
+    "\\frac{\\hat{\\vec{I}}_y}{R_p} M_{y,y} \\frac{\\partial \\vec{I}_y}{\\partial \\vec{\\theta}}\\Bigg|_{GS} \\dot{\\vec{\\theta}}\n",
     "\\end{array}\\right].\n",
     "\n",
     "\\end{aligned}\n",
@@ -698,9 +712,9 @@
     "The process by which the passive modes are chosen is as follows:\n",
     "1. A first \"estimate\" of\n",
     "   $$\n",
-    "   \\lambda_i = \\left\\| \\left( \\frac{\\partial \\mathbf{I}_y}{\\partial \\mathbf{I}_m} \\right)_{:,i} \\right\\|,\n",
+    "   \\lambda_i = \\left\\| \\left( \\frac{\\partial \\vec{I}_y}{\\partial \\vec{I}_m} \\right)_{:,i} \\right\\|,\n",
     "   $$\n",
-    "   is calculated for all passive structure modes (without fully solving the GS equation and calculating the true Jacobian). This is done by perturbing the currents in the modes (passive structures), updating the poloidal flux produced by these currents, and then evaluating the plasma current density $J_p(R,Z, \\psi)$ to obtain the perturbed $\\mathbf{I}_y$. Here, $\\psi$ is the total poloidal flux, i.e. the plasma and the metals contributions combined.\n",
+    "   is calculated for all passive structure modes (without fully solving the GS equation and calculating the true Jacobian). This is done by perturbing the currents in the modes (passive structures), updating the poloidal flux produced by these currents, and then evaluating the plasma current density $J_p(R,Z, \\psi)$ to obtain the perturbed $\\vec{I}_y$. Here, $\\psi$ is the total poloidal flux, i.e. the plasma and the metals contributions combined.\n",
     "\n",
     "2. The $\\lambda_i$ are then sorted into ascending order and with $\\hat{\\lambda} = \\max_i \\{ \\lambda_i\\}$ being the largest (strongest coupling).\n",
     "\n",
@@ -729,7 +743,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "Again, we can plot how these look using similar code as before, except we need to remember that now we've done a mode decomposition, the eigenmode currents $\\bold{I}_d$ need to be transformed back into the original metal currents $\\bold{I}_m$. This is done in the plotting cell below."
+    "Again, we can plot how these look using similar code as before, except we need to remember that now we've done a mode decomposition, the eigenmode currents $\\vec{I}_d$ need to be transformed back into the original metal currents $\\vec{I}_m$. This is done in the plotting cell below."
    ]
   },
   {
@@ -778,7 +792,7 @@
     "print(f\"Timescale = {np.real(timescales[i]):.2e} [s]\")\n",
     "\n",
     "# multiply each metal current (from the eigenvector) with its corresponding Greens matrix and sum\n",
-    "flux = np.sum(mode_currents[:, np.newaxis, np.newaxis] * nonlinear_solver_option_1.vessel_modes_greens, axis=0)\n",
+    "flux = np.sum(mode_currents[:, np.newaxis, np.newaxis] * nonlinear_solver_option_1.vessel_modes_greens[i], axis=0)\n",
     "\n",
     "# plot\n",
     "fig, ax = plt.subplots(1, 1, figsize=(5, 8), dpi=60)\n",
@@ -818,7 +832,7 @@
     "\n",
     "2. As in option 1, calculate a first estimate of  \n",
     "   $$\n",
-    "   \\lambda_i = \\left\\| \\left( \\frac{\\partial \\mathbf{I}_y}{\\partial \\mathbf{I}_m} \\right)_{:,i} \\right\\|\n",
+    "   \\lambda_i = \\left\\| \\left( \\frac{\\partial \\vec{I}_y}{\\partial \\vec{I}_m} \\right)_{:,i} \\right\\|\n",
     "   $$  \n",
     "   for all passive structure modes without fully solving the GS equation for each (this is done for flagged modes too). This is prior to the full Jacobian calculation later on.\n",
     "\n",
@@ -1234,7 +1248,7 @@
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "Python 3",
+   "display_name": "freegsnke",
    "language": "python",
    "name": "python3"
   },
@@ -1248,7 +1262,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.10.16"
+   "version": "3.10.20"
   }
  },
  "nbformat": 4,