高次の電気伝導度の応答計算(理論パート１)

https://qiita.com/hiro949/items/bdb5516fbd9a3b8806b6

では外場の大きさについて１次で近似し、材料の線形応答係数を計算しました。ここではさらに高次の項を計算する方法について紹介します。長いので、まずは理論の概論からスタートします。

# 外場項の取り扱い
システムの１体ハミルトニアンを$\hat{H}_0$としたとき、外場項$\hat{V}_E(t)$を

```math
\hat{V}_E(t) = \sum_{n=1}^{\infty}\frac{1}{n!}\left[ \frac{e}{\hbar} \mathbf{A}(t) \cdot \hat{ \mathbf{D} } \right]^n \hat{H}_0
```

と表します。ここで$\hat{ \mathbf{D} }$が微分で、$\mathbf{A}(t)$はベクトルポテンシャルです。ベクトルポテンシャルは電場と

```math
\begin{align}
\mathbf{E}(t) &= -\frac{\partial \mathbf{A} }{\partial t} \\
\mathbf{E}(\omega) &= -i\omega \mathbf{A}(\omega)
\end{align}
```

の関係があります。2行目はフーリエ変換したものです。よって

```math
\mathbf{A}(t) = \int \frac{d\omega}{2\pi} e^{i\omega t} \mathbf{A}(\omega) =  \int \frac{d\omega}{2\pi} \frac{i}{\omega} e^{-i\omega t} \mathbf{E}(\omega)
```

です。

また、微分は$H_0$に対して作用する演算子です。$H_0$がフェルミオンの生成、消滅演算子を使って

```math
\hat{H}_0 = \int d\mathbf{k} \sum_{pq}h^0_{pq}(\mathbf{k}) \hat{c}_p^\dagger(\mathbf{k}) \hat{c}_q(\mathbf{k})
```

と表されるとき、

```math
\hat{\mathbf{D}}h^0_{pq}(\mathbf{k}) = \frac{\partial{h^0_{pq}}}{\partial \mathbf{k}} 
```

と作用します。

そして全体のハミルトニアンを

```math
\begin{align}
\hat{H}_A(t) &= \hat{H}_0 + \hat{V}_E(t)
= \exp\left\{ \frac{e}{\hbar} \mathbf{A}(t) \cdot \hat{ \mathbf{D} } \right\} \hat{H}_0 \\
&= \hat{W}_\mathbf{A}(t) \hat{H}_0 \\
\hat{W}_\mathbf{A}(t) &= \exp\left\{ \frac{e}{\hbar} \mathbf{A}(t) \cdot \hat{ \mathbf{D} } \right\}
\end{align}
```

とします。
なぜなら、$\mathbf{\Pi}=\mathbf{k} + \frac{e}{\hbar}\mathbf{A}$について、$\mathbf{\Pi}=\mathbf{k}$周りでテイラー展開すると

```math
\hat{H}_\mathbf{A} = \hat{H}_0(\mathbf{k} + \frac{e}{\hbar}\mathbf{A}) = \sum_n \frac{1}{n!} \left(\frac{e}{\hbar}\mathbf{A} \cdot \nabla_\mathbf{k} \right)^n 
\hat{H}_0(\mathbf{k}) = \hat{W}_\mathbf{A} \hat{H}_0
```

となり、$\hat{W}_\mathbf{A}$が$\hat{H}_0(\mathbf{k}) \rightarrow \hat{H}_0(\mathbf{k}+\frac{e}{\hbar}\mathbf{A})$の変換になっているからです。　



# 電流演算子
電流の熱力学的平均値を考えます。まず、分配関数を考えます。

```math
\begin{align}
\hat{S}_c(E) &=  T\exp\left( -i \int dtV_E(t) \right) \\
Z(\mathbf{E}) &= \langle \hat{S}_c(E) \rangle_0 = \frac{Tr[S_c e^{-\beta H_0}]}{Tr[e^{-\beta H_0}]}
\end{align}
```

これを用いて、電流の各成分の熱力学的平均値は(ギリシャ文字)$=x,y,z$として、

```math 
\langle \hat{J}_\mathbf{E}^\mu(t) \rangle = \frac{\langle \hat{J}_\mathbf{E}^\mu(t)\hat{S}_c(E) \rangle_0 }{Z(E)} = \langle \hat{J}_\mathbf{E}^\mu(t) \hat{S}_c \rangle_c
```

となります。ここで、$\langle * \rangle_c$はファインマンダイヤグラムでconnectedダイヤグラムだけ考えることに対応しています。詳しくは以下の記事の補足を参照してください。

https://repository.kulib.kyoto-u.ac.jp/dspace/bitstream/2433/142612/1/KJ00004907236.pdf


具体的には電子の速度演算子に電荷$-e$をかけて

```math
\hat{\mathbf{J}}_\mathbf{E}(t) = -\frac{e}{\hbar} \hat{ \mathbf{D} }\hat{H}_A
= -\frac{e}{\hbar} \hat{ \mathbf{D} } \hat{W}_\mathbf{A}(t) \hat{H}_0 
```
とします。簡単のために以下、$e/\hbar=1$と置くと、

```math
\begin{align}
\hat{W}_\mathbf{A}(t) &= \exp\left\{ \mathbf{A}(t) \cdot \hat{ \mathbf{D} } \right\} \\
\hat{\mathbf{J}}_\mathbf{E}(t) &= - \hat{ \mathbf{D} } \hat{W}_\mathbf{A}(t) \hat{H}_0 
= -\frac{\delta \hat{W}_\mathbf{A}(t)}{\delta \mathbf{A}(t)} \hat{H}_0
= -\frac{\delta \hat{H}_A(t)}{\delta \mathbf{A}(t)} 
= -\frac{\delta \hat{V}_\mathbf{E}(t)}{\delta \mathbf{A}(t)}
\end{align}
```

です。汎関数微分の公式

```math
\frac{\delta F\left[ f,f',\cdots, f^{(n)} \right]}{\delta f} = \sum_{n=0}^\infty (-1)^n \frac{ \partial F }{\partial f^{(n)}}
``` 

を使って具体的な計算を行います。

### 電流の熱力学的平均値

次にこの熱力学的平均値を計算します。平均値は

```math
\langle \hat{\mathbf{J}}_\mathbf{E}(t) \rangle = \langle \hat{\mathbf{J}}_\mathbf{E}(t) \hat{S}_c \rangle_c = - \left\langle \frac{\delta \hat{V}_\mathbf{E}(t)}{\delta \mathbf{A}(t)} \hat{S}_c \right\rangle_c 
= - \left\langle \frac{\delta \hat{S}_c(t)}{\delta \mathbf{A}(t)} \right\rangle_c
= - \frac{\delta}{\delta \mathbf{A}(t)} \left\langle \hat{S}_c(t) \right\rangle_c
```

となります。

# 電気伝導度の計算

## 下準備

まずは下準備の計算を２つします。

### 下準備１
```math
\begin{align}
\frac{\delta}{\delta E_\mu(t) } &= \int d\omega e^{-i\omega t} \frac{\delta}{\delta E_\mu(\omega) } \\
\frac{\delta}{\delta E_\mu(\omega) } &= \int \frac{dt}{2\pi} e^{i\omega t} \frac{\delta}{\delta E_\mu(t) }
\end{align}
```

以下、これらの式を導出します。普通の微分の変数変換と同様にします。

```math
\begin{align}
\frac{\delta}{\delta E_\mu(t) } 
&= \int d\omega \frac{\delta E_\nu(\omega) }{\delta E_\mu(t)} \frac{\delta}{\delta E_\nu(\omega) }
= \int d\omega \int dt' e^{-i\omega t'} \frac{\delta E_\nu(t') }{\delta E_\mu(t)} \frac{\delta}{\delta E_\nu(\omega) } \\
 &= \int d\omega \int dt' e^{-i\omega t'} \delta(t-t') \delta^\mu_\nu \frac{\delta}{\delta E_\nu(\omega) }\\
&= \int d\omega e^{-i\omega t} \frac{\delta}{\delta E_\mu(\omega) }
\end{align}
```

```math
\begin{align}
\frac{\delta}{\delta E_\mu(\omega) } 
&= \int dt \frac{\delta E_\nu(t) }{\delta E_\mu(\omega)} \frac{\delta}{\delta E_\nu(t) }
= \int dt \int \frac{d\omega'}{2\pi} e^{i\omega' t} \frac{\delta E_\nu(\omega') }{\delta E_\mu(\omega)} \frac{\delta}{\delta E_\nu(t) } \\
&= \int \frac{dt}{2\pi} e^{i\omega t} \frac{\delta}{\delta E_\mu(t) }
\end{align}
```

です。

### 下準備２

$\mathbf{E}(t)$を$\mathbf{A}(t)$で汎関数微分すると

```math
\begin{align}
\frac{\delta}{\delta A_\mu(t)} &= \int dt' \frac{\delta E_\nu(t')}{\delta A_\mu(t)} \frac{\delta}{\delta E_\nu(t')} = - \int dt' \frac{d}{dt} \frac{\delta E_\nu(t')}{\delta \dot{A}_\mu(t)} \frac{\delta}{\delta E_\nu(t')} \\
&= \int dt' \frac{d}{dt} \frac{\delta E_\nu(t')}{\delta E_\mu(t)} \frac{\delta}{\delta E_\nu(t')} = \int dt' \frac{d}{dt} \delta(t-t')\frac{\delta}{\delta E_\mu(t')} \\
&= - \int dt' \delta(t-t') \frac{d}{dt} \frac{\delta}{\delta E_\nu(t)} = - \frac{d}{dt} \frac{\delta}{\delta E_\mu(t)} \\
&= - \frac{d}{dt} \int d \omega e^{-i\omega t} \frac{\delta}{\delta E_\mu(\omega)}
\end{align}
```

となります。$d/dt$は$\delta/\delta E_\mu(\omega)$が作用する関数も含めた後ろの積分全体にかかっているので注意してください。よって

```math
\hat{D}^\mu \hat{W}_\mathbf{A}(t) = \frac{\delta \hat{W}_\mathbf{A}(t)}{\delta A_\mu(t)} =  - \frac{d}{dt} \int d \omega e^{-i\omega t} \frac{\delta \hat{W}_\mathbf{A}(t)}{\delta E_\mu(\omega)}
```

```math
\begin{align}
\frac{\delta \hat{W}_\mathbf{A}(t) }{\delta E_\mu(t')} &= \frac{\delta \mathbf{A} (t) }{\delta E_\mu(t')} \cdot \hat{ \mathbf{D} } \exp\left\{ \mathbf{A}(t) \cdot \hat{ \mathbf{D} } \right\} \\
&= \int d\omega \frac{\delta E_\nu(\omega) }{\delta E_\mu(t')} \frac{\delta \mathbf{A}(t) }{\delta E_\nu(\omega)} \cdot \hat{ \mathbf{D} } \exp\left\{ \mathbf{A}(t) \cdot \hat{ \mathbf{D} } \right\} \\
&= \int d\omega e^{-i \omega t'} \frac{i}{2\pi\omega} e^{i\omega t} \hat{D}_\mu \exp\left\{ \mathbf{A}(t) \cdot \hat{ \mathbf{D} } \right\} \\
&= \int \frac{d\omega}{2\pi}  \frac{i e^{i \omega (t-t')}}{\omega} \hat{D}^\mu \hat{W}_\mathbf{A}(t) 
\end{align}
```



## 電気伝導度の計算

以上より

```math
\begin{align}
\sigma(t;t_1\cdots t_n) &= \left( \frac{\delta^n \langle \hat{J}^\mu_\mathbf{E}(t) \rangle }{\delta E_{\alpha_1}(t_1) \cdots \delta E_{\alpha_n}(t_n) } \right)_{\mathbf{E} \rightarrow 0}
= - \left( \frac{\delta^{n+1} \left\langle \hat{S}_c(t) \right\rangle_c }{\delta A_\mu(t) \delta E_{\alpha_1}(t_1) \cdots \delta E_{\alpha_n}(t_n) } \right)_{\mathbf{E} \rightarrow 0} \\
&= \sigma_\mathbf{E}(t;t_1\cdots t_n)|_{\mathbf{E} \rightarrow 0 } 
\end{align}
```

となります。

また、　

$n$次の電気伝導度テンソルの成分を$\sigma^{\mu\alpha_1\cdots \alpha_n}(\omega;\omega_1,\cdots,\omega_n)$とすると、

```math
\langle \hat{J}_\mathbf{E}^\mu(\omega) \rangle = \sum_{n=1}^{\infty} \int \prod_{i=1}^{n} d\omega_i \sigma^{\mu\alpha_1\cdots \alpha_n}(\omega;\omega_1,\cdots,\omega_n) \prod_{i=1}^{n} E_{\alpha_i}(\omega_i)
```

なので、電場について凡関数微分して

```math
\begin{align}
\sigma^{\mu\alpha_1\cdots \alpha_n}(\omega;\omega_1,\cdots,\omega_n) &= \left( \frac{\delta^n \langle \hat{J}_\mathbf{E}^\mu(\omega) \rangle }{\delta E_{\alpha_1}(\omega_1)\cdots \delta E_{\alpha_n}(\omega_n)}\right)_{\mathbf{E}\rightarrow 0} \\
&= \int\cdots \int \frac{dt}{2\pi} e^{i\omega t} \prod_k \frac{dt_k}{2\pi} e^{i\omega_k t_k} \left( \frac{\delta^n \langle \hat{J}^\mu_\mathbf{E}(t) \rangle }{\delta E_{\alpha_1}(t_1) \cdots \delta E_{\alpha_n}(t_n) } \right)_{\mathbf{E} \rightarrow 0} \\
&= \int\cdots \int \frac{dt}{2\pi} e^{i\omega t} \prod_k \frac{dt_k}{2\pi} e^{i\omega_k t_k} \sigma_{\mathbf{E} \rightarrow 0}(t;t_1\cdots t_n)
\end{align}
```

となります。 

具体的に計算するときは次のように式変形します。

```math
\begin{align}
\sigma^{\mu\alpha_1\cdots \alpha_n}(\omega;\omega_1,\cdots,\omega_n) &= \int \frac{dt}{2\pi} e^{i\omega t} \frac{d}{dt} \int d \omega' e^{-i\omega' t} \left( \frac{\delta^{n+1} \left\langle \hat{S}_c(E) \right\rangle_c }{\delta E_\mu(\omega') \delta E_{\alpha_1}(\omega_1) \cdots \delta E_{\alpha_n}(\omega_n) } \right)_{\mathbf{E} \rightarrow 0} \\
&= - \int d \omega' i\omega \int \frac{dt}{2\pi} e^{i(\omega-\omega')t} \left( \frac{\delta^{n+1} \left\langle \hat{S}_c(E) \right\rangle_c }{\delta E_\mu(\omega') \delta E_{\alpha_1}(\omega_1) \cdots \delta E_{\alpha_n}(\omega_n) } \right)_{\mathbf{E} \rightarrow 0} \\ 
&= - i\omega \int d \omega' \delta(\omega-\omega') \left( \frac{\delta^{n+1} \left\langle \hat{S}_c(E) \right\rangle_c }{\delta E_\mu(\omega') \delta E_{\alpha_1}(\omega_1) \cdots \delta E_{\alpha_n}(\omega_n) } \right)_{\mathbf{E} \rightarrow 0} \\
&= - i\omega \left( \frac{\delta^{n+1} \left\langle \hat{S}_c(E) \right\rangle_c }{\delta E_\mu(\omega) \delta E_{\alpha_1}(\omega_1) \cdots \delta E_{\alpha_n}(\omega_n) } \right)_{\mathbf{E} \rightarrow 0}
\end{align}
```

この式は最適化手法でいうところのいわゆる鞍点法という近似になっています。

今回はとりあえずここまでにします。次の記事で具体的な計算の仕方をご紹介します。

https://qiita.com/hiro949/items/175cbfc3654ba3186ce8

# コメント
ご覧の通り、かなり複雑な計算なので計算ミスがあるかもしれません。発見した際はご指摘していただけると幸いです。

# 参考文献
- https://journals.aps.org/prb/abstract/10.1103/PhysRevB.99.045121
- https://journals.aps.org/prb/pdf/10.1103/PhysRevB.61.5337
- https://web.physics.ucsb.edu/~mark/ms-qft-DRAFT.pdf
- http://www-hep.phys.s.u-tokyo.ac.jp/~hama/lectures/lecture.html