S - 東京大学

講義「物理工学実験技法(A)」第13回 (担当 藤堂)
• 計算機実験の技術
• 物性物理におけるシミュレーション手法
• 数値対角化 / モンテカルロ法 / ・・・
• パッケージソフトウェアの利用
• 物質科学シミュレーションのポータルサイト
• MateriApps - http://ma.cms-initiative.jp
• MateriApps Live!
• ALPSシミュレーションパッケージ
• スクリプト言語の利用
• バージョン管理システム
講義資料: http://todo.issp.u-tokyo.ac.jp/ja/lectures
1
東京大学物性研究所計算物質科学研究センター神戸分室
＠理化学研究所 計算科学研究機構 AICS
SPring-8, SACLA
Kobe Port Island
AICS
2
理化学研究所 計算科学研究機構 AICS
2010年5月撮影
3
東京大学物性研究所計算物質科学研究センター神戸分室
AICS R501
Staff：Todo, Kitaura (Kobe Univ.), Sakashita
Yamashita, Matsushita
4
スーパーコンピュータ「京」
• 2011年6月、2011年11月 世界一位。2012年6月完成 @ 神戸
• 計算能力: 毎秒1京回 = 1016 flops = 10 Pflops
5
HPCI戦略プログラム：戦略5分野
6
物質科学分野における大規模シミュレーション
応用
（奔流）
デバイス
エネルギー変換
分子機能
重点；ナノ構造電子機能予測
重点；燃料電池物質
重点；全原子Simによる
ウイルス分子科学
高性能リチウムイオン電池
電子的～機械的性質の
マルチ・シミュレーション
スピントロニクス、
マルチフェロイックス材料
重点；水素・メタンハイドレート
生成、安定
ナノ・生体系反応制御・
化学反応ダイナミクス
太陽電池高効率長寿命化
生体分子高次構造と機能
バイオマス酸素反応解析
光機能分子と非線形外場応答
分子の光物性
光誘起電子ダイナミクスと
光・電子機能量子デバイス
高効率エネルギー
変換構造用材料
新材料探索
ナノ構造体材料の高効率非平
衡エネルギー変換
重点；相関の強い量子系新量子相探求とダイナミックス
電子相間の強い物質
の新機構解明
強相関電子系の
励起ダイナミックス
新しい量子相・量子臨界現象
基礎
12
（源流）
ポリモルフ生起分子集団機能
重点；電子状態・動力学・熱揺らぎ融和と物質理論新
凝縮分子科学系揺らぎとダイナミックス
分子における電子の動的
過程と多体量子動力学
分子の微細量子
構造予測
新量子相・新物質
計算物質科学イニシアティブCMSI: http://cms-initiative.jp
7
物性物理における代表的なシミュレーション手法
• 密度汎関数法 (density functional theory)
• 経験的なパラメタを含まない：「第一原理計算」
• 相関・ゆらぎが強い場合には破綻
• モデル計算
• 有効模型(ハイゼンベルグ模型など)に対する計算手法
• 相関・ゆらぎの効果を正確に取り入れた大規模計算が可能：「第一原理的」
• 数値対角化(厳密対角化) (ED)
• 古典/量子モンテカルロ法 (MC/QMC)
• 密度行列くりこみ群 (DMRG)
• 動的平均場近似 (DMFT)
• 量子化学計算、分子動力学
8
数値対角化
• 有限の大きさの系のハミルトニアンを行列の形で書き下し、数値的に対角化
• 例：ハイゼンベルグ模型
H=J
• スピン数 N
j,k
Sj · Sk = J
Sjz Skz
j,k
1 +
+ (Sj Sk + Sj Sk+ )
2
行列のサイズ (2S+1)N x (2S+1)N
• 例：S=1/2, N=20
行列のサイズ 1048576 x 1048576
メモリ 8TB
• どうやって対角化するか
• 完全対角化：Householder法 (三重対角化) + QR分解など (LAPACKなどを使う)
• 計算量が行列の次元の3乗に比例
• スピンを1つ増やすと、次元は2倍、メモリは4倍、計算時間は8倍!
9
数値対角化における計算上の工夫
• 行列のサイズを小さくする
• 保存量(ハミルトニアンと可換な量)を使う (total momentam, total S, etc)
• ハミルトニアンは保存量毎にブロック対角化。それぞれを独立に対角化
• 例：total Sz
N=20 の系で total Sz=0 の部分空間の次元
20 C10
20!
=
= 184 756
10! 10!
220 = 1 048 576
• 全固有値・固有ベクトルではなく、基底状態と低励起状態だけを求める
• ハミルトニアン行列はかなり「疎」 (1行/1列の非零要素の数 ∼ ボンド数)
• べき乗法の利用： Hn v
0
(絶対値最大の固有ベクトルに収束)
• 行列の非零要素のみ保存、あるいは行列要素をその場で計算
• べき乗法の改良：Lanczos法、CG (共役勾配)法
10
超大規模数値対角化の例 (Nakano et al 2013)
• (歪んだ)三角格子 S=1 反強磁性ハイゼンベルグ模型
次元 712 070 156 203
J. Phys. Soc. Jpn.
メモリ15TB (情報基盤センター FX10 4400ノード)
1
msq
diag
• N=27
0.5
非磁性相
0
120度構造
0
0.2
0.4
0.6
0.8
1
J2/J1
sq
11
モンテカルロ法の種類と応用
• モンテカルロ積分
• 単純サンプリング, 重点的サンプリング
• マルコフ連鎖モンテカルロ (高次元分布の生成, 人工的な「時系列」の生成)
• 平衡状態のシミュレーション: 古典系, 世界線量子モンテカルロ
クラスターアルゴリズム, 拡張アンサンブル法
• 非平衡状態, ダイナミクス, 基底状態の探索(焼き鈍し法)
• 逐次モンテカルロ
• 転送行列モンテカルロ, グリーン関数モンテカルロ
• 格子 QCD, タンパク質の折り畳み, 組み合わせ最適化問題, ベイズ推定, 生態系, 金
融工学, 在庫管理, 交通流, etc
12
イジング模型 Ising Model
• 格子上に規則的に並んだ磁気モーメント(スピン)
• ハミルトニアン Hamiltonian
H=
J
e
= ±1
H
i j
-J
J
H
-H
i
i
i,j⇥
• 分配関数 partition function
Z=
i
E (k)
k
• 状態 k の出現確率 = ボルツマン重み e
E (k)
• 物理量 physical quantities
• ギブス自由エネルギー Gibbs free energy f =
• 内部エネルギー internal energy
E=
1
ln Z
N
⇥
Z=
⇥
k
E
(k)
e
Z
E (k)
13
イジング模型に対するマルコフ連鎖モンテカルロ
メトロポリス法 Metropolis
for (int m = 0; m < total_mcs; ++m) {
// loop over Monte Carlo steps
for (int s = 0; s < num_sites; ++s) {
// loop over lattice sites
double delta = 0;
for (int j = 0; j < num_neighbors; ++j) {
int v = neighbor(s, j);
delta += 2 * J * spin[s] * spin[v]; // calculate energy difference
}
if (random() < exp(-beta * delta))
// accept/reject trial configuation
spin[s] = -1 * spin[s];
}
// measure physical quantities
}
熱浴法 heat bath
if (random() < (1 + tanh(-beta * delta / 2)) / 2)
spin[s] = -1 * spin[s];
14
三次元立方格子強磁性イジング模型
3
比熱
2.5
10000
L=4
8
12
16
24
32
2
L=4
8
12
16
24
32
1000
100
1.5
1
10
0.5
帯磁率
1
0
3
3.5
4
4.5
1
5
5.5
6
3.5
4
4.5
1000
L=4
8
12
16
24
32
0.8
3
5
5.5
6
5.5
6
L=4
8
12
16
24
32
100
0.6
0.4
10
0.2
相関長
(磁化)2
0
3
3.5
4
4.5
5
5.5
6
1
3
3.5
4
4.5
5
15
イジング模型における相転移 Phase Transition
F =E
TS
T=0.995Tc
T=Tc
エネルギー利得
秩序状態
ordered state
T=1.05Tc
エントロピー利得
臨界点 critical point
無秩序状態
disordered state
16
• 量子系の分配関数と物理量の期待値
Z = tre
τ=β
τ7
H
trAe
A⇥ =
Z
b=1
H
τ6
b=2
τ2
b=1
• モンテカルロ法を使うには物理量の期待値を
A(c)W (c)
b=4
τ5
τ4
• 量子系では H, A 等は全て演算子(or 行列)
A⇥ =
τ8
b=4
τ3
b=3
τ1
b=3
τ=0
c
連続虚時間経路積分表示
の形で書きたい (A(c) も W(c) も c-数)
• 量子系の場合の「配位」c と対応する「重み」W(c) は何か?
• 経路積分や高温展開を利用して「配位」や「重み」を定義する
表示・表現 (representation)
• メトロポリス法・熱浴法・拡張アンサンブル法・クラスターアルゴリズムの利用
17
虚時間経路積分表示
古典M粒子 (ポリマー)
量子1粒子
M
4
3
2
1
T
18
量子格子模型 (quantum lattice models)
• 量子スピン模型
!
J xy ! + −
H=
(Si Sj + Si− Sj+ ) + J z
Siz Sjz
2
!i,j"
!i,j"
• ハバード模型 (量子格子気体)
H = −t
!
(c†iσ cjσ
+ h.c.) + U
!i,j"σ
!
i
ni↑ ni↓
τ=β
τ7
• 虚時間経路積分表示が可能 (Suzuki 1976)
b=1
• 粒子がジャンプする場所と(虚)時刻だけを覚えて
おけば良い
• 粒子数は保存
τ8
b=4
τ6
τ5
τ4
b=2
連続虚時間表示 (Beard-Wiese 1996)
τ2
世界線は途中で途切れない
b=1
b=4
τ3
b=3
τ1
b=3
τ=0
19
負符号問題 (negative sign problem)
• フェルミ粒子は粒子の交換に対して、波動関数の符号が入れ替わる
フェルミ統計・フェルミの排他律
• 虚時間経路積分表示した際に粒子が交換すると重みが負になる
• フラストレーションのある量子スピン模型の場合にも負の重みが現れる
• 重みが 0 以下になるので、確率として解釈することができない
• 量子モンテカルロにおける最大の問題 !
3
1
2
3
1
2
20
Spin ladder material Na2Fe2(C2O4)3(H2O)2
J. Phys. Soc. Jpn., Vol. 78, No. 12
ions in octahedral crystal field
O
effective S=1 spins at low T
C
a
H. YAMAGUCHI et al.
b
(c)
c
b
H. YAMAGUCHI et al.
c
Fig. 1. (Color online) (a) Crystal structure of Na2 Fe2 (C2 O4 )3 (H2 O)2 .
(a)
H =sodium
100 Oe
Hydrogen
and
atoms are omitted for clarity. Fe2þ ions are
connected by oxalate molecules and form two-leg ladders along the
0.20 a-axis, as shown by red broken lines. (b) Distorted FeO6 octahedra
H//a
bridged by oxalate ion. (c) Arrangement of the ladders in the bc-plane,
and red broken lines corresponding to rungs of the ladders.
Fig. 2. Tem
(C2 O4 )3 (H
0.15
2.2 Experimental methods
0.10 The magnetic susceptibilities of a single crystal of SIO
were measured at 100 Oe with a superconducting quantum
interference device (SQUID) magnetometer
(Quantum
H⊥a
Design MPMS-XL7). High-field
experimentmagnetizations in pulsed
0.05
calculation
magnetic fields up to about
53 T were measured utilizing a
non-destructive pulse magnet, and the magnetizations were
detected
5 by an10induction
15 method
20 using25a standard
30 pick-up
coil system.
Temperature (K)
ESR measurements on a single crystal of SIO in pulsed
magnetic
(b) fields up to about 53 T were carried out using our
non-destructive
4pulse field ESR apparatus equipped with aH//a
pulse magnet, a far-infrared laser, and several kinds of Gunn
Oscillators for the frequency range from about 50 GHz
Fig. 3. Ma
Honda, Okunishi , Todo, Kindo, Hagiwara (2009)
and H ? a
/ion)
B/ion)
Susceptibility (emu/mol)
the
1623.6 GHz
Susceptibility (emu/mol)
QMC and the DMRG methods, respectively. For
(a) z
these
calculations,
we
use the following parameter values,
0.20
1288.1 GHz
gz ¼ 4:7, gx ¼ 2:0, gy ¼ 2:5, !a ¼ 0:30, and R ¼ 0:1{ 0:5 for
"
977.2 GHzH k a (" ¼ 24" , ’ ¼ 90" ) and H ? a (" ¼ 102" , ’ ¼ H//
45
).
a
y
0.15
847.0 GHzThe principal axis directions and their g-values are deterx
mined from the saturation
magnetizations and the slopes of
730.5 GHzthe ESR resonance
0.10 modes in the following isolated dimer
model analysis. In these field directions, the two sites of the
H⊥ a
H//a ladder are magnetically equivalent. In the calculation, we
0.05
add
the
temperature-independent
Van Vleck paramagnetism
T=1.3 K
H = 100 results
Oe
to both #k and #? . The calculated
do notexperiment
reproduce
ones well, especially the upturn of #k
0
50the experimental0.00
15 Van
20 Vleck
25 para30
at approximately 5 K. In5 Fig. 10
6(b), the
Temperature (K)
manetic contributions are subtracted from both Mk and M? .
)3 (H2 O)2 at 1.3 K for
(b)
Wesignals,
find a small anomaly
at half of the saturation magnetng and weak
z
H// a
Kimura,
ization of R ¼ 0:14 calculated Yamaguchi,
for H k a, which
can
be
considered to originate
from
y a 1/2 magnetization plateau.
24°
The calculated results
also poorly explain the experimental
3
c
Fe
for several theoretical models (ladders,
dimers, etc)
a
b
• Fitting experimental data by QMC results
J. Phys. Soc. Jpn., Vol. 78, No. 12
(b)
(a)
Susceptibility (emu/mol)
•
Fe2+
21
Orbital ordering in eg orbital systems
• Mott insulators with partially filled d-shells
• Non-trivial interplay of charge, spin, and
|x2
orbital degrees of freedom
y2 i
|3z 2
r2 i
• Effective Hamiltonian for orbital degrees of freedom (120o model)
X 1h
z
x
H120 =
Jz Tiz Ti+
+ 3Jx Tix Ti+
4
i, =x,y
i X
p
z
x
z
Jz Tiz Ti+z
± 3Jmix (Tiz Ti+
+ Tix Ti+
)
i
1.5
1
0.5
T x polarized
(xy planes)
T z polarized
0
"
0
"
#
!
#
2
!
4
0
!
0.2
0.4
!
0.6
0.8
!
1
Jx / Jz
1.2
1.4
1.8
1.6
2
Jmix / Jz
orbiton gap !
2
van Rynbach, Todo, Trebst (2010)
!
22
YAMAMOTO,
ated by SSE. System sizes
are from LTODO,
= 8 to AND
L = 16.MIYASHITA
denote linear extrapolations of data points.
0.3
1.4
PHYS
B.1.2
Mean-fi
MOTO, TODO, AND MIYASHITA
PHYSICAL REVIEW B 79, 094503 $
1
We also analyze
the
Therefore we confirm the
the su- Bose-Hubbard Hamiltonian
We existence
analyze theofextended
on a
0.2
SUCCESSIVE PHASE TRANSITIONS
AT FINITE… 22 In order
approximation.
se in the region of
!
"
1
/
2.
cubic
lattice
II. MODEL
B. Mean-field approximation 0.8
sublattice structure which
0.15 0.25
Sπ(L=12)
by the
SSE
We nalso
theNLordering
processes
derTTSof
density.0.6by
Herethe
w
= − t $a†i on
a j +a aia†j % + V
analyze the extended Bose-Hubbard H
Hamiltonian
π
in j analyze
by
SSE
ρ
s
TRANSITIONS
ATsoft-core
FINITE TEMPERATURES
• Interacting
bosons
π(L=10)
"ij#
"ij#
0.1 220.2SIn
order to order
study
the superfluid
solid state,order
we
approximation.
and
attice
0.4
ρs(L=12)
sublattice structure
which is characterized
a staggere
tonian for thisbyMF
is given
ρ
(L=10)
0.15
1
study the ordered† states †at finite temperatures.
s
SF
− !of the
n0.05
$1%NSwe adopt mean
fields0.2
for the
H = − t $ai a j + aia j % + V nin j + U ni$ni − 1%der
i , density. Here
II. MODEL
0.25
Sπ
!
ρs
Supersolid in extended Bose-Hubbard model
!
ρs
Sπ
T/U
!ithe
!i
!
!
ansition between
the
normal
liquid
phase
and
2
=
H
0.1
MF H
"ij#
"ij#
order and superfluid
order
at
sublattices
A
and
B.
The
0
0
SS
s expected to belong to the universality
class
of
†
tonian
for thisoperators
MF is0.1
given
by
0.05
0
0.2
0.3
0.4
1
and
a
are
the
creation
and
annihilation
of
where
a
i
i
odel •because
the
phase
transition
occurs
in
the
Supersolid
co-existence
†of
†
+ U! =
ni$n
1% − !$&a
n
,' = " % and n$1%
!
ia −boson
i
†
(b)
,
a
=
a
T/U+ HA +
i j
ij with i
=
−
zt$a
i ai. The parameter t 0denotes
2 i of the orderi parameter
symmetry breaking
=
H
C,
H
A + a A% & B +
0.04 MF
0.045 A 0.05 B 0.055
0.06
the
hopping
matrix
element,
long-range
(=solid)
and U and V are the on-site and
. On the other
hand,
theorder
phase
transition
y† of Z2diagonal
t/U
!
is
the
nearest-neighbor
repulsions,
respectively,
and
0.3
1.4
ai are thetocreation
and
annihilation
operators of
dai isand
expected
belong
to
the
XY
universality
long-range
order The
(=superfluid)
†
†
FIG.
6.sum
The
t-T phase
diagram for V / UU
= 1 / z and ! / U = 0.7
† over
chemical
potential.
notation
"ij#
means
the
the
aoff-diagonal
'
=
"
%
and
n
=
a
a
.
The
parameter
t
denotes
ne the
$&ai ,order
ij
i
i
=
−
zt$a
+
a
%
&
+
zVn
m
+
1% T−S !n
H
j
i
S
(L=12)
A 0.25
B
A$n
A −state
1.2
parameter has the symmetry of U!1".
A TheAtransition
π A B
obtained
by
temperatures
of n
the
solid
3 SSE.
"
2
,
where
L
is
nearest-neighbor
Theand
system size are
is N
=
L
pping• Experimental
matrix element,realization:
U and V are
thepairs.
on-site
S
(L=10)
optical
lattice?
plotted by solid circles and those
T#s are
π of the superfluid state
n, we study the temperature
dependence
of theseThe order parameter
†
ρ
(L=12)
1
the
length
of
the
system.
of
the
solid
HB the
= −data
zt$a
&A +
-neighbor repulsions, respectively, and ! is the
plotted by open circles. The liness connect
points
for
B%the
B+a
0.2
eters.
ρs(L=10)
guide to the eyes.
state
al potential. The notation
"ij#ismeans the sum over the
0.8
3
-neighbor pairs. The system size is N = L , where L is
0.15
U solid for0.6!t / U
1
successive transitions
of superfluid and
†
iQ·$r
−r
%
j
k
gth A.
of Stochastic
the system.
The expansion
order parameterS#of= the2 !
solid
series
+ aBAs
%&
zVnBinmthe
− !n
e
"nH
#,= − zt$a
$2%
A + seen
A +SSE n
B$nB − 1%here
jn=kB0.045,
!/UB
= 0.7".
was
simulation,
2
N jk
0.1
C =su2zt&
we find again the suppression of the solid order by 0.4
the
show the results obtained by SSE. The simulaperfluid fraction. Namely, S" has a cusp at the superfluid
0.05
where
Q
=
$
#
,
#
,
#
%
is
the
wave
vector
that
represents
thedepictwhere
z$=6%
is and
the
numbe
erformed in the1grandiQ·$r
canonical
ensemble using
0.2comtransition
point. We also
the phase
diagram
−r
%
j k "n n #,
e3.staggered
$2% parameter
As for the order
superfluid
!
# = N2=
j order.
k
pare of
thatthe
to that
of SSE !Fig. 8".
a qualitatively
andshow
mB are
the mean fi
mAThey
12
es N = 103Sand
N jk
0
0!tc ,numbe
&Ais&tation
zVm
Ce.g.,
= 2zt
T c"
good fraction
agreement,
there
the
we adopt
Bose-Einstein
B − tetracritical
AmB ,point
values
of
the
that of condensation
the order parameters state,
of solid,
S#, and
0 temperatures
0.1
0.2
0.4
and the critical
TS" and T0.3
#s are suppressed by
spectively,
= $#
, # , #% is of
thetemperature
wave vector
represents
, as
a function
forthat
various
values
z$=6% isof the
number
of
nearest-neighbor
sites.
sQ
appearance
the
other
order
as
before.
(c)
1 the † where
T/Umentioned
†
$k=0superfluid
= 2 ! "a j akm+Aakand
a j#Let
$3% fields
ed
order.
forthe
the transition
order parameter
of the
http://www.uibk.ac.at/th-physik/qo
5!a",
we As
show
from normal
liquid
studythe
the mean
competition
between
the solid and to
supermBus are
corresponding
the
e
m
=
"n
#
N
A
A
ordersofbythe
analyzing
Ginzburg
Landau
Yamamoto,
Todo,theMiyashita
(2009)
we
condensation
fraction con- jk
tation
number
operators
for
A !GL"
and free
B site
lidadopt
for t /Bose-Einstein
U = 0.02 in which
only S# appears
FIG. fluid
5.values
Temperature
dependence
of order
parameters
obt
$
take
small
value
energy. Since the order parameters
m
and
&normal
Similarly,
A and &
B corre
the same way, the transition
from normal
liquid whilebyinspectively,
SSE
!V
= Uvicinity
/simulation
z". !a"of !tthe/ Utetracritical
= 0.02, $
/ U = 0.7":
liquid-no
in
the
mean-field
analysis,
the
SSE
we
in
the
point,
the
GL
free
energy
is 23
1
パッケージソフトウェアの利用
• 世の中には物性物理のための様々なパッケージソフトウェアが存在する
• 有償：Gaussian、VASP、WIEN2k、、
• 無償：GAMESS、Quantum ESPRESSO、ABINIT、Amber、Gromacs、ALPS、、
• 日本国内でも数多くの高性能なソフトが開発・公開されている
• 小さなグループでこつこつと開発しているソフトが多い
• 知名度は高くない
• ドキュメント、宣伝、ユーザインタフェース、ユーザサポートなどの問題
• 開発者から見ると
• ソフトを開発・公開しただけでは成果にならない (職がない)
• ドキュメント作成やユーザサポートには時間も手間もかかる
• 海外の大規模なソフトに対抗するのはしんどい
24
物質科学シミュレーションのポータルサイト
http://ma.cms-initiative.jp/
MateriApps
25
26
MateriApps Live! -
http://github.com/cmsi/MateriAppsLive
• USBメモリから直接ブートできる Linux システム (Debian Live)
• MateriAppsで紹介している公開ソフト・ツールがあらかじめインストールされている
• Windows、Mac などで利用可
• 2013年7月末 ver. 1 公開予定
• MateriApps サイトで配布
• 学会・アプリ講習会でも
配布予定
27
ALPS プロジェクト
• http://alps.comp-phys.org/
• 量子格子模型(量子スピン系、電子系)のための並列シミュレーションソフトウェアパ
ッケージの開発
• 開発物(ライブラリ・アプリケーション)をフリーウェアとして公開
• 様々な模型に対する様々なシミュレー
Parameter XML Files
Lattice XML File
Model XML File
ションアルゴリズムの入力・計算結果を
Quantum Lattice Model
XML 形式で統一的に入出力可能に
Quantum Monte Carlo
Exact Diagonalization
DMRG
Outputs in XML Format
28
The ALPS project
ALPS = Algorithms and Libraries for Physics Simulations
• International collaboration for developing open-source softwares for
simulation of quantum lattice models, such as quantum spin systems,
electron systems, etc
• ALPS Libraries = collection of generic C++ libraries
• ALPS Applications = collection of application packages using modern
algorithms such as QMC, DMRG, ED, etc
• ALPS Framework = environment for executing large-scale parallel simulations
including XML schemas, tools, scheduler, etc
29
Target audience
• Experimental physicists
• Use “canned codes” to model materials
• Determine microscopic parameters by fitting experimental data to
simulations
• Theoretical physicists
• Quick check of theoretical ideas using many modern algorithms
• Monte Carlo, Diagonalization, DMRG, ...
• Also useful for debugging
• Libraries simplify and accelerate code development
30
ALPS libraries and applications
tools
applications
domain-specific
libraries
numerics
XML manipulation GUI database MC QMC ED DMRG
lattice model observables scheduler
random ublas
iterative eigenvalue solver
Boost library
generic C++
C / Fortran
graph
serialization XML/XSLT
BLAS LAPACK MPI
31
ALPS Lattice XML
periodic chain with length L
<LATTICE name="chain lattice" dimension="1">
<BASIS><VECTOR> 1 </VECTOR></BASIS>
</LATTICE>
<UNITCELL name="simple1d" dimension="1" vertices=”1”>
<EDGE>
<SOURCE vertex="1" offset="0"/><TARGET vertex="1" offset="1"/>
</EDGE>
</UNITCELL>
<LATTICEGRAPH name = "chain lattice">
<FINITELATTICE>
<LATTICE ref="chain lattice"/>
<PARAMETER name=”L”/>
<EXTENT size ="L"/>
<BOUNDARY type="periodic"/>
</FINITELATTICE>
<UNITCELL ref="simple1d"/>
</LATTICEGRAPH>
32
Model XML for describing Hamiltonian
XXZ spin model with two types of bonds (e.g. nearest and next nearest neighbor interactions)
<HAMILTONIAN name=”spin”>
<PARAMETER name=”Jz” default=”J”/>
<PARAMETER name=”Jxy” default=”J”/>
<PARAMETER name=”J” default=”1”/>
<PARAMETER name=”Jz'” default=”J'”/>
<PARAMETER name=”Jxy'” default=”J'”/>
<PARAMETER name=”J'” default=”0”/>
<PARAMETER name=”h” default=”0”/>
<BASIS ref=”spin”/>
<SITETERM site=”i”> -h * Sz(i) </SITETERM>
<BONDTERM type=”0” source=”i” target=”j”>
Jz * Sz(i) * Sz(j) + Jxy * (Splus(x)*Sminus(y)+Sminus(x)*Splus(y)) / 2
</BONDTERM>
<BONDTERM type=”1” source=”i” target=”j”>
Jz' * Sz(i) * Sz(j) + Jxy' * (Splus(x)*Sminus(y)+Sminus(x)*Splus(y)) / 2
</BONDTERM>
</HAMILTONIAN>
!
H=
!i,j"
!
J
xy
(Si+ Sj− + Si− Sj+ )] −
[Jz Siz Sjz +
hSiz
2
i
33
Simulations with ALPS
Parameter XML Files
Lattice XML File
Model XML File
Quantum Lattice Model
Quantum Monte Carlo
Exact Diagonalization
DMRG
Outputs in XML Format
34
ALPS Version 2.1
The ALPS project release 2.0:
open source software for strongly correlated systems
• 2012年5月
• 論文 - JSTAT P05001 (2011)
• http://alps.comp-phys.org/
B. Bauer1 L. D. Carr2 H.G. Evertz3 A. Feiguin4 J. Freire5
S. Fuchs6 L. Gamper1 J. Gukelberger1 E. Gull7 S. Guertler8
A. Hehn1 R. Igarashi9,10 S. Isakov1 D. Koop5 P.N. Ma1
P. Mates1,5 H. Matsuo17 O. Parcollet12 G. Pawlowski13
J.D. Picon14 L. Pollet11,1 T. Pruschke6 E. Santos5
V.W. Scarola15 U. Schollwöck16 C. Silva5 B. Surer1 S. Todo17,10
S. Trebst18 M. Troyer1 ‡ M. L. Wall2 P. Werner1 S. Wessel19,20
Theoretische Physik, ETH Zurich, 8093 Zurich, Switzerland
Department of Physics, Colorado School of Mines, Golden, CO 80401, USA
3
Institut für Theoretische Physik, Technische Universität Graz, A-8010 Graz, Austria
4
Department of Physics and Astronomy, University of Wyoming, Laramie, Wyoming
82071, USA
5
Scientific Computing and Imaging Institute, University of Utah, Salt Lake City,
Utah 84112, USA
6
Institut für Theoretische Physik, Georg-August-Universität Göttingen, Göttingen,
Germany
7
Columbia University, New York, NY 10027, USA
8
Bethe Center for Theoretical Phyics, Universität Bonn, Bonn, Germany
9
Center for Computational Science & e-Systems, Japan Atomic Energy Agency,
110-0015 Tokyo, Japan
10
Core Research for Evolutional Science and Technology, Japan Science and
Technology Agency, 332-0012 Kawaguchi, Japan
11
Physics Department, Harvard University, Cambridge 02138, Massachusetts, USA
12
Institut de Physique Theorique, CEA/DSM/IPhT-CNRS/URA 2306, CEA-Saclay,
F-91191 Gif-sur-Yvette, France
13
Faculty of Physics, A. Mickiewicz University, ul. Umultowska 85, 61-614 Poznan,
Poland
14
Institute of Theoretical Physics, EPF Lausanne, CH-1015 Lausanne, Switzerland
15
Department of Physics, Virginia Tech, Blacksburg, Virginia 24061, USA
16
Arnold Sommerfeld Center for Theoretical Physics and Center for NanoScience,
University of Munich, Theresienstrasse 37, 80333 Munich, Germany
17
Department of Applied Physics, University of Tokyo, 113-8656 Tokyo, Japan
18
Microsoft Research, Station Q, University of California, Santa Barbara, CA 93106,
USA
19
Institute for Solid State Theory, RWTH Aachen University, 52056 Aachen,
Germany
20
Institut für Theoretische Physik III, Universität Stuttgart, Pfaffenwaldring 57,
70550 Stuttgart, Germany
1
2
‡ Corresponding author
http://www.iop.org/EJ/abstract/1742-5468/2011/05/P05001
35
スクリプト言語の利用
• スクリプト言語とは
• コンパイルなしでその場で実行できる(簡易)プログラミング言語
• メジャーなスクリプト言語：Perl、Python、Rubyなど
• Mac OS X や Linux には、あらかじめインストールされている
• Windows では、バイナリパッケージをダウンロード & インストール
• ソースの編集・コンパイル・インストールなどの手間がかからない
• その場で入力
即実行
出力
• 電卓代わりにも
• パッケージ(ライブラリ、モジュール)が充実しているので便利
• CSV/XML/HDF5の読み書き、DBへのアクセス、ネットワーク、CGI、機器の制御、文
字列処理、線形演算・最適化、計算・実験結果の加工、図の作成、GUI、、、
• いろいろなアイデアを手軽に試すことが可能
36
例) Monty Hall 問題
• ゲームのルール
• 三つのドアのうち、一つには車、他の二つにはヤギが隠されている
• 「挑戦者」はまずどれかの一つドアを選ぶ
• 「親切な司会者」は残り二つのうちから、ヤギ (= ハズレ)の入っているドアを開い
て見せてくれる
• 「挑戦者」には、現在選んでいるドアから残りもう一つのドアに変更する権利が
与えられる
• 最初のままが良いのか?
• もう一つのドアに変更すべきか?
37
Pythonによる乱数を使ったシミュレーション
swap
成功回数
import random
success_stay = 0
success_swap = 0
for i in range(1000):
bingo = random.randint(0,2)
choice = random.randint(0,2)
blank = []
for k in range(3):
if k != choice and k != bingo:
blank.append(k)
show = random.choice(blank)
remain = 0 + 1 + 2 - choice - show
stay
0
10
20
30
40
50
試行回数
if choice == bingo:
success_stay += 1
if remain == bingo:
success_swap += 1
print i+1, bingo, choice, show, remain, success_stay, success_swap
38
Buffon の針
• 等間隔 b で線が描かれた床の上に長さ a の針をばらまく (ただし a < b とする)
• 線に針がかかる確率は?
39
Buffon の針 (つづき)
• 最も近い線から、針の中心までの距離 x は 0 ∼ b/2 に均一に分布
• 中心のまわりの針の角度 φ は 0 ∼ π/2 に均一に分布
• x < (a/2) cos
P =
ならば、針は線にかかる Nhits (x, ) = 1
b/2
/2
dx 0 d Nhits (x,
0
b/2
/2
dx 0 d
0
)
a 2
= ·
b
• この確率過程を計算機中でシミュレーションすれば、円周率が求まる!
• モンテカルロ積分、モンテカルロ・サンプリング
40
乱数によるサンプリング − モンテカルロ法
• a = b = 1 のとき
import math
import random
success = 0
for i in range(10000):
x = 0.5 * random.random()
phi = (math.pi/2) * random.random()
if (x < 0.5 * math.cos(phi)):
success += 1
print i, success / float(i+1)
0.7
0.68
2
0.66
0.64
0.62
0.6
0.58
0.56
0.54
0.52
0.5
0
200
400
600
800
1000
41
バージョン管理システム(VCS)の利用
• バージョン管理システムとは？
• ファイルの履歴ををデータベース(リポジトリ)で一括管理す
るシステム
• もともとはプログラムのソースコードのためのシステム
• それ以外のファイル(例えば TeX ファイル)管理にも使える
• 一人で使っても複数人で使っても超便利
• 超優秀な秘書のようなもの
42
バージョン管理システム(VCS)の利用
• 作業者 and/or 作業場所が複数になると、ファイル名や手書
きのログファイルによるバージョン管理はすぐに破綻する
• ネットワーク経由でファイルを check out/check in
• 更新毎に一意なバージョン番号 (リビジョン) を付与
• 任意のバージョン間の比較が容易
• バックアップの代わりにも
• 複数箇所から同時に更新した場合に衝突を回避するしくみ
• ブランチ・マージ・タグ付けなどが可能
43
ありがちなパターン
25/45
1)
2)
A
3)
4)
B
B
A
University of Tokyo
渡辺宙志 CMSI計算科学技術特論A 講義資料より
44
バージョン管理している場合
1)
26/45
2)
3)
5)
6)
B
A
7)
A
4)
(
)
C
University of Tokyo
渡辺宙志 CMSI計算科学技術特論A 講義資料より
45
主なバージョン管理システム
• BitKeeper - かつて Linux のカーネルのソース管理に使われていた
• CVS (Concurrent Versions System) - ネットワークでの利用を考慮とした初めてのバ
ージョン管理システム。以前はよく使われていた
• Git - 現在 Linux の開発に使われている。分散型リポジトリ
• Mercurial - Git のライバル。分散型リポジトリ
• SCCS (Source Code Control System) - 70年代にベル研で開発された世界初のバージ
ョン管理システム。現在は使われない
• Subversion - CVSの改良版として開発された。現在最もポピュラー? Mac OS X や多
くの Linux には最初からインストールされている。Windows 版クライアント
(TortoiseSVN)もある
46
バージョン管理システムの欠点(面倒な点)
• 修正前に最新の状態にアップデートしなければならない
慣れると習慣になります
• 全ての修正を「コミット」しなければならない
慣れると習慣になります
• 衝突(コンフリクト)が発生した時に対処しなければならない
衝突に気づかずに修正してしまうほうが怖いです
• サーバのセットアップが面倒くさい
まずはホスティングサービス(sourceforge, github, bitbucket)を試してみましょう
まわりにいるプロ(?)に相談しましょう
• バージョン管理システムを使うと作業効率が倍以上になる
使わないと人生を半分損する (c) 2013 渡辺宙志
47
Subversion の利用に役立つ資料
• Subversion によるバージョン管理
• http://todo.issp.u-tokyo.ac.jp/archive/subversion/book.html • 「Subversion によるバージョン管理」の読み方
• http://todo.issp.u-tokyo.ac.jp/ja/members/wistaria/log/subversion-intro • 「Subversionって何?」 「リポジトリ」 「バージョン管理モデル」 「最初の
チェックアウト」 「基本的な作業サイクル」 「作業コピーの更新」 「作業コ
ピーに変更を加えること」 「自分の変更点の調査」 「変更点のコミット」
• CVS/Subversionを使ったバージョン管理（前編：バージョン管理の基礎）
• http://sourceforge.jp/magazine/08/09/09/1038233 • CVS/Subversionを使ったバージョン管理（後編：SVNを使ったバージョン管理）
• http://sourceforge.jp/magazine/08/09/24/113215 48