diff --git a/manimlib/camera/camera.py b/manimlib/camera/camera.py
index 1c118623..3dfcd214 100644
--- a/manimlib/camera/camera.py
+++ b/manimlib/camera/camera.py
@@ -30,7 +30,7 @@ class CameraFrame(Mobject):
     }
 
     def init_points(self):
-        self.points = np.array([ORIGIN, LEFT, RIGHT, DOWN, UP])
+        self.set_points([ORIGIN, LEFT, RIGHT, DOWN, UP])
         self.set_width(self.frame_shape[0], stretch=True)
         self.set_height(self.frame_shape[1], stretch=True)
         self.move_to(self.center_point)
@@ -108,14 +108,19 @@ class CameraFrame(Mobject):
         return self.set_rotation(theta=self.euler_angles[2] + dgamma)
 
     def get_shape(self):
-        return (
-            self.points[2, 0] - self.points[1, 0],
-            self.points[4, 1] - self.points[3, 1],
-        )
+        return (self.get_width(), self.get_height())
 
     def get_center(self):
         # Assumes first point is at the center
-        return self.points[0]
+        return self.get_points()[0]
+
+    def get_width(self):
+        points = self.get_points()
+        return points[2, 0] - points[1, 0]
+
+    def get_height(self):
+        points = self.get_points()
+        return points[4, 1] - points[3, 1]
 
     def get_focal_distance(self):
         return self.focal_distance * self.get_height()
@@ -123,7 +128,8 @@ class CameraFrame(Mobject):
     def interpolate(self, frame1, frame2, alpha, path_func):
         self.euler_angles[:] = interpolate(frame1.euler_angles, frame2.euler_angles, alpha)
         self.refresh_camera_rotation_matrix()
-        self.points = interpolate(frame1.points, frame2.points, alpha)
+        # TODO, can probably safely call super
+        self.set_points(interpolate(frame1.get_points(), frame2.get_points(), alpha))
 
 
 class Camera(object):
diff --git a/manimlib/mobject/geometry.py b/manimlib/mobject/geometry.py
index 38e3e531..43496c6a 100644
--- a/manimlib/mobject/geometry.py
+++ b/manimlib/mobject/geometry.py
@@ -165,10 +165,10 @@ class TipableVMobject(VMobject):
         return self.tip_length
 
     def get_first_handle(self):
-        return self.points[1]
+        return self.get_points()[1]
 
     def get_last_handle(self):
-        return self.points[-2]
+        return self.get_points()[-2]
 
     def get_end(self):
         if self.has_tip():
@@ -234,7 +234,7 @@ class Arc(TipableVMobject):
         anchors, and finds their intersection points
         """
         # First two anchors and handles
-        a1, h, a2 = self.points[:3]
+        a1, h, a2 = self.get_points()[:3]
         # Tangent vectors
         t1 = h - a1
         t2 = h - a2
@@ -354,10 +354,10 @@ class AnnularSector(Arc):
             for radius in (self.inner_radius, self.outer_radius)
         ]
         outer_arc.reverse_points()
-        self.append_points(inner_arc.points)
-        self.add_line_to(outer_arc.points[0])
-        self.append_points(outer_arc.points)
-        self.add_line_to(inner_arc.points[0])
+        self.append_points(inner_arc.get_points())
+        self.add_line_to(outer_arc.get_points()[0])
+        self.append_points(outer_arc.get_points())
+        self.add_line_to(inner_arc.get_points()[0])
 
 
 class Sector(AnnularSector):
@@ -382,8 +382,8 @@ class Annulus(Circle):
         outer_circle = Circle(radius=self.outer_radius)
         inner_circle = Circle(radius=self.inner_radius)
         inner_circle.reverse_points()
-        self.append_points(outer_circle.points)
-        self.append_points(inner_circle.points)
+        self.append_points(outer_circle.get_points())
+        self.append_points(inner_circle.get_points())
         self.shift(self.arc_center)
 
 
@@ -534,10 +534,10 @@ class DashedLine(Line):
             return Line.get_end(self)
 
     def get_first_handle(self):
-        return self.submobjects[0].points[1]
+        return self.submobjects[0].get_points()[1]
 
     def get_last_handle(self):
-        return self.submobjects[-1].points[-2]
+        return self.submobjects[-1].get_points()[-2]
 
 
 class TangentLine(Line):
@@ -625,7 +625,7 @@ class Arrow(Line):
         # Tip
         self.add_line_to(tip_width * UP / 2)
         self.add_line_to(tip_length * LEFT)
-        self.tip_index = len(self.points) - 1
+        self.tip_index = len(self.get_points()) - 1
         self.add_line_to(tip_width * DOWN / 2)
         self.add_line_to(points2[0])
         # Close it out
@@ -633,7 +633,7 @@ class Arrow(Line):
         self.add_line_to(points1[0])
 
         if length > 0:
-            self.points *= length / self.get_length()  # Final correction
+            self.scale(length / self.get_length())  # Final correction
 
         self.rotate(angle_of_vector(vect) - self.get_angle())
         self.shift(start - self.get_start())
@@ -645,10 +645,11 @@ class Arrow(Line):
 
     def get_start(self):
         nppc = self.n_points_per_curve
-        return (self.points[0] + self.points[-nppc]) / 2
+        points = self.get_points()
+        return (points[0] + points[-nppc]) / 2
 
     def get_end(self):
-        return self.points[self.tip_index]
+        return self.get_points()[self.tip_index]
 
     def put_start_and_end_on(self, start, end):
         self.set_points_by_ends(start, end, buff=0, path_arc=self.path_arc)
@@ -732,7 +733,7 @@ class Polygon(VMobject):
         # To ensure that we loop through starting with last
         arcs = [arcs[-1], *arcs[:-1]]
         for arc1, arc2 in adjacent_pairs(arcs):
-            self.append_points(arc1.points)
+            self.append_points(arc1.get_points())
             line = Line(arc1.get_end(), arc2.get_start())
             # Make sure anchors are evenly distributed
             len_ratio = line.get_length() / arc1.get_arc_length()
@@ -783,7 +784,7 @@ class ArrowTip(Triangle):
         return self.point_from_proportion(0.5)
 
     def get_tip_point(self):
-        return self.points[0]
+        return self.get_points()[0]
 
     def get_vector(self):
         return self.get_tip_point() - self.get_base()
diff --git a/manimlib/mobject/mobject.py b/manimlib/mobject/mobject.py
index 276efd7d..a81705e7 100644
--- a/manimlib/mobject/mobject.py
+++ b/manimlib/mobject/mobject.py
@@ -63,7 +63,7 @@ class Mobject(object):
         self.family = [self]
 
         self.init_updaters()
-        self.points = np.zeros((0, self.dim))
+        self.init_data()
         self.init_points()
         self.init_colors()
         self.init_shader_data()
@@ -76,9 +76,6 @@ class Mobject(object):
     def __str__(self):
         return self.__class__.__name__
 
-    def reset_points(self):
-        self.points = np.zeros((0, self.dim))
-
     def init_colors(self):
         # For subclasses
         pass
@@ -87,6 +84,37 @@ class Mobject(object):
         # Typically implemented in subclass, unless purposefully left blank
         pass
 
+    # To sort out later
+    def init_data(self):
+        self.data = np.zeros(0, dtype=self.shader_dtype)
+
+    def resize_data(self, new_length):
+        if new_length != len(self.data):
+            self.data = np.resize(self.data, new_length)
+
+    def set_points(self, points):
+        self.resize_data(len(points))
+        self.data["point"] = points
+
+    def get_points(self):
+        return self.data["point"]
+
+    def get_all_point_arrays(self):
+        return [self.data["point"]]
+
+    def get_all_data_arrays(self):
+        return [self.data]
+
+    def get_data_array_attrs(self):
+        return ["data"]
+
+    def clear_points(self):
+        self.resize_data(0)
+
+    def get_num_points(self):
+        return len(self.data)
+    #
+
     # Family matters
     def __getitem__(self, value):
         self_list = self.split()
@@ -102,7 +130,7 @@ class Mobject(object):
         return len(self.split())
 
     def split(self):
-        result = [self] if len(self.points) > 0 else []
+        result = [self] if self.get_num_points() > 0 else []
         return result + self.submobjects
 
     def assemble_family(self):
@@ -117,7 +145,7 @@ class Mobject(object):
         return self.family
 
     def family_members_with_points(self):
-        return [m for m in self.get_family() if m.points.size > 0]
+        return [m for m in self.get_family() if m.has_points()]
 
     def add(self, *mobjects):
         if self in mobjects:
@@ -184,7 +212,8 @@ class Mobject(object):
         copy_mobject = copy.copy(self)
         self.parents = parents
 
-        copy_mobject.points = np.array(self.points)
+        for attr in self.get_data_array_attrs():
+            setattr(copy_mobject, attr, getattr(self, attr).copy())
         copy_mobject.submobjects = []
         copy_mobject.add(*[sm.copy() for sm in self.submobjects])
         copy_mobject.match_updaters(self)
@@ -307,10 +336,6 @@ class Mobject(object):
         return self
 
     # Transforming operations
-    def set_points(self, points):
-        self.points = np.array(points)
-        return self
-
     def apply_to_family(self, func):
         for mob in self.family_members_with_points():
             func(mob)
@@ -318,7 +343,8 @@ class Mobject(object):
     def shift(self, *vectors):
         total_vector = reduce(op.add, vectors)
         for mob in self.get_family():
-            mob.points += total_vector
+            for arr in mob.get_all_point_arrays():
+                arr += total_vector
         return self
 
     def scale(self, scale_factor, **kwargs):
@@ -403,14 +429,14 @@ class Mobject(object):
 
     def wag(self, direction=RIGHT, axis=DOWN, wag_factor=1.0):
         for mob in self.family_members_with_points():
-            alphas = np.dot(mob.points, np.transpose(axis))
+            alphas = np.dot(mob.get_points(), np.transpose(axis))
             alphas -= min(alphas)
             alphas /= max(alphas)
             alphas = alphas**wag_factor
-            mob.points += np.dot(
+            mob.set_points(mob.get_points() + np.dot(
                 alphas.reshape((len(alphas), 1)),
                 np.array(direction).reshape((1, mob.dim))
-            )
+            ))
         return self
 
     def reverse_points(self):
@@ -432,9 +458,8 @@ class Mobject(object):
                 about_edge = ORIGIN
             about_point = self.get_bounding_box_point(about_edge)
         for mob in self.family_members_with_points():
-            mob.points -= about_point
-            mob.points[:] = func(mob.points)
-            mob.points += about_point
+            for arr in mob.get_all_point_arrays():
+                arr[:] = func(arr - about_point) + about_point
         return self
 
     # Positioning methods
@@ -782,22 +807,19 @@ class Mobject(object):
         else:
             return getattr(self, array_attr)
 
-    def get_all_points(self):
+    def get_all_points(self):  # TODO, use get_all_point_arrays?
         if self.submobjects:
             return np.vstack([
-                sm.points for sm in self.get_family()
+                sm.get_points() for sm in self.get_family()
             ])
         else:
-            return self.points
+            return self.get_points()
 
     # Getters
 
     def get_points_defining_boundary(self):
         return self.get_all_points()
 
-    def get_num_points(self):
-        return len(self.points)
-
     def get_bounding_box_point(self, direction):
         result = np.zeros(self.dim)
         bb = self.get_bounding_box()
@@ -896,11 +918,11 @@ class Mobject(object):
 
     def get_start(self):
         self.throw_error_if_no_points()
-        return np.array(self.points[0])
+        return np.array(self.get_points()[0])
 
     def get_end(self):
         self.throw_error_if_no_points()
-        return np.array(self.points[-1])
+        return np.array(self.get_points()[-1])
 
     def get_start_and_end(self):
         return self.get_start(), self.get_end()
@@ -929,7 +951,7 @@ class Mobject(object):
         return z_index_group.get_center()
 
     def has_points(self):
-        return len(self.points) > 0
+        return self.get_num_points() > 0
 
     def has_no_points(self):
         return not self.has_points()
@@ -1166,9 +1188,14 @@ class Mobject(object):
         Turns self into an interpolation between mobject1
         and mobject2.
         """
-        self.points[:] = path_func(mobject1.points, mobject2.points, alpha)
-        self.interpolate_color(mobject1, mobject2, alpha)
-        self.interpolate_light_style(mobject1, mobject2, alpha)
+        mobs = [self, mobject1, mobject2]
+        # for arr, arr1, arr2 in zip(*(m.get_all_data_arrays() for m in mobs)):
+        #     arr[:] = interpolate(arr1, arr2, alpha)
+        # if path_func is not straight_path:
+        for arr, arr1, arr2 in zip(*(m.get_all_point_arrays() for m in mobs)):
+            arr[:] = path_func(arr1, arr2, alpha)
+        # self.interpolate_color(mobject1, mobject2, alpha)
+        self.interpolate_light_style(mobject1, mobject2, alpha)  # TODO, interpolate uniforms instaed
         return self
 
     def interpolate_color(self, mobject1, mobject2, alpha):
@@ -1209,8 +1236,8 @@ class Mobject(object):
         """
         self.align_submobjects(mobject)
         for sm1, sm2 in zip(self.get_family(), mobject.get_family()):
-            sm1.set_points(sm2.points)
-            sm1.interpolate_color(sm1, sm2, 1)
+            for arr1, arr2 in zip(sm1.get_all_data_arrays(), sm2.get_all_data_arrays()):
+                arr1[:] = arr2
         return self
 
     def cleanup_from_animation(self):
@@ -1293,10 +1320,9 @@ class Mobject(object):
 
     # For shader data
     def init_shader_data(self):
-        self.shader_data = np.zeros(len(self.points), dtype=self.shader_dtype)
         self.shader_indices = None
         self.shader_wrapper = ShaderWrapper(
-            vert_data=self.shader_data,
+            vert_data=self.data,
             shader_folder=self.shader_folder,
             texture_paths=self.texture_paths,
             depth_test=self.depth_test,
@@ -1307,18 +1333,18 @@ class Mobject(object):
         self.shader_wrapper.refresh_id()
         return self
 
-    def get_blank_shader_data_array(self, size, name="shader_data"):
-        # If possible, try to populate an existing array, rather
-        # than recreating it each frame
-        arr = getattr(self, name)
-        if arr.size != size:
-            new_arr = np.resize(arr, size)
-            setattr(self, name, new_arr)
-            return new_arr
-        return arr
+    # def get_blank_shader_data_array(self, size, name="data"):
+    #     # If possible, try to populate an existing array, rather
+    #     # than recreating it each frame
+    #     arr = getattr(self, name)
+    #     if arr.size != size:
+    #         new_arr = np.resize(arr, size)
+    #         setattr(self, name, new_arr)
+    #         return new_arr
+    #     return arr
 
     def get_shader_wrapper(self):
-        self.shader_wrapper.vert_data = self.get_shader_data()
+        self.shader_wrapper.vert_data = self.data  # TODO
         self.shader_wrapper.vert_indices = self.get_shader_vert_indices()
         self.shader_wrapper.uniforms = self.get_shader_uniforms()
         self.shader_wrapper.depth_test = self.depth_test
@@ -1348,11 +1374,6 @@ class Mobject(object):
             "shadow": self.shadow,
         }
 
-    def get_shader_data(self):
-        # Typically to be implemented by subclasses
-        # Must return a structured numpy array
-        return self.shader_data
-
     def get_shader_vert_indices(self):
         return self.shader_indices
 
@@ -1390,10 +1411,10 @@ class Point(Mobject):
         return self.artificial_height
 
     def get_location(self):
-        return np.array(self.points[0])
+        return np.array(self.get_points()[0])
 
     def get_bounding_box_point(self, *args, **kwargs):
         return self.get_location()
 
     def set_location(self, new_loc):
-        self.points = np.array(new_loc, ndmin=2, dtype=float)
+        self.set_points(np.array(new_loc, ndmin=2, dtype=float))
diff --git a/manimlib/mobject/svg/svg_mobject.py b/manimlib/mobject/svg/svg_mobject.py
index 69876a8d..09d69f8d 100644
--- a/manimlib/mobject/svg/svg_mobject.py
+++ b/manimlib/mobject/svg/svg_mobject.py
@@ -29,24 +29,25 @@ def check_and_fix_percent_bug(sym):
     # The svg path for percent symbols have a known bug, so this
     # checks if the symbol is (probably) a percentage sign, and
     # splits it so that it's displayed properly.
-    if len(sym.points) not in [315, 324, 468, 483] or len(sym.get_subpaths()) != 4:
+    if len(sym.get_points()) not in [315, 324, 468, 483] or len(sym.get_subpaths()) != 4:
         return
 
     sym = sym.family_members_with_points()[0]
     new_sym = VMobject()
     path_lengths = [len(path) for path in sym.get_subpaths()]
-    if len(sym.points) in [315, 324]:
+    sym_points = sym.get_points()
+    if len(sym_points) in [315, 324]:
         n = sum(path_lengths[:2])
-        p1 = sym.points[:n]
-        p2 = sym.points[n:]
-    elif len(sym.points) in [468, 483]:
+        p1 = sym_points[:n]
+        p2 = sym_points[n:]
+    elif len(sym_points) in [468, 483]:
         p1 = np.vstack([
-            sym.points[:path_lengths[0]],
-            sym.points[-path_lengths[3]:]
+            sym_points[:path_lengths[0]],
+            sym_points[-path_lengths[3]:]
         ])
-        p2 = sym.points[path_lengths[0]:sum(path_lengths[:3])]
-    sym.points = p1
-    new_sym.points = p2
+        p2 = sym_points[path_lengths[0]:sum(path_lengths[:3])]
+    sym.set_points(p1)
+    new_sym.set_points(p2)
     sym.add(new_sym)
     sym.refresh_triangulation()
 
@@ -280,7 +281,7 @@ class SVGMobject(VMobject):
             matrix[:, 1] *= -1
 
             for mob in mobject.family_members_with_points():
-                mob.points = np.dot(mob.points, matrix)
+                mob.set_points(np.dot(mob.get_points(), matrix))
             mobject.shift(x * RIGHT + y * UP)
         except:
             pass
@@ -356,7 +357,7 @@ class VMobjectFromSVGPathstring(VMobject):
         filepath = os.path.join(get_mobject_data_dir(), f"{path_hash}.npy")
 
         if os.path.exists(filepath):
-            self.points = np.load(filepath)
+            self.set_points(np.load(filepath))
         else:
             self.relative_point = np.array(ORIGIN)
             for command, coord_string in self.get_commands_and_coord_strings():
@@ -367,11 +368,11 @@ class VMobjectFromSVGPathstring(VMobject):
                 self.subdivide_sharp_curves()
             if self.should_remove_null_curves:
                 # Get rid of any null curves
-                self.points = self.get_points_without_null_curves()
+                self.set_points(self.get_points_without_null_curves())
             # SVG treats y-coordinate differently
             self.stretch(-1, 1, about_point=ORIGIN)
             # Save to a file for future use
-            np.save(filepath, self.points)
+            np.save(filepath, self.get_points())
         check_and_fix_percent_bug(self)
 
     def get_commands_and_coord_strings(self):
@@ -399,11 +400,11 @@ class VMobjectFromSVGPathstring(VMobject):
                 command = "l"
             if command.islower():
                 leftover_points -= self.relative_point
-                self.relative_point = self.points[-1]
+                self.relative_point = self.get_last_point()
             self.handle_command(command, leftover_points)
         else:
             # Command is over, reset for future relative commands
-            self.relative_point = self.points[-1]
+            self.relative_point = self.get_last_point()
 
     def string_to_points(self, command, coord_string):
         numbers = string_to_numbers(coord_string)
diff --git a/manimlib/mobject/types/vectorized_mobject.py b/manimlib/mobject/types/vectorized_mobject.py
index e94b9e7e..002c32c2 100644
--- a/manimlib/mobject/types/vectorized_mobject.py
+++ b/manimlib/mobject/types/vectorized_mobject.py
@@ -53,14 +53,8 @@ class VMobject(Mobject):
         "long_lines": False,
         # For shaders
         "stroke_shader_folder": "quadratic_bezier_stroke",
-        # "stroke_vert_shader_file": "quadratic_bezier_stroke_vert.glsl",
-        # "stroke_geom_shader_file": "quadratic_bezier_stroke_geom.glsl",
-        # "stroke_frag_shader_file": "quadratic_bezier_stroke_frag.glsl",
         "fill_shader_folder": "quadratic_bezier_fill",
-        # "fill_vert_shader_file": "quadratic_bezier_fill_vert.glsl",
-        # "fill_geom_shader_file": "quadratic_bezier_fill_geom.glsl",
-        # "fill_frag_shader_file": "quadratic_bezier_fill_frag.glsl",
-        # Could also be Bevel, Miter, Round
+        # Could also be "bevel", "miter", "round"
         "joint_type": "auto",
         "flat_stroke": True,
         "render_primitive": moderngl.TRIANGLES,
@@ -69,7 +63,6 @@ class VMobject(Mobject):
             ('point', np.float32, (3,)),
             ('unit_normal', np.float32, (3,)),
             ('color', np.float32, (4,)),
-            # ('fill_all', np.float32, (1,)),
             ('vert_index', np.float32, (1,)),
         ],
         "stroke_dtype": [
@@ -92,6 +85,56 @@ class VMobject(Mobject):
     def get_group_class(self):
         return VGroup
 
+    # To sort out later
+    def init_data(self):
+        self.fill_data = np.zeros(0, dtype=self.fill_dtype)
+        self.stroke_data = np.zeros(0, dtype=self.stroke_dtype)
+
+    def resize_data(self, new_length):
+        self.stroke_data = np.resize(self.stroke_data, new_length)
+        self.fill_data = np.resize(self.fill_data, new_length)
+        self.fill_data["vert_index"][:, 0] = range(new_length)
+
+    def set_points(self, points):
+        if len(points) != len(self.stroke_data):
+            self.resize_data(len(points))
+
+        nppc = self.n_points_per_curve
+        self.stroke_data["point"] = points
+        self.stroke_data["prev_point"][:nppc] = points[-nppc:]
+        self.stroke_data["prev_point"][nppc:] = points[:-nppc]
+        self.stroke_data["next_point"][:-nppc] = points[nppc:]
+        self.stroke_data["next_point"][-nppc:] = points[:nppc]
+
+        self.fill_data["point"] = points
+
+        # # TODO, only do conditionally
+        # unit_normal = self.get_unit_normal()
+        # self.stroke_data["unit_normal"] = unit_normal
+        # self.fill_data["unit_normal"] = unit_normal
+
+        # self.refresh_triangulation()
+
+    def get_points(self):
+        return self.stroke_data["point"]
+
+    def get_all_point_arrays(self):
+        return [
+            self.fill_data["point"],
+            self.stroke_data["point"],
+            self.stroke_data["prev_point"],
+            self.stroke_data["next_point"],
+        ]
+
+    def get_all_data_arrays(self):
+        return [self.fill_data, self.stroke_data]
+
+    def get_data_array_attrs(self):
+        return ["fill_data", "stroke_data"]
+
+    def get_num_points(self):
+        return len(self.stroke_data)
+
     # Colors
     def init_colors(self):
         self.fill_rgbas = np.zeros((1, 4))
@@ -379,36 +422,25 @@ class VMobject(Mobject):
         return result
 
     # Points
-    def set_points(self, points):
-        super().set_points(points)
-        self.refresh_triangulation()
-        return self
-
-    def get_points(self):
-        # TODO, shouldn't points always be a numpy array anyway?
-        return np.array(self.points)
-
     def set_anchors_and_handles(self, anchors1, handles, anchors2):
         assert(len(anchors1) == len(handles) == len(anchors2))
         nppc = self.n_points_per_curve
-        self.points = np.zeros((nppc * len(anchors1), self.dim))
+        new_points = np.zeros((nppc * len(anchors1), self.dim))
         arrays = [anchors1, handles, anchors2]
         for index, array in enumerate(arrays):
-            self.points[index::nppc] = array
+            new_points[index::nppc] = array
+        self.set_points(new_points)
         return self
 
-    def clear_points(self):
-        self.points = np.zeros((0, self.dim))
-
     def append_points(self, new_points):
         # TODO, check that number new points is a multiple of 4?
-        # or else that if len(self.points) % 4 == 1, then
+        # or else that if self.get_num_points() % 4 == 1, then
         # len(new_points) % 4 == 3?
-        self.points = np.vstack([self.points, new_points])
+        self.set_points(np.vstack([self.get_points(), new_points]))
         return self
 
     def start_new_path(self, point):
-        assert(len(self.points) % self.n_points_per_curve == 0)
+        assert(self.get_num_points() % self.n_points_per_curve == 0)
         self.append_points([point])
         return self
 
@@ -422,7 +454,7 @@ class VMobject(Mobject):
         """
         self.throw_error_if_no_points()
         quadratic_approx = get_quadratic_approximation_of_cubic(
-            self.points[-1], handle1, handle2, anchor
+            self.get_last_point(), handle1, handle2, anchor
         )
         if self.has_new_path_started():
             self.append_points(quadratic_approx[1:])
@@ -434,10 +466,10 @@ class VMobject(Mobject):
         if self.has_new_path_started():
             self.append_points([handle, anchor])
         else:
-            self.append_points([self.points[-1], handle, anchor])
+            self.append_points([self.get_last_point(), handle, anchor])
 
     def add_line_to(self, point):
-        end = self.points[-1]
+        end = self.get_points()[-1]
         alphas = np.linspace(0, 1, self.n_points_per_curve)
         if self.long_lines:
             halfway = interpolate(end, point, 0.5)
@@ -473,13 +505,14 @@ class VMobject(Mobject):
         self.add_cubic_bezier_curve_to(new_handle, handle, point)
 
     def has_new_path_started(self):
-        return len(self.points) % self.n_points_per_curve == 1
+        return self.get_num_points() % self.n_points_per_curve == 1
 
     def get_last_point(self):
-        return self.points[-1]
+        return self.get_points()[-1]
 
     def get_reflection_of_last_handle(self):
-        return 2 * self.points[-1] - self.points[-2]
+        points = self.get_points()
+        return 2 * points[-1] - points[-2]
 
     def close_path(self):
         if not self.is_closed():
@@ -487,7 +520,7 @@ class VMobject(Mobject):
 
     def is_closed(self):
         return self.consider_points_equals(
-            self.points[0], self.points[-1]
+            self.get_points()[0], self.get_points()[-1]
         )
 
     def subdivide_sharp_curves(self, angle_threshold=30 * DEGREES, family=True):
@@ -509,7 +542,7 @@ class VMobject(Mobject):
                     ])
                 else:
                     new_points.append(tup)
-            vmob.points = np.vstack(new_points)
+            vmob.set_points(np.vstack(new_points))
         return self
 
     def add_points_as_corners(self, points):
@@ -566,12 +599,12 @@ class VMobject(Mobject):
         self.append_points(points)
 
     def append_vectorized_mobject(self, vectorized_mobject):
-        new_points = list(vectorized_mobject.points)
+        new_points = list(vectorized_mobject.get_points())
 
         if self.has_new_path_started():
             # Remove last point, which is starting
             # a new path
-            self.points = self.points[:-1]
+            self.resize_data(len(self.get_points() - 1))
         self.append_points(new_points)
 
     # TODO, how to be smart about tangents here?
@@ -627,13 +660,13 @@ class VMobject(Mobject):
     def get_nth_curve_points(self, n):
         assert(n < self.get_num_curves())
         nppc = self.n_points_per_curve
-        return self.points[nppc * n:nppc * (n + 1)]
+        return self.get_points()[nppc * n:nppc * (n + 1)]
 
     def get_nth_curve_function(self, n):
         return bezier(self.get_nth_curve_points(n))
 
     def get_num_curves(self):
-        return len(self.points) // self.n_points_per_curve
+        return self.get_num_points() // self.n_points_per_curve
 
     def point_from_proportion(self, alpha):
         num_curves = self.get_num_curves()
@@ -649,21 +682,23 @@ class VMobject(Mobject):
         for any i in range(0, len(anchors1))
         """
         nppc = self.n_points_per_curve
+        points = self.get_points()
         return [
-            self.points[i::nppc]
+            points[i::nppc]
             for i in range(nppc)
         ]
 
     def get_start_anchors(self):
-        return self.points[0::self.n_points_per_curve]
+        return self.get_points()[0::self.n_points_per_curve]
 
     def get_end_anchors(self):
         nppc = self.n_points_per_curve
-        return self.points[nppc - 1::nppc]
+        return self.get_points()[nppc - 1::nppc]
 
     def get_anchors(self):
-        if len(self.points) == 1:
-            return self.points
+        points = self.get_points()
+        if len(points) == 1:
+            return points
         return np.array(list(it.chain(*zip(
             self.get_start_anchors(),
             self.get_end_anchors(),
@@ -671,11 +706,12 @@ class VMobject(Mobject):
 
     def get_points_without_null_curves(self, atol=1e-9):
         nppc = self.n_points_per_curve
+        points = self.get_points()
         distinct_curves = reduce(op.or_, [
-            (abs(self.points[i::nppc] - self.points[0::nppc]) > atol).any(1)
+            (abs(points[i::nppc] - points[0::nppc]) > atol).any(1)
             for i in range(1, nppc)
         ])
-        return self.points[distinct_curves.repeat(nppc)]
+        return points[distinct_curves.repeat(nppc)]
 
     def get_arc_length(self, n_sample_points=None):
         if n_sample_points is None:
@@ -697,8 +733,9 @@ class VMobject(Mobject):
             return np.zeros(3)
 
         nppc = self.n_points_per_curve
-        p0 = self.points[0::nppc]
-        p1 = self.points[nppc - 1::nppc]
+        points = self.get_points()
+        p0 = points[0::nppc]
+        p1 = points[nppc - 1::nppc]
 
         # Each term goes through all edges [(x1, y1, z1), (x2, y2, z2)]
         return 0.5 * np.array([
@@ -711,7 +748,7 @@ class VMobject(Mobject):
         if self.unit_normal_locked:
             return self.saved_unit_normal
 
-        if len(self.points) < 3:
+        if self.get_num_points() < 3:
             return OUT
 
         area_vect = self.get_area_vector()
@@ -719,9 +756,10 @@ class VMobject(Mobject):
         if area > 0:
             return area_vect / area
         else:
+            points = self.get_points()
             return get_unit_normal(
-                self.points[1] - self.points[0],
-                self.points[2] - self.points[1],
+                points[1] - points[0],
+                points[2] - points[1],
             )
 
     def lock_unit_normal(self, family=True):
@@ -747,7 +785,7 @@ class VMobject(Mobject):
     # Alignment
     def align_points(self, vmobject):
         self.align_rgbas(vmobject)
-        if len(self.points) == len(vmobject.points):
+        if self.get_num_points() == len(vmobject.get_points()):
             return
 
         for mob in self, vmobject:
@@ -758,7 +796,7 @@ class VMobject(Mobject):
             # If there's only one point, turn it into
             # a null curve
             if mob.has_new_path_started():
-                mob.add_line_to(mob.points[0])
+                mob.add_line_to(mob.get_points()[0])
 
         # Figure out what the subpaths are, and align
         subpaths1 = self.get_subpaths()
@@ -862,7 +900,7 @@ class VMobject(Mobject):
 
     def pointwise_become_partial(self, vmobject, a, b):
         assert(isinstance(vmobject, VMobject))
-        self.points[:] = vmobject.points[:]
+        self.set_points(vmobject.get_points())
         if a <= 0 and b >= 1:
             return self
         num_curves = self.get_num_curves()
@@ -880,23 +918,26 @@ class VMobject(Mobject):
         i3 = nppc * upper_index
         i4 = nppc * (upper_index + 1)
 
+        points = self.get_points()
+        vm_points = vmobject.get_points()
         if num_curves == 0:
-            self.points[:] = 0
+            points[:] = 0
             return self
         if lower_index == upper_index:
-            tup = partial_quadratic_bezier_points(vmobject.points[i1:i2], lower_residue, upper_residue)
-            self.points[:i1] = tup[0]
-            self.points[i1:i4] = tup
-            self.points[i4:] = tup[2]
-            self.points[nppc:] = self.points[nppc - 1]
+            tup = partial_quadratic_bezier_points(vm_points[i1:i2], lower_residue, upper_residue)
+            points[:i1] = tup[0]
+            points[i1:i4] = tup
+            points[i4:] = tup[2]
+            points[nppc:] = points[nppc - 1]
         else:
-            low_tup = partial_quadratic_bezier_points(vmobject.points[i1:i2], lower_residue, 1)
-            high_tup = partial_quadratic_bezier_points(vmobject.points[i3:i4], 0, upper_residue)
-            self.points[0:i1] = low_tup[0]
-            self.points[i1:i2] = low_tup
+            low_tup = partial_quadratic_bezier_points(vm_points[i1:i2], lower_residue, 1)
+            high_tup = partial_quadratic_bezier_points(vm_points[i3:i4], 0, upper_residue)
+            points[0:i1] = low_tup[0]
+            points[i1:i2] = low_tup
             # Keep points i2:i3 as they are
-            self.points[i3:i4] = high_tup
-            self.points[i4:] = high_tup[2]
+            points[i3:i4] = high_tup
+            points[i4:] = high_tup[2]
+        self.set_points(points)
         return self
 
     def get_subcurve(self, a, b):
@@ -912,8 +953,6 @@ class VMobject(Mobject):
 
     # For shaders
     def init_shader_data(self):
-        self.fill_data = np.zeros(len(self.points), dtype=self.fill_dtype)
-        self.stroke_data = np.zeros(len(self.points), dtype=self.stroke_dtype)
         self.fill_shader_wrapper = ShaderWrapper(
             vert_data=self.fill_data,
             vert_indices=np.zeros(0, dtype='i4'),
@@ -980,6 +1019,7 @@ class VMobject(Mobject):
         return result
 
     def get_stroke_shader_data(self):
+        # TODO, make even simpler after fixing colors
         rgbas = self.get_stroke_rgbas()
         if len(rgbas) > 1:
             rgbas = self.stretched_style_array_matching_points(rgbas)
@@ -988,16 +1028,7 @@ class VMobject(Mobject):
         if len(stroke_width) > 1:
             stroke_width = self.stretched_style_array_matching_points(stroke_width)
 
-        points = self.points
-        nppc = self.n_points_per_curve
-
-        data = self.get_blank_shader_data_array(len(points), "stroke_data")
-        data["point"] = points
-        data["prev_point"][:nppc] = points[-nppc:]
-        data["prev_point"][nppc:] = points[:-nppc]
-        data["next_point"][:-nppc] = points[nppc:]
-        data["next_point"][-nppc:] = points[:nppc]
-        data["unit_normal"] = self.get_unit_normal()
+        data = self.stroke_data
         data["stroke_width"][:, 0] = stroke_width
         data["color"] = rgbas
         return data
@@ -1033,12 +1064,14 @@ class VMobject(Mobject):
         if self.triangulation_locked:
             return self.saved_triangulation
 
-        if len(self.points) <= 1:
+        points = self.get_points()
+
+        if len(points) <= 1:
             return np.zeros(0, dtype='i4')
 
         # Rotate points such that unit normal vector is OUT
         # TODO, 99% of the time this does nothing.  Do a check for that?
-        points = np.dot(self.points, z_to_vector(normal_vector))
+        points = np.dot(points, z_to_vector(normal_vector))
         indices = np.arange(len(points), dtype=int)
 
         b0s = points[0::3]
@@ -1074,18 +1107,10 @@ class VMobject(Mobject):
         return tri_indices
 
     def get_fill_shader_data(self):
-        points = self.points
-        n_points = len(points)
-        unit_normal = self.get_unit_normal()
-
-        # TODO, best way to enable multiple colors?
+        # TODO, make simpler
         rgbas = self.get_fill_rgbas()[:1]
-
-        data = self.get_blank_shader_data_array(n_points, "fill_data")
-        data["point"] = points
-        data["unit_normal"] = unit_normal
+        data = self.fill_data
         data["color"] = rgbas
-        data["vert_index"][:, 0] = range(n_points)
         return data
 
     def get_fill_shader_vert_indices(self):
diff --git a/manimlib/utils/bezier.py b/manimlib/utils/bezier.py
index 03f6f668..50feefa4 100644
--- a/manimlib/utils/bezier.py
+++ b/manimlib/utils/bezier.py
@@ -64,7 +64,14 @@ def partial_quadratic_bezier_points(points, a, b):
 # Linear interpolation variants
 
 def interpolate(start, end, alpha):
-    return (1 - alpha) * start + alpha * end
+    try:
+        return (1 - alpha) * start + alpha * end
+    except TypeError:
+        print(type(start), start.dtype)
+        print(type(end), start.dtype)
+        print(alpha)
+        import sys
+        sys.exit(2)
 
 
 def set_array_by_interpolation(arr, arr1, arr2, alpha):
diff --git a/manimlib/utils/iterables.py b/manimlib/utils/iterables.py
index ac90e562..4118642c 100644
--- a/manimlib/utils/iterables.py
+++ b/manimlib/utils/iterables.py
@@ -81,6 +81,7 @@ def listify(obj):
 
 
 def stretch_array_to_length(nparray, length):
+    # TODO, rename to "resize"?
     indices = np.arange(length) * len(nparray) // length
     return nparray[indices]
 
diff --git a/temp_todo.md b/temp_todo.md
new file mode 100644
index 00000000..b8487426
--- /dev/null
+++ b/temp_todo.md
@@ -0,0 +1,8 @@
+- Write init_data
+- Write set_points
+- Write get_all_points_arrays
+- Write get_primary_points_array
+- Write stretch data array method
+- Replace all .points with .data["point"]
+- Change setting/getting colors based on rgba arrays
+- 
\ No newline at end of file