巻頭言

Vol.18, No.1 ߸ͷൃ‫ʹץ‬౰ͨͬͯ‫౎ژ‬େֶֶज़৘ใϝσΟΞηϯλʔɹ‫ڇ‬ౡ লຊ߸Ͱ͸ɼ‫౎ژ‬େֶֶज़৘ใϝσΟΞηϯλʔʹ͓͍࣮ͯࢪ͞Εͨฏ੒ 30 ೥౓ͷ‫ڞ‬ಉ‫ڀݚ‬ใࠂʹ͍ͭͯಛू͍ͨ͠·͢ɽฏ੒ 30 ೥౓͸ɼʮएखɾঁੑ‫ऀڀݚ‬঑ྭ࿮ʯͰ͸ 14 ݅ɼʮϓϩάϥϜߴ౓Խ‫ڞ‬ಉ‫ڀݚ‬ʯͰ͸ 5 ݅ͷ‫ڞ‬ಉ‫͞ࢪ࣮͕ڀݚ‬Ε·ͨ͠ɽֶज़৘ใϝσΟΞηϯλʔͰ͸ɼεʔύʔίϯϐϡʔλ‫ڞ‬ಉ‫ͯ͠ͱ౓੍ڀݚ‬ɼʮएखɾঁੑ‫ऀڀݚ‬঑ྭ࿮ʯɼ
ʮେ‫ن‬໛‫ࢧࢉܭ‬ԉ࿮ʯɼʮϓϩάϥϜߴ౓Խ‫ڞ‬ಉ‫ڀݚ‬ʯͱ͍͏̏छྨͷ‫ࢧڀݚ‬ԉࣄ‫ۀ‬Λ࣮ࢪ͍ͯ͠·͢ɽ՝୊ืू͸೥౓։࢝લ‫ߦʹޙ‬ΘΕ·͢ɽ·ͨɼ೥౓్த
Ͱ௥Ճืू͕࣮ࢪ͞ΕΔ৔߹͕͋Γ·͢ɽԠื͞Εͨ՝୊ͷ࠾൱͸ɼεʔύʔίϯϐϡʔλ‫ڞ‬ಉ‫اڀݚ‬ըҕһձʹͯ৹ࠪ͞Ε·͢ɽຊ‫ڞ‬ಉ‫΁౓੍ڀݚ‬ͷԠืํ๏΍ਃࠐΈ࣌‫ظ‬౳ͷৄࡉʹ͖ͭ·ͯ͠͸ɼhttp://www.iimc.kyoto-u.ac.jp/ja/services/comp/support/kyodo/Λ͝ࢀরԼ͍͞ɽ
্‫ه‬ͷ‫ڞ‬ಉ‫౓੍ڀݚ‬ͷ͏ͪɼʮएखɾঁੑ‫ऀڀݚ‬঑ྭ࿮ʯͷԠืࢿ֨͸ɼ(1) 40 ࡀະຬͷएख‫ֶ( ऀڀݚ‬ੜΛ‫ؚ‬Ή, ੑผ͸໰Θͳ͍)ɼ͋Δ͍͸ɼ(2) ঁੑ‫( ऀڀݚ‬೥ྸ͸໰Θͳ͍)ɼͱ͞Ε͍ͯ·͢ɽ͜ͷʮएखɾঁੑ‫ऀڀݚ‬঑ྭ࿮ʯͰ͸ɼεʔύʔίϯϐϡʔλΛར༻͢Δ͜ͱͰֶज़తʹΠϯύΫτ͕͋Δ੒ՌΛ૑ग़Ͱ͖Δ՝୊ʹରͯ͠ɼ‫ػࢉܭ‬ར༻ෛ୲ۚͷશֹ·ͨ͸Ұ෦Λຊηϯλʔ͕ෛ୲͍ͯ͠·͢ɽ͞Βʹɼʮएखɾঁੑ‫ऀڀݚ‬঑ྭ࿮ʯ͸ɼJHPCN (ֶࡍେ‫ن‬໛৘ใ‫ج‬൫‫ڞ‬ಉར༻ɾ‫ڞ‬ಉ‫ )఺ڌڀݚ‬ͷ‫׆‬ಈͷҰ‫ͱ؀‬Ґஔ͚ͮΒΕ͓ͯΓɼ༏Εͨ՝୊͸ɼJHPCN ͷ๖ժ‫ڞܕ‬ಉ‫΁ڀݚ‬ਪન͞Ε·͢ɽਪન͞Εͨ՝୊͸ɼ
JHPCN γϯϙδ΢ϜʹͯϙελʔൃදΛߦ͏͜ͱ͕ՄೳͱͳΔͳͲͷಛయ͕͋Γ·͢ͷͰɼੋඇ͜ͷ‫ڞ‬ಉ‫౓੍ڀݚ‬Λ͝‫͍ͩ͘͞༻׆‬ɽຊ߸ʹ‫͞ࡌܝ‬Εͨ 14 ݅ͷʮएखɾঁੑ‫ݚ‬‫ऀڀ‬঑ྭ࿮ʯͷ‫ڀݚ‬੒ՌใࠂͰ͸ɼෳࡶྲྀମ͓Αͼྲྀମߏ଄࿈੒໰୊ɼ΢Ϛͷ‫܈‬Εͷ‫ݸ‬ମؒͷྗֶϞσϧ‫ࢉܭ‬ɼίϯΫϦʔτͷഁյྗֶɼྔࢠԽֶɼ֤छͷ෼ࢠಈྗֶ‫ࢉܭ‬ɼ݁থੜ੒ͳͲʹؔ͢Δଟ༷ͳ‫ߦ͕ڀݚ‬ΘΕɼશࠃ 7 େֶͷएखɾঁੑ‫ऀڀݚ‬ͷ‫ڀݚ‬Λαϙʔτ
͍ͨ͠·ͨ͠ɽͳ͓ɼJHPCN ͷ๖ժ‫ڞܕ‬ಉ‫ʹڀݚ‬͸ 10 ݅ͷ՝୊͕࠾୒͞Ε·ͨ͠ɽ·ͨɼʮϓϩάϥϜߴ౓Խ‫ڞ‬ಉ‫ڀݚ‬ʯͰ͸ɼ‫౎ژ‬େֶͷεʔύʔίϯϐϡʔλΛάϧʔϓίʔε·ͨ͸ઐ༻ΫϥελίʔεͰ͝ར༻ͷ‫( ऀڀݚ‬JHPCN ͓Αͼ HPCI ʹͯ‫౎ژ‬େֶͷεʔύʔίϯϐϡʔλΛར༻͢Δ‫ऀڀݚ‬΋‫ؚ‬Ή) Λର৅ʹɼϓϩάϥϜͷ੍‫ߏޚ‬଄ɾσʔλߏ଄ͷվྑʹΑΔੑೳνϡʔχϯά΍ฒྻԽख๏ͷվྑɼ໰୊෼ׂɾෛՙ෼ࢄํࣜ
ͳͲͷվྑʹΑΔେ‫ن‬໛‫ࢉܭ‬ϓϩάϥϜͷߴ౓ԽɾߴੑೳԽͳͲΛࢧԉ͍ͯ͠·͢ɽฏ੒
30 ೥౓͸ɼ(1) ෼ࢠγϛϡϨʔγϣϯʹΑΔψΫϨΦιʔϜߏ଄มԽͷ໢ཏత୳ࡧɼ(2)‫ڊ‬େ௡೾૎্࣌ͷ໦଄Ո԰ͷ‫נ‬᛽ੜ੒աఔγϛϡϨʔγϣϯɼ(3) ҟํੑ஄ੑ೾ಈ໰୊ʹ
ର͢Δԋࢉࢠੵ෼࣌ؒྖҬ‫ڥ‬քཁૉ๏ͷߴੑೳԽɼ(4) ๞࿨౔ͷେ‫ن‬໛ม‫ܗ‬ɾྲྀಈ‫ࢉܭ‬Λ໨తͱͨ͠‫ݻ‬ӷࠞ߹ MPM ͷ։ൃɼ(5) ‫͍ڧ‬Ի‫ڹ‬ҟํੑΛ༗͢Δ CFRP ʹର͢Δ։‫߹ޱ‬੒๏ͷߴ଎࣮ߦɼͱ͍͏ 5 ݅ͷϓϩάϥϜߴ౓ԽࢧԉΛߦ͍·ͨ͠ɽ͜ΕΒͷ‫ڞ‬ಉ‫Ͱڀݚ‬͸ɼ֤छ‫ڥ؀‬ม਺ͷ‫౓ײ‬෼ੳɼஞ࣍‫ࢉܭ‬ίʔυͷεϨουฒྻԽͱߦྻ‫ࢉܭ‬ϥΠϒϥϦͷ
‫༻׆‬ɼϋΠϒϦουฒྻԋࢉʹ͓͚ΔύϥϝʔλαʔϕΠɼMPI ฒྻίʔυͷϋΠϒϦο
υฒྻԽɼϝϞϦΞΫηεͷޮ཰ԽͳͲ͕࣮ࢪ͞Εɼ10 ഒҎ্ͷߴ଎Խ͕࣮‫ͨ͠ݱ‬՝୊΋͋Γ·͢ɽ͜ͷΑ͏ʹɼଟ෼໺ʹΘͨΔղ๏ͷߴ౓ԽɾߴੑೳԽ͕ߦΘΕ͓ͯΓɼ‫ޮࢉܭ‬཰Λ޲্ͤ͞ΔࡍͷࢀߟʹͳΔͱࢥ͍·͢ɽࠓ‫ޙ‬΋օ༷ͷ‫ڀݚ‬ɼ‫ڭ‬ҭʹ͝‫͚͍ͩͨ༻׆‬ΔΑ͏ʹηϯλʔ‫ڭ‬৬һ΋ਚྗ͍͖ͯ͠·͢ͷͰɼ͝ར༻ɾ͝ࢧԉͷ΄ͲɺΑΖ͓͘͠‫͍ئ‬ਃ্͛͠·͢ɽ

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

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

機械攪拌操作時気泡巻き込みに対する数値解析山本 卓也
東北大学大学院環境科学研究科先端環境創成学専攻

1 緒言

数値シミュレーションでは、Navier-Stoke 式、連続式を解き、気液界面追跡手法としては代数的液体材料を混合する操作は多くの工業プロセス

Volume of Fluid (VOF)法を利用した。本解析条件で用いられ、金属生産プロセス、化学プロセス、

では撹拌 Reynolds 数は 18,000、24,000 であった

食品プロセス等でよく見られる。この混合プロセ

ため、乱流モデルとして LES を利用した。計算

スでは回転する翼を利用した機械的撹拌操作がよ

格子は翼付近、気液界面付近で細分化し、総格子く利用される。

点数はおおよそ 4,000,000 点であった。詳細な数本研究では金属生産プロセスでの機械攪拌操作

値シミュレーション手法に関しては、著者らの既に注目するが、鉄鋼プロセスでは KR（Kanbara

報に記載されている[2]。今回のシミュレーション

Reactor）法と呼ばれる脱硫プロセス、アルミニ

は全てオープンソースであるOpenFOAM を利用

ウム生産プロセスでは脱アルカリ金属等のためにした。

利用されるフラックス処理、脱水素処理、均一な合金生成のための攪拌[1]等に利用される。KR 法
では溶銑表面に散布した CaO 粒子を巻き込み、3 結果と考察

硫黄をスラグ中へと分離する。また、アルミニウ

本研究では数値シミュレーション結果を検証す

ム生産プロセスではアルミニウムは空気中の水蒸るため、高速度ビデオカメラによる撮影と気泡巻

気と反応して不純物となる皮膜を溶融金属表面にき込み条件の測定を行った。中心、偏心のそれぞ形成するため、これを巻き込まずに攪拌する必要

れで気泡巻き込み回転数はそれぞれ 390 rpm,がある。上述のプロセスでは、機械攪拌中での液460 rpm であった。数値シミュレーションでは、

面振動メカニズムが重要となる。

中心位置の 400rpm のみで気泡を巻き込んだ。こ

本研究では、機械撹拌操作中での気液界面巻きれは定性的に一致している。さらに、高速度ビデ込み現象のメカニズムを解明するため、VOF 法を

オカメラにより得られた界面振動位置と数値シミ利用して中心、偏心撹拌時の液面変動メカニズムュレーションによって得られた界面振動位置は良を解明した。

好に一致した。このため、数値シミュレーション
は良好に実験を表していると言える。

2 数値解析手法

回転数 300 rpm, 400 rpm の場合の中心、偏心位置における時間平均化した界面高さ分布を Fig.

内径 192 mm、深さ 400 mm の円筒容器内に水

1 に示す。中心撹拌の方が大きく液面が変形し、を 320 mm まで満たし、6 枚パドル翼を液深 40

軸中心でより大きく凹んだ形状となる。一方で偏mm まで浸漬させた。翼端直径は 60 mm である

心撹拌の場合、液面の平均位置は中心撹拌と比べとし、今回利用したパドル翼は溶融金属系への適ると平坦な分布となる。特に、側壁と撹拌軸の間用を目指し、翼厚が 6 mm とした。翼は 300 rpm,

で液面が盛り上がり、より平坦な液面高さ分布へ

400 rpm で回転させた。翼を円筒容器軸中心位置と変化する。このように平均液面高さが変化したと 60 mm 偏心した位置に配置した。

理由を考察するため、偏心撹拌の場合の Fig. 1 の

2 − Academic Center for Computing and Media Studies, Kyoto University点線上の断面平均速度ベクトルを Fig. 2 に示す。

る。

撹拌翼から吐出された強い流れは側壁と衝突し、上昇流へと変換される。この側壁付近の上昇流が気液界面を持ち上げる。このため、Fig. 1 に示した撹拌軸と側壁間の液面盛り上がりへと繋がり、
気液界面変動幅が小さくなる。

Fig. 3 Snapshots of free-surface height in thecase of on-axis stirring with 400 rpm rotation.

Fig. 1 Time-averaged height of free surface withdifferent impeller location and rotation speed:(a) on-axis 300 rpm, (b) on-axis 400 rpm, (c)eccentric 300 rpm, and (d) eccentric 400 rpm.
Fig. 4 Snapshots of free-surface height in thecase of eccentric stirring with 400 rpm rotation.

4 結言
本研究では機械撹拌時の液面からの気泡巻き込
みメカニズムを解明するため、中心撹拌、偏心撹
拌時の液面変動挙動を解明した。撹拌方法によって液面の平均位置が異なり、動的な変動挙動も全く異なることがわかり、それぞれのメカニズムをFig. 2 Time-averaged cross-sectional velocity
vectors generated by eccentric stirring (a), (c)viewed from point A in Fig. 1, (b), (d) viewed
from point B in Fig. 1. The rotation speeds are(a), (b) 300 rpm and (c), (d) 400 rpm.
中心、偏心撹拌時の動的な界面挙動を Figs. 3, 4に示す。中心撹拌の場合には大きく変動し回転する部分が二箇所存在するものの、偏心撹拌の場合には、そのような大きな変動は存在せず、小さい
ランダムな変動が発生する。これは、中心撹拌ではマクロ不安定性の効果が大きくなり、偏心撹拌
ではマクロ不安定性が発生しにくくなるからであ解明することができた。

引用文献[1] V. S. Warke et al., J. Mater. Process. Technol. 168(2005) 119-126.
[2] T. Yamamoto et al., Chem. Eng. J. 367 (2019)25-36.

謝辞
本研究制度(若手・女性研究者奨励枠)を活用させて頂いたことを、この場を借りて厚く御礼申し上げる。

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

大規模溶解炉における溶融アルミニウム機械撹拌時の流動解明
加藤 賢也
東北大学大学院環境科学研究科先端環境創成学専攻

1

緒言

面に付着した Na+イオンと HCl 水溶液中の H+イオンの交換反応が生じ、水溶液の pH が変動する。

国内アルミニウム産業における最初の生産プロ

pH の変動を Fig. 1 中に示された 4 箇所(A-D 点)でセスは、輸入したアルミニウム地金を溶解炉で溶測定し、pHの変化から物質移動係数を導出した[3]。

かし、その後不要元素を除去し、合金成分を調整

各測定点での物質移動係数の平均値を槽の物質移するための溶湯処理プロセスである。この溶湯処動係数と定義して時計方向 (CW)、半時計方向

理プロセスでは、機械撹拌を行いながら溶融アル(CCW)回転の結果を比較した。

ミニウム中へとハロゲン化物を含む溶融塩、もしくは、塩素ガスが導入され、不純物を浮上分離させる[1]。この溶湯処理は大規模、不透明かつ高温であるため、その実験による観測が困難である。
溶融処理の一つにフラックス処理という不純物除去処理がある。ここでは、フラックスと呼ばれる粉末状の固体粒子を溶融アルミニウム中へと散布し、不純物と固体粒子が反応し、浮上させる。
この反応効率を向上させることが現場技術では求められているが、上記の問題より実験的な観測が困難である。
本研究では、水を用いたモデル実験と数値シミュレーションを活用し、効率の良いフラックス処
理条件を調査した。

2

実験手法、数値解析手法
実機をスケールダウンした水モデル実験装置を作成し、NaOH と HCl の反応を利用した物質移動実験を用いた。Figure 1 に実験装置の概略図を示す。図中に示すような水槽に HCl 水溶液(3.5molL-1)を満たし、攪拌翼を槽の端に配置し、一定速度で回転させた。ここで利用した翼形状は、ア
ルミニウム産業で利用されているものである[2]。
その上に、NaOH 水溶液(2.0 molL-1)に 72 時間浸した後に乾燥させたパーライト粒子(三井金属, 比
重 0.21)を 170 g 散布した。このパーライト粒子表

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

Fig. 1 Schematic representation of water modelexperiment [4].
槽内の流動を可視化するために、同時に数値解析を行った。槽内の流動のみに着目し、Navier-Stokes 式、連続式を解いた。ここで、攪拌Reynolds 数は 30000 (500 rpm), 60000 (1000 rpm)であったため、乱流モデルとして Reynolds Averaged
Navier Stokes (RANS), k-w SST モデルを利用した。
攪拌翼の回転を表現するために、Multi Reference
Frame (MRF)法を利用し、翼付近のみ回転領域であるとした。数値解析手法の詳細は既報[2]と同様である。総計算格子点数は約 125 万点であり、計算は全てオープンソースである OpenFOAM を利用して行われた。

3

結果と考察
水モデル実験によって得られた物質移動係数をFig. 2 に示す。回転数が増大するにつれ物質移動係数は増大し、回転方向が CCW の方が CW より大きい。
このようになった原因を調査するために、

モデル粒子の巻き込みが大きく変化したからであ

数値解析結果より考察する。

る。

Fig. 2 Mass transfer coefficient with different rotation
direction and rotation speed of impeller [4].

Fig. 3 Distribution of turbulent kinetic energy on the

パーライト粒子は比重が小さいため、流れの影
響がなければ自由表面上に浮遊する。これが本実

free surface: (a)-(c) CW, (d)-(f) CCW, (a),(d) 500 rpm,(b), (e) 700 rpm, (c), (f) 1000 rpm [4].

験中では槽内部へと巻き込まれた。このような状
況下で物質移動は、1. 乱流変動による巻き込み、
2. 巻き込まれた粒子表面での物質移動が関与すると考えられる。ここでは、最初に関わる現象である自由表面上での現象に注目する。数値解析により得られた自由表面における乱流エネルギー分
布と速度分布を Figs. 3,4 に示す。CW では乱流エネルギーが小さく、CCW では大きくなる。また、回転数の増大に伴って自由表面上の乱流エネルギ
ーが増大する。速度分布の場合には、回転数の増大に伴って速度が増大するものの、速度の絶対値
としては回転方向に対しては大きく変化しない。
このため、乱流エネルギーが要因で物質移動係数の回転方向依存性が引き起こされたと考えられる。
このため、CCW においては自由表面上で大きい
乱流速度変動によって粒子巻き込みが促進され、槽内部では粒子表面における濃度境界層が薄くなり、物質移動が増大することで、槽全体の物質移
動係数が増大したと考えられる。これらの議論の詳細に関して既報の論文[4]を参照して頂きたい。

4

Fig. 4 Velocity distribution on the free surface: (a)-(c)CW, (d)-(f) CCW, (a),(d) 500 rpm, (b), (e) 700 rpm, (c),(f) 1000 rpm [4].

引用文献
[1] M. E. Schlesinger, Aluminum Recycling, second
edition (CRC press, Florida 2014).
[2] T. Yamamoto, K. Kato, S. V. Komarov, Y. Ueno, M.
Hayashi, Y. Ishiwata, J. Mater. Process. Technol. 259(2018) 409-415.
[3] R. Shiba, M. A. Uddin, Y. Kato, S. Kitamura, ISIJ
Int. 54 (2014) 2754-2760.

結言[4] K. Kato, T. Yamamoto, S. V. Komarov, R.

本研究ではアルミニウムの溶湯処理の一つである機械撹拌操作を伴うフラックス処理を模擬したモデル実験と数値シミュレーションを行った。実験結果によると、撹拌翼の回転方向によって物質移動係数が大きく変化した。これは回転方向の違いによって気液界面での乱流変動が大きく変化し、Taniguchi, Y. Ishiwata, Mater. Trans. 60 (2019)2008-2015.

謝辞
本研究制度(若手・女性研究者奨励枠)を活用させて頂いたことを、この場を借りて厚く御礼申し上げる。

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

界面分光の分子動力学シミュレーション
城塚 達也
茨城大学工学部物質科学工学科

1 はじめにことにより分子メカニズムを解析する。界面構造とダイナミクスとの相関やイオン効果などは実験

固体と液体の界面に代表される埋もれた界面はでも観測されているため、本研究によりそれらの学術的・工業的に極めて重要であるが、その分子分子レベルでの描像を明らかにする。一般に古典レベルでの反応メカニズムには未だに多くの謎がMD シミュレーションでは非分極の分子モデルが存在する。そこで、近年発展が著しい和周波発生

多用されるが分光スペクトルの再現は難しいため、（SFG）分光法を申請者がこれまで開発してきたこの計算では高精度な分極モデルを用いることに手法により解析し、界面現象を明らかにする。特より実験観測量の議論が可能となる。

に、本研究ではシリカ・水界面と酸化チタン・水

このような界面系では系の大きさを十分大きく界面における構造やダイナミクスを調べる。古とらないと実験結果と比較し検討することはでき

典・第一原理分子動力学（MD）シミュレーショ

ないため、全原子数は約数千原子になる。また、ンを用いて界面分光の１つである SFG 分光法を

水分子には分極モデルを採用し、シミュレーショ

解析する。

ン時間は数十 ns になる。計算には高並列化され2 手法
まず SFG スペクトルの計算では申請者が開発

た in-house コードを用いた。

3 結果

した手法を用いた。[1-3] 固体・液体界面を生成

初めに、パルミチン酸単分子膜と水溶液界面のする必要がある場合には、第一原理（DFT）計算

解析を行った。[1] この研究は和周波発生（SFG）

を用いて水との界面を生成し表面構造を実験と比分光法における電気二重層の重要性をアメリカの較する。実験では表面にどれだけの分子が吸着し実験家と理論家との共同研究で初めて明らかにすているか不明な部分が多いため、ここではまず静ることができた。実験は理論的に計算が困難な表的な界面構造を議論する。酸化チタン・水界面に面ポテンシャルを測定し、理論は実験では解明がおける電荷分離状態を電子状態計算により再現し、難しい電気二重層の重要性と界面構造を明らかに

最終的に光触媒や超親水性の分子メカニズムを議することができた。そのため、実験と理論が相補論する。

的な役割を果たすことができた。このように、界

その後、得られた構造で高精度の量子化学計算

面での水溶液の構造・ダイナミクスを解明するこを行い分子動力学（MD）シミュレーションでのとによって、化学反応に与える影響や界面での性

分子力場を構築する。これらの計算では cp2k や質に関する知見を得ることができた。

Gaussian などの汎用計算プログラムを用いる。

また、SFG 分光法における電気二重層の重要性

得られた力場を使い SFG スペクトル、水の拡散を一般的な電解質水溶液に拡張し、実験ではこれ

係数や再配向時間などのダイナミクスを古典 MDまで明らかなっていなかった電気二重層効果のイシミュレーションにより計算し、実験と比較する

オン濃度依存性を始めて明らかにすることができ

6 − Academic Center for Computing and Media Studies, Kyoto Universityた。[2, 3] 本研究の成果は今後 SFG 分光を他の固4 おわりに体と液体界面に応用する際に基礎的なデータ・知
荷電した界面における SFG スペクトルを計算見を提供することができた。
これらに加えて、金属表面の化学反応に関するする手法を開発し、パルミチン酸単分子膜と水溶実験との共同研究も実施した。本研究では主に液界面やシリカ・水界面などに応用した。この研VASP や cp2k などのパッケージプログラムを使究により得られた知見は脂質膜、シリカのような用 し 、 実 験 的 に は 明 ら か で な か っ た CVD

無機化合物、金属酸化物、燃料電池の固体高分子(chemical vapor deposition) における律速段階

膜など幅広い系にとって有益であり、実験の SFGの解析やその分子メカニズムを解明した。[4, 5]の測定はすでになされているのでこの研究のアプ（Figure 1 top）最後に、上記と同様の手法を用ローチが直接適用できる。本研究の解析は荷電しいることにより酸化チタン・水界面における電荷た界面では本質的であり、電気二重層など界面の分離状態の電子状態計算の基盤を構築することが本質に迫る事ができ、分子科学・表面科学としてできた。
（Figure 1 bottom）

も非常に意義深い。また、銅の CVD 法による堆積のシミュレーションを行い実験的に観測された律速段階を明らかにすることができた。今後は、酸化チタン・水界面での水分子の役割や異なる界面における反応性の違いなどを解析する予定である。

謝辞
本研究では、京都大学学術情報メディアセンタ
ーの共同研究制度（若手奨励枠）を活用させていただきました。

参考文献[1] S. K. Reddy, R. Thiraux, B. A. W. Rudd, L.
Lin, T. Adel, T. Joutsuka, F. M. Geiger, H. C.
Allen, A. Morita, and F. Paesani, Chem, 4 (7),1629–1644, 2018.
[2] T. Joutsuka, T. Hirano, M. Sprik, and A.
Morita, Phys. Chem. Chem. Phys., 20 (5),3040–3053, 2018.
[3] T. Joutsuka and A. Morita, J. Phys. Chem. C,
122 (21), 11407–11413, 2018.
Figure 1. (Top) Schematic of CVD processes
considered in Ref. [4]: Adsorption of precursors,migration of adsorbates and desorption of[4] T. Joutsuka and S. Yamauchi (原稿準備中).
[5] T. Nishikawa, K. Horiuchi, T. Joutsuka andS. Yamauchi (投稿中).

by-product. (Bottom) Polaron (excess electron)
simulated in this study at the anatase (101)surface.

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

高効率有機系太陽電池の実現に向けた光機能性分子の構造と電子物性
の相関解明
東野	
  智洋
京都大学大学院工学研究科分子工学専攻

1 緒言

アジアゾール基を用いた例は無かった。そこで本研究では、光捕集能のさらなる向上を狙い、XW17

ポルフィリンは 400~450 nm に Soret 帯と呼ばれの構造にナフトビスチアジアゾール基を導入した

る強い吸収と 550~600 nm に Q 帯と呼ばれる中程

色素 ZnPNTz を設計・合成した（Figure 1）[3]。

度の吸収をもち、増感色素として有望である。しかし、単純な構造のポルフィリン色素では長波長

2

結果と考察

の太陽光を効率よく利用することができないため、
光捕集能を改善させる分子設計が必要である。近年ではプッシュ-プル構造を導入することによって、10%を超える変換効率が実現されてきている。
最近、ベンゾチアジアゾールを用いたプッシュプル型の色素 XW17 を増感色素として用いた系
で 9.5%という高い変換効率が達成されている[1]。
また、有機薄膜太陽電池において、ベンゾチアジ
アゾールよりも電子求引性の高いナフトビスチア
ジアゾールを用いることで光捕集能が改善し、変換効率も向上することが報告されている[2]。しかし、色素増感太陽電池の増感色素にナフトビスチ
C6H13O

OC6H13

C6H13O

OC6H13NN
OC6H13

C6H13OC12H25O

OC6H13

色素増感太陽電池セル作製条件の最適化を行った結果、光電変換効率は XW17 を用いたセルで8.14%となったのに対し、ZnPNTz を用いたセル
では 3.63%となった。光電変換効率は XW17 よりも低くなったが、XW17 を用いたセルで光電流発生が見られたのが 800 nm 程度までであったのに対し、ZnPNTz を用いたセルでは 860 nm 付近まで光電流発生が見られた（Figure 2）。これは、ナフトビスチアジアゾール基の高い電子求引性によってプッシュ−プル特性が強められたことで
ZnPNTz の吸収ピークが長波長化していたためで

C6H13ON

C12H25ON

C12H25ON

Zn
O
C12H25

2.1 新規ポルフィリン色素の物性と太陽電池性能評価N

C12H25
O

N
Zn

N

O
C12H25NN

NO
C12H25NO
C12H25

N

S
S

N

N

S

N

COOH
XW17

COOH
ZnPNTzFigure 1. Molecular structures of porphyrin dyes.

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

Figure 2. Photocurrent action spectra of the DSSCs
based on ZnPNTz (solid line) and XW17 (dashed line)
under the best conditions.

あると考えられる。したがって、ナフトビスチア

よりも ZnPNTz で小さくなっていた。カルボキシ

ジアゾール基を活用することによって、近赤外光

基の LUMO の軌道分布が小さく、励起状態の色を効率よく利用可能できる色素を開発可能である

素からの電子注入の効率が低くなってしまったとと言える。

考えられる。このように、本計算結果は実験結果の理論的解釈の一助となった点で意義がある。

2.2 理論計算
ポルフィリン色素の最安定化構造およびそのフロンティア軌道の電子構造について知見を得るた
めに、Gaussian09 プログラムを用いて密度汎関数
法 （ DFT ） に よ る 理 論 計 算 を 行 っ た（B3LYP/6-31G(d)）
。Figure 3 にポルフィリン色素
ZnPNTz およびXW17 の基底状態での構造および最高被占軌道（HOMO）と最低空軌道（LUMO）
における軌道分布を示す。トリアリールアミン、
ナフトビスチアジアゾール基はどちらもポルフィリンと同一平面を取っており、共役系が効果的に広がっていることがわかった。また、ZnPNTz のLUMO が大きく安定化することによって HOMO–
LUMO ギャップが小さくなっており、吸収ピークの長波長化が見られたことと一致した。一方で、

3

参考文献[1] Tang, Y.; Wang, Y.; Li, X.; Ågren, H.; Zhu,

W.-H.; Xie, Y. ACS Appl. Mater. Interfaces 2015, 7,27976–27985.
[2] (a) Vohra, V.; Kawashima, K.; Kakara, T.;
Koganezawa, T.; Osaka, I.; Takimiya, K.; Murata, H.
Nature Photonics 2015, 9, 403–408. (b) Kawashima,K.; Fukuhara, T.; Suda, Y.; Suzuki, Y.; Koganezawa,
T.; Yoshida, H.; Ohkita, H.; Osaka, I.; Takimiya, K. J.
Am. Chem. Soc. 2016, 138, 10265–10275.
[3] Higashino, T.; Kurumisawa, Y.; Nimura, S.;
Iiyama, H.; Imahori, H. Eur. J. Org. Chem. 2018,
2537–2547.

カルボキシ基の LUMO での軌道分布は、XW17

Figure 3. Selected Kohn−Sham orbitals for (a) ZnPNTz and (b) XW17, obtained by DFT calculations with the
B3LYP/6-31G(d) level. To simplify the calculations, alkoxy groups were replaced with methoxy groups.

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 9無重力下での高プラントル数流体における温度差マランゴニ効果に起因する
液柱内対流場の二次不安定性
小笠原 亨
東京理科大学大学院理工学研究科機械工学専攻

1 緒言
自由界面を有する気液界面において，温度分布が存在する場合，表面張力の温度依存性により表ロケ理論によって Pr = 4 の二次不安定性におけ(/)

る臨界レイノルズ数 Re- を明らかにした．しかし，非線形解析では対流等は明らかにされていな
い．本研究では二次不安定性に関して，直接数値

面張力が不均一となり，対流が生じる．この対流計算および固有直交分解（以下，POD）17)を用いはマランゴニ対流と呼ばれ，微小重力環境下やマて解析を行った．

イクロスケール下で顕在化する．これは自然対流

私が所属している研究室は，国際宇宙ステーシの駆動力の一つとして着目されている．例えば，

ョン日本実験棟「きぼう」
によって行われていた，単結晶材料生成方法の一つである floating-zone（以MEIS ( Marangoni Experiment in Space ) や Dynamic下，FZ）法は，純度の高い単結晶材料を生成する

Surf. と称される流体物理実験に参画していた．そことができる．一方で，材料生産の高効率化のたして，今後，JEREMI ( Japanese European Research

めに，微小重力下で利用するとマランゴニ対流のExperiment on Marangoni Instabilities ) と称される，不安定性によって生じる振動流が材料生成に悪影

流体実験を行う予定である． JEREMI では，気体

1)響を及ぼすことが知られている ．したがって，側に強制対流を付与し，自由界面での熱伝達制御

マランゴニ対流の不安定性に関する理解や制御のを行う予定である．

観点で研究が行われている．
図１は half-zone（以下，HZ）モデルと呼ばれる，FZ 法の半分を模擬した，最も典型的なモデルであ本研究では，マランゴニ対流の不安定化メカニ
ズムの解明による理学的貢献，及び産業分野への波及効果を目指し，研究を実施している．

る．上下壁面間に液柱を形成し，上部を加熱，下部を冷却することで液柱自由表面上に温度勾配を付与し，上から下へマランゴニ対流を発生させる．

z*

T* = T*h

この HZ モデルにおけるマランゴニ対流の不安定
性 に プ ラ ン ト ル 数 （ 以 下 ， Pr = 𝜈𝜈/𝜅𝜅, 𝜈𝜈 ∶動粘性係数，𝜅𝜅 ∶ 熱拡散率 ）が影響を与えることが一般的に知られている 2) - 4)．まず，上下壁間のT* = T*cr*

温度差が小さい時，液柱内対流場は二次元定常流になる．低 Pr 流体の場合，温度差を大きくすると，
対流場は３次元振動流になり（一次不安定性），その後三次元振動流に遷移する（二次不安定性）2) - 11)．
高 Pr 流体においては，二次元定常流から三次元振動流に一気に遷移する（一次不安定性）2), 7), 10)．この対流後の高次不安定性およびカオス化過程においても研究されている 12) - 15). 茂木（2017）16) はフ

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

Fig. 1 Schematics of half-zone geometry

2 計算方法液柱の高さを H，半径を R とし，アスペクト比Γ = 𝐻𝐻⁄𝑅𝑅 と定義する（図１）．液柱は動的・静的な変形を考慮せず，常にストレートな形状を保持しているものとする．計算系全体は無重力環境下を想定している．上下壁間の温度差は ∆𝑇𝑇 = 𝑇𝑇7 −𝑇𝑇- と定義し，表面張力 𝜎𝜎(𝑇𝑇)は温度に依存し，線型的な変化をするものと仮定する．

𝝈𝝈(𝑻𝑻) = 𝝈𝝈(𝑻𝑻𝐜𝐜 ) + 𝝈𝝈𝑻𝑻 (𝑻𝑻 − 𝑻𝑻𝐜𝐜 )(1)

ここで，𝜎𝜎? = 𝜕𝜕𝜕𝜕⁄𝜕𝜕𝜕𝜕 ．流体は非圧縮性ニュート

大きくなるように定義している．また，各成分の温度変動 _𝑇𝑇Z を以下のように定義する．

3 結果

𝑇𝑇XZ = 𝑎𝑎/Z 𝜙𝜙/Z + 𝑎𝑎/Z\] 𝜙𝜙/Z\](7)

周方向波数（以下，m）が m = 3 で周期的な回転振動流であった対流場は Re = 3250 において準ン流体であることを仮定し，Pr が 4 であるものと周期回転流へと遷移した．図２(a) - (c)は FFT によする．無次元化された流体の運動と熱輸送を記述って得られた回転周期を元に準周期回転流を回転する支配方程式は，Navier-Stokes 方程式，連続の座標系に変換し，温度変動の等温面を可視化した式，エネルギー方程式とする．表面張力を除く物

様子を時系列に示している．等温面の形および大性の温度依存性は無視する．

きさが時間によって変形していることが確認できる．したがって，これを二次不安定性後の対流場

𝝏𝝏𝒖𝒖∗
𝟏𝟏 𝟐𝟐 ∗+ (𝒖𝒖∗ ∙ 𝛁𝛁)𝒖𝒖∗ = −𝛁𝛁𝒑𝒑∗ +
𝛁𝛁 𝒖𝒖
𝝏𝝏𝝏𝝏
𝐑𝐑𝐑𝐑

(2)

𝜕𝜕𝑇𝑇 ∗1 / ∗+ (𝒖𝒖∗ ∙ ∇)𝑇𝑇 ∗ =
∇ 𝑇𝑇
𝜕𝜕𝜕𝜕
Ma

(4)

∇ ∙ 𝒖𝒖∗ = 0(3)として，以下に結果を記述する．
二次不安定性前の温度場を POD によって分解した結果，エネルギーの大きい順に成分は m = 3, 6,
9 を示した．これは元の対流場の m および高調波成分である．図２(d) - (f)は二次不安定性後の温度場を POD によって分解した各成分の温度変動の

ここで，u*を無次元速度ベクトル，𝑝𝑝∗ を無次元圧

等温面の瞬時場を可視化した様子である．一番エ

座標系(r, q, z)を用いている．これより以下におい

ネルギー量の大きい成分（k = 1）は二次不安定性

力，𝑇𝑇 ∗ = (𝑇𝑇 − 𝑇𝑇T )⁄∆𝑇𝑇を無次元温度とする．円柱ては*を省略する．無次元数 Re は対流の強さを測るパラメーターであり，以下のように定義する．
Re =|𝜎𝜎? |∆𝑇𝑇𝑇𝑇𝜌𝜌𝜈𝜈 /

前同様 m = 3 であることが確認できる．しかし，次にエネルギーの大きい成分は m = 2, 4 の形を示している．したがって，これらの成分が生じるこ(5)とで二次不安定性が発生したことがわかる．
また，

ここで，r は密度，n は動粘度，Ma をマランゴニ

この対流場からは高調波成分も得られているが，数とすると，Ma = Re Pr という関係がある．

エネルギー値が新しく生じた成分のエネルギーよ

境界条件は上下壁上では，
速度はノンスリップ，りも小さいことが確認されている．

温度は高温壁・低温壁で Th および Tc でそれぞれ一定としている．液柱自由表面上での熱の授受はないものとし，断熱の境界条件を与え，速度に関しては，表面張力とせん断応力の釣り合いから，
マランゴニ効果による流体の駆動を行っている．
計算により得られたデータから時間平均に対する温度変動分 𝑇𝑇Xを算出し，一般的に知られている(b) t = t0 + ⌧ /2(c) t = t0 + ⌧(e) k = 2(m = 2)(f) k = 3(m = 4)

POD16)を用いて解析を行う．POD からは固有値l, 固有ベクトル a, 固有関数 f，が得られる．周期関数を POD した際，同じ固有関数が２つずつ得られることが一般的に知られているため，各固有関数の成分の固有値を足し合わせてエネルギーEk
を以下のように定義する．
𝐸𝐸Z = 𝜆𝜆/Z + 𝜆𝜆/Z\](6)

ここで k は正の整数であり，k が小さいほど Ek は(d) k = 1(m = 3)

Fig. 2 (a)-(c) 3D visualization of thermal-fields in rotating
frame of reference against fundamental frequency ofhydrothermal wave after secondary instability (d)-(f) Largeenergy components obtained by POD

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

図３は二次不安定性後の温度場を POD によっ4 結言て得られたエネルギーの大きい７つの成分の各 m
に対するパワースペクトルを示している．この三Pr = 4 の流体で，無重力下でのストレートな液次元グラフからも二次不安定性前には確認されて柱（Γ= 0.7）を想定して，二次不安定性の遷移過いない m = 3 および高調波以外の m のピークが程を調査するため直接数値計算を行った．二次不存在していることが明確である．

安定性後の対流場は準周期回転流になることを確認した．この対流場を POD によって解析し，二次不安定性前には確認されなかった m = 3 および高調波以外の成分が得られた．この成分のエネル
ギー値の Re 数に対する発達過程から臨界レイノ
ルズ数 Rec(2)を求めた．

5 参考文献Fig. 3 Full spectral decompositions of surface
temperature deviation at 𝒛𝒛 = 𝚪𝚪⁄𝟐𝟐 for each component
after secondary instability

図４に二次不安定性後に生じた m = 2 および 4のエネルギー値 Ek を Re の変化とともに示す．Re1)

A. Cröll et al.: J. Cryst. Growth, 191 (1998) 365.

2)

M. Wanschura et al.: Phys. Fluids, 7 (1995) 912.

3)

I. I. Ryzhkov: Phys. Fluids, 23 (2011) 082103.

4)

K. Fujimura: J. Phys. Soc. Jan., 82 (2013) 074401.

5)

R. Rupp et al.: J. Crystal Growth, 97 (1989) 34.

6)

M. Levenstam and G. Amberg: J. Fluid Mesh., 297 (1995)
357.

7)

J. Leypoldt et al.: J. Fluid Mech., 414 (2000) 285.

8)

M. Levenstam et al.: Phys. Fluids, 13 (2001) 807.

9)N. Imaihi et al.: J. Crystal Growth, 230 (2001) 164.

10)

C. Nienhüser and H. C. Kuhlmann: J. Fluid Mech., 458(2002) 35.

11)

K. Motegi et al.: Phys. Fluids, 29 (2017) 074104.

12)

R. Velten et al.: Phys. Fluids A, 3 (1991) 267.

13)

V. M. Shevtosova et al.: Phys. Review E, 68 (2003)
066311.

14)

I. Ueno et al.: Phys. Fluids, 15 (2003) 408.

15)

T. Matsugase et al.: Int. J. Heat Mass Trans., 89 (2015)
903.

16)

K. Motegi et al.: Annual Meeting of the Japan Society of
Fluid Mechanics 2017 (2017b).

17)

K. Li et al.: J. Crystal Growth, 307 (2007) 155.

が大きくなるにしたがって，それぞれのエネルギ
ーは大きくなることが確認できる．Re < 3200 では
m = 2, 4 の成分の存在を確認することができなかった．このエネルギー値の分布から，二次不安定性における臨界レイノルズ数（Rec(2)）を求める．
臨界値近傍の３点の値を線型近似し，Ek = 0 との交点を Rec(2)と定義した．
30
20
10
0
3200

3400

3600

Fig. 4 Energy of Fig.2 (e) and (f) as a growth of Re

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

Adjoint sensitivity 解析を用いた
パワーデバイス用半導体製造装置の最適設計のための数値解析
堀内 鷹之
大阪大学大学院 基礎工学研究科 物質創成専攻 化学工学領域

1 緒言F
シリコンカーバイド(SiC)は高い熱伝導率や広Γtar  u r  dΓ Γtar  u r  dΓAtar
 Γtar dΓ

(1)いバンドギャップを持つことから、次世代パワー

ここで、式(１)中の積分はtarget領域における面積

デバイス用半導体として期待されている[1,2]。過分を意味し、urは半径方向速度成分である。るつ飽 和 を 駆 動 力 と し た Top-Seeded Solution

ぼ温度微小変化 δTcru に対する目的関数の応答

Growth(TSSG)法によるSiC結晶成長は従来用い（これを局所感度と呼ぶ)、即ち、るつぼ壁の各座られている昇華法に比べて極めて高品質な結晶を標点上の δF/δTcru を計算するために、随伴変数法

作製できる一方で、低炭素溶解度のために成長速を用いて随伴方程式を導出した。以下に解くべき度が極めて遅い。炭素源のるつぼ(crucible)を高温式、及び局所感度を示す。これらの導出の詳細はにし、種結晶(seed)は低温とすることで過飽和度

参考文献[5]を参照されたい。

を高くできると考えられるが、このような場合は、
融液相に大きな温度差が存在する。種結晶付近で
自由表面上の温度差に起因するMarangoni対流が支配的となるため、種結晶面内成長速度が不均一化し結晶の高品質性が失われる恐れがある[3]。

支配方程式:

 u   u  1 p   2 u FE g  T  T0   u 0

 u   T   2T 

そこで、成長速度を低減させることなく u   C D 2C

Marangoni対流の影響のみを低減する制御法が随伴方程式:Q Cp

必要である。るつぼ温度分布条件を変えて制御す  u   u*   u*   u  p*  2 u*  T T *る手法が考えられるが、パラメトリックスタディ  u*  0

では任意の温度分布に対して検討するには非現実  u   T *   2T *  u*  g 的である。そこで本研究では、ノンパラメトリッ

局所感度:

ク感度解析[4]を用いた逆解析を用いて，るつぼ温

F度分布の最適条件を数値的に導出した。

  n T *  r 

 Tcru  r 

(2)

(3)

(4)

感度分布に基づき最急降下法[6]によって、るつ

2 解析手法

ぼ温度Tcruを最適値まで更新した。これらの解析にはオープンソースの有限体積法ソルバーである

解析領域は既報[3]と同一でありるつぼ内溶融液の概略図をFig. 1に示す。図中のtarget領域における種結晶へと向かうMarangoni対流を抑制す

OpenFOAM[7]をカスタマイズして用いた。

3 結果と考察ることを目標とし、次の目的関数Fの最小値問題
を考える。

Figure 2に最適化前後の融液内速度場及び温度
場分布を示す。最適化後には、種結晶に向かう
Marangoni対流による下降流が抑制され、上昇流

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

が見られるようになった。Figure 3に最適化前後の結晶成長速度分布を示す。最適化後は結晶成長速度を大きく損なうことなく、中心から4 mmまでの広範囲にわたって均一な結晶成長速度分布が
見られた。

4 結言
TSSG法によるSiC結晶成長の数値解析に随伴方を用いた逆解析を導入し、結晶成長速度を低下させることなくMarangoni対流を抑制することが可能なるつぼ温度分布をノンパラメトリックにFig. 1 Computational domain of SiC melt.

求めることができ、最適化後の結晶成長速度分布は広範囲にわたって均一化できた。

謝辞
本 研 究の 一部 は、 JSPS 科 研 費 基 盤研究A(18H03839)の助成を受けた。

引用文献
[1]

R. Wei et al., J. Appl. Phys. 113 (2013)053503.

[2]

F. Roccaforte et al., Microelectron. Eng.
187–188 (2018) 66–77.

[3]

T. Yamamoto et al., J. Cryst. Growth 470

Fig. 2 Flow velocity vectors and temperature
distribution of (left) the initial state and(right) the optimized state.

(2017) 75–88.
[4]

K. Momose et al., Heat Tran. Asian Res.
32 (2003) 1–12.

[5]

T. Horiuchi et al., J. Cryst. Growth 517(2019) 59–63.

[6]

D.G. Luenberger, Optimization by Vector
Space Methods, John Wiley & Sons.,
1969.

[7]

OpenFOAM, (see URL:
https://www.openfoam.com/)

Fig. 3 The growth rate profile for the
initial and optimized state.

14 − Academic Center for Computing and Media Studies, Kyoto Universityウマの個体間に作用する力の解明に向けた数値シミュレーション
井上漱太 1、深沢圭一郎 2、平田聡 1
1 : 京都大学野生動物研究センター
2 : 京都大学学術情報メディアセンター

1 はじめにえる。本研究では Attraction と Repulsion の相互作用を、距離を変数にもつ関数によってモデル化し、動物の群れは多様性に富んでいる。近年、計測技術の急激な発展とともに、群れを構成する個体の行

シミュレーションをおこない、野外観察によって得られたデータの再現を狙った。

動ルールに関する研究が高まりを見せている。個体間の単純なインタラクションが群れ全体の維持・協調にどのように寄与しているのだろうか。
Couzin ら

2 野外観察は Attration、Alignment、Repulsion という３種類

野外観察におけるデータ収集は 2016 年 6 月およの力の相互作用により動物の群れにおけるさまざび 2017 年 5-7 月にポルトガルのアルガ山でおこな

まな振る舞いが説明できると提唱した[1]。これによわれた。合計 3 群を終日追跡し、ドローン(Phantom

り、巨大な群れにおいても、近傍の個体との局所的

3、Mavic Pro)により画像データを収集した。それぞなインタラクションが、群れの全体が維持する機構

れの群れの個体数は 7-8 頭だった。画像は 30 分ごであることが理論的に提示された。また、魚類の群とに高度 25-80m 程度で撮影され、一枚の画像が群れにおいて、個体の加速度を定量的に計測すること

れの構成個体全てを含むように高度を調節した。ウ

により、この理論の有効性がすでに検証されている

マの群れは採食、移動、休息という 3 種類の状態に[2]。

大別することができ、本研究では採食場面のみを解

この分野の根底にある一つの大きな疑問は、全て析した。それぞれの群れに関して、合計 60 枚程度の群れに共通する個体間インタラクションルール

の画像を解析に使用した。個体の首の付け根と尾のの存在である。現状では、比較的多くの研究が魚類基部の中点に相当する点を個体の座標として、個体や鳥類を対象にしており、陸上性の哺乳類の群れに間距離の計測に使用した。

おける研究は少ない[3]。特に体サイズの大きな哺乳類の群れにおいて、個体の動きや群れ内での相対位置を定量的に捉えた研究は非常に少ない。そこで、

3 ウマの個体間に働く力の数値シミュレーション

本研究ではポルトガルのアルガ山に生息する野生のウマの群れを対象に、個体間距離を調節している力の数値シミュレーションのために、一個体の周力の推定を目指した。

囲に重力のように”場”を導入し、距離が遠いとウマは極めて社会性の高い動物と一般的に考えられている。野生環境下において、ウマの個体群はAttraction として、距離が近いと Repulsion として働く力を以下のように定めた。

ハレム群とバチェラー群に分かれる。ウマは有蹄類としては珍しく恒常的な群れ構成を維持する。群れのまとまりが安定して観察できることから個体間に働く力を推定するには、非常に適した動物だと言𝐹𝐹 = −𝑎𝑎(𝑟𝑟 − 𝑟𝑟( )− 𝑏𝑏(𝑟𝑟 − 𝑟𝑟( )*+ 𝑑𝑑𝑐𝑐京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 15この関数において r は個体間の距離をあらわし、ra は Attraction と Repulsion の境界値である。本研
究では、各群れの個体数と同数の点を 2 次元平面上にランダムに配置した状態から、F を個体間に作用させ、一定ステップ後に点間の距離を測定した。これを観察で得られた画像の枚数分繰り返した。そして、測定された個体間距離を確率分布とみなし、野外観察により得られた分布とピーク位置および残差平方和を比較し、a、b、c、d、ra の最適パラメー
タセットを探索した。

Figure 1. 個体間距離の確率分布。

0.4

0.5

4 相互作用モデル

0.1

ロジカル距離を使用した。メトリック距離において0.2

Fつかあるが、本研究ではメトリック距離およびトポ0.3ある個体の力の作用範囲を定めるモデルはいく

0.1

距離の近い順に Nt 個体とインタラクションを持つ。
F

タラクションを持つ。トポロジカル距離においては、

0.2

0.0

は、個体は周囲 Rm 内に存在する他個体全てとイン

5 結果

0

10

20
30
Distance (Body length)

40

50

Figure 2. トポロジカル距離を用いたシミュレーションにおける最適パラメータセットから再シミュレーションによる個体間距離の再現は 3 群
中 2 群で成功した (Figure 1)。メトリック距離とト
ポロジカル距離の最適パラメータはそれぞれ、Rm =9、Nt = 3 であった。二つのモデルの比較において、メトリック距離の方が一貫して小さな残差平方和を示したが、劇的な違いは生じなかった。
最適パラメータセットにおける F の挙動は魚群における加速度の振る舞いに近い性質を示した(Figure 2)。具体的には、ごく近い距離のみで、大きな Repulsion の力を示し、7 Body length 付近でAttraction の力に切り替わった。そして、Attractionの力は大きくはないものの、遠くまで作用した。

現された F の挙動。

6 参考文献[1] ID. Couzin, et al. "Collective memory andspatial sorting in animal groups." Journal oftheoretical biology 218.1 (2002): 1-11.
[2] JE. Herbert-Read et al. "Inferring the rules ofinteraction of shoaling fish." PNAS 108.46 (2011):
18726-18731.

[3] S. Inoue, et al. "Spatial positioning ofindividuals in a group of feral horses: a case study

using drone technology." Mammal Research 64.2
(2019): 249-259.

6 考察
本研究ではトポロジカル距離、メトリック距離、両モデルにおいて一個体は近傍のみとインタラク
ションした。インタラクションの作用範囲を広げると結果を全く再現できず、ウマの群れにおいても群れの中での局所的なインタラクションが群れを維持する機構であることが示唆された。

16 − Academic Center for Computing and Media Studies, Kyoto Universityք໘ͷຎࡲ઀৮Λߟྀͨ͠ଛইϞσϧʹΑΔ
మ‫ے‬ίϯΫϦʔτͷ 3 ࣍‫ݩ‬ഁյγϛϡϨʔγϣϯ૬അ ༔ਓ ∗∗

ҵ৓େֶେֶӃཧ޻ֶ‫ڀݚ‬ՊࣾձΠϯϑϥγεςϜՊֶઐ߈1 ͸͡Ίʹ

ࣜΛ࢖༻͢Δɽσ′ = (1 − D)Cε′ࣾձ‫ج‬൫ߏ଄෺ͷ҆શੑධՁͷߴ౓ԽΛਤΔͨΊɼwhen ε′1 > 0 (in opening)਺஋γϛϡϨʔγϣϯʹΑΓɼమ‫ے‬ίϯΫϦʔτͷഁյ‫ڍ‬ಈΛఆྔతʹ೺Ѳ͢Δ͜ͱ΁ͷॏཁੑ͕૿͠
͍ͯΔɽమ‫ے‬ίϯΫϦʔτͷഁյ‫ڍ‬ಈΛਫ਼౓Α͘࠶‫͢ݱ‬ΔͨΊʹ͸ɼమ‫ͱے‬ίϯΫϦʔτͷྗֶ‫ڍ‬ಈͷϞσϧԽʹՃ͑ͯɼమ‫ͱے‬ίϯΫϦʔτؒͷ෇ண‫ڍ‬ಈͷϞσϧԽ͕ղੳਫ਼౓Λࠨӈ͢Δॏཁͳ՝୊ͱͳ
Δɽమ‫ͱے‬ίϯΫϦʔτؒͷ෇ண͸ɼओʹࡐྉք໘ͷ೪ணྗɼࡐྉ͕ؒͣΕΔࡍͷຎࡲྗɼమ‫ے‬ͷઅʹΑΔࢧѹ఍߅ྗʹΑΓ੒Γཱ͍ͬͯΔɽ͔͠͠ͳ͕
Βɼ͜ΕΒͷཁҼͷӨ‫ڹ‬͸ɼమ‫ے‬ͷද໘‫ܗ‬ঢ়΍ɼίϯΫϦʔτͷԠྗঢ়ଶɼమ‫ے‬पΓʹੜ͡ΔͻͼׂΕͳͲͷ༷ʑͳ৚݅ʹΑͬͯมԽ͢ΔͨΊɼ෇ண‫ڍ‬ಈΛ༧ଌ͢Δ͜ͱ͸ࠔ೉ͳ໰୊Ͱ͋Γɼ͜Ε·Ͱʹ෇ண‫ڍ‬ಈΛৄࡉ͔ͭɼ༗ޮʹධՁͰ͖Δํ๏͸ཱ֬͞Ε͍ͯͳ͍ɽ
ͦ͜Ͱຊ‫Ͱڀݚ‬͸ɼ෇ணͷӨ‫ڹ‬ཁҼΛ௚઀൓өͤͨ͞మ‫ے‬ίϯΫϦʔτͷഁյγϛϡϨʔγϣϯख๏Λߏங͢Δɽͦͯ͠ɼ෇ணੑೳͷҟͳΔ‫ͱؙ߯‬ҟ‫ܗ‬మ‫ے‬ͷҾൈ͖ࢼ‫ݧ‬Λ໛ٖͨ͠ഁյγϛϡϨʔγϣϯΛ࣮ࢪ͢Δ͜ͱͰɼ෇ண‫ڍ‬ಈͷ࠶‫ݱ‬ੑΛ‫͢ূݕ‬Δɽ

(1)͜͜Ͱɼσ′ ͱ ε′ ͸ϑΥʔΫτද‫ʹه‬ΑΔ‫࠲ॴہ‬ඪ‫ܥ‬ͷԠྗϕΫτϧͱͻͣΈϕΫτϧɼC ͸஄ੑ܎਺ϚτϦοΫεɼD ͸ 0 ∼ 1 ͷ஋ΛऔΔଛইม਺Ͱ͋Γɼࡐྉ͕݈શͳঢ়ଶͰ͸ 0ɼ‫׬‬શʹഁյͨ͠ঢ়ଶͰ͸1 ͱͳΔɽ
ք໘ʹରͯ͠ਨ௚ํ޲ͷͻͣΈ ε′1 Λ༻͍ɼഁյ໘
ͷ઀৮Λ൑ఆ͢Δɽ઀৮͢Δ৔߹ʹ͸ɼք໘ʹରͯ͠ਨ௚ํ޲ͷॳ‫߶ظ‬ੑΛอ࣋ͤ͞Δ͜ͱͰ઀৮‫ڍ‬ಈΛද‫͢ݱ‬Δɽք໘ʹରͯ͠ਨ௚ํ޲ͷԠྗ σ′1 ͸࣍ࣜͰද͞ΕΔɽ

σ′1 = C11 ε′1 + (1 − D)C12 ε′2 + (1 − D)C13 ε′3when ε′1 ≤ 0(in contact)(2)઀৮ঢ়ଶʹ͓͚Δք໘ͷ‫ݻ‬ண΍͢΂ΓΛߟྀ͢ΔͨΊɼΫʔϩϯͷຎࡲଇΛಋೖ͢Δɽࣜ (2) ʹΑΓಘΒΕͨք໘ʹରͯ͠ਨ௚ํ޲ͷԠྗ σ′1 Λ༻͍Δ
͜ͱʹΑΓɼຎࡲԠྗ τf ͸࣍ࣜͰද͞ΕΔɽ

τf = µσn = µ|σ′1 |(3)͜͜Ͱɼµ ͸ຎࡲ܎਺ɼσn ͸ਨ௚ԠྗͰ͋Δɽք໘ͷ‫ݻ‬ண΍ຎࡲ઀৮Λ൐͏ഁյ‫ڍ‬ಈΛ࠶‫͢ݱ‬Δ

2 ഁյγϛϡϨʔγϣϯख๏ͨΊɼ‫࠲ॴہ‬ඪ‫ܥ‬ͷͤΜஅԠྗ τ′12 ͸࣍ࣜͰද͞ΕΔɽ෇ணͷཁҼͱͳΔք໘Ͱͷഁյ΍ຎࡲྗΛߟྀ͢ΔͨΊɼචऀΒ͕։ൃͨ͠ք໘ͷຎࡲ઀৮Λߟྀͨ͠ଛইϞσϧ [1] Λద༻͢Δɽք໘ʹ௚ަ͢Δ‫ॴہ‬࠲ඪ‫ܥ‬Λߟ͑Δ͜ͱͰɼք໘ͷྗֶ‫ڍ‬ಈΛఆࣜԽ͢Δɽք໘ʹରͯ͠ਨ௚ํ޲ͷͻͣΈ ε′1 Λ༻͍ɼഁյ
໘ͷ։‫ޱ‬Λ൑ఆ͢Δɽ։‫͢ޱ‬Δ৔߹ʹ͸࣍ࣜͷߏ੒′
τ′12 = C44 γ12

when τf ≥ |τ′12 |(not sliding)when τf < |τ′12 |(sliding)′′
τ′12 = (1 − D)C44 γ12)Dτf+ sgn(γ12(4)(5)͜͜Ͱɼsgn ͸ූ߸ؔ਺Ͱ͋Δɽ

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

100

Unit: mm

100
40

100

Displacement
0.5 mmz

Load (kN)

Constraint: xyRound bar
Deformed bar

20

x
0
0

0.4

0.2

Interface

Displacement (mm)

ਤ 2: ‫ؙ͓߯‬Αͼҟ‫ܗ‬మ‫ے‬Ϟσϧʹ͓͚Δՙॏ − มҐؔ܎
Round bar (  16)

Deformed bar (D16)

Cross section

ਤ 1: Ҿൈ͖ࢼ‫ݧ‬ͷ‫ࢼڙ‬ମΛ໛ٖͨ͠༗‫ݶ‬ཁૉϞσϧද 1: ࡐྉύϥϝʔλ
k

Gf

0

0.2

20

0.1 N/mm

0.00010.2

20

0.1 N/mm

0.004

E



Steel

200 GPa

0.3

Concrete

30 GPaInterface

200 MPa

Loading step: 60/100 (6.5 kN)

Loading step: 60/100 (45.0 kN)(a) Round bar(b) Deformed bar

0.5

1.0

0.5
Internal damage variableਤ 3: ‫ؙ͓߯‬Αͼҟ‫ܗ‬మ‫ے‬Ϟσϧʹ͓͚Δ಺෦ͷଛই෼෍

3 ෇ண‫ڍ‬ಈͷ࠶‫ݱ‬ੑ
3.1 ղੳ৚݅

͕ഁյ͍ͯ͠Δͷʹର͠ɼҟ‫ܗ‬మ‫Ͱے‬͸಺෦ͻͼׂΕ΍ɼమ‫ے‬ͷ࣠ʹԊͬͨॎͻͼׂΕΛ࠶‫͖ͨͰݱ‬ɽ෇ண‫ڍ‬ಈͷ࠶‫ݱ‬ੑΛ‫͢ূݕ‬ΔͨΊɼ෇ணੑೳͷҟ
ͳΔ‫ͱؙ߯‬ҟ‫ܗ‬మ‫ے‬ͷҾൈ͖ࢼ‫ݧ‬Λղੳͨ͠ɽ༗‫ݶ‬ཁૉϞσϧΛਤ 1 ʹࣔ͢ɽରশੑΛߟྀ͠ɼ1/4 ϞσϧΛ࢛໘ମཁૉʹΑΓ࡞੒ͨ͠ɽమ‫ے‬ͷද໘‫ܗ‬ঢ়ΛৄࡉʹϞσϧԽ͠ɼද໘ʹް͞ 0.5 mm ͷք໘૬Λઃ͚ͨɽ‫ؙ߯‬Ϟσϧͱҟ‫ܗ‬మ‫ے‬Ϟσϧͷཁૉ਺͸໿60 ສͰ͋ΔɽࡐྉύϥϝʔλΛද 1 ʹࣔ͢ɽమ‫ے‬ͷද໘‫ܗ‬ঢ়ͷΈΛม͑ɼ2 έʔεͷղੳΛ࣮ࢪͨ͠ɽ

3.2 ղੳ݁Ռ
ਤ 2, 3 ʹՙॏ − มҐؔ܎ͱ಺෦ͷଛই෼෍Λࣔ͢ɽઅʹΑΔࢧѹ఍߅ྗͷࠩʹΑΓɼ‫ͱؙ߯‬ҟ‫ܗ‬మ‫Ͱے‬͸࠷େՙॏʹ͕ࠩੜͨ͡ɽ·ͨɼ‫Ͱؙ߯‬͸ք໘ͷΈ18 − Academic Center for Computing and Media Studies, Kyoto University4͓ΘΓʹຊ‫Ͱڀݚ‬͸ɼ‫ͱؙ߯‬ҟ‫ܗ‬మ‫ے‬ͷҾൈ͖ࢼ‫ݧ‬Λ໛ٖͨ͠ղੳΛߦ͏͜ͱͰɼఏҊख๏͕మ‫ے‬ͷද໘‫ܗ‬ঢ়ͷϞσϧԽΛม͑Δ͚ͩͰɼͦΕͧΕͷ෇ண‫ڍ‬ಈͷ܏޲Λৄࡉʹ࠶‫͖Ͱݱ‬Δ͜ͱΛࣔͨ͠ɽࢀߟจ‫ݙ‬[1] ૬അ ༔ਓɼं୩ ຑॹɿຎࡲ઀৮Λߟྀͨ͠ଛইϞσϧʹΑΔ४੬ੑࡐྉͷഁյγϛϡϨʔγϣϯɼ౔໦ֶձ࿦จू A2(Ԡ༻ྗֶ)ɼVol.74, pp.I 233-I 241,2018.

Ｎ結合型糖鎖修飾によるタンパク質の機能制御の関連性
リントゥルオト正美
京都府立大学大学院生命環境科学研究科応用生命科学専攻

N 結合型糖鎖修飾と機能制御の関連性について分子動力学法（MD）法を用いて研究を進めた。
KLK8（ニューロプシン）の基質特異性に大きな影響を与える Loop 99 は KLK ファミリーに共通し N 結合型糖鎖修飾部位を含む。ヒト KLK8 (PDB ID:1npm)に 2 種の糖鎖を結合したモデルと結合なしのモデルを比較した結果、糖鎖修飾がない場合に比べ糖鎖修飾したモデルでは結合部
位周辺のループの揺らぎに変化が表れ、特に Loop 99 の揺らぎが変化した。糖鎖修飾なしでは自由エネルギー曲面上に 3 つの極小点が存在しているが、Loop 99 が主に大きく動いていることがわかった。また、Loop 99 の動きに伴い基質結合部位のポケットの体積が減少し、触媒残基間の距離が大きく変化することがわかった。一方、糖鎖修飾することにより、基質ポケットはさらに広がる動きを見せた。糖鎖修飾によって結合部位の構造が影響を受けることがわかった。

1 諸言に解き明かされてはいない。
KLK ファミリーに属する KLK8（neuropsin）

細胞外のタンパク質や膜たんぱく質の多くが糖は中枢神経系で発現しており、このタンパク質の鎖修飾を受けている。タンパク質の主な糖鎖修飾

機能障害により、記憶障害や精神疾患との関連性

は Ser や Thr 側鎖の酸素原子に糖が結合する Oが示唆されてきた。本研究に先立って KLK8 とそ結合型糖鎖、Asn 側鎖のアミド窒素原子に糖が結の基質である neuregulin-1 の複合体のモデルと合する N 結合型糖鎖が存在する。これらの糖鎖修

して KLK8-ペプチドモデル（糖鎖はなし）を用い、飾はタンパク質の構造変化、機能改善、安定性や分子動力学計算を行った。その結果、KLK8 の基

溶解性の向上、品質管理、分子認識や情報伝達の質特異性は基質のアミノ酸配列に左右され、基質指標として様々な機能を有している。糖鎖修飾はが Loop 99 と効率的に相互作用を行うことが重要

ドラッグデリバリーシステムへの応用、腫瘍マーであることがわかった 2)。N 結合型糖鎖は基本型

カーの開発、新規医薬品開発などで注目されておの糖鎖であっても非常に大きい。Loop 99 への糖り、研究が多くなされてきたが、計算化学を用い鎖修飾は直接結合している Loop 99 の揺らぎはもた例は少ない。

ちろんのことタンパク質全体の揺らぎに影響を与

多様な組織、血漿中で発現し様々な疾病との関えることが予想される。

連性が示唆されている kallikrein（KLK）ファミ

本研究ではKLK8のLoop 99における糖修飾のリーは異なる基質特異性を有するセリンプロテア

反応サイトとその周辺の構造に与える影響に着目

ーゼでこれまでに 15 のタンパク質が報告されてし、分子動力学（MD）計算を用いて研究を行っいる。KLK ファミリーの反応サイトに共通してた。

存在する Loop 99 は酵素によってその長さは異な

2 方法

るが糖鎖結合部位を有しセリンプロテアーゼの基質特異性に関係していると考えられている 1)。しヒト KLK8(PDB ID:1npm) の Loop 99 内かし、糖鎖はそのフレキシビリティーから固体結

Asn95 にアセチルグルコサミン(GlcNAc)、マンノ

晶構造は得られておらず、機能との関連性は完全ース(Man)、ガラクトース(Gal)から成る糖鎖

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

（(GlcNac)2Man(ManGlcNacGal)2）を結合したや BAS では主成分分析の結果、それぞれ 2 つのモデルを基本形として BAN モデル、その末端をエネルギー極小点が存在し、Loop 99 の動きによ

シ ア ル 酸 で 延 長 し た 糖 鎖

る違いであった。WT とは異なり、糖鎖が結合し（(GlcNac)2Man(ManGlcNacGalSia)2）BAS モているモデルでは触媒残基間が水素結合距離にまデル、糖鎖なしの WT モデルを用いた。これらので短くなることはなかった。糖鎖の結合によりモデル周辺 20 Å に水を、生理条件になるように

BAN、BAS では基質ポケットが大きくなり、WT

Na+と Cl-を配置した。

では最大 2605 Å3 出会ったのに対し BAN、BAS

MD 計算には Gromacs3)を用い、アミノ酸残基に は amber ff14SB 力 場

では 5335 Å3、4144 Å3 であった。

4) 、 糖 鎖 に はGLYCAM06-j4)、水分子には TIP3P5)を用いた。

4 まとめ

最急下法によるエネルギー最小化に続いてモデルに位置拘束を課した状態で 300 K まで昇温した。

WT において触媒三残基間の距離が水素結合距

次に 1 ns かけて徐々に拘束を解除しながら、NPT離まで近づくのは 80 ns 以降の短い間であり、エ計算をおこなった。2 ns の平衡化後、100 ns のネルギー曲面上で 2 状態間に平衡関係があるといproduction run を行った。

いきれない。また、糖鎖結合モデルの結果からは

3 結果

N 結合型糖鎖修飾は機能に影響を与えうるということはわかったが、触媒三残基は反応を開始でき
る位置関係にはない。現状では複雑なエネルギー

3 つのモデルともに 15 ns 前後で RMSD はほぼ曲面上の比較的似通った局所的な構造のサンプリ

一定値に到達していた。そこで 20~100 ns のングしかできないことから、さらにサンプリング

trajectory を用いて Cα原子の RMSF の比較を行

空間を広げた計算を行うことを目的として、拡張った。
それぞれのモデルにおける b-factor を Fig.1

アンサンブル法の一つであるレプリカ交換法を用に示す。WT モデルでは Loop 99 に揺らぎがみらいた計算（REMD）を次に行うこととした。300 K

れるが、BAN モデルではその揺らぎが抑えられから450 K までの間の70 個のレプリカを用意し、ていた。一方で BAS モデルでは Loop 99 を含む

REMD を実行中である。さらに計算と解析を続

活性中心周辺のループの揺らぎが大きくなった。

けることによってこれらの糖鎖の基質ポケットや
反応活性点の構造変化と反応制御の関連性について調べたい。

5 参考文献Fig. 1 B-factor coloring of WT and two glycosylated
models.
次に主成分分析を行い自由エネルギーマップを

1)

Skala, Wolfgang et al. J. Biol. Chem. 289,
34267, (2014)2)

リントゥルオト正美他、日本コンピュータ化作成したところ、WT では 3 つの極小点が存在し学会2017 秋季年会精選論文特集号 16 巻
（5ていた。これらの構造は Loop 99 が主に動いてい号）、160-162

ることがわかった。この Loop 99 の揺らぎに伴い

3)(2013).

触媒残基間の距離が大きく変化することがわかった。結晶構造解析結果では２つの触媒残基間の距4)形成距離(1.9~2.3 Å)まで短くなっていた。Loop
99 の動きが活性部位の構造を変えることで活性を制御している可能性が示唆される。
同様に BAN20 − Academic Center for Computing and Media Studies, Kyoto University

J. A. Maier, et al., J. Chem. Theory
Comput., 11, 3696, (2015).

離はプロトンの移動を伴う反応を開始するには非常に大きい(約 4 Å)が、極小点の一つでは水素結合

S. Pronk, et al., Bioinformatics, 29, 845,

5)K. N. Kirschner, et al., J. Comput. Chem.,
29, 622, (2008).

6)

W. L. Jorgensen, et al., J. Chem. Phys., 79,926, (1983).

౎ࢢߏ଄෺ͷ‫ز‬Կతಛ௃͕΋ͨΒ͢େ‫ؾ‬ཚྲྀͷۭؒεέʔϧ΁ͷӨ‫ڹ‬٢ాɹහ࠸ ∗∗

1‫౎ژ‬େֶେֶӃཧֶ‫ڀݚ‬Պ஍‫ٿ‬࿭੕ՊֶՊઐ߈͸͡Ίʹߴ͞ͷඪ४ภࠩͱฏ‫ ͞ߴۉ‬Have ͷൺ Vh ͕ Vh = 0.0౎ࢢԽͷਐߦʹ൐͍ɺ౎ࢢͷߴԹԽʢώʔτΞΠϥϯυ‫ݱ‬৅ʣ΍ਓ޻ഉग़෺ͷ૿େʹΑΔେ‫ؾ‬Ԛછͱ͍ͬͨ౎ࢢಛ༗ͷ‫ڥ؀‬໰୊͕ΑΓਂࠁԽ͍ͯ͠Δɻ͜ΕΒͷ‫ݱ‬৅ͷ༧ଌʹ͸ਫ਼៛ͳ‫ؾ‬৅Ϟσϧͷ࢖༻͕ෆՄܽͰ͋Δɻ͔͠͠ɺ౎ࢢߏ଄෺্ͷ‫ؾ‬৅༧ଌʹ(V00), Vh = 0.5 (V05), Vh = 1.0 (V10) ͱͳΔΑ͏ͳϒϩοΫྻΛ༻͍ͨɻϒϩοΫͷີ౓Λࣔ͢‫΃ݐ‬͍཰ λp ͸ 0.25 ͱ͍ͯ͠Δɻຊ‫Ͱڀݚ‬͸ྲྀೖ‫ڥ‬քʹผ‫ྖࢉܭ‬ҬͰ࡞੒ͨ͠ཚྲྀ‫ڥ‬ք૚ྲྀΛλΠϜεςο
ϓຖʹ༩͍͑ͯΔɻ

͸ɺෳࡶͳ‫ܗ‬ঢ়͓Αͼ഑ஔΛͨ͠ߏ଄෺ʹΑΔେ‫ؾ‬΁ͷӨ‫ڹ‬Λద੾ʹද‫͢ݱ‬Δ͜ͱ͕‫ٻ‬ΊΒΕΔɻಛʹཚྲྀʹΑΔӡಈྔ༌ૹաఔʹ͍ͭͯɺ౎ࢢͷޮՌΛߟྀͯ͠ύϥϝλϥΠζ͢Δ͜ͱ͕ॏཁͰ͋Δɻཚྲྀӡಈྔ༌ૹͷύϥϝλϦθʔγϣϯʹ͸ओʹཚྲྀͷۭؒεέʔϧ͕࢖༻͞ΕΔ͕ɺෳࡶͳ౎ࢢߏ଄෺
ͷӨ‫ڹ‬Λड͚ͨཚྲྀͷۭؒεέʔϧʹؔ͢Δৄࡉͳ

Fig 1: block layouts in the unit area of (a) V00, (b)V05, and (c) V10.

ධՁ͸े෼ߦΘΕ͍ͯͳ͍ɻ‫ط‬ԟ‫ʹڀݚ‬ΑΓɺ࣮ࡏ౎ࢢߏ଄෺ͷॏཁͳಛ௃Ͱ͋Δ‫ݐ‬෺ߴ͞ͷ͹Β͖͕ͭฏ‫ۉ‬෩଎΍ϨΠϊϧζԠྗͱ͍ͬͨ୅දతͳ౷‫ܭ‬ྔʹେ͖͘Ө‫͢ڹ‬Δ͜ͱ͕ࣔ͞Ε͖ͯͨɻͦ͜Ͱຊ‫Ͱڀݚ‬͸ɺ‫ݐ‬෺ߴ͞ͷ͹Β͖͕ͭ‫͢΅ٴ‬ཚ3݁Ռཚྲྀ૊৫ߏ଄Λநग़͢ΔͨΊɺ্ঢӡಈʹରԠ͢ྲྀͷۭؒεέʔϧ΁ͷӨ‫ڹ‬Λ໌Β͔ʹ͢Δ͜ͱΛ໨

Δ ejection ͱԼ߱ӡಈʹରԠ͢Δ sweep Λ൑ผ৚݅

తͱ͢ΔɻͦͷͨΊɺཚྲྀΛཅʹ‫ࢉܭ‬Մೳͳ Largeeddy simulation Λ༻͍ͯɺ౎ࢢߏ଄෺্ͷେ‫ؾ‬ཚྲྀͱͨ͠৚݅෇͖ฏ‫ۉ‬Λߦ͍ɺejection ͓Αͼ sweep पͷ਺஋γϛϡϨʔγϣϯΛ࣮ߦ͠ɺཚྲྀͷۭؒεέʔΑͼ sweep ͸ӡಈྔ༌ૹͷେ෦෼Λ୲͓ͬͯΓɺӡಈϧͷղੳΛߦͬͨɻ

༌ૹΛߦ͏ཚྲྀߏ଄ͷۭؒεέʔϧͷࢉग़ʹదͨ͠Γͷཚྲྀߏ଄ͷۭؒεέʔϧΛࢉग़ͨ͠ɻejection ͓৚݅ͱ͍͑ΔɻFig. 2 ʹओྲྀํ޲ͷ௕͞εέʔϧ Lx

2ͱεύϯํ޲ͷ௕͞εέʔϧ Ly Λࣔ͢ɻૈ౓ϒϩο‫ࢉܭ‬ઃఆ

Ϋߴ্ۭ͞ͷ Have − 5Have ͷߴ౓Ͱ͸ɺsweep Λத‫ྖࢉܭ‬Ҭ͸ओྲྀํ޲ 4 kmɾεύϯํ޲ 2.4 km Ͱɺਫฏղ૾౓͸ 2 m ͱͨ͠ɻ‫ྖࢉܭ‬Ҭ಺ʹ͸ Fig. 1 Ͱࣔ͞ΕͨΑ͏ͳૈ౓ϒϩοΫ͕‫܁‬Γฦ͠ฒ΂ΒΕͯ
͍Δɻ֤ϒϩοΫ͸ਫฏํ޲ʹ 10 ֨ࢠͰղ૾͞Εͯ
͍Δɻ‫ݐ‬෺ߴ͞ͷ͹Β͖ͭͷӨ‫ڹ‬Λௐ΂ΔͨΊɺ‫ݐ‬෺৺ͱͨ͠ཚྲྀ૊৫ߏ଄ͷۭؒεέʔϧͷํ͕ ejectionΛத৺ͱͨ͠௕͞εέʔϧΑΓ௕͘ɺ͜ͷ݁Ռ͸ Vhʹґଘ͍ͯ͠ͳ͍ɻҰํɺૈ౓ϒϩοΫߴ͞ҎԼͰ͸ɺV00 ͰಘΒΕͨ ejection ͓Αͼ sweep पΓͷ Lx
͸ɺV05 ΍ V10 ͷ Lx ΑΓ΋େ͖͘ɺཚྲྀߏ଄͕ओ
ྲྀํ޲ʹҾ͖৳͹͞Ε͍ͯΔ͜ͱ͕෼͔ͬͨɻ

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

Fig 2: Vertical profiles of Lx and Ly in (a, b) V00,(c, d) V05, and (e, f) V10. EJ and SW denote
ejection and sweep, respectively. The black-solid,-dashed and -dotted lines indicate the height ofblocks in V00, V05, and V10, respectively.
ଓ͍ͯɺओྲྀํ޲ͷ௕͞εέʔϧ Lx ͱԣํ޲ͷ௕͞εέʔϧ Ly ͷΞεϖΫτൺʹ͍ͭͯௐ΂ͨʢFig.

3)ɻͦͷ݁Ռɺૈ౓ϒϩοΫߴ্ۭ͞ʹ͓͚Δ ejection ͱ sweep पΓͷཚྲྀۭؒεέʔϧͷΞεϖΫτ
ൺ͸ɺVh ʹΑΒͣ΄΅ҰఆͰ͋ͬͨɻ͜ͷ͜ͱ͸ɺૈ౓ϒϩοΫ্ۭͷཚྲྀ૊৫ߏ଄͸஍ද໘ʹΑΒͣ૬ࣅతͳಛ௃Λ༗͢Δ͜ͱΛ͍ࣔࠦͯ͠ΔɻҰํɺ

Fig 3: Vertical profiles of Lx /Ly in (a) V00, (b)V05, and (c) V10. EJ and SW denote ejection andsweep, respectively. The black-solid, -dashed and-dotted lines indicate the height of blocks in V00,
V05, and V10, respectively.

ૈ౓ϒϩοΫߴ͞ҎԼͷߴ౓Ͱ͸ɺsweep Λத৺ͱྲྀํ޲ʹ৳ͼͨಛ௃Λ΋͍ͬͯΔ͜ͱ͕෼͔Δɻ͢

ߴ͞ͷ͹Β͖͕ͭҟͳΔߏ଄෺‫܈‬Λ࢖༻ͨ͠ Largeeddy simulation Λ࣮ߦͨ͠ɻͦͷͨΊʹɺejectionͳΘͪɺsweep पΓͷ૊৫ߏ଄ͷۭؒ‫ܗ‬ঢ়͸ૈ౓ϒͱ sweep ʹରԠͨ͠ཚྲྀ૊৫ߏ଄Λ৚݅෇͖ฏ‫ʹۉ‬ϩοΫʹΑͬͯม‫͠ܗ‬΍͍͢͜ͱΛ͍ࣔͯ͠Δɻ͜ΑΓநग़ͨ͠ɻͦͷ݁Ռɺૈ౓ϒϩοΫߴ౓্ۭͰ

Ε͸ɺsweep த৺ͷߏ଄͸Լ߱ӡಈΛ൐͏ͨΊૈ౓͸ɺཚྲྀߏ଄ͷۭؒεέʔϧ͸ߴ͞ͷ͹Β͖ͭʹΑϒϩοΫߴ౓ҎԼʹଘࡏ͢Δ͕࣌ؒ௕͘ɺߏ଄෺ʹΒͣ૬ࣅతͳಛ௃Λࣔͨ͠ɻҰํɺૈ౓ϒϩοΫߴ

িಥ͢Δճ਺͕ଟ͘ͳΔ͜ͱ͕‫ݪ‬ҼͰ͋Δͱߟ͑Β౓ҎԼͷߴ౓Ͱ͸ɺૈ౓ϒϩοΫͷӨ‫ڹ‬Λड͚ओྲྀΕΔɻํ޲ʹҾ͖৳͹͞Εͨಛ௃Λࣔͨ͠ɻͨ͠ߏ଄͸ ejection Λத৺ͱͨ͠ߏ଄ʹൺ΂ͯɺओ4·ͱΊຊ‫Ͱڀݚ‬͸‫ݐ‬෺ߴ͞͹Β͖͕͓ͭΑ΅͢ཚྲྀ૊৫

ߏ଄ͷۭؒεέʔϧ΁ͷ‫͢΅ٴ‬Ө‫ڹ‬Λௐ΂ΔͨΊɺ

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

Numerical simulation of InGaSb crystal growth under micro-gravity
onboard the International Space Station
Research on the effects of heating rate for dissolution process
Jin Xin
Osaka University

1 Introduction

axisymmetric

gridsystems.

The

GaSb(feed)/InSb/GaSb(seed) sandwich sample

InGaSb is promising semiconductor that can be

was stacked into a quartz ampoule and sealed

utilized for many thermal-photo-voltaic devices.

with boron Nitride and Carbon sheets. A new

Micro-gravity

phase diagram of GaSb-InSb is utilized in current

environment

onboardthe

International Space Station can suppress all thevolumeadversebetter

OpenFOAM. The simulation is performed under

understanding of the transport phenomenon andzero gravity to save computational cost andcrystal growth mechanism.

further details for numerical methods andeffects

However,andtheoffersus

opportunities

a

for

space

experiments are very rare and require a lot of

averaging

continuum

model

in

procedures can be found in our previous work
(Yamamoto et al., 2016).

significant preparation and money. Therefore, it is

In addition, for different heating rates werenecessary to examine and the determined theadopted for current research: Case A, 3.6K/h; Casemost favorable conditions through numericalB, 7.2 K/h; Case C, 18.0 K/h; Case D, 36.0 K/h.

simulations in advance to shed light on future
space experiments.
Previous simulations have been carried out to
optimize the experimental conditions; however,they predicted a much longer dissolution length
for the feed crystals. The undissolved feed crystals
after the experiments is necessary for the postprocessing and analysis, so the over dissolution
problem of feed crystals should be discussed andprevented for future research. Therefore, a seriesof numerical simulations with different heating
rates are performed.

Figure 1 Schematics of (a) the growth ampoule, (b)its applied temperature and (c) the grid system forcurrent numerical simulation.

2 Numerical methods
Figure 1 shows the schematics of the growth
ampoule, its applied temperature and the 2D

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

Figure 2. Time evolution of the center position of crystal/melt interfaces. Dashed lines represent feed/melt interfaces
and solid lines represents seed (or grown-crystal)/melt interface in the dissolution process in 200,000 seconds.

3 Results and discussion

In the crystal growth process, in all four cases,
we observed the same growth and dissolution
rates as shown in Figure 1 after the heating

Figure 2 shows the computed time evolution ofprocess stopped. At the 200,000th second of thethe crystal/meld interfaces during the dissolution

simulation, in all Cases of B, C and D we observe

process at four different heating rates. The dashed

a complete dissolution of the feed crystal, and thelines represent the positions of the feed/melt

growth process already stopped due to the higher

interface during feed dissolution process, while theheating rates. The dissolution duration was muchsolid lines represent the positions of seed/meld

longer in Case A than that in Cases B, C and D,

interface during the seed dissolution process anddue to the smaller heating rate used in Case A.

the grown-crystal/meld interface during the

The crystal growth was still progressing at the

crystal growth process. It must be mentioned that

200,000th second. Therefore, we recommend thethe crystal growth process begins earlier at theuse of a heating rate of 3.6 K/h or smaller in thehigher heating rate, as seen from Case A to Case

ISS experiments (Inatomi et al., 2015).

D in Figure 1. In addition, at the larger heating。

rate we predict more seed dissolution and less feeddissolution before the beginning of the crystal

4 Conclusion

growth process. We also note that the seed crystal

In summary, current research shows that the

did not dissolve in Case A, being different from theheating rate is an important factor in affecting thecases of B, C and D. This is because in Case A the

dissolution length of the feed and seed crystals.

heating rate was so small that, before even the

Results also suggest that a heating value of 3.6K/

seed started dissolving, the solutes species fromor smaller would be appropriate for the future

the feed crystal diffused through the melt andspace experiments。

reached the seed interface, and accumulated there(which prevents seed dissolution) until the crystal
growth began. However, in Cases B, C and D, thelocal seed interface already dissolved before the5 Reference[1]

solutes from the feed reached the seed interface.
Cases B, C and D had larger seed dissolution
lengths (1.3mm, 2.8mm and 3.3mm, respectively)
than that of Case A.

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

T. Yamamoto et al., Num. Heat Transfer, Part
B: Fundam. 70(5), 441 (2016)[2]

Y. Inatomi et al., NPJ Microgravity. 1(1).
15011 (2015)

3 次元流体変数の予測中井 拳吾
東京大学大学院 数理科学研究科

本稿は斉木 吉隆氏 ( 一橋大学経営管理研究科 )りも圧倒的に計算量が少なくて済むため, 代わりにニューラルネットワークの次元を大きくできる. そとの共同研究に基づくものである.

のためたとえダイナミクスが複雑なふるまいを示す

1 はじめにとしても, それが決定論的である場合にはこの種の機械学習は様々な分野で注目されている. 近年,
機械学習の一種であるリザーバーコンピューティン
グ [1, 2] が時系列データやリャプノフ指数などの予方法が有効である.

3 流体測において有効であることが報告されている [2, 4].

周期境界条件の下で 3 次元非圧縮 Navier–Stokes

リザーバーコンピューティングにおける学習は入力

方程式の直接数値計算によって得られた時系列デー

データから得られたリザーバーベクトルと出力する

タをリザーバーコンピューティングで用いる.

ベクトルを線形にフィッテイングさせることに重点

流体のマクロ変数としてエネルギー関数を考察すが置かれている. この学習はニューラルネットワー

る. 波数 k ∈ N におけるエネルギー E0 (k, t) を次

ク構造を学習しないためフィッティングにかかる計で定義する:

算コストを減らしている. そこで我々はリザーバー
コンピューティングの手法を用いて構成したモデル

E0 (k, t) :=により流体変数の挙動を予測した.

2 リザーバーコンピューティング
dϕ/dt = f (ϕ) で表される力学系の変数, u =h1 (ϕ) ∈ R

Mと s = h2 (ϕ) ∈ R

Pについて考える.

ただし, ある時刻まで (学習時間と呼ぶ) の入力変数 u, 出力変数 s の時系列は既知とする. 出力変数s の時系列を予測する時刻 (予想時間と呼ぶ) において, 入力変数 u の時系列データは未知とする. 学
習時間において u を要素分解することで得た高次元ベクトル r の時系列データに対して出力が s を近似できるように r と s の線形関係を決定する. この決定がリザーバーコンピューティングにおける学習に相当する. (詳細は [3] を参照.)

1
2

∫

3
∑
�
�
�F[v ] (κ, t)�2 dκ,

Dk ζ=1ζ

ただし, Dk := {κ ∈ Z3 |k − 0.5 ≤ |κ| < k + 0.5} とする. F[vζ ] (κ, t) は 速度 v のフーリエ変換を表す.
微細な振動を除くために短時間の時間平均をとったE(k, t) =∑0∗
l=99 E0 (k, t − l∆t )/100 を考察する.

これによりエネルギー関数の本質的な挙動を見ることができる. これ以降時間平均を取った E をエネ
ルギー関数と呼ぶことにする.

4 結果
以下のように時系列データ u(t), s(t) を設定し,リザーバーコンピューティングにより学習を行う.

u = u(t) = (Ẽ(1, t), Ẽ(2, t), · · · , Ẽ(9, t))T ,このように物理的な知見を用いずに時系列データ

s = u(t + ∆t) = (Ẽ(1, t + ∆t), · · · , Ẽ(9, t + ∆t))T ,の学習のみからモデルを構築する. また, ニューラ

ただし, パラメータは [3] の表 I (b) を用いる. また,

ルネットワークの構造自体を学習する機械学習よẼ はエネルギー関数 E を正規化したものである.

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

E(k, t) 各波数 k の時間平均によって得られるエネ
3

ルギースペクトル E(k) = ⟨E(k, t)⟩ について, 時刻

reservoir output
actual1000 < t − T < 2000 の予測した時系列データから

~
E(4)

2
1

再構成することに成功した 図 1 (最下). このこと

0から速度のデータなどの微視的な時系列データを使-1

わずに構成したリザーバーによる力学系は統計量もうまく再現できていることがわかる.

-2
-3
0

50

100

150

200

250

t-T
3

以上のことから, 構成したリザーバーの力学系はエネルギー関数の挙動を表現するカオス力学系としreservoir output
actualて見ることができることが示唆されている.

2

~
E(9)

1

5 まとめ

0
-1
-2

流体のマクロ変数として各波数のエネルギー保有

-3
0

50

100

150

200

250

t-T
10

関数の時系列データの予測ができることが確認でき

ε100

た. また, モデルから得られた時系列データの統計

10-1

10

量としてエネルギースペクトルの一致も確認した.
これらのことから構成したリザーバーモデルがエネ

-2

10-3
0

50

100

150

200

250

t-T
100

E(k)

10-2
-3

10

-4

10-5

ルギー関数の挙動をうまく再現できることが確認で
きた.

6 謝辞

10-1

10

量を表すエネルギー関数の時間発展モデルを機械学習により構成した. 構成したモデルは , エネルギー1

中井は JSPS 科研費 19J12482 の助成を受けている. また, 本研究でおこなった計算の一部は京都大reservoir output
actual
1

10

k

図1

Inference of time-series of macroscopic variables of a fluid flow.

Ẽ(1, ·) から Ẽ(9, ·) の 9 次元の学習によってエ

ネルギー関数の予測に成功した. 図 1 (上 2 つ) にẼ(4, t) と Ẽ(9, t) の時系列データを書き出した. モ
デルから得られた時系列データを赤色の実線, 比較のため Navier–Stokes 方程式の直接数値計算によって得られた時系列データを青色の点線で書いている. 流体のカオス性に由来する初期値鋭敏性により長時間時系列の予測が失敗することが予測される.
実際, 図 1 (下から 2 番目) に書き出したエネルギー関数のエラー ε(t) =

∑ N0

ˆ

k=1 |Ẽ(k, t) − Ẽ(k, t)|/N0(N0 = 9) の成長は t − T ≲ 100 において指数関数的な増加が確認できた. 一方で, エネルギー関数

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

学のスーパーコンピュータ共同研究制度 ( 若手・女
性奨励枠 ) に基づく. ここに感謝の意を表す.

参考文献[1] H. Jaeger, and H. Haas, Scince, 304, (2004), pp.
78-80.
[2] Z. Lu, J. Pathak, B. Hunt, M. Girvan, R. Brockett, and E. Ott, Chaos, 27, (2017), pp. 041102.
[3] K. Nakai, and Y. Saiki, Phys. Rev. E, 98, (2018),
pp. 023111.
[4] J. Pathak, Z. Lu, B. Hunt, M. Girvan, and E.
Ott, Chaos, 27, (2017), pp. 121102.

Numerical simulation of deepwater oil blowout
Turbulent jets and droplet size distribution
Daniel Cardoso Cordeiro
大阪大学大学院 基礎工学研究科

1 Introduction

tension force, FD,k is the drag force calculated by
the Schiller-Naumann model, Fvm,k is the

In 2010, the largest offshore blowout in

virtual mass force and the subscript k indicateshistory happened in the Gulf of Mexico, USA.

the

fluid

phase.

AnLES

The sub-sea injection of chemical dispersantsturbulence model was also used.

Smagorinsky

was used to treat deepwater oil spills. However,The domain was modeled based on thewith only few studies prior to its application,

experimental data [1], i. e., an 1.0 m heightthe overall effectiveness of this method is still

cuboid tank with a 0.3 m x 0.3 m width. Oilbeing questioned as appropriate measures offlows through a round inlet at the center of thethe oil droplets were not performed in situ. The

tank initially filled with water. The DSD waspresent study investigates the droplet size

calculated a posteriori by isolating the droplets

distribution (DSD) and the turbulence featuresin the isosurface  = 0.5, computing theirof the blowout oil treated with chemical

volumedispersants in the area near the oil release

diameter considering them as a sphere.

andcalculating

theirequivalentthrough computational simulations. A hybrid

Six cases were calculated following thevolume-of fluid/Euler-Euler model with a large

experimental Weber number (We = 11 ~ 1.8 x

eddy simulation (LES) turbulent model was106) and Reynolds number (Re = 269 ~ 1076)

used to investigate the turbulent breakup of thegiven by: Re = UjD/j and We = jUj2D /where

oil jet into the water and the effects of surface

Uj is the jet velocity, D is the inlet diameter, j is

tension, inlet diameter and velocity.

the dispersed phase kinematic viscosity, j isthe dispersed phase density and is the surface

2 Numerical Method

tension.

The governing equations solved are thecontinuity (Eq. 1) and Navier-Stokes (Eq. 2)

The volume mean diameter d50 calculated byequations.

 k uk . k  0tk  k uk t

3 Results and Discussions

 k  k uk . uk   k p  .k uk  k  k g  FD,k  FS,k  Fvm,k(1)

the empirical models [2] and the present LES
results is shown in Fig. 1. The black dotted line
represents d50 from the experimental data.

(2)

Prediction by the present LES model agrees
fairly well with the experimental results for all

where u is the velocity,  is the volume fraction,cases, on the other hand, the empirical methods

t is time,  is the density, p is the pressure, g iscannot estimate d50 and the differences are of

the gravity acceleration, FS,k is the surface

several orders of magnitude.

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

average Weber number, which we found to be afunction of the Brownell-Katz number (BK) orof the We, Re and Froude number FrUj2/(gD)as shown in Eq. 3:

We 

max  2

BK  2

We

(3)

ReFr

4 Conclusion
Figure 1: Comparison of the accuracy of
methods for estimation of the volume mean
diameterd50.

Theapplicationofahybridvolume-of-fluid/Euler-Euler model with a large
eddy simulation turbulent model was effective
in modeling the laboratory scale oil blowout in a

Fig. 2a shows a positive isosurface of thewater tank. With a relative mean error of 0.19

second invariant of velocity gradient Q which

for the oil droplet median diameter prediction,

identifies an association of Kelvin-Helmholtz

theinstabilities (initial ‘rings’) and hairpin vortices

empirical models.

model

surpassed

thestate-of-the-art

that ‘peel off’ the jet core into droplets (shown in

However, further investigation is necessary

the white isosurface  = 0.5). This process is

to explain the behavior of the droplet Weberknown as primary atomization. In Fig. 2b, thenumber in order to better understand theinitial orderly hairpins start to become

phenomenon.

increasingly chaotic, enhancing the secondary
breakup of droplets.
The behavior of the average local droplet
Weber number, i.e., the one calculated using theReferences[1] Kujawinski et al., Environ Sci Technol, 45,

droplet diameter and droplet velocity, was also

1298–1306, 2011.

studied. The curves collapse into a similar[2] Socolofsky et al., Proc. of 38th AMOP, 2015.

trend when scaled by the local maximum[3] Miyata et al., JIME., 51, 109-116, 2016.

Figure 2: A positive isosurface Q (colored by the vertical vorticity) and  = 0.5 (white) for Re = 1076 and We= 179. (a) Details of the Kelvin-Helmholtz instabilities near the inlet and the hairpin vortices. (b)Increasing disorder of the turbulent vortices.

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

分子動力学計算による膜貫通型ペプチドとリン脂質二重膜の相互作用ダイナミクス
最上 譲二東北大学大学院工学研究科

本研究では、ペプチドとリン脂質二重膜との相互作用を調べるためにマイクロ秒スケールの分子動力学(MD)シミュレーション及び自由エネルギー計算を行った。
ペプチドはアミノ酸の配列を変えるだけで簡単に分子の特性を制御でき、αヘリックスを形成
させれば疎水性ドメインがリン脂質二重膜に貫通する事で、膜貫通タンパク質のように強固な細胞膜アンカーを実現できる。この細胞膜アンカー型ペプチドを利用する事で、細胞毒性の低い細胞標識ツールへの応用が期待されている。しかし、αヘリックス型ペプチドは横向きに膜表面に吸着する状態なども存在しており、ペプチドが膜に貫通する挙動を明らかにする事が効率的な細胞膜アンカーには重要である。そこで、αヘリックス型ペプチドとリン脂質二重膜との相互作用について MD 計算を通して検討した。

1 緒言

そこで、αヘリックス型ペプチドとリン脂質二重膜の安定性を熱力学的に定量評価するために、より安定に膜と相互作用して脱離しにくいモデ
ル分子として、多くの膜貫通タンパク質の膜ドメ

分子動力学(MD)シミュレーションならびに自由エネルギー計算を行った。

インに共通しているαヘリックスをモチーフにした膜貫通型ペプチド（6K、WALPx など）が研究されている[1]（図 1）。ところが、膜貫通ペプチ

2 計算方法ドはαヘリックスの両末端に親水性残基が存在し、計算対象とする分子は、αヘリックス型ペプチ

安定な膜貫通型をとるようにデザインされている

ドとして比較的よく研究されている 6K ペプチド

ため、ジャイアントリポソーム（LUV: Large(KKAAALAAAAALAAWAALAAAKKKK-NH2)

Unilamellar Vesicle）に膜貫通型ペプチドを作用を用いた[1]。リン脂質二重膜のモデルとしてはしても自発的に膜貫通型に移行しない問題が指摘

POPC 膜を用いて、POPC 膜の外部、内部にペプされている[2]。

チド分子がそれぞれ準安定的に存在する間の全原子のトラジェクトリをサンプリングした。長時間の MD 計算であるがレアイベントを見るのではなく、一定時間内において平衡状態にあると統計的に判断された状態の自由エネルギー計算から、熱力学的安定性を評価した。自由エネルギー計算
は、高速かつ高精度なエネルギー表示法(ER 法)に基づいた手法を用いた。ER 法は Matubayasi
らによって開発された、相互作用エネルギーの分

図 1. リン脂質二重膜へ貫通した-ヘリックス型ペプチド

布関数を汎関数溶液理論に適用した計算手法であ京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 29

り多体の溶媒効果を正確に評価する場合に適している[3]。
MD シミュレーションは以下に示す典型的な条件下で行った。
6K ペプチド1 分子に対してPOPC
は 173 分子、水は 6047 分子を直方体ユニットセルに周期境界条件で配置した。NPT アンサンブル
で 20℃、1 気圧の条件下で 10 ns(サンプリングイ
ンターバル：200 fs)の計算を NAMD パッケージを用いて行った。水分子モデルは TIP3P を用い、
6K ペプチドおよび POPC 分子は Charmm22 お
よび 36 力場を使用した。静電相互作用は PME

図 3. 6K ペプチドの CD スペクトル法を適用し、実空間カットオフを 12Åに設定した。

3 結果と考察
6K ペプチドの溶媒和自由エネルギー計算の結果を図2 に示す。
水中では–704.4 kcal/mol、
POPC
リン脂質膜中では–845.3 kcal/mol と得られた。こ
れより、6K ペプチドのリン脂質膜中への挿入にともなう自由エネルギー変化は DG = –140.9kcal/mol と見積もられた。この値は同じαヘリックス構造をとる 23 残基からなるペプチド
WALP23 と POPC 膜との相互作用を調べた粗視4 結言
以上の結果から、次におこなうべき以下の課題を得た。
(1) 水中における 6K ペプチドの構造を考慮した MD シミュレーション(2) 6K ペプチドの膜中の位置依存性(3) 6K ペプチドの配列を置換し膜透過性をコ
ントロールしたペプチドの MD 計算(3)については、既にアルギニン置換 6K ペプチド化シミュレーションによる自由エネルギー計算のを合成して、HeLa 細胞に暴露した試験により膜値–35 kcal/mol よりも安定性が高く見積もられた透過性が変化する実験結果を得ている（データ省[4]。この原因として、6K ペプチドの水中におけ略）。また、ペプチドとリン脂質二重膜の結合熱のる構造を初期構造のαヘリックスのまま計算して測定よりエンタルピー的寄与も明らかにした（デいる事が考えられる。実際に 6K ペプチドを合成

ータ省略）。このように、膜透過性のコントロール

して得られた水中における 6K と膜中を模倣した
SDS 中の 6K の CD スペクトルを測定したところ、SDS 中においてαヘリックスを示すピークが 2倍程度強調されることが分かった(図 3)。すなわち、水中における 6K ペプチドはαヘリックス構造を維持しているものの、CD スペクトルとして 50%程度緩んだ構造になっている事を示している。

と熱力学的評価を、実験と計算の両面からアプロ
ーチすることで、効率的な膜貫通ペプチドのデザインができると期待される。

5 引用文献[1] Y.-C. Tang, C. M. Deber, Biopolymers 76,
110–118 (2004).
[2] A. Christopher, B. Burkhard, Thermodyn

amics, Kinetics of Dynamic Systems, IntechOpen, DOI:10.5772/1433 (2011).
[3] N. Matubayasi, M. Nakahara, J. Chem.

Phys. 113, 6070 (2000).
[4] W. F. Drew Bennett et al., J. Chem. Phys.
図 2. 6K ペプチドの膜への挿入

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

143, 243127 (2015).

プログラム高度化共同研究報告

分子シミュレーションによるヌクレオソーム構造変化の網羅的探索
生体高分子の並列分子動力学シミュレーションプログラム SCUBA のチューニング
石田 恒
量子科学技術研究開発機構 量子生命科学領域

本共同研究制度課題では、量子科学技術研究開発機構にて開発した生体高分子の並列分子動力学シミュレーションプログラムSCUBAを用いたヌクレオソーム構造変化の網羅的探索を迅速に実行するためのプログラム高度化を実施した。特に、ループ変形、MPI 環境変数、IO 環境変数の最適化を実施することで計算精度を向上することができた。
今後、チューニングされた SCUBA
を用い、ゲノム DNA を収納するヌクレオソーム集合体（染色体）の構造を特徴づけるヌクレオ
ソーム構造変化を特徴づけるヌクレオソーム内相互作用を特定し、染色体上の遺伝子発現メカニ
ズムを原子レベルで理解する。

1 生体高分子の並列分子動力学シミュ
レーションプログラム SCUBA

算するPPPM(Particle-Particle Particle-Mesh)計算法、系のエネルギー、温度、圧力を一定に保つ様々な時間積分アルゴリズム、レプリカ交換、自由エネ

タンパク質や核酸などの生体分子は、その立体構

ルギー計算、エネルギー最小化、基準振動解析、原造を変化させることで機能を発揮している。生体分子結合長を固定することにより時間ステップの増子の立体構造変化を調べるためには、生体分子の大を可能とするSHAKE, RATTLEアルゴリズム、

個々の原子の相互作用を計算しながら立体構造変系の原子分布異方性により引き起こされる並列化化を時々刻々追跡する「分子動力学シミュレーショ

効率の悪化を防ぐための動的ロードバランスなど、ン」は非常に有用な方法である。従来の分子動力学

最新のアルゴリズムを搭載している。

シミュレーションでは計算アルゴリズムや計算機性能の制限により、シミュレーションの対象は単体のタン

2 SCUBA のプログラム高度化

パク質や核酸であった。しかし、生体内ではタンパク質や核酸が単体で働くことはむしろ少数であり、複初めに、プロファイル計測と性能向上に向けた

数のタンパク質、核酸が互いに様々な配向をとりな

調査・検討を実施した。用いた系は小規模系

がら集合体（生体超分子）を形成して機能を発揮す（small、122,862 原子、メモリ使用量：12GB 程る。

度、ディスク使用量：50MB 程度）及び大規模系

SCUBA (Simulation Codes of hUge

（large、2,107110 原子、50GB 程度、ディスク

Biomolecule Assemble) は生体分子（主に核酸、蛋

使用量：500MB 程度）である。実行並列数は小白質）の機能発現メカニズム解析研究を目的として

規 模 系 :1n16p1t(small),

量子科学研究開発機構により開発された、分子動力1n16p4t(small), 1n32p1t(small),1n32p2t(small),学シミュレーションシステムである。用いる運動方

大 規 模 系 :2n64p1t(large),

程式は、古典力学のニュートン運動方程式である。

4n64p1t(large), 4n32p2t(large)とした。
（並列数

SCUBAは空間分割法による並列化により、扱える(XnYpZt)の表記の意味は、X：使用ノード数、Y：

原子数は100万原子以上であり、高い並列化効率（並ノードあたりプロセス数、Z：プロセスあたりス

列化効率99.9%）を有する。

レッド数、X*Y*Z：総コア数）

SCUBAは長距離相互作用を高速かつ高精度に計

1n16p2t(small),
2n32p2t(large),

調査の結果、１．ファンデルワールス力 ellj15

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

の演算がコストの大半を占める ２．プロセス数 2

ジナル版に対して 1.02 倍性能向上した。Cray で倍による性能向上は 1.9 倍、３．スレッド数 2 倍は tune402 が最速となり、オリジナル版に対してによる性能向上は 1.5～1.7 倍程度、４．自動ベク

1.49 倍性能向上した。全体として最速値は Intel

トル化による性能向上は 1.3～1.7 倍程度、５．コ

が Cray に対して、1.21 倍高速であった。

ードジェネレータによりベクトル長やアンロールの制御が可能、６．大規模系では通信(主に 1 対 1通信)比率が増加し、プロセス間のインバランス発生、７．Intel コンパイラでは IO 時間大、バッフ
ァリングを行うための環境変数の指定で性能向上の可能性があることがわかった。
以上の結果を検討、
6 種類
（１．
ベクトル長
（Intel及び Cray コンパイラ）、２．ループ変形＋コンパ
イルオプション（Intel 及び Cray コンパイラ）、
３．OpenMP スケジューリング（Intel 及び Cray
コンパイラ）、４．Hugepage（Intel コンパイラ）、
５．MPI 環境変数（Intel コンパイラ）
、６．IO
環境変数
（Intel コンパイラ））の最適化を実施し、

2.3 OpenMP スケジューリング
OpenMP スケジューリングでは、guided で1-2%性能向上するが、コンパイラや入力データによらずス レッド並列未適用の方が最速であった。

2.4 hugepage
hugepage では、入力データによらず hugepageの使用による性能向上は殆ど無かった。

2.5 MPI 環境変数
MPI 環境変数では、Eager 通信主体の環境変数

性能を評価した。

指定により、2n64p1t(large)でのみ 1.31 倍性能向

2.1 ベクトル長

上した。

ベクトル長においては、VECLEN=16 がほぼ全てのケースで最速であった。ただし、Cray では
4n64p1t(large)で VECLEN=32 で 16 に対して5%性能向上する場合があった。

2.2 ループ変形

2.6 IO 環境変数
IO 環 境 変 数 で は 、 1n161t(small) 、
1n32p1t(small)、2n641t(large)、4n64p1t(large)
でそれぞれ 1.11 倍、1.19 倍、1.15 倍、1.32 倍性能向上した。
（主要な IO は rank0 のみによるリス
タートファイルの read/write で総 IO 量は数ループ変形（(1)最内ループブロッキング(ブロ10MB(small)～数 100MB(large)である。Intel はックサイズ:VECLEN) (2)外側ループフルアンロ

Cray に比べてread/write システムコールのIO 長ール (3)最内ループ内 select 文削除 (4)最内ルーが小さく、呼出回数が read は 7 倍以上、write はプ融合 1(VECLEN の do ループ) (5)最内ループ内

100 倍以上多い。そのため、Intel においては、

if 文削除 (6)最内ループ融合 2(if 文削除に伴うファイル IO のバッファリングを有効にする環

VECLEN の do ループ) (7)最内ループ内除算変更

境(Newton-Raphson 法による除算無効)）とコンパ

FORT_BLOCKSIZE)を指定することで性能向上

イルオプションの組合せにおいては、したと考えられる。
）

2n64p1t(large)では、Intel で tune106(ループ変形(1)(2)(7),-qopenmp コンパイルオプション追加)
が最速となり、オリジナル版に対して 1.05 倍性能

変数(FORT_BUFFERED

、

3 まとめ

向 上 し た 。 Cray で は tune402( ル ー プ 変 形

以上の結果を参考にして、本課題の計算規模に(1)(2)(3)(4),-O3 コンパイル オプション)が最速とおけるループ変形、MPI 環境変数、IO 環境変数

なり、オリジナル版に対して 1.39 倍性能向上した。

の最適化を実施する。そして、分子シミュレーシ

全体として最速値は Intel が Cray に対して 1.07

ョンによるヌクレオソーム構造変化の網羅的探索

倍高速であった。1n16p1t(small)では、Intel でを推進する。

tune006(ループ変形(1)(2)(7))が最速となり、オリ

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

巨大津波遡上時の木造家屋の瓦礫生成過程シミュレーション
浅井 光輝**九州大学大学院 工学研究院 社会基盤部門

1 はじめにる ASI-Gauss 法[2]を採用した．1 部材を 2 つの線
形チモシェンコはり要素だけで分割し，数値積分点

2011 年東北地方太平洋沖地震で生じた大規模なを順応的にシフトすることで部材の弾塑性挙動を

津波により，沿岸域の木造家屋の多くは瓦礫となっ低い計算コストで高精度に表現する手法であり，また．津波と瓦礫が混在して遡上することで，人的・たUpdated Lagrangian 記述を用いた定式化により経済的損失が拡大した．瓦礫の総量は約 2200 万ト

崩壊に至るまでの大変形挙動を解析可能な有限変ンにも及び，災害復旧時にはその膨大な瓦礫の処理

形弾塑性解析方法である．2001 年 9 月 11 日に起こに時間を有し，復旧・復興に遅れが生じた．今後，ったニューヨーク世界貿易センタービル（WTC）の南海トラフなど，同規模の地震・津波被害が生じる

崩壊再現解析[3]や実大 3 次元震動破壊実験装置（E-

危険性が高いことが指摘されている．そのため，津Defense）を用いた実大モデル振動実験との比較[4]

波遡上を事前に把握し災害に備えるだけでなく，同により，本解析コードは性能検証済みである．塑性

時に発生する瓦礫の総量と拡散状況を事前評価し，化の判定は要素両端の断面力を降伏関数に代入す

災害後の復旧活動に備える必要がある．本研究の目ることにより行う．本研究では以下の降伏関数を用

的は，特に木造家屋の倒壊を反映した津波遡上解析

いた．

技術を開発し，遡上解析のさらなる高精度化を実施

すると伴に，津波遡上後の都市全体に発生する瓦礫2

𝑀𝑀𝑦𝑦
𝑀𝑀𝑥𝑥 2𝑁𝑁 2𝑀𝑀𝑧𝑧 2
𝑓𝑓𝑦𝑦 = () +() +( ) +() = 1
𝑀𝑀𝑥𝑥0
𝑀𝑀𝑦𝑦0
𝑁𝑁0
𝑀𝑀𝑧𝑧0

 

総量と拡散・分布状況の事前評価が可能なツールへ

ここで，Mx，My，N，Mz はそれぞれ x，y 軸回りのとの発展を目指すことである．

曲げモーメント，軸力，ねじりモーメントである．

以上の目的から，まずは地震・津波に対する都市

右下添え字の“0”は全断面塑性値であることを示す．

全域の木造家屋倒壊解析から開発を進めており，現部材の破断は要素の断面力を解放することで表現

在以下の 2 点の課題に取り組んでいる．

しており，それにより崩壊現象の解析を可能として(1) 都市全域のモデル化及び地震応答解析(2) 流体構造連成解析手法の開発
本研究グループでは，粒子法の一種である SPH 法

いる．

3 プログラム高度化を用いた高知市の津波遡上解析をすでに実施しており，粒子数 10 億以上の規模の流体解析が可能であることはすでに確認している[1]．一方で，都市全域規模の地震応答解析を実施するにあたって，構造解析手法の計算速度の高速化が必須となった．

3.1 オリジナルコードの問題点
ASI-Gauss 法の従来の使用方法は，1 つの構造物のみを対象として構造解析を行うものであった．そのため，モデルの規模は大きくても節点数・要素数
ともに数万オーダーに抑えられており，計算コストが顕著な問題とはならなかった．しかし本研究で対2 構造解析手法

象とする都市全域規模のモデルは，要素数が 4000

本研究では構造解析手法として有限要素法であ万，節点数が 3000 万程度と大規模である．オリジ
ナルコードは，並列計算に対応しておらず，また要京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 33

素数が数千万オーダーのモデル解析を想定していないため，計算効率の低いアルゴリズムが存在していた．都市全域の大規模地震応答解析をするにあたって，ASI-Gauss 法の並列化およびアルゴリズムの改良が必要となった．

3.2 改良点
ASI-Gauss 法を高速化するため，ソースコードの

図 1 モデル小（要素数：4.4 万，節点数：3.0 万）

スレッド並列化（以下，OpenMP と記述）を段階的に行った．まずは計算コストが高く非効率な場所を計測し，連立一次方程式ソルバーと接触探査アルゴ
リズムの両者から改良を行った．ソルバーについては剛性行列の格納法として CRS 法（Compressed
Row Storage）を採用し，スレッド並列効率の高い
共役勾配法ライブリーを移植した．また接触探査には，特に粒子系の解析における近傍粒子探索法で用いられているバックグラウンドセルを使ったセル検索方法を応用することで，接触探索すべき要素を限定することで大幅な高速化を実現した．

3.3 性能評価

図 2 モデル大（要素数：217 万，節点数：155 万）
その結果，図 1 と同様のモデル・解析条件において，接触探査アルゴリズム改良後のコードでは，解性能評価に使用したモデルは，
図 1 に示す要素数：

析時間が 13 時間 37 分から 2 時間 11 分になり，さ43,934，節点数：30,370 のモデルである．このモデらに 6 倍近く高速化することができた．
（ソルバー

ルを用いて，
1000 ステップの計算に要する時間を測と併せると約 15 倍の高速化を実現．
）

定した．並列化コードの解析環境は，京都大学スー

以上の高速化は，非線形解析までに留まっており，

パーコンピュータ System A，1 ノード，64 スレッ

ASI-Gauss 法一番の特徴である崩壊後の解析までドである．その結果，従来のコードでは 33 時間 32

実施すると，急激に解析速度の低下を招く． 図 1 と分要した計算が，ソルバーの改良後の計算時間は，同様のモデルでは，全体剛性行列構築に 7 秒要した

13 時間37 分になり，
約2.5 倍の高速化を達成した．

のに対し，図 2 に示す要素数：217 万， 節点数 155

特にソルバーのサブルーチン単体に着目すると，35

万オーダーのモデルでは約 5 時間要することが分か

倍高速化（並列化効率：57％）している．この結果った．接触計算では，ギャップ要素を追加し，そのから，ソルバーは大幅に高速化したものの，他の計要素を介して接触要素間の力を伝えている．そして算効率が低いサブルーチンが，解析時間の多くを占全体剛性行列を再び構築している．そのため，現状

めていることが 2.5 倍の高速化に留まっている要因のコードでは接触判定・接触解除判定されるたびに，であることが判明した．

ギャップ要素が増減し，約 5 時間かけて全体剛性行次に実施した高速化は接触探査アルゴリズムの列を一から作り直している．以上，崩壊後の全体剛改良である．従来のコードでは，一要素が，その他性行列のアセンブリングに時間を要していること

すべての要素に対し，接触探査のための計算を行っに起因することが判明した．現在は，崩壊後の解析ていた．無駄な計算を削減するために，解析領域をの高速化に取り組んでいる．

セル単位で分割し，接触探査の候補を自身のセルと近傍のセルのみに限定することで接触探査の高速化を実現した．

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

4 都市全域の地震応答解析
 構造解析手法である ASI-Gauss 法には，改善の余地が依然としてあるものの，都市全域の地震応答解析が実現可能な高速化を達成することができた．

4.1 都市全域のモデル化
都市全域のモデル化には，地理情報システム（GIS）
を用いた[5]．GIS より得られる建物の立体位置情報と数値標高データから，都市全域の建物モデルを構

図 3 高知市の都市モデル

築した（図 3）
．また同データから，津波遡上解析に必要な地表面モデルの作成も可能である．
本研究で対象とする都市は，高知県高知市とした．
前述の通り，本研究グループでは，高知市の津波遡上解析も併行して実施している（図 4）．ASI-Gauss法を用いた都市全域のモデル化及び地震応答解析
を達成後は，津波遡上解析と連成することで，地震・津波に対する都市全域の木造家屋倒壊解析へと発展させる予定である．

4.2 都市全域モデルの領域分割
高知市モデルの概要は次の通りである．領域：7

図 4 高知市の津波遡上解析

km×10 km，建物数：82,916，要素数：43,411,162，
節点数：30,295,639，データサイズ：9.5 GB．単一のパソコンでも解析できないことはないが，数日の間に解析を終了させるためには，複数の計算機を同時に使うことが望ましい．そこで，複数の計算機に計算負荷が一致するように領域分割を行った．その例として，高知市全域モデルを 20 区画に分割した
例を示す（図 5）．
このモデルの平均要素数は 2,170,558（絶対偏差：
0.3 %），平均節点数は 1,514,782（絶対偏差：1.1 %）であった．領域分割の結果，各領域の平均解析時間
は 28 時間 28 分（System A：20 ノード，68 スレッ
ド，1000 ステップ）に抑えることができ，都市全域

図 5 高知市モデルの領域分割の地震動解析が数日の内に完了できるようになった．
図 6 は 20 区画のうち，1 区画を可視化した図であり，カラーコンターはモデルの水平変位を表している．高い建物ほど大きく揺れていることがわかる．
図 7 は図 6 をさらに拡大した図である．建物一つ一つの形状が再現されており，はり要素にてモデル化されていることが確認できる．

図 6 1 区画の地震応答解析結果

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

都市全域モデルにおいては，建物ごとに建築年代を考慮した耐力は設定されておらず，また地震・津波で崩壊の危険性の低い堅牢な建物も同じはり要
素モデル化しているなど，まだ改良の余地があり，引き続きモデルを洗練させていく．

5 まとめ
今回の共同研究により，ソースコードの OpenMP

図 7 1 区画の地震応答解析結果（拡大）化，共益勾配法ライブラリー導入等，特にソルバーを中心に計算速度が高速化した．これにより，都市全域の大規模地震応答解析が実現可能となり，本研

参考文献

究を大きく前進するものとなった．今後，継続的に

計算効率の低いアルゴリズムを改良することで，更[1] 江口史門，浅井光輝，大谷英之，一色正晴：建物群をなる計算高速化を図るとともに，都市全域モデルの含む地表面詳細モデルを用いた粒子法による三次元

修正，流体構造連成解析手法を開発し，巨大津波遡津波遡上解析，土木学会論文集 A1（構造・地震工学），

上時の木造家屋の瓦礫生成過程シミュレーション

Vol. 72，No. 4（地震工学論文集第 35 巻），pp. I_367-

を実施する方針である．

I_377，2016.
[2]

D. Isobe: Progressive Collapse Analysis of

謝辞

Structures: Numerical Codes and Applications,

本研究は，筑波大学・磯部大吾郎教授、田中聖三助

Elsevier, eBook ISBN: 9780128130421, Paperback
ISBN: 9780128129753, 2017.

教との共同研究として実施させていただきました．
コードチューニングをしていただいたクレイ・ジャ[3]

磯部大吾郎，チョウミョウリン：ASI-Gauss 法によ

パン・インクの方々、またサポートしていただきまる世界貿易センタービルの飛行機衝突解析，日本建した京都大学情報メディアセンター関係者各位に

築学会構造系論文集，第 600 号，pp.83-88，2006．

感謝申し上げます．また JSPS 科研費 17H02061 の助[4]

韓元相，磯部大吾郎：ASI-Gauss 法を用いた地震応成を受けて実施しました．ここに感謝の意を表しま

答解析コードの性能検証，計算工学講演会論文集

す．

CD-ROM，第 16 巻，2011．
[5]

大谷英之，陳健，堀宗朗：地震応答解析モデルの堅牢な自動構築のための床形状判読手法の開発，土木学会論文集 A1（構造・地震工学），Vol. 70，No. 4（地震工学論文集第 35 巻），pp. I_1124-I_1131，2014.

36 − Academic Center for Computing and Media Studies, Kyoto Universityҟํੑ஄ੑ೾ಈ໰୊ʹର͢Δ
ԋࢉࢠੵ෼࣌ؒྖҬ‫ڥ‬քཁૉ๏ͷߴੑೳԽࡈ౻ ོହ ∗∗∗

1‫܈‬അେֶ େֶӃཧ޻ֶ෎ ‫ڥ؀‬૑ੜ෦໳ॹ‫ݴ‬͢ɽͳ͓ɼղੳͰ͸ɼ‫౎ژ‬େֶͷεʔύʔίϯϐϡʔλʔͰ OpenMPɼMPI Λ༻͍ͯ‫ࢉܭ‬Λޮ཰తʹߦຊใࠂͰ͸ɼҟํੑ஄ੑ೾ಈ໰୊ʹ͓͚Δԋࢉࢠੵ෼࣌ؒྖҬ‫ڥ‬քཁૉ๏ (CQBEM:Convolution Quadra-ͳ͍ͬͯΔɽ࠷‫ʹޙ‬ɼࠓ‫ޙ‬ͷల๬౳ʹ͍ͭͯ΋ड़΂
Δ͜ͱͱ͢Δɽture time-domain Boundary Element Method) ͷߴੑೳԽʹ͍ͭͯड़΂Δɽۙ೥ɼ୸ૉણҡ‫ڧ‬ԽϓϥενοΫʢCFRPɿCarbon Fiber Reinforced Plastic)΍ΦʔεςφΠτ‫ࡐܥ‬ྉ౳ͷ (Ի‫ )ڹ‬ҟํੑࡐྉ͕ࣗ

2

ղੳϞσϧ

ಈं΍ߤۭɾӉ஦ɼ౔໦෼໺౳ɼ޻ֶͷ༷ʑͳ෼໺Ͱ

ղੳྫͱͯ͠ɼਤ 1 ʹࣔ͢Α͏ͳɼແ‫ݶ‬ҟํੑ஄஫໨ΛूΊ͍ͯΔɽ͜ΕΒҟํੑࡐྉ͸ɼߴ‫౓ڧ‬౳ੑମ V தͷ͖྾‫ʹ܈‬ΑΔ७໘֎೾ͷଟॏࢄཚ໰୊

Λཧ༝ʹ࣮༻Խ͕ਐΜͰ͍Δɽҟํੑࡐྉͷҡ࣋؅ཧʹ͸ɼ௨ৗͷۚଐࡐྉ౳ͱಉ༷ʹɼ௒Ի೾ඇഁյ
‫( ࠪݕ‬UTɿUltrasonic non-destructive Testing) ͷద༻͕ࢼΈΒΕ͍ͯΔɽ͔͠͠ͳ͕ΒɼԻ‫ڹ‬ҟํੑ͕‫ݪ‬ҼͰɼ௒Ի೾͸ෳࡶʹ఻ൖ͢Δ͜ͱ͔Βɼ୳ইਫ਼౓ͷ௿Լ͕‫ݒ‬೦͞Ε͍ͯΔɽͦͷͨΊɼԻ‫ڹ‬ҟํੑΛߟྀͨ͠ɼ஄ੑ೾ಈ (௒Ի೾͸‫ݻ‬ମதͰ஄ੑ೾ಈͷੑ࣭Λࣔ͢) ࢄཚղੳख๏Λ։ൃ͢Δ͜ͱ͸ɼUT ͷਫ਼౓޲্ʹେ͖ͳ‫د‬༩Λ΋ͨΒ͢ͱߟ͑ΒΕΔɽ͜ͷΑ͏ͳཧ༝͔Βɼҟํੑࡐྉதͷ஄ੑ೾ಈղੳ [1] ΋ۙ೥ɼ͍͔ͭ͘ߦͳΘΕ͍ͯΔɽ͔͠͠ͳ͕Βɼ೾ಈղੳʹదͨ͠‫ڥ‬քཁૉ๏ (BEM:Boundary Element

Method) Λ༻͍ͨҟํੑ஄ੑମதͷେ‫ن‬໛஄ੑ೾ಈ
ղੳ͸ɼ΄ͱΜͲߦͳΘΕ͍ͯͳ͍ɽ
ɹͦ͜ͰɼຊใࠂͰ͸ɼCQBEM Λ༻͍ͯɼҟํੑ஄ੑମதͷ͖྾‫ʹ܈‬ର͢Δ஄ੑ೾ಈࢄཚղੳΛߦͳͬͨྫΛ঺հ͢ΔɽCQBEM ͸ɼैདྷͷ࣌ؒྖҬ‫ڥ‬քཁૉ๏ΑΓ΋਺஋҆ఆ [2] Ͱ͋Δ͜ͱͰ஌ΒΕ͓ͯΓɼۙ೥ɼ޻ֶͷ༷ʑͳ࣌ؒྖҬ໰୊ʹద༻͞Εͯ
͍Δɽҟํੑ஄ੑ೾ಈ໰୊ʹ͓͚Δ CQBEM ͷৄࡉઆ໌͸จ‫[ ݙ‬3] ʹৡΓɼҎԼͰ͸ղ͘΂͖໰୊ͷղੳϞσϧʹ͍ͭͯ؆୯ʹઆ໌ͨ͠‫ޙ‬ɼ਺஋ղੳྫΛࣔΛߟ͑Δɽ͖྾͸΋ͪΖΜɼϥϯμϜʹ഑ஔ͞Ε͍ͯͯ΋ߏΘͳ͍ɽ௕͞ 2a ͷ͖྾‫ʹ܈‬ରͯ͠Ԗ௚্
޲͖ʹೖࣹ೾Λૹ৴͢Δͱɼೖࣹ೾ͱ͖྾‫܈‬ͷ૬‫ޓ‬࡞༻ʹΑΓɼࢄཚ೾͕ൃੜ͢Δɽ͜ͷ࣌ɼ໘֎มҐu3 (x, t) ʹରͯ͠ɼҎԼͷ࣌ؒྖҬ‫ڥ‬քੵ෼ํఔ͕ࣜ੒Γཱͭɽu3 (x, t) = uin
3 (x, t) −∫S

T33 (x, y, t) ∗ [u3 (y, t)]dSy(1)ͨͩ͠ɼt ͸࣌ؒ, uin
3 (x, t) ͸ೖࣹ೾ɼT33 (x, y, t) ͸ೋ࣍‫ݩ‬७໘֎ҟํੑ஄ੑ೾ಈ໰୊ʹର͢Δೋॏ૚֩ɼ

[u3 (y, t)] ͸͖྾։‫ޱ‬มҐΛද͢ɽຊ‫Ͱڀݚ‬͸ɼࣜ (1)ͷ৞ΈࠐΈੵ෼ ∗ Λ Lubich ͷԋࢉࢠੵ෼๏ (CQMɿ

Convolution Quadrature Method)[4] ͰධՁ͢Δɽͳ͓ɼCQM Λ༻͍Δ͜ͱʹΑΓࣜ (1) ͷ࣌ؒྖҬੵ෼֩ T33 (x, y, t) Λ௚઀ධՁ͢Δඞཁ͸ͳ͘ͳΓɼϥϓϥεύϥϝʔλ s Λ༻͍ͨ T33 (x, y, s) ΛධՁ͢Δ͜
ͱͱͳΔɽͦͷͨΊɼैདྷͷ࣌ؒྖҬ‫ڥ‬քཁૉ๏͕
ۤखͱͨ͠೪஄ੑ೾ಈ໰୊౳΋ CQBEM Ͱղੳ͢Δ
͜ͱ͕ՄೳͱͳΔ [2].
ɹͯ͞ɼࣜ (1) Λ࣌ؒͱۭؒʹؔͯ͠ద੾ʹ཭ࢄԽ͢Ε͹ɼ࠷ऴతʹ֤࣌ؒεςοϓʹ͓͍ͯ T[u] = uinͷ‫ࣜܗ‬ͷ୅਺ํఔࣜʹ‫ؼ‬ண͞ΕΔɽͦͷ୅਺ํఔࣜ

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

a

crack

2a
a

・ ・ ・・

・・・・

S

・・・・

incident wave

Vscattered wave

ਤ 1: ղ͘΂͖໰୊ɽΛղ͘͜ͱͰɼ֤࣌ؒεςοϓʹ͓͚Δ͖྾։‫ޱ‬มҐΛ‫ٻ‬ΊΔ͜ͱ͕ՄೳͱͳΔɽͨͩ͠ɼT ͸ɼࣜ (1)ͷ‫ڥ‬քੵ෼ʹΑΓ‫·ٻ‬Δ܎਺ߦྻɼ[u] ͸͋Δ࣌ࠁʹਤ 2: ౳ํੑࡐྉʹର͢Δ‫܈‬଎౓‫ۂ‬ઢͷྫɽ͓͚Δ͖྾։‫ޱ‬มҐɼuin ͸ೖࣹ೾͔Β੒ΔϕΫτϧͰ͋ΔɽɹҰൠతʹɼCQBEM ΋‫ؚ‬Ί࣌ؒྖҬ‫ڥ‬քཁૉ๏Ͱ
͸ɼ‫ڥ‬քཁૉ਺Λ M ͱ͢Ε͹ɼT ͸ M × M ఔ౓ͷ܎਺ߦྻͱͳΔɽͦͷͨΊ M ͕େ͖͍৔߹͸ɼ‫ܭ‬ࢉ࣌ؒɾ‫ه‬Ա༰ྔ͕๲େͱͳΔɽ‫ڥ‬քཁૉ๏ͷΞϧΰϦζϜΛվྑͯ͠ߴ଎ԽΛ໨ࢦ͢৔߹͸ɼߴ଎ଟॏ‫ۃ‬๏ (FMM:Fast Multipole Method)[5] ౳ͷํ๏Λద༻͢Δ͜ͱ͕੝ΜʹߦͳΘΕ͖͕ͯͨɼҰൠͷҟํੑ஄ੑ೾ಈ໰୊ͷ৔߹ɼ‫ج‬ຊղͷଟॏ‫ۃ‬ల։Λߦͳ͏͜ͱ͕೉͍͠ [2]ɽͦ͜Ͱຊ‫Ͱڀݚ‬͸ɼ๯಄Ͱड़΂ͨΑ͏ʹɼ·ͣ‫౎ژ‬େֶͷεʔύʔίϯϐϡʔλʔΛ༻͍ͯ Open MP ͱ MPI ͷϋΠϒϦουฒྻԽΛࢪ͢͜ͱͰɼେ‫ن‬໛‫ࢉܭ‬Λޮ཰తʹ࣮ࢪ͢ΔํࡦΛऔͬͨɽ

ਤ 3: Ұํ޲ CFRP ʹର͢Δ‫܈‬଎౓‫ۂ‬ઢͷྫɽ

×1 ͭͷ͖྾͋ͨΓͷ‫ڥ‬քཁૉ਺ = 128 × 20) ͱͨ͠ɽ

਺஋ղੳྫ

3

ҎԼɼ਺஋ղੳྫΛࣔ͢ɽ·ͨɼ࣌ؒ૿෼ cT ∆t/a ͸ɼcT ∆t/a = 0.01ɼ૯࣌ؒεςοϓ਺͸ 256 ͱͨ͠ɽͨͩ͠ɼcT ͸७໘֎೾ͷೖࣹ೾ͷ೾଎Ͱ͋Δɽ͜ͷͱ͖ͷશະ஌਺͸ 655360(શ‫ڥ‬քཁૉ਺ × ૯࣌ؒεςοϓ਺ = 2560 × 256) Ͱ͋

3.1

ղੳ৚݅

Γɼҟํੑ஄ੑ೾ಈ໰୊ʹ͓͚Δ‫ج‬ຊղ͕ดͨ͡‫ܗ‬ࣜͰ༩͑ΒΕͳ͍͜ͱ΋‫ؚ‬Ίͯɼ͔ͳΓͷେ‫ن‬໛೾਺஋ղੳྫͱͯ͠ɼແ‫ݶ‬ҟํੑ஄ੑମ V ΛҰํ޲୸ಈ໰୊ͱͳΔɽҰํɼೖࣹ೾ uin
3 (x, t) ͸ x2 ࣠Λਖ਼ૉણҡ‫ڧ‬ԽϓϥενοΫ (ҎԼɼҰํ޲ CFRPɼਫฏͷํ޲ʹ఻ൖ͢Δฏ໘೾ͱ͠ɼ࣍ͷΑ͏ʹ༩͑ͨɽํ޲Λ୸ૉણҡํ޲) ͱԾఆͨ͠ɽਤ 1 தͷ͖྾‫܈‬͸ɼ128(ਫฏํ޲ × Ԗ௚ํ޲͖྾਺ = 8×16) ‫ݸ‬ͷਫฏ͖྾‫͠ͱ܈‬ɼਫฏํ޲ʷԖ௚ํ޲ = 23.0a×7.0a ͷ௕ํ‫ྖܗ‬Ҭʹ‫ن‬ଇਖ਼͘͠഑ஔͨ͠ɽ֤͖྾͸ 20 ‫ݸ‬ͷ‫ڥ‬քཁૉʹ෼ׂ͓ͯ͠Γɼશ‫ڥ‬քཁૉ਺͸ 2560(͖྾૯਺38 − Academic Center for Computing and Media Studies, Kyoto Universityu0(1 − cos 2πα)2()
 cT t − x + afor (0 ≤ α ≤ 1)
λ
cT
α=
 0
for otherwiseuin
3 (x, t) =

(2)

(3)

͜͜Ͱɼu0 ͸ೖࣹ೾ͷมҐৼ෯ɼλ ͸ೖࣹ೾ͷ೾௕Λද͢ɽͨͩ͠ɼͪΐ͏Ͳ cT t/a = 0.0 Ͱೖࣹ೾͕࠷

crackԼํͷ͖྾ʹ౰ͨΔΑ͏ʹ࣌ؒΛௐ੔͍ͯ͠Δɽ·ͨɼຊใࠂͰ͸ɼOpenMP, MPI ͷϋΠϒϦουฒྻ͕ਖ਼͘͠ߦͳΘΕ͍ͯΔޮՌΛ֬ೝ͢ΔͨΊɼղੳର৅Λҟํੑ஄ੑ೾ಈ໰୊Ͱ࠷΋؆୯ͳ७໘֎೾
ಈ໰୊ͱ͍ͯ͠ΔɽͦͷͨΊɼҰํ޲ CFRP ͷղੳʹඞཁͳ஄ੑఆ਺͸ɼC44 = 1.0ɼC55 = 2.02 ͷΈ

scattered wave(a)30stepͰ͋ΓɼͦΕͧΕ C44 Ͱແ࣍‫ݩ‬Խͨ͠ɽ

3.2

ղੳ݁Ռͷଥ౰ੑΛࣔͨ͢Ίͷ‫܈‬଎౓‫ۂ‬ઢҰൠతʹɼҟํੑ஄ੑ೾ಈ໰୊Ͱ͸ɼղੳʹઌཱͪɼ‫܈‬଎౓‫ۂ‬ઢΛඳ͍ͯ೾ಈ఻ൖํ޲ʹର͢Δ೾଎

(͜͜Ͱ͸‫܈‬଎౓ͷҙຯ) Λ‫ٻ‬Ί͓ͯ͘͜ͱ͕ղੳ݁
Ռͷଥ౰ੑΛࣔ͢ҙຯͰॏཁͰ͋ΔɽࢀߟͷͨΊɼղ(b)90stepੳʹ༻͍ͨ౳ํੑࡐྉͱɼҰํ޲ CFRP ʹର͢Δ‫܈‬

incident wave଎౓‫ۂ‬ઢΛɼͦΕͧΕਤ 2ɼਤ 3 ʹ͓ࣔͯ͘͠ɽຊ‫ݚ‬‫Ͱڀ‬ѻ͏ҟํੑ஄ੑ೾ಈ͸ਤ 3 ʹ͓͚Δ qS2 ೾ (७໘֎೾ɼ੨ઢͰද‫ه‬ɼਤ 2 ͷ౳ํੑࡐྉͷ৔߹͸ SH೾ͱ‫ݺ‬͹ΕΔ) Ͱ͋Δ͜ͱɼ‫܈‬଎౓ϕΫτϧͷ஋͸√ͦΕͧΕͷ c0 = C44 /ρ Ͱແ࣍‫ݩ‬Խ͞Ε͍ͯΔ͜ͱʹ஫ҙ͞Ε͍ͨɽਤ 2 ΑΓɼ౳ํੑࡐྉதͰ͸ɼP

೾ͱ S ೾ͷҟͳΔೋछྨͷ೾ಈ͕ଘࡏ͢Δ͜ͱ͕Θ͔Δɽͨͩ͠ S ೾͸ೋ࣍‫ݩ‬໰୊ͷ৔߹ɼSV ೾ͱ໘(c)150step

֎೾ಈͰ͋Δ SH ೾ʹ෼ྨͰ͖Δ΋ͷͷɼͦΕΒͷ೾ಈ଎౓ࣗମ͸ಉ͡Ͱ͋Δɽ͔͠͠ͳ͕Βɼਤ 3 ͷ

ਤ 4: ౳ํੑࡐྉதͷ͖྾‫ʹ܈‬ΑΔ஄ੑ೾ಈࢄཚղ

Ұํ޲ CFRP ͷ৔߹ɼҟͳΔ೾ಈ଎౓Λ࣋ͭࡾछྨ

ੳ݁Ռ (a)30 (b)90 (c)150 ࣌ؒεςοϓʹ͓͚Δ͖

ͷ೾ಈ͕ଘࡏ͢Δ͜ͱ͕Θ͔Δɽಛʹٖࣅॎ೾Ͱ͋

྾‫܈‬पลͷ໘֎೾ಈ৔ɽ

Δ qP ೾͸ਫฏํ޲ʹૣ͘఻ൖ͢Δ͜ͱ͕Θ͔Δɽ·ͨɼٖࣅԣ೾ͷҰͭͰ͋Δ qS1 ೾͸ 45 ౓ํ޲౳Ͱ೾

ͷ݁Ռ͸ɼͦΕͧΕ 30ɼ90ɼ150 ࣌ؒεςοϓʹ͓໘͕Ϋϩε͠ͳ͕Β఻ൖ͢Δ͜ͱ͕Θ͔Δɽຊใࠂ͚Δɼ͖྾‫܈‬पลͷ७໘֎೾ͷઈର஋ΛՄࢹԽͨ͠

Ͱղੳର৅ͱͨ͠ qS2 ೾ͷ‫܈‬଎౓‫ۂ‬ઢ͸ɼԣ௕ͷପ݁ՌΛ͍ࣔͯ͠Δɽ·ͨɼਤதͷਫฏͷനઢ͸࣮ࡍԁ‫ܗ‬ঢ়Λ͍ࣔͯ͠ΔͨΊɼਫฏํ޲ʹ΍΍଎͍଎౓ͷ͖྾Λ͍ࣔͯ͠Δ͜ͱʹ஫ҙ͞Ε͍ͨɽਤ 4(a) ΑͰ఻ൖ͢Δ͜ͱ͕Θ͔Δɽ͢ͳΘͪɼਤ 3 ΑΓɼԻΓɼ࠷Լํͷ͖྾ʹΑΔࢄཚ೾ͷൃੜΛ֬ೝ͢Δ͜‫ڹ‬ҟํੑͷੑ࣭Λ೺ѲͰ͖ɼղੳ݁ՌͰࣔ͢ࢄཚ೾ͱ͕Ͱ͖Δɽͦͷ‫ޙ‬ɼਤ 4(b),(c) ΑΓɼ֤͖྾ؒͰ

͸ɼਤ 2 ΍ਤ 3 Ͱࣔͨ͠‫܈‬଎౓‫ۂ‬ઢʹैͬͯ఻ൖ͠

ͷଟॏࢄཚ೾΋֬ೝͰ͖Δɽਤ 4(a) ౳ʹ͓͍ͯൃੜ

͍ͯͳ͚Ε͹ਖ਼͍݁͠Ռͱ͸‫͍ͳ͑ݴ‬ɽͨ͠ࢄཚ೾Λ஫ҙਂ͘‫ݟ‬Δͱɼࢄཚ೾ͷ೾໘͸ಉ৺ԁͷ‫ܗ‬ঢ়Λอͬͯ఻ൖ͍ͯ͠Δ͜ͱ͕Θ͔Δɽ͢ͳ

3.3

Θͪɼਤ 4 ͷ݁Ռ͸ɼਤ 2 ͷ౳ํੑࡐྉʹର͢Δ‫܈‬਺஋ղੳ݁Ռ

଎౓‫ۂ‬ઢʹ͕ͨͬͯ͠ɼ౳ํʹ఻ൖ͍ͯ͠Δͱ‫͑ݴ‬

ਤ 4ɼਤ 5 ʹɼͦΕͧΕ౳ํੑࡐྉɼҰํ޲ CFRPʹର͢Δ਺஋ղੳ݁ՌΛࣔ͢ɽ͜ΕΒਤதͷ (a),(b),(c)ΔɽɹҰํɼਤ 5 ͷҰํ޲ CFRP ʹର͢Δ݁Ռʹ஫໨͢

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 39ΊΔͨΊʹཁͨ࣌ؒ͠͸ 4 ࣌ؒ 34 ෼ɼ‫ٻ‬Ίͨ։‫ޱ‬มҐΛ༻͍ͯղੳྖҬ಺෦ͷҟํੑ஄ੑ೾ಈ৔Λ‫ٻ‬ΊcrackΔͨΊʹཁͨ࣌ؒ͠͸ 71 ࣌ؒఔ౓ (಺఺૯਺ 32768
఺) Ͱ͋Δɽ಺఺ͷ‫ࢉܭ‬͸ɼฒྻԽલ͸ 3 िؒఔ౓
Λཁ͢Δ (‫ࠐݟ‬Έ) Ͱ͋ͬͨ͜ͱΛߟ͑Ε͹ɼ‫࣌ࢉܭ‬ؒΛେ෯ʹ࡟‫͢ݮ‬Δ͜ͱ͕Ͱ͖ͨͱ‫͑ݴ‬Δɽ

scattered wave(a)30step4͓ΘΓʹຊߴ౓Խࣄ‫Ͱۀ‬͸ɼ७໘֎೾ಈΛର৅ͱͨ͠ҟํੑ஄ੑ೾ಈ໰୊ʹର͢Δԋࢉࢠੵ෼࣌ؒྖҬ‫ڥ‬քཁૉ๏ͷߴੑೳԽΛߦͳͬͨɽҟํੑ஄ੑ೾ಈ໰୊ͷ‫ڥ‬քཁૉ๏Ͱ͸ɼ‫ج‬ຊղ͕ෳࡶͳͨΊɼଟ͘ͷ‫ࢉܭ‬͕࣌ؒඞཁͰ͋Δɽࠓ‫ޙ‬͸ɼ͞Βʹ FMM Λద༻͢Δ౳ͷํࡦΛऔΓɼߴੑೳԽΛਤΔ༧ఆͰ͋Δɽ·ͨɼೋ࣍‫ݩ‬७໘಺ҟํੑ஄ੑ೾ಈ໰୊΁ͷ֦ு΋ߦ(b)90step
incident waveͳ͏༧ఆͰ͋Δɽँࣙຊ‫ڀݚ‬͸ɼ‫౎ژ‬େֶʹ͓͚Δฏ੒ 30 ೥౓ϓϩάϥϜߴ౓Խࢧԉࣄ‫͞ࢪ࣮Ͱۀ‬Εͨ΋ͷͰ͢ɽ·ͨɼฏ੒ 30 ೥౓ֶࡍେ‫ن‬໛৘ใ‫ج‬൫‫ڞ‬ಉར༻ɾ‫ڞ‬ಉ‫ڌڀݚ‬఺ެื‫ڞܕ‬ಉ‫( ڀݚ‬՝୊൪߸:jh180049, ҟํੑɾඇ‫ۉ‬࣭ࡐྉதΛ఻ൖ͢Δ஄ੑ೾ಈղੳख๏ͷ։ൃͱඇഁ

(c)150stepյ‫΁ࠪݕ‬ͷԠ༻ɼ‫ؔػڀݚ఺ڌ‬:‫౎ژ‬େֶ) ͷࢧԉʹΑΓߦͳΘΕ·ͨ͠ɽਤ 5: Ұํ޲ CFRP தͷ͖྾‫ʹ܈‬ΑΔҟํੑ஄ੑ೾
ಈࢄཚղੳ݁Ռ (a)30 (b)90 (c)150 ࣌ؒεςοϓʹ
͓͚Δ͖྾‫܈‬पลͷ७໘֎೾ಈ৔ɽ
Δͱɼਤ 4 ಉ༷ʹɼਤ 5(a) ΑΓɼ͖྾‫͔܈‬Βࢄཚ೾͕ൃੜ͠ɼਤ 5(b), (c) Ͱଟॏࢄཚ೾ͷൃੜΛ֬ೝ͢
Δ͜ͱ͕Ͱ͖Δɽ͔͠͠ͳ͕Βɼൃੜͨ͠ࢄཚ೾Λ஫ҙਂ͘‫ݟ‬Δͱɼ౳ํੑࡐྉͷ݁Ռʹൺ΂ͯɼࢄཚ೾͸ਫฏํ޲ʹૣ͘఻ൖ͍ͯ͠Δ༷ࢠΛ֬ೝͰ͖Δɽ͢ͳΘͪɼҰํ޲ CFRP ͷ৔߹ͷࢄཚ೾͸ɼਤ 3 ͷ

qS2 ೾ʹ͓͚Δ‫܈‬଎౓‫ۂ‬ઢʹ͕ͨͬͯ͠ɼପԁঢ়ͷ೾໘Λอͬͯ఻ൖ͍ͯ͠Δ༷ࢠΛ֬ೝͰ͖Δ͜ͱ͔
Βɼղੳ݁Ռ͸ଥ౰Ͱ͋Δͱ൑அग़དྷΔɽ
ɹͳ͓ɼղੳʹ͸‫౎ژ‬େֶͷ Cray XC40 Λ༻͍ͨɽฒྻԽ͸ 8 ϊʔυʹ MPI ϓϩηεΛ 1 ׂͭͣͭ౰ͯɼ

1 ϊʔυ 68 ίΞͷ OpenMP ฒྻ‫ࢉܭ‬Λࢪͨ͠ɽ͜
ͷ৔߹ͷղੳʹཁͨ࣌ؒ͠͸ɼશ͖྾։‫ޱ‬มҐΛ‫ٻ‬40 − Academic Center for Computing and Media Studies, Kyoto Universityࢀߟจ‫ݙ‬[1] ࡈ౻ɹོହ,ʮԻ‫ڹ‬ҟํੑࡐྉதͷ௒Ի೾఻ൖγϛϡϨʔγϣϯʯ, ඇഁյ‫ࠪݕ‬, vol.68, No.2, pp.78-83,(2018).
[2] ࡈ౻ɹོହ,ʮ೾ಈղੳͱ࣌ؒྖҬ‫ڥ‬քཁૉ๏ʯ, ‫ܭ‬ࢉ޻ֶ, vol.24, No.3, pp.13-16, (2019).
[3] A. Furukawa, T. Saitoh and S. Hirose, Convolution quadrature time-domain boundary element
method for 2-D and 3-D elastodynamic analyses
in general anisotropic elastic solids, Eng. Anal.
Bound. Elem., vol.39, pp.64-74, (2014).
[4] C. Lubich, Convolution quadrature and discretized operational calculus I , Numer. Math.,vol.52, pp.129-145, (1988).
[5] V. Rokhlin, Rapid solution of integral equationsof classical potential theory, J. Comput. Phys.,vol.60, pp.187-207, (1985).

๞࿨౔ͷେ‫ن‬໛ม‫ܗ‬ɾྲྀಈ‫ࢉܭ‬Λ໨తͱͨ͠
‫ݻ‬ӷࠞ߹ MPM ͷ։ൃࢁ‫ ޱ‬༟໼ ∗∗౦๺େֶࡂ֐Պֶࠃࡍ‫ॴڀݚ‬ʢ˟ 980-8572 ‫ٶ‬৓‫ݝ‬ઋ୆ࢢ੨༿۠ߥ‫ࣈר‬੨༿ 468-1ʣ

1 ͸͡Ίʹ

ρ̄w = θρwۙ೥ଟʑ‫ݟ‬ΒΕΔେ‫ن‬໛౔࠭ࡂ֐ͷඃ֐཈੍ͷͨΊʹ͸ɼෳࡶͳ౔ͷ‫ڍ‬ಈΛද‫͠ݱ‬ಘΔ਺஋γϛϡϨʔγϣϯख๏͕‫ٻ‬ΊΒΕΔɽຊ‫Ͱڀݚ‬͸ɼଟ޸࣭ମཧ
࿦ʹΑΔ‫ݻ‬ӷࠞ߹ମͷࢧ഑ํఔࣜΛ‫طͨ͠ͱࣜૅج‬্ࣜΛ༻͍ͯɼ‫ݻ‬૬ɼӷ૬ͷ࣭ྔอଘଇ͸ɼͦΕͧΕ࣍ࣜʹΑͬͯද͞ΕΔɽρs

D(1 − θ)Dρs
+ (1 − θ)+ (1 − θ)ρs ∇ · vs = 0
Dt
Dt

(3)

Dρw
Dθ
+θ
+ θρw ∇ · vw = 0
Dt
Dt

(4)

ԟͷ‫ݻ‬ӷࠞ߹ MPM[1, 2] Λൃలͤ͞ɼ߱ӍͳͲͷӨ‫ʹڹ‬ΑΔେ‫ن‬໛ͳ஍൫ࡂ֐‫ݱ‬৅ΛΑΓ҆ఆతɾޮ཰
తʹද‫͢ݱ‬Δख๏ΛఏҊ͢Δɽ۩ମతʹ͸ɼैདྷͷख๏͕ӷ૬ʹରͯ͠ऑѹॖੑΛԾఆ͍ͯͨ͜͠ͱʹର͠ɼඇѹॖੑΛԾఆͨ͠ fractional-step ๏ [3] ʹΑΔਫѹͷӄత‫ํࢉܭ‬๏Λಋೖ͢Δ͜ͱͰɼ҆ఆੑͱޮ཰ੑͷվળΛਤΔɽ·ͨɼB-spline ‫ج‬ఈؔ਺Λ࠾
༻͢Δ͜ͱͰɼMPM ͷ௕ॴͷͻͱͭͰ͋ΔྖҬ෼ׂฒྻԽʹΑΔେ‫ن‬໛‫ʹࢉܭ‬ରԠ͢Δͱಉ࣌ʹɼ଎౓ޯ഑͕࿈ଓͱͳΔͨΊɼ‫ࢉܭ‬ͷਫ਼౓ͱ҆ఆੑΛ޲্͕‫ظ‬଴͞ΕΔɽચງɾӡൖɾଯੵͱ͍ͬͨ 3 ࣍‫ݩ‬తʹෳࡶͳ౔ਫࠞ߹ମͷ‫ڍ‬ಈΛ‫ؚ‬Ή໛‫ʹݧ࣮ܕ‬ରͯ͠ఏҊख๏Λద༻͠ɼͦͷද‫ݱ‬ੑೳΛ֬ೝ͢Δɽ

(2)

ρw

͜͜ʹɼv s , vw ͸‫ݻ‬૬ɾӷ૬ͷ଎౓Λද͢ɽ౔͸౔ࠎ֨ͷମੵมԽʹରͯ͠౔ཻࢠͷମੵมԽ͕ඍখͰ͋
Δ͜ͱ͔Βɼ౔ཻࢠͷඇѹॖੑ Dρ s /Dt = 0 ΛԾఆ͢Δɽ·ͨɼਫʹ͍ͭͯ΋ඇѹॖੑ Dρw /Dt = 0 Λಉ༷ʹ༩͑Δɽͦͯ͠ɼؒܺ཰ θ ͷ࿈ଓੑͷԾఆʹΑΓɼࣜ (3), (4) ΑΓ‫ݻ‬ӷࠞ߹ମʹؔ͢Δ࣍ͷ࿈ଓࣜΛಘΔɽ

∇ · [(1 − θ)vs + θvw ] = 0(5)࣍ʹɼ‫ݻ‬૬ͱӷ૬ͷଟ޸࣭ମதͷ෦෼Ԡྗ σ̄ s , σ̄w ͸ɼ֤ʑ࣍ࣜͰ༩͑ΒΕΔɽ

2 ղੳख๏σ̄ s = σ′ − (1 − θ)pw I(6)

2.1 ࢧ഑ํఔࣜ

σ̄w = −θpw I(7)ຊ‫Ͱڀݚ‬͸ଟ޸࣭ମཧ࿦ʹ‫͖ͮج‬ɼ๞࿨౔Λ౔ࠎ͜͜ʹɼσ′ ͸౔ࠎ֨ͷࣔ͢༗ޮԠྗɼ pw ͸ӷମͷ֨ͱؒܺਫͱʹ෼཭͠ɼ‫ݻ‬૬ͱӷ૬͔ΒͳΔ 2 ૬ͷѹྗɼ I ͸ 2 ֊ͷ୯ҐςϯιϧͰ͋Δɽਫͷ೪ੑʹ࿈ଓମͷॏͶ߹Θͤͱͯ͠ߟ͑ΔɽҎԼͰ͸ɼఴ͍͑ͭͯ͸ɼ‫ࢹڊ‬తӡಈྔަ‫׵‬΍ؒܺਫѹͱൺֱͯ͠ࣈ s, w ͸‫ݻ‬૬ͱӷ૬ʹؔ͢Δม਺Λද͢ɽ‫ݻ‬ମͱӷӨ‫͕ڹ‬খ͍ͨ͞ΊɼҰൠతʹ‫ݻ‬ӷࠞ߹ϞσϧͰ͸ಁମͷཻࢠີ౓Λ ρ s , ρw ͱ͢Δͱɼ֤ʑͷଟ޸࣭ମத

ਫ܎਺ʹߟྀ͞Ε͍ͯΔ΋ͷͱΈͳ͢ɽͷ෦෼ີ౓ ρ̄ s , ρ̄w ͸ؒܺ཰ θ Λ༻͍ͯ࣍ࣜʹΑΓද͞ΕΔɽρ̄ s = (1 − θ)ρ s(1)·ͨɼ‫ݻ‬૬ɾӷ૬ؒͷӡಈྔަ‫ ׵‬p̂ ͸ Darcy-ForchheimeଇʹΑΓ༩͑Δɽp̂ = p̂E − pw ∇θ(8)京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 41p̂E =

θ2 µ1.75 ρw θ2 |vw − vs |(vw − vs ) (9)(vw − vs ) + √√k
150 k θ3/2
k=

D250 θ3
150(1 − θ)2

Solid Lagrange point
Water Lagrange point(10)͜͜ʹɼg ͸ॏྗՃ଎౓ɼD50 ͸ฏ‫ؔ͢ʹܘཻۉ‬ΔύϥϝʔλͰ͋Δɽ͜ΕΒͷԠྗͱӡಈྔަ‫׵‬Λ༻͍ͯɼ‫ݻ‬૬ͱӷ૬

Background Euler gridͷӡಈํఔࣜ͸֤ʑҎԼͰ༩͑ΒΕΔɽρ̄ s a s = ∇ · σ̄ s + ρ̄ s b s + p̂(11)ρ̄w aw = ∇ · σ̄w + ρ̄w bw − p̂(12)

ਤ 1: ‫ݻ‬ӷࠞ߹ MPM ʹΑΔ཭ࢄԽͷ໛ࣜਤ͜͜ʹɼ a ͸Ճ଎౓ɼ b ͸෺ମྗͰ͋Δɽ

2.1.1

v ʹؔ͢Δલਐࠩ෼Λద༻͢Δ͜ͱʹΑΓɼࣜ (11),(12) ͸࣍ࣜͰද͞ΕΔɽࡐྉߏ੒ଇ‫ݻ‬૬ͷม‫ޯܗ‬഑ςϯιϧ F ͸ҎԼͷΑ͏ʹ৐ࢉ෼
ղͰ͖ΔͱԾఆ͢Δɽ

F = Fe F p(13)͜͜ʹɼ Fe ͸஄ੑม‫ޯܗ‬഑ɼ F p ͸઼ੑม‫ޯܗ‬഑Ͱ͋Δɽ‫ݻ‬૬ͷ஄ੑԠ౴Λද‫͢ݱ‬ΔͨΊʹɼຊ‫Ͱڀݚ‬͸ Hencky௒஄ੑϞσϧΛ࠾༻͢Δɽ઼ੑม‫ܗ‬ͷද‫ʹݱ‬͸ɼ࣍ʹࣔ͢ Drucker - Prager ͷ߱෬‫ج‬४Λ࠾༻͢Δɽ
√Φ(σ′ , c) = J2 (s(σ′ )) + ηp − ξc(14)͜͜ʹɼs ͸ภࠩԠྗɼ p ͸‫ݻ‬૬ͷ੩ਫѹԠྗɼc ͸

vn+1 − vns= ∇ · σ′n − (1 − θn )∇pn+1(1 − θn )ρs sw∆tn+ (1 − θ )ρs bs + p̂nE

vn+1 − vnwnn= −θn ∇pn+1
θn ρw w
w + θ ρw bw − p̂E
∆t

(16)
(17)

͜͜Ͱɼ‫ݻ‬ӷ֤૬ʹؔͯ͠ҎԼʹࣔ͢தؒ଎౓ v∗s , v∗w
Λಋೖ͢Δɽ

(1 − θn )ρs v∗s = (1 − θn )ρs vns[]+ ∆t ∇ · σ′n + (1 − θn )ρs bs + p̂nE

(18)

= (1 − θn )ρs v∗s − (1 − θn )∆t∇pn+1(19)(1 − θn )ρs vn+1sw[]
θn ρw v∗w = θn ρw vnw + ∆t θn ρw bw − p̂nE೪ணྗɼη, ξ ͸಺෦ຎࡲ֯ ϕ ΑΓఆ·ΔࡐྉύϥϝʔλͰ͋ΔɽDrucker - Prager Ϟσϧʹ͍ͭͯ͸ɼա৒

n
∗
n
n+1
θn ρw vn+1
w = θ ρw vw − θ ∆t∇pw

ͳμΠϨΠλϯγʔΛ཈੍͢ΔͨΊɼҰൠతʹඇؔ

(20)
(21)

࿈ྲྀΕଌ͕༻͍ΒΕΔɽ·ͨɼྲྀಈԽʹ൐͏๲ுม

ࣜ (19), (21) Λ࿈ଓࣜ (5) ʹ୅ೖ͢Δ͜ͱͰɼ࣍ͷѹ‫ʹܗ‬ΑΔ౔ࠎ֨ͷ่յΛද‫͢ݱ‬ΔͨΊʹɼ‫׬‬શ઼ੑྗϙΞιϯํఔ͕ࣜಋ͔ΕΔɽ
)][(1 − θn θn
n+1
n ∗
n ∗
+
∆t∇pw
(22)
∇ · (1 − θ )vs + θ vw +ρsρwΛԾఆ͠ɼ࣍ʹࣔ͢ଛইϞσϧΛಋೖ͢Δɽ

ξ




K if η c − p > 0,
K=



0 otherwise.

(15)

͜͜ʹɼK ͸ମੵ஄ੑ཰Ͱ͋Δɽ

2.2 ‫ݻ‬ӷࠞ߹ MPM
2.2.1

࣌ؒํ޲ͷ཭ࢄԽ

ҎԼͷఆࣜԽʹ͓͍ͯ͸֤ม਺ʹ෇ͨ͠ఴ͑ࣈ

n, n + 1 ͸࣌ࠁ tn , tn+1 ʹ͓͚Δ஋Λҙຯ͢Δɽ଎౓

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

2.2.2ۭؒํ޲ͷ཭ࢄԽ

‫ݻ‬ӷࠞ߹ MPM ͷख๏͸͍͔ͭ͘ఏҊ͞Ε͍ͯΔ
͕ɼຊ‫Ͱڀݚ‬͸ਤ-1 ʹࣔ͢Α͏ʹɼ‫ݻ‬૬ͱӷ૬ΛҟͳΔϥάϥϯδϡཻࢠͰද‫͠ݱ‬ɼ‫ࢠ֨ࢉܭ‬͸ಉҰͱ͢Δख๏Λ༻͍Δɽ͜ΕʹΑΓɼશ෺࣭ྖҬ͸ཻࢠ‫Ͱ܈‬ද‫͞ݱ‬Εɼશ࣭ྔ͸ཻ֤ࢠʹू໿͞ΕΔ͜ͱ͔ρ̄α dv(α = s, w)

(23)

z

0.3 m

Ω0αp

x

Wave maker

Upstream

0.6 mաఔͰҰఆͰ͋Δɽ∫mαp =Downstream
0.3 mΒɼα ૬ͷཻࢠ࣭ྔ mαp ͸ҎԼͰ༩͑ΒΕɼม‫ܗ‬ͷ͜͜ʹɼΩ0p ͸ཻࢠͷॳ‫ظ‬ମੵྖҬͰ͋Δɽ·ͨɼ଎Center point

Initial position ofwater maker౓΍ѹྗɼؒܺ཰ͳͲͷཻࢠ஋͸֨ࢠ఺஋Λ༻͍ͯ

Mobile speed: 0.6 m/s࣍ࣜʹΑΓิؒ͞ΕΔɽ
nn∑

NI (xαp )ΦI(α = s, w)

(24)

I=1

x

y

0.35 m͜͜ʹɼnn ͸શઅ఺਺ɼNI (xαp ) ͸֨ࢠ఺ I ͷ α ૬

Wave maker

Φ(xαp ) =ཻࢠ p ͷҐஔʹ͓͚Δ‫ܗ‬ঢ়ؔ਺ɼΦ ͸ॾʑͷ෺ཧྔΛද͢ɽ‫ܗ‬ঢ়ؔ਺͸༗‫ݶ‬ཁૉ๏ͳͲʹ࢖༻͞ΕΔ΋

0.5 m

ͷͱಉ༷ͷ‫͋Ͱࣜܗ‬Δɽ

0.4 m
0.64 mਤ 2: ࠭ࢁ΁ͷ‫ཱݽ‬೾িಥ࣮‫ݧ‬ͷ૷ஔ໛ࣜਤ

2.2.3 B-spline ‫ج‬ఈؔ਺ຊ‫Ͱڀݚ‬͸ MPM ಛ༗ͷ਺஋ৼಈΛ཈੍͢ΔͨΊʹ

B-spline ‫ج‬ఈؔ਺Λద༻͢ΔɽB-spline ‫ج‬ఈؔ਺͸ύϞʔϧυΛਨ௚ʹҾ্͖͛Δ͜ͱʹΑͬͯ࡞੒ͯ͠

ϥϝʔλۭؒͷ࠲ඪͷฒͼͰ͋ΔϊοτϕΫτϧʹΑ

͓Γɼ਺஋‫ࢉܭ‬Ϟσϧͷ࠭ࢁ͸࣮ࡍͷ࠭ࢁͷ‫ܗ‬ঢ়Λͬͯఆٛ͞ΕɼҰ࣍‫Ͱݩ‬͸ Ξ = {ξ1 , ξ2 , ..., ξn+p , ξn+p+1 }Ϩʔβ‫ʹܭ཭ڑ‬Αͬͯଌఆͨ͠σʔλΛ‫ʹج‬Ϟσϧͱॻ͔ΕΔɽ͜͜ʹɼξi ͸ i ൪໨ͷϊοτΛࣔ͠ɼn‫ܗ‬ঢ়Λܾఆ͢Δɽ·ͨɼ༧Ί 10 ඵؒͷࣗॏੵՙΛߦ͸‫ج‬ఈؔ਺ͷ਺ɼp ͸ଟ߲ࣜͷ࣍਺Λද͢ɽB-spline͏͜ͱʹΑͬͯॳ‫ظ‬৚݅Λ༩͑Δɽ‫ཱݽ‬೾ʹ͍ͭͯ‫ج‬ఈؔ਺͸࣍ͷΑ͏ʹ࠶‫ؼ‬తʹఆٛ͞ΕΔɽ
1 if ξi ≤ ξ < ξi+1
Ni,0 (ξ) = 0 otherwise͸଄೾૷ஔΛ߶ମཻࢠʹΑͬͯද‫͠ݱ‬ɼ࣮‫ͱݧ‬ಉ༷ͷҠಈ଎౓ɾ‫ʹ཭ڑ‬Αͬͯॳ‫ظ‬৚݅Λ࡞੒͢Δɽਫ࿏

(25)

ͷఈ໘͸ϊϯεϦοϓ‫ڥ‬քɼԼྲྀଆ͸ྲྀग़‫ڥ‬քɼͦͷଞͷ໘͸εϦοϓ‫ڥ‬քͱ͠ɼղੳͷ֤छύϥϝʔλ͸ද 1 ʹࣔ͢௨ΓͰ͋Δɽξi+p+1 − ξ
ξ − ξiNi,p−1 (ξ) +Ni+1,p−1 (ξ)
Ni,p (ξ) =ξi+p − ξiξi+p+1 − ξi+1(26)
ຊ‫Ͱڀݚ‬͸࠷ॳͱ࠷‫ޙ‬ͷϊοτ͕ p + 1 ճ‫ݱ‬ΕΔΦʔΓɼ࠭ࢁͷҰ෦͕ྲྀਫʹΑͬͯ࡟ΓऔΒΕɼ࠷ऴతʹϓϯϊοτϕΫτϧΛ༻͍Δ͜ͱͰɼύϥϝʔλۭؒྲྀಈதͷ‫ڍ‬ಈʹண໨͢Δͱɼ೾ͷিಥ௚‫ޙ‬͸͕࠭Լͷ୺෦ʹ͓͍ͯؔ਺஋ͱҰகͤ͞Δ͜ͱʹΑΓɼҰྲྀଆʹେ͖͘޿͕ΔΑ͏ʹ෼෍͢Δ͕ɼঃʑʹྲྀઢൠతͳ MPM ͱಉ༷ʹ֨ࢠ఺ʹ͓͍ͯσΟϦΫϨ‫ڥ‬‫ʹܕ‬มԽ͍ͯ͠Δ͜ͱ͕෼͔Γɼ࣮‫ͱݧ‬ಉ༷ͷಛ௃ք৚݅Λ༩͑Δɽతͳ‫ڍ‬ಈ͕ଊ͑ΒΕ͍ͯΔ͜ͱ͕‫ͯݟ‬औΕΔɽ͜Εղੳ݁Ռͱ࣮‫݁ݧ‬ՌͱͷൺֱΛਤ 3 ʹࣔ͢ɽਤΑԼྲྀଆʹඌΛҾ͘Α͏ʹ෼෍͍ͯ͠Δ͜ͱ͕෼͔ΔɽΑΓɼఏҊ͢Δख๏͸ैདྷͷղੳख๏Ͱ͸ද‫ࠔ͕ݱ‬೉Ͱ͋ͬͨ౔ਫࠞ߹ମͷෳࡶͳେม‫ܗ‬ɾྲྀಈ‫ڍ‬ಈͷ

3 ࣮‫ݧ‬ͷ࠶‫ݱ‬ղੳද‫ݱ‬ੑೳΛ༗͢Δͱߟ͑ΒΕΔɽఏҊख๏ͷද‫ݱ‬ੑೳΛ֬ೝ͢ΔͨΊʹɼ໛‫ݧ࣮ܕ‬ͷ࠶‫ݱ‬ղੳΛ࣮ࢪͨ͠ɽ໛‫ݧ࣮ܕ‬͸ਤ 2 ʹࣔ͢Α͏ʹɼਫ࿏্ʹઃஔ͞Εͨ࠭ࢁʹର͠ɼ଄೾૷ஔʹΑ

4

͓ΘΓʹ

Γൃੜͤͨ͞‫ཱݽ‬೾Λিಥͤ͞Δ΋ͷͰ͋Γɼ࠭ࢁຊ‫Ͱڀݚ‬͸౔ਫࠞ߹ମͷෳࡶ͔ͭେ‫ن‬໛ͳ‫ڍ‬ಈΛͷચງɾӡൖɾଯੵΛ‫ؚ‬Ήɽ࠭ࢁ͸઒࠭Λॆరͨ͠ද‫͢ݱ‬Δ͜ͱΛ໨తͱ͠ɼ2 ૚ཻࢠΛ༻͍ͨ‫ݻ‬ӷࠞ

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

z
0.8s

1.1s

4.5s

x

ਤ 3: ղੳ݁Ռʢ্ஈʣͱ࣮‫ࣸݧ‬ਅʢԼஈʣͷൺֱද 1: ࠶‫ݱ‬ղੳʹ༻͍ͨղੳύϥϝʔλParameter

Value

Particle number / Grid
Density ρ s [kg/m3 ]Young’s modulus E [MPa]Poisson ratio ν
Cohesion c [kPa]Internal friction angle ϕ
Initial porosity θ0Particles’ diameter D50 [m]Particle number / Grid
Density ρw [kg/m3 ]Dynamic viscosity µ [Pa · s]
Gravity force [m/s2 ]Grid size [m]Time increment ∆t [s]

4×4×4
2700
20
0.3
0.05
40◦
0.3
1.0 × 10−3
4×4×4
1000
1.0 × 10−3
−9.8
0.01 × 0.01 × 0.01
1.0 × 10−4

߹ MPM Λൃలͤͨ͞ख๏ΛఏҊ͠ɼ࣮‫ݧ‬ͷ࠶‫ݱ‬ղੳΛߦͬͨɽ‫ূݕ‬ͷ݁Ռɼ҆ఆঢ়ଶʹ͋ͬͨ࠭ࢁ͕
ྲྀਫͷӨ‫ʹڹ‬Αͬͯྲྀಈͨ͠‫ʹޙ‬ଯੵ͢Δ༷ࢠ͕ද‫͞ݱ‬Εɼ࣮‫݁ݧ‬Ռͱͷ੔߹ੑ͕֬ೝ͞Εͨɽࠓ‫ޙ‬͸౔ɾ࠭ͷࡐྉߏ੒ଇʹ͍ͭͯߋʹ‫ݕ‬౼Λߦ͍ɼ࣮ࡍ
ͷ౔࠭ࡂ֐΁ͷద༻ੑΛ‫ݕ‬౼ͯ͠ߦ͘‫ܭ‬ըͰ͋Δɽ

44 − Academic Center for Computing and Media Studies, Kyoto Universityࢀߟจ‫ݙ‬[1] Abe, K., Soga, K., and Bandara, S. : Material point
method for coupled hydromechanical problems. Journal of Geotechnical and Geoenvironmental Engineering, 140(3), 04013033, 2013.
[2] Bandara, S., and Soga, K. : Coupling of soil deformation and pore fluid flow using material point method.
Computers and geotechnics, 63, 199-214, 2015.
[3] Kularathna, S., and Soga, K. : Implicit formulationof material point method for analysis of incompressible materials. Computer Methods in Applied Mechanics and Engineering, 313, 673-686, 2017.

‫͍ڧ‬Ի‫ڹ‬ҟํੑΛ༗͢Δ CFRP ʹର͢Δ։‫߹ޱ‬੒๏ͷߴ଎࣮ߦதാ࿨೭ ∗∗

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

1 ͸͡Ίʹ

ߟྀ͠ͳ͕Β೾‫ܗ‬Λ߹੒͢Δ͜ͱͰ͋Δɽ͜ͷ߹੒୸ૉણҡ‫ڧ‬Խथࢷ (Carbon fiber reinforced plastic:

CFRP) ͸ɼۚଐΑΓ΋ܰ͘ɼ߶ੑɼ଱ർ࿑ੑɼ଱৯
ੑʹ༏Ε͍ͯΔɽ·ͨɼ‫ݱ‬৔Ͱൺֱతࣗࡏʹ੒‫ܗ‬Ճ޻Ͱ͖ΔͨΊɼଟ෼໺Ͱར༻͕Ճ଎͍ͯ͠Δɽ੒‫ܗ‬΍Ճ޻͸༷ʑͳํ๏͕ఏҊ͞Ε͍ͯΔ͕ɼ݁‫ہ‬͸୸ૉણҡͱथࢷͷෳ߹ࡐྉͱͳΔͨΊɼεέʔϧຖʹҟͳΔܽؕ΍ෆྑՕॴ͕ൃੜ͢Δ͜ͱ͕໰୊ͱͳΔɽྫ͑͹ɼϓϦϓϨάΛ༻͍Δ৔߹ɼણҡγʔτؒʹണ཭͕ൃੜ͢Δɽણҡ͏ͶΓ΋ CFRP ͷѹॖ‫౓ڧ‬Λେ෯ʹ௿Լͤ͞Δɽӡ༻‫ޙ‬͸ɼিܸ౳ʹΑΔ಺෦ׂΕΛධՁ͢Δ͜ͱ͕ཁ‫͞ٻ‬ΕΔɽैͬͯɼܽؕͷେ͖͞΍Ґஔʹ߹Θͤͯɼ௒Ի೾๏ɼӔྲྀ୳ই๏ɼX ઢಁաࢼ‫ݧ‬๏౳ͷඇഁյ‫ࠪݕ‬๏͕ద༻͞Ε͍ͯΔɽਂ͞ํ޲ͷ͖ͣͷ৘ใΛఆྔతʹਪఆ͢Δํ๏ͱͯ͠
͸ɼ௒Ի೾๏͕༗ޮͰ͋Δɽ͜͜ͰɼϨʔβʔޫΛର৅ࡐྉʹরࣹ͠ɼরࣹ఺ͷ෺࣭͕೤๲ுΛ‫͜͢͜ى‬ͱͰ௒Ի೾Λൃੜ͢Δํ๏͕͋Δɽൃੜͨ͠௒Ի೾
͸ޫԻ‫ڹ‬೾ͱ‫ݺ‬͹ΕɼϨʔβʔ༠‫ى‬Ԡྗ೾͋Δ͍͸ϑΥτϝΧχΧϧ೾ͱ͸۠ผ͞ΕΔ [1]ɽޫԻ‫ڹ‬೾͸ඇৗʹ޿ଳҬͳ೾ಈ͕ࡐྉද໘͋Δ͍͸಺෦Ͱൃੜ
͢ΔͨΊɼද૚ʹ͍ۙ಺෦Ͱ΋ө૾Խ͕ՄೳͰ͋Δɽॲཧʹ͸਺஋‫ࢉܭ‬Λ༻͍Δ͕ɼըૉຖʹԻ଎ͱ։‫ޱ‬ͷิਖ਼Λߦ͏ͨΊɼϝϞϦΞΫηεͱ෼‫ذ‬Λଟ༻͢Δଟॏϧʔϓߏ଄ʹͳΔɽैͬͯɼ૸ࠪϐον͕ࡉ͔͘ɼө૾Խൣғ͕޿ྖҬʹͳΔ΄Ͳɼେ‫ن‬໛ͳ‫ࢉܭ‬ͱͳΓɼͦͷߴ଎࣮ߦ͕‫ٻ‬ΊΒΕΔɽ͜͜Ͱ͸ɼ։‫ޱ‬߹੒๏ͷίʔυͷνϡʔχϯάΛ‫ڞ‬ಉ‫࣮ͯ͠ͱڀݚ‬ࢪ͠ɼ‫౎ژ‬େֶεʔύʔίϯϐϡʔλ Camphor2(γεςϜ AɼCray XC40) Ͱͦͷੑೳ‫ূݕ‬ΛߦͬͨͷͰɼͦͷ੒ՌΛใࠂ͢Δɽ

2

ޫԻ‫ڹ‬೾ͷ‫ܭ‬ଌͱө૾Խख๏ϨʔβʔޫΛ෺࣭ʹরࣹ͢Δͱɼ෺࣭ͷޫֶಛੑʹैͬͯࢄཚ͋Δ͍͸‫ٵ‬ऩ͞ΕΔɽ‫ٵ‬ऩ͞ΕͨΤωϧΪʔͰޫ೤త࡞༻΍ޫԽֶత࡞༻͕ੜ͡Δɽඇৗʹ୹͍ύϧε෯ͰޫΛ෺࣭ʹরࣹ͢Δͱ୹࣌ؒʹ೤๲ு͕ੜ͡ɼ͜ΕʹΑͬͯԠྗ͕ੜ͡Δɽ͜ͷ‫ݪ‬ཧΛར༻ͯ͠ɼ௒Ի೾ྖҬͷޫԻ‫ڹ‬೾Λൃੜͤ͞ΔɽҎԼͰ͸ɼޫԻ‫ڹ‬೾Λޮ཰Α͘ऩू͢ΔͨΊͷ‫ݦ‬ඍ‫ڸ‬γεςϜͱɼͦͷ೾‫ܗ‬Λ߹੒ͯ͠಺෦Λө૾Խ͢Δ
։‫߹ޱ‬੒๏ͷ֓ཁΛड़΂Δɽ͜Ε·ͰɼஶऀΒ͸ɼޫԻ‫ڹ‬ΠϝʔδϯάͷͨΊͷ‫ݦ‬ඍ‫( ڸ‬Photoacoustic MicroscopyɿPAM) Λ։ൃ͖ͯͨ͠ɽ͜ͷ PAM ͸ɼϨʔβʔޫͷরࣹ఺ͱ௒Ի
೾ϓϩʔϒΛ‫ڞ‬য఺Խ [2] ͨ͠΋ͷͰ͋Γɼরࣹ఺Λ૸ࠪ͠ͳ͕Βɼ֤఺ͰޫԻ‫ڹ‬೾Λड৴͢Δɽ͜ͷޫԻ‫ڹ‬೾ͷ֦͕ΓΛߟྀ͠ͳ͕Βड৴೾Λ߹੒͠ɼ

CFRP ͷ಺෦Λ 3 ࣍‫ݩ‬తʹө૾Խ͢ΔͨΊʹɼ։‫ޱ‬߹੒๏ͷίʔυΛ։ൃ͖ͯͨ͠ɽ͜͜Ͱ͸ɼϞϊελςΟοΫʹө૾ԽΛߦ͏։‫߹ޱ‬੒๏Λϕʔεͱ͓ͯ͠Γɼզʑͷख๏ͷಛ௃͸ɼCFRP ͷԻ‫ڹ‬ҟํੑʹ‫ى‬Ҽ͢ΔԻ଎ͷํ޲ґଘੑͱ೾ಈͷ֦͕Γ (։‫ )ޱ‬Λ

2.1

ޫԻ‫ݦڹ‬ඍ‫( ڸ‬PAM)

ஶऀ͕ CFRP ‫ʹ༻ࠪݕ‬։ൃͨ͠‫ݦ‬ඍ‫ڸ‬͸ɼҩྍ༻ͷ Acoustic-resolution PAM(AR-PAM[2]) ͷ‫ڞ‬য఺‫ݪ‬ཧΛࢀߟʹ͓ͯ͠ΓɼͦͷγεςϜͷ֓ཁΛਤ 1(a)ʹࣔ͢ɽϨʔβʔޫ‫͔ݯ‬Βɼ532nm ͱ 1064nm ͷ 2छྨͷ೾௕ͷϨʔβʔޫΛ੾Γସ͑ͯग़ࣹ͢Δɽ‫܁‬Γฦ͠प೾਺͸ 100Hz Ͱ͋Γɼύϧε෯ 4.0ns ͷ୹ύϧεϨʔβʔΛྭ‫͍ͯ͠ى‬ΔɽϨʔβʔޫ͸ɼϑΝΠόΛ௨ͯ͡‫ݦ‬ඍ‫ڸ‬ϢχοτʹૹΒΕͯɼΞΩγίϯϨϯζʹΑͬͯϦϯάύλʔϯΛ‫ܗ‬੒͢ΔɽϦϯ

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 45άঢ়ͷޫ͸ɼ൓ࣹϓϦζϜΛ༻͍ͯ CFRP ද໘ʹूఄ͢ΔͷͰɼޫԻ‫ڹ‬೾ͷ଎౓෼෍Λ g ͱ͢Δͱɼরޫ͞ΕΔɽൃੜͨ͠ޫԻ‫ڹ‬೾͸ɼਤ 1(b) ʹࣔ͢Α͏ࣹ఺ oi ͱࢄཚ‫ ݯ‬r ·Ͱͷ‫ ཭ڑ‬r(= |r − oi |) ͸ɼ࣍ʹɼϓϦζϜͷத৺ʹઃஔͨ͠ूଋ‫ܕ‬ਫਁ௒Ի೾ϓࣜͱͳΔɽϩʔϒʹΑͬͯ‫ݕ‬ग़͞ΕΔɽਫฏํ޲ʹϥελ૸ࠪ(x-y ฏ໘) ͯ͠ɼ‫ݦ‬ඍ‫ڸ‬ϢχοτΛฏߦҠಈͤ͞Δ
͜ͱʹΑͬͯɼ֤఺ͰޫԻ‫ڹ‬೾Λड৴͢ΔɽϓϩʔϒͰಘΒΕͨ৴߸͸֎෦ΞϯϓͰ૿෯͞Εͨ‫ޙ‬ɼϨγʔόͰड৴͞ΕΔɽ͞ΒʹޫԻ‫ڹ‬৴߸͸ DAQ ͰྔࢠԽ͞ΕͯɼPC ϝϞϦʹอଘ͞Ε͍ͯ͘ɽ

2.2 ։‫߹ޱ‬੒๏ (SAFT)r=|g(r)|T2(1)ࣜ (1) ʹैͬͯɼ࣌ࠁ t = T ͷৼ෯஋ V (oi , T ) ͸Ґஔ r Ͱൃੜͨ͠೾ಈͩͱΈͳ͠ɼ͜ͷৼ෯஋ΛҐஔ r ʹ hi (oi , r) ͱͯ͠Ϛοϐϯά͢ΔɽϨʔβʔͷরࣹ఺ i ΛҠಈ͠ͳ͕Βɼ֤ҐஔͰৼ෯஋ hi ΛॏͶ߹ΘͤΔɽH(r) =

N
∑
i=1

(√
)
)2
(
2
i
i
hi (o , r)w(rx − ox ) + ry − oyiΛߦ͏ɽਤ 2 ʹࣔ͢Α͏ʹɼϨʔβʔޫͷরࣹ఺ͷ

(2)͜͜ͰɼޫԻ‫ڹ‬೾ͷ౸ୡ͢Δൣғ (։‫ )ޱ‬ΛԾఆ͠ɼ࠲ඪΛ oi ͱ͢Δɽi ఺ͰϨʔβʔޫΛরࣹ‫ޙ‬ɼಉҐͦΕʹΑΔৼ෯‫౓ڧ‬ͷิਖ਼ w Λߟྀ͢Δɽ͜ͷ w ͸ɼஔͰड৴ͨ͠ޫԻ‫ڹ‬೾‫ܗ‬Λ V (oi , t) ͱ͢ΔɽCFRP ʹ

Laser͜͜Ͱ͸ SAFT ʹΑΔ CFRP ಺෦ͷ 3 ࣍‫ݩ‬ө૾Խ͓͚ΔޫԻ‫ڹ‬೾͸ਫͱͷ઀৮໘ʢCFRP ද໘ʣͰ୎ӽWatero i ( o xi , o iy , o zi )ͯ͠ൃੜ͢Δ͜ͱ͕֬ೝ͞Ε͍ͯΔɽද໘ͰޫԻ‫ڹ‬೾͕ൃੜͨ࣌͠ࠁΛ t = 0 ͱ͢Δͱɼ࣌ࠁ t ʹ T ͷg೾ಈ͸ɼϨʔβʔরࣹ఺Ͱ͋Δ CFRP ද໘Ͱൃੜ͠
ͨޫԻ‫ڹ‬೾͕ɼ಺෦ͷࢄཚମ͔Β໭ͬͯ͘Δ·Ͱͷd࣌ؒͱΈͳ͢͜ͱ͕Ͱ͖ΔɽCFRP ͸‫͍ڧ‬ҟํੑΛ
Q-sw. gate in

Trig. out
Scanner controller

Switching circuit

Correction of
intensity and jitter

DAQ

Beam splitter

Signal in Amp.

Receiverz

MicroscopeFiberਤ 2: ։‫߹ޱ‬੒๏ʹ͓͚Δର৅ըૉ r ͱরࣹ఺ oi ɽ

Branched lightPDr ( r x, r y, r z)

CFRP

Q-sw Nd:YAG Laser
Trig. outPCr

Photoacoustic signal

y

x

z

(a)

| g| =3.1 km/ s

gz

CFRP

gy

Axicon lens

gx

Amp.

Reflecting prism(condenser lens)Ultrasonic transducer

Confocal area

8

8

6

6

4

4

2

gy 0

-2

-2

-4

-4

-6
-8

(b)

ਤ 1: (a) ޫԻ‫ڹ‬೾‫ܭ‬ଌγεςϜͷ֓ཁਤɼ(b) ‫ݦ‬ඍ‫ڸ‬Ϣχοτɽ

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

2

gz 0

-6
-8

-6

-4

-2

0

gx

2

4

6

8

[km/s]

-8
-8

-6

-4

-2

0

gx

2

4

ਤ 3: ҁ৫ CFRP ͷ‫܈‬଎౓෼෍ɽ

6

8[km/s]

Ψ΢ε෼෍Λ࠾༻͍ͯ͠Δɽ
[w(d) = exp −ৼ෯஋ H ΛελοΫ͢ΔɽҎ্ͷ‫ࢉܭ‬Λશͯͷըૉd2(d0 /2)2]Ͱ౳͘͠ߦ͏ɽैͬͯɼ࣮࣭తʹ͸ɼਤ 4 தͷഁઢ(3)

ਤ 3 ʹɼ‫ࣔ͢ʹޙ‬ҁ৫ CFRP ‫ࢼڙ‬ମͷ‫܈‬଎౓෼෍Λ͓ࣔͯ͘͠ɽx-y ฏ໘Ͱ͸ɼͻ͠‫͍ۙʹܗ‬଎౓෼෍Λࣔ͢ɽͰࣔ͞Εͨ࿮ͷϧʔϓ‫͕ࢉܭ‬ϝΠϯͱͳΓɼ(x, y, z)ํ޲ͷըૉͷ਺ (1 ∼ I, 1 ∼ J, 1 ∼ K) ͚ͩ͜ͷϧʔϓΛ‫܁‬Γฦ͢͜ͱʹͳΔɽ֤ըૉͰͷৼ෯஋ H ͷ‫ࢉܭ‬͸ಠཱͰ͋ΔͷͰɼ͜͜Ͱ͸ର৅ྖҬΛ෼ׂ͠ɼ
MPI ฒྻͰ‫ࢉܭ‬Λߦ͍ͬͯΔɽຊίʔυ͸ଟॏϧʔϓͷ‫ͳͱࢉܭ‬Δ͕ɼͦͷաఔʹ͓͍ͯ‫܈‬଎౓ςʔϒϧΛஞҰࢀর͢Δ͜ͱɼ։‫ޱ‬ͷൣғ͔൱͔Λ൑அ͠ɼ

3 ϓϩάϥϜίʔυͷߴ౓Խ൱ͷ৔߹ϧʔϓΛൈ͚Δඞཁ͕͋ΓɼϝϞϦΞΫη

3.1 ΦϦδφϧίʔυͷϑϩʔͱ໰୊఺
։‫߹ޱ‬੒๏Λ࣮ߦ͢ΔͨΊͷϑϩʔΛਤ 4 ʹࣔ͢ɽղੳ৚݅΍σʔλͳͲΛೖྗ‫ޙ‬ɼଟॏϧʔϓͷ‫ࢉܭ‬εͱ෼‫ذ‬Λଟ༻͢Δ‫͋Ͱࢉܭ‬Δɽ·ͨɼৼ෯஋ H ͷ૯࿨‫ࢉܭ‬͸ OpenMP ʹΑΔϦμΫγϣϯॲཧΛ࠾༻
͍ͯ͠Δ͕ɼॳ‫ظ‬ίʔυͰ͸ɼଟ͘ͷίΞΛ࢖༻͢Δ΄Ͳ‫͕཰ޮࢉܭ‬௿Լ͢Δ͜ͱ͕໰୊Ͱ͋ͬͨɽ

͕࢝·ΔɽϨʔβʔরࣹ఺͔Βө૾Խର৅ͷըૉ r·Ͱͷԟ෮఻ൖ࣌ؒ t = T Λ‫ੵݟ‬΋Δɽ͜ͷͱ͖ɼ‫܈‬଎౓ g ͷςʔϒϧ͔ΒҐஔ (ํ޲) ʹґଘͨ͠଎౓

3.2 վྑ఺

Λ‫ͼݺ‬ग़͠ɼ͞Βʹ։‫ิޱ‬ਖ਼ w Λ‫͢ࢉܭ‬Δɽ͜ͷ࣌෼ੳͷ݁Ռɼຊίʔυͷ‫ؒ࣌ࢉܭ‬ͷ΄ͱΜͲΛ઎

ؒ T ʹରԠ͢Δৼ෯஋ hi Λड৴೾͔Βநग़͢Δɽ͜

ΊΔͷ͸ɼ͜ͷৼ෯஋ H ͷ૯࿨‫͋Ͱࢉܭ‬Γɼߏ଄ͱΕ͸֤Ϩʔβʔরࣹ఺ i ʹ͓͚Δ࡞‫͋Ͱۀ‬Δ͕ɼ͜ͷͯ͠͸ըૉ਺ʹ૬౰͢Δ֎ଆ 3 ॏϧʔϓͱɼরࣹ఺

࡞‫ۀ‬Λ x-y ໘ͷরࣹ఺਺෼ (i = 1 ∼ Nx × Ny ) ߦ͍ɼʹؔ࿈͢Δ಺ଆ 2 ॏϧʔϓͷ‫ ܭ‬5 ॏϧʔϓ͔ΒͳΔɽ

解析条件入力配列確保・変数初期化
波形データ入力͞Βʹɼ࠷಺ଆͷॲཧʹෳ਺ͷ IF จ͕‫·ؚ‬Ε͍ͯΔɽ͜ͷ IF จͷ൑ఆ࣍ୈͰ͸ɼͦͷεςοϓΛଈ࣌தஅ͠ɼ࣍ͷεςοϓʹҠߦ (cycle ॲཧ) ͢ΔΑ͏ʹͳ͍ͬͯΔɽ౰ॳ͸ɼγʔέϯγϟϧͳ 1 ࣍‫ݩ‬ϧʔϓʹม‫ͯ͠׵‬ɼҰ࣍‫ݩ‬഑ྻʹͨ͠ৼ෯஋഑ྻ H Λ༻͍ͯ‫͕͍ͨͯ͠ࢉܭ‬ɼ1 ࣍‫ʹݩ‬ม‫͢׵‬Δඞཁੑ͸௿͘ɼ3

群速度テーブルｇの作成࣍‫ݩ‬഑ྻ΁มߋ͠ɼ3 ॏϧʔϓͰ‫͖Ͱࢉܭ‬ΔΑ͏ʹมߋͨ͠ɽॳ‫ظ‬ίʔυͰ͸಺ଆͷরࣹ఺ͷ 2 ॏϧʔ

映像化対象 r の計算．画素数(1 〜 I, 1 〜 J, 1 〜 K) だけ繰り返す
x-y 面の照射点 oi の計算．照射点数(1 〜 Nx, 1 〜 Ny) 分だけ繰り返す

開口補正 w の計算ϓʹରͯ͠ OpenMP ࢦࣔߦ͕ૠೖ͞Ε͍͕ͯͨɼม
ߋ‫ޙ‬ͷίʔυͰ͸࠷֎ଆͷϧʔϓʹରͯ͠ OpenMPࢦࣔߦΛૠೖ͢ΔΑ͏ʹมߋͨ͠ɽ͞Βʹɼ։‫ޱ‬ͷൣғ಺͔൱͔Λ IF จͰ൑ఆ͍͕ͯͨ͠ɼ։‫ޱ‬ͷൣғ֎͔൱͔ʹมߋΛߦͬͨɽ

群速度 g の参照
H(r) の足し込み

3.3 ੑೳൺֱ
‫౎ژ‬େֶεʔύʔίϯϐϡʔλ Camphor2 (γεH(r) の出力 (MPI-IO)

全ノード計算領域分割 MPI 計算

ਤ 4: ։‫߹ޱ‬੒๏ͷϑϩʔνϟʔτςϜ A, Cray XC40) ͷ 1 ϊʔυΛ༻͍ͯɼվྑίʔυͷੑೳ‫ূݕ‬ΛߦͬͨɽίϯύΠϥ͸ Cray FortranCompiler 8.6.5ɼMPI ϥΠϒϥϦ͸ Cray MPI 7.7.0 Λ༻͍ͨɽγεςϜ A ͷ࣮ߦʹ͍ͭͯ͸ɼશͯ cacheϞʔυΛ࢖༻ͨ͠ɽ౰ॳͷϓϩάϥϜͰɼϓϩηε͋ͨΓ 2 εϨουͱͳΔΑ͏ʹ Hybrid ࣮ߦΛͨ݁͠

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.18 No.1（ Oct. 2019）− 47ੂ਋॥‫ॻش‬੝ଐ॥‫ॻش‬








 




৺  ೅











৺  ೅

8x4x8

8x8x4

32x8x1
16x16x1

16x16x1

4x4x1

1t 2t 4t 1t 2t 4t 1t 2t 4t 1t 1t 2t 1t 2t 1t 2t 1t 1t 1t 1t 1t
8x8x1

16x4x1५ঞॵॻਯ 1t
উট७५ਯ






8x4x4







16x8x1

16x8x1





ਈசக

৺  ೅16x4x1HODSVHGWLPH>V@



8x8x2



ਤ 5: εϨουɾϓϩηε਺Λม͑ͨ৔߹ͷ System A Ͱͷ࣮ߦ࣌ؒͷൺֱՌɼશͯ Flat-MPI ࣮ߦΑΓ΋ੑೳ͕௿Լ͢Δ݁Ռͱϧग़ྗ࣌Ҏ֎ʹ௨৴Λߦ͏͜ͱ͕ແ͍ɽ‫ Ͱ఺࣌ݱ‬I/Oͳ͍ͬͯͨɽ͢ͳΘͪ Hybrid Խ͕ෆ‫׬‬શͰ͋ͬͨɽ

ͷ࣌ؒ͸ඇৗʹখ͍͞ࣄ͔Βɼฒྻ਺Λ૿΍͢ࣄʹ

·ͨɼ౰ॳͷϓϩάϥϜ͓͍ͯɼ1 ϊʔυ಺Ͱ 64 ฒ

ΑΔ MPI-IO ͷ࣌ؒ૿Ճ͸ͦΕ΄Ͳ໰୊ʹ͸ͳΒͳྻɼ128 ฒྻɼ256 ฒྻΛଌఆͨ݁͠Ռɼ64 ฒྻͷੑ͍ͱࢥΘΕΔɽೳΛ‫ج‬४ͱͯ͠ɼ256 ฒྻͰ໿ 2 ഒͷੑೳΛࣔͯ͠
͍ͨɽ256 ฒྻ࣮ߦͷ৔߹ɼIntel KNL ͷϋʔυ΢ΣΞεϨου (HT) Λ࢖ͬͨ৔߹ɼL1 ΩϟογϡΛ 4 ϓ

4ϩηε͕ɼL2 ΩϟογϡΛ 8 ϓϩηε͕‫ڞ‬༗͍ͯ͠ΔࣄΛߟ͑Δͱɼੑೳͱͯ͠͸ѱ͘ͳ͍݁ՌͰ͋Δɽ࠷దԽ‫ޙ‬ͷੑೳධՁʹ͍ͭͯɼਤ 5 ʹ݁ՌΛ·ͱΊΔɽന৭ͷ๮άϥϑ͸౰ॳͷίʔυɼάϨʔͷ๮άϥϑ͕վྑ൛ͷੑೳ݁ՌͰ͋Δɽ̍ϊʔυ಺Ͱɼϓϩηε਺΍εϨου਺Λม࣮͑ͯߦͨ݁͠Ռɼͦͷ಺ͷ࠷଎஋͸ 64 ϓϩηεʷ 4 εϨουͷ৔߹Ͱ͋ΓɼࠇͰ͍ࣔͯ͠Δɽ౰ॳίʔυʹ͓͚Δ࠷଎஋Ͱ͋ͬͨ 256 ϓϩηεฒྻ (Flat MPIɼ16x16x1) ͷ݁Ռͱൺֱ͢Δͱɼվྑ൛ͷ࠷଎஋ (64 ϓϩηεʷ 4 εϨουɼ16x4x1) ͸໿ 1.77 ഒͷੑೳ޲্ͱͳ͍ͬͯ
Δɽ64 ฒྻ (16x4x1) ͷ݁Ռʹ͍ͭͯ஫໨ͯ͠ΈΔͱɼॳ‫ظ‬ίʔυͰ͸ Hybrid ࣮ߦ͸ੑೳ͕௿Լ͢Δ݁
Ռͱͳ͍͕ͬͯͨɼվྑ൛Ͱ͸ 4 εϨου·ͰεέʔϥϏϦςΟ͕‫ݟ‬ΒΕ͓ͯΓɼ໿ 1.81 ഒͷੑೳ޲্ͱ

CFRP ಺෦ͷө૾Խྫ͜͜Ͱ͸ɼCFRP ͷ‫ࢼڙ‬ମʹ͍ͭͯɼޫԻ‫ڹ‬೾Λ༻͍ͨө૾ԽྫΛࣔ͢ɽ‫ࢼڙ‬ମ͸ਤ 6 ʹࣔ͢Α͏ʹɼҁ৫ͷ 8 ૚ߏ଄Ͱ͋Γɼް͞͸ 2.12mm Ͱ͋Δɽલ
અ 2 Ͱࣔͨ͠ޫԻ‫ݦڹ‬ඍ‫ڸ‬Λ༻͍ͯɼվྑίʔυʹ
ΑΔ SAFT ͷө૾ԽΛߦͬͨɽ಺෦ͷ 3 ࣍‫ݩ‬ө૾Խ݁ՌΛਤ 7 ʹࣔ͢ɽՄࢹԽൣғ͸ɼx ํ޲ 20mmɼyํ޲ 10mm ͷൣғͰ͋Γɼ‫ݦ‬ඍ‫ڸ‬ͷεΩϟϯϐον
͸ 0.02mm ͱͨ͠ɽ͜͜Ͱ SAFT Ͱ͸ɼըૉαΠζ͸ 0.025mmɼϏʔϜύϥϝʔλ d0 =0.2mm ͱͯ͠‫ܭ‬ࢉΛߦͬͨɽਤ 7 ͸ௗᛌਤɼಉਤ (b) ͸ද໘ɼಉਤ(c) ͸ x − z அ໘Λද͍ͯ͠Δɽਤ 7 ͷന৭ͷೱ౓͸ɼࣜ (2) ͷ H ͷେ͖͞Λද͍ͯ͠Δɽ͜ͷ݁Ռ͔Βɼද໘ͷ໛༷͸໪࿦ɼ૚ؒߏ଄·Ͱ΋͕઱໌ʹ࠶ߏ੒͞Ε͍ͯΔͷ͕Θ͔Δɽ

ͳ͍ͬͯΔɽ·ͨɼಉ͡ Flat-MPI ಉ࢜Ͱ΋໿ 1.92 ഒͷੑೳ޲্Λࣔͨ͠ɽ64 ϓϩηεͷ৔߹ɼվྑ൛ͷ࠷଎஋͸౰ॳͷ໿ 3.49 ഒͷੑೳ޲্ͱͳͬͨɽ

2 ϊʔυҎ্Λ࢖ͬͨ৔߹ͷੑೳධՁʹ͍ͭͯ͸ɼຊίʔυ͸࣌ؒൃలΛ͠ͳ͍ͷͰɼ‫ج‬ຊతʹϑΝΠ5·ͱΊࠓճͷ‫ڞ‬ಉ‫ʹڀݚ‬Αͬͯɼಛʹ hybrid ฒྻͷੑೳ޲্͕‫ݦ‬ஶͰ͋ͬͨɽ͜ΕʹΑͬͯɼεέʔϥϏϦςΟ͕େ෯ʹվળ͞Εɼ౰ॳͷίʔυΑΓ΋ 3.5 ഒ

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

Thickness: 2.12mm

x

x
y

y

Imaging area

10mm

20mm

ਤ 6: ҁ৫ CFRP ͷࣸਅͱө૾Խର৅ྖҬɽ(a)͍ۙߴ଎Խ͕࣮‫͖ͨͰݱ‬ɽࠓ‫ͳߋޙ‬ΔվྑΛࢼΈΔxz

ͱ͢ΔͳΒ͹ɼҰ൪ͷϗοτεϙοτͱͳ͍ͬͯΔ
ৼ෯஋ H ͷ଍͠ࠐΈ‫ॲࢉܭ‬ཧʹ͍ͭͯ‫ݕ‬౼͢Δඞཁ͕͋Δɽ‫Ͱ఺࣌ݱ‬͸ɼ͜ͷॲཧ͕ϕΫτϧԽ͕ग़དྷ͍ͯͳ͍ࣄ͕ੑೳ޲্ͷ્֐ཁҼͷҰͭͱߟ͑ΒΕyΔɽߏ଄తʹେ෯ͳॻ͖‫͕͑׵‬ඞཁͱͳΔ͕ɼ‫ܧ‬ଓͯ͠‫ݕ‬౼Λߦ͍͍ͨɽ(b)ँࣙɿ ༗ӹͳ͝ॿ‫ݴ‬Λ௖͍ͨʮεʔύʔίϯϐϡʔxλγεςϜ‫ڞ‬ಉ‫اڀݚ‬ըҕһձʯͷօ༷ɼฒͼʹίʔυνϡʔχϯάʹ͍ͭͯ௚઀͝ࢦಋ௖͍ͨΫϨΠɾδϟύϯɾΠϯΫͷླ໦ ޾࿕ࢯʹ‫ँײ‬ਃ্͛͠·͢ɽzຊ‫ؔʹڀݚ‬࿈͢Δ੒Ռ:(1) Kazuyuki Nakahata, Kazuki Karakawa, Keiji Ogi,
Koichi Mizukami, Satoshi Wada, Takeshi Namita, Tsuyoshi
Shiina, Three-dimensional SAFT imaging for anisotropic
materials using photoacoustic microscopy, Ultrasonics,
Vol.98, pp.82-87, 2019.
DOI:10.1016/j.ultras.2019.05.006.
(2) தാ࿨೭, ԫ໦‫ܠ‬ೋ, ਫ্޹Ұ, େฏࠀ‫ݾ‬, ‫ࢁؙ‬ਅ޾, ࿨ాஐ೭, ࿘ా ݈, ௣໊ ‫ؽ‬, ޫԻ‫ڹ‬๏ʹΑΔ୸ૉણҡิ‫ڧ‬थࢷͷද໘௚Լണ཭ͷ 3 ࣍‫ݩ‬Πϝʔδϯά,ి‫ֶؾ‬ձ࿦จࢽ C, Vol.139, No.2, pp.142-148, 2019.

DOI:10.1541/ieejeiss.139.142(c)ਤ 7: ҁ৫ CFRP ͷ಺෦ө૾Խɽ(a) ௗᛌਤɼ(b) ද໘ɼ(c)x − z அ໘ɽࢀߟจ‫ݙ‬[1] ࠤ౻ढ़Ұ, খ‫ݪ‬መ, ޫͱ௒Ի೾ɾѹྗ೾ͷෳ߹త࡞༻Λར༻ͨ͠ҩྍٕज़ͷਐల, ޫֶ, Vol.38, No.6, pp.288-297,
2009.
[2] H. F. Zhang, K. Maslov, G. Stoica, L. V. Wang, Functional photoacoustic microscopy for high-resolution andnoninvasive in vivo imaging, Nature Biotechnology,
Vol.24 (7), pp.848-851, 2006.

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

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

サービス再開日時

障害発生日時

サービス再開日時

2018/11/20

2018/11/20

保守時間[h]2018/10/08

9:00

2018/10/12

16:55

103.92

2018/12/06

9:00

2018/12/07

11:30

26.50

2019/03/27

9:00

2019/03/31

24:00

111.00

15:40

ダウン時間[h]18:50２）サービス状況
ジョブ

サービ
ス時間
[h]

処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼
動
ノード数

10 月

640.0820,785130,153

53,174,864

35,588,686

1574.2

69.8

%

11 月

716.8324,320201,414

67,614,272

48,505,617

1799.3

79.8

%

12 月

717.5034,942

252,345

65,035,363

45,751,924

1794.8

79.1

%

1月

744.0031,000251,81471,854,980

53,598,300

1796.1

78.7

%

2月

672.0017,049

182,801

64,135,474

47,964,320

1799.0

77.5

%

3月

633.0017,619153,333

60,698,374

47,870,627

1736.674.4%

計4,123.41

145,715

1,171,860

382,513,327

279,279,474

1750.0

76.6

%

ノード
利用率

ノード利用率
100%
90%
80%
70%
60%
50%
40%
30%
20%
10%
0%

10月

11月

12月

1月

2月

3月



占有時間 = 合計(経過時間×占有コア数)

平均稼動ノード数 = 電源 ON 状態のノード数の月平均 (10 分間隔のサンプリングデータより算出)



ノード利用率 = 稼動ノードに対するジョブが実行されているノードの割合

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

3.17

システム B 運転状況 ( 2018 年 10 月 ～ 2019 年 3 月 )１）保守作業に伴うサービス休止およびシステムダウン障害発生状況
保守作業に伴うサービス休止
保守開始日時

サービス再開日時

システムダウン障害発生状況
保守時間[h]障害発生日時

サービス再開日時

ダウン時間[h]2018/10/08

9:00

2018/10/10

13:00

52.00

2018/11/20

15:40

2018/11/20

18:50

3.17

2018/12/06

9:00

2018/12/07

10:07

25.12

2019/01/30

23:20

2019/01/31

10:49

11.48

2019/03/27

9:00

2019/03/31

24:00

111.00２）サービス状況
ジョブ

サービ
ス時間
[h]

処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼
動
ノード数

10 月

692.0099,188657,152

13,536,136

11,067,731

830.8

71.5

%

11 月

716.83

198,589765,685

15,102,805

11,910,870

832.6

73.4

%

12 月

718.88

190,644791,775

13,954,791

11,056,737

821.7

69.5

%

1月

732.52

171,032

893,723

15,522,607

11,919,752

839.2

74.2

%

2月

672.00115,571753,756

14,083,765

11,086,523

839.3

72.9

%

3月

633.0052,080624,961

12,167,455

9,999,868

788.3

57.5

%

計

4,165.23

827,104

4,487,052

84,367,559

67,041,481

825.3

69.8

%

ノード
利用率

ノード利用率
100%
90%
80%
70%
60%
50%
40%
30%
20%
10%
0%

10月

11月

12月

1月

2月

3月



占有時間 = 合計(経過時間×占有コア数)

平均稼動ノード数 = 電源 ON 状態のノード数の月平均 (10 分間隔のサンプリングデータより算出)



ノード利用率 = 稼動ノードに対するジョブが実行されているノードの割合

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

システム C 運転状況 ( 2018 年 10 月 ～ 2019 年 3 月 )１）保守作業に伴うサービス休止およびシステムダウン障害発生状況
保守作業に伴うサービス休止
保守開始日時

サービス再開日時

システムダウン障害発生状況
障害発生日時保守時間[h]2018/10/08

9:00

2018/10/10

13:00

52.00

2018/12/06

9:00

2018/12/07

10:07

25.12

2019/03/27

9:00

2019/03/31

24:00

111.002018/11/20

サービス再開日時

15:40

2018/11/20

ダウン時間[h]18:50２）サービス状況
バッチ

サービ
ス時間
[h]

処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼
動
ノード数

10 月

692.003,856

29,887

242,712236,811

15.8

43.3

%

11 月

716.83

2,629

21,859281,109245,698

16.0

47.5

%

12 月

718.88

3,811

24,153

290,756

258,624

15.6

45.8

%

1月

744.0088,968100,266

347,785

292,843

15.9

56.4

%

2月

672.0046,37164,565

318,751

262,690

16.0

50.1

%

3月

633.0027,53061,587

257,685

202,681

15.1

36.0

%

計

4,176.71

173,165

302,317

1,738,798

1,499,347

15.7

46.5

%

ノード
利用率

ノード利用率
100%
90%
80%
70%
60%
50%
40%
30%
20%
10%
0%

10月

11月

12月

1月

2月

3月



占有時間 = 合計(経過時間×占有コア数)

平均稼動ノード数 = 電源 ON 状態のノード数の月平均 (10 分間隔のサンプリングデータより算出)



ノード利用率 = 稼動ノードに対するジョブが実行されているノードの割合

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

3.17

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

資料

別表１ スーパーコンピュータシステム
提供サービス
利用負担額
コース

タイプ

セット

エントリ

タイプA
タイプB
タイプC

基本
基本
基本
基本
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位
最小
追加単位

パーソナル

タイプA1
タイプA2
タイプA3
タイプB1
グループ
タイプB2
タイプB3
タイプC1
タイプC2
タイプA
大規模ジョブ タイプB
タイプC
専用クラスタ-ストレージ容量追加
ライセンスサービス

システム バッチ

12,600 円/年
100,000 円/年
100,000 円/年
100,000 円/年
200,000 円/年
100,000 円/年
240,000 円/年
60,000 円/年
600,000 円/年
300,000 円/年
210,000 円/年
105,000 円/年
252,000 円/年
63,000 円/年
630,000 円/年
315,000 円/年
130,000 円/年
130,000 円/年
156,000 円/年
78,000 円/年
20,000 円/週(7日）
10,000 円/週(7日）
21,000 円/週(7日）
10,500 円/週(7日）
13,000 円/週(7日）
6,500 円/週(7日）
630,000 円/年
315,000 円/年
10,000 円/年
20,000 円/年

システム資源

経過時間 ストレージ
無料
（時間）
（TB)
利用者数

最大1ノード相当((36コア、128GBメモリ)×1)
1
最大4ノード相当((68コア、16+96GBメモリ)×4)
168
最大4ノード相当((36コア、128GBメモリ)×4)
最大1ノード相当((72コア、3072GBメモリ)×1)
4ノード((68コア、16+96GBメモリ)×4)
優先
2ノード((68コア、16+96GBメモリ)×2)
8ノード((68コア、16+96GBメモリ)×8)A
準優先
2ノード((68コア、16+96GBメモリ)×2)
8ノード((68コア、16+96GBメモリ)×8)
占有
4ノード((68コア、16+96GBメモリ)×4)
4ノード((36コア、128GBメモリ)×4)
優先
2ノード((36コア、128GBメモリ)×2)
336
8ノード((36コア、128GBメモリ)×8)
B
準優先
2ノード((36コア、128GBメモリ)×2)
8ノード((36コア、128GBメモリ)×8)
占有
4ノード((36コア、128GBメモリ)×4)
1ノード((72コア、3072GBメモリ)×1)
優先
1ノード((72コア、3072GBメモリ)×1)
C
2ノード((72コア、3072GBメモリ)×2)
準優先
1ノード((72コア、3072GBメモリ)×1)
8ノード((68コア、16+96GBメモリ)×8)A
占有
4ノード((68コア、16+96GBメモリ)×4)
8ノード((36コア、128GBメモリ)×8)B
占有
168
4ノード((36コア、128GBメモリ)×4)
2ノード((72コア、3072GBメモリ)×2)
C
占有
1ノード((72コア、3072GBメモリ)×1)
8ノード((36コア、128GBメモリ)×8)
B
4ノード((36コア、128GBメモリ)×4)
ストレージ容量10TBの追加につき可視化ソフト（AVS,ENVI/IDL)およびプリポストウェアの1ライセンスにつき
B
A
B
C

共有
共有
共有
共有

0.2
3.0
3.0
3.0
24.0
12.0
28.8
7.2
48.0
24.0
24.0
12.0
28.8
7.2
48.0
24.0
24.0
24.0
28.8
14.4

-

48.0
24.0

8
4
16
4
16
8
8
4
16
4
16
8
8
8
16
8

-

16
8

備考
1． 利用負担額は、年度単位で算定している。また、総額表示である。パーソナルコース、グループコース又は専用クラスタコースを、
年度途中から利用を開始する場合及び年度途中で利用を終了する場合の利用負担額は、上記表中の利用負担額を12で除した後、利用月数を乗じて算出するものとし、100円未満に端数が出た場合は、10円単位を四捨五入するものとする。
なお、月途中から利用を開始する場合及び月途中で利用を終了する場合は、それぞれ1月の利用とする。
2． 大型計算機システムの全ての利用者は、上記表のサービスの他、次のサービスを受けることができる。
1） 大判プリンタサービス
2） その他、大型計算機システムが提供するサービス、機器の利用
3． 上記表の大規模ジョブコース、ストレージ容量追加、ライセンスサービスの申請には、スーパーコンピュータシステムの利用者であることが必要である。
4． 「共有」： 当該カテゴリのユーザ間で一定の計算資源を共有するベストエフォートのスケジューリングを行う。
「準優先」： 定常稼働状況において記載値（以上）の計算資源が確保されるように優先スケジューリングを行う。
また、稼働状況によらず記載値の1/4の計算資源が確保されることを保証する。
「優先」： 定常稼働状況において記載値（以上）の計算資源が確保されるように優先スケジューリングを行う。
また、稼働状況によらず記載値の1/2の計算資源が確保されることを保証する。
「占有」： 稼働状況によらず記載値の計算資源が確保されることを保証する。
5． ストレージ容量はバックアップ領域（最大で総容量の1/2）を含む。
6． グループコース及び専用クラスタコースの利用者番号は利用者あたり年額5,000円を負担することで追加できる。
7 ． 機関・部局定額制度
他機関又は学内における部局（『国立大学法人京都大学の組織に関する規程』第3章第2節から第11節で定める組織をいう。）の組織が、その組織単位でグループコースサービスを利用申請する場合の利用負担額は、別表１に規定する
1.5倍の額とする。なお、利用負担額が年額150万円未満の場合は100人、年額150万円を超える場合は、150万円毎に100人までの利用者を認める。ストレージは、1.5倍の容量とする。
8 ． スパコン連携サービス
学術情報メディアセンターのスーパーコンピュータシステムと密な連携により、学内における部局の組織が計算サーバ等を設置する
場合、下記の負担額を支払うものとする。
冷却方式
水冷
空冷

利用負担額
9,800 円/月
11,500 円/月

利用負担額算定単位
水冷冷却方式の計算サーバ等の定格電力 1kWにつき
空冷冷却方式の計算サーバ等の定格電力 1kWにつき

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

別表２ 汎用コンピュータシステム
区

分

利用負担額

仮想サーバホスティングサービス

単位

36,000円／年 1仮想サーバにつき

備考
1． 利用負担額は、総額表示である。
2． 上記表の仮想サーバホスティングサービスを利用するには、スーパーコンピュータシステムの利用者であること。
3． 1仮想サーバに割当てるシステム資源は、CPU：2コア、メモリ：4GB、ディスク：100GBである。
4． 仮想サーバホスティングサービスにおいて、下記の負担額を支払うことによりCPU、メモリ、ディスクを増量することができる。
区

分

利用負担額

単

CPU増量

3,000円／年 2コアにつき(最大8コアまで)

メモリ増量

3,000円／年 4GBにつき(最大64GBまで)

ディスク増量6,000円／年 100GBにつき(最大1,000GBまで)位

5． 利用負担額は、当該年度(4月から翌年3月まで)の利用に対して年額として算定するが、年度途中から利用を開始する場合には月数に応じて減額する。

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

※民間機関利用

別表３ スーパーコンピュータシステム
システム

A

B

C

システム資源

経過時間 ストレージ
無料
（時間）
（TB)
利用者数

利用負担額

8ノード(68コア、16+96GBメモリ)×8)

336

28.8

16

960,000 円/年

12ノード(68コア、16+96GBメモリ)×12)

336

43.2

24

1,440,000 円/年

16ノード(68コア、16+96GBメモリ)×16)

336

57.6

32

1,920,000 円/年

8ノード(36コア、128GBメモリ)×8)

336

28.8

16

1,008,000 円/年

12ノード(36コア、128GBメモリ)×12)

336

43.2

24

1,512,000 円/年

16ノード(36コア、128GBメモリ)×16)

336

57.6

32

2,016,000 円/年

2ノード(72コア、3072GBメモリ)×2)

336

28.8

16

624,000 円/年

3ノード(72コア、3072GBメモリ)×3)

336

43.2

24

936,000 円/年

4ノード(72コア、3072GBメモリ)×4)

336

57.6

32

1,248,000 円/年

備考
1． 利用負担額は、年度単位で算定している。また、総額表示である。パーソナルコース、グループコース又は専用クラスタコース
年度途中から利用を開始する場合及び年度途中で利用を終了する場合の利用負担額は、上記表中の利用負担額を12で除し利用月数を乗じて算出するものとし、100円未満に端数が出た場合は、10円単位を四捨五入するものとする。
なお、月途中から利用を開始する場合及び月途中で利用を終了する場合は、それぞれ1月の利用とする。
2． ストレージ容量はバックアップ領域（最大で総容量の1/2）を含む。

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

全国共同利用版広報・Vol.17(2018)総目次［巻頭言］
Vol.17，No.1 号の発刊にあたって ................................................................................................................1-1
Vol.17，No.2 号の発刊にあたって ................................................................................................................2-1［スーパーコンピュータ共同研究制度（若手・女性研究者奨励枠）研究報告］
高精度視覚質感記憶の心理学的基盤と神経機構の解明...............................................................................1-2
阿蘇山の複雑地形における高精度メッシュを用いた領域大気モデルの感度解析及び性能評価.............1-4
破壊力学に基づく損傷モデルによる鉄筋コンクリートの 3 次元破壊シミュレーション ........................1-6
高効率・高耐久性色素増感太陽電池の実現に向けた新規ポルフィリン色素の開発と電子構造の解明 .1-8
気液混相撹拌操作に対する数値解析 ...........................................................................................................1-10
iPS 細胞懸濁培養の最適化を目的とした撹拌槽内粒子挙動解析..............................................................1-12
パワーデバイス用半導体製造装置設計、最適化に関する数値解析 .........................................................1-14
カスケード選択型分子動力学シミュレーションの開発とタンパク質構造変化予測への適用...............1-16
粒子との接触を伴う液体挙動の直接数値解析 ............................................................................................1-18
無重力下での高プラントル数流体における温度差マランゴニ効果に起因する液柱内対流場の二次不安定性 .............................................................1-21
視聴覚モダリティ間デコーディングによる感覚間協応のメカニズムの解明..........................................1-24
随伴解析を用いた物体表面形状最適化による抵抗低減.............................................................................1-26
アンブレラサンプリングを利用した自由エネルギー反応経路探索法の開発..........................................1-28
タンパク質フォールディング駆動力の確率的解析 ....................................................................................1-30
津波による底泥巻き上げ量の予測と海洋環境変化に関する数値解析 .....................................................1-32
分子動力学シミュレーションを用いた小分子・高分子混合溶液の相溶性の計算手法の確立...............1-34
撹拌時における混合評価に対する数値シミュレーション.........................................................................1-36
1 次元量子スピン鎖の非平衡定常状態における温度勾配の数値的研究 ..................................................1-38
非一様な都市構造物上における大気乱流の組織構造に関する数値解析..................................................1-40
半無限領域のスペクトル法による竜巻を模した渦の数値実験に向けた研究開発 ..................................1-42
高周波分数冪ラプラシアン Navier-Stokes 方程式のエネルギースペクトルの考察 ..............................1-45［プログラム高度化共同研究報告］
動的/静的水～土骨格連成有限変形解析コードの高度化～任意形状に対する領域分割法の適用の検討～ ........................................................................................1-47
高分子分離膜の大規模シミュレーション....................................................................................................1-49
非均質異方性材料評価のための陽的ボクセル有限要素解析の高度化 .....................................................1-52
超小型マイクロ波放電式中和器の電子引き出しを対象とした 3 次元プラズマ粒子計算 ......................1-56［講習会参加報告］
Cray XC40 プログラミング講習会-Scaling DL Training Workloads with the Cray PE Plugin.........2-2［バーストバッファによる高速ファイル I/O］

バーストバッファの利用について～ DataWarp 編～ ................................................................................2-6

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

［サービスの記録・報告］
スーパーコンピュータシステムの稼働状況...................................................................................... 1-60，2-9
センター利用による研究成果（平成 29 年度）..........................................................................................2-12［資料］
大型計算機システム利用負担金 別表...............................................................................................1-64,2-16
全国共同利用版広報・Vol.16(2017)総目次 .................................................................................................1-67
サービス利用のための資料一覧 ..........................................................................................................1-69,2-19［奥付］
奥付 ........................................................................................................................................................1-70,2-20

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

―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

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

奥付京都大学学術情報メディアセンター全国共同利用版広報 Vol．18，No．1
２
０１９年
３１日xx
発行
2019
年 １０月
10 月
日発行
広報編集部会
編集者 京都大学学術情報メディアセンター

深沢 圭一郎（部会長）

全国共同利用版広報編集部会

平石 拓（副部会長）

発行者 〒606-8501 京都市左京区吉田本町

水谷 幸弘

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

熊谷 真由美Academic Center for Computing and Media Studies

尾形 幸亮

Kyoto University
Tel. 075-753-7414
http://www.media.kyoto-u.ac.jp/
印刷所 〒616-8102 京都市右京区太秦森ヶ東町 21-10
株式会社エヌジーピー

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

表紙デザイン：中山 豊（中山商店）