令和４年度高性能汎用計算機高度利用事業「富岳」成果創出加速プログラム「「富岳」が拓く Society 5.0 時代のスマートデザイン」
成果報告書

令和５年５月３０日
国立研究開発法人理化学研究所
坪倉

誠

目次
１．

補助事業の目的 .......................................................................................................................... - 1 -２．令和４年度（報告年度）の実施内容 ................................................................................................ - 1 ２－１．当該年度（令和４年度）の事業実施計画 .............................................................................. - 1 ２－２．実施内容（成果） .................................................................................................................. - 3 ２－３．活動（研究会の活動等） ..................................................................................................... - 64 ２－４．実施体制 ............................................................................................................................... - 69 --2-

補助事業の名称「富岳」成果創出加速プログラム「富岳」が拓く Society 5.0 時代のスマートデザイン

１．補助事業の目的
少子高齢化に伴う労働生産性の低下、確実に来襲する新興ウイルスによるパンデミックや予期せぬ自然災害への対応、さらにはカーボンニュートラル社会の実現等、激変する社会に対応して今後も我が国がものづくり大国としてその国際競争力を維持していくためには、単なる性能試験の代替としてのシミュ
レーションを脱却し、新たなデジタルエンジニアリング技術を創生することが急務である。ここでは、「富岳」の性能を十二分に活用できるアプリケーションソフトウェアによる高性能大規模シミュレーシ
ョンと、多目的最適化、機械学習、データ同化等のデータ科学を融合させることで、サイバー空間とフィジカル空間を高度に結合させた次世代デジタル設計システムを構築する。即ち、深層学習を用いたビッグデータ解析とは異なり、限定されたデータ群に物理シミュレーション技術を効率的に統合することで、演繹的なシミュレーション技術と帰納的な機械学習技術を高度に融合した最適設計システムを開発する。
開発したシステムは、我が国の基幹産業分野を構成する自動車・機械・土木・建築分野に対して以下の四つの課題を選定し、その有用性を実証することで成果創出を狙う。構築したシステムは、我が国の基幹産
業分野を構成する自動車･機械･土木･建築分野の応用課題への展開を目指し、自動車空力設計、都市・建築性能設計、室内環境設計、カーボンフリーガスタービン設計に適用することで、その有用性を実証する。研究開発はそれぞれのターゲットに対応した産学コンソーシアム（既に設置済み）を活用し、産業界での実用化を想定したシステムの設計から実証、さらには産業界メンバーへの技術教育を通した HPC 人
材育成を経て、社会実装までを連携して行う。これにより、我が国を代表する基幹産業分野において、Society 5.0 時代のスマートデザインを実現する。

２．令和４年度（報告年度）の実施内容
２－１．当該年度（令和４年度）の事業実施計画
我が国の基盤産業を構成する機械･材料･土木･建築分野の実用的な課題より、特に経済波及効果の大き
い四つの実証テーマを設定し、大学等の研究者と産業界のエンジニアが強力に連携して、課題を推進する。以下、それぞれのテーマの本年度の実施内容の詳細を示す。
(1) サブ課題 A「意匠空間を考慮した AI 支援多目的最適化による自動車空力デザイン」
性能設計空間と意匠デザイン空間を融合させた新たな空力多目的最適化システムを構築し、実際の自動車空力設計に適用することで、その有用性を実証する。多次元設計空間での最適化を実現するために、遺伝的アルゴリズムに基づく多目的最適化を自動車設計で実用化するフレームワークを構築する。キャパシティコンピューティングによる多目的最適化で問題となる計算負荷低減技術として、機械学習によりサロゲート（代理）モデルやリダクション（縮約）モデルを構築し、多目的設計の大幅な時間短縮を目指す。
本年度は、AI 支援多目的最適化フレームワークの設計、ニューラルネットワークの「富岳」チューニ-1-

ング、基礎流れを対象とした数億セル規模の縮約モデルと実車を対象とした数千万セル規模の縮約モデ
ルの構築技術の開発、多目的最適化実用モデルの構築と解析の実施、及び意匠空間の数理モデル化と空力
性能空間との融合技術開発を行う。
（協力機関：神戸大学、広島大学）(2) サブ課題 B「変容する都市・建築の自然擾乱対応の性能設計」
データサイエンスに基づきメタボリズム（新陳代謝）的に変容する都市・建築モデルをサイバー空間上に再現し、都市・建築の離散的空間での換気性能、グリーンインフラに基づく気温調整機能を定量化する。気象モデルと乱流シミュレーションを融合した広範囲の周波数帯の乱流変動を再現し、極端気象（スーパータイフーン、竜巻、豪雨）下での都市・建築の挙動を把握することで防災性能を予測する。Society
5.0 の理念を踏襲しながら、自然起源の擾乱の都市・建築へのインパクト推定技術の社会貢献をめざし、コンソーシアムを基本として、産業界を巻き込んで実務的展開と不確定性に基づく性能設計の理念と方法を構築する。
本年度は、広域市街地解析に基づく流入条件データベースの大サンプルデータ作成、風・水・構造体の統一解法のためのモデル開発、風・熱・湿度の連成解法の基礎的なモデルの検証、変容都市での環境・防災設計（高い換気性能とレジリエンスの確保）の提案、実市街地を対象とした換気性能評価を行う。
（協力機関：東京工業大学、横浜国立大学、近畿大学、お茶の水女子大学）(3) サブ課題 C「新興ウイルス感染症にロバストで健康・快適・サステイナブルなポストコロナ時代の室
内環境設計」
各種の空気感染性汚染物質の生成から室内空気中での輸送動態、さらに体内動態と免疫システム応答･
生理反応までを包括的に解析し、可視化する統合的システムを開発し、「富岳」に実装する。室内の快適性と感染症に対するロバストネスを両立させる多目的最適化による設計手法を提案する。
本年度は、統合シミュレーション技術の開発（気道内流体解析と粘膜上皮細胞上の粘液輸送解析の連成
解析手法の開発）、汚染物質濃度制御にロバストな室内環境設計法の開発（室内に形成される不均一混合場の定量的評価指標の開発）、実証問題の設定と解析モデル作成、統合シミュレーションの実施、多目的
最適化フレームワークの構築、ニューラルネットワークの「富岳」チューニングを行う。
（協力機関：九州大学、神戸大学）(4) サブ課題 D「マルチコンポーネント統合シミュレーションによるカーボンフリーガスタービン設計」
ガスタービンを構成する圧縮機、燃焼器、タービンを完全に統合して LES で一気に解く、世界初、真のマルチコンポーネント統合シミュレーション技術の開発を行う。本統合シミュレーションが、フラッシ
ュバックや燃焼振動等の異常燃焼、さらには、窒素酸化物（NOx）、一酸化炭素（CO）、およびすす等の様々なエミッションの正確な予測に有用であることを示す。本統合シミュレーション（サイバー空間）結果か
ら、異常燃焼発生時の物理メカニズムを解明し、さらに異常燃焼の兆候を自動で探知するための A I セ
ンシング技術を開発することによって、これらの知見、技術が実機（フィジカル空間）の設計や最適操作
条件選定に貢献しうることを実証する。
本年度は、圧縮機、燃焼器、タービン、各コンポーネントの解析実施（本計算）、統合シミュレーショ

-2-

ンの実施（簡易形状試計算）、2 次元数値シミュレーションによる水素、アンモニア乱流燃焼モデルの精度検証を行う。
（協力機関：京都大学、九州大学）(5) プロジェクトの総合的推進
以上のサブ課題 A～D の進捗を管理し、課題間の連携を促進するために、各課題の代表者を集めて統括する推進委員会を設置する。推進委員会は、３か月に一度程度の頻度で理研が主催して委員会を開催する。また、年に一度の成果報告会や、他の「産業競争力の強化」課題と連携した統合ワークショップを企
画・開催する。
２－２．実施内容（成果）以下の実施内容では、明言しない限り、全て「富岳」を用いて計算している。
(1) サブ課題 A「意匠空間を考慮した AI 支援多目的最適化による自動車空力デザイン」(a) 研究成果の要旨
AI 支援多目的最適化フレームワークの設計、ニューラルネットワークの「富岳」チューニング、基礎
流れを対象とした数億セル規模の縮約モデルと実車を対象とした数千万セル規模の縮約モデルの構築技術の開発、多目的最適化実用モデルの構築と解析の実施、及び意匠空間の数理モデル化と空力性能空間との融合技術開発を行った。この結果、「富岳」上で非定常三次元流れ場に対してニューラルネットワーク
によるモード分解のフレームワークの実装と性能評価を円柱周りの基礎流れでさらに進めることで、自動車の車体形状周辺の流れ場に対するモード分割を実施することができた。これにより、複雑形状周辺の乱流に対する本手法の適用可能性を確認することができた。また、最適化実用モデル構築のため、前年度定義した市販 SUV 車両群の形状特徴量を設計変数として、対象車両形状を変形させるモーフィング技術を確立することができた。さらに、車体外観がヒトに与える印象の特徴量を抽出し、車体形状の特徴量との相関を明らかにすることで意匠空間の数理モデル化を進めることができた。
In this fiscal year, we designed an AI-assisted multi-objective optimization framework, tuned the neural network
application on the supercomputer "Fugaku", developed a technology for constructing reduced models with hundredsof millions of cells for fundamental flow and with tens of millions of cells for real car aerodynamics, constructed andanalyzed a multi-objective optimization practical model, mathematically modeled the aesthetic design space, anddeveloped a fusion technology between the features of aesthetics design and aerodynamic performance. As a result,
implementation and performance evaluation of mode decomposition framework using neural network for unsteadythree-dimensional flow fields on "Fugaku" was further proceeded with the basic flow around a circular cylinder, andthe mode decomposition for flow fields around automobile body was realized. The results confirmed the applicabilityof the method to turbulent flow around complex geometries. In addition, for the purpose of constructing a practicalmulti-optimization problem, we established a morphing technique to deform the body shape of target vehicle using
the geometric features involved in a group of production SUVs as design variables. Furthermore, we extracted thefeature of impression that the exterior view of vehicle gave to humans. Clarifying the correlation with the geometricfeature of vehicle body, we proceeded the mathematical modeling of the design space.

-3-

(b) 本研究の内容
本年度(2022 年度令和 4 年)の研究成果を項目別に示す。
(b)-1 AI 支援多目的最適化フレームワークの設計
前年度に引き続き、意匠空間を考慮した自動車空力形状の多目的最適化問題への適用を想定して、Fig.
2.2.1-1 に示す AI 支援多目的最適化フレームワークを設計した。車体形状を表現する複数のパラメータを設計変数とし（パラメトリック最適化）、最適化手法にはチェビシェフ距離と選好領域と呼ばれる考え方を用いた多目的進化アルゴリズム CHEETAH/R (CHEbyshev-Epsilon optimizer AlgoritHm)を使用する。設
計変数に応じた車体形状は、基準車体形状を市販 CAE ソフトウェア（DEP MeshWorks）の形状モーフィ
ング機能によって変形させることで生成し、形状変更に伴う解析格子の再生成作業を必要としない複雑現象統一解析コード CUBE を用いた数値解析によって目的関数を評価する。そして、設計変数とする車
体形状パラメータを、デザイン的な意図と相関のある形状特性とすることで、意匠特性を考慮した形状最適化の実現を目指す。また、機械学習技術による最適解探索の効率化として、進化計算の過程で実施した
数値解析結果を基にガウス過程回帰モデルによるサロゲートモデルを構築し、目的関数の評価に要する計算コストの低減を図る。
本年度は、多目的進化アルゴリズムのパラメータチューニング、目的関数のハイブリッド評価において数値解析の適用を判定する指標を選定した。そして、既往の多目的最適化課題（「富岳」一般課題 hp210106
において課題設定）を対象に、サロゲートモデルを併用した目的関数のハイブリッド評価による多目的進化計算を実施した。同課題は、車体のプロポーションを構成する部位寸法を変化させる 29 の設計変数によりパラメトリックに車体形状を変形させ、正対風時の空気抵抗、空気揚力、および横風 3°時の空気抵抗増加量についていずれも最小化を目的とする 3 目的の多目的最適化課題である。進化計算の世代サイズは 32 個体/世代とし、各世代中の個体の目的関数（各空力係数）の評価に CFD 解析を適用する割合を50％とした。CFD 解析を適用する個体選定の指標はハイパーボリュームの改善期待値（EHVI）の上位個体とした。目的関数のハイブリッド評価に用いたサロゲートモデルの予測精度評価と、第 3 世代までの計算で得られたパレート解の散布図を Fig. 2.2.1-2 に示す。また、昨年度分析した市販車両の形状特性の主成分を基に、車体形状をパラメトリックに変形させる形状変形技術を構築した（項目(b-4)で後述）
。

START

Multi-objective evolutionary algorithm“CHEETAH/R”
Generate parameter sets
for new generation
Parametric geometry morphing(commercial software “DEP MeshWorks”)

Parameters(text)

Parameters(text)

Shift to the hybrid evaluation of
objective function after some generations.

Morphing geometry
based on the parameters
Geometry (STL file)

Aerodynamics simulation byComplex Unified Simulation Framework“CUBE”

Evaluate objectivefunction by asurrogate modelYes

Is CFD simulationnecessary?

Evaluate objective function
by CFD simulation
Objective functions(text)

No

Reached tothe specified generation?

Yes
END

Fig. 2.2.1-1 Framework of AI supported multi objective optimization.

-4-

No

8.0%

10.0%
CFD (CUBE)

6.0%

Surrogate model

CFD (CUBE)

Surrogate model

8.0%

2.0%

ΔCD

CD'

4.0%

0.0%

6.0%
4.0%

-2.0%
2.0%

-4.0%

-6.0%

0.0%
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
Car No.

1

3
2

5
4

7
6

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
Car No.
9

8

11
10

13
12

15
14

16

Fig. 2.2.1-2 Multi objective optimization of vehicle aerodynamics using the developed optimization framework.
(Top: Accuracy evaluation of the surrogate model, bottom: Scatter plots of the Pareto solutions.)(b)-2 ニューラルネットワークの「富岳」チューニング
「富岳」上で実装した分散深層学習の並列計算性能を示す。Fig. 2.2.1-3 は学習全体及びそれに含まれる各処理のウィークスケーリングを示したものである。順伝搬及び逆伝搬の主要処理 (図中のオレンジ及び赤のプロット) について良好なスケーリングを示すことが確認でき、また学習全体の性能 (青のプロット) についても、ネットワークパラメータ全体のリダクション通信を含むが全体の性能として十分にスケールすることが確認できた。Table 2.2.1-1 は学習全体の演算性能を示している。25,250 CPU (ノード) を使用した場合の学習全体の演算性能は 7.80PFLOPS を示すことが確認できた。

Fig. 2.2.1-3 Week scaling performance of distributed deep learning on "Fugaku"-5-

(dashed line shows ideal scaling).
Table 2.2.1-1 Computational performance of distributed deep learning on "Fugaku".
Number of CPUs

PetaFLOPS

750

0.31

5,250

1.81

10,2503.43

20,2506.28

25,2507.80(b)-3 縮約モデル構築技術の開発
これまで円柱周辺のレイノルズ数 1000 の 3 次元流れ場に対して適用してきたニューラルネットワーク
によるモード分割及び縮約モデル構築手法について、縮約された空間において流れ場を代表する変数(Latent vector)の生成部分を拡張し、これまで固定値として生成していたものを、正規分布に従う乱数として生成するための変更を行った(Fig.2.2.1-4)。また下記の通り、損失関数𝐿に対し、ネットワークの前半部(Encoder)を標準正規分布に近づけるための項𝐷𝐾𝐿 (KL ダイバージェンス)を追加した。
̃) =𝐿(𝑿, 𝑿

̃ 2∑𝑁
𝑖=1(𝑋𝑖 − 𝑋𝑖 )+ 𝐷𝐾𝐿 (𝑞𝜙 (𝒛|𝑿)‖𝑝𝜃 (𝒛))𝑛

1 𝐷
𝐷𝐾𝐿 (𝑞𝜙 (𝒛|𝑿)‖𝑝𝜃 (𝒛))＝ ∑ (1 + 𝑙𝑛 𝜎𝑗2 − 𝜇𝑗2 − 𝜎𝑗2 )2 𝑗=1
̃ はネットここで、𝑿は学習データ、つまり流れ場(流速)のデータをベクトル化したものであり、また、𝑿
ワークが出力したデータ、つまり学習データとなる流れ場データを一旦縮約した後に拡張して再現した
流れ場データである。また、𝒛は縮約された空間での変数をまとめたベクトル(Latent vector)であり、𝑞𝜙 及
び𝑝𝜃 はネットワークの前半部(Encoder)及び後半部(Decoder)に相当する。𝜇及び𝜎 2 は正規分布の平均及び分散を表す。
既存の手法及び上記の新手法を使用し、高精度シミュレーション結果を縮約された変数によって再現した結果を Fig.2.2.1-5 に示す。ここで、POD 及びニューラルネットワークによる手法(MD-CNN-AE、MDCNN-AE+VAE)はいずれも 64 個の変数に縮約した後に再現した流れ場である。高精度シミュレーションの結果(Full-order model)に比べて、POD によって再現した流れ場は、円柱後流の複雑な構造が十分に再現できていないが、一方でニューラルネットワークを使用した手法(MD-CNN-AE)ではより複雑な構造が見えている。更に、上記で説明した新手法(MD-CNN-AE+VAE)についてはより高精度シミュレーション結果に近い複雑な構造を捉えることができており、また鉛直方向(円柱のスパン方向)の流れ(図の最下段)も再現できていることが確認できる。

-6-

�

Decomposed flow fieldsnapshots� �
LatentvectorUnflatten&($! , %!" )$!

3D flow fieldsnapshots
!

�

%!"

Flatten

Input"

CombineCNN +Interpolate$"#

CNN +Pooling

Encoder

FCN%""&($ " , %"" )+Output3D flow fieldsnapshots!
"#�Unflatten
Decomposed flow
field snapshots�

正規分布に従う
乱数を⽣成

DecoderFig. 2.2.1-4 Application of variational auto-encoder on network for mode decomposition of flow field.

Full-order model

POD

MD-CNN-AE

MD-CNN-AE + VAE

Fig. 2.2.1-5 Comparison of flow fields after reduction and reconstruction by the conventional and present
methods.
更に、提案手法について、複雑形状周辺の流れ場への適用性を検討すべく、自動車の車体形状周辺の流れ場を対象とした実験を行った(Fig.2.2.1-6)。高精度シミュレーション実施時の設定を Table 2.2.1-2 に示す。
高精度シミュレーション結果に対し、提案手法を適用して 16 個の変数まで縮約し、その後拡張して元の流れ場の再現を試みた(Fig. 2.2.1-7)。元の流れ場はレイノルズ数 6.5×104 の乱流場であり、16 個の変数を使用して渦の詳細な構造を再現することはできなかったが、大域的な構造は捉えることができた。複雑な流れ場に対し 5 桁ほど計算量を下げた上で最低限の再現精度が担保されることを確認でき、提案手法によって精度と計算量のトレードオフの意味で有効な縮約モデルが構築できる可能性が示された。

-7-

Fig. 2.2.1-6 Numerical simulation of flow around automobile with high accuracy.
Table 2.2.1-2 Numerical conditions for the simulation.
Solver type

Incompressible

Domain size

20m×10m×5m

Number of cubes

1,800

Number of cells per cube

8×8×8 = 512

Total number of cells

1,800×512 = 921,600

Minimum cell size

19.5 mmTime step size

1.0×10-4 sec

Integration time

120 sec (1,200,000 steps)

Reynolds number

6.5×104

Time integration

Crank-Nicolson method

Pressure Poisson

Red-Black SOR

Viscous term

2nd order central difference

Convection term

QUICK

-8-

⾼精度シ ミ ュ レ ーショ ン 結果

縮約後に再現さ れた流れ場

Fig. 2.2.1-7 Reconstructed flow field around automobile from the 16 reduced variables.
(b)-4 多目的最適化実用モデルの構築と解析
最適化実用モデル構築のため、市販 SUV 車両群の形状特徴量を設計変数として、対象車両形状を変形させるモーフィング技術を確立した。このため、前年度実施した市販 SUV 車両群の形状特徴量を、車体形状のモーフィングに用いるパラメータとして再表現することとした。まず、多目的最適化のターゲット
車体に対して、既往のモーフィングパラメータ（「富岳」一般課題 hp210106 における空力多目的最適化
課題で設定）を調整し、各市販 SUV 車両のシルエットやプロポーションを再現して（Fig.2.2.1-8）、その際の各モーフィングパラメータ値を記録した。次に、得られた各市販車両に対応するモーフィングパラ
メータ値群に対して主成分分析を行い、その主成分に基づいて全モーフィングパラメータ値を変化させることとした。これは、その主成分得点をモーフィングによって車体形状を決定するための設計変数として使用することに相当する。結果として、第 1 主成分得点に相当する設計変数値を変化させることで得られる車体形状の変化と、を Fig.2.2.1-9 に示す。図に示されているように、少数の設計変数で車体の大局的な形状変化を与えるモーフィング技術を構築することができた。

Production
vehiclesTarget vehicle
Morphing

MorphedvehiclesFig. 2.2.1-8 Comparison of sideview of production cars (top) and the target vehicle morphed to have similar
proportions and silhouette (bottom)-9-Fig. 2.2.1-9 Parametric shape change based on the developed design variables 1 and 2 (left) and example of body
shape variation based on 10 design variables.
(b-4) 意匠空間の数理モデル化に向けた車体外観がヒトに与える印象の特徴量抽出
意匠空間の数理モデル化のため、車体外観がヒトに与える印象を SD 法（Semantic Differential method）
により評価した。SUV 車体が 360 度回転する動画を被験者に提示し、その動画から受ける印象を 14 の形
容詞対に対する 7 水準で評価回答するアンケート調査を行った。調査には Web 上のアンケートフォームを使用し、輸送工学系の大学生・大学院生および自動車産業界の会社員を対象に 100 名規模の調査とした。評価対象は Fig.2.2.1-10 上図に示す 65 車両とし、被験者 1 名あたりランダムに抽出した 13 車体を評価することとした。また、ある評価対象車両について、被験者に提示した回転動画から代表的フレームを抽出した画像を Fig.2.2.1-10 下図に示す。
主成分分析の結果、車両の外観から受ける印象の第 1、第 2 主成分として、それぞれ「実用的 ⇔ 情緒的」、「フォーマル ⇔ カジュアル」と解釈される主成分が抽出された。さらに、これらの印象の主成分と、前年度抽出した車体プロポーション特性の主成分との主成分得点の関係性を調査したところ、Fig.2.2.1-11 のように両者の間には有意な相関がみられた。このように、アンケート調査で得られた印象

Practical

Emotional

The 1st principal component scoreof impression from the exterior viewの主成分と車体形状の特徴量との相関を明らかにすることにより、意匠空間の数理モデル化を進めた。

The 1st principal component score ofgeometry of SUV vehicles on the market

Fig.2.2.1-10 Evaluated vehicle models(top) and

Fig.2.2.1-11 Correlation between the first principal

typical frames in the rotating movie (bottom) forcomponents of impression from the external view andthe impression survey.

geometry of SUV vehicles on the market.

(2) サブ課題 B「変容する都市・建築の自然擾乱対応の性能設計」(a)研究成果の要旨- 10 -

広域市街地解析に基づく流入条件データベースの大サンプルデータ作成、風・水・構造体の統一解法のためのモデル開発、風・熱・湿度の連成解法の基礎的なモデルの検証、変容都市での環境・防災設計（高い換気性能とレジリエンスの確保）の提案、実市街地を対象とした換気性能評価を行った。この結果、大サンプルデータに基づく評価モデルを提案し、風・水・熱との都市・建築連成問題を解明するモデルの妥当性を示し、変容都市の環境、防災のための具体的な建物配置を提示するとともに、換気性能を明らかにした。
In this fiscal year, this study carried out generation of large sampling data for set up of database of inflow condition
based on analytical results in broad region, development of the coupled simulation for unified solution of wind, water,and structure, validation of fundamental model for coupled solution of wind, heat, and humidity, proposal ofenvironmental and disaster prevention design (ensuring high ventilation performance and resilience) in a transformed
city, and the evaluation of ventilation performance in an actual urban area. As a result, we proposed an evaluation
model based on large sampling data, demonstrated the validity of the model to reveal the coupled problem of wind,
water, and heat for cities and buildings, presented specific building layouts for environment and disaster mitigationin transformed cities, and clarified ventilation performance for the transformation of a city.
(b)本研究の内容(b)-1 広域市街地解析に基づく流入条件データベースの大サンプルデータ作成(b)-1-1 NHM アンサンブルによる強風の極値の評価と極大ケースの解析とサンプルの選択
本研究では、気象擾乱と都市域の乱流構造による不確実性を考慮した実市街地における強風リスクの推定を行うため、東京都心で 42.5m/s の顕著な最大瞬間風速を観測した 2019 年令和 1 年台風 19 号を対象に、気象モデルによる都市レベルのローカルスケール解析、気象場に基づく乱流再生成、および市街地における強風に関する微細構造解析を連続的に継続することでそれぞれアンサンブル解析を実施し、極端台風下の都市を対象とする大サンプルデータの作成を実施した。
台風気象場の気象モデルによる作成は、「富岳」成果創出加速プログラム「防災・減災に資する新時代の大アンサンブル気象・大気環境予測」と協力して実施した。我々の課題では前述の課題において、作成された台風の進路や強度、環境場がそれぞれ異なる 1km 解像度 1000 メンバーの大アンサンブルデータで表現されたメソスケール気象場の特性に基づき、着目市街地での台風通過時の最大風速に基づいて最頻
値・最大値に相当するメンバー、台風経路・強度の再現性の高いメンバー、全体の確率分布を保持するランダム抽出などいくつかの防災上の観点において不確実性の評価が重要となる特徴を持つメンバーを抽出して、ローカルスケールのダウンスケール解析を実施した。
大サンプルのデータに対し、発生する強風の特徴として、平均風速、瞬間風速、最頻風速などを考慮し、最終的に全体で 50 程度のメンバーを選定してデータベース作成を行うことをめざす。ここでは、まず顕著な強風が発生する可能性の高い状況についての知見を蓄積するため、着目市街地での台風通過時の風
速が最も大きい 5 ケースを対象に解析を実施した。
Fig. 2.2.2-1 に気象庁非静力学モデル（NHM）による 1km 解像度 1000 メンバーのアンサンブルデータ
について、東京管区気象台・風速観測点地上 10m 高さにおける、10 分ごとの瞬時値から求めた最大風速のヒストグラムと超過確率分布を示す。台風通過時の最大風速は 18.2m/s から 8.4m/s の間で分布してお- 11 -

り、1m/s ごとに 10 程度の階級に区分し、なめらかに分布が表現されることを確認した。確率分布は、高風速側にテールが伸びる傾向が顕著であり、最頻値は 11m/s から 12m/s の間である一方、約 10%のメン
バーが 15m/s を、約 5%のメンバーが 16m/s を超過した。抽出した 5 ケースについては、風速は 18m/s 付
近の狭い範囲（17.8 m/s から 18.3m/s、相対誤差 3%以内）に収まっており、風速の確率分布の極値付近を見た場合についても一定の稠密性があることが確認された。

Fig. 2.2.2-1

Histogram and exceedance probability distribution of maximum wind speed at

Tokyo (z = 10 m) of the 1000-member NHM ensemble data at 1-km resolution.
Fig. 2.2.2-2 には、今回抽出した東京管区気象台における風速上位 5 ケースについて各地点での最大風速を示す。上位 5 ケースでは、伊豆半島の東岸から三浦半島の間の相模湾内に接線風速最大域があり、大半のケースで東京湾の浦賀水道で縮流が発生していることが確認された。台風の大きさや中心付近の対流構造により、高風速帯の幅が異なっている様子も確認された。
Fig. 2.2.2-3 に、抽出した風速上位 5 ケースと風速最頻値（11.5m/s）1 ケースの時刻歴を示す。上位 5ケースのピーク時刻については、台風進路や環境場の違いにより、ケース間で 1 時間以上の開きが発生していることが確認された。また中心付近の構造と経路の違いによりピーク付近で急激な立ち上がりを示すケース (メンバー398) も見られた。
風速最頻値付近のケースについては、台風の進路や環境場のより大きな不確実性を含んでおり、ピーク
時刻のずれが顕著に大きくなる傾向が見られた。今後の抽出したダウンスケーリングでは、分析対象の現象等によりターゲットを絞った評価の必要性が示唆された。

- 12 -

Fig. 2.2.2-2 Distribution maps of maximum wind speed for the top 5 members with the highest
wind speeds (1-km NHM)

Fig. 2.2.2-3 Time history of the top 5 members with the highest wind speeds and one member near
the wind speed mode.
気象モデル NHM によるローカルスケールのシミュレーションは、前述の解像度 1km でのアンサンブ
ルデータ(D1)に基づき、解像度 250m（D2）および解像度 50m (D3)でのダウンスケーリング計算を実行した。Fig. 2.2.2-4 に計算領域を示す。計算時間範囲は、D2 は 2019 年 10 月 12 日 9:00 UTC からの 6 時間、および D3 は Fig. 2.2.2-2 に示した D1 における各メンバーのピーク時刻の前後 50 分である計 100 分間である。Fig. 2.2.2-5 に気象モデル解析により得られた D2 における気圧と降雨領域、および D3 の高さ約315m と 30m における風速分布を示す。メンバー194・197・398 は、上空に接線風速が最も大きい領域が対象市街地である東京・大手町のほぼ上空にある時刻にピーク風速を生じているのに対し、メンバー336・
597 は接線風速最大域が通過後、より台風中心に近い位置で発生したものであることが確認できる。

- 13 -

Fig. 2.2.2-4 Calculation domains of meteorological model simulation (NHM)

Fig. 2.2.2-5 Simulation results of meteorological model NHM during the peak time. Pressure andprecipitation in D2 (left). Velocity field at z= 315 m and 30 m.

- 14 -

(b)-1-2 NHM アンサンブルに基づく気象モデル・工学 LES 融合解析
気象モデルによる解析結果を市街地 LES の流入条件として利用するためには気象モデルで減衰しやすい高周波変動成分を再生成させることが重要である。気象モデル・工学 LES 融合解析手法[2-1]に基づき、気象場に対する乱流再生成計算を実施した。東京管区気象台のある大手町地区を含む 12km×7km の領域を対象に、水平解像度 25m で計算し、メンバー197・597 については 20 分間、残りのメンバーについては
5 分間の流入風を取得した。
Fig. 2.2.2-6 に乱流再生成後のスパン方向鉛直断面、高さ 300m 水平断面における流れ場を示す。各ケ
ースで 1km ないし 1.5km 程度までの地表近傍において最も高いピーク風速が発生していること、その上では風速が小さくなったのちに上空で再び上昇に転じる傾向が確認された。メンバー194・197・398 は2km から 3km 程度の間隔で並ぶ、流れ方向に軸を持つロール渦構造とそこから生じた乱流構造が確認された。より台風の中心に近い位置の流れ場であるメンバー336・597 では、流速の変化がより大きいロール渦構造が再現され、風向の揺らぎも大きい傾向が確認される。既往文献において報告された台風中心の距離により異なる構造のロール渦[2-2] が発生する状況を再現している可能性があると考えられ、更なる分析が必要である。また鉛直断面を確認すると、後者では、地表近傍のピーク風速が発生する領域の高度が前者よりも低く、上空では風速が大きく低下し、シアが強い状態が発生している。以上より台風の通過時に起こるピーク時の異なる特徴を持つ乱流構造を再現できたことを確認した。

Fig. 2.2.2-6 Result of the turbulence regeneration calculation for meteorological field. Velocity inthe horizontal and spanwise section [unit: m/s].
Fig. 2.2.2-7 に、メンバー197・597 の融合解析前後での高さ 300m での流れ場および、領域の中央位置、高さ 361m での乱流エネルギースペクトルと時刻歴を示す。メンバー197 では、元の流れ場の構造に基づいてより微細な乱流構造が生成されており、スペクトルについても-5/3 乗則に従う風速変動が 0.1Hz 以上の領域まで拡大され、乱流としての適切な性質を持つ風速変動の再生成が確認された。なお、流出境界付
近では非物理的な振動の発達が確認された。一方、メンバー597 では、高周波数域での乱流エネルギーが若干増加したが、顕著な向上は見られなかった。気象場の風速変動のスケールに基づいてフィルタースケ
ールを更に調整することで改善を図ることが有効である可能性が示唆された。

- 15 -

Fig. 2.2.2-7 The velocity fields at horizontal section z = 300m and turbulent energy spectra before
and after the turbulence regeneration calculation in the cases of Members 197 and 597.
(b)-1-3 広域市街地解析に基づく大サンプルデータの作成
市街地で発生するピーク風速や風圧を解析するための気象擾乱を流入条件とする BCM-LES のアンサ
ンブル解析は、理研 R-CCS により開発された CUBE[3]を用いて広域解析(3 km × 6 km × 1.5km、解像度約1.5 m)の領域において、前項にて作成した流入風を東・南・西の 3 面から与え解析を行った。現時点では全領域の 5 分間のデータを 0.5s 間隔の高頻度で取得した。 Fig. 2.2.2-8 にメンバー197 と 597 の大手町地
区の瞬間の風速場（高さ z=100m 水平面および E-W 鉛直断面）を示す。風向はおおむねメンバー197 では
SE、メンバー597 では S であるが、細かな揺らぎを含んでいることが確認される。また鉛直断面からは、
実市街地の不均一な地表面粗度に基づく乱流と台風境界層内の風速変動の重層が確認される。

- 16 -

Member 197

Member 597

Member 197

Member 597

Fig. 2.2.2-8 Instantaneous velocity fields of Otemachi area (horizontal section at z=100m and EW cross section) obtained in the urban LES analysis
Fig. 2.2.2-9 に、計算開始後から 2 分間の 5 秒間隔データに基づいた E-W 鉛直断面における平均風速と最大瞬間風の分布を示す。Member597 は、 Member197 に比べて台風中心に近い位置での境界層構造であり、低高度に風速プロファイルのピークが存在する。気象場の風速変動および鉛直シアと建築物に基づく
地表近傍の乱流の相互作用がより顕著に発生し、上空の高風速な流体が地表近傍により頻繁に輸送される結果、地表近くの平均風速・最大瞬間風速がともに高くなる傾向が確認された。
Fig. 2.2.2-10 に東京管区気象台風速観測点の位置(高さ 35m) における時間平均操作を行っていない瞬時値と 3s 平均の風速時刻歴を示す。台風アンサンブル解析の流入風に基づき平均化時間以下の極めて微細な風速変動を有するデータが得られたことが確認された。

- 17 -

Fig. 2.2.2-9 Mean and maximum wind speed in E-W cross section calculated from the data for
120 s in 5-s intervals.

Fig. 2.2.2-10 Time histories of wind speed at the location of the observational point of Tokyo
Headquarters (z = 35 m).
(b)-2 実建物のファサード形状を再現した大規模 LES に基づく高層建築物の風圧力推定(b)-2-1 複雑な形状のファサード周りに形成される圧力場（赤坂）実建物のファサード形状が複雑化している中で、建築物形状の凹凸については風洞実験では 1/300～
1/600 程度の縮尺に照らし合わせた凹凸サイズが想定されているが、建物表面の部材における局所的なピーク風圧を予測する際には、風洞実験模型では再現できない微細な凹凸形状に基づく風圧力の影響については検討事例が少ない。
本研究ではBCMによる広範囲・高解像度LESにより、複雑形状を有する実建物のファサードの周りに発生する微細な渦構造を明らかにしたうえで高層建築物の風圧力への影響を明らかにした 2-4)。実験モデル- 18 -

の縮尺や解像度の制限から再現が困難なファサードの詳細形状周辺の乱流構造を明らかにし、建物の風圧に与える影響について考察する。対象領域ではFig.2.2.2-11に示される2 つの建物(a)、(b)を含む東京都心部の半径600mの範囲を再現している。建物(a)は、外壁から離れた位置に複数の柱や梁で構成された複雑なファサード形状を有している。(b)の建物では、下部の壁面に凹凸のある日射遮蔽用のルーバーが設置されている。いずれの建物においても20cm程度の部材の幅や凹凸があり、最小空間分解能は建物ファ
サードの近傍で0.25m、総メッシュ数は約5億である。なお、流入条件については筆者らが実施した東京都
心部9kmx9kmの広域解析（流入条件：乱流境界層）によってサンプリングされた流入変動風の時空間デー
タ2-5)を接続している。
Fig.2.2.2-11に建物(a)の隅角部付近の流速の大きさと渦度（z成分）の拡大図を示す。空間解像度25㎝の解析結果では空間解像度50cmの結果2-6)で再現されていなかったファサードを構成する柱の間を流れる流れが再現されている。また、Fig.2.2.2-12には建物(a)のファサードを構成する柱や梁の表面の風圧係数分
布を示している。Fig.2.2.2-12に示す風圧係数分布が柱や梁などを再現しない分布と異なる点として、建物(a)では柱や梁などの部材を再現することで建物隅角部付近のはく離の発生位置が変化しており、その周囲に存在する柱や梁などの部材の表面に作用する風圧力が変化していることが明瞭な差であると考えられる。Fig.2.2.2-12の分布において風上面の周縁部の柱の部材の風圧係数はこのはく離の影響を受け、ほかの風上面の風圧係数より小さく、0付近の値となっている。

Fig.2.2.2-11 Target model（Left）and comparison of spatial resolution (horizontal distribution of vorticityz)(Right)

Fig.2.2.2-12 Distribution of wind pressure coefficient on façade (building (a))- 19 -

Fig.2.2.2-13 には、鉛直のルーバーを有する建物(b)の渦構造と風圧係数の分布を示す。建物(b)ではルー
バーのない上部では直角の隅角部に対応して二次元性の高い渦構造が出ているが、ルーバーのある中層部から下では、隅切状の平面形状に付随して、より 3 次元性の高い細かな渦が発生している。また、風圧係数分布に着目すると、局所的に低圧の構造がルーバー表面に確認される。これは隅角部から剥離した流れが再付着時にルーバーの間に流れ込むことで副次的な細かな渦が形成されるためであると考えれられる。

Fig.2.2.2-13 Occurrence of extreme negative pressure generated on the uneven surface of the louver (building (b))(b)-2-2 広域市街地における乱流構造が風圧力に及ぼす影響（麻布台）
本研究では実建物のファサード形状再現モデルに基づく LES により、実市街地・実建物に存在する様々なスケールの乱流構造が風圧力に及ぼす影響を明らかにした 2-4)。実市街地に存在する乱流場の影響として、市街地に点在する複数の高層建築物から後流の乱流構造が発生する事例について解析を進める。解析
対象領域では都心部において個々の建築物から発生する渦構造の影響を受けるケースを想定し、風洞実験で再現が困難な広範囲な市街地領域（対象地から半径 1.2km 以内）の建物と地形を再現している。対象建物の高さは約 330m であり、150m を超える高層ビルが点在しており、流入面には NNW の風向を想定した乱流境界層（地表面粗度区分 III）を接続している。
まず、Fig.2.2.2-14 に断面 A、B の乱流強度（u 成分）に着目すると、断面 A では建物 F のせん断層の影響が約 800m 下流まで残存し、断面 B では前方建物の剥離せん断層が対象建物の風上面に接近し、正面での高い風圧が確認される。さらに Fig.2.2.2-14 に示すように、対象建物の隅角部において極大な負圧が高さ 250m～300m において出現しており、Fig.2.2.2-15 の時系列に示すように極大な負圧が同時に発生している。Fig.2.2.2-15 に示すように前方にある建物頂部から剥離したせん断層によって対象建物前方の流速が大きくなったことに加えて、隅丸形状の隅角部の前方から発生した剥離せん断層に付随して発生する鉛直の渦によって負圧が発生したと考えられる。一方で、高さ 200m の地点では、250m～300m の負圧と異なる時刻において極大な風圧が発生しているがこれは、風上建物頂部からのせん断層の影響が小さく、前方建物の後流の影響を受け、隅角部の周辺に断続的に高風速が出現するためであると考えられる。

- 20 -

Fig.2.2.2-14 解析対象建物の表面風圧分布と乱流強度

Fig.2.2.2-15 対象建物の隅角部における風圧係数の時系列データ
以上のように、実建物のファサード形状を再現した解析を実施することにより、建物隅角部付近における複雑ファサード形状の再現性やファサードの凹凸に伴って発生する副次的な渦によってファサードの複雑形状に由来する風圧力の分布が発生することを示した。ファサード実形状から発生する圧力分布の局所性・非定常性は外装部材の荷重評価を正確に進める上で把握する必要がある。特に剥離直後の荷重が大きくなる場所でその影響は大きくなることから、特性の系統的吟味は不可欠である。
(b)-3 風・水・構造体の統一解法のためのモデル開発(b)-3-1 マーカー粒子 IB 法により 1 次梁モードで強制振動する角柱の非定常風圧解析
近年、市街地の高層化はますます進行し、国内外で高さ 300m を大きく超える高層建築物が多数建設されている。このような建築物には地震力よりも風力が大きく作用することが知られている。また、空力不
安定性についても、アスペクト比が大きいため線形モードだけでなく曲げモードなどを含めてより精密に検討する必要がある。本研究では CUBE2-3), 2-7)のフレームワークを用いて、直交格子系での constraint IB
法を利用した数値流体解法により、強制振動する剛体角柱周りに生じる空気力の解析を実施した。ロッキングモードと 1 次梁モード振動による振動での負減衰力発生の比較を行った。併せて、振幅の変化に対する変化の検討と、解析時の格子解像度の結果への影響について確認を実施した。なお、CUBE では IB- 21 -

法を適用する領域をラグランジュマーカー粒子により定義する。この粒子の移動を制御することで、物体の変形を表現することができる。特に今回用いた移動境界問題用の CUBE では様々な物理モデルや変形モードに基づいて、ラグランジュマーカー粒子の位置と速度を容易に更新可能である 2-8),2-9)。
解析対象の角柱は、振動周波数を 3.4Hz に固定した 1:1:6 正方角柱である。標準ケースでは角柱頂部の振幅を角柱幅の 9.6%とした。解析の格子解像度は 1/50D とし、時間刻みは 1/25D のケースで 6.13×10-5 s(1/4800 周期)とした。流入風は空気の物性を想定した一様流である。時間積分には Crank Nicolson スキー
ム、空間離散化には QUICK スキームを用いた。
まず各無次元風速にともなうロッキングモードと梁 1 次曲げモードでの振動時の非定常空気力の正減
衰から負減衰力への遷移過程の比較を行った。正・負減衰力の判断は、20 周期以上の風外力の時刻歴データから、三角関数の直交性を利用し、変位に対する風外力の位相を算出した。強制振動による変位に対して、同じ方向の風外力の時刻歴の位相が先行する時、負減衰力が発生する可能性がある。Fig.2.2.2-16に、ロッキングモードおよび梁 1 次曲げモードでの振動時の圧力位相分布を示す。いずれの風速においても角柱低層部での位相が早く、特に前縁から 1/3 程度の位置で最も早くなった。風速の上昇に伴って角柱低層部を中心に位相が大きく前進し正となり、負減衰力を生じる領域が拡大することが確認された。
各時刻の渦構造を可視化した結果、Fig.2.2.2-17 のように無次元風速が小さい場合には、前縁頂部からの渦がより下方まで伸長してより長い時間側面に定在し、下層で先行する渦の剥離による圧力回復が遅れる様子が見られた。
Fig.2.2.2-18 に、ロッキング振動時の、角柱振幅変化に対する圧力位相分布の推移を示す。既往文献 210)と同様に、振幅が小さくなるほど位相の前進が確認された。また、本解析では、低風速時に中間部で見られたピークの位置が下がる傾向が見られた。各時刻の渦構造を可視化した結果、振幅が大きい場合には無次元風速が小さい時ときと同様に、渦の下方への伸長と下層での剥離の遅延が見られた。Fig.2.2.2-19に、ロッキングモードと梁 1 次曲げモードでの角柱全体でのモーメントに基づく正減衰・負減衰の評価を示す。同じ頂部振幅で比較した場合、梁 1 次曲げモードの方が一貫して位相が進んでおりより低い風速で負減衰が発生すること、より大振幅まで成長する可能性があることが確認された。Fig.2.2.2-20 はい
くつかの無次元風速における各振動モードでの位相分布の比較を示す。梁 1 次曲げモード振動時の方が常に全体的な位相が前進していることに加え、負減衰力が大きく働く領域がより高い位置に現れる特徴が見られた。Fig.2.2.2-21 のモード模式図が示すように、梁 1 次曲げモードでは、ロッキングモードよりも低層部での変位が小さくなるが、Fig.2.2.2-17 に示した小振幅で位相が前進する結果と符合する結果となった。
Fig.2.2.2-22 に振動による変位 0、速度が負方向に最大となる瞬間の 20 周期分のアンサンブル平均場
による角柱振動方向面の圧力等値面による渦構造の比較を示す。ロッキングモードの方が前縁の不安定により繰り返し発生する渦が側面の底部の方まで直線的に伸びる渦構造である一方、梁 1 次曲げモード
は底部までは伸びず、中層からやや後方に持ち上げられるように伸びている。その結果、梁 1 次曲げモードでは角柱の下層部で、負圧がロッキングモードより早いタイミングで角柱の後方に到達し、位相が前進していると考えられる。また異なる時間の可視化から、下層での位相の前進は下方で渦の剥離のプロセスが先行して起こることが要因であることが推察されたが、下方まで渦が伸長した場合にはより長い時間
側面付近に定在し剥離が遅れる傾向が確認された。また後流渦の形成の位置やタイミングにも違いが見られ、梁 1 次曲げモードの方が地表近傍の負圧が若干強い傾向が確認された。

- 22 -

梁一次曲げモード

Fig.2.2.2-16

ロッキングモード

Vibration at first-order beam bending and locking mode and lateral pressure
phase distribution at non-dimensional wind speed

Fig.2.2.2-17 Surface contours and 3D isosurfaces of instantaneous pressure during pressure
recovery (first-order bending mode of beam)

Fig.2.2.2-18

Pressure phase distribution on lateral surface of square cylinder at each

amplitude (left) and surface contours of instantaneous pressure during pressure recovery (firstorder bending mode of beam)（Right）

Fig.2.2.2-19

Evaluation of positive and negative damping basedon moments of the entire square cylinder- 23 -

Fig.2.2.2-20 Profile ofdisplacement for each

mode

Fig.2.2.2-21 Comparison of lateral pressure phase distribution during rocking mode
vibration and first-order mode vibration of beam

Fig.2.2.2-22

Comparison of vortex structures by pressure isosurfaces(b)-3-2 渋谷実建物を対象としたガラス破壊による内外圧の連成解析
強風時、局所的な風圧力や飛散物等によりガラスが破壊された場合、室内の圧力が急速に上昇すること
によって、他のガラスや外装材において連鎖的な破壊が進行することが予想される。本研究では、ガラスの破壊が発生した状況を想定した際に他の部材に作用する風圧力の評価を実施することを目的として渋谷地区の実建物を対象とした内外圧連成解析を実施した。
本研究では室内流れの流出入を発生させるため、流出入方向それぞれに開口部を設けた。当初は工場で製作される１枚当たり最大サイズのガラス板を想定した一辺 3m 程度の開口部１つを流入開口として中央に設置し、流出部については住宅の気密性を表す指標から算出した面積の開口部を等間隔に分散させて配置したが、解析上の問題として、 流量バランスの再現性の問題から内部に不自然な分布が生じる問題が確認された。そこで本年度は、モデル前面は全て開口とし、流出開口についても気密性に関わらず、十分大きい寸法を確保したモデルにおける内部流れの再現性を確認することとした。その結果、Fig.2.2.223、Fig.2.2.2-24 に示すように、得られた流速、圧力分布においては基準階内部におけるコア部や柱周りの流れ場の構造が再現された。一方で、モデル躯体からの流れ場の漏れ出しの発生、閉塞した空間の後方
開口部における空気の流出に対する流量バランスの再現性に対する課題が抽出された。

- 24 -

Fig.2.2.2-23 Streamline in indoor space of actual building

Fig.2.2.2-24 Flow and pressure distribution in indoor space of actual building(b)-3-3

LES による強風の影響を受けた降雨の解析

降雨に対する建物の設計では、洪水に対する設計と建物や敷地内の降雨に対する雨水排水設計が行われる。前者では浸水深に応じて敷地のかさ上げ、止水壁や防水扉の設置等の対策を行う。後者では建物や敷地内の雨水の流れを考慮して排水設計を行うが、強風時は雨滴が風で流されるため、壁面・屋根面や地表面への降雨の状況が無風時とは異なる。このため(社)空気調和衛生工学会給排水衛生設備規準・同解説
2-11)では、屋根、庇、ドライエリア等の設計において、壁面積の 1/2 を屋根面積に加算して雨水排水設計を行うよう求めている。また、中高層建物が林立する市街地では、強風による降雨の偏在が下水管路の排水に影響を及ぼすとの報告もある 2-12)。このように強風による降雨の挙動を予測し、雨水排水設計に反映することは重要であるが、その性状を精緻に予測して設計に反映している事例は見られない。
そこで本論文では、強風による降雨の挙動を LES で予測し、Fig 2.2.2-25 の実測結果 2-13)と比較検証を実施した。

- 25 -

(a)全体設置図（１辺 2m の立方体（地面から(b)CUBE２の風上壁面の雨量計

0.17m 浮いている）を 9 個配置）
Fig. 2.2.2-25 Condition of measurement2-13)
計算には、商用版 OpenFOAM の Helyx v3.4 の decoupledEularianSpray のモジュールを用いた。空気は非圧縮粘性流体の支配方程式を PIMPLE 法により解く。雨の挙動に関し Lagrange 法と Euler 法の解法が知られているが、本研究では Euler 法を用い、雨滴の直径毎に、計算セル毎の水の体積分率に関する質量保存則と運動量保存則を時々刻々解くことで雨滴の直径毎の雨滴ベクトルを得る 2-14)。雨滴に作用する空気力の反力に関する空気の NS 式の体積力はここでは無視している。
Fig 2.2.2-26 に計算格子と格子解像度を示す。流入境界には、地表面粗度区分Ⅲ相当の統計量を持つ流入変動風を与えた。実測結果に基づき、立方体頂部での平均風速は 1.7m/s (0:00～0:40)、1m/s (0:40～7:10)とした。降雨量は直径 0.3mm、 1mm、 2mm、 3mm、 4mm、 5mm の雨滴それぞれについて雨ベクトルを計算し、体積分率に掛合せて壁面への付着雨量を算出し、降雨強度および雨滴直径の頻度分布 2-15)に基づいて積算して求めた。

Fig 2.2.2-26 Grid division and spatial resolution.
Fig 2.2.2-27 に平均風速コンターと平均雨滴ベクトルを、Fig 2.2.2-3(c)、(d)に立方体 2、8 の風上壁面
第 1 セルの雨滴ベクトル ux 成分の高さ方向分布を示す。Fig 2.2.2-3 の雨滴ベクトルを見ると、直径 0.3mm
では風の影響が強く、風上壁面の上部で雨ベクトル ux が急に大きくなる。一方、直径 5mm では、上段から下段に向けて線形に減少している。

- 26 -

Fig 2.2.2-27

Average wind speed contour and average raindrop vector

Fig 2.2.2-28 に立方体 2、8 の風上壁面の catch ratio（壁面の降水量／地面の降水量）を実測と計算で比較して示す。上段で地表面の 20%程度、下段で 13%程度であり、給排水衛生設備規準の 50%より小さい。
また、上段から下段に向けて線形に減少する傾向は、Fig 2.2.2-27 (d)と同様であり、直径が大きな雨滴の影響が支配的であると考えらえる。総じて LES と実測の Catch ratio は概ね対応しており、LES の解析手法の妥当性が確認された。

Fig 2.2.2-28 Comparison of catch ratio between field measurement and LES- 27 -

(b)-3-4 CUBE に基づく渋谷地下街の浸水解析
本年度は風・水・構造体の統一解法の豪雨時の被災推定への適用に向け、階層型直交格子 CUBE に基づく渋谷地下街を対象とした浸水解析を実施した。
CUBE 数は 5,145、セル数は 21,073,920、最小セルサイズは 0.19513m である。境界条件は全面滑り壁境
界で、剛体部分は IB 粒子により表現している。Δt は 0.001 秒で 20000 ステップの計算（実時間 20 秒）である。計算はスーパーコンピュータ「富岳」により実行した。各出入り口での流入条件は、雨水が出入口の止水板を越えてくると想定し、全幅堰の式を使用し流量を推定した。流入水の水深はハザードマップにより決定し、これと出入口の幅から断面積を求め、流入速度を設定した。Fig.2.2.2-29 に解析モデル、地下鉄出入口の境界条件を示す。
本研究では、止水と避難誘導を連動させた避難計画や迅速な復旧に向けた排水計画の策定に向けた浸
水経路、浸水深・流速、時間経過の分析を実施した。Fig.2.2.2-30 に 20s における地下街における浸水状
況（水面）を示す。主に 20 秒後には B7、C、B6 各出入り口で流入が始まってから、人の集まる改札口や駅構内の吹抜部分に到達し、Fig.2.2.2-30 に示す範囲で浸水が生じることが観察された。

Fig.2.2.2-29

Analytical model and boundary conditions at subway entrancesFig.2.2.2-30

Isosurface of water based on flooding analysis(b)-4 都市の暑熱環境連成解法の基礎的なモデルの検証(b)-4-1 東京・北の丸モデルを対象とした解析と観測値との比較- 28 -

夏期の猛暑時における都市の暑熱環境連成解法の基礎的なモデルの適用性を精査することを目的とし、2018 年平成 30 年 7 月 23 日に記録的猛暑が観測された東京・北の丸エリア（Fig. 2.2.2-1）を対象として
LES を実施した。同エリアには東京管区気象台（風向風速計、温度計）および一般環境大気測定局（千代田区神田司町）が含まれている。比較的風向が安定し(NNW)且つ気温が上昇過程にある 11：30～12：00
を解析時間帯とした。流入条件は昨年度実施した WRF-LES における風速・温位のデータに基づいて中立を仮定した河合・田村のフィルタリング・リスケーリング手法 2-16）の解析と非等温 LES に基づくドライ
バ解析を実施 2-17)し、流入風を接続している。
東京・北の丸モデルの詳細を Fig. 2.2.2-31 に示す。都市形状の再現範囲は約 2km×1.4km、解析領域は
幅 2.6km×奥行 2.4km×高さ 1.2km を含む広域モデルとし、汎用数値地図「AW3D」に基づく建物・街区形状を反映した大規模計算を実施する。計算格子はテトラ要素を主とする非構造格子系であり、総要素数は約 6000 万である。ここでは、壁面温度境界条件を求めるために、対象時刻より前の 24.5 時間（7/22 11:00～7/23 11:30）について太陽位置の移動を考慮した熱放射・熱伝導解析を実施し、得られた表面温度（同図右下）を初期条件とし、乱流場・放射場の連成解析 2-18)を実施した。

Fig. 2.2.2-31 Detail of target model for Kitanomaru area.
北の丸モデルの解析において、東京管区気象台観測点における風速・温度をサンプリングし、対象日時における観測値との比較を行う。解析結果に基づく平均流れ場において、風速計および温度計の各サイト（図中★）を通過する流線図を Fig. 2.2.2-32 に示す。本計算風向（NNW）下においては、風速計（科学技術館屋上）に関しては流入直後から数百 m に渡り中高層街区が存在しており（当該建物も含め）それらの影響を受けていることが推察される。温度計（北の丸公園露場）は周囲に建物群は存在しないものの、実際の観測サイト周囲に存在する樹木群が解析モデルには反映されていないため、その影響については別途検討を要する。

- 29 -

Fig. 2.2.2-32 Streamlines through observation points for wind speed and temperature
風速および気温について、計算値と観測値を比較した結果を Fig.2.2.2-33 および Fig.2.2.2-34 に示す。
解析値は観測値の統計処理方法と可能な範囲で合わせた「10 秒ごとの平均値」とした。風速変動については計算値が若干過大となっているが、これは本解析では北の丸公園の植栽を再現していない影響だと
平均風速[m/s] 平均風速[m/s]考えられる。温度変動については、計算値・実測値ともに上昇トレンドは確認されるものの、絶対値が解1086析の方が若干過小評価となった。要因として、WRF-LESの地表面近傍温度が若干過小評価されているこ
実測10分値
実測(10分平均値)

4

実測10秒値

前10秒間の平均風速（風程/10秒間)
2 (風速:10秒平均値、気温:10秒瞬時値)とや、市街地スケールでの計算時の温度境界条件の影響だと考えられる。
10

8

6
10
4
8
2
6
0
4

6

4
2

0

解析10秒平均値
FFR(10秒平均値,スカラー風速)
30
35
40
45
50
55
実測10分値
実測(10分平均値)
11時からの経過時間[分]60

実測10秒値
前10秒間の平均風速（風程/10秒間)(風速:10秒平均値、気温:10秒瞬時値)

解析10秒平均値
FFR(10秒平均値,スカラー風速)
30
35
40
45
50
55
11時からの経過時間[分]1060

平均風速[m/s]

平均風速[m/s]
平均風速[m/s]

10

0
8

8
6

実測10分値
実測(10分平均値)

4

実測10秒値
前10秒間の平均風速（風程/10秒

2

0(風速:10秒平均値、気温:10秒瞬時値)

解析10秒平均値
FFR(10秒平均値,スカラー風速)
30
35
40
45
50
55

11時からの経過時間[分]
35
40
45
50
55 12:00
602
10402.2.2-3311時からの経過時間[分]Fig.
Comparison of simulation and observed values for wind velocity0
実測(10分平均値)
838
11: 303540455055 12:00
60
36
6
11時からの経過時間[分]104034
実測(10分平均値)
8
4
38
32
6
実測10分値
実測(10分平均値)
36
30
2
4
実測10秒値
前10秒間の平均風速（風程/10秒
34
11:30
2 (風速:10秒平均値、気温:10秒瞬時値)
35
40
45
50
55 12:00
60
00
32
解析10秒平均値
FFR(10秒平均値,スカラー風速)
11時からの経過時間[分]5030
11: 30
35
40 35
45 40
50 45
55 12:00
60 55
30
実測(10分ごと_前1分平均値)
11時からの経過時間[分]
11:30
35
40
45
50
55 12:00
60 10

平均風速[m/s]

気温[℃] 気温[℃]

平均風速[m/s]

11: 3040

11時からの経過時間[分]気温[℃]
平均風速[m/s]

11時からの経過時間[分]Fig. 2.2.2-34 Comparisonof simulation and observed
temperature.
388 values for air 実測(10分平均値)
実測(10分ごと_前1分平均値)- 30 -

温[℃]366(b)-5 変容都市における高い換気性能とレジリエンスを確保した市街地の提案と換気性能評価344(b)-5-1 最大瞬間風速に基づく耐風性能の評価
322
渋谷駅周囲における地域特性を考慮した変容エリアとして、
300 同駅南側に位置する、旧来からの桜並木の
11:
12:00
35
40
45
50
55 12:00
60
名所「さくら坂」と、広域再開発が進む「桜丘エリア」に挟まれた街区を対象として、中低層建物群によ
11:30
30
35
40
45
50
55
60

40
38
36
34

11時からの経過時間[分]11時からの経過時間[分]
実測(10分ごと_前1分平均値)
実測(10分平均値)

る構成も含め提案された変容パターンを示す（Fig. 2.2.2-35）とともに、各計画概要を以下に述べる。
【C-1：高層集約モデル】対象エリア内の既存建物群を同体積の高層（角柱）建物に置換する。建物寸法
C-0
は、平面 38m×38m、高さ 150m であり、渋谷エリアの超高層建築物と概ね同程度の規模である。周囲はC-1

大規模なオープンスペースが設けられる。
【C-2、3：中低層建物群（コの字型プラン）】中庭空間をコの字平面の建物で囲むプランを基本とし、建物の高低差を設ける。相対的に高い建物について、主風向（南北）方向に直交（壁状に）配置する C2、平行に配置する C-3 を計画した。

C-0

C-2

C-1

C-3

中低層建物群モデル
C-2
C-3
Fig.
2.2.2-35 Various transformation
patterns for Sakurazaka area

本年度は大型台風等による強風・暴風に対する都市の防災（耐風）性能を評価するため、対象エリア内
における最大瞬間風速分布について検討した 2-19)～2-21)。Fig. 2.2.2-36 に各ケースにおける最大瞬間風速（３秒ガスト、基準風速：高さ 220m で 35.9m/s）の空間分布を、Fig. 2.2.2-37 に風速レベルに対する面
中低層建物群モデル
積比率および超過確率を示す。本４ケース中、高層建物に集約された
C1 において強い最大瞬間風速分布が顕在化しており、エリア内のほぼ全域において 20m/s 以上（平均値は現況 C0 の約２倍）となっている。
C1 は一般的な再開発形態と位置付けられるものの、瞬間的な強風が広範囲に発生することが懸念されることから、都市防災（耐風）性能としては適切でない。一方、中低層建物群で構成されるケースのうち、C3 は他ケースに比べ最大瞬間風速が抑制されている。

Fig. 2.2.2-36 Spatial distribution for peak gust speed- 31 -

Fig. 2.2.2-37 Area ratio and cumulative percentage for peak gust speed(b)-5-2 平均風速および Breathability に基づく換気性能の評価
今後の都市の環境性能として、大気中の汚染物質や空気感染性ウイルス等に対する換気性能が要請される。ここでは、平均風速分布および運動エネルギー輸送量に基づいて、対象エリアの換気性能評価を行った 2-20)、2-21)。Fig. 2.2.2-38 に各ケースにおける平均風速（基準風速：高さ 220m で 5.6m/s）の空間分布
を、Fig. 2.2.2-39 に風速レベルに対する面積比率および累積頻度を示す。現状ケース C0 はエリア内に広く低速域がみられ、空気が滞留しやすいものと考えられる。C0 において累積頻度 50%に相当する風速（0.8m/s）を基準とすると、中低層建物群によるケース C2、 C3 は換気性能が大幅に改善されており、当該風速の頻度が 20%まで低減されている。特に C3 は比較的増速域の一様性が高いため、高い換気性能を見込むことができる。高層ケース C1 は広範囲に渡り高い風速レベルとなっているが、前頁の結果より防災（耐風）性能が不十分でないことから変容都市形態としては必ずしも適切でない。

Fig. 2.2.2-38 Spatial distribution for mean wind velocity

Fig. 2.2.2-39 Area ratio and cumulative percentage for mean wind speed
次に、鉛直方向の運動エネルギー輸送量に基づく Breathability2-22)～2-24）により、対象エリアの換気性能評価を行う。鉛直方向の運動エネルギー輸送量 Evt および各コンポーネントを Fig. 2.2.2-40 に示す。
上記の各コンポーネントは、平均運動エネルギー および 乱れの運動エネルギーに関する輸送方程式
において、①②は移流項、③④は拡散項（乱流応力成分）について鉛直成分を抽出したものに相当する。
なお、エネルギー輸送量評価面（水平断面）において算出される各コンポーネントおよび合算値 Evt に関し、正値は上向きの輸送すなわち汚染空気の排出性能、負値は下向きの輸送すなわち上空新鮮空気の導入
性能と位置づけられる。
ここでは、各変容モデルケースごとに換気性能を定量評価するため、評価面高さ（45m、15m）において水平方向の検査範囲を定め、評価領域内の鉛直方向エネルギー輸送量を算出する。Fig. 2.2.2-40 に、上記各コンポーネント別および合算値についてケース間の比較結果を示す。高さ 45m における新鮮空気導入量（負値）に関しては、C1 が顕著に大きく、他３ケースの差はわずかである。全ケースにおいて運動- 32 -

エネルギー輸送の寄与は①成分が主であり、続いて C1・C2 では②、C0・C3 では③であった。同高さにおける汚染空気の排出性能（正値）については、ケース間の差は小さく、①②の寄与が主であった。高さ15m においては、正値・負値ともに現状ケース C0 に比べケース C1～C3 で運動エネルギー輸送量が顕著に改善されている。正味の運動エネルギー輸送量については、ケース C1 を除き、高さ 45m では上向き、
高さ 15m では下向きであることが確認された。

①

②

③

④

コンポーネント

①

②

③

④

コンポーネント

Fig. 2.2.2-40 Case comparison for each component related to vertical kinetic energy transport
Fig. 2.2.2-41 は、breathability について絶対値が大きい領域（青：+1.6m3/s3、赤：-1.6m3/s3）を等値面で描画したものである。現状ケース C0 では対象エリア内に等値面が確認されない一方、高層ケース C1 では吹下しに伴う負領域が広く発生している。中低層ケースについては、C2 では主流方向に直交する中層
建物まわりの渦構造に起因した局所的な等値面が、C3 では風上側の超高層建物群の影響を受けた流れが縦溝状の空間に作用している様子が確認される。

Fig. 2.2.2-41 Isosurfaces for regions with large absolute values of breathability
都市の換気性能に関しては、評価手法・検査空間・物理量等について今後も引き続き検討を行うとともに、防災（耐風）性能、さらには温熱環境等も加味し、多角的視点から最適な都市形態を模索してゆく。
参考文献[2-1] Masaharu Kawaguchi, Tetsuro Tamura, Hidenori Kawai: Analysis of tornado and near-ground turbulence usinga hybrid meteorological model/engineering LES method. Int. J. Heat Fluid Fl., Vol. 80, No. 108464, 2019.
[2-2] Junshi Ito, Tsutao Oizumi, Hiroshi Niino: Near-surface coherent structures explored by large eddy simulationof entire tropical cyclones. Scientific reports, Vol 7.1, pp.1-10, 2017.
[2-3] Niclas Jansson, Rahul Bale, Keiji Onishi, Makoto Tsubokura: CUBE: A scalable framework for large-scale- 33 -

industrial simulation Int J High Perform Comput Appl, Vol 33, pp.678–98, 2018.
[2-4] 河合英徳，田村哲郎：実建物モデルを対象とする広範囲・高解像度 LES -壁面風圧変動を引き起こす乱流構造の抽出，日本流体力学会年会, 2022.9[2-5] 河合英徳，田村哲郎：東京都心部の流入変動風データベースを用いた超高層建築物の風圧力に関する LES, 計算工学講演会論文集，Vol.27 , 2022. 6[2-6] Kondo K., Kawai H., Tamura T. Nakajima K. Wind Pressure Characteristics of High-rise buildings in Middle
and High-height Urban Areas Spread over Local Terrain, WCCM-APCOM, 2022.
[2-7] Rahul Bale, Neelesh A. Patankar, Niclas Jansson, Keiji Onishi, Makoto Tsubokura: Stencil Penalty approachbased constraint immersed boundary method, Computers and Fluids, Vol 200, 104457, 2020.
[2-8] Koji Nishiguchi, Rahule Bale, Shigenobu Okazawa, Makoto Tsubokura: Full Eulerian deformable solid-fluid
interaction scheme based on building-cube method for large-scale parallel computing, Int J Numer Methods
Eng. Vol 117, pp.221–248, 2019.
[2-9] Tokimasa Shimada, Koji Nishiguchi, Rahule Bale, Shigenobu Okazawa, Makoto Tsubokura: Eulerian finite
volume formulation using Lagrangian marker particles for incompressible fluid–structure interaction problems,
Int J Numer Methods Eng. Vol 123, pp.1294–1328, 2022.
[2-10] 丸山勇祐: 高層建物周りの流れ場と空力不安定振動に関する研究, 学位論文, 2000.
[2-11] (社)空気調和衛生工学会，給排水衛生設備規準・同解説 SHASE-S-206, 2009.
[2-12] 永野雄一，伊藤一教，ビル風による地上雨量分布変化が浸水解析に及ぼす影響，土木学会論文集
B1（水工学）, Vol. 75, No.2, I_1273-I_1278, 2019.
[2-13] Kubilay, Aytaç, Numerical simulations and field experiments of wetting of building facades due to winddriven rain in urban areas, A thesis submitted to attain the degree of DOCTOR OF SCIENCES of ETH
ZURICHAYTAÇ, 2014.
[2-14] A. Kubilay, D. Derome b, J. Carmeliet, Analysis of time-resolved wind-driven rain on an array of low-risecubic buildings using large eddy simulation and an Eulerian multiphase model, Building and Environment,
Vol. 114, 2017.
(3) サブ課題 C「新興ウイルス感染症にロバストで健康・快適・サステイナブルなポストコロナ時代の室
内環境設計」(a)研究成果の要旨
令和 4 年度は、統合シミュレーション技術の開発（気道内流体解析と粘膜上皮細胞上の粘液輸送解析の連成解析手法の開発）、汚染物質濃度制御にロバストな室内環境設計法の開発（室内に形成される不均一混合場の定量的評価指標の開発）、実証問題の設定と解析モデル作成、統合シミュレーションの実施、
多目的最適化フレームワークの構築、ニューラルネットワークの「富岳」チューニングを行った。
具体的には、鼻腔･口腔から気管支までの上気道モデルと気管から気管支第 16 分岐(終末細気管支)までの下気道モデルを統合した数値気道モデルの界面境界モデルとして粘液輸送モデルを統合し、気道内気流解析に液膜から液滴形成を再現する Eulerian Wall Film モデルを実装し、気道内 Lagrange 粒子解析と連成解析することで、咳に伴う気道内での飛沫･飛沫核生成と室内環境中への呼出を連続して解析する技術を確立した。また、室内に形成されるスカラー不均一混合場を定量的に評価する指標として新たな換気効- 34 -

率指標(Net Escape Velocity, Net Escape probability)を考案し、その厳密な数学的定義を示すと共に、換気設計への実装可能性を簡易室内モデルを対象とした感度解析にて確認した。
In this fiscal year, we have developed (i) integrated simulation methods for airflow in lumen and mucus transport
on epithelial tissue surfaces, (ii) robust indoor environmental design method for controlling pollutant concentrations
using newly developed ventilation index for directly evaluating heterogeneity of mixing, (iii) multi-objective
optimization framework for "Fugaku" based on neural network.
Specifically, we have integrated a mucus transport model as an interface boundary model for a numerical airway
model that includes an upper airway model from the nasal and oral cavities to the bronchi and a lower airway modelfrom the trachea to the 16th bronchus (terminal bronchus). We also implemented an Eulerian wall film model that
reproduces droplet/droplet nuclei formations from a mucus liquid film for airflow analysis in the lumen and coupledwith Lagrangian particle tracking analysis in the airways. Finally, we have established a series of numerical techniquefor continuous analysis of droplet/droplet nuclei in the airways associated with coughing and expiration into theindoor environment.
In addition, a new ventilation efficiency index (Net Escape Velocity, Net Escape Probability) was developed as anindex to quantitatively evaluate the scalar heterogeneous mixing field formed in an indoor environment, and its
rigorous mathematical definition was presented.
(b)本研究の内容
本年度(2022 年度令和 4 年度)の研究成果を項目別に示す。
(b-1) 細気管支(16 分岐)を含む上気道から下気道までを再現した数値気道モデルの改良
2021 年度に作製した鼻腔・口腔から気管支第 16 分岐(各 lobe 毎に一流路を再現)までを含む上気道･下気道統合型の in silico 肺モデルを改良し、(ボランティア被験者の CT データより忠実に抽出した気管支第8 分岐までの生体モデルを基に)気管支第 8 分岐以降の全ての流路に対して分岐アルゴリズムを適用して気管支 16 分岐まで拡張した改良型の in silico 肺モデルを作製した。概要を Fig. 2.2.3-1 に示す。

- 35 -

Fig. 2.2.3-1 (A) in silico airway from nostrils to 16th generation. (B) The computational interested domain includes
the human respiratory tract, facial features, and external breathing zone.

- 36 -

Fig. 2.2.3-2 Comparison of deposition efficiency in the upper airway (A) and the lower airway (B) and (C) withprevious and experimental data at the three investigated breathing flow rates
上気道から下気道気管支第 16 分岐まで再現した改良型 in silico 肺モデルを対象として定常呼吸を仮定- 37 -

した条件で気道内流れ場解析を実施した上で、鼻孔経由で流入した粒子の気道内粒子沈着確率を解析した結果を Fig. 2.2.3-2 に示す。In vivo ならびに in vitro 結果と十分に一致する結果が得られている。

Fig. 2.2.3-3 Comparison of local deposition distribution for the coarse particles (5, 10 m) at the low breathing
flow rate between the scintigraphy image (A) and CFPD analysis (B, C); at the high breathing flow rate between
the scintigraphy image (D) and CFPD analysis (E, F).

- 38 -

Fig. 2.2.3-3 は下気道で粒子沈着分布の解析結果の一例を in vivo (シンチグラフィーによる可視化画像)
結果を比較した図を示す。定量的な比較は容易では無いが、分岐位置での粒子沈着の傾向、右肺下葉への粒子流入の傾向などは in silico 肺モデルによる解析結果は in vivo 結果と良く一致しており、一定程度の予測精度が確保されていると判断できる。
(b-2) 数値気道モデル界面境界条件としての Eulerian Wall Film モデルの実装
COVID-19 に代表される経気道感染は、咳や発話に伴う呼吸器系での飛沫･飛沫核が汚染発生源であることから、この飛沫･飛沫核発生量の正確な予測は経気道感染予測の基幹である。気道内のウイルス増殖のホットスポットは不均一に分布することも良く知られた事実であることから、飛沫･飛沫核中のウイル
ス含有量は、その発生部位･メカニズムに強く依存する。気道内での飛沫･飛沫核発生の物理的メカニズムに関しては、(i) 気管支収縮に伴う薄膜分裂、(ii) 喉頭領域での粘液剥離、(iii) 口腔内振動、が主要な要因とされている。本研究では、数値気道モデルに Eulerian Wall Film モデルとラグランジュ分散相モデルを連成した数値解析を実施することで、特に咳による気道内での飛沫･飛沫核の発生位置、発生量の定量予測に取り組んだ。
計算負荷を軽減するために、(b-1)で作製した in silico 肺モデルの鼻腔部分と気管支第 3 分岐よりも下流部分を削除することで、口腔から気管支第 2 分岐までの形状を再現した数値気道モデルを作成した。咳
による飛沫の生成は口腔形状および口の開き具合によって変動する口腔内の気流に大きく依存していると考えられるため、このモデルでは歯列形状や口腔開口をより詳細に再現している。

A

GingivaHardpalateOropharynx

Mouth
opening
Teeth Buccal
mucosaLarynxTrachea

Bifurcations

Fig. 2.2.3-4 Numerical respiratory tract (Airway) model from oral cavity to bronchial tube
この数値気道モデルの粘膜界面条件として、Eulerian Wall Film モデルを実装し、咳を想定した気流プロ- 39 -

ファイルを与えた流れ場解析を実施することで、ラグランジュ分散相モデルを連成させた気道内の飛沫拡散解析を行った。粘膜フィルムの初期条件として、気道壁面には 0.33μm、口腔部分には厚さ 68.9μm のフィルム厚を設定した。
Table 2.2.3-1 Model Parameters for EWF model.
Film material

Water

Film thickness (m)3.3E-07 – 6.89E-05 (0.33 μm – 68.9 μm)

Critical shear stress-CSS (Pa)10

Diameter coefficient

0.14

Mass coefficient

0.5

Surface tension (N/m)

0.072

Maximum thickness control (mm)68.9Time step size (s)0.0001

咳に相当する速度プロファイルを気道モデル境界条件として与えた場合の初期 0.1 秒間での気道内での飛沫発生量の時間変化を Fig.2.2.3-5 に示す。粘液フィルムが剥離して粒子生成する部位には明確な特徴があり、口腔内でも特に歯列位置からの粒子生成が過半となる結果が得られている。

- 40 -

Fig. 2.2.3-5 Origin of stripped particles during coughing
Fig. 2.2.3-6 には咳を再現した風速プロファイルを与えた解析の時間進行に伴う特に下気道における剪
断応力分布解析結果を示す。また、下気道での粒子生成位置も併せて示す。

100

t = 0.03s

t = 0.05s

t = 0.08s

Fig. 2.2.3-6 Visualization of stripped particles- 41 -

Fig. 2.2.3-7 Validation with experimental data after optimizing the parameters for the EWF model.
Fig. 2.2.3-7 には EWF モデルのパラメータ調整を行った上で、気管支第二分岐から上流側の上気道モ
デルを対象として解析した気道内発生粒子の粒径分布を示す。Yang らの被験者を対象とした実測結果を精度良く再現する結果が得られており、気道モデル内に EWF モデルを実装することで、飛沫飛沫核生成量の予測が可能となることを確認した。
(b-3) 不完全混合場を定量評価するための新たな換気効率指標の開発と簡易室内モデルを対象とした換
気効率指標の感度解析
室内に形成される不均一濃度場の情報を基に算出した空間平均濃度は、質点系の質量保存則より導かれる完全混合濃度とは異なる値となる。そのため、質点系の質量保存則の完全混合濃度を、空間平均濃度に入れ替えた場合の換気量は、外気導入量とは異なる値となり、この流量は Purging Flow Rate と定義される。また PFR は室空間の全領域から、局所領域の平均濃度を決定する実質的な換気流量という概念へと展開されており、更に、この Local Purging Flow Rate の概念を基に、局所領域への汚染物質の再帰回数
Visitation Frequency や局所領域内の平均滞在時間という新たな換気効率指標も提案され、詳細な数値解析
結果も報告されている。
居住域や呼吸域といった室内の局所領域の平均濃度を環境制御の対象とする場合、対象となる局所領域のサイズは非常に大きな幅を有する。この局所領域を徐々に小さくした場合に、極限として、局所点の濃度を議論することになり、この点の濃度を決定する速度スケールが定義できる。本研究では、この速度スケールを Net Escape Velocity (NEV)と称することとし、NEV が有する物理的なメカニズムと明快な数学的定義を導出した。
室内に形成される不均一流れ場によって、室内に不均一濃度場(スカラー濃度場)が定常的に形成されて- 42 -

いる条件を前提とすれば、このスカラー濃度𝐶 [kg/m3]と濃度場について、スカラーソース𝑆𝜙 [kg/m3s] を内部に含む検査体積を𝑉、𝑉の作る閉曲面ベクトルを𝑺として、その表面におけるスカラーフラックス放出面を𝑆𝑜 、流入面を𝑆𝑖 、𝑺 = 𝑺𝑜 − 𝑺𝑖 と定義する。Net Escape Velocity (NEV)をベクトル表記すれば、次式となる。
𝐶𝑜 ∫ 𝒖𝑵𝑬𝑽 ⋅ 𝑑𝑺 = ∫𝑭 ⋅ 𝑑𝑺 = ∫ 𝑆𝜙 𝑑𝑉 = 𝑚̇ = 𝐶𝑜 𝑞𝑃𝐹𝑅
𝑺𝑜

𝑺

(2.2.3-1)

𝑉

𝐶𝑜 =

1
∫ 𝐶(𝒙)𝑑𝑉
𝑉 𝑉(2.2.3-2)(2.2.3-3)

𝑚̇ = ∫ 𝑆𝜙 𝑑𝑉
𝑉(2.2.3-4)

𝑭 = 𝐶𝒖 − 𝐷∇𝐶
ここで、𝑚̇ [kg/s]はスカラー発生量、𝐶𝑜[kg/m3]は平均濃度、𝒖𝑁𝐸𝑉[m/s]は Net Escape Velocity ベクトル、

𝑞𝑃𝑅𝐹 [m3/s] は Purging Flow Rate、𝑭 [kg/m2s]はスカラーフラックスベクトル、𝐷 [m2/s]は拡散係数、を示す。
NEV は閉曲面ベクトル𝑺で定義される正味のフラックスである。同様に、スカラーフラックス放出面である𝑆𝑜 のみで定義される NEV*をベクトル表記すれば次式となる。
̇𝐶𝑜 ∫ 𝒖∗𝑵𝑬𝑽 ⋅ 𝑑𝑺 = ∫ 𝑭 ⋅ 𝑑𝑺 = ∫ 𝑆𝜙 𝑑𝑉 + ∫ 𝑭 ⋅ 𝑑𝑺 = 𝑚̇ + Δ𝑚

(2.2.3-5)

Δ𝑚̇ = ∫ 𝑭 ⋅ 𝑑𝑺

(2.2.3-6)

𝑺𝑜

𝑺𝒐

𝑉

𝑺𝒊

𝑺𝒊

𝒖∗𝑁𝐸𝑉 =

𝑭𝐶𝑜

(2.2.3-7)

また、検査体積に対するスカラーフラックスの全体の収支式は、(2.2.3-8)

𝑑𝑖𝑣𝑭 = ∇ ∙ 𝑭 = 𝑆𝜙であるから、これを検査体積で積分すれば次式となる。
∫ 𝑑𝑖𝑣𝑭𝑑𝑉 = ∫𝑭 ⋅ 𝑑𝑺 = ∫ 𝑭 ⋅ 𝑑𝑺 − ∫ 𝑭 ⋅ 𝑑𝑺 = ∫ 𝑆𝜙 𝑑𝑉 = 𝑚̇
𝑉

𝑺

𝑺𝒐

𝑺𝒊

(2.2.3-9)

𝑉

この検査体積 CV 内だけにスカラーソースがある場合は、検査体積に対する流入フラックス(Influx)はすべて CV に戻ってくる成分だと考えられることから、この再帰確率 (CV に戻ってくる確率)を𝛼とすると次式が成立する。
∫ 𝑭 ⋅ 𝑑𝑺 = 𝑚̇ ∑ 𝛼 𝑛 =
𝑺𝒊

𝛼
𝑚̇.
1−𝛼

(2.2.3-10)

一方で、全体のバランスを考慮することで再帰確率𝛼が定義できる。
(1 − 𝛼) ∫ 𝑭 ⋅ 𝑑𝑺 = 𝑚̇.

(2.2.3-11)𝑺𝑜

𝛼=∫𝑺 𝑭 ⋅ 𝑑𝑺𝑖∫𝑺 𝑭 ⋅ 𝑑𝑺𝑜=∫𝑺 𝑭 ⋅ 𝑑𝑺𝑖

𝑚̇ + ∫𝑺 𝑭 ⋅ 𝑑𝑺

(2.2.3-12)𝑖

また、検査体積 CV から流出するスカラーフラックスは、I 項の再帰無し(戻ってこない)成分と、II 項- 43 -

の再帰(戻ってくる)成分に分離することができる。結果として、直接排出確率 Net Escape Probability(NEP)が定義できる。
∫ 𝑭 ⋅ 𝑑𝑺 = 𝑚̇ + ∫ 𝑭 ⋅ 𝑑𝑺 = (1 − 𝛼) ∫ 𝑭 ⋅ 𝑑𝑺 + 𝛼 ∫ 𝑭 ⋅ 𝑑𝑺
𝑺𝑜
𝑺𝒊
𝑺𝑜
⏟
⏟ 𝑺𝒐
𝑰

𝑁𝐸𝑃 = 1 − 𝛼(2.2.3-13)

𝑰𝑰(2.2.3-14)

NEV は、局所フラックス𝑭 と、その点の脱出確率 NEP で定義されることから、一部開口を有する閉鎖空間である室内条件において、NEV が物理的には正味の脱出速度 (“Net” Escape Velocity) を意味すると解釈できる。
スカラーフラックス放出面である𝑆𝑜 のみで定義される NEV*ベクトルである𝒖∗𝑁𝐸𝑉 は、物理的にはスカラー放出面𝑺𝑜 上での濃度𝐶𝑜 を持つ空気塊(パーセル)の平均速度を示しており、NEV は、再帰分を取り除いた速度成分となる。
本研究では、自然換気される単純室内モデルを対象として、室内での NEV 分布を解析した。対象とした建物ブロック群と自然換気を想定した建物モデルを Fig. 2.2.3-8 に示す。本解析では主流方向に対して
流れ方向に開口部を設置した STR 条件と流れ方向に対して側面側に開口部を設置した LAT 条件の 2 ケースで解析を行った。
自然換気される室内の空気齢(Age of Air；AoA)分布解析結果を Fig. 2.2.3-9 に、2.2.3-7 式に定義されたNEV の室内分布の解析結果を移流ベクトル分布と共に共に Fig. 2.2.3-10 に示す。空気齢分布は名目換気
時間スケールで無次元化しているため、上流側の自然換気流入口近傍では AoA=0.0 となっており、排気口近傍では AoA=1.0 の値となっているが、室内には不均一分布の形成が確認できる。流入新鮮空気の室
内分配は均一ではなく、局所領域によって非常に大きなバラツキがあり、この不均一分布を考慮した環境制御の重要性が示唆される。

Fig. 2.2.3-8 Schematic diagrams of computational domain. Ventilation building models with (a) streamwise (STR)- 44 -

and (b) lateral (LAT) openings.

Fig. 2.2.3-9 Age of air of (a) STR and (b) LAT cases.
NEV 分布は、移流ベクトルの分布とは異なる分布を示していることから、汚染物質(スカラー)の実質的な排出方向と移流風速の向きは必ずしも一致せず、乱流拡散の影響が相対的に有意となる局所領域では、
移流風速のみによる濃度制御には限界があり、NEV 分布の確認が重要となることが理解できる。

NEV
Flow

Fig. 2.2.3-10 Net Escape Velocity (NEV) of (a) STR and (b) LAT cases.

- 45 -

NEV
Flow

参考文献
[3-1] Naoki Ikegaya, Mats Sandberg and Kazuhide Ito, Rigorous mathematical formulation of net escape velocity
and net escape probability determining a macroscopic concentration, Indoor Air. 2022;32 (7): e13072[3-2] Hana Salati, Mehrdad Khamooshi, Jingliang Dong, Kazuhide Ito, David Fletcher, Sara Vahaji, Kiao Inthavong.
Exhaled Aerosol and Jet Flow Characterization During Nasal Sneezing, Aerosol and Air Quality Research, 22(4),2022, 210338[3-3] Nguyen Dang Khoa, Nguyen Lu Phuong, Ken Takahashi and Kazuhide Ito. Transport and deposition of inhaled
man-made vitreous and asbestos fibers in realistic human respiratory tract models: An in-silico study, Japan
Architectural Review 2022, Vol.5, no.4, pp 592–608[3-4] Haruna Yamasawa, Teruaki Hirayama, Ryota Muta, Kazuki Kuga, Tomohiro Kobayashi, and Kazuhide Ito.
Influence of inlet turbulent condition on the formation mechanism of local scalar concentrations, Japan
Architectural Review, 2022, Vol.5, no.4, pp 691–701[3-5] Sung-Jun Yoo, and Kazuhide Ito, Validation, verification, and quality control of computational fluid dynamics
analysis for indoor environments using a computer-simulated person with respiratory tract, Japan Architectural
Review, 2022, Vol.5, no.4, pp 714–727.
[3-6] Ken Bryan Fernandez, Naoki Ikegaya, Qingyan Chen, and Kazuhide Ito, Determination of age of air, purging
flow rate, and net escape velocity within a cross ventilation model sheltered by urban-like block arrays usinglarge-eddy simulations, Building and Environment, 226 (2022), 109759.
(4) サブ課題 D「マルチコンポーネント統合シミュレーションによるカーボンフリーガスタービン設計」(a)研究成果の要旨
令和 4 年度は、圧縮機、燃焼器、タービン、各コンポーネントの解析実施（本計算）、統合シミュレー
ションの実施（簡易形状試計算）、2 次元数値シミュレーションによる水素、アンモニア乱流燃焼モデルの精度検証を行った。この結果、川崎重工業製ガスタービンを対象に、圧縮機、燃焼器、タービンの簡易形状を対象にした統合シミュレーションを実現した。また、本計算のため、各コンポーネントの形状およびメッシュの再確認・修正を行った。さらに、燃焼モデルとして利用するフレームレット法には非断熱性、選択拡散、火炎伸長の効果を考慮することが燃焼器内乱流燃焼場の高精度予測に有効であることがわかった。
In this fiscal year, a trial integrated simulation of a compressor/combustor/turbine assembly of a gas turbine was
performed, and the validity of flamelet model as a combustion model for hydrogen and ammonia combustion was
examined by two-dimensional (2D) numerical simulations. The trial numerical simulation based on large-eddy
simulations (LES) was successfully carried out, and the 2D numerical simulation showed that the consideration of
the effects of nonadiabaticity, differential diffusion and stretch in the flamelet model is useful for the high-fidelity
simulations.
(b)本研究の内容- 46 -

(b)-1 ガスタービンシミュレーション
本解析では、非構造格子に適用可能な FrontFlow /red large-eddy simulation (LES)を、京都大学と（株）数
値フローデザインにより改良したものを用いた[4-1,4-2]。
(b)-1-1 統合シミュ―ション（簡易形状試計算）(b)-1-1-1 計算対象
Fig. 2.2.4-1 に対象とした川崎重工業株式会社の産業用ガスタービンエンジンを示す。2021 年度令和 3
年度はガスタービンエンジンを構成する圧縮機、燃焼器、タービンをそれぞれ単独で解析した。2022 年
令和 4 年度は 2021 年令和 3 年度の成果に基づき、圧縮機、燃焼器、タービンを統合した解析を行った。
本解析で用いたメッシュの節点数は約 0.9 億、要素数は約 4.1 億である。
(b)-1-1-2 計算方法、モデル、条件
本解析では、予混合燃焼をベースとした Flamelet Generated Manifold (FGM)法[4-3]を用いた LES を実施する。FGM 法の支配方程式は、以下の連続の式、運動量の保存式、混合分率の保存式、反応進行変数の保存式、エンタルピーの保存式、および NOx の生成予測のための輸送方程式で構成される。

質量保存式：̅𝜕𝜌𝜕+ 𝜕𝑥 (𝜌̅ 𝑢̃𝑗 ) = 0
𝜕𝑡𝑗(2.2.4-1)

運動量保存方程式：
̅𝑢
̃𝑢
𝜕𝜌
𝜕𝜎̃
𝜕𝜏
̅𝑢
̃𝑖
𝜕𝜌
𝜕𝑝̅
𝑖 ̃𝑗
+
=−
+ 𝑖𝑗 + 𝑖𝑗
𝜕𝑡
𝜕𝑥𝑗
𝜕𝑥𝑗
𝜕𝑥𝑗
𝜕𝑥𝑗

(2.2.4-2)

混合分率輸送式：
𝜕𝑞𝑍𝑗
̃
̅ 𝑍̃𝑢
̃𝑗
𝜕𝜌
̅ 𝑍̃
𝜕𝜌
𝜕
̃𝑍 𝜕𝑍 ) +
+
=
(𝜌̅
𝐷
𝜕𝑡
𝜕𝑥𝑗
𝜕𝑥𝑗
𝜕𝑥𝑗
𝜕𝑥𝑗

(2.2.4-3)

反応進行度変数輸送式：
𝜕𝑞𝐶𝑗̃
̅ 𝐶̃ 𝑢
̃
𝜕𝜌
̅ 𝐶̃
𝜕𝜌
𝜕
̃𝐶 𝜕𝐶 ) +̃𝐶
+ 𝜕𝑥 𝑗 = 𝜕𝑥 (𝜌̅ 𝐷+ 𝜌̅ 𝑊
𝜕𝑡
𝜕𝑥
𝜕𝑥𝑗
𝑗
𝑗
𝑗

(2.2.4-4)

エンタルピー輸送式：
̃𝑢
̃
̃
̅ℎ
̃𝑗
𝜕𝜌
̅ℎ
𝜕𝜌
𝜕𝑝̅
𝜕𝑝̅
𝜕
𝜆
𝜈𝑡 𝜕ℎ
+
=
+
𝑢
̃
+
+
)
((
)
𝑗
𝜕𝑡
𝜕𝑥𝑗
𝜕𝑡
𝜕𝑥𝑗
𝜕𝑥𝑗
𝐶𝑝
𝑃𝑟𝑡 𝜕𝑥𝑗(2.2.4-5)

NO 輸送式：
𝜕𝑞𝑌𝑁𝑂𝑗
̃𝑢
̃
̃
̅𝑌
̃
𝜕𝜌
̅𝑌
𝜕𝜌
𝜕
𝜕𝑌
𝑁𝑂
𝑁𝑂
̃
̃
+ 𝜕𝑥𝑁𝑂 𝑗 = 𝜕𝑥 (𝜌̅ 𝐷+ 𝜌̅ 𝑊
𝑌𝑁𝑂 𝜕𝑥 ) + 𝜕𝑥
𝑁𝑂
𝜕𝑡
𝑗
𝑗
𝑗
𝑗

(2.2.4-6)

ここで、上付きの―は LES フィルタを表し、～はファーブル平均を意味する。上記輸送方程式より求めた流れ場の混合分率及び反応進行度変数を元に、1 次元予混合燃焼解析から得られた flamelet データベースを参照し、温度や、反応速度、化学種組成を求めることで乱流燃焼場の予測を行う。圧縮機およびタ- 47 -

ービンの動翼を含む領域は非慣性系として扱い、運動量保存式にはコリオリ力および遠心力を含む。慣性系と非慣性系との境界は sliding 境界を設け、sliding 面で慣性系と非慣性系の物理量の交換を各ステップで計算し直す[4-4]。産業用ガスタービンエンジンの解析条件を Table 2.2.4-1 に示す。
(b)-1-1-3 計算成果
産業用ガスタービンエンジンの統合シミュレーションの結果を示す。なお、計算はまだ定常に達しておらず、引き続き継続中である点に注意いただきたい。Fig. 2.2.4-2 と Fig. 2.2.4-3 に、それぞれ、圧縮機およびタービンの回転軸を通る断面の速度、圧力、温度の瞬間分布、および燃焼器の中心を通る断面の速度、混合分率（燃料と酸化剤の比）、速度の瞬間分布を示す。Fig. 2.2.4-2 で見られる断面は空気の取り入
れ口、圧縮機 2 段、燃焼器からのガスがタービンに至り、排気口に至る流路である。図から 2 段ある遠心圧縮機のそれぞれで流れが加速され、圧縮される様子が確認できる。また、圧縮に伴う温度上昇も確認できる。さらに、既燃ガスがタービンの各段を経由することにより降圧していること、圧力の降下に伴い温度が下がっていることが確認できる。一方、Fig. 2.2.4-3 は、圧縮機から出た空気が燃料と混合し、燃焼を経てタービンに向かう過程を示す断面である。圧縮機からの空気と燃料の混合ガスが燃焼器内に流入し、燃焼する様子を確認することができる。
Table 2.2.4-1 Numerical conditions for industrial gas turbine engine.
乱流モデル

LES(標準スマゴリンスキモデル)

移流項の補間法

1 次精度風上

時間積分法

1 次オイラ陰解法

時間刻み

1.0×10－6 秒

壁面の条件

速度：対数則 温度：断熱

流入条件

圧力固定(全圧 0.1MPa)

流出条件

圧力(全圧 0.1MPa)、流量固定(8.2kg/秒)

回転数

22,000rpm

燃料

都市ガス 13A

酸化剤

空気

流入ガス

酸化剤：15℃、0.1MPa、8.0kg/s
燃料：25℃、1MPa
パイロット：1.1×10－3kg/s
メイン：0.1kg/s
追焚き：2.8×10－2kg/s- 48 -

燃焼器

圧縮機

タービン
Fig. 2.2.4-1 Industrial gas turbine.

Fig. 2.2.4-2 Instantaneous distributions of velocity, total pressure and temperature on a cross section with respect
to the axis of rotation of compressor and turbine.

- 49 -

Fig. 2.2.4-3 Instantaneous distributions of velocity, mixture fraction and temperature on a cross section of
combustor.
(b)-1-2 コンポーネントシミュレーション（本計算）(b)-1-2-1 計算対象
本計算は上述の川崎重工業株式会社の産業用ガスタービンエンジンに対して、メッシュ数を 8 倍としている。メッシュ密度は FrontFlow/red の機能の一つである refiner を用いた（次節）。
(b)-1-2-2 計算方法、モデル、条件
Refiner により作成したメッシュの生成方法を示す。FrontFlow/red の refiner は前処理として読み込むメ
ッシュの辺要素の長さを 1/2 倍することにより空間のメッシュ密度を上げるものである。AMR(Adaptive
mesh refinement, [4-5])法は空間の物理量の変化等でメッシュの粗密を調整するが、ここで用いる refiner はメッシャで作成した各セルの辺の中間で分割し、セルの密度を 8 倍にする方法である。四面体要素および六面体要素を代表として refine の前後の様子を Fig. 2.2.4-4 に示す。 Refiner を用いて生成されたメッシュにより作成された産業用ガスタービンエンジンの燃焼器のメッシュ数は約 1.8 億節点、約 11 億要素である。本計算の計算条件を Table 2.2.4-2 に示す。
(b)-1-2-3 計算成果
Fig. 2.2.4-5 に、燃焼器の中心を通る断面の速度、密度、温度、および混合分率（燃料と酸化剤の比）の瞬間分布を示す。図から、圧縮機からの空気と燃料の混合と混合ガスが燃焼器内に流入し、燃焼する様子を確認することができる。本計算は十分なステップの計算が出来ていると考えられる一方、移流項のス
キームにより、解が鈍っていると判断される。今後、移流項のスキームに 2 次中心を用いることにより、詳細な渦構造と燃焼の様子を調査していく予定である。
Table 2.2.4-2 Numerical conditions for combustor.

- 50 -

乱流モデルLES(標準スマゴリンスキモデル)

移流項の補間法

2 次精度風上

時間積分法

1 次オイラ陰解法

時間刻み

1.0×10－6 秒

壁面の条件

速度：対数則 温度：断熱

燃料

都市ガス 13A

酸化剤

空気

流入条件

酸化剤：348℃、1MPa、7.4kg/s
燃料：25℃、1MPa
パイロット：1.1×10－3kg/s
メイン：0.1kg/s
追焚き：2.8×10－2kg/s

Fig. 2.2.4-4 Refined tetrahedral or hexahedral mesh with refiner included FrontFlow/red.

- 51 -

Fig. 2.2.4-5 Instantaneous distributions of velocity, density, temperature and mixture fraction on a cross section of
combustor with a refined mesh.
(b)-1-3 まとめと今後の課題
川崎重工業株式会社の産業用ガスタービンエンジンを対象に、圧縮機、燃焼器、タービンを統合したメッシュ節点数約 0.9 億、要素数約 4.1 億の簡易形状を用いて試計算を実施し、妥当な結果を得ることができた。また、refiner を用いて 8 倍の詳細形状を作成して、そのコンポーネントシミュレーションを実施
し、詳細形状を用いた統合シミュレーションへの見通しを得た。引き続き定常に達するまで計算を継続した後、今後は、移流項に対するスキームの高精度化を試みつつ、詳細形状を用いた統合シミュレーションへと取り掛かる。
(b)-2 燃焼モデル（フレームレット法）の適用性検討（2 次元解析）
本解析では、京都大学開発のインハウスコード FK3 を用いて行った[4-6]。
(b)-2-1 水素燃焼およびアンモニア燃焼(b)-2-1-1 計算対象
予 混 合 燃 焼 を 対 象 と した フ レ ー ム レ ッ ト 法 であ る Flameletgenerated manifold（FGM）法 [4-3]に、選択拡散効果と火炎伸長の影
響を考慮した手法である SPF-PD 法 (Strained Premixed Flamelet
method considering Preferential Diffusion effect)について、その水素燃
焼、およびアンモニア燃焼への適用妥当性を検証するため、円状予
混合伝播火炎の２次元数値シミュレーションを行う。Fig. 2.2.4-6 に計算領域の概略を示す。
Fig.

(b)-2-1-2 計算方法、モデル、条件

2.2.4-6

Schematicof

本研究では、SPF-PD 法を用いた計算に加えて、比較対象として、 computational domain
燃焼モデルを用いない詳細計算法、従来 FGM 法、選択拡散効果のみを考慮した FGM 法である FGM-PD法、火炎伸長のみを考慮した FGM 法である SPF 法を用いた計算を行う。

- 52 -

(b)-2-1-2-1 詳細計算法
詳細計算法では、燃焼モデルを用いず、詳細反応機構に含まれる数百の素反応を直接解く。本研究では、
詳細計算法を用いて得られた結果を真値として取り扱う。詳細計算法では、支配方程式として、以下に示
す質量、運動量、エネルギー、および化学種の質量の保存式を、理想気体の状態方程式とともに解く。
𝜕𝜌+ ∇ ∙ (𝜌𝒖) = 0
𝜕𝑡(2.2.4-7)𝜕𝜌𝒖+ ∇ ∙ (𝜌𝒖𝒖) = −∇𝑃 + ∇ ∙ 𝜏
𝜕𝑡(2.2.4-8)𝜕𝜌ℎ𝜕𝑃+ ∇ ∙ (𝜌ℎ𝒖) =+ 𝒖 ∙ ∇𝑃 + ∇ ∙ [𝜌𝐷ℎ (∇ℎ − ∑ ℎ𝑘 ∇𝑌𝑘 ) − ∑ 𝜌ℎ𝑘 𝑌𝑘 𝑽𝑘 ]𝜕𝑡𝜕𝑡

(2.2.4-9)𝑘𝑘+ 𝜏: ∇𝒖𝜕𝜌𝑌𝑘+ ∇ ∙ (𝜌𝑌𝑘 𝒖) = −∇ ∙ (𝜌𝑌𝑘 𝑽𝑘 ) + 𝜌𝜔̇ 𝑘𝜕𝑡(2.2.4-10)

ここで、𝜌は密度、𝑢は速度ベクトル、𝑝は圧力、𝜏は粘性応力、ℎはエンタルピー、ℎ𝑘 は各化学種𝑘のエンタルピー、𝐷ℎ は熱拡散率、𝑌𝑘 は化学種𝑘の質量分率、𝑽𝑘 は化学種𝑘の拡散速度、𝜔̇ 𝑘 は化学種𝑘の反応率である。𝑽𝑘 は、以下に示す Maxwell-Stefan diffusion により考慮し、温度勾配、圧力勾配、および外力による物質拡散は考慮しない。
∇𝑋𝑘 = ∑ (𝑛

𝑋𝑗 𝑋𝑘
) (𝑽𝑗 − 𝑽𝑘 )𝐷𝑗𝑘

(2.2.4-11)

ここで、𝑋𝑘 は化学種𝑘のモル分率、𝐷𝑗𝑘 は化学種𝑗の化学種𝑘に対する二成分拡散係数である。詳細計算法
においては、式(2.2.4-11)により、化学種ごとの拡散性の違い、すなわち選択拡散効果を考慮する。また、燃焼反応を考慮するため、水素燃焼においては、AramcoMech3.0 [4-7]から C 原子を含む反応を除き、
Glarborg ら [4-8]による N 原子に関連した反応を考慮した 33 化学種、214 素反応から成る詳細反応機構を用いる。アンモニア燃焼においては、32 化学種、204 素反応から成る Otomo 詳細反応機構 [4-9]を用いる。
(b)-2-1-2-2 従来 FGM 法
従来 FGM 法では、燃焼反応を考慮するため、詳細反応機構を直接解く代わりに、事前に作成したデー
タベースを参照する。このデータベースを作成する際には、火炎伸長のない１次元予混合自由伝播火炎の計算を行い、得られた火炎特性を、Control variable と呼ばれるパラメータについてまとめる。本研究では、FlameMaster コード [4-10]を用いて１次元予混合自由伝播火炎の計算を行う。この際、詳細計算法と同様の詳細反応機構を用いる。また、従来 FGM 法においては、Control variable として反応進行変数𝐶を用いる。𝐶は主に燃焼生成物の質量分率で定義され、水素燃焼では𝐶 = 𝑌𝐻2𝑂 、アンモニア燃焼では𝐶 = 𝑌𝑁2 +𝑌𝐻2𝑂 + 𝑌𝐻2により定義する。
２次元数値シミュレーションにおいては、式(2.2.4-7), (2.2.4-8)に加え、式(2.2.4-9), (2.2.4-10)を解く代わりに以下に示す𝐶の保存式を解く。
𝜕𝜌𝐶(2.2.4-12)+ ∇ ⋅ (𝜌𝐶𝒖) = ∇ ⋅ (𝜌𝐷𝑐 ∇𝐶) + 𝜌𝜔̇ 𝑐𝜕𝑡
ここで、𝐷𝑐 は𝐶の拡散係数であり、𝜔̇ 𝑐 は𝐶の反応率である。従来 FGM 法においては、2 次元数値シミュレ- 53 -

ーションにおいては、選択拡散効果は考慮されず、Lewis 数（𝐿𝑒）は１と仮定され、𝐷𝑐 は𝐷𝑐 = 𝐷ℎ で定義される。
(b)-2-1-2-3 FGM-PD 法
従来 FGM 法では、𝐿𝑒 = 1を仮定し、選択拡散効果を考慮できないのに対し、FGM-PD 法 [4-11, 4-12]では選択拡散効果を考慮する。FGM-PD 法における𝐶の保存式は以下の式で与えられる。
𝜕𝜌𝐶+ ∇ ⋅ (𝜌𝐶𝒖) = ∇ ⋅ (𝜌𝐷𝑐𝑐 ∇𝐶) + 𝜌𝜔̇ 𝑐
𝜕𝑡
ただし、拡散係数𝐷𝑐𝑐 は以下の式で定義される。
𝐷𝐶𝐶 =𝜆𝑐𝑘 𝜕𝑊𝑌𝑘∑𝑊𝜌𝑐𝑝
𝐿𝑒𝑘 𝜕𝐶(2.2.4-13)(2.2.4-14)

ここで、λ は熱伝導率、W は混合気の平均モル質量、cp は定圧比熱である。ck は化学種ごとに異なる定数
であり、水素燃焼では H2O のみ 1、アンモニア燃焼では N2 H2O, H2 のみ 1 であり、その他の化学種は 0である。また、各化学種の Lewis 数𝐿𝑒𝑘 は以下の式で定義される。
𝜆𝐿𝑒𝑘 =𝜌𝑐𝑝 𝐷𝑘𝑚
1 − 𝑌𝑘𝐷𝑘𝑚 =∑𝑗≠𝑘 𝑋𝑗 /𝐷𝑗𝑘

(2.2.4-15)(2.2.4-16)

ここで、𝐷𝑘𝑚 は化学種𝑘の混合気に対する拡散係数であり、これにより、選択拡散効果を考慮できる。
なお、FGM-PD 法において、データベース作成時に行う 1 次元計算は従来 FGM 法と同様であり、従来FGM 法のデータベースに加え、𝐷𝑐𝑐 をデータベースにまとめる。
(b)-2-1-2-4 SPF 法
従来 FGM 法では、火炎伸長のない 1 次元予混合自由伝播火炎をもとにデータベースを作成するため、火炎伸長の影響を考慮することができない。一方、SPF 法 [4-13]では、火炎伸長のない 1 次元予混合自由
伝播火炎に加え、火炎伸長のある 1 次元対向流予混合火炎の計算を行い、それらの結果をデータベースにまとめる。この際、従来 FGM 法で用いた𝐶に加え、火炎伸長に対応する Control variable が必要になる。
本研究では、水素燃焼、アンモニア燃焼いずれにおいても、𝐶に加え、混合分率𝑍を Control variable とする。なお、𝑍は Bilger の定義 [4-14]により定義される。そのため、2 次元数値シミュレーションにおいては、式(2.2.4-7), (2.2.4-8), (2.2.4-12)に加え、以下に示す𝑍の保存式を解く。
𝜕𝜌𝑍+ ∇ ⋅ (𝜌𝑍𝒖) = ∇ ⋅ (𝜌𝐷𝑍 ∇𝑍)𝜕𝑡
ここで、𝐷𝑍 は𝑍の拡散係数であり、𝐿𝑒 = 1を仮定し𝐷𝑍 = 𝐷ℎ により与える。

(2.2.4-17)(b)-2-1-2-5 SPF-PD 法
SPF-PD 法では、SPF 法と同様、火炎伸長のない 1 次元予混合自由伝播火炎、および火炎伸長のある 1次元対向流予混合火炎の計算を行い、𝐶, 𝑍についてデータベースをまとめる。2 次元数値シミュレーショ
ンにおいては、式(2.2.4-7), (2.2.4-8)に加え、以下に示す𝐶, 𝑍の保存式を解く。
𝜕𝜌𝐶+ ∇ ⋅ (𝜌𝐶𝒖) = ∇ ⋅ (𝜌𝐷𝑐𝑐 ∇𝐶) + ∇ ⋅ (𝜌𝐷𝐶𝑍 ∇𝑍) + 𝜌𝜔̇ 𝑐
𝜕𝑡- 54 -(2.2.4-18)

𝜕𝜌𝑍+ ∇ ⋅ (𝜌𝑍𝒖) = ∇ ⋅ (𝜌𝐷𝑍𝐶 ∇𝐶) + ∇ ⋅ (𝜌𝐷𝑍𝑍 ∇𝑍)𝜕𝑡
ただし、拡散係数𝐷𝐶𝑍 , 𝐷𝑍𝐶 , 𝐷𝑍𝑍 は以下の式で定義される。
𝜆𝑐𝑘 𝜕𝑊𝑌𝑘∑𝑊𝜌𝑐𝑝
𝐿𝑒𝑘 𝜕𝑍𝜆𝑧𝑘 𝜕𝑊𝑌𝑘
𝐷𝑍𝐶 =
∑
𝑊𝜌𝑐𝑝
𝐿𝑒𝑘 𝜕𝐶𝜆𝑧𝑘 𝜕𝑊𝑌𝑘
𝐷𝑍𝑍 =
∑
𝑊𝜌𝑐𝑝
𝐿𝑒𝑘 𝜕𝑍
𝐷𝐶𝑍 =(2.2.4-19)(2.2.4-20)(2.2.4-21)(2.2.4-22)

ここで、𝑧𝑘 は化学種ごとに異なる定数であり、Bilger の定義 [4-14]から決定される。
(b)-2-1-2-6 計算条件および数値解法
水素燃焼については、当量比 0.45、未燃予混合気温度 300 K、雰囲気圧力 1 atm の条件下において、詳細計算法、従来 FGM 法、FGM-PD 法、SPF 法、SPF-PD 法の５種類の計算手法を用いて数値シミュレー
ションを行う。計算領域は Fig. 2.2.4-6 に示すような正方形領域であり、1 辺は 20.48 mm である。境界条
件として x 方向、y 方向いずれも自由流出境界とする。計算格子は格子間隔 40 µm の等間隔直交スタガー
ド格子であり、火炎厚みを 10 点以上で解像可能である。計算開始時には、中心の着火領域を除いて未燃
予混合気で満たす。中心の着火領域は半径 1 mm であり、火炎伸長なしの 1 次元計算で得られた既燃ガスの温度および組成を与え、火炎を円状に伝播させる。
アンモニア燃焼については、当量比 0.8, 0.9, 1.0, 1.1, 1.2、未燃予混合気温度 673 K、雰囲気圧力 2 MPaの条件下において、詳細計算法、従来 FGM 法、FGM-PD 法、SPF 法、SPF-PD 法の５種類の計算手法を用いて数値シミュレーションを行う。計算領域、および境界条件については、1 辺 10.24 mm、格子間隔10 µm、初期着火半径 0.2 mm という点を除いて、水素燃焼と同様である。
いずれの数値シミュレーションもインハウスコード FK3[4-6]を用いて行った。
(b)-2-1-3 計算成果(b)-2-1-3-1 水素燃焼
Fig. 2.2.4-7 に、火炎半径の時間変化を、5 種類の計算手法について比較して示す。ここで、火炎半径は、火炎面の位置を反応進行変数𝐶 = 0.1の等値線として、中心から火炎面までの距離として定義する。
図より、火炎伸長の影響を考慮しない従来 FGM 法、および FGM-PD 法においては、詳細計算法の結果を過小評価することが分かる。これは、希薄水素/空気予混合燃焼において、火炎伸長を受けると燃焼速度が増加するという傾向[4-15]を、火炎伸長の影響を考慮しない従来 FGM 法、および FGM-PD 法が捉えられないためと考えられる。一方で、選択拡散効果と火炎伸長の影響を考慮した SPF-PD 法は、詳細計算法の結果と良好に一致することが確認できる。しかしながら、火炎伸長の影響のみを考慮した SPF 法は詳
細計算法の火炎伝播速度を過大評価する。この原因を調査するため、Fig. 2.2.4-8 に、火炎半径 4 mm の
時刻における混合分率𝑍の半径方向分布を計算手法について比較して示す。図より、詳細計算法では、火炎面（𝑟 = 4 mm）近傍において、混合分率が減少し、極小値をとることが確認できる。これは選択拡散効果によるものである。また図より、SPF-PD 法、従来 FGM 法、および FGM-PD 法は、この火炎面での混合分率の減少を捉えることが出来るが、SPF 法はこの火炎面近傍での混合分率の減少を再現できず、過大評価することが分かる。ここで、混合分率は当量比に相当する物理量である。そのため、燃料希薄条件- 55 -

において、火炎面での混合分率を過大評価することにより、火炎伝播速度も過小評価すると考えられる。
以上の結果から、火炎伸長を伴う希薄水素/空気予混合燃焼においては、燃焼モデルとして、選択拡散効果と火炎伸長の影響を考慮した SPF-PD 法を用いる必要があることが明らかとなった。

Fig. 2.2.4-7 Comparison of time variations of flame

Fig. 2.2.4-8 Comparison of radial profiles of mixture

radius among detailed, conventional FGM, FGM-PD,

fraction among detailed, conventional FGM, FGM-

SPF, SPF-PD methods.

PD, SPF, SPF-PD methods.

(b)-2-1-3-2 アンモニア燃焼
Fig. 2.2.4-9 に、当量比 0.8 から 1.2 の各条件において、5 種類の計算手法を用いて得られた火炎半径の時間変化を比較して示す。図より、水素火炎と比較するといずれの手法の間にも大きな差が無いが、火炎伸長の影響を考慮しない従来 FGM 法、FGM-PD 法に比べ、火炎伸長の影響を考慮した SPF 法、SPF-PD法の方が詳細計算法に近い結果を示すことが分かる。また、選択拡散効果の考慮の影響は小さく、これはNH3 の Lewis 数が 1 に近いためと考えられる。

(a) 0.8(b) 0.9- 56 -

(c) 1.0(d) 1.1(e) 1.2
Fig. 2.2.4-9 Comparisons of time variation of flame radius among detailed, conventional FGM, FGM-PD,
SPF, SPF-PD methods at equivalence ratios of (a) 0.8, (b) 0.9, (c) 1.0, (d) 1.1, and (e) 1.2.
(b)-2-2 水素－アンモニア混焼
本研究では、水素-アンモニア混焼条件において、Flamelet モデリングの観点から重要と考えられる選
択拡散の影響について、2 次元乱流混合層を対象に直接数値計算（DNS）による検討を行った。
(b)-2-2-1 計算対象
アンモニアは輸送や貯蔵の観点から水素よりもハンドリングが容易であり、水素エネルギーキャリアとして期待される一方、低い着火性や発熱量が課題として挙げられている。そこで本研究では、アンモニ
ア利用技術の形態として、燃焼前アンモニア改質水素の利用を想定し、水素-アンモニア混焼燃焼を対象とした。アンモニア改質後のガス組成が体積分率で水素：窒素が 3：1 になることを想定し、燃料組成を体積分率で、水素：アンモニア：窒素が 0.45：0.40：0.15 と設定した。計算領域は、基礎的な燃焼特性を分析するのに適する 2 次元乱流混合層とした。Fig. 2.2.4-10 に、計算領域を示す。
(b)-2-2-2 計算方法、モデル、条件- 57 -

本研究では、2 次元直接数値計算（DNS）により計算を行う。支配方程式は、質量、運動量、エネルギ
ー、および化学種質量分率の各保存式、状態方程式、Lewis 数 1 近似から導かれる混合分率輸送方程式、および反応進行変数輸送方程式から構成される[4-16]。

(2.2.4-23)(2.2.4-24)(2.2.4-25)(2.2.4-26)(2.2.4-27)(2.2.4-28)(2.2.4-29)
ここで、 は密度、 は流速、 は時間、 は圧力、 は粘性応力テンソル、 は比エンタルピー、 は熱拡散率、 は化学種質量分率、 は拡散速度、 は生成速度、 は混合分率、 は反応進行変数、k は化学種を表す。ただし、反応進行変数の定義は H2 と H2O の質量分率の合計とした。
化学種 k の拡散速度については、次式に基づいて計算する。

(2.2.4-30)(2.2.4-31)

ここで、は化学種 k の体積分率、は j 化学種と k 化学種の 2 成分系拡散係数である。本研究では、式(2.2.4-23)〜(2.2.4-26)の支配方程式の解と、Lewis 数 1 近似のもとで Flamelet progress-variable（FPV）法[4-17, 4-18]に基づいて作成された Flamelet library を式(2.2.4-27)および(2.2.4-28)から求められる検索パラ
メータにより決定されるモデルデータとの比較を行う。Lewis 数の定義は次式による。

(2.2.4-32)

計算領域は、前述の通り 0.02 m x 0.0288 m の矩形領域を 792 x 1,152 の均一な正方格子（格子幅 25μm）- 58 -

で分割した。初期条件として、鉛直方向に 3 分割し、中層に燃料（幅 0.005 m）を、上層と下層に空気を配置し、燃料層には平均的に水平方向に 80 m/s を、空気層には-2 m/s を与え、乱流混合層が形成されるようにした。また、初期流速には、Fig. 2.2.4-11 に示す 3 次元空間に一様等方性乱流場を作成し[4-16,4-19]、
2 次元空間へのプロジェクションを行った。乱流の条件としては、乱れ強さ

m/s、乱流渦の積分mm とした。着火のため、燃料層と空気層の間に Flamelet library における既燃条件のスケール

ごく薄い層（幅 0.0005 m）を配置した。用いる化学反応メカニズムには、Jiang らの 19 化学種 60 素反応で構成されるものを採用した[4-20]。Fig. 2.2.4-12 に、作成した Flamelet library を構成するガス温度の混合分率-反応進行変数空間における分布を示す。
計算条件は前述の通り、式(2.2.4-23)〜(2.2.4-26)を対象に化学反応をアレニウス型式を解いて得られた解（ARF）と、式(2.2.4-27)および(2.2.4-28)の解をもとに Flamelet library を検索した結果（FLM）の 2 ケー
スを比較する。これにより、詳細な拡散流束を考慮する ARF と Lewis 数 1 近似を用いる FLM の差異を議論する。
DNS には、圧力ベースの半陰解法コード FK3[4-6]を用いた。運動量保存式の対流項には保存型 4 次精度中心差分を、スカラー保存式のそれには WENO5 スキームを、各拡散項には 4 次精度中心差分を用いた。時間積分には 3 次精度陽的 TVD Runge-Kutta 法を用いた。化学反応計算の常微分方程式には VODE法を用いた。計算は 576 コア並列計算により実施し、時間刻み 1.0 x 10-8 s にて 100,000 ステップまで約
8 時間の計算時間を要した。

Fig. 2.2.4-10 Initial condition of 2D DNS.

- 59 -

Fig. 2.2.4-11 3D domain for generating isotropic
turbulence.

Fig. 2.2.4-12 Temperature distribution in mixture fraction and progress variable space in Flamelet library.
(b)-2-2-3 計算成果
Fig. 2.2.4-13 に、ARF と FLM のガス温度の瞬間分布を示す。図より、FLM のガス温度は、その分布形状は ARF のそれにおよそ一致するものの、温度レベルが低く捉えられていることがわかる。Lewis 数 1
近似はガス温度を過小評価する可能性があることは、例えば Han らの検討[4-21]でも報告されており、本研究でも改めて確認された。
Fig. 2.2.4-14 に、ARF と FLM の NO 質量分率の瞬間分布を示す。図より、FLM の NO 質量分率は、その分布形状は ARF のそれにおよそ一致するものの、濃度レベルが低く捉えられていることがわかる。通常、FPV 法に基づく燃焼計算において定常 Flamelet 法を用いる場合、NO 濃度は過大に予測される傾向があることが知られている[4-22]。これは、NO 生成反応よりもはるかに短い時間スケールをもつ燃料の酸化反応に関する反応進行変数に基づき Flamelet library が作成されるためであるが、翻って Fig. 2.2.4-14を見てみると、FLM は ARF の予測値を過小評価していることから、上述の FPV 法の特徴を Fig. 2.2.4-13に見られた過小評価された温度レベルが凌駕する結果となっていると考えることができる。
Fig. 2.2.4-15 に、ARF と FLM のガス温度と NO 質量分率について、混合分率空間の散布図を示す。図より、ガス温度のピーク値に 200 K 以上の差異があることがわかる。特に、希薄側に顕著な差異が見られる。希薄側には拡散速度の速い H2 の影響は大きくないものと思われるが、着火やそれに続く連鎖反応において、OH ラジカルや H ラジカル等は重要な役割を果たしていると考えられる。
本研究において明らかとなったガス温度や NO 質量分率の予測値における ARF と FLM の差異は極めて大きく、水素-アンモニアの混焼燃焼場において、Lewis 数 1 近似を用いることの危険性が確認された。
今後は、例えば Yu らの Constant Lewis 数近似[4-23]等、より高精度な手法について、その適用性を検討していく必要があると思われる。

- 60 -

Fig. 2.2.4-13 Instantaneous temperature distribution, (a) ARF, (b) FLM.

Fig. 2.2.4-14 Instantaneous distribution of NO mass fraction, (a) ARF, (b) FLM.

- 61 -

(a)(b)

Fig. 2.2.4-15 Scatter plots of temperature and NO mass fraction, (a) ARF, (b) FLM.
(b)-2-3 まとめと今後の課題
2 次元数値シミュレーションによる水素、アンモニア乱流燃焼モデルの精度検証を行い、燃焼モデルとして利用するフレームレット法には非断熱性、選択拡散、火炎伸長の効果を考慮することが燃焼器内乱流
燃焼場の高精度予測に有効であることを明らかにした。今後は、NO 濃度予測の更なる高精度化を検討するとともに、3 次元数値シミュレーションによる乱流燃焼モデルの妥当性および精度の検証を行う。
参考文献[4-1] H. Moriai et al., J. Eng. Gas Turb. Pow., 135, 091503, (2013).
[4-2] S. Tachibana et al., Combust. Flame, 162, 2621-2637, (2015).
[4-3] J.A. van Oij, and L.P.H. de Goey, Combust. Sci. Tech., 161, 113-137, (2000).
[4-4] M. M. Rai and N. K. Madavant, AIAA paper, No. 88-0361, (1988).
[4-5] M. J. Berger and J. Oliger, J. Comp. Physics, 53, 484-512, (1984).
[4-6] About in-house code FK3, http://www.tse.me.kyoto-u.ac.jp/members/kurose/link.php[4-7] C. W, Zhou et al., Combust. Flame, 197, 423-438 (2018).
[4-8] P. Glarborg et al., Prog. Energy Combust. Sci., 67, 31-68 (2018).
[4-9] A. Otomo et al., Int. J. Hydrog. Energy, 43, 3004-3014, (2018).

- 62 -

[4-10] H. Pitsch, FlameMaster, a C++ computer program for 0D combustion and 1D laminar flame calculations.
[4-11] R. Bastiaans et al., CTR Annual Research Briefs 2007, 195, (2007).
[4-12] J. A. de Swart et al., Combust. Sci. Technol., 161, 113-137 (2010).
[4-13] E. Knudsen et al., Combust. Flame, 160, 2911-2927 (2013).
[4-14] R. W. Bilger, Proc. Combust. Inst., 22, 475-488, (1989).
[4-15] C. K. Law, Proc. Combust. Inst. 22, 1381-1402, (1989).
[4-16] H. Shehab et al., Combust. Flame, 238, 111888, (2022).
[4-17] C. D. Pierce, P. Moin, J. Fluid Mech., 504, 73-97, (2004).
[4-18] P. Yu, R. Kurose, H. Watanabe, Phys. Fluids, 35, 015116, (2023).
[4-19] T. Passot, A. Pouquet, J. Fluid Mech., 181, 441-446, (1987).
[4-20] Y. Jiang et al., Int. J. Energy Res., 44 795-810, (2020).
[4-21] W. Han et al., Proc. Combust. Inst., 38, 2627-2634 (2021).
[4-22] M. Ihme, H. Pitsch, Phys. Fluids, 20, 055110, (2008).
[4-23] P. Yu, H. Watanabe, Proc. JSME FED Conf., Kumamoto (2022).
(5)プロジェクトの総合的推進
以上のサブ課題 A～D の進捗を管理し、課題間の連携を促進するために、各課題の代表者を集めて統括する推進委員会を令和 4 年 5 月、9 月、12 月、3 月の計 4 回、実施した。また、年度末の令和 5 年 3 月 3日に、プロジェクトの成果を報告するための「スマートデザインシンポジウム」を実施した他、令和 4 年
9 月 30 日に、サブ課題間の連携を促進するための「スマートデザインワークショップ」を、他の「産業競争力の強化」課題と連携した統合ワークショップを令和 4 年 9 月 27 日に開催した。

- 63 -

２－３．活動（研究会の活動等）(1) アウトリーチ活動
年 月 日

行

事

R4.9.27「富岳」成果創出加速プログラム第６回 HPC ものづくり統合ワークショップ

R4.9.30

第 1 回スマートデザインプロジェクトワークショップ

R5.3.3

第 1 回スマートデザインプロジェクトシンポジウム

R6.3.7「富岳」成果創出加速プログラム研究交流会(2) プロジェクト推進会議
年 月 日

行

事

R4.5.31

第 4 回「富岳」スマートデザインプロジェクト推進会議

R4.9.12

第 5 回「富岳」スマートデザインプロジェクト推進会議

R4.12.14

第 6 回「富岳」スマートデザインプロジェクト推進会議

R5.3.3

第 7 回「富岳」スマートデザインプロジェクト推進会議(3) コンソーシアム活動
年 月 日
R4.10.14
R5.3.22R4.5.18
R4.12.12
R4.4.22
R4.11.25

行事

都市・建築 CFD コンソーシアム 総会
HPC を活用した自動車用次世代 CAE コンソーシアム総会
燃焼システム用次世代 CAE コンソーシアム総会(4) ワーキンググループ
年 月 日

行

事

R4.5.20R4.8.8R4.9.22
R4.11.1

サブ課題Ａ サブ課題全体研究打合せ（全７回）

R4.12.26
R5.1.27
R5.3.16R4.4.8R4.4.15

サブ課題Ａ 最適化打ち合わせ（全２７回）- 64 -

R4.4.22
R4.5.12R4.6.3R4.6.10R4.6.17R4.6.24
R4.7.22R4.8.5R4.8.17R4.9.9R4.9.16
R4.10.7
R4.10.14
R4.10.21
R4.10.28
R4.11.11
R4.11.18
R4.12.2
R4.12.13
R4.12.20
R5.1.10R5.1.17R5.1.24
R5.1.31R5.2.7R4.4.15
R4.4.22R4.5.13R4.5.30
R4.6.3R4.6.10R4.6.17

サブ課題Ａ デザイン・空力モデル化打ち合わせ（全２７回）

R4.6.24
R4.6.27
R4.7.8R4.7.15
R4.7.29R4.8.5- 65 -

R4.9.9R4.9.16
R4.10.7
R4.10.14
R4.10.21
R4.11.11
R4.11.29
R4.12.6
R4.12.13
R4.12.20
R4.12.27
R5.1.12R5.1.24
R5.2.7R4.4.13R4.5.20
R4.6.27
R4.8.10R4.9.14R4.10.26

サブ課題 B サブ課題全体研究打合せ（研究連絡会）（全 10 回）

R4.12.8
R5.1.11R5.2.22R5.3.24
R4.5.9R4.6.6R4.7.26
R4.9.2R4.10.14
R4.11.22

サブ課題 B 極端気象流入変動風整備 WG（全 10 回）

R4.12.2
R5.1.10R5.2.3R5.3.21R4.5.10R4.6.10R4.7.15

サブ課題 B 台風時の都市・建築における機能性・安全性 WG（全 9 回）R4.9.9- 66 -

R4.10.21
R4.11.22
R4.12.23
R5.1.30
R5.3.10R4.4.27
R4.5.27
R4.6.24
R4.7.22R4.8.26
R4.10.7

サブ課題 B 猛暑時の熱環境推定・安全性指標 WG（全 10 回）

R4.11.11
R4.12.23
R5.1.27
R5.3.10R4.4.8R4.5.17R4.6.17R4.7.25
R4.8.24
R4.9.21

サブ課題 B 早期発展 WG（全 11 回）

R4.10.21
R4.11.22
R4.12.23
R5.2.3R5.3.10R4.4.7R4.4.20
R4.5.10R4.5.25
R4.6.2R4.6.14R4.6.29

サブ課題 C サブ課題全体研究打合せ（全２４回）

R4.7.19R4.8.3R4.8.26
R4.9.9R4.9.22- 67 -

R4.10.6
R4.10.20
R4.11.4
R4.11.22
R4.12.7
R4.12.28
R5.1.13R5.1.30
R5.2.16R5.2.27
R5.3.15
R5.3.28
R4.4.11R4.5.30
R4.9.9R4.10.7

サブ課題 D サブ課題全体研究打合せ（全６回）

R4.11.10
R5.1.20- 68 -

２－４．実施体制
業務項目

担当機関

担当責任者

（１）サブ課題 A「意匠空間理化学研究所計算科学研究センを考慮した AI 支援多目的最ター適化による自動車空力デザイ

神戸大学大学院

教授 坪倉誠広島大学大学院

助教 中島卓司

東京工業大学大学院教授 坂田弘安

チームリーダー 坪倉誠ン」（２）サブ課題 B「変容する都市・建築の自然擾乱対応の性能設計」

研究員 田村哲郎
横浜国立大学大学院

教授 山田貴博

近畿大学

教授 道岡武信

お茶の水女子大学

准教授 河合英徳（３）サブ課題 C「新興ウイ九州大学総合理工学研究院

教授 伊藤一秀

ルス感染症にロバストで健

理化学研究所計算科学研究セン

チームリーダー 坪倉誠

康・快適・サステイナブルなター

ポストコロナ時代の室内環境

神戸大学大学院

教授 坪倉誠（４）サブ課題 D「マルチコ京都大学大学院教授 黒瀬良一

ンポーネント統合シミュレー

九州大学総合理工学研究院

教授 渡邊裕章

設計」

ションによるカーボンフリー
ガスタービン設計」- 69 -

別添１ 学会等発表実績
１．学会等における口頭・ポスター発表
発表した成果（発表題目、
口頭・ポスター発表の別）

発表者氏名発表した場所

発表し 国内・（学会等名）た時期 外の別

日本応用数理学会 ものづ

2022年6

ョン：Society 5.0 時代のくり企業に役立つ応用数月

スマートデザインを目指し

理手法の研究会 第 44 回て」口頭

技術セミナー「燃焼の数値シミュレーシ

黒瀬良一「Data-driven ammonia

M. Saito, J. Nagao,

75th Annual Meeting of

2022年

combustion simulation

J. Xing, R. Kurose

the American Physical

11月

using Neural Ordinary

Society’s Division of Fluid

Differential Equations」口Dynamics (APS DFD)

国内

国外

頭
「A study on FGM method

K. Reo, K. Kinuta,

The 11th Asian Joint

2023年3considering preferential

R. Kurose

Conference on Propulsion

月

diffusion for numerical国外

and Power (AJCPP2023)

simulation of hydrogen/air
premixed flame」口頭
エッサム神田ホール1号館

2022年5

沫・エアロゾル感染リスク（国際呼吸保護学会日本

月

評価のデジタルトランスフ支部2022年度記念講演ォーメーション ～我々は会）

スパコン「富岳」による飛

坪倉 誠

国内

どのようにマスクの感染リ
スク低減効果をシミュレー
ションで評価したか？～,口頭（記念講演）
20th Conference on The

2022年5droplet/aerosol infection

International Society for

月riskassessment andRespiratory Protection

evaluation of(ISRP2022)

Digital transformation of

坪倉 誠

国外

coutnermeasures realizedon the supercomputer“Fugaku” ,口頭
パシフィコ横浜パシフィ

2022年6

沫・飛沫核感染リスク評価

コ

月のデジタルトランスフォー（横浜第37回日本環境感メーション,

染学会総会・学術集会シ

スパコン「富岳」による飛

坪倉 誠- 70 -

国内

口頭

ンポジウム「COVID-19の効果的な予防策をめぐり」）

Evaluation of Infection

Rahul Bale, Chong-Oslo, Norway

2022年6

Risk due to Airborne Virus

Gang Li, Hajime（Interdisciplinary

月

Transmission in a

Fukudome, Saorichallenges towards

Restaurant Environment,

Yumino, Akiyoshi

exascale fluid dynamics

口頭

Iida, and Makoto(MS154), 8th European

TsubokuraCongress on(Keynote Lecture)

Computational Methods

国外

in Applied Sciences andEngineering）
Royal Lahaina Resort,

2022年7and smart design in

Maui, Hawaii

月

manufacturing process（The 11th International

realized on the

Conference ofsupercomputer “Fugaku” ,Computational Fluid

口頭

Dynamics (ICCFD11)）

Digital transformation

坪倉 誠

Osaka, Japan，

2022年7

Simulation on Massively（The twelfth

月

Parallel Environmentsinternational symposium

toward Next-Generation

on Turbulence and Shear

Computer-Aided

Flow Phenomena

Engineering,(TSFP12)）Engineering Tubulence

坪倉 誠

国外

国内

口頭
神戸,

2022年9

感染リスク評価のデジタル（令和4年度空気調和・衛月

トランスフォーメーション

生工学会大会）『富岳』による新型コロナ

坪倉 誠

国内

とウイズコロナ時代の室内
環境デザイン,
口頭
ITフォーラム2023（主2023年2

催：情報処理学会）

月

芝浦工業大学大宮キャン

2023年3

沫感染リスク評価のDXとパス,

月

統合的感染リスク評価シス（2023年電子情報通信学

飛沫・エアロゾル拡散シミ

坪倉 誠ュレーションによる新型コ

国内

ロナ感染リスク評価,
口頭
スパコン「富岳」による飛

坪倉 誠- 71 -

国内

テムの構築,会総合大会, 公開講座「次口頭のパンデミック被害を抑えるには？」）

Multi-Objective

TakujiFahrzeug-Aerodynamik -2022年1

Aerodynamic ShapeNAKASHIMA,

Neue Chancen und

月

Optimization of Simplified

Taichi INUI,

Herausforderungen

Automobile Using HPC,

Yusukedurch Elektromobilität

口頭

NAKAMURA,

und CO2-Gesetzgebung

Keigo Shimizu,(Vehicle Aerodynamics -Takenori Hiraoka,

New opportunities andHiroaki Fukumoto,challenges due to

Akira Oyama,electromobility and CO2

Takahidelegislation)

国外

Nouzawa, Makoto
Tsubokura
Real-World Aerodynamics

Ryo Miyazaki,

横浜

2022年7

Assessment for Road

Akinori（WCCM-APCOM

月31日-

Vehicle Development

Yoshikawa, Rahul

2022）

8月5日

Realized on theBale, KunihikoSupercomputer "Fugaku",

Yoshitake,

口頭

Fortunato Nucera,

国内

Takashi Yoshino,
Makoto Tsubokura“Nonlinear Dimensionality

Kazuto Ando, Keiji

Oslo, Norway,

2022年6

Reduction for Three-

Onishi, Rahul

The 8th European

月

dimensional Flow Field

Bale, Akiyoshi

Congress on

Around Circular Cylinder

Kuroda, andComputational Methods

with Distributed Parallel

Makoto Tsubokura

in Applied Sciences andMachine Learning on

Engineering (ECCOMAS

Fugaku,” ,Congress 2022)

国外

口頭
Investigation of different

Boqi Ren, Chung-13th Asian

2022年

irregular

Gang Li, Makoto

Computational 東北大学

11月

roughness effects on

Tsubokura（Fluid Dynamics
Conference(ACFD)）

laminar naturalconvection flows in cuboid
cavity.,口頭- 72 -

国内

Numerical simulation of

HsuehJui Lu,

13th Asian

2022年

velopharyngeal closure

TsukasaComputational Fluid

10月

and tongue elevation

Yoshinaga,

Dynamics Conference.

during the soundChungGang Li,(ACFD 2022)

production of [s] in /usui/.

Kazunori Nozaki,

国外

Akiyoshi Iida, and
Makoto Tsubokura
Principal Component

Takuji Nakashima,Jeju, Korea.

2022年

Analysis of Flow Field

Kohei Seo,

（The 13th Asian

10月

around Automobile

Takenori Hiraoka,

Computational Fluid

Models under a Crosswind

TakahideDynamics Conference

Condition,

Nouzawa, Keigo(ACFD2022)）

口頭

Shimizu, Yusuke

国外

Nakamura, Akira
Oyama, Gentaro
Hamada and
Makoto Tsubokura
Flow Field Data Mining of

Takuji Nakashima,

Yokohama, Japan.

2022年7

Automobile Models with

Yusuke Nakamura,（15th World Congress

月,

Pareto-Optimal

Keigo Shimizu,

on Computational

Aerodynamic Shape Using

Takenori Hiraoka,

Mechanics (WCCM-XV)

Proper Orthogonal

Gentaro Hamada,and 8th Asian PacificDecomposition,

Takahide

Congress on

口頭

Nouzawa, Akira

Computational

Oyama and Makoto

Mechanics (APCOM-

Tsubokura

VIII)）“Three-dimensional Flow

Kazuto Ando, Keiji

15th World Congress on

2022年7

Field Decomposition with

Onishi, Rahul

Computational

月

Massively Parallel

Bale, Akiyoshi

Mechanics (WCCM-XV)Distributed Learning on

Kuroda, andFugaku,” ,Makoto Tsubokura

国内

国外

口頭
"Voxel topology

Y.Wada,横浜

2022年7

optimization of vehicle

T.Shimada,（WCCM XV & APCOM

月

frame structure subject toK.Nishiguchi,

VIII Book of Abstracts,）

multiple loading using

S.Okazawa, andbuilding cube method

M.Tsubokura

framework," ,- 73 -

国内

口頭“Improvement of

Kazuto Ando, Keiji

33rd International

2022年5

Reduction Performance ofOnishi, Rahul

Conference on Parallel

月

Mode Decomposition for 3-

Bale, Akiyoshi

Computational Fluid

Dimensional Flow Field by

Kuroda, andDynamics (ParCFD2022)

Fugaku,” ,Makoto Tsubokura

国外口頭“富岳での分散並列機械学習

安藤和人, 大西慶治,

計算工学講演会論文集,

2022年

国内

による円柱周り3次元流れ

バレラフール，黒田

Vol.27, (2022)

場の非線形次元縮約,” 口頭明義, 坪倉誠,"車体構造のトポロジー最適

和田有司, 嶋田宗将,

計算工学講演会論文集,2022 年

国内化における板金形成の促進

西口浩司, 岡澤重信,

vol. 27, F-03-06

6月

",,口頭

坪倉誠

“ニューラルネットワークを

安藤和人, 大西慶治,

京都2022 年使用した非線形縮約モデル

黒田明義, バレラフ

（日本流体力学会年会

9月の適用性の検討, ” 口頭,ール, 坪倉誠

2022）“分散深層学習を用いた非線

安藤和人, 大西慶治,

2022 年度統計関連学会連2022 年形モード分割及び縮約モデ

バレラフール, 黒田

合大会

9月

ル構築, ” ,明義, 坪倉誠

国内

国内

口頭
Large eddy simulation of

Alicia Murga,

神戸大学

2022 年

droplet transport andChungGang Li,（空気調和衛生工学会）

9月

inhalatio exposure

Rahul Bale,

講演論文 IS-4

between interacting

Kazuhide Ito,virtual humans.

Makoto Tsubokura

Evaluation of indoor

RahulBale, Chung-Proceedings of the 36th

2022 年

environment safety

GangLi,CG.,Kazuh

Computational Fluid

12 月

through simulation of CO2

ide Ito,Makoto

Dynamics Symposium,

exhalation andTsubokura

E02-1

実走行状態における空気抵

吉川暁典，Rahul

第 36 回数値流体力学シン

2022 年

抗の予測シミュレーション,

Bale，Kim

ポジウム，C02-3

12 月

口頭

Sangwon，吉武 邦

第 36 回数値流体力学シン

2022 年

国内

国内

ventilation,
口頭
国内

彦，Fortunato
Nucera，吉野崇，坪倉 誠自然風を考慮した自動車レ

相本 翔太, Rahul- 74 -

国内

ーンチェンジ運動空 双方

Bale, Kim

ポジウム,

向連成解析,

Sangwon, 吉武邦

C02-2

口頭

彦, Fortunato

12 月

Nucera, 吉野 崇,
坪倉 誠
葺合高等学校高校生国際

2022 年

会議

7月

第 37 回経営に貢献する

2022 年

時代のスマートデザインとCAE セミナー（富士通主7月

感染リスク評価のDX,

催）

スパコン「富岳」による新坪倉誠型コロナとの闘いと「アフ国内

ターコロナ」時代のものづ
くり,口頭「富岳」が拓くSociety 5.0

坪倉誠

国内

口頭
RIKEN の世界 計算科学

2022 年

ロゾル感染デジタルトラン

研究センター ～計算の8月

スフォーメーション,

計算による計算のための口頭

科学～, 神戸市医療センタ「富岳」による飛沫・エア

坪倉誠

国内

ー
「富岳」が拓くSociety 5.0

第 9 回 荏原オープンイ

2022 年

時代のスマートデザイン

ノベーションシンポジウ

10 月～自動車CAE/CFDを例に

ム

坪倉 誠

国内

～,
口頭
「富岳」が拓くSociety 5.0

VINAS Online Users

2022 年

Conference 2022

10 月

帝国ホテル（東京），リー2022 年

テル宴会場の感染リスク評

ガロイヤルホテル（大10 月

価と安心・安全な宴会に向阪）けた対策提案,

安心・安全な立食宴会の口頭体験会（ホテル協会主

坪倉 誠

時代のスマートデザイン 2,

国内

口頭
スパコン「富岳」によるホ

坪倉 誠

国内

催）
日本機械学会関西支部第2022 年のづくり ～自動車CAEを

382 回講演会 実務者のた

11 月

例に～,

めの流体解析技術の基礎

口頭と応用(日本機械学会関西

大規模CFDが拓く新たなも

坪倉 誠

支部主催)

- 75 -

国内

第八波にむけた感染拡大抑

坪倉 誠

東京都モニタリング会議

2022 年

国内

12 月

止のための「富岳」飛沫シ
ミュレーション,口頭
Inflow turbulent database

Hidenori

Kawai，

for urban area based on

Tetruro Tamura

EACWE

2022.9

国外

WCCM-APCOM

2022.8

国内

WCCM-APCOM

2022.8

国内

WCCM-APCOM

2022.8

国内

WCCM-APCOM

2022.8

国内

WCCM-APCOM

2022.8

国内

WCCM-APCOM

2022.8

国内

large-scaled LES
including meteorologicaldisturbance（口頭）A Study of LES Coupling

Maiko Arai,with Thermal Radiation

Hidenori Kawai,for Actual Urban District

Tetsuro Tamura

Local Peak Pressure on

Azusa Ono,

Super High-rise Building

Tsuyoshi Nozu,

in Actual Urban Area

Tetsuro Tamura,
Hidenori Kawai

Study of Aerodynamic

Yuki Nagao,

Characteristics ofNaohiroEccentric Tapered Square

Nakagawa, Yusuke

Cylinder – Analysis of

Maruyama,

Flow field using BCM

Tetsuro Tamura,
Hidenori Kawai,
Masaharu
Kawaguchi,

Large-scale two-phase

Koji Nishiguchi,

flow simulation using

Masashi Morishita,building-cube method for

Tokimasa

urban flooding,

Shimada, Tetsuro
TamuraLES Analysis of

Masaharu

Ventilation Performance

Kawaguchi,

and Wind Gust

Tetsuro Tamura

Occurrence for Strategic
Urban Transformation,
WRF-LES Simulation of

Keigo Nakajima,

Wind Flow over Rough

Hidenori Kawai,

Urban Surface during

Masaharu- 76 -

Typhoon Lan (2017)

Kawaguchi,
Tetsuro Tamura,
Koji Kondo,
Yoshiaki Itoh,Kenji Takagi

LES around a Realistic

Tsuyoshi Nozu,

City Block Designed

Azusa Ono, Tetsuro

Based on a Future City

Tamura andConcept

Hidenori Kawai

Wind Pressure

Koji Kondo,

Characteristics of High-

Hidenori Kawai,

rise buildings in Middle

Tetsuro Tamura,and High-height Urban

Keigo Nakajima

WCCM-APCOM

2022.8

国内

WCCM-APCOM

2022.8

国内

WCCM-APCOM

2022.8

国内

2022.7

国内

2022.7

国内

2022.5

国外

2022.12

国内

Areas Spread over Local
Terrain,
LES on wind pressure

Hidenori Kawai,

acting on high-rise

Tetsuro Tamura,building under strong

Masaharu

wind events of Typhoon,

Kawaguchi

LES for urban boundary

Hidenori Kawai,

12th International

layer using inflow

Tetsuro Tamura,

Symposium on

condition including

Keigo Nakajima,

Turbulence and Shear

meteorological

Flow Phenomena

disturbance of typhoon
LES-based generation ofMasaharu

12th Internationaltime sequential data of

Kawaguchi,

Symposium on

instantaneous urban wind

Tetsuro Tamura

Turbulence and Shearin typhoon -comparison

Flow Phenomena

with observational data,
LES coupling with

Maiko Arai,

14th Americas

thermal radiation for

Hidenori Kawai,

Conference on Wind

actual urban district -Tetsuro Tamura

Engineering

解析条件の異なるLES実用

河合英徳 , 田村哲

第 26 回数値流体力学シンモデルに基づく乱流境界層郎, 畔見徳真, 呉嘉

ポジウム

中の三次元角柱の風圧力の瑞, 長尾悠生:

unsteadiness of turbulent
and thermal field aroundhigh-rise buildings-- 77 -

検証
地下空間浸水時の流体力評

永野雄一, 西口浩司,

第 26 回数値流体力学シン価における三次元乱流解析

田村哲郎

ポジウム

偏心を持つテーパー角柱構

長尾悠生, 丸山勇祐,

第 26 回数値流体力学シン

造物周辺の空力特性 に関田村哲郎, 河合英徳,

ポジウムする研究ーBCMを用いた流

川口真晴

2022.12

国内

2022.12

国内

2022.12

国内

の有用性，

れ場の解析
市街地における台風時最大

田村哲郎，川口真晴

瞬間風速の発生に関する物

第 26 回数値流体力学シン
ポジウム

理機構
河合英徳，田村哲郎

日本流体力学会年会

2022.9国内

川口真晴，田村哲郎

日本流体力学会年会

2022.9

国内

河合英徳，田村哲郎

日本建築学会

2022.9国内

川口真晴，田村哲郎

日本建築学会

2022.9国内2018 年猛暑日の東京都市

新井舞子，河合英日本建築学会

2022.9

国内域における気象擾乱を考慮

徳，田村哲郎河合英徳，田村哲郎

第 27 回計算工学講演会

2022.6

国内

LESによる実在都市の超高

小野 梓，野津剛，

第 27 回計算工学講演会

2022.6

国内層建築物の周辺気流性状と

田村哲郎，河合英徳
第 27 回計算工学講演会

2022.6

国内

実建物モデルを対象とする
広範囲・高解像度LES -壁面風圧変動を引き起こす乱流構造の抽出変容する実市街地のキャノ
ピー上部 ・内部に着目した換気性能評価，
台風時の気象擾乱を考慮した都市広域解析に基づく流入変動風データベースの作
成
台風下の高層建築物に作用するピーク風圧の LES 解析

した流入変動風に関する研究東京都心部の流入変動風データベースを用いた超高層建築物の風圧力に関する
LES

風圧分布特性の検討
都市変容Strategyに向けた

川口真晴，田村哲郎

換気性能と突風発生に関す- 78 -

るLES解析
Building-cube法を用いた

西口浩司，森下将

大規模浸水解析吏，嶋田宗将，田村

第 27 回計算工学講演会

2022.6

国内

発表した場所

発表し

国内・（学会誌・雑誌等名）

た時期

外の別
国外

哲郎

２．学会誌・雑誌等における論文掲載
掲載した論文（発表題目）

発表者氏名

LES flamelet modeling of

R. Kai, T. Tokuoka,

International Journal of

2022年

hydrogen combustion

J. Nagao, A. L.

Hydrogen Energy, 48,

12月

considering preferential

Pillai, R. Kurose

11086-11101 (2023)

Flame self interaction

U. Ahmed, S. P.

Physical Review Fluids,

2023年

during turbulent

Malkeson, A. L.

8, 023202 (2023)

2月

boundary layer flashback

Pillai, N.

of hydrogen-rich premixed

Chakraborty, R.

combustion

KurosePollutant capture

Ryota Muta,

Indoor and Built

2022年

efficiencies in and around

Juyeon Chung,

Environment, Volume

7月

the opening-surface of aCong Li, Sung-Jun

31, Issue 6, July 2022,

fume hood under realistic

Yoo, Kazuhide Ito

Pages 1636-1653.

Breathing zone andKazuki Kuga,

Indoor Air, 2022, 32,

exhaled air re-inhalation

Pawel Wargocki,

e13003rate under transient

Kazuhide Ito

diffusion effect国外

国外

conditions
2022年

国外

2022年

国外

conditions assessed with acomputer-simulated
personThe effects of ventilation

Xiaojun Fan, Huiqi

Building andand temperature on sleep

Shao, Mitsuharu

Environment, 209 (2022)

quality and next-day work

Sakamoto, Kazuki

108666

performance: pilot

Kuga, Li Lan,measurements in a

David P. Wyon,

climate chamber

Kazuhide Ito,
Mariya P.
Bivolarova, Chenxi- 79 -

Liao, and Pawel
Wargocki
CO2 Emission rates from

Mitsuharu

Building andsedentary subjects under

Sakamoto, Kazuki

Environment, 211 (2022)

controlled laboratory

Kuga, Kazuhide

108735conditions

Ito, Gabriel Beko,

2022年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

Mengze Li,
Jonathan Williams
and Pawel
Wargocki
Dispersion characteristics

U Yanagi,

Japan Architecturalof oral microbial

Shinsuke Kato,

Review, 2022, 5 (2), 225-

communities in a built

Hideaki Nagano,

232

environment

Kazuhide Ito,
Toshio Yamanaka,
Yoshihisa Momoi,
Hikaru Kobayashi,
Hirofumi Hayama

Exhaled Aerosol and Jet

Hana Salati,

Aerosol and Air Quality

Flow Characterization

MehrdadResearch, 22 (4),2022,

During Nasal Sneezing

Khamooshi,

210338

Jingliang Dong,
Kazuhide Ito,
David Fletcher,Sara Vahaji, Kiao
Inthavong
SARS-CoV-2 Dynamics in

Hanyu Li, Kazuki

Sustainability, 2022, 14

the Mucus Layer of theKuga, Kazuhide Ito(7), 3896

Multiscale analysis of

Ryota Muta, Sung-

Indoor and Built

material flow and

Jun Yoo, Hyuntae

Environment, Volume

computational fluid

Kim, Toru

31, Issue 9, November

dynamics for predicting

Matsumoto,2022, Pages 2291-2311

individual DEHP exposure

Kazuhide Ito

Human Upper Respiratory
Tract Based on Host-Cell
Dynamics

concentration in indoors
Rigorous mathematical

Naoki Ikegaya,- 80 -

Indoor Air. 2022;32 (7):

formulation of net escapeMats Sandberg andvelocity and net escapeKazuhide Ito

e13072

probability determining amacroscopic concentration
Prediction of exhaled

Kazuki Kuga,

Indoor Air, 2022;32 (8)carbon dioxide

Mitsuharu

e13079concentration using a

Sakamoto, Pawel

computer-simulated

Wargocki,

person that included

Kazuhide Ito2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

alveolar gas exchange
Numerical models for

Hiroshi

Journal of Building

seamlessly predicting

Harashima, Eisaku

Engineering, 57 (2022)

internal diffusion and re-Sumiyoshi, and

104976emission of leaked liquid

Kazuhide Itotoluene from indoor
mortar materials
Transport and deposition

Nguyen Dang

Japan Architecturalof inhaled man-madeKhoa, Nguyen Lu

Review 2022, Vol.5, no.4,

vitreous and asbestos

Phuong, Ken

pp 592–608

fibers in realistic humanTakahashi andrespiratory tract models:

Kazuhide Ito

An in-silico study
Error analysis of human

Yukun Wang,

Building andinhalation exposure

Xiong Shen, Sung-

Environment, 224 (2022)

simulation in industrial

Jun Yoo, Zhengwei

109573workshopLong, Kazuhide ItoInternal diffusion and re-

HiroshiJapan Architecturalemission of leaked liquid

Harashima, Eisaku

Review, 2022, Vol.5, no.4,

ethyl acetate from mortar

Sumiyoshi, andpp 672–681

materials

Kazuhide Ito

Microplastics / Microfibers

Eunsu Lim,

Japan Architecturalin settled indoor house

Hirofumi Tanaka,

Review, 2022, Vol.5, no.4,

dust- Exploratory case

Gen Ni, Yuen Bai,

pp 682–690

study for 10 residential

Kazuhide Itohouses in the Kanto areaof JapanInfluence of inlet

Haruna Yamasawa,

Japan Architectural

turbulent condition on the

Teruaki Hirayama,

Review, 2022, Vol.5, no.4,- 81 -

formation mechanism of

Ryota Muta,

local scalar concentrations

Kazuki Kuga,

pp 691–701

Tomohiro
Kobayashi, andKazuhide Ito
Impact of heat generation

Ryota Muta andJapan Architecturaland use of experimental

Kazuhide Ito

Review, 2022, Vol.5, no.4,

instruments in a fume2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

pp 702–713

hood on pollutant capture
efficiency
Validation, verification,

Sung-Jun Yoo, andJapan Architecturaland quality control of

Kazuhide Ito

Review, 2022, Vol.5, no.4,

computational fluid

pp 714–727

dynamics analysis forindoor environments using
a computer-simulated
person with respiratory
tract
Determination of age ofKen Bryan

Building andair, purging flow rate, andFernandez, Naoki

Environment, 226

net escape velocity within

Ikegaya, Qingyan(2022), 109759

a cross ventilation model

Chen, andsheltered by urban-like

Kazuhide Itoblock arrays using largeeddy simulations
Coupled analysis of

Qianwen Zhang,Physics of Fluids, vol. 34,

transient aerodynamic

Chuqqi Su, Makoto

No. 1, pp.212013, DOI:

and dynamic response of

Tsubokura, Zhen

10.1063/5.0077988

articulated hearvy

Hu, Yiping Wangvehicles under crosswind
Numerical analysis of the

Keiji Onishi,

Physics of Fluids, vol.34,

efficiency of face masks for

Akiyoshi Iida,

No. 3, pp. 033309, DOI:

preventing droplet

Masashi10.1063/5.0083250

airborne infections

Yamakawa,
Makoto Tsubokura

Reference map 法を用いた

嶋田宗将, 西口浩司,

日本計算工学会論文集,

マーカー粒子によるオイラ

Peco Christian, 岡No. 20220002, DOI:

ー型流体-構造連成解析

澤重信, 坪倉誠

10.11421/jsces.2022.2022- 82 -

0002 (2022)
2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

2022 年

国外

ファルマシア（日本薬学202 年 5

国内

室内環境におけるウイルス

会機関紙）, vol.58, No.5,

月

飛沫・エアロゾル感染リス

pp.435-439, DOI:

Adaptively switched time

HsuehJui Lu,

AIP Advances, vol.12,

stepping scheme for direct

ChugGang Li,

035340, DOI:

acoustics computations,

Makoto Tsubokura

10.1063/5.0076657 (2022)

Numerical study ofDeboprasad

International journal ofnatural-convection from

Talkdar, Makoto

heat and mass transfer,horizontal cylinder at

Tsubokuravol. 58, No. 2, pp.849-

eccentric positions with

871, DOI:

charge in aspect ratio of a

10.1007/s00231-021-

cooled square enclosure03145-3 (2022)

Digital Transformation of

Kazuto Ando,

International Journal of

Droplet/Aerosol Infection

Rahul Bale,High Performance

Risk Assessment Realized

ChungGang Li,

Computing Applications,on “Fukuga” for the Fight

Satoshi Matsuoka,

vol.36, No.5-6, pp.568-

against COVID-19, (2022),

Keiji Onishi,

586, DOI:

Makoto Tsubokura

10.1177/10943202211160(alphabetical order56 (2022)

by surname)
Qantifying the COVID 19

Rahul Bale,

Scientific Reports,

infection risk due to

Akiyoshi Iida,

Vol.12, 11186, DOI:

droplet/aerosol inhalation

Masashi10.1038/s41598-022-

Yamakawa,

14862-y (2022)

ChungGang Li,
Makoto Tsubokura
Numerical Simulation of

Eriko Nambu,

Scientific Reports, vol.12,

Air Age in Dental Offices,

Kazunori Nozaki,

14120, DOI:

Makoto Tsubokura,

10.1038/s41598-022-

Mikako Hayashi

18588-9 (2022)

Numerical investigation of

H s u e h J u i L u,

Scientific Reports, vol.12,

effects of tongue

Tsukasa Yoshinaga,C

15361, DOI:

articulation andhungGang Li,

10.1038/s41598-022-

velopharyngeal closure on

Kazunori Nozaki,

18784-7 (2022).

the production of sibilant

Akiyoshi Iida,[s],

Makoto Tsubokura

スパコン「富岳」を用いた

坪倉誠- 83 -

クとリスク低減対策の定量

10.14894/faruawpsj.58.5_

評価

435

コロナ禍から生まれたイノ

石川文博, 坪倉誠,

ファルマシア（日本薬学202 年 5

ベーション

熊本卓哉

会機関紙）, vol.58, No.5,

月

国内

pp.401-404, DOI:
10.14894/faruawpsj.58.5_
401
ACM ゴードン・ベル賞と松岡聡, 坪倉誠

情報処理（情報処理学会

新型コロナ飛沫・エアロゾ誌）, Vol.63, No.7,

ル感染リスク評価のデジタ

pp.328-333 (2022)

2022 年

国内

2022 年

国内

2022 年

国内

2022 年

国内

2022 年

国内

2022 年

国内

2022 年

国内

ルトランスフォーメーション
複雑現象統一的解法フレー

坪倉誠，ラフールバ

ムワーク CUBE による自レ，李崇綱

特集：燃焼科学における大規模数値シミュレー

動車エンジンシミュレーシ

ションの最前線, 燃焼

ョン,

学会

スパコン「富岳」による飛粉体技術, vol.14, No.6,

坪倉誠

pp.119-25 (2022)

沫飛散シミュレーションとマスクの感染拡大防止効果
スパコン「富岳」による室

坪倉誠

病院設備（日本医療福祉

内環境における新型コロナ

設備協会）, vol.64,

ウイルス飛沫感染の予測とNo.4, (2022)

その対策
Numerical simulation of

HsuehJui Lu,

Twelfth International

the sibilant /s/ sound

Sakuya Sugimoto,Symp o s i u m o n

articulation process.

Tsukasa

Turbulence and Shear

Yoshinaga,

Flow Phenomena.

ChungGang Li,(TSFP12 2022)

Kazunori Nozaki,
Akiyoshi Iida, and
Makoto Tsubokura:
Infection risk in a

Rahul Bale,

5TH international

restaurant environment

Chung-Gang Li,

conference on building

due to airborne diseases.

Hajime Fukudome,

energy and environment.

Saori Yumino,(2022)

Akiyoshi Iida,
Makoto Tsubokura
Nonlinear mode

K a z u t o A n d o ,K

- 84 -P ro ce edin gs o fdecomposition andeiji Onishi, Rahul

International Conference

reduced-order modelingBale, Makoto

on High Performancefor three-dimensional

Tsubokura,Computing (ISC2021)

cylinder flow by

Akiyoshi Kuroda,(Jun 24-July 2, Online),

distributed learning on

and Kazuo Minami

Springer, Cham, pp.122–

Fugaku,6月

137 (2021)

Transient flow andAlicia Murga, Seigo

RANS versus LES

2022 年

particle deposition in the

Ohashi, Chung-comparison. Proceedings

6月

respiratory tract

Gang Li, Rahulof Indoor Air 2022,

Bale, Makoto

17th International

Tsubokura,

Conference of theInternational Society of
Indoor Air Quality &Climate.

- 85 -

国内