第4章
電気・機械特性パラメータの算出方法
ファイン MEMS システム化設計プラットフォーム
研究開発プロジェクト
平成 21 年 5 月 25 日

目次
第 4 章 電気・機械特性パラメータの算出方法

4.1

電気機械特性パラメータの算出手順 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2

バネコンポーネントのパラメータ算出 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.3

4.4

4.2.1

バネ定数 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.2.2
4.2.3
4.2.4

等価質量 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

34

直列に繋がれたバネコンポーネントの取り扱い . . . . . . . . . . . . . . . . . . . . .

5892 端子櫛歯コンポーネントのパラメータ算出 . . . . . . . . . . . . . . . . . . . . . . . . . .
4.3.1 容量 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

10
10

4.3.2 機械抵抗 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
音響コンポーネントのパラメータ算出 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.4.1 音響容量 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

13
15
16

音響抵抗 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
イナータンス . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

17
18

4.4.4 音響平行平板アクチュエータ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
その他のコンポネントのパラメータ算出 . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.5.1 磁気駆動梁のパラメータ算出 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

19
21
21

ユニモルフ圧電駆動梁のパラメータ算出 . . . . . . . . . . . . . . . . . . . . . . . .

23
24

4.4.2
4.4.3
4.5

3
3

4.5.2
4.5.3

電気抵抗 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

静電平行平板のパラメータ算出 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2

第 4 章 電気・機械特性パラメータの算出方法
MEMS 等価回路では、各コンポーネントが回路素子で表現される。例えば、バネコンポーネントはバネの機能を表すコンデンサとバネの質量の寄与を表すインダクタを繋げたものとして表現される。これらの素子定数を決定する為には、コンポーネントごとに、その各種パラメータ（形状パラメータ、物性パラメー
タなど）から電気機械特性パラメータを算出する必要がある。等価回路ジェネレータにおいて、これらのパラメータは電気機械特性パラメータ抽出ソフト（以下では、抽出ソフトと表記する）を用いて自動的に計算される。本章では、抽出ソフトに実装されている各種パラメータの算出方法について、概説する。

4.1

電気機械特性パラメータの算出手順

前章までで述べたように、MEMS 等価回路は、Lagrange 方程式を直流動作点において線形化して得られた動作行列に基づいて導出されている。その為、MEMS 等価回路の回路素子の定数には、直流動作点での値を用いる必要がある。従って、電気機械特性パラメータを算出するには、

1. 直流バイアス電圧を印加時の可動部の変位を求める
2. 各コンポーネントの端子の変位を用いて、コンポーネントごとに電気・機械特性パラメータを計算するという手順を踏むことになる。1. の計算では、2 端子櫛歯などの電気機械変換コンポーネントに直流バイ
アス電圧を印加した際に発生する力を可動部に加え、各節点の力の釣り合いから変位を求めることになる。
このようにして求まった変位を用いて、各コンポーネントの電気機械特性パラメータを計算する。

4.2

バネコンポーネントのパラメータ算出

バネコンポーネントは、図 4.1 に示す 4 種類である。バネコンポーネントに算出が必要なパラメータを以下に示す。

1. バネ定数（x, y, z 方向)
2. 等価質量（x, y, z 方向)
3. 長軸方向の抵抗
以下で、それぞれのパラメータの算出方法を説明する。
1

1

2

1

1

y

2

x
z
2

2

Beam

O-shape

Double Folded1

Double Folded2

図 4.1: バネコンポーネントの概略図

3図 4.2: 梁バネ

4.2.1

バネ定数

抽出ソフトでは、直流動作点でのバネの伸縮量、およびバネの非線形性を考慮して、バネ定数の算出を行っている。
図 4.2 に示した櫛歯アクチュエータを支持する梁バネは、変位が梁幅程度になると長軸の伸びの効果に
よって非線形性を示す。このような梁バネのバネ定数は、両持ち梁のバネ定数と等価で、
[ 4
]
[
]
π
π 4 ( x )2
π 4 ( x )2
f = k0+≈ k0 1 +96 128 w
128 w

(4.1)

と表せる [1]。ここで、 k0 は線形なバネ定数（今の場合、両持ち梁なので、k0 = 192EI/L3 ) である。
パラメータ抽出ソフトでは、O 型バネやダブルフォールデッドバネも扱う必要がある。z 方向のバネ定数の算出の際には、曲げだけでなくねじれの寄与も考える必要がある為、ｚ 方向のバネ定数を式 (4.1) に基づくバネの合成によって求めることはできない。そこで、x, y, z 方向のバネ定数算出を統一的に行う為、12自由度の梁モデルを用いる [2]。

12 自由度の梁モデル (図 4.3) は、2 つの端点それぞれに、並進変位 (x, y, z) と角変位 (ϕx , ϕy , ϕz ) を自由度を持たせたモデルで、曲げ、ねじりを一括して扱うことが出来る。この梁モデルでは、力と変位の関係式
が以下のように表わされる。
f⃗ = ([K] + T [KG ])⃗x(4.2)

図 4.3: 12 自由度梁モデル

4ここで、T は長軸方向に働く張力を示しており、変位が大きくなるほど硬くなるというバネの非線形性を引き起こす項である。

O 型バネやダブルフォールデッドバネを上記の梁をつなぎ合わせたものとしてモデル化し、バネの端点に対して各方向に加えた力と変位の関係から、バネ定数を求めている。

図 4.4: 梁の繋ぎあわせ

4.2.2

等価質量

MEMS 等価回路を導出する際、デバイスは集中定数系（バネ-質点系）としてモデル化されている。しかし、バネコンポーネントは質量が分布している為、質量として、実質量（=[密度] × [体積]) ではなく、等価な集中質量系の質量 (等価質量) を用いる必要がある。本小節では、エネルギー法を用いて等価質量を求める方法を示す。

エネルギー法エネルギー法とは、分布質量系と集中質量系の運動エネルギーの最大値を等置することにより、等価質量を求める方法である。具体的には、以下に示す手順で計算する [3]。

1. バネの各部分は同位相で調和振動すると仮定する。
u(ξ, t) = Y (ξ)sin(ωt)(4.3)

ここで、Y (ξ) は振幅分布であり、抽出ソフトでは静撓み曲線を採用している。

2. バネの運動エネルギーの最大値を計算する。微小領域の運動エネルギーを積分することにより、バネ全体の運動エネルギーを計算する。
[∫
{
}2 ]
∫
ρA ∂u(ξ, t)ω2dξ ρA{Y (ξ)}2(4.4)
Tmax = max
dξ
=
2
∂t
2
ここで、ρ は密度、A は断面積である。

3. 等価なバネ-質点系の運動エネルギーの最大値 Tmax,0 を計算する。
4. Tmax,0 , Tmax を等置して、等価質量を定める。

5例えば、図 4.5 に示す先端が回転拘束されたカンチレバーの等価質量は以下のようにして求められる。振幅分布として先端の回転が拘束されたカンチレバーの静たわみ曲線を採用すると、カンチレバーの変位分
布は、

[ ( )
( )3 ]2ξξu(ξ, t) = y0 3−2sin(ωt)ll

ここで、y0 は先端の振幅である。このとき、カンチレバーの運動エネルギーの最大値は、
[∫
{
}2 ]
ρA ∂u(ξ, t)Tmax = max
dξ
2
∂t
[
( )2
( ) 3 ]2∫ρwtω 2 y02 l
ξ
ξ
=
dξ 3
−2
2
l
l
0
)
(
1 13
=
ρwtl ω 2 y02
2 35である。一方、梁の先端にマスが集中すると考えた場合、運動エネルギーの最大値は、
[
}2 ]
{
1
1
∂
= Meq ω 2 y02
Tmax,0 = max Meq(y0 sinωt)
2
∂t
2

(4.5)

(4.6)

(4.7)
(4.8)

(4.9)

式 (4.8) と式 (4.9) を等値すると、先端の回転を拘束されたカンチレバーの等価質量は、

Meq =

13
ρwtl
35

(4.10)と求まる。

図 4.5: 先端を回転拘束したカンチレバー

等価質量の算出手順
梁バネ以外のバネコンポーネント（O 型バネ、ダブルフォールデッド 1, 2）の場合、バネを図 4.3 に示す
梁のつなぎ合わせとして表現し (図 4.4)、エネルギー法を適用する。具体的には、以下の手順を踏むことになる。

1. バネに x, y, z 各方向に力を加えたときの各節点 (梁のつなぎ目) の変位を計算する。
2. 得られた各節点の変位から、それぞれの梁の運動エネルギーの最大値を求める。
3. その表式を用いて、バネ全体の運動エネルギーの最大値を求める。
4. バネに等価な集中質量系の運動エネルギーと等値することにより、等価質量を求める。
以下では、上記 2. 以降の具体的な計算方法について、順に説明していく。

6

梁の運動エネルギーの最大値
まず、梁の左右の端点の変位を xi , yi , zi , ϕxi , ϕyi , ϕzi (i = 1, 2) とする。これらの値は、すでに得られているとして話を進める。
梁の長軸上1 の位置 ξ における各方向の変位は、以下のように表される。
()
ξξ
x(ξ) =
1−
x1 + x2lly(ξ) = f1 (ξ)y1 + f2 (ξ)ϕz1 + f3 (ξ)y2 + f4 (ξ)ϕz2(4.11)= f1 (ξ)z1 − f2 (ξ)ϕy1 + f3 (ξ)z2 − f4 (ξ)ϕy2(4.13)z(ξ)
ここで、形状関数 fi (ξ) はf1 (ξ) =f3 (ξ) =(4.12)

[( )
( )2
( )3
( )2 ( )3 ]
ξ
ξ
ξ
ξ
ξ
1−3
+2
, f2 (ξ) = l
−2
+
l
l
l
l
l
[ ( )
]
( )2
( )3
(
)
2
3
ξ
ξ
ξ
ξ
3
−2
, f4 (ξ) = l −
+
l
l
l
l

梁が x 方向に振動するとき、梁の運動エネルギーの最大値は、
[∫
{
}2 ]
l
ρA
∂
ρω 2
Tx,max = maxdξx(ξ) sin(ωt)=
Th (l, A, x1 , x2 )2∂t20
ここで、Th を以下のように定義した。

[∫ l

Th (l, A, x1 , x2 ) ≡dξA{x(ξ)}2 = lA
0

1 2 1
1
x1 + x1 x2 + x22
3
3
3(4.14)(4.15)(4.16)]

一方、梁が y 方向に振動するとき、梁の運動エネルギーの最大値は、
[∫
{
}2 ]
l
∂
ρω 2
ρA
Ty,max = maxy(ξ) sin(ωt)=
Tv (l, A, y1 , ϕz1 , y2 , ϕz2 )
dξ
2
∂t
2
0
ここで、Tv を以下のように定義した。
∫ l
Tv (l, A, y1 , ϕz1 , y2 , ϕz2 ) ≡dξA{y(ξ)}20[= lA F11 y12 + F22 ϕ2z1 + F33 y22 + F44 ϕ2z2

(4.17)(4.18)(4.19)+2F12 (y1 ϕz1 − y2 ϕz2 ) + 2F13 y1 y2 + 2F14 (y1 ϕz2 − y2 ϕz1 ) + 2F24 ϕz1 ϕz2 ]ここで、Fij は、以下で定義される。

∫ l
Fij ≡

dξfi (ξ)fj (ξ)(4.20)

0l2
11, F12 = −F34 =
l
105
210
9
13
l2
F13 =
, F14 = −F23 = −L, F24 = −
70
430
140
同様にして、梁が z 方向に振動するとき、運動エネルギーの最大値は、
F11 = F33 =

13
,
35

Tz,max==F22 = F44 =ρω 2
Tv (l, A, z1 , −ϕy1 , z2 , −ϕy2 )2ρω 2
Tv (l, A, z2 , ϕy2 , z1 , ϕy1 )2と表すことができる。
1 厳密には、中立軸上。ここでは、梁の断面は長方形なので、中立軸は長軸に一致する

7(4.21)(4.22)(4.23)

O 型バネの x 方向の等価質量
前小節で求めた梁の運動エネルギー最大値の表式を用いて、O 型バネの運動エネルギーの最大値を求めることができる。例えば、x 方向に振動する場合 (図 4.6)、O 型バネの運動エネルギーの最大値は、oshapeTx,max=Ty,max (L0 , W0 t, x0 , ϕz0 , x2 , ϕz2 ) + Ty,max (L2 , W1 t, x3 , ϕz3 , x4 , ϕz4 )+Ty,max (L2 , W1 t, x7 , ϕz7 , x6 , ϕz6 ) + Ty,max (L3 , W0 t, x5 , ϕz5 , x1 , ϕz1 )+Tx,max (L1 , W2 t, x3 , x2 ) + Tx,max (L1 , W2 t, x4 , x5 )+Tx,max (L1 , W2 t, x5 , x6 ) + Tx,max (L1 , W2 t, x2 , x7 )(4.24)と表せる。一方、等価なバネ-マス系の運動エネルギーの最大値は、
0
Tx,max=Mxeﬀ ω 2(x1 − x0 )22(4.25)と書ける。式 (4.24) と (4.25) の右辺を等値すると、O 型バネの x 方向の等価質量

Mxeﬀ =

ρ[Uv (L0 , W0 t, x0 , ϕz0 , x2 , ϕz2 ) + Uv (L2 , W1 t, x3 , ϕz3 , x4 , ϕz4 )(x1 − x0 )2+Uv (L2 , W1 t, x7 , ϕz7 , x6 , ϕz6 ) + Uv (L3 , W0 t, x5 , ϕz5 , x1 , ϕz1 )(4.26)+Uh (L1 , W2 t, x3 , x2 ) + Uh (L1 , W2 t, x4 , x5 )+Uh (L1 , W2 t, x5 , x6 ) + Uh (L1 , W2 t, x2 , x7 )]と求まる。以上と同様にして、y, z 方向の等価質量も計算している。

図 4.6: O 型バネの x 方向運動 (赤線は中央線)

4.2.3

電気抵抗

バネの電気抵抗は、LSI の配線抵抗を見積もる方法に準じて計算を行っている。例えば、ダブルフォール
デッドバネ 2 の抵抗算出では、バネの中央線 (図 4.7 における点線) の長さを L、断面積を A、バネの材質の比抵抗を ρr として、

A
(4.27)
L
から計算を行う。また、O 型バネ、ダブルフォールデッドバネ 1 では中央線の分岐が生じるが、その場合
R = ρrは、抵抗が並列接続していると考えて、全体の抵抗を計算する。

8

図 4.7: ダブルフォールデッドバネ 2 の抵抗算出

4.2.4

直列に繋がれたバネコンポーネントの取り扱い

バネコンポーネントが直列に繋がれた場合、そのバネ定数は個々のバネコンポーネントのバネ定数の単純な重ね合わせ（各バネ定数の逆数の和の逆数）では表せるとは限らない。例えば、図 4.8 に示すミアンダ
梁をダブルフォールデッドバネ 2 を直列につなげて構成した場合、その y 方向のバネ定数は、個々のダブ
ルフォールデッドバネ 2 のｙ方向のバネ定数の重ね合わせでは表現できない。これは、モーメントの効果による。したがって、バネコンポーネントが直列に繋がれた場合のバネ定数を求める為には、個々のバネコ
ンポーネントをつなぎ合わせた全体としての取扱いが必要である。
そこで、抽出ソフトでは、直接連結しているバネをグルーピングし、全体としてのバネ定数を求める。
個々のバネのバネ定数は、バネ形状やサイズパラメータの違いによらず、全体のバネ定数にバネの個数を掛けた値としていることに注意されたい2 。もちろん、このようにして個々のバネ定数を算出しても、MEMS
等価回路を用いたシミュレーション結果には影響を与えない。なぜなら、MEMS 等価回路では、同じバネ
定数のバネを直列つないだ場合、その合成バネ定数は、個々のバネ定数をバネ数で割ったものと等価になるからである。

図 4.8: ミアンダ梁

2 抽出ソフトの処理が非常に煩雑になる為、このようにしている

9

図 4.9: 2 端子櫛歯コンポーネントの概略図

表 4.1: 2 端子櫛歯において算出が必要な電気機械特性パラメータ
パラメータ名単位

質量kg

直流動作点における容量F

直流動作点における容量の 1 階微分
(
)
−
−
−
∂C(→
x ) ∂C(→x ) ∂C(→
x)
,
,
∂x
∂y
∂z

F/m

備考

x, y, z は構造部 2 に対する構造部 1 の相対変位を表す。

2

直流動作点における容量の 2 階微分
( 2 →
)
−
−
−
−
−
∂ C(−
x ) ∂ 2 C(→x ) ∂ 2 C(→x ) ∂ 2 C(→x ) ∂ 2 C(→x ) ∂ 2 C(→
x)
,
,
,
,
,
∂x2
∂x∂y
∂x∂z
∂y 2
∂y∂z
∂z 2

F/m

機械抵抗

Ns/m

(x, y, z 成分)

2 端子櫛歯コンポーネントのパラメータ算出

4.3

図 4.9 に示した 2 端子櫛歯コンポーネントに関して算出が必要なパラメータを、表 4.1 に纏めた。

4.3.1

容量

櫛歯の厚みが等しい場合の容量
櫛歯の厚みが等しい場合、直流動作点での容量は、直流バイアス印加時の櫛歯変位に基づいて、以下の 3
つの寄与の和として算出する (図 4.10)。

1. 重なり容量 (赤字)
重なり長さ部分にある平行平板容量の値を計算した。

2. 櫛歯の先端と対向電極間の容量 (緑字)
一次近似として、櫛歯の先端部と対向部に構成される平行平板容量を計算した。

10

3. 櫛歯の上下面の寄与 (青字)
これはフリンジ電場の寄与に相当する。この寄与を計算するために、Johnson らによって導かれた解析式 [4] を採用する。

[ {(
} (
)2
) (
)]
2ϵxov
w
w
2g
C=
ln
+1 −1 + 1+
ln 1 +πggw(4.28)

図 4.10: 2 端子櫛歯の容量

図 4.11: 解析モデル

櫛歯の厚みが異なる場合の容量
櫛歯の厚みが異なる場合は、垂直駆動櫛歯電極に相当する。垂直駆動櫛形電極は、水平駆動櫛形電極と比較すると、重なり長さが短い、アスペクト比が小さい等の理由でフリンジ電場の寄与が大きい。この為、容量変化を水平駆動櫛形電極に対してしばしば行われるように、平行平板近似を用いて記述すると誤差が大きい。より精度良い結果を得る為には、有限要素法 (FEM) もしくは境界要素法 (BEM) に基づく数値シ
ミュレーションや、等角写像法を用いた解析的な計算を行うことが考えられる。FEM や BEM に基づく数値シミュレーションでは、任意形状に対して容量計算が可能という長所がある反面、計算時間が掛かり、同じ幾何配置でもパラメータが違えばそれぞれ計算が必要という短所がある。一方、等角写像法を用いた解
析的な計算では、容量の解析式が得られるという長所がある反面、非常に単純な形状を除いては閉じた形式の容量式を得られないという欠点がある。パラメトリックな設計に適用することを考えた場合、そして抽出ソフトは Web ライブラリシステム上で動作するアプリケーションで計算に時間を掛けられないという点を勘案した場合、フリンジ電場の寄与を適切に記述しながら、かつ閉じた形式の容量式を得ることは重要11である。そこで、抽出ソフトでは、垂直駆動櫛形電極を対象として、解析的な 2 次元電場計算で得られる
閉じた形式の容量式と 2 次元 FEM 解析結果を組み合わせた半解析的な容量式を別途導出し [5]、この式を用いて、直流動作点における櫛歯容量を算出している。
半解析的な容量式は、以下の仮定に基づいている。

1. 基板の影響は無視する
2. 櫛歯間の電場は x 座標（櫛歯のストローク方向の座標）に依存しない
解析的な 2 次元電場計算による容量式には、Yeh らの導出式 [6] を適用した。Yeh らは、櫛歯の幅が無視できる程に小さいコムドライブを対象に (図 4.11 参照)、2 次元ポアソン方程式を等角写像法で解き、電場に関する閉じた解析式を得た。

E(w) =A(a, b, c) =−jA(a, b, c)√(w − a)(w − b)(w − c)(w + 1)√
V (a − c)(b + 1)
2[F (δ, q) + F (κ, q)](4.29)(4.30)

ここで、V は櫛歯間電圧、F (κ, θ) は第１種不完全楕円関数、δ, κ, q は形状パラメータで決まるパラメータである。式 (4.29) に示す電場分布を用いると、単位長さあたりの櫛歯間容量は、∫
ϵ a
ϵF (π/2, r)
C0 =E(u)du =V b
F (δ, q) + F (κ, q)(4.31)と求まる。この解析では、櫛歯の幅の効果を無視している為、解析式で計算した容量は、2 次元 FEM 解析
結果より小さい (図 4.12)。

図 4.12: 容量の z 方向依存性

その為、櫛歯の幅の効果を考慮する為に、容量式に 2 つの項を加えた。具体的には、ストリップ容量、および 4 次の補正項の和を付加した。ストリップ容量としては、櫛歯幅に等しいストリップの容量式を適用
した。

)
(
πw
K sin 2(w+g))(
Cstrip = ϵπw
2K cos 2(w+g)

ここで、K(q) は第 1 種完全楕円積分である。

12(4.32)

一方、補正項は、(1) 櫛歯対の一方が他方の中央に位置するときに容量が最大になり、(2) この位置に関して容量の z 変位依存性が対称になる点を考慮して設定した。

Ccor = α(z − 2z0 )2 z 2 + β(4.33)

ここで、z0 は可動櫛歯が固定櫛歯の中央に位置する場合の変位である。
補正項には未知のパラメータ α, β が含まれる。これらは、櫛歯幅 1, 2, 3, 4, 5, 6um に対して 2 次元 FEM
解析で得られた容量の z 方向依存性データ（一方が他を包含する範囲において) を C0 + Cstrip + Ccor によって最小二乗フィッティングを行うことで求めた。電極幅を変えた場合 (図 4.13)、およびギャップ幅を変えた
場合 (図 4.14) に対して、得られた半解析的な容量式を z 方向変位の関数としてプロットした。これらのグラフから、得られた容量式は、２次元 FEM 解析とよい一致を示していることが確認できる。

● FEM (w=6µm)－ 半解析的容量式 (w=6µm)■ FEM (w=1µm)－ FEM (w=1µm)
・・・

櫛歯幅ゼロの解析式

図 4.13: 容量の z 方向依存性（幅の値を振った場合)● FEM (g=3µm)－ 半解析的容量式 (g=3µm)■ FEM (g=4µm)
・・・ 半解析的容量式(g=4µm)

図 4.14: 容量の z 方向依存性 (ギャップ幅を振った場合)

4.3.2

機械抵抗

2 端子櫛歯コンポーネントの機械抵抗は、図 4.15 に示す「スライドフィルムダンピング」と「スクイーズフィルムダンピング」の 2 種類の粘性抵抗の解析式に基づいて算出している。そこで、本節では、ダンピングの解析式について説明した後、これらの式を基づいた機械抵抗の算出方法を示す。

13

図 4.15: squeeze ﬁlm damping・sliding ﬁlm damping

粘性抵抗の解析式
スライドフィルムダンピングは、平板が基板に対して平行に移動するときに働く粘性抵抗である。抽出ソ
フトでは、平板-基板間の流体速度の分布が線形な流れ（クエット流）である場合の粘性抵抗を考慮している。この場合、粘性抵抗力は速度勾配 v/G(G は平板-基板間のギャップ距離) に比例し、Fdamping = µAv
G(4.34)と表わされる。
一方、スクイーズフィルムダンピングは、平板が基板に対して垂直に（基板-平板間の流体を押し潰すよ
うな方向に）移動する際に働く粘性抵抗である。抽出ソフトでは、十分長い時間振動させた後、流体が非圧縮性のように振舞うと仮定した場合の解析式を用いている。

Fdamping = µ
ここで、補正係数 β(ϵ) はβ(ϵ) = 1 −β(w/t)w3 t
G3(w < t)∞( nπ )
192ϵ ∑ 1tanh
π 5 n=1,3,··· n52ϵ(4.35)(4.36)と表される。各運動方向の機械抵抗は、櫛歯の各部分に働くダンピングの寄与の和で算出する。詳しくは、次節以降で述べる。

x 方向の機械抵抗
運動方向が x 方向の場合の機械抵抗は、以下の 3 つの寄与の和とする。

1. 櫛歯ギャップの（歯の側面間の）スライドフィルムダンピング
2. 基板との間のスライドフィルムダンピング
3. 歯の先端との間のスクイーズフィルムダンピング

y 方向の機械抵抗
運動方向が y 方向の場合の機械抵抗は、以下の 3 つの寄与の和とする。

1. 櫛歯ギャップの（歯の側面間の）スクイーズフィルムダンピング
2. 基板との間のスライドフィルムダンピング
3. 歯の先端との間のスライドフィルムダンピング
14

z 方向の機械抵抗
運動方向が z 方向の場合の機械抵抗は、以下の 3 つの寄与の和とする。

1. 櫛歯ギャップの（歯の側面間の）スライドフィルムダンピング
2. 基板との間のスクイーズフィルムダンピング
3. 歯の先端との間のスライドフィルムダンピング算出される機械抵抗の値に関する注意
抽出ソフトで算出される機械抵抗の値は、以上で説明したような理想的な設定において得られる解析式に基づいている。しかし、実際のデバイスの運動に伴う雰囲気流体の挙動は非常に複雑であるため、抽出ソ
フトで算出される機械抵抗の値はあくまでも目安とお考えください。

4.4

音響コンポーネントのパラメータ算出

音響コンポーネントは音響回路を構成するコンポーネントである。その一覧を図 4.16 に示した。ここで、音響回路とは、管などを伝播する音波の挙動を電気回路に対応させたものということができる。音波の伝播は、音圧 p と粒子速度 v （音波の進行速度) で特徴づけられる。例えば、体積 W0 の小空洞の表面に断面
積 S の小窓を開けて、そこを目掛けて、音波を入射させたとする。音波は縦波 (進行方向に平行に振動) な
ので、小窓表面を音圧 p で押す。その結果、空洞内の空気が圧縮され、体積減少 δw が生じる。このときの体積減少率は、音圧に比例する [7]。

δw
p
=
W0
K

(4.37)

ここで、K は体積弾性率と呼ばれ、値が大きいほど圧縮しにくいことを表す量である。一方、体積の減った部分 δw には、小窓から外部の空気が流れ込む。その流入総量は減った体積に一致するので、∫
δw = Svdt(4.38)

なお、時間の積分を取ったのは、粒子速度が時々刻々と変動するからである。以上の 2 式をまとめると、
∫
K
p=
udt
(4.39)
W0
ここで、表現の簡略化の為、粒子速度 v の代わりに、それに断面積を掛けた体積速度u ≡ vS(4.40)という量 を用いている（以下も同様に、u を用いる)。一方、電気回路において、容量 C を電圧 V で充電する場合、その挙動は

1
C

∫
Idt

(4.41)

電流 I⇐⇒ 体積速度 u(4.42)

電圧 V⇐⇒ 音圧 p(4.43)

V =と表される。
式 (4.39) と (4.41) を機械的に比較すれば、150という対応があり、空洞は電気回路でいう容量と同じ作用をすることが分かる。その意味で、量 W
K は、音響容量と呼ばれる。同様に、電気回路における式と音響系における式を比較することによって、音響抵抗、
イナータンス (電気回路のインダクタに相当) を定義することができる。音響抵抗は音波が細管を通る際に受ける圧力損失の大きさに関係する量であり、イナータンスは媒質の重さに伴う慣性に関係する量である。
これらについては、以後の節で述べる。
抽出ソフトでは、各々の音響コンポーネントの「音響抵抗」「音響容量」「イナータンス」を算出する。尚、
音響平行平板アクチュエータは、音響-電気-機械の相互変換コンポーネントなので、平板間容量等の電気特性パラメータや平板のバネ定数などの機械特性パラメータも算出される。

図 4.16: 音響コンポーネント

4.4.1

音響容量

前小節において、体積 W0 の小空洞の表面に断面積 S の小窓を開けて、そこを目掛けて、音波を入射させる場合を取り扱った。結果として、音圧 p と音波の体積速度 u の関係は、
∫
K
udt
p=
W0

(4.44)と求められた。この式において、式 (4.42), (4.43) に示す電気系-音響系の対応づけを行うと、小空洞は電気回路における容量と同じ作用を果たし、電気回路の容量に対応する音響容量を定義できる。

CA =

W0
K

(4.45)

ここで、体積弾性率 K は媒質に固有の値で、媒質の密度 ρ と音速 c を用いて、

K = ρc2(4.46)と表せる [7]。尚、音響容量は、空洞だけでなく、両端の開いた管でも定義できる。その場合、式 (4.44) における音圧 p は、両端の音圧差になる。

16

図 4.17: 円筒内の流れ

4.4.2

音響抵抗

音響抵抗の定義
管を音波が通過する場合、空気の流れを妨げる物質 (例えば、繊維) が管に入っていると、摩擦や粘性の効
果によって、音波の圧力の低下が生じる。また、空気の流れを妨げる物質が無い場合でも、空気の内部摩擦の効果によって、音波の圧力が低下する。このような圧力低下は、体積速度 u ≡ vS (S：断面積) に比例し、

p = RA u(4.47)と表すことができる。このとき、RA は音響抵抗と呼ばれている。
抽出ソフトでは、流れを妨げる物質が無い状況での空気の内部摩擦の効果による音響抵抗を算出している（以下の式の導出方法は、次小節を参照）。
円断面の場合： RA
矩形断面の場合： RA

=

8µl
πR4

= f Re(4.48)µl(1 + b/a)2
8ab3(a, b は断面の辺の長さ; a ≥ b)

(4.49)

ここで、R は円断面の半径、a, b (a ≥ b) は矩形断面の縦・横の長さ、f Re は矩形断面のアスペクト比で決まる量、µ は粘度、l は管の長さである。尚、繊維など空気の流れを妨げる物質による音響抵抗を考慮した
いときは、ユーザ自身が別途算出した値を用いることになります。

空気の内部摩擦による音響抵抗
本小節では、流れを妨げる物質が無い状況での空気の内部摩擦の効果による音響抵抗を導出する。
図 4.17 に示したような真っ直ぐ円管内において、入口部で一様に入った流れには、下流に進むにつれて速度勾配が生じるが、そのうち管内の速度分布は変化しなくなる (図 4.17 の十分に発達した領域)。このような定常状態において、管内を平均流速 u で流体が流れているとき、長さ l の部分に生じる圧力損失は、

Darcy-Weisbach 式で表される。
∆p = fl ρu2
Dh 2(4.50)

ここで、f は管摩擦係数、Dh は流体力学半径、ρ は流体の密度である。流体力学半径は、

Dh =4[断面積][周囲長さ]

17(4.51)

と表せる。一方、管摩擦係数は、層流領域3 であれば、以下のように表せる4 。
円断面：f Re=64

(4.52)

矩形断面：f Re=

34.206(b/a)4 − 116.54(b/a)3 + 171.78(b/a)2 − 128.49(b/a) + 95.96(4.53)

ここで、a, b (a ≥ b) は矩形の縦・横の長さ、Re はレイノルズ数であり、粘度 µ 等を用いて、

Re =

ρDh u
µ

(4.54)と表せる。従って、圧力損失は、

8µl
u (R は半径)πR4µl(1 + b/a)2u
矩形断面：∆p = f Re
8ab3
円断面：∆p=(4.55)(4.56)

ここで、管内の流速 v の代わりに、体積速度
円断面：u= πR2 v

(4.57)

矩形断面：u= abv(4.58)

を用いている。式 (4.55)、(4.56) から、円断面・矩形断面を有する管の音響抵抗が求まる。

4.4.3

イナータンス

イナータンスは媒質の重さに伴う慣性に関係する量である。短い管の断面 S に音圧 p が加わると、細管内の媒質に pS という力が作用する。管の長さ l が波長に比べて十分短い場合5 、管内の媒質は一体になって運動すると考えてよい。このとき、一体になって運動する媒質部分の運動方程式は、

dv
dt
ρl du
∴p=
S dt

pS = ρlSp(4.59)

Sl
図 4.18: イナータンス

3 レイノルズ数 Re < 2300 の場合
4 矩形断面の式は、以下の文献記載の f Re のアスペクト比依存性をフィッティングして求めた式です。

F. M. White, Fluid mechanics, McGraw-Hill (1994).
5 波長は [音速 (m/s)]/[周波数 (Hz)] から計算できる。例えば、室温の空気中の音速が 340 m/s、周波数として 100kHz をとると、波長は 3.4 mm となる。多くの場合、MEMS の細管の長さは波長より小さいであろう。

18

表 4.2: 音響平行平板アクチュエータにおいて算出が必要なパラメータ
パラメータ種別
機械特性パラメータ

電気特性パラメータ

パラメータ名単位

振動板の等価質量kg

振動板のバネ定数N/m

振動板の機械抵抗

Ns/m

直流動作点での平板間容量F

直流動作点での平板間容量の 1 階微係数値

F/m

直流動作点での平板間容量の 2 階微係数値

F/m2

音響容量

m3 /Pa

音響抵抗

Ns/m5

イナータンス

kg/m4

音響特性パラメータと表される。この式において、式 (4.42), (4.43) に示す電気系-音響系の対応づけを行うと、右辺の係数は電気回路におけるインダクタに相当し、イナータンスと呼ばれる。

LA =

ρl
S

(4.60)

なお、穴あき管路コンポーネントのように、複数の穴が開いているような場合、イナータンスは上式の断面∑
積 S を断面積の合計 i Si で置き換えて計算を行う。

4.4.4

音響平行平板アクチュエータ

音響平行平板アクチュエータは音響-電気-機械の相互変換コンポーネントなので、上記の音響特性パラメータの他、電気特性パラメータや機械特性パラメータも算出する (表 4.2)。本節では、抽出ソフトにおける音響特性パラメータ以外のパラメータの計算方法を説明する。

振動板のバネ定数
音響平行平板アクチュエータの振動板は円形である。周囲を固定した円形のダイアフラムに一様な圧力∆P を加えたとき、ダイアフラムの撓み分布は、
d(r) =(a2 − r2 )2∆P
64D(4.61)と表せる [8]。ここで、ヤング率を E 、ポアソン比を ν 、ダイアフラムの厚みを t とすると、D は以下の式で定義される。

D=

Et3
12(1 − ν 2 )(4.62)

ダイアフラムの中心の変位に注目したとき、バネ定数は、k=∆P · πa2
64πD=d(r = 0)a2と表すことが出来る。

19(4.63)

図 4.19: 振動板

振動板の等価質量4.2.2 節に示したエネルギー法に基づき、等価質量の算出を行う。前節の撓み分布式を、ダイアフラム中
心での撓み dc [= d(r = 0)] を用いて書き換えると、()
r2
d = dc 1 − 2(4.64)aと掛ける。ダイアフラムが撓み方向に調和振動したとき、運動エネルギーの最大値は、最大値は、
[
]
∫
ω2 a
1 ρπa2 t 2 22
Tmax =
ρ(2πrdr)t · {d(r)} =ω dc(4.65)2 0
25と計算できる。一方、ダイアフラムの中心にマスが集中する考えた場合、運動エネルギーの最大値は、

Tmax,0 =1
Meq ω 2 d2c2(4.66)

ρπa2 t5(4.67)

したがって、ダイアフラムの等価質量は、

Meq =
と求まる。

振動板の機械抵抗
バックプレートに穴が開いているので、振動板に作用する機械抵抗としては、円状の穴が一様に開いた薄い平板の粘性抵抗の式を適用している [9]。

rf =

3µS 2
(4ϕ − ϕ2 − 2 ln ϕ − 3)
2πN G3(4.68)

ここで、S は平板の面積、µ は媒質の粘度、G は平板-基板間の距離、N は穴の総数、ϕ はバックプレートに占める穴の面積の割合である。

直流動作点での平板間容量
まず、ダイアフラム-バックプレート間を等価な平行平板とバネに置き換え、平板間に働く静電引力とバネの復元力の釣り合いから、直流バイアス印加時の平板間距離を求める。尚、平板間に働く静電引力の計算においては、バックプレートの穴の影響は無視している。求めた平板間距離を用いて、平行平板近似により、平板間容量を算出している。

20

図 4.20: 磁気駆動梁

4.5

その他のコンポネントのパラメータ算出

4.5.1

磁気駆動梁のパラメータ算出

磁気駆動梁は、図 4.20 に示すように、カンチレバー上のワンターンの配線を流れる電流と外部磁場によって生じるローレンツ力によって駆動するカンチレバーである。このカンチレバーで算出が必要なパラメー
タを表 4.3 に示す。以下の小節で、各パラメータの算出方法を述べる。
表 4.3: 磁気駆動梁において算出が必要なパラメータ
パラメータ種別

パラメータ名単位

カンチレバー部

バネ定数

N/m

等価質量

kg

機械抵抗

Ns/m

配線部

自己インダクタンス

H

電気抵抗

Ω

カンチレバーのバネ定数
磁場に垂直な配線の一部のうち、根元から遠い線分に働く (根元からの距離 x0 ) 点荷重と撓みの関係から
バネ定数を求める。ここで、磁気駆動カンチレバーは、SiO2(第 1 層目) と Si (第 2 層目) の 2 層構造であることに注意すると、バネ定数は、

k=3(E1 I1 + E2 I2 )x30

21

(4.69)

図 4.21: 磁気駆動梁と表される。ここで、Ei (i = 1, 2) は第 i 層目のヤング率、Ii (i = 1, 2) は第 i 層目の中立軸に関する断面二次モーメントである。中立軸の位置は、∫∑
E1 t21 + E2 t2 (2t1 + t2 )m Em Sm y dSm
∑
ȳN =
=
2(E1 t1 + E2 t2 )m Em Smと決まるので、中立軸に関する断面 2 次モーメントは、
∫ t1
w
3
I1 =
(y − ȳN )wdy = {(t1 − ȳN )3 − ȳN}30∫ t1 +t2w
I2 =(y − ȳN )wdy = {(t − ȳN )3 − (t1 − ȳN )3 }3
t1

(4.70)(4.71)(4.72)

から計算できる。尚、t は厚みの総計、w はカンチレバーの幅である。

カンチレバーの等価質量4.2.2 節に示したエネルギー法に基づいて求めた式により、算出している。
Meq =

33
(ρ1 t1 + ρ2 t2 )wl
140(4.73)

カンチレバーの機械抵抗

4.3.2 節に示したスクイーズフィルムダンピングの解析式により算出している。

配線の自己インダクタンス
ワンターンの配線を 1 回巻きの長方形コイルとみなして導出した自己インダクタンスの式を用いている(図 4.22)。自己インダクタンスの導出にあたっては、以下の仮定をおいている。
1. コイルの断面は半径 r の円形である。半径 r は、長方形断面と同じ断面積を有する円の半径を用いる。
2. 半径 r は、コイルの縦、横の長さに比べて十分小さい (r << a, b)。
3. 電流を流したとき、電流分布は断面で一様である。
これらの仮定のもとで求めた自己インダクタンスは、
[
{
}
]
µ0
µ0
2ab
L =
(a + b) +(a + b) ln− 2 − a ln(a + c) − b ln(b + c) + 2c + ∆ϕ
4π
π
r
[
√
√
√
1
b2 + ξ 2 + ξ
2
2
2
2
∆ϕ =
b + x − b + ξ + ξ ln
2
b
]
√
√
2 + ξ2 + ξ
2xb 2 + x2 + x √ 2r+ r + ξ 2 − x + x ln
− ξ ln
−x ln
b
r
r

22

(4.74)(4.75)

図 4.22: ワンターン配線

ここで、c =√

a2 + b2 は対角線の長さ、∆ϕ はワンターン配線の下辺の一部が欠けていることによるインダクタの減少量、ξ ≡ a − x − r である。式 (4.74) において、第 1 項は配線内部にできる磁束の寄与、第 2 項
は配線外部の磁束で磁力線がコイル面を切るものの寄与である。尚、一部金属 (鉄、コバルト、ニッケルなど) を除く多くの物質で比透磁率がほぼ 1 であるため、透磁率には真空透磁率 µ0 = 4.0 × 10−7 [H/m] を用いることとする。

配線の電気抵抗4.2.3 節に示した算出方法を用いている。

4.5.2

ユニモルフ圧電駆動梁のパラメータ算出

ユニモルフ圧電駆動梁は圧電膜を乗せたカンチレバーであり、圧電膜で電圧を印加することによって撓
む。ユニモルフ圧電駆動梁において計算が必要なパラメータは、梁の等価質量だけである。これは、磁気駆動梁と同様の方法で計算できる。具体的には、式 (4.69) の x0 をカンチレバーの長さ l で置き換えたものになる。

図 4.23: ユニモルフ圧電駆動梁

23表 4.4: 静電平行平板において算出が必要なパラメータ
パラメータ名

4.5.3

単位

平行平板の移動方向

-

平板の質量

kg

直流動作点での平板間容量F

直流動作点での平板間容量の 1 階微係数

F/m

直流動作点での平板間容量の 2 階微係数

F/m2

機械抵抗

Ns/m

静電平行平板のパラメータ算出

静電平行平板において算出が必要なパラメータを表 4.4 に示す。

• 機械抵抗：4.3.2 節に示したスクイーズフィルムダンピングの解析式により算出している。
• 直流動作点での平板間容量：
フリンジ効果を無視した容量式
C(d) = ϵwl
G0 − d

(4.76)

を用いて算出している。ここで、平板間距離の減少量 d には直流動作点解析の結果を用いている。

図 4.24: 静電平行平板の概略図

24

関連図書
[1] S. D. Senturia, Microsystem Design, Kluwer Academic Press (2001) pp.252; W. C. Young and R. G.
Budynas, Roark s Formula for Stress and Strain (seventh edition), McGraw-Hill (2002) TABLE8.10,
pp.245.
[2] J. S. Przemieniecki, Theory of Matrix Structural Analysis, McGraw-Hill (1968) pp.388-391.
[3] 小野左京『メカトロニクス時代の機械力学』培風館 (1999).
[4] W. A. Johnson and L. K. Warne, ”Electrophysics of Mechanical Comb Actuators,” J. Microelectromech. Syst. 4 (1995) pp.49-59.
[5] S. Mochizuki, T. Tsuchiya, G. Hashiguchi, T. koike and C. Mizuta, ”Semi-analytical Formula for
the Capacitance of Vertical Comb Electrodes and its Applications,” proceedings of the 25th Sensor
Symposium (2008) pp.460-463.
[6] J. -L. A. Yeh, C. Y. Hui and N. C. Tien, ”Electrostatic Model for an Asymmetric Combdrive”, J.
Microelectromech. Syst., vol.9, pp.126-135 (2000).
[7] 坂戸健一『基礎音響工学』 コロナ社 (1990) pp.16-18.
[8] W. C. Young and R. G. Budynas, ”Roark’s Formula for Stress and Strain (seventh edition)”, pp.488.
[9] M. Bao and H. Yang, ”Squeezing ﬁlm air damping in MEMS”, Sensors and Actuators A 136 (2007)
pp.3-27.

25