diff --git a/meshmode/mesh/visualization.py b/meshmode/mesh/visualization.py
new file mode 100644
index 0000000000000000000000000000000000000000..7cacd1234fa68e4c794dfe4536bf804a74fffff4
--- /dev/null
+++ b/meshmode/mesh/visualization.py
@@ -0,0 +1,73 @@
+from __future__ import division
+from __future__ import absolute_import
+from six.moves import range
+
+__copyright__ = "Copyright (C) 2014 Andreas Kloeckner"
+
+__license__ = """
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in
+all copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+THE SOFTWARE.
+"""
+
+import numpy as np
+
+
+# {{{ draw_2d_mesh
+
+def draw_2d_mesh(mesh, draw_numbers=True, **kwargs):
+    assert mesh.ambient_dim == 2
+
+    import matplotlib.pyplot as pt
+    import matplotlib.patches as mpatches
+    from matplotlib.path import Path
+
+    for igrp, grp in enumerate(mesh.groups):
+        for iel, el in enumerate(grp.vertex_indices):
+            elverts = mesh.vertices[:, el]
+
+            pathdata = [
+                (Path.MOVETO, (elverts[0, 0], elverts[1, 0])),
+                ]
+            for i in range(1, elverts.shape[1]):
+                pathdata.append(
+                    (Path.LINETO, (elverts[0, i], elverts[1, i]))
+                    )
+            pathdata.append(
+                (Path.CLOSEPOLY, (elverts[0, 0], elverts[1, 0])))
+
+            codes, verts = zip(*pathdata)
+            path = Path(verts, codes)
+            patch = mpatches.PathPatch(path, **kwargs)
+            pt.gca().add_patch(patch)
+
+            if draw_numbers:
+                centroid = (np.sum(elverts, axis=1)
+                        / elverts.shape[1])
+
+                if len(mesh.groups) == 1:
+                    el_label = str(iel)
+                else:
+                    el_label = "%d:%d" % (igrp, iel)
+
+                pt.text(centroid[0], centroid[1], el_label, fontsize=17,
+                        ha="center", va="center",
+                        bbox=dict(facecolor='white', alpha=0.5, lw=0))
+
+# }}}
+
+# vim: foldmethod=marker
diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index e3f55821bd554a477cfe44a81a67e9f1db9bb741..34521057e135ae0e9eae634685f1a4191d44d9ce 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -350,20 +350,10 @@ def test_rect_mesh(do_plot=False):
     mesh = generate_regular_rect_mesh()
 
     if do_plot:
-        cl_ctx = cl.create_some_context()
-        queue = cl.CommandQueue(cl_ctx)
-
-        from meshmode.discretization import Discretization
-        from meshmode.discretization.poly_element import \
-                PolynomialWarpAndBlendGroupFactory
-        discr = Discretization(cl_ctx, mesh,
-                PolynomialWarpAndBlendGroupFactory(order=3))
-
-        from meshmode.discretization.visualization import make_visualizer
-        vis = make_visualizer(queue, discr, vis_order=1)
-        vis.write_vtk_file("rect.vtu", [
-            ("f", discr.nodes()[0]),
-            ])
+        from meshmode.mesh.visualization import draw_2d_mesh
+        draw_2d_mesh(mesh, fill=None)
+        import matplotlib.pyplot as pt
+        pt.show()
 
 
 if __name__ == "__main__":