Simulating a free-fall project with ode45

12 次查看（过去 30 天）
Lukas Menzel
Lukas Menzel 2020-12-23
回答： Mischa Kim 2020-12-23
Hallo everybody, now that I've the basic solution for my ode45-method, I've another three questions.
  1. How do I add differntials of higher or lower orders in one of the other differntial-equation, like in my code. 'rho' is depended on the height h, but that is the first order of the ode, so how do I implement, the first order in the second order?
function dv = Beschleunigung(~, hv)
g = 9.81;
R = 8.31446261;
M = 0.0289586;
m = 120;
c_w = 0.28;
A = 2.7;
p_0 = 101325;
T_EO = 293.15;
rho_0 = (p_0 * M) / (R * T_EO);
v = hv(2); % velocity, positive downwards
dv = [
-v; % dh/dt
g-1/(2*m)*c_w*A*rho_0 * exp(-g*M*h(1)/R*T)*v^2]; % dv/dt
end
2. Since my project needs to show all forms of motion-equations (acceleration, speed, height), how do I plot the acceleration, with is the derivate of the ode45-method
3. And since it is required to name a timespan, and especelly the end, is there a possibility to let the method end, even when you do not exactly know the time needed, to reach a height of 0, in my specific case?

请先登录，再回答此问题。

回答（1 个）

Mischa Kim
Mischa Kim 2020-12-23
Hallo Lukas,
  1. The height is the first component of the state vector. I'll just call the state vector X is in the code below. So you can just back out h = X(1).
  2. Acceleration...you have the expression in the dX vector, the second component. To compute the acceleration as a function of time simply compute this expression in your main function, the function that is calling Beschleunigung().
  3. The best and most elegant way to do this is to use events functions. Run in the command window the built-in demo ballode and check out the code.
function dX = Beschleunigung(~, X) % X is the state vector X = [h; v]
g = 9.81;
R = 8.31446261;
M = 0.0289586;
m = 120;
c_w = 0.28;
A = 2.7;
p_0 = 101325;
T_EO = 293.15;
rho_0 = (p_0 * M) / (R * T_EO);
h = X(1); % position or height
v = X(2); % velocity, positive downwards
dX = [-v; % dh/dt
g-1/(2*m)*c_w*A*rho_0 * exp(-g*M*h/R*T)*v^2]; % dv/dt; T not yet defined?
end

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by