Matlab Solving Differential Equation using Runge Kutta

clc;
h = 50;
x = 0: h: 0.025;
w = 3553.6312;
M = 110.361;
f = 0.02;
A = 24;
Cd = 0.32;
P = 1925.637;
F = 16.80444124;
x(1) = 700; % initial condition
F_yx = @(y, x)(P - f.*w - 0.00118.*Cd * A * y ^ 2) / (M); % function
for i = 1: (length(x) - 1) % calculation loop
k_1 = F_xy(x(i), xy(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
syms y(x)
ode = diff(y, x) - F * y == 0;
ySol(x) = dsolve(ode);
fimplicit(y(x))

답변 (1개)

James Tursa
James Tursa 2019년 10월 25일
In this line:
k_1 = F_xy(x(i), xy(i));
You've got xy(i) as an argument, but there is no xy variable. I think you meant y(i) here.
In this line:
F_yx = @(y, x)(P - f.*w - 0.00118.*Cd * A * y ^ 2) / (M); % function
You've got the order of the inputs as (y,x), but in your downstream code you have the order of the inputs as (x,y). Change your function to use the (x,y) order.
This initial condition:
x(1) = 700; % initial condition
I suspect should apply to y(1), not x(1). You already have the x vector defined above.

카테고리

도움말 센터File Exchange에서 Loops and Conditional Statements에 대해 자세히 알아보기

태그

질문:

2019년 10월 25일

답변:

2019년 10월 25일

Community Treasure Hunt

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

Start Hunting!

Translated by