MATLAB: Modelling behaviour of simple sound source

harmonic oscillationMATLABsimple sound source

I am trying to plot the behaviour of a simple sound source (a harmonically oscillating piston) as the driving frequency passes through the resonant frequency of the piston.
However, the graph I get doesn't look like what I am expecting to see: i.e. it doesn't look like an increase in sound pressure level (dB) up to resonance, and a reduction thereafter. What I get returned in matlab looks more like a graph of the form ln x .
I think I have identified where the problem lies – in the expression for the velocity of the piston (which is then fed into another equation to obtain the air pressure…)
Given that omega (the angular driving frequency) takes a range of values, I would expect the velocity also to display an array of different values, each depending on the particular driving frequency.
However, I appear to be getting just one value for velocity (ut) of the piston, and can't see why. Perhaps it's something to do with the complex numbers in the expression…
Anyway, in the hope that someone can point out the error, this is the relevant code (the expression for velocity has been squeezed down the bottom, along with a couple of other variables):
%————————————————————————– % Parameters %————————————————————————–
mt=0.07; % mass of radiator
At=0.0430; % area of radiator
Rt=3.0; % damping coefficient, t
Kt=68000; % equivalent t stiffness
c2=340^2; % speed of sound squared rho=1.205; % density of air V = 0.015; %Vol of box
pref=50;
kappa=c2*rho/V;
%————————————————————————– % Uncoupled natural frequency %————————————————————————–
omt=sqrt( (Kt + kappa*At^2)/mt);
%————————————————————————– % % frequency, angular frequency % %————————————————————————–
f=(172:0.1:180);
om=2*pi*f; %————————————————————————– % Velocity %————————————————————————– >> F=0.4; >> gamt=Rt/mt; >> ut=F*i*om/(mt*(omt^2-om.^2)+i*gamt*om);

Best Answer

Will, use
ut = F*i*om./(mt*(omt^2-om.^2) + i*gamt*om);
Notice the dot in front of the forward slash.