WO2015162971A1 - 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法 - Google Patents

内燃機関の熱発生率波形算出装置および熱発生率波形算出方法 Download PDF

Info

Publication number
WO2015162971A1
WO2015162971A1 PCT/JP2015/053497 JP2015053497W WO2015162971A1 WO 2015162971 A1 WO2015162971 A1 WO 2015162971A1 JP 2015053497 W JP2015053497 W JP 2015053497W WO 2015162971 A1 WO2015162971 A1 WO 2015162971A1
Authority
WO
WIPO (PCT)
Prior art keywords
ignition
heat generation
ignition timing
delay period
period
Prior art date
Application number
PCT/JP2015/053497
Other languages
English (en)
French (fr)
Inventor
宗矩 今枝
Original Assignee
トヨタ自動車株式会社
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by トヨタ自動車株式会社 filed Critical トヨタ自動車株式会社
Priority to US15/303,331 priority Critical patent/US9885295B2/en
Priority to JP2016514747A priority patent/JP6260692B2/ja
Priority to CN201580020997.XA priority patent/CN106232970B/zh
Priority to EP15782377.4A priority patent/EP3135887B1/en
Publication of WO2015162971A1 publication Critical patent/WO2015162971A1/ja

Links

Images

Classifications

    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02D—CONTROLLING COMBUSTION ENGINES
    • F02D35/00—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for
    • F02D35/02—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for on interior conditions
    • F02D35/023—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for on interior conditions by determining the cylinder pressure
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02D—CONTROLLING COMBUSTION ENGINES
    • F02D35/00—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for
    • F02D35/02—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for on interior conditions
    • F02D35/028—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for on interior conditions by determining the combustion timing or phasing
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02D—CONTROLLING COMBUSTION ENGINES
    • F02D41/00—Electrical control of supply of combustible mixture or its constituents
    • F02D41/02—Circuit arrangements for generating control signals
    • F02D41/14—Introducing closed-loop corrections
    • F02D41/1401—Introducing closed-loop corrections characterised by the control or regulation method
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02D—CONTROLLING COMBUSTION ENGINES
    • F02D45/00—Electrical control not provided for in groups F02D41/00 - F02D43/00
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02D—CONTROLLING COMBUSTION ENGINES
    • F02D41/00—Electrical control of supply of combustible mixture or its constituents
    • F02D41/02—Circuit arrangements for generating control signals
    • F02D41/14—Introducing closed-loop corrections
    • F02D41/1401—Introducing closed-loop corrections characterised by the control or regulation method
    • F02D2041/1413—Controller structures or design
    • F02D2041/1429—Linearisation, i.e. using a feedback law such that the system evolves as a linear one
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02D—CONTROLLING COMBUSTION ENGINES
    • F02D41/00—Electrical control of supply of combustible mixture or its constituents
    • F02D41/02—Circuit arrangements for generating control signals
    • F02D41/14—Introducing closed-loop corrections
    • F02D41/1401—Introducing closed-loop corrections characterised by the control or regulation method
    • F02D2041/1433—Introducing closed-loop corrections characterised by the control or regulation method using a model or simulation of the system
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02D—CONTROLLING COMBUSTION ENGINES
    • F02D35/00—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for
    • F02D35/02—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for on interior conditions
    • F02D35/023—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for on interior conditions by determining the cylinder pressure
    • F02D35/024—Controlling engines, dependent on conditions exterior or interior to engines, not otherwise provided for on interior conditions by determining the cylinder pressure using an estimation
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F02—COMBUSTION ENGINES; HOT-GAS OR COMBUSTION-PRODUCT ENGINE PLANTS
    • F02P—IGNITION, OTHER THAN COMPRESSION IGNITION, FOR INTERNAL-COMBUSTION ENGINES; TESTING OF IGNITION TIMING IN COMPRESSION-IGNITION ENGINES
    • F02P5/00—Advancing or retarding ignition; Control therefor
    • F02P5/04—Advancing or retarding ignition; Control therefor automatically, as a function of the working conditions of the engine or vehicle or of the atmospheric conditions
    • F02P5/145—Advancing or retarding ignition; Control therefor automatically, as a function of the working conditions of the engine or vehicle or of the atmospheric conditions using electrical means
    • F02P5/15—Digital data processing
    • F02P5/1502—Digital data processing using one central computing unit

Definitions

  • the present invention relates to a device for calculating a heat release rate waveform in a spark ignition type internal combustion engine and a method for calculating the same, and in particular, a period from when ignition by a spark plug to ignition of an air-fuel mixture ( In this specification, this period is referred to as “ignition delay period”), and it relates to a technique for obtaining a heat release rate waveform.
  • the heat generation rate in a cylinder is approximated by a Wiebe function.
  • the Weebe function can be used to express the heat generation rate waveform suitably by specifying a plurality of parameters, and is used to estimate the heat generation rate, mass combustion rate, etc. due to combustion of the internal combustion engine. ing.
  • the shape parameter m of the Weibbe function is identified by a predetermined formula based on the combustion ratio at the crank angle at which the heat generation rate is maximum. ing. Then, k, a / ⁇ p m + 1, and identified by respective predetermined formula for other parameters, such as theta b, to fit with high accuracy to the actual heat generation pattern, determining the Wibe function it can.
  • Patent Document 1 a plurality of parameters m, k, a / ⁇ p m + 1 , ⁇ b are identified as described above, and the operation of determining the Weebe function is performed for various operating conditions. It is described that the relationship between these parameters and the operating parameters of the internal combustion engine (load factor, rotational speed, air-fuel ratio, ignition timing, etc.) can be grasped. Further, it is described that if the relationship thus obtained is used, the Weibbe function can be determined for all operating conditions of the internal combustion engine, and the combustion state of the internal combustion engine can be accurately expressed.
  • the parameter m of Wibe function, k, a / ⁇ p m + 1 does not disclose a specific method for specifying the relationship between the operating parameters of the theta b and the internal combustion engine. Therefore, in practice for almost all operating conditions parameters m, k, a / ⁇ p m + 1, identified theta b, must be determined Wibe function for each operating condition. That is, the conventional method leaves room for further reducing the man-hours for creating the heat release rate waveform and reducing the cost.
  • the parameters m, k, a / ⁇ p m + 1 , ⁇ b are identified, and the entire heat release rate waveform can be expressed only after the Weebe function is determined.
  • the combustion state can be evaluated. Therefore, without expressing the entire heat generation rate waveform, for example, a period until the heat generation rate waveform rises after ignition is performed by the spark plug (a period from when ignition by the ignition plug is performed until the air-fuel mixture is ignited) ) Is not possible to estimate and evaluate only the ignition delay period.
  • the present invention has been made in view of such various points, and its object is to create a heat release rate waveform by focusing on the ignition delay period, which is one of the indexes representing the state of the air-fuel mixture in the cylinder.
  • the purpose is to reduce the man-hours for (calculation) and to easily estimate and evaluate, for example, the ignition delay period while ensuring the required accuracy.
  • the inventor of the present invention has a high correlation with the fuel density in the ignition delay period, which is a period from when ignition by the spark plug to ignition of the air-fuel mixture, and the engine load factor with respect to the ignition delay period and The new knowledge that the influence of the ignition timing can be expressed collectively by the fuel density was obtained.
  • the solution principle of the present invention is that the ignition delay period is used as one of the characteristic values of the heat release rate waveform, and the ignition delay period is estimated based on the fuel density. is there.
  • the present invention is directed to a device for calculating a heat release rate waveform in a spark ignition type internal combustion engine, and the period from when ignition by a spark plug to ignition of an air-fuel mixture is performed.
  • the ignition delay period which is one of the characteristic values of the heat release rate waveform, is defined.
  • the ignition delay period is estimated based on the fuel density in the cylinder at the ignition timing, while the air-fuel mixture is
  • the ignition delay period is estimated based on the fuel density in the cylinder at the ignition timing, and the estimated ignition delay period is used.
  • the heat generation rate waveform is calculated.
  • ignition is performed by the spark plug.
  • the ignition delay period until the air-fuel mixture ignites is used. This ignition delay period varies depending on the operating conditions such as the load factor and ignition timing of the internal combustion engine.
  • the engine load factor parameter that defines the fuel injection amount
  • the ignition timing cylinder volume
  • the ignition delay period can be estimated more easily than before while ensuring the required accuracy. And can be evaluated.
  • the ignition timing of the air-fuel mixture is on the advance side of the compression top dead center (assuming that it is on the advance side of the compression top dead center)
  • the ignition timing The ignition delay period is estimated based on the fuel density in the cylinder.
  • the ignition timing of the air-fuel mixture is retarded from the compression top dead center of the piston (assuming that it is retarded from the compression top dead center)
  • the fuel density in the cylinder at the ignition timing is Based on this, the ignition delay period is estimated. This is because, when the ignition timing of the air-fuel mixture is ahead of the compression top dead center of the piston, the cylinder volume decreases after the ignition of the air-fuel mixture, and the fuel density increases accordingly.
  • the ignition timing of the mixture is retarded from the compression top dead center of the piston, the in-cylinder volume increases after the ignition of the mixture, and the fuel density decreases accordingly. It is taken into consideration. In other words, if the ignition timing of the mixture is more advanced than the compression top dead center of the piston, there is a high correlation between the fuel density in the cylinder at the ignition timing and the ignition delay period, and the ignition of the mixture is When the timing is behind the compression top dead center of the piston, based on newly obtained knowledge that there is a high correlation between the fuel density in the cylinder at the ignition timing and the ignition delay period, The method for estimating the ignition delay period is varied depending on the ignition timing.
  • a correction coefficient based on the engine speed (for example, an exponential function of the engine speed) may be multiplied. That is, since the strength of the flow in the cylinder generally changes as the engine speed changes, the ignition delay period changes due to the influence of the disturbance. Therefore, by performing correction based on the engine rotation speed, the ignition delay period can be estimated more accurately.
  • a virtual ignition timing is set, and the estimated ignition delay period estimated according to the virtual ignition timing (e.g., calculated by an arithmetic expression) is calculated from the actual ignition timing to the virtual ignition timing.
  • the estimated ignition delay period estimated according to the virtual ignition timing e.g., calculated by an arithmetic expression
  • the ignition delay period is set based on the fuel density in the cylinder at the ignition timing.
  • the virtual ignition timing is on the retard side with respect to the compression top dead center of the piston, the ignition delay period is estimated based on the fuel density in the cylinder at the ignition timing.
  • the estimated ignition delay period is compared with a virtual ignition delay period that is a period between the actual ignition timing and the virtual ignition timing. Calculated as the true ignition delay period.
  • the heat generation rate waveform is calculated using this true ignition delay period.
  • the ignition timing for estimating the ignition delay period can be accurately obtained by repeated calculation, and the ignition delay period can be calculated with high accuracy.
  • the heat generation rate waveform calculated using the ignition delay period calculated as described above for example, the crank angle period from ignition of the air-fuel mixture to the end of combustion is used as the base, and the heat generation rate at the maximum heat generation rate is A triangular waveform with a vertex is mentioned. If the heat release rate waveform is approximated by this triangular waveform, the period from the ignition timing by the spark plug to the rising timing of the hypotenuse of the triangular waveform is defined as the ignition delay period.
  • the period from the ignition timing of the air-fuel mixture to the maximum heat generation rate period does not depend on at least one of the engine load factor, the air-fuel ratio, the EGR rate, and the oil / water temperature.
  • the predetermined heat generation rate is reached (more specifically, the cylinder volume at the crank angle position at the maximum heat generation rate; a parameter correlated with the turbulence in the cylinder) and the engine speed (also in the cylinder).
  • the triangular waveform is created as determined by a parameter having a correlation with the disturbance.
  • the first half combustion period does not change, and the change in the first half combustion period is an influence of the turbulence in the cylinder. Can be created. By doing so, the man-hours for creating the heat release rate waveform can be further reduced.
  • the present invention relates to a method for calculating a heat release rate waveform in a spark ignition type internal combustion engine. That is, first, a period from when ignition is performed by the spark plug to when the air-fuel mixture is ignited is defined as an ignition delay period that is one of the characteristic values of the heat generation rate waveform.
  • the ignition delay period is estimated based on the fuel density in the cylinder at the ignition timing, while the air-fuel mixture is
  • the ignition delay period is estimated based on the fuel density in the cylinder at the ignition timing, and the estimated ignition delay period is used.
  • a heat release rate waveform is calculated.
  • an ignition delay period which is a period from when ignition by the spark plug is performed until the air-fuel mixture ignites. Since estimation is based on the fuel density in the cylinder, the number of steps required to create the heat release rate waveform can be reduced, and the required accuracy for the ignition delay period can be ensured without creating the entire heat release rate waveform. However, it is possible to estimate and evaluate more easily than in the past.
  • FIG. 9A shows a case where the ignition timing SA is BTDC
  • FIG. 9B shows a case where the ignition timing SA is BTDC
  • FIG. 9B shows a case where the ignition timing SA is ATDC.
  • FIG. A heat generation rate waveform obtained in each engine operating state in which only the load rate is different from each other, and is a diagram in which each heat generation rate waveform in which the ignition timing SA is adjusted so that the maximum heat generation rate timing dQpeakA coincides with each other is displayed. It is.
  • a heat generation rate waveform obtained in each engine operating state in which only the EGR rate differs from each other, and is a diagram in which the respective heat generation rate waveforms in which the ignition timing SA is adjusted so that the maximum heat generation rate timing dQpeakA coincides with each other are displayed in an overlapping manner. It is.
  • FIG. It is a figure which shows the result of having verified the relationship between the prediction first combustion period calculated by Formula (3) with respect to a certain engine, and the measurement first half combustion period measured in the actual machine. It is a figure which shows the result of having verified the relationship between the prediction first half combustion period calculated by Formula (3) with respect to another engine, and the measurement first half combustion period measured in the actual machine.
  • a heat generation rate waveform obtained in each engine operating state in which only the load rate is different from each other is a diagram in which each heat generation rate waveform in which the ignition timing SA is adjusted so that the maximum heat generation rate timing dQpeakA coincides with each other is displayed. It is. It is the figure which displayed in piles the heat release rate waveform obtained in each engine operation state from which only ignition timing SA mutually differs. It is a figure which shows the result of the experiment which investigated the relationship between fuel density ⁇ fuel @ dQpeak at the time of heat release rate, and heat release rate inclination b / a with respect to each different engine speed Ne.
  • FIG. 1 is a diagram illustrating the configuration of the heat release rate waveform calculation device 1 according to the present embodiment and the input / output information of the heat release rate waveform calculation device 1.
  • the heat generation rate waveform calculation apparatus 1 receives various information on the engine state quantity, the control parameter control quantity, and the physical quantity. Examples of such input information include engine speed, load factor, ignition timing, EGR rate, air-fuel ratio, oil / water temperature, intake / exhaust valve open / close timing (valve timing), and the like. Further, the heat generation rate waveform calculation apparatus 1 estimates various characteristic values of the heat generation rate waveform by using the estimation units 2 to 5 storing the following estimation models based on each input information, and uses these various characteristic values. The heat release rate waveform created in this way is output.
  • the heat release rate waveform calculation apparatus 1 includes an ignition delay estimation unit that stores an ignition delay estimation model for estimating an ignition delay, a first half combustion period, a heat release rate gradient, and a heat generation amount as characteristic values of the heat release rate waveform. 2, the first half combustion period estimation unit 3 storing the first half combustion period estimation model, the heat generation rate inclination estimation unit 4 storing the heat generation rate inclination estimation model, and the heat generation amount estimation unit 5 storing the heat generation amount estimation model It has.
  • the ignition delay estimator 2 mixes the spark after the spark is ignited between the spark plug electrodes when the air-fuel mixture is ignited by the engine spark plug (ignition plug). This is a part for estimating the period (hereinafter referred to as the ignition delay period) until the time when the initial flame kernel is formed (hereinafter referred to as the ignition timing) by using the ignition delay estimation model.
  • This ignition delay period is represented by a crank angle [CA].
  • the ignition timing is defined as a time when the heat generation rate (heat generation amount per unit crank angle of crankshaft rotation) reaches 1 [J / CA] after the ignition timing. This value is not limited to this, and can be set as appropriate.
  • a timing when the heat generation amount after the ignition timing reaches a predetermined ratio (for example, 5%) with respect to the total heat generation amount may be set as the ignition timing.
  • the time when the ratio of the heat generation amount to the total heat generation amount reaches a predetermined value for example, the crank angle position when the heat generation amount reaches 10%
  • the time when the ratio of the heat generation amount reaches another predetermined value for example, the ignition timing may be defined based on the crank angle position at the time when 50% is reached.
  • a triangle triangular waveform that approximates the heat generation rate waveform during the period in which the heat generation rate increases is created based on the crank angle position and the rate of heat generation, and ignition is performed based on this triangular waveform.
  • the ignition timing may be defined based on this heat release rate waveform.
  • Each of the values is not limited to this, and can be set as appropriate.
  • the first half combustion period estimation unit 3 generates the heat during the period from the ignition timing to the maximum heat generation rate with the growth of the flame kernel (period from the ignition timing to the combustion end timing).
  • the first half combustion period which is the period until the rate becomes the maximum
  • the maximum heat generation rate timing and the first half combustion period are each represented by a crank angle [CA].
  • the heat generation rate inclination estimation unit 4 calculates an average rate of increase in heat generation rate (slope of heat generation rate) with respect to a change in crank angle in the first half combustion period, that is, a period from the ignition timing to the maximum heat generation rate. This is the part that is estimated using the heat release rate slope estimation model. That is, in this embodiment, as described below with reference to FIG. 2, a triangular waveform that approximates the heat generation rate waveform is created, and the heat generation rate inclination estimation unit 4 ignites the triangular waveform. The slope of the hypotenuse representing the heat generation rate from the time to the maximum heat generation rate is estimated.
  • the unit of inclination of the heat generation rate is represented by [J / CA 2 ].
  • the heat generation amount estimation unit 5 calculates the heat generation amount generated by the combustion of the air-fuel mixture (the heat generation amount generated during the entire combustion period, and the integration of the heat generation rate during the period from the ignition timing to the combustion end timing) Value) using a heat generation amount estimation model.
  • the unit of the heat generation amount is represented by [J].
  • the estimation operation in each of the estimation units 2 to 5 determines the characteristic values of the heat generation rate waveform such as the ignition delay, the first half combustion period, the heat generation rate inclination, and the heat generation amount, and uses these characteristic values.
  • a heat release rate waveform is created.
  • the generated heat generation rate waveform becomes the output of the heat generation rate waveform calculation device 1.
  • the ignition delay estimation unit 2 performs the ignition delay period estimation operation (step ST1), and the first half combustion period estimation unit. 3 in the first half combustion period (step ST2), a heat generation rate inclination estimation operation in the heat generation rate inclination estimation unit 4 (step ST3), and a heat generation amount estimation operation in the heat generation amount estimation unit 5 (step ST4).
  • step ST1 the ignition delay period estimation operation
  • step ST2 the first half combustion period estimation unit. 3 in the first half combustion period
  • step ST3 a heat generation rate inclination estimation operation in the heat generation rate inclination estimation unit 4
  • a heat generation amount estimation operation in the heat generation amount estimation unit 5 step ST4
  • FIG. 2 shows an example of a heat release rate waveform that is created using the characteristic values estimated by the respective estimation units 2 to 5 and is output from the heat release rate waveform calculation device 1.
  • the timing SA in the figure is the ignition timing
  • the timing FA in the figure is the ignition timing.
  • ⁇ in the figure is the ignition delay period.
  • dQpeakA is the maximum heat generation rate time
  • the heat generation rate at this heat generation rate maximum time dQpeakA is b in the figure. That is, this heat generation rate b is the maximum heat generation rate during the combustion period.
  • a in the figure, which is a period from the ignition timing FA to the maximum heat generation rate timing dQpeakA is the first half combustion period.
  • the inclination of the heat release rate in the first half combustion period a is expressed as b / a.
  • c in the figure which is a period from the maximum heat generation rate dQpeakA to the combustion end timing EA, is the second half combustion period.
  • Q1 in the figure is the heat generation amount in the first half combustion period a
  • Q2 is the heat generation amount in the second half combustion period c.
  • the heat generation amount (total heat generation amount Q all ) generated during the entire combustion period is expressed as the sum of the heat generation amount Q1 and the heat generation amount Q2.
  • the heat generation rate waveform calculation apparatus 1 of the present embodiment uses the crank angle period (FA to EA in the figure) from the ignition of the air-fuel mixture to the end of combustion as the base, and the heat generation rate at the maximum heat generation rate dQpeakA.
  • the heat generation rate waveform is approximated by a triangular waveform having b as a vertex.
  • the heat generation rate waveform that is the output of the heat generation rate waveform calculating device 1 is used to examine the system at the time of engine design, the control, and the adaptive value.
  • the ignition delay estimation unit 2 is a part that estimates an ignition delay period ⁇ that is a period from the ignition timing SA to the ignition timing FA.
  • the estimation process of the ignition delay period ⁇ performed in the ignition delay estimation unit 2 is as follows.
  • the ignition delay period ⁇ is estimated using either one of the following formulas (1) and (2) (these formulas correspond to an ignition delay estimation model).
  • ⁇ fuel @ SA is the fuel density in the cylinder at the ignition timing SA (in-cylinder fuel amount [mol] / in-cylinder volume [L] at the ignition timing).
  • ⁇ fuel @ FA is the fuel density in the cylinder at the ignition timing FA (in-cylinder fuel amount [mol] / cylinder volume at the ignition timing [L]).
  • Ne is the engine speed.
  • C 1 , C 2 , ⁇ , ⁇ , ⁇ , ⁇ are coefficients identified based on experiments or the like.
  • the air-fuel ratio is the stoichiometric air-fuel ratio
  • the EGR rate is “0”
  • the engine warm-up operation is complete (the oil-water temperature is equal to or higher than a predetermined value).
  • the equation is established on the condition that the opening / closing timing of the intake valve is fixed.
  • Formula (1) is a formula for calculating the ignition delay period ⁇ when the air-fuel mixture is ignited on the advance side (BTDC) with respect to the timing at which the piston reaches compression top dead center (BTDC) (hereinafter referred to as BTDC ignition).
  • Equation (2) is a calculation of the ignition delay period ⁇ when the air-fuel mixture is ignited on the retard side (ATDC) with respect to the timing at which the piston has reached compression top dead center (ATDC) (hereinafter referred to as ATDC ignition). It is a formula.
  • the ignition delay period ⁇ is calculated by an arithmetic expression using the in-cylinder fuel density ⁇ fuel and the engine rotational speed Ne as variables at a predetermined timing.
  • FIG. 4 is a graph showing the results of experiments measuring changes in the ignition delay period ⁇ with respect to changes in the in-cylinder fuel density ⁇ fuel @ SA at the ignition timing SA in the case of BTDC ignition.
  • the air-fuel ratio is set to the stoichiometric air-fuel ratio
  • the EGR rate is set to “0”
  • the engine warm-up operation is completed (the oil / water temperature is equal to or higher than a predetermined value)
  • the opening / closing timing of the intake valve is fixed. It has been done.
  • the engine rotation speed Ne increases in the order of “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “+”, “ ⁇ ”.
  • ⁇ is 800 rpm
  • ⁇ is 1000 rpm
  • ⁇ is 1200 rpm
  • ⁇ is 1600 rpm
  • ⁇ is 2400 rpm
  • + is 3200 rpm
  • ⁇ is 3600 rpm.
  • the ignition delay period ⁇ is shorter as the fuel density ⁇ fuel @ SA in the cylinder at the ignition timing SA is higher. This is probably because the higher the fuel density ⁇ fuel @ SA , the more fuel molecules around the spark plug, and the rapid growth of the flame kernel after ignition of the spark plug.
  • the engine rotation speed Ne affects the ignition delay period ⁇ . That is, the ignition delay period ⁇ is shorter as the engine speed Ne is higher. This is presumably because the higher the engine speed Ne, the stronger the disturbance of the air-fuel mixture flow in the cylinder (hereinafter simply referred to as turbulence) and the rapid growth of the flame kernel.
  • the in-cylinder fuel density ⁇ fuel @ SA and the engine speed Ne at the ignition timing SA are parameters that affect the ignition delay period ⁇ .
  • FIG. 5 is a graph showing a result of verifying the relationship between the predicted ignition delay period calculated by the equation (1) and the actually measured ignition delay period measured in the actual machine.
  • a prediction formula obtained by identifying the coefficients of C 1 , ⁇ , and ⁇ in equation (1) according to the engine operating conditions is used.
  • the engine rotational speed Ne increases in the order of “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “+”, “ ⁇ ”, “ ⁇ ”.
  • ⁇ is 800 rpm
  • ⁇ is 1000 rpm
  • ⁇ is 1200 rpm
  • ⁇ is 1600 rpm
  • ⁇ is 2000 rpm
  • + is 2400 rpm
  • ⁇ is 3200 rpm
  • ⁇ is 3600 rpm.
  • the predicted ignition delay period substantially coincides with the actually measured ignition delay period, and it can be seen from Equation (1) that the ignition delay period when BTDC ignition is performed is calculated with high accuracy. .
  • FIG. 6 is a graph showing the results of experiments measuring changes in the ignition delay period ⁇ with respect to changes in the fuel density ⁇ fuel @ FA in the cylinder at the ignition timing FA in the case of ATDC ignition.
  • the engine speed is fixed
  • the air-fuel ratio is the stoichiometric air-fuel ratio
  • the EGR rate is set to “0”
  • the engine warm-up operation is complete (the oil water temperature is equal to or higher than a predetermined value)
  • the intake valve This is performed with a fixed opening / closing timing.
  • the load factor increases in the order of “ ⁇ ”, “ ⁇ ”, “+”, and “ ⁇ ”. For example, “ ⁇ ” indicates a load factor of 20%
  • “ ⁇ ” indicates a load factor of 30%
  • “+” indicates a load factor of 40%
  • “ ⁇ ” indicates a load factor of 50%.
  • the correlation between the fuel density ⁇ fuel @ FA in the cylinder at the ignition timing FA and the ignition delay period ⁇ does not depend on the load factor (regardless of the load factor). is there. That is, these correlations can be represented by a single curve.
  • the ignition delay period ⁇ is shorter as the fuel density ⁇ fuel @ FA in the cylinder at the ignition timing FA is higher. This is because, as described above, the higher the fuel density ⁇ fuel @ FA , the greater the number of fuel molecules around the spark plug, and the rapid growth of the flame kernel after ignition of the spark plug. It is done.
  • the fuel density ⁇ fuel @ FA in the cylinder at the ignition timing FA is a parameter that affects the ignition delay period ⁇ .
  • the engine speed Ne is also assumed to be a parameter that affects the ignition delay period ⁇ .
  • FIG. 7 is a graph showing a result of verifying the relationship between the predicted ignition delay period calculated by the equation (2) and the actually measured ignition delay period measured in the actual machine.
  • a prediction formula obtained by identifying each coefficient of C 2 , ⁇ , and ⁇ in the formula (2) according to the engine operating condition is used.
  • the engine rotational speed Ne increases in the order of “ ⁇ ”, “ ⁇ ”, “+”, and “ ⁇ ”. For example, “ ⁇ ” is 800 rpm, “ ⁇ ” is 1200 rpm, “+” is 3600 rpm, and “ ⁇ ” is 4800 rpm.
  • the predicted ignition delay period substantially coincides with the actually measured ignition delay period, and it can be seen from Equation (2) that the ignition delay period when ATDC ignition is performed is calculated with high accuracy. .
  • the inventor of the present invention derived the equations (1) and (2) based on these new findings.
  • the ignition timing SA is also advanced from the timing at which the piston reaches compression top dead center (BTDC). It is. In this case, after reaching the ignition timing SA, the piston moves toward the compression top dead center. That is, the in-cylinder volume decreases, and the fuel density ⁇ fuel increases accordingly.
  • the fuel density [rho Fuel towards the fuel density [rho Fuel @ SA at the ignition timing SA is smaller than the fuel density [rho Fuel @ FA in ignition timing FA.
  • the ignition delay period ⁇ can be obtained with high accuracy.
  • FIG. 9 a diagram showing the ignition timing SA and the heat generation rate waveform
  • the ignition timing SA is advanced (BTDC) from the timing at which the piston reaches compression top dead center. (See FIG. 9 (a)) and on the retard side (ATDC) (see FIG. 9 (b)).
  • BTDC spark-to-distance
  • ATDC retard side
  • the piston moves toward the bottom dead center. That is, the in-cylinder volume increases, and the fuel density ⁇ fuel decreases accordingly. Therefore, as the fuel density ⁇ fuel , there is a high possibility that the fuel density ⁇ fuel @ FA at the ignition timing FA is smaller than the fuel density ⁇ fuel @ SA at the ignition timing SA.
  • the fuel density ⁇ fuel @ FA at the ignition timing FA which is a value correlated with the maximum value of the ignition delay period (the longest ignition delay period among the assumed ignition delay periods), is determined in advance. By multiplying the coefficient, the ignition delay period ⁇ can be obtained with high accuracy.
  • a procedure for determining which of these equations (1) and (2) is used (a procedure for determining whether the ignition timing is BTDC ignition or ATDC ignition) and an ignition delay period (true).
  • the procedure for calculating the ignition delay period) is as follows. First, a virtual ignition timing is set, and the cylinder volume at the virtual ignition timing is obtained. Since this in-cylinder volume can be obtained geometrically by the crank angle position (piston position) corresponding to the virtual ignition timing, the in-cylinder volume is uniquely determined from the virtual ignition timing. Then, the fuel density is obtained from the in-cylinder volume and the fuel injection amount.
  • the estimated ignition delay period is calculated by substituting the fuel density and the engine speed at the virtual ignition timing into Equation (1).
  • the estimated ignition delay period is calculated by substituting the fuel density and the engine speed at the virtual ignition timing into Equation (2). Then, a timing that is advanced by the calculated estimated ignition delay period with respect to the virtual ignition timing is set as the virtual ignition timing.
  • the virtual ignition timing is compared with the actual ignition timing (ignition timing as input information). If the virtual ignition timing does not match the actual ignition timing, the virtual ignition timing is changed. For example, the virtual ignition timing is changed to the retard side. Then, again, the fuel density and the engine speed at the virtual ignition timing are substituted into formula (1) or formula (2) (if the virtual ignition timing is set as BTDC ignition, it is substituted into formula (1).
  • the estimated ignition delay period is calculated by substituting into the equation (2) to obtain the virtual ignition timing, and this and the actual ignition timing (input information) Ignition timing). This operation is repeated, and the virtual ignition timing when the virtual ignition timing and the actual ignition timing coincide with each other is obtained as the true ignition timing.
  • the estimated ignition delay period calculated in Expression (1) or Expression (2) is obtained as the true ignition delay period.
  • the true ignition timing is BTDC (in the case of BTDC ignition)
  • the calculated ignition timing is substituted into the equation (1) to calculate the ignition delay period ⁇ , and the true ignition timing is calculated.
  • the ignition delay period ⁇ may be calculated by substituting the calculated ignition timing into Equation (2).
  • the above operation is as follows. A period between the actual ignition timing and the virtual ignition timing (virtual ignition delay period when the ignition is performed at the virtual ignition timing) and the formula (1) or the formula (2) (estimated) The estimated ignition delay period is compared, and if they do not match, the virtual ignition timing is changed. Then, after again calculating the estimated ignition delay period by the formula (1) or the formula (2), the period between the actual ignition timing and the virtual ignition timing (virtual ignition delay period) and the formula (1) ) Or the estimated ignition delay period calculated by the equation (2). This operation is repeated, and the estimated ignition delay period when these match (when the virtual ignition delay period and the estimated ignition delay period match) is obtained as the true ignition delay period.
  • the ignition delay period ⁇ can be estimated for the entire engine operation area.
  • the ignition timing FA can be obtained by adding the ignition delay period ⁇ to the ignition timing SA.
  • the first half combustion period estimation unit 3 is a part that estimates the first half combustion period a that is a period from the ignition timing FA to the maximum heat release rate dQpeakA.
  • the estimation process of the first half combustion period a performed in the first half combustion period estimation unit 3 is as follows.
  • This first half combustion period a [CA] is estimated using the following equation (3) (this equation corresponds to the first half combustion period estimation model).
  • V @ dQpeak is an in-cylinder volume [L] as a physical quantity at the maximum heat generation rate dQpeakA, and is hereinafter also referred to as a maximum heat generation rate in-cylinder volume.
  • Ne is an engine speed (engine speed).
  • This expression (3) is an expression that is established on the condition that the opening / closing timing of the intake valve is fixed. Further, this equation (3) is established without being affected by the load factor, EGR rate, air-fuel ratio, and oil water temperature. That is, Equation (3) is established based on the fact that the first half combustion period a is not affected by the load factor, EGR rate, air-fuel ratio, and oil water temperature.
  • FIGS. 10 to 13 is a heat generation rate waveform obtained in different engine operating states, and the heat generation rate waveforms obtained by adjusting the ignition timing SA so that the maximum heat generation rate timing dQpeakA coincides with each other. It is displayed.
  • FIG. 10 shows the heat generation rate waveforms obtained in each engine operating state in which only the load factor is different from each other.
  • FIG. 11 shows the heat generation rate waveforms obtained in each engine operating state in which only the EGR rate is different from each other.
  • FIG. 12 shows the heat generation rate waveform obtained in each engine operating state in which only the air-fuel ratio is different from each other.
  • FIG. 13 shows the heat generation rate waveforms obtained when only the oil and water temperatures are different from each other as in the course of engine warm-up operation.
  • the first half combustion period a is maintained constant regardless of any of the load factor, EGR rate, air-fuel ratio, and oil / water temperature. That is, it can be seen that the first half combustion period a is not affected by the load factor, EGR rate, air-fuel ratio, and oil water temperature.
  • FIG. 14 shows the heat generation rate waveforms obtained in each engine operating state in which only the ignition timing SA is different from each other. As can be seen from FIG. 14, the first half combustion period a becomes longer as the ignition timing SA is retarded.
  • FIG. 15 is a heat generation rate waveform obtained in each engine operating state in which only the engine rotation speed Ne is different from each other, and each heat generation in which the ignition timing SA is adjusted so that the maximum heat generation rate dQpeakA coincides with each other.
  • the rate waveform is superimposed and displayed. Since the crank rotation angle [CA] per unit time [ms] increases as the engine speed Ne increases, the first half combustion period a should be longer (longer on the crank angle axis). In the case shown in FIG. 15, the first half combustion period “a” hardly changes even when the engine speed Ne is different. This is considered to be due to the fact that the higher the engine rotation speed Ne, the shorter the first half combustion period a. In other words, the higher the engine speed Ne, the longer the first combustion period a due to the larger crank rotation angle per unit time, and the shorter the first combustion period a due to “other factors”. Can be assumed.
  • the first half combustion period a is affected by the ignition timing SA and the engine speed Ne.
  • the ignition timing SA and the engine rotational speed Ne affect the turbulence in the cylinder.
  • the ignition timing FA and the maximum heat generation rate dQpeakA also shift to the retard side as the ignition timing SA shifts to the retard side.
  • the in-cylinder volume at the maximum heat generation rate dQpeakA increases and the turbulence in the cylinder decreases .
  • the disturbance in a cylinder becomes weak, flame propagation will become slow and the first half combustion period a will become long.
  • the lower the engine speed Ne the lower the flow velocity of the air flowing into the cylinder from the intake system, and the turbulence in the cylinder becomes weaker. And if the disturbance in a cylinder becomes weak, flame propagation will become slow and the first half combustion period a will become long.
  • the higher the engine speed Ne the higher the flow velocity of the air flowing into the cylinder from the intake system, and the turbulence in the cylinder becomes stronger. And if the disturbance in a cylinder becomes strong, flame propagation will become rapid and the first half combustion period a will become short.
  • the above-mentioned “other factors (factors that shorten the first half combustion period a)” is that the flame propagation is rapid due to the fact that the higher the engine speed Ne, the stronger the turbulence in the cylinder. .
  • the inventor of the present invention derived the formula (3) based on this new knowledge.
  • the in-cylinder volume, particularly the maximum in-cylinder volume V @dQpeak is used as a variable as a physical quantity correlated with the ignition timing SA that is the control amount. That is, as described above, as the ignition timing SA shifts to the retard side, the maximum heat generation rate timing dQpeakA also shifts to the retard side, and the in-cylinder volume V @dQpeak increases, which is correlated with the ignition timing SA.
  • the maximum heat release rate in-cylinder volume V @dQpeak is used as a variable.
  • the procedure for obtaining the maximum heat release rate in-cylinder volume V @dQpeak and the procedure for calculating the first half combustion period a are the following variables.
  • a virtual heat generation rate maximum time is set, and the in-cylinder volume at the virtual heat generation rate maximum time is obtained. Since the cylinder volume can be geometrically determined by the crank angle position (piston position) corresponding to the virtual maximum heat generation rate, the cylinder volume is unambiguous from the virtual maximum heat generation rate. It is decided.
  • the estimated first half combustion period is calculated by substituting the in-cylinder volume and the engine rotation speed at the virtual maximum heat generation rate timing into Expression (3). Then, a timing that is advanced by the calculated first half combustion period with respect to the virtual maximum heat generation rate is set as a virtual ignition timing.
  • the ignition timing FA can be calculated by adding the ignition delay period ⁇ to the ignition timing SA.
  • the virtual ignition timing is compared with the calculated ignition timing FA. If the virtual ignition timing does not match the calculated ignition timing FA, the virtual maximum heat generation rate timing is changed. For example, the virtual maximum heat generation rate is changed to the retard side. Then, again, the in-cylinder volume and the engine rotation speed at the virtual maximum heat generation rate time are substituted into the equation (3) to calculate the estimated first half combustion period, and the virtual ignition timing is obtained and calculated.
  • the ignition timing FA ignition timing FA obtained by adding the ignition delay period ⁇ calculated by the ignition delay estimation unit 2 to the ignition timing SA is compared.
  • the above operation is as follows.
  • the period between the ignition timing FA (ignition timing determined according to the actual ignition timing) and the virtual maximum heat generation rate (virtual first half combustion period) was calculated by the formula (3) (estimated) )
  • the estimated first half combustion period (first half combustion period estimated based on the physical quantity at the virtual maximum heat generation rate timing) is compared, and if these do not match, the virtual maximum heat generation rate timing is changed.
  • the period between the ignition timing FA and the virtual maximum heat generation rate (virtual first half combustion period) is calculated using equation (3). Compare the estimated first half combustion period. This operation is repeated, and the estimated first half combustion period when these match (when the virtual first half combustion period and the estimated first half combustion period match) is obtained as the true first half combustion period a.
  • Equation (3) C and ⁇ are obtained by identification based on experiments and the like. Further, ⁇ is a value corresponding to the tumble ratio in the cylinder, and is given as a larger value as the tumble ratio is larger. Note that ⁇ may be set by identification based on experiments or the like. These coefficients can also be identified with respect to changes in the opening / closing timing of the intake valve. As described above, the expression (3) is multiplied by the exponential function (correction coefficient) of the engine rotational speed Ne based on the in-cylinder volume V @dQpeak at the maximum heat generation rate and with the value ⁇ corresponding to the tumble ratio as an exponent. Thus, the first half combustion period a is calculated.
  • FIG. 16 and FIG. 17 are graphs showing the results of verifying the relationship between the predicted first half combustion period calculated by Equation (3) and the first half actual measured combustion period measured in the actual machine for different engines.
  • a prediction formula obtained by identifying the coefficient C in the formula (3) according to the engine operating condition is used.
  • the engine speed Ne increases in the order of “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “+”, and “ ⁇ ”.
  • “ ⁇ ” is 800 rpm
  • “ ⁇ ” is 1000 rpm
  • ⁇ ” is 1200 rpm
  • “ ⁇ ” is 1600 rpm
  • “ ⁇ ” is 2400 rpm
  • “+” is 3200 rpm
  • “ ⁇ ” is 3600 rpm.
  • the engine rotation speed Ne increases in the order of “ ⁇ ”, “ ⁇ ”, “+”, “ ⁇ ”, and “ ⁇ ”. For example, “ ⁇ ” is 800 rpm, “ ⁇ ” is 1200 rpm, “+” is 2400 rpm, “ ⁇ ” is 3600 rpm, and “ ⁇ ” is 4800 rpm.
  • the predicted first half combustion period substantially coincides with the measured first half combustion period, and it can be seen that the first half combustion period a is calculated with high accuracy according to the equation (3).
  • the first half combustion period a is not affected by the load factor, air-fuel ratio, EGR rate, and oil / water temperature, and is based on the maximum heat release rate in-cylinder volume V @dQpeak and the engine speed Ne. Can be estimated.
  • the in-cylinder volume V @dQpeak at the time of maximum heat generation rate and the engine rotation speed Ne are parameters that correlate with the disturbance in the cylinder as described above.
  • the load factor, the air-fuel ratio, the EGR rate, and the oil / water temperature are assumed to have no influence on the first half combustion period a because there is almost no correlation with the disturbance in the cylinder.
  • the maximum heat release rate cylinder volume V @dQpeak and the engine speed Ne which are parameters correlated with the turbulence in the cylinder. Since the first half combustion period a can be estimated based on this, the man-hours for determining the first half combustion period a under various operating conditions of the engine can be greatly reduced.
  • the first half combustion period is not affected by the load factor.
  • This load factor is one of the parameters for controlling the fuel injection amount, and the fuel injection amount is a control parameter that affects the fuel density in the cylinder.
  • the first half combustion period is estimated without depending on the fuel density in the cylinder.
  • the first half combustion period is estimated based on parameters that affect the turbulence in the cylinder, such as the maximum heat generation rate in-cylinder volume V @dQpeak and the engine rotation speed Ne.
  • the heat generation rate gradient is estimated based on the fuel density in the cylinder.
  • the first half combustion period and the heat generation rate gradient estimated in the present embodiment are estimated as values independent of each other (not dependent on each other).
  • the heat generation rate gradient estimation unit 4 is a portion that estimates the heat generation rate gradient b / a (hereinafter referred to as the heat generation rate gradient) in the first half combustion period a.
  • the heat generation rate gradient b / a estimation process performed in the heat generation rate gradient estimation unit 4 is as follows.
  • This heat release rate gradient b / a [J / CA 2 ] is basically estimated using the following equation (4) (this equation corresponds to a heat release rate gradient estimation model).
  • ⁇ fuel @ dQpeak is the fuel density (in-cylinder fuel amount [mol] / in-cylinder volume [L] at the maximum heat generation rate) at the maximum heat generation rate dQpeakA; Also called density.
  • C 3 is a coefficient identified based on experiments or the like.
  • FIGS. 18 (a) to 18 (d) are heat generation rate waveforms obtained in the engine operating states in which only the load rate is different from each other, and the ignition timing SA is set so that the maximum heat generation rate timings dQpeakA coincide with each other.
  • the heat release rate waveforms adjusted for are superimposed and displayed.
  • the ignition timing changes to the retard side in the order of FIGS. 18 (a) to 18 (d), and in each figure, the load factor increases in the order of KL1, KL2, and KL3.
  • KL1 has a load factor of 20%
  • KL2 has a load factor of 30%
  • KL3 has a load factor of 40%.
  • the heat generation rate gradient b / a is affected by the load factor and the ignition timing SA. Specifically, in any of FIGS. 18A to 18D having different ignition timings SA, the heat generation rate gradient b / a increases as the load factor increases. As a factor that the heat generation rate gradient b / a is affected by the load factor in this way, it is conceivable that the fuel density in the cylinder changes according to the load factor. That is, it is considered that the higher the load factor, the greater the in-cylinder fuel amount, and the higher the fuel density in the cylinder and the higher the combustion speed of the air-fuel mixture.
  • 19 (a) and 19 (b) are graphs in which the heat release rate waveforms obtained in each engine operating state in which only the ignition timing SA is different from each other are displayed in an overlapping manner in order to investigate the influence of the change in the ignition timing SA. It is. In FIG. 19A and FIG. 19B, the load factors are different from each other. In either case, the heat generation rate gradient b / a decreases as the ignition timing SA shifts to the retard side. Tend to be.
  • the factor that the heat generation rate gradient b / a is influenced by the ignition timing SA is also considered to be due to the fuel density in the cylinder as in the case of the load factor described above. That is, when the piston is in the vicinity of compression top dead center (TDC), the change in the in-cylinder volume accompanying the change in the crank angle is small, but as it moves away from the TDC in the expansion stroke (for example, from about ATDC 10 ° CA) H) The cylinder volume increases as the retard angle increases, and the fuel density in the cylinder decreases accordingly.
  • TDC compression top dead center
  • the heat generation rate waveform moves to the retard side as a whole with the retard of the ignition timing SA, and the ignition timing FA (waveform of the waveform).
  • the slope of the heat release rate waveform gradually decreases.
  • the slope of the straight line (indicated by the alternate long and short dash line in the figure) connecting from the ignition timing FA (starting point of the waveform) to the heat generation rate b (the peak of the waveform) at the maximum heat generation rate dQpeakA, that is, the heat generation rate.
  • the inclination b / a is also gradually reduced toward the retard side.
  • the influence of the retard angle of the ignition timing SA (that is, the retard angle of the ignition timing FA) on the heat generation rate gradient b / a is that the heat generation rate gradient b / a and the maximum heat generation rate fuel density ⁇ fuel @ This is evident in the relationship with dQpeak . That is, as shown in FIGS. 19 (a) and 19 (b), the maximum heat generation rate dQpeakA shifts to the retard side as the ignition timing SA is retarded, and at this heat generation rate maximum time dQpeakA. As the in-cylinder volume (in-cylinder volume V @dQpeak at the maximum heat generation rate) gradually increases, the fuel density ⁇ fuel @ dQpeak at the maximum heat generation rate decreases accordingly. As a result, the heat generation rate gradient b / a becomes smaller.
  • the present inventor examined how the heat generation rate gradient b / a changes corresponding to the change in the fuel density ⁇ fuel @ dQpeak at the maximum heat generation rate.
  • the results of this experiment are shown in the graphs of FIGS. 20 (a) to 20 (d).
  • the load factor increases in the order of “ ⁇ ”, “ ⁇ ”, “+”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”, “ ⁇ ”. For example, in FIG.
  • ⁇ indicates a load factor of 15%
  • ⁇ indicates a load factor of 20%
  • + indicates a load factor of 25%
  • ⁇ indicates a load factor of 30%
  • ⁇ indicates a load factor of 35%
  • ⁇ indicates a load factor of 40%
  • ⁇ indicates a load factor of 45%
  • ⁇ indicates a load factor of 50%.
  • FIGS. 20 (a) to 20 (d) increases in the order of FIGS. 20 (a) to 20 (d).
  • FIG. 20 (a) is 800 rpm
  • FIG. 20 (b) is 1200 rpm
  • FIG. 20 (c) is 2000 rpm
  • FIG. 20 (d) shows 3200 rpm.
  • the inventor of the present invention derived the formula (4) based on this new knowledge.
  • the maximum heat generation rate fuel density ⁇ fuel @ dQpeak which is a variable of the equation (4), can be obtained by dividing the in-cylinder fuel amount by the maximum heat generation rate in-cylinder volume V @dQpeak as described above. it can.
  • the procedure for obtaining the maximum heat release rate in-cylinder volume V @dQpeak is as described above in the description of the first half combustion period estimation unit 3.
  • the in-cylinder fuel amount is given as input information of the heat release rate waveform calculation device 1.
  • the heat generation rate slope b / a which is one of the characteristic values of the heat generation rate waveform, is basically a linear function of the maximum heat generation rate fuel density ⁇ fuel @ dQpeak (in this embodiment, As an example, it can be calculated as a proportional function.
  • the heat generation rate gradient b / a can be estimated mainly based on the maximum heat generation rate fuel density ⁇ fuel @ dQpeak without considering the load factor and the ignition timing SA, the heat in various operating conditions of the engine can be estimated. The number of steps required to determine the incidence slope b / a is reduced.
  • the heat generation amount estimation unit 5 is a portion that estimates the heat generation amount (total heat generation amount Q all ) generated during the entire combustion period.
  • the heat generation amount Q1 in the first half combustion period a is calculated by the following equation (5).
  • the total heat generation amount Q all can be expressed as in-cylinder fuel amount ⁇ k (combustion efficiency) (this equation corresponds to a heat generation amount estimation model).
  • the combustion efficiency k decreases when the oil / water temperature is low, for example, during warm-up operation, and also changes due to changes in the load factor, engine speed, and the like. Therefore, in this embodiment, a map for determining the value of the combustion efficiency k in association with the oil / water temperature, the load factor, and the engine rotation speed is determined in advance using a database of experimental results of various engines. Then, the total heat generation amount Q all is calculated from the in-cylinder fuel amount using the value of the combustion efficiency k.
  • the heat generation rate b at the maximum heat generation rate dQpeakA is obtained by the following equation (7), and the second half combustion period c is obtained by the following equation (8).
  • the ignition delay estimation unit 2 using the ignition delay estimation model 2 estimates the ignition delay period ⁇
  • the first combustion period estimation unit 3 uses the first combustion period estimation model.
  • Estimation of the first half combustion period a, estimation of the heat generation rate inclination b / a in the heat generation rate inclination estimation unit 4 using the heat generation rate inclination estimation model, and heat generation amount estimation unit 5 using the heat generation amount estimation model The heat generation amount Q all at is estimated, and the maximum heat generation rate b and the second half combustion period c are calculated.
  • a triangular waveform that approximates the heat release rate waveform is created as shown in FIG. 2, and this triangular waveform is output as the heat release rate waveform.
  • the system design, the control, and the adaptive value are examined for engine design.
  • the ignition delay period ⁇ which is one of the characteristic values of the waveform, is calculated from the fuel density ⁇ fuel at a predetermined time. I am trying to calculate. That is, when the ignition timing of the air-fuel mixture is on the advance side (BTDC ignition) from the compression top dead center of the piston, the ignition delay period ⁇ is based on the fuel density ⁇ fuel @ SA in the cylinder at the ignition timing SA. On the other hand, if the ignition timing of the air-fuel mixture is retarded from the compression top dead center (ATDC ignition), the ignition is performed based on the fuel density ⁇ fuel @ FA in the cylinder at the ignition timing FA. The delay period ⁇ is calculated. For this reason, compared with calculating based on both the engine load factor and the ignition timing, man-hours can be reduced.
  • the first half combustion period a which is also one of the characteristic values, is assumed not to depend on any of the engine load factor, air-fuel ratio, EGR rate, and oil water temperature, and the maximum heat release rate in the cylinder.
  • the calculation is made based on the volume V @ dQpeak and the engine speed Ne. Thereby, the man-hour concerning calculation of the first half combustion period a can be reduced significantly.
  • the heat generation rate waveform is created based on the ignition delay period ⁇ , this heat generation rate The waveform is created according to the physical phenomenon in the combustion state in the cylinder.
  • the heat generation rate according to the present embodiment is compared with a method of creating a heat generation rate waveform using a Webe function, in which various parameters such as a shape parameter are mathematically combined to simply match the waveform shape.
  • the heat generation rate waveform created by the waveform calculation device 1 can obtain high reliability.
  • the ignition delay period ⁇ can be calculated from the fuel density ⁇ fuel in the cylinder and the engine rotational speed Ne as described above without creating the entire heat release rate waveform.
  • the ignition delay period ⁇ can be estimated and evaluated more easily than in the past while ensuring the required accuracy.
  • the first half combustion period a and the heat generation rate inclination b / a estimated in the present embodiment are estimated as values independent of each other (not dependent on each other). Therefore, the first half combustion period a is estimated as a value that is not affected by the error when the heat generation rate inclination b / a includes an estimation error, and the heat generation rate inclination b / a is When an estimation error is included in a, it is estimated as a value that is not affected by the error. As a result, high estimation accuracy of these values can be ensured.
  • the heat generation rate waveform calculation method implemented in the heat generation rate waveform calculation apparatus described in the above embodiment is also within the scope of the technical idea of the present invention.
  • the average rate of increase of the heat generation rate in the period from the ignition timing FA of the air-fuel mixture to the maximum heat generation rate dQpeakA is the heat generation rate slope b / a
  • the equation (4) The heat generation rate gradient b / a is calculated as a linear function of the maximum heat generation rate fuel density ⁇ fuel @ dQpeak , but is not limited thereto.
  • the heat generation rate is increasing from the ignition timing FA to the maximum heat generation rate dQpeakA (heat generation rate increase period), for example, a predetermined time slightly before the heat generation rate maximum time dQpeakA from the ignition timing.
  • the rate of increase in the heat generation rate during the period up to the time may be a heat generation rate gradient, and the heat generation rate gradient may be estimated based on the fuel density at the predetermined time.
  • a virtual ignition timing is set and is calculated by repetitive calculation of the formula (1) or the formula (2).
  • the present invention is not limited to this, the ignition timing is sensed in a test with an actual machine, and the ignition timing is set based on this, or a desired ignition timing is input as an input signal to the heat release rate waveform calculation device 1 By doing so, the ignition delay period ⁇ may be obtained.
  • a virtual maximum heat generation rate is set, and the calculation is performed by iterative calculation of Equation (3). I was trying to do it. Not limited to this, the maximum heat generation rate is sensed in a test with an actual machine, and the maximum heat generation rate is set based on this, or the desired heat is input as an input signal to the heat generation rate waveform calculation device 1.
  • the maximum heat generation rate in-cylinder volume V @dQpeak and the first half combustion period a may be obtained by inputting the maximum generation rate timing.
  • the heat release rate waveform calculation apparatus 1 outputs a triangular waveform.
  • the present invention is not limited to this, and a heat generation rate waveform may be generated by performing predetermined filter processing on the generated triangular waveform, and the heat generation rate waveform may be output.
  • the first half combustion period a is calculated as not depending on any of the engine load factor, EGR rate, air-fuel ratio, and oil / water temperature, but at least one of these operating conditions is calculated. You may make it calculate as what does not depend.
  • the present invention it is possible to reduce the man-hours required to create the heat release rate waveform in the spark ignition type internal combustion engine, and it is possible to reduce the cost. Therefore, the invention can be applied to an internal combustion engine for automobiles.

Landscapes

  • Engineering & Computer Science (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Mechanical Engineering (AREA)
  • General Engineering & Computer Science (AREA)
  • Combined Controls Of Internal Combustion Engines (AREA)
  • Electrical Control Of Ignition Timing (AREA)

Abstract

【課題】内燃機関の熱発生率波形を作成するための工数を低減すると共に、例えば着火遅れ期間について、所要の精度を確保しながら従来よりも簡易に推定し、評価できるようにする。 【解決手段】点火栓による点火が行われてから混合気が着火するまでの期間を、熱発生率波形の特性値の一つである着火遅れ期間τとして規定する。混合気の着火時期ＦＡがピストンの圧縮上死点よりも進角側ＢＴＤＣである場合には、点火時期ＳＡにおける気筒内の燃料密度ρfuel@SAに基づいて着火遅れ期間τを推定する一方、混合気の着火時期ＦＡがピストンの圧縮上死点よりも遅角側ＡＴＤＣである場合には、着火時期ＦＡにおける気筒内の燃料密度ρfuel@FAに基づいて前記着火遅れ期間τを推定する。この推定した着火遅れ期間τを用いて熱発生率波形を作成する。

Description

内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
　本発明は、火花点火式の内燃機関における熱発生率波形を算出するための装置、およびその算出方法に関連し、特に、点火栓による点火が行われてから混合気が着火するまでの期間（本明細書では、この期間を「着火遅れ期間」と呼ぶ）に着目して、熱発生率波形を得るための技術に係る。
　従来より、内燃機関の燃焼状態を表現するために、気筒内の熱発生率をウィーベ（Wiebe）関数によって近似することが行われている。ウィーベ関数は、複数のパラメータを特定することによって、熱発生率波形を好適に表現することができるものであって、内燃機関の燃焼による熱発生率や質量燃焼割合などを推定するために利用されている。
　例えば、特許文献１に記載されているウィーベ関数パラメータの決定方法では、熱発生率が最大となるクランク角度での燃焼割合に基づいて、所定の式によりウィーベ関数の形状パラメータｍを同定するようにしている。そして、ｋ、ａ／θp m+1、θbといった他のパラメータについてもそれぞれ所定の式により同定して、実際の熱発生パターンに高い精度で適合するように、ウィーベ関数を決定することができる。
　前記特許文献１には、そのように複数のパラメータｍ，ｋ，ａ／θp m+1，θbを同定し、ウィーベ関数を決定する作業を種々の運転条件に対して行うことにより、これらのパラメータと内燃機関の運転パラメータ（負荷率、回転速度、空燃比、点火時期など）との関係を把握することができる、と記載されている。また、そうして把握した関係を利用すれば、内燃機関のあらゆる運転条件に対してウィーベ関数を決定し、当該内燃機関の燃焼状態を精度良く表すことができると記載されている。
特開２００７－１７７６５４号公報
　しかしながら前記特許文献１には、ウィーベ関数のパラメータｍ、ｋ、ａ／θp m+1、θbと内燃機関の運転パラメータとの関係を特定する具体的な方法について開示されていない。このため、実際には殆ど全ての運転条件についてパラメータｍ、ｋ、ａ／θp m+1、θbを同定し、それぞれの運転条件についてウィーベ関数を決定しなくてはならない。つまり、前記従来の方法には、熱発生率波形を作成するための工数をさらに低減し、コストを削減する余地が残されている。
　また、前記の方法では、パラメータｍ、ｋ、ａ／θp m+1、θbをそれぞれ同定し、ウィーベ関数を決定して初めて熱発生率波形全体を表現することができ、これに基づいて燃焼状態を評価できるものである。よって、熱発生率波形全体を表現することなく、例えば点火プラグによる点火が行われた後に熱発生率波形が立ち上がるまでの期間（点火プラグによる点火が行われてから混合気が着火するまでの期間）である着火遅れ期間だけを簡易な方法で推定し、評価することはできない。
　本発明はかかる諸点に鑑みてなされたものであり、その目的とするところは、気筒内での混合気の状態を表す指標の一つである着火遅れ期間に着目し、熱発生率波形を作成（算出）するための工数を低減すると共に、例えば着火遅れ期間について、所要の精度を確保しながら簡易に推定し、評価できるようにすることにある。
　－発明の解決原理－
　本発明の発明者は、点火栓による点火が行われてから混合気が着火するまでの期間である着火遅れ期間が燃料密度と高い相関を有しており、その着火遅れ期間に対する機関負荷率および点火時期の影響が燃料密度によってまとめて表現できる、という新規な知見を得た。
　かかる新規な知見に基づいて、本発明の解決原理は、熱発生率波形の特性値の一つとして前記着火遅れ期間を用い、この着火遅れ期間を燃料密度に基づいて推定するようにしたことにある。
　－解決手段－
　具体的に、本発明は、火花点火式の内燃機関における熱発生率波形を算出するための装置を対象とするもので、点火栓による点火が行われてから混合気が着火するまでの期間を、前記熱発生率波形の特性値の一つである着火遅れ期間として規定する。そして、前記混合気の着火時期がピストンの圧縮上死点よりも進角側である場合には、前記点火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定する一方、前記混合気の着火時期がピストンの圧縮上死点よりも遅角側である場合には、前記着火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定し、この推定した着火遅れ期間を用いて前記熱発生率波形を算出する構成とする。
　前記の特定事項によると、内燃機関において気筒内の混合気の燃焼による熱発生率の波形を算出するに際して、その熱発生率波形の特性値の一つとして、点火栓による点火が行われてから混合気が着火するまでの着火遅れ期間を用いている。この着火遅れ期間は、内燃機関の負荷率や点火時期等の運転条件によって変化するものであるが、前述したように、機関負荷率（燃料噴射量を規定するパラメータ）および点火時期（筒内体積を規定するパラメータ）の影響については、燃料密度という一つのパラメータによって、まとめて表現することができる。
　そこで、その燃料密度に基づいて着火遅れ期間を推定することにより、機関負荷率および点火時期の両方に基づいて推定するのに比べて、着火遅れ期間の推定にかかる工数を低減できる。そして、こうして推定した着火遅れ期間を用いて熱発生率波形を作成する工数も低減できる。
　しかも、熱発生率波形全体を作成することなく、前記のように燃料密度に基づいて着火遅れ期間だけを推定できるので、この着火遅れ期間について、所要の精度を確保しながら従来よりも簡易に推定し、評価することができる。
　また、着火遅れ期間の推定に際し、混合気の着火時期がピストンの圧縮上死点よりも進角側である場合（圧縮上死点よりも進角側であるとした場合）には、点火時期における気筒内の燃料密度に基づいて着火遅れ期間を推定する。一方、混合気の着火時期がピストンの圧縮上死点よりも遅角側である場合（圧縮上死点よりも遅角側であるとした場合）には、着火時期における気筒内の燃料密度に基づいて着火遅れ期間を推定する。これは、混合気の着火時期がピストンの圧縮上死点よりも進角側である場合、混合気の着火後には、筒内体積が小さくなっていき、これに伴って燃料密度が大きくなっていくのに対し、混合気の着火時期がピストンの圧縮上死点よりも遅角側である場合、混合気の着火後には、筒内体積が大きくなっていき、これに伴って燃料密度が小さくなっていくことを考慮したものである。つまり、混合気の着火時期がピストンの圧縮上死点よりも進角側である場合には、点火時期における気筒内の燃料密度と着火遅れ期間との間に高い相関があり、混合気の着火時期がピストンの圧縮上死点よりも遅角側である場合には、着火時期における気筒内の燃料密度と着火遅れ期間との間に高い相関があるといった新たに得られた知見に基づいて、着火遅れ期間の推定手法を着火時期に応じて異ならせたものである。
　また、着火遅れ期間を推定する際に、機関回転速度に基づく補正係数（例えば機関回転速度の指数関数）を乗算するようにしてもよい。すなわち、機関回転速度が変化すると一般的に気筒内の流動の強さが変化することから、乱れの影響によって着火遅れ期間が変化することになる。そこで、機関回転速度に基づく補正を行うことにより、着火遅れ期間をより精度良く推定することができる。
　前記着火遅れ期間の算出手法としては、仮想の着火時期を設定し、この仮想の着火時期に従って推定された（例えば演算式によって算出された）推定着火遅れ期間が、実際の点火時期から仮想の着火時期までの期間に一致したか否かを、仮想の着火時期を変更しながら繰り返し計算することが挙げられる。つまり、仮想の着火時期を設定し、この仮想の着火時期がピストンの圧縮上死点よりも進角側である場合には、前記点火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定する一方、前記仮想の着火時期がピストンの圧縮上死点よりも遅角側である場合には、前記着火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定する。そして、この推定した着火遅れ期間と、実際の点火時期と前記仮想の着火時期との間の期間である仮想の着火遅れ期間とを比較し、これらが一致した場合における前記推定した着火遅れ期間が真の着火遅れ期間として算出される。そして、この真の着火遅れ期間を用いて前記熱発生率波形を算出する。
　これによれば、着火遅れ期間を推定するための着火時期を繰り返し計算によって正確に求めることが可能になり、着火遅れ期間を高い精度で算出することができる。
　前記のようにして算出した着火遅れ期間を用いて算出する熱発生率波形としては、例えば、混合気の着火から燃焼終了までのクランク角度期間を底辺とし、熱発生率最大時期における熱発生率を頂点とする三角波形が挙げられる。この三角波形によって熱発生率波形を近似すれば、点火栓による点火時期から三角波形の斜辺の立ち上がり時期までの期間が前記着火遅れ期間として規定されることになる。
　その場合に好ましいのは、混合気の着火時期から熱発生率最大時期までの期間（前半燃焼期間と呼ばれる）が、機関負荷率、空燃比、ＥＧＲ率および油水温のうち少なくとも一つに依らず、主に所定の熱発生率に達する時期（より詳しくは、熱発生率最大時期におけるクランク角度位置での筒内体積；気筒内の乱れに相関のあるパラメータ）および機関回転速度（同じく、気筒内の乱れに相関のあるパラメータ）に依って決まるものとして、前記三角波形を作成することである。すなわち、機関負荷率、空燃比、ＥＧＲ率および油水温が変化しても、前半燃焼期間は変化せず、この前半燃焼期間の変化分は気筒内の乱れの影響分であるとして前記三角波形を作成することができる。こうすることで、熱発生率波形を作成する工数のさらなる低減が図られる。
　見方を変えると本発明は、火花点火式の内燃機関における熱発生率波形を算出する方法に係る。すなわち、まず、点火栓による点火が行われてから混合気が着火するまでの期間を、前記熱発生率波形の特性値の一つである着火遅れ期間として規定する。そして、前記混合気の着火時期がピストンの圧縮上死点よりも進角側である場合には、前記点火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定する一方、前記混合気の着火時期がピストンの圧縮上死点よりも遅角側である場合には、前記着火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定し、この推定した着火遅れ期間を用いて熱発生率波形を算出するものである。
　本発明によると、内燃機関における熱発生率波形の特性値の一つとして、点火栓による点火が行われてから混合気が着火するまでの期間である着火遅れ期間を用い、この着火遅れ期間を気筒内の燃料密度に基づいて推定するようにしたので、熱発生率波形の作成にかかる工数を低減できると共に、熱発生率波形の全体を作成することなく、着火遅れ期間について所要の精度を確保しながら、従来よりも簡易に推定し、評価できるようになる。
実施形態に係る熱発生率波形算出装置の構成、および、この熱発生率波形算出装置の入出力情報を表す図である。 熱発生率波形算出装置から出力される熱発生率波形の一例を示す図である。 熱発生率波形算出装置において行われる熱発生率波形の作成手順を示すフローチャート図である。 ＢＴＤＣ着火の場合における、点火時期ＳＡでの筒内の燃料密度ρfuel@SAの変化に対する着火遅れ期間τの変化を実験によって計測した結果を示す図である。 式（１）で算出された予測着火遅れ期間と、実機において計測された実測着火遅れ期間との関係を検証した結果を示す図である。 ＡＴＤＣ着火の場合における、着火時期ＦＡでの筒内の燃料密度ρfuel@FAの変化に対する着火遅れ期間τの変化を実験によって計測した結果を示す図である。 式（２）で算出された予測着火遅れ期間と、実機において計測された実測着火遅れ期間との関係を検証した結果を示す図である。 ＢＴＤＣ着火する場合の点火時期ＳＡと熱発生率波形とを示す図である。 ＡＴＤＣ着火する場合の点火時期ＳＡと熱発生率波形とを示し、図９（ａ）は点火時期ＳＡがＢＴＤＣである場合を、図９（ｂ）は点火時期ＳＡがＡＴＤＣである場合をそれぞれ示す図である。 負荷率のみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示した図である。 ＥＧＲ率のみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示した図である。 空燃比のみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示した図である。 油水温のみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示した図である。 点火時期ＳＡのみが互いに異なる各エンジン運転状態において得られる熱発生率波形を重ねて表示した図である。 エンジン回転速度Ｎｅのみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示した図である。 あるエンジンに対して、式（３）で算出された予測前半燃焼期間と、実機において計測された実測前半燃焼期間との関係を検証した結果を示す図である。 他のエンジンに対して、式（３）で算出された予測前半燃焼期間と、実機において計測された実測前半燃焼期間との関係を検証した結果を示す図である。 負荷率のみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示した図である。 点火時期ＳＡのみが互いに異なる各エンジン運転状態において得られる熱発生率波形を重ねて表示した図である。 熱発生率最大時燃料密度ρfuel@dQpeakと熱発生率傾きｂ／ａとの関係を、互いに異なるエンジン回転速度Ｎｅそれぞれに対して調べた実験の結果を示す図である。
　以下、本発明の実施の形態を図面に基づいて説明する。本実施形態では、自動車用のガソリンエンジン（火花点火機関）を対象とした熱発生率波形を算出（作成）する熱発生率波形算出装置に本発明を適用した場合について説明する。
　図１は、本実施形態に係る熱発生率波形算出装置１の構成、および、この熱発生率波形算出装置１の入出力情報を表す図である。この熱発生率波形算出装置１は、エンジンの状態量、制御パラメータの制御量および物理量の各種情報が入力される。これら入力情報としては、エンジン回転速度、負荷率、点火時期、ＥＧＲ率、空燃比、油水温、吸排気の各バルブの開閉タイミング（バルブタイミング）等が挙げられる。また、熱発生率波形算出装置１は、各入力情報に基づいて以下の各推定モデルを格納した推定部２～５により熱発生率波形の各種特性値を推定して、この各種特性値を利用して作成した熱発生率波形を出力する。
　－熱発生率波形の各特性値の推定部－
　熱発生率波形算出装置１は、熱発生率波形の特性値として着火遅れ、前半燃焼期間、熱発生率傾きおよび熱発生量をそれぞれ推定するために、着火遅れ推定モデルを格納した着火遅れ推定部２、前半燃焼期間推定モデルを格納した前半燃焼期間推定部３、熱発生率傾き推定モデルを格納した熱発生率傾き推定部４、および、熱発生量推定モデルを格納した熱発生量推定部５を備えている。
　着火遅れ推定部２は、エンジンの点火プラグ（点火栓）によって混合気に点火された時期（以下、点火時期という）、即ち点火プラグの電極間に火花放電が行われてから、この火花によって混合気に着火し初期火炎核が形成される時期（以下、着火時期という）までの期間（以下、着火遅れ期間という）を、着火遅れ推定モデルを使用して推定する部分である。この着火遅れ期間はクランク角度［ＣＡ］で表される。なお、前記着火時期の定義として、本実施形態では、点火時期以降に熱発生率（クランクシャフトの回転の単位クランク角度当たりの熱発生量）が１［Ｊ／ＣＡ］に達した時期とする。この値はこれに限らず適宜設定可能であり、例えば、点火時期以降の熱発生量が総熱発生量に対して所定割合（例えば５％）に達した時期を着火時期としてもよい。さらに、総熱発生量に対する熱発生量の割合が所定値に達した時期（例えば１０％に達した時期でのクランク角度位置）と前記熱発生量の割合が他の所定値に達した時期（例えば５０％に達した時期でのクランク角度位置）とに基づいて着火時期を定義するようにしてもよい。つまり、これらのクランク角度位置と熱発生量の割合とによって、熱発生率が上昇していく期間での熱発生率波形に近似した三角形（三角波形）を作成し、この三角波形に基づいて着火時期を定義するものである。また、熱発生率が上昇していく期間での一般的な熱発生率波形の形状を、前記クランク角度位置と熱発生量の割合との関係が成立するように当て嵌めて熱発生率波形を作成し、この熱発生率波形に基づいて着火時期を定義するようにしてもよい。前記各値はこれに限定されることなく適宜設定可能である。
　前半燃焼期間推定部３は、混合気の燃焼期間のうち、前記着火時期から、火炎核の成長に伴って熱発生率が最大となるタイミング（点火時期から燃焼終了時期までの期間中で熱発生率が最大となるタイミング）までの期間である前半燃焼期間を、前半燃焼期間推定モデルを使用して推定する部分である。以下、この熱発生率が最大となるタイミングを熱発生率最大時期という。また、この熱発生率最大時期および前半燃焼期間はそれぞれクランク角度［ＣＡ］で表される。
　熱発生率傾き推定部４は、前記前半燃焼期間、即ち着火時期から熱発生率最大時期までの期間における、クランク角度変化に対する平均的な熱発生率の増加率（熱発生率の傾き）を、熱発生率傾き推定モデルを使用して推定する部分である。すなわち、本実施形態では、図２を参照して以下に述べるように、熱発生率波形に近似した三角波形を作成するようにしており、熱発生率傾き推定部４は、その三角波形における着火時期から熱発生率最大時期までの熱発生率を表す斜辺の傾きを推定するものである。この熱発生率の傾きの単位としては［Ｊ／ＣＡ2］で表される。
　熱発生量推定部５は、混合気の燃焼によって発生した熱発生量（燃焼期間の全期間において発生した熱発生量であって、点火時期から燃焼終了時期までの期間中における熱発生率の積算値）を、熱発生量推定モデルを使用して推定する部分である。この熱発生量の単位としては［Ｊ］で表される。
　各推定部２～５での推定動作によってそれぞれ、着火遅れ、前半燃焼期間、熱発生率傾き、および、熱発生量という熱発生率波形の特性値が求められ、これらの特性値を利用して熱発生率波形が作成される。そして、この作成された熱発生率波形が熱発生率波形算出装置１の出力となる。
　このため、本実施形態に係る熱発生率波形算出装置１にあっては、図３に示すフローチャートの如く、着火遅れ推定部２における着火遅れ期間の推定動作（ステップＳＴ１）、前半燃焼期間推定部３における前半燃焼期間の推定動作（ステップＳＴ２）、熱発生率傾き推定部４における熱発生率傾きの推定動作（ステップＳＴ３）、熱発生量推定部５における熱発生量の推定動作（ステップＳＴ４）が順に行われた後、これら推定された特性値を利用した熱発生率波形の作成動作（ステップＳＴ５）が行われることになる。
　図２は、各推定部２～５で推定された特性値を利用して作成されて、熱発生率波形算出装置１から出力される熱発生率波形の一例を示している。この図２では、図中のタイミングＳＡが点火時期であり、図中のタイミングＦＡが着火時期である。このため、図中のτが着火遅れ期間となる。また、図中においてｄＱｐｅａｋＡが熱発生率最大時期であって、この熱発生率最大時期ｄＱｐｅａｋＡにおける熱発生率が図中のｂである。つまり、この熱発生率ｂが燃焼期間中における最大熱発生率となっている。また、着火時期ＦＡから熱発生率最大時期ｄＱｐｅａｋＡまでの期間である図中のａが前半燃焼期間となっている。このため、前半燃焼期間ａにおける熱発生率の傾きはｂ／ａとして表される。さらに、熱発生率最大時期ｄＱｐｅａｋＡから燃焼終了時期ＥＡまでの期間である図中のｃが後半燃焼期間となっている。また、図中のＱ１は前半燃焼期間ａにおける熱発生量であり、Ｑ２は後半燃焼期間ｃにおける熱発生量である。そして、燃焼期間の全期間において発生する熱発生量（総熱発生量Ｑall）は、これら熱発生量Ｑ１と熱発生量Ｑ２との和として表される。
　言い換えると、本実施形態の熱発生率波形算出装置１は、混合気の着火から燃焼終了までのクランク角度期間（図中のＦＡからＥＡ）を底辺とし、熱発生率最大時期ｄＱｐｅａｋＡにおける熱発生率ｂを頂点とする三角波形によって熱発生率波形を近似するようにしている。本実施形態では、この熱発生率波形算出装置１の出力である熱発生率波形を利用し、エンジン設計に際してのシステムの検討、制御の検討、適合値の検討が行われる。
　以下、各推定部２～５での推定処理について具体的に説明する。
　－着火遅れ推定部－
　着火遅れ推定部２は、前述した如く、点火時期ＳＡから着火時期ＦＡまでの期間である着火遅れ期間τを推定する部分である。
　この着火遅れ推定部２において行われる着火遅れ期間τの推定処理は以下のとおりである。
　この着火遅れ期間τは、以下の式（１）および式（２）の何れか一方を利用して推定される（これらの式が着火遅れ推定モデルに相当する）。
Figure JPOXMLDOC01-appb-M000001
Figure JPOXMLDOC01-appb-M000002
 　ρfuel@SAは点火時期ＳＡにおける筒内の燃料密度（筒内燃料量［ｍｏｌ］／点火時期の筒内体積［Ｌ］）である。ρfuel@FAは着火時期ＦＡにおける筒内の燃料密度（筒内燃料量［ｍｏｌ］／着火時期の筒内体積［Ｌ］）である。Ｎｅはエンジン回転速度である。Ｃ1，Ｃ2，χ，δ，φ，ψはそれぞれ実験等に基づいて同定される係数である。
　これら式（１）および式（２）は、空燃比が理論空燃比であり、ＥＧＲ率が「０」であり、エンジンの暖機運転が完了しており（油水温が所定値以上であり）、吸気バルブの開閉タイミングが固定されていることを条件に成立する式となっている。
　式（１）は、ピストンが圧縮上死点に達するタイミング（ＴＤＣ）よりも進角側（ＢＴＤＣ）で混合気が着火した（以下、ＢＴＤＣ着火という）場合の着火遅れ期間τの算出式である。また、式（２）は、ピストンが圧縮上死点に達したタイミング（ＴＤＣ）よりも遅角側（ＡＴＤＣ）で混合気が着火した（以下、ＡＴＤＣ着火という）場合の着火遅れ期間τの算出式である。
　これらの式に示すように、着火遅れ期間τは、所定タイミングでの筒内の燃料密度ρfuelおよびエンジン回転速度Ｎｅを変数とする演算式によって算出される。
　これら演算式によって着火遅れ期間τが算出可能である根拠について以下に説明する。
　図４は、ＢＴＤＣ着火の場合における、点火時期ＳＡでの筒内の燃料密度ρfuel@SAの変化に対する着火遅れ期間τの変化を実験によって計測した結果を示すグラフである。この実験は、空燃比を理論空燃比とし、ＥＧＲ率を「０」とし、エンジンの暖機運転が完了しており（油水温が所定値以上であり）、吸気バルブの開閉タイミングを固定して行われたものである。また、この図４では、「○」「△」「□」「◇」「×」「＋」「▽」の順でエンジン回転速度Ｎｅが高くなっている。例えば「○」は８００ｒｐｍ、「△」は１０００ｒｐｍ、「□」は１２００ｒｐｍ、「◇」は１６００ｒｐｍ、「×」は２４００ｒｐｍ、「＋」は３２００ｒｐｍ、「▽」は３６００ｒｐｍである。
　この図４に示すように、ＢＴＤＣ着火した場合、点火時期ＳＡにおける筒内の燃料密度ρfuel@SAと着火遅れ期間τとの間にはエンジン回転速度Ｎｅ毎に相関がある。つまり、これらの相関を概ね１本の曲線で表すことができる。図４では、エンジン回転速度Ｎｅが１０００ｒｐｍの場合および２４００ｒｐｍの場合のそれぞれについて、点火時期ＳＡにおける筒内の燃料密度ρfuel@SAと着火遅れ期間τとの相関を１本の曲線で表している。
　図４に示すように、点火時期ＳＡにおける筒内の燃料密度ρfuel@SAが高いほど着火遅れ期間τは短くなっている。これは、燃料密度ρfuel@SAが高いほど、点火プラグ周辺の燃料分子の数が多くなり、点火プラグの点火後における火炎核の成長が急速に行われるためであると考えられる。また、エンジン回転速度Ｎｅは着火遅れ期間τに影響を与えている。つまり、エンジン回転速度Ｎｅが高いほど着火遅れ期間τは短くなっている。これは、エンジン回転速度Ｎｅが高いほど気筒内での混合気流れの乱れ（以下、単に乱れという）が強くなって前記火炎核の成長が急速に行われるためであると考えられる。このように、点火時期ＳＡにおける筒内の燃料密度ρfuel@SAおよびエンジン回転速度Ｎｅは着火遅れ期間τに影響を与えるパラメータとなっている。
　図５は、式（１）で算出された予測着火遅れ期間と、実機において計測された実測着火遅れ期間との関係を検証した結果を示すグラフである。この予測着火遅れ期間を求めるに当たっては、式（１）におけるＣ1，χ，δの各係数をエンジン運転条件に応じて同定することにより得られた予測式を使用している。この図５では、「○」「△」「□」「◇」「×」「＋」「▽」「☆」の順でエンジン回転速度Ｎｅが高くなっている。例えば「○」は８００ｒｐｍ、「△」は１０００ｒｐｍ、「□」は１２００ｒｐｍ、「◇」は１６００ｒｐｍ、「×」は２０００ｒｐｍ、「＋」は２４００ｒｐｍ、「▽」は３２００ｒｐｍ、「☆」は３６００ｒｐｍである。
　この図５から明らかなように、予測着火遅れ期間は実測着火遅れ期間に略一致しており、式（１）によって、ＢＴＤＣ着火した場合の着火遅れ期間が高い精度で算出されていることが解る。
　図６は、ＡＴＤＣ着火の場合における、着火時期ＦＡでの筒内の燃料密度ρfuel@FAの変化に対する着火遅れ期間τの変化を実験によって計測した結果を示すグラフである。この実験は、エンジン回転速度を固定し、空燃比を理論空燃比とし、ＥＧＲ率を「０」とし、エンジンの暖機運転が完了しており（油水温が所定値以上であり）、吸気バルブの開閉タイミングを固定して行われたものである。また、この図６では、「○」「×」「＋」「△」の順で負荷率が高くなっている。例えば「○」は負荷率２０％、「×」は負荷率３０％、「＋」は負荷率４０％、「△」は負荷率５０％である。
　この図６に示すように、ＡＴＤＣ着火した場合、着火時期ＦＡにおける筒内の燃料密度ρfuel@FAと着火遅れ期間τとの間には負荷率に依らず（負荷率に関わりなく）相関がある。つまり、これらの相関を概ね１本の曲線で表すことができる。
　図６に示すように、着火時期ＦＡにおける筒内の燃料密度ρfuel@FAが高いほど着火遅れ期間τは短くなっている。これは、前述したように、燃料密度ρfuel@FAが高いほど、点火プラグ周辺の燃料分子の数が多くなり、点火プラグの点火後における火炎核の成長が急速に行われるためであると考えられる。このように、着火時期ＦＡにおける筒内の燃料密度ρfuel@FAは着火遅れ期間τに影響を与えるパラメータとなっている。また、前述した場合と同様に、エンジン回転速度Ｎｅも着火遅れ期間τに影響を与えるパラメータになっていると想定される。
　図７は、式（２）で算出された予測着火遅れ期間と、実機において計測された実測着火遅れ期間との関係を検証した結果を示すグラフである。この予測着火遅れ期間を求めるに当たっては、式（２）におけるＣ2，φ，ψの各係数をエンジン運転条件に応じて同定することにより得られた予測式を使用している。この図７では、「○」「×」「＋」「△」の順でエンジン回転速度Ｎｅが高くなっている。例えば「○」は８００ｒｐｍ、「×」は１２００ｒｐｍ、「＋」は３６００ｒｐｍ、「△」は４８００ｒｐｍである。
　この図７から明らかなように、予測着火遅れ期間は実測着火遅れ期間に略一致しており、式（２）によって、ＡＴＤＣ着火した場合の着火遅れ期間が高い精度で算出されていることが解る。
　本発明の発明者は、これらの新たな知見に基づいて前記式（１）および式（２）を導き出した。
　以下、着火遅れ期間τの算出に当たり、着火時期で層別する理由について説明する。つまり、ＢＴＤＣ着火した場合と、ＡＴＤＣ着火した場合とを層別し、それぞれに対して異なる演算式（前記式（１）および式（２））を用いて着火遅れ期間τを算出する理由について説明する。
　まず、ＢＴＤＣ着火する場合、図８（点火時期ＳＡと熱発生率波形とを示す図）に示すように、点火時期ＳＡも、ピストンが圧縮上死点に達するタイミングよりも進角側（ＢＴＤＣ）である。この場合、点火時期ＳＡを迎えた後、ピストンは圧縮上死点に向かって移動することになる。つまり、筒内体積が小さくなっていき、これに伴って燃料密度ρfuelは大きくなっていく。このため、燃料密度ρfuelとしては、着火時期ＦＡでの燃料密度ρfuel@FAよりも点火時期ＳＡでの燃料密度ρfuel@SAの方が小さくなっている。そして、着火遅れ期間の最大値（想定される着火遅れ期間のうち最も長い着火遅れ期間）に相関のある値である点火時期ＳＡでの燃料密度ρfuel@SAに対して、予め同定された各種係数を乗算することにより、着火遅れ期間τを高い精度で求めることができる。
　一方、ＡＴＤＣ着火する場合、図９（点火時期ＳＡと熱発生率波形とを示す図）に示すように、点火時期ＳＡは、ピストンが圧縮上死点に達するタイミングよりも進角側（ＢＴＤＣ）である場合（図９（ａ）を参照）と、遅角側（ＡＴＤＣ）である場合（図９（ｂ）を参照）とがある。これらの場合、着火時期ＦＡを迎えた後、ピストンは下死点に向かって移動することになる。つまり、筒内体積が大きくなっていき、これに伴って燃料密度ρfuelは小さくなっていく。このため、燃料密度ρfuelとしては、点火時期ＳＡでの燃料密度ρfuel@SAよりも着火時期ＦＡでの燃料密度ρfuel@FAの方が小さくなっている可能性が高い。そして、着火遅れ期間の最大値（想定される着火遅れ期間のうち最も長い着火遅れ期間）に相関のある値である着火時期ＦＡでの燃料密度ρfuel@FAに対して、予め同定された各種係数を乗算することにより、着火遅れ期間τを高い精度で求めることができる。
　また、これら式（１）および式（２）の何れを使用するかを判断する手順（着火時期がＢＴＤＣ着火であるかＡＴＤＣ着火であるかを判断する手順）および着火遅れ期間（後述する真の着火遅れ期間）を算出する手順は以下のとおりである。まず、仮想の着火時期を設定し、この仮想の着火時期における筒内体積を求める。この筒内体積は仮想の着火時期に対応するクランク角度位置（ピストンの位置）によって幾何学的に求めることが可能であるので、仮想の着火時期から筒内体積は一義的に決まる。そして、この筒内体積と燃料噴射量とから燃料密度を求める。そして、仮想の着火時期をＢＴＤＣ着火として設定した場合には、この仮想の着火時期における燃料密度およびエンジン回転速度を式（１）に代入して推定着火遅れ期間を算出する。一方、仮想の着火時期をＡＴＤＣ着火として設定した場合には、この仮想の着火時期における燃料密度およびエンジン回転速度を式（２）に代入して推定着火遅れ期間を算出する。そして、仮想の着火時期に対し、前記算出された推定着火遅れ期間だけ進角した時期を仮想の点火時期とする。ここで、この仮想の点火時期と実際の点火時期（入力情報である点火時期）とを比較する。仮想の点火時期と実際の点火時期とが一致しない場合、前記仮想の着火時期を変更する。例えば仮想の着火時期を遅角側に変更する。そして、再度、この仮想の着火時期における燃料密度およびエンジン回転速度を式（１）または式（２）に代入して（仮想の着火時期をＢＴＤＣ着火として設定した場合には式（１）に代入し、仮想の着火時期をＡＴＤＣ着火として設定した場合には式（２）に代入する）推定着火遅れ期間を算出して、仮想の点火時期を求め、これと実際の点火時期（入力情報である点火時期）とを比較する。この動作を繰り返し、仮想の点火時期と実際の点火時期とが一致した場合における仮想の着火時期が真の着火時期として得られることになる。そして、この際（真の着火時期が得られた際）、式（１）または式（２）において算出された推定着火遅れ期間が真の着火遅れ期間として得られることになる。また、前記真の着火時期がＢＴＤＣであった場合（ＢＴＤＣ着火であった場合）に前記求められた着火時期を改めて式（１）に代入して着火遅れ期間τの算出を行い、前記真の着火時期がＡＴＤＣであった場合（ＡＴＤＣ着火であった場合）に前記求められた着火時期を改めて式（２）に代入して着火遅れ期間τの算出を行うようにしてもよい。
　以上の動作を言い換えると以下のとおりとなる。前記実際の点火時期と仮想の着火時期との間の期間（仮想の着火時期で着火したとした場合の仮想の着火遅れ期間）と、式（１）または式（２）で算出された（推定された）推定着火遅れ期間とを比較し、これらが一致しない場合には、前記仮想の着火時期を変更する。そして、再度、式（１）または式（２）によって推定着火遅れ期間を算出した後、前記実際の点火時期と仮想の着火時期との間の期間（仮想の着火遅れ期間）と、式（１）または式（２）で算出された推定着火遅れ期間とを比較する。この動作を繰り返し、これらが一致した場合（仮想の着火遅れ期間と推定着火遅れ期間とが一致した場合）における推定着火遅れ期間が真の着火遅れ期間として得られることになる。
　以上のような着火遅れ推定部２による着火遅れ期間τの推定により、エンジン運転領域の全領域に対して着火遅れ期間τの推定が可能になる。
　以上のようにして着火遅れ期間τが算出されると、前記点火時期ＳＡに着火遅れ期間τを加算することで着火時期ＦＡを求めることができる。
　－前半燃焼期間推定部－
　前半燃焼期間推定部３は、前述した如く、着火時期ＦＡから熱発生率最大時期ｄＱｐｅａｋＡまでの期間である前半燃焼期間ａを推定する部分である。
　この前半燃焼期間推定部３において行われる前半燃焼期間ａの推定処理は以下のとおりである。
　この前半燃焼期間ａ［ＣＡ］は、以下の式（３）を利用して推定される（この式が前半燃焼期間推定モデルに相当する）。
Figure JPOXMLDOC01-appb-M000003
 　Ｖ@dQpeakは前記熱発生率最大時期ｄＱｐｅａｋＡにおける物理量としての筒内体積［Ｌ］であり、以下では、熱発生率最大時筒内体積ともいう。Ｎｅはエンジン回転速度（機関回転速度）である。
　この式（３）は、吸気バルブの開閉タイミングが固定されていることを条件に成立する式となっている。また、この式（３）は、負荷率、ＥＧＲ率、空燃比、油水温の影響を受けることなく成立するものとなっている。つまり、式（３）は、前半燃焼期間ａが、負荷率、ＥＧＲ率、空燃比、油水温の影響を受けないことに基づいて成立している。
　この式（３）によって前半燃焼期間ａが算出可能である根拠について以下に説明する。
　図１０～図１３それぞれは、互いに異なるエンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示したものである。図１０は、負荷率のみが互いに異なる各エンジン運転状態において得られる熱発生率波形を重ねて表示したものである。図１１は、ＥＧＲ率のみが互いに異なる各エンジン運転状態において得られる熱発生率波形を重ねて表示したものである。図１２は、空燃比のみが互いに異なる各エンジン運転状態において得られる熱発生率波形を重ねて表示したものである。また、図１３は、エンジンの暖機運転の途中のように油水温のみが互いに異なる場合に得られる熱発生率波形を重ねて表示したものである。
　これら図１０～図１３に示すように、負荷率、ＥＧＲ率、空燃比、油水温の何れが変化しても前半燃焼期間ａは一定に維持されている。つまり、前半燃焼期間ａは、負荷率、ＥＧＲ率、空燃比、油水温の影響を受けないものであることが解る。
　一方、図１４は、点火時期ＳＡのみが互いに異なる各エンジン運転状態において得られる熱発生率波形を重ねて表示したものである。この図１４から解るように、点火時期ＳＡが遅角されるほど前半燃焼期間ａは長くなっている。
　また、図１５は、エンジン回転速度Ｎｅのみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡを互いに一致させるように点火時期ＳＡを調整した各熱発生率波形を重ねて表示したものである。エンジン回転速度Ｎｅが高いほど単位時間［ｍｓ］当たりにおけるクランクの回転角度［ＣＡ］は大きくなるので、その分、前半燃焼期間ａは長くなる（クランク角度軸上で長くなる）はずであるが、図１５に示すものでは、エンジン回転速度Ｎｅが異なっても前半燃焼期間ａは殆ど変化していない。これは、エンジン回転速度Ｎｅが高いほど前半燃焼期間ａが短くなる要因が存在しているものと考えられる。つまり、エンジン回転速度Ｎｅが高いほど単位時間当たりにおけるクランクの回転角度が大きくなることに起因して前半燃焼期間ａが長くなることとは別に、「他の要因」によって前半燃焼期間ａが短くなっていると想定できる。
　このように前半燃焼期間ａは、点火時期ＳＡおよびエンジン回転速度Ｎｅの影響を受けるものであることが解る。
　前半燃焼期間ａが、点火時期ＳＡおよびエンジン回転速度Ｎｅの影響を受ける要因としては、点火時期ＳＡおよびエンジン回転速度Ｎｅが、気筒内の乱れに影響を与えることが考えられる。
　つまり、前記熱発生率最大時期ｄＱｐｅａｋＡがＴＤＣよりも遅角側にある場合を考えると、点火時期ＳＡが遅角側に移行するほど着火時期ＦＡおよび熱発生率最大時期ｄＱｐｅａｋＡも遅角側に移行し、この熱発生率最大時期ｄＱｐｅａｋＡにおける筒内体積（熱発生率最大時筒内体積Ｖ@dQpeak）が大きくなると共に、気筒内の乱れは弱くなる。そして、気筒内の乱れが弱くなると、火炎伝播が緩慢になって前半燃焼期間ａは長くなる。逆に、点火時期ＳＡが進角側に移行するほど着火時期ＦＡおよび熱発生率最大時期ｄＱｐｅａｋＡも進角側に移行し、熱発生率最大時筒内体積Ｖ@dQpeakが小さくなると共に、気筒内の乱れは強くなる。これにより火炎伝播が急速になって、前半燃焼期間ａは短くなる。
　また、エンジン回転速度Ｎｅが低いほど吸気系から気筒内に流入する空気の流速が低くなって気筒内の乱れは弱くなる。そして、気筒内の乱れが弱くなると、火炎伝播が緩慢になって前半燃焼期間ａは長くなる。逆に、エンジン回転速度Ｎｅが高いほど吸気系から気筒内に流入する空気の流速が高くなって気筒内の乱れは強くなる。そして、気筒内の乱れが強くなると、火炎伝播が急速になって前半燃焼期間ａは短くなる。前述した「他の要因（前半燃焼期間ａを短くする要因）」は、このエンジン回転速度Ｎｅが高いほど気筒内の乱れが強くなることに起因して火炎伝播が急速になっていることである。
　本発明の発明者は、この新たな知見に基づいて前記式（３）を導き出した。そして、この式（３）においては、制御量である点火時期ＳＡに相関のある物理量として筒内体積、特に、熱発生率最大時筒内体積Ｖ@dQpeakを変数として用いている。つまり、前述した如く、点火時期ＳＡが遅角側に移行するほど熱発生率最大時期ｄＱｐｅａｋＡも遅角側に移行し、筒内体積Ｖ@dQpeakが大きくなることから、点火時期ＳＡに相関のある物理量として、熱発生率最大時筒内体積Ｖ@dQpeakを変数として用いている。
　前記式（３）の変数である熱発生率最大時筒内体積Ｖ@dQpeakを求める手順および前半燃焼期間ａを算出する手順は以下のとおりである。まず、仮想の熱発生率最大時期を設定し、この仮想の熱発生率最大時期における筒内体積を求める。この筒内体積は仮想の熱発生率最大時期に対応するクランク角度位置（ピストンの位置）によって幾何学的に求めることが可能であるので、仮想の熱発生率最大時期から筒内体積は一義的に決まる。そして、この仮想の熱発生率最大時期における筒内体積およびエンジン回転速度を式（３）に代入して推定前半燃焼期間を算出する。そして、前記仮想の熱発生率最大時期に対し、この算出した推定前半燃焼期間だけ進角した時期を仮想の着火時期とする。前述した着火遅れ推定部２では着火遅れ期間τが算出されているため、前記点火時期ＳＡに着火遅れ期間τを加算することで着火時期ＦＡが算出できる。ここで、前記仮想の着火時期と算出された着火時期ＦＡとを比較する。仮想の着火時期と算出された着火時期ＦＡとが一致しない場合、前記仮想の熱発生率最大時期を変更する。例えば仮想の熱発生率最大時期を遅角側に変更する。そして、再度、この仮想の熱発生率最大時期における筒内体積およびエンジン回転速度を式（３）に代入して推定前半燃焼期間を算出して、仮想の着火時期を求め、これと算出された着火時期ＦＡ（前記点火時期ＳＡに、着火遅れ推定部２で算出された着火遅れ期間τを加算することで得られる着火時期ＦＡ）とを比較する。この動作を繰り返し、仮想の着火時期と算出された着火時期ＦＡとが一致した場合における仮想の熱発生率最大時期が真の熱発生率最大時期ｄＱｐｅａｋＡとして得られることになる。そして、この際（真の熱発生率最大時期ｄＱｐｅａｋＡが得られた際）、式（３）において算出された推定前半燃焼期間が真の前半燃焼期間として得られることになる。また、前記真の熱発生率最大時期ｄＱｐｅａｋＡにおける熱発生率最大時筒内体積Ｖ@dQpeakを幾何学的に求め、これを改めて式（３）に代入して前半燃焼期間ａの算出を行うようにしてもよい。
　以上の動作を言い換えると以下のとおりとなる。前記着火時期ＦＡ（実際の点火時期に従って求められた着火時期）と仮想の熱発生率最大時期との間の期間（仮想の前半燃焼期間）と、式（３）で算出された（推定された）推定前半燃焼期間（仮想の熱発生率最大時期における物理量に基づいて推定された前半燃焼期間）とを比較し、これらが一致しない場合には、前記仮想の熱発生率最大時期を変更する。そして、再度、式（３）によって推定前半燃焼期間を算出した後、前記着火時期ＦＡと仮想の熱発生率最大時期との間の期間（仮想の前半燃焼期間）と、式（３）で算出された推定前半燃焼期間とを比較する。この動作を繰り返し、これらが一致した場合（仮想の前半燃焼期間と推定前半燃焼期間とが一致した場合）における推定前半燃焼期間が真の前半燃焼期間ａとして得られることになる。
　式（３）における各係数について具体的に説明すると、Ｃおよびαは実験等に基づいた同定によって得られる。また、βは気筒内でのタンブル比に応じた値となっており、タンブル比が大きいほど大きな値として与えられる。なお、βは実験等に基づいた同定によって設定するようにしてもよい。また、これら係数は、吸気バルブの開閉タイミングの変化に対して同定することも可能である。このように、式（３）は、熱発生率最大時筒内体積Ｖ@dQpeakに基づき、且つタンブル比に応じた値βを指数とするエンジン回転速度Ｎｅの指数関数（補正係数）を乗算することにより、前半燃焼期間ａを算出するものとなっている。
　図１６および図１７は、互いに異なるエンジンに対して、式（３）で算出された予測前半燃焼期間と、実機において計測された実測前半燃焼期間との関係を検証した結果を示すグラフである。この予測前半燃焼期間を求めるに当たっては、式（３）における係数Ｃをエンジン運転条件に応じて同定することにより得られた予測式を使用している。図１６では、「○」「△」「□」「◇」「×」「＋」「▽」の順でエンジン回転速度Ｎｅが高くなっている。例えば「○」は８００ｒｐｍ、「△」は１０００ｒｐｍ、「□」は１２００ｒｐｍ、「◇」は１６００ｒｐｍ、「×」は２４００ｒｐｍ、「＋」は３２００ｒｐｍ、「▽」は３６００ｒｐｍである。また、図１７では、「○」「×」「＋」「△」「□」の順でエンジン回転速度Ｎｅが高くなっている。例えば「○」は８００ｒｐｍ、「×」は１２００ｒｐｍ、「＋」は２４００ｒｐｍ、「△」は３６００ｒｐｍ、「□」は４８００ｒｐｍである。
　これら図１６および図１７から明らかなように、予測前半燃焼期間は実測前半燃焼期間に略一致しており、式（３）によって前半燃焼期間ａが高い精度で算出されていることが解る。
　以上、述べたように前半燃焼期間ａを、負荷率、空燃比、ＥＧＲ率、油水温の影響を受けないものとして、熱発生率最大時筒内体積Ｖ@dQpeakとエンジン回転速度Ｎｅとに基づいて推定できる。これら熱発生率最大時筒内体積Ｖ@dQpeakとエンジン回転速度Ｎｅは前述した如く気筒内の乱れに相関のあるパラメータとなっている。言い換えると、負荷率、空燃比、ＥＧＲ率、油水温は、気筒内の乱れに殆ど相関が無いことから前半燃焼期間ａに影響を与えないものと想定される。そして、これら負荷率、空燃比、ＥＧＲ率、油水温を考慮することなしに、気筒内の乱れに相関のあるパラメータである熱発生率最大時筒内体積Ｖ@dQpeakとエンジン回転速度Ｎｅとに基づいて前半燃焼期間ａを推定できることから、エンジンの種々の運転条件における前半燃焼期間ａを決定するための工数が大幅に低減できることになる。
　なお、前述したように、前半燃焼期間は負荷率の影響を受けない。この負荷率は、燃料噴射量を制御するためのパラメータの一つであり、燃料噴射量は気筒内の燃料密度に影響を及ぼす制御パラメータである。このため、前半燃焼期間は、気筒内の燃料密度に依ることなく推定される。具体的には、前述したように、熱発生率最大時筒内体積Ｖ@dQpeakおよびエンジン回転速度Ｎｅといった気筒内の乱れに影響を与えるパラメータに基づいて前半燃焼期間は推定される。これに対し、後述するように、熱発生率傾きは気筒内の燃料密度に基づいて推定される。このように、本実施形態で推定される前半燃焼期間と熱発生率傾きとは、互いに独立した（従属関係にない）値として推定されるものとなっている。
　－熱発生率傾き推定部－
　熱発生率傾き推定部４は、前述した如く、前半燃焼期間ａにおける熱発生率の傾きｂ／ａ（以下、熱発生率傾きという）を推定する部分である。
　この熱発生率傾き推定部４において行われる熱発生率傾きｂ／ａの推定処理は以下のとおりである。
　この熱発生率傾きｂ／ａ［Ｊ／ＣＡ2］は、基本的には以下の式（４）を利用して推定される（この式が熱発生率傾き推定モデルに相当する）。
Figure JPOXMLDOC01-appb-M000004
 　ρfuel@dQpeakは、前記熱発生率最大時期ｄＱｐｅａｋＡにおける燃料密度（筒内燃料量［ｍｏｌ］／熱発生率最大時期の筒内体積［Ｌ］）であり、以下では、熱発生率最大時燃料密度ともいう。Ｃ3は実験等に基づいて同定される係数である。
　この式（４）は、エンジン回転速度が固定されており、空燃比が理論空燃比であり、ＥＧＲ率が「０」であり、エンジンの暖機運転が完了しており（油水温が所定値以上であり）、吸気バルブの開閉タイミングが固定されていることを条件に成立する式となっている。なお、エンジン回転速度、空燃比、ＥＧＲ率、エンジンの油水温等の影響については後述する。
　この式（４）によって熱発生率傾きｂ／ａが算出可能である根拠について以下に説明する。
　図１８（ａ）～図１８（ｄ）それぞれは、負荷率のみが互いに異なる各エンジン運転状態において得られる熱発生率波形であって、熱発生率最大時期ｄＱｐｅａｋＡが互いに一致するように点火時期ＳＡを調整した各熱発生率波形を重ねて表示したものである。図１８（ａ）～図１８（ｄ）の順に点火時期は遅角側に変化しており、また、各図において負荷率は、ＫＬ１，ＫＬ２，ＫＬ３の順に大きくなっている。例えば図１８において、ＫＬ１は負荷率２０％、ＫＬ２は負荷率３０％，ＫＬ３は負荷率４０％である。
　これら図１８（ａ）～図１８（ｄ）に表れているように、熱発生率傾きｂ／ａは負荷率および点火時期ＳＡの影響を受けている。具体的には、それぞれ点火時期ＳＡの異なる図１８（ａ）～図１８（ｄ）のいずれにおいても、負荷率が大きいほど熱発生率傾きｂ／ａは大きくなっている。このように熱発生率傾きｂ／ａが負荷率の影響を受ける要因として、負荷率に応じて気筒内の燃料密度が変化することが考えられる。つまり、負荷率が高いほど筒内燃料量が多くなるので、気筒内の燃料密度も高くなって、混合気の燃焼速度が高くなるものと考えられる。
　また、図１８（ａ）～図１８（ｄ）の順に点火時期ＳＡが遅角側に移行するに連れて、熱発生率傾きｂ／ａが小さくなっている。図１９（ａ）および図１９（ｂ）それぞれは、点火時期ＳＡの変化による影響を調べるために、点火時期ＳＡのみが互いに異なる各エンジン運転状態において得られる熱発生率波形を重ねて表示したものである。これら図１９（ａ）および図１９（ｂ）においてはそれぞれ負荷率が異なっているが、そのいずれにおいても熱発生率傾きｂ／ａは、点火時期ＳＡの遅角側への移行に連れて小さくなる傾向がある。
　このように熱発生率傾きｂ／ａが点火時期ＳＡの影響を受ける要因も、前記した負荷率の場合と同じく気筒内の燃料密度によるものと考えられる。すなわち、ピストンが圧縮上死点（ＴＤＣ）付近にあるときは、クランク角度の変化に伴う筒内体積の変化は小さいが、膨張行程においてＴＤＣから離れてゆくに連れて（例えばＡＴＤＣ１０°ＣＡくらいからは）遅角側ほど筒内体積が大きくなってゆき、これに連れて気筒内の燃料密度は低下してゆく。
　そして、前記図１９（ａ）および図１９（ｂ）に表れているように熱発生率波形は、点火時期ＳＡの遅角に伴い全体として遅角側に移動するとともに、着火時期ＦＡ（波形の始点）がＴＤＣ以降になると、徐々に熱発生率波形の傾きも小さくなってゆく。この結果、着火時期ＦＡ（波形の始点）から熱発生率最大時期ｄＱｐｅａｋＡにおける熱発生率ｂ（波形の頂点）までを結んだ直線（図中には一点鎖線で示す）の傾き、即ち熱発生率傾きｂ／ａも遅角側に向かって徐々に小さくなっている。
　そのように点火時期ＳＡの遅角（即ち着火時期ＦＡの遅角）が熱発生率傾きｂ／ａに及ぼす影響は、この熱発生率傾きｂ／ａと熱発生率最大時燃料密度ρfuel@dQpeakとの関係に顕著に表れている。すなわち、図１９（ａ）および図１９（ｂ）に表れているように、点火時期ＳＡの遅角に伴い熱発生率最大時期ｄＱｐｅａｋＡが遅角側に移行し、この熱発生率最大時期ｄＱｐｅａｋＡにおける筒内体積（熱発生率最大時筒内体積Ｖ@dQpeak）が徐々に大きくなってゆくと、これに応じて熱発生率最大時燃料密度ρfuel@dQpeakが小さくなってゆき、これに対応して熱発生率傾きｂ／ａが小さくなってゆく。
　本発明者は、そのように熱発生率最大時燃料密度ρfuel@dQpeakの変化に対応して熱発生率傾きｂ／ａが変化する様子を調べた。この実験の結果を図２０（ａ）～図２０（ｄ）のグラフに示す。これらの各図において負荷率は、「○」「×」「＋」「△」「□」「◇」「▽」「☆」の順に大きくなっている。例えば図２０において、「○」は負荷率１５％、「×」は負荷率２０％、「＋」は負荷率２５％、「△」は負荷率３０％、「□」は負荷率３５％、「◇」は負荷率４０％、「▽」は負荷率４５％、「☆」は負荷率５０％である。
　また、図２０（ａ）～図２０（ｄ）の順にエンジン回転速度Ｎｅが高くなっており、例えば、図２０（ａ）は８００ｒｐｍ、図２０（ｂ）は１２００ｒｐｍ、図２０（ｃ）は２０００ｒｐｍ、図２０（ｄ）は３２００ｒｐｍである。
　図２０（ａ）～図２０（ｄ）それぞれに表れているように、エンジン回転速度を固定すれば、熱発生率最大時燃料密度ρfuel@dQpeakと熱発生率傾きｂ／ａとの関係は、負荷率や点火時期ＳＡが異なっていても概ね１本の直線として表すことができ、両者の間には高い相関（具体的には概略比例関係）のあることが解る。つまり、熱発生率傾きｂ／ａに対するエンジンの負荷率および点火時期ＳＡの影響は、熱発生率最大時燃料密度ρfuel@dQpeakという一つのパラメータによって、まとめて表現されるものである。
　本発明の発明者は、この新たな知見に基づいて前記式（４）を導き出した。
　前記式（４）の変数である熱発生率最大時燃料密度ρfuel@dQpeakは、前述したように筒内燃料量を熱発生率最大時筒内体積Ｖ@dQpeakによって除算して、求めることができる。熱発生率最大時筒内体積Ｖ@dQpeakを求める手順は、前半燃焼期間推定部３の説明において上述したとおりである。また、筒内燃料量は、熱発生率波形算出装置１の入力情報として与えられる。
　このようにして、熱発生率波形の特性値の一つである熱発生率傾きｂ／ａを、基本的には熱発生率最大時燃料密度ρfuel@dQpeakの一次関数（この実施形態では、一例として比例関数）として算出することができる。言い換えると、負荷率および点火時期ＳＡを考慮せず、主に熱発生率最大時燃料密度ρfuel@dQpeakに基づいて熱発生率傾きｂ／ａを推定できることから、エンジンの種々の運転条件における熱発生率傾きｂ／ａの決定にかかる工数が低減される。
　－熱発生量推定部－
　熱発生量推定部５は、前述した如く、燃焼期間の全期間において発生した熱発生量（総熱発生量Ｑall）を推定する部分である。
　以下、この熱発生量推定部における総熱発生量Ｑallの推定動作について説明する。
　まず、前半燃焼期間ａにおける熱発生量Ｑ１は以下の式（５）によって算出される。
Figure JPOXMLDOC01-appb-M000005
 　一方、総熱発生量Ｑallは、筒内燃料量×ｋ（燃焼効率）として表すことができる（この式が熱発生量推定モデルに相当する）。燃焼効率ｋは、例えば暖機運転中のように油水温が低い場合に低下すると共に、負荷率やエンジン回転速度などの変化によっても変化する。そこで、本実施形態では、種々のエンジンの実験結果のデータベースを利用して燃焼効率ｋの値を、油水温、負荷率およびエンジン回転速度に対応づけて決定するマップを予め決定しておく。そして、この燃焼効率ｋの値を用いて、筒内燃料量から総熱発生量Ｑallが算出される。
　図２を参照して上述したように、熱発生率波形を作成するためには、熱発生率最大時期ｄＱｐｅａｋＡにおける熱発生率ｂおよび後半燃焼期間ｃを求める必要がある。この後半燃焼期間ｃにおける熱発生量Ｑ２は以下の式（６）によって算出される。
Figure JPOXMLDOC01-appb-M000006
 　そして、熱発生率最大時期ｄＱｐｅａｋＡにおける熱発生率ｂは以下の式（７）により、また、後半燃焼期間ｃは以下の式（８）によりそれぞれ求められる。
Figure JPOXMLDOC01-appb-M000007
Figure JPOXMLDOC01-appb-M000008
 　以上のようにして熱発生率波形算出装置１では、着火遅れ推定モデルを使用した着火遅れ推定部２での着火遅れ期間τの推定、前半燃焼期間推定モデルを使用した前半燃焼期間推定部３での前半燃焼期間ａの推定、熱発生率傾き推定モデルを使用した熱発生率傾き推定部４での熱発生率傾きｂ／ａの推定、熱発生量推定モデルを使用した熱発生量推定部５での熱発生量Ｑallの推定が行われると共に、最大熱発生率ｂおよび後半燃焼期間ｃが算出される。そして、これらの算出値を用いて、図２のように熱発生率波形を近似する三角波形を作成し、この三角波形を熱発生率波形として出力する。この出力された熱発生率波形を利用して、エンジン設計に際してのシステムの検討、制御の検討、適合値の検討が行われることになる。
　以上説明したように本実施形態では、エンジンの熱発生率波形に近似した三角波形を作成するに際し、その波形の特性値の一つである着火遅れ期間τを、所定時期における燃料密度ρfuelから算出するようにしている。つまり、混合気の着火時期がピストンの圧縮上死点よりも進角側（ＢＴＤＣ着火）である場合には、点火時期ＳＡにおける気筒内の燃料密度ρfuel@SAに基づいて前記着火遅れ期間τを算出する一方、混合気の着火時期がピストンの圧縮上死点よりも遅角側（ＡＴＤＣ着火）である場合には、着火時期ＦＡにおける気筒内の燃料密度ρfuel@FAに基づいて前記着火遅れ期間τを算出するようにしている。このため、エンジンの負荷率および点火時期の両方に基づいて算出するのに比べて、工数の低減が図られる。
　また、本実施形態では、同じく特性値の一つである前半燃焼期間ａを、エンジンの負荷率、空燃比、ＥＧＲ率および油水温のいずれにも依らないものとして、熱発生率最大時筒内体積Ｖ@dQpeakとエンジン回転速度Ｎｅとに基づいて算出するようにしている。これにより前半燃焼期間ａの算出にかかる工数が大幅に低減できる。
　よって、それら着火遅れτや前半燃焼期間ａの算出値を用いて熱発生率波形を作成するための工数を大幅に低減でき、この熱発生率波形を利用して、エンジン設計に際しての種々の検討などを効率的に行うことによって、開発コストを削減できる。
　また、気筒内の燃料密度ρfuelとエンジン回転速度Ｎｅとに基づいて着火遅れ期間τを算出し、この着火遅れ期間τに基づいて熱発生率波形を作成していることから、この熱発生率波形は、気筒内での燃焼状態における物理現象に従って作成されたものとなっている。この点で、単に波形形状を合わせ込むために、形状パラメータなどの各種パラメータを数学的に合わせ込むといったウィーベ関数を利用した熱発生率波形の作成手法に比べて、本実施形態に係る熱発生率波形算出装置１によって作成された熱発生率波形は高い信頼性を得ることができる。
　さらに、本実施形態では、熱発生率波形全体を作成しなくても、前記のように気筒内の燃料密度ρfuelとエンジン回転速度Ｎｅとから着火遅れ期間τを算出することができるので、この着火遅れ期間τについて所要の精度を確保しながら、従来よりも簡易に推定し、評価することができる。
　また、前述した如く、本実施形態で推定される前半燃焼期間ａと熱発生率傾きｂ／ａとは、互いに独立した（従属関係にない）値として推定されるものとなっている。このため、前半燃焼期間ａは、熱発生率傾きｂ／ａに推定誤差が含まれる場合にその誤差の影響を受けることの無い値として推定され、熱発生率傾きｂ／ａは、前半燃焼期間ａに推定誤差が含まれる場合にその誤差の影響を受けることの無い値として推定される。その結果、これら値の推定精度を高く確保することができる。
　－他の実施形態－
　以上説明した実施形態は、自動車用のガソリンエンジンを対象とした熱発生率波形を作成する熱発生率波形算出装置に本発明を適用した場合について説明した。本発明はこれに限らず、自動車用以外の火花点火機関に対しても適用が可能である。また、ガソリンエンジンにも特に限定されるものではなく、例えばガスエンジンに対しても適用が可能である。
　また、前記の実施形態で説明した熱発生率波形算出装置において実施される熱発生率波形算出方法も本発明の技術的思想の範疇である。
　前記の実施形態では、混合気の着火時期ＦＡから熱発生率最大時期ｄＱｐｅａｋＡまでの期間における熱発生率の平均的な増加率を熱発生率傾きｂ／ａとし、式（４）のように、この熱発生率傾きｂ／ａを熱発生率最大時燃料密度ρfuel@dQpeakの一次関数として算出するようにしているが、これにも限定されない。
　すなわち、着火時期ＦＡから熱発生率最大時期ｄＱｐｅａｋＡまでの間で熱発生率が増加している間（熱発生率増加期間）で、例えば、着火時期から熱発生率最大時期ｄＱｐｅａｋＡの少し手前の所定時期までの期間における熱発生率の増加率を熱発生率傾きとし、この熱発生率傾きを、前記所定時期の燃料密度に基づいて推定するようにしてもよい。
　また、前記の実施形態では、着火時期ＦＡおよび着火遅れ期間τを求める手法として、仮想の着火時期を設定し、前記式（１）または式（２）の繰り返し計算によって算出するようにしていた。本発明はこれに限らず、実機での試験において着火時期をセンシングし、これに基づいて着火時期を設定するようにしたり、熱発生率波形算出装置１への入力信号として所望の着火時期を入力したりすることによって、着火遅れ期間τを求めるようにしてもよい。
　また、前記の実施形態では、熱発生率最大時筒内体積Ｖ@dQpeakおよび前半燃焼期間ａを求める手法として、仮想の熱発生率最大時期を設定し、前記式（３）の繰り返し計算によって算出するようにしていた。これに限らず、実機での試験において熱発生率最大時期をセンシングし、これに基づいて熱発生率最大時を設定するようにしたり、熱発生率波形算出装置１への入力信号として所望の熱発生率最大時期を入力したりすることによって、熱発生率最大時筒内体積Ｖ@dQpeakおよび前半燃焼期間ａを求めるようにしてもよい。
　また、前記の実施形態に係る熱発生率波形算出装置１は、三角波形を出力するものであった。本発明はこれに限らず、作成した三角波形に対して所定のフィルタ処理を行って熱発生率波形を作成し、この熱発生率波形を出力するものであってもよい。
　また、前記の実施形態では、前半燃焼期間ａを、エンジンの負荷率、ＥＧＲ率、空燃比および油水温のいずれにも依らないものとして算出していたが、これら運転条件のうち少なくとも一つに依らないものとして算出するようにしてもよい。
　本発明によると、火花点火式の内燃機関における熱発生率波形の作成にかかる工数を低減でき、コストの削減が可能になるので、例えば自動車用の内燃機関に適用可能である。
　１　　　　　　　　熱発生率波形算出装置
　ＳＡ　　　　　　　点火時期
　ＦＡ　　　　　　　混合気の着火時期
　τ　　　　　　　　着火遅れ期間
　ａ　　　　　　　　前半燃焼期間（着火時期から熱発生率最大時期までの期間）
　ｄＱｐｅａｋＡ　　熱発生率最大時期
　ρfuel@SA　　　   点火時期における筒内の燃料密度
　ρfuel@FA　　　   着火時期における筒内の燃料密度
　Ｎｅδ，Ｎｅψ　　　エンジン回転速度に基づく補正係数

Claims (6)

  1. 　火花点火式の内燃機関における熱発生率波形を算出するための装置であって、
    　点火栓による点火が行われてから混合気が着火するまでの期間を、前記熱発生率波形の特性値の一つである着火遅れ期間として規定し、
    　前記混合気の着火時期がピストンの圧縮上死点よりも進角側である場合には、前記点火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定する一方、前記混合気の着火時期がピストンの圧縮上死点よりも遅角側である場合には、前記着火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定し、この推定した着火遅れ期間を用いて前記熱発生率波形を算出する構成となっていることを特徴とする内燃機関の熱発生率波形算出装置。
  2. 　請求項１に記載の内燃機関の熱発生率波形算出装置において、
    　前記着火遅れ期間を、機関回転速度に基づく補正係数を乗算して算出する構成となっている、内燃機関の熱発生率波形算出装置。
  3. 　請求項１または２のいずれかに記載の内燃機関の熱発生率波形算出装置において、
    　前記着火遅れ期間は、
    　仮想の着火時期を設定し、この仮想の着火時期がピストンの圧縮上死点よりも進角側である場合には、前記点火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定する一方、前記仮想の着火時期がピストンの圧縮上死点よりも遅角側である場合には、前記着火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定し、この推定した着火遅れ期間と、実際の点火時期と前記仮想の着火時期との間の期間である仮想の着火遅れ期間とを比較し、これらが一致した場合における前記推定した着火遅れ期間が真の着火遅れ期間として算出され、
    　この真の着火遅れ期間を用いて前記熱発生率波形を算出する構成となっている、内燃機関の熱発生率波形算出装置。
  4. 　請求項１～３のいずれか一つに記載の内燃機関の熱発生率波形算出装置において、
    　混合気の着火から燃焼終了までのクランク角度期間を底辺とし、熱発生率最大時期における熱発生率を頂点とする三角波形によって熱発生率波形を近似し、
    　前記三角波形において、点火栓による点火時期から三角波形の斜辺の立ち上がり時期までの期間を前記着火遅れ期間として規定する構成となっている、内燃機関の熱発生率波形算出装置。
  5. 　請求項４に記載の内燃機関の熱発生率波形算出装置において、
    　前記三角波形における着火時期から熱発生率最大時期までの期間が、機関負荷率、空燃比、ＥＧＲ率および油水温のうち少なくとも一つに依らないものとして、前記三角波形を作成する構成となっている、内燃機関の熱発生率波形算出装置。
  6. 　火花点火式の内燃機関における熱発生率波形を算出する方法であって、
    　点火栓による点火が行われてから混合気が着火するまでの期間を、前記熱発生率波形の特性値の一つである着火遅れ期間として規定し、
    　前記混合気の着火時期がピストンの圧縮上死点よりも進角側である場合には、前記点火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定する一方、前記混合気の着火時期がピストンの圧縮上死点よりも遅角側である場合には、前記着火時期における気筒内の燃料密度に基づいて前記着火遅れ期間を推定し、この推定した着火遅れ期間を用いて前記熱発生率波形を算出することを特徴とする内燃機関の熱発生率波形算出方法。
PCT/JP2015/053497 2014-04-22 2015-02-09 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法 WO2015162971A1 (ja)

Priority Applications (4)

Application Number Priority Date Filing Date Title
US15/303,331 US9885295B2 (en) 2014-04-22 2015-02-09 Heat generation rate waveform calculation device of internal combustion engine and method for calculating heat generation rate waveform
JP2016514747A JP6260692B2 (ja) 2014-04-22 2015-02-09 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
CN201580020997.XA CN106232970B (zh) 2014-04-22 2015-02-09 内燃机的热产生率波形计算装置及热产生率波形计算方法
EP15782377.4A EP3135887B1 (en) 2014-04-22 2015-02-09 Heat generation rate waveform calculation device of internal combustion engine and method for calculating heat generation rate waveform

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2014088035 2014-04-22
JP2014-088035 2014-04-22

Publications (1)

Publication Number Publication Date
WO2015162971A1 true WO2015162971A1 (ja) 2015-10-29

Family

ID=54332145

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2015/053497 WO2015162971A1 (ja) 2014-04-22 2015-02-09 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法

Country Status (5)

Country Link
US (1) US9885295B2 (ja)
EP (1) EP3135887B1 (ja)
JP (1) JP6260692B2 (ja)
CN (1) CN106232970B (ja)
WO (1) WO2015162971A1 (ja)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10001052B2 (en) 2015-02-26 2018-06-19 Toyota Jidosha Kabushiki Kaisha Combustion state estimation method
CN108368794A (zh) * 2015-12-02 2018-08-03 本田技研工业株式会社 内燃机的控制装置

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106255816B (zh) * 2014-04-22 2019-12-17 丰田自动车株式会社 内燃机的热产生率波形计算装置及热产生率波形计算方法
JP6137220B2 (ja) 2015-02-24 2017-05-31 トヨタ自動車株式会社 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
JP7431512B2 (ja) * 2019-05-23 2024-02-15 日立Astemo株式会社 内燃機関制御装置

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2002322927A (ja) * 2001-04-26 2002-11-08 Nissan Motor Co Ltd 自己着火式内燃機関
JP2003278591A (ja) * 2002-03-20 2003-10-02 Denso Corp 内燃機関のトルク制御装置
JP2005048748A (ja) * 2003-07-31 2005-02-24 Nissan Motor Co Ltd 内燃機関の燃焼制御装置
JP2008274905A (ja) * 2007-05-07 2008-11-13 Toyota Motor Corp 多種燃料内燃機関の燃焼制御装置
JP2011080395A (ja) * 2009-10-06 2011-04-21 Cosmo Oil Co Ltd 圧縮自己着火内燃機関の燃焼時期予測方法、該圧縮自己着火内燃機関の制御方法及び圧縮自己着火内燃機関システム

Family Cites Families (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP3823553B2 (ja) * 1998-08-12 2006-09-20 株式会社日立製作所 エンジン燃焼制御装置
US7021286B2 (en) * 2004-02-20 2006-04-04 Nissan Motor Co., Ltd. Ignition timing control for internal combustion engine
JP2005273576A (ja) * 2004-03-25 2005-10-06 Mitsubishi Fuso Truck & Bus Corp エンジンの燃料噴射制御方法及び燃料噴射制御装置
JP4774893B2 (ja) * 2005-09-28 2011-09-14 日産自動車株式会社 エンジンの点火時期制御方法及びエンジンの点火時期制御装置
JP4577211B2 (ja) 2005-12-27 2010-11-10 トヨタ自動車株式会社 Ｗｉｅｂｅ関数パラメータの決定方法および決定装置
WO2011011378A1 (en) * 2009-07-20 2011-01-27 Wayne State University Multi-sensing fuel injection system and method for making the same
JP2011220186A (ja) * 2010-04-08 2011-11-04 Toyota Motor Corp 内燃機関の燃焼制御装置
JP5062340B2 (ja) * 2011-03-11 2012-10-31 株式会社豊田自動織機 燃料噴射装置
JP2013163436A (ja) 2012-02-10 2013-08-22 Daimler Ag モータによるエンジン変動抑制装置
WO2013183163A1 (ja) * 2012-06-08 2013-12-12 トヨタ自動車株式会社 内燃機関の燃焼状態診断装置
DE102012018617B3 (de) * 2012-09-14 2014-03-27 Mtu Friedrichshafen Gmbh Verfahren zur Berechnung motorischer Kenngrößen, Datenverarbeitungssystem und Computerprogrammprodukt
JP6137220B2 (ja) * 2015-02-24 2017-05-31 トヨタ自動車株式会社 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
JP6135695B2 (ja) * 2015-02-26 2017-05-31 トヨタ自動車株式会社 燃焼状態推定方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2002322927A (ja) * 2001-04-26 2002-11-08 Nissan Motor Co Ltd 自己着火式内燃機関
JP2003278591A (ja) * 2002-03-20 2003-10-02 Denso Corp 内燃機関のトルク制御装置
JP2005048748A (ja) * 2003-07-31 2005-02-24 Nissan Motor Co Ltd 内燃機関の燃焼制御装置
JP2008274905A (ja) * 2007-05-07 2008-11-13 Toyota Motor Corp 多種燃料内燃機関の燃焼制御装置
JP2011080395A (ja) * 2009-10-06 2011-04-21 Cosmo Oil Co Ltd 圧縮自己着火内燃機関の燃焼時期予測方法、該圧縮自己着火内燃機関の制御方法及び圧縮自己着火内燃機関システム

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10001052B2 (en) 2015-02-26 2018-06-19 Toyota Jidosha Kabushiki Kaisha Combustion state estimation method
CN108368794A (zh) * 2015-12-02 2018-08-03 本田技研工业株式会社 内燃机的控制装置

Also Published As

Publication number Publication date
US20170037792A1 (en) 2017-02-09
EP3135887A4 (en) 2018-01-24
CN106232970A (zh) 2016-12-14
CN106232970B (zh) 2019-07-19
JP6260692B2 (ja) 2018-01-17
US9885295B2 (en) 2018-02-06
JPWO2015162971A1 (ja) 2017-04-13
EP3135887B1 (en) 2019-04-24
EP3135887A1 (en) 2017-03-01

Similar Documents

Publication Publication Date Title
JP6137220B2 (ja) 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
JP6183546B2 (ja) 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
US9903287B2 (en) Control apparatus for internal combustion engine
JP6260692B2 (ja) 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
JP6020690B2 (ja) 内燃機関の制御装置
JP6256599B2 (ja) 内燃機関の熱発生率波形算出装置および熱発生率波形算出方法
JP2007248119A (ja) Ｗｉｅｂｅ関数パラメータの決定方法および内燃機関の熱発生率推定装置
JP2008255932A (ja) 内燃機関の筒内圧推定方法および筒内圧推定装置
US10001052B2 (en) Combustion state estimation method
JP2008215204A (ja) 内燃機関の熱発生率のシミュレーション方法、内燃機関のトルクモデル作成方法および内燃機関のトルク推定方法
JP3873959B2 (ja) 内燃機関の燃料性状判定装置
JP4803099B2 (ja) 可変圧縮比エンジンのトルク推定装置
JP6332141B2 (ja) 内燃機関の異常診断装置
JP2015218589A (ja) 内燃機関の熱発生量算出装置、熱発生量算出方法、および燃焼モデル作成装置
JP6497202B2 (ja) 点火時期制御装置
송치헌 A Study on 0-D Knock Prediction Model in a Spark-ignited Engine
Yamazaki et al. Boundary Modeling and Identification of Normal Operation for Automobile Engine
JP2015218645A (ja) 内燃機関の開発支援装置
JP2015218704A (ja) 内燃機関の制御装置

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 15782377

Country of ref document: EP

Kind code of ref document: A1

WWE Wipo information: entry into national phase

Ref document number: 15303331

Country of ref document: US

ENP Entry into the national phase

Ref document number: 2016514747

Country of ref document: JP

Kind code of ref document: A

NENP Non-entry into the national phase

Ref country code: DE

REEP Request for entry into the european phase

Ref document number: 2015782377

Country of ref document: EP

WWE Wipo information: entry into national phase

Ref document number: 2015782377

Country of ref document: EP