(57)【要約】【課題】エアフロメータの応答遅れを適切に把握して内燃機関を制御する内燃機関の制御装置を提供する。【解決手段】本制御装置は、機関吸気通路に設けられたエアフロメータ４１と、スロットル通過空気流量を推定するスロットルモデルＭ１０と、スロットルモデルにより算出されたスロットル通過空気流量の推定値に基づいてＡＦＭモデル計算式を用いてエアフロメータの予想出力値を算出するＡＦＭモデルＭ５０とを具備し、エアフロメータの実測値と予想出力値とを用いて内燃機関を制御する。さらに、本制御装置は、エアフロメータの予想出力値と、エアフロメータの実測値との間の誤差を算出する誤差算出手段と、スロットル通過空気流量の増大時に誤差算出手段によって算出された誤差とスロットル通過空気流量の減少時に誤差算出手段によって算出された誤差とに基づいてＡＦＭモデル計算式を修正するモデル計算式修正手段とをさらに具備する。【選択図】図１
2007231840.tif 000002
【特許請求の範囲】
【請求項１】
機関吸気通路に設けられたエアフロメータと、
スロットル通過空気流量を推定するスロットル通過空気流量推定手段と、
上記スロットル通過空気流量推定手段により推定されたスロットル通過空気流量の推定値に基づいてエアフロメータモデル計算式を用いてエアフロメータによって出力されるであろう出力値を予想出力値として算出するエアフロメータモデルとを具備し、
上記エアフロメータの実測値と上記エアフロメータモデルによって算出されたエアフロメータの予想出力値とを用いて内燃機関を制御する、内燃機関の制御装置において、
上記エアフロメータモデルによって算出されたエアフロメータの予想出力値と、上記エアフロメータの実測値との間の誤差を算出する誤差算出手段と、
スロットル通過空気流量が増大しているときに上記誤差算出手段によって算出された誤差とスロットル通過空気流量が減少しているときに上記誤差算出手段によって算出された誤差とに基づいて上記エアフロメータモデル計算式を修正するモデル計算式修正手段とをさらに具備する、内燃機関の制御装置。
【請求項２】
上記モデル計算式修正手段は、上記誤差算出手段によって算出される誤差がスロットル通過空気流量が増大しているときと減少しているときとで等しくなるように上記エアフロメータモデル計算式を修正する、請求項１に記載の内燃機関の制御装置。
【請求項３】
機関吸気通路に設けられたエアフロメータと、
スロットル通過空気流量を推定するスロットル通過空気流量推定手段と、
上記スロットル通過空気流量推定手段により推定されたスロットル通過空気流量の推定値に基づいてエアフロメータモデル計算式を用いてエアフロメータによって出力されるであろう出力値を予想出力値として算出するエアフロメータモデルとを具備し、
上記エアフロメータの実測値と上記エアフロメータモデルによって算出されたエアフロメータの予想出力値とを用いて内燃機関を制御する、内燃機関の制御装置において、
上記エアフロメータモデルによって算出されたエアフロメータの予想出力値と、上記エアフロメータの実測値との間の誤差を算出する誤差算出手段と、
内燃機関が加速しているときに上記誤差算出手段によって算出された誤差と内燃機関が減速しているときに上記誤差算出手段によって算出された誤差とに基づいて上記エアフロメータモデル計算式を修正するモデル計算式修正手段とをさらに具備する、内燃機関の制御装置。
【請求項４】
上記モデル計算式修正手段は、上記誤差算出手段によって算出される誤差が内燃機関が加速しているときと減速しているときとで等しくなるように上記エアフロメータモデル計算式を修正する、請求項３に記載の内燃機関の制御装置。
【請求項５】
上記エアフロメータモデル計算式には実際のスロットル通過空気流量の変化に対するエアフロメータの応答時間に関する時定数が用いられており、上記モデル計算式修正手段は該時定数を修正する、請求項１〜４のいずれか１項に記載の内燃機関の制御装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、内燃機関の制御装置に関する。
【背景技術】
【０００２】
内燃機関の吸気系統を、スロットル弁、吸気管、吸気弁、エアフロメータ等の要素に分けてそれぞれの要素毎にモデル化して数式で表すことにより内燃機関の筒内充填空気量を計算によって求め、求められた筒内充填空気量を利用して内燃機関を制御する内燃機関の制御装置が知られている。
【０００３】
このようなエアモデルを用いた制御装置では、エアフロメータの出力を利用するものも多い。ところが、エアフロメータの出力はエアフロメータが設けられた吸気管内を通過する空気流量（以下、「エアフロ通過空気流量」と称す）に対して遅れが生じてしまう（これは例えば、熱線式のエアフロメータを用いた場合にはエアフロメータ自体のヒートマスによるものである）。このため、内燃機関の定常運転時にはエアフロメータにより比較的正確にエアフロ通過空気流量を検出することができるが、内燃機関の過渡運転時にはエアフロメータではエアフロ通過空気流量を正確に検出することができない。
【０００４】
このようにエアフロメータには実際のエアフロ通過空気流量に対して応答遅れが存在するため、内燃機関の過渡運転時にエアフロメータの出力を実際のエアフロ通過空気流量として利用してエアモデルにより筒内充填空気量を算出すると、算出された筒内充填空気量は実際の筒内充填空気量とは大きく異なった値となってしまい、その結果、内燃機関の運転を最適に制御することができなくなってしまう。
【０００５】
そこで、特許文献１に記載の装置では、吸気系モデルの入力側にエアフロメータの出力の応答遅れを位相進み補償で補償する応答遅れ補償要素を設け、この応答遅れ補償要素の出力を吸気系モデルに入力することとしている。すなわち、時定数を用いた数式を利用してエアフロメータの出力値からエアフロメータに存在する応答遅れ分を補償した値を他の吸気系モデルに入力することとしている。これにより、内燃機関の過渡運転時におけるエアフロメータの応答遅れの問題を改善することとしている。
【０００６】
【特許文献１】特開２００３−３１４３４７
【特許文献２】特開２００２−１３００４２
【特許文献３】特開２０００−３２０３９１
【発明の開示】
【発明が解決しようとする課題】
【０００７】
ところで、エアフロメータの応答遅れ分を適切に補償するためには、補償のために用いられる数式で利用される時定数を最適な値に設定することが必要とされる。ここで、一般に、時定数を最適な値に設定するためには、上記補償のために用いられる数式の入力値（例えば、エアフロメータの出力値）と出力値（例えば、実際のエアフロ通過空気流量）との関係を求めることが必要となる。
【０００８】
ところが、実際のエアフロ通過空気流量を求めることは不可能であるため、上述したようなエアフロメータの出力値と実際のエアフロ通過空気流量との関係を求めることはできず、よってこの関係から時定数を設定することはできない。そこで、特許文献１に記載の装置では、エアフロ通過空気流量の代わりに機関回転数を用いており、よって機関回転数とエアフロメータの出力変化量とに基づいて時定数を設定することとしている。しかしながら、機関回転数はエアフロ通過空気流量を正確に表すものではないため、機関回転数に基づいて時定数を設定しても、時定数を適切に設定することは困難である。
【０００９】
そこで、本発明の目的は、エアフロメータの応答遅れを適切に把握し、把握したエアフロメータの応答遅れを考慮して内燃機関を制御する内燃機関の制御装置を提供することにある。
【課題を解決するための手段】
【００１０】
上記課題を解決するために、第１の発明では、機関吸気通路に設けられたエアフロメータと、スロットル通過空気流量を推定するスロットル通過空気流量推定手段と、上記スロットル通過空気流量推定手段により推定されたスロットル通過空気流量の推定値に基づいてエアフロメータモデル計算式を用いてエアフロメータによって出力されるであろう出力値を予想出力値として算出するエアフロメータモデルとを具備し、上記エアフロメータの実測値と上記エアフロメータモデルによって算出されたエアフロメータの予想出力値とを用いて内燃機関を制御する、内燃機関の制御装置において、上記エアフロメータモデルによって算出されたエアフロメータの予想出力値と、上記エアフロメータの実測値との間の誤差を算出する誤差算出手段と、スロットル通過空気流量が増大しているときに上記誤差算出手段によって算出された誤差とスロットル通過空気流量が減少しているときに上記誤差算出手段によって算出された誤差とに基づいて上記エアフロメータモデル計算式を修正するモデル計算式修正手段とをさらに具備する。
スロットル通過空気流量が増大しているときのエアフロメータの予想出力値と実測値との間の誤差と、スロットル通過空気流量が減少しているときのエアフロメータの予想出力値と実測値との間の誤差とが等しい場合にはエアフロメータモデル計算式は比較的適切に実際のエアフロメータにおける応答遅れを表していると言える。第１の発明によれば、これら誤差に基づいてエアフロメータモデル計算式を修正しているため適切にエアフロメータモデル計算式を修正することができる。
なお、上記エアフロメータモデルによって算出されたエアフロメータの予想出力値は、下記本発明の実施形態における予想出力値に加えて修正出力値を含むものである。
【００１１】
第２の発明では、第１の発明において、上記モデル計算式修正手段は、上記誤差算出手段によって算出される誤差がスロットル通過空気流量が増大しているときと減少しているときとで等しくなるように上記エアフロメータモデル計算式を修正する。
【００１２】
上記課題を解決するために、第３の発明では、機関吸気通路に設けられたエアフロメータと、スロットル通過空気流量を推定するスロットル通過空気流量推定手段と、上記スロットル通過空気流量推定手段により推定されたスロットル通過空気流量の推定値に基づいてエアフロメータモデル計算式を用いてエアフロメータによって出力されるであろう出力値を予想出力値として算出するエアフロメータモデルとを具備し、上記エアフロメータの実測値と上記エアフロメータモデルによって算出されたエアフロメータの予想出力値とを用いて内燃機関を制御する、内燃機関の制御装置において、上記エアフロメータモデルによって算出されたエアフロメータの予想出力値と、上記エアフロメータの実測値との間の誤差を算出する誤差算出手段と、内燃機関が加速しているときに上記誤差算出手段によって算出された誤差と内燃機関が減速しているときに上記誤差算出手段によって算出された誤差とに基づいて上記エアフロメータモデル計算式を修正するモデル計算式修正手段とをさらに具備する。
内燃機関が加速しているときのエアフロメータの予想出力値と実測値との間の誤差と、内燃機関が減速しているときのエアフロメータの予想出力値と実測値との間の誤差とが等しい場合にはエアフロメータモデル計算式は比較的適切に実際のエアフロメータにおける応答遅れを表していると言える。第３の発明によれば、これら誤差に基づいてエアフロメータモデル計算式を修正しているため適切にエアフロメータモデル計算式を修正することができる。
【００１３】
第４の発明では、第３の発明において、上記モデル計算式修正手段は、上記誤差算出手段によって算出される誤差が内燃機関が加速しているときと減速しているときとで等しくなるように上記エアフロメータモデル計算式を修正する。
【００１４】
第５の発明では、第１〜第４のいずれか一つの発明において、上記エアフロメータモデル計算式には実際のスロットル通過空気流量の変化に対するエアフロメータの応答時間に関する時定数が用いられており、上記モデル計算式修正手段は該時定数を修正する。
【発明の効果】
【００１５】
本発明によれば、エアフロメータの応答遅れを表すエアフロメータモデル計算式を適切に修正することができるため、エアフロメータの応答遅れを適切に把握することができるようになる。
【発明を実施するための最良の形態】
【００１６】
以下、図面を参照して本発明の実施形態について説明する。図１に概略的に示した機関本体１は筒内噴射型火花点火式内燃機関を示す。しかしながら、本発明をポート噴射型火花点火式内燃機関等の別の火花点火式内燃機関や圧縮自着火式内燃機関に適用してもよい。
【００１７】
図１に示したように、本発明の実施形態では機関本体１はシリンダブロック２と、シリンダブロック２内で往復動するピストン３と、シリンダブロック２上に固定されたシリンダヘッド４とを具備する。ピストン３とシリンダヘッド４との間には燃焼室５が形成される。シリンダヘッド４には各気筒毎に吸気弁６と、吸気ポート７と、排気弁８と、排気ポート９とが配置される。さらに、図１に示したようにシリンダヘッド４の内壁面の中央部には点火プラグ１０が配置され、シリンダヘッド４内壁面周辺部には燃料噴射弁１１が配置される。また、ピストン３の頂面には燃料噴射弁１１の下方から点火プラグ１０の下方まで延びるキャビティ１２が形成されている。さらに、シリンダヘッド４には、吸気弁６の位相角およびバルブリフト量を連続的に変更可能な吸気弁制御装置１３が設けられている。
【００１８】
各気筒の吸気ポート７は吸気枝管１４を介してサージタンク１５に連結され、サージタンク１５は吸気管１６を介してエアクリーナ１７に連結される。吸気管１６内にはステップモータ１８によって駆動されるスロットル弁１９が配置される。一方、各気筒の排気ポート９は排気管２０に連結され、この排気管２０は排気浄化触媒２１を内蔵したケーシング２２に連結される。なお、以下の説明では、スロットル弁１９から吸気弁６までの吸気枝管１４、サージタンク１５、吸気管１６等の部分を吸気管部分２３と称す。
【００１９】
電子制御ユニット（ＥＣＵ）３１はディジタルコンピュータからなり、双方向性バス３２を介して相互に接続されたＲＡＭ（ランダムアクセスメモリ）３３、ＲＯＭ（リードオンリメモリ）３４、ＣＰＵ（マイクロプロセッサ）３５、入力ポート３６および出力ポート３７を具備する。サージタンク１３には、吸気管部分２３内の空気（吸気ガス）の圧力を検出するための吸気管内圧力センサ４０が設けられており、吸気管内圧力センサ４０は吸気管部分２３内の圧力に比例した出力電圧を発生し、この出力電圧が対応するＡＤ変換器３８を介して入力ポート３６に入力される。
【００２０】
また、スロットル弁１９の上流側の吸気管１６には吸気管１６内を通過して流れる吸入空気流量を検出するエアフロメータ４１が設けられており、本実施形態におけるエアフロメータ４１の構成については後述する。さらに、エアクリーナ１７近傍には吸気温度を検出する吸気温度センサ４２と、大気圧を検出する大気圧センサ４３とが設けられる。スロットル弁１９にはスロットル弁１９の開度を検出するスロットル弁開度センサ４４が設けられており、スロットル弁開度センサ４４はスロットル弁開度に対応する出力信号を発生させる。これらエアフロメータ４１、吸気温度センサ４２、大気圧センサ４３およびスロットル弁開度センサ４４はそれぞれ吸入空気流量（質量流量）、吸気温度（大気温度）、大気圧およびスロットル弁開度に対応する出力信号を発生し、この出力信号が対応するＡＤ変換器３８を介して入力ポート３６に入力される。
【００２１】
また、アクセルペダル４５にはアクセルペダル４５の踏込み量に比例した出力電圧を発生する踏込み量センサ４６が接続され、踏込み量センサ４６の出力電圧は対応するＡＤ変換器３８を介して入力ポート３６に入力される。また、クランク角センサ４７は例えばクランクシャフトが１５度回転する毎に出力パルスを発生し、この出力パルスが入力ポート３６に入力される。ＣＰＵ３５ではこのクランク角センサ４７の出力パルスから機関回転数が計算される。
【００２２】
一方、出力ポート３７は対応する駆動回路３９を介して点火プラグ１０、燃料噴射弁１１、吸気弁制御装置１３およびステップモータ１８に接続される。
【００２３】
図２は本実施形態で用いられるエアフロメータ４１の概略斜視図を、図３は図２に示したエアフロメータ４１の熱線計量部の拡大斜視図である。本実施形態では、エアフロメータ４１は、熱線式流量計とされており、図２に示したように吸気管１６内を流れる空気の一部をバイパスさせるバイパス通路と、このバイパス通路にバイパスされた吸入空気の質量流量を計測する熱線計量部４１ａと、計測された質量流量に応じた電圧を出力する信号処理部４１ｂとを有する。熱線計量部４１ａは図３に示したように白金熱線から成る吸気温度計測用抵抗４１ａ１と、この吸気温度計測用抵抗４１ａ１を信号処理部４１ｂに連結して保持するサポート部４１ａ２と、加熱用抵抗（ボビン部）４１ａ３と、この加熱用抵抗４１ａ３を信号処理部４１ｂに連結して保持するサポート部４１ａ４とを備える。信号処理部４１ｂは、吸気温度計測用抵抗４１ａ１と加熱用抵抗４１ａ３とで構成されるブリッジ回路を有し、このブリッジ回路により吸気温度計測用抵抗４１ａ１と加熱用抵抗４１ａ３との温度差を常に一定に維持するように加熱用抵抗４１ａ３に供給する電力を調整すると共に、この供給する電力を電圧に変換して出力するようになっている。エアフロメータ４１の出力電圧Ｖｇとエアフロメータ４１が配置された吸気管１６内を通過する空気の流量（以下、「エアフロ通過空気流量」と称す）との関係は図４に示した通りである。
【００２４】
ところで、内燃機関の制御装置では、内燃機関の燃焼室５において燃焼される混合気の空燃比を目標空燃比にするために、吸気弁６が閉じたときに燃焼室５内に充填されている空気の量（以下、「筒内充填空気量Ｍｃ」と称す）を推定し、推定された筒内充填空気量Ｍｃに基づいて混合気の空燃比が目標空燃比となるように燃料噴射弁１１によって内燃機関の燃焼室５（または吸気通路）に噴射する燃料の量（以下、「燃料噴射量」と称す）を定めている。したがって、内燃機関の燃焼室５において燃焼される混合気の空燃比を正確に目標空燃比とするためには、筒内充填空気量Ｍｃを正確に推定する必要がある。
【００２５】
通常、筒内充填空気量Ｍｃは、エアフロメータ等の多数のセンサ、およびこれらセンサからの出力値を引数とした多数のマップから推定される。ところが、このようにマップを用いて筒内充填空気量Ｍｃを推定する場合、推定される筒内充填空気量Ｍｃの値をより正確なものにするためには、必要なマップの数およびその引数の数が多くなる。このようにマップの数が多くなると、マップを保存するためのＥＣＵのＲＯＭを記憶容量の大きいものにしなければならず、内燃機関の制御装置の製造コストが高くなってしまう。さらに、各マップを作成するにはマップが用いられる内燃機関の形式毎に適合作業を行わなければならないが、この適合作業における測定点はマップの数およびその引数の数に応じて増大するため、マップの数およびその引数の数が多くなると適合作業の工数も増大してしまう。
【００２６】
そこで、マップを用いずに様々なモデルを用いて、数値計算により筒内充填空気量Ｍｃを算出する内燃機関の制御装置が検討されている。このような制御装置では、数値計算を多用することにより必要なマップの数を極力減らすようにしており、これにより適合作業を行う際の工数を大幅に削減しながらも、筒内充填空気量Ｍｃを正確に算出することができる。このような制御装置のうち本願出願人により提案されたものの一つに、図５に示したエアモデルを搭載した制御装置がある。図示したエアモデルＭ１は内燃機関に適用される最も単純なモデルであり、以下このエアモデルＭ１について説明する。
【００２７】
エアモデルＭ１は、図５に示したようにスロットルモデルＭ１０、吸気管モデルＭ２０、吸気弁モデルＭ３０を備える。スロットルモデルＭ１０には、スロットル弁開度センサ４４によって検出されたスロットル弁１９の開度（スロットル弁開度）θｔと、大気圧センサ４３によって検出された内燃機関周囲の大気圧（または、吸気管１６に吸入される空気の圧力）Ｐａと、吸気温度センサ４２によって検出された内燃機関周囲の大気温度（または、吸気管１６に吸入される空気の温度）Ｔａと、後述する吸気管モデルＭ２０において算出された吸気管部分２３内の圧力（以下、「吸気管内圧力Ｐｍ」と称す）とが入力され、これら入力された各パラメータの値を後述するスロットルモデルＭ１０のモデル計算式に代入することで、単位時間当たりにスロットル弁１９を通過する空気の流量（以下、「スロットル通過空気流量ｍｔ」と称す）が算出される。スロットルモデルＭ１０において算出されたスロットル通過空気流量ｍｔは、吸気管モデルＭ２０へ入力される。
【００２８】
吸気管モデルＭ２０には、スロットルモデルＭ１０において算出されたスロットル通過空気流量ｍｔと、以下で詳述する単位時間当たりに燃焼室５内に流入する空気の流量（以下、「筒内吸入空気流量ｍｃ」と称す。）とが入力され、これら入力された各パラメータの値を後述する吸気管モデルＭ２０のモデル計算式に代入することで、吸気管部分２３内に存在する空気の圧力（以下、「吸気管内圧力Ｐｍ」と称す）と吸気管部分２３内に存在する空気の温度（以下、「吸気管内温度Ｔｍ」と称す）とが算出される。吸気管モデルＭ２０において算出された吸気管内圧力Ｐｍと吸気管内温度Ｔｍは共に吸気弁モデルＭ３０へ入力され、さらに吸気管内圧力ＰｍはスロットルモデルＭ１０にも入力される。
【００２９】
吸気弁モデルＭ３０には、吸気管モデルＭ２０において算出された吸気管内圧力Ｐｍおよび吸気管内温度Ｔｍの他に大気温度Ｔａが入力され、これら入力された各パラメータの値を後述する吸気弁モデルＭ３０のモデル計算式に代入することで、筒内吸入空気流量ｍｃが算出される。算出された筒内吸入空気流量ｍｃは、筒内充填空気量Ｍｃに変換され、この筒内充填空気量Ｍｃに基づいて燃料噴射弁からの燃料噴射量が決定される。また、吸気弁モデルＭ３０において算出された筒内吸入空気流量ｍｃは吸気管モデルＭ２０に入力される。
【００３０】
図５から分かるように、エアモデルＭ１では或るモデルにおいて算出されたパラメータの値が別のモデルへの入力値として利用されるので、エアモデルＭ１全体では、実際に入力される値はスロットル弁開度θｔ、大気圧Ｐａ、および大気温度Ｔａの三つのパラメータのみであり、これら三つのパラメータから筒内充填空気量Ｍｃが算出される。
【００３１】
次に、エアモデルＭ１の各モデルＭ１０〜Ｍ３０について説明する。
スロットルモデルＭ１０では、大気圧Ｐａ、大気温度Ｔａ、吸気管内圧力Ｐｍ、スロットル弁開度θｔから、下記式（１）に基づいてスロットル通過空気流量ｍｔが算出される。ここで、式（１）におけるμｔはスロットル弁における流量係数で、スロットル弁開度θｔの関数であり、よって図６に示したようなマップから定まる。また、Ａｔはスロットル弁の開口断面積を示し、スロットル弁開度θｔの関数であり、図７に示したようなマップから定まる。また、Ｒは気体定数に関する定数であり、実際には気体定数を１mol当たりの気体（空気）の質量Ｍlmolで除算した値である。
【数１】
2007231840.tif 000003
【００３２】
また、Φ（Ｐｍ／Ｐａ）は下記式（２）に示した関数であり、この式（２）におけるκは比熱比（一定値とする）である。この関数Φ（Ｐｍ／Ｐａ）は図８に示したようなグラフに表すことができるので、このようなグラフをマップとしてＥＣＵ３１のＲＯＭ３４に保存し、実際には式（２）を用いて計算するのではなくマップからΦ（Ｐｍ／Ｐａ）の値を求めるようにしてもよい。
【数２】
2007231840.tif 000004
【００３３】
これらスロットルモデルＭ１０の式（１）および式（２）は、スロットル弁１９上流の気体の圧力を大気圧Ｐａ、スロットル弁１９上流の気体の温度を大気温度Ｔａ、スロットル弁１９の下流の気体の圧力を吸気管内圧力Ｐｍとして、図９に示したようなスロットル弁１９のモデルに対して、質量保存則、エネルギ保存則および運動量保存則を適用し、さらに気体の状態方程式、比熱比の定義式、およびマイヤーの関係式を利用することによって得られる。
【００３４】
吸気管モデルＭ２０では、スロットル通過空気流量ｍｔ、筒内吸入空気流量ｍｃ、および大気温度Ｔａから、下記式（３）および式（４）に基づいて吸気管内圧力Ｐｍおよび吸気管内温度Ｔｍが算出される。なお、式（３）および式（４）におけるＶｍはスロットル弁１９から吸気弁６までの吸気枝管１４、サージタンク１５、吸気管１６等の部分（吸気管部分２３）の容積に等しい定数である。
【数３】
2007231840.tif 000005
【数４】
2007231840.tif 000006
【００３５】
ここで、吸気管モデルＭ２０について図１０を参照して説明する。吸気管部分２３の総気体量（総空気量）をＭとすると、総気体量Ｍの時間的変化は、吸気管部分２３に流入する気体の流量、すなわちスロットル通過空気流量ｍｔと、吸気管部分２３から流出する気体の流量、すなわち筒内吸入空気流量ｍｃとの差に等しいため、質量保存則により下記式（５）が得られ、この式（５）および吸気管部分２３における気体の状態方程式（Ｐｍ・Ｖｍ＝Ｍ・Ｒ・Ｔｍ）より、式（３）が得られる。
【数５】
2007231840.tif 000007
【００３６】
また、吸気管部分２３の気体のエネルギＭ・Ｃｖ・Ｔｍの時間的変化量は、吸気管部分２３に流入する気体のエネルギと吸気管部分２３から流出する気体のエネルギとの差に等しい。このため、吸気管部分２３に流入する気体の温度を大気温度Ｔａ、吸気管部分２３から流出する気体の温度を吸気管内温度Ｔｍとすると、エネルギ保存則により下記式（６）が得られ、この式（６）および上記気体の状態方程式より、式（４）が得られる。
【数６】
2007231840.tif 000008
【００３７】
吸気弁モデルＭ３０では、吸気管内圧力Ｐｍ、吸気管内温度Ｔｍ、および大気温度Ｔａから、下記式（７）に基づいて、筒内吸入空気流量ｍｃが算出される。なお、式（７）におけるａ、ｂは、機関回転数Ｎｅから、さらに吸気弁６の位相角（バルブタイミング）および作用角を変更できる可変動弁機構を備えた内燃機関の場合には吸気弁６の位相角、作用角から定まる値である。
【数７】
2007231840.tif 000009
【００３８】
上述した吸気弁モデルＭ３０について図１１を参照して説明する。一般に、吸気弁６が閉じたときに燃焼室５内に充填されている空気の量である筒内充填空気量Ｍｃは、吸気弁６が閉弁するとき（吸気弁閉弁時）に確定し、吸気弁６閉弁時の燃焼室５内の圧力に比例する。また、吸気弁６閉弁時の燃焼室５内の圧力は吸気弁６上流の気体の圧力、すなわち吸気管内圧力Ｐｍと等しいとみなすことができる。したがって、筒内充填空気量Ｍｃは、吸気管内圧力Ｐｍに比例すると近似することができる。
【００３９】
ここで、一定時間（例えば、クランク角７２０°分）当たりに吸気管部分２３から流出する全空気流量を平均化したもの、または一定時間（例えば、クランク角７２０°分）当たりに吸気管部分２３から全ての気筒の燃焼室５に吸入される空気量を上記一定時間で除算したものを筒内吸入空気流量ｍｃ（以下で詳述する）とすると、筒内充填空気量Ｍｃが吸気管内圧力Ｐｍに比例することから、筒内吸入空気流量ｍｃも吸気管内圧力Ｐｍに比例すると考えられる。このことから、理論および経験則に基づいて、上記式（７）が得られる。なお、式（７）における値ａは比例係数であり、値ｂは燃焼室５内に残存していた既燃ガスを表す値（排気弁８閉弁時に燃焼室５内に残る既燃ガス量を後述する時間ΔＴ180°で除算したものに相当）である。また、実際の運転では過渡運転時に吸気管内温度Ｔｍが大きく変化する場合があるため、これに対する補正として理論および経験則に基づいて導かれたＴａ／Ｔｍが乗算されている。
【００４０】
ここで、筒内吸入空気流量ｍｃについて、図１２を参照して内燃機関が４気筒である場合について説明する。なお、図１２は横軸がクランクシャフトの回転角度、縦軸が単位時間当たりに吸気管部分２３から燃焼室５に実際に流入する空気の流量である。図１２に示したように、４気筒の内燃機関では、吸気弁６が例えば１番気筒、３番気筒、４番気筒、２番気筒の順に開弁し、各気筒に対応する吸気弁６の開弁量に応じて吸気管部分２３から各気筒の燃焼室５内へ空気が流入する。例えば、吸気管部分２３から各気筒の燃焼室５内に流入する空気の流量の変位は図１２に破線で示した通りであり、これを総合して吸気管部分２３から全気筒の燃焼室５に流入する空気の流量は図１２に実線で示した通りである。また、例えば１番気筒への筒内充填空気量Ｍｃは図１２に斜線で示した通りである。
【００４１】
これに対して、実線で示した吸気管部分２３から全ての気筒の燃焼室５に流入する空気の流量を平均化したものが筒内吸入空気流量ｍｃであり、図中に一点鎖線で示す。そして、この一点鎖線で示した筒内吸入空気流量ｍｃに、４気筒の場合にはクランクシャフトが１８０°（すなわち、４ストローク式内燃機関において１サイクル中にクランクシャフトが回転する角度７２０°を気筒数で割った角度）回転するのにかかる時間ΔＴ180°を乗算したものが筒内充填空気量Ｍｃとなる。したがって、吸気弁モデルＭ３０で算出された筒内吸入空気流量ｍｃにΔＴ180°を乗算することで、筒内充填空気量Ｍｃが算出される（Ｍｃ＝ｍｃ・ΔＴ180°）。より詳細には、筒内充填空気量Ｍｃが吸気弁閉弁時の圧力に比例することを考慮して、吸気弁閉弁時の筒内吸入空気流量ｍｃにΔＴ180°を乗算したものが筒内充填空気量Ｍｃとされる。
【００４２】
次に、上記エアモデルＭ１を内燃機関の制御装置に実装して、実際に筒内充填空気量Ｍｃを算出する場合について説明する。筒内充填空気量ＭｃはエアモデルＭ１を用いて、上記式（１）、式（３）、式（４）、および式（７）を解くことにより表される。この場合、ＥＣＵ３１で処理するために、これらの式を離散化する必要がある。時刻ｔ、計算間隔Δｔを用いて式（１）、式（３）、式（４）、および式（７）を離散化すると、それぞれ下記式（８）、式（９）、式（１０）、および式（１１）が得られる。なお、吸気管内温度Ｔｍ（ｔ＋Δｔ）は、式（９）および式（１０）によってそれぞれ算出されたＰｍ／Ｔｍ（ｔ＋Δｔ）およびＰｍ（ｔ＋Δｔ）から、式（１２）によって算出される。
【数８】
2007231840.tif 000010
【数９】
2007231840.tif 000011
【数１０】
2007231840.tif 000012
【数１１】
2007231840.tif 000013
【数１２】
2007231840.tif 000014
【００４３】
このようにして実装されたエアモデルＭ１では、スロットルモデルＭ１０の式（８）で算出された時刻ｔにおけるスロットル通過空気流量ｍｔ（ｔ）と、吸気弁モデルＭ３０の式（１１）で算出された時刻ｔにおける筒内吸入空気流量ｍｃ（ｔ）とが、吸気管モデルＭ２０の式（９）および式（１０）に代入され、これにより時刻ｔ＋Δｔにおける吸気管内圧力Ｐｍ（ｔ＋Δｔ）および吸気管内温度Ｔｍ（ｔ＋Δｔ）が算出される。次いで、算出されたＰｍ（ｔ＋Δｔ）およびＴｍ（ｔ＋Δｔ）は、スロットルモデルＭ１０および吸気弁モデルＭ３０の式（８）および式（１１）に代入され、これにより時刻ｔ＋Δｔにおけるスロットル通過空気流量ｍｔ（ｔ＋Δｔ）および筒内吸入空気流量ｍｃ（ｔ＋Δｔ）が算出される。そして、このような計算を繰り返すことによって、スロットル弁開度θｔ、大気圧Ｐａ、および大気温度Ｔａから、任意の時刻ｔにおける筒内吸入空気流量ｍｃが算出され、算出された筒内吸入空気流量ｍｃに上記時間ΔＴ180°を乗算することで、任意の時刻ｔにおける筒内充填空気量Ｍｃが算出される。
【００４４】
なお、内燃機関の始動時には、すなわち時刻ｔ＝０においては、吸気管内圧力Ｐｍは大気圧と等しい（Ｐｍ（０）＝Ｐａ）とされ、吸気管内温度Ｔｍは大気温度と等しい（Ｔｍ（０）＝Ｔａ）とされて、各モデルＭ１０〜Ｍ３０における計算が開始される。
【００４５】
なお、上記エアモデルＭ１では、大気温度Ｔａおよび大気圧Ｐａが一定であるとしているが、時刻によって変化する値としてもよく、例えば、大気温度を検出するための吸気温度センサによって時刻ｔにおいて検出された値を大気温度Ｔａ（ｔ）、大気圧を検出するための大気圧センサによって時刻ｔにおいて検出された値を大気圧Ｐａ（ｔ）として上記式（８）、式（１０）、および式（１１）に代入するようにしてもよい。
【００４６】
ところで、上述したように上記エアモデルＭ１は内燃機関に適用される最も単純なモデルであり、上記エアモデルＭ１によって算出される筒内充填空気量Ｍｃは実際の筒内充填空気量に対して誤差が生じ易い。そこで、図１３に示した後述するエアモデルＭ２では、エアモデルによって筒内充填空気量をより高い精度で算出すべく、上記スロットル弁開度θｔ、大気圧Ｐａ、および大気温度Ｔａの三つのパラメータに加えて、エアフロメータ４１の出力値ＡＦＭに基づいて筒内充填空気量Ｍｃを算出することとしている。以下、このエアモデルＭ２について説明する。
【００４７】
エアモデルＭ２は、上記スロットルモデルＭ１０、吸気管モデルＭ２０、および吸気弁モデルＭ３０に加えて、電子制御スロットルモデルＭ４０、およびエアフロメータモデル（ＡＦＭモデル）Ｍ５０を備える。電子制御スロットルモデルＭ４０には、アクセルペダル４５の踏込み量に比例した出力電圧を発生する踏込み量センサ４６の出力値Ａｃｃｐが入力され、これに基づいて現時刻から所定時間Ｔ０先の時刻においてスロットル弁１９が到達すると予想されるスロットル弁開度（以下、「所定時間後のスロットル弁開度」と称す）θｔｆが算出される。
【００４８】
また、ＡＦＭモデルＭ５０には、スロットルモデルＭ１０によって算出されたスロットル通過空気流量（以下、「スロットル通過空気流量の推定値」と称す）ｍｔが入力され、この入力された値を後述するＡＦＭモデルＭ５０のモデル計算式に代入することで、実際にスロットル弁１９の配置された吸気管１６内を上記スロットル通過空気流量の推定値ｍｔだけ空気が流れていると仮定した場合にエアフロメータ４１が出力するであろう出力値（以下、「予想出力値」と称す）ＡＦＭｍｔが算出される。なお、以下の説明では、各モデルＭ１０等によって算出された値については「推定値」と表記することとする。
【００４９】
すなわち、エアフロメータ４１の出力は、実際のエアフロ通過空気流量（エアフロ通過空気流量とスロットル通過空気流量はほぼ同一であると考えられるため、以下ではスロットル通過空気流量として説明する）に対して固有の応答特性に基づく応答遅れを有しており、したがって現在の実際のスロットル通過空気流量のみに対応した値となっているわけではない。ＡＦＭモデルＭ５０は、この応答特性をシュミレートしたモデルであり、上記エアフロメータ４１の応答遅れを考慮した上で、予想出力値を算出する。
【００５０】
ここで、エアモデルＭ２は、図１３に示したように、三つのモデルブロックＭ２'、Ｍ２''、Ｍ２'''に分けて考えることができる。モデルブロックＭ２'では、踏込み量センサ４６の出力値Ａｃｃｐに基づいて電子制御スロットルモデルＭ４０により所定時間後のスロットル弁開度の推定値θｔｆが算出されると共に、算出された所定時間後のスロットル弁開度の推定値θｔｆに基づいてスロットルモデルＭ１０'により所定時間後のスロットル通過空気流量の推定値ｍｔｆが算出される。そして、算出された所定時間後のスロットル通過空気流量の推定値ｍｔｆに基づいて吸気弁モデルＭ２０'により所定時間後の吸気管内圧力の推定値Ｐｍｆが算出される。算出された所定時間後の吸気管内圧力の推定値ＰｍｆはスロットルモデルＭ１０'及び吸気弁モデルＭ３０'に入力される。吸気弁モデルＭ３０'では、所定時間後の吸気管内圧力の推定値Ｐｍｆに基づいて所定時間後の筒内吸入空気流量の推定値ｍｃｆが算出され、吸気管モデルＭ２０'に入力される。
【００５１】
したがって、エアモデルＭ２のモデルブロックＭ２'全体では、所定時間後のスロットル弁開度θｔｆを入力することにより、所定時間後の筒内吸入空気流量の推定値ｍｃｆおよび所定時間後の吸気管内圧力の推定値Ｐｍｆが算出される。このように、所定時間後の筒内吸入空気流量の推定値ｍｃｆ等を算出することにより、実際の空燃比を目標空燃比に一致させ易くなる。すなわち、内燃機関の過渡運転時には筒内吸入空気流量は刻々と変化しており、よって現在の筒内吸入空気流量（或いは、現在の筒内吸入空気流量から算出された筒内充填空気量）に基づいて燃料噴射量を算出しても実際に噴射するときには既に筒内吸入空気流量が変化してしまっている場合があり、現在の筒内吸入空気流量を正確に算出することができても結果的に実際の空燃比が目標空燃比と異なったものとなってしまうことがある。これに対して、所定時間後の筒内吸入空気流量に基づいて燃料噴射量を算出すれば、実際に噴射するときにその筒内吸入空気流量となっており、よって実際の空燃比を目標空燃比に一致させ易くなる。
【００５２】
ただし、モデルブロックＭ２'では電子制御スロットルモデルＭ４０が存在する点を除いて上記エアモデルＭ１と同様な計算が行われているため、モデルブロックＭ２'で算出された所定時間後の筒内吸入空気流量の推定値ｍｃｆに基づいて筒内充填空気量Ｍｃを算出すると、上述したように実際の筒内充填空気量に対して誤差が生じ易い。
【００５３】
そこで、エアモデルＭ２では、エアフロメータ４１の出力に基づいて、モデルブロックＭ２'において算出された吸気管内圧力の推定値Ｐｍｆを修正し、修正した吸気管内圧力Ｐｍｆｉｎに基づいて筒内充填空気量Ｍｃを算出することとしている。
【００５４】
すなわち、モデルブロックＭ２''において、スロットル弁開度センサ４４によって検出された現在のスロットル弁開度θｔｃに基づいて上記エアモデルＭ１と同様な計算が行われる。すなわち、スロットルモデルＭ１０''には現在のスロットル弁開度θｔｃが入力されると共に現在のスロットル通過空気流量の推定値ｍｔｃが出力され、吸気管モデルＭ２０''には現在のスロットル通過空気流量の推定値ｍｔｃが入力されると共に現在の吸気管内圧力の推定値Ｐｍｃが出力され、この吸気管内圧力ＰｍｃがスロットルモデルＭ１０''及び吸気弁モデルＭ３０''に入力される。吸気弁モデルＭ３０''では現在の筒内吸入空気流量の推定値ｍｃｃが算出され、吸気管モデルＭ２０''に入力される。
【００５５】
そして、上記計算において算出された現在のスロットル通過空気流量の推定値ｍｔｃに基づいてＡＦＭモデルＭ５０によりエアフロメータ４１の応答遅れを考慮した予想出力値ＡＦＭｍｔが算出され、このようにして算出された予想出力値ＡＦＭｍｔがスロットル通過空気流量として吸気管モデルＭ２０'''に入力される。そして、吸気管モデルＭ２０'''および吸気弁モデルＭ３０'''により上述したのと同様な方法で吸気管内圧力の推定値Ｐｍｍｄｌが算出される。このようにして算出された吸気管内圧力の推定値Ｐｍｍｄｌは、現在のスロットル弁開度θｔｃに基づいてモデルブロックＭ２'と同様なエアモデルを利用して算出されたエアフロメータ４１によって出力されるであろう出力値にスロットル通過空気流量が等しいと仮定した場合における吸気管内圧力を示している（以下、「ＡＦＭモデルを用いた吸気管内圧力の推定値」と称す）。
【００５６】
一方、モデルブロックＭ２'''では、実際のエアフロメータ４１の出力値ＡＦＭがスロットル通過空気流量として吸気管モデルＭ２０''''に入力され、吸気管モデルＭ２０''''により吸気管内圧力の推定値Ｐｍａｆｍが出力される。吸気弁モデルＭ３０''''には吸気管内圧力の推定値Ｐｍａｆｍが入力され、筒内吸入空気流量ｍｃａｆｍが出力され、吸気管モデルＭ２０''''に再び入力される。このようにして算出された吸気管内圧力の推定値Ｐｍａｆｍは、スロットル通過空気流量がエアフロメータ４１の出力値に等しいと仮定した場合における吸気管内圧力を示している（以下、「エアフロメータの出力に基づく吸気管内圧力の推定値」と称す）。
【００５７】
したがって、ＡＦＭモデルを用いた吸気管内圧力の推定値Ｐｍｍｄｌとエアフロメータの出力に基づく吸気管内圧力の推定値Ｐｍａｆｍとの差分は、エアモデルによって算出された現在の吸気管内圧力とエアフロメータ４１に基づいて算出された吸気管内圧力との誤差を表している。したがって、この差分をモデルブロックＭ２'において算出された吸気管内圧力の推定値Ｐｍｆに加算した吸気管内圧力の推定値Ｐｍｆｉｎ（＝Ｐｍｆ＋Ｐｍａｆｍ−Ｐｍｍｄｌ）は実際の吸気管内圧力とほぼ一致する。そして、この吸気管内圧力の推定値Ｐｍｆｉｎが吸気管モデルＭ３０に入力され、吸気管モデルＭ３０によって算出された筒内吸入空気流量の推定値ｍｃｆｉｎに基づいて筒内充填空気量の推定値Ｍｃが算出され、この筒内充填空気量の推定値Ｍｃに基づいて燃料噴射量が算出される。
【００５８】
次に、電子制御スロットルモデルＭ４０及びＡＦＭモデルＭ５０について説明する。
電子制御スロットルモデルＭ４０は、アクセルペダル４５の踏込み量に基づいて現時刻から所定時間Ｔ０後のスロットル弁開度θｔｆを算出するモデルである。本実施形態においては、踏込み量センサ４６によって検出されたアクセルペダル踏込み量Ａｃｃｐと、図１４に示したマップとに基づいて暫定的な目標スロットル弁開度θｒ１が求められ、この暫定的な目標スロットル弁開度θｒ１を所定時間Ｔ（例えば、６４ｍｓｅｃ）だけ遅延させた値が最終的な目標スロットル弁開度θｒとして決定される。そして、ＥＣＵ３１は、実際のスロットル弁開度が目標スロットル弁開度θｒとなるようにステップモータ１８に対して駆動信号を送出する。
【００５９】
このように、現時刻の目標スロットル弁開度θｒは、現時刻から所定時間Ｔだけ前の時刻におけるアクセルペダル踏込み量Ａｃｃｐに応じて決定された暫定的な目標スロットル弁開度θｒ１と等しいから、現時刻から所定時間Ｔ０（０≦Ｔ０≦Ｔ）だけ後の時刻における目標スロットル弁開度θｒは現時刻から時間（Ｔ−Ｔ０）前における暫定的な目標スロットル弁開度θｒ１と等しい。また、ステップモータ１８の作動遅れ時間を無視すれば、目標スロットル弁開度θｒは実際のスロットル弁開度と等しい。このような考えに基づき、電子制御スロットルモデルＭ４０では、現時刻から時間（Ｔ−Ｔ０）前における暫定的な目標スロットル弁開度θｒ１を現時刻から所定時間Ｔ０だけ後の時刻ｔにおけるスロットル弁開度（所定時間後のスロットル弁開度）θｔｆとして推定する。なお、ステップモータ１８の作動遅れ時間を考慮に加えて、所定時間後のスロットル弁開度θｔｆを設定してもよい。
【００６０】
次に、ＡＦＭモデルＭ５０について具体的に説明する。エアフロメータ４１では、上述したように吸気温度計測用抵抗４１ａ１と加熱用抵抗（ボビン部）４１ａ３との温度差を常に一定に維持するように加熱用抵抗４１ａ３に供給する電力を調整するとともに、その供給電力に基づいてエアフロ通過空気流量を算出するようにしている。ここで、この供給電力はボビン部４１ａ３およびサポート部４１ａ４から吸気管１６内を通過する空気への放熱量を示すことから、エアフロメータ４１はボビン部４１ａ３およびサポート部４１ａ４からの放熱量に基づいてエアフロ通過空気流量を算出すると言い換えることができる。
【００６１】
ここで、ボビン部４１ａ３は、より詳細には円筒状のセラミックスボビンに白金熱線を巻回し、その外周にガラスをコーティングすることにより形成される。このため、白金熱線から周囲の空気への放熱は白金熱線と周囲の空気との間にガラス層が介在することにより遅れてしまう。したがって、吸気管１６内を通過する空気の流量が急激に増大したような場合であっても、ボビン部４１ａ３からの単位時間当たりの放熱量（以下、単に「放熱量」と称す）は直ぐには増大せず、或る程度の遅れをもって増大することになる。換言すると、ボビン部４１ａ３からの放熱量は実際のエアフロ通過空気流量に対して応答遅れが存在する。また、同様なことがサポート部４１ａ４からの放熱量についても言え、サポート部４１ａ４からの放熱量は実際のエアフロ通過空気流量に対して応答遅れが存在する。
【００６２】
この応答遅れは、一次遅れに近似することができることがわかっており、ボビン部４１ａ３の放熱量の応答遅れは次式（１３）で表される。ここで、式（１３）におけるωbはエアフロメータ４１の出力値から換算したボビン部４１ａ３の放熱量、すなわち放熱遅れが生じた結果ボビン部４１ａ３の白金熱線から実際に放熱される放熱量（以下、「遅れ放熱量」と称す）を示している。また、式（１３）におけるＷbは、応答遅れを補償した放熱量、すなわち放熱遅れが生じないと仮定した場合におけるボビン部４１ａ３の白金熱線から放熱される放熱量（以下、「完全放熱量」と称す）を示している。すなわち、完全放熱量Ｗbは、内燃機関が定常運転を行っているときにおける放熱量に等しく、基本的にエアフロメータ４１近傍の吸気管１６内を通過する空気流量のみの関数である。エアフロメータ４１近傍の吸気管１６内を通過する空気流量はスロットル通過空気流とほぼ等しいため、完全放熱量Ｗbはスロットル通過空気流量の関数と考えることができ、スロットル通過空気流量とボビン部４１ａ３からの完全放熱量Ｗbとの関係は図１５のように表すことができる。さらに、式（１３）におけるτbは、ボビン部４１ａ３からの放熱における一次遅れの時定数であり、その算出方法については後述する。
【数１３】
2007231840.tif 000015
【００６３】
同様に、サポート部４１ａ４の放熱量の応答遅れは次式（１４）で表される。サポート部４１ａ４からの完全放熱量Ｗbもスロットル通過空気流量の関数と考えることができ、スロットル通過空気流量とサポート部４１ａ４からの完全放熱量Ｗsとの関係は図１５のように表すことができる。また、式（１４）におけるτsは、サポート部４１ａ４からの放熱における一次遅れの時定数である。
【数１４】
2007231840.tif 000016
【００６４】
式（１３）および式（１４）によりボビン部４１ａ３およびサポート部４１ａ４からの遅れ放熱量ωb、ωsが算出される。エアフロメータ４１では、放熱遅れが生じた結果ボビン部４１ａ３およびサポート部４１ａ４から実際に放熱される放熱量、すなわち遅れ放熱量に応じて出力電圧が変化するため、式（１３）および式（１４）によって算出された遅れ放熱量の和（ωb＋ωs）に基づいてエアフロメータ４１が出力するであろう出力値（予想出力値）ＡＦＭｍｔが算出される。本実施形態では、エアフロメータ４１における遅れ放熱量の和（ωb＋ωs）とエアフロメータ４１の予想出力値ＡＦＭｍｔとの関係を予め実験的にまたは計算によって求め、図１６に示したようなマップとしてＥＣＵ３１のＲＯＭ３４に記憶させておく。そして、式（１３）および式（１４）から算出された遅れ放熱量の和ωb＋ωsに基づいて上記マップを用いてエアフロメータ４１の予想出力値ＡＦＭｍｔが算出される。なお、上記マップを用いた遅れ放熱量の和ωb＋ωsに基づくエアフロメータ４１の予想出力値ＡＦＭｍｔの算出は下記式（１５）のように表すことができる。
【数１５】
2007231840.tif 000017
【００６５】
なお、ボビン部４１ａ３の一次遅れの時定数τbは下記式（１６）によって、またサポート部４１ａ４の一次遅れの時定数τsは下記式（１７）によって算出される。式（１６）および式（１７）において、ｕは、エアフロメータ４１の検出部における流路、すなわちエアフロメータ４１のバイパス通路における単位断面積当たりの空気流量である。単位断面積当たりの空気流量ｕはエアフロメータ４１の出力値ＡＦＭに基づいてマップによりまたは所定の計算式により算出される。また、ｋb、ｋs、ｍb、ｍsは予め実験または計算によって求められる定数であり、ｋb、ｍbがボビン部４１ａ３についての定数、ｋs、ｍsがサポート部４１ａ４についての定数をそれぞれ示している。ボビン部４１ａ３とサポート部４１ａ４とでは応答遅れの度合いが異なるので、ボビン部４１ａ３とサポート部４１ａ４とを分離して時定数を設定することによってスロットル通過空気流量から予想出力値を算出するための算出精度を向上させることとしている。
【数１６】
2007231840.tif 000018
【数１７】
2007231840.tif 000019
【００６６】
次に、上記ＡＦＭモデルＭ５０を内燃機関の制御装置に実装して、実際にスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔからエアフロメータ４１の予想出力値ＡＦＭｍｔを算出する場合について説明する。エアフロメータ４１の予想出力値ＡＦＭｍｔはＡＦＭモデルＭ５０を用いて、上記式（１３）、式（１４）を解くことにより表される。この場合、ＥＣＵ３１で処理するために、これらの式を離散化する必要がある。時刻ｔ、計算間隔Δｔを用いて式（１３）、式（１４）を離散化すると、それぞれ下記式（１８）、式（１９）が得られる。
【数１８】
2007231840.tif 000020
【数１９】
2007231840.tif 000021
【００６７】
このようにして実装されたＡＦＭモデルＭ５０では、時刻ｔにおけるボビン部４１ａ３からの遅れ放熱量ωb（ｔ）、サポート部４１ａ４からの遅れ放熱量ωs（ｔ）およびスロットルモデルＭ１０の式（８）で算出された時刻ｔにおけるスロットル通過空気流量の推定値ｍｔ（ｔ）が式（１８）および式（１９）に代入され、これにより時刻ｔ＋Δｔにおけるボビン部４１ａ３からの遅れ放熱量ωb（ｔ＋Δｔ）およびサポート部４１ａ４からの遅れ放熱量ωs（ｔ＋Δｔ）が算出される。そして、これら遅れ放熱量ωb（ｔ＋Δｔ）およびωs（ｔ＋Δｔ）を用いて式（１５）により時刻ｔ＋Δｔにおけるエアフロメータ４１の予想出力値ＡＦＭｍｔ（ｔ＋Δｔ）が算出される。
【００６８】
ところで、上記したように図５および図１３に記載の装置ではスロットル弁開度θｔ、大気温度Ｔａおよび大気圧Ｐａ等を入力するモデル計算式を用いて機関始動時からの逐次計算により筒内充填吸気量が算出される。このため、エアモデルの各要素モデルにモデル化誤差が生じていると、算出される筒内充填吸気量は不正確になり燃料噴射量の正確な算出ができなくなる。
【００６９】
一方、エアモデルの各要素モデル（スロットルモデルＭ１０、吸気管モデルＭ２０、吸気弁モデルＭ３０、電子制御スロットルモデルＭ４０、ＡＦＭモデルＭ５０）を考えると、モデル化そのものには誤差がない場合であっても、使用による特性の変化や製作公差によるばらつきなどによりモデル化誤差が生じている場合がある。そこで、本実施形態では各要素モデルについて修正を行うこととしており、特に本実施形態ではスロットルモデルＭ１０とＡＦＭモデルＭ５０とを修正するようにしている。
【００７０】
以下、本実施形態のモデル修正方法について説明する。本実施形態では、現在の実際のエアフロ通過空気流量、スロットル弁開度、吸気管内圧力を、それぞれエアフロメータ４１、スロットル弁開度センサ４４、吸気圧センサ４０で実測し、これらの値を用いて図５の筒内充填吸気量計算とは独立してスロットルモデルＭ１０、ＡＦＭモデルＭ５０等についての修正操作を行う。
【００７１】
まず、スロットルモデルＭ１０の修正操作について説明する。図１７は、スロットルモデルＭ１０の修正操作を説明するブロック図である。図１７に示したように、スロットルモデルの修正操作を行う際には、ＡＦＭモデルＭ５０が用いられる。
【００７２】
本実施形態では、吸気圧センサ４０とスロットル弁開度センサ４４とで実測した現在の吸気管内圧力Ｐｍとスロットル弁開度θｔ（および大気圧Ｐａ、大気温度Ｔａ）とを用いてスロットルモデルＭ１０（図５）によって算出されるスロットル通過空気流量の推定値ｍｔが、エアフロメータ４１の出力から算出した現在のスロットル通過空気流量に一致するようにスロットルモデルＭ１０の式（１）を修正する。
【００７３】
すなわち、図１７において、実測値Ｐｍとスロットル弁開度θｔとを用いてスロットルモデルＭ１０の式（１）によって算出されたスロットル通過空気流量の推定値ｍｔをＡＦＭモデルＭ５０に入力すると、実際のスロットル通過空気流量が上記スロットル通過空気流量の推定値ｍｔとなっていたと仮定した場合におけるエアフロメータ４１の予想出力値ＡＦＭｍｔに変換される。
【００７４】
したがって、スロットルモデルＭ１０にモデル化誤差が含まれていない場合、すなわち実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとが一致している場合には、実際のエアフロメータ４１の出力値ＡＦＭとＡＦＭモデルＭ５０の予想出力値ＡＦＭｍｔとが一致するはずであり、一致していない場合には実際のスロットル通過空気流量とスロットル通過空気流量の推定値ｍｔとが一致していないと考えられ、よってスロットルモデルＭ１０にモデル化誤差が生じていると考えられる。
【００７５】
したがって、出力値ＡＦＭと予想出力値ＡＦＭｍｔとが一致しない場合には、スロットルモデルＭ１０に生じているモデル化誤差を補償するために、予想出力値ＡＦＭｍｔがエアフロメータ４１の出力値（実測値）ＡＦＭに一致するようにスロットルモデルＭ１０のモデル計算式を修正する必要がある。
【００７６】
ここで、本実施形態では、上述したようにスロットルモデルＭ１０の計算式として、下記式（２０）を用いている。
【数２０】
2007231840.tif 000022
上記式（２０）のパラメータのうち流量係数μｔは製品毎のばらつきや使用による特性変化の影響を受けやすい。したがって、本実施形態では、スロットルモデルＭ１０の式（２０）により算出されたスロットル通過空気流量の推定値ｍｔをＡＦＭモデルＭ５０の式（１３）および式（１４）に代入して算出された予想出力値ＡＦＭｍｔがエアフロメータ４１の出力値ＡＦＭに一致するように流量係数μｔの値を修正することとしている。これにより、スロットルモデルＭ１０が修正され、スロットルモデルＭ１０を用いた計算結果をエアフロメータ４１の出力値と一致させることができるようになり、よってスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔを実際のスロットル通過空気流量に一致させることができるようになる。なお、当然ながらスロットルモデルＭ１０の上記修正結果は図５筒内充填吸気量Ｍｃの計算にも同時に反映され、筒内充填吸気量Ｍｃの推定精度が向上するようになる。
【００７７】
このような流量係数μｔの値の修正方法としては、例えば、下記式（２１）によるものが考えられる。なお、式（２１）におけるμｔ’は修正後の流量係数である。
【数２１】
2007231840.tif 000023
【００７８】
ところが、式（２１）により流量係数μｔの値を修正すると、内燃機関が定常運転を行っているときには適切に修正することができるが、内燃機関が過渡運転を行っているときには適切に修正することができない。以下、このことについて図１８を参照して簡単に説明する。
【００７９】
図１８は、スロットル弁開度θｔと空気流量に関する各種パラメータの値とのタイムチャートである。図中、時刻ｔ0はスロットル弁開度θｔがθｔ0となったとき、時刻ｔ1はスロットル弁開度θｔがθｔ1となったときを示している。また、図１８は、実際のスロットル通過空気流量（図１８中の実線）と、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔ（図１８中の一点鎖線）とが時刻ｔ0前まではずれているが、時刻ｔ0以降はほぼ一致している場合のタイムチャートを示している。すなわち、図１８は、時刻ｔ0前、すなわちスロットル弁開度θｔがθｔ0よりも小さいときには、流量係数μｔの値に誤差が生じていてスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔが実際のスロットル通過空気流量からずれてしまっていると共に、時刻ｔ0以降、すなわちスロットル弁開度θｔがθｔ0以上であるときには、流量係数μｔの値に誤差がほとんどなく、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔが実際のスロットル通過空気流量にほぼ一致している場合のタイムチャートを示している。
【００８０】
上述したように、エアフロメータ４１には放熱遅れに基づく応答遅れがあるため、図１８に示したように実際のスロットル通過空気流量が変化している場合に、エアフロメータ４１の出力値ＡＦＭ（図１８中の破線）には実際のスロットル通過空気流量に対して遅れが生じる。特に、エアフロメータ４１の出力値ＡＦＭを実際のスロットル通過空気流量に対して一次遅れで近似することができることを考えると、エアフロメータ４１の現在の出力値ＡＦＭは現在の実際のスロットル通過空気流量のみに応じて変化するわけではなく、過去の実際のスロットル通過空気流量の履歴にも応じて変化する。
【００８１】
また、上述したようにＡＦＭモデルＭ５０ではエアフロメータ４１の放熱遅れに基づく応答遅れが考慮されているため、図１８に示したように実際のスロットル通過空気流量が変化している場合には、ＡＦＭモデルＭ５０によって算出される予想出力値ＡＦＭｍｔにもスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔに対して遅れが生じている。
【００８２】
このため、図１８に示したように、時刻ｔ0以降においてスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔが実際のスロットル通過空気流量にほぼ一致しているにも関わらず、時刻ｔ0以降においてエアフロメータ４１の出力値ＡＦＭとＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔとは一致していない。したがって、例えば、出力値ＡＦＭと予想出力値ＡＦＭｍｔとを用いて式（２１）によってスロットルモデルＭ１０の流量係数μｔを修正しようとすると、時刻ｔ0以降において流量係数μｔには誤差が無いにも関わらず、流量係数μｔが修正されてしまうことになる。
【００８３】
或いは、流量係数μｔの値の修正方法としては、例えば、下記式（２２）によるものが考えられる。このように、積分を用いることにより、エアフロメータ４１の出力値ＡＦＭに存在するノイズの影響を抑制することができる。
【数２２】
2007231840.tif 000024
【００８４】
ところが、式（２２）により流量係数μｔの値を修正する場合であっても、内燃機関が定常運転を行っているときには適切に修正することができるが、内燃機関が過渡運転を行っているときには適切に修正することができない。以下、このことについて簡単に説明する。
【００８５】
ここで、流量係数μｔのマップは、図６では連続的に示されているが、実際には一定の機関運転領域に区切られており、流量係数μｔのマップ値は機関運転領域毎に設定される。具体的には、流量係数μｔのマップは図１９に示したようにスロットル弁開度θｔを複数のマップ領域に分割し、スロットル弁開度θｔの各マップ領域毎に一つの流量係数μｔのマップ値が設定される。したがって、流量係数μｔを修正する場合にはスロットル弁開度θｔの各マップ領域毎にそのマップ領域に対応する流量係数μｔの値を修正することになる。
【００８６】
図１８に示した例を参照すると、スロットル弁開度θｔの一つのマップ領域がθｔ0からθｔ1である場合、スロットル弁開度θｔがそのマップ領域内にある期間は時刻ｔ0から時刻ｔ1までとなっている。したがって、このスロットル弁開度θｔがこのマップ領域にある場合における流量係数μｔの値の誤差率を、時刻ｔ0から時刻ｔ1までのエアフロメータ４１の出力値ＡＦＭの平均値とＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔの平均値との比に基づいて算出すると、比較的正確に流量係数μｔの値の誤差を補償することができると考えられ、斯かる考えに基づいて式（２２）による流量係数μｔの値の修正が提案されている。より詳細には、式（２２）では、エアフロメータ４１の出力値ＡＦＭの平均値とその予想出力値ＡＦＭｍｔの平均値との比は、エアフロメータ４１の出力値ＡＦＭの積分値とその予想出力値ＡＦＭｍｔの積分値との比に等しいため、式（２２）では現在マップに保存されている流量係数μｔ（θｔ）に、エアフロメータ４１の出力値ＡＦＭの積分値（図１８中の面積Ｓ1＋Ｓ2）を予想出力値ＡＦＭｍｔの積分値（図１８中の面積Ｓ1）で除算したものを乗算することによって（θｔ・（Ｓ1＋Ｓ2）／Ｓ1）、流量係数μｔ（θｔ）を修正した値μｔ’（θｔ）を算出している。
【００８７】
ところが、上述したように、図１８に示した例では、時刻ｔ0以降は実際のスロットル通過空気流量とスロットルモデルＭ１０により算出されたスロットル通過空気流量の推定値ｍｔとがほぼ一致しているにも関わらず、エアフロメータ４１に応答遅れが存在するため、エアフロメータ４１の出力値ＡＦＭの時刻ｔ0から時刻ｔ1までの積分値（Ｓ1＋Ｓ2）とＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔの時刻ｔ0から時刻ｔ1までの積分値（Ｓ1）とが一致しない。このため、スロットル弁開度θｔがθｔ1からθｔ2であるマップ領域における流量係数μｔには誤差が無いにも関わらず流量係数μｔが修正されてしまうことになる。
【００８８】
このように、式（２１）および式（２２）により流量係数μｔの値を修正すると、内燃機関が定常運転を行っているときには的確に修正することができるが、内燃機関が過渡運転を行っているときには的確に修正することができない。そこで、本発明では、内燃機関が過渡運転を行っているときにおいても的確に流量係数μｔの値を修正することができるような内燃機関の制御装置を提供する。
【００８９】
図２０は、図１８と同様に、スロットル弁開度θｔと空気流量に関する各種パラメータの値とのタイムチャートである。以下では、図２０を参照して、スロットル弁開度θｔがθｔ0からθｔ1であるマップ領域Ａにおける、流量係数μｔの修正について考える。式（１３）および式（１４）に示したように、ＡＦＭモデルＭ５０では、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔに対して一次遅れが存在するとしてエアフロメータ４１の予想出力値ＡＦＭｍｔを算出している。このため、ＡＦＭモデルＭ５０の予想出力値ＡＦＭｍｔは、スロットルモデルＭ１０によって算出された現在のスロットル通過空気流量の推定値ｍｔだけでなく所定時間前のスロットル通過空気流量の推定値ｍｔに基づいて算出されることになる。
【００９０】
例えば、図２０に示したように時刻ｔ0においてスロットル弁開度θｔがマップ領域Ａに侵入したとき、時刻ｔ0以降においてＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔは、時刻ｔ0以降においてスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔだけでなく、時刻ｔ0以前にスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔを用いて算出されることになる。すなわち、ＡＦＭモデルＭ５０によって算出される予想出力値ＡＦＭｍｔは、スロットル弁開度θｔがマップ領域Ａ以外のマップ領域にあるときのスロットル通過空気流量の推定値ｍｔに基づいて算出されることになる。
【００９１】
ここで、上述したように図２０に示した例では、スロットル弁開度θｔがθｔ0からθｔ1であるマップ領域Ａに対応する流量係数μｔの修正を行おうとしており、流量係数μｔの修正を行うに当たって実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとの比をエアフロメータ４１の出力値ＡＦＭとＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔとの比から算出しようとしている。ところが、スロットル弁開度θｔがマップ領域Ａにある場合においてＡＦＭモデルＭ５０によって算出される予想出力値ＡＦＭｍｔは、マップ領域Ａ以外のマップ領域におけるスロットル通過空気流量の推定値ｍｔの影響を受けており、よってスロットル弁開度θｔがマップ領域Ａにある場合におけるエアフロメータ４１の出力値ＡＦＭとＡＦＭモデルＭ５０によって算出される予想出力値ＡＦＭｍｔとの比は、スロットル弁開度θｔがマップ領域Ａにある場合における実際のスロットル通過空気流量とスロットル通過空気流量の推定値ｍｔとの比を表したものではないものとなっている。
【００９２】
そこで、本実施形態では、スロットル弁開度θｔがマップ領域Ａにある場合における流量係数μｔを修正するときには、ＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔを直接利用せずに、この予想出力値ＡＦＭｍｔからスロットル弁開度θｔがマップ領域Ａに侵入する前（図２０に示した例では、スロットル弁開度θｔがθｔ0よりも小さいとき）のスロットル通過空気流量の推定値ｍｔの影響を除去した値（以下、「修正出力値」と称す）ＡＦＭｍｔｅを用いることとする。
【００９３】
本実施形態では、修正出力値ＡＦＭｍｔｅの算出は基本的にＡＦＭモデルＭ５０における予想出力値ＡＦＭｍｔの算出と同様に行われる。すなわち、下記式（２３）および式（２４）によりボビン部４１ａ３およびサポート部４１ａ４の遅れ放熱量ωbe、ωseが算出され、算出された遅れ放熱量ωbe、ωseに基づいて式（２５）により図１８に示したマップと同様なマップから修正出力値ＡＦＭｍｔｅが算出される。式（２３）および式（２４）は、ＥＣＵ３１で処理するために時刻ｔ、計算間隔Δｔを用いて下記式（２６）および式（２７）のように離散化されて用いられる。
【数２３】
2007231840.tif 000025
【数２４】
2007231840.tif 000026
【数２５】
2007231840.tif 000027
【数２６】
2007231840.tif 000028
【数２７】
2007231840.tif 000029
【００９４】
ただし、ＡＦＭモデルＭ５０における予想出力値ＡＦＭｍｔの算出と異なり、修正出力値ＡＦＭｍｔｅの算出はマップ領域毎に行われる。すなわち、スロットル弁開度θｔが或るマップ領域に侵入したとき、図２０に示した例ではスロットル弁開度θｔがマップ領域Ａに侵入したとき（すなわちスロットル弁開度θｔがθｔ0となったとき）にそのマップ領域における修正出力値ＡＦＭｍｔｅの算出が開始されると共に、スロットル弁開度θｔがそのマップ領域から離脱したとき、図２０に示した例ではスロットル弁開度θｔがマップ領域Ａから離脱したとき（すなわちスロットル弁開度θｔがθｔ1となったとき）にそのマップ領域における修正出力値ＡＦＭｍｔｅの算出が終了せしめられる。
【００９５】
そして、修正出力値ＡＦＭｍｔｅの算出にあたっては、スロットル弁開度θｔがマップ領域Ａに侵入する前にスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔ（以下、「過去のスロットル通過空気流量の推定値ｍｔ」と称す）の影響を除去するために、スロットル弁開度θｔがマップ領域Ａに侵入したときの修正出力値ＡＦＭｍｔｅの値（以下、「初期値」と称す）を、エアフロメータ４１の出力値ＡＦＭとほぼ同様な値になるようにしている。
【００９６】
すなわち、上述したように修正出力値ＡＦＭｍｔｅはＥＣＵ３１において上記離散化された式（２６）および式（２７）に基づいて算出されるが、これら式（２６）および式（２７）には初期値以外に過去の（スロットル弁開度θｔがマップ領域Ａに侵入する前の）スロットル通過空気流量の推定値ｍｔが影響するような項は存在しない。逆に言うと、初期値を過去のスロットル通過空気流量の推定値ｍｔの影響を受けていない値に設定しさえすれば、修正出力値ＡＦＭｍｔｅに対する過去のスロットル通過空気流量の推定値ｍｔの影響を最小限に抑制することができる。
【００９７】
本実施形態では、スロットル弁開度θｔがマップ領域Ａに侵入したときの修正出力値ＡＦＭｍｔｅの初期値をそのときのエアフロメータ４１の出力値ＡＦＭとほぼ同様な値としているため、修正出力値ＡＦＭｍｔｅの初期値は過去のスロットル通過空気流量の推定値ｍｔに基づいて定まる値ではなく、スロットル弁開度θｔがマップ領域Ａに侵入したときのエアフロメータ４１の出力値ＡＦＭに基づいて定まる値であり、よって修正出力値ＡＦＭｍｔｅから過去のスロットル通過空気流量の推定値ｍｔの影響が除去される。
【００９８】
このように、修正出力値ＡＦＭｍｔｅの初期値をエアフロメータ４１の出力値ＡＦＭとほぼ同様な値として式（２３）および式（２４）によってその後の修正出力値ＡＦＭｍｔを算出するためには、時刻ｔ0におけるボビン部４１ａ３の遅れ放熱量およびサポート部４１ａ４の遅れ放熱量を時刻ｔ0におけるボビン部４１ａ３からの実際の放熱量およびサポート部４１ａ４からの実際の放熱量に一致させることが必要である。ところが、時刻ｔ0におけるエアフロメータ４１の出力値ＡＦＭから時刻ｔ0におけるボビン部４１ａ３およびサポート部４１ａ４からの実際の放熱量を算出することは困難である。すなわち、図１８に示したマップ等によってエアフロメータ４１の出力値ＡＦＭからボビン部４１ａ３およびサポート部４１ａ４からの実際の放熱量を合計した値を求めることはできるが、ボビン部４１ａ３およびサポー部４１ａ４それぞれからの実際の放熱量を求めることはできない。
【００９９】
そこで、本実施形態では、エアフロメータ４１のボビン部４１ａ３およびサポート部４１ａ４それぞれからの実際の放熱量を比較的正確に算出することができるオブザーバを用いるようにしている。このオブザーバは、時刻やスロットル弁開度θｔのマップ領域とは無関係に、常にエアフロメータ４１のボビン部４１ａ３およびサポート部４１ａ４それぞれからの実際の放熱量の近似値を算出している。
【０１００】
具体的には、オブザーバは、下記式（２８）および式（２９）によってボビン部４１ａ３からの実際の放熱量の近似値（以下、「近似放熱量」と称す）ωboおよびサポート部４１ａ４からの実際の放熱量の近似値（以下、「近似放熱量」と称す）ωsoをそれぞれ算出している。式（２８）および式（２９）において、ＫIは積分制御等に用いられる積分係数と同様な係数である。また、ＡＦＭｍｔｏはオブザーバによって算出された近似放熱量ωbo、ωsoを合計した値に基づいて下記式（３０）により（すなわち、図１８のマップにより）算出された値である。
【数２８】
2007231840.tif 000030
【数２９】
2007231840.tif 000031
【数３０】
2007231840.tif 000032
【０１０１】
したがって、オブザーバでは、エアフロメータ４１の出力値ＡＦＭとオブザーバによって算出された近似放熱量ωbo、ωsoに基づいて算出されるエアフロメータ４１の出力値（以下、「近似出力値」と称す）ＡＦＭｍｔｏとの偏差（ＡＦＭ−ＡＦＭｍｔｏ）を時間的に累積する。この出力値ＡＦＭと近似出力値ＡＦＭｍｔｏとの偏差の累積値は、実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとの誤差に基づいて変化し、この誤差の絶対値が大きいほど上記偏差の累積値の絶対値が大きくなる。したがって、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔに、出力値ＡＦＭとオブザーバによって算出された近似出力値ＡＦＭｍｔｏとの偏差の累積値に係数ＫIを乗じた値を加算することで、実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとの誤差を或る程度補償することができる。これにより、オブザーバでは、実際のスロットル通過空気流量とほぼ等しいスロットル通過空気流量に基づいて近似放熱量ωbo、ωsoが算出されるため、オブザーバによって算出された近似放熱量ωbo、ωsoはボビン部４１ａ３およびサポート部４１ａ４からの実際の放熱量とほぼ等しい値となっている。
【０１０２】
なお、上記オブザーバを内燃機関の制御装置に実装して、スロットル通過空気流量の推定値ｍｔおよびエアフロメータ４１の出力値ＡＦＭから近似放熱量ωbo、ωsoを算出するためには、式（２８）および式（２９）を離散化する必要がある。時刻ｔ、計算間隔Δｔを用いて式（２８）および式（２９）を離散化すると、それぞれ下記式（３１）および式（３２）が得られる。なお、式（３１）および式（３２）中のΔＡＦＭＩＮＴ（ｔ）は、時刻ｔにおけるＡＦＭ（ｔ）−ＡＦＭｍｔｏ（ｔ）の累積値であり、下記式（３３）により表される。
【数３１】
2007231840.tif 000033
【数３２】
2007231840.tif 000034
【数３３】
2007231840.tif 000035
【０１０３】
このようにして実装されたオブザーバでは、時刻ｔにおけるボビン部４１ａ３およびサポート部４１ａ４からの近似放熱量ωbo（ｔ）およびωso（ｔ）、スロットルモデルＭ１０の式（８）によって算出された時刻ｔにおけるスロットル通過空気流量の推定値ｍｔ（ｔ）、時刻ｔにおけるエアフロメータ４１の出力値ＡＦＭ（ｔ）、および時刻ｔにおいてオブザーバによって算出された近似出力値ＡＦＭｍｔｏ（ｔ）が式（３１）、式（３２）および式（３３）に代入され、これにより時刻ｔ＋Δｔにおける近似放熱量ωbo（ｔ＋Δｔ）およびωso（ｔ＋Δｔ）が算出される。そして、これら近似放熱量ωbo（ｔ＋Δｔ）およびωso（ｔ＋Δｔ）に基づいて式（３０）により時刻ｔ＋Δｔにおける近似出力値ＡＦＭｍｔｏ（ｔ＋Δｔ）が算出され、この時刻ｔ＋Δｔにおける近似出力値ＡＦＭｍｔｏ（ｔ＋Δｔ）はその後の近似放熱量ωboおよびωsoの算出に利用される。そして、このような計算を繰り返すことによって、スロットル通過空気流量の推定値ｍｔ、エアフロメータ４１の出力値ＡＦＭから、任意の時刻ｔにおける近似放熱量ωbo、ωsoが算出される。
【０１０４】
図２１は、図２０の時刻ｔ0〜ｔ1における流量に関するパラメータの推移を拡大して示したタイムチャートである。図２１から分かるように、オブザーバによって算出された近似出力値ＡＦＭｍｔｏ（図２１中の短い破線）は、エアフロメータ４１の出力値ＡＦＭとほぼ同様な値となっている。上述したように修正出力値ＡＦＭｍｔｅの算出にあたっては、スロットル弁開度θｔがマップ領域Ａに侵入したときに、すなわち時刻ｔ0においてオブザーバによって算出されたボビン部４１ａ３の近似放熱量ωboおよびサポート部４１ａ４の近似放熱量ωsoが用いられる。したがって、修正出力値ＡＦＭｍｔｅの算出においては、遅れ放熱量ωbe、ωseの初期値、すなわち時刻ｔ0における遅れ放熱量ωbe（ｔ0）、ωse（ｔ0）はオブザーバによって算出された時刻ｔ0における近似放熱量ωbo（ｔ0）、ωse（ｔ0）とされる（ωbe（ｔ0）＝ωbo（ｔ0）、ωse（ｔ0）＝ωso（ｔ0））。このため、時刻ｔ0における修正出力値ＡＦＭｍｔｅは、時刻ｔ0におけるオブザーバの近似出力値ＡＦＭｍｔｏと同一の値となる。
【０１０５】
そして、修正出力値ＡＦＭｍｔｅの算出に当たっては、遅れ放熱量ωbe、ωseの初期値、すなわち時刻ｔ0における遅れ放熱量ωbe（ｔ0）、ωse（ｔ0）が算出された後に、上述した式（２６）および式（２７）によって時刻ｔ0以降における遅れ放熱量ωbe（ｔ）、ωse（ｔ）が算出され、そしてこれら時刻ｔ0以降の時刻ｔにおける遅れ放熱量ωbe（ｔ）、ωse（ｔ）に基づいて上記式（２５）を用いて時刻ｔ0以降の時刻ｔにおける修正出力値ＡＦＭｍｔｅ（ｔ）が算出される。このようにして算出された修正出力値ＡＦＭｍｔｅは図２１に点線で示したように推移する。すなわち、修正出力値ＡＦＭｍｔｅは時刻ｔ0においてはオブザーバによって算出された近似出力値ＡＦＭｍｔｏと同一の値となっており、その後徐々にＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔに近づくように推移する。
【０１０６】
このようにして算出された修正出力値ＡＦＭｍｔｅは、上述したようにスロットル弁開度θｔがマップ領域Ａに侵入する前のスロットル通過空気流量の推定値ｍｔの影響が除去された値となっている。このため、マップ領域Ａにおけるエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの比は、マップ領域Ａにおける実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量との比を比較的正確に示しており、すなわちマップ領域ＡにおけるスロットルモデルＭ１０の流量係数μｔの誤差率を示している。したがって、上記式（２２）と同様に、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとに基づいて下記式（３４）により誤差率Ｅｒｒが算出される。
【数３４】
2007231840.tif 000036
【０１０７】
すなわち、誤差率Ｅｒｒは、スロットル弁開度θｔが特定のマップ領域（図２０、２１に示した例ではマップ領域Ａ）にある間に亘ってエアフロメータ４１の出力値ＡＦＭを積分した値、すなわち時刻ｔ0からｔ1に亘って出力値ＡＦＭを積分した値（図２１中の面積Ｓ3＋Ｓ4）を、スロットル弁開度θｔが上記特定のマップ領域（図２０、２１に示した例ではマップ領域Ａ）にある間に亘って修正出力値ＡＦＭｍｔｅを積分した値、すなわち修正出力値ＡＦＭｍｔｅを時刻ｔ0からｔ1に亘って積分した値（図２１中の面積Ｓ3）で除算した値（（Ｓ3＋Ｓ4）／Ｓ3）に等しい。このようにして算出された誤差率Ｅｒｒを上記特定のマップ領域（図２０、２１に示した例ではマップ領域Ａ）の流量係数μｔのマップ値に乗算することで、流量係数μｔのマップ値が的確に修正される。
【０１０８】
なお、上記実施形態では、スロットル弁開度θｔがθｔ0〜θｔ1であるマップ領域Ａに対応する流量係数μｔのマップ値を修正する場合について示したが、他のマップ領域においても同様にそのマップ領域に対応する流量係数μｔのマップ値が修正せしめられる。これにより、内燃機関の過渡運転中であっても各マップ領域に対応する流量係数μｔのマップ値が的確に修正せしめられる。
【０１０９】
ところで、上述したように、修正出力値ＡＦＭｍｔｅを用いることでスロットル弁開度θｔがマップ領域Ａに侵入する前のスロットル通過空気流量の推定値ｍｔの影響が除去される。ところが、エアフロメータ４１の出力値ＡＦＭもスロットル弁開度θｔがマップ領域Ａに侵入する前の実際のスロットル通過空気流量の影響を受けており、上述したように修正出力値ＡＦＭｍｔｅを用いたとしてもエアフロメータ４１の出力値ＡＦＭから過去の実際のスロットル通過空気流量の影響を除去することができない。このため、上述したように、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとに基づいて算出されたスロットルモデルＭ１０の流量係数μｔの誤差率は、実際の流量係数μｔの誤差率に対して僅かに異なった値となる。
【０１１０】
そこで、上記実施形態の変更例では、エアフロメータ４１の出力値ＡＦＭから過去の実際のスロットル通過空気流量の影響を除去することとしている。以下、エアフロメータ４１の出力値ＡＦＭから過去の実際のスロットル通過空気流量の影響を除去してスロットルモデルＭ１０の流量係数μｔの誤差率を正確に算出する方法について説明する。
【０１１１】
ところで、一般に、エアフロメータ４１等のセンサの応答は、初期値応答とゼロ応答とに分離することができる。ここで、初期値応答とは、或る時刻ｔ以降の状態量がゼロであると仮定した場合における時刻ｔ以降のセンサの出力値の推移であり、ゼロ応答とは、或る時刻ｔまで状態量がゼロであると仮定した場合における時刻ｔ以降のセンサの出力値の推移である。
【０１１２】
図２２は、エアフロメータ４１の出力値ＡＦＭについての初期値応答曲線とゼロ応答曲線とを示す図である。ここで、Φ1は時刻０を基準とした初期値応答曲線（破線）を示しており、時刻０以降に実際のスロットル通過空気流量がゼロであると仮定した場合におけるエアフロメータ４１の出力値の推移を示している。また、Φ2は時刻０を基準としたゼロ応答曲線（一点鎖線）を示しており、時刻０までの実際のスロットル通過空気流量がゼロであると仮定した場合、すなわち時刻０からエアフロメータ４１での計測が開始された場合におけるエアフロメータ４１の出力値の推移を示している。
【０１１３】
すなわち、初期値応答曲線は時刻０におけるエアフロメータ４１の出力値に応じて、すなわち時刻０よりも前の実際のスロットル通過空気流量に応じて変化するのに対して、ゼロ応答曲線は時刻０以降における実際のスロットル通過空気流量に応じて変化する。そして、時刻０以降の各時刻ｔにおけるエアフロメータ４１の出力値ＡＦＭ（ｔ）は、初期値応答曲線の時刻ｔにおける値Φ1（ｔ）とゼロ応答曲線の時刻ｔにおける値Φ2（ｔ）とを合計した値に等しい（ＡＦＭ（ｔ）＝Φ1（ｔ）＋Φ2（ｔ））。そして、このようなエアフロメータ４１の出力値ＡＦＭと、初期値応答曲線の各時刻における値Φ1と、ゼロ応答曲線の各時刻における値Φ2との関係は、基準となる時刻（図２２に示した例では時刻０）を如何なる時刻にしても成立する。
【０１１４】
上述したようにゼロ応答曲線は所定時刻以降における実際のスロットル通過空気流量に応じて変化するため、特定のマップ領域に侵入した時刻を基準としたエアフロメータ４１の出力値ＡＦＭのゼロ応答曲線からは、その特定のマップ領域に侵入する前の実際のスロットル通過空気流量の影響が除去されている。したがって、斯かるゼロ応答曲線を求め、これと修正出力値ＡＦＭｍｔｅについてのゼロ応答曲線（正確には、修正出力値ＡＦＭｍｔｅについてのゼロ応答曲線に対応する曲線）とを後述するように比較することにより、上記特定のマップ領域に侵入する前の実際のスロットル通過空気流量の影響を除去して、流量係数μの誤差率を的確に算出することができる。
【０１１５】
ただし、時刻ｔ0を基準としたゼロ応答を算出するのは困難であることから、本実施形態では、エアフロメータ４１の出力値ＡＦＭから時刻ｔ0を基準とした初期値応答曲線の各時刻における値Φ1を減算することでゼロ応答を算出することとしている（Φ2（ｔ）＝ＡＦＭ（ｔ）−Φ1（ｔ））。
【０１１６】
ここで、時刻ｔ0を基準とした初期値応答曲線の各時刻における値Φ1の算出方法について説明する。上述したように初期値応答曲線の各時刻における値Φ1は、時刻ｔ0以降に実際のスロットル通過空気流量がゼロであると仮定した場合におけるエアフロメータ４１の出力値の推移を示しているため、上記式（１３）および式（１４）のスロットル通過空気流量の推定値ｍｔをゼロとすることにより求めることができる。すなわち、初期値応答曲線についての放熱量ωbi、ωsiは下記式（３５）および式（３６）により表される。この場合、時刻ｔ0における放熱量ωbi（ｔ0）、ωsi（ｔ0）は、時刻ｔ0における実際の放熱量の近似値である近似放熱量ωbo（ｔ0）、ωso（ｔ0）が用いられる。
【数３５】
2007231840.tif 000037
【数３６】
2007231840.tif 000038
【０１１７】
そして、式（３５）、（３６）により算出された放熱量ωbi、ωsiに基づいて、式（３７）により図１８に示したマップと同様なマップから初期値応答曲線の各時刻における値Φ1が算出される。なお、式（３５）および式（３６）は、ＥＣＵ３１で処理するために離散化されて用いられる。
【数３７】
2007231840.tif 000039
【０１１８】
このようにして算出された時刻ｔ0を基準とした初期値応答曲線の各時刻における値Φ1をエアフロメータ４１の出力値ＡＦＭから減算することにより、エアフロメータ４１の出力値についての時刻ｔ0を基準としたゼロ応答曲線の各時刻における値Φ2が算出される。
【０１１９】
また、修正出力値ＡＦＭｍｔｅも時刻ｔ0において近似放熱量ωbo、ωsoを利用しているため、修正出力値ＡＦＭｍｔｅから上述した初期値応答曲線の各時刻における値Φ1を減算することにより、修正出力値ＡＦＭｍｔｅについてのゼロ応答曲線の各時刻における値が算出される。このため、本実施形態では、スロットルモデルＭ１０の流量係数μｔの誤差率Ｅｒｒを算出するにあたって、エアフロメータ４１の出力値ＡＦＭと、修正出力値ＡＦＭｍｔｅと、エアフロメータ４１の出力値についての時刻ｔ0を基準とした初期値応答曲線の各時刻における値Φ1とを用いて下記式（３８）により算出される。
【数３８】
2007231840.tif 000040
【０１２０】
このように、流量係数μｔの誤差率を算出するにあたって、修正出力値ＡＦＭｍｔｅを用いることにより、ＡＦＭモデルＭ５０によって算出された予想出力値ＡＦＭｍｔに存在する過去のスロットル通過空気流量の推定値ｍｔの影響が除去され、すなわち所定時間前のスロットル通過空気流量の影響で予想出力値ＡＦＭｍｔに含まれている成分を予想出力値ＡＦＭｍｔから除去される。さらに、エアフロメータ４１の出力値ＡＦＭおよび修正出力値ＡＦＭｍｔｅから初期値応答曲線の各時刻における値Φ1を減算することにより、エアフロメータ４１の出力値ＡＦＭ等に存在する過去の実際のスロットル通過空気流量の影響が除去され、すなわち所定時間前のスロットル通過空気流量の影響でエアフロメータの出力値ＡＦＭに含まれている成分をこのエアフロメータの出力値から除去することができる。
【０１２１】
このため、上記式（３８）により算出された誤差率Ｅｒｒはスロットル弁開度θｔが所定のマップ領域内にあるときの実際のスロットル通過空気流量およびスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値に基づいて算出されることになるため、スロットル弁開度θｔが所定のマップ領域内にあるときの流量係数の誤差率を正確に算出することができる。
【０１２２】
次に、ＡＦＭモデルＭ５０の修正操作について説明する。ＡＦＭモデルＭ５０の修正操作について説明するにあたり、まずＡＦＭモデルＭ５０のモデル計算式（１３）、（１４）で用いられるボビン部４１ａ３およびサポート部４１ａ４からの放熱における一次遅れの時定数τBおよびτSを修正する必要性について説明する。
【０１２３】
内燃機関の運転に伴ってエアフロメータ４１のボビン部４１ａ３およびサポート部４１ａ４にはデポジットが付着する。このようにデポジットが付着すると、ボビン部４１ａ３およびサポート部４１ａ４の熱容量が大きくなり、ボビン部４１ａ３およびサポート部４１ａ４からの放熱に遅れが生じる。このため、ボビン部４１ａ３およびサポート部４１a４へのデポジットの付着に伴って、ＡＦＭモデルＭ５０の式（１３）、（１４）で用いられる時定数τBおよびτSに誤差が生じる。このように時定数τBおよびτSに誤差が生じると、上記エアモデルＭ２によって算出される筒内充填空気量Ｍｃに誤差が生じる。
【０１２４】
図２３は、モデルブロックＭ２'、Ｍ２''、Ｍ２'''によって算出される吸気管内圧力のタイムチャートである。図１３を用いて説明したように、エアモデルＭ２では、モデルブロックＭ２'で算出された所定時間後の吸気管内圧力Ｐｍｆ（図中の実線）に、モデルブロックＭ２'''で算出されたエアフロメータの出力に基づく吸気管内圧力の推定値Ｐｍａｆｍ（図中の破線）とモデルブロックＭ２''で算出されたＡＦＭモデルＭ５０を用いた吸気管内圧力の推定値Ｐｍｍｄｌ（図中の一点鎖線）との差分（図中のα）を加算することとしている。
【０１２５】
ここで、例えばＡＦＭモデルＭ５０の時定数τbおよびτsが実際のエアフロメータ４１における一次遅れに対応していなくて大きすぎる場合、モデルブロックＭ２''で算出されたＡＦＭモデルＭ５０を用いた吸気管内圧力の推定値は図２３にＰｍｍｄｌ’で示した二点鎖線のように推移してしまう。すなわち、上述したように差分αは、ＡＦＭモデルＭ５０の時定数τbおよびτsが適切な値になっていればエアモデルによって算出された現在の吸気管内圧力とエアフロメータ４１に基づいて算出された吸気管内圧力との誤差を的確に表しているが、ＡＦＭモデルＭ５０の時定数τbおよびτsが適切な値になっていないと斯かる誤差を的確に表していない。したがって、ＡＦＭモデルＭ５０の時定数τbおよびτsが適切でない場合に所定時間後の吸気管内圧力Ｐｍｆに差分αを加えたものに基づいて筒内充填空気量Ｍｃを算出すると、筒内充填空気量Ｍｃを適切に算出することができない。このため、エアモデルＭ２において筒内充填空気量Ｍｃを適切に算出するためには、ＡＦＭモデルＭ５０において用いられる時定数τbおよびτsを実際のエアフロメータ４１における一次遅れに対応するような適切な値に設定する必要がある。
【０１２６】
次に、ＡＦＭモデルＭ５０のモデル計算式、特に時定数τbおよびτsの修正についての基本概念について説明する。上述したように、スロットルモデルＭ１０の修正操作において用いられた修正出力値ＡＦＭｍｔｅは、ＡＦＭモデルＭ５０を用いて算出された値であり、したがって当然時定数τbおよびτsの影響を受ける。
【０１２７】
図２４は、空気流量に関する各種パラメータの値のタイムチャートである。図２４は、特に、実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されるスロットル通過空気流量の推定値ｍｔとの間に誤差があり、実際のスロットル通過空気流量よりもスロットル通過空気流量の推定値ｍｔの方が小さい場合の各種パラメータの値の推移を示している。
【０１２８】
図２４（ａ）は、実際のスロットル通過空気流量が徐々に増大している場合（すなわち、内燃機関の加速時）、図２４（ｂ）は、実際のスロットル通過空気流量が徐々に減少している場合（すなわち、内燃機関の減速時）をそれぞれ示している。また、図中のＡＦＭはエアフロメータ４１の出力値、ＡＦＭｍｔｅはＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが実際の一次遅れに対応した適切な値となっている場合に算出される修正出力値、ＡＦＭｍｔｅ'はＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも小さい値となっている場合に算出される修正出力値、ＡＦＭｍｔｅ''はＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも大きな値となっている場合に算出された修正出力値をそれぞれ示している。
【０１２９】
図２４（ａ）に示したようにスロットル通過空気流量が徐々に増大している場合、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも小さな値になっていると、修正出力値ＡＦＭｍｔｅが大きくなる傾向にある。この理由は以下の通りである。すなわち、修正出力値ＡＦＭｍｔｅは、基本的に、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔに基づいてＡＦＭモデルＭ５０によって算出される。そして、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも小さいと、エアフロメータ４１の実際の応答遅れよりも小さい応答遅れが生じているものとして修正出力値ＡＦＭｍｔｅが算出されることになる。実際のスロットル通過空気流量が増大している場合には、このようにエアフロメータ４１の応答遅れが実際よりも小さいものとして修正出力値ＡＦＭｍｔｅが算出されると、増大しているスロットル通過空気流量の推定値ｍｔに迅速に応答して修正出力値ＡＦＭｍｔｅは直ぐに上昇し、図２４（ａ）に一点鎖線で示したＡＦＭｍｔｅ'のように推移し、大きい値となる。特に、図２４（ａ）に示した例では、修正出力値がＡＦＭｍｔｅ'のように大きい値となると、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅ'との間の差分ΔＡＦＭは小さいものとなる。
【０１３０】
一方、図２４（ａ）に示したようにスロットル通過空気流量が徐々に増大している場合、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも大きな値となっていると、修正出力値ＡＦＭｍｔｅが小さくなる傾向にある。すなわち、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが大きいと、エアフロメータ４１の実際の応答遅れよりも大きい応答遅れが生じているものとして修正出力値ＡＦＭｍｔｅが算出されることになる。実際のスロットル通過空気流量が増大している場合には、このようにエアフロメータ４１の応答遅れが実際よりも大きいものとして修正出力値ＡＦＭｍｔｅが算出されると、増大しているスロットル通過空気流量の推定値ｍｔにゆっくり応答して修正出力値ＡＦＭｍｔｅは大きく遅れて上昇し、図２４（ａ）に二点鎖線で示したＡＦＭｍｔｅ''のように推移し、小さい値となる。特に、図２４（ａ）に示した例では、修正出力値がＡＦＭｍｔｅ''のように小さい値となると、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅ''との間の差分ΔＡＦＭは大きいものとなる。
【０１３１】
また、図２４（ｂ）に示したように、スロットル通過空気流量が徐々に減少している場合、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも小さな値となっていると、修正出力値ＡＦＭｍｔｅが小さくなる傾向にある。すなわち、上述したように時定数τbおよびτsが小さいとエアフロメータ４１の実際の応答遅れよりも小さな応答遅れが生じているものとして修正出力値ＡＦＭｍｔｅが算出されるため、減少しているスロットル通過空気流量の推定値ｍｔに迅速に応答して修正出力値ＡＦＭｍｔｅは直ぐに下降し、図２４（ｂ）に一点鎖線で示したＡＦＭｍｔｅ'のように推移し、小さい値となる。特に、図２４（ｂ）に示した例では、修正出力値がＡＦＭｍｔｅ’のように小さい値となると、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅ'との間の差分ΔＡＦＭは大きいものとなる。
【０１３２】
一方、図２４（ｂ）に示したように、スロットル通過空気流量が徐々に減少している場合、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも大きな値となっていると、修正出力値ＡＦＭｍｔｅが大きくなる傾向にある。すなわち、上述したように時定数τbおよびτsが大きいとエアフロメータ４１の実際の応答遅れよりも大きな応答遅れが生じているものとして修正出力値ＡＦＭｍｔｅが算出されるため、減少しているスロットル通過空気流量の推定値ｍｔにゆっくりと応答して修正出力値ＡＦＭｍｔｅは大きく遅れて下降し、図２４（ｂ）に二点鎖線で示したＡＦＭｍｔｅ''のように推移し、大きい値となる。特に、図２４（ｂ）に示した例では、修正出力値がＡＦＭｍｔｅのように大きい値となると、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅ''との間の差分ΔＡＦＭは小さいものとなる。
【０１３３】
上記説明をまとめると、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが適切な値よりも小さい場合、スロットル通過空気流量が増大しているときには修正出力値ＡＦＭｍｔｅが大きい値となる（図２４（ａ）のＡＦＭｍｔｅ'参照）のに対して、スロットル通過空気流量が減少しているときには修正出力値ＡＦＭｍｔｅが小さい値となる（図２４（ｂ）のＡＦＭｍｔｅ'参照）。逆に、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが最適な値よりも大きい場合、スロットル通過空気流量が増大しているときには修正出力値ＡＦＭｍｔｅが小さい値となる（図２４（ａ）のＡＦＭｍｔｅ''参照）のに対して、スロットル通過空気流量が減少しているときには修正出力値ＡＦＭｍｔｅが大きい値となる（図２４（ｂ）のＡＦＭｍｔｅ''参照）。
【０１３４】
特に、図２４に示したように、実際のスロットル通過空気流量よりもスロットルモデルＭ１０によって算出されるスロットル通過空気流量の推定値ｍｔの方が小さい場合、時定数τbおよびτsが最適な値よりも小さいと、スロットル通過空気流量が増大しているときにはエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分ΔＡＦＭが小さくなり、逆にスロットル通過空気流量が減少しているときにはこの差分ΔＡＦＭが大きくなる。一方、時定数τbおよびτsが最適な値よりも大きいと、スロットル通過空気流量が増大しているときにはエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分ΔＡＦＭが大きくなり、逆にスロットル通過空気流量が増大しているときにはこの差分ΔＡＦＭが小さくなる。さらに、時定数τbおよびτsが最適な値となっている場合、スロットル通過空気流量が増大しているときであってもスロットル通過空気流量が減少しているときであってもエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分ΔＡＦＭはほぼ同じ値となる。
【０１３５】
そこで、本発明では、ＡＦＭモデルＭ５０で用いられる時定数τbおよびτsが最適な値となるようにすべく、スロットル通過空気流量が増大しているときとスロットル通過空気流量が減少しているときとでエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分ΔＡＦＭが異なる値となる場合には、差分ΔＡＦＭが同じ値となるように時定数τbおよびτsを修正することとしている。
【０１３６】
特に、本実施形態では、上述したように時定数τbおよびτsは、式（１６）および式（１７）によってそれぞれ算出されるが、このうち定数ｋbおよびｋsを修正することとしている。また、この定数ｋbおよびｋsは、一定の機関運転領域毎に定められて、図１９に示したマップと同様なマップとしてＥＣＵ３１のＲＯＭ３４に保存されている。具体的には、定数ｋbおよびｋsのマップは、スロットル弁開度θｔを複数のマップ領域に分割し、スロットル弁開度θｔの各マップ領域毎に一つの定数ｋbおよびｋsのマップ値が設定される。したがって、定数ｋbおよびｋsを修正する場合にはスロットル弁開度θｔの各マップ領域毎にそのマップ領域に対応する定数ｋbおよびｋsのマップ値を修正することになる。
【０１３７】
したがって、本実施形態では、スロットル弁開度θｔが或るマップ領域にあるときのエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとに基づいて、そのマップ領域に対応する定数ｋbおよびｋsのマップ値を修正することとしている。
【０１３８】
具体的には、時定数τbおよびτsの定数ｋbおよびｋsは下記のように修正される。まず、スロットル通過空気流量が増大しているときに、スロットル弁開度θｔが或るマップ領域に侵入した時刻をｔ0、離脱した時刻をｔ1とすると、時刻ｔ0から時刻ｔ1までのエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分ΔＡＦＭの積分値（すなわち、図２５において斜線で示された出力値ＡＦＭの実線、修正出力値ＡＦＭｍｔｅの破線、時刻ｔ0の直線および時刻ｔ1の直線で囲まれた部分の面積）を、時刻ｔ0から時刻ｔ1までの時間で除算したもの（以下、これを「平均誤差」と称す）ΔＡＦＭａｖｅが下記式（３９）によって算出される。
【数３９】
2007231840.tif 000041
【０１３９】
そして、スロットル通過空気流量が減少しているときにも同様にして、スロットル弁開度θｔが上記或るマップ領域に侵入してから離脱するまでのエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分ΔＡＦＭの積分値を、スロットル弁開度θｔが上記或るマップ領域に侵入してから離脱するまでの時間で除算した平均誤差ΔＡＦＭａｖｅが上記式（３９）によって算出される。
【０１４０】
このようにして、スロットル通過空気流量が増大しているときおよび減少しているときそれぞれについて少なくとも一回、或るマップ領域に対応する平均誤差ΔＡＦＭａｖｅが算出される。そして、平均誤差ΔＡＦＭａｖｅが正である場合、すなわちエアフロメータ４１の出力値が修正出力値ＡＦＭｍｔｅよりも大きい場合、スロットル通過空気流量が増大しているときに算出された平均誤差（以下、「増大時の平均誤差ΔＡＦＭａｖｅｉ」と称す）がスロットル通過空気流量が減少しているときに算出された平均誤差（以下、「減少時の平均誤差ΔＡＦＭａｖｅｄ」と称す）よりも小さいときには時定数τbおよびτsが大きくなるように定数ｋbおよびｋsの値を修正し、増大時の平均誤差ΔＡＦＭａｖｅｉが減少時の平均誤差ΔＡＦＭａｖｅｄよりも大きいときには時定数τbおよびτsが小さくなるように定数ｋbおよびｋsの値を修正する。
【０１４１】
一方、平均誤差ΔＡＦＭａｖｅが負である場合、すなわちエアフロメータ４１の出力値が修正出力値ＡＦＭｍｔｅよりも小さい場合、増大時の平均誤差ΔＡＦＭａｖｅｉが減少時の平均誤差ΔＡＦＭａｖｅｄよりも大きいときには時定数τbおよびτsが大きくなるように定数ｋbおよびｋsの値を修正し、増大時の平均誤差ΔＡＦＭａｖｅｉが減少時の平均誤差ΔＡＦＭａｖｅｄよりも小さいときには時定数τbおよびτsが小さくなるように定数ｋbおよびｋsの値を修正する。
【０１４２】
さらに、増大時の平均誤差ΔＡＦＭａｖｅｉと減少時の平均誤差ΔＡＦＭａｖｅｄとがほぼ同一である場合には、ＡＦＭモデルＭ５０の時定数τbおよびτsは適切な値となっていると考えられるため、これら時定数τbおよびτsの修正を行わない。
【０１４３】
また、増大時の平均誤差ΔＡＦＭａｖｅｉと減少時の平均誤差ΔＡＦＭａｖｅｄとの間の差分ΔｉｄはＡＦＭモデルＭ５０で用いられる時定数τbおよびτsのずれの大きさを表しており、この差分Δｉｄが大きいほど適切な値からの時定数τbおよびτsのずれが大きい。
【０１４４】
そこで、本実施形態では、時定数τbおよびτsを修正する際には、すなわち定数ｋbおよびｋsの値を修正する際には、増大時の平均誤差ΔＡＦＭａｖｅｉと減少時の平均誤差ΔＡＦＭａｖｅｄとの間の差分Δｉｄに応じて定数ｋbおよびｋsの修正率Ｍｄｆを定めるようにしており、この差分Δｉｄが大きくなるほど定数ｋbおよびｋsの修正割合が大きくなるようにしている。なお、平均誤差ΔＡＦＭａｖｅの差分Δｉｄと時定数τbおよびτsの修正量との関係はシミュレーションによって計算することができるため、或いは予め実験によって算出することができるため、これに基づいて平均誤差ΔＡＦＭａｖｅの差分Δｉｄと修正率Ｍｄｆとの関係をマップとしてＥＣＵ３１のＲＯＭ３４に保存するようにしている。そして、機関運転時において上述したように算出された平均誤差ΔＡＦＭａｖｅの差分Δｉｄに基づいてこのマップにより修正率Ｍｄｆが算出されると共に、現在の定数ｋbおよびｋsにこの修正率Ｍｄｆを乗算したものがそれぞれ新たな定数ｋbおよびｋsとされる（ｋb＝Ｍｄｆ・ｋb、ｋs＝Ｍｄｆ・ｋs）。
【０１４５】
このようにして、ＡＦＭモデルＭ５０の時定数τbおよびτsの定数ｋbおよびｋsを修正することにより、時定数τbおよびτsを実際のエアフロメータ４１における一次遅れに対応するような適切な値へと修正することができる。
【０１４６】
なお、スロットル通過空気流量が増大しているときまたは減少しているときであるか否かは、例えば、スロットル弁開度θｔが或るマップ領域に侵入してから離脱するまでの間におけるエアフロメータ４１の出力値ＡＦＭまたはスロットル通過空気流量ｍｔの増減によって判断される。また、スロットル弁開度θｔが或るマップ領域に侵入してから離脱するまでの時間が所定時間よりも長い場合には単位時間当たりのスロットル通過空気流量の増大または減少の程度が小さいため、平均誤差ΔＡＦＭａｖｅの算出を行わないようにしてもよい。
【０１４７】
図２６は、増大時の平均誤差ΔＡＦＭａｖｅｉまたは減少時の平均誤差ΔＡＦＭａｖｅｄの算出・保存制御の制御ルーチンを示すフローチャートである。図２６に示した制御ルーチンは一定時間間隔の割り込みによって行われる。
【０１４８】
まず、ステップ１０１では、スロットル弁開度θｔが所定のマップ領域にあるか否かが判定される。スロットル弁開度θｔが上記所定のマップ領域に侵入したときには、ステップ１０１においてスロットル弁開度θｔが所定のマップ領域にあると判定されてステップ１０２およびステップ１０３へと進む。ステップ１０２およびステップ１０３では、エアフロメータ４１の出力値ＡＦＭが検出されると共に修正出力値ＡＦＭｍｔｅが算出される。次いで、ステップ１０４では、侵入フラグＸが１とされる。なお、この侵入フラグＸは、スロットル弁開度θｔが所定のマップ領域にある場合に１とされ、それ以外の場合に０とされるフラグである。そしてスロットル弁開度θｔが上記所定のマップ領域内にある間はこれらステップ１０１〜１０４が繰り返される。
【０１４９】
その後、スロットル弁開度θｔが上記所定のマップ領域から離脱すると、ステップ１０１からステップ１０５へと進む。ステップ１０５では、侵入フラグＸが１であるか否かが判定される。スロットル弁開度θｔが上記所定のマップ領域から離脱した直後には侵入フラグＸが１のままとなっており、この場合にはステップ１０６へと進む。ステップ１０６では、ステップ１０２およびステップ１０３において検出および算出されたエアフロメータの出力値ＡＦＭおよび修正出力値ＡＦＭｍｔｅに基づいて式（３９）により平均誤差ΔＡＦＭａｖｅが算出される。
【０１５０】
次いで、ステップ１０７では、ステップ１０６において算出された平均誤差ΔＡＦＭａｖｅが所定値Ｐ以上であるか否かが判定される。所定値Ｐは、比較的小さい値であり、それよりも平均誤差ΔＡＦＭａｖｅが小さいとエアフロメータ４１のノイズの影響が大きくなってしまうような値である。ステップ１０７において、平均誤差ΔＡＦＭａｖｅが所定値Ｐよりも小さいと判定された場合には制御ルーチンが終了せしめられる。一方、平均誤差ΔＡＦＭａｖｅが所定値Ｐ以上であると判定された場合にはステップ１０８へと進む。ステップ１０８では、スロットル弁開度θｔが上記所定のマップ領域内にある間にエアフロメータ４１の出力値ＡＦＭが増大したかまたは減少したかが判定される。
【０１５１】
次いで、ステップ１０９では、ステップ１０８においてスロットル弁開度θｔが上記所定のマップ領域内にある間にエアフロメータ４１の出力値ＡＦＭが増大したと判定された場合には、ステップ１０６で算出された平均誤差ΔＡＦＭａｖｅを増大時の平均誤差ΔＡＦＭａｖｅｉとしてＥＣＵ３１のＲＯＭ３４に保存し、一方、ステップ１０８においてスロットル弁開度θｔが上記所定のマップ領域内にある間にエアフロメータ４１の出力値ＡＦＭが減少したと判定された場合には、ステップ１０６で算出された平均誤差ΔＡＦＭａｖｅを減少時の平均誤差ΔＡＦＭａｖｅｉとしてＥＣＵ３１のＲＯＭ３４に保存する。次いで、ステップ１１０において侵入フラグが０とされて制御ルーチンが終了せしめられる。そして、その次のルーチンでは、ステップ１０５において侵入フラグが１でないと判定されて、制御ルーチンが終了せしめられる。
【０１５２】
図２７は、図２６に示した制御ルーチンによって算出・保存された増大時の平均誤差ΔＡＦＭａｖｅｉおよび減少時の平均誤差ΔＡＦＭａｖｅｄに基づいて時定数τbおよびτsの定数ｋbおよびｋsを修正する制御ルーチンを示すフローチャートである。図２７に示した制御ルーチンは、前回本制御ルーチンを実行してからΔＡＦＭａｖｅｉおよびΔＡＦＭａｖｅｄが少なくともそれぞれ一つずつ図２６に示した制御ルーチンによって算出された後に開始される。
【０１５３】
まず、ステップ１２１では、図２６に示したΔＡＦＭａｖｅｄの算出・保存制御のステップ１０９で算出された増大時の平均誤差ΔＡＦＭａｖｅｉおよび減少時の平均誤差ΔＡＦＭａｖｅｄが取得される。次いで、ステップ１２２では、ステップ１２１で取得された増大時の平均誤差ΔＡＦＭａｖｅｉ（複数取得された場合にはその平均値）からステップ１２１で取得された減少時の平均誤差ΔＡＦＭａｖｅｄを減算して差分Δｉｄが算出される（Δｉｄ＝ΔＡＦＭａｖｅｉ−ΔＡＦＭａｖｅｄ）。次いで、ステップ１２３では、ステップ１２２で算出されたΔｉｄに基づいてマップにより修正率Ｍｄｆが算出される。次いで、ステップ１２４では、ＡＦＭモデルＭ５０の時定数τbおよびτsの現在の定数ｋbおよびｋsにステップ１２３で算出された修正率Ｍｄｆを乗算したものが新たな定数ｋbおよびｋsとされる（ｋb＝Ｍｄｆ・ｋb、ｋs＝Ｍｄｆ・ｋs）。
【０１５４】
なお、上記実施形態では、平均誤差ΔＡＦＭａｖｅをエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分に基づいた値としているが、エアフロメータ４１の出力値ＡＦＭと予想出力値ＡＦＭｍｔとの差分に基づいた値としてもよい。
【０１５５】
また、上記実施形態では、スロットル通過空気流量が増大しているときとスロットル通過空気流量が減少しているときとのエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分ΔＡＦＭに基づいて時定数τbおよびτsを修正することとしているが、これに限らず、例えば内燃機関の加速時と内燃機関の減速時とのエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの差分等、スロットル通過空気流量が変化するようにスロットル通過空気流量以外の運転パラメータが変化しているときのこれら出力値の差分に基づいて時定数を修正するようにしてもよい。
【図面の簡単な説明】
【０１５６】
【図１】本発明の制御装置を搭載した内燃機関全体を示す図である。
【図２】図１に示したエアフロメータの概略斜視図である。
【図３】図２に示したエアフロメータの熱線計量部の拡大斜視図である。
【図４】エアフロメータの出力電圧とエアフロ通過空気流量との関係を示す図である。
【図５】エアモデルのブロック図である。
【図６】スロットル弁開度と流量係数との関係を示す図である。
【図７】スロットル弁開度と開口断面積との関係を示す図である。
【図８】関数Φ（Ｐｍ／Ｐａ）を示す図である。
【図９】スロットルモデルの基本概念を示す図である。
【図１０】吸気管モデルの基本概念を示す図である。
【図１１】吸気弁モデルの基本概念を示す図である。
【図１２】筒内充填ガス量および筒内吸気ガス量の定義に関する図である。
【図１３】図５に示したエアモデルとは別のエアモデルのブロック図である。
【図１４】アクセルペダル踏込み量と目標スロットル弁開度との関係を示す図である。
【図１５】スロットル通過空気流量と完全放熱量との関係を示す図である。
【図１６】遅れ放熱量の和とエアフロメータの予想出力値との関係を示す図である。
【図１７】スロットルモデルの修正操作を示す図である。
【図１８】スロットル弁開度と空気流量に関する各種パラメータの値とのタイムチャートである。
【図１９】スロットル弁開度に基づいて定められた流量係数のマップを示すである。
【図２０】スロットル弁開度と空気流量に関する各種パラメータの値とのタイムチャートである。
【図２１】図２０の一部を拡大した図である。
【図２２】初期値応答とゼロ応答とを説明するための図である。
【図２３】吸気管内圧力のタイムチャートである。
【図２４】エアフロメータの出力値と修正出力値とのタイムチャートである。
【図２５】空気流量に関する各種パラメータの値とのタイムチャートである。
【図２６】増大時の平均誤差または減少時の平均誤差の算出・保存制御の制御ルーチンを示すフローチャートである。
【図２７】時定数τbおよびτsの定数ｋbおよびｋsを修正する制御ルーチンを示すフローチャートである。
【符号の説明】
【０１５７】
１    機関本体
５    燃焼室
６    吸気弁
７    吸気ポート
８    排気弁
１１    燃料噴射弁
１４    吸気枝管
１５    サージタンク
１６    吸気管
１９    スロットル弁
２３    吸気管部分
【図１】
2007231840.tif 000042
【図２】
2007231840.tif 000043
【図３】
2007231840.tif 000044
【図４】
2007231840.tif 000045
【図５】
2007231840.tif 000046
【図６】
2007231840.tif 000047
【図７】
2007231840.tif 000048
【図８】
2007231840.tif 000049
【図９】
2007231840.tif 000050
【図１０】
2007231840.tif 000051
【図１１】
2007231840.tif 000052
【図１２】
2007231840.tif 000053
【図１３】
2007231840.tif 000054
【図１４】
2007231840.tif 000055
【図１５】
2007231840.tif 000056
【図１６】
2007231840.tif 000057
【図１７】
2007231840.tif 000058
【図１８】
2007231840.tif 000059
【図１９】
2007231840.tif 000060
【図２０】
2007231840.tif 000061
【図２１】
2007231840.tif 000062
【図２２】
2007231840.tif 000063
【図２３】
2007231840.tif 000064
【図２４】
2007231840.tif 000065
【図２５】
2007231840.tif 000066
【図２６】
2007231840.tif 000067
【図２７】
2007231840.tif 000068