MATHEMATICA IN ACTION for Issue 8-2 of Mathematica in Education and Research

The Costa Surface, in Snow and in Mathematica

The ability to finely control 3-dimensional polygons allows us to create pleasing views and instructive animations of the surprising and mysterious Costa surface.

Dedicated to the Memory of Al Gray

by  Dan Schwalbe and Stan Wagon

It came as a great surprise when, in 1983, Celso Costa, a graduate student in Brazil discovered an entirely new shape that met all the rigorous conditions of being a minimal surface. Further work by Hoffman, Hoffman, and Meeks led to stunning visualizations and new results about the surface, such as the surprising amount of symmetry it has. Thus when S.W. proposed to sculptor Helaman Ferguson that he captain a team for the Breckenridge International Snow Sculpture Championships in Jan. 1999, the Costa design was a natural. D.S. and Macalester senior Tamas Nemeth filled out the team and we spent four fabulous days in the rarefied air of Breckenridge scupting a 12-foot-high Costs surface. This article amplifies earlier work on the Costa surface in Mathematica to obtain an esthetic view of the surface by cutting along a sphere, and also an animation that shows that the Costa surface is a torus in disguise.

.  Minimal Surfaces

The discovery in 1983 of a new minimal surface by Celso Costa, then a Brazilian graduate student, was remarkable. Before then only three such surfaces (complete and embeddable in (ℝ^3) with no self-intersections) were known: the plane, the catenoid, and the helicoid. A minimal surface is one such that any small perturbation at any point increases its area. We will restrict attention here to examples that are embeddable (there are examples that can be visualized in  (ℝ^3) but only with self-intersections, analogous to the familiar situation with the Klein bottle) and complete (have no borders: the borders that do appear in the images of this article should be viewed merely as signposts en route to infinity).

It is clear that the plane works, and easy to believe that the helicoid, with its flat surfaces, works. To understand why a catenoid (the surface of revolution of a catenary) is a minimal surface, one would have to delve into the differential geometry of the surface. Note the different topologies: the plane has the topology of a sphere with one point removed (genus 0, one puncture); the other two have the topology of a sphere with two points removed (genus 0, two punctures). The Costa surface is striking for its complex topology: it is a torus with three points removed; a movie showing a transformation from one to the other is presented in section 4. For more information on minimal surfaces and the important role played by the Costa discovery, see [Hoffman, 1987]. While it is perhaps stretching things to say that the Costa surface generalizes to 3-space the familiar notion of a line, the analogy is worth mentioning: perturb a line anywhere and its length goes up; the same is true for this remarkable family of surfaces.

One striking consequence of being a minimal surface is that every point has negative Gaussian curvature. Think of the familiar saddle point: it curves positively in one direction and negatively in another. This is not the case for a bump or a pit on a plane; such features have positive curvature. Thus every point on a minimal surface is a saddle point. This adds structural stability to the shape and is one reason why we chose it for the snow sculpture (Figure ). Because of the anticipated strength, we felt we could make the surface very thin. And indeed we did. Invisible Handshake (so-called because the Costa surface is the complementary space of two hands coming together but not touching) is 12 feet tall, but less than 4 inches thick. For more information on the Breckenridge Snow Sculpture Championships and the diverse physical and mental challenges required to sculpt such an object out of a 20-ton block of hard-packed snow in only a few days, see [Ferguson et al, 1999].

Figures  and  show the catenoid and helicoid. As noted, they are topologically equivalent, and one can make a movie in Mathematica that shows one being transformed to the other [Ogawa 1992].

Needs["Graphics`SurfaceOfRevolution`"]

SurfaceOfRevolution[Cosh[u], {u, -2.6, 2.6}, RevolutionAxis -> {1, 0, 0}, 
   BoxRatios -> {1, 1, 1}, ViewPoint -> {0.1, 2, 0}];

SurfaceOfRevolution[Cosh[u], {u, -2.399999999999999, 2}, 
   RevolutionAxis -> {1, 0, 0}, BoxRatios -> {1, 1, 1}, 
   ViewPoint -> {2, 0.4, 0.3}];

Code for Figure

Block[{$DisplayFunction = Identity},
	im1=SurfaceOfRevolution[Cosh[u], {u,  -2.6,  2.6},
RevolutionAxis -> {1, 0, 0}, Axes -> None, ColorOutput -> GrayLevel,BoxRatios -> {1, 1, 1},Boxed -> False, ViewPoint -> {0.1, 2, 0}];
\(im2 = SurfaceOfRevolution[Cosh[u], \ {u, \ \(-2.4\), \ 2. }, 

                RevolutionAxis -> {1, 0, 0}, \ Axes -> None, 
                ColorOutput -> GrayLevel, 
			BoxRatios -> {1, 1, 1}, 

                			\ Boxed -> False, \ ViewPoint -> {2, 0.4, 0.3}]\)];
Show[GraphicsArray[{im1,  im2}]];

Show[Block[{$DisplayFunction = Identity}, 
    {ParametricPlot3D[{u*Cos[t], u*Sin[t], t}, {t, 0, 4*Pi}, {u, -1, 1}, 
      PlotPoints -> {100, 20}], 
     ParametricPlot3D[{Cos[t], Sin[t], t, Thickness[0.02]}, {t, 0, 4*Pi}, 
      PlotPoints -> 100], ParametricPlot3D[{-Cos[t], -Sin[t], t, 
       {Thickness[0.02], Hue[0]}}, {t, 0, 4*Pi}, PlotPoints -> 100]}], 
   BoxRatios -> {1, 1, 2}, ViewPoint -> {0, -2, 1}, Boxed -> False, 
   Axes -> None];

Code for Figure

\(\(Show[
      Block[{$DisplayFunction = Identity}, 

        		{ParametricPlot3D[{u\ Cos[t], u\ Sin[t], t}, {t, 0, 4\ π}, {
              u, \(-1\), 1}, PlotPoints -> {100, 20}], 
			
          ParametricPlot3D[{Cos[t], Sin[t], t, Thickness[0.02]}, {t, 0, 
              4\ π}, PlotPoints -> 100], 
			
          ParametricPlot3D[{\(-Cos[t]\), \(-Sin[t]\), 
              t, {Thickness[0.02], GrayLevel[ .7]}}, {t, 0, 4\ π}, 
            PlotPoints -> 100]}], BoxRatios -> {1, 1, 2}, 
      ViewPoint -> {0, \(-2\), 1}, Boxed -> False, 
      Axes -> None, ColorOutput -> GrayLevel, 
      BoxRatios -> {1, 1, 2}]; \)\)

.  A Parametrization of the Costa Surface

The late Al Gray found some relatively simple parametric equations that allow the surface to be generated by a single ParametricPlot3D command in Mathematica. More details of his work can be found in [Ferguson et al, 1996]. These equations were critical to the speedy generation of the surface and allowed Ferguson to use his computer-driven drill to sculpt a Costa surface out of stone. The main tools are the special functions known as Weierstrass functions. First we define two constants and an abbreviation for a Weierstrass functions; this is followed by the Gray formula.

c = WeierstrassInvariants[{0.5, 0.5*I}][[1]]; e = WeierstrassP[0.5, {c, 0}]; 
  Z = WeierstrassZeta[#1, {c, 0}] & ; 
  GrayCosta[z_] := 
   With[{a = z + (Z[z - 1/2] - Z[z - I/2])/(2*e), 
     b = Z[z]/Pi - (Pi*(1 - I))/(4*e), d = WeierstrassP[z, {c, 0}]}, 
    Pi/2*{Re[a - b], Im[a + b], Log[Abs[(d - e)/(d + e)]]/Sqrt[2*Pi]}]

This function is periodic with period 1 in each direction. And the function is infinite at (0, 0), ((1/2), 0), and (0, (1/2)); these three spikes correspond to the three ends of the surface. The fat columns in Figure  correspond to a single asymptote, and this is the equatorial plane in the Costa images that we will see shortly.

norm[v_] := Sqrt[v . v]

Show[Graphics3D[EdgeForm[]], Plot3D[norm[GrayCosta[u + I*v]], 
    {u, -0.3, 2.299999999999999}, {v, -0.3, 1.6}, PlotPoints -> 120, 
    DisplayFunction -> Identity], PlotRange -> {0, 3.5}];

Code for Figure

\(\(Show[Graphics3D[EdgeForm[]], 
	\ \ 
      Plot3D[norm[GrayCosta[u + I\ v]], {u, \(-0.3\), 2.3}, {v, \(-0.3\), 
          1.6}, PlotPoints -> 120, DisplayFunction -> Identity], 
      PlotRange -> {0, 3.5}, \ ColorOutput -> GrayLevel, \ 
      PlotRegion\  -> \ {{0, 1}, \ {\(- .3\), 1.2}}]; \)\)

A straightforward plot leads to difficulties because some of the polygons will extend to vast distances. So, following [Ferguson 1996], we generate the polygons and then remove those that lie entirely outside the sphere of radius 5.5. Our goal is to get a view of that part of the surface lying inside the sphere of radius 2.7.

We first set some options and a useful viewpoint.

\(br = {6.6, 6.4, 4.9}; 

    pr = {{\(-3.4\), 3.2}, {\(-3.2\), 3.2}, {\(-2.4\), 2.5}}; 

    vp = {\(-2\), 0.2, 1.4}; 

    SetOptions[Graphics3D, Boxed -> False, BoxRatios -> br, 
      PlotRange -> pr]; \)

We next define a function to produce the set of polygons we want: those lying inside the sphere of radius 5.5.

\(sphereQ[pts_, r_] := And@@\((\(# . # < r\^2\ &\)\ /@\ pts)\)
\)\(\(costaPolys[res_, r_] := 	
      Cases[
		
        ParametricPlot3D[
          Evaluate[GrayCosta[u + \ I\ v]], {u, 0.001, 1.001}, {v, 0.001, 
            1.001}, PlotPoints -> res, 
			
          DisplayFunction\  -> Identity], 
			
        Polygon[pts_] /; sphereQ[pts, r], ∞]; \)\)

We generate and store a large collection of these polygons using an 100×100 grid.

\(CostaPolygons = costaPolys[100, 5.5]; 
Length[CostaPolygons]\)

And we now compare our polygon set to the large sphere. It takes some trial and error to come up with suitable sphere-radii. Figure  shows, at the left, the general Costa shape. The plane forming an equator is what we must modify to get our spherical cutoff. The right-hand image shows the sphere of radius 2.7. We can see that our polygons adequately extend through the sphere; thus we should be able to modify the protruding polygons so that they stop at the sphere. A trick similar to one mentioned by Allan Hayes in his column in the last issue are used to generate a sphere of solid color while the surface parts use the default illumination. We use Hayes’s idea of overloading RGBColor with arguments beyond 1 — Mathematica hiccups at this, but does enhance the color appropriately — and an undocumented feature that a vector can be added as second argument to Polygon that plays the role of a fake normal vector. By adding the same vector — {1, 0, 0} — to all the polygons of the sphere, we get a yellowish sphere. The howmany function counts the number of points lyiung outside a sphere.

\(rad = 2.7; 
howmany[pts_, r_] := Count[pts, p_ /; p . p > r\^2]\)

Needs["Graphics`Shapes`"]

Remove@Sphere

Show[
      Graphics3D[CostaPolygons, ViewVertical -> {0, -1, 0.6}, 
        ViewPoint -> vp]];

\(Off[Graphics::realu]; 

    Show[Graphics3D[{{SurfaceColor[RGBColor[2, 2, 0.6]], 
            Sphere[rad, 20, 20] /. \ 
			
              Polygon[x_]\  :> \ Polygon[x, {1, 0, 1}]}, 
          DeleteCases[CostaPolygons, \ 
            Polygon[pts_] /; \ howmany[pts, rad] == 0, ∞]}, 
        ViewPoint -> {0.4, \(-2\), 0.38}, \ 
        PlotRange -> {{\(-3\), 3}, {\(-3\), 3}, {\(-3\), 3}}, \ 
        BoxRatios -> {1, 1, 1}]]; \)

Code for Figure

\(\(Show[
      GraphicsArray[{
		
          Graphics3D[CostaPolygons, ViewVertical -> {0, \(-1\), 0.6}, \ 
            ViewPoint -> vp, \ ColorOutput -> GrayLevel], 
		
          Graphics3D[{{SurfaceColor[RGBColor[2, 2, 0.6]], 
                Sphere[rad, 20, 20] /. 
                  Polygon[x_] :> Polygon[x, {1, 0, 0}]}, 
              DeleteCases[CostaPolygons, \ 
                Polygon[pts_] /; \ howmany[pts, rad] == 0, \ ∞]}, 
            ViewPoint -> {0.4, \(-2\), 0.38}, \ 
            ColorOutput -> GrayLevel, 
			\ 
            PlotRange -> {{\(-3\), 3}, {\(-3\), 3}, {\(-3\), 3}}, \ 
            BoxRatios -> {1, 1, 1}]}], \ ColorOutput -> GrayLevel]; 
    \)\)

.  The Spherical Cutoff

Now we turn to the manipulation of the polygons to get a nice-looking  spherical cutoff. We locate each polygon that intersects the sphere of radius 2.7 and fix it depending on how many of its four vertices lie outside the sphere. If 2 vertices lie inside and 2 outside, we just move the outside ones to the intersection of the protruding edge with the sphere. If three points lie outside, we do the same to each of the three points. The hardest case is where only one point lies outside the sphere; then we take edge intersections to get three new points. This gives a hexagon that we break into two quadrilaterals.

The intersection of PQ with a sphere of radius r is easily found by setting X to t P+(1-t)Q and solving X·X==(r^2). This yields two values of t, and we choose the one that lies in [0, 1]. The resulting formulas for the two values of t that work are not too bad and are encoded in the compiled functions f1 and f2 below. Then intersect finds the point of intersection by using the correct value of t.

\(f1\  = \ 
      Compile[{x, y, z, x1, y1, z1, r}, 
		
        \((2\ x\^2 - 2\ x\ x1 + 2\ y\^2 - 2\ y\ y1 + 2\ z\^2 - 2\ z\ z1 + 
              √\((
                  \((\(-2\)\ x\^2 + 2\ x\ x1 - 2\ y\^2 + 2\ y\ y1 - 2\ z\^2 + 
                          2\ z\ z1)\)\^2 - 
                    4\ \((\(-r\^2\) + x\^2 + y\^2 + z\^2)\)\ 
                      \((x\^2 - 2\ x\ x1 + x1\^2 + y\^2 - 2\ y\ y1 + y1\^2 + 
                          z\^2 - 2\ z\ z1 + z1\^2)\))\))\)/
          \((2\ \((
                x\^2 - 2\ x\ x1 + x1\^2 + y\^2 - 2\ y\ y1 + y1\^2 + z\^2 - 
                  2\ z\ z1 + z1\^2)\))\)]; 


    f2\  = \ Compile[{x, y, z, x1, y1, z1, r}, 
		
        \((2\ x\^2 - 2\ x\ x1 + 2\ y\^2 - 2\ y\ y1 + 2\ z\^2 - 2\ z\ z1 - 
              √\((
                  \((\(-2\)\ x\^2 + 2\ x\ x1 - 2\ y\^2 + 2\ y\ y1 - 2\ z\^2 + 
                          2\ z\ z1)\)\^2 - 
                    4\ \((\(-r\^2\) + x\^2 + y\^2 + z\^2)\)\ 
                      \((x\^2 - 2\ x\ x1 + x1\^2 + y\^2 - 2\ y\ y1 + y1\^2 + 
                          z\^2 - 2\ z\ z1 + z1\^2)\))\))\)/
          \((2\ \((
                x\^2 - 2\ x\ x1 + x1\^2 + y\^2 - 2\ y\ y1 + y1\^2 + z\^2 - 
                  2\ z\ z1 + z1\^2)\))\)]; 


    intersect[{{x_, y_, z_}, {u_, v_, w_}}, r_] := 
      Module[
		{t = f1[x, y, z, u, v, w, r]}, 
		
        If[\(! \((0 ≤ t ≤ 1)\)\), 
          t = f2[x, y, z, u, v, w, r]]; 
		
        \((1 - t)\)\ {x, y, z} + t\ {u, v, w}]\)

The next step is to fix a protruding polygon. The functions fix[i], i = 1, 2, 3 do that. The cutoff function cuts off a polygonal set of points according to which of the five cases it falls in (if 4 points lie outside, it is deleted, if no points lie outside it is retained, and the other cases are sent to fix).

\(\(fix[1]\)[pts_, \ r_]\  := 
      Module[{v, p, q, s}, 
\ \ \ \ \ 
        v\  = \ \(Select[pts, \ # . # > r\^2\ &, 1]\)[[1
            ]]; 	
		
        vv = \(Position[pts, v]\)[[1, 1
            ]]; 

        		{p, q, s}\  = \ 
          pts[[Mod[vv + {1, 2, 3}, 4]\  /. \ 
              0 -> 4]]; 

        \ \ 	{pv, qv, sv} = 
          \((\(intersect[{#, v}, \ r]&\)\ /@\ {p, q, s})\); 

        \ \ \ \ {{p, q, qv, pv}, {q, s, sv, qv}\ }]
		\)\(\(fix[2]\)[pts_, \ r_]\  := 
      Module[{p, q}, 
		{p, q}\  = \ Select[pts, \ # . # > r\^2&]; 

        \ \ \ {pp, qq} = 
          \(\(Position[pts, \ #]\)[[1, 1
                ]]&\)\ /@\ {p, q}; 
\ \ \ \ 
        pts\  /. 
          \ {p\  -> \ 
              intersect[{p, \ 
                  pts[[pp
                      \  /. {1 | 3 -> 
                          \((qq /. {2 -> 4, 4 -> 2})\), 

                        											\ \ \ \ \ \ \ 
                        2 -> \((qq /. {1 -> 3, 3 -> 1})\)}
                        ]]}, \ r], 
				\ 
            q\  -> \ 
              intersect[{q, \ 
                  pts[[qq\  /. \ 

                      										{
                        2 | 4 -> \((pp /. {1 -> 3, 3 -> 1})\), 
                        
											\ \ \ \ \ \ \ \ \ 
                        3 -> \((pp /. {2 -> 4, 4 -> 3})\)}
                        ]]}, \ r]}]
\)\(\(fix[3]\)[pts_, \ r_]\  := 
      Module[{v, p, q, s}, 
\ \ \ \ 
        v\  = \ \(Select[pts, \ # . # < r\^2&, 1]\)[[1
            ]]; 
\ \ \ \ 
        vv = \((\(Position[pts, v]\)[[1, 1
              ]])\); 

        		{p, q, s}\  = \ 
          pts[[Mod[vv + {1, 2, 3}, 4]\  /. \ 0 -> 4\ 
            ]]; 
\ \ 	
        pts\  /. 
          \ {p\  -> \ intersect[{p, v}, r], 

            				\ \ \ \ \ \ \ \ \ \ \ \ 
            q\  -> \ intersect[{q, v}, r], 

            			\ \ \ \ \ \ \ \ \ \ \ \ \ \ 
            s\  -> \ intersect[{s, v}, \ r]}]
\)\(cutoff[pts_, \ r_]\  := \ 
      \((h = howmany[pts, \ r]; 
		
        Which[h == 0, \ Polygon[pts], 
		\ \ \ h == 1, 
          Polygon\ /@\((\ \(fix[1]\)[pts, r])\), 
			h == 2\  || h == 3, 
          Polygon[\(fix[h]\)[pts, r]], 
	\ \ 	h == 4, \ {}])\)
\)SphericalCutoff[a_,  r_]  :=  
      a  /.  Polygon[pts_]  :>  cutoff[pts, r];

Code for Figure   DO NOT PRINT

\(\(
WireSphere[r_, \ m_, \ n_]\  := \ 
      Sphere[r, m, n]\ 	 /. 
        Polygon[{x_, y__}] :> Line[{x, y, x}]; 


    Show[Graphics3D[{{GrayLevel[0.2], WireSphere[1.2, 40, 40]}, 

          GrayLevel[0.95], 		
          Polygon[\((
                p = {{ .3,  .3, 0}, \ {1,  .3, 0}, \ {1, 1, 0}, 
                    \ { .3, 1, 0}})\)\  /. 
              \ {x_, y_, z_} :> \ {x, y, z -  .001}], 
			
          new = \ \ \ \(fix[1]\)[p, 1.2]; 
GrayLevel[ .3], 	
          Polygon[new[\([1]\)]], \ 	GrayLevel[0.7], 	
          Polygon[new[\([2]\)]]}], 

      PlotRange -> {{\(- .1\), 1}, \ {\(- .4\), 1.3}, \ {\(- .4\),  .5}}, \ 
      BoxRatios -> \ {1.1, 1.7,  .9}, 
Lighting -> False, \ 
      ViewPoint -> {1, \(-1.5\), 2}]; \)\)

And now we can see the spherical cutoff. This was the shape that we sculpted out of snow in Breckenridge, with all surfaces being about four inches thick.

finalPolys = {EdgeForm[{}], SphericalCutoff[CostaPolygons, rad]};

Show[Graphics3D[finalPolys], ViewPoint -> vp, 
      ViewVertical -> {0, -1, 0.6}];

Code for Figure

Show[
      Graphics3D[
        finalPolys = {EdgeForm[{}], SphericalCutoff[CostaPolygons, rad]}], 
      ViewPoint -> vp, ViewVertical -> {0, -1, 0.6},  
      ColorOutput -> GrayLevel];

.  Adding the Ends

An interesting enhancement is to find the points on the ends and add them to the surface. Certainly an argument ca be made that the borderless view in Figure  is preferable. But the point here is to show how to locate the ends and draw lines there if we wish. First we find the points on both Costa and the sphere (in three cases) and sort them for ease of curve-drawing.

\(angleSort[pts_]\  := \ 
      Sort[pts, 
        ArcTan[#1[[1]], 
              \ #1[[2]]] < 
            ArcTan[#2[[1]], 
              \ #2[[2]]]\ &]
\)\(threeEnds = 
      angleSort/@

        		{Cases[
            finalPolys, {x_?NumericQ, y_, z_} /; 
              Abs[{x, y, z} . {x, y, z} - rad\^2] ≤ 0.1 && z > 1, 
            ∞], 
          Cases[finalPolys, {x_?NumericQ, y_, z_} /; 
              Abs[{x, y, z} . {x, y, z} - rad\^2] ≤ 0.1 && 
                z < \(-1\), ∞], 
			
          Cases[finalPolys, {x_?NumericQ, y_, z_} /; 
              Abs[{x, y, z} . {x, y, z} - rad\^2] ≤ 0.1 && 
                Abs[z] < 0.9, ∞]}; 


    ends = Graphics3D[{Thickness[0.005], 
          \((Line[Append[#, First[#]]]&)\)/@\ threeEnds}]; \)

Next we sperimpose the ends on the sphere. This shows us that we have properly found and sorted the points. We are rewarded with a nice view of the baseball curve (it is reminiscent of the seams on a baseball) that forms the boundary of the central end.

\(WireSphere[r_, m_, n_] := 
      Sphere[r, m, n] /. Polygon[{x_, y__}] :> Line[{x, y, x}]; 

    
Show[Graphics3D[{GrayLevel[0.6], WireSphere[rad, 20, 20]}], 
			
      ends /. Thickness[_] :> Thickness[0.018], 
      BoxRatios -> {1, 1, 1}, ViewPoint -> {1, \(-2\), 0.1}, 
      PlotRange -> All]; \)

And we now put it all together, to get the desired spehrical cutoff of the Costa surface (Figure ). This brings a common problem to the fore: the lines forming the ends might be hidden from view by some polygons. To fix this we define a function, moveToEye, that moves a point 2% closer to the viewer’s eye. This function needs to know the viewpoint, box ratios, and plot range to do its job. We use the formula worked out by Steven Wilkinson in his article in this issue.

moveToEye[a_, vp_, br_, pr_, lf_:  0.02] := 
      a /. {x_?NumericQ, y_?NumericQ, z_?NumericQ} :> 
          (vpConverted =  
              Apply[Plus, pr, {1}]/2 - 
                vp   Apply[Subtract, pr, {1}] Max[br]/br; {x, y, z} + 
              ((lf (vpConverted - {x, y, z}))/(norm[
                  vpConverted - {x, y, z}] )))

\(\(Show[Graphics3D[{finalPolys}], 
	moveToEye[ends, vp, br, pr], 
      Prolog -> {{RGBColor[1, 1,  .4], Disk[{0.5, 0.5}, 0.31]}, 
          Circle[{0.5, 0.5}, 0.31]}, ViewVertical -> {0, \(-1\), 0.6}, 
      ViewPoint -> vp]; \)\)

Code for Figure

\(\(Show[Graphics3D[{finalPolys}], moveToEye[\ ends, \ vp, \ br, \ pr], \ 
      
	Prolog -> {{RGBColor[1, 1,  .4], Disk[{0.5, 0.5}, \ 0.31]}, \ 

          			\ Circle[{0.5, 0.5}, \ 0.31]}, 
	
      ViewVertical -> {0, \(-1\), 0.6}, \ \ ViewPoint -> vp, 
\ 
      PlotRegion -> {{\(-0.14\), \ 1.21}, \ {\(-0.16\), \ 1.11}}, \ 
      ColorOutput -> GrayLevel]; \)\)

Figure  shows a different projection of the surface. The two elliptical openings are orthogonal to each other and this view shows the horizontal opening.

\(\(Show[Graphics3D[{finalPolys}], 
	
      Prolog -> {{RGBColor[1, 1,  .4], Disk[{0.5, 0.5}, 0.31]}, 
          Circle[{0.5, 0.5}, 0.31]}, ViewVertical -> {0, \(-1\), 0.6}, 
      ViewPoint -> {\(-0.2\), \ \(-0.38\), \(-1\)}]; \)\)

Code for Figure

\(\(Show[Graphics3D[{finalPolys}], 
	
      Prolog -> {{RGBColor[1, 1,  .4], Disk[{0.5, 0.5}, 0.31]}, 
          Circle[{0.5, 0.5}, 0.31]}, ViewVertical -> {0, \(-1\), 0.6}, 
      ViewPoint -> {\(-0.2\), \ \(-0.38\), \(-1\)}, \ 
      ColorOutput -> GrayLevel]; \)\)

We conclude with two other views of the surface. First we cut it in half perpendicular to the baseball curve.

\(pr1 = {{0, 3.2}, \ {\(-3.2\), \ 3.2}, \ {\(-2.4\), \ 2.5}}; 

    br1 = {3.2, 6.4, 4.9}; 
vp1 = {\(-2\), \ \(-1.4\), \  .7}; 


    Show[Graphics3D[{Thickness[0.008], finalPolys\ }], 

      moveToEye[\ ends, \ vp1, \ br1, \ pr1], 
\ \ ViewPoint -> vp1, \ 
      BoxRatios -> br1, \ PlotRange -> pr1]; \)

Code for Figure  (do not print)

\(\(
pr2 = {{0, 3.2}, \ {\(-3.2\), \ 3.2}, \ {\(-2.4\), \ 2.5}}; 

    br2 = {3.2, 6.4, 4.9}; 
vp2 = {\(-2\), \ \(-1.4\), \  .7}; 
			
    Show[Graphics3D[{Thickness[0.008], 
			finalPolys\ }], 

      moveToEye[\ ends, \ vp2, \ br2, \ pr2], ColorOutput -> GrayLevel, 
\ \ 
      ViewPoint -> vp2, \ BoxRatios -> br2, \ PlotRange -> pr2]; \)\)

And next we slice it along the equator, thus revealing half of the baseball curve.

\(vp2 = {1.3, \(-2.4\), 2}; 

    pr2 = {{\(-3.2\), 3.2}, {\(-3.2\), 3.2}, {\(-2.4\), 0}}; 

    br2 = {6.4, 6.4, 2.4}; 

    Show[Graphics3D[{EdgeForm[], finalPolys}], 
	\ \ \ \ \ 
      moveToEye[ends, vp2, br2, pr2], 
	ViewPoint -> vp2, 
      BoxRatios -> br2, PlotRange -> pr2]; \)

Code for Figure

\(vp1 = {1.3, \(-2.4\), 2}; 

    pr1 = {{\(-3.2\), 3.2}, \ {\(-3.2\), \ 3.2}, \ {\(-2.4\), \ 0}}; 

    br1 = {6.4, 6.4, 2.4}; 
	

    Show[Graphics3D[{EdgeForm[], finalPolys}], 
	
      moveToEye[\ ends, \ vp1, \ br1, \ pr1], \ 
\ \ ViewPoint -> vp1, \ 
      BoxRatios -> br1, \ PlotRange -> pr1, \ ColorOutput -> GrayLevel]; \)

Figure  shows the amazing symmetry of the Costa surface: amazing because it was not at all suspected when the surface was first discovered. The symmetry group of the surface turns out to be (Subscript[D, 4]), the symmetry group of a square. One can take one quarter of the image in Figure  and use it to generate the entire figure, and then flip and rotate to get the rest. Knowing this allows one to speed up the computation of the polygons by about a factor of 8, by focusing just on the small patch forming one eighth of the unit square in the u-v domain, and then making the necessary transformations. We leave that as an exercise.

.  Toroidal Trousers

A pair of trousers is topologically the same as a sphere with three punctures. We mentioned earlier that the Costa surface is equivalent to a triply-punctured torus. Looking at Costa, it is by no means clear where the torus is hiding. And there is also the issue that the Invisible Handshake interpretation indicates that Costa has six entry/exit points: three for each hand. But, at first glance, a triply-punctured torus has fewer than six such points. This is an illusion, and by the end of this section the reader should be able to see where the six entry/exit points are in the punctured torus.

Helaman Ferguson described to us how to go from a torus to Costa. Consider the largest circle on the torus going around the entire surface and the smaller inside circle that goes around the visible hole. Stand the torus on end and insert a puncture at the very top and bottom of the large circle and a third puncture at the top of the small inner circle. The topmost hole is then stretched out, down, and around the torus to form the skirt (the baseball curve). The top inside puncture is stretched and pushed up to form the top ellipse and the bottom puncture is stretch out to form the bottom ellipse. This is difficult to picture at first and an animation can help a lot. On the other hand, figuring out formulas to do all the stretching seems like a near-impossible task. We shall see that, thanks to Al Gray’s parameterization, it is almost trivial! The resulting animation is shown on the cover of this issue.

Gray’s parametrization of the Costa surface is amazing. Watching the surface form as the parameters vary immediately gives a feeling for the
standard parametrization of a torus.  The Gray parametrization uses 0 ≤ u ≤ 1, 0 ≤ v  ≤ 1.  Figure  shows frames from a side-by-side movie of the Costa surface and the torus forming as the parameters vary across the square, shown in the middle.

The code for doing all this is in the electronic supplement. To get an idea of the main steps, here is a redefinition of the GrayCosta routine so that it caches the values, along with the definition of the torus. A similar modification must be made to intersect, since that is where new points on the surface are obtained.

\(\(torus = 
      Compile[{u, v}, {
			
          Sin[2\ π\ \((u + 0.5)\)]\ \((2 + Cos[v\ 2\ π])\), 
			
          Sin[\(-2\)\ π\ v], 
			
          \(-Cos[2\ π\ \((u + 0.5)\)]\)\ \((2 + Cos[2\ π\ v])\)}]; 
    \)\)

\(GrayCosta[u_, v_] := 
      \((z = u + v\ I; 
        ans\  = 
\ \ \ \ \ \ 
          With[{a = z + \(Z[z - 1\/2] - Z[z - I\/2]\)\/\(2  e\), 

              \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 
              b = Z[z]\/π - \(π\ \((1\  - \ I)\)\)\/\(4  e\), 

              \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 
              d = WeierstrassP[z, {c, 0}]}, 

\ \ \ \ \ \ \ 
            \(π\/2\) {Re[a - b], Im[a + b], 
                \(1\/\@\(2  π\)\) Log[Abs[\(d - e\)\/\(d + e\)]]}]; 
\ \ 
        torCache[ans] = torus[u, v]; \ 
\ \ uvCache[ans] = \ {u, v, 0}; 

        \ \ \ GrayCosta[u, v] = ans)\)\)

Code for Figure

The code is in the separate file for generating the cover images, since both require the same preliminaries

Notice the eight cut-outs near the edges of the square. These correspond to the skirt and ellipses in the Costa surface. What we do is define an inverse Costa map. After doing the spherical cutoff, the inverse Costa function maps the Costa surface back onto the square, but with the deleted polygons leading to missing pieces in the square. We get the inverse function simply by building an array indexed by the vertices of the Costa polygons. We build an array of torus values at the same time, thus getting the punctures in the torus that correspond exactly to the three ends in the spherically cutoff Costa surface.

The simplest approach to try to generate a movie showing the torus deforming into a Costa surface would be just a linear deformation of the form k torus + (1-k) costa. After watching the Costa surface parametrization and seeing how closely it resembles the familiar torus parametrization, this amazingly simple approach seemed like it might work. And indeed it does! The simple linear combination accomplishes the stretching and pulling as if the virtual torus were made from the proverbial topologist’s rubber. The final result (shown on the cover) makes it easy for anyone to see what mathematical sculptors naturally visualize in their imagination.

A final point: the initial punctured torus that arises in our method has a rather large puncture, so we added several frames that show a small puncture expanding. These were constructed simply by restricting v (via Min and Max) with two or three circular functions of u when u was within 0.03 of 0 or (1/2) or 1. The lighting was enhanced by using FaceForm to get different interior and exterior shadings, a point that is crucial to a visualization such as this.

.  Next Year

We plan to enter the snow sculpture competition again in January, 2000. Although we do not yet know what we will sculpt, we are sure that Mathematica will play a role in the development of the final vision, and the plan to get there.

REFERENCES

H. Ferguson, A. Gray, and S. Markvorsen, Costa’s minimal surface via Mathematica, Mathematica in Education and Research 5:1 (1996) 5–10.
H. Ferguson, C. Ferguson, T. Nemeth, D. Schwalbe, and S. Wagon, Invisible Handshake, The Mathematical Intelligencer (to appear).
D. Hoffman, The computer-aided discovery of new embedded minimal surfaces, The Mathematical Intelligencer 9:3 (1987) 8–21.
A. Ogawa, Helicatenoid, The Mathematica Journal 2:2 (1992) 21.

About the Authors

Dan Schwalbe has been a Mathematica user since coming to Macalester College in 1994.  Before that he co-authored the Maple Flight Manual and taught MAA Interactive Mathematics Text Workshops.  He is currently developing course materials with Mathematica and Java for an Introduction to Computational Science course.

Dan Schwalbe
Department of Mathematics and Computer Science
Macalester College
St. Paul, MN 55105
schwalbe@macalester.edu

Stan Wagon uses Mathematica extensively in his research, teaching, and exposition, and is especially appreciative of the new ways Mathematica allows us to look at mathematical objects. He is author or coauthor of several books, including the just-released second edition of Mathematica in Action, Animating Calculus, VisualDSolve (with Dan Schwalbe), and Which Way Did the Bicycle Go?. He also teaches a high-altitude course, Rocky Mountain Mathematica, every summer in the mountains of Colorado.

Stan Wagon
Department of Mathematics and Computer Science
Macalester College
St. Paul, MN 55105
wagon@macalester.edu