第4章: 電子-フォノン結合

レベル: 中級読了時間: 約40分 Version 1.0

学習目標

電子-フォノン相互作用ハミルトニアンの理解と導出
フレーリッヒハミルトニアンと変形ポテンシャル結合の理解
ポーラロン概念（大きなポーラロン・小さなポーラロン）の習得
電子自己エネルギーと質量繰り込みの計算方法
Eliashberg関数α²F(ω)と電子-フォノン結合定数λの理解
超伝導理論（BCS理論）への応用とMcMillan方程式
第一原理計算による電子-フォノン結合の計算手法

1. 電子-フォノン相互作用ハミルトニアン

1.1 基本的な枠組み

電子-フォノン相互作用は、結晶中の電子が格子振動（フォノン）と結合する現象を記述します。全ハミルトニアンは次のように書けます：

$$\hat{H} = \hat{H}_\text{el} + \hat{H}_\text{ph} + \hat{H}_\text{el-ph}$$

ここで：

$\hat{H}_\text{el}$: 電子系のハミルトニアン
$\hat{H}_\text{ph}$: フォノン系のハミルトニアン
$\hat{H}_\text{el-ph}$: 電子-フォノン相互作用ハミルトニアン

1.2 電子ハミルトニアン

電子系は通常、Bloch状態で記述されます：

$$\hat{H}_\text{el} = \sum_{n\mathbf{k}} \varepsilon_{n\mathbf{k}} c_{n\mathbf{k}}^\dagger c_{n\mathbf{k}}$$

ここで、$c_{n\mathbf{k}}^\dagger$、$c_{n\mathbf{k}}$はバンド$n$、波数$\mathbf{k}$の電子の生成・消滅演算子、 $\varepsilon_{n\mathbf{k}}$はバンドエネルギーです。

1.3 フォノンハミルトニアン

調和近似でのフォノン系は：

$$\hat{H}_\text{ph} = \sum_{\mathbf{q}\nu} \hbar\omega_{\mathbf{q}\nu} \left(b_{\mathbf{q}\nu}^\dagger b_{\mathbf{q}\nu} + \frac{1}{2}\right)$$

$b_{\mathbf{q}\nu}^\dagger$、$b_{\mathbf{q}\nu}$は波数$\mathbf{q}$、ブランチ$\nu$のフォノンの生成・消滅演算子です。

1.4 電子-フォノン相互作用ハミルトニアン

電子-フォノン相互作用は、原子の変位により電子ポテンシャルが変化することで生じます。最も一般的な形は：

$$\hat{H}_\text{el-ph} = \sum_{nm\mathbf{k}\mathbf{q}\nu} g_{nm\nu}(\mathbf{k},\mathbf{q}) c_{n\mathbf{k}+\mathbf{q}}^\dagger c_{m\mathbf{k}} (b_{\mathbf{q}\nu} + b_{-\mathbf{q}\nu}^\dagger)$$

ここで、$g_{nm\nu}(\mathbf{k},\mathbf{q})$は電子-フォノン結合行列要素（matrix element）です。この項は電子が波数$\mathbf{k}$から$\mathbf{k}+\mathbf{q}$に遷移する際に、フォノン$\mathbf{q}\nu$を吸収または放出することを表します。

物理的解釈

電子-フォノン相互作用の物理的意味：

散乱過程: 電子がフォノンによって散乱される
運動量保存: $\mathbf{k}' = \mathbf{k} + \mathbf{q}$（または$\mathbf{k} - \mathbf{q}$）
エネルギー交換: 電子とフォノン間でエネルギーが交換される
抵抗の起源: 電気抵抗や熱抵抗の微視的機構

1.5 結合行列要素の導出

結合行列要素は、原子変位によるポテンシャルの変化から導出されます：

$$g_{nm\nu}(\mathbf{k},\mathbf{q}) = \sqrt{\frac{\hbar}{2N M \omega_{\mathbf{q}\nu}}} \langle n\mathbf{k}+\mathbf{q}| \frac{\partial V}{\partial u_{\mathbf{q}\nu}} |m\mathbf{k}\rangle$$

ここで：

$N$: 単位胞の数
$M$: 原子質量
$V$: 電子に作用するポテンシャル
$u_{\mathbf{q}\nu}$: フォノンモード$\mathbf{q}\nu$の変位

2. フレーリッヒハミルトニアンと変形ポテンシャル

2.1 フレーリッヒハミルトニアン（極性結合）

極性結晶（イオン結晶）では、縦光学（LO）フォノンによる長距離クーロン相互作用が重要です。フレーリッヒ（Fröhlich）ハミルトニアンは：

$$\hat{H}_\text{Fr} = \sum_{\mathbf{k}\mathbf{q}} V_\text{Fr}(\mathbf{q}) c_{\mathbf{k}+\mathbf{q}}^\dagger c_{\mathbf{k}} (b_{\mathbf{q}} + b_{-\mathbf{q}}^\dagger)$$

結合の大きさは：

$$V_\text{Fr}(\mathbf{q}) = -i \left(\frac{2\pi e^2 \hbar\omega_\text{LO}}{V}\right)^{1/2} \left(\frac{1}{\varepsilon_\infty} - \frac{1}{\varepsilon_0}\right)^{1/2} \frac{1}{q}$$

ここで：

$\varepsilon_\infty$: 高周波誘電率（光学的誘電率）
$\varepsilon_0$: 静的誘電率
$\omega_\text{LO}$: LO フォノン周波数
$e$: 電子電荷

特徴的な性質

$1/q$依存性: 長距離相互作用（クーロン的）
極性依存性: $1/\varepsilon_\infty - 1/\varepsilon_0$に比例
LOフォノン: 縦光学モードのみが寄与
無次元結合定数: $\alpha = e^2/(2\hbar\omega_\text{LO}a_B)(\varepsilon_\infty^{-1} - \varepsilon_0^{-1})$

2.2 変形ポテンシャル結合

共有結合性結晶では、格子変形によるバンド端のシフトが主要な結合機構です。変形ポテンシャル（deformation potential）理論では：

$$\hat{H}_\text{DP} = \sum_{\mathbf{k}\mathbf{q}} D_\text{ac} \sqrt{\frac{\hbar}{2\rho V \omega_{\mathbf{q}}}} q \, c_{\mathbf{k}+\mathbf{q}}^\dagger c_{\mathbf{k}} (b_{\mathbf{q}} + b_{-\mathbf{q}}^\dagger)$$

$D_\text{ac}$は変形ポテンシャル定数で、バンド端のひずみ依存性：

$$D_\text{ac} = \frac{\partial E_\text{band}}{\partial \epsilon}$$

ここで、$\epsilon$は体積ひずみです。

2.3 両者の比較

特性	フレーリッヒ結合	変形ポテンシャル
結晶タイプ	極性（イオン性）	非極性（共有結合性）
相互作用範囲	長距離（$1/q$）	短距離（$q$に比例）
主要フォノン	縦光学（LO）	音響フォノン
結合の強さ	小$q$で強い	大$q$で強い
典型例	GaAs, ZnO, LiF	Si, Ge, Diamond
無次元結合	$\alpha = e^2/(2\hbar\omega a_B)\Delta(1/\varepsilon)$	$\lambda = D^2/(2\rho v_s^2 E_F)$

2.4 実際の材料での例

材料パラメータの比較（Python）

import numpy as np
import matplotlib.pyplot as plt

# 材料パラメータ
materials = {
    'GaAs': {
        'type': 'polar',
        'eps_inf': 10.9,
        'eps_0': 12.9,
        'omega_LO': 36.2e12,  # rad/s
        'a_B': 10e-9,  # m (effective Bohr radius)
    },
    'Si': {
        'type': 'nonpolar',
        'D_ac': 9.0e9,  # Pa (deformation potential)
        'rho': 2330,  # kg/m^3
        'v_s': 9000,  # m/s (sound velocity)
        'E_F': 5.0 * 1.6e-19,  # J
    }
}

# フレーリッヒ結合定数
def frohlich_alpha(mat):
    e = 1.6e-19  # C
    hbar = 1.054e-34  # J·s
    eps_0_SI = 8.854e-12  # F/m

    omega = mat['omega_LO']
    a_B = mat['a_B']
    delta_eps_inv = 1/mat['eps_inf'] - 1/mat['eps_0']

    alpha = (e**2 / (4*np.pi*eps_0_SI)) / (2*hbar*omega*a_B) * delta_eps_inv
    return alpha

# 変形ポテンシャル結合定数
def deformation_lambda(mat):
    D = mat['D_ac']
    rho = mat['rho']
    v_s = mat['v_s']
    E_F = mat['E_F']

    lam = D**2 / (2 * rho * v_s**2 * E_F)
    return lam

# 計算と表示
print("電子-フォノン結合定数の比較\n")
print(f"GaAs (極性結晶):")
alpha_GaAs = frohlich_alpha(materials['GaAs'])
print(f"  フレーリッヒ α = {alpha_GaAs:.3f}")
print(f"  結合の強さ: {'中程度' if alpha_GaAs < 6 else '強い'}\n")

print(f"Si (非極性結晶):")
lambda_Si = deformation_lambda(materials['Si'])
print(f"  変形ポテンシャル λ = {lambda_Si:.3f}")
print(f"  結合の強さ: {'弱い' if lambda_Si < 0.3 else '中程度'}\n")

# 結合の波数依存性をプロット
q_values = np.linspace(1e8, 1e10, 100)  # m^-1

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))

# フレーリッヒ結合（1/q依存性）
V_frohlich = 1 / q_values  # 相対値
ax1.plot(q_values/1e9, V_frohlich/V_frohlich[0], 'b-', linewidth=2, label='Fröhlich (極性)')
ax1.set_xlabel('波数 q (nm⁻¹)', fontsize=12)
ax1.set_ylabel('相対結合強度', fontsize=12)
ax1.set_title('フレーリッヒ結合（極性結晶）', fontsize=13, fontweight='bold')
ax1.grid(True, alpha=0.3)
ax1.legend(fontsize=11)
ax1.set_yscale('log')

# 変形ポテンシャル結合（q依存性）
V_deformation = q_values  # 相対値
ax2.plot(q_values/1e9, V_deformation/V_deformation[-1], 'r-', linewidth=2, label='変形ポテンシャル（非極性）')
ax2.set_xlabel('波数 q (nm⁻¹)', fontsize=12)
ax2.set_ylabel('相対結合強度', fontsize=12)
ax2.set_title('変形ポテンシャル結合（非極性結晶）', fontsize=13, fontweight='bold')
ax2.grid(True, alpha=0.3)
ax2.legend(fontsize=11)

plt.tight_layout()
plt.savefig('coupling_comparison.png', dpi=300, bbox_inches='tight')
plt.show()

print("プロット保存: coupling_comparison.png")

3. ポーラロンの物理

3.1 ポーラロンとは

ポーラロン（polaron）は、電子-フォノン相互作用により、電子が周囲の格子変形を伴って運動する準粒子です。電子が格子を極性化（polarize）することからこの名前がつきました。

graph LR A[自由電子] -->|電子-フォノン結合| B[ポーラロン] B --> C[電子 + 格子変形雲] C --> D[有効質量増加] C --> E[移動度低下] C --> F[バンドギャップ繰り込み] style A fill:#e3f2fd style B fill:#fff9c4 style C fill:#f3e5f5

3.2 大きなポーラロン（Large Polaron）

格子変形が多数の単位胞に広がる場合、大きなポーラロンが形成されます。条件：$\alpha < 6$（フレーリッヒ結合定数が小さい）

$$E_\text{polaron} = -\alpha \hbar\omega_\text{LO} \quad \text{（結合エネルギー）}$$

$$m^* = m_\text{band} \left(1 + \frac{\alpha}{6}\right) \quad \text{（有効質量）}$$

特徴：

格子変形は広範囲に分布（数十～数百Å）
摂動論的取り扱いが可能
有効質量の増加は比較的小さい
バンド的な輸送特性を保持

3.3 小さなポーラロン（Small Polaron）

強い電子-フォノン結合（$\alpha > 6$）では、電子が1つの原子サイトに局在化し、小さなポーラロンを形成します。

$$E_\text{polaron} \approx -\frac{g^2}{\hbar\omega} \quad \text{（strong coupling limit）}$$

ホッピング輸送：

$$\mu = \frac{ea^2\nu}{k_B T} \exp\left(-\frac{E_a}{k_B T}\right)$$

ここで、$E_a \approx \hbar\omega/2$は活性化エネルギーです。

特徴：

格子変形は1～数単位胞に局在
電子は局在状態間をホッピング
活性化型の輸送（温度とともに移動度増加）
非摂動的取り扱いが必要

3.4 ポーラロン転移

性質	大きなポーラロン	小さなポーラロン
結合定数	$\alpha < 6$	$\alpha > 6$
サイズ	多数の単位胞	1～数単位胞
有効質量	$m^* \approx m(1 + \alpha/6)$	$m^* \gg m$（発散的）
輸送機構	バンド輸送	ホッピング輸送
温度依存性	$\mu \propto T^{-n}$	$\mu \propto \exp(-E_a/k_B T)$
理論手法	摂動論	変分法、経路積分

注意：中間領域

$\alpha \sim 6$付近では、大きなポーラロンから小さなポーラロンへの転移が起こります。この領域では：

摂動論も変分法も不十分
数値的手法（量子モンテカルロなど）が必要
輸送特性が複雑
温度依存性が非単調になることも

3.5 典型的な材料例

ポーラロン材料の分類（Python）

import matplotlib.pyplot as plt
import numpy as np

# 材料とそのフレーリッヒ結合定数
materials_polaron = {
    'InSb': {'alpha': 0.023, 'type': 'large'},
    'GaAs': {'alpha': 0.068, 'type': 'large'},
    'CdTe': {'alpha': 0.39, 'type': 'large'},
    'AgBr': {'alpha': 1.53, 'type': 'large'},
    'SrTiO3': {'alpha': 3.77, 'type': 'intermediate'},
    'KCl': {'alpha': 5.0, 'type': 'intermediate'},
    'TlBr': {'alpha': 5.8, 'type': 'intermediate'},
    'NaCl': {'alpha': 6.8, 'type': 'small'},
    'KI': {'alpha': 7.5, 'type': 'small'},
    'Transition metal oxides': {'alpha': 10, 'type': 'small'},
}

# プロット
fig, ax = plt.subplots(figsize=(12, 6))

colors = {'large': 'blue', 'intermediate': 'orange', 'small': 'red'}
labels = {'large': '大きなポーラロン', 'intermediate': '中間領域', 'small': '小さなポーラロン'}

for i, (name, data) in enumerate(materials_polaron.items()):
    color = colors[data['type']]
    label = labels[data['type']] if i == 0 or data['type'] != list(materials_polaron.values())[i-1]['type'] else None
    ax.scatter(data['alpha'], 0, s=200, c=color, alpha=0.7, edgecolors='black', linewidth=1.5, label=label, zorder=3)
    ax.text(data['alpha'], 0.05, name, rotation=45, ha='left', va='bottom', fontsize=9)

# 転移境界
ax.axvline(x=6, color='gray', linestyle='--', linewidth=2, alpha=0.7, label='転移境界 (α ≈ 6)')
ax.axvspan(0, 6, alpha=0.1, color='blue', label='大きなポーラロン領域')
ax.axvspan(6, 12, alpha=0.1, color='red', label='小さなポーラロン領域')

ax.set_xlabel('フレーリッヒ結合定数 α', fontsize=13, fontweight='bold')
ax.set_xlim(-0.5, 12)
ax.set_ylim(-0.1, 0.2)
ax.set_yticks([])
ax.set_title('材料のポーラロン分類', fontsize=14, fontweight='bold')
ax.legend(loc='upper right', fontsize=10)
ax.grid(True, axis='x', alpha=0.3)

plt.tight_layout()
plt.savefig('polaron_materials.png', dpi=300, bbox_inches='tight')
plt.show()

print("\n材料別ポーラロンタイプ:")
for name, data in materials_polaron.items():
    print(f"{name:25s}: α = {data['alpha']:5.2f}  →  {data['type']}")

4. 電子自己エネルギーと質量繰り込み

4.1 電子自己エネルギー

電子-フォノン相互作用により、電子の自己エネルギー（self-energy）$\Sigma$が生じます。これはGreen関数形式で記述されます：

$$\Sigma_{n\mathbf{k}}(\omega) = \sum_{m\mathbf{q}\nu} \frac{|g_{nm\nu}(\mathbf{k},\mathbf{q})|^2}{\hbar\omega - \varepsilon_{m,\mathbf{k}+\mathbf{q}} + \hbar\omega_{\mathbf{q}\nu} + i\delta}$$ $$+ \sum_{m\mathbf{q}\nu} \frac{|g_{nm\nu}(\mathbf{k},\mathbf{q})|^2}{\hbar\omega - \varepsilon_{m,\mathbf{k}+\mathbf{q}} - \hbar\omega_{\mathbf{q}\nu} + i\delta}$$

自己エネルギーは複素数で：

$$\Sigma = \Sigma' + i\Sigma''$$

$\Sigma'$: 実部 - エネルギー繰り込み（バンドギャップ、有効質量の変化）
$\Sigma''$: 虚部 - 準粒子寿命（散乱率）

4.2 質量繰り込み（Mass Renormalization）

電子-フォノン相互作用により、電子の有効質量が変化します。質量増強因子（mass enhancement factor）は：

$$\frac{m^*}{m_\text{band}} = 1 + \lambda$$

ここで、$\lambda$は電子-フォノン結合定数です。自己エネルギーの実部から：

$$\lambda = -\left.\frac{\partial \Sigma'(\omega)}{\partial \omega}\right|_{\omega = \varepsilon_F}$$

4.3 準粒子寿命と散乱率

自己エネルギーの虚部は準粒子寿命に関連します：

$$\frac{1}{\tau} = -\frac{2\Sigma''(\varepsilon_F)}{\hbar}$$

これは電気抵抗率と直接関係します：

$$\rho = \frac{m^*}{ne^2\tau}$$

4.4 温度依存性

散乱率の温度依存性は、フォノン占有数を通じて現れます：

$$\frac{1}{\tau(T)} \propto \int d\omega \, \alpha^2F(\omega) \left[2n_\text{B}(\omega, T) + 1\right]$$

ここで、$n_\text{B}(\omega, T) = 1/(e^{\hbar\omega/k_B T} - 1)$はBose-Einstein分布です。

温度領域による振る舞い：

高温（$T \gg \Theta_D$）: $1/\tau \propto T$（線形）
低温（$T \ll \Theta_D$）: $1/\tau \propto T^5$（Bloch-Grüneisen）

散乱率の温度依存性（Python）

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad

# Bose-Einstein分布
def bose_einstein(omega, T):
    """ボーズ・アインシュタイン分布"""
    hbar = 1.054e-34  # J·s
    kB = 1.381e-23  # J/K
    if T == 0:
        return 0
    x = hbar * omega / (kB * T)
    if x > 100:  # オーバーフロー回避
        return 0
    return 1.0 / (np.exp(x) - 1)

# Debyeモデルのフォノン状態密度
def debye_dos(omega, omega_D):
    """デバイモデル状態密度"""
    if omega <= omega_D:
        return 3 * omega**2 / omega_D**3
    else:
        return 0

# 散乱率計算（簡略版）
def scattering_rate(T, omega_D, alpha2F_constant=1.0):
    """
    電子-フォノン散乱率

    Parameters:
    -----------
    T : float
        温度 (K)
    omega_D : float
        デバイ周波数 (rad/s)
    alpha2F_constant : float
        α²F(ω)の定数因子
    """
    def integrand(omega):
        n_B = bose_einstein(omega, T)
        dos = debye_dos(omega, omega_D)
        return dos * (2*n_B + 1)

    result, _ = quad(integrand, 0, omega_D)
    return alpha2F_constant * result

# 計算
theta_D = 300  # K (Debye温度)
hbar = 1.054e-34
kB = 1.381e-23
omega_D = kB * theta_D / hbar  # rad/s

temperatures = np.logspace(-1, 3, 100)  # 0.1 K to 1000 K
rates = np.array([scattering_rate(T, omega_D) for T in temperatures])

# 理論的な温度依存性
T_high = temperatures[temperatures > theta_D]
T_low = temperatures[temperatures < theta_D/5]
rate_linear = T_high / theta_D  # 高温での線形
rate_T5 = (T_low / theta_D)**5  # 低温でのT^5

# プロット
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))

# 線形スケール
ax1.plot(temperatures, rates/rates[len(rates)//2], 'b-', linewidth=2, label='数値計算')
ax1.plot(T_high, rate_linear/rate_linear[0] * rates[len(rates)//2]/rates[len(rates)//2],
         'r--', linewidth=2, label='高温極限 (∝ T)')
ax1.axvline(theta_D, color='gray', linestyle=':', linewidth=2, alpha=0.7, label=f'Θ_D = {theta_D} K')
ax1.set_xlabel('温度 (K)', fontsize=12)
ax1.set_ylabel('相対散乱率', fontsize=12)
ax1.set_title('電子-フォノン散乱率（線形スケール）', fontsize=13, fontweight='bold')
ax1.legend(fontsize=11)
ax1.grid(True, alpha=0.3)
ax1.set_xlim(0, 1000)

# 対数スケール
ax2.loglog(temperatures, rates/rates[len(rates)//2], 'b-', linewidth=2, label='数値計算')
ax2.loglog(T_high, rate_linear/rate_linear[0] * rates[len(rates)//2]/rates[len(rates)//2],
           'r--', linewidth=2, label='高温 (∝ T)')
ax2.loglog(T_low, rate_T5/rate_T5[-1] * rates[10]/rates[len(rates)//2],
           'g--', linewidth=2, label='低温 (∝ T⁵)')
ax2.axvline(theta_D, color='gray', linestyle=':', linewidth=2, alpha=0.7, label=f'Θ_D = {theta_D} K')
ax2.set_xlabel('温度 (K)', fontsize=12)
ax2.set_ylabel('相対散乱率', fontsize=12)
ax2.set_title('電子-フォノン散乱率（対数スケール）', fontsize=13, fontweight='bold')
ax2.legend(fontsize=11)
ax2.grid(True, alpha=0.3, which='both')

plt.tight_layout()
plt.savefig('scattering_rate_temperature.png', dpi=300, bbox_inches='tight')
plt.show()

print(f"デバイ温度: {theta_D} K")
print(f"高温領域 (T > {theta_D} K): 散乱率 ∝ T")
print(f"低温領域 (T < {theta_D/5} K): 散乱率 ∝ T^5")

5. Eliashberg関数と電子-フォノン結合定数

5.1 Eliashberg関数α²F(ω)

Eliashberg関数は、電子-フォノン相互作用の周波数依存性を記述する重要な物理量です：

$$\alpha^2F(\omega) = \frac{1}{N(E_F)} \sum_{nm\mathbf{k}\mathbf{q}\nu} |g_{nm\nu}(\mathbf{k},\mathbf{q})|^2 \delta(\varepsilon_{n\mathbf{k}} - E_F) \delta(\varepsilon_{m,\mathbf{k}+\mathbf{q}} - E_F) \delta(\omega - \omega_{\mathbf{q}\nu})$$

物理的意味：

$\alpha^2$: 電子-フォノン結合の強さの2乗
$F(\omega)$: フォノン状態密度に類似
フェルミ面での電子状態とフォノンの結合を反映

5.2 電子-フォノン結合定数λ

Eliashberg関数から、総電子-フォノン結合定数が定義されます：

$$\lambda = 2\int_0^\infty \frac{\alpha^2F(\omega)}{\omega} d\omega$$

この$\lambda$は以下の物理量に直接関係します：

有効質量増強: $m^*/m = 1 + \lambda$
超伝導転移温度（後述）
電気抵抗率
比熱の電子項の増強

5.3 フォノンモード別の寄与

各フォノンモードの寄与を評価できます：

$$\lambda_\nu = 2\int_0^\infty \frac{\alpha^2F_\nu(\omega)}{\omega} d\omega$$

全結合定数は：

$$\lambda = \sum_\nu \lambda_\nu$$

5.4 典型的なα²F(ω)の形状

Eliashberg関数の計算と可視化（Python）

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import simps

def alpha2F_debye(omega, omega_D, lambda_total):
    """
    デバイモデルのα²F(ω)

    Parameters:
    -----------
    omega : array
        角周波数 (rad/s)
    omega_D : float
        デバイ周波数 (rad/s)
    lambda_total : float
        総結合定数
    """
    result = np.zeros_like(omega)
    mask = (omega > 0) & (omega <= omega_D)
    # 規格化条件: ∫(α²F/ω)dω = λ/2
    # デバイモデル: α²F ∝ ω²
    result[mask] = (3 * lambda_total / (2 * omega_D**2)) * omega[mask]**2
    return result

def alpha2F_einstein(omega, omega_E, lambda_E, width):
    """
    アインシュタインモデル（ローレンツ型）のα²F(ω)

    Parameters:
    -----------
    omega : array
        角周波数 (rad/s)
    omega_E : float
        アインシュタイン周波数 (rad/s)
    lambda_E : float
        このモードの結合定数
    width : float
        ローレンツ幅 (rad/s)
    """
    # ローレンツ型
    lorentzian = (width / (2*np.pi)) / ((omega - omega_E)**2 + (width/2)**2)
    # 規格化
    norm = simps(lorentzian / omega_E, omega)
    return lambda_E / (2*norm) * lorentzian

def alpha2F_realistic(omega, params):
    """
    現実的な多モードα²F(ω)
    音響+光学フォノンの組み合わせ
    """
    omega_D = params['omega_D']
    lambda_ac = params['lambda_acoustic']
    omega_opt = params['omega_optical']
    lambda_opt = params['lambda_optical']
    width_opt = params['width_optical']

    # 音響フォノン（デバイ型）
    alpha2F_ac = alpha2F_debye(omega, omega_D, lambda_ac)

    # 光学フォノン（ローレンツ型）
    alpha2F_opt = alpha2F_einstein(omega, omega_opt, lambda_opt, width_opt)

    return alpha2F_ac + alpha2F_opt

# パラメータ設定
hbar = 1.054e-34
kB = 1.381e-23
theta_D = 400  # K
omega_D = kB * theta_D / hbar

# 周波数軸
omega = np.linspace(0, 3*omega_D, 1000)
freq_THz = omega * hbar / (2*np.pi*1e12)  # THz単位

# ケース1: 単純金属（アルミニウム型）
params_simple = {
    'omega_D': omega_D,
    'lambda_acoustic': 0.4,
    'omega_optical': 2*omega_D,
    'lambda_optical': 0.0,
    'width_optical': 0.1*omega_D
}

alpha2F_simple = alpha2F_realistic(omega, params_simple)

# ケース2: 複雑な材料（MgB2型）
params_complex = {
    'omega_D': omega_D,
    'lambda_acoustic': 0.3,
    'omega_optical': 2.5*omega_D,
    'lambda_optical': 0.6,
    'width_optical': 0.2*omega_D
}

alpha2F_complex = alpha2F_realistic(omega, params_complex)

# 結合定数を計算
lambda_simple = 2 * simps(alpha2F_simple / omega, omega)
lambda_complex = 2 * simps(alpha2F_complex / omega, omega)

# プロット
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(12, 10))

# ケース1
ax1.fill_between(freq_THz, 0, alpha2F_simple*1e-12, alpha=0.3, color='blue', label='Total')
ax1.plot(freq_THz, alpha2F_simple*1e-12, 'b-', linewidth=2)
ax1.axvline(kB*theta_D/(2*np.pi*1e12*hbar), color='gray', linestyle='--',
            linewidth=1.5, alpha=0.7, label=f'Θ_D = {theta_D} K')
ax1.set_xlabel('周波数 (THz)', fontsize=12)
ax1.set_ylabel('α²F(ω) (THz⁻¹)', fontsize=12)
ax1.set_title(f'単純金属型 Eliashberg関数 (λ = {lambda_simple:.3f})',
              fontsize=13, fontweight='bold')
ax1.legend(fontsize=11)
ax1.grid(True, alpha=0.3)
ax1.set_xlim(0, freq_THz[-1])

# ケース2
alpha2F_ac_only = alpha2F_debye(omega, omega_D, params_complex['lambda_acoustic'])
alpha2F_opt_only = alpha2F_einstein(omega, params_complex['omega_optical'],
                                     params_complex['lambda_optical'],
                                     params_complex['width_optical'])

ax2.fill_between(freq_THz, 0, alpha2F_complex*1e-12, alpha=0.3, color='purple', label='Total')
ax2.plot(freq_THz, alpha2F_complex*1e-12, 'purple', linewidth=2.5, label=f'Total (λ = {lambda_complex:.3f})')
ax2.plot(freq_THz, alpha2F_ac_only*1e-12, 'b--', linewidth=1.5, alpha=0.7,
         label=f'Acoustic (λ = {params_complex["lambda_acoustic"]:.3f})')
ax2.plot(freq_THz, alpha2F_opt_only*1e-12, 'r--', linewidth=1.5, alpha=0.7,
         label=f'Optical (λ = {params_complex["lambda_optical"]:.3f})')
ax2.set_xlabel('周波数 (THz)', fontsize=12)
ax2.set_ylabel('α²F(ω) (THz⁻¹)', fontsize=12)
ax2.set_title('複雑な材料型 Eliashberg関数（音響+光学）', fontsize=13, fontweight='bold')
ax2.legend(fontsize=11)
ax2.grid(True, alpha=0.3)
ax2.set_xlim(0, freq_THz[-1])

plt.tight_layout()
plt.savefig('eliashberg_function.png', dpi=300, bbox_inches='tight')
plt.show()

print("\n=== Eliashberg関数と結合定数 ===")
print(f"\nケース1（単純金属型）:")
print(f"  λ_total = {lambda_simple:.3f}")
print(f"  有効質量増強: m*/m = {1 + lambda_simple:.3f}")

print(f"\nケース2（複雑な材料型）:")
print(f"  λ_acoustic = {params_complex['lambda_acoustic']:.3f}")
print(f"  λ_optical = {params_complex['lambda_optical']:.3f}")
print(f"  λ_total = {lambda_complex:.3f}")
print(f"  有効質量増強: m*/m = {1 + lambda_complex:.3f}")

5.5 対数平均フォノン周波数

超伝導理論で重要な対数平均周波数は：

$$\omega_\text{log} = \exp\left[\frac{2}{\lambda}\int_0^\infty \frac{d\omega}{\omega} \alpha^2F(\omega) \ln\omega\right]$$

これは超伝導転移温度の見積もりに使用されます。

6. 超伝導への応用

6.1 BCS理論の基礎

BCS（Bardeen-Cooper-Schrieffer）理論では、フォノンを媒介とした電子間の引力的相互作用により、クーパー対が形成され、超伝導が発現します。

graph TD A[電子1がフォノンを放出] --> B[格子が変形] B --> C[電子2がフォノンを吸収] C --> D[実効的な引力相互作用] D --> E[クーパー対形成] E --> F[超伝導状態] G[電子-フォノン結合 λ] --> D H[クーロン斥力 μ*] -.阻害.-> E style A fill:#e3f2fd style D fill:#fff9c4 style F fill:#c8e6c9

6.2 McMillan方程式

超伝導転移温度$T_c$は、McMillan方程式で近似的に与えられます：

$$T_c = \frac{\omega_\text{log}}{1.20} \exp\left[-\frac{1.04(1+\lambda)}{\lambda - \mu^*(1+0.62\lambda)}\right]$$

ここで：

$\omega_\text{log}$: 対数平均フォノン周波数
$\lambda$: 電子-フォノン結合定数
$\mu^*$: 繰り込みクーロン斥力パラメータ（典型的に0.1～0.15）

6.3 強結合理論（Eliashberg理論）

$\lambda > 1$の強結合領域では、より精密なEliashberg方程式を解く必要があります：

$$\Delta(\omega_n) = \pi T \sum_m \int_0^\infty d\omega' \frac{\alpha^2F(\omega')}{\omega' - i(\omega_n - \omega_m)} \frac{\Delta(\omega_m)}{\sqrt{\omega_m^2 + \Delta^2(\omega_m)}}$$

ここで、$\omega_n = (2n+1)\pi T$は松原振動数です。

6.4 Allen-Dynes修正式

McMillan方程式の改良版として、Allen-Dynes式があります：

$$T_c = \frac{\omega_\text{log}}{1.20} \exp\left[-\frac{1.04(1+\lambda)}{\lambda - \mu^*(1+0.62\lambda)}\right] f_1 f_2$$

ここで、$f_1$、$f_2$は$\lambda$と$\omega_\text{log}$に依存する補正因子です。

Tc予測の計算（Python）

import numpy as np
import matplotlib.pyplot as plt

def mcmillan_tc(omega_log, lambda_ep, mu_star=0.13):
    """
    McMillan方程式によるTc計算

    Parameters:
    -----------
    omega_log : float
        対数平均フォノン周波数 (K)
    lambda_ep : float
        電子-フォノン結合定数
    mu_star : float
        繰り込みクーロン斥力パラメータ

    Returns:
    --------
    Tc : float
        超伝導転移温度 (K)
    """
    numerator = 1.04 * (1 + lambda_ep)
    denominator = lambda_ep - mu_star * (1 + 0.62 * lambda_ep)

    if denominator <= 0:
        return 0  # 超伝導にならない

    Tc = (omega_log / 1.20) * np.exp(-numerator / denominator)
    return Tc

def allen_dynes_tc(omega_log, lambda_ep, omega_2=None, mu_star=0.13):
    """
    Allen-Dynes修正式によるTc計算

    Parameters:
    -----------
    omega_log : float
        対数平均フォノン周波数 (K)
    lambda_ep : float
        電子-フォノン結合定数
    omega_2 : float, optional
        二次モーメント周波数 (K)
    mu_star : float
        繰り込みクーロン斥力パラメータ
    """
    if omega_2 is None:
        omega_2 = omega_log * 1.2  # 典型的な近似

    # 補正因子f1
    f1 = (1 + (lambda_ep / 2.46)**1.5)**(1/3)

    # 補正因子f2
    f2 = 1 + ((omega_2 / omega_log - 1) * lambda_ep**2) / (lambda_ep**2 + 1.82)

    # McMillan基本式
    Tc_base = mcmillan_tc(omega_log, lambda_ep, mu_star)

    return Tc_base * f1 * f2

# 実際の超伝導材料データ
materials_sc = {
    'Al': {'omega_log': 375, 'lambda': 0.43, 'Tc_exp': 1.2},
    'Pb': {'omega_log': 108, 'lambda': 1.55, 'Tc_exp': 7.2},
    'Nb': {'omega_log': 275, 'lambda': 1.04, 'Tc_exp': 9.3},
    'MgB2': {'omega_log': 600, 'lambda': 0.87, 'Tc_exp': 39},
    'Nb3Sn': {'omega_log': 220, 'lambda': 1.50, 'Tc_exp': 18},
}

print("=== 超伝導転移温度の理論計算 ===\n")
print(f"{'材料':<8} {'λ':>6} {'ω_log':>8} {'Tc(McM)':>10} {'Tc(AD)':>10} {'Tc(exp)':>10}")
print("-" * 60)

for name, params in materials_sc.items():
    Tc_mcm = mcmillan_tc(params['omega_log'], params['lambda'])
    Tc_ad = allen_dynes_tc(params['omega_log'], params['lambda'])
    Tc_exp = params['Tc_exp']

    print(f"{name:<8} {params['lambda']:>6.2f} {params['omega_log']:>6.0f} K "
          f"{Tc_mcm:>9.2f} K {Tc_ad:>9.2f} K {Tc_exp:>9.2f} K")

# λとTcの関係をプロット
lambda_range = np.linspace(0.1, 2.5, 100)
omega_log_values = [200, 400, 600]  # K

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))

# McMillan方程式
for omega_log in omega_log_values:
    Tc_values = [mcmillan_tc(omega_log, lam) for lam in lambda_range]
    ax1.plot(lambda_range, Tc_values, linewidth=2, label=f'ω_log = {omega_log} K')

# 実験データをプロット
for name, params in materials_sc.items():
    ax1.plot(params['lambda'], params['Tc_exp'], 'o', markersize=10,
             label=f'{name} (exp)', zorder=10)

ax1.set_xlabel('電子-フォノン結合定数 λ', fontsize=12)
ax1.set_ylabel('Tc (K)', fontsize=12)
ax1.set_title('McMillan方程式: Tc vs λ', fontsize=13, fontweight='bold')
ax1.legend(fontsize=10, ncol=2)
ax1.grid(True, alpha=0.3)
ax1.set_xlim(0, 2.5)
ax1.set_ylim(0, 50)

# 理論vs実験の比較
materials_list = list(materials_sc.keys())
Tc_exp_list = [materials_sc[m]['Tc_exp'] for m in materials_list]
Tc_mcm_list = [mcmillan_tc(materials_sc[m]['omega_log'], materials_sc[m]['lambda'])
               for m in materials_list]
Tc_ad_list = [allen_dynes_tc(materials_sc[m]['omega_log'], materials_sc[m]['lambda'])
              for m in materials_list]

x = np.arange(len(materials_list))
width = 0.25

bars1 = ax2.bar(x - width, Tc_exp_list, width, label='実験値', color='green', alpha=0.7)
bars2 = ax2.bar(x, Tc_mcm_list, width, label='McMillan式', color='blue', alpha=0.7)
bars3 = ax2.bar(x + width, Tc_ad_list, width, label='Allen-Dynes式', color='red', alpha=0.7)

ax2.set_xlabel('材料', fontsize=12)
ax2.set_ylabel('Tc (K)', fontsize=12)
ax2.set_title('超伝導転移温度: 理論 vs 実験', fontsize=13, fontweight='bold')
ax2.set_xticks(x)
ax2.set_xticklabels(materials_list)
ax2.legend(fontsize=11)
ax2.grid(True, alpha=0.3, axis='y')

plt.tight_layout()
plt.savefig('superconductivity_tc.png', dpi=300, bbox_inches='tight')
plt.show()

6.5 高温超伝導への拡張

銅酸化物やFeベース超伝導体などの高温超伝導体では、フォノン以外の媒介も重要です：

スピン揺らぎ（spin fluctuations）
軌道揺らぎ（orbital fluctuations）
電荷揺らぎ（charge fluctuations）

一般化されたEliashberg関数では、これらの寄与も含めることができます。

7. 計算手法とPython実装

7.1 第一原理計算の枠組み

電子-フォノン結合定数は、密度汎関数摂動理論（DFPT）により第一原理計算できます。主要なソフトウェア：

Quantum ESPRESSO: ph.x（フォノン）、EPW（電子-フォノン）
ABINIT: 電子-フォノン結合計算機能
VASP: DFPT計算

7.2 計算ワークフロー

graph TD A[結晶構造] --> B[DFT基底状態計算] B --> C[フォノン計算 DFPT] C --> D[電子-フォノン行列要素] D --> E[α²F ω の計算] E --> F1[λの計算] E --> F2[ω_logの計算] F1 --> G[Tc予測] F2 --> G style A fill:#e3f2fd style D fill:#fff9c4 style G fill:#c8e6c9

7.3 簡略化モデルによる計算例

電子-フォノン結合の計算ツール（Python）

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import simps
from scipy.optimize import brentq

class ElectronPhononCalculator:
    """電子-フォノン結合計算クラス"""

    def __init__(self, omega_D, N_EF, g_avg):
        """
        Parameters:
        -----------
        omega_D : float
            デバイ周波数 (rad/s)
        N_EF : float
            フェルミ面での状態密度 (states/eV/atom)
        g_avg : float
            平均電子-フォノン結合行列要素 (eV)
        """
        self.omega_D = omega_D
        self.N_EF = N_EF
        self.g_avg = g_avg
        self.hbar = 1.054e-34  # J·s
        self.kB = 1.381e-23  # J/K
        self.eV = 1.602e-19  # J

    def phonon_dos(self, omega):
        """
        デバイモデルフォノン状態密度

        F(ω) = 3ω²/ω_D³ for ω ≤ ω_D
        """
        if isinstance(omega, np.ndarray):
            F = np.zeros_like(omega)
            mask = (omega > 0) & (omega <= self.omega_D)
            F[mask] = 3 * omega[mask]**2 / self.omega_D**3
            return F
        else:
            if 0 < omega <= self.omega_D:
                return 3 * omega**2 / self.omega_D**3
            else:
                return 0

    def alpha2F(self, omega):
        """
        Eliashberg関数 α²F(ω)

        簡略化: α²F(ω) = N(EF) *  * F(ω)
        """
        g_squared = (self.g_avg * self.eV)**2  # J²
        F = self.phonon_dos(omega)
        # 単位: (states/J) * J² * (rad/s)⁻¹ = J·rad/s
        alpha2F = self.N_EF / self.eV * g_squared * F
        return alpha2F

    def calculate_lambda(self, omega_max=None):
        """
        電子-フォノン結合定数λを計算

        λ = 2 ∫ [α²F(ω)/ω] dω
        """
        if omega_max is None:
            omega_max = self.omega_D

        omega = np.linspace(1e-6, omega_max, 5000)
        integrand = self.alpha2F(omega) / omega
        lambda_ep = 2 * simps(integrand, omega)

        return lambda_ep

    def calculate_omega_log(self, lambda_ep=None, omega_max=None):
        """
        対数平均フォノン周波数を計算

        ω_log = exp[(2/λ) ∫ (α²F(ω)/ω) ln(ω) dω]
        """
        if lambda_ep is None:
            lambda_ep = self.calculate_lambda(omega_max)

        if omega_max is None:
            omega_max = self.omega_D

        omega = np.linspace(1e-6, omega_max, 5000)
        integrand = (self.alpha2F(omega) / omega) * np.log(omega)
        integral = 2 * simps(integrand, omega)

        omega_log = np.exp(integral / lambda_ep)

        # ケルビンに変換
        omega_log_K = omega_log * self.hbar / self.kB

        return omega_log_K

    def plot_alpha2F(self, omega_max=None):
        """α²F(ω)をプロット"""
        if omega_max is None:
            omega_max = self.omega_D * 1.2

        omega = np.linspace(0, omega_max, 1000)
        alpha2F_values = np.array([self.alpha2F(w) for w in omega])

        # THz単位に変換
        freq_THz = omega * self.hbar / (2*np.pi*1e12)

        fig, ax = plt.subplots(figsize=(10, 6))
        ax.fill_between(freq_THz, 0, alpha2F_values, alpha=0.3, color='blue')
        ax.plot(freq_THz, alpha2F_values, 'b-', linewidth=2)
        ax.axvline(self.omega_D * self.hbar / (2*np.pi*1e12),
                   color='red', linestyle='--', linewidth=2,
                   label=f'Debye frequency')
        ax.set_xlabel('周波数 (THz)', fontsize=12)
        ax.set_ylabel('α²F(ω)', fontsize=12)
        ax.set_title('Eliashberg関数', fontsize=13, fontweight='bold')
        ax.legend(fontsize=11)
        ax.grid(True, alpha=0.3)
        plt.tight_layout()
        return fig

# 使用例: アルミニウムの簡略計算
print("=== 電子-フォノン結合計算例（アルミニウム）===\n")

# パラメータ（文献値に基づく近似）
theta_D_Al = 428  # K
hbar = 1.054e-34
kB = 1.381e-23
omega_D_Al = kB * theta_D_Al / hbar  # rad/s

N_EF_Al = 0.32  # states/eV/atom (フェルミ面状態密度)
g_avg_Al = 0.05  # eV (平均結合行列要素)

# 計算
calc = ElectronPhononCalculator(omega_D_Al, N_EF_Al, g_avg_Al)

lambda_Al = calc.calculate_lambda()
omega_log_Al = calc.calculate_omega_log(lambda_Al)

print(f"デバイ温度: {theta_D_Al} K")
print(f"フェルミ面状態密度: {N_EF_Al} states/eV/atom")
print(f"平均結合行列要素: {g_avg_Al} eV")
print(f"\n計算結果:")
print(f"  λ = {lambda_Al:.3f}")
print(f"  ω_log = {omega_log_Al:.1f} K")
print(f"  有効質量増強: m*/m = {1 + lambda_Al:.3f}")

# Tc予測
def mcmillan_tc(omega_log, lambda_ep, mu_star=0.13):
    numerator = 1.04 * (1 + lambda_ep)
    denominator = lambda_ep - mu_star * (1 + 0.62 * lambda_ep)
    if denominator <= 0:
        return 0
    Tc = (omega_log / 1.20) * np.exp(-numerator / denominator)
    return Tc

Tc_pred = mcmillan_tc(omega_log_Al, lambda_Al)
Tc_exp_Al = 1.2  # K

print(f"\n超伝導転移温度:")
print(f"  Tc (計算): {Tc_pred:.2f} K")
print(f"  Tc (実験): {Tc_exp_Al:.2f} K")
print(f"  誤差: {abs(Tc_pred - Tc_exp_Al)/Tc_exp_Al*100:.1f}%")

# α²F(ω)プロット
fig = calc.plot_alpha2F()
plt.savefig('alpha2F_aluminum.png', dpi=300, bbox_inches='tight')
plt.show()

# パラメータスキャン: g_avgの影響
print("\n=== パラメータスキャン: 結合行列要素の影響 ===")
g_values = np.linspace(0.02, 0.10, 20)
lambda_values = []
Tc_values = []

for g in g_values:
    calc_temp = ElectronPhononCalculator(omega_D_Al, N_EF_Al, g)
    lam = calc_temp.calculate_lambda()
    omega_log = calc_temp.calculate_omega_log(lam)
    Tc = mcmillan_tc(omega_log, lam)

    lambda_values.append(lam)
    Tc_values.append(Tc)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))

ax1.plot(g_values, lambda_values, 'b-', linewidth=2)
ax1.axhline(lambda_Al, color='red', linestyle='--', linewidth=1.5,
            alpha=0.7, label=f'Al (λ = {lambda_Al:.3f})')
ax1.set_xlabel('平均結合行列要素 g (eV)', fontsize=12)
ax1.set_ylabel('電子-フォノン結合定数 λ', fontsize=12)
ax1.set_title('g vs λ', fontsize=13, fontweight='bold')
ax1.legend(fontsize=11)
ax1.grid(True, alpha=0.3)

ax2.plot(lambda_values, Tc_values, 'r-', linewidth=2)
ax2.axvline(lambda_Al, color='blue', linestyle='--', linewidth=1.5,
            alpha=0.7, label=f'Al (Tc = {Tc_exp_Al:.2f} K)')
ax2.set_xlabel('電子-フォノン結合定数 λ', fontsize=12)
ax2.set_ylabel('Tc (K)', fontsize=12)
ax2.set_title('λ vs Tc', fontsize=13, fontweight='bold')
ax2.legend(fontsize=11)
ax2.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('parameter_scan.png', dpi=300, bbox_inches='tight')
plt.show()

7.4 実践的なワークフロー

Quantum ESPRESSOを使った実際の計算ステップ：

Quantum ESPRESSO計算スクリプト例

#!/bin/bash
# 電子-フォノン結合計算ワークフロー

# 1. SCF計算（自己無撞着場）
pw.x < scf.in > scf.out

# 2. フォノン計算（DFPT）
ph.x < ph.in > ph.out

# 3. 電子-フォノン結合行列要素の計算（粗いqメッシュ）
ph.x < elph.in > elph.out

# 4. Wannier補間でqメッシュを高密度化
epw.x < epw.in > epw.out

# 5. α²F(ω)とλの計算
# EPWの出力からα²F(ω)を抽出
python analyze_epw_output.py

# 6. Tc予測
python predict_tc.py --lambda 0.43 --omega_log 375

8. 実験的測定手法

8.1 角度分解光電子分光（ARPES）

ARPESは電子-フォノン結合を直接観測できる強力な手法です。

観測される量：

キンク（kink）: フォノンエネルギーでのバンド分散の折れ曲がり
自己エネルギー: スペクトル関数から$\Sigma'(\omega)$と$\Sigma''(\omega)$を抽出
質量増強: $m^*/m = 1 - \partial\Sigma'/\partial\omega$

8.2 トンネル分光

超伝導トンネル接合のコンダクタンスから、$\alpha^2F(\omega)$を直接得られます（McMillan-Rowell法）：

$$\alpha^2F(\omega) = \frac{1}{\pi} \left|\frac{d^2I}{dV^2}\right| / \left|\frac{dI}{dV}\right|$$

8.3 中性子散乱

非弾性中性子散乱により、フォノン線幅の波数・エネルギー依存性を測定：

$$\Gamma_{\mathbf{q}\nu} = 2\pi\omega_{\mathbf{q}\nu} \lambda_{\mathbf{q}\nu}$$

ここで、$\Gamma_{\mathbf{q}\nu}$は測定されるフォノン線幅です。

8.4 輸送測定

電気抵抗率の温度依存性から、電子-フォノン散乱率を抽出：

$$\rho(T) = \rho_0 + A_\text{ep} \int_0^{\omega_\text{max}} d\omega \, \alpha^2F(\omega) \left[\frac{2n_B(\omega,T) + 1}{\omega}\right]$$

8.5 比熱測定

電子比熱の増強から$\lambda$を見積もれます：

$$\gamma_\text{obs} = \gamma_\text{band}(1 + \lambda)$$

測定手法	得られる情報	利点	制限
ARPES	自己エネルギー、λ(k,ω)	波数・エネルギー分解	表面敏感
トンネル分光	α²F(ω)、λ	直接測定	超伝導体のみ
中性子散乱	λ(q,ν)	バルク敏感	大型試料必要
電気抵抗	λ（積分値）	簡便	間接的
比熱	λ（積分値）	簡便	間接的

まとめ

本章で学んだこと

1. 電子-フォノン相互作用の基礎

相互作用ハミルトニアンの導出と物理的意味
結合行列要素$g_{nm\nu}(\mathbf{k},\mathbf{q})$の役割
運動量・エネルギー保存則

2. 結合機構の分類

フレーリッヒ結合: 極性結晶、長距離相互作用（$\propto 1/q$）
変形ポテンシャル: 非極性結晶、短距離相互作用（$\propto q$）
材料の性質による使い分け

3. ポーラロン形成

大きなポーラロン: $\alpha < 6$、摂動論的、バンド輸送
小さなポーラロン: $\alpha > 6$、強結合、ホッピング輸送
有効質量の増強と移動度への影響

4. 自己エネルギーと繰り込み

複素自己エネルギー: $\Sigma = \Sigma' + i\Sigma''$
実部からの質量繰り込み: $m^*/m = 1 + \lambda$
虚部からの準粒子寿命: $1/\tau = -2\Sigma''/\hbar$

5. Eliashberg関数と結合定数

$\alpha^2F(\omega)$: 周波数依存性を記述
$\lambda = 2\int [\alpha^2F(\omega)/\omega] d\omega$: 総結合定数
対数平均周波数$\omega_\text{log}$の重要性

6. 超伝導への応用

フォノン媒介のクーパー対形成
McMillan方程式による$T_c$予測
強結合Eliashberg理論

7. 計算と測定手法

DFPT（密度汎関数摂動理論）による第一原理計算
ARPES、トンネル分光、中性子散乱による実験測定
輸送・熱物性からの間接的評価

重要な公式のまとめ

電子-フォノン相互作用ハミルトニアン：

$$\hat{H}_\text{el-ph} = \sum_{nm\mathbf{k}\mathbf{q}\nu} g_{nm\nu}(\mathbf{k},\mathbf{q}) c_{n\mathbf{k}+\mathbf{q}}^\dagger c_{m\mathbf{k}} (b_{\mathbf{q}\nu} + b_{-\mathbf{q}\nu}^\dagger)$$

Eliashberg関数：

$$\alpha^2F(\omega) = \frac{1}{N(E_F)} \sum |g|^2 \delta(\varepsilon - E_F) \delta(\omega - \omega_{\mathbf{q}\nu})$$

電子-フォノン結合定数：

$$\lambda = 2\int_0^\infty \frac{\alpha^2F(\omega)}{\omega} d\omega$$

質量繰り込み：

$$\frac{m^*}{m} = 1 + \lambda$$

McMillan方程式：

$$T_c = \frac{\omega_\text{log}}{1.20} \exp\left[-\frac{1.04(1+\lambda)}{\lambda - \mu^*(1+0.62\lambda)}\right]$$

次のステップ

第5章では、フォノン輸送現象を学びます：

熱伝導率の微視的理論
フォノン-フォノン散乱（Umklapp過程）
境界散乱とサイズ効果
格子熱伝導率の計算手法

演習問題

問題1: フレーリッヒ結合定数の計算

GaAsの電子-フォノン結合を評価します。以下のパラメータを使用してください：

$\varepsilon_\infty = 10.9$
$\varepsilon_0 = 12.9$
$\hbar\omega_\text{LO} = 36$ meV
有効ボーア半径 $a_B^* = 10$ nm

(a) フレーリッヒ結合定数$\alpha$を計算してください。
(b) GaAsはどちらのポーラロンタイプに分類されますか？
(c) 有効質量増強を見積もってください。

ヒント

フレーリッヒ結合定数は： $$\alpha = \frac{e^2}{4\pi\varepsilon_0} \frac{1}{2\hbar\omega_\text{LO} a_B^*} \left(\frac{1}{\varepsilon_\infty} - \frac{1}{\varepsilon_0}\right)$$

問題2: 散乱率の温度依存性

デバイ温度$\Theta_D = 400$ Kの金属における電子-フォノン散乱率を考えます。

(a) $T = 10$ K（低温）と$T = 500$ K（高温）での散乱率の比$\tau^{-1}(500K)/\tau^{-1}(10K)$を見積もってください。
(b) 電気抵抗率の温度依存性を議論してください。
(c) Pythonコードを用いて、$T = 1$ K～1000 Kの範囲で散乱率をプロットしてください。

問題3: Eliashberg関数とλの計算

簡略化されたEliashberg関数が以下で与えられるとします：

$$\alpha^2F(\omega) = \begin{cases} C\omega^2 & 0 < \omega \leq \omega_D \\ 0 & \omega > \omega_D \end{cases}$$

(a) 規格化条件$\lambda = 2\int_0^\infty [\alpha^2F(\omega)/\omega] d\omega = 0.5$から、定数$C$を求めてください。
(b) 対数平均周波数$\omega_\text{log}$を$\omega_D$の関数として求めてください。
(c) $\omega_D = 300$ K、$\mu^* = 0.13$として、McMillan方程式から$T_c$を計算してください。

問題4: 超伝導転移温度の材料設計

より高い$T_c$を持つ超伝導体を設計したいと考えています。

(a) McMillan方程式から、$T_c$を最大化するための戦略を3つ提案してください。
(b) $\lambda$を0.5から1.5まで変化させたときの$T_c$の変化を、$\omega_\text{log} = 400$ Kとして計算し、プロットしてください。
(c) なぜ$\lambda$を無限に大きくしても$T_c$は無限に大きくならないのか、物理的に説明してください。

問題5: 第一原理計算データの解析

Quantum ESPRESSOから得られた架空の$\alpha^2F(\omega)$データ（CSV形式）を解析します。

(a) Pythonを用いて、データから$\lambda$と$\omega_\text{log}$を数値積分で計算するコードを書いてください。
(b) フォノンモードを音響と光学に分離し、それぞれの$\lambda$への寄与を評価してください。
(c) 計算した$\lambda$と$\omega_\text{log}$から$T_c$を予測し、実験値と比較してください。

サンプルデータ構造

# omega (THz), alpha2F (THz^-1)
0.0, 0.0
1.0, 0.05
2.0, 0.12
...
20.0, 0.0

問題6: ARPESデータからの自己エネルギー抽出（発展）

ARPES測定から得られたスペクトル関数$A(\mathbf{k}, \omega)$から、自己エネルギーの実部$\Sigma'(\omega)$と虚部$\Sigma''(\omega)$を抽出する方法を調べてください。

(a) スペクトル関数とGreen関数の関係を示してください。
(b) 「キンク」が現れるエネルギーとフォノンエネルギーの関係を説明してください。
(c) $\Sigma'(\omega)$の傾きから$\lambda$を見積もる方法を示してください。

参考文献

Grimvall, G. "The Electron-Phonon Interaction in Metals" (North-Holland, 1981)
Ziman, J. M. "Electrons and Phonons" (Oxford University Press, 2001)
Alexandrov, A. S. & Devreese, J. T. "Advances in Polaron Physics" (Springer, 2010)
Eliashberg, G. M. "Interactions between electrons and lattice vibrations in a superconductor" Sov. Phys. JETP 11, 696 (1960)
McMillan, W. L. "Transition temperature of strong-coupled superconductors" Phys. Rev. 167, 331 (1968)
Allen, P. B. & Dynes, R. C. "Transition temperature of strong-coupled superconductors reanalyzed" Phys. Rev. B 12, 905 (1975)
Giustino, F. "Electron-phonon interactions from first principles" Rev. Mod. Phys. 89, 015003 (2017)
Baroni, S., de Gironcoli, S., Dal Corso, A. & Giannozzi, P. "Phonons and related crystal properties from density-functional perturbation theory" Rev. Mod. Phys. 73, 515 (2001)
Poncé, S., Margine, E. R., Verdi, C. & Giustino, F. "EPW: Electron-phonon coupling, transport and superconducting properties using maximally localized Wannier functions" Comput. Phys. Commun. 209, 116 (2016)
Damascelli, A., Hussain, Z. & Shen, Z.-X. "Angle-resolved photoemission studies of the cuprate superconductors" Rev. Mod. Phys. 75, 473 (2003)

性質	大きなポーラロン	小さなポーラロン
結合定数	\(\alpha < 6\)	\(\alpha > 6\)
サイズ	多数の単位胞	1～数単位胞
有効質量	\(m^* \approx m(1 + \alpha/6)\)	\(m^* \gg m\)（発散的）
輸送機構	バンド輸送	ホッピング輸送
温度依存性	\(\mu \propto T^{-n}\)	\(\mu \propto \exp(-E_a/k_B T)\)
理論手法	摂動論	変分法、経路積分

中級フォノン物理学

第4章: 電子-フォノン結合

学習目標

1. 電子-フォノン相互作用ハミルトニアン

1.1 基本的な枠組み

1.2 電子ハミルトニアン

1.3 フォノンハミルトニアン

1.4 電子-フォノン相互作用ハミルトニアン

物理的解釈

1.5 結合行列要素の導出

2. フレーリッヒハミルトニアンと変形ポテンシャル

2.1 フレーリッヒハミルトニアン（極性結合）

特徴的な性質

2.2 変形ポテンシャル結合

2.3 両者の比較

2.4 実際の材料での例

3. ポーラロンの物理

3.1 ポーラロンとは

3.2 大きなポーラロン（Large Polaron）

3.3 小さなポーラロン（Small Polaron）

3.4 ポーラロン転移

注意：中間領域

3.5 典型的な材料例

4. 電子自己エネルギーと質量繰り込み

4.1 電子自己エネルギー

4.2 質量繰り込み（Mass Renormalization）

4.3 準粒子寿命と散乱率

4.4 温度依存性

5. Eliashberg関数と電子-フォノン結合定数

5.1 Eliashberg関数α²F(ω)

5.2 電子-フォノン結合定数λ

5.3 フォノンモード別の寄与

5.4 典型的なα²F(ω)の形状

5.5 対数平均フォノン周波数

6. 超伝導への応用

6.1 BCS理論の基礎

6.2 McMillan方程式

6.3 強結合理論（Eliashberg理論）

6.4 Allen-Dynes修正式

6.5 高温超伝導への拡張

7. 計算手法とPython実装

7.1 第一原理計算の枠組み

7.2 計算ワークフロー

7.3 簡略化モデルによる計算例

7.4 実践的なワークフロー

8. 実験的測定手法

8.1 角度分解光電子分光（ARPES）

8.2 トンネル分光

8.3 中性子散乱

8.4 輸送測定

8.5 比熱測定

まとめ

本章で学んだこと

1. 電子-フォノン相互作用の基礎

2. 結合機構の分類

3. ポーラロン形成

4. 自己エネルギーと繰り込み

5. Eliashberg関数と結合定数

6. 超伝導への応用

7. 計算と測定手法

重要な公式のまとめ

次のステップ

演習問題

問題1: フレーリッヒ結合定数の計算

問題2: 散乱率の温度依存性

問題3: Eliashberg関数とλの計算

問題4: 超伝導転移温度の材料設計

問題5: 第一原理計算データの解析

問題6: ARPESデータからの自己エネルギー抽出（発展）

参考文献