從MathWorks上找了些Runge-Kutta算法及應(yīng)用的源碼，轉(zhuǎn)到這里，希望對要解常微分方程（組）的蟲友有借鑒。

1、Runge Kutta 4th order ode（用四階RK法解ODE）

https://www.mathworks.com/matlabcentral/fileexchange/29851-runge-kutta-4th-order-ode

by Judah S
29 Dec 2010
Code covered by the BSD License

CODE:

%It calculates ODE using Runge-Kutta 4th order method
%Author Ido Schwartz

clc; %Clears the screen
clear all;

h=1.5; %step size
x = 0:h:3; % Calculates upto y(3)
y = zeros(1,length(x));
y(1) = 5;
F_xy = @(t,r) 3.*exp(-t)-0.4*r; %change the function as you desire

for i=1length(x)-1)   %calculation loop
    k_1 = F_xy(x(i),y(i));
    k_2 = F_xy(x(i)+0.5*h,y(i)+0.5*h*k_1);
    k_3 = F_xy((x(i)+0.5*h),(y(i)+0.5*h*k_2));
    k_4 = F_xy((x(i)+h),(y(i)+k_3*h));

    y(i+1) = y(i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*h;  %main equation
end

2.  sixth order Runge-Kutta ODE（It calculates 6th order RK ODE using Fehlberg formulation）

https://www.mathworks.com/matlabcentral/fileexchange/29931-6th-order-runge-kutta-ode

by Judah S
06 Jan 2011
Code covered by the BSD License

CODE:

clc
clf
clear all

% % Author
%  Judah S
% Sixth Order RK-Fehlberg Equation
% Function y' + 2xy = 0; initial condition y(0)=1; where 0 <= x <= 10
% You can easily change the function and rest condition according to your
% requirment.
   fcnstr='-2*x*y' ;
   f=inline(fcnstr) ;

   x0=0 ;
   y0=1 ;

% xf, x location at where you wish to see the solution to the ODE

   xf=10 ;

% n, number of steps to take

   n=5 ;

format long g
h=(xf-x0)/n ;
xa(1)=x0 ;
ya(1)=y0 ;

for i=1:n
    % Adding Step Size
  xa(i+1)=xa(i)+h ;

% Calculating k1, k2, k3, k4, K5 and K6
  k1 = f(xa(i),ya(i)) ;
  k2 = f(xa(i)+(1/4)*h,ya(i)+(1/4)*k1*h) ;
  k3 = f(xa(i)+(3/8)*h,ya(i)+(3/32)*h*(k1+3*k2)) ;
  k4 = f(xa(i)+(12/13)*h,ya(i)+(12/2197)*h*(161*k1-600*k2+608*k3)) ;
  k5 = f(xa(i)+h,ya(i)+(1/4104)*h*(8341*k1-32832*k2+29440*k3-845*k4)) ;
  k6 = f(xa(i)+(0.5)*h,ya(i)+h*(-(8/27)*k1+2*k2-(3544/2565)*k3+(1859/4104)*k4-(11/40)*k5)) ;
% Using 6th Order Runge-Kutta formula
  ya(i+1)=ya(i)+1/5*((16/27)*k1+(6656/2565)*k3+(28561/11286)*k4-(9/10)*k5+(2/11)*k6)*h ;

  end


disp(sprintf('\n\n********************************Results**********************************'))

% The following finds what is called the 'Exact' solution
xspan = [x0 xf];
[x,y]=ode45(f,xspan,y0);
[yfi dummy]=size(y);
yf=y(yfi);

% Plotting the Exact and Approximate solution of the ODE.
hold on
xlabel('x');ylabel('y');
title('Exact and Approximate Solution of the ODE by the 4th Order Runge-Kutta Method');
plot(x,y,'--','LineWidth',2,'Color',[0 0 1]);            
plot(xa,ya,'-','LineWidth',2,'Color',[0 1 0]);
legend('Exact','Approximation');


disp(sprintf('\n   Approximate = %g',ya(n+1)))
disp(sprintf('   Exact       = %g',yf))
disp(sprintf('\n   True Error = Exact - Approximate'))
disp(sprintf('              = %g - %g',yf,ya(n+1)))
disp(sprintf('              = %g',yf-ya(n+1)))
disp(sprintf('\n   Absolute Relative True Error Percentage'))
disp(sprintf('              = | ( Exact - Approximate ) / Exact | * 100'))
disp(sprintf('              = | %g / %g | * 100',yf-ya(n+1),yf))
disp(sprintf('              = %g',abs( (yf-ya(n+1))/yf )*100))

3. RKN1210 - A 12th/10th order Runge-Kutta-Nystrom integrator（Integrator for second-order ODE's with very stringent error tolerances）

https://www.mathworks.com/matlabcentral/fileexchange/25291-rkn1210-a-12th10th-order-runge-kutta-nystrom-integrator

by Rody Oldenhuis
11 Sep 2009 (Updated 28 Oct 2010)
Code covered by the BSD License
源碼及示例見上面鏈接。

4. Accurate Fast Marching（Multistencils second order Fast Marching 2D and 3D including rk4 shortest path and skeletonize）

https://www.mathworks.com/matlabcentral/fileexchange/24531-accurate-fast-marching

by Dirk-Jan Kroon
23 Jun 2009 (Updated 21 Dec 2010)
Code covered by the BSD License
源碼比較多，感興趣就請從上面的鏈接下吧。

5. double pendulum（the function is GUI that simulate double pendulum，the function is GUI that simulate double pendulum by using 4th-order Runge-Kutta algorithm for the differential equations）

https://www.mathworks.com/matlabcentral/fileexchange/28490-double-pendulum

by Moshe Lindner
16 Aug 2010
Code covered by the BSD License
源碼見以上連接。

6. Weighted Kinematic Control of PPRR Manipulator（PPRR機(jī)械手運(yùn)動(dòng)控制）

https://www.mathworks.com/matlabcentral/fileexchange/23910-weighted-kinematic-control-of-pprr-manipulator

by Hrishi Shah
26 Apr 2009 (Updated 16 Feb 2010)
Code covered by the BSD License
源碼見以上鏈接。

7. Program to Solve initial value problems by various methods（求解初值問題）

https://www.mathworks.com/matlabcentral/fileexchange/24614-program-to-solve-initial-value-problems-by-various-methods

by Daniel Klawitter
01 Jul 2009 (Updated 01 Jul 2009)
Code covered by the BSD License
很實(shí)用的初值問題求解代碼，下載見以上鏈接。

8. Mackey-Glass time series generator（generates a Mackey-Glass time series using the 4-th order Runge-Kutta method）

https://www.mathworks.com/matlabcentral/fileexchange/24390-mackey-glass-time-series-generator

by Marco Cococcioni
09 Jun 2009
Code covered by the BSD License

9. CHAO-AND-HEART-DATABASE（Each volume is composed of figures of simulations with matlab R2007a, with runge kutta 4,with comments.）

https://www.mathworks.com/matlabcentral/fileexchange/24272-chao-and-heart-database

by Dit Papa Lamine ndao
27 May 2009
Code covered by the BSD License

10. Cardiac Chao with modele pacemaker VI1（The model is electrical stimulated with current Io=constant or Io=Asinwt.The model is integrated with runge kutta and R2007a.）

https://www.mathworks.com/matlabcentral/fileexchange/21531-cardiac-chao-with-modele-pacemaker-vi1

by Dit Papa Lamine ndao
24 Sep 2008 (Updated 24 Sep 2008)
No BSD License

11. Runge Kutta Ballistic Propagation（一個(gè)彈道的計(jì)算程序）

https://www.mathworks.com/matlabcentral/fileexchange/16626-runge-kutta-ballistic-propagation

by Luke
27 Sep 2007 (Updated 28 Sep 2007)
No SD License

12. M-File Comput N°1 Of Model Pacemaker VI1 Of All The Cardiac Cells With R2007a（computing with runge kutta 4 model pacemaker VI1 for all the cardiac cells）

https://www.mathworks.com/matlabcentral/fileexchange/16039-m-file-comput-n%C2%B01-of-model-pacemaker-vi1-of-all-the-cardiac-cells-with-r2007a

by Dit Papa Lamine NDao
21 Aug 2007 (Updated 22 Aug 2007)
No BSD License

13. Designing Algorith and Programation Of Cardiac Software Guide with Mat Lab（algorith and programation of a guide computing the symbolic model pacemake）

https://www.mathworks.com/matlabcentral/fileexchange/15840-designing-algorith-and-programation-of-cardiac-software-guide-with-mat-lab

by Dit Papa Lamine NDao
06 Aug 2007 (Updated 06 Aug 2007)
No BSD License

14. M-File Application OF Model Pcemaker VI1 of Cardiac Cells With Runge Kutta 4

https://www.mathworks.com/matlabcentral/fileexchange/15795-m-file-application-of-model-pcemaker-vi1-of-cardiac-cells-with-runge-kutta-4

by Dit Papa Lamine NDao
02 Aug 2007 (Updated 02 Aug 2007)
No BSD License

15. ROCK4（explicit solver for stiff equations）

https://www.mathworks.com/matlabcentral/fileexchange/12129-rock4

by Gerd Steinebach
01 Sep 2006 (Updated 01 Sep 2006)
No BSD License

16. Symplectic Integrators（Symplectic integrators for separable and nonseparable Hamiltonians）

https://www.mathworks.com/matlabcentral/fileexchange/7686-symplectic-integrators

by Francisco Beron-Vera
18 May 2005 (Updated 21 Jul 2006)
Code covered by the BSD License

17. TriStream（Compute streamlines on triangular mesh）

https://www.mathworks.com/matlabcentral/fileexchange/11278-tristream

by Matthew Wolinsky
02 Jun 2006 (Updated 06 Jun 2006)
No BSD License

18. equationGUI（A simple GUI for ODE (van der pol-prey predator - logistic model) and PDE (vibrating string)）

https://www.mathworks.com/matlabcentral/fileexchange/10784-equationgui

by giancarlo zaccone
18 Apr 2006 (Updated 20 Apr 2006)
Code covered by the BSD License

19. fieldlines3D（fieldlines3D(X,Y,Z,U,V,W,x0,y0,z0) plots a field line with the initial point P(x0,y0,z0).The matrices X,Y,Z,U,V,and W must all be the same size and contain corresponding position and field components.The function employs the 4th order Runge-Kutta method for solving numerically differential equation y'=f(x,y).）

https://www.mathworks.com/matlabcentral/fileexchange/7768-fieldlines3d

by Avni Pllana
02 Jun 2005 (Updated 06 Jun 2005)
Code covered by the BSD License

20. fieldlines（Plots a field line with the initial point P(x0,y0) for a 2D-field.）

https://www.mathworks.com/matlabcentral/fileexchange/5710-fieldlines

by Avni Pllana
14 Aug 2004 (Updated 18 Aug 2004)
Code covered by the BSD License

21. Dynamics of Some Classical System Models（Animated response of some classical systems）

https://www.mathworks.com/matlabcentral/fileexchange/6558-dynamics-of-some-classical-system-models

by Howard Wilson
15 Dec 2004 (Updated 20 Dec 2004)
No BSD License

22. rkn86（Solves numerically a system of ODEs of the form d^2 y/ dx^2=f(x,y).）

https://www.mathworks.com/matlabcentral/fileexchange/3153-rkn86

by Charalampos Tsitouras
17 Mar 2003 (Updated 20 Feb 2004)
No BSD License

23. ode87 Integrator（integrates ode system with high accuracy）

https://www.mathworks.com/matlabcentral/fileexchange/3616-ode87-integrator

by Vasiliy Govorukhin
18 Jun 2003 (Updated 19 Jun 2003)
Code covered by the BSD License

24. ode86（integrates a system of ordinary differential equations at stringent tolerances）

https://www.mathworks.com/matlabcentral/fileexchange/2911-ode86

by Charalampos Tsitouras
06 Jan 2003 (Updated 07 Jan 2003)
No BSD License 返回小木蟲查看更多