UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Electroreduction of oxygen to hydrogen peroxide of particulate electrodes Oloman, Colin 1974

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
831-UBC_1974_A6_7 M47.pdf [ 5.39MB ]
Metadata
JSON: 831-1.0058928.json
JSON-LD: 831-1.0058928-ld.json
RDF/XML (Pretty): 831-1.0058928-rdf.xml
RDF/JSON: 831-1.0058928-rdf.json
Turtle: 831-1.0058928-turtle.txt
N-Triples: 831-1.0058928-rdf-ntriples.txt
Original Record: 831-1.0058928-source.json
Full Text
831-1.0058928-fulltext.txt
Citation
831-1.0058928.ris

Full Text

INTERACTIVE SPLINE APPROXIMATION  by  MARIAN MERCHANT B.Sc,  Simon F r a s e r U n i v e r s i t y , 1971  A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF  Master o f S c i e n c e i n t h e Department of Computer  Science  We accept t h i s t h e s i s as conforming tp t h e r e q u i r e d s t a n d a r d  THE UNIVERSITY  OF BRITISH COLUMBIA  January, 1974  In  presenting  an  advanced  the I  Library  further  for  his  of  this  written  thesis  degree shall  agree  scholarly  by  this  at  purposes  for  it  freely  permission  representatives. thesis  may It  Depa rtment  The U n i v e r s i t y o f B r V a n c o u v e r 8, Canada  is  fulfilment  of  of  Columbia,  British  available  for  for extensive  be g r a n t e d  financial  permission.  partial  the U n i v e r s i t y  make  that  in  by  shall  that  not  requirements I  agree  r e f e r e n c e and copying  t h e Head o f  understood  gain  the  of  this  be a l l o w e d  or  that  study. thesis  my D e p a r t m e n t  copying  for  or  publication  without  my  - i i ABSTRACT  The approximation  use o f s p l i n e b a s i s f u n c t i o n s i n s o l v i n g l e a s t problems i s i n v e s t i g a t e d .  The  squares  q u e s t i o n as to which a r e  a p p r o p r i a t e b a s i s f u n c t i o n s to use i s d i s c u s s e d a l o n g w i t h the why  the f i n a l c h o i c e was  made.  The Householder t r a n s f o r m a t i o n method  f o r s o l v i n g the f i x e d knot s p l i n e approximation D e s c r i p t i o n s of both an automatic and  reasons  problem i s examined.  procedure u s i n g f u n c t i o n m i n i m i z a t i o n  an i n t e r a c t i v e procedure u s i n g a g r a p h i c s t e r m i n a l f o r s o l v i n g the  v a r i a b l e knot s p l i n e approximation numerical  problem a r e g i v e n .  In c o n c l u s i o n ,  r e s u l t s u s i n g the i n t e r a c t i v e system are presented  and  analyzed.  - iii  -  TABLE OF CONTENTS  Section  1  S p l i n e s i n I n t e r a c t i v e Approximation  1  Section  2  The  3  2.1  Introduction  3  2.2  D e f i n i t i o n of a S p l i n e F u n c t i o n  4  2.3  P i e c e w i s e Continuous P o l y n o m i a l S p l i n e Representation  5  2.4  Mathematical S p l i n e Representation  6  2.5  B-Spline  9  2.6  Example o f B - S p l i n e  2.7  Equivalence  3  The  3.1  Introduction  18  3.2  D e f i n i t i o n o f the L i n e a r Problem  19  3.3  S o l u t i o n of the L i n e a r Problem  20  3.4  S o l u t i o n of the F i x e d Knot Problem  2 7  4  The V a r i a b l e Knot Approximation Problem  30  4.1  Introduction  30  4.2  D e f i n i t i o n of the V a r i a b l e Knot Problem  31  4.3  S o l u t i o n o f the V a r i a b l e Knot Problem  32  4.4  Knot O p t i m i z a t i o n  34  5  Numerical Results  37  5.1  Introduction  37  5.2  Use  39  Section  Section  Section  Spline Representation  Problem  Representation Representation  of Spline Representations  L i n e a r Approximation Problem  of the System  12 14  18  - iv -  Section  5.3  T i t a n i u m Heat Data  52  5.4  The Bug  62  6  Summary and F u t u r e P o s s i b i l i t i e s  Bibliography Appendix  669,  72  A  Program L i s t i n g s  B  U s e r s ' Guide  74 102  - v TABLES  I  T e s t Data - 2 Non-uniform Knots  44  II  T e s t Data -,2 Knots Optimized  45  III  T e s t Data - 3 U n i f o r m Knots  47  IV  T e s t Data - 4 Non-uniform Knots  50  V  T i t a n i u m Heat Data - 5 U n i f o r m Knots  53  VI  T i t a n i u m Heat Data - 5 Non-uniform Knots  56  VII  T i t a n i u m Heat Data - 5 Knots Optimized  59  VIII  The Bug  64  - vi -  FIGURES  1  Elementary C u b i c S p l i n e B a s i s F u n c t i o n s  8  2.  Cubic B - s p l i n e B a s i s F u n c t i o n s  15  3.  T e s t Data - 2 U n i f o r m Knots  43  4.  T e s t Data - 2 Knots Optimized  46  5.  T e s t D a t a - 3 U n i f o r m Knots  48  6.  T e s t Data - 4 U n i f o r m Knots  49  7.  T e s t .Data - 4 Non-uniform Knots  51  8.  T i t a n i u m Heat Data - 5 U n i f o r m Knots  55  9.  T i t a n i u m Heat Data - 5 Non-uniform Knots  58  10.  T i t a n i u m Heat Data - 5 Knots Optimized  61  11.  The Dented Bug  6  12.  The Bug  68  3  - vii-  ACKNOWLEDGEMENT  There a r e many p e o p l e who a i d e d i n t h e development o f t h i s thesis,  S i n c e I cannot hope t o i n c l u d e a l l i n d i v i d u a l s who p r o v i d e d  h e l p ; I s h a l l mention only; those who p r o v i d e d major a s s i s t a n c e and a p o l o g i z e t o anyone who may f e e l s l i g h t e d by b e i n g  omitted.  I would l i k e t o thank Dr. J . M. V a r a h who p r o v i d e d a s s i s t a n c e through our. d i s c u s s i o n s o f t h e t o p i c . the f i n a n c i a l a i d which he p r o v i d e d  academic  A l s o I am g r a t e f u l f o r  i n t h e form o f a r e a s e a r c h a s s i s t a n t s h i p .  I would l i k e t o thank t h e Department o f Computer S c i e n c e and t h e Computing Centre, f o r t h e a s s i s t a n t s h i p s they p r o v i d e d .  In p a r t i c u l a r I  a p p r e c i a t e d t h e IBM f e l l o w s h i p w h i c h I r e c e i v e d f o r r e s e a r c h i n I n t e r a c t i v e Numerical A n a l y s i s .  Finally,  This provided  the i n i t i a l  ideas f o r t h i s  thesis.  I would l i k e t o thank my husband P e t e r f o r a l l t h e  d i s h e s he d i d d u r i n g the w r i t i n g o f t h i s  thesis.  - viii  -  NOTATION  The  f o l l o w i n g common n o t a t i o n  [n],  [n,pp.  x €  [a,b]  ]  i s used throughout the t h e s i s :  r e f e r s to reference  C(a,b)  a < x < ,b;  s  CC [a,b]  s  s(x) dx  i n the b i b l i o g r a p h y ;  a <_ x <^ b;  x  m  n  has  i s the  m  continuous d e r i v a t i v e s on  j - t h d e r i v a t i v e of  s  [a,b];  with respect  to  x  J  evaluated  a) (x) r  {(?£,y£):£=1,2,...,n} £i  6;  i s t h e f i r s t d e r i v a t i v e o f the f u n c t i o n respect  A = {6.:i=l,2,...,k}  at  to  to w i t h  x;  i s the s e t A  w i t h elements  i s the s e t o f orderd. p a i r s i s the v e c t o r  a;  6,6,...,6 ; (x^,y^);  that i s ,  a =  m i s the l e a s t squares norm o f aj min a.  f (a^)  max {c,d}  i s t h e minimum o f the f u n c t i o n  f  with respect to  a; i s the maximum o f  c  and  d;  i s the t r a n s p o s e o f the m a t r i x  Q .  - 1 -  Section 1  SPLINES IN INTERACTIVE APPROXIMATION  P o l y n o m i a l s a r e f r e q u e n t l y d e s i r e d as a s e t o f b a s i s f u n c t i o n s for approximation.  Problems i n o b t a i n i n g a c c u r a t e r e s u l t s  2 standard set of b a s i s f u n c t i o n s  with,the  IB.  {1  }  l e a d to the development  of o r t h o g o n a l p o l y n o m i a l s as the s e t o f b a s i s f u n c t i o n s .  The use o f  o r t h o g o n a l p o l y n o m i a l s i s a w e l l - e s t a b l i s h e d technique i n a p p r o x i m a t i o n . Consequently,  they have been s t u d i e d i n depth and a r e not d i s c u s s e d  further here. P o l y n o m i a l s do have one disadvantage is,  i n approximation.  t h e i r n a t u r e over the e n t i r e r e g i o n o f approximation  by t h e i r b e h a v i o r i n o n l y a s m a l l area o f t h i s r e g i o n .  That  i s determined Higher-order  p o l y n o m i a l s do not a l l e v i a t e t h i s problem but merely impose more o s c i l l a t o r y b e h a v i o r on the a p p r o x i m a t i o n .  However, p o l y n o m i a l s p l i n e s  can  counteract t h i s r e s t r i c t i v e nature of polynomials. P o l y n o m i a l s p l i n e s , c o n s i s t of p i e c e w i s e p o l y n o m i a l s at p o i n t s known as knots over the r e g i o n o f a p p r o x i m a t i o n . p o l y n o m i a l determines  the shape o f the a p p r o x i m a t i o n  r e l a t i v e l y i n d e p e n d e n t l y of s u r r o u n d i n g a r e a s . i s determined  by imposing  i n a small area  a t the k n o t s .  o f these knots t h a t makes s p l i n e s  i d e a l l y s u i t e d f o r i n t e r a c t i v e approximation. cedures  Each p i e c e w i s e  The amount of dependence  c o n t i n u i t y requirements  I t i s p r e c i s e l y the involvement  connected  P r e v i o u s l y , automatic  i n v o l v i n g the m i n i m i z a t i o n o f the l e a s t squares  pro-  e r r o r were used  - 2 to  determine the b e s t l o c a t i o n s f o r the k n o t s .  linear ly,  appropriate i n i t i a l  As t h e knots occur non-  guesses had to be made and t h e r e e x i s t s t h e  p o s s i b i l i t y t h a t t h e procedure would converge t o a l o c a l minimum r a t h e r than the g l o b a l minimum.  Such a l o c a l minimum might n o t have been a  s u i t a b l e approximation.  G r a p h i c a l i n t e r a c t i o n a l l o w s immediate procedure.  c o n t a c t w i t h the s o l u t i o n  The l o c a t i o n s o f t h e knots can be chosen v i s u a l l y and t h e  a p p r o x i m a t i o n attempted w i t h t h i s s e t can be observed. procedure enables poor i n i t i a l eliminated.  e s t i m a t e s f o r t h e knot l o c a t i o n s t o be  A l s o i t a l l o w s f o r the m a n i p u l a t i o n o f t h e knots u n t i l a  s a t i s f a c t o r y approximation i s obtained. the  The i n t e r a c t i v e  i n t e r a c t i v e procedure as i n i t i a l  U s i n g the knots o b t a i n e d from  guesses, an automatic procedure should  converge r a p i d l y to a s u i t a b l e minimum.  - 3 -  THE SPLINE REPRESENTATION PROBLEM  2.1  Introduction  There i s no u n i v e r s a l l y accepted basis functions.  representation f o r spline  There a r e s e v e r a l k n o w n ^ r e p r e s e n t a t i o n s each w i t h i t s  own p a r t i c u l a r advantages and d i s a d v a n t a g e s .  Carasso and L a u r e n t 14] d i s c u s s t h r e e methods o f n u m e r i c a l  con-  s t r u c t i o n o f s p l i n e s - a p r o j e c t i o n method, a method o f d i r e c t r e s o l u t i o n and  a method u s i n g a b a s i s .  involving a basis.  Of t h e s e ,  they recommend  With the c h o i c e o f a r e a s o n a b l e  the use of a method  b a s i s , Carasso and  Laurent conclude t h a t t h i s method p r o v i d e s more a c c u r a t e . r e s u l t s than the p r o j e c t i o n method and t h r e e times l e s s computation than the method o f direct resolution.  G r e v i l l e [8] p r o v i d e s for splines.  truncated  one u s i n g B ^ s p l i n e s  power f u n c t i o n s  continuous p o l y n o m i a l s  Schultz Schultz describes  (discussed i n Section  ( d i s c u s s e d i n S e c t i o n 2.5).  summarize these r e p r e s e n t a t i o n s piecewise  functions  From the d e f i n i t i o n o f a s p l i n e f u n c t i o n , he develops a  representation-using and  a comprehensive overview o f b a s i s  2.4)  de Boor and R i c e [5]  and a l s o i n c l u d e a r e p r e s e n t a t i o n i n v o l v i n g ( d i s c u s s e d i n S e c t i o n 2.3).  [17] g i v e s a g e n e r a l b a s i s f o r B^-splines. the r e p r e s e n t a t i o n  f o r cubic B-splines  In [ 1 8 ] ,  i n more d e t a i l .  The b a s i s f u n c t i o n r e s u l t i n g from a p p l y i n g t h e s e t of c u b i c B-rsplines t o the s p e c i a l case of u n i f o r m l y  spaced knots i s s t a t e d .  t h i s r e s u l t i s given i n d e t a i l i n Section  2.6.  The d e r i v a t i o n o f  - 4 -  2,2  D e f i n i t i o n , of a S p l i n e Function  A l t h o u g h s p l i n e s e x i s t i n e n g i n e e r i n g and d r a f t i n g as a d e v i c e f o r curve smoothing, the'.basic mathematical  formulationcof a s p l i n e  f u n c t i o n comes from p i e c e w i s e c o n t i n u o u s p o l y n o m i a l s .  The  mathematical  d e f i n i t i o n f o r m a l i z e s the e n g i n e e r i n g concept.  D e f i n i t i o n 1;'  A  ( p o l y n o m i a l ) s p l i n e f u n c t i o n o f degree  i s a p o l y n o m i a l o f degree,  m  m  which Is i n  on C  m-1  [a,b] [a,b] .  A l t h o u g h t h i s d e f i n i t i o n i n c o r p o r a t e s the b a s i c n o t i o n o f a s p l i n e f u n c t i o n , i t does not p r o v i d e the e s s e n t i a l components needed f o r the use o f s p l i n e s i n n u m e r i c a l p r o b l e m - s o l v i n g .  For t h i s purpose,  the  f o l l o w i n g , more c o n s t r u c t i v e d e f i n i t i o n o f s p l i n e s i s b e t t e r .  D e f i n i t i o n ,2;' —•• ——  Given a p a r t i t i o n  a =  8  n  0  < 6_ < 1  ....'< 6.. < 6\ ,. = b k k+1  then a (polynomial) s p l i n e f u n c t i o n o f degree k  i n t e r n a l knots  S-^j  function  S(x)  w i t h the f o l l o w i n g  1.  S(x)  i s a p o l y n o m i a l of degree  2,  S(x)  and  A = {S^  ; i = 0,  S(x) = ( s ^ ( x ) ; i = 0, 1, s.(x) x  is in  [6., 6..-] i ' x+1  n  f »^-' a  ^  with s  a  properties; m  or l e s s i n  1, 2,  m^l  everywhere.  1,.  ..., k} for  °  i t s d e r i v a t i v e s of o r d e r s  are continuous  Let  ' '"' ^k  m  k+1}  be the s e t o f knots  be the s e t o f p o l y n o m i a l s such  i = 0,1,  • '  k .  The s e t  s a t i s f y the f o l l o w i n g c o n t i n u i t y c o n d i t i o n s a t the k n o t s :  S(x)  and that must  - 5 -  =  7j  -  Sl  . .."k;  l(x)  — .  dx  6. l  j = 0, 1,  3  s. 1  (x) '6\  l  for  i = 1, 2,  ..., m - 1 .  2.3  P i e c e w i s e Continuous P o l y n o m i a l S p l i n e R e p r e s e n t a t i o n  The p i e c e w i s e continuous p o l y n o m i a l d e f i n i t i o n f o r m u l a t e d i n t o an approximation problem. [6\,6^ ^]'+  for  i = 0, 1,  nomial of degree  m  or l e s s .  {c_^_. : i = 0, 1, c. . ' s  k;  k;  (2) can be  Suppose t h a t on each  the d a t a i s approximated  interval  by a  poly-  Given the s e t of c o e f f i c i e n t s  j =1,  1,  .. . , m}  the problem  i s to f i n d  the  where m  S(x) for  6 . < x < <$.,.. .1 — — l+l  =  where  i = 0, 1,  T h i s system r e s u l t s i n c.. ij  which i s  mk  representation. {c_  : i = 0  and  J  s.(x) =  c.Cx-S.)  3  .... k .  (m + 1 ) ( k + 1) = mk + m + k + 1  more than a r e r e q u i r e d  unknowns  f o r a non-redundant s p l i n e  T h e r e f o r e i t i s o n l y n e c e s s a r y to compute the s e t j = 0, 1,  ..., m  and  i = 1, 2,  and  j = m}  The remaining c o e f f i c i e n t s can be computed from the c o n s t r a i n t s  derived  from the c o n t i n u i t y c o n d i t i o n s ; t h a t i s  1 d —, — . dx J  c. . 1 J  =  j !  , . s . . (x)  6.  3  l  for  i = 1, 2,  k;  j = :0, 1,  m-1  ..., k  .  - 6 T h i s system i s u s e f u l f o r a p p r o x i m a t i o n as i t g i v e s an r e p r e s e n t a t i o n f o r each p i e c e w i s e  polynomial  explicit  between each knot p a i r .  But  d e s p i t e the f a c t t h a t the b a s i s f u n c t i o n s are d e f i n e d over a p a r t i c u l a r knot i n t e r v a l , they must be computed over the e n t i r e i n t e r v a l . the system o f e q u a t i o n s formed tend r e s u l t i n g s o l u t i o n i s not  Also,  to be i l l - c o n d i t i o n e d ; t h a t i s , the  accurate.  Since t h i s r e p r e s e n t a t i o n i s  i d e n t i c a l to the mathematical r e p r e s e n t a t i o n between each knot p a i r , i l l - c o n d i t i o n i n g occurs  2.4  f o r the same r e a s o n (as d e s c r i b e d  Mathematical S p l i n e  The splines  standard  power f u n c t i o n s ) .  mainly i n mathematical a n a l y s i s .  manipulate  elementary  This representation  i s used  Most theorems i n v o l v i n g s p l i n e f u n c t i o n s  proved u s i n g elementary s p l i n e s as they are easy to  analytically.  D e f i n i t i o n 3:  An  elementary s p l i n e f u n c t i o n of degree  defined  m,  y™  ,  is  by m y y  +  =  7  \  y™  for  0  for  y > 0  J  y < 0  Elementary s p l i n e s g i v e r i s e to a s e t of b a s i s f u n c t i o n s splines.  2.4).  Representation  r e p r e s e n t a t i o n of s p l i n e s i s t h a t of  ( a l i a s truncated  are d e r i v e d and  i n Section  In p a r t i c u l a r the  {1, x, x , 2  for  set  x , m  (x-6^)™ ,  (x-6 >™  forms a s e t of b a s i s f u n c t i o n s f o r a s p l i n e of degree of these b a s i s f u n c t i o n s f o r a cubic s p l i n e with i n t e r n a l knots i s g i v e n i n F i g u r e 1 .  }  k  m  .  An  four uniformly  example spaced  - 7 These b a s i s f u n c t i o n s can a l s o be problem.  Given the s e t of c o e f f i c i e n t s  the problem i s the d e t e r m i n a t i o n  m+1  S(x)  =  T  ,  L  1=1  T h i s system r e s u l t s  in  n  1  a  j_'  s  2,  approximation  . .., m + k +  1}  where  I • a.j.-Cx-fi'.)"  .  x+m+1  L  1=1  m + k + 1  : i = 1,  i n t o an  k  +  L  l  f o r a unique r e p r e s e n t a t i o n .  {a^  of the  ._.  a.x  formulated  n  l +  .  unknowns - e x a c t l y the number needed  Therefore  a l l the c o e f f i c i e n t s must  be  computed.  Despite  i t s s i m p l i c i t y the mathematical r e p r e s e n t a t i o n  s p l i n e s should never be used f o r c o m p u t a t i o n a l purposes. by  t h i s method w i l l produce i l l - c o n d i t i o n e d  m + k + 1  increases.  of  S p l i n e s computed  systems of e q u a t i o n s  as  I n t u i t i v e l y , a r e a s o n f o r t h i s can be seen from  the example p l o t t e d i n F i g u r e 1.  Notice  t h a t the l a s t  basis  function  3 (x-6^)  +  i s zero n e a r l y everywhere as w e l l as b e i n g  r e l a t i v e to the o t h e r b a s i s f u n c t i o n s .  extremely  small  Consequently i t i s p o s s i b l e to  produce a l i n e a r combination o f these b a s i s f u n c t i o n s which i s almost z e r o ; t h a t i s , the s e t o f b a s i s f u n c t i o n s i s almost l i n e a r l y dependent. Using  t h i s s e t w i l l produce a system of l i n e a r e q u a t i o n s f o r the  i m a t i o n problem whose c o r r e s p o n d i n g m a t r i x w i l l be  ill-conditioned  matrix  i s nearly singular.  i n most c a s e s .  are more s t a b l e  computationally.  This  Hence i t i s b e t t e r to  choose b a s i s f u n c t i o n s which are more d i f f i c u l t to c o n c e i v e but  approx-  analytically  ELEMENTARY CUBIC SPLINE BASIS FUNCTIONS  FIGURE 1  - 9  2.5  B-spline  Representation  In d e r i v i n g a s e t o f b a s i s conditioned  functions  which produce w e l l ^ •  systems; one would l i k e t o s a t i s f y these two c r i t e r i a ;  1.  That t h e support  (region  i n which the f u n c t i o n a l value  i s nonszero) o f t h e s p l i n e s i s f i n i t e ; 2.  That t h e number o f knot i n t e r v a l s i n v o l v e d i s minimal  To  (as s m a l l  t h i s end a s e t o f b a s i s  divided  as p o s s i b l e ) .  functions  i s derived  on t h e concept o f  differences.  Definition^;  The f u n c t i o n a l  Df  =  f(6  ± F  «  D  d e f i n e d by  6  1 + 1 >  • l+l  i + m + 1  )  i+m+1 i+m+1  is  The  i n t h e support  the d i v i d e d  difference of  d i v i d e d d i f f e r e n c e depends l i n e a r l y on  important, t h e d i v i d e d d i f f e r e n c e o f o r d e r p o l y n o m i a l o f degree  l  l+m  i f_ o f  f(x) . m + 1  order  Also,  m + 1 .  and more  i s zero f o r any  m .  Now, f o r any f u n c t i o n  g,  by t h e Lagrange i n t e r p o l a t i o n formula,  -  where  co(x) = v^-^) * ( x - S ^ ) • • • ( x - 6 ^ ^ ^ )  .  i +  I n order t o define  g(o"., ) b e t h e ( m + 1 ) s t i+n.  B-splines l e t the function of  10 -  C^^^^x)™ . S u b s t i t u t i n g f o r  g , t h e Lagrange  divided difference interpolation formula  gives C6.  m+1 gC6 6 i ?  where  i + 1  ,...,6  w(x) i sas previously Letting  i + m + 1  )  I n=0  =  + n  -x)  m  T~  i+n  defined.  i range over  -m t o k  g i v e s e x a c t l y t h e number  of B - s p l i n e s r e q u i r e d t oform a s e t o fb a s i s f u n c t i o n s ; that i s m + k + 1  functions.  : i = -m,-m + 1 ,  Thus, t h es e t (s^(x)  k}  where S .  m+1  r \ (x)  (6  V  =  >  ,  F A  w ( x ) = (x-6 .) • (x-6 . .. ) • • • (x-6 l  B-spline basis  i + I  ) i s precisely t h es e t o f  XTTirrl  functions.  To c o m p l e t e 2m  t —  i+n  n=0  with  -x)?  i + n +  supplementary  the definition, this s e trequires t h eaddition o f  knots t ot h e o r i g i n a l knot s e t A =  '^1'"'"'^k+1^ "  These k n o t s must b e e x t e r n a l t o t h e o r i g i n a l knot s e tw i t h less than  6^ a n d m  knots greater than  $^+1 * ^  of choosing these e x t r a knots i st h ef o l l o w i n g :  Vol: *-i 6  1+k+1  for  i = 1 , 2 , ..., m  o  kTT~  6  =  =  6  k+1  1  +  %+I^DI* 1  k+  1  n e  P  o  s  s  l  D  m  knots  le  method  - 11 -  For  B - s p l i n e s i t can-Be shown t h a t :  1.  s . Cx) i  2.  i s strictly positive in  The s u p p o r t of interval  3>  s  [6\. ,8.., I  J  ^(- )  is finite  x  ...]:•  x+m+1 '  and r e s t r i c t e d  to the  1 ^ , 6 . ^ ] ;  Any s p l i n e  S(x)  can be u n i q u e l y r e p r e s e n t e d as a l i n e a r  combination o f  ( s ^ ( x ) : 1 = *-m,  -m + 1, ..., k} .  I t i s s i m p l e to f o r m u l a t e an a p p r o x i m a t i o n problem from B-spline functions. the  Given the s e t o f c o e f f i c i e n t s  problem i s to f i n d the  a  ^'  s  S(x)  where for  {s_^(x) : i = -m,  B-splines.  =  k}  n  k T x=-m  -m + 1,  k}  T h i s system r e s u l t s i n  mathematical system. the  i  {a^ : i = -m  a.s.(x)  i s the s e t o f b a s i s m + k + 1  functions  unknowns as i n the  There a r e no redundant parameters and hence a l l  c o e f f i c i e n t s must be computed.  The system of e q u a t i o n s d e r i v e d u s i n g B - s p l i n e s remains wel'lc o n d i t i o n e d as  m + k + 1  increases.  In f a c t , one can produce n u m e r i c a l  upper bounds on the c o n d i t i o n number of the m a t r i x of normal e q u a t i o n s for  a u n i f o r m l y spaced knot s e t f o l l o w i n g the method d e s c r i b e d i n  Schultz  [18, pp. 70-72] .  A l s o , because the b a s i s f u n c t i o n s g i v e minimal  support the systems produced a r e banded w i t h the band-width dependent the  degree o f the s p l i n e .  on  - 12 -  2.6  Example o f B - s p l l n e  Representation  To g i v e a c o n c r e t e look l i k e , partition. (s^(x)  consider  representation  : i = -3, -2,  k}  Given a u n i f o r m  for  i = -3, -2,  f o r the b a s i s  functions  p a r t i t i o n , the mesh l e n g t h i s  and t h e r e f o r e the  k + 3, k + 4 .  ±  i = -3, -2,...., k  uniform  can be developed i n the f o l l o w i n g manner:  s (x)  for  functions  the p a r t i c u l a r case o f c u b i c s p l i n e s on a  An e x p l i c i t  h = T~r k+1  example o f what B - s p l i n e b a s i s  n=0  i - t h knot i s  6. = i  Thus  »'<W  with  4 co'(x)  =  n  (x  j=o  S u b s t i t u t i n g f o r the knots  gives  n=0  j n=0  i+n - x) :fk+l i+n CJ' ^k+lj  (i+n -  3~ +  (k+l)x)^  (k+1)'  i+n k+1  ^ k+1  - 13 -  Now 4  (i+n-i-j) k+1  n  j=0  4  n j=o k+1  1  •  O+ij*  4  ^  C n _ j ) 0  Jr^V S u b s t i t u t i n g back i n t o  s  ^Cx)  4 s.(x)  (i+ri-(k+l)x)'  I  =  n=0  (k + l )  •(t'+l) 4 n (n-j) 3=0  3  4  =  (k+1)  (i+n-(k+1)x)' 4 n=0 (n-j) n  I  j=o  Expanding f o r n  s  .  ( x )  and .  j :  ( k + 1 )  {  Ci-(k+l)x) , +  Cl+2-(k+l)x)  ,  (i+4-(k+l)x)  +  for  i = -3, -2, ..., k .  _ (i+l-(k+l)x)  3  i  ^ -  r  —  3  + _ —( i—+ 3-- ( k + l ) x+) 3  3  +  ,. i  3  - 14 -  Letting  x' = (k+l)x - i -• 2  and d e f i n i n g  S(x*) = s ^ k + D x  r  S(x')  gives  5  - i - 2)  x' <_ -2  o q+x') 6  3  (x') 4  3  (x*) q-x') 4 ~ 6 3  d-x') ~ 6  3  (2-x') 24  3  , (2-x') 24  3  -2 <_ x' £ -1  3  -1 £ x' £ 0  (k+l)< (1-x*) 6  3  , (2-x') 24  3  0 <_ x* £ 1  (2-x') . 24  1 £ x' £ 2  0'  2 < x'  3  When t h e e x p l i c i t r e p r e s e n t a t i o n i s used on a. u n i f o r m of  [0,1] w i t h f o u r i n t e r n a l k n o t s ; t h e s e t  {s^(x)  .  partition,  : i = -3, -2, ..., 4}  i s as shown i n F i g u r e 2.  Equivalence of Spline Representations  A l t h o u g h B - s p l i r i e s p r o v i d e a s u i t a b l e method f o r s o l v i n g approximation  problems t h e c o e f f i c i e n t s . o b t a i n e d a r e n o t extremely  spline useful.  In p a r t i c u l a r , . i t i s p r e f e r a b l e t o know t h e c o e f f i c i e n t s o f t h e p i e c e w i s e continuous p o l y n o m i a l s between a d j a c e n t knots than t o know t h e  CUBIC B-SPLINE BASIS FUNCTIONS  FIGURE 2  - 16 -  coefficients  of B - s p l i n e b a s i s f u n c t i o n s .  To t h i s end, i t i s a p p r o p r i a t e  to d e r i v e t h e unknown s e t o f t h e p i e c e w i s e continuous p o l y n o m i a l {c..  : i - 0,1,  coefficients  ..., k, j = 0, 1, ..., m}. from the known v a l u e s o f the  B-spline coefficients  { a ^ • 1 " 1> 2,  m + k +1}  .  C o n s i d e r t h e f u n c t i o n a l v a l u e o f the s p l i n e f o r each knot; that i s , at ' Case I ' :  x = 6. . x  F o r p i e c e w i s e continuous  for  Case II:  d dx  "  J  1  c. .  polynomials  1X  dx  J <  i = 1, 2,..., k + 1 ;  j =0,  J  1,  m  For B-splines  d  - S(x) 6\ dx x  j  v  4  dx  J  £=-m  J  d  k  I £=-m for  i - l ,  k + 1;  2,  * ^  j  a„ — dx  3  s„(x) ^  6. x  j = 0, 1, ..., m  Now  o£  —, s„(x) dx  J  6.  dx j  %+n - + X) m  *f  riio  U  ' W C  x  6". X  ,J _ •- d " .m io-^W. d^V' + m+1. t  i  1  1  C  X)  Differentiating with respect to x  ~i  dx  -e+n - +  C6  x)  to)  X  !  C ^ - < -  J  0  By equating the terms of the derivatives:  i  for  i = 1, 2,  k + 1;  *  j = 0/1,  £  m  - 18 -  Section 3  THE LINEAR APPROXIMATION PROBLEM 3.1  Introduction  S e v e r a l methods a r e known f o r s o l v i n g l i n e a r problems.  These methods can be a p p l i e d t o problems  sets of basis functions.  approximation  involving general  S p l i n e a p p r o x i m a t i o n w i t h a f i x e d knot s e t i s  a p a r t i c u l a r a p p l i c a t i o n o f the g e n e r a l problem/  de Boor and R i c e [5] d e s c r i b e an a p p r o x i m a t i o n method  involving  o r t h o g o n a l p r o j e c t i o n . ' The b a s i c i d e a i s t o minimize t h e e r r o r  ||y-u||  of a p p r o x i m a t i n g  y .  Py  y  by  u  by t h e o r t h o g o n a l p r o j e c t i o n  i s b e s t c a l c u l a t e d u s i n g an orthonormal b a s i s .  Py  of  Therefore, given a  g e n e r a l s e t o f b a s i s f u n c t i o n s f o r t h e a p p r o x i m a t i o n , an orthonormal s e t of b a s i s f u n c t i o n s must be d e r i v e d ,  de Boor and R i c e use a m o d i f i e d  Gram-Schmidt p r o c e s s to generate such a s e t o f b a s i s  functions.  The most common t e c h n i q u e used f o r . s o l v i n g l i n e a r problems  approximation  i s the method o f normal e q u a t i o n s ( d e s c r i b e d i n S e c t i o n 3.3).  [13] d i s c u s s e s t h e g e n e r a l l i n e a r l e a s t squares problem and l i n e a r  least  squares problem u s i n g s p l i n e s i n d e t a i l g i v i n g r e s u l t s c o n c e r n i n g t h e uniqueness o f t h e s o l u t i o n and the symmetric  and p o s i t i v e d e f i n i t e  e r t i e s o f t h e a s s o c i a t e d l e a s t squares m a t r i x .  prop-  Patent a l s o i n c l u d e s a  program s o l v i n g the s p l i n e a p p r o x i m t i o n problem w i t h f i x e d k n o t s .  The  b a s i s f u n c t i o n s used i n g e n e r a t i n g t h e system,of normal e q u a t i o n s were B-splines.  Patent  - 19 -  Smith [19, pp. 110-119] develops the method of noraal equations for the p a r t i c u l a r case of s p l i n e approximation.  The p a r t i c u l a r set of  basis functions used are the mathematical representation.  The.least  squares matrix derived using these basis functions i s given i n f u l l . Golub [7] develops a method using Householder transformations Cdescribed i n Section 3.3). An A l g o l program based on t h i s procedure f o r a general set of b a s i s functions i s given i n Businger and Golub [ 3 ] . Although a s u i t a b l e s e t of b a s i s functions i s a v a i l a b l e which prevents i l l - c o n d i t i o n i n g as tbe number of knots increases, there i s s t i l l the problem of preventing i l l - c o n d i t i o n i n g as the knots become nonuniformly spaced.  The orthogonal p r o j e c t i o n method counteracts t h i s  problem because the basis functions are orthonormalized before  solution.  The method of normal equations frequently produces i l l - c o n d i t i o n e d systems as shown i n the example c i t e d i n Golub [ 7 ] .  Solving the s p l i n e approx-  imation problem using normal equations on a non-uniformly a prime example of t h i s i l l - c o n d i t i o n i n g .  knot set i s  However, the method of  Householder transformations counteracts t h i s problem because of the orthogonality of the transformations.  Consequently, t h i s method i s  necessary f o r a s t a b l e s o l u t i o n to the s p l i n e approximation  3.2  problem.  D e f i n i t i o n of the Linear Problem A l l aspects of the l i n e a r problem are incorporated i n the  following d e f i n i t i o n :  - 20 -  Definition  5:  Given a d i s c r e t e s e t of d a t a  { (x£,y_£) : £ ~ 1, 2, . .. , n}  and a f u n c t i o n  S(x)  where  =  m T i=l  a.s.(x)  {s^(x) : i = 1, 2, ..., m} • i s a s e t o f b a s i s  f u n c t i o n s and  {a^ : i = 1, 2, ..., m}  unknown c o e f f i c i e n t s  i s a s e t of  occurring l i n e a r l y ;  the l i n e a r  a p p r o x i m a t i o n problem i s to determine v a l u e s f o r the a^'s  to produce the " b e s t f i t " o f  S(x)  to the s e t  of d a t a .  3.3  S o l u t i o n o f the L i n e a r  Problem  Most approximation problems least icular  squares e r r o r as s a t i s f y i n g form o f the l e a s t  o f the l e a s t  Definition  6:  c o n s i d e r the m i n i m i z a t i o n of the  the " b e s t f i t " c r i t e r i a .  squares e r r o r used i n t h i s case i s the square  squares norm where:  Given a v e c t o r  the  least  The p a r t -  v  =>  squares norm o f  v,  | v  , is  Thus, g i v e n the f u n c t i o n • S(x)  and  t h e data s e t  { C x ^ y ^ ) .: £.=  1, 2,  the s o l u t i o n of the l i n e a r problem becomes the. m i n i m i z a t i o n of t h e squares e r r o r by the a p p r o p r i a t e c h o i c e o f the unknowns  least  ^ a ^ : i = 1,  That i s , b y . f i n d i n g  n  I [y  min {a }  m £ as i - l  -  p  l«X- *-  ±  1  2 (x„)]  : .  1  The u s u a l method o f s o l u t i o n , t h a t o f normal e q u a t i o n s , i s developed  i n the f o l l o w i n g  way:.  In o r d e r to f i n d the minimum  {a } ±  2  n  m  £=1  i=l  1  *-  1  d i f f e r e n t i a t e the summation w i t h r e s p e c t to each o f the parameters {a. l  : i = 1, 2,  m}  „ 9 a  VM: j  and s e t to z e r o .  n  i=i  = -2  Thus  m £-  [y  X i=i  *-  n I,  1=1  [y» L  9 a  i  s  i ^ > ]  m - I . i - l  }  a s (x„).] • s (x») * *-  1  1  3  = 0  for  j = 1, 2,  m .  L=l ^ J l  Rearranging  the terms g i v e s  'J £ S  ( X  ) ]  =  J  ^ J £ S  x  ( X  }  ..  2,  - 22 -  Using  the f o l l o w i n g changes:  1.  Denote by. S  the  S  for 2.  1=1,  Denote by  m x m  n  =  S.. :=  2, y_  matrix  m the  and  j = 1, 2,..., m ;  m-dimension v e c t o r  n  for 3.  j = l , 2, ..., m ;  Denote by  &_  the  and  m-dimension v e c t o r of unknown c o e f -  i ficients  a_ = {a_^ : i = 1, 2, ..., m} ;  the problem becomes one o f f i n d i n g the s o l u t i o n t o the system of  Sa =  equations  y_ . Another method o f s o l v i n g  n I  min  m [y  {a.} 1=1 l i s by u s i n g o r t h o g o n a l  ~  p  2  I  *-  a.s.(xo)]  i=l  1  1  transformations.  X  '  T h i s r e q u i r e s the f o l l o w i n g  changes:  1.  Denote by  S  the  n x m  matrix of f u n c t i o n values;  is, S  for  £ = 1, 2,  =  s  l±  =  n  s  i  and  ( x  £  }  i = l,  2  m;  that  2.  Denote by  y_ 3.  y_ the  =  {y£ •  Denote by  n-dimension vector of ordinates  &  a. the  1?  =  n) ;  and  m-dimension vector of unknown coef-  ficients a,  (a^ : i - = 1, 2,  =  m} ;  then the problem becomes to find  min a.  2  | |_y_ - SaJ |  Consider multiplying the previous equation by an orthogonal T matrixs Q  .  Because multiplying by an orthogonal matrix does not  change the norm; the l i n e a r least squares problem remains the same.  Thus  the problem becomes to find min a_  |  Q Sa| | T  .  2  T Now  consider  Q  to be a series of orthogonal transformations  T which transforms  Q S  into an upper triangular matrix  R .  I f such a  series can be found then the l i n e a r least squares problem reduces to finding min a. Since the zero part of to. solve the system  R  | | Q,y_ - Ra | | . T  2  i s independent of  Ra = b_ where  a,;  i t i s only necessary •  T b. = (Q j ) f o r  i = 1, 2,  m .  - 24 -  The remainder o f  Q y_ c o n t a i n s the l e a s t squares e r r o r ;  that i s ,  n  IZ " S a l  I (Q y)' i=m+l T  There e x i s t s a s e r i e s o f o r t h o g o n a l t r a n s f o r m a t i o n s w i l l reduce  S  formations.  They can be c o n s t r u c t e d as f o l l o w s :  Q  t o an upper t r i a n g u l a r m a t r i x known as Householder  Given a v e c t o r matrix  P  v  which trans-  c o n s t r u c t a symmetric o r t h o g o n a l  such t h a t  vector.whose f i r s t  Pv = w  where  element i s  w  ±||v||  i s a unit and whose  remaining elements a r e z e r o .  Householder showed.that  t h e r e e x i s t s a symmetric w = Pv . [1,  f o r any two v e c t o r s  orthogonal matrix  The symmetry and o r t h o g o n a l i t y o f  v,  w  with T  T T v v = w w  P = I - 2uu . such t h a t P  i s proven i n A c t o n  P. 327] . The problem now i s t o determine t h e r e q u i r e d v e c t o r  P .  The method i s d e s c r i b e d i n Acton [ 1 , pp. 324-329]  m o d i f i c a t i o n s and t h e d e r i v e d  _u i n  with s l i g h t  u is  v  where  K  "2 2 K = 2 . v | ± 2v^| v|  n  Two c o m p u t a t i o n a l c o n s i d e r a t i o n s come i n t o e f f e c t when u s i n g Householder t r a n s f o r m a t i o n s .  The f i r s t  i s which s i g n t o choose i n  - 25 -  computing  _u .  K  The c h o i c e i s to p i c k  max{2| |v| | + ^ 1 v .  =  2  2  =  choose  such t h a t  |v| | , /  2| |v| | s  K,  max!  i n order to avoid c a n c e l l a t i o n . v^ < 0  K  2  - 2v  ±  | | v| | } •  K„  Thus i f  v^ >^ 0  choose  K^;  if  IC, .  The second c o n s i d e r a t i o n i s the computation o f  Pv .  Rewriting Pv. =  (I -  2uu )v  Pv  =  T I v - 2 mi v  =  T v - _u(2_u v)  T  as  T the s c a l a r  2u_ v  Hence the m a t r i x f a r more e f f i c i e n t  i s computed f i r s t P  need never be formed e x p l i c i t l y . than forming  P  To manipulate Householder t r i a n g u l a r matrix consider applying e q u i v a l e n t to a p p l y i n g . P P^  such t h a t  P^S  f o l l o w e d by the v e c t o r s u b t r a c t i o n .  and p e r f o r m i n g a m a t r i x m u l t i p l i c a t i o n . transformations P  reduces column 1 o f  to form the upper  to the m a t r i x  t o each column i n  w„  T h i s method i s  S,  S . w^,  S .  This i s  That i s , t h e r e i s a to  - 26 -  where  w|  i s the transformed column 1 o f  columns o f  S  S,  w_2,  Note that, a l l the o t h e r  a r e a l t e r e d by t h i s t r a n s f o r m a t i o n .  S i m i l a r l y , there i s a of  S .  . such that.  P^S  reduces  column 2  to  •±11*2 I  where  w^  of  are a l t e r e d a l s o .  S  i s , t h e transformed column 2 o f  S .  However a l l o t h e r columns  In p a r t i c u l a r , column 1 r e v e r t s back to non-  zero s t a t u s which i s not d e s i r a b l e .  T h e r e f o r e , i n o r d e r to p r e s e r v e the zeroes i n column 1 l e t P_  be the Householder  t r a n s f o r m a t i o n which reduces  column 2 o f  S,  w ,  ±11*2 I  where  w^  element o f  i s the transformed column 2 o f w' , "2  This.transformation w i l l  w i l l a l t e r a l l the remaining columns of m  times;.  S  S  can be reduced  and  w^  i s the f i r s t  l e a v e column 1 unchanged but  S .  C o n t i n u i n g the p r o c e s s  to an upper t r i a n g u l a r m a t r i x o f the form  - 27 -  ±11^1I1I • !w +n  2  w  w  m  w„  w P  m  x  :  %  P.P-S 2 1  0  =  m-1 m  ±||w |  0  Thus i n the l i n e a r l e a s t squares problem u s i n g applying  3.4  Q  T  to  S,  Q  T  reduces  S  Q  =  p m  **'  P P  to the upper t r i a n g u l a r  2  and  1  matrix  R .  S o l u t i o n o f the F i x e d Knot Problem  The f i x e d knot problem f o r s p l i n e s of the g e n e r a l  is a particular  l i n e a r l e a s t squares problem.  the l i n e a r problem a l l o w the f u n c t i o n s f u n c t i o n s w i t h a f i x e d knot s e t .  In the d e f i n i t i o n 5 f o r  s^(x)  This creates  squares a p p r o x i m a t i o n problem f o r s p l i n e s .  application  t o be the B - s p l i n e the f i x e d knot  least  T h i s problem can be  e x a c t l y l i k e the g e n e r a l problem u s i n g Householder  basis  solved  transformations.  However because a p a r t i c u l a r b a s i s s e t i s b e i n g used a few c o m p u t a t i o n a l considerations  come i n t o  effect.  The  first  of these i s t h a t o f b a s i s  Because the b a s i s f u n c t i o n s a r e e v a l u a t e d  function  o f t e n i t i s important t h a t  e f f i c i e n t method of computation be a v a i l a b l e . the v a l u e s  o f the b a s i s f u n c t i o n denominator vw  throughout a l l c o m p u t a t i o n s oo'(x)  initially  and  Therefore  s  r e t a i n the v a l u e s  evaluation. an  In the f i x e d knot case (x)  1  remain"'.constant  i t i s advantageous to c a l c u l a t e for future  use.  A l s o , because the b a s i s f u n c t i o n s have m i n i m a l support i t i s o n l y n e c e s s a r y to compute a b a s i s f u n c t i o n a l v a l u e range of s u p p o r t . to c a l c u l a t e the of s u p p o r t .  Otherwise the f u n c t i o n a l v a l u e  T h i s has  two  disadvantages.  i t i s more time-consuming Second because of  i s the r e d u c t i o n o f  S  S . has  Because the b a s i s f u n c t i o n s a banded s t r u c t u r e .  I f the d a t a p o i n t  [6., 6..  ... ]  x^  T h i s banded  then the  Hence the m a t r i x w i l l have b l o c k s  s».  matrices  the l o c a t i o n v o f the  i s outside  abscissa  the support r e g i o n ,  element of  S  w i l l be  o f non-zero elements a l o n g  where  6-. < x„ < 6., ,.. . i — Z — x+m+1  the  T h i s s p e c i a l banded s t r u c t u r e can r  that  0' .  where each non^zero b l o c k of the m a t r i x i s formed from the' v a l u e s Xn Z  the  have  the f a m i l a r handedness u s u a l l y a s s o c i a t e d w i t h  r a t h e r a s p e c i a l s t r u c t u r e dependent on  x» £  region  second c o m p u t a t i o n a l c o n s i d e r a t i o n  s t r u c t u r e i s not  .  It i s possible  o f machine a c c u r a c y  l i m i t e d support the m a t r i x  x^  .  the o r d e r  r e s u l t i n g l e a s t squares m a t r i x  point  0  the  zero.  The  is,  First  than to t e s t f o r the range.  r o u n d - o f f e r r o r , the summation w i l l be  but  is  f u n c t i o n at a l l p o i n t s w i t h o u t t e s t i n g f o r the  to c a l c u l a t e a v a l u e  r a t h e r than  i f i t i s within  diagonal of be  - 29 -  used i n t h e f i x e d knot a p p r o x i m a t i o n . by m a n i p u l a t i o n by Householder reduce t h e banded p a r t o f of  computational  S .  S i n c e the zeroes a r e unchanged  transformations i t i s only necessary to T h i s causes a c o n s i d e r a b l e r e d u c t i o n  effort.  The t h i r d c o m p u t a t i o n a l c o n s i d e r a t i o n i s t h e c o n d i t i o n number of  the matrix.  ficult  Since  S  is initially  t o s t a t e a n y t h i n g about  a non-square m a t r i x i t i s d i f -  i t s c o n d i t i o n number.  However t h e o r t h o -  gonal t r a n s f o r m a t i o n s do not a f f e c t t h e norm (and hence t h e c o n d i t i o n number) o f  S  i n any way.  Thus, by d e t e r m i n i n g t h e c o n d i t i o n number  of  t h e square upper t r i a n g u l a r p a r t o f t h e m a t r i x  S  has e f f e c t i v e l y been determined.  the c o n d i t i o n of  The s t i p u l a t i o n f o r a w e l l - c o n d i t i o n e d  m a t r i x i s t h a t t h e c o n d i t i o n number be l e s s than the c o n d i t i o n number f o r a g i v e n problem the l o c a t i o n o f t h e k n o t s .  R  in + k + 1  s h o u l d remain  and t h a t  independent  of  - 30 -  Section 4  THE VARIABLE KNOT APPROXIMATION PROBLEM  4.1  Introduction  The f i x e d knot problem has been i n v e s t i g a t e d t h o r o u g h l y and adequate methods e x i s t f o r i t s s o l u t i o n .  The v a r i a b l e knot  problem  i n v o l v e s c h o o s i n g l o c a t i o n s f o r t h e knots so as t o p r o v i d e t h e " b e s t approximation" p o s s i b l e . problem  Because t h e knots occur n o n - 1 i n e a r l y , . t h i s  i s more complex.  One a l t e r n a t i v e t o s o l v i n g the . v a r i a b l e knot problem minimize  t h e l e a s t squares e r r o r w i t h r e s p e c t  to the knots.  i s to  Another  a l t e r n a t i v e , i s t o p o s i t i o n t h e knots v i s u a l l y t o a r r i v e a t an a p p r o x i m a t i o n which i s r e a s o n a b l e to t h e eye a l t h o u g h n o t n e c e s s a r i l y t h e " b e s t " i n the l e a s t squares sense. approximation problem.  This requires  g r a p h i c a l i n t e r a c t i o n w i t h the s p l i n e  By i n t e r a c t i o n a r e a s o n a b l e a p p r o x i m a t i o n can be  d e r i v e d manually, and good i n i t i a l guesses  f o r an automatic m i n i m i z a t i o n  technique can be o b t a i n e d .  de Boor and R i c e [6] s o l v e t h e v a r i a b l e knot problem by an a u t o matic t e c h n i q u e .  They minimize  t h e l e a s t squares e r r o r i n i n t e g r a l form  fb ., \ . { [y - S ( x , A ) ] V ^a 2  over a l l s p l i n e s o f degree  ,m  with  k  knots.  The t r a p e z o i d a l r u l e i s  used t o o b t a i n an approximation t o t h i s , i n t e g r a l .  A d i s c r e t e Newton's  method i s a p p l i e d t o minimize each knot i n d i v i d u a l l y w h i l e t h e r e s t o f  - 31 -  the  knots remain s t a t i o n a r y .  The knots a r e o p t i m i z e d by sweeping  the  knot s e t from r i g h t to l e f t and r e - e v a l u a t i n g the f i x e d knot  through problem  each time.  Smith [19, p. 110^119] p r e s e n t s t h e u s e o f s p l i n e s i n i n t e r a c t i v e data f i t t i n g .  A l t h o u g h he does n o t a l l o w f o r t h e p o s s i b i l i t y o f  r e s p e c i f y i n g c e r t a i n knot l o c a t i o n s ; he does a l l o w the p o s s i b i l i t y o f r e s p e c i f y i n g t h e e n t i r e knot s e t and a t t e m p t i n g the f i x e d knot  problem  again.  4.2  D e f i n i t i o n of the V a r i a b l e Knot  Problem  Whereas t h e l i n e a r l e a s t squares problem can be g e n e r a l i z e d to any s e t o f b a s i s f u n c t i o n s ; the d e f i n i t i o n o f the n o n - l i n e a r a p p r o x i m a t i o n problem i s r e s t r i c t e d a b l e knot  to s p l i n e f u n c t i o n s and i s r e f e r r e d t o as the v a r i -  problem.  D e f i n i t i o n 7:  Given a d i s c r e t e s e t o f d a t a  {(x^,y^)  : I •= 1, 2  n} ;  a s e t o f knots i n s t r i c t l y  A  =  : 1=1,  increasing order  2, ..., k> ;  a s e t o f s p l i n e b a s t s functions:  {s.(x) : i =» -m,  -m .+ 1,  k} >  - 32 -  and a s e t o f l i n e a r  coefficients  {a.^ : i =" -m,  l , ••... j k} ' k  with  S(x,A)  =  y  a.s.Cx.A.)  •1' « i x x--m  the v a r i a b l e knot problem i s t o determine v a l u e s f o r the s e t  A  to produce t h e " b e s t f i t " o f  S  to the  set of data.  Note t h a t i n t h i s case the unknown parameters a r e t h e knots - not t h e linear coefficients  {a^  :  i = -m,  -m + 1 ,  ..., k} .  To be c o m p l e t e l y  c o r r e c t b o t h t h e c o e f f i c i e n t s and the knots should be e v a l u a t e d s i m u l t a n e o u s l y t o produce t h e " b e s t f i t " .  T h i s problem i s much more d i f f i c u l t .  However, a f e a s i b l e a l t e r n a t i v e i s f o r each knot s e t coefficients  {a^ : i = -m,_-m + 1,  k}  A,  t o p i c k the  by s o l v i n g t h e f i x e d knot  problem.  4.3  S o l u t i o n o f the V a r i a b l e Knot  Problem  In the i n t r o d u c t i o n to t h i s s e c t i o n two a l t e r n a t i v e s were proposed f o r s o l v i n g t h e v a r i a b l e knot problem: an automatic . procedure and an i n t e r a c t i v e p r o c e d u r e . order s i n c e the f i r s t  I t i s b e s t to d i s c u s s t h e s e i n the r e v e r s e  f o l l o w s i n h e r e n t l y from t h e second.  The use o f i n t e r a c t i o n f o r s o l v i n g s p l i n e a p p r o x i m a t i o n problems can be summarized i n t h e f o l l o w i n g s t e p s :  - 33 -  1.  Read i n the d a t a s e t and graph i t on a g r a p h i c s t e r m i n a l .  2.  S p e c i f y an i n i t i a l  s e t of knots and o v e r l a y t h e i r  location  on the x - a x i s . 3.  Solve the f i x e d knot problem u s i n g t h i s knot s e t and lay  over-  the r e s u l t i n g s p l i n e c u r v e .  4.  A l l o w r e s p e c i f i c a t i o n of the l o c a t i o n o f any of the k n o t s .  5.  R e c a l c u l a t e the f i x e d knot problem u s i n g the new  knot s e t  and graph the r e s u l t i n g s p l i n e a p p r o x i m a t i o n curve i f i t i s wanted. 6.  A l l o w the o p t i o n s of r e s p e c i f y i n g k n o t s , changing the number of k n o t s , or o p t i m i z i n g the knot s e t .  T h i s i n t e r a c t i v e procedure produces a r e a s o n a b l e f i t much f a s t e r than an automatic procedure.  For example, i f the o r i g i n a l knot s e t  g i v e s a poor a p p r o x i m a t i o n , the s i t u a t i o n can immediately be  remedied  by m a n i p u l a t i n g the k n o t . s e t d r a s t i c a l l y as opposed to the more c a u t i o u s procedures of automatic t e c h n i q u e s . fast i n i t i a l  T h i s approach a l l o w s an extremely  approach t o a good f i t .  Then automatic r e f i n e m e n t c o u l d  use the r e s u l t i n g knot s e t t o r e a c h an o p t i m a l knot s e t q u i c k l y .  There a r e s e v e r a l p o s s i b l e approaches  f o r p o s i t i o n i n g the k n o t s ,  de Boor and R i c e [6, p. 12-18] p r e s e n t some o f these f o r an automatic procedure but v a r i a t i o n s seem s u i t a b l e f o r i n t e r a c t i v e  placement.  One p o s s i b l e approach suggested i s t h a t a d d i t i o n a l knots be p l a c e d near the l o c a t i o n of the maximum e r r o r . initially area.  as i t i s u s u a l l y d e s i r e d  T h i s seems r e a s o n a b l e  to produce a b e t t e r f i t i n t h a t  E v e n t u a l l y , ihowever, the d a t a i n t h a t r e g i o n w i l l become i n t e r -  - 34 -  p o l a t e d which i s not the d e s i r e d phenomenon.  A l s o t h i s procedure  not compensate f o r the a p p r o p r i a t e placement o f  does  fewer knots over t h e  e n t i r e range r a t h e r than the c o n c e n t r a t i o n o f knots  i n one p l a c e .  Another p o s s i b i l i t y i s to p l a c e knots a t p o s i t i o n s o f r a p i d change i n the d a t a .  T h i s a l l o w s the p o l y n o m i a l t o determine  i t s own  shape i n the i n t e r v a l n e a r l y i n d e p e n d e n t l y o f the s u r r o u n d i n g  interval.  T h i s i s because a t p o s i t i o n s o f r a p i d change, the h i g h e s t o r d e r term o f the p i e c e w i s e p o l y n o m i a l dominates.  I t i s p r e c i s e l y t h i s term t h a t i s  not i n c l u d e d i n the. c o n t i n u i t y c o n s t r a i n t s .  T h i s h e l p s overcome one o f  the b a s i c problems w i t h .polynomials - t h a t t h e i r o s c i l l a t o r y n a t u r e makes it  difficult  f o r .them to adequately  With a graph  approximate d a t a .  o f t h e d a t a w i t h i n reach i t i s p o s s i b l e t o make  r e a s o n a b l e p r e d i c t i o n s about.knot  placement.  This i s the greatest value  of g r a p h i c a l i n t e r a c t i o n - the d a t a and t h e a b i l i t y t o manipulate t h e approximation  4.4  a r e d i r e c t l y a t hand.  Knot O p t i m i z a t i o n  Despite the fact  t h a t r e a s o n a b l e approximations  can be made t o  a s e t o f d a t a i n t e r a c t i v e l y , o f t e n a more f o r m a l f i t c r i t e r i a i s d e s i r e d . T h i s can be a c h i e v e d by a u t o m a t i c a l l y r e f i n i n g the e x i s t i n g knot by m i n i m i z i n g the l e a s t squares the d a t a s e t  {( £> £) x  v  error.  : Z - 1,- 2,  Given t h e f u n c t i o n  .. ., n}  That  S(x,A)  e r r o r over t h e  i s , find  min min A {a } f  n £ ' [y„ -  Z=l  and  the s o l u t i o n o f the v a r i a b l e  knot problem becomes the m i n i m i z a t i o n o f t h e l e a s t squares knot s e t A .  set locally  k £ a s (x„,A)] i=-m 1  1  L  - 35 -  I t i s not  too c r u c i a l to o b t a i n the g l o b a l minimum i n the  i n t e r a c t i v e case as a r e a s o n a b l e e s t i m a t e o f the knot s e t a l r e a d y The  purpose of o p t i m i z a t i o n  The tially  minimization  i s merely to r e f i n e t h i s  method chosen i s known as COMPLEX which e s s e n -  discussed  here but  are adequately d e s c r i b e d  f o r c h o o s i n g COMPLEX i n v o l v e the and  estimate.  i n v o l v e s ! r e f l e c t i n g the f u n c t i o n around i t s c e n t r o i d .  are not  that i t w i l l  i n Box  f a c t t h a t i t does not  The d e t a i l s  [2],  Reasons  require derivatives  converge f a i r l y r a p i d l y towards the minimum.  COMPLEX c o n t a i n s s t r a i n t s can be  exists.  one  additional feature.  imposed on the  function.  This i s that  These c o n s t r a i n t s can  conbe  e i t h e r e x p l i c i t - meaning t h a t the independent v a r i a b l e can be bounded by some f u n c t i o n or c o n s t a n t ; value  or i m p l i c i t - meaning t h a t the  can be bounded by some f u n c t i o n or  constant.  In a p p r o x i m a t i o n u s i n g s p l i n e f u n c t i o n s have e x p l i c i t  i t i s o n l y n e c e s s a r y to  c o n s t r a i n t s to p r e v e n t the knots from c o a l e s c i n g .  c o n s t r a i n t s i n v o l v e k e e p i n g the knots separated How  functional  this distance  by  These  a certain distance.  i s determined i s p a r t i a l l y dependent on the machine  p r e c i s i o n and hence the m a t r i x used f o r s o l v i n g the f i x e d knot problem. Because o f machine p r e c i s i o n the knots must be as l e a s t as great m a t r i x w i l l be  as the machine a c c u r a c y .  separated  by  a  distance  Otherwise the l e a s t squares  singular.  More important, however, i s the c o n d i t i o n number of the m a t r i x used i n the f i x e d knot problem. the two  As  two  knots converge towards each  c o r r e s p o n d i n g rows o f the l e a s t . s q u a r e s  dependent c a u s i n g  the c o n d i t i o n number to r i s e .  other,  m a t r i x become more l i n e a r l y Therefore  adequate  con-  - 36 s t r a i n t s must be put on the knots to prevent them from c a u s i n g n u m e r i c a l instabilities.  F o r these r e a s o n s , the f o l l o w i n g c o n s t r a i n t was p l a c e d  on each knot:  set  h  =  6, ,- - 6_  k+1  0  and c o n s t r a i n each knot  6. . + .0001-h < 6. < i-l  for  i =  1, 2,  k .  —  l —  S. l  - .0001-h . l+l  by  - 37 -  NUMERICAL RESULTS  5.1  Introduction  There a r e t h r e e main areas where l e a s t squares can be used.  These a r e :  approximation  the a p p r o x i m a t i o n o f mathematical  functions;  the approximation of e x p e r i m e n t a l d a t a ; and computer-aided, d e s i g n . the f i r s t  o f these s p l i n e s do not p r o v i d e s u f f i c i e n t a c c u r a c y to warrant  t h e i r use as f u n c t i o n a l approximations. good r e s u l t s . of  In the. second arear.splines g i v e  In the t h i r d a r e a s p l i n e s a r e a b l e to f i t the c o n t o u r s  a d e s i g n extremely w e l l because o f t h e i r p i e c e w i s e n a t u r e .  In  o r d e r t o demonstrate  the p o s s i b i l i t i e s of i n t e r a c t i v e  s p l i n e approximation t h r e e examples a r e g i v e n .  The f i r s t e x a m p l e . f i t s  c u b i c s p l i n e s to an i n t e r e s t i n g s e t of d a t a m a i n l y t o demonstrate p o s s i b i l i t i e s of the method. to  The second example f i t s  the  a c u b i c s p l i n e curve  a d a t a s e t g i v e n i n de Boor and R i c e [ 5 ] , [6] i n v o l v i n g d a t a from a  T i t a n i u m heat experiment. o u t l i n e of a Volkswagen.  The f i n a l example i s the a p p r o x i m a t i o n of the T h i s r e s u l t i s merely  t h e . p o s s i b i l i t i e s o f s p l i n e s i n computer-aided any p r a c t i c a l  All to  In  an IBM  intended to  demonstrate  d e s i g n r a t h e r than h a v i n g  importance.  examples were run on an Adage G r a p h i c s T e r m i n a l connected  360/67 duplex o p e r a t i n g under MTS  (Michigan T e r m i n a l System)  l o c a t e d a t the U n i v e r s i t y o f B r i t i s h Columbia. o b t a i n e d on a Calcomp p l o t t e r . p r e s e n t e d i n Appendices  A  Hardcopy p l o t s were,  Program l i s t i n g s and a u s e r ' s guide are  and B r e s p e c t i v e l y .  - 38 -  The i n t e r a c t i v e system produces  1.  the f o l l o w i n g o u t p u t :  Questions r e g a r d i n g what o p t i o n the u s e r would l i k e p r i n t e d on a c o n v e r s a t i o n a l t e r m i n a l .  are  An example o f one  such s e s s i o n w i t h the i n t e r a c t i v e system i s g i v e n i n Appendix 2.  B.  A graph o f the d a t a p o i n t s , knots and f i t t e d produced, on the g r a p h i c s t e r m i n a l . the  hardcopy  that can be produced  curve are  This i s i d e n t i c a l to from i t a s , f o r example,  t h a t o f F i g u r e 3. 3.  A hardcopy  p l o t o f the graph can be produced  T h i s p l o t i s i d e n t i f i e d by a t i t l e a run number  'n'  which  the p l o t .  n-th  session.  A hardcopy p r i n t o u t as g i v e n i n T a b l e I which to  as i n p u t . a n d  i n d i c a t e s t h a t i t i s the  hardcopy o f the c u r r e n t t e r m i n a l 4.  specified  upon.request.  corresponds  I t can be matched t o the p l o t by the  and run number.  The hardcopy  printout  title  c o n t a i n s the  following information: the  a b s c i s s a e and o r d i n a t e s of the d a t a p o i n t s  (the  o r i g i n a l i n p u t t o the system); the  fitted  o r d i n a t e s (the a p p r o x i m a t i o n to the  o r d i n a t e s by the system); the  residuals  and.the the of the  fitted  (the d i f f e r e n c e between the o r d i n a t e s ordinates);  l e a s t squares e r r o r  (the square r o o t of the  sum  the squares o f the r e s i d u a l s ) ; l o c a t i o n o f the knots  6 , 6 , .. ., <S  and  the c o e f f i c i e n t s o f t h e p i e c e w i s e continuous (s^x)  : i = 0, 1, ..., k}  between knot p a i r  polynomial  between each knot.  0^,6^]  That i s ,  the p o l y n o m i a l c o e f f i c i e n t  i s the v a l u e f o r c . . i n m .(x)  5.2  =  I  j=0 A —A  c  (x - 6 ) J  j  .  -L  Use o f t h e System  To demonstrate  the use o f the i n t e r a c t i v e system a simple s e t  of 11 d a t a p o i n t s was chosen. (ignoring  As c a n be seen from t h e p l o t i n F i g u r e 3  t h e f i t t e d c u r v e f o r t h e moment) t h e eye tends to approximate t h e  data w i t h the f o l l o w i n g  curve:  .9  .7 .6 .5 .4 .3 .2 .1 0  One would l i k e  .1  .2  .3  t o manipulate  data w i t h the above curve.  .4  splines  ,9  so t h a t  1,  they a l s o approximate  this  -  The f i r s t  40 -  attempt was made to approximate  two u n i f o r m l y spaced i n t e r n a l knots over the  '[0,1] .  p l o t i n F i g u r e 3 and. l e a s t squares e r r o r o f  the d a t a w i t h As can be seen  .1104,  from  the r e s u l t  wasn.not s u i t a b l e .  A second attempt was made by moving the two i n t e r n a l knots f u r t h e r apart to  .25  and  .75  respectively.  a p p r o x i m a t i o n g i v e n i n F i g u r e 5, and T a b l e I .  This r e s u l t e d i n the The r e s u l t s a r e somewhat  worse than the p r e v i o u s a p p r o x i m a t i o n (the l e a s t squares e r r o r was -.1574  as compared t o  .1104).  Consequently t h i s a p p r o x i m a t i o n was  e l i m i n a t e d and t h e p r e v i o u s u n i f o r m l y spaced knot s e t r e t a i n e d .  F u r t h e r attempts were made t o produce a b e t t e r a p p r o x i m a t i o n by moving t h e two knots c l o s e r t o g e t h e r . the  When i t became apparent  that  l e a s t squares e r r o r had been reduced adequately w i t h t h e knots  located at  14  and  .6  respectively  ( p l o t and r e s u l t s n o t g i v e n ) , these  v a l u e s were g i v e n t o the m i n i m i z a t i o n procedure t o f i n d t h e o p t i m a l knot locations.  The r e s u l t s a r e g i v e n i n F i g u r e 4 and T a b l e I I .  z a t i o n was t e r m i n a t e d by the c o n s t r a i n t s on the k n o t s . l e a s t s q u a r e s , e r r o r was.reduced s i g n i f i c a n t l y .0544).  (the f i n a l  The o p t i m i -  However, the e r r o r was  Presumably - t h i s i s the b e s t a p p r o x i m a t i o n p o s s i b l e w i t h two  knots.  The next s t e p i n the procedure would then be t o i n c r e a s e the number o f knots to t h r e e .  T h i s was done and the i n i t i a l  results of  t h r e e u n i f o r m l y spaced i n t e r n a l knots a r e shown i n . F i g u r e .5. and T a b l e I I . The i n t e r e s t i n g t h i n g t o note i s that the l e a s t squares e r r o r and approxi m a t i o n a r e i d e n t i c a l to t h a t i n T a b l e I .  T h i s i s because t h e t h i r d  -  41 -  knot l o c a t e d a t the p o i n t of symmetry Is i n e r t and does not to the a p p r o x i m a t i o n at a l l . between  .5  and  .75  contribute  Consequently the p i e c e w i s e p o l y n o m i a l  i s merely the p o l y n o m i a l between , .25  and  .75  shifted.  S i n c e i t was was  increased to four.  internal the  u s e l e s s to c o n t i n u e w i t h t h r e e knots the number The i n i t i a l  knots gave the r e s u l t  i n F i g u r e 6.  o p t i m i z e d r e s u l t w i t h 2 knots  compared .with  .0544).  approximation with four uniform This r e s u l t  i s better  ( l e a s t squares e r r o r o f  .0247  than as  But the curve i n the end r e g i o n s , a l t h o u g h s m a l l e r ,  c o n t a i n s more o s c i l l a t i o n s .  The next s t e p i s to v a r y the knots i n t e r a c t i v e l y . adjustment  of the f i r s t , and l a s t knots towards the b o u n d a r i e s  s i g n i f i c a n t l y b e t t e r r e s u l t s , which t e r m i n a t e d around Next, adjustment continuously. .49999  and  The movement of the two knots was  .50001  as i t was  and remaining f a i r l y  are  spectacular.  and  produced .75 . results  f i n a l l y terminated at  f e l t t h a t they were c o a l e s c i n g t o o much  constant).  squares m a t r i x was  The r e s u l t s  of t h i s  still  only  approximation  As can be seen from T a b l e VI and F i g u r e 7 the l e a s t  squares e r r o r was  (  .25  o f the second and t h i r d knots produced b e t t e r  (although the c o n d i t i o n number o f the l e a s t 5.96  First,  .0000043  It i s interesting  i n F i g u r e s 5 and 7.  and the p l o t  resembles  the one  expected.  t o note the s i m i l a r i t y of the knot  Although F i g u r e 7 r e p r e s e n t s two  locations  n e a r l y e q u a l knots  at the c e n t e r the d i f f e r e n c e i n the:approximations o b t a i n e d i n d i c a t e s t h a t i t i s not n e c e s s a r i l y a good s t r a t e g y t o r e p l a c e two by one knot and r e d u c i n g the o r d e r o f the system by  one.  c o a l e s c i n g knots  0  - 42 -  T h i s example demonstrates to approximate  data.  the power o f the i n t e r a c t i v e  system  In p a r t i c u l a r , the r e s u l t s o f the i n t e r a c t i v e  procedure on the f o u r knot approximation produced such i m p r e s s i v e r e s u l t s t h a t automatic refinement was  unnecessary.  FIGURE 3  - 44 TEST DATA - 2 NON- UNIFORM KNOTS ABSCISSAE ORDINATES FITTED OROINATES RESIDUALS 0.0 0.0 -7.871840E- 03 7.8718U0E-03 1.000000E-01 0.0 3. 829566E-02 -3.829566E-02 2.000000E-01 0.0 -4.988915E- 02 4.988915.E-02 3.000000E-01 0.0 -2. 5950 19E-02 2.595019E-02 4.000000E-01 L000000E-01 1. 877502E-01 -8.775020E-02 5.000000E-01 5.000000E-01 5.000004E- 01 -4.619360E-07 6.000000E-01 9.000000E-01 8.122501E- 01 8.774978E-02 7.000000E-01 1.000000E 00 L025949E 00 -2. 594873E-02 8.000000E-01 LOOOOOOE 00 L049888E 00 -4.988807E-02 9.000000E-01 LOOOOOOE 00 9. 617022E- 01 3. 829774E-02 LOOOOOOE 00 LOOOOOOE 00 L007872E 00 -7.872522E-03 THE LEAST SQUARES ERROR IS L574225E-01 KNOT LOCATION POLYNOMIAL POWER POLYNOMIAL COEFFICIENT 0 0.0 0 -7.871840E-03 1 L979496E 00 2 -1.940852E 01 3 4.230307E 01 2,500000E-01 0 -6.504732E-02 1 2.070669E-01 2 L231876E 01 3 -1.642503E 01 7.500000E-01 0 L065044E 00 1 2.070355E-01 2 -1.231879E 01 3 4.230357E 01 LOOOOOOE 00  TEST DATA - 2 NON-UNIFORM KNOTS  Table I  - -45 TEST DATA - 2 KNOTS OPTIMIZED ABSCISSAE ORDINATES FITTED ORDINATES RESIDUALS 0.0 Oo 0 -7.379014E-03 7.379014E-03 0.0 1. 000000E-01 2.207708E-02 •2. 207708E-02 0.0 2. 000000E-01 1.452189E-02 1.452189E-02 0.0 3. OCOOOOE-01 1. 511128E-02 1. 51 1 128E-02 1.000000E-01 4. 000000E-01 1.223733E-01 •2.237328E-02 5.000000E-01 5. 000000E-01 4.999992E-01 8.030709E-07 9.000000E-01 6. O00O0OE-O1 8.776275E-01 2.237248E-02 7. 000000E-01 000000E 00 015110E 00 •1. 51 1062E-02 8. 000000E-01 •000O00E 00 014520E 00 •1 .452056E-02 9. 00000OE-O1 000000E 00 779216E-01 2.207839E-02 1. 000000E 00 .00O000E 00 007378E 00 •7.378042E-03 5. i,43563E-02 THE L EAST SQUARES ERROR IS LOCATION POLYNOMIAL POWER KNOT POLYNOMIAL COEFFICIENT 0.0 0 0 -7.379014E-03 1 9.650481E-01 2 -8.405960E 00 3 1.701073E 01 4.999000E-01 0 4.994677E-01 1 5.313705E 00 2 1.710503E 01 3 -5.700724E 04 5.001000E-01 0 5.005385E-01 1  1.000000E 00  2 3  TEST DATA - 2 KNOTS OPTIMIZED  Table I I  '\ s r.1UU r\ ry 3.3~\ *IOOCOI1  r-  "A  -1.710493E 1.701070E  01 01  FIGURE 4  - 47 TEST DATA - 3 UNIFORM KNOTS AESCISSAE ORDINATES FITTED ORDINATES RESIDUALS 0. 0 0 .0 - 7 . 872656E- 03 7. 872656E- 03 0 .0 1. 000000E- 01 3. 829633E- 02 - 3 . 829633E- 02 2. OOOOOOE- 01 0 .0 - 4 . 988962E- 02 4. 988962E- 02 3. OOOOOOE- 01 0 .0 -2. 595067E- 02 2. 595067E- 02 4. OOOOOOE- 01 1.OOOOOOE- 01 1. 877483E- 01 - 8 . 774823E- 02 5. OOOOOOE- 01 5 .OOOOOOE- 01 4. 999984E- 01 1. 564622E- 06 6. OOOOOOE- 01 9oOOOOOOE- 01 8. 122493E- 01 8. 775061E- 02 7. OOOOOOE- 01 1•OOOOOOE 00 1. 025949E 00 -2. 594928E- 02 8. OOOOOOE- 01 1 .OOOOOOE 00 -4. 988800E- 02 1. 049888E 00 9. OOOOOOE- 01 1.OOOOOOE 00 9. 617013E- 0 1 3„ 829866E- 02 1- OOOOOOE 00 1.OOOOOOE 00 1 .007872E 00 - 7 . 871866E- 03 THE LEAST SQUARES ERROR IS 1.574225E-01 KNOT LOCATION P O L Y N O M I A L POWER POLYNOMIAL COEFFICIENT 0 0.0 0 -7.872656E-03 1 1.979532E 00 2 -1.940877E 01 3 4.230345E 01 2.500000E-01 0 -6.504667E-02 1 2.070408E-01 2 1.231880E 01 3 -1.642502E 01 5.000000E-01 0 4.999984E-01 1 3.286754E 00 2 5.340576E-05 3 -1.642529E 01 7.500000E-01 0 1.065044E 00 1 2.070377E-01 2 -1.231889E 01 3 4.230394E 01 1.OOOOOOE 00  TEST DATA - 3 UNIFORM KNOTS  Table I I I  FIGURE 5  FIGURE 6  - 50 -  TEST DATA - 4 NON-UNIFORM KNOTS ABSCISSAE ORDINATES FITTED ORDINATES RESIDUALS 0.0 0.0 3.96103 1E- 08 -3.961031E-08 1. 000000E- 01 0.0 -1.974404E- 07 1.974404E-07 2. 000000E- 01 0.0 2. 697052E-07 -2. 697052E-07 3. OOOOOOE-01 0.0 2.083834E- 07 -2.083834E-07 4. OOOOOOE-01 1.000000E-01 9.999895E- 02 1.032065E-06 5. O000OOE- 01 5.000000E-01 5.000019E- 01 -1.941880E-06 6. OOOOOOE-01 9.000000E-01 8. 999965E-01 3.502471E-06 7.000000E- 01 1.000000E 00 1.OOOOOOE 00 -4.593458E-07' 8.000000E- 01 1.000000E 00 9.99999 1E- 01 8.458737E-07 9. 0 00000E-01 1.OOOOOOE 00 9.999999E- 01 7.450581E-08 1. 0O0OOOE 00 1.000000E 00 9.999999E- 01 1.192093E-07 THE LEAST SQUARES ERROR IS 4 266889E-06 KNOT LOCATION POLYNOMIAL POWER POLYNOMIAL COEFFICIENT 0 0.0 0 3.961031E-08 1 -1.450448E-02 2 2.175139E-01 3 -7.249289E-01 2.500000E-01 0 -1.357395E-03 1 -4.168792E-02 2 -3.261279E-01 3 3.405853E 01 4.999900E-01 0 4.999400E-01 1 6. 180718E 00 2 2.521675E 01 3 -8.394170E 05 5.000100E-01 0 5.000489E-01 1 6.180779E 00 2 -2.521672E 01 3 3.405827E 01 7.500000E-01 0 1.001356E 00 1 -4.164546E-02 2 3.259857E-01 3 -7.245331E-01 1.OOOOOOE 00  TEST DATA - 4 NON-UNIF0RM  T a b l e IV  KNOTS'  FIGURE 7  - 52 -  5.3  T i t a n i u m Heat. Data  de Boor and R i c e [ 5 ] , [6] p r e s e n t e d a s e t o f d a t a o b t a i n e d from a heat experiment i n v o l v i n g the  f i x e d knot program  Titanium.  T h i s d a t a was  used to  demonstrate  u s i n g a knot s e t w i t h f i v e e q u a l l y spaced knots  [5, p. 18] and the v a r i a b l e knot program which i n v o l v e d the o p t i m i z a t i o n of t h i s f i x e d knot s e t .  F u r t h e r s t u d i e s i n v o l v e d the o p t i m i z a t i o n o f  a non-uniform knot s e t c o n t a i n i n g f i v e k n o t s .  One o f the b e s t ways to t e s t a system i s to t r y i t on p r e v i o u s l y done r e s u l t s and check the answers. was  Initially  the i n t e r a c t i v e  system  done w i t h the u n i f o r m knot s e t s p e c i f i e d i n de Boor and R i c e [5, p. 18].  The r e s u l t s and p l o t are g i v e n i n T a b l e V and F i g u r e 8.  These  results  compare f a v o r a b l e w i t h those i n de Boor and R i c e (the r e s i d u a l was as compared t o t h e i r the  1.24).  1.16  However, as can be seen from the p l o t ,  r e s u l t i n g a p p r o x i m a t i o n i s none too s a t i s f a c t o r y .  It effect.  i s here t h a t the power o f the i n t e r a c t i v e system comes i n t o  A f t e r moving the knots so t h a t the p l o t produced i s more a c c u r a t e ,  b e t t e r r e s u l t s are o b t a i n e d . indicate initial  substantial  improvement.  The r e s u l t i n g knot s e t was  g u e s s s f o r knot o p t i m i z a t i o n g i v i n g  and F i g u r e 10. (.092  The r e s u l t s g i v e n i n T a b l e VI and F i g u r e 9  final  r e s u l t shown i n T a b l e V I I  The r e d u c t i o n i n the l e a s t squares e r r o r was  as compared t o  1.16)  ment a l l o w e d f o r the d e r i v i n g  used i n an  substantial  l a r g e l y because the i n t e r a c t i v e knot of accurate s t a r t i n g values.  place-  -.53 " TITANIUM HEAT DATA - 5 UNIFORM KNOTS 1 ABSCISSAE ORDINATES FITTED ORDINATES RESIDUALS 6. 440000E- 01 5. 950000E 02 6. 235025E- 01 2. 049747E- 02 6. 050000E 02 6.220000E- 01 6. 425120E- 01 -2. 051 199E- 02 6. 150000E 02 6.3 800 00E-01 6. 516960E- 01 -1. 369604E- 02 6. 250000E 02 6. 490000E- 01 6. 534591E- 0 1 -4. 459172E- 03 6. 350000E 02 6. 520000E- 01 6. 502057E- 01 1. 794243E- 03 6. 450000E 02 6. 390000E- 01 6„ 443403E- 01 -5. 340301E- 03 6. 550000E 02 6. 460000E- 01 6. 382672E- 01 7.732764E- 03 6. 570000E- 01 6. 343911E- 0 1 6. 650000E 02 2. 260887E- 02 6. 750000E 02 6. 520000E- 01 6. 351163E- 01 1. 688367E- 02 6. 850000E 02 6. 550000E- 01 6. 420718E- 01 1. 292809E- 02 6. 950000E 02 6. 640000E- 01 6. 537848E- 01 1. 021518E- 02 7.050000E 02 6. 630000E- 01 6. 680065E- 0 1 -5. 006507E- 03 7. 150000E 02 6. 630000E- 01 6. 824884E- 01 -1. 948840E- 02 7.250000E 02 6. 680000.E-01 6. 949821E- 01 -2. 698214E- 02 7.350000E 02 6. 760000E- 01 7. 032389E- 01 -2. 723893E- 02 7.450000E 02 6. 760000E- 01 7. 050105E- 01 -2. 901048E- 02 7.550000E 02 6. 860000E- 01 6. 980482E- 01 -1. 204824E- 02 7. 650000E 02 6. 790000E- 01 6. 815894E- 01 -2. 589412E- 03 7. 750000E 02 6.780000E- 01 6. 608155E- 01 1. 718443E- 02 7. 850000E 02 6. 830000E- 01 6. 4 2 39 3 3E01 4. 060671E- 02 7,950000E 02 6.940000E- 01 6. 329899E- 01 6. 101005E- 02 8. 050000E 02 6. 990000E- 01 6. 392726E- 0 1 5. 972740E- 02 8. 150000E 02 4.209192E- 02 7. 100000E- 01 6. 679080E- 01 8.250000E 02 7.300000E- 01 4 436404E- 03 7. 255636E- 01 8o 350000E 02 7. 630000E- 01 8. 189063E- 01 -5o 590630E- 02 8.450000E 02 8. 120000E- 01 9. 505556E- 01 -1. 385556E- 01 8. 550000E 02 9. 070000E- 01 1. 106939E 00 -1. 999393E- 01 8.650000E 02 1. 270436E 00 - 2 , 264375E- 01 1. 0440O0E 00 8. 750000E 02 1. 336000E 00 1. 423430E 00 -8. 743048E- 02 8. 850000E 02 3. 326998E- 01 1.881000E 00 1. 548299E 00 8.950000E 02 2. 169000E 00 1. 627424E 00 5. 415753E- 01 9. 050000E 02 2. 075000E 00 4. 318160E- 01 1. 643184E 00 9. 150000E 02 1. 598000E 00 1. 583930E 00 1. 406908E- 02 9. 250000E 02 1. 461898E 00 -2. 508975E- 01 1. 21 1000E 00 9. 350000E 02 9. 160000E- 01 1. 295290E 00 -3. 792901E- 01 9. 450000E 02 7.460000E- 01 1. 10231 1E 00 -3. 563114E- 01 9. 550000E 02 6. 720000E- 01 9. 011649E- 01 -2. 291650E- 01 9, 650000E 02 6.270000E- 01 7. 100576E- 01 -8. 305758E- 02 6„ 150000E- 01 9. 750000E 02 5. 471910E- 01 6.780899E- 02 9. 850000E 02 4. 307705E- 01 6.070000E- 01 1. 762294E- 01 9, 950000E 02 6. 060000E- 01 3. 790006E- 01 2, 2&9993E- 01 6.090000E- 01 4. 022449E- 01 1. 005000E 03 2. 067550E- 01 4. 795058E- 01 6. 030000E- 01 1. 015000E 03 1. 234941E- 01 6.010000E- 01 5. 819451E- 0 1 1. 025000E 03 1. 905493E- 02 6.030000E- 01 6. 807239E- 01 -7. 772392E- 02 1. 035000E 03 6.010000E- 01 7. 470052E- 01 -1. 460052E- 01 1. 045000E 03 6. 1 10000E-01 1. 055000E 03 7. 519506E- 01 -1. 409506E- 01 6.010000E- 01 6. 667216E- 01 -6o 572163E- 02 1. 065000E 03 6.080000E- 01 4.624799E- 01 . 1. 455200E- 01 1. 075000E 03 THE LEAST SQUARES ERROR IS 1.157334E 00 KNOT LOCATION POLYNOMIAL POWER POLYNOMIAL COEFFICIENT 0 5.950000E 02 0  - 54 -  6.750000E 02  7.550000E 02  8.350000E 02  9.050000E 02  9.950000E 02  0 1 2 3  6.235024E-01 2.472371E-03 -6.1 14945E-05 4.007443E-07  0 1 2 3  6.351163E-01 3.827438E-04 3.502920E-05 -3.747538E-07  0 1 2 3  6. 980482E-01 -1.207871E-03 -5.49116;E-05 1.1 11178E-06  0 1 2 3  8.189063E-01 1. 134087E-02 2.117710E-04 -2.936617E-06  0 1 2 3  1.643183E 00 -2.179519E-03 -4.049190E-04 3.034045E-06  0 1 2 3  3.790006E-01 -1.337662E-03 4.142732E-04 -4.806359E-06  1.075000E 03  TITANIUM HEAT DATA - 5 UNIFORM KNOTS  Table V  - 56 TITANIUM HEAT DATA 5 NON-UNIFORM KNOTS ABSCISSAE ORDINATES FITTED ORDINATES RESIDUALS 5.950000E 02 6 .440000E -01 6 .252043E -01 1.879563E -02 6.050000E 02 6 .220000E -01 6 . 338547E-01 -1. 185482E -02 6.150000E 02 6 .380000E -01 6 .407921E -01 -2. 792177E--03 6.250000E 02 6 .490000E -01 6 .462172E -01 2.782739E -03 6.350000E 02 6 .520000E -01 6 .503309E -0 1 1. 669 100E--03 6.450000E 02 6 . 390000E-01 6 .533335E -01 -1.433358E -02 6.550000E 02 6 .460000E -01 6 .554263E -0 1 -9. 426299E--03 6.650000E 02 6 .570000E -01 6 . 568095E-01 1.904927E -04 6.750000E 02 6 . 520000E'-01 6 -576843E -01 -5. 684260E--03 6.850000E 02 6 .550000E--01 6 .582511E -01 -3.251 180E--03 6.950000E 02 6 .640000E--01 6 .587108E -01 5.2891 15E--03 7.050000E 02 6 . 630000E--01 6,. 592643E--01 3. 735680E--03 7.150000E 02 6 -630000E--01 6 .601120E--01 2. 888002E--03 7.250000E 02 6 . 680000E--01 6,. 614549E--01 6. 545052E--03 7.350000E 02 6 .760000E--01 6 .634936E--01 1.250640E--02 7.450000E 02 6,.760000E--01 6,» 664289E-01 9. 571020E--03 7.550000E 02 6 . 860000E--01 6,. 704616E--01 1.553836E--02 7.650000E 02 6.. 790000E--01 6,. 757923E-01 3.207672E--03 7.750000E 02 6..780000E--01 6..826219E--01 -4.621979E--03 7.850000E 02 6,• 830000E-01 6,. 9 115 11E-01 -8. 151092E--03 7.950000E 02 6,.940000E--01 7,.015807E--01 -7.580712E--03 8.050000E 02 6,.990000E--01 7., 14 1112E--01 - 1. 511123E-02 8.150000E 02 7,, 100000E--01 7..289435E--01 -1.894355E--02 8.250000E 02 7,. 300000E-01 7,• 462784E-01 -1. 627839E-•02 8.350000E 02 -A 1 . A -1 7,> 630000E-01 7».6 6 3165E• ~J I u J *f I U " vj J v/ 1 8.450000E 02 8.. 120000E•01 7o 908105E--01 2.118939E-•02 8.550000E 02 9..070000E-•01 8,•572057E--01 4.979425E-•02 8.650000E 02 1..044000E 00 1. 038639E 00 5.359702E-•03 8.750000E 02 1.336000E 00 1.• U02930E 00 -6.693059E-•02 8.850000E 02 1. 881000E 00 1. 859838E 00 2.116160E- 02 8.950000E 02 2.. 169000E 00 2. 161066E 00 7.933423E-•03 9.050000E 02 2. 075000E 00 2. 064366E 00 1.063279E- 02 9.150000E 02 1. 598000E 00 624727E 00 -2.672801E-•02 r. 9.250000E 02 1. 211000E 00 1. 185819E 00 2.518102E- 02 9.350000E 02 9. 160000E- 01 9. 078093E-•01 . 8. 190691E-•03 9.450000E 02 7. 460000E- 01 7. 544307E-•01 -8, 430697E- 03 9.550000E 02 6. 720000E- 01 6. 808758E-•01 -8.875843E- 03 9.650000E 02 6. 270000E- 01 6. 432318E- 01 - 1. 623188E-02 9.750000E 02 6 o 150000E- 01 6. 181573E- 01 -3.157331E- 03 9.850000E 02 6. 070000E- 01 6. 028832E-•Q 1 4.116789E- 03 9.950000E 02 6- 060000E- 01 5. 955343E- 01 1.046569E- 02 1.005000E 03 6. 090000E- 01 5. 942356E- 0 1 1. 476442E-02 1.015000E 03 6. 030000E- 01 5. 971119E- 01 5.888034E03 1.025000E 03 6, 0 10000E-01 6. 022885E- 01 -1. 288563E- 03 1.035000E 03 6. 030000E- 01 6. 078902E- 01 -4.890207E- 03 1.045000E 03 6. 010000E- 01 6. 120420E- 01 - 1. 104198E-02 1.055000E 03 6. 110000E- 01 6. 128688E- 01 -1.868859E- 03 1.065000E 03 6. 010000E- 01 6. 084959E- 01 -7.495925E- 03 1.075000E 03 6. 080000E- 01 5. 970479E- 01 1.095206E- 02 THE LEAST SQUARES ERROR IS 1.142650E-01 KNOT LOCATION POLYNOMIAL POWER POLYNOMIAL COEFFICIENT 0 5.950000E 02  I {  - 57 -  8.400000E  8.700000E  9.OOOOOOE  1.075000E  6.252043E-01 9.573922E-04 -9.568968E-06 3.345709E-08  0 1 2 3  7.774121E-01 2.293389E-03 1.502197E-05 1.244829E-05  0 1 2 3  1 . 1 9 5 8 3 7 E 00 3.680510E-02 1o135369E-03 -4.252815E-05  0 1 2 3  2 . 1 7 3 5 5 8 E 00 -9.898979E-03 -2.692169E-03 6.085630E-05  0 1 2 3  1.385562E 00 -4.455817E-02 9.592120E-04 -7.467897E-06  0 1 2 3  6.600300E-01 -3.667146E-03 6.306361E-05 -3.125027E-07  02  02  02  9.200000E 02  9.600000E  0 1 2 3  02  03  TITANIUM HEAT DATA - 5 NON-UNIFORM KNOTS  T a b l e VI  - 58 -  O  - 59 TITANIUM HEAT DATA - 5 KNOTS OPTIMIZED ABSCISSAE ORDINATES FITTED ORDINATES RESIDUALS 5..950000E 02 6,.440000E--01 225190E--0 1 6. 2. 148103E--02 6., 050000E 02 6,.220000E--01 6..327376E--01 -1. 073763E-•02 6., 150000E 02 6., 380000E-01 6. 408629E--01 -2. 862922E-•03 6., 250000E 02 6..490000E--01 6..471329E--01 1. 867094E-•03 6..350000E 02 6.• 520000E-01 6. 517856E--01 2. 144140E-•04 6..450000E 02 6..390000E--01 6,.550590E--01 -1. 605903E-•02 6. 550000E 02 6., 460000E-01 6. 571912E--01 -1. 1 191 16E-•02 6.. 650000E 02 6., 570000E-01 6. 584202E--01 -1. 420241E-•03 6. 750000E 02 6. 520000E--01 6. 58984 1E--01 -6. 984055E- 03 6. 850000E 02 6., 550000E-01 6. 591210E--01 -4. 121054E-•03 6. 950000E 02 6. 640000E-•01 6. 590686E-•01 4. 931286E- 03 7. 050000E 02 6..630000E--01 6. 590654E-•01 3. 934622E- 03 7. 150000E 02 6. 630000E--01 6. 593490E-•0 1 3. 65094 1E-03 7. 250000E 02 6. 680000E-•01 6. 601578E-•01 7. 842168E-•03 7. 350000E 02 6. 760000E-•01 6. 617296E-•01 1. 427035E- 02 7. 450000E 02 6. 760000E-•01 6. 643026E-•01 1. 169730E-•02 7. 550000E 02 6. 860000E-•01 6. 681147E-•01 1. 788527E- 02 7o 650000E 02 6. 790000E-•01 6. 734042E-•01 5. 595822E-•03 7. 750000E 02 6. 780000E-•01 6. 804088E-•01 -2. 408907E- 03 7o 850000E 02 6. 830000E-•01 6. 893667E-•01 -6 =366715E-•03 7. 950000E 02 6. 940000E-•01 7. 005159E-•01 -6. 515928E- 03 8. 050000E 02 6. 990000E-•01 7. 1409U6E-•01 -1. 509461E- 02 8. 150000E 02 7. 100000E-•01 7. 303407E-•01 -2. 034071E- 02 8. 250000E 02 7. 300000E-•01 7. 494920E-•01 -1. 949206E- 02 8, 350000E 02 7 630000E-•0 1 7^ 717870E-•0 1 -8. 786 9 95E-03 8. 450000E 02 8. 120000E-•01 7. 992001E-•01 1. 279991E- 02 8. 550000E 02 9. 070000E- 01 8. 665996E-•01 4. 040042E- 02 8. 650000E 02 1. 044000E 00 1. 038002E 00 5. 996864E- 03 8. 750000E 02 1. 336000E 00 1. 378159E 00 -4. 215827E- 02 8. 850000E 02 1. 881000E 00 1 .851600E 00 2. 939950E- 02 8. 950000E 02 2. 169000E 00 2. 178109E 00 -9. 109914E- 03 9, 050000E 02 2. 075000E 00 2. 067245E 00 7. 754855E- 03 9. 150000E 02 1. 598000E 00 1. 616635E 00 -1. 863577E- 02 9. 250000E 02 1. 21 1000E 00 1. 193341E 00 1. 765861E- 02 9. 350000E 02 9. 160000E- 01 9. 150418E- 01 9. 581815E- 04 9. 450000E 02 7. 460000E- 01 7. 509134E- 01 -4. 913419E- 03 9. 550000E 02 6. 720000E- 01 6. 683230E- 01 3. 676936E- 03 9. 650000E 02 6. 270000E- 01 6. 346373E- 01 -7. 637367E- 03 9. 750000E 02 6. 150000E- 01 6. 188952E- 0 1 -3. 895170E- 03 9. 850000E 02 6. 070000E- 01 6. 086131E- 01 -1. 613097E- 03 9. 950000E 02 6. 060000E- 01 6. 027005E- 01 3. 299463E- 03 6. 090000E- 01 1. 005000E 03 6. 002449E- 01 8. 755121E- 03 6. 030000E- 01 1. 015000E 03 6. 00 3 335E-01 2. 666444E- 03 1. 025000E 03 6. 010000E- 01 6, 020537E- 01 -1. 053706E- 03 1. 035000E 03 6. 030000E- 01 6. 044927E- 01 -1. 492694E- 03 1. 045000E 03 6. 010000E- 01 6. 067377E- 01 -5. 737711E- 03 6. 1 10000E-01 1. 055000E 03 6. 078761E- 01 3. 123831E- 03 1. 065000E 03 6. 010000E- 01 6. 069952E- 01 -5. 995244E- 03 1. 075000E 03 6. 0800O0E- 01 6. 031824E- 01 4. 817545E- 03 THE LEAST SQUARES ERROR IS 9.286332E-02 KNOT LOCATION POLYNOMIAL POWER POLYNOMIAL COEFFICIENT 0 5.950000E 02 S  - 60 -  8.395486E 02  8.733201E 02  8.989514E 02  9.179270E 02  9.681765E 02  1.075000E 03  0 1 2 3  6.225189E-01 1.134473E-03 -1.165708E-05 3.967317E-08  0 1 2 3  7.830166E-01 2.551379E-03 1 .744409E-05 1.084204E-05  0 1 2 3  1.306674E 00 4.082611E-02 1.115898E-03 -5.281300E-05  0 1 2 3  2.196890E 00 -6.059237E-03 -2.945114E-03 6.665658E-05  0 1 2 3  1.476898E 00 -4.5826 13E-02 8.494323E-04 -5.438899E-06  0 1 2 3  6.288878E-01 -1.658810E-03 2.952565E-05 -1. 521178E-07  TITANIUM HEAT DATA - 5 KNOTS OPTIMIZED  Table VIII  FIGURE 10  - 62 -  5.4  The Bug  The use o f s p l i n e a p p r o x i m a t i o n i n computer-aided i s a r e l a t i v e l y unexplored a r e a . s e c t i o n demonstrates  design  The Volkswagen d a t a p r e s e n t e d i n t h i s  p o s s i b l e a p p l i c a t i o n s of s p l i n e s i n the area.  the a p p r o x i m a t i o n demonstrates  Also  the use o f a l a r g e knot s e t .  The. i n i t i a l a p p r o x i m a t i o n o f 18 e q u a l l y spaced knots  produced  a Volkswagen w i t h a dented hood as can be seen i n F i g u r e 11.  Inter-  a c t i v e m a n i p u l a t i o n o f t h e knots produced  resemblance  to a Volkswagen shown i n  the more r e a s o n a b l e  Figure 12 and T a b l e VIII..  S i n c e the d e s i r e d  r e s u l t s were based on r e p r e s e n t i n g the d a t a a c c u r a t e l y by s i g h t  rather  than m i n i m i z i n g the l e a s t squares e r r o r the knots were not o p t i m i z e d .  - 64 THE BUG ABSCISSAE 3. 500000E-•01 3. 525000E- 01 3. 550000E-•01 4. 300000E-•01 4. 700000E-•01 5. 000000E-•01 5. 500000E-•01 6.OOOOOOE-•01 7. OOOOOOE-•01 8. OOOOOOE-•01 9. OOOOOOE-•01 1. OOOOOOE 00 1. 100000E 00 1. 200000E 00 1, 300000E 00 1. 400000E 00 1. 500000E 00 1. 600000E 00 1 .700000E 00 1. 800000E 00 1. 900000E 00 2. OOOOOOE 00 2. 100000E 00 2. 200000E 00 2. 300000E 00 2, 4U0O0OE 00 2. 500000E 00 2. 600000E 00 2. 700000E 00 2. 800000E 00 2. 9C0000E 00 3. OOOOOOE 00 3. 100000E 00 3. 200000E 00 3. 350000E 00 3. 400000E 00 3. 500000E 00 3. 600000E 00 3. 700000E 00 3. 800000E 00 3. 900000E 00 4. OOOOOOE 00 4. 100000E 00 4. 2C0000E 00 4. 300000E 00 4. 400000E 00 4. 500000E 00 4. 6O0000E 00 4. 700000E 00 4. 800000E 00 4. 900000E 00 5. OOOOOOE 00  ORDINATES 3 .OOOOOOE-•01 4.OOOOOOE-•01 5.500000E-•01 6. 500000E-•01 8.OOOOOOE--01 9.500000E-•01 1.060000E 00 1.200000E 00 1.250000E 00 1•290000E 00 1.3.00000E 00 1„430000E 00 1.500000E 00 1.560000E 00 1.620000E 00 1.675000E 00 1•720000E 00 1.750000E 00 1.780000E 00 1.810000E 00 1.840000E 00 1.860000E 00 1„890000E 00 1.915000E 00 1.930000E 00 1.950000E 00 1.960000E 00 1.970000E 00 1.980000E 00 1 .990000E 00 2 .OOOOOOE 00 2 .005000E 00 2 .005000E 00 2 .005000E 00 2 .005000E 00 2 . 100000E 00 2.250000E -00 2 .430000E 00 2 .600000E 00 2 .750000E 00 2 .890000E 00 2 .920000E 00 2 .950000E 00 2 .970000E 00 3 .OOOOOOE 00 3.020000E 00 3 .040000E 00 3 .050000E 00 3 .055000E 00 3 .065000E 00 3.075000E 00 ' 3.080000E 00  ORDINATES RESIDUALS 4.035096E-•01 -1. 035095E- 01 4.113406E-•0 1 -1. 134065E- 02 4.192352E-•01 1. 307648E-•01 6.758505E-•01 -2. 585047E- 02 8.169750E-•01 -1. 697499E- 02 9.183491E-•01 3. 165080E- 02 1 .068680E 00 -8. 680243E- 03 1. 18203 1E 00 1. 796837E- 02 1.275474E 00 -2. 547405E- 02 1.276578E 00 1. 342188E- 02 1 .304118E- 00 -4. 118513E-03 1.421022E 00 8. 977752E- 03 1.509212E 00 -9. 212483E- 03 1.564281E 00 -4. 281554E- 03 1.616894E 00 3. 105875E- 03 1 .668297E 00 6.702900E- 03 1.713498E 00 6.501570E- 03 1.752367E 00 -2. 367944E- 03 1.785924E 00 -5. 924519E- 03 1•815188E 00 -5. 188584E-03 1.841183E 00 -1, 182909E- 03 1 .864781E 00 -4. 782557E-03 1.887129E 00 2. 871351E- 03 1.908055E 00 6.944049E- 03 1.927378E 00 2. 621699E- 03 1 .944908E 00 5. 091667E-•03 1.960458E 00 -4. 580617E- 04 1 .973843E 00 -3. 842760E- 03 1 .984876E 00 -4. 87637 1E-03 1.993369E 00 -3. 369596E- 03 1.999138E 00 8. 617891E- 04 2 .002012E 00 2. 987819E- 03 2 .002145E 00 2. 855293E- 03 2 .001301E 00 3. 698572E-03 2 .03421 1E 00 -2. 921 108E- 02 2 .079496E 00 2. 050392E- 02 2 .238708E 00 1. 129068E- 02 2 .428134E 00 1. 865786E- 03 2 .612069E 00 -1. 206875E- 02 2 .767555E 00 -1. 755604E- 02 2 .871636E 00 1. 836338E- 02 2 .916960E 00 3. 039550E-03 2 .946946E 00. 3. 053293E- 03 2 .973553E 00 -3. 552493E- 03 2 .996966E 00 3. 032722E- 03 3 .017381E 00 2. 619695E- 03 3 .034980E 00 5. 019724E- 03 3 .049955E 00 4. 446507E- 05 3 .062496E 00 -7. 496823E- 03 3 .072790E 00 -7. 791314E- 03 3.081028E 00 -6. 028723E- 03 3 .087397E 00 -7. 396754E-03  - 65 5. 100000E 5. 200000E 5. 300000E 5. 400000E 5. 500000E 5. 600000E 5. 700000E 5. 800000E 5. 900000E 6. OOOOOOE 6. 100000E 6. 200000E 6. 300000E 6. 400000E 6. 500000E 6, 600000E 6. 700000E 6. 800000E 6. 900000E 7. OOOOOOE 7. 100000E 7. 200000E 7, 300000E 7. 400000E 7, 500000E 7. 600000E 7. 700000E  00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 "7 00 / « Q A A A A A T? 7. 900000E 00 VJ O W V W S J .LJ 00 8, OOOOOOE 8. 100000E 00 8.200000E 00 8. 300000E 00 8. 400000E 00 So 500000E 00 8.600000E 00 8. 700000E 00 8. 800000E 00 8.900000E 00 9. OOOOOOE 00 9. 100000E 00 9. 200000E 00 9. 300000E 00 9, 400000E 00 9. 500000E 00 9. 550000E 00  3. 090000E 3. 100000E 3. 10000OE 3. 100000E 3. 100000E 3o 100000E 3. 100000E 3. 090000E 3. 080000E 3. 070000E 3. 060000E 3. 050000E 3. 040000E 3. 030000E 3. 010000E 3. OOOOOOE 2. 980000E 2. 950000E 2. 920000E 2. 890000E 2. 850000E 2. 800000E 2. 750000E 2. 690000E 2. 630000E 2. 560000E 2. 490000E  00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00 00  O  AA  fi  T  nnnnw  L « ~T t- W \J \J \J U  2. 340000E 00 2. 250000E 00 2. 160000E 00 2 ©080000E 00 2. OOOOOOE 00 1. 900000E 00 1. 805000E 00 1. 7050C0E 00 1. 605000E 00 1. 500000E 00 1. 375000E 00 1. 240000E 00 1. 030000E 00 8. 300000E--01 6. 400000E--01 4. 450000E--01 2. 500000E--01 2. 000000E--01  -2. 069191E- 03 3. 092069E 00 4. 853774E- 03 3. 095146E 00 3. 283732E- 03 3. 096716E 00 3. 132608E- 03 3. 096868E 00 4. 319567E- 03 3. 095680E 00 6. 755099E- 03 3. 093245E 00 1. 035092E- 02 3. 089649E 00 5. 025774E- 03 3. 084974E 00 6. 905058E- 04 3. 079309E 00 -2. 739966E- 03 3. 072740E 00 - 5 . 259581E- 03 3. 065260E 00 -6. 500497E- 03 3. 056500E 00 - 5 . 999859E- 03 3. 046000E 00 - 3 . 295366E- 03 3. 033295E 00 - 7 . 928025E- 03 3. 017928E 00 2. 999435E 00 ' 5.641840E- 04 2. 645336E- 03 2. 977354E 00 - 1 . 224987E- 03 2. 951224E 00 - 5 . 867698E- 04 2. 920587E 00 5. 003192E-•03 2. 884996E 00 5. 948193E- 03 2. 844051E 00 2. 174579E-•03 2. 797825E 00 3. 533684E- 03 2. 746466E 00 - 1 . 264103E-•04 2. 690125E 00 1. 045644E- 03 2. 628954E 00 - 3 . 100801E-•03 2. 563101E 00 -2„ 719902E- 03 2. 492720E 00 -A O * r\ it i n n c o_vy ~J £. o/i n n c c u nr\ \J T - T 4. \J _ / i _ 1. 036749E- 03 2. 338963E 00 - 5 . 881384E-•03 2. 255831E 00 - 9 . 028815E-•03 2. 169028E 00 7. 202886E-•04 2. 079279E 00 1. 233252E-•02 1, 987667E 00 4. 775889E-•03 1. 895224E 00 2. 018948E-•03 1. 802980E 00 - 5 . 694207E-•03 1. 710693E 00 - 8 . 008420E-•03 1. 613008E 00 - 3 . 283981E--03 1. 503284E 00 374889E 00 1. 102020E-•04 1. 881186E--02 00 1. 221188E 1. 940164E-•03 8 . 038939E 00 1. 487869E--03 -8. 8. 384878E--01 6. 435510E-•03 6. 335645E--01 7, 090468E--03 4. 379095E--01 - 1 . 525517E--02 2. 652552E--01 -03 8., 153468E1..918465E--01  THE LEAST SQUARES ERROR IS 1.897547E-01 KNOT LOCATION POLYNOMIAL POWER 0 3.500000E-01 0 1 2 3 1 6.000000E-01  POLYNOMIAL COEFFICIENT 4.035096E-01 3.119434E 00 5.245779E 00 -2.106841E 01  - 66 -  2  9.000000E-01  3  1.OOOOOOE OO  4  1.2000OOE OO  5  1.400000E OO  6  2.OOOOOOE OO  7  3.OOOOOOE OO  8  3.200000E OO  9  3.400000E OO  10  3.500000E 00  11  3.900000E 00  0 1 2 3  1.182031E 1.792031E -1.055562E 1.979568E  00 00 01 01  0 1 2 3  1.304118E 00 8.034888E-01 7.260439E 00 -3.605061E 01  0 1 2 3  1.421021E 1.174047E -3.554706E 6.331032E  0 1 2 3  1.564281E 00 5.119148E-01 2.439298E-01 -1.019521E 00  0 1 2 3  1.668296E 00 4.870993E-01 -3.677540E-01 1.702099E-01  0 1 3  1.864781E 00 2.299367E-01 -6. 161 657E-02 -3.110615E-02  0 1 2 3  2.002012E 00 1.324511E-02 -1.546146E-01 3.537191E-01  0 1 2 3  2.001300E 00 -6.130829E-03 5.760336E-02 9.639282E 00  0 1 2 3  2.079494E 1.173615E 5.841224E -1.655806E  0 1 2 3  2.238708E 00 1.845140E 00 8.737149E-01 -3.827115E 00  0 1 2 3  2.871635E 00 7.071088E-01 -3.718792E 00 1.180118E 01  00 00 00 00  00 00 00 01  - 67 -  12  4.OOOOOOE OO  13  5.OOOOOOE OO  14  6.OOOOOOE OO  15  7.OOOOOOE OO  16  8.OOOOOOE OO  17  8.500000E OO  0 1 2 3  2.916960E OO 3. 173848E-01 -1.784239E-01 3.148174E-02  0 1 2 3  3.087396E OO 5.497074E-02 -8.398247E-02 1.435216E-02  0 1 2 3  3.072739E OO -6.993294E-02 -4.092169E-02 -7.690954E-02  0 1 2 3  2.884996E OO -3.825288E-01 -2.716408E-01 2.506399E-02  0 1 2 3  2.255880E OO -8.505980E-01 -1.964577E-01 1.721067E-01  0 2 3  1.802979E 00 -9. 179758E-0 I 6. 170338E-02 -1.105849E 00  0 1 2 3  1.221187E -1.685660E -1.597072E 2.289328E  1  18  9.OOOOOOE OO  19  9.550000E OO  THE BUG  Table VIII  00 00 00 00  - 69  -  Section 6  SUMMARY AND  The  FUTURE POSSIBILITIES  r e s u l t s shown i n the p r e v i o u s  section provide  a brief  out-  l o o k i n t o the p o s s i b i l i t i e s of - i n t e r a c t i v e l e a s t squares a p p r o x i m a t i o n u s i n g s p l i n e f u n c t i o n s . There remain, however, many problems and sibilities  for further research  With regards  in this  suggests choosing  6.., l+l  l+l  6. l  spacing and  i n S e c t i o n 2.  i s different.  first  1  - <S  0  -m  6, ,, - 6, k+1 k  For the case o f a u n i f o r m described  the e x t e r n a l k n o t s .  the e x t e r n a l knots so  6  =.6. l -  area.  to the b a s i s f u n c t i o n s ; a p o s s i b i l i t y i s to  more a p p r o p r i a t e methods of choosing p. 73]  pos-  Schultz  find [18,  that  < i < -1  ;  k + l < i < k — —  + m  knot s e t t h i s i s e q u i v a l e n t , to the method However, f o r the non-uniform knot s e t the knot  S c h u l t z ' s c h o i c e i s poor when the boundary knot  i n t e r i o r knot are c o a l e s c i n g as i t causes h i g h peaks i n the  f u n c t i o n s i n v o l v i n g the e x t e r n a l and boundary k n o t s .  T h i s causes a sharp  i n c r e a s e i n t h e . c o n d i t i o n number of the l e a s t squares m a t r i x which should be  The  - a situation  avoided.  second p o s s i b i l i t y  f i n d optimal weights.for  (interconnected with  the f i r s t ) i s to  the b a s i s f u n c t i o n s i n v o l v i n g the  knots so as to reduce the c o n d i t i o n number.  At p r e s e n t  l a s t b a s i s f u n c t i o n are a r b i t r a r i l y m u l t i p l i e d by  external  the f i r s t  (m + 1)  so as  and to  o b t a i n an adequate p o r t i o n of these b a s i s f u n c t i o n s when d e r i v i n g the l e a s t squares m a t r i x .  An  basis  i n t e r e s t i n g sutdy c o u l d be made of  the  - 70c-  theoretical  and p r a c t i c a l importance o f such weights f o r p r o d u c i n g  good bounds f o r the norm of the l e a s t squares m a t r i x .  A t h i r d p o s s i b i l i t y i s to produce a more e f f i c i e n t method of computing  the B - s p l i n e b a s i s f u n c t i o n s u s i n g r e c u r s i v e formulae.  These  methods a r e d e s c r i b e d i n Lyche and Schumaker [TO]. They p r o v i d e subs t a n t i a l r e d u c t i o n s i n the time i n v o l v e d i n e v a l u a t i n g the B - s p l i n e s as i d e n t i c a l terms are not recomputed.  S i n c e the B - s p l i n e b a s i s  function  e v a l u a t i o n uses the most s i g n i f i c a n t amount of time i n the system; derivation  of an improved, more e f f i c i e n t a l g o r i t h m f o r s p l i n e  f u n c t i o n e v a l u a t i o n would be h i g h l y  basis  justified.  With r e g a r d s to the f i x e d knot l e a s t squares problem;  a  b e t t e r e s t i m a t e o f the l e a s t squares norm might a i d i n the a p p r o x i m a t i o n process.  Both de Booraand R i c e and P a t e n t use an i n t e g r a l norm i n t h e i r  a p p r o x i m a t i o n as opposed  t o the d i s c r e t i z e d norm (sum o f the squares  of the r e s i d u a l s ) used i n t h i s system,  de Boor and R i c e [ 5 ] , [6]  e s t i m a t e the norm by the t r a p e z o i d a l r u l e ; hence, i n e f f e c t , a linear interpolation  between f i t t e d d a t a p o i n t s .  F i l o n q u a d r a t u r e (an i n t e r p o l a r y discussed i n h i s thesis.  obtaining  P a t e n t [13] used  q u a d r a t u r e ) , the t h e o r y o f which i s  S i n c e the main i n t e r e s t of t h i s system  was  i n t e r a c t i v e a p p r o x i m a t i o n r a t h e r than l e a s t squares a p p r o x i m a t i o n , the d i s c r e t i z e d norm was  adequate.  With r e g a r d s to the v a r i a b l e knot problemj b e t t e r a l g o r i t h m s are needed  f o r the n o n - l i n e a r l e a s t squares problem.  In p a r t i c u l a r ,  a l g o r i t h m f o r the s p e c i f i c case of knot o p t i m i z a t i o n would be  useful.  an  - 71  -  A l s o the p o s s i b i l i t y of k e e p i n g knots s t a t i o n a r y w h i l e o t h e r knots c o u l d be u s e f u l i n the i n t e r a c t i v e  Interactiveetechniques tation.  For too l o n g , n u m e r i c a l  box* a l g o r i t h m s  case.  are a v i t a l a s s e t to n u m e r i c a l a n a l y s t s have been d e v e l o p i n g  i n which numbers are read i n and  answers are p r i n t e d out.  By  optimizing  i n t e r a c t i n g with  compu'black  (hopefully) correct  the procedure, what i s  happening to the s o l u t i o n can be observed immediately. g i v e s the power to i n t e r r u p t the s o l u t i o n p r o c e s s  and  Interaction change the  s t a r t i n g c o n d i t i o n s of the problem.  P a r t of the f u t u r e of n u m e r i c a l development of i n t e r a c t i v e methods.  H o p e f u l l y , the use of an  a c t i v e approach r a t h e r than a 'black box' s o l v i n g many n u m e r i c a l  problems.  computation depends on  the  inter-  approach w i l l be v a l u a b l e i n  BIBLIOGRAPHY  Acton, F. S.,  Numerical Methods t h a t  Harper and Box,  M.  J.,  Row, "A New  pp.  Method of C o n s t r a i n e d Other Methods".  and  Golub, G.,  8.  Squares S o l u t i o n s  Numerische Mathematik 7,  and L a u r e n t , P.  the P r a c t i c a l Use  J.,  "On  the Numerical  of I n t e r p o l a t i n g S p l i n e  P r o c e e d i n g s of the IFIP Congress, 1968, de Boor C.  and  imation  TR 20,  imation  "Least  by  1965,  TR  A6-A9.  Squares C u b i c S p l i n e ApproxT e c h n i c a l Report  "Least  Squares C u b i c S p l i n e Approx-  Purdue U n i v e r s i t y :  Technical  21, 1968..  "Numerical Methods f o r S o l v i n g L i n e a r L e a s t  Problems".  Numerische Mathematik .7,, 1965,  T. N.  E.,  T. N.  York:  Lyche, T. and  E.  (ed.),  206-216.  1969,  pp.  Schumaker, L.,  Theory  1-35.  Theory and A p p l i c a t i o n o f S p l i n e  Academic P r e s s ,  polating Natural  pp.  Squares  " I n t r o d u c t i o n to S p l i n e F u n c t i o n s " .  . and A p p l i c a t i o n of S p l i n e F u n c t i o n s , Greville,  Functions".  Purdue U n i v e r s i t y :  I I - V a r i a b l e Knots".  Golub, G.,  pp.  Construction  1968.  and R i c e , J . R.,  Report CSD  Greville,  R i c e , J . R.,  I - F i x e d Knots".  de Boor, C.  New  a  269-276.  Carasso, C.  CSD  and  Computer J o u r n a l  "Linear Least  Householder T r a n s f o r m a t i o n s " .  and  The  Optimization  42-52.  B u s i n g e r , P.  pp.  York:  1970.  Comparison w i t h 1965,  ( U s u a l l y ) Work, New  Functions,  1969. "Computation of Smoothing and  S p l i n e s v i a L o c a l Bases".  Center f o r Numerical A n a l y s i s ,  1971.  Inter-  U n i v e r s i t y o f Texas  -  [11]  O ' R e i l l y , D.,  73 ~  ' " I n t r o d u c t i o n to N o n - l i n e a r  Function  U n i v e r s i t y of B r i t i s h Columbia Computing Centre, [12]  O ^ R e i l l y , D.,  1  "Non-linear  Method of M. Centre, [13]  Patent,  J . Box".  Function  Optimization  P. D.,  R i c e , J...R. ,  "Least  Addison-Wesley, and  Schoenberg, I . J . ( e d l ) ,  S c h u l t z , M., Problems". Functions,  [18]  S c h u l t z , M.,  [19]  Smith, L., Problems". 1969.  t  Ph.D.  Dissertion,  New  "The  Comm. A.CM.  Complex Method f o r 8,  Approximations w i t h pp.  Use  pp.  487-489.  S p e c i a l Emphasis  and  1969. Elliptic  S p e c i a l Emphasis on  Spline  279-347.  S p l i n e A n a l y s i s , New The  1973,  Academic P r e s s ,  " M u l t i v a r i a t e Spline Functions  1969,  - V o l . I I , Reading,  Approximations w i t h York:  1972.  1969.  K u e s t e r , J . L.,  Optimization".  on S p l i n e F u n c t i o n s , [17]  Complex  Square P o l y n o m i a l S p l i n e A p p r o x i m a t i o n " .  The Approximation o f F u n c t i o n s  R i c h a r d s o n , J . A. Constrained  [16]  - The  1971.  Massachusetts; [15]  1971.  U n i v e r s i t y o f B r i t i s h Columbia Computing  C a l i f o r n i a I n s t i t u t e of Technology: [14]  Optimization".  York:  Prentice-Hall,  1973.  o f Man-Machine I n t e r a c t i o n i n D a t a - F i t t i n g  Stanford U n i v e r s i t y :  T e c h n i c a l Report CS  131,  - 74 -  Appendix A  .  The  PROGRAM LISTINGS  f o l l o w i n g pages c o n t a i n l i s t i n g s o f t h e m a j o r i t y o f t h e  program used i n the i n t e r a c t i v e s p l i n e approximation.  The s u b r o u t i n e s  c o n s i s t o f a main c o n t r o l r o u t i n e ..which c a l l s t h e v a r i o u s i n the proper  As  subroutines  order.  i s the. problem w i t h most programs w r i t t e n f o r s p e c i f i c  hardware c o n f i g u r a t i o n s most o f the i n p u t / o u t p u t machine dependent.  T h i s program i s no e x c e p t i o n .  should n o t be d i f f i c u l t  subroutines are However the s u b r o u t i n e s  t o change to a system c o n s i s t i n g o f  correspond-  i n g hardware; t h a t i s , a c o n v e r s a t i o n a l t e r m i n a l , a g r a p h i c s and  a plotter.  1.  terminal  The r o u t i n e s i n q u e s t i o n a r e :  INFREE - which i s a f r e e format i n p u t r o u t i n e .  A l l calls  to INFREE can be r e p l a c e d by FORTRAN READ and FORMAT statements. 2.  The p l o t s u b r o u t i n e s  PLOT, SCALE, AXIS, SYMBOL and PLOTND  need o n l y minor m o d i f i c a t i o n s from one p l o t system to another. subroutine  Note, however, t h a t PSEND i s a s p e c i f i c f o r the graphics  currently available plot 3.  t e r m i n a l which s i m p l y p l o t s a l l  information.  AGPLOT c r e a t e s a hardcopy o f a p l o t c u r r e n t l y b e i n g p l a y e d on a g r a p h i c s  terminal.  I f such a r o u t i n e i s  u n a v a i l a b l e t h i s r o u t i n e can be o m i t t e d .  However t h i s  causes the l o s s o f t h e a b i l i t y to keep the p l o t permanently.  dis-  information  -  75 "  There i s a l s o one m a t r i x m a n i p u l a t i o n s u b r o u t i n e which, i s not g i v e n i n the program.  The s u b r o u t i n e INVERT i s used t o o b t a i n t h e  c o n d i t i o n number o f the m a t r i x . a b l e almost  However a s i m i l a r s u b r o u t i n e i s a v a i l -  everywhere.  The  code f o r the f u n c t i o n m i n i m i z a t i o n i s a l s o n o t g i v e n as  i t was not w r i t t e n by the author.  However a s i m i l a r s u b r o u t i n e i s  a v a i l a b l e i n R i c h a r d s o n and Kuester [ 1 5 ] .  c* C*  C O N T R O L  R O U T I N E  C *  C* DIMENSION X ( 2 0 0 ) , Y ( 2 0 0 ) , YF (200) , RES (200) DIMENSION A ( 2 0 0 ) , DELTA (50) INTEGER YES /'YES'/, NO /'NO'/, ANSWER /» '/ C  C *** ENTER AND PLOT THE DATA POINTS 207 CALL DATAIN (X, Y, N) 205 CALL DATAPL (X, Y, XMIN, XMN, XMAX, EX, YMIN, EY, N) C  C * * * ENTER AND PLOT THE ORIGINAL KNCT SET 206 ERRORP = 1000. CALL KNOTIN (X, DELTA, M, N, MK, S206) CALL KNOTPL (DELTA, XMIN, DX, YMIN, DY, M, HK) C C * * * COMPUTE THE I N I T I A L SPLINE APPROXIMATION CALL FIXED (X, Y, YF, RES, A, DELTA, ERROR, ERRCBF, M, N, MK) C c  ***  C C ***  PLOT THE SPLINE APPROXIMATION CALL SPLNPL (A, DELTA, XMIN, XMN,  XMAX, DX, YMIN, DY, M,  ASK I F HARDCOPY OUTPUT I S WANTED CALL SPLOUT (X, Y, YF, RES, A, DELTA, ERROR, N, M,  MK)  C *** ASK WHETHER OR NOT TO CONTINUE ITERATING 202 WRITE (6,600) 600 FORMAT ('SDO YOU WISH TO VARY THE KNOTS') CALL INFREE (2058, ANSWER, 4) IF (ANSWER.EQ.YES) GO TO 200 I F (ANSWER.EQ.NO) GO TO 201 WRITE (6,601) 601 FORMAT (* SPLEASE ANSWER YES OR NO') GO TO 202 C  C *** VARY THE KNOTS 200 CALL CONTRL (X, Y, YF, RES, A, DELTA, ERROR, N, M, GO TO 202 C C *** ASK WHETHER OR NOT TO OPTIMIZE THE KNOT SET 201 WRITE (6,602) 602 FORMAT ('SDO YOU WISH TO OPTIMIZE THE KNOT SET') CALL INFREE (2058, ANSWER, 4) IF (ANSWER.EQ.YES) GO TO 203 I F (ANSWER.EQ.NO) GO TO 204 WRITE (6,601) GO TO 201  MK)  C  C * * * OPTIMIZE THE KNOT SET 203 CALL OPT (X, Y, YF, RES, A, DELTA, ERROR, N, H , C  HK)  MK)  - 77 -  C * * * ASK WHETHER OR NOT TO RESTART WITH A NEW KNOT SET 204 WRITE (6,603) 603 FORMAT ('&DO YOU WISH TO RESTART WITH A NEW KNOT SET') CALL INFREE (2058, ANSWER, 4) IF (ANSWER.EQ.YES) GO TO 205 IF (ANSWER.EQ.NO) GO TO 209 WRITE (6,601) GO TO 204 C C * * * ASK WHETHER OF NOT TO RESTART 209 WRITE (6,604) 604 FORMAT (»6D0 YOU WISH TO RESTART WITH A NEW DATA SET') CALL INFREE (2058, ANSWER, 4) IF (ANSWER.EQ.YES) GO TO 207 I F (ANSWER.EQ.NO) GO TO 208 WRITE (6,601) GO TO 209 208 CALL PLOTND STOP END  SUBROUTINE DATAIN (X, I, H ) C £***«**************************** ************************************** C* C* I N P U T OF D A T A P O I N T S C* C* C********************************************************************** DIMENSION X (1) , Y (1) C C * * * READ IN THE NUMBER OF DATA POINTS CALL INFREE (4, N) C C *** READ IN THE ABSCISSA AND ORDINATE VALUES DO 100 L=1,N CALL INFREE (20, X (L) , Y (L) ) 100 CONTINUE C C / / / DEBUGGING INFORMATION WRITE (3,301) N 301 FORMAT {•OTHE ', 1 3 , * DATA POINTS ARE: ) DO 400 L=1,N WRITE (3,300) X ( L ) , Y (L) 400 CONTINUE 300 FORMAT (' 1P2E15.6) RETURN END 1  - 79 -  SUBROUTINE KNOTIN  (X, DELTA, M, N, HK, *)  C c *  C* C*  I N P U T  DIMENSION C C *** 600 C C *** 601 C C ***  OF  K N O T  S E T  X ( 1 ) , DELTA (1)  REQUEST THE DEGREE OF THE SPLINE WRITE (6,600) FORMAT ('8DEGREE OF THE S P L I N E » ) CALL INFREE (10, M) READ IN THE NUMBER OF KNOTS WRITE (6,601) FORMAT ('SNUMBER OF KNOTS') CALL INFREE (10, K) CALCULATE THE TOTAL NUMBER OF MK = M + K + 1  KNOTS  C C *** 603  CHECK FOR AN UNDERDETERMINED SYSTEM IF (MK.LE.N) GO TO 200 WRITE (6,603) FORMAT ('OSYSTEM I S UNDERDETERMINED:'/ * ' REDUCE DEGREE OF SPLINE OR NUMBER OF KNOTS') RETURN 1  C C *** READ IN THE POSITION CF THE KNOTS 200 DO 100 1=1,K IM = I • H + 1 WRITE (6, 602) I CALL INFREE (26, DELTA (IM)) 100 CONTINUE 602 FORMAT ('8P0SITI0N OF KNOT (', 1 2 , » ) • ) C C *** CALCULATE THE POSITION OF THE SUPPLEMENTARY H = (X(N) - X ( 1 ) ) / (K + 1) M1 = H + 1 DO 101 1=1,M1 DELTA (I) = X(1) + (I - M1) * H DELTA (MK + I) = X (N) + ( 1 - 1 ) * H 101 CONTINUE RETURN END  KNOTS  -80-  SUBROUTINE  DATAPL  (X, Y, XMIN,  XMN , XMAX, EX, YMIN,  DY,  N)  C  Q*************************************************  c*  C* P L O T OF D A T A P O I N T S C* ' • C* Q**************************** ************************************* DIMENSION X (1) , Y (1) DIMENSION XS (200) , YS (200) C C *** RESET THE PLOTTER CALL PLOT ( 1 5 . , 0., -3) C C *** PROTECT THE ORIGINAL DATA POINTS XMN = X (1) XMAX = X(N) DO 100 L=1,N XS(L) = X(L) YS(L) = Y (L) 100 CONTINUE C C *.** SCALE THE ABSCISSA AND ORDINATE VALUES CALL SCALE (XS, N, 10., XMIN, DX, 1) CALL SCALE (YS, N, 10., YMIN, EY, 1) C *** PLOT THE DATA POINTS DO 101 L=1,N CALL SYMBOL ( X S ( L ) - . 0 3 5 , Y S ( I ) - . 0 7 , .14, • + •, 0., 1) 101 CONTINUE * * * PLOT THE AXES ORY = -YMIN / DY IF (ORY.LT.O.) ORY = 0. CALL AXIS ( 0 . , ORY, • - 1 , 10., 0., XHIN, DX) ORX = -XMIN / DX I F (ORX.LT.O.) ORX = 0. CALL AXIS (ORX, 0., » 1, 10., 9 0 . YMIN, DY) c  C  f  •  C  C  ***  DISPLAY THE PLOT CALL PSEND RETURN END  - 81 -  SUBROUTINE SPLNPL #  (ft, DELTA, XMIN, XMN, DY, M, MK)  XMAX, DX,  YMIN,  C c *  C* C* C*  P L O T  OF  S P L I N E  C U R V E  DIMENSION A ( 1 ) , DELTA (1) COMMON /DW/ DOMEGA(50,5) C C *** C C ***  C C ***  CALCULATE THE SPLINE FUNCTION CALL OMEGA (DELTA, M, MK) RAISE THE PEN TO REPOSITION IPEN =3 X = XMN XSTEP = .02 * DX NP = 50 * INT ( (XMAX - XMN) DO 100 J=1,NP  DENOMINATOR  IT  / DX  + .5)  CALCULATE THE SPLINE FUNCTION VALUE AT THE CURRENT X VALUE Y = 0. DO 101 1=1,MK Y = Y + A (I) * SPLINE (DELTA, X, M, MK, I , 0) CONTINUE  101 C C *** SCALE THE X AND Y VALUES XP = (X - XMIN) / DX YP = (Y - YMIN) / DY C *** I F THE POINT IS WITHIN THE RANGE I F ( (YP.GT. 10.5) .OR. (YP.LT.-.5) ) GO TO C *** THEN PLOT I T CALL PLOT (XP, YP, IPEN) IPEN =2 GO TO 201 C C *** OTHERWISE RAISE THE PEN 200 IPEN = 3 201 X = XMN + J * XSTEP 100 CONTINUE C C *** DISPLAY THE SPLINE CURVE CALL PSEND RETURN END c  C  200  - 82SUBROUTINE KNOTPL  (DELTA, XMIN, DX, YMIN, DY , M, MK)  C C* C* C* C*  P L O T  O F  K N O T  S E T  Q***************************************** * * * * * * * * * * * * * * * * * * * * * * * * * * * * *  DIMENSION  DELTA(1)  C  C ***  C C *** C C *** 100 C C ***  SCALE THE KNOT ORDINATE ONTO THE X-AXIS YP = -YMIN / DY IF (YP.LT.0.0) YP = 0.0 M2 = M + 2 DO 100 I=M2,HK SCALE  THE KNOT XP = (DELTA (I) - XMIN) / DX  PLOT THE KNOT CALL SYMBOL CONTINUE DISPLAY THE KNOTS CALL PSEND RETURN END  (XP-.0525, YP-.115,  . 2 1 , «X», 0., 1)  - 83 -  SUBROUTINE FIXED (X, Y, YF, RES, A, DELTA, ERROR, ERRORP, M, N, MK C C********************************************************************** c *  C* C* C*  F I X E D K N O T L E A S T S Q U A R E S  L I N E A R A P P R O X I M A T I O N  c *  Q**********************************************************************  DIMENSION DIMENSION  X(1), Y(1), S (200,50)  YF(1),  RES{1), A ( 1 ) ,  E ELT A (1)  C C * * * CALCULATE THE SPLINE FUNCTION DENOMINATOR CALL OMEGA (DELTA, M, MK) C C * * * FORM THE MATRIX OF SPLINE BASIS FUNCTION VALUES DO 100 1=1,MK  DO 101 L=1,N S ( L , I ) = SPLINE CONTINUE CONTINUE  (DELTA, X (L) , M, MK, I , 101 100 C C / / / DEBUGGING INFORMATION WRITE (3,300) 300 FORMAT (»0AT », 13X, "THE HOUSEHOLDER MATRIX I S : * ) DO 400 L=1,N WRITE (3,301) X (L) WRITE (3,302) ( S ( L , I ) , 1=1,MK) 400 CONTINUE 301 FORMAT (« », 1P8E15.6) 302 FORMAT (» «, 15X, 1P7E15.6) C C * * * TRANSFER Y INTO A DO 102 L=1,N A (L) = Y ( L ) 102 CONTINUE C C / / / DEBUGGING INFORMATION WRITE (3,303) 303 FORMAT (•OTHE HOUSEHOLDER VECTOR I S : » ) WRITE (3,301) (A (L) , L=1,N) C C * * * PERFORM THE HOUSEHOLDER TRANSFORMATIONS CALL TRANS (S, A, N, MK) C C * * * CALCULATE THE RESIDUALS AND LEAST SQUARES ERROR ERROR = RESERR (X, Y, Y F , RES, A, DELTA, N, M, MK) CALL ITER (ERRORP, ERROR) RETURN END  0)  - 84 -  FUNCTION RESERR  (X, Y, YF, RES, A, DELTA, N, M, HK)  C  C*******************************  C* C* C* C* C*  **************** * ^ « ^ « 4 0 H 4 4 * * 4 < > » * * 4 * * *  R E S I D U A L S L E A S T S Q U A R E S  A N E E R R O R  c********************************************************************** DIMENSION X ( 1 ) , Y ( 1 ) , Y F ( 1 ) , R E S ( 1 ) , DOUBLE PRECISION DY, DRES, DERR COMMON /DW/ DOMEGA(50,5) C C ***  C C ***  101 C C ***  C C *** 100 C C /// 300  CALCULATE THE SPLINE FUNCTION CALL OMEGA (DELTA, M, MK) DERR = 0.D0 DO 100 L=1,N  A ( 1 ) , DELTA (1)  DENOMINATOR  CALCULATE THE SPLINE FUNCTION VALUES AT T E E DATA POINTS DY = 0.D0 DO 101 1=1,MK DY = DY + A (I) * SPLINE (DELTA, X (L) , M, MK, I , 0) CONTINUE YF (L) = DY CALCULATE THE RESIDUALS DRES = Y(L) - DY RES (L) = DRES CALCULATE THE LEAST SQUARES ERROR DERR = DERR + DRES * DRES CONTINUE RESERR = DSQRT (DERR) DEBUGGING INFORMATION WRITE (3,300) RESEBR FORMAT ('OTHE LEAST SQUARES ERROR RETURN END  IS •, 1 P E 1 5 . 6 )  - 85 -  SUBROUTINE CONTRL  (X, Y, YF, RES, A, DELTA, ERRORP, N, M, MK)  C c * c *  C* C* C*  V A R I A B L E  K N O T  C O N T R O L  DIMENSION X ( 1 ) , Y ( 1 ) , Y F ( 1 ) , R E S ( 1 ) , A ( 1 ) , D ELT A (1) DIMENSION AP{200), DELTAP (50) , RESP(200), Y F P ( 2 0 0 ) INTEGER YES /'YES'/* NC /'NO'// ANSWER / • •/ MMK = MK • M + 1 C C ***  TRANSFER THE PREVIOUS KNOT SET DO 100 1=1,MMK DELTAP (I) = DELTA (I) CONTINUE FORMAT ('SPLEASE ANSWER YES OR NO')  100 601 C C *** OBTAIN WHICH KNOT TO VARY 200 CALL VARYIN (I) C C *** REPOSITION THAT KNOT IM = I + M + 1 CALL VARYKT (DELTAP (IM), I) C C *** ASK WHETHER OR NOT TO REPOSITION OTHER KNOTS 204 WRITE (6,602) 602 FORMAT ('SDO YOU WISH TO VARY ANOTHER KNOT') CALL INFREE (2058, ANSWER, 4) IF (ANSWER.EQ.YES) GO TO 200 I F (ANSWER.EQ.NO) GO TO 203 WRITE (6,601) GO TO 204 203 CONTINUE C C / / / DEBUGGING INFORMATION WRITE (3,300) 300 FORMAT ('OTHE NEW KNOT SET I S : ' ) WRITE (3,301) (DELTAP ( I ) , 1=1,MMK) 301 FORMAT (• «, 1P8E15.6) C C *** RECALCULATE THE SPLINE APPROXIMATION CALL FIXED (X, Y, YFP, RESP, AP, DELTAP, ERROR, ERRORP, M, N,*MK) C C * * * ASK WHETHER OR NOT TO PLOT THE NEW APPROXIMATION 207 WRITE (6,603) 603 FORMAT ('SDO YOU WISH TO SEE THE NEW APPROXIMATION') CALL INFREE (2058, -ANSWER, 4) I F (ANSWER. EQ. YES) .GO TO 205 IF (ANSWER.EQ.NO) GO TO 201 WRITE (6,601) GO TO 207  - 86 205  C  CALL CALL  DATAPL KNOTPL  (X, Y, XHIN, XHN, XMAX, EX, YMIN, EY, (DELTAP, XMIN, DX, YMIN, DY, M , MK)  CALL  SPLNPL  (AP,  DELTAP, OR  NOT  TO  XMIN, KEEP  XMN, XMAX,  C *** 206  DECIDE WHETHER WRITE (6,604)  NEW  604  FORMAT («SDO YOU WISH TO C O N T I N U E CALL INFREE (2058, ANSWER, 4) IF (ANSWER.EQ.YES) GO TO 2 0 8 IF (ANSWER.EQ.NO) GO T O 2 0 9 WRITE (6,601) GO T O 2 0 6  DX,  N)  YMIN,  EY,  M,  MK)  APPROXIMATION WITH  THE  NEW  APPROXIMATION•)  C C *** 208 101  IF SO, RETAIN DO 1 0 1 1=1,MK A (I) CONTINUE DO  102  103 C C  ***  102  =  T H E NEW  VALUES  AP(I)  1=1,MMK  DELTA (I) CONTINUE DO 1 0 3 L=1,N  =  DELTAP  (I)  YF (L) = YFP(L) RES (L) = RESP(L) CONTINUE ASK  IF  A  HARDCOPY  CALL SPLOUT (X, ERRORP = ERROR GO  TO  Y,  IS YF,  WANTED RES,  A,  DELTA,  ERROR,  N ,  M , MK)  201  ,  •  C C  ***  209  201  REPLOT CALL CALL  THE OLD  DATAPL KNOTPL  CALL SPLNPL RETURN END  SPLINE  (X, Y, XMIN, XMN, XMAX, (DELTA, XMIN, DX, YMIN, (A,  DELTA,  XMIN,  DX, EY,  XMN, XMAX,  YMIN, DY, M , MK)  N)  DX,  DY,  YMIN,  M,  MK)  - 87 SUBROUTINE  C* C* C* C*  600  VARYIN  I N P U T  O F  WRITE (6,600) FORMAT ( ' S K N O T TO BE CALL INFREE (10, I) RETURN END  SUBROUTINE  C* C* C* C*  600  (I)  V  A  VARYKT  R  I  A  WRITE (6,600) FORMAT ('SNEW CALL  IN F R E E  RETURN END  B  V  I  A  E  I  E  K N O T  1  (DELTAI,  L  R  VARIED )  E  I POSITION  (26,  A  I)  K N O T  FOR  DELTAI)  KNOT  V A L U E  (',  12,  F  )  ')  - 88 SUBROUTINE ITER  (ERRORP, ERROR)  C  C*******************«*******4$****^ c *  C* C* C* 600 601  I T E R A T I O N  O U T P U T  WRITE (6,600) ERRORP FORMAT (•OTHE LEAST SQUARES ERROR OF THE PREVIOUS «, * 'APPROXIMATION WAS 1PE15.6) WRITE (6,601) ERROR FORMAT (• THE LEAST SQUARES ERROR OF THE NEW «, * 'APPROXIMATION IS •, 1PE15.6) RETURN END  - 89 -  SUBROUTINE SPLOUT  (X, Y, YF, RES, A, DELTA, ERROR, N, M, HK)  C c *  C* * C*  O U T P U T  C  DIMENSION X{1), Y ( 1 ) , Y F ( 1 ) R E S ( 1 ) , A ( 1 ) , DELTA (1) INTEGER YES /'YES'/* NO /'NO'/, ANSWER / • •/ INTEGER T I T L E ( 1 0 ) , BLANK / ' •/, WUM / 1 / f  C C *** ASK WHETHER OR NOT TO RETAIN RESULTS 200 WRITE (6,600) 600 FORMAT ('SDO YOU WISH TO RETAIN THE APPROXIMATION RESULTS*) CALL INFREE (2058, ANSWER, 4) IF .(ANSWER. EQ. YES) GO TO 201 I F (ANSWER.EQ.NO) GO TO 202 WRITE (6,601) 601 FORMAT ('&PLEASE ANSWER YES OR NO') GO TO 200 C C *** BLANK OUT THE T I T L E 201 DO 101 1=1,10 T I T L E (I) = BLANK 101 CONTINUE C C * * * ASK FOR THE T I T L E OF THE PLOT WRITE (6,602) 602 FORMAT (•STITLE FOR PLOT') CALL INFREE (2058, T I T L E , 40) C C * * * HARDCOPY PRINTOUT WRITE (8, 603) T I T L E , NUM 603 FORMAT ('1', 10A4, 20X, 15) WRITE (8, 604) 604 FORMAT ('0', 5X, 'ABSCISSAE', 5X, 'ORDINATES', 3X, • 'FITTED ORDINATES', 3X, 'RESIDUALS') C C *** PRINT OUT THE DATA, F I T T E D RESULTS AND RESIDUALS DO 100 L=1,N WRITE (8, 605) X (L) , Y (L) , YF (L) , RES (L) 100 CONTINUE 605 FORMAT (» ', 1P5E15.6) C C * * * PRINT OUT THE LEAST SQUARES ERROR WRITE (8,606) ERROR 606 FORMAT ('OTHE LEAST SQUARES ERROR I S ' , 1PE15.6) C C *** CALCULATE THE POLYNOMIAL COEFFICIENTS CALL COEFF (X, Y, A-, DELTA, N, M, M K )  C C ***  PLOT A HARDCOPY OF THE APPROXIMATION CALL SYMBOL ( 1 . , 10., .14, T I T L E , 0., 40)  - 90 -  C  FNUM  =  CALL  NUMBER  CALL  PSEND  CALL NUM  202  NUM  AGPLOT =  RETURN END  NUM  (9.,  10.,  (10.0, + 1  .14,  6000)  FNUM,  0 . ,  -1)  - 91 -  SUBROUTINE NORM  (S, MK)  C C* C* M A T R I X C O N D I T I O N N U M B E R C* C* C«**************************«*^ DIMENSION S ( 2 0 0 , 5 0 ) , ST(50,50) C C *** SAVE THE ORIGINAL MATRIX DO 100 1=1,MK DO 100 L=1,MK ST(L,I) = S(L,I) 100 CONTINUE C C *** INVERT THE MATRIX AND CALCULATE ITS CONDITION NUKEER CALL INVERT (ST, MK, 50, DET, COND) C C / / / DEBUGGING INFORMATION WRITE (3,300) 300 FORMAT ('OTHE INVERTED TRANSFORMED HOUSEHOLDER MATRIX I S : ' ) DO 400 1=1,MK WRITE (3,301) ( S T ( L , I ) , L=1,MK) WRITE (3,302) 400 CONTINUE 301 FORMAT (• ', 1P8E15.6) 302 FORMAT (' ») C C *** PRINT OUT THE CONDITION NUMBER WRITE (6,600) COND 600 FORMAT {•OTHE CONDITION NUMBER I S ' , 1PE15.6) RETURN END  r- 92 -  SUBROUTINE TRANS  (S, B  f  N, MK)  C  c* C* H O U S E H O L D E R T R A N S F O R M A T I O N C* C O N T R O L R O U T I N E C* C* C************************************************************4********* DIMENSION S ( 2 0 0 , 5 0 ) , V (200) DIMENSION B(1) REAL KSQ C C * * * AUGMENT THE MATRIX S BY B MK1 = MK + 1 DO 100 L=1,N S (L,MK 1) = B (L) 100 CONTINUE C * * * PERFORM THE HOUSEHOLDER TRANSFORMS DO 101 1=1,MK C C * * * TRANSFER CURRENT COLUMN OF HOUSEHOLDER MATRIX DO 102 L=1,N V(L) = S ( L , I ) 102 CONTINUE C C * * * DERIVE HOUSEHOLDER TRANSFORMATION VECTOR CALL UVEC (V, N, I , KSQ) C C * * * PERFORM HOUSEHOLDER TRANSFORMATION ON THE VECTOR CALL HOUSE (S, V, N, MK1, I , KSQ) 101 CONTINUE C C / / / DEBUGGING INFORMATION WRITE (3,300) 300 FORMAT ('OTHE TRANSFORMED HOUSEHOLDER MATRIX I S : ' ) DO 400 L=1,N WRITE (3,301) ( S ( L , I ) , 1=1,MK) WRITE (3,302) 400 CONTINUE 301 FORMAT (« • , 1P8E15.6) 302 FORMAT (» ') C C * * * CALCULATE THE CONDITION NUMBER OF THE MATRIX CALL NORM ( S , MK) C C * * * PERFORM THE BACKWARDS SUBSTITUTION CALL SOLVE (S, B, N, MK) RETURN END C  SUBROUTINE UVEC (V, N I , KSQ) C Q************************************************* #  C* C* H O U S E H O L D E R V E C T O R C* C* Q************************ ******************************************** DIMENSION V(1) DOUBLE PRECISION DNORH REAL KSQ C C *** COMPUTE THE NORM OF V DNORM = 0.D0 DO 100 L=I,N DNORM = DNORM + V (L) * V (L) 100 CONTINUE DNORM = DSQRT (DNORM) C C * * * CALCULATE THE COEFFICIENT NSGN = 1 I F ( V ( I ) . L T . O . ) NSGN = -1 KSQ = 2 * DNORM * (DNORM • NSGN * V ( I ) ) C C *** CREATE THE HOUSEHOLDER VECTOR V ( I ) = V ( I ) + NSGN * DNORM RETURN END  - 94 -  SUBROUTINE HOUSE  (S, U, N, HK1, I , KSQ)  C C* C* C*  H O U S E H O L D E R  T R A N S F C R M A T  c*  Q*************************************** ******************  DIMENSION REAL KSQ  I C N  *******  S ( 2 0 0 , S 0 ) , U(1)  C ***  C  PERFORM HOUSEHOLDER DO 100 J=I MK1  TRANSFORMATION  ON REQUIRED COLUMNS OF MATRIX  f  C C *** C  101 C C ***  102 100  i  MULTIPLY THE TRANSPOSE OF THE VECTOR BY THE COLUMN OF THE MATRIX UT = 0. DO 101 L=I,N UT = UT + U (L) * S (L, J) CONTINUE UT = 2. * UT / KSQ PERFORM THE HOUSEHOLDER TRANSFORMATION DO 102 1=1,H S (L, J) = S (L,J) - UT * U (L) CONTINUE CONTINUE RETURN END  - 95 -  SUBROUTINE SOLVE (S, B, N, HK) C C***************4*******4***4**4**** C* C* B A C K W A R D S S U B S T I T U T I O N C* C* , C C ***  101  DIMENSION S (200, 1) , B (1) MK1 = MK + 1 COMPUTE THE FINAL ELEMENT B (MK) = S(MK,MK1) / S (MK , MK) I F (MK.EQ.1) GO TO 200 MKM1 = MK - 1 DO 100 I=1,MKM1 MKI = MK - I B(MKI) = 0. MKI1 = MKI • 1 DO 101 J=MKI1,MK B (MKI) = B (MKI) + B (J) * S(MKI,J) CONTINUE B (MKI) = (S(MKI,MK1) - B (MKI) ) / S(MKI,MKI) CONTINUE  100 C C * * * TRANSFER THE REMAINING 200 I F (N.EQ.MK) GO TO 201 DO 102 L= MK1 , N B(L) = S(L,MK1) 102 CONTINUE 201 RETURN END  ELEMENTS INTO B  - 96 SUBROUTINE OPT (X, Y, YF, RES, A, DELTA, EBRCRP, N, M, MK) C C**********************************************************************  c* C* C*  O P T I M I Z A T I O N  OF  K N O T S  c* Q********************************************************************** DIMENSION X {1) , Y ( 1 ) , Y F ( 1 ) , R E S ( 1 ) , A ( 1 ) , BELT A (1) DIMENSION A P ( 2 0 0 ) , DELTAP (50) , Y F P ( 2 0 0 ) , RESP (200) DIMENSION P (50,75) INTEGER YES / • Y E S » / » NO / ' K O V , ANSWER / • »/ EXTERNAL RESERR, GDELTA, HDELTA C C ***  TRANSFER KNOTS FOR OPTIMIZATION MMK = MK + M + 1 DO 300 1=1,MMK DELTAP (I) = DELTA (I) CONTINUE  300 C C ***  TRANSFER PARAMETERS DO 301 1=1,MK A P ( I ) = A (I) CONTINUE  301 C C ***  TRANSFER THE KNOTS TO THE SIMPLEX DO 102 J=1,MMK P ( J , 1) a DELTAP (J) CONTINUE  102 C C *** C C ***  FOR OPTIMIZATION  COMPUTE THE NUMBER OF POINTS IN THE SIMPLEX I F (MMK.GE.10) NPLUS = (3 * MMK + 1) / 2 I F (MMK.LT.10) NPLUS = 2 * MMK MINIMIZE THE LEAST SQUARES ERROR CALL COMPLX (X, Y, YFP, RESP, AP, DELTAP, EBRORP, # P, 50, MMK, NPLUS, MMK, M, N, MK, 1.3, - 1 , # 100, 100, 1, - . 0 0 0 1 , RESERR, GDELTA, # HDELTA, G.DELTA, 5200, S200)  C C *** C C *** 200 801 C C *** 207 603  RECALCULATE THE FIXED KNOT APPROXIMATION CALL FIXED (X, Y, YFP, RESP, AP, DELTAP, ERROR, EBRORP, GO TO 207 I F KNOT OPTIMIZATION FAILED, PRINT AN ERROR WRITE (6,801) FORMAT (» 0OPTIMIZATION FAILED. *)  MESSAGE  ASK WHETHER OR NOT TO PLOT THE NEW APPROXIMATION WRITE (6,603)' FORMAT (•SDO YOU WISH TO SEE THE NEW APPROXIMATION') CALL INFREE (2058, ANSWER, U) I F (ANSWER.EQ.YES) GO TO 205 IF (ANSWER.EQ.NO) GO TO 201  M, N,  MK)  - 97 -  601 205  WRITE (6,601) FORMAT (•SPLEASE ANSWER YES OR NO') GO TO 207 CALL DATAPL (X, Y, XMIN, XMN, XMAX, DX, YMIN, DY, N) CALL KNOTPL (DELTAP, XMIN, DX, YMIN, DY, M, MK) CALL SPLNPL (AP, DELTAP, XMIN, XMN, XMAX, DX, YMIN, DY, M,  MK) C C *** DECIDE WHETHER OR NOT TO KEEP NEW APPROXIMATION 206 WRITE (6,604) 604 FORMAT ('SDO YOD WISH TO CONTINUE WITH THE NEW APPROXIMATION') CALL INFREE (2058, ANSWER, 4) I F (ANSWER.EQ.YES) GO TO 208 IF (ANSWER.EQ.NO) GO TO 209 WRITE (6,601) GO TO 206 C *** I F SO, RETAIN THE NEW VALUES 208 DO 101 1=1,MK A (I) = A P ( I ) 101 CONTINUE DO 302 1=1,MMK DELTA (I) = DELTAP (I) 302 CONTINUE DO 103 L=1,N YF (L) = YFP(L) RES (L) = RESP (L) 103 CONTINUE C • . C *** ASK I F A HARDCOPY I S WANTED CALL SPLOUT (X, Y, YF, RES, A, DELTA, ERROR, N, M, MK) ERRORP = ERROR GO TO 201 C C *** REPLOT THE OLD SPLINE 209 CALL DATAPL (X, Y, XMIN, XMN, XMAX, EX, YMIN, DY, N) CALL KNOTPL (DELTA, XMIN, DX, YMIN, DY, M, MK) CALL SPLNPL (A, DELTA, XMIN, XMN, XMAX, EX, YMIN, DY, M, MK) 201 RETURN C  END  '  - 98 -  FUNCTION GDELTA  C  (T, H, MMK,  J)  C***********************4*4*44*****************************************  C* C* C*  L O W E R  B O U N D  C O N S T R A I N T  c *  £********************************************************************** DIMENSION T ( 1 ) C C *** SET THE LOWER BOUNDS ON THE KNOTS IF ( ( J . L E . M+1) .OR. (J.GT. (MMK-(M + 1) )) ) G O T O 100 H = T (MMK- (M + 1) ) - T (M + 1) GDELTA = T ( J - 1 ) + .0001 * H GO TO 999 100 GDELTA = T (J) 999 RETURN END  FUNCTION HDELTA (T, M, MMK, J) C C********************************************************************** C* C* U P P E R B O U N D C O N S T R A I N T C* C*  c * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *  DIMENSION  T (1)  C c  ***  100 999  SET THE UPPER BOUNDS ON THE KNOTS IF ( (J.LE.M + 1) .OR. (J.GT. (MMK-(M + 1)))) H = T (MMK- (M+1) ) - T (M + 1) HDELTA = T ( J + 1) - .0001 * H GO TO 999 HDELTA = T ( J ) RETURN END  GOTO  100  - 99 -  SUBROUTINE COEFF  C* C* C* C*  (X, Y, A, DELTA, N, M,  P O L Y N O M I A L  MK)  C O E F F I C I E  T S  **** DIMENSION X ( 1 ) , Y ( 1 ) A ( 1 ) , DELTA (1) DIMENSION C(50,4) COMMON /DW/ DOMEGA(50,5) r  C C ***  CALCULATE THE SPLINE DENOMINATOR CALL OMEGA (DELTA, M, MK)  C ***  C  DETERMINE THE NUMBER OF PIECEWISE POLYNOMIALS K = MK - M M1 = M + 1  C c  #**  102 101  DETERMINE THE COEFFICIENTS OF THE PIECEWISE POLYNOMIALS JFACT = 1 DO 100 J=1,M1 J1 = J - 1 DO 101 1=1,K IM = I + M C ( I , J ) = 0.0 DO 102 L=1,MK C ( I , J ) = C ( I , J ) + A (L) * * SPLINE (DELTA, DELTA (IM) , H, MK, L,\ J 1 ) CONTINUE C ( I , J ) = C ( I , J ) / JFACT CONTINUE JFACT = JFACT * J CONTINUE  100 C C ***  800  104 103  801 802  PRINT OUT THE COEFFICIENTS WRITE (8,800) FORMAT («0 KNOT', 5X, 'LOCATION', 5X, 'POLYNOMIAL POWER *, * 5X, 'POLYNOMIAL COEFFICIENT') DO 103 1=1,K 11 = I - 1 IM = I + M WRITE (8,801) 11, DELTA (IM) DO 104 J=1,M1 J1 = J - 1 WRITE (8,802) J 1 , C ( I , J ) CONTINUE CONTINUE MK1 = MK + 1 WRITE (8,801) K, DELTA (MK1) FORMAT (' *, 14, 2X, 1PE15.6) FORMAT (• ', 30X, 11, 15X, 1PE15.6) RETURN END  - 100.FUNCTION  SPLINE  (DELTA,  X ,  H,  C********************************  MK, I ,  J )  :  c * S P L A N D  C* C* C* C*  I  N  E D  E  R  DIMENSION D E L T A (1) COMMON / D W / D O M E G A ( 5 0 , 5 ) DOUBLE PRECISION DSPLN, DX, MCON MJ  ***  N -  C c  ***  100  ***  I F  THE POINT  IS  I N THE REGION  EVALUATE DX  =  X  M2  =  M +  DO  100  THE SPLINE  SUPPORT  AT  THE GIVEN  "  LL = I + L 1 DDELTA = DELTA (LL) DY = D D E L T A - DX IF (DY.LE.O.DO) GO TO  100  DSPLN  MJ /  =  DSPLN  +  DY  1) . O R . ( I . E Q . M K ) THE CONSTANT  (J.EQ.O) =  MJ +  101  ** )  1  L=MJ1,M MCON  SPLINE END  TERM  GO T O 2 0 0  =  -MCON  *  CONTINUE RETURN  POINT  2  COMPUTE  DO  GO T O 2 0 0  L=1,M2  CONTINUE IF ( (I.EQ.  MJ1  200  OF  = M • 1 ( ( X . L T . D E L T A ( I ) ) . O R . ( X . G T . D E L T A ( I + M1) ) )  IF  101  BY  J  C c  DDELTA,  = 0 . D 0  CHECK M1 IF  I O N V E S  = 1  =  DSPLN C c  F U N C T I V A T I  =  MCON  *  DSPLN  L  DSPLN  DOMEGA  ( I , L )  =  DSPLN  M1  *  - 101 -  SUBROUTINE C* C* C* C*  S  OMEGA  S P L I N E  (DELTA, M,  MK)  F U N C T I O N  D E N O M I N A T O R  DIMENSION DELTA (1) DOUBLE PRECISION DDELTA COMMON /DW/ DOMEGA(50,5) C C *** *** C  100 C C /// 300 400 301  CALCULATE THE MATRIX OF VALUES FOR THE DENOMINATOR OF THE SPLINE FUNCTION M2 = M + 2 DO 100 1=1,MK DO 100 J=1,M2 DOMEGA(I,J) = 1. JJ = I + J - 1 DO 100 1=1,M2 I F (L.EQ.J) GO TO 100 LL = I + L - 1 DDELTA = DELTA (JJ) - DELTA (LL) DOMEGA(I,J) = DOMEGA(I,J) * DDELTA CONTINUE DEBUGGING INFORMATION WRITE (3,300) FORMAT (*OTHE BASIS FUNCTION DENOMINATOR VALUES DO 400 1=1,MK WRITE (3,301) (DOMEGA ( I , J ) , J=1,M2) CONTINUE FORMAT (' •, 1P5E15.6) RETURN END  ARE:»)  - 102 -  Appendix B  INTERACTIVE SPLINE  APPROXIMATION  Purpose  T h i s system enables the u s e . t o p e r f o r m l e a s t squares imations u s i n g s p l i n e f u n c t i o n s .  approx-  These approximations a r e performed  i n t e r a c t i v e l y w i t h the a i d o f a g r a p h i c s t e r m i n a l .  The approximations  are d i s p l a y e d immediately a f t e r computation and theuuser can r e s p e c i f y knot l o c a t i o n s and r e c a l c u l a t e the f i t .  Features are included to  o p t i m i z e the knot s e t and change the number o f knots i n t h e knot s e t .  Type o f Routine  T h i s i s a s e l f - c o n t a i n e d program w r i t t e n i n FORTRAN IV.  How to Use  To r u n t h i s program under MTS a t t h e Adage G r a p h i c s T e r m i n a l ; e n t e r t h e command:  $RUN IRAM:SPLINE+AGT:BASIC  3=debugfile 4 = d a t a f i l e 8 = p r i n t f i l e 9 - p l o t f i l e  where debugfile  c o n t a i n s the debugging  information.  U n l e s s the system  runs i n t o problems w i t h t h e a p p r o x i m a t i o n t h i s s h o u l d be s e t t o *DUMMY* . datafile  i s the f i l e c o n t a i n i n g t h e i n p u t d a t a . d e s c r i b e d i n S e c t i o n A:  The format i s  Data Input Format.  - 103 -  printflle  i s the f i l e . c o n t a i n i n g the p r i n t e d output from an approxi m a t i o n c o r r e s p o n d i n g to a p l o t .  plotfile  i s the f i l e . c o n t a i n i n g the p l o t i n f o r m a t i o n to be r e t a i n e d for  a  hardcopy.  Upon c o m p l e t i o n o f t h e program a hardcopy :  o f the p r i n t o u t s  and p l o t s which were saved can be o b t a i n e d by i s s u i n g the f o l l o w i n g commands: $C0PY p r i n t f i l e  *PRINT*  , $R PL0T:Q E A R = p l o t f i l e  where ' p r i n t f i l e '  c o n t a i n s the p r i n t output d e s c r i b e d above and  'plotfile'  c o n t a i n s the p l o t i n f o r m a t i o n .  Description A.  Data Input Format  The d a t a f i l e has the f o l l o w i n g The  first  structure:  l i n e s t a t e s the number o f p o i n t s i n the f i r s t  data set  f o l l o w e d by d a t a l i n e s w i t h the sequence  ABSCISSA  ORDINATE.  The r e m a i n i n g d a t a s e t s f o l l o w w i t h an i d e n t i c a l format. itself  i s i n f r e e - f o r m a t p r o v i d i n g t h a t a t l e a s t one b l a n k  each e n t r y .  The d a t a delimits  - 104  B.  -  Displays  Because of the i n t e r a c t i v e c a p a b i l i t i e s o f the system, p l a y an i n t e g r a l p a r t i n the s t r u c t u r e .  displays  Inoorder t o b e s t d e s c r i b e the  f l o w o f c o n t r o l through the system a sample run f o l l o w s w i t h comments i n s e r t e d to augment the d i s p l a y i n f o r m a t i o n .  Example;  Demonstration  $SIG IRAM 'DEMONSTRATION  run  RUN  1,  PASSWORD Signon i n f o r m a t i o n $C0MMENT LOAD GRAPH IF NECESSARY $C0PY AGT:GRAPH > AGTI $COMMENT SPLINE APPROXIMATION PROGRAM $RUN IRAM:SPLINE+AGT;BASIC 3=*DUMMY* 4=DATA 8=-PRINT 9=-PL0T EXECUTION BEGINS (One d a t a s e t from u n i t 4 i s immediately read i n .  A p l o t of t h i s s e t  appears on the g r a p h i c s d i s p l a y . ) DEGREE OF THE  SPLINE? 3  (The degree of the p i e c e w i s e p o l y n o m i a l s wanted i s r e q u e s t e d , most f r e q u e n t l y used degree i s 3; t h a t i s , a c u b i c s p l i n e  approximation.)  NUMBER OF KNOTS? 4 (The number o f i n t e r n a l knots wanted i s r e q u e s t e d .  Boundary knots  and supplementary k n o t s are computed by the program.)  - 105 -  POSITION OF KNOT ( 1)? .2 POSITION OF KNOT ( 2)1 .4 POSITION OF KNOT ( 3)? .6 POSITION OF KNOT ( 4 ) ? .8 (The l o c a t i o n o f t h e a b s c i s s a o f each knot i s r e q u e s t e d . be e n t e r e d i n a s t r i c t l y  i n c r e a s i n g sequence.  The knots must  Immediately  f o l l o w i n g the  i n p u t of, a l l t h e knots t h e knot s e t i s p l o t t e d a l o n g t h e x - a x i s . ) (A f i x e d knot l e a s t squares approximation i s now computed.  The r e s u l t i n g  s p l i n e f u n c t i o n i s o v e r l a i d on t h e g r a p h i c s d i s p l a y ) . THE LEAST SQUARES ERROR OF THE PREVIOUS APPROXIMATION WAS (The i n i t i a l l e a s t squares e r r o r i s a r b i t r a r i l y  1.OOOOOOE 03  s e t t o 1000 s i n c e no  p r e v i o u s a p p r o x i m a t i o n was done.) THE LEAST SQUARES ERROR OF THE NEW. APPROXIMATION IS DO YOU WISH TO RETAIN THE APPROXIMATION RESULTS?  1.537876E-02  NO  (A message a s k i n g whether o r n o t t o keep a hardcopy  of the present  approximation s e t i s p r i n t e d . ) DO YOU WISH TO VARY THE KNOT SET?  YES  (A message a s k i n g whether o r n o t t o r e p o s i t i o n any o f t h e knots i s printed.) KNOT TO BE VARIED?  1  NEW POSITI'ONFFOR KNOT ( 1)? .3 DO YOU WISH TO VARY ANOTHER KNOT? KNOT TO BE VARIED?  YES  4  NEW POSITION FOR KNOT ( 4)? .7 DO YOU WISH TO VARY ANOTHER KNOT?  NO  (A f i x e d knot l e a s t squares, a p p r o x i m a t i o n ( i s computed u s i n g t h e new knot  - 106  set.  The r e s u l t i n g l e a s t squares  -  error information i s printed.)  THE LEAST SQUARES.ERROR OF THE  PREVIOUS APPROXIMATION WAS  THE  NEW  LEAST SQUARES ERROR OF THE  DO YOU  WISH TO SEE THE  NEW  APPROXIMATION.IS  APPROXIMATION?  9;576172E-03  YES  ( I f the answer t o t h i s q u e s t i o n i s a f f i r m a t i v e the new s p l i n e curve r e p l a c e s the one  1.537876E-03  knot  s e t and  c u r r e n t l y b e i n g d i s p l a y e d on the g r a p h i c s  terminal.) DO YOU  WISH TO CONTINUE WITH. THE  (The v a l u e s of the new  NEW' APPROXIMATION?  approximation  YES  are t r a n s f e r r e d t o the  proper  arrays.) DO YOU  WISH TO RETAIN THE  TITLE FOR (A t i t l e  APPROXIMATION RESULTS?  YES  PLOT? T E S T DATA'' i s p r i n t e d on each graph and p r i n t o u t .  These are a l s o  labelled  w i t h a number which i s the number o f h a r d c o p i e s o b t a i n e d thus f a r d u r i n g the run.  T h i s method enables  p a r t i c u l a r run to be u n i q u e l y DO YOU  WISH TO VARY THE  DO YOU  WISH TO OPTIMIZE THE  c o r r e s p o n d i n g p l o t s and p r i n t o u t s from a identified.)  KNOT SET?  NO  KNOT SET?  NO  (A r e q u e s t i s p r i n t e d as to whether oranot t o o p t i m i z e the c u r r e n t knot set.  WARNING:  o p t i m i z a t i o n i s time-consuming and  expensive.  I t is.  b e s t done n o n - i n t e r a c t i v e l y . ) DO YOU  WISH TO RESTART WITH A NEW  KNOT SET?  NO  (A request, i s i s s u e d whether o r not to r e s t a r t w i t h an expanded or c o n t r a c t e d knot s e t . r e p l o t t e d and DO YOU  I f the answer i s a f f i r m a t i v e the d a t a s e t i s  the program r e s t a r t s from the  WISH TO RESTART WITH A NEW  DATA SET?  beginning.) NO  ( I f the answer i s a f f i r m a t i v e , the program reads another d a t a s e t  from  - 107 -  unit  4 and r e t u r n s t o t h e b e g i n n i n g o f the program).  PLOTTING WILL TAKE APPROX. 2:MINS  32 SECONDS  STOP 0 EXECUTION  TERMINATED  $C0PY -PRINT *PRINT-* $RUN PL0T:Q PAR=^PL0T $SIG  C.  Output  Approximation a file.  i n f o r m a t i o n c o r r e s p o n d i n g t o a p l o t i s put i n  The format o f each approximation p r i n t o u t i s as f o l l o w s :  TITLE  NUMBER  ABSCISSAE  ORDINATES  FITTED ORDINATES  RESIDUALS  f o r one d a t a s e t  THE LEAST SQUARES, ERROR IS KNOT  LOCATION  This prints  POLYNOMIAL POWER  POLYNOMIAL  the c o e f f i c i e n t s o f every power f o r each p i e c e -  wise p o l y n o m i a l between each knot p a i r . c o r r e s p o n d t o t h i s , output.  P l o t s o b t a i n e d f o r the p l o t t e r  They can be matched by the t i t l e and l a b e l  number.  Restrictions  The number o f d a t a p o i n t s must be l e s s The degree  COEFFICIENT  that  o f t h e s p l i n e msut be l e s s than  The number o f knots must be l e s s t h a t  42 .  200 . 3 .  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0058928/manifest

Comment

Related Items