MATLAB: Bike model with and without air resistance using matlab

bicycle model; air resistance

Hi I'm trying to follow along in a computational physics book (Computational Physics; Nicholas Giordano 2nd edition) that attempts to model a bike with and without air resistance. I did everything that the book did but I got a completely different plot at the end.This is my code for the bike with air resistance.

function bike(v0,dt,tf)

t = 0:dt:tf;

v(1) = v0;

P = 400;

m = 70;

p = 1.225

A = 0.33

for i = 1:length(t)-1

    v(i+1) = v(i) + (P/m*v(i)-(p*A*v(i)^2/2*m))*dt;
    t(i+1) = t(i) + dt;
end

plot(t,v,'b');

title('Comparison of Euler approximation to actual solution')

xlabel('time')

ylabel('v')

disp(v(end));

P stands for power, m is mass, p is density of air, and A is the frontal area of the rider.

This is my code for the bike without air resistance.

function bike(v0,dt,tf)

t = 0:dt:tf;

v(1) = v0;

P = 400;

m = 70;

for i = 1:length(t)-1

    v(i+1) = v(i) + (P/m*v(i))*dt;
    t(i+1) = t(i) + dt;
end

plot(t,v,'b');

title('Comparison of Euler approximation to actual solution')

xlabel('time')

ylabel('v')

disp(v(end));

If someone could please help me out, I would really appreciate it.

Best Answer

The dimensions in your formulas do not look right. I think it is a matter of operator precedence. 1/2*3 is not the same as 1/(2*3).

Also, note that you build the time vector twice: once at the beginning

 t = 0:dt:tf

and then in the loop

t(i+1) = t(i) + dt;

However this should be no problem.

Related Solutions

MATLAB: Variable Usage in Function File

1) Return vel as an output of your function. E.g.,

function [vel,vend] = velocity1(dt, ti,tf,vi)

2) Have your driver code and function code in separate files. E.g.,

% A file called velocity1_driver.m (or some other name of your choosing)
dt = 0.5;
ti = 0;
tf = 12;
vi = 0;
[vel,vend] = velocity1(dt, ti, tf, vi);
  :
 etc

And then a function file

% A file called velocity1.m
function [vel, vend] = velocity1(dt,ti,tf,vi)
  :
 etc

MATLAB: Keep on getting ‘Error using plot, Vectors must be the same length’. How to plot it correctly

When the while loop breaks, you have defined T to have length of i_max, and t has length i_max+1. This is because of the line

t(i+1) = t(i) + t_step

which causes the vector t to be one element longer than T. To be able to plot correctly, you could either drop the last element in t (since it doesn't appear to affect the last element of T) or compute a final element of T. It looks like you already thought of my first idea, so if you remove the lines

t(end) = [];
h(end) = [];

your code should work.

Edit I made some errors in my previous answer, this is corrected now.

Related Question