シミュレーションの例として，電気抵抗の計算を行ってみたいと思います。フォノンによる電気抵抗はグリューナイゼンの式によると，下記のように表されます。
デバイ温度を300Kとして，1～300Kまでの電気抵抗の温度変化を計算したのが下記のスクリプトです。
//------------------------------------------------------------------------//
//--- グリューナイゼンの式のシミュレーション
clear;
DT = 300;//デバイ温度。
C = 10000;//適当な定数。
temp = [1:1:300]';//温度の設定
DebR = (C/DT)*((temp./DT)^5).*integrate('x^5/((exp(x)-1)*(1-exp(-x)))','x',0.0,DT./temp);
//------------------------ 図の描画とデータ保存 ------------------------//
clf();//表示している図のクリア
plot2d(temp,DebR,5)//3項目の5は色を示します。5→赤
xtitle('Resistivity simulation','Temperature (K)','R_phonon');
//下記の項目を有効にすると計算結果を保存します。
//TvR=[temp, DebR];
//savematfile ResisSimulation.txt TvR -ascii -double
うまくいけば，下図の様なグラフが表示されるはずです。計算のポイントは前に書いたようにループは使用しないということです。
