スポンサーサイト

上記の広告は1ヶ月以上更新のないブログに表示されています。
新しい記事を書く事で広告が消せます。
スポンサー広告 | --/--/-- --:--

ドルーデ模型の計算例

各模型での誘電率や屈折率の式をイメージしたいので、
実際に計算して、グラフ化してみることにしました。

まずは、ドルーデ模型
文献 [1] の計算例にあるパラメータを使って、計算してみました。

パラメータ
$\hbar\omega_p = \sqrt{30} {\rm eV} \simeq 5.48 {\rm eV}$
$\hbar \gamma = 0.02 {\rm eV}$

具体的に何の金属を想定した値なのかはわかりませんが、
プラズマ周波数はほぼ可視~紫外域にあるようなので、
紫外の値(波長で 225 nm)になってますね。

γ の方は通常、逆数 τ がほぼ 10 ~ 100 fs 程度のようですが、
このパラメータは τ = 33 fs に対応してますね。

計算した結果はこちら。

誘電率
Drude-e1e2-hw.png

プラズマ周波数でゼロになり、
プラズマ周波数以下で実部は大きく負になっていく。

屈折率
Drude-n1n2-hw.png

高周波数極限で、実部は1に近づいていく。
虚部は小さくなっていく。
電子が追随できなくなるので、光がそのまま透過していくことに対応している。

低周波領域では、虚部が大きくなり、光が侵入できなくなる。
これは実部の誘電率が大きく負になることに起因している。

反射率
Drude-R-hw.png

プラズマ周波数より低周波でほぼ全反射になる。
逆に高周波領域では、ほとんんど反射しなくなる。

グラフにしてみると、分かりやすいですね!

参考文献
[1] 工藤恵栄 「光物性基礎」 (オーム社)
スポンサーサイト
物理>物性・固体物理 | コメント(0) | 2015/10/06 12:17

プラズマ周波数 (3) 解決編

以前にプラズマ周波数の意味が分からないと言って、
怪しげな記事を書きましたが、ようやく理解できました!
(と言っても普通の人は当然理解しているような話ですが・・・汗)

参考文献 [1] に挙げた室蘭工大の矢野先生の講義ノート
ほんとに分かりやすかったです。
学生向けに作られていて(講義ノートだからそりゃそうですが・・・)、
とにかく計算過程が懇切丁寧で素晴らしいです。
講義ノートと言っても、これほど丁寧に計算してくれるノートは
あまり見かけませんね。

どこで思考が混乱していたかというと、
外部電場にとらわれていたことです。

プラズマ周波数というのは、外部電場がない状態で、
自由電子の集合が自らの作る反電場を復元力として
振動する時の固有振動数
のことだったんですね。

以下、[1] を参考にして、理解したことをまとめてみたいと思います。

簡単のために、一次元で考える。
平衡状態(電気的に中性)から自由電子全体が x 軸正の方向に x だけシフトした状態を考える。
この時、外場がなければどのような運動をするかということを考える。

電子が x だけシフトすると、表面に分極電荷が現れる。
x 軸正負両側の表面に単位面積あたり現れる電荷は $\mp Nex$ であるから、
この表面電荷によりコンデンサが構成され、内部には\[
E = \frac{Nex}{\varepsilon_0}
\tag{1}
\]の電場が生成されるので、これが復元力となり、振動が起きる。

自由電子の運動方程式は、\[
m \ddot{x} = -eE = -\frac{Ne^2}{\varepsilon_0} x
\tag{2}
\]となり、固有振動数\[
\omega_p = \frac{Ne^2}{m\varepsilon_0}
\tag{3}
\]の振動が起きることが分かる。

これがプラズマ周波数の正体だ!

・・・って、大げさに言うほどのことでもなくて、
分かってしまえば、たったこれだけの単純な話だったんですね(^^;

この時、分極は $P = -Nex$ であるから、\[
\varepsilon = \varepsilon_0 E + P = 0
\tag{4}
\]となるから、誘電率はゼロである。

文献 [2] によると、
縦波のモードが存在しうるのは、誘電率がゼロの時のみである
ことが分かります。

マックスウェル方程式によると、\[
\nabla \cdot {\bf D} = \nabla \cdot ({\varepsilon {\bf E}}) = 0
\tag{5}
\]電場を\[
{\bf E} = {\bf E}_0 e^{i(\omega t - {\bf k}\cdot{\bf r})}
\tag{6}
\]として、εやE0 は空間的に一様であるとすると、(5) より\[
\varepsilon {\bf k} \cdot {\bf E}_0 = 0
\tag{7}
\]なる関係式が得られる。
これが成立する条件は、
${\bf k} \cdot {\bf E}_0 = 0$ (横波)であるか、または $\varepsilon = 0$ である。

まとめると、

縦波が存在するためには、誘電率がゼロでなければならない。
その時は、分極電荷による反電場を復元力とした固有振動が起きていることになる。


[1] 矢野隆治 「応用光学講義ノート」
http://www.muroran-it.ac.jp/rikagaku/ap/zyuken/labo/yano/opt_found_for_Student_20150326.pdf
[2] 佐藤勝昭 「プラズモンの基礎」
http://home.sato-gallery.com/research/principles_of_plasmons.pdf
物理>物性・固体物理 | コメント(0) | 2015/10/02 18:52

屈折率と反射率・透過率の関係

今後のために必要なので、
反射率・透過率が屈折率でどのように表されるか
を導出します(直入射の場合)。

z 軸正方向に媒質Aから媒質Bへ光が入射したとする。
z = 0 に界面があり、
z < 0 が複素屈折率 $n_a$ の媒質 A、
z > 0 が複素屈折率 $n_b$ の媒質 B とする。

入射波を i 、透過波を t、反射波を r の添え字で表すこととすると、
各々の電場は、\[
E_{ix} = E_i \exp[i(n_a kz-\omega t)] \\
E_{tx} = E_t \exp[i(n_b kz-\omega t)] \\
E_{rx} = E_r \exp[-i(n_a kz+\omega t)]
\tag{1}
\]と表せる(ただし、k は真空中の波数)。
付随する磁場は、\[
\nabla\times{\bf E} = -\frac{\partial {\bf B}}{\partial t}
\tag{2}
\]を用いて、計算できる。
E の空間微分のうち、ゼロでないのは $\partial E_x/\partial z$ のみであるから、
磁場は y 成分のみゼロでない。
これを用いて各々の磁場の y 成分を計算すると、\[
B_{iy} = \frac{n_a E_i}{c} \exp [i(n_a kz-\omega t)] \\
B_{ty} = \frac{n_b E_t}{c} \exp [i(n_b kz-\omega t)] \\
B_{ry} = -\frac{n_a E_r}{c} \exp [i(n_a kz+\omega t)]
\tag{3}
\]となる。
境界(z=0) における電場の連続性より\[
E_i + E_r = E_t
\tag{4}
\]
磁性体ではないと仮定して、境界両側で透磁率に差はないとすると、
磁場 H の連続性より\[
n_a E_i - n_a E_r = n_b E_t
\tag{5}
\]
(4) (5) 両式を解いて、振幅反射率振幅透過率を求めると、\[
r = \frac{E_r}{E_i} = \frac{n_a - n_b}{n_a + n_b}
\tag{6}
\]\[
t = \frac{E_t}{E_i} = \frac{2n_a}{n_a + n_b}
\tag{7}
\]と表せる。

え~と、ここからは我流なのですが、
まず、屈折率の虚部(つまり減衰)がないとし、n がすべて実数であると仮定すると、
光強度は $n|E|^2$ に比例するから、強度反射率強度透過率は\[
R = |r|^2 = \frac{(n_a - n_b)^2}{(n_a + n_b)^2}
\tag{8}
\]\[
T = \frac{n_b}{n_a}|t|^2 = \frac{4n_a n_b}{(n_a + n_b)^2}
\tag{9}
\]となり、\[
R + T = 1
\tag{10}
\]というエネルギー保存則が導かれる。
虚部が存在すると、どうも計算がうまくいきません。
減衰があっても、境界面上では成立しているように思うのですが・・・

よくある系として、空気から媒質へ入射する場合を考えておきます。
媒質の複素屈折率を n として、$n_a = 1$、$n_b = n$ であるから、\[
r = \frac{1 - n}{1 + n}
\tag{11}
\]\[
t = \frac{2}{1 + n}
\tag{12}
\]となる。
この場合、反射に関しては、減衰がないので、強度反射率を考えると、\[
R = |r|^2 = \frac{(n_1 - 1)^2 + n_2^2}{(n_1 + 1)^2 + n_2^2}
\tag{13}
\]となる。(ここで、$n_1$、$n_2$ は屈折率の実部と虚部)

参考文献
[1] 花村榮一 「固体物理学」(裳華房)
ジャンル:[学問・文化・芸術]  テーマ:[自然科学
物理>物性・固体物理 | コメント(0) | 2015/10/01 08:18

ローレンツ模型

ドルーデ模型に引き続き、ローレンツ模型(Lorentz model)を扱います。

ドルーデ模型は、自由電子の散乱のみを考慮して、
金属の分散をうまく説明するモデルであるが、
ローレンツ模型は、さらにバンド間遷移を考慮して、
復元力のあるバネに束縛された電子の振動子として扱うモデルで、
誘電体の分散をうまく説明できる。

復元力の項を取り入れると、運動方程式は\[
m\ddot{x} + m\gamma \dot{x} + m\omega_0^2 x = -eE
\tag{1}
\]となる。これをいつものように振動解として解くと、\[
x_0 = \frac{eE_0}{m} (\omega^2 - \omega_0^2 - i\gamma\omega)^{-1}
\tag{2}
\]となり、振動子の密度を N とすると、分極は\[
P_0 = -Nex_0 = -\frac{Ne^2E_0}{m} (\omega^2 - \omega_0^2 - i\gamma\omega)^{-1}
\tag{3}
\]となる。これより、誘電率と比誘電率は\[
\varepsilon(\omega) = \varepsilon_0 - \frac{Ne^2}{m} (\omega^2 - \omega_0^2 - i\gamma\omega)^{-1}
\tag{4}
\]\[
\varepsilon_r(\omega) = 1 - \frac{Ne^2}{m\varepsilon_0} (\omega^2 - \omega_0^2 - i\gamma\omega)^{-1}
\tag{5}
\]と表される。
実部と虚部はそれぞれ、\[
\varepsilon_{1r} = 1 + \frac{Ne^2}{m\varepsilon_0}
\frac{\omega_0^2 - \omega^2}{(\omega^2 - \omega_0^2)^2 + \gamma^2\omega^2}
\tag{6.1}
\]\[
\varepsilon_{2r} = \frac{Ne^2}{m\varepsilon_0}
\frac{\gamma\omega}{(\omega^2 - \omega_0^2)^2 + \gamma^2\omega^2}
\tag{6.2}
\]となる。

参考文献
[1] 小林浩一 「光物性入門」 (裳華房)
[2] 佐藤勝昭 「プラズモンの基礎」
http://home.sato-gallery.com/research/principles_of_plasmons.pdf
物理>物性・固体物理 | コメント(0) | 2015/09/29 07:22

ドルーデ模型

一週間、夏休みでした。
勉強もさぞかし進んだだろうと思いきや、
最近入会したスポーツジムとピアノ三昧で勉強はほとんどゼロです^^;
唯一、都内に用事があった時の往復の電車で代数の本を数行読んだぐらい(笑)
勉強は、他にやることがない環境に置かれないと難しいですね・・・汗

さて、プラズマ周波数の理解は置いておいて、先に進みます。

自由電子模型に続いて、今度はドルーデ模型 (Drude model) を扱います。

金属の状態をうまく説明するモデルのようです。

今度は、完全に自由な電子ではなく、衝突散乱を受ける場合を考える。
散乱の部分をダンピング項として取り入れて、\[
m \ddot{x} + m\gamma \dot{x} = -eE
\tag{1}
\]とする。ここで、取り入れたダンピング項の係数\[
\gamma = 1/\tau
\tag{2}
\]は、電子同士や電子とイオンの衝突に関する平均自由時間の逆数と思えばよい。

前記事と同様に、ωで振動する成分を考えることにすると、振幅は\[
x_0 = \frac{eE_0}{m\omega(\omega - i\gamma)}
\tag{3}
\]となり、分極は、\[
P_0 = -Nex_0 = -\frac{Ne^2E_0}{m\omega(\omega-i\gamma)}
\tag{4}
\]となる。これより、誘電率と比誘電率は、\[
\varepsilon(\omega) = \varepsilon_0 - \frac{Ne^2}{m\omega(\omega-i\gamma)}
\tag{5}
\]\[
\varepsilon_r(\omega) = 1 - \frac{Ne^2}{m\varepsilon_0\omega(\omega-i\gamma)}
= 1 - \frac{\omega_p^2}{\omega(\omega-i\gamma)}
\tag{6}
\]と表される。

実部と虚部はそれぞれ、\[
\varepsilon_{1r}(\omega) = 1 - \frac{\omega_p^2}{\omega^2 + \gamma^2}
\tag{7.1}
\]\[
\varepsilon_{2r}(\omega) = \frac{\omega_p^2 \gamma} {\omega(\omega^2 + \gamma^2)}
\tag{7.2}
\]となる。

低周波極限(ω→0)で、実部は負の無限大に、虚部は正の無限大にそれぞれ発散する。

誘電率が負になると、屈折率の虚部が大きくなり、
光は内部に入り込めなくなり、強い反射が起きる。


これが金属が光を反射するメカニズムのようですが、
前述のとおり、誘電率のことがよくわかってないので、
数式上は理解できても、いまいちイメージがつかめませんね。

参考文献
[1] 小林浩一 「光物性入門」 (裳華房)
[2] 佐藤勝昭 「プラズモンの基礎」
http://home.sato-gallery.com/research/principles_of_plasmons.pdf
ジャンル:[学問・文化・芸術]  テーマ:[自然科学
物理>物性・固体物理 | コメント(0) | 2015/08/31 08:00

プラズマ周波数 (2)

プラズマ周波数についてまだ理解できないままですが、
とりあえず、まったくの我流で考察をすすめてみます。

プラズマ周波数\[
\omega_p^2 \equiv \frac{Ne^2}{m\varepsilon_0}
\tag{1}
\]を用いて、自由電子模型の記事の (5) 式を書き直すと、\[
P_0 = -\varepsilon_0 E_0 \frac{\omega_p^2}{\omega^2}
\tag{2}
\]となる。

まず、低周波(DC)極限 ω→0 を考えると、分極は無限大に発散する。

しかし、分極は有限のはずだから、これはむしろ、
内部の電場が限りなくゼロになっていると考えるべきでしょう。

では、どうして、ゼロになるのか?

完全なDC電場の場合は、
電子はどんどん一方向に動いて行き、分極が作られ、
その分極による逆向きの電場により、外場が打ち消されて、
やがて内部電場はゼロになる。
いわゆる静電遮蔽の効果です。

しかし、振動電場の場合は事情が変わってきます。
(2) 式にマイナスがついていることから分かる通り、
電場の振動と分極の振動は逆位相になっています。

これは、振動する外力による強制振動の定常解を考えると、
確かにそうなって当然ですね。
なぜなら、変位が正の領域に入っている時には既に
減速フェーズに入っていなければならないので、
外力の向きは負の向きになっているはずだから。

ということは、静電遮蔽の場合と異なり、
分極電場の向きは外部電場と同じ方向となり、
打ち消しあうのではなく、増強されることになります。

ただし、ここで見えている電場は外部電場と分極電場の
足しあわされた後の内部電場であり、
分極が応答するのはこの内部電場に対してだから、
内部電場と分極電場の位相が同位相というだけで、
実際に外部電場の位相との比較でどうなっているかは分かりません。

実際のところ、電子の振動振幅が大きくなってくると、
電場は電子に仕事をしているわけだから、
電場のエネルギーは失われていって、弱まっていく気はします。
ただ、電磁気学としてどのような機構でゼロになるのか
よく分かりません。

ところで、光物性の本 [1] によると、
誘電率は、入射光の電場と誘起された双極子からの双極子放射光電場の足し合わせ
考えるとよいようです。
[1] の第一章の初めに、非常にシンプルで教育的な計算例が載っていて、
一度読んでなるほど・・・と思った記憶はあるのですが、
今回はめんどくさくて割愛してしまいました。
やはり、もう一度おさらいした方がよさそうですね^^;

この [1] の計算例は、ローレンツ模型について行われているのですが、
自由電子模型についても、この考え方を適用して、
分極からの双極子放射を考えれば、
おそらく、双極子放射電場と入射電場がうまく打ち消しあうように働いている
という結果が出るのではないかと考えています。

さて、次に、高周波極限ω→∞を考えると、
これは簡単で、振動周期が速くなるにつれ、
電子が一方向に加速される時間が短くなるため、
電子の移動振幅が小さくなり、分極が小さくなっていく。
そのため、入射光に与える影響は少なくなっていき、
極限的には、真空中と同様に光が透過していくようになる。

このイメージはおそらく正しいでしょう。

この間に位置する臨界点値的なものがプラズマ周波数だと思うのですが、
その時にたとえば何と何がバランスした状態になっているかとか
そういうイメージが描き出せません。

分かる方、是非、ご教授ください。
非常にシンプルなモデルなのに、考え出すと、難解ですね!


参考文献
[1] 小林浩一 「光物性入門」 (裳華房)
ジャンル:[学問・文化・芸術]  テーマ:[自然科学
物理>物性・固体物理 | コメント(0) | 2015/08/17 12:54

プラズマ周波数 (1)

先日、少しの間、「プラズマ周波数」についての記事をアップしましたが、
その後、恥ずかしながら、とんでもない勘違いをしていたことが分かり、
取り下げました。

結局、プラズマ周波数とは何なのか?
誘電率がゼロとか負とは、どういうことなのか?

さっぱりイメージできなくなりました・・・

まあ、ここは理解できていなくても、
とりあえず、数式だけで先へ進むことはできるのですが、
それもどうかと・・・
物理>物性・固体物理 | コメント(0) | 2015/08/17 00:06

自由電子模型

いよいよ、物質中の電子の運動を古典的にモデル化して、
誘電率の分散を求めていきます。
普通はドルーデ模型から入るのですが、その前に、
完全な自由電子の模型を考えてみようと思います。

簡単のために、一次元で考える。
衝突散乱もなく、復元力もない完全に自由な電子の運動方程式は、\[
m\ddot{x} = -eE
\tag{1}
\]となる。振動数ωの外場を与えると、x が振動数ωで強制振動されるから、\[
x = x_0 e^{i\omega t}
\tag{2.1}
\]\[
E = E_0 e^{i\omega t}
\tag{2.2}
\]を代入すると、\[
x_0 = \frac{eE_0}{m\omega^2}
\tag{3}
\]となる。
先日、議論したように、巨視的分極 P は数密度 N とすると、\[
P = Np = -Nex
\tag{4}
\]と表されるから、分極のω振動成分は\[
P_0 = -Nex_0 = - \frac{Ne^2E_0}{m\omega^2}
\tag{5}
\]となる。分極と誘電率の関係\[
D = \varepsilon_0 E + P = \varepsilon E
\tag{6}
\]から、\[
\varepsilon(\omega) = \varepsilon_0 - \frac{Ne^2}{m\omega^2}
\tag{7}
\]\[
\varepsilon_r(\omega) = 1 - \frac{Ne^2}{m\varepsilon_0\omega^2}
\tag{8}
\]という誘電率の分散式が得られる。\[
\omega_p^2 \equiv \frac{Ne^2}{m\varepsilon_0}
\tag{9}
\]と定義すると、\[
\varepsilon_r(\omega) = 1 - \frac{\omega_p^2}{\omega^2}
\tag{10}
\]と表される。$\omega_p$ をプラズマ周波数と呼ぶ。


参考文献
[1] 小林浩一 「光物性入門」 (裳華房)
[2] 佐藤勝昭 「プラズモンの基礎」
http://home.sato-gallery.com/research/principles_of_plasmons.pdf
[3] 矢野隆治 「応用光学講義ノート」
http://www.muroran-it.ac.jp/rikagaku/ap/zyuken/labo/yano/opt_found_for_Student_20150326.pdf
[4] 花村栄一 「固体物理学」 (裳華房)
物理>物性・固体物理 | コメント(0) | 2015/08/13 12:55

屈折率の虚部の意味

波動方程式の z 方向に進行する平面波解 ( v = c/n) \[
{\bf E} = {\bf E}_0 \exp \{ i\omega (nz/c - t) \}
\tag{1}
\]に複素屈折率\[
n = n_1 + i n_2
\tag{2}
\]を代入してみる。
\[
{\bf E} = {\bf E}_0 \exp (-\omega n_2z/c) \exp \{ i\omega (n_1z/c - t) \}
\tag{3}
\]電磁波強度(単位時間に単位面積を通過するエネルギー)のサイクル平均を I とすると、
I は $|{\bf E}|^2$ に比例するから、\[
I = I_0 \exp (-2\omega n_2z/c)
\tag{4}
\]と表される。
これを見ると、屈折率の虚部は吸収を表すことが分かる。
吸収係数 αを\[
I = I_0 e^{-\alpha z}
\tag{5}
\]で定義すると、\[
\alpha = 2\omega n_2 / c
\tag{6}
\]となる。


参考文献
[1] 小林浩一「光物性入門」(裳華房)
ジャンル:[学問・文化・芸術]  テーマ:[自然科学
物理>物性・固体物理 | コメント(0) | 2015/08/12 12:23

誘電率と屈折率

このところ、誘電分極のことを理解しようとしていたのは、
光物性やプラズモニクスについて勉強したくて、
基礎的なところを整理したかったというのが動機です。

まずは、誘電率屈折率の関係について、整理しておきます。
単位系は引き続き、SI です。

電場がそれほど強くないときは、分極はおおむね電場に比例すると考えてよい。\[
{\bf P} = \varepsilon_0 \chi {\bf E}
\tag{1}
\]この時の比例係数 $\chi$ を電気感受率と呼ぶ。
電束密度は、\[
{\bf D} = \varepsilon_0 {\bf E} + {\bf P} = \varepsilon_0 (1 + \chi) {\bf E}
\tag{2}
\]で表されるから、\[
\varepsilon = \varepsilon_0 (1 + \chi)
\tag{3}
\]を誘電率と定義すると、\[
{\bf D} = \varepsilon {\bf E}
\tag{4}
\]と表せる。
また、真空の誘電率との比\[
\varepsilon_r = \varepsilon / \varepsilon_0
\tag{5}
\]を比誘電率と呼ぶ。

ここまでが定義。
マックスウェル方程式では、$\varepsilon_0$ が $\varepsilon$ に置き変わるので、
そこから得られる電磁波の速度を v とすると、\[
v^2 = \frac{1}{\varepsilon \mu_0} = \frac{c^2}{\varepsilon_r}
\tag{6}
\]となる。
一方、屈折率を n とすると、\[
v = c / n
\tag{7}
\]であるから、比誘電率と屈折率の間には、\[
n^2 = \varepsilon_r
\tag{8}
\]の関係が導かれる。

(4) で位相遅れまでを考慮すると、誘電率は複素数になる。
それに応じて、屈折率も複素数となる。
これらを実部と虚部に分けて、\[
\varepsilon_r = \varepsilon_{1r} + i\varepsilon_{2r}
\tag{9}
\]\[
\varepsilon = \varepsilon_1 + i\varepsilon_2
\tag{10}
\]\[
n = n_1 + i n_2
\tag{11}
\]と書くことにする。
(8) に代入して、実部と虚部を等置すると、\[
\varepsilon_{1r} = n_1^2 - n_2^2
\tag{12.1}
\]\[
\varepsilon_{2r} = 2n_1n_2
\tag{12.2}
\]を得る。

逆に、n について解く。\[
n_1^2 + n_2^2 = \sqrt{(n_1^2 - n_2^2)^2 + 4n_1^2n_2^2} = \sqrt{\varepsilon_{1r}^2 + \varepsilon_{2r}^2}
\tag{13.1}
\]\[
n_1^2 n_2^2 = \varepsilon_{2r}^2 / 4
\tag{13.2}
\]であるから、$n_1^2$, $n_2^2$ は2次方程式\[
x^2 - \sqrt{\varepsilon_{1r}^2 + \varepsilon_{2r}^2} x + \varepsilon_{2r}^2 / 4 = 0
\tag{14}
\]の解であり、\[
n_1^2 = \frac{1}{2} \left( \varepsilon_{1r} + \sqrt{\varepsilon_{1r}^2 + \varepsilon_{2r}^2} \right)
\tag{15.1}
\]\[
n_2^2 = \frac{1}{2} \left( -\varepsilon_{1r} + \sqrt{\varepsilon_{1r}^2 + \varepsilon_{2r}^2} \right)
\tag{15.2}
\]と解くことができる。

参考文献
[1] 小林浩一「光物性入門」(裳華房)
ジャンル:[学問・文化・芸術]  テーマ:[自然科学
物理>物性・固体物理 | コメント(0) | 2015/08/11 12:48
 | HOME | Next »

上記広告は1ヶ月以上更新のないブログに表示されています。新しい記事を書くことで広告を消せます。