![[PukiWiki]](image/c60.png "[PukiWiki]")
[2018/08 started by Maruoka and Maeda]
Linux環境でPythonを用いてMEEPを使えるようすることを目標とします。
WindowsではMEEPがLinuxかMacでしか動作しないのでLinuxのWSL(Windows Subsystem for Linux)のインストールをしてください。
やり方はこちら。
WindowsでLinuxが使えるようになったらLinuxのと同じ手順で環境が構築できます。
MacではAnacondaをコマンドラインインストーラーを使ってインストールすればLinuxの時と同様にできると思います。(未確認)
まずはAnacondaを用いてPythonの環境を構築します。(参考)
AnacondaのサイトにてPython3.x系のインストーラーをダウンロードしてください。
インストーラーは.shのファイルなので以下のように実行してください。
$ bash [ダウンロードフォルダ]/Anaconda3-5.2.0-Linux-x86_64.sh
次にPyMeepのインストールをします。(参考)
以下のコマンドを実行してください。
$ conda create -n mp -c chogan -c conda-forge pymeep
これで必要なセットアップは終わりました。
以下のコマンドでアクティベートしてください。
$ source activate mp
"-bash: activate: No such file or directory" とエラーが出たときは
export PATH=~/anaconda3/bin:$PATH
を実行してください。 flexやtubeでは
#!/bin/sh _CONDA_ROOT="/home/students/maru/anaconda3" \. "$_CONDA_ROOT/etc/profile.d/conda.sh" || return $? _conda_activate "$@"
というエラーが出ると思います。これはログイン先のシェルが、tcshなのに、shのスクリプトをsourceしてしまっているのが原因です。 tcshからはtcshスクリプトしかsource出来ません。
bash
と打って、ログインシェルをbashに変更してから実行してみてください。 実行はhogehoge.pyの部分を任意のファイルに置き換えてください。
(mp) $ python hogehoge.py
サンプルソースはGithub等で入手できるのでこのソースコードで試してみてください。
無効かするには以下のコマンドを実行してください。
(mp) $ source deactivate
並列処理をサポートしたpymeepをインストールするには上記の作業をすべて行ったうえで次のコマンドを実行してください。
$ conda create -n pmp -c chogan -c conda-forge pymeep-parallel
$ source activate pmp
にて並列処理をサポートしたpymeepをアクティベートしたのち、
$ mpirun -np 4 python <script_name>.py
を実行してください。-np 4とすると4プロセッサで並列処理されます。数字は適当に変更してください。
meep <file-name>.ctl &
とすると、バックグランドでプログラムが実行される。可変パラメータを外部から渡すことで、多数のプログラムを同時に動かせる。例えば、次のようなシェルスクリプトを著者は使用している。入射角を0から85度まで5度ずつ変えたプログラムを同時に走らせている。
for((i=0;i<90;i+=5))
do
meep incident-angle=${i} Main.ctl &
done
MEEPで用いられるMaxwell方程式は
など通常のMaxwell方程式と異なります。
MEEPで用いられるMaxwell方程式は次のような式変形の最後の形です。
（ドキュメント参考）
通常の以下のMaxwell方程式について、
以下のように変形します。
以下のようにx',y',z',t',ω',E',H',σ'を定めます。
このとき、Maxwell方程式は次のようになります。
この形のMaxwell方程式を使っています。
このとき、真空の誘電率、透磁率はどちらも1です。
また、単位長さはaです。
単位時間に真空中の光は単位長さaだけ進みます。
周波数fの真空中の光の波長は1/fです。
apt-get install python-h5py丸岡のbash on windowsではpip install h5pyでは、インストールはできても実行できなかった（ImportError: No module named h5py）。
次のコードで左回りの円偏光になります。 amplitudeによって、ソースの位相を変化することができます。
(set! sources
(list
(make source
(src (make continuous-src (frequency (/ 1 7))))
(component Ex)
(center 0 0 8)
(size 50 50 0)
)
(make source
(src (make continuous-src (frequency (/ 1 7))))
(component Ey)
(center 0 0 8)
(size 50 50 0)
(amplitude (exp (* 0+1i 1.570796327)))
)
)
)
次のコードで右回りの円偏光になります。 amplitudeによって、ソースの位相を変化することができます。
(set! sources
(list
(make source
(src (make continuous-src (frequency (/ 1 7))))
(component Ey)
(center 0 0 8)
(size 50 50 0)
)
(make source
(src (make continuous-src (frequency (/ 1 7))))
(component Ex)
(center 0 0 8)
(size 50 50 0)
(amplitude (exp (* 0+1i 1.570796327)))
)
)
)
x方向に直線偏光した平面波とy方向に直線偏光した平面波を 中心(0,0,8)、大きさ50×50×0の平面から、xy平面について45度方向((x,y)=(1,1)方向)に直線偏光した平面波を出力するコード。
(set! sources
(list
(make source
(src (make continuous-src (frequency (/ 1 7))))
(component Ex)
(center 0 0 8)
(size 50 50 0)
(amplitude (exp (* 0+1i 1.570796327)))
)
(make source
(src (make continuous-src (frequency (/ 1 7))))
(component Ey)
(center 0 0 8)
(size 50 50 0)
(amplitude (exp (* 0+1i 1.570796327)))
)
)
)
斜め入射のシミュレーションの方法は3つあります。
(define (pw-amp k x0) (lambda (x)
(exp (* 0+1i (vector3-dot k (vector3+ x x0))))))
(define-param theta-deg 48)
(define theta-rad (deg->rad theta-deg))
(define epsilon1 2.25)
(define epsilon2 1.25)
(define-param fcen 0.8) ; pulse center frequency
(define-param df 0.02) ; turn-on bandwidth
(define-param kdir (vector3 (sin theta-rad) 0 (cos theta-rad))) ; direction of k (length is irrelevant)
(define k (vector3-scale (* 2 pi fcen (sqrt epsilon1))
(unit-vector3 kdir))) ; k with correct length
(set! pml-layers (list (make pml (thickness 2)(direction Z))))
(set! k-point (vector3* (* fcen (sqrt epsilon1) ) (vector3 (sin theta-rad) 0 (cos theta-rad))))
(set! geometry
(list
(make block
(center 0 0 -2.5)
(size 10 infinity 5)
(material (make dielectric(epsilon epsilon1)))
)
(make block
(center 0 0 2.5)
(size 10 infinity 5)
(material (make dielectric(epsilon epsilon2)))
)
)
)
(set! sources
(list
(make source
(src (make continuous-src (frequency fcen) (fwidth df)))
(component Ez) (center 0 0 -3) (size 10 0 0)
(amp-func (pw-amp k (vector3 0 0 -3))))))
(define-param T 40) ; run time
(run-until T
(to-appended "ez" (at-every 0.1 output-efield-z))
(to-appended "ex" (at-every 0.1 output-efield-x))
)
Ezのプロット
Exのプロット
Ezのプロット
光源が一つしかない場合のEzのプロット
MEEPのbuilt-in関数には電荷分布を求める機能がないので、自作する必要がある。
ガウスの法則を差分化すればよい。
前述の全反射直前のTM波(ブロッホ周期境界条件)のシミュレーションについて例を示す。
電荷分布(+電場ベクトル)のプロット
コード
PMLをある方向にだけ設置するにはオプションdirectionを設定すればよい。例えば、Y方向だけに設置するには以下のようにすればよい。
pml_layers = [mp.PML(thickness=1.0, direction=mp.Y)]
PMLが設置されていない方向は固定端の全反射が起きる。 真ん中にwaveguideを置いて、Y方向にだけPMLを設置したときの電場は次のようになる。 X方向にPMLがないので全反射が起きて定在波が生じている。
コードは次の通り。
import meep as mp
cell = mp.Vector3(16,8,0)
geometry = [mp.Block(mp.Vector3(1e20,1,1e20),
center=mp.Vector3(0,0),
material=mp.Medium(epsilon=9))]
sources = [mp.Source(mp.ContinuousSource(frequency=0.1),
component=mp.Ez,
center=mp.Vector3(-7,0))]
pml_layers = [mp.PML(thickness=1.0, direction=mp.Y)]
resolution = 10
sim = mp.Simulation(cell_size=cell,
boundary_layers=pml_layers,
geometry=geometry,
sources=sources,
resolution=resolution)
dT = 0.2
sim.run(
mp.at_beginning(mp.output_epsilon),
mp.to_appended("ez",
mp.at_every(dT, mp.output_efield_z)),
until=100
)
物性値のライブラリが公開されています。
materials-library.scmをダウンロードしてください。
ライブラリを使用するには、プログラムの最初に
(include "/path/to/materials-library.scm")
を追加してライブラリをロードしてください。
"/path/to/materials-library.scm"は自分がmaterials-library.scmを設置した箇所に書き換えてください。
丸岡がダウンロードしたmaterials.scmを使用する場合は次のパスを使用してください。
(include "/home/students/maru/meep/scheme/materials.scm")
ただし最新のmaterials.scmではない可能性があります。
このライブラリはtube61, 62, 63でしか動きません（flex, tube60では動きません）。
詳しい使い方は https://meep.readthedocs.io/en/latest/Materials/ を参照してください。
初期設定ではMEEPの単位長さが1μmであるものとしてmaterials.scmが設計されています。
異なる単位長さを使用したい場合、特別に指定する必要があります。
例えば、0.1µmを単位長さにしたいときは
(set! um-scale (* um-scale 0.1))
としてください。
物性値のライブラリmaterials.pyをロードするのに特別な操作をする必要はありません。
初期設定ではMEEPの単位長さが1μmであるものとしてmaterials.pyが設計されています。
異なる単位長さを使用したい場合、materials.pyを直接書き換える必要があります。
um_scaleが宣言されている行以降でum_scaleを目的の単位長さに書き換えてください。
例えば、0.1µmを単位長さにしたいときは
um_scale=um_scale*0.1
としてください。
materials.pyの場所は次のようにして調べてください。
(pmp) maru@tube60:~$ python Python 3.6.6 | packaged by conda-forge | (default, Jul 26 2018, 09:53:17) [GCC 4.8.2 20140120 (Red Hat 4.8.2-15)] on linux Type "help", "copyright", "credits" or "license" for more information. >>> import meep Using MPI version 3.1, 1 processes >>> meep.__file__ '/home/students/maru/anaconda3/envs/pmp/lib/python3.6/site-packages/meep/__init__.py'
この場合、materials.pyの場所は/home/students/maru/anaconda3/envs/pmp/lib/python3.6/site-packages/meep/materials.pyです。
このmaterials.pyの書き換えは(pmp)の環境でだけ有効であることに注意してください。
((mp)でmaterials.pyの単位長さを変える場合は(mp)のmaterials.pyを書き換えてください）
pythonによるhdf5ファイルの取り扱いを説明する。
import h5py #以降のプログラムのpukiwikiでの掲載ではこの行は省略する
f = h5py.File("file-name.h5", "r") #rとすると読み込み専用で開く
print(f)
とすると、<HDF5 file "file-name" (mode r)>と表示される。
print(list(f).keys())
とすると、 h5fileのkeyが表示される。例えば、['ex']と表示されたする。
Ex = f['ex'][:]
とすると、キー'ex'に紐づいた配列をExに取り出せる。注意として、f['ex'][i]のようにして、何回もアクセスしないこと。このようなアクセスは非常に時間がかかる。
print(Ex) #print(h5file['ex'].value)も同様の動作をする。
とすると、Exの配列の中身が表示される。要素が多すぎるときは・・・と表示されて省略される。 Exが二次元配列でEx[i][j]と参照するとしたとき、iを第一インデックス、jを第二インデックスと呼ぶことにする。 第一インデックスの長さ(iの動く長さは)
print(len(Ex))
で分かる。第二インデックスの長さは
print(len(Ex[0]))
で分かる。
以下のスクリプトを.bash_profileに貼り付けて使用すると、h5ファイルを動画に自動で変換できて便利です。
function h5mp4(){
mkdir $1
cp $1.h5 $1
cd $1
h5topng $1.h5 -Z -c dkbluered -R -t 0:$2
convert $1*png $1.gif
ffmpeg -i $1.gif -qscale:v 1 $1.mp4
}
例えば時間の区切り回数が1000回あるようなhogehoge.h5を動画にしたいときは
h5mp4 hogehoge 999
とするとhogehoge.gifとhogehoge.mp4とhogehoge-<数字>.pngが自動で生成されます。 時間の区切り回数が知りたいときは
h5ls hogehoge.h5
としてください。
ez Dataset{100, 200, 1000/Inf}
と表示される数字のうち、1000/Infが区切り回数です。
import h5py
import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import matplotlib.animation as animation
from PIL import Image
import moviepy.editor as edit
from mpl_toolkits.axes_grid1 import make_axes_locatable
plt.rcParams["font.size"] = 18
#----------input----------------
surfaceExFilename ="After-ex_+0.5.h5"
surfaceEyFilename = "After-ey_+0.5.h5"
surfacePxFilename = "Pre-ex_+0.5.h5"
surfacePyFilename = "Pre-ey_+0.5.h5"
savename = "current-line-withPreE-takingMiddleCurrent.png"
startTime = 0
endTime = 1000
waveLength = 6
timeResolution = 10
xid=16 #array[xid][yid][:]の電場をプロット
yid=46
#--------------------------------
def plotf(arr,label_):
time=np.array([])
amplitude=np.array([])
for i in range(startTime, endTime):
x=i * 1./(waveLength * timeResolution)
y=arr[xid][yid][i]
time=np.append(time,[x])
amplitude=np.append(amplitude,[y])
plt.plot(time,amplitude,label=label_)
plt.legend(bbox_to_anchor=(1.05, 1), loc="upper left")
h5file_A_x = h5py.File(surfaceExFilename, "r")
h5file_A_y = h5py.File(surfaceEyFilename, "r")
h5file_P_x = h5py.File(surfacePxFilename, "r")
h5file_P_y = h5py.File(surfacePyFilename, "r")
Ax = h5file_A_x['ex'][:]
Ay = h5file_A_y['ey'][:]
Px = h5file_P_x['ex'][:]
Py = h5file_P_y['ey'][:]
Ellは2次元配列でxyに応じて楕円率が入っていると想定
ヒストグラムの正規化はググると「normed=True」が出てくるがうまくいかなかった
stackoverflowが役に立った
日本語解説ブログ
ecのキーワードでヒストグラムの縦線の色を指定できる(参考)
Ell_ = Ell.reshape(-1,)#plt.histは1次元配列を対象とするのでEllを1次元化
weights_ = np.ones(len(Ell_))/float(len(Ell_))#ヒストグラムを正規化するための処理
fig = plt.figure()
ax = fig.add_subplot(111)
#binsはヒストグラムの縦棒の数
ax.hist(Ell_, weights=weights_, bins=21, range=(-1,0,1.0),ec='black')
ax.set_xlim(-1.0,1.0)
ax.set_ylim(0.0,0.5)
fig.savefig("histgram.png")
array.reshape(-1)の挙動
>>> a = np.linspace(1,8,num=8)
>>>
>>> a
array([1., 2., 3., 4., 5., 6., 7., 8.])
>>> b = a.reshape((2,4))
>>> b
array([[1., 2., 3., 4.],
[5., 6., 7., 8.]])
>>> c = a.reshape((2,2,2))
>>> c
array([[[1., 2.],
[3., 4.]],
[[5., 6.],
[7., 8.]]])
>>> d = c.reshape(-1,)
>>> d
array([1., 2., 3., 4., 5., 6., 7., 8.])
代わりにweightsオプションを指定すればこの問題は解決できる. これはdataの1つの値の重みを指定するパラメータで, デフォルトでは実質的に weights = np.ones(len(data)) を指定していることになる. いま総和を1にしたいのだから, ひとつのデータの重みは1/len(data)である.
ndarrayの場合は配列 / scalarが可能。
>>> x = [1,2,3] >>> x [1, 2, 3] >>> x /2 Traceback (most recent call last): File "<stdin>", line 1, in <module> TypeError: unsupported operand type(s) for /: 'list' and 'int' >>> y = np.array([1,2,3]) >>> y / 2 array([0.5, 1. , 1.5])
楕円率のヒストグラムの作成の際の hist()のキーワードに
range=(-1,1)
を入れなないと空気を読んで，最大値と最小値の間をbinsの数だけ分割するようです
ax.set_xlim(-1., 0, 1., 0)
だと，-1から0までにしかプロットされないので
ax.set_xlim(-1.,1.)
に修正
setenv GUILE_WARN_DEPRECATED no
「unbound variable: 」と、エラーが出てるにもかかわらず、unbound variableの名前が表示されないときは、全角空白文字がどこかに含まれている。全角空白を削除すればよい。
meep hoehoge.ctl var=1
とすると、var=1は無視される。
meep var=1 hogehoge.ctl
とする。
e1,e2,e3[vector3]により、方向を定められる。 例えば、xy平面内でx軸からΘだけ傾けて置きたいならば、
(make block (size a b c) (e1 (vector3 (cos theta) (sin theta) 0)) (e2 (vector3 (- (sin theta)) (cos theta) 0) (e3 (vector3 0 0 1)) )
X, Y, Z方向をそれぞれ厚さa, b, cにしたいなら、
(set! pml-layers (list (make pml (thickness a)(direction X)))) (set! pml-layers (list (make pml (thickness b)(direction Y)))) (set! pml-layers (list (make pml (thickness c)(direction Z))))
とする。
負の透過率、反射率、吸収率は電磁波が発散しきっていないのにシミュレーションを中断したことが原因のことが丸岡の場合多い。
(run-sources+ T
(stop-when-fields-decayed dt Ex (vector3 x y z) 1e-3))
とすると、T秒間シミュレーションをとりあえず回した後に、(x,y,z)におけるdt秒間のExの最大値が、全ての時間を通しての最大値と比較して1e-3倍以下になるまで回してくれる。 しかしこの関数は落とし穴があって、非常に長時間共振するような構造の場合は、まだ、電磁波の発振が続いているのにシミュレーションを中断する場合がある。 電磁波が発散しきっているかどうかの判断は、h5topngなどで電磁波を画像化してを直接確認するしかない。
下のリンクの方法で取り組んでいました。理論値を複数のローレンツの重ね合わせで合わせるのはよくなさそうです。 https://www.mail-archive.com/meep-discuss@ab-initio.mit.edu/msg05147.html
マニュアルのトップページにschemeがobsoleteでpythonに変えろと書いてあったことと、メーリングリストで（バグはまだ多いですが）pythonに機能が追加されていることもあって、先週からpythonに移行のために動作検証しているところです。
PMLに垂直に入射するときは半波長。 以下、MEEP作成者のPMLに関するドキュメン(http://math.mit.edu/~stevenj/18.369/pml.pdf)からの引用
With PML, however, the constant factor is very good to start with, so experience shows that a simple quadratic or cubic turn-on of the PML absorption usually produces negligible refections for a PML layer of only half a wavelength or thinner [1,21]. (Increasing the resolution also increases the resolution alos increase the effectiveness of the PML, because it approaches the exact wave equation.)
斜め入射の時はPMLの吸収率が落ちるので、厚みを増す必要がある。 反射率は入射角をXとすると、ある定数δを用いて e^{-2*k*d*cos(X)*δ} と書ける。dはPMLの厚さ。kは波数。(詳細はhttps://ieeexplore.ieee.org/stamp/stamp.jsp?arnumber=477075を参照してください)。
高専ではmaterials.pyが使用できるが、tube,flexでは使用できていない(2019/01/20:使用できるようになった)。次のようなエラーが出る。
tube60:~$ python Python 3.6.7 | packaged by conda-forge | (default, Nov 21 2018, 03:09:43) [GCC 7.3.0] on linux Type "help", "copyright", "credits" or "license" for more information. >>> import meep as mp >>> from meep.materials import Au Traceback (most recent call last): File "<stdin>", line 1, in <module> ModuleNotFoundError: No module named 'meep.materials'
高専の環境では~/anaconda3/envs/mp/lib/python3.6/site-packages/meep/にmaterials.pyが存在するが、flex,tubeには存在せずこのことが原因の可能性がある。 flex,tubeの~/anaconda3/pkgs/pymeep-1.5-py36_nomkl_2が高専では~/anaconda3/pkgs/pymeep-1.6-py36_nomkl_2になっており、バージョンを1.5から1.6にすることで使用できるようになる可能性がある。
ネット上からmaterias.pyをダウンロードしてホームディレクトリに置いたのちに実行すると、実行はできるものの次のようなワーニングが出てくる。
>>> from materials import Au /home/guest/maeda/.conda/envs/mp/lib/python3.6/site-packages/meep/geom.py:175: RuntimeWarning: Epsilon < 1 may require adjusting the Courant parameter. See the 'Numerical Stability' entry under the 'Materials' section of the documentation warnings.warn(eps_warning, RuntimeWarning)
高専ではワーニングなしに実行できる。
conda list
とすると、
pymeep 1.5
となり、pymeepのバージョンが1.5であることが分かった。
オプションをつけてpymeepを1.5から1.6にupdateすると実行できるようになった。
(mp) maru@flex:~$ conda update -c chogan -c conda-forge pymeep
Solving environment: done
## Package Plan ##
environment location: /home/students/maru/anaconda3/envs/mp
added / updated specs:
- pymeep
The following packages will be downloaded:
package | build
---------------------------|-----------------
libctl-4.1.4 | 1 97 KB chogan
libgdsii-0.2.dev | 0 800 KB chogan
python-3.6.6 | hd21baee_1003 29.0 MB conda-forge
openssl-1.0.2p | h14c3975_1002 3.1 MB conda-forge
mpb-1.7.0 | nomkl_1 74 KB chogan
pymeep-1.7.0 | py36_nomkl_1 3.3 MB chogan
openblas-0.3.5 | h9ac9557_1000 15.8 MB conda-forge
nomkl-3.0 | 0 48 KB
ca-certificates-2018.11.29 | ha4d7672_0 143 KB conda-forge
certifi-2018.11.29 | py36_1000 145 KB conda-forge
------------------------------------------------------------
Total: 52.5 MB
The following NEW packages will be INSTALLED:
libgdsii: 0.2.dev-0 chogan
nomkl: 3.0-0
The following packages will be UPDATED:
ca-certificates: 2018.03.07-0 --> 2018.11.29-ha4d7672_0 conda-forge
certifi: 2018.11.29-py36_0 --> 2018.11.29-py36_1000 conda-forge
libctl: 4.1.2-0 chogan --> 4.1.4-1 chogan
mpb: 1.6.2-4 chogan --> 1.7.0-nomkl_1 chogan [nomkl]
openblas: 0.2.20-8 conda-forge --> 0.3.5-h9ac9557_1000 conda-forge
pymeep: 1.5-py36_nomkl_2 chogan [nomkl] --> 1.7.0-py36_nomkl_1 chogan [nomkl]
The following packages will be DOWNGRADED:
openssl: 1.1.1a-h7b6447c_0 --> 1.0.2p-h14c3975_1002 conda-forge
python: 3.6.8-h0371630_0 --> 3.6.6-hd21baee_1003 conda-forge
Proceed ([y]/n)? y
Downloading and Extracting Packages
libctl-4.1.4 | 97 KB | ########################################################################################### | 100%
libgdsii-0.2.dev | 800 KB | ########################################################################################### | 100%
python-3.6.6 | 29.0 MB | ########################################################################################### | 100%
openssl-1.0.2p | 3.1 MB | ########################################################################################### | 100%
mpb-1.7.0 | 74 KB | ########################################################################################### | 100%
pymeep-1.7.0 | 3.3 MB | ########################################################################################### | 100%
openblas-0.3.5 | 15.8 MB | ########################################################################################### | 100%
nomkl-3.0 | 48 KB | ########################################################################################### | 100%
ca-certificates-2018 | 143 KB | ########################################################################################### | 100%
certifi-2018.11.29 | 145 KB | ########################################################################################### | 100%
Preparing transaction: done
Verifying transaction: done
Executing transaction: done
pymeep-parallelにおいてmaterials.pyを使用したいときは、以下のように1.7.0-py36_nomklh39e3cac_1にアップデートしてください。
(pmp) maru@tube60:~$ conda update -c chogan -c conda-forge pymeep-parallel
Solving environment: done
## Package Plan ##
environment location: /home/students/maru/anaconda3/envs/pmp
added / updated specs:
- pymeep-parallel
The following packages will be downloaded:
package | build
---------------------------|-----------------
libgcc-ng-7.3.0 | hdf63c60_0 6.1 MB conda-forge
pymeep-parallel-1.7.0 |py36_nomklh39e3cac_1 3.3 MB chogan
------------------------------------------------------------
Total: 9.4 MB
The following NEW packages will be INSTALLED:
libgdsii: 0.2.dev-0 chogan
nomkl: 3.0-0
The following packages will be UPDATED:
ca-certificates: 2018.8.24-ha4d7672_0 conda-forge --> 2018.11.29-ha4d7672_0 conda-forge
certifi: 2018.8.24-py36_1 conda-forge --> 2018.11.29-py36_1000 conda-forge
libctl: 4.1.2-0 chogan --> 4.1.4-1 chogan
libgcc-ng: 7.2.0-hdf63c60_3 conda-forge --> 7.3.0-hdf63c60_0 conda-forge
mpb: 1.6.2-4 chogan --> 1.7.0-nomkl_1 chogan [nomkl]
openssl: 1.0.2p-h470a237_0 conda-forge --> 1.0.2p-h14c3975_1002 conda-forge
pymeep-parallel: 1.5-py36_nomklh39e3cac_4 chogan [nomkl] --> 1.7.0-py36_nomklh39e3cac_1 chogan [nomkl]
Proceed ([y]/n)? y
Downloading and Extracting Packages
libgcc-ng-7.3.0 | 6.1 MB |
################################################################################################################################################################################################## | 100%
pymeep-parallel-1.7. | 3.3 MB |
##################################################################################################################################################### ############################################# | 100%
Preparing transaction: done
Verifying transaction: done
Executing transaction: done
以下のようなエラーが出た場合は、計算を行っているディレクトリの書き込み・変更を許可してください。
Working in 3D dimensions.
Computational cell is 20 x 20 x 20 with resolution 10
block, center = (0,0,0)
size (1,3,1)
axes (1,0,0), (0,1,0), (0,0,1)
dielectric constant epsilon diagonal = (1,1,1)
time for set_epsilon = 11.7491 s
-----------
creating output file "./hogehoge.h5"...
on time step 1 (time=0.05), 14.4708 s/step
HDF5-DIAG: Error detected in HDF5 (1.10.2) MPI-process 0:
#000: H5F.c line 445 in H5Fcreate(): unable to create file
major: File accessibilty
minor: Unable to open file
meep: error on line 497 of h5file.cpp: error opening HDF5 output file
application called MPI_Abort(MPI_COMM_WORLD, 1) - process 1
flexのGLIBCのバージョンが古いので実行できていない。t70に移動して実行するとよい。
(pmp) maeda@tube60:/abinitio2/maeda/pymeep/Rectangle/drude/Omote$ mpirun -np 5 python after-Rec-drude-Au.py 45
Working in 3D dimensions.
Computational cell is 20 x 20 x 20 with resolution 10
block, center = (0,0,0)
size (1,3,1)
axes (1,0,0), (0,1,0), (0,0,1)
dielectric constant epsilon diagonal = (1,1,1)
time for set_epsilon = 5.05642 s
lorentzian susceptibility: frequency=3.47141, gamma=2.01155
lorentzian susceptibility: frequency=2.39466, gamma=0.701702
lorentzian susceptibility: frequency=0.66944, gamma=0.278261
lorentzian susceptibility: frequency=0.33472, gamma=0.19438
drude susceptibility: frequency=1e-10, gamma=0.0427474
-----------
creating output file "./after-Rec-drude-Au-ex-surface-45deg.h5"...
creating output file "./after-Rec-drude-Au-ey-surface-45deg.h5"...
creating output file "./after-Rec-drude-Au-ez-inner-45deg.h5"...
creating output file "./after-Rec-drude-Au-ez-outer-45deg.h5"...
Fatal error in PMPI_Waitall: Other MPI error, error stack:
PMPI_Waitall(405)...............: MPI_Waitall(count=92, req_array=0x9adaad0, status_array=0xc9e3670) failed
MPIR_Waitall_impl(221)..........: fail failed
PMPIDI_CH3I_Progress(623).......: fail failed
pkt_RTS_handler(317)............: fail failed
do_cts(662).....................: fail failed
MPID_nem_lmt_dcp_start_recv(302): fail failed
dcp_recv(165)...................: Internal MPI error! Cannot read from remote process
Two workarounds have been identified for this issue:
1) Enable ptrace for non-root users with:
echo 0 | sudo tee /proc/sys/kernel/yama/ptrace_scope
2) Or, use:
I_MPI_SHM_LMT=shm
Fatal error in PMPI_Waitall: Other MPI error, error stack:
PMPI_Waitall(405)...............: MPI_Waitall(count=85, req_array=0xa328890, status_array=0xd21c7f0) failed
MPIR_Waitall_impl(221)..........: fail failed
PMPIDI_CH3I_Progress(658).......: fail failed
MPID_nem_handle_pkt(1439).......: fail failed
pkt_RTS_handler(317)............: fail failed
do_cts(662).....................: fail failed
MPID_nem_lmt_dcp_start_recv(302): fail failed
dcp_recv(165)...................: Internal MPI error! Cannot read from remote process
Two workarounds have been identified for this issue:
1) Enable ptrace for non-root users with:
echo 0 | sudo tee /proc/sys/kernel/yama/ptrace_scope
2) Or, use:
I_MPI_SHM_LMT=shm
Fatal error in PMPI_Waitall: Other MPI error, error stack:
PMPI_Waitall(405)...............: MPI_Waitall(count=85, req_array=0x9c2b8f0, status_array=0xcc31ed0) failed
MPIR_Waitall_impl(221)..........: fail failed
PMPIDI_CH3I_Progress(658).......: fail failed
MPID_nem_handle_pkt(1439).......: fail failed
pkt_RTS_handler(317)............: fail failed
do_cts(662).....................: fail failed
MPID_nem_lmt_dcp_start_recv(302): fail failed
dcp_recv(165)...................: Internal MPI error! Cannot read from remote process
Two workarounds have been identified for this issue:
1) Enable ptrace for non-root users with:
echo 0 | sudo tee /proc/sys/kernel/yama/ptrace_scope
2) Or, use:
I_MPI_SHM_LMT=shm
解決策 エラーメッセージの2つめの方法で解決した
export I_MPI_SHM_LMT=shm
を実行させてから
(pmp)$mpirun -np N python hogehoge.py
で並列計算してくれた
git clone https://github.com/stevengj/meep ./autogen.sh make
とするとコードの完全体が得られる。 libctlというschemeの自作ライブラリが読み込まれているため読解が困難になっている。
Autotools で処理されるべきファイルに .in という拡張子が使われます。 つまり、その段階ではまだ Scheme のコードのテンプレートであり、処理されることで Scheme のコードとして妥当なものになります。ビルドプロセスによりますが、通常は configure が済んだ段階では処理済みのもの (この場合は meep.scm) が出来ていると思います。 もちろん make が済んだ時点なら確実に出来ているでしょう。
(define-class A B (define-property C D E))とすると、Aという名前のクラスが定義される。Bはオーバーライドするクラスで、BのプロパティをAは引き継ぐ。Bをno-parentにすると、オーバーライドは起きない。 Aに対してCというプロパティが設定され、その初期値はDとなる。またCの型はEとなる。 以下にmeep.scm.inに含まれるクラスの定義の例を示す。
(define-class susceptibility no-parent (define-property sigma-offdiag (vector3 0 0 0) 'vector3) (define-property sigma-diag no-default 'vector3)) (define-class lorentzian-susceptibility susceptibility (define-property frequency no-default 'number)
株式会社EEMが無料で公開している、OpenFDTDとOpenMOMはOpenFDTD DLサイト、OpenMOM DLサイトからダウンロードすることができます。
使い方については、DLサイトならびに、研究室が所蔵しているRFワールド No.39 2017年 08月号 に記載されています。
書籍の方は現在（2018/6/12)丸岡が借りています。
