巻頭言

Vol.12 No.2 号の発刊にあたって
京都大学学術情報メディアセンター
岩下武史

本号では，学術情報メディアセンターの共同研究に関する報告を掲載いたします．昨年度に実施されました若手研究者の支援事業に関して 10 件，プログラム高度化支援事業について 3 件，
大規模計算の支援事業について 1 件の報告をそれぞれ掲載しています．このような報告を特に
若手研究者や学生の方にとっての新しい出会いの場として活用頂ければ幸いです．
これまでにセンターサービスにおける広報発行の意義に関し，WEB を中心としたインターネットによる情報提供と紙媒体としての広報誌の役割の差別化について継続的に検討を行ってきました．その中で，紙媒体の良さの一つは「自分にとって一見必要でない情報」に出会う機会が多いということだと考えています．即ち，自分が知りたいと思う情報について調べる場合，イン
ターネットは大変便利ですが，例えば自身の研究分野と異なる分野の研究についてネット上で調べたりすることは稀で，望まない限りそのような機会に乏しいといえます．一方，紙媒体の場合，ぱらぱらとめくるだけでも何からの情報が入ってくることになります．そこに，興味深い図や式，
単語を見出すことがあれば，その先に何か新しい繋がりが出てくるかもしれません．例えば，今号の中畑氏（愛媛大学）の原稿では，非破壊検査という応用の中で，有限積分法の GPU における効率的な実装法について述べています．本解析におけるプログラムは原稿中の式(4)～(6)に見られるように，反復型のステンシル計算となります．一方，反復型のステンシル計算には，他に拡散方程式の差分解析や電磁波解析等の波動解析で用いられる FDTD 法等があります．つまり，応用分野や解法が異なっている場合においても，反復型ステンシル計算を用いている研究では，本原稿で述べられている GPU 向けの実装法が有用となります．このように，異なる応用分野においても同一の解法やプログラム構造が用いられている場合には，これらの研究間において知見を相互に活用することが可能となります．今回の広報における研究報告がなんらかの形でユーザの皆様にとって有用な情報を提供できれば幸いです．
次に，現在のスパコンサービスについて述べさせていただきますと，今号で掲載されている運転状況に見られるように，主たるシステムである Camphor（システム A），Laurel（システム B）
について 60%を超える高いノード稼働率が続いています．結果として，ユーザの皆様にとりましてはなかなかジョブが走らないというご不満がある状況かと思います．こうした状況への対応策の一つとして，学術情報メディアセンターでは HPCI への資源提供を強化することを主目的として，新たな並列演算サーバを 2014 年度に導入，運用開始予定です．当該システムはメニー
コアプロセッサ（Intel Xeon Phi）を搭載した計算ノードをビルディングブロックとするクラス
タシステムとなります．旧来のシステムに加え，新システムについてもユーザの皆様にご活用頂
ければ幸いです．センター広報としましても最新のアーキテクチャに基づくプロセッサであるIntel Xeon Phi を使ったプログラミングや研究事例の紹介に取り組んでいきたいと考えております．今後も皆様の研究・教育活動に益するところがありますように取り組んで参りますので，センターシステムの利用を含め，どうぞよろしくお願いいたします．

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 1

スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

色素増感太陽電池を指向した新奇ポルフィリン系色素の構造と	
  	
  	
 
電子構造の解明
梅山有和
京都大学大学院工学研究科分子工学専攻

1 緒言

キサリン環を持たないポルフィリン（ZnP）と比較して、優れた光捕集能と太陽電池性能（η= 6.3%）

現在、化石燃料の枯渇と環境へ与える影響が深を実現した 2,3。本研究では、更なる性能の向上を刻な問題になってきており、太陽電池活用の拡充

目指し、電子供与性基のトリアリールアミンを有が期待されている。しかしシリコンに代表されるするイミダゾールを ZnPQ に縮環させた ZnPQI を

無機太陽電池は、製造コストが高く、広範に普及設計し、その合成と物性評価を行った 4。

させることに障害が生じている。一方、有機太陽
電池は低コスト、着色による意匠性、柔軟性、軽2 結果と考察量という利点が指摘されている。中でも色素増感
太陽電池はエネルギー変換効率（η）が高く、次世代のエネルギー源として期待されている。現在
までにルテニウム錯体と多孔性酸化チタンを用いた系において 11%を越える高い η 値が報告されている。しかし、ルテニウムは高価で、資源制約もある。そのため、ルテニウム錯体に変わる安価で高性能な色素の開発が求められている。そこで
我々は、安価な金属を利用できるポルフィリン錯体に注目し、新規ポルフィリン系色素の開発を行ってきた 1。特に、キノキサリンが縮環したポル

2.1 新規ポルフィリン色素の諸物性と色
素増感太陽電池への応用
ZnPQI および ZnPQ のジクロロメタン溶液中の吸収スペクトルを Figure 2 に示す。ZnPQI は ZnPQに比べて幅広い範囲で吸収を示し、吸収端は 680

nm に達することがわかった。また、電気化学測定により求めた ZnPQI の第一酸化電位（Eox = 0.92

V vs. NHE）は ZnPQ（Eox = 0.98 vs. NHE）と比較して低くなった。

フィリン色素（ZnPQ、Figure 1）を開発し、キノ
C6H13

Ar

N

Ar

C6H13

N

Ar

Ar

ZnP

N

Ar =
Ar
N

N Zn N

N

Ar
COOH

COOH

N
Ar

ZnPQI

N Zn NAr

N

N Zn N

N

Ar

Ar

N

N

ZnPQ

N

COOH

Ar

Figure 1. Structures of porphyrin compounds.

2 − Academic Center for Computing and Media Studies, Kyoto University

Figure 2. UV/Vis absorption spectra of ZnPQI (solidline) and ZnPQ (dashed line) measured indichloromethane.

ZnPQI を酸化チタンに吸着させ、色素増感太陽

密度汎関数法（DFT）による理論計算（B3LYP）（JSC = 13.9 mA cm 、VOC = 0.68 V、ff = 0.71）であGaussian03 プログラムを用いた。分子軌道の可視電池を作製して性能評価を行ったところ、η= 6.8%を行った。基底関数系に 3-21G(d)を選択し、った。これは、ZnPQ の場合（η= 6.3%、JSC = 13.2

化には Molstudio3.0 プログラムを用いた。Figure 4–2mA cm 、VOC = 0.71 V、ff = 0.67）と比べて高い値–2に ZnPQI および ZnPQ の最高被占軌道（HOMO）となっている。最適条件下での作用スペクトルをと最低空軌道（LUMO）における電子密度分布を

Figure 3 に示す。吸収の長波長化を反映し、650-750示す。一般に、アンカー部位の LUMO の電子密nm で ZnPQI の IPCE 値は ZnPQ に比べて高くなっている。

度が大きいほど、励起状態の色素と酸化チタンの3d 軌道との間の電子カップリングが大きくなる

ことが知られている。ZnPQI では ZnPQ と同様に、アンカー部位の LUMO の電子密度が大きくなっている。そのため、色素から酸化チタンの伝導帯への電子注入の効率が、ZnPQ と同様に ZnPQI を用いた系でも高くなっており、かつ吸収がブロード化したことが、JSC 値が向上したことの主な原因であると考えられる。このように、本計算結果は実験結果の理論的解釈の一助となった点で意義がある。
Figure 3. Action spectra of ZnPQI cell (solid line) andZnPBQ cell (dashed line)

3 参考文献[1] H. Imahori, T. Umeyama, K. Kurotobi and Y.
Takano,2.2 理論計算
ポルフィリンの再安定化構造およびそのフロン
ティア軌道の電子構造について知見を得るために、“Self-Assembling

Porphyrins

andPhthalocyanines for Photoinduced Charge Separation

and Charge Transport”, Chem. Commun., 48,
4032-4045 (2012).

[2] A. Kira, Y. Matsubara, H. Iijima, T. Umeyama, Y.

Matano, S. Ito, M. Niemi, N. V. Tkachenko, H.

(a)

Lemmetyinen and H. Imahori, “Effects of π-Elongation
and Fused Position of Quinoxaline-Fused Porphyrins as

Sensitizers in Dye-Sensitized Solar Cells on Optical,

Electrochemical, and Photovoltaic Properties”, J. Phys.
Chem. C, 114, 11293-11304 (2010).

[3] H. Imahori, H. Iijima, H. Hayashi, Y. Toude, T.

(b)

Umeyama,Y.

Matano

and

S.

Ito,

“Bisquinoxaline-Fused Porphyrins for Dye-SensitizedSolar Cells”, ChemSusChem, 4, 797-805 (2011).

[4] H. Hayashi, A. S. Touchy, Y. Kinjo, Y. Toude, K.

Kurotobi, T. Umeyama, Y. Matano, and H. Imahori,
HOMO“Triarylamine-Substituted

LUMO

Figure 4. Sets of molecular orbital diagrams for (a) ZnPQI and(b) ZnPQ obtained by DFT calculations with
B3LYP/6-31G(d).

Imidazole-andQuinoxaline-Fused Push-Pull Porphyrin for High
Performance

Dye-Sensitized

ChemSusChem, 6, 508-517 (2013).

SolarCell”,

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 3スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

後周期遷移金属錯体触媒を用いる置換アセチレンの重合機構の解析
Jesus Rodriguez Castanon, 三田 文雄*
京都大学大学院工学研究科, *関西大学化学生命工学部

1 緒言

スーパーコンピューター上の Gaussian 09 プログ
ラムを活用し，higher layer (Pd, フェニル基を除く

置換ポリアセチレンは剛直な共役主鎖構造に基

dppf 基，CH≡CCH2OCONH，CH≡CCH2OCOO) にづく光・電気特性，液晶形成能，気体透過性，刺密度汎関数法 (DFT，B3LYP/LANL2DZ) を，lower

激応答性などの様々な有用な物性を示すことから，

layer に PM6 を用いた ONIOM により行った。

機能性高分子材料として活発に研究されている。
置換ポリアセチレンは対応する置換アセチレンモ

3 結果と考察

ノマーの遷移金属触媒重合により合成される。触媒の中心金属と配位子の設計により，重合活性，

Scheme 2 に，dppf，p-シアノフェニル基 (R)，

生成ポリマーの主鎖の幾何構造，分子量を制御可

CH3CN に配位された PdII 種 I への 1 番目と 2 番目能である。近年は，高い官能基耐性と，立体規則の置換アセチレンモノマーCH≡CH2OCOXC6H13

性の観点から，ロジウム，パラジウムをはじめと(X = NH, O) の配位–挿入における各中間体の自する後周期遷移金属錯体触媒が広く用いられてい由エネルギー G を示す。I は，AgOTf との反応にる。筆者らはこれまで置換アセチレン重合触媒のより，Pd 錯体から Br が脱離し，平面四角形の空構造を，
核磁気共鳴分光法，
単結晶X線構造解析，

配位座に CH3CN が配位して生成する。最初のス

元素分析，質量分析等により検討してきた。しかテップは，Pd に配位している CH3CN とモノマーしながら，触媒によってはこれらの分析が困難で，との交換である。モノマーの三重結合部が Pd 中

十分な知見を得ることができていなかった。
また，心に配位し，中間体 II が生成する。II は，モノマ重合反応の中間体の観測は極めて困難で，反応機ーのカルボニル基が垂直方向から Pd に配位した構の解析は不十分であった。本研究では，パラジ中間体 II’に転換する。次のステップは，Pd 中心

ウム錯体触媒による，カルバメートまたはカーボと R 基の間へのモノマーの 1,2-型挿入である。遷

ネート基含有アセチレンモノマーの重合 (Scheme

移状態はメタラシクロブタン TSI と考えられる。

1) における，配位–挿入反応の各種中間体のエネ

モノマー挿入により III を生成する反応は吸熱反

1)ルギーを計算し，反応機構について考察した。

応であるが，III に CH3CN が配位した中間体
III-CH3CN の生成は発熱反応である。CH3CN は中間体の安定化に大きく寄与していると考えられる。
III あるいは III-CH3CN が生成すると，カルボニル基は立体的な制約で，Pd に配位できない。すな
わち，大きな歪みのため III’は生成しない。
2 番目のモノマーの挿入は，Pd の空配位座へのScheme 1. Polymerization of acetylene monomer by[(dppf)PdBr(C6H4-p-CN)/AgOTf catalyst. Tf =CF3SO2.

2 計算手法計算には京都大学学術情報メディアセンターのモノマーの配位による IV の生成から始まる。1番目のモノマーの挿入と同様に，IV はモノマーのカルボニル基の垂直方向からの配位により安定化された IV’に異性化する。メタラシクロブタン遷
移状態 TSII を経由して V が生成する。IV’から V
への ΔG は，II’から III の生成，あるいは II’から
III-CH3CN 生成の ΔG よりも負にずっと大きい。

4 − Academic Center for Computing and Media Studies, Kyoto University

Scheme 2. Reaction pathway for 1,2-insertion of CH≡CH2OCOXC6H13 (X = NH, O) into[(dppf)PdII(NCCH3)(C6H4-p-CN)]+; dppf and C6H4-p-CN are shown as an arc and R, respectively. The values below thecompound numbers represent their relative G values (kJ/mol), based on I. The counteranion (TfO–) is omitted forclarity.
このことより，TSII 生成の活性化エネルギーは，ェニルアセチレンモノマーの重合機構についても

TSI のそれよりも低いと推測される。これは共役

DFT 計算を行い，この場合は本系とは逆に 2,1-型が長くなる効果と考えられ，カルバメートとカー挿入が有利であることを明らかにしている 2)。こボネートモノマーのいずれも同様の結果であった。

れは，挿入したモノマーのフェニル基の π-電子がアセチレンモノマーのパラジウム錯体への配位ロジウムに配位安定化する効果が，立体障害によ−挿入機構として，上述の 1,2-型に加えて，2,1-型がある (Scheme 3)。2,1-型についても同様に各中間体のエネルギーを計算したところ，2,1-型挿入る効果よりも大きいためである。

4 結言

は，1,2-型挿入に比べて配位−挿入時の立体障害がより大きく，エネルギー的に不利であることが確

本研究により，これまで明らかにされていなか認された。筆者らは，ロジウム錯体触媒によるフった，パラジウム錯体触媒によるアセチレン重合の反応機構に関する有用な知見が得られた。本研究成果の活用により，新しい触媒の分子設計と，重合挙動の解析における一層の進展が期待される。
文献	 1) J. Rodriguez Castanon, K. Kuwata, M. Shiotsuki, F.
Sanda, Chem. Eur. J. 18, 14085–14093 (2012).

Scheme 3. 1,2- and 2,1-insertions.

2) N. Onishi, M. Shiotsuki, T. Masuda, N. Sano, F.
Sanda, Organometallics 32, 846–853 (2013).

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 5スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

メニスカスを有する液膜内マランゴニ対流の流動特性
山本卓也
大阪大学大学院基礎工学研究科物質創成専攻

1 緒言

体領域の非定常熱伝導は近似式を用いて境界条件として与えた。
本研究で用いた近似を以下に示す。

気液界面で温度差や濃度差が発生すると局所の(1) 重力の影響は無視できる。

表面張力が変化し、界面上で表面張力差が発生す(2) 液膜から周囲気体への熱移動は無い。

る。表面張力差によって駆動する流れはマランゴ(3) 表面張力を除く物性値は温度に依存しない。

ニ対流として知られ、温度差が発生すると一般に(4) 自由表面は流れによって変形しない。

低温部に向かって流れる。マランゴニ対流は浮力	  解析領域の概略図を Fig. 1 に示す。無重力、固対流と共存し、浮力による対流の影響の方が大き

液界面が平滑、固液気体の接触点ではピン止めさい場合が多い。そこで、浮力対流のほとんど発生れているとし、Young-Laplace 式を用い自由表面形

しない微小重力環境においてのマランゴニ対流の状を推定すると、自由表面形状は球の一部と同じ研究が数多く報告されている[1]。

形状になる。そして、液膜を作製するのに用いる

宇宙飛行士の D. Pettit 博士は 2003 年に国際宇宙液体積（平坦な液膜の液体積に対する体積比

ス テ ー シ ョ ン (ISS) の 中 で ’Saturday Morning(V/V0)）によって自由表面形状を決めた[4]。

Science’と呼ばれる一連の簡単な実験を行った[2]。

	  支配方程式は Navier-Stokes 式、連続式、エネル

その実験の一つに温度差マランゴニ対流の観察がギー式を用いた。境界条件として上下の自由表面ある。ステンレスのワイヤを用いて直径 50 mm 程度の円形枠を作製しその中に水の液膜を張り、スには以下の温度差マランゴニ効果の条件を与える。
(1)µ∇ nu = −σ T ∇ tT

テンレス枠に加熱されたハンダごてを近づけるこ

ここで、u は速度ベクトル、T は温度、µは粘度、とにより液膜の一部を加熱し、マランゴニ対流を

σT は表面張力の温度勾配係数であり、添字 t, n は発生させた。その結果、熱源に近づく流れが発生

それぞれ自由表面に対して接線、法線方向を表す。

した(Outward flow と呼ぶ)。熱源から離れる流れ	  また、自由表面上では断熱であるとする。固体(Inward flow)が予測されたが、それとは逆に向かう

壁では、滑りなし条件と以下の近似化された非定流れが観測された。この実験について現在様々な常熱伝導を表す条件を用いる[4]。

議論がされてきたが[3]、宇宙実験での流れのメカ
ニズムについては判っていない。

q
T = T0 + ⋅A( 6αt − r θ )2

2k 6α t(2)

そこで、本研究では数値シミュレーションを用

ここで、T0 は初期温度、q/A は熱流束、αは熱拡散

いて、宇宙実験で Outward flow が発生したメカニ

率、t は時間、r は液膜半径、θは加熱点からの角ズムを解明することを目的とする。特に流れ方向

度、k は熱伝導率であり、熱流速 q/A の値は実験[2]に及ぼす液膜形状に着目し、様々な液膜の形状をより熱伝導と輻射の影響を 考慮した概算値

用いて非定常三次元数値計算を行った。

1.0×104 W/m2 を用いた。

2 数値解析手法

	  数値解析にはフリーソースである OpenFOAM
を用い、空間に対して二次精度中心差分、時間に対して陰的オイラー法、速度場と圧力場の連成に本研究では液相領域のみの流体計算を行い、固6 − Academic Center for Computing and Media Studies, Kyoto University

は PISO アルゴリズムを用いた。

3 結果と考察“Let’s enjoy the microgravity science with Don Pettit
convection experiment in liquid thin film with thin

3.1 平坦な液膜(体積比 V/V0 = 1)

metallic wire ring”, J. Jpn. Soc. Microgravity Appl., 23,

平坦な液膜の場合の加熱 3 秒後の流跡線を Fig.
2 に示す。上から見た図(a)、加熱部を横から見た図(b)をそれぞれ示している。Figure 2 より液膜内に上下に対称な渦対が発生していることが判る。
しかし、加熱部近傍の渦対の他に流れは発生しなかった。

157-160 (2006).
[4] T. Yamamoto, Y. Takagi, Y. Okano, and S. Dost,“Numerical investigation for the effect of the liquidfilm volume on thermocapillary flow direction in a thincircular liquid film”, Phys. Fluids, 25, 082108 (2013).

謝辞

3.2 歪んでいる液膜(体積比 V/V0 ≠ 1)	  本共同研究制度（若手研究者奨励枠）を活用さ歪んでいる液膜の場合には平坦な液膜とは異なせて頂きましたことを、この場を借りて深謝致し

り、加熱近傍での渦対の他に流れが発生した。

ます。また、本研究の一部は科学研究費補助金・Figure 3 に液膜が膨らみ形状(V/V0 > 1)、凹み形状

基盤 B（課題番号 23360343, 25289087）の援助を(V/V0 < 1)の場合の加熱 3 秒後の流跡線を示す。ど受けた。

ちらの場合においても加熱部近傍で流跡線が上下(a)に動いている、つまり、平坦な場合と同様に渦対が発生していたが、液膜が膨らんでいる場合には(b)

加熱部に向かう液膜表面に平行な二次元流が発生した(outward flow)。その一方、液膜が凹んでいる
場合には、加熱部から液膜中心に向かう流れが発

Fig. 1 Schematics of numerical domain: (a) convex

生した(inward flow)。

water film (V/V0 > 1), and (b) concave water film (V/V0
< 1).

4 結言

(a)

Heat (b)Heat

円形液膜内マランゴニ対流の流れの方向性を解明することを目的として液膜形状に着目して数値解析を行った。数値解析によって液膜が凹んでいる場合(V/V0 < 1)、実験と逆向きの inward flow が発生し、液膜が膨らんでいる場合(V/V0 > 1)、実験と同方向の outward flow が発生した。これより、
D. Pettit 博士が行った実験では液膜が膨らんでいFig. 2 Path lines in a flat water film at 3 s: (a) top view,

たことにより加熱部に向かう流れが発達し、予測

and (b) side view.

と逆向きの流れが発生したと考えられる。

(a)

Heat

(b)

Heat

引用文献
[1] M. Lappa, Thermal Convection: Pattern,
Evolution and Stability, Elsevier Ltd., Oxford (2004).
[2] D. Pettit, Saturday Morning Science Videos,
2003, see http://mix.msfc.nasa.gov/IMAGES/QTVR/
0601211.mov.
[3] H. Kawamura, K. Iwamoto, and H. Harashima,

Fig. 3 Path lines in a: (a) convex water film (V/V0 > 1),and (b) concave water film (V/V0 < 1) at 3 s.

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 7スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告␩ ⍬␨ ␯ ␀ ‫⎯ذ‬೑Α଑࠵ՃՕஞԝ⎯␉ ⏿ ␄ ␫ ⍬⏫ ۚਏॵ෱ խງ∗λֹ੒Ӂ੒Ӂώ࢙Ӂ٥֋ѕ

1

⎰⎙ ⏂⎬

␉ ⏿ ␄ ␫ ⍬⏫ ⏓ ๆ⎅⎠⏧ ⏦ ⏵ ѽॾ3⏧ ⏦⏵ ⏞ ␄ ␥ ⏫ ⏻ ⎬ ⎰ਪॅ⎯೑Α଑࠵ՃՕஞ⎍ඝ⏂‫ܛ‬⎔ ⎔ ⎨⎰ろ ߹‫ ␩⎬ݘ‬⍬␨ ␯ ␀ ൘ଟߧ⏓ ๆ⎅⎧ ␉⏿ ␄ ␫ ⍬⎿ ⏍⎧ ⎋ ⏋ [1], ߹‫ ⎬ݘ‬ろ ⎞ ⏍⏊ ⎯೑Α଑࠵ՃՕஞ⎰⏧ ⏦⏫ ш⏓ ۱⎤ ⎠ٌї⏓ ߟ⎚ゎ ඕ٥֋⎨⎰ろ ੒‫ޔظؠ‬Ն⏓ ๆ⎅⏵ຢӁ‫⎫์࡚ ⎧⎅⎋⎬ذ‬ฉӍ⏓ Ж⎙ ⎧ ⎅⏌ ߏ⎍โゖ ⎫٥֋⎧ ෭๚ଯ⎬٠ࡥ⎖ ⏍⎠␙⏞␯ ⏧ ␨ ࠃਉ࠵Ճ⎨ 14 ࠵Ճ⎿ ⎨⎌⏊ ෗⏊ ⎌⎬⎫⎤ ⎧ ⎅⏌ ゎ ⎘ ⎌⎘ ⎫ ⎍⏊ ろ ⎪ ⎯⏉ ⎇ ⎫ ೑Α ⎯೑Α଑࠵ՃՕஞ⏓ ๆ⎅⎧ ѽॾ⎚⏌ ゎ ␙⏞ ␯ ⏧ ␨ ࠃਉ⏓଑࠵ՃՕஞ⎍⏧ ⏦⏵ ⎯ஂ‫ظ‬ຟ⏓ ຜ⎐ ‫⎌⎯ ⏌⎚ߐפ‬ろ ⎿ ⎠ ろਪॅ⎯೑Α଑࠵ՃՕஞ⎯ય⎌⏊ ຜ⎐ ‫ ⏌⎚ߐפ‬೑Α଑࠵Ճ଑զ⎚⏌ ⎠ ⏂⎬ ๆ⎅⎠ ␙⏞ ␯ ⏧ ␨ ઐ෣⎰ z = r − 1 ⎌⎥ż > 0 ⎨⎃⏌ ゎ ␉ ⏿ ␄ ␫ ⍬⏫ ш⏓ ۱⎇ ‫ॅ␧␒ ␯⏱ ␯⏞⎯ݘ‬⎫⎅ゎ ඕ٥֋⎨⎰⎔ ⏍⏊ ⎯դ෻⎬୿⎉ ⏌ ⎠⏂ろ ␩ ⍬␨ ␯ ␀(T ) ⎬ ԟ⎘ ⎧ ⎰ろ ඕۣ⎨⎰ T = 106 ⏓ Ά⎇ ゎ ⎥⎿ ⏋ ろ 106Ѿ⎯߇۱⎨ু⎧ ⎯␦ ␯ ⏫ ۚਏ⏓ ‫ى‬଑⎚⏌ ゎ T = 2 × 106 ⎩‫ ⎬ذ‬ඝ⏂‫ ⎠⏍ ⎿ܛ‬λ଑࠵Ճ⎿ ⎨⎯೑Α଑࠵ՃՕஞ⏓ ੒‫ؠ‬⎘ ⎧ ␉⏿ ␄ ␫ ⍬⏫ ⏓ ۚढ़⎘ ⎠ࣟ‫ ⏃܄‬ろ ଑ड़ଯ⎫ٌї⎰ബ⏐‫ޔظ‬Ն⏓ ๆ⎅⎧ ෭๚ଯ⎬٠ࡥ⎘ , ⎞⏍⏊ ⎯ԝ⎯ԟ‫ ⏓ؙ‬೾‫ށ‬⏊ ⎫ ⎅⎔ ⎩ ⎰ҷ௨‫ ⏌ ⎃⎨⏀݌‬ゎ ⏧ ⏦⏵ Օஞ⎌⏊ ண⏊ ⏍⎠␉ ⏿ ␄ ␫ ⍬⏫ [2] ⎯Ԙ୅⎨۪‫ ⏌⎚ݺ‬.

࠵ՃՕஞԝ␉⏿ ␄ ␫ ⍬⏫ ⎬਻⎘ ⎧ ೾‫ ␫ ␄ ⏿␉ށ‬⍬⏫ ѽॾ⏓Օஞ⏓ ষ⎴ࡥ⎚⎠ ⏂⎯൘൚⎫ ⎪ ⎰⎃ ⎿ ⏋ ຜ⎐ ગ⏊ ⏍⎧ ⎅۱⎇ ߏ⎨⏧ ⏦⏵ ⎯ஂ‫ظ‬ຟ⎩ ೾‫ ␫ ␄ ⏿ ␉ށ‬⍬⏫ ⎬ ⎋ ⎅⎧ ࡚2೑Α଑࠵ՃՕஞԝ⎯␉⏿ ␄ ␫ ⍬⏫ ш์⎫ ຟ⎩ ⎯ԟ‫ ⏓ؙ‬෗⏊ ⎌⎬ ⎚⏌ ゎΙч⎯ࠡਞ⎎ ⎬ࡓ⎤ ⎧ , ⏧ ⏦⏵Օஞ⎯೑Ց਒সΫѬଟ⏓

1000

100

଑࠵ՃՕஞ⏓ า⏂ਪॅ‫ڃ‬ֆ⏂⏌ ゎ (2) ଱୻⎫ ␙⏞␯ ⏧ ␨ ઐ෣⏓ ঎଑⎘ ろ า⏂ֆ⏂⎠ ೑Α଑࠵ՃՕஞ⎯␙⏞ ␯ ⏧ ␨ ઐ

100

distribution␫ ⍬⏫ ⎩ ⎘ ⎧ ॅ๲␞␃␧ш⎚⏌ : (1) Ά⎇ ⏳⏵ ␂ ␜ ⎯೑Αdistribution‫ ⎬ذ‬ඝ⏂‫ ⎠⏍ ⎿ܛ‬೑Α଑࠵ՃՕஞ⏓ ␊ ⍬␅ ⎩ ⎘ ⎠ ␉ ⏿ ␄1010෣ࣘ⎬ ⎋ ⎒ ⏌ Λટ⏓ Ւм⎘ ろ ⎞ ⏍⏊ ⎯೑Α଑࠵ՃՕஞ⎬਻⎘ ⎧ ෕ঽ⏓ ೒⎒ ⎧ ⎋⎐ ゎ (3) ␙⏞␯ ⏧ ␨ ઐ෣⎬⎋⎒ ⏌ ⏧1

1
1⏦⏵ ⏞ ␄ ␥ ⏫ ⏻ ࣘ⎯଱୻⎫ λ୅⏓ ࡾՃِࣥ⎩ ⎘ ⎧ ろ ݁⎴10out degree

100

1

10
in degree100␙⏞ ␯ ⏧ ␨ ઐ෣⏓ ૳Ѭ⎚⏌ ⎿ ⎨⏧ ⏦⏵ Օஞ⎯ߗԝ౎ଽ⏓‫ ⏌⎚ޔظ‬ゎ (4) ⏧ ⏦⏵Օஞ⎍␙⏞␯ ⏧ ␨ ઐ෣⏓ ૳Ѭ⎚⏌ ‫ݘ‬य 1: ೑Α଑࠵ՃՕஞ⎯ࡥߘॅ (out-degree) ഇ೘ (ܳ). ௡⎬ろ ⏧⏦⏵Օஞ⎩ ␢⍬⏫ ␦ ⏿ ␅ ֜๸⎍݂⏃ ‫⎅פ‬೑Α଑࠵Ճߘॅ (in-degree) ഇ೘ (ϓ).

Օஞ⏓ ষ⎷ゎ ⎞⎯೑Α଑࠵ՃՕஞ⎯෕ঽ⏓ Ւм⎚⏌ (⎔ ⎔⎨⎰ i ⎨⎃ ⏌ ⎩ ⎚⏌ )ゎ (5) ⎖ ⏊ ⎬ ⏃ ⎇ λ୞␙⏞ ␯ ⏧ ␨ ઐ෣⏓ ૳Ѭ⎚⏌ ⎿ ⎨ߗԝ౎ଽ⏓ ਞ⎒ ⏌ ゎ (6) ૳Ѭ⎘ ⎠⏊ ろ ݁⎴⏧ ⏦⏵ Օஞ⎩ ␢⍬⏫ ␦ ⏿ ␅ ֜๸⎍݂⏃ ‫⎅פ‬೑Α଑࠵ՃՕஞ⏓ ষ⎷ゎ ⎞⎯೑Α଑࠵ՃՕஞ⎯෕ঽ⏓ Ւм⎚⏌ (⎔ ⎔⎨⎰ j ⎨⎃⏌ ⎩ ⎚⏌ )ゎ (7) i ⎩ ⎅⎇ ೑Α଑࠵ՃՕஞ⎩ j ⎩⎅⎇ ೑Α଑࠵ՃՕஞ⎯ԝ⎬ ⎰␦ ␯ ⏫ ⎍ਣ‫⎞ ⎧ ⎘ ⎩ ⏌⎚ݚ‬⏍⏊ ⎯ԝ⎬ ย‫ ⏫ ␯ ␦ڽ‬i → j ⏓ ิ⎉ ⏌ ゎ (8) (3) ⎌⏊ (7)⏓ T (≫ 1) Ѿ،⏋ റ⎚ߏ⎨೑Α଑࠵ՃՕஞ⏓ ␊ ⍬␅ ⎩ ⎘⎠ ␉ ⏿ ␄ ␫ ⍬⏫ ⏓ ۚढ़⎚⏌ ゎ∗ ඕ٥֋⎰ୱඅ੒Ӂ AIMR ⎯ࢡ຦ӻ޺⎩ ⎯֩ஔ٥֋⎬ Է⎦⎐ . Email:yoshi.saiki@r.hit-u.ac.jp

8 − Academic Center for Computing and Media Studies, Kyoto Universityय 1(ܳ) ⎰ࣘՒ⎯൘൚⎨ண⏊ ⏍⎠␉ ⏿ ␄ ␫ ⍬⏫ ⎬⎋⎒ ⏌ࡥߘॅഇ೘⏓ ಺⎘ ⎧ ⎋⏋ ろ ഇ೘⎰␕⏩ଯ⎨⎃⏌ ⎔ ⎩ ⎍٪⎧ࠟ⏍⏌ ゎ ⎔ ⏍⎰ろ ⏧ ⏦⏵Օஞ⎌⏊ ண⏊ ⏍⎠␉ ⏿ ␄ ␫ ⍬⏫ ⎍೾‫ ␫ ␄ ⏿ ␉ށ‬⍬⏫ ⎯஦૎⎯λ⎥⎨⎃ ⏌ ”⏵ ⏭ ⍬␧␑ ␦ ⍬ ”
ड़⏓ ย⎘ ⎧ ⎅⏌ ⎔ ⎩ ⏓ ߟܱ⎘ ⎧ ⎅⏌ ゎ λ൘ろ य 1(ϓ) ⎰௡ߘॅഇ೘⏓ ಺⎘ ⎧ ⎅⏌ ⎍ろ ࡥߘॅ⎯ࣟ‫⎰ ⎩܄‬Ϊ⎫⏋ ろ ௡ߘॅഇ೘⎰␕⏩ଯ⎩ ⎰⎫⎤ ⎧ ⎅⎫ ⎅ゎ ⎚⎫ ⏐⎢ , ‫׎‬અ⎬੒⎎⎅௡ߘॅ⏓ ߖ⎤ ⎠ ೑Α଑࠵ՃՕஞ⎍ਣ‫ ⏓ ⎩ ⎔⎅ ⎫ ⎘ݚ‬Τය⎚⏌ ゎ
ߘ⎬ろ ೾‫ ␫ ␄ ⏿␉ށ‬⍬⏫ ⎯ഇฅ⎨ຜ⎐ ‫ ⏌⎅ ⎧⏍⏐ޡ‬યअड़⎯Ԙ୅⎌⏊ ⏧ ⏦⏵ ⎯ஂ‫ظ‬ຟ⏓ ຜ⎐ ‫ ␊ ⎫์࡚ ⏌⎚ߐפ‬⍬␅ (೑Α଑࠵ՃՕஞ) ⏓ ஦૎⎦⎒ ⏌ ⎔ ⎩ ⏓ ۱⎇ ゎ યअड़⎩

50

0.4

-50

0.2

⎍ , ඕൎ‫ॅߘ⎰⎨܍‬યअड़⎯⏀⏓ Ά⎇ ゎ
130

λ1upo-λ1chaosॅયअड़⏅టѻયअड़⎫ ⎪ โゖ ⎫ ࠦຳ⎍ગ⏊ ⏍⎧ ⎅⏌ [2]<z3>upo-<z3>chaos⎰⎪ ⎯␊ ⍬␅ ⎍␉⏿ ␄ ␫ ⍬⏫ ⎯์⎌⏓ ஦૎⎦⎒⏌ ຟ⎨ろ ߘ

-150
-250

0
-0.2

0.5

-350

10
0

y

<Z>upo-<Z>chaos

20

0
-10

20

40
60
out degree

80

100

-30
-20

0

100 200 300 400 500 600
out degreeय 3: ⏧ ⏦⏵ Օஞ⎬ Е⎤ ⎠ z 3 ⎯ߗԝച‫ ⎩ה‬ү೑Α଑࠵Ճ-20
0

100 200 300 400 500 600
out degree

-0.5
-1

-0.4
0

Օஞ⎬ Е⎤ ⎠ z 3 ⎯ߗԝച‫�( ܴ⎯ה‬z 3 �chaos − �z 3 �upo ) ⎩
-10

0
x

10

20

೑Α଑࠵ՃՕஞ⎯ࡥߘॅ⎯ԟ‫)ܳ( ؙ‬ゎ ⏧ ⏦⏵ Օஞ⎯␦ ⏞␓␊ ␑ ‫ ॅޱ‬λ1chaos ⎩ ү೑Α଑࠵ՃՕஞ⎯␦ ⏞ ␓␊ ␑ ‫ॅޱ‬

य 2: ⏧ ⏦⏵Օஞ⎬Е⎤ ⎠ z ⎯ߗԝച‫ ⎩ה‬ү೑Α଑࠵ՃՕ

λ1upo ⎯ܴ (λ1chaos − λ1upo ) ⎩ ࡥߘॅ⎯ԟ‫( ؙ‬ϓ)ゎ

ஞ⎬Е⎤ ⎠ z ⎯ߗԝച‫�( ܴ⎯ה‬z�chaos − �z�upo ) ⎩ ೑Α଑࠵ՃՕஞ⎯ࡥߘॅ⎯ԟ‫)ܳ( ؙ‬ゎ ܳय⎨Ά⎤ ⎠ ೑Α଑࠵ՃՕஞ⎯ય⎨݂⏃ ੒⎎ ⎫ યअड़⏓ ߖ⎤ ⎠೑Α଑࠵ՃՕஞ (࠵Ճ 9.7) (ϓ)ゎ⎘ ⎧ ⏃ ࡥߘॅ⎩ ⎯৷ԟ⎰ֲ⎅⎩ ‫ ⏌ ⎉ٿ‬ゎ4⎿⎩⏂

य 2(ܳ) ⎰ബॅ z ⎯ߗԝച‫⏵⏦ ⏧ ⏓ה‬Օஞ⎬Е⎤ ⎧ ۱⎤ඕ٥֋⎨⎰ろ ⏵ ␄ ␨ ␯ ⏴⏞ ␄ ␥ ⏫ ⏻ ⎬ ඝ⏂‫ ⎠⏍ ⎿ܛ‬೑⎠ખ⎩ ೑Α଑࠵ՃՕஞ⎬Е⎤ ⎧ ۱⎤ ⎠ખ⎩ ⎯⎛⏍⎯ࡥߘॅΑ଑࠵ՃՕஞ⎯ԝ⏓ ⏧ ⏦⏵ Օஞ⎍೑Ց਒⎬ সΫ⎚⏌ ‫࣐ټ‬Μਣड़⏓ ߟ⎘ ⎠ゎ ࡥߘॅ⎍੒⎎ ⎅⎼⎪ ろ ຑࠋ⎯⎛⏍⎰ࢡ⎖⏓ ろ ␉⏿ ␄ ␫ ⍬⏫ ⎯Ԙ୅⎌⏊ ॅ๲␞␃␧ш⎚⏌ ൘൚⏓ ଚΕ⎐ ⎫ ⏌ ‫ ⏌⎌⏐⎍ ⎩ ⎔ ⏌⎃⎍ڽؚ‬ゎ ࡓ⎤ ⎧ ろ z ⎯ച‫ה‬ખ⎬ԟ⎘ ⎠ゎ ⎖ ⏊ ⎬ろ ۚढ़⎖ ⏍⎠␉ ⏿ ␄ ␫ ⍬⏫ ⎌⏊ ண⏊ ⏍⏌ ஦૎⎘ ⎧ ⎰ろ ࡥߘॅ⎍੒⎎ ⎅೑Α଑࠵ՃՕஞ⎼⎪ ⏧ ⏦⏵ ஂ‫ظ‬ຟ⏓ ๆ⎅⎧ ろ ⏧ ⏦⏵⎯ѽॾ⏓ ۱⎤ ⎠ゎ ⎩ ⎐ ⎬ろ ␉ ⏿ ␄ ␫ ⍬⏓ ຜ⎐ ‫ ⏌⎚ߐפ‬೑Α଑࠵ՃՕஞ⎬ ⎫ ⏌ ‫⎇⎅ ⎩ ⏌ ⎃ ⎬ڽؚ‬⏫ ⎬ ⎋ ⎒ ⏌ યअड़⎩ ろ ⏧ ⏦⏵ ⎯ஂ‫ظ‬ຟ⏓ ຜ⎐ ‫ ⏌⎚ߐפ‬೑⎔ ⎩ ⎨⎃⏌ ゎ य 2(ϓ) ⎰ろ ⎔ ⎔ ⎨Ά⎤ ⎠ય⎨ҙߘॅ⎍݂⏃Α଑࠵ՃՕஞ⎩ ⎯ԟ‫ ⏓ؙ‬૛⎺⎠ゎ ٌї⎩ ⎘ ⎧ ろ ࡥߘॅયअ੒⎎ ⎅೑Α଑࠵ՃՕஞ⎯Օஞ⎯ x-y ച෣⎹⎯ࠄϳ⎨⎃ ⏌ ゎड़⎯‫⎅ۿ‬೑Α଑࠵ՃՕஞ⎍ろ ⏧ ⏦⏵ ⎯ஂ‫ ⏓ظ‬ຜ⎐ ‫⎚ߐפ‬⎔ ⎯೑Α଑࠵ՃՕஞ⎰ಒ࣢⎬ ຜ⎐ ⏧ ⏦⏵ ⎯ஂ‫⎚ߐפ ⏓ظ‬⏌ ‫ٌ ⎇⎅ ⎩ ⏌ ⎃ ⎬ڽؚ‬ї⏓ ண⎠ ゎ⏌ Օஞ⎩ ⎫ ⎤ ⎧ ⎅⏌ . ⎫ ⎋ , Ιࣘ⎯ѽॾ⎰␙⏞ ␯ ⏧ ␨ ࠃਉ࠵Ճ 12 ⎿ ⎨⎯೑Α଑࠵ՃՕஞ⏓ ๆ⎅⎠ ٌї⎨⎃ ⏌ .
ߘ⎬ろ z ⎯ച‫ה‬ખ⎡⎒ ⎨⎫ ⎐ ろ ਩⎯⏃ ⎇ ࢢ⎘ ‫ظஂ⎯ߘۿ‬ඕൎ‫⎰⎨܍‬ろ ߘॅયअड़⎯⏀⎬ࢺ୅⏓ ୻⎧ ⎠⎍ろ ⎞⎯਩⎯યअड़⎍⏧ ⏦⏵ ⎯⎪ ⎯ड़߸⎩ ٌ⎴⎥⎅⎧ ⎅⏌ ⎌⎩ ⎅⎇Ѩ੕⏃ ಒ࣢⎬ ࡚์⎨⎃ ⏌ ゎຟ⎬ ԟ⎘ ⎧ ⏃ ஔโ⎯ٌї⎍ண⏊ ⏍⏌ ⎌⏓ ҷ⎌⏂⏌ ゎ Ιч⎰␙⏞␯ ⏧ ␨ ࠃਉ࠵Ճ 14 ⎿ ⎨⎯೑Α଑࠵ՃՕஞ⏓ ๆ⎅⎠ٌї⎨⎃⏌ . य 3(ܳ) ⎰ബॅ z ⎯ 3 ߘ␞⍬␝ ␯ ␄ ⏓ ⏧ ⏦⏵Օஞ⎬ Е⎤ ⎧ ߗԝച‫ ⎠ ⎘ה‬ખ⎩ ೑Α଑࠵ՃՕஞ⎬ Е⎤ ⎧ߗԝച‫ ⎠ ⎘ה‬ખ⎩ ⎯⎛⏍⎯ࡥߘॅΜਣड़⏓ ߟ⎘ ⎧ ⎅⏌ ゎ z⎯ച‫ה‬ખ⎯ࣟ‫ ⎩܄‬ஔโ⎬ろ z ⎯ 3 ߘ␞⍬␝ ␯ ␄ ⎬ԟ⎘ ⎧ ⏃ࡥߘॅ⎍੒⎎ ⎅೑Α଑࠵ՃՕஞ⎍⏧ ⏦⏵ ⎯ஂ‫ ⏓ظ‬ຜ⎐ ‫פ‬ߐ⎚⏌ ೑Α଑࠵ՃՕஞ⎨⎃⏌ ⎩ ‫ ⏌ ⎉ٿ‬ゎ य 3(ϓ) ⎰␦ ⏞␓‫۪ދ‬ഓ٤
[1] E. Ott, Chaos in Dynamical Systems, Cambridge
Univ. Press,Cambridge, 1993.
[2] ਊ୉૥Տ, ‫ܞ‬ฅՏฯ, ೾‫ ␫ ␄ ⏿ ␉ށ‬⍬⏫ -Է৐⎌⏊ Ыๆ⎿ ⎨ , ‫פ‬੐ѕӁࠉ, 2010.

␊ ␑ ‫ ⏵⏦ ⏧ ⏓ॅޱ‬Օஞ⎬ Е⎤ ⎧ ച‫ ⎠ ⎘ה‬ખ⎩ ೑Α଑࠵ՃՕஞ⎬ Е⎤ ⎧ ച‫ ⎠ ⎘ה‬ખ⎩ ⎯⎛⏍⎯ࡥߘॅΜਣड़⏓ ߟ⎘⎧ ⎅⏌ ゎ ࣘՒ⎬ߟ⎘ ⎠਩⎯ஂ‫ظ‬ຟ⎩ ஔโ⎬ろ ࡥߘॅ⎍੒⎎⎅೑Α଑࠵ՃՕஞ⎍⏧ ⏦⏵ ⎯␦ ⏞ ␓␊ ␑ ‫ ⏓ॅޱ‬ຜ⎐ ‫ߐפ‬⎚⏌ ‫ ⏌⎃⎬ڽؚ‬ゎ z ⎯ച‫ה‬ખຟ⎯⏉ ⎇ ⎫ ଊߘ⎯ஂ‫ظ‬ຟ⎡⎒⎨⎫ ⎐ ろ z 3 ⎯⏉ ⎇ ⎫ ⏉ ⏋ ‫ظஂ ⎫ߘۿ‬ຟ⏅೑Α଑‫ ⎬ॅޱ‬ԟ
京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 9スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

フラグメント分子軌道法を用いたタンパク質中電子移動の高精度経路解析法の開発鬼頭(西岡)宏任 a), 安藤耕司
京都大学大学院理学研究科 ( a 現所属: 名古屋大学大学院理学研究科)

1 INTRODUCTION
蛋白質中でおこる電子移動反応は、光合成や内路の結果をどの程度改善するか検証を行った。

2 MODEL & METHODS

呼吸などの生体エネルギー変換に重要な役割を果たしている。光合成反応中心蛋白質や、チトクロ

本研究では、トリプトファン(Trp)とプロリン

ムやアズリンなどの金属蛋白質に対する分光研究(Pro)からなる Trp-(Pro)3-Trp ペプチド(図 1)をモから、蛋白質環境下のトンネル効果によって、電デル系として採用し、Trp 間のホール(正孔)の長子は 10〜30Åもの長距離を移動していることが分距離移動に対して我々の理論を適用した。FMO

かっている[1]。そこで我々は、蛋白質環境が電子計算では、炭素原子間の共有結合を切断(図中の破移動を媒介する詳細なメカニズムを明らかにする

線)した五つのフラグメントを用いた。 (Pro)3 ブ

ために、大規模系のための量子化学計算手法であリッジが電子移動を媒介する強度は、一般化マリるフラグメント分子軌道(FMO)法[2]を利用した

ケン・ハッシュ(GMH)法または架橋グリーン関数

電子移動経路解析法の開発を行った。

(bridge GF)法で求めた電子的相互作用(TDA)から

FMO 法は、蛋白質などの巨大分子を小さなフ評価した[6,7]。反応速度はこの TDA の二乗に比例

ラグメントに分割し、フラグメント単量体とフラする[1]。 (Pro)3 ブリッジ中の電子移動経路は、トグメント二量体を他のフラグメントの静電ポテン

ンネル電流解析法から解析した[6,7]。

シャル下で解いたエネルギーから、分子全体のエネルギーを計算する[2]。このオリジナルの FMO
法(FMO2)に対して、計算コストは増えるが、フPro
Donor FMO

ラグメント三量体補正まで考慮した FMO3 法を

Trp

Pro

用いると、全系エネルギーに対する計算精度は大きく改善される[3]。FMO 法には全系の波動関数

Trp

Pro

を求めることが出来ないという欠点があったが、常行らのフラグメント分子軌道線形結合(FMO-LCMO)法の開発によって解決された[4]。

Acceptor FMO

図 1 Trp-(Pro)3-Trp ペプチド

同グループは、FMO3 結果に基づく LCMO 計算

GMH 法には、クープマンスの定理から、全系(FMO3-LCMO)から、三量体補正が全系の波動関の荷電状態を中性として得られる最高占有軌道数の精度も大きく改善することを示した[5]。

(HOMO)と、その一つ下のエネルギーの軌道

これまで我々は、FMO2-LCMO 法を適用して(HOMO-1) を使用した[7]。Bridge GF 法とトン得た分子軌道関数を基に、電子移動の経路を解析

ネル電流解析法では、図１中に示した 左端の Trpする手法を開発し[6,7]、光合成系への適用を行っフラグメント単量体の HOMO にあった電子が、てきた[7,8]。今回は、解析手法を FMO3-LCMO

右端の Trp フラグメントの HOMO へ移動するとを使ったものに拡張し、三量体補正が電子移動経して計算を行った[6,7]。
（従って、ホールは右端

10 − Academic Center for Computing and Media Studies, Kyoto UniversityTrp から左端 Trp に移動するとみなす。
）

現実的な蛋白質系での経路解析を行う予定である。

本研究では全て 6-31G(d)基底関数を用いた。今(a)RHF

回のように、いわゆるダブルゼータ以上の基底関

0.16

2.37

数を使用し、共有結合の切断を伴う FMO-LCMO

2.6 0.38 0.22

計算を行う場合、切断部分の基底関数が両側のフラグメントでダブルカウントされることから生じる線形従属性を取り除く必要がある[5]。そこで我々は、フラグメント分子軌道線形結合を (価電

0.85

1.37

(b)FMO3-LC(VC)MO子+ 内殻電子) 軌道空間に制限する

0.17

2.29

FMO-LC(VC)MO 計算を行うことで、基底関数の2.61 0.29

次元を減らし、線形従属性の問題を回避した。この手続きは、電子移動解析の計算コストの大幅な削減にもなっている [6,7]。

0.17
0.15

0.92

1.40

(c)FMO2-LC(VC)MO0.41

3 RESULTS

0.96

0.98

通常の全電子 Hartree-Fock(RHF)計算から求1.43

めた全エネルギーに対する FMO 計算の誤差は、FMO2: 2.26 milli-Hartree (1,42 kcal/mol)
FMO3: 0.0349milli-Hartree (0.0219 kcal/mol)となり、三体補正を考慮することによって大きく

0.37

図 2	
  フラグメント間トンネル電流マップ

5 ACKNOWLDGEMENTS

改善される。
次に TDA の計算結果を下記の表に示す。参照と本研究は、スーパーコンピュータ共同研究制度する RHF では、全波動関数を(価電子+内殻電子)(若手奨励枠)と科研費新学術領域研究
「π空間」
(領

軌道空間に射影した[7]。表から分かるように、TDA

域番号 2007)による。

の計算結果は計算方法に大きく依存せず参照値を[1] 垣谷俊昭(1998)『光・物質・生命と反応(下)』良く再現している。

20 章と 23 章, 丸善株式会社

GMH

Bridge GF(×10−2cm−1)(×10−2cm−1)

701 (1999).

RHF

10.2

10.0

[3] D.G. Fedorov, K. Kitaura, J. Chem. Phys.

FMO3-LC(VC)MO

11.0

10.2

120, 6832 (2004).

FMO2-LC(VC)MO

9.86

9.21

[4] S. Tsuneyuki, et al., Chem. Phys. Lett. 476,

Method[2] K. Kitaura, et al., Chem. Phys. Lett. 313,

	  最後に、トンネル電流を図２にプロットした。

104 (2009).

FMO2-LC(VC)MO の電子移動経路(c)と比べて、[5] T. Kobori, et al., J. Chem. Phys. 139,

FMO3-LC(VC)MO(b)は、マイナーなカレント分

094113 (2013).

布も含め RHF の結果(a)を良く再現している。

[6] H. Nishioka and K. Ando, J. Chem. Phys.
134, 204109 (2011).

4 SUMMARY[7] H. Kitoh-Nishioka and K. Ando, J. Phys.

Chem. B, 116, 12933 (2012).

Trp-(Pro)3-Trp モデル系に対する解析から、共[8] 鬼頭-西岡 宏任、安藤耕司 “光合成反応中心に有結合の切断を伴う FMO 計算では、三量体補正

おける長距離電子移動の経路解析” 第４章１節,	
 が電子移動経路の結果を大きく改善することを明

高次π空間の創発と機能開発(監修: 赤阪健、大須らかにした。今回の解析結果をふまえて、今後、

賀篤弘、福住俊一、神取秀樹) シーエムシー出版.

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 11スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

分子軌道計算を用いた水溶性ポリアセチレンの高次構造解析
寺田 佳世
奈良先端科学技術大学院 大学物質創成科学研究科

1緒言
Rh 触媒を用いた重合によって得られるシス立

体規則性の高い置換ポリアセチレンは、側鎖に光学活性な置換基を導入することで、巻き方向の片寄ったらせん構造を形成することが知られている。
このらせん構造は、熱、pH、極性溶媒などの外部刺激に応答して、らせん反転やピッチ変化、ランダムコイルへの転移などの構造変化を起こす。
我々はこれまでに、置換ポリアセチレンの刺激応
答性を活用した新規機能性生体材料の開発を目的として、水溶性および生体適合性を付与する置換基である双生イオン型のベタイン基をポリアセチ
レン側鎖に導入し、CD スペクトルなどを用いた高Figure 1. poly(1)の水中におけるCD および UV-vis スペクトル次構造の検討を行ってきた。その結果、例えばカルボキシベタイン基を側鎖に有するポリ（プロパルギルアミド）の poly(1)（Chart）は、有機溶媒

2

実験

中と水中でとり得る構造が異なること、温度変化
や pH 変化に応答して構造が変化することなどがまず spartan という分子モデリングソフトを用

明らかとなった。しかし、具体的にどのような二いて、分子力学計算(MMFF94、水中)により 1 の 12次構造変化が起きているのかについて、詳しい情量体の主鎖単結合の二面角（）を 10°ずつ変化

報は得られていない。そこで本研究では、分子モさせて分子力学計算によりエネルギーを求め、コデリング、分子シミュレーション計算によって、ンフォメーションの最適化を行った。ここで得らポリマー構造の最適化や溶媒依存性などに関するれたコンフォメーションの座標を用いて、解析を行うことを目的とした。

Gaussian 09、ZINDO/S 法によって Poly(1)の励起状態計算を行い、velocity モードにおける遷移状
態におけるパラメータである rotatory strength(Rvel)および oscillator strength (f)を低エネル
ギー側からそれぞれ 20 個算出した。
これらの値に半値幅を 14 nm としたガウス分布に従う波形の肉付けを行うことにより、CD および UV–vis スペク
トルをシミュレーションし、実測値との比較を試

Chart. Poly(1)の構造

12 − Academic Center for Computing and Media Studies, Kyoto University

みた。

3

結果•考察
1 の 12 量体の主鎖二重結合–単結合間の二面角()を 10°ずつ変化させてエネルギーをプロットした（Figure 2）。Poly(1)は、が-80°のときに最も安定であり、右巻きらせんより左巻きらせんの方が 9 kJ/mol・unit 安定であった。

Figure 2. 1‐12 量体の主鎖単結合の二面角に対するエネルギープロット

Figure 3. poly(1)の CD およびUV-vis スペクトルシミュレーション
結果ォメーションの予測は行えたが、有機溶媒中や高再安定となったが-80°のときの座標に対し

温での解析はまだ行えていない。MM 計算による最て、poly(1)の水中における CD および UV-vis スペ

安定構造の予測は、水以外の溶媒中においては困

クトルのシミュレーションを行った。その結果を難であるため、有機溶媒中のコンフォメーション

Figure 3 に示す。CD スペクトルにおいて 300–340

解析は他の計算手法を用いて構造最適化を行う必

nm 付近に長波長側から、負、正の分裂型のコット

要がある。分子軌道法(MO)法、密度汎関数(DFT)

ン効果が、UV–vis スペクトルにおいて 333 nm を法では、溶媒分子を連続誘電体として考慮したmax とする単峰性のピークが観測された。水中でCOSMO、SCRF、PCM 法や、溶媒の分子としての性質観測された実測の CD および UV-vis スペクトルを反映させた QM/MM 法を用いることができるため、（Figure 1）は、このシミュレーションとある程

有機溶媒中における構造解析が可能になる。また

度一致し、Poly(1)は水中で左巻きのらせん構造を

MO、DFT 法では異なる温度下での構造計算も可能

形成していることが示唆された。

である。そこで今後ソフトウェア MOPAC や
Gaussian を用いて MO 法や DFT 法で構造最適化を4行うことで、今回達成できなかった有機溶媒中や結言•今後の展望

高温において Poly(1)が形成する構造に関する詳
今回、ソフトウェア Spartan を使用し、分子力細な情報を得ることができると期待される。構造学(MM)計算による poly(1)のコンフォメーション

解析により得られた最安定構造の整合性についての最適化を行った。そして最安定となった座標には、今回と同様の手法によるスペクトルシミュレ

対して、ソフトウェア Gaussian 09 を使用し、ーションにより検討可能であると考えられる。まZINDO/S 法による励起状態計算を行い Rvel および f

た、今回は ZINDO/S 法のみを用いたが、TD–DFT 法を求めた。これらのパラメータを用いて、poly(1)による励起状態計算も合わせて行うことにより、の CD および UV–vis スペクトルシミュレーション

スペクトルシミュレーションの精度をより増すこを行った。その結果、MM 計算より再安定となったとも可能であると推測される。

1 の 12 量体の構造が、実際のポリマーが水中で形成する構造と類似していることを確認することに

5

謝辞

成功した。
今回の研究期間内に Poly(1)の水中でのコンフ

本共同研究制度(若手奨励枠)を活用させて頂きましたことを，この場を借りて感謝致します。

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 13

スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

核融合プラズマにおける第一原理に基づいた乱流輸送シミュレーション
今寺 賢志
京都大学大学院 エネルギー科学研究科 エネルギー基礎科学専攻

磁場閉じ込め核融合プラズマでは、密度勾配や温度勾配を自由エネルギー源として励起される不安定性の非線形発展の結果、系が乱流状態に移行し、複雑な輸送現象が現出する。このような現象を解明することを目的に、我々は第一原理に基づいた乱流輸送コードの開発を行ってきた。本稿では、平成 24 年度スーパーコンピュータ共同研究制度（若手研究者奨励枠）の基で行った(A)開放系プラズマ中の過渡的な変化に対する輸送の応答解析、(B)MHD-イオン系乱流間の多階層相互作用に関する解析、(C)トロイダルコードへの拡張について報告する。

1

研究背景

や有限体積法を基礎とした解法が広く用いられてきた。それらは、空間微分や数値流束を評価する

核融合プラズマでは、急峻な温度勾配を自由エ上で、物理的意味とは独立に、離散化したグリッネルギー源としてミクロな乱流が発生し、閉じ込

ド値をいかに「滑らか」に結ぶかを第一義の拘束めを劣化させる一方、非線形効果による自己組織

条件としている。一方、我々のグループでは、そ化過程を通して、マクロな渦構造や層流（帯状流）のようなシミュレーションに適用する数値計算手が生成され、閉じ込めが改善される。このような法として、支配方程式に基づいて物理的に意味のミクロとマクロが混在した複雑な揺らぎの微細構ある拘束条件の下で離散化を行う「多重拘束マル造やダイナミックスが、現象を大きく左右する。

チモーメント概念」を提案してきた [2]。この概

宇宙・天体プラズマでも、銀河・星形成時におけ

念では、複数の物理モーメント値を参照することる降着円盤内の乱流現象や粒子加速現象などで、で精度を確保し、更にそれらを支配方程式のモー同様の過程が重要な役割を果たす。

メントに沿って独立に時間積分することで、物理

このようなプラズマのダイナミックスをより正確に再現する方法として、6 次元の位相空間に格的に意味のある多重の拘束条件を課している（表1）。

子を設けて、速度分布関数を連続媒質として数値的に解く Vlasov シミュレーション（核融合プラ
ズマでは、磁力線による旋回運動を平均化した位
相空間5 次元ジャイロ運動論的Vlasov シミュレー
ション [1]）が世界的に広く用いられている。

Vlasov シミュレーションは、粒子シミュレーショ
ンに比べて数値ノイズが低い半面、高次元の位相空間を直接扱うため、格子点数が膨大となり、位相空間の微細構造を保持したシミュレーションは容易ではない。特に、位相空間の微細渦構造やセパラトリックス近傍での引き伸ばしや畳み込みが、数値散逸（数値エントロピー）を増大させる要因

表 1：従来の差分概念と保存型 IDO 法、多重拘束
マルチモーメント法の概念の比較となり、大きな障害となっている。
このような Vlasov シミュレーションの数値計算手法として、これまでの歴史では、有限差分法

14 − Academic Center for Computing and Media Studies, Kyoto University

我々は、この概念を用いて、第一原理に基づいた乱流輸送コード（以下 GKMUST コード）を開発し [3]、核融合プラズマ中の多階層乱流輸送に早く伝播しているのに対し、power を切ると遅く

関する解析を行ってきた [4]。本稿では、平成 24

なり、結果として径電場が強く形成されたものと年度スーパーコンピュータ共同研究制度（若手研考えられる。

究者奨励枠）の基で行った「開放系プラズマ中の過渡的な変化に対する輸送の応答解析」、「MHDイオン系乱流間の多階層相互作用に関する解析」、「トロイダルコードへの拡張」の 3 点について報告する。

2 開放系プラズマ中の過渡的な変化に対する輸送の応答
核融合プラズマ中の乱流輸送に対するこれまでの解析アプローチは、ある空間点の流束と勾配の一意的関係に準拠した「局所輸送」の考え方が基本であった。しかし、過渡的な変化に対する輸送の応答を調べる研究が進展するにつれ、単純な局所輸送モデルでは説明できない現象が数多く観測されるようになった。例えば、L（低閉じ込め）モードから H（高閉じ込め）モードへの遷移過程における広範囲での温度上昇の同時性や、ペレッ
トや不純物入射によって生じるコールドパルスの弾道的伝搬がそれに該当し、プラズマ輸送が質的に非局所・非拡散であることが示唆されている。
そこで本研究では、これまで行ってきた多階層
乱流シミュレーションの結果から、ミクロ（～mm）スケールの乱流と帯状流に代表されるメゾ（～cm）スケールの揺動間の非線形結合が非局

図 1: (a)乱流輸送に伴う熱流束、および(b)径電場を
縦軸、
温度勾配を横軸としたリサージュ曲線。Case(A)は t<4×105 で外部からの入力が on、t>4×105 でoff、Case (B)は常に off。

所・非拡散的なダイナミックスに深く関与しているとの着想に基づき、外部からの熱の入出力を考慮したシミュレーションにおいて、途中で input

power を切った場合にどのように乱流、および径

3 MHD-イオン系乱流間の多階層相互作用

電場が応答するかについて解析を行った。
図 1(a)は熱流束を縦軸、温度勾配を横軸とした典型的な磁場閉じ込め核融合プラズマでは、磁リサージュ曲線を示している。power を切ると温気流体力学的(MHD)不安定性の結果、イオンジャ度分布の緩和が起きるものの、power を入れない

イロ半径の数 10 倍程度の磁気島が形成される一ケースに較べて高い温度勾配で定常状態が維持さ方、温度勾配を駆動源として、イオン温度勾配

れていることがわかる。これは臨界温度勾配の上(ITG)不安定性が励起され、イオンジャイロ半径程

方遷移に対応するが、原因としては、図 1(b)から度の乱流構造が形成される。これまでは、単階層わかるように、power を切ると径電場が大きく生シミュレーションによってこれらは個別に解析さ成されたためと考えられる。その要因としては、れてきたが、我々はその両者の効果を考慮したジ

power を入れている間は、径電場が励起される特性時間と比較して弾道的な輸送が半径方向により

ャイロ運動シミュレーションを行うことで、それらの多階層間相互作用を解析した [5]。

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 15

図 2 は、各磁気島幅に対する ITG モードの成長

張し、まず線形においてベンチマークテストを行率のポロイダル波数依存性を示している。磁気島

った。用いたジャイロ運動論的 Vlasov 方程式を以によるモード間結合の結果、異なる波数の成長率下に示す。

が等しくなり、
高波数における Landau 減衰や有限
ラーマ－半径効果(FLR)による散逸によって全体の成長率が減少することがわかった ( w = 6,12 ) 。
一方で、流体モデルでの結果 [6]と比較すると、その減少は小さく、さらに、磁気島が十分大きい場合 ( w = 18 ) 、逆に不安定化されることが確認された。これは、磁気島が増大するとそれに伴って新たな有理面が形成され、モード幅が狭い短波長
ITG モードがその影響をより強く受けるためであることが新たに明らかとなった(図 3 参照)。

∂f dr ∂f dθ ∂f dϕ ∂f dv ∂f++++= 0 (1)∂t dt ∂r dt ∂θ dt ∂ϕ dt ∂v
1 ∂ Φ φ v + µ B
dr
sin θ
=−
−
dt
rB ∂θ
R0
2

2
1 ∂ Φ φ v|| v + µ B
dθ
cos θ
=
+
−
dt rB ∂rqR
rR0

dϕ v||=
dt R
dv||
1 ∂ Φφ 1∂ Φφr sin θ
=−
−
− µB
dt
qR ∂θ
R ∂ϕ
qR
+

∂ Φ φ cos θ v  ∂ Φ φsin θ +∂θR0  ∂rr −

1 ∂n0 ( r )  ∂ 2Φ 1 ∂ 2Φ
∂ 2Φ  1−−+∂r 2  r n0 ( r ) ∂r  ∂r 2 r 2 ∂θ 2

+Φ =図 2: 異なる磁気島の幅 ( w = 0,6,12,18 ) に対す
るイオン系乱流の分散関係

1
f φ dv||d µ − n0 ( r )n0 ( r ) ∫∫

(

)

(2)

( r, θ , ϕ , v , µ ) はジャイロ中心座標、f はイオン分
布関数、 < Φ >φ はジャイロ平均された静電ポテ
ン シ ャ ル 、 n0 ( r ) は イ オ ン の 平 衡 密 度 、
B = B0 R0 / ( R0 + r cos θ ) は磁場をそれぞれ表している。5 次元の位相空間に格子を設けて速度分布関数を連続媒質として数値的に解く Vlasov シ
ミュレーションの場合、数値流体力学と同様のア
プローチが適用可能であるが、問題点としては、


高次元の位相空間を直接扱うため、格子点数が膨大となり、高効率並列化が必要



位相空間の微細渦構造やセパラトリックス
近傍での引き伸ばしや畳み込みが、数値散逸

図 3: 異なる磁気島の幅 ( w = 0,6,12,18 ) に対す
る静電ポテンシャルの構造（数値エントロピー）を増大させる要因となるため、適切な数値計算手法が必要が挙げられる。そこで本研究では、まず前者の対4策として、 ( r, θ , µ ) の 3 次元領域分割法を採用しトロイダルコードへの拡張
2 章、および 3 章はスラブ配位における解析結た。計算のフローチャートは以下の通りである。
(1)f ( r, θ , ϕ , v , µ ) を Fourier 変 換 し 、
f ( r, kθ , kϕ , v , µ ) を計算。その際に座標を集果であったが、より現実に近いトロイダル配位の団通信(MPI_ALLtoALL)で転置させることを解析を行うことは極めて重要である。そこで本研利用して、 ( r, θ , µ ) 領域分割を ( kθ , kϕ , µ ) 領

究では、GKMUST コードをトロイダル配位に拡

域分割に変更。

16 − Academic Center for Computing and Media Studies, Kyoto University

(2) 領域分割されていない r 方向に対して有限差分 法 を 用 い る こ と で 式 (2) を 計 算 し 、

Φ ( r, kθ , kϕ ) を求める。

(3) 逆 Fourier 変換を施すことで Φ ( r, θ , ϕ ) を求める。その際に座標を集団通信で転置させること を 利 用 し て 、 ( kθ , kϕ , µ ) 領 域 分 割 を( r, θ , µ ) 領域分割に戻す。
(4) ( r, θ ) 方向に関する Φ ( r, θ , ϕ ) の裾値を、隣接する CPU に一対一通信(MPI_SENDRECV)で送る。
(5) 式 (1) を 有 限 差 分 法 を 用 い て 計 算 し 、f ( r, θ , ϕ , v , µ ) を時間発展させる。
(6) ( r, θ ) 方向に関する f ( r, θ , ϕ , v , µ ) の裾値を

図 5：初期値で規格化された径電場の時間発展。
実赤線が数値計算結果、青点線が理論計算で得られた収束値 [9]に対応。

隣接する CPU に一対一通信で送る。

後者の対策としては、前述の多重拘束マルチモーメント法を用いる予定であるが、初期のベンチマ
ークテストのため、本計算では 4 次精度中心差分法を用いた。
図4 (a)はCyclone DIII-D base parameter [8]と呼ばれる典型的な物理パラメータ－を用いて行った際のイオン温度勾配乱流(ITG)モードの線形成長率を、図 4 (b)はポロイダル断面における最も不安定なモードの静電ポテンシャルを示しており、他の文献とよい一致を見せている。また図 5 は、帯状流の無衝突減衰テスト [9]の結果であり、他の文

5

総括と今後の課題
本研究では、スラブ配位における開放系プラズ

マ中の過渡的な変化に対する輸送の応答解析と、
MHD-イオン系乱流間の多階層相互作用に関する
解析、さらにはスラブコードのトロイダル配位への拡張を行った。
今後は、トロイダル GKMUST コードを非線形に拡張し、スラブ配位において明らかとなった過渡的な変化に対する応答や、多階層相互作用を軸に解析を進めていく予定である。

献に加え、理論計算で得られた収束値ともよい一致を見せている。
(a)

Reference(b)[1] Y. Idomura, L. Villard and T.-H. Watanabe, Nucl.
Fusion 50, 043002 (2010).

[2] T. Kawano, K. Imadera, J. Q. Li and Y. Kishimoto,Plasma Fusion Res. 6, 2401097 (2011).

[3] J. Q. Li, K. Imadera, Z. X. Wang, P. Hilscher and Y.
Kishimot, Progr. Nucl. Sci. Tech. 2, 64 (2011).

[4] K. Imadera, Y. Kishimoto and J. Q. Li, Plasma
Fusion Res. 5, 019 (2010).

図 4：Cyclone DIII-D base parameter を用いた際のITG モードの(a)線形成長率と(b)最も不安定なモードのポロイダル構造[6] P. P. Hilscher, K. Imadera, J. Q. Li and Y. Kishimoto,Plasma Fusion Res. 8, 2403040 (2013).

[7] Z. X. Wang, J. Q. Li, Y. Kishimoto and J. Q. Dong,

Phys. Plasmas 16, 060703 (2009).

[8] A. M. Dimits, et al., Phys. Plasmas 7, 969 (2000).

[9] F. L. Hinton and M. N. Rosenbluth, Plasma Phys.
Control. Fusion 41, 653 (1999).

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 17スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

実河川に繁茂する植生の揺動の数値解析
岡本 隆明
京都大学大学院工学研究科社会基盤工学専攻水理環境ダイナミクス分野

1 研究目的を利用することで計算格子の数を大幅に増やすことが可能となり(図 1)，１つ１つの柔軟な植生要素

河道内植生群落は河道の流下能力に大きな影響

を計算領域上に正確に再現した(図 2) ．柔軟植生を及ぼしている．環境保全上，河道内の全ての植は水流によって大きく変形するが，変形形状は植生群落を除去することは不可能であるため，流水生群落内での位置によって大きく変化することがの阻害となる植生群落のみを正確に評価して伐採わかった．すなわち，植生群落の上流端では植生することが求められる．しかしながら，実河川のは水流によって大きくたわみ，
倒伏する．
しかし，

植生の多くは柔軟性を有し流れによって変形する

上流側の植生によって流れが遮蔽されるため，下ため，植生の伐採位置・量による流下能力の改善

流にいくにしたがって植生のたわみは小さくなる効果については十分な知見は得られていない．

(図 3)．

そこで本研究では植生群落内外の流れを水路実

乱流構造の変化に着目して研究し，柔軟植生の

験及び数値解析することで，植生群落の伐採によ

変形は大規模組織渦の発達過程とも密接な関係がる抵抗特性・粗度係数の変化を定量評価し，河道あることがわかった．すなわち，植生群落の上流

内樹木群落の維持管理のための知見と手法を提供

側では植生の揺れは小さいが，群落内を下流にいする．

くにしたがって植生は起き上がり，組織的な揺動

植生にかかる抗力を数値解析し，植生変形形状をするようになる．これは植生先端部には組織渦を把握することで正確に評価することで植生群落が発生したが下流にいくにしたがって発達していの伐採位置・量による抵抗特性・粗度係数の量的ったためである(図 4)．

に知ることができ，植生群落を保護した上で河川

また計算領域を大きくとることが可能となった

管理するための手法を提供する．また本研究ではので，横断方向の計算格子数を増やして数値解析

植生にかかる抗力を数値解析し，植生変形形状をし大規模渦の３次元構造について考察した．数値

正確に評価することで数値計算モデルも精度を向解析結果から大規模渦は横断方向に位相差をもち，上させる．これにより樹木群落による流況変化をかぶり水深の低いケースほど組織構造の横断方向

正確に数値シュミレーションすることが可能とし，スケールが小さくなることを解明した(図-5)．植生

樹木群落による河床の地形変化・群落の拡大状況

流れの組織乱流が横断方向に非一様な構造をもつを将来予測するなど有意義なデータを提供する．

ことは既往研究(Ghisalberti & Nepf (2005))の可視化実験によって確認されているが，詳細な実験デ2 計算結果ータは得られていない．今後植生乱流の研究を進めていく上で，今回スーパーコンピュータを用い

本研究では Large Eddy Simulation (LES)を用いて，柔軟植生流れを数値計算した．植生変形によって抵抗係数は大きく変化することが知られているため，植生の変形形状を詳細に調べた．表 1 に本研究の計算条件を示す．スーパーコンピュータ

18 − Academic Center for Computing and Media Studies, Kyoto University

た数値計算結果から大変有用なデータが得られた．

表-1 計算条件
λ = ahL v (cm)

B v (cm)

H (cm)h (cm)H/h

U m (cm)

B-10

0.78

1.6

3.2

15.0

5.0

3.0

12.0

A-10

1.56

1.6

1.6

15.0

5.0

3.0

Fr

12.0

0.10

0.10

C-10

0.39

3.2

3.2

15.0

5.0

3.0

12.0

0.10

C-21

0.39

3.2

3.2

6.3

5.0

1.3

10.0

0.13

C-20
C-22
C-23

C-24
C-25
C-26

0.39

3.2

0.39

3.2

3.2

0.39

3.2

3.2

0.39

3.2

3.2

0.39

3.2

3.2

0.39

3.2

3.2

3.2

5.0

5.0

7.5

1.0

5.0

10.0

20.0

0.10

10.0

3.0

5.0

0.12

10.0

2.5

5.0

0.14

10.0

2.0

5.0

15.0

10.0

1.5

5.0

12.5

hd 0.9
h 0.8
1

Case

0.09

10.0

4.0

0.08

10.0

0.07

Re
18000
18000

0.7

18000

0.6

5000

0.5

6300

0.4

7500
10000

hd min h
hd max h

0.3

12500

0.2

15000

hd h

0.1

20000

0

0

10 20 30 40 50 60 70 80 90 100 110 120 130 140 150

x (cm)

図-3 柔軟植生変位の流下方向変化
Free surface
(a)

(a)

t=0.0(s)3Vegetation
elements

y/h

2

2

p1
0

p

1

Ejection(A)
0

5

10

(b)

15

t=3.4(s)30x / Lv 20

Bv

H

h

y

Flow

V

x

p
Sweep(B)0510~ (cm/s)
u

-1

0

+

q

0

20

Sweep(C)

0

0

5

12

Sweep(B)q
Ejection(D)0x / Lv 10

5~ (cm/s)ux / Lv

t=3.4(s)z / Bv-1z

50

Ejection(A)4x / Lv 10

5

t=3.4(s)80Vegetation
element- uv (cm/s)2

t=0.0(s)8Sweep(C)4

900 +1200

10z / Bv
q

12

600

5
-50

z / Bv

300

x / Lv

Sweep(B)010Ejection(A)

8

1500

500

0

12

4

1000

15

1

x / Lv 20z / Bvq015

t=0.0(s)12yx / Lv 20

布（縦断面）

W

図-1 計算領域

15

図-4 瞬間流速ベクトルと瞬間レイノルズ応力分

z

B

10

t=3.4(s)

2

p1

U

5

(b)

y/h

2

0

Ejection(A)
0

3

y/h

Lv

t=0.0(s)

3

y/h

10

Sweep(B)840Ejection(D)
0

10

5
-50

- uv (cm/s)21050

x / Lv

図-5 瞬間流速ベクトルと瞬間レイノルズ応力分
布（水平面）

図-2 計算メッシュ

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 19

スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

流水中の自然石に働く流体力の評価
吉田 圭介 1 , 牛島 省 2 , 宮木 伸 3, 竹内 章人 1
1 岡山大学大学院 環境生命科学研究科, 2 京都大学 学術情報メディアセンター, 3 京都大学 工学研究科

1 諸言
近年，河川のインフラ整備では構造物の堅牢さに加えて自然環境への調和が要求されており，現地で調達した自然の素材（石や木など）を用いた河川構
造物が河道内に設置される場合がある．例えば，適度な空隙を有する水制工や根固工（図 1）や水深寸法以上の巨石を河道に配置した巨石投入工法は，川の流れに多様性をもたせて河川景観を創造したり，水生生物の生息環境を保全するなどの効果がある．
しかし，コンクリート構造物とは異なり，自然石を用いた工法では材料の要素が画一ではなく，構造物の流水に対する力学的応答（治水機能）が明確ではない場合がある．現在，実務では技術者の経験や既往の実験例（マニュアル[1]など）からこうした構
造物の安全度を推定するのが通例であるが，これには労力や費用が多くかかる上，射流などの急変流場
では計測そのものが困難であるため，設計の信頼性
は必ずしも十分ではないと指摘されている[2]．
そこで，本研究では固気液多相場の数値解析法[3]を利用して，一様な開水路等流場の底面に設置された自然石に働く流体力を解析し，併せて行った室内
模型実験との比較から数値解析法の妥当性を検討することとした．本研究の最終目標は，流水中の自然石に働く流体力を精度良く評価できる解析手法を確立し，水工学の実務に資することである．

図 1 根固工の例（椹野川）
勾配 1/1000

0.4m

z
0.16m
10m

フラップゲート

y

x

ピット

1m

5m

0.6m

図 2 実験水路，計測ピット及び座標系
表 1 実験条件表
Case
1
2
3

流量Q [ l /s ]
40
60
80

断面平均流速V [m/s]
0.33
0.50
0.67

2 実験
本研究では長さ 16m，幅 0.6m，高さ 0.4m，勾配
1/1000 の可変勾配型循環水路（図 2 参照）
を用いた．
座標系は図に示す通りである．実験では自然石をア
クリル製の水路底面に接着して通水し，流水中の流体力と周辺の流れ場を測定した．河川の自然石は通常，上流から下流に沿って徐々に扁平な形状を呈す
るが，本研究では図 3 に示すように表面の凹凸や粗さが明瞭な石を用いている．

20 − Academic Center for Computing and Media Studies, Kyoto University

図 3 実験で用いた自然石（左：石 A, 右：石 B ）流体力の計測では水路のピット内中央部に四分
力計(東京計測社製，Y116M2 型)を設置し，自然石に働く抗力・揚力を 30Hz の時間解像度で測定した．
水流の流速計測には 2 次元電磁流速計(KENEK 社
製，VM-802H 型)を用い，開水路の底面境界層が十分に発達して，主流速の鉛直分布が対数則分布を示

図 4 自然石の 4 面体分割（左：石 A, 右：石 B ，
図 3 参照）
表 2 数値解析条件表
時間刻み Dt [s]計算時間T [s]空間格子幅 (Dx , Dy, Dz ) [mm]5.0×10-4
25.0(5.0, 4.8, 3.8)

計算領域 (Lx , Ly , Lz ) [m](1.0, 0.125, 0.1)

流入条件
流出条件
側壁・水面条件
底面条件

断面平均流速自由流出
Slip 条件
Nonslip 条件

すことを確認した．また，水深の計測にはポイント
ゲージを使用し，等流場であることを確認した．実験条件は表 1 に示す通りである．

図 5 自然石周辺の流れ場の様子（水路中央の鉛直断面，流速ベクトルと渦度コンター，石 A）

図 6 自然石周辺の流れ場の様子（底面近くの水平断面，渦度コンター，石 A）
0.1

3 数値解析
本研究では牛島ら[3]の数値解析手法を用いて，自然石に働く流体力を解析した．牛島らの手法では固気液多相場を一流体モデルとして扱い，混合体に対する連続式と運動方程式を非圧縮性条件の下でDNS（直接数値計算）により解く．また，任意形状の固体は四面体要素の集合として表現され，固体領域に働く流体力は抗力係数などの経験的なパラメータからではなく，直接的に適切に算定される．
本研究では実験で用いた自然石（図 3 参照）の 3
次元形状測定を行って，数値解析に必要なデータを作成した（図 4 参照）．数値解析の条件は表 2 に示
す通りである．解析では自然石を剛体とみなした．
また，解析では直方体格子を用い， x , y, z 各方向に均一に格子分割を行っている．なお，負荷軽減のため，解析では領域分割による並列計算を行っている．

4 結果と考察

抗力D(N)

0.08
実験値:A

0.06

解析値:A
実験値:B

0.04

解析値:B

0.02
0
20

40

60

80

100

流量(l/s)

図 7 流量と抗力の関係（石 A と石 B）の抗力の算定値を示す．実験値と解析値を比較すると，定性的には流量と抗力には比例関係が認められる．一方，定量的には解析値の方が小さい値を示す．
この要因として，境界条件や，自然石と水路底面の接着条件など，詳細な条件の相違が考えられる．今後は条件の詳細を検討するなどして，流体力評価の再現性を向上させたい．

5 参考文献[1] (財)土木研究センター：護岸ブロックの水理特性試験法
マニュアル（第 2 版）
，2003.

図 5，図 6 は数値解析結果の一例であり，水路中央の鉛直断面における流速場と渦度コンター，および底面近くの水平面内の渦度コンターを各々示した．自然石周辺では周期的な渦構造が生成され，それらが移流された後に拡散する様子が観察される．
図 7 には水理実験と数値解析から求めた，流量毎[2] 忰熊公子，前野詩朗，吉田圭介，高田大資，山村明：
射流場の護岸・根固めブロックに作用する流体力評価，
水工論文集，第 57 巻，I_691-I_696, 2013.
[3] 牛島省，福谷彰，牧野統師：3 次元自由水面流中の接触を伴う任意形状物体運動に対する数値解法,土木学会
論文集 B，Vol.64，pp.128-138，2008.

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 21

スーパーコンピュータ共同研究制度（若手研究者奨励枠）研究報告

水～土骨格連成有限変形問題を対象とした土構造物・地盤系の固有振動解析
山田正太郎 1)

野田利弘 2)

1) 名古屋大学大院 工学研究科 社会基盤工学専攻
2) 名古屋大学 減災連携研究センター

1 はじめに
本研究では、水～土骨格連成有限変形解析コー
ド GEOASIA1) により、地盤とその上部に構築された盛土を対象に土構造物・地盤系全体の地震直前の固有振動数および固有振動モードを算出し、地盤と土構造物がそれぞれ大きく変形するモード

図 1 有限要素メッシュと境界条件を抽出する。さらに、各モードに対応する固有振動数付近に卓越振動数を有する 2 種類の地震波を入力波として、地震応答解析を行い、卓越する固有振動モードの違いが、土構造物・地盤系の変形
破壊挙動に大きく影響を与えることを示す。

(a) 地盤が大きく変形するモード(MODE 1)

なお、同解析手法から得られる固有値問題は複素固有値・複素固有ベクトルを解とする大規模なものとなるため、スーパーコンピュータを利用して、上記の問題を解くことにした。

2 解析条件計算は二次元平面ひずみ条件で行った。水平地
盤上に盛土を施工した場合を考えた。計算に用いた有限要素メッシュと境界条件を図 1 に示す。この図に示す土構造物・地盤系に対し、固有振動解
析と地震応答解析を行った。

3 固有振動解析結果
固有振動解析によって得られた初期(圧密放置後)の固有振動モードのうち、地盤が大きく変形するモード(MODE 1)と盛土が大きく変形するモード(MODE 2)を図 2 に示す。なお、本論文では固有振動モードが複素固有ベクトルの実部と虚部の線形結合で表されるが、ここでは、複素固有ベク
トルの虚部で表されるモードのみを示している。

22 − Academic Center for Computing and Media Studies, Kyoto University(b) 盛土が大きく変形するモード(MODE 2)

図 2 固有振動数および固有振動モード(虚部)4 地震応答解析結果
入力地震波の加速度履歴およびフーリエ振幅ス
ペクトルを図 3 に示す。
WAVE 1 は MODE 1 に、WAVE 2 は MODE 2 に対応する固有振動数に卓越振動数を合わせた地震波である。地盤底面の全節点の水平方向に課せられた粘性境界に、これらの地震動を入力した。
図 4 に地震中から地震後にかけてのせん断ひず
み分布の推移を示す。WAVE 1，WAVE 2 どちらを入力した場合でも地震後に遅れ破壊を生じている様子が見て取れる。せん断ひずみの進展過程に着目すると、WAVE 1 では、地震中に砕石層上端から粘土層中央へとせん断ひずみが進展していることが分かる。その後、時間の経過とともに地盤から盛土へとひずみが進展し、遅れ破壊を生じている。一方、WAVE 2 では、地震中に砕石層上端

にひずみが集中するものの、盛土、地盤ともに大感謝いたします。

きな変形は生じていない。しかし、時間が経過す
るにつれて、まず盛土に滑り面が発生し、その後地盤内へ滑り線が進展していることが分かる。つまり、
地震波が有する卓越振動数の違いによって、遅れ破壊過程に大きな違いが見られることが分かった。

Acceleration a (cm/sec2)

300
200
100
0
-100
-200

Fourier amplitude spectrum (gal・sec)

-300
0

10

Time t (sec)

20

30

10

MODE 10.699(Hz) 卓越振動数：0.708(Hz)86(a) WAVE 1 を入力した場合（地盤から盛土へ破壊が進展）

4
2
0
0

1

2

3

Frequency f (Hz)

4

5

6

(a) WAVE 1

2

Acceleration a (cm/sec )

300
200
100
0
-100
-200

Fourier amplitude spectrum (gal・sec)

-300
0

10

Time t (sec)

20

30

10
8

MODE 2
1.935(Hz)

卓越振動数：1.929(Hz)

6
4
2
0
0

1

2

3

Frequency f (Hz)

4

5

6

(b) WAVE 2 を入力した場合（盛土から地盤へ破壊が進展）(b) WAVE2

図 3 入力地震波
図 4 せん断ひずみ分布

5 結論
本研究では、地盤の初期値・境界値問題に対して、固有振動数および固有振動モードを計算した
上で、土構造物・地盤系に対して地震応答解析を行った結果、地震動の卓越振動数が異なれば、異なる固有振動モードが卓越して現れ、地震中、地震後の破壊を含む地盤変形の進行はその固有振動
モードの影響を強く受けることを示した。

6 参考文献1) Noda et al. (2008): Soil-water coupled finite
deformation analysis based on a rate-type
equation of motion incorporating the SYS
Cam-slay model, Soils and Foundations, 45( 6),pp. 771-790.

最後に、スーパーコンピュータ共同利用制度に京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 23

プログラム高度化支援事業研究報告

コンクリート材料の物質拡散・非線形力学を連成した経年劣化シミュレータ
浅井 光輝
九州大学大学院工学研究院 社会基盤部門

コンクリートの経年劣化は，中性化・アルカリシリカ反応・塩害等さまざまな要因がある．その多くは力学的要因と化学的要因が複雑に連成しており，劣化メカニズムを詳細に把握したうえ
で合理的な対策することは困難である．そこで，コンクリート内部での物質の浸透・拡散現象と同時に，亀裂・はく離損傷等の現象を再現可能なシミュレータを開発した．
同シミュレータを用い，コンクリートアリカリシリカ反応時に観測される自由表面での複雑なひび割れが定性的に再現することができた．材料劣化の促進実験結果との比較検証することで計算の信頼性を確保するためにも大規模計算が不可欠であり，現在は大規模並列計算への対応を実施中である．

1 緒言

そこで本研究では，まずコンクリート材料のさ
まざまな劣化要因の中から ASR を取りあげ，連コンクリート構造物の経年劣化の要因は，中性

成解析を実施することで，ASR 現象の特徴的な損化・アルカリシリカ反応(ASR)・塩害等のさまざ傷である亀甲状のひび割れパターン（Fig.1 を参まな現象が考えられる．その多くは力学的要因と照）を再現することを目標とし，コンクリート材

化学的要因が複雑に連成しており，劣化メカニズ

料をモルタルと骨材，そして両者の界面の 3 相にムを詳細に把握したうえで合理的な対策を行うこ

分けたメゾコンクリート材料モデルを用いた数値とは非常に困難である．コンクリート内部での物解析モデルの確立を試みた．

質の浸透・拡散現象は，材料内部の膨張あるいは劣化の起点となり，微視的な亀裂，材料界面での2 拡散膨張・不連続面進展の連成解析

はく離等の損傷を招く．また逆に，微視的な損傷を受けた領域は，拡散現象において高速に物質がASR とはコンクリート材料の母材であるモル

拡散する経路となり，物質が局所的な領域に拡散

タル内部のアルカリイオンと介在物である骨材がしやすくなる．以上のように，拡散問題と力学の反応し，骨材表面にゲル状の物質が生成されるこ損傷問題は両者が双方向に連成する現象であり，とで内部膨張が生じ，この内部膨張に伴い材料界

材料劣化の微視的なメカニズムを解明していくに

面のはく離，あるいは母材内部に微視的なひび割は，連成現象までを再現可能なシミュレータが必れが進展していくものと考えられている．この一要となる．

連の流れをモデル化するものとした．
コンクリート中へのイオンの浸透には拡散，移流，吸着などの現象が考えられるが，本研究においては拡散現象のみを考慮するものとした．基本
的には，非定常拡散解析と内在物膨張に伴う不連続面進展解析とも有限要素法で解析を行い，この両者を組み合わせることで，拡散膨張を考慮した
不連続面進展解析手法の構築を行うことにした．
写真 1 ASR による損傷事例

24 − Academic Center for Computing and Media Studies, Kyoto University

このコンクリートのアルカリ骨材反応を想定した

連成解析の手順を図 1 に示す．

な破壊状態を表す．損傷に伴う微小空隙の発達は解析手順としては，まず非定常拡散問題を解き拡
散物質の空間分布を予測する．次に，濃度の空間分

材料の剛性低下を引き起こすため本研究では弾性係数の低下によって損傷状態を表わす．

布の結果を基に膨張力を決定する．この際，膨張力

D=によって引張り応力が基準値を超え不連続面と判

E0 − EE0(2)

定される領域（要素）は，拡散係数の高い仮想空隙

E0 は非損傷状態の弾性係数，E は損傷後の弾性領域に置換し，再び拡散問題を解く．以上の手順を係数である．Ce を一般的な等方弾性テンソルとす

繰り返すことで，浸透・拡散に伴う不連続面進展解ると，損傷を考慮した材料構成則は以下のように析を実施する．

規定される．

NON-STATIONARYDIFFUSION ANALYSIS

σ = (1 − D )Ce : ε econcentration of ion

コンクリートは引張り応力に弱い材料であり，COUPLING ANALYSIS ①
next step

圧縮強度に比べて引張強度は，約 1/10 程度である．

Estimation of expansion force
from concentration of ion

そこで，圧縮と引張強度の相対的な影響度を考慮

expansion force

した相当ひずみを用いて破壊損傷を定義すること

DISCONTINUOUS SURFACEDEVELOPMENT ANALYSISが望ましい．本研究では Peerlings ら 1)が提案した圧縮域と引張域の強度比を考慮した相当ひずみ

crack propagation

図1

COUPLING ANALYSIS ②(4)の両者を用いることにした．

Element determined as crack
is displaced for porous material

ε eq =

拡散膨張・不連続面進展の連成解析フロー

2.1 ASR 反応に伴う骨材膨張のモデル化
カリ反応性骨材と微細空隙中のアルカリ溶液が化学反応することにより生成するアルカリシリカゲ
ルが吸水，膨張し骨材周辺からひび割れが進展する．
本研究ではアルカリシリカゲルの吸水による体積膨張を骨材の体積膨張として表現した．具体的には，骨材要素のみに次式に示すように濃度に比例した
膨張ひずみを与えることにした．
c

c
x

c
ε cy ε zc γ xy
γ cyz

= {αc αc αc 0 0 0}T

}

c T
γ zx

m −1
I1
2m(1 − 2ν )2

1  m −1 
2m+
I1  +J
2m  1 − 2ν 
(1 + ν )2 2

(4)

ここで，I1 はひずみの 1 次不変量，J2 は偏差ひず

アルカリシリカ反応では，コンクリート内のアル{ε }= {ε(3)

みの 2 次不変量であり，m はコンクリートの圧縮
強度と引張強度の比，νはコンクリートのポアソン比を示す．
この相当ひずみを用い，損傷の進展則に従うも
のとした．
D =1−

κ0
[1 − A + A exp{− B(κ − κ 0 )}]κ

(5)

ここで A，B は損傷の進展を表すパラメータ，κ
は材料が過去に受けた最大の相当ひずみである．

(1)

ここで，αは膨張係数，cはイオン濃度を表す．この膨張ひずみをもとに膨張力を評価し，力学計算においてはこの力を駆動力とした非線形問題を解く．

図 2 に示すように，κがκ 0 となると損傷が開始するものとし，対応して損傷変数 D は非ゼロの値を持ち始め，最終的には完全な損傷状態を示す 1.0に漸近する．
D
1.0

2.2 損傷モデルによる不連続面の表現
損傷モデルでは連続体損傷力学を導入し，損傷
変数 D を用い，要素剛性（材料割線係数）を段階的に減少させる．損傷変数は 0≦D≦1 であり，
D=0 のときには健全な状態を示し，D=1 は最終的

0 κ0

κ

図 2 損傷パラメータ D の発展則

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 25

2.3 損傷係数による拡散係数の変動

3.2 解析結果

先に説明した損傷モデルを用いれば，微視的なク解析結果として，ここでは本研究で採用した損ラックが成長することで段階的に剛性を低下して傷モデルを用いた場合と，簡易的に完全に脆性的いく過程を表現できる．その損傷係数を参照しなが材料として扱い瞬時に応力解放するモデル（完全ら，拡散係数も空隙に相当した大きな値c Bへと段階

脆性モデル）を採用した場合の結果の両者を示す．

的に変化させることにした．具体的には，損傷係数

損傷モデルを用いた場合の結果として, 図 4 にはを用いて以下の式により損傷後の拡散係数cDを与えアルカリイオン分布，応力分布，損傷領域(損傷定るものとした．

数が 1 の領域)を示し，図 5 は 2500 日後の自由表

cD = c + D(cB − c )(6)

面上の損傷パラメータの分布図を示す．また，図6 には完全脆性モデル採用時の結果を示す．

3 解析事例

1.0g/mm33.1 解析条件
上記の解析手法を用い，3 次元コンクリートモ
デルを解析した結果を報告する．図 3 に示す解析モデルは 100mm 立法のコンクリートを研磨し，1mm 間隔で断面をスキャンした画像から作成し0.0g/mm3(a) アルカリイオン分布
4.0MPaた 2)．画素を 1mm に調整し，ボクセル要素によりモデルを作成したため，総要素数は 100 万であ
0.0MPaる．

(b) 最大主応力分布yxz(c) 損傷領域分布

図 3 3D コンクリート解析モデル

250 日後解析に用いた物性値は表-1 にまとめて示す．そ

図4の他，損傷モデルのパラメータは A=0.99，B=100，

1,000 日後

2,500 日後損傷モデルによる解析結果

m=10 とした．境界条件としては，yz 平面以外，

全面の変位を拘束し，非拘束面である yz 平面(x=1)から拡散物質を浸透させるものとした．
表-1 材料パラメータ一覧
弾性定数

ポアソン

κ0(損傷開始相

(GPa)

比

当ひずみ)

モルタル

20

0.25

2.0×10-4

界面

20

0.25

5.0×10-5

骨材

60

0.25

―

26 − Academic Center for Computing and Media Studies, Kyoto University

図52500 日後の自由表面内の損傷定数分布

3 次元問題のひび割れの進展を簡易的に表現でき

1.0g/mm3

るボクセル FEM をベースとした．
ボクセル FEM
では，要素サイズが小さいことを前提として，各要素は完全に脆性的な材料モデルを用いた例も多0.0g/mm3

いが，実際には有限の要素サイズを有しており，(a) アルカリイオン分布

またそのサイズは完全に微細なひび割れ間隔とは一致しないため，本研究では損傷理論を導入する

4.0MPaことでその矛盾の解消を試みた．
両者を連成させる基本的なフローを提案し，内部膨張により不連続面が進展する様子を解析した．

0.0MPa

この際，ひび割れ進展のモデル化方法として損傷(b) 最大主応力分布モデルを導入し，提案手法の妥当性について検討した．
損傷モデルを導入した連成解析手法では，損傷部分の剛性を段階的に減少させることで，モルタ
ル部分にも損傷が段階的に進展する様子を再現す(c) 損傷領域分布
250days

1,000days

ることが可能となった．また，解析に用いるパラ2,500day

メータの値を変更すれば，モルタルへのひび割れ
発生時間，およびひび割れ損傷幅などの損傷進展

図6

完全脆性モデルによる解析結果

挙動を調節可能であることを確認している．
今後は，解析対象を拡大させ，実際に供用時に発生するコンクリート板での亀甲状のひび割れの図 6 に示す完全脆性モデルを採用した場合には，再現等を試み，本手法の妥当性を検証していく予

アルカリイオンが骨材に到達すると応力が高くな定であり,それに向けた大規模並列化を実施中でり，骨材周辺のみに不連続面が進展していく．

ある．これに際して，京都大学学術メディアセン

2,500 日後以降には，ひび割れの進展が止まり，

ター平成 24 年度プログラム高度化共同研究支援

モルタルへとひび割れが進展せず，表面での損傷事業に採用いただいた成果として,現在では 1000

は観測できなかった．

万要素以上を用いた解析が可能となった．また，一方，図4に示す提案手法である損傷モデルを導平成 25 年の同支援事業にも採択いただき，現在入した結果では，段階的に損傷が進展しながら，骨はメモリ分散化を実施しており，数千万から億単

材周辺での膨張力が継続して外部へと伝達するこ位の要素を用いた解析を目指している.

とで，ひび割れがモルタル中へと進展する過程を解析できた．結果として図5では，アルカリシリカ反応の特徴である表面上の亀甲状のひび割れパター

参考文献

ンを定性的に表現することができた．

1) R. H. J. Peerlings, R. de Borst, W. A. M.
Brekelmans

4 結語

and

M.

G.

D.

Geers,

Gradient-enhanced damage modeling of
concrete

本研究では，物質の濃度に併せて介在物表面がfracture,

Mech.

Cohes.-Frict.

Mater., 3, pp.323-342, 1998.

膨張し，かつ膨張時には物質を消費しないものと

2) 浅井光輝，山岸道弘，寺田賢二郎，永井学志，

仮定した上で，非定常拡散解析と内在物膨張に伴非局所型有限要素法の開発とその破壊挙動解

う不連続面進展解析とを組み合わせた連成解析手

析 へ の 適 用 ， 土 木 学 会 論 文 集 , No.759,

法を開発した．解析手法としては，非均質材料のpp.233-245, 2004.

微視構造のモデル化に優れ，かつ要素除去により京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 27

プログラム高度化支援事業研究報告

GPU ίϯϐϡʔςΟϯάʹΑΔΠϝʔδϕʔε FIT ͷޮ཰Խதാ ࿨೭ ∗∗

ѪඤେֶେֶӃཧ޻ֶ‫ڀݚ‬Պ

1 ͸͡Ίʹ
‫ܦ‬೥ྼԽͨ͠ߏ଄෦ࡐͷ݈શ౓ΛඇഁյతʹධՁ
͢ΔͨΊʹɼ௒Ի೾΍ి࣓೾Λ༻͍ͨ‫ߦ͕ࠪݕ‬ΘΕ͍ͯΔɽ͜ΕΒͷ೾ಈ͸ෆՄࢹͰ͋ΔͨΊɼγϛϡϨʔγϣϯΛԉ༻ͯ͠਺஋తʹ఻ൖ‫ܦ‬࿏Λਪఆ͢Δ
͜ͱ͕ߦΘΕ͍ͯΔɽஶऀΒ͸ɼ͜Ε·Ͱʹ༗‫ੵݶ‬෼๏ [1](Finite Integration Technique: FIT) ͱΠϝʔδϕʔεॲཧ [2] Λ૊Έ߹Θͤͨ஄ੑ೾ಈͷγϛϡϨʔγϣϯख๏ [3] ʹ͍ͭͯ‫ڀݚ‬Λߦ͖ͬͯͨɽFIT͸ɼࢧ഑ํఔࣜ (ภඍ෼ํఔࣜ) Ληϧͱ‫ݺ‬͹ΕΔඍখਖ਼ํ‫( ܗ‬2 ࣍‫)ݩ‬or ཱํମ (3 ࣍‫ྖͰ )ݩ‬Ҭੵ෼ͨ͠‫ޙ‬ɼͦͷੵ෼ํఔࣜΛ཭ࢄԽ͢Δํ๏Ͱ͋Δɽ࣌ؒ

͕Ϛ΢ϯτ͞Ε͓ͯΓɼ࠷େͰ 64 ‫ج‬ͷ GPU ͕ฒྻͰར༻Ͱ͖ΔɽGPU1 ‫ͨ͋ج‬Γ 512 ‫ݸ‬ͷฒྻॲཧίΞ͕౥ࡌ͞Ε͍ͯΔͨΊɼ‫ࢉܭ‬ͷΞϧΰϦζϜʹΑͬͯ͸ CPU ฒྻ‫ࢉܭ‬ΑΓ΋ GPU ‫ࢉܭ‬ͷํ͕ɼFLOPS͕େ͖͘ͳΔ͜ͱ͕ใࠂ͞Ε͍ͯΔ [5]ɽฏ੒ 24 ೥౓ͷϓϩάϥϜߴ౓Խࣄ‫Ͱۀ‬͸ɼΠϝʔδϕʔε FIT Λ GPU Ͱฒྻ‫ͨ͠ࢉܭ‬৔߹ͷޮ཰վળʹ͍ͭͯΞυόΠεΛ௖͍ͨɽߴ౓ԽલʹɼCUDA

Fortran Λ༻͍ͨஶऀࣗ਎ͷίʔυʹΑͬͯɼ͋Δఔ౓ͷੑೳ͕ग़‫͕͍͖ͨͯͰݱ‬ɼ෼ੳґཔͷ݁ՌɼҎԼͷΑ͏ͳੑೳ޲্ͷ༨஍͕͋Δ͜ͱ͕෼͔ͬͨɽ

ɾɹ GPU ίʔυͷ࣮ߦ‫( ࣜܗ‬Χʔωϧ) ͓Αͼσόํ޲ʹ͸ɼཅతʹղͷߋ৽Λߦ͏ɽΠϝʔδϕʔεΠεϝϞϦ (GPU ଆͷϝϞϦ) ͷར༻ʹ͍ͭͯॲཧʹΑͬͯɼ2 ࣍‫ݩ‬ղੳͰ͸ϕΫλʔਤ΍σδΧ࠷దԽͷ༨஍͕͋Δ͜ͱɽϝ౳ͷը૾͔Βɼ3 ࣍‫ݩ‬ղੳͰ͸ CADɼCT ࣸਅɼද໘Ԝತσʔλ౳ͷσδλϧը૾͔Β਺஋ϞσϧΛ࡞ɾɹෳ਺ͷ GPU Λ࢖༻͢Δ৔߹ɼMPI ʹΑͬͯσόΠεؒͷϝϞϦసૹΛ͍ͯ͠Δ͕ɼͦͷࡍͷ੒Ͱ͖ɼϓϦϓϩηεͷ࿑ྗΛ‫ݮ‬Β͢͜ͱʹΑͬͯɼ௨৴Λ࠷దԽ͢Δ༨஍͕͋Δ͜ͱɽ

ଟ༷ͳର৅ʹରͯ͠༰қʹγϛϡϨʔγϣϯΛߦ͏͜ͱ͕Ͱ͖Δɽฏ੒ 22 ೥౓ͷϓϩάϥϜߴ౓Խࢧԉɾɹ೾ಈ఻ൖΛՄࢹԽ͢ΔͨΊʹྖҬ಺ͷ෺ཧྔΛίΞʹΑΔେ‫ن‬໛ͳ‫͕ࢉܭ‬ՄೳͱͳͬͨɽͦͷੑೳϝϞϦ͔ΒϗετϝϞϦ (CPU ଆͷϝϞϦʣʹ͸ֶࡍେ‫ن‬໛৘ใ‫ج‬൫‫ڞ‬ಉར༻ɾ‫ڞ‬ಉ‫఺ڌڀݚ‬ͷฏσʔλΛసૹ͢Δඞཁ͕͋Γɼ͜ͷϓϩηεΛ

੒ 23 ೥౓‫ڞ‬ಉ‫͞ূݕ͍͓ͯʹڀݚ‬Εɼ1000 ίΞརӅณ͢ΔͨΊͷํࡦ͕ඞཁͰ͋Δ͜ͱɽ༻࣌Ͱ΋ฒྻԽੑೳ͕ϦχΞʹൃ‫͢ݱ‬Δ͜ͱΛࣔ͜͜͠Ͱ͸ɼҎ্ͷΑ͏ͳཁҼʹର͢Δվળࡦɼ͓Αࣄ‫Ͱۀ‬Πϝʔδϕʔε FIT Λվྑͯ͠௖͖ɼϚϧνͨ [4]ɽ

ҰํɼάϥϑΟοΫϘʔυ্ͷूੵճ࿏Ͱ͋Δ Graph-

ics Processing Unit(GPU) ΛՊֶٕज़‫ʹࢉܭ‬Ԡ༻͢Δ
ಈ͖͕೥ʑՃ଎͍ͯ͠ΔɽGPU Λը૾ॲཧઐ༻ͷϓϩηοαͱ͚ͯͩ͠Ͱ͸ͳ͘ɼฒྻ‫ࢉܭ‬Ϣχοτͱͯ͠ར༻͢ΔࢼΈ͸ General purpose computing on GPU(GPGPU) ͱ‫ݺ‬͹Ε͍ͯΔɽ‫౎ژ‬େֶʹ͓͍ͯ΋ GPUͷ‫ʹͨ৽͕ػࢉܭ‬੔උ͞ΕɼαϒγεςϜ B Ͱ͸ 1ϊʔυ͋ͨΓ 1 ‫ج‬ͷ GPU(NVIDIA ࣾͷ Tesla M2090)

28 − Academic Center for Computing and Media Studies, Kyoto Universityఆ‫ظ‬తʹग़ྗ͍ͯ͠Δ͕ɼͦͷͱ͖ʹσόΠε

ͼͦͷ੒Ռʹ͍ͭͯใࠂ͢Δɽ

2 ‫ػࢉܭ‬ͷ࢓༷ͱ։ൃ‫ڥ؀‬
2.1‫ ͍ͨ༻ʹূݕ‬GPU ‫ػࢉܭ‬ຊ‫Ͱڀݚ‬͸ɼαϒγεςϜ B Λར༻ͨ͠ɽGPUϘʔυͱͯ͠ NVIDIA ࣾͷ Tesla M2090 ͕ 1 ϊʔυ͋ͨΓ 1 ‫ج‬Ϛ΢ϯτ͞Ε͍ͯΔɽM2090 ʹ͸ɼ512‫ݸ‬ͷ CUDA ίΞͱΑ͹ΕΔ࠷খ୯ҐͷԋࢉॲཧϢχοτ͕಺ଂ͞Ε͍ͯΔɽ͜ͷ CUDA ίΞ͕ 32 ‫ݸ‬ू·ͬͯ 1 ͭͷετϦʔϛϯάɾϚϧνϓϩηοα

䜾䝸䝑䝗(SM) Λ‫ܗ‬੒͢ΔɽSM ಺෦ʹ͸ɼߴ଎ϝϞϦͰ͋Δx2 ᪉ྥ

ϨδελʹՃ͑ɼγΣΞʔυϝϞϦɼL1 Ωϟογϡ

1౳ͷΦϯνοϓϝϞϦ (‫ ܭ‬64KB) ͕͋ΔɽGPU ͷϝ

․

‣

‧‧‥‧

x1 ᪉ྥ

65535
2

‣

ϞϦ֊૚͸ɼL2 ΩϟογϡΛ‫ܦ‬༝ͯ͠ɼάϩʔόϧ

․

x3 ᪉ྥ

ϝϞϦɼϩʔΧϧϝϞϦɼςΫενϟɼίϯελϯ
τϝϞϦ౳ͷσʔλΛ‫ه‬Ա͢Δ DRAM(GDDR5) ʹ‧‧‥‧Ͱ͋ΓɼαϒγεςϜ B ্ͷ CPU(DDR3) ͷͦΕͱ䝤䝻䝑䜽‫͍͕ͯͬܨ‬ΔɽϝϞϦόϯυ෯͸ 177GB/s(ެশ஋)ൺ΂Δͱ 1.7 ഒఔ౓େ͖͍ɽैͬͯɼ͜ͷେ͖͍ό

x2 ᪉ྥ
1

ϯυ෯Λ༗ޮར༻͢Δ͜ͱ͕ GPU ‫ࢉܭ‬ͷߴ଎Խͷ 1

2(1)

ͭͷΞϓϩʔνͱͳΔɽ͔͠͠ɼGPU1 ‫ͨ͋ج‬Γͷ

(2)

DRAM ͸ 6GB Ͱ͋Γɼେ‫ن‬໛ͳ‫ࢉܭ‬Λߦ͏৔߹ʹ͸ɼෳ਺ͷ GPU Λ࢖༻͢Δ͜ͱʹͳΔɽຊ‫Ͱڀݚ‬x3 ᪉ྥ͸ɼ࠷େ 8 ‫ج‬ͷ GPU Λ༻͍ͯ‫ূݕ‬Λߦͬͨɽ1

2.2 ։ൃ‫ڥ؀‬

2.2.1 CUDA Fortran

2

䝇䝺䝑䝗
x1 ᪉ྥ

ਤ 1: CUDA Fortran ʹ͓͚Δ‫ࢉܭ‬ͷ֊૚ɽ͜͜Ͱ͸ɼϓϩάϥϛϯά‫ ͯ͠ͱڥ؀‬CUDA For-

tran Λ༻͍ͨɽCUDA Fortran Λఏ‫͢ڙ‬Δͷ͸ɼ‫ࡏݱ‬ॲཧΛ࣮ߦ͢Δ͜ͱʹͳΔ ∗ ɽ

tor Fortran Ver.12.8 for Linux(CUDA Toolkit 4.1) Λ༻

3

FIT ΁ͷ GPU ‫ࢉܭ‬ͷಋೖ

3.1

༗‫ੵݶ‬෼๏ (FIT)ͷͱ͜Ζ PGI ࣾͷΈͰ͋Γɼ͜͜Ͱ͸ PGI Accelera͍ͨɽCUDA Fortran ͷϓϩάϥϜ͸ɼCPU Λಈ࡞ͤ͞Δϗετ༻ίʔυͱ GPU Λಈ࡞ͤ͞ΔσόΠε༻ίʔυͷ 2 ͭͷίʔυ͔Β੒Γཱͭɽϗετ༻ίʔυʹ͸ɼ௨ৗͷϓϩάϥϛϯάͷଞʹ CUDA ͷϥϯλΠϜ APIɼBuilt-in ͷϕΫτϧ‫ܕ‬એ‫ݴ‬ɼ͓ΑͼσόΠε༻ίʔυΛ‫ Ϳݺ‬call จ౳͕‫·ؚ‬Ε͍ͯΔɽ

2.2.2 CUDA ͷ֊૚ߏ଄ਤ 1 ʹࣔ͢Α͏ʹɼCUDA ͸άϦουɼϒϩοΫɼεϨουͱ͍͏ฒྻ‫ࢉܭ‬ͷཻ౓୯ҐΛ͍࣋ͬͯΔɽԼ޲͖ͷ೾ઢ໼ҹͰࣔͨ͠εϨου͕‫ࢉܭ‬Λߦ͏࠷খ୯ҐͰ͋ΓɼͦͷεϨουͷू·Γ͕ϒϩοΫɼϒϩοΫͷू·Γ͕άϦουͱͳΔɽϒϩοΫͱεϨο
υ͸ 3 ࣍‫ݩ‬తͳߏ଄Λ͓ͯ͠ΓɼͦΕͧΕ࠷େ਺͕ܾΊΒΕ͍ͯΔɽਤ 1 ্ଆʹࣔ͢Α͏ʹɼ1 ͭͷάϦου಺ͷ 1 ํ޲ͷϒϩοΫͷ࠷େ਺͸ 65535 ‫Ͱݸ‬͋Δɽ·ͨɼ1 ϒϩοΫ಺Ͱ࢖༻Ͱ͖ΔεϨουͷ࠷Πϝʔδϕʔε FIT ͷৄࡉͳఆࣜԽʹ͍ͭͯ͸ɼա‫ڈ‬ͷ࿦จ [3] Ͱ‫͞ࡌه‬Ε͍ͯΔͷͰɼ͜͜Ͱ͸ུ֓ʹ͍ͭͯड़΂Δɽຊ࿦จͰ͸ɼ2 ࣍‫ݩ‬஄ੑ೾ಈ৔

Λର৅ͱ͠ɼ೾ಈͷ఻ൖํ޲͸ (x1 –x2 )ɼภ޲ํ޲

͸໘֎ (x3 ) ͱ͢Δɽ͜ͷ৔߹ͷ஄ੑ೾͸ SH(Shear

Horizontal) ೾ͱ‫ݺ‬͹ΕΔɽx3 ࣠ํ޲ͷཻࢠ଎౓Λv3 ɼͤΜஅԠྗΛ τ31 ɼτ32 ͱ͓͍ͨͱ͖ɼSH ೾ͷ఻ൖΛࢧ഑͢Δӡಈํఔ͓ࣜΑͼߏ੒ࣜ͸ҎԼͱͳΔɽ∂τ31 (x, t) ∂τ32 (x, t)
+
∂x1
∂x2
∂v3 (x, t)τ̇31 (x, t)
=
µ(x)
∂x1
∂v3 (x, t)τ̇32 (x, t)
=
µ(x)
∂x2

ρ(x)v̇3 (x, t) =

(1)
(2)
(3)

େ਺͸ɼx1 ํ޲͕ 1024 ‫ݸ‬ɼx2 ํ޲͕ 1024 ‫ݸ‬ɼx3

˙ ͸࣌ؒ t ʹؔ͢Δภඍ෼ (∂/∂t) Λද͢ɽ্
͜͜Ͱɼ()

Ϩου਺͸ 1024 ‫͕͋ݶ੍͏͍ͱͰ·ݸ‬Δɽͳ͓ɼϒࣜͰɼρ ͸ࡐྉͷີ౓ɼµ ͸ͤΜஅ஄ੑ܎਺Λද͠ɼ౳√ํੑࡐྉͷ৔߹ʹ͸ԣ೾Ի଎ͱͷؔ܎ࣜ cT = µ/ρํ޲͕ 64 ‫͋Ͱݸ‬Δɽͨͩ͠ɼ1 ϒϩοΫ಺ͷ߹‫ܭ‬εϩοΫ΍εϨουͷ ID ͸ߏ଄ମͰ؅ཧ͞ΕΔɽ

͕੒ཱ͢Δɽ

ϓϩάϥϜ࣮ߦ࣌ʹɼάϦου্ʹΧʔωϧΛੜ੒͠ɼεϨουͦΕͧΕ͕ͦͷΧʔωϧΛ࣮ߦ͢Δɽ
ैͬͯɼGPU ্ͷ֤ϓϩηοα͸ SIMD ‫ܕ‬ͷσʔλ∗ CUDA ʹ͓͚ΔΧʔωϧ͸ if จ౳ͷ෼‫ྩ໋ذ‬΋࣮ߦͰ͖ΔͨΊɼ‫ ͳີݫ‬SIMD ॲཧͰ͸ͳ͘ɼNVIDIA ࣾ͸ SIMT(Single
Instruction Multiple Thread) ͱ‫ݺ‬ΜͰ۠ผͯ͠ද‫͍ͯ͠ه‬Δɽ

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 29๪㡿ᇦ 3཭ࢄԽ‫ޙ‬͸ɼཅతʹղΛߋ৽͢ΔΞϧΰϦζϜͱͳΔɽx1 ํ޲ͷ഑ྻཁૉΛ iɼx2 ํ޲ͷͦΕΛ j ͱ
͢Δͱɼ্ࣜͷํఔࣜ͸࣍ͷΑ͏ʹղ͘͜ͱʹͳΔɽ

V3(1,3)V3 (i, j) ⇐ V3 (i, j) + δ1 [T31 (i + 1, j) − T31 (i, j)]

V3(1,3)x2

V3(2,2)
MPI 通信V3(4,4)V3(i,j)V3(i,j)V3(1,2)

¯∆t
µ̄
ε2 =
∆x

(7)

Ͱ͋Δɽ͜͜Ͱɼ⇐ ͸ӈลͷ‫݁ࢉܭ‬ՌΛࠨลʹ୅ೖx1

V3(3,2)V3(i,j)V3(2,4)V3(3,4)V3(1,1) V3(2,1) V3(3,1) V3(4,1)

µ̄∆t
ε1 =
,
∆x

V3(3,1)V3(1,2)

V3(1,4)͜͜Ͱɼ∆t,
δ1 = δ2 =ρij ∆xV3(i,j)V3(4,1)

MPI 通信(5)(6)

V3(1,1)V3(1,1) V3(2,1) V3(3,1) V3(4,1)

T31 (i, j) ⇐ T31 (i, j) + ε1 [V3 (i, j) − V3 (i − 1, j)]䝯䝰䝸䛾୪䜃

V3(2,1)V3(1,2)(4)

T32 (i, j) ⇐ T32 (i, j) + ε2 [V3 (i, j) − V3 (i, j − 1)]

V3(i,j)

MPI 通信+ δ2 [T32 (i, j + 1) − T32 (i, j)]

MPI 通信V3(1,4)

๪㡿ᇦ 4

๪㡿ᇦ 1

๪㡿ᇦ 2

V3(4,4)

ਤ 2: ࠨɿFIT Ͱ༻͍ΒΕΔ 2 ࣍‫ݩ‬഑ྻͷΠϯσοΫ

͢Δ͜ͱΛࣔ͢ɽ্ࣜͰɼ∆x ͱ ∆t ͸ɼͦΕͧΕηεɼ͓Αͼ෭ྖҬͱକ௨৴ɽӈɿFortran ʹ͓͚Δ 2

ϧ௕͓Αͼ࣌ؒεςοϓࠁΈͰ͋Γɼ·ͨɼµ̄ ͸ x1࣍‫ݩ‬഑ྻͷϝϞϦ্ͷฒͼɽ¯ ͸ x2 ํ޲ʹྡΓ߹͏ηϧʹํ޲ʹྡΓ߹͏ηϧɼµ̄
͓͍ͯ µ ͷฏ‫ۉ‬ԽॲཧΛ࣮ߦͨ͠஋Ͱ͋Δɽ

3.2 εϨουͱ഑ྻ

ਤ 2 ʹࣔ͢Α͏ʹɼ‫ྖࢉܭ‬Ҭ͕େ͖͍৔߹ʹ͸഑ྻΛ෭ྖҬʹ෼഑͢Δ † ɽਤ 1 Ͱ͸ɼϒϩοΫͱεϨου‫ ʹڞ‬3 ࣍‫ݩ‬഑ྻͰઆ໌Λ͕ͨ͠ɼSH ೾఻ൖ໰୊ͷ FIT ͸ 2 ࣍‫ݩ‬ղੳͰ͋ΔͷͰɼ͜͜Ͱ͸ϒϩοΫͱεϨου͸ 2 ࣍‫ݩ‬഑ྻͱ͢Δ (x3 ํ޲ͷϒϩοΫɼεϨου਺͸‫ ʹڞ‬1 ͱ͢Δ)ɽ

ਤ 2 ʹ഑ྻͷׂΓ౰ͯͷ֓ཁΛࣔ͢ɽ͜ͷྫͰ͸ɼ‫ྖࢉܭ‬ҬΛ 4 ͭͷ෭ྖҬʹ෼ׂ͢Δɽ෭ྖҬΛ

1 ͭͷάϦουͱ͠ɼ෭ྖҬ⃝1 Ͱ͸఺ઢͰғΜ্ͩԼͷϒϩοΫʹ෼͔Ε͍ͯΔͱ͢ΔɽάϦου಺෦ͷϒϩοΫͷେ͖͞͸ CUDA Fortran ͷߏ଄ମද‫ه‬Λ༻͍ͯද͞ΕΔɼ͜ͷ৔߹ɼgridDim%x=1, grid-

Dim%y=2 ͱͳΔɽ·ͨɼ1 ͭͷϒϩοΫʹ͸ 8 ͭͷ഑ྻ͕͋Δɽຊ‫Ͱڀݚ‬͸ɼ1 ͭͷ഑ྻͷ‫͑ྫ( ࢉܭ‬͹ࣜ (4) ʹ͓͚Δ V3 (1, 1) ͷߋ৽) ʹ 1 ͭͷεϨο

1 ʹ͓͚Δ 1 ͭͷϒϩουΛׂΓ౰ͯΔɽਤ 2 ͷ⃝

Fortran ͷಛ௃ͱͯ͠ɼਤ 2 ͷӈଆʹࣔ͢Α͏ʹɼ഑ྻ͸ i ͷΠϯσοΫε͕ઌʹಈ͘ॱ൪ͰϝϞϦ্ʹฒͿɽैͬͯɼϒϩοΫ಺ͰϝϞϦ͕࿈ଓ͢Δͷ͸ɼx1 ํ޲ʹฒͿ৔߹Ͱ͋Δ͜ͱʹ஫ҙ͢Δɽ͜Ε͸ɼCUDA Fortran Ͱ΋ಉ༷Ͱ͋Δɽ·ͨɼ࣮ࡍɼεϨουॲཧ͸ 1 ϫʔϓ (32 εϨου) ୯ҐͰߦΘΕΔ
ͨΊɼεϨου਺͸ 32 ͷഒ਺͕ਪ঑͞Ε͍ͯΔɽ

3.3

Ϛϧν GPU ‫ࢉܭ‬͜͜Ͱ͸ɼ1 ͭͷ෭ྖҬͷ‫ࢉܭ‬Λ 1 ͭͷσόΠε

ʹλΠΞοϓ͢Δɽ͢ͳΘͪɼϚϧνϓϩηεʹ෼ׂ‫ޙ‬ɼ1 ͭͷ CPU ϓϩηε͕ 1 ͭͷ෭ྖҬΛ୲౰͢
Δɽͦͷ 1 ͭͷϓϩηεʹ 1 ͭͷ GPU σόΠεΛ

ׂ౰ͯͯɼϚϧν GPU ‫ࢉܭ‬Λߦ͏ɽຊ‫Ͱڀݚ‬͸ɼ2࣍‫ ݩ‬MPI ௨৴Λ༻͍ͯɼϓϩηεؒͷσʔλ௨৴Λߦ͍ͬͯΔɽ2 ࣍‫ ݩ‬MPI ௨৴࣌ʹ͸ɼྡ઀͢Δ෭ྖҬͷؒͰକ෦ (ਤ 2 ʹࣔ͢৭͚ͮͨ͠ଳ෦෼) Λަ‫׵‬͢Δඞཁ͕͋Δɽຊ‫Ͱڀݚ‬͸ɼ্ड़ͷΑ͏ʹ 1 ϓϩηεͱ 1 σόΠεΛλΠΞοϓ͍ͤͯ͞ΔͨΊɼ·ΫͰ͸ɼx1 ํ޲ͷεϨου਺͸ 4ɼx2 ํ޲ͷεϨͣϗετଆͷϝϞϦʹσόΠεଆ͔ΒσʔλΛసૹ

Dim%y=2 ͱ‫ه‬ड़͞ΕΔɽैͬͯɼ1 ͭͷάϦουͰ͸ɼx1 ํ޲ʹ (gridDim%x)×(blockDim%x) ‫ݸ‬ɼx2 ํόΠεϝϞϦʹૹΔͱ͍͏ஈ֊Λ౿ΜͰ͍Δɽ͜͜ου਺͸ 2 ͱͳΔɽ͜Ε͸ɼblockDim%x=4, block-͠ɼ࣍ʹϗετϝϞϦΛհͯ͠ަ‫ͨ͠׵‬σʔλΛσ޲ʹ (gridDim%y)×(blockDim%y) ‫ݸ‬ͷཻ౓ͰεϨο෼ͷΈΛ௨৴ʹ༻͍͍ͯΔɽͳ͓ɼ௨৴ϥΠϒϥϦ

υ͕࡞ΒΕΔ͜ͱʹͳΔɽ

(Ver. 1.8.1) Λ༻͍ͨɽਤ 3 ʹϓϩάϥϜߴ౓Խલͷ† GPU ͸ϝϞϦόϯυ෯͕େ͖͍ͷͰɼ1 ͭͷσόΠεͰ‫ࢉܭ‬ՄೳͳϝϞϦྔͰ͋Δ৔߹͸ɼ෭ྖҬʹ෼ׂ͢Δ͜ͱ͸ආ͚ͨํ͕Α͍ɽ

30 − Academic Center for Computing and Media Studies, Kyoto UniversityͰ͸ɼͰ͖Δ͚ͩ࠷খ‫ݶ‬ͷ௨৴ͱͳΔΑ͏ʹɼV3 ੒
ͱͯ͠ MPI-1 ͱ MPI-2 Λαϙʔτ͢Δ MVAPICH2

ϑϩʔνϟʔτΛࣔ͢ɽͳ͓ɼV3 ͸਺஋ੵ෼ʹΑͬετؒͰίϐʔ͍ͯͨ͠ɽ͜ͷཧ༝͸ɼ෭ྖҬ⃝
1 ͱ

㛤ጞ

⃝
3 ͷ௨৴ͷ৔߹ʹ͸ɼକ෦ͷϝϞϦ͕࿈ଓ͢ΔͨΊσόΠε-ϗετؒͷ௨৴͸ίΞϨευɾΞΫηε͕

V3ࡢィ⟬

1 ͱ⃝
2 ͷ௨৴Ͱ͸ϝϞϦ͕ෆ
ՄೳͰ͋Δ͕ɼ෭ྖҬ⃝࿈ଓͱͳͬͯ͠·͏ͨΊɼ͜ΕΛආ͚ΔͨΊʹ୯७

MPI㏻ಙฎ⌮ʹશ෦Λίϐʔ͍ͯͨ͠ɽ͔͠͠ɼࠓճɼGPU ্ͷෆ࿈ଓϝϞϦΛసૹ༻όοϑΝʹߴ଎ʹ٧ΊࠐΉ৽

V3ࡢ
ࢹࣂ࢖ࢫ-࣍ࢫࢺ㌿㏦

T31ࡢィ⟬

͍ͯɼ8GPU ͷฒྻ࣌Ͱ໿ 1.5 ഒͷߴ଎Խ͕ୡ੒Ͱ
͖ͨɽ

U3ࡢィ⟬
T32ࡢィ⟬

͍͠ΧʔωϧΛ࡞੒ͨ͠ɽ͜ͷ݁ՌɼMPI సૹʹ͓

4.3݁Ռग़ྗ࣌ͷσόΠεϝϞϦ͔ΒϗετϝϞ
Ϧ΁ͷసૹ࣌ؒͷӅณࣇ࢓࢖ࣝฟຊ

4.1 અͷΧʔωϧ࠷దԽʹΑͬͯɼ཯଎ஈ֊͸ΧʔωϧͰ͸ͳ͘ͳΓɼϗετॲཧͰ͋ΔϑΝΠϧग़ྗYesͷํ͕࣌ؒΛཁ͢ΔΑ͏ʹͳͬͯ͠·ͬͨɽߴ౓Խલ

8SGDWH㸽͸ɼਤ 3 ʹࣔ͢Α͏ʹɼϑΝΠϧग़ྗͷલʹ V3 ͸ U3

ࠉNo
⤊஢
GPUฎ⌮

GPU-CPU㏻ಙ

ʹม‫͓ͯ͠׵‬Γɼ͜ΕΒ͸͢΂ͯϗετॲཧ͞Εͯ
CPUฎ⌮͍ͨɽ͜͜Ͱ͸ɼਤ 4 ʹࣔ͢Α͏ʹ U3 ͷม‫׵‬ΛΧʔωϧ্Ͱߦ͍ɼϑΝΠϧग़ྗΛ࣍ͷ࣌ؒεςοϓʹ

ਤ 3: ϓϩάϥϜߴ౓ԽલͷϑϩʔνϟʔτɽͯมҐ U3 ʹม‫ͯ͠׵‬ɼϑΝΠϧग़ྗ͍ͯ͠Δɽ

͢Δ͜ͱͰɼV3 ͷ‫ͱࢉܭ‬Φʔόʔϥοϓ͢ΔΑ͏ʹͨ͠ɽ͜ͷ݁Ռɼߴ౓Խલͷίʔυʹൺ΂ͯɼ1GPU‫Ͱ࣌ࢉܭ‬໿ 3.2 ഒͷߴ଎Խ͕ୡ੒Ͱ͖ͨɽ
㛤ጞ

4 ߴ଎Խɾޮ཰ԽͷͨΊͷվྑ
4.1 Χʔωϧͷ࠷దԽ͓ΑͼσόΠεϝϞϦͷར༻ํ๏ΦϦδφϧͷίʔυͰ͸‫ڞ‬༗ϝϞϦΛଟ༻͍ͯ͠

๓ࡢ᫬㛫ࢫࢸࢵࣉࡢ
ࣇ࢓࢖ࣝฟຊ

V3ࡢィ⟬͕ͨɼNVIDIA ͷ Fermi ΞʔΩςΫνϟͰ͸ɼΞΫ

MPI㏻ಙฎ⌮ηεύλʔϯΛେ෯ʹվળͨ͠ΓɼσʔλΛ࠶ར༻͢ΔͷͰͳ͚Ε͹ɼ‫ڞ‬༗ϝϞϦΛ࢖༻͢ΔϝϦοτ
͸গͳ͍ɽैͬͯɼ͜͜Ͱ͸‫ڞ‬༗ϝϞϦͷ࢖༻Λഇࢭͨ͠ɽ·ͨɼഇࢭʹ൐͏σόΠεϝϞϦΞΫηε

T31ࡢィ⟬

U3ࡢィ⟬

T32ࡢィ⟬

V3ࡢ
ࢹࣂ࢖ࢫ-࣍ࢫࢺ㌿㏦ͷόϯΫίϯϑϦΫτ͸ൃੜ͍ͯ͠ͳ͍͜ͱΛ֬ೝͨ͠ɽҰํɼ‫ڞ‬༗ϝϞϦ্ʹอ༗͍ͯͨ͠ഔ࣭ͷࡐྉఆ਺౳͸ɼίϯελϯτϝϞϦ্ʹ഑ஔͨ͠ɽ·ͨɼFortran ͸Ϣʔβఆٛͷߏ଄ମ͕࢖༻Ͱ͖Δ͕ɼ

YesίΞϨευɾΞΫηεΛ๦͍͛ͯͨͨΊɼ͜ΕΛഇࢭͨ͠ɽ͜ͷ݁ՌɼΦϦδφϧίʔυʹൺ΂ͯ໿ 1.5

8SGDWH㸽ࠉNo⤊஢ഒͷΧʔωϧͷߴ଎ԽΛୡ੒Ͱ͖ͨɽ

4.2 Ϛϧν GPU ࣌ͷσόΠεؒͷϝϞϦసૹਤ 2 ʹࣔ͢Α͏ʹɼFIT ͷ৔߹ɼMPI సૹͷର৅

ͱͳΔ഑ྻ͸෭ྖҬͷକ෦ 1 ૚͚ͩͰྑ͍ɽ͔͠͠ɼ

ߴ౓ԽલͷϓϩάϥϜͰ͸ɼ഑ྻશମΛσόΠε-ϗ

GPUฎ⌮GPU-CPU㏻ಙ

CPUฎ⌮ਤ 4: ϓϩάϥϜߴ౓Խ‫ޙ‬ͷϑϩʔνϟʔτ (ϑΝΠϧग़ྗ࣌ؒΛӅณॲཧʹΑͬͯ୹ॖ͍ͯ͠Δ)ɽ

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 31

5 ඇഁյ‫΁ࠪݕ‬ͷԠ༻

6 ·ͱΊνϡʔχϯάͯ͠௖͍ͨίʔυΛ༻͍ͯɼඇഁյ

ࠓճɼΧʔωϧͷ࠷దԽɼMPI ௨৴ͷޮ཰Խʹͭ‫ࠪݕ‬ͷܽؕ࠶ߏ੒γϛϡϨʔγϣϯΛ࣮ࢪͨ͠ɽ͜Ε͍ͯ͝ࢦಋΛ௖͍ͨɽվળͷաఔͰɼ཯଎ஈ֊͕Χʔ͸શ೾‫ܗ‬αϯϓϦϯάॲཧ (FSAP) ํࣜ [6] ʹΑΔ։ωϧͰ͸ͳ͘ɼϗετॲཧʢϑΝΠϧग़ྗʣ΁ͱγ‫߹ޱ‬੒๏Λ໛ٖͨ͠΋ͷͰ͋ΓɼΞϨΠ୳৮ࢠΛ༻͍ͯܽؕΛ࠶ߏ੒͢Δɽਤ 5 ͷ্෦ʹࣔ͢Α͏ͳܽؕ഑ஔΛ૝ఆ͠ɼࢼ‫ݧ‬ମ্෦ʹ 32 ૉࢠͷΞϨΠ୳৮ࢠΛઃஔͨ͠ɽFSAP ํࣜͷγϛϡϨʔγϣϯͰ͸ɼϑτͨ͠ɽैͬͯɼFIT ͷ‫ࢉܭ‬աఔΛ‫ݟ‬௚͠ɼϑΝΠϧग़ྗΛӅณ͢Δ͜ͱΛࢼΈͨɽ͜ΕΒͷ݁Ռɼ౰ॳͷίʔυʹൺ΂ͯେ෯ͳߴ଎Խ͕࣮‫͖ͨͰݱ‬ɽ͜ͷ੒ՌΛඇഁյ‫ʹࠪݕ‬Ԡ༻͠ɼܽؕͷ࠶ߏ੒γϛϡૉࢠ਺෼ (32 ճ) ͷ௒Ի೾ͷૹ৴Λߦ͏ඞཁ͕͋ΔɽϨʔγϣϯ͕࣮࣌ؒͰ࣮ࢪͰ͖Δ͜ͱΛࣔͨ͠ɽεύ͖ 5000 εςοϓͰ͋Δɽղੳର৅͕ PC Ͱ΋े෼࣮͖ͨ͜ͱ͸ඇৗʹҙ͕ٛ͋Δɽࠓ‫ޙ‬ͷ՝୊͸ɼ MPIϞσϧͷେ͖͞͸ 200 ສཁૉɼ࣌ؒߋ৽͸֤ճʹͭߦՄೳͳେ͖͞Ͱ͋ͬͨͷͰɼ͜͜Ͱ͸ Windows ൛ίϯͷߴ౓ԽͰಘΒΕٕͨज़͕ɼPC ϨϕϧʹҠసͰసૹ࣌ؒͷӅณɼແ൓ࣹ‫ڥ‬ք૚ (PML) ͷޮ཰Խɼ͞ͷ PGI ίϯύΠϥͱ NVIDIA GeForce GTX680 Λ༻Βʹ͸ 3 ࣍‫ݩ‬େ‫ن‬໛໰୊΁ͷԠ༻ΛࢼΈ͍͖͍ͯͨɽ

ՌΛࣔ͢ɽܽؕͷҐஔ͓Αͼେ͖͞Λਫ਼౓ྑ͘࠶‫ݱ‬ँࣙɿ ༗ӹͳ͝ॿ‫ݴ‬Λ௖͍ͨʮεʔύʔίϯϐϡʔ

Ͱ͖͍ͯΔɽͳ͓ɼ࠶ߏ੒ʹཁͨ͠‫ؒ࣌ࢉܭ‬͸໿ 10λγεςϜ‫ڞ‬ಉ‫اڀݚ‬ըҕһձʯͷօ༷ɼฒͼʹɼ෼Ͱ͋Γɼ͜Ε·Ͱͷඇഁյ‫ࠪݕ‬γϛϡϨʔγϣϯνϡʔχϯάίʔυΛఏ‫ͯ͠ڙ‬௖͍ͨΫϨΠɾδϟπʔϧͱͯ͠͸ྨΛ‫ߴͲ΄͍ͳݟ‬଎ԽΛ࣮‫͖ͨͰݱ‬ɽύϯɾΠϯΫͷ෢ాେีࢯʹ‫ँײ‬Λਃ্͛͠·͢ɽ͍ͯ‫ূݕ‬Λߦͬͨɽਤ 5 ͷԼ෦ʹɼܽؕͷ࠶ߏ੒݁ࣇ࢙࣮ࢬࢻ࢔ࣞ࢖᥈ゐᏊࢀߟจ‫ݙ‬ ⣲Ꮚ㸪ࣆࢵࢳ PP[1] P. Fellinger, R. Marklein, K.J. Langenberg and S. Klaholz, “Numerical modeling of elastic wave propagationand scattering with EFIT –elastodynamic ﬁnite integration technique”, Wave Motion, 21, pp.47-66 (1995).

32mm
Ḟ㝗

2mm

42mm

[2] K. Terada, T. Miura and N. Kikuchi, “Digital imagebased modeling applied to the homogenization analysis of composite materials”, Computational Mechanics,
20(4), pp.331-346 (1997).
[3] K. Nakahata, J. Tokunaga, K. Kimoto and S. Hirose,“Large scale simulation of ultrasonic wave propagation in
concrete using parallelized EFIT”, Journal of Solid Mechanics and Materials Engineering, 2(11), pp.1462-1469(2008).
[4] ֽࢁ࿨உɼ‫ڇ‬ౡলɼࣉా‫ݡ‬ೋ࿠ɼԬᖒॏ৴ɼ໦ଜҰ࿠ɼதാ࿨೭ɼઙҪޫًɼদຊ७Ұɼ‫ؠ‬Լ෢࢙ɼখࢁాߞ
ೋ: ϚϧνϑΟδοΫε͓Αͼ࠷దԽ໰୊ʹ޲͚ͨϋΠύϑΥʔϚϯε‫ֶྗࢉܭ‬, ֶࡍେ‫ن‬໛৘ใ‫ج‬൫‫ڞ‬ಉར༻ɾ‫ڞ‬ಉ‫ ఺ڌڀݚ‬ฏ੒ 23 ೥౓‫ڞ‬ಉ‫ ڀݚ‬தؒใࠂॻ,11-NA04 (2011).
[5] T. Okamoto, H. Takenaka, T. Nakamura and T. Aoki,“Accelerating large-scale simulation of seismic wave
propagation by multi-GPUs and three-dimensional domain decomposition”, Earth, Planets and Space, 62,
pp.939-942 (2010).

ਤ 5: FSAP ํࣜʹΑΔ 5 ͭ݀ͷ࠶ߏ੒݁Ռ (FIT ʹ

ΑΔ‫ࢉܭ‬೾‫ܗ‬Λར༻)ɽ

32 − Academic Center for Computing and Media Studies, Kyoto University[6] தാ࿨೭, ஛಺༟޾, ໦ຊ࿨ࢤ, ࣼ֯୳ই༻ΞϨΠ୳৮ࢠΛ༻͍ͨ FSAP ํࣜʹΑΔ಺෦ܽؕͷ௒Ի೾ը૾Խ,೔ຊ‫ػ‬ցֶձ࿦จू, A ฤ, 77 ‫ר‬, 779 ߸, pp.106-116(2011).

プログラム高度化支援事業研究報告

人工島／埋立地等の耐震性評価を目的とした
動的／静的水～土骨格連成有限変形解析コードの高度化
野田利弘
名古屋大学 減災連携研究センター

1 研究背景と本事業の目的

しい状況となっていた。
このような中、京都大学学術情報メディアセン

日本の戦後の経済発展は、東京、名古屋、阪神

ターの平成 24 年度プログラム高度化支援事業の地域などの大規模太平洋沿岸地域に生産活動に依追加募集に応募し、10 月に採択をいただいた。当拠してきたところが大きい。この間、当該沿岸地

該事業（共同研究）では、第一段階として、上記域では、軟弱な沖積地盤上に浚渫土などを利用しの課題を有する当該解析コードを、MPI またはて埋立地や人工島を大規模に造成してきた。とこOpenMP と MPI とのハイブリッド並列化を図る

ろがこのような人工地盤は、兵庫県南部地震や東とともに、他のソルバの利用も検討に入れて、大北太平洋沖地震などの地震のほかには、巨大地震

規模・高速化計算が実現できるよう高度化・高性に見舞われたことがない。

能化を図ることであった。このため主として、当名古屋大学地盤力学グループでは、これまでに

該解析コードにおける、接線剛性マトリクス作成

開発してきた動的／静的水～土骨格連成有限変形

ルーチンにおける OpenMP と MPI とのハイブリ○
R

解析コード（GEOASIA [1], [2]）を用いて、上ッド並列化を実現するとともに、ソルバについて記の粘土や砂からなる自然堆積地盤とその間の中直接法として MUMPS [3]を、反復法として

間土からなる人工島・埋立地や、地盤と土構造物

PETSc[4] などを組み込み、スーパーコンピューシステムについて、地震中・地震後の二次元・三タを用いてそれらの効果を、主に三次元解析と有次元地震応答解析を実施し、巨大地震発生時の耐限要素数が多い二次元解析について検証した。

震性の再評価や問題点・課題の抽出等を行ってきた。具体的には、自然堆積（沖積・洪積）粘土～
砂互層地盤、および人工島をモデル化した後、予想される巨大地震を入力して、自然・人工地盤とその上の護岸や盛土、建物基礎、各種埋設施設の数秒から数分までの地震中およびその後の数日から数十年までの地震後の耐震性を解析的に評価してきた（例えば、図１）。
しかしながら、当該解析コードは高速化を実現するための並列化が一部で OpenMP を用いて施されているだけで、また（パソコン／WS での計算実現を意図していたため）使用できるソルバがPARDISO に限定しているなど、チューニングが十分でなかった。このため、大規模な解析領域を

図１ 軟弱地盤上の構造物の三次元応答解析事例

対象にした高速計算を常態的に実現することが難（せん断ひずみ、右上：中央付近の拡大図）

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 33

2 動的／静的水～土骨格連成有限変形
解析コードの概要

方程式を解くことになるが、材料的・幾何学的の両非線形であるためイタラティブに処理する。このとき、これまでは、計算機はパソコン／WS で、

一般に、土は土粒子が構成する土骨格とその間ソルバには PARDISO を用いて解を求めている。

隙に水と空気が存在する、すなわち固体、液体お⑤地震中の計算時間間隔 1/1000 程度～地震後

よび気体からなる三相混合体である。間隙に空気

100 年程度までの圧密変形のための計算時間間隔がない場合、あるいは実質上空気がないとみなさ

1 年程度と広範な時間帯を対象に計算する。

れる場合は、間隙に水だけが存在することになるが、このような土を飽和土と呼び、空気がある場
合を不飽和土と呼ぶ。当該解析コードは、固液二相混合体理論（例えば、[5], [6]）に基づいて、こ3 直接法ソルバ MUMPS または反復法
ソルバ PETSc を用いたハイブリッド並
列化の検討の飽和土からなる地盤が盛土や地震外力などの荷重を受けた場合に生じる変形や破壊の力学挙動をソルバの PARDISO の部分を、プロセス並列が求める。当該解析コードには、土骨格の構成式と可能な直接法ソルバ MUMPS およびに反復法ソして報告者らが開発してきた弾塑性構成式（SYS

ルバ PETSc に置き換え、係数行列作成サブルー

Cam-clay model [7]）を搭載し、同一の理論（弾チン（以下 KGRMAT と書く）の OpenMP と塑性理論）のもとで、砂から粘土までその間に稠

MPI のハイブリッド並列化を行い、速度向上性を密に存在するあらゆる土をすべて計算対象にする

含め計算結果がオリジナルのコードと一致するか（All Soils）。また、変形を通じて破壊へ至る進行

どうかの動作検証を行った。ここで、直接法ソル

性破壊や液状化していた砂地盤がやがて示す圧密バの中から MUMPS を用いたのは非対称行列を変形など、地盤の安定状態から不安定状態へ、そ対象にした MPI 並列が可能で、継続利用性の観して不安定状態からまたつぎの安定状態へ、時間点から入手が比較的容易であったためである。ま的にも空間的にも連続して解く（All States）。さた、PETSc を用いたのは、MPI 並列可能性など

らにここでは、外力形態として静的問題と動的問

加えて、非常に多くのユーザコミュニティがあり、題を区別せずに、どちらの問題も解析を連続してドキュメント、例題なども非常に充実しているか

進める（All Round）
。

らである。

当該解析コードの定式化または求解における主

図２は、図１とほぼ同等の有限要素メッシュをな特徴は、次の点である。①混合体である飽和土用いて検証した三次元解析の結果を示す。なお、の運動方程式を、飽和土からの間隙水の出入りに

二次元解析についても検証を行っているが、省略

伴う「体積変化拘束式」の下で解く必要があるたする。計算結果はオリジナルの解析コードの結果

め、一相系材料を扱う場合よりも解くべき元が多と一致し、まずスレッド並列だけでなくプロセスい。②土骨格の弾塑性的特性、つまり材料非線形並列が可能であることを確認した。全体の計算速性を考慮するとともに、地盤や土構造物の刻々の度においては、
シリアル（p1t1）
比から見て PETSc

幾何的形状変化、つまり幾何的非線形性をによって MUMPS よりもやや高速化が図られる

updated Lagrangian 法を用いて考慮するため、

一方、KGRMAT においては OpenMP 並列によ運動方程式を時間微分し、躍度（＝加加速度＝加る計算時間短縮（ただし、PETSc のスレッド並列

速度の一階時間微分）を有する速度型の運動方程効果はない）が図られ、並列度が 16（p1t16）程式を用いる。③空間的には有限要素法を用いて展度であればPARDISOと同様の計算時間になるこ

開するとともに、時間積分においては躍度が線形とが分かった。また、PETSc を用いた場合におい的に変化すると仮定しつつ、Wilson のθ法[8]をては、全体の計算時間のうち KGRMAT と PETSc

拡大的に用いて、位置ベクトルの時間に関する三以外の「その他」の部分が占める比率が、並列度

階の連立常微分方程式を連立一次方程式に展開し、128（p8t16）で 80/96 であることがわかり、高速①の体積変化拘束式と連立して解く。④最終的に計算を実施する上で、
KGRMAT などの MPI 通信

は、非対称疎行列を係数行列とするこの連立一次量削減や、ソルバと KGRMAT 以外の箇所の並列

34 − Academic Center for Computing and Media Studies, Kyoto University

化等によるチューンが必要なことも明らかになっ5 参考文献た。
なお、他の解析事例から MUMPS で可能な計算が PETSc では実現できない場合もあり、行列[1] Asaoka, A. and Noda, T. (2007): All soils all states

all round geo-analysis integration, Int. Workshop on

Constitutive Modelling - Development, Implementation,の性質に関する課題も明らかになった。

Evaluation, and Application, Hong Kong, China,
11-27.

[2] Noda, T, Asaoka, A. and Nakano, M.(2008):
Soil-water coupled finite deformation analysis based ona rate-type equation of motion incorporating the SYS
Cam-clay model, Soils and Foundations, 48(.6),

771-790.

[3] MUMPS: MUltifrontal Massively Parallel sparse
direct Solver, http://graal.ens-lyon.fr/MUMPS/index.
php[4] PETSc: Portable, Extensible Toolkit for Scientific
Computation, http://www.mcs.anl.gov/petsc/

図２ 動作検証結果[5] de Boer, R.(1998): Theory of porous media -past

4 おわりにand present-, ZAMM - Journal of Applied Mathematics

本事業の採択から半年足らずの間に、OpenMPと MPI によるハイブリッド並列化が可能である
ことがわかる一方、課題がいくつかあることがわ
かってきた。大規模計算を実現するためには、その解決とともに領域分割法の採用などな抜本的な検討も必要である。幸いなことに、H25 年度の支援事業にも採択いただき、ソルバ・KGRMAT 以外の部分の並列化なども継続的に進んでいる。
今後、大規模計算が可能となれば、巨大地震を経験していない戦後の埋立地・人工島や数多くのand Mechanics, 78, 441-466.

[6] 西村直志（1999）：第 3 章土質力学, 地盤工学
ハンドブック, 地盤工学会, 51-64.

[7] Asaoka, A., Noda, T., Yamada, E., Kaneda, K. andNakano, M. (2002): An elasto-plastic description of two

distinct volume change mechanisms of soils, Soils andFoundations, 42(5), 47-57.

[8] Wilson, E. L., Farhoomand, I. and Bathe, K. -J.
(1973): Nonlinear dynamic analysis of complex

structures, Earthquake Engineering and Structural
Dynamics, 1, 241-252.

土構造物などにおける巨大地震発生時の問題点・
課題の抽出や耐震・防災対策の代替案の検討に繋がるだけでなく、当該学問分野の高度化にも貢献できるものと信じている。
謝辞 京都大学学術情報メディアセンターの牛島省教授には本事業をご紹介いただくとともに、本事業に直接携わっていただいた同センターの中島浩教授・岩下武史准教授とクレイ・ジャパン・イ
ンクの方々には本報には現れてこない事項も含め、
今もなお、数多く検討いただいている。ここに謝意を表します。
なお、スパコン使用にあたっては、JSPS 科研
費 20126012 の助成を受けた。

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 35

スーパーコンピュータ共同研究制度（大規模計算支援枠）研究報告

パーティクルフィルタを用いた分布型降雨流出シミュレーションの高度化と不確実性評価Seong Jin NOH１, 立川康人 2
Korea Institute of Construction Technology, 2 京都大学大学院工学研究科 社会基盤工学専攻

1

高精度な洪水予測は、洪水被害を軽減するために不可欠である。この研究では、パーティクル
フィルタという、モデルの形式と水文量の確率分布の種類に限らずに適用できるデータ同化手法
を流出モデルに導入して水文予測の高精度化を図るとともに、流出モデルの不確実性を評価する。
特に多数のモデルを同時に実行する計算機能力が要求されるパーティクルフィルタを適用するため、超高速、かつ並列処理が可能なスーパーコンピュータを駆使する。パーティクルフィルタで分布型水文モデルの状態量を更新する方法を、桂川流域(面積：887 km2)の洪水予測に適用し、予測精度の改善を確認した。また、不確かさを評価するために、レーダー雨量の推定に用いるパラメータ値と水文モデルのパラメータ値の両方を、水文モデルの状態量とともに同時に同化し、統合的に入力データと予測モデルの不確実性を推定した。スーパーコンピューティング環境を駆使した新たな降雨流出予測システムの開発の適用性を確認することができた。今後、分布型流出モデルおよびパーティクルフィルタを基本とする高精度の実時間流出予測が期待される。

1

概要

して、モデルの形式と水文量の確率分布形によらず適用できる手法で、観測情報や予測システムの治水施設を設計したり、リアルタイムで河川流不確かさを評価しつつ精度の高い予測情報を提供

量や水位を予測したりするためには、降雨量からする極めて有望な手法であるが、多数のモデルを

河川の流量や水位を求める必要がある。それに関同時に実行する計算機能力が要求される。
そこで、する水循環の数値シミュレーションシステムでは、この研究では、超高速、かつ並列処理が可能なス

自然現象の単純化、
降雨・河川流量の測定誤差等、ーパーコンピュータを駆使することで、従来にな様々な要因に基づく不確実性が避けられず、予測い分布型降雨流出シミュレーションシステムの高情報の提供においては、予測値とともにその不確度化を実現し不確実性を評価する。パーティクル

実性を評価することが極めて重要である。特に洪フィルタで分布型水文モデルの状態量を更新する

水は世界で最も一般的な自然災害であるが、様々方法を、桂川流域（面積：887 km2）の洪水予測な不確実性のために、必要な予測時間のもとでのに適用して予測精度の改善が確認する。分布型流

正確な洪水予測は十分になされていない。

出モデルにおいて、水文プロセスの異なる時間ス

データ同化に基づくアンサンブル予報技術は、ケールに起因する遅延応答に対処するために、ラ

明示的に水文モデルによる予測と観測の両方の不グ(lag)フィルタリングアプローチを適用する。ま

確実性を定量化し、 新しい観測値を得るたびに、た、不確かさを評価するため、レーダー雨量の推水文モデルの状態量（ときにはパラメータ値も同定と水文モデルの両方のパラメータを、水文モデ時に）を逐次に更新する方法である。パーティク

ルの状態量とともに同時に同化し、統合的に入力

ルフィルタは非線形システムのデータ同化手法とデータと予測モデルの不確実性を推定する。

36 − Academic Center for Computing and Media Studies, Kyoto University

2

同化手法

こで、逐次インポータンスリサンプリング（SIR）
パーティクルフィルタでは、各時刻、重みに比例

2.1 ベイジアンフィルタリング

する確率で粒子の復元抽出を行い、現時刻のフィルタ分析に従う粒子群を得る。

次のような動的システムを考える。

xk  f ( xk 1, uk )  k k ~ N (0,Wk ) (1)yk  h( xk )  kvk ~ N (0,Vk )(2)
ここで、 xk と yk は時刻 k におけるシステムと観測状態を示す nx 次元のベクトルである。 f は、
システム方程式、 h は観測方程式である。 k と
 k はそれぞれモデル誤差と観測誤差、Wk とVk

2.3 ラグ(lag)フィルタリングアプローチ
流域では、降雨、浸透、地表面水、地下水、河川の流れなど、異なる水文プロセスが、異なる時間スケールで動いている。また、河川流量は降水が地表面に到達して浸透し、地中流・地表面流がは、誤差の共分散行列である。ベイジアンフィル

発生して流下する時間ラグを持つ状態量である。

タリング理論では、システムおよび観測モデルが本研究では、水文プロセスの異なる時間スケール

非線形、非ガウスの場合、解析的に観測ベクトルに起因する遅延応答に対処するための枠組みとし

y1:k = {y i , i = 1, , k} における現在の状態 xk のて、ラグ(lag)フィルタリングアプローチを次の手

distribution)を求めることが不可能である。この

1) 現時刻 k の予測を得るため、過去の時刻
k  j から計算を開始する。

事後確率密度分布(posterior probability density

場合には、パーティクルフィルタなどのシミュレ順で適用する。

2) 時刻 k  j  1 で、システム状態変数を保存

ーションベースのアンサンブル方法を使用して、事後分布を近似的に求めることができる。

する。

2.2 パーティクルフィルタ

各粒子の重みをパーティクルフィルタの方法で推3) 粒子が現在の時間ステップ k に到着したら、定する。

パーティクルフィルタおよび逐次モンテカルロ
法とは動的システムの状態やパラメータ推定問題
を状態空間の中の多数の粒子の数値計算により近似的に解く方法である。現時刻までの状態系列の事後分布は、各粒子に適切な重み値を与えること
で重み付き粒子群が目的分布の近似となるようにする方法である。
p ( x k | y1:k ) 

ここで、 xnw  (x  x )
i
i
k

k

i
k

(3)

wk は、それぞれ i 番目の状態（"
粒子"）と、その重みで、 () は、ディラックデル
i
k と

i 1

4) 推定した重みに従って、粒子をリサンプリン
グする。リサンプリング中に、時刻 k  j  1 の粒子を現時刻 k の重みに従ってリサンプリンする。
更新された時刻 k  j  1 の状態量は次の時刻

k  1 の予測値を計算するとき（ステップ 1)）に、予測計算の初期条件として使用する。

3

状態量更新による洪水予測

3.1 シミュレーションの概要

タ関数である。重みは時間ごとに更新され、前時刻の重み wki 1 と現時刻の観測データ yk を観測モデルに代入して各粒子の尤度 p( yk | xki ) に比例す

る重み値を得る。

wki  wki 1 p( yk | xki )開発したデータ同化法を桂川流域(面積：887km2)の洪水予測に適用した。日吉ダムからの放流
量データは対象流域の上端からの入力データとし(4)て用い、水文モデルではダム下流の降雨流出を予上述した逐次インポータンスサンプリング

測の対象とする。使用した水文モデルは物理的分（SIS）アルゴリズムがパーティクルフィルタの布型流出モデルである WEP (water and energy

基礎をなす方法である。SIS アルゴリズムには、transter processes)モデルである。WEP の状態量

縮退現象(sample degeneracy)という、数回の反復は、土壌水分量、表面流出、地下水位、熱フラッ

計算後に、重みが 1 つの粒子に集中され、他の粒クス成分、斜面や河川の流出などである。水文モ

子の重みが 0 になる問題点が指摘されている。そデルの空間分解能は 250m であり、時間更新間隔
1 時間とした。アンサンブルシミュレーションは京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 37

オープンメッセージパッシング·インタフェース（openMPI）の並列計算技術を適用して京都大学

4 入力データとモデルパラメータの不確実性評価のスーパーコンピュータで実行した。

4.1 シミュレーションの概要
3.2 シミュレーションの結果

ここでは、レーダー雨量変換方程式と水文モデ

表 1 はシミュレーションの期間と観測結果を、

ルのパラメータを、パーティクルフィルタにより

図 1 は 6 時間先の流量予測結果を示している。降同時に推定する。適用した流域は、前節で説明し

雨データは観測値を予測値として用い、完全に降た桂川流域の日吉ダムの上流域で、六甲山に位置

雨予測がなされた状態を仮定している。データ同する空間分解能 250m、観測半径 80km の X-MP化を適用していない場合(Open loop)は、大幅に河レーダーのデータを使用する(図 3)。このレーダー

川流量を過小評価しているが、パーティクルフィ

は、水平および垂直方向の 2 つの偏波（マルチパ

ルタの場合、予測値が観測値とほぼ一致し、デー

ラメータ, MP）を使用して、次の式で降雨量を推タ同化によって、水文モデルの状態量が適切に更定する。

新されたことが分かった。図 2 は予測時間ごとの予測シミュレーションの二乗誤差(RMSE)を示している。パーティクルフィルタによる更新の効果

R  B  Z HH b

KDP < KDP_thres(5)

R  C  K DP c

KDP > KDP_thres(6)が 24 時間先まで継続していることが確認できる。

ここで、 Z HH と K DP は、X-MP レーダーから観表１．シミュレーションの期間と観測結果(6)で降雨強度に変換される。B、b、C、c は、通測される観測値であり、 K DP_thres によって式(5)、観測

シミュレーションの期

最大流量

降水量

常、地上観測によって較正されるパラメータであ地点

間(m3s-1)(mm)

る。
本研究では、
４つのレーダパラメータの中で、

桂2003 年 8 月 1－31 日

361

258

B と C の２つのパラメータ値を、水文モデルのパラメータと同時に推定する。水文モデルは、前のセクションで使用した WEP モデルより簡単な構
造の分布型水文モデルを適用し、そのモデルが持つ不確実性の高い４つのパラメータを推定する。

図１．6 時間先の流量予測結果と観測の比較

図 3． 日吉ダム流域と六甲山 X-MP レーダー
図２． 予測時間ごとの二乗誤差(RMSE)

38 − Academic Center for Computing and Media Studies, Kyoto University

4.2 シミュレーションの結果
図 4 はレーダパラメータ、B、C とモデルパラ

メータ K DP をパーティクルフィルタにより推定したものである、
その時間的な変化を示している。
各パラメータが観測値とのデータ同化によって不確実性の範囲が減少していることが確認できる。
図 5 と 6 はレーダパラメータのフィルタリング適

図 5．レーダパラメータをフィルタリングする場合用の有無による流量計算結果の違いを示している。
レーダとモデルのパラメータの同時フィルタリン
グで流量計算の精度が大きく向上することが確認できる。

図 6．レーダパラメータをフィルタリングしない場合

5

まとめ
本研究では、パーティクルフィルタを、流出モデルに導入して水文予測の高精度化を図るとともに、流出モデルの不確実性を評価した。スーパーコンピューティング環境を駆使した新たな降雨流
出予測システムの開発の適用性を確認し、分布型
流出モデルを基本とする高精度の実時間流出予測が期待される。

6

参考文献

S. J. Noh, Y. Tachikawa, M. Shiiba, and S. Kim,
“Applying sequential Monte Carlo methods
into a distributed hydrologic model: lagged
particlefilteringapproachwithregularization,” Hydrol. Earth Syst. Sci., vol.
15, pp. 3237-3251, 2011.
S. J. Noh, Y. Tachikawa, M. Shiiba, and S. Kim,
“Ensemble Kalman filtering and particle
filtering in a lag-time window for short-term
streamflow forecasting with a distributedhydrologic model,” J. Hydrol. Eng., ASCE, inpress.

図 4．パラメータ分布の変化

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 39

研究紹介

SIMD 演算を考慮した fill-in 戦略を用いた IC 分解前処理
武村 直和 ∗ 岩下 武史 † ，中島 浩 †∗†1

京都大学大学院情報学研究科

京都大学学術情報メデ ィアセンター

はじめに
近年，多くのプロセッサが演算性能向上のために一命令で複数のデータを処理する SIMD(Single In-

struction Multiple Data) 演算器を備えている．例えば ，x86 系のプロセッサで利用できる SIMD 拡張命
令セットとして SSE(Streaming SIMD Extensions)や Intel AVX(Advanced Vector eXtensions) が知られている．また，これらの SIMD 命令において，一命令で処理できるデータ量（ SIMD 幅）も増加の傾向にある．こうした背景から，数値シミュレーショ
ンにおいて，プロセッサの演算性能を十分に利用するためには SIMD 命令の活用が重要となってきている [1][2]．そこで，本研究では，有限要素解析等により生ずる連立一次方程式の求解法としてよく用いられる ICCG 法（ より正確には，IC(0) 分解前処理付
き共役勾配法）に対して，SIMD 演算を利用することによる性能改善を試みる．

本稿では，上記の検討の結果として，前処理行列をブロック化し ，ブロック内での fill-in を考慮する

IC 分解前処理（ IC(B) 前処理と表記する）を提案す
る．本手法では，IC(0) 前処理と比べ，分解中に許容される fill-in が増加することによる前処理行列の近似度の向上，即ち ICCG 法における収束性の向上が期待できる．一方，提案手法の導入により，一般に一反復あたりの計算量は増加する．これは，非零ブ
ロック内の要素はすべて非零要素として扱われことによる．従って，提案する前処理が IC(0) 前処理と比べて有利となるかど うかは，収束性向上（ 反復数の低減）と一反復あたりの計算時間の増加のトレード オフにより決まる．しかし ，提案手法を用いた場
合，反復中の主要な計算核である代入計算と行列ベクトル積はブロック毎の密行列計算により構成され，本計算に対しては SIMD 演算を活用することが期待できる．つまり，SIMD 演算により，ブロック化に伴う計算量増加の影響を相殺できれば ，収束性向上

ICCG 法は線形反復法の一種で，反復計算における主要な計算核は前処理部分の前進後退代入計算とに伴って全体の計算時間を削減することが期待でき

共役勾配（ CG ）法に含まれる疎行列ベクトル積であ利用した数値実験より，提案する前処理の有効性にる．そこで，本稿ではこれらの二つの計算核，特に従来からベクトル化・並列化が困難であるとされてる．本稿では，フロリダ大学の行列データベースをついて，IC(0) 前処理との比較から検証する．

きた前進後退代入計算を中心に，SIMD 演算を利用

2

疎行列計算の SIMD 化のスカラ型プロセッサを用いる．疎行列演算におい2.1

SIMD 演算て，このようなスカラ型のプロセッサにおける SIMD

SIMD は一命令で複数のデータを処理する技術であり，x86 CPU 向けの SIMD 拡張命令セットとしてした性能改善法について述べる．
本研究では，計算プラットフォームとして x86 系

演算機構を活用する方法として，SIMD 幅に応じたサイズを持つ密行列演算を導入する手法がある．そこで，著者らは ICCG 法において小規模な密行列演算を導入することによる性能改善手法について検討を行った．

40 − Academic Center for Computing and Media Studies, Kyoto University

SSE や Intel AVX がある．これらの命令セットを用いて効率的に SIMD 演算を行うには以下の条件がある．プロセッサ上での演算は 1. メモリからレジスタ
へのデータのロード → 2. レジスタ上のデータによる演算 → 3. レジスタからメモリへのデータのスト⎛

ア の順で行われる．ここで，1. のロード と 3. のス

a11
⎜
⎜ 0
⎜
⎜ 0
⎝
0

トアを効率的に行うには処理するデータが メモリ上で連続に存在することが望ましい．これは，非連続のデータを扱う際には gather や scatter といった低効率な処理が必要となるためである．また，2. の演
算を SIMD 化するには一度に行う演算に依存性がな

0
a22

a13
0

0
0

a33
0

⎞
0
⎟
a24 ⎟
⎟
a34 ⎟
⎠
a44

図 1: サンプル行列

いことが必要となる．

val={a11 ,a13 ,a22 ,a24 ,a33 ,a34 ,a44 }col ind={0,2,1,3,2,3,3}

ICCG 法は非構造メッシュによる有限要素解析等で用いられることが多いため，通常の実装手法では、

row ptr={0,2,4,6,7}

代入計算や疎行列ベクトル積において間接参照が必要となる．上記のとおり，間接参照はメモリに対し

図 2: CRS 形式て非連続なデータアクセスとなることが多く，SIMD
演算の効率性の点では好ましくない．そこで，本研
究では次小節で述べる疎行列のブロック化を行う．

れる．この場合，乗算ベクトル pは col ind 配列を用いて間接的に参照されるため，非連続なメモリアク2.2

セスを要し ，SIMD 演算の有効利用を阻害する．

疎行列のブロック化

次に，係数行列をブロック化し，Block CRS 形式

行列全体をブロックサイズ nl × nw のブロックに分け，非零要素を含む非零ブロックのみをブロック毎にメモリに格納する．この時，非零ブロックに含まれている零要素も格納することで，各非零ブロッ
クは密行列として扱われる．疎行列ベクトル積や代入計算をブロック単位で行うことで，ブロック内でSIMD 演算を利用可能となる [3]．
CRS 方式は，疎行列のデータ格納形式として最もよく使われる方式で，係数行列の非零要素の情報のみを格納する．それに対し，本研究では非零ブロックをを用いて格納した場合，密行列として扱われる非零ブロックに対して SIMD 演算が利用可能となる．本
研究で用いる x86 系プロセッサの場合，SIMD 演算を用いた実装法として，デ ィレクティブによる場合と C++組み込み関数による場合の二種類がある．ここで，後者の場合について，BCRS(2×2) 形式で行列を格納した場合の SIMD 化された疎行列・ベクト
ル積を図 5 に示す．

単位としたデータ格納形式である Block CRS 形式を

3IC(B) 分解前処理用いる．本形式では，非零ブロックの行・列に関する

3.1

前処理付き共役勾配法

位置情報については CRS 形式を用い，各要素の値については非零ブロックを単位として順に一次元配列

共役勾配 (CG) 法は連立一次方程式に格納する．図 1 の様なサンプル行列を CRS 形式でAx=b

格納する場合，図 2 に示すような 3 つの一次元配列を(1)用いる．また，同じ行列を BCRS(2×2) 形式で格納すを解くための Krylov 部分空間法の一つであり，係数る場合，図 3 のようになる．ここで，BCRS(nl × nw )

行列A が正定値かつ対称である場合に用いることが形式はブロックサイズ nl × nw の Block CRS 形式をできる [4, 5]．

表す．

2.3

ブロック化された疎行列による行列・
ベクト ル積

係数行列を CRS 形式で格納した場合，疎行列・ベクトル積Ap=q の計算は，一般に図 4 の様に記述さ

bval={a11 ,0,0,a22 ,a13 ,0,0,a24 ,a33 ,0,a34 ,a44 }bcol ind={0,1,1}
brow ptr={0,2,3}
図 3: BCRS(2×2) 形式

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 41次のように表せる．
⎧
i−1
�
⎪
⎪
⎨ ai,i −u2k,i dk,kdi,j =
k=1⎪⎪⎩0// Compute Ap=qfor(i=0; i<N; i++){q[i]=0.0
for(j=row_ptr[i]; j<row_ptr[i+1]; j++){q[i]+=val[j] × p[col_ind[j]]}}(i = j)(3)(i �= j)

⎧
i−1
�
⎪
⎪
⎪
(
a
−
uk,i uk,j dk,k )/di,i
⎪
i,j
⎪
⎪
⎪
k=1
⎪
⎨
(i < j ∩ fi,j = 1) (4)ui,j =⎪⎪⎪⎪1(i = j)⎪⎪⎪⎪⎩0(otherwise)

図 4: CRS 形式を用いた疎行列ベクトル積// Compute Ap=q__m128d XA0, XA1, XP0, XP1, XQ0;
Iofs=0;Vofs=0;for(k=0; k<Nd; k++){XQ0 = _mm_setzero_pd();for(j=brow_ptr[k] ;j<brow_ptr[k+1];j++){Jofs = bcol_ind[j]*2;
XP0 = _mm_load1_pd(&p[Jofs]);
XP1 = _mm_load1_pd(&p[Jofs+1]);
XA0 = _mm_load_pd(&bval[Vofs]);
XA1 = _mm_load_pd(&bval[Vofs+2]);
XA0 = _mm_mul_pd(XA0,XP0);
XQ0 = _mm_add_pd(XQ0,XA0);
XA1 = _mm_mul_pd(XA1,XP1);
XQ0 = _mm_add_pd(XQ0,XA1);
Vofs+=4;}_mm_store_pd(&q[Iofs],XQ0);
Iofs+=2;}

ここで，行列を IC 分解する際に，ai,j = 0であっても，ui,j �= 0 となることを fill-in と呼ぶ．従って，行列F = {fi,j } は fill-in を制御するために用いられ，これにより前処理行列の非零パターンを変える事ができる．

F を次のように定めた時，fill-in は全て棄却され Aの上三角部分と U の非零パターンは等しくなり，これを IC(0) 分解と呼ぶ．IC(0) 分解によって得られる行列による前処理を IC(0) 分解前処理と呼び，単に ICCG 法という場合には IC(0)-CG 法を指すと考えるのが一般的である．

⎧⎨1 (a �= 0)i,jfi,j =⎩0 (a = 0)

図 5: BCRS(2×2) 形式を用いた疎行列ベクトル積

CG 法の収束性を改善する方法として前処理があTる．これは A を近似した前処理行列M = ML ML

を用いて前処理を行い，解く方程式を(ML−1 AML−T )MLT x=ML−1 b(2)とすることで収束性を向上させる方法である．

(5)

i,j

3.3

IC(B) 前処理

IC(0) 分解により生成される前処理行列は係数行列と同じ非零パターンを有しており，前処理過程である代入計算において効果的に SIMD 演算を利用することは難しい．そこで本稿では，ブロック化した

3.2

疎行列に対する IC 分解の手法として IC(B,nl ×nw )IC(0) 分解前処理

分解を提案する．

CG 法に対してよく用いられる前処理として，係数行列A= {ai,j } を不完全コレスキー分解 (IC 分解)

することで得られるLDLT (≈ A) を前処理行列と提案手法では，まず係数行列をブロックサイズ nl ×nw

でブロック化し，非零ブロックを BCRS 形式で格納する．次に，IC 分解において，各ブロック内での fill-

する IC 分解前処理がある．IC 分解による求まる対in を許容して分解を行う．即ち，本手法は，非零ブ角行列D= {di,j }と上三角行列U = {ui,j } (=L ) は

ロックに含まれる要素の行・列番号の集合を B と表T

42 − Academic Center for Computing and Media Studies, Kyoto University

した場合，F として次式を用いて分解することと定義できる．

fi,j =⎧⎨1⎩0

表 2: 疎行列ベクトル積の比較((i, j) ∈ B)

(6)(otherwise)

計算時間 [s]本手法では，非零ブロック内の fill-in が許される

parabolicfem

thermal2

G3circuitため，fill-in を一切考慮しない IC(0) 分解前処理と比

CRS

0.79

2.75

2.14

べ収束性の向上が期待できる．また，この場合，前進・後退代入計算は，ブロック内で密行列計算とし

BCRS(1×2)
BCRS(2×1)

0.90
0.73

3.03
2.31

2.15
1.59

て扱われる．その結果，非対角ブロックに対する代BCRS(2×2)

1.01

2.88

2.23

入計算は図 5 に示した疎行列・ベクトル積と同様にSIMD 化することができる．なお，対角ブロックにBCRS(1×4)
BCRS(4×1)

1.17
1.03

3.34
2.69

2.62
1.98

ついては処理に逐次性があるため，SIMD 化は行わ

BCRS(4×4)

2.08

5.49

3.75

ない．

4

数値実験

4.1

実験条件

いる．従って，提案手法である IC(B) 前処理を用いた場合においても係数行列の格納形式として CRS 形式を利用することや，IC(B) 前処理と異なるブロッ
クサイズを用いた BCRS 形式を使用することは可

本稿では，フロリダ大学の Matrix Collection[6] の行列データを用いた数値実験結果について報告する．

Matlix Collection データベースから，問題サイズが大きく，右辺ベクトルが提供されている対称行列問
題を選び，実験の対象とした．数値実験に用いた計算機環境を表 1 に示す．今回の実験では，1 コアによる実行とし，SIMD の効果について検証を行う．ICCG
法の収束判定基準は，相対残差が 10−7 以下となった時点とした．ここで相対残差とは初期残差に対する残差の 2 ノルムの比 ||r||2 /||r0 ||2 である．IC(B) 前処理と BCRS 形式におけるブロックサイズには，(1×2),(2×1), (2×2), (1×4), (4×1), (4×4) を用いる．

用いた疎行列ベクトル積 100 回分の計算時間を示す．
いずれの数値実験においても，BCRS(2×1) 形式を利用した場合が最良の結果となった．CRS 形式との比
較では，扱われる行列要素の数は増加しているもののブロック内の計算に SIMD 計算が導入されたことによる効果が大きく，結果として CRS 形式を上回る性能が得られた．今回扱った行列は対称行列である
ため，ブロックサイズ (1×2) と (2×1) および (1×4)と (4×1) はそれぞれ計算量が等しいが，ブロックサ
イズ (2×1) と (4×1) の方が計算時間が短い．従って，
クサイズ (1×2) や (1×4) の利用は避けるのが良いと考えられる．

(1 コアのみ使用)
コンパイラ

BCRS 形式の有効性について検証を行った．表 2 に
parabolic fem,thermal2,G3 circuit の行列データをブロック内でリダ クション演算が必要となるブロッ

表 1: 実験環境
プロセッサ
Intel Xeon E5-2670

次に，表 3 に IC(0) 分解前処理と提案手法である

64GB

メモリ

能である．そこで，まず疎行列ベクトル積におけるIC(B) 分解前処理の比較を示す．ここでは疎行列ベク

Intel Compiler 13.1.0

トル積のための係数行列の格納形式としては，最も結果が良かったブロックサイズ (2×1) の BCRS 形式を用いる．3 つの数値実験の全てにおいて，適切なブ

4.2

ロックサイズを設定することにより提案手法は IC(0)-

数値実験結果IC 分解前処理のための代入計算には分解により得られた前処理行列が用いられ，これは CG 法中の疎行列ベクトル積で用いられる係数行列とは異なって

CG 法を上回る性能を実現している．いずれのデータ
セットの場合でも，ブロックのサイズが大きい方がfill-in 増加の効果が大きく，収束性が向上しているが，一反復あたりの計算時間も増加している．G3 circuit

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 43

では，ブロックサイズを (4×4) とした場合の収束性の向上幅が大きく，最も大きいブロックを使用した場合が最良となった．一方，parabolic fem と thermal2
では，それぞれブロックサイズ (2×1)，(2×2) の時に最も計算時間が最短となる結果となった．
表 3: IC(B) 前処理と IC(0) 前処理の比較

5

おわりに

ICCG 法において，係数行列をブロック化し ，ブロック内の fill-in を許容する IC(B) 分解前処理を提案した．本手法では，各ブロックを密行列として扱うことができ，前進・後退代入計算をブロック単位
で SIMD 化することができる．また，IC(0) 前処理と比べて許容される fill-in が増加するため，前処理効果の改善（ 反復数の低減）が期待できる．フロリ

parabolic fem
反復回数

計算時間 [s]時間/反復 [s]ダ大学の行列データベースを用いた数値実験によっ

IC(0)

678

12.8

1.88E-2て，適切なブロックサイズを設定することで提案手IC(B,1×2)
IC(B,2×1)

629
620

12.3
12.1

1.95E-2
1.96E-2

時間の短縮が可能であることを示した．

IC(B,2×2)

588

12.2

2.08E-2

IC(B,1×4)IC(B,4×1)

590
560

13.1
13.0

2.22E-2
2.31E-2

IC(B,4×4)

463

14.5

3.13E-2

thermal2法により，IC(0) 前処理と比べて求解に必要な計算

参考文献[1] 中島 研吾, 前処理付きマルチスレッド 並列疎行列ソル
バー情報処理学会研究報告, 2013-HPC-139, 6, pp.16 (2013).
[2] 西田 晃, 複数アーキテクチャ上での疎行列ベクトル
積の性能最適化手法, 情報処理学会計算機アーキテク
チャ研究会報告, 2009-ARC-186, 7, pp.1-8 (2009).

反復回数

計算時間 [s]時間/反復 [s]

IC(0)

1657

114

6.90E-2

IC(B,1×2)
IC(B,2×1)

1666
1629

104
101

6.25E-2
6.20E-2

IC(B,2×2)

1604

96

5.98E-2

IC(B,1×4)

1622

106

6.51E-2

IC(B,4×1)

1525

101

6.62E-2

IC(B,4×4)

1415

121

8.54E-2[4] Richard Barrett et al., Templates for the Solutionof Linear System : Building Blocks for Iterative
Methods, Society for Industrial and Applied Mathematics, 1987.
[5] Yousef Saad, Iterative Methods for Sparse Linear
Systems(2nd), Society for Industrial and Applied
Mathematics, 2003.

G3 circuit
反復回数

計算時間 [s]時間/反復 [s]

IC(0)

184

9.82

5.34E-2

IC(B,1×2)

156

8.06

5.17E-2

IC(B,2×1)
IC(B,2×2)

154
93

8.07
4.67

5.24E-2
5.03E-2

IC(B,1×4)

122

6.60

5.41E-2

IC(B,4×1)
IC(B,4×4)

121
51

6.50
3.36

5.37E-2
6.59E-2

44 − Academic Center for Computing and Media Studies, Kyoto University[3] Samuel Williams, Leonid Oliker, Richard Vuduc,John Shalf, Katherine Yelick and James Demmel, Optimization of sparse matrix-vector multiplication on emerging multicore platforms, Parallel Computing, 35, (2009), pp. 178-194.

[6] The University of Florida Sparse Matrix Collection, http://www.cise.ufl.edu/research/sparse/matrices/.

サービスの記録・報告システム A 運転状況 ( 2013 年 4 月 ～ 2013 年 9 月 )１）保守作業に伴うサービス休止およびシステムダウン障害発生状況
保守作業に伴うサービス休止
保守開始日時システムダウン障害発生状況

サービス再開日時

保守時間[h]障害発生日時

サービス再開日時

ダウン時間[h]2013/04/01

0:00

2013/04/02

9:00

33.00

2013/07/03

23:25

2013/07/04

5:40

6.25

2013/06/11

9:00

2013/06/11

20:00

11.00

2013/07/20

18:14

2013/07/20

18:32

0.30

2013/08/06

9:00

2013/08/07

9:00

24.00

2013/07/23

10:30

2013/07/23

19:00

8.50

２）サービス状況
バッチ

サービ
ス時間
[h]

処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼動
ノード数

4月

687.00105,507

102,169

12,428,500

9,982,010

939.9

60

%

5月

744.0040,155

121,387

15,411,300

11,949,000

939.9

68

%

6月

709.0017,78998,261

12,937,900

10,320,600

939.8

62

%

7月

728.9521,005

120,975

14,529,200

11,096,500

939.7

63

%

8月

720.0019,082

124,330

13,918,300

10,852,000

939.8

64

%

9月

720.0019,038

121,665

13,183,600

10,468,200

939.9

62

%

計

4308.95

222,576

688,787

82,408,800

64,668,310

939.8

63

%

ノード
利用率

ノード利用率(棒グラフ)＋コア利用率(折れ線グラフ)
100%
90%
80%
70%
60%
50%

ノード利用率

40%

コア利用率

30%
20%
10%
0%



4月

5月

6月

7月

8月

9月

占有時間 = 合計(経過時間×占有コア数)

平均稼動ノード数 = 電源 ON 状態のノード数の月平均 (10 分間隔のサンプリングデータより算出)



ノード利用率 = 稼動ノードに対するジョブが実行されているノードの割合

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 45

システム B 運転状況 ( 2013 年 4 月 ～ 2013 年 9 月 )１）保守作業に伴うサービス休止およびシステムダウン障害発生状況
保守作業に伴うサービス休止
保守開始日時システムダウン障害発生状況

サービス再開日時

保守時間[h]障害発生日時

サービス再開日時

ダウン時
間[h]

2013/04/01

0:00

2013/04/02

9:00

33.00

2013/07/17

9:40

2013/07/17

10:10

0.50

2013/06/11

9:00

2013/06/11

20:00

11.00

2013/07/23

10:30

2013/07/23

19:00

8.50

2013/08/06

9:00

2013/08/07

9:00

24.00

サービ
ス時間
[h]

処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼動
ノード数

4月

687.0066,103263,934

3,648,490

300,156

556.9

63

%

5月

744.0083,587

174,377

4,058,320

326,492

557.0

67

%

6月

709.0060,252434,553

4,154,320

333,411

551.1

68

%

7月

735.00121,759398,277

3,998,060

316,426

556.9

65

%

8月

720.00129,387

279,907

4,005,990

327,144

543.7

63

%

9月

720.00387,864

251,454

3,671,620

305,383

551.1

64

%

計

4315.00

848,952

1,802,502

23,536,800

1,909,012

552.8

65

%

２）サービス状況
バッチ
ノード
利用率

ノード利用率(棒グラフ)＋コア利用率(折れ線グラフ)
100%
90%
80%
70%
60%
50%

ノード利用率

40%

コア利用率

30%
20%
10%
0%

4月

5月

6月

7月

8月

9月



占有時間 = 合計(経過時間×占有コア数)

平均稼動ノード数 = 電源 ON 状態のノード数の月平均 (10 分間隔のサンプリングデータより算出)



ノード利用率 = 稼動ノードに対するジョブが実行されているノードの割合

46 − Academic Center for Computing and Media Studies, Kyoto University

システム C 運転状況 ( 2013 年 4 月 ～ 2013 年 9 月 )１）保守作業に伴うサービス休止およびシステムダウン障害発生状況
保守作業に伴うサービス休止
保守開始日時システムダウン障害発生状況

サービス再開日時

保守時間[h]障害発生日時

サービス再開日時

ダウン時
間[h]

2013/04/01

0:00

2013/04/02

9:00

33.00

2013/07/17

9:40

2013/07/17

10:10

0.50

2013/06/11

9:00

2013/06/11

20:00

11.00

2013/07/23

10:30

2013/07/23

19:00

8.50

2013/08/06

9:00

2013/08/07

9:00

24.00

サービ
ス時間
[h]

処理
件数

経過
時間[h]

4月

687.00913

3,112

66,072

46,270

16.0

25

%

5月

744.00

1,761

7,326

111,538

64,946

16.0

35

%

6月

709.00

1,295

9,854

108,761

76,927

15.8

45

%

7月

735.00

1,172

7,183

166,717

122,342

16.0

52

%

8月

720.00928

7,425

209,175

124,625

15.7

61

%

9月

720.00536

8,769

144,123

99,357

16.0

54

%

計

4315.00

6,605

43,669806,386

534,467

15.9

45

%

２）サービス状況
バッチ
占有
時間[h]

CPU
時間[h]

平均稼動
ノード数

ノード
利用率

ノード利用率(棒グラフ)＋コア利用率(折れ線グラフ)
100%
90%
80%
70%
60%
50%

ノード利用率

40%

コア利用率

30%
20%
10%
0%

4月

5月

6月

7月

8月

9月



占有時間 = 合計(経過時間×占有コア数)

平均稼動ノード数 = 電源 ON 状態のノード数の月平均 (10 分間隔のサンプリングデータより算出)



ノード利用率 = 稼動ノードに対するジョブが実行されているノードの割合

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 47




  


 











































!  "  #!  "  #!  "  #!  "  #$%&

$ '(!  "  #)*+&          !"# $%&'  %()*+,-./

48 − Academic Center for Computing and Media Studies, Kyoto University!  "  #サービスの記録・報告

センター利用による研究成果―平成 24 年度計算機利用結果報告書から―【1000 情報学】

Nakashima, Osni Marque：Parallel Smoother1.

Takeshi Iwashita, Yu Hirotani, Takeshi

Based on Block Red-Black Ordering forMifune, Toshio Murayama & Hideki Ohtani：

Multigrid Poisson Solver：Proc. Intl. Mtng.

Large-scale time-harmonic electromagnetic

High-Performancefield analysis using a multigrid solver on a

Computational Science, 2012

distributed memory parallel computer ：

2.

4.

2012

Simulations on Stray Photoelectron Current

河合直聡, 岩下武史, 中島 浩：ブロック化赤‐Flows Around a Spacecraft：J. Geophysics

黒順序付け法に基づく並列マルチグリッドポア

Rearch, Vol.117, No.09210, pp.1 – 13, 2012
Tasuku Hiraishi, Masaru Ueno, Tatsuya Abe,
Takeshi Iwashita, Hiroshi Nakashima ：10, 2012

Xcrypt in Lisp：A Scripting System for Job

野瀬田裕樹, 河合直聡, 岩下武史, 高橋康人, 美Level Parallel Programming in Lisp: Proc.

舩 健, 中島 浩：並列境界要素解析フレーム

Intl. Lisp Conf., 2012

ワークの設計と実装: 先進的計算基盤システム

10. Yasuhito Takahashi, Tadashi Tokumasu,

シンポジウム論文集（SACSIS2012）, pp.69 –Masafumi Fujita, Takeshi Iwashita, Hiroshi

70, 2012

Nakashima, Tadashi Tokumasu, Masafumi

Yasuhito

Takahashi,

Takeshi

Iwashita,

Fujita,Shinj

Wakao,

Koji

Fujiwara,

Hiroshi Nakashima, Tadashi Tokumasu,

Yoshiyuki Ishihara：Parallel Time-Periodic

Masafumi Fujita, Shinj Wakao, Koji Fujiwara,

Finite-Element Method for Fast Magnetic

Yoshiyuki Ishihara：Parallel Time-Periodic

Field Analysis of Induction Motors：Proc. 15th

Finite-Element Method for Steady-State

Biennial

Analysis of Rotating Machines：IEEE Trans.

Compuation, 2012Vol.48,

No.2,

pp.1019–1022,

Conf.

Electromagneti

Field

11. Kazuki Semba, Koji Tani, Takashi Yamada,

2012

TakeshiTakeshi Iwashita, Yasuto Takahashi, Hiroshi

Hiroshi Nakashima：Parallel Performance ofNakashima ： Algebraic Block Multi-Color

Multi-threaded

Ordering Method for Parallel Multi-Threaded

Algebraic Block Multi-color Ordering In

Sparse Triangular Solver in ICCG Method：

FiniteProc.

Analyses ： Proc.

Intl.

ParallelandDistributed

Iwashita,

Yasuhito

ICCG Solver

Element

Takahashi,
Based onElectromagnetic

Field

15th

Conf.

Biennial

Electromagneti Field Compuation, 2012

Tasuku Hiraishi, Tatsuya Abe, Takeshi

12. 南 武志，岩下武史，中島 浩：冗長な計算を

Iwashita, Hiroshi Nakashima：Xcrypt: a Perl

伴わない 3 次元 FDTD 法の時空間タイリング：

Extension

情報処理学会論文誌 コンピューティングシス

Programming

for

Job
：

Proc.

Level
Intl.

Parallel
WS.

High-performance Infrastructure for Scalable
7.

9.

ティングシステム（ACS）, Vol 5, No.3, pp.1 –Processing Symp., 20126.

Yohei Miyake, Hideyuki Usui, Hirotsugu
Kojima, Hiroshi Nakashima：Plasma Particle

Magnetics,
5.

for

Parallel Computing, Vol.38, pp.485 – 500,

ソンソルバ: 情報処理学会論文誌：コンピュー

3.

8.

Computing

テム(ACS), Vol.6, No.1, pp.56 – 65, 2013
13. 仙波和樹, 山田 隆, 岩下武史, 中島 浩, 高橋

Tools, 2012

康人, 美舩健：実応用有限要素電磁場解析におMasatoshi Kawai, Takeshi Iwashita, Hiroshi

ける折り畳み前処理を用いた反復解法の有効京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 49

性：電気学会論文誌, Vol.133, No.2D, pp.203 –22. 塘三千代, 里村雄彦：火星をイメージした地形
風の数値シミュレーション：日本気象学会 2012

213, 2013
14. Takeshi Iwashita, Akihiro Ida, Masatoshi
Kawai, Hiroshi Nakashima ： Performance

年度秋季大会予稿集, Vol.102, pp.P162 – P162,
2012

Evaluation of Multi-threaded Iterative Solver

23. Ogohara, K. & Satomura, T.：Climatology ofon Recent Processors：Proc. 2013 SIAM Conf.

Dust Haze Expansion on Mars：AOGS-AGU

Computational Science and Engineering,(WPGM) Joint Assembly 2012
24. Nattapon Mahavik, Takehiko Satomura, &2013
15. Tasuku Hiraishi, Masaru Ueno, Tatsuya Abe,

Somchai Baimuang：Z-R parameter variation

Motoharu Hibino, Takeshi Iwashita, Hiroshi

using conventional weatehr radar in the

Nakashima ： Xcrypt on Lisp: A Scripting

middle of Indochina peninsular：日本気象学会

System for Job Level Parallel Programming2012 年度秋季大会予稿集, Vol.102, pp.P383

in Lisp：International Lisp Conference, pp.

-P383, 2012
25. Nattapon Mahavik & Takehiko Satomura：

107 – 114, 2012
16. Tasuku Hiraishi, Tatsuya Abe, Takeshi

Calibration of weather radar for analyzing

Iwashita, Hiroshi Nakashima：Xcrypt: A Perlprecipitation distribution in Thailand：日本気

Extension

象学会 2012 年度秋季大会予稿集, Vol.102, pp.

for

Job

Programming：

Workshop

on

Level

SecondParallelInternational

High-performance

P143 - P143, 2012
26. Yuri

Komori&T.

Satomura ： The

Infrastructure for Scalable Tools (WHIST

precipitation peak off the western coast of2012), pp. 1 – 9, 2012

Indochina

17. 上野 優，平石 拓，日比野元春，岩下武史，
中島 浩：ジョブ並列スクリプト言語 Xcrypt

Peninsula ： Proc.

2nd

Int.

Workshop on Nonhydrostatic Numerical
Models, pp. 95 – 96, 2012における RPC に基づく多言語対応：情報処理【2300 ゲノム科学】

学会第 92 回プログラミング研究会, pp.1 – 14,

27. 大和勝幸，河内孝之：見えてきたゼニゴケゲノ

2013

ム：植物科学最前線, Vol.3, pp.71 – 83, 2012

18. 松井 健，平石 拓，八杉昌宏，馬谷誠二：高28. Shigeyuki Tsukamoto, Satoshi Hirao, Lixy

速版 Barnes-Hut 多体シミュレーションの並列

Yamada, Hitoshi Sawada, Kimitsune Ishizaki,

実装：先進的計算基盤システムシンポジウム

Takayuki Kohchi, Katsuyuki T. Yamato ：（SACSIS2012）, pp.1 – 8, 2012

Transcriptomic and proteomic analyses of

19. 松井 健，平石 拓，八杉昌宏，馬谷誠二：ワ

sperm of Marchantia polymorpha L.：第 54 回

ークスティーリングフレームワークにおける集日本植物生理学会年会 要旨集, pp.356 – 356,

団通信機能：並列/分散/協調処理に関するサマ2013

ー・ワークショップ（SWoPP2012）, pp.1 – 8,【3900 心理学】

2012

29. 川島康裕, 山城博幸, 山本洋紀, 村瀬智一, 梅田

20. 日比野元春，南 武志，平石 拓，岩下武史，雅宏, 樋口敏宏：照明強度のヒト大脳視覚野活中島 浩：Xcrypt を用いた 3 次元 FDTD 法プ動への影響：fMRI 研究：平成 24 年度照明学会

ログラムの自動チューニング：日本応用数理学

第 45 回全国大会講演要旨集, pp.8-17 - 8-17,

会年会（JSIAM2012）, pp.277 – 278, 20122012【2000 環境学】

【4100 数学】

21. 小森友里, 里村雄彦：領域気象モデル（WRF）

30. 新谷誠：5 元体上のいくつかの Optimalを用いたインドシナ半島西岸沖の降水極大の再Linear Code について：2012 年度応用数学合

現実験：日本気象学会 2012 年度秋季大会予稿

同研究集会, pp.62 – 63, 2012

集,Vol.102, pp.P387‐P387, 2012

50 − Academic Center for Computing and Media Studies, Kyoto University

31. Kiori Obuse, Shin-ichi Takehiro & Michio

Yamada：Weak interaction between zonal jets

39. 鵜沼 昂, 竹見哲也：2012 年 5 月 6 日につくばon a β plane：Japan Journal of Industrial

市付近で発生した竜巻の領域気象モデルによるand Applied Mathematics, Vol.30, No.1,

再現実験：平成 24 年度京都大学防災研究所研pp.111 – 127, 2013

究発表講演会講演要旨 C03, pp.1 – 1, 2013
40. Y.

【4200 天文学】

Miyake,

H.

Usui,

H.

Kojima,H.

32. Kenji Hasegawa, Benoit Semelin ： The

Nakashima：Plasma particle simulations on

impacts of ultraviolet radiation feedback onstray photoelectron current flows around agalaxies during the epoch of reionization：

spacecraft：Journal of Geophysical Research,

Monthly Notices of the Royal Astronomical

Vol.117, pp. A09210:1 – 13, 2012

Society, Vol.428, No.1, pp.154 – 166, 2013

41. 石井翔大, 内山雄介, 宮澤泰正：黒潮続流フロン

33. Kenji Hasegawa, Benoit Semelin：Impacts of

ト周辺におけるサブメソスケール渦の発生機

ultraviolet radiation feedback on the cosmic

構：2013 年度日本海洋学会春季大会・講演要旨

reionization

集, 2013

history ： AIP

conference

proceedings, Vol.1480, pp.267 – 270, 2012

42. 甲斐田秀樹, 内山雄介, J.C. McWilliams：超高34. Masayuki Umemura, Hajime Susa, Kenji

解像度 3 次元海洋モデルによる砕波帯‐陸棚相

Hasegawa, Tamon Suwa, Benoit Semelin：

互作用の解析：2013 年度日本海洋学会春季大Formation and Radiative Feedback of First

会・講演要旨集, 2013

Objects and First Galaxies ： Progress ofTheoreticalandExperimental

Physics,

43. 栗山貴生, 内山雄介, 宮澤泰正：流路変動に伴う黒潮系水塊の瀬戸内海への波及について：2013年度日本海洋学会春季大会・講演要旨集, 2013

pp.306 – 328, 2012

44. Uchiyama, Y., Ishii, T., Tsumune, D. &【4400 地球惑星科学】
35. Kei Nishina, Yoichi Ishikawa & Toshiyuki

Miyazawa,

Y ： Coastal

dispersal

ofAwaji：Impact of simultaneous assimilation ofradionuclides released from the Fukushima

intermediate velocities derived from Argo

Nuclear Power Plant：AGU Fall Meeting

float trajectories together with temperature

2012, San Francisco, CA, USA, 2012and salinity profiles in a high-resolution

45. Ishii, S., Uchiyama, Y. & Miyazawa, Y.：

4-dimensional variational data assimilation

Submesoscale Frontal Dynamics in thesystem：4th Argo Science Workshop, Venice,Kuroshio Extension Jet：AGU Fall Meeting

Italy, 20122012, San Francisco, CA, USA, 2012

36. Kei Nishina, Yoichi Ishikawa & Toshiyuki

46. Kuriyama, T., Uchiyama, Y. & Miyazawa,Awaji ：Impact of surface drifter data for

Y.：Impact of the Kurishio paths on oceanic“ocean weather forecast”- an example ofand estuarine circulations in and around Seto

mesoscalevariability

Inland Sea, Japan：AGU Fall Meeting 2012,

Extension

Region ： 4thInternatiional

Workshop

on

the

in

Modelliing

theKuroshio
Ocean,San Francisco, CA, USA, 2012
47. 宮地哲朗, 榎本剛, 太田洋一郎, 三好建正,

Hann-Ming Henry Juang：NCEP-GFS を用い

Yokohama, Japan, 201237. 仁科 慧, 淡路敏之, 石川洋一：中規模変動の再た Atmospheric River の予報実験：研究集会「週

現性向上に対する Argo による水温・塩分観測間及び 1 か月予報における顕著現象の予測可能値と中層流速推定値の同時同化の有効性：日本性」, pp.76 – 80, 2012

海洋学会秋季大会, 静岡, 2012

48. Kei

Nishina,Yoichi

Ishikawa,

Toshiyuki

38. 鵜沼 昂, 竹見哲也：高分解能計算による 2012

Awaji：Simultaneous assimilation of surface

年 5 月 6 日につくば市付近で発生した竜巻の再drifter data, satellite and in-situ observations

現実験：日本流体力学会 2012 年度年会拡張要

for

旨集, USB, 15, pp.1 – 3, 2012variability in the Kuroshio Extension Region：

improved

estimatesof

meso-scale

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 51

Journal of Oceanography, Vol.69, No.2,

SWCNT Growth: Armchair versus Zigzag

pp.175 – 190, 2013

Nanotubes：J. Am. Chem. Soc, Vol.134, pp.

49. 鹿倉洋介, 深畑幸俊, 平原和朗：プレート運動と9311 – 9319, 2012

地震発生に伴う西南日本内陸活断層の長期的な57. M. A. Addicoat, A. J. Page, Z. E. Brain, L.

破壊関数の変化：地球惑星科学連合 2012 年大Flack, K. Morokuma & S. Irle：Optimization

会, No.SSS38-11, 2012

of50. 鹿倉洋介, 深畑幸俊, 平原和朗：昭和東南海・南aGenetic

Algorithm

for

the

Functionalization of Fullerenes ：J. Chem.

海地震後のプレート間固着分布を考慮した西南

Theo. Comp, Vol.8, pp.1841 – 1851, 2012日本内陸活断層のクーロン破壊関数の長期変58. T. Kawatsu , J. Hasegawa ： Sequentially

化：日本地震学会 2012 年秋季大会, No.P1-67,

Coupled Hole-Electron Transfer Pathways for2012

Bridge-Mediated Triplet Excitation Energy
Transfer：J. Phys. Chem. C, Vol.116, No.44,【4500 プラズマ科学】
51. Masafumi

Shoji,

Yoshiharu

Omura

&

pp.23252 – 23256, 2012

Multidimensional

59. J. Hasegawa, K. J. Fujimoto & T. Kawatsu：A

Nonlinear Mirror-mode Structures in the

Configuration Interaction Picture for a

Earth ’ s Magnetosheath ： Journal of

Molecular Environment Using Localized

Geophysical Research, Vol.117, pp.A08208 -Molecular Orbitals：The Excited States ofA08208, 2012

Retinal Proteins: J. Chem. Theo. Comp., Vol.8,

Lou-Chuang

Lee：

52. Masafumi Shoji , Yoshiharu Omura ：Precipitation of Highly Energetic Protons by

No.11, pp.4452 – 4461, 2012
60. T. Kawatsu, K. Matsuda & J. Hasegawa：Ion

Singlet Excitation Energy Transfer MediatedCyclotron Triggered Emissions：Journal of

by Local Exciton Bridges：J. Phys. Chem. C,

HeliumBranchElectromagnetic

Geophysical Research, Vol.117, pp.A12210 -Vol.116, No.26, pp.13865 – 13876, 2012
61. T. Kawatsu , J. Hasegawa：Excitation Energy

A12210, 2012
【4600 基礎化学】

Transfer in GFP-X-CFP Model Peptides

53. A. J. Page, C.-P. Chou, B. Q. Pham, H. Witek,(X=amino acids): direct vs. through-bridge

S. Irle & K. Morokuma：Quantum Chemicalenergy transfers：Intern. J. Quantum. Chem,

Investigation of Epoxide and Ether Groups in

Vol.113, No.4, pp.563 – 568, 2013

Graphene Oxide and Their Vibrational

62. J. Hasegawa , K. Matsuda：Theoretical study

Spectra：Phys. Chem. Chem. Phys., Vol.15, pp.

of the excited states and the redox potentials

3725 – 3735, 2013

of

54. L. Ding, L. W. Chung, K. Morokuma ：Reaction

Mechanism

Decarboxylation
Green

of

the

Fluorescentofunusually

distortedβ-trifluoromethylporphycene ： Theor. Chem.

Photo-induced

Acc., Vol.130, No.2-3, pp.175 – 185, 2011

Photoactivatable

63. T. Kawatsu, K. Matsuda & J. Hasegawa：Protein:An

Bridge-Mediated Excitation Energy Transfer

ONIOM(QM:MM) Study：J. Phys. Chem. B,

Pathways through Protein Media: a Slater

Vol.117, pp.1075 – 1078, 2013

Determinant-Based

55. A. J. Page, T. Isomoto, J. Knaup, S. Irle & K.

Calculation

Electronic

Combinedwith

Coupling
Localized

Morokuma：Effects of Molecular Dynamics

Molecular Orbitals：J. Phys. Chem. A, Vol.115,

Thermostats on Descriptions of Chemical

No.39, pp.10814-10822, 2011

Non-Equilibrium：J. Chem. Theo. Comp, Vol.
8, pp.4019 – 4028, 2012
56. J. Kim, A. J. Page, S. Irle & K. Morokuma：
The Dynamics of Local Chirality during

52 − Academic Center for Computing and Media Studies, Kyoto University

64. Y. Kitagawa, K. Matsuda & J. Hasegawa：
Theoretical Study of the Excited States of the
Photosynthetic
Photosystem

Reaction
II:

Electronic

Centerin

Structure,

Interactions, and Their Origin ： Biophys.

71. Yoshino, J., Nakamura, Y., Kunitomo, S.,

Chem., Vol.159, No.2-3, pp.227 – 236, 2012

Hayashi, N. & Higuchi, H. ：Synthesis of

65. J. Hasegawa, T. Kawatsu, K. Toyoya & K.

tris(4-amino-2,6-dimethylphenyl)borane andMatsuda ： Chemical-intuition based LMO

facile extension of its pi-conjugated system bytransformation simplifies excited-state wave

utilizing the reactivity of the amino groups：

functions of peptides：Chem. Phys. Letters,

Tetrahedron Lett., in press , 2013
72. Kano, N., Furuta, A., Kambe, T., Yoshino, J.,

Vol.508, No.1-3, pp.171 – 176, 2011
66. Yasuhito Karino , Nobuyuki Matubayasi：

Shibata, Y., Kawashima, T., Mizorogi, N. &Interaction-component analysis of the urea

Nagase, S. ： 2,2'-Diborylazobenzenes witheffect on amino acid analogs：Phys.Chem.

doubleChem. Phys., Vol.15, No.12, pp.4377 – 4391,

fluorescent properties by substituents and2013

redox reactions：Eur. J. Inorg. Chem., pp.1584

67. Tomonori Kawakami, Isamu Shigemoto &N-B

coordination:

controlof– 1587, 2012

Nobuyuki Matubayasi：Free-energy analysis

73. 吉野惇郎, 中村優里, 國友静香, 林 直人, 樋口

of water affinity in polymer studied by

弘行：トリ(p-アミノフェニル)ボラン誘導体の合atomistic molecular simulation combined

成と性質：日本化学会近畿支部・高分子学会北with the theory of solutions in the energy

陸支部合同福井大会(平成 24 年度北陸地区講演

representation：J. Chem. Phys., Vol.137, No.

会と研究発表会), 2012
74. 吉野惇郎, 中村優里, 國友静香, 林 直人, 樋口

23, pp.234903-1 - 234903-9, 2012
68. Kazuhiro

Takemura,

Hao

Guo,

Shun

弘行：トリ(p-アミノフェニル)ボラン誘導体の合Sakuraba, Nobuyuki Matubayasi & Akio

成，構造，および性質：日本化学会第 93 春季

Kitao：Evaluation of protein-protein docking

年会, 2013

model structures using all-atom molecular【5000 機械工学】

dynamics simulations combined with the

75. Naoto Yanagawa, Tetsuaki Takeda, Shumpei

energyFunatani：Numerical Analysis of Control ofrepresentation ： J. Chem. Phys., Vol.137,

Natural Circulation by Injection of Helium

No.21, pp.215105-1 - 215105-10, 2012

Gas：6th International Topical Meeting on

solution

theoryinthe

69. Hideaki Takahashi, Atsushi Omi, Akihiro
Morita & Nobuyuki Matubayasi：Simple andexact approach to the electronic polarization
effectontheFormulation

solvationformolecularmechanicalfree

Vol.136,

Technology

HTR2012, pp.1 – 5, 201276. 柳川直澄，武田哲明，舩谷俊平，野村真司：Hequantum-mechanical/

析：日本原子力学会関東・甲越支部第１１回若system手研究者発表討論会: pp.8 – 9, 2012

and

its

No.21,

pp.214503-1-

77. 上田 桂, 西川雅章, 北條正樹, 武田一朗：繊維よれを含む繊維強化複合材料の圧縮強度に関する座屈分岐解析：第４回日本複合材料合同会議

70. Munetaka Takeuchi, Nobuyuki Matubayasi,
Yasuo Kameda, Babak Minofar, Shin-ichi&

Reactor

ガスを用いた自然循環流の制御に関する数値解

214503-12, 2012

Ishiguro,Temperatureenergy:applications to aqueous solutions：J. Chem.
Phys.,High

Yasuhiro

講演論文, pp.(1A-05)1 - (1A-05)2, 201378. 加藤幹, 西川雅章, 田渕雄一郎, 北條正樹：

Umebayashi ：PEFC 電極マイクロポーラス層用薄膜構造体のFree-Energy and Structural Analysis of Ion

細孔構造による押込み変形への影響の評価：電Solvation and Contact Ion-Pair Formation of

気化学会・第 53 回電池討論会講演論文集,

Li+ with BF4 － and PF6 － in Water andpp.363 – 363, 2012

Carbonate Solvents：J. Phys. Chem B, Vol.116,
No.22, pp.6476 – 6487, 2012

79. 津村祐介, 新玉重貴, 西川雅章, 北條正樹：
CFRP 積層板初期破壊に及ぼすボイドの影響の京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 53

評価：日本機械学会・第 20 回機械材料・材料

of Spouted Bed：Proceedings of 5th Asian

加工技術講演会講演論文集, pp.(614)1 - (614)6,

Particle Technology Symposium, No.078,

2012

2012

80. 津村祐介, 新玉重貴, 西川雅章, 北條正樹：画像

90. S. Miyai, T. Katsuo, T. Tsuji, T. Takayama &解析による繊維位置探索を用いた CFRP90°

T. Tanaka：DEM Modelling of the digging

層内部のひずみ定量化：第 37 回複合材料シン

process of gravel: Influence of particle

ポジウム講演論文集, pp.(A-19)1 - (A-19)2,

roundness ： Discrete element modeling of2012

particulate media, pp.141 – 148, 2012

81. 後藤 聡, 西川雅章, 北條正樹：マイクロボンド【5100 電気電子工学】

法による炭素繊維への熱硬化性樹脂の接着・硬91. Hiromitsu Awano, Hiroshi Tsutsui, Hiroyuki

化過程の観察：第 37 回複合材料シンポジウム

Ochi & Takashi Sato：Bayesian Estimation of

講演論文集, pp.(A-28)1 - (A-28)2, 2012

Multi-Trap RTN Parameters using Markov

82. 松田慎平，西川雅章，北條正樹：高分子複合材

Chain

Monte

Carlo

Method ： IEICE

料の緩和時間の温度依存性解析：日本機械学会

Transactions on Fundamentals of Electronics,

M&M2012 材料力学カンファレンス講演論文

Communications and Computer Sciences, Vol.

集, pp.(OS0902)1 - (OS0902)3, 2012

E95-A, No.12, pp.2272 – 2283, 2012

83. 上田 桂，西川雅章，北條正樹，武田一朗：繊

92. Hiromitsu Awano, Hiroshi Tsutsui, Hiroyuki

維未含浸部を有する複合材料の圧縮強度に関す

Ochi & Takashi Sato ： Multi-Trap RTN

る有限要素解析：第 54 回構造強度に関する講

Parameter Extraction Based on Bayesian

演会講演論文集, pp.288 – 290, 2012

Inference ： International Symposium on

84. Masaaki NISHIKAWA, Masaki HOJO &Satoshi GOTO：Finite Element Simulation

Quality Electrical Design (ISQED) (Santa
Clara,CA), 2013

for Curing and Adhesion Process of Epoxy【5200 土木工学】

Resin to Carbon Fiber：Proc. of International

93. 伊藤靖晃，白土博通，野澤剛二郎：等方性乱流

ConferenceonComposite

Interfaces(Interface 21), pp.(OR41)1 - (OR41)4, 201285. 藤井勧人：二次元振動噴流の発振制御に関する
数値解析：日本機械学会 2012 年度年次大会
DVD 論文集, 2012
86. 松田直樹，琵琶志朗，林高弘：材料非線形性の高感度測定のためのラム波モードの考察: 第 20
回超音波による非破壊評価シンポジウム講演論
文集, pp.29 – 30, 201387. 藤原 忍，辻 拓也，田中敏嗣：高濃度粒子が存在する流れ場に壁面が与える影響について（粒子が移動を伴う場合）：混相流, Vol.26,
pp.489 – 497, 2013
88. K. Higashida, Y. Okuyama, T. Tsuji & T.

場における構造基本断面の空気力の空間相関特性：日本流体力学会年会 2012, pp.229-1 - 229-2,
2012
94. 伊藤靖晃，白土博通，野澤剛二郎：LES による
正方形角柱の変動空気力の空間相関評価に関する検討: 数値流体力学シンポジウム,
pp.E10-2-1 - E10-2-6, 2012
95. 伊藤靖晃，白土博通，野澤剛二郎：LES による
等方性乱流場における変動空気力の空間相関評
価に関する検討：構造工学論文集, Vol.59A,
pp.616 – 626, 2013
96. 鳥生大祐，牛島 省：多相場モデルによる自由
水面流中の砂面変形の数値解法：計算工学講演
会論文集, Vol.17, No.E-1-1, 2012

Tanaka：A numerical model for the motion of

97. 鳥生大祐，牛島 省：多相場モデルを用いた多large object in gas-fluidized bed：Proceedings

孔質体内の自然対流の３次元計算：土木学会論

of 5th Asian Particle Technology Symposium,

文集 B1(水工学), Vol.69, No.4, pp.I_613 - I_618,

No.163, 20122013

89. Z.A. Ali, H. Yada, T. Tsuji & T. Tanaka：

98. Kishida, K., Sawada, A., Yasuhara, H. &Microscopic DEM-CFD Coupling Simulation

Hosoda, T. ： Estimation of Fracture flow

54 − Academic Center for Computing and Media Studies, Kyoto University

considering the inhomogeneous structure of

108.Nurmalia, N. Nakamura, H. Ogi, M. Hirao &single rock fractures：Soils and Foundations,K. Nakahata：Mode conversion behavior of

Vol.53, No.1, pp.105 – 116, 2013

SH guided wave in a tapered plate：NDT &99. 崔 瑛, 野々村政一, 井浦智美，岸田 潔，木村International, Vol.45, pp.156 – 161, 2012

亮：小土被り土砂地山でのトンネル掘削におけ【5300 建築学】

る事前地山改良工の効果に関する解析的検討：109.石川浩一郎、松下千裕：屋根型単層ラチスシェ

土木学会論文集 F1（トンネル工学）特集号,

ルの地震荷重伝達能力に限界づけられるトラス

Vol.68, No.3, pp.I_27 - I_37, 2012

アーチ架構の動的挙動と保有耐力評価：日本建100.藤岡 奨, 牛島 省：収束計算を伴う界面捕捉

築 学 会 構 造 系 論 文 集 , Vol.78, No.683,

手法による３次元 PLIC のアルゴリズム：土木

pp.101~110, 2013

学 会 論 文 集 A2( 応 用 力 学 ), Vol.68, No.2,【5400 材料工学】

pp.I_207~I_216, 2012

110.Y. Hinuma, F. Oba, Y. Kumagai, & I. Tanaka：

101.山下 遼, 藤岡 奨, 牛島 省：自由水面流中のIonization potentials of (112) and (11-2) facet

接触を伴う複雑形状物体運動の並列計算法：土surfaces of CuInSe2 and CuGaSe2：Physical

木学会論文集 A2( 応用力学), Vol.68, No.2,

Review

pp.I_171~ I_178, 2012

245433-7, 2012

102.藤岡 奨，牛島 省：界面の平面近似によるフB,

Vol.86,

No.24,

pp.245433-1-

111.M. Choi, F. Oba, Y. Kumagai, & I. Tanaka：

ラックス評価を用いた VOF 関数移流計算：17

Antiferrodistortive-like

回計算工学講演会講演論文集, pp.E-2-3, 2012rotation induced by the oxygen vacancy in

103.山下 遼, 藤岡 奨, 牛島 省：流体との力学的

cubic SrTiO3：Advanced Materials, Vol.25,連成を考慮した接触を伴う複雑形状物体運動の並列計算法：第 17 回計算工学講演会講演論文集, pp.B-4-3, 2012

oxygen-octahedron

No.1, pp.86 – 90, 2013
112.R. Ishikawa, N. Shibata, F. Oba, T. Taniguchi,
S. D. Findlay, I. Tanaka, & Y. Ikuhara：

104.中畑和之, 中 民矢, 飯田洋志, 石井 武, 松田

Functional complex point-defect structure in

吉人：EMFIT を用いた電磁波の伝搬解析とア

a huge-size-mismatch system ： Physical

レイアンテナを用いた鉄筋の映像化シミュレー

Review Letters, Vol.110, No.6, pp.065504-1-

ションへの応用：コンクリート構造物の非破壊検査論文集, Vol.4, pp.253 – 262, 2012065504-5, 2013
113.Mary

Clare

Sison

Escano ： Water

105.中畑和之，田中貴之：GPU 計算を用いたイメ

co-adsorption and electric field effects on

ージベース EFIT の高速化と超音波計測による

borohydride structures on Os(111) by first

解析妥当性の検証：土木学会論文集 A2(応用力

principles calculations：Journal of alloys and学), Vol.68, No.2, pp.I227 - I236, 2012

compounds, in press, 2013

106.K. Nakahata, K. Kimoto：Reconstruction of

114.大場史康，チェ・ミンソク，熊谷 悠，田中 功：

using

SrTiO3 における酸素空孔の特異な原子・電子

image-based FIT and time reversal approach:

構造：日本金属学会 2013 年春期(第 152 回)講Journal of Solid Mechanics and Materials

演大会概要, pp.270 – 270, 2013

flaws

in

heterogeneous

media

Engineering, Vol.6, No.6, pp.771 – 781, 2012

115.横山智康，藤田達也，林 博之，大場史康，田107.K. Nakahata, K. Terada, T. Kyoya, M.

中 功：第一原理計算による Cu2O 基中間バン

Tsukino & K. Ishii：Simulation of ultrasonic

ド型太陽電池の光吸収層の設計：日本金属学会

and electromagnetic wave propagation for2013 年春期(第 152 回)講演大会概要, pp.271 –nondestructive testing of concrete using

271, 2013

image-based FIT：Journal of Computational【5500 プロセス工学】

Science and Technology, Vol.6, No.1, pp.28 –116.山本卓也、高木洋平、伴 貴彦、岡野泰則、水37, 2012

口 尚：液膜内マランゴニ対流に及ぼすメニス

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 55

カス形状の影響に関する数値解析：第 44 回化
学工学秋季大会, pp.N108 - N108, 2012

125.建部勝利, 藤本 仁, 浜 孝之, 宅田 裕彦：移動平板上の水膜流へ衝突する複数棒状水噴流の117.T. Yamamoto, Y. Takagi, T. Ban, Y. Okano, H.

流動特性：CAMP-ISIJ, pp.1047 – 1047, 2012

Dost ： Numerical

126.齊藤翔造, 吹場活佳, 丸 祐介, 室伏洋希：微小investigation for the effect of meniscus

重力実験機の姿勢制御ジェットと超音速流の干

geometry on thermocapillary flow in a thin

渉に関する風洞実験：平成 24 年度宇宙航行のliquid film：23rd International Symposium on

力学シンポジウム講演集, 2012

Minakuchi,&

S.

Transport Phenomena, pp.43 – 43, 2012

【5700 基礎生物学】

118.山本卓也, 高木洋平, 伴 貴彦, 岡野泰則, 水口

127.Yutaka Satou, Tadasu Shin-i, Yuji Kohara,

尚：リング状液膜内マランゴニ対流に及ぼすメ

Nori Satoh & Shota Chiba ： A genomic

ニスカス形状の影響に関する数値解析：日本マ

overview of short genetic variations in a basal

イクログラビティ応用学会第 26 回学術講演会,

chordate, Ciona intestinalis：BMC Genomics,

pp.21A06 - 21A06, 2012

Vol.13, pp.208 – 208, 2012

119.山本卓也, 高木洋平, 伴 貴彦, 岡野泰則, 水口【6800 薬学】

尚：円形液膜内マランゴニ対流の三次元シミュ

128.東 巧, 小林祐輔, 竹本佳司：β－ターンミミッ

レーション：第 26 回数値流体力学シンポジウ

クを用いたアミノウレア触媒反応の機構解析：

ム, pp.C04-1 - C04-1, 2012

反応と合成のシンポシンポジウム, Vol.38,

120.山本卓也, 高木洋平, 岡野泰則：液膜上サーモキ

pp.208 – 209, 2012ャピラリー流れ方向に及ぼす液膜形状の影響に【9999 その他】

関する数値解析：化学工学会第７８年会講演要129.Takahiro Murashima, Masatoshi Toda, &Toshihiro Kawakatsu：Multiscale simulation

旨集, pp.F102 - F102, 2013
121.T. Yamamoto: A three-dimensional simulation

for soft matter: Application to wormlikeof thermocapillary flow in a circular thinmicellarliquid film under zero gravity：Proceedings ofProceedings, Vol.1518, pp.436 – 439, 2012

the 29th International Symposium on Space

130.Takahiro Murashima, Masatoshi Toda, andTechnology

andScience(ISTS-29),

Vol.2013-s-84-e, 2013
simulationofthermocapillary flow in a circular thinconcave and convex liquid film under
micro-gravity

：

：

AIP

Conference

Toshihiro Kawakatsu ： Multiscale Fluid
Dynamics Simulation Applied to Micellar

122.T. Yamamoto, Y. Takagi, & Y. Okano ：Three-dimensional

solution

Proceedings

of

8th

International OpenFOAM Workshop, 2013
123.T. Yamamoto, Y. Takagi, Y. Okano, & S.
Dost：Numerical investigation for the effect of
the liquid film volume on thermocapillary
flow direction in a thin circular liquid film：
Physics of Fluids, Vol . 25, pp.082108 –
082108, 2013
124.山本卓也、高木洋平、岡野泰則：静止液膜形状の影響を考慮した円形液膜内温度差マランゴニ
対流に関する数値解析：化学工学第４５回秋季
大会要旨集, Vol.K119, 2013
【5600 総合工学】

56 − Academic Center for Computing and Media Studies, Kyoto University

Solution：NIC Series, Vol.46, pp.187 – 192,
2013

資料

大型計算機システム利用負担金
別表１ スーパーコンピュータシステム
コース

タイプ

セット

タイプA
パーソナル タイプB
タイプC

基本
基本
基本
基本
最小

エントリ

タイプA1
タイプA2
タイプA3
タイプB1
グループ タイプB2
タイプB3
タイプC1
タイプC2
タイプG1
タイプA
大規模ジョブ タイプB
タイプC
専用クラスタ

―

追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

最小
追加単位

ライセンスサービス

提供サービス

利用負担額システム バッチ

12,600 円/年
100,000 円/年
100,000 円/年
100,000 円/年
200,000 円/年
200,000 円/年
240,000 円/年
120,000 円/年
600,000 円/年
300,000 円/年
250,000 円/年
250,000 円/年
300,000 円/年
150,000 円/年
750,000 円/年
375,000 円/年
400,000 円/年
200,000 円/年
240,000 円/年
120,000 円/年
250,000 円/年
250,000 円/年
20,000 円/週(7日）
5,000 円/週(7日）
24,000 円/週(7日）
6,000 円/週(7日）
20,000 円/週(7日）
10,000 円/週(7日）
750,000 円/年
375,000 円/年
20,000 円/年

システム資源

経過時間
（時間）

ディスク
（GB)

共有 最大1ノード相当((16コア、64GBメモリ)×1)
1
60
共有 最大4ノード相当((32コア、64GBメモリ)×4)
168
1,000
共有 最大4ノード相当((16コア、64GBメモリ)×4)
168
1,000
共有 最大2ソケット相当((8コア、384GBメモリ)×2)
168
1,000
4ノード((32コア、64GBメモリ)×4)
336
8,000
優先
4ノード((32コア、64GBメモリ)×4)
8,000
8ノード((32コア、64GBメモリ)×8)
336
9,600
A
準優先
4ノード((32コア、64GBメモリ)×4)
4,800
8ノード((32コア、64GBメモリ)×8)
336 16,000占有
4ノード((32コア、64GBメモリ)×4)
8,000
4ノード((16コア、64GBメモリ)×4)
336
8,000
優先
4ノード((16コア、64GBメモリ)×4)
8,000
8ノード((16コア、64GBメモリ)×8)
336
9,600
B
準優先
4ノード((16コア、64GBメモリ)×4)
4,800
8ノード((16コア、64GBメモリ)×8)
336 16,000占有
4ノード((16コア、64GBメモリ)×4)
8,000
4ソケット((8コア、384GBメモリ)×4)
336
8,000
優先
2ソケット((8コア、384GBメモリ)×2)
4,000
C
4ソケット((8コア、384GBメモリ)×4)
336
4,800
準優先
2ソケット((8コア、384GBメモリ)×2)
2,400
2ノード((16コア、64GBメモリ + 1GPU)×2)
336
4,000
B
優先
(GPU)
2ノード((16コア、64GBメモリ + 1GPU)×2)
4,000
8ノード((32コア、64GBメモリ)×8)A
占有
2ノード((32コア、64GBメモリ)×2)
8ノード((16コア、64GBメモリ)×8)
B
占有
2ノード((16コア、64GBメモリ)×2)
4ソケット((8コア、384GBメモリ)×4)
C
占有
2ソケット((8コア、384GBメモリ)×2)
8ノード((16コア、64GBメモリ)×8)
16,000B
4ノード((16コア、64GBメモリ)×4)
8,000
可視化ソフト（AVS,ENVI/IDL)およびプリポストウェアの1ライセンスにつき
B
A
B
C

利用者
番号

8
8
16
8
16
8
8
8
16
8
16
8
16
8
16
8
8
8
16
8

備考
1． 利用負担額は、年度単位で算定している。また、総額表示である。
2． 大型計算機システムの全ての利用者は、上記表のサービスの他、次のサービスを受けることができる。
1） 大判プリンタサービス
2） その他、大型計算機システムが提供するサービス、機器の利用
3． 上記表の大規模ジョブコース、ライセンスサービスの申請には、大型計算機システムの利用者であることが必要である。
4． 「共有」： 当該カテゴリのユーザ間で一定の計算資源を共有するベストエフォートのスケジューリングを行う。
「準優先」： 定常稼働状況において記載値（以上）の計算資源が確保されるように優先スケジューリングを行う。
また、稼働状況によらず記載値の1/4の計算資源が確保されることを保証する。
「優先」： 定常稼働状況において記載値（以上）の計算資源が確保されるように優先スケジューリングを行う。
また、稼働状況によらず記載値の1/2の計算資源が確保されることを保証する。
「占有」： 稼働状況によらず記載値（以上）の計算資源が確保されることを保証する。
5． ディスク容量はバックアップ領域（最大で総容量の1/2）を含む。
6． グループコース及び専用クラスタコースのシステム資源は、下記の負担額を支払うことにより増量することができる。
なお増量は各月１日に実施し、増量した資源は当該年度末までの期間にわたって利用されるものとする。
コース

タイプ
タイプA1
タイプA2
タイプA3
タイプB1

グループ

タイプB2
タイプB3
タイプC1
タイプC2
タイプG1

専用クラスタ―

追加負担金額（増量単位あたり）

20,000 円/月
12,000 円/月
30,000 円/月
25,000 円/月
15,000 円/月
37,500 円/月
20,000 円/月
12,000 円/月
25,000 円/月
37,500 円/月

システム資源増量単位
4ノード((32コア、64GBメモリ)×4)
4ノード((32コア、64GBメモリ)×4)
4ノード((32コア、64GBメモリ)×4)
4ノード((16コア、64GBメモリ)×4)
4ノード((16コア、64GBメモリ)×4)
4ノード((16コア、64GBメモリ)×4)
2ソケット((8コア、384GBメモリ)×2)
2ソケット((8コア、384GBメモリ)×2)
2ノード((16コア、64GBメモリ + 1GPU)×2)
4ノード((16コア、64GBメモリ)×4)

ディスク増量
（ＧＢ）

8,000
4,800
8,000
8,000
4,800
8,000
4,000
2,400
4,000
8,000

7． グループコース及び専用クラスタコースを通年でなく利用する場合には、下記の負担額を支払うものとする。
ただし、利用期間は当該年度内に限るものとする。
利用期間
最小
追加単位
最小
タイプA2
追加単位
最小
タイプA3
追加単位
最小
タイプB1
追加単位
最小
グループ
タイプB2
コース
追加単位
最小
タイプB3
追加単位
最小
タイプC1
追加単位
最小
タイプC2
追加単位
最小
タイプG1
追加単位
最小
専用クラス
―
タ コース
追加単位

タイプA1

3ヶ月
80,000 円
80,000 円
96,000 円
48,000 円
240,000 円
120,000 円
100,000 円
100,000 円
120,000 円
60,000 円
300,000 円
150,000 円
160,000 円
80,000 円
96,000 円
48,000 円
100,000 円
100,000 円
300,000 円
150,000 円

6ヶ月
120,000 円
120,000 円
144,000 円
72,000 円
360,000 円
180,000 円
150,000 円
150,000 円
180,000 円
90,000 円
450,000 円
225,000 円
240,000 円
120,000 円
144,000 円
72,000 円
150,000 円
150,000 円
450,000 円
225,000 円

9ヶ月
180,000 円
180,000 円
216,000 円
108,000 円
540,000 円
270,000 円
225,000 円
225,000 円
270,000 円
135,000 円
675,000 円
337,500 円
360,000 円
180,000 円
216,000 円
108,000 円
225,000 円
225,000 円
675,000 円
337,500 円

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 57

8． グループコース及び専用クラスタコースの利用者番号は利用者あたり年額5,000円を負担することで追加できる。
9 ． 機関・部局定額制度
他機関又は学内における部局（『国立大学法人京都大学の組織に関する規程』第3章第2節から第11節で定める組織をいう。）の組織が、その組織単位でグループコースサービス（年間）の利用を申請する場合、料金表(年間）に掲載額の1.5倍を利用負担金とする。なお、利用負担金額が150万円未満の場合は100人、150万円を超える場合は、150万円毎に100人までの利用者を認める。

別表２（汎用コンピュータシステム）
区分
VMホスティングサービス
ホームページサービス
ストリーミングサービス

利用負担額
72,000円／年
6,000円／年
6,000円／年

単位
1仮想マシンにつき1ドメイン名につき1申請につき

備考
１．利用負担額は、総額表示である。
２．上記表の汎用コンピュータシステムのサービスを利用するためには、大型計算機システムの利用者であることが必要である。
３．VM ホスティングサービスにおいて、下記の負担額を支払うことにより CPU、メモリ、ディスクを増量することができる。
区分
CPU増量
メモリ増量
ディスク増量

利用負担額
18,000円／年
18,000円／年
18,000円／年

単位
2コアにつき（最大8コアまで）
8GBにつき（最大64GBまで）
200GBにつき（最大1,000GBまで）

４．VM ホスティングサービスにおいて VMware を用いる場合は、下記の負担額を支払うことにより VMware の利用及び CPU、メモリ、ディスクを増量することができる。ただし、システム資源が限られているためサービスの提供を限定することがある。
区分
VMware利用
CPU増量
メモリ増量
ディスク増量

利用負担額
72,000円／年
36,000円／年
36,000円／年
18,000円／年

単位
1仮想マシンにつき2コアにつき（最大8コアまで）
8GBにつき（最大64GBまで）
200GBにつき（最大1,000GBまで）

５．ホームページサービス及びストリーミングサービスにおいて、下記の負担額を支払うことにより公開スペー
スの上限を拡大することができる。
区分
公開スペース20GBプラン
公開スペース50GBプラン

利用負担額
3,000円／年
9,000円／年

６．利用負担額は、当該年度（4月から翌年3月まで）の利用に対して年額として算定するが、年度途中から利用を開始する場合には月数に応じて減額する。

別表３ スーパーコンピュータシステム（民間機関利用）
システム資源システム
AB

8ノード(32コア、64GBメモリ)×8)
12ノード(32コア、64GBメモリ)×12)
16ノード(32コア、64GBメモリ)×16)
8ノード(16コア、64GBメモリ)×8)
12ノード(16コア、64GBメモリ)×12)
16ノード(16コア、64GBメモリ)×16)

経過時間
（時間）

ディスク
（GB)

利用者
番号

336
336
336
336
336
336

9,600
14,40019,2009,600
14,40019,200

16
24
32
16
24
32

利用負担額
960,000 円/年
1,440,000 円/年
1,920,000 円/年
1,200,000 円/年
1,800,000 円/年
2,400,000 円/年

備考
1． 利用負担額は、年度単位で算定している。また、総額表示である。
2． ディスク容量はバックアップ領域（最大で総容量の1/2）を含む。
3． 通年でなく利用する場合には、下記の負担額を支払うものとする。
ただし、利用期間は当該年度内に限るものとする。
システム

A

B

システム資源
8ノード
12ノード
16ノード
8ノード
12ノード
16ノード

3ヶ月
240,000 円
360,000 円
480,000 円
300,000 円
450,000 円
600,000 円

利用期間
6ヶ月
480,000 円
720,000 円
960,000 円
600,000 円
900,000 円
1,200,000 円

58 − Academic Center for Computing and Media Studies, Kyoto University

9ヶ月
720,000 円
1,080,000 円
1,440,000 円
900,000 円
1,350,000 円
1,800,000 円

資料

―
1.

サービス利用のための資料一覧―

スーパーコンピュータシステム・ホスト一覧
システム A：camphor.kudpc.kyoto-u.ac.jp



システム B・C：laurel. kudpc.kyoto-u.ac.jp


システム B（SAS 利用時）
：sas.kudpc.kyoto-u.ac.jp※ ホストへの接続は SSH(Secure SHell) 鍵認証のみ、パスワード認証は不可

2.

問い合わせ先 & リンク集


情報環境機構のホームページ
http://www.iimc.kyoto-u.ac.jp/



学術情報メディアセンターのホームページ
http://www.media.kyoto-u.ac.jp/



スーパーコンピュータシステムに関する問合わせ先

利用申請などに関する問合わせ先【共同利用掛】E-mail : zenkoku-kyo@media.kyoto-u.ac.jp / Tel : 075-753-7424
URL: http://www.iimc.kyoto-u.ac.jp/ja/services/comp/

システムの利用など技術的な問い合わせ【研究支援グループ】
E-mail : consult@kudpc.kyoto-u.ac.jp / Tel : 075-753-7426
URL: http://www.iimc.kyoto-u.ac.jp/ja/services/comp/contact.html


ホームページ・ホスティングサービス関する問合わせ先【情報環境支援グループ】
E-mail : whs-qa@media.kyoto-u.ac.jp / Tel : 075-753-7494
URL: http://www.iimc.kyoto-u.ac.jp/ja/services/whs/

コンテンツ作成支援サービスに関する問合わせ先【コンテンツ作成室】
E-mail : cpt@media.kyoto-u.ac.jp / Tel : 075-753-9012
URL: http://www.iimc.kyoto-u.ac.jp/ja/services/content/

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.12 No.2（ Dec. 2013）− 59

編集後記

編 集 後 記狭いながらも新居に住み始め 4 年目となりました。1 年目に縁があって甲斐犬の雌をペットショ
ップで見つけ飼い始めたのですが、翌年偶然同じ店で同犬種の雄を見つけ、我が家には２頭の甲斐犬が同居しています。甲斐犬は秋田犬紀州犬柴犬等と同じく日本古来から本土に住む日本犬です。日本犬は洋犬と対象的に保守的で飼い主以外の人間には慣れにくいとのことですが個体差もあるようです。身体の小さめな雌の方は警戒心が強く散歩のときもよく人に吠えますが、雄の方は対照的に人好きで初対面の人でも犬好きな人には馴れ馴れしく抱きついたりします。
先日犬の散歩中に公園で段ボール箱に捨てられた猫(スコティッシュ)を拾い飼うことになり、現在犬２頭猫１匹と同居の動物数が増えつつあります。
動物家族2013 年も 12 月になり、2014 年が近づいています。来年はサッカーファンにとってのビッグイベ
ント W 杯がブラジルで行われます。先日、組み合わせ抽選会が行われ、コロンビア、ギリシャ、
コートジボアールと同じ組になりました。早くも楽な組み合わせだとか、いやいや実際にはかな
り厳しいぞといった様々な意見がメディアやネット上で取り交わされています。確かに個人的な第一印象は、グループリーグ突破の可能性は高いぞと楽観的なものでしたが、よくよく考えてみると欧州の中でギリシャはフィジカルも強く、守備も堅いため、どちらかといえば日本が苦手とするタイプであり、コートジボアール、コロンビアとも圧倒的なフィジカルとスピードを有する
世界的選手が存在しており、とても楽観視はできないと思いなおした所です。日本も海外で活躍する選手を含め、技術とスタミナの両方に長けた選手を擁しており、他国からすると侮れない力をもっています。来年の W 杯での日本の活躍、期待したいです。
ロマニスタ

京都大学学術情報メディアセンター全国共同利用版広報 Vol．12，No．2
2013 年 12 月 20 日発行
広報編集部会
編集者 京都大学学術情報メディアセンター
広報教育委員会・全国共同利用版広報編集部会
発行者 〒606-8501 京都市左京区吉田本町

岩下 武史（部会長）
平石 拓 （副部会長）
秋田 祐哉

小西 満

京都大学学術情報メディアセンター

小林 寿

斎藤 紀恵

Academic Center for Computing and Media Studies

高見 好男

元木 環

Kyoto University
Tel. 075-753-7400
http://www.media.kyoto-u.ac.jp/
印刷所 〒616-8102 京都市右京区太秦森ヶ東町 21-10
株式会社エヌジーピー

60 − Academic Center for Computing and Media Studies, Kyoto University

表紙デザイン：谷 卓司（ティアンドティ･デザインラボ）