Implicit Runge Kutta Order 4.

4 次查看(过去 30 天)
Tom Craven
Tom Craven 2015-6-1
评论: Jan 2015-6-1
I am trying to write an Solver for the implicit runge kutta method of order 4. at the moment my function is not working im not sure if it because of the unappropriate initial value for solving the system or I made some mistakes when making the function.
Any suggestions for help or sources of helpful information would be appreciated.
function [ tout,yout ] = IRK4Solver( f,t,y0 ) % INPUT: f(t,y) is an anonymous function that defines
% the right-hand side of the ODE ydot = f(t,y)
% t =[t0 t1 ... tfinal] is a vector of grid points
% with length N
% y0=[a b c] is a column vector that contain the
% initial values x(0)=a, y(0)=b, z(0)=c.
% OUTPUT:tout is a column vector of grid points.
% yout is an 3 x N matrix containing the solution
% at different grid points.
N=numel(t);
yout=zeros(3,N);
yout(:,1)=y0;
for n=2:N
h=t(n)-t(n-1); f=@(x1,x2)([yout(:,n-1)+0.25*h.*f((t(n-1)+(0.5-sqrt(3)/6)*h),x1)+(0.25-sqrt(3)/6)*h.*f((t(n-1)+(0.5+sqrt(3)/6)*h),x2)-x1;yout(:,n-1)+(0.25+sqrt(3)/6)*h.*f((t(n-1)+(0.5-sqrt(3)/6)*h),x1)+0.25*h.*f((t(n-1)+(0.5+sqrt(3)/6)*h),x2)-x2]);
fsolve(f,[1,1,1;1,1,1]);
yout(:,n)=yout(:,n-1)+0.5*h.*f((t(n-1)+(0.5-sqrt(3)/6)*h),x1)+0.5*h.*f((t(n-1)+(0.5+sqrt(3)/6)*h),x2);
end
tout=t(n);
end
  1 个评论
Jan
Jan 2015-6-1
Please format the code properly, when you want others to read it. Remove the empty lines between each lines of the code. Use the "{} Code" button. Then explain "does not work" with any Details.

请先登录,再进行评论。

回答(0 个)

标签

尚未输入任何标签。

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by