home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 10.ddi / CONTROL.DI$ / INITIAL.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  2.2 KB  |  77 lines

  1. function [yout,x,t] = initial(a,b,c,d,x0,t)
  2. %INITIAL Initial condition response of continuous-time linear systems.
  3. %    INITIAL(A,B,C,D,X0) plots the time response of the linear system
  4. %        .
  5. %        x = Ax + Bu
  6. %        y = Cx + Du
  7. %
  8. %    due to an initial condition on the states.  The time vector is 
  9. %    automatically determined based on the system poles and zeros.  
  10. %
  11. %    INITIAL(A,B,C,D,X0,T) plots the initial condition response for the
  12. %    times specified in the vector T.  The time vector must be 
  13. %    regularly spaced.  When invoked with left hand arguments:
  14. %    
  15. %        [Y,X,T] = INITIAL(A,B,C,D,X0,...)
  16. %
  17. %    returns the output and state responses (Y and X), and the time 
  18. %    vector (T).  No plot is drawn on the screen.  The matrix Y has as
  19. %    many columns as outputs and one row for element in T.  Similarly,
  20. %    the matrix X has as many columns as states and length(T) rows.
  21. %    
  22. %    See also: IMPULSE,STEP,LSIM, and DINITIAL.
  23.  
  24. %    Clay M. Thompson  7-6-90
  25. %    Revised ACWG 6-21-92
  26. %    Copyright (c) 1986-93 by the MathWorks, Inc.
  27.  
  28. if nargin==0, eval('exresp(''initial'',1)'), return, end
  29.  
  30. error(nargchk(4,6,nargin));
  31. error(abcdchk(a,b,c,d));
  32.  
  33. [ny,nu] = size(d);
  34. if (nu*ny==0)|isempty(a),
  35.   x = []; t = []; if nargout~=0, yout=[]; end, return, 
  36. end
  37.  
  38. [nx,na] = size(a);
  39. if nargin==4        % No x0 specified, use random initial condition.
  40.   x0 = randn(nx,1);
  41.  
  42. else
  43.   x0 = x0(:);        % Make sure x0 is a column vector
  44.   if length(x0)~=nx,
  45.     error('The length of X0 must be equal to the number of states.');
  46.   end
  47. end
  48.  
  49. if (nargin==4)|(nargin==5)  % Workout time vector if not supplied.
  50.   % The next two constants control the precision of the plot
  51.   % and the time interval of the plot.
  52.   st=0.005; % Set settling time bound  = 0.5%
  53.   precision=30; % Show approx 30 points for simple graph
  54.   t=timvec(a,b,c,x0,st,precision);
  55. end
  56.  
  57. % --- Simulation ---
  58. dt = t(2)-t(1);
  59. [aa,bb] = c2d(a,b,t(2)-t(1));
  60. n = length(t);
  61. [nx,nu] = size(b);
  62. x = ltitr(aa,bb,zeros(n,nu),x0(:));
  63. y=x*c.';
  64. tout = t;
  65.  
  66. if nargout==0,        % If no output arguments, plot graph
  67.   plot(t, y, [t(1),t(length(t))],[0;0],'w:')
  68.   xlabel('Time (secs)'), ylabel('Amplitude')
  69.   return % Suppress output
  70. end
  71. yout = y; 
  72.  
  73.  
  74.  
  75.  
  76.  
  77.