B-spline Curve and Surface Evaluation

Abstract

Using the Taylor expansion, a B-spline segment can be efficiently represented in the power basis. Consequently, we can use the Horner algorithm to evaluate points and derivatives with respect to a given set of parameters. Our approach is a generalization of W. Böhm's algorithm.

1  Introduction

A vector-valued rational curve may be represented as R(t)=[(P(t))/W(t)], where P:®Âd and W:® are polynomials. Using the Leibnitz formula, we can derive that the mth derivative of R(t) is given by

R(m)(t)=
P(m)(t)- m
å
 i=1 
 CmiR(m-i)(t)W(i)(t)
W(t)
 ,
where, Cmi=[m!/(i!(m-i)!)]. Similarly, differentiating a rational surface, S(u,v)=[(P(u,v))/W(u,v)], m times in u and l times in v yields:

S(m+l)umvl(u,v)=
P(m+l)umvl(u,v)- (m,l)
å
 (i,j) ¹ (0,0) 
 CljCmiSum-ivl-j(l-j+m-i)(u,v) W(j+i)uivj(u,v)
W(u,v)
 .
As it is seen, the evaluation of a rational curve/surface is actually based on evaluation of a polynomial curve/surface. Therefore, we limit our discussion to polynomial curves and surfaces in this paper.

A piecewise polynomial B-spline curve of order k (or degree n=k-1) may be given by

P(t)=
å
 i 
 di Ni,k(t),
(1-1)
where, di Î Âd are the de Boor (control) points and Ni,k(t) the B-spline basis functions.3 In CAGD applications, we often need to evaluate points and derivatives on a certain B-spline segment with respect to a given set of parameters. In this case, we may represent this segment in the power basis and evaluate it efficiently using the Horner algorithm.

One way to represent a B-spline segment in the power basis is to convert the segment first into the Bézier basis using knot insertion technique, then into the power basis. Böhm proposed an alternative approach in which he computes the Taylor expansion of a B-spline segment.2 His approach is more efficient than the first one. In this paper, we extend his algorithm to B-spline surfaces and give details for implementation.

2  Taylor expansion of B-spline curves

A Taylor expansion of a B-spline segment was proposed by Böhm.2 For completeness, we briefly review his method in this section and give our pseudo-code that minimizes the use of computer memory.

Let d0,0i=di with di being the de Boor points in (1-1). Then, the mth derivative of P(t) is given by:4,5

P(m)(t)=(k-1)(k-2)¼(k-m)
å
 i 
 dim,mNi,k-m(t)
where,

dim,m= dim-1,m-1-di-1m-1,m-1
ti+k-m-ti
 = dim-1,m-1-di-1m-1,m-1
gim
     (m=1,¼,n).
(2-1)
Using the de Boor's knot insertion method, we may represent P(m)(t) in a lower degree form:

P(m)(t)=(k-1)(k-2)¼(k-m)
å
 i 
 dip,mNi,k-p(t)
where,2,3

dp,mi=dp-1,mi-1+(t-ti)dp,m+1i = dp-1,mi-1+bidp,m+1i     (p=m+1,¼,n).
(2-2)
For t Î (tj,tj+1), Ni,1(t)=di,j (di,j=1 if i=j, otherwise di,j=0). Therefore, when p=n, we have

P(m)(t)=(k-1)(k-2)¼(k-m)djn,m    (m=0,1,¼, n).
Accordingly, the evaluation of a point and m (m £ n) derivatives at t Î (tj,tj+1) reduces to computing dn,0j,dn,1j,¼,dn,mj. When m > n, P(m)(t)=0 since P(t) is a polynomial curve of degree n.

For t* Î (tj,tj+1), it is known Ni,k(t*) > 0 if i=j-n,j-n+1,¼,j and Ni,k(t*)=0 elsewhere. Therefore, to evaluate a point and derivatives at t* Î (tj,tj+1) requires:

In consideration of the above facts, we present the following pseudo-code which minimizes the number of divisions and the amount of memory space used:

(1) Compute gpi and bp:
       Let [^t]i=tj-n+i, where i=1,¼,2×n.
       do p=1, n
             do i=p, n
                   gpi = 1 / ([^t]i+k-p -[^t]i)
             enddo
             bp=t*-[^t]p
       enddo
 
(2) Compute dm,mi using equation (2-1):
       Let d0,0i=dj-n+i, where, i=0,1,¼,n.
       do m=1, n
             do i=m, n
                   dim,m=gmi(dim-1,m-1 -dm-1,m-1i-1)
             enddo
       enddo
 

(3) From dm,mj-n+m we compute dp,mj-n+p using equation (2-2):
       do p=1, n
             do m=p-1, 0, -1
                   dp,mm=dp-1,mm +bpdp,m+1m+1
             enddo
       enddo

For di Î Âd, it is seen that the computation of dn,mj (m=0,1,¼,n) requires [(n(n+1))/2] divisions and d×n×(n+1) multiplications. (Evaluation of a point using the knot insertion algorithm requires [(n(n+1))/2] divisions and [(d×n×(n+1))/2] multiplications.) After obtaining dn,ij at t* Î (tj,tj+1), the Taylor expansion of the B-spline segment corresponding to t Î [tj,tj+1] is given by
P(t)= n
å
 i=0 
 P(i)(t*)
i!
 (t-t*)i = n
å
 i=0 
 dn,ij
Cni
 (t-t*)i.
We have assumed that t* Î (tj,tj+1). However, if the knot sequence has a multiplicity r with respect to t*, we may choose j such that t* Î (tj,tj+1] or t* Î [tj,tj+1). This depends on whether the ``left'' or ``right'' segment is evaluated.

3  Taylor expansion of B-spline surfaces

A B-spline surface of order ku in u and kv in v is given by
S(u,v)=
å
 i 
å
 j 
 Ni,ku(u)Nj,kv(v)di,j.
We may write a B-spline surface as
S(u,v)=
å
 j 
 Nj,kv(v)dj(u),
(3-1)
where
dj(u)=
å
 i 
 Ni,ku(u)di,j.
Assume (u*,v*) Î (uju,uju+1)×(vjv,ujv+1). By considering (3-1) as a B-spline curve with parameter v, we can evaluate S(l)vl(u,v) (l=0,1,¼,nv) at v* using the same approach for curve and obtain:

    S(0)(u,v*)=djvnv,0(u) = åi Ni,ku(u)di,jvnv,0

    S(1)v(u,v*)=(kv-1)djvnv,1(u) = åi Ni,ku(u)di,jvnv,1

    ¼

    S(l)vl(u,v*)=(kv-1)¼(kv-l)djvnv,l(u) = åi Ni,ku(u)di,jvnv,l

Further evaluation of S(l)vl(u,v*) at u* Î (uju,uju+1) requires only ku de Boor points whose i indices range from ju-nu to ju. Therefore, we should implement the algorithm so that it minimizes the use of memory space and requires only [(kv×(kv-1))/2] divisions and d×kv×(kv-1)×ku multiplications to compute all S(l)vl(u,v*). Having obtained S(l)vl(u,v*), we can further evaluate S(m+l)umvl(u,v*) at u*, which is equivalent to evaluating kv B-spline curves.

The above approach evaluates S(u,v) first in v then in u. However, it should be noted that if ku < kv, we will save computation time by evaluating S(u,v) first in u then in v. Since the alternative approach is symmetric, we shall not discuss it again.

Given all partial derivatives of a polynomial surface S(u,v) at (u*,v*), we can represent the surface in the power basis by the Taylor expansion:

S(u,v)=S(u*,v*)+ æ
ç
è 		x
u
 + y
v
 ö
÷
ø 		S(u*,v*)+¼+ 1
n!
 æ
ç
è 		x
u
 +y
v
 ö
÷
ø 		n

 
 S(u*,v*)
where n=nu+nv, x=u-u*, y=v-v*, and

æ
ç
è 		x
u
 +y
v
 ö
÷
ø 		k

 
 S(u*,v*)= k
å
 i=0 
 Cki k S(u*,v*)
uk-ivi
 xk-iyi   , k=1,2,¼,n
For a polynomial surface, it is noted that, if k-i > nu or i > nv,

k S(u*,v*)
uk-ivi
 =0.
Therefore,

1
k!
 æ
ç
è 		x
u
 +y
v
 ö
÷
ø 		k

 
 S= k
å
 i=0 
 Cki
k!
 k S
uk-ivi
 xk-iyi = min{k,nv}
å
 i=max{0,k-nu} 
 Cki
k!
 k S
uk-ivi
 xk-iyi.
To speed up computation, we may create a coefficient table for [(Cki)/k!]=[1/(i!(k-i)!)].

4  Horner algorithm

Given a cubic curve P(t)=a0+a1t+a2t2+a3t3 with ai Î Âd, we may evaluate the curve efficiently using the well-known Horner method:

P(t) = t[t(a3t+a2)+a1]+a0.
As it is seen, only three multiplications are involved. To evaluate the first derivative, we may compute first P¢(t)=a1+2a2t+3a3t2, then the value of the first derivative using the Horner method. In this way, it is readily checked that d(3n-2) multiplications are required to compute a point and the first derivative. An improved algorithm in computing a point and the first derivative can be found in Acton's book.1 Writing the first derivative as

P¢(t)=t[a3t+(a3t+a2)]+t(a3t+a2)+a1,
reveals the relation between the evaluation of a polynomial and its first derivative. This is summarized in the following pseudo-code: 

   point = a(n) * t + a(n-1)
   deriv = a(n)
   do j = n-2, 0, -1
      deriv = deriv * t + pt
      point = point * t + a(j)
   enddo
It is seen that only d(2n-1) number of multiplications involved. Let eval(j) denote the jth derivative. We extend the improved algorithm to the evaluation of m derivatives. 
   eval(0) = a(n) * t + a(n-1)
   eval(1) = a(n)
   do i=n-2, 0 , step=-1
      k = min(m,n-i)
      do j=k, 1, step=-1
         eval(j) = eval(j) * t + eval(j-1)
      enddo
   enddo
c  multiply derivatives by factorial constant:
   k = 1
   do i=2, m
      k = k * i
      eval(i) = eval(i) * k
   enddo
It is readily checked that only d([(n(n+1))/2] + 3(n-1)) multiplications are required to evaluate a point and all n derivatives.

We can extend the above method to a surface. Given a vector-valued surface

S(u,v)= m
å
 j=0 
 n
å
 i=0 
 aijuivj,
we may represent it in the following form:

S(u,v)= m
å
 j=0 
 aj(u) vj,
(4-1)
where aj(u)=åi=0naijui. By considering (4-1) as a polynomial curve in v, we can use the algorithm in evaluating a curve to obtain S(u,v),Sv(1)(u,v),¼. Then, with respect to S(u,v),S(1)v(u,v) and so on, we compute derivatives in terms of u.

5  Conclusions

In CAD/CAM systems, curves and surfaces are usually represented in the B-spline form. Evaluation of a B-spline curve/surface may be carried out using the knot insertion technique. However, if we want to evaluate a set of points and derivatives on certain B-spline segments, algorithms based on knot insertion are very costly. In this case, a better solution is to represent these B-spline segments in the power basis using the Taylor expansion. Then, evaluation of a curve/surface can be done efficiently using the Horner method.

6  References

[1] Acton, Forman S.: Numerical Methods That Work, Harper and Row, New York, 1970.

[2] Böhm, W.: Efficient evaluation of splines, Computing (33), 171-177 (1984)

[3] de Boor, C.: A Practical Guide to B-spline, Springer-Verlag, New York, 1978.

[4] Hoscheck, J. and Lasser, D.: Fundamentals of Computer Aided Geometric Design, A K Peter, Ltd, 1993.

[5] Farin, G.: Curves and Surfaces for Computer Aided Geometric Design, Academic Press, Inc, 1988.

RETURN