Horn Loudspeaker Simulation part 3: Multiple segments and more T-matrices

In part 2 we combined the horn simulation from part 1 with a simple driver model with front and rear chambers. The method used for that relied on more or less manually finding impedances at certain points in the circuit, and converting back and forth between the different domains when required. 

A more efficient method is to utilize matrix algebra in all stages. In many cases it will be possible to find a composite matrix that describe the system from input to output. And this approach also makes it easy to find the total T-matrix of a horn consisting of several segments. 

Matlab and Octave are great for this, since they have all the matrix functions built-in. Except, since we have a 2D matrix per frequency, we must either do a loop over frequencies, which makes the code a bit messy, or we can overload the matrix functions so that they can accept 3D matrices and do the loop automatically. That simplifies the code considerably. I have added the code for this to the repo. (Unfortunately, function overloading isn't as easy in Octave as in Matlab, and seems a bit clunky. So instead I added separate functions for the Octave case, and separate demo files. If anyone knows a better and more elegant way, please update the code and send a pull request on GitHub.)

Multiple Segments

With multiple segments, we can easily form a composite horn matrix by multiplying together the matrices from one end to the other. Using the 12-type matrices, we would start from the throat end and multiply the next matrix on the right hand side until we reach the mouth:

MatMult12

We can demonstrate this by simulating the throat impedance of a two-segment horn, with a conical first segment and an exponential second segment. 

The T-matrix for a conical horn is defined by

ConicalHornEqs

where

con x1x2

Again we can implement this as a function:

function [a,b,c,d] = conicalHornMatrix(k,Zrc,S1,S2,L)
kL = k*L;
x1 = L/(sqrt(S2/S1)-1);
x2 = x1+L;
kx1 = k*x1;
kx2 = k*x2;
sinkL = sin(kL);
coskL = cos(kL);
a = x2/x1 * coskL - sinkL./kx1;
b = Zrc*1i/S2*x2/x1*sinkL;
c = 1i/Zrc* S1*((x2/x1 + 1./(kx1.*kx1)).*sinkL - coskL.*L./(kx1*x1));
d = x1/x2.*coskL + sinkL./kx2;

 

The function in the repo has a few more bells and whistles. First, if S1=S2, we get a division by zero when calculating x1. So in this case we can either call the function for an exponential horn, which handles this directly, or the function for a straight duct. I've chosen the latter, since it's faster (fewer calculations), and this function can also be used for calculating front and rear chambers.

Second, I've added the option to use type-21 matrices. These are very similar to the type-12 matrices (actually, they are the inverse of each other), and sometimes the analysis can be a bit simpler using one over the other.

What is left now is to calculate the composite matrix:

S1 = 80e-4;
S2 = 350e-4;
L12 = 60e-2;
S3 = 5000e-4;
L23 = 75e-2;
Zrc = rho*c;
% Calculate submatrices
Me = expoHornMatrix(k,Zrc,S2,S3,L23);
Mc = conicalHornMatrix(k,Zrc,S1,S2,L12);
% Calculate composite horn matrix
Mh = Mc*Me;

 

Simple as that. Now we can calculate the radiation impedance and throat impedance as in part 1, and we get the following results:

conExpHorn

The Hornresp results are overlaid with dashed lines. 

Horn Speaker

We can now apply the T-matrix method to the full horn speaker system. First we need the driver matrix. It consists of multiplying together the matrices of the electrical, mechanical and transduction parts:

DriverMatrix

where

DriverTe

DriverBl

DriverTm

DriverSd

Additionally, we need to add the matrices for the front and rear chambers. Looking at our equivalent diagram, repeated below, the rear chamber is a series impedance, and the front chamber is a parallel impedance.

HornSpeakerEquivalentAcoustic

 The matrix for the series impedance is similar to the matrices for TE and TM above. For a parallel impedance,

ParallelZ

We can now set up the total system matrix as

SystemTmatrix

Given that we have the radiation impedance, we can find the electrical input impedance directly from this matrix, using

ThroatZ12

How do we find the rest of the quantities? Previously we calculated the Thevening equivalent of the driver, which we also could do from the driver matrix, but here's another approach: With the electrical impedance, we can find the input current, given the input voltage, and then we can use the matrices to find the output pressure and volume velocity:

InputCurrent

PmUmMatrix

The output power is calculated as before, or we can use the mouth pressure and volume velocity directly. 

The diaphragm displacement is calculated from the driver volume velocity, which is found by using the above process on the driver matrix only. In code:

 Iin = eg./Ze.';
egv = ones(1,length(freq)) * eg;
inputs = [egv; Iin];
outputs = inv(T) * inputs;
driverOut = inv(TD) * inputs;
pm = outputs(1,:);
Um = outputs(2,:);
% Power into the load
Pa = real(pm.*conj(Um));
% Convert power to SPL
I = Pa / (2*pi);
prad = sqrt(I*rho*c);
% Diaphragm displacement
UaL = driverOut(2,:);
x = UaL./w / Sd * sqrt(2);

 

The results:

blog3spl

blog3displ

blog3ze

Again, we have a good match with the hornresp results.

The code for this blog post, including the utility functions and overloaded matrix operations, can be found here.