JP2007239484A

JP2007239484A - 内燃機関の吸入空気量推定装置

Info

Publication number: JP2007239484A
Application number: JP2006059328A
Authority: JP
Inventors: Satoshi Tanaka; 聡田中
Original assignee: Toyota Motor Corp
Current assignee: Toyota Motor Corp
Priority date: 2006-03-06
Filing date: 2006-03-06
Publication date: 2007-09-20

Abstract

【課題】エアモデルを構成する吸気管モデルに生じた誤差を常時正確に修正することができる吸入空気量推定装置を提供する。
【解決手段】吸入空気量推定装置は、スロットル通過空気流量に基づいて吸気管モデル計算式を用いて吸気管内圧力を算出する吸気管モデルＭ２０を具備し、吸気管モデルによって算出された吸気管内圧力の推定値Ｐｍに基づいて筒内充填空気量を推定する。また、この装置は、実際の吸気管内圧力Ｐｍａｃｔを検出する吸気管内圧力検出手段をさらに具備し、吸気管モデルの吸気管モデル計算式においてパラメータとして吸気管容積Ｖｍが用いられ、吸気管モデルによって吸気管内圧力の推定値を算出すると共に吸気管内圧力検出手段によって吸気管内圧力の実測値を検出し、これら吸気管内圧力の推定値と実測値とに基づいて吸気管モデルで用いられる吸気管容積の値を修正するようにした。
【選択図】図１４

Description

本発明は、内燃機関の吸入空気量推定装置に関する。

内燃機関の吸気系統を、スロットル弁、吸気管、吸気弁等の要素に分けてそれぞれの要素毎にモデル化して数式で表すことにより機関の吸入空気量（筒内充填空気量）を計算によって求める、いわゆるエアモデルを用いた内燃機関の吸入空気量推定装置が知られている。

このようなエアモデルを用いた吸入空気量推定装置では、例えば、後述するように大気圧、大気温度以外には機関回転数とスロットル弁開度とのみによって筒内充填空気量を算出することが可能となる。また、モデル計算式を用いた計算により筒内充填空気量が算出できるため、スロットル弁開度等の変化速度が大きい過渡運転時にも応答性良好に吸入空気量を算出することが可能となっている。

このようなエアモデルを用いた吸入空気量推定装置では、吸気系をモデル化した計算式に基づいて筒内充填吸気量を算出するため、吸気系のモデル化の際に誤差があると正確な筒内充填空気量を算出することができない問題がある。

また、仮に初期状態では実際の吸気系を正確にモデル化していた場合であっても、例えば、スロットル弁等は長期間の間に徐々に作動特性や流量特性が変化するような場合がある。また、初期状態においてもスロットル弁等の特性は公差内でばらつきを生じており、必ずしも正確にモデルと一致していない場合がある。したがって、モデル計算式によって正確な筒内充填空気量を算出するためには、例えば機関運転中に或る運転パラメータの実測値とモデル計算式によって算出されたこの運転パラメータの推定値とを比較する等してモデル計算式の誤差を修正することが必要となる。

そこで、特許文献１に記載の吸入空気量推定装置では、エアフロメータ及び吸気管内圧力センサによって検出された空気流量及び吸気管内圧力に基づいてエアモデルを構成する各要素モデル毎にモデル計算式を修正することとしている。特に、上記吸入空気量推定装置では、スロットル弁モデル及び吸気管モデルのモデル計算式をそれぞれ個別に修正することとしており、スロットル弁モデルのモデル計算式をエアフロメータによって検出された実測値に基づいて修正し、吸気管モデルのモデル計算式を吸気管内圧力センサによって検出された実測値に基づいて修正することとしている。このように、各要素モデル毎にモデル計算式を実測値に合わせて修正することにより、エアモデル全体が実際の吸気系に合わせて的確に修正されるようになり、吸入空気流量の推定精度が向上するとされている。

特開２００４−２１１５９０号公報特表２００４−５２８５０４号公報特開平１１−９４０５７号公報

ところで、特許文献１に記載の吸入空気量推定装置には、エアモデルの一つの要素モデルとして吸気管モデルが含まれている。吸気管モデルはスロットル弁モデルによって算出されたスロットル通過空気流量と吸気弁モデルによって算出された筒内吸入空気流量とに基づいて吸気管内圧力を算出するモデルであり、この吸気管モデルにも誤差が生じる場合がある。ところが、特許文献１に記載の装置では、吸気管モデルのモデル計算式の修正については全く考慮していない。

そこで、本発明の目的は、エアモデルを構成する吸気管モデルに生じた誤差を常時正確に修正することにより吸入空気量の計算精度を向上させることが可能な内燃機関の吸入空気量推定装置を提供することにある。

上記課題を解決するために、第１の発明では、スロットル通過空気流量に基づいて吸気管モデル計算式を用いて吸気管内圧力を算出する吸気管モデルを具備し、該吸気管モデルによって算出された吸気管内圧力の推定値に基づいて筒内充填空気量を推定する内燃機関の吸入空気量推定装置において、実際の吸気管内圧力を検出する吸気管内圧力検出手段をさらに具備し、上記吸気管モデルの吸気管モデル計算式においてパラメータとして吸気管容積が用いられ、上記吸気管モデルによって吸気管内圧力の推定値を算出すると共に上記吸気管内圧力検出手段によって吸気管内圧力の実測値を検出し、これら吸気管内圧力の推定値と実測値とに基づいて上記吸気管モデルで用いられる吸気管容積の値を修正するようにした。
吸気管モデルは吸気管容積（吸気管の一部の容積。下記実施形態では、スロットル弁から吸気弁までの間の吸気管の容積）のみが問題となるモデルであり、吸気管モデルで用いられる吸気管容積の値が実際の吸気管容積に一致していれば吸気管モデルには誤差がほとんど生じない。第１の発明によれば、吸気管内圧力の推定値と実測値とに基づいて誤差の生じる可能性がある吸気管容積の値を修正することで、吸気管モデルを常時適切に修正することができる。

第２の発明では、第１の発明において、上記吸気管モデルで用いられる吸気管容積の値の修正は、吸気管内圧力の変化率が基準変化率以上であるときの吸気管内圧力の推定値と実測値とに基づいて行われる。
吸気管モデルで用いられる吸気管容積の値に誤差が生じている場合、吸気管内圧力の変化率が高いときに吸気管内圧力の推定値と実測値との間に誤差が生じ易い。第２の発明では、吸気管内圧力の変化率が基準変化率以上であるときの吸気管内圧力の推定値と実測値との誤差に基づいて吸気管容積の値の誤差が修正されるため、吸気管容積の値の誤差を的確に修正することができる。

第３の発明では、第２の発明において、上記吸気管モデルで用いられる吸気管容積の値は、吸気管内圧力の変化率が基準変化率以上であるときの吸気管内圧力の推定値と実測値とが等しくなるように修正される。

第４の発明では、第１〜第３のいずれか一つの発明において、上記吸気管モデルは、スロットル通過空気流量に加えて筒内吸入空気流量に基づいて吸気管内圧力を算出し、上記スロットル通過空気流量は、スロットル弁開度に基づいてスロットルモデル計算式を用いてスロットル通過空気流量を算出するスロットルモデルによって算出され、上記筒内吸入空気流量は、吸気管内圧力に基づいて吸気弁モデル計算式を用いて筒内吸入空気流量を算出する吸気弁モデルによって算出され、上記吸気管内圧力検出手段によって検出された吸気管内圧力に基づいて上記スロットルモデル計算式又は吸気弁モデル計算式を修正するスロットル・吸気弁モデル修正手段をさらに具備し、上記吸気管モデルの吸気管容積の値の修正は、スロットル・吸気弁モデル修正手段によってスロットルモデル計算式又は吸気弁モデル計算式を修正した後に行われる。
吸気管モデルではスロットルモデル及び吸気弁モデルの出力値を利用しているため、これらスロットルモデル及び吸気弁モデルに誤差が生じていると、吸気管モデルで用いられる吸気管容積の値に誤差が生じていなくても吸気管内圧力の推定値と実測値との間で誤差が生じてしまう。第４の発明によれば、予めスロットルモデル及び吸気弁モデルの修正が行われるため、吸気管容積の値に誤差が生じている場合に限り吸気管内圧力の推定値と実測値との間に誤差が生じるようになり、適切に吸気管容積の値の誤差を修正することができるようになる。

本発明によれば、吸気管モデルに生じた誤差を常時正確に修正することにより吸入空気量の計算精度を向上させることが可能な内燃機関の吸入空気量推定装置が提供される。

以下、図面を参照して本発明の実施形態について説明する。図１に概略的に示した機関本体１は筒内噴射型火花点火式内燃機関を示す。しかしながら、本発明をポート噴射型火花点火式内燃機関等の別の火花点火式内燃機関や圧縮自着火式内燃機関に適用してもよい。

図１に示したように、本発明の実施形態では機関本体１はシリンダブロック２と、シリンダブロック２内で往復動するピストン３と、シリンダブロック２上に固定されたシリンダヘッド４とを具備する。ピストン３とシリンダヘッド４との間には燃焼室５が形成される。シリンダヘッド４には各気筒毎に吸気弁６と、吸気ポート７と、排気弁８と、排気ポート９とが配置される。さらに、図１に示したようにシリンダヘッド４の内壁面の中央部には点火プラグ１０が配置され、シリンダヘッド４内壁面周辺部には燃料噴射弁１１が配置される。また、ピストン３の頂面には燃料噴射弁１１の下方から点火プラグ１０の下方まで延びるキャビティ１２が形成されている。さらに、シリンダヘッド４には、吸気弁６の位相角およびバルブリフト量を連続的に変更可能な吸気弁制御装置１３が設けられている。

各気筒の吸気ポート７は吸気枝管１４を介してサージタンク１５に連結され、サージタンク１５は吸気管１６を介してエアクリーナ１７に連結される。吸気管１６内にはステップモータ１８によって駆動されるスロットル弁１９が配置される。一方、各気筒の排気ポート９は排気管２０に連結され、この排気管２０は排気浄化触媒２１を内蔵したケーシング２２に連結される。なお、以下の説明では、スロットル弁１９から吸気弁６までの吸気枝管１４、サージタンク１５、吸気管１６等の部分を吸気管部分２３と称す。

電子制御ユニット（ＥＣＵ）３１はディジタルコンピュータからなり、双方向性バス３２を介して相互に接続されたＲＡＭ（ランダムアクセスメモリ）３３、ＲＯＭ（リードオンリメモリ）３４、ＣＰＵ（マイクロプロセッサ）３５、入力ポート３６および出力ポート３７を具備する。サージタンク１３には、吸気管部分２３内の空気（吸気ガス）の圧力を検出するための吸気管内圧力センサ４０が設けられており、吸気管内圧力センサ４０は吸気管部分２３内の圧力に比例した出力電圧を発生し、この出力電圧が対応するＡＤ変換器３８を介して入力ポート３６に入力される。

また、スロットル弁１９の上流側の吸気管１６には吸気管１６内を通過して流れる吸入空気流量を検出するエアフロメータ４１が設けられており、本実施形態におけるエアフロメータ４１の構成については後述する。さらに、エアクリーナ１７近傍には吸気温度を検出する吸気温度センサ４２と、大気圧を検出する大気圧センサ４３とが設けられる。スロットル弁１９にはスロットル弁１９の開度を検出するスロットル弁開度センサ４４が設けられており、スロットル弁開度センサ４４はスロットル弁開度に対応する出力信号を発生させる。これらエアフロメータ４１、吸気温度センサ４２、大気圧センサ４３およびスロットル弁開度センサ４４はそれぞれ吸入空気流量（質量流量）、吸気温度（大気温度）、大気圧およびスロットル弁開度に対応する出力信号を発生し、この出力信号が対応するＡＤ変換器３８を介して入力ポート３６に入力される。

また、アクセルペダル４５にはアクセルペダル４５の踏込み量に比例した出力電圧を発生する負荷センサ４６が接続され、負荷センサ４６の出力電圧は対応するＡＤ変換器３８を介して入力ポート３６に入力される。また、クランク角センサ４７は例えばクランクシャフトが１５度回転する毎に出力パルスを発生し、この出力パルスが入力ポート３６に入力される。ＣＰＵ３５ではこのクランク角センサ４５の出力パルスから機関回転数が計算される。

一方、出力ポート３７は対応する駆動回路３９を介して点火プラグ１０、燃料噴射弁１１、吸気弁制御装置１３およびステップモータ１８に接続される。

図２は本実施形態で用いられるエアフロメータ４１の概略斜視図を、図３は図２に示したエアフロメータ４１の熱線計量部の拡大斜視図である。本実施形態では、エアフロメータ４１は、熱線式流量計とされており、図２に示したように吸気管１６内を流れる空気の一部をバイパスさせるバイパス通路と、このバイパス通路にバイパスされた吸入空気の質量流量を計測する熱線計量部４１ａと、計測された質量流量に応じた電圧を出力する信号処理部４１ｂとを有する。熱線計量部４１ａは図３に示したように白金熱線から成る吸気温度計測用抵抗４１ａ１と、この吸気温度計測用抵抗４１ａ１を信号処理部４１ｂに連結して保持するサポート部４１ａ２と、加熱用抵抗（ボビン部）４１ａ３と、この加熱用抵抗４１ａ３を信号処理部４１ｂに連結して保持するサポート部４１ａ４とを備える。信号処理部４１ｂは、吸気温度計測用抵抗４１ａ１と加熱用抵抗４１ａ３とで構成されるブリッジ回路を有し、このブリッジ回路により吸気温度計測用抵抗４１ａ１と加熱用抵抗４１ａ３との温度差を常に一定に維持するように加熱用抵抗４１ａ３に供給する電力を調整すると共に、この供給する電力を電圧に変換して出力するようになっている。エアフロメータ４１の出力電圧Ｖｇとエアフロメータ４１が配置された吸気管１６内を通過する空気の流量（以下、「エアフロ通過空気流量」と称す）との関係は図４に示した通りである。

ところで、内燃機関の制御装置では、内燃機関の燃焼室５において燃焼される混合気の空燃比を目標空燃比にするために、吸気弁６が閉じたときに燃焼室５内に充填されている空気の量（以下、「筒内充填空気量Ｍｃ」と称す）を推定し、推定された筒内充填空気量Ｍｃに基づいて混合気の空燃比が目標空燃比となるように燃料噴射弁１１によって内燃機関の燃焼室５（または吸気通路）に噴射する燃料の量（以下、「燃料噴射量」と称す）を定めている。したがって、内燃機関の燃焼室５において燃焼される混合気の空燃比を正確に目標空燃比とするためには、筒内充填空気量Ｍｃを正確に推定する必要がある。

通常、筒内充填空気量Ｍｃは、エアフロメータ等の多数のセンサ、およびこれらセンサからの出力値を引数とした多数のマップから推定される。ところが、このようにマップを用いて筒内充填空気量Ｍｃを推定する場合、推定される筒内充填空気量Ｍｃの値をより正確なものにするためには、必要なマップの数およびその引数の数が多くなる。このようにマップの数が多くなると、マップを保存するためのＥＣＵのＲＯＭを記憶容量の大きいものにしなければならず、内燃機関の制御装置の製造コストが高くなってしまう。さらに、各マップを作成するにはマップが用いられる内燃機関の形式毎に適合作業を行わなければならないが、この適合作業における測定点はマップの数およびその引数の数に応じて増大するため、マップの数およびその引数の数が多くなると適合作業の工数も増大してしまう。

そこで、マップを用いずに様々なモデルを用いて、数値計算により筒内充填空気量Ｍｃを算出する内燃機関の制御装置が検討されている。このような制御装置では、数値計算を多用することにより必要なマップの数を極力減らすようにしており、これにより適合作業を行う際の工数を大幅に削減しながらも、筒内充填空気量Ｍｃを正確に算出することができる。このような制御装置のうち本願出願人により提案されたものの一つに、図５に示したエアモデルを搭載した制御装置がある。図示したエアモデルは内燃機関に適用される最も単純なモデルであり、以下このエアモデルについて説明する。

エアモデルは、図５に示したようにスロットルモデルＭ１０、吸気管モデルＭ２０、吸気弁モデルＭ３０を備える。スロットルモデルＭ１０には、スロットル弁開度センサ４４によって検出されたスロットル弁１９の開度（スロットル弁開度）θｔと、大気圧センサ４３によって検出された内燃機関周囲の大気圧（または、吸気管１６に吸入される空気の圧力）Ｐａと、吸気温度センサ４２によって検出された内燃機関周囲の大気温度（または、吸気管１６に吸入される空気の温度）Ｔａと、後述する吸気管モデルＭ２０において算出された吸気管部分２３内の圧力（以下、「吸気管内圧力」と称す）Ｐｍとが入力され、これら入力された各パラメータの値を後述するスロットルモデルＭ１０のモデル式に代入することで、単位時間当たりにスロットル弁１９を通過する空気の流量（以下、「スロットル通過空気流量の推定値ｍｔ」と称す）が算出される。スロットルモデルＭ１０において算出されたスロットル通過空気流量の推定値ｍｔは、吸気管モデルＭ２０へ入力される。

吸気管モデルＭ２０には、スロットルモデルＭ１０において算出されたスロットル通過空気流量の推定値ｍｔと、以下で詳述する単位時間当たりに燃焼室５内に流入する空気の流量（以下、「筒内吸入空気流量ｍｃ」と称す。なお、筒内吸入空気流量ｍｃの定義については、吸気弁モデルＭ３０の説明において詳述する）とが入力され、これら入力された各パラメータの値を後述する吸気管モデルＭ２０のモデル式に代入することで、吸気管部分２３内に存在する空気の圧力（以下、「吸気管内圧力Ｐｍ」と称す）と吸気管部分２３内に存在する空気の温度（以下、「吸気管内温度Ｔｍ」と称す）とが算出される。吸気管モデルＭ２０において算出された吸気管内圧力Ｐｍと吸気管内温度Ｔｍは共に吸気弁モデルＭ３０へ入力され、さらに吸気管内圧力ＰｍはスロットルモデルＭ１０にも入力される。

吸気弁モデルＭ３０には、吸気管モデルＭ２０において算出された吸気管内圧力Ｐｍおよび吸気管内温度Ｔｍの他に大気温度Ｔａが入力され、これら入力された各パラメータの値を後述する吸気弁モデルＭ３０のモデル式に代入することで、筒内吸入空気流量ｍｃが算出される。算出された筒内吸入空気流量ｍｃは、筒内充填空気量Ｍｃに変換され、この筒内充填空気量Ｍｃに基づいて燃料噴射弁からの燃料噴射量が決定される。また、吸気弁モデルＭ３０において算出された筒内吸入空気流量ｍｃは吸気管モデルＭ２０に入力される。

図５から分かるように、エアモデルでは或るモデルにおいて算出されたパラメータの値が別のモデルへの入力値として利用されるので、エアモデル全体では、実際に入力される値はスロットル弁開度θｔ、大気圧Ｐａ、および大気温度Ｔａの三つのパラメータのみであり、これら三つのパラメータから筒内充填空気量Ｍｃが算出される。

次に、エアモデルの各モデルＭ１０〜Ｍ３０について説明する。
スロットルモデルＭ１０では、大気圧Ｐａ、大気温度Ｔａ、吸気内圧力Ｐｍ、スロットル弁開度θｔから、下記式（１）に基づいてスロットル通過空気流量の推定値ｍｔが算出される。ここで、式（１）におけるμｔはスロットル弁における流量係数で、スロットル弁開度θｔの関数であり、よって図６に示したようなマップから定まる。また、Ａｔはスロットル弁の開口断面積を示し、スロットル弁開度θｔの関数であり、図７に示したようなマップから定まる。なお、これら流量係数μｔおよび開口断面積Ａｔをまとめたμｔ・Ａｔをスロットル弁開度θｔから一つのマップで求まるようにしてもよい。また、Ｒは気体定数に関する定数であり、実際には気体定数を１mol当たりの気体（空気）の質量Ｍlmolで除算した値である。

また、Φ（Ｐｍ／Ｐａ）は下記式（２）に示した関数であり、この式（２）におけるκは比熱比（一定値とする）である。この関数Φ（Ｐｍ／Ｐａ）は図８に示したようなグラフに表すことができるので、このようなグラフをマップとしてＥＣＵ３１のＲＯＭ３４に保存し、実際には式（２）を用いて計算するのではなくマップからΦ（Ｐｍ／Ｐａ）の値を求めるようにしてもよい。

これらスロットルモデルＭ１０の式（１）および式（２）は、スロットル弁１９上流の気体の圧力を大気圧Ｐａ、スロットル弁１９上流の気体の温度を大気温度Ｔａ、スロットル弁１９の下流の気体の圧力を吸気管内圧力Ｐｍとして、図９に示したようなスロットル弁１９のモデルに対して、質量保存則、エネルギ保存則および運動量保存則を適用し、さらに気体の状態方程式、比熱比の定義式、およびマイヤーの関係式を利用することによって得られる。

吸気管モデルＭ２０では、スロットル通過空気流量の推定値ｍｔ、筒内吸入空気流量ｍｃ、および大気温度Ｔａから、下記式（３）および式（４）に基づいて吸気管内圧力Ｐｍおよび吸気管内温度Ｔｍが算出される。なお、式（３）および式（４）におけるＶｍはスロットル弁１９から吸気弁６までの吸気枝管１４、サージタンク１５、吸気管１６等の部分（吸気管部分２３）の容積に等しい定数である。

ここで、吸気管モデルＭ２０について図１０を参照して説明する。吸気管部分２３の総気体量（総空気量）をＭとすると、総気体量Ｍの時間的変化は、吸気管部分に流入する気体の流量、すなわちスロットル通過空気流量の推定値ｍｔと、吸気管部分から流出する気体の流量、すなわち筒内吸入空気流量ｍｃとの差に等しいため、質量保存則により下記式（５）が得られ、この式（５）および吸気管部分２３における気体の状態方程式（Ｐｍ・Ｖｍ＝Ｍ・Ｒ・Ｔｍ）より、式（３）が得られる。

また、吸気管部分２３の気体のエネルギＭ・Ｃｖ・Ｔｍの時間的変化量は、吸気管部分２３に流入する気体のエネルギと吸気管部分２３から流出する気体のエネルギとの差に等しい。このため、吸気管部分２３に流入する気体の温度を大気温度Ｔａ、吸気管部分２３から流出する気体の温度を吸気管内温度Ｔｍとすると、エネルギ保存則により下記式（６）が得られ、この式（６）および上記気体の状態方程式より、式（４）が得られる。

吸気弁モデルＭ３０では、吸気管内圧力Ｐｍ、吸気管内温度Ｔｍ、および大気温度Ｔａから、下記式（７）に基づいて、筒内吸入空気流量ｍｃが算出される。なお、式（７）におけるａ、ｂは、機関回転数Ｎｅから、さらに吸気弁６の位相角（バルブタイミング）および作用角を変更できる可変動弁機構を備えた内燃機関の場合には吸気弁６の位相角、作用角から定まる値である。

上述した吸気弁モデルＭ３０について図１１を参照して説明する。一般に、吸気弁６が閉じたときに燃焼室５内に充填されている空気の量である筒内充填空気量Ｍｃは、吸気弁６が閉弁するとき（吸気弁閉弁時）に確定し、吸気弁６閉弁時の燃焼室５内の圧力に比例する。また、吸気弁６閉弁時の燃焼室５内の圧力は吸気弁６上流の気体の圧力、すなわち吸気管内圧力Ｐｍと等しいとみなすことができる。したがって、筒内充填空気量Ｍｃは、吸気管内圧力Ｐｍに比例すると近似することができる。

ここで、一定時間（例えば、クランク角７２０°分）当たりに吸気管部分２３から流出する全空気流量を平均化したもの、または一定時間（例えば、クランク角７２０°分）当たりに吸気管部分２３から全ての気筒の燃焼室５に吸入される空気量を上記一定時間で除算したものを筒内吸入空気流量ｍｃ（以下で詳述する）とすると、筒内充填空気量Ｍｃが吸気管内圧力Ｐｍに比例することから、筒内吸入空気流量ｍｃも吸気管内圧力Ｐｍに比例すると考えられる。このことから、理論および経験則に基づいて、上記式（７）が得られる。なお、式（７）におけるａは比例係数であり、ｂは燃焼室５内に残存していた既燃ガスを表す適合値（排気弁８閉弁時に燃焼室５内に残る既燃ガス量を後述する時間ΔＴ180°で除算したものと考えられる）である。また、実際の運転では過渡運転時に吸気管内温度Ｔｍが大きく変化する場合があるため、これに対する補正として理論および経験則に基づいて導かれたＴａ／Ｔｍが乗算されている。

ここで、筒内吸入空気流量ｍｃについて、図１２を参照して内燃機関が４気筒である場合について説明する。なお、図１２は横軸がクランクシャフトの回転角度、縦軸が単位時間当たりに吸気管部分２３から燃焼室５に実際に流入する空気の流量である。図１２に示したように、４気筒の内燃機関では、吸気弁６が例えば１番気筒、３番気筒、４番気筒、２番気筒の順に開弁し、各気筒に対応する吸気弁６の開弁量に応じて吸気管部分２３から各気筒の燃焼室５内へ空気が流入する。例えば、吸気管部分２３から各気筒の燃焼室５内に流入する空気の流量の変位は図１２に破線で示した通りであり、これを総合して吸気管部分２３から全気筒の燃焼室５に流入する空気の流量は図１２に実線で示した通りである。また、例えば１番気筒への筒内充填空気量Ｍｃは図１２に斜線で示した通りである。

これに対して、実線で示した吸気管部分２３から全ての気筒の燃焼室５に流入する空気の流量を平均化したものが筒内吸入空気流量ｍｃであり、図中に一点鎖線で示す。そして、この一点鎖線で示した筒内吸入空気流量ｍｃに、４気筒の場合にはクランクシャフトが１８０°（すなわち、４ストローク式内燃機関において１サイクル中にクランクシャフトが回転する角度７２０°を気筒数で割った角度）回転するのにかかる時間ΔＴ180°を乗算したものが筒内充填空気量Ｍｃとなる。したがって、吸気弁モデルＭ３０で算出された筒内吸入空気流量ｍｃにΔＴ180°を乗算することで、筒内充填空気量Ｍｃが算出される（Ｍｃ＝ｍｃ・ΔＴ180°）。より詳細には、筒内充填空気量Ｍｃが吸気弁閉弁時の圧力に比例することを考慮して、吸気弁閉弁時の筒内吸入空気流量ｍｃにΔＴ180°を乗算したものが筒内充填空気量Ｍｃとされる。

次に、上記エアモデルを内燃機関の制御装置に実装して、実際に筒内充填空気量Ｍｃを算出する場合について説明する。筒内充填空気量Ｍｃはエアモデルを用いて、上記式（１）、式（３）、式（４）、および式（７）を解くことにより表される。この場合、ＥＣＵ３１で処理するために、これらの式を離散化する必要がある。時刻ｔ、計算間隔Δｔを用いて式（１）、式（３）、式（４）、および式（７）を離散化すると、それぞれ下記式（８）、式（９）、式（１０）、および式（１１）が得られる。なお、吸気管内温度Ｔｍ（ｔ＋Δｔ）は、式（９）および式（１０）によってそれぞれ算出されたＰｍ／Ｔｍ（ｔ＋Δｔ）およびＰｍ（ｔ＋Δｔ）から、式（１２）によって算出される。

このようにして実装されたエアモデルでは、スロットルモデルＭ１０の式（８）で算出された時刻ｔにおけるスロットル通過空気流量の推定値ｍｔ（ｔ）と、吸気弁モデルＭ３０の式（１１）で算出された時刻ｔにおける筒内吸入空気流量ｍｃ（ｔ）とが、吸気管モデルＭ２０の式（９）および式（１０）に代入され、これにより時刻ｔ＋Δｔにおける吸気管内圧力Ｐｍ（ｔ＋Δｔ）および吸気管内温度Ｔｍ（ｔ＋Δｔ）が算出される。次いで、算出されたＰｍ（ｔ＋Δｔ）およびＴｍ（ｔ＋Δｔ）は、スロットルモデルＭ１０および吸気弁モデルＭ３０の式（８）および式（１１）に代入され、これにより時刻ｔ＋Δｔにおけるスロットル通過空気流量の推定値ｍｔ（ｔ＋Δｔ）および筒内吸入空気流量ｍｃ（ｔ＋Δｔ）が算出される。そして、このような計算を繰り返すことによって、スロットル弁開度θｔ、大気圧Ｐａ、および大気温度Ｔａから、任意の時刻ｔにおける筒内吸入空気流量ｍｃが算出され、算出された筒内吸入空気流量ｍｃに上記時間ΔＴ180°を乗算することで、任意の時刻ｔにおける筒内充填空気量Ｍｃが算出される。

なお、内燃機関の始動時には、すなわち時刻ｔ＝０においては、吸気管内圧力Ｐｍは大気圧と等しい（Ｐｍ（０）＝Ｐａ）とされ、吸気管内温度Ｔｍは大気温度と等しい（Ｔｍ（０）＝Ｔａ）とされて、各モデルＭ１０〜Ｍ３０における計算が開始される。

なお、上記エアモデルでは、大気温度Ｔａおよび大気圧Ｐａが一定であるとしているが、時刻によって変化する値としてもよく、例えば、大気温度を検出するための吸気温度センサによって時刻ｔにおいて検出された値を大気温度Ｔａ（ｔ）、大気圧を検出するための大気圧センサによって時刻ｔにおいて検出された値を大気圧Ｐａ（ｔ）として上記式（８）、式（１０）、および式（１１）に代入するようにしてもよい。

ところで、上記したように図５の装置ではスロットル弁開度θｔ、大気温度Ｔａおよび大気圧Ｐａのみを入力するモデル式を用いて機関始動時からの逐次計算により筒内充填吸気量が算出される。このため、エアモデルの各要素モデルにモデル化誤差が生じていると、算出される筒内充填吸気量は不正確になり燃料噴射量の正確な算出ができなくなる。

一方、エアモデルの各要素モデル（スロットルモデルＭ１０、吸気管モデルＭ２０および吸気弁モデルＭ３０）を考えると、モデル化そのものには誤差がない場合であっても、使用による特性の変化や製作公差によるばらつきなどによりモデル化誤差が生じている場合がある。このため、エアモデルにより求めた筒内充填吸気量等の値をエアフロメータや吸気管内圧力から求めた吸気量と比較し、モデル式から計算した吸気量をエアフロメータ等により求めた実測値に合致するように修正することが行われている。

しかし、このような従来の修正方法は各要素モデルの集合であるエアモデルの全体を一律に修正するものであるのに対して、実際には各モデルに発生する誤差は一律ではなく、それぞれのモデル毎に相違している。このため、全体を一律に修正する方法では正確なモデルの修正を行うことはできず、例えば内燃機関が定常運転を行っているときに筒内充填吸気量を実測値に合致するように補正すると、内燃機関が過渡運転を行っているときなどでは吸気圧力の計算値が実測値から大きく離れてしまい、過渡運転を行っているときの筒内充填吸気量に大きな誤差が生じるような場合がある。

そこで、本実施形態では、実測値に基づいてモデル修正するものの、エアモデル全体（Ｍ１０〜Ｍ３０）を一律に修正するのではなく、それぞれの要素モデルに応じた補正を別個に行う。

以下、本実施形態のモデル修正方法について説明する。本実施形態では、現在の実際のエアフロ通過空気流量、スロットル弁開度、吸気管内圧力を、それぞれエアフロメータ４１、スロットル弁開度センサ４４、吸気圧センサ４０で実測し、これらの値を用いて図５の筒内充填吸気量計算とは独立してスロットルモデルＭ１０、吸気管モデルＭ２０、吸気弁モデルＭ３０等についての修正操作を行う。

まず、吸気管モデルＭ２０の修正操作について説明する。吸気管モデルＭ２０の式（３）及び式（４）から分かるように、吸気管モデルＭ２０では、スロットル通過空気流量ｍｔ及び筒内吸入空気流量ｍｃ等を正確に算出することができれば、吸気管部分２３の容積（吸気管容積）Ｖｍのみが問題となる。逆に言うと、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分２３の容積に正確に一致していれば、吸気管モデルＭ２０にはモデル化誤差がほとんど生じない。

ここで、吸気管部分２３の容積は、基本的に内燃機関の設計時において算出することができるため、通常、吸気管モデルＭ２０においてはこのようにして算出された吸気管部分２３の容積がパラメータとして用いられる。ところが、実際の吸気管部分２３の容積は、製作公差によるばらつきが生じることや、吸気管部分２３内にデポジットが付着することにより、設計時に算出された容積とは異なる値となってしまう場合がある。吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍを内燃機関の設計時において算出される吸気管部分２３の容積とすると、このような場合には吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍは実際の吸気管部分２３の容積とは異なった値となってしまう。

そこで、本実施形態では、吸気管モデルＭ２０を修正するにあたり、吸気管モデルＭ２０の式（３）及び（４）で用いられる吸気管部分２３の容積Ｖｍに生じた誤差を修正することとしている。

吸気管モデルＭ２０に生じた誤差の修正方法について図１３を参照して説明する。図１３は、吸気管内圧力Ｐｍのタイムチャートであり、特にスロットル弁開度が時刻ｔ₀においてステップ的に大きくなった場合又は小さくなった場合における吸気管内圧力の時間推移を示したものである。図中の実線は吸気管部分２３の容積が大きい場合、破線は吸気管部分２３の容積Ｖｍが小さい場合の吸気管内圧力の時間推移を示している。

図１３（ａ）から分かるように、時刻ｔ₀以前、すなわちスロットル弁開度がステップ的に大きくなる前であって内燃機関が定常運転を行っているときには、吸気管部分２３の容積の大小以外の条件が同一であれば、吸気管部分２３の容積の大小が異なっていたとしても、吸気管内圧力は同一の一定の値となる。すなわち、吸気管部分２３の容積が大きい場合（すなわち、図１３（ａ）に実線で示した場合）と吸気管部分２３の容積が小さい場合（すなわち、図１３（ａ）に破線で示した場合）との間で吸気管部分２３に流入する空気の流量も吸気管部分２３から流出する空気の流量も等しいため吸気管部分２３の容積の大小に関わらず吸気管内圧力は同一の値となり、且つ吸気管部分２３に流入する空気の流量と吸気管部分２３から流出する空気の流量とは等しいため吸気管内圧力は一定の値となる。

その後、時刻ｔ₀においてスロットル弁開度が急激に大きくなると、それに伴って吸気管部分２３の容積に関わらず吸気管内圧力が上昇する。このとき、吸気管部分２３の容積が大きい場合と吸気管部分２３の容積が小さい場合とでは吸気管内圧力の上昇速度が異なり、吸気管部分２３の容積が小さい場合（図１３（ａ）中の破線）の方が吸気管内圧力の上昇速度が速い。なぜなら、スロットル弁開度が大きくなってから、吸気管部分２３内に流入する空気の流量は吸気管部分２３の容積が大きい場合と小さい場合とで同じように増大すると共に、同量の空気が吸気管部分２３に充填されると吸気管部分２３の容積が小さい場合の方が吸気管内圧力が高くなるためである。この結果、吸気管内圧力の上昇中においては、吸気管部分２３の容積が小さい場合の方が大きい場合よりも吸気管内圧力が大きな値となる。

一方、図１３（ｂ）から分かるように、時刻ｔ₀においてスロットル弁開度が急激に小さくなると、それに伴って吸気管部分２３の容積に関わらず吸気管内圧力が降下する。このとき、吸気管部分２３の容積が大きい場合と吸気管部分の容積が小さい場合とでは吸気管内圧力の降下速度が異なり、吸気管部分２３の容積が小さい場合の方が吸気管内圧力の降下速度が速い。なぜなら、スロットル弁開度が小さくなってから吸気管部分２３から各気筒内へ流出する空気の流量は吸気管部分２３の容積が大きい場合と小さい場合とでほぼ等しいものとなっていると共に、同量の空気が吸気管部分２３から流出すると吸気管部分２３に残っている空気の量は吸気管部分２３の容積が小さい場合の方が少ないためである。この結果、吸気管内圧力の降下中においては、吸気管部分２３の容積が小さい場合の方が大きい場合よりも吸気管内圧力が小さな値となる。

以上の観点より、吸気管モデルＭ２０の式（３）及び（４）で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分２３の容積よりも小さい値となっている場合、スロットル弁開度が急激に大きくなったときに吸気管モデルＭ２０によって算出される吸気管内圧力Ｐｍは実際の吸気管内圧力よりも大きく、逆にスロットル弁開度が急激に小さくなったときに吸気管モデルＭ２０によって算出される吸気管内圧力Ｐｍは実際の吸気管内圧力よりも小さい。一方、吸気管モデルＭ２０の式（３）及び（４）で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分２３の容積よりも大きい値となっている場合、スロットル弁開度が急激に大きくなったときに吸気管モデルＭ２０によって算出される吸気管内圧力Ｐｍは実際の吸気管内圧力よりも小さく、逆にスロットル弁開度が急激に小さくなったときに吸気管モデルＭ２０によって算出される吸気管内圧力Ｐｍは実際の吸気管内圧力よりも大きい。そして、吸気管モデルＭ２０の式（３）及び（４）で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分２３の容積にほぼ一致している場合、スロットル弁開度が急激に大きくなったときでも急激に小さくなったときでも吸気管モデルＭ２０によって算出される吸気管内圧力Ｐｍは実際の吸気管内圧力にほぼ一致している。

そこで、本発明では、吸気管内圧力センサ４０により吸気管内圧力の実測値（実際の吸気管内圧力に対応する値）Ｐｍａｃｔを検出すると共に、吸気管モデルＭ２０によって吸気管内圧力の推定値Ｐｍを算出し、これら吸気管内圧力の実測値Ｐｍａｃｔと推定値Ｐｍとに基づいて、吸気管モデルＭ２０の式（３）及び（４）で用いられる吸気管部分２３の容積Ｖｍの値を修正することとしている。特に、本実施形態では、スロットル弁開度が急激に大きくなったとき又は急激に小さくなったときに吸気管モデルＭ２０によって算出される吸気管内圧力Ｐｍが実際の吸気管内圧力とほぼ一致するように吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍの値が修正される。

具体的には、本実施形態では、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍの値の修正は下記のようにして行われる。以下の説明では、スロットル弁開度が急激に大きくなった場合における修正方法を示しているが、スロットル弁開度が急激に小さくなった場合においても同様にして修正が行われる。

図１４は、吸気管内圧力及びスロットル弁開度のタイムチャートである。図示した例では、時刻ｔ₀までスロットル弁開度θｔがアイドリング状態での開度（すなわち、アクセルペダルの踏込み量がゼロであるときに対応する開度）θｔｉになっていると共に、時刻ｔ₀においてスロットル弁開度が急激に大きくされ、その後、大きい開度θｔｌで維持される。吸気管内圧力センサ４０によって検出される吸気管内圧力の実測値Ｐｍａｃｔ及び吸気管モデルＭ２０によって算出された吸気管内圧力の推定値Ｐｍは、時刻ｔ₀においてスロットル弁開度が急激に大きくされると、それに伴って上昇し始める。その後、或る程度時間が経過すると、大きいスロットル弁開度θｔｌに対応して吸気管内圧力の実測値Ｐｍａｃｔ及び推定値Ｐｍ共に一定の値に収束する。

図１４から分かるように、単位時間当たりの吸気管内圧力の実測値Ｐｍａｃｔの変化量（以下、「変化率」と称す）ＲＣは、時刻ｔ₀において急激に大きくなると共に時間の経過に伴って徐々に小さくなっている。ここで、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣについて後述するような基準変化率ＲＣｒｅｆを設定すると、図１４に示した例では、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣは時刻ｔ₀から時刻ｔ₁まで基準変化率ＲＣｒｅｆ以上となる。すなわち、時刻ｔ₀からｔ₁までに亘って吸気管内圧力が大きく上昇している。なお、上記基準変化率ＲＣｒｅｆは、スロットル弁開度が急激に変化した場合に限り吸気管内圧力の実測値の変化率が取り得る値であって、スロットル弁開度が変化していない場合又はスロットル弁開度が小さくしか変化していない場合には吸気管内圧力の実測値の変化率がほとんど取り得ない程度に大きい値である。

本実施形態では、このように吸気管内圧力の変化率が大きい期間において、吸気管内圧力センサ４０によって吸気管内圧力の実測値Ｐｍａｃｔを検出すると共に、吸気管モデルＭ２０によって吸気管内圧力の推定値Ｐｍを算出する。そして、下記式（１３）により、吸気管内圧力の変化率が大きい期間における吸気管内圧力の実測値Ｐｍａｃｔと推定値Ｐｍとの平均比率Ｐｒａｖｅを算出する。

上述したように、吸気管モデルＭ２０の式（３）及び（４）で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分２３の容積にほぼ一致している場合には、吸気管内圧力の実測値Ｐｍａｃｔと推定値Ｐｍとはほぼ一致するため、上記平均比率Ｐｒａｖｅはほぼ１となる。逆に、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分の容積に一致していない場合には上記平均比率Ｐｒａｖｅは１より大きい値又は１より小さい値となると共に、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍと実際の吸気管部分２３の容積とのずれが大きいほど上記平均比率Ｐｒａｖｅの値が１から離れた値となる。そこで、本実施形態では、上記平均比率Ｐｒａｖｅの１からのずれ量が大きいときほど吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍの修正率ηを大きいものとすることとしている。

また、スロットル弁開度が急激に大きくなったときに、上述したように算出された平均比率Ｐｒａｖｅが１よりも小さい場合（図１４に示した例のようになっている場合）、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分２３の容積よりも小さい値となっている。従ってこの場合、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍの値を大きくするように修正する。一方、スロットル弁開度が急激に大きくなったときに、上述したように算出された平均比率Ｐｒａｖｅが１よりも大きい場合、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍが実際の吸気管部分２３の容積よりも大きい値となっていることを表している。従ってこの場合、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍの値を小さくするように修正する。

特に、本実施形態では、吸気管部分２３の容積Ｖｍの誤差率ηと平均比率Ｐｒａｖｅとの関係を予め実験により又は計算により求め、図１５（ａ）に示したようなマップとしてＥＣＵ３１のＲＯＭ３４に保存しておく。そして、機関運転時において上述したように平均比率Ｐｒａｖｅを算出すると共に、算出された平均比率Ｐｒａｖｅと図１５（ａ）に示したマップとに基づいて吸気管部分２３の容積Ｖｍの誤差率ηを算出すると共に、算出された誤差率ηを現在吸気管モデルＭ２０で用いられている吸気管部分２３の容積Ｖｍに乗算することとしている（Ｖｍ＝η・Ｖｍ）。なお、図１５（ａ）に示したマップはスロットル弁開度が急激に大きくなった場合において用いられるマップを示しており、スロットル弁開度が急激に小さくなった場合には例えば図１５（ｂ）に示したようなマップが用いられる。

また、このように吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍの値を修正するに当たっては、吸気管モデルＭ２０に入力されるパラメータの値、すなわちスロットル通過空気流量の推定値ｍｔ及び筒内吸入空気流量ｍｃが実際のスロットル通過空気流量及び筒内吸入空気流量を正確に表している必要がある。すなわち、吸気管モデルＭ２０に入力されるパラメータの値が実際の値に対して誤差を生じていると、例え吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍが正確であっても、吸気管モデルＭ２０によって算出される吸気管内圧力の推定値Ｐｍが吸気管内圧力センサ４０によって検出された吸気管内圧力と異なる値となってしまう。その結果、吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍが不必要に修正されてしまうことになる。

そこで、本実施形態では、上述したように吸気管モデルＭ２０で用いられる吸気管部分２３の容積Ｖｍの値を修正する前に、後述するような方法でスロットルモデルＭ１０及び吸気弁モデルＭ３０のモデル計算式の修正を行うこととしている。このように、スロットルモデルＭ１０及び吸気弁モデルＭ３０のモデル計算式の修正を予め行っておくことで、スロットルモデルＭ１０によって算出されるスロットル弁通過空気流量ｍｔが実際のスロットル通過空気流量を正確に表すようになり、また吸気弁モデルＭ３０によって算出される筒内吸入空気流量ｍｃが実際の筒内吸入空気流量を正確に表すようになる。このため、修正後のスロットルモデルＭ１０及び吸気弁モデルＭ３０によって算出されたスロットル通過空気流量ｍｔ及び筒内吸入空気流量ｍｃに基づいて吸気管モデルＭ２０を用いて算出された吸気管内圧力の推定値Ｐｍを用いて吸気管モデルＭ２０に用いられる吸気管部分２３の容積Ｖｍの値を修正することで、吸気管部分２３の容積Ｖｍが不必要に修正されてしまうことを防止することができる。

図１６は、吸気管モデルＭ２０の式（３）及び（４）に用いられる吸気管部分２３の容積Ｖｍの修正制御のフローチャートである。図示した制御ルーチンは、一定時間間隔の割り込みによって行われる。

図１６を参照すると、ステップ１０１では、変化率フラグｆｌａｇが０であるか否かが判定される。変化率フラグｆｌａｇは、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆ以上である場合、すなわちスロットル弁開度が急激に大きくなった場合に１とされ、それ以外の場合に０とされるフラグである。ステップ１０１において、変化率フラグｆｌａｇが０であるであると判定された場合（ｆｌａｇ＝０）には、ステップ１０２へと進む。

ステップ１０２では、スロットルモデルＭ１０及び吸気弁モデルＭ３０の修正が完了しているか否かが判定される。具体的には、スロットルモデルＭ１０の修正に関しては後述するスロットルモデルＭ１０の修正操作における誤差率Ｅｒｒが予め定められた所定値（比較的小さい値）以下である場合、吸気弁モデルＭ３０の修正に関しては後述する誤差率αが予め定められた所定値（比較的小さい値）以下である場合にこれらモデルＭ１０及びＭ３０の修正が完了していると判定され、それ以外の場合には修正が完了していないと判定される。ステップ１０２において、スロットルモデルＭ１０及び吸気弁モデルＭ３０の修正が完了していないと判定された場合には、制御ルーチンが終了せしめられる。従って、これらモデルＭ１０及びＭ３０の修正が完了するまでは吸気管部分２３の容積Ｖｍの修正が行われないことになる。

一方、ステップ１０２において、スロットルモデルＭ１０及び吸気弁モデルＭ３０の修正が完了していると判定された場合には、ステップ１０３へと進む。ステップ１０３では、吸気管内圧力センサ４０によって検出された吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆ以上であるか否かが判定される。ステップ１０３において、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆよりも小さいと判定された場合には、制御ルーチンが終了せしめられる。従って、変化率ＲＣが基準変化率ＲＣｒｅｆよりも小さい場合には基準変化率ＲＣｒｅｆ以上になるまで待つことになる。

次いで、ステップ１０４において時間カウンタｃｎｔ及び積算値ｓｕｍが０にリセットされる。ここで、時間カウンタｃｎｔは吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆ以上になってからの時間を表しており、積算値ｓｕｍは吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆ以上になってからの制御ルーチンの実行回数分だけＰｍａｃｔ／Ｐｍの値を積算した値を表している。次いで、ステップ１０５において変化率フラグが１にセットされる。

このように変化率フラグが１にセットされると、次の制御ルーチンではステップ１０１において変化率フラグが０でないと判定されて、ステップ１０６へと進む。ステップ１０６では、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆ以上であるか否かが判定され、変化率ＲＣが基準変化率ＲＣｒｅｆ以上のままになっている場合にはステップ１０７へと進む。ステップ１０７では、吸気管内圧力センサ４０によって検出された現在の吸気管内圧力の実測値Ｐｍａｃｔを吸気管モデルＭ２０によって算出された現在の吸気管内圧力の推定値Ｐｍで除算した値を積算値ｓｕｍに加算したものが新たな積算値ｓｕｍとされる（ｓｕｍ＝ｓｕｍ＋（Ｐｍａｃｔ／Ｐｍ））。次いで、ステップ１０８では、現在の時間カウンタｃｎｔの値に１を加算したものが新たな時間カウンタｃｎｔの値とされる（ｃｎｔ＝ｃｎｔ＋１）。そして、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆ以上となっている間、ステップ１０６及びステップ１０７が繰り返し実行される。

その後、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆよりも小さくなった場合には、ステップ１０６からステップ１０９へと進む。ステップ１０９では、ステップ１０７を繰り返すことによって算出された積算値ｓｕｍを時間カウンタｃｎｔで除算した値が平均比率Ｐｒａｖｅとされる（Ｐｒａｖｅ＝ｓｕｍ／ｃｎｔ）。次いで、ステップ１１０では、ステップ１０９で算出された平均比率Ｐｒａｖｅに基づいて図１５（ａ）に示したマップを用いて吸気管部分２３の容積Ｖｍの誤差率ηが算出される。ステップ１１１では、吸気管モデルＭ２０の式（３）及び（４）で現在用いられている吸気管部分２３の容積Ｖｍの値に誤差率ηを乗算したものが新たな吸気管部分２３の容積Ｖｍとされる（Ｖｍ＝η・Ｖｍ）。そして、ステップ１１２では、変化率フラグｆｌａｇが０にリセットされ、制御ルーチンが終了せしめられる。

なお、上記実施形態では、スロットル弁開度が急激に大きくなった場合として、吸気管内圧力の実測値Ｐｍａｃｔの変化率ＲＣが基準変化率ＲＣｒｅｆ以上である場合を用い、斯かる場合の吸気管内圧力の実測値Ｐｍａｃｔ及び推定値Ｐｍに基づいて吸気管モデルＭ２０を修正している。しかしながら、スロットル弁開度が急激に大きくなった場合としては、斯かる場合に限らず、例えば吸気管モデルＭ２０によって算出された吸気管内圧力の推定値Ｐｍの変化率が基準変化率以上である場合や、スロットル弁開度の変化量又はアクセルペダルの踏込み量の変化量が基準変化量以上である場合を用いても良い。

次に、スロットルモデルＭ１０の修正操作について説明する。図１７は、スロットルモデルＭ１０の修正操作を説明するブロック図である。図１７に示したように、スロットルモデルの修正操作を行う際には、エアフロメータモデルＭ４０が用いられる。エアフロメータモデルＭ４０には、スロットルモデルＭ１０によって算出されたスロットル通過空気流量（以下、「スロットル通過空気流量の推定値」と称す）ｍｔが入力され、この入力されたパラメータの値を後述するエアフロメータモデルＭ４０のモデル式に代入することで、実際にスロットル弁１９を通過する空気の流量（以下、「実際のスロットル通過空気流量」と称す）がこの値のときにエアフロメータ６１が出力するであろう出力値が算出される。

すなわち、実際のエアフロメータ４１の出力は、実際のエアフロ通過空気流量（エアフロ通過空気流量とスロットル通過空気流量はほぼ同一であると考えられるため、以下ではスロットル通過空気流量として説明する）に対して固有の応答特性に基づく応答遅れを有しており、したがってそのときの実際のスロットル通過空気流量のみに対応した値となっているわけではない。エアフロメータモデルＭ４０は、この応答特性をシュミレートしたモデルであり、上記エアフロメータ４１の応答遅れを考慮した上で、スロットル弁１９をスロットル通過空気流量の推定値ｍｔだけ空気が流れている場合にエアフロメータ４１が出力するであろう出力値（以下、「予想出力値」と称す）を算出する。

ここで、エアフロメータモデルＭ４０について具体的に説明する。エアフロメータ４１では、上述したように吸気温度計測用抵抗４１ａ１と加熱用抵抗（ボビン部）４１ａ３との温度差を常に一定に維持するように加熱用抵抗４１ａ３に供給する電力を調整するとともに、その供給電力に基づいてエアフロ通過空気流量を算出するようにしている。ここで、この供給電力はボビン部４１ａ３およびサポート部４１ａ４から吸気管１６内を通過する空気への放熱量を示すことから、エアフロメータ４１はボビン部４１ａ３およびサポート部４１ａ４からの放熱量に基づいてエアフロ通過空気流量を算出すると言い換えることができる。

ここで、ボビン部４１ａ３は、より詳細には円筒状のセラミックスボビンに白金熱線を巻回し、その外周にガラスをコーティングすることにより形成される。このため、白金熱線から周囲の空気への放熱は白金熱線と周囲の空気との間にガラス層が介在することにより遅れてしまう。したがって、吸気管１６内を通過する空気の流量が急激に増大したような場合であっても、ボビン部４１ａ３からの単位時間当たりの放熱量（以下、単に「放熱量」と称す）は直ぐには増大せず、或る程度の遅れをもって増大することになる。換言すると、ボビン部４１ａ３からの放熱量は実際のエアフロ通過空気流量に対して応答遅れが存在する。また、同様なことがサポート部４１ａ４からの放熱量についても言え、サポート部４１ａ４からの放熱量は実際のエアフロ通過空気流量に対して応答遅れが存在する。

この応答遅れは、一次遅れに近似することができることがわかっており、ボビン部４１ａ３の放熱量の応答遅れは次の式（１４）で表される。ここで、式（１４）におけるω_bはエアフロメータ４１の出力値から換算したボビン部４１ａ３の放熱量、すなわち放熱遅れが生じた結果ボビン部４１ａ３の白金熱線から実際に放熱される放熱量（以下、「遅れ放熱量」と称す）を示している。また、式（１４）におけるＷ_bは、応答遅れを補償した放熱量、すなわち放熱遅れが生じないと仮定した場合におけるボビン部４１ａ３の白金熱線から放熱される放熱量（以下、「完全放熱量」と称す）を示している。すなわち、完全放熱量Ｗ_bは、内燃機関が定常運転を行っているときにおける放熱量に等しく、基本的にエアフロメータ４１近傍の吸気管１６内を通過する空気流量のみの関数である。エアフロメータ４１近傍の吸気管１６内を通過する空気流量はスロットル通過空気流とほぼ等しいため、完全放熱量Ｗ_bはスロットル通過空気流量の関数と考えることができ、スロットル通過空気流量とボビン部４１ａ３からの完全放熱量Ｗ_bとの関係は図１８のように表すことができる。さらに、式（１４）におけるτ_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（ｔ＋Δｔ）を用いて式（１６）により時刻ｔ＋Δｔにおけるエアフロメータ４１の予想出力値ＡＦＭｍｔ（ｔ＋Δｔ）が算出される。

再びスロットルモデルの修正操作について話を戻すと、本実施形態では、吸気圧センサ４０とスロットル弁開度センサ４４とで実測した現在の吸気管内圧力Ｐｍとスロットル弁開度θｔ（および大気圧Ｐａ、大気温度Ｔａ）とを用いてスロットルモデルＭ１０（図５）によって算出されるスロットル通過空気流量の推定値ｍｔが、エアフロメータ４１の出力から算出した現在のスロットル通過空気流量に一致するようにスロットルモデルの式（１）を修正する。

すなわち、図１７において、実測値Ｐｍとスロットル弁開度θｔとを用いてスロットルモデルＭ１０の式（１）によって算出されたスロットル通過空気流量の推定値ｍｔをエアフロメータモデルＭ４０に入力すると、実際のスロットル通過空気流量が上記スロットル通過空気流量の推定値ｍｔとなっていたと仮定した場合におけるエアフロメータ４１の予想出力値ＡＦＭｍｔに変換される。

したがって、スロットルモデルＭ１０にモデル化誤差が含まれていない場合、すなわち実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとが一致している場合には、実際のエアフロメータ４１の出力値ＡＦＭとエアフロメータモデルＭ４０の予想出力値ＡＦＭｍｔとが一致するはずであり、一致していない場合には実際のスロットル通過空気流量とスロットル通過空気流量の推定値ｍｔとが一致していないと考えられ、よってスロットルモデルＭ１０にモデル化誤差が生じていると考えられる。

したがって、出力値ＡＦＭと予想出力値ＡＦＭｍｔとが一致しない場合には、スロットルモデルＭ１０に生じているモデル化誤差を補償するために、予想出力値ＡＦＭｍｔがエアフロメータ４１の出力値（実測値）ＡＦＭに一致するようにスロットルモデルＭ１０の計算式を修正する必要がある。

ここで、本実施形態では、上述したようにスロットルモデルＭ１０の計算式として、下記式（２１）を用いている。

上記式（２１）のパラメータのうち流量係数μｔは製品毎のばらつきや使用による特性変化の影響を受けやすい。したがって、本実施形態では、スロットルモデルＭ１０の式（２１）により算出されたスロットル通過空気流量の推定値ｍｔをエアフロメータモデルＭ４０の式（１４）および式（１５）に代入して算出された予想出力値ＡＦＭｍｔがエアフロメータ４１の出力値ＡＦＭに一致するように流量係数μｔの値を修正することとしている。これにより、スロットルモデルＭ１０が修正され、スロットルモデルＭ１０を用いた計算結果をエアフロメータ４１の出力値と一致させることができるようになり、よってスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔを実際のスロットル通過空気流量に一致させることができるようになる。なお、当然ながらスロットルモデルＭ１０の上記修正結果は図５筒内充填吸気量Ｍｃの計算にも同時に反映され、筒内充填吸気量Ｍｃの推定精度が向上するようになる。

このような流量係数μｔの値の修正方法としては、例えば、下記式（２２）によるものが考えられる。なお、式（２２）におけるμｔ’は修正後の流量係数である。

ところが、式（２２）により流量係数μｔの値を修正すると、内燃機関が定常運転を行っているときには適切に修正することができるが、内燃機関が過渡運転を行っているときには適切に修正することができない。以下、このことについて図２０を参照して簡単に説明する。

図２０は、スロットル弁開度θｔと空気流量に関する各種パラメータの値とのタイムチャートである。図中、時刻ｔ₀はスロットル弁開度θｔがθｔ₀となったとき、時刻ｔ₁はスロットル弁開度θｔがθｔ₁となったときを示している。また、図２０は、実際のスロットル通過空気流量（図２０中の実線）と、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔ（図２０中の一点鎖線）とが時刻ｔ₀前まではずれているが、時刻ｔ₀以降はほぼ一致している場合のタイムチャートを示している。すなわち、図２０は、時刻ｔ₀前、すなわちスロットル弁開度θｔがθｔ₀よりも小さいときには、流量係数μｔの値に誤差が生じていてスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔが実際のスロットル通過空気流量からずれてしまっていると共に、時刻ｔ₀以降、すなわちスロットル弁開度θｔがθｔ₀以上であるときには、流量係数μｔの値に誤差がほとんどなく、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔが実際のスロットル通過空気流量にほぼ一致している場合のタイムチャートを示している。

上述したように、エアフロメータ４１には放熱遅れに基づく応答遅れがあるため、図２０に示したように実際のスロットル通過空気流量が変化している場合に、エアフロメータ４１の出力値ＡＦＭ（図２０中の破線）には実際のスロットル通過空気流量に対して遅れが生じる。特に、エアフロメータ４１の出力値ＡＦＭを実際のスロットル通過空気流量に対して一次遅れで近似することができることを考えると、エアフロメータ４１の現在の出力値ＡＦＭは現在の実際のスロットル通過空気流量のみに応じて変化するわけではなく、過去の実際のスロットル通過空気流量の履歴にも応じて変化する。

また、上述したようにエアフロメータモデルＭ４０ではエアフロメータ４１の放熱遅れに基づく応答遅れが考慮されているため、図２０に示したように実際のスロットル通過空気流量が変化している場合には、エアフロメータモデルＭ４０によって算出される予想出力値ＡＦＭｍｔにもスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔに対して遅れが生じている。

このため、図２０に示したように、時刻ｔ₀以降においてスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔが実際のスロットル通過空気流量にほぼ一致しているにも関わらず、時刻ｔ₀以降においてエアフロメータ４１の出力値ＡＦＭとエアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔとは一致していない。したがって、例えば、出力値ＡＦＭと予想出力値ＡＦＭｍｔとを用いて式（２２）によってスロットルモデルＭ１０の流量係数μｔを修正しようとすると、時刻ｔ₀以降において流量係数μｔには誤差が無いにも関わらず、流量係数μｔが修正されてしまうことになる。

或いは、流量係数μｔの値の修正方法としては、例えば、下記式（２３）によるものが考えられる。このように、積分を用いることにより、エアフロメータ４１の出力値ＡＦＭに存在するノイズの影響を抑制することができる。

ところが、式（２３）により流量係数μｔの値を修正する場合であっても、内燃機関が定常運転を行っているときには適切に修正することができるが、内燃機関が過渡運転を行っているときには適切に修正することができない。以下、このことについて簡単に説明する。

ここで、流量係数μｔのマップは一定の機関運転領域に区切られており、流量係数μｔのマップ値は機関運転領域毎に設定される。具体的には、流量係数μｔのマップは図２１に示したようにスロットル弁開度θｔを複数のマップ領域に分割し、スロットル弁開度θｔの各マップ領域毎に一つの流量係数μｔのマップ値が設定される。したがって、流量係数μｔを修正する場合にはスロットル弁開度θｔの各マップ領域毎にそのマップ領域に対応する流量係数μｔの値を修正することになる。

図２０に示した例を参照すると、スロットル弁開度θｔの一つのマップ領域がθｔ₀からθｔ₁である場合、スロットル弁開度θｔがそのマップ領域内にある期間は時刻ｔ₀から時刻ｔ₁までとなっている。したがって、このスロットル弁開度θｔがこのマップ領域にある場合における流量係数μｔの値の誤差率を、時刻ｔ₀から時刻ｔ₁までのエアフロメータ４１の出力値ＡＦＭの平均値とエアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔの平均値との比に基づいて算出すると、比較的正確に流量係数μｔの値の誤差を補償することができると考えられ、斯かる考えに基づいて式（２３）による流量係数μｔの値の修正が提案されている。より詳細には、式（２３）では、エアフロメータ４１の出力値ＡＦＭの平均値とその予想出力値ＡＦＭｍｔの平均値との比は、エアフロメータ４１の出力値ＡＦＭの積分値とその予想出力値ＡＦＭｍｔの積分値との比に等しいため、式（２３）では現在マップに保存されている流量係数μｔ（θｔ）に、エアフロメータ４１の出力値ＡＦＭの積分値（図２０中の面積Ｓ₁＋Ｓ₂）を予想出力値ＡＦＭｍｔの積分値（図２０中の面積Ｓ₁）で除算したものを乗算することによって（θｔ・（Ｓ₁＋Ｓ₂）／Ｓ₁）、流量係数μｔ（θｔ）を修正した値μｔ’（θｔ）を算出している。

ところが、上述したように、図２０に示した例では、時刻ｔ₀以降は実際のスロットル通過空気流量とスロットルモデルＭ４０により算出されたスロットル通過空気流量の推定値ｍｔとがほぼ一致しているにも関わらず、エアフロメータ４１に応答遅れが存在するため、エアフロメータ４１の出力値ＡＦＭの時刻ｔ₀から時刻ｔ₁までの積分値（Ｓ₁＋Ｓ₂）とエアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔの時刻ｔ₀から時刻ｔ₁までの積分値（Ｓ₁）とが一致しない。このため、スロットル弁開度θｔがθｔ₁からθｔ₂であるマップ領域における流量係数μｔには誤差が無いにも関わらず流量係数μｔが修正されてしまうことになる。

このように、式（２２）および式（２３）により流量係数μｔの値を修正すると、内燃機関が定常運転を行っているときには的確に修正することができるが、内燃機関が過渡運転を行っているときには的確に修正することができない。そこで、本発明では、内燃機関が過渡運転を行っているときにおいても的確に流量係数μｔの値を修正することができるような内燃機関の吸入空気量推定装置を提供する。

図２２は、図２０と同様に、スロットル弁開度θｔと空気流量に関する各種パラメータの値とのタイムチャートである。以下では、図２２を参照して、スロットル弁開度θｔがθｔ₀からθｔ₁であるマップ領域Ａにおける、流量係数μｔの修正について考える。式（１４）および式（１５）に示したように、エアフロメータモデルＭ４０では、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔに対して一次遅れが存在するとしてエアフロメータ４１の予想出力値ＡＦＭｍｔを算出している。このため、エアフロメータモデルＭ４０の予想出力値ＡＦＭｍｔは、スロットルモデルＭ１０によって算出された現在のスロットル通過空気流量の推定値ｍｔだけでなく所定時間前のスロットル通過空気流量の推定値ｍｔに基づいて算出されることになる。

例えば、図２２に示したように時刻ｔ₀においてスロットル弁開度θｔがマップ領域Ａに侵入したとき、時刻ｔ₀以降においてエアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔは、時刻ｔ₀以降においてスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔだけでなく、時刻ｔ₀以前にスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔを用いて算出されることになる。すなわち、エアフロメータモデルＭ４０によって算出される予想出力値ＡＦＭｍｔは、スロットル弁開度θｔがマップ領域Ａ以外のマップ領域にあるときのスロットル通過空気流量の推定値ｍｔに基づいて算出されることになる。

ここで、上述したように図２２に示した例では、スロットル弁開度θｔがθｔ₀からθｔ₁であるマップ領域Ａに対応する流量係数μｔの修正を行おうとしており、流量係数μｔの修正を行うに当たって実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとの比をエアフロメータ４１の出力値ＡＦＭとエアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔとの比から算出しようとしている。ところが、スロットル弁開度θｔがマップ領域Ａにある場合においてエアフロメータモデルＭ４０によって算出される予想出力値ＡＦＭｍｔは、マップ領域Ａ以外のマップ領域におけるスロットル通過空気流量の推定値ｍｔの影響を受けており、よってスロットル弁開度θｔがマップ領域Ａにある場合におけるエアフロメータ４１の出力値ＡＦＭとエアフロメータモデルＭ４０によって算出される予想出力値ＡＦＭｍｔとの比は、スロットル弁開度θｔがマップ領域Ａにある場合における実際のスロットル通過空気流量とスロットル通過空気流量の推定値ｍｔとの比を表したものではないものとなっている。

そこで、本実施形態では、スロットル弁開度θｔがマップ領域Ａにある場合における流量係数μｔを修正するときには、エアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔを直接利用せずに、この予想出力値ＡＦＭｍｔからスロットル弁開度θｔがマップ領域Ａに侵入する前（図２２に示した例では、スロットル弁開度θｔがθｔ₀よりも小さいとき）のスロットル通過空気流量の推定値ｍｔの影響を除去した値（以下、「修正出力値」と称す）ＡＦＭｍｔｅを用いることとする。

本実施形態では、修正出力値ＡＦＭｍｔｅの算出は基本的にエアフロメータモデルＭ４０における予想出力値ＡＦＭｍｔの算出と同様に行われる。すなわち、下記式（２４）および式（２５）によりボビン部４１ａ３およびサポート部４１ａ４の遅れ放熱量ω_be、ω_seが算出され、算出された遅れ放熱量ω_be、ω_seに基づいて式（２６）により図１９に示したマップと同様なマップから修正出力値ＡＦＭｍｔｅが算出される。式（２４）および式（２５）は、ＥＣＵ３１で処理するために時刻ｔ、計算間隔Δｔを用いて下記式（２７）および式（２８）のように離散化されて用いられる。

ただし、エアフロメータモデルＭ４０における予想出力値ＡＦＭｍｔの算出と異なり、修正出力値ＡＦＭｍｔｅの算出はマップ領域毎に行われる。すなわち、スロットル弁開度θｔが或るマップ領域に侵入したとき、図２２に示した例ではスロットル弁開度θｔがマップ領域Ａに侵入したとき（すなわちスロットル弁開度θｔがθｔ₀となったとき）にそのマップ領域における修正出力値ＡＦＭｍｔｅの算出が開始されると共に、スロットル弁開度θｔがそのマップ領域から離脱したとき、図２２に示した例ではスロットル弁開度θｔがマップ領域Ａから離脱したとき（すなわちスロットル弁開度θｔがθｔ₁となったとき）にそのマップ領域における修正出力値ＡＦＭｍｔｅの算出が終了せしめられる。

そして、修正出力値ＡＦＭｍｔｅの算出にあたっては、スロットル弁開度θｔがマップ領域Ａに侵入する前にスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔ（以下、「過去のスロットル通過空気流量の推定値ｍｔ」と称す）の影響を除去するために、スロットル弁開度θｔがマップ領域Ａに侵入したときの修正出力値ＡＦＭｍｔｅの値（以下、「初期値」と称す）を、エアフロメータ４１の出力値ＡＦＭとほぼ同様な値になるようにしている。

すなわち、上述したように修正出力値ＡＦＭｍｔｅはＥＣＵ３１において上記離散化された式（２７）および式（２８）に基づいて算出されるが、これら式（２７）および式（２８）には初期値以外に過去の（スロットル弁開度θｔがマップ領域Ａに侵入する前の）スロットル通過空気流量の推定値ｍｔが影響するような項は存在しない。逆に言うと、初期値を過去のスロットル通過空気流量の推定値ｍｔの影響を受けていない値に設定しさえすれば、修正出力値ＡＦＭｍｔｅに対する過去のスロットル通過空気流量の推定値ｍｔの影響を最小限に抑制することができる。

本実施形態では、スロットル弁開度θｔがマップ領域Ａに侵入したときの修正出力値ＡＦＭｍｔｅの初期値をそのときのエアフロメータ４１の出力値ＡＦＭとほぼ同様な値としているため、修正出力値ＡＦＭｍｔｅの初期値は過去のスロットル通過空気流量の推定値ｍｔに基づいて定まる値ではなく、スロットル弁開度θｔがマップ領域Ａに侵入したときのエアフロメータ４１の出力値ＡＦＭに基づいて定まる値であり、よって修正出力値ＡＦＭｍｔｅから過去のスロットル通過空気流量の推定値ｍｔの影響が除去される。

このように、修正出力値ＡＦＭｍｔｅの初期値をエアフロメータ４１の出力値ＡＦＭとほぼ同様な値として式（２４）および式（２５）によってその後の修正出力値ＡＦＭｍｔを算出するためには、時刻ｔ₀におけるボビン部４１ａ３の遅れ放熱量およびサポート部４１ａ４の遅れ放熱量を時刻ｔ₀におけるボビン部４１ａ３からの実際の放熱量およびサポート部４１ａ４からの実際の放熱量に一致させることが必要である。ところが、時刻ｔ₀におけるエアフロメータ４１の出力値ＡＦＭから時刻ｔ₀におけるボビン部４１ａ３およびサポート部４１ａ４からの実際の放熱量を算出することは困難である。すなわち、図１９に示したマップ等によってエアフロメータ４１の出力値ＡＦＭからボビン部４１ａ３およびサポート部４１ａ４からの実際の放熱量を合計した値を求めることはできるが、ボビン部４１ａ３およびサポー部４１ａ４それぞれからの実際の放熱量を求めることはできない。

そこで、本実施形態では、エアフロメータ４１のボビン部４１ａ３およびサポート部４１ａ４それぞれからの実際の放熱量を比較的正確に算出することができるオブザーバを用いるようにしている。このオブザーバは、時刻やスロットル弁開度θｔのマップ領域とは無関係に、常にエアフロメータ４１のボビン部４１ａ３およびサポート部４１ａ４それぞれからの実際の放熱量の近似値を算出している。

具体的には、オブザーバは、下記式（２９）および式（３０）によってボビン部４１ａ３からの実際の放熱量の近似値（以下、「近似放熱量」と称す）ω_boおよびサポート部４１ａ４からの実際の放熱量の近似値（以下、「近似放熱量」と称す）ω_soをそれぞれ算出している。式（２９）および式（３０）において、Ｋ_Iは積分制御等に用いられる積分係数と同様な係数である。また、ＡＦＭｍｔｏはオブザーバによって算出された近似放熱量ω_bo、ω_soを合計した値に基づいて下記式（３１）により（すなわち、図１９のマップにより）算出された値である。

したがって、オブザーバでは、エアフロメータ４１の出力値ＡＦＭとオブザーバによって算出された近似放熱量ω_bo、ω_soに基づいて算出されるエアフロメータ４１の出力値（以下、「近似出力値」と称す）ＡＦＭｍｔｏとの偏差（ＡＦＭ−ＡＦＭｍｔｏ）を時間的に累積する。この出力値ＡＦＭと近似出力値ＡＦＭｍｔｏとの偏差の累積値は、実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとの誤差に基づいて変化し、この誤差の絶対値が大きいほど上記偏差の累積値の絶対値が大きくなる。したがって、スロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔに、出力値ＡＦＭとオブザーバによって算出された近似出力値ＡＦＭｍｔｏとの偏差の累積値に係数Ｋ_Iを乗じた値を加算することで、実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値ｍｔとの誤差を或る程度補償することができる。これにより、オブザーバでは、実際のスロットル通過空気流量とほぼ等しいスロットル通過空気流量に基づいて近似放熱量ω_bo、ω_soが算出されるため、オブザーバによって算出された近似放熱量ω_bo、ω_soはボビン部４１ａ３およびサポート部４１ａ４からの実際の放熱量とほぼ等しい値となっている。

なお、上記オブザーバを内燃機関の制御装置に実装して、スロットル通過空気流量の推定値ｍｔおよびエアフロメータ４１の出力値ＡＦＭから近似放熱量ω_bo、ω_soを算出するためには、式（２９）および式（３０）を離散化する必要がある。時刻ｔ、計算間隔Δｔを用いて式（２９）および式（３０）を離散化すると、それぞれ下記式（３２）および式（３３）が得られる。なお、式（３２）および式（３３）中のΔＡＦＭＩＮＴ（ｔ）は、時刻ｔにおけるＡＦＭ（ｔ）−ＡＦＭｍｔｏ（ｔ）の累積値であり、下記式（３４）により表される。

このようにして実装されたオブザーバでは、時刻ｔにおけるボビン部４１ａ３およびサポート部４１ａ４からの近似放熱量ω_bo（ｔ）およびω_so（ｔ）、スロットルモデルＭ１０の式（８）によって算出された時刻ｔにおけるスロットル通過空気流量の推定値ｍｔ（ｔ）、時刻ｔにおけるエアフロメータ４１の出力値ＡＦＭ（ｔ）、および時刻ｔにおいてオブザーバによって算出された近似出力値ＡＦＭｍｔｏ（ｔ）が式（３２）、式（３３）および式（３４）に代入され、これにより時刻ｔ＋Δｔにおける近似放熱量ω_bo（ｔ＋Δｔ）およびω_so（ｔ＋Δｔ）が算出される。そして、これら近似放熱量ω_bo（ｔ＋Δｔ）およびω_so（ｔ＋Δｔ）に基づいて式（３１）により時刻ｔ＋Δｔにおける近似出力値ＡＦＭｍｔｏ（ｔ＋Δｔ）が算出され、この時刻ｔ＋Δｔにおける近似出力値ＡＦＭｍｔｏ（ｔ＋Δｔ）はその後の近似放熱量ω_boおよびω_soの算出に利用される。そして、このような計算を繰り返すことによって、スロットル通過空気流量の推定値ｍｔ、エアフロメータ４１の出力値ＡＦＭから、任意の時刻ｔにおける近似放熱量ω_bo、ω_soが算出される。

図２３は、図２２の時刻ｔ₀〜ｔ₁における流量に関するパラメータの推移を拡大して示したタイムチャートである。図２３から分かるように、オブザーバによって算出された近似出力値ＡＦＭｍｔｏ（図２３中の短い破線）は、エアフロメータ４１の出力値ＡＦＭとほぼ同様な値となっている。上述したように修正出力値ＡＦＭｍｔｅの算出にあたっては、スロットル弁開度θｔがマップ領域Ａに侵入したときに、すなわち時刻ｔ₀においてオブザーバによって算出されたボビン部４１ａ３の近似放熱量ω_boおよびサポート部４１ａ４の近似放熱量ω_soが用いられる。したがって、修正出力値ＡＦＭｍｔｅの算出においては、遅れ放熱量ω_be、ω_seの初期値、すなわち時刻ｔ₀における遅れ放熱量ω_be（ｔ₀）、ω_se（ｔ₀）はオブザーバによって算出された時刻ｔ₀における近似放熱量ω_bo（ｔ₀）、ω_se（ｔ₀）とされる（ω_be（ｔ₀）＝ω_bo（ｔ₀）、ω_se（ｔ₀）＝ω_so（ｔ₀））。このため、時刻ｔ₀における修正出力値ＡＦＭｍｔｅは、時刻ｔ₀におけるオブザーバの近似出力値ＡＦＭｍｔｏと同一の値となる。

そして、修正出力値ＡＦＭｍｔｅの算出に当たっては、遅れ放熱量ω_be、ω_seの初期値、すなわち時刻ｔ₀における遅れ放熱量ω_be（ｔ₀）、ω_se（ｔ₀）が算出された後に、上述した式（２７）および式（２８）によって時刻ｔ₀以降における遅れ放熱量ω_be（ｔ）、ω_se（ｔ）が算出され、そしてこれら時刻ｔ₀以降の時刻ｔにおける遅れ放熱量ω_be（ｔ）、ω_se（ｔ）に基づいて上記式（２６）を用いて時刻ｔ₀以降の時刻ｔにおける修正出力値ＡＦＭｍｔｅ（ｔ）が算出される。このようにして算出された修正出力値ＡＦＭｍｔｅは図２３に点線で示したように推移する。すなわち、修正出力値ＡＦＭｍｔｅは時刻ｔ₀においてはオブザーバによって算出された近似出力値ＡＦＭｍｔｏと同一の値となっており、その後徐々にエアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔに近づくように推移する。

このようにして算出された修正出力値ＡＦＭｍｔｅは、上述したようにスロットル弁開度θｔがマップ領域Ａに侵入する前のスロットル通過空気流量の推定値ｍｔの影響が除去された値となっている。このため、マップ領域Ａにおけるエアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとの比は、マップ領域Ａにおける実際のスロットル通過空気流量とスロットルモデルＭ１０によって算出されたスロットル通過空気流量との比を比較的正確に示しており、すなわちマップ領域ＡにおけるスロットルモデルＭ１０の流量係数μｔの誤差率を示している。したがって、上記式（２３）と同様に、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとに基づいて下記式（３５）により誤差率Ｅｒｒが算出される。

すなわち、誤差率Ｅｒｒは、スロットル弁開度θｔが特定のマップ領域（図２２、２３に示した例ではマップ領域Ａ）にある間に亘ってエアフロメータ４１の出力値ＡＦＭを積分した値、すなわち時刻ｔ₀からｔ₁に亘って出力値ＡＦＭを積分した値（図２３中の面積Ｓ₃＋Ｓ₄）を、スロットル弁開度θｔが上記特定のマップ領域（図２２、２３に示した例ではマップ領域Ａ）にある間に亘って修正出力値ＡＦＭｍｔｅを積分した値、すなわち修正出力値ＡＦＭｍｔｅを時刻ｔ₀からｔ₁に亘って積分した値（図２３中の面積Ｓ₃）で除算した値（（Ｓ₃＋Ｓ₄）／Ｓ₃）に等しい。このようにして算出された誤差率Ｅｒｒを上記特定のマップ領域（図２２、２３に示した例ではマップ領域Ａ）の流量係数μｔのマップ値に乗算することで、流量係数μｔのマップ値が的確に修正される。

なお、上記実施形態では、スロットル弁開度θｔがθｔ₀〜θｔ₁であるマップ領域Ａに対応する流量係数μｔのマップ値を修正する場合について示したが、他のマップ領域においても同様にそのマップ領域に対応する流量係数μｔのマップ値が修正せしめられる。これにより、内燃機関の過渡運転中であっても各マップ領域に対応する流量係数μｔのマップ値が的確に修正せしめられる。

ところで、上述したように、修正出力値ＡＦＭｍｔｅを用いることでスロットル弁開度θｔがマップ領域Ａに侵入する前のスロットル通過空気流量の推定値ｍｔの影響が除去される。ところが、エアフロメータ４１の出力値ＡＦＭもスロットル弁開度θｔがマップ領域Ａに侵入する前の実際のスロットル通過空気流量の影響を受けており、上述したように修正出力値ＡＦＭｍｔｅを用いたとしてもエアフロメータ４１の出力値ＡＦＭから過去の実際のスロットル通過空気流量の影響を除去することができない。このため、上述したように、エアフロメータ４１の出力値ＡＦＭと修正出力値ＡＦＭｍｔｅとに基づいて算出されたスロットルモデルＭ１０の流量係数μｔの誤差率は、実際の流量係数μｔの誤差率に対して僅かに異なった値となる。

そこで、上記実施形態の変更例では、エアフロメータ４１の出力値ＡＦＭから過去の実際のスロットル通過空気流量の影響を除去することとしている。以下、エアフロメータ４１の出力値ＡＦＭから過去の実際のスロットル通過空気流量の影響を除去してスロットルモデルＭ１０の流量係数μｔの誤差率を正確に算出する方法について説明する。

ところで、一般に、エアフロメータ４１等のセンサの応答は、初期値応答とゼロ応答とに分離することができる。ここで、初期値応答とは、或る時刻ｔ以降の状態量がゼロであると仮定した場合における時刻ｔ以降のセンサの出力値の推移であり、ゼロ応答とは、或る時刻ｔまで状態量がゼロであると仮定した場合における時刻ｔ以降のセンサの出力値の推移である。

図２４は、エアフロメータ４１の出力値ＡＦＭについての初期値応答曲線とゼロ応答曲線とを示す図である。ここで、Φ₁は時刻０を基準とした初期値応答曲線（破線）を示しており、時刻０以降に実際のスロットル通過空気流量がゼロであると仮定した場合におけるエアフロメータ４１の出力値の推移を示している。また、Φ₂は時刻０を基準としたゼロ応答曲線（一点鎖線）を示しており、時刻０までの実際のスロットル通過空気流量がゼロであると仮定した場合、すなわち時刻０からエアフロメータ４１での計測が開始された場合におけるエアフロメータ４１の出力値の推移を示している。

すなわち、初期値応答曲線は時刻０におけるエアフロメータ４１の出力値に応じて、すなわち時刻０よりも前の実際のスロットル通過空気流量に応じて変化するのに対して、ゼロ応答曲線は時刻０以降における実際のスロットル通過空気流量に応じて変化する。そして、時刻０以降の各時刻ｔにおけるエアフロメータ４１の出力値ＡＦＭ（ｔ）は、初期値応答曲線の時刻ｔにおける値Φ₁（ｔ）とゼロ応答曲線の時刻ｔにおける値Φ₂（ｔ）とを合計した値に等しい（ＡＦＭ（ｔ）＝Φ₁（ｔ）＋Φ₂（ｔ））。そして、このようなエアフロメータ４１の出力値ＡＦＭと、初期値応答曲線の各時刻における値Φ₁と、ゼロ応答曲線の各時刻における値Φ₂との関係は、基準となる時刻（図２４に示した例では時刻０）を如何なる時刻にしても成立する。

上述したようにゼロ応答曲線は所定時刻以降における実際のスロットル通過空気流量に応じて変化するため、特定のマップ領域に侵入した時刻を基準としたエアフロメータ４１の出力値ＡＦＭのゼロ応答曲線からは、その特定のマップ領域に侵入する前の実際のスロットル通過空気流量の影響が除去されている。したがって、斯かるゼロ応答曲線を求め、これと修正出力値ＡＦＭｍｔｅについてのゼロ応答曲線（正確には、修正出力値ＡＦＭｍｔｅについてのゼロ応答曲線に対応する曲線）とを後述するように比較することにより、上記特定のマップ領域に侵入する前の実際のスロットル通過空気流量の影響を除去して、流量係数μの誤差率を的確に算出することができる。

ただし、時刻ｔ₀を基準としたゼロ応答を算出するのは困難であることから、本実施形態では、エアフロメータ４１の出力値ＡＦＭから時刻ｔ₀を基準とした初期値応答曲線の各時刻における値Φ₁を減算することでゼロ応答を算出することとしている（Φ₂（ｔ）＝ＡＦＭ（ｔ）−Φ₁（ｔ））。

ここで、時刻ｔ₀を基準とした初期値応答曲線の各時刻における値Φ₁の算出方法について説明する。上述したように初期値応答曲線の各時刻における値Φ₁は、時刻ｔ₀以降に実際のスロットル通過空気流量がゼロであると仮定した場合におけるエアフロメータ４１の出力値の推移を示しているため、上記式（１４）および式（１５）のスロットル通過空気流量の推定値ｍｔをゼロとすることにより求めることができる。すなわち、初期値応答曲線についての放熱量ω_bi、ω_siは下記式（３６）および式（３７）により表される。この場合、時刻ｔ₀における放熱量ω_bi（ｔ₀）、ω_si（ｔ₀）は、時刻ｔ₀における実際の放熱量の近似値である近似放熱量ω_bo（ｔ₀）、ω_so（ｔ₀）が用いられる。

そして、式（３６）、（３７）により算出された放熱量ω_bi、ω_siに基づいて、式（３８）により図１９に示したマップと同様なマップから初期値応答曲線の各時刻における値Φ₁が算出される。なお、式（３６）および式（３７）は、ＥＣＵ３１で処理するために離散化されて用いられる。

このようにして算出された時刻ｔ₀を基準とした初期値応答曲線の各時刻における値Φ₁をエアフロメータ４１の出力値ＡＦＭから減算することにより、エアフロメータ４１の出力値についての時刻ｔ₀を基準としたゼロ応答曲線の各時刻における値Φ₂が算出される。

また、修正出力値ＡＦＭｍｔｅも時刻ｔ₀において近似放熱量ω_bo、ω_soを利用しているため、修正出力値ＡＦＭｍｔｅから上述した初期値応答曲線の各時刻における値Φ₁を減算することにより、修正出力値ＡＦＭｍｔｅについてのゼロ応答曲線の各時刻における値が算出される。このため、本実施形態では、スロットルモデルＭ１０の流量係数μｔの誤差率Ｅｒｒを算出するにあたって、エアフロメータ４１の出力値ＡＦＭと、修正出力値ＡＦＭｍｔｅと、エアフロメータ４１の出力値についての時刻ｔ₀を基準とした初期値応答曲線の各時刻における値Φ₁とを用いて下記式（３９）により算出される。

このように、流量係数μｔの誤差率を算出するにあたって、修正出力値ＡＦＭｍｔｅを用いることにより、エアフロメータモデルＭ４０によって算出された予想出力値ＡＦＭｍｔに存在する過去のスロットル通過空気流量の推定値ｍｔの影響が除去され、すなわち所定時間前のスロットル通過空気流量の影響で予想出力値ＡＦＭｍｔに含まれている成分を予想出力値ＡＦＭｍｔから除去される。さらに、エアフロメータ４１の出力値ＡＦＭおよび修正出力値ＡＦＭｍｔｅから初期値応答曲線の各時刻における値Φ₁を減算することにより、エアフロメータ４１の出力値ＡＦＭ等に存在する過去の実際のスロットル通過空気流量の影響が除去され、すなわち所定時間前のスロットル通過空気流量の影響でエアフロメータの出力値ＡＦＭに含まれている成分をこのエアフロメータの出力値から除去することができる。

このため、上記式（３９）により算出された誤差率Ｅｒｒはスロットル弁開度θｔが所定のマップ領域内にあるときの実際のスロットル通過空気流量およびスロットルモデルＭ１０によって算出されたスロットル通過空気流量の推定値に基づいて算出されることになるため、スロットル弁開度θｔが所定のマップ領域内にあるときの流量係数の誤差率を正確に算出することができる。

次に、吸気弁モデルＭ３０の修正操作について説明する。上述したように、吸気弁モデルＭ３０の式（７）において用いられる適合値ａは比例係数であり、また適合値ｂは燃焼室５内に残存していた既燃ガスを表す値である。これら適合値ａ、ｂは、機関回転数ＮＥ、吸気弁の位相角および作用角（以下では、吸気弁の開弁時期ＶＶＴを例にとって説明する）に応じて変化する。このため、機関回転数ＮＥおよび吸気弁の開弁時期ＶＶＴとこれら適合値ａ、ｂとの関係を予め実験的にまたは計算により算出し、図２５（ａ）および図２５（ｂ）に示したようにマップとしてＥＣＵ３１に保存すると共に、吸気弁モデルＭ３０において上記式（７）の計算をするにあたってはクランク角センサ４７によって検出された機関回転数ＮＥおよび吸気弁制御装置１３への吸気弁開弁時期ＶＶＴを用いて上記マップによりこれら適合値ａ、ｂを算出している。なお、図２５（ａ）は適合値ａのマップを、図２５（ｂ）は適合値ｂのマップをそれぞれ示しており、図中のｘ軸は機関回転数ＮＥを、ｙ軸は吸気弁開弁時期ＶＶＴ（ｙ軸が大きくなるほど進角）をそれぞれ示している。

ところが、これら適合値ａ、ｂには、上述したように使用による特性の変化や製作公差によるばらつきなどにより誤差が生じる場合がある。このように、適合値ａ、ｂに誤差が生じている場合には、吸気弁モデルＭ３０のモデル化誤差を補償するためにこれら適合値ａ、ｂのマップ値を修正していく必要がある。ここで、本実施形態では、適合値ａ、ｂは、機関運転領域毎に定められている。図２５に示した例では、適合値ａ、ｂは、マップ領域（例えば、図２５中の斜線部分が一つのマップ領域を表している）毎に、すなわち一定間隔の機関回転数毎および一定間隔の吸気弁開弁時期毎に定められている。したがって適合値ａ、ｂの修正は、各マップ領域毎に行われなければならない。

そこで、本実施形態では、吸気弁モデルＭ３０の式（１１）で用いられる適合値ａ、ｂをマップ領域毎に修正することとしている。以下では、例として、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴがマップ領域Ｂ（すなわち、機関回転数ＮＥがＮＥ₁からＮＥ₂の間であって且つ吸気弁開弁時期ＶＶＴがＶＶＴ₁からＶＶＴ₂の間であるマップ領域）にある場合における適合値ａ、ｂのマップ値ａ_ij、ｂ_ijの修正方法について説明する。

図２６は、内燃機関の運転中における吸気管内圧力Ｐｍ、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴのタイムチャートである。図示した例では、時刻ｔ₀以降において機関回転数ＮＥがＮＥ₀以上となっており、また時刻ｔ₁以降において機関回転数ＮＥがＮＥ₁以上となっている。すなわち、機関回転数ＮＥは時刻ｔ₀から時刻ｔ₁までに亘ってＮＥ₀とＮＥ₁との間の値となっており、それ以外の時刻においてはＮＥ₀とＮＥ₁との間の値となっていない。一方、吸気弁開弁時期ＶＶＴは時刻ｔ₀から時刻ｔ₁までに亘ってＶＶＴ₁とＶＶＴ₂との間の値となっている。

したがって、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴの値は、時刻ｔ₀以前はマップ領域Ｂ内にはなく、時刻ｔ₀においてマップ領域Ｂに侵入すると共に時刻ｔ₀から時刻ｔ₁までに亘ってマップ領域Ｂ内にある。その後、時刻ｔ₁においてマップ領域Ｂから離脱して、時刻ｔ₁以降はマップ領域Ｂ内にはない。そして、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴの値がマップ領域Ｂ内にあるとき、すなわち時刻ｔ₀から時刻ｔ₁までの間において、吸気弁モデルＭ３０ではマップ領域Ｂに対応した適合値ａ_ij、ｂ_ijが用いられる。

ここで、上記式（６）について再び考えると、式（６）の左辺は吸気管部分２３の気体の内部エネルギの時間的変化量を表しており、式（６）の右辺のＣｐ・ｍｔ・Ｔａは吸気管部分２３に流入する気体のエネルギを、Ｃｐ・ｍｃ・Ｔｍは吸気管部分２３から流出した気体のエネルギをそれぞれ表している。ここで、吸気管部分２３における気体の状態方程式（Ｐｍ・Ｖｍ＝Ｍ・Ｒ・Ｔｍ）を考慮すると、式（６）の左辺は下記式（４０）のように表せる。

したがって、時刻ｔ₀から時刻ｔ₁までの間の内部エネルギの変化量は下記式（４１）のように表せる。特に、吸気管内圧力センサ４０により実際の吸気管内圧力Ｐｍｒが検出されているため、式（４１）において実際の吸気管内圧力Ｐｍｒを用いることによって、時刻ｔ₀から時刻ｔ₁までの間の実際の内部エネルギの変化量を正確に算出することができる。なお、式（４１）においてＰｍｒ（ｔ₀）は時刻ｔ₀における実際の吸気管内圧力を、Ｐｍｒ（ｔ₁）は時刻ｔ₁における実際の吸気管内圧力をそれぞれ示している。

一方、吸気管部分２３に流入する気体のエネルギ（すなわち、式（６）の左辺のＣｐ・ｍｔ・Ｔａ）に関しては、上述したようにスロットルモデルＭ１０の流量係数μｔのマップ値が逐次修正されるため、スロットルモデルＭ１０の流量係数μｔは適切な値となっている。したがって、スロットルモデルＭ１０の式（１）の吸気管内圧力Ｐｍとして吸気管内圧力センサ４０によって検出された実際の吸気管内圧力Ｐｍｒを用いると、吸気管部分２３に流入する気体のエネルギを正確に算出することができる。特に、時刻ｔ₀から時刻ｔ₁までに吸気管部分２３に流入する気体のエネルギは下記式（４２）により算出することができる。

このように時刻ｔ₀から時刻ｔ₁までの間の内部エネルギの変化量が上記式（４１）で、時刻ｔ₀から時刻ｔ₁までに吸気管部分２３に流入する気体のエネルギが上記式（４２）でそれぞれ正確に算出されることから、式（４２）から式（４１）を減算した下記式（４３）により時刻ｔ₀から時刻ｔ₁までに吸気管部分２３から流出した気体のエネルギを正確に算出することができる。

ここで、上述したように吸気管部分２３から流出した気体のエネルギをＣｐ・ｍｃ・Ｔｍで表すことができ、またこのｍｃとして吸気弁モデルＭ３０によって算出された筒内吸入空気流量ｍｃ（すなわち、式（７））を用いると、時刻ｔ₀から時刻ｔ₁までに吸気管部分２３から流出した気体のエネルギは下記式（４４）のように表せる。そして、式（４３）及び式（４４）はいずれも吸気管部分２３から流出した気体のエネルギを表しているため、本来であればこれら式によって算出された値は一致するはずである。

ところが、吸気弁モデルＭ３０のモデル式（７）に誤差があると、式（４３）と式（４４）とによって算出された値は一致しない。特に上記式（４４）では吸気管内圧力として吸気管内圧力センサ４０により検出された実際の吸気管内圧力Ｐｍｒが用いられていることを考慮すると、式（４３）と式（４４）とが一致しない場合は式（７）で用いられる適合値ａ、ｂのマップ値に誤差がある場合であると考えられる。

そこで、本実施形態では、実際の筒内吸入空気流量に対する吸気弁モデルＭ３０によって算出された筒内吸入空気流量ｍｃの誤差率（実際の筒内吸入空気流量／吸気弁モデルＭ３０によって算出された筒内吸入空気流量）をαとする。この誤差率αを式（４４）に乗算して算出された値と式（４３）によって算出された値は等しくなるため、下記式（４５）が成立する。そして、式（４５）を誤差率αについて解くと、誤差率αは下記式（４６）のように表せる。

このようにして算出された誤差率αは、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴがマップ領域Ｂにある場合における実際の筒内吸入空気流量と吸気弁モデルＭ３０によって算出された筒内吸入空気流量の推定値ｍｃとの誤差率を表している。そして、この筒内吸入空気流量における誤差を補償するためには、本実施形態では、マップ領域Ｂに対応するマップ値ａ_ij、ｂ_ijに誤差率αを乗算したものを新たな適合値ａ、ｂのマップ領域Ｂにおけるマップ値ａ_ij、ｂ_ijとすることとしている（ａ_ij＝α・ａ_ij、ｂ_ij＝α・ｂ_ij）。このように吸気弁モデルＭ３０の式（７）に用いられる適合値ａ、ｂのマップ値を修正することにより、吸気弁モデルＭ３０によって算出される筒内吸入空気流量の推定値ｍｃを実際の筒内吸入空気流量とほぼ等しい値とすることができる。

このように、本実施形態では、機関運転状態（例えば、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴ等）が特定の機関運転領域（例えば、上記マップ領域Ｂ）に侵入してから離脱するまでに吸気管部分２３から流出した気体のエネルギの実際の値（すなわち、吸気管部分２３のエネルギ保存則を用いて算出された値）と、斯かる気体のエネルギの吸気弁モデルＭ３０を用いて算出された値とが求められ、実際の値と吸気弁モデルＭ３０を用いて算出された値との比率に基づいて上記特定の機関運転領域に対応する適合値ａ、ｂのマップ値を修正することにより、内燃機関が過渡運転を行っている場合であっても適合値ａ、ｂのマップ値を的確に修正することができる。

ここで、吸気管モデルＭ２０によって算出された吸気管内圧力の推定値Ｐｍは所定時間前に吸気弁モデルＭ３０によって算出された筒内吸入空気流量の推定値ｍｃを用いて算出されている。このため、例えば吸気管内圧力の実測値と吸気弁モデルＭ３０によって算出された筒内吸入空気流量の推定値ｍｃを用いて吸気管モデルＭ２０によって算出された吸気管内圧力の推定値との偏差に基づいて適合値ａ、ｂのマップ値を修正しようとすると、吸気管内圧力の実測値が現在の実際の筒内吸入空気流量に対応した値であるのに対して吸気管モデルＭ２０によって算出された吸気管内圧力の推定値は吸気弁モデルＭ３０によって算出された過去の筒内吸入空気流量の推定値ｍｃに対応した値となっているため、内燃機関が定常運転を行っているときには比較的的確に吸気弁モデルを修正することができるが、内燃機関が過渡運転を行っているときには的確に吸気弁モデルを修正することができない。

一方、本実施形態では、機関運転状態（例えば、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴ等）が特定の機関運転領域（例えば、上記マップ領域Ｂ）に侵入してから離脱するまでに吸気管部分２３から流出した気体のエネルギの実際の値（すなわち、吸気管部分２３のエネルギ保存則を用いて算出された値）と、斯かる気体のエネルギの吸気弁モデルＭ３０を用いて算出された値とを比較し、これら値が一致するように適合値ａ、ｂのマップ値を修正するようにしている。このように機関運転状態が特定の機関運転領域にある間に吸気管部分２３から流出した気体のエネルギの総量同士を比較することにより、すなわち斯かる間の単位時間当たりのエネルギ流出量の積分値同士を比較することにより、内燃機関が過渡運転を行っている場合であっても、この特定の機関運転領域に対応する適合値ａ、ｂのマップ値を的確に修正することができるようになる。

次に、上記式（４６）を内燃機関の制御装置に実装して、実際にスロットル通過空気流量の推定値ｍｔおよび吸気管内圧力Ｐｍの実測値から誤差率αを算出する場合について説明する。ＥＣＵ３１では式（４６）中の積分項をそのまま算出することはできない。このため、式（４６）中のスロットル通過空気流量の推定値ｍｔについての積分項は時刻ｔ、計算間隔Δｔを用いて離散化された下記式（４７）によって算出され、吸気管内圧力Ｐｍについての積分項は時刻ｔ、計算間隔Δｔを用いて離散化された下記式（４８）によって算出される。

そして、スロットル通過空気流量の推定値ｍｔについての積分項に関しては、時刻ｔにおけるスロットル通過空気流量の推定値ｍｔ（ｔ）が式（４７）に代入され、これにより時刻ｔにおけるスロットル通過空気流量の推定値ｍｔについての積分項の値ｍｔｉｎｔ（ｔ）が算出される。また、吸気管内圧力Ｐｍについての積分項に関しては、時刻ｔにおける吸気管内圧力Ｐｍ（ｔ）が式（４８）に代入され、これにより時刻ｔにおける吸気管内圧力Ｐｍについての積分項の値ｈｉｎｔ（ｔ）が算出される。なお、いずれの積分項の算出にあたっても、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが対象となるマップ領域に侵入したとき、すなわち時刻ｔ₀におけるスロットル通過空気流量の推定値ｍｔ（ｔ₀）および吸気管内圧力Ｐｍ（ｔ₀）を初期値として上記式（４７）および式（４８）の計算が行われる。

そして、式（４６）による誤差率αの算出は、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが対象となるマップ領域から離脱した後、すなわち時刻ｔ₁以降に行われる。式（４７）についての初期値が時刻ｔ₀におけるスロットル通過空気流量の推定値ｍｔ（ｔ₀）であることから、時刻ｔ₁におけるスロットル通過空気流量の推定値ｍｔについての積分項の値ｍｔｉｎｔ（ｔ₁）は時刻ｔ₀から時刻ｔ₁までのｍｔの積分値を表しており、また、式（４８）についての初期値が時刻ｔ₀における吸気管内圧力Ｐｍ（ｔ₀）であることから、時刻ｔ₁における吸気管内圧力Ｐｍについての積分項の値ｈｉｎｔ（ｔ₁）は時刻ｔ₀から時刻ｔ₁までのａ・Ｐｍ＋ｂの積分値を表している。このため、上記マップ領域における誤差率αは下記式（４９）により算出せしめられる。

ところで、上記説明では、吸気管内圧力センサ４０によって実測された吸気管内圧力を用いて式（４６）により誤差率αを算出することとしている。ところが、実際の吸気管部分２３内の圧力は吸気脈動により大きく変動しており、吸気管内圧力センサ４０によって実測される吸気管内圧力も図２７に実線で示したように大きく変動する。一方、上述した吸気管モデルＭ２０や吸気弁モデルＭ３０では、吸気管部分２３内に生じる圧力脈動を無視して吸気管部分２３内の実際の圧力の１サイクル（本実施形態では、クランク角７２０°）当たりの平均値（以下、「サイクル間平均値」と称す。）を吸気管内圧力Ｐｍとして用いている。また、吸気管内圧力センサ４０の出力にはノイズが含まれてしまうことも多い。

このため、吸気管内圧力センサ４０によって実測された時刻ｔ₀および時刻ｔ₁における吸気管内圧力Ｐｍ（ｔ₀）およびＰｍ（ｔ₁）を直接用いて誤差率αを算出すると、すなわち機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが或るマップ領域に侵入した時およびこのマップ領域から離脱した時に吸気管内圧力センサ４０によって実測された吸気管内圧力Ｐｍｒを直接用いて誤差率αを算出すると、算出された誤差率αの値が不正確なものになってしまう場合がある。

そこで、本実施形態では、上述したように吸気管内圧力センサ４０によって実測された吸気管内圧力Ｐｍを直接用いて誤差率αを算出せずに、吸気管内圧力センサ４０によって実測された吸気管内圧力のサイクル間平均値（図２７中の破線）を用いて誤差率αを算出することとしている。すなわち、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが或るマップ領域に侵入した時を含むサイクルにおいて吸気管内圧力センサ４０によって実測された吸気管内圧力のサイクル間平均値をこのマップ侵入時の吸気管内圧力Ｐｍｒとして用い、このマップ領域から離脱した時を含むサイクルにおいて吸気管内圧力センサ４０によって実測された吸気管内圧力のサイクル間平均値をこのマップ離脱時の吸気管内圧力Ｐｍｒとして用いて、上記式（４６）により誤差率αを算出する。

このように、式（４６）に用いる吸気管内圧力Ｐｍｒとして吸気管内圧力センサ４０によって検出された吸気管内圧力のサイクル間平均値を用いることにより、式（４６）に入力される吸気管内圧力Ｐｍｒが、吸気管モデルＭ２０や吸気弁モデルＭ３０で用いられている吸気管内圧力Ｐｍと同様に、吸気管部分２３内に生じる圧力脈動の影響を排除した値となっており、よって式（４６）によって正確に誤差率αを算出することができるようになる。また、誤差率αを算出するにあたって吸気管内圧力センサ４０の出力に含まれているノイズの影響を抑制することができる。

次に、上記実施形態の変更例について説明する。上記実施形態では、誤差率αを算出した後、この誤差率αを適合値ａおよび適合値ｂに乗算して適合値ａ、ｂの両マップ値を同じ割合で修正している。しかしながら、適合値ａおよび適合値ｂの両マップ値に必ずしも同程度の誤差が生じている場合ばかりではなく、適合値ａのみの誤差が大きい場合や適合値ｂのみの誤差が大きい場合も存在する。したがって、適合値ａおよび適合値ｂのマップ値を正確に修正するためには、適合値ａおよび適合値ｂを同率で修正せずに、適合値ａおよび適合値ｂを別々に修正する必要がある。

そこで、本変更例では、上述したように算出された誤差率αに基づいて適合値ａおよび適合値ｂを別々に修正することとしている。本変更例では、まず、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが或るマップ領域にある場合において、上記式（４６）に基づいて誤差率αを算出すると共に、誤差率αの算出における積分期間、すなわち機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが上記マップ領域にある期間中における吸気管内圧力Ｐｍの平均値Ｐｍａｖｅを算出する。便宜上、このようにして算出された誤差率および吸気管内圧力の平均値をそれぞれα₁、Ｐｍ₁とする。その後、再び機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが上記マップ領域にある場合において、上述したように誤差率αと吸気管内圧力の平均値Ｐｍａｖｅを算出する。便宜上、このようにして算出された誤差率および吸気管内圧力の平均値をそれぞれα₂、Ｐｍ₂とする。

図２８は、機関回転数ＮＥおよび吸気弁開弁時期ＶＶＴが上記マップ領域にある場合における吸気管内圧力Ｐｍと筒内吸入空気流量ｍｃとの関係を示した図である。適合値ａ、ｂのマップ値が修正される前に吸気弁モデルＭ３０によって算出される筒内吸入空気流量ｍｃと吸気管内圧力Ｐｍとの関係は、例えば、図２８に実線で示された関係になっている。

ここで、上述したようにして算出された誤差率α₁および吸気管内圧力Ｐｍ₁を用いて算出されたα₁（ａ・Ｐｍ₁＋ｂ）は吸気管内圧力がＰｍ₁であるときの実際の筒内吸入空気流量を正確に表しており、また、誤差率α₂および吸気管内圧力Ｐｍ₂を用いて算出されたα₂（ａ・Ｐｍ₂＋ｂ）は吸気管内圧力がＰｍ₂であるときの実際の筒内吸入空気流量を正確に表している。このことから、点（Ｐｍ₁、α₁（ａ・Ｐｍ₁＋ｂ））および点（Ｐｍ₂、α₂（ａ・Ｐｍ₂＋ｂ））の二点を通る直線（図２８中の破線）は、吸気管内圧力と実際の筒内吸入空気流量との関係を比較的正確に表していると考えられる。

したがって、適合値ａ、ｂのマップ値を修正することで吸気管内圧力Ｐｍと吸気弁モデルＭ３０によって算出される筒内吸入空気流量ｍｃとの関係を、上記二点を通る直線とすることで、吸気弁モデルＭ３０によって算出される筒内吸入空気流量の推定値ｍｃを実際の筒内吸入空気流量とほぼ等しい値にすることができる。

そこで、本実施形態では、上述した誤差率α₁、α₂および吸気管内圧力Ｐｍ₁、Ｐｍ₂を用いることで、適合値ａ、ｂのマップ値を下記式（５０）、（５１）により修正することとしている。なお、式（５０）、（５１）において、適合値ａ’および適合値ｂ’はそれぞれ修正後のマップ値を示している。

このように適合値ａ、ｂのマップ値を修正することで、修正後の適合値ａ’および適合値ｂ’のマップ値を用いた吸気弁モデルＭ３０によって算出される筒内吸入空気流量ｍｃは実際の筒内空乳空気流量とほぼ等しい値となっている。したがって、本変更例によれば、適合値ａおよび適合値ｂのマップ値を別々に修正することができ、これにより適合値ａおよび適合値ｂのマップ値を正確に修正することができるようになる。

本発明の吸入空気量推定装置を搭載した内燃機関全体を示す図である。図１に示したエアフロメータの概略斜視図である。図２に示したエアフロメータの熱線計量部の拡大斜視図である。エアフロメータの出力電圧とエアフロ通過空気流量との関係を示す図である。エアモデルのブロック図である。スロットル弁開度と流量係数との関係を示す図である。スロットル弁開度と開口断面積との関係を示す図である。関数Φ（Ｐｍ／Ｐａ）を示す図である。スロットルモデルの基本概念を示す図である。吸気管モデルの基本概念を示す図である。吸気弁モデルの基本概念を示す図である。筒内充填ガス量および筒内吸気ガス量の定義に関する図である。吸気管モデルに生じた誤差の修正方法についての基本概念を示す図である。吸気管内圧力及びスロットル弁開度のタイムチャートである。吸気管部分の容積の誤差率と平均比率との関係を表すマップを示す図である。吸気管モデルに用いられる吸気管部分の容積の値の修正制御のフローチャートである。スロットルモデルの修正操作を示す図である。スロットル通過空気流量と完全放熱量との関係を示す図である。遅れ放熱量の和とエアフロメータの予想出力値との関係を示す図である。スロットル弁開度と空気流量に関する各種パラメータの値とのタイムチャートである。スロットル弁開度に基づいて定められた流量係数のマップを示すである。スロットル弁開度と空気流量に関する各種パラメータの値とのタイムチャートである。図２２の一部を拡大した図である。初期値応答とゼロ応答とを説明するための図である。機関回転数および吸気弁開弁時期に基づいて定められた値ａ、ｂのマップを示す図である。吸気管内圧力、機関回転数、吸気弁開弁時期のタイムチャートである。吸気管内圧力のタイムチャートである。吸気管内圧力と筒内吸入空気流量との関係を示した図である。

符号の説明

１機関本体
５燃焼室
６吸気弁
７吸気ポート
８排気弁
１１燃料噴射弁
１４吸気枝管
１５サージタンク
１６吸気管
１９スロットル弁
２３吸気管部分

Claims

スロットル通過空気流量に基づいて吸気管モデル計算式を用いて吸気管内圧力を算出する吸気管モデルを具備し、該吸気管モデルによって算出された吸気管内圧力の推定値に基づいて筒内充填空気量を推定する内燃機関の吸入空気量推定装置において、
実際の吸気管内圧力を検出する吸気管内圧力検出手段をさらに具備し、
上記吸気管モデルの吸気管モデル計算式においてパラメータとして吸気管容積が用いられ、
上記吸気管モデルによって吸気管内圧力の推定値を算出すると共に上記吸気管内圧力検出手段によって吸気管内圧力の実測値を検出し、これら吸気管内圧力の推定値と実測値とに基づいて上記吸気管モデルで用いられる吸気管容積の値を修正するようにした、内燃機関の吸入空気量推定装置。
上記吸気管モデルで用いられる吸気管容積の値の修正は、吸気管内圧力の変化率が基準変化率以上であるときの吸気管内圧力の推定値と実測値とに基づいて行われる、請求項１に記載の内燃機関の吸入空気量推定装置。
上記吸気管モデルで用いられる吸気管容積の値は、吸気管内圧力の変化率が基準変化率以上であるときの吸気管内圧力の推定値と実測値とが等しくなるように修正される、請求項２に記載の内燃機関の吸入空気量推定装置。
上記吸気管モデルは、スロットル通過空気流量に加えて筒内吸入空気流量に基づいて吸気管内圧力を算出し、上記スロットル通過空気流量は、スロットル弁開度に基づいてスロットルモデル計算式を用いてスロットル通過空気流量を算出するスロットルモデルによって算出され、上記筒内吸入空気流量は、吸気管内圧力に基づいて吸気弁モデル計算式を用いて筒内吸入空気流量を算出する吸気弁モデルによって算出され、
上記吸気管内圧力検出手段によって検出された吸気管内圧力に基づいて上記スロットルモデル計算式又は吸気弁モデル計算式を修正するスロットル・吸気弁モデル修正手段をさらに具備し、
上記吸気管モデルの吸気管容積の値の修正は、スロットル・吸気弁モデル修正手段によってスロットルモデル計算式又は吸気弁モデル計算式を修正した後に行われる、請求項１〜３のいずれか１項に記載の内燃機関の吸入空気量推定装置。