diff --git a/depixel/bspline.py b/depixel/bspline.py
index ac7ab0e..68053a8 100644
--- a/depixel/bspline.py
+++ b/depixel/bspline.py
@@ -119,12 +119,19 @@ def __init__(self, knots, points, degree=None):
         self.degree = degree
         self._reset_cache()
 
-    def _reset_cache(self):
-        self._cache = {}
+    def _reset_cache(self, i=None):
+        ps = {}
+        if i is not None:
+            ps = self._cache['ps']
+            for ps_u in ps.values():
+                for j in range(len(ps_u)):
+                    [ps_u[j].pop(i - k, None) for k in range(j)]
+
+        self._cache = {'ps': ps}
 
     def move_point(self, i, value):
         self._points[i] = value
-        self._reset_cache()
+        self._reset_cache(i)
 
     def __str__(self):
         return "<%s degree=%s, points=%s, knots=%s>" % (
@@ -147,6 +154,23 @@ def points(self):
     def useful_points(self):
         return self.points
 
+    def _get_p(self, u, k, s):
+        ps = self._cache['ps'].setdefault(
+            u, [{} for _ in range(self.degree - s + 1)])
+
+        for i in range(k - self.degree, k - s + 1):
+            ps[0][i] = self.points[i]
+
+        def get_p(r, i):
+            if i not in ps[r]:
+                a = (u - self.knots[i]) / (self.knots[i + self.degree - r + 1]
+                                           - self.knots[i])
+                ps[r][i] = (1 - a) * get_p(r - 1, i - 1) + a * get_p(r - 1, i)
+
+            return ps[r][i]
+
+        return get_p(self.degree - 1, k - s)
+
     def __call__(self, u):
         """
         De Boor's Algorithm. Made out of more maths.
@@ -161,16 +185,8 @@ def __call__(self, u):
             if k == len(self.points):
                 k -= 1
             return self.points[k]
-        ps = [dict(zip(range(k - self.degree, k - s + 1),
-                       self.points[k - self.degree:k - s + 1]))]
 
-        for r in range(1, self.degree - s + 1):
-            ps.append({})
-            for i in range(k - self.degree + r, k - s + 1):
-                a = (u - self.knots[i]) / (self.knots[i + self.degree - r + 1]
-                                           - self.knots[i])
-                ps[r][i] = (1 - a) * ps[r - 1][i - 1] + a * ps[r - 1][i]
-        return ps[-1][k - s]
+        return self._get_p(u, k, s)
 
     def quadratic_bezier_segments(self):
         """