東京工業大学 TSUBAME 共同利用
令和４年度利用終了課題 利用成果報告書集

東京工業大学 学術国際情報センター
共同利用推進室https://www.gsic.titech.ac.jp/tsubame本報告書集および個別の報告書の PDF ファイルは、以下の URL にあります。
令和４年度採択実績および利用終了課題報告書
https://www.gsic.titech.ac.jp/kyodou/kadai_r4■令和４年度産業利用 利用成果報告書 一覧
申請課題名
所属機関／利用課題責任者頁

3D 超音波 CT における腫瘍の検出と判別手法の開発
1
株式会社 Lily MedTech／陣内佑■令和４年度学術利用 利用成果報告書別紙の提出免除課題 一覧（利用成果を論文/学会等にて発表した要旨等の提出により提出免除）
申請課題名
所属機関／利用課題責任者
論文/学会等における発表済利用成果の情報
固体電解質セルによる CO2 電気分解反応の第一原理計算
立命館大学理工学部／渡部弘達[1] Hirotatsu Watanabe, Teppei Ogura, International Journal of Hydrogen Energy, Vol. 47, pp. 29027-29036,
2022, Pages 29027-29036(査読あり) DOI:10.1016/j.ijhydene.2022.06.243[2] Hirotatsu Watanabe, S. Higashidani, A. Tawa, I. Ban. Teppei Ogura, Oxidation process of Ni cathode in CO2electrolysis in SOEC: X-Ray and DFT study, ECS Transactions, Vol. 111, 2023, Pages 1319-1326
メタサーフェスレンズのマルチスケール解析
東京農工大学／岩見健太郎[1] Satoshi Ikezawa, Ryota Yamada, Kosuke Takaki, Chikara Ogawa, and Kentaro Iwami, “Micro-Optical Line
Generator Metalens for a Visible Wavelength Based on Octagonal Nanopillars Made of Single-Crystalline
Silicon”, IEEE Sensors Journal, 22 (15) 14851-14861 (2022) 査読有 DOI: 10.1109/JSEN.2022.3186060[2] Naoyuki Yamada, Hiroki Saito, Satoshi Ikezawa, and Kentaro Iwami, “Demonstration of a multicolor
metasurface holographic movie based on a cinematographic approach,” Optics Express, 30 (10), 17591–
17603, (2022) 査読有[3] ( 招 待講演 ) Kentaro Iwami, “Multicolor metasurface holographic movie based on a cinematographic
approach”, SPIE/COS Photonics Asia 108: PA 109: Holography, Diffractive Optics and Applications XI,
Beijing, 12318-47, 5-11 December 2022 (Online)[4] (招待講演) Kentaro Iwami, “Dielectric metasurfaces for holography and focusing at visible wavelengths”,
to be presented at the 13th International Conference on Optics-Photonics Design and Fabrication(ODF ’22),OThA4A-01, Sapporo, 3-5, Aug. 2022
データ同化援用フェーズフィールド計算による高精度多結晶組織予測
東京農工大学大学院工学研究院／三好英輔[1] E. Miyoshi, M. Ohno, Y. Shibuta, A. Yamanaka, T. Takaki, Validating a mean-field theory via large-scale
phase- field simulations for abnormal grain growth induced by nonuniform grain boundary properties, J.
Mater. Sci. 57 (2022/09) 16690-16709. [peer-reviewed] DOI:10.1007/s10853-022-07660-4[2] E. Miyoshi, M. Ohno, Y. Shibuta, A. Yamanaka, T. Takaki, "Abnormal Grain Growth Induced by Nonuniform
Grain Boundary Properties: Large-scale Phase-field Simulation vs. Mean-field Theory", WCCMAPXOM2022, Yokohama, Japan, July 15–Aug. 18, 2022.
[3] E. Miyoshi, M. Ohno, Y. Shibuta, A. Yamanaka, T. Takaki, "Combined Bayesian inference and phase-fieldmodelling for evaluating triple-junction drag on grain boundary migration", MMM2022, Baltimore, USA, Oct.
2–7, 2022.

申請課題名所属機関／利用課題責任者
論文/学会等における発表済利用成果の情報
分子動力学シミュレーションによる新興再興感染症に関する研究
国立感染症研究所病原体ゲノム解析研究センター／横山勝[1] Masuda, T.; Kotani, O.; Yokoyama, M.; Abe, Y.; Kawai, G.; Sato, H. Cis-Allosteric Regulation of HIV-1 Reverse
Transcriptase by Integrase. Viruses 2023, 15, 31. DOI:10.3390/v15010031

■令和４年度学術利用 利用成果報告書 一覧
申請課題名
所属機関／利用課題責任者

頁

MPI アプリケーションにおけるプロファイルおよびトレース予測手法の評価
5
電気通信大学／三輪忍
計算化学による薬物とトリメチル化シクロデキストリンの相互作用評価
15
千葉大学大学院薬学研究院製剤工学研究室／東顕二郎
アレスチン分子の活性化機構の分子動力学シミュレーション
17
東北大学大学院薬学研究科／井上飛鳥
新規 BRCAness 誘導阻害剤開発に向けた理論研究
19
名古屋大学大学院医学系研究科／小関準
揺動電磁気シミュレーションを用いた光放射現象に関する研究23
新潟大学工学部工学科／櫻井篤
臨床情報統合データベースの機械学習解析25
日本大学医学部生体機能医学系薬理学分野／浅井聰
海洋底探査を加速する AI とシミュレーション技術の開発
29
大阪公立大学大学院工学研究科／橋本博公
継続的ベイズ推論の改善
31
理化学研究所革新知能統合研究センター／モハマド エムティヤーズ カーン
ヒト血清アルブミン-環状ペプチド複合体の相互作用解析
35
長岡工業高等専門学校／和久井直樹
廃棄物最終処分場における間隙内流体挙動の数値解析
37
埼玉県環境科学国際センター／鈴木和将
LRnLA アルゴリズムを用いた物理シミュレーション
41
法政大学情報科学部／善甫康成
発達期多元自発活動の数理モデルとその学習理論の構築
45
日本大学医学部生体機能医学系薬理学分野／早川隆
深層学習を用いた分子動力学シミュレーションの高速化
49
名古屋大学大学院医学系研究科／林周斗

申請課題名所属機関／利用課題責任者

頁

電気コイル隙間内に流れ込む冷却液挙動に関する大規模数値解析
53
大阪公立大学工学研究科／金田昌之
深層学習を用いた分子動力学法ベースの水和熱力学量計算法の高速化
55
東北大学大学院工学研究科／吉留崇
粘性の温度依存性を考慮した乱流熱流動の直接数値解析
57
大阪公立大学工学研究科／須賀一彦
GPU クラスタを用いたミリ波帯大規模広帯域電波伝搬シミュレーション
61
国立研究開発法人情報通信研究機構／チャカロタイ ジェドヴィスノプ
ベイナイト変態中の準安定オーステナイトにおける固溶元素と炭素原子の相関の解明
67
茨城大学大学院理工学研究科／永野隆敏
人工画像を用いた大規模事前学習
69
産業技術総合研究所 デジタルアーキテクチャ研究センター／横田理央
人工画像データ及び事前学習モデルのリポジトリ化による AI ハブ構築
73
産業技術総合研究所 デジタルアーキテクチャ研究センター／谷村勇輔
集風レンズ付き風車とそのマルチロータシステムの流体シミュレーション
77
九州大学応用力学研究所／胡長洪■令和４年度学術利用 利用成果報告書（公開延期）
申請課題名
所属機関／利用課題責任者

頁

機能性ペプチド提示エクソソームの創製
公開延期

京都大学化学研究所／川口祥正

TSUBAME 共同利用 令和 4 年度 産業利用 成果報告書

利用課題名 3D 超音波 CT における腫瘍の検出と判別手法の開発
英文：Tumor Detection and Classification for 3D Ultrasound Computed Tomography
利用課題責任者
Yuu Jinnai
株式会社 Lily MedTech
Lily MedTech Inc.
https://www.lilymedtech.com/
邦文抄録（300 字程度）

リングエコーを読影するためには大きな三次元ボリュームのすべてを読影し、その中から小さな病変部位がな
いかを確認する必要があり、診断のための労力が大きい。この問題を解決するための手段として、自動診断支援 (CAD) ソフトウェアによってあらかじめ病変候補領域を提示することで読影負担を軽減する方法が考えられる。本研究の目的はリングエコー画像から腫瘍症例を提示する CAD ソフトウェアを開発し、その有用性を評価することである。各乳房画像中に腫瘍があるか否かを推定する機械学習モデルを開発し、ソフト単体による判定精度を ROC AUC で評価する。
英文抄録（100 words 程度）

Ultrasound Computed Tomography (USCT) is an emerging imaging technology that has received
significant attention for its promising capabilities for breast cancer screening. To reduce the effort
spent by radiologists we propose a deep learning model to predict the presence of breast tumors in a
USCT volume. We propose a convolutional neural network (CNN)-based classification model for
predicting if a USCT volume contains a tumor. We evaluated the model by an area under thereceiver operating characteristic curve (AUROC) in detecting USCT volumes containing a tumor.
The result suggests that the proposed model has a potential to reduce the workload of radiologists
in a clinical setting.
Keywords: Deep Learning, Computer-Aided Diagnosis, Ultrasound Computed Tomography, Machine
Learning

背景と目的

リュームを撮像する。一方、大きな三次元ボリュームマンモグラフィ (MG) は乳がんの診断のために最ものすべてを確認し、その中から小さな病変部位を検広く使われている画像診断装置である。しかしながら出するのは労力が増大する。この問題を解決するた

MG は乳腺組織が多い高濃度乳房の症例には感度めの手段として、自動診断支援 (CAD) ソフトウェア

が低くなってしまう問題が知られている。

によってあらかじめ病変候補領域を提示する方法

検査の見落としを防ぐ方法として、超音波 (US) を補助検査として併用することができる。しかし US の問題点として、技師のスキルの依存度が高いこと。
また、画像は技師によって撮像位置が決定されるため、同じ撮像位置で同じ画像を再現することは困難である。
超音波 CT (USCT) は技師のスキルに依存しない、再現性の高い三次元撮像を行う新しい超音波画像診断装置である (図 1)。USCT はコロナル断面で撮図 1が考えられる。

られた二次元断層画像のシリーズからなる三次元ボ1

本研究の目的は USCT の画像から腫瘍症例を提示

大値をその三次元ボリュームの腫瘍の有無の判定

する CAD ソフトウェアを開発し、有用性を評価するこスコアとした。また、Grad-CAM を用い、腫瘍の位置とである。USCT の画像中に腫瘍があるか否かを推推定を行い、腫瘍領域とのヒット率を評価した。

定する機械学習モデルを開発し、ソフト単体による精度を AUROC と AUPRC で評価する。
概要

本研究の目的は USCT のための CAD ソフトウェ
アを実現するため、画像から腫瘍症例を提示する
機械学習モデルを開発し、その性能を評価することである。データ並列化によって複数 GPU を利用した深層学習を行い、USCT の画像中の病変候補領域を提示する検出タスクと、提示された領域

図 2が悪性腫瘍であるかを判定する判別タスクを解決

学習を終えたモデルを用いて評価データセットでするためのモデルを学習させる。MRI・針生検な

ROC AUC をしたところ、断層画像の AUROC はどの診断で悪性と判定された症例を真値として、0.86 であった。

判別モデルの感度・特異度・AUROC を評価する。
結果および考察

２つのサイトから得られた 186 名の症例を対象に後ろ向き研究を行った。各乳房の画像は 24 から 91 枚の断層画像シリーズからなる三次元ボリュームであ図 3る。各症例は左右の乳房の画像があるため、計 372

図 3 は顕著性マップの出力結果である。正しく腫瘍個の三次元ボリュームからなる。各断層画像の画素が検出された断層画像の顕著性マップと腫瘍領域と数は 512x512 ピクセルである。各ボリュームは腫瘍のヒット率の平均は 0.844 であった。(a)の画像はの有無がラベル付けされており、医師による腫瘍領

USCT の画像である。画像上方に黒く見えているの域のアノテーションが付けられている。

が腫瘍領域である。(b)の画像は簡単なデータ拡張

実験には NVIDIA Tesla V100 GPU を用いた。実のみを用いて学習したモデルによるマップであり、(c)験プログラムはすべて Python 3.7 で実装した。機械の画像は CutMix によるデータ拡張を用いて学習し学習モデルは PyTorch を用いて実装した。

たモデルによるマップである。赤い領域が診断結果

図 2 は三次元ボリューム中に腫瘍の有無を判定するに最も影響を与えたと考えられる領域である。簡単

モデルの構成である。ImageNet で事前学習したなデータ拡張のみのモデルだと腫瘍とは関係のない

DenseNet-121 モデルを用いて各断層画像に対して領域が診断結果に影響を与えてしまっているが、よ腫瘍の有無を推定した。画像はモデルに合わせてり効率的なデータ拡張を用いることで腫瘍領域のみ

448x448 の解像度にサイズ変更した。ピクセルの輝が診断結果に大きな影響を与えている。データ拡張

度値は平均 0、標準偏差 1 に正規化した。モデルのを加えることによってモデルがよりノイズに対して頑過適合を防ぐために画像の反転・切り抜き・回転・ガ強になっていると考えられる。

ウスノイズ・CutMix のデータ拡張を行った。
各断層画像に対する腫瘍の有無の判定スコアの最2

まとめ、今後の課題

learning. Journal of Big Data, 2019.

本研究は USCT の腫瘍の自動検出のための機械

Ramprasaath R Selvaraju et al. Grad-cam:

学習モデルを評価した。モデルの精度は高く、自動Visual explanations from deep networks via

診断支援プログラムとして有効である可能性を示唆

gradient-based localization. In Proceedings ofしている。このプログラムを用いて病変のある症例の

the IEEE international conference on computer

候補を医師に提示することで見落とし防止と診断のvision, 2017.

効率化の効果が期待される。

Jianming

本研究では腫瘍の悪性・良性の区別なく、病変候補

attention by excitation backprop. International

領域として検出することを目的とした。しかしながら、

Journal of Computer Vision, 2018.

実際には提示された腫瘍の悪性度も診断支援には重要な情報となる。検出された領域が悪性である確率を推定する判別モデルを開発することで更に有益なプログラムになると考えられる。
参考文献

Katharina Holland et al. Quantification ofmasking risk in screening mammography withvolumetric breast density maps. Breast cancer
research and treatment, 2017.
Wendie A Berg et al. Combined screening withultrasoundandmammography

vsmammography alone in women at elevated riskof breast cancer. Jama, 2008.
RachelF Brem et al. Screening breast

ultrasound:

past,

present,

andfuture.

American Journal of Roentgenology, 2015.
Nebojsa Duric et al. Detection of breast cancer
with ultrasound tomography: First results withthe computed ultrasound risk evaluation (cure)prototype. 2007.
Jia Deng et al. Imagenet: A large-scale
hierarchical image database. In 2009 IEEE
conference on computer vision and patternrecognition, 2009.
Gao

Huang

et

al.

Densely

connected

convolutional networks. In Proceedings of theIEEE conference on computer vision andpattern recognition, 2017.
Connor Shorten and Taghi M Khoshgoftaar. A
survey on image data augmentation for deep

3

Zhang

et

al.

Top-downneural4

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

利用課題名 MPI アプリケーションにおけるプロファイルおよびトレース予測手法の評価
英文： Evaluation of a Method for Predicting Profiles and Traces of MPI Applications
利用課題責任者 三輪忍First name Surname Shinobu Miwa所属 電気通信大学
Affiliation The University of Electro-Communications
URL http://www.hpc.is.uec.ac.jp/miwa_lab/index.html
邦文抄録（300 字程度）
本研究では，小規模実行によって得た性能メトリクスを用いて大規模実行時の並列アプリケーションの性能メトリクスを予測するツールである Extra-P を関数コール回数予測に応用した場合の評価を行う．我々の実験結果によると，Extra-P は 5 つの並列アプリケーションに対して関数コール回数を 2.23%以下の誤差で予測できた．また，我々は，関数コール回数予測の応用として，関数コール予測の結果をも言いて関数の実行時間を予測する手法を提案する．我々の実験結果によると，提案手法はさまざまな関数の総実行時間の予測精度を平均 1.41%改善した．
英文抄録（100 words 程度）
In this study, we provide the first evaluation of Extra-P, which is a tool for predicting the performance
metrics of a parallel application executed at a large scale using the performance metrics of theapplication executed at multiple small scales, in a use case of function call count prediction. Our
experimental results demonstrated that Extra-P predicted the function call counts for five parallel
applications within an error of 2.23%. Additionally, we propose a technique to predict the total
execution time of a function using the call count predicted as a use case of function call countprediction. Our experimental results demonstrated that the proposed technique improved theprediction accuracy of the total execution time of various functions by 1.41% on average.

Keywords: Extra-P, parallel applications, function call count, evaluation

1. Introduction

In this context, various methods for modeling

The analysis of the dynamic behavior of parallel

the performance of parallel applications haveapplications

bybeen developed. These methods enable us tosupercomputing users for various purposes (e.g.,collect the various performance metrics (e.g.,

performance tuning). Because an application's

execution time) of a given application executed at

behavior depends on the scale of execution (i.e.,a large scale with reasonable computing time andprocess count and input data), the collection of

resources, by creating their scalability models

performance metrics (e.g., function call, andempirically or analytically. Specifically, Extra-P,

memory and cache access counts) usually requires

which

the execution of a targeted application combined

performance analysis of parallel applications,with a profiling tool (e.g., TAU) at a targeted scale.

automatically generates a model to provide the

However,theextrapolated performance from a subset of metricperformance metrics of large-scale applications

data of a given application executed at small

consumes a large amount of computing time andscales. In previous work, the authors reported

resources. A more cost-effective approach tothat Extra-P is very useful for estimating somecollecting performance metrics is required for

performance

metricslarge-scale applications.

numberfloating-point

this

is

widelyapproachperformedforcollecting

5

is

the

of

most

powerful(e.g.,

tool

for

the

execution

time,

instructions,andcommunication

given

targeted application. Such metrics include theapplication executed at a large scale; however, but

number of function calls and the total execution

it is still unknown whether Extra-P can be used

time of functions. To collect performance metrics,

for scalability prediction for the other metrics. In

an instrumentation code is inserted into theparticular, the function call count is an important

targeted application using a compiler or dynamic

metric for understanding the application behavior

library, and the application is then executed on aand is widely used for the performance analysis oftargeted system. Profiling tools, such as TAU andparallel applications; however, to the best of our

Score-P, are available for the insertion of theknowledge, researchers have not reported the

instrumentation

effectiveness of Extra-P in function call countperformance metrics collected through profilingprediction.

are measured values. We believe that the costs

In this paper, we show the first evaluation of(computing

Extra-P

prediction.

performance-metric collection through profiling

Therefore, as a use case for function call countare high because this approach requires aprediction, we propose a method that predicts thetargeted application to be executed on a targeted

total execution time of function using the number

system. Many techniques, such as selectiveof predicted function calls. Our experimental

instrumentation, have been developed to reduce

results demonstrated that Extra-P predicted thethenumber of function calls with sufficiently highhowever, performance is limited by the execution

accuracy for various functions and our proposedof the application plus the required collection

method predicted the total execution time ofoverhead.

various functions with higher accuracy than the

2.2. Extra-P

conventional method.

Extra-P is a performance analysis tool for

We summarize the main contributions of thisparallel applications that can be used for the

study as follows:

quick analysis of parallel applications. Extra-P







for

messagefunction

size)

call

of

count

a

costs

code.

time

of

We

and

note

that

the

resources)

performance-metricofcollection;

We demonstrate that Extra-P is also usefulextrapolates a performance metric of a targeted

for creating models of function call counts.

application executed at a large scale from those of

We demonstrate that combining function callthe application executed at multiple small scales.

count prediction is effective in predicting theTototal execution time of functions.

scalability model for the performance metric

We demonstrate that Extra-P greatly reducesusing regression. The independent variables of

the cost of collecting the total execution timethe model are the process count and problem size,of functions.

and the dependent variable is the performanceaccomplish

this,

Extra-P

generatesametric to be predicted.
2. BackgroundThe general form of the model used in Extra-P is2.1. Profilingas follows:

Profiling parallel applications is widely used for
performance analysis and tuning in the field of
where m is the number of independent variables,

high-performance computing. Profiling involves
the collection of the performance metrics of ais an independent variable, and

6

is a

regression parameter. n is the number of product

Table 1. System configuration ofterms, each of which has a different combination

TSUBAME3.0

of. The exponents

andare called the hypotheses and characterize theform of the equation. Extra-P users can specify
the sets ofTable 2. Node configuration of TSUBAME3.0(referred to as I and J,

andrespectively) used for the regression from thecommand line interface.
Using

Extra-P,wecancollect

variousperformance metrics of a targeted application

static code analysis. Thus, to obtain the exact

without executing the application at a targeted

number of function calls, we require a targeted

scale; however, the performance metrics reported

application to be executed at a targeted scale.

by Extra-P are predicted values. Additionally, the3.2. Experimental Methodologies

effectiveness of Extra-P has been proven for a few

In this subsection, we evaluate the accuracy of

performance metrics (e.g., total execution time ofExtra-Pa function and communication message size). It is

Specifically, we first collect sample data (i.e., callstill unknown whether Extra-P can also be usedcounts per function) by executing a targeted

for function call count prediction.

application at various small scales and then input

infunction

call

count

prediction.

them into Extra-P. Next, Extra-P automatically3.

Function Call Count Prediction using Extra-P

3.1. Function

CallCount andgenerates an optimal model to predict the

Performance

scalability of the call count for each function.

TuningFinally, we predict the function call count at a

The function call count is helpful for tuning

large scale with using the model generated andthe performance of an application. Because

compare the result with the measured call count.

calling a function has a performance overhead,

To enable the generation of multi-variable

a reduction in the function call count may

models, we use multi-parameter for the modeler

improve overall application performance. For

option in Extra-P. Moreover, we use I = {-1, 0, 1, 2,

example, inlining a function can remove the3} for the exponents of the polynomial and use J =calling overhead from an application, but it{0, 1} for the exponents of the logarithm inincreases

Equationtheapplication

code

size.

The

(1).

We

function call count is fundamental information

force_combination_exponentsinallow_negative_exponents to True.

performancetuning;therefore,

many

also

set
and

profiling tools support functionality to reportWethe call counts of the functions executed.

supercomputer TSUBAME3.0. It consists of 540

The function call count is a dynamic attributecompute nodes, each of which has two CPUsof an application. It is usually affected by(Intel Xeon E5-2680 V4). The system and nodevarious factors of running applications (e.g.,

configurations of TSUBAME3.0 are shown in

inputs, process count, and results of branch

Tables 1 and 2, respectively. We used TAU as a

instructions) and cannot be computed using

profiling tool when collecting the actual number

7

conducted

ourexperiment

on

the

Table 3. Benchmark programs and scales

Table 5. Fitting errors in function call countused for model generation. Abbreviations in

prediction

italics in the Scales column represent the as
following: pc=process count, ps=problem size,

gs=grid size, it=iterations, kc=key count, andmv=max value.

the scales listed in Tables 3 and 4. Our
experimental results exclude a few functions
executed only at a specific scale.
We used the following two metrics to evaluate
prediction accuracy.
MAPE: The accuracy of a model is often
evaluated using the mean absolute percent

Table 4. Scales used for prediction

error (MAPE), which is the average of theabsolute error rate between the model's
output values ( ) and the observed values ().

It can be expressed by the following formula:of function calls.
We tested five HPC applications shown on
TSUBAME3.0. Four out of the five applications

where N is the number of data points.

(i.e., EP, FT, IS, and MG), which we selected from

Weighted MAPE: In the MAPE, the errors in

the NAS Parallel Benchmark suite, are simple

each data point contribute to the overall error

and stress particular properties. The remaining

equally; however, the more precise prediction

application, LULESH, is more complex and weof a large function call count is helpful for

selected it from the CORAL-2 Benchmark suite.

performance tuning. Therefore, we introduce

The scales used for model generation are

the weighted MAPE as the other metric forsummarizedanprediction accuracy. The weighted MAPE is aapplication for every combination of the scaleweighted average of the absolute error rate byparameters shown in the table and then entered

the number of function calls and can be

all the resulting function call counts into Extra-P

expressed by the following formula:

in

Table

3.

We

executed

to generate the models. For example, we executed
EP at 30 scales (6 types of process counts

5types of problem sizes) in total.

where

Using the models generated by Extra-P, wefunction calls measured for the function

predicted the function call counts of each$t$ and the total number of function calls

application for the three types of process counts

measured during the execution of a targeted

and input data shown in Table 4. We note thatapplication at a given scale, respectively.

andrepresent the number ofthese targeted scales are larger than the data

3.3. Experimental Results

collection scales shown in Table 3. In the next

Table 5 summarizes the fitting errors of thesection, we focus on the functions executed at allfunction-call-count models generated by Extra-P.

8

(a) MAPE(b) MAPE(b) Weighted MAPE(b) Weighted MAPE

Figure 1. Accuracy of function call countFigure 2. Accuracy of function call countprediction for various process counts using

prediction for various input data using large

large input dataprocess counts

As shown in the table, Extra-P produced models

applications. We used large input data for thiswell -fitted to the collected data, in many cases. In

experiment.

particular,modelsFigure 1 (a) shows that two out of the fivegenerated by Extra-P had an error of 1.58% on

applications (i.e., IS, and LULESH) had small

average for FT. A few applications, such as EP

errors (within 10%). As shown in Table 5, Extra-P

and IS, had large MAPEs because they included

produced function-call-count models well -fitted tosome functions that had irregular patterns in thethese applications so that the generated models

function call count. However, the numbers of calls

had high prediction accuracy at the targeted

for such functions were relatively small. Extra-P

scales.

provided very accurate models for many functions

By contrast, FT and MG had large MAPEs

that were important to overall performance, evenbecause of the inaccuracy of the call-count models

in such applications.

for some functions. However, this inaccuracy may

Figure 1 shows the prediction accuracy of the

not become a problem in many use cases because

generated models for various process counts. The

such functions have relatively small call counts

x-axis represents the process count, and theand are therefore less important than the othery-axis represents MAPE or weighted MAPE. We

functions. Figure 1 (b) exemplifies this. The

note that the y-axis uses the logarithmic scale.

weighted MAPE was up to 2.23%.

The five bars on each x-label represent the five

Figure 2 shows the prediction accuracy of thethefunction-call-count9i-th total execution time of function F.
Generally, Extra-P predicts the total execution
time of function directly. Specifically, we firstcollect the total execution time of each function() while executing a targeted application at

various small scales. Next, we enter the collected
data into Extra-P and then create a total
execution time of function model for each function.

Figure 3. Overview of prediction of total

Finally, we predict the total execution time of

execution time of functions

function executed at a targeted scale using thegenerated models for various input data. We usedgenerated model (the direct method in Figure 3).

large process counts for this experiment. Similar

Hereafter, we call this the direct approach

to Figure 1, the weighted MAPEs of the models

because we need to distinguish between thisgenerated were very small for all combinations ofapproach and our approach proposed in the next

applications and input data, whereas the MAPEs

section.

were large for a few applications. Thus, we

4.2. Combined Approachconclude that Extra-P produced models sufficient

Theoretically,for practical use in function call count prediction.

function can be expressed as follows.

4.

Combining Function Call Count Prediction and

where

Total Execution Time of Functions Prediction

function per call.

the

total

execution

time

of

is the average execution time of

As described in the previous section, Extra-P

Because the execution scales (i.e., process counts

produced accurate function-call-count models, inand input data) have different influences on the

many cases. In this section, as a use case forcall count (predicting function call counts, we combinefunction per call ( ), the direct modeling offunction call count prediction with the scalability

difficult for some functions. For such functions,prediction of the total execution time of aprediction

function.

developing individual models for4.1. Predicting the Total Execution Time ofthen combining the predicted results using the) and total execution time ofaccuracy

can

be

is

improved
and

by

, and

models (the combined method in Figure 3). We

Functions
The prediction of the total execution time of acall this the combined approach and use Extra-P

function

to model both

is

often

performed

toidentify

a

and

.

scalability bug in a targeted application. As

4.3. Experimental Methodologies

shown in previous work, it is the main use case of

Many profiling tools report two types of total

Extra-P.

execution time of function. One is inclusive time,

The total execution time of function (which includes the time taken by all callees. The) can beexpressed as follows.

other is exclusive time, which represents the time
taken by the targeted function itself. In this study,
we performed the prediction for both inclusive

where

is the function call count andand exclusive time.

is the10

Table 6. Fitting errors in the total execution time of function prediction. In "x (y-z)," x is the
MAPE, and y and z are the minimum and maximum MAPE per function, respectively

We

used

the

same

experimental

systemdescribed in Section 3.2. We also used thebenchmark programs and execution scales listed
in Tables 3 and 4.
We evaluated the cost of data collection, in
addition to fitting errors and prediction accuracy.
We define the cost of collecting data ( C) for an
application (app) executed on M cores as follows:(c) MAPE

Because many supercomputing services require
users to buy points, which are consumed based on
the product of the computing time and number ofcomputing resources.
Both direct and combined approaches require
model generation and prediction processes in
addition to data collection. These two processes
also consume computing time and resources;(b) Weighted MAPE

however, they can be executed very quickly, evenFigure 4. Accuracy of prediction of total

on a single compute node. Because the costs ofexecution time of functions for various

model generation and prediction are relatively

process counts using large input data

small compared with the cost of data collection,
we ignore them in this study.

the direct approach, the combined approach

4.4. Experimental Results

reduced the fitting error by 20,302% for FT. ThisThe fitting errors for the total execution time ofis because, with respect to the exclusive time,function models are shown in Table 6. Each entry

Extra-P can create more precise models forrepresents the MAPE of all models generated byandan approach for an application, and the range

Figure 4 shows the accuracy in predicting therepresents the minimum and maximum values intotal execution time of function. The x-axis

the case of computing the MAPE for each function.

represents the process count, and the y-axis

The table shows that the combined approach

represents the MAPE or weighted MAPE. We

produced models with higher accuracy than thenote that the y-axis uses logarithmic scale. We

direct approach when predicting the exclusive

used large input data.

time of functions. In particular, compared withAs shown in Figure 4 (a), Extra-P on average

11

than

.

large process counts. Similar to Figure 4, theweighted MAPEs were small, whereas the
MAPEs were large. The models generated had aweighted MAPE of 1.41%, on average. Thus, weconclude that Extra-P produced models sufficientfor practical use in total execution time offunction prediction.
Figure 6 shows the cost of collecting the total
execution time of functions. The x-axes represent(d) MAPE

the number of execution scales to be collected,
while the y-axes represent costs. The y-axis uses
the logarithmic scale. The blue lines represent thecase of predicting the total execution time offunction and the red lines represent the case ofcollecting the total execution time of function
based on actual runs (i.e., using profilers). We
note that the data collection time for the direct(b) Weighted MAPE

approach was the same as that for the combined

Figure 5. Accuracy of prediction of total

approach.

execution time of functions for various input

The figure shows that the cost of predictively

data using large process counts

collecting the total execution time of functions

produced inaccurate models in the total execution

was completely constant across the number oftime of function prediction. The smallest MAPE

predicted scales. This is because we ignored

was 58.95+%, even when we used the combined

Extra-P's execution cost, which was relatively

approach. In particular, FT had a MAPE of

small compared with the cost of data collection, as

10,000,000+%

count.

described in Section 4.3. By contrast, the cost of

all-to-allcollecting the total execution time of function

communication functions and they were calledbased on actual runs increased gradually as thefrom several locations within the code, it was

numberdifficult for Extra-P to precisely extrapolate theparticular, in the case of collecting the total

performance of such functions.

execution time of function at 27 scales, theBy contrast, we observed different properties ofpredictive approach reduced the cost of data

the models produced by Extra-P in view of thecollection by 99% for EP compared with theweighted MAPE. Figure 4 (b) shows that theapproach of actually-running the application.

models achieved very high accuracy (within 1% ofThus, the predictive approach was effective,

the weighted MAPE) for four out of five

particularly for collecting the total execution time

applications. Therefore, we consider that theof function at multiple scales.

Because

FT

for

the

large

executedprocessmanyofpredicted

scalesincreased.

Inmodels achieved sufficient accuracy in practice.
Figure 5 shows the prediction accuracy of the

5.

generated models for various input data. We used

In this study, we evaluated Extra-P in function

12

Conclusions

Figure 6. Cost of prediction of total execution time of functions
call count prediction. Our experimental results

The predicted function call counts may be usefulshowed that Extra-P produced highly accurate

for the prediction of performance metrics, withscalability models for various functions in thethe exception of the total execution time ofcase of function call count prediction. Additionally,functions.

we showed that the predicted function call counts

In future work, we will extend the use case ofwerefunction call count prediction to the prediction of

helpfulfor

improvingtheprediction

accuracy of the total execution time of functions.

other metrics.

13

14

TSUBAME 共同利用 令和 4 年度 TSUBAME 利用 成果報告書

利用課題名 計算化学による薬物とトリメチル化シクロデキストリンの相互作用評価
英文：Evaluation of Interaction between Drug and Tri-Methyl Cyclodextrin
by Computational Chemistry
利用課題責任者 東 顕二郎
First name Surname Kenjirou Higashi
所属 千葉大学大学院薬学研究院製剤工学研究室
Affiliation Laboratory of Pharmaceutical Technology, Graduate School of Pharmaceutical Sciences,
Chiba University
URL http://www.p.chiba-u.jp/lab/seizai/index.html
邦文抄録（300 字程度）
難水溶性薬物の溶解性改善手法として、過飽和溶液の調製が広く研究されている。一方、過飽和溶液は熱力学的に不安定であるため、その安定化を目的として結晶化抑制剤が用いられている。シクロデキストリン類は薬物をその空洞内に包接することで、みかけ溶解度を改善する可溶化剤として知られている。近年我々はトリメチル化シクロデ
キストリンが可溶化作用に加えて、強い結晶化抑制作用を有することを見出した。本課題では、トリメチル化シクロデ
キストリンの結晶化抑制メカニズムを明らかにすることを目的として、計算化学（主に分子動力学シミュレーション）を用いて、過飽和溶液中の薬物とトリメチル化シクロデキストリンの相互作用を評価した。
英文抄録（100 words 程度）
The drug-supersaturation formulations are developed to overcome the solubility issue of poorly
water-soluble drugs. However, the supersaturated drugs are subject to crystallization due to their
thermodynamic instability. The addition of excipient that can inhibit drug crystallization and stabilize
drug supersaturation is of high importance. Cyclodextrin (CD) and its derivatives are widely employed
as solubilizers to increase drug apparent solubility through inclusion complex formation. Our group
found a new function of permethylated CD as effective crystallization inhibitor. In this study, we aimed
to investigate the crystallization inhibition mechanism of permethylated CD. The MD simulation was
conducted to unravel the interaction mode between drug and permethylated CD.

Keywords: シクロデキストリン、過飽和製剤、結晶化抑制作用、MD シミュレーション
背景と目的れている。一方、CD 誘導体による薬物の結晶化抑制

難溶性薬物の溶解性や経口吸収性改善技術として、作用に関する研究は未だ十分な検討が行われていな

薬物過飽和製剤が広く研究されている 。薬物過飽和製い。我々の先行研究において、12 種類の CD 誘導体に剤は、経口投与時に、消化管において薬物過飽和状態

よる薬物結晶化抑制作用を評価した結果、疎水性の高を形成することにより、難溶性薬物の経口吸収性を改いメチル化 CD がより強い薬物結晶化抑制作用を示す

善する。しかし、過飽和状態は熱力学的に不安定であことを見出した (Fig. 1)。しかしながら、その結晶化抑り、薬物の結晶化により、時間経過に伴い過飽和状態

制機構は不明であった。

は失われる。そのため、薬物結晶化を抑制または遅延させ、より長時間薬物過飽和状態を維持することが、過飽和製剤を利用した吸収性改善には重要となる。
近年の研究において、環状オリゴ糖であるシクロデ
キストリン(CD)の誘導体の一部が薬物結晶化抑制作
用を示すことが報告されている。CD は、生体適合性を有する添加剤として古くから医薬品、食品、化粧品など
様々な分野で応用されており、内部の疎水性空洞内にFig. 1. Comparison of solubilization and crystallization

薬物を包接することで可溶化作用を示すことが報告さ

inhibition ability of CD derivatives.

15そこで、本研究では、モデル薬物としてクロルタリドン

結果および考察(CLT)を用い、分子動力学(MD)シミュレーションによる

(a)

(b)

CD 誘導体の結晶化抑制メカニズム解明を試みた。

概要
MD シミュレーションによる薬物と CD 誘導体の相互作用評価
Fig. 4. Snapshots of (a) β-CD and (b) TM-β-CD with

6 分子の CLT(Fig.2a)及び、1 分子の β-CD またはトリ

CLT molecules from MD simulation at 50 ns.

メチル-β-CD (TM-β-CD)(Fig. b)を 81.0 Å3 の体積の立1.

方体セルに配置し、水分子を 18555 分子加え、エネル

可溶化作用

ギーを極小化したものを初期構造とした(Fig. 3)。全て

MD シミュレーションの結果、β-CD は CLT 分子を空の分子について第一原理計算を行い、各原子の電荷洞内に包接したのに対して、TM-β-CD は CLT 分子を空を補正した。AMBER14: EHT 力場を用い、25ºC (298

洞内に包接しなかった (Fig. 4)。実験的手法による可K)で 50 ns の平衡過程を設定し、NAMD プログラムに

溶化作用試験の結果では、β-CD は TM-β-CD と比較しより MD シミュレーションを実行した。

て強い CLT 可溶化作用を示しており、これは各 CD に
よる CLT 包接能の違いによることが明らかとなっている。

O(a)すなわち、MD シミュレーションの結果は実験的に求め

NH
HO

NH2
O
S
O

た可溶化作用の強さを支持しするものであり、CD 空洞内における薬物との相互作用が可溶化作用の主な駆

Cl

(b)

動力であることが確認された。

RO

2.
O
OR
O
OR

結晶化抑制作用MD シミュレーションの結果、β-CD 及び TM-β-CD は

β-CD: R=H
TM-β-CD: R=CH3

いずれも外側の官能基を介して CLT と相互作用を形成した(Fig. 4)。実験的手法による結晶化抑制作用試験のFig. 2. Chemical structure of (a) CLT and (b) β-CD or

結果では、TM-β-CD は β-CD よりも強い CLT 結晶化抑TM-β-CD.

制作用を示した。これは各 CD 空洞外の置換基と CLTの相互作用の強さの違いを反映しており、疎水性の高いメチル基が CLT と強く相互作用を形成するためと考えられている。すなわち、今回の MD シミュレーションの結果と実験的に求めた結晶化抑制作用の強さは一致しておらず、本 MD シミュレーション条件では各 CD による CLT 結晶化抑制作用を適切に評価できていないと考えられた。

CLT: Blue β-CD: Pink

CLT: Blue

まとめ、今後の課題

TM-β-CD: Pink

MD シミュレーションにより、各 CD の可溶化作用評

Fig. 3. Initial structures of CLT-supersaturated

価は達成できたが、結晶化作用評価は困難であった。

solution containing 1 molecule of β-CD (left) or

今後は、薬物結晶化機構のうち核成長に焦点を絞り、TM-β-CD (right). The water molecules are hidden for

MD シミュレーションを用いた検討を行う。具体的には、

the convenience of observation.

CLT 結晶界面への各 CD の吸着挙動を評価することで、各 CD の CLT 核成長抑制機構の解明を行う。

16

TSUBAME 共同利用 令和４年度 学術利用 成果報告書

利用課題名 アレスチン分子の活性化機構の分子動力学シミュレーション
英文：Molecular dynamics simulations of the arrestin activation mechanism
井上 飛鳥
Asuka Inoue
東北大学大学院薬学研究科
Graduate School of Pharmaceutical Sciences, Tohoku University
http://www.pharm.tohoku.ac.jp
細胞質シグナル制御タンパク質であるβアレスチンは G タンパク質共役型受容体 (GPCR) により活性化を受ける。この活性化にはリン酸化された GPCR の C 末端との結合が重要だと考えられてきたが、最近の研究によって GPCR の膜貫通ヘリックスから成る疎水性コアや膜脂質であるホスファチジルイノシトール 4,5 ビスリン酸 (PIP2) も寄与していることが分かってきた。そこで、リン酸化 C 末端、疎水性コア及び PIP2 の β アレスチンへの活性化の寄与を分子動力学シミュレーションによって検証し、β アレスチンの活性化機構を調べた。その結果、それぞれの結合構造では β ア
レスチンの活性化の程度が異なり、リン酸化 C 末端は細胞質側コアよりも活性化への寄与が大きいことが分かった。また、PIP2 によって β アレスチンは活性化構造に近づくことが分かった。
Activation of β-arrestins is triggered by binding to G protein-coupled receptors (GPCRs) on the plasma
membrane. Although binding to the phosphorylated C-terminal of GPCRs has been considered important for thisactivation, recent studies have revealed that β-arrestin activation is also influenced by the intracellular hydrophobic
core of the GPCR and phosphatidylinositol 4,5-bisphosphate (PIP2) in the plasma membrane. Here, by using
molecular dynamics simulation, we examined the contribution of the phosphorylated C-terminus, the intracellularcore, and PIP2 to β-arrestin activation. Our results show that the degree of β-arrestin activation differs for eachbinding mode, with the phosphorylated C-terminus contributing more to activation than the intracellular core. In
addition, we find that PIP2 molecules shift β-arrestin to the active conformation.
Keywords: β アレスチン、PIP2、GPCR、MD シミュレーション

背景と目的性化に対する寄与を調べることを目的とした。

β アレスチンは細胞膜上の G タンパク質共役型受容概要

体 (GPCR) に結合することで活性化する。活性化した

始めに、バソプレシン V2 受容体 (V2R) と β アレスチ

β アレスチンは GPCR の内在化を引き起こすことで G タンの複合体モデルをホモロジーモデリングによって構築

ンパク質シグナルを終結させるほかに、多様なシグナル

した（図１）。この複合体立体構造は GPCR のリン酸化 C

分子の活性化の足場となることが知られている。そのた末端と細胞質側コアの両方に結合している。モデリングめ、β アレスチンによる GPCR シグナル制御機構はの鋳型構造には V2R の C 末端を融合させたキメラ M2

GPCR を標的とした創薬において注目されている。これ

ムスカリン受容体 (M2R) と β アレスチンの複合体構造までの研究において、β アレスチンの活性化にはリン酸(PDBID: 6U1N) ならびに V2R と G タンパク質の複合体化された GPCR の C 末端との結合が重要だと考えられ

構造 (PDBID: 7KH0) を用いた。リン酸化 C 末端と細胞てきたが、最近の研究によって GPCR の細胞質側のコ質側コアによる活性化の寄与を調べるために、①リン酸

アや膜脂質であるホスファチジルイノシトール 4,5 ビス

化 C 末端と細胞質側コアの両方に結合した GPCR-βリン酸 (PIP2) によっても活性化が制御されることが分か

アレスチン複合体モデル、②リン酸化 C 末端のみと結ってきた。

合したモデル、③細胞質側コアのみと結合したモデル

本研究では、リン酸化 C 末端、細胞質側の疎水性コの３種を構築した。これらのモデルを 1-パルミトイル-2-

ア及び PIP2 による β アレスチンの構造変化を分子動力

オレイル-sn-グリセロ-3-ホスホコリン (POPC) 膜に埋学 (MD) シミュレーションによってとらえ、それぞれの活め込んだ（図２）。また、PIP2 による β アレスチンの活性17化の寄与を検証するために、PIP2 を含んだ POPC 膜に③のモデルを埋め込み、100 ns の MD シミュレーションを行った。
MD シミュレーション中における β アレスチンの活性化の評価には β アレスチンの 2 つのドメイン間のねじれの角度を用いた（図３）。N 末端側の N ドメインと C 末端側の C ドメイン間のねじれの角度が大きいほど活性化構造に近く、ねじれの角度が小さいと不活性化構造に近いことを示している。
結果および考察
POPC 膜におけるシミュレーションの結果、活性化を
図１ β アレスチンのモデル構築

示すドメイン間のねじれの角度は① > ② > ③ の順になっていた（図 4a）。この結果からリン酸化 C 末端と細胞質側コアが活性化に与える程度は異なり、リン酸化 C
末端のほうが細胞質側コアよりも活性化に対する寄与が大きいと考えられる。
また、PIP2 を含有した POPC 膜におけるシミュレーショ
ンでは、PIP2 含有 POPC 膜の β アレスチンは POPC のみの膜よりもドメイン間のねじれの角度が大きくなっていた（図 4b）。この結果から、細胞質側コアとの結合ではほ

図２ 脂質膜への埋め込み

とんど活性化を示さない β アレスチンについて、PIP2 が存在する脂質二重膜環境では活性化構造に近づくことが示唆された。
まとめ、今後の課題
今回の MD シミュレーションの結果から GPCR のリン酸
化 C 末端、細胞質側コア、PIP2 による β アレスチンの活性化構造変化への寄与が分かった。今後はこれらの活性化の程度の違いが β アレスチンのシグナル伝達分子との相互作用に及ぼす影響の解析を進める。
図３ β アレスチンのドメイン間のねじれ

図４ βアレスチンの活性化(a) POPC 膜における 3 つのモデルの比較(b) POPC 膜と PIP2 含有 POPC 膜の比較
18

TSUBAME 共同利用 令和４年度 学術利用 成果報告書

利用課題名 新規 BRCAness 誘導阻害剤開発に向けた理論研究
英文：Theoretical study for the development of novel BRCAness inhibitor
利用課題責任者
小関 準
所属

国立大学法人東海国立大学機構 名古屋大学大学院医学系研究科

邦文抄録（300 字程度）
近年、合成致死の原理を利⽤した PARP 阻害剤による BRCA 変異陽性がんの治療法が再脚光を浴びている。
PARP 阻害剤の有効性は様々な癌腫においても確認されているため、幅広くリバース・イノベーションを起こしうる可能性がある。そこで、未だ詳細が明らかではない「BRCAness」に関わる分子を標的とした、革新な新規誘導阻害剤の研究開発を行うことを目的とし、研究を推進した。CRISPR スクリーニングから得られた創薬標的候補（38 遺伝子）について創薬標的としてのランク付けを行い、ランクが高いものから順にタンパク質の構造学的および熱力学的解析を実施した。機能阻害剤となるリード化合物探索を行った。
英文抄録（100 words 程度）
Recently, PARP inhibitors based on the principle of synthetic lethality have received renewed
attention as a treatment for BRCA mutation-positive cancers, and their efficacy has been confirmed ina variety of carcinomas, leading to broad reverse innovation. 38 genes obtained by CRISPR screening
are drug discovery targets, ranked from highest to lowest in terms of protein structure and molecularbiology. Protein structures and thermodynamic analyses were performed to find candidate compounds
for functional inhibitors.
Keywords: PARP inhibitor, BRCAness, in silico drug design, thermodynamical analysis, structural
analysis

Alphafold2 により得られたタンパク質構造の構造学的背景と目的
近年、合成致死の原理を利用した PARP 阻害剤によ解析および熱力学的解析からタンパク質内の有効標的

る BRCA 変異陽性がんの治療法が、再注目されてき部位の同定と機能阻害のための候補化合物探索を実ている。国内だけではなく海外の臨床試験からも、様々施した。その結果、複数の標的に対するリード化合物をながん種において有効性が確認されてきている。2016

同定することができた。

年に行われた転移性前立腺がん患者の大規模解析において、PARP 阻害剤の標的となりうる DNA 修復不全

概要

“BRCAness”に関わる候補として 15 個の遺伝子に変

Domainfocused CRISPR スクリーニングから得ら異が同定された。しかしながら BRCA 変異陽性でも全くれた合成致死性を示す機能的ドメイン近傍のタンパク効果を認めない症例の存在や PARP 阻害剤に対する質構造を予測し、当該領域に対して独⾃の in silico

有意なレスポンスが認められない変異患者の存在が報

創薬デザイン戦略を適⽤することで、ドメイン活性特異告されている。そこでこれらの問題点を解決するために、的な新規阻害剤を設計した。申請者は Alphafold2 を

CRISPR スクリーニングから得られた“BRCAness”に用いて予測した様々な候補タンパク質の構造学的解析

関わる分子を標的とした、革新的な新規誘導阻害剤のと、分子動力学シミュレーションによる熱力学的解析を開発を目指した。

実施することで、独自の指標によって創薬標的の可能

本プロジェクトでは、実験的に“BRCAness”に関わる性をランク付けした。その結果を用いて、可能性の高い遺伝子を抽出することで従来の問題克服を試み、

タンパク質から順に阻害剤候補化合物の探索を実施し

19

た。

結果および考察CRISPR スクリーニングから得られた共通トップヒッ
トのうち、親株と耐性株に共通する 2 つのタンパク質、および親株のみでのヒットと耐性化株のみでのヒットした計 38 遺伝子ついて、インシリコ創薬標的に適しているかどうかを検討した。その結果、親株と耐性株の共通
ヒットと親株のみでのヒット因子についてはインシリコ設
図 2 in silico 創薬標的としてのランク付け計が容易ではないことが判明した。具体的には有効な阻害部位を絞り込むことが困難、もしくは特定の 2 次構造を有していないためコンピュータ設計が容易ではな

標的部位の同定のために、α-ヘリックスや β-シート等

いことが想定された。また特定の 2 次構造を有してはいの特定の 2 次構造を有している中で、構造に柔軟性をるものの、標的部位の候補の絞り込みが困難であるた有している箇所の有無を判定した。これは柔軟性があめ、Domain focused CIRPSR screening 等の標的部ることで標的タンパク質の界面が結合対象タンパク質位絞り込み実験から阻害部位が同定されない限り設計に対して Fit できる可能性を有するためである。一方で、が困難であるタンパク質パターンが存在することがわ

変性領域部位に関しては、PPI 結合に関与する可能性かった。

はあるが、コンピュータ上で取り扱う上で計算時間が発散するため（現在の計算限界）、インシリコ設計では除外している。この指標をもとにして、RANK-A に 1 分子、RANK-B に 9 分子抽出することができた。特に創薬標的としての可能性が高いタンパク質の例として、低分子化合物による阻害が見込めそうなもの（図 3 左）や、タ
ンパク質の構造的柔軟性と表面構造からタンパク質間相互作用阻害が見込めそうなもの（図 3 右）のようなものがピックアップされてた。現在までに複数のタンパク質に対する阻害剤候補化合物の抽出に成功している（図 4）。現在もなお、順次可能性が高そうなタンパク質に対して引き続きスクリーニングを実施している。

図 1 CRISPR スクリーニングから得られた候補分子のタンパク質構造パターン

そこで耐性株に特異的なヒット遺伝子群にフォーカスし詳細な解析を進めることとした。まず、mkDSSP を用
いてタンパク質全体に有する特定を有している割合、部分構造が集中しているのか否か、標的部位の同定もしくは推定できるかを指標として 6 つにランク付けした。
下記が創薬標的としてのランキングで用いた指標の詳細である（図 2）。
図 3 創薬標的としての可能性が高いタンパク質の例20

図 4 in silico 手法で抽出したリード化合物候補

まとめ、今後の課題
CRISPR スクリーニングから得られた創薬標的候補（38 遺伝子）について創薬標的のランク付けを行い、可能性が高いものから構造学・熱力学的解析を実施することで、機能阻害剤となるリード化合物探索を行った。
Domain-focused CRISPR screening と in silico
drug screening の統合解析から得られた候補化合物の検証を引き続き行い、PARP 阻害剤耐性細胞に対して強い BRCAness 誘導効果を示す化合物を in vitro、

in vivo の両面から選定する。さらに、その機序の解明に 向 け 、 Frataxin complex に 関 し て は 代 謝 、
JMJD6/BRD2 についてはエピゲノムの変化を阻害剤
およびノックアウト細胞株を用いて詳細に検証する。特に、「PARP 阻害剤耐性細胞は Frataxin の働きによって TCA 回路からエネルギーを得ることができ、それにより活性酸素種による DNA 損傷を抑止している」可能性や「PARP 阻害剤耐性細胞では DNA 修復反応に重要なクロマチンリモデリングが JMJD6/BRD2 を介して活性化している」可能性等の仮説について解析を進める。また、今回の in silico 解析では酵素阻害剤や PPIの設計のみ可能となっており、degrader については検討できていない。そのため、今後創薬 標的について
degrader の設計が可能かどうかを関連ベンチャー企業等との連携の上で検討を進めていく予定である。

21

22

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

揺動電磁気シミュレーションを用いた光放射現象に関する研究
A study of light emission phenomena using fluctuational electromagnetic simulation
櫻井 篤
Atsushi Sakurai
新潟大学工学部
Faculty of Engineering, Niigata University
http://www.eng.niigata-u.ac.jp/~rad/index.html
本研究では，半導体発光素子（light emitting diode ; LED）を用いた 熱フォトニクス(Thermophotonics ; TPX)
発電システムに注目した．これはエミッタに LED を用いることが特徴である．LED では，電圧を印加して非平衡状態にすることで，プランクの法則による黒体放射とは異なったエレクトロ・ルミネセンスと呼ばれる発光を生じる．そのため，遠方場においても高効率なエネルギー輸送ができる可能性がある．TPX システムは，大きな電力密度を有する発電手法である．また，バイアスのかけ方を工夫することで冷却機としても応用できる．
In this study, we focused on a thermophotonics (TPX) power generation system using semiconductor lightemitting diodes (LEDs). This is characterized by using LEDs as the emitter. In LEDs, by applying a voltage tobring them into an unbalanced state, light emission called electroluminescence, which is different from blackbody
radiation according to Planck's law, is generated. Therefore, there is a possibility that highly efficient energy
transportation can be performed even in a distant field. The TPX system is a power generation method with a large
power density. It can also be applied as a cooler by devising how to apply bias.
Keywords: thermophotonics, light emitting diodes, electroluminescence, power generation(Thermophotonics ; TPX)発電システムに注目した．

背景と目的
現在，世界規模での急速な工業化に伴い，エネルギ

こ れ は NTPV 発 電 シ ス テ ム と 異 な り ， エ ミ ッ タ にー消費の増加，化石燃料の枯渇や地球温暖化が大きLED を用いることが特徴である．LED では，電圧をな問題となっている．更に東日本大震災の発生などに印加して非平衡状態にすることで，プランクの法則によより，脱原子力エネルギーの潮流がある．このような現る黒体放射とは異なったエレクトロ・ルミネセンス(EL)と状から，エネルギー利用の高効率化や自然エネルギー呼ばれる発光を生じる．そのため，遠方場においてもによって生み出される再生可能エネルギーの利用に注

高効率なエネルギー輸送ができる可能性がある．TPX

目が集まっている．

システムは，大きな電力密度を有する発電手法である．

これらの問題の解決策の一つとして，近接場熱光起

また，バイアスのかけ方を工夫することで冷却機として電力発電（Near-field Thermophotovoltaic ; NTPV）

も応用できる．

システムがある．NTPV システムでは，太陽光だけでなく工場の排熱など多様な熱源が利用可能である．発結果および考察

電原理は，太陽光発電と大きな違いはないが，波長選TPX システムの LED-PV セル間のふく射熱伝達を

択的エミッタによって PV セルの高感度波長域に整合計算するために，揺動電磁気学シミュレーションを用いした熱ふく射光を放射することで高効率な発電ができる．

た．一般的な電磁波解析手法では，有限温度を持つ物

しかし，エミッタと PV セル間の距離をナノメートルスケ体から生じる電磁波を表現できない．揺動電磁気学のールで維持する必要があり，電力密度が小さいという理論では，物体内部でランダムに熱揺動する格子振動

問題点がある．

を揺動電流としてマクスウェル方程式の電流密度の項

そ こ で 本 研 究 で は ， 半 導 体 発 光 素 子 （ lightに導入される．従来の電磁波解析手法に揺動電磁気

emitting diode ; LED ） を 用 い た 熱 フ ォ ト ニ ク ス

学の理論を組み入れることによって，2 物体間の熱ふく23

いて TPX 発電システムのエネルギー輸送メカニズムに射解析が可能となる．
そこで本研究では，TPX 発電システムの LED と PVついて検討した．今後は異なる半導体材料の活用も視セルに用いる半導体材料として，黒リンを用いる．黒リ

野に入れて新たな光と物質の相互作用に起因するエネ

ンは中赤外領域にバンドギャップを持つ二次元材料でルギー輸送メカニズムの解明を目指す．

ある。図 1 に TPX 発電システムの概要を示す。TPX
発電では、半導体による非平衡熱ふく射を光起電力(PV)セルに照射し、光電変換を行うものである。PV セルで生成された電力のうち一部を非平衡熱ふく射エミッターの駆動電力として使用する。
実際、理想的な黒リン半導体が作製できたと仮定して TPX 発電システムに応用した場合の揺動電磁気学
シミュレーションを行った結果を以下に示す。高温側熱源（450 K）に接するエミッターと低温側熱源(300 K)に接する PV セル間におけるエネルギーフラックスに対して、化学ポテンシャルが 0～0.15 eV の範囲にある場合

Fig. 1 TPX システムの概要図のエネルギーフラックスの計算結果を図 2 に示す。比較のため、黒体のスペクトルエネルギーフラックスも示している。黒リンエミッターに印加する化学ポテンシャル
を大きくすると、バンドギャップより短い波長域の遠方場エネルギーフラックスが増加し、黒体限界を超えることが示された。温度が同じであれば、化学ポテンシャル
が正の光子は、化学ポテンシャルのない光子（熱放射）よりも光子状態あたりのエネルギーが高くなり、黒体放射を超えるエネルギーフラックスとなることがわかった。
次に、ギャップ距離の影響を調べた結果が図 3 である。通常、ギャップ距離 10 m 以上あれば、近接場効

Fig. 2 化学ポテンシャルによる遠方場光放射エネル

果の影響はほとんど無いと考えられるが、10 mm 離しギーの変化た場合でもエネルギーフラックスの大きさが維持される
ことが示されている。この結果はギャップ距離を何 cmと離しても結果は変わらないと考えられるため、かなり遠方場までエネルギーを伝達できる可能性がある。従来技術である、熱電発電や NTPV 発電と比較して、熱源から離れた状況で発電ができれば産業応用上の自由度が格段に上がる。ただし、エミッター・PV セル間の距離を離すとそれらの形態係数が関わってくるため、将来的にはこの非平衡熱ふく射エミッターに指向性を持たせることが必須である
Fig. 3 エネルギーフラックスの距離依存性

まとめ、今後の課題

本研究では，揺動電磁気学シミュレーションを用

24

TSUBAME 共同利用 令和４年度 学術利用 成果報告書

臨床情報統合データベースの機械学習解析
Machine Learning for Integrated Database of Clinical Information
浅井 聰
Satoshi Asai
日本大学 医学部 生体機能医学系 薬理学分野
Department of Pharmacology, School of Medicine, Nihon University
https://www.med.nihon-u.ac.jp/department/pharmacology/
臨床情報統合データベースに含まれる薬物投与・過去の診断や検査結果の履歴から疾患の発症・合併症進行リスクを予測するベイズモデルを構築しようとした。この目的のため、ベイズモデルにおいて少数のリスク因子を非線形効果まで含めて精密に評価するために、アニールド重点標本可逆ジャンプリーマン多様体上ハミルトンモンテカルロ
法の大規模実装を開発した。効率化の工夫により従来の可逆ジャンプモンテカルロ法と異なり、大規模データのベイ
ズマルチカーネル学習に実用可能なものが得られた。さらにこれまで開発した手法を実際の臨床データに適用し、既存の医療統計手法では検出できなかったバイアス効果を明らかにする等の結果を得た。
We developed an efficient large-scale implementation of Bayesian multiple-kernel learning based onannealed-importance-sampling-Riemannian-manifold-Hamiltonian Monte Carlo, with the goal of
estimating risks of disease onset and progression based on clinical information recorded in a large
database of electronic medical records. Developing and applying different techniques for efficientcomputation, we obtained an implementation of the algorithm that can deal with a Bayesian
multiple-kernel model with large numbers of kernels and samples, unlike previously reported
reversible-jump algorithms.

Keywords: electronic medical records, pharmacoepidemiology, Bayesian analysis, multiple-kernel
learning, Monte-Carlo method
背景と目的

実際に疾患進行に影響する変数は少数であり、変数選

本研究の背景及び目的は昨年度の我々の同名の課択を効率的に行う必要がある。また臨床医が判断する題と概ね同じであり、まずこれを以下に引用する。（引重症度のような潜在変数の情報をアルゴリズムによっ

用開始）近年、電子カルテデータなどを含む大規模なて推定できることも必要だ。これらを変数間の非線形関

医療データベースから薬剤の効果・疾患リスクなどの係も考慮しながら行うことは，潜在変数付きノンパラメト

医学的知見を機械学習・人工知能技術を用いて抽出すリックベイズモデルでスパースな変数選択を行う問題にる試みが注目を集めている。実際に、世界各国の大学帰着される(次ページ図参照)。上記のような潜在変数

病院や地域中核病院の電子カルテデータを用いた解付きのノンパラメトリックベイズの枠組みにおいて、医療析が行われはじめており[1,2]、さらに適切な倫理的なデータベースのような大規模データ上でスパースかつ枠組みのもと、電子カルテデータに患者から採取した

効率的な推定に成功した研究はまだ存在しない。先行遺伝•生化学的測定データを組み合わせ、これを患者研究のほとんどで潜在変数のモデリングは巧妙に避け

個別医療に役立てようとする試みがはじまっている。

られている。実際、医療データに限らずノンパラメトリッ

しかしながら、近年の機械学習•人工知能技術の発

クベイズモデルにおける潜在変数の推定やスパース変展にもかかわらず、上記のような大規模医療データか数選択は、機械学習一般の問題として依然として困難

ら有益な医学的知見を抽出するためには解決すべきな問題の一つであり、効率的な解法が模索されている

技術的課題が複数存在する。例えば、医療データベー[3]。申請者のグループは、下にも述べるように、この問スの解析では薬剤投与・検査値推移・疾患発症など数題を解決しうる新奇なアルゴリズムを得たので、本課題

千以上の予測変数から疾患進行を予測しようとするが、において TSUBAME を利用した大規模実装を行いこ25

の手法の有用性を示そうとする。（引用終了）

た。

本年度では上記に加えて、これまでに開発した機械

また、昨年度開発した、非ベイズ設定におけるマル

学習アルゴリズムを具体的にデータベースの解析に適

チカーネル法を実際にデータベースに適用し、新たな用していくという目的を設定した。そして、具体的に解析

医学的知見を得るとともに、開発中のアルゴリズムを実を行ったところ、従来の医療統計の手法の限界を超え際にデータベースに適用するためのインターフェースをて新たな知見を得られたため、その内容を以下にまと開発した。具体的には、近年論争の題材となっていた、める。また、新奇アルゴリズムの開発についても一定のベンゾジアゼピン系薬物・抗コリン作用を持つ薬物の長成果が得られ、論文報告の準備中であり、内容を以下期投与による認知症リスクに関する従来の研究結果ににまとめる。

多分にバイアスが含まれている可能性が高いことを機械学習手法を用いて明らかにした。

概要
臨床情報統合データベースから医学的知見を引き結果および考察出すための新奇アルゴリズムの大規模実装が部分的

ベイズマルチカーネル法の大規模実装に完成した。昨年度から引き続きガウス過程に基づい実装すべきアルゴリズムは下図 A のように、薬物投たノンパラメトリックベイズの大規模実装を可能にする与量に対して疾患・合併症発症などのような予測すべく、高効率のモンテカルロ法を開発した。リーマン多べき対象のリスク関数を設定し、それらを統計的に様体上のハミルトンモンテカルロ法とアニールド重点サ

推定することである。このとき、疾患発症の予測に役ンプリングとよばれる手法を用いつつ独自の実装上の立つ変数の数は、比較的少数と考えられ、これらの効率化を施すことで高効率のベイズマルチカーネル法

少数の変数が何であるかを見抜きたい。また、疾患の事後過程サンプリングアルゴリズムの開発に成功しリスクは薬物投与量などに対して非線形に変化する

26

ことが通常であり、非線形推定の手法を用いたい。

体 上 の ハ ミ ル ト ン モ ン テ カ ル ロ 法 (RMHMC) は将来的に潜在変数を導入することからベイズの枠組

Girolami らにより 2011 年に提案されているが、そみでアルゴリズムを構築することを考慮すると、ベイの実行のためには、各モンテカルロ試行で毎回修正

ズマルチカーネル法が適しているという結論になっヘッセ行列の全ての 3 回偏微分や逆行列を計算すた。ベイズマルチカーネル法はこれまでも理論的にる必要があり、計算コスト上の難点があった。しかしは論じられてきたが[4]、技術的な困難から実用上はながら我々は、修正ヘッセ行列の 3 回微分をもとのあまり用いられてこなかった。議論は省略するもの

ヘッセ行列の微分値と特異値分解を用いて高速にの、いくつかのアプローチの中で我々はノンパラメト

計算しつつ、前状態のヘッセ行列の特異値分解の結リックベイズの枠組みで変数選択をする上で根本的果を出発点として現在のヘッセ行列の特異値分解にに問題を解決する可能性が最も高いと思われる可2 次収束する Jacobi 法をベースとするアルゴリズム

逆ジャンプモンテカルロ法を用いた。この方法はを用いて計算することでこの計算量を大幅に低減さ1995 年に Green らによって提唱されて以来[5]、用せ、実用に足るアルゴリズムを得た。

いる変数・カーネルの集合を変更する際の可逆ジャ

このアルゴリズムを、機械学習の標準ベンチマーク

ンプの採択率の悪さが問題とされ、様々な改良案がデータである UCI データセットに適用した結果、及び提案されていた。我々は、図 B に示すように、あるリ

以前我々が用いた日本大学医学部における認知症

スク関数の集合から、別のリスク関数の集合にジャリスク解析のためのデータセットに適用した結果と合ンプする際に、それらのリスク関数集合に対応する

わせて論文公表を準備中である。

モデルを人工的に内挿した中間モデルを多数用意しこれらの間を可逆的にジャンプすると、計算時間を

非ベイズマルチカーネル法の大規模実装の日本大増やせば増やすほど理想的な採択率に近づくという学医学部臨床情報統合データベースへの応用

アルゴリズムを定式化したが、調査の結果これはア

令和３年度から令和４年度初頭に我々はベイズマル

ニールド重点標本可逆ジャンプモンテカルロ法

チカーネル法の実装に先立って非ベイズのマルチカ(AIS-RJMC)として 2013 年に一般的•抽象的な形でーネル法の TSUBAME3.0 における複数 GPU を利Karagiannis と Andrieu により提案されている方法

用した大規模実装を構築していた。そこで、日本大[6]であった。しかしながら、この方法は実応用上で学医学部倫理委員会の承認のもと、この大規模実の高性能にまだ結びついていないためか、あまり応装の臨床情報統合データベースへの適用を開始し用例が報告されていなかった。なぜ応用上の高性能た。特に令和 4 年度後半は、まず認知症リスクの解に結びついていないかは明らかで、理想的な採択率析に注力した。認知症のリスクは多数知られているの実現のためには可逆ジャンプの一回一回においが、特に長期投与されている薬物によるリスクが論てモンテカルロ法を収束させたまま一つのモデルか

争を引き起こしていた。従来法による医療統計解析

ら別のモデルにゆっくりと移行させる必要があったが、によりベンゾジアゼピン系薬物や抗コリン作用を持

通常の AIS-RJMC の実装ではそのようなことが容つ薬物の長期投与が認知症リスクを高めるとする複易でなかったからである(図 B)。しかしながら我々は、数の研究結果が報告されていた一方、これは統計

ベイズマルチカーネルの問題設定において応用上解析のデザインに起因する逆因果バイアス・交絡バ十分な精度でこれを達成することに成功した。具体

イアスではないかという疑いを唱える論文も発表され的には、モンテカルロ法としてスケール不変性を持ていた。そこで、我々は逆因果やその他の交絡を表つエネルギーの修正ヘッセ行列(ヘッセ行列の負固す予測変数・カーネルをこれらの薬物の効果を表す

有値を正に固有値に置き換えたもの)を計量とするリ予測変数・カーネルと同時に用いて統計解析を行っーマン多様体上のハミルトンモンテカルロ法を組みたところ、逆因果・その他の交絡の影響が明確に示込んだ AIS-RJ-RMHMC を設計した。リーマン多様される結果となった。そしてこれらのバイアスの調整

27

下では(少なくとも我々のデータセットからは)これらの薬物の長期投与による認知症発症リスクの上昇は有意でないという結果を得た。ベンゾジアゼピン系
薬物の効果の解析結果については論文投稿後、出版社から小修正後に受理する連絡を受けており、抗コリン作用を持つ薬物の影響については論文投稿準備中である。
まとめ、今後の課題完成した大規模実装によって効率的な変数選択を行いながら疾患リスクや薬物効果を判定できるようになった。しかしながら、潜在変数を組み込んだモデルのための実装はまだ行っていない。臨床情報統合データベースのデータのような時系列データの潜在変数の推定のためには、祖先サンプリング付き粒子ギブス法という確立した手法があり、これを上記のAIS-RJ-RMHMC に組み込むことで達成できると考えられるが、データが大規模であるだけにさらなる入念な調整が必要と思われる。令和 5 年度はこの点においてアルゴリズムを向上させつつ大規模応用の応用範囲を広げていくことを展望としている。
参考文献[1]S. Lee and H.-S. Kim J Lipid, Atheroscler.
(2021) 10(3):282-290[2]M. Chowdhury et al., Front. Psychiatry(2021) 738466[3]M. Gönen, Proceedings of ICML (2012)[4]T. Suzuki, CoLT 2012, pp.8.1-8.20[5]P. Green, Biometrika (1995) 82(4):711-732[6]G. Karagiannis & C. Andrieu, J. Comp.
Graph. Stat. (2013) 22(3): 623-648[7]M. Girolami & B. Calderhead, J. R. Stat. Soc.
B (2011) 73(2): 123-214

28

TSUBAME 共同利用 令和４年度 学術利用 成果報告書

利用課題名 海洋底探査を加速する AI とシミュレーション技術の開発
英文：Development of AI and simulation technologies to accelerate ocean-bottom exploration
橋本 博公
Hirotada Hashimoto
大阪公立大学 大学院工学研究科
Osaka Metropolitan University
https://kyoiku-kenkyudb.omu.ac.jp/html/100002465_ja.html
海溝型地震や熱水活動観測による資源評価など，調査船を用いた海洋底探査の重要性が増している。外乱が複雑に変化する状況では，緻密な操船制御を長時間にわたり求められるため，その代替として深層強化学
習にもとづく自律操船 AI の開発を行った。また，海底での採掘・揚鉱装置の開発を加速させるために，個別要素法と MPS 法を連成させた混相流解析用の GPU シミュレーションコードを開発し，その検証を行った。
The importance of seabed exploration using research vessels is increasing. A careful and precise
ship control is required in a situation where natural disturbances complicatedly change in time. Asan alternative, we developed AI for autonomous ship maneuvering based on deep reinforcement
learning. In addition, in order to accelerate the development of mining and mining equipment, a
GPGPU DEM-MPS code for simulating multi-phase flows reproducing seabed environment wasdeveloped.

Keywords: Seabed exploration, AI, DEM,-MPS, GPU

背景と目的

成シミュレーションが求められる。本研究では，個別要

海溝型地震や海底火山噴火の予測，熱水活動観測

素法と粒子法を連成させて解く固液混相流シミュレーシ

による資源評価など，調査船を用いた海洋底探査の重ョンを開発し，スパコンを用いた大規模数値解析が実要性が増している。気象海象が複雑に変化する状況下行可能な環境の構築を図る。

において，海洋底探査に課せられた高度ミッションを達成するためには，革新的な船舶制御が求められる。本

概要

研究では，既存の制御理論では取り扱いが困難な海

今年度は，深層強化学習のひとつである DDPG ベ

洋底探査における操船問題の新たな解として，深層強ースの自動航路追従 AI の開発を行った。ニューラルネ

化学習にもとづく自律操船 AI の開発を行う。スパコン

ットワークの入力には予め与えられる航路情報と船位・を利用した学習環境を構築し，学習スキームの改良や方位の偏差を用いる。追従精度の検証とロバスト性を学習用パラメータの最適化を行うことで，海洋底探査の確認するため，海洋調査船を用いて航路追従に関する

高度化と効率化に資する船舶制御を実現する。開発し検証実験を行った。その結果，海技者と同等の精度にた自律操船 AI は，模型船を用いた水槽試験による検て追従が可能であるものの，操船指示の左右非対称

証に加えて，実船を用いた実海域での実証実験を実施

性，指令値の振動，大出力値の頻発という課題が確認することにより精度とロバスト性の評価を行う。これらのされた。そこで，これらの課題を解決するための提案を

研究成果をもとに，今後の海洋底探査を加速させるた行い，改めて実船実験を実施することでその有効性をめの基盤的な船舶制御技術として確立することを目的

実証した。

とする。

海底環境の模擬シミュレーションについては，自前でまた，探査・開発の対象となる海底環境は，砂，泥，開発した個別要素法と MPS 法を組み合わせた固液混

礫などの離散体で構成されており，採掘・揚鉱装置の相流解析用の GPGPU コードを開発し，ダムブレイク問開発を加速させるためには，連続体である海水との連題でシミュレーションモデルの妥当性を検証した。

29

結果および考察
開発した自動操船 AI について，海洋調査船を用いて航路追従を行った結果を Fig.1 に示す。自然外乱が存在する実海域においても，目標とする直線と円弧を組み合わせた航路を高精度に追従可能であること，オ
リジナルの DDPG の課題であった指令値の振動や大出力値を回避できることを確認した。この結果から，開発した AI は人間と同程度の操船能力を有するものと判断できる。

Fig.2 Multiphase flow analysis of dam-breaking

液相，固相それぞれの単相計算に比較して，液相はダムブレイク直後に固相粒子に引っかかることで進行が遅れていることや,それによって速度が低下したことで水しぶきの高さが低くなることが観測された。その後，アーチ状の巻波が形成されている様子も見られた。これは，固体粒子によって形成された斜面との摩擦で流体が減速し，水面付近との速度差が発生することによるものと考えられる。
Fig.1

Ship path-tracking by AI captain(upper: ship trajectory, lower: rudder command)

まとめ，今後の課題
今年度は，深層強化学習ベースの自動航路追従 AI開発した DEM-MPS の GPGPU 連成解析コードのの開発と GPU シミュレーション用の DEM-MPS 連成

妥当性を検証するため，ダム崩壊解析を行った。その解析コードの開発を行った。前者については，実船を用結果を Fig.2 に示す。

いてバリデーションを実施した。オリジナルの DDPG では色々な課題が確認されたが，それらの解決策を提案
し，その効果として，開発した自動操船 AI は自然外乱
下でも航路を高精度に追従しつつ，操舵の量と頻度を抑制できることを確認した。今後は，より外乱の厳しい条件での検証が必要となる。
DEM-MPS の混相流計算コードについては，ダムブ
レイク解析により連成影響が考慮されていることが確認できた。今後は検証実験を実施することで精度を明らかにするとともに，大規模計算への対応を進めていく。

30

TSUBAME 共同利用 令和年度 学術利用 成果報告書

利用課題名 22IAJ - Improving Sequential Bayesian Inference
利用課題責任者
Mohammad Emtiyaz Khan所属
Approximate Bayesian Inference Team, RIKEN Center for Advanced Intelligence Project
http://team-approx-bayes.github.io
Bayesian principles offer a mathematically rigorous way to fundamentally improve and rethinkcontinual learning and adaptive AI systems. However, performing approximate Bayesian inference
for deep learning and large neural networks is still a challenging open question and requires large
amounts of compute resources. In this project, we employed TSUBAME3.0 to develop newstate-of-the-art Bayesian deep learning methods. These lay the foundation for our future work on
developing robust practical AI systems which can continually learn and update their beliefs byperforming sequential Bayesian inference.

Keywords: Approximate Bayesian Inference, Continual Learning, Deep Learning, Optimization,
Numerical Algorithms

背景と目的

Most of modern deep learning uses optimizers
to train deep neural networks. The following
equations show the update-rule of the popular
RMSprop / Adam deep learning optimizer:

Current machine-learning and AI systems arevery rigid and struggle to adapt to newsituations. Sequential Bayesian inference
offers a principled way to update the current
beliefs based on newly observed information.
However, performing exact Bayesian inference
for large neural network models is
computationally intractable. The goal of thisproject was to develop new state-of-the-art
methods for approximate Bayesian inference.
These algorithms are expected to enable AI
systems
tobecomemorerobusttoperturbations and continually adapt likehumans and animals.

To enable deep learning systems to adapt to
new situations and continually learn, it isessential for these models to know what theyknow. However, these models should also
understand what they don’t know yet andshould not be overconfident in such situations.
Starting from a rigorous mathematical theory,
we have developed a modification of RMSpropwhich can estimate uncertainty in deeplearning by computing an approximate
Bayesian posterior. The algorithm is shown in
the following equations, where the changes to
RMSprop/Adam are highlighted in red.

概要

The report is structured in three sections. Eachsection corresponds to an algorithm forBayesian deep learning we researched anddeveloped using TSUBAME3.0.
1. Second-Order Optimization via Bayes(SOBA)
2. The Lie-Group Bayesian Learning Rule
(LieGroup-BLR)
3. Sharpness-Aware Minimization as an
Optimal Relaxation of Bayes (bSAM)
結果および考察

1. Second-Order(SOBA)

Optimization

viaThis method allows the deep learning model tobe uncertain about its decisions in regions

Bayes

31

where there is no data, as shown in thefollowing Figure 1.

Figure 2. In approximate Bayesian inference, aprobability distribution (called the Bayesianposterior) is estimated. SOBA is restricted toGaussian distributions, but using the Lie-groupformalism, we were able to generalize to other
distributions. The figure shows Uniform,
Rayleigh and Laplace distributions.

Figure 1. When training neural networks, using
our SOBA optimizer instead of RMSprop allows
the model to be uncertain in situations where
there is little observed data. This uncertainty isessential to quickly adapt to new beliefs once
more information is observed.

In Figure 3, we compare the LieGroup-BLR to
SGD and the SOBA optimizer (here called“iVON”).
We observe that the additional
flexibility to choose a non-Gaussian noise
distribution can sometimes improve the results.

The SOBA optimizer offers state-of-the-art
performance for Bayesian deep learning. A
previous version was part of our winning solution
to the NeurIPS 2021 Challenge on Approximate
Inferencein
Bayesian
Deep
Learning(https://izmailovpavel.github.io/neurips_bdl_comp
etition/).
TSUBAME3.0 was essential in developing andbenchmarking this method. Large computational
resources are required when training neuralnetworks on big datasets such as ImageNet or
simply to train modern neural networkarchitectures. In future projects on sequential
Bayesian inference for large neural networks,
SOBA will be one of the key ingredients to obtainaccurate models with well-calibrated uncertainty.

Figure3.
Thedeveloped
Lie-Group
optimization methods improve the performanceover traditional SGD optimization in terms ofaccuracy and uncertainty on image classification
datasets with neural network models.

2. The Lie-Group Bayesian Learning Rule
(LieGroup-BLR)
Due to its mathematical formulation, the SOBA
optimizer appears to be restricted to perturb
the weights with Gaussian noise (see the
Algorithm before Figure 1). In a follow-up work,
we found a way to generalize SOBA fromGaussian noise to arbitrary noise distributions.
This is achieved by using the mathematicalframework of Lie groups (see the following
Figure 2).

3. Sharpness-Aware Minimization
Optimal Relaxation of Bayesasan

Sharpness-Aware Minimization (SAM) is a
recent method proposed by Google researchers to
improve generalization and robustness in deeplearning. The method has been very influential
over the last years, accumulating over 500
citations in a short timespan.

32まとめ、今後の課題

In this project, we have developed newstate-of-the-art methods for approximate
Bayesian inference (SOBA, LieGroup-BLR andbSAM). Our next step is to use thesealgorithms to improve sequential and continual
learning of AI systems, making them more
adaptive and robust. Furthermore, future effort
will be spent on applying the approximate
Bayesian inference algorithms developed inthis project on larger problems and machinelearning models. For instance, we want to showthat our methods can be useful in large-scalesettings like large language models.

Figure 4. Sharpness-Aware Minimizationconsiders the worst-case perturbation in aneighborhood, whereas Bayesian approaches
average over many random perturbations.
In our recent work we discovered a connectionof SAM to Bayesian approaches, as illustrated in
the above Figure 4. This connection enabled us topropose a Bayesian-SAM method (called bSAM)
which extends SAM to obtain uncertainty
estimates.
We illustrate the uncertainty estimation in
Figure 5. The decision boundary obtained by
bSAM is more “blurred out”, especially in regions
where there is no data. In these regions, themodel is uncertain about the correct prediction.

Figure 5. bSAM extends SAM to obtain
uncertainty estimates without computational
overhead.

Figure 6. The proposed bSAM method improves
accuracy and uncertainty estimation whencompared to SAM.
Using TSUBAME3.0 was essential indeveloping and researching our bSAM method onlarge deep learning models.

33

34

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

利用課題名 ヒト血清アルブミン-環状ペプチド複合体の相互作用解析
英文：Interaction analysis of human serum albumin-cyclic peptide complex
利用課題責任者 和久井 直樹
First name Surname Naoki Wakui
所属 独立行政法人国立高等専門学校機構 長岡工業高等専門学校
Affiliation National Institute of Technology (KOSEN), Nagaoka College
URL http://www.nagaoka-ct.ac.jp/
邦文抄録（300 字程度）
医薬品として近年注目を集めている環状ペプチドであるが、腎排出によって急速に体外へと排出されてしまう
課題を抱えている。この課題を克服するためには輸送体であるヒト血清アルブミンとの結合が重要である。環状ペプチドとアルブミンの複合体構造の相互作用を解析することで、アルブミンとの結合に重要な相互作用を明らかにし、環状ペプチド設計に役立てることが出来る。結晶構造解析によって得られたヒト血清アルブミンとコリスチンの複合体構造に対して分子動力学シミュレーションと相互作用解析を実施した。相互作用解析の結果、コリスチンの直鎖構造および環構造のどちらにおいても主鎖を構成する原子との水素結合が主であることが明らかとなった。
英文抄録（100 words 程度）
Cyclic peptides have attracted much attention in recent years as pharmaceuticals, but they face
the problem of being rapidly eliminated from the body by renal efflux. In order to overcome thisproblem, binding to the transporter, human serum albumin, is important. Analysis of the interaction
between cyclic peptides and albumin in complex structures will reveal the interactions that areimportant for binding to albumin and will be useful for cyclic peptide design. Molecular dynamics
simulations and interaction analysis were performed on the complex structure of human serum
albumin and colistin. Interaction analysis revealed that hydrogen bonding with atoms

constituting the main chain is the predominant mechanism in both the linear and ring
structures of colistin.

Keywords: Cyclic peptides, human serum albumin, molecular dynamics simulation

概要

背景と目的
低分子医薬品、抗体医薬品に続く第 3 の医薬品とし

結晶構造解析によって得られた、ヒト血清アルブミンて環状ペプチド医薬品が注目を集めている。しかし、環とコリスチンの複合体構造を初期構造とし、1 µs の分子状ペプチドは腎臓でろ過され、体外へと急速に排出さ動力学シミュレーションを 8 回実施した。また、ヒト血清れてしまうという課題がある。この課題を克服するため

アルブミンとコリスチンの複合体構造からコリスチンをには血漿中に存在する輸送体であるヒト血清アルブミ

削除したヒト血清アルブミンに対しても 1 µs の分子動力

ンとの結合が重要になってくる。ヒト血清アルブミンと環学シミュレーションを 8 回実施した。分子動力学シミュレ状ペプチドの複合体構造の報告数は限られており、ヒトーションによって得られたトラジェクトリに対して相互作

血清アルブミンへの結合メカニズムを解明するにはより用解析を実施し、ヒト血清アルブミンとコリスチンの間に多くの複合体構造に対して相互作用解析を行う必要が形成されている重要な相互作用の同定を行った。また、ある。結晶構造解析によって新たに得られたヒト血清ア複合体構造およびコリスチンを削除したヒト血清アルブ

ルブミンと環状ペプチドの複合体構造に対して分子動ミンの構造変化の様子を調査した。

力学シミュレーションおよび相互作用解析を実施し、ヒ
結果および考察ト血清アルブミンとの結合に重要となる相互作用を明ら

相互作用解析の結果、複数の分子動力学シミュレ

かにする。

35

ーションにおいて全シミュレーション時間中の 30%以上にわたって相互作用を形成していたアミノ酸はChain A では LYS212 と ASP324 であった。Chain
B では LYS12 と ASP13 であった。Chain A のLYS212 はコリスチンの直鎖構造における主鎖酸素と水を介した水素結合を形成していた。Chain A のASP324 はコリスチンの直鎖構造における主鎖 NHとの水素結合を形成していた。Chain B の LYS12はコリスチンの環構造における主鎖 NH との水素結(b)コリスチンを削除したヒト血清アルブミンの rmsd 変合を形成していた。Chain B の ASP13 はコリスチン

化の環構造における主鎖酸素と水を介した水素結合、

図 1. 分子動力学シミュレーションによる構造変化

環構造における側鎖 NH2 と水素結合を形成していた。興味深いことにコリスチンの環構造および直鎖
まとめ、今後の課題

構造どちらも主鎖を用いた相互作用が多く見られヒト血清アルブミンとコリスチンの複合体構造はヒた。
ヒト血清アルブミンとコリスチン複合体構造およびト血清アルブミンとダルババンシンの複合体構造とコリスチンを削除したヒト血清アルブミンの構造変化

大きく異なっていることから、２つの複合体の形成メの様子を rmsd (root mean square deviation)を用

カニズムについて詳細な解析が必要である。

いて可視化した。ヒト血清アルブミンとコリスチン複合体構造では 1 µs の分子動力学シミュレーションで

参考文献

5.7～9.7 Å 程度の構造変化が観測された（図 1(a)）。

Ito, Sho, et al. "Structural basis for the binding

一方、コリスチンを削除したヒト血清アルブミンでは

mechanism

5.9～22.9 Å 程度の構造変化が観測された（図complexed with cyclic peptide dalbavancin."1(b)）。以上のことから、コリスチンの結合によるヒト

Journal of medicinal chemistry 63.22 (2020):

血清アルブミン構造の安定化が示唆された。

14045-14053.

(a) ヒト血清アルブミンとコリスチン複合体構造 rmsd
変化

36

of

human

serum

albuminTSUBAME 共同利用 令和４年度 学術利用 成果報告書

利用課題名 廃棄物最終処分場における間隙内流体挙動の数値解析
英文：Numerical analyses of fluid dynamics in pores of landfill waste layer
利用課題責任者 鈴木和将
Kazuyuki Suzuki所属 埼玉県環境科学国際センター
Affiliation Center for Environmental Science in Saitama
URL http://www.pref.saitama.lg.jp/cess/index.html

邦文抄録
本研究では、廃棄物最終処分場内部の水やガスの流れ問題の高品質な計算スキームの開発を目的として、数値流体解析の検討を行った。具体的には支配方程式に Navier‐Stokes 方程式を用いて SUPG/PSPG 法に基づく安定
化有限要素法により離散化した。さらに,GPGPU コンピューティングによる計算の高速化を試みた。連立一次方程式の求解には GPBi‐CG 法を導入し、マルチ GPU を利用して高速並列計算によるシミュレーションを行った。
英文抄録
The purpose of this study is to establish a numerical simulation model of fluid flow in a landfill layer
with high quality and high precision. A 3D finite element method is a powerful tool for flows having
complex geometry such as porous media in landfill. However, it leads a huge amount of computation
cost. In this study, we examined to accelerate the 3D FEM by using the Graphics Processing Unit as a
general-purpose use (GPGPU).

Keywords: landfill, Numerical simulation, FEM, Navier-Stokes equations, GPGPU

ル化した流れ場において流体シミュレーションを行い、背景と目的

近年、資源開発や地球科学の分野では、岩石の流体特性と間隙形状の関係を明らかにすることを目ような多孔質媒体における流体特性を明らかにする的とする。

ため、デジタル岩石物理学の研究が進んでいる 1)。
デジタル岩石物理学は、岩石をデジタル化し、それ

概要

に対し数値シミュレーションを適用することにより、不試料には、一般廃棄物焼却施設から採取した主均質な岩石の間隙形状を考慮して水理特性といった灰等の廃棄物試料を用いた。

物理特性を求めることが可能となる。この手法は、複これらの試料は、マイクロフォーカス X 線 CT 装置

雑な廃棄物層間隙形状を有する廃棄物最終処分場を用いて撮影した。その後、試料の CT 画像を画像の流体挙動の解析においても有用なツールとなるも処理ソフトウェアに読み込み、モデル間隙形状の作のと期待される。廃棄物層内の流体挙動は、間隙形成を行った。さらに、作成した間隙形状のデータを状に強く依存すると考えられるものの、従来の方法

STL ファイルへ変換し、この STL ファイルを基に、自では、定量的な間隙形状分布や間隙に関する幾何

動メッシュ生成ソフトウェア Gmsh2)を用いて有限要学的な位置情報を表すことができなかったため、間素メッシュの生成を行った。

隙形状と流体挙動との相互作用は十分に解明され数値シミュレーションの支配方程式には、下記ていない。そこで、本研究では、マイクロフォーカス X(1),(2)式に示す非圧縮性 Navier-Stokes 方程式、

線 CT を用いて、廃棄物材料の間隙に関する 3 次元

連続の式を用いた。

ボリュームデータを取得し、このデータから定量的な

u i
u i
p
1   u i u j 
in , (1)
+u j
=−
+
+
t
x j
xi Re x j  x j xi 間隙の幾何学的な情報の抽出を試みるともに、モデ

37

u i
= 0 inx i

. (2)

ここで、 u は流速、 p は圧力、Re は Reynolds 数、
 は計算領域である。
上記(1),(2)の支配方程式に対して、安定化有限要素法（SUPG/PSPG 法）3,4)を適用すると以下のよ
うに弱形式が導かれる。
∫ 𝑤𝑖 (
𝛺

𝜕𝑢𝑖
𝜕𝑢𝑖
𝜕𝑤𝑖
+ 𝑢̅𝑗
) 𝑑𝛺 − ∫
𝜕𝑡
𝜕𝑥𝑗
𝛺 𝜕𝑥𝑖
1 𝜕𝑤𝑖 𝜕𝑢𝑖 𝜕𝑢𝑗
(
+
)
𝜕𝑥𝑖
𝛺 𝑅𝑒 𝜕𝑥𝑗 𝜕𝑥𝑗

+∫
𝑛

𝑒𝑙
+ ∑𝑒=1̅̅̅𝑘∫𝛺 (𝜏𝑢

𝜕𝑢

𝜕𝑤𝑖
𝜕𝑥𝑘

)(

𝜕𝑢𝑖
𝜕𝑡

+ 𝑢̅𝑗
𝜕𝑞

𝑛

𝜕𝑢𝑖
𝜕𝑥𝑗
𝜕𝑢

+

𝜕𝑝
𝜕𝑥𝑖

図 1 焼却灰試料における数値シミュレーション結果) 𝑑𝛺 = 0 ,(3)

𝜕𝑢

𝜕𝑝

𝑗

𝑖

そこで、これらの試料に対して定量的な間隙の構𝑒𝑙∫𝛺 𝑞 𝜕𝑥𝑖 𝑑𝛺 + ∑𝑒=1
∫𝛺 (𝜏 𝜕𝑥 ) ( 𝜕𝑡𝑖 + 𝑢̅𝑗 𝜕𝑥 𝑖 + 𝜕𝑥 ) 𝑑𝛺 =
𝑖

𝑖

0,

造解析を行った。
図 2 に、焼却灰試料の間隙のネットワーク構造を(4)

ここで、 wi , q は、それぞれ式(1),(2)に対する細線化した結果を示す（中心線分析）。焼却灰は、多Galerkin 項の重み関数である。また、 u i は移流速度

くのノードが存在しており、それらが複雑に結びつくを表し、  は安定化パラメータを表している。これらことにより多数のパスを形成している。これは、流体(3),(4)式に対して Crank-Nicolson 法により時間方向移動が可能となるルートがより複雑になることを示唆の離散化を行うことで、有限要素方程式が導かれしており、数値シミュレーション結果と同様の傾向でる。導かれた連立一次方程式の解法には GPBi-CGあった。

法 5)を適用した。
本研究では、計算の高速化を達成するために、
GPGPU を用いた並列計算を取り入れ、NVIDIA 社の統合開発環境である CUDA とそれらのライブラリ(cuBLAS, cuSPARSE)を用いて GPBi-CG 法を適用した。
また、X 線 CT 画像データからの三次元モデルの作製及び間隙構造解析は、Simpleware 社製の
ScanIP を用いて行った。
結果および考察図 2 焼却灰のセンターラインネットワーク

上部の面より垂直に流入する一様流を与え、流れの数値シミュレーションを実施した。一例として焼却
灰の間隙内速度ベクトルを図 1 に示す。焼却灰は、他の廃棄物試料と比較して、粒径が小さく、密に充まとめ、今後の課題

填されているため、より複雑なパスと通り水が流れて今後、間隙内流体挙動と間隙形状との関係性についる様子が確認できた。

いて評価を行い、廃棄物最終処分場の水分移動現
象のメカニズム解明につなげていきたいと考えてい
る。

38

参考文献
1)
2)

3)

4)

5)

辻健:貯留層マネージメントに向けたデジタル岩石の利用、石油技術協会誌,84(6) 403-410 (2019)
Geuzaine, C. and Remacle, J.F.: Gmsh: a three-dimensional
finite element mesh generator with built-in pre- and postprocessing facilities, International Journal for Numerical
Methods in Engineering, 79(11), 1309-1331(2009)
Tezduyar, T.E.: Stabilized finite element formulations forincompressible flow computations, Advanced in Applied
Mechanics, 28, 1-44(1991)
Tezduyar, T.E., Mittal, S., Ray, S.E. ＆ Shih, R.:
Incompressible flow computations with stabilized bilinear
and linear equal-order-interpolation velocity-pressure
elements, Computer Methods in Applied Mechanics andEngineering, 95, 221-242(1992)
Zhang, S.L.: GPBi-CG: Generalized product-type methods
based on Bi-CG for solving nonsymmetric linear system,
SIAM J. Sci. Comput., 18, 537-551(1997)

39

40

TSUBAME 共同利用

令和 4 年度 学術利用

成果報告書

利用課題名 LRnLA アルゴリズムを用いた物理シミュレーション
英文：Simulation of Physical Processes with LRｎLA Algorithm
善甫 康成
Yasunari Zempo
法政大学 情報科学部
Computer and Information Sciences,Hosei University
http://cis.k.hosei.ac.jp/
邦文抄録（300 字程度）
計算の中で CPU と GPU を同時に使い並列処理を行うヘテロな計算を行えば、利用可能なすべての処
理能力を活用することができる．もちろん CPU のピーク性能とメモリ帯域幅が GPU と比べ 20 倍程度の低いので、シミュレーションに CPU と GPU の両方を同時に使用することは非効率的な場合が多い．
そこで CPU の並列化効率を向上させ、GPU と同等の性能を発揮できる新しい時間ブロックアルゴリズ
ムを考案した．空間と時間にまたがるサブタスク間でデータを交換するようにデータ構造を工夫し、利用することとした．我々がこれまで開発してきた LRnLA 型の 1 つの有用なアルゴリズムとなる．
英文抄録（100 words 程度）
A heterogeneous computation can engage both CPU and GPU parallelism at the same time, andtake advantage of all available processing power. But CPU peak performance and memory bandwidth
are at least 20 times lower than that of GPU, so simultaneous use of both CPU and GPU for
simulation may be inefficient. We propose a new temporal blocking algorithm of the LRnLA family
that can increase the efficiency of CPU parallelization so that the performance is on par with GPU.
The algorithm uses the data structure designed specifically to exchange data between subtasks
spanning space and time.

Keywords: LRnLA ・temporal blocking・LBM・hwloc・many-core
概要
Finite Difference と Lattice Boltzmann などの数値計算スキームを CPU で効率的に実行するために、実装したアルゴリズムは ConeTorre (CT)というLRnLA(Locally
RecursivenonLocally
Asynchronous) である[1, 3]．並列ではない計算でも、メモリバンド幅が性能ボトルネックになる場合、CT は L1, L2, L3 キャッシュでデータを局在化させることができ、データアクセス時間を減少させることができるという特徴を持っている．CT は𝑁𝑁と𝑁𝑁𝑇𝑇という引数がある関数かサブルーチンとして理解す
るとわかりやすい．図 1 に示したように、時空間次元𝐷𝐷の計算で、𝐶𝐶(𝑁𝑁)を計算領域とする．つまり各辺に 𝑁𝑁 個 の セ ル を 持 つ 立 方 体 と す る ． ConeTorre
𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑁𝑁𝑇𝑇 ) という関数で、𝐶𝐶(𝑁𝑁)の中のスキーム計算
実行していく．時間軸が1から𝑁𝑁𝑇𝑇 まで分割されている．各軸で１セルずつずれて𝑁𝑁𝑇𝑇 回繰り返し、𝐶𝐶(𝑁𝑁)の計算を実行する．
さらに𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑁𝑁𝑇𝑇 )は、𝐶𝐶𝐶𝐶(𝑛𝑛, 𝑛𝑛 𝑇𝑇 ) (𝑛𝑛 < 𝑁𝑁, 𝑛𝑛 𝑇𝑇 < 𝑁𝑁𝑇𝑇 )というサブタスクに分解できる(図 1)．𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑁𝑁𝑇𝑇 )で計算すべき時空間を覆い、全計算を実行する．𝑁𝑁と𝑁𝑁𝑇𝑇 というパレメターがあるため、𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑁𝑁𝑇𝑇 )というサブタスクのデータ量を調整できる．分解方法を上

背景と目的
GPU を搭載したマシンにおいて、GPU ばかりでなく CPU の計算処理能力も利用したヘテロな計算が行えれば、すべての計算処理能力を使うことができる．しかし奈良が、GPU と CPU のピーク性能やメモリーバンド幅には大きな差がある．大規模な数値計算を GPU を用いて計算を行うときでも、CPUも当然利用している．この CPU をどれほど効率的に行うのが、今回の研究テーマであった．
前回までの報告で、CPU 制御下の RAM に保存されているデータを GPU で処理しても計算速度が落ちない計算手法の開発の報告を行った[1,2]. その手法を使うと、GPU グロバールメモリに入れるほどの計算ばかりではなく、CPU 制御下の RAM やSSD に入れなければならないほどの数 TB のサイズの計算までできるようになった．さらに計算性能を向上させるてめには、CPU の処理能力まで利用するが望ましい．さて TSUBAME3.0 で設置されている CPU の性能は、GPU に比べるとかなり低く、その性能は約 20 分の 1 程度である．この環境で CPU性能を利用し、CPU でも GPU 並みのパフォーマンスで同じ計算できるような新しいアルゴリズムを開発することが目的である．
41

タを L1 キャシュに入れるという方法である．

手 く 選 択 す れ ば 、 L2 キ ャ ッ シ ュ に 入 れ る
𝐶𝐶𝐶𝐶(𝑛𝑛𝐿𝐿2 , 𝑛𝑛 𝑇𝑇 𝐿𝐿2 )が得られる．これをさらに分解すれば L1 に入れるタスクが得られる．これによりキャッシュでのデータの局在性を十分活用することができ、計算のパーフォマンスを向上させることができ
る．

図 3:

FArShFold の並列化

図 4．FArShFold の並列化とデータ依存の様子

図 1. 𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑁𝑁)の数値計算領域と更に分解した𝐶𝐶𝐶𝐶の様
子．矢印はデータ依存を示す．

現在の CPU ではマルチコアが一般的である．従って、マルチコアによる並列処理を行う場合、コアあたりの L2, L3 キャッシュ量が減少するため、データの局在化によるパーフォマンス向上が期待できない．そこで、この課題を解決するため、新しい
FArShFold という並列法を開発した．
FArShFold というアルゴリズムを導入するため、
FArSh (Functionally Arranged Shadow) というデータ配置が必須である[4]． FArSh データは、隣にある CT の間の交換のために考案されたデータ配列である．CT 計算を行うための元の計算領域の𝐶𝐶(𝑁𝑁)のデータは、主メモリストレージ（RAM か SSD)からロードし、処理後同じのメモリストレージに保存される．CT 計算が行われる時空間での斜めのスロ
ープのデータは、FArSh データ配列からロードしCT 計算の後、同じ FArSh データ配列に上書きされる (図 2)．

その結果、𝐶𝐶𝐶𝐶の間のデータ交換では、𝐶𝐶𝐶𝐶スロー
プのデータは各コアの L2 および L3 キャッシュに残る．コア間の交換には、𝐶𝐶(2) のデータしか送らないことになるので、並列化の効率は飛躍的に高まる.
TSUBAME3.0 ノードの場合は Intel Xeon E52680 V4 が 2 台ある． 各 CPU の L2 キャッシュは256KB であって、L3 は L2 より 10 倍大きい．コア
数は 14 で(5 + 5 + 4)の 3 コアクラスター構成である．これらのコアは L3 キャッシュを共有している．
そこで𝑛𝑛’ = 3𝑛𝑛𝑠𝑠 として𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑛𝑛𝑠𝑠 )を5(4)スレッドで処理すれば、その5(4)スレッドが L3 キャッシュに局
在化することになるので、FArSh を共有できるはずである(図 4)．
結果および考察考 案 し た FArShFold ア ル ゴ リ ズ ム を 実 際 の
FDTD コードに実装した．Hwloc ライブラリを使用
してコアあたりの処理を分ける.AVX ベクトル演算を使用も実装した．数値スキームとして、LBM を選択した．LBM で、セルあたりのデータとセル更新あたりの浮動小数点演算（FLOP)を詳しく調べることができる．具体的には LBM 法のステンシル（Q)と衝突演算子を変更する．それによりセルデータ量によってパーフォマンスの変化を調べることができる．

図 2. FArSh データ配列の使い方

マルチコアの並列化で L2, L3 上のデータ局在性を保つため、次のような並列化法が考えられる．まず、𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑁𝑁)で計算領域を確保する(図 1)．次に、𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑁𝑁)を𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑛𝑛’)に分解する(𝑛𝑛’ < 𝑁𝑁)．𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑛𝑛’)をさらに𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑛𝑛𝑠𝑠 ) に分解する (𝑛𝑛𝑠𝑠 < 𝑛𝑛’) (図 3)．こ
こで𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑛𝑛𝑠𝑠 )それぞれ処理はスレッド並列により実施する．その結果、空間並列化ではなく、時間軸で並列化を行うことが可能になる．その時 FArSh
データ配列はスレッド数に合わせて分解し、1 コアの L2 キ ャ ッ シ ュ 分 を 各 コ ア に 割 り 付 け る ．
𝐶𝐶𝐶𝐶(𝑁𝑁, 𝑛𝑛𝑠𝑠 )を更に𝐶𝐶𝐶𝐶(2, 𝑛𝑛𝑠𝑠 )に分解する．これに対応する元の計算領域（𝐶𝐶(2)という立方体）にあるデー

図 5．パフォーマンスベンチマーク

42基本的なベンチマークは Ryzen R9 5950X で行った ( 図 5) ． 測 定 単 位 は Pups (LBM population
update per second)を用いた．これは 1 秒あたりのセル更新数 (LUps - Lattice cell Update per second)に𝑄𝑄を掛けたものである．
時空間ブロッキングを使用していない LRnLA の場合、パフォーマンスピークは~6 Pups であるが、CT を用いた場合、その 5 倍約 30Pups となりかな
りの性能向上が期待できる(図 5)．またこのベンチマークで𝑛𝑛’の最適値を得ることもできる．最適パー
フォマンスは１GLUps を超えることも分かった．
まとめ、今後の課題
新しいアルゴリズム ConeTorre を考案した．このアルゴリズムを FDTD コードに実装して、１個のマルチコアの CPU において GPU と同等のパー
フォマンスを得た．これはマルチコアの並列化において各コアに対応する L2, L3 にデータを局在化させることができたためである．
今回は基礎検討の段階であるが、GPU でも計算処理を実行する場合、同じ CT の時空間データ領域を分解できる．今回考案したアルゴリズムを用いれば GPU と CPU 組み合わせたヘテロ計算を実行できるはずであり、GPU だけの計算のパーフォマン
スより更に向上させるものと期待できる．
参考文献[1] Zakirov, A, V Levchenko, A Perepelkina, andYasunari Zempo. "High performance FDTD
algorithm for GPGPU supercomputers." In
Journal of Physics: Conference Series, vol. 759,
no. 1, p. 012100. IOP Publishing, 2016.
[2] Shimokawabe, T., Endo, T., Onodera, N.,
Aoki, T.: A stencil framework to realize largescale computations beyond device memory
capacity on GPU supercomputers. In: 2017
IEEE International Conference on Cluster
Computing (CLUSTER), pp. 525–529 (2017).
IEEE[3] Levchenko, V., Perepelkina, A. Locally
Recursive
Non-LocallyAsynchronous
Algorithms
forStencilComputation.
Lobachevskii J Math 39, 552–561 (2018).
https://doi.org/10.1134/S1995080218040108[4] Perepelkina, A., Levchenko, V.D. (2021).
Functionally Arranged Data for Algorithms
with Space-Time Wavefront. In: Parallel
Computational Technologies. PCT 2021. CCIS,
vol 1437.

4344TSUBAME 共同利用 令和４年度 学術利用 成果報告書

発達期多元自発活動の数理モデルとその学習理論の構築
The mathematical model and learning theory for pluralistic brain activity in early development
早川 隆
Hayakawa Takashi
日本大学 医学部 生体機能医学系 薬理学分野
Department of Pharmacology, School of Medicine, Nihon University
https://www.med.nihon-u.ac.jp/department/pharmacology/
本研究では、課題責任者の科研費プロジェクト「発達期多元自発活動の数理モデルとその学習理論の構築」の遂行のために、大脳皮質局所神経回路モデルの数値シミュレーション及びその神経回路を用いた機械学習を GPU アク
セラレーターを用いて大規模かつ高速に実行するシミュレーターを作成した。これを用いて発達期に形成されるとされるモジュール構造を持つ局所神経回路モデルが、学習後に実験的に観察されるモジュールごとの同期的なダイナ
ミクスを再現すること、そしてこの振る舞いが新たに構築した統計力学理論によって記述できることを示した。
In the present study, I developed a generic program code for simulation of mammalian neocortical-circuit
models and machine learning using these models, for the goal set out in a KAKENHI grant project
entitled “the mathematical model and learning theory for pluralistic brain activity in early
development.” This fast simulator that takes full advantage of GPU accelerators allowed me to show thata cortical-circuit model with module structures exhibits module-wise synchronous dynamics after
learning as do neocortical circuits of animals, and that these dynamics can be described by a newly
developed statistical mechanical theory.

Keywords:

mammalian neocortex, computational neuroscience, machine learning, statisticalmechanics, large-scale simulation

背景と目的な大脳皮質の神経細胞集団が示す、同期的かつ不規

本課題は研究課題責任者の早川が同名の学術変革

則/規則的な振動現象の起源を記述できる可能性のあ領域(B)科研費課題を遂行するためのものである。以る新しい統計力学理論を世界に先駆けて構築していた下に科研費課題の背景と目的を述べるが、本共同利[3]。また、学術変革領域(B)の共同研究者の千葉は、

用課題は、以下の科研費課題の目的を達成するため

集団同期現象を取り扱う際に広く用いることができる数に、TSUBAME の持つ複数 GPU を利用した大規模計学理論を世界に先駆けて構築していた[4]。

算環境を利用し、脳神経回路モデルの高速シミュレー

このような領域の研究者達の先行研究と並行して、ションによって研究を加速しようとするものである。

近年、神経回路の機械学習の研究者達が、統計力学

本研究の背景としてはまず、近年になって発達初期の理論を用いて深層学習の振る舞いを解析しており、一脳回路の活動様式が実験的に観察できるようになり、定の初期条件から学習を開始すると、その後の学習性

この時期の動物脳が従来想像されていたものよりダイ能が高くなることを理論的に示していた[5]。

ナミックで普遍的な多元自発活動を示すことが明らか

以上の背景を踏まえて、本課題で研究代表者は、動になってきたことがある。学術変革領域(B)の共同研究

物脳の初期状態である発達初期の脳活動は、リカレン

者の上阪・水野らはこの時期の脳活動が、発達と共に、ト神経回路である大脳皮質回路が何らかの意味で最同期的な活動から非同期的な活動に切り替わっていく

適な状態に調整された結果として現れているのか、と様子を世界に先駆けて観察することに成功していたいう根本的な疑問に対して、共同研究を通して答えよう[1,2]。

とする。より具体的には、先行研究を踏まえて大脳皮質

一方、課題責任者の早川は、脳波などに見られるようの神経回路モデルを構築し、シナプス結合強度分布な45

どのパラメータを発達初期の各段階の動物脳における脳皮質の神経回路モデルを立てれば良いのかを追究

実験結果と合致するものに設定した際に、領域の共同した。共同研究者との議論を通じて、先行研究で課題研究者らがこれまでに観察してきた現象を再現できる

責任者が調べてきた無構造の局所神経回路モデルにかどうかを確認する。そして再現された同期活動から非則り解析・議論を進めるだけでは、観察される現象を十同期活動への転移が、神経回路モデルにその後機械分に説明しきれないという結論に至った。そこで、先行学習を行わせた際の学習性能にどのように影響するの研究において構築した神経回路モデルと統計力学理

かを、課題責任者とその他のグループが先行研究で構論を、実際の大脳皮質神経回路と同様にモジュール構

築した統計力学理論及び数値シミュレーションの技術造を持つ場合に拡張した。というのも、成体の脳では神を用いて調べる。これによって、脳の非常に柔軟な情経回路はこのようなモジュールごとに同期活動を示して報処理能力の根源を探ると同時に、フィードフォワード

おり、それが確率的な情報を表現していることが最近神経回路に比べてあまり理論的に明らかにされてこなの研究で明らかになってきているからである[6,7]。従っかったリカレント神経回路の機械学習理論を打ち立てて、本研究で問うべき課題は、無構造の大脳皮質からることを目的とする。

学習に伴ってこのモジュール構造がどのようにして現概要れるかである。

発達を通して大脳皮質局所回路に形成されるモジュー

そこで課題責任者は、課題責任者自身が先行研究で

ル構造を組み込んだ神経回路モデルとその学習の数用いた個々のニューロンの神経発火率(0~2 の連続値

値シミュレーションを GPU アクセラレーターを利用してをとる)が常微分方程式で表される単純モデル、本研究

大規模かつ高速に実行できるシミュレーターを構築した。

課題の主な対象である大脳皮質感覚野の神経細胞の神経発火頻度が常微分方程式に従い連続的に時間変発火機構をより精緻にモデル化したリアリスティックモ

化するモデル及び、現実の脳神経細胞の発火活動機

デルの両方をベースに、神経結合が無構造及びモジュ構を取り入れた精緻モデルの 2 通りのシミュレーターをール構造有りの両方の場合のシミュレーションを実行

作成し、これらのシミュレーション結果を裏付ける統計

でき、さらに単純モデルの場合は機械学習課題も実行

力学理論も構築した。これらを用いて大脳皮質神経回できるような数値シミュレーターを C++・OpenACC を用

路の多元的自発活動と学習の一端が解明された。

いて構築することに成功した。例えば、下図に示すよう

結果および考察にモジュール構造を持つ単純モデルにおいて、モジュ令和４年度は背景と目的で述べた本研究における根ールの秩序変数が特定のダイナミクスを示すよう学習本的な疑問に答えようとするにあたって、どのような大する数値シミュレーションを実行できる。発達過程を想46

定し様々なモデルの予備的なシミュレーションを行なっ

参考文献

たが、いずれの場合も GPU アクセラレーターによって[1] Ikezoe et al. (2018) Cell Rep.

大幅な計算速度の向上が得られた。

[2] Uesaka et al. (unpublished data)

単純モデルをベースとした場合、モジュール構造を持[3] Hayakawa & Fukai (2020) Phys. Rev. Res.

つ場合に神経細胞集団がモジュールごとに同期活動す[4] Chiba (2015) Adv. in Mathる様子を統計力学理論を用いて記述することにも成功[5] Bahri et al. (2020) Ann.Rev.Cond.Matt.Phys.

し、このシミュレーターの結果を理論的に裏付けること[6] Ko et al. (2011) Natureにも成功した。また、モジュール構造を持った場合、モ[7] Shi et al. (2022) Nat. Comm.

ジュールの同期による確率表現を変更する学習は、新[8] Bernardi et al. (2021) PLoS Comput. Biol.

たな表現を学習し直す場合に比べてはるかに高速に行うことができることも、数値的・理論的に示した。また単純モデルの範囲で、成体の情報表現として知られる持続的・漸増的な集団発火活動を機械学習の結果として再現することにも成功した。これらの結果は過去の実
験観察事実と定性的にはよく整合し、単純モデルであることによる限界はあるが、成体脳の学習の本質を捉えている可能性が示唆された。
次に、正則発火ニューロン・高頻度発火ニューロン・
低域値発火ニューロンの 3 種類から構成される、現実の大脳皮質神経細胞の発火活動機構をある程度忠実にモデル化した精密モデルのシミュレーターを先行研
究[8]を参考に構築した。このシミュレーターでも、単純モデルの場合と同様に、モジュールごとの集団同期活
動を GPU アクセラレータを用いた大規模かつ高速な数値シミュレーションによって再現し、これを統計力学理論によって裏付けることに成功した。

まとめ、今後の課題
令和 4 年度は研究プロジェクトの 1 年目ということもあり、基盤となるシミュレーター及びそれに対する理論の適用を確認する作業が主となった。令和 5 年度はいよい
よ、このシミュレーターのパラメータや構造を発達初期のものに連続的に変更していった際に、どのように活動が変化していくのかを調べていく。また、無構造な発達初期の神経回路モデルに機械学習を行った際にモジュール構造が現れ、確率的情報表現が可能となるのか、その際の最初の神経回路モデルのパラメータはどのよ
うな値に取るべきか、などを数値シミュレーション・理論
解析を駆使して調べていく。

47

48

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

深層学習を用いた分子動力学シミュレーションの高速化
Acceleration of molecular dynamics simulation using deep learning
林 周斗
Shuto Hayashi
東京医科歯科大学 難治疾患研究所
Medical Research Institute, Tokyo Medical and Dental University
https://www.shimamlab.info
タンパク質が持つ機能を解析するため、分子動力学シミュレーションによりタンパク質立体構造のダイナミクスを解析する研究が広く行われている。しかしながら、分子動力学シミュレーションには大きな計算コストを必要とするという問題点があり、長時間・大規模なシミュレーションを行うにはスーパーコンピュータを用いても膨大な時間を要する。
そこで本課題では、分子動力学シミュレーションを高速かつ高精度に模倣する深層ニューラルネットワークを開発した。ベンチマークデータセットを用いた検証により、我々の手法が任意のタンパク質に対して高速かつ高精度にシミュレーションを実行できることを示した。
To analyze the functions of proteins, molecular dynamics simulations are widely used to study the
dynamics of protein structures. However, molecular dynamics simulations have the drawback ofrequiring significant computational resources, and even with the use of supercomputers, conducting
long and large-scale simulations can be time-consuming. In this project, we developed a deep neuralnetwork that can mimic molecular dynamics simulations with high speed and accuracy. Throughvalidation using benchmark datasets, we demonstrated that our method can perform fast and accurate
simulations for any given protein.

Keywords: molecular dynamics simulation, deep learning, E(3)-equivariant neural network, graphneural network, protein structure

背景と目的

きるのは長くともマイクロ秒〜ミリ秒スケールのシミュレ

タンパク質が持つ機能はその立体構造によって決定

ーションまでである。

づけられており、タンパク質の立体構造を精査すること

そこで本課題では深層学習をベースとしたデータ駆はその機能を知る上で重要な手がかりとなる。これまで

動型のアプローチを用いることで分子動力学シミュレー

X 線結晶構造解析、核磁気共鳴分光法、電子顕微鏡

ションを高速に行う手法を開発することを目的とするなどの実験手法により、さまざまなタンパク質の立体構（図１）。

造情報が解析されてきた。しかしながら、これらの手法
分⼦動⼒学シミュレーション

によって得られる立体構造は一時刻点における静止像であるため、タンパク質立体構造の動的変化を捉える⊿t⊿t

ことはできない。実際、薬学分野においては、薬・標的

・・・

⊿t

膨⼤な計算コストタンパク質間の相互作用を見極めるために立体構造の動的変化を解析することが重要視されている。
深層学習

これまでタンパク質立体構造の動的変化を解析する

深層学習により１ステップで予測できるか？

研究が数多くなされてきたが、それらの中でも特に広く用いられているのが分子動力学シミュレーションである。
しかしながら、分子動力学シミュレーションには大きな

図１：本研究の背景

計算コストを必要とするという問題点があり、スーパーコンピュータを用いたとしても現実的な時間内で計算で49の比をプロットした（図４）。平均して約 20 倍の高速化を

概要

達成した。

本課題ではタンパク質の全原子座標、速度、力、原子種を入力として、一定時間後の全原子座標、速度、H

力を高速・高精度に予測する深層ニューラルネットワー

C

ク（LAMDA）を学習する。学習データセットとして、実際のタンパク質分子動力学シミュレーションデータを用いる。ネットワークの構築、および学習には PyTorch を用いる。
図３：LAMDA の予測誤差結果および考察タンパク質の全原子座標、速度、力、原子種を入力として、64 フェムト秒後の全原子座標、速度、力を予測するモデル（LAMDA）を学習した。LAMDA を学習させるため、AlphaFold2[1,2] によって推定されたヒトタン
パク質立体構造を初期構造として MD シミュレーションを実行し、そのトラジェクトリデータを利用した。
性能評価のため、ベンチマークデータセットに対して
分子動力学シミュレーションを適用して得られたトラジェ
クトリと LAMDA の予測を比較した。図２はあるタンパ
ク質のある水素原子とある炭素原子におけるトラジェク
トリ（黒丸）と LAMDA の予測（赤丸）である。これらの原子は非線形な動きをしているにも関わらず、
LAMDA が 64 フェムト秒後の座標を正確に予測できて

図４：MD シミュレーションと LAMDAいることが分かる。

の実行時間比HC

まとめ、今後の課題
本課題では分子動力学シミュレーションを模倣する
ニューラルネットワーク（LAMDA）を学習し、64 フェム
ト秒後の原子座標を高速かつ高精度に計算することが
error = 1.01×10!#Å

できた。長時間のシミュレーションを行うには LAMDA

error = 3.25×10!"Å

を繰り返し適用する必要があるが、安定して実行するこ

図２：MD によるトラジェクトリ（黒）ととは現状困難である。今後は安定した長時間シミュレー

LAMDA による予測（赤）ションを実現するための技術開発に取り組んでいく予

また、ベンチマークデータセット中の全 2172 タンパク定である。

質の全水素原子と全炭素原子に対して、タンパク質の全原子数と予測誤差をプロットした（図３）。ほとんどの参考文献

原子に関しては 10-2Åスケールの誤差であり、誤差の1.

大きさはタンパク質の大きさに依存しないことが分か

Jumper, J et al. Highly accurate protein structure
prediction with AlphaFold. Nature (2021).

る。

2.

最後に、MD シミュレーションと LAMDA の実行時間

Varadi, M et al. AlphaFold Protein Structure
Database: massively expanding the structural

50

coverageofprotein-sequence

space

with

high-accuracy models. Nucleic Acids Research(2021).

51

52

TSUBAME 共同利用 令和４年度 学術利用 成果報告書

電気コイル隙間内に流れ込む冷却液挙動に関する大規模数値解析
英文：Large-scale computing of pouring coolant on the electric coil
利用課題責任者 金田 昌之
First name Surname Masayuki Kaneda
所属 大阪公立大学
Affiliation Osaka Metropolitan University
URL https://www.omu.ac.jp/eng/htlab/staff/kaneda.html
邦文抄録（300 字程度）
電動自動車などに搭載されているモータの電気コイルを冷却するために，冷却液を直接流下する手法がとられている．液体の濡れ広がりと除熱の関係を明らかにするため，これまで開発してきたフェーズフィールド LBM 大規模二相流数値解析手法に温度場解析の実装を試みた．その結果解析不安定や非物理的な温度場を呈することがわかったため，LBM の時間発展方程式を改良し，安定解析できる手法へと発展させた．本手法を用いて，簡易化ステー
タコイル構造を対象とし，流下温度及び流量が除熱に及ぼす影響を定性的に明らかにした．さらに電動モータのロー
タ・ステータ間流れに着目し，大規模乱流解析を行うことで構造が及ぼす影響について明らかにした．
英文抄録（100 words 程度）
A direct liquid cooling is used to cool the electric coils of motors installed in electric vehicles. To clarify
the relationship between the liquid wetting spreading on the coil and heat removal, we implemented atemperature field analysis into the phase-field LBM large-scale two-phase flow numerical code wehave developed. The results showed that the analysis was unstable and the temperature field wasunphysical, so the time evolution equations of the LBM were improved and the method was developed
into a stable analysis method. Using this method, the effects of flow temperature and flow rate on heatremoval were qualitatively clarified for a simplified stator coil structure.

Keywords: Cooling of electric devices, two-phase flow, phase-field lattice Boltzmann method,
thermal simulation, Taylor-Couette flow
背景と目的検討は，複雑構造を取り扱える手法に限られ，さらに二近年の温室効果ガス削減目標に向けて，内燃機関

相流の解析スキル及び伝熱工学の知識が必要となる．

の電動化が促進されており，そのパワートレインである

一般商用ソフトである程度解析することは可能である

電動モータの小型高効率化への要請が高まっている．

が，解析安定性を重視しており，複雑構造二相熱流動

これはモータのコイル部分からの発熱密度の上昇を意

解析の精度がどの程度担保されているかは不明であ味しており，その効果的な除熱が課題となっている．実る．一方で我々は，複雑構造を通過する二相流を解析際の自動車用のモータでは，非電気伝導性の冷却液できる数値解析コードの開発を行ってきた．

を発熱する電気コイル部分に直接流下することで冷却一方で自動車用モータの高速回転化が模索されてしている．この場合，冷却液がまんべんなくコイル内をおり，特にロータ・ステータ間で生じる流体摩擦の低減

浸潤して流下することが求められるのはいうまでもない．

を議論する段階に来ている．しかしながら電磁場の最しかしながら冷却液の濡れと交換熱量の関連は明らか

適化のために設けられている内外円筒の溝構造が流にされておらず，それを予測するためのツールもないのれ場に及ぼす影響までは検討されていない．

が現状である．冷却液の挙動と除熱の関連が明らかに

そこで本プロジェクトでは，前者についてはこれまで

なれば，今後の機器設計のヒントとなる．

に開発してきたコードを改良し，温度場解析も可能とす

しかしながら，実際のコイル形状は大変複雑な構造ることを目的とした．フェーズフィールド格子ボルツマンをしているため，実験的に明らかにするには実機を用法（PELBM）に温度場解析用の分布関数を追加し，解いて経験的なアプローチしかできない．数値解析による析安定性を担保できる手法を開発した．後者について

53

は直接数値解析を実施することで溝が回転トルクに及結果および考察ぼす影響を議論した．

二相流熱流動場解析の一例として，解析領域を上面から見た濡れおよび液体の温度分布を示す．温度は構概要

造体近傍の水平面の分布のため，流体が流下してくぼ

これまで簡易モデルを対象とし，その流動特性をフェんでいる箇所は表示していない．詳細は省略するが，ーズフィールド格子ボルツマン法（PFLBM）による二相流速に応じて妥当な温度場が形成されていることが確

流解析により解明してきた．今年度は継続課題として，認できた簡易化ステータコイル構造の熱問題解析手法を構築した．そのため温度場解析用の LBM の分布関数を追加したものの，既存の方法では解析が発散したり，非物理的な温度場を呈するなどした．そこで高密度比解析
で問題となる，速度場の時間発展方程式に着目して検討したところ，速度場の解析の不安定さが非物理的な温度場を誘発していると考えた．文献調査を行い安定解析が報告されている組み合わせを実装することで，最終的に速度場温度場ともに安定的に解析できる手法へとたどり着いた．解析対象を図１に示す．

図３二相熱流動解析結果↓g

ロータ・ステータ間流れの速度分布を図４に示す．溝
内部の速度分布も精度良くとらえることができた．トル
クについては溝無しの場合と比較して４０％以上上昇
することが確認できた．

図１二相熱流動解析の解析対象

1.0

またロータ・ステータ間の流れ場については，上記のTC Re=1050
TC Re=3500
GTC Re=1344
GTC Re=4480

0.8<U>z /riLBM を改良して効率的に解析できる手法へと発展させた．ロータ・ステータにはそれぞれ溝が設けられており，そこで誘起される流れが回転トルクへ影響を及ぼす

0.6
0.4
0.2

と考え，乱流単相直接数値解析を実施し，既往の溝無0.0

しの場合と比較した．ステータのみに溝を設置した構造

0.0

0.2

0.4

0.6

0.8

1.0

1.2

1.4

1.6

1.8

(r-ri)/dを図２に示す．

図４溝付きロータ・ステータ間流れ解析の速度分布
まとめ、今後の課題

二相流熱流動解析を安定に実施できるコードを作成
し，大規模解析を行うことで交換熱量を予測することができた．また溝付きのロータ・ステータ間流れの直図２溝付きロータ・ステータ間流れ構造

接乱流解析により，トルクに及ぼす影響まで明らかとなった．今後は各種パラメータでの解析及び超高
速回転を考慮した解析を展開していく．

54

TSUBAME 共同利用 令和 ４年度 学術利用 成果報告書

利用課題名 深層学習を用いた分子動力学法ベースの水和熱力学量計算法の高速化
英文：A Fast Computation Method for Solvation Thermodynamic Quantities Using Deep Learning
吉留 崇
Takashi Yoshidome
東北大学 大学院工学研究科 応用物理学専攻
Department of Applied Physics, Graduate School of Engineering, Tohoku University
Grid Inhomogeneous Solvation Theory (GIST)は、水和エネルギーと水和エントロピーの空間的な分布を計算することで、水和を詳細に表すことが出来るものの、分子動力学法を用いるため大きな計算コストがかかることが問題である。そこで、本研究では水和自由エネルギー分布を予測する深層学習モデルの開発を行った。我々の深層学習モデルは、水和自由エネルギー分布を約 1 分程度で計算でき、かつ GIST の計算結果を定量的に再現することに成功した。
Grid Inhomogeneous Solvation Theory (GIST) can describe hydrations around proteins using the
spatial distribution of the hydration energy and hydration entropy. However, huge computation isrequired for the GIST because the GIST requires the trajectory obtained using molecular dynamics
simulations. In the present study, a deep learning (DL) model to predict the hydration free energy
distribution was proposed. It was found that our DL model can quantitatively predict the hydration
free energy distribution in less than 1 minute.

Keywords: タンパク質、水和、深層学習、グリッド不均一溶液理論

背景と目的

く水和分布予測法[2]と同じ深層学習モデルを採用した。

タンパク質における水和は、折り畳みやリガンドとの深層学習の学習データとテストデータとして、以前、分結合などの様々な過程において重要な役割を果たす。

子動力学法と GIST を用いて計算した 76 のタンパク質

Grid Inhomogeneous Solvation Theory (GIST)[1]の水和自由エネルギー分布を使用した。ハイパーパラ

は、水和エネルギーと水和エントロピーの空間的な分メータを最適化するために、128 通りのハイパーパラメ

布を計算することで、水和を詳細に表すことに成功した。

ータを用意し、それぞれに対し 2-fold cross validation

しかし、分子動力学法を用いるため大きな計算コストがを実行した。実行には 68 個のタンパク質を用いた。誤かかる。このため、GIST のドッキングアルゴリズムなど

差関数を最小とするハイパーパラメータを選び、残り 10への応用例は未だに少ない。

個でテスト計算を行った。さらに、誤差関数が小さいも本プロジェクトでは、GIST 計算を高速化することをのから 40 個選び、水和エントロピーと水和エネルギー目的とし、タンパク質立体構造からその周りの水和自を 予 測 す る 深 層 学 習 モ デ ル の 構 築 と 2-fold cross

由エネルギーの分布を予測する深層学習モデルの開validation を実行した。

発を行った。2-fold cross validation を用いて深層学
習モデルの最適なハイパーパラメータを決定し、深層結果および考察学習モデルを構築した。リガンド結合に伴い、タンパク

最適なハイパーパラメータを決定し、深層学習モデ

質表面から移動する水分子の自由エネルギーを計算しルを構築した後、そのパフォーマンスを調べるためにた所、GIST で計算したものと良く一致した。

10 個のタンパク質に適用した。具体的には、リガンドが結合する際にタンパク質から移動する水分子の自由エ
ネルギーを計算した。構築した深層学習モデルと GIST

概要
本研究では、以前我々が提案した、深層学習に基づの両方で行い、結果を比較した。その結果、深層学習

55

モデルで得た水分子の自由エネルギー（∆𝐺!""# ）のGIST 計算結果（∆𝐺$%&' ）からの誤差率（#$∆𝐺!""# −∆𝐺$%&' &/∆𝐺$%&' #）は、約 4%であった。この結果は、
我々の深層学習が高精度で GIST の計算結果を再現できることを意味する。また、我々の深層学習モデルは、
約 1 分程度で水和自由エネルギーの分布を計算でき、
この計算時間は GIST 計算（2 日程度）の約 1/1000 であった。よって、我々の深層学習モデルは、高速にかつ
正確に水和自由エネルギー分布を計算できることが分かった。

まとめ、今後の課題
本プロジェクトでは、高速にかつ正確に水和自由エ
ネルギー分布を計算できる深層学習モデルを開発した。
今後は、開発した深層学習とドッキングアルゴリズムとの融合を行う。

[1] C.N. Nguyen, T.K. Young, and M.K. Gilson,
J. Chem. Phys., 137, 044101 (2012).
[2] K. Kawama, Y. Fukushima, M. Ikeguchi, M.
Ohta, and T. Yoshidome, J. Chem. Inf. Model.,62, 4460 (2022).

56

TSUBAME 共同利用 令和４年度 学術利用 成果報告書

利用課題名 粘性の温度依存性を考慮した乱流熱流動の直接数値解析
英文：Direct numerical simulations of turbulent heat transfer with temperature-dependent viscosity
利用課題責任者
須賀 一彦
所属
大阪公立大学 工学研究科 機械系専攻
https://www.omu.ac.jp/eng/htlab/
流体の温度変化に伴う物性値の変化が乱流熱流動場に与える影響を理解する為に，平行平板間乱流を対象として温度による物性値変化を考慮した直接数値解析を実施した．流体には水を想定し，チャネル壁の上壁は 293[K],下壁は 353[K]の等温壁とし，流体の密度・熱拡散率は一定で粘性係数のみが温度によって変化する乱流熱流動場を解析した．解析は格子ボルツマン法を用いて行い，摩擦レイノルズ数 300, 650, 1000 の３条件で解析を行った．その結果，粘性一定条件と比較して，高温壁側（低温壁）では乱れが減少（増加）し，乱流摩擦抵抗が減少（増大）する
ことが分かった．また，粘性が変化する熱流動場においても，対数速度・温度領域を確認することができた．
To understand the effects of the temperature-dependent fluid properties on turbulent flows, we
performed direct numerical simulations (DNSs) of turbulent channel flows with temperaturedependent fluid properties where isothermal conditions were specified for the top (293K) and bottom(353K) walls. We considered a water fluid with variable viscosity but constant fluid density and thermal
diffusivity. The DNSs at the friction Reynolds number of 300, 650, and 1000 were performed by thelattice Boltzmann method. The results show that turbulence reduction (enhancement) occurs near thehot (cold) walls, leading to a decreased (an increased) frictional resistance. It is also found that thelogarithmic region is still retained for the mean velocity and temperature of variable-viscosity turbulent
flows.

Keywords: Lattice Boltzmann method, Turbulent channel flows, Heat transfer, variable-viscosity flow,
logarithmic law
解析を行った．流体には水を想定し，粘性係数の温度

背景と目的
一般的な数値シミュレーションでは，温度はパッシブ依存性のみを考慮した． 本研究では，粘性係数𝜇(𝑇)スカラーとして取り扱われ，温度が流動場に与える影は流体温度T[K]の関数としたモデル𝜇(𝑇) = 𝐴 × 10𝑇−𝐶

響は無視されることが多い．しかし，燃焼など急激な温を 用 い ， そ れ ぞ れ の モ デ ル 係 数 は A = 2.414 ×

度変化を伴う熱流動場では，物性値が局所の流体温10−5 [Pa s], B = 247.8[K]，C = 140[K]を用いた （１）．そ

度に応じて変化するため，パッシブスカラーとして取りの他の物性値は参照温度𝑇𝑟𝑒𝑓 = (𝑇𝑡𝑜𝑝 + 𝑇𝑏𝑜𝑡 )/2にお扱うことができない．しかし，このようなアクティブスカラ

ける値を使用した．流体は主流方向に圧力差を課してーの乱流熱輸送に関する知見は少なく，対応する乱流駆動し，摩擦レイノルズ数𝑅𝑒𝜏 = 300, 650, 1000の 3 条

モデルや壁モデルも少ない．本研究では，温度によっ

件で解析を行った．計算領域は，Reτ = 300のケースでて粘性係数が大きく変化する流体として水を対象とし，は ， 主 流 方 向 に 6𝐻 , ス パ ン 方 向 に 3𝐻 と し ， Reτ =アクティブスカラーの直接数値解析を実施する．得られ650, 1000のケースでは主流方向に3𝐻, スパン方向にた結果より，粘性係数の変化が乱流場に与える影響を1.5𝐻とした．直接数値解析は格子ボルツマン法を用い理解するとともに，乱流 モデルの構築に必要な基礎て行い，速度場は３次元２７方向多緩和時間格子ボル的なデータベースを作成することを目的とする．

ツマン法（２）を用いて解析し，温度場には３次元１９方向

概要

速度正規化格子ボルツマン法 （３）を用いた．また，不均

本研究は，上 壁𝑇𝑡𝑜𝑝 = 293[K], 下壁𝑇𝑏𝑜𝑡 = 293 [K]

衡修正の局所細密格子法 （４）を用いて，物理量が急激の等温壁とした平行平板間乱流を対象とした直接数値に変化する壁近傍領域に細かい計算格子を配置した．

𝐵

57

本研究では，壁面近傍では内層スケールで無次元化された格子幅がΔ+ = 1.5となるように設定した．総格子
点数は𝑅𝑒𝜏 = 300, 650, 1000のケースにおいて，それぞれ 4.5 億点，7.3 億点，22.5 億点となった．プログラ
ムコードは CUDA Fortran によって記述されており，複数の GPU を用いた領域分割法によって並列化した．
結果および考察
図 １ に 壁 面 近 傍 𝑦 + = 𝑦𝑢𝜏 /𝜈𝑟𝑒𝑓 ≃ 10 に お け る ，𝑅𝑒𝜏 = 1000のケースの主流方向変動速度のスナップ
ショットを示す．ここで，𝑦はそれぞれの壁からの距離を

図２ 𝑅𝑒𝜏 = 1000における主流方向・壁面垂直方向

表し，𝑢𝜏 はそれぞれの壁における摩擦速度を表し，𝜈𝑟𝑒𝑓

の乱流強度分布．乱流強度𝑢𝑟𝑚𝑠 , 𝑣𝑟𝑚𝑠 はそれぞれ

は参照温度における動粘性係数を表している．図より，の壁面の摩擦応力の平均値を基にした摩擦速度で粘性係数一定条件の流れで見られる乱流構造と同様

無次元化を行っている．

に，高温・低温壁の近傍では変動速度の高速・低速領域が主流方向に伸びるストリーク構造が確認される．し次 に ， 乱 流 統 計 量 に つ い て 議 論 を 行 う ． 図２ に ，かし，高温壁・低温壁の近傍の乱流構造には違いが見𝑅𝑒𝜏 = 1000における主流方向・壁面垂直方向の乱流られる．高温壁の近傍に見られるストリーク構造は，低強度を摩擦速度で無次元化した結果を示す．図より，温壁近傍の構造と比較してスパン方向幅が狭く，変動

𝑦 + < 10の壁面の近傍では，高温壁側の乱流強度が強度も強くなっていることが確認できる．

低温壁側の結果を上回っており，主流方向・垂直方向の乱れ強度ともに鋭い立ち上がりを見せている．いっぽうで，壁面から離れた領域（𝑦 + > 10）では，乱れ強度の最大値は低温壁側のほうが大きく，チャネル中央
まで強い乱れが維持されている様子が分かる．なお，
ここでは𝑅𝑒𝜏 = 1000の結果のみを示しているが，この傾向は他の摩擦レイノルズ数条件でも見られる．
最後に，図３に対数速度分布，対数温度分布を示す．
なお，温度には，壁面温度からの偏差Θを示している．
また，比較の為に粘性係数一定条件の結果として速度分布には直接数値解析結果 （５），温度分布には経験式（６）を示している．図より，得られた平均速度分布の傾き

は，粘性係数一定条件の結果とおおむね一致しており，粘性係数が温度に応じて変化する流れ場においても，対数領域が存在することが示唆される．いっぽう，平均温度に関しては，𝑅𝑒𝜏 = 650, 1000 のケースにおいて，
壁面遠方で分布が一致する領域（50 < 𝑦 + < 200）が
図１ 高温，低温壁近傍の主流方向変動速度．変動見られる．この領域の平均温度の傾きは，粘性係数一速度 𝑢′ は 摩 擦 速度 𝑢𝜏 で 無 次元 化さ れ， −3 < 𝑢′ /定条件の経験式の傾きとおおむね一致しており，対数𝑢𝜏 < 3の範囲で可視化されている．赤は高速，青は域に対応していると考えられる．この対数領域において，低速領域を示している．上図は低温壁（上壁）近傍，

平均速度分布は低温壁側では下方にシフトし，高温壁

下図は高温壁（下壁）近傍の結果を示す．

側では上方にシフトする．これは，低温壁側の壁面摩58

擦係数が増加し，高温壁側で増加したことを意味する．
いっぽうで，平均温度に関しては，逆の傾向が見られ，高温壁側の平均温度は下方にシフトし，低温壁側の平均温度は上方にシフトする．つまり，高温側で壁面熱伝達率が増加し，低温側で増加したことが示唆される．実際に，壁面摩擦係数は，高温壁側では低温壁側と比べて 20％ほど小さいが，ヌセルト数は高温壁側では低温
壁側と比べて 7％ほど大きく，粘性係数の変化が速度
場・温度場に与える影響が異なっていることが確認された．

まとめ、今後の課題

粘性係数の温度依存性が乱流熱流動に与える影響を議論した．高温・低温壁の近傍では粘性変化が乱流渦構造，乱流強度に影響を及ぼし，対数速度・温度分布の上下シフトを導くことが分かった．対数領域の分布において，高温壁側の分布は上方シフトし，低温壁側の分布は下方シフトするが，温度分布では逆の傾向が見られることが分かった．今後は，さらに詳細な統計量を取得し，乱流輸送現象を明らかにする．
参考文献

図４ 平均速度・温度分布の比較．上図には速度[1] C.O. Popiel, J. Wojtkowiak, Heat Transfer

分布，下図には温度分布を示している．比較のEngineering 19 (3) (1998), 87–10.

為，粘性係数一定条件の結果として，平均速度分[2] K. Suga et al., Computers & Mathematics with

布の直接数値解析結果（５），平均温度分布の経験

Applications 69 (6) (2015), 518-529.

式（６）を示している．

[3] K. Suga et al., International Journal of Heat andFluid Flow 68 (2017), 225-236.
[4] Y. Kuwata, K. Suga, Journal of Computational
Physics 311 (2016), 348-362.
[5] K. Iwamoto, K. et al., Database of Fully Developed
Channel Flow-THTLAB Internal Report, Rapport
technique, THTLAB, Dept. of Mech. Engng., The Univ.
of Tokyo.
[6] B.A. Kader, International Journal of Heat and Mass
Transfer 24.9 (1981), 1541-1544.

59

60

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

利用課題名 GPU クラスタを用いたミリ波帯大規模広帯域電波伝搬シミュレーション
英文：Large-Scale Propagation Simulations of Ultra-Wideband Electromagnetic Fields in Millimeters and Terahertz
Wave Frequencies Using GPU Cluster
チャカロタイ ジェドヴィスノプ
Jerdvisanop Chakarothai
国立研究開発法人情報通信研究機構
National Institute of Information and Communications Technology
https://www.nict.go.jp
ミリ波帯・テラヘルツ帯における電波は他の周波数帯と比較して、電波の直進性が大きく、見通し伝搬の形態で利用されることが多い。しかし、ひとたび電波伝搬経路のフレネルゾーン内に吸収体・散乱体・反射体などの障害物が存在すると、受信電力の予測手法の複雑さが増大し、不確かさが大きくなる。特に屋内伝搬環境は、壁、天井、机、棚や椅子等の様々な障害物が存在し、伝搬状況を正しく予測することは困難である。そこで、本研究では、オフィス
等の室内環境においてミリ波帯（特に 28 GHz 帯）電波の受信電力を高精度に予測するために、大規模電磁界シミュ
レーション技術を開発した。特に実際の構造物の離散化及び機能的な電波散乱シート（入射された電波を様々な方向へ散乱させるための構造）を従来のプログラムに組み込むために、改良・拡張した。(1)屋内電波環境の数値解析のために、床・壁のみではなく、より実際のオフィスの構造に近い解析モデルにするために、3 次元構造モデル（テーブル、椅子、パソコンなど）を解析空間の解像度に合わせて離散化して組み込むことを可能にした。 (2)電波散乱シートを解析空間に組み込むために、完全磁気導体及び完全導体を同一のメッシュで取り扱えるようにするために、
FDTD 法の定式化を行い、三次元解析による定式化の妥当性検証を行った。実際に散乱板を解析モデルに組み込
み、FDTD 法によって解析した結果とモーメント法によるものと比較した。今後、散乱板を組み込んだ大規模電磁界解析を行う予定である。
Millimeter-wave and terahertz electromagnetic (EM) waves are often used in line-of-sight propagation. However, if
there are absorbers, scatterers, reflectors, or other obstacles within the Fresnel zone of the wave propagation path, theuncertainty of the received power prediction method increases significantly. Especially in indoor propagation
environments, there are various obstacles such as walls, ceilings, desks, shelves, and chairs, making accurate predictionof propagation conditions difficult. Therefore, in this study, we developed large-scale EM simulation techniques for
accurate prediction of the received power of millimeter-wave (especially 28 GHz band) waves in indoor environments.
In particular, we improved and extended conventional programs to discretize actual structures and incorporate
functional structures such as EM scattering sheet. (1) To perform numerical analysis of indoor radio wave environments,
in addition to floors and walls, we discretize 3D structural models (tables, chairs, computers, etc.) to match resolutionof the analysis space. (2) To incorporate EM scattering sheets into the analysis region, we incorporate a surface
impedance into the FDTD method to handle perfect magnetic conductors and perfect conductors on the same surface
plane and demonstrated the validity of the formulation through 3D numerical analysis. We compared the results by the
FDTD analysis which includes a scattering sheet with those obtained by the moment method. In the future, we plan toperform large-scale electromagnetic field analysis incorporating scattering sheets.
Keywords: millimeter-wave, FDTD, indoor propagation, large-scale simulation, electromagnetic scattering sheet

内電波環境においては、様々な障害物によって電波が背景と目的
5G 無線通信や Beyond 5G の次世代無線通信技術

散乱され、複雑な電波環境を形成するため、受信電力

は主にミリ波帯電波が利用される。しかしながら、ミリ波を正確に予測することが難しい。そこで、本研究では、

帯電波の伝搬経路のフレネルゾーン内に吸収体・散乱

様々な伝搬環境(特に室内環境)におけるミリ波帯電波体・反射体などの障害物が存在すると、受信電力の予の受信電力を高精度に予測するために、超大規模電測手法の複雑さが増大し、不確かさが大きくなる。また、磁界シミュレーション技術を新たに開発すると共に、機波長が非常に短いため、物体の表面形状等による散能的なの構造物(リフレクトアレーや電波散乱シート等)乱が生じ、従来のレイトレーシング法や光学近似法にに対する数値解析法を提案し、その妥当性を示す．

よる高精度な解析が困難である。特にオフィス等の屋

61

有する．

概要
本研究では、オフィス等の室内環境においてミリ波帯

まず，一次元 FDTD 法における媒質境界周辺での電（特に 5G 無線通信に使われる 28 GHz 帯）電波の受信電力を高精度に予測するために、TSUBAME3.0 GPU

界および磁界の配置を示す．電界の接線成分𝐸𝐸𝑦𝑦𝑡𝑡 を次クラスタを用いた大規模並列化電磁界シミュレーション

によって表すことができる[11]．

式のように複素周波数領域での表面インピーダンス𝑍𝑍𝑠𝑠技術を開発した。解析に用いた計算手法は、電磁界解

𝑍𝑍𝑠𝑠 (𝑠𝑠)� 𝑠𝑠𝐻𝐻𝑧𝑧𝑡𝑡 (𝑠𝑠)𝐸𝐸𝑦𝑦𝑡𝑡 (𝑠𝑠) = �𝑠𝑠析 で よ く 使 わ れ て い る 時 間 領 域 有 限 差 分 (Finitedifference time-domain, FDTD)法である[1]。

(1)

ただし，𝐻𝐻𝑧𝑧𝑡𝑡 は境界における磁界の接線成分である．式(1)を時間領域で表すと，

本研究に用いたプログラムは文献[2]-[7]にも用いられた人体ばく露評価のためのプログラムを改良・拡張し𝐸𝐸𝑦𝑦𝑡𝑡 (𝑡𝑡) = 𝑋𝑋(𝑡𝑡) ∗たものである。昨年度に引き続き、物体形状を忠実に

𝜕𝜕𝐻𝐻𝑧𝑧𝑡𝑡 (𝑡𝑡)
= 𝑋𝑋(𝑡𝑡) ∗ 𝐻𝐻𝑧𝑧′ (𝑡𝑡)
𝜕𝜕𝜕𝜕

(2)

再現できるようにするために、きめ細かくセルの各辺にとなり，X(t) = L-1{X(s)},細かく設定することができるように変更した他、解析モで表される．ここで，Zs(s)は以下のような式で与えられデルのデータをハードディスクに保存せずに、直接プロ

る．

グラムの中に組み込めるようにした。加えて、新たな機X(s) = Zs(s)/s, 𝐻𝐻𝑧𝑧′ = ∂𝐻𝐻𝑧𝑧𝑡𝑡 (t)/∂t

𝑍𝑍𝑠𝑠 (𝑠𝑠) = �

能的な構造物（リフレクトアレー）を解析空間に組み込

𝑠𝑠𝜇𝜇0 𝜇𝜇𝑟𝑟
𝑠𝑠𝜀𝜀0 𝜀𝜀𝑟𝑟 + 𝜎𝜎

(3)

むために、表面インピーダンスに対して、FDTD 法の定σ, εr, µr, ε0, µ0 はそれぞれ媒質の導電率，比誘電率，比式化を行い、提案する手法の妥当性を検討した。

透磁率，真空の誘電率および透磁率である．次のステ

本研究の解析対象としては、屋内での電波環境を可

ップでは，X(s)を FILT によって時間領域波形へ変換し視化することを目的とした。実際にテーブルや椅子等た後，Prony 法を適用すると，z 領域において X(z)を次のオフィス器具を含めた解析モデルは、大きさ 3×2.4式のように表すことができる[2],[3]．
𝑁𝑁𝑙𝑙×2.6 m であり、必要な計算資源は、約 193 億セルの3𝑋𝑋(𝑧𝑧) = �

約 720GB の総使用メモリを要した。大規模 GPU クラス

𝑙𝑙=1

𝐴𝐴𝑙𝑙
,
1 − 𝑝𝑝𝑙𝑙 𝑧𝑧 −1(4)

タにより高速化を行った結果、16 ノード(64 GPUs)によるNl は Prony 法によって得られた実数の極に対応してお計算時間は 1000 ステップ当たり、約 32 分で終えること

り，Al および pl は文献[12],[13]の手順に従って Prony 法

ができ、実用的な時間で解析が可能であることを確認

から抽出された係数である．なお，通常 Prony 法を適用した。

すると，複素数の極も得られるが，式(3)で表される表
面インピーダンスは時間領域において単調的に減衰していく関数で求まるため，実数の極のみで表現すること

表面インピーダンスの定式化ができる．式(4)を z 領域での式(2)に代入すると，本研究において，高伝導性媒質のための表面インピーダンスが適用される．表面インピーダンス法に関して

𝑁𝑁𝑙𝑙

𝑁𝑁𝑙𝑙

𝑙𝑙=1

𝑙𝑙=1

𝑡𝑡
𝐸𝐸𝑦𝑦𝑡𝑡 (𝑧𝑧) = � 𝐸𝐸𝑦𝑦,𝑙𝑙(𝑧𝑧) = �は，これまで様々な定式化手法が提案されてきたが𝐴𝐴𝑙𝑙
𝐻𝐻 ′ (𝑧𝑧)1 − 𝑝𝑝𝑙𝑙 𝑧𝑧 −1 𝑧𝑧(5)[8],[9]，本提案手法では，PEC および PMC を同一境界が得られる．式(5)より各実数の極に対応した電界の離面上で定義できることが特徴である．また，従来手法で散的な時間 n ステップ目における更新式は以下のようは，時間領域での厳密解や周波数領域でのベクトルフに求まる．

ィッティング(Vector Fitting, VF)法[10]を用いて定式化が

𝑡𝑡,𝑛𝑛
𝑡𝑡,𝑛𝑛−1
𝐸𝐸𝑦𝑦,𝑙𝑙= 𝑝𝑝𝑙𝑙 𝐸𝐸𝑦𝑦,𝑙𝑙+

行われたが，本手法は周波数領域での表面インピータ

1
1
𝐴𝐴𝑙𝑙
𝑡𝑡,𝑛𝑛+
𝑡𝑡,𝑛𝑛−
�𝐻𝐻𝑧𝑧 2 − 𝐻𝐻𝑧𝑧 2 �
∆𝑡𝑡

(6)

Prony 法を用いて時間領域でフィッティングを行うことで，𝑡𝑡,𝑛𝑛の総和で計算できる．
従って，電界の接線成分は𝐸𝐸𝑦𝑦,𝑙𝑙

高精度な解を与えながらも定式化が簡素になる利点を

以下のような更新方程式が導出される．

ンスを FILT によって時間領域波形に変換してから一方では，磁界の接線成分は前進差分を適用すると，

62

1𝑡𝑡,𝑛𝑛+

1
𝑡𝑡,𝑛𝑛−

𝐻𝐻𝑧𝑧 (𝑖𝑖+21) = 𝐻𝐻𝑧𝑧 (𝑖𝑖+21) −
2

2

2∆𝑡𝑡
�𝐸𝐸 𝑡𝑡,𝑛𝑛 1 − 𝐸𝐸𝑦𝑦𝑛𝑛 (𝑖𝑖) �
𝜇𝜇0 ∆𝑥𝑥 𝑦𝑦 (𝑖𝑖+2)

600

(7)

400

最終的に，式(5)の𝐸𝐸𝑦𝑦𝑡𝑡,𝑛𝑛(𝑖𝑖+1) を式(7)の中に代入すれば，
2371

媒質境界面における磁界の更新方程式が得られる．
1𝑡𝑡,𝑛𝑛+

1
𝑡𝑡,𝑛𝑛−

𝐻𝐻𝑧𝑧 (𝑖𝑖+21) = 𝐻𝐻𝑧𝑧 (𝑖𝑖+21) −
2

2

2∆𝑡𝑡

𝑁𝑁

𝑙𝑙
𝐴𝐴𝑙𝑙
𝜇𝜇0 ∆𝑥𝑥 + 2 ∑𝑙𝑙=1

𝑁𝑁𝑙𝑙

𝑡𝑡,𝑛𝑛−1
�𝐸𝐸𝑦𝑦𝑛𝑛 (𝑖𝑖) − � 𝑝𝑝𝑙𝑙 𝐸𝐸𝑦𝑦,𝑙𝑙
�

(8)

𝑙𝑙=1

単位: mm

これによって従来の FDTD 法の磁界の配置において，

図 1 三次元モデルのメッシュデータからボクセル化（離散化）して得られた結果

完全電気導体を表面インピーダンスによって表すことができる．本提案法の計算フローチャートを図 2 に示す．

テーブル

なお，本提案法は，高導電性媒質に限らずミリ波帯における生体組織へ適用することが可能であることが示されている[13]．一方，PMC の場合は，単純に境界面における磁界の接線成分をゼロとすることで実現でき椅子る．従って，本提案法を用いて，同一面内で PEC または PMC をモデルすることが可能である．
結果および考察
まず、図 1 のようなテーブルおよび椅子を 3 × 2.4 ×

図 2 模擬オフィス部屋内にテーブルおよび椅子を配置したときのミリ波帯(28 GHz)屋内電磁界分布2.6 m の模擬部屋内に配置して解析を行った。扉の材
3

質は木材とし、比誘電率および導電率は、それぞれ 2.0

次に、機能的な構造を有する電波散乱シートに対し及び 0.001 S/m とした。施工ボードはεr = 2.4, σ = 0.001て解析を行う。図 3 に異なる高さを有する 2 種類の金属

S/m、床面は、εr = 13, σ = 0.01 S/m の厚さ 10 cm のコン板で構成されている EMSS を示す。基準面における反クリートとした。解析モデルの解像度は、1 mm で、8 層

射位相は、PMC が 0o で，PEC が-180o となっている。0oの完全整合層 (Perfectly matched layers)を含めた解析の反射位相は、金属板の高さをλ/4 分だけ変えることで領域は、3,024 × 2,424 × 2,624 セルである。解析周実現される。また、図 3 には、実際の EMSS の構造(モ

波数は 28 GHz である。また、入射はモーメント法によっーメント法における解析モデル)と FDTD シミュレーショ

て求めたホーンアンテナの開口分布を波源とした。ホーンで使用される解析モデルを示している。解析周波数

ンアンテナを用いた場合、開口面分布は磁流のみで表を 30 GHz と し 、 EMSS の 全 体 サ イ ズ は 20×20

現できるが、任意波源の場合は電流及び磁流の両方

cm2(20λ×20λ)であり、金属板の幅は 2 cm (2λである．

をモデリングする必要がある。図 2 には、テーブルや椅

基準面での反射位相が-180o の金属板は、導電率σ =子を含む模擬部屋の外に波源（K バンドの標準ホーン

58×1010 S/m を有する金属の表面インピーダンスで置きアンテナ）を設置したときの電界分布を示す。図 2 より換えられ、反射位相がゼロの板は完全磁気導体(PMC)

28 GHz 帯電波はスタッドの間を透過することが分かり、で置き換えられている（磁界の接線成分がゼロ）。解析

ホーンアンテナからは鋭いビームが形成されていること

領域のサイズは 241×241×61，セルサイズは∆x = ∆y =を確認できた。結果から、波源が 1.5 m の高さに配置さ∆z = 1 mm、時間ステップ間隔は 1.67 ps とした。入射波れているため、高さ約 37 cm のテーブルおよび椅子か形は連続正弦波とし、20λ×20λ の大きさを有する EMSS

らの散乱は少ないことが確認された。

に対する散乱断面積(Bistatic Scattering Cross Section,
BSCS)を導出するための計算時間はわずか 23 秒であった。図 6 は垂直入射時(θinc = 0o)の表面インピーダン

63

ある。

Ex
Hy

2 cm

謝辞 本研究を進めるあたり，表面インピーダンス法にk関して検討して頂いた法政大学大学院修了生の竹谷和真氏に感謝の意を表す．
20 cm
z
20 cm

PEC

Hy

x

Ex
k

参考文献

y

[1] K. S. Yee, “Numerical solution of initial boundary

PMC

value problems involving Maxwell’s equations inisotropic media,” IEEE Trans. Antennas Propag., vol. AP14, no. 3, pp. 302–307, May 1966.
[2] J. Chakarothai, K. Wake, S. Watanabe, “Scalable GPUparallelized FDTD method for analysis of large-scale

図 3 電波散乱シート(EMSS)の実際の解析モデルおよび表面イン
ピーダンスを用いた場合の FDTD 解析モデル[14]= ∆y = ∆z = 1 mm
= 1.67 ps
FDTD analysis

0oEincelectromagnetic

dosimetryproblems,”

AppliedComputational Electromagnetics (ACES) Journal, vol. 31,
no. 6, pp. 661-668, June 2016.

MoM
FDTD

[3] J. Chakarothai, “Novel FDTD scheme for analysis
of frequency-dependent medium using fast inverseθLaplace transform and Prony’s method,” IEEE Trans.
Antennas Propagat., vol. 67, no. 9, pp. 6076-6089, Sep.
2019.

90o
-20 -10 0 10 20 [dBsm]
図 4 垂直入射時の電波散乱シートの散乱特性(実線：MoM, 破線:
FDTD)[14]90o[4] J. Chakarothai, S. Watanabe, K. Wake, “Numerical
dosimetry of electromagnetic pulse exposures using
FDTD method,” IEEE Trans. Antennas Propagat., vol.

スを使った FDTD 解析で得られた BSCS およびモーメ

66, no. 10, pp. 5397-5408, Oct. 2018.

[5] チャカロタイジェドヴィスノプ, 和氣加奈子,

ント法による電磁界解析で得られた BSCS を示している。

渡辺聡一, 陳強, 澤谷邦男, “超広帯域電磁界解析

両者は非常によく一致しており、提案手法の妥当性がのための周波数依存性 FDTD 法,” 電子情報通信確認された。

学会 和文論文 C, vol. J102-C, no. 5, pp. 102-113,
May 2019.

まとめ、今後の課題[6] J. Chakarothai, K. Fujii, “A unified approach for

本研究では、28 GHz における室内電波伝搬環境を

treatment of frequency-dependent materials in FDTD

高精度に解析できる時間領域有限差分法を開発した。

method,” ISAP 2019, Xi’an, China, Oct. 30, 2019.

本手法では、これまでのプログラムに対して、(1)より実

[7]

J.

Chakarothai,

K.

Wake,

S.

Watanabe,

際のオフィスの構造に近い解析モデルにするために、3“Convergence of a single-frequency FDTD solution in次元構造モデル（テーブル、椅子、パソコンなど）を解析

numerical dosimetry,” IEEE Trans. Microw. Theory

空間の解像度に合わせて離散化して組み込むことを可Techn., vol. 64, no. 3, pp. 707-714, March 2016.
[8] J. G. Maloney and G. S. Smith, “The use of surface

能にした。(2) 機能的な構造体（電波散乱シート等）をimpedance concepts in the finite-difference time-

解析空間に組み込むために、FDTD 法の定式化を行い、domain method,” IEEE Trans. Antennas Propag., vol.

三次元解析による定式化の妥当性検証を行った。モー40, no. 1, pp. 38-48, 1992.

メント法による解析結果と比較して、同等な解析精度が[9] J. H. Beggs, R. J. Luebbers, K. S. Yee, and K. S.
Kunz, “Finite-difference time-domain implementation

得られることが示された。

of surface impedance boundary conditions,” IEEE

今後は、電波散乱シートやリフレクトアレー等を実際Trans. Antennas Propag., vol. 40, pp. 49-56, 1992.

に組み込んだ大規模電波伝搬解析へ拡張する予定で

64

[10] B. Gustavsen and A. Semlyen, “Rational
approximation of frequency domain responses by
vector fitting,” IEEE Trans. Power Delivery, vol. 14, no.
3, pp. 1052-1061, July 1999.
[11] M. A. Leontovich, “On the approximate boundary
conditions for electromagnetic fields on the surface of
wellconducting

bodies,”in

Invest.

Propag.

Radiowaves, B. A. Vedensky, Ed. Moscow: Academyof Science USSR, pp. 5-20, 1948.
[12] J. Chakarothai, K. Fujii, Y. Suzuki, J. Shibayama,
and K. Wake, “Analyses of transient energy depositionin biological bodies exposed to electromagnetic pulses
using parameter extraction method,” IEICE Trans.
Commun., vol. E105-B, no. 6, pp. 694-706, June 2022.
[13] K. Takeya, J. Chakarothai, J. Shibayama, Y.
Suzuki, and K. Fujii, “Novel surface impedance
formulation via FILT and Prony method for FDTD
analyses of lossy media,” IEEE Antennas and WirelessPropagation Letters, vol. 23, no. 5, pp. 1593-1597, May
2024.
[14] J. Chakarothai, K. Takeya, J. Shibayama, Y.
Suzuki, K. Fujii, “FDTD Analysis of Electromagnetic
Scattering Sheet Using Surface Impedance Method,”2024 IEEE International Workshop on Antenna
Technology (iWAT2024), April 15-18 2024 (Sendai,
Japan).

65

66

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

利用課題名 ベイナイト変態中の準安定オーステナイトにおける固溶元素と炭素原子の相関の解明
英文：Correlation between solid solution elements and carbon atoms in metastable austenite during
bainitic transformation
利用課題責任者 永野 隆敏
Takatoshi Nagano
所属 茨城大学大学院理工学研究科 工学野物質科学工学領域

邦文抄録（300 字程度）

Fe-Si-C 合金中の炭素の挙動と、炭素の濃度変化による結晶構造への影響を第一原理計算の VASP によって算出する。その結果によって Fe-Si-C 合金中の準安定オーステナイト相の炭素濃度が高濃度化するメカ
ニズムを明らかにした。また炭素濃度の違いによって変態の優位性が変わることも明らかにした。オーステナイト
中の炭素濃度が異なる理由として、オーステナイト中にケイ素が固溶しているとき炭素濃度が 0.4 wt%付近まではケイ素が固溶していないときよりも安定することが挙げられる。このことからオーステナイトに固溶しているケイ素近傍に炭素がトラップされ、炭素濃度の差が生じたと考えられる。変態の優位性の差については高炭素濃度であるオーステナイトはケイ素が固溶している可能性が高く、このケイ素が変態を阻害したと考えることで説明できる。
英文抄録（100 words 程度）

The behavior of carbon in Fe-Si-C alloys and the effect of carbon concentration changes on the
crystal structure are calculated by first-principles VASP calculations. The results reveal themechanism of the high carbon concentration in the metastable austenite phase in Fe-Si-C alloys.
It is also clarified that the predominance of transformation changes with carbon concentration.
The reason for the difference in carbon concentration in austenite is that when silicon is solidly
soluble in austenite, the carbon concentration is more stable up to around 0.4 wt% than whensilicon is not solidly soluble in austenite. This suggests that carbon is trapped in the vicinity ofsilicon solid solution in austenite, resulting in the difference in carbon concentration. The
difference in the predominance of transformation can be explained by considering that austenitewith a high carbon concentration is likely to have silicon in solid solution, and that this siliconinhibits transformation.
Keywords: ５つ程度
熱処理 高炭素濃度オーステナイト ベイナイト変態 セメンタイト 鉄中のケイ素固溶 up-hill diffusion
背景と目的
概要

次世代の強度-延性を兼ね備えた鉄鋼材料として、
準安定オーステナイト（FCC）相を含む複相組織が有力

超高強度、高延性鋼板を開発することは建物や自動視されている。その一つにベイナイト組織を利用した鋼車をはじめ鉄鋼を用いた工業製品の性能を向上させるの開発がある。この組織形成には炭素濃度コントロール

ために重要である。そのために鉄鋼中炭素位置と固溶が重要であるがことがわかっているが、鉄鋼中炭素位置

元素の相関を解明する必要がある。

と固溶元素の相関を解明するには実験的手法だけではFe-Si-Mn-C 鋼を用いたこれまでの研究で、全量γ

不明な点が残ってしまう。そこで第一原理計算である

相となる温度からの急冷によるオーステンパー処理中

VASP によって、炭素との相関を示す元素とその程度をでのベイナイト変態進行を観察し、処理温度を 573-

比較しメカニズムを明らかにし、その金属組織のコントロ

773K としたすべての条件でγ中炭素濃度の二極化がールをより容易にすることを目的とする。

観測された。この炭素濃度二極化のメカニズム解明とし67

て、固溶元素と炭素原子との相関を第一原理計算によ多い別の物質についても同様の評価を行う必要があって明らかにする。

る。

結果および考察

図 1 オーステナイト中の C 濃度変化による安定性の変化
シミュレーションによって以下のことがわかった。
・オーステナイト中を炭素が拡散するとき、ケイ素が固溶していると炭素の経路に指向性が生じる。
・オーステナイト中の炭素濃度が変化するとき、炭素濃度が 一定値までケイ素が固溶しているときの方が安定している(図 1)。
・セメンタイトはケイ素の固溶が少ないほど安定する。
・セメンタイト中、炭素濃度が小さいときはケイ素が固溶しているときの方が安定する。
オーステナイト中の炭素濃度が異なる理由として、オーステナイト中にケイ素と炭素の相関が考えられる。
低炭素濃度オーステナイトから先に変態が行われ、その後高炭素濃度オーステナイトが変態を行うという現象について、このケイ素が変態を阻害したと考えられる。
まとめ、今後の課題

第一原理計算の VASP によって Fe-Si-C 合金
中の準安定オーステナイト相の炭素濃度が高濃度化するメカニズムを明らかにした。また炭素濃度の違いによって変態の優位性が変わることも明らかにした。
今後はケイ素だけでなく鉄鋼材料に含有する機会が68

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

人工画像を用いた大規模事前学習
Large Scale Pretraining Using Synthetic Images
横田 理央
産業技術総合研究所 デジタルアーキテクチャ研究センター

邦文抄録（300 字程度）
これまでの研究では学習中に tar ファイルから直接画像を読む方法により inode 数を 1/1000 に削減、パラメータ
分散・再計算・オフロードの機構を導入することでメモリ消費をノードに反比例して削減する方法を開発してきた。
本課題では、これらの技術を用いて初めて可能になる超巨大な vision transformer の事前学習を実施した。ただし、このような事前学習に必要な画像データは JFT-300M/3B 規模であり、それらのデータセットは非公開である。
そこで、本申請課題では Visual Atom や Newton Fractal などの人工画像を用いて事前学習を行った。その結果、人工画像を用いた事前学習を行うことで ImageNet-21k を超える事前学習性能を達成した。
英文抄録（100 words 程度）
Previous studies have developed methods to reduce the number of inodes to 1/1000 by reading images
directly from tar files during training, and to reduce memory consumption inversely proportional tonodes by introducing parameter distribution, recalculation, and offloading mechanisms. In thisassignment, we conducted pre-training of a very large vision transformer, which is possible only byusing these techniques. However, the image data required for such pre-training is on the scale of
JFT-300M/3B, and these data sets are not publicly available. Therefore, in this proposal, pre-training
was conducted using artificial images such as Visual Atom and Newton Fractal. As a result, weachieved a pre-training performance exceeding ImageNet-21k by pre-training with artificial images.

Keywords: 深層学習、事前学習、Vision Transformer、人工画像、大規模データ

背景と目的

表的な深層ニューラルネットの変遷を示す。画像処理

深層学習技術の中で格段に高い精度を発揮している分野では LeNet が畳み込みニューラルネットを採用すのが超巨大な transformer の事前学習である。本課

ることで精度を向上させ、AlexNet は GPU を用いるこ

題はこのような巨大な transformer の事前学習に必要とで高速にこれを処理できるようにし、ResNet ではスになる膨大なデータを数式から人工的に作った画像でキップ接続とバッチ正規化により学習を容易にした。ま代替する革新的な基盤技術を提案する。既にた 、 MobileNet の squeeze and excite 機 構 やImageNet-21k に関しては同等のデータ量で同じ精度

EfficientNet の neural architecture search を用いてを達成しており、数式から無限に生成できる人工画像

効率的なニューラルネットが設計できるようになり、精で自然画像と同程度の事前学習効果が得られるという度を維持しつつ小型化できるようになった。一方、自然

事実は、深層学習分野に革命をもたらすものであると

言語処理分野では LSTM などの再帰的ニューラルネッ

予想される。

ト が 長 く 用 い ら れ て い た が 、 2017 年 に 登 場 し たtransformer により大幅な性能向上が実現され、現在

概要

では transformer が主流なニューラルネットとなってい本課題で用いる計算モデルは深層ニューラルネットのる。この transformer を画像処理分野で使えるようにし一種である vision transformer である。深層ニューラたものが本課題で用いる vision transformer である。

ルネットはこの 20 年間で目覚ましい進歩を遂げてきた。

LeNet から EfficientNet に至るまで画像処理分野で

図１に初期の LeNet や LSTM などの深層ニューラル

用いられてきたニューラルネットは、全て畳み込みニュ

ネットから現在の Vision Transformer に至るまでの代

ーラルネットであったが、vision transformer は畳み込69

みのような機構を予め人手で組み込むことなくデータかの更新を行わずに勾配を累積するサブステップを設けらそのような構造をも学習することができるため、大量ている。左図に示す更新ステップあたりの計算時間がノのデータがある場合には優位になる。ただし、最も大きード数に反比例して減少しているのはこのサブステップ

い実画像データセット JFT を Google が非公開としていの数がノード数に反比例して減少しているためである。

ることが画像処理分野全体の発展にとって大きな障壁

Original は最適化前の計算時間で Optimized は最適となっている。本課題では、自動的に生成・ラベル付け

化後の計算時間を表す。右図に示すように 1024 ノード

で き る フ ラ ク タ ル な ど の 人 工 画 像 を 用 い て visionを用いたときの最適化後の演算性能は 368.6TFLOPs

transformer の事前学習を行い、ImageNet-21k やとなっている。「富岳」のノードあたりの単精度の理論演

JFT-300M などの大規模自然画像での事前学習効果

算性能が 6.8 TFLOPs であることを考えると改善の余と比較する。

地が大いにあることが分かる。Transformer の演算の9 割以上は行列積となっているため、高い理論演算性
能比が期待できるが、現状では大きな正方行列でないと高い演算性能がでないため、transformer が必要とする行列の次元で性能最適化を行う必要がある。また、Python から行列積のカーネルを呼び出すオーバーヘ
ッドも極力小さくする必要があるため、バッチ行列積のカーネルを用いたり、カーネルフュージョンを行う必要がある。

図１ 主要な深層ニューラルネットモデルの変遷

結果および考察
本課題では、vision transformer に対して４種類の異なる分散並列化を行った。１つ目は、データをそれぞれのプロセスに分散させるデータ並列。２つ目は、層内

図２ データ並列による強スケーリングのテンソルを分散させるテンソル並列。３つ目は、パラメータを層ごとに分散し層間でパイプライン化を行うパイプライン並列。４つ目は、勾配を AllReduce する代わまとめ、今後の課題

りに ReduceScatter しておき、forward や backward

本課題では、これまでの継続元課題で培ってきた「富で必要になったときだけ AllGather する手法の ZeRO岳」上で大規模深層学習を行うための超分散並列・省である。これらの手法を「富岳」上で動作させることは継メモリ・省 I/O の技術を集積し、画像処理分野の最大の続元課題で既に行っており、今年度はその技術を応用課題となっているデータセットの不足・寡占の問題を、して vision transformer による大規模事前学習を行っ人工画像データセットによって解決する方法を提案するた。

ことができた。これまでは、インターネットからスクレー

図２にデータ並列数のみを変化させたときの強スケーピングされた写真などの画像を用いて学習することがリングの様子を示す。左図は更新ステップあたりの計

最 も 有 効 と さ れ て き た が 、 本 課 題 で は vision

算時間、右図はそのときの総 TFLOPs を示す。この実transformer のような最新のニューラルネットでは自然

験では、ノード数によらずグローバルバッチサイズを固画像が必ずしも必要でないことを、このような大きな規定するために、1024 ノード以外の実験ではパラメータ

模で世界で初めて示すことができた。ただし、図２の結70

果からも分かるように並列化効率は良いものの、ノード単体性能は「富岳」の理論演算性能の 1/20 程度の実効性能しかでておらず、今後これを向上させることが最重要課題である。

71

72

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

人工画像データ及び事前学習モデルのリポジトリ化による AI ハブ構築
AI Hub of Synthetic Image Datasets and Pre-trained Models
谷村 勇輔
産業技術総合研究所 デジタルアーキテクチャ研究センター

邦文抄録（300 字程度）
大規模学習を高速に実行可能にする次世代システムの設計に向けて、アクセラレータの効果的な利用や大規模並列実行における通信の最適化、ストレージ I/O の効率化等に関する手法の研究開発を行っており、TSUBAME
にてその評価検証を行った。特に、学習中に tar ファイルから直接画像を読む方法により inode 数を 1/1000 に削減した。

英文抄録（100 words 程度）
Toward the design of next-generation systems that enable high-speed execution of large-scale learning,
we have been conducting research and development of methods for effective use of accelerators,
optimization of communication in large-scale parallel execution, and efficiency improvement of storage
I/O, etc. We evaluated and verified these methods in TSUBAME. In particular, the number of inodes
was reduced to 1/1000 by reading images directly from tar files during training.

Keywords: 深層学習、事前学習、Vision Transformer、人工画像、大規模データ

背景と目的ることで精度を向上させ、AlexNet は GPU を用いるこ

深層学習技術の中で格段に高い精度を発揮しているとで高速にこれを処理できるようにし、ResNet ではスのが超巨大な transformer の事前学習である。本課

キップ接続とバッチ正規化により学習を容易にした。ま題はこのような巨大な transformer の事前学習に必要た 、 MobileNet の squeeze and excite 機 構 やになる膨大なデータを数式から人工的に作った画像でEfficientNet の neural architecture search を用いて代替する革新的な基盤技術を提案する。既に効率的なニューラルネットが設計できるようになり、精ImageNet-21k に関しては同等のデータ量で同じ精度度を維持しつつ小型化できるようになった。一方、自然を達成しており、数式から無限に生成できる人工画像

言語処理分野では LSTM などの再帰的ニューラルネッ

で自然画像と同程度の事前学習効果が得られるというト が 長 く 用 い ら れ て い た が 、 2017 年 に 登 場 し た

事実は、深層学習分野に革命をもたらすものであると

transformer により大幅な性能向上が実現され、現在予想される。

では transformer が主流なニューラルネットとなっている。この transformer を画像処理分野で使えるようにし概要たものが本課題で用いる vision transformer である。

本課題で用いる計算モデルは深層ニューラルネットの

LeNet から EfficientNet に至るまで画像処理分野で一種である vision transformer である。深層ニューラ

用いられてきたニューラルネットは、全て畳み込みニュ

ルネットはこの 20 年間で目覚ましい進歩を遂げてきた。

ーラルネットであったが、vision transformer は畳み込

図１に初期の LeNet や LSTM などの深層ニューラル

みのような機構を予め人手で組み込むことなくデータか

ネットから現在の Vision Transformer に至るまでの代らそのような構造をも学習することができるため、大量

表的な深層ニューラルネットの変遷を示す。画像処理のデータがある場合には優位になる。ただし、最も大き分野では LeNet が畳み込みニューラルネットを採用す

い実画像データセット JFT を Google が非公開としてい73

ることが画像処理分野全体の発展にとって大きな障壁で と し た 。 デ ー タ セ ッ ト は FractalDB-100k ととなっている。本課題では、大規模学習を高速に実行

FractalDB-300k を中心に事前学習し、余ったノードで可能にする次世代システムの設計に向けて、アクセラ

より小規模なデータセットの学習を行った。また、レータの効果的な利用や大規模並列実行における通FractalDB-100k では 10 epoch、FractalDB-300k で信の最適化、ストレージ I/O の効率化等に関する手法

は 3 epoch だけの学習を行った。これらのデータセットの研究開発を行っており、TSUBAME にてその評価

での事前学習は初めてであったことと、事前の検証を行った。特に、学習中に tar ファイルから直接画

FractalDB-50k までの profile では Webdatasets を

像を読む方法により inode 数を 1/1000 に削減した。

用いなくとも I/O の時間がほとんど影響しないという計測 結 果 が 出 て い た た め 、 安 全 を と っ てFractalDB-100k や FractalDB-300k の事前学習を
Webdatasets を使わずにまずは行った。しかし、3 時間
待っても FractralDB-100k の事前学習すら全く始まらず最初の画像読込に膨大な時間がかかっていた。後にtorch.profile の計測の仕方が間違っていることが判明
し FractalDB-10k などでも I/O の影響は非常に大きかったことが分かった。本実施中にはこのことは判明して
い な か っ た が 、 Webdatasets か ら FractalDB100k/300k を読み込む方式に切替えて事前学習を行

図１ 主要な深層ニューラルネットモデルの変遷

った。

結果および考察

まとめ、今後の課題

JFT-300M 規模に相当する FractalDB-100k (1 億画

少なくとも ViT-Tiny/Small/Base の 3 種類のモデル像)及び FractalDB-300k (3 億画像)のデータセットをと FractalDB-10k/21k/50k と ImageNet-21k の 4 種生成する必要があるが、Fractal 画像を生成するには類のデータセットの全ての組み合わせにおける事前学

ディスプレイに接続されている GPU を用いるのが最も習を十分な epoch まで完了することができた。これまで高速であるため、研究室のサーバで全ての画像を生成の FractalDB を用いた ViT の事前学習の研究例では、

し TSUBAME に 転 送 す る こ と に し た 。 た だ し 、ImageNet-1k で の 事 前 学 習 と CIFAR10 で のFractalDB-300k を転送のために tar ファイルに固め

fine-tuning にとどまっていたが、本チャレンジではそるのに 2 週間ほどかかった。さらに、PyTorch の標準のの 20 倍程度の規模となる ImageNet-21k による事前

DataLoader から読み込む場合、3 億画像をファイルシ

学習と ImageNet-1k による fine-tuning を行うことがステムからアクセスすることになり、I/O がボトルネックできた。Fine-tuning は ImageNet-1k だけでなく、となる。Webdatasets を用いることで、少数の巨大な

CIFAR100 や CIFAR10 でも行った。興味深い結果とtar ファイルに画像を固めて、tar から直接読むことがでし て は 、 CIFAR10 や CIFAR100 な ど で はきる。ただし、Webdatasets を用いた場合、ノード間シ

FractralDB-21k

ャッフルは起きず、ノード内のシャッフルのみが可能でImageNet-21k で事前学習したものよりも精度が低かある 。こ の ことが 学 習に ど のよ う な影 響が 出る かをったのに対して、最も困難とされる ImageNet-1k のImageNet-1k の事前学習を用いて検証したところ、学

fine-tuning

習曲線が一致したため、問題ないと判断した。

ImageNet-21k よ り も 高 精 度 に な っ た 。 さ ら に 、以上の様々な問題により、本実施では当初計画を縮ViT-Tiny や ViT-Small な ど の 小 さ い モ デ ル で は小し、モデルサイズは ViT-Large ではなく ViT-Base まImageNet-21k の方が FractalDB-21k よりも精度が74

で 事 前 学 習 し た も の がで はFractalDB-21kの 方 がよ く な っ た が 、 最 も 大 き い ViT-Base で はFractalDB-21k の方が ImageNet-1k の fine-tuning
において高い精度を達成した。つまり、fine-tuning タ
スクが困難であるほど、また、ViT の規模が大きいほど
FractalDB の ImageNet に対する優位性は向上するという興味深い結果となった。

75

76

TSUBAME 共同利用 令和 4 年度 学術利用 成果報告書

利用課題名 集風レンズ付き風車とそのマルチロータシステムの流体シミュレーション
英文：CFD Simulations of Diffuser Augmented Wind Turbines and Multirotor Systems
胡 長洪Changhong Hu
九州大学 応用力学研究所
Research Institute for Applied Mechanics, Kyushu University
URL https://www.riam.kyushu-u.ac.jp/
邦文抄録（300 字程度）
高い発電性能を持つ集風レンズ付き風車が注目されており，定格出力 200kW の中型レンズ風車の技術開発を進めている．本利用課題では，数値流体シミュレーションにより集風レンズ付き風車の発電性能と風抵抗を評価し，中型機用の集風レンズ形状の空力設計を行った．実機スケールのシミュレーション手法を開発し，模型スケールと実機スケールの風車の発電性能を評価した．また，集風レンズ付き風車を同一構造上に 2 基配置するマルチロータシス
テムに対し流体シミュレーションを行い，風車配置が発電性能に与える影響を調べた．
英文抄録（100 words 程度）
A wind lens turbine (WLT) has been attracting attention for its high-power output efficiency, and weare developing a medium-sized WLT with a rated output of 200 kW. In this project, we evaluated thepower and drag coefficients of WLTs by CFD and aerodynamically designed the diffuser shape for the
200kW WLT. A full-scale simulation method was developed to evaluate the scale effect on the power
coefficient of the WLT. CFD was performed on a multi-rotor system with two WLTs to investigate therelationship between turbine spacing and power gain.

Keywords: Diffuser augmented wind turbine, Aero dynamics, Multi rotor system, Lattice Boltzmann
method, GPU

背景と目的

概要

高い発電性能を持つ風車として集風レンズ付き風車

本研究課題の 3 つの研究目的に対し，以下の内容をが注目されており，九州大学応用力学研究所では定格実施した．レンズ風車の流体計算には GPU と大規模

出力 200kW の中型レンズ風車の開発を進めている．

計算に適した格子ボルツマン法を利用した．

レンズ風車の発電性能と風抵抗は，集風レンズの形状① ベイズ最適化のソフトウェア PHYSBO を用いた多に大きく影響される．レンズ風車の大型化，複数のロー目的最適化により，高発電性能で低抵抗となる集タを同一構造上に配置するマルチロータシステム化で風レンズ形状の探索を行った．CFD により発電性

は，集風レンズに掛かる大きな風荷重が課題となる．

能と風抵抗を評価した．

本研究課題では，② 実機スケールレンズ風車のブレードに対する空力① 数値流体解析により高い発電性能を維持したまま特性を QBlade により算出した．そのデータを利用

低抵抗となる集風レンズ形状を提案することしフルスケールの CFD を行い，実験機スケールと② 風洞試験が実施できない実機スケールに対する流結果を比較した．
③ 2 基のレンズ風車の隙間間隔を変えたシミュレーシ

体シミュレーション手法を確立し，レイノルズ数の影
響を明らかにすること

ョンを行い，単基の場合からの発電量の増加率と③ 2 基レンズ風車から構成するマルチロータシステム

離間距離の簡易を調べた．

の風車離間間隔が発電性能に影響を調べ，適切な風車配置を CFD により検討することを目的とする．
77

表 1 模型スケールと実機スケールのパワー係数 Cp

結果および考察
① 集風レンズ形状の最適化
発電性能 Cp*（レンズ直径基準）の最大化と抗力係

Ci タイプCi タイプCii タイプ

Cii タイプ

つば 3%

つば 5%

つば 3%

つば 5%

数 Cd の最小化を目的とし，集風レンズ形状のベイズ最

模型

0.683

0.709

0.739

0.794

適化を行った．入力変数は，入口位置，出口位置，曲

実機

0.731

0.764

0.786

0.826

率，つば高さの 4 つとした．最適化で提案された 64 個のレンズ形状（赤：BO）を用いた場合の Cp と Cd のプ③ 2 基マルチロータシステムのシミュレーション

ロットを図 1 に示し，従来型レンズ形状

従来型のレンズ風車 Ci タイプを横に 2 基並べた場合（青:Conventional）と初期データセットに用いた乱数をの発電性能 Cp を CFD により調べた．風車間の隙間は，用いて生成した形状（緑:Random）の結果と比較する．

集風レンズ直径の 0%, 10%, 20%, 50%, 100%, 200%,

プロットエリアの右下に行くほど高性能で低抵抗となり，300%の 7 ケースとした．単基風車と比較した場合の各最適化で提案された形状は従来型よりも良い結果を示

風車の Cp の増加率を図 2 に示す．2 基の風車を近づした．

けて（離間距離が直径の 50%ほど）配置した場合，各風車の性能は単基よりも 1%ほど向上することが分かった．2 基風車が離れるほどマルチロータの相乗効果は小さくなった．

図 2 2 基マルチロータの発電量の増加率
図 1 ベイズ最適化で得られた解（評価値：パワー係数

まとめ、今後の課題

Cp*と抗力係数 Cd）のプロット

本研究課題では，レンズ風車に対して①集風レンズ
形状の最適化，②実機スケールのシミュレーション，③② 実機スケールのシミュレーション

2 基マルチロータシステムのシミュレーションを実施し格子ボルツマン法による CFD を用い，従来型のレンた．

ズ風車 Ci タイプと Cii タイプに対して，模型スケール（ロベイズ最適化と CFD を組み合わせた空力設計によーター直径 1m）と実機スケール（ローター直径 21m）でり，従来型のレンズ風車よりも高性能で低抵抗となる集発電性能 Cp の評価した．その結果を表 1 に示し，実機風レンズ形状を提案した．実機スケールでは模型スケスケールでは模型スケールよりも発電性能が高くなる

ールよりも発電性能が高くなる傾向が確認された．レン

ことが分かった．

ズ風車のマルチロータシステムでは，2 基を近づけて配置することで発電性能が向上することが分かった．
今後，実機スケールレンズ風車と 2 基マルチロータシ
ステムに対する空力弾性解析手法の開発を行う．

78

東京工業大学 TSUBAME 共同利用
令和４年度利用終了課題 利用成果報告書集

発行 ： 令和６年 ７月
国立大学法人 東京工業大学 学術国際情報センター 共同利用推進室
住所 ： 〒152-8550 東京都目黒区大岡山 2-12-1 E2-6
E-mail: kyoyo@gsic.titech.ac.jp
URL: https://www.gsic.titech.ac.jp/tsubame
本書に記載の記事・写真等の二次利用を禁じます。これらの情報は著作権法上認められた「私的利用」または「引用」の条件をみたした場合を除いて、著作権者に無断で転載、複製、
放送、公衆送信、翻訳、販売、賃与等の利用を禁じます。