JP2010019105A - 制御装置 - Google Patents

Abstract

【課題】むだ時間の変化を逐次推定して高精度な制御を実施可能にしつつ、プラントモデルのパラメータが同定することに伴い振動又は発散することを抑制した制御装置を提供する。
【解決手段】所定のサンプリング周期ｄｔで離散化された数式により制御対象を表したプラントモデルのパラメータを逐次同定し、同定された同定パラメータを用いてプラントモデルに含まれている余むだ時間を推定し、推定した余むだ時間に基づきむだサンプリング回数の値を繰り上げ／繰り下げ更新し、更新されたむだサンプリング回数に基づきむだ時間を推定する。そして、同定手段は、前記同定パラメータを存在領域内に制限しつつ実行される。
【選択図】 図１７

Description

本発明は、内燃機関等の制御対象を数式で表したプラントモデルを用い、該プラントモデルのパラメータを逐次同定することで制御対象の動特性を自動的に適応させるようにした制御装置に関する発明である。

従来から、制御対象の動特性を逐次同定する制御装置が提案されている。この制御装置は、制御対象を数式で表したプラントモデルを持ち、実際の制御入力をプラントモデルに入力した時のプラントモデル出力と実際の制御対象の出力との誤差をゼロに近づけるようにプラントモデルのパラメータを逐次推定するものであるが、制御対象がむだ時間を含む系の場合、同定に用いるプラントモデルは、実際の制御対象が持つむだ時間を考慮する必要がある。そこで特許文献１，２では、以下の同定手法によりむだ時間の推定を可能にしている。

すなわち、所定のサンプリング周期（ｄｔ）で離散化された数式により制御対象を表したプラントモデルについて、そのプラントモデルのむだ時間（Ｌ）をサンプリング周期で除算した商をむだサンプリング回数（ｄ）、余りを余むだ時間（Ｌ１）とする。そして、逐次同定されたパラメータを用いれば、プラントモデルに含まれている余むだ時間を推定することができることを利用して、以下に説明する図８のタイムチャートに例示される如くむだ時間を推定する。

図８（ａ）に示すように、ある制御入力に対するノミナルむだ時間と実むだ時間との間に図示のような誤差が生じている場合を想定する。なお、サンプリングタイミング（ｄ・ｄｔ）を超えた部分が余むだ時間Ｌ１であり、ノミナルむだ時間は、予め設定されたむだ時間のノミナル値である。

この場合、同定が逐次行われることにより、図８（ｂ）に示すように、余むだ時間の推定値Ｌ１＿hat（以下、ｘ＿hatはｘの同定値もしくは推定値を表す）が変化し、ｄｔ近傍まで増加する。そして、余むだ時間の推定値Ｌ１＿hatがｄｔ近傍で所定時間以上とどまっている場合には、図８（ｃ）に示すように、むだサンプリング回数の推定値ｄ＿hatが１つ繰り上げられ（ｄ＿hat←ｄ＋１）、更に、余むだ時間の推定値Ｌ１＿hatがゼロ近傍の値に更新される。その後、同定が逐次行われることにより、図８（ｄ）に示すように、余むだ時間の推定値Ｌ１＿hatが増加し、それに伴いむだ時間の推定値が実むだ時間に収束する。

このようにして更新されたむだサンプリング回数ｄ＿hat及び余むだ時間Ｌ１＿hatを用いれば、Ｌ（むだ時間）＝ｄ（むだサンプリング回数）×ｄｔ（サンプリング周期）＋Ｌ１（余むだ時間）、との算出式に基づきむだ時間Ｌを推定することができる。

なお、図８では、余むだ時間の推定値Ｌ１＿hatが増加する場合を例示し、それに伴い、むだサンプリング回数の推定値ｄ＿hatが１つ繰り上げられる様子を説明したが、これとは逆に、余むだ時間の推定値Ｌ１＿hatが減少してゼロ近傍で所定時間以上とどまっている場合には、むだサンプリング回数の推定値ｄ＿hatが１つ繰り下げられると共に、余むだ時間の推定値Ｌ１＿hatがｄｔ近傍の値に更新される。
特開２００６−１１８４２９号公報 特開２００７−２８５１６７号公報

しかしながら、上述の如く推定したむだ時間と実際のむだ時間との差が大きいと、プラントモデルのパラメータ（例えば応答遅れを示すパラメータや定常ゲイン）を逐次同定するにあたり、パラメータの同定値が振動又は発散するおそれのあることが分かった。

以下、図１９を用いて上記発散及び振動について説明する。図１９は、縦軸をむだ時間推定値とし、横軸をプラントモデルのパラメータ（ここでは時定数をパラメータの一例として説明）とした図であり、むだ時間推定値及び時定数が、逐次同定されることにより真値に向けて変化する際の挙動を示す。また、図中右側の（ａ）（ｂ）（ｃ）欄はプラントモデルへのステップ入力に対する出力を縦軸とし、横軸を経過時間とした場合における前記出力の時間変化を示す図であり、入力に対する出力の応答特性を表している。そして、（ａ）（ｂ）（ｃ）欄中の点線は、制御対象による実際の出力変化特性を示し、一点鎖線は、時定数について同定を開始した時点におけるプラントモデルの出力変化特性を示し、実線は、時定数について逐次同定を行った後における出力変化特性を示す。

（ｃ）欄の一点鎖線に示すように推定むだ時間と実際のむだ時間との差が小さい場合には、時定数の値は逐次同定により少し小さくなるよう補正され、その結果、実線に示すようにプラントモデルの出力応答特性は真値に近づくこととなる。

しかしながら、（ｂ）欄の一点鎖線に示すように推定むだ時間と実際のむだ時間との差が大きい場合には、真値に近づくように時定数の値を小さく補正した結果、実線に示すようにプラントモデルの出力応答特性は過剰に立てられ、真値に対してオーバーシュートする（つまり振動する）特性となってしまう。これは、同定値（時定数及びむだ時間）の同定挙動も振動的になってしまうことを意味する。さらに（ａ）欄の一点鎖線に示すように推定むだ時間と実際のむだ時間との差が更に大きい場合には、真値に近づくように時定数の値を小さく補正した結果、実線に示すようにプラントモデルの出力応答特性はさらに過剰に立てられ、出力値が無限大となるよう発散する特性となってしまう。これは、同定値（時定数及びむだ時間）の同定挙動も発散してしまうことを意味する。

本発明は、上記課題を解決するためになされたものであり、その目的は、むだ時間の変化を逐次推定して高精度な制御を実施可能にしつつ、プラントモデルのパラメータが同定することに伴い振動又は発散することを抑制した制御装置を提供することにある。

以下、上記課題を解決するための手段、及びその作用効果について記載する。

請求項１記載の発明では、
所定のサンプリング周期（ｄｔ）で離散化された数式により制御対象を表したプラントモデルを用い、前記制御対象への入力を前記プラントモデルに加えた際のプラントモデル出力と前記制御対象の実出力との差である同定誤差をゼロに近づけるように、前記プラントモデルのパラメータを逐次同定する同定手段と、
前記制御対象のむだ時間（Ｌ）を前記サンプリング周期で除算した商をむだサンプリング回数（ｄ）、余りを余むだ時間（Ｌ１）としてモデル化した場合に、前記同定手段により同定された同定パラメータを用いて前記プラントモデルに含まれている前記余むだ時間を推定する余むだ時間推定手段と、
前記余むだ時間推定手段により推定された推定余むだ時間が、所定の増加条件を満たす程度に増加した場合には前記むだサンプリング回数の推定値を一つ繰り上げるよう更新し、所定の減少条件を満たす程度に減少した場合には前記むだサンプリング回数の推定値を一つ繰り下げるよう更新する更新手段と、
前記むだサンプリング回数更新手段により更新された前記むだサンプリング回数の推定値に基づき前記むだ時間を推定するむだ時間推定手段と、
を備え、
前記同定手段は、前記同定パラメータの値が予め設定された制約領域内に制限されつつ実行されることを特徴とする。

これによれば、離散化されたプラントモデルのパラメータを逐次同定し（同定手段）、同定されたパラメータを用いてプラントモデルに含まれている余むだ時間を推定し（余むだ時間推定手段）、推定された余むだ時間がサンプリング周期近傍へ変化した場合にはむだサンプリング回数の推定値を一つ繰り上げるよう更新し、ゼロ近傍へ変化した場合には一つ繰り下げるよう更新し（むだサンプリング回数更新手段）、更新されたむだサンプリング回数推定値に基づきむだ時間を推定する（むだ時間推定手段）。よって、むだ時間の変化を逐次推定して高精度な制御を実施できる。

しかも、「前記同定手段は、前記同定パラメータの値が予め設定された制約領域内に制限されつつ実行される」との構成を備えるので、図１９に例示される如くパラメータ（図１９の例では時定数）が振動又は発散してしまうことを抑制できる。

請求項２記載の発明では、前記同定パラメータは、射影アルゴリズムにより設定された制約領域内に制限されることを特徴とする。このように射影アルゴリズムを用いれば、パラメータの値を予め設定された制約領域内に制限しつつ逐次同定することを容易に実現でき、好適である。

射影アルゴリズムとは、例えば「現代制御シリーズ７.アダプティブコントロール（コロナ社）」において公知のアルゴリズムであるが、図１８を参照しながら以下にその概略を説明する。可変パラメータの存在領域について事前知識がある場合において、そのパラメータを存在領域に閉じ込めたままで同定処理するものである。例えば、存在領域の境界上における更新ベクトルが、存在領域から外れる向きのベクトルである場合には、その更新ベクトルを境界上の接線方向、又はその接線方向よりも存在領域側に向かう向きのベクトル（つまり射影ベクトル）に補正して、逐次同定による更新を行う。これにより、更新ベクトルが存在領域から外れない向きのベクトルに制限され、その結果、パラメータは存在領域内に強制的に留められることとなる。

請求項３記載の発明では、前記同定パラメータを制限する前記制約領域は、前記プラントモデルの時定数に相当する値が所定範囲内となるように設定されることを特徴とする。

このように時定数に相当する値が所定範囲内となるように、同定パラメータの制約領域を設定することにより、図１９に例示されるようなパラメータ振動又は発散を抑止でき好適である。ちなみに、前記所定範囲は例えば、プラントモデルが安定系となるように設定するとよい（請求項４）。また例えば、実際の制御対象の応答特性が過制動であると予め分かっている場合であれば、オーバーシュートしない安定範囲に設定するとなおよい。また、実際の制御対象のばらつき，経年劣化も含めた応答特性の取りうる範囲が予め分かっている場合であれば、前記所定範囲をその範囲に設定するとよい。ちなみに、「プラントモデルが安定系」とは、有界な大きさの任意の入力に対してその出力も有界であることを意味する。

請求項５記載の発明では、前記同定パラメータを制限する前記制約領域は、前記プラントモデルの定常ゲインが、前記制御対象の定常ゲインと符号が同一で所定範囲内となるように設定されることを特徴とする。これによれば、例えば実際の制御対象の定常ゲインが正である場合に、同定手段が誤って負であると同定してしまうことを回避しながら同定でき、推定むだ時間が振動又は発散することを抑制でき好適である。また、実際の制御対象のばらつき，経年劣化も含めた定常ゲインの取りうる範囲が予め分かっている場合であれば、前記所定範囲をその範囲に設定するとなおよい。

請求項６記載の発明では、前記同定パラメータを制限する前記制約領域は、前記プラントモデルのすべての極が安定であるように設定されることを特徴とする。これによれば、制御対象を安定系として捉えている場合において、同定手段が誤って不安定系であると同定してしまうことを回避しながら同定でき、推定むだ時間が振動又は発散することを抑制でき好適である。ちなみに、「極が安定である」とは、例えば「システム工学シリーズ３.フィードバック制御入門（コロナ社）」において定義されているように、例えば離散系においては、伝達関数の分母多項式の根の絶対値が１より小さいことを意味する。

請求項７記載の発明では、
前記制御対象は、内燃機関の排気通路に配設された排ガス浄化用の触媒と、前記排気通路において前記触媒の上流側で特定排ガス濃度を検出する上流側排ガスセンサと、前記排気通路において前記触媒の下流側で特定排ガス濃度を検出する下流側排ガスセンサと、を有する空燃比制御システムであり、
前記上流側排ガスセンサで検出した空燃比を前記触媒の上流側目標空燃比に一致させるようにフィードバックして、供給空燃比を前記制御対象への制御入力として演算するメインフィードバック制御手段と、
前記下流側排ガスセンサの検出した空燃比を前記触媒の下流側目標空燃比に一致させるようにフィードバックして、前記上流側目標空燃比を前記制御入力の算出に用いる値として演算するサブフィードバック制御手段と、
を備えることを特徴とする。

以下、本発明を内燃機関の空燃比制御システムに適用して具体化した一実施形態について図面を用いて説明する。

まず、図１に基づいてエンジン制御システム全体の概略構成を説明する。内燃機関であるエンジン１１の吸気管１２の最上流部には、エアクリーナ１３が設けられ、このエアクリーナ１３の下流側に、吸入空気量を検出するエアフローメータ１４が設けられている。このエアフローメータ１４の下流側には、モータ１５によって開度調節されるスロットルバルブ１６と、このスロットルバルブ１６の開度（スロットル開度）を検出するスロットル開度センサ１７とが設けられている。

更に、スロットルバルブ１６の下流側には、サージタンク１８が設けられ、このサージタンク１８には、吸気管圧力を検出する吸気管圧力センサ１９が設けられている。また、サージタンク１８には、エンジン１１の各気筒に空気を導入する吸気マニホールド２０が設けられ、各気筒の吸気マニホールド２０の吸気ポート近傍に、それぞれ燃料を噴射する燃料噴射弁２１が取り付けられている。また、エンジン１１のシリンダヘッドには、各気筒毎に点火プラグ２２が取り付けられ、各点火プラグ２２の火花放電によって各気筒内の混合気に着火される。

また、エンジン１１のシリンダブロックには、冷却水温を検出する冷却水温センサ２６や、エンジン１１のクランク軸２７が所定クランク角回転する毎にクランク角信号（パルス信号）を出力するクランク角センサ２８が取り付けられている。このクランク角センサ２８のクランク角信号に基づいてクランク角やエンジン回転速度が検出される。

一方、エンジン１１の排気管２３には、排気ガスを浄化するための２つの触媒２５，３０が直列に設けられている。各触媒２５，３０は、例えば、三元触媒、ＮＯｘ吸蔵型三元触媒等により構成され、上流側の触媒２５の上流側と下流側には、それぞれ特定排気ガス濃度（例えば酸素濃度、空燃比等）を検出する上流側排気ガスセンサ３１と下流側排気ガスセンサ３２とが設けられている。本実施形態では、上流側排気ガスセンサ３１として空燃比センサを用い、下流側排気ガスセンサ３２として酸素センサ（Ｏ_２センサ）を用いているが、この構成に限定されないことは言うまでもない。

これら各種センサの出力は、エンジン制御回路（以下「ＥＣＵ」と表記する）２９に入力される。このＥＣＵ２９は、マイクロコンピュータを主体として構成され、内蔵されたＲＯＭ（記憶媒体）に記憶された各種のエンジン制御プログラムを実行することで、エンジン運転状態に応じて燃料噴射弁２１の燃料噴射量や点火プラグ２２の点火時期を制御する。

また、ＥＣＵ２９は、図２に示すように、上流側排気ガスセンサ３１で検出した空燃比を上流側触媒２５上流側の目標空燃比に一致させるように供給空燃比（燃料噴射量）をフィードバック補正するメインフィードバック制御（以下の説明では「フィードバック」を「Ｆ／Ｂ」で表記する）を行うメインＦ／ＢコントローラＭ１１と、下流側排気ガスセンサ３２の検出電圧（検出空燃比）を目標電圧（上流側触媒２５下流側の目標空燃比）に一致させるように上流側触媒２５上流側の目標空燃比をＦ／Ｂ補正するサブＦ／Ｂ制御を行うサブＦ／ＢコントローラＭ２１として機能する。

メインＦ／ＢコントローラＭ１１は、事前にモデル化して適合したプラントモデル（１次遅れ系＋むだ時間のモデル）を基に設計されたものであり、このメインＦ／ＢコントローラＭ１１により上流側排気ガスセンサ３１の検出空燃比と上流側触媒２５上流側の目標空燃比との偏差を小さくするように空燃比補正係数が演算され、最適な空燃比Ｆ／Ｂ制御が実現される。

一方、サブＦ／ＢコントローラＭ２１は、下流側排気ガスセンサ３２の検出電圧（検出空燃比）を目標電圧（上流側触媒２５下流側の目標空燃比）に一致させるように上流側触媒２５上流側の目標空燃比をＦ／Ｂ補正するサブＦ／Ｂ制御を行うように、事前にモデル化して適合したプラントモデル（２次遅れ系＋むだ時間のモデル）を基に設計されている。

しかしながら、このサブＦ／ＢコントローラＭ２１の実際の制御対象（上流側触媒２５、下流側排気ガスセンサ３２）の個体差や劣化等によってサブＦ／Ｂ制御誤差が生じる。そこで、本実施形態では、適応制御と称される制御方式を用い、サブＦ／ＢコントローラＭ２１におけるＦ／Ｂゲインを制御対象（プラント）の現時点の動特性に自動的に適応させ、制御系の性能を常に最良の状態に保持するようにしている。すなわち、制御対象を数式モデルで表したプラントモデルを用い、制御対象の入力（上流側触媒２５の流入過不足酸素量）を該プラントモデルに入力した時のプラントモデル出力と制御対象の出力（下流側排気ガスセンサ３２の検出電圧）との誤差ｅをゼロに近づけるようにプラントモデルの可変パラメータを逐次同定する（この機能が特許請求の範囲でいう同定手段に相当する）。

ここで、プラントモデルは、離散時間で表した離散プラントモデルであり、むだ時間を持つ２次遅れ系で近似すると共に、図３に示すように、２次遅れ系は、２つの１次遅れ系に分割して離散化した後結合して近似するようにしている。

プラントモデル（上流側触媒２５＋下流側排気ガスセンサ３２のモデル）の離散化は、オフラインで以下のように行う。本実施形態では、プラントモデルをむだ時間を持つ２次遅れ系で近似する。

連続プラントモデルは、出力をＹ（ｓ）、入力をＵ（ｓ）として次式で表される。

ここで、ω、ζ、Ｋ、Ｌは連続モデルパラメータで、ωは固有角振動数、ζは減衰係数、Ｋは定常ゲイン、Ｌはむだ時間である。

また、むだ時間Ｌは、次式で表される。

Ｌ＝ｄ×ｄｔ＋Ｌ１
ここで、ｄはむだサンプリング回数、ｄｔはサンプリング周期、Ｌ１は余むだ時間である。つまり、むだ時間Ｌを離散化のサンプリング周期ｄｔで除算した商ｄをむだサンプリング回数とし、余りＬ１を余むだ時間とするものである。むだ時間Ｌをサンプリング周期ｄｔで割り切れない場合は、余むだ時間Ｌ１は、後述するむだ時間更新処理により０〜ｄｔの正の値となる。

上記［数１］式で表される連続プラントモデルの伝達関数Ｇ（ｓ）は次式で表される。

この連続プラントモデルを離散化（拡張ｚ変換）する際に、部分分数展開を利用して、次のように１次遅れ系の和に分解する。ここで、Z［＊］は＊の拡張ｚ変換を表す。拡張ｚ変換は、例えば「ディジタル制御システム−解析と設計（日刊工業新聞社）」において公知である。

従って、離散プラントモデルは、次式で表される。

ここで、ｐ1，ｐ2，q1，q2，q3は離散モデルパラメータであり，q1，q2，q3の中に余むだ時間Ｌ１が含まれる。

以上説明した離散化式をまとめると、図４のように表される。

この離散モデルパラメータを連続化すると、図５に示すように、連続モデルパラメータω、ζ、Ｋ、Ｌ１が導き出される。図５中のＺｒｅ，Ｚｉｍは図６に示す算出式により求められる。

上記各式に関して、実際はオンボードで演算可能にするため近似演算を用いてもよい。例えば三角関数、指数関数、対数関数は、２次までテイラー展開し、誤差大の範囲はテーブルでもつものとする。

さて、図４に示すように、拡張ｚ変換により離散化すると、離散モデルパラメータq1，q2，q3に余むだ時間Ｌ１の情報（これが「むだ時間の一部の情報」に相当）が含まれる。このことが拡張ｚ変換の特徴であり、一般的なｚ変換による離散化では現れないものである。同定システムにより同定するのは、離散モデルパラメータθ（θ＝［ｐ1，ｐ2，q1，q2，q3］^Ｔ:上付きＴは転置を表す）であり、同定値θ＿hat（ｘ＿hatはｘの同定値もしくは推定値を表すものとする、以下同様）を用いて、図５の連続化式により連続モデルパラメータの推定値ω＿hat，ζ＿hat，Ｋ＿hat，Ｌ１＿hatを算出できるようになっている。この機能が特許請求の範囲でいう余むだ時間推定手段に相当する。そして、この余むだ時間の推定値Ｌ１＿hatを用いることで、むだ時間Ｌ（むだサンプリング回数ｄ）を推定することが可能となる。この機能が特許請求の範囲でいうむだ時間推定手段に相当する。

以下、離散モデルパラメータを同定する同定処理について説明する。本同定処理では、同定誤差ｅ（制御対象の実出力と離散プラントモデル出力との偏差）をゼロにするように離散モデルパラメータを逐次最小二乗法により推定する。以下、図２を基に詳細に説明する。

離散プラントモデルは上記［数４］式で表される数式モデルであり、離散モデルパラメータθ（θ＝［ｐ1，ｐ2，q1，q2，q3］^Ｔ）は適応機構により逐次推定される。

図２に示すように、適応機構Ｍ３２は、離散プラントモデルの予測出力と制御対象の出力（下流側排気ガスセンサ３２の検出電圧）との同定誤差ｅをゼロに近づけるように離散モデルパラメータθを推定するものである。同定誤差ｅは、離散プラントモデルに実際の制御入力ｕを入力した時のモデル出力ｙ＿hatと制御対象の実出力ｙとの誤差である。尚、離散プラントモデルに入力される制御入力ｕや制御対象の実出力ｙ（下流側排気ガスセンサ３２の検出電圧）には、不要な直流成分やノイズ成分などが含まれる。それ故に、制御入力ｕは、低周波成分除去手段としてのＨＰＦ（ハイパスフィルタ）により直流成分が除去される。

また、制御対象の実出力ｙにおいて、Ｆ／Ｂ制御用の出力ｙは、高周波成分除去手段としてのＬＰＦ（ローパスフィルタ）によりノイズ成分が除去され、パラメータ同定用の出力ｙは、ＨＰＦとＬＰＦとを組み合わせた（バンドパスフィルタ）により直流成分とノイズ成分とが除去される。このとき、空燃比Ｆ／Ｂ用の出力ｙと同定用の出力ｙとには各々別のＬＰＦを設定できるようになっている。

また、同定誤差ｅには必要に応じてフィルタ処理（ＬＰＦ）や不感帯処理が施される。このフィルタ処理により、同定パラメータの振動が抑制され、不感帯処理により過同定が抑制される。尚、このフィルタ処理や不感帯処理は無くても良い。

ここで、前記［数４］式を、入力にむだ時間を有する自己回帰モデルとして表現しなおすと、次のように表される。

上式を線形パラメトリックな自己回帰モデルとして表現しなおすと、次のように表される。

ここで、θ＿hat(k-1)＝［ｐ1＿hat(k-1)，ｐ2＿hat(k-1)，ｑ１＿hat(k-1)，ｑ２＿hat(k-1)，ｑ３＿hat(k-1)］^Ｔ，ζ(k-1)＝［ｙ(k-1)，ｙ(k-2)，ｕ(k-d)，ｕ(k-d-1)，ｕ(k-d-2)］^Ｔである。同定誤差ｅは次式で求められる。

ここで，同定誤差eに対してノイズ除去や不感帯処理を施すことにより，同定パラメータの振動や過同定を避ける構成としてもよい．
次に、適応機構Ｍ３２において、同定誤差ｅをゼロにするように離散モデルパラメータの推定値θ＿hatを算出するパラメータ調整則は、本実施形態では重み付き最小二乗法の原理に基づいて導出される。次式に示す同定誤差ｅの２乗和を評価関数として考える。

ここで、λは忘却係数とも呼ばれる重み係数である。上記評価関数が最小となるようなパラメータ調整則は次式のように与えられる。

上式において、同定誤差ｅに乗算する行列ゲインΓ及びスカラゲインγは次式のようになる。

そして、上式で求められた離散モデルパラメータの推定値が離散プラントモデルに反映される。

尚、パラメータ調整則は上記以外に、例えば、固定ゲイン則、漸減ゲイン則、固定忘却ゲイン則，可変ゲイン則、固定トレースゲイン則，上限トレースゲイン則，下限トレースゲイン則，両限トレースゲイン則等を用いても良い。

適応機構Ｍ３２で同定された離散モデルパラメータは、学習判定処理により同定誤差ｅが所定範囲内であるか否かが判定され、同定誤差ｅが所定範囲内であれば、むだ時間Ｌの学習が許可される。

同定された離散モデルパラメータｐ1，ｐ2，ｑ1，ｑ2，ｑ3は、図５に示す連続化式により連続化され、連続モデルパラメータω、ζ、Ｋ、Ｌ１に変換される。そして、後述する方法で、余むだ時間Ｌ１からむだ時間Ｌが学習される。この学習値は、ＥＣＵ２９のバックアップＲＡＭ等の書き換え可能な不揮発性メモリに記憶され、電源オフ中でも、学習値の記憶を保持できるようになっている。

サブＦ／ＢコントローラＭ２１に対しては、連続モデルパラメータω、ζ、Ｋ、Ｌ１の学習値が離散化されて入力される。サブＦ／ＢコントローラＭ２１は、学習値から算出された離散モデルパラメータを用いて制御対象の入力（入力酸素量）を算出する。この際、学習値（同定値）は逐次反映しても良いし、学習値の収束後に反映しても良い。また、学習値の範囲によって、予め算出した離散モデルパラメータ又は制御ゲインそのものを切り替えるようにしても良い。

次に、むだ時間更新処理について説明する。むだ時間Ｌを持つ連続系は、次式で表される。ただし，説明の簡略化のため１次遅れ系を例にとり説明するが，高次の系においても考え方は同様である．

そして、上式を拡張ｚ変換により離散化することで、次式が得られる。

但し、上記［数１１］式、［数１２］式において、ｙ（ｋ）＝ｘ（ｋ）である。

ここで、図７は、ｋｄｔ≦ｔ＜（ｋ＋１）ｄｔにおける出力ｙ（ｋ＋１）に影響を与える制御入力ｕ（むだ時間分だけ正の方向にずらした制御入力ｕ）の変化を表すタイムチャートである。図７において、むだ時間Ｌがサンプリング周期ｄｔで割り切れない場合、入力ｕ（τ−Ｌ）は、τ＝ｋｄｔ〜（ｋ＋１）ｄｔの期間内で一度だけ値が変わり、当該期間では、ｄ回前の入力ｕの影響とｄ＋１回前の入力ｕの影響を受ける。なお、上記［数１２］式では、右辺の第２項によりｄ回前の入力ｕの影響を表し、同第３項によりｄ＋１回前の入力ｕの影響を表している。

この場合、制御対象の実むだ時間は未知であるが、拡張ｚ変換により離散化した離散モデルパラメータには余むだ時間Ｌ１の情報が含まれており、この余むだ時間Ｌ１を用いることによりむだ時間Ｌの推定が可能となる。本実施形態では、離散モデルパラメータの推定値θ＿hatの連続化により算出される余むだ時間の推定値Ｌ１＿hatの変化に基づき、離散プラントモデルのむだ時間を制御対象の実むだ時間に近づけるように、むだ時間（むだサンプリング回数ｄ）を更新し、更新したむだ時間（むだサンプリング回数ｄ）を離散プラントモデルでのプラントモデル出力の算出に反映させるようにしている。

より具体的には、余むだ時間の推定値Ｌ１＿hatが、サンプリング周期ｄｔで規定される上限値近傍又は下限値近傍にあるか否かを判定する。ここで、推定値Ｌ１＿hatの上限値はｄｔ、下限値は０であり、微小な正の定数ε2により上限値近傍を「ｄｔ−ε2〜ｄｔの範囲」、下限値近傍を「０〜０＋ε2の範囲」として設定する。この場合、余むだ時間の推定値Ｌ１＿hatが上限値近傍（ｄｔ−ε2〜ｄｔ）に変化してその状態が所定時間継続した場合に、むだサンプリング回数ｄを１つ繰り上げると共に、該余むだ時間の推定値Ｌ１＿hatを下限値近傍の所定値（但し下限値＝０以上）とする。また、余むだ時間の推定値Ｌ１＿hatが下限値近傍（０〜０＋ε2）に変化してその状態が所定時間継続した場合に、むだサンプリング回数ｄを１つ繰り下げると共に、該余むだ時間の推定値Ｌ１＿hatを上限値近傍の所定値（但し上限値＝ｄｔ以下）とする。

その概要を図８のタイムチャートを用いて説明する。図８（ａ）に示すように、ある入力に対するノミナルむだ時間（むだ時間のノミナル値）と実むだ時間との間に図示のような誤差が生じている場合を想定する。ノミナルむだ時間は、予め設定したノミナルモデルスケジューラを用いて算出される。サンプリングタイミング（ｄ・ｄｔ）を超えた部分が余むだ時間Ｌ１である。

この場合、同定が逐次行われることにより、図８（ｂ）に示すように、余むだ時間の推定値Ｌ１＿hatが変化し、ｄｔ近傍（上述したｄｔ−ε2〜ｄｔの範囲）まで増加する。そして、余むだ時間の推定値Ｌ１＿hatがｄｔ近傍で所定時間以上とどまっている場合には、図８（ｃ）に示すように、むだサンプリング回数の推定値ｄ＿hatが１つ繰り上げられ（ｄ＿hat←ｄ＋１）、更に、余むだ時間の推定値Ｌ１＿hatがゼロ近傍の値に更新される。

その後、同定が逐次行われることにより、図８（ｄ）に示すように、余むだ時間の推定値Ｌ１＿hatが増加し、それに伴いむだ時間の推定値が実むだ時間に収束する。この場合、ノミナルむだ時間と実むだ時間との誤差がサンプリング周期の２倍以上であっても、上記更新を繰り返すことにより、むだ時間の推定値が実むだ時間に収束する。また、実むだ時間に収束後、劣化等により実むだ時間の変動がある場合にも、上記更新を実行することにより常にむだ時間の推定値を実むだ時間に収束させることができる。

なお、図８では、余むだ時間の推定値Ｌ１＿hatが増加する場合を例示し、それに伴い、むだサンプリング回数の推定値ｄ＿hatが１つ繰り上げられる様子を説明したが、これとは逆に、余むだ時間の推定値Ｌ１＿hatが減少する場合には、それに伴いむだサンプリング回数の推定値ｄ＿hatが１つ繰り下げられると共に、余むだ時間の推定値Ｌ１＿hatがｄｔ近傍の値に更新される。

以上説明した適応制御処理等は、ＥＣＵ２９によって図９乃至図１５の各プログラムに従って実行される。以下、図９乃至図１５の各プログラムの処理内容を説明する。
［燃料噴射制御プログラム］
図９の燃料噴射制御プログラムは、エンジン運転中に所定周期（例えば３０℃Ａ周期）で実行される。本プログラムが起動されると、まずステップ１０１で、例えば基本噴射量マップ等を用いて、その時点のエンジン回転速度や負荷等の運転状態パラメータに基づいて基本噴射量ＴＰを算出する。

この後、ステップ１０２に進み、空燃比Ｆ／Ｂ制御実行条件が成立しているか否かを判定する。ここで、空燃比Ｆ／Ｂ制御実行条件は、例えば、次の条件(1)〜(4)によって判定される。

(1)上流側排気ガスセンサ３１と下流側排気ガスセンサ３２が活性化していること
(2)空燃比制御システムが正常に動作すること（フェイル判定が行われていないこと）
(3)エンジン冷却水温が所定温度以上（例えば７０℃以上）であること、つまりエンジン１１の暖機後であること
(4)エンジン運転状態がＦ／Ｂ実行領域であること
これらの条件(1)〜(4)の中で、１つでも満たさない条件があれば、空燃比Ｆ／Ｂ制御実行条件が不成立となり、ステップ１０３に進み、空燃比補正係数ＦＡＦを「１．０」に設定する。この場合は、空燃比Ｆ／Ｂ制御は行われないこととなる。

これに対して、上記条件(1)〜(4)を全て満たせば、空燃比Ｆ／Ｂ制御実行条件が成立して、ステップ１０４に進み、後述する図１０の目標空燃比算出プログラムを実行して、上流側触媒２５の上流側の目標空燃比λＴＧを算出した後、ステップ１０５に進み、上流側排気ガスセンサ３１で検出した空燃比と目標空燃比λＴＧとの偏差を小さくするように空燃比補正係数ＦＡＦを算出する。

以上のようにして、ステップ１０３又は１０５で、空燃比補正係数ＦＡＦを算出した後、ステップ１０６に進み、空燃比補正係数ＦＡＦやその他各種の補正係数（例えば冷却水温補正係数、学習補正係数、加減速時の補正係数等）により基本噴射量ＴＰを補正して、要求燃料噴射量ＴＡＵを算出する。

ＴＡＵ＝ＴＰ×ＦＡＦ×各種補正係数
［目標空燃比算出プログラム］
図１０の目標空燃比算出プログラムは、上記図９の燃料噴射制御プログラムのステップ１０４で実行されるサブルーチンである。本プログラムが起動されると、まずステップ２０１で、後述する図１１の同定実行条件判定プログラムを実行して、同定実行条件の成否を判定し、その判定結果に応じて同定実行フラグをセット／リセット（ＯＮ／ＯＦＦ）する。

その後、ステップ２０２に進み、図１１の同定実行条件判定プログラムの実行結果に基づいて同定実行フラグがＯＮであるか否かを判定し、同定実行フラグがＯＮであれば、同定実行条件が成立していると判断して、ステップ２０３に進み、下流側排気ガスセンサ３２の目標電圧を設定する。この際、元々の目標電圧に例えば所定振幅を加えた電圧を目標電圧として設定してもよい。

この後、ステップ２０４に進み、後述する図１３の同定処理プログラムを実行して、制御対象の入力（上流側触媒２５の流入過不足酸素量）を離散プラントモデルに入力した時のプラントモデル出力と制御対象の出力（下流側排気ガスセンサ３２の検出電圧）との誤差ｅをゼロに近づけるように離散プラントモデルの可変パラメータを同定する。この後、ステップ２０５の処理に進む。

一方、同定実行フラグがＯＦＦであれば、同定実行条件が不成立であると判断して、上記同定処理（ステップ２０３、２０４）を行わずに、ステップ２０５に進む。

このステップ２０５では、サブＦ／Ｂ制御入力算出処理を実行する。このサブＦ／Ｂ制御入力算出処理では、同定処理により学習したパラメータに基づいて、下流側排気ガスセンサ３２の検出電圧を目標電圧に近付けるように、モデルベース制御則（例えば現代制御やロバスト制御）によりサブＦ／Ｂ制御入力となる入力酸素量（上流側触媒２５の流入過不足酸素量）を算出する。この際、ＰＩ制御又はＰＩＤ制御のような古典制御を用いてサブＦ／Ｂ制御入力を算出するようにしても良い。

本実施形態では、入力酸素量＝目標酸素量と見なすようにしているが、入力酸素量＝検出酸素量としても良いし、空燃比補正係数から入力酸素量を算出するようにしても良い。

この後、ステップ２０６に進み、目標酸素量（入力酸素量）より目標空燃比λＴＧを次のようにして算出する。まず、今回演算時の目標燃料過剰率φTGを次式により算出する。

φTG＝Ｏ(k)／（Ｋ・Ｗ・ｄｔ）
Ｏ(k)：今回の目標酸素量
Ｋ：標準空気酸素質量割合
Ｗ：上流側触媒２５の流入排気ガス流量［ｇ／ｓ］
ｄｔ：演算周期
ここで、Ｗは、上流側触媒２５の流入排気ガス流量の代わりに、エンジン吸入空気量を用いても良い。また、演算周期ｄｔは、時間同期でも良いし、エンジン回転速度同期でも良い。エンジン回転速度同期の場合は、演算周期ｄｔはエンジン回転速度Ｎeの関数となる。また、ここでφTGは正確には演算周期当たりの目標燃料過剰量を表すが、単位時間当たりの目標燃料過剰量を用いてもよい。その場合、φTGは次式となる。

φTG＝Ｏ(k)／（Ｋ・Ｗ）
この目標燃料過剰率φＴＧの逆数が目標空燃比λＴＧとなる。

λＴＧ＝１／φＴＧ
この後、ステップ２０７に進み、目標空燃比λＴＧにベース目標空燃比を加算（又は乗算）して最終的な目標空燃比λＴＧを求める。この際、ベース目標空燃比は、例えばエンジン運転条件をパラメータとするマップ等により算出される。
［同定実行条件判定プログラム］
図１１の同定実行条件判定プログラムは、図１０の目標空燃比算出プログラムのステップ２０１で実行されるサブルーチンである。本プログラムが起動されると、まずステップ３０１で、触媒劣化検出時であるか否かを判定し、触媒劣化検出時であれば、ステップ３０５に進み、同定実行フラグをＯＮにセットする。

一方、上記ステップ３０１で、触媒劣化検出時でないと判定されれば、ステップ３０２に進み、後述する図１２の運転状態判定プログラムを実行して、現在のエンジン運転状態が定常か過渡かを判定する。この後、ステップ３０３に進み、図１２の運転状態判定プログラムの実行結果に基づいて現在のエンジン運転状態が過渡であるか否かを判定し、過渡であれば、ステップ３０４に進み、同定実行フラグをＯＦＦにリセットし、定常であれば、ステップ３０５に進み、同定実行フラグをＯＮにセットする。

尚、同定実行条件は、触媒劣化検出時であることのみとしたり、逆に触媒劣化検出時であることを同定実行条件から外したりすることも可能である。また、空燃比Ｆ／Ｂ制御が実行される全運転領域で同定実行フラグをＯＮにセットするようにしても良い。
［運転状態判定プログラム］
図１２の運転状態判定プログラムは、図１１の同定実行条件判定プログラムのステップ３０２で実行されるサブルーチンであり、運転状態判定手段としての役割を果たす。本プログラムが起動されると、まずステップ３０１で、エンジン回転速度変化量（変化率の絶対値）が所定値以下であるか否かを判定し、エンジン回転速度変化量が所定値以下でなければ、ステップ４０４に進み、現在のエンジン運転状態が過渡であると判定する。

一方、上記ステップ４０１で、エンジン回転速度変化量が所定値以下と判定されれば、ステップ４０２に進み、負荷変化量（負荷変化率の絶対値）が所定値以下であるか否かを判定し、負荷変化量が所定値以下でなければ、ステップ４０４に進み、現在のエンジン運転状態が過渡であると判定する。

上記ステップ４０１と４０２でいずれも「Ｙｅｓ」と判定された場合、つまりエンジン回転速度変化量が所定値以下で且つ負荷変化量が所定値以下である場合には、ステップ４０３に進み、現在のエンジン運転状態が定常であると判定する。

尚、エンジン運転中に、所定時間以上の期間、定常と判定されない場合は、同定の実行頻度を確保するために、エンジン回転速度変化量と負荷変化量に対する定常判定値（ステップ４０１と４０２の所定値）を大きい値に変更するようにしても良い。
［同定処理プログラム］
図１３の同定処理プログラムは、図１０の目標空燃比算出プログラムのステップ２０４で実行されるサブルーチンであり、特許請求の範囲でいう同定手段としての役割を果たす。本プログラムが起動されると、まずステップ５０１で、入力酸素量（上流側触媒２５の流入過不足酸素量）を次の(1)、(2)のいずれかの方法で算出する。

(1)目標酸素量をベース目標空燃比変化分に応じて補正して入力酸素量を算出する。

(2)上流側触媒２５の流入ガス空燃比を、上流側排気ガスセンサ３１の検出空燃比又は燃料噴射補正係数に基づく値とし、次式により入力酸素量を算出する。

入力酸素量＝Ｋ・Ｗ・（理論当量比−触媒流入ガス当量比）・ｄｔ
Ｋ：標準空気酸素質量割合
Ｗ：触媒流入排気ガス流量［ｇ／ｓ］（エンジン吸入空気量でもよい）
ｄｔ：演算周期
ここで、Ｗは、上流側触媒２５の流入排気ガス流量の代わりに、エンジン吸入空気量を用いても良い。また、演算周期ｄｔは、時間同期でも良いし、エンジン回転速度同期でも良い。エンジン回転速度同期の場合は、演算周期ｄｔはエンジン回転速度Ｎeの関数となる。また、ここで入力酸素量は正確には演算周期当たりの入力酸素量を表すが、単位時間当たりの入力酸素量を用いてもよい。その場合、入力酸素量は次式により算出される。

入力酸素量＝Ｋ・Ｗ・（理論当量比−触媒流入ガス当量比）
入力酸素量の算出後、ステップ５０２に進み、離散プラントモデルの出力を算出する。この際、離散プラントモデルの出力と検出出力に同一の所定信号処理を施すものとする。この所定信号処理は、例えば、ＬＰＦ、ＨＰＦ、ＢＰＦのいずれかであるが、他の処理であっても良い。

この後、ステップ５０３に進み、制御対象の実出力と離散プラントモデルの出力とから同定誤差ｅを算出する。そして、次のステップ５０４で、前述したパラメータ調整則に従ってパラメータ適応処理を実行する。これにより、制御対象の実出力と離散プラントモデルの出力との誤差をゼロに近づけるように離散モデルパラメータの推定値θ＿hatを算出する。

この後、ステップ５０５に進み、学習判定処理により同定誤差ｅが所定範囲内であるか否かを判定し、同定誤差ｅが所定範囲内であれば、むだ時間Ｌの学習が許可される。この後、ステップ５０６に進み、前記算出した離散モデルパラメータの推定値θ＿hatに対して図５の連続化式を用いて連続化の処理を実行し、連続モデルパラメータ（余むだ時間の推定値Ｌ１＿hat）を算出する。尚、上記ステップ５０５の学習判定処理は、ステップ５０６の連続化処理後に行うようにしても良い。

そして、次のステップ５０７（更新手段）で、後述する図１４及び図１５のむだ時間更新処理プログラムを実行して、余むだ時間Ｌ１を学習する。この後、ステップ５０８に進み、図４の離散化式を用いて連続モデルパラメータを離散化処理して離散プラントモデルのパラメータを更新する。

尚、パラメータ適応処理のアルゴリズムは、本実施形態に限定されず、他のアルゴリズムを用いても良い。
［むだ時間更新処理のプログラム］
図１４及び図１５のむだ時間更新処理プログラムは、図１３の同定処理プログラムのステップ５０７で、当該プログラムの演算周期ｄｔ（例えば１８０℃Ａ毎又は所定時間毎）で繰り返し実行されるサブルーチンであり、むだ時間推定手段及び学習手段としての役割を果たす。本プログラムで用いる微小な正の定数ε1，ε2，δは、ε1＜ε2＜δの関係に設定されている。

本プログラムが起動されると、まずステップ６０１で、前記図１３のステップ５０６で算出した余むだ時間の推定値Ｌ１＿hatが所定の判定値「０＋ε2」よりも小さいか否かを判別する。すなわち、余むだ時間の推定値Ｌ１＿hatが余むだ時間の下限値近傍（０〜０＋ε2）にあるか否かを判別する。そして、Ｌ１＿hat≧０＋ε2であれば、ステップ６０３に進み、カウンタｃｎｔ＿ｍｉｎｓを０にクリアする。

一方、上記ステップ６０１で、Ｌ１＿hat＜０＋ε2と判定されば（つまりＬ１＿hatが下限値近傍にあれば）、ステップ６０２に進み、カウンタｃｎｔ＿ｍｉｎｓを１インクリメントする。その後、ステップ６０４に進み、カウンタｃｎｔ＿ｍｉｎｓの値が所定値Ｎ１以上で且つ余むだ時間の推定値Ｌ１＿hatが微小な負の定数（−δ）よりも小さいこと、又は、カウンタｃｎｔ＿ｍｉｎｓの値が所定値Ｎ２以上であることのいずれかが成立するか否かを判別する。ここで、Ｎ１＞Ｎ２である。なお、このステップ６０４の判定条件を、カウンタｃｎｔ＿ｍｉｎｓの値が所定値Ｎ１以上で且つ余むだ時間の推定値Ｌ１＿hatが微小な負の定数（−δ）よりも小さいことのみとしたり、或は、カウンタｃｎｔ＿ｍｉｎｓの値が所定値以上であることのみとしたりすることも可能である。

そして、このステップ６０４で「Ｎｏ」と判定されれば、ステップ６０６に進み、余むだ時間Ｌ１＿ａｄｐを「０」にする。

一方、上記ステップ６０４で「Ｙｅｓ」と判定されれば、ステップ６０５に進み、余むだ時間Ｌ１＿ａｄｐ、むだサンプリング回数変化量Δｄ＿ａｄｐ、カウンタｃｎｔ＿ｍｉｎｓをそれぞれ更新する。

Ｌ１＿ａｄｐ＝ｄｔ−δ
Δｄ＿ａｄｐ＝−１
ｃｎｔ＿ｍｉｎｓ＝０
以上のようにして、ステップ６０３又は６０５又は６０６の処理を実行した後、図１５のステップ６０７に進み、余むだ時間の推定値Ｌ１＿hatが所定の判定値「ｄｔ−ε2」以上であるか否かを判定する。すなわち、余むだ時間の推定値Ｌ１＿hatが余むだ時間の上限値近傍（ｄｔ−ε2〜ｄｔ）にあるか否かを判定する。

もし、Ｌ１＿hat＜ｄｔ−ε2であれば、ステップ６０８に進み、カウンタｃｎｔ＿ｐｕｌｕｓを０にクリアする。

一方、上記ステップ６０７で、Ｌ１＿hat≧ｄｔ−ε2と判定されれば（つまりＬ１＿hatが上限値近傍にあれば）、ステップ６０９に進み、カウンタｃｎｔ＿ｐｕｌｓを１インクリメントする。

その後、ステップ６１０に進み、カウンタｃｎｔ＿ｐｕｌｓの値が所定値Ｎ１以上で且つ余むだ時間の推定値Ｌ１＿hatが上限値から微小な正の定数δだけ大きい値（ｄｔ＋δ）よりも小さいこと、又は、カウンタｃｎｔ＿ｐｕｌｓの値が所定値Ｎ２以上であることの何れかが成立するか否かを判定する。尚、このステップ６１０の判定条件を、カウンタｃｎｔ＿ｐｕｌｓの値が所定値Ｎ１以上で且つ余むだ時間の推定値Ｌ１＿hatが（ｄｔ＋δ）よりも小さいことのみとしたり、カウンタｃｎｔ＿ｐｕｌｓの値が所定値Ｎ２以上であることのみとしたりすることも可能である。

このステップ６１０で「Ｎｏ」と判定されれば、ステップ６１１に進み、余むだ時間Ｌ１＿ａｄｐを「ｄｔ−ε1」として更新する。

一方、上記ステップ６１０で「Ｙｅｓ」と判定されれば、ステップ６１２に進み、余むだ時間Ｌ１＿ａｄｐ、むだサンプリング回数変化量Δｄ＿ａｄｐ、カウンタｃｎｔ＿ｐｕｌｓをそれぞれ更新する。

Ｌ１＿ａｄｐ＝δ
Δｄ＿ａｄｐ＝＋１
ｃｎｔ＿ｐｕｌｓ＝０
以上のようにして、ステップ６０８又は６１１又は６１２の処理を実行した後、図１４のステップ６１３に進み、学習フラグが学習値の更新を意味するＯＮであるか否かを判定し、学習フラグがＯＮであれば、ステップ６１４に進み、前回までのむだサンプリング回数変化量の積算値Σ（Δｄ＿ａｄｐ）をクリアし、学習フラグがＯＦＦであれば、積算値Σ（Δｄ＿ａｄｐ）のクリアは行わない。

この後、ステップ６１５に進み、むだサンプリング回数学習値ｄ＿ｌｒｎに積算値Σ（Δｄ＿ａｄｐ）を加算して、むだサンプリング回数ｄ＿ａｄｐを更新する。

ｄ＿ａｄｐ＝ｄ＿ｌｒｎ＋Σ（Δｄ＿ａｄｐ）
図１６は、図２に記載の離散プラントモデルＭ３１、適応機構Ｍ３２及びむだ時間学習部Ｍ４０の概略を示す図であり、図１６中の（ａ）（ｂ）欄は、プラントモデルの入力（酸素量）に対する出力（酸素量）を縦軸とし、横軸を経過時間とした場合における前記出力の時間変化を示す図であり、入力に対する出力の応答特性を表している。そして、（ａ）（ｂ）欄中の点線は、制御対象による実際の出力変化を示し、実線は、プラントモデルのパラメータについて同定を開始した時点におけるプラントモデルの出力変化を示す。

適応機構Ｍ３２では、先述の如く誤差ｅをゼロに近づけるように、プラントモデルパラメータのうちむだ時間以外のパラメータ、つまり応答遅れ及び定常ゲイン等を同定する（図１６（ａ）参照）。このように応答遅れを示す時定数を同定により調整することで、応答特性の傾きが点線に示す実際の特性に近づくように調整される。

むだ時間学習部Ｍ４０では、先述の如く誤差ｅをゼロに近づけるように、プラントモデルパラメータのうちむだ時間を同定する（図１６（ｂ）参照）。このようにむだ時間を同定により調整することで、応答特性の時間軸方向のずれ量が点線に示す実際の特性に近づくように調整される。

ここで、むだ時間学習部Ｍ４０にて学習されたむだ時間の学習初期値が、実際のむだ時間に対して大きくずれていると、適応機構Ｍ３２においてむだ時間以外のプラントモデルパラメータ（例えば時定数や定常ゲイン）を逐次同定するにあたり、時定数や定常ゲインの同定値が振動又は発散することが懸念される（図１９参照）。

そこで、このような懸念を解消すべく本実施形態では、適応機構Ｍ３２においてプラントモデルパラメータ（時定数や定常ゲイン）を同定するにあたり、これらのパラメータの値を予め設定された範囲内に制限する射影アルゴリズムを用いて同定することで、上記振動及び発散の抑制を図っている。

以下、適応機構Ｍ３２においてパラメータの更新範囲を制限して同定するにあたり、その制限手法について説明する。本実施形態では以下の制約条件１，２を満たすよう制限している。

はじめに、［数１］〜［数４］式等にて先述した離散プラントモデルＭ３１及び同定パラメータθ＿hatは、次のように表すことができる。

＜制約条件１について＞
ここで、制御対象は、有界な大きさの任意の入力に対してその出力も有界となる安定系であると考えられるため、本実施形態では、離散プラントモデルＭ３１が、安定系となるように、適応機構Ｍ３２にて同定されるパラメータθ＿hatに制約条件を与える。具体的には、以下の［数１５］［数１６］式を満たすよう制限することを制約条件１とする。

上記ω、ζは連続モデルパラメータであり、ωは固有角振動数、ζは減衰係数である。そして、ω×ζは２次遅れ系における広義の時定数に相当する。また、ＴＭは予め設定された所定の安定度制限値である。この安定度制限値は、大きくなるほどより安定度が高いことを示す値であり、例えばＴＭ＝０とすると単に系が安定であるように制限され，正の所定値とするとより安定側（系の応答が振動的でなくなるよう）に制限される。なお、上記［数１５］式を離散モデルパラメータθ＿hat（同定値）で表現すると下記［数１６］式となる。

＜制約条件２について＞
本実施形態では、制御対象はその定常ゲインの符号が時不変であることを前提としているが、前記同定手段における同定過程において定常ゲイン推定値の符号が誤って切り替わってしまう虞があり，この符号変化は応答が逆位相となることを意味するため、同定手段及びむだ時間推定手段の同定値及び推定値の振動又は発散を招いてしまう。そこで、例えば定常ゲインＫ＿hatが正となるようプラントモデルの入出力を定義している場合には、定常ゲインが負にはならないよう制限する。具体的には、以下の［数１７］［数１８］式を満たすよう制限することを制約条件２とする。

ここで、制約条件１で述べたプラントモデルが安定系であるという条件満たすためには［数１７］式中の分母が正である必要があることが分かっている。そこで、［数１７］式の条件は次の［数１８］式で表せる。

ここで、εｎ，εｄは予め設定される正の微小な所定値である。以上により、［数１６］［数１８］式が同定パラメータθ_hatに対する制約条件１，２である。

このように導出される制約条件１，２を満たしながら同定パラメータθ_hatを更新する方法としては、例えば、図１８を用いて先述した射影アルゴリズムを用いることが挙げられる。また、これら２つの制約条件１，２のすべてを満たすよう構成してもよいし、いずれか一方の制約条件を満たす構成としてもよい。また、［数１５］式において時定数に相当するω×ζを所定上限値で制限する条件を追加する構成としてもよいし、［数１７］式において定常ゲインＫを所定上限値で制限する条件を追加する構成としてもよい。

図１７は、上記制約条件１を満たすよう時定数を同定することによる、振動、発散の抑制効果を示す図である。そして、図１７（ａ）（ｂ）（ｃ）は図１９（ａ）（ｂ）（ｃ）に対応するものであり、図１７（ｂ）（ｃ）欄の一点鎖線に示すようにむだ時間の初期値が大きい場合であっても、制約条件１により時定数の値は［数１６］式を満たすよう制限されるので、プラントモデルの入力に対する出力の応答特性を立てすぎないよう（時定数が過小にならないよう）抑制でき、（ｂ）欄の実線に示すように同定できる。よって、むだ時間及び時定数を同定するにあたり、オーバーシュート領域にて振動すること或いは不安定領域にて発散することが回避される。

以上詳述した本実施形態によれば、以下の効果が得られるようになる。

（１）適応機構Ｍ３２及びむだ時間学習部Ｍ４０により離散プラントモデルを同定するにあたり、［数１６］［数１８］式に示す制約条件１，２の少なくとも一方を満たすように同定するので、適応機構Ｍ３２及びむだ時間学習部Ｍ４０にて同定するパラメータ（図１７の例ではむだ時間及び時定数）が振動、発散することを回避でき、むだ時間及び時定数を真値に近づけるように同定する同定速度の向上を図ることができる。

（２）制約条件１，２を満たしながら同定パラメータθ_hatを更新する手法として射影アルゴリズムを用いるので、同定パラメータθ_hatの値を予め設定された範囲内に制限しつつ逐次同定することを容易に実現できる。

（３）制御対象は安定系であり、離散プラントモデルＭ３１の同定対象となるパラメータθ＿hatも安定系となる。この点に着目した本実施形態では、２次遅れ系における広義の時定数「ω＿hat×ζ＿hat」が安定度制限値ＴＭより大きくなるよう制限するので（［数１６］式参照）、時定数を存在領域内に制限しつつ各種同定パラメータθ_hat（図１７の例ではむだ時間及び時定数）を逐次同定することを容易に実現できる。

（４）同定パラメータθ_hatのうち定常ゲインＫ＿hatは、プラントモデルにおいて正の値となるよう定義されている。この点に着目した本実施形態では、定常ゲインＫ＿hatが定義に反して負の値とならないよう制限するので（［数１８］式参照）、定常ゲインＫ＿hatを存在領域内に制限しつつ各種同定パラメータθ_hatを逐次同定することを容易に実現できる。

（５）本実施形態では、拡張ｚ変換により離散化を行って離散モデルパラメータに余むだ時間Ｌ１の情報を含ませると共に、該離散モデルパラメータから算出される余むだ時間Ｌ１の推定値を基に、離散プラントモデルのむだ時間を制御対象の実むだ時間に近づけるように推定し、さらに、離散化のむだサンプリング回数ｄを増加又は減少させることで、推定むだ時間を実むだ時間の変化に追従させるようにしたので、応答遅れだけでなくむだ時間も逐次推定できることとなり、制御対象のばらつきや経時変化によらず制御対象の動特性を十分把握することができ、高い制御精度を維持することができる。また、むだ時間の情報を種々に活用することもできる。

（６）本実施形態では、同定誤差に基づいて逐次推定された離散モデルパラメータを連続モデルパラメータに変換し、その連続モデルパラメータから余むだ時間の推定値Ｌ１＿hatを算出するようにしたので、離散モデルパラメータから演算を行う場合に比べて、余むだ時間の推定値Ｌ１＿hatの演算を容易に行うことができる。

（７）推定むだ時間を学習する構成としたため、制御対象の個体差や経時変化等により実むだ時間のずれが生じている場合にも、それに追従した制御を継続的に実施できる。このとき、プラントモデル出力と制御対象の実出力との誤差が所定範囲内である場合に学習値を更新するようにしているので、適正な同定が行われていることを条件に学習処理が行われ、学習値の信頼性が高められる利点もある。

（他の実施形態）
上記実施形態は、以下のように変更して実施してもよい。また、本発明は上記実施形態の記載内容に限定されず、以下に説明する各実施形態の特徴的構成をそれぞれ任意に組み合わせるようにしてもよい。

・上記実施形態では、［数１５］［数１６］式で用いられる安定度制限値ＴＭを固定値としているが、エンジン運転状態や制御対象（プラント）の運転状態に応じて可変設定してもよい。前記エンジン運転状態の具体例としては、エンジン出力軸（クランク軸）の回転速度、エンジン負荷等が挙げられる。前記制御対象の運転状態の具体例としては、触媒２５及び下流側排気ガスセンサ３２の劣化度合や温度等が挙げられる。

・上記実施形態では、応答遅れを示すパラメータを存在領域に制限しながら同定するにあたり、固有角振動数ω及び減衰係数ζにより定義される２次遅れの値（ω×ζ）を応答遅れパラメータとしているが、１次遅れにより定義される時定数を応答遅れパラメータとしてもよいし、３次以上の高次遅れにより定義される値を応答遅れパラメータとしてもよい。

・上記実施形態では、離散モデルパラメータの推定値θ＿hatを連続モデルパラメータに変換し、その連続モデルパラメータに含まれる余むだ時間の推定値Ｌ１＿hatを用いてむだ時間Ｌ（むだサンプリング回数ｄ）を推定する構成としたが、これに代えて、離散パラメータの連続化を行わず、余むだ時間Ｌ１を含む離散パラメータの推定値を用いてむだ時間Ｌ（むだサンプリング回数ｄ）を推定することも可能である。

・上記実施形態では、制御対象（上流側触媒２５と下流側排気ガスセンサ３２）をむだ時間を含む２次遅れ系で近似するようにしたが、本発明は、制御対象をむだ時間を含む３次以上の高次遅れ系で近似するようにしても良い。この場合、３次以上の高次遅れ系は、３つ以上の１次遅れ系に分割して離散化した後結合して近似するようにすれば良い。

・本発明の適用範囲は、上流側触媒２５と下流側排気ガスセンサ３２を含むサブＦ／Ｂ系（空燃比制御システム）に限定されず、エンジン１１と上流側排気ガスセンサ３１を含むメインＦ／Ｂ系（内燃機関制御システム）に本発明を適用しても良い。この場合、同定手段は、エンジン運転条件に応じて設定されたノミナルモデルパラメータ（モデルパラメータのノミナル値）からの偏差を同定し、エンジン運転条件を表す信号の高周波振動成分を抑制するためのフィルタ手段を設けるようにすると良い。ノミナルモデルパラメータは、予め設定したノミナルモデルスケジューラを用いて算出される。

・上記実施形態では、排気管２３に２つの触媒２５，３０を設けたが、触媒を１つのみ設けたシステムにも本発明を適用して実施できる。その他、本発明は、サブＦ／Ｂ系、メインＦ／Ｂ系に限定されず、他の制御対象に適用しても良い。

本発明の一実施形態におけるエンジン制御システムの概略を示す構成図である。 ＥＣＵで実現するメインＦ／Ｂ制御機能とサブＦ／Ｂ制御機能の概要を示すブロック図である。 ２次遅れ＋むだ時間で近似する連続プラントモデルを、１次遅れ系に分割して離散化した後結合して離散プラントモデルを作成する手法を説明する図である。 連続モデルパラメータを離散化して離散モデルパラメータを演算する方法を説明する図である。 離散モデルパラメータ（同定値）を連続化して連続モデルパラメータを演算する方法を説明する図である。 図５中のＺｒｅ，Ｚｉｍを算出する方法を説明する図である。 制御入力ｕの変化を表すタイムチャートである。 むだ時間更新処理の概要を具体的に示すタイムチャートである。 燃料噴射制御プログラムの処理の流れを示すフローチャートである。 目標空燃比算出プログラムの処理の流れを示すフローチャートである。 同定実行条件判定プログラムの処理の流れを示すフローチャートである。 運転状態判定プログラムの処理の流れを示すフローチャートである。 同定処理プログラムの処理の流れを示すフローチャートである。 むだ時間更新処理プログラムの処理の流れを示すフローチャートである（その１）。 むだ時間更新処理プログラムの処理の流れを示すフローチャートである（その２）。 （ａ）は適応機構Ｍ３２の同定による調整機能を示す図、（ｂ）はむだ時間学習部Ｍ４０の同定による調整機能を示す図である。 同実施形態によりモデルパラメータを制限した場合における、同定値の挙動を示す図である。 同実施形態が採用する射影アルゴリズムの概念を説明する図である。 同実施形態に反してモデルパラメータを制限しなかった場合における、同定値の挙動を示す図である。

符号の説明

１１…エンジン（内燃機関）、２３…排気管、２５…上流側の触媒、２９…ＥＣＵ（同定手段，余むだ時間推定手段，むだ時間推定手段、更新手段）、３０…下流側の触媒、３１…上流側排気ガスセンサ、３２…下流側排気ガスセンサ、Ｍ３１…離散プラントモデル、Ｍ３２…適応機構、Ｍ４０…むだ時間学習部。

Claims (7)

1. 所定のサンプリング周期（ｄｔ）で離散化された数式により制御対象を表したプラントモデルを用い、前記制御対象への入力を前記プラントモデルに加えた際のプラントモデル出力と前記制御対象の実出力との差である同定誤差をゼロに近づけるように、前記プラントモデルのパラメータを逐次同定する同定手段と、
   前記制御対象のむだ時間（Ｌ）を前記サンプリング周期で除算した商をむだサンプリング回数（ｄ）、余りを余むだ時間（Ｌ１）としてモデル化した場合に、前記同定手段により同定された同定パラメータを用いて前記プラントモデルに含まれている前記余むだ時間を推定する余むだ時間推定手段と、
   前記余むだ時間推定手段により推定された推定余むだ時間が、所定の増加条件を満たす程度に増加した場合には前記むだサンプリング回数の推定値を一つ繰り上げるよう更新し、所定の減少条件を満たす程度に減少した場合には前記むだサンプリング回数の推定値を一つ繰り下げるよう更新する更新手段と、
   前記むだサンプリング回数更新手段により更新された前記むだサンプリング回数の推定値に基づき前記むだ時間を推定するむだ時間推定手段と、
   を備え、
   前記同定手段は、前記同定パラメータの値が予め設定された制約領域内に制限されつつ実行されることを特徴とする制御装置。
2. 前記同定パラメータは、射影アルゴリズムにより設定された制約領域内に制限されることを特徴とする請求項１に記載の制御装置。
3. 前記同定パラメータを制限する前記制約領域は、前記プラントモデルの時定数に相当する値が所定範囲内となるように設定されることを特徴とする請求項１又は２に記載の制御装置。
4. 前記プラントモデルの時定数に相当する値の前記所定範囲は、前記プラントモデルが安定系であるという条件より設定されることを特徴とする請求項３に記載の制御装置。
5. 前記同定パラメータを制限する前記制約領域は、前記プラントモデルの定常ゲインの符号が前記制御対象の定常ゲインの符号に一致するように設定されることを特徴とする請求項１〜４のいずれか１つに記載の制御装置。
6. 前記同定パラメータを制限する前記制約領域は、前記プラントモデルの極が安定であるように設定されることを特徴とする請求項１〜５のいずれか１つに記載の制御装置。
7. 前記制御対象は、内燃機関の排気通路に配設された排ガス浄化用の触媒と、前記排気通路において前記触媒の上流側で特定排ガス濃度を検出する上流側排ガスセンサと、前記排気通路において前記触媒の下流側で特定排ガス濃度を検出する下流側排ガスセンサと、を有する空燃比制御システムであり、
   前記上流側排ガスセンサで検出した空燃比を前記触媒の上流側目標空燃比に一致させるようにフィードバックして、供給空燃比を前記制御対象への制御入力として演算するメインフィードバック制御手段と、
   前記下流側排ガスセンサの検出した空燃比を前記触媒の下流側目標空燃比に一致させるようにフィードバックして、前記上流側目標空燃比を前記制御入力の算出に用いる値として演算するサブフィードバック制御手段と、
   を備えることを特徴とする請求項１〜６のいずれか１つに記載の制御装置。