JP4872383B2 - エンジンの点火時期制御方法及びエンジンの点火時期制御装置 - Google Patents
エンジンの点火時期制御方法及びエンジンの点火時期制御装置 Download PDFInfo
- Publication number
- JP4872383B2 JP4872383B2 JP2006059318A JP2006059318A JP4872383B2 JP 4872383 B2 JP4872383 B2 JP 4872383B2 JP 2006059318 A JP2006059318 A JP 2006059318A JP 2006059318 A JP2006059318 A JP 2006059318A JP 4872383 B2 JP4872383 B2 JP 4872383B2
- Authority
- JP
- Japan
- Prior art keywords
- combustion
- combustion chamber
- mass ratio
- ignition timing
- ignition
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Images
Classifications
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T10/00—Road transport of goods or passengers
- Y02T10/10—Internal combustion engine [ICE] based vehicles
- Y02T10/40—Engine management systems
Landscapes
- Electrical Control Of Ignition Timing (AREA)
- Combined Controls Of Internal Combustion Engines (AREA)
Description
本発明は、エンジン（内燃機関）の点火時期制御方法及びエンジンの点火時期制御装置に関する。
ＭＢＴ（最大トルクの得られる最小進角値）で混合気に点火した場合に、混合気の燃焼圧力が最大値Ｐｍａｘとなるクランク角θｐｍａｘを基準クランク角とすると、基準クランク角は燃焼方式によらずほぼ一定である。また、燃焼室内における燃焼解析によれば、燃焼室に供給された燃料に対する燃焼質量の比率を表す燃焼質量割合は点火時に０％であり、完全燃焼によって１００％に達する。そして、基準クランク角における燃焼質量割合は一定で約６０％であるとして、そのときの運転条件（エンジンの負荷と回転速度）より燃焼速度を求め、この燃焼速度に基づいて燃焼期間ＢＵＲＮを算出し、基準クランク角よりこの燃焼期間ＢＵＲＮと着火遅れ時間相当角ＩＧＮＤＥＡＤとの合計のクランク角区間だけ進角側のクランク角位置をＭＢＴの得られる基本点火時期として算出する点火時期制御方法を提案している（特許文献１参照）。
特開２００４−３３２６４７号公報
ところで、特許文献１でいう着火遅れ時間相当角ＩＧＮＤＥＡＤは、点火時期から火炎の形成が開始されるタイミング（火炎形成開始時期）までの時間である点火無駄時間ＤＥＡＤＴＩＭＥをクランク角区間へと換算した値である。これを燃焼質量割合でいうと、点火タイミングを始期とし、燃焼質量割合が０％より上昇を開始するタイミングを終期とする時間をクランク角区間へと換算した値である。従って、火炎形成開始時期は燃焼質量割合が０％より上昇を開始するタイミングに相当する。
この場合に、火炎形成開始時期は、一般的に燃焼室内圧力に基づいて計測している。すなわち、火炎形成開始前後で燃焼室内圧力が変化するはずであるから、その燃焼室内圧力の変化タイミングをもって火炎形成開始時期としている。
しかしながら、実際には火炎形成開始時期を精度良く計測することができていないのが現状である。というのも、可視化した燃焼室内を目視すると火炎が形成されているのに、火炎形成開始直後の燃焼室内圧力の上昇（変化）が少なすぎるためか、燃焼解析装置では燃焼質量割合が０％のままである。つまり、燃焼解析装置によれば、燃焼が始まっていないように見えるのである。このように、現状の燃焼解析装置では火炎形成開始時期そのものを精度よく計測することが難しく、実際よりも遅れて火炎形成開始時期が計測されてしまう。
そして、火炎形成開始時期の計測に誤差が生じているということは点火無駄時間ＤＥＡＤＴＩＭＥの適合に誤差が生じることを意味する。点火無駄時間ＤＥＡＤＴＩＭＥの適合に誤差が生じると、最終的にＭＢＴの得られる点火時期の推定精度が悪化してしまう。
そこで本発明は、着火遅れ時間の算出方法及び算出手段を改良した点火時期制御方法及び点火時期制御装置を提供することを目的とする。
本発明は、点火の瞬間から火炎が形成されるとみなし、この点火直後の火炎の速度を、層流燃焼速度と燃焼室内ガスの運動エネルギーと燃焼室内ガスの熱膨張率に基づいて算出し、初期燃焼時の燃焼質量割合を基準燃焼質量割合として、この基準燃焼質量割合になったときの燃焼室内容積に基づいて基準燃焼質量割合になったときの火炎半径を算出し、点火から燃焼質量割合が前記基準燃焼質量割合に上昇するまでの時間を着火遅れ時間とみなして、この着火遅れ時間を前記点火直後の火炎速度と前記基準燃焼質量割合になったときの火炎半径とに基づいて算出し、燃焼質量割合が前記基準燃焼質量割合から基準クランク角での燃焼質量割合になるまでの時間である燃焼時間を算出し、これら着火遅れ時間と燃焼時間とを加算した値をクランク角区間に変換し、前記基準クランク角よりこの変換されたクランク角区間だけ進角側の値をＭＢＴの得られる基本点火時期を算出するように構成する。
火炎形成開始直後に燃焼室内圧力の上昇（変化）が十分となれば、現状の燃焼解析装置でも燃焼質量割合の計測結果が精度の良いものとなることに着目し、燃焼室内圧力の上昇（変化）が十分となる初期燃焼時の燃焼質量割合を基準燃焼質量割合として定める。この基準燃焼質量割合として例えば２〜１０％を設定する。２〜１０％の範囲でよいのは、似たような初期燃焼の範囲内にあるためである。こうして、点火時期から燃焼質量割合がこの基準燃焼質量割合に上昇するまでの時間を燃焼解析装置により計測するのであれば、この時間計測は現状の燃焼解析装置によっても精度の良いものとなる。そして、点火時期から燃焼質量割合が基準燃焼質量割合に上昇するまでのこの時間に基づいてＭＢＴの得られる基本点火時期を算出したとき、十分な精度が得られることを確認している。
このように、本発明によれば、点火から燃焼質量割合が０％より上昇を開始するタイミングまでの時間（コンベンショナルエンジンに適用している先行点火時期制御方法における点火無駄時間ＤＥＡＤＴＩＭＥ）ではなく、点火から燃焼質量割合が基準燃焼質量割合に上昇するまでの時間を着火遅れ時間とみなして算出（推定）し、この算出した着火遅れ時間に基づいてＭＢＴの得られる基本点火時期を算出するので、点火から燃焼質量割合が０％より上昇を開始するタイミングまでの時間に基づいてＭＢＴの得られる基本点火時期を算出する場合よりも、ＭＢＴの得られる基本点火時期（ＭＢＴＣＡＬ）の予測精度を向上させることができる。
以下、図面に基づき本発明の実施形態について説明する。
図１は、エンジンの点火時期制御方法の実施に直接使用するエンジンの点火時期制御装置の概略構成を示している。
空気は吸気コレクタ２に蓄えられた後、吸気マニホールド３を介して各気筒の燃焼室５に導入される。燃料は各気筒の吸気ポート４に配置された燃料インジェクタ２１より噴射供給される。空気中に噴射された燃料は気化しつつ空気と混合してガス（混合気）を作り、燃焼室５に流入する。この混合気は吸気弁１５が閉じることで燃焼室５内に閉じこめられ、ピストン６の上昇によって圧縮される。
この圧縮混合気に対して高圧火花により点火を行うため、パワートランジスタ内蔵の点火コイルを各気筒に配した電子配電システムの点火装置１１を備える。すなわち、点火装置１１は、バッテリからの電気エネルギーを蓄える点火コイル１３と、点火コイル１３の一次側への通電、遮断を行うパワートランジスタ（図示しない）と、燃焼室５の天井に設けられ点火コイル１３の一次電流の遮断によって点火コイル１３の二次側に発生する高電圧を受けて、火花放電を行う点火プラグ１４とからなっている。
圧縮上死点より少し手前で点火プラグ１４により火花が飛ばされ圧縮混合気に着火されると、火炎が広がりやがて爆発的に燃焼し、この燃焼によるガス圧がピストン６を押し下げる仕事を行う。この仕事はクランクシャフト７の回転力として取り出される。燃焼後のガス（排気）は排気弁１６が開いたときに排気通路８へと排出される。
排気通路８には三元触媒９、１０を備える。三元触媒９、１０は排気の空燃比が理論空燃比を中心とした狭い範囲（ウインドウ）にあるとき、排気中に含まれるＨＣ、ＣＯ、ＮＯｘといった有害三成分を同時に効率よく除去できる。空燃比は吸入空気量と燃料量の比であるので、エンジンの１サイクル（４サイクルエンジンではクランク角で７２０°区間）当たりに燃焼室５に導入される吸入空気量と、燃料インジェクタ２１からの燃料噴射量との比が理論空燃比となるように、エンジンコントローラ３１ではエアフローメータ３２からの吸入空気流量の信号とクランク角センサ（３３、３４）からの信号に基づいて燃料インジェクタ２１からの燃料噴射量を定めると共に、三元触媒９の上流に設けたＯ2センサ３５からの信号に基づいて空燃比をフィードバック制御している。
吸気コレクタ２の上流には絞り弁２３がスロットルモータ２４により駆動される、いわゆる電子制御スロットル２２を備える。運転者が要求するトルクはアクセルペダル４１の踏み込み量（アクセル開度）に現れるので、エンジンコントローラ３１ではアクセルセンサ４２からの信号に基づいて目標トルクを定め、この目標トルクを実現するための目標空気量を定め、この目標空気量が得られるようにスロットルモータ２４を介して絞り弁２３の開度を制御する。
吸気弁１５、排気弁１６は、クランクシャフト７を動力源として、各々吸気側カムシャフト２５及び排気側カムシャフト２６に設けられたカムの動作により開閉駆動される。吸気側には、吸気弁１５のバルブリフト量及び作動角を連続的に可変制御する多節リンク状の機構で構成される可変バルブ機構（以下、「ＶＥＬ機構」という。）２８を備える。このＶＥＬ機構２８には吸気弁１５のバルブリフト量及び作動角を検出する作動角センサ（図示しない）が併設されている。
同じく吸気側には、クランクシャフト７と吸気側カムシャフト２５との回転位相差を連続的に可変制御して、吸気弁１５のバルブタイミングを進遅角する可変バルブタイミング機構（以下、「ＶＴＣ機構」という。）２７を備える。また、吸気側カムシャフト２５の他端には吸気側カムシャフト２５の回転位置を検出するためのカム角センサ３４が併設されている。
これらＶＥＬ機構２８及びＶＴＣ機構２７（可変動弁装置）の具体的な構成は特開２００３−３８７２号公報により公知であるので、その詳しい説明は省略する。
ＶＥＬ機構２８、ＶＴＣ機構２７の各アクチュエータに指令して、吸気弁１５のリフト特性（バルブタイミング（開閉時期）や吸気弁１５のバルブリフト量）を変えると燃焼室５に残留する不活性ガスの量が変化する。燃焼室５内の不活性ガスの量が増えるほどポンピングロスが減って燃費がよくなるので、運転条件によりどのくらいの不活性ガスが燃焼室５内に残留したらよいかを目標吸気弁閉時期や目標バルブリフト量にして予め定めており、エンジンコントローラ３１ではそのときの運転条件（エンジンの負荷と回転速度）より目標吸気弁閉時期と目標バルブリフト量とを定め、それら目標値が得られるようにＶＴＣ機構２７及びＶＥＬ機構２８の各アクチュエータを介して吸気弁１５の閉時期とバルブリフト量とを制御する。
吸気温度センサ４３からの吸気温度の信号、吸気圧力センサ４４からの吸気圧力の信号、排気温度センサ４５からの排気温度の信号も入力されるエンジンコントローラ３１では、パワートランジスタ１３を介して点火プラグ１４の一次側電流の遮断時期である点火時期を制御する。
さて、ＭＢＴ（最大トルクの得られる最小進角値）で混合気に点火した場合に、混合気の燃焼圧力が最大値Ｐｍａｘとなるクランク角θｐｍａｘを基準クランク角とすると、基準クランク角は燃焼方式によらずほぼ一定である。また、燃焼室内における燃焼解析によれば、燃焼室５に供給された燃料に対する燃焼質量の比率を表す燃焼質量割合は点火時に０％であり、完全燃焼によって１００％に達する。そして、基準クランク角における燃焼質量割合は一定で約６０％であるとして、そのときの運転条件（エンジンの負荷と回転速度）より燃焼速度を求め、この燃焼速度に基づいて燃焼期間ＢＵＲＮを算出し、基準クランク角よりこの燃焼期間ＢＵＲＮと着火遅れ時間相当角ＩＧＮＤＥＡＤとの合計のクランク角区間だけ進角側のクランク角位置を基本点火時期ＭＢＴＣＡＬとして算出する点火時期制御方法を提案している（特開２００４−３３２６４７号公報参照）。
燃焼解析に基づくこうした点火時期制御方法（この点火時期制御方法を、以下「先行点火時期制御方法」という。）は吸気弁１５のバルブリフト量及び吸気弁１５の開閉タイミングが変化しないエンジン（コンベンショナルエンジン）を対象に本出願人が開発してきた経過があるので、現在でもコンベンショナルエンジンに適用する限り先行点火時期制御方法でなんら問題ないのであるが、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンに対しても、この先行点火時期制御方法を適用したとき、上記の基本点火時期ＭＢＴＣＡＬがＭＢＴの得られる点火時期と合わないことが判明している。その理由は、ＶＥＬ機構２８及びＶＴＣ機構２７が非作動状態であるときに点火時期制御に用いる各種の値を適合していれば、当然ながら、ＶＥＬ機構２８及びＶＴＣ機構２７が非作動状態であるときに基本点火時期ＭＢＴＣＡＬはＭＢＴの得られる点火時期と一致する。
この状態つまり同じ運転条件においてＶＥＬ機構２８やＶＴＣ機構２７を作動させてリフト特性（吸気弁１５のバルブリフト量や吸気弁閉時期）を変化させると、燃焼室５内のガス流動であるタンブルやスワールの各強度が変化し、これに伴って燃焼室内ガスの乱流状態での燃焼速度である乱流燃焼速度が変化する。すると、ＭＢＴの得られる点火時期も変化する。しかしながら、ＶＥＬ機構２８やＶＴＣ機構２７を作動させたからといって運転条件が同じであるため、基本点火時期ＭＢＴＣＡＬは変化しない。従って、ＶＥＬ機構２８やＶＴＣ機構２７を作動させた途端に基本点火時期ＭＢＴＣＡＬがＭＢＴの得られる点火時期からずれて燃費が悪くなる。
例えば、ＶＥＬ機構２８の作動で吸気弁１５のバルブリフト量が大きい状態から小さい状態に切換えられたとすると、切換前よりタンブルやスワールが形成されにくくなり乱流燃焼速度が遅くなって燃焼期間が長引く。従って、同じ運転条件で考えると、バルブリフト量が小さい状態であるときには、バルブリフト量が大きい状態のときよりＭＢＴの得られる点火時期が進角側に移動するのであるが、バルブリフト量が小さい状態のときにも、バルブリフト量が大きい状態のときと同じ燃焼期間ＢＵＲＮを算出したのでは、実際より短い燃焼期間を算出してしまうことになり、基本点火時期ＭＢＴＣＡＬがＭＢＴの得られる点火時期よりも遅すぎることとなる。
ＶＴＣ機構２７の作動で吸気弁閉時期が遅い状態から早い状態へと進角される場合も同様である。すなわち、吸気弁閉時期が遅い状態から早い状態へと進角されると、進角される前よりタンブルやスワールが形成されにくいために乱流燃焼速度が遅くなって燃焼期間が長くなる。従って、同じ運転条件で考えると、吸気弁閉時期が早い状態であるときには、吸気弁閉時期が遅い状態のときよりＭＢＴの得られる点火時期が進角側に移動するのであるが、吸気弁閉時期が早いときにも吸気弁閉時期が遅い状態のときと同じ燃焼期間ＢＵＲＮを算出したのでは、実際より短い燃焼期間を算出してしまうことになり、基本点火時期ＭＢＴＣＡＬがＭＢＴの得られる点火時期より遅すぎてしまう。
そこで本実施形態では、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンを対象としていても、物理モデルに従った制御構造をもち、可能な限り適合の不要な制御とするため、コンベンショナルエンジンに適用している先行点火時期制御方法に対して〈１〉燃焼速度の算出方法、〈２〉基準クランク角の設定方法、〈３〉着火遅れ時間の算出方法を変更する。
ここで、上記〈１〉だけでなく、〈２〉と〈３〉をも追加しているのは、制御方法を再検討した結果である。まず、上記〈２〉をも追加しているのは次の理由による。すなわち、コンベンショナルエンジンに適用している先行点火時期制御方法に対してこれまで正しいとされてきた、基準クランク角と燃焼質量割合が６０％のときのクランク角位置との関係を、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンについて実験してみると、基準クランク角と燃焼質量割合が６０％のときのクランク角位置との間に大きなずれ（クランク角差）があり、そのずれがそのまま点火時期算出値（ＭＢＴＣＡＬ）の推定誤差になってしまうため、基準クランク角位置の設定方法を改める必要が生じたためである。
上記〈３〉をも追加しているのは次の理由による。〈３〉の着火遅れ時間とは、点火タイミングを起点として燃焼質量割合が０％である間の時間のことであるが、この着火遅れ時間の推定誤差が大きい。その理由は、燃焼解析装置によって燃焼質量割合が０％より上昇を開始するタイミングを計測したときの計測誤差が大きいために、この計測した値を正にして先行点火時期制御方法における点火無駄時間ＤＥＡＤＴＩＭＥの適合を行うと、後述する基本点火時期ＭＢＴＡＤＶの推定誤差にこの計測誤差が含まれてしまうためである。従って、燃焼解析装置による計測誤差が入ってこないように、着火遅れ時間の算出方法を改めるようにしたものである。
以下、項を分けて説明する。
〈１〉燃焼速度の算出方法
コンベンショナルエンジンに適用している先行点火時期制御方法では、燃焼速度Ｓｂを燃焼室内ガスの層流状態での燃焼速度である層流燃焼速度ＳＬと、燃焼室内ガスの乱流状態での燃焼速度である乱流燃焼速度ＳＴの和、つまり次の（補１）式により燃焼速度Ｓｂを求めている。
〈１〉燃焼速度の算出方法
コンベンショナルエンジンに適用している先行点火時期制御方法では、燃焼速度Ｓｂを燃焼室内ガスの層流状態での燃焼速度である層流燃焼速度ＳＬと、燃焼室内ガスの乱流状態での燃焼速度である乱流燃焼速度ＳＴの和、つまり次の（補１）式により燃焼速度Ｓｂを求めている。
Ｓｂ＝ＳＬ＋ＳＴ …（補１）
（補１）式の層流燃焼速度ＳＬ、乱流燃焼速度ＳＴはそれぞれ次の（補２）式、（補３）式により与えている。
（補１）式の層流燃焼速度ＳＬ、乱流燃焼速度ＳＴはそれぞれ次の（補２）式、（補３）式により与えている。
ＳＬ＝ＳＬ０×（Ｔ／２９８）＾ｎ×（ｐ／１０１．３２５）＾ｄ
×（１−ｂ×ＭＲＥＳＲ＾ｋ） …（補２）
ＳＴ＝Ｓｔ・Ｎｅ …（補３）
ただし、ＳＬ０ ：標準状態での層流燃焼速度［ｍ／ｓｅｃ］、
Ｔ ：未燃ガス温度［Ｋ］、
ｐ ：燃焼室内圧力［ｋＰａ］、
ＭＲＥＳＲ ：内部不活性ガス率［％］、
Ｓｔ ：乱流燃焼速度係数、
ｎ、ｄ、ｂ、ｋ：係数、
（補３）式の乱流燃焼速度ＳＴは、燃焼室５内におけるガスの乱れ強さによって決まる。また、（補３）式のように、乱流燃焼速度Ｓｔをエンジン回転速度Ｎｅの関数としたのは、燃焼室５内におけるガスの乱れ強さは吸気弁１５のバルブリフト量、開閉タイミングが同じ場合、ピストンスピード（エンジン回転速度）に比例するため、乱流燃焼速度Ｓｔの推定にはエンジン回転速度Ｎｅの感度を持たせれば十分であると考えたものである。
×（１−ｂ×ＭＲＥＳＲ＾ｋ） …（補２）
ＳＴ＝Ｓｔ・Ｎｅ …（補３）
ただし、ＳＬ０ ：標準状態での層流燃焼速度［ｍ／ｓｅｃ］、
Ｔ ：未燃ガス温度［Ｋ］、
ｐ ：燃焼室内圧力［ｋＰａ］、
ＭＲＥＳＲ ：内部不活性ガス率［％］、
Ｓｔ ：乱流燃焼速度係数、
ｎ、ｄ、ｂ、ｋ：係数、
（補３）式の乱流燃焼速度ＳＴは、燃焼室５内におけるガスの乱れ強さによって決まる。また、（補３）式のように、乱流燃焼速度Ｓｔをエンジン回転速度Ｎｅの関数としたのは、燃焼室５内におけるガスの乱れ強さは吸気弁１５のバルブリフト量、開閉タイミングが同じ場合、ピストンスピード（エンジン回転速度）に比例するため、乱流燃焼速度Ｓｔの推定にはエンジン回転速度Ｎｅの感度を持たせれば十分であると考えたものである。
しかしながら、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンでは、コンベンショナルエンジンと相違して、吸気弁１５のバルブリフト量や開閉タイミングが運転条件によって大きく変わり、それらバルブリフト量や開閉タイミングの変化による燃焼室５内におけるガスの乱れ強さの変化を無視できないので、本実施形態では、上記（補３）式に代え、次のようにして乱流燃焼速度を算出する。すなわち、ＶＥＬ機構２８の作動に伴う吸気弁１５のバルブリフト量の変化と、ＶＴＣ機構２７の作動に伴う吸気弁閉時期の変化とに対応するため、特にタンブルの影響を燃焼速度の推定式に加える。
まず、燃焼速度Ｓｂ［ｍ／ｓｅｃ］は次の（１）式により表される。（１）式右辺第２項は、乱流燃焼速度ＳＴを与える式で、このように乱流燃焼速度ＳＴを、ｕ／ＳＬを変数としてこの変数の指数関数で与えるとする点は公知である。
Ｓｂ＝ＳＬ＋ｂ・（ｕ／ＳＬ）＾ａ …（１）
ただし、ＳＬ ：層流燃焼速度、
ｕ ：乱れ強さ、
ａ，ｂ：適合係数、
（１）式の、燃焼室内のガス流動に伴う乱れのうちの１つの特性値である乱れ強さｕ［ｍ／ｓｅｃ］は次の（２）式により表わされる。これは、燃焼室内ガスの乱れ強さｕは回転速度Ｎｅと比例の関係を持つことが、また、その比例定数は、スワール強度、タンブル強度と相関があることが知られているので、これを式に表したものである。
ただし、ＳＬ ：層流燃焼速度、
ｕ ：乱れ強さ、
ａ，ｂ：適合係数、
（１）式の、燃焼室内のガス流動に伴う乱れのうちの１つの特性値である乱れ強さｕ［ｍ／ｓｅｃ］は次の（２）式により表わされる。これは、燃焼室内ガスの乱れ強さｕは回転速度Ｎｅと比例の関係を持つことが、また、その比例定数は、スワール強度、タンブル強度と相関があることが知られているので、これを式に表したものである。
ｕ＝ｃ・Ｉｔ・Ｎｅ …（２）
ただし、Ｉｔ：タンブル強度［無名数］
Ｎｅ：エンジン回転速度［ｒｐｍ］、
ｃ ：適合係数、
ここで、本実施形態では（２）式にスワール強度は入れていない。これは、今回対象としているエンジンにおいては、ピストン６の中心軸に対して周方向に旋回する流れであるスワールは、ピストン６の中心軸に直交する軸に対して周方向に旋回する流れであるタンブルに比べ、乱流燃焼速度に与える影響が少ないと判断し、今回は考慮しないためである。
ただし、Ｉｔ：タンブル強度［無名数］
Ｎｅ：エンジン回転速度［ｒｐｍ］、
ｃ ：適合係数、
ここで、本実施形態では（２）式にスワール強度は入れていない。これは、今回対象としているエンジンにおいては、ピストン６の中心軸に対して周方向に旋回する流れであるスワールは、ピストン６の中心軸に直交する軸に対して周方向に旋回する流れであるタンブルに比べ、乱流燃焼速度に与える影響が少ないと判断し、今回は考慮しないためである。
ただし、本発明ではスワール強度は対象外というのではない。対象とするエンジンにおいてタンブルよりスワールのほうが乱流燃焼速度に与える影響が大きいときには、スワール強度Ｉｓを、後述するタンブル強度Ｉｔの推定方法と同様の推定方法により推定すればよい。
次に、（２）式のタンブル強度Ｉｔ（燃焼室内のガス流動）の推定方法を説明する。まず、図２（ａ）〜図２（ｅ）は燃焼室５内におけるタンブルの生成から消滅までをピストン６の動きに合わせて図解したものである。順に説明すると、図２（ａ）は吸気弁開時期ＩＶＯ、つまり吸気弁１５が開いた瞬間で、燃焼室５内においてさまざまな方向にガスが流れ込んでいる。図２（ｂ）はピストン６が下降している吸気上死点後４５ｄｅｇ付近（４５ｄｅｇＡＴＤＣ）を示し、この４５ｄｅｇＡＴＤＣ付近からタンブルが形成され始める。図２（ｃ）は吸気弁閉時期ＩＶＣ付近を示し、吸気弁閉時期ＩＶＣ付近までタンブル強度が増加していく。したがって、吸気弁閉時期ＩＶＣが遅いほどタンブルが強くなる。
図２（ｄ）は吸気上死点後３１５ｄｅｇ付近まで、つまりピストン６が上昇する圧縮行程中を示し、圧縮工程中に徐徐にタンブルが小渦（乱れ）に変換されていく。図２（ｅ）は吸気上死点後３１５ｄｅｇ付近を示し、この３１５ｄｅｇＡＴＤＣ付近でタンブルは消え全て小渦に変換される。
一方、図２（ｆ）と図２（ｇ）はＶＥＬ機構２８の作動、非作動に伴う吸気弁１５のバルブリフト量の相違でタンブルの形成がどのように違うのかを示している。このうち、左側に示す図２（ｆ）は吸気弁１５のバルブリフト量が小さい場合で、バルブリフト量が小さいときには、吸気がシリンダ壁側とシリンダ中心側の両方に分かれて流れ、お互いに流れを打ち消しあうため、タンブルが形成されにくい。これに対して右側に示す図２（ｇ）は吸気弁１５のバルブリフト量が大きい場合で、バルブリフト量が大きいときには、ほとんどの吸気がシリンダ中心側へと流れるため、タンブルを形成しやすい。
さて、タンブル強度Ｉｔ［無名数］は次式により小渦の有する角運動量の総和として算出されることがＳＡＥペーパーに記載されている（ＳＡＥ９８１０４８）。
Ｉｔ＝ΣＭiｒi×ｖi …（３）
ただし、ｍi：ｉ番目の小渦におけるガスの質量、
ｒi：シリンダの中心からｉ番目の小渦までの距離、
ｖi：ｉ番目の小渦でのガス流速、
しかしながら、燃焼室内の全ての小渦についてこれらの値ｍi、ｒi、ｖiをエンジンコントローラ３１においてオンラインで算出することは現在のところ不可能であるので、上記（３）式に代えて次の（４）式でタンブル強度Ｉｔ［無名数］を近似する。
ただし、ｍi：ｉ番目の小渦におけるガスの質量、
ｒi：シリンダの中心からｉ番目の小渦までの距離、
ｖi：ｉ番目の小渦でのガス流速、
しかしながら、燃焼室内の全ての小渦についてこれらの値ｍi、ｒi、ｖiをエンジンコントローラ３１においてオンラインで算出することは現在のところ不可能であるので、上記（３）式に代えて次の（４）式でタンブル強度Ｉｔ［無名数］を近似する。
Ｉｔ＝ａ・（ｍ＾ｂ・Ｎｅ＾ｃ） …（４）
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ａ ：係数、
ｂ、ｃ：適合係数、
（４）式のｍ＾ｂはタンブル強度に対する燃焼室内ガス質量ｍの、またＮｅ＾ｃはタンブル強度に対する回転速度Ｎｅの各影響を反映させたものである。ここで、燃焼室内ガス質量ｍは吸気弁閉時期ＩＶＣに燃焼室５に入っているガス質量のことである。
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ａ ：係数、
ｂ、ｃ：適合係数、
（４）式のｍ＾ｂはタンブル強度に対する燃焼室内ガス質量ｍの、またＮｅ＾ｃはタンブル強度に対する回転速度Ｎｅの各影響を反映させたものである。ここで、燃焼室内ガス質量ｍは吸気弁閉時期ＩＶＣに燃焼室５に入っているガス質量のことである。
さらに、タンブルは上記図２（ａ）〜図２（ｅ）で示したように、吸気弁閉時期ＩＶＣが４５ｄｅｇＡＴＤＣ付近より遅いほど強くなり、また、図２（ｆ）、図２（ｇ）で示したように吸気弁１５のバルブリフト量が大きいほうが強くなると考え、（４）式の係数ａを新たに次の（５）式で表す。
ａ＝ｋＶＥＬ・（ＩＶＣ−θ０ｔ） …（５）
ただし、ｋＶＥＬ：タンブル強度のバルブリフト量補正係数［１／ｄｅｇ］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θ０ｔ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
従って、（５）式を（４）式に代入した次の（６）式によりタンブル強度Ｉｔを近似する。
ただし、ｋＶＥＬ：タンブル強度のバルブリフト量補正係数［１／ｄｅｇ］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θ０ｔ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
従って、（５）式を（４）式に代入した次の（６）式によりタンブル強度Ｉｔを近似する。
Ｉｔ＝ｋＶＥＬ・（ＩＶＣ−θ０ｔ）・ｍ＾ｂ・Ｎｅ＾ｃ …（６）
ただし、ｋＶＥＬ：タンブル強度のバルブリフト量補正係数［１／ｄｅｇ］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θ０ｔ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｂ、ｃ ：適合係数、
（６）式の吸気弁閉時期ＩＶＣの単位は吸気上死点を起点として遅角側に計測する値［ｄｅｇＡＴＤＣ］であるので、吸気弁閉時期ＩＶＣが吸気上死点より遅れるほどＩＶＣの値が大きくなり、従って（６）式よりタンブル強度Ｉｔが大きくなる。また、（６）式のバルブリフト量補正係数ｋＶＥＬはバルブリフト量の関数であり、バルブリフト量が大きいときには、バルブリフト量が小さいときより大きくなり、従って（６）式よりタンブル強度Ｉｔが大きくなる。
ただし、ｋＶＥＬ：タンブル強度のバルブリフト量補正係数［１／ｄｅｇ］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θ０ｔ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｂ、ｃ ：適合係数、
（６）式の吸気弁閉時期ＩＶＣの単位は吸気上死点を起点として遅角側に計測する値［ｄｅｇＡＴＤＣ］であるので、吸気弁閉時期ＩＶＣが吸気上死点より遅れるほどＩＶＣの値が大きくなり、従って（６）式よりタンブル強度Ｉｔが大きくなる。また、（６）式のバルブリフト量補正係数ｋＶＥＬはバルブリフト量の関数であり、バルブリフト量が大きいときには、バルブリフト量が小さいときより大きくなり、従って（６）式よりタンブル強度Ｉｔが大きくなる。
（６）式によりタンブル強度Ｉｔを近似した新しい考え方が妥当かどうかを確認するため、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンについて、ＶＥＬ機構２８は非作動状態つまり吸気弁１５のバルブリフト量を一定として、ＶＴＣ機構２７のみを作動し吸気弁閉時期ＩＶＣを変化させる実験を行ったところ、図３に示すようにタンブル強度と吸気弁閉時期ＩＶＣとの関係を示す結果が得られた。図３より吸気弁１５のバルブリフト量一定の条件で、吸気弁閉時期ＩＶＣが大きくなる（つまりＩＶＣが吸気上死点より遅れる）ほどタンブル強度Ｉｔが強くなることが確かめられた。
上記（６）式のタンブル形成開始角θ０ｔは一定値（例えば４５ｄｅｇＡＴＤＣ）である。タンブル形成開始角θ０ｔはエンジン仕様に依存し、ＶＥＬ機構２８、ＶＴＣ機構２７の有無には依存しない。また、図２より、吸気弁閉時期ＩＶＣはタンブル形成開始角θ０ｔより必ず遅角側の値である。従って、上記（６）式の（ＩＶＣ−θ０ｔ）の値は必ず正の値になると考えている。
上記（６）式のタンブル強度のバルブリフト量補正係数ｋＶＥＬはバルブリフト量Ｌｉｆｔ［ｍ］（あるいはバルブ作動角）に応じて定めている。ここでは簡単のため、例えば、ＶＥＬ機構２８の作動、非作動により吸気弁１５のバルブリフト量が大小の２段に切換可能であり、通常の運転条件（運転条件１とする）ではＶＥＬ機構２８を非作動状態とし、運転条件１より所定の運転条件（運転条件２とする）に移行したときにはＶＥＬ機構２８を作動させ、図４下段に示したように実線のバルブリフト特性から破線のバルブリフト特性へとバルブリフト量を小さくしているとして説明すると、ＶＥＬ機構２８の非作動状態でバルブリフト量が第１バルブリフト量Ｌｉｆｔ１［ｍ］となり、ＶＥＬ機構２８の作動状態でバルブリフト量が第２バルブリフト量Ｌｉｆｔ２［ｍ］となる（Ｌｉｆｔ２＜Ｌｉｆｔ１）。従って、図４上段に示したように、ＶＥＬ機構２８の非作動時（あるいは運転条件１のとき）には第１バルブリフト量Ｌｉｆｔ１に対応する第１バルブリフト量補正係数ｋＶＥＬ１を、これに対してＶＥＬ機構２８の作動時（あるいは運転条件２のとき）には第２バルブリフト量Ｌｉｆｔ２に対応する第２バルブリフト量補正係数ｋＶＥＬ２をバルブリフト量補正係数ｋＶＥＬとして設定する。このように、ＶＥＬ機構２８の作動で吸気弁１５のバルブリフト量が第１バルブリフト量Ｌｉｆｔ１より第２バルブリフト量Ｌｉｆｔ２へと小さくなったときに、バルブリフト量補正係数ｋＶＥＬを第１バルブリフト量補正係数ｋＶＥＬ１より第２バルブリフト量補正係数ｋＶＥＬ２へと小さくするのは、バルブリフト量が小さいときのほうがタンブルが弱くなる（従って、タンブル強度Ｉｔが小さくなる）からである。
このように、（６）式によれば、タンブル強度Ｉｔを、新たに、吸気弁１５のバルブリフト量、吸気弁閉時期ＩＶＣ、燃焼室内ガス質量ｍ、エンジン回転速度Ｎｅの関数としている。
〈２〉基準クランク角の設定方法
コンベンショナルエンジンに適用している先行点火時期制御方法では、次の（補４）式のように、燃焼圧が最大となるクランク角θｐｍａｘを、燃焼質量割合が６０％となるときのクランク角位置に設定すると共に、この燃焼圧が最大となるクランク角θｐｍａｘを基準クランク角としている。なお、燃焼質量割合がｘ％となるときのクランク角位置を、以下、「θｍｂ ｘ％」で表す。例えば、燃焼質量割合が６０％となるときのクランク角位置は「θｍｂ６０％」である。また、θｍｂ ｘ％の起点は圧縮上死点とする。
〈２〉基準クランク角の設定方法
コンベンショナルエンジンに適用している先行点火時期制御方法では、次の（補４）式のように、燃焼圧が最大となるクランク角θｐｍａｘを、燃焼質量割合が６０％となるときのクランク角位置に設定すると共に、この燃焼圧が最大となるクランク角θｐｍａｘを基準クランク角としている。なお、燃焼質量割合がｘ％となるときのクランク角位置を、以下、「θｍｂ ｘ％」で表す。例えば、燃焼質量割合が６０％となるときのクランク角位置は「θｍｂ６０％」である。また、θｍｂ ｘ％の起点は圧縮上死点とする。
θｍｂ６０％＝θｐｍａｘ …（補４）
また、コンベンショナルエンジンに適用している先行点火時期制御方法では、燃焼圧が最大となるクランク角θｐｍａｘを次の（補５）式のようにエンジン回転速度Ｎｅの関数で与えている。
また、コンベンショナルエンジンに適用している先行点火時期制御方法では、燃焼圧が最大となるクランク角θｐｍａｘを次の（補５）式のようにエンジン回転速度Ｎｅの関数で与えている。
θｐｍａｘ＝ｈ・（Ｎｅ）＾ｉ …（補５）
ただし、ｈ、ｉ：適合係数、
一方、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンを対象として実験してみたところ、図５に示したように、回転速度Ｎｅが一定の条件のもとで、燃焼圧が最大となるクランク角θｐｍａｘとθｍｂ ｘ％との関係を表す実験結果が得られた。ただし、図５にはθｍｂ ｘ％として燃焼質量割合が０％となるときのクランク角位置であるθｍｂ０％、燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％、燃焼質量割合が６０％となるときのクランク角位置であるθｍｂ６０％の３つの場合だけを示しており、図５において左はθｐｍａｘとθｍｂ０％の、中央はθｐｍａｘとθｍｂ１０％の、右はθｐｍａｘとθｍｂ６０％の各関係を整理したものである。実際には燃焼質量割合が０％、１０％、６０％以外の値（つまり２０％、３０％、４０％、５０％）となるときのクランク角位置での実験結果も得ている（図７参照）。
ただし、ｈ、ｉ：適合係数、
一方、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンを対象として実験してみたところ、図５に示したように、回転速度Ｎｅが一定の条件のもとで、燃焼圧が最大となるクランク角θｐｍａｘとθｍｂ ｘ％との関係を表す実験結果が得られた。ただし、図５にはθｍｂ ｘ％として燃焼質量割合が０％となるときのクランク角位置であるθｍｂ０％、燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％、燃焼質量割合が６０％となるときのクランク角位置であるθｍｂ６０％の３つの場合だけを示しており、図５において左はθｐｍａｘとθｍｂ０％の、中央はθｐｍａｘとθｍｂ１０％の、右はθｐｍａｘとθｍｂ６０％の各関係を整理したものである。実際には燃焼質量割合が０％、１０％、６０％以外の値（つまり２０％、３０％、４０％、５０％）となるときのクランク角位置での実験結果も得ている（図７参照）。
燃焼圧が最大となるクランク角θｐｍａｘと、燃焼質量割合が６０％となるときのクランク角位置であるθｍｂ６０％とでは、図５右側に示したように最大誤差が２．６ｄｅｇもあり、θｐｍａｘとθｍｂ０％との間の誤差を示す図５左側と同様であることがわかる。
次に、図６は、吸入空気量（エンジン負荷）が一定の条件のもとでの、燃焼圧が最大となるクランク角θｐｍａｘとエンジン回転速度Ｎｅとの関係を示したものである。図示のように、傾向としては回転速度Ｎｅが大きくなるほどθｐｍａｘが進角側にずれているが、回転速度Ｎｅに対するθｐｍａｘのバラツキが大きく、燃焼圧が最大となるクランク角θｐｍａｘは上記（補５）式のようにエンジン回転速度Ｎｅだけでは表しきれないことがわかる。
図７は燃焼質量割合を横軸に、燃焼質量割合がｘ％となるときのクランク角位置であるθｍｂ ｘ％をθｐｍａｘで直線近似した場合のＲの２乗値を縦軸に採り、図５に示したθｍｂ０％時、θｍｂ１０％時、θｍｂ６０％時の３つの実験結果及び図示しない図５と同様のθｍｂ２０％時、θｍｂ３０％時、θｍｂ４０％時、θｍｂ５０％時の残りの実験結果とを整理したものである。縦軸のＲの２乗値は、この値が１．０に近づくほどθｐｍａｘとθｍｂ ｘ％（ｘ＝０、１０、２０、３０、４０、５０、６０）との間に相関があることを表すので、図７によれば、燃焼圧が最大となるクランク角θｐｍａｘはθｍｂ１０％付近で一番相関があることがわかる。
なお、図７にはＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンでもエンジン仕様の異なる２種類のエンジンの特性を示している。２つのエンジンの違いは主に燃焼室形状にあり、第２エンジンは第１エンジンよりも燃焼室が平べったい形状のものである。第１エンジンでは燃焼圧が最大となるクランク角θｐｍａｘは燃焼質量割合が１０％となるときのクランク角位置から５０％となるときのクランク角位置との間で広く相関があるのに対して、第２エンジンになると燃焼圧が最大となるクランク角θｐｍａｘは燃焼質量割合が１０％となるときのクランク角位置付近でだけ相関がある。従って、２つのエンジンに共通して燃焼圧が最大となるクランク角θｐｍａｘが高い相関を有するのは燃焼質量割合が１０％となるときのクランク角位置付近にあるときだけである。このように、エンジン機種が相違しても、燃焼圧が最大となるクランク角θｐｍａｘはθｍｂ１０％となるときのクランク角位置付近で一番相関がある、という結果が初めて得られた。
こうした実験結果を受けて本実施形態では、燃焼圧が最大となるクランク角θｐｍａｘに基づいて、燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％を算出し、このθｍｂ１０％を基準クランク角として設定する、つまり次の（７）式によりθｍｂ１０％［ｄｅｇＡＴＤＣ］を算出する。
θｍｂ１０％＝ａ・θｐｍａｘ−ｂ …（７）
ただし、ａ：適合係数［無名数］、
ｂ：適合クランク角［ｄｅｇ］、
（７）式は燃焼圧が最大となるクランク角θｐｍａｘをｂ（正の値）のクランク角だけ進角側にシフトしてθｍｂ１０％を求めるとする式である。この結果、図８に示したように本実施形態での基準クランク角（＝θｍｂ１０％）はコンベンショナルエンジンに適用している先行点火時期制御方法における基準クランク角（＝θｍｂ６０％）よりも進角側にくることとなる。
ただし、ａ：適合係数［無名数］、
ｂ：適合クランク角［ｄｅｇ］、
（７）式は燃焼圧が最大となるクランク角θｐｍａｘをｂ（正の値）のクランク角だけ進角側にシフトしてθｍｂ１０％を求めるとする式である。この結果、図８に示したように本実施形態での基準クランク角（＝θｍｂ１０％）はコンベンショナルエンジンに適用している先行点火時期制御方法における基準クランク角（＝θｍｂ６０％）よりも進角側にくることとなる。
そして、ＶＥＬ機構２８及びＶＴＣ機構２７を備えるエンジンでは、燃焼圧が最大となるクランク角θｐｍａｘは、吸気弁閉時期ＩＶＣが遅いほど、また吸気弁１５のバルブリフト量が大きいほど圧縮上死点より遅れると考え、また、燃焼室内ガス質量（負荷）によって影響されると判断し、さらにＦＴＡ（Ｆａｕｌｔ Ｔｒｅｅ Ａｎａｌｙｓｉｓ）より燃焼室内のガス流速（タンブル、スワール）の影響をもつと判断し、次の（８）式により燃焼圧が最大となるクランク角θｐｍａｘ［ｄｅｇＡＴＤＣ］を算出することとした。従って、（８）式と上記（補５）式とを比較すれば、ｋ１ＶＥＬ・（ＩＶＣ−θｔ０）・ｍ＾ｄを新たに導入したものである。
θｐｍａｘ＝ｋ１ＶＥＬ・（ＩＶＣ−θｔ０）・ｍ＾ｄ・Ｎｅ＾ｅ
…（８）
ただし、ｋ１ＶＥＬ：流速係数［無名数］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ｄ、ｅ ：適合定数［無名数］、
（８）式に示す燃焼圧が最大となるクランク角θｐｍａｘに対するこうした新しい考え方が妥当かどうかを確認するため、吸気弁１５のバルブリフト量及び吸気弁閉時期ＩＶＣを一定として、燃焼圧が最大となるクランク角θｐｍａｘと吸入空気量（燃焼室内ガス質量ｍの代表値）、エンジン回転速度Ｎｅ（ガス流速の代表値）の関係を実験したところ、図９に示す結果が得られた。図９によれば、吸気弁１５のバルブリフト量及び吸気弁閉時期ＩＶＣが一定の条件で、吸入空気量が同じであればエンジン回転速度Ｎｅが大きいときのほうがエンジン回転速度Ｎｅが小さいときより燃焼圧が最大となるクランク角θｐｍａｘが大きくなる（圧縮上死点よりの遅角量が大きくなる）こと、またエンジン回転速度Ｎｅが同じであれば吸入空気量が多くなるほど燃焼圧が最大となるクランク角θｐｍａｘが大きくなる（圧縮上死点よりの遅角量が大きくなる）ことが確かめられた。
…（８）
ただし、ｋ１ＶＥＬ：流速係数［無名数］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ｄ、ｅ ：適合定数［無名数］、
（８）式に示す燃焼圧が最大となるクランク角θｐｍａｘに対するこうした新しい考え方が妥当かどうかを確認するため、吸気弁１５のバルブリフト量及び吸気弁閉時期ＩＶＣを一定として、燃焼圧が最大となるクランク角θｐｍａｘと吸入空気量（燃焼室内ガス質量ｍの代表値）、エンジン回転速度Ｎｅ（ガス流速の代表値）の関係を実験したところ、図９に示す結果が得られた。図９によれば、吸気弁１５のバルブリフト量及び吸気弁閉時期ＩＶＣが一定の条件で、吸入空気量が同じであればエンジン回転速度Ｎｅが大きいときのほうがエンジン回転速度Ｎｅが小さいときより燃焼圧が最大となるクランク角θｐｍａｘが大きくなる（圧縮上死点よりの遅角量が大きくなる）こと、またエンジン回転速度Ｎｅが同じであれば吸入空気量が多くなるほど燃焼圧が最大となるクランク角θｐｍａｘが大きくなる（圧縮上死点よりの遅角量が大きくなる）ことが確かめられた。
なお、図９の縦軸の単位［ｄｅｇＡＴＤＣ］は圧縮上死点を起点とするクランク角であり、吸気上死点を起点とするものでない。
上記（８）式のタンブル形成開始角θ０ｔについては上記（６）式のところで前述したところと同じである。すなわち、タンブル形成開始角θ０ｔは一定値（例えば４５ｄｅｇＡＴＤＣ）である。タンブル形成開始角θ０ｔはエンジンの仕様に依存し、ＶＥＬ機構２８、ＶＴＣ機構２７の有無には依存しない。図２より、吸気弁閉時期ＩＶＣはタンブル形成開始角θ０ｔより必ず遅角側の値である。従って、（８）式の（ＩＶＣ−θ０ｔ）の値は必ず正の値になる。
上記（８）式の流速係数ｋ１ＶＥＬは吸気弁１５のバルブリフト量Ｌｉｆｔ（あるいはバルブ作動角）に応じて定めている。ここでも簡単のため、前述したようにＶＥＬ機構２８の作動、非作動により吸気弁１５のバルブリフト量が大小の２段に切換可能であり、運転条件１ではＶＥＬ機構２８を非作動状態とし、運転条件１より運転条件２に移行したときにはＶＥＬ機構２８を作動させ、図１０下段に示したように実線のバルブリフト特性から破線のバルブリフト特性へとバルブリフト量を小さくしているとして説明すると、図１０上段に示したように、ＶＥＬ機構２８の非作動時（あるいは運転条件１のとき）には第１バルブリフト量Ｌｉｆｔ１に対応する第１流速係数ｋ１ＶＥＬ１を、これに対してＶＥＬ機構２８の作動時（あるいは運転条件２のとき）には第２バルブリフト量Ｌｉｆｔ２に対応する第２流速係数ｋ１ＶＥＬ２を流速係数ｋ１ＶＥＬとして設定する。このように、バルブリフト量が小さくなったときに流速係数ｋ１ＶＥＬを第１流速係数ｋ１ＶＥＬ１より第２流速係数ｋ１ＶＥＬ２へと小さくするのは、バルブリフト量が小さいときのほうがタンブルが弱くなる（従って、θｐｍａｘがより遅角側に移動する）からである。
また、上記（８）式のｍ＾ｄはθｐｍａｘに対する筒内ガス質量ｍの、Ｎｅ＾ｅはθｐｍａｘに対する回転速度Ｎｅの各影響を反映させている。
このように、（８）式によれば、燃焼圧が最大となるクランク角θｐｍａｘを、エンジン回転速度Ｎｅに加えて、燃焼室内ガス質量ｍ、吸気弁１５のバルブリフト量及び吸気弁閉時期ＩＶＣの関数でもあるとして新たに構成している。
〈３〉着火遅れ時間の算出方法
コンベンショナルエンジンに適用している先行点火時期制御方法では、点火タイミングを起点とし火炎形成開始時期までの時間を点火無駄時間ＤＥＡＤＴＩＭＥとして適合している。これを燃焼質量割合でいうと、点火タイミングを始期とし、燃焼質量割合が０％より上昇を開始するタイミングを終期とする時間である。
〈３〉着火遅れ時間の算出方法
コンベンショナルエンジンに適用している先行点火時期制御方法では、点火タイミングを起点とし火炎形成開始時期までの時間を点火無駄時間ＤＥＡＤＴＩＭＥとして適合している。これを燃焼質量割合でいうと、点火タイミングを始期とし、燃焼質量割合が０％より上昇を開始するタイミングを終期とする時間である。
この場合に、火炎形成開始時期は、一般的に燃焼室内圧力に基づいて計測している。すなわち、火炎形成開始前後で燃焼室内圧力が変化するはずであるから、その燃焼室内圧力の変化タイミングをもって火炎形成開始時期としている。
しかしながら、燃焼質量割合が０％より上昇を開始するタイミングは冷却損失の算出誤差やガス漏れの影響を受けやすく、燃焼解析装置で正確に火炎形成開始時期を計測することが難しい。そのため、燃焼解析装置で計測した値を正として点火無駄時間ＤＥＡＤＴＩＭＥの適合を行うと、後述する基本点火時期ＭＢＴＣＡＬの推定誤差にこの点火無駄時間ＤＥＡＤＴＩＭＥの計測誤差が含まれてしまう。
これに対して本実施形態では、点火タイミングを起点として燃焼質量割合が２％に上昇するまでの時間を着火遅れ時間τであるとみなす。すなわち、本発明で導入しているこの着火遅れ時間τは物理的な意味での着火遅れ時間ではなく、点火タイミングから燃料が２％燃えるまでの時間である。このため、本実施形態では燃焼質量割合が２％より１０％までの時間が燃焼時間となるが、これはコンベンショナルエンジンに適用している先行点火時期制御方法では主燃焼期間ＢＵＲＮ２に相当しているので、本実施形態では、燃焼時間の推定を主燃焼期間のみとすることとなる。
ここで、２％という値は初期燃焼時の燃焼質量割合である。初期燃焼時とは、火炎形成開始時期よりも時間的に後のタイミングのことである。初期燃焼時の燃焼質量割合であるこの２％を、以下では「基準燃焼質量割合」ともいう。なお、本発明は、基準燃焼質量割合を２％とする場合に限定するものでなく、基準燃焼質量割合として２〜１０％を設定すればよい。２〜１０％の範囲で設定してよいとした理由は、似たような初期燃焼の範囲内にあるためである。
また、燃焼解析装置によっては、例えば燃焼質量割合が２％に上昇して以降、精度良く計測できるものと、燃焼質量割合が５％にまで上昇しないと精度良く計測できないものとがあったとして、燃焼質量割合が２％に上昇して以降、精度良く計測できる燃焼解析装置を用いるのであれば、基準燃焼質量割合として２％を設定すればよいし、燃焼質量割合が５％にまで上昇しないと精度良く計測できない燃焼解析装置を用いるのであれば、基準燃焼質量割合として５％を設定すればよい。
このように、本実施形態において点火タイミングから燃焼質量割合が２％になるまでの時間を着火遅れ時間τとして設定することとすると、コンベンショナルエンジンに適用している先行点火時期制御方法において提案しているアレニウスの式（自己着火の式）を用いることができないため、本実施形態では、着火遅れ時間τ［ｓｅｃ］を新たに次の（９）式により算出することとする。
τ＝（Ｄkernel−Ｄ０）
／２［（Ｔａｄ／Ｔ）・ＳＬ＋｛（２／３）・ｋ｝＾（１／２）］
…（９）
ただし、Ｄkernel：θｍｂ２％時の火炎直径［ｍ］、
Ｄ０ ：点火タイミング直後の火炎直径［ｍ］、
Ｔａｄ ：火炎温度［Ｋ］、
Ｔ ：未燃ガス温度［Ｋ］、
ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｋ ：燃焼室内ガスの運動エネルギー、
この（９）式は点火タイミングの瞬間から火炎が形成されるとみなし、点火タイミング直後の火炎速度を算出する式を用いて算出するものである。すなわち、点火タイミング直後の火炎速度ＳＦＬＡＭＥ［ｍ／ｓｅｃ］を算出する式とは、
ＳＦＬＡＭＥ＝Ｒkernel／ｔ …（追１）
ただし、Ｒkernel：燃焼質量割合が２％になったときの火炎半径［ｍ］、
ｔ ：燃焼質量割合が２％になるまでの時間［ｓｅｃ］、
であり、この（追１）式を時間について解くと次式が得られる。
／２［（Ｔａｄ／Ｔ）・ＳＬ＋｛（２／３）・ｋ｝＾（１／２）］
…（９）
ただし、Ｄkernel：θｍｂ２％時の火炎直径［ｍ］、
Ｄ０ ：点火タイミング直後の火炎直径［ｍ］、
Ｔａｄ ：火炎温度［Ｋ］、
Ｔ ：未燃ガス温度［Ｋ］、
ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｋ ：燃焼室内ガスの運動エネルギー、
この（９）式は点火タイミングの瞬間から火炎が形成されるとみなし、点火タイミング直後の火炎速度を算出する式を用いて算出するものである。すなわち、点火タイミング直後の火炎速度ＳＦＬＡＭＥ［ｍ／ｓｅｃ］を算出する式とは、
ＳＦＬＡＭＥ＝Ｒkernel／ｔ …（追１）
ただし、Ｒkernel：燃焼質量割合が２％になったときの火炎半径［ｍ］、
ｔ ：燃焼質量割合が２％になるまでの時間［ｓｅｃ］、
であり、この（追１）式を時間について解くと次式が得られる。
ｔ＝Ｒkernel／ＳＦＬＡＭＥ …（追２）
そして、この（追２）式の時間ｔを着火遅れ時間τで置き換える、つまり次の（追３）式により着火遅れ時間τ［ｓｅｃ］を算出するものである。
そして、この（追２）式の時間ｔを着火遅れ時間τで置き換える、つまり次の（追３）式により着火遅れ時間τ［ｓｅｃ］を算出するものである。
τ＝Ｒkernel／ＳＦＬＡＭＥ …（追３）
ただし、Ｒkernel ：燃焼質量割合が２％のときの火炎半径［ｍ］、
ＳＦＬＡＭＥ：点火タイミング直後の火炎速度［ｍ／ｓｅｃ］、
言い換えると、（追３）式は、点火タイミング直後に種火（火炎核）より火炎が球状に拡がっていく（成長する）ことに着目し、このときの火炎速度をＳＦＬＡＭＥとし、点火タイミングから時間τが経過して基準燃焼質量割合になったときに火炎半径がＲkernelになっているとしてこれらＳＦＬＡＭＥ、τ、Ｒkernelの関係を表した式である。このようにして、点火タイミング直後の火炎速度ＳＦＬＡＭＥと基準燃焼質量割合になったときの火炎半径Ｒkernelとに基づいて着火遅れ時間τを算出する。
ただし、Ｒkernel ：燃焼質量割合が２％のときの火炎半径［ｍ］、
ＳＦＬＡＭＥ：点火タイミング直後の火炎速度［ｍ／ｓｅｃ］、
言い換えると、（追３）式は、点火タイミング直後に種火（火炎核）より火炎が球状に拡がっていく（成長する）ことに着目し、このときの火炎速度をＳＦＬＡＭＥとし、点火タイミングから時間τが経過して基準燃焼質量割合になったときに火炎半径がＲkernelになっているとしてこれらＳＦＬＡＭＥ、τ、Ｒkernelの関係を表した式である。このようにして、点火タイミング直後の火炎速度ＳＦＬＡＭＥと基準燃焼質量割合になったときの火炎半径Ｒkernelとに基づいて着火遅れ時間τを算出する。
（追３）式の点火タイミング直後の火炎速度ＳＦＬＡＭＥは次の式により与えられる。
ＳＦＬＡＭＥ＝２［（Ｔａｄ／Ｔ）・ＳＬ＋｛（２／３）・ｋ｝＾（１／２）］
…（追４）
ただし、Ｔａｄ ：火炎温度［Ｋ］、
Ｔ ：未燃ガス温度［Ｋ］、
ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｋ ：燃焼室内ガスの運動エネルギー、
（追４）式そのものは公知で、この（追４）式は点火タイミング直後の火炎速度ＳＦＬＡＭＥが点火タイミング直後の層流燃焼速度ＳＬと、燃焼室内ガスの運動エネルギーｋとに依存することを表している。ここで、未燃ガス温度Ｔは点火タイミングから燃焼質量割合が２％に上昇するまでの時間（着火遅れ時間τ）における平均の温度である。Ｔａｄ／Ｔは燃焼室内ガスの熱膨張率を表している。
…（追４）
ただし、Ｔａｄ ：火炎温度［Ｋ］、
Ｔ ：未燃ガス温度［Ｋ］、
ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｋ ：燃焼室内ガスの運動エネルギー、
（追４）式そのものは公知で、この（追４）式は点火タイミング直後の火炎速度ＳＦＬＡＭＥが点火タイミング直後の層流燃焼速度ＳＬと、燃焼室内ガスの運動エネルギーｋとに依存することを表している。ここで、未燃ガス温度Ｔは点火タイミングから燃焼質量割合が２％に上昇するまでの時間（着火遅れ時間τ）における平均の温度である。Ｔａｄ／Ｔは燃焼室内ガスの熱膨張率を表している。
また、火炎半径Ｒkernelと火炎直径Ｄkernelとの間には、
Ｒkernel＝Ｄkernel／２ …（追５）
の関係があるので、この（追５）式と（追４）式とを（追３）式に代入すれば、上記（９）式が得られる。
Ｒkernel＝Ｄkernel／２ …（追５）
の関係があるので、この（追５）式と（追４）式とを（追３）式に代入すれば、上記（９）式が得られる。
（９）式の火炎直径Ｄkernelは燃焼質量割合が２％となったときの火炎の直径であるので、次の（１０）式により新たに算出する。
Ｄkernel＝｛（Ｖｃｙｌ・［１−（１−ｘｂ）／（３・ｘｂ＋１）＾（１／κ）］
・６）／π｝＾（１／３） …（１０）
ただし、Ｖｃｙｌ：基準燃焼質量割合になったときの燃焼室内容積［ｍ3］、
ｘｂ ：基準燃焼質量割合（＝０．０２）、
κ ：ポリトロープ指数（＝１．３４）、
（１０）式は次のようにして導いたものである。点火タイミング直後の火炎の形状は球形に近いので、火炎直径Ｄkernelを与える式は次の（１１ａ）式である。
・６）／π｝＾（１／３） …（１０）
ただし、Ｖｃｙｌ：基準燃焼質量割合になったときの燃焼室内容積［ｍ3］、
ｘｂ ：基準燃焼質量割合（＝０．０２）、
κ ：ポリトロープ指数（＝１．３４）、
（１０）式は次のようにして導いたものである。点火タイミング直後の火炎の形状は球形に近いので、火炎直径Ｄkernelを与える式は次の（１１ａ）式である。
Ｄkernel＝｛（Ｖｂ・６）／π｝＾（１／３） …（１１ａ）
ただし、Ｖｂ：基準燃焼質量割合になったときの火炎体積、
熱力学の公式より体積燃焼割合の式は次の（１１ｂ）式で表される（公知）。
ただし、Ｖｂ：基準燃焼質量割合になったときの火炎体積、
熱力学の公式より体積燃焼割合の式は次の（１１ｂ）式で表される（公知）。
体積燃焼割合＝１−（１−ｘｂ）／（３・ｘｂ＋１）＾（１／κ）
…（１１ｂ）
ただし、ｘｂ：基準燃焼質量割合（＝０．０２）、
κ ：ポリトロープ指数（＝１．３４）、
この（１１ｂ）式と、
Ｖｂ＝Ｖｃｙｌ・体積燃焼割合 …（１１ｃ）
ただし、Ｖｃｙｌ：基準燃焼質量割合になったときの燃焼室内容積［ｍ3］、
の関係式とを利用すれば基準燃焼質量割合になったときの火炎体積Ｖｂを求める次の（１１ｄ）式が得られる。
…（１１ｂ）
ただし、ｘｂ：基準燃焼質量割合（＝０．０２）、
κ ：ポリトロープ指数（＝１．３４）、
この（１１ｂ）式と、
Ｖｂ＝Ｖｃｙｌ・体積燃焼割合 …（１１ｃ）
ただし、Ｖｃｙｌ：基準燃焼質量割合になったときの燃焼室内容積［ｍ3］、
の関係式とを利用すれば基準燃焼質量割合になったときの火炎体積Ｖｂを求める次の（１１ｄ）式が得られる。
Ｖｂ＝Ｖｃｙｌ・［１−（１−ｘｂ）／（３・ｘｂ＋１）＾（１／κ）］
…（１１ｄ）
この（１１ｄ）式を（１１ａ）式に代入すれば上記（１０）式が得られる。
…（１１ｄ）
この（１１ｄ）式を（１１ａ）式に代入すれば上記（１０）式が得られる。
上記（９）式の点火タイミング直後の火炎の直径Ｄ０は１ｍｍとする。火炎温度Ｔａｄと未燃ガス温度Ｔの比であるＴａｄ／Ｔは適合定数とする。層流燃焼速度ＳＬの算出方法はコンベンショナルエンジンに適用している先行点火時期制御方法と同じとする。
上記（９）式の、着火遅れ時間τにおける燃焼室内ガスの運動エネルギーｋは次式により算出する。
ｋ＝（１／２）ｍ（ｆ・Ｎｅ）＾２ …（１２Ａ）
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ：エンジン回転速度［ｒｐｍ］、
ｆ ：適合定数、
あるいは、ピストン速度を用いて次式により算出する。
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ：エンジン回転速度［ｒｐｍ］、
ｆ ：適合定数、
あるいは、ピストン速度を用いて次式により算出する。
ｋ＝（１／２）ｍ（ｆ０・ピストン速度）＾２ …（１２Ｂ）
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
ｆ０：適合定数、
（１２Ａ）式、（１５Ｂ）式は次のようにして導いたものである。着火遅れ時間τにおける燃焼室内ガスの平均ガス流速をＶＧＡＳとすれば、燃焼室内ガスの運動エネルギーｋは次の式により算出することができる。
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
ｆ０：適合定数、
（１２Ａ）式、（１５Ｂ）式は次のようにして導いたものである。着火遅れ時間τにおける燃焼室内ガスの平均ガス流速をＶＧＡＳとすれば、燃焼室内ガスの運動エネルギーｋは次の式により算出することができる。
ｋ＝０．５×ｍ×ＶＧＡＳ＾２ …（追６）
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
ＶＧＡＳ：平均ガス流速［ｍ／ｓｅｃ］、
ここで、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳはピストン速度に依存し、ピストン速度はエンジン回転速度Ｎｅに比例するとして次の２つ式を導入する。
ただし、ｍ ：燃焼室内ガス質量［ｇ］、
ＶＧＡＳ：平均ガス流速［ｍ／ｓｅｃ］、
ここで、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳはピストン速度に依存し、ピストン速度はエンジン回転速度Ｎｅに比例するとして次の２つ式を導入する。
ＶＧＡＳ＝ｆ０・ピストン速度 …（追７）
ＶＧＡＳ＝ｆ・Ｎｅ …（追８）
ただし、ｆ０、ｆ：適合係数、
このうち（追８）式を（追６）式に導入すれば、上記（１２Ａ）式が、また（追７）式を（追６）式に導入すれば、上記（１２Ｂ）式が得られる。
ＶＧＡＳ＝ｆ・Ｎｅ …（追８）
ただし、ｆ０、ｆ：適合係数、
このうち（追８）式を（追６）式に導入すれば、上記（１２Ａ）式が、また（追７）式を（追６）式に導入すれば、上記（１２Ｂ）式が得られる。
これで、項分け説明を終了する。
次に、上記〈１〉〜〈３〉の変更を行った後の新しい点火時期制御方法（この新しい点火時期制御方法を以下「変更後点火時期制御方法」という。）を、以下にまとめる。ただし、適合係数については改めて取り直している。
ＭＢＴの得られる点火時期（この点火時期を「基本点火時期」という。）ＭＢＴＣＡＬ［ｄｅｇＢＴＤＣ］を次式により算出する。
ＭＢＴＣＡＬ＝−｜（τ＋ＢＴ）・Ｎｅ・６−θｍｂ１０％｜…（１３）
ただし、ＢＴ ：燃焼時間［ｓｅｃ］、
τ ：着火遅れ時間［ｓｅｃ］、
θｍｂ１０％：基準クランク角［ｄｅｇＡＴＤＣ］、
これを図解したのが図１１である。変更後火時期制御方法では、点火タイミングから燃焼質量割合が基準燃焼質量割合（＝２％）に上昇するまでの時間である着火遅れ時間τ［ｓｅｃ］と、燃焼質量割合が２％から１０％になるまでの時間である燃焼時間ＢＴ［ｓｅｃ］とを加算し、この加算値［ｓｅｃ］にＮｅ［ｒｐｍ］・６を乗算することによってクランク角区間［ｄｅｇ］に変換する。この換算されたクランク角区間（（τ＋ＢＴ）・Ｎｅ・６）は点火タイミングから基準クランク角（＝θｐｍａｘ）であるθｍｂ１０％までのクランク角区間である。従って、基準クランク角であるθｍｂ１０％［ｄｅｇＡＴＤＣ］よりこの換算されたクランク角区間だけ進角側の値を基本点火時期ＭＢＴＣＡＬとして算出する。
ただし、ＢＴ ：燃焼時間［ｓｅｃ］、
τ ：着火遅れ時間［ｓｅｃ］、
θｍｂ１０％：基準クランク角［ｄｅｇＡＴＤＣ］、
これを図解したのが図１１である。変更後火時期制御方法では、点火タイミングから燃焼質量割合が基準燃焼質量割合（＝２％）に上昇するまでの時間である着火遅れ時間τ［ｓｅｃ］と、燃焼質量割合が２％から１０％になるまでの時間である燃焼時間ＢＴ［ｓｅｃ］とを加算し、この加算値［ｓｅｃ］にＮｅ［ｒｐｍ］・６を乗算することによってクランク角区間［ｄｅｇ］に変換する。この換算されたクランク角区間（（τ＋ＢＴ）・Ｎｅ・６）は点火タイミングから基準クランク角（＝θｐｍａｘ）であるθｍｂ１０％までのクランク角区間である。従って、基準クランク角であるθｍｂ１０％［ｄｅｇＡＴＤＣ］よりこの換算されたクランク角区間だけ進角側の値を基本点火時期ＭＢＴＣＡＬとして算出する。
なお、（１３）式において絶対値をとりマイナスの符号を付けているのは、θｍｂ１０％の単位は圧縮上死点より遅角側に計測するクランク角［ｄｅｇＡＴＤＣ］であるのに対して、ＭＢＴＣＡＬの単位は圧縮上死点より進角側に計測するクランク角［ｄｅｇＢＴＤＣ］であるため、θｍｂ１０％の単位をＭＢＴＣＡＬの単位へと変換してやる必要があるためである。
ここで、基準クランク角であるθｍｂ１０％［ｄｅｇＡＴＤＣ］は、次の（１４）式により算出する。
θｍｂ１０％＝ｋ１ＶＥＬ・（ＩＶＣ−θｔ０）・ｍ＾ｃ１・Ｎｅ＾ｃ２＋ｃ３
…（１４）
ただし、ｋ１ＶＥＬ：流速係数［無名数］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ｃ１〜ｃ３：適合係数、
前述したように、吸気弁閉時期ＩＶＣとタンブル形成開始角θｔ０の単位は吸気上死点を起点としており、一方、θｍｂ１０％は圧縮上死点を起点としている。従って、吸気上死点を起点とする単位より圧縮上死点を起点とする単位への変換は適合係数ｃ３により行うこととなる。例えば４気筒エンジンでは吸気上死点と圧縮上死点の間に１８０ｄｅｇのずれがあるので、適合係数ｃ３に１８０ｄｅｇを入れてやればよい。
…（１４）
ただし、ｋ１ＶＥＬ：流速係数［無名数］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ｃ１〜ｃ３：適合係数、
前述したように、吸気弁閉時期ＩＶＣとタンブル形成開始角θｔ０の単位は吸気上死点を起点としており、一方、θｍｂ１０％は圧縮上死点を起点としている。従って、吸気上死点を起点とする単位より圧縮上死点を起点とする単位への変換は適合係数ｃ３により行うこととなる。例えば４気筒エンジンでは吸気上死点と圧縮上死点の間に１８０ｄｅｇのずれがあるので、適合係数ｃ３に１８０ｄｅｇを入れてやればよい。
（１３）式の燃焼時間ＢＴ［ｓｅｃ］は次式により算出する。
ＢＴ＝０．０８・Ｖｃｙｌ／｛（（３・ｘｂ＋１）＾（１／κ））・Ａｂ・Ｓｂ｝
＝０．０７０７・Ｖｃｙｌ／（Ａｂ・Ｓｂ） …（１５）
ただし、Ｖｃｙｌ：燃焼室内容積［ｍ3］、
ｘｂ ：燃焼時間算出時の平均燃焼質量割合（＝０．０６）、
κ ：ポリトロープ指数、
Ａｂ ：火炎表面積［ｍ2］、
Ｓｂ ：燃焼速度［ｍ／ｓｅｃ］、
（１５）式はコンベンショナルエンジンに適用している先行点火時期制御方法における式と基本的に同じ式である。すなわち、（１５）式においてＶｃｙｌ／｛（３・ｘｂ＋１）＾（１／κ）｝は既燃ガスの質量割合であり、燃焼時間ＢＴはこの既燃ガスの質量割合に比例し、燃焼速度Ｓｂに反比例するとする式である。（１５）式の火炎表面積Ａｂはエンジンの仕様から算出することができる。
＝０．０７０７・Ｖｃｙｌ／（Ａｂ・Ｓｂ） …（１５）
ただし、Ｖｃｙｌ：燃焼室内容積［ｍ3］、
ｘｂ ：燃焼時間算出時の平均燃焼質量割合（＝０．０６）、
κ ：ポリトロープ指数、
Ａｂ ：火炎表面積［ｍ2］、
Ｓｂ ：燃焼速度［ｍ／ｓｅｃ］、
（１５）式はコンベンショナルエンジンに適用している先行点火時期制御方法における式と基本的に同じ式である。すなわち、（１５）式においてＶｃｙｌ／｛（３・ｘｂ＋１）＾（１／κ）｝は既燃ガスの質量割合であり、燃焼時間ＢＴはこの既燃ガスの質量割合に比例し、燃焼速度Ｓｂに反比例するとする式である。（１５）式の火炎表面積Ａｂはエンジンの仕様から算出することができる。
なお、（１５）式の燃焼質量割合ｘｂとしては、θｍｂ２％からθｍｂ１０％までのクランク角区間（燃焼期間）における値を用いる必要がある。この場合に、燃焼質量割合はθｍｂ２％からθｍｂ１０％までのクランク角区間で２％から１０％へと変化するので、ここでは、θｍｂ２％とθｍｂ１０％とを平均した値であるθｍｂ６％時の燃焼質量割合の値つまり６％を用いる。
上記（１５）式の燃焼速度Ｓｂ［ｍ／ｓｅｃ］は次式により算出する。
Ｓｂ＝ＳＬ＋［｛ｋＶＥＬ・（ＩＶＣ−θｔ０）・ｍ＾ｃ４・Ｎｅ＾ｃ５｝
／ＳＬ］＾ｃ６ …（１６）
ただし、ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｋＶＥＬ ：タンブル強度のバルブリフト量補正係数［１／ｄｅｇ］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θ０ｔ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｃ４、ｃ５：適合係数、
ｃ６ ：適合係数、
（１６）式右辺第２項が乱流燃焼速度［ｍ／ｓｅｃ］を与える式であり、乱流燃焼速度を（１６）式右辺第２項で近似している点が新しい。つまり、本実施形態では、燃焼室内ガスの乱れ強さｕ［ｍ／ｓｅｃ］を次の（１６−１）式により算出し、この燃焼室内ガスの乱れ強さｕに基づいて次の（１６−２）式により乱流燃焼速度ＳＴ［ｍ／ｓｅｃ］を算出している。
／ＳＬ］＾ｃ６ …（１６）
ただし、ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｋＶＥＬ ：タンブル強度のバルブリフト量補正係数［１／ｄｅｇ］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θ０ｔ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｃ４、ｃ５：適合係数、
ｃ６ ：適合係数、
（１６）式右辺第２項が乱流燃焼速度［ｍ／ｓｅｃ］を与える式であり、乱流燃焼速度を（１６）式右辺第２項で近似している点が新しい。つまり、本実施形態では、燃焼室内ガスの乱れ強さｕ［ｍ／ｓｅｃ］を次の（１６−１）式により算出し、この燃焼室内ガスの乱れ強さｕに基づいて次の（１６−２）式により乱流燃焼速度ＳＴ［ｍ／ｓｅｃ］を算出している。
ｕ＝ｋＶＥＬ・（ＩＶＣ−θｔ０）・ｍ＾ｃ４・Ｎｅ＾ｃ５ …（１６−１）
ＳＴ＝［ｕ／ＳＬ］＾ｃ６ …（１６−２）
これに対して、（１６）式の層流燃焼速度ＳＬ［ｍ／ｓｅｃ］としては、コンベンショナルエンジンに適用している先行点火時期制御方法と同じでよく、従って次の（１７）式により算出する。
ＳＴ＝［ｕ／ＳＬ］＾ｃ６ …（１６−２）
これに対して、（１６）式の層流燃焼速度ＳＬ［ｍ／ｓｅｃ］としては、コンベンショナルエンジンに適用している先行点火時期制御方法と同じでよく、従って次の（１７）式により算出する。
ＳＬ＝ＳＬ０・（Ｔ／２９８）＾ａ・（ｐ／１０１．３２５）＾ｂ
・（１−２．１×ＭＲＥＳＲ＾ｃ） …（１７）
ただし、ＳＬ０ ：標準状態での層流燃焼速度［ｍ／ｓｅｃ］、
Ｔ ：未燃ガス温度［Ｋ］、
ｐ ：燃焼室内圧力［ｋＰａ］、
ＭＲＥＳＲ：内部不活性ガス率［％］、
ａ、ｂ、ｃ：係数、
一方、上記（１３）式の着火遅れ時間τ［ｓｅｃ］は次の（１８）式により新たに算出する。
・（１−２．１×ＭＲＥＳＲ＾ｃ） …（１７）
ただし、ＳＬ０ ：標準状態での層流燃焼速度［ｍ／ｓｅｃ］、
Ｔ ：未燃ガス温度［Ｋ］、
ｐ ：燃焼室内圧力［ｋＰａ］、
ＭＲＥＳＲ：内部不活性ガス率［％］、
ａ、ｂ、ｃ：係数、
一方、上記（１３）式の着火遅れ時間τ［ｓｅｃ］は次の（１８）式により新たに算出する。
τ＝［｛｛Ｖｃｙｌ・［１−（１−ｘｂ）／（３・ｘｂ＋１）＾（１／κ）］・６｝ ／π｝＾（１／３）−０．００１］
／｛ｃ７・ＳＬ＋ｃ８・Ｎｅ・ｍ＾（１／２）｝
＝｛０．４９０・Ｖｃｙｌ＾（１／３）｝
／｛ｃ７・ＳＬ＋ｃ８・Ｎｅ・ｍ＾（１／２）｝ …（１８）
ただし、ｘｂ ：基準燃焼質量割合（＝０．０２）、
κ ：ポリトロープ指数（＝１．３４）、
ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
Ｖｃｙｌ ：基準燃焼質量割合になったときの燃焼室内容積［ｍ3］、
ｃ７、ｃ８：適合定数、
上記（１８）式の層流燃焼速度ＳＬ［ｍ／ｓｅｃ］としては、次の（１９）式により算出する。
／｛ｃ７・ＳＬ＋ｃ８・Ｎｅ・ｍ＾（１／２）｝
＝｛０．４９０・Ｖｃｙｌ＾（１／３）｝
／｛ｃ７・ＳＬ＋ｃ８・Ｎｅ・ｍ＾（１／２）｝ …（１８）
ただし、ｘｂ ：基準燃焼質量割合（＝０．０２）、
κ ：ポリトロープ指数（＝１．３４）、
ＳＬ ：層流燃焼速度［ｍ／ｓｅｃ］、
ｍ ：燃焼室内ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
Ｖｃｙｌ ：基準燃焼質量割合になったときの燃焼室内容積［ｍ3］、
ｃ７、ｃ８：適合定数、
上記（１８）式の層流燃焼速度ＳＬ［ｍ／ｓｅｃ］としては、次の（１９）式により算出する。
ＳＬ＝ＳＬ０・（Ｔ／２９８）＾ａ・（ｐ／１０１．３２５）＾ｂ
・（１−２．１×ＭＲＥＳＲ＾ｃ） …（１９）
ただし、ＳＬ０ ：標準状態での層流燃焼速度［ｍ／ｓｅｃ］、
Ｔ ：未燃ガス温度［Ｋ］、
ｐ ：燃焼室内圧力［ｋＰａ］、
ＭＲＥＳＲ：内部不活性ガス率［％］、
ａ、ｂ、ｃ：係数、
ここで、着火遅れ時間τを算出するのに用いる（１９）式の層流燃焼速度ＳＬと、燃焼時間ＢＴを算出するのに用いる上記（１７）式の層流燃焼速度ＳＬとでは、式そのものは変わらず、また標準状態での層流燃焼速度ＳＬ０、内部不活性ガス率ＭＲＥＳＲの値も変わらないが、未燃ガス温度Ｔ、燃焼室内圧力ｐとして代入する値が（１９）式と（１７）式とで、後述するように相違することとなる。
・（１−２．１×ＭＲＥＳＲ＾ｃ） …（１９）
ただし、ＳＬ０ ：標準状態での層流燃焼速度［ｍ／ｓｅｃ］、
Ｔ ：未燃ガス温度［Ｋ］、
ｐ ：燃焼室内圧力［ｋＰａ］、
ＭＲＥＳＲ：内部不活性ガス率［％］、
ａ、ｂ、ｃ：係数、
ここで、着火遅れ時間τを算出するのに用いる（１９）式の層流燃焼速度ＳＬと、燃焼時間ＢＴを算出するのに用いる上記（１７）式の層流燃焼速度ＳＬとでは、式そのものは変わらず、また標準状態での層流燃焼速度ＳＬ０、内部不活性ガス率ＭＲＥＳＲの値も変わらないが、未燃ガス温度Ｔ、燃焼室内圧力ｐとして代入する値が（１９）式と（１７）式とで、後述するように相違することとなる。
上記（１７）式及び（１９）式の標準状態での層流燃焼速度ＳＬ０［ｍ／ｓｅｃ］、係数ａ、ｂ、ｃ［無名数］は次の（２０）式〜（２３）式により算出する。
ＳＬ０＝（０．２６３２−０．８４７２／（φ−１．１３）＾２）
…（２０）
ａ＝２．１８−０．８０・（φ−１） …（２１）
ｂ＝−０．１６＋０．２２・（φ−１） …（２２）
ｃ＝１ …（２３）
ただし、φ：当量比［無名数］、
ここで、（２０）式〜（２３）式についてはＳＡＥペーパーにより公知である（ＳＡＥ１９９９１０１７５参照）。
…（２０）
ａ＝２．１８−０．８０・（φ−１） …（２１）
ｂ＝−０．１６＋０．２２・（φ−１） …（２２）
ｃ＝１ …（２３）
ただし、φ：当量比［無名数］、
ここで、（２０）式〜（２３）式についてはＳＡＥペーパーにより公知である（ＳＡＥ１９９９１０１７５参照）。
これで変更後点火時期制御方法をまとめたものの説明を終了する。
次に、上記（１１ｂ）式に示した体積燃焼割合の求め方と、上記（１７）式、（１９）式に示した未燃ガス温度の算出方法とについて補足説明を行う。
まず、上記（１１ｂ）式に示した体積燃焼割合の求め方を簡単に説明する。
既燃ガスの質量燃焼割合をｘ、体積燃焼割合をｙとして、燃焼室内の燃焼が定容燃焼同様に、
ｙ＝ｆ（ｘ）
＝１＋（ｘ−１）／（１＋ｘ（ｋ１−１）＾（１／κ））
…（２４）
ただし、ｋ１≒４〜５、
の式により書けると仮定すると、（２４）式のｋ１に４を代入して計算することにより、上記（１１ｂ）式が容易に得られる。
ｙ＝ｆ（ｘ）
＝１＋（ｘ−１）／（１＋ｘ（ｋ１−１）＾（１／κ））
…（２４）
ただし、ｋ１≒４〜５、
の式により書けると仮定すると、（２４）式のｋ１に４を代入して計算することにより、上記（１１ｂ）式が容易に得られる。
上記（１７）式、（１９）式に示した未燃ガス温度Ｔの算出方法をまとめて説明すると、断熱変化と仮定し、熱力学の次の（２５）式を用いて燃焼室５のθｍｂ ｘ％時の未燃ガス温度Ｔｍｂ ｘ％を算出する。
Ｔｍｂ ｘ％＝Ｔｉｖｃ・（Ｐｍｂ ｘ％／Ｐｉｖｃ）＾｛（κ−１）／κ｝ …（２５）
ただし、Ｔｍｂ ｘ％：θｍｂ ｘ％時の未燃ガス温度［Ｋ］、
Ｔｉｖｃ ：ＩＶＣ時の燃焼室内温度［Ｋ］、
Ｐｍｂ ｘ％：θｍｂ ｘ％時の燃焼室内圧力［ｋＰａ］、
Ｐｉｖｃ ：ＩＶＣ時の燃焼室内圧力［ｋＰａ］、
κ ：比熱比（固定値で１．３〜１．４）、
ここで、（２５）式の燃焼室５のθｍｂ ｘ％時の未燃ガス圧力Ｐｍｂ ｘ％の算出方法を説明すると、この燃焼室５のθｍｂ ｘ％時の未燃ガス圧力Ｐｍｂ ｘ％は次の（２６）式により与えられる。
ただし、Ｔｍｂ ｘ％：θｍｂ ｘ％時の未燃ガス温度［Ｋ］、
Ｔｉｖｃ ：ＩＶＣ時の燃焼室内温度［Ｋ］、
Ｐｍｂ ｘ％：θｍｂ ｘ％時の燃焼室内圧力［ｋＰａ］、
Ｐｉｖｃ ：ＩＶＣ時の燃焼室内圧力［ｋＰａ］、
κ ：比熱比（固定値で１．３〜１．４）、
ここで、（２５）式の燃焼室５のθｍｂ ｘ％時の未燃ガス圧力Ｐｍｂ ｘ％の算出方法を説明すると、この燃焼室５のθｍｂ ｘ％時の未燃ガス圧力Ｐｍｂ ｘ％は次の（２６）式により与えられる。
Ｐｍｂ ｘ％＝Ｐｉｖｃ＋Δｐ …（２６）
ただし、Δｐ：燃焼室内未燃ガス圧力の変化量、
ここで、燃焼室内未燃ガス圧力の変化量Δｐは、容積変化による圧力変化分Δｐｖと、燃焼による圧力変化分Δｐｃとに分けて算出する。つまり、次の（２７）式により燃焼室内未燃ガス圧力の変化量Δｐを算出する。
ただし、Δｐ：燃焼室内未燃ガス圧力の変化量、
ここで、燃焼室内未燃ガス圧力の変化量Δｐは、容積変化による圧力変化分Δｐｖと、燃焼による圧力変化分Δｐｃとに分けて算出する。つまり、次の（２７）式により燃焼室内未燃ガス圧力の変化量Δｐを算出する。
Δｐ＝Δｐｖ＋Δｐｃ …（２７）
このうち、（２７）式の容積変化による圧力変化分Δｐｖは熱力学の式である次の（２８）式により算出することができる。
このうち、（２７）式の容積変化による圧力変化分Δｐｖは熱力学の式である次の（２８）式により算出することができる。
Δｐｖ＝Ｐｉｖｃ・｛（Ｖｉｖｃ／Ｖｍｂ ｘ％）＾κ−１｝ …（２８）
ただし、Ｖｉｖｃ ：ＩＶＣ時の燃焼室内容積、
Ｖｍｂ ｘ％：θｍｂ ｘ％時の燃焼室内容積、
ここで、燃焼室５の吸気弁閉時期における容積Ｖｉｖｃはピストン６の位置から求めることができる。燃焼室５のθｍｂ ｘ％時における未燃ガス容積Ｖｍｂ ｘ％については、図１２に示したように、θｍｂ ｘ％とＶｍｂ ｘ％の関係を表す特性を予め作成しておき、θｍｂ ｘ％から当該特性を用いて、燃焼室５のθｍｂ ｘ％時における未燃ガス容積Ｖｍｂ ｘ％を求めればよい。図１２は要するにクランク角に対する燃焼室容積の特性である。
ただし、Ｖｉｖｃ ：ＩＶＣ時の燃焼室内容積、
Ｖｍｂ ｘ％：θｍｂ ｘ％時の燃焼室内容積、
ここで、燃焼室５の吸気弁閉時期における容積Ｖｉｖｃはピストン６の位置から求めることができる。燃焼室５のθｍｂ ｘ％時における未燃ガス容積Ｖｍｂ ｘ％については、図１２に示したように、θｍｂ ｘ％とＶｍｂ ｘ％の関係を表す特性を予め作成しておき、θｍｂ ｘ％から当該特性を用いて、燃焼室５のθｍｂ ｘ％時における未燃ガス容積Ｖｍｂ ｘ％を求めればよい。図１２は要するにクランク角に対する燃焼室容積の特性である。
一方、上記（２７）式の燃焼による圧力変化分Δｐｃは、次の（２９）式により簡易に算出すればよい。
Δｐｃ＝（ｘ／１００）×Δｐｃｔｏｔａｌ …（２９）
ただし、Δｐｃｔｏｔａｌ：燃焼によるトータルの圧力上昇分、
ここで、燃焼によるトータルの圧力上昇分Δｐｃｔｏｔａｌについては、図１３に示したように、燃料噴射パルス幅Ｔｉ［ｍｓｅｃ］（または燃料噴射量）とΔｐｃｔｏｔａｌの関係を表す特性を予め作成しておき、燃料噴射パルス幅Ｔｉ（または燃料噴射量）から当該特性を用いて燃焼によるトータルの圧力上昇分Δｐｃｔｏｔａｌを求める。
ただし、Δｐｃｔｏｔａｌ：燃焼によるトータルの圧力上昇分、
ここで、燃焼によるトータルの圧力上昇分Δｐｃｔｏｔａｌについては、図１３に示したように、燃料噴射パルス幅Ｔｉ［ｍｓｅｃ］（または燃料噴射量）とΔｐｃｔｏｔａｌの関係を表す特性を予め作成しておき、燃料噴射パルス幅Ｔｉ（または燃料噴射量）から当該特性を用いて燃焼によるトータルの圧力上昇分Δｐｃｔｏｔａｌを求める。
これで補足説明を終える。
次に、エンジンコントローラ３１で実行される上記変更後点火時期制御方法における基本点火時期ＭＢＴＣＡＬの算出方法を図１４、図１６、図１７、図１８、図１９のフローチャートを参照しながら詳述する。
図１４は基本点火時期ＭＢＴＣＡＬの算出に必要な各種の物理量を算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。
まずステップ１では、吸気弁閉時期ＩＶＣ［ｄｅｇＡＴＤＣ］、バルブリフト量Ｌｉｆｔ［ｍ］、温度センサ４３により検出されるコレクタ内温度ＴＣＯＬ［Ｋ］、圧力センサ４４により検出されるコレクタ内圧力ＰＣＯＬ［Ｐａ］、温度センサ４５により検出される排気温度ＴＥＸＨ［Ｋ］、内部不活性ガス率ＭＲＥＳＲ［％］、燃料噴射パルス幅Ｔｉ［ｍｓｅｃ］、クランク角センサ（３３、３４）により検出されるエンジン回転速度Ｎｅ［ｒｐｍ］、総ガス質量ＭＧＡＳ［ｇ］を読み込む。
ここで、吸気弁閉時期ＩＶＣはＶＴＣ機構２７に与える指令値から既知である。あるいはカム角センサ３４により実際の吸気弁閉時期を検出してもかまわない。
吸気弁１５のバルブリフト量ＬｉｆｔはＶＥＬ機構２８に与える指令値から既知である。つまり、図４下段、図１０下段に示したように、ＶＥＬ機構２８の非作動時には第１バルブリフト量Ｌｉｆｔ１が、これに対してＶＥＬ機構２８の作動時には第２バルブリフト量Ｌｉｆｔ２がバルブリフト量Ｌｉｆｔである。
燃料噴射パルス幅Ｔｉは図示しない燃料噴射パルス幅の算出ルーチンにおいて算出されている。例えば、シーケンシャル噴射時に燃料インジェクタ２１に与える燃料噴射パルス幅Ｔｉ［ｍｓｅｃ］の算出式は次のようなものである。
Ｔｉ＝Ｔｐ×Ｔｆｂｙａ×（α＋αｍ−１）×２＋Ｔｓ …（３０）
ただし、Ｔｐ ：基本噴射パルス幅［ｍｓｅｃ］、
Ｔｆｂｙａ：［無名数］、
α ：空燃比フィードバック補正係数［無名数］、
αｍ ：空燃比学習値［無名数］、
Ｔｓ ：無効パルス幅［ｍｓｅｃ］、
（３０）式の目標当量比Ｔｆｂｙａは無名数であり、理論空燃比を１４．７とすると、次の（３１）式により表される値である。
ただし、Ｔｐ ：基本噴射パルス幅［ｍｓｅｃ］、
Ｔｆｂｙａ：［無名数］、
α ：空燃比フィードバック補正係数［無名数］、
αｍ ：空燃比学習値［無名数］、
Ｔｓ ：無効パルス幅［ｍｓｅｃ］、
（３０）式の目標当量比Ｔｆｂｙａは無名数であり、理論空燃比を１４．７とすると、次の（３１）式により表される値である。
Ｔｆｂｙａ＝１４．７／目標空燃比 …（３１）
例えば（３１）式より目標空燃比が理論空燃比（１４．７）のときＴｆｂｙａ＝１．０となり、目標空燃比が例えば２２．０といったリーン側の値であるとき、Ｔｆｂｙａは１．０未満の正の値となる。
例えば（３１）式より目標空燃比が理論空燃比（１４．７）のときＴｆｂｙａ＝１．０となり、目標空燃比が例えば２２．０といったリーン側の値であるとき、Ｔｆｂｙａは１．０未満の正の値となる。
クランク角センサはクランクシャフト７のポジションを検出するポジションセンサ３３と、吸気用カムシャフト２５ポジションを検出するフェーズセンサ（＝カム角センサ）３４とからなり、これら２つのセンサ３３、３４からの信号に基づいてエンジン回転速度Ｎｅ［ｒｐｍ］が算出されている。
内部不活性ガス率ＭＲＥＳＦＲは燃焼室５内に残留する不活性ガス量を燃焼室５内の総ガス量で除した値で、その算出については総ガス質量ＭＧＡＳと共に図１６により後述する。
ステップ２では、燃焼室５の吸気弁閉時期ＩＶＣにおける容積（圧縮開始時期の容積）Ｖｉｖｃ［ｍ3］を算出する。燃焼室５の吸気弁閉時期における容積Ｖｉｖｃは、ピストン６のストローク位置によって決まる。ピストン６のストローク位置はエンジンのクランク角位置によって決まる。
図１５を参照して、エンジンのクランクシャフト７１の回転中心７２がシリンダの中心軸７３からオフセットしている場合を考える。コネクティングロッド７４、コネクティングロッド７４とクランクシャフト７１との結節点７５、コネクティングロッド７４とピストンをつなぐピストンピン７６が図に示す関係にあるとする。このときの、燃焼室５の吸気弁閉時期における容積Ｖｉｖｃは次の（４１）式〜（４５）式で表すことができる。
Ｖｉｖｃ＝ｆ１（θivc）＝Ｖｃ＋（π／４）Ｄ2・Ｈivc
…（４１）
Ｖｃ＝（π／４）Ｄ2・Ｈx／（ε−１） …（４２）
Ｈivc＝｛（ＣＮＤ＋ＳＴ2／２）−（ＣＲoff−ＰＩＳoff）2｝1/2
−｛（ＳＴ／２）・cos（θivc＋θoff）｝
＋（ＣＮＤ2−Ｘ2）1/2 …（４３）
Ｘ ＝（ＳＴ／２）・sin（θivc＋θoff）−ＣＲoff＋ＰＩＳoff
…（４４）
θoff＝arcsin｛（ＣＲoff−ＰＩＳoff）／（ＣＮＤ・（ＳＴ／２））｝
…（４５）
ただし、Ｖｃ ：隙間容積［ｍ3］、
ε ：圧縮比、
Ｄ ：シリンダボア径［ｍ］、
ＳＴ ：ピストンの全ストローク［ｍ］、
Ｈivc ：吸気弁閉時期におけるピストンピン７６の
ＴＤＣからの距離［ｍ］、
Ｈx ：ピストンピン７６のＴＤＣからの距離の最大値と最小値の差 ［ｍ］、
ＣＮＤ ：コネクティングロッド７４の長さ［ｍ］、
ＣＲoff ：結節点７５のシリンダ中心軸７３からのオフセット距離
［ｍ］、
ＰＩＳoff：クランクシャフト回転中心７２のシリンダ中心軸７３から のオフセット距離［ｍ］、
θivc ：吸気弁閉時期のクランク角［ｄｅｇＡＴＤＣ］、
θoff ：ピストンピン７６とクランクシャフト回転中心７２とを結ぶ 線がＴＤＣにおいて垂直線となす角度［ｄｅｇ］、
Ｘ ：結節点７５とピストンピン７６との水平距離［ｍ］、
吸気弁閉時期のクランク角θivcは前述のように、エンジンコントローラ３１からＶＴＣ機構２７への指令信号によって決まるので、既知である。（４１）式〜（４５）式にこのときのクランク角θivc（＝ＩＶＣ）を代入すれば、燃焼室５の吸気弁閉時期における容積Ｖｉｖｃを算出することができる。したがって、実用上は燃焼室５の吸気弁閉時期における容積Ｖｉｖｃは吸気弁閉時期ＩＶＣをパラメータとするテーブルで設定したものを用いる。
…（４１）
Ｖｃ＝（π／４）Ｄ2・Ｈx／（ε−１） …（４２）
Ｈivc＝｛（ＣＮＤ＋ＳＴ2／２）−（ＣＲoff−ＰＩＳoff）2｝1/2
−｛（ＳＴ／２）・cos（θivc＋θoff）｝
＋（ＣＮＤ2−Ｘ2）1/2 …（４３）
Ｘ ＝（ＳＴ／２）・sin（θivc＋θoff）−ＣＲoff＋ＰＩＳoff
…（４４）
θoff＝arcsin｛（ＣＲoff−ＰＩＳoff）／（ＣＮＤ・（ＳＴ／２））｝
…（４５）
ただし、Ｖｃ ：隙間容積［ｍ3］、
ε ：圧縮比、
Ｄ ：シリンダボア径［ｍ］、
ＳＴ ：ピストンの全ストローク［ｍ］、
Ｈivc ：吸気弁閉時期におけるピストンピン７６の
ＴＤＣからの距離［ｍ］、
Ｈx ：ピストンピン７６のＴＤＣからの距離の最大値と最小値の差 ［ｍ］、
ＣＮＤ ：コネクティングロッド７４の長さ［ｍ］、
ＣＲoff ：結節点７５のシリンダ中心軸７３からのオフセット距離
［ｍ］、
ＰＩＳoff：クランクシャフト回転中心７２のシリンダ中心軸７３から のオフセット距離［ｍ］、
θivc ：吸気弁閉時期のクランク角［ｄｅｇＡＴＤＣ］、
θoff ：ピストンピン７６とクランクシャフト回転中心７２とを結ぶ 線がＴＤＣにおいて垂直線となす角度［ｄｅｇ］、
Ｘ ：結節点７５とピストンピン７６との水平距離［ｍ］、
吸気弁閉時期のクランク角θivcは前述のように、エンジンコントローラ３１からＶＴＣ機構２７への指令信号によって決まるので、既知である。（４１）式〜（４５）式にこのときのクランク角θivc（＝ＩＶＣ）を代入すれば、燃焼室５の吸気弁閉時期における容積Ｖｉｖｃを算出することができる。したがって、実用上は燃焼室５の吸気弁閉時期における容積Ｖｉｖｃは吸気弁閉時期ＩＶＣをパラメータとするテーブルで設定したものを用いる。
ステップ３では、燃焼室５の吸気弁閉時期ＩＶＣにおける温度（圧縮開始時期温度）Ｔｉｖｃ［Ｋ］を算出する。燃焼室５に流入するガスの温度は、燃焼室５に流入する新気と燃焼室５に残留する不活性ガスとが混じったガスの温度であり、燃焼室５に流入する新気の温度は吸気コレクタ２内の新気温度ＴＣＯＬに等しく、また燃焼室５内に残留する不活性ガスの温度は排気ポート部近傍の排気温度ＴＥＸＨで近似できるので、燃焼室５の吸気弁閉時期ＩＶＣにおける温度Ｔｉｖｃは吸気弁閉時期ＩＶＣになったタイミングでの、吸気コレクタ２内の新気温度ＴＣＯＬ、排気温度ＴＥＸＨ、燃焼室５内に残留する不活性ガスの割合である内部不活性ガス率ＭＲＥＳＲから次の（４６）式により求めることができる。
Ｔｉｖｃ＝ＴＥＸＨ×ＭＲＥＳＲ＋ＴＣＯＬ×（１−ＭＲＥＳＲ）
…（４６）
ステップ４では燃焼室５の吸気弁閉時期ＩＶＣにおける圧力（圧縮開始時期圧力）Ｐｉｖｃ［ｋＰａ］を算出する。すなわち、吸気弁閉時期ＩＶＣになったタイミングでのコレクタ内圧力ＰＣＯＬを吸気弁閉時期ＩＶＣにおける圧力Ｐｉｖｃとして取り込む。
…（４６）
ステップ４では燃焼室５の吸気弁閉時期ＩＶＣにおける圧力（圧縮開始時期圧力）Ｐｉｖｃ［ｋＰａ］を算出する。すなわち、吸気弁閉時期ＩＶＣになったタイミングでのコレクタ内圧力ＰＣＯＬを吸気弁閉時期ＩＶＣにおける圧力Ｐｉｖｃとして取り込む。
ステップ５では、燃料噴射パルス幅Ｔｉ（または燃料噴射量）から図１３を内容とするテーブルを検索することにより、トータルの燃焼による圧力上昇分Δｐｔｏｔａｌ［Ｐａ］を算出する。
ステップ６では、バルブリフト量Ｌｉｆｔから図４上段を内容とするテーブルを検索することにより、タンブル強度のバルブリフト量補正係数ｋＶＥＬ［１／ｄｅｇ］を算出する。すなわち、ＶＥＬ機構２８の非作動時には第１バルブリフト量補正係数ｋＶＥＬ１を、これに対してＶＥＬ機構２８の作動時には第２バルブリフト量補正係数ｋＶＥＬ２を、バルブリフト量補正係数ｋＶＥＬとして算出する。
ステップ７では、同じくバルブリフト量Ｌｉｆｔから図１０上段を内容とするテーブルを検索することにより、流速係数ｋ１ＶＥＬ［無名数］を算出する。すなわち、ＶＥＬ機構２８の非作動時には第１流速係数ｋ１ＶＥＬ１を、これに対してＶＥＬ機構２８の作動時には第２流速係数ｋ１ＶＥＬ２を、流速係数ｋ１ＶＥＬとして算出する。
ステップ８では、このうちの流速係数ｋ１ＶＥＬ、吸気弁閉時期ＩＶＣ、総ガス質量ＭＧＡＳ、エンジン回転速度Ｎｅを用いて、次の（４７）式により燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％［ｄｅｇＡＴＤＣ］を算出する。
θｍｂ１０％＝ｋ１ＶＥＬ・（ＩＶＣ−θｔ０）・ＭＧＡＳ＾ｃ１・Ｎｅ＾ｃ２
＋ｃ３ …（４７）
ただし、ｋ１ＶＥＬ：流速係数［無名数］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ＭＧＡＳ ：総ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ｃ１〜ｃ３：適合係数、
（４７）式は、上記（１４）式において、燃焼室内ガス質量ｍとして総ガス質量ＭＧＡＳを用いたもので、基本的に上記（１４）式と変わらない。
＋ｃ３ …（４７）
ただし、ｋ１ＶＥＬ：流速係数［無名数］、
ＩＶＣ ：吸気弁閉時期［ｄｅｇＡＴＤＣ］、
θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ＭＧＡＳ ：総ガス質量［ｇ］、
Ｎｅ ：エンジン回転速度［ｒｐｍ］、
ｃ１〜ｃ３：適合係数、
（４７）式は、上記（１４）式において、燃焼室内ガス質量ｍとして総ガス質量ＭＧＡＳを用いたもので、基本的に上記（１４）式と変わらない。
ここで、燃焼圧が最大となるクランク角θｐｍａｘは（４７）式右辺第１項、つまり次の（４７−１）式により与えていることになる。
θｐｍａｘ＝ｋ１ＶＥＬ・（ＩＶＣ−θｔ０）・ＭＧＡＳ＾ｃ１・Ｎｅ＾ｃ２
…（４７−１）
ステップ９では燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％から単純にクランク角で４ｄｅｇを差し引いた値を燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％［ｄｅｇＡＴＤＣ］として求める。これはθｍｂ２％よりθｍｂ１０％までのクランク角区間で燃焼質量割合が直線的に変化するとみなしてθｍｂ６％を求めるようにしたものある。
…（４７−１）
ステップ９では燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％から単純にクランク角で４ｄｅｇを差し引いた値を燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％［ｄｅｇＡＴＤＣ］として求める。これはθｍｂ２％よりθｍｂ１０％までのクランク角区間で燃焼質量割合が直線的に変化するとみなしてθｍｂ６％を求めるようにしたものある。
ステップ１０では、このθｍｂ６％から図１２を内容とするテーブルを検索することにより、燃焼室５のθｍｂ６％時における容積Ｖｍｂ６％を算出する。
ステップ１１、１２はステップ９、１０と同様である。 ステップ１１では燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％から単純にクランク角で８ｄｅｇを差し引いた値を燃焼質量割合が２％となるときのクランク角位置であるθｍｂ２％［ｄｅｇＡＴＤＣ］として求める。ステップ１２ではこのθｍｂ２％から図１２を内容とするテーブルを検索することにより、燃焼室５のθｍｂ２％時における容積Ｖｍｂ２％を算出する。
図１６は燃焼室５内の内部不活性ガス率ＭＲＥＳＲ［％］を算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。このフローは上記図１４のフローに先立って実行する。
ステップ２１ではエアフローメータ３２の出力と目標当量比Ｔｆｂｙａを読み込む。ステップ２２ではエアフロメータ３２の出力に基づいて、燃焼室５に流入する新気量（シリンダ新気量）ＭＡＣＹＬ［ｇ］を算出する。このシリンダ新気量ＭＡＣＹＬの算出方法については公知の方法を用いればよい（特開２００１−５００９１号公報参照）。
ステップ２３では、燃焼室５内の内部不活性ガス量ＭＲＥＳ［ｇ］を算出する。この内部不活性ガス量ＭＲＥＳの算出についても公知の方法を用いればよい（特開２００５−１７１８５６号公報参照）。
ステップ２４では、この内部不活性ガス量ＭＲＥＳ、シリンダ新気量ＭＡＣＹＬ、目標当量比Ｔｆｂｙａから次の（４８）式により燃焼室５の総ガス質量ＭＧＡＳ［ｇ］を算出する。
ＭＧＡＳ＝ＭＡＣＹＬ×(１＋Ｔｆｂｙａ／１４．７)＋ＭＲＥＳ
…（４８）
（４８）式右辺の括弧内の「１」は新気分、「Ｔｆｂｙａ／１４．７」は燃料分である。
…（４８）
（４８）式右辺の括弧内の「１」は新気分、「Ｔｆｂｙａ／１４．７」は燃料分である。
ステップ２５では、この総ガス質量ＭＧＡＳ、内部不活性ガス量ＭＲＥＳを用いて、次の（４９）式により内部不活性ガス率ＭＲＥＳＲ（燃焼室５内の総ガス質量に対する内部不活性ガス量の割合）［％］を算出する。
ＭＲＥＳＲ＝ＭＲＥＳ／ＭＧＡＳ …（４９）
図１７は燃焼時間ＢＴ［ｓｅｃ］を算出するためのもの、図１８は着火遅れ時間τ［ｓｅｃ］を算出するためのもので、それぞれ一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。図１７、図１８は図１４に続けて実行する。図１７、図１８はどちらを先に実行してもかまわない。
図１７は燃焼時間ＢＴ［ｓｅｃ］を算出するためのもの、図１８は着火遅れ時間τ［ｓｅｃ］を算出するためのもので、それぞれ一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。図１７、図１８は図１４に続けて実行する。図１７、図１８はどちらを先に実行してもかまわない。
ここでは図１７から先に説明すると、ステップ３１では、吸気弁閉時期ＩＶＣ［ｄｅｇＡＴＤＣ］、図１４のステップ２で算出されている燃焼室５の吸気弁閉時期における容積Ｖｉｖｃ［ｍ3］、図１４のステップ３で算出されている燃焼室５の吸気弁閉時期における温度Ｔｉｖｃ［Ｋ］、図１４のステップ４で算出されている燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ［Ｐａ］、図１４のステップ５で算出されているトータルの燃焼による圧力上昇分Δｐｃｔｏｔａｌ［ｋＰａ］、図１４のステップ６で算出されているタンブル強度のバルブリフト量補正係数ｋＶＥＬ［１／ｄｅｇ］、図１４のステップ１０で算出されている燃焼室５のθｍｂ６％時における容積Ｖｍｂ６％［ｄｅｇＡＴＤＣ］、図１６のステップ２４で算出されている総ガス質量ＭＧＡＳ［ｇ］、図１６のステップ２５で算出されている内部不活性ガス率ＭＲＥＳＲ［％］、燃料噴射パルス幅Ｔｉ［ｍｓｅｃ］、目標当量比Ｔｆｂｙａ［無名数］、エンジン回転速度Ｎｅ［ｒｐｍ］を読み込む。
ステップ３２では、燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ、燃焼室５の吸気弁閉時期における容積Ｖｉｖｃ、燃焼室５のθｍｂ６％時における容積Ｖｍｂ６％から、次の（５０）式により、容積変化による圧力変化分Δｐｖ［ｋＰ］を算出する。
Δｐｖ＝Ｐｉｖｃ・｛（Ｖｉｖｃ／Ｖｍｂ６％）＾κ−１｝ …（５０）
ただし、κ：比熱比（固定値で１．３〜１．４）、
（５０）式は上記（２８）式と同じものである。ここでの容積変化による圧力上昇分Δｐｖはθｍｂ２％よりθｍｂ１０％までのクランク角区間（燃焼期間）における圧力上昇分であるので、θｍｂ２％とθｍｂ１０％の平均値であるθｍｂ６％を用いている。
ただし、κ：比熱比（固定値で１．３〜１．４）、
（５０）式は上記（２８）式と同じものである。ここでの容積変化による圧力上昇分Δｐｖはθｍｂ２％よりθｍｂ１０％までのクランク角区間（燃焼期間）における圧力上昇分であるので、θｍｂ２％とθｍｂ１０％の平均値であるθｍｂ６％を用いている。
ステップ３３では、燃焼によるトータルの圧力上昇分Δｐｃｔｏｔａｌから、次の（５１）式により、燃焼による圧力上昇分Δｐｃ［ｋＰ］を算出する。
Δｐｃ＝０．０６×Δｐｃｔｏｔａｌ …（５１）
（５１）式は上記（２９）式において、燃焼質量割合ｘとして６％を用いたものである。
（５１）式は上記（２９）式において、燃焼質量割合ｘとして６％を用いたものである。
ステップ３４では、これら燃焼による圧力上昇分Δｐｃ、容積変化による圧力上昇分Δｐｖを燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃに加算して、つまり次の（５２）式により燃焼室５のθｍｂ６％時における未燃ガス圧力Ｐｍｂ６％［ｋＰ］を算出する。
Ｐｍｂ６％＝Ｐｉｖｃ＋Δｐｖ＋Δｐｃ …（５２）
ステップ３５では、この燃焼室５のθｍｂ６％時における未燃ガス圧力Ｐｍｂ６％、燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ、燃焼室５の吸気弁閉時期における温度Ｔｉｖｃを用いて、次の（５３）式により燃焼室５のθｍｂ６％時における未燃ガス温度Ｔｍｂ６％［Ｋ］を算出する。
ステップ３５では、この燃焼室５のθｍｂ６％時における未燃ガス圧力Ｐｍｂ６％、燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ、燃焼室５の吸気弁閉時期における温度Ｔｉｖｃを用いて、次の（５３）式により燃焼室５のθｍｂ６％時における未燃ガス温度Ｔｍｂ６％［Ｋ］を算出する。
Ｔｍｂ６％＝Ｔｉｖｃ・（Ｐｍｂ６％／Ｐｉｖｃ）＾｛（κ−１）／κ｝ …（５３）
ただし、κ：比熱比（固定値で１．３〜１．４）、
（５３）式は上記（２５）式において燃焼質量割合ｘが６％のときの値である。
ただし、κ：比熱比（固定値で１．３〜１．４）、
（５３）式は上記（２５）式において燃焼質量割合ｘが６％のときの値である。
ステップ３６では、目標当量比Ｔｆｂｙａから次の（５４）式〜（５６）式により標準状態での層流燃焼速度ＳＬ０［ｍ／ｓｅｃ］、係数ａ、ｂ［無名数］を求め、ステップ３７で、これら標準状態での層流燃焼速度ＳＬ０、係数ａ、ｂと、燃焼室５のθｍｂ６％時における未燃ガス温度Ｔｍｂ６％、燃焼室５のθｍｂ６％時における未燃ガス圧力Ｐｍｂ６％、内部不活性ガス率ＭＲＥＳＲを用いて、次の（５７）式によりθｍｂ６％時の層流燃焼速度ＳＬ［ｍ／ｓ］を算出する。
ＳＬ０＝（０．２６３２−０．８４７２／（Ｔｆｂｙａ−１．１３）＾２）
…（５４）
ａ＝２．１８−０．８０・（Ｔｆｂｙａ−１） …（５５）
ｂ＝−０．１６＋０．２２・（Ｔｆｂｙａ−１） …（５６）
ＳＬ＝ＳＬ０・（Ｔｍｂ６％／２９８）＾ａ・（Ｐｍｂ６％／１０１．３２５）＾ｂ
・（１−２．１×ＭＲＥＳＲ） …（５７）
（５４）式〜（５６）式は、上記（２０）式〜（２２）式において、当量比φとして目標当量比Ｔｆｂｙａを用いたもので、基本的に（２０）式〜（２２）式と変わらない。同様に、（５７）式は、上記（１７）式において、未燃ガス温度Ｔ、燃焼室内圧力ｐとしてＴｍｂ６％、Ｐｍｂ６％を用いたもので、基本的に（１７）式と変わらない。
…（５４）
ａ＝２．１８−０．８０・（Ｔｆｂｙａ−１） …（５５）
ｂ＝−０．１６＋０．２２・（Ｔｆｂｙａ−１） …（５６）
ＳＬ＝ＳＬ０・（Ｔｍｂ６％／２９８）＾ａ・（Ｐｍｂ６％／１０１．３２５）＾ｂ
・（１−２．１×ＭＲＥＳＲ） …（５７）
（５４）式〜（５６）式は、上記（２０）式〜（２２）式において、当量比φとして目標当量比Ｔｆｂｙａを用いたもので、基本的に（２０）式〜（２２）式と変わらない。同様に、（５７）式は、上記（１７）式において、未燃ガス温度Ｔ、燃焼室内圧力ｐとしてＴｍｂ６％、Ｐｍｂ６％を用いたもので、基本的に（１７）式と変わらない。
ステップ３８では、このθｍｂ６％時の層流燃焼速度ＳＬ、タンブル強度のバルブリフト量補正係数ｋＶＥＬ、吸気弁閉時期ＩＶＣ、総ガス質量ＭＧＡＳ、エンジン回転速度Ｎｅを用いて、次の（５８）式により、燃焼室内のガス流動としてのタンブル強度Ｉｔ［無名数］を算出する。
Ｉｔ＝ｋＶＥＬ・（ＩＶＣ−θｔ０）・ＭＧＡＳ＾ｃ４・Ｎｅ＾（ｃ５−１）
…（５８）
ただし、θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｃ４、ｃ５：適合係数、
ＶＴＣ機構２７の作動で吸気弁閉時期ＩＶＣが遅角側に遅らされるほど、つまりタンブル開始角を起点としてＩＶＣまでのクランク角区間（ＩＶＣ−θｔ０）が大きくなるほど実際のタンブル強度が大きくなるが、（５８）式よれば、吸気弁閉時期ＩＶＣが遅れるとき（ＩＶＣの値が大きくなるとき）タンブル強度Ｉｔも大きくなるのであり、算出値としてのＩｔは実際のタンブル強度と良く一致することとなる。
…（５８）
ただし、θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｃ４、ｃ５：適合係数、
ＶＴＣ機構２７の作動で吸気弁閉時期ＩＶＣが遅角側に遅らされるほど、つまりタンブル開始角を起点としてＩＶＣまでのクランク角区間（ＩＶＣ−θｔ０）が大きくなるほど実際のタンブル強度が大きくなるが、（５８）式よれば、吸気弁閉時期ＩＶＣが遅れるとき（ＩＶＣの値が大きくなるとき）タンブル強度Ｉｔも大きくなるのであり、算出値としてのＩｔは実際のタンブル強度と良く一致することとなる。
また、ＶＥＬ機構２８の作動で吸気弁のバルブリフト量が大きくなるほど実際のタンブル強度が大きくなるが、図４上段によればバルブリフト量Ｌｉｆｔが大きくなるときバルブリフト量補正係数ｋＶＥＬが大きくなり、（５８）式よれば、バルブリフト量補正係数ｋＶＥＬが大きくなるときタンブル強度Ｉｔも大きくなるのであり、算出値としてのＩｔは実際のタンブル強度と良く一致することとなる。
さらに、（５８）式よれば、総ガス質量ＭＧＡＳが大きくなるほどタンブル強度Ｉｔが大きくなる。
ステップ３９では、このタンブル強度Ｉｔと回転速度Ｎｅから次の（５９）式により燃焼室内ガスの乱れ強さｕ［ｍ／ｓｅｃ］を算出する。
ｕ＝Ｉｔ・Ｎｅ …（５９）
これは、タンブルが圧縮工程中に小渦に変換されることにより、燃焼室内ガスの乱れが強くなるので、タンブル強度Ｉｔが強いほど燃焼室内ガスの乱れが強くなると考え、タンブル強度Ｉｔ（燃焼室内のガス流動）に基づいて燃焼室内ガスの乱れ強さｕを推定するものである。
これは、タンブルが圧縮工程中に小渦に変換されることにより、燃焼室内ガスの乱れが強くなるので、タンブル強度Ｉｔが強いほど燃焼室内ガスの乱れが強くなると考え、タンブル強度Ｉｔ（燃焼室内のガス流動）に基づいて燃焼室内ガスの乱れ強さｕを推定するものである。
ステップ４０では、この燃焼室内ガスの乱れ強さｕから次の（６０）式により乱流燃焼速度ＳＴ［ｍ／ｓｅｃ］を算出する。
ＳＴ＝［ｕ／ＳＬ］＾ｃ６ …（６０）
ただし、ｃ６：適合係数、
ステップ４１では、このようにして求めた乱流燃焼速度ＳＴと、ステップ３７で算出している、燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％時の層流燃焼速度ＳＬとから、次の（６１）式により燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％時の燃焼速度Ｓｂ［ｍ／ｓ］を算出する。
ただし、ｃ６：適合係数、
ステップ４１では、このようにして求めた乱流燃焼速度ＳＴと、ステップ３７で算出している、燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％時の層流燃焼速度ＳＬとから、次の（６１）式により燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％時の燃焼速度Ｓｂ［ｍ／ｓ］を算出する。
Ｓｂ＝ＳＬ＋ＳＴ …（６１）
なお、ステップ３８〜４１に代えて、燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％時の燃焼速度Ｓｂ［ｍ／ｓ］を次の一つの式により算出してもかまわない。
なお、ステップ３８〜４１に代えて、燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％時の燃焼速度Ｓｂ［ｍ／ｓ］を次の一つの式により算出してもかまわない。
Ｓｂ＝ＳＬ＋［｛ｋＶＥＬ・（ＩＶＣ−θｔ０）・ＭＧＡＳ＾ｃ４・Ｎｅ＾ｃ５｝
／ＳＬ］＾ｃ６ …（６２）
ただし、θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｃ４、ｃ５、ｃ６：適合係数、
（６２）式は、上記（１６）式において、燃焼室内ガス質量ｍとして総ガス質量ＭＧＡＳを用いたもので、基本的に（１６）式と変わらない。
／ＳＬ］＾ｃ６ …（６２）
ただし、θｔ０ ：タンブル形成開始角［ｄｅｇＡＴＤＣ］、
ｃ４、ｃ５、ｃ６：適合係数、
（６２）式は、上記（１６）式において、燃焼室内ガス質量ｍとして総ガス質量ＭＧＡＳを用いたもので、基本的に（１６）式と変わらない。
ステップ４２では、この燃焼質量割合が６％となるときのクランク角位置であるθｍｂ６％時の燃焼速度Ｓｂ、燃焼室５のθｍｂ６％時におけ容積Ｖｍｂ６％を用いて、次の（６３）式により燃焼時間ＢＴ［ｓ］を算出する。
ＢＴ＝０．０７０７・Ｖｍｂ６％／（Ａｂ・Ｓｂ） …（６３）
ただし、Ａｂ：火炎表面積［ｍ2］、
（６３）式は、上記（１５）式において、燃焼室内容積ＶｃｙｌとしてＶｍｂ６％を用いたもので、基本的に（１５）式と変わらない。
ただし、Ａｂ：火炎表面積［ｍ2］、
（６３）式は、上記（１５）式において、燃焼室内容積ＶｃｙｌとしてＶｍｂ６％を用いたもので、基本的に（１５）式と変わらない。
図１８は着火遅れ時間τ［ｓｅｃ］を算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。
ステップ５１〜５７は図１６のステップ３１〜３７と同様である。すなわち、ステップ５１では、吸気弁閉時期ＩＶＣ［ｄｅｇＡＴＤＣ］、図１４のステップ２で算出されている燃焼室５の吸気弁閉時期における容積Ｖｉｖｃ［ｍ3］、図１４のステップ３で算出されている燃焼室５の吸気弁閉時期における温度Ｔｉｖｃ［Ｋ］、図１４のステップ４で算出されている燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ［Ｐａ］、図１４のステップ５で算出されているトータルの燃焼による圧力上昇分Δｐｃｔｏｔａｌ［Ｐａ］、図１４のステップ６で算出されているタンブル強度のバルブリフト量補正係数ｋＶＥＬ［１／ｄｅｇ］、図１４のステップ１２で算出されている燃焼室５のθｍｂ２％時における容積Ｖｍｂ２％［ｄｅｇＡＴＤＣ］、図１６のステップ２４で算出されている総ガス質量ＭＧＡＳ［ｇ］、図１６のステップ２５で算出されている内部不活性ガス率ＭＲＥＳＲ［％］、燃料噴射パルス幅Ｔｉ［ｍｓ］、目標当量比Ｔｆｂｙａ［無名数］、エンジン回転速度Ｎｅ［ｒｐｍ］を読み込む。
ステップ５２では、燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ、燃焼室５の吸気弁閉時期における容積Ｖｉｖｃ、燃焼室５のθｍｂ２％時における容積Ｖｍｂ２％から、次の（６４）式により、容積変化による圧力変化分Δｐｖ［ｋＰ］を算出する。
Δｐｖ＝Ｐｉｖｃ・｛（Ｖｉｖｃ／Ｖｍｂ２％）＾κ−１｝ …（６４）
ただし、κ：比熱比（固定値で１．３〜１．４）、
（６４）式は上記（２８）式と同じものである。ここでの容積変化による圧力上昇分Δｐｖはθｍｂ０％よりθｍｂ２％までのクランク角区間（着火遅れクランク角区間）における圧力上昇分であるが、θｍｂ０％とθｍｂ２％の平均値であるθｍｂ１％は用いず、着火遅れクランク角区間の終期の値であるθｍｂ２％を用いている。
ただし、κ：比熱比（固定値で１．３〜１．４）、
（６４）式は上記（２８）式と同じものである。ここでの容積変化による圧力上昇分Δｐｖはθｍｂ０％よりθｍｂ２％までのクランク角区間（着火遅れクランク角区間）における圧力上昇分であるが、θｍｂ０％とθｍｂ２％の平均値であるθｍｂ１％は用いず、着火遅れクランク角区間の終期の値であるθｍｂ２％を用いている。
ステップ５３では、燃焼によるトータルの圧力上昇分Δｐｃｔｏｔａｌから、次の（６５）式により、燃焼による圧力上昇分Δｐｃ［ｋＰ］を算出する。
Δｐｃ＝０．０２×Δｐｃｔｏｔａｌ …（６５）
（６５）式は上記の（２９）式において、燃焼質量割合ｘとして２％を用いたものである。
（６５）式は上記の（２９）式において、燃焼質量割合ｘとして２％を用いたものである。
ステップ５４では、これら燃焼による圧力上昇分Δｐｃ、容積変化による圧力上昇分Δｐｖを燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃに加算して、つまり次の（６６）式により燃焼室５のθｍｂ２％時における未燃ガス圧力Ｐｍｂ２％［ｋＰ］を算出する。
Ｐｍｂ２％＝Ｐｉｖｃ＋Δｐｖ＋Δｐｃ …（６６）
ステップ５５では、この燃焼室５のθｍｂ２％時における未燃ガス圧力Ｐｍｂ２％、燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ、燃焼室５の吸気弁閉時期における温度Ｔｉｖｃを用いて、次の（６７）式により燃焼室５のθｍｂ２％時における未燃ガス温度Ｔｍｂ２％［Ｋ］を算出する。
ステップ５５では、この燃焼室５のθｍｂ２％時における未燃ガス圧力Ｐｍｂ２％、燃焼室５の吸気弁閉時期における圧力Ｐｉｖｃ、燃焼室５の吸気弁閉時期における温度Ｔｉｖｃを用いて、次の（６７）式により燃焼室５のθｍｂ２％時における未燃ガス温度Ｔｍｂ２％［Ｋ］を算出する。
Ｔｍｂ２％＝Ｔｉｖｃ・（Ｐｍｂ２％／Ｐｉｖｃ）＾｛（κ−１）／κ｝ …（６７）
ただし、κ：比熱比（固定値で１．３〜１．４）、
（６７）式は上記（２５）式において燃焼質量割合ｘが２％のときの値である。
ただし、κ：比熱比（固定値で１．３〜１．４）、
（６７）式は上記（２５）式において燃焼質量割合ｘが２％のときの値である。
ステップ５６では、目標当量比Ｔｆｂｙａから次の（６８）式〜（７０）式により標準状態での層流燃焼速度ＳＬ０［ｍ／ｓｅｃ］、係数ａ、ｂを求め、ステップ５７で、これら標準状態での層流燃焼速度ＳＬ０、係数ａ、ｂと、燃焼室５のθｍｂ２％時における未燃ガス温度Ｔｍｂ２％、燃焼室５のθｍｂ２％時における未燃ガス圧力Ｐｍｂ２％、内部不活性ガス率ＭＲＥＳＲを用いて、次の（７１）式によりθｍｂ２％時の層流燃焼速度ＳＬ［ｍ／ｓ］を算出する。
ＳＬ０＝（０．２６３２−０．８４７２／（Ｔｆｂｙａ−１．１３）＾２）
…（６８）
ａ＝２．１８−０．８０・（Ｔｆｂｙａ−１） …（６９）
ｂ＝−０．１６＋０．２２・（Ｔｆｂｙａ−１） …（７０）
ＳＬ＝ＳＬ０・（Ｔｍｂ２％／２９８）＾ａ・（Ｐｍｂ２％／１０１．３２５）＾ｂ
・（１−２．１×ＭＲＥＳＲ） …（７１）
（６８）式〜（７０）式は、上記（２０）式〜（２２）式において、当量比φとして目標当量比Ｔｆｂｙａを用いたもので、基本的に（２０）式〜（２２）式と変わらない。同様に、（７１）式は、上記（１９）式において、未燃ガス温度Ｔ、燃焼室内圧力ｐとしてＴｍｂ２％、Ｐｍｂ２％を用いたもので、基本的に（１９）式と変わらない。
…（６８）
ａ＝２．１８−０．８０・（Ｔｆｂｙａ−１） …（６９）
ｂ＝−０．１６＋０．２２・（Ｔｆｂｙａ−１） …（７０）
ＳＬ＝ＳＬ０・（Ｔｍｂ２％／２９８）＾ａ・（Ｐｍｂ２％／１０１．３２５）＾ｂ
・（１−２．１×ＭＲＥＳＲ） …（７１）
（６８）式〜（７０）式は、上記（２０）式〜（２２）式において、当量比φとして目標当量比Ｔｆｂｙａを用いたもので、基本的に（２０）式〜（２２）式と変わらない。同様に、（７１）式は、上記（１９）式において、未燃ガス温度Ｔ、燃焼室内圧力ｐとしてＴｍｂ２％、Ｐｍｂ２％を用いたもので、基本的に（１９）式と変わらない。
ステップ５８では、このθｍｂ２％時の層流燃焼速度ＳＬ、燃焼室５のθｍｂ２％時における容積Ｖｍｂ２％、エンジン回転速度Ｎｅ、総ガス質量ＭＧＡＳを用いて、次の（７２）式により着火遅れ時間τ［ｓｅｃ］を算出する。
τ＝｛０．４９０・Ｖｍｂ２％＾（１／３）｝
／｛ｃ７・ＳＬ＋ｃ８・Ｎｅ・ＭＧＡＳ＾（１／２）｝
…（７２）
ただし、ｃ７、ｃ８：適合係数、
（７２）式は、上記（１８）式において、燃焼室内容積Ｖｃｙｌ、燃焼室内ガス質量ｍとして、それぞれ燃焼室５のθｍｂ２％時における容積Ｖｍｂ２％、総ガス質量ＭＧＡＳを用いたもので、基本的に（１８）式と変わらない。
／｛ｃ７・ＳＬ＋ｃ８・Ｎｅ・ＭＧＡＳ＾（１／２）｝
…（７２）
ただし、ｃ７、ｃ８：適合係数、
（７２）式は、上記（１８）式において、燃焼室内容積Ｖｃｙｌ、燃焼室内ガス質量ｍとして、それぞれ燃焼室５のθｍｂ２％時における容積Ｖｍｂ２％、総ガス質量ＭＧＡＳを用いたもので、基本的に（１８）式と変わらない。
図１９は基本点火時期ＭＢＴＣＡＬ［ｄｅｇＢＴＤＣ］を算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。図１７、図１８のうち遅く実行されるフローに続けて実行する。
ステップ６１では、図１４のステップ８で算出されている燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％［ｄｅｇＡＴＤＣ］、図１７のステップ４２で算出されている燃焼時間ＢＴ［ｓｅｃ］、図１８のステップ５８で算出されている着火遅れ時間τ［ｓｅｃ］、エンジン回転速度Ｎｅ［ｒｐｍ］を読み込む。
ステップ６２では、着火遅れ時間τ、燃焼期間ＢＴ、燃焼質量割合が１０％となるときのクランク角位置であるθｍｂ１０％、エンジン回転速度Ｎｅを用いて、次の（７３）式により基本点火時期ＭＢＴＣＡＬ［ｄｅｇＡＴＤＣ］を算出し、ステップ６３で次の（７４）式によりその絶対値を採り、マイナスの符号を付けたものを改めて基本点火時期ＭＢＴＣＡＬ［ｄｅｇＢＴＤＣ］として算出する。
ＭＢＴＣＡＬ＝（τ＋ＢＴ）・Ｎｅ・６−θｍｂ１０％ …（７３）
ＭＢＴＣＡＬ＝−｜ＭＢＴＣＡＬ｜ …（７４）
（７３）式、（７４）式は、上記（１３）式と同じ式である。
ＭＢＴＣＡＬ＝−｜ＭＢＴＣＡＬ｜ …（７４）
（７３）式、（７４）式は、上記（１３）式と同じ式である。
このようにして算出した基本点火時期ＭＢＴＣＡＬは、点火時期指令値として点火レジスタに移され、実際のクランク角がこの点火時期指令値と一致したタイミングでエンジンコントローラ３１より一次電流を遮断する点火信号が点火コイル１３に出力される。
ここで、本実施形態の作用効果を説明する。ここでは、上記〈１〉の燃焼速度の算出方法及び上記〈２〉の基準クランク角の設定方法による作用効果をまず説明し、その後で上記〈３〉の着火遅れ時間の算出方法による作用効果を説明する。
〔１〕上記〈１〉の燃焼速度の算出方法及び上記〈２〉の基準クランク角の設定方法による作用効果：
タンブル強度Ｉｔ（燃焼室内のガス流動）は、燃焼室内のさまざまな場所に分布する小渦の角運動量の総和で表されるので、例えばシリンダ（燃焼室５）の中心からｉ番目の小渦までの距離をｒi、ｉ番目の小渦のガス質量をｍi、ｉ番目の小渦のガス流速をｖiとすると、上記（３）式によりタンブル強度Ｉｔ（燃焼室内のガス流動）が表されるのであるが、燃焼室内の全ての小渦についてこれらの値ｍi、ｒi、ｖiをエンジンコントローラ３１においてオンラインで算出することは現在のところ不可能である。
〔１〕上記〈１〉の燃焼速度の算出方法及び上記〈２〉の基準クランク角の設定方法による作用効果：
タンブル強度Ｉｔ（燃焼室内のガス流動）は、燃焼室内のさまざまな場所に分布する小渦の角運動量の総和で表されるので、例えばシリンダ（燃焼室５）の中心からｉ番目の小渦までの距離をｒi、ｉ番目の小渦のガス質量をｍi、ｉ番目の小渦のガス流速をｖiとすると、上記（３）式によりタンブル強度Ｉｔ（燃焼室内のガス流動）が表されるのであるが、燃焼室内の全ての小渦についてこれらの値ｍi、ｒi、ｖiをエンジンコントローラ３１においてオンラインで算出することは現在のところ不可能である。
そこで、本発明者は、タンブルが圧縮工程中に小渦に変換されることにより、燃焼室内ガスの乱れｕが強くなるため、タンブル強度Ｉｔが強いほど燃焼室内ガスの乱れが強くなると考え、タンブル強度Ｉｔ（燃焼室内のガス流動）に基づいて燃焼室内ガスの乱れ強さｕを推定するものとした。すなわち、本実施形態によれば、タンブル強度Ｉｔ（燃焼室内のガス流動）に基づいて燃焼室内ガスの乱れ強さｕを算出し（図１７のステップ３８、３９参照）、この燃焼室内ガスの乱れ強さｕに基づいて乱流燃焼速度ＳＴを算出し（図１７のステップ４０参照）、この乱流燃焼速度ＳＴに基づいてＭＢＴの得られる基本点火時期を算出する（図１７のステップ４１、４２、図１９のステップ６１〜６３参照）ので、タンブル強度Ｉｔ（燃焼室内のガス流動）の影響を受けて乱流燃焼速度ＳＴが変動しても、その変動する乱流燃焼速度ＳＴを精度よく算出できることから、基本点火時期ＭＢＴＣＡＬ（ＭＢＴの得られる基本点火時期）の予測精度を向上させることができる。
上記のように、燃焼室内の全ての小渦についてｍi、ｒi、ｖiをエンジンコントローラ３１においてオンラインで算出することは現在のところ不可能であるので、タンブル強度Ｉｔを求めるには何らかの近似が必要である。そこで、本発明者は燃焼室内のガス流動を新たに吸気弁１５のリフト特性とエンジン回転速度Ｎｅとに基づいて推定するものとした。ここで、吸気弁１５のリフト特性とは、吸気弁のバルブリフト量（あるいは吸気弁のバルブ作動角）または吸気弁の開閉時期のことである。
すなわち、本実施形態によれば、タンブル強度Ｉｔ（燃焼室内のガス流動）を吸気弁１５のバルブリフト量Ｌｉｆｔ及び吸気弁閉時期ＩＶＣ（リフト特性）とエンジン回転速度Ｎｅとに基づいて推定するので（図１４のステップ６、図１７のステップ３８参照）、簡易な構成でタンブル強度Ｉｔ（燃焼室内のガス流動）を精度よく算出できることになった。
ＶＥＬ機構２８を備えるときには、ＶＥＬ機構２８の作動状態で吸気弁１５のリフト量が小さいときと、ＶＥＬ機構２８の非作動状態で吸気弁１５のリフト量が大きいときでは燃焼室内のガス流動の向きが変わるため（図２（ｆ）と図２（ｇ）を参照）、燃焼室内ガスの乱れ強さｕも変わってしまうのであるが、本実施形態によれば、乱流燃焼速度ＳＴを少なくとも吸気弁１５のバルブリフト量Ｌｉｆｔ（リフト特性）に基づいて算出するので（図１４のステップ６、図１７のステップ３８、３９、４０参照）、ＶＥＬ機構２８の作動、非作動に伴う吸気弁１５のバルブリフト量Ｌｉｆｔの大小に関係なく乱流燃焼速度ＳＴの推定精度を向上できる。
ＶＴＣ機構２７の働きにより、吸気弁閉時期ＩＶＣが遅いほど燃焼室内ガスの乱れを形成する期間が長くなり、燃焼室内ガスの乱れ強さｕが強くなることに対応し、本実施形態によれば、乱流燃焼速度ＳＴを少なくとも吸気弁閉時期ＩＶＣに基づいて算出するので（図１７のステップ３８、３９、４０参照）、吸気弁閉時期ＩＶＣの進角、遅角に関係なく乱流燃焼速度ＳＴの推定精度を向上できる。
燃焼室５の総ガス質量ＭＧＡＳ（燃焼室内のガス量）が多いほどガスの持つ運動エネルギーが多くなり燃焼室内ガスの乱れ強さｕも強くなり燃焼が早くなることに対応し、本実施形態によれば、乱流燃焼速度ＳＴを燃焼室５の総ガス質量ＭＧＡＳ（燃焼室内のガス量）に基づいても算出するので（図１７のステップ３８、３９、４０参照）、燃焼室５の総ガス質量ＭＧＡＳ（燃焼室内のガス量）の多少に関係なく乱流燃焼速度ＳＴの推定精度を向上できる。
燃焼室内ガスの乱れが少ない状態では層流燃焼速度ＳＬが乱流燃焼速度ＳＴに与える影響が大きいことに対応し、本実施形態によれば、層流燃焼速度ＳＴを算出する層流燃焼速度算出処理手順を含み（図１７のステップ３７参照）、乱流燃焼速度ＳＴをこの層流燃焼速度ＳＬに基づいても算出するので（図１７のステップ４０参照）、燃焼室内ガスの乱れが少ない状態においても乱流燃焼速度ＳＴの推定精度を向上できる。
ところで、実施形態にいう吸気弁閉時期ＩＶＣは、実際に吸気弁１５が閉じる瞬間のことである。実際には、吸気管内圧力の脈動による影響で、吸気弁閉時期ＩＶＣよりも早く燃焼室内圧力が上昇し始めると、吸気弁は開いているのにそれ以上空気が燃焼室５へと入らないことがある。こうした場合には、燃焼室内のガス流動（燃焼室内ガスの乱れ）に影響するのは、吸気弁閉時期ＩＶＣではなく、吸気弁閉時期ＩＶＣの前に燃焼室内圧力が実際に上昇し始める時期である。ここで、吸気弁閉時期ＩＶＣの前に燃焼室内圧力が実際に上昇し始める時期を「圧縮行程開始時期」として新たに定義すれば、吸気弁閉時期ＩＶＣと圧縮行程開始時期とは異なる値となる。
そこで、吸気弁閉時期ＩＶＣに代え、この新たに定義した圧縮行程開始時期に基づいて乱流燃焼速度ＳＴを算出することで、吸気管内圧力の脈動による影響により、吸気弁閉時期ＩＶＣよりも早く燃焼室内圧力が上昇し始め、吸気弁１５は開いているのにそれ以上空気が燃焼室内に入らない場合においても、乱流燃焼速度ＳＴの推定精度を向上できることとなる。
圧縮行程開始時期を求める方法には、次のように検出する方法と推定する方法の２つが考えられる。
検出する方法：燃焼室に臨んで燃焼室内の圧力を検出するセンサ（圧縮行程開始時期検出装置）を設けておき、このセンサにより検出される燃料室内圧力に基づいて、吸気弁閉時期付近で、燃焼室内圧力が上昇を開始するときのクランク角を圧縮行程開始時期として計測（検出）する。
推定する方法：吸気弁閉時期ＩＶＣ［ｄｅｇＡＴＤＣ］と圧縮行程開始時期［ｄｅｇＡＴＤＣ］の間のクランク角度差ＤＣＡ［ｄｅｇ］はエンジン回転速度Ｎｅの感度を持っているので、次の式により圧縮行程開始時期を算出（推定）することができる。
圧縮行程開始時期＝ＩＶＣ−ＤＣＡ …（７５）
（７５）式のクランク角度差ＤＣＡは、そのときのエンジン回転速度Ｎｅから所定のテーブル（エンジン回転速度Ｎｅ大きくなるほどクランク角度差ＤＣＡが大きくなる特性）を検索することにより求める。
〔２〕上記〈３〉の着火遅れ時間の算出方法による作用効果：
図２０はコンベンショナルエンジンに適用している先行点火時期制御方法における点火無駄時間相当クランク角ＩＧＮＤＥＡＤと、本実施形態の変更後点火時期制御方法における着火遅れクランク角区間（＝τ・Ｎｅ・６）の違いを示している。なお、図２０は点火タイミング直後の燃焼質量割合の変化を拡大して示すもので、縦軸、横軸の各スケールは図８と同じでない。
（７５）式のクランク角度差ＤＣＡは、そのときのエンジン回転速度Ｎｅから所定のテーブル（エンジン回転速度Ｎｅ大きくなるほどクランク角度差ＤＣＡが大きくなる特性）を検索することにより求める。
〔２〕上記〈３〉の着火遅れ時間の算出方法による作用効果：
図２０はコンベンショナルエンジンに適用している先行点火時期制御方法における点火無駄時間相当クランク角ＩＧＮＤＥＡＤと、本実施形態の変更後点火時期制御方法における着火遅れクランク角区間（＝τ・Ｎｅ・６）の違いを示している。なお、図２０は点火タイミング直後の燃焼質量割合の変化を拡大して示すもので、縦軸、横軸の各スケールは図８と同じでない。
図２０より、本実施形態における着火遅れクランク角区間（＝τ・Ｎｅ・６）は、物理的な意味での着火遅れクランク角区間（つまり先行点火時期制御方法における点火無駄時間相当クランク角ＩＧＮＤＥＡＤ）ではなく、点火タイミングから燃料が２％燃えるまでのクランク角区間である（第１の仮定）。すなわち、火炎形成開始直後に燃焼室内圧力の上昇（変化）が十分となれば、現状の燃焼解析装置でも燃焼質量割合の計測結果が精度の良いものとなることに着目し、本実施形態では燃焼室内圧力の上昇（変化）が十分となる初期燃焼時（火炎形成開始時期よりも時間的に後のタイミング）の燃焼質量割合を基準燃焼質量割合として定めている。この基準燃焼質量割合は本実施形態では２％である。そして、点火時期から燃焼質量割合がこの２％に上昇するまでの時間を計測するのであれば、この時間計測は現状の燃焼解析装置によっても精度の良いものとなる。
しかも、点火タイミングの瞬間から火炎が形成されるとみなし（第２の仮定）、点火タイミング直後の火炎速度と基準燃焼質量割合になったときの火炎半径とに基づいて着火遅れ時間τを算出している。
本実施形態ではこうした２つの仮定を入れることによって、実機により検証してみたところ、基本点火時期ＭＢＴＣＡＬの算出に十分な精度が得られることを確認している。
このように、本実施形態（請求項１、６に記載の発明）によれば、点火タイミングから燃焼質量割合が０％より上昇を開始するタイミングまでの時間（コンベンショナルエンジンに適用している先行点火時期制御方法における点火無駄時間ＤＥＡＤＴＩＭＥ）ではなく、点火時期から燃焼質量割合が初期燃焼時の燃焼質量割合である基準燃焼質量割合に上昇するまでの時間である着火遅れ時間τを算出（推定）し、この着火遅れ時間τに基づいて基本点火時期ＭＢＴＣＡＬ（ＭＢＴの得られる基本点火時期）を算出するので、点火タイミングから燃焼質量割合が０％より上昇を開始するタイミングまでの時間に基づいて基本点火時期ＭＢＴＣＡＬを算出する場合よりも、基本点火時期ＭＢＴＣＡＬの予測精度を向上させることができる。
点火時期から燃焼質量割合が基準燃焼割合に上昇するまでの時間は、燃焼シミュレータを用いれば算出可能であるのに対して、車載のエンジンコントローラ３１で算出するには演算時間がかかりすぎることとなる。しかしながら、本実施形態（請求項１、６に記載の発明）によれば、点火タイミングの瞬間から火炎が形成されるとみなし、点火タイミング直後の火炎速度ＳＦＬＡＭＥと基準燃焼質量割合になったときの火炎半径Ｒkernelとに基づいて、着火遅れ時間τ（点火時期から燃焼質量割合が基準燃焼質量割合に上昇するまでの時間）を算出するので（上記（追３）式を参照）、演算負荷が飛躍的に低減することとなり、車載のエンジンコントローラ３１上でも、演算時間がかかりすぎることなく着火遅れ時間τを算出（演算）させることができる。
点火タイミング直後の火炎速度は層流燃焼速度と燃焼室内ガスの運動エネルギーの影響を多く受けることから、本実施形態（請求項１、６に記載の発明）によれば、火炎速度ＳＦＬＡＭＥを、少なくとも層流燃焼速度ＳＬと燃焼室内ガスの運動エネルギーｋとに基づいて算出するようにしているので（上記（追４）式を参照）、点火タイミング直後の層流燃焼速度ＳＬや燃焼室内ガスの運動エネルギーｋが相違しても、点火タイミング直後の火炎速度ＳＦＬＡＭＥを精度良く算出することができる。
本実施形態（請求項１、６に記載の発明）によれば、点火タイミング直後の火炎速度ＳＦＬＡＭＥを、さらに燃焼室内ガスの熱膨張率（Ｔａｄ／Ｔ）に基づいても算出するので（上記（追４）式を参照）、さらに点火タイミング直後の火炎速度ＳＦＬＡＭＥを精度良く推定できる。
燃焼室内ガスの運動エネルギーを演算するには複雑なガス流動シュミレーションを行う必要があり、車載のエンジンコントローラ３１ではこの演算は演算時間上不可能である。しかしながら、本実施形態（請求項２、７に記載の発明）によれば、燃焼室内ガスの運動エネルギーｋを、少なくとも燃焼室内ガス質量ｍとピストン速度に基づいて算出するので（上記（１２Ｂ）式を参照）、点火タイミング直後の燃焼室内ガスの運動エネルギーの演算時間が低減することとなり、車載のエンジンコントローラ３１上でも点火タイミング直後の燃焼室内ガスの運動エネルギーを算出（演算）させることができる。
基準燃焼質量割合になったときの火炎半径を正確に演算（推定）するには３次元の燃焼シュミレーションを行う必要があり、車載のエンジンコントローラ３１ではこの演算は演算時間上不可能である。しかしながら、本実施形態（請求項１、６に記載の発明）によれば、基準燃焼質量割合になったときの火炎半径Ｒkernelを基準燃焼質量割合になったｔきの燃焼室容積Ｖｃｙｌとに基づいて算出するので（上記（追５）式、（１０）式を参照）、車載のエンジンコントローラ３１上において大幅に精度を落とすことなく基準燃焼質量割合になったときの火炎半径Ｒkernelを演算することが可能になり、高価なＣＰＵを使う必要がなくなり、コストを低減できる。
点火タイミング直後には、火炎の形状は球形に近いため、本実施形態（請求項５、１０に記載の発明）によれば、火炎の形状を球形として基準燃焼質量割合での火炎半径Ｒkernelを算出する（上記（追５）式、（１１ａ）式を参照）。これにより、３次元の燃焼シュミレーションを使わずとも車載のエンジンコントローラ３１上において精度よく基準燃焼質量割合になったときの火炎半径Ｒkernelを算出（推定）できることとなり、高価なＣＰＵを使う必要が無くなり、コストを低減できる。
図２１は第２実施形態の着火遅れ時間τ［ｓｅｃ］を算出するためのもので、第１実施形態の図１８と置き換わるものである。図２１も一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。図２１において図１８と同一の部分には同一のステップ番号を付けている。
第２実施形態では、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳの算出方法が第１実施形態と相違するものである。このため、図２１において第１実施形態の図１８と相違する部分はステップ７１、７２である。すなわち、ステップ７１では、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳを算出する。この平均ガス流速ＶＧＡＳの算出については、図２２のフローにより説明する。
図２２（図２１のステップ７１のサブルーチン）においてステップ８１では充填効率ＩＴＡＣ［無名数］、エンジン回転速度Ｎｅ［ｒｐｍ］、吸気弁１５のバルブリフト量Ｌｉｆｔ［ｍ］、吸気弁閉時期［ｄｅｇＡＴＤＣ］を読み込む。
ここで、充填効率ＩＴＡＣの算出方法としては、公知の算出方法を用いればよい。また、エアフローメータ３２により検出される吸入空気流量を充填効率ＩＴＡＣとして用いてもかまわない。
ステップ８２では充填効率ＩＴＡＣとエンジン回転速度Ｎｅから図２３を内容とするマップを検索することにより、着火遅れ時間τにおける燃焼室内ガスの基本平均ガス流速ＶＧＡＳ０［ｍ／ｓｅｃ］を算出する。図２３に示したように、基本平均ガス流速ＶＧＡＳ０は充填効率ＩＴＡＣが小さくかつ高回転速度の領域において最も大きくなり、この領域より離れるほど小さくなる値である。
ステップ８３では吸気弁１５のバルブリフト量Ｌｉｆｔから図２４上段を内容とするテーブルを検索することにより、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速のバルブリフト量補正係数ｋ２ＶＥＬ［無名数］を算出する。すなわち、ＶＥＬ機構２８の非作動時には第１バルブリフト量補正係数ｋ２ＶＥＬ１を、これに対してＶＥＬ機構２８の作動時には第２バルブリフト量補正係数ｋ２ＶＥＬ２を、燃焼室内ガスの平均ガス流速のバルブリフト量補正係数ｋ２ＶＥＬとして算出する。このように、ＶＥＬ機構２８の作動時のほうを、ＶＥＬ機構２８の非作動時より燃焼室内ガスの平均ガス流速のバルブリフト量補正係数ｋ２ＶＥＬを小さくしているのは、吸気弁１５のバルブリフト量が小さくなるＶＥＬ機構２８の作動時のほうが燃焼室内ガスの平均ガス流速が小さくなるためである。
なお、図２４下段に示したように、ＶＥＬ機構２８の非作動時に第１バルブリフト量Ｌｉｆｔ１が、これに対してＶＥＬ機構２８の作動時に第２バルブリフト量Ｌｉｆｔ２がバルブリフト量Ｌｉｆｔになる。
ステップ８４では、吸気弁閉時期ＩＶＣから図２５を内容とするテーブルを検索することにより、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速の吸気弁閉時期補正係数ｋＩＶＣ［無名数］を算出する。図２５に示したように、燃焼室内ガスの平均ガス流速の吸気弁閉時期補正係数ｋＩＶＣは吸気弁閉時期ＩＶＣが遅くなるほど大きくなる値である。これは、吸気弁閉時期ＩＶＣが遅くなるほど燃焼室内ガスの平均ガス流速が大きくなるためである。
ステップ８５では、このようにして算出した基本平均ガス流速ＶＧＡＳ０、バルブリフト量補正係数ｋ２ＶＥＬ、吸気弁閉時期補正係数ｋＩＶＣの３つの値を乗算することにより、つまり次式により、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳ［ｍ／ｓｅｃ］を算出する。
ＶＧＡＳ＝ＶＧＡＳ０・ｋ２ＶＥＬ・ｋＩＶＣ …（７６）
（７６）式は次のようにして導いたものである。すなわち、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳは、ピストン速度、吸気ポート４と燃焼室５との圧力差（圧力比）、吸気弁開口面積、吸気弁１５のバルブリフト量Ｌｉｆｔ、吸気弁閉時期ＩＶＣにより変化する。これを式に表せば、次式になる。
（７６）式は次のようにして導いたものである。すなわち、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳは、ピストン速度、吸気ポート４と燃焼室５との圧力差（圧力比）、吸気弁開口面積、吸気弁１５のバルブリフト量Ｌｉｆｔ、吸気弁閉時期ＩＶＣにより変化する。これを式に表せば、次式になる。
ＶＧＡＳ＝ｇ（ピストン速度、吸気ポート圧力、燃焼室圧力、吸気弁開口面積）
…（７７）
ただし、ｇ：関数、
ここで、ピストン速度はエンジン回転速度Ｎｅにより決まる。吸気ポート４圧力は充填効率ＩＴＡＣと吸気弁１５のバルブリフト量Ｌｉｆｔで決まる。燃焼室内圧力は充填効率ＩＴＡＣで決まる。吸気弁開口面積は吸気弁１５のバルブリフト量Ｌｉｆｔで決まる。従って、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳは、（７７）式に代えて、次の（７８）式のようにエンジン回転速度Ｎｅ、充填効率ＩＴＡＣ、吸気弁１５のバルブリフト量Ｌｉｆｔ、吸気弁閉時期ＩＶＣの関数ｆ１として与えればよいことになる。
…（７７）
ただし、ｇ：関数、
ここで、ピストン速度はエンジン回転速度Ｎｅにより決まる。吸気ポート４圧力は充填効率ＩＴＡＣと吸気弁１５のバルブリフト量Ｌｉｆｔで決まる。燃焼室内圧力は充填効率ＩＴＡＣで決まる。吸気弁開口面積は吸気弁１５のバルブリフト量Ｌｉｆｔで決まる。従って、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳは、（７７）式に代えて、次の（７８）式のようにエンジン回転速度Ｎｅ、充填効率ＩＴＡＣ、吸気弁１５のバルブリフト量Ｌｉｆｔ、吸気弁閉時期ＩＶＣの関数ｆ１として与えればよいことになる。
ＶＧＡＳ＝ｇ（Ｎｅ、ＩＴＡＣ、Ｌｉｆｔ、ＩＶＣ） …（７８）
ただし、ｇ：関数、
しかしながら、この（７８）式だと燃焼室内ガスの平均ガス流速ＶＧＡＳが４つの変数をパラメータとすることになり適合工数が膨大なものとなってしまうので、（７８）式を次の（７９）式で近似する。
ただし、ｇ：関数、
しかしながら、この（７８）式だと燃焼室内ガスの平均ガス流速ＶＧＡＳが４つの変数をパラメータとすることになり適合工数が膨大なものとなってしまうので、（７８）式を次の（７９）式で近似する。
ＶＧＡＳ＝ｈ（ＩＴＡＣ、Ｎｅ）×ｊ（Ｌｉｆｔ）×ｎ（ＩＶＣ）
…（７９）
ただし、ｈ：ＩＴＡＣとＮｅを変数とする関数、
ｊ：Ｌｉｆｔを変数とする関数、
ｎ：ＩＶＣを変数とする関数、
（７９）式は、１つの関数ｇを３つの関数ｈ、ｊ、ｎの積で表したものである。（７９）式の第１の関数ｈ（ＩＴＡＣ、Ｎｅ）を燃焼室内ガスの基本平均ガス流速ＶＧＡＳ０、第２の関数ｊ（Ｌｉｆｔ）を燃焼室内ガスの平均ガス流速のバルブリフト量補正係数ｋ２ＶＥＬ、第３の関数ｎ（ＩＶＣ）を燃焼室内ガスの平均ガス流速の吸気弁閉時期補正係数ｋＩＶＣとおけば、（７９）式は次のようになり、上記（７６）式が得られることとなる。
…（７９）
ただし、ｈ：ＩＴＡＣとＮｅを変数とする関数、
ｊ：Ｌｉｆｔを変数とする関数、
ｎ：ＩＶＣを変数とする関数、
（７９）式は、１つの関数ｇを３つの関数ｈ、ｊ、ｎの積で表したものである。（７９）式の第１の関数ｈ（ＩＴＡＣ、Ｎｅ）を燃焼室内ガスの基本平均ガス流速ＶＧＡＳ０、第２の関数ｊ（Ｌｉｆｔ）を燃焼室内ガスの平均ガス流速のバルブリフト量補正係数ｋ２ＶＥＬ、第３の関数ｎ（ＩＶＣ）を燃焼室内ガスの平均ガス流速の吸気弁閉時期補正係数ｋＩＶＣとおけば、（７９）式は次のようになり、上記（７６）式が得られることとなる。
ＶＧＡＳ＝ｈ（ＩＴＡＣ、Ｎｅ）×ｊ（Ｌｉｆｔ）×ｎ（ＩＶＣ）
＝ＶＧＡＳ０・ｋ２ＶＥＬ・ｋＩＶＣ
上記図２３、図２４上段、図２５の各特性は、実機において計測し求めておくかまたはガス流動シミュレーションを行って求めておく。
＝ＶＧＡＳ０・ｋ２ＶＥＬ・ｋＩＶＣ
上記図２３、図２４上段、図２５の各特性は、実機において計測し求めておくかまたはガス流動シミュレーションを行って求めておく。
このようにして着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳを算出したら図２１に戻りステップ７２で次式により着火遅れ時間τ［ｓｅｃ］を算出する。
τ＝｛０．４９０・Ｖｍｂ２％＾（１／３）｝
／｛ｃ７・ＳＬ＋ｃ９・ＶＧＡＳ・ＭＧＡＳ＾（１／２）｝
…（８０）
ただし、ｃ７、ｃ９：適合係数、
ここで、第１実施形態と相違する第２実施形形態に固有の作用効果を説明する。
／｛ｃ７・ＳＬ＋ｃ９・ＶＧＡＳ・ＭＧＡＳ＾（１／２）｝
…（８０）
ただし、ｃ７、ｃ９：適合係数、
ここで、第１実施形態と相違する第２実施形形態に固有の作用効果を説明する。
吸気弁開口面積の大小、吸気ポート４と燃焼室５の圧力差（圧力比）の大小によって点火タイミング直後の燃焼室内ガス流速が大きく相違するため、これらの影響を考慮しないとすれば、点火タイミング直後の燃焼室内ガスの運動エネルギーの算出精度が低下してしまうのであるが、本実施形態（請求項３、８に記載の発明）によれば、燃焼室内ガスの運動エネルギーｋを、少なくとも、吸気ポート圧力と、燃焼室内圧力と、吸気弁開口面積とに基づいて算出するので（上記（７７）式を参照）、吸気弁開口面積や吸気ポート４と燃焼室５の圧力差（圧力比）が相違しても、点火タイミング直後の燃焼室内ガスの運動エネルギーｋを精度良く算出することができる。
吸気弁１５のプロファイル、具体的には吸気弁１５のバルブリフト量Ｌｉｆｔや吸気弁閉時期ＩＶＣに応じて点火タイミング直後の燃焼室内ガスの流れる方向が変わる。流れの方向によっては、燃焼室内で流れが打ち消され点火タイミング直後の燃焼室内ガスの流速が小さくなるケースがあり、点火タイミング直後の燃焼室内ガスの運動エネルギーに影響する。本実施形態（請求項４、９に記載の発明）によれば、点火タイミング直後の燃焼室内ガスの運動エネルギーｋを、少なくとも吸気弁１５のプロファイル（吸気弁１５のバルブリフト量Ｌｉｆｔと吸気弁閉時期ＩＶＣの少なくとも一方）に基づいて算出するので（上記（７９）式を参照）、吸気弁１５のプロファイルが相違しても、点火タイミング直後の燃焼室内ガスの運動エネルギーｋを精度良く算出することができる。
第２実施形態に対して、さらに、スワールコントロールバルブ、タンブルコントロールバルブのように、燃焼室５への吸気流動を変化させ得る吸気流動可変装置を備えさせることが考えられる。
このときには、タンブルコントロールバルブ（スワールコントロールバルブ）の作動時と非作動時とで、着火遅れ時間τにおける燃焼室内ガスの平均ガス流速ＶＧＡＳが相違してくるので、図２３を内容とする基本平均ガス流速ＶＧＡＳのマップ、図２４上段を内容とするバルブリフト量補正係数ｋ２ＶＥＬのテーブル、図２５を内容とする吸気弁閉時期補正係数ｋＩＶＣのテーブルは、タンブルコントロールバルブ（スワールコントロールバルブ）の作動時と、タンブルコントロールバルブ（スワールコントロールバルブ）の非作動時とに別個に備えさせる。例えば、タンブルコントロールバルブ（スワールコントロールバルブ）の作動時のほうが、タンブルコントロールバルブ（スワールコントロールバルブ）の非作動時よりも燃焼室内のガス流速が小さくなるため、図２３において同じ充填効率ＩＴＡＣとエンジン回転速度Ｎｅでも、タンブルコントロールバルブ（スワールコントロールバルブ）の作動時のほうが、タンブルコントロールバルブ（スワールコントロールバルブ）の非作動時よりも基本平均ガス流速ＶＧＡＳ０の値が小さくなるように設定する。同様にして、図２４上段において同じバルブリフト量Ｌｉｆｔでも、タンブルコントロールバルブ（スワールコントロールバルブ）の作動時のほうが、タンブルコントロールバルブ（スワールコントロールバルブ）の非作動時よりもバルブリフト量補正係数ｋ２ＶＥＬの値が小さくなるように、また図２５において同じ吸気弁閉時期ＩＶＣでも、タンブルコントロールバルブ（スワールコントロールバルブ）の作動時のほうが、タンブルコントロールバルブ（スワールコントロールバルブ）の非作動時よりも吸気弁閉時期補正係数ｋＩＶＣの値が小さくなるように設定する。
吸気流動可変装置はタンブルコントロールバルブ、スワールコントロールバルブに限られるものでなく、吸気管長を可変に制御可能な装置（吸気可変装置）が含まれる。この吸気可変装置を備える場合には、吸気可変装置の作動時と非作動時とで燃焼室内ガスのガス流速がどうなるかに基づいて、基本平均ガス流速ＶＧＡＳ０のマップ、バルブリフト量補正係数ｋ２ＶＥＬのテーブル、吸気弁閉時期補正係数ｋＩＶＣのテーブルを個別に設定してやればよい。例えば、吸気可変装置の作動時のほうが、吸気可変装置の非作動時よりも燃焼室内のガス流速が小さくなるのであれば、図２３において同じ充填効率ＩＴＡＣとエンジン回転速度Ｎｅでも、吸気可変装置の作動時のほうが、吸気可変装置の非作動時よりも基本平均ガス流速ＶＧＡＳ０の値が小さくなるように設定する。同様にして、図２４上段において同じバルブリフト量Ｌｉｆｔでも、吸気可変装置の作動時のほうが、吸気可変装置の非作動時よりもバルブリフト量補正係数ｋ２ＶＥＬの値が小さくなるように、また図２５において同じ吸気弁閉時期ＩＶＣでも、吸気可変装置の作動時のほうが、吸気可変装置の非作動時よりも吸気弁閉時期補正係数ｋＩＶＣの値が小さくなるように設定する。
実施形態では、ＶＥＬ機構２８及びＶＴＣ機構２７を備える場合で説明したが、本発明は、ＶＥＬ機構２８、ＶＴＣ機構２７のいずれかを備える場合あるいはいずれも備えない場合にも適用がある。ＶＥＬ機構２８、ＶＴＣ機構２７のいずれも備えない場合には、バルブリフト量補正係数ｋＶＥＬ及び吸気弁閉時期ＶＴＣが固定値となるだけある。
請求項１に記載の発明において、時間算出処理手順は図１８のステップ５８により、基本点火時期算出処理手順は図１９のステップ６１、６２、６３によりそれぞれ果たされている。
請求項１０に記載の発明において、時間算出手段の機能は図１８のステップ５８により、基本点火時期算出手段の機能は図１９のステップ６１、６２、６３によりそれぞれ果たされている。
１ エンジン
５ 燃焼室
１１ 点火装置
１５ 吸気弁
２１ 燃料インジェクタ
２７ ＶＴＣ機構（可変動弁装置）
２８ ＶＥＬ機構（可変動弁装置）
３１ エンジンコントローラ
３３、３４ クランク角センサ
５ 燃焼室
１１ 点火装置
１５ 吸気弁
２１ 燃料インジェクタ
２７ ＶＴＣ機構（可変動弁装置）
２８ ＶＥＬ機構（可変動弁装置）
３１ エンジンコントローラ
３３、３４ クランク角センサ
Claims (10)
- 点火の瞬間から火炎が形成されるとみなし、この点火直後の火炎の速度を、層流燃焼速度と燃焼室内ガスの運動エネルギーと燃焼室内ガスの熱膨張率に基づいて算出する火炎速度算出処理手順と、
初期燃焼時の燃焼質量割合を基準燃焼質量割合として、この基準燃焼質量割合になったときの燃焼室内容積に基づいて基準燃焼質量割合になったときの火炎半径を算出する火炎半径算出処理手順と、
点火から燃焼質量割合が前記基準燃焼質量割合に上昇するまでの時間を着火遅れ時間とみなして、この着火遅れ時間を前記点火直後の火炎速度と前記基準燃焼質量割合になったときの火炎半径とに基づいて算出する着火遅れ時間算出処理手順と、
燃焼質量割合が前記基準燃焼質量割合から基準クランク角での燃焼質量割合になるまでの時間である燃焼時間を算出する燃焼時間算出処理手順と、
これら着火遅れ時間と燃焼時間とを加算した値をクランク角区間に変換するクランク角区間変換処理処理手順と、
前記基準クランク角よりこの変換されたクランク角区間だけ進角側の値をＭＢＴの得られる基本点火時期として算出する基本点火時期算出処理手順と
を含むことを特徴とするエンジンの点火時期制御方法。 - 前記燃焼室内ガスの運動エネルギーを、少なくとも燃焼室内ガス質量とピストン速度に基づいて算出することを特徴とする請求項１に記載のエンジンの点火時期制御方法。
- 前記燃焼室内ガスの運動エネルギーを、少なくとも、吸気ポート圧力と、燃焼室内圧力と、吸気弁開口面積とに基づいて算出することを特徴とする請求項１に記載のエンジンの点火時期制御方法。
- 前記燃焼室内ガスの運動エネルギーを、少なくとも吸気弁のプロファイルに基づいて算出することを特徴とする請求項１に記載のエンジンの点火時期制御方法。
- 前記火炎の形状を球形として前記火炎半径を算出することを特徴とする請求項１に記載のエンジンの点火時期制御方法。
- 点火の瞬間から火炎が形成されるとみなし、この点火直後の火炎の速度を、層流燃焼速度と燃焼室内ガスの運動エネルギーと燃焼室内ガスの熱膨張率に基づいて算出する火炎速度算出手段と、
初期燃焼時の燃焼質量割合を基準燃焼質量割合として、この基準燃焼質量割合になったときの燃焼室内容積に基づいて基準燃焼質量割合になったときの火炎半径を算出する火炎半径算出手段と、
点火から燃焼質量割合が前記基準燃焼質量割合に上昇するまでの時間を着火遅れ時間とみなして、この着火遅れ時間を前記点火直後の火炎速度と前記基準燃焼質量割合になったときの火炎半径とに基づいて算出する着火遅れ時間算出手段と、
燃焼質量割合が前記基準燃焼質量割合から基準クランク角での燃焼質量割合になるまでの時間である燃焼時間を算出する燃焼時間算出手段と、
これら着火遅れ時間と燃焼時間とを加算した値をクランク角区間に変換するクランク角区間変換手段と、
前記基準クランク角よりこの変換されたクランク角区間だけ進角側の値をＭＢＴの得られる基本点火時期として算出する基本点火時期算出手段と
を含むことを特徴とするエンジンの点火時期制御装置。 - 前記燃焼室内ガスの運動エネルギーを、少なくとも燃焼室内ガス質量とピストン速度に基づいて算出することを特徴とする請求項６に記載のエンジンの点火時期制御装置。
- 前記燃焼室内ガスの運動エネルギーを、少なくとも、吸気ポート圧力と、燃焼室内圧力と、吸気弁開口面積とに基づいて算出することを特徴とする請求項６に記載のエンジンの点火時期制御装置。
- 前記燃焼室内ガスの運動エネルギーを、少なくとも吸気弁のプロファイルに基づいて算出することを特徴とする請求項６に記載のエンジンの点火時期制御装置。
- 前記火炎の形状を球形として前記火炎半径を算出することを特徴とする請求項６に記載のエンジンの点火時期制御装置。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2006059318A JP4872383B2 (ja) | 2006-03-06 | 2006-03-06 | エンジンの点火時期制御方法及びエンジンの点火時期制御装置 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2006059318A JP4872383B2 (ja) | 2006-03-06 | 2006-03-06 | エンジンの点火時期制御方法及びエンジンの点火時期制御装置 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2007239483A JP2007239483A (ja) | 2007-09-20 |
| JP4872383B2 true JP4872383B2 (ja) | 2012-02-08 |
Family
ID=38585274
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2006059318A Active JP4872383B2 (ja) | 2006-03-06 | 2006-03-06 | エンジンの点火時期制御方法及びエンジンの点火時期制御装置 |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP4872383B2 (ja) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP5237314B2 (ja) * | 2010-02-19 | 2013-07-17 | ヤンマー株式会社 | ディーゼルエンジン |
| JP5696568B2 (ja) * | 2011-04-04 | 2015-04-08 | トヨタ自動車株式会社 | 内燃機関の制御装置 |
Family Cites Families (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2003148316A (ja) * | 2001-11-09 | 2003-05-21 | Nissan Motor Co Ltd | 内燃機関の点火時期制御装置 |
| JP4366992B2 (ja) * | 2003-05-09 | 2009-11-18 | 日産自動車株式会社 | 内燃機関の点火時期制御装置 |
| JP4075862B2 (ja) * | 2004-06-04 | 2008-04-16 | 日産自動車株式会社 | エンジンの点火時期制御装置 |
-
2006
- 2006-03-06 JP JP2006059318A patent/JP4872383B2/ja active Active
Also Published As
| Publication number | Publication date |
|---|---|
| JP2007239483A (ja) | 2007-09-20 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP3975936B2 (ja) | ノッキング指標値算出装置 | |
| JP4036138B2 (ja) | 火花点火式内燃機関の燃焼制御装置 | |
| JP4774893B2 (ja) | エンジンの点火時期制御方法及びエンジンの点火時期制御装置 | |
| JP4604903B2 (ja) | エンジンの点火時期制御装置 | |
| JP4872383B2 (ja) | エンジンの点火時期制御方法及びエンジンの点火時期制御装置 | |
| JP4158747B2 (ja) | 内燃機関の点火時期制御装置 | |
| JP5696568B2 (ja) | 内燃機関の制御装置 | |
| JP4747763B2 (ja) | エンジンの点火時期制御方法及びエンジンの点火時期制御装置 | |
| JP4120620B2 (ja) | エンジンの制御装置 | |
| JP4241511B2 (ja) | エンジンのノック制御装置 | |
| JP6536613B2 (ja) | 内燃機関の制御装置 | |
| JP4075862B2 (ja) | エンジンの点火時期制御装置 | |
| JP4251069B2 (ja) | エンジンのノック検出装置及びノック制御装置 | |
| JP4055647B2 (ja) | エンジンの点火時期制御装置 | |
| JP4254605B2 (ja) | エンジンの点火時期制御装置 | |
| JP4244654B2 (ja) | エンジンの点火時期制御装置 | |
| JP4239607B2 (ja) | エンジンの点火制御装置 | |
| Lahti et al. | Dynamic engine control for HCCI combustion | |
| JP4135655B2 (ja) | 内燃機関の点火時期制御装置 | |
| JP2005226481A (ja) | エンジンのデポジット量検出装置及びノック制御装置 | |
| JP4241510B2 (ja) | エンジンのノック制御装置及びエンジンの点火時期制御装置。 | |
| JP2008196506A (ja) | エンジンのノック制御装置及びエンジンの点火時期制御装置 | |
| JP4151605B2 (ja) | エンジンのノック制御装置 | |
| Kaiadi et al. | Closed-loop combustion control for a 6-cylinder port-injected natural-gas engine | |
| JP4941280B2 (ja) | エンジンの残留ガス量推定装置及び残留ガス量推定方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20090225 |
|
| A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20110203 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20110405 |
|
| A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20110606 |
|
| TRDD | Decision of grant or rejection written | ||
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20111025 |
|
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20111107 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20141202 Year of fee payment: 3 |
|
| R150 | Certificate of patent or registration of utility model |
Ref document number: 4872383 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 Free format text: JAPANESE INTERMEDIATE CODE: R150 |