Kohlausch-Williams-Watts(KWW)式は、高分子やガラスの誘電緩和スペクトルの解析に使われている。デバイ型の緩和では、双極子の相関が時間領域では指数関数1つで表されるのに対し、KWW式では、引き延ばされた指数型(stretched exponential)で表される。すなわち、φ(t)=exp{-(t/τ)β}の形になる。これは、緩和過程が単一の緩和時間ではあらわせず、緩和時間が分布しているという効果を現象論的に簡単な形でとりいれたものである。
デバイの式は、周波数領域で簡単な解析的な形に書けるので、誘電緩和スペクトルのフィッティングに広く用いられてきたが、KWW式の周波数領域での形は積分や級数の形になる。このためかどうか知らないが、KWW式とHN式のパラメータの関係について調べた論文があったり、KWWで表せるといっておきながら実際のフィッティングはHN式だったりということがある。
私はこれまでにIgorProをデータ解析に使ってきたが、今回、KWW式の周波数領域の関数を使って誘電緩和スペクトルのフィッティングをできるようにしたので、ソースも含めて公開する。
ダウンロードしたあと、KWW XFUNCフォルダの中にある、KWWXFUNC PPCをIgorProフォルダの中のIgor Extensionsフォルダに入れて、IgorProを起動する。使い方は、IgorProのMiscメニューから、Function Helpを選び、External functionのリストからnKWWを選ぶとヘルプが表示される。また、これを使ってIgorProで計算を行ったサンプルをつけてあるのでそちらも参考にしてほしい。
具体的になにをやっているかは、Documentationフォルダの中のKWW.pdfを見ること。必要な参考文献もこれに書いてある。
KWW XFUNCを使うには、Power Macintosh、OSバージョン7以降、Igor Proバージョン3以降が必要。
KWW XFUNCの改造およびその他の付属ファイルを改造するには、Igor Pro XOP ToolkitとCodeWarrior Proが必要。
一応、68KMac用のプロジェクトファイルも作ってあるが、環境がないため動作確認はできていない。
この分野の研究者ならたいてい知っている誘電緩和のテキスト
C.J.F Bottcher and P. Bordewijk, "THEORY OF ELECTRIC POLARIZATION" 2nd. ed. volume II, Elsevier (1978), ISBN 0-444-41579-3
は誘電緩和の緩和関数型や現象論に詳しいが、この本のKWWの部分の記述は間違っている。パラメータの範囲と使うべき計算式が食い違っていて、この本の通りにプログラムを作ってもうまく動かない。私はこれで1週間程悩んだ。WilliamsとWattsの元論文にあたること。こちらが正しい。
IgorProのFuncFitを使ってフィッティングを行うときに、βの初期値を1にすると、明らかに緩和時間が分布しているようなデータであるにも関わらず、イテレーションでβが変化せず1のままである。実装したときに、βが負の値になったり1を越えたりしたときは、数値計算エラーを避けるために、ガンマ関数などの計算は行わず十分大きな数値を返すようにしたのが原因かもしれない。初期値を1以外の、実際のデータに近いスペクトルになるような値から始めれば問題は起きない。
フォルダは朱で示した。
Documentation KWW.pdf ドキュメント。必ず読むこと Igor Procedure KWWproc フィッティングをするためのプロシージャ。HN式のときとほぼ同じ。 Test experiment フィッティング例のIgor experimentファイル Integration 実装のために行った数値積分のためのルーチンと結果。 MIDEXP.C MIDINF.C MIDPNT.C MIDSQL.C MIDSQU.C NRUTIL.C nrutil.h POLINT.C QROMB.C QROMO.C QSIMP.C romb1.results ロンバーグ積分による結果 rimb1_50_06 romb1_50_01 romb1_50_02 romb1_50_03 romb1_50_04 romb1_50_05 romb1_50_07 romb1_50_08 romb1_50_09 romb1_50_10 romb1_5_01 romb1_5_02 romb1_5_03 romb1_5_04 romb1_5_05 romb1_5_06 romb1_5_07 romb1_5_08 romb1_5_09 romb1_5_10 romb2.results ロンバーグ積分による結果(その2) romb2_50_01 romb2_50_02 romb2_50_03 romb2_50_04 romb2_50_05 romb2_50_06 romb2_50_07 romb2_50_08 romb2_50_09 romb2_50_10 romb2_5_01 romb2_5_02 romb2_5_03 romb2_5_04 romb2_5_05 romb2_5_06 romb2_5_07 romb2_5_08 romb2_5_09 romb2_5_10 romb2_1e-6 ロンバーグ積分による結果(その2)、収束条件変更。 romb2e6_50_02 romb2e6_50_03 romb2e6_50_04 romb2e6_50_05 romb2e6_50_06 romb2e6_50_07 romb2e6_50_08 romb2e6_50_09 romb2e6_50_10 romberg1 romberg1.c romberg2 romberg2.c romberg2.prj romberg2.prj Data CW Settings.stm romberg2.tdm simp simp.c simp.prj simp.prj Data CW Settings.stm simp.tdm simp.results シンプソン法による結果 simp50_02 simp50_03 simp50_04 simp50_05 simp50_06 simp50_07 simp50_08 simp50_09 simp50_10 simp5_01 simp5_02 simp5_03 simp5_04 simp5_05 simp5_06 simp5_07 simp5_08 simp5_09 simp5_10 TRAPZD.C xmidont xmidont.prj xmidont.prj Data CW Settings.stm xmidont.tdm XMIDPNT.C xpolint XPOLINT.C xpolint.prj xpolint.prj Data CW Settings.stm xpolint.tdm xqromb XQROMB.C xqromb.prj xqromb.prj Data CW Settings.stm xqromb.tdm xqromo XQROMO.C xqromo.prj xqromo.prj Data CW Settings.stm xqromo.tdm KWW XFUNC KWW XFUNCのソースコードとプロジェクト calcKWW.c GAMMLN.C KWWXFUNC 68K.prj KWWXFUNC Help KWWXFUNC Help Resources KWWXFUNC PPC KWWXFUNC PPC.prj KWWXFUNC PPC.prj Data CW Settings.stm KWWXFUNC PPC.tdm KWWXFUNC Resources KWWXFUNC.c KWWXFUNC.h KWWXFUNC.r POLINT.C QROMB.C TRAPZD.C KWW XFUNC test KWW XFUNCのテストのIgor Experimentファイル test0 test1 test2 test3 test4 tokai.data 東海大学八木原研・新屋敷さんよりいただいた各βの値での計算値。 計算のチェックに使用した。 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 simulation.exp 積分方法の検討の結果をグラフ表示したもの。IgorPro experiment file。 history1 KWW.org1 KWW.org2 kww02.exp kww03.exp kww04.exp kww05.exp kww06.exp kww07.exp kww08.exp kww09.exp kww10.exp kwwseries.results 級数を使ったときの計算結果。IgorPro wave file。 kww01.awav kww02.awav kww03.awav kww04.awav kww05.awav kww06.awav kww07.awav kww08.awav kww09.awav kww10.awav romberg2conv