e851.m

matlab算法集 matlab算法集

%------------------------------------------------------------------
% Example 8.5.1: Predictor-Corrector Method                       
%------------------------------------------------------------------
   clc
   clear
   n = 2;
   m = 11;
   alpha = 0;	
   beta = 2;
   h = (beta-alpha)/(m-1);
   x0 = [0 1]';
   x1 = zeros(n,1);
   x = zeros(n,1);
   e = zeros(m,1);
   E = zeros(m,5);
   X = zeros (m,n);
   X0 = zeros (4,n);

% Exact solution 

   fprintf ('Example 8.5.1: Predictor-Corrector Method\n');
   t = linspace (alpha,beta,m)';
   for k = 1 : m
      X(k,1) = t(k)*exp(-t(k));
      X(k,2) = (1-t(k))*exp(-t(k));
   end
   graphxy (t,X,'Exact Solution','t','x');  
 
% Euler solution 
   
   x = x0;
   for k = 1 : m-1
      x1 = funf851(t(k),x);
      for j = 1 : n
         x(j) = x(j) + h*x1(j);
         e(j) = x(j) - X(k+1,j);
      end
      E(k,1) = norm (e,inf);
   end

% Runge-Kutta solution 
   
   x = x0;
   for k = 1 : m-1
      [x1,E(k,2)] = rkfstep (x,t(k),t(k+1),'funf851');
      x = x1;
      for j = 1 : n
         e(j) = x(j) - X(k+1,j);
      end
      E(k,3) = norm (e,inf);
   end

% Adams predictor-corrector solution 
   
   for k = 1 : 4
      j = 5 - k;
      X0(k,1) = t(j)*exp(-t(j));
      X0(k,2) = (1-t(j))*exp(-t(j));
   end

   for k = 4 : m-1 
      [X0,E(k,4)] = predcorr (X0,t(k),t(k+1),'funf851');
      for j = 1 : n
         e(j) = X0(1,j) - X(k+1,j);
      end
      E(k,5) = norm (e,inf);
   end

% Display results

   fprintf ('\n t     Euler     RKF(local)    RKF(global)');
   fprintf ('     ABM(local)    ABM(global)\n');
   fprintf ('--------------------------------------------');
   fprintf ('-----------------------------\n');             
   for k = 1 : m-1
      fprintf ('%2.1f %9.6f  %12.6e  %12.6e  %12.6e  %12.6e\n',...
               t(k+1),E(k,:));
   end
   fprintf ('--------------------------------------------');
   fprintf ('-----------------------------\n');             
   wait
%------------------------------------------------------------------