/* [wxMaxima batch file version 1] [ DO NOT EDIT BY HAND! ]*/
/* [ Created with wxMaxima version 0.8.2 ] */

/* [wxMaxima: title   start ]
speaker.wxm Speaker Design Program
   [wxMaxima: title   end   ] */

/* [wxMaxima: section start ]
Introduction
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
Written by Francis Deck, May 19, 2009

Distributed under GNU Public License. Details at: 

    http://www.gnu.org/copyleft/gpl.html

The purpose of this workbook is to provide a speaker design "program"
that is both multi-platform and open-source. Rather than mess with yet
another programming language and user interface design tool, I decided
to use a computer algebra system (CAS) with powerful computation and
graphing functions built-in. In addition, it makes the nuts and bolts
visible to everybody, and shows how simple it really is.

There is more documentation at the bottom of the workbook.
   [wxMaxima: comment end   ] */

/* [wxMaxima: section start ]
Installation
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
If you are reading this for the first time, you are probably reading a
PDF printout of the workbook. You can look at it, but it won't do
anything. Here are the installation steps for the actual workbook:

1. Download WxMaxima for free from SourceForge, and install it

2. Download the workbook, speaker.wxm, from my site

3. Run WxMaxima, and open the workbook
   [wxMaxima: comment end   ] */

/* [wxMaxima: section start ]
Brief Instructions
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
1. Modify the user input section to reflect your design choices

2. Type Ctrl-R to re-evaluate the entire workbook

3. Scroll down, to find your results

4. Save the workbook, preferably with a new name

5. You can run more than one copy of the workbook at a time, to
   compare different designs.
   [wxMaxima: comment end   ] */

/* [wxMaxima: section start ]
User input section
Modify for your design
   [wxMaxima: section end   ] */

/* [wxMaxima: input   start ] */
/* Thiele-Small Parameters of driver 
Eminence 2512-ii 
Edit these parameters for your driver. */

"Driver resonant frequency in Hz = ";         Fs: 37 ;
"Coil resistance in Ohms = ";                 Re: 5.04 ;
"Electrical Q factor, dimensionless = ";      Qes: 0.44 ;
"Mechanical Q factor, dimensionless = ";      Qms: 3.13 ;
"Compliance equivalent volume in liters = ";  Vas_liters: 147 ;
"Cone area in cm^2 = ";                       Sd_cm2: 519.5 ;
"Excursion limit in mm = ";                   Xmax_mm: 4.9 ;
"Coil inductance in mH = ";                   Le_mH: 0.46 ;
"Nominal driver impedance in Ohms = ";        Znom: 8 ;
"Thermal power rating of driver in Watts = "; Pe: 250 ;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Box design 
My little 12" neo system
Edit these parameters for your design. */

"Box volume in liters = "; Vbox_liters: 32 ;
"Number of drivers = "; Ndrivers: 1 ;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Port design

Note that the program lets you specify a round or rectangular port.
You can enter dimensions for both. The program will only use the
dimensions that it needs.

shapePort   = 0 for circular, 1 for rectangular
endPort     = 0 for both ends free
            = 1 for one end flanged, the other end free (typical)
            = 2 for both ends flanged */

shapePort: 0 ;
endPort: 1 ;

"Port tuning frequency in Hz = ";               Fport: 45 ;
"Diameter of port in cm = ";                    Dport_cm: 7.5 ;
"One dimension of rectangular port in cm = ";   Aport_cm: 8 ;
"Other dimension of rectangular port in cm = "; Bport_cm: 5 ;
"Number of ports = ";                           nPorts: 1 ;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Electrical inputs for graphs
These are typical values, unlikely to need changing. */

"Input power for sensitivity curve in Watts = "; PinSens: 1 ;
"Input power for excursion curve in Watts = ";   PinExc: Pe ;
"Distance to listener in meters = ";             Dist: 1 ;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Computation section
Modify at your own risk
   [wxMaxima: section end   ] */

/* [wxMaxima: input   start ] */
/* Physical constants */

"Atmospheric pressure in Pa = "; Patm: 101300;
"Adiabatic constant (dimensionless) = "; b: 1.4;
"Product of these in Pa = "; bpatm: b*Patm;
"Density of air in kg/m^3 = "; rho: 1.18;
"Speed of sound in air in m/s = "; c: 345;
"Reference pressure for SPL in Pa = "; splref: 2e-5;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Convert TS parameters into SI units */

"Vas in m^3 = "; Vas: Vas_liters/1000, numer;
"Sd in m^2 = "; Sd: Sd_cm2/10000, numer;
"Xmax in m = "; Xmax: Xmax_mm/1000, numer;
"Le in H = "; Le: Le_mH/1000, numer;
"Vbox in m^3 = "; Vbox: Vbox_liters/1000, numer;

"Port area in m^2 = "; if shapePort = 0 then (
    Rport: Dport_cm/200,
    Sport: nPorts * %pi * Rport^2)
else (
    Aport: aPort_cm/100,
    Bport: bPort_cm/100,
    if Aport < Bport then Rport: Aport/2 else Rport: Bport/2,
    Sport: nPorts*Aport*Bport), numer;

"End correction factor = "; if endPort = 0 then endCorrect: 0.85
    else if endPort = 1 then endCorrect: 0.732
        else endCorrect: 0.614;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
Ndrivers;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Electromechanical parameters */

"Driver resonant angular frequency in 1/s = "; ws: Fs*2*%pi, numer;
"Port angular frequency in 1/s = "; wport: Fport*2*%pi, numer;
"Driver spring constant in N/m = "; Km: bpatm*Sd^2/Vas, numer;
"Driver damping factor in N s / m = "; Rm: Km/Qms/ws, numer;
"Cone mass in kg = "; Mm: Km/ws^2, numer;
"Field length factor in T m = "; BL: sqrt(Re*Km/Qes/ws), numer;
"Box spring constant in N/m = "; Kbox: Sd^2*bpatm*Ndrivers/Vbox, numer;
"Voltage input for sensitivity curve in V RMS = "; VinSens: sqrt(PinSens/Ndrivers*Znom), numer;
"Voltage input for excursion curve in V peak = "; VinExc: sqrt(2*PinExc/Ndrivers*Znom), numer;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Port correction function */

portFact(w) := ev(if w = wport then 0.00001 else w^2/(w^2 - wport^2), numer)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* General purpose excursion formula */

X(f,V) := block(
    w: f*2*%pi,
    Zcoil: Re + %i*w*Le,
    Keff: (Km+Kbox)*portFact(w),
    ev(BL*V/Mm/Zcoil/(Keff/Mm + %i*w/Mm*(BL^2/Zcoil + Rm) - w^2),numer))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Sensitivity formula */

P(f,V) := block(
    w: f*2*%pi,
    ev(20*log(cabs(X(f,V)*w^2*rho*Sd*Ndrivers*portFact(w)/2/Dist/%pi/splref))/log(10),numer))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Waveform based excursion formula */

Xwavf(f,V) := ev((cabs(X(f,V)) + cabs(X(2*f,V)) + cabs(X(3*f,V)))/3, numer)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Port air speed formula */

Vport(f,V) := block(
    w: f*2*%pi,
    ev(cabs(w*X(f,V)*Sd*Ndrivers/Sport/abs(w^2/wport^2-1)),numer))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
/* Impedance formula */

Impedance(f, V) := block(
    w = f*2*%pi,
    ev(Re*V/cabs(V - BL*%i*w*X(f,V)),numer))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Output graphs
   [wxMaxima: section end   ] */

/* [wxMaxima: input   start ] */
/* Sensitivity Plot */

wxplot2d([P(x,VinSens)], [x,30,250],

[xlabel, "Frequency (Hz)"], 
[ylabel, "Sensitivity (dB SPL)"],
[gnuplot_preamble, "set title 'Sensitivity Plot';
set logscale x;
set xtics add (30, 40, 60, 80, 120, 140, 200)"])$

/* Excursion plot */

wxplot2d([1000*cabs(X(x,VinExc)), Xmax_mm], [x,30,250],

[xlabel, "Frequency (Hz)"], 
[ylabel, "Excursion (mm)"],
[legend, "Excursion", "Xmax"],
[gnuplot_preamble, "set title 'Excursion Plot';
set logscale x;
set xtics add (30, 40, 60, 80, 120, 140, 200)"])$

/* Waveform based excursion plot */

wxplot2d([1000*Xwavf(x,VinExc), Xmax_mm], [x,30,250],

[xlabel, "Frequency (Hz)"], 
[ylabel, "Excursion (mm)"],
[legend, "Excursion", "Xmax"],
[gnuplot_preamble, "set title 'Waveform Based Excursion Plot';
set logscale x;
set xtics add (30, 40, 60, 80, 120, 140, 200)"])$

/* Port air speed plot */

wxplot2d([Vport(x,VinExc)/c], [x,30,250],

[xlabel, "Frequency (Hz)"], 
[ylabel, "Air speed (Mach)"],
[gnuplot_preamble, "set title 'Port air speed plot';
set logscale x;
set xtics add (30, 40, 60, 80, 120, 140, 200)"])$

/* Impedance plot */

wxplot2d([Impedance(x, VinSens)], [x,30,250],

[xlabel, "Frequency (Hz)"], 
[ylabel, "Impedance (Ohms)"],
[gnuplot_preamble, "set title 'Impedance plot';
set logscale x;
set xtics add (30, 40, 60, 80, 120, 140, 200)"])$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Useful statistics
   [wxMaxima: section end   ] */

/* [wxMaxima: input   start ] */
ground(x,m):=floor(x*m)/m$
ft3: (0.254*12)^3$
in3: (0.254)^3$
"Box volume in cubic feet = "; vBox_ft3: Vbox_liters/ft3, numer;
"Box volume in cubic inches = "; vBox_in3 = Vbox_liters/in3, numer;
"Port length in m = "; lPort: Sport*bpatm/rho/Vbox/wport^2 - Rport*2*endCorrect, numer;
"Port length in cm = "; lPort_cm: lPort*100;
"Port length in inches = "; lPort_in: lPort_cm/2.54;
"Port diameter in inches = "; Dport_cm/2.54;
"Port width in inches = "; Aport_cm/2.54;
"Port height in inches = "; Bport_cm/2.54;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Additional Documentation
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
This project started many years ago with my own effort to understand
the physics behind speaker design. I ended up deriving the electromagnetic
theory from scratch, with input from various online references. The
report on that derivation is at:

    http://personalpages.tds.net/~fdeck/bass/speaker.pdf

I wrote a speaker design "program" as an Excel spreadsheet to see if I
could get the same results as the popular WinISD program. Then I noticed
that people occasionally ask for a speaker design program that runs on
an Apple Mac. Here it is.

I like using a computer algebra system (CAS) for computation, for a 
number of reasons. It handles complex numbers automatically, allowing my
formulas to be a lot cleaner. The workbook combines source code, user
input, results, and documentation, all in one place. The CAS is free, and
runs on multiple platforms.

The "Waveform Based Excursion Plot" is experimental. The standard excursion
plot assumes a sinewave input, but a bass guitar waveform contains a
mixture of harmonics, as does a typical music program. I am using as a
model, a simple scalar average of the first three harmonics to predict
peak cone excursion. I assume that the amplifier power rating determines
the amplitude of the input signal.
   [wxMaxima: comment end   ] */

/* Maxima can't load/batch files which end with a comment! */
"Created with wxMaxima"$