int1.m

一个时滞系统的工具箱

function rez=int1(sb,se,t,tt,xt,initfun,F,num);
%INT1  calculates the integrals of function F(t,s)*y(s) 
%      from -tau(t) to 0 over variable s, included in the right
%      hand of a system of functional differential equations
%      to be solved numerically.
%      For calculations INT1 uses initial function as well arrays
%      tt, xt used for construction of degenerate cubic splines.
%      INT1 implements Simpson's method of the integration.
%      INT1 is used in special functions where user describes
%      the system in order to solve it by the program dde45.
%      For details see descriptions.

if nargin < 8 
  num=1:size(xt,2);
end;

% Initialization
M=20;           % The number of partitions
h=2*(se-sb)/M;
j=0:M; 
tj=sb+0.5*h*j;  % Array of internal division of the segment [sb, se]
x1=zeros(length(num), M+1);

% Calculation of integrand at the internal points
% of the segment [sb, se].
k=find(t+tj<=tt(1));
len=length(k);
if ~isempty(k) 
   for i=1:len
       s=t+tj(k(i));
       temp1=eval(initfun);
       x1(:,k(i))=temp1(num,:);
   end;
end;
k=(len+1):(M+1);
if ~isempty(k), 
  x1(:,k)=Spline3(tt,xt(:,num),tj(k)+t)'; 
end;

s=tj(1);
F1=eval(F);
x2=zeros(size(F1,1), M+1);
x2(:,1)=F1*x1(:,1);
for i=2:M+1
  s=tj(i);
  F1=eval(F);
  x2(:,i)=F1*x1(:,i);
end;

% The final summation
summ=x2(:,1)+4*sum(x2(:,2:2:M)')'+2*sum(x2(:,3:2:M-1)')'+x2(:,M+1);
rez=summ*h/6;