▶ 独立行政法人海上技術安全研究所の特許一覧 ▶ 三井造船株式会社の特許一覧
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】特許公報(B2)
(11)【特許番号】
(24)【登録日】2023-07-06
(45)【発行日】2023-07-14
(54)【発明の名称】燃焼機関の実時間性能予測方法及び実時間性能予測プログラム
(51)【国際特許分類】
F02D 45/00 20060101AFI20230707BHJP
G06F 30/20 20200101ALI20230707BHJP
G06F 30/10 20200101ALI20230707BHJP
【ＦＩ】
F02D45/00 372
F02D45/00 370
G06F30/20
G06F30/10
【請求項の数】
16
(21)【出願番号】P 2019064402
(22)【出願日】2019-03-28
(65)【公開番号】P2020165341
(43)【公開日】2020-10-08
【審査請求日】2022-02-21
(73)【特許権者】
【識別番号】501204525
【氏名又は名称】国立研究開発法人 海上・港湾・航空技術研究所
(73)【特許権者】
【識別番号】000005902
【氏名又は名称】株式会社三井Ｅ＆Ｓ
(74)【代理人】
【識別番号】100098545
【弁理士】
【氏名又は名称】阿部 伸一
(74)【代理人】
【識別番号】100087745
【弁理士】
【氏名又は名称】清水 善廣
(74)【代理人】
【識別番号】100106611
【弁理士】
【氏名又は名称】辻田 幸史
(74)【代理人】
【識別番号】100189717
【弁理士】
【氏名又は名称】太田 貴章
(72)【発明者】
【氏名】ボンダレンコ オレクシー
(72)【発明者】
【氏名】福田 哲吾
(72)【発明者】
【氏名】北川 泰士
(72)【発明者】
【氏名】小川 靖之
(72)【発明者】
【氏名】福島 健史
(72)【発明者】
【氏名】宮川 修二郎
【審査官】北村 亮
(56)【参考文献】
【文献】特開２００４－２３９１２９（ＪＰ，Ａ）
【文献】特開２０１８－１９３９９４（ＪＰ，Ａ）
【文献】特開２０１１－１０６４５６（ＪＰ，Ａ）
【文献】特開２００８－１５１０５１（ＪＰ，Ａ）
【文献】特開２００８－００８２３６（ＪＰ，Ａ）
【文献】特開２００８－０９０４８８（ＪＰ，Ａ）
【文献】特開２０１８－１７８８７０（ＪＰ，Ａ）
(58)【調査した分野】(Int.Cl.，ＤＢ名)
Ｆ０２Ｄ ４５／００
Ｇ０６Ｆ ３０／２０
Ｇ０６Ｆ ３０／１０
(57)【特許請求の範囲】
【請求項1】
初期条件に基づき、燃焼機関の掃気系、排気系をモデル化した過給機を含んだサイクル平均値（ＣＭＶ）モデルによる掃排気状態の予測と、前記燃焼機関の複数の気筒の燃焼挙動をモデル化した現象モデルによる燃焼状態の予測とを統合して用い、前記燃焼機関の性能に関連したパラメータを実時間で予測して導出するにあたり、前記燃焼機関を通過する掃気と排気を連続した流れとして、前記気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するとともに、前記掃気量、前記排気量、及び前記燃焼挙動に伴う熱バランスから掃気圧力と排気圧力とを導出し、複数の前記気筒に対して圧力を同一に扱うことを特徴とする燃焼機関の実時間性能予測方法。
【請求項2】
前記サイクル平均値（ＣＭＶ）モデルは、前記掃気系の前記気筒の直前の掃気レシーバにおける前記掃気圧力と、前記排気系の前記気筒の直後の排気レシーバにおける前記排気圧力とを、前記掃気量及び前記排気量と前記燃焼挙動に伴う熱バランスから導出することを特徴とする請求項１に記載の燃焼機関の実時間性能予測方法。
【請求項3】
前記現象モデルは、前記気筒の排気弁の閉成から前記排気弁の開成までのクランク角度毎に前記気筒における燃料の燃焼に伴う物理量の変化を計算するものであることを特徴とする請求項１又は請求項２に記載の燃焼機関の実時間性能予測方法。
【請求項4】
前記クランク角度毎の前記燃焼に伴う前記物理量の変化として、前記燃焼に伴う温度上昇の変化を用い、前記温度上昇の前記計算に当たっては、前記燃焼の熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いることを特徴とする請求項３に記載の燃焼機関の実時間性能予測方法。
【請求項5】
前記クランク角度毎の前記温度上昇の前記計算は、前記クランク角度毎の前記温度上昇の差分を求める計算であることを特徴とする請求項４に記載の燃焼機関の実時間性能予測方法。
【請求項6】
前記温度上昇の差分を求める前記計算は、前記差分として未知の温度増加分に対する計算を、二次方程式として式（３８）及び式（３９）に基づいて行うことを特徴とする請求項５に記載の燃焼機関の実時間性能予測方法。【数38】
【数39】
【請求項7】
前記クランク角度毎の前記物理量の変化の前記計算は、前記気筒の前記排気弁の閉成から前記排気弁の開成までの間において導出することを特徴とする請求項３から請求項６のいずれか１項に記載の燃焼機関の実時間性能予測方法。
【請求項8】
前記燃焼機関の性能に関連した前記パラメータは、前記燃焼機関への前記掃気系、前記排気系の各部の圧力、温度、流量、前記気筒からの排気エネルギー、前記気筒の気筒内圧力、気筒内温度、及び前記燃焼機関のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることを特徴とする請求項１から請求項７のいずれか１項に記載の燃焼機関の実時間性能予測方法。
【請求項9】
導出された前記燃焼機関の性能に関連した前記パラメータに基づいて、前記燃焼機関の状態表示、前記燃焼機関の状態判断、前記燃焼機関の制御の少なくとも１つを行うことを特徴とする請求項１から請求項８のいずれか１項に記載の燃焼機関の実時間性能予測方法。
【請求項10】
コンピュータに、燃焼機関の掃気系、排気系をモデル化した過給機を含んだサイクル平均値（ＣＭＶ）モデル、及び前記燃焼機関の複数の気筒の燃焼挙動をモデル化した現象モデルを設定するモデル設定ステップと、
初期条件を入力する初期条件入力ステップと、
入力された前記初期条件に基づいて前記サイクル平均値（ＣＭＶ）モデルで掃排気状態の計算を行うサイクル平均値（ＣＭＶ）モデル計算ステップと、
前記初期条件及び前記サイクル平均値（ＣＭＶ）モデル計算ステップの計算結果に基づいて前記現象モデルで燃焼状態の計算を行う現象モデル計算ステップと、
前記サイクル平均値（ＣＭＶ）モデル計算ステップと前記現象モデル計算ステップの計算結果から導出される前記燃焼機関の性能に関連したパラメータを出力する出力ステップと
を実行させるにあたり、前記サイクル平均値（ＣＭＶ）モデル計算ステップで、前記燃焼機関を通過する掃気と排気を連続した流れとして、前記気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するとともに、前記掃気量、前記排気量、及び前記燃焼挙動に伴う熱バランスから掃気圧力と排気圧力とを導出して、複数の前記気筒に対して圧力を同一に扱うことを特徴とする燃焼機関の実時間性能予測プログラム。
【請求項11】
前記初期条件入力ステップで入力する前記初期条件として、前記燃焼機関の回転数、トルク、燃料投入量 、掃気圧、掃気温度、及び排ガス温度のうちの少なくとも１つを含むことを特徴とする請求項１０に記載の燃焼機関の実時間性能予測プログラム。
【請求項12】
前記サイクル平均値（ＣＭＶ）モデル計算ステップにおいて、サイクルごとに前記掃排気状態の計算を繰り返すことを特徴とする請求項１０又は請求項１１に記載の燃焼機関の実時間性能予測プログラム。
【請求項13】
前記現象モデル計算ステップにおいて、サイクルごとに前記燃焼状態の計算を繰り返すことを特徴とする請求項１０から請求項１２のいずれか１項に記載の燃焼機関の実時間性能予測プログラム。
【請求項14】
前記現象モデル計算ステップにおける前記燃焼状態の計算は、前記気筒の前記燃焼挙動を前記気筒の排気弁の閉成から前記排気弁の開成までの前記気筒のクランク角度毎に、燃焼に伴う熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いて行なうことを特徴とする請求項１０から請求項１３のいずれか１項に記載の燃焼機関の実時間性能予測プログラム。
【請求項15】
前記クランク角度毎の前記燃焼状態の計算は、前記クランク角度毎の温度上昇の差分を、二次方程式として式（３８）及び式（３９）に基づいて行うことを特徴とする請求項１４に記載の燃焼機関の実時間性能予測プログラム。【数38】
【数39】
【請求項16】
前記出力ステップにおける出力は、前記パラメータとして前記燃焼機関への前記掃気系、前記排気系の各部の圧力、温度、流量、前記気筒からの排気エネルギー、前記気筒の気筒内圧力、気筒内温度、及び前記燃焼機関のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることを特徴とする請求項１０から請求項１５のいずれか１項に記載の燃焼機関の実時間性能予測プログラム。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、舶用エンジン等、燃焼機関の性能を精度よく予測する燃焼機関の実時間性能予測方法及び実時間性能予測プログラムに関する。
【背景技術】
【0002】
貨物やエネルギー等の海上輸送の増加は、環境問題への深刻な懸念を引き起こし、これまで以上に効率的でクリーンとなる商用船が求められている。２０１３年に国際海事機関（ＩＭＯ）は、船舶からのＣＯ２排出量を制限するエネルギー効率設計指標（ＥＥＤＩ）を導入し、２０２５年以降に建造された船舶については最終的に３０％の削減を要求している。更に、ＭＥＰＣ７２（２０１８年４月）では、２００８年比で輸送１回あたりのＣＯ２排出量を２０３０年までに４０％削減し、２０５０年までの７０％削減に向けて努力することが合意された。また、２０５０年までにＧＨＧ（温室効果ガス）を少なくとも半分に削減し、今世紀内にできるだけ早く排出量をゼロにするというビジョンと目標を含むＩＭＯ ＧＨＧ戦略も決められた。これに対応して、推進システムの設計は、絶えず複雑さを増すという課題に直面している。ただし、エネルギー効率は設計上の問題だけでなく、運用運航時にも維持考慮する必要があるものである。
したがって、常にすべての運転点で推進システムの最適動作を追求する方向において、状態の推定及び実際の効率の追跡は非常に重要である。この点で、運用中に共存する推進システムの仮想モデル（デジタルツイン）を構築して、運用の予測と洞察を提供することは、実現可能な解決策の１つである。デジタルツインは推進システムの現在の状態をモデル化するために使用される物理学ベースのモデルの組織的な修正版である。デジタルツインの重要な要件は、すべての運転条件下で、ほぼリアルタイムに性能が物理的な対応物に確実に一致するようにし、推進プラント固有の特性を正確に反映するモデルであることである。
ディーゼルエンジンは、燃料の化学エネルギーを機械エネルギーに効率的に変換するため、船舶推進システムにとってなくてはならないものであり続けているため、デジタルツインの中核部分と見なされており、ディーゼルエンジンのモデリングは最も重要である。
ディーゼルエンジンのモデリングは、コンピュータシミュレーションの開発以来長年にわたって発展してきており、その複雑さの程度に応じて、伝達関数モデル、準定常平均値モデル、及び充填－空化現象論モデルなど、さまざまなモデルタイプを区別することができる。特定のモデルの選択は、前述のデジタルツインの要件（リアルタイムの実行時間とエンジンの作業プロセスについての洞察）によって決まる。推進システムの分野においては、サイクル平均値（ＣＭＶ）エンジンモデリングアプローチは、単純さと詳細さとの間の妥協において、エンジン定常性能及び過渡応答の評価に広く使用されている。サイクル平均値（ＣＭＶ）モデルで採用されている主な仮定は、エンジンはスロットルで直列に接続されたボリュームと見なされ、コンプレッサーによって押し出された空気と排気ガスはエンジンの動作周期に関係なく連続的に流れているとみなすことである。この点において、サイクル平均値（ＣＭＶ）エンジンモデルは、エンジンの極めて重要な要素である燃焼挙動の予測を欠いており、その結果として、エンジンの効率と同様に、正味燃料消費率（ｂｓｆｃ）に対する異なる設定の影響は、無視される。この欠点を克服するために多くの試みがなされている。例えば、経験的マップ又は人工ニューラルネットワーク（ＡＮＮ）を利用して燃焼サイクルを模擬するような、古典的なサイクル平均値（ＣＭＶ）モデルへの様々な拡張が導入されている（非特許文献１－３）。
一方、燃焼関連パラメータは、シリンダープロセスをいくつかの離散的な段階に分割し、温度、圧力と仕事の平均値を使用して代数関数のみを利用して計算する、Seiligerサイクルアプローチで推定することもできる（非特許文献４）。このアプローチは、非特許文献５で小型タンカーの完全な推進システムを構築するために使用され、非特許文献６では、排気ガス再循環システムと可変ノズル面積過給機タービンを備えたエンジンの制御指向モデルを開発するために使用された。非特許文献７では、サイクル平均値（ＣＭＶ）モデル化手法に対するさらなる改良が提案され、そこでは現象論的燃焼モデルがエンジンサイクルの閉部分（特にＥＶＣ（排気弁閉成タイミング：exhaust valve close timing, deg）からＥＶＯ（排気弁開成タイミング：exhaust valve open timing, deg）まで）を表すために使用され、サイクル平均値（ＣＭＶ）モデルは空気及び排気ガスのシミュレーションに使用され、他のエンジン部品の計算に使用された。サイクルモデルは各時間ステップで呼び出される必要があるクランク角分解保存則と現象論的燃焼モデルで構成されているが、微分形式で表され、それらの解はリアルタイム実行制約を満たさない。
したがって、常にすべての運転点で推進システムの最適動作を追求する方向において、状態の推定及び実際の効率の追跡は非常に重要である。この点で、運用中に共存する推進システムの仮想モデル（デジタルツイン）を構築して、運用の予測と洞察を提供することは、実現可能な解決策の１つである。デジタルツインは推進システムの現在の状態をモデル化するために使用される物理学ベースのモデルの組織的な修正版である。デジタルツインの重要な要件は、すべての運転条件下で、ほぼリアルタイムに性能が物理的な対応物に確実に一致するようにし、推進プラント固有の特性を正確に反映するモデルであることである。
ディーゼルエンジンは、燃料の化学エネルギーを機械エネルギーに効率的に変換するため、船舶推進システムにとってなくてはならないものであり続けているため、デジタルツインの中核部分と見なされており、ディーゼルエンジンのモデリングは最も重要である。
ディーゼルエンジンのモデリングは、コンピュータシミュレーションの開発以来長年にわたって発展してきており、その複雑さの程度に応じて、伝達関数モデル、準定常平均値モデル、及び充填－空化現象論モデルなど、さまざまなモデルタイプを区別することができる。特定のモデルの選択は、前述のデジタルツインの要件（リアルタイムの実行時間とエンジンの作業プロセスについての洞察）によって決まる。推進システムの分野においては、サイクル平均値（ＣＭＶ）エンジンモデリングアプローチは、単純さと詳細さとの間の妥協において、エンジン定常性能及び過渡応答の評価に広く使用されている。サイクル平均値（ＣＭＶ）モデルで採用されている主な仮定は、エンジンはスロットルで直列に接続されたボリュームと見なされ、コンプレッサーによって押し出された空気と排気ガスはエンジンの動作周期に関係なく連続的に流れているとみなすことである。この点において、サイクル平均値（ＣＭＶ）エンジンモデルは、エンジンの極めて重要な要素である燃焼挙動の予測を欠いており、その結果として、エンジンの効率と同様に、正味燃料消費率（ｂｓｆｃ）に対する異なる設定の影響は、無視される。この欠点を克服するために多くの試みがなされている。例えば、経験的マップ又は人工ニューラルネットワーク（ＡＮＮ）を利用して燃焼サイクルを模擬するような、古典的なサイクル平均値（ＣＭＶ）モデルへの様々な拡張が導入されている（非特許文献１－３）。
一方、燃焼関連パラメータは、シリンダープロセスをいくつかの離散的な段階に分割し、温度、圧力と仕事の平均値を使用して代数関数のみを利用して計算する、Seiligerサイクルアプローチで推定することもできる（非特許文献４）。このアプローチは、非特許文献５で小型タンカーの完全な推進システムを構築するために使用され、非特許文献６では、排気ガス再循環システムと可変ノズル面積過給機タービンを備えたエンジンの制御指向モデルを開発するために使用された。非特許文献７では、サイクル平均値（ＣＭＶ）モデル化手法に対するさらなる改良が提案され、そこでは現象論的燃焼モデルがエンジンサイクルの閉部分（特にＥＶＣ（排気弁閉成タイミング：exhaust valve close timing, deg）からＥＶＯ（排気弁開成タイミング：exhaust valve open timing, deg）まで）を表すために使用され、サイクル平均値（ＣＭＶ）モデルは空気及び排気ガスのシミュレーションに使用され、他のエンジン部品の計算に使用された。サイクルモデルは各時間ステップで呼び出される必要があるクランク角分解保存則と現象論的燃焼モデルで構成されているが、微分形式で表され、それらの解はリアルタイム実行制約を満たさない。
【0003】
ここで、特許文献１には、エンジンの状態量に関わる第１パラメータの値を検出するセンサの検出値を構成要素とする第１観測ベクトルを取得し、エンジンの状態量に関わる第２パラメータの仮想観測値を構成要素とする第２観測ベクトルを演算し、第１観測ベクトル及び第２観測ベクトルを用いたカルマンフィルター理論を適用してエンジンの状態量に関わる第３パラメータの予測推定値を構成要素とする予測推定ベクトルを演算し、その予測推定ベクトルに基づいて制御対象への制御指示値を演算するエンジンの制御装置が開示されている。
また、特許文献２には、実際のエンジンの吸気通路に、上流側から下流側へ順次、エアフローセンサ、スロットル弁、排気ターボ式過給機、インタークーラ、電動式過給機、吸気圧力センサが配設され、エンジン制御系として、吸気状態について影響を与える実際の機器類の特性と同一の特性を有するように設定された仮想機器類を複数種組み合わせて構成された同定モデルが設定されるエンジンの制御装置が開示されている。
また、特許文献３には、４サイクル多気筒ガスエンジンの各気筒のサイクルを吸気行程、圧縮行程、膨張行程、及び排気行程に区分して、サイクルタイム毎の１気筒当りの性能を計算するシリンダモデルを用い、４サイクル多気筒ガスエンジンの全気筒について、各気筒の行程順序に対応させて、シフト関数を用いてサイクルタイム分ずらして同時に各気筒分のエンジン性能を計算し、これを統合して４サイクル多気筒ガスエンジンのエンジン性能を計算する４サイクル多気筒ガスエンジンのシミュレーション方法が開示されている。
また、特許文献２には、実際のエンジンの吸気通路に、上流側から下流側へ順次、エアフローセンサ、スロットル弁、排気ターボ式過給機、インタークーラ、電動式過給機、吸気圧力センサが配設され、エンジン制御系として、吸気状態について影響を与える実際の機器類の特性と同一の特性を有するように設定された仮想機器類を複数種組み合わせて構成された同定モデルが設定されるエンジンの制御装置が開示されている。
また、特許文献３には、４サイクル多気筒ガスエンジンの各気筒のサイクルを吸気行程、圧縮行程、膨張行程、及び排気行程に区分して、サイクルタイム毎の１気筒当りの性能を計算するシリンダモデルを用い、４サイクル多気筒ガスエンジンの全気筒について、各気筒の行程順序に対応させて、シフト関数を用いてサイクルタイム分ずらして同時に各気筒分のエンジン性能を計算し、これを統合して４サイクル多気筒ガスエンジンのエンジン性能を計算する４サイクル多気筒ガスエンジンのシミュレーション方法が開示されている。
【先行技術文献】
【特許文献】
【0004】
【文献】特開２０１８－１７８８７０号公報
【文献】特開２００８－１５１０５１号公報
【文献】特開２００８－８２３６号公報
【非特許文献】
【0005】
【文献】Zhu J., Modeling and Simulation of Container Ship’s Main Diesel Engine. In: Proceedings of the International Multi Conference of Engineers and Computer Scientists. (IMECS 2008, Hong Kong); 2008, vol. 2, pp.1980-1983.
【文献】Livanos G., Papalambrou G., Kyrtatos NP., Electronic Engine Control for Ice Operation of Tankers. In: Proceedings of CIMAC Congress, Paper No.: 44, Viena, Austria; 2007
【文献】Nikzadfar K., Samekhi AH., An extended mean value model (MVEM) for control-oriented modelling of diesel engines transient performance and emissions. Fuel 2015; 154:275-292.
【文献】Ding Y. Characterising Combustion in Diesel Engines, using parametrised finite stage cylinder process models [dissertation]. Delft Technical University; 2011.
【文献】Sui C.-B., Ding Y., Propulsion System Simulation of Small Tanker Based on Cyclic Engine Model. In: Proceedings of the International Symposium on Marine Engineering. (ISME 2014, Harbin, China), Paper PID155, 2014.
【文献】Hansen J.M., Zander C.-G., Pedersen N., Blanke M., Vejlgaard-L.M., Modelling for Control of Exhaust Gas Recirculation on Large Diesel Engines. In: Proceedings of the 9th IFAC Conference on Control Applications in Marine Systems. (CAMS 2013, Osaka, Japan), 2013.
【文献】Baldi F, Theotokatos G, Andersson K. Development of a combined mean value-zero dimensional model and application for a large marine four-stroke diesel engine simulation. Appl Energy 2015;154:402-15
【発明の概要】
【発明が解決しようとする課題】
【0006】
特許文献１及び特許文献２は、エンジンの吸排気系の状態に基づいて制御しようとするものであり、気筒の燃焼挙動については考慮されていない。
特許文献３には、システム全体のサイクルタイムを決め、各シリンダ（気筒）の行程をシフト関数によりサイクルタイム分ずらして同時に各気筒分のエンジン性能を計算することが記載されているが、当該サイクルタイムはクランク軸の半回転に相当するものであり、精度が制御に使えるほど十分とはいえない。
また、エンジンの性能を精度よく予測するプログラムは実用されているが、いずれも計算時間が長く、実時間（オンライン）では使用できない。
そこで本発明は、精度が制御に使えるほど十分高く、しかも計算時間が短くリアルタイムで使用可能な燃焼機関の実時間性能予測方法及び実時間性能予測プログラムを提供することを目的とする。
特許文献３には、システム全体のサイクルタイムを決め、各シリンダ（気筒）の行程をシフト関数によりサイクルタイム分ずらして同時に各気筒分のエンジン性能を計算することが記載されているが、当該サイクルタイムはクランク軸の半回転に相当するものであり、精度が制御に使えるほど十分とはいえない。
また、エンジンの性能を精度よく予測するプログラムは実用されているが、いずれも計算時間が長く、実時間（オンライン）では使用できない。
そこで本発明は、精度が制御に使えるほど十分高く、しかも計算時間が短くリアルタイムで使用可能な燃焼機関の実時間性能予測方法及び実時間性能予測プログラムを提供することを目的とする。
【課題を解決するための手段】
【0007】
請求項１記載に対応した燃焼機関の実時間性能予測方法においては、初期条件に基づき、燃焼機関の掃気系、排気系をモデル化した過給機を含んだサイクル平均値（ＣＭＶ）モデルによる掃排気状態の予測と、燃焼機関の複数の気筒の燃焼挙動をモデル化した現象モデルによる燃焼状態の予測とを統合して用い、燃焼機関の性能に関連したパラメータを実時間で予測して導出するにあたり、燃焼機関を通過する掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するとともに、掃気量、排気量、及び燃焼挙動に伴う熱バランスから掃気圧力と排気圧力とを導出し、複数の気筒に対して圧力を同一に扱うことを特徴とする。
請求項１に記載の本発明によれば、掃排気状態の予測と燃焼挙動による予測とを合理的に統合して、制御に使用できるほどの十分高い精度と計算時間をもって燃焼機関の性能予測を実時間で行うことができる。また、掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
請求項１に記載の本発明によれば、掃排気状態の予測と燃焼挙動による予測とを合理的に統合して、制御に使用できるほどの十分高い精度と計算時間をもって燃焼機関の性能予測を実時間で行うことができる。また、掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
【0008】
請求項２記載の本発明は、サイクル平均値（ＣＭＶ）モデルは、掃気系の気筒の直前の掃気レシーバにおける掃気圧力と、排気系の気筒の直後の排気レシーバにおける排気圧力とを、掃気量及び排気量と燃焼挙動に伴う熱バランスから導出することを特徴とする。
請求項２に記載の本発明によれば、掃気レシーバと排気レシーバの圧力を複数の気筒に対して同一の圧力をもって扱うことができる。
請求項２に記載の本発明によれば、掃気レシーバと排気レシーバの圧力を複数の気筒に対して同一の圧力をもって扱うことができる。
【0009】
請求項３記載の本発明は、現象モデルは、気筒の排気弁の閉成から排気弁の開成までのクランク角度毎に気筒における燃料の燃焼に伴う物理量の変化を計算するものであることを特徴とする。
請求項３に記載の本発明によれば、排気弁の閉成から開成までという区間に限定することで、実際に気筒内での燃焼に関連した計算に限定ができ、気筒の初期状態が明確となり、微分方程式ではなく差分法を用いて燃焼に伴う物理量の変化を計算することができるため、迅速に解を求めクランク角度毎に現象モデルによる燃焼状態の計算を行うことができる。
請求項３に記載の本発明によれば、排気弁の閉成から開成までという区間に限定することで、実際に気筒内での燃焼に関連した計算に限定ができ、気筒の初期状態が明確となり、微分方程式ではなく差分法を用いて燃焼に伴う物理量の変化を計算することができるため、迅速に解を求めクランク角度毎に現象モデルによる燃焼状態の計算を行うことができる。
【0010】
請求項４記載の本発明は、クランク角度毎の燃焼に伴う物理量の変化として、燃焼に伴う温度上昇の変化を用い、温度上昇の計算に当たっては、燃焼の熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いることを特徴とする。
請求項４に記載の本発明によれば、現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して迅速に行うことができる。
請求項４に記載の本発明によれば、現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して迅速に行うことができる。
【0011】
請求項５記載の本発明は、クランク角度毎の温度上昇の計算は、クランク角度毎の温度上昇の差分を求める計算であることを特徴とする。
請求項５に記載の本発明によれば、微分方程式を解く場合等と比較して現象モデルによる燃焼状態の計算を迅速に行うことができる。
請求項５に記載の本発明によれば、微分方程式を解く場合等と比較して現象モデルによる燃焼状態の計算を迅速に行うことができる。
【0012】
請求項６記載の本発明は、温度上昇の差分を求める計算は、差分として未知の温度増加分に対する計算を、二次方程式として式（３８）及び式（３９）に基づいて行うことを特徴とする。
請求項６に記載の本発明によれば、微分方程式を解いたり近似解を求めるのではなく差分法で二次方程式化することで、現象モデルによる燃焼状態の予測を迅速に行うことができる。
請求項６に記載の本発明によれば、微分方程式を解いたり近似解を求めるのではなく差分法で二次方程式化することで、現象モデルによる燃焼状態の予測を迅速に行うことができる。
【0013】
請求項７記載の本発明は、クランク角度毎の物理量の変化の計算は、気筒の排気弁の閉成から排気弁の開成までの間において導出することを特徴とする。
請求項７に記載の本発明によれば、現象モデルによる燃焼状態の計算を排気弁の閉成から開成までの間で迅速に完了することができる。
請求項７に記載の本発明によれば、現象モデルによる燃焼状態の計算を排気弁の閉成から開成までの間で迅速に完了することができる。
【0014】
請求項８記載の本発明は、燃焼機関の性能に関連したパラメータは、燃焼機関への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒の気筒内圧力、気筒内温度、及び燃焼機関のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることを特徴とする。
請求項８に記載の本発明によれば、燃焼機関の性能評価に必要なこれらのパラメータに関連した性能予測やモデルの改善を行うことができる。
請求項８に記載の本発明によれば、燃焼機関の性能評価に必要なこれらのパラメータに関連した性能予測やモデルの改善を行うことができる。
【0015】
請求項９記載の本発明は、導出された燃焼機関の性能に関連したパラメータに基づいて、燃焼機関の状態表示、燃焼機関の状態判断、燃焼機関の制御の少なくとも１つを行うことを特徴とする。
請求項９に記載の本発明によれば、得られたパラメータを有効に活用することができる。
請求項９に記載の本発明によれば、得られたパラメータを有効に活用することができる。
【0016】
請求項１０記載に対応した燃焼機関の実時間性能予測プログラムにおいては、コンピュータに、燃焼機関の掃気系、排気系をモデル化した過給機を含んだサイクル平均値（ＣＭＶ）モデル及び燃焼機関の複数の気筒の燃焼挙動をモデル化した現象モデルを設定するモデル設定ステップと、初期条件を入力する初期条件入力ステップと、入力された初期条件に基づいてサイクル平均値（ＣＭＶ）モデルで掃排気状態の計算を行うサイクル平均値（ＣＭＶ）モデル計算ステップと、初期条件及びサイクル平均値（ＣＭＶ）モデル計算ステップの計算結果に基づいて現象モデルで燃焼状態の計算を行う現象モデル計算ステップと、サイクル平均値（ＣＭＶ）モデル計算ステップと現象モデル計算ステップの計算結果から導出される燃焼機関の性能に関連したパラメータを出力する出力ステップとを実行させるにあたり、サイクル平均値（ＣＭＶ）モデル計算ステップで、燃焼機関を通過する掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するとともに、掃気量、排気量、及び燃焼挙動に伴う熱バランスから掃気圧力と排気圧力とを導出して、複数の気筒に対して圧力を同一に扱うことを特徴とする。
請求項１０に記載の本発明によれば、掃排気状態の予測と燃焼挙動による予測とを合理的に統合して、制御に使用できるほどの十分高い精度と計算時間をもって性能予測を実時間で行うことができる。また、掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
請求項１０に記載の本発明によれば、掃排気状態の予測と燃焼挙動による予測とを合理的に統合して、制御に使用できるほどの十分高い精度と計算時間をもって性能予測を実時間で行うことができる。また、掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
【0017】
請求項１１記載の本発明は、初期条件入力ステップで入力する初期条件として、燃焼機関の回転数、トルク、燃料投入量 、掃気圧、掃気温度、及び排ガス温度のうちの少なくとも１つを含むことを特徴とする。
請求項１１に記載の本発明によれば、掃排気状態と燃焼状態の予測に必要な燃焼機関の回転数、トルク、燃料投入量、掃気圧、排気圧、排気温度、又は排ガス温度を考慮した性能予測を行うことができる。
請求項１１に記載の本発明によれば、掃排気状態と燃焼状態の予測に必要な燃焼機関の回転数、トルク、燃料投入量、掃気圧、排気圧、排気温度、又は排ガス温度を考慮した性能予測を行うことができる。
【0018】
請求項１２記載の本発明は、サイクル平均値（ＣＭＶ）モデル計算ステップにおいて、サイクルごとに掃排気状態の計算を繰り返すことを特徴とする。
請求項１２に記載の本発明によれば、実時間性能予測の精度を向上させることができる。
請求項１２に記載の本発明によれば、実時間性能予測の精度を向上させることができる。
【0019】
請求項１３記載の本発明は、現象モデル計算ステップにおいて、サイクルごとに燃焼状態の計算を繰り返すことを特徴とする。
請求項１３に記載の本発明によれば、実時間性能予測の精度を向上させることができる。
請求項１３に記載の本発明によれば、実時間性能予測の精度を向上させることができる。
【0020】
請求項１４記載の本発明は、現象モデル計算ステップにおける燃焼状態の計算は、気筒の燃焼挙動を気筒の排気弁の閉成から排気弁の開成までの気筒のクランク角度毎に、燃焼に伴う熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いて行なうことを特徴とする。
請求項１４に記載の本発明によれば、クランク角度毎に現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して行うことができる。
請求項１４に記載の本発明によれば、クランク角度毎に現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して行うことができる。
【0021】
請求項１５記載の本発明は、クランク角度毎の燃焼状態の計算は、クランク角度毎の温度上昇の差分を、二次方程式として式（３８）及び式（３９）に基づいて行うことを特徴とする。
請求項１５に記載の本発明によれば、微分方程式を解いたり近似解を求めるのではなく差分法で二次方程式化することで、現象モデル計算ステップにおける燃焼状態の計算を迅速に行なうことができる。
請求項１５に記載の本発明によれば、微分方程式を解いたり近似解を求めるのではなく差分法で二次方程式化することで、現象モデル計算ステップにおける燃焼状態の計算を迅速に行なうことができる。
【0022】
請求項１６記載の本発明は、出力ステップにおける出力は、パラメータとして燃焼機関への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒の気筒内圧力、気筒内温度、及び燃焼機関のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることを特徴とする。
請求項１６に記載の本発明によれば、これらのパラメータに関連した性能予測やモデルの改善を行うことができる。
請求項１６に記載の本発明によれば、これらのパラメータに関連した性能予測やモデルの改善を行うことができる。
【発明の効果】
【0023】
本発明の燃焼機関の実時間性能予測方法によれば、掃排気状態の予測と燃焼挙動による予測とを合理的に統合して、制御に使用できるほどの十分高い精度と計算時間をもって燃焼機関の性能予測を実時間で行うことができる。また、掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
【0024】
また、サイクル平均値（ＣＭＶ）モデルは、掃気系の気筒の直前の掃気レシーバにおける掃気圧力と、排気系の気筒の直後の排気レシーバにおける排気圧力とを、掃気量及び排気量と燃焼挙動に伴う熱バランスから導出する場合には、掃気レシーバと排気レシーバの圧力を複数の気筒に対して同一の圧力をもって扱うことができる。
【0025】
また、現象モデルは、気筒の排気弁の閉成から排気弁の開成までのクランク角度毎に気筒における燃料の燃焼に伴う物理量の変化を計算するものである場合には、排気弁の閉成から開成までという区間に限定することで、実際に気筒内での燃焼に関連した計算に限定ができ、気筒の初期状態が明確となり、微分方程式ではなく差分法を用いて燃焼に伴う物理量の変化を計算することができるため、迅速に解を求めクランク角度毎に現象モデルによる燃焼状態の計算を行うことができる。
【0026】
また、クランク角度毎の燃焼に伴う物理量の変化として、燃焼に伴う温度上昇の変化を用い、温度上昇の計算に当たっては、燃焼の熱発生パターンを近似するウィーベ (Ｗｉｅｂｅ) 関数を用いる場合には、現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して迅速に行うことができる。
【0027】
また、クランク角度毎の温度上昇の計算は、クランク角度毎の温度上昇の差分を求める計算である場合には、微分方程式を解く場合等と比較して現象モデルによる燃焼状態の計算を迅速に行うことができる。
【0028】
また、温度上昇の差分を求める計算は、差分として未知の温度増加分に対する計算を、二次方程式として式（３８）及び式（３９）に基づいて行う場合には、微分方程式を解いたり近似解を求めるのではなく差分法で二次方程式化することで、現象モデルによる燃焼状態の予測を迅速に行うことができる。
【0029】
また、クランク角度毎の物理量の変化の計算は、気筒の排気弁の閉成から排気弁の開成までの間において導出する場合には、現象モデルによる燃焼状態の計算を排気弁の閉成から開成までの間で迅速に完了することができる。
【0030】
また、燃焼機関の性能に関連したパラメータは、燃焼機関への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒の気筒内圧力、気筒内温度、及び燃焼機関のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものである場合には、燃焼機関の性能評価に必要なこれらのパラメータに関連した性能予測やモデルの改善を行うことができる。
【0031】
また、導出された燃焼機関の性能に関連したパラメータに基づいて、燃焼機関の状態表示、燃焼機関の状態判断、燃焼機関の制御の少なくとも１つを行う場合には、得られたパラメータを有効に活用することができる。
【0032】
また、本発明の燃焼機関の実時間性能予測プログラムによれば、掃排気状態の予測と燃焼挙動による予測とを合理的に統合して、制御に使用できるほどの十分高い精度と計算時間をもって性能予測を実時間で行うことができる。また、掃気と排気を連続した流れとして、気筒の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
【0033】
また、初期条件入力ステップで入力する初期条件として、燃焼機関の回転数、トルク、燃料投入量 、掃気圧、掃気温度、及び排ガス温度のうちの少なくとも１つを含む場合には、掃排気状態と燃焼状態の予測に必要な燃焼機関の回転数、トルク、燃料投入量、掃気圧、排気圧、排気温度、又は排ガス温度を考慮した性能予測を行うことができる。
【0034】
また、サイクル平均値（ＣＭＶ）モデル計算ステップにおいて、サイクルごとに掃排気状態の計算を繰り返す場合には、実時間性能予測の精度を向上させることができる。
【0035】
また、現象モデル計算ステップにおいて、サイクルごとに燃焼状態の計算を繰り返す場合には、実時間性能予測の精度を向上させることができる。
【0036】
また、現象モデル計算ステップにおける燃焼状態の計算は、気筒の燃焼挙動を気筒の排気弁の閉成から排気弁の開成までの気筒のクランク角度毎に燃焼に伴う熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いて行なう場合には、クランク角度毎に現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して行うことができる。
【0037】
また、クランク角度毎の燃焼状態の計算は、クランク角度毎の温度上昇の差分を、二次方程式として式（３８）及び式（３９）に基づいて行う場合には、微分方程式を解いたり近似解を求めるのではなく差分法で二次方程式化することで、現象モデル計算ステップにおける燃焼状態の計算を迅速に行なうことができる。
【0038】
また、出力ステップにおける出力は、パラメータとして燃焼機関への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒の気筒内圧力、気筒内温度、及び燃焼機関のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものである場合には、これらのパラメータに関連した性能予測やモデルの改善を行うことができる。
【図面の簡単な説明】
【0039】
【図1】本実施形態による燃焼機関の実時間性能予測を適用したＨｉＬテストベンチを示す図
【図2】同燃焼機関の実時間性能予測を適用する船舶の一般的な推進プラントモデルを示す図
【図3】同燃焼機関の実時間性能予測の概念図
【図4】同燃焼機関の実時間性能予測の説明図
【図6】同燃焼機関の実時間性能予測のフローチャート
【図7】本実施例による定常状態シミュレーション結果を示す図
【図8】同モデルの枠組みにおける燃焼性能を示す図
【図9】同計算手順を評価した試験の結果を示す図
【発明を実施するための形態】
【0040】
以下に、本発明の実施形態による燃焼機関の実時間性能予測方法及び実時間性能予測プログラムについて説明する。
【0041】
図１は本実施形態による燃焼機関の実時間性能予測を適用したハードウェアインザループ（ＨｉＬ：Hardware in the Loop）テストベンチを示す図である。
仮想プラント（サイクル平均値（ＣＭＶ）モデル、現象モデル）２は実プラント（燃焼機関）１をモデル化したものであり、仮想プラント２を用いて実時間性能予測を行う。また、予測した結果の信号を変換するアナログ信号インターフェース３、燃焼機関の制御等を行うマイクロコントローラ４．及び燃焼機関の状態予測結果や状態判断結果を表示するモニタを備えたコンピュータ５の間ではデータ交換が行われる。
仮想プラント（サイクル平均値（ＣＭＶ）モデル、現象モデル）２は実プラント（燃焼機関）１をモデル化したものであり、仮想プラント２を用いて実時間性能予測を行う。また、予測した結果の信号を変換するアナログ信号インターフェース３、燃焼機関の制御等を行うマイクロコントローラ４．及び燃焼機関の状態予測結果や状態判断結果を表示するモニタを備えたコンピュータ５の間ではデータ交換が行われる。
【0042】
図２は本実施形態による燃焼機関の実時間性能予測を適用する船舶の一般的な推進プラントモデルを示す図である。
船舶１０に設置されている指令器１１からの指令に基づきガバナー１２が燃料ラックを制御して過給機１３を備えた燃焼機関（エンジン）１４へ燃料を供給する。燃焼機関１４にはシャフト１５を介してプロペラ１６が接続されている。
船体ダイナミクスは下式（１）で表される。
また、燃焼機関ダイナミクスは下式（２）で表される。
また、ガバナーは下式（３）で表される。
ここで、ｖ：船速、ne：エンジン回転速度、ｈp：燃料ポンプインデックス、Ｔp：プロペラスラスト、Ｑp：プロペラトルク、Ｒtot：船体抵抗、Ｆw：波により加わる抵抗（外乱）、Ｑe：エンジントルク、Ｋp,Ｋi,τsm：ガバナーのパラメータである。
エンジントルクＱeはブレーキ平均有効圧力Ｐb（Ｑe∝Ｐb）に比例し、これはエンジンシリンダ内の燃料燃焼の結果である。
船舶１０に設置されている指令器１１からの指令に基づきガバナー１２が燃料ラックを制御して過給機１３を備えた燃焼機関（エンジン）１４へ燃料を供給する。燃焼機関１４にはシャフト１５を介してプロペラ１６が接続されている。
船体ダイナミクスは下式（１）で表される。
【数1】
【数2】
【数3】
エンジントルクＱeはブレーキ平均有効圧力Ｐb（Ｑe∝Ｐb）に比例し、これはエンジンシリンダ内の燃料燃焼の結果である。
【0043】
図３は本実施形態による燃焼機関の実時間性能予測の概念図である。
従来、詳細に応じて、発生トルクに関してエンジン性能をシミュレートする２つの主な方法がある。一つは排気系をモデル化したサイクル平均値（ＣＭＶ：Cycle-Mean Value Engine）モデルを使用するものであり、もう一つは燃焼機関１４の気筒の燃焼挙動モデルをモデル化した現象モデルを使用するものである。現象モデル（Phenomenological Model）としては、充填－空化モデル（Filling-Emptying Model）が挙げられる。
サイクル平均値（ＣＭＶ）モデルによる予測は、時間領域における連続微分方程式の結果を高速で得られるが、燃焼機関１４の性能に関する限られた情報である。
また、現象モデルによる予測は、燃焼機関１４の性能に関する豊富な情報が得られるが、計算が遅いためクランクシャフト角度領域における連続微分方程式の計算結果を迅速に得られない。
そこで、本実施形態では、燃焼機関１４の掃気系、排気系をモデル化したサイクル平均値（ＣＭＶ）モデルによる掃排気状態の予測と、燃焼機関１４の気筒の燃焼挙動をモデル化した現象モデルによる燃焼状態の予測とを統合して用いる。
従来、詳細に応じて、発生トルクに関してエンジン性能をシミュレートする２つの主な方法がある。一つは排気系をモデル化したサイクル平均値（ＣＭＶ：Cycle-Mean Value Engine）モデルを使用するものであり、もう一つは燃焼機関１４の気筒の燃焼挙動モデルをモデル化した現象モデルを使用するものである。現象モデル（Phenomenological Model）としては、充填－空化モデル（Filling-Emptying Model）が挙げられる。
サイクル平均値（ＣＭＶ）モデルによる予測は、時間領域における連続微分方程式の結果を高速で得られるが、燃焼機関１４の性能に関する限られた情報である。
また、現象モデルによる予測は、燃焼機関１４の性能に関する豊富な情報が得られるが、計算が遅いためクランクシャフト角度領域における連続微分方程式の計算結果を迅速に得られない。
そこで、本実施形態では、燃焼機関１４の掃気系、排気系をモデル化したサイクル平均値（ＣＭＶ）モデルによる掃排気状態の予測と、燃焼機関１４の気筒の燃焼挙動をモデル化した現象モデルによる燃焼状態の予測とを統合して用いる。
【0044】
図４は本実施形態による燃焼機関の実時間性能予測の説明図であり、図５は図４を簡易化した計算フロー図である。
図４には、空気圧縮機１７、空気冷却器１８、掃気レシーバ１９、気筒（シリンダ）２０、排気レシーバ２１、燃料ポンプ２２、タービン２３、コンピュータ５を示している。また、図５では、空気圧縮機１７、掃気レシーバ１９、気筒２０、排気レシーバ２１、タービン２３を示している。
図４には、空気圧縮機１７、空気冷却器１８、掃気レシーバ１９、気筒（シリンダ）２０、排気レシーバ２１、燃料ポンプ２２、タービン２３、コンピュータ５を示している。また、図５では、空気圧縮機１７、掃気レシーバ１９、気筒２０、排気レシーバ２１、タービン２３を示している。
【0045】
図６は、本発明の実施形態による燃焼機関の実時間性能予測のフローチャートである。
まず、初期条件として燃焼機関１４等に関する定数やパラメータを入力する（Ｓ１：初期条件入力ステップ）。初期条件としては、燃焼機関１４の回転数、トルク、燃料投入量、掃気圧、排気圧、排気温度、又は排ガス温度のうちの少なくとも１つを含むことが好ましい。これにより、掃排気状態と燃焼状態の予測に必要な燃焼機関１４の回転数、トルク、燃料投入量、掃気圧、排気圧、排気温度、又は排ガス温度を考慮した性能予測を行うことができる。なお、初期条件として入力するパラメータには、実船モニタリングにより取得した実データを用いることもできる。
また、使用するサイクル平均値（ＣＭＶ）モデルと、現象モデルを設定する（Ｓ２：モデル設定ステップ）。
まず、初期条件として燃焼機関１４等に関する定数やパラメータを入力する（Ｓ１：初期条件入力ステップ）。初期条件としては、燃焼機関１４の回転数、トルク、燃料投入量、掃気圧、排気圧、排気温度、又は排ガス温度のうちの少なくとも１つを含むことが好ましい。これにより、掃排気状態と燃焼状態の予測に必要な燃焼機関１４の回転数、トルク、燃料投入量、掃気圧、排気圧、排気温度、又は排ガス温度を考慮した性能予測を行うことができる。なお、初期条件として入力するパラメータには、実船モニタリングにより取得した実データを用いることもできる。
また、使用するサイクル平均値（ＣＭＶ）モデルと、現象モデルを設定する（Ｓ２：モデル設定ステップ）。
【0046】
次に、入力された初期条件に基づいてサイクル平均値（ＣＭＶ）モデルで掃排気状態の計算を行う（Ｓ３：サイクル平均値（ＣＭＶ）モデル計算ステップ）。このとき、サイクル平均値（ＣＭＶ）モデルは過給機１３を含み、掃排気状態の計算は、燃焼機関１４を通過する掃気と排気を連続した流れとして、気筒２０の１サイクル当たりの平均値に基づいた掃気量及び排気量の計算であることが好ましい。これにより、サイクル平均値（ＣＭＶ）モデルによる掃排気状態の予測と、現象モデルによる燃焼状態の予測とを統合して用いるに当り、掃気と排気を連続した流れとして、気筒２０の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒２０に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
【0047】
次に、初期条件と、サイクル平均値（ＣＭＶ）モデル計算ステップＳ３の計算結果に基づいて現象モデルで燃焼状態の計算を行う（Ｓ４：現象モデル計算ステップ）。現象モデル計算ステップＳ４における燃焼状態の計算は、気筒２０の燃焼挙動を気筒２０の排気弁の閉成から排気弁の開成までの気筒２０のクランク角度毎に、燃焼に伴う熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いて行なう。これにより、クランク角度毎に現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して行うことができる。
現象モデル計算ステップＳ４の後、排気弁が開成か否かを判断する（Ｓ５：排気弁開成確認ステップ）。
排気弁開成確認ステップＳ５において排気弁が開成でないと判断した場合は、ｉ＝ｉ＋１として現象モデル計算ステップＳ４に戻る。現象モデル計算ステップＳ４において、サイクルごとに燃焼状態の計算を繰り返すことで、実時間性能予測の精度を向上させることができる。
排気弁開成確認ステップＳ５において排気弁が開成であると判断した場合は、現象モデル計算ステップＳ４における計算結果を出力する（Ｓ６：現象モデル計算結果出力ステップ）。
現象モデル計算ステップＳ４の後、排気弁が開成か否かを判断する（Ｓ５：排気弁開成確認ステップ）。
排気弁開成確認ステップＳ５において排気弁が開成でないと判断した場合は、ｉ＝ｉ＋１として現象モデル計算ステップＳ４に戻る。現象モデル計算ステップＳ４において、サイクルごとに燃焼状態の計算を繰り返すことで、実時間性能予測の精度を向上させることができる。
排気弁開成確認ステップＳ５において排気弁が開成であると判断した場合は、現象モデル計算ステップＳ４における計算結果を出力する（Ｓ６：現象モデル計算結果出力ステップ）。
【0048】
次に、サイクル平均値（ＣＭＶ）モデル計算ステップＳ３における計算結果と、現象モデル計算結果出力ステップＳ６で出力された計算結果とから導出される燃焼機関１４の性能に関連したパラメータを更新して出力する（Ｓ７：更新ステップ）。
更新、出力するパラメータとしては、燃焼機関１４への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒２０の気筒内圧力、気筒内温度、及び燃焼機関１４のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることが好ましい。これにより、これらのパラメータに関連した性能予測やモデルの改善を行うことができる。なお、図６では更新、出力するパラメータとして、コンプレッサーパワーＷc、空気圧Ｐs、排気温度Ｔexh、排気レシーバ２１内の排ガス質量Ｍexh、燃焼機関１４の回転速度ｎeを例示している。
更新、出力するパラメータとしては、燃焼機関１４への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒２０の気筒内圧力、気筒内温度、及び燃焼機関１４のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることが好ましい。これにより、これらのパラメータに関連した性能予測やモデルの改善を行うことができる。なお、図６では更新、出力するパラメータとして、コンプレッサーパワーＷc、空気圧Ｐs、排気温度Ｔexh、排気レシーバ２１内の排ガス質量Ｍexh、燃焼機関１４の回転速度ｎeを例示している。
【0049】
更新ステップＳ７の後、終了時刻Ｔendに達したか否かを判断する（Ｓ８：終了判断ステップ）。
終了判断ステップＳ８において終了時刻Ｔendに達していないと判断した場合は、ｊ＝ｊ＋１としてサイクル平均値（ＣＭＶ）モデル計算ステップＳ３に戻る。サイクル平均値（ＣＭＶ）モデル計算ステップＳ３において、サイクルごとに掃排気状態の計算を繰り返すことで、モデルが実際の燃焼機関１４に益々近似し、実時間性能予測の精度を向上させることができる。
終了判断ステップＳ８において終了時刻Ｔendに達したと判断した場合は、サイクル平均値（ＣＭＶ）モデル計算ステップＳ３と現象モデル計算ステップＳ４の計算結果から導出される燃焼機関１４の性能に関連したパラメータを出力し（Ｓ９：出力ステップ）、実時間性能予測を終了する。
終了判断ステップＳ８において終了時刻Ｔendに達していないと判断した場合は、ｊ＝ｊ＋１としてサイクル平均値（ＣＭＶ）モデル計算ステップＳ３に戻る。サイクル平均値（ＣＭＶ）モデル計算ステップＳ３において、サイクルごとに掃排気状態の計算を繰り返すことで、モデルが実際の燃焼機関１４に益々近似し、実時間性能予測の精度を向上させることができる。
終了判断ステップＳ８において終了時刻Ｔendに達したと判断した場合は、サイクル平均値（ＣＭＶ）モデル計算ステップＳ３と現象モデル計算ステップＳ４の計算結果から導出される燃焼機関１４の性能に関連したパラメータを出力し（Ｓ９：出力ステップ）、実時間性能予測を終了する。
【0050】
このように、初期条件に基づき、燃焼機関１４の掃気系、排気系をモデル化したサイクル平均値（ＣＭＶ）モデルによる掃排気状態の予測と、燃焼機関１４の気筒２０の燃焼挙動をモデル化した現象モデルによる燃焼状態の予測とを統合して用い、燃焼機関１４の性能に関連したパラメータを実時間で予測して導出することで、掃排気状態の予測と燃焼挙動による予測とを合理的に統合して、制御に使用できるほどの十分高い精度と計算時間をもって性能予測を実時間で行うことができる。
また、サイクル平均値（ＣＭＶ）モデルは過給機を含み、燃焼機関１４を通過する掃気と排気を連続した流れとして、気筒２０の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出することで、掃気と排気を連続した流れとして、気筒２０の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒２０に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
また、サイクル平均値（ＣＭＶ）モデルは、掃気系の気筒２０の直前の掃気レシーバ１９における掃気圧力と、排気系の気筒２０の直後の排気レシーバ２１における排気圧力とを、掃気量及び排気量と燃焼挙動に伴う熱バランスから導出することで、掃気レシーバ１９と排気レシーバ２１の圧力を複数の気筒２０に対して同一の圧力をもって扱うことができる。また、掃気系と排気系は気体の系であり、爆発等の急激な圧力変化がないところ、サイクル平均値（ＣＭＶ）を用いても複数の気筒２０に対して圧力を同一に扱うことができる。このことは、現象モデルの条件として排気弁の閉成時と開成時に各気筒２０に作用する掃気圧と排気圧を同一に扱うことができ、計算が単純化できる。
また、現象モデルは、気筒２０の排気弁の閉成から排気弁の開成までのクランク角度毎に気筒２０における燃料の燃焼に伴う物理量の変化を計算することで、排気弁の閉成から開成までという区間に限定して実際に気筒２０内での燃焼に関連した計算に限定ができ、気筒２０の初期状態が明確となり、微分方程式ではなく差分法を用いて燃焼に伴う物理量の変化を計算することができるため、迅速に解を求めクランク角度毎に現象モデルによる燃焼状態の計算を行うことができる。
また、クランク角度毎の燃焼に伴う物理量の変化として、燃焼に伴う温度上昇の変化を用い、温度上昇の計算に当たっては、燃焼の熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いることで、現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して迅速に行うことができる。
また、クランク角度毎の温度上昇の計算は、クランク角度毎の温度上昇の差分を求める計算とすることで、微分方程式を解く場合等と比較して現象モデルによる燃焼状態の計算を迅速に行うことができる。
また、クランク角度毎の物理量の変化の計算は、気筒２０の排気弁の閉成から排気弁の開成までの間において導出することで、現象モデルによる燃焼状態の計算を排気弁の閉成から開成までの間で迅速に完了することができる。
また、燃焼機関１４の性能に関連したパラメータは、燃焼機関１４への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒２０の気筒内圧力、気筒内温度、及び燃焼機関１４のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることで、燃焼機関の性能評価に必要なこれらのパラメータに関連した性能予測やモデルの改善を行うことができる。
また、導出された燃焼機関の性能に関連したパラメータに基づいて、燃焼機関１４の状態表示、燃焼機関１４の状態判断、燃焼機関１４の制御の少なくとも１つを行うことで、得られたパラメータを有効に活用することができる。
また、サイクル平均値（ＣＭＶ）モデルは過給機を含み、燃焼機関１４を通過する掃気と排気を連続した流れとして、気筒２０の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出することで、掃気と排気を連続した流れとして、気筒２０の１サイクル当たりの平均値に基づいた掃気量及び排気量を導出するため、複数ある気筒２０に対してサイクル平均値（ＣＭＶ）モデルを用いても誤差を少なくすることができる。
また、サイクル平均値（ＣＭＶ）モデルは、掃気系の気筒２０の直前の掃気レシーバ１９における掃気圧力と、排気系の気筒２０の直後の排気レシーバ２１における排気圧力とを、掃気量及び排気量と燃焼挙動に伴う熱バランスから導出することで、掃気レシーバ１９と排気レシーバ２１の圧力を複数の気筒２０に対して同一の圧力をもって扱うことができる。また、掃気系と排気系は気体の系であり、爆発等の急激な圧力変化がないところ、サイクル平均値（ＣＭＶ）を用いても複数の気筒２０に対して圧力を同一に扱うことができる。このことは、現象モデルの条件として排気弁の閉成時と開成時に各気筒２０に作用する掃気圧と排気圧を同一に扱うことができ、計算が単純化できる。
また、現象モデルは、気筒２０の排気弁の閉成から排気弁の開成までのクランク角度毎に気筒２０における燃料の燃焼に伴う物理量の変化を計算することで、排気弁の閉成から開成までという区間に限定して実際に気筒２０内での燃焼に関連した計算に限定ができ、気筒２０の初期状態が明確となり、微分方程式ではなく差分法を用いて燃焼に伴う物理量の変化を計算することができるため、迅速に解を求めクランク角度毎に現象モデルによる燃焼状態の計算を行うことができる。
また、クランク角度毎の燃焼に伴う物理量の変化として、燃焼に伴う温度上昇の変化を用い、温度上昇の計算に当たっては、燃焼の熱発生パターンを近似するウィーベ(Ｗｉｅｂｅ) 関数を用いることで、現象モデルによる燃焼状態の計算を燃焼の熱発生パターンを近似して迅速に行うことができる。
また、クランク角度毎の温度上昇の計算は、クランク角度毎の温度上昇の差分を求める計算とすることで、微分方程式を解く場合等と比較して現象モデルによる燃焼状態の計算を迅速に行うことができる。
また、クランク角度毎の物理量の変化の計算は、気筒２０の排気弁の閉成から排気弁の開成までの間において導出することで、現象モデルによる燃焼状態の計算を排気弁の閉成から開成までの間で迅速に完了することができる。
また、燃焼機関１４の性能に関連したパラメータは、燃焼機関１４への掃気系、排気系の各部の圧力、温度、流量、気筒からの排気エネルギー、気筒２０の気筒内圧力、気筒内温度、及び燃焼機関１４のトルク、回転数、ガバナー状態、燃料消費量、負荷の少なくとも１つを含むものであることで、燃焼機関の性能評価に必要なこれらのパラメータに関連した性能予測やモデルの改善を行うことができる。
また、導出された燃焼機関の性能に関連したパラメータに基づいて、燃焼機関１４の状態表示、燃焼機関１４の状態判断、燃焼機関１４の制御の少なくとも１つを行うことで、得られたパラメータを有効に活用することができる。
【0051】
次に、本実施形態による燃焼機関の実時間性能予測の計算方法について説明する。
１．モデルの構成
［１．１ 推進プラントのシステム解析］
システム分析の方法（非特許文献８：Koz’minykh A.V. Fundamentals of ship propulsion plant system analysis. (in Russian). Odessa National Maritime Academy; 2000.）に従って、検討中のシステムは階層的にいくつかの下位レベルの部位に分解される。次に、設計情報等と物理変数を決定して、各部の相互接続を確立する。さらに、すべての部位は、入出力関係に関して一般的かつ再構成可能な数学モデルによって記述される有限数の構成要素に分解される。分解の深さは、必要な詳細レベルと与えられる情報量によって異なる。
このように、従来の推進プラントは、下式（４）に従って、軸の回転運動を介して連動するプロペラとエンジン（燃焼機関）という２つの主要な部位で構成されていると見なされる。
ここで、θeとθpはそれぞれエンジンとプロペラのトルクであり、Ｉeはエンジン、プロペラ、シャフトライン及び追加された水の質量を含む総慣性モーメントであり、ｎeはエンジン回転速度である。
エンジントルクθeは、１サイクルの間にシリンダ内に発生した正味平均有効圧力（ＢＭＥＰ）Ｐbの結果であり下式（５）で表される。
次のステップは、エンジンとプロペラの部位をシステム解析して、式（４）で必要なトルクを求めることである。
１．モデルの構成
［１．１ 推進プラントのシステム解析］
システム分析の方法（非特許文献８：Koz’minykh A.V. Fundamentals of ship propulsion plant system analysis. (in Russian). Odessa National Maritime Academy; 2000.）に従って、検討中のシステムは階層的にいくつかの下位レベルの部位に分解される。次に、設計情報等と物理変数を決定して、各部の相互接続を確立する。さらに、すべての部位は、入出力関係に関して一般的かつ再構成可能な数学モデルによって記述される有限数の構成要素に分解される。分解の深さは、必要な詳細レベルと与えられる情報量によって異なる。
このように、従来の推進プラントは、下式（４）に従って、軸の回転運動を介して連動するプロペラとエンジン（燃焼機関）という２つの主要な部位で構成されていると見なされる。
【数4】
エンジントルクθeは、１サイクルの間にシリンダ内に発生した正味平均有効圧力（ＢＭＥＰ）Ｐbの結果であり下式（５）で表される。
【数5】
【0052】
［１．２ プロペラのシステム解析］
船舶推進シミュレーションに関しては、一般に、トルクに関するプロペラ性能のシミュレーションは準定常的アプローチが採用されている。これはプロペラへの水の流入速度の変動部分も考慮した、プロペラ単独特性のマップの形でのトルク特性の表現に基づいている（非特許文献９：Taskar B., Yum K.K., Pedersen E., Steen S. Dynamics of a Marine Propulsion System With a Diesel Engine and a Propeller Subject to Waves. In: Proceedings of the ASME 2015 34th International Conference on Ocean, Offshore and Arctic Engineering. (OMAE 2015, Newfoundland, Canada).）。本実施形態ではプロペラの詳細な考察は対象外であるので、下式（６）のようにエンジン最大連続定格（ＭＣＲ）点を通る単純なプロペラ二乗則を代わりに使い、次のように考察する。
ここで、ｃp＝θe0／ｎ2
e0はＭＣＲポイントで決定される公称プロペラ係数である。
船舶推進シミュレーションに関しては、一般に、トルクに関するプロペラ性能のシミュレーションは準定常的アプローチが採用されている。これはプロペラへの水の流入速度の変動部分も考慮した、プロペラ単独特性のマップの形でのトルク特性の表現に基づいている（非特許文献９：Taskar B., Yum K.K., Pedersen E., Steen S. Dynamics of a Marine Propulsion System With a Diesel Engine and a Propeller Subject to Waves. In: Proceedings of the ASME 2015 34th International Conference on Ocean, Offshore and Arctic Engineering. (OMAE 2015, Newfoundland, Canada).）。本実施形態ではプロペラの詳細な考察は対象外であるので、下式（６）のようにエンジン最大連続定格（ＭＣＲ）点を通る単純なプロペラ二乗則を代わりに使い、次のように考察する。
【数6】
【0053】
［１．３ エンジンのシステム解析］
商船の大多数は、原動機として低速２ストローク舶用ディーゼルエンジンを利用している。エンジンモデルの目的は発生する正味平均有効圧力（ＢＭＥＰ）に関するエンジンの外部特性を表すことであり、それは一般に回転速度、空気質量流量、燃料質量流量などのエンジン状態の関数である（非特許文献１０：Xiros N. Robust control of diesel ship propulsion. Springer; 2002.）。必要な状態は、後述するようにサイクル平均値（ＣＭＶ）モデルにおいて生成される。
サイクル平均値（ＣＭＶ）モデルで使用される仮定は、集中定数モデルを形成する有限数の制御体積及び抵抗要素にエンジンを分解することができるということである。図４に示すこれらの構成要素は、シリンダ１８、掃気レシーバ１９、排気レシーバ２１、空気冷却器１８、過給機の空気圧縮機１７及タービン２３である。
エンジンの状態を記述するのに必要な基本方程式は、上述の構成要素に適用される熱力学的法則から得ることができる。したがって、掃気レシーバ１９内の空気圧は、次の式（７）の形式で質量流量のバランスから求められる。
空気冷却器１８が掃気レシーバ１９内の空気の温度Ｔsを制御する。後者の概念の当然の結果である仮定は、レシーバ内空気温度がゆっくり変化しているので、定常状態の値と見なされるということである。次に、空気冷却器１８の出口温度は、冷却器の有効性と相対空気質量流量を考慮して次の式（８）のように計算される。
ここで、Ｔwは空気冷却器１８への冷却水入口の温度を示し、ηacは一定の空気冷却器効率である。
次の式（９）に従って、空気圧縮機１７を出る空気の温度Ｔcは、コンプレッサーの等エントロピー効率の定義を使用して評価される（非特許文献１１：Slobodyanyuk L.I., Polyakov V.I. Ship’s steam and gas turbines principles of operation. (in Russian). Shipbuilding, Leningrad; 1983.）。
商船の大多数は、原動機として低速２ストローク舶用ディーゼルエンジンを利用している。エンジンモデルの目的は発生する正味平均有効圧力（ＢＭＥＰ）に関するエンジンの外部特性を表すことであり、それは一般に回転速度、空気質量流量、燃料質量流量などのエンジン状態の関数である（非特許文献１０：Xiros N. Robust control of diesel ship propulsion. Springer; 2002.）。必要な状態は、後述するようにサイクル平均値（ＣＭＶ）モデルにおいて生成される。
サイクル平均値（ＣＭＶ）モデルで使用される仮定は、集中定数モデルを形成する有限数の制御体積及び抵抗要素にエンジンを分解することができるということである。図４に示すこれらの構成要素は、シリンダ１８、掃気レシーバ１９、排気レシーバ２１、空気冷却器１８、過給機の空気圧縮機１７及タービン２３である。
エンジンの状態を記述するのに必要な基本方程式は、上述の構成要素に適用される熱力学的法則から得ることができる。したがって、掃気レシーバ１９内の空気圧は、次の式（７）の形式で質量流量のバランスから求められる。
【数7】
【数8】
次の式（９）に従って、空気圧縮機１７を出る空気の温度Ｔcは、コンプレッサーの等エントロピー効率の定義を使用して評価される（非特許文献１１：Slobodyanyuk L.I., Polyakov V.I. Ship’s steam and gas turbines principles of operation. (in Russian). Shipbuilding, Leningrad; 1983.）。
【数9】
【0054】
上述のように、サイクル平均値（ＣＭＶ）モデルにおける中心的な仮定は、エンジンシリンダが同等のオリフィスによって表すことができ、それが所与の圧力比に対して同じ質量流量を生み出すことができるということである。これにより、オリフィスを通る準一次元流方程式を以下の式（１０）の形式で採用することが可能になる。
ここで、Ａv（Ａは上部に「～」付）は吸気ポートと排気弁を含むエンジンシリンダの有効等価面積である。それは吸気ポートと排気弁の瞬時面積（非特許文献１０）を使って、あるいは利用可能な実験データから見積もることができる。
掃気レシーバ１９と同様に、排気レシーバ２１内のガス状態は、次の式（１１）のように質量とエネルギーのバランスと理想ガスの法則を使用して計算される。
【数10】
掃気レシーバ１９と同様に、排気レシーバ２１内のガス状態は、次の式（１１）のように質量とエネルギーのバランスと理想ガスの法則を使用して計算される。
【数11】
【0055】
燃料質量流量は、サイクル及びシリンダ当たりに噴射される燃料の量がエンジン速度及び燃料ポンプラック位置の一次関数であると仮定して下式（１２）のように評価される。
【数12】
【0056】
ここで、古典的なサイクル平均値（ＣＭＶ）手法とより高度な本実施形態による統合サイクル平均値（ＣＭＶ）モデルとの間の著しい相違、すなわちシリンダから流れる排気ガスのエネルギー率に論点が到達した。
古典的なアプローチ（非特許文献１０、非特許文献１２：Theotokatos G. On the cycle mean value modelling of a large two-stroke marine diesel engine. Proc. IMechE, vol. 224 Part M: J. Engineering for the Maritime Environment, p.193-205, 2010.）では、排気ガスのエネルギー率は燃焼による掃気のエネルギー率の増分と考えられ、次の式（１３）のように与えられた。
反対に、シリンダを出るエネルギー率は、以下の式（１４）に説明するように、静止閉鎖系（エンジンシリンダ）に適用されるエネルギー保存則から直接得ることができる。
ここで、Ｗiは１サイクルのエンジンシリンダ図示パワーであり、Ｑwは１サイクルの熱損失率である。これら２つの量はエンジンサイクル計算の結果であり、後に詳細に説明する。
古典的なアプローチ（非特許文献１０、非特許文献１２：Theotokatos G. On the cycle mean value modelling of a large two-stroke marine diesel engine. Proc. IMechE, vol. 224 Part M: J. Engineering for the Maritime Environment, p.193-205, 2010.）では、排気ガスのエネルギー率は燃焼による掃気のエネルギー率の増分と考えられ、次の式（１３）のように与えられた。
【数13】
【数14】
【0057】
上記の２つの方法（式（１３）及び式（１４））の類似性、すなわち古典的なサイクル平均値（ＣＭＶ）方法におけるサイクルシミュレーションは、排気ガス中に残っている燃料化学エネルギーの割合を表す係数ζaによって置換されることに気付くであろう。この係数は、次の式（１５）の形式の線形関数を使用して正味平均有効圧力（ＢＭＥＰ）と相関している（非特許文献１０）。
【数15】
【0058】
最後に重要なことを述べると、少なくとも２つのモデルに共通して説明されるモデルの一部は、過給機のタービンと圧縮機の性能である。過給機は通常、圧縮機及びタービン性能マップを使用してモデル化されている。ただし、この設計独自データは通常第三者には利用できない。代わりに、解析モデルと経験的相関の組み合わせが過給機の性能を表すために使用される。
式（７）及び式（８）で必要とされる圧縮機の空気質量流量Ｇcは、下式（１６）のように空気圧縮に必要な等エントロピー仕事の方程式から計算することができる。
式（７）及び式（８）で必要とされる圧縮機の空気質量流量Ｇcは、下式（１６）のように空気圧縮に必要な等エントロピー仕事の方程式から計算することができる。
【数16】
【0059】
次に、圧縮機の動力Ｗcはタービンと圧縮機間の動的動力伝達としてモデル化され、次の式（１７）の形式の一次システムで近似される（非特許文献１３：Samokhin S., Sarjovaara T., Zenger K., Larmi M. Modelling and Control of Diesel Engines with a High-Pressure Exhaust Gas Recirculation System. In: Proceedings of the 19th IFAC World Congress. Cape Town, Africa, p. 3006-3011, 2014.）。
ここで、ＷTはタービン内のガス膨張による出力で次の式（１８）のように定義する。
【数17】
【数18】
【0060】
タービンを通過する排気ガスの質量流量は、エンジンのシリンダと同様に計算される。等価オリフィスを考慮し、次の式（１９）に従って準１次元アプローチを適用する。
【数19】
【0061】
幾何学的面積ＡTを修正するタービン流量係数μは、下式（２０）のようにタービン圧力比との関数関係を保持する（非特許文献１０、１１）。
【数20】
【0062】
通常、定常状態下での圧縮機動作点は圧縮機マップ上の単一の曲線上にあり、したがって圧縮機等エントロピー効率ηiCは一定とすることができる。タービンの等エントロピー効率ηiTはタービンの特性速度Ｕt/Ｃsに強く依存するが（非特許文献１０、１１）、次の経験式である式（２１）を採用することができる（非特許文献１４：Ray A. Dynamic Modelling of Power Plant Turbines for Controller Design. Applied Mathematical Modelling, vol. 4, issue 2; 1980.）。
【数21】
【0063】
式（２１）は過給機回転速度についての情報を必要とするが、提示されたサイクル平均値（ＣＭＶ）モデルではそれは明確に定義されていない。それは、圧縮機の相対等エントロピー仕事量と相対回転速度との間の二乗則関係から、次の式（２２）の形で推定することができる（非特許文献１５：Bondarenko O., Fukuda T. Development of Diesel Engine Simulator for Use with Self-Propulsion Model. J. Japan Institute of Marine Engineers. Vol. 48, p. 98-105; 2013.）。
【数22】
【0064】
最後に、エンジンの正味平均有効圧力（ＢＭＥＰ）は指示平均有効圧力（ＩＭＥＰ）Ｐiと摩擦平均有効圧力（ＦＭＥＰ）Ｐfとの間の差として計算される。ここで言及しておく価値があるのは、ＩＭＥＰ計算に関連した、従来のＣＭＶと本実施形態による統合ＣＭＶアプローチのもう１つの重要な違いである。従来のサイクル平均値（ＣＭＶ）モデルでは、指示平均有効圧力（ＩＭＥＰ）は燃料ポンプラック位置Ｆpに比例すると考えられ、下式（２３）のように不完全燃焼を考慮した燃焼効率ηcで修正される（非特許文献１０）。
顕著なターボチャージャーラグ効果を伴う大きな過渡的遷移のみが１を下回る燃焼効率の低下を引き起こす可能性があるため、これは残りのエンジン状態計算との弱い関係を意味する（非特許文献１６：Rakopoulos C.D., Giakoumis E.G. Diesel Engine Transient Operation. Springer; 2009.）。
【数23】
【0065】
本実施形態による燃焼機関の実時間性能予測（統合サイクル平均値（ＣＭＶ）手法）では指示平均有効圧力（ＩＭＥＰ）がサイクルシミュレーションの結果であり、サイクルシミュレーションは過給機の性能によって直接決定されるシリンダ内の初期条件に依存する。
ＦＭＥＰ計算では、ChenとFlynnの摩擦相関（非特許文献１７：Chen S.K., Flynn P. Development of a compression ignition research engine. SAE Paper No. 650733; 1965.）が選択され、実験的測定値と一致する摩擦損失が得られるように倍率が使用された。
ＦＭＥＰ計算では、ChenとFlynnの摩擦相関（非特許文献１７：Chen S.K., Flynn P. Development of a compression ignition research engine. SAE Paper No. 650733; 1965.）が選択され、実験的測定値と一致する摩擦損失が得られるように倍率が使用された。
【0066】
２．エンジンシリンダーモデリング
［２．１ クローズドサイクルの説明］
全エンジンサイクルは、ガス交換（新規空気導入及び燃焼ガス排出）、圧縮、燃焼及び膨張のようないくつかの段階からなる。本実施形態による燃焼機関の実時間性能予測では、サイクルのガス交換部分を特徴付けるパラメータは、（前述のように）１サイクルにわたって平均化された連続変数としてモデル化される。サイクルの残りのパラメータは、定常開放系に適用されるゼロ次元熱力学的アプローチを考慮して計算される（非特許文献７）。このような場合、質量とエネルギーに関する保存則のみが、作動媒体の状態が空間的に均質であり時間と共に変化する理想気体であると仮定して考慮される。
［２．１ クローズドサイクルの説明］
全エンジンサイクルは、ガス交換（新規空気導入及び燃焼ガス排出）、圧縮、燃焼及び膨張のようないくつかの段階からなる。本実施形態による燃焼機関の実時間性能予測では、サイクルのガス交換部分を特徴付けるパラメータは、（前述のように）１サイクルにわたって平均化された連続変数としてモデル化される。サイクルの残りのパラメータは、定常開放系に適用されるゼロ次元熱力学的アプローチを考慮して計算される（非特許文献７）。このような場合、質量とエネルギーに関する保存則のみが、作動媒体の状態が空間的に均質であり時間と共に変化する理想気体であると仮定して考慮される。
【0067】
下式（２４）のようにエネルギー方程式又は熱力学の第一法則は、運動エネルギーを無視してシリンダ内の内部エネルギーの変化を提供する。
【数24】
【0068】
下式（２５）のように作動媒体とシリンダ壁との間の熱流ｄＱwは、一定の平均シリンダ壁温度を仮定して、対流熱伝達のための標準式を用いて考慮される。
主な関心はサイクル中に伝達される熱であるので、平均熱伝達係数αが使用される。計算には、Woschni又はHohenbergによって提供された経験的相関を使用することができる（非特許文献１９：Merker GP, Schwarz C, Stiesch G, Otto F. Simulating combustion. Springer; 2006.）。
燃料燃焼による熱の流れｄＱfは次の式（２６）の通りである。
ここで、熱発生率dqx＝mf.cdxcは、その単純さと汎用性のためにゼロ次元アプローチで広く使用されているウィーベ（Ｗｉｅｂｅ）関数に従ってモデル化されている（非特許文献２０：Ghojel J.I. Review of the development and applications of the Wiebe function: a tribute to the contribution of Ivan Wiebe to engine research. Int. J. Engine Res. Vol. 11, p. 297-312; 2010）。無次元燃焼速度dxcに対するウィーベ関数は次の式（２７）のように書くことができる。
ここで、Φ＝（φ－φsoi－φid）/φzは正規化燃焼期間を表す。
これらのパラメータＣ，ｍとφzは燃焼速度の形状を特徴付けるものであり、これはエンジンのあらゆる動作点に適合させるべきである。しかしながら、実際的な考察のために、エンジン運転条件との上述のパラメータの様々な相関関係が導入されており（非特許文献４、非特許文献１８：Medica V., Simulation of turbocharged diesel engine driving electrical generator under dynamic working conditions [dissertation]. Rijeka, Croatia: University of Rijeka; 1988.、非特許文献１９）、それらの相関関係における定数が較正パラメータとして考慮されている。さらに、解離過程の損失及び不完全な燃料燃焼を考慮に入れるために、燃焼効率ηcが式（２６）に導入される。
【数25】
燃料燃焼による熱の流れｄＱfは次の式（２６）の通りである。
【数26】
【数27】
これらのパラメータＣ，ｍとφzは燃焼速度の形状を特徴付けるものであり、これはエンジンのあらゆる動作点に適合させるべきである。しかしながら、実際的な考察のために、エンジン運転条件との上述のパラメータの様々な相関関係が導入されており（非特許文献４、非特許文献１８：Medica V., Simulation of turbocharged diesel engine driving electrical generator under dynamic working conditions [dissertation]. Rijeka, Croatia: University of Rijeka; 1988.、非特許文献１９）、それらの相関関係における定数が較正パラメータとして考慮されている。さらに、解離過程の損失及び不完全な燃料燃焼を考慮に入れるために、燃焼効率ηcが式（２６）に導入される。
【0069】
式（２７）で導入された着火遅れφidは下式（２８）の修正トルストフの式から決定される（非特許文献２１：Kuleshov A.S. Multi-Zone DI Diesel Spray Combustion Model for Thermodynamic Simulation of Engine with PCCI and High EGR Level. SAE Paper No. 2009-01-1956; 2009.）。
ここで、Ｅa＝２３０００Ｋ ２８０００kJ/kmoleは自己着火プロセスの見かけの活性化エネルギーである。
式（２８）で導入された補正係数ＣTは、燃料噴射がＴＤＣの後に行われる場合には負であり得る遅延期間中の温度勾配を説明する。
シリンダ内のガスの膨張又は圧縮によってピストンに伝達される仕事は、次の式（２９）ように評価される。
シリンダの体積変化ｄＶはピストンの運動学から計算される（非特許文献１９）。
【数28】
式（２８）で導入された補正係数ＣTは、燃料噴射がＴＤＣの後に行われる場合には負であり得る遅延期間中の温度勾配を説明する。
シリンダ内のガスの膨張又は圧縮によってピストンに伝達される仕事は、次の式（２９）ように評価される。
【数29】
【0070】
前述した式の系及び質量保存則で補足された理想気体の式は、ピストンの圧力、温度及び運動の間の特定の関係を設定せず、代わりに、時間ステップｄｔ又はクランク角度ステップｄφのいずれかにおけるパラメータの相対変化を評価することを可能にする。しかしながら、解の許容可能な精度を得るためには、小さい時間ステップを有するルンゲクッタ法又はオイラー法のいずれかが必要であり、両方とも時間がかかる。次に、方程式の分解能を管理する別の方法について説明する。
【0071】
［２．２ 微分方程式］
下式（３０）のように２つの有限状態間の内部エネルギーの変化を考慮して、式（２４）の熱力学の第１法則を積分形式で書き換えることから検討を始める。
内部エネルギーＵは、平均熱容量Ｃv（Ｃは上部に「￣」付）から温度の関数として評価することができる。さらに、作動媒体は給気と化学量論的燃焼生成物との均質混合物と考えられる。よって下式（３１）で表される。
ここで、空気及び燃焼生成物の比熱容量は、表にしたデータの回帰関数で表すことができ（非特許文献２２：Grigor’eva V.A., Zorina V.M., editors. Theoretical foundations of thermal engineering. (in Russian). Energoatomizdat, Moscow, Vol. 2; 1988.）、ガス混合物の熱容量は相加率として評価される。よって下式（３２）で表される。
ここで、ｒは燃焼生成物の比率である。
下式（３０）のように２つの有限状態間の内部エネルギーの変化を考慮して、式（２４）の熱力学の第１法則を積分形式で書き換えることから検討を始める。
【数30】
【数31】
【数32】
【0072】
初期状態での熱容量は、気筒内のガスの状態によって明確に決定されるが、遷移終了時のガスの状態は不明である。遷移過程におけるガス状態の小さな変化を仮定すると、高次項を無視して、熱容量に対するテイラー展開は次の式（３３）のようになる。
同じ仮定の下で、状態遷移の結果としてピストンに伝達された仕事は次の式（３４）のように評価されることができる。
考慮された間隔の中央である平均圧力ｐc（ｐは上部に「￣」付）は、下式（３５）のように平均パラメータを同様に考慮して理想気体の法則から評価できる。
式（３５）を式（３４）に代入して追加の変数を導入すると、下式（３６）のようにシリンダ内のガスの仕事のための式が得られる。
式（３５）の混合物のガス定数Ｒは、式（３２）と同様に加成率として評価される。
ここで、式（３２）と式（３３）で定義された熱容量を使って、下式（３７）のように遷移の両端での内部エネルギーを定義する。
見て分かるように、遷移及び増分の始めの状態のみが式（３６）及び式（３７）に導入される。
【数33】
【数34】
【数35】
【数36】
ここで、式（３２）と式（３３）で定義された熱容量を使って、下式（３７）のように遷移の両端での内部エネルギーを定義する。
【数37】
【0073】
最後に、式（３６）と式（３７）を式（３０）に代入すると、下式（３８）のように変換により未知の温度増分ΔＴに関して二次方程式を生成する（非特許文献２３：Mizernyuk G.N., Kuleshov A.S. Computational evaluation of internal combustion engine working process. (in Russian). In: News of Higher Education, Mashinostroenie, Moscow; 1986.）。
ここで、
従って、ＥＶＣにおけるガスの初期状態から出発してクランクシャフトを角度Δφだけ前進させると、下式（４０）のように作動媒体の新しい状態は非線形代数方程式の式（３８）のみを解くことによって評価することができる。
最後に重要なことは、提示された計算方法はオイラー法と比較して比較的大きな計算ステップで高い精度を維持することである。さらに、計算手順は計算ステップの柔軟な調整を可能にする。したがって、サイクルのさまざまな段階を、精度を損なうことなく可変ステップで計算することができる。
【数38】
【数39】
【数40】
【実施例】
【0074】
上記説明からわかるように、完全なモデルには、特定のエンジンに合うように調整する必要がある一連の経験的パラメータが含まれている。必要なデータは、エンジンプロジェクトガイド、エンジンショップテストデータ、及び海上試験データから取得できる。本発明のために、データは株式会社三井Ｅ＆Ｓマシナリーによって提供され、そこでテストエンジンは様々な負荷で運転された。試験エンジンの仕様を表１に示す。
【表1】
【0075】
モデルチューニングはＣＭＶ部と燃焼部を別々に行った。またＣＭＶモデルの較正も２段階で行われた。これは、モデルに過給機のタービンの未知の特性、すなわちタービン流量係数と等エントロピー効率が含まれているためである。モデル予測値と実際のデータとの間の残差を最小化することによって、必要な変数、すなわち排気ガス圧力、タービン有効流路面積及び効率がエンジンのあらゆる動作点で見出された。その後、回帰関数の定数は、ニュートン - ラプソン法のファミリーを使用して全ての利用可能なデータにわたって調整された。
【0076】
式（２７）の燃焼モデルは、エンジンの全動作点においてウィーベ関数の形状を決定する３つのパラメータを必要とする。非特許文献１９で提案された経験的相関のパラメータは、エンジンのすべての動作点に対して燃焼モデルの良好な適合性を提供し、またｂｓｆｃとの良好なトレードオフを提供するように調整された。パラメータ最適化プロセスは、非特許文献２４で導入された修正を伴うParticle Swarm Optimizationアルゴリズムの助けを借りて達成された。
【0077】
実験データに対してプロットした統合ＣＭＶモデル（本実施例）の定常状態シミュレーション結果を図７に示す。図７では、実験結果を「〇」印で、本実施例による予測結果を「◆」印を結ぶ線で示している。また、図７（ａ）は、縦軸が燃料消費量（正味燃料消費率：ｂｓｆｃ）であり、横軸が負荷である。図７（ｂ）は、縦軸が空気質量流量、横軸が負荷である。図７（ｃ）は、縦軸が掃気圧、横軸が負荷である。図７（ｄ）は、縦軸が排ガス温度、横軸が負荷である。さらに、図８は相対エラーが報告されている下表２と組み合わせて、統合ＣＭＶモデルの枠組みにおける燃焼性能を示す。図８では、縦軸を気筒内圧力、横軸をクランク角度とし、実験結果を点線で示し、本実施例による予測結果を実線で示している。図８（ａ）は負荷が１００％の場合、図８（ｂ）は負荷が８５％の場合、図８（ｃ）は負荷が７５％の場合、図８（ｄ）は負荷が６０％の場合である。図から分かるように、よく調整されたモデルは、モデルのＣＭＶ部分と燃焼部分の両方について、実際のエンジンとかなりよく一致している（非特許文献７、１２）。ここで再び言及する価値があるのは、このモデルには過給機の詳細な特性が同時に含まれていないことで、空気質量流量や排気ガス温度などの関連する変数を正しく予測することができないということである。
【表2】
【0078】
最後に、本発明の重要な目的は、エンジン全体のリアルタイムシミュレーションに適した計算性能の向上である。本実施例の計算手順を評価するために、簡単なテストを設定した。テストエンジンの幾何学形状と運動学を用いて、ルンゲクッタ（Runge－Kutta）法と差分方程式の提案した方法を比較して、断熱圧縮－膨張をシミュレートした。燃焼と熱伝達のモデルは、関連する不確実性を避けるために無視されている。計算ステップ、特にクランク角は、ＴＤＣで一貫した精度の圧力を提供する０．１°の微小角度から最大の最大角度まで変化した。シミュレーションは科学的シミュレーションソフトウェアScilab（https://www.scilab.org/ [accessed 11.02.2019]）（Matlabと同様）で行い、アルゴリズムの実行時間は組み込み関数tic（）/ toc（）を使って測定した。結果を図９に縦軸を対数目盛で示す。図９では、縦軸を計算時間、横軸をクランク角度ステップとし、ルンゲクッタ法の結果を「●」印を結ぶ線で示し、本実施例による予測結果を「◆」印を結ぶ線で示している。見て分かるように、本実施例による燃焼機関の実時間性能予測法の計算速度はルンゲクッタ法よりもおよそ一桁優れている。さらに、シリンダ内サイクルの透明な計算方式により、サイクルの各段階に対する可変ステップは計算の加速をもたらす。
【産業上の利用可能性】
【0079】
本発明による燃焼機関の実時間性能予測方法及び実時間性能予測プログラムは、船舶の実海域における船舶の燃焼機関の状態予測、状態表示、状態判断、又制御等に用いることができる。また、船舶の燃焼機関のみならず一般的な燃焼機関にも展開が可能なものである。
【符号の説明】
【0080】
１３ 過給機
１４ 燃焼機関
１９ 掃気レシーバ
２０ 気筒
２１ 排気レシーバ
Ｓ１ 初期条件入力ステップ
Ｓ２ モデル設定ステップ
Ｓ３ サイクル平均値（ＣＭＶ）モデル計算ステップ
Ｓ４ 現象モデル計算ステップ
Ｓ９ 出力ステップ
１４ 燃焼機関
１９ 掃気レシーバ
２０ 気筒
２１ 排気レシーバ
Ｓ１ 初期条件入力ステップ
Ｓ２ モデル設定ステップ
Ｓ３ サイクル平均値（ＣＭＶ）モデル計算ステップ
Ｓ４ 現象モデル計算ステップ
Ｓ９ 出力ステップ
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
