%
%   CHAIN.M   Plotting and Function Evaluation with MATLAB    
%
%   This demo gives a simple example of the use of plots in MATLAB.      
%   It also shows how to use MATLAB to perform simple functional evaluations    
%   using vectorized operations and loops. See the MATLAB HELP files for more 
%   description.
%                                                                               
%   This is an example of plotting a simple 2-d function, T(a,b), versus a and b.
%   In this case, T represents the time it takes a chain to slide off a frictionless 
%   table.  The variables a and b represent the total length of the chain and the
%   initial length of the overhang, respectively.  The function T(b) is drawn for
%   several values of the variable a (giving a family of curves).
%   Solving the differential equation that describes this system gives the following 
%   result:  
%           T(a,b) = sqrt(a/g) * ln((a+sqrt(a*a-b*b))/b)
% 
%   The goal here is simply to evaluate and plot this function for various values of 
%   a and b.
%
%   File prepared by J. R. White, UMass-Lowell (July 2003)
%

%
%   getting started
      clear all, close all, nfig = 0;
%
%   parameter definition 
%     a = length of chain
%     b = initial length of overhang                                             
%                                                                               
      g = 9.8;                 % gravitational accel  (m/sec^2)
%
      ao = 5;                  % start value for a (m)
      af = 25;                 % end value for a (m)
      Na = 5;                  % number of values of a to use 
      a = linspace(ao,af,Na);  % actual vector of a values
%
      bo = 0.1;                % start value for b (fraction of a)
      bf = 0.9;                % end value for b (fraction of a)
      Nb = 10;                 % number of values of b to use for each a value
% 
      T = zeros(Nb,Na);        % initialize array for storage of results 
      b = zeros(Nb,Na);        % initialize array to keep Nb b values for Na a values
      txt = cell(Na,1);        % cell array to store text for legend on plot
%
%   loop over all values of a
      for n = 1:Na
        xa = a(n);
        b(:,n) = linspace(bo*xa,bf*xa,Nb)';
        T(:,n) = sqrt(xa/g)*log((xa+sqrt(xa*xa-b(:,n).*b(:,n)))./b(:,n));
        txt(n) = {['a',' = ',num2str(a(n)),' m']};  
      end
%
%   set color and marker code for creating plots
      Ncm = 6;
      scm = ['r-';      % red solid
             'g:';      % green dotted
             'b-';      % blue solid
             'm:';      % magenta dotted
             'c-';      % cyan solid
             'y:'];     % yellow dotted
%                                                                               
%   plot family of curves (one for each value of a)                                       
      nfig = nfig+1;   figure(nfig)
      for j = 1:Na
        plot(b(:,j),T(:,j),scm(j,:),'LineWidth',2), grid on, hold on
      end 
      title('Chain:  Time to Fall Off Table for Several Chain Lengths')                               
      xlabel('b-value (m)'),ylabel('T(a,b)-value (sec)')
      legend(txt)
%
%   place some additional text on the plot
      v = axis;   xt = v(1)+0.52*(v(2)-v(1));
      yt1 = v(3) + 0.07*(v(4)-v(3)); 
      yt2 = v(3) + 0.03*(v(4)-v(3));
      text(xt,yt1,'a = length of chain (m)')
      text(xt,yt2,'b = length of initial overhang (m)')
      hold off
%
%  end simulation