diff --git a/doc/tutorials/02-charged_system/02-charged_system-1.ipynb b/doc/tutorials/02-charged_system/02-charged_system-1.ipynb
index 9b3dc777b46..a9bb5e41d6d 100644
--- a/doc/tutorials/02-charged_system/02-charged_system-1.ipynb
+++ b/doc/tutorials/02-charged_system/02-charged_system-1.ipynb
@@ -19,7 +19,7 @@
     "#define EXTERNAL_FORCES\n",
     "#define MASS\n",
     "#define ELECTROSTATICS\n",
-    "#define LENNARD_JONES\n",
+    "#define WCA\n",
     "```"
    ]
   },
@@ -46,7 +46,7 @@
     "plt.ion()\n",
     "\n",
     "# Print enabled features\n",
-    "required_features = [\"EXTERNAL_FORCES\", \"MASS\", \"ELECTROSTATICS\", \"LENNARD_JONES\"]\n",
+    "required_features = [\"EXTERNAL_FORCES\", \"MASS\", \"ELECTROSTATICS\", \"WCA\"]\n",
     "espressomd.assert_features(required_features)\n",
     "print(espressomd.features())\n",
     "\n",
@@ -67,12 +67,8 @@
     "types       = {\"Anion\":          0, \"Cation\": 1}\n",
     "numbers     = {\"Anion\": n_ionpairs, \"Cation\": n_ionpairs}\n",
     "charges     = {\"Anion\":       -1.0, \"Cation\": 1.0}\n",
-    "lj_sigmas   = {\"Anion\":        1.0, \"Cation\": 1.0}\n",
-    "lj_epsilons = {\"Anion\":        1.0, \"Cation\": 1.0}\n",
-    "\n",
-    "WCA_cut = 2.**(1. / 6.)\n",
-    "lj_cuts     = {\"Anion\":  WCA_cut * lj_sigmas[\"Anion\"],\n",
-    "               \"Cation\": WCA_cut * lj_sigmas[\"Cation\"]}"
+    "wca_sigmas   = {\"Anion\":        1.0, \"Cation\": 1.0}\n",
+    "wca_epsilons = {\"Anion\":        1.0, \"Cation\": 1.0}"
    ]
   },
   {
@@ -139,12 +135,11 @@
    "metadata": {},
    "source": [
     "Before we can really start the simulation, we have to specify the\n",
-    "interactions between our particles. We already defined the Lennard-Jones parameters at the beginning,\n",
+    "interactions between our particles. We already defined the WCA parameters at the beginning,\n",
     "what is left is to specify the combination rule and to iterate over particle type pairs. For simplicity, \n",
     "we implement only the *Lorentz-Berthelot* rules. \n",
     "We pass our interaction pair to <tt>system.non_bonded_inter[\\*,\\*]</tt> and set the \n",
-    "pre-calculated LJ parameters <tt>epsilon</tt>, <tt>sigma</tt> and <tt>cutoff</tt>. With <tt>shift=\"auto\"</tt>,\n",
-    "we shift the interaction potential to the cutoff so that $U_\\mathrm{LJ}(r_\\mathrm{cutoff})=0$."
+    "pre-calculated WCA parameters <tt>epsilon</tt> and <tt>sigma</tt>."
    ]
   },
   {
@@ -165,14 +160,13 @@
     "    else:\n",
     "        return ValueError(\"No combination rule defined\")\n",
     "\n",
-    "# Lennard-Jones interactions parameters\n",
+    "# WCA interactions parameters\n",
     "for s in [[\"Anion\", \"Cation\"], [\"Anion\", \"Anion\"], [\"Cation\", \"Cation\"]]:\n",
-    "    lj_sig = combination_rule_sigma(\"Berthelot\", lj_sigmas[s[0]], lj_sigmas[s[1]])\n",
-    "    lj_cut = combination_rule_sigma(\"Berthelot\", lj_cuts[s[0]], lj_cuts[s[1]])\n",
-    "    lj_eps = combination_rule_epsilon(\"Lorentz\", lj_epsilons[s[0]], lj_epsilons[s[1]])\n",
+    "    wca_sig = combination_rule_sigma(\"Berthelot\", wca_sigmas[s[0]], wca_sigmas[s[1]])\n",
+    "    wca_eps = combination_rule_epsilon(\"Lorentz\", wca_epsilons[s[0]], wca_epsilons[s[1]])\n",
     "\n",
-    "    system.non_bonded_inter[types[s[0]], types[s[1]]].lennard_jones.set_params(\n",
-    "        epsilon=lj_eps, sigma=lj_sig, cutoff=lj_cut, shift=\"auto\")"
+    "    system.non_bonded_inter[types[s[0]], types[s[1]]].wca.set_params(\n",
+    "        epsilon=wca_eps, sigma=wca_sig)"
    ]
   },
   {
@@ -182,7 +176,7 @@
     "## 3 Equilibration\n",
     "\n",
     "With randomly positioned particles, we most likely have huge overlap and the strong repulsion will\n",
-    "cause the simulation to crash. The next step in our script therefore is a suitable LJ equilibration.\n",
+    "cause the simulation to crash. The next step in our script therefore is a suitable WCA equilibration.\n",
     "This is known to be a tricky part of a simulation and several approaches exist to reduce the particle overlap.\n",
     "Here, we use the steepest descent algorithm and cap the maximal particle displacement per integration step\n",
     "to 1% of sigma.\n",
@@ -197,8 +191,8 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "# Lennard-Jones Equilibration\n",
-    "max_sigma = max(lj_sigmas.values())\n",
+    "# WCA Equilibration\n",
+    "max_sigma = max(wca_sigmas.values())\n",
     "min_dist = 0.0\n",
     "system.minimize_energy.init(f_max=0, gamma=10.0, max_steps=10,\n",
     "                            max_displacement=max_sigma * 0.01)\n",
@@ -423,7 +417,7 @@
     "\n",
     "To make your life easier, don't try to equilibrate randomly positioned particles,\n",
     "but set them up in a crystal structure close to equilibrium. If you do it right,\n",
-    "you don't even need the Lennard-Jones equilibration. \n",
+    "you don't even need the WCA equilibration. \n",
     "To speed things up, don't go further than 1000 particles and use a P$^3$M accuracy of $10^{-2}$.\n",
     "Your RDF should look like the plot in figure 2. If you get stuck,\n",
     "you can look at the solution script <tt>/doc/tutorials/02-charged_system/scripts/nacl_units.py</tt> (or <tt>nacl_units_vis.py</tt> with visualization)."