第 24 回数値流体力学シンポジウム
F6-1

産業機械における POD を用いた形状最適化手法
Shape Optimization Method for Industrial Machineries using POD°

米倉一男，IHI

解析技術部，横浜市磯子区新中原町 1 番地，E-mail：kazuo yonekura@ihi.co.jp

市東素明，IHI

解析技術部，横浜市磯子区新中原町 1 番地，E-mail：motoaki shito@ihi.co.jp

渡邊修，IHI

解析技術部，横浜市磯子区新中原町 1 番地，E-mail：osamu watanabe@ihi.co.jp

久保世志，IHI

解析技術部，横浜市磯子区新中原町 1 番地，E-mail：seiji kubo@ihi.co.jp

Kazuo Yonekura, IHI, 1, Shin-Nakahara-Cho, Isogo-ku, Yokohama
Motoaki Shito, IHI, 1, Shin-Nakahara-Cho, Isogo-ku, Yokohama
Osamu Watanabe, IHI, 1, Shin-Nakahara-Cho, Isogo-ku, Yokohama
Seiji Kubo, IHI, 1, Shin-Nakahara-Cho, Isogo-ku, YokohamaThis paper investigates on optimization method using POD for machineries associated with ﬂuids. We apply POD
for extracting the critical geometries of machineries. If we limit the design space using the critical geometries, wecan eﬃciently search high performance geometries. We conﬁrmed an advantage of the POD based method via
two numerical examples. The ﬁrst example is optimization of wing-section geometry. In this case, we can extract
critical geometries from NACA wing-sections and obtained high performance solutions. The second example iscompressor inlet valve. In this example, we compared the POD based method to genetic algorithm and conﬁrmedthat POD based method is superior to naive genetic algorithm in view of accuracy of optimal solution.

1.

としては主成分分析の一手法である固有直交分解 (POD:

はじめに翼や圧縮機などの流体が関係する産業機械において，その形状は製品の性能を大きく左右する．より高性能な機械を設計するために流体が関連する部分の形状最適化が広く行われている1 ．流体の支配方程式は非線形性が強いため発見的な最適化手法が多く用いられる．発見的手法の例としては例えば遺伝的アルゴリズム（GA: Ge-netic Algorithm）2 を用いた手法や応答曲面法（RSM:Response Surface Method）3 を用いた手法が挙げられる．他にも逐次二次計画法を用いた局所最適化も行われている4 ．GA は様々なアルゴリズムが提案され，多くの成果を挙げている．また RSM も実用の分野に浸透し，広く用いられる手法の一つになっている．しかしながら発見的な最適化手法を用いた場合一般に探索空間が広い場合には最適解の探索に時間がかかり，その精度も悪くなるという特徴がある．また最適化の過程がブラックボックスになっていて設計者に最適化の過程が見えにくいという欠点もある．そこで本稿ではこれらの欠点を克服するため高性能な解を含む探索空間をあらかじめ抽出してから最適化を行う手法を検討する．

Proper Orthogonal Decomposition) が有効である．PODは与えられた複数のベクトルを直交する基底ベクトルに分解する．特によく似た形状については高次の基底関数の影響が小さくなることが期待され，少ない基底関数から新しい形状を作成することが可能である．POD を用いた最適化手法として，Toal ら5 は POD を用いて形状をふるいにかけてから最適化を実行する手法を提案した．
本稿で検討する POD を用いた最適解の探索手法はFig. 1 のとおりである．まず事前にいくつかの性能の良い形状を得ておく．これは経験や知見に従って得るか，あるいはランダムサンプリングにより得ることができる．次にそれらの形状について POD を実行する．POD により得られた基底関数のうち低次の基底関数だけを用いてパラメトリックに新しい形状を生成する．そしてパラメータを変数として RSM を用いて最適解を探索する．PODを用いることでパラメータは例えば 3 個程度の少ない数に限定できるので，探索領域が大幅に限定される．また探索領域は高性能な形状から生成するため，高性能な形状を含んでいることが期待できる．従って RSM や GAを用いたときに，より良い解を探索できることが期待できる．さらに設計者が特徴的な形状を実際に見ることができ，最適化の過程がわかりやすいという利点もある．

高性能な形状にはいくつかの特徴的な形状が存在することが多い．その特徴的な形状を抽出して組み合わせ，別の形状を生成した場合，それらの中には高性能な形状が本稿では次の表記を用いる．すなわちベクトルは v な
どの太字で表記し，転置ベクトルを v > で表すこととす

多く含まれている可能性がある．従ってこのように特徴形状を利用すれば，効率的に探索できる可能性がある．例
その知見を利用して特徴形状を抽出することができればる．またベクトルのノルムはすべてユークリッドノルム√kvk = v > v とする．なお，扱う値はすべて実数である．

形状最適化を効率よく行うことができる．またこのとき

2.

設計者は高性能な形状に共通する形状を知ることができ，固有直交分解 (POD: Proper Orthogonal Decompoition) は主成分分析の一手法で，相関があることが期

えば二次元翼形状に関しては多くの知見が得られており，最適化の過程が把握しやすい．特徴形状を抽出する手法

1

固有直交分解c 2010 by JSFM
Copyright °第 24 回数値流体力学シンポジウム
F6-1と定義すると，aij ，uj は行列 V > V の固有値分解を用いPrepare high performance geometries.
(Via randam sampling, knowledge or else.)て求めることができる．式 (1) を用いると V > V は次の↓ように記述できる．

V >V

Calculate basis vectors via POD.

↓=nX

kuj k2 aj a>j .

j=1

Cut oﬀ high-order basis vectors.

ここで両辺に左から a1 > ，右から a1 をかけると，↓>Generate new geometries from basis vectors.

>a1 V V a1

↓

=nX

¡¢2
kuj k2 a1 > aj

j=1を得る．従って a1 を固定したとき，ku1 k を最大化するFind the optimal geometry byResponce Surface Method.

aj は a1 と直交する．さらにこのとき ku1 k を最大化す
る a1 は V > V の最大固有値に対応する正規固有ベクトルである．同様にして aｊ はすべて V > V の j 番目に大きい固有値に対応する正規固有ベクトルであることが従う．
最後に aｊ が直交しているから，
Xaij v i = ui

Fig. 1: Flow chart of the algorithm.

待される複数のベクトルを，直交する基底ベクトルに分解する．本稿では Sirovich 6 による snapshotPOD を用iいる．POD の概略は Fig. 2 のとおりである．
以下ではアルゴリズムを簡単に述べる．主成分分析にが従う．以上のように行列の固有値分解を利用して POD関してはラグランジュの未定乗数法を用いた説明がなさがを行うことができる．本稿の数値実験では行列のサイれることが多いが，本稿ではアルゴリズムと関連付けるズが比較的小さいこともあり，固有値分解の方法としてために線形代数を用いて説明する．今回はあらかじめ n実装が簡単なべき乗法を用いた．

個の形状の幾何情報がベクトル v̂ i ∈ Rm (i = 1, . . . , n)

得られた直交分解から新しい形状を生成することがでの形で与えられているものとする．このときベクトル v̂ i

きる．基底関数 uj は n 個得られるが，j が大きくなるにの j 番目の値を (v̂i )j で表すと約束すると，(v̂i )j は各形つれてそのノルムは小さくなる．特に最初に与えた形状状に共通する x 座標 xj に対する y 座標値である．相関を求めるために，まず与えられたベクトル v̂ i を平均成分の個数が多く，さらに相関が大きければ高次の基底関数のノルムは無視できるほど小さくなることが期待される．とそれ以外の成分に次のように分解する．

そこで例えば基底ベクトルとして，そのノルムが閾値以上であるノルムを採用する．すなわち例えば R = 0.01 にv̂ i = v̄ + v i ,v̄ =n1X

n i=1

対して

½v̂ i .

N = maxj

固有直交分解はベクトル v i を直交基底 uj の線形和に分¯¾¯ kuj k¯≥Rj ¯
ku1 k(3)などとして j ≤ N であるような基底関数を用いて新しい解する．

形状を生成する．このとき新しい形状を表すベクトル pvi =nX

はaij uj .

(1)

j=1

た だ し uj は 直 交 基 底 で あ る ．ま た ，aj³´>とすると，
a1j a2j · · · anj

kaj k = 1 (∀j ∈ {1, . . . , n})p = v̄ +=NX

αj uj

(4)

j=1で表される．ただし α ∈ RN は任意に設定するパラ メ ー タ で あ り，与 え 方 の 自 由 度 は 大 き い ．今 回 は(2)max {|αi | | i = 1, . . . , N } ≤ 2.5 の範囲で適当な格子点を選択した．ただし制約条件が課されている場合にはその制約条件を満たす範囲で α を選択する．
式 (4) により超平面が与えられる．これに拘束条件を考慮に入れることで超平面上の部分空間が設計空間となる．この設計空間は最初に与えられた形状ベクトルで張を満たす．POD の目的は条件 (1)，(2) を満たし，添え字 j が小さい順に kuj k が最大となるような分解を求める．uj と aij を求めることである．ここで行列 V を³´
V = v1 v2 · · · vn ,2c 2010 by JSFM
Copyright °第 24 回数値流体力学シンポジウム
F6-1
1POD-design-space
original points
feasible reagion

Given: Vectors v̂ i (i = 1, 2, . . . , n).
Find: Orthogonal basis vectors uj
Pnwhere v̂ i = v̄ + j=1 aij uj .
STEP0. Extract averagevector; v̂ i =´ v̄ + v i .
³

STEP1. Set V = v 1 v 2 · · · v n .
STEP2. Calculate normal eigen vectors ai of V > V .
STEP3. Calculate basis vectors ui
Pby uj = i aij v i .

0.5

0Fig. 2: Abstract of POD.

00.5

1Fig. 3: POD disign space
られる空間の部分空間でもある．この設計空間と最初に与える形状との関係は Fig. 3 のようになる．ここでは簡単のために 2 次元平面上に描いているが，実際には次元はもっと大きい．例えば第 3 節で扱う翼形状最適化の例

目的関数

CL → 最大化,
CD → 最小化,

制約条件

コード長 = 1,

題においては設計変数は 200 個あり自由度は 200 自由度
であるが，POD により 3 次元超平面を抽出している．

翼厚比 = r.

式 (3) における R を適切に選ぶことで POD により得られる設計空間は Fig. 3 のように最初に与える幾何情報の点集合の十分近くを通る平面にできる．仮に元の幾何情報の間の相関が小さい場合であっても R = 0 とするこ

一般に翼厚比が大きいほど CL が大きく CD は小さくとで式 (4) の超平面が元の幾何情報ベクトルをすべて含なる傾向があることが知られている．POD を用いた場むようにできる．ここで R を変化させると N が変化し，合，式 (4) においてパラメータ α を大きくすると翼厚比POD により得られる設計空間の次元が変化することになる．

が大きくなり，それだけで揚力係数が大きくなることがありえる．すなわち，新規に高性能な形状を生成できたとしても，それが POD の効果であるのか翼厚比を大き

第 1 節で述べたようにこの超平面上の目的関数値の集くしたことによる効果であるのか判断が困難である．そ合は比較的良い値を含んでいる可能性が高い．したがっこで本節ではすべての翼について翼厚比が一定の値になて POD を用いない場合には Fig. 3 の feasible reagion のるよう縮尺を変える．

ような設計空間全体を探索しなければならないのに対して，POD を用いた場合には Fig. 3 の POD-disign-spaceを探索すれば十分である可能性が高い．すなわち PODにより高性能な解の集合を抽出しようとしていると考えればよい．特に目的関数が連続でかつ，R を正として打ち切った基底ベクトルの影響が無視できるならば，超平面上の目的関数値は，少なくとも元形状に相当する目的関数値の最大値と最小値に挟まれた区間を含む．このとき，この超平面上を探索すれば元形状と同じかそれより高性能な解が求まる．また設計空間が狭められたことで応答曲面を生成する場合に曲面の精度が向上して，最適解の精度が向上し探索時間が短縮することが期待できる．

3.

２次元翼の形状最適化
本節では実際に二次元翼の形状最適化を行う．翼形状

設計において重要な因子は揚力係数 CL と抗力係数 CDである．そこで次のような多目的最適化問題を考える．

3

Abbot7 を 参 考 に 揚 力 係 数 が 比 較 的 大 き く 抗 力係 数 が 比 較 的 小 さ い 翼 形 状 と し て ，NACA2415，
NACA4415，NACA23012，NACA631 -412，NACA641 412，NACA651 -212 を元形状として選択する．ただし翼厚比 r はこれらの翼の平均翼厚比である 0.129730 とした．元形状と，それぞれの形状に対する CL および CDの値を Fig. 4 に示した．x 軸を 100 分割し，翼正圧側と負圧側のそれぞれに対して各 x 座標に対応する y 座標を並べたベクトルを最後につなげて幾何情報ベクトルとする．このベクトルを入力として POD を実行し固有ベクトルと基底ベクトル（Fig. 5）を得た．基底ベクトルのノルムは Fig. 6 に示したように 4 次までの基底ベクトルが卓越している．今回はパラメータ数を少なくするため3 次までの基底ベクトルを用いる．いくつかの係数ベク
トル α ∈ R3 に対して式 (4) に従い新規形状を 108 形状生成する．パラメータ α は第 2 節で述べたとおりに設定
c 2010 by JSFM
Copyright °第 24 回数値流体力学シンポジウム
F6-1

した．具体的には

α1 ∈ {1.0, 1.5, 2.0, 2.5},0.080.04
0-0.0400.20.40.60.8

α2 ∈ {−2.5, −2.0, −1.5, −1.0, 0.5, 0.0, 0.5, 1.0, 1.5, 2.0},

1

α3 ∈ {−1.5, −1.0, 0.0, 1.0},(a) NACA2415, CL = 0.7022, CD = 0.00731

からそれぞれ選択した．元の形状と新規形状について CL

0.080.040-0.04

および CD を計算した結果を Fig. 7 に示す．CL ，CD00.2

0.4

0.60.8の分布は元形状の値の周囲に広範囲に分布しており，元1形状よりも高性能な形状が多く存在していることがわか(b) NACA4415, CL = 0.9182, CD = 0.00919

る．従ってこの例では設計空間をうまく限定することが0.080.040-0.04

できたといえる．ただし CL および CD の計算にはパネル法を実装したソルバを用いた．計算条件はレイノルズ
00.20.40.60.8

数 Re = 3.0 × 106 ，迎角 5 度とした．

1

新規形状の情報を元に応答曲面を生成し最適解を探索(c) NACA23012, CL = 0.9011, CD = 0.00795

する．応答曲面の生成には ModeFRONTIER ver.4.2.1 80.080.040-0.04を用いる．最終的に CL /CD が最大となる翼を最適解として選択した．選択解の形状を Fig. 8 (a) に示す．比較のた
00.20.40.60.8

1め勝井と冨田4 より，揚抗比の最大化を目的とした形状最(d) NACA631 -412, CL = 0.9681, CD = 0.0077

適化の結果を，本稿に合わせて修正した結果を Fig. 8(b)0.080.040-0.04に示す．勝井と冨田4 はベジェ曲線を用いて翼を表現し，00.2

0.4

0.60.8NACA0012 翼を元にして逐次二次計画法による形状最適化を行った．引用するにあたって，POD による結果との比較を行うために翼厚比を本稿の条件に合わせて変え，パネル法を用いて再計算した．そのため引用元の論文とは形状と揚力係数，抗力係数，揚坑比が異なっている．レイノルズ数と迎角の計算条件は引用元においても本稿の条件と同じである．揚力係数，抗力係数，揚抗比の比較を Tab. 1 に示す．最適解の幾何形状は大きく異なっている．また揚力係数はどちらもほぼ同じであるのに対して抗力係数は POD の方が小さい．さらに POD を持ちいた手法の方が揚抗比が大きく，高性能な解を得ることができていると言える．
最後に最適解と基底ベクトルの関係について述べる．最適解を生成する係数パラメータは α = (1.5, −2, −1)> であった．元形状として採用した翼はすべてキャンバーが比較的小さいが，選択解はキャンバーが大きいという特徴がある．また翼正圧側後縁が凹な形状をしている．ここで 1 次の基底ベクトルを見るとキャンバーが大きい．
また 2 次の基底ベクトルは正圧側後縁で凸になっている．
POD を用いることでキャンバーを持つ成分と正圧側後縁が凸な成分を独立に抜き出し，これにパラメータを乗じることでキャンバーが大きく後縁側で凹な形状をした選択解を生成することができたことがわかる．また平板翼よりもキャンバーが大きく後縁側が凹な翼の方が揚力係数が大きくなることが知られており，選択解のような翼を選択することは妥当である．このように設計者は基底ベクトルを確認することで選択解の妥当性を確認することができる．

1(e) NACA641 -412, CL = 0.8755, CD = 0.008
0.080.04
0-0.0400.20.40.60.8

1(f) NACA651 -212, CL = 0.6989, CD = 0.00747
Fig. 4: Original geometries.

0.1
ave.
mode1
mode2
mode3
mode4

0.08
0.06
0.04
0.02
0-0.02
-0.04
-0.06
00.2

0.4

0.6

0.8

1

Fig. 5: Average vector and 1–4th basis vectors.

4c 2010 by JSFM
Copyright °第 24 回数値流体力学シンポジウム
F6-1

0.6

0.12
0.08
0.04
0-0.04

0.5

00.2

0.4

0.6

0.8

1

0.4
Power

(a) Selected optimal solution.
0.3
0.12
0.08
0.04
0-0.04

0.2
0.1

00.2

0.4

0.6

0.8

1(b) Katsui and Tomita wing-section (adjusted).

01234
Mode56

Fig. 8: Geometries of optimal solutions.

Fig. 6: Wing section; Norms of basis vectors.

Tab. 1: Comparison of optimal solutions
Generated
Original

0.014

CD

0.012

POD

Katsui and Tomita (adjusted)

CL

1.28591.2944CD

0.007070.01302

CL /CD

181.90

99.42

0.01

した．具体的には

0.008

α1 , α2 ∈ {0.1, 0.8, 1.5, 2.0}

0.006
0.5

0.6

0.7

0.8

0.9
1CL

1.1

1.2

1.3

α3 ∈ {0.1, 0.8, 1.5}

1.4

からそれぞれ選択した．次に 48 個の設計変数のうち，制Fig. 7: Plots of CL and CD .

約条件を満たさないものを取り除き，残ったものを元に形状をつくり，Fluent を用いて全圧損失の値を求める．
得られた全圧損失の値とパラメータ α の情報を元に応答4.

曲面を生成し，応答曲面のみを用いて最適解を探索する．

吸吐弁の形状最適化

探索された最適解は応答曲面の生成点であったため，再前節で扱った例題は形状が座標値で与えられている場計算は行う必要がなかった．なお，応答曲面を生成する合であった．しかし実際の設計においては部材の肉厚など際は ModeFrontier ver. 4.2.1 8 を用いた．Tab. 2 に最適の決められた部分の寸法を最適化したい場合も多い．そ解の全圧損失係数の値を示す．POD を用いた手法の最適のような場合を想定してレシプロ圧縮機吸吐弁の形状最解と併せて，POD を用いずに GA を用いて最適化を行っ適化を行う．圧縮機吸吐弁は圧縮機への吸気を制御するた結果も併せて示す．ζ が全圧損失係数である．GA を実ための弁である．弁が開いているときの形状と 10 個の行したときの条件を Tab. 3 に示す．GA を用いた場合は設計変数を Fig. 9 に示す．各設計変数には制約条件が課されている．流路が合流する部分でバックステップ流れ全圧損失係数が ζ = 27.3 である一方で，POD を用いた場合は全圧損失係数が ζ = 25.6 となっている．POD をがあり全圧損失が大きい．そこで流路形状に関して最適用いることでより高性能な解が探索できたことが確認で化を行い全圧損失を小さくすることを目標とする．本節では設計変数の値を入力ベクトルとして POD を実行すきる．以上のように限られた設計変数を元に形状が決まる．いま，タグチメソッドを用いて圧力損失が比較的小性能な解を求めることができる．

さい形状が 5 形状与えられている．これらに対して POD

5.

る場合に関しても POD を用いた最適化手法を用いて高結論を実行して基底ベクトルを得る．基底ベクトルのノルム

流体関連機械の形状最適化における POD を用いた高は Fig. 10 に示したように 3 次までの基底ベクトルが卓性能な解の探索手法を検討した．検討した手法では，まず越している．そこで 3 次までの基底ベクトルを用いて式

流体性能の良い構造物の幾何情報の集合を元にして POD(4) に従って新たにパラメトリックに 48 個の設計変数を生成する．パラメータ α は第 2 節で述べたとおりに設定によりそれらを基底ベクトルに分解する．次に低次の基底ベクトルからパラメトリックに新規形状を表現し，RSM

5c 2010 by JSFM
Copyright °第 24 回数値流体力学シンポジウム
F6-1

Tab. 2: Pressure loss coeﬃcient.
POD + RSM GA + RSM (without POD)
ζ25.6

27.3

Tab. 3: Property of GA.

Fig. 9: Geometry model and parameters.

36

世代数20

交叉率0.5

突然変異確率0.1

RSM 使用率0.5

総計算数7200.7

参考文献

0.6

1. 吉清水宗，下山幸治，鄭信圭，大林茂，横野泰之: 最適解探索能力向上に向けた翼形状表現方法に関する検討, 第 8 回最適化シンポジウム 2008 講演論文集, 4,175–183 (2006).

0.5Power

初期個体数

0.4
0.3
0.2

2. Deb, K.: Multi-objective optimization using evolutionary algorithm, John Wiley & Sons, New York,
2001.

0.1
01

23
Mode

45

3. Myers, R. H. and Montgomery, D. C.: Responce Surface Methodology: process and product optimization using designed experiments, second edition, A
Willey-Interscience publication, New York, 2002.

Fig. 10: Valve: Norms of basis vectors.などの非線形最適化手法を用いて最適解を探索する．この手法を用いることで性能の良い形状に共通する特徴を持った形状が得られる．またパラメータを減らせるため，最適化を行う際の探索領域が狭くなり，発見的手法を用いた場合に，より良い解が求められることが期待できる．
実際に二次元翼形状の最適化問題と，吸吐弁の形状最適化問題について最適化を行った．どちらの場合も性能の良い解を得ることができた．さらに低次の主成分を分析することで，高性能な形状の特徴を知ることができる．このためブラックボックス化しがちな最適設計の過程が設計者にわかりやすく，ノウハウの蓄積にも繋がる．また逆に，これまでに得られている高性能な形状を用いて新たな形状を生成することで知識の有効利用が可能である．

4. 勝井辰博，冨田高嗣: CAD/CFD 統合型高効率最適化システムの構築とその検証 –2 次元翼型の最適化による最適化システムの検証–, 日本船舶海洋工学会論文集, 4, 175–183 (2006).
5. Toal, D. J. J., Bressloﬀ, N. W. and Keane, A. J.:
Geometric Filtration Using POD for Aerodynamic
Design Optimization, Collection of Technical Papers- 26th AIAA Applied Aerodynamics Conference, 2,859–871 (2008).
6. Sirovich, L.: Turbulence and dynamics of coherent
structures part 1: Coherent structures, Quarterly of
applied mathematics, 45, (3), 561–571 (1987).
7. Abbot, I. H. and von Doenhoﬀ, A. E.: Theory ofwing sections including a summary of airfoil data,
Dover publications, INC., New York, 1959.

POD を用いることの最大の利点は探索領域が小さくなることである．例えば第 2 節の翼形状の最適化の問題
では設計変数を 3 個のパラメータに減少させている．この利点は特に設計変数の多い大規模な問題において効果が大きい．今後の課題はより大規模な問題に対して手法の有効性を確認することである．

8. 多目的ロバスト設計最適化支援ツール modeFRONTIER, http://www.cdaj.co.jp/product/020000mod
efrontier/index.html.
6c 2010 by JSFM
Copyright °