want to solve (d^2 y)/(dx^2 )+dy/dx-6y=0 using 4th order Runge-Kutta method with y(0) = 3 and y’(0) = 1
조회 수: 33 (최근 30일)
이전 댓글 표시
My code is :
clc
clear all
h = 0.5;
x = 1:h:5;
y = zeros(2,length(x)); % y vector declaration
y(1,1) = 3; % y value
y(2,1) = 1; % y' value
f = @(x) (dy^2)/(dx^2 )+dy/dx-6*y;
for i=1:(length(x)-1) % calculating y values
z1 = f( x(i) , y(:,i));
z2 = f( x(i)+0.5*h, y(:,i)+0.5*h*z1);
z3 = f( x(i)+0.5*h, y(:,i)+0.5*h*z2);
z4 = f( x(i)+h, y(:,i)+h*z3);
y(:,i+1) = y(:,i) + (1/6)*(z1 + 2*z2 + 2*z3 + z4)*h;
end
%plotting results
Y = x.^2 - x - 6;
figure
hold on
plot(x,Y,'b'); % analytic solution plot--- blue
plot(x,y(1,:),'r'); % runga kutta solution --- Red
grid on
legend('Analytical solution','RK4');
xlabel('time (T)')
ylabel('y')
title('Runga gutta 4th order');
But My output is :
Too many input arguments.
Error in Untitled3 (line 11)
z1 = f( x(i) , y(:,i));
How can i fix this problem sir?
댓글 수: 0
채택된 답변
James Tursa
2021년 2월 7일
편집: James Tursa
2021년 2월 7일
Change your function handle from this
f = @(x) (dy^2)/(dx^2 )+dy/dx-6*y;
to this
f = @(x,y) [y(2);-y(2)+6*y(1)];
The last part of that function handle is simply solving this
(d^2 y)/(dx^2 )+dy/dx-6y=0
for the 2nd derivative and then plugging in y(2) and y(1) for the other parts.
추가 답변 (1개)
Riccardo Bonomelli
2021년 2월 6일
Hi, the problem seems to be in the definition of the function f, you defined f as a single variable function or f(x) while, according to your code, f should be a function of two variables, like f(x,y). Furthermore, the dy and dx inside the function f in your definition must be defined in some way, say dx= something and dy= something. I guess you are trying to express a derivative, but you must define it in some way.
댓글 수: 3
Riccardo Bonomelli
2021년 2월 6일
The definition of the derivatives depends on the equation you are trying to solve, in your case you are solving a second order differential equation which requires a system of two first order differential equations to be solved using 4th order Runge Kutta. I'm no expert on the subject, so I am not comfortable explaining the details of the implementation.
참고 항목
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!