吸気弁閉弁時期が遅角されていて吸気がスロットル弁を亜音速で通過しているときに、スロットル弁の開度を推定するスロットル弁開度推定方法であって、前記吸気弁閉弁時期を進角して吸気が前記スロットル弁を音速で通過するようにする第一段階と、エアフローメータにより吸入空気量を検出する第二段階と、前記エアフローメータにより検出された前記吸入空気量に基づき前記スロットル弁の開度を推定する第三段階とを有することを特徴とするスロットル弁開度推定方法。前述のようにしてスロットル弁の開度を推定するときに、吸気弁の閉弁時期が吸気下死点より大きく遅角されていて、圧縮行程の吸気弁開弁中に比較的多量の吸気が気筒内からスロットル弁下流側の吸気管へ戻されて吸気管内の圧力が比較的高くなっていると、スロットル弁を通過する吸気量はスロットル弁の開度だけでなく、吸気管内の圧力によっても変化することとなり、正確なスロットル弁の開度を推定することが困難となる。従って、本発明の目的は、機関定常時においてエアフローメータにより検出される吸入空気量に基づき、正確にスロットル弁の開度を推定することができるスロットル弁開度推定方法を提供することである。本発明による請求項1に記載のスロットル弁開度推定方法は、吸気弁閉弁時期が遅角されていて吸気がスロットル弁を亜音速で通過しているときに、スロットル弁の開度を推定するスロットル弁開度推定方法であって、前記吸気弁閉弁時期を進角して吸気が前記スロットル弁を音速で通過するようにする第一段階と、エアフローメータにより吸入空気量を検出する第二段階と、前記エアフローメータにより検出された前記吸入空気量に基づき前記スロットル弁の開度を推定する第三段階とを有することを特徴とする。吸気弁閉弁時期が遅角されていて圧縮行程の吸気弁開弁中に比較的多量の吸気が気筒内からスロットル弁下流側の吸気管へ戻されると、吸気管内の圧力が比較的高くなって、吸気はスロットル弁を亜音速で通過することとなる。このときには、スロットル弁を通過する吸気量はスロットル弁の開度だけでなく、吸気管内の圧力によっても変化することとなるために、正確なスロットル弁の開度を推定することが困難である。本発明による請求項1に記載のスロットル弁開度推定方法によれば、第一段階として、吸気弁閉弁時期を進角して圧縮行程の吸気弁開弁中に気筒内から吸気管へ戻される吸気を少なくして吸気管内の圧力を低下させ、吸気がスロットル弁を音速で通過するようにし、第二段階として、エアフローメータにより吸入空気量を検出し、第三段階として、エアフローメータにより検出された吸入空気量に基づきスロットル弁の開度を推定するようになっている。それにより、第二段階においてエアフローメータにより検出される吸入空気量は、第一段階によって吸気がスロットル弁を音速で通過しているときの吸気量となり、スロットル弁の開度だけにより変化する値となるために、第三段階において、エアフローメータにより検出された吸入空気量に基づき正確にスロットル弁の開度を推定することが可能となる。図1は、本発明によるスロットル弁開度推定方法が適用される内燃機関を示す概略図である。同図において、1は機関本体であり、2は各気筒共通のサージタンクである。3はサージタンク2と各気筒とを連通する吸気管シリンダヘッドに形成された吸気ポートを含むであり、4はサージタンク2の上流側の吸気通路である。各吸気管3には燃料噴射弁5が配置され、吸気通路4におけるサージタンク2の直上流側にはスロットル弁6が配置されている。スロットル弁6は、アクセルペダルに連動するものではなく、ステップモータ等のアクチュエータによって駆動されるものであり、自由に開度設定可能となっている。7は吸気通路4のスロットル弁6より上流側において吸気流量を検出するエアフローメータである。機関本体1において、8は吸気弁であり、9は排気弁であり、10はピストンである。11はスロットル弁6の開度を検出するためのスロットル弁開度センサである。12はエアクリーナである。吸気弁8は可変バルブタイミング機構図示せずにより少なくとも閉弁時期が可変とされている。可変バルブタイミング機構は、例えば、吸気弁用カムの位相を変化させるものでも良いし、吸気弁を任意の時期に開閉させることができるアクチュエータとしても良い。内燃機関1における燃焼空燃比を、例えば、理論空燃比等の所望空燃比にするためには、機関過渡時を含めて気筒内へ流入した吸入空気量を正確に推定することが必要とされる。エアフローメータ7は、機関定常時においては、比較的正確に吸入空気量を測定することができる。しかしながら、機関過渡時においては、急激に変化する吸入空気量に対してエアフローメータ7の出力が直ぐに応答せず、正確な吸入空気量の測定は不可能である。それにより、機関過渡時においても、正確な吸入空気量を把握することを可能とするために、機関吸気系をモデル化して電子制御装置図示せずによって吸入空気量を推定するようにしている。先ずは、スロットル弁6をモデル化する。吸気がスロットル弁6を通過する際のエネルギ保存則、運動量保存則、及び、状態方程式を使用して、今回のスロットル弁通過空気量mtgsecが、次式1によって表される。以下の式を含めて、スロットル弁通過空気量等の変数の添え字iは今回を示し、i-1は前回を示している。mt=μA12ΦPa1ここで、μは流量係数であり、Aはスロットル弁6の開口面積m3である。もちろん、機関吸気系にアイドルスピードコントロールバルブISC弁が設けられている時には、Aには、ISC弁の開口面積が加えられる。流量係数及びスロットル弁の開口面積は、それぞれがスロットル弁開度TA度の関数となっており、図2及び3には、それぞれのスロットル弁開度TAに対するマップが図示されている。Rは気体定数であり、Taはスロットル弁上流側の吸気温度Kであり、Paはスロットル弁上流側の吸気通路圧力kPa、すなわち、大気圧必要に応じて変数として大気圧センサにより測定するようにしても良いであり、Pmはスロットル弁下流側の吸気管圧力kPaである。また、関数ΦPaは、比熱比kを使用して次式2及び2’によって表されるものであり、図4にはPm/Paに対するマップが図示されている。PmPa<=1の場合ΦPa=122PmPa>1の場合ΦPa={2kPaPmPaPa}122’次いで、吸気弁をモデル化する。気筒内へ供給される吸入空気量mcgsecは、吸気管圧力Pmに基づきほぼ線形に変化するものであるために、次式3によって表すことができる。mc=TaTmaPm3ここで、Tmはスロットル弁下流側の吸気温度Kであり、aは、可変バルブタイミング機構により吸気弁8の閉弁時期が遅角されるほど、気筒内から吸気系へ戻される吸気量が増加して結果的に吸気量が減少し、また、機関回転数Nが高いほど、単位時間当たりの吸気量が減少するために、小さくなる係数であり、実機を使用する適合試験等によって吸気弁8の閉弁時期と機関回転数Nとに基づきマップ化しておくことが好ましい。次いで、吸気管をモデル化する。吸気管内に存在する吸気の質量保存則、エネルギ保存則、及び、状態方程式を使用して、吸気管圧力Pmとスロットル弁下流側の吸気温度Tmとの比における時間変化率は次式4によって表され、また、吸気管圧力Pmの時間変化率は次式5によって表される。ここで、Vは吸気管の容積m3であり、具体的には、サージタンク2と吸気管3との合計容積である。ddt=RV4dPmdt=kRV5式4及び式5は離散化され、それぞれ、次式6及び7が得られ、式7によって今回の吸気管圧力Pmが得られれば、式6によって今回の吸気管内の吸気温度Tmを得ることができる。式6及び7において、離散時間Δtは、現在の吸入空気量mcを算出する実行間隔とされ、例えば8msである。PmTm=PmTmΔtRVmc6Pm=PmΔtkRVTamcTm7こうして、吸気量を推定する際には、先ず、式7を使用して吸気管圧力Pmが算出される。式7は、前回の吸気管圧力Pmと、前回のスロットル弁通過空気量mtと、前回のスロットル弁上流側の吸気温度Taと、前回の吸入空気量mcと、前回の吸気管内の吸気温度Tmとに基づき、今回の吸気管圧力Pmを算出するようになっている。これらの初期値として、Pmには大気圧Paが、TaとTmにはスロットル弁上流側の吸気温度Taがそれぞれ実測されて使用され、mtには、これらのPm及びTaを使用して式1から算出された値が使用され、また、mcには、これらのPm、Ta及びTmを使用して式3により算出された値が使用される。次いで、式6を使用して今回の吸気管内の吸気温度Tmが算出される。次いで、式1を使用して今回のスロットル弁通過空気量mtが算出される。この式1を使用するスロットル弁通過空気量mtの算出において、現在の流量係数μ及び現在のスロットル弁6の開口面積Aを決定するための現在のスロットル弁開度TAは、スロットル弁開度センサ11の出力が使用される。次いで、式3を使用して今回の吸入空気量mcが算出される。その後は、今回の吸気管圧力Pmは前回の吸気管圧力Pmとされ、今回の吸気管内の吸気温度Tmは前回の吸気管内の吸気温度Tmとされ、今回のスロットル弁通過空気量mtは前回のスロットル弁通過空気量mtとされ、今回の吸入空気量mcは前回の吸入空気量mcとされる。こうして、吸入量mcは、機関始動完了と同時に逐次算出される吸気管圧力Pmに基づき、逐次推定されることとなる。ところで、スロットル弁開度センサ11が故障したときには、現在のスロットル弁6の開度を推定することが必要となる。本実施例では、スロットル弁開度センサ11が故障したときには、電子制御装置により図5に示すフローチャートに従ってスロットル弁6の開度を推定するようになっている。先ず、ステップ101において、スロットル弁開度センサ11が故障したか否かが判断される。この判断が否定されるときには、スロットル弁6の開度の推定は必要なく、そのまま終了する。しかしながら、例えば、スロットル弁開度センサ11が出力信号を発しないときや、アクチュエータを作動してスロットル弁6の開度を変更しているにも係わらずにスロットル弁開度センサ11の出力信号が変化しないときには、スロットル弁開度センサ11は故障していると判断され、ステップ101が肯定されてステップ102へ進む。ステップ102では、式7により算出されるスロットル弁下流側の現在の吸気管圧力と大気圧との比PmPaが1より大きいか否かが判断される。ここでkは前述したように比熱比であり、ステップ102の判断が否定されるときには、吸気はスロットル弁6を音速で通過しており、式2及び図4に示すように、関数Φの値がスロットル弁下流側の現在の吸気管圧力と大気圧との比PmPaに対して一定値となる。それにより、式1で示すスロットル弁通過空気量mtは、スロットル弁6の開口面積A図2及び3に示すように、流量係数μと共にスロットル弁開度TAから一義的に定まる値であり、また、スロットル弁上流側の吸気温度Taは正確に測定可能であるだけにより変化する値となる。こうして、ステップ104においてエアフローメータ7により現在の吸入吸気量Qを検出し、次いで、ステップ105において、吸入空気量Qが一定となる機関定常時には、エアフローメータ7により検出される吸入空気量Qとスロットル弁通過空気量mtとは等しくなることを利用して、吸入空気量Qに基づき式1を使用してスロットル弁7の開口面積Aを逆算し、こうして算出された開口面積Aからスロットル弁6の開度TAを正確に推定する。一方、現在の吸気管圧力と大気圧との比PmPaが1より大きいときには、吸気はスロットル弁6を亜音速で通過しており、式2’及び図4に示すように、関数Φの値は吸気管圧力Pmによって変化する。それにより、式1で示すスロットル弁通過空気量mtは、スロットル弁6の開口面積Aだけでなく、吸気管圧力Pmによっても変化する値となり、式1を使用してスロットル弁7の開口面積Aを正確に逆算することが困難であり、こうして算出された開口面積Aからスロットル弁6の開度TAを正確に推定することも難しい。そのために、本フローチャートでは、ステップ102の判断が肯定されるときには、ステップ103において、スロットル弁6の開度を変化させることなく、可変バルブタイミング機構によって吸気弁8の閉弁時期を進角するようになっている。それにより、圧縮行程の吸気弁開弁中に、気筒内から吸気管へ戻される吸気を減少させることにより、吸気管圧力Pmを低下させる。このような吸気弁8の閉弁時期の進角によってエアフローメータ7により検出される吸入空気量Qは増加することとなる。こうして、吸気管圧力Pmを十分に低下させれば、吸入空気量Qが吸気管圧力の低下により増加せずに一定となり、ステップ102の判断が否定され、吸気がスロットル弁6を音速で通過するようになり、ステップ104において吸入空気量Qが一定となる機関定常時にエアフローメータ7により現在の吸入空気量Qを検出し、次いで、ステップ105において現在の吸入空気量Qに基づき正確にスロットル弁6の開度TAを推定することが可能となる。ステップ102の判断は、吸気弁8の閉弁時期を遅角してもエアフローメータ7により検出される吸入空気量Qが増加しなくなったときに否定されるようにしても良い。もし、スロットル弁開度センサ11の故障時にスロットル弁6の開度が比較的大きいと、吸気弁8の閉弁時期を進角しても吸気管圧力Pmは十分に低くならず、ステップ102の判断が否定されないこともある。この場合には、スロットル弁6の開度を正確に推定することはできない。