第1章: カイラルフォノンの理論的基礎

1.1 カイラルフォノンの導入と歴史的背景

カイラルフォノンとは

カイラルフォノン（chiral phonon）は、格子振動が特定の回転方向を持つフォノンモードです。電磁波の円偏光に類似して、フォノンも左回り（L: left-handed）または右回り（R: right-handed）の角運動量を持つことができます。

歴史的発展

カイラルフォノンの理論は以下のように発展してきました：

物理的重要性

カイラルフォノン研究の重要性は以下の点にあります：

• 角運動量輸送: 熱流に伴う角運動量の輸送が可能
• 円二色性: 光学的活性と円偏光選択則の起源
• スピン-格子結合: 電子スピンとの相互作用（スピントロニクス応用）
• トポロジカル物性: ベリー曲率に基づく幾何学的性質
• 量子情報: フォノン量子状態の制御手段

1.2 フォノン角運動量（PAM）の定義

古典的定義

フォノン角運動量（Phonon Angular Momentum, PAM）は、格子振動における各原子の軌道角運動量の総和として定義されます：

ここで：

• \( m_i \): 原子 \( i \) の質量
• \( \mathbf{u}_i \): 原子 \( i \) の変位ベクトル
• \( \dot{\mathbf{u}}_i \): 原子 \( i \) の速度ベクトル

調和振動子展開

調和近似では、変位ベクトルを正規モード展開します：

ここで：

• \( \lambda \): フォノンモード（波数 \(\mathbf{q}\) とバンド指数）
• \( \omega_\lambda \): フォノン角周波数
• \( a_\lambda, a_\lambda^\dagger \): 生成・消滅演算子
• \( \mathbf{e}_{i\lambda} \): 固有ベクトル（偏極ベクトル）

量子化されたPAM

時間微分を考慮すると、角運動量の \( z \) 成分は：

ここで \( \sigma_\lambda \) はカイラル偏極度（chiral polarization）：

円偏光基底

2次元系（xy平面）では、円偏光基底が有用です：

線形基底 \((\hat{\mathbf{x}}, \hat{\mathbf{y}})\) と円偏光基底 \((\mathbf{e}_L, \mathbf{e}_R)\) の関係：

具体例: 単純2原子系

質量 \( m \) の2つの原子が円運動する場合を考えます：

このとき、各原子の角運動量は：

全角運動量：

この例は、2つの原子が反対方向に動いても、同じ回転方向であれば角運動量が打ち消されないことを示しています。

1.3 対称性要件とカイラル性

空間反転対称性

カイラルフォノンが存在するための必要条件は、空間反転対称性の破れです。

時間反転対称性

時間反転 \(\mathcal{T}\) の下で：

時間反転対称性がある系では、左回りと右回りのフォノンがエネルギー縮退します（Kramers縮退）。しかし、カイラル性は依然として定義可能です。これは電磁波の円偏光と同様です。

カイラルフォノンを支持する結晶構造

空間反転対称性が破れている主な結晶クラス：

対称性と選択則

点群の既約表現により、どのフォノンモードがカイラルになるかが決まります。例えば、C_3v対称性を持つ単層TMDCsでは：

• \( A_1 \): z方向振動（非カイラル）
• \( E' \): 面内振動（カイラル可能、低周波数）
• \( E'' \): 面内振動（カイラル可能、高周波数）

対称性破れの起源

実験的に空間反転対称性を破る方法：

1. 材料選択: 本質的に非対称な結晶構造（TMDCs, ペロブスカイトABOH₃）
2. 基板効果: 単層材料を基板上に配置（グラフェン/SiC）
3. 電場印加: 外部電場による対称性破れ
4. 歪み: 機械的歪みによる構造変化
5. 積層: 異種材料のvan der Waals積層

1.4 ベリー位相とトポロジカル性質

ベリー接続とベリー曲率

フォノン固有状態 \(|\lambda(\mathbf{q})\rangle\) のパラメータ（波数 \(\mathbf{q}\)）依存性から、ベリー接続（Berry connection）が定義されます：

ベリー曲率（Berry curvature）は：

または、直接計算では：

ベリー位相

閉じた経路 \(\mathcal{C}\) に沿ったベリー位相：

（Stokesの定理により、曲率の面積分に等しい）

カイラル性とベリー曲率の関係

カイラル偏極度 \(\sigma_\lambda\) とベリー曲率 \(\Omega_\lambda^z\) の間には関連があります。高対称点（例: \(\Gamma\)点）では：

トポロジカルフォノンとの関連

ベリー曲率の積分であるChern数は、トポロジカル不変量です：

非ゼロのChern数を持つフォノンバンドはトポロジカルフォノンと呼ばれます。

ワイル点とディラック点

3次元系では、バンド交差点でベリー曲率が特異性を持ちます：

• ワイル点: 2つのバンドが線形に交差、Chern電荷 \(\pm 1\)
• ディラック点: 4重縮退、2つのワイル点の重なり

例: グラフェンのK点

単層グラフェンの \( K \) と \( K' \) 点では、音響フォノンとTA/LAバンドがディラック点を形成します。K点周りの有効ハミルトニアンは：

ここで \( v_s \) は音速、\(\sigma_{x,y}\) はパウリ行列です。ベリー曲率は：

（符号は \( K \) と \( K' \) で反対）

1.5 群論による選択則

点群とフォノンモード

結晶の点群 \( G \) の既約表現により、フォノンモードの対称性が分類されます。\(\Gamma\)点（波数 \(\mathbf{q} = 0\)）でのフォノンは、点群の既約表現に従います。

C_3v点群の解析

C_3v対称性（TMDCs）の既約表現：

E表現は2次元表現で、2つの基底関数 \((E_x, E_y)\) を持ちます。これらは円偏光基底に変換できます：

D_3h点群の解析

D_3h対称性（グラフェン）は空間反転対称性を持ちます。\(\Gamma\)点では：

• \( A'_{2u} \): out-of-plane音響モード（非カイラル）
• \( E' \): in-plane光学モード（カイラル可能だが \(\Gamma\) では縮退）
• \( A''_{2u} \): out-of-plane光学モード（非カイラル）
• \( E'' \): in-plane音響モード（カイラル可能だが \(\Gamma\) では縮退）

重要: \(\Gamma\)点では空間反転対称性により左回り・右回りが縮退しますが、K点では対称性が下がり、カイラル性が現れます。

指標表（Character Table）の利用

C_3vの指標表：

3回回転 \( C_3 \) の下で、E表現の基底は以下のように変換されます：

円偏光基底では：

この位相因子の符号が異なることが、左回りと右回りの区別の起源です。

選択則の導出

ラマン散乱の遷移確率は、ラマンテンソル \( \mathcal{R} \) により決まります：

ここで \(\mathbf{e}_i, \mathbf{e}_s\) は入射・散乱光の偏光ベクトルです。

円偏光 \(\sigma^+ = (1, i, 0)/\sqrt{2}\) と \(\sigma^- = (1, -i, 0)/\sqrt{2}\) を使うと：

• \( \sigma^+ \to \sigma^+ \): 左回りフォノン励起
• \( \sigma^- \to \sigma^- \): 右回りフォノン励起
• \( \sigma^+ \to \sigma^- \) または \( \sigma^- \to \sigma^+ \): 禁制（選択則）

1.6 フォノン円偏光の数学的定式化

複素固有ベクトル表示

カイラルフォノンモードの固有ベクトルは複素数です。2原子からなる単位胞の場合：

カイラル偏極度の計算：

簡略化すると：

正規化条件

固有ベクトルの正規化：

ここで \( \alpha = x, y, z \) は座標成分です。

直交性条件

異なるモード間の直交性：

位相自由度と測定可能量

固有ベクトルには全体位相の自由度 \( e^{i\phi} \) がありますが、カイラル偏極度 \(\sigma\) はこの位相に依存しません（gauge不変）：

時間発展

時間依存性を含めた変位場：

瞬間的な軌跡を描くと、カイラルフォノンでは各原子が楕円軌道を描きます。

例: MoS₂のE'モード

単層MoS₂の \(\Gamma\) 点における E' モード（約383 cm^-1）の固有ベクトル（概略）：

（Mo原子のみが面内で円運動、S原子は静止）

カイラル偏極度：

1.7 実装: PAM計算プログラム

概要

ここでは、フォノン固有ベクトルからカイラル偏極度（PAM）を計算するPythonプログラムを実装します。VASP, Quantum ESPRESSO, Phonopy等の第一原理計算ソフトから得られた固有ベクトルデータを解析できます。

実装仕様

• 入力: 固有ベクトル（複素数配列）、原子質量
• 出力: カイラル偏極度 \(\sigma\)、角運動量 \(L_z\)
• 可視化: 原子軌道のアニメーション、\(\sigma\) vs 周波数プロット

コード実装

import numpy as np
import matplotlib.pyplot as plt
from typing import Tuple, List
from dataclasses import dataclass

@dataclass
class PhononMode:
    """フォノンモードのデータ構造"""
    frequency: float  # cm^-1 or THz
    eigenvector: np.ndarray  # shape: (n_atoms, 3), complex
    masses: np.ndarray  # shape: (n_atoms,)

    def __post_init__(self):
        """初期化後の検証"""
        assert self.eigenvector.shape[0] == len(self.masses), \
            "原子数と質量配列のサイズが不一致"
        assert self.eigenvector.shape[1] == 3, \
            "固有ベクトルは3次元（x,y,z）でなければならない"

def normalize_eigenvector(eigenvector: np.ndarray,
                          masses: np.ndarray) -> np.ndarray:
    """
    質量加重規格化

    Parameters:
    -----------
    eigenvector : ndarray, shape (n_atoms, 3)
        複素固有ベクトル
    masses : ndarray, shape (n_atoms,)
        原子質量（amu）

    Returns:
    --------
    normalized : ndarray, shape (n_atoms, 3)
        規格化された固有ベクトル
    """
    norm_sq = np.sum(masses[:, None] * np.abs(eigenvector)**2)
    return eigenvector / np.sqrt(norm_sq)

def compute_chiral_polarization(eigenvector: np.ndarray,
                                 masses: np.ndarray,
                                 direction: str = 'z') -> float:
    """
    カイラル偏極度を計算

    σ = -i Σᵢ (eᵢ* × eᵢ) · ẑ

    Parameters:
    -----------
    eigenvector : ndarray, shape (n_atoms, 3)
        複素固有ベクトル（正規化済み）
    masses : ndarray, shape (n_atoms,)
        原子質量
    direction : str
        角運動量の方向 ('x', 'y', 'z')

    Returns:
    --------
    sigma : float
        カイラル偏極度（-1: 右回り, +1: 左回り, 0: 線偏光）
    """
    # 正規化
    e_norm = normalize_eigenvector(eigenvector, masses)

    # ベクトル積 e* × e の計算
    cross_product = np.zeros(len(masses), dtype=complex)

    if direction == 'z':
        # (e* × e)_z = e_x* e_y - e_y* e_x
        cross_product = (np.conj(e_norm[:, 0]) * e_norm[:, 1] -
                         np.conj(e_norm[:, 1]) * e_norm[:, 0])
    elif direction == 'x':
        cross_product = (np.conj(e_norm[:, 1]) * e_norm[:, 2] -
                         np.conj(e_norm[:, 2]) * e_norm[:, 1])
    elif direction == 'y':
        cross_product = (np.conj(e_norm[:, 2]) * e_norm[:, 0] -
                         np.conj(e_norm[:, 0]) * e_norm[:, 2])
    else:
        raise ValueError(f"Unknown direction: {direction}")

    # -i を掛けて虚部を取る → 実部を取る
    sigma = -1j * np.sum(cross_product)

    return sigma.real

def compute_angular_momentum(eigenvector: np.ndarray,
                             masses: np.ndarray,
                             amplitude: float = 1.0,
                             frequency: float = 1.0) -> float:
    """
    古典的角運動量を計算

    L_z = Σᵢ mᵢ (uᵢ × u̇ᵢ)_z

    Parameters:
    -----------
    eigenvector : ndarray
        固有ベクトル
    masses : ndarray
        原子質量
    amplitude : float
        振動振幅（Å）
    frequency : float
        角周波数（rad/s or a.u.）

    Returns:
    --------
    L_z : float
        z方向角運動量
    """
    sigma = compute_chiral_polarization(eigenvector, masses, 'z')

    # L_z = σ * ℏ * n_phonon （量子論）
    # 古典的には L_z ∝ A² ω σ
    return sigma * amplitude**2 * frequency * np.sum(masses)

def create_circular_mode(n_atoms: int = 1,
                         chirality: str = 'L') -> Tuple[np.ndarray, np.ndarray]:
    """
    理想的な円偏光モードを生成（テスト用）

    Parameters:
    -----------
    n_atoms : int
        原子数
    chirality : str
        'L' (左回り) or 'R' (右回り)

    Returns:
    --------
    eigenvector : ndarray, shape (n_atoms, 3)
    masses : ndarray, shape (n_atoms,)
    """
    masses = np.ones(n_atoms)
    eigenvector = np.zeros((n_atoms, 3), dtype=complex)

    if chirality == 'L':
        eigenvector[:, 0] = 1.0
        eigenvector[:, 1] = 1j
    elif chirality == 'R':
        eigenvector[:, 0] = 1.0
        eigenvector[:, 1] = -1j
    else:
        raise ValueError("chirality must be 'L' or 'R'")

    return eigenvector, masses

def visualize_trajectory(eigenvector: np.ndarray,
                         masses: np.ndarray,
                         atom_index: int = 0,
                         n_steps: int = 100) -> None:
    """
    原子軌道を可視化

    Parameters:
    -----------
    eigenvector : ndarray
        固有ベクトル
    masses : ndarray
        原子質量
    atom_index : int
        可視化する原子のインデックス
    n_steps : int
        時間ステップ数
    """
    # 時間配列
    t = np.linspace(0, 2*np.pi, n_steps)

    # 複素固有ベクトルから実空間軌道を計算
    e = eigenvector[atom_index]
    trajectory = np.real(e[None, :] * np.exp(-1j * t[:, None]))

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

    # xy平面
    ax1.plot(trajectory[:, 0], trajectory[:, 1], 'b-', linewidth=2)
    ax1.plot(trajectory[0, 0], trajectory[0, 1], 'go', markersize=10, label='Start')
    ax1.plot(trajectory[-1, 0], trajectory[-1, 1], 'ro', markersize=10, label='End')
    ax1.set_xlabel('x displacement')
    ax1.set_ylabel('y displacement')
    ax1.set_title(f'Atom {atom_index} trajectory (xy plane)')
    ax1.axis('equal')
    ax1.grid(True)
    ax1.legend()

    # 時間発展
    ax2.plot(t, trajectory[:, 0], label='x(t)')
    ax2.plot(t, trajectory[:, 1], label='y(t)')
    ax2.plot(t, trajectory[:, 2], label='z(t)')
    ax2.set_xlabel('Time (arbitrary units)')
    ax2.set_ylabel('Displacement')
    ax2.set_title('Time evolution')
    ax2.legend()
    ax2.grid(True)

    plt.tight_layout()
    plt.savefig('phonon_trajectory.png', dpi=150, bbox_inches='tight')
    plt.show()

def analyze_phonon_modes(modes: List[PhononMode]) -> None:
    """
    複数のフォノンモードを解析

    Parameters:
    -----------
    modes : list of PhononMode
        解析するモードのリスト
    """
    frequencies = []
    sigmas = []

    for mode in modes:
        sigma = compute_chiral_polarization(mode.eigenvector, mode.masses)
        frequencies.append(mode.frequency)
        sigmas.append(sigma)

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

    colors = ['red' if s < -0.5 else 'blue' if s > 0.5 else 'gray'
              for s in sigmas]

    ax.scatter(frequencies, sigmas, c=colors, s=100, alpha=0.7)
    ax.axhline(y=0, color='k', linestyle='--', linewidth=0.5)
    ax.set_xlabel('Frequency (cm⁻¹)')
    ax.set_ylabel('Chiral Polarization σ')
    ax.set_title('Phonon Chirality Spectrum')
    ax.grid(True, alpha=0.3)

    # 凡例
    from matplotlib.patches import Patch
    legend_elements = [
        Patch(facecolor='blue', label='Left-handed (σ > 0.5)'),
        Patch(facecolor='red', label='Right-handed (σ < -0.5)'),
        Patch(facecolor='gray', label='Non-chiral (|σ| < 0.5)')
    ]
    ax.legend(handles=legend_elements)

    plt.tight_layout()
    plt.savefig('chirality_spectrum.png', dpi=150, bbox_inches='tight')
    plt.show()

# ===== 使用例 =====
if __name__ == "__main__":
    print("=== カイラルフォノン解析プログラム ===\n")

    # テスト1: 理想的な左回りモード
    print("Test 1: 理想的な左回りモード")
    e_L, m_L = create_circular_mode(n_atoms=1, chirality='L')
    sigma_L = compute_chiral_polarization(e_L, m_L)
    print(f"  カイラル偏極度: {sigma_L:.6f}")
    print(f"  期待値: +1.0")
    print(f"  判定: {'PASS' if abs(sigma_L - 1.0) < 1e-6 else 'FAIL'}\n")

    # テスト2: 理想的な右回りモード
    print("Test 2: 理想的な右回りモード")
    e_R, m_R = create_circular_mode(n_atoms=1, chirality='R')
    sigma_R = compute_chiral_polarization(e_R, m_R)
    print(f"  カイラル偏極度: {sigma_R:.6f}")
    print(f"  期待値: -1.0")
    print(f"  判定: {'PASS' if abs(sigma_R + 1.0) < 1e-6 else 'FAIL'}\n")

    # テスト3: 線偏光モード（x方向のみ）
    print("Test 3: 線偏光モード")
    e_linear = np.array([[1.0, 0.0, 0.0]], dtype=complex)
    m_linear = np.array([1.0])
    sigma_linear = compute_chiral_polarization(e_linear, m_linear)
    print(f"  カイラル偏極度: {sigma_linear:.6f}")
    print(f"  期待値: 0.0")
    print(f"  判定: {'PASS' if abs(sigma_linear) < 1e-6 else 'FAIL'}\n")

    # テスト4: 2原子系（MoS2風）
    print("Test 4: 2原子系（Mo + S）")
    # Mo原子が円運動、S原子が静止
    e_MoS2 = np.array([
        [1.0, 1j, 0.0],   # Mo原子（左回り）
        [0.0, 0.0, 0.0]   # S原子（静止）
    ], dtype=complex)
    m_MoS2 = np.array([95.95, 32.06])  # Mo, S の質量（amu）
    sigma_MoS2 = compute_chiral_polarization(e_MoS2, m_MoS2)
    print(f"  カイラル偏極度: {sigma_MoS2:.6f}")
    print(f"  Mo原子のみが寄与\n")

    # 可視化
    print("軌道可視化を生成中...")
    visualize_trajectory(e_L, m_L, atom_index=0)

    # 複数モードの解析例
    print("\n複数モードの解析例...")
    modes = [
        PhononMode(100.0, create_circular_mode(2, 'L')[0],
                   create_circular_mode(2, 'L')[1]),
        PhononMode(200.0, create_circular_mode(2, 'R')[0],
                   create_circular_mode(2, 'R')[1]),
        PhononMode(300.0, np.array([[1, 0, 0], [0, 1, 0]], dtype=complex),
                   np.ones(2)),
    ]
    analyze_phonon_modes(modes)

    print("\n解析完了！")

出力例

=== カイラルフォノン解析プログラム ===

Test 1: 理想的な左回りモード
  カイラル偏極度: 1.000000
  期待値: +1.0
  判定: PASS

Test 2: 理想的な右回りモード
  カイラル偏極度: -1.000000
  期待値: -1.0
  判定: PASS

Test 3: 線偏光モード
  カイラル偏極度: 0.000000
  期待値: 0.0
  判定: PASS

Test 4: 2原子系（Mo + S）
  カイラル偏極度: 0.968254
  Mo原子のみが寄与

軌道可視化を生成中...
複数モードの解析例...

解析完了！

実データへの適用

VASP計算結果（OUTCAR）からのデータ読み込み例：

def read_vasp_eigenvectors(outcar_file: str) -> List[PhononMode]:
    """
    VASPのOUTCARからフォノン固有ベクトルを読み込み

    Returns:
    --------
    modes : list of PhononMode
    """
    modes = []
    # ファイル解析ロジック（省略）
    # VASP出力形式に応じて実装
    return modes

# 実データ解析
modes = read_vasp_eigenvectors('OUTCAR')
analyze_phonon_modes(modes)

演習問題