%程序1
- arg1 = 2;
- arg2 = 1;
- [T,Y] = ode45('vdp1000',[0 10],[2 0], [], arg1, arg2);
- plot(T,Y(:,1),'-o');
%%程序2
- function dy = vdp1000(t, y, flag, arg1, arg2)
- dy = zeros(2,1); % a column vector
- dy(1) = y(2);
- dy(2) = arg1*(arg2 - y(1)^2)*y(2) - y(1);
%%ode5
- function Y = ode5(odefun,tspan,y0,varargin)
- %ODE5 Solve differential equations with a non-adaptive method of order 5.
- % Y = ODE5(ODEFUN,TSPAN,Y0) with TSPAN = [T1, T2, T3, ... TN] integrates
- % the system of differential equations y' = f(t,y) by stepping from T0 to
- % T1 to TN. Function ODEFUN(T,Y) must return f(t,y) in a column vector.
- % The vector Y0 is the initial conditions at T0. Each row in the solution
- % array Y corresponds to a time specified in TSPAN.
- %
- % Y = ODE5(ODEFUN,TSPAN,Y0,P1,P2...) passes the additional parameters
- % P1,P2... to the derivative function as ODEFUN(T,Y,P1,P2...).
- % This is a non-adaptive solver. The step sequence is determined by TSPAN
- % but the derivative function ODEFUN is evaluated multiple times per step.
- % The solver implements the Dormand-Prince method of order 5 in a general
- % framework of explicit Runge-Kutta methods.
- %
- % Example
- % tspan = 0:0.1:20;
- % y = ode5(@vdp1,tspan,[2 0]);
- % plot(tspan,y(:,1));
- % solves the system y' = vdp1(t,y) with a constant step size of 0.1,
- % and plots the first component of the solution.
- if ~isnumeric(tspan)
-
- error('TSPAN should be a vector of integration steps.');
-
- end
- if ~isnumeric(y0)
-
- error('Y0 should be a vector of initial conditions.');
-
- end
- h = diff(tspan);
- if any(sign(h(1))*h <= 0)
-
- error('Entries of TSPAN are not in order.')
-
- end
- try
-
- f0 = feval_r(odefun,tspan(1),y0,varargin{:});
-
- catch
-
- msg = ['Unable to evaluate the ODEFUN at t0,y0. ',lasterr];
-
- error(msg);
-
- end
- y0 = y0(:); % Make a column vector.
- if ~isequal(size(y0),size(f0))
-
- error('Inconsistent sizes of Y0 and f(t0,y0).');
-
- end
- neq = length(y0);
- N = length(tspan);
- Y = zeros(neq,N);
- % Method coefficients -- Butcher's tableau
- %
- % C | A
- % --+---
- % | B
- C = [1/5; 3/10; 4/5; 8/9; 1];
- A = [ 1/5, 0, 0, 0, 0
-
- 3/40, 9/40, 0, 0, 0
- 44/45 -56/15, 32/9, 0, 0
- 19372/6561, -25360/2187, 64448/6561, -212/729, 0
- 9017/3168, -355/33, 46732/5247, 49/176, -5103/18656];
- B = [35/384, 0, 500/1113, 125/192, -2187/6784, 11/84];
- % More convenient storage
- A = A.';
- B = B(:);
- nstages = length(B);
- F = zeros(neq,nstages);
- Y(:,1) = y0;
- for i = 2:N
-
- ti = tspan(i-1);
-
- hi = h(i-1);
-
- yi = Y(:,i-1);
-
- % General explicit Runge-Kutta framework
-
- F(:,1) = feval_r(odefun,ti,yi,varargin{:});
-
- for stage = 2:nstages
-
- tstage = ti + C(stage-1)*hi;
-
- ystage = yi + F(:,1:stage-1)*(hi*A(1:stage-1,stage-1));
-
- F(:,stage) = feval_r(odefun,tstage,ystage,varargin{:});
-
- end
-
- Y(:,i) = yi + F*(hi*B);
-
- end
- Y = Y.';
来自：jlxiaohuo的博客 |
