matlab ode45 函数传自定义参数用法及定步长ode5解算函数

sorry

%程序1

1. arg1 = 2;
   
2. arg2 = 1;
   
3. [T,Y] = ode45('vdp1000',[0 10],[2 0], [], arg1, arg2);
   
4. plot(T,Y(:,1),'-o');

复制代码

%%程序2

1. function dy = vdp1000(t, y, flag, arg1, arg2)
   
2. dy = zeros(2,1);    % a column vector
   
3. dy(1) = y(2);
   
4. dy(2) = arg1*(arg2 - y(1)^2)*y(2) - y(1);

复制代码

%%ode5

1. function Y = ode5(odefun,tspan,y0,varargin)
   
2. %ODE5 Solve differential equations with a non-adaptive method of order 5.
   
3. % Y = ODE5(ODEFUN,TSPAN,Y0) with TSPAN = [T1, T2, T3, ... TN] integrates
   
4. % the system of differential equations y' = f(t,y) by stepping from T0 to
   
5. % T1 to TN. Function ODEFUN(T,Y) must return f(t,y) in a column vector.
   
6. % The vector Y0 is the initial conditions at T0. Each row in the solution
   
7. % array Y corresponds to a time specified in TSPAN.
   
8. %
   
9. % Y = ODE5(ODEFUN,TSPAN,Y0,P1,P2...) passes the additional parameters
   
10. % P1,P2... to the derivative function as ODEFUN(T,Y,P1,P2...).
    
11. % This is a non-adaptive solver. The step sequence is determined by TSPAN
    
12. % but the derivative function ODEFUN is evaluated multiple times per step.
    
13. % The solver implements the Dormand-Prince method of order 5 in a general
    
14. % framework of explicit Runge-Kutta methods.
    
15. %
    
16. % Example
    
17. % tspan = 0:0.1:20;
    
18. % y = ode5(@vdp1,tspan,[2 0]);
    
19. % plot(tspan,y(:,1));
    
20. % solves the system y' = vdp1(t,y) with a constant step size of 0.1,
    
21. % and plots the first component of the solution.
    
22. if ~isnumeric(tspan)
    
23.    
    
24.     error('TSPAN should be a vector of integration steps.');
    
25.    
    
26. end
    
27. if ~isnumeric(y0)
    
28.    
    
29.     error('Y0 should be a vector of initial conditions.');
    
30.    
    
31. end
    
32. h = diff(tspan);
    
33. if any(sign(h(1))*h <= 0)
    
34.    
    
35.     error('Entries of TSPAN are not in order.')
    
36.    
    
37. end
    
38. try
    
39.    
    
40.     f0 = feval_r(odefun,tspan(1),y0,varargin{:});
    
41.    
    
42. catch
    
43.    
    
44.     msg = ['Unable to evaluate the ODEFUN at t0,y0. ',lasterr];
    
45.    
    
46.     error(msg);
    
47.    
    
48. end
    
49. y0 = y0(:); % Make a column vector.
    
50. if ~isequal(size(y0),size(f0))
    
51.    
    
52.     error('Inconsistent sizes of Y0 and f(t0,y0).');
    
53.    
    
54. end
    
55. neq = length(y0);
    
56. N = length(tspan);
    
57. Y = zeros(neq,N);
    
58. % Method coefficients -- Butcher's tableau
    
59. %
    
60. % C | A
    
61. % --+---
    
62. % | B
    
63. C = [1/5; 3/10; 4/5; 8/9; 1];
    
64. A = [ 1/5, 0, 0, 0, 0
    
65.    
    
66. 3/40, 9/40, 0, 0, 0
    
67. 44/45 -56/15, 32/9, 0, 0
    
68. 19372/6561, -25360/2187, 64448/6561, -212/729, 0
    
69. 9017/3168, -355/33, 46732/5247, 49/176, -5103/18656];
    
70. B = [35/384, 0, 500/1113, 125/192, -2187/6784, 11/84];
    
71. % More convenient storage
    
72. A = A.';
    
73. B = B(:);
    
74. nstages = length(B);
    
75. F = zeros(neq,nstages);
    
76. Y(:,1) = y0;
    
77. for i = 2:N
    
78.    
    
79.     ti = tspan(i-1);
    
80.    
    
81.     hi = h(i-1);
    
82.    
    
83.     yi = Y(:,i-1);
    
84.    
    
85.     % General explicit Runge-Kutta framework
    
86.    
    
87.     F(:,1) = feval_r(odefun,ti,yi,varargin{:});
    
88.    
    
89.     for stage = 2:nstages
    
90.       
    
91.         tstage = ti + C(stage-1)*hi;
    
92.       
    
93.         ystage = yi + F(:,1:stage-1)*(hi*A(1:stage-1,stage-1));
    
94.       
    
95.         F(:,stage) = feval_r(odefun,tstage,ystage,varargin{:});
    
96.       
    
97.     end
    
98.    
    
99.     Y(:,i) = yi + F*(hi*B);
    
100.    
     
101. end
     
102. Y = Y.';

复制代码

来自：jlxiaohuo的博客