巻頭言

Vol.21, No.1 号の発刊に当たって京都大学学術情報メディアセンター

牛島 省

本号では，京都大学学術情報メディアセンターにおいて実施された 2021 年度の共同研
究報告について特集いたします．2021 年度は，12 件の「若手・女性研究者奨励枠」と 2 件の「プログラム高度化共同研究」による共同研究が実施されました．
学術情報メディアセンターでは，スーパーコンピュータ共同研究制度として，「若手・女性研究者奨励枠」，「大規模計算支援枠」，「プログラム高度化共同研究」という３種類の研究支援事業を実施しています．課題募集は年度開始前後に行われます．また，年度途中で追加募集が実施される場合があります．応募された課題の採否は，スーパーコンピュータ共同研究企画委員会にて審査されます．本共同研究制度への応募方法や申込み時期等の詳細につきましては，http://www.iimc.kyoto-u.ac.jp/ja/services/comp/support/kyodo/
をご参照下さい．
上記の共同研究制度のうち，「若手・女性研究者奨励枠」の応募資格は，(1) 40 歳未満の若手研究者 (学生を含む, 性別は問わない)，あるいは，(2) 女性研究者 (年齢は問わない)，とされています．この「若手・女性研究者奨励枠」では，スーパーコンピュータを利用することで学術的にインパクトがある成果を創出できる課題に対して，計算機利用負担金の全額または一部を本センターが負担しています．さらに，「若手・女性研究者奨励
枠」は，JHPCN (学際大規模情報基盤共同利用・共同研究拠点) の活動の一環と位置づけられており，優れた課題は，JHPCN の萌芽型共同研究へ推薦されます．推薦された課題
は，JHPCN シンポジウムにてポスター発表を行うことが可能となるなどの特典がありますので，是非この共同研究制度をご活用ください．
本号に掲載された 12 件の「若手・女性研究者奨励枠」の研究成果報告では，複雑流体
および流体構造連成問題，プラズマ計算，各種の分子構造・分子動力学計算などに加え，ニューラルネットワークや機械学習に関する課題など多様な研究が行われ，全国 7 大学の若手・女性研究者の研究をサポートいたしました．なお，JHPCN の萌芽型共同研究には9 件の課題が採択されました．
「大規模計算支援枠」は，一定規模の計算資源を 2 週間程度ご利用頂くことが可能な共同研究制度で，利用負担金は本センターが全額負担するものです．2021 年度は応募がありませんでしたが，大規模計算の実施を希望される方はご検討ください．
また，「プログラム高度化共同研究」では，京都大学のスーパーコンピュータをグルー
プコースまたは専用クラスタコースでご利用の研究者 (JHPCN および HPCI にて京都大学のスーパーコンピュータを利用する研究者も含む) を対象に，プログラムの制御構造・
データ構造の改良による性能チューニングや並列化手法の改良，問題分割・負荷分散方式などの改良による大規模計算プログラムの高度化・高性能化などを支援しています．2021
年度は，(1) 高密度領域まで適用可能なモンテカルロ法の開発と有限密度 2 カラー QCDの相図の決定，(2) 低メモリ・低コスト通信乱流大規模直接数値計算コードの高度化，という 2 件のプログラム高度化支援を行いました．これらの共同研究では，プロセス並列の通信の効率化やスレッド並列計算の高度化などの検討を加え，応募者に満足いただける性
能向上を実現しました．
今後も皆様の研究，教育にご活用いただけるようにセンター教職員も尽力していきます
ので，ご利用・ご支援のほど、よろしくお願い申し上げます．

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

特集

新スーパーコンピュータシステムの紹介
深沢 圭一郎京都大学学術情報メディアセンター

京都大学では、1969 年に FACOM230 をスーパ

5.8PFLOPS 程度ということになっています。見かーコンピュータという名称ではなく汎用コンピュ

け上の性能はあまり前システム A と比べて向上しータとして導入したことが、大規模計算機システムていませんが、HBM2e の効果により、これまでよ

導入の始まりになっています。当初はその時代のス

り高い実行性能となるアプリケーションが多いとーパーコンピュータの代表的な構成ということも考えています。

あり、ベクトル型 CPU を搭載したシステムが導入システム B（Laurel3）ではシステム A と同じ新されていましたが、2000 年に HPC2500 というス

型 Xeon ですが、DDR5 を 512GB/ノード搭載した

カラ型 CPU（SPARC64V）を初めて導入し、2007

システムです。システム A の 4 倍のメモリサイズ、

年にはいわゆる T2K スパコンとして 4 ソケット構

DDR 系ということで低メモリアクセスレイテンシ

成の x86 型 CPU を導入しています。その次のリプとなっています。システム全体では 370 ノードとな

レイスでは、現在の構成に続くメニーコア型のシス

り、前システム B と比べてノード数が半分以下とな

テム A、汎用型のシステム B、大規模共有メモリ型っていますが、ノード当たりの性能が AVX-512 対のシステム C が始まっています。今回のリプレイス

応やコア数増加により向上しており、2PFLOPS とにおいても、この考え方を踏襲し、メニーコア型かいう性能となっています。

つ高メモリバンド幅のシステム A、汎用型かつ高メシステム C（Cinnamon3）では、システム B のモリサイズのシステム B、大規模共有メモリ型のシ構成にメモリを 2TB/ノードと増加させ、大規模共有

ステム C の導入を計画しました。更に機械学習向け

メモリシステムとしています。こちらは前回と同じとして小規模ですが、GPU といったアクセラレー

で 16 ノードとなっています。

タを搭載したシステムの導入も計画しました。スト

今回新しく導入するシステム G（Gardenia）としレージシステムでは、今回はバーストバッファではて、NVIDIA A100（80GB）を 4 枚/ノード搭載し、なく、SSD で構成された高速アクセスが可能なストCPU は Milan 世代の EPYC×2 と DDR4 512GB

レージの導入も計画しました。このようなシステムとなっています。全体として 16 ノードと小規模での導入を目指して、調達を行っていましたが、新型

すが、主として 1 ノード利用を想定しており、ノーコロナウイルス感染症に影響された半導体不足、ウ

ド内 GPU メモリが 320GB（ノード内メインメモリ

クライナ危機による物価高や円安などにより、調達は 512GB）と多くあることで、ノード内で閉じる機

時期が遅れてしまい、利用者の方には大変ご迷惑を械学習向けシステムとしています。

おかけしています。この記事を執筆時点（10 月末）

ストレージでは HDD による 40PB と SSD によ

では、新スパコンの契約が終わり、稼働に向けて全る 4PB の容量を持つシステムを導入し、計算によ力で準備をしているところです。

るデータ出力だけでなく、様々なデータを集める場

この新スパコンですが、まとめると図 1 のような所としての利用も想定しています。データ収集にはシステムになりました。システム A（Camphor3）

VM ホスティングサーバを利用することができます。

には新型 Xeon と高メモリバンド幅を持つ HBM2e

今回の新システムとVMホスティングサーバとの間を搭載したノードとなっており、ノード当たりの性

は、物理的に直結するネットワーク構成を踏襲して

能は 6TFLOPS に近い性能となっています。システ

おり、これまでと同様に効率的にご利用いただける

ム全体では 1120 ノードとなっており理論性能は環境となっています。

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

図 1 新スーパーコンピュータシステム概要
利用開始時期に関しては、情報環境機構の次期ス

Slurm は世界で多く利用されているスケジューラ

パコンの Web ページ（https://www.iimc.kyoto-となっています。これまでの京都大学のスパコンジ

u.ac.jp/ja/services/comp/supercomputer/other/com

ョブスケジューラで行えていたことは、Slurm でも

ing_2022.html）を確認いただけると最新情報が公実現できるように準備していますので、一部スクリ開されていますが、現在はシステム A と SSD によ

プト記述書式の変更が出るかもしれませんが、これる高速ストレージ以外を 2023 年 1 月頃に試験運用までと同様の機能はカバーされています。

を開始し、すべてのシステムの運用開始を 2023 年

一方でアプリケーションですが、近年のライセン

7 月に予定しています。このように計算資源が無い

ス料の高騰に伴い、スパコン導入アプリケーション

時期が長く続いてしまうため、クラウド資源の利用の整理を行いました。できる限り削減をしないためを準備しています。新スパコンではクラウドバースに、一部アプリケーションでは永続ライセンスへのティングと呼ばれる計算資源のクラウド環境への変更（サポートがなく、バージョンアップに制限あ

拡張機能を持たせていました。これを利用し、2022り）し、SAS は導入を見送りました。SAS に関して年の 11 月中には本年度 7 月まで利用いただいたスは個人利用ライセンスのご購入や他センターのごパコンユーザのみなさまへクラウド環境による計利用をお願いします。

算機資源の提供を行う予定です。これ以外にも何ら

今回のリプレイスでは、世界的な様々な状況によかの計算資源の提供や新システムのできる限りのり導入システムの資源削減や導入遅れなど現在も早期提供を考えています。

ご迷惑をおかけしております。本来であれば、一部

これまではハードウェアの紹介でしたが、ジョブシステムの性能評価などをご報告したいところでスケジューラやアプリケーションに関しても新ス

したが、実機が無いため、性能に関しては次号へのパコンでは変更があります。ジョブスケジューラで掲載を予定しています。利用者の期待に応えられるは、前システムでは PBS、その前は NQS と変更がシステムを全力で準備していますので、運用開始ま続いていますが、今回も変更があり Slurm と呼ばれでお待ちください。リプレイスへのご理解ご協力どるオープンソースのスケジューラになりました。

うぞよろしくお願いいたします。

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

シンポジウム開催報告

α×SC2022K 群れとスーパーコンピュータに関するシンポジウム
開催報告
深沢 圭一郎京都大学学術情報メディアセンター

1 α×SC とは
京都大学学術情報メディアセンターと九州大学
情報基盤研究開発センターでは「なにか×スパコン（SC）」というテーマを掲げることにより、その「なにか」の分野とスパコンの分野で新しい共同研究が生まれることを期待したシンポジウムを年二回開催しています。
2018 年から開催しており今回が 8 回目の開催になります。様々な分野にコンピュータが利用されている現在では、スパコンセンターのあり方もこれまでとは変わる必要があります。そこで、これまでのスパコン利用者だけでなく、新しい分野にスパコンを使ってもらい、より良い成果やこれまでにない成果の創出に繋がるきっかけ作りが必要となると考え、シンポジウムを開催しています。
今回は、2022 年度の京都大学担当回として、「なにか」を「群れ」としたα×SC2022K 群れとスー
パーコンピュータに関するシンポジウムを 2022 年
9 月 1 日 13 時 30 分～18 時 15 分に開催しました（2022K の K は Kyoto の K です）。新型コロナウイルス感染症への対策として 2021 年度は一部講演者のみ現地参加とし、基本的にはリモートでの参加という形でしたが、2022 年度は現地参加も可能なハ
イブリッド開催としました。
今回テーマを「群れ」とした理由は、筆者が研究代表者として、京都大学研究支援 SPIRITS 2021-

図 1 α×SC2022K のポスター
1 のようなポスターを作成することができ、広報に活用できました。参加者は現地 12 名、リモートでは 48 名となり、現地参加者は少ないですが、この 2年程度の間にリモートでの参加が便利であり、会議などへの通常の参加方法になってきている結果で2022 年度「プラズマ粒子シミュレーションの応用にはないかと思われます。

よる野生ウマの行動数値実験モデルの確立」という2 シンポジウムの講演

プロジェクトに採択されており、群れに関する研究者と繋がりを持てたことによります（このプロジェ
クトとの共催としています）。また、九州大学情報基盤研究開発センターとの共催でもあります。
SPIRITS との共催のおかげで、昨年度に引き続き図4 − Academic Center for Computing and Media Studies, Kyoto University

シンポジウムでは、4 名の研究者に「群れ」に関連するご講演をいただき、最後に全体を通した総合討論を行いました。シンポジウムの講演資料などは

図 2 リモートでの開催風景（上津原さんのご講演時）
Web ページで公開していますので、ご興味のある方た。動物や人の集団の動きを時空間の単位要素へ分は ご 覧 く だ さ い 。 Web ペ ー ジ の URL は割する動的モード分解を利用した映像分析手法のhttps://sites.google.com/kyoto-u.ac.jp/sc2022k/ にご紹介がありました。

なります。リモートの開催風景は図 2 のような様子

総合討論では、シンポジウムの講演内容に関する追加の議論やシンポジウムのテーマであるスパコでした。
シンポジウムでは、まず、村上 久（京都工芸繊維

ンをどう活用するか、自分の研究にどういう機能が大学情報工学・人間科学系）助教に、「互いに動きをあれば使いやすいのかなどの議論を行いました。分

読み合うこと＝相互予期」に基づく群れの自己組織野が異なる群れに関係する研究者が集まった場で化に関してご講演をいただきました。村上先生は、したので、お互いの研究分野との比較など、群れにイグノーベル賞を 2021 年に受賞されており、相互関してこれまでに余りない視点からの議論もあり

予期の研究はその受賞内容でもあります。続いて、ました。実際の群れを観測した後、解析する流れがSPIRITS プロジェクトの共同研究者でもある前田

一般的ですが、近年のドローンやセンサなどを活用

玉青（京都大学大学院理学研究科生物科学専攻 博

し、定量的にデジタル解析が進んでいるために、今士後期課程）さんが、ドローンによる観察が明らか

後スパコンを利用した解析が力を発揮しそうな一にするウマの重層社会という題目で講演をしまし

方で、使いやすさ（利用申請含む）には改善がまだた。重層社会とは、ヒト社会のように、小さな群れ

まだ必要だとコメントをいただきました。

が集合して高次の集団を作る社会であり、その様子
をドローンでの観測と定量的解析を行っていると

3 シンポジウムを終えていうご講演でした。
休憩を挟んだ後、後半では衛星の群れということ

今回は SPIRITS プロジェクトとの共催によるシで、小型 SAR 衛星コンステレーションのご講演をンポジウムの宣伝協力やハイブリッド開催により、

上津原 正彦（株式会社 QPS 研究所）開発部部長に

例年より参加者が多くなりました。他のシンポジウ行っていただきました。QPS 研究所では 2025 年以

ムでも感じますが、総合討論や、質疑では会場にい降に合成開口レーダー(SAR)を有する小型衛星 36

る関係者からの発言が多く、議論を活発に行いたい

機でコンステレーションを構築計画しており、その場合は現地開催を強化する必要があると感じまし

運用に関してご紹介いただきました。次に、高見 利

た。また、現地参加者間で、共同研究が始まったと也（大分大学理工学部）教授による、動的モード分の報告もあり、スパコンを基盤に異分野を上手く連解による群れの運動分析というご講演がありまし

携させることを今後も目指していきます。

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

Numerical simulation of Marangoni convection in a shallow rectangularcavity with a linear solutal boundary condition
ZHANG Jiangao
Graduate School of Engineering Science, Osaka University

1 Introduction

values, Ch and Cl, are set at the four boundary
corners, and the linear solutal boundary condition

Marangoni convection usually occurs in a liquid
layer with a free surface driven by the variation ofis applied, as shown in Fig. 1(b).
The governing equations of the fluid flow in thefree surface tension due to temperature gradients

cavity are written in dimensionless forms:

or concentration gradient. This convection has

Continuity: V = 0(1)

Momentum:V+ V  V = −P +  2V


(2)


1
+ V   = −  2 

Sc

(3)

been of interest to researchers due to its rich
dynamical features and its wide existence innatureandindustrial

processes,

such

as

oceanography, droplets, thin-film coatings, andcrystal growth.

Mass transport:

The dimensionless concentrations is defined as on

= (C - Cl)/(Ch - Cl),  and P are the dimensionless

Marangoni convection have been reported [1-3].

time and pressure. Sc = D/ is the Schmidt number,

However, most of the studies only considered thewhere  is the dynamic viscosity and D is thepure thermal Marangoni effect, while the solutal

diffusion coefficient of working fluid.

Marangoni convection plays an important role instudy, calculation has been performed for the fluid

coating and crystal growth. To the best of ourwith Sc = 10.

knowledge, there are few works focused on the(a) 3D view

Numerous

fruitful

research

findings

solutal Marangoni flow. It would be an interesting

In thisphenomenon to investigate since it occurs in someprocesses. In this work, we perform a series ofthree-dimensional numerical simulations to study
pure solutal Marangoni convection in a shallow
rectangular cavity to shed further light on flow
characteristics and flow pattern transition.

(b) top view

2 Numerical method
We consider a cavity filled with fluid, as shownin Fig. 1. The top boundary is the free surface,
while the other boundaries are the cavity walls.
The Cartesian coordinate system is adopted. The
length, width, and depth of the cavity are L, L, andFig. 1 The computational domain and the

0.1L, respectively. High and low concentration

boundary conditions.

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

3 Results and discussion

corners, as seen in Fig. 3(a). Also, as the solutal
Marangoni

numberbecomeslarger,the

When the Marangoni number is relatively small,concentration drops increase, and the solutal

the solutal Marangoni convection is steady state,

Marangoni effect enhances, resulting in a stronger

which is called the basic flow hereafter. Figure 2

fluid flow. Furthermore, every peak of the curves

shows the streamlines and the iso-concentrationin Fig. 3(b) coincides with each vortex-center

lines of the typical basic flow at Sc = 10. The

position.

surface fluid flows from the high concentration(a)spot (A) to the low concentration spot (C) with acertain curvature, and then returns back near thebottom due to the mass conservation, as shown in
the internal field of Fig. 2(a). The Schmidt numberof the working fluid in the present study is muchlarger than unity. It means that the relative
contribution of the advective mechanism is more
dominant than that of the diffusion. Thus, the(b)mass transfer is greatly sensitive to the flow region
and the iso-concentration lines are drasticallydistorted on the free surface, as shown in Fig. 2(b).
A

A

Internal field

D



z

z

C
y

y

x

x

y
x

(a) streamlines of the internal field (left) andAD

4 Conclusion



z

x

y

velocity U (b) along the diagonals on the freesurface.

boundaries (right)l field

C Fig. 3 Distributions
Bof concentration  (a) andThe basic flow is steady state and the flowyxCstructure remains symmetric along the diagonals.

B

Large concentration drops appear in the vicinity of(b) iso-concentration lines on the top free surface

the boundary corners so that the effect of solutal

Fig. 2 A typical snapshot of the streamlines (a) andMarangoni in the middle part of the free surface isconcentration field (b) of the basic flow.

weaker than at the corners，which greatly affectFigure 3 shows the concentration and velocity
distributions along the AB and CD diagonals at
different

solutalMarangoni

numbers.

The

positions of AB and CD lines are shown in Fig. 2(b).
Large concentration drops appear near theboundary corners due to the low mass diffusivity
so that the concentration gradient in the middle
part of the cavity is smaller than that of thethe flow structure in the whole system.

5 Reference[1] M. K. Smith, Phys. Fluids, 29 (1986) 3182.
[2] L. Peng, et al., Int. J. Heat Mass Transf., 50(2007) 872.
[3] T. Yano, et al., Adv. Space Res., 66 (2020) 2047.

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

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

1 緒言
C8H17OnBu
nBu

ポルフィリンは 400~450 nm に Soret 帯と呼ばれる強い吸収と 550~600 nm に Q 帯と呼ばれる中程

S

度の吸収をもち、増感色素として有望である。特に、ポルフィリンにドナー・アクセプターを導入

OC8H17
N
N

Zn

tBu

C8H17OnBu
nBu

することでも光捕集性能を向上させることができ、で中程度の変換効率が達成されてきた[2]。さらに
我々は、メチレン鎖を介してチオフェン環を縮環させた色素 T-iPr を開発し、高いエネルギー変換
効率を実現することに成功している（Figure 1）[3]。

COOH

T-Ph

を超える変換効率が達成されている[1]。また、ポ実際に芳香環を直接縮環させた色素を用いること

N

OC8H17

集能の改善に伴って光電変換効率が向上し、10%
ルフィリンに芳香環を縮環させて π 共役系を拡張

S

N

C8H17Oしたプッシュ−プル型色素は長波長領域での光捕NNS
C8H17OtBu

OC8H17
N
N

Zn

N

N

S

N
COOH

N
OC8H17
T-tBuPh

Figure 2. Molecular structures of porphyrin dyes.

また、ポルフィリン色素の末端に嵩高い置換基を導入して電解液中のレドックス対が酸化チタン表

2 結果と考察

面に接近することを抑制することで、酸化チタン−
レドックス対間の逆電子移動を抑制でき、エネル
ギー変換効率の向上が可能であることを見出している。そこで本研究では、ポルフィリン色素に嵩高い置換基を導入した色素 T-Ph および T-tBuPhを設計・合成した（Figure 2）。

2.1 ポルフィリン色素の物性と太陽電池性
能評価
ポルフィリン色素 T-Ph, T-tBuPh の THF 中での吸収スペクトルを Figure 3 に示す。これら 2 つの色素間にはほとんど差がなく、色素末端に導入した置換基の嵩高さは光捕集特性に影響を与えない

C8H17OnBu
nBu

S
C8H17O

OC8H17
N
N

Zn

N

N

S

ことがわかった。次に、色素増感太陽電池を作製

N
COOH

N
OC8H17

Figure 1. Methylene-bridged fused porphyrin dye T-iPr.

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

し、セル性能評価を行ったところ、色素としてT-Ph と T-tBuPh を用いたセルでの光電変換効率はそれぞれ 7.36%、7.43%となった。これは、置換基の嵩高さが増大することで酸化チタン−レドッ
クス対間の逆電子移動が抑制され、光電変換効率の向上につながったと考えられる。

ε / 105 M–1 cm–1















Wavelength / nmFigure 3. UV/vis absorption spectra of T-Ph (solid line)
and T-tBuPh (dashed line) in THF.
現在、光電変換効率の向上を目指し、色素増感
太陽電池セル作製条件の最適化を進めている。また、嵩高さを増大させて逆電子移動を効果的に抑制できるポルフィリン色素を設計・合成すること
で、さらなるエネルギー変換効率の向上が実現できると期待できる。
Figure 4. Selected Kohn−Sham orbitals for (a) T-Ph, and2.2 理論計算(b) T-tBuPh obtained by DFT calculations with the

ポルフィリン色素の最安定化構造およびそのフロンティア軌道の電子構造について知見を得るた
めに、Gaussian09 プログラムを用いて密度汎関数
法 （ DFT ） に よ る 理 論計 算 を 行 っ た（B3LYP/6-31G(d)）。ポルフィリン色素 T-Ph,
T-tBuPh の HOMO、LUMO の準位は同じであり、それぞれ–4.57 eV, –2.67 eV であった。末端に導入したフェニル基や t-ブチルフェニル基にはほとん
ど HOMO や LUMO の分布が見られず、ポルフィリンの電子状態には影響を与えないことがわかった（Figure 4）。このため、光捕集特性に違いが見られなかったと考えられる
また、T-tBuPh の t-ブチルフェニル基は T-Ph におけるフェニル基より嵩高く、ポルフィリン色素の横幅と同程度にまで張り出しているため、酸化チタン上に吸着している色素の隙間を通ってレドックス対が酸化チタン表面に近づくことを抑制で

B3LYP/6-31G(d) level. To simplify the calculations,
octyl and octyloxy groups on the meso-aryl groups were
replaced with methyl and methoxy groups, respectively.

3 参考文献[1] A. Yella, C.-L. Mai, S. M. Zakeeruddin, S.-N.
Chang, C.-H. Hsieh, C.-Y. Yeh, M. Grätzel, Angew.
Chem. Int. Ed. 2014, 53, 2973–2977.
[2] H. Hayashi, A. S. Touchy, Y. Kinjo, K. Kurotobi, Y.
Toude, S. Ito, H. Saarenpää, N. V. Tkachenko, H.
Lemmetyinen, H. Imahori, ChemSusChem 2013, 6,
508–517; T. Higashino, H. Imahori, Dalton Trans.
2015, 44, 448–463.
[3] Y. Kurumisawa, T. Higashino, S. Nimura, Y. Tsuji,
H. Iiyama, H. Imahori, J. Am. Chem. Soc. 2019, 141,
9910–9919.

きると考えられる。その結果、酸化チタン−レドッ
クス対間での逆電子移動を抑制でき、T-tBuPh を用いた場合に太陽電池性能の向上に繋がったと考えられる。このように、本計算結果は実験結果の理論的解釈の一助となった点で意義がある。

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

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

1 はじめにン間距離を反応座標に採用する．[1, 2] この座標は均一系だけでなく不均一系でも利用可能である．

埋もれた界面における電荷移動などの化学反応プログラムパッケージは cp2k を用いて，水中のは生体機能や触媒反応など様々な分野で重要であ

水分子の自己解離反応をシミュレーションした．

る．近年分光などの実験技術の向上により，その自由エネルギー曲面の計算ではアンブレラサンプ

分子レベルでの反応メカニズムが明らかになってリングを使用し、1 次元および 2 次元の平均力ポきた．理論計算の高精度化が求められている．しテンシャルを計算した．水の構造を再現するにはかし，実験を解析できる高精度な計算手法は発展

DFT 汎関数の選定も重要である．そこで、本研究途上である．そこで本研究では，第一原理計算を

では水 256 分子と revPBE-D3 汎関数を用いた．

用いて荷電した界面のモデリング手法を開発し、この汎関数では純水の動径分布関数をよく再現す

その構造やダイナミクスを解析する．特に，本研ることが知られており，水の pKa（pKw）もよく

究では純水中のベンチマーク計算を行ったのち、再現することが期待されるためである．[3, 4] 続

アモルファスシリカ・水界面などにおける構造やいて，不均一系への応用例としてアモルファスシ

ダイナミクスを調べる．

リカ・水界面に本手法を応用し、界面分光などの2 手法

実験で観測された pKa と比較した．

3 結果

まず、第一原理計算により水中の水分子の自己
解離のシミュレーションでプロトン移動の反応座

配位数とイオン間距離を反応座標に取った 2 次標を検討した．先行研究ではプロトン移動反応を

元平均力ポテンシャル曲面を計算したところ，計記述する反応座標は結合長や配位数などに限られ算された水分子の自己解離の pKa は 13.7 となりていたが、実験の pKa の定量的な再現だけでなく実験値の 14.0 と良い一致を示した．本手法により遷移状態やプロトン移動後のイオン対の状態を定高精度（実験と定量的に比較できる）に反応自由

性的にも再現することが困難であった．pKa だけ

エネルギー曲面を計算できることがわかった．一でなく活性化自由エネルギーも実験と定量的に比方で，活性化自由エネルギーは 20.9 kcal/mol と較できる計算手法が必要であるが，配位数はそのなり実験値の 23.9 kcal/mol を約 3 kcal/mol 過小

シンプルさ・扱いやすさから幅広い計算プログラ評価した．これは一般化勾配近似によるものだと

ムパッケージで利用可能であるため，多くの計算いうことが分かっており，ハートリーフォック交例がある．また，配位数にはパラメーター依存性

換を加えて revPBE0-D3 汎関数などを用いるこなどの問題があるため新たな反応座標の開発が必となどよって制度を改善できることが示せた，化要であった．加えて，界面へなどの不均一系にお学精度（約 1 kcal/mol）が求められる計算の場合，

いても応用できるものが望ましい．

配位数だけでなくイオン間距離も加えて 2 次元のそこで，本研究ではアニオンとカチオンのイオ

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

平均力ポテンシャルを用いる必要があるとわかった．
 続いて，アモルファスシリカ・水界面での計算
では，表面のシラノール基の密度が 4.9 OH/nm2となり実験値と一致した．しかし，約 1.5 nm2 の表面積を持つスラブではシラノール基の密度の標準偏差が非常に大きいため，1 つのスラブを用いる先行研究とは異なりいくつかのスラブの平均を用いる必要があることも判明した．シラノール基の界面となる角度は界面とほぼ平行となることも判明し，先行研究と多少異なることもわかり今後の界面分光で注意する必要があることが分かった．

参考文献[1] T. Joutsuka, J. Phys. Chem. B, 126 (24),4565–4571, 2022.
[2]

T.

Joutsuka,

ChemRxic,

DOI:

10.26434/chemrxiv-2022-nfl1j.
[3] T. Joutsuka and K. Ando, J. Phys. Chem.
B, 124 (38), 8323–8330, 2020.
[4] T. Joutsuka and K. Ando, Chem. Lett., 50(7), 1325–1328, 2021.
[5] S. K. Reddy, R. Thiraux, B. A. W. Rudd, L.

計算で得られたシラノール基の pKa は 9.4（ケイLin, T. Adel, T. Joutsuka, F. M. Geiger, H. C.

酸 水 溶 液 ）、 7.7 (isolated silanol), 6.9-10.4

Allen, A. Morita, and F. Paesani, Chem, 4 (7),(geminal silanol), 5.4-12.6 (vicinal silanol) と

1629–1644, 2018.

なり実験と良い一致を示した．本研究でも 2 次元

[6] T. Joutsuka, T. Hirano, M. Sprik, and A.

の平均力ポテンシャルを用いており実験を定量的

Morita, Phys. Chem. Chem. Phys., 20 (5),に再現したことから本手法の有効性を示すことが3040–3053, 2018.
[7] T. Joutsuka and A. Morita, J. Phys. Chem.

できた．

C, 122 (21), 11407–11413, 2018.
[8] T. Joutsuka, H. Yoshinari, and S.

4 おわりにYamauchi, Bull. Chem. Soc. Jpn., 94 (1), 106–
111, 2021.

界面での脱プロトン化をシミュレーションする[9] T. Joutsuka, Ryu Hamamura, Kakeru

ため配位数とイオン間距離による 2 次元の平均力

Fujiwara, Tetsuo Honma, Masahiko Nishijima,

ポテンシャルを用いてpKa を計算する手法を構築

and Shohei Tada, Int. J. Hydrog. Energy, 47,

した．今後は，ディープラーニングを組み合わせ

21369–21374, 2022.

ることにより計算コストを削減する予定である．

[10] Shohei Tada, Nagomu Ochiai, Hiroka

以上で計算した界面について，現在界面分光スペ

Kinoshita,

Mitsuhiro

Yoshida,Natsumiクトルの計算法の構築をしている．申請者の先行

Shimada,TatsuyaJoutsuka,

Masahiko

研究[5-8]をさらに拡張し，分光スペクトルからわ

Nishijima, Tetsuo Honma, Noriko Yamauchi,

かる界面構造と光触媒反応メカニズムを議論する

Yoshio Kobayashi, Kenta Iyoki, ACS Catal., 12,

予定である．最終的に，第一原理分子動力学（MD）

7748–7759, 2022.

シミュレーションを用いて固液界面をシミュレー
ションし，界面分光の１つである和周波発生（SFG）分光法を解析する．
以上に加えて，表面における二酸化炭素の水素化触媒のモデリング[9-10]も進めている．このような系に対しても本研究で開発した手法を応用していく予定である．
謝辞
本研究では、京都大学学術情報メディアセンタ
ーの共同研究制度（若手奨励枠）を活用させていただきましたので、御礼申し上げます。

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.21 No.1（ Dec. 2022）− 11三次元的な MHD プラズマダイナミクスの解明

橋本 翼
九州大学大学院理学府地球惑星科学専攻

1 はじめにとするためスーパーコンピュータを用い計算した。
ジョブはハイブリッド並列を用い高速化を図って地球を取り巻く宇宙環境は、太陽風ー磁気圏ーいる。

電離圏（S-M-I）といった大きさや性質の異なる領域が結合されている。それぞれの領域は独立したものではなく、磁場や電流を介して結ばれ、諸3 結果現象があった際はもう片方の系へとその情報が伝達されることで、物理ダイナミクスが進展する。

シミュレーション解析により、IMF の方向不連

そのため、諸現象を包括的に理解することが重要

続面は磁気圏の昼側から夜側に移動するまでに、となる。

複数回の磁力線つなぎ代わり現象（リコネクショ

 本研究は、太陽風―磁気圏(S-M)領域において、ン）を繰り返し、磁気圏の磁場構造やプラズマ分

磁場の向きが北向きから南向きに不連続に変化す

布を大きく変化させているということが分かった。

る方向不連続面を持つ惑星間空間磁場

太陽風にともない IMF が北向きから南向きに不(Interplanetary Magnetic Field; IMF)が到達し連続に変化することは、衛星磁場観測からも確認たときに生じる諸現象の再現と、そこで生じる物

されている。その IMF が地球磁気圏に到達した

理ダイナミクスの解明を目的にした。IMF が実際とき、地球磁場による減速から北向き IMF と南にこのような状況になることは観測的にも確認さ

向き IMF は急速に近づきリコネクションを起これている。

す（図１上）。先行研究で、北向き IMF と南向き

S-M を複合的に再現できるのが数値シミュレ

IMF がリコネクションを起こすことは示唆されーションであるが、広い宇宙環境において各物理ていた（Maynard et al.,2002[2]; Samsonov et al.,量をグリッドごとに短いタイムスケールで計算す

2018[3]）が、本研究のシミュレーションによってるためには大規模な計算処理能力が必要になる。

それを再現した。さらに本研究では、このリコネ

そのため、スーパーコンピュータを利用し研究をクションにより変化した磁場構造が、夜側まで移実施した。

動する過程で複数回リコネクションが発生していることを明らかにした（図１下）。すなわち、リコ
ネクションによる磁力線のつなぎ代わりによって、2 研究手法IMF や地球起源の磁力線が複雑化し、大きく磁気圏構造の変化をもたらすことが分かった（図２）。

シミュレーションに関して、本研究では太陽風

このような磁力線構造変化は過去のシミュレーシから磁気圏までを磁気流体とし、磁気圏内部境界

ョンや地上・衛星観測では指摘されたことがない。

と し て 薄 層 電 離 圏 を 採 用 し た

また本研究では、昼側におけるリコネクションや

REPPU(Reproduce Plasma Universe) コ ー ド

高緯度側で発生するリコネクションがそれぞれ８(Tanaka, 2015[1])を用いた。このコードでは 50

分間というタイムスケールで行われていることがもの物理量を解像度約 300×300×300 で計算、分かった。

および再現するため大規模な計算処理能力を必要

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

4 おわりに本研究では、太陽風―磁気圏(S-M)領域において、IMF の方向不連続面が磁気圏周辺に到達した
際に複数回のリコネクションが生じることを、シミュレーションにより明らかにした。今後はこのリコネクションにともない放出されるリコネクシ
ョンフローや磁気張力・磁気圧の時間変化についてその結果を可視化することで、発動しているダ
イナミクスの解明を行う。また、リコネクションの継続時間が 8 分間というタイムスケールであったことに対し、この数字がどの IMF 状況下でも一

般的なものであるかを検証するため、IMFBy 成分
や動圧といったリコネクション効率を変化させう
る物理パラメータの値を変え、細かくシミュレー
ションを行う。

5 謝辞
 本研究の計算の一部は、京都大学のスーパーコ
ンピュータ共同研究制度（若手・女性奨励枠）に
図１：IMF 不連続面にともなうリコネクション、

基づくものであり、ここに感謝の意を表す。

ラインは磁力線を示し、背景の色はリコネクションの発生の目安となる磁気拡散を示す。
（上図：IMF 不連続面間で生じる、北向き IMFと南向き IMF のリコネクション、

6 参考文献下図：IMF 不連続面と地球ローブ磁力線）[1] Tanaka, T. (2015). In Auroral Dynamics
and Space Weather (pp. 177–190).
Hoboken, NJ: John Wiley & Sons, Inc.
[2] Maynard, N. C., Sonnerup, B. U. Ö,
Siscoe, G. L., Weimer, D. R., Siebert, K. D.,
Erickson, G. M., et al. (2002). Journal ofGeophysical Research, 107(A12), 1456.
[3] Samsonov, A.A., Sibeck, D.G.,
Dmitrieva, N.P., Semenov, Geophys. Res.
Lett. (2017)

図 2：図１の概念図、ラインは磁力線を示す。

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

Numerical investigation for hydrodynamic characterization of orbitally
shaken bioreactor in a suspension cell culture
WANG Liya
Graduate School of Engineering Science, Osaka University

1 Introductionof cells behavior, shaking velocity and shaking
radius still need for clarification. To this end, theOrbitally shaken bioreactor (OSR) is one of thenumerical investigations on the OSR with amost frequently utilized bioreactor for thehollowcultivation of stem cells, because they are

investigating the oxygen supply ability, cells

manageable and easy-to-operate. However, thesuspension amount and shear stress at various

common OSR hardly meets the anticipatedoperational conditions.

demands since cells migrate to the center and2 Numerical methodbottom of the tank due to the “Einstein's tea leaf
paradox” [1,2]. It is well-known that the selectionof a suitable bioreactor is especially important forbiomedical and biotechnological development.
Thus,

a

novel

structureandgeometry

improvement of OSR is urgently developed.
Nowadays, many researchers have made great
efforts to facilitate the innovation of culturebioreactors for laboratory studies and industrial
applications. Horiguchi et al. developed a novel
disposable bioreactor with an O-shaped dish andan O-shaped bag, neither of which contains a

wall

are

furthercarriedoutbyThe bioreactor rotates in an orbital manner.
Therefore, the coordinates of the hollow vessel
center is expressed as Eq. (1). The bioreactor
which does not contain a central region is
illustrated in Figure 1. The outer diameter of thetank is 84 mm and the tank height is 42 mm. The
culture medium is assumed as water and its
height keeps a fixed value 15 mm. The ratio of theinner diameter and outer diameter, shaking

central region. The cells, therefore, are prevented

velocity (ω), and shaking radius (R) are the crucialfrom gathering at the bottom center region. In

variables, which play a key role in the hollow

their experiment study, the cell aggregates canorbitally shaken bioreactor.

grow in the newly designed tank with similar
diameter to aggregates without the evidence of
necrotic cores. Because the O-shaped bioreactor
could provide sufficient oxygen to the cellaggregates[3].

Furthermore,

Zhu

et

al.

investigated the effect the ratio between the innerdiameter and outer diameter on the flow field,
mass transfer and shear stress in the hollow OSR
using numerical simulation. They suggested thatthe diameter ratio should be kept lower to achieve
high oxygen transfer rate [4]. However, the effects

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

Fig. 1 Numerical simulation domain and shaking
trajectory.
𝑥𝑥
𝑅𝑅 cos(𝜔𝜔𝜔𝜔)
𝒓𝒓 = (𝑦𝑦) = (−𝑅𝑅 sin(𝜔𝜔𝜔𝜔))
𝑧𝑧
0

(1)

The governing equations for the liquid phase(air and water) are the Naiver Stokes equations.
In the multi-phase configuration, continuous
species transfer (CST) method that is analogous tothe VOF model, can be used to catch the oxygen
concentration (C) in each phase [5], written as:
𝜕𝜕𝜕𝜕
𝜕𝜕𝜕𝜕

+ ∇ ∙ (𝐶𝐶𝒗𝒗) = ∇ ∙ (𝐷𝐷∇C) − ∇ ∙ (𝐷𝐷

𝐻𝐻−1

𝐻𝐻𝐻𝐻+(1−𝛼𝛼)

(2)

𝐶𝐶𝐶𝐶)

The percentage of afloat particles can be

Fig. 2 Oxygen concentration distribution at aexpressed by a nondimensional parameter. Here,vertical cross section for the hollow bioreactor

we introduce the suspension ratio PS defined as:(upper panel) and conventional bioreactor (lower

PS=NS/NA(3)

where NS is the number of floating particles, andNA is the total number of particles.
The shear stress acting on the particles is an
essential factor to qualify the condition of cellsuspension

culture.

Applying

Stokes’approximation for flow around a sphere, themagnitude of shear stress over the surface of aspherical particle can be estimated in the followingform:|𝜏𝜏| =(4)

𝑑𝑑p

new bioreactor, the suspension ratio (PS) and shear
stress (|τ|) are further quantitatively evaluated.
It is found that the suspension ratio is greatly
increased without the increasing of the shear
stress, which is highly expected.

4 Conclusion

suspension ability in the cell culture processing.

The oxygen mass transfer capability including isa crucial factor in the design, optimization andscale-up of the cell culture process. The results
Equation

here would offer a better process environment forsustaining adequate oxygen supply as well as

3 Results and discussion

using

To get a comprehensive understanding of theGenerally, the new bioreactor design developed

2𝜇𝜇f |𝒗𝒗f −𝒗𝒗p |calculated

panel)

(2)

for

two

configurations are given in Figure 2. Resultspresent a conspicuous distinction between thesetwo cases. It is evident that the highest oxygen
concentration is observed at the gas-liquidThis confirms that the bioreactor designed without
an inner wall can be considered as a better
candidate for iPS cell cultivation.

5 Reference[1] Yamamoto T, et al., J. Chem. Eng. Jpn., 51(2018) 423–430.
[2] Horiguchi I, et al., the Society of Chemicalinterface in both bioreactors. However, it also

Engineers, Japan 81st Annual Meeting, E204

should be noted that the concentration is relatively(2016).

poorer (blue region) at the central region for theconventional bioreactor. This is because theabsence of the inner wall makes it easier for theflow to become turbulent. Consequently, it can bestated that the new bioreactor design introduced
here optimizes the culture condition, particularly[3] Horiguchi I, et al., J. Vis. Exp., e57922 (2019).
[4] Zhu, L, et al., J. Chem. Technol. Biotechnol.,

94 (2019) 2212–2218.
[5] Marschall H, et al., Chem Eng Sci 78 (2012)
111–127.

enhancing the oxygen mass transfer.

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

水-土二相連成解析とマルチスケールモデリングによる相転移型地盤崩壊シミュレータの構築辻 勲平
九州大学大学院工学府建設システム工学専攻，日本学術振興会特別研究員 DC2

1 はじめに係数を示す．式(1)は𝜀𝜀𝜀𝜀をパラメータとして，地表流と浸透流の両方に対応する式であり，影響範囲内に土津波による防波堤の浸透・洗掘破壊や，豪雨により発生した土石流や斜面崩壊，河川堤防の決壊など，
水と土の混相流による自然災害が多発している．水と地盤に関する混相流計算には，FEM や MPM を用いて地盤を連続体として扱うものが多い．一方，
浸透流によって地盤内で局所的に生じる内部侵食の解析には，土粒子の挙動を直接表現する個別要素
法 DEM を用いた数値解析が実施されている．本研
究では，豪雨による堤防決壊，津波による港湾構造物の被災など，水と土砂礫が入り混じる固液混相流
による自然災害を再現し得る数値シミュレータの開発を目指した．

粒子がない場合（𝜀𝜀𝜀𝜀 = 1）は通常の Navier-stokes 方程式，土粒子が存在する場合（𝜀𝜀𝜀𝜀 < 1）は，右辺第 4 項以降の抵抗力項が作用する拡張ダルシー則を解く．

2.2 DEM による地盤の変形解析
土粒子の挙動は DEM で計算し，球もしくは複数の球を剛体結合した要素で土粒子の形状を表現す
る．各粒子・計算ステップ毎に運動方程式を解くことで捨石マウンドの変形を表現する．
𝑚𝑚𝑚𝑚𝑠𝑠𝑠𝑠

𝑑𝑑𝑑𝑑𝒗𝒗𝒗𝒗𝑠𝑠𝑠𝑠= 𝑚𝑚𝑚𝑚𝑠𝑠𝑠𝑠 𝒈𝒈𝒈𝒈 − 𝛻𝛻𝛻𝛻𝛻𝛻𝛻𝛻𝑉𝑉𝑉𝑉𝑠𝑠𝑠𝑠 + 𝒇𝒇𝒇𝒇𝑑𝑑𝑑𝑑 + � 𝒇𝒇𝒇𝒇𝑐𝑐𝑐𝑐
𝑑𝑑𝑑𝑑𝐷𝐷𝐷𝐷
𝐼𝐼𝐼𝐼

2 解析手法𝑑𝑑𝑑𝑑𝝎𝝎𝝎𝝎𝒔𝒔𝒔𝒔= � 𝒎𝒎𝒎𝒎𝑐𝑐𝑐𝑐 − 𝒎𝒎𝒎𝒎𝑟𝑟𝑟𝑟
𝑑𝑑𝑑𝑑𝐷𝐷𝐷𝐷

ここで，𝑚𝑚𝑚𝑚𝑠𝑠𝑠𝑠 , 𝑉𝑉𝑉𝑉𝑠𝑠𝑠𝑠 , −𝛻𝛻𝛻𝛻𝛻𝛻𝛻𝛻𝑉𝑉𝑉𝑉𝑠𝑠𝑠𝑠 , 𝒇𝒇𝒇𝒇𝑑𝑑𝑑𝑑 , 𝒇𝒇𝒇𝒇𝑐𝑐𝑐𝑐 はそれぞれ DEM 粒子の質量，体積，浮力，抗力，土粒子同士の接触力を

2.1 ISPH 法による流体解析

示す．接触力には，一般の DEM 解析に用いられる

水は非圧縮性流体として扱い，地盤の間隙率𝜀𝜀𝜀𝜀を介して，地表流と浸透流を統一的に記述した Darcy-

バネ・ダッシュポットモデルを採用しているが詳細は割愛する．流体との相互作用である抗力𝒇𝒇𝒇𝒇𝑑𝑑𝑑𝑑 に関し

Brinkman 型の統一方程式を採用し，これを ISPHて，比較的大きな土粒子に対しては前報と同様に法によって解く．

Darcy 流速を用いた抗力を算出し，流体の抵抗力の�𝑓𝑓𝑓𝑓
𝐶𝐶𝐶𝐶𝑟𝑟𝑟𝑟 (𝜀𝜀𝜀𝜀) 𝐷𝐷𝐷𝐷𝒗𝒗𝒗𝒗
1
�𝑓𝑓𝑓𝑓
= − 𝛻𝛻𝛻𝛻𝛻𝛻𝛻𝛻 + 𝒈𝒈𝒈𝒈 + 𝜈𝜈𝜈𝜈𝐸𝐸𝐸𝐸 (𝜀𝜀𝜀𝜀)𝛻𝛻𝛻𝛻 2 𝒗𝒗𝒗𝒗
𝜌𝜌𝜌𝜌𝑓𝑓𝑓𝑓𝜀𝜀𝜀𝜀 𝐷𝐷𝐷𝐷𝐷𝐷𝐷𝐷

反作用として計算する．実際の現象では，大きな土粒子同士の間隙に存在する小さな粒子は，その流失−𝑎𝑎𝑎𝑎(𝜀𝜀𝜀𝜀)𝜀𝜀𝜀𝜀𝒗𝒗𝒗𝒗𝑟𝑟𝑟𝑟 − 𝑏𝑏𝑏𝑏(𝜀𝜀𝜀𝜀)𝜀𝜀𝜀𝜀 2 |𝒗𝒗𝒗𝒗𝑟𝑟𝑟𝑟 |𝒗𝒗𝒗𝒗𝑟𝑟𝑟𝑟 (𝜀𝜀𝜀𝜀 < 0.8)

しやすさが想定されるため，実流速を用いた抗力を

ここで𝛻𝛻𝛻𝛻, 𝒈𝒈𝒈𝒈, 𝜌𝜌𝜌𝜌𝑓𝑓𝑓𝑓 , 𝒗𝒗𝒗𝒗𝑓𝑓𝑓𝑓 , 𝒗𝒗𝒗𝒗𝑟𝑟𝑟𝑟 , 𝐶𝐶𝐶𝐶𝑟𝑟𝑟𝑟 はそれぞれ，圧力，重力加した．土粒子の回転運動について，𝐼𝐼𝐼𝐼, 𝝎𝝎𝝎𝝎𝒔𝒔𝒔𝒔 , 𝒎𝒎𝒎𝒎𝑐𝑐𝑐𝑐 , 𝒎𝒎𝒎𝒎𝑟𝑟𝑟𝑟 は�−𝑐𝑐𝑐𝑐(𝜀𝜀𝜀𝜀)|𝒗𝒗𝒗𝒗𝑟𝑟𝑟𝑟 |𝒗𝒗𝒗𝒗𝑟𝑟𝑟𝑟(𝜀𝜀𝜀𝜀 ≥ 0.8)

適用して，簡易的な半解像型連成を実施する方針と速度，流体密度，流体速度，仮想質量係数，有効粘慣性モーメント，角速度，接線方向接触力によるト性係数を表す．ここで，𝜀𝜀𝜀𝜀は地盤の間隙率を示し，計

ルク，回転摩擦を示す．

算対象の SPH 粒子の影響範囲に入る DEM 粒子の総体積から数値的に求める．間隙率を用いて，ダルシ�𝑓𝑓𝑓𝑓 = 𝜀𝜀𝜀𝜀𝒗𝒗𝒗𝒗𝑓𝑓𝑓𝑓 及び地盤内の見かけの流体密度
ー流速 𝒗𝒗𝒗𝒗
𝜌𝜌𝜌𝜌̅𝑓𝑓𝑓𝑓 = 𝜀𝜀𝜀𝜀𝜌𝜌𝜌𝜌𝑓𝑓𝑓𝑓 を得る．また，𝑎𝑎𝑎𝑎, 𝑏𝑏𝑏𝑏, 𝑐𝑐𝑐𝑐, 𝜈𝜈𝜈𝜈𝐸𝐸𝐸𝐸 は抵抗係数，粘性16 − Academic Center for Computing and Media Studies, Kyoto University

3 非解像型連成による防波堤崩壊解析
Coarse gravel→ Resolved coupling

上述の非解像型連成に基づいて，ISPH-DEM に

Fine particles→ Un-resolved coupling

よる防波堤マウンドの浸透崩壊シミュレーションを実施した．図-1 に示す解析結果のように，非解像型連成の場合は，捨石マウンドの変形を表現したが，ケーソンを倒壊に至らしめる局所的な噴砂を再現[Central cross-section]

できなかった．そこで，Terzhagi の限界動水勾配に基づく破壊判定を導入し，破壊された DEM 粒子を流失したと仮定し，その粒子の質量を減少させた場
合，ケーソン近傍で噴砂が生じ，ケーソンブロック
が転倒，防波堤マウンドの崩壊に至る，実験[1]と同じ傾向を再現できた[2]（破壊判定を導入した計算につ:Velocity vector

いては現在投稿中）．
[A] 粒子質量一定[B] 質量補正あり

0s0s

図-2 ハイブリッド連成による水-土粒子連成解析

5 まとめと今後の課題
2.0 s

2.0 s

本研究では，ISPH 法による流体解析と DEM による地盤の変形解析を，非解像型，解像型，そして
4.0 s

4.0 s

それらハイブリッド法に基づく連成を行い，水-土混
相流シミュレータを開発した．現在，ハイブリッド

0.075

質量分布[g]1.10.0

流速分布[cm/s] 10.0

図-1 非解像型連成による防波堤浸透崩壊解析

連成の高精度化と精度検証に取り組んでおり，今後，これら連成手法の使い分けにより，地盤の内部侵食を効率的に計算できる手法を目指す．また，均質化

4 解像・非解像型ハイブリッド連成法等の平均化手法と組み合わせることで，内部侵食時の地盤の力学的特性を把握するマルチスケール

地盤の局所的な破壊を精緻に再現するには，間隙モデリングへの展開を目指す．

流を直接表現する解像型連成に基づく流体解析が有効だが，現状の手法では計算規模，精度の点から

謝辞：本研究はJSPS 科研費21J12338・20H02418，学際大困難を極める．そこで，間隙内の複雑な流れを表現

規模情報基盤共同利用・共同研究拠点（jh170016-NAH，する手法へと改良すべく，まずは SPH 法の従来のEX21604）の支援を受けた．また，九州大学情報基盤研究

勾配モデルを見直し，１次コンシステンシーを有す

開発センター研究用計算機システムの一般利用，京都大学

る勾配モデルを適用することで，従来手法で不可能

学術メディアセンターのスーパーコンピュータを利用しだった負圧計算と，Karman 渦列のような複雑な流た．ここに記して謝意を表する．

れを表現出来る手法へと改良した．しかし，上述の参考文献：

間隙内部の流れを直接表現する解像型連成解析は，[1] Kasama, K., Zen, K., Nakagawa, Y., and Furukawa, Z. :微視的なスケールの解析には有用であるものの，計Instability evaluation of rubble mound for breakwaters

算対象領域が大きくなるに連れて莫大な計算コス

subjected to tsunami-induced overflow and seepage flow.

トが生じることが予想される．そのため，大小の粒

Soils and Foundations, Vol.60 (6), pp.1532–1545, 2020.

径の礫地盤内の局所的に大きな流速を持つ流れを[2] 辻勲平，浅井光輝，笠間清伸：ISPH-DEMによるマウ

効率的に計算する「解像・非解像型ハイブリッド連

ンド変形・ケーソン滑動を考慮した防波堤浸透崩壊シ成手法」の開発を進めた．

ミュレーション，土木学会論文集A2（応用力学），
Vol.77，No.2, pp.I_105-I_116，2021.

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

超音波キャビテーションマルチバブル動力学の探索
山本卓也
東北大学大学院工学研究科金属フロンティア工学専攻

1 緒言

計算領域は一辺 400 µm からなる立方体領域であるとし、計算領域中央に直径 20 µm の気泡を配液体中に超音波を照射すると、音響キャビテー

置し、その気泡の周囲に正多面体の頂点位置に気ションが生じ、キャビテーション気泡が大きく膨泡を配置した。本計算では、正四面体、正六面体、張、圧縮するのを繰り返す。圧縮時には、気泡内正八面体の三種類の配置で、気液の物性は空気とぶが 100 気圧以上、5000 K 以上の温度になると水であるとした。中心気泡-周囲気泡間距離を 100,

され、熱分解に伴う化学反応が進行することで材

150 µm とし、音圧振幅を 0.5 – 0.9 atm.として計算料合成技術が開発されてきた[1]。また、このキャを行った。境界面で 20 kHz の超音波を照射し、気ビテーションによって、合成された無機化合物の泡振動を調査した。

形態が変化することが知られている。このような超音波処理プロセスはキャビテーション気泡の運

3 結果と考察

動や化学反応、音響流等が同時に発生することにより、複雑な現象となるため未解明な現象が数多

Figure 1 に様々な気泡配置の場合における中心

く見受けられる。特に、複数気泡存在する条件で気泡半径の時間変化を示す。音圧振幅が 0.7 atm.

は、クエンチングと呼ばれる反応が停止する現象の場合の結果のみを示す。
図中で Tetra は正四面体、[2]が知られており、複数存在する気泡の運動によHexa は正六面体、Okta は正八面体の配置であり、って反応が停止すると考えられているが、その気それぞれ周囲気泡数は、4, 8, 6 個である。また、R

泡運動までを詳細に解明できていない。

は中心気泡、周囲気泡間距離を表す。全ての条件

本研究では、マルチバブル系における音響キャビ

において、時間と共に気泡は膨張し、圧縮時に非テーション気泡の運動を解明するため、複数気泡存

線形な振動を行う。また、周囲気泡数が増加する

在する条件において超音波を照射し、数値解析によに伴い、気泡膨張時の最大半径が増加し、非線形ってその気泡運動を調査した。

振動時の振動振幅が増幅される。つまり、マルチ
バブル系においては、気泡の振動振幅が増幅され

2 数値解析手法

る。
気泡の非線形振動の圧縮強度は最大気泡径に依支配方程式として、圧縮性 Navier-Stokes 式、連存することが報告されている[7]。このため、気泡

続式、energy 式を用い、気液混相流モデルの一つ径が最大となる現象が重要となる。気泡径が最大である Volume of Fluid (VOF)法と連成して解いた。

になる時刻における気泡の時間変化を Fig. 2 に示VOF 法には代数的に界面再構築を行う方法[3]を利す。周囲気泡が先に圧縮し、中心気泡の圧縮が遅用 し 、 計 算 負 荷 低 減 の た め Adaptive Meshれる。この圧縮の遅れによって気泡振動と超音波

Refinement (AMR)と組み合わせた。表面張力モ

振動の位相差がより大きくなり、マルチバブル系

デルとして Continuum Surface Force (CSF)モデ

において最大気泡径がより大きくなったと考えら

ル[4]を利用し、計算精度向上のために表面張力項れる。

計算時には Laplacian フィルタ[5]を導入した。詳細な数値モデルに関しては既報[6]の通りである。

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

気泡振動振幅増幅を定量化するために、新しい変数を導入する。周囲気泡の影になることで音波伝播遅れが生じるので、カバー率と名付けた中心気泡[4] J. U. Brackbill et al., J. Comput. Phys. 100

に対する投影面積と中心気泡の表面積の比を用い(1992) 335-354.

て定量化する。カバー率 rcov は以下のように定義し[5] D. A. Hoang et al., Comput. Fluids 86 (2013)

た。

28-36.

(1)((4𝜋𝜋𝜋𝜋𝑅𝑅$%$$&'
𝑛𝑛𝑅𝑅$%$$&'
𝑟𝑟!"# =
=
(
4𝜋𝜋𝑅𝑅
𝑅𝑅(

[6] T. Yamamoto, S. V. Komarov, J. Appl. Phys.
128 (2020) 044702.

ここで、n は周囲気泡数、Rbubble は周囲気泡の半径である。以下の結果では、周囲気泡配置、周囲気泡数[7] W. Lauterborn, T. Kurz, Rep. Prog. Phys. 73(2010) 106501.

が異なる場合においても、全てカバー率のみから評価を行った。
Figure 3 に最大気泡径とカバー率の関係性を示す。Figure 3 (b)には単一気泡の場合との差を示す。
最大気泡径は音圧振幅と共に増加するが、カバー率にも比例して増加する。このことから、気泡配
置や気泡間距離に寄らず、カバー率のみから気泡
振動振幅増幅を定量化することができた。

4 結言
Figure 1 Time variation of equivalent bubble radius at
本研究では数値解析を利用してマルチバブル系

0.7 atm. in single- and multi-bubble systems.

における気泡振動を調査し、以下の知見を得ることができた。
・

マルチバブル系では気泡振動振幅が大きくなる。

・

周囲気泡の個数が大きくなり、周囲気泡との距離が短くなると、気泡の振動振幅が大きくなる。

・

周囲気泡の投影面積比であるカバー率を利用すると気泡振動増幅を整理でき、超音波振動
中の最大気泡径はカバー率に比例する。

Figure 2 Time variation of cavitation bubble shapeduring their contraction half-period under the pressure
amplitude of 0.9 atm. The distance between the center
bubble and surrounding bubbles is 100 µm, and the
bubble arrangement is regular hexahedron. Time is

謝辞

本研究制度(若手・女性研究者奨励枠)を活用さ

counted from the beginning of simulation.

せて頂いたことを、この場を借りて厚く御礼申し
上げる。
参考文献

[1] 日本音響学会、音響バブルとソノケミストリー（コロナ社、2012 年）[2] S. I. Hatanaka et al., Jpn. J. Appl. Phys. 40(2001) 3856-3860.
[3] H. Rusche, Ph. D. Thesis of Imperial Collageof Science, Technology and Medicine, London(2002).

Figure 3 Relationship between the maximum
equivalent radius of bubble and cover ratio atdifferent pressure amplitudes: (a) minimumpressure and (b) difference in the minimumpressures between the multi-bubble and thesingle-bubble systems.

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.21 No.1（ Dec. 2022）− 19ࡾ࣍‫਺ݩ‬஋ྲྀମղੳʹΑΔւ‫ྛ؛‬ͷ௡೾๷ࡂޮՌ‫ߏػ‬ͷղ໌໺ଜྯՂ ∗∗1౦๺େֶࡂ֐Պֶࠃࡍ‫ॴڀݚ‬ং࿦ࣗવੜଶ‫ܥ‬Λ‫ͨ͠༻׆‬αεςφϒϧͳ๷ࡂɾ‫ࡂݮ‬ઃඋͱͯ͠ɼ৿ྛͷଘࡏʹ‫ظ‬଴͕‫ͤد‬ΒΕ͍ͯΔɽઋ
୆๷ࡂ࿮૊ʢࠃ࿈๷ࡂੈքձٞʣͰ΋ɼͦͷ‫ޮࡂݮ‬͞Ε͍ͨɽ͜Ε͸ɼຊ‫ط͕ڀݚ‬ใ [1] ಉ༷ʹɼϚϧνεέʔϧϞσϧʹΑΔւ‫ྲྀ಺ྛ؛‬ΕͷΑΓ߹ཧతɾਖ਼֬ͳγϛϡϨʔγϣϯΛ࠷ऴ໨ඪͱͯ͠‫͍ͯ͛ܝ‬Δ͜ͱʹΑΔɽՌʹ͍ͭͯ‫͞ٴݴ‬ΕΔͳͲɼͦͷ‫ظ‬଴͸ࠃࡍతʹ΋
ߴ·͍ͬͯΔɽ๷ࡂޮՌ΍‫ػ‬ೳ‫ݶ‬քͷਖ਼֬ͳΞηεϝϯτʹ͸ɼਁਫγϛϡϨʔγϣϯΛߦ͏͜ͱ͕༗
ޮ͕ͩɼ‫ط‬ԟͷਁਫγϛϡϨʔγϣϯͷଟ͘͸ɼա౓ʹ୯७Խͨ͠ւ‫ྛ؛‬Ϟσϧͷར༻Λओྲྀͱ͓ͯ͠
ΓաখɾաେͳੑೳධՁͷՄೳੑ͕͋ͬͨɽෳࡶͳ
थ໦‫ܗ‬ঢ়ΛϞσϧԽ͠ɼथ໦ͷ‫ݮ‬ਰޮՌΛධՁ͢ΔͨΊʹ͸ɼࡾ࣍‫਺ݩ‬஋ྲྀମղੳͷ࣮ࢪ͕ඞਢͰ͋Δɽຊ‫Ͱڀݚ‬͸ɼ͜Ε·ͰաখɾաେͳධՁ͕͞Ε͍ͯͨथ໦‫܈‬ͷ఍߅ޮՌΛΑΓ߹ཧతʹධՁ͢Δ͜ͱΛ໨తʹɼա౓ʹ୯७Խ͞Ε͍ͯͨւ‫ྛ؛‬ΛৄࡉʹϞσϧԽͨ͠ࡾ࣍‫਺ݩ‬஋ྲྀΕࢼ‫ݧ‬Λ࣮ࢪ͠ɼथ໦ͷෳࡶͳ‫ܗ‬ঢ়ʹ‫ى‬Ҽͨ͠఍߅ಛੑ͕ϚΫϩϨϕϧͰͷ‫ݮ‬ਰʹͲΕ΄ͲޮՌΛ‫͢΅ٴ‬ͷ͔Λௐࠪ͢Δɽ۩ମྲྀΕͷࢧ഑ํఔࣜ

2.1

࣍ʹࣔ͞ΕΔࡾ࣍‫ ݩ‬Navier-Stokes ͷӡಈํఔࣜ‫ͼٴ‬࿈ଓͷࣜΛ҆ఆԽ༗‫ݶ‬ཁૉ๏ʢSUPG/PSPGʣ

[2, 3] ʹΑΓ཭ࢄԽ͠ɼղੳ͢Δ͜ͱͰࡾ࣍‫਺ݩ‬஋ྲྀମղੳΛ࣮ࢪͨ͠ɽ

)
(
∂u

+ u · ∇u = ∇ · σf + b 
ρf
∂t


∇·u=0

in

Yf ,

(1)
͜͜Ͱɼu ͸ྲྀ଎ϕΫτϧΛ͓ࣔͯ͠ΓɼԠྗσf ͸(( ∂u )T )∂uѹྗ P ‫ͤͼٴ‬ΜஅԠྗϕΫτϧ τ = µf ∂x + ∂xͱͱ΋ʹɽ࣍ͷࣜͰఆΊΒΕΔɽతʹ͸‫ط‬ԟͷ‫ڀݚ‬ͷଟ͕͘࠾༻͍ͯͨ͠ԁபͷΑ͏ͳγϯϓϧͳϞσϧͱൺֱ࣮‫ݧ‬Λߦ͏͜ͱͰ (1) ৄ

σf = −P I + τ ,ࡉͳथ໦Ϟσϧͳ͠ʹ͸աখɾաେͳධՁ͕‫͜ى‬Δ

͜͜Ͱ I ͸ೋ֊ͷ߃౳ςϯιϧͰ͋Δɽྲྀମ෺ੑ஋

ྲྀΕɼ(2) ‫ط‬ԟͷ୯७ϞσϧͰ߹ཧతͳධՁ͕Մೳͳ

Ͱ͋Δີ౓ ρf ‫ͼٴ‬೪ੑ܎਺ µf ʹ͍ͭͯ͸ͦΕͧྲྀΕɼͷೋͭΛ໌֬Խ͢Δ͜ͱΛ໨ࢦ͢ɽΕ 25C◦ Թ౓ԼͰͷਫͷ෺ੑ஋Λ༩͑Δɽͳ͓ɼ࣌

(2)ؒํ޲ͷ཭ࢄԽʹ͸ΫϥϯΫχίϧιϯ๏Λ༻͍ͨɽ2਺஋ྲྀΕࢼ‫ݧ‬ࣗ༝ද໘ͷัଊʹ͸ Phase-Field ๏ [4, 5] Λ༻͍͍ͯΔɽຊ‫Ͱڀݚ‬͸ۣ‫ܗ‬ਫ࿏Ϟσϧ಺ʹઃஔͨ͠ԁபϞσϧɼथ໦Ϟσϧͷ૒ํʹରͯ͠ࡾ࣍‫਺ݩ‬஋ྲྀମղੳ
Λ࣮ࢪ͠ɼ྆ऀͷղੳ݁Ռ͔ΒͦΕͧΕ‫ݮ‬ਰޮՌʢྫ

3݁Ռ͑͹ѹྗɾਫਂͷ௿‫޲܏ݮ‬ʣΛଌఆ͢Δ͜ͱͰ֤Ϟং࿦Ͱड़΂ͨ (1) ‫( ͼٴ‬2) ͷ֤ྲྀΕΛ໌Β͔ʹ͢σϧͷ๷ࡂޮՌΛධՁ͢Δɼ“਺஋ྲྀΕࢼ‫ ”ݧ‬ͷखଓ

ΔͨΊɼෳࡶͳथ໦Ϟσϧɾ؆ૉͳԁபϞσϧͷ͖྆Λ౿Ήɽ਺஋ྲྀΕࢼ‫ݧ‬ͷ֓೦΍֤छ‫ڥ‬ք৚݅ͷৄํΛ༻͍ͨࡾ࣍‫਺ݩ‬஋ྲྀମղੳΛߦͬͨɽΑΓ۩ମࡉʹ͍ͭͯ͸ಉ༷ͷઃఆΛߦ͍ͬͯΔ‫ط‬ใ [1] Λࢀরతʹ͸ɼ(i) थ໦Ϟσϧ͕୯ಠͰ഑ஔ͞Εͨ৔߹ɼ(ii)20 − Academic Center for Computing and Media Studies, Kyoto UniversityԁபϞσϧ͕୯ಠͰ഑ஔ͞Εͨ৔߹ɼ(iii) थ໦Ϟσaϧ͕֨ࢠঢ়ʹ‫܈‬഑ஔ͞Εͨ৔߹ɼ(iv) ԁபϞσϧ͕֨ࢠঢ়ʹ‫܈‬഑ஔ͞Εͨ৔߹ɼͷ̐ύλʔϯͷղੳΛ࣮ࢪͨ͠ɽथ໦ϞσϧɼԁபϞσϧͷ૒ํ͸ਤ 1 aʹࣔ͢ɽҰ࿈ͷղੳ݁Ռͷ͏ͪɼ࣍ࣜʹࣔࣜ͢ΑΓࢉఆͨ͠ྲྀମͷ୯Ґମੵ౰ͨΓӡಈΤωϧΪʔʢ=ಈѹʣΛ
TreeՄࢹԽͨ͠ྫΛਤ 1 b ʹࣔ͢ɽ

Dynamic pressure =1
ρf u · u [Pa].
2

(3)

b

0.0

10

Cylinder20

Dynamic pressure [Pa]
30
40
50

60

75+

ਤ 1 ͔Β෼͔ΔΑ͏ʹɼෳࡶͳ‫ܗ‬ঢ়Λ༗͢Δथ໦ϞσϧͰ͸ɼ؆қͳԁபϞσϧͱൺֱͯ͠ɼྲྀΕͷ‫ݮ‬ਰޮՌ͕޿ൣ͔ͭྲྀԼํ޲ʹ௕͘‫ܧ‬ଓ͢Δɽਫਂɾzྲྀ଎Λมߋͨ͠‫ ܭ‬24 έʔεͷ਺஋ྲྀΕࢼ‫݁ݧ‬Ռ͔

d-1
x

Tree

x

CylinderΒ·ͱΊͨͱ͜Ζɼ࣍ͷΑ͏ͳ஌‫ݟ‬ΛಘΔ͜ͱ͕Ͱ͖ͨɽz• थ໦Ϟσϧͷෳࡶ‫ܗ‬ঢ়͕ਫ຅͍ͯ͠Δ৔߹Ͱd-1

͋ͬͯ΋؆қͳԁபϞσϧͱಉఔ౓ͷ‫ݮ‬ਰޮՌͱͳΔྲྀΕ৚݅ͷଘࡏ͕֬ೝ͞Εͨɽ• ԁபϞσϧͱथ໦Ϟσϧͷ‫ݮ‬ਰޮՌ͸ɼෳࡶ‫ܗ‬ঢ়͕‫׬‬શʹਫ຅͠ɼߴ͍ྲྀ଎஋͕༩͑ΒΕͨέʔεʹ͓͍ͯಛʹ͔͍཭͕ܹ͘͠ͳͬͨɽ• ैདྷͷ߅ྗࢉఆʹ༻͍ΒΕ͍ͯͨ౤Ө໘ੵͷΑ͏ͳೋ࣍‫ݩ‬తࢦඪͰ͸ͳ͘ɼ຅ਫ෦෼ͷମੵ઎༗཰ͳͲͷࡾ࣍‫ݩ‬తࢦඪΛ༻͍Δํ͕ਖ਼֬ͳ‫ݮ‬ਰޮՌධՁ͕ߦ͑Δɽ

4

݁࿦ຊ‫Ͱڀݚ‬͸ɼ(1) ৄࡉͳथ໦Ϟσϧͳ͠ʹ͸աখɾաେͳධՁ͕‫͜ى‬ΔྲྀΕɼ(2) ‫ط‬ԟͷ୯७ϞσϧͰ߹ཧతͳධՁ͕ՄೳͳྲྀΕɼͷೋͭΛ໌֬Խ͢Δ͜ͱΛ໨ࢦ͠ɼྲྀ଎΍ਫਂΛมߋͨ͠ෳ਺ͷ਺஋ྲྀΕࢼ‫ݧ‬ΛɼԁபϞσϧͱथ໦Ϟσϧͷ૒ํʹର࣮ͯ͠ࢪͨ͠ɽҰ࿈ͷ݁Ռ͔Βɼ(1) ‫( ͼٴ‬2) ͷྲྀΕͷ৚݅Λఆੑతʹ·ͱΊΔ͜ͱ͕Ͱ͖ͨɽैདྷख๏Ͱ༻͍
ΒΕ͖ͯͨ౤Ө໘ੵͷΑ͏ͳೋ࣍‫ݩ‬తͳࢦඪΛආ͚ɼෳࡶ‫ܗ‬ঢ়ͷମੵ઎༗཰ͳͲΛར༻ͨ͠‫ݮ‬ਰੑೳධՁ
͕ॏཁͰ͋Δ͜ͱ͕໌Β͔ʹͳͬͨɽΑΓଟ͘ͷ਺஋ྲྀΕࢼ‫ݧ‬Λߦ͍ɼࡾ࣍‫ݩ‬తͳࢦඪΛ༻͍ͨ‫ݮ‬ਰޮਤ 1: a ਺஋ྲྀΕࢼ‫͍ͨ༻ʹݧ‬थ໦ϞσϧɾԁபϞσϧɼb ֤Ϟσϧ‫ޙ‬എͰͷಈѹ (3) ෼෍ͷ༷ࢠʢx-zฏ໘ɼྲྀԼํ޲͸ӈखଆʣࢀߟจ‫ݙ‬[1] Nomura, R., Takase, S., Moriguchi, S., Terada, K.,
and LeVeque, R. J. “Multiscale evaluation methodof the drag eﬀect on shallow water ﬂow throughcoastal forests based on 3D numerical simulations”,
Int J Numer Meth Fl, 94(1), 32-58, 2022.
[2] Brooks, A. N. and Hughes, T. J., “Streamline
upwind/Petrov-Galerkin formulations for convection dominated ﬂows with particular emphasis on
the incompressible Navier-Stokes equations”, Comput Method Appl M, 32(1-3), 199-259, 1982.
[3] Bova, S. W. and Carey, G. F., “A symmetric formulation and SUPG scheme for the shallow-waterequations”, Adv Water Resour, 19(3), 123-131,
1996.
[4] Chiu, P. H. and Lin, Y. T., “A conservative phaseﬁeld method for solving incompressible two-phaseﬂows”, J Comput Phys, 230(1), 185-204, 2011.
[5] Takada, N., Matsumoto, J. and Matsumoto, S.,“Phase-ﬁeld model-based simulation of motions ofa two-phase ﬂuid on solid surface”, J Comput Sci
Tech, 7(2), 322-337, 2013.

ՌੑೳΛߦ͏͜ͱͰɼԁபϞσϧͱथ໦Ϟσϧͷ‫ݮ‬ਰੑೳʹ‫ݦ‬ஶͳ͔͍཭͕‫ݱ‬ΕΔ͖͍͠஋Λ໌Β͔ʹ͍ͯ͘͜͠ͱΛࠓ‫ޙ‬ͷ՝୊ͱ͍ͯ͠Δɽ

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

代数幾何的手法とランダム行列理論を用いた再帰的ニューラル
ネットワークによる時系列スパースモデリング中野 直人 ∗

背景

アトラクタ

観測可能
40

20

人工ニューラルネットワークは脳神経科学のネッ

10

トワークを模擬するために設計され，現代の機械学

20
10

習技術の重要な基盤となっている．深層学習で用いいるが，人工ニューロンの相互的な情報処理を真似た

x(t)

°20

全体を復元

られる前向きニューラネットワークがよく知られて構造の再帰的ニューラルネットワークは時系列解析

40
30

z(t)

1

明治大学大学院先端数理科学研究科

°10

0
0

10

0

°10

観測可能

20 °20

観測

10
0

10

20
t

30

20

30

20

30

t

0

°20
0

10

20

推定可能

°10

20

推定可能

y(t)

∗

0

10
t

と相性が良いとされている．その中でも，人工ニュー図 1: Lorenz 方程式の数値解．左上の (x, y, z) のを用いたモデルは，集団的カオスダイナミクスやノよってアトラクタを再構成して (y, z) を推定する．

ロン同士の結合を学習しないランダムネットワーク
イズによるカオス抑制といった豊かな動力学的特性を持つ．さらに共通入力同期と呼ばれ，過去の履歴

3 変数の系から右下の x 変数のみを観測し，RC に知られている [2]．本研究では RC の普遍的な性質に関わらず同じ入力の下で内部ノード状態の任意のを導くためにランダムネットワークを経路積分の方

これらの特徴からランダムネットワークは，学習

トワークを確率分布のパラメータ座標ごとにそれぞ軌道が収束する現象が知られている．

の負荷が低いこともあいまって，リアルタイムでの法を用いて解析した [3]. サイズ 500 のランダムネッ
れ 100 アンサンブル構成し，RC による時系列推定

情報処理へ応用されており，そのような枠組みの 1

タスクを学習する．性能評価はテストデータに対す

る．ここ数十年，計算性能の向上や動作原理の解明いる．相空間において確率分布のパラメータを変えつとして reservoir computing (RC) が用いられていのために RC は多くの拡張がなされ，情報処理能力

もさまざま尺度で評価されてきた．しかし RC の計

算性能を最適化するための設計指針やそもそも RCがうまく機能する理由について，今のところ十分な理解は得られておらず，RC 自体の普遍的な記述が必要とされている．

2

手法

RC のようにノード間の相互作用によって情報処理を行うシステムに対しては，同じく多体間相互作用に対する統計物理による記述が重要であることが22 − Academic Center for Computing and Media Studies, Kyoto University

る平均二乗誤差 (MSE) の中央値と四分位範囲を用て多アンサンブルの計算を行うため，スーパーコン
ピュータによる並列計算が威力を発揮する．

3

データ

ここではカオス時系列の代表例である Lorenz 方

程式 [1] の時系列データを用いる．Lorenz 方程式は

(x, y, z) の 3 変数の方程式で記述されるが（図 1 左上），観測可能な変数は x だけとする（図 1 左下）．

観測が部分的な状況で，観測データから力学全体をRC で再構成し，その他の変数（図 1 右）を推定するタスクを行い，RC の性能を評価する．

1/J

2.0

0

1.5

−1
−2
−3

1.0

−4

0.5

−5

0.0
0.0

0.5

1.0

1.5

2.0

−6

1/J

J0 /J
2.0

0

1.5

−1
−2
−3

1.0

−4

0.5

−5

0.0
0.0

0.5

1.0

1.5

2.0

−6

J0 /J
図 2: 十分時間発展させたガウシンアンネットワー

クのノード値に対する (上) ノード間平均と (下) 分

散のアンサンブル平均．

4

結果

図 2 は十分時間発展させたガウシンアンネット

ワークのノード値に対する (上) ノード間平均と (下)分散のアンサンブル平均の相図である．各パラメー
タ座標で N(J0 /N, J 2 /N ) に従うガウス分布に従う

ランダムネットワークを用いた．横軸は J0 /J, 縦軸

は 1/J である．これにより，ランダムネットワーク

は，ノード平均が 0, 分散が 0 の秩序相，平均が正，

分散が 0 の偏極秩序相，平均が 0，分散が正のカオ
ス相に分かれることが示された．

さらに，時系列タスクの性能に対しては相転移が確認され，相図で現れるネットワークの性質が大きく異なる相分離曲線の近傍で最適な学習性能が得られることがわかった．これを他の変数の推定タスク
や確率分布を対象とすることで（未掲），RC のカオ

スの縁に対する新たな知見を得ることができた．

図 3: ガウシアンネットワークを用いた RC による

時系列推定タスクのテストデータに対する平均二乗
誤差の (上) アンサンブル中央値と (下) アンサンブ

ル四分位範囲．

謝辞

本研究は京都大学スーパーコンピュータ共同研究制度(若手・女性奨励枠) の活用によるものである．また，本研
究は春名純一氏, 兎子尾理貴氏 (京都大学理学研究科) との共同研究を含み，科学研究費基盤研究 (C)(20K03747) の助成を受けた．ここに関係各位に対して感謝の意を表する．

参考文献[1] E. N. Lorenz, Deterministic Nonperiodic Flow, J.
Atmos. Sci., 20, (1963).
[2] H. Sompolinsky, A. Crisanti, and H. J. Sommers,
Chaos in Random Neural Networks, Phys. Rev.
Lett. 61, 259–262 (1988).
[3] J. Haruna, R. Toshio and N. Nakano, Universality
in reservoir computing and higher-order statistics,
arXiv:2112.01886 (2021)

京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.21 No.1（ Dec. 2022）− 23౷‫ܭ‬෺ཧֶ໛‫ͱܕ‬ͷରԠΛར༻ͨ͠‫ػ‬ցֶशख๏ͷੑೳධՁय़໊७Ұ ∗∗‫౎ژ‬େֶେֶӃཧֶ‫ڀݚ‬Պ෺ཧֶɾӉ஦෺ཧֶઐ߈

ຊߘ͸த໺௚ਓࢯ (໌࣏େֶେֶӃ ઌ୺਺ཧՊֶ
‫ڀݚ‬Պ ಛ೚।‫ڭ‬त) ͓Αͼీࢠඌཧ‫౎ژ( ࢯو‬େֶେֶӃཧֶ‫ڀݚ‬Պ) ͱͷ‫ڞ‬ಉ‫[ ڀݚ‬1] ʹ‫ͮ͘ج‬΋ͷͰ
͋Δɽ

1

͸͡Ίʹ
ਓ޻χϡʔϥϧωοτϫʔΫ͸‫ݱ‬୅ͷ‫ػ‬ցֶशཧ࿦ʹ͓͚Δຊ࣭తͳ‫͋Ͱૅج‬Δɽಛʹ࠶‫ؼ‬తχϡʔϥϧωοτϫʔΫ (RNN) ͸๛͔ͳಈྗֶతಛੑ [2]Λࣔ͠ɼϦΞϧλΠϜͰͷ৘ใॲཧʹԠ༻͞Ε͍ͯ
ΔɽͦͷΑ͏ͳ৘ใॲཧͷ࿮૊ΈͷҰͭͱͯ͠ɼϦβόʔ‫( ػࢉܭ‬RC) ͕‫͛ڍ‬ΒΕΔɽRC ͸ϩʔϨϯπΧΦεͳͲͷ࣌‫ྻܥ‬༧ଌʹԠ༻͞Ε੒ޭΛऩΊ͍ͯΔɽҰํͰɼϦβόʔ‫ʹػࢉܭ‬ΑΓ࣌‫ྻܥ‬༧ଌ͕͏·͍͘͘ཧ༝΍‫ࢉܭ‬ੑೳΛ࠷దԽ͢ΔͨΊͷઃ‫ࢦܭ‬਑ʹ͍ͭͯɼ౷Ұతͳཧղ͸ಘΒΕ͍ͯͳ͍ɽզʑ͸͜ΕΒͷղ໌Λ໨ࢦͯ͠ɼ‫ܦ‬࿏ੵ෼දࣔͰͷੜ੒ؔ਺ͷํ๏ [3–5] Λ༻͍ɼϊʔυ਺͕े෼େ͖͍৔߹
ͷ RC ͷಈྗֶʹ͍ͭͯཧ࿦తʹղੳͨ͠ɽ͜ͷ݁
ՌɼRC ͷಈྗֶ͸݁߹ఆ਺ͷ෼෍ؔ਺ʹԠͯ͡෼ྨ͞ΕΔ͜ͱ͕൑໌ͨ͠ɽ͜Ε͸ɼҟͳΔ෼෍ؔ਺Λ͍࣋ͬͯͯ΋ɼಉ͡ྨ (“ීวྨ”) ʹଐ͢ΔωοτϫʔΫ͸ϊʔυ਺͕େ͖͍‫͍͓ͯʹݶۃ‬ɼಉ͡ಈྗ
ֶΛࣔ͢͜ͱΛҙຯ͢Δɽզʑ͸͜ͷීวྨͷ‫͢ূݕ‬ΔͨΊʹɼ࣌‫ྻܥ‬༧ଌλεΫͷେ‫ن‬໛‫ࢉܭ‬Λߦͬͨɽ۩ମతʹ͸ɼීวྨͷٞ࿦͔Βɼ݁߹ߦྻ Jij Λਖ਼‫ن‬෼෍ͱΨϯϚ෼෍͔ΒαϯϓϦϯάͨ͠ͱ͖ɼϦβόʔ‫ػࢉܭ‬͸ҟͳΔಈྗֶΛࣔ͢͜ͱ͕༧‫͞ݴ‬ΕΔɽ͜ΕΒͦΕͧΕͷ৔߹ʹɼ‫ࢉܭ‬ੑೳ͸ͲͷΑ͏ʹมΘΔ͔ɼ͓Αͼɼ‫ࢉܭ‬ੑೳ͕޲্͢Δύϥϝʔλͷબ୒ࢦ਑ΛಘΔͨΊɼϩʔϨϯπΧΦεͷ࣌‫ྻܥ‬༧ଌΛࢼΈͨɽ

24 − Academic Center for Computing and Media Studies, Kyoto University2ηοτΞοϓ
RC Ͱ͸ҎԼͷ (཭ࢄԽ͞Εͨ) ࣌ؒൃలํఔࣜΛߟ͑Δɽ

ri (t + 1) = −(1 − α)ri (t)


Nin
N
∑
∑
+ α
Jij tanh(ri ) +Wijin xj + bj 
i=j

(1)

j=1

ri (t) ͸ N ੒෼࣮਺ϕΫτϧɼα ͸ֶश཰ͱ‫ݺ‬͹ΕΔ࣮਺ɼJij ͸ N × N ߦྻɼbi ͸όΠΞεͱ‫ݺ‬͹Ε

Δ N ੒෼࣮ϕΫτϧɼW in ͸ N × Nin ߦྻɼxi (t)͸ೖྗ࣌‫ Ͱྻܥ‬Nin ੒෼࣮ϕΫτϧͰ͋Δɽ

۩ମతʹ͸ϩʔϨϯπΧΦε࣌‫ྻܥ‬ͷ X ࠲ඪͷೖྗ࣌‫͠ͱྻܥ‬ɼಉ࣌ࠁͷ Y ࠲ඪΛ༧ଌ͢Δ͜ͱΛࢼΈͨɽ͜͜Ͱ͸ϩʔϨϯπΧΦεͷ࣌ؒൃలํఔࣜ ′x (t) = −σ(x + y),

y ′ (t) = x(ρ − z) − y, ′z (t) = xy − βzʹͯॳ‫ظ‬஋Λ x(0) = y(0) = z(0) = 1 ͱͯ࣌ؒ͠ൃ
లͤͨ࣌͞‫ྻܥ‬Λѻͬͨɽύϥϝʔλ͸ σ = 10, ρ =28, β = 8/3 ͱͨ͠ɽ۩ମతʹ͸ ri (t = 0) Λ [0, 1] ͷҰ༷෼෍͔ΒϥϯμϜʹ༩͑ɼͦͷॳ‫ظ‬஋Λ༻͍ͯࣜ (1) Λ t = 0
͔Β t = 200000 ·Ͱ‫͠ࢉܭ‬ɼri (t) ͷ࣌‫ྻܥ‬Λಘͨɽ͜ͷͱ͖݁߹ߦྻ Jij ͸༩͑ΒΕͨ J0 ͱ J ʹର͠ɼฏ‫ ͕ۉ‬J0 /N ɼ෼ࢄ͕ J 2 /N ʹͳΔΑ͏ʹύϥϝʔλΛௐઅͨ͠ਖ਼‫ن‬෼෍ͱΨϯϚ෼෍͔ΒαϯϓϦϯ
άͨ͠ɽ͜ͷ্ͰɼϩʔϨϯπΧΦεͷ Y ࠲ඪͷ∑༧ଌ஋Λ ŷ(t) =j (Wout )j rj (t) ͱ͍͏ઢ‫ܗ‬ม‫׵‬ͰఆΊΔɽri (t) ͷ t = 5000 ͔Β t = 100000 ·Ͱͷ࣌‫ྻܥ‬Λֶशσʔλͱ͠ɼWout Λฏ‫ۉ‬ೋ৐‫ࠩޡ‬∑100000
M SE := t=5000 (y(t) − ŷ(t))2 Λ࠷খʹ͢Δ΋ͷͱͯ͠બͿɽ͜͏ͯ͠ಘΒΕͨ Wout Λ༻͍ͯ ŷ(t) (t ∈[105000, 200000]) ΛϩʔϨϯπΧΦεͷ Y ࠲ඪͷ༧ଌ஋ͱͨ͠ɽ
͜͏ͯ͠ಘΒΕֶͨशσʔλͱ༧ଌσʔλͷฏ‫ۉ‬ೋ৐‫ࠩޡ‬ΛɼJ0 /J ͱ 1/J ΛͦΕͧΕ [0.1, 2.0] ·Ͱ

0.02 ࠁΈͰมԽͤ͞ɼΧϥʔϚοϓΛ࡞੒ͨ͠ɽ·ͨ૬ਤͷҰ఺ (J0 /J, 1/J) ʹରͯ͠ɼJ0 ͱ J Λ‫ݻ‬ఆ͠ ri (t = 0) ʹ͍ͭͯ 100 Ξϯαϯϒϧ‫ࢉܭ‬Λߦ্ͬͨͰɼͦ͜ͰಘΒΕֶͨशσʔλͱ༧ଌσʔλͷฏ‫ۉ‬ೋ৐‫ ࠩޡ‬100 ‫ͭͣݸ‬ͷฏ‫( ۉ‬Ξϯαϯϒϧฏ‫ )ۉ‬ΛΧϥʔϚοϓʹදࣔͨ͠ɽ

ਤ 2: ݁߹ఆ਺ΛΨϯϚ෼෍͔Βαϯϓϧͨ͠ͱ͖ͳ͓ϊʔυ਺ N ͸ 500ɼֶश཰ α ͸ 0.2 ͱ͠ɼWin͸ [−1, 1] ͷҰ༷෼෍͔Βαϯϓϧ͠ɼbj ͸ 0 ͱͨ͠ɽ3

ͷֶश‫ࠩޡ‬ͷΧϥʔϚοϓ (N = 500, 100 Ξϯαϯϒϧ)݁Ռ͓Αͼ݁࿦
ࢴ໘ͷ౎߹্ɼֶश‫ࠩޡ‬ͷΧϥʔϚοϓͷΈΛ঺

հ͢Δɽਤ 1 ͕ਖ਼‫ن‬෼෍ͷ৔߹ɼਤ 2 ͕ΨϯϚ෼෍ͷ৔߹Ͱ͋Δɽಉ࣌ʹɼ༧උ‫Ͱࢉܭ‬ఆৗঢ়ଶͷϦβόʔঢ়ଶͷ૬ਤʹ͍ͭͯ΋‫͓ͯ͠ࢉܭ‬Γɼಛʹ‫ऀޙ‬ͷ৔߹Λਤ 3 ʹࣔ͢ɽ(౷‫ֶྗܭ‬తͳྨਪͰ͸Ϧβόʔঢ়ଶ ri (t) ͸εϐϯʹରԠ͠ɼ͜ͷ૬ਤ͸ͦͷฏ‫ͱۉ‬෼ࢄΛ΋ͱʹඳ͍ͨ΋ͷͰ͋Δɽ)ֶश‫ࠩޡ‬ͷΧϥʔϚοϓͱ౷‫ܭ‬෺ཧֶ໛‫ܕ‬ͷ૬ਤΛൺֱ͢Δͱɼ͜ΕΒ͕ؔ࿈͍ͯ͠Δ͜ͱ͕ࣔࠦ͞
ΕΔɽಛʹ౷‫ܭ‬෺ཧ໛‫ͯʹܕ‬૬సҠ͕‫͜ى‬Δύϥϝʔλ෇ۙ (֤৭ͷྖҬͷ‫ڥ‬ք) Ͱɼֶश‫͕ࠩޡ‬গͳ͘ͳΔͱ͍͏ࣔࠦΛಘͨɽࠓ‫ޙ‬͸ RC ͱ౷‫ܭ‬෺ཧֶ໛‫ܕ‬ͱͷؔ࿈ʹ͍ͭͯɼ‫ܦ‬࿏ੵ෼දࣔͰͷੜ੒ؔ਺ͷํ๏͔Β୳͍ͬͯͭ͘΋ΓͰ͋Δɽਤ 3: ݁߹ఆ਺ΛΨϯϚ෼෍͔Βαϯϓϧͨ͠ͱ͖ͷఆৗϦβόʔঢ়ଶ૬ਤ

ࢀߟจ‫ݙ‬[1] J. Haruna, R. Toshio, and N. Nakano, arXiv
preprint arXiv:2112.01886 (2021).
[2] H. Sompolinsky, A. Crisanti, and H. J. Sommers, Phys. Rev. Lett. 61, 259 (1988).
[3] P. C. Martin, E. D. Siggia, and H. A. Rose,Phys.
Rev. A 8, 423 (1973).
[4] C. De Dominicis and L. Peliti, Phys. Rev. B 18,
353 (1978).
[5] H. K. Janssen, Zeitschrift fr Physik B Condensed Matter 23, 377 (1976)ਤ 1: ݁߹ఆ਺Λਖ਼‫ن‬෼෍͔Βαϯϓϧͨ͠ͱ͖ͷֶश‫ࠩޡ‬ͷΧϥʔϚοϓ (N = 500, 100 Ξϯαϯϒϧ)京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.21 No.1（ Dec. 2022）− 25‫˿ٿߣ܊‬ƪ˾̅˖‫ػ‬ʍՑҗӌࡌ˴˙˼̅˂઺π ٩‫∗ ڭ‬∗୼֮Ҙ๎੝ӌ ӌࡱٰ֖ϙ ຌ૾࣮൙‫ۑ‬ӌ೼ฉච‫ۮ‬ʎࢬັ Ԇ߅ (ງ९੝ӌ‫ݗؼ‬ӌ೼) , ঀ෼ ոຒ߅ (φׄ੝ӌ‫ؼ‬ϼԛ๽ٰ֖Ѡ ) ʇʍִடٰ֖ʊՂʄ

ʪʝɣʱߪɸʇɶʅʡ, ɼʫɫٔଜ໼଺ʆɡʪ࣪܏
ʊʎɲʍ࠱ʍӌࡌൣ൥ɫอ‫ۄ‬ʆɡʪ.

ɮʡʍʆɡʪ.

1ʎɷʠʊ.
‫ׯ‬௻, Ցҗӌࡌʍφ࠱ʆ, ˼ʽƪ˾̅˚˟˕˚́ƪˁʱʡʂɾ˼ˈƪˢƪ˅̅˦˷ƪ˘ʵ̅˂ [1, 2] ɫߢ‫ػ‬໑ʣ˼˵˩ˠ˧޼ॐʉʈʍ฽ਢʊɩɣʅอ‫ۄ‬ʆ

3

˴˙˽ʍঙଜʇٗѢ.

3.1

ঙଜ.

N ߣ‫ٿ‬ʍ˿ƪ˾̅˖ൣପ߲ (Lorenz96 system) [6]ɡʪɲʇɫ൙‫ܘ‬ɴʫʅɣʪ. ѼƧʎɲʍՑҗӌࡌࠬdxk= xk−1 (xk+1 −xk−2 )−xk +f, for k = 1, . . . , N,
dt൥ʱຌ੄ʍʺ˟˽ˀƪഷॐʊ଼๑ɶ, ߢԨౙୈʍ฽

ʱ۵ɧʪ. ɾɿɶ, xN +k = xk (k = 1, . . . , N ) ʇɸ৳உʱ२‫ۃ‬ɴɺʅɣʪ [3, 4].

ʪ. ɲɲʆʎ೜Μଜߣ‫( ٿ‬९ʍ˼ʴ˩ˠ˧޼ॐ) ɫ 2

φൣʆ, ຌ੄ອӌʱ޳ʠʇɶɾߢ‫ػ‬໑˙ƪˑʎ, ࢉʆɡʪ (N, f ) = (8, 6) ʍ࣪܏ʱ۵ɧʪ. ɲʍʴ˚˻Վʍ‫ܿڷ‬ɫ޼ॐԪॐ଺ʊਕљɸʪʽʼˋ१ʊʧʩ૫

ˁˑʎ, ೜Μଜߣ‫ٿ‬ɫ 1 ʍࡀՎՠனʇ೜Μଜߣ‫ٿ‬ɫ 2ߢԨ฽ਢɫʆɬʉɣɲʇɫઢʨʫʅɣʪ. ɲʍɾʠ,ʍࡀՎՠனʍɣɹʫʡɫᘰ෋ʊඨʠ‫ܦ‬ʝʫʅɣʪ˪

ӌࡌʊʧʩமɾߢԨౙୈ˴˙˽ʍ҈উʊʎອӌ‫ۥػ‬˘˿ʽʼˋ [7] ʆɡʪʇ۵ɧʨʫʅɣʪ [8].

ਚʍԣ୐ɪʨʍ۵‫ޅ‬ɫ಴๗ʆɡʪ. ࠄ‫ݣ‬, ѼƧʎ೜ஞൣପ߲ʍߢԨʊʧʩமʨʫɾഷॐ (x1 (t), . . . , x8 (t))୐ , ࡀՎ୐ʍʧɥʉອӌ‫ػ‬ʍ೜ഷࡘ܏ʇɣʂɾҺ࠱ʍߢ‫ػ‬໑˙ƪˑʱ {u} ʇɶʅ˼ˈƪˢƪ˅̅˦˷ƪອӌ‫଺ػ‬१ࠃɫ‫ڇ݌‬ɶɥʪɲʇ෢ʨɪʊɶɾ [5]. ච

˘ʵ̅˂ʊʧʩӌࡌɴɺʅமʨʫɾߢԨౙୈ˴˙˽ٰ֖ʆʎ, ࠄ‫ݣ‬ʍ‫ࣛڇ‬ʆɶʏɶʏ‫ٵ‬ʨʫʪອӌ‫ػ‬ʍʱʈʍପ୩‫ڇ݌‬ɸʪɪ҈෢ɸʪ. ˪˘˿ʽʼˋ१ʎ

২‫ט‬१ʱగ઒ɴɺʪ‫ۥ‬ਚ೜Μଜʉອӌ‫ۥػ‬ਚʍ‫ڇ݌‬ອӌ‫ػ‬ʍ২‫ט‬१ʱగ઒ɴɺʪ୅‫଺ث‬ʉ‫ۥ‬ਚʆɡʩ, ɼ१ʍ۵‫ޅ‬ʊ‫ۈ‬ɰʅ, ɼʍʧɥʉ‫ۥ‬ਚ೜Μଜʉອӌ‫ػ‬ʍ‫ۥ‬ਚ೜Μଜ१ɪʨՑҗӌࡌ˴˙˼̅˂ʍ‫ܪ‬௟ɴɫ‫ۥ‬ਚʱʡʃອӌ‫ػ‬ʍՑҗӌࡌ˴˙˽ʍ‫ۥ‬યʱ‫ۼ‬ɥ.

฽Վɴʫɾ. றʊචߓʆʎએߢ‫ػ‬໑˙ƪˑ฽ਢʣ૫ߢ

2˼ˈƪˢƪ˅̅˦˷ƪ˘ʵ̅˂.

dφ/dt = f (φ) ʆ೅ɴʫʪອӌ‫ػ‬ʍഷॐ, u = h1 (φ)∈ RM ʊʃɣʅ۵ɧʪ. ɡʪߢ‫ܗ‬ʝʆʍߢ‫ػ‬໑˙ƪˑɪʨ௬ອഷॐ u(t) ʊ੆ɶʅࡰອɫ u(t + ∆t) ʇʉʪʧɥʉߢԨౙୈ˴˙˽ʱՑҗӌࡌʊʧʂʅ‫ۥ‬યɸʪ.

Ԩౙୈɴɺʅமʨʫɾߢ‫ػ‬໑ʍ੐‫ߢݥ‬Ԩഒೣʇɣʂ
ɾ೜ഷࡘ܏ʍ‫ڇ݌‬१ʊુ෾ɸʪ.

3.2

ٗѢ.

ʽʼˋ१ʉʈʊʧʩ૫ߢԨʍߢ‫ػ‬໑˙ƪˑ฽ਢʎ˼ˈƪˢƪ˅̅˦˷ƪ˘ʵ̅˂ʍற૙ʍφʃʎ˝‫ܪ‬௟ʆɡʪɫ, மʨʫɾՑҗӌࡌ˴˙˽ʊʧʂʅߢ

˷ƪ˻˽˟˕˚́ƪˁʍ‫ۥ‬ਚ߭੄ʱӌࡌɶʉɣɾʠ

‫ػ‬໑˙ƪˑʱॴ୩ວɮ฽ਢʆɬʅɣʪɲʇɫӂ௳ɴΏ୭଺ʊ‫ޟك‬ສɫࢭʉɮʅ‫ݗ‬ʞ, ਜ਼ʮʩʊ˝˷ƪ˻˽ʫɾ (ऺ 1). ʝɾ, Ցҗӌࡌ˴˙˽ʍ૫ߢԨౙୈʱ‫ٵ‬˟˕˚́ƪˁʍߣ‫ٿ‬ʱ੝ɬɮɶʅɣʪ୐ʆɡʪ (ࣚʪɲʇʊʧʂʅ९҈ʍ˿ƪ˾̅˖ൣପ߲ʍʴ˚˻ˁ‫ݟ‬ʎ [5]). ɲʍɲʇʊʧʩ˒ʶ˜˱ˁˋɫഉ‫ތ‬ʉʔˑʱ‫ڇ݌‬ɸʪɲʇʡӂ௳ɴʫɾ (ऺ 2).

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

3

0.2

reservoir output
actual

2

0.15
PDF

1
0

0.1

-1

0.05

-2

0

-3

reservoir
actual

0

5

10

15

-6

-3

0

3

6

9

x1

20

t

ऺ 2: ഷॐʍ੐‫ߢݥ‬Ԩഒೣ. ‫ۥ‬२ɶɾՑҗӌࡌ˴˙

ऺ 1: ߢ‫ػ‬໑˙ƪˑ฽ਢ. ɡʪߢ‫ ܗ‬t = 0 ɪʨ‫ۥ‬२ɶ˽ʧʂʅ฽ਢɴʫɾഷॐ x1 ʍߢ‫ػ‬໑˙ƪˑɪʨɼɾՑҗӌࡌ˴˙˽ʧʂʅ฽ਢɴʫɾഷॐ x1 ʍߢ‫ػ‬ʍ੐‫ߢݥ‬Ԩഒೣʱ‫ޟك‬ɶ˩˿˕˚ɶɾ. ʝɾ, ੆жɸ໑˙ƪˑ (ঐः) ʇ੆жɸʪ˿ƪ˾̅˖ൣପ߲ʍߢԨ

ʪ९҈ʇɶʅ, ˿ƪ˾̅˖ൣପ߲ʍߢԨౙୈʊʧʩ

ౙୈʊʧʩமʨʫɾഷॐ x1 ʍߢ‫ػ‬໑˙ƪˑ (ॾः)மʨʫɾഷॐ x1 ʍߢ‫ػ‬໑˙ƪˑɪʨ੐‫ߢݥ‬Ԩഒೣʱ࢑ɬࡰɶɾ.

ʡ‫ޟك‬ɶ܏ʮɺʅ˩˿˕˚ɶɾ.

4ʝʇʠ.
‫ۥ‬ਚ೜Μଜʉອӌ‫ۥػ‬ਚʍ୅‫ث‬ແʍφʃʍ˪˘˿ʽʼˋʆɡʪ˿ƪ˾̅˖ൣପ߲ʊʃɣʅ, ɼʍߢ‫ػ‬໑˙ƪˑʍʞʍӌࡌʊʧʂʅՑҗӌࡌߢԨౙୈ˴˙˽ʍ‫ۥ‬યʱ‫ۼ‬ʂɾ. ɲʍ˴˙˽ʊʧʩߢ‫ػ‬໑˙ƪˑ฽ਢʣ೜ഷ෋୩ഒೣʉʈɫ‫ڇ݌‬ɴʫʪɲʇʱ෢ʨɪ
ʊɶɾ. ʝɾ, ˣ˻˳ƪˑ (N, f ) = (8, 5.8), (8, 6.5)
ʉʈʍ࣪܏ʍ˿ƪ˾̅˖ൣପ߲ʊ੆ɶʅʡடํʊ˴
˙˼̅˂ɶ, ฽ਢɫɥʝɮ‫ۼ‬ɮɲʇʡӂ௳ɶɾ.

5[3] K. Nakai, and Y. Saiki, Physical Review E 98,(2018), 023111:1-6.
[4] K. Nakai, and Y. Saiki, Discrete and Continuous Dynamical Systems Series S, (2021),
14:1079-1092.
[5] M. Kobayashi, K. Nakai, Y. Saiki, and N.
Tsutsumi, Physical Review E 104, (2021),
044215:1-7.
[6] E. Lorenz, and K. Emanuel. Journal of the Atmospheric Sciences 45, (1998), 399414.

ࠗ߯઺πʎ, ӌ૫‫ݡ‬ສ‫ؼ‬ಛ-੝ӌҔӋƔՑఉֽѓஉीदߚ‫ז‬- ʍ࢕२ʱ࠷ɰʅɣʪ. ʝɾ, චٰ֖ʆɩɲʉʂ
ɾ‫ޟك‬ʍφ೼ʎ֮୤੝ӌʍˋƪˣƪ˅̅˦˷ƪˑִடٰ֖ॣ୩ ( ࠥࠬƔࢗ१ࢧໃ༄ ) ʊՂʄɮ. ɲɲʊԈࠗʍίʱ೅ɸ.

[7] Y. Saiki, H. Takahasi, and J. A. Yorke, Nonlinearity 34, (2021), 5744.
[8] Y. Saiki, M. A. F. Sanjun, and J. A. Yorke.
Chaos 28, (2018), 103110.

‫ޖ‬۵ഞ‫ٯ‬[1] H. Jaeger, and H. Haas, Scince, 304, (2004),
78-80.
[2] Z. Lu, J. Pathak, B. Hunt, M. Girvan, R.
Brockett, and E. Ott, Chaos 27, (2017),
041102.

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

計算化学的手法による様々な界面に対するカテコール構造の接着能解明
曽川 洋光
関西大学化学生命工学部

1 緒言

PVC 基板と PMMA 基板では，カルボン酸部位を有する PMMA が ALGDA のカテコール部分とより強近年，環境負荷低減の観点から，バイオマスを活く相互作用することが予想されたが，実際の結果は用した高性能材料の開発が強く求められている。申これに反したものとなった。
そこで本研究では，PVC

請者の研究グループでは，天然のイガイ接着タンパ

や PMMA に対するカテコールの相互作用の大きさ

1)

ク質の接着機構に学び ，海藻由来の多糖類であるを見積もることで，実測との比較検討を行なった。

アルギン酸（ALG）にドーパミン（DA）を化学修飾したカテコール含有アルギン酸（ALGDA）を合成し，

2 実験

これが乾燥条件下で良好な接着能を示すことを見出している（図 1）2)。ALGDA はとりわけマイカ基

第一に，PVC および PMMA 六量体モデルに対し板で最も良好な接着性を示し，次いでポリ塩化ビニて Gaussian 16 を利用し，DFT 法を用いた構造最適

ル（PVC）基板，アクリル（PMMA）基板の順に接化を行なった。PMMA 六量体モデルの立体規則性に着強度が低下することも明らかとなっている。

ついては，計算を容易にするためシンジオタクチッ

ALGDA に置換したカテコール基とマイカ基板上のク構造を採用した。また，DA との相互作用を算出

水酸基は，水素結合を介して強く相互作用すること

す る た め ， N-(3,4-dihydroxyphenethyl)acetamide

で，接着強度が大きく向上したと考えられる。一方，（DAamide）を PVC および PMMA 六量体モデル近傍に配置した構造最適化を行なった。各ポリマーモデ
ルに対し，DAamide の配置位置をいくつか変えたもの
を初期座標として計算を実施し，最も安定なエネル
ギー値をとったものを最適化後構造とした。いずれの計算も汎関数および基底関数にはそれぞれ
ωB97X-D，6-31G(d,p)を使用した。

3 結果・考察DFT 結果より算出したエネルギー値を表 1 に示す。DAamide と各ポリマー六量体モデルのエネルギー値の和を基準値とし，DAamide が各ポリマー鎖近傍に配置された際のエネルギー値を求めたところ，PVC
および PMMA の各ポリマー六量体モデルで，エネルギーの値はそれぞれ−94.9 kJ/mol および−63.3
Figure 1. Chemical structure of ALGDA and adhesive

kJ/mol 安定化されることが明らかとなった。
つまり，strength evaluated by tensile tests (experimental values).

ALGDA にアミド結合を介して導入されたカテコー

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

Table 1. Relative Energies estimated by DFT method
combinationa

ΔE (kJ/mol)

DAamide + 6-mer model of PVC

0.0

DAamide on 6-mer model of PVC–94.9

DAamide + 6-mer model of PMMA

0.0

DAamide on 6-mer model of PMMA

–63.3

Calculated by the DFT method (ωB97X-D/6-31G*).

The ΔE were relative values using DAamide + 6-mermodel of PVC, and DAamide + 6-mer model of PMMA
as standards for estimation the interaction of catechol
unit, respectively.

状態で再安定となった構造を図 2 に示す。PMMA では，主鎖上に位置するメチル基が立体障害として作用することで DAamide の近接を防ぎ，結果としてより相互作用しにくくなったものと推測される。加えて，
PVC ではカテコール部分だけでなく，アミド結合部
分もポリマー近傍に近いコンフォメーションを取っており，これが構造安定化に寄与している可能性も示唆された。

4 まとめと今後の展望
本研究では，ALGDA の界面に依存した接着能に関する知見を得るため，DAamide と PVC および
PMMA 六量体モデルの相互作用を DFT 計算から算出した。その結果，DAamide は PVC とより相互作用
しやすいことが明らかとなり，実測の傾向を反映した結果が得られた。このことは，基板モデルと相互作用する部位の構造安定化の度合いを見積もることで，接着強度の大きさをある程度予測できることを示しているとも言える。ただし，界面での分子のふるまい・挙動は複雑なものであることから，今回示したような単純な一分子同士の相互作用モデル
ではなく，分子動力学計算を用いた大規模シミュレーションによる相互作用予測が必要であるとも考えられる。これらの計算を実施することで，より正確な界面接着強度の計算科学的予測が達成可能になると期待される。

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

6 引用文献1. Y. Liu, K. Ai, L. Lu, Chem. Rev. 114, 5057–5115
Figure 2. Optimized conformation of DAamide on a 6-mermodel of PVC (top) and PMMA (bottom). Gray: carbon;
white: hydrogen; red: oxygen; blue: nitrogen, and green:(2014).
2. S. Inata, H. Sogawa, F. Sanda, Polym. Prepr. Jpn.
70(2), 1Pc53 (2021).

chloride.
ルユニットは，PMMA 鎖よりも PVC 鎖とより強く相互作用することが計算科学的にも示された。
DAamide が各ポリマー六量体モデルと相互作用した

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

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

高密度領域まで適用可能なモンテカルロ法の開発と有限密度２カラーQCD の相図の決定
伊藤 悦子
理化学研究所 数理創造プログラム(iTHEMS)
スーパーコンピュータ共同研究制度の「2021 年度プログラム高度化」により、格子 QCD の分野でよく用いられるラショナ
ルハイブリッドモンテカルロ (RHMC)法のコードの高速化を行なった。
演算・通信の改良により12 桁の精度を保ったまま、
約 2.4 倍の高速化に成功した。

1 研究概要高速化が急務となり、今回の共同研究に至った。

本共同研究の背景には、関連研究として 2018年度から学際大規模情報基盤共同利用・共同研究
拠点(JHPCN)の公募型共同研究課題 「高密度領域まで適用可能なモンテカルロ法の開発と有限密度２カラーQCD の相図の決定」でシステム A を利用していた事がある。
JHPCN での研究課題では、原子核の中の物理を記述するミクロな理論である量子色力学(QCD)のゲージ対称性の自由度を一つ落とした「２カラーQCD」の高密度領域で生じる超流動相の性質を解明することが目標となっている。これは、現在全くと言っていいほどわかっていない「中性子星の内部の現象に新たな知見を得る」と期待されて

図 1: 2 カラーQCD 相図 (T=0.79Tc と 0.39Tc
を調べた)

いる研究テーマである。
これまでに、図 1 のように 2 カラーQCD の相

2 計算コードについて図の決定をおこなってきた。主な成果としては超流動性の秩序変数であるダイクォーク凝縮を測定

本研究は、いわゆる「格子 QCD」分野の第一し、超流動相の直接証拠を得た事や、その超流動原理計算である。その計算手順は大きく２つであ

相は 100MeV 以下の低温で出現する事、またクォる。

ークの閉じ込めは BCS 相でも引き続き起こることなどが示された(文献[1],[2],[3])。つまり摂動論的な描像ではない超流動相が広がっていることを示しており、非摂動的な第一原理計算の重要性がStep 1: モンテカルロ法を使い「ゲージ配位」
をサンプリングする
Step 2: 生成した配位を使って物理量を測定する高まったと言える。今後、超流動相でのハドロン

Step 1 で一旦生成された配位は保存しておくのスペクトル、相互作用ポテンシャル、状態方程事で、いつでも新しい物理量を測定するのに活用

式などの物理現象を調べるにあたって、コードのでき、分野の財産となっている。実際に、国際的30 − Academic Center for Computing and Media Studies, Kyoto University

4 高速化の実行に公開・共有されている配位も存在する。
我々の研究においてもJHPCN で獲得した計算
時間のほとんどを Step 1 の配位生成に利用して

4.1 共同研究の組織おり、例えば 32^4 の格子サイズの配位生成は足ヒューレットパッカード(HPE)から伊藤泰善、

掛け３年に渡り行っている。このため、この部分

原山卓也。JHPCN での研究組織から、飯田圭の計算コードを少しでも高速化することで、研究(高知大) 石黒克也(高知大)、村上耕太郎(京に要する時間を大幅に縮減できる。

都大) 、伊藤悦子(理研)が参加した。基本的配位生成のアルゴリズムは、ハイブリッドモンには、HPE の二人が中心になって進めて、

テカルロアルゴリズムにラショナル近似を使った

JHPCN 側からの参加者は、時々進捗等を伺う「ラショナルハイブリッドモンテカルロ法程度だった。

(RHMC)」を用いている。これは、現実の QCD
でも奇数個のクォークが存在する際に使われる手法であり、非常に汎用性が高い。オリジナルのコ
ードは Fortran で書かれ、創業数 6587 行、並列化手法は分散メモリ型、MPI 並列化されていた。

3 プログラムの概要と高速化の問題点
高速化の対象となる RHMC アルゴリズムのコー

4.2 プログラムの構造の特徴
HPE の方によるプログラム分析結果では、(i)サブルーチンコールのネストの深さが大きい(最大 15 段)、(ii)ループ構造の末端で呼ばれるサブルーチン(ベクトル関連演算、MPI関連)の粒度が小さい、(iii)MPI 関数は全て同期型関数を使用している点が挙げられた。

4.3 高速化の検討事項

ドは、以下で成り立っている。

(4)のマルチシフトソルバーの演算はインラ(1) ゲージ配位(バイナリーデータ)の読み込み

イン展開/ループ融合/スレッド並列化/BLAS(2) 乱数の生成ライブラリ置換、通信はインライン展開/非(3) 行列の足し算・掛け算

同期関数化/複数要素化の組合せで高速化の(4) マルチシフトソルバーによる逆行列の計算可能性を検討した。また、(3)のベクトル(行(5) ゲージ配位(バイナリーデータ)の書き出し列)の演算は、インライン展開/スレッド並列

また、MPI 並列に関しては、現在対象としている

化、通信はインライン展開/非同期関数化/演

計算は時空点が 32^4 の格子サイズであり、1core

算・通信オーバーラップの組合せで高速化のあたり 8x4x4x16 の格子サイズの 512MPI 並列に可能性がある。そのほか、コンパイルオプシ

分割している。演算のほとんどが各 core の計算で

ョン、Hugepage、 Cray MPI 環境変数、rank

あり、時々通信を行なっている状況である。

配置等での高速化の可能性がある。

計算で最も時間のかかっている部分は、(4)のマ

一方で、並列化については、並列 I/O の適

ルチシフトソルバーの部分である。また、この時用による性能向上は期待できず、既存の I/Oに何度も(3)の行列演算のライブラリーが呼ばれを含むサブルーチンのデータコピー部分をている。また、(1),(2),(5)の部分は並列化されておスレッド並列化することが、最小限の修正でらず効率が悪くなっているため、この辺りを MPI

最大の効果をもたらすことがわかった。

化する改善も、当初は見込まれた。
基本的に物理学の研究者で構成されている
JHPCN の共同研究のメンバーでは、これ以上のチューニングが難しかったため、このプログラム

4.4 演算・通信高速化の実施
高速化にあたっての実施条件は以下のとおりである。

高速化の共同研究制度を利用して、ヒューレッド
パッカードの方を紹介していただき、共同研究を行うこととなった。

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

この元で、まずは演算部分の高速化として、６

明示的指定)オプション」「複数要素化」「アルゴリ

種類の演算高速化手法を組み合わせて、オリジナ

ズムの変更」「プロセスデコンポジション」「演算・

ル版(に PPN オプションを実行時に付与したも通信オーバーラップ」「サブルーチンマージ」を行の)に対して、どの程度性能が向上するかを調べたなっている。効果として最も効いたのは、「複素要結果が下の図２である。

素化」の部分で、ここで 8.1%の性能向上が見られた。
これら６種類の通信高速化のみを行うと、合計
として 1.265 倍の性能向上があり、4.23[h]かかっていた計算が 3.32[h]となった。
最後に、演算高速化と通信高速化の両方を組み合わせた結果が以下の図 4 である。

図 2: 演算高速化性能評価結果のまとめ(HPE からの報告書より抜粋)
具体的には、左から「インライン展開」「ループ
融合」「スレッド並列化(omp2)」「演算・通信オー
バーラップ」「ループスケジューリング」「SIMD化」をおこなっている。
効果として最も効いたのは「演算・通信オーバ

図 4: 演算・通信高速化性能評価結果のまとめーラップ」の部分で、ここで 13.8%の性能向上が見られた。

具体的には、左から「PPN（通信）」「インライ

６種類の演算高速化を行うことで、合計するとン展開（演算）」「ループ融合（演算）」「SIMD 化1.622 倍の性能向上があり、4.02[h]かかっていた（通信）」「複素要素化＋アルゴリズムの変更（通計算が 2.48[h]となった。

信）」「プロセスデコンポジション（通信）」「演算
通信オーバーラップ」「サブルーチンマージ（通信）」

次に、通信部分の高速化についてまとめると、を施している。全体を組み合わせた際は、SIMD

こちらも６種類の高速化手法の組み合わせを行い、化が大きく性能向上に影響し、56.2%の効果を上オリジナル版(に PPN オプションを実行時に付与

げた。

したもの)に対して、どの程度性能が向上するかを調べた結果が下の図 3 である。

これら全てを取り入れると、オリジナル版（にPPN 適用）した base コードと比較して、4.23[h]から 2.02[h]の実行時間の削減となり、2.088 倍の性能向上となった。

5 まとめと今後の展望
最終的に、本当のオリジナル版(PPN 適用なし)のコード実行から比較して、2.424 倍の性能向上が得られた。一方で実行結果は 12 桁以上の精度
図 3: 通信高速化性能評価結果のまとめで一致しており、計算精度をほとんど落とすことなく高速化ができた事になる。約 2.4 倍の高速化

具体的には、左から「PPN(プロセッサー数の32 − Academic Center for Computing and Media Studies, Kyoto Universityというのは驚くべき成果である。おそらく物理の

研究者のみではここまでのチューニングを行うことはできなかったため、本共同研究ならではの成果と言えるだろう。
2022 年度は早速この高速化コードをシステム
A で利用し、研究を進めている。
また、
この RHMC
コードは、格子 QCD 分野で標準的なコードの一つとなっており、他のグループとも情報シェアが可能なため、ここでの報告の情報が広く伝わることを期待する。
参考文献[1]K.Iida, E.Itou, T.-G. Lee:
PTEP 2021 (2021) 1, 013B05[2] K.Iida, E.Itou, T.-G. Lee:JHEP2001(2020)18[3] K.Ishiguro, K.Iida, E.Itou,
PoS (LATTICE2021) 063

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

低メモリ・低コスト通信乱流大規模直接数値計算コードの高度化
山本 義暢
山梨大学大学院総合研究部

ノイマン型スーパーコンピュータの基礎を作ったフォン・ノイマンはその活用先として、乱流現象解明へのブレイクスルーを想定していた。実際彼の予測通り、1972 年には初めての乱流直接数値計算(Direct Numerical Simulation, DNS)が行われ [1]、計算機性能の発展に伴いその対象レイノルズ数の拡張が試みられている(例えば、 [2], [3])。

乱流の DNS においては周期境界条件が適用できる形状においては、フーリエスペクトル法の適用が(FX シリーズ;FX1, FX10, FX100 及び FX1000)を導入
しているのに対し、京都大学は x86 タイプのスカラ並列機(T2K, Cray XE6, XC40, 及び 2023 年からはDELL C6620)が導入されている。一方東北大学はベクトル並列機(SX シリーズ; SX-9, SX-ACE, SX-Aurora
TSUBASA)が導入さているのが特徴である。
103

SX-Aurora TSUBASAGFLOS/core

1 研究背景SX-9

102

SX-ACE

101

Kyoto U.
Nagoya U.
Tohoku U.

計算精度・効率の両面から最高のパフォーマンスを与え、標準手法として確立している。このスペクト

100
2000

ル法を流体計算に適用する場合、Navier-Stokes 式のル法が適用される。疑スペクトル法では非線形項を物理空間で計算し、フーリエ空間に戻して計算を行う。そのため疑スペクトル法の計算コストは高速フーリエ変換(Fast Fourier Transform, FFT)のコストと同等となる。しかし FFT は並列効果を作るためには

101

ば、 [4], [3]）。ここで FFT の並列化効率は、MPI プ

2025

FX1
HPC2500

10-1
2000

FX10
T2K

100

Kyoto U.
Nagoya U.
Tohoku U.

2005

は通信コストが極めて高く、大規模化つまりは対象
Re 数の拡張において深刻な問題となっている（例え

2020

102

領域分割軸の転置(all-to-all 通信を伴う)が必要となり、近年のノード間結合が疎結合の超並列計算機で2010 2015year(a) 理論演算速度/core
GB/core

非線形項が主要演算部となり、一般には疑スペクト

2005

XE6

XC40

FX100

2010 2015
year

Dell
FX10002020

2025(b)メモリ量/core

図 1 京大、名大、東北大におけるスーパーコン

ロセスあたりのメモリ量が大きいほど有利となる。

ピュータにおける理論演算速度・メモリ量の変化,

図 1 は京都大学、名古屋大学及び東北大学の全国共(a)演算速度/core、(b)メモリ量/core

同利用施設として提供されているスーパーコンピ
ューにおける flat MPI を想定した場合の 1MPI プロ
セスあたりの理論演算速度(GFLOPS)とメモリ量(GB)の経年変化を示している。京都大学と名古屋大学は 2000 年初めまでは同機種(富士通 HPC2500)であったが、名古屋大学は富士通独自のスカラ並列機

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

興味深いことに、図 1(a)に示すようにコア当たりの演算速度は増加し続けているに対し、FX シリーズ
及び SX シリーズのメモリ量は減少傾向にあることがわかる(図 1(b)参照)。一方京都大学のシステムは演算性能を向上させつつ、メモリ量は一定に保たれていることが確認できる。つまり乱流の大規模 DNS

ル法を適用し、壁垂直方向(z)に 2 次精度中心差分法

る。

Adams-Bashforth 法を適用する。波数空間での圧力においては、その優位性が FX 及び SX シリーズに対して x86 タイプとの差が縮小していることを示唆す

そこで本研究では京都大学学術情報メディアセン

ターのシステム A(Cray XC40)を対象とし、筆者らがベクトル並列計算機場上で開発してきた壁面乱流場の大規模 DNS コード [5] [6]の移植評価及びその高速化を試みる共同研究を実施した。

2 DNS コードの概要を適用する。
時間積分はFractional step法に基づき、
圧力項に Euler 陰解法、その他の項に 2 次精度
ポアソン方程式は、3 重対角行列解法(Tri-Diagonal

Matrix Algorithm, TDMA)による直接解法が適用される。従って主要演算部は、x, y 方向への 2 次元実フーリエ変換/逆変換と z 方向への TDMA となる。FFTについて FFTW3.0 の逐次版を用い、コンパイラは

Cray コンパイラを使用した。

2.3 領域分割方法

2.1 対象とする流動場

図 2 チャンネル流とその座標系

・従来方法

図 3 2 次元領域分割

領域分割方法として、超並列計算に本研究では図 2 に示す十分に発達したチャンネル対応できるように図 3(a)に示す y 及び z 方向への 2

2h, h:チャンネル半値幅)において no-slip 条件が適用

従来の方法 [4], [8]では並列効果を作るために y 方流を対象とする。本流動場では、主流方向(x)及びス
パン方向(y)に周期境界条件が適用でき、壁面(z=0,される。流れ場の支配方程式は非圧縮性流体のNavier-Stokes 式と連続式であり、基本物理量は x, y,z 方向の速度𝑢𝑢�, ��, ��及び一定密度で除した圧力𝑝𝑝� となる。ここに上付き添え字~は瞬間値であることを示している。流れ場の駆動力としては x 方向に一定の圧力勾配(Fx=u2/h, u:壁面摩擦速度)を課す。従って流動場の支配パラメータであるレイノルズ数は、次元領域分割を適用する。しかし y 方向の FFT 及びz 方向への TDMA は領域分割方向への演算となり、
向 FFT 時には y→x への分割軸の変更(図 3(b)参照)

及び z 方向 TDMA 時には z→y への分割軸の変更(図3(c))を行う必要がある。つまり速度 3 成分、非線形

項6 項とすると1step あたり、
x→y への転置を9 回、
z→y, y→z の転置を 1 回、y→x への転置を 9 回行うことになる。この転置は各領域分割コミュニケータ

内での all-to-all 通信(a2a)により行われる。従って乱壁面摩擦速度(u)とチャンネル半値幅(h)によって規流の DNS における計算時間の大半はこの転置によ流のカノニカル流として知られており、1987 年に初向の TDMA に対し、shift 通信と x, y 方向へのオーバ

ネル流の DNS は一様等方性の DNS [1], [2], [3]に比べ置は必要となるが、この分割方向(xy_communicator、

2.2 離散化手法

Aries+Dragonfly ネットワークシステムの場合は、4定される(Re=uh/, :動粘性係数)。

本流動場は実際の工学現象で頻出する非一様乱

めての DNS [7]が行われた(x, y 方向にフーリエスペ

クトル法、z 方向に Chebyshev-tau 法)。このチャン

時間積分長が大きくなり、計算コストはさらに増大する。

る通信時間が占めることになる [4], [8]。
・本手法

この転置コストを低減するために、
z方ーラップによる並列化手法 [9]を適用し、z→y, y→zの転置を回避する。一方 FFT 時の x→y, y→x への転図 3 参照)を直接結合ノードに割り当てることにより通信負荷を低減させる。Cray XC40 の Cray
ノードが直接結合ノードとなる。

離散化手法としては、x, y 方向にフーリエスペクト

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

.2.4 a2a 通信の低減化

テムA の場合は実際のノード配置を指定することが本 DNS コードでは、FFT 時の x→y, y→x への計 18

回の a2a が発生する。これを低減させるために、演できないため、この 4 ノードが直接結合ノードに配置される場合とそうでない場合がランダムに発生し、性能評価において大きな負担となった。

算と通信のオーバラップを導入している。具体的には x, y 方向FFT は z 方向へ独立に行われることから、
z 方向 1 断面の a2a とその z 方向前後の FFT を同時に実行し a2a 通信の低減を図っている。この実装はOpenMP 機能を用いている。

4 性能評価結果
表 2 簡易体系での性能評価結果：camphor2

FFT/a2a オ

ーバラップ

3 性能評価条件ありなし

3.1 解析実体系
超大規模計算として想定している計算条件は x, y, z
方向に 14400, 13824, 5760 格子であり、必要メモリ量は約 150TB を要する。システム A(camphor2)では

全 1920 ノード(メモリ量 160TB)で動作する規模であり、xy_communicator に 88 分割、z 方向 480 分

割、各 MPI プロセスに 3 つの OpenMP スレッドを配置する（88x480x3=126,720 コア）。但しシステムA ではこの体系での計測はできないため、簡易解析「 富 岳 」 (Fugaku) 東 京 大 学 情 報 基 盤 セ ン タ ー
Oakforest-PACS(OFP)での実行条件を表 1 に示す。
表 1 実体系での実行条件

system

#node

#mpi

#mpi

OFP

1920

88

480

Camphor2

FugakuFugaku

1920

5760
8640

3.2 簡易解析対象

_xy
88

192
144

#

_z

AP

480

3

480

720

3

3

4

#core

126,720

126,720

276,480
414,720システム A(Cray XC40 8 ノード/528 コア)を対象として、x, y, z 方向に 14400, 13824, 24 格子を適用する。領域分割は xy_communicator に 88 分割、z 方向に 2 分割とし各 MPI プロセスに 3 つの OpenMP ス

レッドを配置した(88x2x3=528 コア)。この場合、

xy_communicator 内の a2a は 4 ノード間となり、直接結合ノードでの実行を想定している。但し、シス36 − Academic Center for Computing and Media Studies, Kyoto University

#mpi

#mpi

8

88

2

8

_xy
88

_z
2

#

elapse

3

2056

AP
3

[s]

1842

表 2 に簡易体系を用いて 100step 分の時間積分に要した経過時間(wall-time)を示す。FFT と a2a のオーバラップを用いた場合、18.4[s/step]を要する。また FFT と a2a のオーバラップにより約 200 秒の時

間短縮効果が得られていることが確認できる。

表 3 に理化学研究所「富岳」(Fugaku)東京大学情報基盤センターOakforest-PACS(OFP)を用いた実体系での測定結果を示す。

対象として z 方向に縮小した条件を設定した。なお参考のため、本体系を対象とできる理化学研究所#node

表 3 実体系での測定結果

system

#node

#mpi

#mpi

Fugaku

5760

192

480

OFP

Fugaku

1920

8640

_xy
88

144

#

s/step

3

4.8

_z

AP

720

4

480

3

18.1
3.2

OFP とシステム A ではネットワークシステムがことなるため超並列計算時の単純な比較には注意を要するが、表 2 の簡易体系で見積もった 18[s/step]とほぼ同等のであることが確認できる。理化学研究所の「富岳」においては、8640 ノードで 3.2[s/step]

演算速度：516[TFLOPS]での計算が可能であり、5760 ノードとの比較により算出した並列化効率は99%を超えている。
「富岳」での実行効率は理論性能の 1.8%程度であり、これを基に OFP での演算速度を見積もると、約 92[TFLOPS]となる。これは理論性

能の 1.6%に相当する。
なお OFP では大規模 HPC チ
ャレンジにより 8192 ノードを使用しての評価も予定していたが、新型コロナウイルス感染症拡大防止によりキャンセルとなってしまった。

5 高速化の検討成を受けて遂行した。記して謝意を表する。

参照文献
本 DNS コードのさらなる高速化を目指し、非同期型 a2a による a2a 時間の短縮、MPI プロセス配置の手動設定等を検討したが、現段階で有効な結果は得られなかった。この理由は実体系では z 方向の MPIプロセスを 480 程度(480 ノード)用いるため、TDMAの評価も重要となるが、簡易体系では高々2 ノードであるため TDMA 部分について簡易体系となっていないことが大きいと考えられる。またこれを評価するにはかなりノード数を増やすことが必要であり、今回の検討規模（8 ノード）では困難であった。

6 まとめと今後の課題
本研究では、壁面乱流場の超並列大規模直接数値
コードにおける京都大学システムＡ(camphor2)での性能評価を実施した。従来乱流の DNS においては、

NEC/SX シリーズあるいは富士通/FX シリーズの適用が有用であったが、コア当たりのメモリ量が小さ
く、コア当たりのメモリ量及び演算速度のバランスに優れる Xeon Phi ベースの Camphor2 との差が短縮することが想定された(図 1 参照)。実際実コードを用いた評価においても Camphor2 と同じ Xeon Phi を用いる東京大学情報基盤センターOakforest-PACS に

おいては、理化学研究所「富岳」(FX1000 相当)に匹敵する実行効率が得られることが確認できてた(表3 参照)。京都大学の次期システム（DELL C6620）に
おいては Camphor2 と同等のコア当たりのメモリ量が確保される予定であるため、本 DNS コードが有効に動作する可能性が高い。

一方京都大学システムでは最大実行可能ノード数がかなり小さく抑えられているため実際のPRODUCTION RUN は実行できずまさに絵に描いた餅となってしまうことが残念である(Camphor2 では全 1920 ノードに対し、最大ジョブは高々128 ノードしか利用できない)。今後は本評価に基づき、さらなる大規模化と高速化を実現するための基本戦
略の策定、さらには実際の PRODUCTION RUN を実行できる環境の模索も重要と考えられる。

謝辞：本研究は京都大学学術情報メディアセンターとの共同研究（高度化研究）により実施した。理化学研究所「富岳」を用いた評価は、HPCI 公募課題(ID:hp210125)により提供された資源を利用して行った。
また本研究は栢森情報科学振興財団による研究助[1] S. A. Orszag and G. S. Patterson, "Numerical
Simulation of Three-Dimensional Homogeneous

Isotropic Turbulence," Phys. Rev. Lett., vol. 28, no. 2,
pp. 76-79, 1972.

[2] Y. Kaneda, T. Ishihara, M. Yokokawa, K. Itakura and
A. Uno, "Energy dissipation rate and energy

spectrum in high resolution direct numericalsimulations of turbulence in a periodic box," Phys.
Fluids, vol. 15, no. 2, pp. L21-L24, 2003.

[3] P. K. Yeung and K. Ravikumar, "Advancing

understanding of turbulence through extremescale computation: Intermittency and simulationsat large problem sizes," Phys. Rev. Fluids, vol. 5, no.
11, p. 110517, 2020.

[4] J. Xu, "Benchmarks on tera-scalable models for

DNS of turbulent channel flow," Parallel
Computing, vol. 33, no. 12, pp. 780-794, 2007.

[5] Y. Yamamoto and T. Kunugi, "Direct numericalsimulation of a high-Froude-number turbulent

open-channel flow," Phys. Fluids, vol. 23, no. 12, p.

125108, 2011.

[6] Y. Yamamoto and T. Kunugi, "MHD effects on
turbulent dissipation process in channel flows
with an imposed wall-normal magnetic field,"

Fusion Eng. Des., vol. 109, no. 111, pp. 1137-1142,
2016.

[7] J. Kim, P. Moin and R. Moser, "Turbulence statistics
in fully developed channel flow at low Reynolds
number," J. Fluid Mech., vol. 177, pp. 133-166,

1987.

[8] M. Lee, N. Malaya and R. D. Moser, "Petascale direct

numerical simulation of turbulent channel flow on
up to 786k cores," in Proceedings of theInternational Conference on High Performance
Computing, Networking, Storage and Analysis,
2013.

[9] Y. Yamamoto and Y. Tsuji, "Numerical evidence oflogarithmic regions in channel flow at Reτ= 8000,"
Phys. Rev. Fluids, vol. 3, no. 1, p. 012602, 2018.

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

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

サービス再開日時

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

9:00

2022/12/02

14:00

5.00

2022/03/05

6:00

2022/03/07

9:00

51.00

2022/03/30

9:00

2021/03/31

24:00

39.00

サービス再開日時

ダウン時間[h]なし２）サービス状況
サービス
時間[h]

ジョブ
処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼動
ノード数

ノード
利用率

10 月

744.0066,087

378,923

74,112,416

61,047,462

1794.3

80.3

%

11 月

720.0058,052

340,178

66,744,907

60,079,010

1797.0

78.3

%

12 月

744.00104,339361,725

75,090,388

63,162,220

1798.1

82.8

%

1月

744.0087,189320,878

73,046,361

54,796,703

1799.4

79.3

%

2月

672.0088,311

271,140

64,812,314

47,165,309

1798.1

79.6

%

3月

654.0097,807

207,969

62,517,789

45,850,338

1705.070.9%

計4,278.00

501,785

1,880,812

416,324,175

332,101,042

1782.078.5%

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

10月

11月

12月

1月

2月

3月

l

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

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

l

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

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

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

サービス再開日時

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

9:00

2022/12/02

14:00

5.00

2022/03/05

6:00

2022/03/07

9:00

51.00

2022/03/30

9:00

2021/03/31

24:00

39.00

サービス再開日時

ダウン時間[h]なし２）サービス状況
サービス
時間[h]

ジョブ
処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼動
ノード数

ノード
利用率

10 月

744.00109,110509,050

13,358,964

10,408,289

807.9

66.7

%

11 月

720.00135,180

921,764

13,768,690

11,345,109

807.8

67.3

%

12 月

744.0095,017652,116

14,258,013

11,249,362

807.8

70.3

%

1月

744.00136,636

809,456

14,688,547

11,116,104

807.9

71.8

%

2月

672.00153,039

634,854

13,198,992

10,890,750

808.0

70.3

%

3月

654.00101,477355,201

12,944,393

10,920,315

787.9

66.1

%

計

4,278.00

730,459

3,882,440

82,217,598

65,929,928

804.6

68.8

%

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

10月

11月

12月

1月

2月

3月

l

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

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

l

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

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

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

サービス再開日時

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

9:00

2022/12/02

14:00

5.00

2022/03/05

6:00

2022/03/07

9:00

51.00

2022/03/30

9:00

2021/03/31

24:00

39.00

サービス再開日時

ダウン時間[h]なし２）サービス状況
サービス
時間[h]

ジョブ
処理
件数

経過
時間[h]

占有
時間[h]

CPU
時間[h]

平均稼動
ノード数

ノード
利用率

10 月

744.006,017

22,211

368,371

207,114

16.0

60.8

%

11 月

720.008,581

28,686383,286

262,904

16.0

67.3

%

12 月

744.002,550

43,474341,437

213,482

16.0

58.2

%

1月

744.002,554

39,770250,922

153,604

16.0

53.8

%

2月

672.001,403

31,093426,882

300,515

16.0

59.6

%

3月

654.003,936

31,659321,845

167,326

15.8

55.9

%

計

4,278.00

25,041

196,8932,092,742

1,304,944

16.0

59.2

%

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

10月

11月

12月

1月

2月

3月

l

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

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

l

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

40 − Academic Center for Computing and Media Studies, Kyoto University京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.21 No.1（ Dec. 2022）− 41

資料

別表１

スーパーコンピュータシステム
区分

提供サービス
利用負担額

コース
エントリ

パーソナル

グループ

タイプ

セット

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

-

基本

12,600 円/年

B

共有

最大0.5ノード相当((112コア、512GBメモリ)×0.5)

1

0.2

基本

100,000 円/年

A

共有

最大2ノード相当((112コア、128GB高速メモリ)×2)

168

8.0

-

タイプB

基本

100,000 円/年

B

共有

最大2ノード相当((112コア、512GBメモリ)×2)

168

8.0

-

タイプC

基本

100,000 円/年

C

共有

最大1ノード相当((112コア、2048GBメモリ)×1)

168

8.0

-

タイプG

基本

100,000 円/年

G

共有

最大1GPU相当((16コア、128GBメモリ＋1GPU)×1)

168

8.0

-

タイプA0

最小/追加

80,000 円/年

準々優先

1ノード((112コア、128GB高速メモリ)×1)

168

6.4

2

タイプA1

最小/追加

200,000 円/年

優先

1ノード((112コア、128GB高速メモリ)×1)

336

16.0

4

タイプA2

最小/追加

240,000 円/年

準優先

2ノード((112コア、128GB高速メモリ)×2)

336

19.2

6

タイプA3

最小/追加

300,000 円/年

占有

1ノード((112コア、128GB高速メモリ)×1)

336

16.0

4

タイプB0

最小/追加

88,000 円/年

準々優先

1ノード((112コア、512GBメモリ)×1)

168

6.4

2

タイプB1

最小/追加

220,000 円/年

優先

1ノード((112コア、512GBメモリ)×1)

336

16.0

4

タイプB2

最小/追加

264,000 円/年

準優先

2ノード((112コア、512GBメモリ)×2)

336

19.2

6

タイプB3

最小/追加

330,000 円/年

占有

1ノード((112コア、512GBメモリ)×1)

336

16.0

4

タイプC0

最小/追加

96,000 円/年

準々優先

1ノード((112コア、2048GBメモリ)×1)

168

6.4

2

タイプC1

最小/追加

240,000 円/年

4

タイプG0

最小/追加

64,000 円/年

タイプG1

最小/追加

320,000 円/年

大規模ジョブ
タイプB

ストレージ

システム資源タイプAタイプA専用クラスタシステム バッチ

最小80,000 円/週(7日）

追加単位

20,000 円/週(7日）

最小

88,000 円/週(7日）

追加単位

22,000 円/週(7日）

A

B

C
G
A
B

優先

1ノード((112コア、2048GBメモリ)×1)

336

16.0

準々優先

1GPU((16コア、128GBメモリ＋1GPU)×1)

168

6.4

2

優先

2GPU((16コア、128GBメモリ＋1GPU)×2)

336

32.0

8

占有

8ノード((112コア、128GB高速メモリ)×8)
2ノード((112コア、128GB高速メモリ)×2)

占有

8ノード((112コア、512GBメモリ)×8)
2ノード((112コア、512GBメモリ)×2)

最小

660,000 円/年

追加単位

330,000 円/年

大容量ストレージ

10,000 円/年

大容量ストレージ容量10TBの追加につき高速ストレージ

10,000 円/年

高速ストレージ容量2TBの追加につきネットワークストレージ

5,000 円/年

ネットワークストレージ容量1TBの追加につき20,000 円/年

可視化ソフトおよびプリポストウェアの1ライセンスにつき

タイプB

ライセンスサービスB--

2ノード((112コア、512GBメモリ)×2)
1ノード((112コア、512GBメモリ)×1)

168
168
-

-

-

-

-

-

-

-

-

32.0

8

16.0

4

備考
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 ．スパコン連携サービス
学術情報メディアセンターのスーパーコンピュータシステムと密な連携により、学内における部局の組織が計算サーバ等を設置する場
下記の負担額を支払うものとする。
区分
スパコン連携サービス

冷却方式

利用負担額

水冷

16,800 円/月

水冷冷却方式の計算サーバ等の定格電力 1kWにつき利用負担額算定単位

空冷

20,000 円/月

空冷冷却方式の計算サーバ等の定格電力 1kWにつき

42 − 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.21 No.1（ Dec. 2022）− 43

別表３
システム

A

B

スーパーコンピュータシステム（民間機関利用）
経過時間
（時間）

ストレージ
（TB)

無料
利用者数

2ノード((112コア、128GB高速メモリ)×2)

336

19.2

6

960,000 円/年

3ノード((112コア、128GB高速メモリ)×3)

336

28.8

9

1,440,000 円/年

4ノード((112コア、128GB高速メモリ)×4)

336

38.4

12

1,920,000 円/年

2ノード((112コア、512GBメモリ)×2)

336

19.2

6

1,056,000 円/年

3ノード((112コア、512GBメモリ)×3)

336

28.8

9

1,584,000 円/年

4ノード((112コア、512GBメモリ)×4)

336

38.4

12

2,112,000 円/年

システム資源利用負担額

備考
1．利用負担額は、年度単位で算定している。また、総額表示である。
年度途中から利用を開始する場合及び年度途中で利用を終了する場合の利用負担額は、上記表中の利用負担額を12で除した後、利用月数を乗じて算出するものとし、100円未満に端数が出た場合は、10円単位を四捨五入するものとする。
なお、月途中から利用を開始する場合及び月途中で利用を終了する場合は、それぞれ1月の利用とする。
2．ストレージ容量はバックアップ領域（最大で総容量の1/2）を含む。

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

全国共同利用版広報・Vol.20(2021)総目次［巻頭言］
Vol.20，No.1 号の発刊に当たって ................................................................................................................1-1
Vol.20，No.2 号の発刊に当たって ................................................................................................................2-1［スーパーコンピュータ共同研究制度（若手・女性研究者奨励枠）研究報告］
機械攪拌中浮遊油膜巻き込みメカニズムの解明 ..........................................................................................1-2
高効率有機系太陽電池の実現に向けた光機能性分子の構造と電子物性の相関解明 ................................1-4
時系列リモートセンシングデータによる全球陸域統合環境モニタリング................................................1-6
自己集合性 1,3,5-Tri(benzoimidazolyl)benzene 誘導体の構造解析 .........................................................1-8
界面分光の分子動力学シミュレーション....................................................................................................1-10
Numerical simulation of thermal-solutal Marangoni convection in a shallow rectangular modelwith mutually perpendicular temperature and concentration gradients........................................1-12
浮遊培養法による iPS 細胞作製技術の最適化を目的とした流動解析 .....................................................1-14
数値シミュレーションを用いたウマの個体間に働く力の推定 .................................................................1-16
ローレンツ方程式に対する機械学習によるモデリング.............................................................................1-18
力学系と代数幾何の理論に基づくデータ駆動的スパースモデリング .....................................................1-20
第一原理計算による BaSi2 太陽電池の接合パートナー探索....................................................................1-22
統計物理学模型との対応を利用した機械学習手法の性能評価.................................................................1-24
遷移金属による協奏的結合組み換えを利用した効率的分子変換 .............................................................1-26［プログラム高度化共同研究報告］
高周波加熱と磁気ノズルを利用した無電極プラズマ推進機の 2 次元軸対称プラズマ粒子計算...........1-28［中島浩先生追悼特集］
中島先生を偲んで ............................................................................................................................................2-2
師曰く................................................................................................................................................................2-3
中島浩先生を偲んで ........................................................................................................................................2-4
中島浩先生を偲んで ........................................................................................................................................2-5［シンポジウム開催報告］
α×SC2021K 行動学とスーパーコンピュータに関するシンポジウム 開催報告 ..................................2-6［サービスの記録・報告］
スーパーコンピュータシステムの稼働状況...................................................................................... 1-32，2-8
センター利用による研究成果（2020 年度） ..............................................................................................2-11
データ消失障害報告 ......................................................................................................................................2-15［資料］
大型計算機システム利用負担金 別表...............................................................................................1-36,2-16
全国共同利用版広報・Vol.19(2020)総目次 .................................................................................................1-39
サービス利用のための資料一覧 ..........................................................................................................1-41,2-19［訃報］
訃報 .................................................................................................................................................................1-42
京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.21 No.1（ Dec. 2022）− 45―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
URL: http://www.iimc.kyoto-u.ac.jp/ja/services/comp/contact.html

46 − Academic Center for Computing and Media Studies, Kyoto University京都大学 学術情報メディアセンター 全国共同利用版［広報］Vol.21 No.1（ Dec. 2022）− 47

京都大学学術情報メディアセンター全国共同利用版広報 Vol．21，No．1
２
０２２年
９日 発行
2022
年 １２月
12 月
xx 日発行
広報編集部会
編集者 京都大学学術情報メディアセンター

深沢 圭一郎（部会長）

全国共同利用版広報編集部会

鳥生 大祐（副部会長）

発行者 〒606-8501 京都市左京区吉田本町

中西 健史

京都大学学術情報メディアセンター當山 達也

Academic Center for Computing and Media Studies

熊谷 真由美Kyoto University
Tel. 075-753-7407
http://www.media.kyoto-u.ac.jp/
印刷所 〒616-8102 京都市右京区太秦森ヶ東町 21-10
株式会社エヌジーピー

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

表紙デザイン：中山 豊（中山商店）