Algorithme de la Méthode de Hamming pour résoudre les équations différentielles vectorielles en Matlab

Algorithme de la Méthode de Hamming pour résoudre les équations différentielles vectorielles en Matlab :

function [t,y] = ode_Ham(f,tspan,y0,N,KC,varargin) % Hamming method to solve vector d.e. y’(t) = f(t,y(t)) % for tspan = [t0,tf] and with the initial value y0 and N time steps % using the modifier based on the error estimate depending on KC = 1/0 if nargin < 5, KC = 1; end %with modifier by default if nargin < 4 | N <= 0, N = 100; end %default maximum number of iterations if nargin < 3, y0 = 0; end %default initial value y0 = y0(:)’; end %make it a row vector h = (tspan(2)-tspan(1))/N; %step size tspan0 = tspan(1)+[0 3]*h; [t,y] = ode_RK4(f,tspan0,y0,3,varargin{:}); %Initialize by Runge-Kutta t = [t(1:3)’ t(4):h:tspan(2)]’; for k = 2:4, F(k - 1,:) = feval(f,t(k),y(k,:),varargin{:}); end p = y(4,:); c = y(4,:); h34 = h/3*4; KC11 = KC*112/121; KC91 = KC*9/121; h312 = 3*h*[-1 2 1]; for k = 4:N p1 = y(k - 3,:) + h34*(2*(F(1,:) + F(3,:)) - F(2,:)); %Eq.(6.4.9a) m1 = p1 + KC11*(c - p); %Eq.(6.4.9b) c1 = (-y(k - 2,:) + 9*y(k,:) +... h312*[F(2:3,:); feval(f,t(k + 1),m1,varargin{:})])/8; %Eq.(6.4.9c) y(k+1,:) = c1 - KC91*(c1 - p1); %Eq.(6.4.9d) p = p1; c = c1; %update the predicted/corrected values F = [F(2:3,:); feval(f,t(k + 1),y(k + 1,:),varargin{:})]; end

Algorithme de la Méthode d'Adams-Bashforth-Moulton pour résoudre les équations différentielles vectorielles en Matlab

Algorithme de la Méthode d'Adams-Bashforth-Moulton pour résoudre les équations différentielles vectorielles en Matlab :

function [t,y] = ode_ABM(f,tspan,y0,N,KC,varargin) %Adams-Bashforth-Moulton method to solve vector d.e. y’(t) = f(t,y(t)) % for tspan = [t0,tf] and with the initial value y0 and N time steps % using the modifier based on the error estimate depending on KC = 1/0 if nargin < 5, KC = 1; end %with modifier by default if nargin < 4 | N < = 0, N = 100; end %default maximum number of iterations y0 = y0(:)’; %make it a row vector h = (tspan(2) - tspan(1))/N; %step size tspan0 = tspan(1)+[0 3]*h; [t,y] = rk4(f,tspan0,y0,3,varargin{:}); %initialize by Runge-Kutta t = [t(1:3)’ t(4):h:tspan(2)]’; for k = 1:4, F(k,:) = feval(f,t(k),y(k,:),varargin{:}); end p = y(4,:); c = y(4,:); KC22 = KC*251/270; KC12 = KC*19/270; h24 = h/24; h241 = h24*[1 -5 19 9]; h249 = h24*[-9 37 -59 55]; for k = 4:N p1 = y(k,:) +h249*F; %Eq.(6.4.8a) m1 = pk1 + KC22*(c-p); %Eq.(6.4.8b) c1 = y(k,:)+ ... h241*[F(2:4,:); feval(f,t(k + 1),m1,varargin{:})]; %Eq.(6.4.8c) y(k + 1,:) = c1 - KC12*(c1 - p1); %Eq.(6.4.8d) p = p1; c = c1; %update the predicted/corrected values F = [F(2:4,:); feval(f,t(k + 1),y(k + 1,:),varargin{:})]; End

Algorithme de la Méthode de Runge-Kutta pour résoudre les équations différentielles vectorielles en Matlab

Algorithme de la Méthode de Runge-Kutta pour résoudre les équations différentielles vectorielles en Matlab :

function [t,y] = ode_RK4(f,tspan,y0,N,varargin) %Runge-Kutta method to solve vector differential eqn y’(t) = f(t,y(t)) % for tspan = [t0,tf] and with the initial value y0 and N time steps if nargin < 4 | N <= 0, N = 100; end if nargin < 3, y0 = 0; end y(1,:) = y0(:)’; %make it a row vector h = (tspan(2) - tspan(1))/N; t = tspan(1)+[0:N]’*h; for k = 1:N f1 = h*feval(f,t(k),y(k,:),varargin{:}); f1 = f1(:)’; %(6.3.2a) f2 = h*feval(f,t(k) + h/2,y(k,:) + f1/2,varargin{:}); f2 = f2(:)’;%(6.3.2b) f3 = h*feval(f,t(k) + h/2,y(k,:) + f2/2,varargin{:}); f3 = f3(:)’;%(6.3.2c) f4 = h*feval(f,t(k) + h,y(k,:) + f3,varargin{:}); f4 = f4(:)’; %(6.3.2d) y(k + 1,:) = y(k,:) + (f1 + 2*(f2 + f3) + f4)/6; %Eq.(6.3.1) end

Algorithme de la Méthode de Heun pour résoudre les équations différentielles vectorielles en Matlab

Algorithme de la Méthode de Heun pour résoudre les équations différentielles vectorielles en Matlab :

function [t,y] = ode_Heun(f,tspan,y0,N) %Heun method to solve vector differential equation y’(t) = f(t,y(t)) % for tspan = [t0,tf] and with the initial value y0 and N time steps if nargin<4 | N <= 0, N = 100; end if nargin<3, y0 = 0; end h = (tspan(2) - tspan(1))/N; %stepsize t = tspan(1)+[0:N]’*h; %time vector y(1,:) = y0(:)’; %always make the initial value a row vector for k = 1:N fk = feval(f,t(k),y(k,:)); y(k+1,:) = y(k,:)+h*fk; %Eq.(6.2.3) y(k+1,:) = y(k,:) +h/2*(fk +feval(f,t(k+1),y(k+1,:))); %Eq.(6.2.4) end

Algorithme de la Méthode d'Euler pour résoudre les équations différentielles de premier ordre en Matlab

Algorithme de la Méthode d'Euler pour résoudre les équations différentielles de premier ordre en Matlab :

%nm610: Euler method to solve a 1st-order differential equation clear, clf a = 1; r = 1; y0 = 0; tf = 2; t = [0:0.01:tf]; yt=1- exp(-a*t); %Eq.(6.1.5): true analytical solution plot(t,yt,’k’), hold on klasts = [8 4 2]; hs = tf./klasts; y(1) = y0; for itr = 1:3 %with various step size h = 1/8,1/4,1/2 klast = klasts(itr); h = hs(itr); y(1)=y0; for k = 1:klast y(k + 1) = (1 - a*h)*y(k) +h*r; %Eq.(6.1.3): plot([k - 1 k]*h,[y(k) y(k+1)],’b’, k*h,y(k+1),’ro’) if k < 4, pause; end end end

Algorithme de la Méthode de la Sécante en Matlab

Algorithme de la Méthode de la Sécante en Matlab :

function [x,fx,xx] = secant(f,x0,TolX,MaxIter,varargin) % solve f(x) = 0 by using the secant method. %input : f = ftn to be given as a string ’f’ if defined in an M-file % x0 = the initial guess of the solution % TolX = the upper limit of |x(k) - x(k - 1)| % MaxIter = the maximum # of iteration %output: x = the point which the algorithm has reached % fx = f(x(last)), xx = the history of x h = 1e-4; h2 = 2*h; TolFun=eps; xx(1) = x0; fx = feval(f,x0,varargin{:}); for k = 1: MaxIter if k <= 1, dfdx = (feval(f,xx(k) + h,varargin{:})-... feval(f,xx(k) - h,varargin{:}))/h2; else dfdx = (fx - fx0)/dx; end dx = -fx/dfdx; xx(k + 1) = xx(k) + dx; %Eq.(4.5.2) fx0 = fx; fx = feval(f,xx(k+1)); if abs(fx) < TolFun | abs(dx) < TolX, break; end end x = xx(k + 1); if k == MaxIter, fprintf(’The best in %d iterations\n’,MaxIter), end

Algorithme de la Méthode de la Position Fausse en Matlab

Algorithme de la Méthode de la Position Fausse en Matlab :

function [x,err,xx] = falsp(f,a,b,TolX,MaxIter) %bisct.m to solve f(x)=0 by using the false position method. %input : f = ftn to be given as a string ’f’ if defined in an M-file % a/b = initial left/right point of the solution interval % TolX = upperbound of error(max(|x(k)-a|,|b-x(k)|)) % MaxIter = maximum # of iterations %output: x = point which the algorithm has reached % err = max(x(last)-a|,|b-x(last)|) % xx = history of x TolFun = eps; fa = feval(f,a); fb=feval(f,b); if fa*fb > 0, error(’We must have f(a)f(b)<0!’); end for k = 1: MaxIter xx(k) = (a*fb-b*fa)/(fb-fa); %Eq.(4.3.1) fx = feval(f,xx(k)); err = max(abs(xx(k) - a),abs(b - xx(k))); if abs(fx) < TolFun | err<TolX, break; elseif fx*fa > 0, a = xx(k); fa = fx; else b = xx(k); fb = fx; end end x = xx(k); if k == MaxIter, fprintf(’The best in %d iterations\n’,MaxIter), end