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
19 次查看(过去 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!