特許 6167554 流路形状最適化方法及び流路形状最適化装置

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】6167554

(24)【登録日】2017年7月7日

(45)【発行日】2017年7月26日

(54)【発明の名称】流路形状最適化方法及び流路形状最適化装置

(51)【国際特許分類】

   G06F 17/50 20060101AFI20170713BHJP

   G06F 19/00 20110101ALI20170713BHJP

【ＦＩ】

   G06F17/50 604A

   G06F17/50 612G

   G06F19/00 110

【請求項の数】4

【全頁数】16

(21)【出願番号】特願2013-31184(P2013-31184)

(22)【出願日】2013年2月20日

(65)【公開番号】特開2014-160398(P2014-160398A)

(43)【公開日】2014年9月4日

【審査請求日】2015年12月21日

【新規性喪失の例外の表示】特許法第３０条第２項適用 ２０１２年８月２１日 ウェブサイト ｈｔｔｐ：／／ｊｓｉａｍ２０１２．ｋａｓｈｉ．ｉｎｆｏ．ｗａｓｅｄａ．ａｃ．ｊｐ／〜ｋａｓｈｉ／ｄｏｗｎｌｏａｄ／で公開された日本応用数理学会２０１２年度年会講演予稿集に公開

(73)【特許権者】

【識別番号】000000099

【氏名又は名称】株式会社ＩＨＩ

(74)【代理人】

【識別番号】100175802

【弁理士】

【氏名又は名称】寺本 光生

(74)【代理人】

【識別番号】100064908

【弁理士】

【氏名又は名称】志賀 正武

(74)【代理人】

【識別番号】100167553

【弁理士】

【氏名又は名称】高橋 久典

(72)【発明者】

【氏名】米倉 一男

【審査官】 松浦 功

(56)【参考文献】

【文献】 特開２０１０−１６５２４０（ＪＰ，Ａ）

【文献】 特開２００４−１０２９６２（ＪＰ，Ａ）

【文献】 特開２０１２−１９８６１０（ＪＰ，Ａ）

【文献】 特開２０１３−０３７４３４（ＪＰ，Ａ）

(58)【調査した分野】（Int.Cl.，ＤＢ名）

Ｇ０６Ｆ １７／５０

Ｇ０６Ｆ １９／００

Ｇ０６Ｑ ５０／０４

Ｇｏｏｇｌｅ Ｓｃｈｏｌａｒ

(57)【特許請求の範囲】

【請求項1】

コンピュータが、
流体解析により流れ場を求める流れ場計算工程と、
流れ場計算工程で求められた流れ場の形状を変更し、当該変更による前記流れ場のエネルギー散逸を示す目的関数に対する感度を計算する感度計算工程と、
前記感度に基づいて流れ場の形状を更新する形状更新工程と、
を実行し、
前記コンピュータが、前記流体解析のタイムステップ数を制御して非定常状態で前記流れ場の形状を更新しながら、前記流れ場計算工程、前記感度計算工程及び前記形状更新工程を繰り返し、前記流れ場の形状が更新されなくなる状態となるまで収束計算する
ことを特徴とする流路形状最適化方法。

【請求項2】

前記コンピュータは、前記流体解析に格子ボルツマン法を用いていることを特徴とする請求項１記載の流路形状最適化方法。

【請求項3】

前記コンピュータは、各位置における空隙率を設定可能な多孔質媒体モデルを用い、当該多孔質媒体モデルの内部に流れを形成し、前記感度に基づいて前記空隙率を更新することによって流れ場の形状を更新することを特徴とする請求項１または２記載の流路形状最適化方法。

【請求項4】

流体解析により流れ場を求める流れ場計算工程と、
流れ場計算工程で求められた流れ場の形状を変更し、当該変更による前記流れ場のエネルギー散逸を示す目的関数に対する感度を計算する感度計算工程と、
前記感度に基づいて流れ場の形状を更新する形状更新工程と、
を実行する演算処理手段を有し、
前記演算処理手段が、前記流体解析のタイムステップ数を制御して非定常状態で前記流れ場の形状を更新しながら、前記流れ場計算工程、前記感度計算工程及び前記形状更新工程を繰り返し、前記流れ場の形状が更新されなくなる状態となるまで収束計算する
ことを特徴とする流路形状最適化装置。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、流路形状最適化方法及び流路形状最適化装置に関するものである。

【背景技術】

【0002】

従来、流体解析による流路形状の最適化は、初期設定された流れ場において定常状態となるまで収束計算を行い、その後、感度計算に基づいて流れ場の形状を変化させ、新たな流れ場において再び定常状態となるまで収束計算を行うことで最適解を探索している（非特許文献１）。また、例えば、非特許文献２においては、流体解析の１つである格子ボルツマン法を用いた構造の最適化について開示されている。

【先行技術文献】

【非特許文献】

【0003】

【非特許文献1】日本機械学会論文集（Ｂ編） ７６巻７６８号（２０１０−８） ｐ．３２−４１ 論文Ｎｏ．１０−０１１２ 「低レイノルズ数流れにおける最適物体形状の探索（トポロジー最適化法による抗力最小化と揚力最大化）」 近藤継男 松森唯益 川本敦史

【非特許文献2】日本機械学会 第２２回設計工学・システム部門講演会・２０１２．９．２６−２８ ｐ．４８０−４８７ 「格子ボルツマン法を用いた構造最適化における一考察」 矢地謙太郎 山田祟恭 吉野正人 松本敏郎 泉井一浩 西脇眞二

【発明の概要】

【発明が解決しようとする課題】

【0004】

しかしながら、従来の流体解析による流路形状の最適化では、流れ場の形状を変化させるごとに流体解析の計算を行うことから計算量が膨大であり、最適解を見つけ出すために非常に多くの時間が必要となっている。このため、１回の収束計算に特に多くの時間が必要となるような大規模な乱流を含むような流れ場の流路形状の最適化には適用することは現実的ではない。

【0005】

本発明は、上述する問題点に鑑みてなされたもので、流体解析による流路形状の最適化に要する時間を短縮することを目的とする。

【課題を解決するための手段】

【0006】

本発明は、上記課題を解決するための手段として、以下の構成を採用する。

【0007】

第１の発明の流路形状最適化方法であって、流体解析により流れ場を求める流れ場計算工程と、流れ場計算工程で求められた流れ場の形状を変更し、当該変更による感度を計算する感度計算工程と、上記感度に基づいて流れ場の形状を更新する形状更新工程と、を有し、上記流れ場が非定常状態の間に上記流れ場計算工程、上記感度計算工程及び上記形状更新工程を繰り返しながら、上記流れ場の形状が更新されなくなる状態となるまで収束計算するという構成を採用する。

【0008】

第２の発明は、上記第１の発明において、上記流体解析が格子ボルツマン法を用いているという構成を採用する。

【0009】

第３の発明は、上記第１または第２の発明において、各位置における空隙率を設定可能な多孔質媒体モデルを用い、当該多孔質媒体モデルの内部に流れを形成し、上記感度に基づいて上記空隙率を更新することによって流れ場の形状を更新するという構成を採用する。

【0010】

第４の発明は、流路形状最適化装置であって、流れ場が非定常状態の間に、流れ場の計算と、流れ場の形状を変更した場合の感度計算と、感度に基づく流れ場の形状の更新とを繰り返し、さらに上記流れ場が定常状態となるまで収束計算する演算処理手段を備えるという構成を採用する。

【発明の効果】

【0011】

本発明によれば、流れ場が非定常状態のうちに流れ場の計算と、流れ場の形状を変更した場合の感度計算と、感度に基づく流れ場の形状の更新とを繰り返す。このため、最終的に流れ場の形状が更新されなくなる状態となったときには、最適な流れ場の形状（すなわち流路形状）を得ることができる。このような本発明によれば、１度の収束計算によって流路の最適形状を得られるため、流体解析による流路形状の最適化に要する時間を極めて短縮することができる。

【図面の簡単な説明】

【0012】

【図1】本発明の考え方を説明するためのフローチャートである。

【図2】本発明の一実施形態における流路形状最適化装置のハードウェア構成を概略的に示すブロック図である。

【図3】（ａ）が本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第１実施例で用いた多孔質媒体モデルの模式図であり、（ｂ）が第１実施例におけるタイムステップと流体損失との関係を示すグラフである。

【図4】本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第１実施例における各タイムステップにおける流路形状を可視化した図である。

【図5】（ａ）が本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第２実施例で用いた多孔質媒体モデルの模式図であり、（ｂ）が第２実施例におけるタイムステップと流体損失との関係を示すグラフである。

【図6】本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第２実施例における各タイムステップにおける流路形状を可視化した図である。

【図7】（ａ）が本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第３実施例で用いた多孔質媒体モデルの模式図であり、（ｂ）が第３実施例におけるタイムステップと流体損失との関係を示すグラフである。

【図8】本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第３実施例における各タイムステップにおける流路形状を可視化した図である。

【図9】（ａ）が本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第４実施例で用いた多孔質媒体モデルの模式図であり、（ｂ）が第４実施例におけるタイムステップと流体損失との関係を示すグラフである。

【図10】本発明の一実施形態における流路形状最適化装置及び流路形状最適化方法の第４実施例における各タイムステップにおける流路形状を可視化した図である。

【発明を実施するための形態】

【0013】

本発明に係る流体最適化方法及び流体最適化装置の一実施形態についての説明に先立ち、本発明の考え方について説明する。

【0014】

まず、基礎方程式の定式化について説明する。本発明は、流路の入口から供給されて出口に至るまでの流体のエネルギ損失が最小限となる流路形状（最適な流路形状）を短時間で得ることを目的とするものである。本発明では、最適な流路形状を得るために、流体解析の手法として、格子ボルツマン法を用いる。格子ボルツマン法は、ナビエ−ストークス方程式の近似工程式としてボルツマン方程式を格子上で離散化して解く手法である。このような格子ボルツマン法は、計算の過程で収束計算を必要とせず、単純な四則演算のみで数値計算を行うため、ナビエ−ストークス方程式を用いた計算と比較して収束性が良くかつ計算時間も短い。

【0015】

また、本発明では、最適化の過程で流体と固体とを区別するために多孔質媒体モデルを用い、設計変数は多孔質媒体モデルの各位置の空隙率とする。また、本発明では、格子ボルツマン法による非定常計算を１タイムステップ計算するごとに感度解析（感度の計算）を行い、この感度解析に基づいて空隙率を更新する。このため、比較的短時間で最適化計算を終えることができる。なお、非定常計算中に空隙率を更新すると一般的な流体解析では収束性が悪くなることが予測されるが、格子ボルツマン法を用いることで安定して計算を行うことができる。

【0016】

多孔質媒体モデルとしては、例えばブリンクマンモデルを使用する。この場合の流体の支配方程式は、［数１］に示す通りである。なお、［数１］において、εは空隙率、Ｆは流体に働く外力、ｐは圧力、ｕは流体の速度ベクトル、νは粘性係数、αは自由に設定できる定数を示している。

【0017】

【数1】

【0018】

これを、格子ボルツマン法を用いて解く場合、［数２］の式にしただって計算すれば良い。なお、［数２］において、ｆｉはｉ方向の速度分布関数、ｆ^ｅｑ_ｉは平衡分布関数、ｆ^ｃｏｌ_ｉは便宜上使用している変数、Ｆ_ｉはｉ方向の外力項、δ_ｔは時間ステップ幅、ｅ_ｉはｉ方向の格子速度ベクトル、ｗ_ｉは重み係数、τは時間緩和係数、ｃは格子速度（定数）を示している。

【0019】

【数2】

【0020】

なお、多孔質媒体モデルを用いた格子ボルツマン法については、下記文献に詳しく記載されており、ここでの詳細な説明は省略する。

【0021】

Z. Guo and T. S. Zhao, Lattice Boltzmann model for incompressible flows through porous media, PHISICAL REVIEW E, 66 (2002), 036304.

【0022】

続いて、最適化問題の定式化について説明する。本発明では、流れにおけるエネルギ散逸を目的関数として最適化を行う。ここでの最適化問題は、［数３］〜［数６］の式に示すようになる。なお、［数３］に示す目的関数は、［数４］の式及び［数５］の式によって定義される。ここで、［数４］の式すなわち［数３］の式の第一項は粘性項に起因する粘性散逸を示し、［数５］の式すなわち［数３］の式の第二項は多孔質媒体の抵抗に起因する摩擦損失を示す。なお、［数３］において、Ｄは粘性散逸関数、Ｃは摩擦損失関数、Ｔは計算終了時刻を示す。また、［数４］において、ｐは自由に設定できる定数（例えばｐ＝２〜５程度、Ωは設計領域を示す。また、［数６］において、Ｖは設計者が指定する規定体積を示す。

【0023】

【数3】

【0024】

【数4】

【0025】

【数5】

【0026】

【数6】

【0027】

ここで、最適化問題におけるフローチャートを図１に示す。図１に示すように、本発明では、初期化の工程と（ステップＳ１、ステップＳ２及びステップＳ３）と、最適化の工程（ステップＳ４、ステップＳ５、ステップＳ６、ステップＳ７及びステップＳ８）とに大別される。

【0028】

まずステップＳ１においては、初期条件の設定を行う。ここでは、流体の流量や流速の設定、流体の物性値（粘性等）の設定、多孔質媒体モデルの定義（例えば流体の入口位置や出口位置の設定）等を行う。続いて、ステップＳ２では、ステップＳ１で設定した初期条件に基づいて格子ボルツマン法を用いて計算を行う。ステップＳ３では、ステップＳ２による計算が収束したか否かを判定し、収束していない場合には引き続きステップＳ２を行う。なお、ステップＳ３では、ステップＳ２の計算結果が周期的な流れを示す状態である場合においても収束したと判定する。これらのステップＳ１、ステップＳ２及びステップＳ３とで構成される初期化工程により、多孔質媒体モデルに入口から出口を繋ぐ流路が形成される。

【0029】

なお、本発明において初期化工程は、必ずしも必須の工程ではない。例えば、予め流路が形成された多孔質媒体モデルを利用することによって初期化工程を省略しても良い。また、初期化工程では、格子ボルツマン法以外の手法によって収束計算を行っても良い。

【0030】

上述のステップＳ３で計算が収束したと判定された場合には、最適化工程に移行する。最適化工程では、初期工程によって流路が形成された多孔質媒体モデルに対して流体が供給されるところから計算が開始される。まず、ステップＳ４において格子ボルツマン法で非定常流れを１タイムステップのみ計算を行う（流れ場計算工程）。すなわち、ステップＳ４では、Ｔ＝１において、［数４］〜［数６］の式に基づいて上述の［数３］の式を解くことにより、流体解析により流れ場を求める。

【0031】

続いてステップＳ５では感度を計算する（感度計算工程）。ここでは、流れ場の形状を変更（すなわち空隙率を仮に変更）し、［数７］の式を用いて感度を計算する。なお、［数７］の式のうち「数８」で示す部分は感度を示し、［数９］で示す係数は制約条件を満たすための修正係数であり、［数１０］で示す部分は制約条件を満たすための修正項である。

【0032】

【数7】

【0033】

【数8】

【0034】

【数9】

【0035】

【数10】

【0036】

なお、［数７］を解くにあたり、［数１１］で示す式を用いる。

【0037】

【数11】

【0038】

また、本発明のように格子ボルツマン法を用いる場合には［数１２］となる。なお、［数１２］において、ｆ^ｅｑは平衡分布関数、ρは密度、ｕ（太字）は速度ベクトル、ｅ（太字）は格子点での粒子速度、ｃはｃ＝ｄｘ／ｄｔで、ｄｘは格子間距離、ｄｔは時間刻み幅である。αは多孔質媒体モデルでの抵抗に相当し、εの関数である。τは緩和時間係数と呼ばれ、レイノルズ数に依存する定数である。また、βは便宜上使用する変数である。

【0039】

【数12】

【0040】

なお、［数１１］及び［数１２］には、［数１３］で示す式が含まれている。

【0041】

【数13】

【0042】

この［数１３］で示す式は、ナビエ−ストークスの式を用いた場合には陽に解けないので随伴変数法が必要になるのに対して、格子ボルツマン法を用いた場合には陽に解けるので、随伴変数法は不要となる。

【0043】

続いて、ステップＳ６では空隙率の更新を行う（形状更新工程）。ここでは、［数１４］を用いて、εをステップＳ５で求めたｄεだけ更新する。なお、ステップＳ６では、ステップＳ５で感度計算した結果、最も感度の高いものに基づいて更新を行う。

【0044】

【数14】

【0045】

続いて、ステップＳ７では、収束判定を行う。ここでは、例えば、求められた目的関数の変動両が予め記憶された基準値を下回った場合に、収束したものとする。ここで、収束していないと判定された場合には、タイムステップを単位時間だけ進め、タイムステップの更新（ステップＳ８）を行った後、再びステップＳ４を行う。一方、ステップＳ７で収束していると判定された場合には、最適化が完了したものとして終了する。

【0046】

次に、上述したような本発明の考え方を実現するための本実施形態の流路形状最適化装置１０について説明する。図２は、本実施形態の流路形状最適化装置１０のハードウェア構成を概略的に示すブロック図である。

【0047】

本実施形態の流路形状最適化装置１０は、パーソナルコンピュータやワークステーション等のコンピュータによって構成されるものであり、ＣＰＵ１（演算処理手段）、記憶装置２、リムーバブルメディアドライブ３、入力装置４、出力装置５、及び通信装置６を備えている。

【0048】

ＣＰＵ１は、記憶装置２、リムーバブルメディアドライブ３、入力装置４、出力装置５、及び通信装置６と電気的に接続されており、これらの各種装置から入力される信号を処理すると共に、処理結果を出力するものである。このＣＰＵ１は、流路形状最適化プログラムＰに基づいて、流れ場の形状が更新されなくなる状態となるまで収束計算を行う。なお、流れ場の形状が更新されなくなる状態とは、例えば、流れ場が定常状態となった場合や、流れ場の形状変化の程度が予め設定された基準値よりも小さくなった場合を示す。

【0049】

記憶装置２は、メモリ等の内部記憶装置及びハードディスクドライブ等の外部記憶装置によって構成されており、ＣＰＵ１から入力される情報を記憶すると共にＣＰＵ１から入力される指令に基づいて記憶した情報を出力するものである。この記憶装置２は、プログラム記憶部２ａとデータ記憶部２ｂとを備えている。

【0050】

プログラム記憶部２ａは、流路形状最適化プログラムＰを記憶している。この流路形状最適化プログラムＰは、所定のＯＳにおいて実行されるアプリケーションプログラムであり、コンピュータから構成される本実施形態の流路形状最適化装置１０に、図１で示すフローチャートに沿った処理（流路形状最適化方法）を実行させる。

【0051】

データ記憶部２ｂは、多孔質媒体モデルデータＤ１、目的関数データＤ２、計算条件データＤ３、初期値データＤ４、計算結果データＤ５等を記憶する。多孔質媒体モデルデータＤ１は、上述した多孔質媒体モデルを定義するデータが含まれており、例えば、多孔質媒体モデルの形状データや、特性データが含まれている。目的関数データＤ２は、例えば、［数３］〜［数５］の式を示すデータである。計算条件データＤ３は、［数２］、［数６］、［数７］、［数１２］及び［数１４］の式を示すデータである。初期値データＤ４は、流体の流量や流速、流体の物性値（粘性等）を示すデータである。これらの多孔質媒体モデルデータＤ１、目的関数データＤ２、計算条件データＤ３及び初期値データＤ４は、例えば、リムーバブルメディアドライブ３、入力装置４あるいはネットワークＮ等を介して入力される。なお、計算結果データＤ５は、ＣＰＵ１が計算した結果を示すデータであり、最終的な演算結果やその途中で生成される中間データを含む。

【0052】

リムーバブルメディアドライブ３は、ＤＶＤ（Digital Versatile Disc）やＵＳＢ（Universal Serial Bus）メモリ等のリムーバブルメディアを取り込み可能あるいは接続可能に構成されており、ＣＰＵ１から入力される指令に基づいて、リムーバブルメディアに記憶されるデータを出力するものである。例えば、リムーバブルメディアに流路形状最適化プログラムＰが記憶されている場合には、リムーバブルメディアドライブ３は、ＣＰＵ１から入力される指令に基づいて、リムーバブルメディアに記憶される流路形状最適化プログラムＰを出力する。なお、流路形状最適化プログラムＰは、必ずしもリムーバルディスクに格納されている必要はない。例えば、ネットワークＮを介して流路形状最適化プログラムＰの取得を行い、この取得した流路形状最適化プログラムＰを記憶装置２に記憶させるようにしても良い。

【0053】

入力装置４は、本実施形態の流路形状最適化装置１０と作業者とのマンマシンインターフェイスであり、ポインティングデバイスであるキーボード４ａやマウス４ｂを備えている。出力装置５は、ＣＰＵ１から入力される信号を可視化して出力するものであり、ディスプレイ５ａ及びプリンタ５ｂを備えている。

【0054】

通信装置６は、本実施形態の流路形状最適化装置１０とネットワークＮとを電気的に接続し、ネットワークＮに接続された他の機器との間においてデータの受け渡しを行うものである。なお、ネットワークＮとしては、社内ＬＡＮ（Local Area Network）やインターネット等が考えられる。

【0055】

このような本実施形態の流路形状最適化装置１０は、流路形状最適化プログラムＰに基づいて、図１で示すフローチャートに沿った処理（流路形状最適化方法）を実行する。すなわち、ＣＰＵ１は、流路形状最適化プログラムＰに基づいて、ステップＳ１〜ステップＳ８の処理を実行することによって、最適化された流路形状を求める。

【0056】

まずＣＰＵ１は、初期条件の設定を行う。ここでは、ＣＰＵ１は、リムーバブルメディアドライブ３、入力装置４あるいは通信装置６から入力される目的関数、流体の流量や流速の、流体の物性値、多孔質媒体モデルデータ等を、多孔質媒体モデルデータＤ１、目的関数データＤ２、計算条件データＤ３及び初期値データＤ４としてデータ記憶部２ｂに記憶させる。続いて、ＣＰＵ１は、これらのデータを用い、流路形状最適化プログラムＰに基づいて、格子ボルツマン法を用いた流体解析を行い、計算が収束あるいは周期的な流れが形成された時点で計算を終了する。これらの工程は、上述した初期化工程に相当し、この工程によって多孔質媒体モデルに対して入口と出口とを繋ぐ流路が形成される。

【0057】

ＣＰＵ１は、上述のように多孔質媒体モデルに対して流路が形成されると、この流路形状の最適化を行う。ここでは、ＣＰＵ１は、まず流路形状最適化プログラムＰに基づいて、格子ボルツマン法を用いた流体解析で非定常流れを１タイムステップのみ流れ場の計算を行う。ＣＰＵ１は、この計算を、多孔質媒体モデルデータＤ１、目的関数データＤ２、計算条件データＤ３及び初期値データＤ４を用いて行ない、その計算結果を計算結果データＤ５としてデータ記憶部２ｂに記憶させる。さらにＣＰＵ１は、流路形状最適化プログラムＰに基づいて、上述のステップＳ５で説明した感度の計算及びステップＳ６で説明した空隙率の更新を行うことで流路形状（流れ場の形状）を変化させ、解が収束するまでタイムステップの進行を行いながら、格子ボルツマン法による流体解析、感度の計算及び空隙率の更新を続ける。そして、ＣＰＵ１は、解が収束すると、その流路形状（流れ場の形状）を最適解とする。

【0058】

以上のような本実施形態の流路形状最適化装置１０及び流路形状最適化方法によれば、流れ場が非定常状態のうちに流れ場の計算と、流れ場の形状を変更した場合の感度計算と、感度に基づく流れ場の形状の更新とを繰り返す。このため、最終的に解が収束して流れ場の形状が更新されなくなる状態となったときには、最適な流れ場の形状（すなわち流路形状）を得ることができる。このような本実施形態の流路形状最適化装置１０及び流路形状最適化方法によれば、最適化工程において１度の収束計算によって流路の最適形状を得られるため、流体解析による流路形状の最適化に要する時間を極めて短縮することができる。

【0059】

また、本実施形態の流路形状最適化装置１０及び流路形状最適化方法によれば、格子ボルツマン法を用いた流体解析を行うことによって流れ場を計算している。格子ボルツマン法を用いた流体解析は、単純な四則演算のみで数値計算をおこなうため、ナビエ−ストークス方程式を用いた流体解析と比較して収束性が良くかつ計算時間も短い。このため、本実施形態の流路形状最適化装置１０及び流路形状最適化方法によれば、安定してかつ短時間で流路形状の最適化を行うことが可能となる。

【0060】

また、本実施形態の流路形状最適化装置１０及び流路形状最適化方法によれば、各位置における空隙率を設定可能な多孔質媒体モデルを用い、この多孔質媒体モデルの内部に流れを形成し、感度に基づいて空隙率を更新することによって流れ場の形状を更新している。このため、空隙率の更新により容易に流れ場の形状を変化させることができる。

【0061】

なお、図３〜図１０を参照して、上記実施形態に即した本発明の実施例について説明する。

【0062】

（第１実施例）
図３（ａ）は、本実施例で用いた二次元多孔質媒体モデルを示す模式図である。この図に示すように、本実施例においては、左辺全域が流体の入口とされ、右辺の中央が流体の出口とされた多孔質媒体モデルを用いている。また、本実施例では、レイノルズ数を１０としている。図３（ｂ）は、本実施例におけるタイムステップと流体損失（目的関数の計算値）との関係を示すグラフである。図３（ｂ）において、「Ｄ」は粘性散逸関数に基づくグラフであり、「Ｃ」は摩擦損失関数に基づくグラフであり、「ｔｏｔａｌ」は粘性散逸関数と摩擦損失関数とに基づく（すなわち目的関数に基づく）グラフである。また、図４は、タイムステップ（図中Ｓｔｅｐ）が０、１００、２００、３００、４００、５００、７００、１０００のときの流路形状を可視化した図である。なお、図４においては、黒く示している部分が流路である。

【0063】

本実施例においては、図３（ｂ）及び図４から分かるように、４５０ステップまでは物理的に意味のある解を探索し、その後細部の最適化に移行している。また、本実施例において計算が破綻することはない。このような本実施例における計算時間は３０．１秒であった。一方、従来手法（流路形状を変更するたびに収束計算を行う手法）で同様の最適化をおこなったところ、計算結果は同じであったが、計算時間は約５分であった。このことから、本発明が極めて短時間で最適化を完了できることが確認された。

【0064】

（第２実施例）
図５（ａ）は、本実施例で用いた二次元多孔質媒体モデルを示す模式図である。この図に示すように、本実施例においては、左辺の中央及び下辺の左寄りの部位が流体の入口とされ、上辺の左寄りの部位が流体の出口とされた多孔質媒体モデルを用いている。また、本実施例では、レイノルズ数を１０としている。図５（ｂ）は、本実施例におけるタイムステップと流体損失（目的関数の計算値）との関係を示すグラフである。図５（ｂ）において、「Ｄ」は粘性散逸関数に基づくグラフであり、「Ｃ」は摩擦損失関数に基づくグラフであり、「ｔｏｔａｌ」は粘性散逸関数と摩擦損失関数とに基づく（すなわち目的関数に基づく）グラフである。また、図６は、タイムステップ（図中Ｓｔｅｐ）が０、１０００、２０００、３０００、４０００、５０００、６０００、６８００のときの流路形状を可視化した図である。なお、図６においては、黒く示している部分が流路である。

【0065】

本実施例においては、図５（ｂ）及び図６から分かるように、３０００ステップまでは物理的に意味のある解を探索し、その後細部の最適化に移行している。また、本実施例においては部分的に乱流になり、目的関数が振動した。このような本実施例における計算時間は１９５．５秒であった。

【0066】

（第３実施例）
図７（ａ）は、本実施例で用いた二次元多孔質媒体モデルを示す模式図である。なお、図７（ａ）における矢印が流れのベクトルを示しており、本実施例では多孔質媒体モデルの左辺全域が入口、右辺全域が出口とされ、上辺と下辺にも境界条件として流れ（外部流れ）が与えられている。つまり、この図に示すように、本実施例においては、二次元多孔質媒体モデルが流れの中に晒されている。また、本実施例では、レイノルズ数を０．１としている。図７（ｂ）は、本実施例におけるタイムステップと流体損失（目的関数の計算値）との関係を示すグラフである。図７（ｂ）において、「Ｄ」は粘性散逸関数に基づくグラフであり、「Ｃ」は摩擦損失関数に基づくグラフであり、「ｔｏｔａｌ」は粘性散逸関数と摩擦損失関数とに基づく（すなわち目的関数に基づく）グラフである。また、図８は、タイムステップ（図中Ｓｔｅｐ）が０、３００、５００、７００、１０００、１３００、１５００、２５００のときの流路形状を可視化した図である。なお、本実施例では、図８において黒く示している部分が流路であり、図８（ａ）に描かれている円形の物体が、形状の初期状態として多孔質媒体モデルとして与えられており、この物体の形状が計算の過程で変化する。

【0067】

本実施例においては、６００ステップまでは体積制約を満たす解を探索し、その後細部の最適化に移行している。このような本実施例における計算時間は７２．６秒であった。

【0068】

（第４実施例）
図９（ａ）は、本実施例で用いた二次元多孔質媒体モデルを示す模式図である。この図に示すように、本実施例においては、左辺の中央が流体の入口とされ、右辺の中央が流体の出口とされ、さらに中央に円部を有する多孔質媒体モデルを用いている。なお、本実施例では、円部は回転される円柱等を模擬しており、回転速度Ｕｒｏｔで時計回りに回転されている。また、本実施例では、レイノルズ数を１００としている。図９（ｂ）は、本実施例におけるタイムステップと流体損失（目的関数の計算値）との関係を示すグラフである。図９（ｂ）において、「Ｄ」は粘性散逸関数に基づくグラフであり、「Ｃ」は摩擦損失関数に基づくグラフであり、「ｔｏｔａｌ」は粘性散逸関数と摩擦損失関数とに基づく（すなわち目的関数に基づく）グラフである。また、図１０は、タイムステップ（図中Ｓｔｅｐ）が０、１０００、２０００、３０００、４０００、５０００、６０００、９５００のときの流路形状を可視化した図である。なお、図１０においては、黒く示している部分が流路である。

【0069】

本実施例においては、８０００ステップまでは全領域で乱流化している。また、流体領域と目的関数とが振動し、最終的にはほぼ層流になって安定した。このような本実施例における計算時間は２６６．７秒であった。

【0070】

以上、添付図面を参照しながら本発明の好適な実施形態について説明したが、本発明は、上記実施形態に限定されないことは言うまでもない。上述した実施形態において示した各構成部材の諸形状や組み合わせ等は一例であって、本発明の趣旨から逸脱しない範囲において設計要求等に基づき種々変更可能である。

【0071】

例えば、上記実施形態においては、最適化工程において格子ボルツマン法を用いた流体解析を行う構成について説明した。しかしながら、本発明はこれに限定されるものではなく、収束性の低下について許容できるのであれば、最適化工程において他の流体解析手法を用いても良い。

【0072】

また、上記実施形態においては、空隙率を変化させることによって流路の形状を変化させる構成について説明した。しかしながら、本発明はこれに限定されるものではなく、他の方法にて流路の形状を変化させる構成を採用することも可能である。

【符号の説明】

【0073】

１ ＣＰＵ（演算処理手段）、２ 記憶装置、２ａ プログラム記憶部、２ｂ データ記憶部、３ リムーバブルメディアドライブ、４ 入力装置、４ａ キーボード、４ｂ マウス、５ 出力装置、５ａ ディスプレイ、５ｂ プリンタ、６ 通信装置、１０ 流路形状最適化装置、Ｄ１ 多孔質媒体モデルデータ、Ｄ２ 目的関数データ、Ｄ３ 計算条件データ、Ｄ４ 初期値データ、Ｄ５ 計算結果データ、Ｎ ネットワーク、Ｐ 流路形状最適化プログラム

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】