以下、本発明による内燃機関(ディーゼル機関)の混合気状態推定装置を含んだエミッション発生量推定装置の各実施形態について図面を参照しつつ説明する。
(第1実施形態)
図1は、本発明の第1実施形態に係る内燃機関のエミッション発生量推定装置を4気筒内燃機関(ディーゼル機関)10に適用したシステム全体の概略構成を示している。このシステムは、燃料供給系統を含むエンジン本体20、エンジン本体20の各気筒の燃焼室(筒内)にガスを導入するための吸気系統30、エンジン本体20からの排ガスを放出するための排気系統40、排気還流を行うためのEGR装置50、及び電気制御装置60を含んでいる。
エンジン本体20の各気筒の上部には燃料噴射弁(噴射弁、インジェクタ)21が配設されている。各燃料噴射弁21は、図示しない燃料タンクと接続された燃料噴射用ポンプ22に燃料配管23を介して接続されている。燃料噴射用ポンプ22は、電気制御装置60と電気的に接続されていて、同電気制御装置60からの駆動信号(後述する指令最終燃料噴射圧力Pcrfinに応じた指令信号)により燃料の実際の噴射圧力(吐出圧力)が同指令最終燃料噴射圧力Pcrfinになるように同燃料を昇圧するようになっている。
これにより、燃料噴射弁21には、燃料噴射用ポンプ22から前記指令最終燃料噴射圧力Pcrfinまで昇圧された燃料が供給されるようになっている。また、燃料噴射弁21は、電気制御装置60と電気的に接続されていて、同電気制御装置60からの駆動信号(指令燃料噴射量(質量)Qfinに応じた指令信号)により噴射期間TAUだけ開弁し、これにより各気筒の燃焼室内に前記指令最終燃料噴射圧力Pcrfinにまで昇圧された燃料を前記指令燃料噴射量Qfinだけ直接噴射するようになっている。
吸気系統30は、エンジン本体20の各気筒の燃焼室にそれぞれ接続された吸気マニホールド31、吸気マニホールド31の上流側集合部に接続され同吸気マニホールド31とともに吸気通路を構成する吸気管32、吸気管32内に回動可能に保持されたスロットル弁33、電気制御装置60からの駆動信号に応答してスロットル弁33を回転駆動するスロットル弁アクチュエータ33a、スロットル弁33の上流において吸気管32に順に介装されたインタクーラー34と過給機35のコンプレッサ35a、及び吸気管32の先端部に配設されたエアクリーナ36とを含んでいる。
排気系統40は、エンジン本体20の各気筒にそれぞれ接続された排気マニホールド41、排気マニホールド41の下流側集合部に接続された排気管42、排気管42に配設された過給機35のタービン35b、及び排気管42に介装されたディーゼルパティキュレートフィルタ(以下、「DPNR」と称呼する。)43を含んでいる。排気マニホールド41及び排気管42は排気通路を構成している。
EGR装置50は、排気ガスを還流させる通路(EGR通路)を構成する排気還流管51と、排気還流管51に介装されたEGR制御弁52と、EGRクーラー53とを備えている。排気還流管51はタービン35bの上流側排気通路(排気マニホールド41)とスロットル弁33の下流側吸気通路(吸気マニホールド31)を連通している。EGR制御弁52は電気制御装置60からの駆動信号に応答し、再循環される排気ガス量(排気還流量、EGRガス流量)を変更し得るようになっている。
電気制御装置60は、互いにバスで接続されたCPU61、CPU61が実行するプログラム、テーブル(ルックアップテーブル、マップ)、及び定数等を予め記憶したROM62、CPU61が必要に応じてデータを一時的に格納するRAM63、電源が投入された状態でデータを格納するとともに同格納したデータを電源が遮断されている間も保持するバックアップRAM64、並びにADコンバータを含むインターフェース65等からなるマイクロコンピュータである。
インターフェース65は、吸気管32に配置された熱線式エアフローメータ71、スロットル弁33の下流であって排気還流管51が接続された部位よりも下流の吸気通路に設けられた吸気温センサ72、スロットル弁33の下流であって排気還流管51が接続された部位よりも下流の吸気通路に配設された吸気管圧力センサ73、クランクポジションセンサ74、アクセル開度センサ75、燃料噴射用ポンプ22の吐出口の近傍の燃料配管23に配設された燃料温度センサ76、気筒毎にそれぞれ配置された着火時期取得手段としての筒内圧力センサ77、及び、スロットル弁33の下流であって排気還流管51が接続された部位よりも下流の吸気通路に配設された吸気酸素濃度センサ78と接続されていて、これらのセンサからの信号をCPU61に供給するようになっている。また、インターフェース65は、燃料噴射弁21、燃料噴射用ポンプ22、スロットル弁アクチュエータ33a、及びEGR制御弁52と接続されていて、CPU61の指示に応じてこれらに駆動信号を送出するようになっている。
熱線式エアフローメータ71は、吸気通路内を通過する吸入空気の質量流量(単位時間当りの吸入空気量、単位時間あたりの新気量)を計測し、同質量流量Ga(空気流量Ga)を表す信号を発生するようになっている。吸気温センサ72は、エンジン10のシリンダ(即ち、燃焼室、筒内)に吸入されるガスの温度(即ち、吸気温度)を検出し、同吸気温度Tbを表す信号を発生するようになっている。吸気管圧力センサ73は、エンジン10のシリンダに吸入されるガスの圧力(即ち、吸気管圧力)を検出し、同吸気管圧力Pbを表す信号を発生するようになっている。
クランクポジションセンサ74は、各気筒の絶対クランク角度を検出し、実クランク角度CAactを表すとともにエンジン10の回転速度であるエンジン回転速度NEをも表す信号を発生するようになっている。アクセル開度センサ75は、アクセルペダルAPの操作量を検出し、アクセル操作量Accpを表す信号を発生するようになっている。燃料温度センサ76は、燃料配管23を通過する燃料の温度を検出し、燃料温度Tcrを表す信号を発生するようになっている。各筒内圧力センサ77は、燃焼室内のガスの圧力(従って、上記筒内ガスの圧力)を検出し、筒内ガス圧力Paを表す信号を発生するようになっている。この筒内圧力センサ77は、後述するように着火時点を検出するためにのみ使用される。吸気酸素濃度センサ78は、吸気中の酸素濃度を検出し、吸気酸素濃度RO2inを表す信号を発生するようになっている。
(混合気状態の推定方法の概要)
次に、上記のように構成された混合気状態推定装置を含んだエミッション発生量推定装置(以下、「本装置」と云う。)による混合気状態の推定方法について説明する。
図2は、或る一つの気筒のシリンダ内(筒内、燃焼室内)に吸気マニホールド31からガスが吸入され、燃焼室内に吸入されたガスが排気マニホールド41へ排出される様子を模式的に示した図である。図2に示したように、燃焼室内に吸入されるガス(従って、筒内ガス)には、吸気管32の先端部からスロットル弁33を介して吸入された新気と、排気還流管51からEGR制御弁52を介して吸入されたEGRガスが含まれる。
吸入される新気量(質量)と吸入されるEGRガス量(質量)の和に対するEGRガス量の割合(即ち、EGR率)は、運転状態に応じて電気制御装置60(CPU61)により適宜制御されるスロットル弁33の開度、及びEGR制御弁52の開度に応じて変化する。
かかる新気、及びEGRガスは、吸気行程において開弁している吸気弁Vinを介してピストンの下降に伴って燃焼室内に吸入されて筒内ガスとなる。筒内ガスは、ピストンが圧縮下死点に達する時点近傍で吸気弁Vinが閉弁することにより燃焼室内に密閉され、その後の圧縮行程においてピストンの上昇に伴って圧縮される。
そして、ピストンが圧縮上死点近傍に達すると(具体的には、後述する燃料噴射開始時期(クランク角度)CAinjが到来すると)、本装置は、前記指令燃料噴射量Qfinに応じた噴射期間TAUだけ燃料噴射弁21を開弁することで燃料を燃焼室内に直接噴射する。この結果、燃料噴射弁21の噴孔から噴射された(液体の)燃料は、圧縮により高温になっている筒内ガスから受ける熱により直ちに燃料蒸気になるとともに、時間の経過に伴って同筒内ガスを取り込みながら混合気となって燃焼室内において円錐状に拡散していく。
上述したように、指令燃料噴射量Qfinの燃料は、実際には、上記燃料噴射開始時期CAinjから噴射期間TAUだけ連続して噴射されるが、以下においては便宜上、指令燃料噴射量Qfinの燃料が燃料噴射開始時期CAinjにて一時(瞬時)に噴射されるものとして説明を続ける。
図3(a)は、燃料噴射弁21の噴孔から、噴射期間TAUに相当する指令燃料噴射量(質量)Qfinの燃料が一時に噴射された時点(即ち、噴射後経過時間t=0)での質量Qfinの燃料(燃料蒸気)の様子を模式的に示した図である。図3(b)は、その後の或る時点(任意の噴射後経過時間t)での図3(a)に示した質量Qfinの燃料蒸気の様子を模式的に示した図である。
図3(b)に示すように、質量Qfinの燃料蒸気は、燃料噴射開始時期CAinj(即ち、噴射後経過時間t=0)において噴射された後、噴霧角θをもって円錐状に拡散しながら筒内ガスを順次取り込んでいく。そして、質量Qfinの燃料蒸気は、任意の噴射後経過時間tにおいて、噴射後経過時間tの関数である質量Gの筒内ガス(以下、「混合気形成筒内ガス」と云うこともある。)と混ざり合って質量(Qfin+G)の混合気となっているものと仮定する。
本装置は、この混合気の任意の噴射後経過時間tにおける状態を推定する。混合気の状態としては、後述するエミッションの発生量の推定に必要となる混合気温度Tmix、混合気内の燃料濃度[Fuel]mix、混合気内の酸素濃度[O2]mix、混合気内の窒素濃度[N2]mixが推定される。以下、先ず、係る混合気の状態推定に必要となる、任意の噴射後経過時間tにおける空気過剰率λの取得方法について説明する。
<空気過剰率λの取得>
噴射後経過時間tにおける空気過剰率λは、下記(1)式に示すように定義される。下記(1)式において、stoichは、単位質量の燃料の燃焼に必要な筒内ガスの質量(以下、「筒内ガス理論空燃比stoich」と呼ぶ。)である。筒内ガス理論空燃比stoichの値は、吸気中の酸素濃度に応じて変化すると考えられるから、上記吸気酸素濃度RO2inを引数とする所定の関数に従って取得され得る。
このように定義される空気過剰率λは、例えば、日本機械学会論文集 25-156(1959年),820ページ 「ディーゼル機関の噴霧到達距離に関する研究」 和栗雄太郎,藤井勝,網谷竜夫,恒屋礼次郎 (以下、「非特許文献1」と称呼する。)にて紹介された実験式である下記(2)式、及び下記(3)式に基づいて噴射後経過時間tの関数として求めることができる。
上記(3)式において、tは上記噴射後経過時間であり、dλ/dtは噴射後経過時間tの関数である燃料希釈率である。また、cは収縮係数、dは燃料噴射弁21の噴孔径、ρfは(液体の)燃料密度、Lは論理希釈ガス量であって、これらの各値は全て定数である。
上記(3)式において、ΔPは有効噴射圧力であって、上記最終燃料噴射圧力Pcrfinから噴射開始時点(即ち、噴射後経過時間t=0)での筒内ガス圧力Pg0を減じた値である。筒内ガス圧力Pg0は、圧縮行程(及び膨張行程)における筒内ガスの状態が吸気弁Vinの閉弁時(即ち、筒内ガスが密閉された時点。以下、「IVC」と呼ぶ。)以降断熱変化するとの仮定のもと、下記(4)式に従って求めることができる。
上記(4)式において、PgivcはIVCにおける筒内ガス圧力である。上述したように、IVCは圧縮下死点近傍であるから、IVCにおいて筒内ガス圧力は吸気管圧力Pbと略等しいと考えられる。従って、IVCにおいて吸気管圧力センサ73により検出される吸気管圧力PbがPgivcとして使用され得る。Vg(CAivc)はIVCにおけるクランク角度CAに対応する筒内容積であり、Vg(CAinj)は噴射後経過時間t=0におけるクランク角度CAに対応する筒内容積である。筒内容積Vgは機関10の設計諸元に基づいてクランク角度CAの関数Vg(CA)として取得することができるから、Vg(Caivc),Vg(CAinj)も取得することができる。κは筒内ガスの比熱比(本例では、一定)である。
また、上記(3)式において、θは図3(b)に示した噴霧角である。噴霧角θは、噴射開始時点(即ち、噴射後経過時間t=0)における筒内ガスの密度ρg0、及び上記有効噴射圧力ΔPに応じて変化すると考えられるから、筒内ガスの密度ρg0、及び有効噴射圧力ΔPと噴霧角θとの関係を予め規定したテーブルMapθに基づいて取得することができる。筒内ガスの密度ρg0は、筒内ガスの全質量Mgを、噴射後経過時間t=0における上記筒内容積Vg(CAinj)で除することで取得することができる。筒内ガスの全質量Mgは、IVCにおける気体の状態方程式に基づく下記(5)式に従って取得され得る。下記(5)式において、TgivcはIVCにおける筒内ガス温度である。IVCは圧縮下死点近傍であるから、IVCにおいて筒内ガス温度は吸気温度Tbと略等しいと考えられる。従って、IVCにおいて吸気温センサ72により検出される吸気温度TbがTgivcとして使用され得る。Rは筒内ガスのガス定数(本例では、一定)である。
また、上記(3)式において、ρgは噴射後経過時間tにおける筒内ガス密度であって、前記筒内ガスの全質量Mgを、噴射後経過時間tにおける上記筒内容積Vg(CA)で除することで、噴射後経過時間tの関数として取得することができる。
以上のことから、有効噴射圧力ΔPと噴霧角θとを上述のようにして求めれば、噴射後経過時間tの値と同噴射後経過時間tの関数である筒内ガス密度ρgの値とを使用して、上記(3)式に従って燃料希釈率dλ/dtが噴射後経過時間tの関数として求められる。そして、噴射後経過時間t=0から微小時間Δt(例えば、0.1msec)毎に求めた燃料希釈率dλ/dtの値を上記(2)式に従って時間で積分(積算)していくことで噴射後経過時間tにおける空気過剰率λを噴射後経過時間t=0から微小時間Δt毎に取得することができる。
なお、上記(3)式から取得される燃料希釈率dλ/dtの値は常に正の値となることから上記(2)式から取得される空気過剰率λの値は噴射後経過時間tの増大に従って増加していく。そうすると、上記(1)式から理解できるように、混合気形成筒内ガスの質量Gが噴射後経過時間tの増大に従って増加していく。このことは、噴射後の燃料蒸気が円錐状に拡散していくことに伴って燃料蒸気と混ざり合う(燃料蒸気が取り込む)筒内ガス(従って、混合気形成筒内ガス)の量が増大していくことに対応している。
<混合気温度Tmixの取得>
次に、上述のように取得される空気過剰率λの値を利用して任意の噴射後経過時間tにおける混合気温度Tmixを取得する方法について説明する。一般に、混合気の熱エネルギー(エンタルピ)Hmixは、混合気温度Tmixを用いて下記(6)式に従って表すことができる。
上記(6)式において、Mmixは混合気の総質量(混合気質量)、Cmixは混合気の定圧比熱である。従って、混合気のエンタルピHmix、混合気質量Mmix、及び混合気の定圧比熱Cmixをそれぞれ噴射後経過時間t=0から微小時間Δt毎に逐次求めていく(更新していく)ことで、下記(7)式に従って混合気温度Tmixを噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができる。以下、先ず、混合気質量Mmixの求め方について説明する。
<<混合気質量Mmix>>
上述したように、質量Qfinの燃料蒸気は、任意の噴射後経過時間tにおいて質量Gの混合気形成筒内ガスと混ざり合って質量(Qfin+G)の混合気となっているから、任意の噴射後経過時間tにおける混合気質量Mmixは(Qfin+G)である。ここで、上記(1)式より「G=stoich・λ・Qfin」と表すことができるから、混合気質量Mmixは、空気過剰率λを用いて下記(8)式にて表すことができる。
よって、上述したように噴射後経過時間t=0から微小時間Δt毎に取得され得る空気過剰率λの値を上記(8)式に順次適用していくことで、混合気質量Mmixを噴射後経過時間t=0から微小時間Δt毎に取得することができる。
<<混合気の定圧比熱Cmix>>
次に、混合気の定圧比熱Cmixの求め方について説明する。一般に、混合気の定圧比熱Cmixは、同混合気内の酸素濃度[O2]mix、及び混合気温度Tmixに大きく依存すると考えられる。ここで、混合気内の酸素濃度[O2]mixは、後述するように、噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができる。従って、混合気温度Tmixを噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができれば、下記(9)式に従って、混合気の定圧比熱Cmixを微小時間Δt毎に求めることができる。
上記(9)式において、funcCmixは、混合気の酸素濃度[O2]mix、及び混合気温度Tmixを引数とする混合気の定圧比熱Cmixを求めるための関数である。なお、上記(9)式を使用して混合気の定圧比熱Cmixを微小時間Δt毎に求めていく際における[O2]mix,Tmixの引数値としては、それぞれ現時点(即ち、噴射後経過時間t)よりも微小時間Δt前の値が使用される。
<<混合気のエンタルピHmix>>
次に、混合気のエンタルピHmixの求め方について説明する。いま、噴射後経過時間(t−Δt)における混合気のエンタルピHmix(t−Δt)が既知である場合において、同噴射後経過時間(t−Δt)から噴射後経過時間tまでの微小時間Δtの間における混合気のエンタルピの増加分ΔHmixについて考える。この混合気のエンタルピの増加分ΔHmixとしては、微小時間Δtの間に混合気に新たに取り込まれる筒内ガスの熱エネルギーΔHgと、同微小時間Δtの間において混合気内で発生する化学反応による反応熱Hrとが挙げられる。
先ず、上記筒内ガスの熱エネルギーΔHgは、下記(10)式により表すことができる。ここにおいて、gは微小時間Δtの間に混合気に新たに取り込まれる筒内ガスの質量である。この質量gは、噴射後経過時間tにおける混合気形成筒内ガスの質量から噴射後経過時間(t−Δt)における混合気形成筒内ガスの質量を減じた値である。従って、上述した関係「G=stoich・λ・Qfin」を利用して下記(11)式により求めることができる。(11)式において、λ(t),λ(t−Δt)はそれぞれ、噴射後経過時間t,(t−Δt)における空気過剰率であり、上記(2)式、(3)式から求めることができる。
また、上記(10)式において、Tgは、噴射後経過時間tにおける筒内ガスの温度であり、筒内ガスの状態がIVC以降断熱変化するとの仮定のもと、下記(12)式に従って求めることができる。下記(12)式において、上述したように、TgivcはIVCにおける筒内ガス温度であり、Vg(CAivc)はIVCにおけるクランク角度CAに対応する筒内容積である。また、Vg(CA)は現時点(即ち、噴射後経過時間t)における上記筒内容積Vg(CA)である。
また、上記(10)式において、Cgは、噴射後経過時間tにおける筒内ガスの定圧比熱であり、混合気の定圧比熱Cmixを求める上記(9)式と同様、下記(13)式に従って求めることができる。下記(13)式において、funcCgは、吸気中の酸素濃度[O2]in、及び筒内ガス温度Tgを引数とする筒内ガスの定圧比熱Cgを求めるための関数である。
なお、上記(13)式を使用して筒内ガスの定圧比熱Cgを微小時間Δt毎に求めていく際における[O2]inの引数値としては、吸気酸素濃度センサ78により検出される上記吸気酸素濃度RO2inが使用される。また、上記筒内ガス温度Tgの引数値としては、現時点(即ち、噴射後経過時間t)における値が使用される。以上により、上記(10)式の右辺の項の全てを求めることができるから、同(10)式に従って上記筒内ガスの熱エネルギーΔHgを求めることができる。
次に、上記微小時間Δtの間において混合気内で発生する化学反応による反応熱Hrは、下記(14)式で表すことができる。下記(14)式において、Hfは所定の定数であり、qrは上記微小時間Δtの間において混合気内で発生する化学反応による燃料消費量である。
上記燃料消費量qrの対象となる化学反応としては、着火反応(熱炎反応)、低温酸化反応(冷炎反応)のみならず、その他の種々の化学反応が含まれる。この燃料消費量qrは、混合気内の酸素濃度[O2]mix、混合気内の燃料濃度[Fuel]mix、混合気温度Tmixに大きく依存すると考えられるから、下記(15)式に従って表すことができる。
上記(15)式において、funcqrは、混合気内の酸素濃度[O2]mix、混合気内の燃料濃度[Fuel]mix、及び混合気温度Tmixを引数とする上記燃料消費量qrを求めるための関数である。混合気内の燃料濃度[Fuel]mixも、上記混合気内の酸素濃度[O2]mixと同様、後述するように、噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができる。上記(15)式を使用して上記燃料消費量qrを微小時間Δt毎に求めていく際における[O2]mix,[Fuel]mixの引数値としては、ぞれぞれ現時点(即ち、噴射後経過時間t)よりも微小時間Δt前の値が使用される。
また、上記(15)式の混合気温度Tmixの引数値としては化学反応前の混合気温度Tpreが使用される。この化学反応前の混合気温度Tpreは、上記(11)式にて算出される質量gの筒内ガスが混合気に新たに取り込まれた後であって、且つ、噴射後経過時間(t−Δt)からの上記微小時間Δtの間における化学反応が発生する前の段階における混合気温度であり、下記(16)式に従って求めることができる。
上記(16)式において、Mmix(t−Δt),Cmix(t−Δt)はそれぞれ、噴射後経過時間(t−Δt)における混合気質量、及び混合気の定圧比熱であって、上記(8)式、及び上記(9)式によりそれぞれ取得され得る。更には、噴射後経過時間(t−Δt)における混合気のエンタルピHmix(t−Δt)は既知である。以上により、上記化学反応前の混合気温度Tpreを求めることができる。従って、上記(15)式の右辺の引数値を全て取得することができるから、(15)式,(14)式に従って上記化学反応による反応熱Hrを求めることができる。
以上のことから、噴射後経過時間(t−Δt)における混合気のエンタルピHmix(t−Δt)が既知である場合において、同噴射後経過時間(t−Δt)から噴射後経過時間tまでの微小時間Δtの間における混合気のエンタルピの増加分ΔHmix(=ΔHg+Hr)が求められるから、噴射後経過時間tにおける混合気のエンタルピHmix(t)(=Hmix(t−Δt)+ΔHmix)を求めることができる。
更には、噴射後経過時間t=0における混合気は、筒内ガスが取り込まれる前の状態(即ち、燃料蒸気のみ)であるから(図3(a)を参照)、この時点での混合気のエンタルピHmix(0)は、下記(17)式にて求めることができる。ここにおいて、Cfは燃料(蒸気)の定圧比熱(ここでは、定数)である。
また、Tfは、燃料蒸気そのものの温度であり、液体の燃料が噴射直後に燃料蒸気に変化する際の単位質量当たりの潜熱Qvaporを考慮して下記(18)式に従って求めることができる。下記(18)式において、Tcrは噴射後経過時間t=0において燃料温度センサ76により検出される液体の燃料温度である。αcrは燃料が燃料噴射用ポンプ22の吐出口近傍から燃料噴射弁21までの燃料配管23を通過する際の熱損失分を考慮するための補正係数である。
従って、噴射後経過時間t=0における混合気のエンタルピHmix(0)も求めることができる。以上より、混合気のエンタルピHmixを、噴射後経過時間t=0から微小時間Δt毎に取得することができる。
このようにして、混合気のエンタルピHmix、混合気質量Mmix、及び混合気の定圧比熱Cmixがそれぞれ噴射後経過時間t=0から微小時間Δt毎に逐次求められるから、上記(7)式に従って混合気温度Tmixを噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができる。
<混合気内の燃料濃度[Fuel]mixの取得>
次に、混合気内の燃料濃度(質量濃度)[Fuel]mixを取得する方法について説明する。噴射後経過時間tにおける混合気内の燃料濃度[Fuel]mixは、上記(8)式により取得される噴射後経過時間tにおける混合気質量Mmixに対する、「噴射後経過時間tにおいて混合気内に存在する燃料の質量」の割合である。
ここで、「噴射後経過時間tにおいて混合気内に存在する燃料の質量」は、噴射後経過時間t=0にて噴射された燃料量(指令燃料噴射量Qfin)から、噴射から現時点(噴射後経過時間t)までの間に化学反応で消費された燃料分を減じた値である。従って、噴射後経過時間tにおける混合気内の燃料濃度[Fuel]mixは、下記(19)式で表すことができる。
上記(19)式において、「Σqr」は、噴射から現時点(噴射後経過時間t)まで微小時間Δt毎に上記(15)式に従って逐次取得・更新されていくそれぞれの燃料消費量qrの和である。このように、燃料消費量qr、及び混合気質量Mmixをそれぞれ噴射後経過時間t=0から微小時間Δt毎に逐次求めることで、上記(19)式に従って混合気内の燃料濃度[Fuel]mixを噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができる。
<混合気内の酸素濃度[O2]mixの取得>
次に、混合気内の酸素濃度(質量濃度)[O2]mixを取得する方法について説明する。噴射後経過時間tにおける混合気内の酸素濃度[O2]mixは、噴射後経過時間tにおける上記混合気質量Mmixに対する、「噴射後経過時間tにおいて混合気内に存在する酸素の質量」の割合である。
ここで、「噴射後経過時間tにおいて混合気内に存在する筒内ガスの質量」は、上述した噴射後経過時間tにおける混合気形成筒内ガスの質量Gから、噴射から現時点(噴射後経過時間t)までの間に化学反応で消費された筒内ガス分を減じた値である。上記微小時間Δtにおける燃料消費量qrの燃料と反応して同微小時間Δtの間において混合気内で化学反応により消費される筒内ガスの消費量grは、下記(20)式にて表すことができる。
従って、「噴射後経過時間tにおいて混合気内に存在する筒内ガスの質量」は、「G−Σgr」と表すことができる。ここで、「Σgr」は、噴射から現時点(噴射後経過時間t)まで微小時間Δt毎に上記(20)式に従って逐次取得・更新されていくそれぞれの筒内ガス消費量grの和である。
「噴射後経過時間tにおいて混合気内に存在する酸素の質量」は、上記「噴射後経過時間tにおいて混合気内に存在する筒内ガスの質量」に、筒内ガス内の酸素濃度(従って、上記吸気中の酸素濃度[O2]in)を乗じることで求めることができる。以上のことから、噴射後経過時間tにおける混合気内の酸素濃度[O2]mixは、下記(21)式で表すことができる。
このように、混合気形成筒内ガスの質量G、筒内ガス消費量gr、及び混合気質量Mmixをそれぞれ噴射後経過時間t=0から微小時間Δt毎に逐次求めることで、上記(21)式に従って混合気内の酸素濃度[O2]mixを噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができる。
<混合気内の窒素濃度[N2]mixの取得>
次に、混合気内の窒素濃度(質量濃度)[N2]mixを取得する方法について説明する。噴射後経過時間tにおける混合気内の窒素濃度[N2]mixは、噴射後経過時間tにおける上記混合気質量Mmixに対する、「噴射後経過時間tにおいて混合気内に存在する窒素の質量」の割合である。
ここで、筒内ガス中の窒素は不活性ガスであって噴射から現時点(噴射後経過時間t)までの間に混合気内で化学反応により消費されることがない。従って、「噴射後経過時間tにおいて混合気内に存在する窒素の質量」は、噴射後経過時間tにおける混合気形成筒内ガスの質量Gに、筒内ガス内の窒素濃度(即ち、吸気中の窒素濃度[N2]in)を乗じることで求めることができる。以上のことから、噴射後経過時間tにおける混合気内の窒素濃度[N2]mixは、下記(22)式で表すことができる。
このように、混合気形成筒内ガスの質量G、混合気質量Mmixをそれぞれ噴射後経過時間t=0から微小時間Δt毎に逐次求めることで、上記(22)式に従って混合気内の窒素濃度[N2]mixを噴射後経過時間t=0から微小時間Δt毎に逐次求めていくことができる。以上、エミッション発生量の推定に必要となる混合気温度Tmix、混合気内の燃料濃度[Fuel]mix、混合気内の酸素濃度[O2]mix、及び混合気内の窒素濃度[N2]mixの取得方法について説明した。
(エミッション発生量の推定方法)
次に、本装置によるエミッション発生量の推定方法について説明する。エミッション発生量としては、本例では、スート(Soot)発生量と、NO発生量とが推定される。以下、先ず、スート発生量の推定について説明する。
<スート発生量>
スート発生量は、例えば、日本機械学会論文集(B編) 48巻432号「直接噴射式ディーゼル機関の燃焼モデルと性能予測」(以下、「非特許文献2」と称呼する。)にて紹介された実験式である下記(23)式から得られる、混合気内のスート濃度[Soot]mixの変化速度(以下、「スート発生速度d[Soot]mix/dt」と称呼する。)を利用して求めることができる。
上記(23)式において、dmsf/dtは混合気内でのスート生成に伴う混合気内のスート濃度[Soot]mixの増加速度(以下、「スートの生成速度」と称呼する。)であり、dmso/dtは混合気内で生成されたスートの酸化に伴う混合気内のスート濃度[Soot]mixの減少速度(以下、「スートの酸化速度」と称呼する。)である。これらはそれぞれ、下記(24)式、(25)式に従って表される。
上記(24)式、(25)式において、Af,Aoは定数、Esf,Esoは活性化エネルギー(ここでは、定数)である。また、Pgは、筒内ガスの圧力であり、筒内ガスの状態がIVC以降断熱変化するとの仮定のもと、下記(26)式に従って求めることができる。下記(26)式において、上述したように、PgivcはIVCにおける筒内ガス圧力であり、Vg(CAivc)はIVCにおけるクランク角度CAに対応する筒内容積である。また、Vg(CA)は現時点(即ち、噴射後経過時間t)における上記筒内容積Vg(CA)である。
上記(24)式は、混合気内の燃料濃度[Fuel]mixが大きいほど、混合気温度Tmixが高いほど、筒内ガス圧力Pgが高いほど、スートの生成速度が大きくなる(従って、スートが生成され易くなる)ことを示している。また、上記(25)式は、混合気内のスート濃度[Soot]mixが大きいほど、混合気内の酸素濃度[O2]mixが大きいほど、混合気温度Tmixが高いほど、筒内ガス圧力Pgが高いほど、スートの酸化速度が大きくなる(従って、スートが酸化して消滅し易くなる)ことを示している。
上記(23)式〜(25)式に、上述したように噴射後経過時間t=0から微小時間Δt毎に逐次取得され得る混合気温度Tmix、混合気内の燃料濃度[Fuel]mix、及び混合気内の酸素濃度[O2]mixを同微小時間Δt毎に順次適用していくことで、スート発生速度d[Soot]mix/dtを噴射後経過時間t=0から微小時間Δt毎に求めていくことができる。
従って、噴射後経過時間t=0から微小時間Δt毎に求められるスート発生速度d[Soot]mix/dtのそれぞれの値を時間で積分(積算)していくことで噴射後経過時間tにおける混合気内のスート濃度[Soot]mixを噴射後経過時間t=0から微小時間Δt毎に更新していくことができる。このようにして、噴射後経過時間tにおける混合気内のスート濃度[Soot]mixを求めることができれば、この混合気内のスート濃度[Soot]mixの値に混合気質量Mmixの値を乗じることで噴射後経過時間tにおけるスート発生量を求めることができる。
<NO発生量>
次に、NO発生量の推定について説明する。NO発生量は、例えば、森北出版「燃焼工学」水谷幸夫著(以下、「非特許文献3」と称呼する。)にて紹介されている拡大ゼルドヴィッチ(Zel’dovich)機構に基づく式である下記(27)式から得られる、混合気内のNO濃度[NO]mixの変化速度(以下、「NO発生速度d[NO]mix/dt」と称呼する。)を利用して求めることができる。
上記(27)式において、Anは定数、Enは活性化エネルギー(ここでは、定数)である。また、Kfo(Tmix)は平衡定数であって混合気温度Tmixに基づいて決定される。上記(27)式は、混合気内の酸素濃度[O2]mixが大きいほど、混合気内の窒素濃度[N2]mixが大きいほど、混合気温度Tmixが高いほど、NOの発生速度が大きくなる(従って、NOが生成され易くなる)ことを示している。
上記(27)式に、上述したように噴射後経過時間t=0から微小時間Δt毎に逐次取得され得る混合気温度Tmix、混合気内の窒素濃度[N2]mix、及び混合気内の酸素濃度[O2]mixを同微小時間Δt毎に順次適用していくことで、NO発生速度d[NO]mix/dtを噴射後経過時間t=0から微小時間Δt毎に求めていくことができる。
従って、噴射後経過時間t=0から微小時間Δt毎に求められるNO発生速度d[NO]mix/dtのそれぞれの値を時間で積分(積算)していくことで噴射後経過時間tにおける混合気内のNO濃度[NO]mixを噴射後経過時間t=0から微小時間Δt毎に更新していくことができる。このようにして、噴射後経過時間tにおける混合気内のNO濃度[NO]mixを求めることができれば、この混合気内のNO濃度[NO]mixの値に混合気質量Mmixの値を乗じることで噴射後経過時間tにおけるNO発生量を求めることができる。
(噴射燃料の区分、混合気の区分)
以上のように取得される混合気の状態(温度Tmix等)、及びエミッション発生量は、図3(a),(b)に示すように、指令燃料噴射量Qfinの燃料が燃料噴射開始時期CAinjにて一時(瞬時)に噴射されるとの仮定のもとで得られる値である。即ち、質量Qfinの燃料の全てが、燃料噴射開始時期CAinj以降、燃料噴射開始時期CAinjでの有効噴射圧力ΔP、同燃料噴射開始時期CAinjでの筒内ガス密度ρg0に基づいて決定される上記噴霧角θ(図3(b)を参照)をもって「ひとつの混合気」となって円錐状に拡散していくものとして扱われる。
しかしながら、上述したように、指令燃料噴射量Qfinの燃料は、実際には、上記燃料噴射開始時期CAinjから噴射期間TAUだけ連続して噴射される。従って、燃焼室内にて進行していく混合気を微視的に見れば、混合気の末尾部分により近い部分に含まれる燃料ほど上記燃料噴射開始時期CAinjから遅れて噴射されることになる。換言すれば、混合気が占める領域内の位置によって、その位置に対応する燃料の部分の噴射時点が異なる。
他方、噴射時点が異なれば、噴射後の経過時間に対する筒内ガスの温度Tg、密度ρg等の状態が異なることに加え、噴射時点での有効噴射圧力ΔP、及び噴射時点での筒内ガス密度ρg0(従って、噴霧角θ)も異なる。換言すれば、燃焼室内にて進行していく混合気の状態は、同混合気に含まれる燃料の噴射時点によって異なる。
以上のことから、噴射後の経過時間に対する混合気の状態(温度等)は、同混合気が占める領域内の位置に応じて異なる。従って、混合気内でのエミッション発生量も同混合気が占める領域内の位置に応じて異なることになる。即ち、燃焼室内にて進行していく混合気の不均一性に起因して混合気内でのエミッション発生量も不均一になる。
係る混合気の不均一性、エミッション発生量の不均一性を考慮して混合気の状態、及びエミッション発生量を精度良く推定するためには、噴射期間TAUを複数の期間に区分して指令燃料噴射量Qfinの燃料を同区分された各期間に対応して順次噴射されていくそれぞれの部分に区分することが考えられる。そして、上記区分された燃料の各部分に基づいて形成されていくそれぞれの混合気の状態を個別に推定することが考えられる。
そこで、本例では、本装置は、図4に示すように、噴射期間TAUを「前期1/3TAU」「中期1/3TAU」「後期1/3TAU」の3つの期間に均等に区分し、前期1/3TAU、中期1/3TAU、後期1/3TAUに対応してそれぞれ質量Q(1),Q(2),Q(3)の燃料が順次個別に噴射されるものとする。
より具体的には、1番目噴射に係わる噴射期間「前期1/3TAU」分の質量Q(1)の燃料が燃料噴射開始時期CAinjにて一時に噴射され、1番目噴射から1/3TAU経過した後に2番目噴射に係わる噴射期間「中期1/3TAU」分の質量Q(2)の燃料が一時に噴射され、2番目噴射から1/3TAU経過した後に3番目噴射に係わる噴射期間「後期1/3TAU」分の質量Q(3)の燃料が一時に噴射されるものとする。なお、「Q(1)+Q(2)+Q(3)=Qfin」の関係があるが、噴射期間TAUに亘る噴射圧力Pcrが一定に維持され得ないこと等に起因してQ(1),Q(2),Q(3)は互いに等しい値とはならない。
本装置は、1番目噴射に基づいて形成されていく混合気(1番目混合気)、2番目噴射に基づいて形成されていく混合気(2番目混合気)、及び3番目噴射に基づいて形成されていく混合気(3番目混合気)の「3つの混合気」をそれぞれ個別に扱い、混合気毎に、上述した手法により混合気の状態(温度Tmix等)及びエミッション発生量を個別に推定する。
これにより、混合気間での、噴射後の経過時間に対する筒内ガスの温度Tg、密度ρg等の状態の相違、噴射時点での有効噴射圧力ΔP及び筒内ガス密度ρg0(従って、噴霧角θ)の相違が考慮されて混合気の状態及びエミッション発生量が個別に推定され得る。
そして、本装置は、混合気毎に推定されたエミッション発生量を合計することでエミッションの総発生量(具体的には、スート総発生量Soot、NO総発生量NO)を推定する。これにより、上述したエミッション発生量の不均一性が考慮されてエミッションの総発生量が精度良く推定され得る。
(燃料噴射制御の概要)
本装置は、上述した混合気の状態の推定、及びエミッション発生量の推定に係わる計算を筒内ガスの量が確定するIVCの直後に開始し、燃料噴射開始時期CAinjの到来前にエミッション総発生量(具体的には、スート総発生量Soot、及びNO総発生量NO)の推定を完了する。
そして、本装置は、目標エミッション発生量(スート目標発生量Sootter、NO目標発生量NOter)を機関の運転状態から求め、上記推定されたスート総発生量Soot、及びNO総発生量NOの一方が対応する目標発生量に対して十分に大きい場合、その一方の総発生量が小さくなるように燃料噴射圧力をフィードバック制御する。
具体的には、本装置は、スート総発生量Sootからスート目標発生量Sootterを減じた値が所定量よりも大きいとき、燃料噴射開始圧力を基本燃料噴射圧力Pcrbaseよりも所定量だけ高くする。これにより、スート発生量が小さくなる方向に制御される。一方、本装置は、NO総発生量NOからNO目標発生量NOterを減じた値が所定量よりも大きいとき、燃料噴射開始圧力を基本燃料噴射圧力Pcrbaseよりも所定量だけ低くする。これにより、NO発生量が小さくなる方向に制御される。以上が、燃料噴射制御の概要である。
(実際の作動)
次に、上記のように構成された内燃機関のエミッション発生量推定装置の実際の作動について説明する。
<混合気温度等及びエミッション量の算出>
CPU61は、図5〜図9に一連のフローチャートにより示した混合気温度等、及びエミッション量の算出を行うためのルーチンを所定時間の経過毎に、気筒毎に、繰り返し実行するようになっている。従って、所定のタイミングになると、CPU61はステップ500から処理を開始し、ステップ505に進んで吸気弁Vinが開状態から閉状態へと変化したか否か(即ち、IVCが到来したか否か)を判定し、「No」と判定する場合、ステップ595に直ちに進んで本ルーチンを一旦終了する。
いま、或る気筒においてIVCが到来したものとすると、CPU61はステップ505に進んだとき「Yes」と判定してステップ510に進み、IVC時クランク角度CAivcをクランクポジションセンサ74から取得される現時点での実クランク角度CAactの値に設定し、IVC時筒内ガス圧力Pgivcを吸気管圧力センサ73から得られる現時点での吸気管圧力Pbの値に設定し、IVC時筒内ガス温度Tgivcを吸気温センサ72から得られる現時点での吸気温度Tbの値に設定し、吸気酸素濃度[O2]inを吸気酸素濃度センサ78から得られる現時点での吸気酸素濃度RO2inの値に設定する。
続いて、CPU61はステップ515に進んで、上記設定されたIVC時筒内ガス圧力Pgivcと、上記設定されたIVC時筒内ガス温度Tgivcと、上記(5)式とに基づいて筒内ガスの全質量Mgを求める。
次いで、CPU61はステップ520に進み、アクセル開度センサ75により得られる現時点でのアクセル開度Accp、クランクポジションセンサ74から取得される現時点でのエンジン回転速度NE、及び図10に示したテーブル(マップ)MapQfinから指令燃料噴射量Qfin(即ち、燃料噴射期間TAU)を求める。テーブルMapQfinは、アクセル開度Accp及びエンジン回転速度NEと指令燃料噴射量Qfinとの関係を規定するテーブルであり、ROM62内に格納されている。
次に、CPU61はステップ525に進み、指令燃料噴射量Qfin、エンジン回転速度NE、及び図11に示したテーブルMapCAinjから燃料噴射開始時期CAinjを決定する。テーブルMapCAinjは、指令燃料噴射量Qfin及びエンジン回転速度NEと燃料噴射開始時期CAinjとの関係を規定するテーブルであり、ROM62内に格納されている。
続いて、CPU61はステップ530に進んで、指令燃料噴射量Qfin、エンジン回転速度NE、及び図12に示したテーブルMapPcrbaseから基本燃料噴射圧力Pcrbaseを決定する。テーブルMapPcrbaseは、指令燃料噴射量Qfin及びエンジン回転速度NEと基本燃料噴射圧力Pcrbaseとの関係を規定するテーブルであり、ROM62内に格納されている。
次に、CPU61はステップ535に進み、上記求めた燃料噴射期間TAUを3で除した値「TAU/3」と、上記求めた基本燃料噴射圧力Pcrbaseと、関数funcPcrとから、1番目噴射、2番目噴射、及び3番目噴射のそれぞれの燃料噴射圧力Pcr(1),Pcr(2),Pcr(3)を求める。上述のごとく、1番目噴射(質量Q(1))は燃料噴射開始時期CAinjにて一時に実行され、2番目噴射(質量Q(2))は1番目噴射から1/3TAU経過した後に一時に実行され、3番目噴射(質量Q(3))は2番目噴射から1/3TAU経過した後に一時に実行されるものとする。
関数funcPcrは、噴射圧力が燃料噴射開始時期CAinjにおいて上記求められた基本燃料噴射圧力Pcrbaseに調整されている状態にて、燃料噴射開始時期CAinjから燃料噴射期間TAUに亘って燃料が連続して噴射される場合における噴射圧力の変動(低下)を考慮して、燃料噴射圧力Pcr(1),Pcr(2),Pcr(3)を決定する関数である。具体的には、これにより、1番目噴射の燃料噴射圧力Pcr(1)は基本燃料噴射圧力Pcrbaseと等しい値に設定され、2,3番目噴射圧力Pcr(2),Pcr(3)は、基本燃料噴射圧力Pcrbaseよりも所定量だけ低い値に設定される。
続いて、CPU61はステップ540に進み、上記求めた燃料噴射期間TAU、値「TAU/3」と、上記求めた基本燃料噴射圧力Pcrbaseと、関数funcQとから、1番目噴射、2番目噴射、及び3番目噴射のそれぞれの燃料量(質量)Q(1),Q(2),Q(3)を求める。このステップ540が噴射燃料区分手段に相当する。
関数funcQは、噴射圧力が燃料噴射開始時期CAinjにおいて上記求められた基本燃料噴射圧力Pcrbaseに調整されている状態にて、燃料噴射開始時期CAinjから燃料噴射期間TAUに亘って燃料が連続して噴射される場合における噴射圧力の変動(低下)を考慮して、燃料量Q(1),Q(2),Q(3)を決定する関数である。これにより、燃料量(質量)Q(1),Q(2),Q(3)は、上述のごとく、「Qfin=Q(1)+Q(2)+Q(3)」の関係が成立するように設定される一方で、互いに等しい値とはならない。
次に、CPU61はステップ545に進んで上記値「TAU/3」と、現時点でのエンジン回転速度NEと、上記決定された燃料噴射開始時期CAinjと、関数funcCAinjとから、1番目噴射、2番目噴射、及び3番目噴射のそれぞれの噴射時点でのクランク角度CAinj(1),CAinj(2),CAinj(3)を求める。
関数funcCAinjは、エンジン回転速度NEが現時点での値で一定であるものとしてクランク角度CAinj(1),CAinj(2),CAinj(3)を求める関数である。具体的には、これにより、1番目噴射時クランク角度CAinj(1)は上記燃料噴射開始時期CAinjに対応する値に設定され、2番目噴射時クランク角度CAinj(2)は燃料噴射開始時期CAinj(従って、1番目噴射時)から期間「TAU/3」だけ遅れた時期に対応する値に設定され、3番目噴射時クランク角度CAinj(3)は2番目噴射時から期間「TAU/3」だけ遅れた時期に対応する値に設定される。
次いで、CPU61はステップ550に進み、燃料温度センサ76から得られる現時点での燃料温度Tcrと、上記(18)式とに基づいて燃料蒸気温度Tfを求める。続いて、CPU61はステップ555に進んで、上記吸気酸素濃度[O2]inと、[O2]inを引数とする筒内ガス理論空燃比stoichを求めるための関数funcstoichとに基づいて筒内ガス理論空燃比stoichを求める。
次に、CPU61はステップ560に進み、上記吸気酸素濃度[O2]inと、[O2]inを引数とする筒内ガス中の窒素濃度を求めるための関数func[N2]とに基づいて筒内ガス中の窒素濃度(即ち、吸気窒素濃度[N2]in)を求める。
続いて、CPU61はステップ565に進んで、現時点でのエンジン回転速度NEと、微小時間Δt(例えば、0.1msec)と、NE,Δtを引数とする微小クランク角度ΔCAを求めるための関数funcΔCAとに基づいて、同微小時間Δtに相当するクランク角度である微小クランク角度ΔCAを求める。この微小クランク角度ΔCAは、エンジン回転速度NEが現時点(即ち、IVC直後)での値である場合における、微小時間Δtに相当するクランク角度である。
そして、CPU61はステップ570に進み、スート総発生量Sootを初期値Soot0に設定するとともにNO総発生量NOを初期値NO0に設定し、続くステップ575にて変数iの値を「0」にクリアする。ここで、初期値Soot0は、燃料噴射前の時点から筒内ガス中に含まれている(従って、EGRガス中に予め含まれている)スート量に相当し、初期値NO0は、燃料噴射前の時点から筒内ガス中に含まれている(従って、EGRガス中に予め含まれている)NO量に相当する。
また、変数iの値は、何番目の噴射かを識別するための値であり、「i=1」は1番目噴射(従って、1番目混合気)を表し、「i=2」は2番目噴射(従って、2番目混合気)を表し、「i=3」は3番目噴射(従って、3番目混合気)を表す。
次に、CPU61は図6のステップ605へと進み、変数iの値をその時点での値(現時点では「0」)に「1」を加えた値に設定する。これにより、変数iは現時点では「1」になる。従って、以降、「i=1」である限りにおいて1番目噴射に係わる算出がなされていく。具体的には、先ず、以下のようにして、1番目噴射に関する各種初期値が決定されていく。
先ず、CPU61はステップ610に進んで、先のステップ515にて求めた筒内ガスの全質量Mgを、先のステップ545にて求めたi番目(1番目)噴射時クランク角度CAinj(i)から得られる1番目噴射時筒内容積Vg(CAinj(1))で除することで(1番目噴射時の)筒内ガス密度ρg0を求める。
続いて、CPU61はステップ615に進み、先のステップ510にて求めたIVC時筒内ガス圧力Pgivcと、上記IVC時筒内容積Vg(CAivc)と、上記i番目(1番目)噴射時筒内容積Vg(CAinj(i))と、上記(4)式に相当する式とに基づいて(1番目噴射時の)筒内ガス圧力Pg0を求める。
次に、CPU61はステップ620に進んで、先のステップ535にて求めたi番目(1番目)噴射圧力Pcr(i)から上記筒内ガス圧力Pg0を減じることで(1番目噴射時の)有効噴射圧力ΔPを求め、続くステップ625にて、上記求めた有効噴射圧力ΔPと、筒内ガス密度ρg0と、テーブルMapθとに基づいて(1番目噴射に係わる)噴霧角θを求める。これにより、噴霧角θは、i番目(1番目)噴射時点(即ち、クランク角度CAinj(1))での有効噴射圧力ΔP及び筒内ガス密度ρg0に基づいて決定される。
次いで、CPU61はステップ630に進み、空気過剰率前回値λbを初期値「0」に設定し、続くステップ635にて(1番目混合気に係わる)混合気形成筒内ガス質量Gの値を初期値「0」に設定するとともに、続くステップ640にて(1番目混合気に係わる)燃料消費量積算値sumqr、及び筒内ガス消費量積算値sumgrを共に初期値「0」に設定する。
続いて、CPU61はステップ645に進んで、(1番目混合気に係わる)混合気のエンタルピHmixを、上記(17)式に相当する式に従って、初期値(即ち、先のステップ540にて求めたi番目(1番目)噴射量Q(i)と、燃料の定圧比熱Cfと、先のステップ550にて求めた燃料蒸気温度Tfの積)に設定する。
次に、CPU61はステップ650に進み、(1番目混合気に係わる)混合気の定圧比熱Cmixを初期値である上記燃料の定圧比熱Cfに設定し、続くステップ655にて(1番目混合気に係わる)混合気質量Mmixを初期値であるi番目(1番目)噴射量Q(i)に設定する。
次いで、CPU61はステップ660に進んで、(1番目混合気に係わる)混合気内のNO濃度[NO]mix、スート濃度[Soot]mix、酸素濃度[O2]mixをそれぞれ初期値「0」に設定するとともに、(1番目混合気に係わる)混合気内の燃料濃度[Fuel]mixを初期値「1」に設定する。
続いて、CPU61はステップ665に進み、(1番目混合気に係わる)噴射後経過時間tを初期値「0」に設定するとともに、(1番目噴射に係わる)クランク角度CAを初期値であるi番目(1番目)噴射時クランク角度CAinj(i)に設定する。これにより、i番目(1番目)混合気に係わる噴射後経過時間tがi番目(1番目)噴射時からカウントされることになる。
次に、CPU61はステップステップ670に進んで、フラグENDsootの値、フラグENDnoの値を共に「0」に初期化する。ここで、フラグENDsootは、その値が「1」のときスート濃度[Soot]mixの値が更新中であることを示し、その値が「0」のときスート濃度[Soot]mixの更新が終了していることを示す。また、フラグENDnoは、その値が「1」のときNO濃度[NO]mixの値が更新中であることを示し、その値が「0」のときNO濃度[NO]mixの更新が終了していることを示す。このようにして、i番目(1番目)噴射に関する各種初期値が決定される。
次に、CPU61は図7のルーチンに進み、i番目(1番目)噴射に関する混合気温度の算出のための処理を開始する。具体的には、CPU61は先ずステップ705に進み、(1番目噴射に係わる)噴射後経過時間tの値(現時点では、「0」)を微小時間Δtだけ増大・更新させるとともに、(1番目噴射に係わる)クランク角度CAの値(現時点では、「CAinj(1)」)を微小クランク角度ΔCAだけ増大・更新させる。このように、クランク角度CAの値が噴射後経過時間tに対応する値に維持されていく。これにより、以降、(1番目噴射に係わる)噴射後経過時間t=Δtとなり、(1番目噴射に係わる)クランク角度CA=CAinj(1)+ΔCAとなる。
続いて、CPU61はステップ710に進んで、先のステップ515にて求めた筒内ガスの全質量Mgを、上記ステップ705にて更新されたクランク角度CAに対応する筒内容積Vg(CA)で除することで噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)筒内ガス密度ρgを求める。
次に、CPU61はステップ715に進んで、上記ステップ510にて求められたIVC時筒内ガス圧力Pgivcと、上記IVC時筒内容積Vg(CAivc)と、上記クランク角度CAに対応する筒内容積Vg(CA)と、上記(26)式とに基づいて噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)筒内ガス圧力Pgを求める。
次いで、CPU61はステップ720に進み、上記ステップ510にて求められたIVC時筒内ガス温度Tgivcと、上記IVC時筒内容積Vg(CAivc)と、上記クランク角度CAに対応する筒内容積Vg(CA)と、上記(12)式とに基づいて噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)筒内ガス温度Tgを求める。
次に、CPU61はステップ725に進んで、上記ステップ510にて求められた上記吸気酸素濃度[O2]inと、ステップ720にて求めた筒内ガス温度Tgと、[O2]in,Tgを引数とする筒内ガスの定圧比熱Cgを求めるための関数funcCgと、上記(13)式とに基づいて噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)筒内ガスの定圧比熱Cgを求める。
続いて、CPU61はステップ730に進み、先のステップ710にて求めた筒内ガス密度ρgと、先のステップ625にて求めた噴霧角θと、先のステップ620にて求めた有効噴射圧力ΔPと、先のステップ705にて更新したi番目(1番目)噴射に係わる噴射後経過時間tと、上記(3)式とに基づいてi番目(1番目)噴射に係わる燃料希釈率dλ/dtを求める。
次いで、CPU61はステップ735に進んで、上記(2)式に従って、i番目(1番目)噴射に係わる空気過剰率λを、その時点での空気過剰率前回値λb(現時点では、ステップ630の処理により「0」)に、上記求めた燃料希釈率dλ/dtに微小時間Δtを乗じた値「dλ/dt・Δt」を加えた値に更新する。これにより、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)空気過剰率λが求められる。
次に、CPU61はステップ740に進み、先のステップ555にて求めた筒内ガス理論空燃比stoichと、ステップ735にて求めた空気過剰率λと、空気過剰率前回値λb(現時点では、ステップ630の処理により「0」。次回からは、後述するステップ785にて設定されている値)と、ステップ540にて設定されたi番目(1番目)噴射量Q(i)と、上記(11)式に相当する式とに基づいて(1番目混合気に係わる)微小時間Δt(噴射後経過時間(t−Δt)〜tの間)において混合気に新たに取り込まれた筒内ガス質量gを求める。
続いて、CPU61はステップ745に進んで、混合気形成筒内ガス質量Gを、その時点での値(現時点では、ステップ635の処理により「0」)に上記新たに取り込まれた筒内ガス質量gを加えた値に更新する。これにより、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)混合気形成筒内ガス質量Gが求められる。
次に、CPU61はステップ750に進み、混合気質量Mmixを、その時点での値(現時点では、ステップ655の処理により1番目噴射量Q(1))に上記新たに取り込まれた筒内ガス質量gを加えた値に更新する。これにより、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)混合気質量Mmixが求められる。
続いて、CPU61はステップ755に進んで、化学反応前の混合気のエンタルピHpreを、その時点での混合気のエンタルピHmix(現時点では、ステップ645の処理により値「Q(1)・Cf・Tf」)に、上記(10)式に従って求められる「上記新たに取り込まれた筒内ガスの熱エネルギーΔHg=g・Cg・Tg」を加えた値に設定する。
次に、CPU61はステップ760に進み、上記求めた化学反応前の混合気のエンタルピHpreを、上記ステップ750にて求めた混合気質量Mmixにその時点での混合気の定圧比熱Cmix(現時点では、ステップ650の処理により燃料の定圧比熱Cf。次回からは、後述するステップ830にて設定されている値)を乗じた値で除することにより上記(16)式に相当する式に従って化学反応前の混合気の温度Tpreを求める。
次いで、CPU61はステップ765に進んで、その時点での混合気内の酸素濃度[O2]min(現時点では、ステップ660の処理により「0」。次回からは、後述するステップ820にて設定されている値)と、燃料濃度[Fuel]mix(現時点では、ステップ660の処理により「0」。次回からは、後述するステップ815にて設定されている値)と、上記求めた化学反応前の混合気の温度Tpreと、上記(15)式とに基づいて微小時間Δt(噴射後経過時間(t−Δt)〜tの間)において混合気内で発生する化学反応による燃料消費量qrを求める。
続いて、CPU61はステップ770に進み、上記求めた消費燃料量qrと、上記(14)式とに基づいて微小時間Δt(噴射後経過時間(t−Δt)〜tの間)において混合気内で発生する化学反応による反応熱Hrを求め、続くステップ775にて混合気のエンタルピHmixを、上記求めた化学反応前の混合気のエンタルピHpreに上記求めた反応熱Hrを加えた値に設定・更新する。これにより、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)混合気のエンタルピHmixが求められる。
そして、CPU61はステップ780に進んで、上記ステップ775にて求めた混合気のエンタルピHmixと、上記ステップ750にて求めた混合気質量Mmixと、その時点での混合気の定圧比熱Cmix(現時点では、ステップ650の処理により燃料の定圧比熱Cf。次回からは、後述するステップ830にて設定されている値)と、上記(7)式とに基づいて混合気温度Tmixを算出する。これにより、噴射後経過時間t=Δt(従って、クランク角度CA=CAinj(1)+ΔCA)における(1番目噴射に係わる)混合気温度Tmixが求められる。このステップ780は混合気状態推定手段に相当する。
次に、CPU61はステップ785に進んで、次回の計算の準備のため、空気過剰率前回値λbを上記ステップ735にて求めた空気過剰率λの値に設定する。この値は、以降、上述したステップ735にて使用される。このようにして、噴射後経過時間tにおけるi番目(1番目)噴射に関する混合気温度Tmixが算出される。
次に、CPU61は図8のルーチンに進み、i番目(1番目)噴射に関する各種濃度の算出のための処理を開始する。具体的には、CPU61は先ずステップ805に進み、上記ステップ765にて求めた微小時間Δt(噴射後経過時間(t−Δt)〜tの間)において混合気内で発生する化学反応による燃料消費量qrと、ステップ555にて求めた筒内ガス理論空燃比stoichと、上記(20)式とに基づいて、微小時間Δt(噴射後経過時間(t−Δt)〜tの間)において混合気内で発生する化学反応による筒内ガス消費量grを求める。
続いて、CPU61はステップ810に進んで、(1番目混合気に係わる)燃料消費量積算値sumqrを、その時点での値(現時点では、ステップ640の処理により「0」)にステップ675にて求めた上記燃料消費量qrを加えた値に設定・更新し、(1番目混合気に係わる)筒内ガス消費量積算値sumgrを、その時点での値(現時点では、ステップ640の処理により「0」)にステップ805にて求めた上記筒内ガス消費量grを加えた値に設定・更新する。これにより、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)燃料消費量積算値sumqr、及び筒内ガス消費量積算値sumgrが求められる。
次いで、CPU61はステップ815に進み、ステップ540にて求めたi番目(1番目)噴射量Q(i)と、上記求めた燃料消費量積算値sumqrと、ステップ750にて求めた混合気質量Mmixと、上記(19)式に相当する式とに基づいて、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)混合気内の燃料濃度[Fuel]mixを求める。
次に、CPU61はステップ820に進んで、ステップ745にて求めた混合気形成筒内ガス質量Gと、上記求めた筒内ガス消費量積算値sumgrと、ステップ510にて設定された吸気酸素濃度[O2]inと、ステップ750にて求めた混合気質量Mmixと、上記(21)式とに基づいて、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)混合気内の酸素濃度[O2]mixを求める。
続いて、CPU61はステップ825に進み、ステップ745にて求めた混合気形成筒内ガス質量Gと、ステップ560にて設定された吸気窒素濃度[N2]inと、ステップ750にて求めた混合気質量Mmixと、上記(22)式とに基づいて、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)混合気内の窒素濃度[N2]mixを求める。このステップ815〜825も混合気状態推定手段に相当する。
次に、CPU61はステップ830に進んで、上記ステップ820にて求めた混合気内の酸素濃度[O2]mixと、ステップ780にて求めた混合気温度Tmixと、上記(9)式とに基づいて噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)混合気の定圧比熱Cmixを求める。この値は、以降、ステップ760、780にて使用される。
次いで、CPU61はステップ835に進み、フラグENDsootの値が「0」であるか否かを判定し、「No」と判定する場合、後述する混合気内のスート濃度[Soot]mixの更新(スート発生速度d[soot]mix/dtの時間積分処理。ステップ855)を行うことなくステップ870に直ちに進む。
現時点では、先のステップ670の処理によりフラグENDsootの値は「0」になっている。従って、CPU61はステップ835にて「Yes」と判定してステップ840に進み、ステップ815にて求めた混合気内の燃料濃度[Fuel]mixと、ステップ715にて求めた筒内ガス圧力Pgと、ステップ780にて求めた混合気温度Tmixと、上記(24)式とに基づいてスートの生成速度dmsf/dtを求める。
続いて、CPU61はステップ845に進んで、その時点での混合気内のスート濃度[Soot]mix(現時点では、ステップ660の処理により「0」。次回からは、後述するステップ855にて設定されている値)と、ステップ820にて求めた混合気内の酸素濃度[O2]mixと、ステップ715にて求めた筒内ガス圧力Pgと、ステップ780にて求めた混合気温度Tmixと、上記(25)式とに基づいてスートの酸化速度dmso/dtを求める。
次に、CPU61はステップ850に進み、上記求めたスートの生成速度dmsf/dtと、上記求めたスートの酸化速度dmso/dtと、上記(23)式とに基づいてスート発生速度d[soot]mix/dtを求め、続くステップ855にて、混合気内のスート濃度[Soot]mixを、その時点での値(現時点では、ステップ660の処理により「0」)に上記求めたスート発生速度d[soot]mix/dtに微小時間Δtを乗じた値「d[Soot]mix/dt・Δt」を加えた値に設定・更新する。これにより、噴射後経過時間t=Δt(従って、クランク角度CA=CAinj(1)+ΔCA)における(1番目噴射に係わる)混合気内のスート濃度[Soot]mixが求められる。このステップ855がエミッション発生量推定手段に相当する。
続いて、CPU61はステップ860に進んで、クランク角度CAが圧縮上死点(以下、「TDC」と称呼する。)以降であって、且つ、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)ステップ780で求められた混合気温度Tmixがスートの反応限界温度TminSootより低いか否かを判定する。なお、スートの反応限界温度TminSootとは、混合気温度Tmixがその温度未満になるとスートが殆ど発生しなくなる温度である。
ここで、「Yes」と判定する場合、CPU61はステップ865に進んでフラグENDsootの値を「0」から「1」に変更する。これにより、以降、CPU61はステップ835に進んだとき「No」と判定するようになり、上述したように、混合気内のスート濃度[Soot]mixの更新(ステップ855)が実行されなくなる。これにより、スート濃度[Soot]mixの更新に係わる不必要な計算が省略でき、CPU61の計算負荷を低減することができる。
現時点はIVCの直後であるから、TDC前である。従って、CPU61はステップ860にて「No」と判定してステップ870に直ちに進む。即ち、フラグENDsootの値が「0」に維持される。
CPU61はステップ870に進むと、フラグENDnoの値が「0」であるか否かを判定する。ここで、「No」と判定する場合、CPU61は後述する混合気内のNO濃度[NO]mixの更新(NO発生速度d[NO]mix/dtの時間積分処理。ステップ880)を行うことなく図9のステップ905に直ちに進む。
現時点では、先のステップ670の処理によりフラグENDnoの値は「0」になっている。従って、CPU61はステップ870にて「Yes」と判定してステップ875に進み、ステップ825にて求めた混合気内の窒素濃度[N2]mixと、ステップ820にて求めた混合気内の酸素濃度[O2]mixと、ステップ780にて求めた混合気温度Tmixと、上記(27)式とに基づいてNO発生速度d[NO]mix/dtを求める。
続いて、CPU61はステップ880に進み、混合気内のNO濃度[Soot]mixを、その時点での値(現時点では、ステップ660の処理により「0」)に上記求めたNO発生速度d[NO]mix/dtに微小時間Δtを乗じた値「d[NO]mix/dt・Δt」を加えた値に設定・更新する。これにより、噴射後経過時間t=Δt(従って、クランク角度CA=CAinj(1)+ΔCA)における(1番目噴射に係わる)混合気内のNO濃度[NO]mixが求められる。このステップ880もエミッション発生量推定手段に相当する。
続いて、CPU61はステップ885に進んで、クランク角度CAがTDC以降であって、且つ、噴射後経過時間t(従って、クランク角度CA)における(1番目噴射に係わる)ステップ780で求められた混合気温度TmixがNOの反応限界温度TminNOより低いか否かを判定する。なお、NOの反応限界温度TminNOとは、混合気温度Tmixがその温度未満になるとNOが殆ど発生しなくなる温度である。
ここで、「Yes」と判定する場合、CPU61はステップ890に進んでフラグENDnoの値を「0」から「1」に変更する。これにより、以降、CPU61はステップ870に進んだとき「No」と判定するようになり、上述したように、混合気内のNO濃度[NO]mixの更新(ステップ880)が実行されなくなる。これにより、NO濃度[NO]mixの更新に係わる不必要な計算が省略でき、CPU61の計算負荷を低減することができる。
上述と同様、現時点はIVCの直後であるから、TDC前である。従って、現時点では、CPU61はステップ885にて「No」と判定して図9のルーチンに直ちに進む。即ち、フラグENDnoの値が「0」に維持される。以上により、変数i=1についての(従って、1番目噴射に係わる)噴射後経過時間t=Δt(クランク角度CA=CAinj(1)+ΔCA)における、混合気(1番目混合気)の空気過剰率λ、1番目混合気の状態(温度Tmix等)、及び1番目混合気に係わるエミッション濃度([Soot]mix,[NO]mix)が算出される。
CPU61は図9のステップ905に進むと、フラグENDsootの値、及びフラグENDnoの値が共に「1」であるか、又は、クランク角度CAがTDC以降の所定の終了判定クランク角度CAendに一致したか否かを判定する。現時点では、上述のごとく、フラグENDsootの値、及びフラグENDnoの値が共に「0」であり、且つ、クランク角度CAは、1番目噴射時クランク角度CAinj(1)に微小クランク角度ΔCAを加えた値(従って、TDC前)であるから上記終了判定クランク角度CAendに達していない。
従って、現時点では、CPU61はステップ905にて「No」と判定して図7のステップ705に戻り、(1番目噴射に係わる)噴射後経過時間t(現時点では、「1・Δt」)を微小時間Δtだけ増大・更新させるとともに、(1番目噴射に係わる)クランク角度CA(現時点では、「CAinj(1)+ΔCA」)を微小クランク角度ΔCAだけ増大・更新させた後、上述した図7のステップ710〜図9のステップ905の処理を再び実行する。
これにより、変数i=1についての(従って、1番目噴射に係わる)噴射後経過時間t=2・Δt(クランク角度CA=CAinj(1)+2・ΔCA)における、混合気(1番目混合気)の空気過剰率λ(ステップ735を参照。)、1番目混合気の状態(温度Tmix等。ステップ780を参照。)、及び(フラグENDsoot=フラグENDno=0である限りにおいて)1番目混合気に係わるエミッション濃度([Soot]mix,[NO]mix。ステップ855、880を参照。)が算出される。
そして、図9のステップ905の判定において「No」と判定される毎に、図7のステップ705〜図9のステップ905の処理が繰り返し実行されていく。これにより、図9のステップ905の判定において「No」と判定される限りにおいて、変数i=1に維持された状態で、1番目混合気の空気過剰率λ、1番目混合気の状態(温度Tmix等)、及び(フラグENDsoot=フラグENDno=0である限りにおいて)1番目混合気に係わるエミッション濃度([Soot]mix,[NO]mix)が1番目噴射時から微小時間Δt毎に(即ち、CAinj(1)から微小クランク角度ΔCA毎に)求められていく。
以降、TDC後の膨張行程における筒内容積の増大等に伴って混合気温度Tmixが低下することで上述したステップ860の条件が成立した場合、フラグENDsootの値が「0」から「1」に変更されて、それ以降、上述のごとく、スート濃度[Soot]mixの更新(ステップ855)が実行されなくなる。また、膨張行程における筒内容積の増大等に伴って上述したステップ885の条件が成立した場合、フラグENDnoの値が「0」から「1」に変更されて、それ以降、上述のごとく、NO濃度[NO]mixの更新(ステップ880)が実行されなくなる。
そして、フラグENDsootの値とフラグENDnoの値が共に「1」となった場合、又は、(フラグENDsootの値とフラグENDnoの値が共に「1」となっていなくても)クランク角度CAが上記終了判定クランク角度CAendに達した場合、CPU61は図9のステップ905に進んだとき「Yes」と判定してステップ910以降に進み、変数i=1についての(従って、1番目噴射に係わる)計算を終了するための処理を行う。
即ち、CPU61はステップ910に進むと、先のステップ855の処理にて更新されている現時点での混合気内のスート濃度[Soot]mixの値に、先のステップ750にて更新されている現時点での混合気質量Mmixを乗じることでi番目(1番目)混合気内のスート発生量Soot(i)を求めるとともに、先のステップ880の処理にて更新されている現時点での混合気内のNO濃度[NO]mixの値に、同混合気質量Mmixを乗じることでi番目(1番目)混合気内のNO発生量NO(i)を求める。
次いで、CPU61はステップ915に進み、スート総発生量Sootを、その時点での値(現時点では、ステップ570の処理により初期値Soot0)に上記求めたi番目(1番目)混合気内のスート発生量Soot(i)を加えた値(Soot0+Soot(1))に更新するとともに、NO総発生量NOを、その時点での値(現時点では、ステップ570の処理により初期値NO0)に上記求めたi番目(1番目)混合気内のNO発生量NO(i)を加えた値(NO0+NO(1))に更新する。以上により、1番目噴射(従って、1番目混合気)に関する計算が終了する。
そして、CPU61はステップ920に進んで変数iの値が「3」であるか否かを判定する。現時点では、変数i=1であるから、CPU61はステップ920にて「No」と判定して図6のステップ605に戻り、変数iの値を「1」だけ増大させる。これにより、以降、変数iが「2」に設定されるから、2番目噴射(従って、2番目混合気)に係わる算出がなされていく。
即ち、先ず、上述した図6のステップ610〜670にて、2番目噴射に関する各種初期値が設定される。具体的には、ステップ610にて、先のステップ545にて求めた2番目噴射時クランク角度CAinj(2)から得られる2番目噴射時筒内容積Vg(CAinj(2))を使用して(2番目噴射時の)筒内ガス密度ρg0が求められる。
続いて、ステップ615にて、上記2番目噴射時筒内容積Vg(CAinj(2))を使用して(2番目噴射時の)筒内ガス圧力Pg0が求められる。この結果、ステップ620にて、2番目噴射圧力Pcr(2)から上記筒内ガス圧力Pg0を減じることで(2番目噴射時の)有効噴射圧力ΔPが求められ、続くステップ625にて、上記求めた有効噴射圧力ΔPと、筒内ガス密度ρg0と、テーブルMapθとに基づいて(2番目噴射に係わる)噴霧角θが求められる。これにより、噴霧角θは、2番目噴射時点(即ち、クランク角度CAinj(2))での有効噴射圧力ΔP及び筒内ガス密度ρg0に基づいて決定される。
また、ステップ645にて、(2番目混合気に係わる)混合気のエンタルピHmixが、初期値(即ち、先のステップ540にて求めた2番目噴射量Q(2)と、燃料の定圧比熱Cfと、先のステップ550にて求めた燃料蒸気温度Tfの積)に設定され、ステップ655にて、(2番目混合気に係わる)混合気質量Mmixが初期値である2番目噴射量Q(2)に設定される。
そして、ステップ665にて、(2番目混合気に係わる)噴射後経過時間tが初期値「0」に設定されるとともに、(2番目噴射に係わる)クランク角度CAが初期値である2番目噴射時クランク角度CAinj(2)に設定される。これにより、2番目混合気に係わる噴射後経過時間tが2番目噴射時からカウントされることになる。
このようにして、図6のステップ610〜670にて2番目噴射に関する各種初期値が設定されると、次に、上述した図7のステップ705〜図9のステップ905の処理が実行される。これにより、図9のステップ905の判定において「No」と判定される限りにおいて、変数i=2に維持された状態で、2番目混合気の空気過剰率λ、2番目混合気の状態(温度Tmix等)、及び(フラグENDsoot=フラグENDno=0である限りにおいて)2番目混合気に係わるエミッション濃度([Soot]mix,[NO]mix)が2番目噴射時から微小時間Δt毎に(即ち、CAinj(2)から微小クランク角度ΔCA毎)に求められていく。
そして、上述したステップ905の条件が成立すると、ステップ910にて2番目混合気内のスート発生量Soot(2)とNO発生量NO(2)が求められ、ステップ915にて、スート総発生量Sootが、現時点での値(初期値Soot0+Soot(1))に上記求めた2番目混合気内のスート発生量Soot(2)を加えた値(Soot0+Soot(1)+Soot(2))に更新されるとともに、NO総発生量NOが、現時点での値(初期値NO0+NO(1))に上記求めた2番目混合気内のNO発生量NO(2)を加えた値(NO0+NO(1)+NO(2))に更新される。以上により、2番目噴射(従って、2番目混合気)に関する計算が終了する。
そして、ステップ920にて変数iの値が「3」であるか否かが判定される。現時点では、変数i=2であるから、「No」と判定されて図6のステップ605に再び戻り、変数iの値が「1」だけ増大させられる。これにより、以降、変数iが「3」に設定されるから、3番目噴射(従って、3番目混合気)に係わる算出がなされていく。
即ち、先ず、上述した図6のステップ610〜670にて、3番目噴射に関する各種初期値が設定される。具体的には、ステップ610〜625にて(3番目噴射に係わる)噴霧角θが、3番目噴射時点(即ち、クランク角度CAinj(3))での有効噴射圧力ΔP及び筒内ガス密度ρg0に基づいて決定される。
また、ステップ645にて、(3番目混合気に係わる)混合気のエンタルピHmixが、初期値(即ち、先のステップ540にて求めた3番目噴射量Q(3)と、燃料の定圧比熱Cfと、先のステップ550にて求めた燃料蒸気温度Tfの積)に設定され、ステップ655にて、(3番目混合気に係わる)混合気質量Mmixが初期値である3番目噴射量Q(3)に設定される。
そして、ステップ665にて、(3番目混合気に係わる)噴射後経過時間tが初期値「0」に設定されるとともに、(3番目噴射に係わる)クランク角度CAが初期値である3番目)噴射時クランク角度CAinj(3)に設定される。これにより、3番目混合気に係わる噴射後経過時間tが3番目噴射時からカウントされることになる。
このようにして、図6のステップ610〜670にて3番目噴射に関する各種初期値が設定されると、次に、上述した図7のステップ705〜図9のステップ905の処理が実行される。これにより、図9のステップ905の判定において「No」と判定される限りにおいて、変数i=3に維持された状態で、3番目混合気の空気過剰率λ、3番目混合気の状態(温度Tmix等)、及び(フラグENDsoot=フラグENDno=0である限りにおいて)3番目混合気に係わるエミッション濃度([Soot]mix,[NO]mix)が3番目噴射時から微小時間Δt毎に(即ち、CAinj(3)から微小クランク角度ΔCA毎)に求められていく。
そして、上述したステップ905の条件が成立すると、ステップ910にて3番目混合気内のスート発生量Soot(3)とNO発生量NO(3)が求められ、ステップ915にて、スート総発生量Sootが、現時点での値(初期値Soot0+Soot(1)+Soot(2))に上記求めた3番目混合気内のスート発生量Soot(3)を加えた値(Soot0+Soot(1)+Soot(2)+Soot(3))に更新されるとともに、NO総発生量NOが、現時点での値(初期値NO0+NO(1)+NO(2))に上記求めた3番目混合気内のNO発生量NO(3)を加えた値(NO0+NO(1)+NO(2)+NO(3))に更新される。
以上により、3番目噴射(従って、3番目混合気)に関する計算が終了する。また、これにより、スート総発生量Sootが各混合気内のスート発生量の和として求められるとともに、NO総発生量NOが各混合気内のNO発生量の和として求められる。即ち、エミッション発生量の推定が終了する。CPU61は、以上説明した混合気状態の推定、及びエミッション発生量の推定をIVCの直後(即ち、燃料噴射開始時期CAinjの前の時点)において完了する。
そして、ステップ920にて変数iの値が「3」であるか否かが判定される。現時点では、変数i=3であるから、「Yes」と判定される。この場合、CPU61はステップ925に進み、現時点(即ち、IVCの直後の時点)におけるエンジン回転速度NEと、上記ステップ520にて決定されている指令燃料噴射量Qfinと、テーブルMapSootter,MapNOterとに基づいてスート目標発生量Sootter、及びNO目標発生量NOterを求める。
次いで、CPU61はステップ930に進んで、スート発生量偏差Δsootを、ステップ915にて求められているスート総発生量Sootから上記求めたスート目標発生量Sootterを減じた値に設定するとともに、NO発生量偏差ΔNOを、ステップ915にて求められているNO総発生量NOから上記求めたNO目標発生量NOterを減じた値に設定する。
続いて、CPU61はステップ935に進み、上記求めたスート発生量偏差Δsootが基準値Sootrefより大きいか否かを判定し、「Yes」と判定する場合、ステップ940に進んで最終燃料噴射圧力Pcrfinを、上記ステップ530にて決定されている基本燃料噴射圧力Pcrbaseに所定値ΔPcrを加えた値に設定する。これにより、スート発生量が小さくなる方向に燃料噴射圧力が補正されることになる。
一方、ステップ935にて「No」と判定する場合、CPU61はステップ945に進み、上記求めたNO発生量偏差ΔNOが基準値NOrefより大きいか否かを判定し、「Yes」と判定する場合、ステップ950に進んで最終燃料噴射圧力Pcrfinを、上記ステップ530にて決定されている基本燃料噴射圧力Pcrbaseから上記所定値ΔPcrを減じた値に設定する。これにより、NO発生量が小さくなる方向に燃料噴射圧力が補正されることになる。
他方、ステップ945にて「No」と判定する場合(即ち、ΔSoot≦Sootref、且つ、ΔNO≦NOrefの場合)、CPU61はステップ955に進んで、最終燃料噴射圧力Pcrfinを上記ステップ530にて決定されている基本燃料噴射圧力Pcrbaseと等しい値に設定する。即ち、この場合、燃料噴射圧力が補正されない。
そして、CPU61はステップ960に進むと、燃料噴射圧力が上記設定された最終燃料噴射圧力Pcrfinになるように燃料噴射用ポンプ22(の駆動回路)に対して制御指示を行い、ステップ595に進んで図5〜図9の一連の本ルーチンを一旦終了する。以降、CPU61は、次のIVCが到来するまでの間、ステップ505に進む毎に「No」と判定し続ける。
この結果、本ルーチンの実行により、IVCが到来する毎に、燃料噴射形態(噴射量、噴射圧力、噴射時期)が決定されるとともに混合気の状態、及びエミッション総発生量が直ちに推定され、係るエミッション総発生量の推定結果に基づいて噴射圧力が補正されていく。
また、CPU61は、図13にフローチャートにより示した燃料噴射制御を行うためのルーチンを所定時間の経過毎に、気筒毎に、繰り返し実行するようになっている。従って、所定のタイミングになると、CPU61はステップ1300から処理を開始し、ステップ1305に進んで実クランク角度CAactが先のステップ525にて決定されている燃料噴射開始時期CAinjに一致したか否かを判定し、「No」と判定する場合、ステップ1395に直ちに進んで本ルーチンを一旦終了する。
いま、実クランク角度CAactが上記燃料噴射開始時期CAinjに一致したものとすると、CPU61はステップ1310に進んで、対応する燃料噴射弁21に対してステップ520にて決定されている指令燃料噴射量Qfinの燃料の噴射指示(具体的には、燃料噴射期間TAUに亘る開弁指示)を行い、ステップ1395に進んで本ルーチンを一旦終了する。これにより、指令燃料噴射量Qfinの燃料が先のステップ940、950、955の何れかにて設定されている最終燃料噴射圧力Pcrfinをもって噴射される。
以上、説明したように、本発明による混合気状態推定装置、及びエミッション発生量推定装置の第1実施形態によれば、噴射期間TAUを「前期1/3TAU」「中期1/3TAU」「後期1/3TAU」の3つの期間に均等に区分し、「前期1/3TAU」に対応する1番目噴射(質量Q(1))が燃料噴射開始時期CAinjにて一時に実行され、1番目噴射から1/3TAU経過した後に「中期1/3TAU」に対応する2番目噴射(質量Q(2))が一時に実行され、2番目噴射から1/3TAU経過した後に「後期1/3TAU」に対応する3番目噴射(質量Q(3))が一時に実行されるものとする。そして、1番目噴射に基づく1番目混合気、2番目噴射に基づく2番目混合気、及び3番目噴射に基づく3番目混合気の「3つの混合気」をそれぞれ個別に扱い、混合気毎に、噴射後経過時間tに対する混合気の空気過剰率λ(筒内ガスが混ざり合っていく程度を表す値)、混合気の状態(温度Tmix等)及び混合気内でのエミッション発生量(スート発生量、及びNO発生量)を個別に推定する。
これにより、混合気間での、噴射後経過時間tに対する空気過剰率λの相違、噴射後経過時間tに対する筒内ガスの温度Tg、密度ρg等の状態の相違、噴射時点での有効噴射圧力ΔP及び筒内ガス密度ρg0(従って、噴霧角θ)の相違が考慮されて混合気の状態及びエミッション発生量が個別に推定され得る。そして、混合気毎に推定されたエミッション発生量を合計することでエミッションの総発生量(具体的には、スート総発生量Soot、NO総発生量NO)が推定される。これにより、上述した混合気の不均一性、エミッション発生量の不均一性が考慮されてエミッションの総発生量が精度良く推定され得る。
本発明は上記第1実施形態に限定されることはなく、本発明の範囲内において種々の変形例を採用することができる。例えば、上記第1実施形態では、指令燃料噴射量Qfin(従って、燃料噴射期間TAU)にかかわらず、噴射燃料(従って、混合気)を3つに区分しているが、指令燃料噴射量Qfinに応じて区分数を変更してもよい。この場合、指令燃料噴射量Qfinが大きいほど区分数を大きくすることが好適である。
また、上記第1実施形態においては、噴射期間TAUを複数の期間に均等に区分しているが、区分されたそれぞれの期間に対応する噴射燃料量が等しくなるように噴射期間TAUを複数の期間に区分してもよい。この場合、複数の期間は不均一になる。
また、上記第1実施形態においては、上記(15)式(図7のステップ765)にて算出される燃料消費量qrの対象となる化学反応として、着火反応(熱炎反応)、及び低温酸化反応(冷炎反応)のみならず、その他の種々の化学反応が含まれているが、着火反応、及び低温酸化反応に比してその他の種々の化学反応に係わる燃料消費量は十分に小さいことを考慮して、上記燃料消費量qrの対象となる化学反応として着火反応、及び低温酸化反応のみを扱うように構成してもよい。
この場合、燃料消費量qrを求めるための関数funcqrは、化学反応前の混合気温度Tpreが低温酸化反応が発生する温度範囲内にある場合に低温酸化反応により消費される燃料量の値を出力し、同化学反応前の混合気温度Tpreが着火反応が発生する温度範囲内にある場合には着火反応により消費される燃料量の値を出力し、同化学反応前の混合気温度Tpreがいずれの温度範囲内にもない場合には「0」を出力するように構成される。
また、上記第1実施形態においては、クランク角度CAがTDC以降であって、且つ、混合気温度Tmixがエミッションの反応限界温度より低くなった時点以降、CPU61の計算負荷低減のためにエミッション発生速度の時間積分処理が終了せしめられるように構成されているが、これに加え、クランク角度CAがTDC以前であって、且つ、混合気温度Tmixがエミッションの反応限界温度を超えるまでの間も、エミッション発生速度の時間積分処理を行わないように構成してもよい。これにより、エミッション発生量の算出に係わる不必要な計算が省略でき、より一層CPU61の計算負荷を低減することができる。
(第2実施形態)
次に、本発明の第2実施形態に係る内燃機関の混合気状態推定装置を含んだエミッション発生量推定装置について説明する。この第2実施形態は、燃料噴射期間TAUの途中において混合気が着火する場合を想定し、噴射燃料を混合気の着火前に噴射される部分と、同混合気の着火後に噴射される部分とに区分する点のみで上記第1実施形態と異なっている。従って、以下、係る相違点を中心に説明する。
(第2実施形態における噴射燃料の区分、混合気の区分)
燃料噴射期間TAUが比較的長い場合、燃料噴射期間TAUの途中において(即ち、燃料噴射継続中において)燃焼室内を進行していく混合気が着火する場合が多い。この場合、着火前に噴射された燃料に基づいて形成されていく混合気(以下、「着火前噴射に基づく混合気」と称呼する。)については、噴射から着火までの期間が比較的長いから、着火時点では既に混合気が十分に分散している。従って、予混合圧縮着火燃焼に類似した予混合的な燃焼が支配的になる場合が多い。
一方、着火後に噴射された燃料に基づいて形成されていく混合気(以下、「着火後噴射に基づく混合気」と称呼する。)については、噴射後直ちに着火することになるから、混合気が十分に分散していない状態で拡散しながら着火する。従って、拡散燃焼に類似した拡散的な燃焼が支配的になる場合が多い。
他方、予混合的燃焼と拡散的燃焼とでは燃料の反応速度(燃焼速度)が異なる。即ち、予混合的燃焼では、燃料と酸素とが十分に混ざり合った状態で着火が発生するから、燃料と容易に化学反応し得る酸素が十分に存在し得る。この結果、燃料の反応速度(燃焼速度)が比較的大きくなる。
これに対し、拡散的燃焼では、燃料と酸素が十分に混ざり合っていない状態で着火が発生するから、燃料と容易に化学反応し得る酸素が不足する。この結果、燃料の反応速度(燃焼速度)が比較的小さくなる。
更には、燃料の反応速度(燃焼速度)が異なると、上記(15)式にて表される微小時間Δtの間において混合気内で発生する化学反応による燃料消費量qr(従って、反応熱Hr)も異なる。即ち、着火前噴射に基づく混合気と、着火後噴射に基づく混合気とでは、上記(15)式の関数funcqrの引数である、混合気内の酸素濃度[O2]mix、混合気内の燃料濃度[Fuel]mix、混合気温度Tmixが同じであっても、燃料消費量qrが異なることになる。
係る着火前噴射に基づく混合気と着火後噴射に基づく混合気の間の燃料消費量qr(従って、反応熱Hr)の相違を考慮して混合気の状態、及びエミッション発生量を精度良く推定するためには、上記(15)式に相当する関数funcqrを、着火前噴射に基づく混合気用(関数funcqrpre)と着火後噴射に基づく混合気用(関数funcqrpost)の2種類準備し、噴射燃料を混合気の着火前に噴射される部分と同混合気の着火後に噴射される部分とに区分することが考えられる。そして、着火前噴射に基づく混合気の状態を上記着火前噴射に基づく混合気用の関数funcqrpreを利用して得られる燃料消費量qrに基づいて推定するとともに、着火後噴射に基づく混合気の状態を上記着火後噴射に基づく混合気用の関数funcqrpostを利用して得られる燃料消費量qrに基づいて推定することが考えられる。
そこで、第2実施形態は、着火遅れ時間(噴射開始から着火までの時間)Tdelayを推定し、図14に示すように、噴射期間TAUを「前期Tdelay」「後期(TAU−Tdelay)」の2つの期間に区分するとともに、前期Tdelay、後期(TAU−Tdelay)に対応してそれぞれ質量Q(1),Q(2)の燃料が順次個別に噴射されるものとする。着火遅れ時間Tdelayの推定については後述する。
より具体的には、着火前噴射(1番目噴射)に係わる噴射期間「前期Tdelay」分の質量Q(1)の燃料が燃料噴射開始時期CAinjにて一時に噴射され、着火前噴射からTdelay経過した後に着火後噴射(2番目噴射)に係わる噴射期間「後期(TAU−Tdelay)」分の質量Q(2)の燃料が一時に噴射されるものとする。なお、「Q(1)+Q(2)=Qfin」の関係がある。
第2実施形態は、着火前噴射に基づく混合気(1番目混合気)の状態を求める際には上記着火前噴射に基づく混合気用の関数funcqrpreを利用して得られる燃料消費量qrを使用し、着火後噴射に基づく混合気(2番目混合気)の状態を求める際には上記着火後噴射に基づく混合気用の関数funcqrpostを利用して得られる燃料消費量qrを使用する。
そして、第2実施形態は、着火前噴射に基づく混合気と着火後噴射に基づく混合気とをそれぞれ個別に扱い、混合気毎に、上述した第1実施形態と同様の手法により混合気の状態(温度Tmix等)及びエミッション発生量を個別に推定する。そして、第2実施形態は、混合気毎に推定されたエミッション発生量を合計することでエミッションの総発生量(具体的には、スート総発生量Soot、NO総発生量NO)を推定する。
これにより、混合気間での、噴射後の経過時間に対する筒内ガスの温度Tg、密度ρg等の状態の相違、噴射時点での有効噴射圧力ΔP及び筒内ガス密度ρg0(従って、噴霧角θ)の相違に加え、混合気間での上述した燃料消費量qr(従って、反応熱Hr)の相違が考慮されて混合気の状態及びエミッション発生量が個別に推定され得る。即ち、上述した混合気の不均一性、エミッション発生量の不均一性に加えて上述した燃料消費量qr(従って、反応熱Hr)の相違が考慮されてエミッションの総発生量が精度良く推定され得る。
(第2実施形態の実際の作動)
以下、第2実施形態に係る混合気状態推定装置を含んだエミッション発生量推定装置の実際の作動について説明する。この装置のCPU61は、第1実施形態のCPU61が実行する図5〜図9、及び図13に示したルーチンのうち、図6、図8、及び図13に示したルーチンをそのまま実行するとともに、図5、図7、及び図9のルーチンをそれぞれ以下のように一部だけ変更して実行する。以下、主として図5、図7、及び図9のルーチンの変更部分について説明する。
第2実施形態のCPU61は、図5に示したルーチンのステップ535〜545に代えて図15に示したステップ1505〜1520を実行する。即ち、CPU61はステップ530からステップ1505に進むと、現時点から最も近い過去において燃焼が発生した気筒の筒内圧力センサ77の出力経緯からその燃焼に係わる着火遅れ時間を取得し、同取得した着火遅れ時間を今回の燃焼についての着火遅れ時間Tdelayとする。
ここで、着火遅れ時間(噴射開始から着火までの時間)の取得は、着火発生時点にて筒内圧力が急増することに着目して筒内圧力センサ77の出力経緯から着火時期を特定することで取得し得る。このようにして取得された着火遅れ時間を今回の燃焼についての着火遅れ時間Tdelayとするのは、今回の燃焼についての着火遅れ時間Tdelayは現時点から最も近い過去における燃焼に係わる着火遅れ時間に近い値になるであろうとの予測に基づく。
次に、CPU61はステップ1510に進んで、上記求めた着火遅れ時間Tdelayと、図5のステップ530にて求めた基本燃料噴射圧力Pcrbaseと、上記関数funcPcrとから、1番目噴射(着火前噴射)、及び2番目噴射(着火後噴射)のそれぞれの燃料噴射圧力Pcr(1),Pcr(2)を求める。上述のごとく、着火前噴射(質量Q(1))は燃料噴射開始時期CAinjにて一時に実行され、着火後噴射(質量Q(2))は着火前噴射から着火遅れ時間Tdelay経過した後に一時に実行されるものとする。これにより、第1実施形態と同様、着火前噴射の燃料噴射圧力Pcr(1)は基本燃料噴射圧力Pcrbaseと等しい値に設定され、着火後噴射の燃料噴射圧力Pcr(2)は、基本燃料噴射圧力Pcrbaseよりも所定量だけ低い値に設定される。
続いて、CPU61はステップ1515に進み、図5のステップ520にて求めた燃料噴射期間TAU、着火遅れ時間Tdelayと、上記基本燃料噴射圧力Pcrbaseと、上記関数funcQとから、着火前噴射、及び着火後噴射のそれぞれの燃料量(質量)Q(1),Q(2)を求める。これにより、燃料量(質量)Q(1),Q(2)は、上述のごとく、「Qfin=Q(1)+Q(2)」の関係が成立するように設定される。このステップ1515が噴射燃料区分手段に相当する。
次いで、CPU61はステップ1520に進み、上記着火遅れ時間Tdelayと、現時点でのエンジン回転速度NEと、図5のステップ525にて決定された燃料噴射開始時期CAinjと、上記関数funcCAinjとから、着火前噴射、及び着火後噴射のそれぞれの噴射時点でのクランク角度CAinj(1),CAinj(2)を求める。これにより、着火前噴射時クランク角度CAinj(1)は上記燃料噴射開始時期CAinjに対応する値に設定され、着火後噴射時クランク角度CAinj(2)は燃料噴射開始時期CAinj(従って、着火前噴射時)から着火遅れ時間Tdelayだけ遅れた時期に対応する値に設定される。そして、CPU61は図5のステップ550へと進む。
以上、ステップ1505〜1520の処理は、変数i=1が着火前噴射(従って、着火前噴射に基づく混合気)に対応し、変数i=2が着火後噴射(従って、着火後噴射に基づく混合気)に対応することを意味している。
また、第2実施形態のCPU61は、図7に示したルーチンのステップ765に代えて図16に示したステップ1605〜1615を実行する。即ち、CPU61はステップ760からステップ1605に進むと、変数iの値が「1」であるか否かを判定する。
いま、変数i=1であるものとすると(即ち、着火前噴射に係わる算出がなされているものとすると)、CPU61はステップ1605にて「Yes」と判定してステップ1610に進み、その時点での混合気内の酸素濃度[O2]minと、燃料濃度[Fuel]mixと、図7のステップ760にて求めた化学反応前の混合気の温度Tpreと、上記(15)式に相当する着火前噴射に基づく混合気用の関数Funcqrpreとに基づいて微小時間Δt(噴射後経過時間(t−Δt)〜tの間)において混合気内で発生する化学反応による燃料消費量qrを求める。
一方、変数i=2であるものとすると(即ち、着火後噴射に係わる算出がなされているものとすると)、CPU61はステップ1605にて「No」と判定してステップ1615に進み、その時点での混合気内の酸素濃度[O2]minと、燃料濃度[Fuel]mixと、図7のステップ760にて求めた化学反応前の混合気の温度Tpreと、上記(15)式に相当する着火後噴射に基づく混合気用の関数Funcqrpostとに基づいて微小時間Δt(噴射後経過時間(t−Δt)〜tの間)において混合気内で発生する化学反応による燃料消費量qrを求める。そして、CPU61は図7のステップ770へと進む。
また、第2実施形態のCPU61は、図9に示したルーチンのステップ920に代えて図17に示したステップ1705を実行する。即ち、CPU61はステップ915からステップ1705に進むと、変数iが「2」であるか否かを判定し、変数iが「2」でない場合(即ち、i=1である場合)、図6のステップ605に戻る一方、変数iが「2」である場合、ステップ925以降へと進む。
これにより、ステップ1705にて「Yes」と判定された時点でエミッション発生量の推定が完了し、スート総発生量Sootが、初期値Soot0に、着火前噴射に基づく混合気内のスート発生量Soot(1)と、着火後噴射に基づく混合気内のスート発生量Soot(2)とを加えた値として求められる。同様に、NO総発生量NOが、初期値NO0に、着火前噴射に基づく混合気内のNO発生量NO(1)と、着火後噴射に基づく混合気内のNO発生量NO(2)とを加えた値として求められる。
以上、説明したように、本発明による混合気状態推定装置、及びエミッション発生量推定装置の第2実施形態によれば、噴射期間TAUを「前期Tdelay」「後期(TAU−Tdelay)」の2つの期間に区分し、「前期Tdelay」に対応する着火前噴射(質量Q(1))が燃料噴射開始時期CAinjにて一時に実行され、着火前噴射から着火遅れ時間Tdelay経過した後に「後期(TAU−Tdelay)」に対応する着火後噴射(質量Q(2))が一時に実行されるものとする。
そして、着火前噴射に基づく1番目混合気、及び着火後噴射に基づく2番目混合気の「2つの混合気」をそれぞれ個別に扱い、混合気毎に、噴射後経過時間tに対する混合気の空気過剰率λ、混合気の状態(温度Tmix等)及び混合気内でのエミッション発生量(スート発生量、及びNO発生量)を個別に推定する。そして、混合気毎に推定されたエミッション発生量を合計することでエミッションの総発生量が推定される。
加えて、着火前噴射に基づく混合気(1番目混合気)の状態を求める際には上記着火前噴射に基づく混合気用の関数funcqrpreを利用して得られる燃料消費量qrを使用し、着火後噴射に基づく混合気(2番目混合気)の状態を求める際には上記着火後噴射に基づく混合気用の関数funcqrpostを利用して得られる燃料消費量qrを使用する。
これにより、上述した混合気の不均一性、エミッション発生量の不均一性のみならず、上述した「着火前噴射に基づく混合気と着火後噴射に基づく混合気の間の燃料消費量qrの相違」が考慮されてエミッションの総発生量が精度良く推定され得る。
本発明は上記第2実施形態に限定されることはなく、本発明の範囲内において種々の変形例を採用することができる。例えば、上記第2実施形態では、噴射燃料を、混合気の着火前に噴射される部分と、同混合気の着火後に噴射される部分との「2つの部分」に区分するように構成されているが、混合気の着火前に噴射される部分を更に複数の(例えば、M個の)部分に区分し、混合気の着火後に噴射される部分も更に複数の(例えば、N個の)部分に区分することで噴射燃料を合計(M+N)個の部分に区分するように構成してもよい。
この場合、上記M個の部分に基づくM個の混合気の状態を求める際には、上記着火前噴射に基づく混合気用の関数funcqrpreを利用して得られる燃料消費量qrが使用され、上記N個の部分に基づくN個の混合気の状態を求める際には、上記着火後噴射に基づく混合気用の関数funcqrpostを利用して得られる燃料消費量qrが使用される。
これによると、上述した「着火前噴射に基づく混合気と着火後噴射に基づく混合気の間の燃料消費量qrの相違」が考慮されつつ、上述した混合気の不均一性、エミッション発生量の不均一性がより一層明確に考慮されてエミッションの総発生量がより一層精度良く推定され得る。
(第3実施形態)
次に、本発明の第3実施形態に係る内燃機関の混合気状態推定装置を含んだエミッション発生量推定装置について説明する。この第3実施形態は、燃料噴射期間TAUを微小時間Δt(例えば、0.1msec)毎に(多数の(n個の)期間に)区分して噴射燃料を同区分された各期間に対応して噴射されるそれぞれの部分(n個の部分)に区分する点、並びに、2番目以降に噴射される部分に基づいて形成されていく混合気の空気過剰率を上記(2)式、(3)式を利用することなく1番目に噴射される部分に基づいて形成されていく混合気の空気過剰率に基づいて決定する点で、上記第1、第2実施形態と異なっている。以下、係る相違点を中心に説明する。なお、以下、i番目(i:n以下の自然数)に噴射される部分に基づいて形成されていく混合気を「i番目混合気」と称呼する。
上記第1、第2実施形態では、(2つ、又は3つに区分された)混合気毎に、且つ、微小時間Δt(例えば、0.1msec)毎に、上記(2)式、(3)式を利用して空気過剰率λが算出されている(図7のステップ730、735を参照)。これは、混合気毎に、且つ、微小時間Δt毎に、混合気に新たに取り込まれる筒内ガス量gを求める必要があるからである(図7のステップ740を参照)。
ここで、(3)式の右辺には、時々刻々と変化する変数ρg,tに関するべき演算が含まれている。べき演算には多大な計算負荷が伴う。従って、上記(2)式、(3)式を利用した空気過剰率λの算出には多大な計算負荷が伴うから、係る算出回数をできるだけ減らすことが望ましい。
しかしながら、噴射燃料を多数の(n個の)部分に区分するこの第3実施形態で上記第1、第2実施形態と同様に、混合気毎、且つ、微小時間Δt毎に上記(2)式、(3)式を利用した空気過剰率λの算出を行うものとすると、係る算出回数が非常に大きな値となってしまい、計算負荷が膨大となる。
そこで、第3実施形態では、2番目以降に噴射される部分に基づいて形成されていくそれぞれの混合気の噴射後経過時間t=0からの微小時間Δt毎の空気過剰率は、1番目混合気の噴射後経過時間t=0からの微小時間Δt毎の空気過剰率と等しいと仮定される。係る仮定により、2番目以降に噴射される部分に基づいて形成されていくそれぞれの混合気について、微小時間Δt毎に混合気に新たに取り込まれる筒内ガス量gを上記(2)式、(3)式を利用することなく算出することができる。この結果、2番目以降に噴射される部分に基づいて形成されていくそれぞれの混合気についての上記(2)式、(3)式を利用した空気過剰率λの算出が一切不要となる。以下、この点について図18を参照しながら説明する。
図18は、n個に区分された噴射燃料のそれぞれの部分(q(1),q(2),・・・,q(n−1),q(n))に基づいて形成されていくn個の混合気の時間経過に対する様子(噴射開始から噴射終了までの様子)を模式的に示した図である。
図18に示すように、第3実施形態では、1番目噴射に係わる噴射期間「1番目の微小時間Δt」分の質量q(1)の燃料が燃料噴射開始時期CAinjにて一時に噴射され、2番目噴射に係わる噴射期間「2番目の微小時間Δt」分の質量q(2)の燃料が燃料噴射開始時期CAinjからΔt経過した後に一時に噴射され、・・・、n番目噴射に係わる噴射期間「n番目の微小時間Δt」分の質量q(n)の燃料が燃料噴射開始時期CAinjから(n−1)・Δt経過した後に一時に噴射されるものとする。なお、指令燃料噴射量Qfinとq(i)(i:n以下の自然数)には、下記(28)式に示す関係が成立するが、燃料噴射期間TAUに亘る噴射圧力の変動等に起因してq(i)(i:n以下の自然数)は互いに等しい値とはならない。
更には、図18に示すように、i番目混合気についてi番目混合気に係わる噴射後経過時間t=(k−1)・Δt
〜 k・Δtの間にi番目混合気に新たに取り込まれる筒内ガス量を「g(i,k)」と表すものとする(i:n以下の自然数、k:自然数。以下、同様)。
先ず、1番目噴射(従って、1番目混合気)のみについて考える。1番目混合気に係わる噴射後経過時間t=k・Δtにおける1番目混合気の空気過剰率λ(k)(第1番目部分混合指標値)は、上記第1、第2実施形態と同様、上記(2)式、(3)式を利用して順次求めることができる。
従って、1番目混合気に新たに取り込まれる筒内ガス量g(1,k)(図18において斜線で示した部分に対応する)は、上記(11)式に相当する下記(29)式に従って求めることができる。なお、λ(0)=0とする。即ち、1番目混合気に新たに取り込まれる筒内ガス量g(1,k)だけは、上記第1、第2実施形態と同様、上記(2)式、(3)式を利用して求められる。
ここで、図18、及び上記(29)式から理解できるように、1番目混合気に係わる噴射後経過時間t=k・Δtの時点で1番目混合気に取り込まれている筒内ガスの総量sumg(k)は、下記(30)式にて表される。
例えば、1番目混合気に係わる噴射後経過時間t=2・Δtの時点で1番目混合気に取り込まれている筒内ガスの総量sumg(2)=g(1,1)+g(1,2)であり、t=3・Δtの時点ではsumg(3)=g(1,1)+g(1,2)+g(1,3)となる。従って、1番目混合気に係わる噴射後経過時間t=k・Δtの空気過剰率λ(k)は、下記(31)式に従って表すことができる。
従って、上記(31)から下記(32)式が導出され得るから、値「λ(k)−λ(k−1)」を、既知の値である「g(1,k)」と「q(1)」とを利用して求めることができる。
ここで、上記仮定によれば、2番目以降のi番目混合気(i:2以上n以下の自然数)に係わる噴射後経過時間t=k・Δtにおけるi番目混合気の空気過剰率は、上記求められた「1番目混合気に係わる噴射後経過時間t=k・Δtにおける1番目混合気の空気過剰率λ(k)」と等しい。従って、2番目以降のi番目混合気についてi番目混合気に係わる噴射後経過時間t=(k−1)・Δt
〜 k・Δtの間にi番目混合気に新たに取り込まれる筒内ガス量g(i,k)は、上述した仮定を利用すれば、上述した「1番目混合気に新たに取り込まれる筒内ガス量g(1,k)」を表す上記(29)式と同様に、下記(33)式にて表すことができる。
この上記(33)式に上記(32)式を代入すると、下記(34)式を得ることができる。この下記(34)式によれば、2番目以降のi番目混合気に新たに取り込まれる筒内ガス量g(i,k)を、計算負荷の大きい上記(2)式、(3)式を利用することなく、既知の値である「g(1,k)」と「q(i)」と「q(1)」とを利用して簡易に求めることができる。即ち、2番目以降のi番目混合気に新たに取り込まれる筒内ガス量g(i,k)を求めるために、2番目以降のi番目混合気についての上記(2)式、(3)式を利用した空気過剰率λを算出する必要がなくなる。
このように、第3実施形態では、2番目以降のi番目の混合気の噴射後経過時間t=0からの微小時間Δt毎の空気過剰率が、1番目混合気の噴射後経過時間t=0からの微小時間Δt毎の空気過剰率λ(k)と等しいと仮定することで、上記(2)式、(3)式を利用した空気過剰率λの算出を必要とする混合気を1番目混合気のみとすることができる。従って、計算負荷の大きい上記(2)式、(3)式を利用した空気過剰率λの算出回数を減らすことができ、CPU61の計算負荷を減らすことができる。
なお、この第3実施形態では、第1、第2実施形態と同様に混合気間での噴射後経過時間tに対する筒内ガスの状態(温度Tg、圧力Pg等)の相違は考慮され得る一方、混合気間での噴射時点での有効噴射圧力ΔP及び筒内ガス密度ρg0(従って、噴霧角θ)の相違、即ち、混合気間での噴射後経過時間tに対する空気過剰率λの相違が考慮されることなく、混合気の状態及びエミッション発生量が個別に推定されることになる。従って、第1、第2実施形態に比して上述した混合気の不均一性、及びエミッション発生量の不均一性が考慮される程度は低くなる。
(第3実施形態の実際の作動)
以下、第3実施形態に係る混合気状態推定装置を含んだエミッション発生量推定装置の実際の作動について説明する。この装置のCPU61は、第1実施形態のCPU61が実行する図5〜図9に示した一連のルーチン、及び図13に示したルーチンのうち、図13に示したルーチンのみをそのまま実行するとともに、図5〜図9に示した一連のルーチンに代えて図5〜図9のルーチンにそれぞれ対応する図19〜図23にフローチャートにより示した一連のルーチンを実行する。加えて、この装置のCPU61は、図24〜図27にフローチャートにより示した一連のルーチンを追加的に実行する。
なお、図19〜図27に示したルーチンのステップであって図5〜図9に示したルーチンのステップと同一のものについては、図5〜図9に示したルーチンの対応するステップの番号と同一の番号を付すことでそれらの説明を省略する。以下、第3実施形態に特有の図19〜図27に示した各ルーチンについて説明する。
第3実施形態のCPU61は、図5〜図9に示した一連のルーチンに対応する図19〜図23に示した一連のルーチンを所定時間の経過毎に繰り返し実行するようになっている。従って、所定のタイミングになると、CPU61は図19のステップ1900から処理を開始し、ステップ505に進んで、「Yes」と判定する場合(即ち、IVCが到来した場合)、ステップ510〜530の処理を順に実行し、その後、ステップ1905に進んで、ステップ520にて求めた燃料噴射期間TAUを微小時間Δt(例えば、0.1msec)で除することで区分数nを求める。
次いで、CPU61はステップ1910に進み、図5のステップ540と同様にして、上記求めた燃料噴射期間TAU、微小時間Δtと、上記求めた基本燃料噴射圧力Pcrbaseと、上記関数funcQとから、1番目噴射からn番目噴射までのそれぞれの燃料量(質量)q(1),q(2),・・・,q(n)を求める。これにより、燃料量(質量)q(1),q(2),・・・,q(n)は、上述のごとく、「Qfin=q(1)+q(2)+・・・+q(n)」の関係が成立するように設定される一方で、互いに等しい値とはならない。このステップ1910は噴射燃料区分手段に相当する。
続いて、CPU61はステップ1915に進んで、図5のステップ545と同様にして、上記微小時間Δtと、現時点でのエンジン回転速度NEと、ステップ525にて決定された燃料噴射開始時期CAinjと、上記関数funcCAinjとから、1番目噴射からn番目噴射までのそれぞれの噴射時点でのクランク角度CAinj(1),CAinj(2),・・・,CAinj(n)を求める。
これにより、1番目噴射時クランク角度CAinj(1)は上記燃料噴射開始時期CAinjに対応する値に設定され、2番目噴射時クランク角度CAinj(2)は燃料噴射開始時期CAinjからΔtだけ遅れた時期に対応する値に設定され、・・・、n番目噴射時クランク角度CAinj(n)は燃料噴射開始時期CAinjから(n−1)・Δtだけ遅れた時期に対応する値に設定される。
次に、CPU61はステップ550〜570の処理を順に実行し、その後、図20のルーチンに進んで1番目噴射のみに関する各種初期値の決定を行うための処理を行う。具体的には、CPU61はステップ2005に進んで、図6のステップ610と同様にして、ステップ515にて求めた筒内ガスの全質量Mgを、先のステップ1915にて求めた1番目噴射時クランク角度CAinj(1)から得られる1番目噴射時筒内容積Vg(CAinj(1))で除することで、1番目噴射時の筒内ガス密度ρg0を求める。
続いて、CPU61はステップ2010に進み、図6のステップ615と同様にして、ステップ510にて求めたIVC時筒内ガス圧力Pgivcと、上記IVC時筒内容積Vg(CAivc)と、上記1番目噴射時筒内容積Vg(CAinj(1))と、上記(4)式に相当する式とに基づいて1番目噴射時の筒内ガス圧力Pg0を求める。
次に、CPU61はステップ2015に進んで、ステップ530にて求めた基本燃料噴射圧力Pcrbaseから上記1番目噴射時の筒内ガス圧力Pg0を減じることで1番目噴射時の有効噴射圧力ΔPを求め、続くステップ625にて、上記求めた有効噴射圧力ΔPと、筒内ガス密度ρg0と、上記テーブルMapθとに基づいて1番目噴射に係わる噴霧角θを求める。これにより、噴霧角θは、1番目噴射時点(即ち、クランク角度CAinj(1))での有効噴射圧力ΔP及び筒内ガス密度ρg0に基づいて決定される。
続いて、CPU61はステップ2020に進み、上述したように、計算の便宜上、1番目混合気に係わる噴射後経過時間t=0における1番目混合気の空気過剰率λ(0)(この値は、後述する図21のステップ2110にて使用される。)を「0」に設定する。
次に、CPU61はステップ635、640の処理を順に実行し、その後、ステップ2025に進んで、図6のステップ645と同様にして、1番目混合気に係わる混合気のエンタルピHmixを、上記(17)式に相当する式に従って、初期値(即ち、先のステップ1910にて求めた1番目噴射量q(1)と、燃料の定圧比熱Cfと、ステップ550にて求めた燃料蒸気温度Tfの積)に設定する。
次いで、CPU61はステップ650の処理を実行し、その後、ステップ2030に進み、1番目混合気に係わる混合気質量Mmixを初期値である1番目噴射量q(1)に設定する。続いて、CPU61はステップ660の処理を実行し、その後、ステップ2035に進んで、図6のステップ665と同様、1番目混合気に係わる噴射後経過時間tを初期値「0」に設定する。また、CPU61は、1番目噴射に係わるクランク角度CAを初期値である1番目噴射時クランク角度CAinj(1)に設定する。これにより、1番目混合気に係わる噴射後経過時間tが1番目噴射時からカウントされることになる。加えて、このステップ2035では、CPU61は変数kの値を「0」に設定する。この変数kの値は、噴射後経過時間tが「k・Δt」であることを表す。
次に、CPU61はステップ670の処理を実行し、その後、図21のルーチンに進んで、1番目噴射に関する混合気温度の算出のための処理を開始する。具体的には、CPU61は、先ず、ステップ2105に進み、図7のステップ705と同様に、1番目混合気に係わる噴射後経過時間tをΔtだけ進めるとともに、1番目噴射に係わるクランク角度CAを図19のステップ565で求めたΔCAだけ進める。
加えて、このステップ2105では、CPU61は変数kの値を「1」だけインクリメントする。これにより、1番目噴射に係わるクランク角度CAの値、及び変数kの値が、1番目混合気に係わる噴射後経過時間tに対応する値に維持されていく。
次に、CPU61はステップ710〜730の処理を順に実行し、その後、ステップ2110に進み、図7のステップ735と同様、1番目混合気に係わる噴射後経過時間t=k・Δtにおける1番目混合気の空気過剰率λ(k)(第1番目部分混合気指標値)を、上記(2)式、(3)式に従って、1番目混合気に係わる噴射後経過時間t=(k−1)・Δtにおける1番目混合気の空気過剰率λ(k−1)(k=1の場合、λ(k−1)=λ(0)=0)に、ステップ730にて求めた燃料希釈率dλ/dtにΔtを乗じた値「dλ/dt・Δt」を加えることで求める。このステップ2110が混合指標値取得手段に相当する。
続いて、CPU61はステップ2115に進んで、ステップ2110にて既に求められているλ(k),λ(k-1)と、図19のステップ1910に求められているq(1)と、上記(29)式とに従って、1番目混合気に係わる噴射後経過時間t=(k−1)・Δt
〜 k・Δtの間に1番目混合気に新たに取り込まれる筒内ガス量g(1,k)を求める。
次いで、CPU61はステップ2120に進み、図7のステップ745と同様、混合気形成筒内ガス質量Gを、その時点での値(初期値は、図20のステップ635により「0」に設定されている。)に上記求めた1番目混合気に新たに取り込まれる筒内ガス質量g(1,k)を加えた値に更新する。これにより、噴射後経過時間t=k・Δtにおける1番目混合気に係わる混合気形成筒内ガス質量Gが求められる。
次に、CPU61はステップ2125に進んで、図7のステップ750と同様、混合気質量Mmixを、その時点での値(初期値は、図20のステップ2030により「q(1)」に設定されている。)に上記1番目混合気に新たに取り込まれる筒内ガス質量g(1,k)を加えた値に更新する。これにより、噴射後経過時間t=k・Δtにおける1番目混合気に係わる混合気質量Mmixが求められる。
続いて、CPU61はステップ2130に進み、図7のステップ755と同様、化学反応前の混合気のエンタルピHpreを、その時点での混合気のエンタルピHmix(初期値は、図20のステップ2025の処理により値「q(1)・Cf・Tf」)に、「上記1番目混合気に新たに取り込まれる筒内ガスの熱エネルギーΔHg=g(1,k)・Cg・Tg」を加えた値に設定する。
次に、CPU61はステップ760〜780の処理を順に実行する。これにより、ステップ780において、噴射後経過時間t=k・Δt(従って、クランク角度CA=CAinj(1)+k・ΔCA)における1番目混合気の温度Tmixが求められる。
次いで、CPU61は図22のルーチンに進み、1番目噴射に関する各種濃度の算出のための処理を開始する。この図22のルーチンは、図8のルーチンのステップ815をステップ2205に代えた点のみにおいて図8のルーチンと異なる。このステップ2205では、CPU61は、図19のステップ1910にて求めた1番目噴射量q(1)と、ステップ810にて求めた燃料消費量積算値sumqrと、図21のステップ2125にて求めた混合気質量Mmixと、上記(19)式に相当する式とに基づいて、噴射後経過時間t=k・Δtにおける1番目混合気内の燃料濃度[Fuel]mixを求める。
そして、CPU61は図23のルーチンに進み、ステップ905の判定において「No」と判定する毎に、図21のステップ2105〜図23のステップ905の処理を繰り返し実行する。これにより、ステップ905の判定において「No」と判定される毎に、図21のステップ2105にて、変数kの値が「1」だけインクリメントされ、1番目混合気に係わる噴射後経過時間tがΔtだけ進められるとともに、1番目噴射に係わるクランク角度CAがΔCAだけ進められる。
即ち、ステップ905の判定において「No」と判定される限りにおいて、1番目混合気の空気過剰率λ(k)、1番目混合気に新たに取り込まれる筒内ガス質量g(1,k)、1番目混合気の状態(温度Tmix等)、及び(フラグENDsoot=フラグENDno=0である限りにおいて)1番目混合気に係わるエミッション濃度([Soot]mix,[NO]mix)が1番目噴射時から微小時間Δt毎に(即ち、CAinj(1)から微小クランク角度ΔCA毎に)求められていく。
そして、ステップ905の条件が成立すると、CPU61は図23のステップ905に進んだとき「Yes」と判定してステップ2305以降に進み、1番目混合気に係わる計算を終了するための処理を行う。
即ち、CPU61はステップ2305に進むと、図9のステップ910と同様に、1番目混合気内のスート発生量Soot(1)、及び1番目混合気内のNO発生量NO(1)を求める。次いで、CPU61はステップ2310に進み、図9のステップ915と同様、スート総発生量Sootを、その時点での値(現時点では、図19のステップ570の処理により初期値Soot0)に上記求めた1番目混合気内のスート発生量Soot(1)を加えた値(Soot0+Soot(1))に更新するとともに、NO総発生量NOを、その時点での値(現時点では、図19のステップ570の処理により初期値NO0)に上記求めた1番目混合気内のNO発生量NO(1)を加えた値(NO0+NO(1))に更新する。以上により、1番目噴射(従って、1番目混合気)に関する計算が終了する。
次に、CPU61はステップ2315に進んで、変数iを「1」に設定する。ここで、変数iは、第1、第2実施形態と同様、何番目の噴射か(従って、何番目の混合気か)を識別するための値である。
続いて、CPU61はステップ2320に進み、変数iの値を「1」だけインクリメントし、続くステップ2325を経由して後述する図24〜図27に示した一連の「i(2≦i≦n)番目噴射に関するエミッション量の算出」ルーチンを実行することにより、i番目混合気内のエミッション発生量である、スート発生量Soot(i)、及びNO発生量NO(i)を求める(2≦i≦n)。ここで、値nは、図19のステップ1905にて求められた噴射燃料の区分数である。
次いで、CPU61はステップ2330に進んで、スート総発生量Sootを、その時点での値(i=2の場合、Soot0+Soot(1))に上記求めたi番目混合気内のスート発生量Soot(i)を加えた値に更新するとともに、NO総発生量NOを、その時点での値(i=2の場合、NO0+NO(1))に上記求めたi番目混合気内のNO発生量NO(i)を加えた値に更新する。これにより、スート総発生量Soot=Soot0+Soot(1)+・・・+Soot(i)となり、NO総発生量NO=NO0+NO(1)+・・・+NO(i)となる。
そして、CPU61はステップ2335に進んで、変数iの値が上記区分数nに一致したか否かを判定し、「No」と判定する場合、ステップ2320に戻る。即ち、ステップ2320の繰り返し実行により変数iが区分数nに達するまでの間、ステップ2320〜2330の処理が繰り返し実行される。この結果、変数iの値が「1」ずつインクリメントされていくとともに、ステップ2330にてスート総発生量Sootの値、及びNO総発生量NOが更新されていく。
そして、変数iの値が区分数nに達すると、CPU61はステップ2335に進んだとき「Yes」と判定してステップ925に進む。この時点で、エミッション発生量の推定が終了して、スート総発生量Sootが「Soot0+Soot(1)+・・・+Soot(n)」に確定し、NO総発生量NOが「NO0+NO(1)+・・・+NO(n)」に確定する。CPU61は、以上説明した混合気状態の推定、及びエミッション発生量の推定をIVCの直後(即ち、燃料噴射開始時期CAinjの前の時点)において完了する。
以下、図24〜図27に示した一連の「i(2≦i≦n)番目噴射に関するエミッション量の算出」ルーチンについて説明する。このルーチンは、i(2≦i≦n)番目混合気のスート発生量Soot(i)、及びNO発生量NO(i)を求めるルーチンであって、図24〜図27のルーチンはそれぞれ、図6〜図9のルーチンに対応している。
CPU61は図23のステップ2320を実行した後、ステップ2325を経由して、先ず、i(2≦i≦n)番目噴射に関する各種初期値を決定するために図24のルーチンから実行する。具体的には、CPU61は、先ず、ステップ635、640の処理を順に実行し、その後、ステップ2405に進んで、図20のステップ2025と同様にして、i番目混合気に係わる混合気のエンタルピHmixを、初期値(即ち、先のステップ1910にて求めたi番目噴射量q(i)と、燃料の定圧比熱Cfと、ステップ550にて求めた燃料蒸気温度Tfの積)に設定する。
次いで、CPU61はステップ650の処理を実行し、その後、ステップ2410に進み、i番目混合気に係わる混合気質量Mmixを初期値である上記i番目噴射量q(i)に設定する。続いて、CPU61はステップ660の処理を実行し、その後、ステップ2415に進んで、図20のステップ2035と同様、i番目噴射に係わるクランク角度CAを初期値であるi番目噴射時クランク角度CAinj(i)に設定する。加えて、このステップ2415では、CPU61は変数kの値を「0」に設定する。この変数kの値は、i番目噴射からの経過時間が「k・Δt」であること(即ち、i番目噴射に係わるクランク角度CAが「CAinj(i)+k・ΔCA」であること)を表す。
次に、CPU61はステップ670の処理を実行し、その後、図25のルーチンに進んで、i(2≦i≦n)番目噴射に関する混合気温度の算出のための処理を開始する。具体的には、CPU61は、先ず、ステップ2505に進み、図21のステップ2105と同様に、i番目噴射に係わるクランク角度CAを図19のステップ565で求めたΔCAだけ進める。加えて、このステップ2505では、CPU61は変数kの値を「1」だけインクリメントする。これにより、i番目噴射に係わるクランク角度CAの値と変数kの値とが、互いに対応する値に維持されていく。
次に、CPU61はステップ715〜725の処理を順に実行し、その後、ステップ2510に進み、図21のステップ2115にて求めた「1番目混合気に新たに取り込まれる筒内ガス量g(1,k)」と、図19のステップ1910にて求めたq(i)と、同ステップ1910にて求めたq(1)と、上記(34)式とに基づいて、上記(2)式、(3)式を利用することなく、i(2≦i≦n)番目混合気に係わる噴射後経過時間t=(k−1)・Δt
〜 k・Δtの間にi番目混合気に新たに取り込まれる筒内ガス量g(i,k)を求める。
次いで、CPU61はステップ2515に進み、図21のステップ2120と同様、混合気形成筒内ガス質量Gを、その時点での値(初期値は、図24のステップ635により「0」に設定されている。)に上記求めたi番目混合気に新たに取り込まれる筒内ガス質量g(i,k)を加えた値に更新する。これにより、噴射後経過時間t=k・Δtにおけるi番目混合気に係わる混合気形成筒内ガス質量Gが求められる。
次に、CPU61はステップ2520に進んで、図21のステップ2125と同様、混合気質量Mmixを、その時点での値(初期値は、図24のステップ2410により「q(i)」に設定されている。)に上記i番目混合気に新たに取り込まれる筒内ガス質量g(i,k)を加えた値に更新する。これにより、噴射後経過時間t=k・Δtにおけるi番目混合気に係わる混合気質量Mmixが求められる。
続いて、CPU61はステップ2525に進み、図21のステップ2130と同様、化学反応前の混合気のエンタルピHpreを、その時点での混合気のエンタルピHmix(初期値は、図24のステップ2405の処理により値「q(i)・Cf・Tf」)に、「上記i番目混合気に新たに取り込まれる筒内ガスの熱エネルギーΔHg=g(i,k)・Cg・Tg」を加えた値に設定する。
次に、CPU61はステップ760〜780の処理を順に実行する。これにより、ステップ780において、噴射後経過時間t=k・Δt(従って、クランク角度CA=CAinj(i)+k・ΔCA)におけるi(2≦i≦n)番目混合気の温度Tmixが求められる。
次いで、CPU61は図26のルーチンに進み、i番目噴射に関する各種濃度の算出のための処理を開始する。この図26のルーチンは、図8のルーチンのステップ815をステップ2605に代えた点のみにおいて図8のルーチンと異なる。このステップ2605では、CPU61は、図19のステップ1910にて求めたi番目噴射量q(i)と、ステップ810にて求めた燃料消費量積算値sumqrと、図25のステップ2520にて求めた混合気質量Mmixと、上記(19)式に相当する式とに基づいて、噴射後経過時間t=k・Δtにおけるi番目混合気内の燃料濃度[Fuel]mixを求める。
そして、CPU61は図27のルーチンに進み、ステップ905の判定において「No」と判定する毎に、図25のステップ2505〜図27のステップ905の処理を繰り返し実行する。これにより、ステップ905の判定において「No」と判定される毎に、図25のステップ2505にて、変数kの値が「1」だけインクリメントされ、i番目噴射に係わるクランク角度CAがΔCAだけ進められる。
即ち、ステップ905の判定において「No」と判定される限りにおいて、i番目混合気に新たに取り込まれる筒内ガス質量g(i,k)、i番目混合気の状態(温度Tmix等)、及び(フラグENDsoot=フラグENDno=0である限りにおいて)i番目混合気に係わるエミッション濃度([Soot]mix,[NO]mix)がi番目噴射時から微小時間Δt毎に(即ち、CAinj(i)から微小クランク角度ΔCA毎に)求められていく。
そして、ステップ905の条件が成立すると、CPU61は図27のステップ905に進んだとき「Yes」と判定してステップ910に進んで、i番目混合気内のスート発生量Soot(i)、及びi番目混合気内のNO発生量NO(i)を求める。i番目混合気内のスート発生量Soot(i)は、図26のステップ855の処理にて更新されている現時点でのi番目混合気内のスート濃度[Soot]mixの値に、図25のステップ2520にて更新されている現時点での混合気質量Mmixを乗じることで求められる。また、i番目混合気内のNO発生量NO(i)は、図26のステップ880の処理にて更新されている現時点でのi番目混合気内のNO濃度[NO]mixの値に、上記混合気質量Mmixを乗じることで求められる。
そして、CPU61はステップ2495を経由して(即ち、図24〜図27に示した一連のルーチンの実行を終了して)上述した図23のステップ2330に進む。このようにして、上述した図23のステップ2320〜2330の処理が繰り返し実行される毎に、変数i(2≦i≦n)の値が「1」ずつインクリメントされていくとともに、図27のステップ910にてi番目混合気内のスート発生量Soot(i)、及びi番目混合気内のNO発生量NO(i)が順次求められていく。以上、図24〜図27に示した一連のルーチンについて説明した。
上述したように、変数iの値が区分数nに達すると、CPU61は図23のステップ2335に進んだとき「Yes」と判定してステップ925〜960を実行し、図23のステップ2330にて求められているスート総発生量Sootの値、及びNO総発生量NOの値(従って、エミッション総発生量の推定結果)に基づいて噴射圧力を補正し、ステップ1995に進んで図19〜図23の一連の本ルーチンを一旦終了する。以降、CPU61は、次のIVCが到来するまでの間、図19のステップ505に進む毎に「No」と判定し続ける。
以上、説明したように、本発明による混合気状態推定装置、及びエミッション発生量推定装置の第3実施形態によれば、噴射期間TAUが多数のn(=TAU/Δt)個の期間に区分され、i(1≦i≦n)番目噴射に係わる噴射期間「i番目の微小時間Δt」分の質量q(i)の燃料が燃料噴射開始時期CAinjから(i−1)・Δt経過した後に一時に噴射されるものとする。1番目混合気に係わる噴射後経過時間t=k・Δtにおける1番目混合気の空気過剰率λ(k)(第1番目部分混合指標値)は、上記第1、第2実施形態と同様、上記(2)式、(3)式を利用して求められ、1番目混合気の状態(温度Tmix等)は係るλ(k)を利用して推定される。
そして、2番目以降のi番目混合気(i:2以上n以下の自然数)に係わる噴射後経過時間t=k・Δtにおけるi番目混合気の空気過剰率は、上記「1番目混合気に係わる噴射後経過時間t=k・Δtにおける1番目混合気の空気過剰率λ(k)」と等しいとの仮定のもと、上記(2)式、(3)式を利用することなく、2番目以降のi番目混合気の状態(温度Tmix等)が推定される。
従って、上記(2)式、(3)式を利用した空気過剰率λの算出を必要とする混合気を1番目混合気のみとすることができる。この結果、計算負荷の大きい上記(2)式、(3)式を利用した空気過剰率λの算出回数を減らすことができるから、CPU61の計算負荷を減らすことができる。
本発明は上記第3実施形態に限定されることはなく、本発明の範囲内において種々の変形例を採用することができる。例えば、上記第3実施形態では、2番目以降のi番目混合気(i:2以上n以下の自然数)に係わる噴射後経過時間t=k・Δtにおけるi番目混合気の空気過剰率は、上記「1番目混合気に係わる噴射後経過時間t=k・Δtにおける1番目混合気の空気過剰率λ(k)」と等しいと仮定されているが、同空気過剰率λ(k)に所定の係数を乗じた値になると仮定してもよい。
即ち、2番目以降のi番目混合気に係わる噴射後経過時間t=k・Δtにおけるi番目混合気の空気過剰率をλ(i,k)(i:2以上n以下の自然数)とすると、例えば、「λ(i,k)=h(i)・λ(k)」と設定される。ここで、h(i)は変数iに応じて決定される係数である。
この場合、i(2≦i≦n)番目混合気に係わる噴射後経過時間t=(k−1)・Δt 〜
k・Δtの間にi番目混合気に新たに取り込まれる筒内ガス量g(i,k)は、図25のステップ2510に記載された上記(34)式ではなく、上記(33)式に相当する式である「g(i,k)=stoich・(λ(i,k)−λ(i,k−1))・q(i)」に従って取得され得る。
(第4実施形態)
次に、本発明の第4実施形態に係る内燃機関のエミッション発生量推定装置について説明する。この第4実施形態は、混合気の着火後において燃焼室内で所謂定常火炎が発生している場合を想定し、係る定常火炎に特有の特性を利用して混合気の状態(温度等)、及びエミッション発生量を推定する点で、上記第1〜第3実施形態と異なっている。以下、係る相違点を中心に説明する。
噴射開始時点以降における噴射された燃料(従って、混合気)の燃料噴射弁21の噴孔からの到達距離(以下、「混合気到達距離X」と称呼する。)は、例えば、上記非特許文献1にて紹介された実験式である下記(35)式、(36)式に従って噴射後経過時間tの関数として求めることができる。下記(36)式において、dX/dtは噴射後経過時間tの関数である混合気移動速度である。なお、下記(36)式の右辺に示される各種値は、上記(3)式の右辺に示されるものと同一である。
ここで、上記(3)式の両辺を上記(36)式の両辺で除することにより、dλ/dXを、ρgとtanθと各種定数とを用いて表すことができる。いま、ρgとtanθが噴射後において一定であるとすると、dλ/dXは或る定数(正の値)になる。更に、空気過剰率λは混合気到達距離X=0(従って、噴射時点)において「0」である。
従って、この場合、混合気到達距離Xと空気過剰率λとの関係は図28に示すような線形関係となる。即ち、混合気到達距離Xが「0」から増加するにつれて空気過剰率λが「0」から増加し、X=X0のとき、λ=1となる。
ところで、燃料の噴射期間TAUが比較的長いと、図28に示すように、混合気の着火後において燃焼室内にて所謂定常火炎(或いは、定常火炎に非常に近い形態を有する火炎)が発生する場合がある。
定常火炎が発生している場合、定常火炎内における燃料が過剰な領域(即ち、空気過剰率λ<1となる領域。以下、「リッチ領域」と称呼する。)では、混合気内に取り込まれた筒内ガス中の酸素が不足していることにより酸素が全て燃焼により消費されている。即ち、定常火炎のリッチ領域内(図28では、0≦X≦X0に対応する。)では、定常状態における酸素濃度(以下、「定常酸素濃度[O2]mixsteady」と称呼する。)が「0」になっている。なお、定常火炎のリッチ領域内では、酸素が全て燃焼により消費されても燃料が残存しているから、定常状態における燃料濃度(以下、「定常燃料濃度[Fuel]mixsteady」と称呼する。)は「0」より大きい値になっている。
一方、定常火炎が発生している場合、定常火炎内における酸素が過剰な領域(即ち、空気過剰率λ>1となる領域。以下、「リーン領域」と称呼する。)では、燃料が不足していることにより燃料が全て燃焼により消費されている。即ち、定常火炎のリーン領域内(図28では、X0≦Xに対応する。)では、定常燃料濃度[Fuel]mixsteadyが「0」になっている。なお、定常火炎のリーン領域内では、燃料が全て燃焼により消費されても酸素が残存しているから、定常酸素濃度[O2]mixsteadyは「0」より大きい値になっている。
他方、上述したように、混合気内のスート発生速度d[Soot]mix/dtを求める上記(23)式には、上記(24)式に示すように混合気内の燃料濃度[Fuel]mixの値を因数に含むスートの生成速度dmsf/dtを求める項と、上記(25)式に示すように混合気内の酸素濃度[O2]mixの値を因数に含むスートの酸化速度dmso/dtを求める項とが存在する。
従って、定常火炎のリッチ領域内におけるスート発生速度d[Soot]mix/dtを上記(23)式に従って計算する場合、上記(25)式に示すスートの酸化速度dmso/dtを求める項の値が常に「0」に維持される。換言すれば、この場合、上記スートの酸化速度dmso/dtを求める項の計算を省略することができる。
同様に、定常火炎のリーン領域内におけるスート発生速度d[Soot]mix/dtを上記(23)式に従って計算する場合、上記(24)式に示すスートの生成速度dmsf/dtを求める項の値が常に「0」に維持される。換言すれば、この場合、上記スートの生成速度dmsf/dtを求める項の計算を省略することができる。
加えて、上記(24)式、(25)式の右辺には、時々刻々と変化する変数Pg,Tmixに関するべき演算が含まれている。べき演算には多大な計算負荷が伴う。従って、(24)式、(25)式の算出回数をできるだけ減らすことが望ましい。
そこで、定常火炎が占める領域内(即ち、混合気内)におけるスート発生速度d[Soot]mix/dtを上記(23)式に基づいて推定する第4実施形態は、定常火炎のリッチ領域内(λ<1)では、上記(23)式においてdmso/dtの項を省略することで得られる式である「d[Soot]mix/dt=dmsf/dt」に従ってスート発生速度d[Soot]mix/dtを求める。
同様に、第4実施形態は、定常火炎のリーン領域内(λ>1。実際には、λ≧1)では、上記(23)式においてdmsf/dtの項を省略することで得られる式である「d[Soot]mix/dt=−dmso/dt」に従ってスート発生速度d[Soot]mix/dtを求める。これらにより、(24)式、(25)式の算出回数を減らすことができ、この結果、CPU61の計算負荷を減らすことができる。
以下、定常火炎が占める領域内(即ち、混合気内)におけるスート発生速度d[Soot]mix/dtを上記(23)式を利用して求める際に必要となる、混合気温度(具体的には、定常状態における燃焼による混合気温度上昇量。以下、「定常混合気温度上昇量ΔTmixsteady」と称呼する。)、上記定常燃料濃度[Fuel]mixsteady、及び、上記定常酸素濃度[O2]mixsteady、並びに、定常火炎が占める領域内(即ち、混合気内)におけるNO発生速度d[NO]mix/dtを上記(27)式を利用して求める際に必要となる、定常状態における窒素濃度(以下、「定常窒素濃度[N2]mixsteady」と称呼する。)を求める手法について順に説明する。
以下、係る説明に必要となる空気過剰率λを、上記筒内ガス理論空燃比stoichと、燃料量(質量)Qと、混合気形成筒内ガス質量Gとを用いて下記(37)式に従って定義する。
<定常混合気温度上昇量ΔTmixsteady>
定常混合気温度上昇量ΔTmixsteadyは、下記(38)式にて表すことができる。下記(38)式において、Qreacは混合気内で燃焼により消費される燃料量(質量)であり、Cgは筒内ガスの定圧比熱であり、Cfは燃料の定圧比熱である。HfはQreacを反応熱Hr(=Hf・Qreac)に変換するための係数である。この(38)式は、「混合気内の質量Qreacの燃料の燃焼により発生した反応熱Hrは、質量Qの燃料と質量Gの筒内ガスとを含む混合気の温度を定常混合気温度上昇量ΔTmixsteadyだけ上昇させるために消費される」との仮定に基づくものである。
ここで、定常火炎のリッチ領域内(λ<1)では、混合気内で燃焼により消費される燃料量Qreacは、混合気内における筒内ガス中の全ての酸素と反応して消費される分の燃料量となるから、下記(39)式にて求めることができる。
以上より、上記(37)式〜(39)式からQ,G,Qreacを消去して整理すると、定常火炎のリッチ領域内(λ<1)における定常混合気温度上昇量ΔTmixsteadyは、空気過剰率λの関数として下記(40)式にて表すことができる。
他方、定常火炎のリーン領域内(λ>1。実際には、λ≧1)では、混合気内の燃料が燃焼により全て消費されるから、混合気内で燃焼により消費される燃料量Qreacは「Q」と等しくなる。以上より、「Qreac=Q」の関係と、上記(37)式及び(38)式とからQ,G,Qreacを消去して整理すると、定常火炎のリーン領域内(λ≧1)における定常混合気温度上昇量ΔTmixsteadyは、空気過剰率λの関数として下記(41)式にて表すことができる。
図29は、上記(40)式、及び(41)式にて表される空気過剰率λと定常混合気温度上昇量ΔTmixsteadyとの関係を示したグラフである。なお、混合気到達距離Xが図28に示すように空気過剰率λと線形関係にあることを考慮すると、図29は、混合気到達距離Xと定常混合気温度上昇量ΔTmixsteadyとの関係を示したグラフであるということもできる。図29に示すように、定常混合気温度上昇量ΔTmixsteadyは、λ=1(即ち、X=X0。図28を参照。)にて最大値ΔTmix0=Hf/(stoich・Cg+Cf)を採る。
<定常燃料濃度[Fuel]mixsteady>
上述したように、定常火炎のリーン領域内では、定常燃料濃度[Fuel]mixsteadyは、空気過剰率λ(従って、混合気到達距離X)にかかわらず「0」に維持される。一方、定常火炎のリッチ領域内(λ<1)では、混合気内に残存している燃料量(質量)は、上記(39)式を考慮すると、(Q−(G/stoich))と表すことができる。従って、混合気質量(Q+G)に対する「混合気内に残存している燃料量(質量)」の割合である定常燃料濃度[Fuel]mixsteadyは、下記(42)式にて表すことができる。
従って、上記(37)式、及び(42)式からQ,Gを消去して整理すると、定常火炎のリッチ領域内(λ<1)における定常燃料濃度[Fuel]mixsteadyは、空気過剰率λの関数として下記(43)式にて表すことができる。
図30は、上記(43)式にて表される空気過剰率λと定常燃料濃度[Fuel]mixsteadyとの関係を示したグラフである。なお、混合気到達距離Xが空気過剰率λと線形関係にあることを考慮すると、図30は、混合気到達距離Xと定常燃料濃度[Fuel]mixsteadyとの関係を示したグラフであるということもできる。
<定常酸素濃度[O2]mixsteady>
上述したように、定常火炎のリッチ領域内では、定常酸素濃度[O2]mixsteadyは、空気過剰率λ(従って、混合気到達距離X)にかかわらず「0」に維持される。一方、定常火炎のリーン領域内(λ>1。実際には、λ≧1)では、混合気内に残存している酸素量(質量)は、上記「Qreac=Q」の関係式を考慮すると、(G−(Q・stoich))・[O2]in
と表すことができる。ここで、[02]inは吸気酸素濃度(従って、筒内ガス内の酸素濃度)である。従って、混合気質量(Q+G)に対する「混合気内に残存している酸素量(質量)」の割合である定常酸素濃度[02]mixsteadyは、下記(44)式にて表すことができる。
従って、上記(37)式、及び(44)式からQ,Gを消去して整理すると、定常火炎のリーン領域内(λ≧1)における定常酸素濃度[O2]mixsteadyは、空気過剰率λの関数として下記(45)式にて表すことができる。
図31は、上記(45)式にて表される空気過剰率λと定常酸素濃度[O2]mixsteadyとの関係を示したグラフである。なお、混合気到達距離Xが空気過剰率λと線形関係にあることを考慮すると、図31は、混合気到達距離Xと定常酸素濃度[O2]mixsteadyとの関係を示したグラフであるということもできる。
<定常窒素濃度[N2]mixsteady>
上述したように、筒内ガス中の窒素は不活性ガスであるから混合気内で化学反応により消費されることがない。従って、混合気内に残存している窒素量(質量)は、G・[N2]in と表すことができる。ここで、[N2]inは吸気窒素濃度(従って、筒内ガス内の窒素濃度)である。従って、混合気質量(Q+G)に対する「混合気内に残存している窒素量(質量)」の割合である定常窒素濃度[N2]mixsteadyは、下記(46)式にて表すことができる。
従って、上記(37)式、及び(46)式からQ,Gを消去して整理すると、定常窒素濃度[N2]mixsteadyは、空気過剰率λの関数として下記(47)式にて表すことができる。
図32は、上記(47)式にて表される空気過剰率λと定常窒素濃度[N2]mixsteadyとの関係を示したグラフである。なお、混合気到達距離Xが空気過剰率λと線形関係にあることを考慮すると、図32は、混合気到達距離Xと定常窒素濃度[N2]mixsteadyとの関係を示したグラフであるということもできる。
以上説明したように、定常火炎が占める領域内(即ち、混合気内)においては、定常混合気温度上昇量ΔTmixsteady、定常燃料濃度[Fuel]mixsteady、定常酸素濃度[O2]mixsteady、及び定常窒素濃度[N2]mixsteadyを、空気過剰率λの関数でそれぞれ表すことができる。
他方、空気過剰率λは、上記(2)式、(3)式により噴射後経過時間tの関数funcλ(t)として表すことができる。図33は、関数funcλ(t)から得られる噴射後経過時間tと空気過剰率λとの関係を示したグラフである。なお、t=t0のとき、λ=1となるものとする。
以上のことから、図29〜図32に示したそれぞれの関係に対して上記関数funcλ(t)から得られる図33に示した関係を適用することにより、定常火炎が占める領域内(即ち、混合気内)においては、定常混合気温度上昇量ΔTmixsteady、定常燃料濃度[Fuel]mixsteady、定常酸素濃度[O2]mixsteady、及び定常窒素濃度[N2]mixsteadyをそれぞれ、噴射後経過時間tの関数である、funcΔTmixsteady(t),func[Fuel]mixsteady(t),func[02]mixsteady(t),func[N2]mixsteady(t)で表すことができる。
図34〜図37は、funcΔTmixsteady(t),func[Fuel]mixsteady(t),func[02]mixsteady(t),func[N2]mixsteady(t)で表される関係をそれぞれ示したグラフである。第4実施形態は、係る噴射後経過時間tの関数であるfuncΔTmixsteady(t),func[Fuel]mixsteady(t),func[02]mixsteady(t),func[N2]mixsteady(t)を利用して、定常火炎が占める領域内(即ち、混合気内)におけるスート発生速度d[Soot]mix/dt、及びNO発生速度d[NO]mix/dtを求める。
(第4実施形態の実際の作動)
以下、第4実施形態に係るエミッション発生量推定装置の実際の作動について説明する。この装置のCPU61は、第1実施形態のCPU61が実行する図5〜図9に示した一連のルーチン、及び図13に示したルーチンのうち、図13に示したルーチンのみをそのまま実行するとともに、図5〜図9に示した一連のルーチンに代えて図38〜図41にフローチャートにより示した一連のルーチンを実行する。
なお、図38〜図41に示したルーチンのステップであって図5〜図9に示したルーチンのステップと同一のものについては、図5〜図9に示したルーチンの対応するステップの番号と同一の番号を付すことでそれらの説明を省略する。以下、第4実施形態に特有の図38〜図41に示した各ルーチンについて説明する。なお、図38のルーチンは図5のルーチンに、図39のルーチンは図6のルーチンに、図40のルーチンは図7、及び図8のルーチンに、図41のルーチンは図9のルーチンに対応している。
第4実施形態のCPU61は、図5〜図9に示した一連のルーチンに対応する図38〜図41に示した一連のルーチンを所定時間の経過毎に繰り返し実行するようになっている。従って、所定のタイミングになると、CPU61は図38のステップ3800から処理を開始し、ステップ505に進んで、「Yes」と判定する場合(即ち、IVCが到来した場合)、ステップ510〜530の処理を順に実行し、その後、ステップ550〜560の処理を順に実行する。
次に、CPU61は、図39のルーチンに進んで各種関数の決定を行うための処理を行う。具体的には、CPU61はステップ3905に進んで、定常火炎が発生している間の筒内ガス密度ρgの平均値である定常火炎時平均筒内ガス密度ρgaveを求める。定常火炎時平均筒内ガス密度ρgaveは、筒内ガスの全質量Mgに依存して決定されるから、図38のステップ515にて求めた筒内ガスの全質量Mgと、Mgを引数とするρgaveを求める関数funcρgaveとから求められる。
続いて、CPU61はステップ3910に進み、定常火炎が発生している間の筒内ガス圧力Pgの平均値である定常火炎時平均筒内ガス圧力Pgaveを求める。定常火炎時平均筒内ガス圧力Pgaveは、IVC時筒内ガス圧力Pgivc、IVC時クランク角度CAivc、指令燃料噴射量Qfinに依存して決定されるから、図38のステップ510にて求めたIVC時筒内ガス圧力Pgivc及びIVC時クランク角度CAivcと、図38のステップ520にて求めた指令燃料噴射量Qfinと、Pgivc,CAivc,Qfinを引数とするPgaveを求める関数funcPgaveとから求められる。
次いで、CPU61はステップ3915に進んで、図38のステップ530にて求めた基本燃料噴射圧力Pcrbaseから上記定常火炎時平均筒内ガス圧力Pgaveを減じることで定常火炎時平均有効噴射圧力ΔPaveを求め、続くステップ3920にて、上記求めた定常火炎時平均有効噴射圧力ΔPaveと、上記求めた定常火炎時平均筒内ガス密度ρgaveと、上記テーブルMapθとに基づいて定常火炎時平均噴霧角θaveを求める。これにより、定常火炎時平均噴霧角θaveは、ΔPave及びρgaveに基づいて決定される。
続いて、CPU61はステップ3925に進み、上記(3)式の有効噴射圧力ΔP、筒内ガス密度ρg、及び噴霧角θとして、上記求めた定常火炎時平均有効噴射圧力ΔPave、定常火炎時平均筒内ガス密度ρgave、及び定常火炎時平均噴霧角θaveをそれぞれ使用することで、上記(2)式、(3)式を利用して、噴射後経過時間tと空気過剰率λとの関係を規定する上記関数funcλ(t)(図33を参照)を決定する。
次に、CPU61はステップ3930に進んで、上記求めた関数funcλ(t)と、上記(40)式、及び(41)式にて表される空気過剰率λを引数とするΔTmixsteadyを求めるための関数funcΔTmixsteady(λ)とから、噴射後経過時間tを引数とする定常混合気温度上昇量ΔTmixsteadyを求めるための関数funcΔTmixsteady(t)を決定する。
続いて、CPU61はステップ3935に進んで、上記求めた関数funcλ(t)と、上記(43)式にて表される空気過剰率λを引数とする[Fuel]mixsteadyを求めるための関数func[Fuel]mixsteady(λ)とから、噴射後経過時間tを引数とする定常燃料濃度[Fuel]mixsteadyを求めるための関数func[Fuel]mixsteady(t)を決定する。
次いで、CPU61はステップ3940に進んで、上記求めた関数funcλ(t)と、上記(45)式にて表される空気過剰率λを引数とする[O2]mixsteadyを求めるための関数func[O2]mixsteady(λ)とから、噴射後経過時間tを引数とする定常酸素濃度[O2]mixsteadyを求めるための関数func[O2]mixsteady(t)を決定する。
次に、CPU61はステップ3945に進んで、上記求めた関数funcλ(t)と、上記(47)式にて表される空気過剰率λを引数とする[N2]mixsteadyを求めるための関数func[N2]mixsteady(λ)とから、噴射後経過時間tを引数とする定常窒素濃度[N2]mixsteadyを求めるための関数func[N2]mixsteady(t)を決定する。
続いて、CPU61はステップ3950に進んで、図6のステップ660と同様、混合気内のNO濃度[NO]mix、及びスート濃度[Soot]mixをそれぞれ初期値「0」に設定し、続くステップ3955にて、図6のステップ665と同様、噴射後経過時間tを初期値「0」に設定する。
次に、CPU61はステップ670の処理を実行し、その後、図40のルーチンに進み、定常火炎が占める領域内(従って、混合気内)における混合気温度、各種濃度等の算出のための処理を開始する。具体的には、CPU61は先ずステップ4002に進み、噴射後経過時間tの値(初期値は、図39のステップ3955の処理により「0」)を微小時間Δt(例えば、0.1msec)だけ進め、続くステップ4004にて、現時点での噴射後経過時間tの値と、図39のステップ3925にて決定した関数funcλ(t)とに基づいて噴射後経過時間tに対応する空気過剰率λを求める。
次いで、CPU61はステップ4006に進み、現時点での噴射後経過時間tの値と、図39のステップ3930にて決定した関数funcΔTmixsteady(t)とに基づいて噴射後経過時間tに対応する定常混合気温度上昇量ΔTmixsteadyを求め、同求めたΔTmixsteadyに図38のステップ550にて求めた燃料蒸気温度Tfを加えて、噴射後経過時間tに対応する定常混合気温度Tmixsteadyを求める。
続いて、CPU61はステップ4008に進んで、現時点での噴射後経過時間tの値と、図39のステップ3935にて決定した関数func[Fuel]mixsteady(t)とに基づいて噴射後経過時間tに対応する定常燃料濃度[Fuel]mixsteadyを求める。
次に、CPU61はステップ4010に進んで、現時点での噴射後経過時間tの値と、図39のステップ3940にて決定した関数func[O2]mixsteady(t)とに基づいて噴射後経過時間tに対応する定常酸素濃度[O2]mixsteadyを求める。
次いで、CPU61はステップ4012に進んで、現時点での噴射後経過時間tの値と、図39のステップ3945にて決定した関数func[N2]mixsteady(t)とに基づいて噴射後経過時間tに対応する定常窒素濃度[N2]mixsteadyを求める。
次に、CPU61はステップ835の判定を行う。いま、「Yes」と判定するものとすると、CPU61はステップ4016に進み、ステップ4004にて更新された現時点での噴射後経過時間tに対応する空気過剰率λが「1」より小さいか否か(即ち、リッチ領域か否か)を判定する。
ここで、「Yes」と判定する場合(即ち、リッチ領域である場合)、CPU61はステップ4018に進んで、ステップ4008にて求めた現時点での噴射後経過時間tに対応する定常燃料濃度[Fuel]mixsteadyと、図39のステップ3910にて求めた定常火炎時平均筒内ガス圧力Pgaveと、ステップ4006にて求めた現時点での噴射後経過時間tに対応する定常混合気温度Tmixsteadyと、上記(24)式とに基づいてスートの生成速度dmsf/dtを求める。続いて、CPU61はステップ4020に進み、上記(23)式から「dmso/dt」の項を削除したステップ4020内に記載の式と、上記求めたスートの生成速度dmsf/dtとから、現時点での噴射後経過時間tに対応するスート発生速度d[soot]mix/dtを求める。即ち、この場合、上記(25)式の計算が省略される。
一方、ステップ4016にて「No」と判定する場合(即ち、リーン領域である場合)、CPU61はステップ4022に進んで、その時点での混合気内のスート濃度[Soot]mix(初期値は、図39のステップ3950の処理により「0」)と、ステップ4010にて求めた現時点での噴射後経過時間tに対応する定常酸素濃度[O2]mixsteadyと、図39のステップ3910にて求めた定常火炎時平均筒内ガス圧力Pgaveと、ステップ4006にて求めた現時点での噴射後経過時間tに対応する定常混合気温度Tmixsteadyと、上記(25)式とに基づいてスートの酸化速度dmso/dtを求める。続いて、CPU61はステップ4024に進み、上記(23)式から「dmsf/dt」の項を削除したステップ4024内に記載の式と、上記求めたスートの酸化速度dmso/dtとから、現時点での噴射後経過時間tに対応するスート発生速度d[soot]mix/dtを求める。即ち、この場合、上記(24)式の計算が省略される。
上記ステップ4020、或いはステップ4024にて現時点での噴射後経過時間tに対応するスート発生速度d[soot]mix/dtが求められると、CPU61はステップ855の処理を実行し、その後、ステップ4028に進んで、現時点での噴射後経過時間tに対応する空気過剰率λが「1」より大きく、且つ、上記求めた現時点での噴射後経過時間tに対応する定常混合気温度Tmixsteadyが上記スートの反応限界温度TminSootより低いか否かを判定する。
そして、図8のステップ860の判定と同様、CPU61は、ステップ4028にて「Yes」と判定する場合、ステップ865の処理を実行してからステップ870に進み、ステップ4028にて「No」と判定する場合、ステップ870に直ちに進む。
CPU61は、ステップ870に進むと、同ステップ870の判定を行う。いま、「Yes」と判定するものとすると、CPU61はステップ4034に進み、ステップ4006にて求めた現時点での噴射後経過時間tに対応する定常混合気温度Tmixsteadyと、ステップ4010にて求めた現時点での噴射後経過時間tに対応する定常酸素濃度[O2]mixsteadyと、ステップ4012にて求めた現時点での噴射後経過時間tに対応する定常窒素濃度[N2]mixsteadyと、上記(27)式とに基づいて現時点での噴射後経過時間tに対応するNO発生速度d[NO]mix/dtを求める。
続いて、CPU61はステップ880の処理を実行し、その後、ステップ4038に進んで、現時点での噴射後経過時間tに対応する空気過剰率λが「1」より大きく、且つ、上記求めた現時点での噴射後経過時間tに対応する定常混合気温度Tmixsteadyが上記NOの反応限界温度TminNOより低いか否かを判定する。
そして、図8のステップ885の判定と同様、CPU61は、ステップ4038にて「Yes」と判定する場合、ステップ890の処理を実行してから図41のルーチンに進み、ステップ4038にて「No」と判定する場合、図41のルーチンに直ちに進む。
CPU61は図41のルーチンに進むと、図9のステップ905に対応するステップ4105の判定を行う。CPU61は、ステップ4105の判定において「No」と判定する毎に、図40のステップ4002〜図41のステップ4105の処理を繰り返し実行する。これにより、ステップ4105の判定において「No」と判定される毎に、図40のステップ4002にて、噴射後経過時間tがΔtだけ進められる。
即ち、ステップ4105の判定において「No」と判定される限りにおいて、空気過剰率λ、及び(フラグENDsoot=フラグENDno=0である限りにおいて)エミッション濃度([Soot]mix,[NO]mix)が、噴射後経過時間tの微小時間Δt毎の進行に対応して求められていく。
そして、ステップ4105の条件が成立すると、CPU61は図41のステップ4105に進んだとき「Yes」と判定してステップ4110以降に進み、エミッション発生量に係わる計算を終了するための処理を行う。
即ち、CPU61はステップ4110に進むと、現時点(即ち、ステップ4105の条件が成立した時点)での空気過剰率λと、図38のステップ520にて求めた指令燃料噴射量Qfinと、ステップ4110内に記載の式とに基づいて、指令燃料噴射量Qfinの燃料を含んだ定常火炎に係る混合気であって空気過剰率が値λとなる混合気の質量(定常混合気質量Mmixsteady)を求める。
次に、CPU61はステップ4115に進んで、図40のステップ855にて更新されている現時点での混合気内のスート濃度[Soot]mixの値に上記定常混合気質量Mmixsteadyを乗じることで定常火炎が占める領域内でのスート発生量を求め、同求めたスート発生量に初期値Soot0を加えることでスート総発生量Sootを求める。
同様に、CPU61はステップ4115にて、図40のステップ880にて更新されている現時点での混合気内のNO濃度[NO]mixの値に上記定常混合気質量Mmixsteadyを乗じることで定常火炎が占める領域内でのNO発生量を求め、同求めたNO発生量に初期値NO0を加えることでNO総発生量NOを求める。
そして、CPU61はステップ925〜960を順に実行し、図41のステップ4115にて求められているスート総発生量Sootの値、及びNO総発生量NOの値(従って、エミッション総発生量の推定結果)に基づいて噴射圧力を補正し、ステップ3895に進んで図38〜図41の一連の本ルーチンを一旦終了する。以降、CPU61は、次のIVCが到来するまでの間、図38のステップ505に進む毎に「No」と判定し続ける。
以上、説明したように、本発明によるエミッション発生量推定装置の第4実施形態は、混合気の着火後において燃焼室内で所謂定常火炎が発生している場合を想定している。そして、第4実施形態は、定常火炎のリッチ領域内(λ<1)では、定常状態における酸素濃度(定常酸素濃度[O2]mixsteady)が「0」になることを利用して、上記(23)式においてdmso/dtの項を省略することで得られる式である「d[Soot]mix/dt=dmsf/dt」に従ってスート発生速度d[Soot]mix/dtを求める。
同様に、第4実施形態は、定常火炎のリーン領域内(λ≧1)では、定常状態における燃料濃度(定常燃料濃度[Fuel]mixsteady)が「0」になることを利用して、上記(23)式においてdmsf/dtの項を省略することで得られる式である「d[Soot]mix/dt=−dmso/dt」に従ってスート発生速度d[Soot]mix/dtを求める。これらにより、(24)式、(25)式の算出回数を減らすことができ、この結果、CPU61の計算負荷を減らすことができる。
21…燃料噴射弁、60…電気制御装置、61…CPU、72…吸気温センサ、73…吸気管圧力センサ、74…クランクポジションセンサ、76…燃料温度センサ、77…筒内圧力センサ、78…吸気酸素濃度センサ