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