% PPALLDEM	Demonstrate form and usage of piecewise polynomials.
% Copyright (c) 1990-92 by Carl de Boor and The MathWorks, Inc.
                                            % latest change: September 28, '90
clc;clg;echo on;home                             
%       A piecewise polynomial, or   pp   for short, is characterized by its
%  b r e a k  sequence,   break   say, and the   c o e f f i c i e n t 
%  array,   coefs   say, of the local power form of its polynomial pieces.
%
%     The break sequence is assumed to be strictly increasing,
%               break(1) < break(2) < ... < break(l+1),
%  with  l  the number of polynomial pieces which make up  pp .
%
%      While these polynomials may be of varying degrees, they are
%  all recorded as polynomials of the same  o r d e r   k , i.e., the
%  coefficient array  c  is of size  (l,k) , with  coefs(j,:)  containing
%  the  k  coefficients in the local power form for the j-th polynomial
%  piece.
pause                    % Touch any key to continue
%
%       With this, the precise description of  pp  in terms of the break
%   sequence  break  and the coefficient array  coefs  is
% 
%   pp(t) = polyval(coefs(j,:),t-break(j))  for   break(j) <= t < break(j+1) ,
%
%  Here,  polyval(a,x)  is the MATLAB function; it returns the number 
%     sum{a(j)x^(k-j):j=1,..,k} = a(1)*x^(k-1) + a(2)*x^(k-2) + ... + a(k)*x^0 .
%
%  For  t not in [break(1),break(l+1)),  pp(t)  is defined by extending the 
% first, resp. last, polynomial piece.  
%
pause                    % Touch any key to continue
%    A pp is usually constructed in some M-file, through a process of
%  interpolation or approximation, or conversion from some other form (e.g.,
%  from a spline), and is output as a single (row) vector. But it is also
%  possible to make one up from scratch, using the statement

%           pp=ppmak(breaks,coefs).

pause                    % Touch any key to continue

%    Let's try that.

pp=ppmak([-5:-1],[-22:-11])
echo off;
ttype=10;td=1;tl=4;tbreak=[-5:-1];tk=3;tcoef=[-22:-20;-19:-17;-16:-14;-13:-11];
zcheck(pp-[ttype td tl tbreak tk tcoef(:)']);
% Check also integration and differentiation. 
% Evaluation and taking a piece are checked below.
zcheck(pp-fnder(fnint(pp)));echo on
pause                    % Touch any key to continue

%    Back comes this row vector in which you could actually find all those
%  numbers input, since they are all negative and distinct. But you thereby
%  also see some other numbers. In fact, it is not worthwhile to wonder just
%  how the pp is described in this vector.

pause                    % Touch any key to continue

%    If you really want to know the constituent parts of pp, use
%           [breaks,coefs,l,k]=ppbrk(pp)

pause(5);

[breaks,coefs,l,k]=ppbrk(pp)
pause                    % ... and there they are. Touch any key to continue
echo off
zcheck(breaks-tbreak);zcheck(coefs-tcoef);zcheck([l k]-[tl tk]);echo on
%    Here are some operations you can perform on  pp.

%       v=fnval(pp,x)     evaluates,

%       dpp=fnder(pp)     differentiates,

%       dpp=fnint(pp)     integrates,

%       pj=pppce(pp,j)   pulls out the j-th polynomial piece.

%       pc=ppcut(pp,[a b])  pulls out the part relevant for the interval
%                           [a,b] .

%       (values = ) fnplt(pp)  plots pp between its extreme breaks.

%     For example, here is a plot of the particular pp we just made up.

pause                     % Touch any key to continue

x=[-55:-5]/10; plot(x,fnval(pp,x),'-.');pause


%   Here, added to the plot, are the breaklines

pause                     % Touch any key to continue

yy=axis;hold on;
for j=1:l+1;plot([breaks(j) breaks(j)],yy(3:4));end;pause

%   And here is, added to that plot, the plot of its third piece

pause                     % Touch any key to continue

plot(x,fnval(pppce(pp,3),x));pause
hold off; echo off
zcheck(-14+(x+3).*(-15 + (x+3)*(-16))-fnval(pppce(pp,3),x));echo on

%     While the pp-representation of a pp is efficient for evaluation, the
%  c o n s t r u c t i o n  of a pp from some data is usually more
%  efficiently handled by determining first its
%   B - r e p r e s e n t a t i o n , i.e., its representation as a linear
%  combination of B-splines.	For this, look at      spalldemo  .

pause                     % Touch any key to finish
