diff --git a/doc/conf.py b/doc/conf.py
index dae53133f996479359e6a2fd5e4149ad603db899..ef68d84b995dcaa539ea683bb059cb39fb89ed40 100644
--- a/doc/conf.py
+++ b/doc/conf.py
@@ -48,7 +48,7 @@ copyright = u"2008, Andreas Klöckner"
 # The short X.Y version.
 def get_version():
     conf = {}
-    execfile("../meshpy/__init__.py", conf)
+    exec(compile(open("../meshpy/__init__.py").read(), "../meshpy/__init__.py", 'exec'), conf)
     return conf["version"]
 
 version = get_version()
diff --git a/examples/airfoil3d.py b/examples/airfoil3d.py
index d5872d25434e9b901e436baced703ad2e1e1b103..d1307184359efe9176dc7fca1f91514651909966 100644
--- a/examples/airfoil3d.py
+++ b/examples/airfoil3d.py
@@ -1,3 +1,6 @@
+from __future__ import absolute_import
+from __future__ import print_function
+from six.moves import zip
 def main():
     import numpy
     #from math import pi, cos, sin
@@ -59,7 +62,7 @@ def main():
     mesh_info.set_holes([(0, 0, 0.5)])
 
     mesh = build(mesh_info)
-    print "%d elements" % len(mesh.elements)
+    print("%d elements" % len(mesh.elements))
     mesh.write_vtk("airfoil3d.vtk")
 
 
diff --git a/examples/box-in-box.py b/examples/box-in-box.py
index 63e61db997c3f019f0901b0b855ae141c252e27a..028fe18860ad0a982e65b999b454933121b521fc 100644
--- a/examples/box-in-box.py
+++ b/examples/box-in-box.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from __future__ import print_function
 def main():
     import numpy
     from meshpy.tet import MeshInfo, build
@@ -33,7 +35,7 @@ def main():
 
     mesh = build(mesh_info, max_volume=0.06,
             volume_constraints=True, attributes=True)
-    print ("%d elements" % len(mesh.elements))
+    print(("%d elements" % len(mesh.elements)))
     mesh.write_vtk("box-in-box.vtk")
 
 
diff --git a/examples/demo.py b/examples/demo.py
index 289fddf5a7538565f8e991ffcc529f7c416c770d..d4271c6514bdf976ed30977ef8d4ed2ec53d588e 100644
--- a/examples/demo.py
+++ b/examples/demo.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from __future__ import print_function
 from meshpy.tet import MeshInfo, build
 
 mesh_info = MeshInfo()
@@ -14,11 +16,11 @@ mesh_info.set_facets([
     [3,7,4,0],
     ])
 mesh = build(mesh_info)
-print "Mesh Points:"
+print("Mesh Points:")
 for i, p in enumerate(mesh.points):
-    print i, p
-print "Point numbers in tetrahedra:"
+    print(i, p)
+print("Point numbers in tetrahedra:")
 for i, t in enumerate(mesh.elements):
-    print i, t
+    print(i, t)
 mesh.write_vtk("test.vtk")
 
diff --git a/examples/mesh_ply.py b/examples/mesh_ply.py
index 602431a300ba595626a423c1598dcda82adcdfec..4c785385afd1af2527813aa78f3e1f93d90bd417 100644
--- a/examples/mesh_ply.py
+++ b/examples/mesh_ply.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from __future__ import print_function
 def main():
     from ply import parse_ply
     import sys
@@ -16,7 +18,7 @@ def main():
     builder.set(mi)
     mi.set_holes([builder.center()])
     mesh = build(mi)
-    print "%d elements" % len(mesh.elements)
+    print("%d elements" % len(mesh.elements))
     mesh.write_vtk("out.vtk")
 
 
diff --git a/examples/nico_mesh.py b/examples/nico_mesh.py
index b1d0bd92d9aec97caddef8d7434af8c585ddc99f..72ade6ef4e0e001254ab77ef5d7b8dca6431a5d0 100644
--- a/examples/nico_mesh.py
+++ b/examples/nico_mesh.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from six.moves import range
 def main():
     import meshpy.triangle as triangle
     import math
diff --git a/examples/test_ball.py b/examples/test_ball.py
index 5c564ba1cb9a39f997cbcca6ca8db13534dfa56e..86f6466826a6a540b07a41d37a578d6f48c8eeb4 100644
--- a/examples/test_ball.py
+++ b/examples/test_ball.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from six.moves import range
 def main():
     from math import pi, cos, sin
     from meshpy.tet import MeshInfo, build
diff --git a/examples/test_cylinder.py b/examples/test_cylinder.py
index 2b4040de9a5543e2001d3adf8eaa2cebf5bdc2ce..b3624dd63f5d708e27d8a4aaca88a9b403a36cb5 100644
--- a/examples/test_cylinder.py
+++ b/examples/test_cylinder.py
@@ -19,7 +19,7 @@ def main():
 
     mesh = build(mesh_info, max_volume=0.01)
     mesh.write_vtk("cylinder.vtk")
-    mesh.write_neu(file("cylinder.neu", "w"), {
+    mesh.write_neu(open("cylinder.neu", "w"), {
         1: ("minus_z", 1),
         2: ("outer", 2),
         3: ("plus_z", 3),
diff --git a/examples/test_tet_torus.py b/examples/test_tet_torus.py
index ab45c4d8d49bbdcf421296098617937eab8eac06..51a0ccefb0cc607c9e980c93fe55bbb5f77bd264 100644
--- a/examples/test_tet_torus.py
+++ b/examples/test_tet_torus.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from six.moves import range
 def main():
     from math import pi, cos, sin
     from meshpy.tet import MeshInfo, build
@@ -27,7 +29,7 @@ def main():
     mesh.save_elements("torus_mesh")
     mesh.save_nodes("torus_mesh")
 
-    mesh.write_neu(file("torus.neu", "w"),
+    mesh.write_neu(open("torus.neu", "w"),
             {1: ("pec", 0)})
 
 
diff --git a/examples/test_tri_pml.py b/examples/test_tri_pml.py
index 62d59566b71fafdc9ceca074b4a42bb1750b8a0b..c206181dbbe1bd4404e3251257fb034903a2c25b 100644
--- a/examples/test_tri_pml.py
+++ b/examples/test_tri_pml.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from six.moves import range
 def main():
     import meshpy.triangle as triangle
     import math
diff --git a/examples/test_tri_quadratic.py b/examples/test_tri_quadratic.py
index 7683b2c969bb9f7eec17d0cf954775ae4083a2fc..14fc7b6cdb97fa2979ea0397eaad81519adcf538 100644
--- a/examples/test_tri_quadratic.py
+++ b/examples/test_tri_quadratic.py
@@ -7,11 +7,13 @@ from meshpy.triangle import MeshInfo, build
 # Utility function to create lists of the form [(1,2), (2,3), (3,4),
 #(4,1)], given two numbers 1 and 4
 from itertools import islice, cycle
-loop = lambda a, b: zip(range(a, b), islice(cycle(range(a, b)), 1, None))
+from six.moves import range
+from six.moves import zip
+loop = lambda a, b: list(zip(list(range(a, b)), islice(cycle(list(range(a, b))), 1, None)))
 
 info = MeshInfo()
 info.set_points([(0,0), (1,0), (1,1), (0,1), (2,0), (3,0), (3,1), (2,1)])
-info.set_facets(loop(0,4) + loop(4,8), range(1,9)) # Create 8 facets and apply markers 1-8 on them
+info.set_facets(loop(0,4) + loop(4,8), list(range(1,9))) # Create 8 facets and apply markers 1-8 on them
 info.regions.resize(2)
 info.regions[0] = [0.5, 0.5, 1, 0.1] # Fourth item specifies maximum area of triangles as a region attribute
 info.regions[1] = [2.5, 0.5, 2, 0.1] # Replace 0.1 by a smaller value to produce a finer mesh
diff --git a/examples/test_tri_simple_square.py b/examples/test_tri_simple_square.py
index 6607017f755845d9097622f2a87e3bdeb78a2eb8..cbb2bc3c70e82076e9e3f83faa8266b24657a775 100644
--- a/examples/test_tri_simple_square.py
+++ b/examples/test_tri_simple_square.py
@@ -1,3 +1,6 @@
+from __future__ import absolute_import
+from __future__ import print_function
+from six.moves import range
 def main():
     import meshpy.triangle as triangle
 
@@ -16,7 +19,7 @@ def main():
 
     mesh = triangle.build(info, max_volume=1e-3, min_angle=25)
 
-    print "A"
+    print("A")
     triangle.write_gnuplot_mesh("triangles.dat", mesh)
 
 if __name__ == "__main__":
diff --git a/examples/test_triangle.py b/examples/test_triangle.py
index 33b25b02f050457f8e8a8650fc161cc91a534e36..d88e7c81675bf53268fc9a96e1c75f0948370cec 100644
--- a/examples/test_triangle.py
+++ b/examples/test_triangle.py
@@ -1,8 +1,10 @@
 from __future__ import division
+from __future__ import absolute_import
 
 import meshpy.triangle as triangle
 import numpy as np
 import numpy.linalg as la
+from six.moves import range
 
 
 def round_trip_connect(start, end):
diff --git a/examples/test_triangle_refine.py b/examples/test_triangle_refine.py
index 2f845f86b72153d17f03922059cc02ba6b5692f5..16183f892f6a19489adff6af245b152dcf098276 100644
--- a/examples/test_triangle_refine.py
+++ b/examples/test_triangle_refine.py
@@ -1,3 +1,6 @@
+from __future__ import absolute_import
+from __future__ import print_function
+from six.moves import range
 def main():
     import meshpy.triangle as triangle
     import math
@@ -36,7 +39,7 @@ def main():
             )
 
     triangle.write_gnuplot_mesh("triangles-unrefined.dat", mesh)
-    print len(mesh.elements)
+    print(len(mesh.elements))
 
     mesh.element_volumes.setup()
 
@@ -46,7 +49,7 @@ def main():
         mesh.element_volumes[i] = 1e-8
 
     mesh = triangle.refine(mesh)
-    print len(mesh.elements)
+    print(len(mesh.elements))
 
     triangle.write_gnuplot_mesh("triangles.dat", mesh)
 
diff --git a/examples/write_dolfin.py b/examples/write_dolfin.py
index a390c8edbe746b510fb01f77774442cb075f3e5e..ec8c9317a00a0bc73e28caab4ca338c617c868b8 100644
--- a/examples/write_dolfin.py
+++ b/examples/write_dolfin.py
@@ -1,3 +1,5 @@
+from __future__ import absolute_import
+from __future__ import print_function
 def main():
     import meshpy.triangle as triangle
 
@@ -7,32 +9,32 @@ def main():
 
     mesh = triangle.build(info, max_volume=1e-3, min_angle=25)
 
-    print """
+    print("""
         <?xml version="1.0" encoding="UTF-8"?>
 
         <dolfin xmlns:dolfin="http://www.fenics.org/dolfin/">
           <mesh celltype="triangle" dim="2">
             <vertices size="%d">
-        """ % len(mesh.points)
+        """ % len(mesh.points))
 
     for i, pt in enumerate(mesh.points):
-      print '<vertex index="%d" x="%g" y="%g"/>' % (
-              i, pt[0], pt[1])
+      print('<vertex index="%d" x="%g" y="%g"/>' % (
+              i, pt[0], pt[1]))
 
-    print """
+    print("""
         </vertices>
         <cells size="%d">
-        """ % len(mesh.elements)
+        """ % len(mesh.elements))
 
     for i, element in enumerate(mesh.elements):
-      print '<triangle index="%d" v0="%d" v1="%d" v2="%d"/>' % (
-              i, element[0], element[1], element[2])
+      print('<triangle index="%d" v0="%d" v1="%d" v2="%d"/>' % (
+              i, element[0], element[1], element[2]))
 
-    print """
+    print("""
             </cells>
           </mesh>
         </dolfin>
-        """
+        """)
 
 if __name__ == "__main__":
     main()
diff --git a/meshpy/common.py b/meshpy/common.py
index e620021e3eac8829ca7675d76769f9a83ab2425e..5715a6fd35f074d7a23e0763e17468d614ccfb05 100644
--- a/meshpy/common.py
+++ b/meshpy/common.py
@@ -1,3 +1,7 @@
+from __future__ import absolute_import
+from __future__ import print_function
+from six.moves import range
+from six.moves import zip
 class _Table:
     def __init__(self):
         self.Rows = []
@@ -101,7 +105,7 @@ class MeshInfoBase:
         outfile.write("PROGRAM: MeshPy VERSION: %s\n" % version)
         outfile.write("%s\n" % datetime.now().ctime())
 
-        bc_markers = bc.keys()
+        bc_markers = list(bc.keys())
         if periodicity:
             periodic_marker, periods = periodicity
             bc_markers.append(periodic_marker)
@@ -149,7 +153,7 @@ class MeshInfoBase:
         outfile.write("ELEMENT GROUP 1.3.0\n")
         # FIXME
         i = 0
-        grp_elements = range(len(self.elements))
+        grp_elements = list(range(len(self.elements)))
         material = 1
         flags = 0
         outfile.write("GROUP:%11d ELEMENTS:%11d MATERIAL:%11s NFLAGS: %11d\n"
@@ -190,7 +194,7 @@ class MeshInfoBase:
                     face2el.setdefault(face, []).append((ti, fi+1))
 
         else:
-            raise ValueError, "invalid number of dimensions (%d)" % dim
+            raise ValueError("invalid number of dimensions (%d)" % dim)
 
         # actually output bc sections
         if not self.faces.allocated:
@@ -257,7 +261,7 @@ class MeshInfoBase:
 
 
 def dump_array(name, array):
-    print "array %s: %d elements, %d values per element" % (name, len(array), array.unit)
+    print("array %s: %d elements, %d values per element" % (name, len(array), array.unit))
 
     if len(array) == 0 or array.unit == 0:
         return
@@ -265,12 +269,12 @@ def dump_array(name, array):
     try:
         array[0]
     except RuntimeError:
-        print "  not allocated"
+        print("  not allocated")
         return
 
     for i, entry in enumerate(array):
         if isinstance(entry, list):
-            print "  %d: %s" % (i, ",".join(str(sub) for sub in entry))
+            print("  %d: %s" % (i, ",".join(str(sub) for sub in entry)))
         else:
-            print "  %d: %s" % (i, entry)
+            print("  %d: %s" % (i, entry))
 
diff --git a/meshpy/geometry.py b/meshpy/geometry.py
index c6e6abb9e6853e607e99dd22fc9c093dbbaf3082..fcf1a3ba3c23d5863029f0502a30b8bb74dd942e 100644
--- a/meshpy/geometry.py
+++ b/meshpy/geometry.py
@@ -1,5 +1,8 @@
 from __future__ import division
+from __future__ import absolute_import
 import numpy as np
+from six.moves import range
+from six.moves import zip
 
 __doc__ = """
 
@@ -290,7 +293,7 @@ def make_circle(r, center=(0, 0), subdivisions=40, marker=Marker.SHELL):
     y = r*np.sin(phi) + cy
 
     return ([np.array(pt) for pt in zip(x, y)],
-            round_trip_connect(range(subdivisions)),
+            round_trip_connect(list(range(subdivisions))),
             None,
             subdivisions*[marker])
 
@@ -414,7 +417,7 @@ def generate_extrusion(rz_points, base_shape, closure=EXT_OPEN,
             p_indices = (first_idx,)
             points.append((0, 0, z))
         else:
-            p_indices = tuple(xrange(first_idx, first_idx+len(base_shape)))
+            p_indices = tuple(range(first_idx, first_idx+len(base_shape)))
             points.extend([(x*r, y*r, z) for (x, y) in base_shape])
 
         rings[ring_idx] = p_indices
diff --git a/meshpy/gmsh_reader.py b/meshpy/gmsh_reader.py
index ac67ca5b21c7f69b054acb3222e17ce0bb2db4d9..b4d0ef2e7da5f52ea5156c948ba754762993ab98 100644
--- a/meshpy/gmsh_reader.py
+++ b/meshpy/gmsh_reader.py
@@ -1,6 +1,9 @@
 """Reader for the GMSH file format."""
 
 from __future__ import division
+from __future__ import absolute_import
+from six.moves import range
+from functools import reduce
 
 __copyright__ = "Copyright (C) 2009 Xueyu Zhu, Andreas Kloeckner"
 
@@ -98,7 +101,7 @@ class LineFeeder:
             return True
 
         try:
-            self.next_line = self.line_iterable.next()
+            self.next_line = next(self.line_iterable)
         except StopIteration:
             return False
         else:
@@ -111,7 +114,7 @@ class LineFeeder:
             return nl.strip()
 
         try:
-            nl = self.line_iterable.next()
+            nl = next(self.line_iterable)
         except StopIteration:
             raise GmshFileFormatError("unexpected end of file")
         else:
diff --git a/meshpy/naca.py b/meshpy/naca.py
index deb72937b5b642d2a48df1d03c6b8bf9c39d52ab..1f55d65d0c64851667ff236efdbc57dbe3e78484 100644
--- a/meshpy/naca.py
+++ b/meshpy/naca.py
@@ -1,6 +1,9 @@
 from __future__ import division
+from __future__ import absolute_import
+from __future__ import print_function
 
 import numpy
+from six.moves import range
 
 
 class FourDigitsSymmetric:
@@ -130,7 +133,7 @@ def get_naca_points(naca_digits, number_of_points=100,
 
     if verbose:
         def explain(*s):
-            print " ".join(str(s_i) for s_i in s)
+            print(" ".join(str(s_i) for s_i in s))
     else:
         def explain(*s):
             pass
@@ -146,7 +149,7 @@ def get_naca_points(naca_digits, number_of_points=100,
 
     raw_abscissae = numpy.linspace(0, 1, number_of_points, endpoint=True)
     abscissae = numpy.empty_like(raw_abscissae)
-    for i in xrange(number_of_points):
+    for i in range(number_of_points):
         abscissae[i] = abscissa_map(raw_abscissae[i])
 
     digits_int = int(naca_digits)
@@ -227,7 +230,7 @@ def get_naca_points(naca_digits, number_of_points=100,
 def write_points(points, filename):
     file = open(filename, "w")
     for pt in points:
-        print >>file, "\t".join(repr(p_comp) for p_comp in pt)
+        print("\t".join(repr(p_comp) for p_comp in pt), file=file)
 
 
 def main():
@@ -263,7 +266,7 @@ def main():
     if options.output is None:
         options.output = "naca-%s.dat" % digits
 
-    print "Output file:", options.output
+    print("Output file:", options.output)
     write_points(points, options.output)
 
 
diff --git a/meshpy/ply.py b/meshpy/ply.py
index e85edfe683d5f61e2d5ae65071a1be40965a7e32..872620c0ccd0ea2f67564ac08578ac76bd1094e7 100644
--- a/meshpy/ply.py
+++ b/meshpy/ply.py
@@ -1,7 +1,10 @@
+from __future__ import absolute_import
+from six.moves import range
+from six.moves import zip
 def parse_int(it):
-    return int(it.next())
+    return int(next(it))
 def parse_float(it):
-    return float(it.next())
+    return float(next(it))
 class ListParser:
     def __init__(self, len_parser, item_parser):
         self.len_parser = len_parser
@@ -11,7 +14,7 @@ class ListParser:
         return [self.item_parser(it) for i in range(self.len_parser(it))]
 
 def make_parser(it):
-    tp = it.next()
+    tp = next(it)
     if tp == "list":
         len_parser = make_parser(it)
         item_parser = make_parser(it)
@@ -49,7 +52,7 @@ def parse_ply(name):
         elif words[0] in ["comment", "created"]:
             i += 1
         else:
-            raise ValueError, "invalid header field"
+            raise ValueError("invalid header field")
     i += 1 # skip end_header
 
     result = {}
@@ -66,7 +69,7 @@ def parse_ply(name):
         pass
 
     for name, line_count, props in data_queue:
-        prop_names, parsers = zip(*props)
+        prop_names, parsers = list(zip(*props))
         result[name] = DataBlock(
                 properties=prop_names, 
                 data=[parse_line(parsers, l) for l in lines[i:i+line_count]])
diff --git a/meshpy/tet.py b/meshpy/tet.py
index 934f803effdabdf5fbf2e27e92aca9af16c850d1..113aed600c8a4651cb2bc0af5af209264c3a9519 100644
--- a/meshpy/tet.py
+++ b/meshpy/tet.py
@@ -1,5 +1,9 @@
+from __future__ import absolute_import
+from __future__ import print_function
 from meshpy.common import MeshInfoBase, dump_array
 import meshpy._tetgen as internals
+import six
+from six.moves import range
 
 
 class MeshInfo(internals.MeshInfo, MeshInfoBase):
@@ -93,10 +97,10 @@ class MeshInfo(internals.MeshInfo, MeshInfoBase):
         for name in ["points"]:
             dump_array(name, getattr(self, name))
         for ifacet, facet in enumerate(self.faces):
-            print "facet %d:" % ifacet
+            print("facet %d:" % ifacet)
             for ipolygon, polygon in enumerate(facet.polygons):
-                print "  polygon %d: vertices [%s]" % \
-                        (ipolygon, ",".join(str(vi) for vi in polygon.vertices))
+                print("  polygon %d: vertices [%s]" % \
+                        (ipolygon, ",".join(str(vi) for vi in polygon.vertices)))
 
     def write_vtk(self, filename):
         import pyvtk
@@ -116,7 +120,7 @@ class MeshInfo(internals.MeshInfo, MeshInfoBase):
     def set_element_constraints(self, element_constraints):
         self.element_volumes.setup()
 
-        for i in xrange(len(self.element_volumes)):
+        for i in range(len(self.element_volumes)):
             if i in element_constraints:
                 self.element_volumes[i] = element_constraints[i]
             else:
@@ -132,7 +136,7 @@ class Options(internals.Options):
         self.parse_switches(switches)
         self.quiet = 1
 
-        for k, v in kwargs.iteritems():
+        for k, v in six.iteritems(kwargs):
             try:
                 getattr(self, k)
             except AttributeError:
diff --git a/meshpy/tools.py b/meshpy/tools.py
index c58ee95fc0a9f9cf5d4d2ed54f2d3ac5f0edb5f5..837384f786ec6c2a2b8c9c5d9d261bd8a31a7a17 100644
--- a/meshpy/tools.py
+++ b/meshpy/tools.py
@@ -1,3 +1,6 @@
+from __future__ import absolute_import
+import six
+from six.moves import range
 def uniform_refine_triangles(points, elements, factor=2):
     new_points = points[:]
     new_elements = []
@@ -103,13 +106,13 @@ def make_swizzle_matrix(spec):
 
     assert set(mapping.keys()) == set(axes), \
             "axis mapping not complete"
-    assert set(axis.lstrip("-") for axis in mapping.itervalues()) == set(axes), \
+    assert set(axis.lstrip("-") for axis in six.itervalues(mapping)) == set(axes), \
             "Axis mapping not onto"
 
     n = len(axes)
     result = numpy.zeros((n, n), dtype=int)
 
-    for imp_axis, final_axis in mapping.iteritems():
+    for imp_axis, final_axis in six.iteritems(mapping):
         imp_axis = axes.index(imp_axis)
 
         sign = 1
diff --git a/meshpy/triangle.py b/meshpy/triangle.py
index 1da0b8570ba23cabac97858952a30d672ccb8bf9..3bf887fc48105d91a0b08fbbd19d7bb83e655c0a 100644
--- a/meshpy/triangle.py
+++ b/meshpy/triangle.py
@@ -1,6 +1,9 @@
 from __future__ import division
+from __future__ import absolute_import
 from meshpy.common import MeshInfoBase, dump_array
 import meshpy._triangle as internals
+from six.moves import range
+from six.moves import zip
 
 
 
@@ -22,12 +25,13 @@ class MeshInfo(internals.MeshInfo, MeshInfoBase):
                self.number_of_element_attributes, \
                [(name, getattr(self, name)) for name in self._constituents]
 
-    def __setstate__(self, (p_attr_count, e_attr_count, state)):
+    def __setstate__(self, xxx_todo_changeme):
+        (p_attr_count, e_attr_count, state) = xxx_todo_changeme
         self.number_of_point_attributes = p_attr_count
         self.number_of_element_attributes = e_attr_count
         for name, array in state:
             if name not in self._constituents:
-                raise RuntimeError, "Unknown constituent during unpickling"
+                raise RuntimeError("Unknown constituent during unpickling")
 
             dest_array = getattr(self, name)
 
@@ -151,7 +155,7 @@ def build(mesh_info, verbose=False, refinement_func=None, attributes=False,
     if volume_constraints:
         opts += "a"
         if max_volume:
-            raise ValueError, "cannot specify both volume_constraints and max_area"
+            raise ValueError("cannot specify both volume_constraints and max_area")
     elif max_volume:
         opts += "a%.20f" % max_volume
 
diff --git a/setup.py b/setup.py
index 89e20c520e3203338e9d5f4482e18014df67ff9a..16af323efe2d31585a9c9c2c42eda7ea3e57b825 100644
--- a/setup.py
+++ b/setup.py
@@ -49,12 +49,6 @@ def main():
     init_filename = "meshpy/__init__.py"
     exec(compile(open(init_filename, "r").read(), init_filename, "exec"), conf)
 
-    try:
-        from distutils.command.build_py import build_py_2to3 as build_py
-    except ImportError:
-        # 2.x
-        from distutils.command.build_py import build_py
-
     import codecs
     setup(name="MeshPy",
           version=conf["version"],
@@ -74,7 +68,12 @@ def main():
               'Natural Language :: English',
               'Programming Language :: C++',
               'Programming Language :: Python',
+              'Programming Language :: Python :: 2.6',
+              'Programming Language :: Python :: 2.7',
               'Programming Language :: Python :: 3',
+              'Programming Language :: Python :: 3.2',
+              'Programming Language :: Python :: 3.3',
+              'Programming Language :: Python :: 3.4',
               'Topic :: Multimedia :: Graphics :: 3D Modeling',
               'Topic :: Scientific/Engineering',
               'Topic :: Scientific/Engineering :: Mathematics',
@@ -115,11 +114,7 @@ def main():
                   extra_compile_args=conf["CXXFLAGS"],
                   extra_link_args=conf["LDFLAGS"],
                   ),
-              ],
-
-          # 2to3 invocation
-          cmdclass={'build_py': build_py},
-          )
+              ])
 
 
 if __name__ == '__main__':
diff --git a/test/test_meshpy.py b/test/test_meshpy.py
index 1a1e8e149e1901d6b31618cba624d0668b269a36..1d7b2864cd82fe26793a27b0c95f9a61d2fe1411 100644
--- a/test/test_meshpy.py
+++ b/test/test_meshpy.py
@@ -1,4 +1,6 @@
 from __future__ import division
+from __future__ import absolute_import
+from six.moves import range
 
 __copyright__ = "Copyright (C) 2013 Andreas Kloeckner"