Driven Damped Harmonic Oscillator

General Instruction

All italicized text must be deleted before submitting your lab report.

the symbol indicates code that needs to be completed.

All major headings must also be filled in.

Introduction

In this section introduce the subject matter and provide the necessary math and physics relevant to this experiment.

Experimental strategy

In this section explain your methods of determining relevant values as well as the Euler's loop you created. Also, explain briefly why each variable is necessary.

Mathematical Modeling

Mass-Spring oscillator differential equation:

alternatively we can use instead of

Here we have to be very careful to realize that the friction is always in the opposite direction compared to the velocity. If the velocity is positive, the air resistance is in the negative direction and vice versa. We can accomplish that by using the function "sign(v)" which returns "-1" if the velocity is in the negative direction and "+1" if the velocity is in the positive direction.Solving the differential equation numerically using a modified Euler's method.

Our equations approximating the velocities and positions are:

clear

Now that we have the necessary algorithms, we can program the model.

Numerical Model (MATLAB code)

First we declare our known constants:

Some explanations:

• P is the amplitude of the "driver position". It can be estimated by measuring how much the second cart (the which is attached directly to the rotating driver arm) moves back and forth.
• w is the radial frequency of the oscillator and thus the radial frequency with which the second cart (the which is attached directly to the rotating driver arm) moves back and forth.

% code example

m=....; % mass of cart in kg

k_1=....; % Spring constant 1 in N/m

k_2=....; % Spring constant 2 in N/m

b=....; % magnetic damping coefficient in Ns/m

f=....; % rolling resistance & axle friction in N

k = k_1+k_2; % calculated total k in N/m

P=....; % Driver position amplitude in m

w=....; % Radial frequency of driver position in rad/s

F=0.085; % Amplitude of driving force in N

Computing some frequencies for convencience

% code example

fd=w/2/pi() % Frequency of driving force in Hz

wo=sqrt(k/m) % undamped resonance frequency (in rad/s)

fo=wo/2/pi() % undamped resonance frequency (in Hz)

w_damped=sqrt(k/m-(b/2*m)^2) % damped resonance frequency in rad/s

f_damped=w_damped/2/pi() % damped resonance frequency in Hz

Set time stepsize and number for Euler's method calculation

% code example

Total_Time=20; % total time in seconds

dt=.00002; % Set step size (seconds)

N=round(Total_Time/dt); % Set the amount of steps to nearest integer

Preallocate arrays (for speed):

% code example

x=zeros(N,1);

v=zeros(N,1);

t=zeros(N,1);

PD=zeros(N,1); % array for driving force values

Set initial conditions:

% code example

x(1)=0; % Set initial condition for position in m

v(1)=0; % Set initial condition for velocity in m/s

t(1)=0; % Set initial time (0).

The Euler's method loops:

Some explanations:

• w_a = .... is used to allow you to change the frequency of the driving force within the for-loop in once convenient place.
• You can also use w_a=w_damped or w_a=w_damped/2 etc. to run the simulation at particular frequencies related to the system properties.

Remember (from above):

% code example

w_a =....; % set active w for simulation

for n=1:N % Set the number of loops

v(n+1)=....; % calculate next velocity

x(n+1)=....; % calculate next position

t(n+1)=....; % calculate next time value

end

PD=P*sin(w_a*t); % Position of driver as function of time aray

Plot the results

% code example

plot(t,x)

title('Position and driving force')

xlabel('time in s')

ylabel('Position in m')

grid 'on'

grid 'minor'

hold on

plot(t,PD)

legend('position oscillator','position driver')

hold off

[David Hill]

What is your question?