%
%   PSERIES_EX0.M   Evaluate the Taylor Series Expansion
%                        for y(x) = ln(1+x)
%
%   This demo simply illustrates how to evaluate an infinite series
%   expansion.  It also illustrates the importance of the radius of
%   convergence of a particular series.  The series used here is
%     ln(1+x) = x - (1/2)x^2 + (1/3)x^3 - ...
%             = sum (-1)^m * x^(m+1 )/ (m+1)  for m = 0 to inf
%   and a theoretical analysis shows that the radius of convergence, 
%   R, is unity (meaning that the series will converge for |x| < 1).
%
%   File prepared by J. R. White, UMass-Lowell (Oct. 2003)
%

%
%   getting started
      clear all,  close all,   nfig = 0;
%
%   define range of x values 
      iopt = menu('Choose Option?','Input single x value', ...
                  'Use default range 0 <  x < 1.2');
      disp('  ')
      if iopt == 1,  
        x = input(' Enter x-value:  ');
      else
        x = 0:0.05:1.2;  Nx = length(x);
      end 
%
%   evaluate y(x) for single value
      maxT = 101;  tol = 1e-6;     
%
%   Single X Value: evaluate y1(x) solution using infinite series expansion
      if iopt == 1
        m = 0;  mrerr = 1.0;
        T = x;  y1 = T; 
        while mrerr > tol  &  m < maxT
          r = - ((m+1)/(m+2))*x;   T = r*T;  
          y1 = y1 + T;   
          m = m+1;  
          if y1 ~= 0,   mrerr = abs(T/y1);   end
        end
        NT1 = m;
%
        fprintf(1,'\n Results for Series Evaluation of ln(1 + %5.2f)\n',x);
        fprintf(1,'    exact result                 = %8.5f \n',log(1+x));
        fprintf(1,'    approx result (for NT = %3d) = %8.5f \n',NT1,y1);
      end
%
%
%   Vector of X Values: evaluate y2(x) solution using infinite series expansion
      if iopt == 2
        m = 0;  mrerr = 1.0;
        T = x;  y2 = T; 
        while mrerr > tol  &  m < maxT
          r = - ((m+1)/(m+2))*x;   T = r.*T;  
          y2 = y2 + T;   
          m = m+1;  
          ii = find(y2);  rerr = T(ii)./y2(ii);   mrerr = max(abs(rerr));
        end
        NT2 = m;
%
        ye = log(1+x);
        fprintf(1,'\n Results for Series Evaluation of ln(1 + x)\n');
        fprintf(1,'     x         exact         approx          ratio \n')
        for i = 2:Nx
          fprintf(1,'  %5.3f   %12.3e   %12.3e   %12.3e\n', ...
                   x(i),ye(i),y2(i),ye(i)/y2(i));
        end
%
      nfig = nfig+1;   figure(nfig)
      plot(x,ye,'r--',x,y2,'b-','LineWidth',2), grid
      v = [0 1.2 -2 2];  axis(v)
      title('PSeries\_Ex0:  Power Series Example -- y(x) = ln(1+x)'); 
      xlabel('x-values'),ylabel('y-values')
      legend('Exact Soln',['Power Series (NT = ',num2str(NT2),')'])
      end
%
%   end of program
