Wolfram Computation Meets Knowledge

Splines Come to Mathematica

One of the areas I contributed to Mathematica 7 was support for splines. The word “spline” originated from the term used by ship builders referring to thin wood pieces.

Over the last 40 years, splines have become very popular in computer graphics, computer animation and computer-aided design fields. From containers for household goods to state-of-the-art airplanes, it is hard to find any industrial product without spline surfaces. Also, they are widely used in other mathematical studies, such as interpolation and approximation.

Through its integration of numerics, symbolics and graphics, Mathematica has the opportunity to go much further with splines than has ever been possible before. Mathematica has had basic spline packages for a long time. But in Mathematica 7 we decided to make highly general spline support a core feature of the system.

Splines give another way to represent classes of functions. For decades, mathematicians had been using polynomials for numerical analysis. In the early 20th century, with advances in approximation theory, spline functions were beginning to emerge. The basic idea is simple. In essence, they consist of piecewise polynomials with local supports.

Since Version 5.1, Mathematica has offered general support for piecewise functions, both numerically and symbolically. In Mathematica 7, the B-spline functions can be expanded using PiecewiseExpand. For example, a uniform cubic B-spline basis function can be expanded to the following.

In[1]:= PiecewiseExpand[BSplineBasis[3,x]]

Splines have many attractive properties for interpolation. For instance, you can achieve good smoothness without relying on high degrees. Also, it is easier to avoid undesirable effects, such as Runge’s phenomenon. Ever since Mathematica 2, we’ve had the concept of InterpolatingFunction, which allows us to have a symbolic representation for an approximated function defined by interpolation. In Mathematica 7, the spline method is fully integrated in Interpolation. You can easily compare the results of the old "Hermite" method and the new "Spline" method.

In[2]:= data=Table[{i,RandomReal[]},{i,0,10}]; f=Interpolation[data,InterpolationOrder->3,Method->

In the late 1950s, two French mathematicians, Pierre Bézier from Renault and Paul de Casteljau from Citroën, independently developed spline-based computer-aided design systems for the automakers. Afterward, the popularity of splines literally exploded, and the rest is history.

Before we get into the fun side of splines, let’s discuss basis functions of B-splines, from which the B of B-splines comes. With interactive capability introduced in Mathematica 6, Mathematica is an incredible tool to explore the properties of the spline functions. For instance, a basis function can be created by adding of two lower-degree splines, as illustrated here using Manipulate.

In[6]:= f=BSplineBasis[{3,Range[5]},0,x]; Subscript[g, 0][d_]:=((x+d)-1)/3 BSplineBasis[{2,Range[5]},0,x+d]; Subscript[g,1][d_]:=(5-(x-d))/3 BSplineBasis[{2,Range[5]},1,x-d]; Manipulate[Plot[{f,Subscript[g,0][d],Subscript[g,1][d],Subscript[g,0][d]+Subscript[g,1][d]},{x,1,5},PlotStyle->{Automatic,Automatic,Automatic,Directive[Thick,Dashed]}],{d,0,1},SaveDefinitions->True]

In the preceding examples, each piecewise polynomial is defined on an interval with the same length and thus forms a “uniform” basis. However, by using intervals with different lengths, one can generate a family of nonuniform basis functions. The unions of the endpoints of intervals are called “knots.” The following example shows how the knot change affects the cubic basis functions.

In[10]:= Manipulate[Plot[Evaluate[Table[BSplineBasis[{3,Sort[knots[[All,1]]]},n,x],{n,0,5}]],{x,0,9},Filling->Axis,PlotRange->{All,{0,1}}],{{knots,Tuples[{Range[0,9],{0}}]},{0,0},{9,0},Locator}]

From these univariate functions, we can build up multidimensional functions using the tensor products. For instance, a basis function in 2D can be defined by multiplying two univariate functions.

In[12]:= basis[u_,v_]:=BSplineBasis[3,u]BSplineBasis[3,v]; Plot3D[basis[u,v],{u,0,1},{v,0,1}]

Or you can generate a family of 2D functions using different bases.

In[14]:= knots={0,0,0,1,1,1};GraphicsGrid[Table[Plot3D[BSplineBasis[{2,knots},i,u]BSplineBasis[{2,knots},j,v],{u,0,1},{v,0,1},PlotRange->All,Ticks->None],{i,0,2},{j,0,2}]]

The Wolfram Demonstrations Project already has a couple of interesting Demonstrations regarding B-spline basis. “B-Spline Basis Functions” is one such example.

B-Spline Basis Functions
B-Spline Basis Functions

Now, let’s move toward more interesting topics. We will start with curves. B-spline curves can be formulated as a linear summation of points and univariate basis functions.

curve[t_]:=\!\(\*UnderoverscriptBox[\(\[Sum]\),\(i=0\),\(n\)]\(pts[[i]]\BSplineBasis[{d,k},i-1,t]\)\)

And the formula can be further simplified by using Dot. The points are traditionally called “control points,” and they define the shape of the curve in the range space. This is a typical example of a cubic B-spline curve in 2D.

In[16]:= knots={0,0,0,0,1,1,1,1}; pts={{0,0},{1,2},{3,0},{4,2}}; c[t_]:=Dot[Table[BSplineBasis[{3,knots},i,t],{i,0,3}],pts]; g1=ParametricPlot[c[t],{t,0,1}];

But all those indices and dot products can be frustrating. As an alternative, Mathematica 7 supplies a new function, called BSplineFunction, which makes the construction much simpler.

In[20]:= f=BSplineFunction[pts,SplineDegree->3,SplineKnots->knots];<br /> g2=ParametricPlot[f[t],{t,0,1}];

It is not just about the convenience. Internally, when BSplineFunction is called, Mathematica prepares important parameters that make overall computation much faster in subsequent calls.

But we haven’t stopped here. Mathematica 7 also provides new B-spline graphics primitives, such as BSplineCurve and BSplineSurface. The preceding example can be written purely in graphics language.

In[22]:= g3=Graphics[BSplineCurve[pts,SplineDegree->3,SplineKnots->knots],Axes->True];

And you can see that all three approaches result in identical curves.

In[23]:= GraphicsRow[{g1,g2,g3}]

The new features in Mathematica provide very detailed controls over B-splines. We can deal with arbitrary degrees and arbitrary knot configurations. I made a curve editor based on a Demonstration I wrote before Mathematica 7. This example allows you to control every aspect of B-spline curves, as well as the ability to inspect their basis functions.

B-Spline Curve Editor with Knot Control
B-Spline Curve Editor with Knot Control—Click to download notebook

In addition, SplineWeights are provided for rational B-splines, or NURBS. With rational splines, one can create shapes that are not possible to represent with polynomial—for instance, a perfect circle:

In[24]:= pts={{.5,0},{1,0},{1,1},{.5,1},{0,1},{0,0},{.5,0}}; w={1,.5,.5,1,.5,.5,1}; k={0,0,0,1/4,1/2,1/2,3/4,1,1,1}; Graphics[{BSplineCurve[pts,SplineWeights->w,SplineKnots->k,SplineClosed->False],Gray,Dashed,Line[pts],Red,Point[pts]}]

B-spline surfaces can be created in the same fashion, except that now we need to use tensor products. Instead of a list of points, we need to use a matrix of points, which is sometimes called a control mesh. The following example shows a simple B-spline surface using BSplineSurface.

In[28]:= pts={{{0,0,2},{1,0,1},{3,0,2},{4,0,0}},{{0,1,3},{1,1,2},{3,1,3},{4,1,1}},{{0,2,3},{1,2,2},{3,2,3},{4,2,1}},{{0,3,2},{1,3,1},{3,3,2},{4,3,0}}}; Graphics3D[{Yellow,Specularity[White,50],BSplineSurface[pts],Red,PointSize[Medium],Point/@pts,Dashed,Blue,Line[pts],Line[Transpose[pts]]}]

Of course, one can use BSplineFunction with ParametricPlot3D instead. Naturally, all the options for ParametricPlot3D work perfectly well with B-spline surfaces.

In[30]:= ParametricPlot3D[BSplineFunction[pts][u,v],{u,0,1},{v,0,1},BoundaryStyle->Black,ColorFunction->(Hue[ArcTan[#4-1/2,#5-1/2]/(2 \[Pi])]&),Mesh->5,RegionFunction->((#4-1/2)^2+(#5-1/2)^2>1/8+1/16 Sin[10 ArcTan[#4-1/2,#5-1/2]]&)]

As I mentioned, splines surfaces have been used in car design for almost half a century. I wrote this little Demonstration to illustrate how a boxy design can be transformed into a futuristic body, by increasing spline degrees. Of course, this Demonstration is just for fun, and does not reflect real-world car design (which is an incredibly elaborate task, by the way).

Designing a Car Body with Splines
Designing a Car Body with Splines

The new spline functionality in Mathematica 7 is very rich and general. For instance, one can easily make a 4D spline manifold (or any number of dimensions, for that matter) using BSplineFunction.

In[31]:= f=BSplineFunction[RandomReal[1,{7,7,7}],3]; ContourPlot3D[f[x,y,z],{x,0,1},{y,0,1},{z,0,1},Mesh->None,MaxRecursion->1,Contours->5,ContourStyle->ColorData[27,

I will continue the story about splines in Mathematica 7 in a later post. I would like to discuss a few Demonstrations that I wrote specifically for the B-splines, and some other features that I haven’t mentioned here. It will be more interesting (and practical), so please stay tuned.