JP2018021693A - 温度分布推定装置および温度分布推定方法 - Google Patents
温度分布推定装置および温度分布推定方法 Download PDFInfo
- Publication number
- JP2018021693A JP2018021693A JP2016152212A JP2016152212A JP2018021693A JP 2018021693 A JP2018021693 A JP 2018021693A JP 2016152212 A JP2016152212 A JP 2016152212A JP 2016152212 A JP2016152212 A JP 2016152212A JP 2018021693 A JP2018021693 A JP 2018021693A
- Authority
- JP
- Japan
- Prior art keywords
- temperature distribution
- time
- initial
- furnace
- physical model
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Images
Landscapes
- Vertical, Hearth, Or Arc Furnaces (AREA)
- Waste-Gas Treatment And Other Accessory Devices For Furnaces (AREA)
- Coke Industry (AREA)
Abstract
【課題】限られた数の温度センサ情報から精度良く炉1内温度分布を推定する。【解決手段】向流熱交換型の反応炉1内の温度分布を推定する。設定した過去の設定時刻での初期温度分布を初期値として物理モデルを時間順方向に解いて現在の温度分布を推定し、検出した検出温度の時系列情報に基づき時間逆方向に設定時刻までアジョイントモデルを解いて、設定時刻でのアジョイント変数の値を補正情報として求めて初期温度分布を修正する。修正した後の初期温度分布を初期値として物理モデルを時間順方向に解いて設定時刻よりも後の時刻での炉1内の温度分布を推定する。初期温度分布の修正を誤差が収束するまで繰り返す。【選択図】 図2
Description
本発明は、高炉、コークス乾式消火設備(CDQ)及びガス化溶融炉などで例示される反応炉のような、炉の上側から固形分を主要成分とした原料が装入されると共に、炉の下側から原料と熱交換されるガスが供給されるシャフト型の反応炉に代表される、向流熱交換型の反応炉内の温度分布を推定するための温度分布推定の技術に関する。
製鉄業の施設においては、向流熱交換型の反応炉を数多く使用している。向流熱交換型の反応炉では、上側から固形分を主要成分とする原料が装入され、その上側からの原料と下側から供給されるガスとの間で向流熱交換が行われる。そのような向流熱交換型の反応炉においては、炉内部における多数箇所での炉内温度を計測することが困難、つまり炉内温度分布のオンライン測定が困難な場合が多い。
例えば、高炉においては、上部から原料を装入し、予熱、還元、溶解の各過程を経て、炉の下部から製品である溶銑が排出される。この際に、炉内での原料の温度履歴はその製品の品質に大きな影響を与えるため、炉の高さ方向の炉内温度分布を適正に保つことが望まれる。このとき、適正な温度分布が実現されているかどうかを監視するには、炉の高さ方向に沿って多数点で温度を観測して上下方向の温度分布を実測することが好ましい。しかしながら、経済的な理由などから温度計の設置場所は限られているのが実情であることから、限られた温度センサ情報から炉全体の温度分布を精度良く推定することが望まれている。
そこで近年、炉の限られた位置での温度センサ情報と物理モデルを組み合わせ、物理モデルのモデル化誤差を補償することで、温度分布の推定精度を保つ技術の開発が提案されている。しかし、従来技術は、線形近似に基づいたカルマンフィルタ(非特許文献1)や、集中乗数系近似により簡略化されたモデルに基づいた粒子フィルタ(非特許文献2,3)にとどまっている。また、本来は、非線形、分布乗数系として扱うべきである流動、反応、伝熱現象を含んだ物理モデルの特徴を保ちつつ、温度センサ情報を融合させてプラント(反応炉)の状態推定を行った例はない。
また、高炉などの反応炉では、炉内に固形成分が充填された状態で操業を行うためにプロセス全体の熱容量が大きい。このため、向流熱交換型の反応炉では、アクションに対する応答の時定数が長いという特徴がある。また、炉上部で装入された原料が炉下部に降下するまでには数時間オーダーの無駄時間が存在する。そのため、過去に生じたプロセス外乱の影響が一定時間後になって初めて観測されるという場合が多い。
J. Brannbacka and H. Saxen: ISIJ Int., 41(2001), 1131.
S. Sonoda, N. Murata, H. Hino, H. Kitada, and M. Kano: ISIJ Int., 52(2012), 1086.
N. Kaneko, S. Matsuzaki, M. Ito, H. Oogai, and K. Uchida: ISIJ Int., 50(2010), 939.
上記で述べたように、高炉などの反応炉では、炉内の温度分布を精度良く検出することが反応炉の操業にとって重要であるが、炉内を観測するために設置可能な温度センサの数が限られている。
また、温度分布の乱れの原因となる外乱は、外乱が発生してから所定時間後に初めてその影響が実測される場合が多い。このため、現在の温度センサ情報では現在の炉の状況を精度良く推定出来ない場合も想定されている。
また、温度分布の乱れの原因となる外乱は、外乱が発生してから所定時間後に初めてその影響が実測される場合が多い。このため、現在の温度センサ情報では現在の炉の状況を精度良く推定出来ない場合も想定されている。
本発明は、上記のような点に着目してなされたもので、向流熱交換型の反応炉において限られた数の温度センサ情報から、精度良く炉内温度分布を推定可能な技術を提供することを目的としている。
本発明では、限られた数の温度センサ情報からの温度分布推定であっても、アジョイント法を用いて物理モデルに初期設定される初期温度分布の精度を向上させることで、現在や将来などの初期設定時刻よりも将来の炉内温度分布を推定する。
すなわち、課題を解決するために、本発明の一態様は、向流熱交換型の反応炉内の温度分布を推定するに際し、設定した炉内の温度分布である初期温度分布から、その初期温度分布の設定時刻よりも後の時刻における炉内の温度分布を推定するための物理モデルと、その物理モデルの逆変換モデルであるアジョイントモデルとを有し、且つ炉内の温度を温度センサで順次検出し、設定した過去の設定時刻での初期温度分布を初期値として上記物理モデルを時間順方向に解いて、上記設定した過去の設定時刻よりも後の所定時刻における炉内の温度分布を推定し、検出した検出温度の時系列情報に基づき時間逆方向に上記設定時刻まで上記アジョイントモデルを解いて、設定時刻でのアジョイント変数の値を、上記初期温度分布を補正するための補正情報として求め、その求めた補正情報によって設定した過去の初期温度分布を修正し、修正した後の初期温度分布を初期値として上記物理モデルを時間順方向に解いて上記設定時刻よりも後の時刻での炉内の温度分布を推定する。
すなわち、課題を解決するために、本発明の一態様は、向流熱交換型の反応炉内の温度分布を推定するに際し、設定した炉内の温度分布である初期温度分布から、その初期温度分布の設定時刻よりも後の時刻における炉内の温度分布を推定するための物理モデルと、その物理モデルの逆変換モデルであるアジョイントモデルとを有し、且つ炉内の温度を温度センサで順次検出し、設定した過去の設定時刻での初期温度分布を初期値として上記物理モデルを時間順方向に解いて、上記設定した過去の設定時刻よりも後の所定時刻における炉内の温度分布を推定し、検出した検出温度の時系列情報に基づき時間逆方向に上記設定時刻まで上記アジョイントモデルを解いて、設定時刻でのアジョイント変数の値を、上記初期温度分布を補正するための補正情報として求め、その求めた補正情報によって設定した過去の初期温度分布を修正し、修正した後の初期温度分布を初期値として上記物理モデルを時間順方向に解いて上記設定時刻よりも後の時刻での炉内の温度分布を推定する。
上記のアジョイント法によって初期温度分布を修正し、その修正後の初期温度分布を初期値として物理モデルを解いて炉内の温度分布を再推定する再推定処理は、予め設定した所定回数、若しくは、物理モデルによる解析値と検出温度との差が収束したと評価されるまで繰り返すことが好ましい。
本発明の一態様によれば、アジョイント法を用いて、実測した温度センサ情報の履歴と、炉内の温度分布を推定するための物理モデルの計算値との誤差情報を時間逆方向に伝搬させることにより、誤差が最小化するように初期温度分布を修正する。そして、その修正した初期温度分布を初期値として、現在や将来の温度分布を再推定する。
このように、本発明の一態様によれば、過去一定区間内に測定された温度センサ情報に物理モデルの解析値が合致する方向に初期温度分布を修正することで、初期温度分布を精度良く求めることができる。この結果、限られた温度センサ情報から、現在や将来などの所望時刻での温度分布の推定精度が向上する。
このように、本発明の一態様によれば、過去一定区間内に測定された温度センサ情報に物理モデルの解析値が合致する方向に初期温度分布を修正することで、初期温度分布を精度良く求めることができる。この結果、限られた温度センサ情報から、現在や将来などの所望時刻での温度分布の推定精度が向上する。
以下、本発明に基づく実施形態について図面を参照して説明する。
<構成>
本実施形態では、向流熱交換型のシャフト型の反応炉として高炉を想定して説明する。本発明は、コークス乾式消火設備(CDQ)、ガス化溶融炉などの向流熱交換型の反応炉であっても適用可能である。
<構成>
本実施形態では、向流熱交換型のシャフト型の反応炉として高炉を想定して説明する。本発明は、コークス乾式消火設備(CDQ)、ガス化溶融炉などの向流熱交換型の反応炉であっても適用可能である。
なお、本実施形態では、原料及びガスの移動方向である、炉の上下方向の温度分布を推定する場合を例に説明する。温度分布としては、熱交換する原料(固形分)の温度分布とガスの温度分布を推定することとする。
適用する炉には、図1に示すように、その炉1の温度を測定する複数の温度センサが設定されている。複数の温度センサは、炉1の上下方向に離隔して設定されている。なお、高炉内は、上下方向に沿って向流移動層、融着層、及び液滴下層の複数の層に機能的に区画されるが、複数の温度センサは、例えば向流移動層に配設される。
適用する炉には、図1に示すように、その炉1の温度を測定する複数の温度センサが設定されている。複数の温度センサは、炉1の上下方向に離隔して設定されている。なお、高炉内は、上下方向に沿って向流移動層、融着層、及び液滴下層の複数の層に機能的に区画されるが、複数の温度センサは、例えば向流移動層に配設される。
各温度センサはそれぞれ、設定された計測サンプリングタイム毎に設定位置での温度を検出し、温度分布推定装置2に出力する。
温度分布推定装置2は、図2に示すように、温度記憶部2A、制御部本体2Bと、温度分布推定部2Cと、補正情報演算部2Dと、初期値修正部2Eと、評価部2Fとを備える。
温度分布推定装置2は、図2に示すように、温度記憶部2A、制御部本体2Bと、温度分布推定部2Cと、補正情報演算部2Dと、初期値修正部2Eと、評価部2Fとを備える。
温度記憶部2Aは、入力した検出温度を検出時刻に紐付けて記憶部3に順次記憶する。
制御部本体2Bは、温度分布推定部2C、補正情報演算部2D、初期値修正部2Eの順番に処理を実行させる。制御部本体2Bは、更に、初期値修正部2Eで初期温度分布の修正が完了したら、再度、修正後の初期温度分布を初期値として温度分布推定部2Cを作動させる。その後、制御部本体2Bは、評価部2Fが繰り返し要と判定した場合には、補正情報演算部2D、初期値修正部2E及び温度分布推定部2Cを上述のように作動させる。ここで、温度分布を推定したい時刻が現在時刻よりも将来の時刻であれば、制御部本体2Bは、評価部2Fが繰り返し完了と判定されたときに、分布推定する所定時刻をその求めたい将来の時刻に設定して、修正後の温度分布を初期値として上記温度分布推定部2Cを作動させることで、その将来時刻での温度分布を求める。
制御部本体2Bは、温度分布推定部2C、補正情報演算部2D、初期値修正部2Eの順番に処理を実行させる。制御部本体2Bは、更に、初期値修正部2Eで初期温度分布の修正が完了したら、再度、修正後の初期温度分布を初期値として温度分布推定部2Cを作動させる。その後、制御部本体2Bは、評価部2Fが繰り返し要と判定した場合には、補正情報演算部2D、初期値修正部2E及び温度分布推定部2Cを上述のように作動させる。ここで、温度分布を推定したい時刻が現在時刻よりも将来の時刻であれば、制御部本体2Bは、評価部2Fが繰り返し完了と判定されたときに、分布推定する所定時刻をその求めたい将来の時刻に設定して、修正後の温度分布を初期値として上記温度分布推定部2Cを作動させることで、その将来時刻での温度分布を求める。
温度分布推定部2Cは、設定した過去の設定時刻での初期温度分布を初期値とし、炉1の操業条件に応じたパラメータとに基づき物理モデルを上記の設定時刻から所定時刻まで時間順方向に解いて、初期温度分布の設定時刻よりも後の所定時刻における炉1内の温度分布を推定する処理を行う。本実施形態では所定時刻を現在の時刻とする。所定時刻は、現在の時刻よりも過去の時刻であっても良い。
温度分布推定部2Cで使用する物理モデルは、設定した炉1内の温度分布である初期温度分布を初期値として、その初期温度分布の設定時刻よりも後の時刻における炉1内の温度分布を推定するための物理モデルである。具体的には、物理モデルは、炉1内部の物理収支、運動量収支、及び熱収支を考慮した、上側から装入される燃料を含む固形成分を主要成分とした原料と、下部から吹き込まれる燃焼用空気などのガスとの間の向流熱交換プロセスを表現した向流熱交換モデルからなるモデル計算式である。
この物理モデルに入力される炉1の操業条件としては、例えば原料の投入温度、原料の切り出し速度や、供給されるガスの温度、流量などである。
本実施形態の温度分布推定部2Cでは、初期温度分布の設定時刻を、過去の設定時刻の初期データとして物理モデルに設定して、物理モデルを、時間順方向に現在の時刻まで逐次解くことで、現在の炉1の温度分布を推定する。なお、過去の設定時刻から現在までに時間発展を辿るようにモデル計算するため、その間に行われた炉1の操業条件は決まっている。その操業条件に応じて、物理モデルの各パラメータ値を求めて設定すればよい。すなわち、原料の投入温度、原料の切り出し速度や、供給されるガスの温度、流量はそれぞれ一定値に設定しても良いが、現在までの操業条件は分かっているので、時間の経過に応じて値を設定変更するようにしてもよい。
本実施形態の温度分布推定部2Cでは、初期温度分布の設定時刻を、過去の設定時刻の初期データとして物理モデルに設定して、物理モデルを、時間順方向に現在の時刻まで逐次解くことで、現在の炉1の温度分布を推定する。なお、過去の設定時刻から現在までに時間発展を辿るようにモデル計算するため、その間に行われた炉1の操業条件は決まっている。その操業条件に応じて、物理モデルの各パラメータ値を求めて設定すればよい。すなわち、原料の投入温度、原料の切り出し速度や、供給されるガスの温度、流量はそれぞれ一定値に設定しても良いが、現在までの操業条件は分かっているので、時間の経過に応じて値を設定変更するようにしてもよい。
上記の向流熱交換モデルからなる物理モデルは、公知の物理モデルを使用すればよい。
上下方向の1次元の向流熱交換モデルからなる物理モデルの例を(1)式及び(2)式に示す。
上下方向の1次元の向流熱交換モデルからなる物理モデルの例を(1)式及び(2)式に示す。
本実施形態で使用する物理モデルは、温度センサが検出した温度センサ情報をパラメータに持たないモデル式を採用している。その代わり、本実施形態では、アジョイント法(4次元変分法)によって、記憶部3に記憶されている検出温度のうち、設定時刻から現在の時刻までの検出温度の履歴情報を参照し、アジョイントモデルを用いて時間を遡ることで、物理モデルの解析値と温度センサ情報(実測値)との差を最小にするような設定時刻でのアジョイント変数を求めて、初期温度分布の修正(更新)を行う。
補正情報演算部2Dは、初期温度分布の設定時刻(時刻t=0)から現時刻までに、温度センサが検出した検出温度の時系列情報に基づき、アジョイント法で、現時刻から、初期温度分布を設定した設定時刻まで時間逆方向にアジョイントモデルを解いて、その解である、時刻t0でのアジョイント変数λ(0)を、初期温度分布を補正するための補正情報として求める。
アジョイントモデルは、物理モデルの逆変換モデルである。アジョイントモデルは、公知のアジョイント手法によって、物理モデルを逆変換することで求めることが出来る。アジョイントモデルは、データミスフィットの情報を時間逆方向へ逆伝搬させる役割を有し、t=0のアジョイント変数λ(0)は、現時刻Tのアジョイント変数λ(T)=0を初期値としてアジョイントモデルを時刻t=Tからt=0まで時間逆方向に積分することで求めることが出来る。このアジョイント変数λ(t)は、温度分布の初期値に対する評価関数(モデルによる解析値と実測値との誤差二乗値の積算)の勾配情報を時間逆方向に伝える役割の変数である。またアジョイント変数は、物理モデルに設定される温度分布を表現する要素数と同じ要素数となっている。
初期値修正部2Eは、補正情報演算部2Dが求めた補正情報によって温度分布推定部2Cに設定した上記初期温度分布を修正する。
例えば、現在の初期温度分布x(0)を、上記補正情報をλ(0)とした場合に、下記(3)式によって、現在の初期温度分布値の修正を実施する。
X(0) ← X(0)+α・λ(0) ・・・(3)
ここで、αは、係数であって、例えば0.005〜0.05のような範囲の小さな定数である。このαによって1回の初期値修正による補正情報の反映度合いが設定される。αを小さな値とすることで、徐々に初期温度分布を真値に向けて変更することになる。
例えば、現在の初期温度分布x(0)を、上記補正情報をλ(0)とした場合に、下記(3)式によって、現在の初期温度分布値の修正を実施する。
X(0) ← X(0)+α・λ(0) ・・・(3)
ここで、αは、係数であって、例えば0.005〜0.05のような範囲の小さな定数である。このαによって1回の初期値修正による補正情報の反映度合いが設定される。αを小さな値とすることで、徐々に初期温度分布を真値に向けて変更することになる。
そして、制御部本体2Bは、温度分布推定部2Cを再作動させて、初期値修正部2Eが修正した後の初期温度分布を初期値として物理モデルを再度時間順方向に解いて現在(所定時刻)での炉1内の温度分布を推定する。この処理が温度分布再推定部を構成する。
評価部2Fは、上記の初期温度分布の修正(更新)が予め設定した回数実行された否かを判定し、設定回数未満の場合には、繰り返し要とする。例えば繰り返し数を100回に設定する。この場合には、上記の初期温度分布の修正(更新)し、その修正後の初期温度分布で再度、物理モデルによって現在の温度分布を推定する処理が予め設定した回数だけ繰り返されることになる。
評価部2Fは、上記の初期温度分布の修正(更新)が予め設定した回数実行された否かを判定し、設定回数未満の場合には、繰り返し要とする。例えば繰り返し数を100回に設定する。この場合には、上記の初期温度分布の修正(更新)し、その修正後の初期温度分布で再度、物理モデルによって現在の温度分布を推定する処理が予め設定した回数だけ繰り返されることになる。
または評価部2Fは、予め設定した最大回数(例えば100回)だけ繰り返したか、設定した評価関数による評価に基づき所定誤差まで収束したか否かを判定する。評価関数は、例えば、初期値を設定した設定時刻から現在時刻までの間における予め設定した時間範囲(同化ウインドウ)内にある、複数の温度センサ情報と、その温度センサ情報に対応する物理モデルを解いた際の解析値との乖離度合いを評価する関数とする。そして、評価部2Fは、その乖離度合いが所定位置以下となるまで、若しくは初期温度分布の修正の繰り返しに応じた乖離度の変化の勾配が所定勾配に収まるまで、初期温度分布の修正(更新)が必要と判定する。
ここで、初期温度分布を設定する設定時刻は、過去の任意の時刻を設定することが出来る。設定時刻は、炉の操業途中の時刻であってもよい。
ここで、初期温度分布を設定する設定時刻は、過去の任意の時刻を設定することが出来る。設定時刻は、炉の操業途中の時刻であってもよい。
<動作その他>
以上の温度分布推定装置2によれば、アジョイント法を用いて、実測した温度センサ情報の履歴と、炉1内の温度分布を推定するための物理モデルの解析値との誤差情報を時間逆方向に伝搬させることにより、誤差を最小化するための補正量として設定時刻でのアジョイント変数λ(0)を求める。そしてその補正量に基づき初期温度分布を真値に近づく方向に修正する。この結果、その最適化した過去の初期温度分布を初期値として再度、物理モデルを時間順方向に解くことで、現在や将来などの任意の時刻での温度分布を精度良く推定することができる。
このように、本実施形態によれば、過去一定区間内に測定された温度センサ情報に合致若しくは近似する初期温度分布を精度良く求めることができる。この結果、本実施形態では、限られた温度センサ情報から、現在や将来などの所望で時刻での温度分布の推定精度が向上する。
以上の温度分布推定装置2によれば、アジョイント法を用いて、実測した温度センサ情報の履歴と、炉1内の温度分布を推定するための物理モデルの解析値との誤差情報を時間逆方向に伝搬させることにより、誤差を最小化するための補正量として設定時刻でのアジョイント変数λ(0)を求める。そしてその補正量に基づき初期温度分布を真値に近づく方向に修正する。この結果、その最適化した過去の初期温度分布を初期値として再度、物理モデルを時間順方向に解くことで、現在や将来などの任意の時刻での温度分布を精度良く推定することができる。
このように、本実施形態によれば、過去一定区間内に測定された温度センサ情報に合致若しくは近似する初期温度分布を精度良く求めることができる。この結果、本実施形態では、限られた温度センサ情報から、現在や将来などの所望で時刻での温度分布の推定精度が向上する。
ここで、時間の次元も考慮したデータ同化手法として、アジョイント法の他に、カルマンフィルタを用いたデータ同化手法もある。しかし、上述のように高炉1などの反応炉1にあっては、温度分布の乱れの原因となる外乱は、外乱が発生してから所定の遅れ時間後に初めてその影響が実測される場合が多い。このため、カルマンフィルタを用いた従来のデータ同化手法では、アジョイント法を用いたデータ同化に比べて収束性が悪いと考えられる。
すなわち、本発明を採用することで、目的の時刻での温度分布の推定精度が向上すると考えられる。
また、上記の実施形態では、上下方向の温度分布を推定する場合で説明しているが、上下方向の温度分布推定の共に径方向の温度分布を推定するようにしても良い。
すなわち、本発明を採用することで、目的の時刻での温度分布の推定精度が向上すると考えられる。
また、上記の実施形態では、上下方向の温度分布を推定する場合で説明しているが、上下方向の温度分布推定の共に径方向の温度分布を推定するようにしても良い。
<本実施形態の妥当性について>
本実施形態の妥当性について、図3に示すような、真値シミュレーションを行って確認した。
各モデルには、上述のように向流熱交換モデルからなる物理モデルを設定する。この物理モデルでの向流熱交換モデルは、図1に示すような、炉1の上下方向の1次元方向で向流熱交換を行うモデルとして設定する。そのモデルの計算式として、上述の(1)式及び(2)式を採用する。
本実施形態の妥当性について、図3に示すような、真値シミュレーションを行って確認した。
各モデルには、上述のように向流熱交換モデルからなる物理モデルを設定する。この物理モデルでの向流熱交換モデルは、図1に示すような、炉1の上下方向の1次元方向で向流熱交換を行うモデルとして設定する。そのモデルの計算式として、上述の(1)式及び(2)式を採用する。
図3に記載の3つのモデル4,5,6は同一の物理モデルである。そして、各物理モデル4,5,6に対し、初期値として、過去の設定時刻における所定の初期温度分布を設定して、各物理モデルを時間順方向に逐次解くことで、3つのモデル4,5,6で所定時刻での温度分布を推定するシミュレーションを実施する。
この例では、各温度分布を20要素で表現する場合とする。なお、固体側の固体温度分布とガス側のガス温度分布とがあるため、分布の行列であるx(t)は40要素で表現されている。
この例では、各温度分布を20要素で表現する場合とする。なお、固体側の固体温度分布とガス側のガス温度分布とがあるため、分布の行列であるx(t)は40要素で表現されている。
物理モデル(真値)4には、真値の温度分布を初期温度分布として設定されて、時間順方向に逐次、解いて推定した温度分布を出力する。
物理モデル(補正あり)5及び物理モデル(補正なし)6には、上記の真値の温度分布とは異なる仮定の温度分布を初期温度分布として設定する。物理モデル(補正あり)5及び物理モデル(補正なし)6に設定する初期の初期温度分布は同じ分布値とする。
物理モデル(補正あり)5及び物理モデル(補正なし)6には、上記の真値の温度分布とは異なる仮定の温度分布を初期温度分布として設定する。物理モデル(補正あり)5及び物理モデル(補正なし)6に設定する初期の初期温度分布は同じ分布値とする。
但し、物理モデル(補正あり)5では、物理モデル(真値)4が解いた温度分布のうちの一部の温度値を、温度センサ情報(実測値)と見なす。そして、その逐次計算された一部の温度値を用いて、仮定の初期温度分布(過去値)を補正する補正情報をアジョイント法で求めて、当該仮定の初期温度分布(過去値)の初期値を逐次修正する。この仮定の初期温度分布(過去値)のを修正する処理のことを、データと物理モデルとを同化させるという意味でデータ同化と記す。データ同化は、一般に、物理モデルで表現される数値シミュレーションに実測データを埋め込み、馴染ませることを言う。
次に、この物理モデル(補正あり)5による同化した温度分布が、同化なしの温度分布(物理モデル(補正なし)6による推定値に比較して、真値シミュレーションにおける温度分布(物理モデル(真値)4による推定値)に近づくことを確かめ、本発明に基づく実施形態が妥当性のあることを示す。
ここで、温度分布を左右する操業条件としては、図1に併記した通り、ガス流量、ガス入り側温度、原料降下速度を考慮している。
ここで、温度分布を左右する操業条件としては、図1に併記した通り、ガス流量、ガス入り側温度、原料降下速度を考慮している。
また、温度測定が可能なセンサは図1中の温度センサに示すように2箇所とした。すなわち、以下の解析では、物理モデル(真値)4が解いた温度分布を表現する40要素のうちの2つの要素の値を測定した温度センサ情報と見なした。
物理モデルとして設定する熱分布の計算モデルとしては、上述の通り、微小区間内における熱収支を考慮したモデルである。この例におけるデータ同化とは、センサ位置におけるモデルの解析値と実測値との解離を使用して、過去の初期温度分布の推定を行うことをいう。
物理モデルとして設定する熱分布の計算モデルとしては、上述の通り、微小区間内における熱収支を考慮したモデルである。この例におけるデータ同化とは、センサ位置におけるモデルの解析値と実測値との解離を使用して、過去の初期温度分布の推定を行うことをいう。
次に、アジョイント法によるデータ同化を行う手段について述べる。本実施形態で採用するアジョイント法によるデータ同化の手法は、非特許文献(データ同化:淡路ら著:京大出版)に従った。ここで、以下の説明では、説明を簡易にするために、物理モデルMを線形化しているが、非線形の物理モデルからアジョイントモデルを求めるようにしても良い。
この例におけるアジョイントモデルを使用したデータ同化の処理は、具体的には次の通りである。
まず、物理モデルMを線形化したモデルをM_linとする。また、タイムステップtにおける計算した温度分布をx(t)とする。温度分布x(t)は固体およびガスの温度分布を纏めてベクトル化したものである。すなわち、本例の温度分布x(t)は、40要素の縦ベクトルで表現出来る。すなわち、ガス温度の温度分布が20要素で表現され、固体温度の温度分布が20要素で表現されている。
まず、物理モデルMを線形化したモデルをM_linとする。また、タイムステップtにおける計算した温度分布をx(t)とする。温度分布x(t)は固体およびガスの温度分布を纏めてベクトル化したものである。すなわち、本例の温度分布x(t)は、40要素の縦ベクトルで表現出来る。すなわち、ガス温度の温度分布が20要素で表現され、固体温度の温度分布が20要素で表現されている。
ここで、線形化されたモデルM_linにより、xに微小な擾乱dxが生じた場合、その時間発展は dx(t+1) = M_lin [ dx(t) ] で記述される。
また、M_linの転置をとったものをMa_linと記す。
一方、温度分布の初期値に対する評価関数(物理モデルと実測値との誤差二乗値の積算)の勾配情報を時間逆方向に伝える役割の変数 λを定義する。λはx(t)と同じ要素数である。
また、M_linの転置をとったものをMa_linと記す。
一方、温度分布の初期値に対する評価関数(物理モデルと実測値との誤差二乗値の積算)の勾配情報を時間逆方向に伝える役割の変数 λを定義する。λはx(t)と同じ要素数である。
そして、次のステップの処理によって、アジョイント法によるデータ同化の処理を行う。
(1)ステップ1:
過去の設定時刻における初期温度分布x(0)を初期値として物理モデルに設定する。
(2)ステップ2:
物理モデル x(t+1) = M[x(t)]によって、時間順方向に逐次計算することで、物理量の時間発展を計算する。
(1)ステップ1:
過去の設定時刻における初期温度分布x(0)を初期値として物理モデルに設定する。
(2)ステップ2:
物理モデル x(t+1) = M[x(t)]によって、時間順方向に逐次計算することで、物理量の時間発展を計算する。
(3)ステップ3:
上記の物理モデルの逆変換モデルであるアジョイントモデルによって、時間逆方向に逐次計算することで、初期温度分布を設定した設定時刻t=0におけるアジョイント変数λ(0)を計算する。
ここでのアジョイントモデルは次のように表される。
λ(t−1) = Ma_lin[λ(t)]
+ H’* inv(R) * (H*x(t)− y(t))
上記の物理モデルの逆変換モデルであるアジョイントモデルによって、時間逆方向に逐次計算することで、初期温度分布を設定した設定時刻t=0におけるアジョイント変数λ(0)を計算する。
ここでのアジョイントモデルは次のように表される。
λ(t−1) = Ma_lin[λ(t)]
+ H’* inv(R) * (H*x(t)− y(t))
このとき、アジョイントモデルにおける第二項については、実測値(温度センサ情報)が得られるタイムステップのみについて加える。
y(t)は、時刻tでの物理モデル(真値)4の出力値のうち、実測値とみなした検出温度を有する行列である。
上記の式におけるHとRは、下記のように表現される行列である。本例では、温度分布の要素のうちの7番目と12番目とを実測値(センサ情報に対応)としているため、Hの40×2の行列のうち7行目と12行目にそれぞれ「1」が設定され、その他に「0」を設定している。
y(t)は、時刻tでの物理モデル(真値)4の出力値のうち、実測値とみなした検出温度を有する行列である。
上記の式におけるHとRは、下記のように表現される行列である。本例では、温度分布の要素のうちの7番目と12番目とを実測値(センサ情報に対応)としているため、Hの40×2の行列のうち7行目と12行目にそれぞれ「1」が設定され、その他に「0」を設定している。
(4)ステップ4:
ステップ3で求めたアジョイント変数λ(0)に基づき、初期温度分布x(0)を、下記式のように修正(更新)する。
x(0) ← x(0) + α* λ(0)
ここでは、αは小さな定数で有り、この例では0.01とした。
ステップ3で求めたアジョイント変数λ(0)に基づき、初期温度分布x(0)を、下記式のように修正(更新)する。
x(0) ← x(0) + α* λ(0)
ここでは、αは小さな定数で有り、この例では0.01とした。
(5)ステップ5:
ステップ2〜4をx(0)が予め設定した収束量以下に収束するまで繰り返す。
本例では簡易のため、収束の評価を行うことなく、上記のデータ同化による分布修正を100回繰り返した場合である。
このように物理モデル(補正あり)5においては、このようなアジョイント法を用いて向流熱交換反応炉1における過去の温度分布を修正していく。
ステップ2〜4をx(0)が予め設定した収束量以下に収束するまで繰り返す。
本例では簡易のため、収束の評価を行うことなく、上記のデータ同化による分布修正を100回繰り返した場合である。
このように物理モデル(補正あり)5においては、このようなアジョイント法を用いて向流熱交換反応炉1における過去の温度分布を修正していく。
ここで、物理モデル(真値)4による真値シミュレーションの初期温度分布は、ここでは、図4(a)のように設定した。この真値の初期温度分布が分からない状態で、物理モデル(補正あり5)において初期温度分布が精度良く求めることができるかを実証した。
物理モデル(補正あり)5及び物理モデル(補正なし)6の最初の初期温度分布として、図4(b)の分布を設定した。すなわち、真値の初期温度分布はピークがあるのに対し、仮定として設定した初期の初期温度分布(図4(b))はピークなしの状態から推定を開始した。
物理モデル(補正あり)5及び物理モデル(補正なし)6の最初の初期温度分布として、図4(b)の分布を設定した。すなわち、真値の初期温度分布はピークがあるのに対し、仮定として設定した初期の初期温度分布(図4(b))はピークなしの状態から推定を開始した。
ここで、その他、物理モデルの条件は下記の通りした。
・物理モデルの要素数(温度分布上の温度推定値の数)は20とした。
・シミュレーション開始タイムステップはt=0、最終タイムステップはt=100とした。つまり、現在のタイムステップを100に設定した。
・固体温度分布の20の要素のうち7番目と12番目の2要素のみ、タイムステップ10ごとに観測可能とした。すなわち、温度センサで、全部で10回の観測が行われた場合の例とした。
・物理モデルの要素数(温度分布上の温度推定値の数)は20とした。
・シミュレーション開始タイムステップはt=0、最終タイムステップはt=100とした。つまり、現在のタイムステップを100に設定した。
・固体温度分布の20の要素のうち7番目と12番目の2要素のみ、タイムステップ10ごとに観測可能とした。すなわち、温度センサで、全部で10回の観測が行われた場合の例とした。
以上の部分的な温度測定から初期温度分布をどの程度まで推定できるか、コンピュータによってシミュレーションを実施して確認した。
図4(c)が物理モデル(補正あり)の演算結果である。図4から分かるように、データ同化無しの物理モデル(補正なし)の場合に比べて、物理モデル(補正あり)による推定値である初期温度分布は、真値の図4左に近づくように修正されていることが分かる。
図4(c)が物理モデル(補正あり)の演算結果である。図4から分かるように、データ同化無しの物理モデル(補正なし)の場合に比べて、物理モデル(補正あり)による推定値である初期温度分布は、真値の図4左に近づくように修正されていることが分かる。
またこのときの、3つモデル4,5,6の演算による、現在(タイムステップ=100)での温度分布の推定値を図5に示す。この図5から分かるように、データ同化なし(図5(b))の場合と比較して、データ同化あり(図5(c))の場合は、真値(図5(a))に近い温度分布を演算出来ることが分かる。なお、図4及び5の結果は、上記のステップ2〜4のデータ同化処理を、タイムステップ100の時刻で、100回繰り返した後の収束結果である。実際に収束の評価を行うことで、モデル(補正あり)5では、より真値の温度分布に接近するものと考えられる
また、3つモデル4,5,6による、今回の解析区間(タイムステップ1〜100)における温度分布の推移についても演算してみた。その結果を図6に示す。図6では横軸に時間、縦軸に高さ位置を取り、温度分布の推移の等高線図を示したものである。
また、3つモデル4,5,6による、今回の解析区間(タイムステップ1〜100)における温度分布の推移についても演算してみた。その結果を図6に示す。図6では横軸に時間、縦軸に高さ位置を取り、温度分布の推移の等高線図を示したものである。
ここで、図6及び図7は、1タイムステップ進める度に、上記のデータ処理を100回繰り返した場合の例である。
図6から分かるように、データ同化を行うことで、温度が流れに乗って拡散していく推移も真値に近づくことも確認できた。
また、物理モデル(補正あり)5と物理モデル(補正なし)6の2つモデルによる、誤差(真値と推定値との差異の全高さ位置の平均値)の推移についても演算してみた。その誤差を定量的に示したものが図7である。破線が物理モデル(補正あり)5の場合を示している。
図6から分かるように、データ同化を行うことで、温度が流れに乗って拡散していく推移も真値に近づくことも確認できた。
また、物理モデル(補正あり)5と物理モデル(補正なし)6の2つモデルによる、誤差(真値と推定値との差異の全高さ位置の平均値)の推移についても演算してみた。その誤差を定量的に示したものが図7である。破線が物理モデル(補正あり)5の場合を示している。
図7から分かるように、現在(タイムステップ=100)における誤差のみならず、過去の温度分布についても誤差を低減できていることが分かる。
このように本発明に基づく、アジョイント法によるデータ同化処理を採用して初期温度分布を更新することで、初期値として設定した初期温度分布が真値に近づくように修正される結果、推定する温度分布の精度が向上することが分かる。
このように本発明に基づく、アジョイント法によるデータ同化処理を採用して初期温度分布を更新することで、初期値として設定した初期温度分布が真値に近づくように修正される結果、推定する温度分布の精度が向上することが分かる。
1 炉
2 温度分布推定装置
2A 温度記憶部
2B 制御部本体
2C 温度分布推定部(温度分布再推定部)
2D 補正情報演算部
2E 初期値修正部
2F 評価部
3 記憶部
4,5,6 物理モデル
λ(t) アジョイント変数
2 温度分布推定装置
2A 温度記憶部
2B 制御部本体
2C 温度分布推定部(温度分布再推定部)
2D 補正情報演算部
2E 初期値修正部
2F 評価部
3 記憶部
4,5,6 物理モデル
λ(t) アジョイント変数
Claims (2)
- 向流熱交換型の反応炉内の温度分布を推定する温度分布推定装置であって、
設定した炉内の温度分布である初期温度分布から、その初期温度分布の設定時刻よりも後の時刻における炉内の温度分布を推定するための物理モデルと、上記物理モデルの逆変換モデルであるアジョイントモデルとを有し、
設定した過去の設定時刻での初期温度分布を初期値として上記物理モデルを時間順方向に解いて、上記過去の設定時刻よりも後の所定時刻における炉内の温度分布を推定する温度分布推定部と、
上記設定時刻から上記所定時刻までにおける、上記反応炉内の温度を検出する1又は2以上の温度センサが検出した炉内温度を記憶部に記憶する温度記憶部と、
上記記憶部に記憶されている検出した炉内温度を用いて時間逆方向に上記所定時刻から上記設定時刻まで上記アジョイントモデルを解いて、上記初期温度分布を補正するための補正情報を求める補正情報演算部と、
上記補正情報演算部が求めた補正情報に基づき上記温度分布推定部に設定した上記初期温度分布を修正する初期値修正部と、
上記初期値修正部が修正した後の初期温度分布を初期値として上記物理モデルを時間順方向に解いて炉内の温度分布を推定する温度分布再推定部と、
を備えることを特徴とする温度分布推定装置。 - 向流熱交換型の反応炉内の温度分布を推定する温度分布推定方法であって、
設定した炉内の温度分布である初期温度分布から、その初期温度分布の設定時刻よりも後の時刻における炉内の温度分布を推定するための物理モデルと、その物理モデルの逆変換モデルであるアジョイントモデルとを有し、
炉内の温度を温度センサで順次検出して記憶し、
設定した過去の設定時刻での初期温度分布を初期値として上記物理モデルを時間順方向に解いて、上記設定した過去の設定時刻よりも後の所定時刻における炉内の温度分布を推定し、
上記記憶した検出温度の時系列情報に基づき時間逆方向に上記所定時刻から上記設定時刻まで上記アジョイントモデルを解いて、上記設定時刻でのアジョイント変数の値を、上記初期温度分布を補正するための補正情報として求め、その求めた補正情報に基づき上記設定した初期温度分布を修正し、
上記修正した後の初期温度分布を初期値として上記物理モデルを時間順方向に解いて上記設定時刻よりも後の時刻での炉内の温度分布を推定することを特徴とする温度分布推定方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2016152212A JP2018021693A (ja) | 2016-08-02 | 2016-08-02 | 温度分布推定装置および温度分布推定方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2016152212A JP2018021693A (ja) | 2016-08-02 | 2016-08-02 | 温度分布推定装置および温度分布推定方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
JP2018021693A true JP2018021693A (ja) | 2018-02-08 |
Family
ID=61164904
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2016152212A Pending JP2018021693A (ja) | 2016-08-02 | 2016-08-02 | 温度分布推定装置および温度分布推定方法 |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP2018021693A (ja) |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH1161284A (ja) * | 1997-08-15 | 1999-03-05 | Nkk Corp | 焼結鉱の還元粉化性の評価試験方法 |
WO2007122677A1 (ja) * | 2006-04-13 | 2007-11-01 | Osaka University | 熱対流場・物質拡散場の設計支援方法、設計支援システム、および設計支援用プログラム |
US20140278314A1 (en) * | 2013-03-13 | 2014-09-18 | The Government Of The United States Of America, As Represented By The Secretary Of The Navy | System and method for correcting a model-derived vertical structure of ocean temperature and ocean salinity based on velocity observations |
WO2016051292A1 (en) * | 2014-10-03 | 2016-04-07 | International Business Machines Corporation | Tunable miniaturized physical subsurface model for simulation and inversion |
-
2016
- 2016-08-02 JP JP2016152212A patent/JP2018021693A/ja active Pending
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH1161284A (ja) * | 1997-08-15 | 1999-03-05 | Nkk Corp | 焼結鉱の還元粉化性の評価試験方法 |
WO2007122677A1 (ja) * | 2006-04-13 | 2007-11-01 | Osaka University | 熱対流場・物質拡散場の設計支援方法、設計支援システム、および設計支援用プログラム |
US20140278314A1 (en) * | 2013-03-13 | 2014-09-18 | The Government Of The United States Of America, As Represented By The Secretary Of The Navy | System and method for correcting a model-derived vertical structure of ocean temperature and ocean salinity based on velocity observations |
WO2016051292A1 (en) * | 2014-10-03 | 2016-04-07 | International Business Machines Corporation | Tunable miniaturized physical subsurface model for simulation and inversion |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US20190019096A1 (en) | Estimator, estimation method, program and storage medium where program stored for model parameter estimation and model parameter estimation system | |
CN109583585B (zh) | 一种电站锅炉壁温预测神经网络模型的构建方法 | |
CN107526927B (zh) | 一种高炉铁水质量在线鲁棒软测量方法 | |
CN109935280B (zh) | 一种基于集成学习的高炉铁水质量预测系统及方法 | |
CN102176221B (zh) | 基于动态工况的焦炉加热燃烧过程焦炉温度预测方法 | |
CN101211383A (zh) | 一种高炉铁水硅含量的特征分析预报方法 | |
JP7027536B2 (ja) | 解析システム及び解析方法 | |
CN103048058B (zh) | 一种焦炉火道温度在线检测方法 | |
CN107299170B (zh) | 一种高炉铁水质量鲁棒软测量方法 | |
CN106709197A (zh) | 基于滑动窗口t‑s模糊神经网络模型的铁水硅含量预测方法 | |
CN113821984A (zh) | 一种基于时域卷积模型的加热炉钢坯温度计算方法 | |
JP6024718B2 (ja) | 高炉炉熱予測装置及び高炉炉熱予測方法 | |
JP6003909B2 (ja) | 高炉通気性予測装置及び高炉通気性予測方法 | |
CN110569567B (zh) | 焦炉火道温度的测量方法及系统 | |
JP6834295B2 (ja) | 状態予測装置 | |
JP5854171B2 (ja) | 補正装置、補正方法及び鉄鋼精錬方法 | |
JP5900026B2 (ja) | 炉内温度分布の推定方法および推定装置 | |
CN117930780A (zh) | 石英产品生产质量控制方法及系统 | |
JP6094127B2 (ja) | 温度分布推定方法及び温度分布推定装置 | |
JP5900386B2 (ja) | 乾留炉の制御方法および制御装置 | |
JP2018021693A (ja) | 温度分布推定装置および温度分布推定方法 | |
JP5900025B2 (ja) | 炉内温度分布の推定方法および推定装置 | |
CN115164115A (zh) | 基于物理模型驱动机器学习的燃气管道泄漏速率预测方法 | |
CN103245481B (zh) | 基于变频技术的变负荷大型换热器气阻特性的检测方法 | |
CN113609739A (zh) | 材料热处理工艺与微观组织、性能关系数据库的构建方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20180323 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20190108 |
|
A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20190227 |
|
A02 | Decision of refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A02 Effective date: 20190625 |