EXCELでスピーカーシミュレーション（１０）

１０．シミュレーションに使う定数の求め方

１０．シミュレーションに使う定数の求め方

　スピーカーの音圧をシミュレーションするに当たって必要になる種々の定数について、その求め方を説明します。

１０－１．ドライバーユニットの振動系質量　 $m$

　ドライバーユニットの振動系質量は、振動板をはじめ、エッジやダンパーの一部、ボビンやボイスコイル等の質量を合計することによって求めることができます。しかし、これらを測ろうとすると、せっかくのドライバーユニットを破壊しなければ測定できません。
　ところが、破壊しなくても、インピーダンス特性を測定することによって、質量を間接的に知ることができます。正しくは、インピーダンス特性から分かる共振周波数 $f_0$ を用います。
　まず、共振周波数 $f_0$ についてですが、ドライバーユニットは、質量 $m$ 、機械抵抗 $r$ 、スチフネス $s$ の直列回路でした。その合成インピーダンスは、

$j \omega m + r + \dfrac{s}{j \omega}$
で表すことができます。このうち、 $j \omega m$ と $\frac{s}{j \omega}$ は、角周波数 $\omega$ の関数ですので、周波数によってインピーダンスが変化します。共振とは、これらのインピーダンスが打ち消し合うことでインピーダンスが限りなく小さくなる（厳密には機械抵抗のみになる）ことを言います。
　つまり、これらのインピーダンスを足し合わせた結果ゼロになる周波数が、共振周波数ということになりますので、
$\begin{align} j \omega m + \dfrac{s}{j \omega} &= 0 \\ j \omega m &= - \dfrac{s}{j \omega} \\ j \omega^2 m &= - \dfrac{s}{j} \\ \omega^2 &= \dfrac{s}{m} \\ \omega &= \sqrt{\dfrac{s}{m}} \end{align}$
$\therefore f_0 = \dfrac{1}{2 \pi} \sqrt{\dfrac{s}{m}}$
　このように、共振周波数は、スチフネスと質量が分かれば計算によって求めることができます。
　この式において、質量 $m$ は分母に入っています。つまり、質量が増加すれば、共振周波数は低くなる方向に変化するということです。
　そこで、ドライバーユニット単体のときの共振角周波数 $\omega_1$ と、振動板に何か"おもり"をくっつけて質量 $m$ を $m+\Delta m$ としたときの共振角周波数 $\omega_2$ を、それぞれ測定します。すると、
$\begin{align} \omega_1^2 &= \dfrac{s}{m} \\ \omega_2^2 &= \dfrac{s}{m +\Delta m} \end{align}$
となります。ここで、 $s$ を消去するために、
$\begin{align} \omega_1^2 &= \dfrac{s}{m} \\ s &= \omega_1^2 m \end{align}$
を、 $\omega_2$ の式に代入すると、
$\begin{align} \omega_2^2 &= \dfrac{\omega_1^2 m}{m+\Delta m} \\ \omega_2^2 m+ \omega_2^2 \Delta m &= \omega_1^2 m \\ m (\omega_1^2-\omega_2^2) &= \omega_2^2 \Delta m \\ \therefore m &= \dfrac{\omega_2^2 \Delta m}{\omega_1^2 - \omega_2^2} \end{align}$
となり、質量 $m$ が求まります。ちなみに、角周波数を周波数に置き換えてもＯＫです。
$m = \dfrac{f_2^2 \Delta m}{f_1^2 - f_2^2}$
　つまり、測定した共振周波数と、くっつけた"おもり"の質量が分かれば、ドライバーユニットの振動系の質量 $m$ を計算によって求めることができるのです。
　ちなみに、インピーダンス特性は、ドライバーユニットに下記のような測定回路を接続して、発信器の周波数を変化させながら交流電圧計で電圧を測定することで求めることができます。

　この回路はドライバーユニットに負荷抵抗 $R$ を直列に接続しています。負荷抵抗 $R$ の抵抗値は、ドライバーユニットの定格インピーダンスよりも十分小さいものを選びます。この回路全体に流れる電流 $I$ は、ドライバーユニットのインピーダンスを $Z$ とすると、
$I = \dfrac{V_t}{Z+R}$
となります。直流回路ですので、同じ電流がドライバーユニットにも流れるため、ドライバーユニットに印加される電圧を $V_{sp}$ とすると、
$I = \dfrac{V_{sp}}{Z} = \dfrac{V_t}{Z+R}$
という等式が成り立ちます。これを、 $Z$ について解くと、
$\begin{align} V_{sp}Z+V_{sp}R &= V_tZ \\ Z(V_t-V_{sp}) &= V_{sp}R \\ \therefore Z &= \dfrac{V_{sp}}{V_t-V_{sp}}R \end{align}$
となりますので、回路全体の電圧と、ドライバーユニットの電圧を交流電圧計で測定すれば、その周波数のドライバーユニットのインピーダンスを知ることができます。
　これを周波数毎に測定することで、インピーダンス特性を知ることができます。

１０－２．ボイスコイルのインダクタンス　 $L$

　ボイスコイルのインダクタンスは、ソレノイドのインダクタンスの求め方を応用して算出します。ボイスコイルは、それ自体は空芯ソレノイドになりますが、ドライバーユニットとして組み込まれると、磁器回路のポールが磁芯のようになります。ここでは長岡係数を用いた方法で計算してみることにします。
　長岡係数を用いたソレノイドのインダクタンスは、

$L = \dfrac{k_n \mu_0 \mu_e \dfrac{\pi d^2}{4} n^2}{l}$
で求めることができます。
ここで、
$k_n$ は長岡係数
$\mu_0$ は真空の透磁率 $=4 \pi×10^{-7}$
$\mu_e$ はポールの実効比透磁率
$d$ はボイルコイルの巻きの直径
$n$ はボイスコイルの巻き数
$l$ はボイスコイルの素線の長さ
です。
　 $\mu_e$ はポールの直径に対する高さの比がたいていは1前後であると思われるので、一桁くらいの数字になるのではないかと思われます。
　 $d$ 、 $n$ 、 $l$ はドライバーユニットを分解してみないことには分かりませんので、だいたいの値を入れてからインピーダンス特性を測定してフィッティングしてみてください。
　そして長岡係数 $k_n$ ですが、下記の式で求めることができます。
$k_n= \dfrac {4}{3 \pi \sqrt {1-k^2}} \left( \dfrac {1-k^2}{k^2} K(k) - \dfrac { 1-2k^2}{k^2} E(k)-k \right)$
　ここで、 $k$ は、
$\dfrac{2r}{l} = \dfrac{k}{\sqrt{1-k^2}}$
で定義されますので、
$\begin{align} 2r \sqrt{1-k^2} &= kl \\ 4r^2(1-k^2) &= k^2l^2 \\ 4r^2-4r^2k^2 &= k^2l^2 \\ 4r^2 &= k^2(l^2+4r^2) \\ k^2 &= \dfrac{4r^2}{l^2+4r^2} \\ k &= \dfrac{2r}{ \sqrt{l^2+4r^2}} \\ &= \dfrac{1}{\sqrt{\frac{l^2}{4r^2}+\frac{4r^2}{4r^2} } } \\ \therefore k &= \dfrac{1}{\sqrt{1+(\frac{l}{2r})^2}} \end{align}$
となります。※ $r=\frac{d}{2}$ です。
　また、 $K(k)$ は $k$ についての第一種完全楕円積分で、 $E(k)$ は $k$ についての第二種完全楕円積分です。どうやらEXCELでは楕円積分の関数が用意されていないようです。そのため、下記の近似式を用いて計算するしかないようです。
第一種完全楕円積分
$\begin{align} K(k) &= \int_{0}^{\frac{\pi}{2}}(1-ksin^2 \theta )^{\frac{1}{2}} d \theta \\ &= \int_{0}^{\frac{\pi}{2}} \sum_{n=0}^\infty \dfrac{(2n-1)!!}{(2n)!!} k^n sin^{2n} \theta d \theta \\ &= \sum_{n=0}^\infty \dfrac{(2n-1)!!}{(2n)!!} \int_{0}^{\frac{\pi}{2}} sin^{2n} \theta d \theta \\ &= \dfrac{\pi}{2} \left\{ \sum_{n=0}^\infty \left( \dfrac{(2n-1)!!}{(2n)!!} \right)^2 k^n \right\} \end{align}$
第二種完全楕円積分
$\begin{align} E(k) &= \int_{0}^{\frac{\pi}{2}} \sqrt{1-ksin^2 \theta} \\ &= \dfrac{\pi}{2} \left\{ 1-\sum_{n=1}^\infty \left( \dfrac{(2n-1)!!}{(2n)!!} \right)^2 \dfrac{k^n}{2n-1} \right\} \end{align}$
　実際には $n$ は $\infty$ までは計算できませんので、50くらいで計算してみますと、こんな感じです。

１０－３．放射インピーダンス　 $m_a, r_a$

　以前の記事で、放射インピーダンスにおける放射質量 $m_a$ は、ka<<1の周波数では近似式として $\frac{8}{3}\rho_0(\frac{S_d}{\pi})^{\frac{3}{2}}$ を用いることができ、放射抵抗は近似式として $\rho_0c\frac{k^2S_d^2}{2\pi}$ を用いることができるという説明をしました。しかし、実際にはこれらのインピーダンスは周波数の関数として表すことができます。
　振動板を半径aの円盤と仮定したときの放射インピーダンスは、導出は省略しますが、下記の式で表すことができます。

$m_a = S_d \rho_0 \dfrac{1}{k} \left( \dfrac{H_1(2ka)}{ka} \right)$
$r_a = S_d \rho_0 c \left( 1- \dfrac{1}{ka}\right)J_1(2ka)$
　ここで、 $H_1$ はストルーブ関数、 $J_1$ は第一種ベッセル関数です。 $k$ は波定数 $\frac{\omega}{c}$ です。
　第一種ベッセル関数はEXCELでもサポートしており、
=BESSELJ(*,1)
です。*のところは、 $2ka$ を計算したセルを指定しましょう。次数は1です。
　ところが、ストルーブ関数はEXCELではサポートしていないようですので、下記の近似式を使うしかないようです。
$H_1(2ka) = \dfrac{2}{\pi}-J_0(2ka)+\left( \dfrac{16}{\pi}-5 \right) \left( \dfrac{sin(2ka)}{2ka} \right)+\left( 12-\dfrac{36}{\pi} \right) \left( \dfrac{1-cos(2ka)}{(2ka)^2} \right)$
　 $J_0$ は前記と同じ第一種ベッセル関数ですが、次数は0です。
　ちなみに、上記の式は機械インピーダンス表示ですので、比音響インピーダンスを用いる場合は $S_d$ で割り、音響インピーダンスを用いる場合は $S_d^2$ で割る必要があることをお忘れなく。

前へ　９．ダブルバスレフ型エンクロージャ
 最初に戻る　１．波動方程式を求める