Initial population.

1 files changed, 1200 insertions, 0 deletions

diff --git a/src/algebra/draw.spad.pamphlet b/src/algebra/draw.spad.pamphlet
new file mode 100644
index 00000000..de5dc3e3
--- /dev/null
+++ b/src/algebra/draw.spad.pamphlet
@@ -0,0 +1,1200 @@
+\documentclass{article}
+\usepackage{axiom}
+\begin{document}
+\title{\$SPAD/src/algebra draw.spad}
+\author{Clifton J. Williamson, Scott Morrison, Jon Steinbach, Mike Dewar}
+\maketitle
+\begin{abstract}
+\end{abstract}
+\eject
+\tableofcontents
+\eject
+\section{package DRAWCFUN TopLevelDrawFunctionsForCompiledFunctions}
+<<package DRAWCFUN TopLevelDrawFunctionsForCompiledFunctions>>=
+)abbrev package DRAWCFUN TopLevelDrawFunctionsForCompiledFunctions
+++ Author: Clifton J. Williamson
+++ Date Created: 22 June 1990
+++ Date Last Updated: January 1992 by Scott Morrison
+++ Basic Operations: draw, recolor
+++ Related Constructors:
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ Description: TopLevelDrawFunctionsForCompiledFunctions provides top level 
+++ functions for drawing graphics of expressions.
+TopLevelDrawFunctionsForCompiledFunctions():
+ Exports == Implementation where
+  ANY1 ==> AnyFunctions1
+  B    ==> Boolean
+  F    ==> Float
+  L    ==> List
+  SEG  ==> Segment Float
+  SF   ==> DoubleFloat
+  DROP ==> DrawOption
+  PLOT ==> Plot
+  PPC  ==> ParametricPlaneCurve(SF -> SF)
+  PSC  ==> ParametricSpaceCurve(SF -> SF)
+  PSF  ==> ParametricSurface((SF,SF) -> SF)
+  Pt   ==> Point SF
+  PSFUN ==> (SF, SF) -> Pt
+  PCFUN ==> SF -> Pt
+  SPACE3 ==> ThreeSpace(SF)
+  VIEW2 ==> TwoDimensionalViewport
+  VIEW3 ==> ThreeDimensionalViewport
+
+  Exports ==> with
+
+--% Two Dimensional Function Plots
+
+    draw: (SF -> SF,SEG,L DROP) -> VIEW2
+      ++ draw(f,a..b,l) draws the graph of \spad{y = f(x)} as x
+      ++ ranges from \spad{min(a,b)} to \spad{max(a,b)}.
+      ++ The options contained in the list l of
+      ++ the domain \spad{DrawOption} are applied.
+    draw: (SF -> SF,SEG) -> VIEW2
+      ++ draw(f,a..b) draws the graph of \spad{y = f(x)} as x
+      ++ ranges from \spad{min(a,b)} to \spad{max(a,b)}.
+
+--% Parametric Plane Curves
+
+    draw: (PPC,SEG,L DROP) -> VIEW2
+      ++ draw(curve(f,g),a..b,l) draws the graph of the parametric
+      ++ curve \spad{x = f(t), y = g(t)} as t ranges from \spad{min(a,b)} to 
+      ++ \spad{max(a,b)}.
+      ++ The options contained in the list l of the domain \spad{DrawOption}
+      ++ are applied.
+    draw: (PPC,SEG) -> VIEW2
+      ++ draw(curve(f,g),a..b) draws the graph of the parametric
+      ++ curve \spad{x = f(t), y = g(t)} as t ranges from \spad{min(a,b)} to 
+      ++ \spad{max(a,b)}.
+
+--% Parametric Space Curves
+
+    draw: (PSC,SEG,L DROP) -> VIEW3
+      ++ draw(curve(f,g,h),a..b,l) draws the graph of the parametric
+      ++ curve \spad{x = f(t), y = g(t), z = h(t)} as t ranges from 
+      ++ \spad{min(a,b)} to \spad{max(a,b)}.
+      ++ The options contained in the list l of the domain
+      ++ \spad{DrawOption} are applied.
+    draw: (PSC,SEG) -> VIEW3
+      ++ draw(curve(f,g,h),a..b,l) draws the graph of the parametric
+      ++ curve \spad{x = f(t), y = g(t), z = h(t)} as t ranges from 
+      ++ \spad{min(a,b)} to \spad{max(a,b)}.
+    draw: (PCFUN,SEG,L DROP) -> VIEW3
+      ++ draw(f,a..b,l) draws the graph of the parametric
+      ++ curve \spad{f} as t ranges from 
+      ++ \spad{min(a,b)} to \spad{max(a,b)}.
+      ++ The options contained in the list l of the domain
+      ++ \spad{DrawOption} are applied.
+    draw: (PCFUN,SEG) -> VIEW3
+      ++ draw(f,a..b,l) draws the graph of the parametric
+      ++ curve \spad{f} as t ranges from 
+      ++ \spad{min(a,b)} to \spad{max(a,b)}.
+
+    makeObject: (PSC,SEG,L DROP) -> SPACE3
+      ++ makeObject(curve(f,g,h),a..b,l) returns a space of the
+      ++ domain \spadtype{ThreeSpace} which contains the graph of the
+      ++ parametric curve \spad{x = f(t), y = g(t), z = h(t)} as t ranges from 
+      ++ \spad{min(a,b)} to \spad{max(a,b)};
+      ++ The options contained in the list l of the domain
+      ++ \spad{DrawOption} are applied.
+    makeObject: (PSC,SEG) -> SPACE3
+      ++ makeObject(sp,curve(f,g,h),a..b) returns the space \spad{sp}
+      ++ of the domain \spadtype{ThreeSpace} with the addition of the graph
+      ++ of the parametric curve \spad{x = f(t), y = g(t), z = h(t)} as t
+      ++ ranges from \spad{min(a,b)} to \spad{max(a,b)}.
+    makeObject: (PCFUN,SEG,L DROP) -> SPACE3
+      ++ makeObject(curve(f,g,h),a..b,l) returns a space of the
+      ++ domain \spadtype{ThreeSpace} which contains the graph of the
+      ++ parametric curve \spad{x = f(t), y = g(t), z = h(t)} as t ranges from 
+      ++ \spad{min(a,b)} to \spad{max(a,b)}.
+      ++ The options contained in the list l of the domain
+      ++ \spad{DrawOption} are applied.
+    makeObject: (PCFUN,SEG) -> SPACE3
+      ++ makeObject(sp,curve(f,g,h),a..b) returns the space \spad{sp}
+      ++ of the domain \spadtype{ThreeSpace} with the addition of the graph
+      ++ of the parametric curve \spad{x = f(t), y = g(t), z = h(t)} as t
+      ++ ranges from \spad{min(a,b)} to \spad{max(a,b)}.
+
+--% Three Dimensional Function Plots
+
+    draw: ((SF,SF) -> SF,SEG,SEG,L DROP) -> VIEW3
+      ++ draw(f,a..b,c..d,l) draws the graph of \spad{z = f(x,y)}
+      ++ as x ranges from \spad{min(a,b)} to \spad{max(a,b)} and y ranges from
+      ++ \spad{min(c,d)} to \spad{max(c,d)}.
+      ++ and the options contained in the list l of the domain
+      ++ \spad{DrawOption} are applied.
+    draw: ((SF,SF) -> SF,SEG,SEG) -> VIEW3
+      ++ draw(f,a..b,c..d) draws the graph of \spad{z = f(x,y)}
+      ++ as x ranges from \spad{min(a,b)} to \spad{max(a,b)} and y ranges from
+      ++ \spad{min(c,d)} to \spad{max(c,d)}.
+    makeObject: ((SF,SF) -> SF,SEG,SEG,L DROP) -> SPACE3
+      ++ makeObject(f,a..b,c..d,l) returns a space of the domain
+      ++ \spadtype{ThreeSpace} which contains the graph of \spad{z = f(x,y)}
+      ++ as x ranges from \spad{min(a,b)} to \spad{max(a,b)} and y ranges from
+      ++ \spad{min(c,d)} to \spad{max(c,d)}, and the options contained in the
+      ++ list l of the domain \spad{DrawOption} are applied.
+    makeObject: ((SF,SF) -> SF,SEG,SEG) -> SPACE3
+      ++ makeObject(f,a..b,c..d) returns a space of the domain
+      ++ \spadtype{ThreeSpace} which contains the graph of \spad{z = f(x,y)}
+      ++ as x ranges from \spad{min(a,b)} to \spad{max(a,b)} and y ranges from
+      ++ \spad{min(c,d)} to \spad{max(c,d)}.
+
+--% Parametric Surfaces
+
+    draw: (PSFUN, SEG, SEG, L DROP) -> VIEW3
+      ++ draw(f,a..b,c..d) draws the
+      ++ graph of the parametric surface \spad{f(u,v)}
+      ++ as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)}.
+      ++ The options contained in the
+      ++ list l of the domain \spad{DrawOption} are applied.
+    draw: (PSFUN, SEG, SEG) -> VIEW3
+      ++ draw(f,a..b,c..d) draws the
+      ++ graph of the parametric surface \spad{f(u,v)}
+      ++ as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)}
+      ++ The options contained in the list
+      ++ l of the domain \spad{DrawOption} are applied.
+    makeObject: (PSFUN, SEG, SEG, L DROP) -> SPACE3
+      ++ makeObject(f,a..b,c..d,l) returns a
+      ++ space of the domain \spadtype{ThreeSpace} which contains the
+      ++ graph of the parametric surface \spad{f(u,v)}
+      ++ as u ranges from \spad{min(a,b)} to
+      ++ \spad{max(a,b)} and v ranges from \spad{min(c,d)} to \spad{max(c,d)};
+      ++ The options contained in the
+      ++ list l of the domain \spad{DrawOption} are applied.
+    makeObject: (PSFUN, SEG, SEG) -> SPACE3
+      ++ makeObject(f,a..b,c..d,l) returns a
+      ++ space of the domain \spadtype{ThreeSpace} which contains the
+      ++ graph of the parametric surface \spad{f(u,v)}
+      ++ as u ranges from \spad{min(a,b)} to
+      ++ \spad{max(a,b)} and v ranges from \spad{min(c,d)} to \spad{max(c,d)}.
+    draw: (PSF,SEG,SEG,L DROP) -> VIEW3
+      ++ draw(surface(f,g,h),a..b,c..d) draws the
+      ++ graph of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)};
+      ++ The options contained in the
+      ++ list l of the domain \spad{DrawOption} are applied.
+    draw: (PSF,SEG,SEG) -> VIEW3
+      ++ draw(surface(f,g,h),a..b,c..d) draws the
+      ++ graph of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)};
+    makeObject: (PSF,SEG,SEG,L DROP) -> SPACE3
+      ++ makeObject(surface(f,g,h),a..b,c..d,l) returns a
+      ++ space of the domain \spadtype{ThreeSpace} which contains the
+      ++ graph of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to
+      ++ \spad{max(a,b)} and v ranges from \spad{min(c,d)} to \spad{max(c,d)}.
+      ++ The options contained in the
+      ++ list l of the domain \spad{DrawOption} are applied.
+    makeObject: (PSF,SEG,SEG) -> SPACE3
+      ++ makeObject(surface(f,g,h),a..b,c..d,l) returns a
+      ++ space of the domain \spadtype{ThreeSpace} which contains the
+      ++ graph of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to
+      ++ \spad{max(a,b)} and v ranges from \spad{min(c,d)} to \spad{max(c,d)}.
+    recolor: ((SF,SF) -> Pt,(SF,SF,SF) -> SF) -> ((SF,SF) -> Pt)
+      ++ recolor(), uninteresting to top level user; exported in order to 
+      ++ compile package.
+
+  Implementation ==> add
+    --!!  I have had to work my way around the following bug in the compiler:
+    --!!  When a local variable is given a mapping as a value, e.g.
+    --!!  foo : SF -> SF := makeFloatFunction(f,t),
+    --!!  the compiler cannot distinguish that local variable from a local
+    --!!  function defined elsewhere in the package.  Thus, when 'foo' is
+    --!!  passed to a function, e.g.
+    --!!  bird := fcn(foo),
+    --!!  foo will often be compiled as |DRAW;foo| rather than |foo|. This,
+    --!!  of course, causes a run-time error.
+    --!!  To avoid this problem, local variables are not given mappings as
+    --!!  values, but rather (singleton) lists of mappings.  The first element
+    --!!  of the list can always be extracted and everything goes through
+    --!!  as before.  There is no major loss in efficiency, as the computation
+    --!!  of points will always dominate the computation time.
+    --!!                                     - cjw,  22 June MCMXC
+
+    import PLOT
+    import TwoDimensionalPlotClipping
+    import GraphicsDefaults
+    import ViewportPackage
+    import ThreeDimensionalViewport
+    import DrawOptionFunctions0
+    import MakeFloatCompiledFunction(Ex)
+    import MeshCreationRoutinesForThreeDimensions
+    import SegmentFunctions2(SF,Float)
+    import ViewDefaultsPackage
+    import AnyFunctions1(Pt -> Pt)
+    import AnyFunctions1((SF,SF,SF) -> SF)
+    import DrawOptionFunctions0
+    import SPACE3
+
+    EXTOVARERROR : String := _
+      "draw: when specifying function, left hand side must be a variable"
+    SMALLRANGEERROR : String := _
+      "draw: range is in interval with only one point"
+    DEPVARERROR : String := _
+      "draw: independent variable appears on lhs of function definition"
+
+------------------------------------------------------------------------
+--                     2D - draw's  
+------------------------------------------------------------------------
+
+    drawToScaleRanges: (Segment SF,Segment SF) -> L SEG
+    drawToScaleRanges(xVals,yVals) ==
+      -- warning: assumes window is square
+      xHi := convert(hi xVals)@Float; xLo := convert(lo xVals)@Float
+      yHi := convert(hi yVals)@Float; yLo := convert(lo yVals)@Float
+      xDiff := xHi - xLo; yDiff := yHi - yLo
+      pad := abs(yDiff - xDiff)/2
+      yDiff > xDiff =>
+        [segment(xLo - pad,xHi + pad),map(convert(#1)@Float,yVals)]
+      [map(convert(#1)@Float,xVals),segment(yLo - pad,yHi + pad)]
+
+    drawPlot: (PLOT,L DROP) -> VIEW2
+    drawPlot(plot,l) ==
+      branches := listBranches plot
+      xRange := xRange plot; yRange := yRange plot
+      -- process clipping information
+      if (cl := option(l,"clipSegment" :: Symbol)) case "failed" then
+        if clipBoolean(l,clipPointsDefault()) then
+          clipInfo :=
+            parametric? plot => clipParametric plot
+            clip plot
+          branches := clipInfo.brans
+          xRange := clipInfo.xValues; yRange := clipInfo.yValues
+        else
+          "No explicit user-specified clipping"
+      else
+        segList := retract(cl :: Any)$ANY1(L SEG)
+        empty? segList =>
+          error "draw: you may specify at least 1 segment for 2D clipping"
+        more?(segList,2) =>
+          error "draw: you may specify at most 2 segments for 2D clipping"
+        xLo : SF := 0; xHi : SF := 0; yLo : SF := 0; yHi : SF := 0
+        if empty? rest segList then
+          xLo := lo xRange; xHi := hi xRange
+          yRangeF := first segList
+          yLo := convert(lo yRangeF)@SF; yHi := convert(hi yRangeF)@SF
+        else
+          xRangeF := first segList
+          xLo := convert(lo xRangeF)@SF; xHi := convert(hi xRangeF)@SF
+          yRangeF := second segList
+          yLo := convert(lo yRangeF)@SF; yHi := convert(hi yRangeF)@SF
+        clipInfo := clipWithRanges(branches,xLo,xHi,yLo,yHi)
+        branches := clipInfo.brans
+        xRange := clipInfo.xValues; yRange := clipInfo.yValues
+      -- process scaling information
+      if toScale(l,drawToScale()) then
+        scaledRanges := drawToScaleRanges(xRange,yRange)
+        -- add scaled ranges to list of options
+        l := concat(ranges scaledRanges,l)
+      else
+        xRangeFloat : SEG := map(convert(#1)@Float,xRange)
+        yRangeFloat : SEG := map(convert(#1)@Float,yRange)
+        -- add ranges to list of options
+        l := concat(ranges(ll : L SEG := [xRangeFloat,yRangeFloat]),l)
+      -- process color information
+      ptCol := pointColorPalette(l,pointColorDefault())
+      crCol := curveColorPalette(l,lineColorDefault())
+      -- draw
+      drawCurves(branches,ptCol,crCol,pointSizeDefault(),l)
+
+    normalize: SEG -> Segment SF
+    normalize seg ==
+      -- normalize [a,b]:
+      -- error if a = b, returns [a,b] if a < b, returns [b,a] if b > a
+      a := convert(lo seg)@SF; b := convert(hi seg)@SF
+      a = b => error SMALLRANGEERROR
+      a < b => segment(a,b)
+      segment(b,a)
+
+--% functions for creation of maps SF -> Point SF (two dimensional)
+
+    myTrap1: (SF-> SF, SF) -> SF
+    myTrap1(ff:SF-> SF, f:SF):SF ==
+      s := trapNumericErrors(ff(f))$Lisp :: Union(SF, "failed")
+      s case "failed" => _$NaNvalue$Lisp
+      r:=s::SF
+      r >max()$SF or r < min()$SF => _$NaNvalue$Lisp
+      r
+
+    makePt2: (SF,SF) -> Point SF
+    makePt2(x,y) == point(l : List SF := [x,y])
+
+--% Two Dimensional Function Plots
+ 
+    draw(f:SF -> SF,seg:SEG,l:L DROP) ==
+      -- set adaptive plotting off or on
+      oldAdaptive := adaptive?()$PLOT
+      setAdaptive(adaptive(l,oldAdaptive))$PLOT
+      -- create function SF -> Point SF
+      ff : L(SF -> Point SF) := [makePt2(myTrap1(f,#1),#1)]
+      -- process change of coordinates
+      if (c := option(l,"coordinates" :: Symbol)) case "failed" then
+        -- default coordinate transformation
+        ff := [makePt2(#1,myTrap1(f,#1))]
+      else
+        cc : L(Pt -> Pt) := [retract(c :: Any)$ANY1(Pt -> Pt)]
+        ff := [(first cc)((first ff)(#1))]
+      -- create PLOT
+      pl := pointPlot(first ff,normalize seg)
+      -- reset adaptive plotting
+      setAdaptive(oldAdaptive)$PLOT
+      -- draw
+      drawPlot(pl,l)
+ 
+    draw(f:SF -> SF,seg:SEG) == draw(f,seg,nil())
+ 
+--% Parametric Plane Curves
+
+    draw(ppc:PPC,seg:SEG,l:L DROP) ==
+      -- set adaptive plotting off or on
+      oldAdaptive := adaptive?()$PLOT
+      setAdaptive(adaptive(l,oldAdaptive))$PLOT
+      -- create function SF -> Point SF
+      f := coordinate(ppc,1); g := coordinate(ppc,2)
+      fcn : L(SF -> Pt) := [makePt2(myTrap1(f,#1),myTrap1(g,#1))]
+      -- process change of coordinates
+      if not (c := option(l,"coordinates" :: Symbol)) case "failed" then
+        cc : L(Pt -> Pt) := [retract(c :: Any)$ANY1(Pt -> Pt)]
+        fcn := [(first cc)((first fcn)(#1))]
+      -- create PLOT
+      pl := pointPlot(first fcn,normalize seg)
+      -- reset adaptive plotting
+      setAdaptive(oldAdaptive)$PLOT
+      -- draw
+      drawPlot(pl,l)
+ 
+    draw(ppc:PPC,seg:SEG) == draw(ppc,seg,nil())
+
+------------------------------------------------------------------------
+--                     3D - Curves  
+------------------------------------------------------------------------
+
+--% functions for creation of maps SF -> Point SF (three dimensional)
+
+    makePt4: (SF,SF,SF,SF) -> Point SF
+    makePt4(x,y,z,c) == point(l : List SF := [x,y,z,c])
+
+--% Parametric Space Curves
+
+    id: SF -> SF
+    id x == x
+
+    zCoord: (SF,SF,SF) -> SF
+    zCoord(x,y,z) == z
+
+    colorPoints: (List List Pt,(SF,SF,SF) -> SF) -> List List Pt
+    colorPoints(llp,func) ==
+      for lp in llp repeat for p in lp repeat
+        p.4 := func(p.1,p.2,p.3)
+      llp
+
+    makeObject(psc:PSC,seg:SEG,l:L DROP) ==
+      sp := space l
+      -- obtain dependent variable and coordinate functions
+      f := coordinate(psc,1); g := coordinate(psc,2); h := coordinate(psc,3)
+      -- create function SF -> Point SF with default or user-specified
+      -- color function
+      fcn : L(SF -> Pt) := [makePt4(myTrap1(f,#1),myTrap1(g,#1),myTrap1(h,#1),_
+                            myTrap1(id,#1))]
+      pointsColored? : Boolean := false
+      if not (c1 := option(l,"colorFunction1" :: Symbol)) case "failed" then
+        pointsColored? := true
+        fcn := [makePt4(myTrap1(f,#1),myTrap1(g,#1),myTrap1(h,#1),_
+                retract(c1 :: Any)$ANY1(SF -> SF)(#1))]
+      -- process change of coordinates
+      if not (c := option(l,"coordinates" :: Symbol)) case "failed" then
+        cc : L(Pt -> Pt) := [retract(c :: Any)$ANY1(Pt -> Pt)]
+        fcn := [(first cc)((first fcn)(#1))]
+      -- create PLOT
+      pl := pointPlot(first fcn,normalize seg)$Plot3D
+      -- create ThreeSpace
+      s := sp
+      -- draw Tube
+--      print(pl::OutputForm)
+      option?(l,"tubeRadius" :: Symbol) =>
+        pts := tubePoints(l,8)
+        rad := convert(tubeRadius(l,0.25))@DoubleFloat
+        tub := tube(pl,rad,pts)$NumericTubePlot(Plot3D)
+        loops := listLoops tub
+        -- color points if this has not been done already
+        if not pointsColored? then
+          if (c3 := option(l,"colorFunction3" :: Symbol)) case "failed"
+            then colorPoints(loops,zCoord)  -- default color function
+            else colorPoints(loops,retract(c3 :: Any)$ANY1((SF,SF,SF) -> SF))
+        mesh(s,loops,false,false)
+        s
+      -- draw curve
+      br := listBranches pl
+      for b in br repeat curve(s,b)
+      s
+
+    makeObject(psc:PCFUN,seg:SEG,l:L DROP) ==
+      sp := space l
+      -- create function SF -> Point SF with default or user-specified
+      -- color function
+      fcn : L(SF -> Pt) := [psc]
+      pointsColored? : Boolean := false
+      if not (c1 := option(l,"colorFunction1" :: Symbol)) case "failed" then
+        pointsColored? := true
+        fcn := [concat(psc(#1), retract(c1 :: Any)$ANY1(SF -> SF)(#1))]
+      -- process change of coordinates
+      if not (c := option(l,"coordinates" :: Symbol)) case "failed" then
+        cc : L(Pt -> Pt) := [retract(c :: Any)$ANY1(Pt -> Pt)]
+        fcn := [(first cc)((first fcn)(#1))]
+      -- create PLOT
+      pl := pointPlot(first fcn,normalize seg)$Plot3D
+      -- create ThreeSpace
+      s := sp
+      -- draw Tube
+      option?(l,"tubeRadius" :: Symbol) =>
+        pts := tubePoints(l,8)
+        rad := convert(tubeRadius(l,0.25))@DoubleFloat
+        tub := tube(pl,rad,pts)$NumericTubePlot(Plot3D)
+        loops := listLoops tub
+        -- color points if this has not been done already
+        mesh(s,loops,false,false)
+        s
+      -- draw curve
+      br := listBranches pl
+      for b in br repeat curve(s,b)
+      s
+
+    makeObject(psc:PSC,seg:SEG) ==
+      makeObject(psc,seg,nil())
+
+    makeObject(psc:PCFUN,seg:SEG) ==
+      makeObject(psc,seg,nil())
+
+    draw(psc:PSC,seg:SEG,l:L DROP) ==
+      sp := makeObject(psc,seg,l)
+      makeViewport3D(sp, l)
+
+    draw(psc:PSC,seg:SEG) ==
+      draw(psc,seg,nil())
+
+    draw(psc:PCFUN,seg:SEG,l:L DROP) ==
+      sp := makeObject(psc,seg,l)
+      makeViewport3D(sp, l)
+
+    draw(psc:PCFUN,seg:SEG) ==
+      draw(psc,seg,nil())
+
+------------------------------------------------------------------------
+--                     3D - Surfaces  
+------------------------------------------------------------------------
+
+    myTrap2: ((SF, SF) -> SF, SF, SF) -> SF
+    myTrap2(ff:(SF, SF) -> SF, u:SF, v:SF):SF ==
+      s := trapNumericErrors(ff(u, v))$Lisp :: Union(SF, "failed")
+      s case "failed" =>  _$NaNvalue$Lisp
+      r:SF := s::SF
+      r >max()$SF or r < min()$SF => _$NaNvalue$Lisp
+      r
+
+    recolor(ptFunc,colFunc) ==
+      pt := ptFunc(#1,#2)
+      pt.4 := colFunc(pt.1,pt.2,pt.3)
+      pt
+
+    xCoord: (SF,SF) -> SF
+    xCoord(x,y) == x
+
+--% Three Dimensional Function Plots
+
+    makeObject(f:(SF,SF) -> SF,xSeg:SEG,ySeg:SEG,l:L DROP) ==
+      sp := space l
+      -- process color function of two variables
+      col2 : L((SF,SF) -> SF) := [xCoord]     -- dummy color function
+      pointsColored? : Boolean := false
+      if not (c2 := option(l,"colorFunction2" :: Symbol)) case "failed" then
+        pointsColored? := true
+        col2 := [retract(c2 :: Any)$ANY1((SF,SF) -> SF)]
+      fcn : L((SF,SF) -> Pt) :=
+        [makePt4(myTrap2(f,#1,#2),#1,#2,(first col2)(#1,#2))]
+      -- process change of coordinates
+      if (c := option(l,"coordinates" :: Symbol)) case "failed" then
+        -- default coordinate transformation
+        fcn := [makePt4(#1,#2,myTrap2(f,#1,#2),(first col2)(#1,#2))]
+      else
+        cc : L(Pt -> Pt) := [retract(c :: Any)$ANY1(Pt -> Pt)]
+        fcn := [(first cc)((first fcn)(#1,#2))]
+      -- process color function of three variables, if there was no
+      -- color function of two variables
+      if not pointsColored? then
+        c := option(l,"colorFunction3" :: Symbol)
+        fcn := 
+          c case "failed" => [recolor((first fcn),zCoord)]
+          [recolor((first fcn),retract(c :: Any)$ANY1((SF,SF,SF) -> SF))]
+      -- create mesh
+      mesh := meshPar2Var(sp,first fcn,normalize xSeg,normalize ySeg,l)
+      mesh
+
+    makeObject(f:(SF,SF) -> SF,xSeg:SEG,ySeg:SEG) ==
+      makeObject(f,xSeg,ySeg,nil())
+
+    draw(f:(SF,SF) -> SF,xSeg:SEG,ySeg:SEG,l:L DROP) ==
+      sp := makeObject(f, xSeg, ySeg, l)
+      makeViewport3D(sp, l)
+
+    draw(f:(SF,SF) -> SF,xSeg:SEG,ySeg:SEG) ==
+      draw(f,xSeg,ySeg,nil())
+
+--% parametric surface
+
+    makeObject(s:PSF,uSeg:SEG,vSeg:SEG,l:L DROP) ==
+      sp := space l
+      -- create functions from expressions
+      f : L((SF,SF) -> SF) := [coordinate(s,1)]
+      g : L((SF,SF) -> SF) := [coordinate(s,2)]
+      h : L((SF,SF) -> SF) := [coordinate(s,3)]
+      -- process color function of two variables
+      col2 : L((SF,SF) -> SF) := [xCoord]     -- dummy color function
+      pointsColored? : Boolean := false
+      if not (c2 := option(l,"colorFunction2" :: Symbol)) case "failed" then
+        pointsColored? := true
+        col2 := [retract(c2 :: Any)$ANY1((SF,SF) -> SF)]
+      fcn : L((SF,SF) -> Pt) := 
+        [makePt4(myTrap2((first f),#1,#2),myTrap2((first g),#1,#2),myTrap2((first h),#1,#2),_
+                 myTrap2((first col2),#1,#2))]
+      -- process change of coordinates
+      if not (c := option(l,"coordinates" :: Symbol)) case "failed" then
+        cc : L(Pt -> Pt) := [retract(c :: Any)$ANY1(Pt -> Pt)]
+        fcn := [(first cc)((first fcn)(#1,#2))]
+      -- process color function of three variables, if there was no
+      -- color function of two variables
+      if not pointsColored? then
+        col3 : L((SF,SF,SF) -> SF) := [zCoord]  -- default color function
+        if not (c := option(l,"colorFunction3" :: Symbol)) case "failed" then 
+          col3 := [retract(c :: Any)$ANY1((SF,SF,SF) -> SF)]
+        fcn := [recolor((first fcn),(first col3))]
+      -- create mesh
+      mesh := meshPar2Var(sp,first fcn,normalize uSeg,normalize vSeg,l)
+      mesh
+
+    makeObject(s:PSFUN,uSeg:SEG,vSeg:SEG,l:L DROP) ==
+      sp := space l
+      -- process color function of two variables
+      col2 : L((SF,SF) -> SF) := [xCoord]     -- dummy color function
+      pointsColored? : Boolean := false
+      if not (c2 := option(l,"colorFunction2" :: Symbol)) case "failed" then
+        pointsColored? := true
+        col2 := [retract(c2 :: Any)$ANY1((SF,SF) -> SF)]
+      fcn : L((SF,SF) -> Pt) := 
+        pointsColored? => [concat(s(#1, #2), (first col2)(#1, #2))]
+        [s]
+      -- process change of coordinates
+      if not (c := option(l,"coordinates" :: Symbol)) case "failed" then
+        cc : L(Pt -> Pt) := [retract(c :: Any)$ANY1(Pt -> Pt)]
+        fcn := [(first cc)((first fcn)(#1,#2))]
+      -- create mesh
+      mesh := meshPar2Var(sp,first fcn,normalize uSeg,normalize vSeg,l)
+      mesh
+
+    makeObject(s:PSF,uSeg:SEG,vSeg:SEG) ==
+      makeObject(s,uSeg,vSeg,nil())
+
+    draw(s:PSF,uSeg:SEG,vSeg:SEG,l:L DROP) ==
+      -- draw
+      mesh := makeObject(s,uSeg,vSeg,l)
+      makeViewport3D(mesh,l)
+
+    draw(s:PSF,uSeg:SEG,vSeg:SEG) ==
+      draw(s,uSeg,vSeg,nil())
+ 
+    makeObject(s:PSFUN,uSeg:SEG,vSeg:SEG) ==
+      makeObject(s,uSeg,vSeg,nil())
+
+    draw(s:PSFUN,uSeg:SEG,vSeg:SEG,l:L DROP) ==
+      -- draw
+      mesh := makeObject(s,uSeg,vSeg,l)
+      makeViewport3D(mesh,l)
+
+    draw(s:PSFUN,uSeg:SEG,vSeg:SEG) ==
+      draw(s,uSeg,vSeg,nil())
+ 
+@
+\section{package DRAW TopLevelDrawFunctions}
+<<package DRAW TopLevelDrawFunctions>>=
+)abbrev package DRAW TopLevelDrawFunctions
+++ Author: Clifton J. Williamson
+++ Date Created: 23 January 1990
+++ Date Last Updated: October 1991 by Jon Steinbach
+++ Basic Operations: draw
+++ Related Constructors:
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ Description: TopLevelDrawFunctions provides top level functions for 
+++ drawing graphics of expressions.
+TopLevelDrawFunctions(Ex:Join(ConvertibleTo InputForm,SetCategory)):
+ Exports == Implementation where
+  B    ==> Boolean
+  BIND ==> SegmentBinding Float
+  L    ==> List
+  SF   ==> DoubleFloat
+  DROP ==> DrawOption
+
+  PPC  ==> ParametricPlaneCurve Ex
+  PPCF ==> ParametricPlaneCurve(SF -> SF)
+  PSC  ==> ParametricSpaceCurve Ex
+  PSCF ==> ParametricSpaceCurve(SF -> SF)
+  PSF  ==> ParametricSurface Ex
+  PSFF ==> ParametricSurface((SF,SF) -> SF)
+  SPACE3 ==> ThreeSpace(SF)
+  VIEW2 ==> TwoDimensionalViewport
+  VIEW3 ==> ThreeDimensionalViewport
+
+  Exports ==> with
+
+--% Two Dimensional Function Plots
+
+    draw: (Ex,BIND,L DROP) -> VIEW2
+      ++ draw(f(x),x = a..b,l) draws the graph of \spad{y = f(x)} as x
+      ++ ranges from \spad{min(a,b)} to \spad{max(a,b)}; \spad{f(x)} is the 
+      ++ default title, and the options contained in the list l of
+      ++ the domain \spad{DrawOption} are applied.
+    draw: (Ex,BIND) -> VIEW2
+      ++ draw(f(x),x = a..b) draws the graph of \spad{y = f(x)} as x
+      ++ ranges from \spad{min(a,b)} to \spad{max(a,b)}; \spad{f(x)} appears 
+      ++ in the title bar.
+
+--% Parametric Plane Curves
+
+    draw: (PPC,BIND,L DROP) -> VIEW2
+      ++ draw(curve(f(t),g(t)),t = a..b,l) draws the graph of the parametric
+      ++ curve \spad{x = f(t), y = g(t)} as t ranges from \spad{min(a,b)} to 
+      ++ \spad{max(a,b)}; \spad{(f(t),g(t))} is the default title, and the
+      ++ options contained in the list l of the domain \spad{DrawOption}
+      ++ are applied.
+    draw: (PPC,BIND) -> VIEW2
+      ++ draw(curve(f(t),g(t)),t = a..b) draws the graph of the parametric
+      ++ curve \spad{x = f(t), y = g(t)} as t ranges from \spad{min(a,b)} to 
+      ++ \spad{max(a,b)}; \spad{(f(t),g(t))} appears in the title bar.
+
+--% Parametric Space Curves
+
+    draw: (PSC,BIND,L DROP) -> VIEW3
+      ++ draw(curve(f(t),g(t),h(t)),t = a..b,l) draws the graph of the
+      ++ parametric curve \spad{x = f(t)}, \spad{y = g(t)}, \spad{z = h(t)}
+      ++ as t ranges from \spad{min(a,b)} to \spad{max(a,b)}; \spad{h(t)}
+      ++ is the default title, and the options contained in the list l of
+      ++ the domain \spad{DrawOption} are applied.
+    draw: (PSC,BIND) -> VIEW3
+      ++ draw(curve(f(t),g(t),h(t)),t = a..b) draws the graph of the parametric
+      ++ curve \spad{x = f(t)}, \spad{y = g(t)}, \spad{z = h(t)} as t ranges
+      ++ from \spad{min(a,b)} to \spad{max(a,b)}; \spad{h(t)} is the default
+      ++ title.
+    makeObject: (PSC,BIND,L DROP) -> SPACE3
+      ++ makeObject(curve(f(t),g(t),h(t)),t = a..b,l) returns a space of
+      ++ the domain \spadtype{ThreeSpace} which contains the graph of the
+      ++ parametric curve \spad{x = f(t)}, \spad{y = g(t)}, \spad{z = h(t)}
+      ++ as t ranges from \spad{min(a,b)} to \spad{max(a,b)}; \spad{h(t)}
+      ++ is the default title, and the options contained in the list l of
+      ++ the domain \spad{DrawOption} are applied.
+    makeObject: (PSC,BIND) -> SPACE3
+      ++ makeObject(curve(f(t),g(t),h(t)),t = a..b) returns a space of the
+      ++ domain \spadtype{ThreeSpace} which contains the graph of the
+      ++ parametric curve \spad{x = f(t)}, \spad{y = g(t)}, \spad{z = h(t)}
+      ++ as t ranges from \spad{min(a,b)} to \spad{max(a,b)}; \spad{h(t)} is
+      ++ the default title.
+
+--% Three Dimensional Function Plots
+
+    draw: (Ex,BIND,BIND,L DROP) -> VIEW3
+      ++ draw(f(x,y),x = a..b,y = c..d,l) draws the graph of \spad{z = f(x,y)}
+      ++ as x ranges from \spad{min(a,b)} to \spad{max(a,b)} and y ranges from
+      ++ \spad{min(c,d)} to \spad{max(c,d)}; \spad{f(x,y)} is the default
+      ++ title, and the options contained in the list l of the domain
+      ++ \spad{DrawOption} are applied.
+    draw: (Ex,BIND,BIND) -> VIEW3
+      ++ draw(f(x,y),x = a..b,y = c..d) draws the graph of \spad{z = f(x,y)}
+      ++ as x ranges from \spad{min(a,b)} to \spad{max(a,b)} and y ranges from
+      ++ \spad{min(c,d)} to \spad{max(c,d)}; \spad{f(x,y)} appears in the title bar.
+    makeObject: (Ex,BIND,BIND,L DROP) -> SPACE3
+      ++ makeObject(f(x,y),x = a..b,y = c..d,l) returns a space of the
+      ++ domain \spadtype{ThreeSpace} which contains the graph of
+      ++ \spad{z = f(x,y)} as x ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and y ranges from \spad{min(c,d)} to \spad{max(c,d)}; \spad{f(x,y)}
+      ++ is the default title, and the options contained in the list l of the
+      ++ domain \spad{DrawOption} are applied.
+    makeObject: (Ex,BIND,BIND) -> SPACE3
+      ++ makeObject(f(x,y),x = a..b,y = c..d) returns a space of the domain
+      ++ \spadtype{ThreeSpace} which contains the graph of \spad{z = f(x,y)}
+      ++ as x ranges from \spad{min(a,b)} to \spad{max(a,b)} and y ranges from
+      ++ \spad{min(c,d)} to \spad{max(c,d)}; \spad{f(x,y)} appears as the
+      ++ default title.
+
+--% Parametric Surfaces
+
+    draw: (PSF,BIND,BIND,L DROP) -> VIEW3
+      ++ draw(surface(f(u,v),g(u,v),h(u,v)),u = a..b,v = c..d,l) draws the
+      ++ graph of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)}; \spad{h(t)}
+      ++ is the default title, and the options contained in the list l of
+      ++ the domain \spad{DrawOption} are applied.
+    draw: (PSF,BIND,BIND) -> VIEW3
+      ++ draw(surface(f(u,v),g(u,v),h(u,v)),u = a..b,v = c..d) draws the
+      ++ graph of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)}; \spad{h(t)} is
+      ++ the default title.
+    makeObject: (PSF,BIND,BIND,L DROP) -> SPACE3
+      ++ makeObject(surface(f(u,v),g(u,v),h(u,v)),u = a..b,v = c..d,l) returns
+      ++ a space of the domain \spadtype{ThreeSpace} which contains the graph
+      ++ of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)}; \spad{h(t)} is
+      ++ the default title, and the options contained in the list l of
+      ++ the domain \spad{DrawOption} are applied.
+    makeObject: (PSF,BIND,BIND) -> SPACE3
+      ++ makeObject(surface(f(u,v),g(u,v),h(u,v)),u = a..b,v = c..d) returns
+      ++ a space of the domain \spadtype{ThreeSpace} which contains the
+      ++ graph of the parametric surface \spad{x = f(u,v)}, \spad{y = g(u,v)},
+      ++ \spad{z = h(u,v)} as u ranges from \spad{min(a,b)} to \spad{max(a,b)}
+      ++ and v ranges from \spad{min(c,d)} to \spad{max(c,d)}; \spad{h(t)} is
+      ++ the default title.
+
+  Implementation ==> add
+    import TopLevelDrawFunctionsForCompiledFunctions
+    import MakeFloatCompiledFunction(Ex)
+    import ParametricPlaneCurve(SF -> SF)
+    import ParametricSpaceCurve(SF -> SF)
+    import ParametricSurface((SF,SF) -> SF)
+    import ThreeSpace(SF)
+
+------------------------------------------------------------------------
+--                     2D - draw's  (given by formulae)
+------------------------------------------------------------------------
+
+--% Two Dimensional Function Plots
+ 
+    draw(f:Ex,bind:BIND,l:L DROP) ==
+      -- create title if necessary
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM2D",l)
+        else l := concat(title s,l)
+      -- call 'draw'
+      draw(makeFloatFunction(f,variable bind),segment bind,l)
+ 
+    draw(f:Ex,bind:BIND) == draw(f,bind,nil())
+ 
+--% Parametric Plane Curves
+
+    draw(ppc:PPC,bind:BIND,l:L DROP) ==
+      f := coordinate(ppc,1); g := coordinate(ppc,2)
+      -- create title if necessary
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM2D",l)
+        else l := concat(title s,l)
+      -- create curve with functions as coordinates
+      curve : PPCF := curve(makeFloatFunction(f,variable bind),_
+                            makeFloatFunction(g,variable bind))$PPCF
+      -- call 'draw'
+      draw(curve,segment bind,l)
+ 
+    draw(ppc:PPC,bind:BIND) == draw(ppc,bind,nil())
+
+------------------------------------------------------------------------
+--                     3D - Curves  (given by formulas)
+------------------------------------------------------------------------
+
+    makeObject(psc:PSC,tBind:BIND,l:L DROP) ==
+      -- obtain dependent variable and coordinate functions
+      t := variable tBind; tSeg := segment tBind
+      f := coordinate(psc,1); g := coordinate(psc,2); h := coordinate(psc,3)
+      -- create title if necessary
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM3D",l)
+        else l := concat(title s,l)
+      -- indicate draw style if necessary
+      if not option?(l,"style" :: Symbol) then
+        l := concat(style unparse(convert(f)@InputForm),l)
+      -- create curve with functions as coordinates
+      curve : PSCF := curve(makeFloatFunction(f,t),_
+                            makeFloatFunction(g,t),_
+                            makeFloatFunction(h,t))
+      -- call 'draw'
+      makeObject(curve,tSeg,l)
+
+    makeObject(psc:PSC,tBind:BIND) ==
+      makeObject(psc,tBind,nil())
+
+    draw(psc:PSC,tBind:BIND,l:L DROP) ==
+      -- obtain dependent variable and coordinate functions
+      t := variable tBind; tSeg := segment tBind
+      f := coordinate(psc,1); g := coordinate(psc,2); h := coordinate(psc,3)
+      -- create title if necessary
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM3D",l)
+        else l := concat(title s,l)
+      -- indicate draw style if necessary
+      if not option?(l,"style" :: Symbol) then
+        l := concat(style unparse(convert(f)@InputForm),l)
+      -- create curve with functions as coordinates
+      curve : PSCF := curve(makeFloatFunction(f,t),_
+                            makeFloatFunction(g,t),_
+                            makeFloatFunction(h,t))
+      -- call 'draw'
+      draw(curve,tSeg,l)
+
+    draw(psc:PSC,tBind:BIND) ==
+      draw(psc,tBind,nil())
+
+------------------------------------------------------------------------
+--                     3D - Surfaces  (given by formulas)
+------------------------------------------------------------------------
+
+--% Three Dimensional Function Plots
+
+    makeObject(f:Ex,xBind:BIND,yBind:BIND,l:L DROP) ==
+      -- create title if necessary
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM3D",l)
+        else l := concat(title s,l)
+      -- indicate draw style if necessary
+      if not option?(l,"style" :: Symbol) then
+        l := concat(style unparse(convert(f)@InputForm),l)
+      -- obtain dependent variables and their ranges
+      x := variable xBind; xSeg := segment xBind
+      y := variable yBind; ySeg := segment yBind
+      -- call 'draw'
+      makeObject(makeFloatFunction(f,x,y),xSeg,ySeg,l)
+
+    makeObject(f:Ex,xBind:BIND,yBind:BIND) ==
+      makeObject(f,xBind,yBind,nil())
+
+    draw(f:Ex,xBind:BIND,yBind:BIND,l:L DROP) ==
+      -- create title if necessary
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM3D",l)
+        else l := concat(title s,l)
+      -- indicate draw style if necessary
+      if not option?(l,"style" :: Symbol) then
+        l := concat(style unparse(convert(f)@InputForm),l)
+      -- obtain dependent variables and their ranges
+      x := variable xBind; xSeg := segment xBind
+      y := variable yBind; ySeg := segment yBind
+      -- call 'draw'
+      draw(makeFloatFunction(f,x,y),xSeg,ySeg,l)
+
+    draw(f:Ex,xBind:BIND,yBind:BIND) ==
+      draw(f,xBind,yBind,nil())
+
+--% parametric surface
+
+    makeObject(s:PSF,uBind:BIND,vBind:BIND,l:L DROP) ==
+      f := coordinate(s,1); g := coordinate(s,2); h := coordinate(s,3)
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM3D",l)
+        else l := concat(title s,l)
+      if not option?(l,"style" :: Symbol) then
+        l := concat(style unparse(convert(f)@InputForm),l)
+      u := variable uBind; uSeg := segment uBind
+      v := variable vBind; vSeg := segment vBind
+      surf : PSFF := surface(makeFloatFunction(f,u,v),_
+                             makeFloatFunction(g,u,v),_
+                             makeFloatFunction(h,u,v))
+      makeObject(surf,uSeg,vSeg,l)
+
+    makeObject(s:PSF,uBind:BIND,vBind:BIND) ==
+      makeObject(s,uBind,vBind,nil())
+
+    draw(s:PSF,uBind:BIND,vBind:BIND,l:L DROP) ==
+      f := coordinate(s,1); g := coordinate(s,2); h := coordinate(s,3)
+      -- create title if necessary
+      if not option?(l,"title" :: Symbol) then
+        s:String := unparse(convert(f)@InputForm)
+        if sayLength(s)$DisplayPackage > 50 then
+          l := concat(title "AXIOM3D",l)
+        else l := concat(title s,l)
+      -- indicate draw style if necessary
+      if not option?(l,"style" :: Symbol) then
+        l := concat(style unparse(convert(f)@InputForm),l)
+      -- obtain dependent variables and their ranges
+      u := variable uBind; uSeg := segment uBind
+      v := variable vBind; vSeg := segment vBind
+      -- create surface with functions as coordinates
+      surf : PSFF := surface(makeFloatFunction(f,u,v),_
+                             makeFloatFunction(g,u,v),_
+                             makeFloatFunction(h,u,v))
+      -- call 'draw'
+      draw(surf,uSeg,vSeg,l)
+
+    draw(s:PSF,uBind:BIND,vBind:BIND) ==
+      draw(s,uBind,vBind,nil())
+
+@
+\section{package DRAWCURV TopLevelDrawFunctionsForAlgebraicCurves}
+<<package DRAWCURV TopLevelDrawFunctionsForAlgebraicCurves>>=
+)abbrev package DRAWCURV TopLevelDrawFunctionsForAlgebraicCurves
+++ Author: Clifton J. Williamson
+++ Date Created: 26 June 1990
+++ Date Last Updated:  October 1991 by Jon Steinbach
+++ Basic Operations: draw
+++ Related Constructors:
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ Description: TopLevelDrawFunctionsForAlgebraicCurves provides top level 
+++ functions for drawing non-singular algebraic curves.
+
+TopLevelDrawFunctionsForAlgebraicCurves(R,Ex): Exports == Implementation where
+  R  : Join(IntegralDomain, OrderedSet, RetractableTo Integer)
+  Ex : FunctionSpace(R)
+
+  ANY1  ==> AnyFunctions1
+  DROP  ==> DrawOption
+  EQ    ==> Equation
+  F     ==> Float
+  FRAC  ==> Fraction
+  I     ==> Integer
+  L     ==> List
+  P     ==> Polynomial
+  RN    ==> Fraction Integer
+  SEG   ==> Segment
+  SY    ==> Symbol
+  VIEW2 ==> TwoDimensionalViewport
+
+  Exports ==> with
+
+    draw: (EQ Ex,SY,SY,L DROP) -> VIEW2
+      ++ draw(f(x,y) = g(x,y),x,y,l) draws the graph of a polynomial
+      ++ equation.  The list l of draw options must specify a region
+      ++ in the plane in which the curve is to sketched.
+
+  Implementation ==> add
+    import ViewportPackage
+    import PlaneAlgebraicCurvePlot
+    import ViewDefaultsPackage
+    import GraphicsDefaults
+    import DrawOptionFunctions0
+    import SegmentFunctions2(RN,F)
+    import SegmentFunctions2(F,RN)
+    import AnyFunctions1(L SEG RN)
+
+    drawToScaleRanges: (SEG F,SEG F) -> L SEG F
+    drawToScaleRanges(xVals,yVals) ==
+      -- warning: assumes window is square
+      xHi := hi xVals; xLo := lo xVals
+      yHi := hi yVals; yLo := lo yVals
+      xDiff := xHi - xLo; yDiff := yHi - yLo
+      pad := abs(yDiff - xDiff)/2
+      yDiff > xDiff =>
+        [segment(xLo - pad,xHi + pad),yVals]
+      [xVals,segment(yLo - pad,yHi + pad)]
+
+    intConvert: R -> I
+    intConvert r ==
+      (nn := retractIfCan(r)@Union(I,"failed")) case "failed" =>
+        error "draw: polynomial must have rational coefficients"
+      nn :: I
+
+    polyEquation: EQ Ex -> P I
+    polyEquation eq ==
+      ff := lhs(eq) - rhs(eq)
+      (r := retractIfCan(ff)@Union(FRAC P R,"failed")) case "failed" =>
+        error "draw: not a polynomial equation"
+      rat := r :: FRAC P R
+      retractIfCan(denom rat)@Union(R,"failed") case "failed" =>
+        error "draw: non-constant denominator"
+      map(intConvert,numer rat)$PolynomialFunctions2(R,I)
+
+    draw(eq,x,y,l) ==
+      -- obtain polynomial equation
+      p := polyEquation eq
+      -- extract ranges from option list
+      floatRange := option(l,"rangeFloat" :: Symbol)
+      ratRange := option(l,"rangeRat" :: Symbol)
+      (floatRange case "failed") and (ratRange case "failed") =>
+        error "draw: you must specify ranges for an implicit plot"
+      ranges : L SEG RN := nil()             -- dummy value
+      floatRanges : L SEG F := nil()         -- dummy value
+      xRange : SEG RN := segment(0,0)        -- dummy value
+      yRange : SEG RN := segment(0,0)        -- dummy value
+      xRangeFloat : SEG F := segment(0,0)    -- dummy value
+      yRangeFloat : SEG F := segment(0,0)    -- dummy value
+      if not ratRange case "failed" then
+        ranges := retract(ratRange :: Any)$ANY1(L SEG RN)
+        not size?(ranges,2) => error "draw: you must specify two ranges"
+        xRange := first ranges; yRange := second ranges
+        xRangeFloat := map(convert(#1)@Float,xRange)@(SEG F)
+        yRangeFloat := map(convert(#1)@Float,yRange)@(SEG F)
+        floatRanges := [xRangeFloat,yRangeFloat]
+      else
+        floatRanges := retract(floatRange :: Any)$ANY1(L SEG F)
+        not size?(floatRanges,2) =>
+          error "draw: you must specify two ranges"
+        xRangeFloat := first floatRanges
+        yRangeFloat := second floatRanges
+        xRange := map(retract(#1)@RN,xRangeFloat)@(SEG RN)
+        yRange := map(retract(#1)@RN,yRangeFloat)@(SEG RN)
+        ranges := [xRange,yRange]
+      -- create curve plot
+      acplot := makeSketch(p,x,y,xRange,yRange)
+      -- process scaling information
+      if toScale(l,drawToScale()) then
+        scaledRanges := drawToScaleRanges(xRangeFloat,yRangeFloat)
+        -- add scaled ranges to list of options
+        l := concat(ranges scaledRanges,l)
+      else
+        -- add ranges to list of options
+        l := concat(ranges floatRanges,l)
+      -- process color information
+      ptCol := pointColorPalette(l,pointColorDefault())
+      crCol := curveColorPalette(l,lineColorDefault())
+      -- draw
+      drawCurves(listBranches acplot,ptCol,crCol,pointSizeDefault(),l)
+
+@
+\section{package DRAWPT TopLevelDrawFunctionsForPoints}
+<<package DRAWPT TopLevelDrawFunctionsForPoints>>=
+)abbrev package DRAWPT TopLevelDrawFunctionsForPoints
+++ Author: Mike Dewar
+++ Date Created: 24 May 1995
+++ Date Last Updated: 25 November 1996
+++ Basic Operations: draw
+++ Related Constructors:
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ Description: TopLevelDrawFunctionsForPoints provides top level functions for 
+++ drawing curves and surfaces described by sets of points.
+ 
+TopLevelDrawFunctionsForPoints(): Exports == Implementation where
+
+  DROP  ==> DrawOption
+  L     ==> List
+  SF    ==> DoubleFloat
+  Pt    ==> Point SF
+  VIEW2 ==> TwoDimensionalViewport
+  VIEW3 ==> ThreeDimensionalViewport
+
+  Exports ==> with
+    draw: (L SF,L SF) -> VIEW2
+      ++ draw(lx,ly) plots the curve constructed of points (x,y) for x
+      ++ in \spad{lx} for y in \spad{ly}.
+    draw: (L SF,L SF,L DROP) -> VIEW2
+      ++ draw(lx,ly,l) plots the curve constructed of points (x,y) for x
+      ++ in \spad{lx} for y in \spad{ly}.
+      ++ The options contained in the list l of
+      ++ the domain \spad{DrawOption} are applied.
+    draw: (L Pt) -> VIEW2
+      ++ draw(lp) plots the curve constructed from the list of points lp.
+    draw: (L Pt,L DROP) -> VIEW2
+      ++ draw(lp,l) plots the curve constructed from the list of points lp.
+      ++ The options contained in the list l of the domain \spad{DrawOption}
+      ++ are applied.
+    draw: (L SF, L SF, L SF) -> VIEW3
+      ++ draw(lx,ly,lz) draws the surface constructed by projecting the values
+      ++ in the \axiom{lz} list onto the rectangular grid formed by the 
+      ++ \axiom{lx X ly}.
+    draw: (L SF, L SF, L SF, L DROP) -> VIEW3
+      ++ draw(lx,ly,lz,l) draws the surface constructed by projecting the values
+      ++ in the \axiom{lz} list onto the rectangular grid formed by the 
+      ++ The options contained in the list l of the domain \spad{DrawOption}
+      ++ are applied.
+
+  Implementation ==> add
+
+    draw(lp:L Pt,l:L DROP):VIEW2 ==
+      makeViewport2D(makeGraphImage([lp])$GraphImage,l)$VIEW2
+
+    draw(lp:L Pt):VIEW2 == draw(lp,[])
+
+    draw(lx: L SF, ly: L SF, l:L DROP):VIEW2 ==
+      draw([point([x,y])$Pt for x in lx for y in ly],l)
+
+    draw(lx: L SF, ly: L SF):VIEW2 == draw(lx,ly,[])
+
+    draw(x:L SF,y:L SF,z:L SF):VIEW3 == draw(x,y,z,[])
+
+    draw(x:L SF,y:L SF,z:L SF,l:L DROP):VIEW3 ==
+      m  : Integer := #x
+      zero? m => error "No X values"
+      n  : Integer := #y
+      zero? n => error "No Y values"
+      zLen : Integer := #z
+      zLen ~= (m*n) => 
+        zLen > (m*n) => error "Too many Z-values to fit grid"
+        error "Not enough Z-values to fit grid"
+      points : L L Pt := []
+      for j in n..1 by -1 repeat
+        row : L Pt := []
+        for i in m..1 by -1 repeat
+          zval := (j-1)*m+i
+          row := cons(point([x.i,y.j,z.zval,z.zval]),row)
+        points := cons(row,points)
+      makeViewport3D(mesh points,l)
+
+@
+\section{License}
+<<license>>=
+--Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd.
+--All rights reserved.
+--
+--Redistribution and use in source and binary forms, with or without
+--modification, are permitted provided that the following conditions are
+--met:
+--
+--    - Redistributions of source code must retain the above copyright
+--      notice, this list of conditions and the following disclaimer.
+--
+--    - Redistributions in binary form must reproduce the above copyright
+--      notice, this list of conditions and the following disclaimer in
+--      the documentation and/or other materials provided with the
+--      distribution.
+--
+--    - Neither the name of The Numerical ALgorithms Group Ltd. nor the
+--      names of its contributors may be used to endorse or promote products
+--      derived from this software without specific prior written permission.
+--
+--THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+--IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+--TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+--PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
+--OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+--EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+--PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+--PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+--LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+--NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+--SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+@
+<<*>>=
+<<license>>
+
+<<package DRAWCFUN TopLevelDrawFunctionsForCompiledFunctions>>
+<<package DRAW TopLevelDrawFunctions>>
+<<package DRAWCURV TopLevelDrawFunctionsForAlgebraicCurves>>
+<<package DRAWPT TopLevelDrawFunctionsForPoints>>
+@
+\eject
+\begin{thebibliography}{99}
+\bibitem{1} nothing
+\end{thebibliography}
+\end{document}