diff --git a/FloryHugg_DiffUnbleached.ipynb b/FloryHugg_DiffUnbleached.ipynb
index 08b6006c1506659f0febddab18b0b9abaf835273..f9d2427fbc75ba1e75a53285192b6dca219d1883 100644
--- a/FloryHugg_DiffUnbleached.ipynb
+++ b/FloryHugg_DiffUnbleached.ipynb
@@ -16,7 +16,7 @@
     "# domain = ms.Sphere(df.Point(0, 0, 0), 1.0)\n",
     "# mesh = ms.generate_mesh(domain, 50)\n",
     "mesh = df.UnitIntervalMesh(10000)\n",
-    "dt = 0.000001\n",
+    "dt = 0.0001\n",
     "\n",
     "F = df.FunctionSpace(mesh, 'CG', 1)"
    ]
@@ -40,21 +40,21 @@
     "#              tc*df.inner(df.grad(c_tot), df.grad((1-c_tot)/c_tot*c/Ga0))) * df.dx\n",
     "    \n",
     "#     # Weak form radial symmetry:\n",
-    "#     form = ((df.inner((c-c0)/dt, tc*X[0]*X[0]) +\n",
-    "#              df.inner(df.grad(c), df.grad((1-c_tot)/Ga0*tc*X[0]*X[0]))) -\n",
-    "#              df.inner(df.grad(c_tot), df.grad((1-c_tot)/Ga0/c_tot*c*tc*X[0]*X[0]))-\n",
-    "#              tc*df.inner(df.grad(c), df.grad((1-c_tot)/Ga0*X[0]*X[0]))+\n",
-    "#              tc*df.inner(df.grad(c_tot), df.grad((1-c_tot)/c_tot*c/Ga0*X[0]*X[0]))-\n",
-    "#              (1-c_tot)/Ga0*2*X[0]*c.dx(0)*tc+\n",
-    "#              (1-c_tot)/Ga0/c_tot*c*2*X[0]*c_tot.dx(0)*tc) * df.dx\n",
-    "  # Weak form radial symmetry:\n",
     "    form = ((df.inner((c-c0)/dt, tc*X[0]*X[0]) +\n",
-    "             df.inner(df.grad(c), df.grad((1-c_tot+a*c_tot*c_tot)*tc*X[0]*X[0]))) -\n",
-    "             df.inner(df.grad(c_tot), df.grad((1-c_tot+a*c_tot*c_tot)/c_tot*c*tc*X[0]*X[0]))-\n",
-    "             tc*df.inner(df.grad(c), df.grad((1-c_tot+a*c_tot*c_tot)*X[0]*X[0]))+\n",
-    "             tc*df.inner(df.grad(c_tot), df.grad((1-c_tot+a*c_tot*c_tot)/c_tot*c*X[0]*X[0]))-\n",
-    "             (1-c_tot+a*c_tot*c_tot)*2*X[0]*c.dx(0)*tc+\n",
-    "             (1-c_tot+a*c_tot*c_tot)/c_tot*c*2*X[0]*c_tot.dx(0)*tc) * df.dx\n",
+    "             df.inner(df.grad(c), df.grad((1-c_tot)*Ga0*tc*X[0]*X[0]))) -\n",
+    "             df.inner(df.grad(c_tot), df.grad((1-c_tot)*Ga0/c_tot*c*tc*X[0]*X[0]))-\n",
+    "             tc*df.inner(df.grad(c), df.grad((1-c_tot)*Ga0*X[0]*X[0]))+\n",
+    "             tc*df.inner(df.grad(c_tot), df.grad((1-c_tot)/c_tot*c*Ga0*X[0]*X[0]))-\n",
+    "             (1-c_tot)*Ga0*2*X[0]*c.dx(0)*tc+\n",
+    "             (1-c_tot)*Ga0/c_tot*c*2*X[0]*c_tot.dx(0)*tc) * df.dx\n",
+    "  # Weak form radial symmetry:\n",
+    "#     form = ((df.inner((c-c0)/dt, tc*X[0]*X[0]) +\n",
+    "#              df.inner(df.grad(c), df.grad((1-c_tot+a*c_tot*c_tot)*tc*X[0]*X[0]))) -\n",
+    "#              df.inner(df.grad(c_tot), df.grad((1-c_tot+a*c_tot*c_tot)/c_tot*c*tc*X[0]*X[0]))-\n",
+    "#              tc*df.inner(df.grad(c), df.grad((1-c_tot+a*c_tot*c_tot)*X[0]*X[0]))+\n",
+    "#              tc*df.inner(df.grad(c_tot), df.grad((1-c_tot+a*c_tot*c_tot)/c_tot*c*X[0]*X[0]))-\n",
+    "#              (1-c_tot+a*c_tot*c_tot)*2*X[0]*c.dx(0)*tc+\n",
+    "#              (1-c_tot+a*c_tot*c_tot)/c_tot*c*2*X[0]*c_tot.dx(0)*tc) * df.dx\n",
     "    \n",
     "    t = 0\n",
     "    # Solve in time\n",
@@ -119,9 +119,11 @@
     "# c0_1 = calc_sim(c0_1, c_tot_1, Ga0_1)\n",
     "# c0_9 = calc_sim(c0_9, c_tot_9, Ga0_9)\n",
     "\n",
+    "# c0_1 = calc_sim(c0_1, c_tot1, 0)\n",
+    "# c0_2 = calc_sim(c0_2, c_tot2, a)\n",
     "\n",
-    "c0_1 = calc_sim(c0_1, c_tot1, 0)\n",
-    "c0_2 = calc_sim(c0_2, c_tot2, a)"
+    "c0_1 = calc_sim(c0_1, ct_1, g_1)\n",
+    "c0_2 = calc_sim(c0_2, ct_2, g_2)"
    ]
   },
   {
@@ -131,9 +133,9 @@
    "outputs": [],
    "source": [
     "# 1D:\n",
-    "plt.plot(np.linspace(0, 1, 10000), [c0_1([x]) for x in np.linspace(0, 1, 10000)])\n",
+    "plt.plot(np.linspace(0, 1, 10000), [1.31*c0_1([x]) for x in np.linspace(0, 1, 10000)])\n",
     "plt.plot(np.linspace(0, 1, 10000), [c0_2([x]) for x in np.linspace(0, 1, 10000)])\n",
-    "# plt.xlim(0.06, 0.125)\n",
+    "plt.xlim(0.06, 0.125)\n",
     "# plt.ylim(0, 0.3)\n",
     "# 3D:\n",
     "# plt.plot(np.linspace(0, 0.5, 1000), [c0([x, 0, 0]) for x in np.linspace(0, 0.5, 1000)])"
@@ -177,14 +179,16 @@
     "ct_2 = df.Function(F)\n",
     "c0_1 = df.Function(F)\n",
     "c0_2 = df.Function(F)\n",
+    "g_1 = df.Function(F)\n",
+    "g_2 = df.Function(F)\n",
     "\n",
     "ct_1.interpolate(df.Expression(p_tot(p1_i, p1_o), degree=1))\n",
     "ct_2.interpolate(df.Expression(p_tot(p2_i, p2_o), degree=1))\n",
     "P1 = c_tot1(0)/c_tot1(1)\n",
     "P2 = c_tot2(0)/c_tot2(1)\n",
     "\n",
-    "g_2 = df.Function(F)\n",
-    "g_2.interpolate(df.Expression(p_tot(p2)))\n",
+    "g_1.interpolate(df.Expression('1', degree=1))\n",
+    "g_2.interpolate(df.Expression(p_tot(p1_o/p2_o, p1_i/p2_i), degree=1))\n",
     "D_out1 = 1-c_tot1(1)\n",
     "a = (P2*D_out1/P1-1+p2_o)/p2_o**2\n",
     "\n",
@@ -199,6 +203,10 @@
    "outputs": [],
    "source": [
     "plt.plot(np.linspace(0, 1, 1000), [ct_1([x]) for x in np.linspace(0, 1, 1000)])\n",
+    "plt.plot(np.linspace(0, 1, 1000), [ct_2([x]) for x in np.linspace(0, 1, 1000)])\n",
+    "plt.show()\n",
+    "\n",
+    "plt.plot(np.linspace(0, 1, 1000), [g_2([x]) for x in np.linspace(0, 1, 1000)])\n",
     "# plt.plot(np.linspace(0, 1, 1000), [(1-c_tot1([x])-a*c_tot1([x])**2) for x in np.linspace(0, 1, 1000)])"
    ]
   },