Getting linear results from ODE45 instead of exponential

7 次查看(过去 30 天)
I aam very new to Matlab and have this project to figure out the velocity vs time and displacement vs time of a sphere falling from rest. With my code I am able to get the displacement to work which is r(1) but my velocity vs time plot is linear when it should be exponential. I've tried retyping the whole equation for the code and inputting the squared term in all the ways I could figure out with no success.
function Projecto
clear all
tr=[0 15]; %seconds
initv=[0 0]; %starting point
[t,y]=ode45(@sphere, tr, initv);
plot(t,y(:,1))
ylabel('x (m)')
xlabel('time(s)')
figure
plot(t,y(:,2))
ylabel('velocity (m/s)')
xlabel('time(s)')
function r=sphere(t,y)
mass=167.6;
g=9.81;
rho=1.21; %kilogram per meter cubed
d=.01; %meters
r=zeros(2,1);
r(1)=y(2);
r(2)=(mass*g-rho*g*(pi/6)*d^3-0.5*rho*y(2)*abs(y(2))*(pi/4)*d^2*0.5)/(mass+(pi/12)*d^3*rho
  2 个评论
Ameer Hamza
Ameer Hamza 2020-9-30
Why do you expect it to be exponential? For a sphere falling under constant gravity, the velocity is supposed to increase linearly.
Connor Jack
Connor Jack 2020-9-30
I have to consider drag in the equation so I have a velocity squared term. In the equation r(2) the y(2) is my velocity variable, this is the equation that I started with:
m(dv/dt)=mg-rho*g*(pi/6 * d^3)-0.5rho*v^2*(pi/4 d^2)Cd-(pi/12)d^3*rho*(dv/dt)
rho*g*(pi/6 * d^3) - being the buoyant force
0.5rho*v^2*(pi/4 d^2)Cd - being the drag force
(pi/12)d^3*rho*(dv/dt) - being the force on accelerating body

请先登录,再进行评论。

回答(1 个)

Bjorn Gustavsson
Bjorn Gustavsson 2020-9-30
When you look at the result you have to start looking at how big the different components of the forces are along the trajectory. So write a function where you can calculate the different forces - that way you can judge how big the drag is relative to your force of gravity.
One more thing you can do is to extend the time of your fall - I increased it by a factor of 10 and then you start to see a reduction of the acceleration.
One further thing to try is to modify the mass and density (why are you using both?) to give you a lighter object that should have a lower terminal velocity, when I reduced that by a factor of 10 the terminal velocity can be clearly seen within 150 s.
HTH

类别

Help CenterFile 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