step_response.m

一个matlab的将军模型

function xfinal = step_response(A,Ainv,b,tf)

% Compute the `step response` of the affine system.
%
% Syntax:
%   "xfinal = step_response(A,Ainv,b,tf)"
%
% Description:
%   Compute the final value of the state vector "x" at the final time "tf"
%   for the affine dynamics "dx/dt = A*x + b". The inputs are
%
%   * "A": the system matrix
%
%   * "Ainv": the inverse of "A" if it exists, otherwsie it should be
%     "[]"
%
%   * "b": constant input vector for the affine dynamics
%
%   * "tf": final time for the step response
%
%   The final value of "x" is computed by evaluating the expression
%
%
%
%   "x(tf) = e^{A*tf}*x(0) + e^{A*tf} * integral_{s=0}^{s=tf} e^{-A*s}*b ds"
%
%
%
% Implementation:
%   If "A" is invertible, then the above integral reduces to
%
%
%
%   "x(tf) = (e^{A*tf}-I)*Ainv*b".
%   
%   
%
%   Otherwise, the value of "x(tf)" can be obtained by calling the ODE
%   solver to integrate the system equation "dx/dt = A*x + b" using the `ode
%   file` "affine.m" from time "0" to "tf".

global GLOBAL_ODE_PAR

if (tf == 0) | all(b == 0)
  xfinal = zeros(size(b));
else
  if ~isempty(Ainv)
    % if A is invertible, then use the closed form solution
    xfinal = (expm(A*tf)-eye(size(A)))*Ainv*b;
  else
    % otherwise, perform numerical integration to get the final value of
    % the step response
    % Alternative method, use ode45 to integrate the step response
  
    [T,Y] = ode45('affine',[0 tf],zeros(size(A,1),1),GLOBAL_ODE_PAR,A,b);
    xfinal = Y(length(T),:)';
  end
end
return