大規模並列計算による多電子系非線形ダイナミクスの量子的記述―分子から固体まで―
乙部智仁 1、矢花一浩 2、岩田潤一 2、中務孝 3、山極満 1、G.F.Bertsch4
１

日本原子力機構、2 筑波大計算センター、3 理化学研究所、4 ワシントン大
メールアドレス、otobe.tomohito@jaea.go.jp

概要： これまでの解析的手法では分からなかった電子ダイナミクスが並列計算による量子シミュレーションを行うことで明らかにできる可能性がある。我々は時間依存密度汎関数法（TD-DFT）の実時間実空間解法による大規模計算を行うことで強レーザー場中にある分子のトンネルイオン化と固体の絶縁破壊過程の第一原理シミュレーションを行った。

1

はじめに科学シミュレーションにおいて多粒子系の記ギー項を計算する。また時間発展はテイラー展開の 4 次で近似する[2]。

述は大きな問題である。特に量子多体系の精密トンネルイオン化に対しては外場を静電場とな記述は現在不可能であり何らかの近似が必要することで静的な Kohn-Sham 方程式を用い、である。密度汎関数法（DFT）は多電子系の基トンネルイオン化率を外向き波の境界条件を満底状態を高い精度で再現する理論の1つである。

たす状態の複素固有値の虚数部分から計算する

近年DFTを時間変化する系に拡張した時間依存[3,4]。電子の軌道について並列化したコードを密度汎関数法（TD-DFT）[1]による電子ダイナ作成し、効率のよい計算を可能とした。外向きミクスを記述する研究が注目されている。

波の境界条件は吸収境界条件により擬似的に満TD-DFTはこれまで線形領域の記述において成たす。このため分子の外側の非常に大きな空間果を上げており、現在非線形領域への応用とそが必要となる。の有用性が研究されている。

固体電子は単位胞のみについて実空間法を用

また近年強レーザー場中にある多電子系の非いBlochの定理から無限周期系を考える。また空線形現象の研究が盛んである。その理解には量間的に一様なベクトルポテンシャルにより結晶子的シミュレーションが必要とされているが、中の電場を扱った（式2）。

計算方法の確立と非常に大きな計算資源が必要であるため解析的手法や1，2電子系の記述が主である。我々はTD-DFTの基礎方程式である時間依存Khon-Sham方程式を実時間実空間法を用いK∂ψ i , k ( r , t )2⎡ 1 ⎛e K ⎞K= ⎢−⎜ =∇ + A(t ) ⎟ + Vion ( r ) +c∂t⎠⎢⎣ 2m ⎝KK ρ (r , t )K ⎤Ke 2 ∫ dr ′ K K + Vxc [ ρ ( r , t )]⎥ψ i.k ( r , t )r − r′⎦

i=(2)にある分子のトンネルイオン化率とダイアモン

K
Vion (r ) は イ オ ン か ら の ポ テ ン シ ャ ル 、KVxc [ ρ (r , t )] は相関交換ポテンシャルである。波ドの絶縁破壊過程について大規模並列計算を行

動関数の添え字 k は Bloch の位相を示す。このうことで調べた。

ときのベクトルポテンシャル A(t ) はレーザー場

計算方法( Alaser (t ) )と共に表面電荷による電場( Aind (t ) )を含て計算するアプローチにより、強レーザー場中

2KKKKKK

計算は時間依存密度汎関数法の基礎方程式であ。計算は Blochんでいる（ A(t ) = Alaser (t ) + Aind (t ) ）る時間依存 Khon-Sham（TD-KS）方程式を実の位相のサンプリングｋ（1152 点）について並時間実空間法を用いて計算する。具体的にはハ列化し計算を行った。

ミルトニアン及び波動関数を 3 次元メッシュ上3の値として近似し、9 点差分により運動エネル

3.1 分子のトンネルイオン化率

計算結果

図1に最外殻分子軌道（HOMO）とHOMO-1

ザー強度ではレーザー通過後も電子振動が長くの静電場中の電子密度分布図を示した。本来電のこる結果となった。この振動は伝導帯へ励起場と逆方向に分極するはずの電子密度が電場方した電子が自由電子として振舞うことによるプ向に動くことが分かる。この計算から分子軌道ラズマ振動数であり、レーザー振動数と励起電毎の分極の特徴は単純ではなく多体効果の重要子のプラズマ振動数が近い時に誘起されること性を示す結果となった。またイオン化率の角度が分かった。このようなダイナミクスの詳細は依存性やσ軌道であるHOMO-1の影響を明ら微視的シミュレーションを行うことで初めて明かにした[4]。これはこれまでのHOMOだけを考らかとなった。

慮した計算が不十分であることを示している。

4 まとめ
強レーザー場中にある多電子系の非線形ダイナミクスについて実時間実空間法を用いた T計算法とその結果を報告した。我々の計算は多体効果を取り入れた第一原理計算であり、多電子系のフェムト秒領域のダイナミクスをシミュレートすることが可能である。
原子分子、クラスターといった孤立系から固体等の無限系までを統一的に記述可能であり、今までに無い数値計算法である。次世代の計算機によりX線自由電子レーザーによる生体分子

図1 ベンゼンの最外殻分子軌道(HOMO：左)とさらに深く束縛された軌道(HOMO-1：右)の静電場中での電子密度分布。

3.2 透明素材の絶縁破壊過程
レーザーによる誘電破壊過程は電子励起の確率や電子―電子衝突過程の重要性等不明な点が多い。量子シミュレーションにより微視的な高速現象の解明が可能となる。の損傷シミュレーションや表面や薄膜、アモルファスなど複雑でより実際に即した物質の加工・制御の最適化研究の可能性があり、現在準備的研究を進めている。
このような大規模な量子計算の可能性はより高性能な計算機の登場によりさらに高まっていくものと思われる。
参考文献[1] E.Runge and E.K.U. Gross,[Density-Functional
Theory for Time-Dependent System]、Phys.
Rev. Lett., 52, pp997、1984 年[2] K.Yabana and G.F.Bertsch, [Time-Dependent

時間発展
図 2 パルスレーザー場中にあるダイアモンドの電子ダイナミクス。基底状態からの密度のずれをしめした。パルス初期（左）、レーザーのピーク（中）、レーザー通過後（右）。

local-density approximation in real time], Phys.
Rev. B, 54, pp4484, 1996 年[3] T.Otobe, et al, [First-Principle
calculations for the tunnel ionization rate ofatoms and molecules] Phys.Rev. A, 69,

図2に時間幅16fsのレーザーをダイアモンド

pp053404, 2004 年に照射した時の電子密度の規定状態からのずれ[4] T.Otobe and K.Yabana, [Density-Functionalをしめした。時間と共に励起電子が増加し結晶

calculation for the tunnel ionization rate of内に広がっていく様子が分かる。またあるレーhydrocarbon molecules] Phys. Rev. A, 75,pp062507, 2007 年