(57)【要約】【課題】この発明は、内燃機関の気筒内の熱発生を精度良くモデル化することのできるWiebe関数パラメータの決定方法を提供することを目的とする。また、この発明は、内燃機関の気筒内の熱発生率を精度良く推定することのできる内燃機関の熱発生率推定装置を提供することを他の目的とする。【解決手段】クランク角度毎の実熱発生率を求め、実熱発生率が最大となった時点でのクランク角度である熱発生率最大クランク角度ＣＡｍをWiebe関数パラメータの一つとして決定する。実熱発生率が最大となった時点での実燃焼割合αに基づいて、形状パラメータｍの値を決定する。熱発生率最大クランク角度ＣＡｍ、形状パラメータｍ等に基づいて、Wiebe関数の熱発生開始点のクランク角度をＣＡｓを定める。【選択図】図２
2007248119.tif 000002
【特許請求の範囲】
【請求項１】
内燃機関の気筒内の燃焼による熱発生をモデル化するためのWiebe関数のパラメータの値を決定する方法であって、
前記内燃機関の気筒内におけるクランク角度毎の実熱発生率を求めるステップと、
前記実熱発生率が最大となった時点でのクランク角度である熱発生率最大クランク角度をWiebe関数パラメータの一つとして決定するステップと、
前記実熱発生率が最大となった時点での実燃焼割合を求めるステップと、
前記実燃焼割合に基づいて、Wiebe関数パラメータの一つである形状パラメータｍの値を決定するステップと、
を備え、
Wiebe関数の熱発生開始点のクランク角度を、前記熱発生率最大クランク角度の値と、前記形状パラメータｍの値と、Wiebe関数パラメータの一つである燃焼期間の値と、Wiebe関数パラメータの一つであるパラメータａの値と、に基づいて設定することを特徴とするWiebe関数パラメータの決定方法。
【請求項２】
前記形状パラメータｍの値を決定するステップでは、前記実燃焼割合の値をαとしたとき、下記（Ｉ）式により形状パラメータｍの値が決定され、
前記熱発生率最大クランク角度をＣＡｍ［°ATDC］、Wiebe関数の熱発生開始点のクランク角度をＣＡｓ［°BTDC］、前記燃焼期間をθｐ［°］、βを下記（II）式で表す値、としたとき、下記（III）式が成り立つことを特徴とする請求項１記載のWiebe関数パラメータの決定方法。
【数１】
2007248119.tif 000011
ＣＡｓ＝β×θｐ−ＣＡｍ      ・・・（III）
【請求項３】
前記決定された熱発生率最大クランク角度および形状パラメータｍの値を用いてWiebe関数モデルにより算出されるクランク角度毎のモデル熱発生率と、前記実熱発生率との誤差の大きさを表す誤差評価値が最小となるように、残りのWiebe関数パラメータの値を決定するステップを更に備えることを特徴とする請求項１または２記載のWiebe関数パラメータの決定方法。
【請求項４】
前記残りのWiebe関数パラメータは、前記燃焼期間θｐと、気筒内に供給された燃料の熱量が気筒内での熱発生量に変換される効率ｋとであることを特徴とする請求項３記載のWiebe関数パラメータの決定方法。
【請求項５】
内燃機関の運転条件と、気筒内の実熱発生率が最大となるクランク角度である熱発生率最大クランク角度との関係を記憶した関係記憶手段と、
Wiebe関数モデルに基づいて、気筒内の熱発生率の推定値であるモデル熱発生率をクランク角度毎に算出する熱発生率推定手段と、
を備え、
前記熱発生率推定手段は、前記関係記憶手段に記憶された関係に従って運転条件に応じて取得される熱発生率最大クランク角度に、前記モデル熱発生率が最大となるクランク角度が一致するように、前記モデル熱発生率を算出することを特徴とする内燃機関の熱発生率推定装置。
【請求項６】
前記関係記憶手段は、前記内燃機関の運転条件と、Wiebe関数パラメータである形状パラメータおよび燃焼期間との関係を更に記憶しており、
前記熱発生率推定手段は、前記関係記憶手段に記憶された関係に従って運転条件に応じて取得される熱発生率最大クランク角度をＣＡｍ［°ATDC］、形状パラメータをｍ、燃焼期間をθｐ［°］とし、Wiebe関数パラメータの一つであるａを所定の定数とし、βを下記（IV）式で表す値としたとき、Wiebe関数の熱発生開始点のクランク角度ＣＡｓ［°BTDC］が下記（Ｖ）式で表す値となるように、前記モデル熱発生率を算出することを特徴とする請求項５記載の内燃機関の熱発生率推定装置。
【数２】
2007248119.tif 000012
ＣＡｓ＝β×θｐ−ＣＡｍ      ・・・（Ｖ）
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、Wiebe関数パラメータの決定方法および内燃機関の熱発生率推定装置に係り、特に、内燃機関の気筒内の熱発生をモデル化するためのWiebe関数パラメータの決定方法およびWiebe関数モデルを利用した内燃機関の熱発生率推定装置に関する。
【背景技術】
【０００２】
従来より、内燃機関の気筒内の熱発生をWiebe関数によってモデル化する試みがなされている。Wiebe関数によれば、クランク角度毎の熱発生率や燃焼割合等を算出することができる。そして、クランク角度毎の熱発生率の値が分かれば、筒内圧や図示トルク等を算出することができる。ゆえに、気筒内の熱発生をWiebe関数によって正確にモデル化することができれば、内燃機関の開発段階や、内燃機関の制御装置などにおいて、Wiebe関数を有効に利用することが可能となる。
【０００３】
気筒内の熱発生をWiebe関数によって正確にモデル化するためには、Wiebe関数に含まれるいくつかのパラメータの値を、対象とする内燃機関の機種に応じて、精度良く決定する必要がある。また、それらのパラメータの値は、空燃比、機関回転数、負荷率、点火時期などの運転条件によっても異なる。
【０００４】
特開２００４−３３２６５９号公報には、Wiebe関数を利用して点火時期をＭＢＴ（Minimum advance for the Best Torque）付近に制御する技術が開示されている。この公報では、同一エンジンのＭＢＴにおいては、燃焼期間の長さにかかわらず熱発生率の形状（波形）が同じになり、かつ熱発生率が最大となるときのクランク角度および燃焼割合も一定になることを仮定している。そして、Wiebe関数パラメータのうち、熱発生率形状に関与する形状パラメータｍを定数とし、それらの値を実験で求めることとしている（上記公報の段落番号００３０参照）。なお、上記公報でのパラメータ「ｎ」が本明細書における「ｍ」に相当する。
【０００５】
【特許文献１】特開２００４−３３２６５９号公報
【発明の開示】
【発明が解決しようとする課題】
【０００６】
しかしながら、上記公報には、形状パラメータｍその他のWiebe関数パラメータを実験データに基づいて決定する具体的な方法については何ら開示されていない。また、本発明者の知見によれば、同一エンジンのＭＢＴにおいても、形状パラメータｍの値は運転条件によって変化し、一定とはならない。また、種々の制御が行われる現代のエンジンにおいては、必ずしも点火時期をＭＢＴに合わせればよいとは限らない。このため、ＭＢＴ以外の点火時期の下でのトルク等を推定する要請も存在し、そのためには点火時期の変化によるWiebe関数パラメータの値の変化傾向を把握する必要がある。
【０００７】
このように、内燃機関の開発や内燃機関の制御装置においてWiebe関数を有効に利用するには、様々な運転条件の下でWiebe関数パラメータを精度良く決定し、それらの値が運転条件に応じてどのように変化するかを正確に把握することが必要である。しかしながら、従来、Wiebe関数パラメータを精度良く決定する手法は確立されていないのが実状である。
【０００８】
この発明は、上述のような課題を解決するためになされたもので、内燃機関の気筒内の熱発生を精度良くモデル化することのできるWiebe関数パラメータの決定方法を提供することを目的とする。また、この発明は、内燃機関の気筒内の熱発生率を精度良く推定することのできる内燃機関の熱発生率推定装置を提供することを他の目的とする。
【課題を解決するための手段】
【０００９】
第１の発明は、上記の目的を達成するため、Wiebe関数パラメータの決定方法であって、
前記内燃機関の気筒内におけるクランク角度毎の実熱発生率を求めるステップと、
前記実熱発生率が最大となった時点でのクランク角度である熱発生率最大クランク角度をWiebe関数パラメータの一つとして決定するステップと、
前記実熱発生率が最大となった時点での実燃焼割合を求めるステップと、
前記実燃焼割合に基づいて、Wiebe関数パラメータの一つである形状パラメータｍの値を決定するステップと、
を備え、
Wiebe関数の熱発生開始点のクランク角度を、前記熱発生率最大クランク角度の値と、前記形状パラメータｍの値と、Wiebe関数パラメータの一つである燃焼期間の値と、Wiebe関数パラメータの一つであるパラメータａの値と、に基づいて設定することを特徴とする。
【００１０】
また、第２の発明は、第１の発明において、
前記形状パラメータｍの値を決定するステップでは、前記実燃焼割合の値をαとしたとき、下記（Ｉ）式により形状パラメータｍの値が決定され、
前記熱発生率最大クランク角度をＣＡｍ［°ATDC］、Wiebe関数の熱発生開始点のクランク角度をＣＡｓ［°BTDC］、前記燃焼期間をθｐ［°］、βを下記（II）式で表す値、としたとき、下記（III）式が成り立つことを特徴とする。
【数１】
2007248119.tif 000003
ＣＡｓ＝β×θｐ−ＣＡｍ      ・・・（III）
【００１１】
また、第３の発明は、第１または第２の発明において、
前記決定された熱発生率最大クランク角度および形状パラメータｍの値を用いてWiebe関数モデルにより算出されるクランク角度毎のモデル熱発生率と、前記実熱発生率との誤差の大きさを表す誤差評価値が最小となるように、残りのWiebe関数パラメータの値を決定するステップを更に備えることを特徴とする。
【００１２】
また、第４の発明は、第３の発明において、
前記残りのWiebe関数パラメータは、前記燃焼期間θｐと、気筒内に供給された燃料の熱量が気筒内での熱発生量に変換される効率ｋとであることを特徴とする。
【００１３】
また、第５の発明は、内燃機関の熱発生率推定装置であって、
内燃機関の運転条件と、気筒内の実熱発生率が最大となるクランク角度である熱発生率最大クランク角度との関係を記憶した関係記憶手段と、
Wiebe関数モデルに基づいて、気筒内の熱発生率の推定値であるモデル熱発生率をクランク角度毎に算出する熱発生率推定手段と、
を備え、
前記熱発生率推定手段は、前記関係記憶手段に記憶された関係に従って運転条件に応じて取得される熱発生率最大クランク角度に、前記モデル熱発生率が最大となるクランク角度が一致するように、前記モデル熱発生率を算出することを特徴とする。
【００１４】
また、第６の発明は、第５の発明において、
前記関係記憶手段は、前記内燃機関の運転条件と、Wiebe関数パラメータである形状パラメータおよび燃焼期間との関係を更に記憶しており、
前記熱発生率推定手段は、前記関係記憶手段に記憶された関係に従って運転条件に応じて取得される熱発生率最大クランク角度をＣＡｍ［°ATDC］、形状パラメータをｍ、燃焼期間をθｐ［°］とし、Wiebe関数パラメータの一つであるａを所定の定数とし、βを下記（IV）式で表す値としたとき、Wiebe関数の熱発生開始点のクランク角度ＣＡｓ［°BTDC］が下記（Ｖ）式で表す値となるように、前記モデル熱発生率を算出することを特徴とする。
【数２】
2007248119.tif 000004
ＣＡｓ＝β×θｐ−ＣＡｍ      ・・・（Ｖ）
【発明の効果】
【００１５】
第１の発明によれば、内燃機関の気筒内の熱発生をWiebe関数でモデル化する上で、気筒内の実熱発生率が最大となった時点でのクランク角度である熱発生率最大クランク角度をWiebe関数パラメータの一つとして決定することができる。熱発生率最大クランク角度によれば、Wiebe関数の熱発生開始点（燃焼開始点）の絶対的なクランク角度を定めることができる。そして、熱発生率最大クランク角度には、点火時期等の運転条件の影響が精度良く反映される。更に、熱発生率最大クランク角度は、実熱発生率をクランク角毎に取得することにより一義的に定まるので、その値を正確かつ容易に決定することができる。また、第１の発明によれば、実熱発生率が最大となった時点での実燃焼割合に基づいて形状パラメータｍの値を決定することにより、形状パラメータｍを高精度に決定することができる。このようなことから、第１の発明によれば、あらゆる運転条件において気筒内の熱発生を精度良く推定することのできる、高精度なWiebe関数モデルを構築することが可能となる。
【００１６】
第２の発明によれば、上記（Ｉ）式を用いることにより、形状パラメータｍの値を高精度かつ容易に決定することができる。また、上記（II）式および（III）式を用いることにより、Wiebe関数の熱発生開始点の絶対的なクランク角度ＣＡｓを高精度かつ容易に設定することができる。
【００１７】
第３の発明によれば、決定済みの熱発生率最大クランク角度および形状パラメータｍの値を用いてWiebe関数モデルにより算出されるモデル熱発生率と、実熱発生率との誤差が最小となるように、残りのWiebe関数パラメータの値を精度良く決定することができる。
【００１８】
第４の発明によれば、残りのWiebe関数パラメータとして、燃焼期間θｐおよび効率ｋの値を精度良く決定することができる。
【００１９】
第５の発明によれば、Wiebe関数モデルを用いて内燃機関の気筒内の熱発生率を推定する上で、気筒内の実熱発生率が最大となった時点でのクランク角度である熱発生率最大クランク角度をWiebe関数パラメータの一つとして使用することができる。すなわち、予め記憶された関係に基づいて運転条件に応じて取得される熱発生率最大クランク角度に、熱発生率の推定値が最大となるクランク角度が一致するように、熱発生率を推定することができる。熱発生率最大クランク角度には、点火時期等の運転条件の影響が精度良く反映される。このため、第５の発明によれば、熱発生率最大クランク角度をWiebe関数パラメータとして使用することにより、種々の運転条件の下での熱発生率を精度良く推定することができる。
【００２０】
第６の発明によれば、クランク角度毎の熱発生率を推定する上で、熱発生開始点のクランク角度を、熱発生率最大クランク角度に基づいて精度良く算出することができる。
【発明を実施するための最良の形態】
【００２１】
実施の形態１．
［システム構成の説明］
図１は、本発明の実施の形態１で用いられるシステム構成を説明するための図である。図１に示すように、本実施形態のシステムは、火花点火式の内燃機関１０を備えている。内燃機関１０の各気筒には、吸気弁１４、排気弁１６および点火プラグ３０が設けられている。
【００２２】
また、内燃機関１０には、クランク軸の回転位置（クランク角度）を検出するクランク角センサ１２が設けられている。クランク角センサ１２は、クランク軸が所定回転角だけ回転する毎に、Hi出力とLo出力を反転させるセンサである。クランク角センサ１２の出力によれば、機関回転数NEを検知することもできる。
【００２３】
また、内燃機関１０には、筒内圧センサ１８が組み込まれている。筒内圧センサ１８によれば、気筒内（燃焼室内）に生ずる圧力を検出することができる。
【００２４】
内燃機関１０の吸気通路１９の途中には、サージタンク２０が設けられている。サージタンク２０には、その内部の圧力、すなわち吸気管圧力を検出する吸気圧センサ２１が設置されている。また、吸気通路１９には、その内部を流通する吸入空気量GAを検出するエアフローメータ２２が配置されている。吸気圧センサ２１の出力によれば、内燃機関１０の負荷率KL［％］などを取得することができる。なお、負荷率KLは、吸入空気量GAと機関回転数NEとに基づいて算出してもよい。
【００２５】
エアフローメータ２２の下流には、スロットル弁２４が配置されている。スロットル弁２４の近傍には、スロットル開度TAを検出するスロットル開度センサ２６が組み付けられている。
【００２６】
内燃機関１０の吸気ポートには、ガソリン等の燃料を噴射するためのインジェクタ２８が配設されている。なお、本発明における内燃機関は、このようなポート噴射式の内燃機関に限定されるものではなく、筒内直接噴射式の内燃機関や、ポート噴射と筒内噴射とを併用する内燃機関であってもよい。
【００２７】
内燃機関１０の排気通路３２には、排気ガスの空燃比を検出する空燃比センサ３３が設置されている。また、排気通路３２には、排気ガスを浄化するための触媒３４が組み込まれている。
【００２８】
本実施形態のシステムは、制御装置として、ECU(Electronic Control Unit)５０を備えている。ECU５０には、上述した各種のセンサからセンサ信号が供給されている。ECU５０は、それらのセンサ信号に基づいて、スロットル弁２４や、インジェクタ２８、点火プラグ３０などの各種アクチュエータを制御することができる。また、後述するように、ECU５０は、クランク角センサ１２、筒内圧センサ１８等のセンサ信号に基づいて、Wiebe関数パラメータの決定処理を行うことができる。更に、ECU５０は、予め記憶されたWiebe関数モデルに基づいて、気筒内の燃焼による熱発生率を推定する処理を行うことができる。
【００２９】
（Wiebe関数）
本実施形態のWiebe関数は、下記（１）式あるいは下記（２）式で表される。
【数３】
2007248119.tif 000005
【００３０】
上記（１）および（２）式中の各記号の意味は、次の通りである。
Ｑ：内燃機関１０の気筒内の発熱量［J］
Ｑtotal：気筒内に供給された燃料の熱量［J］
θ：熱発生開始点を基準とするクランク角度［°］
θｐ：燃焼期間［°］
ｍ：形状パラメータ
ｋ：効率
【００３１】
上記（１）式は、内燃機関１０の気筒内における熱発生率dＱ/dθを算出するためのWiebe関数である。上記（２）式の左辺のＱ/ｋＱtotalは、後述するように、内燃機関１０の気筒内における燃焼割合を表す。つまり、上記（２）式は、燃焼割合Ｑ/ｋＱtotalを算出するためのWiebe関数である。
【００３２】
図２は、横軸にクランク角度を、縦軸に熱発生率dＱ/dθをとったグラフである。上記（１）式のWiebe関数によれば、クランク角度毎の熱発生率dＱ/dθを算出することができる。以下の説明では、このWiebe関数に基づいて算出される熱発生率dＱ/dθを「モデル熱発生率」と称する。図２中の太い実線は、モデル熱発生率を表す。
【００３３】
一方、図１に示すシステムでは、筒内圧センサ１８により、筒内圧をクランク角度毎に計測することができる。以下の説明では、筒内圧センサ１８による計測値を「実筒内圧」と称する。後述するように、実筒内圧からは、エネルギー保存則の式を用いて、クランク角度毎の熱発生率を算出することができる。以下の説明では、実筒内圧に基づいて算出される熱発生率を「実熱発生率」と称する。図２中の太い破線は、実熱発生率を表す。
【００３４】
なお、本明細書において、「クランク角度毎」と言ったとき、そのクランク角度の間隔は、特に限定されるものではなく、例えば１°間隔程度とすることができる。
【００３５】
上記（１）式のWiebe関数によれば、θ＝０のときには、モデル熱発生率dＱ/dθ＝０と算出される。そして、θが０を超えた時点からdＱ/dθ＞０となり、熱発生が開始することを意味する。つまり、上記（１）式のWiebe関数では、θ＝０が熱発生開始点（燃焼開始点）となるので、θは上記のように熱発生開始点を基準とするクランク角度（熱発生開始後の経過クランク角度）を意味することになる。
【００３６】
本明細書では、Wiebe関数によってクランク角度毎のモデル熱発生率を算出する上で値を決定する必要のある複数のパラメータを「Wiebe関数パラメータ」と称する。実機（内燃機関１０）の気筒内の熱発生を、上記Wiebe関数を用いて正確にモデル化する上では、実熱発生率とモデル熱発生率との誤差がなるべく小さくなるように、各Wiebe関数パラメータの値を精度良く決定することが必要となる。以下、Wiebe関数の各パラメータについて順次説明する。
【００３７】
（燃焼期間θｐ）
図２に示すように、燃焼期間θｐは、燃焼による熱の発生が継続する期間、つまり熱発生開始点（燃焼開始点）から熱発生終了点（燃焼終了点）までの期間をクランク角度で表したものとしての物理的意味を有している。Wiebe関数は、熱発生が継続している状態をモデル化するものであるので、その定義域は０≦θ≦θｐとなる。
【００３８】
（形状パラメータｍ）
形状パラメータｍは、図２のようなグラフにおけるモデル熱発生率の形状に大きく関与するパラメータであり、特に、モデル熱発生率が最大値をとるクランク角度に大きく関与するパラメータである。
【００３９】
（効率ｋ）
上記（１）および（２）式中、気筒内に供給された燃料の熱量Ｑtotalは、筒内に供給された燃料量に、その燃料の低位発熱量を乗ずることにより算出することができる。なお、低位発熱量は、真発熱量とも呼ばれる物性値である。低位発熱量とは、単位量の燃料が完全燃焼したときに発生する熱量から、燃料中に含まれる水分および燃焼によって生じる水分を蒸発させるのに必要な熱量（潜熱）を差し引いた残りの熱量を意味する。
【００４０】
内燃機関における燃焼では、冷却損失や燃料の燃え残りなどに起因して、何らかの熱損失を伴うのが普通である。このため、気筒内に供給された燃料の熱量Ｑtotalが全部そのまま気筒内の発熱量Ｑに変換されることは実際上はあり得ない。本実施形態では、このことをWiebe関数に反映させるパラメータとして、効率ｋを導入している。すなわち、効率ｋは、気筒内に供給された燃料の熱量Ｑtotalが発熱量Ｑに変換される効率としての物理的意味を有しており、０＜ｋ＜１なる範囲の数である。
【００４１】
効率ｋを用いると、気筒内での総発熱量は、Ｑtotalにｋを乗じたｋＱtotalとして表すことができる。総発熱量ｋＱtotalは、図２のようなグラフにおいて、熱発生率dＱ/dθの曲線と、dＱ/dθ＝０の直線とが囲む面積で表される。実機における総発熱量ｋＱtotalの値は、実熱発生率を積分することによって算出することが可能である。
【００４２】
上記のように総発熱量がｋＱtotalとして表されるので、あるクランク角度θにおける燃焼割合は、そのクランク角度θまでの発熱量Ｑと、総発熱量ｋＱtotalとの比Ｑ/ｋＱtotalとして表すことができる。この燃焼割合Ｑ/ｋＱtotalは、Wiebe関数によれば、上記（２）式で表される。上記（２）式をθで微分すると、上記（１）式が得られる。つまり、上記（２）式は、上記（１）式と等価である。
【００４３】
（パラメータａ）
燃焼割合の本来の定義からすると、燃焼割合の値は、燃焼終了時、つまりθ＝θｐにおいて１（１００％）となるべきである。一方、Wiebe関数では、燃焼終了時の燃焼割合は、上記（２）式にθ＝θｐを代入することにより、下記式で表される。
Ｑ/ｋＱtotal（燃焼終了時）＝１−exp（−ａ）      ・・・（３）
【００４４】
上記（３）式が示すように、Wiebe関数で算出される燃焼終了時の燃焼割合Ｑ/ｋＱtotalの値は、パラメータａのみの関数となる。そして、この値は１よりも「−exp（−ａ）」だけ小さくなることになり、１になり得ない。そこで、本実施形態では、Wiebe関数で算出される燃焼終了時の燃焼割合が１近傍の値（ここでは０．９９９とする）になると仮定することとした。この仮定によれば、下記式が得られる。
Ｑ/ｋＱtotal（燃焼終了時）＝１−exp（−ａ）＝０．９９９      ・・・（４）
上記（４）式を方程式として解くと、ａ≒６．９となる。そこで、本実施形態では、パラメータａを定数６．９として扱うものとする。
【００４５】
以上、Wiebe関数パラメータとして、燃焼期間θｐ、形状パラメータｍ、効率ｋ、パラメータａの４つについて説明した。Wiebe関数によってクランク角度毎のモデル熱発生率を算出する上では、次の理由から、もう一つのWiebe関数パラメータを設定する必要がある。前述したように、上記（１）式のWiebe関数中のクランク角度θは、Wiebe関数の熱発生開始点を基準としている。このため、クランク角度毎のモデル熱発生率を算出するためには、Wiebe関数のクランク角度θと、絶対的なクランク角度（例えば上死点基準クランク角度）との関係を定めるためのWiebe関数パラメータが更に必要となる。このようなWiebe関数パラメータとして、一つには、着火遅れθｂを想定することができる。
【００４６】
点火プラグ３０により筒内の混合気に点火された後、燃焼が開始して実熱発生率が立ち上がるまでには、若干の遅れが生ずる。このため、図２に示すように、Wiebe関数によるモデル熱発生率の波形を実熱発生率の波形に精度良く一致させるためには、Wiebe関数の熱発生開始点を点火時期SAよりも遅らせた位置にする必要がある。着火遅れθｂは、この要請に応えるためのパラメータであり、Wiebe関数の熱発生開始点の、点火時期SAからの遅れをクランク角度で表した値として定義される。
【００４７】
上記着火遅れθｂをWiebe関数パラメータとして追加すれば、Wiebe関数の熱発生開始点の絶対的なクランク角度を定めることができるので、クランク角度毎のモデル熱発生率を算出することが可能となる。しかしながら、図２に示すように、実熱発生率は、立ち上がり部分がなだらかな形状になり易いのに対し、モデル熱発生率は、熱発生開始点から一気に立ち上がるような形状になり易い。このため、実熱発生率とモデル熱発生率とは、熱発生開始直後の形状が一致しない場合がある。このような場合には、着火遅れθｂの値を一義的に決めることができず、着火遅れθｂを精度良く決定することが困難になり易い。このようなことから、着火遅れθｂをWiebe関数パラメータとして使用することとすると、着火遅れθｂの値を十分な精度で決定しにくくなり、その結果、正確なモデル化が困難となり易い。
【００４８】
そこで、本実施形態では、着火遅れθｂをWiebe関数パラメータとして使用することはせず、代わりに、実熱発生率が最大となった時点でのクランク角度（以下、「熱発生率最大クランク角度」と称する）をWiebe関数パラメータとして使用することとした。以下の説明では、熱発生率最大クランク角度を、圧縮上死点を基準とするクランク角度で表すものとし、記号をＣＡｍ［°ATDC］と表記する。
【００４９】
熱発生率最大クランク角度ＣＡｍは、実熱発生率がピークになったときのクランク角度であるので、実熱発生率の波形から一義的に決めることが容易である。このため、熱発生率最大クランク角度ＣＡｍをWiebe関数パラメータとして使用することにすれば、高い精度で容易にその値を決定することができる。また、後述するように、熱発生率最大クランク角度ＣＡｍを用いれば、着火遅れθｂを用いなくても、Wiebe関数の熱発生開始点および終了点を定めることができるので、クランク角度毎のモデル熱発生率を算出することが可能となる。
【００５０】
［実施の形態１におけるWiebe関数パラメータ決定方法］
以下、実施の形態１における具体的処理について説明する。図３は、実機データに基づいてWiebe関数パラメータの値を決定するために、本実施形態においてECU５０が実行するルーチンのフローチャートである。
【００５１】
図３に示すルーチンによれば、まず、内燃機関１０の運転中に、クランク角センサ１２および筒内圧センサ１８の出力に基づいて、クランク角度毎の実筒内圧Ｐdata［Pa］が計測される（ステップ１００）。このステップ１００において、実筒内圧Ｐdataの計測は、全クランク角度範囲に渡って行わなくてもよく、少なくとも実際の燃焼期間が含まれる範囲で行えばよい。
【００５２】
次に、上記ステップ１００で実測されたクランク角度毎の実筒内圧Ｐdataに基づいて、次のような方法により、クランク角度毎の実熱発生率が算出される（ステップ１０２）。次式は、熱発生率dＱ/dθと、筒内圧Ｐと、筒内容積Ｖ［m3］と、クランク角度θとを用いて、熱力学上のエネルギー保存則を表した式である。なお、次式中のクランク角度θは、上記（１）式等のWiebe関数におけるθとは異なり、上死点基準の絶対的なクランク角度を表すものとする。
【数４】
2007248119.tif 000006
【００５３】
上記（５）式中、κは比熱比であり、燃焼ガスの組成などに基づいて定まる既知の値である。また、筒内容積Ｖおよびその変化率dＶ/dθは、クランク角度θに応じて幾何学的に決まる値である。よって、これらの既知の値と、上記ステップ１００で実測されたクランク角度毎の実筒内圧Ｐdataの値とを上記（５）式に代入することにより、図２中の太い破線で示すようなクランク角度毎の実熱発生率dＱ/dθを算出することができる。
【００５４】
なお、上記ステップ１０２の実熱発生率の算出に際しては、例えば、内燃機関１０の冷却水温等の温度情報や排気ガス成分の情報等を加味して、それらに基づく補正を施して算出するようにしてもよい。
【００５５】
続いて、上記ステップ１０２で算出されたクランク角度毎の実熱発生率のうちの最大値が探索される。そして、実熱発生率が最大値をとった時点でのクランク角度の値が取得され、その値が熱発生率最大クランク角度ＣＡｍとして決定される（ステップ１０４）。
【００５６】
次いで、上記ステップ１０２で算出されたクランク角度毎の実熱発生率に基づいて、総発熱量ｋＱtotalが算出される（ステップ１０６）。具体的には、上記ステップ１０２で算出されたクランク角度毎の実熱発生率を積分（積算）することにより、総発熱量ｋＱtotalが算出される。
【００５７】
続いて、実熱発生率が最大となった時点での燃焼割合、つまり熱発生率最大クランク角度ＣＡｍでの燃焼割合が算出される（ステップ１０８）。以下の説明では、熱発生率最大クランク角度ＣＡｍでの燃焼割合を「実燃焼割合」と称し、その値を記号αで表す。このステップ１０８においては、具体的には、まず、上記ステップ１０２で算出されたクランク角度毎の実熱発生率を熱発生率最大クランク角度ＣＡｍまで積分（積算）することにより、熱発生率最大クランク角度ＣＡｍまでの発熱量の値が算出され、次に、その値を上記ステップ１０６で算出された総発熱量ｋＱtotalで除算した値が、実燃焼割合αとして算出される。
【００５８】
次に、上記ステップ１０８で算出された実燃焼割合αの値を用いて、形状パラメータｍの値を決定する処理が行われる（ステップ１１０）。この処理は、以下に説明するような考え方に基づいて行われる。
【００５９】
実燃焼割合αに相当するWiebe関数上の値、つまり、モデル熱発生率が最大となるときの燃焼割合の値は、図２のモデル熱発生率のグラフが囲む面積のうち、モデル熱発生率が最大となったクランク角度までの面積（図２中のハッチング部分）を、全部の面積で除した値として算出することができる。
【００６０】
この値を算出するには、まず、上記（１）式をθで微分し、d２Ｑ/dθ２＝０とおくことにより、モデル熱発生率dＱ/dθが最大となるクランク角度θ＊を求める。下記の式展開を行うことにより、クランク角度θ＊を表す式が得られる。
【数５】
2007248119.tif 000007
【００６１】
上記（６）式で表されるクランク角度θ＊を上記（２）式に代入することにより、モデル熱発生率dＱ/dθが最大となるクランク角度θ＊での燃焼割合は、下記（７）式で表すことができる。
【数６】
2007248119.tif 000008
【００６２】
上記（７）式に示されるように、モデル熱発生率dＱ/dθが最大となるクランク角度θ＊での燃焼割合の値は、形状パラメータｍのみの関数として表される。この値、つまり上記（７）式の右辺と、対応する実機データである実燃焼割合αとを等しいと置いた式をｍについて解くと、下記（８）式が得られる。
【数７】
2007248119.tif 000009
【００６３】
上記ステップ１１０においては、以上説明した考え方に基づき、上記ステップ１０８で算出された実燃焼割合αの値を上記（８）式に代入して算出された値が、形状パラメータｍの値として決定される。
【００６４】
ところで、本実施形態では、上記のようにして形状パラメータｍの値が決定されると、以下に説明するように、Wiebe関数の熱発生開始点（燃焼開始点）および熱発生終了点（燃焼終了点）の絶対的なクランク角度を決定することが可能となる。以下では、Wiebe関数の熱発生開始点の、圧縮上死点を基準とした絶対的なクランク角度を記号ＣＡｓ［°BTDC］で表す。同様に、Wiebe関数の熱発生終了点の、圧縮上死点を基準とした絶対的なクランク角度を記号ＣＡｅ［°ATDC］で表す。
【００６５】
まず、上記（６）式の両辺を燃焼期間θｐで除した値をβとおく。つまり、βを下記式で定義する。
【数８】
2007248119.tif 000010
【００６６】
上記（９）式が示すように、βは、パラメータａおよび形状パラメータｍのみで表される。本実施形態では、前述した通りａ＝６．９として扱うことにしている。また、形状パラメータｍの値は、上記ステップ１１０の処理で決定されている。よって、上記βの値は、ここまでの処理によって決定済みの値として扱うことができる。そして、上記（９）式より、下記式が成り立つ。
θ＊＝β×θｐ      ・・・（１０）
【００６７】
また、図２に示すように、モデル熱発生率が最大となるクランク角度θ＊から、熱発生終了点までの期間をθresとおくと、θresの定義および上記（１０）式より、下記式が成り立つ。
θres＝θｐ−θ＊＝（１−β）×θｐ      ・・・（１１）
【００６８】
本実施形態では、実熱発生率が最大となる熱発生率最大クランク角度ＣＡｍと、モデル熱発生率が最大となるクランク角度θ＊とが一致するように、Wiebe関数パラメータを決定する。圧縮上死点基準の熱発生率最大クランク角度ＣＡｍと、熱発生開始点基準のクランク角度θ＊とを一致させるようにすると、ＣＡｍ、θ＊、θres、熱発生開始点ＣＡｓ、熱発生終了点ＣＡｅの間に図２に示す関係が成り立つ。この関係と、上記（１０）式および（１１）式とから、下記式を導くことができる。
ＣＡｓ＝θ＊−ＣＡｍ＝β×θｐ−ＣＡｍ      ・・・（１２）
ＣＡｅ＝ＣＡｍ＋θres＝ＣＡｍ＋（１−β）×θｐ      ・・・（１３）
【００６９】
上記（１２）式および（１３）式が示すとおり、Wiebe関数の燃焼開始点および燃焼終了点の絶対的なクランク角度ＣＡｓおよびＣＡｅは、熱発生率最大クランク角度ＣＡｍと、燃焼期間θｐと、決定済みの値であるβとで表すことができる。このように、熱発生率最大クランク角度ＣＡｍをWiebe関数パラメータとして使用すれば、Wiebe関数の熱発生開始点ＣＡｓを決定することができる。よって、着火遅れθｂを使用せずに、上記（１）式のWiebe関数によってクランク角度毎のモデル熱発生率を算出することが可能となる。
【００７０】
なお、定義より、ＣＡｓ＋ＣＡｅ＝θｐである。この関係によれば、上記（１２）式および（１３）式は互いに等価な式である。
【００７１】
さて、ここまでの処理により、Wiebe関数パラメータとして、熱発生率最大クランク角度ＣＡｍと形状パラメータｍとが決定された。残りの決定すべきWiebe関数パラメータは、効率ｋおよび燃焼期間θｐである。本実施形態では、上記ステップ１１０の処理に続いて、効率ｋおよび燃焼期間θｐの値を最小二乗法により決定する処理が行われる（ステップ１１２）。
【００７２】
具体的には、まず、上記（１）式のWiebe関数の効率ｋおよび燃焼期間θｐに仮の値を代入し、パラメータａ、熱発生率最大クランク角度ＣＡｍおよび形状パラメータｍには上記の決定済みの値を使用して、クランク角度毎のモデル熱発生率が算出される。次いで、このクランク角度毎のモデル熱発生率と、クランク角度毎の実熱発生率（上記ステップ１０２で算出されたもの）との誤差（偏差）の二乗を熱発生開始点ＣＡｓから熱発生終了点ＣＡｅまで足し合わせた値（誤差二乗和）が算出される。この誤差二乗和の算出処理が、効率ｋおよび燃焼期間θｐの仮の値を変化させながら、繰り返し行われる。そして、誤差二乗和を最小とするような効率ｋおよび燃焼期間θｐの仮の値の組み合わせが、効率ｋおよび燃焼期間θｐの最終的な値として決定される。このような処理によれば、実熱発生率とモデル熱発生率との誤差が最小となるように、効率ｋおよび燃焼期間θｐの値を決定することができる。
【００７３】
以上説明した通り、本実施形態では、着火遅れθｂに代えて、熱発生率最大クランク角度ＣＡｍをWiebe関数パラメータとして使用することとした。熱発生率最大クランク角度ＣＡｍは、実熱発生率の波形から容易に一義的に決定することができるので、着火遅れθｂに比して、その値を高精度に定めることができる。このため、上述した手法によれば、内燃機関１０の気筒内の熱発生をモデル化するためのWiebe関数モデルを、高精度かつ容易に構築することができる。
【００７４】
また、本発明者の知見によれば、熱発生率最大クランク角度ＣＡｍは、点火時期SAを進角するにつれて（ＭＢＴに近づけるにつれて）、連続的に進角するという性質がある。つまり、熱発生率最大クランク角度ＣＡｍには、内燃機関１０の運転条件の影響、特に点火時期SAの影響が精度良く反映されると言える。このため、熱発生率最大クランク角度ＣＡｍをWiebe関数パラメータとすることは、理に適っていると言え、高精度なWiebe関数モデルを構築する上で有利に作用する。
【００７５】
また、本実施形態では、上記ステップ１１０の処理により、実燃焼割合αの値のみから形状パラメータｍの値を直接に決定することができる。本発明者の知見によれば、実熱発生率が最大となった時点での実燃焼割合αは、熱発生率最大クランク角度ＣＡｍと同様に、点火時期SAの変化に対して連続的に変化する。つまり、実燃焼割合αの値には、内燃機関１０の運転条件の影響、特に点火時期SAの影響が精度良く反映されると言える。このため、本実施形態のように、実燃焼割合αから形状パラメータｍの値を決定することは、理に適っていると言え、高精度なWiebe関数モデルを構築する上で有利に作用する。
【００７６】
なお、図３に示すルーチンのステップ１１２においては、誤差二乗和を誤差評価値とする最小二乗法によって効率ｋおよび燃焼期間θｐの値を決定しているが、最小二乗法以外の最適化手法によってそれらの値を決定するようにしてもよい。また、本発明では、効率ｋおよび燃焼期間θｐの値は、上記ステップ１１２のような最適化手法以外の手法で決定することとしてもよい。例えば、燃料の低位発熱量と、インジェクタ２８からの燃料噴射量等に基づいて、気筒内に供給された燃料の熱量Ｑtotalを算出し、上記ステップ１０６で算出された総発熱量ｋＱtotalをこのＱtotalで除することによって効率ｋの値を決定してもよい。
【００７７】
［実施の形態１における内燃機関の熱発生率推定装置］
図３に示すルーチンの処理により決定される各Wiebe関数パラメータの値は、内燃機関１０の運転条件（本実施形態では、点火時期SA、機関回転数NE、空燃比A/F、負荷率KLの４つとする）に応じて変化する。そこで、点火時期SA、機関回転数NE、空燃比A/F、負荷率KLをそれぞれ変化させ、各々の運転条件の下で図３に示すルーチンの処理を行ってデータを収集および解析することにより、各Wiebe関数パラメータと内燃機関１０の運転条件との関係を把握することができる。
【００７８】
そのようにして把握された関係は、関数等の演算式あるいはマップの形にまとめることができる。ここでは、熱発生率最大クランク角度ＣＡｍ、形状パラメータｍ、効率ｋおよび燃焼期間θｐの各々と運転条件との関係が、点火時期SA、機関回転数NE、空燃比A/Fおよび負荷率KLを変数とする下記の関数ｆ１〜ｆ４としてまとめられたものとする。
ＣＡｍ＝ｆ１（SA，NE，A/F，KL）
ｍ＝ｆ２（SA，NE，A/F，KL）
ｋ＝ｆ３（SA，NE，A/F，KL）
θｐ＝ｆ４（SA，NE，A/F，KL）
【００７９】
上記関数ｆ１〜ｆ４のような関係をECU５０に記憶させておけば、あらゆる運転条件について、クランク角度毎のモデル熱発生率をWiebe関数モデルに基づいて算出可能な熱発生率推定装置を構成することができる。
【００８０】
図４は、そのように構成された熱発生率推定装置において、ECU５０が、熱発生率の推定値であるモデル熱発生率を算出するために実行するルーチンのフローチャートである。図４に示すルーチンによれば、まず、熱発生率を推定すべき運転条件の点火時期SA，機関回転数NE，空燃比A/F，負荷率KLの値を、ECU５０に記憶された上記の関数ｆ１〜ｆ４に代入することにより、その運転条件に対応する熱発生率最大クランク角度ＣＡｍ、形状パラメータｍ、効率ｋおよび燃焼期間θｐの各々の値が取得される（ステップ１２０）。
【００８１】
続いて、Wiebe関数の熱発生開始点のクランク角度ＣＡｓが算出される（ステップ１２２）。具体的には、まず、上記ステップ１２０で取得された形状パラメータｍの値に基づき、ａ＝６．９として、上記（９）式に従ってβの値が算出される。次いで、このβの値と、上記ステップ１２０で取得された熱発生率最大クランク角度ＣＡｍおよび燃焼期間θｐの値に基づき、上記（１２）式に従って、熱発生開始点のクランク角度ＣＡｓが算出される。
【００８２】
そして、上記ステップ１２０および１２２で取得された値を用いることにより、熱発生率の推定値として、クランク角度毎のモデル熱発生率が算出される（ステップ１２４）。具体的には、上記ステップ１２０で算出された形状パラメータｍ、効率ｋおよび燃焼期間θｐの値を上記（１）式に代入し、熱発生開始点のクランク角度が上記ステップ１２２で算出されたＣＡｓになるようにして、クランク角度毎のモデル熱発生率が算出される。なお、上記（１２）式と上記（１３）式とは等価なので、熱発生終了点がＣＡｓになるようにすれば、熱発生終了点は上記（１３）式のＣＡｅに必然的に一致する。
【００８３】
また、上記ステップ１２４の処理において、モデル熱発生率の熱発生開始点をＣＡｓとすることにより、ＣＡｍ、θ＊、ＣＡｓの間に図２に示す関係が成り立つ。このため、上記ステップ１２４で算出されるモデル熱発生率が最大となるクランク角度θ＊は、上記ステップ１２０で取得された熱発生率最大クランク角度ＣＡｍに一致する。このため、図４に示すルーチンの処理によれば、推定された熱発生率が最大となるクランク角度を、推定対象の運転条件にかかわらず、実機データから得られた傾向に常に一致させることができる。よって、あらゆる運転条件について、気筒内の熱発生率を精度良く推定することができる。
【００８４】
ところで、クランク角度毎の熱発生率の推定値が得られれば、これと上記（５）式のエネルギー保存則とを組み合わせることにより、クランク角度毎の筒内圧の推定値を算出することができる。更に、クランク角度毎の筒内圧の値からは、内燃機関１０の図示トルクを算出することができる。本発明によれば、図４に示すルーチンの処理によって熱発生率を精度良く推定することができるので、その推定された熱発生率を用いて、あらゆる運転条件の下での筒内圧やトルクを精度良く推定することが可能となる。
【００８５】
上記のようにして得られた内燃機関１０の筒内圧やトルクの推定結果は、内燃機関１０の開発段階などで有効に利用することができる。すなわち、従来より、内燃機関１０の開発段階において制御ロジックの検討などを行う場合には、様々な運転条件の下で筒内圧やトルクなどを実測する作業が必要である。これに対し、本発明によれば、Wiebe関数パラメータの値と運転条件との関係を正確に把握して高精度なWiebe関数モデルを構築することが可能である。そして、そのWiebe関数モデルを記憶させた上述の熱発生率推定装置を用いれば、筒内圧やトルクの実測作業をWiebe関数モデルに基づく予測で置き換えることが可能となる。このため、実機データの計測点数を大幅に削減することができ、開発工数、開発コストを大幅に削減することができる。
【００８６】
また、そのようなWiebe関数モデルは、車載ECUに記憶させておくこともできる。Wiebe関数モデルを車載ECUに記憶させておけば、車両に搭載された内燃機関１０について、種々の運転条件下における筒内圧やトルクを、精度良く予測することができる。そして、その予測結果を利用すれば、例えばトルクデマンド制御やドライバビリティ改善制御などの各種制御を車載ECUが精度良く行うことが可能となる。
【００８７】
また、図３に示すWiebe関数パラメータを決定するためのルーチンを車載ECUに実行させることも可能である。図３に示すルーチンを車載ECUで実行することによってWiebe関数パラメータを算出（決定）し、その算出されたWiebe関数パラメータの値を現在の運転状態と関係付けて学習することとすれば、車載ECUに記憶されたWiebe関数モデルを現在の内燃機関１０の状態に応じて修整することが可能となる。つまり、内燃機関１０の特性の経時変化や個体差を反映させて、車載ECUに記憶されたWiebe関数モデルをより高精度に修整することができる。
【００８８】
なお、以上の説明では、前記のような各種の数式を利用したが、本発明では、それらの数式は前記のものと完全に同じである必要はなく、それらと等価な数式を利用してもよい。
【００８９】
また、上述した実施の形態１においては、ECU５０が前記第４および第５の発明における「関係記憶手段」に相当している。また、ECU５０が、図４に示すルーチンの処理を実行することにより、前記第４および第５の発明における「熱発生率推定手段」が実現されている。
【図面の簡単な説明】
【００９０】
【図１】本発明の実施の形態１のシステム構成を説明するための図である。
【図２】横軸にクランク角度を、縦軸に熱発生率dＱ/dθをとったグラフである。
【図３】本発明の実施の形態１において実行されるルーチンのフローチャートである。
【図４】本発明の実施の形態１において実行されるルーチンのフローチャートである。
【符号の説明】
【００９１】
１０  内燃機関
１２  クランク角センサ
１４  吸気弁
１６  排気弁
１８  筒内圧センサ
２１  吸気圧センサ
３０  点火プラグ
３３  空燃比センサ
３４  触媒
５０  ECU(Electronic Control Unit)
【図１】
2007248119.tif 000013
【図２】
2007248119.tif 000014
【図３】
2007248119.tif 000015
【図４】
2007248119.tif 000016