角速度和角加速度计算装置、转矩推算装置、燃烧状态推算装置

技术领域

[0001]

本发明是涉及角速度和角加速度计算装置、转矩推算装置、燃烧状态推算装置的发明，特别是涉及内燃机的角速度和角加速度计算装置、转矩推算装置、燃烧状态推算装置。

背景技术

[0002]

为了高精度地控制作为内燃机的发动机，采用测量发动机的燃烧状态、反馈该测量结果而进行控制的方法比较有效。为此，正确地测量发动机的燃烧状态很重要。作为测量发动机的燃烧状态的方法，具有由缸内压力传感器测量汽缸缸内压力而计算转矩的方法和由曲柄角度传感器信号推算气压转矩(ガス圧トルク)的方法。

[0003]

由曲柄角度传感器信号推算气体压力转矩的方法例如公开于专利文献1中。在专利文献1中，根据由曲柄角度传感器测量的曲柄角度，求出规定的曲柄角度范围中的曲轴的平均角加速度，由该平均角加速度、预先求出的绕曲轴的惯性矩、预先存储的摩擦转矩，推算气压转矩的平均值。

[0004]

[专利文献1]日本特开2005-61371号公报

[0005]

然而，在由缸内压力传感器测量汽缸缸内压力而计算转矩的方法中，耐久性、成本成为问题。在专利文献1公开的方法中，不计算曲轴的角加速度，而是计算规定的曲柄角度范围内的平均角加速度，仅能推算气体压力转矩的平均值。另外，在专利文献1公开的方法中，仅是由曲柄角度对角速度进行差分近似，由规定的2个曲柄角度的角速度差计算出平均角加速度，所以，不能计算出角加速度。

[0006]

另外，在专利文献1公开的方法中，关于对气体压力转矩的推算所需要的外力转矩，仅将作为该外力转矩的一部分的摩擦转矩作为与发动机转速和冷却水温相对应的映像预先存储并利用。为此，在专利文献1公开的方法中，仅在得知摩擦转矩以外的外力转矩为0的情况下或能够由倾斜传感器等其它测量装置推算负荷转矩的情况下，才能够测定气体压力转矩，适用范围受到限制。

发明内容

[0007]

因此，本发明的目的在于提供一种由曲柄角度传感器信号计算曲柄角速度和角加速度的角速度及角加速度计算装置、推算气体压力转矩和外力转矩的转矩推算装置、及推算燃烧状态的燃烧状态推算装置。

[0008]

作为解决手段，本发明的内燃机中的曲轴的角速度和角加速度计算装置、转矩推算装置、燃烧状态推算装置至少具有：曲柄角度检测单元，该曲柄角度检测单元检测内燃机的曲柄角度和处于该曲柄角度的时刻；曲柄角度修正单元，该曲柄角度修正单元根据规定的曲柄角度间隔，修正由曲柄角度检测单元检测出的曲柄角度和时刻；角速度和角加速度计算单元，该角速度和角加速度计算单元将与由曲柄角度修正单元修正了的曲柄角度对应的时刻的函数变换成可看作周期函数的形式，并对函数利用离散傅里叶变换进行运算，从而计算曲柄的角速度和角加速度。

[0009]

本发明中记载的角速度和角加速度计算装置、转矩推算装置、燃烧状态推算装置，将与由曲柄角度修正单元修正的曲柄角度对应的时刻的函数变换成可看作周期函数的形式，对该函数利用离散傅里叶变换进行运算，从而计算曲柄的角速度和角加速度，所以，能够提供由曲柄角度传感器信号计算出精度良好的曲柄角速度和角加速度的角速度和角加速度计算装置、推算气体压力转矩和外力转矩的转矩推算装置、及推算燃烧状态的燃烧状态推算装置。

附图说明

[0084]

图1为包含本发明实施方式的燃烧状态推算装置的内燃机的示意图。

图2为本发明实施方式的燃烧状态推算装置的框图。

图3为表示本发明实施方式的燃烧状态推算装置的动作的流程图。

图4为表示由本发明实施方式的燃烧状态推算装置计算的曲柄角加速度的图。

图5为用于说明在本发明实施方式的燃烧状态推算装置中使用的运动方程式的图。

图6为表示用本发明实施方式的燃烧状态推算装置推算的气体压力转矩的图。

具体实施方式

[0010]

(实施方式)

图1为包含本实施方式的燃烧状态推算装置的内燃机的示意图。图1所示的内燃机以4冲程直列4汽缸的发动机为例进行说明。首先，图1所示的发动机10在各汽缸内具有活塞20，在活塞20的上方形成燃烧室23。

[0011]

另外，在发动机10中，进气通道30和排气通道31相连接。吸入的空气通过进气通道30的上游端的空气滤清器32，经由气流传感器52和节气门33，进入稳压箱34。在节气门33的近旁设置检测节气门开度的油门传感器53。进入了稳压箱34的空气从进气通道30的下游端的进气门25被吸入到各汽缸的燃烧室23。在稳压箱34的下游配置喷射燃烧的燃料喷射阀27。

[0012]

被吸入到燃烧室23的空气和从燃料喷射阀27喷射的燃料由火花塞24进行火花点火，废气从排气门26排出。废气通过排气通道31，在由排气净化催化剂36净化后，排出到外部。从燃烧室23排出的废气的一部分通过EGR(废气再循环)阀35，被引导至进气通道30。在排气通道31中设置检测排气中的氧浓度的空燃比传感器54。

[0013]

对于各汽缸的活塞20，为了将其往复运动变换为旋转运动，通过连杆21连接于曲轴22。在曲轴22上，为了检测其旋转角度，设有曲柄角度传感器50。另外，在各汽缸上设置有检测冷却水温的水温传感器51。

[0014]

另外，发动机10具有ECU(发动机控制单元)11，向该ECU输入气流传感器52、油门传感器53、空燃比传感器54、曲柄角度传感器50、水温传感器51、及设置于加速踏板40的加速踏板传感器55等各传感器的输出，该ECU对应于运行状态控制吸入空气量、燃料喷射量、点火时刻等。此外，在ECU11中，曲柄角度检测方法、曲柄角度修正方法、角速度和角加速度计算方法、轴转矩计算方法、外力转矩推算方法、气体压力转矩计算方法、气体压力做功量计算方法、燃烧状态推算方法作为程序被保存并执行。

[0015]

将本实施方式的燃烧状态推算装置作为功能块表示的图表示在图2中。在图2所示的燃烧状态推算装置中，具有曲柄角度检测单元1、曲柄角度修正单元2、及角速度和角加速度计算单元3；该曲柄角度检测单元1根据从曲柄角度传感器50输入的信号，检测曲柄角度和与之对应的时刻；该曲柄角度修正单元2修正由曲柄角度检测单元1检测出的曲柄角度和时刻；该角速度和角加速度计算单元3根据曲柄角度修正单元2的修正结果，计算曲柄的角速度和角加速度。此外，由图2所示的曲柄角度检测单元1、曲柄角度修正单元2、及角速度和角加速度计算单元3构成角速度和角加速度计算装置4。

[0016]

另外，在图2所示的燃烧状态推算装置中，具有轴转矩计算单元5、外力转矩推算单元6、及气体压力转矩计算单元7；该轴转矩计算单元5根据角速度和角加速度计算单元3的计算结果，由包含活塞20、连杆21、及曲轴22的系统的运动方程式计算轴转矩；该外力转矩推算单元6根据特异曲柄角度处的轴转矩推算外力转矩；该气体压力转矩计算单元7根据轴转矩计算单元5和外力转矩推算单元6的结果，计算气体压力转矩。此外，通过组合图2所示的角速度和角加速度计算装置4、轴转矩计算单元5、外力转矩推算单元6、及气体压力转矩计算单元7，构成转矩推算装置8。

[0017]

另外，在图2所示的燃烧状态推算装置中，具有气体压力做功量计算单元9a和燃烧状态推算单元9b；该气体压力做功量计算单元9a根据气体压力转矩计算单元7的计算结果，计算气体压力做功量；该燃烧状态推算单元9b根据角速度和角加速度计算单元3、气体压力转矩计算单元7、及气体压力做功量计算单元9a各自的计算结果中的至少任意一个，推算内燃机的燃烧状态，输出其推算结果。

[0018]

下面，具体地说明图1和图2所示的本实施方式的燃烧状态推算装置实施曲柄角速度和角加速度的计算、从发动机输出的转矩的推算、发动机燃烧状态的推算的方法。具体地说，图3示出了表示本实施方式的燃烧状态推算装置的动作的流程图。

[0019]

首先，在图3所示的步骤S1中，曲柄角度检测单元1根据曲柄角度传感器50的信号，检测与曲柄角度θ_n(n＝1，2，...，N_sample)相对应的时刻t_n。N_sample为样本数，曲柄角度θ_n表示累积角度。在曲柄角度传感器50测量的转子上没有缺齿、误差也小的情况下，曲柄角度间隔θ_n+1-θ_n为固定值Δθ＝2π/N_pulse(N_pulse为包含缺齿的转子的齿数)。然而，在实际的转子上，由于存在缺齿或制造误差等，所以，预先使曲轴22转1圈，对于缺齿以外的齿评价并求出曲柄角度θ_n。即使在不能预先求出误差的情况下，也预先计算出与缺齿以外的齿对应的曲柄角度θ_n。

[0020]

然后，在图3所示的步骤S2中，曲柄角度修正单元2由曲柄角度θ_n和其时刻t_n，计算出曲柄角度间隔Δθ的曲柄角度<hat>θ_j(在文章中，将在θ_j上加了“^”的文字记为<hat>θ_j)的时刻<hat>t_j(在文章中，将在t_j上加了“^”的文字记为<hat>t_j)。时刻<hat>t_j的计算可使用例如公知的3次样条插补法。

[0021]

然后，在图3所示的步骤S3中，角速度和角加速度计算单元3由规定的曲柄角度范围的曲柄角度<hat>θ_j和其时刻<hat>t_j，计算出曲轴的曲柄角速度和角加速度。在发动机为4冲程发动机并且其汽缸数为N_cylinder的情况下，从发动机行程的周期性考虑，规定的曲柄角度范围最好为4π/N_cylinder(弧度)的整数倍。设规定的曲柄角度范围的样本数为N_scope，为了标注方便，将规定的曲柄角度范围的曲柄角度<hat>θ_j和其时刻<hat>t_j的下标设为j＝1，2，...，N_scope。但是，在以后的运算中，也可利用j＝N_scope+1的曲柄角度<hat>θ_j和其时刻<hat>t_j。

[0022]

下面，具体地说明计算曲轴的曲柄角速度和角加速度的方法。首先，在以曲柄角度为横轴、时刻为纵轴的情况下，考虑通过曲柄角度<hat>θ_j和其时刻<hat>t_j的点的函数t(θ)，使用数学式1的α定义数2。此外，从发动机行程的周期性考虑，数学式2可看作是以规定的曲柄角度范围为周期的周期函数。

[0023]

[数学式1]

[0024]

[数学式2]

x(θ)＝t(θ)-αθ

[0025]

若用θ对数学式2进行微分并整理，则作为曲柄角度θ的函数，可求出角速度<dot>θ(θ)(在文章中，将在θ上加了“·”的文字表示为<dot>θ)。

[0026]

[数学式3]

$\stackrel{·}{\theta }\left(\theta \right)=\frac{\mathrm{d\theta }}{\mathrm{dt}}\left(\theta \right)=\frac{1}{\frac{\mathrm{dx}\left(\theta \right)}{\mathrm{d\theta }}+\alpha }$

[0027]

在实际的计算中，由x(θ)的样本点x_j＝<hat>t_j-α<hat>θ_j(j＝1，2，...，N_scope)求出dx(<hat>θ_j)/dθ，并计算<dot>θ(<hat>θ_j)。此外，在设x_j的离散傅里叶变换为X(ω)(|ω|≤π/Δθ)时，dx(<hat>θ_j)/dθ通过对iωX(ω)进行逆离散傅里叶变换而求出。i为虚数单位。

[0028]

曲柄角度θ处的角加速度<2dot>θ(θ)(在文章中，将在θ上加了“¨”的文字表示为<2dot>θ)可由数学式4求出。此外，d²t(<hat>θ_j)/dθ²通过对(iω)²X(ω)进行逆离散傅里叶变换而求出。

[0029]

[数学式4]

$\stackrel{··}{\theta }\left(\theta \right)=\frac{d^2\theta }{{\mathrm{dt}}^2}\left(\theta \right)=-\stackrel{·}{\theta }{\left(\theta \right)}^3\frac{d^{2}t\left(\theta \right)}{{\mathrm{d\theta }}^2}$

[0030]

在使用数学式3和数学式4计算曲柄角速度和角加速度的情况下，也可对X(ω)的高频成分进行频带限制。在对X(ω)的比规定的频率高的频率成分(高频成分)进行频带限制的情况下，能够减小曲柄角度传感器50的误差和噪声成分。在这里，规定的频率是指，不对曲柄角速度和角加速度的计算产生影响的频率成分中的较高一侧的频率。

[0031]

另外，在上述计算方法中，为了计算规定的曲柄角度范围的曲柄角速度和角加速度，使用数学式2的周期函数，但本发明并不仅限定于此，也可利用窗函数计算。具体地说，例如在计算j＝1，2，...，N_scope的曲柄角速度和角加速度时，利用j＝-M+1，...，N_scope+M中的<hat>θ_j，<hat>t_j。M为适当的自然数。

[0032]

然后，使用数学式5的x′(θ)代替x(θ)，用由数学式6定义的窗函数乘以数学式5，由数学式3和数学式4计算出j＝1，2，...，N_scope的曲柄角速度和角加速度。此外，数学式6的a和d为适当的常数，例如可以设a＝10，d＝4。通过利用数学式6所示的窗函数计算出曲柄角速度和角加速度，即使在发动机的转速急剧变化的情况下，也能够正确地计算曲柄角速度和角加速度。另外，利用窗函数计算曲柄角速度和角加速度的方法并不限于数学式5和数学式6。

[0033]

[数学式5]

[0034]

[数学式6]

$w_j=\left(\begin{array}{cc}\exp (-a{\left|\frac{-j+1}{M}\right|}^d)& \mathrm{ifj}\in \{-M+1,···,0\}\\ 1& \mathrm{ifj}\in \{1,···,N_{\mathrm{scope}}\}\\ \exp (-a{\left|\frac{j-N_{\mathrm{scope}}}{M}\right|}^d)& \mathrm{ifj}\in \{N_{\mathrm{scope}}+1,···,N_{\mathrm{scope}}+M\}\end{array}\right)$

[0035]

由角速度和角加速度计算单元3计算出的曲柄角加速度的例子表示于图4的曲线图。图4为重叠绘出跨越多个点火循环的曲柄角加速度的曲线图，横轴表示曲柄角度，纵轴表示曲柄角加速度。由图4可以得知，在曲柄角度为40度(degree)附近，存在由燃烧形成的曲柄角加速度的峰值。

[0036]

然后，在图3所示的步骤S4中，轴转矩计算单元5由包含活塞20、连杆21、及曲轴22的系统的运动方程式计算出轴转矩。这样，能够以良好精度推算曲轴22周围的轴转矩。

[0037]

下面，具体地说明在轴转矩计算单元5中使用的运动方程式。首先，为了使问题简单，考虑1个汽缸的运动方程式。图6表示活塞20、连杆21及曲轴22的示意图。在图6中，曲轴22水平地设置，在其垂直线上存在活塞20的上止点。另外，在图6中，将水平面内的与曲轴22直交的方向设为x轴，将垂直方向设为y轴。另外，在图6中，将曲柄角度表示为θ，将连杆的角度表示为φ。

[0038]

另外，设活塞20的质量为m_p，连杆21的长度为l，质量为m_l，惯性矩为I_l，重心为G_l，从重心G_l到曲柄销28的长度为ρ_ll。另外，设曲轴22的曲柄半径为r，质量为m_r，惯性矩为I_r，重心为G_r，从曲轴22的中心到曲柄重心G_r的长度为ρ_rr。另外，λ＝l/r，重力加速度为g。

[0039]

设包含图6所示的活塞20、连杆21及曲轴22的系统整体的动能为K、势能为U，并将它们算出，由拉格朗日算子L＝K-U计算数学式7所示的拉格朗日运动方程式，从而获得该系统的运动方程式。其中，Q表示气体压力转矩，F表示摩擦转矩等外力转矩。

[0040]

[数学式7]

$\frac{d}{\mathrm{dt}}\left(\frac{\partial L}{\partial \stackrel{·}{\theta }}\right)-\frac{\partial L}{\partial \theta }=Q+F$

[0041]

在设系统整体的绕曲轴22的惯性矩为I(θ)、重力产生的势能为U_g(θ)时，系统整体的动能K、势能U及运动方程式分别可表示为数学式8、数学式9、数学式10。

[0042]

[数学式8]

$K=\frac{1}{2}I\left(\theta \right){\stackrel{·}{\theta }}^2$

[0043]

[数学式9]

U＝U_g(θ)

[0044]

[数学式10]

$I\left(\theta \right)\stackrel{··}{\theta }+\frac{1}{2}{I}^{\prime }\left(\theta \right){\stackrel{·}{\theta }}^2+G\left(\theta \right)=Q+F$

[0045]

其中，$G\left(\theta \right)=\partial U_g\left(\theta \right)/\partial \left(\theta \right),$I′(θ)＝dI(θ)/dθ。

[0046]

数学式10所示的运动方程式的I(θ)、I′(θ)、G(θ)分别整理为数学式11、数学式12、数学式13。

[0047]

[数学式11]

$I\left(\theta \right)=I_r+m_r{\left({\rho }_{r}r\right)}^2+I_l{\left(\frac{\cos \theta }{\lambda \cos \phi }\right)}^2+m_l{r}^2{(1-{\rho }_l)}^2{\cos }^2\theta$

$+m_p{r}^2{(\sin \theta +\tan \phi \cos \theta )}^2+m_l{r}^2{(\sin \theta +{\rho }_l\tan \phi \cos \theta )}^2$

[0048]

[数学式12]

$I^{\prime }\left(\theta \right)=2I_l(\frac{{\cos }^2\theta }{{\lambda }^2{\cos }^2\phi }-1)\frac{\sin \theta \cos \theta }{{\lambda }^2{\cos }^2\phi }-2m_l{r}^2{(1-{\rho }_l)}^2\sin \theta \cos \theta$

$+2m_p{r}^2(\sin \theta +\tan \phi \cos \theta )(\cos \theta -\tan \phi \sin \theta +\frac{{\cos }^2\theta }{\lambda {\cos }^3\phi })$

$+2m_l{r}^2(\sin \theta +{\rho }_l\tan \phi \cos \theta )(\cos \theta -{\rho }_l\tan \phi \sin \theta +{\rho }_l\frac{{\cos }^2\theta }{\lambda {\cos }^3\phi })$

[0049]

[数学式13]

G(θ)＝-gr{(m_p+m_l+ρ_rm_r)sinθ+(m_p+ρ_lm_l)tanφcosθ}

[0050]

其中，φ＝sin^-1((1/λ)sinθ)。

[0051]

另外，气体压力转矩Q用数学式14表示。

[0052]

[数学式14]

Q＝q(t)r(sinθ+tanφcosθ)

[0053]

其中，q(t)为通过气体压力推压活塞20的力。

[0054]

以曲轴22为刚体的4冲程直列4汽缸发动机的I(θ)、I′(θ)、G(θ)、Q可利用数学式11到数学式14如数学式15到数学式18那样导出。

[0055]

[数学式15]

$I\left(\theta \right)=I_w+4I_r+4I_l{\left(\frac{\cos \theta }{\lambda \cos \phi }\right)}^2+4m_r{\left({\rho }_{r}r\right)}^2+4m_l{r}^2{(1-{\rho }_l)}^2{\cos }^2\theta$

$+4m_p{r}^2({\sin }^2\theta +{\tan }^2\phi {\cos }^2\theta )+4m_l{r}^2({\sin }^2\theta +{{\rho }_l}^2{\tan }^2\phi {\cos }^2\theta )$

[0056]

[数学式16]

$I^{\prime }\left(\theta \right)=4r^2\{\frac{I_l}{l^2{\cos }^2\phi }(\frac{{\cos }^2\theta }{{\lambda }^2{\cos }^2\phi }-1)+m_p(\frac{{\cos }^2\theta }{{\lambda }^2{\cos }^4\phi }-{\tan }^2\phi +1)$

$+m_l{{\rho }_l}^2(\frac{{\cos }^2\theta }{{\lambda }^2{\cos }^4\phi }-{\tan }^2\phi +\frac{2}{{\rho }_l}-1)\}\sin 2\theta$

[0057]

[数学式17]

G(θ)＝-4gr(m_p+ρ_lm_l)tanφcosθ

[0058]

[数学式18]

Q(θ)＝r(q₁(t)+q₄(t))(sinθ+tanφcosθ)

+r(q₂(t)+q₃(t))(-sinθ+tanφcosθ)

[0059]

其中，I_w为连接在曲轴22上的飞轮的惯性矩，q_k(t)为通过气体压力推压第k汽缸的活塞20的力。

[0060]

上述运动方程式是一个例子，在轴转矩计算单元5中使用的运动方程式并不限定于上述运动方程式。例如，可考虑将图6所示的连杆21近似为仅在曲柄侧端部(大端部)和活塞侧端部(小端部)的2点存在质量的2质量点模型的运动方程式。在将运动方程式近似为2质点模型的情况下，数学式15到数学式17可简化为数学式19到数学式21所示的式子。

[0061]

[数学式19]

I(θ)＝{I_w+4Ir+4m_r(ρ_rr)²+4m_l1r²}

+4r²(mp+m_l2)(sin²θ+tan²φcos²θ)

[0062]

[数学式20]

$I^{\prime }\left(\theta \right)=4r^2(m_p+m_{l2})(\frac{{\cos }^2\theta }{{\lambda }^2{\cos }^4\phi }-{\tan }^2\phi +1)\sin 2\theta$

[0063]

[数学式21]

G(θ)＝-4gr(m_p+m_l2)tanφcosθ

[0064]

其中，m_l1表示曲柄侧端部(大端部)的质量，m_l2表示活塞侧端部(小端部)的质量。

[0065]

在轴转矩计算单元5中，将作为曲柄角度、曲柄角速度和角加速度的函数的数学式10的左边设为轴转矩T_r(θ，<dot>θ，<2dot>θ)，使用由角速度和角加速度计算单元3计算出的曲柄角速度和角加速度，可计算出轴转矩T_r(θ，<dot>θ，<2dot>θ)。

[0066]

下面，在图3所示的步骤S5中，在外力转矩推算单元6中，由特异曲柄角度的轴转矩推算外力转矩。这样，能够以良好的精度推算摩擦转矩、运动负荷等外力转矩。

[0067]

下面，具体地说明外力转矩的推算方法。在外力转矩推算单元6中，利用了这样的情况，即，在作为活塞的上下止点的曲柄角度θ＝0，π方面，数学式16、数学式17、数学式18式成为0，数学式10的运动方程式成为数学式22。

[0068]

[数学式22]

$I\left(\theta \right)\stackrel{··}{\theta }=F$

[0069]

然后，外力转矩推算单元6通过在数学式22中使用由角速度和角加速度计算单元3计算出的曲柄角加速度，从而推算外力转矩F。作为外力转矩F的推算方法，可假设外力转矩F为一定，由多个活塞20的上下止点的曲柄角加速度计算外力转矩，求出这些外力转矩的平均值作为外力转矩F，也可对由活塞20的上下止点的曲柄角加速度计算出的外力转矩使用公知的卡尔曼滤波等，将由此获得的值作为外力转矩F。另外，作为外力转矩F的其它的推断方法，也可由特定汽缸的活塞20处于上止点的曲柄角加速度计算外力转矩，将该外力转矩的值作为外力转矩F。

[0070]

外力转矩推算单元6能够适用的发动机并不仅限于4冲程4汽缸的发动机。具体地说，设关于任意汽缸数的发动机的拉格朗日算子L为曲柄角度θ和曲柄角速度<dot>θ的函数L(θ，<dot>θ)。当对于某一曲柄角度θ^＊的微小变量，拉格朗日算子L(θ，<dot>θ)具有对称性$L({\theta }^{*}+\partial \theta ,<\mathrm{dot}>\theta )=L({\theta }^{*}-\partial \theta ,<\mathrm{dot}>\theta )$时，绕曲轴22的惯性矩I(θ)和重力势能U_g(θ)都具有对称性$I({\theta }^{*}+\partial \theta )=I({\theta }^{*}-\partial \theta ),$$U_g({\theta }^{*}+\partial \theta )=U_g({\theta }^{*}-\partial \theta ).$这在成为对象的系统(在上述例子中，包含活塞20、连杆21及曲轴22的系统)的机构相对于从到的变换对称时明显得到满足。当在曲柄角度θ^＊方面能够满足对称性条件时，绕曲轴22的惯性矩I(θ)和重力势能U_g(θ)的与θ相对应的导出函数I′(θ)和G(θ)分别在曲柄角度θ^＊方面成为0。

[0071]

4冲程1汽缸的发动机中的曲柄角度θ^＊＝0，π(活塞20的上下止点)。另外，在活塞20的上下止点，气体压力转矩Q也成为0。因此，上下止点的运动方程式成为数学式22。

[0072]

4冲程直列4汽缸的发动机也与4冲程1汽缸的发动机相同。这样，运动方程式记述为惯性矩I(θ)、曲柄角加速度<2dot>θ、及外力转矩F的关系式的曲柄角度θ^＊被称为特异曲柄角度。由外力转矩推算单元6，在该特异曲柄角度θ^＊方面利用数学式22推算外力转矩F。

[0073]

另外，即使在数学式22不严格地成立的情况下，当数学式22近似地成立时，也能够推算外力转矩F。例如，在4冲程直列6汽缸的发动机的情况下，成为各活塞20的上止点的曲柄角度θ＝0、±2π/3，I′(θ)＝0，G(θ)＝0，气体压力转矩Q也成为足够小的值，所以，数学式22近似地成立。数学式22近似地成立的曲柄角度在本发明的外力转矩推算单元6中也称为特异曲柄角度。

[0074]

然后，在图3所示的步骤S6中，气体压力转矩计算单元7根据由轴转矩计算单元5计算出的轴转矩T_r(θ，<dot>θ，<2dot>θ)以及由外力转矩推算单元6推算的外力转矩<hat>F计算气体压力转矩<hat>Q。具体地说，气体压力转矩<hat>Q通过从轴转矩T_r(θ，<dot>θ，<2dot>θ)减去推算的外力转矩<hat>F而求出。

[0075]

由气体压力转矩计算单元7计算出的气体压力转矩的一个例子表示在图5的曲线图中。在图5中，横轴表示曲柄角度，纵轴表示转矩。另外，在图5中，在表示有由气体压力转矩计算单元7计算出的气体压力转矩Q1的同时，为了比较，还表示有由与曲柄角度一起测量的缸内压力计算出的气体压力转矩Q2的值。从图5可以看出，由缸内压力计算出的气体压力转矩Q2与由气体压力转矩计算单元7计算出的气体压力转矩Q1精度良好地相符。

[0076]

然后，在图3所示的步骤S7中，气体压力做功量计算单元9a关于规定的曲柄角度范围的曲柄角度θ对由气体压力转矩计算单元7计算出的气体压力转矩<hat>Q进行积分，从而计算出气体压力做功量W_Q。这样，能够以良好的精度计算出反映了规定的曲柄角度范围内的燃烧状态的气体压力做功量W_Q。

[0077]

然后，在图3所示的步骤S8中，燃烧状态推算单元9b根据包含以下量中的至少任意一个的值推算内燃机的燃烧状态，所述量是，由角速度和角加速度计算单元3计算出的曲柄角加速度<2dot>θ、由气体压力转矩计算单元7计算出的气体压力转矩<hat>Q、由气体压力做功量计算单元9a计算出的气体压力做功量W_Q。

[0078]

具体地说，燃烧状态推算单元9b能够由气体压力做功量W_Q推算内燃机是否由燃烧产生必要的做功量，是否没有不发火。另外，燃烧状态推算单元9b比较各汽缸的气体压力做功量W_Q，从而也能够推算汽缸间偏差。另外，燃烧状态推算单元9b通过计算出气体压力做功量W_Q的COV值(将气体压力做功量W_Q的标准偏差除以平均值后获得的值)，也能够推算内燃机的燃烧偏差的大小。

[0079]

另外，燃烧状态推算单元9b通过由气体压力转矩<hat>Q计算气体压力峰值位置及其大小等，能够推算点火定时是否适当，是否没有不发火，是否没有发生异常燃烧(爆震、提前点火)。另外，燃烧状态推算单元9b通过计算气体压力转矩<hat>Q的变动量，也能够推算异常燃烧。

[0080]

另外，燃烧状态推算单元9b利用曲柄角加速度<2dot>θ代替气体压力转矩<hat>Q，从而能够简单地推算燃烧状态。

[0081]

如上述那样，燃烧状态推算单元9b根据曲柄角加速度<2dot>θ、气体压力转矩<hat>Q、气体压力做功量W_Q中的至少任意一个的值，或将它们组合起来的值，能够把握内燃机的状态。

[0082]

如上述那样，本实施方式的燃烧状态推算装置能够通过具有曲柄角度检测单元1、曲柄角度修正单元2、及角速度和角加速度计算单元3的角速度和角加速度计算装置4计算出高精度的曲柄角速度和角加速度。此外，在本实施方式的燃烧状态推算装置中，不使用缸内压力传感器，而是使用曲柄角度传感器，以良好的精度计算曲柄角速度和角加速度，所以，耐久性高，寿命长。

[0083]

另外，本实施方式的燃烧状态推算装置能够通过具有角速度和角加速度计算装置4、轴转矩计算单元5、外力转矩推算单元6、及气体压力转矩计算单元7的转矩推算装置8推算高精度的轴转矩T_r(θ，<dot>θ，<2dot>θ)、外力转矩F、及气体压力转矩Q。另外，本实施方式的燃烧状态推算装置具有转矩推算装置8、气体压力做功量计算单元9a、及燃烧状态推算单元9b，能够以良好精度把握内燃机的状态，能够实现内燃机的节能化、低废气排出化。