FoxとSuperflipを使った粉末X線構造解析の手順

  • 2014 2/19: fox2sflip.pyを少し改良して、EDMA用の入力ファイルも同時に出力するようにしました。その関係で使い方も少し変化しました(標準出力ではなくファイルに出力、xplorファイル名自動的に設定など)。以下はそれに応じて書き直しました。
  • charge flipping法(Oszlanyi & Suto, 2004)は比較的最近に提案された結晶構造を解くための方法である。最近の論文で単結晶データで直接法とcharge flipping法が比較されていたが(van der Lee, 2013)、構造を解く能力に両者でほとんど差がないことが示されている。これはcharge flipping 法が今後広く使われていくことを示唆している。charge flipping法を実装したコードの中でもSuperflip(Palatinus & Chapuis, 2007)は上記比較論文でも使われていて、粉末回折データにも対応しており、各種OSで動作する実行ファイルが提供されている。
  • 新たにSuperflipを評価する必要があったので、粉末データを自分で測定して、解析まで行ったので、その手順について記録しておく。なお、Superflipには各回折ピークの強度を入力する必要があるが、粉末パターンからの強度抽出する機能はもってない。そこでここでは強度抽出のためにFoxプログラムを使う。したがってFoxの出力データからSuperflipとEDMAへの入力データを作る必要があるが、その手間を省くためにPythonで入力データを自動作成するプログラムを作成した。このプログラムと使った入出力データも最後のところに置いておく。この解析ではMacを使っているが、Superflip, EDMA, Fox, VestaはWindows, Linux版も提供されており、基本的に操作は同じである。

charge flipping法の原理と粉末回折への応用

  • charge flipping法では、最初ランダムな位相と測定強度から求めた構造因子の絶対値部分を使って仮の構造因子をつくり、それをフーリエ変換してセル内の電子密度マップを得る。この電子密度マップで小さい正の値(δ)よりも低い電子密度部分の符号を反転させる。これが名前の由来である。電子密度をフーリエ変換して、新しい構造因子を得る。この構造因子の位相部分だけは更新されたものを使い、構造因子の絶対値部分は実測値に直す。後はこの手順を繰り返すと位相と電子密度分布が改善されていき、最終的に構造が解かれたことになる。
  • 粉末回折パターンではピークの重なりが不可避なため、重なったピークについては正確な強度が得られない。Superflipでは分離していないピークをグループとして扱い、histogram matching法で強度の再分配を行う。histogram matching法では、上記のループの途中で、得られた電子密度分布を期待される電子密度のヒストグラムになるように変換する。ヒストグラムは類似組成の結果を使ってもいいが、そう都合よくあるわけでもないので、ほとんどの場合、セル内にあるはずの原子から推定する。そのため化学組成を与えてやる(compositionのオプション)。この新しい電子分布をフーリエ変換して、新しい構造因子を得る。グループに属するピークについては、ここで得られた構造因子の絶対値部分を再分配する。こうして得られた新たな構造因子を使って、charge flipping本来のループに戻る。

使用するデータ

  • 今回は例としてforsterite (Mg2SiO4)の粉末X線回折データを使っている。構造は良く知られていて、リチウム電池の電極として研究されているLiFePO4もこの構造である。測定は本センターのリガクSmartlab粉末X線回折計を使った。この装置は入射側にモノクロメーターがあり、Kalpha2をカットできる。検出器はD/TeX1次元検出器である。試料は標準のガラス板に詰めて、2thetaで5~140˚を測定した。わずかに不純物によるピークが含まれているが、特に考慮する必要はない。データはリガクのras形式でセーブしている。その測定データをエディタで開き、強度データ行以外を削除する。そのようにして作ったデータがforsterite.intである。Foxで読める形式になっている。

Foxによる強度データの抽出

  • ここからはFoxを使ったピーク強度の抽出である。Foxは直接空間法を使った結晶構造を解くプログラムであるが(Favre-Nicolin & Cerny, 2002)、Le Bail法が使えるのでここでは強度抽出のために使っている。もちろんFoxで構造を解くことも可能であるが、今回はSuperflipを使う。Foxを使った構造解析は既に別に書いた方を参照。
  • まずFoxを起動したら、ObjectsメニューからNew Crystal, New Powderpatternを各1つ作成する。Powder DiffractionのDataタブから"Import 2Theta-Obs-Sigma pattern"でデータを読み込む。正常に読み込まれると粉末パターンが別ウィンドウに表示される。表示されない時は入力データがおかしいか、選択したデータ形式が間違っている。正常に読み込まれたら、Radiationタブで"X-ray Tube Cu Ka1"を選ぶ。Foxでは指数付けと空間群探索の機能もあるが、今回は格子常数と空間群は分かっているのでそれらはスキップ。140˚までの全パターンをLe Bail fittingするため、Max Sin(theta)/lamdaを0.61くらいにする。PhasesタブからAdd Backgroundを選ぶ。Max Sin(theta)/lamdaは普通にFoxで構造を解く場合は0.25くらいで十分で、大きくしても計算時間が長くなるだけである。今回は全パターンでのLe Bail fittingのために特別に大きくしている。PhasesタブからAdd Crystalline Phaseを選び、作成したCrystalline objectを選ぶ。Crystalsのところで、格子常数と空間群を入力する。格子常数はa=4.75, b=10.20, c=5.98を入力、空間群にはPbnmを入力(Pbnmは非標準なセッティングであるが、よく使うのでsflip-symmetry.txtに追加してある)。Powder patternのウィンドウから右クリックで、Le Bail fittingを行う。Rwpがほとんど下がらなくなるまで、何度か繰り返して行う。ここでメインのメニューからデータをxml.gzでセーブする。今回の場合、Foxはここで終了。

Superflip用入力データの作成

  • Foxでセーブしたデータはgzipで圧縮されている。ダブルクリックすると解凍され、現れたxmlファイルはエディタで読める。その中の<HKLIobssigmaweightlist>以下にLe Bail fittingで求められた強度データが保存されているので、これをSuperflipに使う。そのままコピー&ペーストする場合は、強度の後に余分なものがついているが、dataformatでdummyと置けば、その部分は読み飛ばしてくれる。
  • histogram matching法を使う場合には近接したピークをグループ化する必要がある。そのためにどのピークがグループに属するかをこの強度リストで示してやる必要がある。具体的には各行で、最初は指数(hkl),次に強度、最後にグループを示す整数を付ける。独立したピークは0で、グループに属するピークは同じ整数をつける。最初のグループは1で、次のグループが2という風に増えていく。また、キーワードdataformatでグループのコラムに対応した位置にgroupを追加する。
  • これを手動で変換するのは面倒なので、この変換も含めて、Superflip用の入力ファイルのテンプレートを作成するプログラムをPythonで作った(fox2sflip.py)。最近手直しして、EDMA用の入力ファイルも同時に作成できるようにした。これらの入力ファイルで必要な対称操作部分も自動で作成するために、最後に置いてあるsflip-symmetry.txtを同じディレクトリーに置いておく必要がある(スクリプト内で指定しているので、置き場所はそこで指定してもよい)。なおFoxで使う空間群の名称に対応するためにsflip-symmetry.txtの内容は更新されており、また誤りも直した(底面心格子で軸を取り替えた場合のcentersが正しくなかった)。Foxで実測パターンと計算された強度を見比べたところ、強度の補正は行われているように見受けられるので、多重度とLP因子の補正は行っていない。このプログラムはターミナルで使う。2つのパラメーターを取り、最初はxmlファイル名、2番目は重なりを評価するためのfwhm(角度)である。ピークがfwhm/2以下しか離れていない時はグループとして扱う。fwhmは省略すると0.1が使われる。経験的には実際のピークのfwhmよりもさらに小さくした方がいいようである。以前出力は標準出力に出していたが、今は使ったxmlファイルのベース名を使った名前のファイルを自動的に作るようにした。実行はターミナル上で以下のようにする。
    > ./fox2sflip.py forsterite.xml 0.05
    これでできるforsterite.inflipファイルはまだ完全ではないが、とりあえずはそのままでSuperflipで読み込むことができるようになっている。直す部分としては組成とZが分かっているなら、キーワードcompositionを正しく与えてやる必要がある。これにはセル内にある元素数を指定する必要がある。forsteriteの場合は、Mg2SiO4でZ=4なので、"composition Mg8 Si4 O16"となる。また先頭の#を削除する。さらにhistogram compositionの先頭の#を削除する。もし組成自体が不確定な場合はcompositionとhistogram行は#を先頭につけたままにする。このテンプレートは標準的なセッティングなので、必要に応じてパラメーターなどを変える必要がある。
  • Foxの空間群探索では非標準な空間群を含めて候補として示される。International Table Aに載っている軸と座標原点の取り方はsflip-symmetry.txtでもカバーしているが、それ以外は標準に戻す必要がある。Foxで空間群の指定と格子常数を直して、Le Bail fittingを再度行えば良い。また、立方晶でFoxの格子常数がa = b= cに正確になっていないことがある。これは正方晶、六方晶のa,bについても起こる。このような場合もエラーとなるので編集する必要がある。もしテンプレートにsymmetry - endcentersが書き込まれてない場合は非標準な空間群を使ってないか、格子常数がブラベ格子と一致しているかをチェックする。Pythonはインタープリターなので、もし不具合があればコード自体の手直しもユーザーが簡単にできる。特にSuperflipのinflipファイルのデフォルトのパラメータは利用者の要望に合わせて変えた方がよい。

Superflipの実行

  • 入力ファイルが出来たら、Superflipを実行するが、ターミナルで使う必要がある。Superflipのバイナリーがあるディレクトリにいて、入力ファイルも同じところにあるなら、ターミナルから
    > ./superflip forsterite.inflip
    で計算を実行させることができる。ここでよくあるトラブルはhistogram compositionが指定されているが、compositonが与えられていない、格子常数がおかしいなどである。計算が終了したらSuccessがあるかチェック。Successがゼロでも実際には正解の電子密度が得られていることもある。計算で出来たxplorファイルをダブルクリックして、Vestaで電子密度を表示させる。下にVestaで電子密度を表示したものを示してある。forsteriteの場合、Z=4なので、セル中には原子が28個あり、28個の電子密度の高い領域がきれいに確認でき、構造は正しく得られている(文献のデータと比較するとb軸方向で座標が1/2シフトしているが)。また、電子密度の高い部分が四面体や八面体を作っているのが見て取れる。もし28個見えない場合は、VestaのObjects->Properties->Isosurfaces...の設定からIsosurface levelを下げてやると見えるかもしれない。重い元素がある場合、酸素などの軽元素は見えないことも多い。なお、このデータの場合はhistogram matching法を使わなくても、ほぼ同じ結果が得られた。
    forsterite-map.png

電子密度の解析

  • EDMAプログラムはSuperflipで得られた電子密度を解析して、原子を同定して、等価な原子位置を整理してcifファイルに書き出してくれる。EDMAはSuperflipと同じサイトからダウンロードできる。EDMAの入力ファイルforsterite.edmaはfox2sflip.py実行時に既に作成されている。こちらは化学組成を定義しないと実行できないので、composition行で組成を入力する。定義はSuperflipの入力ファイルと同じである。Zが最初分かっていない場合もあろうが、先にVestaで表示した電子密度からZを推定できるかもしれない、EDMAの実行もターミナルから
    > EDMA forsterite.edma
    で実行する。cifファイルが作られるので、それをVestaで読む。EDMAの結果では28個が6つの非等価なサイトに正しく整理されていた。Mg, Siは原子番号近いので区別できていない場合もある。原子種の間違いなどはVesta上で編集すればよい。
  • Vestaでも電子密度ピークの位置などを求めることができる。やってみるとリストの最初の28個までは電子密度ピークが有意に高い。そのうち後半の16個は密度が前半のものより低く、これらが酸素原子だと分かる。位置から等価な席を割り出し、6つの非等価な席を得ることは簡単にできる。ここは結合距離などを見て判断する。Siは酸素4つに囲まれている方の席である。今回のデータの場合、Macbook Airでも最速数秒で構造が解かれており、charge flipping法の威力がよく分かる。Vestaで得られる位置は、特殊位置に近いと特殊位置と解釈してしまう傾向があり、逆にEDMAは特殊位置に非常に近くても特殊位置とはしない傾向がある。

その他のデータ

  • Foxで以前解析した自前データ等を10個くらい試してみた。当時もSuperflipを試したが、それらしい電子密度がうまく得られなかった。今回は空間群が既に確定しているところが違うが、多くで正しい電子密度が得られた。AlPO4-moganiteの例を下に示す。全ての原子が見えている。回折データが既にあれば、10分程度の作業。重い元素を含む酸化物やデータが悪いものでは酸素が見えないこともあった。また、低対称性でかなり複雑な場合には全く意味ある電子密度が得られなかった。
    AlPO4-mog-superflip.png

結果がよくない場合の対応

  • 試してみる価値のある調整は、2thetaの範囲(Fox上で), fwhmを変えてみる、delta, weakratio, Bisoを変えてみる、試行回数(repeatmode)を増やすなど。
  • もし重い原子だけ見えているなら、それらの原子座標をVesta又はEDMAで求めて、それをFoxを使った実空間法で残りの元素の位置を求めることができる。FoxのCrystal objectを定義して、Superflipで見えた原子位置の非等価な座標のみを入力し、座標は固定する。一度Le Bail fittingをしてから、差フーリエマップを描かせて軽い元素が見えるかチェックする。PやSiの場合は配位多面体を定義してよい。さらに見えてない原子をCrystalタブにatomとして追加してglobal optimizationを行い、そららの原子位置を求める。
  • 全く見込みがない場合は、Foxで実空間探索法の計算を行う、または直接法(EXPO2013など)を試してみる。

空間群が不明な場合

  • Foxの空間群探索では正解がRwpトップ3位以内に来ていることが多いが、全くダメな場合もある。私の経験では小さいピークが空間群を左右する鍵を握っているのに、それらの小さいピークを無視して、より少ないピーク数で強度の強いピークをフィットした方がRwpが低くなるケースがあった。この場合正解からはほど遠い空間群が探索リスト上位に出て来ることになる。そのような場合は、粉末パターンを拡大してよく観察して、小さいピークがちゃんとその格子で説明されるかどうかを見る必要がある。指数を見て、消滅則から空間群を制約する。International Table Aには消滅則が載っているが、桜井敏雄の「X線結晶解析の手引き」裳華房の付表1の「空間群決定の表」が便利である。
  • 空間群が不明な場合、ありそうな空間群をテストする。この場合はまずFoxに戻って、そこで使用する空間群を指定してLe Bail fittingを行わなくてはならない。その空間群を使って、Superflipを使う。
  • charge flipping法は空間群を知らずに(P1で)解析できると論文等に書かれている。これは正しいが、粉末データの場合は強度抽出時に空間群が必ず必要となる。
  • Superflipで得られた電子密度から空間群をチェックしたい場合は.inflip中のderivesymmetryを有効にする(#を削除する)。強度抽出時に空間群を使っていて、.inflipでも空間群を使っているのでそれにバイアスされないか少し疑問。

構造の精密化

  • 上記のように正しそうな初期構造を得たら、EDMAを実行してcifファイルを出力させる。このファイルをVestaで読んで、必要なら修正を行う。原子種の取り違えが多い。VestaでRIETAN-FPのinsファイルを作成する。insファイルを修正して、RIETAN-FPでRietveld解析を行う。

作業をshell scriptで制御

  • 一連の作業を行う時にはshell scriptでプログラムを逐次動作させる方が便利。以下は簡単なscriptであるが、Foxの起動、xml.gzファイルの解凍, 入力ファイルの作成、Jeditの起動, Superflipの実行、Vestaで電子密度を見る、EDMAの実行, Vestaでcifを見る、を行う。
#!/bin/bash

base1="./"
# directory for programs, change for your system
base2="./"
# directory for data files, change for your system

flip="superflip"
edma="EDMA"
f2sf="fox2sflip.py"
file=$1

open -a "fox"

echo -n "Input file name> "
read file
echo -n "Input fwhm (default=0.1) > "
read fwhm

gunzip $base2${file}.xml.gz
 
$base1$f2sf $base2${file}.xml $fwhm
echo -n "To edit inflip file, hit return key! > "
read x
open -a "Jedit X" $base2$file.inflip
echo -n "To start Superflip, hit return key! > "
read x
$base1$flip $base2${file}.inflip
echo -n "To view xplor file by VESTA, hit return key! > "
read x
open -a "VESTA" $base2${file}.xplor
echo -n "To edit edma file, hit return key! > "
read x
open -a "Jedit X" $base2$file.edma
echo -n "To start EDMA, hit return key! > "
read x
$base1$edma $base2${file}.edma
echo -n "To view cif file by VESTA, hit return key! > "
read x
open -a "VESTA" $base2${file}.cif
 

文献:

  • Favre-Nicolin, V. and Cerny, R. (2002) J. Appl. Cryst. 35, 734-743
  • Momma, K.; Izumi, F. (2011) J. Appl. Cryst. 44, 1272−1276.
  • Oszlanyi, G. & Suto, A. (2004) Acta Cryst. A60, 134–141
  • Palatinus L., & Chapuis G. (2007) J. Appl. Cryst. 40, 786-790.
  • van der Lee, A. (2013) J. Appl. Cryst., 46, 1306-1315

プログラムの入手先:

ここで使った入出力ファイルとコード等

もう少し難しい構造解析の例

SDPD Round Robin 2のsample2(Sr5V3(O/OH/F)22))

  • SDPD Round Robin 2のsample2(Sr5V3(O/OH/F)22))をSuperflipで解いてみた。指数付け、空間群推定は全く問題なし。SuperflipでSr,Vはもちろん、酸素もほぼ見えたので、酸素位置をFoxで見つける必要はなかった。EXPO2013も試したが、こちらでも問題なく初期構造が得られた。酸素位置の細かい違いを無視すれば同じ構造が得られた。原子は全て一般位置にあるので自由度は90になる。もしこれをFox単独で解析を行う場合には、結構な困難が予想される。SrとV位置は比較的簡単に決まるだろう。しかし酸素位置はVOx多面体として解析するのがいいだろうが、Vの配位数は4,5,6いずれの可能性も考慮しないといけないだろう。また組成からはVにはつながってない酸素も存在するはずなので、そちらは単体の原子として追加する必要がある。Sr,V位置を固定すれば、自由度は21~39(Vの配位数に依存)に下げられる。

SDPD Round Robin 3のsample2 (LaWO)

  • SDPD Round Robin 3のsample2 (LaWO)をFox+Superflipで解いてみた。情報としてhexagonalであり、大まかな格子常数が与えられている。X線と中性子のデータが与えられているが、まずX線データを使う。Foxの指数付けではAdvancedで、この例ではc軸がかなり長いので、Length maxとVolume maxをデフォルト設定より増やしておく必要がある。hexagonalを指定する。ただそれだけでは解は得られず。不純物ピークを1にすると(Nb spurious lines)、scoreの高い解が1つ得られた。これは与えられている値と大体合う。これを使って空間群の推定。ユニークには決まらず5つが同じRwp, GoFで並ぶ。Superflipを各空間群で実行すると、それぞれそれなりの電子密度が得られる。なおinflip, edmaファイルで格子常数がa=bになってないので、エラーが出る。実行前に直しておく必要がある。Superflipでderivesymmetryオプションを付けて、得られた電子密度から空間群をチェックした。全てでP-62cとなり、これ自体は5つの候補の1つであったので、以下P-62cで解析。Vestaで電子密度を見て、ピーク位置を求めてみると、58個の明瞭なピークがあり、これらがWとLaであろう。そのうち、最初の18個までで一区切りなので、これらがWで、残り40個がLaであろうと予想される。La/Wの比率が与えられている組成(2つの可能性)とずれているのでLaとWを取り違えている部分があるのだろう。とりあえずその数をcompositionに入れて、EDMAを実行するか、Vestaでピーク位置を求めて、その位置から非等価な席を決めて行く。EDMAではWが28個同定されていたが、これは多すぎるので、一般位置にあるWが多分Laであると考えられるので訂正する。また、よく見ると4f特殊位置にあるLaは隣の同じく4f位置のLaとの距離が2.2 Aくらいと非常に短く、これは明らかに何かおかしいが、電子密度は明瞭に見えている。1/2の部分占有が想定される。また、その近くのWとの距離が2.5 Aくらいであり、これもLa-Wの距離としては短すぎる。W-Wでも短いが、八面体の面共有ならちょうどいい距離。従って、このLaは実際にはWである。しかし、こういった結晶化学的考察をしないといけないのでは、経験のない人には解析は難しい。結局、Laは非等価な席が4席、Wは6席あって、Wの1つは部分占有(1/2)と考えざるを得ない。組成はLa36W20O114となり、与えられていた組成の最初のものに近い。Wは普通6配位であろうから、FoxではW. Laを得られた位置に固定しておいて、WO6の回転のみを考えればよいはず。多分Wに結合してない酸素はないはず。X線では酸素位置はあまり決まらなそうなので、与えられている中性子回折データでFoxの計算を行うことになろう。そこまでまだやってない。

トップ   編集 凍結解除 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS
Last-modified: 2014-02-19 (水) 22:21:38 (1370d)