    function ShowPoly()
  % Examines the Chebychev, Legendre, and Hermite Polynomials
  clc
  close all
  % Check out EvalPoly, PlotPoly, ScalePoly, and DerPoly...
  P = MakePoly([-6 5 5 -5 1])
  PlotPoly(P,-1.2,3.2)
  Three_times_P = ScalePoly(3,P)
  Derivative_of_P = DerPoly(P)
  % Check out SubPoly and MultPoly...
  P = MakePoly([10 20 30 40])
  Q = MakePoly([1 2 3 4 5 6 7])
  PminusQ = SubPoly(P,Q)
  QminusP = SubPoly(Q,P)
  PtimesQ = MultPoly(P,Q)
  QtimesP = MultPoly(Q,P)
  [PtimesQ.a' QtimesP.a']
  
  % Check out MakeChebychev...
  figure
  fprintf('\n\nChebychev Polynomial Coefficients...\n')
  T = MakeChebychev(10);
  for k=1:8
      fprintf('\nk = %1d:  ',k)
      fprintf('%8.3f  ',T(k).a)
      PlotPoly(T(k),-1,1)
      title(sprintf('Chebychev Polynomial: Degree = %1d',k),'Fontsize',16)
      pause(.5)
  end
  
  % Check out MakeLegendre...
  figure
  fprintf('\n\nLegendre Polynomial Coefficients...\n')
  P = MakeLegendre(10);
  for k=1:8
      fprintf('\nk = %1d:  ',k)
      fprintf('%8.3f  ',P(k).a)
      PlotPoly(P(k),-1,1)
      title(sprintf('Legendre Polynomial: Degree = %1d',k),'Fontsize',16)
      pause(.3)
  end

  % Check out MakeHermite...
  figure
  fprintf('\n\nHermite Polynomial Coefficients...\n')
  H = MakeHermite(10);
  for k=1:8
      fprintf('\nk = %1d:  ',k)
      fprintf('%8.3f  ',H(k).a)
      PlotPoly(H(k),-3,3)
      title(sprintf('Hermite Polynomial: Degree = %1d',k),'Fontsize',16)
      pause(.3)
  end
 
  function P = MakePoly(a)
% a is a length n vector.
% P is a two-field structure that encodes a polynomial p(x).
% P.a specifies the coefficients:
%   p(x) = P.a(1) + P.a(2)x + a(3)x^2 + ... + P.a(n)x^(n-1)
% P.deg is the degree of the polynomial (i.e.,P.deg == n-1)
deg = length(a)-1;
P = struct('a',a,'deg',deg);

  function PlotPoly(P,L,R)
% P encodes a polynomial p(x) and L<R
% Opens a figure window and displays p(x) across the interval [L,R] 
x = linspace(L,R)';
y = EvalPoly(P,x);
plot(x,y,[L,R],[0 0],'r--')

  function yVals = EvalPoly(P,xVals)
% P encodes a polynomial p(x)
% xVals is a column n-vector
% yVals is a column n-vector with yVals(k) = p(xVals(k)), k=1:n
n = length(xVals);
z = ones(n,1);
yVals = P.a(1)*z;
for k=2:P.deg+1
    z = xVals.*z;
    yVals = yVals + P.a(k)*z;
end

  function R = ScalePoly(alfa,P)
% alfa is a scalar and P encodes a polynomial p(x)
% R encodes the polynimal alfa*p(x)
R = P;
R.a = alfa*R.a;

  function R = DerPoly(P)
% P encodes a polynomial p(x) 
% R encodes its derivative
if P.deg==0
    R = MakePoly(0);
else
    R = MakePoly(zeros(1,P.deg));
    for k=1:P.deg
        R.a(k) = P.a(k+1)*k;
    end
end

  function R = SubPoly(P,Q)
% P encodes the polynomial p(x)
% Q encodes the polynomial q(x)
% R encodes the polynomial r(x) = p(x)-q(x)
d = max(P.deg,Q.deg);
c = zeros(1,d+1);
for k=1:d+1
    if k<=P.deg+1
        c(k) = c(k) + P.a(k);
    end
    if k<=Q.deg+1
        c(k) = c(k) - Q.a(k);
    end
end
R = MakePoly(c);

  function R = MultPoly(P,Q)
% P encodes the polynomial p(x)
% Q encodes the polynomial q(x)
% R encodes the polynomial r(x) = p(x)*q(x)
d = P.deg + Q.deg;
R = MakePoly(zeros(1,d+1));
for k=1:d+1
    for j= max(1,k-Q.deg):min(k,P.deg+1)    
        R.a(k) = R.a(k) + P.a(j)*Q.a(k-j+1);
    end
end

    function T = MakeChebychev(n)
 % T is a length-n structure array of polynomials
 % T(k) encodes the k-th Chebychev polynomial    
 T(1)  = MakePoly([0 1]);
 T(2)  = MakePoly([-1 0 2]);
 Two_x = MakePoly([0 2]);
 for k=2:n-1
     T(k+1) = SubPoly(MultPoly(Two_x,T(k)),(T(k-1)));
 end
 
   function P = MakeLegendre(n) 
 % P is a length-n structure array of polynomials
 % P(k) encodes the k-th Legendre polynomial 
 P(1)  = MakePoly([0 1]);
 P(2)  = MakePoly([-1/2 0 3/2]);
 for k=2:n-1
     Z = ScalePoly((2*k+1)/(k+1),MakePoly([0 1]));
     P(k+1) = SubPoly(MultPoly(Z,P(k)),ScalePoly(k/(k+1),P(k-1)));
 end
  
   function H = MakeHermite(n) 
 % H is a length-n structure array of polynomials
 % H(k) encodes the k-th Hermite polynomial 
 H(1) = MakePoly([0 1]);
 H(2) = MakePoly([-1 0 1]);
 x = MakePoly([0 1]);
 for k=2:n-1
    H(k+1) = SubPoly(MultPoly(x,H(k)),DerPoly(H(k)));
 end

    

