はじめに
最近は「メタマテリアル」というワードを良く聞きますよね。
私自身もこの分野に非常に可能性を感じており、自己研鑽も兼ねて勉強しています。
メタマテリアルや周期構造をがっつり取り組んでいる方からするともの足りない内容かと思いますが、私と同レベルの初心者向けにわかりやすく説明できればと思っています。
メタマテリアルや周期構造は分散関係を用いてバンドギャップを持つ構造を指します。
普段、振動モードや共振などのモード次元や周波数次元でメタマテリアルや周期構造を理解しようとするとうまくいきません。
メタマテリアルや周期構造を理解するためには、「波動」的な考え方が必要です。
本日はその導入部分について説明します。
格子振動
最初の導入として、「格子振動」を考えます。
ちなみに、下記文献を参考にしています。
私の所有しているものは、古いハードカバーのもので、下記のように日焼けしてしまっています。
理論的なことをつらつらと書くのは好きではないので、できるだけ理論式を用いないで説明したいと思います。
1原子1次元格子
まず、図1のような質量mとバネkが無限に連続するモデルを考えます。
図1
振動工学を勉強したことある方なら、下記のように考えるのではないでしょうか？
（Nが無限の）N自由度モデルと考えると、共振周波数fnは無限に存在するし、自由度が増えるにつれて共振周波数が無限に近づくんじゃない？
上記の考えは誤りです。
実は、共振周波数fnは無限に存在するのですが、モデルの自由度が増えてると上限の共振周波数fmaxに漸近していきます。
図2のように質量mが距離aで等間隔に並んでいるとして、そのときのxnの運動方程式を考えると、下式になります。
図2
$$
m \frac{d^2 x_n}{dt^2}= k(x_{n+1} – x_{n}) + k(x_{n-1} – x_{n})
$$
ここで、xnを波数q（=ω/v）を用いて下記とします。（vは速度、qは空間周波数とも呼ぶ）
この角振動数ωを波数qの関数を分散関係と呼びます。周期構造やメタマテリアルではこの分散関係を解いてバンドギャップを求めることになります。
$$
x_n= \frac{1}{\sqrt{m}}u(q)e^{i(qan-ωt)}
$$
上記を運動方程式に代入して、ωについて整理すると下式になります。
$$
ω= ±2\sqrt{\frac{k}{m}}sin{\frac{qa}{2}}
$$
これを図示すると図3のようになります。（m=1, k=1000）
図3
clear all;clc;close all
m=1;
k=10^3;
qa=-pi:0.1:pi;
figure
plot(qa,+2*sqrt(k/m)*sin(qa/2)/(2*pi))
hold on
plot(qa,-2*sqrt(k/m)*sin(qa/2)/(2*pi))
grid on
ylim([0 20])
xlim([-pi pi])
ylabel('共振周波数 fn')
xlabel('単位格子内の波数 qa')
2原子1次元格子
1原子1次元格子と同様の考えで、図4の2原子1次元格子を考えます。
図4
\(\left\{\begin{array}{l}m_1 \frac{d^2 S_{n,1}}{dt^2}= k(S_{n,2} – S_{n,1}) + k(S_{n+1,2} – S_{n,1})\\m_2 \frac{d^2 S_{n,2}}{dt^2}= k(S_{n-1,1} – S_{n,2}) + k(S_{n,1} – S_{n,2})\end{array}\right.\)
上記を運動方程式に代入して、ωについて整理すると下式になります。
$$
ω^2=k (\frac{1}{m_1}+\frac{1}{m_2}) ± k \sqrt{(\frac{1}{m_1}+\frac{1}{m_2})^2 – \frac{4}{m_1 m_2} sin^2 (\frac{qa}{2})}
$$
これを図示すると図5のようになります。（m1=1, m2=2, k=1000）
図5
m1とm2をkで無限につなぐと、図5に示すように、ω^2=1000～2000に波を伝搬しないバンドギャップが表れます。
周期構造やメタマテリアルでは、このバンドギャップをうまくデザインしましょうというのが根っこのコンセプトになっています。
clear all;clc;close all
m1=1;
m2=2; % m1<m2
k=10^3;
qa=-pi:0.1:pi;
w2_p=k*(1/m1+1/m2)+k*sqrt( (1/m1+1/m2)^2 - 4/(m1*m2)*(sin(qa/2)).^2 );
w2_m=k*(1/m1+1/m2)-k*sqrt( (1/m1+1/m2)^2 - 4/(m1*m2)*(sin(qa/2)).^2 );
figure
plot(qa,w2_p)
hold on
plot(qa,w2_m)
grid on
% ylim([0 20])
xlim([-pi pi])
ylabel('ω^2')
xlabel('単位格子内の波数 qa')
text(pi,2*k/m2,'2k/m2')
text(pi,2*k/m1,'2k/m1')
ちなみにqが-π/a～π/aの区間のみ扱う考え方は”第1ブリルアン・ゾーン”と呼びます。
先ほど参考文献で載せたL. Brillouinさんが発案したんですねー。
次回は2次元に拡張したいと思います。
コメント