diff --git a/mathics/builtin/drawing/graphics3d.py b/mathics/builtin/drawing/graphics3d.py
index efbf2c4cbf..a9df10705f 100644
--- a/mathics/builtin/drawing/graphics3d.py
+++ b/mathics/builtin/drawing/graphics3d.py
@@ -4,34 +4,36 @@
 Three-Dimensional Graphics
 """
 
+from mathics.version import __version__  # noqa used in loading to check consistency.
 
+import html
+import json
 import numbers
-from mathics.version import __version__  # noqa used in loading to check consistency.
+
 from mathics.core.expression import (
     Expression,
     from_python,
     system_symbols_dict,
     SymbolList,
 )
+from mathics.core.formatter import lookup_method
+
 from mathics.builtin.base import BoxConstructError, Builtin, InstanceableBuiltin
 from mathics.builtin.graphics import (
     Graphics,
     GraphicsBox,
     PolygonBox,
-    create_pens,
     _Color,
     LineBox,
     PointBox,
     Style,
     RGBColor,
     get_class,
-    asy_number,
     CoordinatesError,
     _GraphicsElements,
 )
 
-import json
-import html
+from mathics.formatter.asy_fns import asy_create_pens, asy_number
 
 
 def coords3D(value):
@@ -443,7 +445,11 @@ def boxes_to_tex(self, leaves=None, **options):
 
         elements._apply_boxscaling(boxscale)
 
-        asy = elements.to_asy()
+        format_fn = lookup_method(elements, "asy")
+        if format_fn is not None:
+            asy = format_fn(elements)
+        else:
+            asy = elements.to_asy()
 
         xmin, xmax, ymin, ymax, zmin, zmax, boxscale = calc_dimensions()
 
@@ -466,11 +472,11 @@ def boxes_to_tex(self, leaves=None, **options):
 
         for i, line in enumerate(boundbox_lines):
             if i in axes_indices:
-                pen = create_pens(
+                pen = asy_create_pens(
                     edge_color=RGBColor(components=(0, 0, 0, 1)), stroke_width=1.5
                 )
             else:
-                pen = create_pens(
+                pen = asy_create_pens(
                     edge_color=RGBColor(components=(0.4, 0.4, 0.4, 1)), stroke_width=1
                 )
 
@@ -483,7 +489,7 @@ def boxes_to_tex(self, leaves=None, **options):
 
         # Draw axes ticks
         ticklength = 0.05 * max([xmax - xmin, ymax - ymin, zmax - zmin])
-        pen = create_pens(
+        pen = asy_create_pens(
             edge_color=RGBColor(components=(0, 0, 0, 1)), stroke_width=1.2
         )
         for xi in axes_indices:
@@ -778,9 +784,6 @@ def __init__(self, content, evaluation, neg_y=False):
     def extent(self, completely_visible_only=False):
         return total_extent_3d([element.extent() for element in self.elements])
 
-    def to_asy(self):
-        return "\n".join([element.to_asy() for element in self.elements])
-
     def _apply_boxscaling(self, boxscale):
         for element in self.elements:
             element._apply_boxscaling(boxscale)
@@ -821,22 +824,6 @@ def to_json(self):
             )
         return data
 
-    def to_asy(self):
-        face_color = self.face_color
-
-        # Tempoary bug fix: default Point color should be black not white
-        if list(face_color.to_rgba()[:3]) == [1, 1, 1]:
-            face_color = RGBColor(components=(0, 0, 0, face_color.to_rgba()[3]))
-
-        pen = create_pens(face_color=face_color, is_face_element=False)
-
-        return "".join(
-            "path3 g={0}--cycle;dot(g, {1});".format(
-                "--".join("(%.5g,%.5g,%.5g)" % coords.pos()[0] for coords in line), pen
-            )
-            for line in self.lines
-        )
-
     def extent(self):
         result = []
         for line in self.lines:
@@ -871,18 +858,6 @@ def to_json(self):
             )
         return data
 
-    def to_asy(self):
-        # l = self.style.get_line_width(face_element=False)
-        pen = create_pens(edge_color=self.edge_color, stroke_width=1)
-
-        return "".join(
-            "draw({0}, {1});".format(
-                "--".join("({0},{1},{2})".format(*coords.pos()[0]) for coords in line),
-                pen,
-            )
-            for line in self.lines
-        )
-
     def extent(self):
         result = []
         for line in self.lines:
@@ -930,29 +905,6 @@ def to_json(self):
             )
         return data
 
-    def to_asy(self):
-        l = self.style.get_line_width(face_element=True)
-        if self.vertex_colors is None:
-            face_color = self.face_color
-        else:
-            face_color = None
-        pen = create_pens(
-            edge_color=self.edge_color,
-            face_color=face_color,
-            stroke_width=l,
-            is_face_element=True,
-        )
-
-        asy = ""
-        for line in self.lines:
-            asy += (
-                "path3 g="
-                + "--".join(["(%.5g,%.5g,%.5g)" % coords.pos()[0] for coords in line])
-                + "--cycle;"
-            )
-            asy += "draw(surface(g), %s);" % (pen)
-        return asy
-
     def extent(self):
         result = []
         for line in self.lines:
@@ -1146,21 +1098,6 @@ def init(self, graphics, style, item):
         self.points = [Coords3D(graphics, pos=point) for point in points]
         self.radius = item.leaves[1].to_python()
 
-    def to_asy(self):
-        # l = self.style.get_line_width(face_element=True)
-
-        if self.face_color is None:
-            face_color = (1, 1, 1)
-        else:
-            face_color = self.face_color.to_js()
-
-        return "".join(
-            "draw(surface(sphere({0}, {1})), rgb({2},{3},{4}));".format(
-                tuple(coord.pos()[0]), self.radius, *face_color[:3]
-            )
-            for coord in self.points
-        )
-
     def to_json(self):
         face_color = self.face_color
         if face_color is not None:
diff --git a/mathics/builtin/graphics.py b/mathics/builtin/graphics.py
index cef531b68f..d1911f70c4 100644
--- a/mathics/builtin/graphics.py
+++ b/mathics/builtin/graphics.py
@@ -8,7 +8,6 @@
 
 from math import floor, ceil, log10, sin, cos, pi, sqrt, atan2, degrees, radians, exp
 import base64
-from itertools import chain
 
 from mathics.version import __version__  # noqa used in loading to check consistency.
 from mathics.builtin.base import (
@@ -32,10 +31,12 @@
     system_symbols_dict,
     from_python,
 )
-from mathics.core.formatter import lookup_method
 
 from mathics.builtin.drawing.colors import convert as convert_color
+from mathics.core.formatter import lookup_method
 from mathics.core.numbers import machine_epsilon
+from mathics.formatter.asy_fns import asy_bezier, asy_color, asy_number
+
 
 GRAPHICS_OPTIONS = {
     "AspectRatio": "Automatic",
@@ -50,7 +51,6 @@
     "$OptionSyntax": "Ignore",
 }
 
-
 class CoordinatesError(BoxConstructError):
     pass
 
@@ -144,10 +144,6 @@ def create_css(
     return "; ".join(css)
 
 
-def asy_number(value):
-    return "%.5g" % value
-
-
 def _to_float(x):
     x = x.round_to_float()
     if x is None:
@@ -155,29 +151,6 @@ def _to_float(x):
     return x
 
 
-def create_pens(
-    edge_color=None, face_color=None, stroke_width=None, is_face_element=False
-):
-    result = []
-    if face_color is not None:
-        brush, opacity = face_color.to_asy()
-        if opacity != 1:
-            brush += "+opacity(%s)" % asy_number(opacity)
-        result.append(brush)
-    elif is_face_element:
-        result.append("nullpen")
-    if edge_color is not None:
-        pen, opacity = edge_color.to_asy()
-        if opacity != 1:
-            pen += "+opacity(%s)" % asy_number(opacity)
-        if stroke_width is not None:
-            pen += "+linewidth(%s)" % asy_number(stroke_width)
-        result.append(pen)
-    elif is_face_element:
-        result.append("nullpen")
-    return ", ".join(result)
-
-
 def _data_and_options(leaves, defined_options):
     data = []
     options = defined_options.copy()
@@ -321,55 +294,15 @@ def _extract_graphics(graphics, format, evaluation):
 
     # generate code for svg or asy.
 
-    if format == "asy":
-        code = "\n".join(element.to_asy() for element in elements.elements)
-    elif format == "svg":
-        format_fn = lookup_method(elements, "svg")
-        if format_fn is not None:
-            code = format_fn(elements)
-        else:
-            code = elements.to_svg()
+    if format in ("asy", "svg"):
+        format_fn = lookup_method(elements, format)
+        code = format_fn(elements)
     else:
         raise NotImplementedError
 
     return xmin, xmax, ymin, ymax, ox, oy, ex, ey, code
 
 
-class _ASYTransform:
-    _template = """
-    add(%s * (new picture() {
-        picture saved = currentpicture;
-        picture transformed = new picture;
-        currentpicture = transformed;
-        %s
-        currentpicture = saved;
-        return transformed;
-    })());
-    """
-
-    def __init__(self):
-        self.transforms = []
-
-    def matrix(self, a, b, c, d, e, f):
-        # a c e
-        # b d f
-        # 0 0 1
-        # see http://asymptote.sourceforge.net/doc/Transforms.html#Transforms
-        self.transforms.append("(%f, %f, %f, %f, %f, %f)" % (e, f, a, c, b, d))
-
-    def translate(self, x, y):
-        self.transforms.append("shift(%f, %f)" % (x, y))
-
-    def scale(self, x, y):
-        self.transforms.append("scale(%f, %f)" % (x, y))
-
-    def rotate(self, x):
-        self.transforms.append("rotate(%f)" % x)
-
-    def apply(self, asy):
-        return self._template % (" * ".join(self.transforms), asy)
-
-
 class Show(Builtin):
     """
     <dl>
@@ -596,15 +529,6 @@ def to_css(self):
             alpha,
         )
 
-    def to_asy(self):
-        rgba = self.to_rgba()
-        alpha = rgba[3] if len(rgba) > 3 else 1.0
-        return (
-            r"rgb(%s, %s, %s)"
-            % (asy_number(rgba[0]), asy_number(rgba[1]), asy_number(rgba[2])),
-            alpha,
-        )
-
     def to_js(self):
         return self.to_rgba()
 
@@ -845,12 +769,18 @@ class ColorDistance(Builtin):
 
     options = {"DistanceFunction": "Automatic"}
 
+    requires = ("numpy",)
+
     messages = {
         "invdist": "`1` is not Automatic or a valid distance specification.",
         "invarg": "`1` and `2` should be two colors or a color and a lists of colors or "
         + "two lists of colors of the same length.",
     }
 
+    # If numpy is not installed, 100 * c1.to_color_space returns
+    # a list of 100 x 3 elements, instead of doing elementwise multiplication
+    requires = ("numpy",)
+
     # the docs say LABColor's colorspace corresponds to the CIE 1976 L^* a^* b^* color space
     # with {l,a,b}={L^*,a^*,b^*}/100. Corrections factors are put accordingly.
 
@@ -887,13 +817,6 @@ class ColorDistance(Builtin):
     def apply(self, c1, c2, evaluation, options):
         "ColorDistance[c1_, c2_, OptionsPattern[ColorDistance]]"
 
-        # If numpy is not installed, 100 * c1.to_color_space returns
-        # a list of 100 x 3 elements, instead of doing elementwise multiplication
-        try:
-            import numpy as np
-        except:
-            raise RuntimeError("NumPy needs to be installed for ColorDistance")
-
         distance_function = options.get("System`DistanceFunction")
         compute = None
         if isinstance(distance_function, String):
@@ -1219,24 +1142,6 @@ def extent(self):
             )
         return result
 
-    def to_asy(self):
-        l = self.style.get_line_width(face_element=True)
-        x1, y1 = self.p1.pos()
-        x2, y2 = self.p2.pos()
-        pens = create_pens(self.edge_color, self.face_color, l, is_face_element=True)
-        x1, x2, y1, y2 = asy_number(x1), asy_number(x2), asy_number(y1), asy_number(y2)
-        return "filldraw((%s,%s)--(%s,%s)--(%s,%s)--(%s,%s)--cycle, %s);" % (
-            x1,
-            y1,
-            x2,
-            y1,
-            x2,
-            y2,
-            x1,
-            y2,
-            pens,
-        )
-
 
 class _RoundBox(_GraphicsElement):
     face_element = None
@@ -1270,28 +1175,6 @@ def extent(self):
         ry += l
         return [(x - rx, y - ry), (x - rx, y + ry), (x + rx, y - ry), (x + rx, y + ry)]
 
-    def to_asy(self):
-        x, y = self.c.pos()
-        rx, ry = self.r.pos()
-        rx -= x
-        ry -= y
-        l = self.style.get_line_width(face_element=self.face_element)
-        pen = create_pens(
-            edge_color=self.edge_color,
-            face_color=self.face_color,
-            stroke_width=l,
-            is_face_element=self.face_element,
-        )
-        cmd = "filldraw" if self.face_element else "draw"
-        return "%s(ellipse((%s,%s),%s,%s), %s);" % (
-            cmd,
-            asy_number(x),
-            asy_number(y),
-            asy_number(rx),
-            asy_number(ry),
-            pen,
-        )
-
 
 class _ArcBox(_RoundBox):
     def init(self, graphics, style, item):
@@ -1342,6 +1225,7 @@ def _arc_params(self):
 
         return x, y, abs(rx), abs(ry), sx, sy, ex, ey, large_arc
 
+
 class DiskBox(_ArcBox):
     face_element = True
 
@@ -1423,16 +1307,6 @@ def init(self, graphics, style, item=None):
         else:
             raise BoxConstructError
 
-    def to_asy(self):
-        pen = create_pens(face_color=self.face_color, is_face_element=False)
-
-        asy = ""
-        for line in self.lines:
-            for coords in line:
-                asy += "dot(%s, %s);" % (coords.pos(), pen)
-
-        return asy
-
 
 class Line(Builtin):
     """
@@ -1467,15 +1341,6 @@ def init(self, graphics, style, item=None, lines=None):
         else:
             raise BoxConstructError
 
-    def to_asy(self):
-        l = self.style.get_line_width(face_element=False)
-        pen = create_pens(edge_color=self.edge_color, stroke_width=l)
-        asy = ""
-        for line in self.lines:
-            path = "--".join(["(%.5g,%5g)" % coords.pos() for coords in line])
-            asy += "draw(%s, %s);" % (path, pen)
-        return asy
-
 
 def _svg_bezier(*segments):
     # see https://www.w3.org/TR/SVG/paths.html#PathDataCubicBezierCommands
@@ -1509,65 +1374,6 @@ def path(max_degree, p):
             yield s
 
 
-def _asy_bezier(*segments):
-    # see http://asymptote.sourceforge.net/doc/Bezier-curves.html#Bezier-curves
-
-    while segments and not segments[0][1]:
-        segments = segments[1:]
-
-    if not segments:
-        return
-
-    def cubic(p0, p1, p2, p3):
-        return "..controls(%.5g,%.5g) and (%.5g,%.5g)..(%.5g,%.5g)" % tuple(
-            list(chain(p1, p2, p3))
-        )
-
-    def quadratric(qp0, qp1, qp2):
-        # asymptote only supports cubic beziers, so we convert this quadratic
-        # bezier to a cubic bezier, see http://fontforge.github.io/bezier.html
-
-        # CP0 = QP0
-        # CP3 = QP2
-        # CP1 = QP0 + 2 / 3 * (QP1 - QP0)
-        # CP2 = QP2 + 2 / 3 * (QP1 - QP2)
-
-        qp0x, qp0y = qp0
-        qp1x, qp1y = qp1
-        qp2x, qp2y = qp2
-
-        t = 2.0 / 3.0
-        cp0 = qp0
-        cp1 = (qp0x + t * (qp1x - qp0x), qp0y + t * (qp1y - qp0y))
-        cp2 = (qp2x + t * (qp1x - qp2x), qp2y + t * (qp1y - qp2y))
-        cp3 = qp2
-
-        return cubic(cp0, cp1, cp2, cp3)
-
-    def linear(p0, p1):
-        return "--(%.5g,%.5g)" % p1
-
-    forms = (linear, quadratric, cubic)
-
-    def path(max_degree, p):
-        max_degree = min(max_degree, len(forms))
-        while p:
-            n = min(max_degree, len(p) - 1)  # 1, 2, or 3
-            if n < 1:
-                break
-            yield forms[n - 1](*p[: n + 1])
-            p = p[n:]
-
-    k, p = segments[0]
-    yield "(%.5g,%.5g)" % p[0]
-
-    connect = []
-    for k, p in segments:
-        for s in path(k, list(chain(connect, p))):
-            yield s
-        connect = p[-1:]
-
-
 class BernsteinBasis(Builtin):
     attributes = ("Listable", "NumericFunction", "Protected")
     rules = {
@@ -1611,18 +1417,6 @@ def init(self, graphics, style, item, options):
             raise BoxConstructError
         self.spline_degree = spline_degree.get_int_value()
 
-    def to_asy(self):
-        l = self.style.get_line_width(face_element=False)
-        pen = create_pens(edge_color=self.edge_color, stroke_width=l)
-
-        asy = ""
-        for line in self.lines:
-            for path in _asy_bezier((self.spline_degree, [xy.pos() for xy in line])):
-                if path[:2] == "..":
-                    path = "(0.,0.)" + path
-                asy += "draw(%s, %s);" % (path, pen)
-        return asy
-
 
 class FilledCurve(Builtin):
     """
@@ -1688,20 +1482,6 @@ def parse_component(segments):
         else:
             raise BoxConstructError
 
-    def to_asy(self):
-        l = self.style.get_line_width(face_element=False)
-        pen = create_pens(edge_color=self.edge_color, stroke_width=l)
-
-        if not pen:
-            pen = "currentpen"
-
-        def components():
-            for component in self.components:
-                transformed = [(k, [xy.pos() for xy in p]) for k, p in component]
-                yield "fill(%s--cycle, %s);" % ("".join(_asy_bezier(*transformed)), pen)
-
-        return "".join(components())
-
     def extent(self):
         l = self.style.get_line_width(face_element=False)
         result = []
@@ -1777,52 +1557,6 @@ def process_option(self, name, value):
         else:
             raise BoxConstructError
 
-    def to_asy(self):
-        l = self.style.get_line_width(face_element=True)
-        if self.vertex_colors is None:
-            face_color = self.face_color
-        else:
-            face_color = None
-        pens = create_pens(
-            edge_color=self.edge_color,
-            face_color=face_color,
-            stroke_width=l,
-            is_face_element=True,
-        )
-        asy = ""
-        if self.vertex_colors is not None:
-            paths = []
-            colors = []
-            edges = []
-            for index, line in enumerate(self.lines):
-                paths.append(
-                    "--".join(["(%.5g,%.5g)" % coords.pos() for coords in line])
-                    + "--cycle"
-                )
-
-                # ignore opacity
-                colors.append(
-                    ",".join([color.to_asy()[0] for color in self.vertex_colors[index]])
-                )
-
-                edges.append(
-                    ",".join(["0"] + ["1"] * (len(self.vertex_colors[index]) - 1))
-                )
-
-            asy += "gouraudshade(%s, new pen[] {%s}, new int[] {%s});" % (
-                "^^".join(paths),
-                ",".join(colors),
-                ",".join(edges),
-            )
-        if pens and pens != "nullpen":
-            for line in self.lines:
-                path = (
-                    "--".join(["(%.5g,%.5g)" % coords.pos() for coords in line])
-                    + "--cycle"
-                )
-                asy += "filldraw(%s, %s);" % (path, pens)
-        return asy
-
 
 class RegularPolygon(Builtin):
     """
@@ -2129,7 +1863,7 @@ def draw(points):
 
     def make_draw_asy(self, pen):
         def draw(points):
-            for path in _asy_bezier((self.spline_degree, points)):
+            for path in asy_bezier((self.spline_degree, points)):
                 yield "draw(%s, %s);" % (path, pen)
 
         return draw
@@ -2332,23 +2066,6 @@ def draw(px, py, vx, vy, t1, s):
 
         return make
 
-    def to_asy(self):
-        width = self.style.get_line_width(face_element=False)
-        pen = create_pens(edge_color=self.edge_color, stroke_width=width)
-        polyline = self.curve.make_draw_asy(pen)
-
-        arrow_pen = create_pens(face_color=self.edge_color, stroke_width=width)
-
-        def polygon(points):
-            yield "filldraw("
-            yield "--".join(["(%.5g,%5g)" % xy for xy in points])
-            yield "--cycle, % s);" % arrow_pen
-
-        extent = self.graphics.view_width or 0
-        default_arrow = self._default_arrow(polygon)
-        custom_arrow = self._custom_arrow("asy", _ASYTransform)
-        return "".join(self._draw(polyline, default_arrow, custom_arrow, extent))
-
     def extent(self):
         width = self.style.get_line_width(face_element=False)
 
@@ -2418,20 +2135,6 @@ def extent(self):
         y = p[1] - h / 2.0 + opos[1] * h / 2.0
         return [(x, y), (x + w, y + h)]
 
-    def to_asy(self):
-        x, y = self.pos.pos()
-        content = self.content.boxes_to_tex(evaluation=self.graphics.evaluation)
-        pen = create_pens(edge_color=self.color)
-        asy = 'label("$%s$", (%s,%s), (%s,%s), %s);' % (
-            content,
-            x,
-            y,
-            -self.opos[0],
-            -self.opos[1],
-            pen,
-        )
-        return asy
-
 
 def total_extent(extents):
     xmin = xmax = ymin = ymax = None
@@ -2731,9 +2434,6 @@ def extent(self, completely_visible_only=False):
             ymax *= 2
         return xmin, xmax, ymin, ymax
 
-    def to_asy(self):
-        return "\n".join(element.to_asy() for element in self.elements)
-
     def set_size(
         self, xmin, ymin, extent_width, extent_height, pixel_width, pixel_height
     ):
@@ -2974,13 +2674,13 @@ def boxes_to_tex(self, leaves=None, **options):
         elements.view_width = w
 
         asy_completely_visible = "\n".join(
-            element.to_asy()
+            lookup_method(element, "asy")(element)
             for element in elements.elements
             if element.is_completely_visible
         )
 
         asy_regular = "\n".join(
-            element.to_asy()
+            lookup_method(element, "asy")(element)
             for element in elements.elements
             if not element.is_completely_visible
         )
@@ -2993,7 +2693,7 @@ def boxes_to_tex(self, leaves=None, **options):
         )
 
         if self.background_color is not None:
-            color, opacity = self.background_color.to_asy()
+            color, opacity = asy_color(self.background_color)
             asy_background = "filldraw(%s, %s);" % (asy_box, color)
         else:
             asy_background = ""
diff --git a/mathics/core/definitions.py b/mathics/core/definitions.py
index 57d0228be8..75cb5adb58 100644
--- a/mathics/core/definitions.py
+++ b/mathics/core/definitions.py
@@ -125,7 +125,8 @@ def __init__(
             self.builtin.update(self.user)
             self.user = {}
             self.clear_cache()
-        import mathics.formatter.svg
+        import mathics.formatter.svg  # noqa
+        import mathics.formatter.asy  # noqa
 
     def load_pymathics_module(self, module, remove_on_quit=True):
         """
diff --git a/mathics/core/formatter.py b/mathics/core/formatter.py
index c39b72f69d..b9514b52db 100644
--- a/mathics/core/formatter.py
+++ b/mathics/core/formatter.py
@@ -5,7 +5,7 @@
 format2fn = {}
 
 
-def lookup_method(self, format: str) -> Callable:
+def lookup_method(self, format: str, module_fn_name=None) -> Callable:
     """
     Find a conversion method for `format` in self's class method resolution order.
     """
@@ -19,7 +19,7 @@ def lookup_method(self, format: str) -> Callable:
     )
 
 
-def add_conversion_fn(cls) -> None:
+def add_conversion_fn(cls, module_fn_name=None) -> None:
     """Add to `format2fn` a mapping from a conversion type and builtin-class
     to a conversion method.
 
@@ -41,8 +41,13 @@ def add_conversion_fn(cls) -> None:
     # The last part of the module name is expected to be the conversion routine.
     conversion_type = module_dict["__name__"].split(".")[-1]
 
-    # Derive the conversion function from the passed-in class argument.
-    module_fn_name = cls.__name__.lower()
+    # Derive the conversion function from the passed-in class argument,
+    # unless it is already set.
+    if module_fn_name is None:
+        module_fn_name = cls.__name__.lower()
+    elif hasattr(module_fn_name, "__name__"):
+        module_fn_name = module_fn_name.__name__
+
 
     # Finally register the mapping: (Builtin-class, conversion name) -> conversion_function.
     format2fn[(conversion_type, cls)] = module_dict[module_fn_name]
diff --git a/mathics/formatter/asy.py b/mathics/formatter/asy.py
new file mode 100644
index 0000000000..a2cffe1992
--- /dev/null
+++ b/mathics/formatter/asy.py
@@ -0,0 +1,374 @@
+# -*- coding: utf-8 -*-
+
+"""
+Format a Mathics object as an Aymptote string
+"""
+
+from mathics.builtin.graphics import (
+    _Color,
+    _ArcBox,
+    ArrowBox,
+    BezierCurveBox,
+    FilledCurveBox,
+    GraphicsBox,
+    GraphicsElements,
+    InsetBox,
+    LineBox,
+    PointBox,
+    PolygonBox,
+    RGBColor,
+    RectangleBox,
+    _RoundBox,
+)
+
+from mathics.builtin.drawing.graphics3d import (
+    Graphics3DElements,
+    Line3DBox,
+    Point3DBox,
+    Polygon3DBox,
+    Sphere3DBox,
+)
+
+from mathics.core.formatter import lookup_method, add_conversion_fn
+from mathics.formatter.asy_fns import asy_bezier, asy_color as _color, asy_create_pens, asy_number
+
+
+class _ASYTransform:
+    _template = """
+    add(%s * (new picture() {
+        picture saved = currentpicture;
+        picture transformed = new picture;
+        currentpicture = transformed;
+        %s
+        currentpicture = saved;
+        return transformed;
+    })());
+    """
+
+    def __init__(self):
+        self.transforms = []
+
+    def matrix(self, a, b, c, d, e, f):
+        # a c e
+        # b d f
+        # 0 0 1
+        # see http://asymptote.sourceforge.net/doc/Transforms.html#Transforms
+        self.transforms.append("(%f, %f, %f, %f, %f, %f)" % (e, f, a, c, b, d))
+
+    def translate(self, x, y):
+        self.transforms.append("shift(%f, %f)" % (x, y))
+
+    def scale(self, x, y):
+        self.transforms.append("scale(%f, %f)" % (x, y))
+
+    def rotate(self, x):
+        self.transforms.append("rotate(%f)" % x)
+
+    def apply(self, asy):
+        return self._template % (" * ".join(self.transforms), asy)
+
+
+def arrow_box(self) -> str:
+    width = self.style.get_line_width(face_element=False)
+    pen = asy_create_pens(edge_color=self.edge_color, stroke_width=width)
+    polyline = self.curve.make_draw_asy(pen)
+
+    arrow_pen = asy_create_pens(face_color=self.edge_color, stroke_width=width)
+
+    def polygon(points):
+        yield "filldraw("
+        yield "--".join(["(%.5g,%5g)" % xy for xy in points])
+        yield "--cycle, % s);" % arrow_pen
+
+    extent = self.graphics.view_width or 0
+    default_arrow = self._default_arrow(polygon)
+    custom_arrow = self._custom_arrow("asy", _ASYTransform)
+    return "".join(self._draw(polyline, default_arrow, custom_arrow, extent))
+
+
+add_conversion_fn(ArrowBox, arrow_box)
+
+
+def bezier_curve_box(self) -> str:
+    line_width = self.style.get_line_width(face_element=False)
+    pen = asy_create_pens(edge_color=self.edge_color, stroke_width=line_width)
+
+    asy = ""
+    for line in self.lines:
+        for path in asy_bezier((self.spline_degree, [xy.pos() for xy in line])):
+            if path[:2] == "..":
+                path = "(0.,0.)" + path
+            asy += "draw(%s, %s);" % (path, pen)
+    return asy
+    if self.arc is None:
+        return _roundbox(self)
+
+    x, y, rx, ry, sx, sy, ex, ey, large_arc = self._arc_params()
+
+add_conversion_fn(BezierCurveBox, bezier_curve_box)
+
+
+def filled_curve_box(self) -> str:
+    line_width = self.style.get_line_width(face_element=False)
+    pen = asy_create_pens(edge_color=self.edge_color, stroke_width=line_width)
+
+    if not pen:
+        pen = "currentpen"
+
+    def components():
+        for component in self.components:
+            transformed = [(k, [xy.pos() for xy in p]) for k, p in component]
+            yield "fill(%s--cycle, %s);" % ("".join(asy_bezier(*transformed)), pen)
+
+    return "".join(components())
+
+
+add_conversion_fn(FilledCurveBox, filled_curve_box)
+
+
+def graphics_elements(self) -> str:
+    result = []
+    for element in self.elements:
+        format_fn = lookup_method(element, "asy")
+        if format_fn is None:
+            result.append(element.to_asy(offset))
+        else:
+            result.append(format_fn(element))
+
+    return "\n".join(result)
+
+
+add_conversion_fn(GraphicsElements, graphics_elements)
+graphics3delements = graphics_elements
+
+
+add_conversion_fn(Graphics3DElements)
+
+
+def insetbox(self) -> str:
+    x, y = self.pos.pos()
+    content = self.content.boxes_to_tex(evaluation=self.graphics.evaluation)
+    pen = asy_create_pens(edge_color=self.color)
+    asy = 'label("$%s$", (%s,%s), (%s,%s), %s);' % (
+        content,
+        x,
+        y,
+        -self.opos[0],
+        -self.opos[1],
+        pen,
+    )
+    return asy
+
+
+add_conversion_fn(InsetBox)
+
+
+def line3dbox(self):
+    # l = self.style.get_line_width(face_element=False)
+    pen = asy_create_pens(edge_color=self.edge_color, stroke_width=1)
+
+    return "".join(
+        "draw({0}, {1});".format(
+            "--".join("({0},{1},{2})".format(*coords.pos()[0]) for coords in line),
+            pen,
+        )
+        for line in self.lines
+    )
+
+
+add_conversion_fn(Line3DBox)
+
+
+def linebox(self) -> str:
+    line_width = self.style.get_line_width(face_element=False)
+    pen = asy_create_pens(edge_color=self.edge_color, stroke_width=line_width)
+    asy = ""
+    for line in self.lines:
+        path = "--".join(["(%.5g,%5g)" % coords.pos() for coords in line])
+        asy += "draw(%s, %s);" % (path, pen)
+    # print("### linebox", asy)
+    return asy
+
+
+add_conversion_fn(LineBox)
+
+
+def point3dbox(self) -> str:
+    face_color = self.face_color
+
+    # Tempoary bug fix: default Point color should be black not white
+    if list(face_color.to_rgba()[:3]) == [1, 1, 1]:
+        face_color = RGBColor(components=(0, 0, 0, face_color.to_rgba()[3]))
+
+    pen = asy_create_pens(face_color=face_color, is_face_element=False)
+
+    return "".join(
+        "path3 g={0}--cycle;dot(g, {1});".format(
+            "--".join("(%.5g,%.5g,%.5g)" % coords.pos()[0] for coords in line), pen
+        )
+        for line in self.lines
+    )
+
+
+add_conversion_fn(Point3DBox)
+
+
+def pointbox(self) -> str:
+    pen = asy_create_pens(face_color=self.face_color, is_face_element=False)
+
+    asy = ""
+    for line in self.lines:
+        for coords in line:
+            asy += "dot(%s, %s);" % (coords.pos(), pen)
+
+    # print("### pointbox", asy)
+    return asy
+
+
+add_conversion_fn(PointBox)
+
+
+def polygon3dbox(self) -> str:
+    l = self.style.get_line_width(face_element=True)
+    if self.vertex_colors is None:
+        face_color = self.face_color
+    else:
+        face_color = None
+    pen = asy_create_pens(
+        edge_color=self.edge_color,
+        face_color=face_color,
+        stroke_width=l,
+        is_face_element=True,
+    )
+
+    asy = ""
+    for line in self.lines:
+        asy += (
+            "path3 g="
+            + "--".join(["(%.5g,%.5g,%.5g)" % coords.pos()[0] for coords in line])
+            + "--cycle;"
+        )
+        asy += "draw(surface(g), %s);" % (pen)
+
+    # print("### polygon3dbox", asy)
+    return asy
+
+
+add_conversion_fn(Polygon3DBox)
+
+
+def polygonbox(self):
+    line_width = self.style.get_line_width(face_element=True)
+    if self.vertex_colors is None:
+        face_color = self.face_color
+    else:
+        face_color = None
+    pens = asy_create_pens(
+        edge_color=self.edge_color,
+        face_color=face_color,
+        stroke_width=line_width,
+        is_face_element=True,
+    )
+    asy = ""
+    if self.vertex_colors is not None:
+        paths = []
+        colors = []
+        edges = []
+        for index, line in enumerate(self.lines):
+            paths.append(
+                "--".join(["(%.5g,%.5g)" % coords.pos() for coords in line]) + "--cycle"
+            )
+
+            # ignore opacity
+            colors.append(
+                ",".join([_color(color)[0] for color in self.vertex_colors[index]])
+            )
+
+            edges.append(",".join(["0"] + ["1"] * (len(self.vertex_colors[index]) - 1)))
+
+        asy += "gouraudshade(%s, new pen[] {%s}, new int[] {%s});" % (
+            "^^".join(paths),
+            ",".join(colors),
+            ",".join(edges),
+        )
+    if pens and pens != "nullpen":
+        for line in self.lines:
+            path = (
+                "--".join(["(%.5g,%.5g)" % coords.pos() for coords in line]) + "--cycle"
+            )
+            asy += "filldraw(%s, %s);" % (path, pens)
+
+    # print("### polygonbox", asy)
+    return asy
+
+
+add_conversion_fn(PolygonBox)
+
+
+def rectanglebox(self) -> str:
+    line_width = self.style.get_line_width(face_element=True)
+    x1, y1 = self.p1.pos()
+    x2, y2 = self.p2.pos()
+    pens = asy_create_pens(
+        self.edge_color, self.face_color, line_width, is_face_element=True
+    )
+    x1, x2, y1, y2 = asy_number(x1), asy_number(x2), asy_number(y1), asy_number(y2)
+    return "filldraw((%s,%s)--(%s,%s)--(%s,%s)--(%s,%s)--cycle, %s);" % (
+        x1,
+        y1,
+        x2,
+        y1,
+        x2,
+        y2,
+        x1,
+        y2,
+        pens,
+    )
+
+
+add_conversion_fn(RectangleBox)
+
+
+def _roundbox(self):
+    x, y = self.c.pos()
+    rx, ry = self.r.pos()
+    rx -= x
+    ry -= y
+    line_width = self.style.get_line_width(face_element=self.face_element)
+    pen = asy_create_pens(
+        edge_color=self.edge_color,
+        face_color=self.face_color,
+        stroke_width=line_width,
+        is_face_element=self.face_element,
+    )
+    cmd = "filldraw" if self.face_element else "draw"
+    return "%s(ellipse((%s,%s),%s,%s), %s);" % (
+        cmd,
+        asy_number(x),
+        asy_number(y),
+        asy_number(rx),
+        asy_number(ry),
+        pen,
+    )
+
+
+add_conversion_fn(_RoundBox)
+
+
+def sphere3dbox(self) -> str:
+    # l = self.style.get_line_width(face_element=True)
+
+    if self.face_color is None:
+        face_color = (1, 1, 1)
+    else:
+        face_color = self.face_color.to_js()
+
+    return "".join(
+        "draw(surface(sphere({0}, {1})), rgb({2},{3},{4}));".format(
+            tuple(coord.pos()[0]), self.radius, *face_color[:3]
+        )
+        for coord in self.points
+    )
+
+
+add_conversion_fn(Sphere3DBox)
diff --git a/mathics/formatter/asy_fns.py b/mathics/formatter/asy_fns.py
new file mode 100644
index 0000000000..ca319d0e1e
--- /dev/null
+++ b/mathics/formatter/asy_fns.py
@@ -0,0 +1,105 @@
+# -*- coding: utf-8 -*-
+
+""" Asymptote-related functions"""
+
+from itertools import chain
+
+def asy_bezier(*segments):
+    # see http://asymptote.sourceforge.net/doc/Bezier-curves.html#Bezier-curves
+
+    while segments and not segments[0][1]:
+        segments = segments[1:]
+
+    if not segments:
+        return
+
+    def cubic(p0, p1, p2, p3):
+        return "..controls(%.5g,%.5g) and (%.5g,%.5g)..(%.5g,%.5g)" % tuple(
+            list(chain(p1, p2, p3))
+        )
+
+    def quadratric(qp0, qp1, qp2):
+        # asymptote only supports cubic beziers, so we convert this quadratic
+        # bezier to a cubic bezier, see http://fontforge.github.io/bezier.html
+
+        # CP0 = QP0
+        # CP3 = QP2
+        # CP1 = QP0 + 2 / 3 * (QP1 - QP0)
+        # CP2 = QP2 + 2 / 3 * (QP1 - QP2)
+
+        qp0x, qp0y = qp0
+        qp1x, qp1y = qp1
+        qp2x, qp2y = qp2
+
+        t = 2.0 / 3.0
+        cp0 = qp0
+        cp1 = (qp0x + t * (qp1x - qp0x), qp0y + t * (qp1y - qp0y))
+        cp2 = (qp2x + t * (qp1x - qp2x), qp2y + t * (qp1y - qp2y))
+        cp3 = qp2
+
+        return cubic(cp0, cp1, cp2, cp3)
+
+    def linear(p0, p1):
+        return "--(%.5g,%.5g)" % p1
+
+    forms = (linear, quadratric, cubic)
+
+    def path(max_degree, p):
+        max_degree = min(max_degree, len(forms))
+        while p:
+            n = min(max_degree, len(p) - 1)  # 1, 2, or 3
+            if n < 1:
+                break
+            yield forms[n - 1](*p[: n + 1])
+            p = p[n:]
+
+    k, p = segments[0]
+    yield "(%.5g,%.5g)" % p[0]
+
+    connect = []
+    for k, p in segments:
+        for s in path(k, list(chain(connect, p))):
+            yield s
+        connect = p[-1:]
+
+
+def asy_color(self) -> str:
+    """Return an asymptote string fragment for object's RGB or RGBA value"""
+    rgba = self.to_rgba()
+    alpha = rgba[3] if len(rgba) > 3 else 1.0
+    return (
+        r"rgb(%s, %s, %s)"
+        % (asy_number(rgba[0]), asy_number(rgba[1]), asy_number(rgba[2])),
+        alpha,
+    )
+
+
+def asy_create_pens(
+    edge_color=None, face_color=None, stroke_width=None, is_face_element=False
+) -> str:
+    """
+    Return an asymptote string fragment that creates a drawing pen.
+    """
+    result = []
+    if face_color is not None:
+        brush, opacity = asy_color(face_color)
+        if opacity != 1:
+            brush += "+opacity(%s)" % asy_number(opacity)
+        result.append(brush)
+    elif is_face_element:
+        result.append("nullpen")
+    if edge_color is not None:
+        pen, opacity = asy_color(edge_color)
+        if opacity != 1:
+            pen += "+opacity(%s)" % asy_number(opacity)
+        if stroke_width is not None:
+            pen += "+linewidth(%s)" % asy_number(stroke_width)
+        result.append(pen)
+    elif is_face_element:
+        result.append("nullpen")
+    return ", ".join(result)
+
+
+def asy_number(value) -> str:
+    """Format an asymptote number"""
+    return "%.5g" % value
diff --git a/mathics/formatter/svg.py b/mathics/formatter/svg.py
index 7a741efb4b..835a0e234d 100644
--- a/mathics/formatter/svg.py
+++ b/mathics/formatter/svg.py
@@ -27,7 +27,6 @@
 
 class _SVGTransform:
     def __init__(self):
-        from trepan.api import debug; debug()
         self.transforms = []
 
     def matrix(self, a, b, c, d, e, f):
@@ -77,7 +76,7 @@ def create_css(
     return "; ".join(css)
 
 
-def arrowbox(self, offset=None):
+def arrow_box(self, offset=None):
     width = self.style.get_line_width(face_element=False)
     style = create_css(edge_color=self.edge_color, stroke_width=width)
     polyline = self.curve.make_draw_svg(style)
@@ -95,7 +94,7 @@ def polygon(points):
     return "".join(self._draw(polyline, default_arrow, custom_arrow, extent))
 
 
-add_conversion_fn(ArrowBox)
+add_conversion_fn(ArrowBox, arrow_box)
 
 
 def beziercurvebox(self, offset=None):
@@ -113,7 +112,7 @@ def beziercurvebox(self, offset=None):
 add_conversion_fn(BezierCurveBox)
 
 
-def filledcurvebox(self, offset=None):
+def filled_curve_box(self, offset=None):
     line_width = self.style.get_line_width(face_element=False)
     style = create_css(
         edge_color=self.edge_color, face_color=self.face_color, stroke_width=line_width
@@ -131,9 +130,9 @@ def components():
     )
 
 
-add_conversion_fn(FilledCurveBox)
+add_conversion_fn(FilledCurveBox, filled_curve_box)
 
-def graphicsbox(self, leaves=None, **options) -> str:
+def graphics_box(self, leaves=None, **options) -> str:
         if not leaves:
             leaves = self._leaves
 
@@ -186,10 +185,10 @@ def graphicsbox(self, leaves=None, **options) -> str:
         return svg_main  # , width, height
 
 
-add_conversion_fn(GraphicsBox)
+add_conversion_fn(GraphicsBox, graphics_box)
 
 
-def graphicselements(self, offset=None):
+def graphics_elements(self, offset=None):
     result = []
     for element in self.elements:
         format_fn = lookup_method(element, "svg")
@@ -201,13 +200,13 @@ def graphicselements(self, offset=None):
     return "\n".join(result)
 
 
-add_conversion_fn(GraphicsElements)
-graphics3delements = graphicselements
+add_conversion_fn(GraphicsElements, graphics_elements)
+graphics3delements = graphics_elements
 
 add_conversion_fn(Graphics3DElements)
 
 
-def insetbox(self, offset=None):
+def inset_box(self, offset=None):
     x, y = self.pos.pos()
     if offset:
         x = x + offset[0]
@@ -238,10 +237,10 @@ def insetbox(self, offset=None):
     return svg
 
 
-add_conversion_fn(InsetBox)
+add_conversion_fn(InsetBox, inset_box)
 
 
-def linebox(self, offset=None):
+def line_box(self, offset=None):
     line_width = self.style.get_line_width(face_element=False)
     style = create_css(edge_color=self.edge_color, stroke_width=line_width)
     svg = ""
@@ -254,7 +253,7 @@ def linebox(self, offset=None):
     return svg
 
 
-add_conversion_fn(LineBox)
+add_conversion_fn(LineBox, line_box)
 
 
 def pointbox(self, offset=None):
diff --git a/test/test_formatter/test_asy.py b/test/test_formatter/test_asy.py
new file mode 100644
index 0000000000..a996de1e7c
--- /dev/null
+++ b/test/test_formatter/test_asy.py
@@ -0,0 +1,116 @@
+import re
+from mathics.core.expression import Symbol, Integer0, Integer1, Expression
+from mathics.core.evaluation import Evaluation
+from mathics.session import MathicsSession
+from mathics.builtin.inout import MakeBoxes
+from mathics.core.formatter import lookup_method
+
+session = MathicsSession(add_builtin=True, catch_interrupt=False)
+evaluation = Evaluation(session.definitions)
+
+GraphicsSymbol = Symbol("Graphics")
+ListSymbol = Symbol("List")
+
+asy_wrapper_pat = r"""^\s*
+\s*\\begin{asy}
+\s*usepackage\("amsmath"\);
+\s*size\(.+\);
+\s*
+"""
+
+
+def extract_asy_body(asy):
+    matches = re.match(asy_wrapper_pat, asy)
+    body = asy[len(matches.group(0)) :]
+    assert matches
+    print(body)
+    return body
+
+
+def get_asy(expression):
+    options = {}
+    boxes = MakeBoxes(expression).evaluate(evaluation)
+    return boxes.boxes_to_tex()
+
+
+def test_asy_circle():
+    expression = Expression(
+        GraphicsSymbol,
+        Expression("Circle", Expression(ListSymbol, Integer0, Integer0)),
+    )
+
+    asy = get_asy(expression)
+    inner_asy = extract_asy_body(asy)
+
+    # Circles are implemented as ellipses with equal major and minor axes.
+    # Check for that.
+    matches = re.match(
+        r'^draw\(ellipse\(\((.+),\s*(.+)\),(.*),(.*)\), .*', inner_asy
+    )
+
+    assert matches
+    # Check that center point is centered and
+    # major and minor axes are the same
+    assert matches.group(1) == matches.group(2)
+    assert matches.group(3) == matches.group(4)
+
+
+def test_asy_point():
+    expression = Expression(
+        GraphicsSymbol,
+        Expression("Point", Expression(ListSymbol, Integer0, Integer0)),
+    )
+
+    asy = get_asy(expression)
+    inner_asy = extract_asy_body(asy)
+
+    matches = re.match(r'^dot\(\((.+), (.+)\), .+\);.*', inner_asy)
+    assert matches
+    # Since the x,y pont is the same, we'll check that whatever this
+    # coordinate mapped to, it is the same.
+    assert matches.group(1) == matches.group(2)
+
+
+def test_asy_arrowbox():
+    expression = Expression(
+        GraphicsSymbol,
+        Expression(
+            "Arrow",
+            Expression(
+                ListSymbol,
+                Expression(ListSymbol, Integer0, Integer0),
+                Expression(ListSymbol, Integer1, Integer1),
+            ),
+        ),
+    )
+    asy = get_asy(expression)
+    inner_asy = extract_asy_body(asy)
+
+    matches = re.match(r'^draw\(.*\)', inner_asy)
+    # TODO: Match line and arrowbox
+    assert matches
+
+
+def test_asy_bezier_curve():
+
+    expression = Expression(
+        GraphicsSymbol,
+        Expression(
+            "BezierCurve",
+            Expression(
+                ListSymbol,
+                Expression(ListSymbol, Integer0, Integer0),
+                Expression(ListSymbol, Integer1, Integer1),
+            ),
+        ),
+    )
+    asy = get_asy(expression)
+    inner_asy = extract_asy_body(asy)
+
+    matches = re.match(r'^draw\(.*\)', inner_asy)
+    # TODO: Match line and arrowbox
+    assert matches
+
+
+if __name__ == "__main__":
+    test_asy_bezier_curve()
diff --git a/test/test_formatter.py b/test/test_formatter/test_svg.py
similarity index 62%
rename from test/test_formatter.py
rename to test/test_formatter/test_svg.py
index 998141a9b0..ed251ce20f 100644
--- a/test/test_formatter.py
+++ b/test/test_formatter/test_svg.py
@@ -1,5 +1,5 @@
 import re
-from mathics.core.expression import Symbol, Integer0, Expression
+from mathics.core.expression import Symbol, Integer0, Integer1, Expression
 from mathics.core.evaluation import Evaluation
 from mathics.session import MathicsSession
 from mathics.builtin.inout import MakeBoxes
@@ -63,6 +63,7 @@ def test_svg_circle():
     assert matches
     assert matches.group(1) == matches.group(2) == matches.group(3)
 
+
 def test_svg_point():
     expression = Expression(
         GraphicsSymbol,
@@ -75,12 +76,53 @@ def test_svg_point():
     # Circles are implemented as ellipses with equal major and minor axes.
     # Check for that.
     print(inner_svg)
-    matches = re.match(
-        r'^<circle cx="(\S+)" cy="(\S+)" r="(\S+)" .*/>', inner_svg
-    )
+    matches = re.match(r'^<circle cx="(\S+)" cy="(\S+)" r="(\S+)" .*/>', inner_svg)
     assert matches
     assert matches.group(1) == matches.group(2)
 
 
+def test_svg_arrowbox():
+    expression = Expression(
+        GraphicsSymbol,
+        Expression(
+            "Arrow",
+            Expression(
+                ListSymbol,
+                Expression(ListSymbol, Integer0, Integer0),
+                Expression(ListSymbol, Integer1, Integer1),
+            ),
+        ),
+    )
+    svg = get_svg(expression)
+    inner_svg = extract_svg_body(svg)
+
+    matches = re.match(r'^<polyline points=".+"\s+style=".*"\s*/><', inner_svg)
+    # TODO: Could pick endpoint of this line and match with beginnign of arrow polygon below
+    assert matches
+    arrow_polygon = inner_svg[len(matches.group(0)) - 1 :]
+    matches = re.match(r'^<polygon points=".+"\s+style=".*"\s*/>', arrow_polygon)
+    assert matches
+
+
+def test_svg_bezier_curve():
+
+    expression = Expression(
+        GraphicsSymbol,
+        Expression(
+            "BezierCurve",
+            Expression(
+                ListSymbol,
+                Expression(ListSymbol, Integer0, Integer0),
+                Expression(ListSymbol, Integer1, Integer1),
+            ),
+        ),
+    )
+    svg = get_svg(expression)
+    inner_svg = extract_svg_body(svg)
+
+    matches = re.match(r'^<path d=".+"\s*/>', inner_svg)
+    assert matches
+
+
 if __name__ == "__main__":
-    test_svg_point()
+    test_svg_bezier_curve()