JP4075862B2 - エンジンの点火時期制御装置 - Google Patents

Description

本発明は、エンジン（内燃機関）の点火時期制御装置、特に最大の軸トルクを発生するのに必要な最小点火進角値（いわゆるＭＢＴ）となるように点火時期を制御するものに関する。

燃焼室内で混合気が燃焼する燃焼期間における燃焼速度を複数に分割し、その分割した各燃焼期間における乱流火炎速度を運転条件に基づいて算出し、これら複数の乱流火炎速度に基づいてＭＢＴの得られる点火時期を基本点火時期として算出するものがある（特許文献１参照）。
特開２００３−１４８２３６号公報

ところで、特許文献１に記載の技術によれば、気体の流れがない状態での火炎の伝播速度である層流火炎速度ＳＬ１、ＳＬ２に対して、燃焼期間における平均乱れ強さＳＴ１、ＳＴ２を乗算することによって気体の流れがある状態での火炎の伝播速度である乱流火炎速度ＦＬＡＭＥ１、ＦＬＡＭＥ２を算出している。

この場合に、燃焼期間における平均乱れ強さＳＴ１、ＳＴ２を算出するに際してはエンジンの回転速度ＮＲＰＭに応じた値とするだけの構成であり、特にスキッシュの影響そのものを扱ってはいない。

ここで、スキッシュとは、ピストン頭部とシリンダヘッド下面との間の狭い隙間に存在する混合気が、圧縮行程の終わり近くで燃焼室の主要空間に押し込まれることにより生じるシリンダ半径方向のガス流動のことで、このスキッシュは燃焼期間における平均乱れ強さに影響する。また、この影響はエンジンの機種毎に相違する。

このため、上記の特許文献１に記載の技術によれば、エンジンの機種が相違する毎に燃焼期間における平均乱れ強さを適合する必要があり、適合工数が増大することになっている。

このため本発明者らはエンジンの機種が相違してもなんとか適合工数を減らすことができないかとスキッシュと燃焼期間における平均乱れ強さの関係について各種実験を行い、新たな知見を得た。

そこで本発明は、この新たな知見に基づき、燃焼期間における平均乱れ強さを算出するに際して適合工数を最小限にした装置を提供することを目的とするものである。

本発明は、燃焼室内で混合気が燃焼する燃焼期間を複数に分割し、その分割した各燃焼期間（ＢＵＲＮ１、ＢＵＲＮ２）における燃焼速度（ＦＬＡＭＥ１、ＦＬＡＭＥ２）を運転条件に基づいて算出し、これら複数の燃焼速度に基づいてＭＢＴの得られる点火時期を基本点火時期として算出し、この基本点火時期で火花点火を行うようにしたエンジンの点火時期制御装置において、前記いずれか一つの燃焼期間がスキッシュの生じる期間を含み、このスキッシュの生じる期間を含む燃焼期間における燃焼速度がこのスキッシュの生じる期間を含む燃焼期間における乱流火炎速度であり、スキッシュによる火炎速度補正量（Ｋｑ）を少なくともエンジンの仕様に基づいて算出し、前記スキッシュの生じる期間を含む燃焼期間における乱流火炎速度を、このスキッシュの生じる期間を含む燃焼期間における層流火炎速度とエンジンの運転条件のうち少なくともいずれか一つと、前記スキッシュによる火炎速度補正量（Ｋｑ）とに基づいて算出するように構成する。

また、本発明は、燃焼室内で混合気が燃焼する燃焼期間を初期燃焼期間（ＢＵＲＮ１）と主燃焼期間（ＢＵＲＮ２）の複数に分割し、その分割した各燃焼期間（ＢＵＲＮ１、ＢＵＲＮ２）における燃焼速度（ＦＬＡＭＥ１、ＦＬＡＭＥ２）を運転条件に基づいて算出し、これら複数の燃焼速度（ＦＬＡＭＥ１、ＦＬＡＭＥ２）に基づいてＭＢＴの得られる点火時期を基本点火時期ＭＢＴＣＡＬとして算出し、この基本点火時期ＭＢＴＣＡＬで火花点火を行うようにしたエンジンの点火時期制御装置において、前記主燃焼期間がスキッシュの生じる期間を含み、前記主燃焼期間における燃焼速度が前記主燃焼期間における乱流火炎速度であり、スキッシュによる火炎速度補正量（Ｋｑ）を少なくともエンジンの仕様に基づいて算出し、前記主燃焼期間における乱流火炎速度（ＦＬＡＭＥ２）を、前記主燃焼期間における層流火炎速度とエンジンの運転条件のうち少なくともいずれか一つと、前記スキッシュによる火炎速度補正量（Ｋｑ）とに基づいて算出するように構成する。

スキッシュの形状によりスキッシュが燃焼速度に影響する程度が定まる、という実験結果を新たに得ている。ということは、スキッシュの形状はエンジンの仕様により定まるのであるから、エンジンの仕様によりスキッシュがスキッシュの生じる期間を含む燃焼期間における乱流火炎速度や主燃焼期間における乱流火炎速度に影響する程度をほぼ推定できることを意味する。すなわち、本発明によれば、いずれか一つの燃焼期間がスキッシュの生じる期間を含み、スキッシュによる火炎速度補正量（Ｋｑ）を少なくともエンジンの仕様に基づいて算出し、前記スキッシュの生じる期間を含む燃焼期間における乱流火炎速度（ＦＬＡＭＥ２）を、このスキッシュの生じる期間を含む燃焼期間における層流火炎速度とエンジンの運転条件のうち少なくともいずれか一つと、前記スキッシュによる火炎速度補正量（Ｋｑ）に基づいて算出するように、また本発明によれば主燃焼期間がスキッシュの生じる期間を含み、スキッシュによる火炎速度補正量（Ｋｑ）を少なくともエンジンの仕様に基づいて算出し、前記主燃焼期間における乱流火炎速度（ＦＬＡＭＥ２）を、前記主燃焼期間における層流火炎速度とエンジンの運転条件のうち少なくともいずれか一つと、前記スキッシュによる火炎速度補正量（Ｋｑ）に基づいて算出するように構成するので、エンジン機種が異なる毎に適合し直す必要が無くなり、適合工数を低減できる。

以下、図面に基づき本発明の実施形態について説明する。

図１は、本発明のシステムを説明するための概略図である。

空気は吸気コレクタ２に蓄えられた後、吸気マニホールド３を介して各気筒の燃焼室５に導入される。燃料は各気筒の吸気ポート４に配置された燃料インジェクタ２１より噴射供給される。空気中に噴射された燃料は気化しつつ空気と混合してガス（混合気）を作り、燃焼室５に流入する。この混合気は吸気弁１５が閉じることで燃焼室５内に閉じこめられ、ピストン６の上昇によって圧縮される。

この圧縮混合気に対して高圧火花により点火を行うため、パワートランジスタ内蔵の点火コイルを各気筒に配した電子配電システムの点火装置１１を備える。すなわち、点火装置１１は、バッテリからの電気エネルギーを蓄える点火コイル１３と、点火コイル１３の一次側への通電、遮断を行うパワートランジスタと、燃焼室５の天井に設けられ点火コイル１３の一次電流の遮断によって点火コイル１３の二次側に発生する高電圧を受けて、火花放電を行う点火プラグ１４とからなっている。

圧縮上死点より少し手前で点火プラグ１４により火花が飛ばされ圧縮混合気に着火されると、火炎が広がりやがて爆発的に燃焼し、この燃焼によるガス圧がピストン６を押し下げる仕事を行う。この仕事はクランクシャフト７の回転力として取り出される。燃焼後のガス（排気）は排気弁１６が開いたとき排気通路８へと排出される。

排気通路８には三元触媒９を備える。三元触媒９は排気の空燃比が理論空燃比を中心とした狭い範囲（ウインドウ）にあるとき、排気に含まれるＨＣ、ＣＯ、ＮＯｘといった有害三成分を同時に効率よく除去できる。空燃比は吸入空気量と燃料量の比であるので、エンジンの１サイクル（４サイクルエンジンではクランク角で７２０°区間）当たりに燃焼室５に導入される吸入空気量と、燃料インジェクタ２１からの燃料噴射量との比が理論空燃比となるように、エンジンコントローラ３１ではエアフローメータ３２からの吸入空気流量の信号とクランク角センサ（３３、３４）からの信号に基づいて燃料インジェクタ２１からの燃料噴射量を定めると共に、三元触媒９の上流に設けたＯ₂センサ３５からの信号に基づいて空燃比をフィードバック制御している。

吸気コレクタ２の上流には絞り弁２３がスロットルモータ２４により駆動される、いわゆる電子制御スロットル２２を備える。運転者が要求するトルクはアクセルペダル４１の踏み込み量（アクセル開度）に現れるので、エンジンコントローラ３１ではアクセルセンサ４２からの信号に基づいて目標トルクを定め、この目標トルクを実現するための目標空気量を定め、この目標空気量が得られるようにスロットルモータ２４を介して絞り弁２３の開度を制御する。

吸気弁用カムシャフト２５、排気弁用カムシャフト２６及びクランクシャフト７の各前部にはそれぞれカムスプロケット、クランクスプロケットが取り付けられ、これらスプロケットにタイミングチェーン（図示しない）を掛け回すことで、カムシャフト２５、２６がエンジンのクランクシャフト７により駆動されるのであるが、このカムスプロケットと吸気弁用カムシャフト２５との間に介在して、作動角一定のまま吸気弁用カムの位相を連続的に制御し得る吸気バルブタイミングコントロール機構（以下、「吸気ＶＴＣ機構」という。）２７と、カムスプロケットと排気弁用カムシャフト２６との間に介在して、作動角一定のまま排気弁用カムの位相を連続的に制御し得る排気バルブタイミングコントロール機構（以下、「排気ＶＴＣ機構」という。）２８とを備える。吸気弁１５の開閉時期や排気弁１６の開閉時期を変えると燃焼室５に残留する不活性ガスの量が変化する。燃焼室５内の不活性ガスの量が増えるほどポンピングロスが減って燃費がよくなるので、運転条件によりどのくらいの不活性ガスが燃焼室５内に残留したらよいかを目標吸気弁閉時期や目標排気弁閉時期にして予め定めており、エンジンコントローラ３１ではそのときの運転条件（エンジンの負荷と回転速度）より目標吸気弁閉時期と目標排気弁閉時期を定め、それら目標値が得られるように吸気ＶＴＣ機構２７、排気ＶＴＣ機構２８の各アクチュエータを介して吸気弁閉時期と排気弁閉時期を制御する。

さらに、吸気弁１５の弁リフト量及び作動角を連続的に可変制御する多節リンク状の機構で構成される可変バルブ機構（以下「ＶＥＬ機構」という。）５１を備える。このＶＥＬ機構５１及び上記吸気ＶＴＣ機構２７の具体的な構成は特開２００３−３８７２号公報により公知であるのでその詳しい説明は省略する。

弁リフト量を小さくすると吸気弁１５の周囲を通過して燃焼室５に流入する吸入空気の流速が増大するためこの吸入空気の流速増大で噴射燃料の霧化が改善されることから、エンジンコントローラ３１ではエンジンより排出されるＨＣが増える傾向にある始動直後からしばらくは噴射燃料の霧化を促進してＨＣを低減するため、ＶＥＬ機構５１のアクチュエータを介しＶＥＬ機構５１を作動させて吸気弁１５の弁リフト量を小さくする。そして、その後はＶＥＬ機構５１を非作動として吸気弁１５の弁リフト量を元に戻し大量の吸入空気を燃焼室５へと流入させる。

吸気温度センサ４３からの吸気温度の信号、吸気圧力センサ４４からの吸気圧力の信号、排気温度センサ４５からの排気温度の信号、排気圧力センサ４６からの排気圧力の信号が、水温センサ３７からの冷却水温の信号と共に入力されるエンジンコントローラ３１では、パワートランジスタ１３を介して点火プラグ１４の一次側電流の遮断時期である点火時期を制御する。

図２はエンジンコントローラ３１内で行われる点火時期制御のブロック図で、大きくは点火時期演算部５１と点火時期制御部６１とからなる。点火時期演算部５１はさらに初期燃焼期間算出部５２、主燃焼期間算出部５３、燃焼期間算出部５４、基本点火時期算出部５５及び前回燃焼開始時期算出部５６からなる。

初期燃焼期間算出部５２では、混合気が着火してから火炎核が形成されるまでの期間を初期燃焼期間ＢＵＲＮ１として算出する。主燃焼期間算出部５３では、火炎核が形成されてから燃焼圧力が最大値Ｐｍａｘに達するまでの期間を主燃焼期間ＢＵＲＮ２として算出する。燃焼期間算出部５４では、これら初期燃焼期間ＢＵＲＮ１と主燃焼期間ＢＵＲＮ２との合計を、燃焼開始より最大燃焼圧力Ｐｍａｘに至るまでの燃焼期間ＢＵＲＮとして算出する。基本点火時期算出部５５では、この燃焼期間ＢＵＲＮに基づいてＭＢＴの得られる点火時期（この点火時期を「基本点火時期」という。）ＭＢＴＣＡＬを算出する。

点火時期制御部６１ではこのようにして算出された基本点火時期を点火時期指令値とし、この指令値で点火プラグ１４が燃焼室５内の混合気に対して着火するように、イグニッションコイル１３への通電角と非通電角を制御する。

上記のように燃焼期間ＢＵＲＮを初期燃焼期間ＢＵＲＮ１と主燃焼期間ＢＵＲＮ２に分けて算出し、燃焼期間ＢＵＲＮに応じて基本点火時期ＭＢＴＣＡＬを求めるようにしたのは、燃焼解析より得られた結果に基づくものである。以下、燃焼解析に基づくこの点火時期制御をさらに説明する。

図３に示すようにＭＢＴ（最大トルクの得られる最小進角値）で混合気に点火した場合に混合気の燃焼圧力が最大値Ｐｍａｘとなるクランク角を基準クランク角θＰＭＡＸ［ｄｅｇＡＴＤＣ］とする。基準クランク角θＰＭＡＸは燃焼方式によらずほぼ一定であり、一般に圧縮上死点後１２〜１５度、最大で圧縮上死点後１０〜２０度の範囲にある。

図４に火花点火エンジンにおける燃焼室内の燃焼解析により得られた燃焼質量割合ＢＲ（燃焼ガス質量割合）の変化を示す。燃焼室に供給された燃料に対する燃焼質量の比率を表す燃焼質量割合ＢＲは、点火時に０％であり、完全燃焼によって１００％に達する。基準クランク角θＰＭＡＸにおける燃焼質量割合は一定で約６０％であることが実験により確かめられている。

燃焼質量割合ＢＲが０％から基準クランク角θＰＭＡＸ相当の約６０％に達するまでの変化代に相当する燃焼期間は、燃焼開始直後で燃焼質量割合にも燃焼圧力にもほとんど変化のない期間である初期燃焼期間と、燃焼質量割合と燃焼圧力が急激に増加する主燃焼期間とに分けられる。初期燃焼期間は、燃焼開始から火炎核が形成されるまでの段階であり、火炎核が形成されるのは燃焼質量割合が０％から２％〜１０％まで変化したときである。この初期燃焼期間中は、燃焼圧力や燃焼温度の上昇速度が小さく、燃焼質量割合の変化に対して初期燃焼期間は長い。初期燃焼期間の長さは燃焼室内の温度や圧力の変化の影響を受けやすい。

一方、主燃焼期間においては、火炎核から外側へと火炎が伝播するのであり、その火炎速度（つまり燃焼速度）が急上昇する。そのため、主燃焼期間の燃焼質量割合の変化は初期燃焼期間の燃焼質量割合の変化に比べて大きい。

エンジンコントローラ３１では、燃焼質量割合が２％に達する（変化する）までを初期燃焼期間ＢＵＲＮ１［ｄｅｇ］とし、初期燃焼期間ＢＵＲＮ１の終了後、基準クランク角θＰＭＡＸに至るまでの区間（燃焼室量割合でいえば２％より約６０％に達するまでの間）を主燃焼期間ＢＵＲＮ２［ｄｅｇ］として区別する。そして、初期燃焼期間ＢＵＲＮ１に主燃焼期間ＢＵＲＮ２を加えた合計である燃焼期間ＢＵＲＮ［ｄｅｇ］を算出し、この燃焼期間ＢＵＲＮから基準クランク角θＰＭＡＸ［ｄｅｇＡＴＤＣ］を差し引き、さらに後述する点火無駄時間相当クランク角ＩＧＮＤＥＡＤ［ｄｅｇ］を加えたクランク角位置を、ＭＢＴの得られる点火時期である基本点火時期ＭＢＴＣＡＬ［ｄｅｇＢＴＤＣ］として設定する。

火炎核の形成される初期燃焼期間での燃焼室５内の圧力、温度は、点火時の圧力、温度とほぼ等価になるが、これから点火時期を算出しようとしているのに、最初から正確な点火時期を設定することはできない。そこで、図２に示したように前回燃焼開始時期算出部５６で基本点火時期の前回値を前回燃焼開始時期ＭＢＴＣＹＣＬ［ｄｅｇＢＴＤＣ］として算出し、この値を初期燃焼期間算出部５２に対して与えるようにし、初期燃焼期間算出部５２において初期燃焼期間の算出をサイクリックに繰り返すことで、精度の高い結果を時間遅れなしに出すようにしている。

次に、エンジンコントローラ３１で実行される上記の基本点火時期ＭＢＴＣＡＬの算出を以下のフローチャートを参照しながら詳述する。

図５は点火時期の算出に必要な各種の物理量を算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。

まずステップ１１では、吸気弁閉時期ＩＶＣ［ｄｅｇＢＴＤＣ］、温度センサ４３により検出されるコレクタ内温度ＴＣＯＬ［Ｋ］、圧力センサ４４により検出されるコレクタ内圧力ＰＣＯＬ［Ｐａ］、温度センサ４５により検出される排気温度ＴＥＸＨ［Ｋ］、内部不活性ガス率ＭＲＥＳＦＲ［％］、温度センサ３７により検出される冷却水温ＴＷＫ［Ｋ］、目標当量比ＴＦＢＹＡ、クランク角センサにより検出されるエンジン回転速度ＮＲＰＭ［ｒｐｍ］、点火無駄時間ＤＥＡＤＴＩＭＥ［μｓｅｃ］を読み込む。

ここで、クランク角センサはクランクシャフト７のポジションを検出するポジションセンサ３３と、吸気用カムシャフト２５ポジションを検出するフェーズセンサ３４とからなり、これら２つのセンサ３３、３４からの信号に基づいてエンジン回転速度ＮＲＰＭ［ｒｐｍ］が算出されている。

吸気弁閉時期ＩＶＣは吸気ＶＴＣ機構２７に与える指令値から既知である。あるいはフェーズセンサ３４により実際の吸気弁閉時期を検出してもかまわない。

内部不活性ガス率ＭＲＥＳＦＲは燃焼室内に残留する不活性ガス量を燃焼室内の総ガス量で除した値で、その算出については後述する。点火無駄時間ＤＥＡＤＴＩＭＥは一定値である。

目標当量比ＴＦＢＹＡは図示しない燃料噴射量の算出フローにおいて算出されている。目標当量比ＴＦＢＹＡは無名数であり、理論空燃比を１４．７とすると、次式により表される値である。

ＴＦＢＹＡ＝１４．７／目標空燃比 …（１）
例えば（１）式より目標空燃比が理論空燃比のときＴＦＢＹＡ＝１．０となり、目標空燃比が例えば２２．０といったリーン側の値であるとき、ＴＦＢＹＡは１．０未満の正の値である。

ステップ１２では、燃焼室５の吸気弁閉時期ＩＶＣにおける容積（つまり圧縮開始時期での容積）ＶＩＶＣ［ｍ³］を算出する。燃焼室５の吸気弁閉時期における容積ＶＩＶＣは、ピストン６のストローク位置によって決まる。ピストン６のストローク位置はエンジンのクランク角位置によって決まる。

図６を参照して、エンジンのクランクシャフト７１の回転中心７２がシリンダの中心軸７３からオフセットしている場合を考える。コネクティングロッド７４、コネクティングロッド７４とクランクシャフト７１との結節点７５、コネクティングロッド７４とピストンをつなぐピストンピン７６が図に示す関係にあるとする。このときの、燃焼室５の吸気弁閉時期における容積ＶＩＶＣは次式（２）〜（６）で表すことができる。

ＶＩＶＣ＝ｆ１（θivc）
＝Ｖｃ＋（π／４）Ｄ²・Ｈivc …（２）
Ｖｃ＝（π／４）Ｄ²・Ｈx／（ε−１） …（３）
Ｈivc＝｛（ＣＮＤ＋ＳＴ²／２）−（ＣＲoff−ＰＩＳoff）²｝^1/2
−｛（ＳＴ／２）・cos（θivc＋θoff）｝＋（ＣＮＤ²−Ｘ²）^1/2
…（４）
Ｘ ＝（ＳＴ／２）・sin（θivc＋θoff）−ＣＲoff＋ＰＩＳoff
…（５）
θoff＝arcsin｛（ＣＲoff−ＰＩＳoff）／（ＣＮＤ・（ＳＴ／２））｝
…（６）
ただし、Ｖｃ：隙間容積［ｍ³］、
ε ：圧縮比、
Ｄ ：シリンダボア径［ｍ］、
ＳＴ ：ピストンの全ストローク［ｍ］、
Ｈivc ：吸気弁閉時期におけるピストンピン７６のＴＤＣからの
距離［ｍ］、
Ｈx ：ピストンピン７６のＴＤＣからの距離の最大値と最小値
の差［ｍ］、
ＣＮＤ ：コネクティングロッド７４の長さ［ｍ］、
ＣＲoff ：結節点７５のシリンダ中心軸７３からのオフセット距離 ［ｍ］、
ＰＩＳoff：クランクシャフト回転中心７２のシリンダ中心軸７３から のオフセット距離［ｍ］、
θivc ：吸気弁閉時期のクランク角［ｄｅｇＡＴＤＣ］、
θoff ：ピストンピン７６とクランクシャフト回転中心７２と
を結ぶ線がＴＤＣにおいて垂直線となす角度［ｄｅｇ］、
Ｘ ：結節点７５とピストンピン７６との水平距離［ｍ］、
吸気弁閉時期のクランク角θivcは前述のように、エンジンコントローラ３１から吸気ＶＴＣ機構２７への指令信号によって決まるので、既知である。式（２）〜（６）にこのときのクランク角θivc（＝ＩＶＣ）を代入すれば、燃焼室５の吸気弁閉時期における容積ＶＩＶＣを算出することができる。したがって、実用上は燃焼室５の吸気弁閉時期における容積ＶＩＶＣは吸気弁閉時期ＩＶＣをパラメータとするテーブルで設定したものを用いる。吸気ＶＴＣ機構２７を備えないときには定数で与えることができる。

ステップ１３では、燃焼室５の吸気弁閉時期ＩＶＣにおける温度（つまり圧縮開始時期温度）ＴＩＮＩ［Ｋ］を算出する。燃焼室５に流入するガスの温度は、燃焼室５に流入する新気と燃焼室５に残留する不活性ガスとが混じったガスの温度であり、燃焼室５に流入する新気の温度は吸気コレクタ２内の新気温度ＴＣＯＬに等しく、また燃焼室５内に残留する不活性ガスの温度は排気ポート部近傍の排気温度ＴＥＸＨで近似できるので、燃焼室５の吸気弁閉時期ＩＶＣにおける温度ＴＩＮＩは吸気弁閉時期ＩＶＣになったタイミングでの、吸気コレクタ２内の新気温度ＴＣＯＬ、排気温度ＴＥＸＨ、燃焼室５内に残留する不活性ガスの割合である内部不活性ガス率ＭＲＥＳＦＲから次式により求めることができる。

ＴＩＮＩ＝ＴＥＸＨ×ＭＲＥＳＦＲ＋ＴＣＯＬ×（１−ＭＲＥＳＦＲ）
…（７）
ステップ１４では燃焼室５の吸気弁閉時期ＩＶＣにおける圧力（つまり圧縮開始時期圧力）ＰＩＮＩ［Ｐａ］を算出する。すなわち、吸気弁閉時期ＩＶＣになったタイミングでのコレクタ内圧力ＰＣＯＬを吸気弁閉時期ＩＶＣにおける圧力ＰＩＮＩとして取り込む。

ステップ１５では、燃焼室５内の混合気の燃えやすさを表す反応確率ＲＰＲＯＢＡ［％］を算出する。反応確率ＲＰＲＯＢＡは無次元の値であり、残留不活性ガス率ＭＲＥＳＦＲ、冷却水温ＴＷＫ［Ｋ］、目標当量比ＴＦＢＹＡの３つのパラメータに依存するので、次式により表すことができる。

ＲＰＲＯＢＡ＝ｆ３（ＭＲＥＳＦＲ、ＴＷＫ、ＴＦＢＹＡ） …（８）
具体的に説明すると、ＭＲＥＳＦＲ、ＴＷＫ、ＴＦＢＹＡの３つのパラメータの組み合わせによって得られる反応確率の最大値を１００％とし、これらのパラメータと反応確率ＲＰＲＯＢＡの関係を実験的に求め、求めた反応確率ＲＰＲＯＢＡをパラメータに応じたテーブルとしてエンジンコントローラ３１のメモリに予め格納しておく。ステップ１５ではパラメータに応じてこのテーブルを検索することにより反応確率ＲＰＲＯＢＡを求める。

具体的には、冷却水温ＴＷＫに応じて図７に示すような特性を有する水温補正係数のテーブルと、同様に設定された内部不活性ガス率補正係数のテーブル（図示しない）と、目標当量比ＴＦＢＹＡに応じて図８に示すような特性を有する当量比補正係数のテーブルを予めメモリに格納しておく。各補正係数の最大値はそれぞれ１．０であり、３種類の補正係数の積に反応確率の最大値１００％を掛け合わせることで、反応確率ＲＰＲＯＢＡを算出する。

各テーブルを説明すると、図７に示す水温補正係数は冷却水温ＴＷＫが高いほど大きく、冷却水温ＴＷＫが８０℃以上では１．０になる。図８に示す当量比補正係数は目標当量比ＴＦＢＹＡが１．０のとき、つまり理論空燃比のときに最大値の１．０となり、目標当量比が１．０より大きくても小さくても当量比補正係数は減少する。内部不活性ガス率補正係数は図示しないが、内部不活性ガス率ＭＲＥＳＦＲがゼロの場合に１．０となる。

ステップ１６では、基準クランク角θＰＭＡＸ［ｄｅｇＡＴＤＣ］を算出する。前述のように基準クランク角θＰＭＡＸはあまり変動しないが、それでもエンジン回転速度ＮＲＰＭの上昇に応じて進角する傾向があるため、基準クランク角θＰＭＡＸはエンジン回転速度ＮＲＰＭの関数として次式で表すことができる。

θＰＭＡＸ＝ｆ４（ＮＲＰＭ） …（９）
具体的にはエンジン回転速度ＮＲＰＭから、エンジンコントローラ３１のメモリに予め格納された図９に示す特性のテーブルを検索することにより基準クランク角θＰＭＡＸを求める。算出を容易にするために、基準クランク角θＰＭＡＸを一定とみなすことも可能である。

最後にステップ１７では、点火無駄時間相当クランク角ＩＧＮＤＥＡＤ［ｄｅｇ］を算出する。点火無駄時間相当クランク角ＩＧＮＤＥＡＤは、エンジンコントローラ３１から点火コイル１３の一次電流を遮断する信号を出力したタイミングから点火プラグ１４が実際に点火するまでのクランク角区間で、次式により表すことができる。

ＩＧＮＤＥＡＤ＝ｆ５（ＤＥＡＤＴＩＭＥ、ＮＲＰＭ） …（１０）
ここでは、点火無駄時間ＤＥＡＤＴＩＭＥを２００μｓｅｃとする。（１０）式は、エンジン回転速度ＮＲＰＭから点火無駄時間ＤＥＡＤＴＩＭＥに相当するクランク角である点火無駄時間相当クランク角ＩＧＮＤＥＡＤを算出するためのものである。

図１０は初期燃焼期間ＢＵＲＮ１［ｄｅｇ］を算出するためのもの、また図１２は主燃焼期間ＢＵＲＮ２［ｄｅｇ］を算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。図１０、図１２は図５に続けて実行する。図１０、図１２はどちらを先に実行してもかまわない。

まず図１０から説明すると、ステップ１６１では、前回燃焼開始時期ＭＢＴＣＹＣＬ［ｄｅｇＢＴＤＣ］、図５のステップ１２で算出されている燃焼室５の吸気弁閉時期における容積ＶＩＶＣ［ｍ³］、図５のステップ１３で算出されている燃焼室５の吸気弁閉時期における温度ＴＩＮＩ［Ｋ］、図５のステップ１４で算出されている燃焼室５の吸気弁閉時期における圧力ＰＩＮＩ［Ｐａ］、エンジン回転速度ＮＲＰＭ［ｒｐｍ］、図５のステップ１５で算出されている反応確率ＲＰＲＯＢＡ［％］を読み込む。

ここで、前回燃焼開始時期ＭＢＴＣＹＣＬは、基本点火時期ＭＢＴＣＡＬの［ｄｅｇＢＴＤＣ］の１サイクル前の値であり、その算出については図１３により後述する。

ステップ１６２では燃焼室５の燃焼開始時期における容積Ｖ０［ｍ³］を算出する。前述したように、ここでの点火時期（燃焼開始時期）は今回のサイクルで演算する基本点火時期ＭＢＴＣＡＬではなく基本点火時期の１サイクル前の値である。すなわち、基本点火時期の１サイクル前の値であるＭＢＴＣＹＣＬから次式により燃焼室５の燃焼開始時期における容積Ｖ０を算出する。

Ｖ０＝ｆ６（ＭＢＴＣＹＣＬ） …（１１）
具体的には前回燃焼開始時期ＭＢＴＣＹＣＬにおけるピストン６のストローク位置と、燃焼室５のボア径から、燃焼室５のＭＢＴＣＹＣＬにおける容積Ｖ０を算出する。図５のステップ１２では、燃焼室５の吸気弁閉時期ＩＶＣにおける容積ＶＩＶＣを、吸気弁閉時期をパラメータとする吸気弁閉時期容積のテーブルを検索することにより求めたが、ここではＭＢＴＣＹＣＬをパラメータとする前回燃焼開始時期容積のテーブルを検索することにより、燃焼室５の前回燃焼開始時期ＭＢＴＣＹＣＬにおける容積Ｖ０を求めればよい。

ステップ１６３では燃焼開始時期における有効圧縮比Ｅｃを算出する。有効圧縮比Ｅｃは無次元の値であり、次式に示すように燃焼室５の燃焼開始時期における容積Ｖ０を燃焼室５の吸気弁閉時期における容積ＶＩＶＣで除した値である。

Ｅｃ＝ｆ７（Ｖ０、ＶＩＶＣ）＝Ｖ０／ＶＩＶＣ …（１２）
ステップ１６４では吸気弁閉時期ＩＶＣから燃焼開始時期に至る間の燃焼室５内の温度上昇率ＴＣＯＭＰを次式に示すように有効圧縮比Ｅｃに基づいて算出する。

ＴＣＯＭＰ＝ｆ８（Ｅｃ）＝Ｅｃ＾（κ−１） …（１３）
ただし、κ：比熱比、
（１３）式は断熱圧縮されるガスの温度上昇率の式である。なお、（１３）式右辺の「＾」は累乗計算を表している。この記号は後述する式でも使用する。

κは断熱圧縮されるガスの定圧比熱を定容比熱で除した値で、断熱圧縮されるガスが空気であればκ＝１．４であり、簡単にはこの値を用いればよい。ただし、混合気に対してκの値を実験的に求めることで、一層の算出精度の向上が可能である。

図１１は（１３）式を図示したものである。従って、このような特性のテーブルを予めエンジンコントローラ３１のメモリに格納しておき、有効圧縮比Ｅｃに基づき当該テーブルを検索することにより温度上昇率ＴＣＯＭＰを求めることも可能である。

ステップ１６５では、燃焼室５の燃焼開始時期における温度Ｔ０［Ｋ］を、燃焼室５の吸気弁閉時期における温度ＴＩＮＩに温度上昇率ＴＣＯＭＰを乗じることで、つまり
Ｔ０＝ＴＩＮＩ×ＴＣＯＭＰ …（１４）
の式により算出する。

ステップ１６６、１６７はステップ１６４、１６５と同様である。すなわち、ステップ１６６では吸気弁閉時期ＩＶＣから燃焼開始時期に至る間の燃焼室５内の圧力上昇率ＰＣＯＭＰを次式に示すように有効圧縮比Ｅｃに基づいて算出する。

ＰＣＯＭＰ＝ｆ９（Ｅｃ）＝Ｅｃ＾κ …（４１）
ただし、κ：比熱比、
（４１）式も（１３）式と同じに断熱圧縮されるガスの圧力上昇率の式である。（４１）式右辺の「＾」も（１３）式と同じに累乗計算を表している。

κは上記（１３）式で用いている値と同じで、断熱圧縮されるガスが空気であればκ＝１．４であり、簡単にはこの値を用いればよい。ただし、混合気に対してその組成、温度からκの値を求めることで、一層の算出精度の向上が可能である。

図１１と同様の特性のテーブルを予めエンジンコントローラ３１のメモリに格納しておき、有効圧縮比Ｅｃに基づき当該テーブルを検索することにより圧力上昇率ＰＣＯＭＰを求めることも可能である。

ステップ１６７では、燃焼室５の燃焼開始時期における圧力Ｐ０［Ｐａ］を、燃焼室５の吸気弁閉時期における圧力ＰＩＮＩに圧力上昇率ＰＣＯＭＰを乗じることで、つまり
Ｐ０＝ＰＩＮＩ×ＰＣＯＭＰ …（４２）
の式により算出する。

ステップ１６８では、初期燃焼期間における層流火炎速度ＳＬ１［ｍ／ｓｅｃ］を次式（公知）により算出する。

ＳＬ１＝ｆ１０（Ｔ０、Ｐ０）
＝ＳＬstd×｛（Ｔ０×Ｔstd）＾２．１８｝
×｛（Ｐ０／Ｐstd）＾（−０．１６）｝ …（１５）
ただし、Ｔstd ：基準温度［Ｋ］、
Ｐstd ：基準圧力［Ｐａ］、
ＳＬstd：基準温度Ｔstdと基準圧力Ｐstdにおける基準層流火炎速度
［ｍ／ｓｅｃ］、
Ｔ０ ：燃焼室５の燃焼開始時期における温度［Ｋ］、
Ｐ０ ：燃焼室５の燃焼開始時期における圧力［Ｐａ］、
層流火炎速度は気体の流れがない状態での火炎の伝播速度のことであり、燃焼室５内の圧縮速度、燃焼室５内の吸気流速に因らず、燃焼室５の温度及び圧力の関数となることが知られていることから、初期燃焼期間における層流火炎速度を燃焼開始時温度Ｔ０と燃焼開始時圧力Ｐ０の関数として、また後述するように主燃焼期における層流火炎速度を圧縮上死点時温度ＴＴＤＣと圧縮上死点圧力ＰＴＤＣの関数としている。これは、層流火炎速度は一般的に、エンジン負荷、燃焼室５内の不活性ガス率、吸気弁閉時期、比熱比、吸気温度により変化するのであるが、これらは燃焼室５内の温度Ｔと圧力Ｐに影響する因子であるので、層流火炎速度は最終的に燃焼室５内の温度Ｔと圧力Ｐにより規定できるとするものである。

上記の（１５）式において基準温度Ｔstdと基準圧力Ｐstdと基準層流火炎速度ＳＬstdは実験により予め定められる値である。

燃焼室５の通常の圧力である２ｂａｒ以上の圧力下では、（１５）式の圧力項（Ｐ０／Ｐstd）＾（−０．１６）は小さな値となる。従って、圧力項（Ｐ０／Ｐstd）＾（−０．１６）を一定値として、基準層流火炎速度ＳＬstdを基準温度Ｔstdのみで規定することも可能である。

従って、基準温度Ｔstdが５５０［Ｋ］で、基準層流火炎速度ＳＬstdが１．０［ｍ／ｓｅｃ］で、圧力項が０．７である場合の燃焼開始時期における温度Ｔ０と層流火炎速度ＳＬ１との関係は近似的に次式で定義することができる。

ＳＬ１＝ｆ１１（Ｔ０）
＝１．０×０．７×（Ｔ０／５５０）＾２．１８ …（１６）
ステップ１６９では、初期燃焼期間におけるガス流動の平均乱れ強さＳＴ１を算出する。このガス流動の平均乱れ強さＳＴ１は無次元の値であり、詳細には後述する。

ステップ１７０では初期燃焼期間における層流火炎速度Ｓ１とこの初期燃焼期間におけるガス流動の平均乱れ強さＳＴ１から、初期燃焼期間におけるガスの乱流火炎速度ＦＬＡＭＥ１［ｍ／ｓｅｃ］（燃焼速度）を次式により算出する。

ＦＬＡＭＥ１＝ＳＬ１×ＳＴ１ …（１８）
燃焼室５内にガス乱れがあるとガスの火炎速度が変化する。（１８）式はこのガス乱れに伴う火炎速度への寄与（影響）を考慮して、気体の流れがある状態での火炎速度である乱流火炎速度を算出するようにしたものである。

ステップ１７１では、このようにして算出した乱流火炎速度ＦＬＡＭＥ１に基づいて次式により初期燃焼期間ＢＵＲＮ１［ｄｅｇ］を算出する。

ＢＵＲＮ１＝｛（ＮＲＰＭ×６）×（ＢＲ１×Ｖ０）｝
／（ＲＰＲＯＢＡ×ＡＦ１×ＦＬＡＭＥ１） …（１９）
ただし、ＡＦ１：火炎核の反応面積（固定値）［ｍ²］、
この（１９）式および後述する（２２）式は、燃焼ガス質量を乱流火炎速度（燃焼速度）で割ると燃焼期間が得られるとする次の基本式より導いたものであるが、（１９）式、後述する（２２）式右辺の分子、分母ががただちに燃焼ガス質量、燃焼速度を表すものではない。

燃焼期間［ｓｅｃ］＝シリンダ内総質量［ｇ］
／（未燃ガス密度［ｇ／ｍ³］
×火炎表面積［ｍ²］×乱流火炎速度［ｍ／ｓｅｃ］）
…（補１）
（補１）式右辺分母の未燃ガス密度は、未燃ガス質量［ｇ］を未燃ガス体積［ｍ³］で割った値であるので、従来装置（特開平１０−３０５３５号公報）のように質量に相当する充填効率ＩＴＡＣのみの関数では未燃ガス密度を正確に計算できているとはいえない。そこで、（補１）式に対して実験結果とを照らし合わせつつ所定の近似を導入して初めて得られたのが上記（１９）式及び後述する（２２）式に示す実験式である。

ここで、（１９）式右辺のＢＲ１は燃焼開始時期より初期燃焼期間ＢＵＲＮ１の終了時期までの燃焼質量割合の変化代であり、ここではＢＲ１＝２％に設定している。（１９）式右辺の（ＮＲＰＭ×６）は単位をｒｐｍからクランク角（ｄｅｇ）に変換するための処理である。火炎核の反応面積ＡＦ１は実験的に設定される。

また、初期燃焼期間中はほぼ燃焼室容積は変わらないとみなすことができる。従って、初期燃焼期間ＢＵＲＮ１を算出するに際して最初の燃焼室容積である燃焼開始時の燃焼室容積Ｖ０を採用している。

次に図１２のフローに移ると、ステップ１８１では図１０のステップ１６１と同様に、図５のステップ１２で算出されている燃焼室５の吸気弁閉時期における容積ＶＩＶＣ［ｍ³］、図５のステップ１３で算出されている燃焼室５の吸気弁閉時期における温度ＴＩＮＩ［Ｋ］、図５のステップ１４で算出されている燃焼室５の吸気弁閉時期における圧力ＰＩＮＩ［Ｐａ］、エンジン回転速度ＮＲＰＭ［ｒｐｍ］、図５のステップ１５で算出されている反応確率ＲＰＲＯＢＡ［％］を読み込み、さらにシリンダ新気量ＭＡＣＹＬ［ｇ]、目標当量比ＴＦＢＹＡ、内部不活性ガス量ＭＲＥＳ［ｇ］、外部不活性ガス量ＭＥＧＲ［ｇ］を読み込む。

ここで、図１には外部ＥＧＲ装置は示していないが、図１２に関する限り外部ＥＧＲ装置を備えているエンジンを前提として説明する。この場合に、外部不活性ガス量ＭＥＧＲは例えば公知の手法（特開平１０−１４１１５０号公報参照）を用いて算出すればよい。なお、図１に示す本実施形態のように外部ＥＧＲ装置を備えていないエンジンを対象とするときには外部不活性ガス量ＭＥＧＲ＝０で扱えば足りる。シリンダ新気量ＭＡＣＹＬ、内部不活性ガス量ＭＲＥＳの算出については図１４以降で後述する。

ステップ１８２、１８３は図１０のステップ１６３、１６４と同様である。すなわち、ステップ１８２で圧縮上死点時期における有効圧縮比Ｅｃ ２を算出する。有効圧縮比Ｅｃ ２も上記（１２）式の有効圧縮比Ｅｃと同様に無次元の値であり、次式に示すように燃焼室５の圧縮上死点時における容積ＶＴＤＣを燃焼室５の吸気弁閉時期における容積ＶＩＶＣで除した値である。

Ｅｃ ２＝ｆ１３（ＶＴＤＣ、ＶＩＶＣ）＝ＶＴＤＣ／ＶＩＶＣ
…（４３）
（４３）式において燃焼室５の圧縮上死点時における容積ＶＴＤＣは運転条件によらず一定であり、予めエンジンコントローラ３１のメモリに格納しておけばよい。

ステップ１８３では吸気弁閉時期ＩＶＣから圧縮上死点に至る間の燃焼室５内の断熱圧縮による温度上昇率ＴＣＯＭＰ ２を次式に示すように有効圧縮比Ｅｃ ２に基づいて算出する。

ＴＣＯＭＰ ２＝ｆ１４（Ｅｃ ２）
＝Ｅｃ ２＾（κ−１）…（４４）
ただし、κ：比熱比、
図１１と同様の特性のテーブルを予めエンジンコントローラ３１のメモリに格納しておき、有効圧縮比Ｅｃ ２から当該テーブルを検索することにより温度上昇率ＴＣＯＭＰ ２を求めることも可能である。

ステップ１８４ではシリンダ新気量ＭＡＣＹＬ、目標当量比ＴＦＢＹＡ、内部不活性ガス量ＭＲＥＳ、外部不活性ガス量ＭＥＧＲから次式により燃焼室５の総ガス質量ＭＧＡＳ［ｇ］を算出する。

ＭＧＡＳ＝ＭＡＣＹＬ×(１＋ＴＦＢＹＡ／１４．７)＋ＭＲＥＳ＋ＭＥＧＲ
…（４５）
（４５）式右辺の括弧内の「１」は新気分、「ＴＦＢＹＡ／１４．７」は燃料分である。

ステップ１８５ではこの燃焼室５の総ガス質量ＭＧＡＳと、シリンダ新気量ＭＡＣＹＬ、目標当量比ＴＦＢＹＡを用い、次式により混合気の燃焼による温度上昇量（燃焼上昇温度）ＴＢＵＲＮ［Ｋ］を算出する。

ＴＢＵＲＮ＝｛ＭＡＣＹＬ×ＴＦＢＹＡ／１４．７×ＢＲｋ×Ｑ｝
／（Ｃｖ×ＭＧＡＳ）
…（４６）
ただし、Ｑ ：燃料の定発熱量、
ＢＲｋ：シリンダ内燃料の燃焼質量割合、
Ｃｖ ：定積比熱、
（４６）式右辺の分子はシリンダ内燃料による発生総熱量［Ｊ］、分母は単位発生熱量当たりの温度上昇率［Ｊ／Ｋ］を意味している。すなわち、（４６）式は熱力学の公式に当てはめた近似式である。

ここで、シリンダ内燃料の燃焼質量割合ＢＲｋとしては予め実験等で適合しておく。簡易的には例えば６０％／２＝３０％を設定する。これは、本実施形態では燃焼質量割合ＢＲが約６０％に達するまでを燃焼期間として扱うので、そのちょうど中間の３０％をＢＲｋとして設定するものである。

燃料の定発熱量Ｑは燃料の種類により異なる値であるので、燃料の種類に応じ予め実験等で求めておく。定積比熱Ｃｖは２〜３の値であり予め実験等で代表値を適合しておく。ただし、混合気に対してその組成、温度から定積比熱Ｃｖの値を求めることで、一層の算出精度の向上が可能である。

ステップ１８６では、燃焼室５の圧縮上死点における温度ＴＴＤＣ［Ｋ］を、燃焼室５の吸気弁閉時期における温度ＴＩＮＩに圧縮上死点までの温度上昇率ＴＣＯＭＰ ２を乗じその乗算値に上記の燃焼上昇温度ＴＢＵＲＮを加算することで、つまり次式により算出する。

ＴＴＤＣ＝ＴＩＮＩ×ＴＣＯＭＰ ２＋ＴＢＵＲＮ
…（４７）
ステップ１８７では、この燃焼室５の圧縮上死点における温度ＴＴＤＣと容積ＶＴＤＣ及び燃焼室５の吸気弁閉時期における圧力ＰＩＮＩ、容積ＶＩＶＣ及び温度ＴＩＮＩから次式により燃焼室５の圧縮上死点における圧力ＰＴＤＣ［Ｋ］を算出する。

ＰＴＤＣ＝ＰＩＮＩ×ＶＩＶＣ×ＴＴＤＣ／（ＶＴＤＣ×ＴＩＮＩ）
…（４８）
（４８）式は状態方程式を用いて得たものである。すなわち、吸気弁閉時期における圧力、容積及び温度（ＰＩＮＩ、ＶＩＶＣ、ＴＩＮＩ）を用いて次の状態方程式が成立する。

ＰＩＮＩ×ＶＩＶＣ＝ｎ・Ｒ・ＴＩＮＩ…（補２）
ただし、ｎ：モル数、
Ｒ：ガス定数、
圧縮上死点近傍では容積はほぼ等しいので、圧縮上死点での圧力、容積及び温度（ＰＴＤＣ、ＶＴＤＣ、ＴＴＤＣ）を用いて次の状態方程式が成立する。

ＰＴＤＣ×ＶＴＤＣ＝ｎ・Ｒ・ＴＴＤＣ…（補３）
この（補３）式と上記（補２）との両式からｎ・Ｒを消去しＰＴＤＣについて解くと、上記（４８）式が得られる。

ステップ１８８では図１０のステップ１６８と同様にして、次式（公知）により、主燃焼期間における層流火炎速度ＳＬ２［ｍ／ｓｅｃ］を算出する。

ＳＬ２＝ｆ１５（ＴＴＤＣ、ＰＴＤＣ）
＝ＳＬstd×｛（ＴＴＤＣ×Ｔstd）＾２．１８｝
×｛（ＰＴＤＣ／Ｐstd）＾（−０．１６）｝
…（４９）
ただし、Ｔstd ：基準温度［Ｋ］、
Ｐstd ：基準圧力［Ｐａ］、
ＳＬstd：基準温度Ｔstdと基準圧力Ｐstdにおける基準層流火炎速度
［ｍ／ｓｅｃ］、
ＴＴＤＣ：燃焼室５の圧縮上死点における温度［Ｋ］、
ＰＴＤＣ：燃焼室５の圧縮上死点における圧力［Ｐａ］、
（４９）式の解説は上記（１６）式と同様ある。すなわち、（４９）式の基準温度Ｔstdと基準圧力Ｐstdと基準層流火炎速度ＳＬstdは実験により予め定められる値である。燃焼室５の通常の圧力である２ｂａｒ以上の圧力下では、（４９）式の圧力項（ＰＴＤＣ／Ｐstd）＾（−０．１６）は小さな値となる。従って、圧力項（ＰＴＤＣ／Ｐstd）＾（−０．１６）を一定値として、基準層流火炎速度ＳＬstdを基準温度Ｔstdのみで規定することも可能である。よって、基準温度Ｔstdが５５０［Ｋ］で、基準層流火炎速度ＳＬstdが１．０［ｍ／ｓｅｃ］で、圧力項が０．７である場合の圧縮上死点における温度ＴＴＤＣと層流火炎速度ＳＬ２との関係は近似的に次式で定義することができる。

ＳＬ２＝ｆ１６（ＴＴＤＣ）
＝１．０×０．７×（ＴＴＤＣ／５５０）＾２．１８
…（５０）
ステップ１８９では主燃焼期間におけるガス流動の平均乱れ強さＳＴ２を算出する。この主燃焼期間におけるガス流動の平均乱れ強さＳＴ２も初期燃焼期間におけるガス流動の平均乱れ強さＳＴ１と同様に、詳細は後述する。

ステップ１９０では、主燃焼期間における層流火炎速度ＳＬ２［ｍ／ｓｅｃ］と主燃焼期間におけるガス流動の平均乱れ強さＳＴ２とから、主燃焼期間における乱流火炎速度ＦＬＡＭＥ２［ｍ／ｓｅｃ］（燃焼速度）を次式により算出する。

ＦＬＡＭＥ２＝ＳＬ２×ＳＴ２ …（２１）
ただし、ＳＬ２：層流火炎速度［ｍ／ｓｅｃ］、
（２１）式は（１８）式と同様、ガス乱れに伴う乱流火炎速度への寄与を考慮して、気体の流れがある状態での火炎速度である乱流火炎速度を算出するようにしたものである。

ステップ１９１では、このようにして算出した主燃焼期間における乱流火炎速度ＦＬＡＭＥ２に基づいて主燃焼期間ＢＵＲＮ２［ｄｅｇ］を上記の（１９）式に類似した次式で算出する。

ＢＵＲＮ２＝｛（ＮＲＰＭ×６）×（ＢＲ２×ＶＴＤＣ）｝
／（ＲＰＲＯＢＡ×ＡＦ２×ＦＬＡＭＥ２） …（２２）
ただし、ＡＦ２：火炎核の反応面積［ｍ²］
ここで、（２２）式右辺のＢＲ２は主燃焼期間の開始時期より終了時期までの燃焼質量割合の変化代である。初期燃焼期間の終了時期に燃焼質量割合ＢＲが２％になり、その後、主燃焼期間が開始し、燃焼質量割合ＢＲが６０％に達して主燃焼期間が終了すると考えているので、ＢＲ２＝６０％−２％＝５８％を設定している。ＡＦ２は火炎核の成長行程における平均の反応面積であり、（１９）式のＡＦ１と同様に、予め実験的に定めた固定値とする。

主燃焼期間では圧縮上死点を挟んで燃焼室容積が変化する。つまり、主燃焼期間の開始時期と、主燃焼期間の終了時期のほぼ中央に圧縮上死点位置が存在するとみなすことができる。また、圧縮上死点付近ではクランク角が変化しても燃焼室容積があまり変化しない。そこで主燃焼期間での燃焼室容積としてはこの圧縮上死点での燃焼室容積ＶＴＤＣで代表させることとしている。

図１３は基本点火時期ＭＢＴＣＡＬ［ｄｅｇＢＴＤＣ］を算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。図１０、図１２のうち遅く実行されるフローに続けて実行する。

ステップ４１では、図１０のステップ１７１で算出されている初期燃焼期間ＢＵＲＮ１、図１２のステップ１９１で算出されている主燃焼期間ＢＵＲＮ２、図５のステップ１７で算出されている点火時期無駄時間相当クランク角ＩＧＮＤＥＡＤ、図５のステップ１６で算出されている基準クランク角θＰＭＡＸを読み込む。

ステップ４２では、初期燃焼期間ＢＵＲＮ１と主燃焼期間ＢＵＲＮ２の合計を燃焼期間ＢＵＲＮ［ｄｅｇ］として算出する。

ステップ４３では次式により基本点火時期ＭＢＴＣＡＬ［ｄｅｇＢＴＤＣ］を算出する。

ＭＢＴＣＡＬ＝ＢＵＲＮ−θＰＭＡＸ＋ＩＧＮＤＥＡＤ …（２３）
ステップ４４では、この基本点火時期ＭＢＴＣＡＬから点火無駄時間相当クランク角ＩＧＮＤＥＡＤを差し引いた値を前回燃焼開始時期ＭＢＴＣＹＣＬ［ｄｅｇＢＴＤＣ］として算出する。

このようにして算出した基本点火時期ＭＢＴＣＡＬは、点火時期指令値として点火レジスタに移され、実際のクランク角がこの点火時期指令値と一致したタイミングでエンジンコントローラ３１より一次電流を遮断する点火信号が点火コイル１３に出力される。

また、今サイクルの点火時期指令値としてステップ４３で算出された基本点火時期ＭＢＴＣＡＬが用いられたとすると、次サイクルの点火時期になるまでの間、ステップ４４で算出された前回燃焼開始時期ＭＢＴＣＹＣＬが図１０のステップ１６２において用いられる。

次に、図１４は燃焼室５内の内部不活性ガス率ＭＲＥＳＦＲを算出するためのもので、一定時間毎（例えば１０ｍｓｅｃ毎）に実行する。このフローは上記図５のフローに先立って実行する。

ステップ５１ではエアフローメータ３２の出力と目標当量比ＴＦＢＹＡを読み込む。ステップ５２ではエアフロメータ３２の出力に基づいて、燃焼室５に流入する新気量（シリンダ新気量）ＭＡＣＹＬを算出する。このシリンダ新気量ＭＡＣＹＬの算出方法については公知の方法を用いればよい（特開２００１−５００９１公報参照）。

ステップ５３では、燃焼室５内の内部不活性ガス量ＭＲＥＳを算出する。この内部不活性ガス量ＭＲＥＳの算出については、図１５のフローにより説明する。

図１５（図１４ステップ５３のサブルーチン）においてステップ６１では、燃焼室５内の排気弁閉時期ＥＶＣにおける不活性ガス量ＭＲＥＳＣＹＬを算出する。この不活性ガス量ＭＲＥＳＣＹＬの算出についてはさらに図１６のフローにより説明する。

図１６（図１５ステップ６１のサブルーチン）においてステップ７１では、排気弁閉時期ＥＶＣ［ｄｅｇＢＴＤＣ］、温度センサ４５により検出される排気温度ＴＥＸＨ［Ｋ］、圧力センサ４６により検出される排気圧力ＰＥＸＨ［ｋＰａ］を読み込む。

ここで、吸気弁閉時期ＩＶＣが吸気ＶＴＣ機構２７に与える指令値から既知であったように、排気弁閉時期ＥＶＣも排気ＶＴＣ機構２８に与える指令値から既知である。

ステップ７２では燃焼室５の排気弁閉時期ＥＶＣにおける容積ＶＥＶＣを算出する。これは吸気弁閉時期ＩＶＣにおける容積ＶＩＶＣと同様に、排気弁閉時期をパラメータとするテーブルを検索することにより求めればよい。すなわち、排気弁ＶＴＣ機構２８を備える場合には、排気弁閉時期ＥＶＣから図２３に示すテーブルを検索することにより、燃焼室５の排気弁閉時期ＥＶＣにおける容積ＶＥＶＣを求めればよい。排気ＶＴＣ機構２８を備えないときには定数で与えることができる。

また、図示しないが圧縮比を変化させる機構を有する場合には、圧縮比の変化量に応じた排気弁閉時期における燃焼室容積ＶＥＶＣをテーブルから求める。排気ＶＴＣ機構２８に加えて圧縮比を変化させる機構をも有する場合には、排気弁閉時期と圧縮比変化量とに応じたマップを検索することにより排気弁閉時期における燃焼室容積を求める。

ステップ７３では、目標当量比ＴＦＢＹＡから図２４に示すテーブルを検索することにより、燃焼室５内の不活性ガスのガス定数ＲＥＸを求める。図２４に示すように、不活性ガスのガス定数ＲＥＸは目標当量比ＴＦＢＹＡが１．０のとき、つまり理論空燃比のとき最も小さく、これより大きくても小さくても大きくなる。

ステップ７４では、排気温度ＴＥＸＨに基づいて燃焼室５の排気弁閉時期ＥＶＣにおける温度ＴＥＶＣを推定する。簡単には排気温度ＴＥＸＨをそのままＴＥＶＣとおけばよい。なお、燃焼室５の排気弁閉時期における温度ＴＥＶＣは、インジェクタ２１の燃料噴射量に応じた熱量により変化するため、このような特性をも加味すれば、ＴＥＶＣの算出精度が向上する。

ステップ７５では、排気圧力ＰＥＸＨに基づいて燃焼室５の排気弁閉時期ＥＶＣにおける圧力ＰＥＶＣを算出する。簡単には排気圧力ＰＥＸＨをＰＥＶＣと置けばよい。

ステップ７６では、燃焼室５の排気弁閉時期ＥＶＣにおける容積ＶＥＶＣ、排気弁閉時期ＥＶＣにおける温度ＴＥＶＣ、排気弁閉時期ＥＶＣにおける圧力ＰＥＶＣ及び不活性ガスのガス定数ＲＥＸから、燃焼室５の排気弁閉時期ＥＶＣにおける不活性ガス量ＭＲＥＳＣＹＬを次式により算出する。

ＭＲＥＳＣＹＬ＝（ＰＥＶＣ×ＶＥＶＣ）／（ＲＥＸ×ＴＥＶＣ） …（２４）
このようにして燃焼室５の排気弁閉時期ＥＶＣにおける不活性ガス量ＭＲＥＳＣＹＬの算出を終了したら図１５に戻り、ステップ６２で吸排気弁１５、１６のオーバーラップ（図では「Ｏ／Ｌ」と略記する）中に排気側から吸気側へ吹き返す不活性ガス量であるオーバーラップ中吹き返し不活性ガス量ＭＲＥＳＯＬを算出する。

この不活性ガス量ＭＲＥＳＯＬの算出については図１７のフローにより説明する。

図１７（図１５ステップ６２のサブルーチン）においてステップ８１では、吸気弁開時期ＩＶＯ［ｄｅｇＢＴＤＣ］と、排気弁閉時期ＥＶＣ［ｄｅｇＢＴＤＣ］、図１６のステップ７４で算出されている燃焼室５の排気弁閉時期ＥＶＣにおける温度ＴＥＶＣを読み込む。

ここで、吸気弁開時期ＩＶＯは、吸気弁閉時期ＩＶＣより吸気弁１５の開き角だけ前の時期となるので、吸気弁閉時期ＩＶＣより吸気弁１５の開き角（予め分かっている）とから求めることができる。

ステップ８２では吸気弁開時期ＩＶＯと排気弁閉時期ＥＶＣとから、吸排気弁のオーバーラップ量ＶＴＣＯＬ［ｄｅｇ］を次式により算出する。

ＶＴＣＯＬ＝ＩＶＯ＋ＥＶＣ …（２５）
例えば、吸気ＶＴＣ機構２７用アクチュエータへの非通電時に吸気弁開時期ＩＶＯが吸気上死点位置にあり、吸気ＶＴＣ機構２７用アクチュエータへの通電時に吸気弁開時期が吸気上死点より進角する特性であり、かつ排気ＶＴＣ機構２８用アクチュエータへの非通電時に排気弁閉時期ＥＶＣが排気上死点にあり、排気弁ＶＴＣ機構２８用アクチュエータへの通電時に排気弁閉時期ＥＶＣが排気上死点より進角する特性である場合には、ＩＶＯとＥＶＣの合計が吸排気弁のオーバーラップ量ＶＴＣＯＬとなる。

ステップ８３では、吸排気弁のオーバーラップ量ＶＴＣＯＬから、図２５に示すテーブルを検索することによりオーバーラップ中の積算有効面積ＡＳＵＭＯＬを算出する。図２５に示すようにオーバーラップ中の積算有効面積ＡＳＵＭＯＬは吸排気弁のオーバーラップ量ＶＴＣＯＬが大きくなるほど大きくなる値である。

ここで、図２６は、吸排気弁のオーバーラップ中の積算有効面積ＡＳＵＭＯＬの説明図であり、横軸はクランク角、縦軸は吸気弁１２と排気弁１５とのそれぞれの開口面積を示している。オーバーラップ中の任意の時点における有効開口面積は、排気弁開口面積と吸気弁開口面積とのうち小さい方とする。オーバーラップ中の全期間における積算有効面積ＡＳＵＭＯＬは、吸気弁１５及び排気弁１６が開いている期間の積分値（図中の斜線部）である。

このようにオーバーラップ中積算有効面積ＡＳＵＭＯＬを算出することで、吸気弁１５と排気弁１６とのオーバーラップ量を１つのオリフィス（流出孔）であると近似することができ、排気系の状態と吸気系の状態とからこの仮想オリフィスを通過するガス流量を簡略的に算出し得る。

ステップ８４では、目標当量比ＴＦＢＹＡと、燃焼室５の排気弁閉時期ＥＶＣにおける温度ＴＥＶＣとから、図２７に示すマップを検索することにより、燃焼室５に残留する不活性ガスの比熱比ＳＨＥＡＴＲを算出する。図２７に示したように、燃焼室に残留する不活性ガスの比熱比ＳＨＥＡＴＲは目標当量比ＴＦＢＹＡが１．０の近傍にあるときが最も小さくなり、それより大きくても小さくても大きくなる。また、目標当量比ＴＦＢＹＡが一定の条件では、燃焼室５の排気弁閉時期ＥＶＣにおける温度ＴＥＶＣが高くなるほど小さくなる。

ステップ８５では過給判定フラグＴＢＣＲＧ及びチョーク判定フラグＣＨＯＫＥを設定する。この過給判定フラグＴＢＣＲＧ及びチョーク判定フラグＣＨＯＫＥの設定については図１８のフローにより説明する。

図１８（図１７ステップ８５のサブルーチン）においてステップ１０１では、吸気圧力センサ４４により検出される吸気圧力ＰＩＮと、図１６のステップ７５で算出されている燃焼室５の排気弁閉時期ＥＶＣにおける圧力ＰＥＶＣを読み込む。

ステップ１０２では、吸気圧力ＰＩＮと、燃焼室５の排気弁閉時期ＥＶＣにおける圧力ＰＥＶＣとから、次式により吸気排気圧力比ＰＩＮＢＹＥＸを算出する。

ＰＩＮＢＹＥＸ＝ＰＩＮ／ＰＥＶＣ …（２６）
この吸気排気圧力比ＰＩＮＢＹＥＸは無名数であり、これと１をステップ１０３で比較する。吸気排気圧力比ＰＩＮＢＹＥＸが１以下の場合には過給無しと判断し、ステップ１０４に進んで過給判定フラグＴＢＣＲＧ（ゼロに初期設定）＝０とする。

吸気排気圧力比ＰＩＮＢＹＥＸが１より大きい場合には過給有りと判断し、ステップ１０５へ進んで過給判定フラグＴＢＣＲＧ＝１とする。

ステップ１０６では、図１４のステップ５１で読み込まれている目標当量比ＴＦＢＹＡから図２８に示すテーブルを検索することにより、混合気の比熱比ＭＩＸＡＩＲＳＨＲを求め、これをステップ１０７で不活性ガスの比熱比ＳＨＥＡＴＲと入れ換える。図２８に示したように、混合気の比熱比ＭＩＸＡＩＲＳＨＲは、目標当量比ＴＦＢＹＡが小さくなるほど大きくなる値である。

ステップ１０６、１０７において、不活性ガスの比熱比ＳＨＥＡＴＲを混合気の比熱比ＭＩＸＡＩＲＳＨＲに置き換えるのは、ターボ過給や慣性過給等の過給時を考慮したものである。すなわち、過給時には吸排気弁のオーバーラップ中のガス流れが吸気系から排気系へ向かう（吹き抜ける）ので、この場合においては、上記の仮想オリフィスを通過するガスの比熱比を不活性ガスの比熱比から混合気の比熱比に変更することで、吹き抜けるガス量を精度良く推定し、内部不活性ガス量を精度良く算出するためである。

ステップ１０８では、図１７のステップ８４または図１８のステップ１０６、１０７で算出している不活性ガスの比熱比ＳＨＥＡＴＲに基づき、最小と最大とのチョーク判定しきい値ＳＬＣＨＯＫＥＬ、ＳＬＣＨＯＫＥＨを次式により算出する。

ＳＬＣＨＯＫＥＬ＝｛２／（ＳＨＥＡＴＲ＋１）｝
＾｛ＳＨＥＡＴＲ／（ＳＨＥＡＴＲ−１）｝
…（２７ａ）
ＳＬＣＨＯＫＥＨ＝｛−２／（ＳＨＥＡＴＲ+１）｝
＾｛−ＳＨＥＡＴＲ／（ＳＨＥＡＴＲ−１）｝
…（２７ｂ）
これらのチョーク判定しきい値ＳＬＣＨＯＫＥＬ、ＳＬＣＨＯＫＥＨは、チョークする限界値を算出している。

ステップ１０８において、（２７ａ）右辺、（２７ｂ）右辺の各累乗計算が困難な場合には、（２７ａ）、（２７ｂ）式の算出結果を、最小チョーク判定しきい値ＳＬＣＨＯＫＥＬのテーブルと最大チョーク判定しきい値ＳＬＣＨＯＫＥＨのテーブルとしてそれぞれエンジンコントローラ３１のメモリに予め記憶しておき、不活性ガスの比熱比ＳＨＥＡＴＲから当該テーブルを検索することにより求めてもよい。

テップ１０９では、吸気排気圧力比ＰＩＮＢＹＥＸが、最小チョーク判定しきい値ＳＬＣＨＯＫＥＬ以上でかつ最大チョーク判定しきい値ＳＬＣＨＯＫＥＨ以下の範囲内にあるか否か、すなわちチョーク状態にないか否かを判定する。吸気排気圧力比ＰＩＮＢＹＥＸが範囲内にある場合にはチョーク無しと判断し、ステップ１１０に進んでチョーク判定フラグＣＨＯＫＥ（ゼロに初期設定）＝０とする。

吸気排気圧力比Ｐ１ＮＢＹＥＸが範囲内にない場合にはチョーク有りと判断し、ステップ１１１に進んでチョーク判定フラグＣＨＯＫＥ＝１とする。

このようにして過給判定フラグとチョーク判定フラグの設定を終了したら図１７に戻り、ステップ８６〜８８で次の４つの場合分けを行う。

〈１〉過給判定フラグＴＢＣＲＧ＝０かつチョーク判定フラグＣＨＯＫＥ＝０のとき
〈２〉過給判定フラグＴＢＣＲＧ＝０かつチョーク判定フラグＣＨＯＫＥ＝０のとき
〈３〉過給判定フラグＴＢＣＲＧ＝０かつチョーク判定フラグＣＨＯＫＥ＝１のとき
〈４〉過給判定フラグＴＢＣＲＧ＝１かつチョーク判定フラグＣＨＯＫＥ＝０のとき
そして、上記〈１〉のときにはステップ８９に進んで、過給無しかつチョーク無し時のオーバーラップ中の平均吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ１を、上記〈２〉のときにはステップ９０に進んで過給無しかつチョーク有り時のオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ２を、上記〈３〉のときにはステップ９１に進んで過給有りかつチョーク無し時のオーバーラップ中の平均吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ３を、上記〈４〉のときにはステップ９２に進んで過給有りかつチョーク有り時の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ４をそれぞれ算出し、算出結果をオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐに移す。

ここで、過給無しかつチョーク無し時のオーバーラップ中吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ１の算出について図１９のフローにより説明する
図１９（図１７ステップ８９のサブルーチン）においてステップ１２１では、図１６のステップ７３、７５で算出されている不活性ガスのガス定数ＲＥＸ、燃焼室５の排気弁閉時期における圧力ＰＥＶＣを読み込む。

ステップ１２２では、不活性ガスのガス定数ＲＥＸと、図１７のステップ８１で読み込まれている燃焼室５の排気弁閉時期における温度ＴＥＶＣとに基づき、後述するガス流量の算出式に用いる密度項ＭＲＳＯＬＤを次式により算出する。

ＭＲＳＯＬＤ＝ＳＱＲＴ｛１／（ＲＥＸ×ＴＥＶＣ）｝ …（２８）
ここで、（２８）式右辺の「ＳＱＲＴ」はすぐ右のカッコ内の値の平方根を計算させる関数である。

なお、密度項ＭＲＳＯＬＤの平方根計算が困難な場合は、（２８）式の算出結果をマップとしてエンジンコントローラ３１のメモリに予め記憶しておき、ガス定数ＲＥＸと燃焼室５の排気弁閉時期における温度ＴＥＶＣとからそのマップを検索することにより求めてもよい。

ステップ１２３では、図１７のステップ８４で算出されている不活性ガスの比熱比ＳＨＥＡＴＲと、図１８のステップ１０２で算出されている吸気排気圧力比ＰＩＮＢＹＥＸとに基づき、後述するガス流量の算出式に用いる圧力差項ＭＲＳＯＬＰを次式により算出する。

ＭＲＳＯＬＰ＝ＳＱＲＴ［ＳＨＥＡＴＲ／（ＳＨＥＡＴＲ−１）
×｛ＰＴＮＢＹＥＸ＾（２／ＳＨＥＡＴＲ）
−ＰＴＮＢＹＥＸ＾（（ＳＨＥＡＴＲ＋１）／ＳＨＥＡＴＲ）｝］ …（２９）
ステップ１２４では、これら密度項ＭＲＳＯＬＤ、圧力差項ＭＲＳＯＬＰと、燃焼室５の排気弁閉時期における圧力ＰＥＶＣとから、過給無しかつチョーク無し時のオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ１を次式（ガス流量の算出式）により算出し、その算出値をステップ１２５でオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐに移す。

ＭＲＥＳＯＬｔｍｐ１＝１．４×ＰＥＶＣ×ＭＲＳＯＬＤ×ＭＲＳＯＬＰ
…（３０）
次に、過給無しかつチョーク有り時の吹き返し不活性ガス流量の算出について図２０のフローにより説明する
図２０（図１７ステップ９０のサブルーチン）においてステップ１３１、１３２では、図１９のステップ１２１、１２２と同様にして、不活性ガスのガス定数ＲＥＸ、燃焼室５の排気弁閉時期における圧力ＰＥＶＣを読み込み、これらから前述の（２８）式により密度項ＭＲＳＯＬＤを算出する。

ステップ１３３では、図１７のステップ８４で算出されている不活性ガスの比熱比ＳＨＥＡＴＲに基づき、チョーク時圧力差項ＭＲＳＯＬＰＣを次式により算出する。

ＭＲＳＯＬＰＣ＝ＳＱＲＴ［ＳＨＥＡＴＲ×｛２／（ＳＨＥＡＴＲ＋１）｝ ＾｛（ＳＨＥＡＴＲ＋１）／〔ＳＨＥＡＴＲ−１）｝］
…（３１）
なお、（３１）式の累乗計算と平方根計算とが困難な場合には、（３１）式の算出結果を、チョーク時圧力差項ＭＲＳＯＬＰＣのテーブルとしてエンジンコントローラ３１のメモリに予めに記憶しておき、不活性ガスの比熱比ＳＨＥＡＴＲからそのテーブルを検索することにより求めてもよい。

ステップ１３４では、これら密度項ＭＲＳＯＬＤ、チョーク時圧力差項ＭＲＳＯＬＰＣと、燃焼室５の排気弁閉時期における圧力ＰＥＶＣとから、過給無しかつチョーク有り時のオーバーラップ中吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ２を次式により算出し、その算出値をステップ１３５でオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐに移す。

ＭＲＥＳＯＬｔｍｐ２＝ＰＥＶＣ×ＭＲＳＯＬＤ×ＭＲＳＯＬＰＣ
…（３２）
次に、過給有りかつチョーク無し時の吹き返しガス流量の算出について図２１のフローにより説明する
図２１（図１７ステップ９１のサブルーチン）においてステップ１４１では、吸気圧力センサ４４により検出される吸気圧力ＰＩＮを読み込む。

ステップ１４２では、図１８のステップ１０６、１０７で算出されている不活性ガスの比熱比ＳＨＥＡＴＲと、図１８のステップ１０２で算出されている吸気排気圧力比ＰＩＮＢＹＥＸとから、過給時圧力差項ＭＲＳＯＬＰＴを次式により算出する。

ＭＲＳＯＬＰＴ＝ＳＱＲＴ［ＳＨＥＡＴＲ／（ＳＨＥＡＴＲ−１）
×｛ＰＩＮＢＹＥＸ＾（−２／ＳＨＥＡＴＲ）
−ＰＩＮＢＹＥＸ＾（−（ＳＨＥＡＴＲ＋１）／ＳＨＥＡＴＲ）｝］ …（３３）
なお、（３３）式の累乗計算と平方根計算とが困難な場合は、（３３）式の算出結果を、過給時圧力差項ＭＲＳＯＬＰＴのマップとしてエンジンコントローラ３１のメモリに予め記憶しておき、不活性ガスの比熱比ＳＨＥＡＴＲと吸気排気圧力比ＰＩＮＢＹＥＸとからそのマップを検索することにより求めてもよい。

ステップ１４３では、この過給時圧力差項ＭＲＳＯＬＰＴと吸気圧力ＰＩＮとに基づいて、過給有りかつチョーク無し時のオーバーラップ中吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ３を次式により算出し、その算出値をステップ１４４でオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐに移す。

ＭＲＥＳＯＬｔｍｐ３＝−０．１５２×ＰＩＮ×ＭＲＳＯＬＰＴ …（３４）
ここで、（３４）式の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ３は負の値とすることで、オーバーラップ中に吸気系から排気系へ吹き抜ける混合気のガス流量を表すことができる。

次に、過給有りかつチョーク有り時のオーバーラップ中吹き返し不活性ガス流量の算出について図２２のフローにより説明する
図２２（図１７ステップ９２のサブルーチン）においてステップ１５１、１５２では、図２１のステップ１４１と同じく吸気圧力センサ４４により検出される吸気圧力ＰＩＮを読み込むと共に、図２０のステップ１３２と同じくチョーク時圧力差項ＭＲＳＯＬＰＣを前述の（３１）式により算出する。

ステップ１５３では、このチョーク時圧力差項ＭＲＳＯＬＰＣと吸気圧力ＰＩＮとに基づいて、過給有りかつチョーク有り時のオーバーラップ中吹き返しガス流量ＭＲＥＳＯＬｔｍｐ４を次式により算出し、その算出値をステップ１５４でオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐに移す。

ＭＲＥＳＯＬｔｍｐ４＝−０．１０８×ＰＩＮ×ＭＲＳＯＬＰＣ …（３５）
ここで、（３５）式の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐ４も、ＭＲＥＳＯＬｔｍｐ３と同様、負の値とすることで、オーバーラップ中に吸気側から排気側へ吹き抜ける混合気のガス流量を表すことができる。

このようにして、過給の有無とチョークの有無との組み合わせにより場合分けした、オーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐの算出を終了したら図１７に戻り、ステップ９３においてこのオーバーラップ中の吹き返し不活性ガス流量ＭＲＥＳＯＬｔｍｐとオーバーラップ期間中の積算有効面積ＡＳＵＭＯＬとから、次式によりオーバーラップ中の吹き返し不活性ガス量ＭＲＥＳＯＬを算出する。

ＭＲＥＳＯＬ＝（ＭＲＥＳＯＬｔｍＰ×ＡＳＵＭＯＬ×６０）
／（ＮＲＰＭ×３６０） …（３６）
このようにしてオーバーラップ中の吹き返し不活性ガス量ＭＲＥＳＯＬの算出を終了したら図１５に戻り、ステップ６３において燃焼室５内の排気弁閉時期ＥＶＣにおける不活性ガス量ＭＲＥＳＣＹＬと、このオーバーラップ中吹き返しガス量ＭＲＥＳＯＬとを加算して、つまり次式により内部不活性ガス量ＭＲＥＳを算出する。

ＭＲＥＳ＝ＭＲＥＳＣＹＬ＋ＭＲＥＳＯＬ …（３７）
前述のように、過給有り時にはオーバーラップ中吹き返し不活性ガス流量（ＭＲＥＳＯＬｔｍｐ３、ＭＲＥＳＯＬｔｍｐ４）が負となるため、上記（３６）式のオーバーラップ中の吹き返し不活性ガス量ＭＲＥＳＯＬも負となり、このとき（３７）式によれば、オーバーラップ中の吹き返し不活性ガス量ＭＲＥＳＯＬの分だけ内部不活性ガス量が減じられる。

このようにして内部不活性ガス量ＭＲＥＳの算出を終了したら図１４に戻り、ステップ５４においてこの内部不活性ガス量ＭＲＥＳと、目標当量比ＴＦＢＹＡとを用いて、次式により内部不活性ガス率ＭＲＥＳＦＲ（燃焼室５内の総ガス量に対する内部不活性ガス量の割合）を算出する。

ＭＲＥＳＦＲ＝ＭＲＥＳ
／｛ＭＲＥＳ＋ＭＡＣＹＬ×（１＋ＴＦＢＹＡ／１４．７）｝
…（３８）
これで内部不活性ガス率ＭＲＥＳＦＲの算出を総て終了する。

このように本実施形態によれば、内部不活性ガス量ＭＲＥＳを、燃焼室５の排気弁閉時期における不活性ガス量ＭＲＥＳＣＹＬと、吸排気弁のオーバーラップ中の吹き返しガス量ＭＲＥＳＯＬとで構成し（図１５のステップ６３参照）、この場合に、燃焼室５の排気弁閉時期における温度ＴＥＶ及び圧力ＰＥＶＣを算出し（図１６のステップ７４、７５）、これら温度ＴＥＶＣ、圧力ＰＥＶＣと不活性ガスのガス定数ＲＥＸとに基づいて状態方程式（上記（２４）式）により燃焼室５の排気弁閉時期における不活性ガス量ＭＲＥＳＣＹＬを算出する（図１６のステップ７６参照）ようにしたので、特に、燃焼室５内部の状態量（ＰＥＶＣ、ＶＥＶＣ、ＴＥＶＣ）が刻々と変化する過渡運転時においても、運転条件に関わらず精度良く燃焼室５の排気弁閉時期における不活性ガス量ＭＲＥＳＣＹＬを算出（推定）できる。

また、燃焼室５の排気弁閉時期における温度ＴＥＶＣ及び圧力ＰＥＶＣ、不活性ガスのガス定数ＲＥＸ及び比熱比ＳＨＥＡＴＲ、吸気圧力ＰＩＮに基づいてオーバーラップ中の吹き返し不活性ガス流量（ＭＲＥＳＯＬｔｍｐ１、ＭＲＥＳＯＬｔｍｐ２）を算出し（図１９、図２０参照）、このガス流量にオーバーラップ中の積算有効面積ＡＳＵＭＯＬを乗算して、オーバーラップ中の吹き返しガス量ＭＲＥＳＯＬを算出する（図１７のステップ９３参照）ようにしたので、精度良くオーバーラップ中吹き返しガス量ＭＲＥＳＯＬを算出（推定）できる。

このように、燃焼室５の排気弁閉時期における不活性ガス量ＭＲＥＳＣＹＬ、オーバーラップ中吹き返しガス量ＭＲＥＳＯＬとも精度良く算出（推定）できると、これらの和である内部不活性ガス量ＭＲＥＳも精度良く算出（推定）できることになり、この精度良く推定することが可能となった内部不活性ガス量ＭＲＥＳに基づいて算出される内部不活性ガス率ＭＲＥＳＦＲを、点火時期の算出に用いる燃焼室５内の吸気弁閉時期ＩＶＣにおける温度ＴＩＮＩに活かすことで（図５のステップ１３参照）、燃焼室５内の吸気弁閉時期ＩＶＣにおける温度ＴＩＮＩを精度良く算出できる。また、精度良く推定することが可能となった内部不活性ガス量ＭＲＥＳを、燃料噴射量、バルブ開閉タイミング（オーバーラップ量）などにも活かすことで、エンジンを適切に制御することが可能である。

また、不活性ガスのガス定数ＲＥＸや不活性ガスの比熱比ＳＨＥＡＴＲは、目標当量比ＴＦＢＹＡに応じた値としているので（図２４、図２７参照）、理論空燃比を外れた空燃比での運転時（例えば理論空燃比よりもリーンな空燃比で運転を行うリーン運転時、冷間始動時のようにエンジンが元々不安定な状態を安定させるために理論空燃比の空燃比よりもリッチ側の空燃比で運転するエンジン始動直後、同じく大きな出力が要求されるために理論空燃比の空燃比よりもリッチ側の空燃比で運転する全負荷運転時）にも、燃焼室５の排気弁閉時期における不活性ガス量ＭＲＥＳＣＹＬ、オーバーラップ中吹き返しガス量ＭＲＥＳＯＬ、これらの合計である内部不活性ガス量ＭＲＥＳ、これに基づく内部不活性ガス率ＭＲＥＳＦＲを精度良く算出できる。

また、オーバーラップ期間の積算有効面積ＡＳＵＭＯＬを仮想オリフィスの面積とし、この仮想オリフィスを排気が燃焼室５から吸気系へと吹き抜けると仮定しているので、オーバーラップ中の吹き返し不活性ガス量ＭＲＥＳＯＬの算出が簡略化されている。

次に、図２９は初期燃焼期間における平均乱れ強さＳＴ１を算出するためのものである。

図２９（図１０のステップ１６９のサブルーチン）においてステップ２０１ではエンジン回転速度ＮＲＰＭ［ｒｐｍ］を読み込む。

ステップ２０２、２０３では吸気弁開時期ＩＶＯの状態での吸気弁通過平均流速Ｖｏ、吸気弁閉時期ＩＶＣの状態での吸気弁通過平均流速Ｖｃをそれぞれ算出する。吸気弁開時期ＩＶＯの状態での吸気弁通過平均流速Ｖｏの算出については図３０のフローにより、吸気弁閉時期ＩＶＣの状態での吸気弁通過平均流速Ｖｃの算出については図３４のフローにより説明する。

まず、一方の吸気弁通過平均流速Ｖｏから説明すると、図３０（図２９のステップ２０２のサブルーチン）においてステップ２１１では吸気弁開時期ＩＶＯの吸気ポート圧Ｐｏ［Ｐａ］、吸気弁開時期ＩＶＯの状態での開弁クランク角区間θｏ［ｄｅｇ］、吸気弁開時期ＩＶＯの状態での吸気弁リフト量Ｌｏ［ｍｍ］、エンジン回転速度ＮＲＰＭ［ｒｐｍ］を読み込む。

ここで、吸気弁開時期の吸気ポート圧Ｐｏとしては吸気弁開時期ＩＶＯにおいて圧力センサ４４により検出されるコレクタ内圧力ＰＣＯＬをサンプリングすればよい。

ＶＥＬ機構５１を備えるエンジンでは、吸気弁開時期ＩＶＯの状態での吸気弁リフト量Ｌｏ及び開弁クランク角区間θｏと、吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃ及び開弁クランク角区間θｃとが相違することがあるので、これらθｏ、Ｌｏ、θｃ、Ｌｃを次のようにして算出する。すなわち、図３１は吸気弁開時期ＩＶＯにはＶＥＬ機構５１が作動しており従って吸気弁１５の弁リフトが小さい状態（実線参照）にあったものが、吸気弁開時期の直後にＶＥＬ機構５１が非作動となり吸気弁閉時期ＩＶＣの前に吸気弁１５の弁リフトが大きい状態（一点鎖線参照）へと移行している場合を示している。

ここでいう吸気弁リフト量は、開弁期間中の平均のリフト量のことであり、吸気弁開時期ＩＶＯの状態（ＶＥＬ機構５１の作動時）にあっては吸気弁リフト量Ｌｏは小さな側のリフト特性の開弁期間中の平均のリフト量（つまり図示のＬａ）、吸気弁閉時期ＩＶＣの状態（ＶＥＬ機構５１の非作動時）にあっては吸気弁リフト量Ｌｃは大きな側のリフト特性の開弁期間中の平均のリフト量（つまり図示のＬｂ）である。

また、吸気弁開時期ＩＶＯの状態（ＶＥＬ機構５１の作動時）にあっては開弁クランク角区間θｏは小さな側のリフト特性の開弁クランク角区間（つまり図示のθａ）、吸気弁閉時期ＩＶＣの状態（ＶＥＬ機構５１の非作動時）にあっては開弁クランク角区間θｃは大きな側のリフト特性の開弁クランク角区間（つまり図示のθｂ）である。

図示の大小二つの吸気弁リフト特性はＶＥＬ機構５１と吸気弁用カムの仕様により定まるので、図示の値Ｌａ、Ｌｂ、θａ、θｂは予め図上で計算することができる。従って、吸気弁開時期ＩＶＯの状態での吸気弁リフト量Ｌｏ、吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃ、吸気弁開時期ＩＶＯの状態での開弁クランク角区間θｏ、吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃとしては、ＶＥＬ機構５１の作動、非作動の各状態をみて、ＶＥＬ機構５１の作動時にはＬａとθａを、またＶＥＬ機構５１の非作動時にはＬｂとθｂを用いればよい（Ｌａ、Ｌｂ、θａ、θｂはいずれも一定値）。

あるいは、計算によって吸気弁開時期の状態での吸気弁リフト量、吸気弁閉時期の状態での吸気弁リフト量、吸気弁開時期の状態での開弁クランク角区間、吸気弁閉時期の状態での開弁クランク角区間を求めることもできる。例えば、吸気弁１５が開いているクランク角区間にわたって瞬時吸気弁リフト量を積算し、その積算値を開弁クランク角区間で除算して平均値を求め、その平均値を吸気弁開時期の状態での吸気弁リフト量Ｌｏや吸気弁開時期の状態での吸気弁リフト量Ｌｃとすればよい。また、吸気弁開時期のクランク角と、吸気弁閉時期のクランク角を検出し、その間のクランク角区間を吸気弁開時期の状態での開弁クランク角区間θｏや吸気弁閉時期の状態での開弁クランク角区間θｃとすればよい。

一方、吸気弁リフト特性が図３１のいずれか一方の特性にあって吸気弁開時期ＩＶＯより吸気弁閉時期ＩＶＣまでのクランク角区間の間にＶＥＬ機構５１の作動状態が変わらないとき、例えば吸気弁開時期ＩＶＯより吸気弁閉時期ＩＶＣまでのクランク角区間の間中ずっとＶＥＬ機構５１が作動していれば、吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃと、吸気弁開時期ＩＶＯの状態での吸気弁リフト量Ｌｏとは共に図３１に示されるＬａに等しく、吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃと、吸気弁開時期ＩＶＯの状態での開弁クランク角区間θｏとは共に図３１に示されるθａに等しい。この反対に、吸気弁開時期ＩＶＯより吸気弁閉時期ＩＶＣまでのクランク角区間の間中ずっとＶＥＬ機構５１が非作動であるときには、吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃと、吸気弁開時期ＩＶＯの状態での吸気弁リフト量Ｌｏとは共に図３１に示されるＬｂに等しく、吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃと、吸気弁開時期ＩＶＯの状態での開弁クランク角区間θｏとは共に図３１に示されるθｂに等しい。

ステップ２１２では吸気行程中の吸気ポート圧が吸気弁開時期の吸気ポート圧Ｐｏを維持すると仮定したときに燃焼室５へと流入する吸入空気量（燃焼室５へと流入する吸入空気量を以下「シリンダ吸入空気量」という。）Ｑｏ［ｍｍ³］を吸気弁開時期の吸気ポート圧Ｐｏに比例させて、つまり次式により算出する。

Ｑｏ＝Ｐｏ×ｆ２１（θｏ、ＮＲＰＭ） …（５１）
ここで、（５１）式の補正係数ｆ２１（θｏ、ＮＲＰＭ）は吸気弁開時期の状態での吸気弁リフト量θｏと回転速度ＮＲＰＭの関数であることを表している。この補正係数ｆ２１（θｏ、ＮＲＰＭ）は吸気弁開時期の吸気ポート圧Ｐｏが同じでも、吸気弁開時期の状態での開弁クランク角区間θｏや回転速度ＮＲＰＭの相違によってシリンダ吸入空気量Ｑｏが変化するので、これを考慮するものである。補正係数ｆ２１（θｏ、ＮＲＰＭ）を、関数で与えるのではなくマップ検索により求めるようにしてもかまわない。図３２はこの場合のマップの特性であり、図示のように補正係数ｆ２１（θｏ、ＮＲＰＭ）は吸気弁開時期の状態での開弁クランク角区間θｏを大きくするほど、また回転速度ＮＲＰＭが高くなるほど大きくなる値である。

ステップ２１３では吸気弁開時期ＩＶＯの状態での吸気弁開時期ＩＶＯより吸気弁閉時期ＩＶＣまでの平均吸気弁開口面積Ｓｏ［ｍｍ²］を吸気弁開時期の状態での吸気弁リフト量Ｌｏの関数として、つまり次式により算出する。

Ｓｏ＝ｆ２２（Ｌｏ） …（５２）
ここで、吸気弁開時期の状態での平均吸気弁開口面積Ｓｏを、関数で与えるのではなくテーブル検索により求めるようにしてもかまわない。図３３はこの場合のテーブル特性であり、吸気弁開時期の状態での平均吸気弁開口面積Ｓｏは図示のように吸気弁開時期の状態での吸気弁リフト量Ｌｏが大きくなるほど大きくなる値である。

ステップ２１４では吸気弁開時期ＩＶＯの状態での吸気弁通過平均流速Ｖｏ［ｍｍ／ｓ］を次式により算出する。

Ｖｏ＝Ｑｏ／（Ｓｏ・ＮＲＰＭ） …（５３）
すなわち、吸気弁開時期の状態での吸気弁通過平均流速Ｖｏは、吸気弁開時期の状態でのシリンダ吸入空気量Ｑｏに比例して大きくなり、吸気弁開時期の状態での平均吸気弁開口面積Ｓｏ及び回転速度ＮＲＰＭに反比例して小さくなる値である。

このようにして吸気弁開時期の状態での吸気弁通過平均流速Ｖｏの算出を終了したら図２９のステップ２０３に戻り、吸気弁閉時期ＩＶＣの状態での吸気弁通過平均流速Ｖｃを算出する。この吸気弁通過平均流速Ｖｃの算出については図３４のフローにより説明する。吸気弁閉時期の状態での吸気弁通過平均流速Ｖｃの算出方法は、図３０により前述した吸気弁開時期の状態での吸気弁通過平均流速Ｖｏの算出方法と同様である。

図３４（図２９のステップ２０３のサブルーチン）においてステップ２２１では吸気弁閉時期ＩＶＣの吸気ポート圧Ｐｃ［Ｐａ］、吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃ［ｄｅｇ］、吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃ［ｍｍ］、エンジン回転速度ＮＲＰＭ［ｒｐｍ］を読み込む。

ここで、吸気弁閉時期の吸気ポート圧Ｐｃとしては吸気弁閉時期ＩＶＣにおいて圧力センサ４４により検出されるコレクタ内圧力ＰＣＯＬをサンプリングすればよい。

吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃ及び吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃとしては、図３１で前述したように、吸気弁閉時期ＩＶＣの前にＶＥＬ機構５１が非作動となり、吸気弁閉時期ＩＶＣにおいて吸気弁リフト特性が大きい側（一点鎖線参照）に移行していれば、図示のＬｂ及びθｂを用いる。また、吸気弁開時期ＩＶＯより吸気弁閉時期ＩＶＣまでのクランク角区間の間中ずっとＶＥＬ機構５１が作動しているときには、吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃは図示のＬａに、吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃは図示のθａになる。この反対に、吸気弁開時期ＩＶＯより吸気弁閉時期ＩＶＣまでのクランク角区間の間中ずっとＶＥＬ機構５１が非作動であるときには、吸気弁閉時期ＩＶｃの状態での吸気弁リフト量Ｌｃは図示のＬｂに、吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃは図示のθｂになる。

ステップ２２２では吸気行程中の吸気ポート圧が吸気弁閉時期の吸気ポート圧Ｐｃを維持すると仮定したときのシリンダ吸入空気量Ｑｃを吸気弁閉時期の吸気ポート圧Ｐｃに比例させて、つまり次式により算出する。

Ｑｃ＝Ｐｃ×ｆ２１（θｃ、ＮＲＰＭ） …（５４）
ここで、（５４）式の補正係数ｆ２１（θｃ、ＮＲＰＭ）は上記（５１）式と同じ関数である。補正係数ｆ２１（θｃ、ＮＲＰＭ）を、関数で与えるのではなくマップ検索により求めるようにしてもかまわない。図３２はこの場合のマップの特性であり、図示のように補正係数ｆ２１（θｃ、ＮＲＰＭ）は吸気弁閉時期の状態での開弁クランク角区間θｃを大きくするほど、また回転速度ＮＲＰＭが高くなるほど大きくなる値である。

ステップ２２３では吸気弁閉時期ＩＶＣの状態での吸気弁開時期より吸気弁閉時期までの平均吸気弁開口面積Ｓｃ［ｍｍ²］を吸気弁閉時期の状態での吸気弁リフト量Ｌｃの関数として、つまり次式により算出する。

Ｓｃ＝ｆ２２（Ｌｃ） …（５５）
ここで、吸気弁閉時期の状態での平均吸気弁開口面積Ｓｃを、関数で与えるのではなくテーブル検索により求めるようにしてもかまわない。図３３はこの場合のテーブル特性であり、図示のように吸気弁閉時期の状態での吸気弁リフト量Ｌｃが大きくなるほど大きくなる値である。

ステップ２２４では吸気弁閉時期ＩＶＣの状態での吸気弁通過平均流速Ｖｃ［ｍｍ／ｓ］を次式により算出する。

Ｖｃ＝Ｑｃ／（Ｓｃ・ＮＲＰＭ） …（５６）
すなわち、吸気弁閉時期の状態での吸気弁通過流速Ｖｃは、吸気弁閉時期の状態でのシリンダ吸入空気量Ｑｃに比例して大きくなり、吸気弁閉時期の状態での平均吸気弁開口面積Ｓｃ及び回転速度ＮＲＰＭに反比例して小さくなる値である。

このようにして吸気弁閉時期の状態での吸気弁通過平均流速Ｖｃの算出を終了したら図２９に戻りステップ２０４、２０５で吸気弁開時期ＩＶＯの状態での初期燃焼期間における平均乱れ強さＵｏ１、吸気弁閉時期ＩＶＣの状態での初期燃焼期間における平均乱れ強さＵｃ１をそれぞれ次式により算出する。

Ｕｏ１＝ｆ２３（Ｖｏ、ＮＲＰＭ）×Ｋｃ×Ｋｔ …（５７）
Ｕｃ１＝ｆ２３（Ｖｃ、ＮＲＰＭ）×Ｋｃ×Ｋｔ …（５８）
ただし、Ｋｃ：スワールコントロールバルブ補正係数、
Ｋｔ：タンブルコントロールバルブ補正係数、
（５７）、（５８）式の初期燃焼期間における平均乱れ強さ基本値ｆ２３（Ｖｏ、ＮＲＰＭ）、ｆ２３（Ｖｃ、ＮＲＰＭ）は吸気弁通過平均流速（Ｖｏ、Ｖｃ）と回転速度ＮＲＰＭの関数であることを表している。この初期燃焼期間における平均乱れ強さ基本値ｆ２３（Ｖｏ、ＮＲＰＭ）、ｆ２３（Ｖｃ、ＮＲＰＭ）は吸気弁通過平均流速Ｖｏ、Ｖｃや回転速度ＮＲＰＭの相違によって初期燃焼期間における平均乱れ強さが変化するので、これを考慮するものである。初期燃焼期間における平均乱れ強さ基本値ｆ２３（Ｖｏ、ＮＲＰＭ）、ｆ２３（Ｖｃ、ＮＲＰＭ）を、関数で与えるのではなくマップ検索により求めるようにしてもかまわない。図３５はこの場合のマップの特性であり、図示のように初期燃焼期間における平均乱れ強さ基本値ｆ２３（Ｖｏ、ＮＲＰＭ）、ｆ２３（Ｖｃ、ＮＲＰＭ）は吸気弁通過平均流速Ｖｏ、Ｖｃが大きくなるほど、また回転速度ＮＲＰＭが高くなるほど大きくなる値である。

ここで、図１には図示しなかったが、低負荷状態での燃焼状態の改善ため、燃焼室５のシリンダ軸を中心として周方向（シリンダ壁に沿った方向）にガス流動（スワール）を与えるためのスワールコントロールバルブや、シリンダ軸に直交する軸を中心として旋回するガス流動（タンブル）を与えるためのタンブルコントロールバルブが吸気ポート４に設けられることがあり、上記（５７）、（５８）式のスワールコントロールバルブ補正係数Ｋｃ、タンブルコントロールバルブ補正係数Ｋｔはこれら付加的なガス流動を考慮するものである。すなわち、低負荷域でスワールコントロールバルブ、タンブルコントロールバルブが作動して燃焼室５内のガス流動を強化するときには、初期燃焼期間における平均乱れ強さが大きくなるので、補正係数Ｋｃ、Ｋｔとして１．０を超える値を与えるのである。実際には補正係数Ｋｃ、Ｋｔは適合値である。

なお、制御上は汎用性を持たせるためスワールコントロールバルブ、タンブルコントロールバルブの両方を扱っているが、実際のエンジンではこれらスワールコントロールバルブ、タンブルコントロールバルブのうちいずれか一方しか備えられない。いずれが備えられるのかはエンジン仕様によりわかるので、不要なほうの補正係数は１．０とする。

ステップ２０６では加重平均係数ｋ１を算出する。この加重平均係数ｋ１の算出については図３６のフローにより説明する。

図３６（図２９のステップ２０６のサブルーチン）においてステップ２３１では吸気弁開時期ＩＶＯの吸気ポート圧Ｐｏ［Ｐａ］、吸気弁開時期ＩＶＯの状態での開弁クランク角区間θｏ［ｄｅｇ］、吸気弁開時期ＩＶＯの状態での吸気弁リフト量Ｌｏ［ｍｍ］、吸気弁閉時期ＩＶＣの吸気ポート圧Ｐｃ［Ｐａ］、吸気弁閉時期ＩＶＣの状態での開弁クランク角区間θｃ［ｄｅｇ］、吸気弁閉時期ＩＶＣの状態での吸気弁リフト量Ｌｃ［ｍｍ］、平均吸気ポート圧Ｐａｖｅ［Ｐａ］、平均開弁クランク角区間θａｖｅ［ｄｅｇＡＴＤＣ］、平均吸気弁リフト量Ｌａｖｅ［ｍｍ］を読み込む。

ここで、平均吸気ポート圧ａｖｅ、平均吸気弁リフト量Ｌａｖｅ、平均開弁クランク角区間θａｖｅは次式により算出している。

Ｐａｖｅ＝ΣＰ（ｔ）／（ｔｃ−ｔｏ） …（５９）
Ｌａｖｅ＝ΣＬ（ｔ）／（ｔｃ−ｔｏ） …（６０）
θａｖｅ＝Σθ（ｔ）／（ｔｃ−ｔｏ） …（６１）
ただし、Ｐ（ｔ）：瞬時吸気ポート圧、
Ｌ（ｔ）：瞬時吸気弁リフト量、
θ（ｔ）：瞬時開弁クランク角区間、
ｔｃ ：吸気弁閉時期の時刻、
ｔｏ ：吸気弁開時期の時刻、
（５９）式の瞬時吸気ポート圧Ｐ（ｔ）としては圧力センサ４４により検出されるコレクタ内圧力ＰＣＯＬを用いる。

ＶＥＬ機構５１ではＶＥＬ機構アクチュエータに与える制御量に対して吸気弁リフト量と開弁クランク角区間とが一義的に定まり、時間と共にＶＥＬ機構アクチュエータに与える制御量が変化するので、ＶＥＬ機構アクチュエータに与える制御量に対する吸気弁リフト量と開弁クランク角区間の各テーブルを予め用意しておき、そのときにＶＥＬ機構アクチュエータに与える制御量から各テーブルを検索して求めた値をそれぞれ（６０）式の瞬時吸気弁リフト量Ｌ（ｔ）、（６１）式の瞬時開弁クランク角区間θ（ｔ）として用いる。（６０）式の瞬時吸気弁リフト量Ｌ（ｔ）としてはこの瞬時吸気弁リフト量Ｌ（ｔ）を検出するリフトセンサを設けておき、このリフトセンサからの信号を用いることでもかまわない。

ステップ２３２では吸気弁閉時期の吸気ポート圧Ｐｃと、吸気弁開時期の吸気ポート圧Ｐｏを比較する。吸気弁閉時期の吸気ポート圧Ｐｃと、吸気弁開時期の吸気ポート圧Ｐｏとが等しくないときにはステップ２３３に進み吸気ポート圧に関する加重平均係数（第１係数）Ｋｐを次式により算出する。

Ｋｐ＝（Ｐａｖｅ−Ｐｏ）／（Ｐｃ−Ｐｏ） …（６２）
一方、吸気弁閉時期の吸気ポート圧Ｐｃと、吸気弁開時期の吸気ポート圧Ｐｏとが等しいときにはステップ２３２よりステップ２３４に進み、Ｋｐ＝１／２とする。

ステップ２３５では吸気弁閉時期の状態での開弁クランク角区間θｃと、吸気弁開時期の状態での開弁クランク角区間θｏとを比較する。吸気弁閉時期の状態での開弁クランク角区間θｃと、吸気弁開時期の状態での開弁クランク角区間θｏとが等しくないときにはステップ２３６に進み次式により開弁クランク角区間に関する加重平均係数（第２係数）Ｋθを次式により算出する。

Ｋθ＝（θａｖｅ−θｏ）／（θｃ−θｏ） …（６３）
一方、吸気弁閉時期の状態での開弁クランク角区間θｃと吸気弁開時期の状態での開弁クランク角区間θｏとが等しいときにはステップ２３５よりステップ２３７に進み、Ｋθ＝１／２とする。

ステップ２３８では吸気弁閉時期の状態での吸気弁リフト量Ｌｃと、吸気弁開時期の状態での吸気弁リフト量Ｌｏとを比較する。吸気弁閉時期の状態での吸気弁リフト量Ｌｃと吸気弁開時期の状態での吸気弁リフト量Ｌｏとが等しくないときにはステップ２３９に進み次式により吸気弁リフト量に関する加重平均係数（第３係数）ＫＬを次式により算出する。

ＫＬ＝（Ｌａｖｅ−Ｌｏ）／（Ｌｃ−Ｌｏ） …（６４）
一方、吸気弁閉時期の状態での吸気弁リフト量Ｌｃと、吸気弁開時期の状態での吸気弁リフト量Ｌｏとが等しいときにはステップ２３８よりステップ２４０に進み、ＫＬ＝１／２とする。

ステップ２４１ではこのようにして算出した３つの係数Ｋｐ、Ｋθ、ＫＬを乗算してつまり次式により加重平均係数ｋ１を算出する。

ｋ１＝４×Ｋｐ×Ｋθ×ＫＬ …（６５）
このようにして加重平均係数ｋ１の算出を終了したら図２９のステップ２０７に戻りこの加重平均係数ｋ１を用いて吸気弁開時期の状態での初期燃焼期間における平均乱れ強さＵｏ１と、吸気弁閉時期の状態での初期燃焼期間における平均乱れ強さＵｃ１とを加重平均した値に回転速度ＮＲＰＭを乗算した値を初期燃焼期間における平均乱れ強さＳＴ１とする。つまり次式により初期燃焼期間における平均乱れ強さＳＴ１を算出する。

ＳＴ１＝｛ｋ１×Ｕｏ１＋（１−ｋ１）×Ｕｃ１｝×ＮＲＰＭ …（６６）
（６６）式において回転速度ＮＲＰＭを露わにし、この露わにした回転速度ＮＲＰＭに比例させて初期燃焼期間における平均乱れ強さＳＴ１を求めるようにしたのは実験結果によるものである。もちろん、（６６）式に代えて次式としてもかまわない。

ＳＴ１＝｛ｋ１×Ｕｏ１＋（１−ｋ１）×Ｕｃ１｝ …（６７）
次に、図３７は主燃焼期間における平均乱れ強さＳＴ２を算出するためのものである。主燃焼期間における平均乱れ強さＳＴ２の算出方法そのものは、前述した初期燃焼期間における平均乱れ強さＳＴ２の算出方法と同様である。

図３７（図１２のステップ１８９のサブルーチン）においてステップ２５１ではエンジン回転速度ＮＲＰＭ［ｒｐｍ］を読み込む。

ステップ２５２、２５３では図２９のステップ２０２、２０３と同じに吸気弁開時期ＩＶＯの状態での吸気弁通過平均流速Ｖｏ、吸気弁閉時期ＩＶＣの状態での吸気弁通過平均流速Ｖｃをそれぞれ算出する。吸気弁開時期ＩＶＯの状態での吸気弁通過平均流速Ｖｏの算出については図３０のフローにより、吸気弁閉時期ＩＶＣの状態での吸気弁通過平均流速Ｖｃの算出については図３４のフローにより前述した。

ステップ２５４では、図２９では算出しなかったスキッシュ補正係数Ｋｑ（スキッシュによる火炎速度補正量）を算出する。

ここで、スキッシュとは、ピストン６頭部とシリンダヘッド下面との間の狭い隙間に存在する混合気が、圧縮行程の終わり近くで燃焼室５の主要空間に押し込まれることにより生じるシリンダ半径方向のガス流動のことで、燃焼室５の形状によってその大きさや発生する位置が決まっている。そして、燃焼室５内の乱れ強さに影響する。

そこで、本発明者が各種の実験を行って調べたところスワール、タンブル、スキッシュが燃焼室５内の乱れ強さに与える影響がそれぞれ違っていることが分かってきたので、本発明はこの実験結果を初期燃焼期間における平均乱れ強さＳＴ１及び主燃焼期間における平均乱れ強さＳＴ２の算出にそれぞれ反映させたものである。

すなわち、スワール、タンブルは初期燃焼期間における平均乱れ強さＳＴ１には影響するものの、主燃焼期間における平均乱れ強さＳＴ２にはほとんど影響しないことが判明したので、スワールコントロールバルブ補正係数Ｋｃ及びタンブルコントロールバルブ補正係数Ｋｔは上記（５７）、（５８）式で前述したように初期燃焼期間における平均乱れ強さＳＴ１の算出にのみ導入し、後述するように主燃焼期間における平均乱れ強さＳＴ２の算出には導入しない。

一方、スキッシュはスワール、タンブルと異なり、初期燃焼期間における平均乱れ強さＳＴ１には影響せず、主燃焼期間における平均乱れ強さＳＴ２に影響することが判明したので、スキッシュ補正係数Ｋｑは主燃焼期間における平均乱れ強さＳＴ２の算出にのみ導入し、上記（５７）、（５８）式で前述したように初期燃焼期間における平均乱れ強さＳＴ１の算出には導入していない。これは、スキッシュが生じるのは圧縮行程後半つまり主燃焼期間であり、初期燃焼期間においてはスキッシュが生じないためである。

また、本発明ではスキッシュ補正係数Ｋｑをスキッシュエリア体積やスキッシュエリアからのガス吹き出し方向といったエンジンの仕様に少なくとも基づいて算出する。

これに対して、スキッシュ補正係数Ｋｑそのものは先願装置（特願２００３−９５７５１参照）において既に開示しているものの、その先願装置の出願時点ではまだ詳しい実験を行っていなかったので、スキッシュ補正係数Ｋｑとしてはスワールコントロールバルブ補正係数Ｋｃ、タンブルコントロールバルブ補正係数Ｋｔと同列に扱い、Ｋｃ、Ｋｔと共に、このスキッシュ係数Ｋｑを初期燃焼期間における平均乱れ強さＳＴ１、主燃焼期間における平均乱れ強さＳＴ２の算出のいずれにも導入している。また、このときのスキッシュ補正係数Ｋｑは適合値（一定値）で与えている。

さて、ここではエンジン仕様の具体例を挙げて説明する。図３９は１気筒当たり２つの吸気弁１５及び２つの排気弁１６を備えている４弁エンジンについて１気筒分のみを拡大した概略平面図である。こうした構成の燃焼室形状に対して、本実施形態では、スキッシュエリアを吸気弁１５、排気弁１６と干渉する部分を除いて４つの部位に分割する。これら４つのスキッシュエリアを区別するため図で真上より時計回りに番号をふり、第１スキッシュエリアＡ１、第２スキッシュエリアＡ２、第３スキッシュエリアＡ３、第４スキッシュエリアＡ４とする。この場合に、各スキッシュエリアＡ１〜Ａ４の体積は図面により定まるので、図面により定まるこれら第１スキッシュエリア体積、第２スキッシュエリア体積、第３スキッシュエリア体積、第４スキッシュエリア体積をそれぞれｓｑｖ１、ｓｑｖ２、ｓｑｖ３、ｓｑｖ４（いずれも一定値）とする。

こうした複数のスキッシュエリアＡ１〜Ａ４を有する燃焼室形状を備えるエンジンを対象として構成したのが図３８のフローである。

図３８（図３７のステップ２５４のサブルーチン）においてステップ２６１ではエンジン回転速度ＮＲＰＭ［ｒｐｍ］を読み込む。

ステップ２６２では変数ｉに１を入れる。この変数ｉはスキッキュエリア番号を指示するものである。すなわち、ｉ＝１のとき第１スキッシュエリアＡ１を、ｉ＝２のとき第２スキッシュエリアＡ２を、ｉ＝３のとき第３スキッシュエリアＡ３を、ｉ＝４のとき第４スキッシュエリアＡ４をそれぞれ指示する。

この変数ｉとスキッシュエリアの総数である４とをステップ２６３で比較する。このとき変数ｉには１が入っているので、ステップ２６４に進みｉ＝１であるときの部位別スキッシュエリア体積ｓｑｖ（ｉ）を算出する。例えば、図４０に示したように変数ｉと第１から第４までの部位別スキッシュエリア体積ａｑｖ１〜ｓｑｖ４とを対照できる表を予め作成しておき、この表から変数ｉの指示するスキッシュエリア体積を拾ってくればよい。このときにはｉ＝１であるので、部位別スキッシュエリア体積ｓｑｖ（１）＝第１スキッシュエリア体積ｓｑｖ１である。

ステップ２６５ではスキッシュエリアからのガス吹き出し方向角とエンジン回転速度ＮＲＰＭとから図４１を内容とするマップを検索することにより部位別ガス吹き出し方向補正係数ｓｑｔ（ｉ）を算出する。

ここで、スキッシュエリアからのガス吹き出し方向には図４２に示したように２種類あることが知られており（上段の（ａ）はガス吹き出し方向が水平でないもの、下段の（ｂ）はガス吹き出し方向が水平方向である）、エンジンの仕様により予めいずれかに定まっている。また、図４２上段の（ａ）の場合にガス吹き出し方向はエンジン回転速度ＮＲＰＭが基準回転速度より外れたときに変化する。そこで、エンジンの仕様（基準回転速度を含む）により予め定まるガス吹き出し方向角と、エンジン回転速度ＮＲＰＭとをパラメータとして図４１に示すような特性を予め実験等により求めておけば、エンジンの仕様により定まっているガス吹き出し方向角と、そのときのエンジン回転速度ＮＲＰＭとから図４１を内容とするマップを検索することにより部位別ガス吹き出し方向補正係数ｓｑｔ（ｉ）を求めることができる。

図４１に示したように部位別ガス吹き出し方向補正係数ｓｑｔ（ｉ）は回転速度が一定の条件においてガス吹き出し方向角が大きくなるほど大きくなってピークを採りその後は小さくなる特性である。これは、ピークを採る位置の吹き出し方向角のときちょうど点火プラグ６に向かってスキッシュエリアよりガスが吹き出しスキッシュが最大となるためである。

図４１には第１回転速度（例えば低回転速度）の場合、第２回転速度（例えば中回転速度）の場合、第３回転速度（例えば高回転速度）の場合の代表的な３つの特性を示しているが、さらに特性の数を増やしてもかまわない。また、図４１において実際の回転速度が３つの代表的な回転速度（第１、第２、第３の回転速度）の間にあるときには補間計算により部位別ガス吹き出し方向補正係数ｓｑｔ（ｉ）を求めればよい。

ステップ２６６では部位別スキッシュ補正係数ｋｓｑ（ｉ）ここでは第１スキッシュエリアＡ１についてのスキッシュ補正係数を次式により算出する。

ｋｓｑ（ｉ）＝ｓｑｖ（ｉ）×ｓｑｔ（ｉ）×ｋ２ …（６８）
ただし、ｋ２：係数、
（６８）式はスキッシュエリア体積に比例させて、かつ吹き出し方向角及び回転速度ＮＲＰＭにも考慮してスキッシュ補正係数ｋｓｑを算出するようにしたものである。（６８）式の係数ｋ２はスキッシュエリア体積を乱れ強さへと変換するための係数で、適用対象のエンジンで予め適合しその係数を記憶させておく。

ステップ２６７ではこの第１スキッシュエリアＡ１についての部位別スキッシュ補正係数ｋｓｑ（１）をメモリに保存したあと、ステップ２６８で変数ｉを１だけ増やし、ステップ２６３に戻る。このときｉ＝２であることよりステップ２６４〜２６６の操作を実行すると、第２スキッシュエリアＡ２についての部位別スキッシュ補正係数ｋｓｑ（２）が求まるので、これをステップ２６７でメモリに保存したあとステップ２６８で変数ｉを１だけ増やしステップ２６３に戻る。このときｉ＝３であることよりステップ２６４〜２６６の操作を実行すると、第３スキッシュエリアＡ３についての部位別スキッシュ補正係数ｋｓｑ（３）が求まるので、これをステップ２６７でメモリに保存したあと、ステップ２６８で変数ｉを１だけ増やしステップ２６３に戻る。このときｉ＝４であることよりステップ２６４〜２６６の操作を実行すると、第４スキッシュエリアＡ４についての部位別スキッシュ補正係数ｋｓｑ（４）が求まるので、これをステップ２６７でメモリに保存したあとステップ２６８で変数ｉを１だけ増やしステップ２６３に戻る。このときにはｉ＝５つまり全てのスキッシュエリアについてのスキッシュ補正係数の算出を終了しているので、ステップ２６３よりステップ２６９に進み、全てのスキッシュエリアについての補正係数を合算した値をスキッシュ補正係数Ｋｓｑとして、つまり次式により一気筒全体についてのスキッシュ補正係数Ｋｓｑを算出する。

Ｋｓｑ＝ｋｓｑ（１）＋ｋｓｑ（２）＋ｋｓｑ（３）＋ｋｓｑ（４）…（６９）
なお、ステップ２６３〜２６８のループ操作は一瞬にして終了するものである。

このようにしてスキッシュ補正係数Ｋｓｑの算出を終了したら図３７に戻り、ステップ２５５、２５６で吸気弁開時期ＩＶＯの状態での主燃焼期間における平均乱れ強さＵｏ２、吸気弁閉時期ＩＶＣの状態での主燃焼期間における平均乱れ強さＵｃ２をそれぞれ次式により算出する。

Ｕｏ２＝ｆ２４（Ｖｏ、ＮＲＰＭ）×Ｋｑ …（７０）
Ｕｃ２＝ｆ２４（Ｖｃ、ＮＲＰＭ）×Ｋｑ …（７１）
（７０）、（７１）式の主燃焼期間における平均乱れ強さ基本値ｆ２４（Ｖｏ、ＮＲＰＭ）、ｆ２４（Ｖｃ、ＮＲＰＭ）も吸気弁通過平均流速（ＶｏあるいはＶｃ）と回転速度ＮＲＰＭの関数であることを表している。この平均乱れ強さ基本値ｆ２４（Ｖｏ、ＮＲＰＭ）、ｆ２４（Ｖｃ、ＮＲＰＭ）は吸気弁通過平均流速Ｖｏ、Ｖｃや回転速度ＮＲＰＭの相違によって主燃焼期間における平均乱れ強さが変化するので、これを考慮するものである。主燃焼期間における平均乱れ強さ基本値ｆ２４（Ｖｏ、ＮＲＰＭ）、ｆ２４（Ｖｃ、ＮＲＰＭ）を、関数で与えるのではなくマップ検索により求めるようにしてもかまわない。図４３はこの場合のマップの特性であり、図示のように主燃焼期間における平均乱れ強さ基本値ｆ２４（Ｖｏ、ＮＲＰＭ）、ｆ２４（Ｖｃ、ＮＲＰＭ）は吸気弁通過平均流速Ｖｏ、Ｖｃが大きくなるほど、また回転速度ＮＲＰＭが高くなるほど大きくなる値である。

ステップ２５７では図２９のステップ２０６と同様にして加重平均係数ｋ１の算出し、ステップ２５８でこの加重平均係数ｋ１を用いて吸気弁開時期の状態での主燃焼期間における平均乱れ強さＵｏ２と、吸気弁閉時期の状態での主燃焼期間における平均乱れ強さＵｃ２とを加重平均した値に回転速度ＮＲＰＭを乗算した値を主燃焼期間における平均乱れ強さＳＴ２とする。つまり次式により主燃焼期間における平均乱れ強さＳＴ２を算出する。

ＳＴ２＝｛ｋ１×Ｕｏ２＋（１−ｋ１）×Ｕｃ２｝×ＮＲＰＭ …（７２）
（７２）式において回転速度ＮＲＰＭを露わにし、この露わにした回転速度ＮＲＰＭに比例させて主燃焼期間における平均乱れ強さＳＴ２を求めるようにしたのは実験結果によるものである。もちろん、（７２）式に代えて次式としてもかまわない。

ＳＴ２＝｛ｋ１×Ｕｏ２＋（１−ｋ１）×Ｕｃ２｝ …（７３）
次に、図４４は、上より吸気ポート圧、シリンダ吸入空気量、吸気弁通過平均流速、燃焼期間における平均乱れ強さ、燃焼期間における乱流火炎速度が低負荷状態からの加速時にどのように変化するのかを示す波形図である。ただし、簡単のためＶＥＬ機構５１は、加速中ずっと作動状態か非作動状態のいずれかにあって加速途中に作動状態から非作動状態へとあるいはこの逆へと変化しないものとする。

アクセルペダルを踏み込んでの急加速により特にＮｏ．２気筒の吸気行程において吸気ポート圧が小さい値から大きい値へと急激に変化している。このため、Ｎｏ．２気筒では吸気弁開時期ＩＶＯの状態と吸気弁閉時期ＩＶＣの状態とで吸気ポート圧が大きく異なることから、吸気弁開時期ＩＶＯの状態での吸入空気量Ｑｏと吸気弁閉時期ＩＶＣの状態での吸入空気量Ｑｃとが異なり、図示のように吸気弁開時期ＩＶＯの状態での吸入空気量Ｑｏ（第３段目の一点鎖線参照）よりも吸気弁閉時期ＩＶＣの状態での吸入空気量Ｑｃ（第３段目の実線参照）のほうが多くなる。

吸気弁通過平均流速、平均乱れ強さについても、図示のように吸気弁開時期ＩＶＯの状態での吸気弁通過平均流速Ｖｏ（第４段目の一点鎖線参照）、平均乱れ強さＵｏ１、Ｕｏ２（第５段目の一点鎖線参照）よりも吸気弁閉時期ＩＶＯの状態での吸気弁通過平均流速Ｖｃ（第４段目の実線参照）、平均乱れ強さＵｃ１、Ｕｃ２（第５段目の実線参照）のほうが大きくなる。

このため、吸気弁開時期ＩＶＯの状態での初期燃焼期間における平均乱れ強さＵｏ１を初期燃焼期間における層流火炎速度ＳＬ１に乗算して得られる初期燃焼期間における乱流火炎速度ＦＬＡＭＥ１ｏ（＝Ｕｏ１×ＳＴ１）（第６段目の二点鎖線参照）よりも吸気弁閉時期ＩＶＣの状態での初期燃焼期間における平均乱れ強さＵｃ１を初期燃焼期間における層流火炎速度ＳＬ１に乗算して得られる初期燃焼期間における乱流火炎速度ＦＬＡＭＥ１ｃ（＝Ｕｃ１×ＳＴ１）（第６段目の実線参照）のほうが大きく、かつ吸気弁開時期ＩＶＯの状態での主燃焼期間における平均乱れ強さＵｏ２を、主燃焼期間における層流火炎速度ＳＬ２に乗算して得られる主燃焼期間における乱流火炎速度ＦＬＡＭＥ２ｏ（＝Ｕｏ２×ＳＴ２）（最下段の二点鎖線参照）よりも吸気弁閉時期ＩＶＣの状態での初期燃焼期間における平均乱れ強さＵｃ２を主燃焼期間における層流火炎速度ＳＬ２に乗算して得られる主燃焼期間における乱流火炎速度ＦＬＡＭＥ２ｃ（＝Ｕｃ２×ＳＴ２）（最下段の実線参照）のほうが大きくなる。

この場合に、吸気弁開時期ＩＶＯでの状態と吸気弁閉時期ＩＶＣの状態とが異なっているＮｏ．２気筒について初期燃焼期間及び主燃焼期間を算出するための乱流火炎速度として吸気弁開時期ＩＶＯの状態での上記乱流火炎速度ＦＬＡＭＥ１ｏ及びＦＬＡＭＥ２ｏを用いたのでは、実際の初期燃焼期間、実際の主燃焼期間よりいずれも短くなるし、この反対に吸気弁閉時期ＩＶＣの状態での上記乱流火炎速度ＦＬＡＭＥ１ｃ及びＦＬＡＭＥ２ｃを用いたのでは、実際の初期燃焼期間、実際の主燃焼期間よりいずれも長くなってしまう。

これに対して本実施形態によれば、吸気弁開時期ＩＶＯの状態での初期燃焼期間における平均乱れ強さＵｏ１と、吸気弁閉時期ＩＶＣの状態での初期燃焼期間における平均乱れ強さＵｃ１とを加重平均した値を初期燃焼期間における平均乱れ強さＳＴ１として算出するので（図２９のステップ２０７参照）、この加重平均値としての初期燃焼期間における平均乱れ強さＳＴ１を初期燃焼期間における層流火炎速度ＳＬ１に乗算して得られる初期燃焼期間における乱流火炎速度ＦＬＡＭＥ１（＝ＳＬ１×ＳＴ１）は吸気弁開時期ＩＶＯの状態での上記乱流火炎速度ＦＬＡＭＥ１ｏと、吸気弁閉時期ＩＶＣの状態での上記乱流火炎速度ＦＬＡＭＥ１ｃの間の値を採り（第６段目の破線参照）、また吸気弁開時期ＩＶＯの状態での主燃焼期間における平均乱れ強さＵｏ２と、吸気弁閉時期ＩＶＣの状態での主燃焼期間における平均乱れ強さＵｃ２とを加重平均した値を主燃焼期間における平均乱れ強さＳＴ２として算出するので（図３７のステップ２５８参照）、この加重平均値としての主燃焼期間における平均乱れ強さＳＴ２を主燃焼期間における層流火炎速度ＳＬ２に乗算して得られる主燃焼期間における乱流火炎速度ＦＬＡＭＥ２（＝ＳＬ２×ＳＴ２）は吸気弁開時期ＩＶＯの状態での上記乱流火炎速度ＦＬＡＭＥ２ｏと、吸気弁閉時期ＩＶＣの状態での上記乱流火炎速度ＦＬＡＭＥ２ｃの間の値を採る（最下段の破線参照）。

これにより、Ｎｏ．２気筒のように、吸気行程中に吸気ポート圧が変化して、吸気弁開時期の状態と吸気弁閉時期の状態とが異なる気筒が存在する場合においても、燃焼期間（初期燃焼期間、主燃焼期間）を精度よく算出することができる。なお、減速時にも同様の作用効果が得られる。

図４５は燃焼室内圧力の変化を示している。

（Ａ）、（Ｂ）はそれぞれ低負荷状態、高負荷状態（いずれも定常時）のもので、高負荷状態のほうが低負荷状態より吸気ポート圧が高い分、基準クランク角θＰＭＡＸでの最大値が高くなっている。なお、いずれの状態でも基本点火時期ＭＢＴＣＡＬで点火するので、基準クランクθＰｍａｘの位置で燃焼室内圧力が最大となっている。

これら２つの定常状態では吸気弁開時期の状態と、吸気弁閉時期の状態とが変わらないので、基本点火時期ＭＢＴＣＡＬとしては、燃焼が急激に進行する高負荷状態での基本点火時期ＭＢＴＣＡＬのほうが燃焼が緩やかな低負荷状態での基本点火時期ＭＢＴＣＡＬより遅角側にくる。

一方、（Ｄ）は低負荷状態から高負荷状態へと過渡的に変化し、吸気行程中に吸気ポート圧が上昇した、つまり吸気弁開時期の状態での吸気ポート圧より吸気弁閉時期の状態での吸気ポート圧が高くなっている気筒（上記のＮｏ．２気筒）について、吸気弁開時期の状態での吸気ポート圧に基づいて算出した基本点火時期で点火した場合の燃焼室内圧力の変化で、ノックが発生している。これは、当該気筒では吸気弁開時期後に吸気ポート圧が上昇するのに伴って、実際の平均乱れ強さは、吸気弁開時期の状態での平均乱れ強さより上昇するため、実際の乱流火炎速度（燃焼速度）も吸気弁開時期の状態での平均乱れ強さより算出される乱流火炎速度より上昇する。従って、吸気弁開時期の状態での平均乱れ強さを用いて算出した基本点火時期ＭＢＴＣＡＬで点火した場合には、実際の要求点火時期より早過ぎることになり燃焼が急激に進行してノックを起こすことになるのである。

これに対して、（Ｃ）は、本実施形態により吸気弁開時期の状態での平均乱れ強さＵｏ１、Ｕｏ２と、吸気弁閉時期の状態での平均乱れ強さＵｃ１、Ｕｃ２とを加重平均して得られる値ＳＴ１、ＳＴ２に基づいて基本点火時期ＭＢＴＣＡＬを算出した場合の燃焼室内圧力の変化を示したものである。（Ｃ）によれば乱流火炎速度と燃焼室内圧力の最大値とは（Ａ）と（Ｂ）の間に位置し、この場合にも基準クランクθＰｍａｘで燃焼室内圧力が最大となるのであり、ノックは発生しない。

上記の図４５は上記（Ａ）〜（Ｄ）の場合について燃焼室内圧力の変化で示したのに対して、図４６は上記（Ａ）〜（Ｄ）の場合について燃焼質量割合の変化で示したものである。

図４６において（Ａ）、（Ｂ）はそれぞれ低負荷状態、高負荷状態（いずれも定常時）のもので、高負荷状態のほうが低負荷状態より吸気ポート圧が高い分燃焼状態が良好であるため、燃焼質量割合の立ち上がりが急激となっている。なお、いずれの定常状態でも基本点火時期ＭＢＴＣＡＬで点火するので、基準クランク角θＰｍａｘでの燃焼質量割合が一定値Ｒｍａｘ（＝６０％）となっている。

これら２つの定常状態では吸気弁開時期の状態と、吸気弁閉時期の状態とが変わらないので、基本点火時期ＭＢＴＣＡＬとしては、燃焼が急激に進行する高負荷状態での基本点火時期ＭＢＴＣＡＬのほうが燃焼が緩やかな低負荷状態での基本点火時期ＭＢＴＣＡＬより遅角側にくる。

一方、（Ｄ）は低負荷状態から高負荷状態へと過渡的に変化し、吸気行程中に吸気ポート圧が上昇した、つまり吸気弁開時期の状態での吸気ポート圧より吸気弁閉時期の状態での吸気ポート圧が高くなっている気筒（上記のＮｏ．２気筒）について、吸気弁開時期の状態での吸気ポート圧に基づいて算出した基本点火時期ＭＢＴＣＡＬで点火した場合の燃焼質量割合の変化で、ノックが発生している。これは、当該気筒では吸気弁開時期後に吸気ポート圧が上昇するのに伴って、実際の平均乱れ強さは、吸気弁開時期の状態での平均乱れ強さより上昇するため、実際の乱流火炎速度（燃焼速度）も吸気弁開時期の状態での平均乱れ強さより算出される乱流火炎速度より上昇する。従って、吸気弁開時期の状態での平均乱れ強さを用いて算出した基本点火時期ＭＢＴＣＡＬで点火した場合には、実際の要求点火時期より早過ぎることになり、このときの基準クランク角θＰｍａｘでの燃焼質量割合が上記の一定値Ｒｍａｘより大きくなり、ノックが発生する。

これに対して、（Ｃ）は、本実施形態により吸気弁開時期の状態での平均乱れ強さＵｏ１、Ｕｏ２と、吸気弁閉時期の状態での平均乱れ強さＵｃ１、Ｕｃ２とを加重平均した値に基づいて基本点火時期ＭＢＴＣＡＬを算出した場合の燃焼質量割合の変化を示したものである。（Ｃ）によれば基準クランク角θＰｍａｘでの燃焼質量割合が上記の一定値Ｒｍａｘへと一致することになっておりノックは発生しない。

ここで、本実施形態の作用効果を説明する。

スキッシュの形状によりスキッシュが主燃焼期間における乱流火炎速度に影響する程度が定まる、という実験結果を新たに得ている。ということは、スキッシュの形状はエンジンの仕様により定まるのであるから、エンジンの仕様によりスキッシュが主燃焼期間における乱流火炎速度に影響する程度をほぼ推定できることを意味する。すなわち、本実施形態（請求項１、２に記載の発明）によれば主燃焼期間ＢＵＲＮ２がスキッシュの生じる期間を含み、スキッシュ補正係数Ｋｑ（スキッシュによる火炎速度補正量）をスキッシュエリア体積やスキッシュエリアからのガス吹き出し方向といったエンジンの仕様に少なくとも基づいて算出し（図３８のステップ２６４、２６５参照）、主燃焼期間における燃焼乱流火炎速度ＦＬＡＭＥ２を主燃焼期間における層流火炎速度ＳＬ２と、このスキッシュ補正係数Ｋｑとに基づいて算出するので（図１２のステップ１８９、１９０、図３７のステップ２５５、２５６、２５８参照）、エンジン機種が異なる毎に適合し直す必要が無くなり、適合工数を低減できる。

実験結果によれば、初期燃焼期間にはスキッシュの生じる期間を含まず、従ってスキッシュは初期燃焼間における乱流火炎速度ＦＬＡＭＥ１に影響しないことに対応して、本実施形態（請求項３に記載の発明）によれば、初期燃焼期間における乱流火炎速度ＦＬＡＭＥ１をスキッシュ補正係数Ｋｑ（スキッシュによる火炎速度補正量）に基づいては算出しないので、無駄な算出をしないで済む。

スキッシュエリア体積が大きいほど主燃焼期間における平均乱れ強さが強くなり主燃焼期間における乱流火炎速度ＦＬＡＭＥ２が速くなることに対応して、本実施形態（請求項４に記載の発明）によれば、スキッシュ補正係数Ｋｑ（スキッシュによる火炎速度補正量）をエンジンの仕様であるスキッシュエリア体積（ｓｑｖ１〜ｓｑｖ４）に基づいて算出するので（図３８のステップ２６４参照）、スキッシュエリア体積の大小に関係なく主燃焼期間における乱流火炎速度ＦＬＡＭＥ２を正確に算出できる。

スキッシュエリアからのガス吹き出し方向も主燃焼期間における平均乱れ強さを変化させ主燃焼期間における乱流火炎速度ＦＬＡＭＥ２が変化することに対応して、本実施形態（請求項５に記載の発明）によれば、スキッシュ補正係数Ｋｑ（スキッシュによる火炎速度補正量）をエンジンの仕様であるスキッシュエリアからのガス吹き出し方向に基づいて算出するので（図３８のステップ２６５参照）、スキッシュエリアからのガス吹き出し方向に関係なく主燃焼期間における乱流火炎速度ＦＬＡＭＥ２を正確に算出できる。

スキッシュエリアの形状は一気筒の中でも部位毎に異なることに対応して、本実施形態（請求項６に記載の発明）によれば、一気筒の中でスキッシュエリアを複数の部位に分割し、その分割した部位別スキッシュエリア毎に部位別スキッシュ補正係数ｋｓｑ（ｉ）を算出し、これらの部位別スキッシュ補正係数ｋｓｑ（ｉ）を合計した値を一気筒全体についてのスキッシュ補係数Ｋｑ（スキッシュによる火炎速度補正量）とするので（図３８のステップ２６９参照）、一気筒の中で部位毎にスキッシュエリアの形状が異なる場合においても、主燃焼期間における乱流火炎速度ＦＬＡＭＥ２を正確に算出できる。

本実施形態（請求項７に記載の発明）によれば、吸気弁開時期ＩＶＯから吸気弁閉時期ＩＶＣまでの間にＶＥＬ機構５１（可変バルブ機構）が作動から非作動へと変化する場合に、吸気弁開時期ＩＶＯの状態での吸気弁通過平均流速Ｖｏとスキッシュ補正係数Ｋｑ（スキッシュによる火炎速度補正量）とに基づいて吸気弁開時期ＩＶＯの状態での主燃焼期間における平均乱れ強さＵｏ２を算出する手段（図３７のステップ２５５参照）と、吸気弁閉時期ＩＶＣの状態での吸気弁通過平均流速Ｖｃとスキッシュ補正係数Ｋｑ（スキッシュによる火炎速度補正量）とに基づいて吸気弁閉時期ＩＶＣの状態での主燃焼期間における平均乱れ強さＵｃ２を算出する手段（図３７のステップ２５６参照）と、これらの加重平均値を主燃焼期間における平均乱れ強さＳＴ２として算出する手段（図３７のステップ２５８参照）と、この主燃焼期間における平均乱れ強さＳＴ２と、主燃焼期間における層流火炎速度ＳＬ２とから主燃焼期間における乱流火炎速度ＦＬＡＭＥ２を算出する手段（図１２のステップ１９０参照）とを含んでいるので、吸気弁開時期ＩＶＯから吸気弁閉時期ＩＶＣまでの間にＶＥＬ機構５１が作動から非作動へと変化する場合においても、主燃焼期間における平均乱れ強さＳＴ２を精度よく算出できる。

主燃焼期間における平均乱れ強さＳＴ２がエンジン回転速度ＮＲＰＭに比例するという実験結果を得ていることに対応して、本実施形態（請求項８に記載の発明）によれば、吸気弁開時期ＩＶＯの状態での主燃焼期間における平均乱れ強さＵｏ２と、吸気弁閉時期ＩＶＣの状態での主燃焼期間における平均乱れ強さＵｃ２との加重平均値にエンジン回転速度ＮＲＰＭをさらに乗算した値を主燃焼期間における平均乱れ強さＳＴ２として算出するので（図３７のステップ２５８参照）、主燃焼期間における平均乱れ強さＳＴ２を簡易に算出できる。

本実施形態（請求項９に記載の発明）によれば、吸気弁開時期の状態での吸気弁通過平均流速Ｖｏとスキッシュ補正係数Ｋｑ（スキッシュによる燃焼速度補正量）とに基づいて吸気弁開時期の状態での主燃焼期間における平均乱れ強さＵｏ２を算出する手段（図３７のステップ２５５参照）と、吸気弁閉時期の状態での吸気弁通過平均流速Ｖｃとスキッシュ補正係数Ｋｑ（スキッシュによる燃焼速度補正量）とに基づいて吸気弁閉時期の状態での主燃焼期間における平均乱れ強さＵｃ２を算出する手段（図３７のステップ２５６参照）と、これらの加重平均値を主燃焼期間における平均乱れ強さＳＴ２として算出する手段（図３７のステップ２５８参照）と、この主燃焼期間における平均乱れ強さＳＴ２と、前記主燃焼期間における層流火炎速度ＳＬ２とから前記主燃焼期間における乱流火炎速度ＦＬＡＭＥ２を算出する手段（図１２のステップ１９０参照）とを含んでいるので、ＶＥＬ機構５１が作動状態や非作動状態のいずれかを保っている場合において吸気弁開時期ＩＶＯから吸気弁閉時期ＩＶＣまでの間に吸気ポート圧が大きく変化する加速時や減速時においても、主燃焼期間における平均乱れ強さＳＴ２を精度よく算出できる。

実施形態では、初期燃焼期間を燃焼質量割合の変化代としてゼロから２％まで（つまりＢＲ１＝２％）、主燃焼期間を燃焼質量割合の変化代として２〜６０％まで（つまりＢＲ２＝５８％）と規定したが、本発明は必ずしもこの数値に限定されるものでない。

実施形態では、燃焼ガス質量割合で説明したが、燃焼ガス質量そのものを用いてもかまわない。

実施形態では燃焼期間をスキッシュの生じない初期燃焼期間と、スキッシュの生じる主燃焼期間との２つに分割する場合で説明したが、３以上に分割する場合でもかまわない（請求項１に記載の発明）。例えば、燃焼開始より所定クランク角までに燃焼室内で燃焼する燃焼ガス質量の変化代を３以上の複数に分割し、その分割されたそれぞれの燃焼ガス質量の変化代に対応する分割燃焼期間を算出し、それら全ての分割燃焼期間を合計した値を燃焼開始から所定クランク角までの燃焼期間として算出してもかまわない（請求項１、２に記載の発明）。

実施形態では、吸気弁開時期から吸気弁閉時期までの間にＶＥＬ機構５１（可変バルブ機構）が作動から非作動へと変化する場合で説明したが、吸気弁開時期から吸気弁閉時期までの間にＶＥＬ機構５１が非作動から作動へと変化する場合にも適用することができる（請求項７に記載の発明）。

実施形態では、ＶＥＬ機構５１を備える場合で説明したが、ＶＥＬ機構５１を備えない場合にも適用があることはいうまでもない。

請求項１に記載の発明において、燃焼速度算出手段の機能は図１０のステップ１７０、図１２のステップ１９０により、基本点火時期算出手段の機能は図１３のステップ４３により、火炎速度補正量算出手段の機能は図３８のステップ２６４、２６５により、乱流火炎速度算出手段の機能は図１２のステップ１８９、１９０、図３７のステップ２５４〜２５８によりそれぞれ果たされている。

請求項２に記載の発明において、燃焼速度算出手段の機能は図１０のステップ１７０、図１２のステップ１９０により、基本点火時期算出手段の機能は図１３のステップ４３により、火炎速度補正量算出手段の機能は図３８のステップ２６４、２６５により、乱流火炎速度算出手段の機能は図１２のステップ１８９、１９０、図３７のステップ２５４〜２５８によりそれぞれ果たされている。

一実施形態のエンジンの制御システム図。 エンジンコントローラで実行される点火時期制御のブロック図。 燃焼室の圧力変化図。 燃焼質量割合の変化を説明する特性図。 物理量の算出を説明するためのフローチャート。 エンジンのクランクシャフトとコネクティングロッドの位置関係を説明するダイアグラム。 水温補正係数の特性図。 当量比補正係数の特性図。 基準クランク角の特性図。 初期燃焼期間の算出を説明するためのフローチャート。 温度上昇率の特性図。 主燃焼期間の算出を説明するためのフローチャート。 基本点火時期の算出を説明するためのフローチャート。 内部不活性ガス率の算出を説明するためのフローチャート。 内部不活性ガス量の算出を説明するためのフローチャート。 ＥＶＣ時不活性ガス量の算出を説明するためのフローチャート。 オーバーラップ中吹き返し不活性ガス量の算出を説明するためのフローチャート。 過給判定フラグ、チョーク判定フラグの設定を説明するためのフローチャート。 過給無しかつチョーク無し時のオーバーラップ中吹き返し不活性ガス流量の算出を説明するためのフローチャート。 過給無しかつチョーク有り時のオーバーラップ中吹き返し不活性ガス流量の算出を説明するためのフローチャート。 過給有りかつチョーク無し時のオーバーラップ中吹き返し不活性ガス流量の算出を説明するためのフローチャート。 過給有りかつチョーク有り時のオーバーラップ中吹き返し不活性ガス流量の算出を説明するためのフローチャート。 排気弁閉時期における燃焼室容積の特性図。 不活性ガスのガス定数の特性図。 オーバーラップ中の積算有効面積の特性図。 オーバーラップ中の積算有効面積の説明図。 不活性ガスの比熱比の特性図。 混合気の比熱比の特性図。 初期燃焼期間における平均乱れ強さの算出を説明するためのフローチャート。 吸気弁開時期の状態での吸気弁通過平均流速の算出を説明するためのフローチャート。 吸気行程中にＶＥＬ機構が非作動へと変化した場合の、吸気弁開時期の状態での開弁クランク角区間θａ、吸気弁開時期の状態での吸気弁リフト量Ｌａ、吸気弁閉時期の状態での開弁クランク角区間θｂ、吸気弁閉時期の状態での吸気弁リフト量Ｌｂを説明するための波形図。 補正係数の特性図。 平均吸気弁開口面積の特性図。 吸気弁閉時期の状態での吸気弁通過平均流速の算出を説明するためのフローチャート。 初期燃焼期間における平均乱れ強さ基本値の特性図。 加重平均係数の算出を説明するためのフローチャート。 主燃焼期間における平均乱れ強さの算出を説明するためのフローチャート。 スキッシュ補正係数の算出を説明するためのフローチャート。 一気筒中のスキッシュエリアの分割を説明するためのエンジンの概略平面図。 部位別スキッシュエリア体積の表図。 部位別ガス吹き出し方向補正係数の特性図。 スキッシュエリアからのガス吹出し方向を説明するための概略断面図。 主燃焼期間における平均乱れ強さ基本値の特性図。 加速時の代表的な値の変化を示す波形図。 燃焼室内圧力の変化を示す波形図。 燃焼質量割合の変化を示す波形図。

符号の説明

１ エンジン
５ 燃焼室
１１ 点火装置（火花点火手段）
１５ 吸気弁
３１ エンジンコントローラ
５１ ＶＥＬ機構（可変バルブ機構）

Claims (9)

1. 燃焼室内で混合気が燃焼する燃焼期間を複数に分割し、その分割した各燃焼期間における燃焼速度を運転条件に基づいて算出する燃焼速度算出手段と、
   これら複数の燃焼速度に基づいてＭＢＴの得られる点火時期を基本点火時期として算出する基本点火時期算出手段と、
   この基本点火時期で火花点火を行う火花点火手段と
   を備えるエンジンの点火時期制御装置において、
   前記いずれか一つの燃焼期間がスキッシュの生じる期間を含み、
   このスキッシュの生じる期間を含む燃焼期間における燃焼速度がこのスキッシュの生じる期間を含む燃焼期間における乱流火炎速度であり、
   スキッシュによる火炎速度補正量を少なくともエンジンの仕様に基づいて算出する火炎速度補正量算出手段と、
   前記スキッシュの生じる期間を含む燃焼期間における乱流火炎速度を、このスキッシュの生じる期間を含む燃焼期間における層流火炎速度とエンジンの運転条件のうち少なくともいずれか一つと、前記スキッシュによる火炎速度補正量とに基づいて算出する乱流火炎速度算出手段と
   を備えることを特徴とするエンジンの点火時期制御装置。
2. 燃焼室内で混合気が燃焼する燃焼期間を初期燃焼期間と主燃焼期間の複数に分割し、その分割した各燃焼期間における燃焼速度を運転条件に基づいて算出する燃焼速度算出手段と、
   これら複数の燃焼速度に基づいてＭＢＴの得られる点火時期を基本点火時期として算出する基本点火時期算出手段と、
   この基本点火時期で火花点火を行う火花点火手段と
   を備えるエンジンの点火時期制御装置において、
   前記主燃焼期間がスキッシュの生じる期間を含み、
   前記主燃焼期間における燃焼速度が前記主燃焼期間における乱流火炎速度であり、
   スキッシュによる火炎速度補正量を少なくともエンジンの仕様に基づいて算出する火炎速度補正量算出手段と、
   前記主燃焼期間における乱流火炎速度を、前記主燃焼期間における層流火炎速度とエンジンの運転条件のうち少なくともいずれか一つと、前記スキッシュによる火炎速度補正量とに基づいて算出する乱流火炎速度算出手段と
   を備えることを特徴とするエンジンの点火時期制御装置。
3. 前記初期燃焼期間における乱流火炎速度を前記スキッシュによる火炎速度補正量に基づいて算出しないことを特徴とする請求項２に記載のエンジンの点火時期制御装置。
4. 前記スキッシュによる火炎速度補正量を前記エンジンの仕様であるスキッシュエリア体積に基づいて算出することを特徴とする請求項１または２に記載のエンジンの点火時期制御装置。
5. 前記スキッシュによる乱流火炎速度補正量を前記エンジンの仕様であるスキッシュエリアからのガス吹き出し方向に基づいて算出することを特徴とする請求項１または２に記載のエンジンの点火時期制御装置。
6. 一気筒の中で前記スキッシュエリアを複数の部位に分割し、その分割した部位別スキッシュエリア毎にスキッシュによる燃焼速度補正量を算出し、これらの部位別スキッシュエリア毎のスキッシュによる火炎速度補正量を合計した値を一気筒全体についてのスキッシュによる火炎速度補正量とすることを特徴とする請求項４または５に記載のエンジンの点火時期制御装置。
7. 吸気弁の弁リフト量及び作動角を連続的に可変制御する可変バルブ機構を備え、
   吸気弁開時期から吸気弁閉時期までの間にこの可変バルブ機構が作動から非作動へとまたはこの逆へと変化する場合に、
   前記吸気弁開時期の状態での吸気弁通過平均流速と前記スキッシュによる燃焼速度補正量とに基づいて吸気弁開時期の状態での主燃焼期間における平均乱れ強さを算出する手段と、
   前記吸気弁閉時期の状態での吸気弁通過平均流速と前記スキッシュによる燃焼速度補正量とに基づいて吸気弁閉時期の状態での主燃焼期間における平均乱れ強さを算出する手段と、
   これらの加重平均値を主燃焼期間における平均乱れ強さとして算出する手段と、
   この主燃焼期間における平均乱れ強さと、前記主燃焼期間における層流火炎速度とから前記主燃焼期間における乱流火炎速度を算出する手段と
   を含んでいることを特徴とする請求項２に記載のエンジンの点火時期制御装置。
8. 前記加重平均値にエンジン回転速度をさらに乗算した値を前記主燃焼期間における平均乱れ強さとして算出することを特徴とする請求項７に記載のエンジンの点火時期制御装置。
9. 前記吸気弁開時期の状態での吸気弁通過平均流速と前記スキッシュによる燃焼速度補正量とに基づいて吸気弁開時期の状態での主燃焼期間における平均乱れ強さを算出する手段と、
   前記吸気弁閉時期の状態での吸気弁通過平均流速と前記スキッシュによる燃焼速度補正量とに基づいて吸気弁閉時期の状態での主燃焼期間における平均乱れ強さを算出する手段と、
   これらの加重平均値を主燃焼期間における平均乱れ強さとして算出する手段と、
   この主燃焼期間における平均乱れ強さと、前記主燃焼期間における層流火炎速度とから前記主燃焼期間における乱流火炎速度を算出する手段と
   を含んでいることを特徴とする請求項２に記載のエンジンの点火時期制御装置。

Images