Journal of Computer Chemistry, Japan
Online ISSN : 1347-3824
Print ISSN : 1347-1767
ISSN-L : 1347-1767
研究論文
計算科学的手法による物質中におけるミュオン位置同定
渡邊 功雄Sari Dita PuspitaRamadhan RedoAdiperdana BudiSulaiman Shukri
著者情報
ジャーナル フリー HTML

2020 年 19 巻 3 号 p. 115-123

詳細
Abstract

素粒子であるミュオンを物質科学研究に利用するミュオンスピン緩和法(μSR)が開発されてから半世紀が過ぎようとしている.これまで数多くの研究成果が報告されるとともに幅広い研究分野に応用されている.近年,有機分子性物質にける磁性体・超伝導体が数多く開発されるに伴い,μSRの化学分野における利用もポピュラーになりつつある.一方で,μSRの利用が広がるとともに,そこからより詳しい情報を得たいという要求も増大し,計算科学に関するアクティビティも生まれつつある.本稿では,その活動の一環であるミュオン位置計算に関して,最近の著者のスーパーコンピュータを用いたミュオン位置同定に関する研究活動を概観する.

ミュオンは素粒子の一つである.この素粒子が物質科学研究への応用が可能あることが示されて半世紀がたとうとしている.ミュオンは目に見えないものであるが,素粒子物理学の知識よりその素性は明らかになっており,おおまかではあるが,水素の約1/9,しかしながら電子の約200倍の質量をもった水素のアバター的粒子である.このミュオンを測定試料外部より打ち込み,物質の電子状態を探る微視的研究手法をミュオンスピン緩和法(μSR)という [1].ミュオンは自然界には安易(安定)に存在しておらず,地球上では毎秒手のひらに一個程度のミュオンが天空から降り注ぐ程度である.そして寿命2.2マイクロ秒で消失する.このため,実際の物質科学研究にミュオンを用いる際には人工的にミュオンを生成する加速器が必要であり [2],結果,やや敷居が高い測定手法と思われているかもしれない.しかしながら,その応用範囲はとても広く,固体物理学を中心として化学分野や生物学分野においてもその活用が広まり,特に近年では化学分野における応用が注目を集めている.その主たるターゲットは有機分子性磁性体や超伝導体である.これら一連の有機分子性物質は,化学的合成手法によって高純度な試料の合成ができることに加え,自在な分子修飾や元素置換により細かく電子状態を制御する分子設計を可能にする.複雑な分子構造を持つにも関わらず,モデル化した単純な電子状態を用いて学理課題が議論できる分子性物質は,μSRを応用するには格好のターゲットである.

ミュオンの応用において最たる問題は,その物質中の位置である.ミュオンは物質に停止した後に周囲の電子状態を反映した物理情報を与える.例えば,周囲の磁気スピンがどのように秩序化するか,ミュオンがいかに化学的分子性結合を行うか,という情報である.どちらも原子レベルの情報であり,その情報を正確にひも解くことがミュオンの応用として極めて重要である.もし物質中におけるミュオン位置情報を正しく得ることができれば,μSR実験データと突き合わせてより詳細に物質の電子状態を議論することができる.しかしながら実際には,目で見えないミュオンを試料外部から投げ入れてミュオンに任せて好ましい位置で自発的に停止してもらうため,その位置を知ることは決して容易なことではない.そこで検討されているのが計算科学的研究手法である.

本稿においては,あくまでも筆者の研究の限られた研究の範囲内ではあるが,ミュオン位置決定をどのように行うことができるか,また,その手法のためにはどのような計算機科学的努力がなされてきたいかということに関しての報告を行う.これまでの我々の研究アクティビティを外観し,どのような研究手法が求められているか,また,それによりどのような効果を得ることができるか.という点が読者の皆様に伝われば本望であるといえる.

ミュオン位置を求めることに関して,筆者の知りえる範囲内では確立した手法はまだない.ある程度限られた条件下での議論が多く,計算的手法からはもっともな回答が得られたように思えるが,その解答が実際の実験結果をどれほどうまく説明したかどうかというと,結論としてはなかなか難しい.それほどにミュオンの位置を同定することはこれまで困難な作業であった.ミュオン位置が明確になったら物質科学研究上どのような有用な情報を得られるかというと,

① 磁気秩序状態におけるスピン構造が求められる.

② ミュオン近傍における微細な電子状態の時間的変化,空間的変化の情報を知ることができる.と大きく2点があげられる.これらを目的として研究においてはμSRとミュオン位置計算を通じて幅広く有益な成果を獲得することが期待できる.例えば,有機分子性磁性体はπ電子のスピンが持つ磁気秩序状態が主とした磁気構造をとることが考えられるが,π電子のほかに磁性原子を含む物質群がさらに磁気原因を複雑化させること,また,微量しか合成できないことや多くの水素を含む軽元素によって構築されるソフトな結晶構造を有するなど,様々な要因が研究を複雑化させる.磁気構造の同定には中性子散乱実験が極めて有効であることはよく知られているが,空間的に低いスピン密度と,水素原子による中性子散乱の高いバックグランドが小さい磁気散乱信号をマスクしてしまうというジレンマがある [3, 4].一方,μSRはそのような制限に関わることなく磁気状態に関する情報を提供するケースが多い.ここでもしミュオン位置が分かれば,双極子-双極子相互作用に基づいてどのような磁気スピン構造が発生しているかを議論することが可能になる.また,ミュオンが停止している近隣の電子伝達 [5, 6]やイオン拡散 [7, 8]などの動的振る舞いを定量的に解き明かすことも可能にする.これらの情報は他の測定手法では得難い場合も多く,μSRとスーパーコンピュータの活用を組み合わせた新しい研究手法としてμSR研究分野では近年重要な研究課題とされている.

最大の問題は,このミュオン位置計算の手法がまだ確立していないことにある.ミュオンは正電荷を持つ素粒子であり(負電荷のものもあるが,本稿では割愛する),基本的には結晶中で電場ポテンシャル的にもっとも安定な位置に停止するものと考えられている [1].このため,シンプルな電場ポテンシャル計算研究は過去に数多くチャレンジされている.筆者の知る限り一つの試料系でもっともミュオン位置計算が行われたものとして,La系高温超伝導酸化物がある.この物質は化学的元素置換によって試料中の電荷濃度を連続的に変化させることにより,電子状態が絶縁体反強磁性状態から超伝導状態へと変化する [9].酸化物高温超伝導体が発見されて30年あまり,その研究の歴史からミュオン位置計算は行われいてる [10,11,12,13,14,15,16,17,18,19].長距離磁気秩序状態におけるμSR測定からミュオン位置における自発的内部磁場の値が同定できるので,その値を再現できるミュオン位置を求めるということが基本的計算方針である.この方針は決して特殊なものではなく,現在においてもミュオン位置計算に対して第一義的に検討する事項となる.

Figure 1に,最も基本的なLa系酸化物高温超伝導体の母物質であるLa2CuO4に関し,計算的手法を用いて提案されたミュオン位置を示す [10,11,12,13,14,15].各研究はそれぞれ異なった計算方法があり,結果として提案されているミュオン位置は大変広範囲にばらついてしまい,どの結果も実験値をうまく説明できていないこの状態は現在まで続いており,ミュオン研究分野における一つの解決課題として残されている.これら一連の計算結果には二つの問題点がある.一つは,計算モデルが実際のミュオン実験の状況を再現していないことである.ミュオンは結晶中では大変希薄な磁性不純物として振る舞うため,ユニットセル1つにミュオンを仮想導入した単純なモデルでは試料全体のユニットセルにミュオンがあるという仮定になり,実際の実験条件を再現できない.このため,希薄不純物を含んだ大きなセル模型(スーパーセル)による大規模計算が必要である [19].二つ目の問題は,電子スピンとミュオン自身のもつ量子的効果が考慮されていないことにある.物質中の電子スピンは実際には空間的に広がる電子軌道を持つ [20].また,この電子軌道も周囲の電子雲との相互作用を通じてその形を大きく変える.同様にミュオン自体もその素粒子としての性格よりゼロ点振動による空間的広がりを持つ [21].これらの空間的広がりの効果をどのように取り込むかが問題となる.さらには,ミュオンは軽いとはいえ,水素原子のアバターでありミュオン近傍では自らの質量によって結晶構造を局所的に歪ませることが指摘されている [15].結果,ミュオン位置計算は次の命題として捉えなおすことができる.「希薄磁性不純物による局所的電子状態の変化を求める大規模計算」.このような計算は,CPUが向上したとはいえ我々が通常用いるパソコンやワークステーション程度の計算機ではたちうちできず,スーパーコンピュータを活用した大規模計算が必要になる.

Figure 1.

 Muon positions in La2CuO4 suggested in previous studies.

我々は,上記のような方針で理化学研究所情報システム本部の管理するスーパーコンピュータ(HOKUSAI:Figure 2)を活用し,La2CuO4に関するミュオン位置計算を進めている.HOKUSAIの性能等に関する情報は理化学研究所のホームページからご参照いただきたい(http://i.riken.jp/supercom/overview/).我々は自力開発のプログラムを持たないため,Vienna ab-initio Simulation Package (VASP)による密度汎関数法(DFT)を用いた計算を行っている.計算の手順として,まずユニットセルを用いて静電ポテンシャル計算からポテンシャルが極小になる点を求める.その点をミュオンが入射直後に停止する初期位置と設定する.ミュオンを計算中に取り込む際に軽い水素として扱い,質量を1/9に設定した水素の擬ポテンシャルを用いる.実際の実験条件である結晶中の超希薄不純物としてのミュオンの状態を実現するために,32ユニットセルを持つスーパーセルを設定し,その中心のユニットセル1つに水素の擬ポテンシャルとしてのミュオンを一つ配して,ミュオン位置の緩和や結晶構造の局所的変化,および影響を受ける周囲の電子状態の変化を探っている.最終的な結果からミュオン周囲での電場ポテンシャルを用いて数値的にミュオンのゼロ点振動による空間的広がりも計算し,それらを総合してμSRの実験結果を説明することを試みている.この計算過程においてミュオンのスピンの影響は考慮していない.基本的に希薄不純物であることに加え,結晶の原子点間に停止して周囲の電子と化学的結合を行わないかぎり,磁気双極子相互作用を主とした低エネルギーな磁気相互作用しか持たないためと考えている.ミュオンスピンの影響にかかる詳細な議論はなされていないように思えるが,実際の実験結果を見る限り高温超伝導酸化物中においてミュオンスピンが周囲の電子状態を強く変化させたという証拠はない.

Figure 2.

 RIKEN supercomputing system (HOKUSAI).

Figure 3に近年得られた成果の一部を示す.これは静電ポテンシャル計算から得られたLa2CuO4に期待できる3つのミュオン位置(図中の黒マーク)に関して,ミュオンが存在した場合の最近接にあるCuの電子状態の変化を求めたものである.図はCuO6面におけるスピン密度の相対変化を表しており,ミュオンがない場合からミュオンによってどのようにスピン密度に変化が生じるかという結果を,予想されるミュオン位置であるM1,M2およびM3について個別に表示したものである.図からわかるように,ミュオンがあるとないとではCuの電子状態が大きく変わることがわかる.この変化は極めて局所的であり,ミュオンから離れるにしたがってユニットセル単位で急速に変化が失われる.計算の収束状態などの詳細は拙著をご参照いただきたい [17,18,19].ちなみに電子相関関数はGGA-PW91を用い,Cu位置に電子間相互作用のクーロン反発力を加えている.この電子相関関数は,Cuなどの3d電子を含む絶縁体などの電子状態をよく再現できる関数として知られている.当然ではあるが,どの電子相関関数を選ぶかということは,計算規模に関わらず実際の電子状態を説明する上で重要である.近年,より多くの電子相関関数が提案されており,その関数がターゲットする試料系にあっているかどうかは,関数形を変えた同様な計算の繰り返しによる検証が必要になる.このような計算の繰り返しのためにも大規模な高速計算の繰り返しが必要になり,高速・大容量のスーパーコンピュータの利用が必須となる.我々の計算例ではあるが,この大規模計算一回に対して,192コアと約90 GBメモリを活用し,1計算を大体48∼72時間程度で完了させている.CPUはIntel XEON Gold 6148である.当然ながら計算条件を変えることによってこの時間は変化する.しかしながら,実際の研究を進める上では十分容認できる時間要素であるが,逆にHOKUSAIを活用せずして本研究は成り立たない.この銅酸化物に関する研究成果を取りまとめた論文は現在査読中であるが,最終的な結果として我々の計算結果がμSR測定から得られた時間スペクトルをよく再現することが確認された [22].

Figure 3.

 Our suggested muon positions in La2CuO4. The upper panel shows distributions of the spin density of the Cu atom within the CuO6 plane around the most possible three candidate positions for the muon. Relative changes in the spin density are described by the color contour. Red and blue colors indicate the high- and low-spin densities, respectively. The left Figure shows the case without the muon and others from the left to right indicate cases with the muon at M1, M2 and M3. The lower panel shows the muon positions. Each Figure corresponds to the one above in the upper panel. The CuO6 unit which compose the basic structure of La2CuO4 is shown. Red and blue balls indicate O and Cu, respectively [22].

このLa系酸化物高温超伝導体の経験を基に,次のステップとして有機分子性磁性体・超伝導体にこの計算手法活用してその電子状態を解決することに取り組んでいる.上述したように,有機分子性物質は最近のμSR研究の重要なターゲットである.特に,磁気スピンの起源や磁気秩序状態のスピン構造が当面の主体研究課題である.これらの問題解決に関しては成功した例が極めて少なく,μSRとDFT計算を併用した磁気スピン状態へのアプローチが強く求められている.しかし,有機分子性物質には酸化物とは異なった特有の問題がある.第一に,計算対象となる試料の分子量がとても多いこと.これにより,結晶構造や電子状態の最適化に膨大な計算量が必要になる.また,酸化物では可能であったスーパーセルの取り扱いが難しくなる.第二は,初期条件としての電子スピンの取り扱いである.酸化物の場合には局在化した電子スピンを格子点に配すれば必要十分な初期条件となる.一方,有機分子性物質は空間的に大きく広がったπ電子に起因する磁気スピンが磁性の主役になることが多く,その計算上の取り扱いが難しい.また基本構造を担う分子がダイマー化することにより,スピン密度中心をどこするかという問題がさらに複雑になる.VASPにおいては,初期条件はあくまで格子点での電子スピンの設置になる一方,最終的な結果は分子上に大きく広がった電子雲状態で全体のスピン状態が議論されることが望ましい.VASPがどこまでこの条件を満たしてくるかどうかという点に関しては,まだまだ検証が必要である.もともとVASPは有機分子性物質向けに開発されたプログラムではないため,HOKUSAIに限らず他のスーパーコンピュータ上でも安定に稼働する他のオープンソースを活用した計算による補完作業が今後求められるかもしれない.

これらの問題をかかえつつも,我々はまずもって有機分子性物質においてもミュオン位置を決定するプロセスは一緒と考え,第一ステップとしてのユニットセルを用いた静電ポテンシャル計算からスタートし,ミュオンが停止すると期待される位置の初期状態を求めた.我々グループは,この方針のもとに典型的な有機磁性体であるλ-(BEDT-STF)2FeCl4 [23]に関してのμSRとDFT研究を継続している.この有機磁性体はダイマーを作るBEDT-STFの2分子上に広がるπスピンに起因する磁気モーメントとFeの持つ磁気モーメントの2種類の磁気スピンがある [23,24,25,26,27].λ-(BEDT-STF)2FeCl4では約15K近傍で自発的な磁気秩序化が発生するが,どちらのスピンがどのような磁気構造をもって秩序化するかは未だ自明ではない.また,隣接するπ電子とFe電子間に期待される磁気相互作用とその磁気秩序状態に対する役割など,新しい磁性状態を含む大変興味ある有機磁性体である [23,24,25,26,27].さらに,高圧下で磁気秩序状態が失われて約4K以下で超伝導が発現することも報告されており [28],磁性と超伝導の相関を調べる上からも注目を集めている.

Figure 4に今回研究のターゲットであるλ-(BEDT-STF)2FeCl4の分子構造を示す [23].BEDT-STF分子がλ型に積層することによって2次元性の強い結晶構造を持つ.ここで,Seの原子位置に関してはBEDT-STFの一分子内における配位にランダム性があり,Sとの入れ替えが起こることが示されている [23].2つのBEDT-STF分子がダイマーを形成し,お互いのπ電子軌道を共有し2つのBEDT-STF分子で1つの電子を共有する.つまり,2分子に広がる広大な軌道全体がπ電子スピンの存在領域となる.この効果をDFT計算でどのように取り入れるか,もしくはVASPの計算アルゴリズムに任せるのかは計算する側の熟考が必要である.一方,ダイマー間にはFeイオンが配位し,π電子軌道状の磁気スピンに加えてもう一つのFeの磁気スピンが局在する.結果として,λ-(BEDT-STF)2FeCl4には2種類の磁気スピンが存在し,それぞれの秩序化状態とお互いの磁気相関を解き明かすことが電子状態を総合的に理解する上で重要な課題である.我々はλ-(BEDT-STF)2FeCl4に関するμSR測定を行い,長距離秩序状態が発現することでミュオン位置に発生する周囲の磁気スピンからの内部磁場を決定した.この内部磁場がπ電子もしくはFeのどちらに起因するのかはμSRデータのみでは判別ができず,そこでDFT計算の出番となる.

Figure 4.

 (left) Molecular structure of λ-(BEDT-STF)2FeCl4. Two BEDT-STF molecules stack and form a dimer state (red area). The FeCl4 molecule attaches to the dimer. (right) BEDT-STF molecule, indicating C (brown), Se (green), S (yellow) and H (blue).

まず静電ポテンシャル計算を実行する.静電ポテンシャル計算としてはマーデルングポテンシャル計算が知られている.しかし,過去の経験ではこのポテンシャル計算でミュオン位置の推察に成功した例はなく,擬ポテンシャルを活用した静電計算を行わなければならない.Figure 5にVASPによる静電ポテンシャル計算結果を示す.ユニットセルを細かいグリッドに分割し,各点における静電ポテンシャルを繰り返し計算する.擬ポテンシャルを用いることに伴って指定する電子相関関数はGGAを用いた.6×5×2のk点サンプリングを行い,500meVに計算のカットオフをかけている.それでも計算規模は小さく,HOKUSAIを活用すれば30分以下の計算で完了する.ミュオンは正電荷を持つために,基本的には負電荷イオンを持つ元素に近いところを好む.計算結果では,2次元面上のBEDT-STF分子間のFeCl4分子近傍に一つ極小ポテンシャルを確認した.面白いことに,さらにもう1点,BEDT-STF分子の端の方のH位置付近にも同様な極小ポテンシャルを持つ点が確認された.静電ポテンシャル的観点から見ると,λ-(BEDT-STF)2FeCl4には2つのミュオン位置が可能である.前者はよりFeに近い電子情報を,後者はダイマー上のπ電子の情報を見やすいことが推察される.

Figure 5.

 Map of the static potential obtained on λ-(BEDT-STF)2FeCl4 by using VASP. The isosurface area shows the 100 meV higher region from the local minimum potential. (left) The whole view of the λ-(BEDT-STF)2FeCl4 crystal. The local minimum static potential region can be seen near FeCl4 and H. (right) Expanded views of the minimum potential regions.

これらを始点として次の計算を進めていくことになる.我々が他の有機磁性体における静電ポテンシャルを調べる限り,ミュオンがダイマー間,もしくはダイマーの端の方に付着しやすいという傾向は,似たような結晶系を持つ有機磁性体・超伝導体には共通の現象である.一方,他の有機磁分子性物質において,DFT計算は用いずに実験結果を解析することでミュオンが分子面内に停止する傾向があるという研究結果もあり [29],今後の我々の研究成果との比較に興味が持たれる.

我々はこの時点まで計算が完了している.本稿が読者によって読まれる頃にはもう少し計算が進められているかと期待する.ただし,冒頭にも書いたように,有機分子性物質は多くの原子によるユニットセルを持つために,どこまで大型のスーパーセルによるミュオン位置計算を進めることができるかどうか,という点が今後の重要な課題になる.スーパーセルを用いた希薄不純物計算という方針は,実際のμSRの実験条件を再現するという点から避けがたい条件である.できればより多くのスーパーセルを用いた計算が望まれるが,自分がアクセスできる計算機性能によって大きく律速される.いかに大容量・高速のスーパーコンピュータを利用できるかどうかが鍵となる.

今後の計算方針は次のような手順になる.

① 静電ポテンシャル計算から得られた初期停止にミュオンを置き,結晶全体の最適化を行う.

② 磁気スピンを仮定して電子状態の最適化を行う.

③ ミュオンのゼロ点振動を求める.

④ 双極子-双極子計算から,磁気スピンの構造の同定を行う.

プロセスを列挙することは簡単であるが,実際の計算はなかなかうまくいかない.かくいう筆者グループも試験的にミュオンを含んだλ-(BEDT-STF)2FeCl4の結晶構造最適化の計算を実行してみた.しかしながら,結晶自体がソフトな有機分子性物質の特性により,結晶が大きく歪む場合もあり,なかなかうまく最適化できないのが現状である.また,2×2×2のスーパーセルを用いた計算も実行してみたが,原子数が多すぎてHOKUSAIを用いてもメモリ不足になってしまう.より高性能なスーパーコンピュータを用いた計算が望まれるが,そのようなリソースはなかなか活用が難しい.筆者が邪推するに,そのあたりのジレンマをどううまく折り合いをつけていくかが有機分子性物質におけるミュオン位置計算の現時点での落としどころになるのかもしれない.

ミュオン位置が分かればより多くの情報をμSR研究から得ることができる.今は磁性体を中心とした手法開発ではあるが,この手法が確立されてきたら次は有機超伝導体や有機スピン液体物質における計算を展開し,その特有な電子状態を実験と計算の両面から解き明かしていくことが可能になるかもしれない.μSRの研究者人口は世界でも多い方ではなく,あまつさえ計算ができる研究者は大変数が限られている.μSR測定結果のみで研究分野を切り開く時代ではなく,μSRを知らずともμSR研究者と共同して計算を実施し,そこから新しい研究展開が得られていくことを強く期待している.本稿によって,多くの計算科学研究者がミュオン科学にも興味を持ち,ともに新しい研究成果を獲得するチームを作り上げられることを祈ってやまない.

謝辞

本研究は,東北大の小池洋二教授,是恒隆准教授と川股隆行博士,上智大学の足立匡教授,東京理科大の南館隆亮博士,チョクロアミノト=パロポ大学(インドネシア)のIrwan Ramli博士,インドネシア大学(インドネシア)のBudhy Kurniawan教授,との共同研究による成果である.本研究の計算には理化学研究所情報システム本部の管理運営するHOKUSAIを活用している(課題番号G19007).また,理化学研究所のInternational Program Associate (IPA) 制度および日本学術振興会の科研費基盤B (20H04463)の支援を受けて実施されている.チーム全体を代表して深く感謝申し上げる.

参考文献
 
© 2020 日本コンピュータ化学会
feedback
Top