デバイ緩和の式からHavriliak-Negamiの式まで

背景

誘電体は、絶縁体の一種であり、印加される電圧の周波数に応じて電気的特性が変わる。電子回路のパーツではコンデンサに使用されていたりする。どの程度の周波周波数の電圧で、どの程度、電気的特性が変動するかは物体に依存する。物体内部の特性を議論する際には誘電率という物性値を計測して議論される。単純に回路素子の特性として議論したい場合はインピーダンスもしくは容量と表現するが、試料寸法で割れば誘電率が算出される。

その誘電率でFittingに採用されがちなデバイ緩和の式、Havriliak-Negamiの式を順に紹介する。

デバイ緩和

誘電体の中でも有機物の高分子で構成される誘電体は、単純なモデルとしてデバイ緩和というモデルが使用される。デバイ緩和の式は以下の通りである。

$$
\varepsilon^* = \varepsilon_{\infty}+ \frac{\varepsilon_s-\varepsilon_\infty}{\left(1+i\omega\tau\right)^\beta}
$$

上式の通り誘電率$\varepsilon^*$は複素数の形で表現さる。また、$\varepsilon_\infty$は高周波における誘電率であり、$\varepsilon_s$は低周波の誘電率である。よって第二項は $\varepsilon_s-\varepsilon_\infty$は高周波と低周波の誘電率の差で決まる。$\omega$はrad/sの周波数であり、$\tau$は緩和時間である。グラフ化すると以下の通りである。

import numpy as np
import matplotlib.pyplot as plt
import math
# パラメータ
epsilon_inf = 4.0  # 高周波領域での誘電率
epsilon_s = 10.0   # 低周波領域での誘電率

n = 1              # 緩和の指数

f_m=100
tau = 1/(2*math.pi*f_m)  #緩和時間

# 周波数の範囲を設定
f =  np.logspace(-1, 6, 10000) # [Hz]
omega = 2*math.pi*f


# デバイ緩和の式
epsilon = epsilon_inf + (epsilon_s - epsilon_inf) / (1 + (1j * omega * tau)**n)

# 実部と虚部を取得
epsilon_real = epsilon.real
epsilon_imag = -1*epsilon.imag
# print(epsilon_imag)
# プロット
plt.figure(figsize=(8, 6))
plt.plot(f, epsilon_real, label='Real part')
plt.plot(f, epsilon_imag, label='Imaginary part')
plt.xscale('log')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Dielectric Constant')
plt.legend()
plt.grid(True)
plt.show()

>>>出力結果

Havriliak-Negamiの式

実際の物質では、より複雑な誘電率スペクトルが現れる。デバイ緩和の式であると対応できないためデバイ緩和の式を発展させた、以下のHavriliak-Negamiの式を用いる。

$$
\varepsilon^* = \varepsilon_{\infty}+ \frac{\varepsilon_s-\varepsilon_\infty}{\left(1+(i\omega\tau\right)^\alpha)^\beta}
$$
pythonで実装すると以下の通りである。今回はパラメータの大半を1にしているため、デバイ緩和の式と一致する。

import numpy as np
import matplotlib.pyplot as plt
import math
# パラメータ
epsilon_inf = 4.0  # 高周波領域での誘電率
epsilon_s = 10.0   # 静電誘電率と無限周波誘電率の差

f_m=100
tau = 1/(2*math.pi*f_m)  #緩和時間

beta = 1              # 緩和の指数
alpha = 1

# 周波数の範囲を設定

f =  np.logspace(-1, 6, 10000) # [Hz]
omega = 2*math.pi*f

# HavriliakとNegamiの 式
epsilon = epsilon_inf + (epsilon_s - epsilon_inf) / (1 + (1j * omega * tau)**beta)**alpha

# 実部と虚部を取得
epsilon_real = epsilon.real
epsilon_imag = -1*epsilon.imag

# プロット
plt.figure(figsize=(8, 6))
plt.plot(f, epsilon_real, label='Real part')
plt.plot(f, epsilon_imag, label='Imaginary part')
plt.xscale('log')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Dielectric Constant')
plt.legend()
plt.grid(True)
plt.show()

>>>出力結果

複数個の誘電緩和が確認される場合の式

実際の物質では、複数の誘電緩和が生じ、計測結果に現れる。
2個のHavriliak-Negamiの式が重なっている様子は以下の通りである。実際の物質と表現していますが、かなり誇大表現です。

import numpy as np
import matplotlib.pyplot as plt
import math
# パラメータ
epsilon_inf = 4.0  # 高周波領域での誘電率

delta_1 = 1
delta_2 =1

f_1=10
tau_1 = 1/(2*math.pi*f_1)        # 緩和時間

f_2=1000
tau_2 = 1/(2*math.pi*f_2)       # 緩和時間

beta_1 = 1              # 緩和の指数
alpha_1 = 1

beta_2 = 1              # 緩和の指数
alpha_2 = 1

# 周波数の範囲を設定
f =  np.logspace(-1, 6, 10000) # [Hz]
omega = 2*math.pi*f

# HavriliakとNegamiの 式
epsilon_1 = (delta_1) / (1 + (1j * omega * tau_1)**beta_1)**alpha_1
epsilon_2 = (delta_2) / (1 + (1j * omega * tau_2)**beta_2)**alpha_2

epsilon = epsilon_inf + epsilon_1 + epsilon_2
# 実部と虚部を取得
epsilon_real = epsilon.real
epsilon_imag = -1*epsilon.imag

# プロット
plt.figure(figsize=(8, 6))
plt.plot(f, epsilon_real, label='Real part')
plt.plot(f, epsilon_imag, label='Imaginary part')
plt.xscale('log')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Dielectric Constant')
plt.legend()
plt.grid(True)
plt.show()

>>>出力結果

まとめ

誘電率のfitting式として使われるデバイ緩和とHavriliak-Negamiの式を紹介した。
最後にHavriliak-Negamiの緩和が2つ現れる場合のコードを示した。
物理系、材料系の学生はこの式に最小2乗法でも実装するなり、手打ちでパラメータの数値を変えて実験データに対してFittingしてみてほしい。