diff --git a/doc/source/models/PCSAFT.ipynb b/doc/source/models/PCSAFT.ipynb
index 973f43bd867d0fefdfb2d68d33553a5f96c9a32a..57dfb230b10335568c8628141519e79913e8c5d6 100644
--- a/doc/source/models/PCSAFT.ipynb
+++ b/doc/source/models/PCSAFT.ipynb
@@ -104,7 +104,51 @@
    "source": [
     "## Superancillary\n",
     "\n",
-    "The superancillary equation for PC-SAFT has been developed, and is much more involved than that of the cubic EOS. As a consequence, the superancillary equation has been provided as a separate package rather than integrating it into to teqp to minimize the binary size of teqp. It can be installed from PYPI with: ``pip install PCSAFTsuperanc``"
+    "The superancillary equation for PC-SAFT has been developed, and is much more involved than that of the cubic EOS. As a consequence, the superancillary equation has been provided as a separate package rather than integrating it into to teqp to minimize the binary size of teqp. It can be installed from PYPI with: ``pip install PCSAFTsuperanc``\n",
+    "\n",
+    "The scaling in the superancillaries uses reduced variables:\n",
+    "\n",
+    "$$ \\tilde T = T/(\\epsilon/k_{\\rm B}) $$\n",
+    "$$ \\tilde\\rho = \\rho_{\\rm N}\\sigma^3 $$\n",
+    "\n",
+    "where $\\rho_{\\rm N}$ is the number density, and the other parameters are from the PC-SAFT model"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "f6e3b8d2",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import PCSAFTsuperanc\n",
+    "\n",
+    "sigma_m = 3e-10 # [meter]\n",
+    "e_over_k = 150.0 # [K]\n",
+    "m = 5\n",
+    "\n",
+    "# The saturation temperature\n",
+    "T = 300\n",
+    "\n",
+    "[Ttilde_crit, Ttilde_min] = PCSAFTsuperanc.get_Ttilde_crit_min(m=m)\n",
+    "print('Ttilde crit:', Ttilde_crit)\n",
+    "\n",
+    "# Get the scaled densities for liquid and vapor phases\n",
+    "[tilderhoL, tilderhoV] = PCSAFTsuperanc.PCSAFTsuperanc_rhoLV(Ttilde=T/e_over_k, m=m)\n",
+    "# Convert back to molar densities\n",
+    "N_A = PCSAFTsuperanc.N_A # The value of Avogadro's constant used in superancillaries\n",
+    "rhoL, rhoV = [tilderho/(N_A*sigma_m**3) for tilderho in [tilderhoL, tilderhoV]]\n",
+    "\n",
+    "# As a sanity check, confirm that we got the same pressure in both phases\n",
+    "c = teqp.SAFTCoeffs()\n",
+    "c.sigma_Angstrom = sigma_m*1e10\n",
+    "c.epsilon_over_k = e_over_k \n",
+    "c.m = m\n",
+    "model = teqp.PCSAFTEOS([c])\n",
+    "z = np.array([1.0])\n",
+    "pL = rhoL*model.get_R(z)*T*(1+model.get_Ar01(T, rhoL, z))\n",
+    "pV = rhoV*model.get_R(z)*T*(1+model.get_Ar01(T, rhoV, z))\n",
+    "print('Pressures are:', pL, pV, 'Pa')"
    ]
   },
   {