計測自動制御学会東北支部 第165回研究集会（1卵6．12．13）
資料番号165−13

流動囁のオブザーバに関する数値実験
NumericalExperimentonObseTYerOfF］owField
O早瀬 敏幸＊，林 叡＊

O ToshiyukiHAYASE＊andSatoruHAYASIⅡ＊＊東北大学 流体科学研究所＊Institut巳OfFluidScience，TohokuUniverslty
キーワード：計算機援用流動場制御（Computer−Aid巳dFlowFieldControl），数値流体力学（ComputadonalFluidDynamics），オブザーバ（Observer），分布定数系PistributedParameterSystem），非線形系（NonlinearThTnamicSystem）

連絡先：〒980−77仙台市青葉区片平 東北大学流体科学研究所 流動場制御研究部門
Tel：（022）217，5253，Fax：（022）217−5253，E−mail：hayase＠血tohoku．acjp

1．轄への拡張に関する研究がさかんに行われているが（ヱ〕，流れ場の制御を含めて，非線形の分布

本研究では，制御系を構成する高性能の計定数系の制御に関する研究は少ない．これを算機内部に流れ場の数学モデルを構成し，流実現するためには，理論的な面からの検討以

動場を非線形かつ分布定数系の制御対象とし外にも，空間に連続的に分布した状態量を測て制御するための基礎的嘩討を行う・流動場定するためのセンサ，および制御入力を生成の制御は現在広く行われているが，通常は流するためのアクチュエータの開発など，解決れ場を比較的低次元のシステムでモデル化し，すべき多くの問題がある．

少数のセンサとアクチュエータで制御系を構

近年めざましい性能向上を遂げたパーソナ成するものがほとんどである．たとえば気流

ルコンピュータは，低次元の制御系に関して中の航空機の飛行制御などのように，適切なは十分な性能を有しており，複雑な行列演算モデル化を行えば，比較的低次元の制御問題等を含む最近の制御系の構成要素として幅広とく用いられている．本研究は，コンピュータ方，流れ場申乱流変動等を直接制御しようとの性能が現在も急速に向上していることを踏する場合には，流れ場本来の非線形かつ分布まえて，現状のスーパーコンピュータを制御定数系として平取り扱いが必要となる・現屯系の基本構成要素として用いることにより，線形系を対象として，非線形・分布定数系としで流れ場を制御する基本的な構造を図1（a）に示す．流れ場を制御ことを目標とし，そのための基礎的検討を行対象とし，それに作用する境界条件と，体積うものである．

力を入力と考える．非圧縮粘性流体の流れ場

本論文ではその第1段階として，制御を行う際に必要となる流動場の状態量の推定の問は，ナビエ・ストークス式

題を取り上げる．系の数学モデルと，出力の＝−gradp−（u・grad）u・∇2u・f，（1）測定値から内部の状態量を推定するオブザーバ（状態観測器）の理論（3）を流れ場に適用し，有限体積法に基づく差分近似数学モデルと，流れ場の有限個の測定値より，状態量の推定を行う．最初に，流動場のオブザーバの構造について検討した接，正方形管埠申乱流の状態量を推定する問題について数値シミュレー

連続式

divu＝0，（2）および初期条件，境界条件によって支配される．これを固2に示すような無限次元の動的システムとしてとらえると，状態量は領域内ションを行い，提案したオブザーパの有効性を検討する．の速度ベクトル場u（rカであり，この速度場l巾，王）と対応する圧力場紳，Ⅹ）の関数である有限個の観測量（例えば，流れ場中町物体表面

コ．流動場の制御の圧力およびせん断応力の分布，あるいはそ制御理論の立場から見た流動場の制御系のれらの積分値である抗九出力が定義される．この出力信号は；図1（a）に示すように基準信号と比較され，その偏差
Sbl亡F虎dba慮巾，タ）が制御装置に入力される．制御装置は，制御規則とアクチュエータよりなる∴制御規則は，偏差信号と流れ場の状態量よりテクチュエータへの入力を決定する．アクチュエータは，流れ場の境界条件，体積力のいずれか，ある（a）Fundarnental免edbackcontroIsystem

5【加F亡edbad【hp）C心彗erV可・（b）FeedbackcontroIsystemwithstateobserver

Fig．1BlockdiagramofflowfieldcontroIsystem

Fig・2 Blockdiagramoftheflowsystem■1■

いはその両方を発生する．この制御系の構成際の流れ場に収束するが，乱流場やカルマンは，形式的には通常のフィードバック制御系

渦を伴うような不安定な流れ場に対しては，と同一であるが，系が強い非線形性を有する

たとえ同一の初期条件から計算を開始しても，ことと，状態変数とアクチュエータの自由度変動成分を含めて実際の流れ堤に収束する解を無限大としていることが本質的に異なってを得ることは困難である．
これに対して，オブザーバでは，実際の流いる．

本論文では，まず系の状態量に注目し，無れ場の観測出力と，シミュレーションの対応限自由度の状態量の代わりに十分高次の離散する計算値との偏差を求め，この信号を適切ベクトル場を与えるオブザーバを構成する．に数学モデルにフィードバックすれば，不安

図1（叫こオブザーパを含んだ流れ場の制御系

定な流れ場の場合でもシミュレーションの結の構成を示す．オブザーパは，流れ場に与え果を実際の流れ場に収束させ，また安定な流るのと同一の制御入力を数学モデルに与え，れ場に対しては定常解への収束を加速するこ

その出力を実際の観測出力と比較し，偏差をとが可能となる．さらに，数値解析手法に存数学モデルにフィードバックすることにより，速やかに数学モデル内部の状態量を実際の状態量に収束させるものである．なお，図1仲）ではオブザーパの動作に注目するため，制御装置とオブザーバのブロックを別々に示したが，実際の制御系では制御規則とオブザーパは同一のコンピュータ内部にブログラムされる．

3．オブザーパの構成（a）Flowsimulation

以下に，流れの数値シミュレーションのアルゴリズムを系の数学モデルとして用いた場合の，オブザーバの構成について述べる．通常の流れのシミュレーションとの基本的な差

異は，出力誤差信号のフィードバックの有無にある・すなわち，図3（a）に示すように，通常の流れのシミュレーションでは，境界条件と体積力を実際の流れに一致させた場合，もし対象とする流れが大域的に安定ならば，初期条件が異なっていても，シミュレーションの結果は数値計算手法の精度の範囲内で，実（b）Stateobserver

Fig．3Comparisonbetweenflowsimulationandstate
observer在する誤差（有限の格子間隔による打ち切りの有効性を検討する．すなわち，制御対象で誤差など）を低減させることも期待できる．

ある実際の流れ場の代わりに，前もって計算

有限次元の時不変線形システムにおいては，した発達乱流の数値解（以下「基準解」と呼出力誤差信号の線形関数としてモデルへの入ぷ）を用い，制御対象の数学モデルとして全力を与えることにより，任意の収束性をもつての速度成分を0とおいた初期値からの数値オブザーバを設計することができるが（l），非線シミュレーションを行って，これに適当なフ形系に対するオブザーパの一般的な理論は完ィードバックを施すことにより，数学モデル成されていない（ユ）．本論文では，流れ場の数学の解を基準解に収束させようとするものであモデルとして，有限体積法に基づく流れ場のる．

シミュレーションスキーム（4）を用いてオブザ

本論文で取り扱う対象領域および座標系をーバを構成する．流れ場・を対象とした場合，図4に示す．計算領域は正方形断面管内の流系の可観測性を理論的に検討することは困難

れ場で，管路長Jは正方形管の1辺の長さの4であり，さらに入力である境界条件および体倍とする・非圧嘩・粘性流体を仮定し，ナビ積力を，出力誤差の関数として生成する一般

エ・スト一夕ス式と連続式を基礎式にとる．

的手法は確立され七いないので，ここでは流速度の境界条件革も 壁面上では滑りなしの条れ場のある境界面上での速度誤差を，その点件，上流と下流の断面には周期境界条件を課の圧力境界条件にフィードバックする最も単した．圧力については，上流断面と下流断面純な場合について数値解析により検討する．の対応する各点間に圧力差を与えた．
数値解析手法については，既報〔5）で詳しく述

4．数値計算例べているので，ここでは計算の概略を述べるにとどめる．等間隔の矩形スタガード格子系

著者らは既報で，正方形管内の発達乱流に関する数値シミュレーションを行った（5）．以下
でⅠも この数値計算結果を用いてオブザーバを用い，有限体積法に基づいて離散化された基礎方程式群を，SIMPLER法に類似の手法により解いた．対流項の離散化には，物理的考察に基づき再定式化されたqⅥ亡兄スキーム

TablelComputalional亡Onditions．

Pipelen昌仙LJ

4Pressure血飽renc∈ ∠b0．0649

Reynoldsnumber Res

9000

Gridpoints 弗朝AN320xlOxlO

GridspacL喝：hlh2h30．2xO，1xO．1

Timeincrement hL0．05を用いた（句．また時間微分項の離散化には，2〟1如＋む

Total residual at canvergence0．01
CPUtimefbronetimestep 5］ 2Fig．4I）omainandcoordinatesystem

する際の出力および入力（図1参照）は，オブザーパの出力および入力に，制御に必要な信号を追加して得られる】．各位置において，軸速度成分の基準解からの誤差を求め，それに適当なフィードバック係数を掛けてその位置での圧力境界条件に加える．その結果，断面上の各点で，流体は軸方向の速度誤差を減少させる方向に加速あるいは減速される．なお，入力の各位置での圧力は，一般のオブザ
Fig．5 StruCtureOfthestateobserver

ーバでは断面内の全ての出力変数（軸速度）の関数として与えられるべきものであるが，次精度の陸解法を用いた．

主な計算条件を表1に示す．なお，諸量は正方形管の1辺の長さ占，平均軸速度以mおよび流体の密度pを用いて無次元化されている．
基準解として，レイノルズ数粘（＝〟m占／巧＝9000の発達乱流解が得られるように圧力差車（p打m2で無次元化）を設定した．ここで用いた格子系は，乱流場を精度よく再現するのに十分な格子解像度をもっているとは言い難いが，既報でも述べたように，得られた解は平均軸速度分布を良好に再現し〔5），また比較的振幅の大きな低周波領域の変動成分を良好に表している．このことは，本論文の基準解が実際の乱流場にローバスフィルタを通過させたものに対応していることを意味する．本研究の目的である乱流場の制御には，エネルギの大きな低周波乱れ成分の抑制が有効であると予想されるので，乱流場の低周波成分を保存した
本基準解への収束性によりオブザーバの有効

ここでは対応する位置における軸速度のみの1次阻数として与え，さらに，その係数も断面上で全て一定であるとした最も単純な場合を取り扱った．
フィードバックゲインた（図5参照）を種々に変化させた計算の結果を図6（a），（b）に示す．
図は，エl＝0の断面中心における軸速度成分の時間変化を示したものである．図中の太線は，基準解（制御対象）の発達乱流を示している．
他の結果は，時刻声0において全ての速度成分を0として非定常数値計算を開始したオブザーバの出力である．図6（a）の破線はゲイン丘＝0の場合で，フィードバックを行わない通常の数値シミュレーションに対応する．この場合，与えられた圧力差4）（一定）による流れ場の加速は非常に緩やかである．時定数て（軸速度の変化を1次遅れ系とLたときに，定常値の郎％まで増加するのに要する時間）の借を破線の傾きより見積もるとて＝120（無次元）性を評価することは妥当であると考えられる．となり，フィードバックを行わない場合には，オブザーバの構成を図5に示す．本論文では，出力として，上流端断面上での軸方向速度成分をとり，入力としてその断面に作用する圧力をとった［なお，流動場の制御系を構成発達乱流解に達するのに非常に長い時間を要することがわかる．さらにこの場合には，発達乱流解に落ちついた接でも，変動成分が基準解のそれに収束することは期待できない．

ることにより，上流断面中心の軸速度成分が，このことは，フィードバックを施さないシミ速やかに基準解のそれに収束することが示さュレーション結果を用いて乱流れた．

ることが，本質的に困難であることを示している．

次に，フィードバックゲイン拒2の場合に，上流端断面中心における叱打5速度成分と圧力

一方，適当なゲインを用いてフィードバッpの収束状況を囲7に示す．図7（a），桝より，クを行った場合には，速やかに基準解に収束

速度成分に関してオブザーバの結果はPlOのしている．例えば，図6（a）の庵＝ヱにおける時範囲で基準解の変動成分によく追従している．

定数はT可となり，フィードバックを行わな

ただし，叫3速度成分については，本計算条件い場合に比べて60倍の収束性の改善が見られの範囲内でみられる基準解のバイアス成分がる．さらに変動成分についても，オブザーパ十分再現されていない．また圧力に関しては，の結果は基準解によく追従していることがわ

fく10の範囲で，オブザーパの結果はフィードかる．これは，断面上の各点において，基準バックによる流体の加速に対応して大きな値解からの誤差に応じて，これを減少させるよをとるが，流れ場が基準解に収束するにつれうに圧力差に修正を加えるフィードバックのて，圧力も基準解のそれに追従している．

効果である．囲6（孔）では，フィードバックゲイ

これまでは，上流断面の中心における各速ンの増加とともに，収束性，追従性が向上し度成分と圧力の収束状況を調べたが，流れ場ているが，更にゲインを増すと，やがて系は全体の収束状況につt不安定となる．図印小も過大なフィードバ

は，ゲイン拒2の場合に，3つの異なる断面ックゲインによりオブザーパが不安定化した上（∫l＝0，2，ヨ）の，4つの異なる代表位置につ例である．以上より，適当なゲインを選択すいて，それぞれ軸速度成分の収束状況を示し

0 ユ 4 6 8 10 12 1ヰ 16 18 ユO

0 2 4 石 亀 10 1ヱ 14 16 18 20tt（a）Convergentresults（b）Unstableresult鈷rlarge免edbackgain

Fig・6 Variationoftheaxialvelocityatthecenteroftheupstreamboundary（Comparisonfbrtheftedbackgain）

たものである．図8（a）−（c）は，断面の対角線に下流方向への劣化の程度は小さく，一方，図沿って，g（亡），（d）に見るように，位置が管壁に近づき，中心からコーナーに向かう位置，国旨（d）は，対角線上にない壁面近傍の位置におけかつ対角線から離れるほど劣化の程度は大きる結果である（中列の回申にそれぞれの断面くなる．

内の位置を示した）．図の左列の上流断面￠1＝0）

上で述べたように．オブザーパにより得らの結果では，軸速度の誤差が圧力境界条件にれる結果の基準解に対する追従性は，領域内直接フィードバックされるため，何れの位置の位置によってかなり異なる．そこで，あるにおいてもオブザーパの結果は基準解に良〈領域全体におけるオブザーバの平均的な追従追従している．図の中列，右列と下流断面に性を定量的に評価す．るために，速度場に関す向かうに従って，オブザーバの結果の追従性る以下のノルムを定義する．

は劣化しているが，その程度は，断面内の位置によって異なる．すなわち，図呂（a），仲）のよ

1 1【u（叫し＝【い小串仰げ】盲

うに，断面中心に近い対角線上の位置では，（ヨ）

ここで，積分は対象領域Ⅴにわたる体積積分を表す．また，rは，領域の体積を表す．
上記のノルムを用いれば，オブザーパにより得られた速度場u（ちエ）の基準解u＊（ちⅩ）に関する誤差ノルムが，時間の関数として次式○ ユ 4 6 葛 10 1ユ 14 16 柑 ユO

で評価できる．

t（a）zLユVelocityCOmPOnentu（痛一u●（叫し1 ＝出（叫一顧・（甜ユー琉㌧（町軒軒ザ（4）

0 コ 4 6 箕 10 1ユ 14 1古 川 加tO））u3VelocitycoInpOnent

式（4）で与えられる誤差ノルムの時間変化を，ゲイン腱ユおよび后4の場合について示したものが図9である．図よりこ初期時刻巨0においては，オブザーパの初期値として全ての速Stateobserver

度成分を0としたため，誤差ノルムの値はほStandard solution

ぼ1となっている．その後FlO程度までは，（con加11edsystem）
0 2 4 d ‡ 10 1ユ 14 15 1琶 之Ot（c）Pressure

Fi苫．7 Convergenceofst血Observeratthe
CenterOftheupstreamboundary（k＝2）

基準解への収束に伴って誤差ノルムは単調に減少し，それ以降はほぼ一定値となる．2種類のゲインに対する結果を比較すると，収束性は，ゲインの大きな場合のほうが優っているが，誤差ノルムの定常値についてはほとんXl＝0

xl＝3Ⅹl＝ユE−■−−−■一■−一一一一■■−■・・一一一一■＿一■一−−r■・ け−−−・20020t（a）氾＝0．45，粘⇒〕J520020t匝）xユ三0．25，粘＝0．苫

t20020tt20t

2I3（c）池＝0．05，Ⅹ戸0．05

t20020t（d）粘＝0．25，Ⅹ戸0．05

Fig，8 Variationoftheaxialvelocitycomponentatseverallocationsofthedomain（k＝2）

ど差は見られない．

また図9中の実線は，初期条件の異なる2し，声2以降ではほぼ一定値（0．2）となる．このことは，先に述べたように，フィードバックつの発達乱流解の間の誤差ノルムを．式（4）中を行わない数値シミュレーションにより発達の u（r刃 を u¶0刃 で 置 き 換 え た

乱流解が得られたとしても，基準解との間にlu■（0，小u●（叫ivで評価したものである・

は，常に0．2程度の誤差ノルムが存在すること

時刻戸0においては，この値は当然0となるが，を示している．一方，先に得られたフィードバ

どの増加とともに，2つの解の間の差異は増大ックを施した場合の定常誤差ノルムの平均値

′h︶
4

コ績ヲコニ

0∧U＞

00 2 4 6 8 10 12 14 16 柑 200123

ヰtⅩlFig．9 Normoftheestimationerrorinthe
domainasafunctionoftime

は0．1ユで，フィードバックにより速度場全体

Fig・10 Normofthesteadyestimationerror
evaluatedateach肌iallocationで平均した結果（図9）では，慮＝2と丘＝4の誤の誤差ノルムが半分程度にまで減少すること

差ノルムの定常値に差は認められなかったが，を意味している．

図10では，特に上流端近傍で，ゲイン鬼が大上で述べた結果は，対象領域全体にわたるきいほど誤差ノルムが減少しており，腱4の場オブザーパの，基準解への追従性を表したも合には，上流端で，誤差ノルムがフィードバのであるが，既に述べたように，この追従性ックを施さない場合の1／4程度にまで減少しは流れ方向の各断面によって異なっていたている．

（図8参照）．そこで，式（4）の積分領域Ⅴを，以上，上流断面の各点において軸速度成分∫l軸に垂直な断面A（れ）で置き換えて，軸方向の誤差を圧力境界条件にフィードバックするの各位置における誤差ノルムを求め，これが単純なオブザーバによっても，基準解への良十分収束した時間の範囲（10＜fくユ0）で平均したものを囲10に示す．なお，図中の実線は，上好な追従が行われることを，速度場のノルムを用いて定量的に示した．にも述べたようにulO，Ⅹ）とu竹，Ⅹ）より求めた誤差ノルムであり，独立な2つの発達乱流解5．結に対する誤差ノルムを与えている．図中の3種類のゲイン慮＝1，ユ，4に対する結果をみると，計算機を援用した流動場の制御を実現する

ともに上流端でもっとも誤差が小さく，下流ための基礎的研究として，本論文では，オブに行くに従って誤差が増大し，フィードバッザーバーの構成に関する検討を行った．流れクを施さない場合の結果（実線）に漸近して場の数値シミュレーションのアルゴリズムをいる．なお，エ】≧3．6でみられる不自然な値の数学モデルとして，流れ場の観測量と計算結減少は，周期境界条件の影響であり，管路長J果との誤差をフィードバックすることにより，をさらに増せば単調に実線に漸近するものと状態量の推定を行うオブザーパーの基本構成考えられる．また，誤差ノルムを流れ場全体を示すとともに，基本的な正方形管路内の発達乱流盲対象として，オブザーバーの有効性参考文献を数値解析的に検証した．すなわち実際の流れ場を，前もって計算した発達乱流の数値解で置き換え，初期値を全て0とおいた数値計（1）R．E，Skelton：DynamicSystemsControl，JohlWiley＆sons，（198即（2）X，LiandJ．Yong：OptimalControITheoryfbr

算を行って，上流断面上の軸速度の誤差を圧

hfiniteDimensionalSystems，Birkhauser，力境界条件にフィードバックした．適当なゲ（3）E，A．MisawaandJ．K．Hedrick：Nonlinear

インを選択することにより，フィードバック

ObserYerS−AState−Of」the−ArtSurvey，J，を行わない場合に比べて，発達乱流解への収束性が60倍程度加速され，定常乱流状態の誤（1995）

Dyn．Sys．andMea5．，TranSASME，111−3，344／352（1989）（4）T．Hayase，J．A．C．HumphreyandR．Greif：

差ノルムも流れ場全体では1仏上流端近傍でMimiManualhrROTFLO2，Dept．Mech．Eng．
Rep・，＃FM−90−1，Univ．Calif．Berkeley，は1／4程度に減少させることが出来た．

（1990）

今後の課題としては，比碗ゲインの増加に伴う系の不安定性に関して理論的検討を行うとともに，本論文では断面内で一定とした比（5）早瀬，末松，菊池：正方形断面管路内乱
流の直接数値シミュレーション，日本機械学会論文集，5邑・54占，B，3糾／370（1992）（6）T．Hayase，J．A．C．HumphreyandR．GreiftA

例ゲインを，各点毎に異なるとした場合，さConsistentlyFormulatedQUICKSchemefbr
FastandStableConvergenceUsingFinite−

らに，各点の入力（圧力境界条件）が，周囲VolumeIterativeCalculationProcedure，J．の各点の軸速度の線形関数で与えられるような，より一般的な場合についても検討したい．

Comput．Phys．，98−1，108／11i（1992）