JP7376886B2 - 分子反応モデリング方法 - Google Patents

Description

特許法第３０条第２項適用 １．集会名 平成３０年度ＪＰＥＣフォーラム ２．開催日 平成３０年５月９日 ３．公開者 高橋祐樹、辻浩二及び三谷尚洋

特許法第３０条第２項適用 １．掲載アドレス ｈｔｔｐ：／／ｗｗｗ．ｐｅｃｊ．ｏｒ．ｊｐ／ｊａｐａｎｅｓｅ／ｒｅｐｏｒｔ／２０１８ｒｅｐｏｒｔ／ｉｎｄｅｘ．ｐｄｆ ２．掲載日 平成３０年５月９日 ３．公開者 高橋祐樹、辻浩二及び三谷尚洋

特許法第３０条第２項適用 １．刊行物名 平成２９年度 高効率な石油精製技術の基礎となる石油の構造分析・反応解析等に係る研究開発事業 事業報告書 ２．発行日 平成３０年３月３０日 ３．公開者 一般財団法人石油エネルギー技術センター

特許法第３０条第２項適用 １．刊行物名 平成３０年度 高効率な石油精製技術の基礎となる石油の構造分析・反応解析等に係る研究開発事業 事業報告書 ２．発行日 平成３１年３月２９日 ３．公開者 一般財団法人石油エネルギー技術センター

特許法第３０条第２項適用 １．刊行物名 石油学会創立６０周年記念東京大会（特別講演，招待講演，第４８回石油・石油化学討論会）（講演要旨）、第１９１頁及び第２０４頁，石油学会 ２．発行日 平成３０年１０月１７日 ３．公開者 高橋祐樹、辻浩二及び三谷尚洋

特許法第３０条第２項適用 １．集会名 石油学会創立６０周年記念東京大会（特別講演、招待講演、第４８回石油・石油化学討論会） ２．開催日 平成３０年１０月１８日 ３．公開者 高橋祐樹、辻浩二及び三谷尚洋

特許法第３０条第２項適用 １．刊行物名 ペトロテック，第４２巻，第１号，第２３～２８頁，公益社団法人 石油学会 ２．発行日 平成３１年１月１日 ３．公開者 三谷尚洋

特許法第３０条第２項適用 １．掲載アドレス ｈｔｔｐｓ：／／ｗｗｗ．ｊｓｔａｇｅ．ｊｓｔ．ｇｏ．ｊｐ／ｂｒｏｗｓｅ／ｓｅｋｉｙｕ／ｌｉｓｔ／ｃｈａｒ／ｊａ ２．掲載日 平成３０年１２月１７日 ３．公開者 高橋祐樹、辻浩二及び三谷尚洋

特許法第３０条第２項適用 １．集会名 ２０１８ ＡＩＣｈＥ Ａｎｎｕａｌ Ｍｅｅｔｉｎｇ ２．開催日 平成３０年１０月２９日 ３．公開者 グェン・テイ・ホン・テュイ、片岡祥、高橋祐樹、辻浩二及び田中隆三

特許法第３０条第２項適用 １．掲載アドレス ｈｔｔｐｓ：／／ｗｗｗ．ａｉｃｈｅ．ｏｒｇ／ｃｏｎｆｅｒｅｎｃｅｓ／ａｉｃｈｅ－ａｎｎｕａｌ－ｍｅｅｔｉｎｇ／２０１８ ２．掲載日 平成３０年１０月 ３．公開者 グェン・テイ・ホン・テュイ、片岡祥、高橋祐樹、辻浩二及び田中隆三

本発明は、ランプ反応モデルに基づく石油の分子反応モデリング技術に関する。

石油精製に関する諸装置の運転においては、通常、比重や粘度、蒸留性状(沸点)などの全体の物理的性状に基づいて原料油を分析し、過去の類似のデータを有する油種の運転実績を参考にして運転条件を決めるという手法がとられている。しかしながら、昨今では、輸入原油種が多様化しており、類似する過去のデータを探すことは容易ではない。さらに運転効率の向上や環境保護という面からも、単純に過去の運転実績を踏襲すればよいというものではなくなっている。

そこで、比重や粘度、蒸留性状というような石油全体を一括りにした観点で捉えるのではなく、石油を構成している炭化水素分子というレベルでその化学構造や存在割合を把握し、それにより得られた推定物性値等の知見に基づいて運転条件を設定することができれば、より客観性に基づいた効率的な運転ができると考えられてきた。そのため、石油業界においては、石油を分子レベルで把握する技術の出現が待ち望まれてきた。

一方、石油は膨大数の炭化水素分子からなる混合物であり、特に重質油は分子量が大きく、かつ複雑な化学構造を有する分子が極めて多種類存在する。そのため、分子の一つ一つについて化学構造を特定し、それらの存在割合をも特定するというのは、非常に困難なことであった。

そのため、石油を分子レベルで分析し化学構造を解析するにあたっては、高分解能質量分析装置であるフーリエ変換イオンサイクロトロン共鳴方式による質量分析計（ＦＴ－ＩＣＲ－ＭＳ）を用いて高精度に分子量を計測する技術が用いられてきた（特許文献１、２等）。また、石油精製プロセスの反応モデルでは、沸点等の一般性状で括られた擬似成分で反応を扱うモデルが主として用いられてきた。これはランピングモデルと呼ばれており常圧残油（ＡＲ）のような残渣油等の原料油を分子組成ではなく、沸点ごとに括られたランプ（仮想的に括った単一成分）として取り扱い、ランプ間の反応速度から物質収支を推定するモデルである。ランピングモデルは、製油所でも取得可能な一般性状を用いるため取り扱いが容易であるという利点がある（非特許文献１、２等）。

特表２０１４－５００５０６号公報 特表２０１４－５０３８１６号公報

Ｄ．Ｍ．Ｃａｍｐｂｅｌｌ，Ｃ．Ｂｅｎｎｅｔｔ，Ｚ，Ｈｏｕ，Ｍ．Ｔ．Ｋｌｅｉｎ，Ｉｎｄ．Ｅｎｇ．Ｃｈｅｍ Ｒｅｓ．，４８，１６８３（２００９） Ｓ．Ｒ．Ｈｏｒｔｏｎ，Ｌ．Ｚｈａｎｇ，Ｚ，Ｈｏｕ，Ｃ．Ｂｅｎｎｅｔｔ，Ｍ．Ｔ．Ｋｌｅｉｎ，Ｓ．Ｚｈａｏ，Ｉｎｄ．Ｅｎｇ．ＣｈｅｍＲｅｓ．，５４，４２２６（２０１５）

例えば、常圧残油水素化脱硫（ＲＤＳ）は、常圧残油（ＡＲ）中の硫黄元素を除去する脱硫反応が主たる目的の装置であるが、脱硫（ＨＤＳ）触媒の被毒となるバナジウム（Ｖ）やニッケル（Ｎｉ）の除去を行うと共に、芳香環の核水添反応と核水添反応を受けた含窒素化合物からの脱窒素反応が並行して起こる。現状では、実装置のように触媒劣化に基づく反応性の変化が生じた場合には、都度、反応試験を行い、反応速度定数を求める必要がある。この点を改良するため、触媒の劣化に伴い、反応性がどのように変化するか、その際触媒が受けるダメージはどのようなものかという、連続運転における触媒が受ける影響と反応性の相関について解析を行い、同一触媒であれば、運転時間の経過に伴う反応性の変化を組込んだ分子反応モデリング技術が必要とされている。

実際にＲＤＳ原料油・生成油のＦＴ－ＩＣＲ－ＭＳを測定すると、数万から数十万のピークが検出され、膨大な数の分子が含まれていることがわかる。しかも、これらの分子は芳香環、ナフテン環、ヘテロ環で構成されるコアと、主にアルキル基で構成される側鎖や架橋が結合した複雑な構造を有していると考えられる。そのため分子及び反応経路の数が膨大となり、分子組成をベースとした分子反応モデリングは非常に困難である。

そこで、本発明は、ＦＴ－ＩＣＲ－ＭＳ分析結果とランピングモデルとを組み合わせ、先ずコアのランピング及びランピングした反応ネットワークを定義し、コンピュータを用いて定義された反応ネットワークを簡略化し、それらの結果に基づいて分子の反応を描写し、構造に基づき生成物を推定する速度モデルを効率よく構築する新たな手法を提供することを目的とするものである。

上記の目的を達成するため、本発明者らは、以下の本発明を創出した。即ち、本発明によるコンピュータによる石油の分子反応モデリング方法は、
（１）基本コア、ヘテログループ及びナフテン環数から成るコア因子に基づいて、基本コアの構造を比較した後に、ヘテログループの構造を比較することにより、前記構造に基づくランプへ、コアのランピングを行うステップと、
（２）前記ステップ（１）で得られたランプに含まれるコアの反応を追跡することにより、前記ランプに基づいたコアの反応ネットワークを構築するステップと、
（３）前記ステップ（２）で得られたランプに基づいたコアの反応の重要度の判定を、
（Ａ）前記反応モデリングに供する原料のコアについて、該コアが含まれる全反応における個々の反応の影響度を、下式（I）：
（式中、
Ｃ_ｉ，ｔ：は、反応時間ｔにおいて、反応ｉに存在する原料のモル濃度であり、
ｋ_ｉは、反応ｉの反応速度定数であり、
Ｎは、反応ネットワークに存在する反応の総数、である）
により算出し、各反応のＧＬＯ_ｔを比較し、ＧＬＯ_ｔが最も小さい反応を除外し、
（Ｂ）前記工程（Ａ）の影響度の算出を、下式（II）：
（式中、
ｋ_ｒは、反応ネットワークから除外された反応速度定数であり、
Ｃ_ｒは、反応ネットワークから除外された原料コアの濃度であり、
αは、０～１である）
により表される条件を満たさなくなるまで、繰り返し、
（Ｃ）反応時間ｔ_ｑ及びｔ_ｑ－１における各原料のモル流量の差分を算出し、前記ステップ（Ａ）において除外された反応の重要度を、下式（III）：
（式中、
Ｃ_{ｉ，ｔｑ－１}は、反応時間t_q-1において、反応ｉに存在する原料のモル濃度であり、
Ｃ_ｉ，ｔｑは、反応時間t_qにおいて、反応ｉに存在する原料のモル濃度であり、
βは、任意の値である）
により算出したＭＯＬがβよりも大きい場合には、重要なコアとして保持し、
（Ｄ）前記ステップ（Ｃ）において得られた重要なコアを、下式（IV）：
（式中、
Ｃ_ｉは、コアｉのモル濃度であり、
ｋ_ｉは、原料としてのコアｉの反応速度定数であり、
ｋ’_ｉは、生成物としてのコアｉの反応速度定数であり、
ｍ及びｎは、原料及び生成物としてのコアｉを含む反応数
χは、任意の値である）
により算出したＬＯＣがχよりも大きい場合には、前記コアの含まれる反応は、重要であると判定すること、
により行う、ステップと、
（４）前記ステップ（２）で得られた反応ネットワークから、前記各ランプの速度変化を表す速度モデルを構築するステップと、
（５）前記ステップ（４）で得られた各ランプの速度モデルにもとづいた反応速度定数を、最小二乗法を用いて算出するステップと、
を含んでなることを特徴とするものである。

また、本発明の別の実施態様においては、石油の分子反応モデリング装置、並びに石油の分子反応モデリング方法や石油の分子反応モデリング装置を実行するためのコンピュータプログラムも提供される。

本発明によれば、先ずコアのランピング及びランピングした反応ネットワークを定義し、コンピュータを用いて定義された反応ネットワークを簡略化し、それらの結果に基づいて分子の反応を描写し、構造に基づき生成物を推定する速度モデルを効率よく構築する新たな手法を提供することができる。

石油に含まれる分子の反応解析の基本概念を示した概念図である。 石油の分子反応モデリング方法手法の基本フローである。 コアのランピング手法の基本フローである。 類似した構造を有する基本コアの反応パス。 反応ネットワークに含まれる全基本コアのそれぞれの関係。 コアに含まれるヘテログループの反応パス。 ランピングに基づいて構築された反応ネットワーク。 反応の重要度判定の手法を示すフロー図。 ３種の触媒が充填されたＲＤＳ反応器。 ２７のコア及び３０の反応から構築した簡略ランプ。 ２３のコア及び２５の反応から構築した簡略ランプ。 １９のコア及び２１の反応から構築した簡略ランプ。 １４のコア及び１５の反応から構築した簡略ランプ。 反応速度モデル構築のための反応モデル。 モル流量の参考値及び推定値の比較結果（ＨＤＭ＿１０）。 モル流量の参考値及び推定値の比較結果（ＨＤＭ＿１１）。 モル流量の参考値及び推定値の比較結果（ＨＤＳ＿３）。 参考値及び構築したモデルより推定したモル流量の比較結果（Ｌ１～Ｌ３）。 参考値及び構築したモデルより推定したモル流量の比較結果（Ｌ４～Ｌ６）。 参考値及び構築したモデルより推定したモル流量の比較結果（Ｌ７～Ｌ９）。 参考値及び構築したモデルより推定したモル流量の比較結果（Ｌ１０～Ｌ１２）。 参考値及び構築したモデルより推定したモル流量の比較結果（Ｌ１３～Ｌ１４）。

本発明の実施形態を説明するにあたり、先ず、本明細書にて使用する用語ないし表現について説明する。

＜定義＞
「石油」
本明細書において、「石油」とは、原油、並びに原油を蒸留して得られる諸留分及び諸留分に改質や分解等の二次装置による処理を加えて得られる留分等をも含む総称的な概念をいう。或いは、原油を蒸留して得られたある留分について、さらに飽和炭化水素や芳香族炭化水素等の成分に分画した分画物を指すこともある。

「重質油」
本明細書において、「重質油」とは、原油を蒸留して得られる諸留分のうち、常圧残油以上の沸点を有する石油留分を意味する。

＜石油に関する装置＞
本明細書において、「石油に関する装置」とは、蒸留装置や抽出装置をはじめ、改質装置、水素添加反応装置、脱硫装置等の化学反応を伴う装置等、石油の処理に関する装置をすべて含む。「石油に関する装置」を総じて、「石油精製装置」ともいう。

「物性値」
「物性値」とは、上記の方法により特定された分子構造及びその存在割合に基づいて得られる値であって、物質の物理的又は化学的な性質や性状、特性を表現するものであれば、名称の如何に拘わらず、「物性値」に含まれる。本明細書において、「物性値」とは、これらに限定されるものではないが、例えば、融点、ハンセン溶解度指数値、生成ギブス自由エネルギー、イオン化ポテンシャル、分極率、誘電率、蒸気圧、液体密度、API度、気体粘度、液体粘度、表面張力、沸点、臨界温度、臨界圧力、臨界体積、生成熱、熱容量、双極子モーメント、エンタルピー、エントロピー等である。

「成分」
「成分」とは、「混合物をある特定の物理的又は化学的性状を基準として括った塊」、即ち、「ある特定の物理的又は化学的性状を基準として分画された分画物（フラクション）」を意味する。特定の物理的又は化学的性状を基準として括る方法としては、例えば、蒸留試験における沸点範囲を特定して、その温度範囲にあるものを一つの成分として分画する方法等が挙げられる。この場合、混合物は「分画物（フラクション）の集合体」ということになる。或いは、「成分」を、多成分混合物を構成する一つ一つの構成員であって、「同一の分子種に属すると認められる分子の集合体」と捉えてもよい。ここで、「同一の」とは、「分子構造を完璧に特定し、その上で同一である」、或いは、「分子構造上の異性体（分子式は同じであるが構造が異なるもの）同士は同一のものとする」という意味と捉えてもよく、例えば、後述する「ＪＡＣＤのような方式で特定された構造において同一である」という意味と捉えてもよい。さらには、広く「任意に定めた基準に基づいて一括りにした分子の集合体」という意味と捉えてもよい。

「分子」
上記「成分」における「分子」に関し、分子が持つ構造に関する何等かの情報を特定するという行為であれば、あらゆる行為を包含するものである。目的及び必要性に応じて、その度合い、表示の方式を適宜選択すればよい。分子全体の構造を特定するという行為のみならず、分子の一部分についての構造に関する情報を組み込んでもよい。例えば、コア部分の構造のみを特定し、側鎖部分や架橋部分については構造は特定せず分子式のままにしておいてもよい。
本明細書において、好ましくは、後述するＪＡＣＤで分子構造を特定する。「ＪＡＣＤ」で構造が特定された分子というのは、後述するアトリビュートの結合位置の違いによる異性体をすべて含む概念である。本明細書において、「分子」は、異性体をすべて含む概念と捉えてもよい。

「コア」
「コア」とは、後述の「ＪＡＣＤ」の項で記載する「アトリビュート」の一種であって、具体的には、芳香環又はナフテン環そのもの、芳香環とナフテン環が架橋ではなく直接結合しているもの、芳香環又はナフテン環にヘテロ環が架橋ではなく直接結合しているものである。架橋又は側鎖は、コアとは別のアトリビュートであるため、「コア」とは、架橋又は側鎖を一切有しないものを意味している。
一方、「シングルコア」とは、上記コアを１個だけ有する分子を指す概念である。分子を指す概念であるため、コアに側鎖が結合しているものも包含している。上記コアの２個以上が架橋してなる分子を「マルチコア」という。「マルチコア」も分子を意味するため、コアに側鎖が結合しているものも包含している。２個のコアが架橋してなる分子を「ダブルコア」という。例えばナフタレン分子は、１個の芳香環からなるものであるため「シングルコア」であり、ベンゼン環２個からなるダブルコアではない。

「ＪＡＣＤ」
「ＪＡＣＤ（Juxtaposed Attributes for Chemical-structure Description）」とは、分子構造に関する新規な表示方式であって、分子の構造を、アトリビュートの種類及びアトリビュートの数により表示するものである。アトリビュートが他のアトリビュートのいずれの位置において結合しているかについては表示しない。
上記において、「アトリビュート」とは、分子を構成している化学構造上の部品（パーツ）を指す概念である。芳香族化合物においては、具体的には、前述の「コア」、「架橋」及び「側鎖」を指す。
この表示方式によると、石油を構成する膨大数の分子の各々に関し、それらの構造を、必要かつ十分な程度に特定することができる。 以下の化学式で表された分子を例にとって説明する。

この化合物をＪＡＣＤで表すと、以下の表１のようになる。

ＪＡＣＤで表示され、構造が特定された分子とは、アトリビュートの結合位置の違いによる異性体をすべて含む概念である。

「ＪＫＭＴ」
「ＪＫＭＴ」とは、デラウェア大Klein教授のグループが開発した分子反応モデリング技術ＫＭＴ（Kinetic Modelar’s Toolbox）を一般財団法人石油エネルギー技術センター向け仕様に改変したものであり、ＣＭＥ（Composition Model Editor）、ＩＮＧｅｎ（Interactive Network GENerator）及びＫＭＥ（Kinetic Model Editor）の各モジュールからなるものである。ＣＭＥはバルク物性推算と構造・組成モデリングのためのモジュールであり、ＣｏｒｅＧｅｎ、ＦｏｏｔＧｅｎ、ＣｏｍｐＧｅｎ、ＰｒｏｐＧｅｎと呼ばれる４つのサブモジュールからなり、バルクの分析データから分子組成をモデリング化するツールである。ＩＮＧｅｎは反応ネットワーク構築のためのモジュールであり、ＫＭＥはユーザーインターフェースである。反応種の構造タイプと反応ファミリーの組み合わせをグラフィカルユーザーインタフェースで選択し、それに基づいて反応経路とネットワーク、及び生成分子が自動的に算出され、それを構造式作画ソフトＣｈｅｍＤｒａｗ（登録商標）にて確認できるようになっている。

＜石油の分子反応モデリング＞
本発明は、ＲＤＳ原料油等の石油の分子組成を入力データとして石油分子の反応解析を行うための定量的な反応モデルを提供しようとするものである。以下、具体例としてＲＤＳの場合を説明する。具体的には図１に示すように、本発明では、各々の反応の中心が十分に離れているという仮定のもと、分子単位ではなく構造属性の単位（コア、側鎖、架橋）毎の反応モデルを作成し、最後に連成して分子反応モデリング、即ち反応速度解析を行うものである。これにより、反応モデルで取り扱う成分数や反応経路数を大幅に減らすことができる。

本発明の石油の分子反応モデリング方法は、図２のフローチャートに示すように、ＦＴ－ＩＣＲ－ＭＳ分析によりＲＤＳ原料油及び生成油に含まれるコアに関する構造解析を行い、コアのランピング及びランピングした反応ネットワークの定義を行い、定義されたネットワークをコンピュータを用いて簡略化し、その結果に基づいて反応速度モデルを構築するとともに反応速度定数を算出するものである。各工程について、以下詳述する。

（１）コアのランピング
先ず、各コアを、基本コア、ヘテログループ、及びナフテン環数の３つの因子により分類する。例えば、ＸＸを基本コア、Ｙをナフテン環数、ＺＺをヘテログループとすると、Ｍ０ＸＸＹＺＺで示されるコードによって認識可能となり、例えば、後記するようにＭＡＴＬＡＢ（登録商標）の計算ソフトを使用する場合にコード情報として使用することができる。

上記した３つの因子のうち、基本コア及びヘテログループは主にコアの反応パスを決める因子である。図３は、コアのランピングの基本的なフローチャートを示したものである。初めに基本コアの構造を比較し、次にヘテログループの孔層を比較する。これらが構造的に類似していれば同じ反応パスを有するとして同じランプに分類する。基本コアが異なる場合、それぞれの基本コアの反応パスを比較し、これらが類似しており、かつヘテログループの構造が類似していれば同じランプに分類する。ヘテログループの構造が異なる場合は、それぞれのヘテログループの反応パスを比較し、類似していれば同じランプに分類し、異なる場合は異なるランプに分類する。この手法に基づき原料のコアをグループ分けし、ランプを作成する。原料のランプに含まれるコアの反応を追跡することにより生成物のランプを作成することができる。原料及び生成物のランプした反応が構築され、これらの反応を組み合わせることでランプによる反応ネットワークを構築することができる。

（２）コアの反応ネットワークの構築
上記のようにしてコアを同一のランプ及び異なるランプにグループ分けし、同一のランプに属するコアの反応を追跡する。ベンゼン環の減少に伴いナフテン環が増加するという仮定に基づくと、ＲＤＳにおける基本コアの反応パスでは、図４に示すように例えば基本コア３８が有する３つのベンゼン環は全て飽和しナフテン環になり、基本コア４２が生成する。

上記したそれぞれの基本コアと反応パスとの関係は、図５に示すように、多くの基本コアは構造が異なっていても類似した反応パスを有していることがわかる。図５においては、便宜的に基本コアを７つのグループに分類し、それぞれのグループをＡ～Ｇとした。基本コアのうち、重要な７つのグループと非重要なコアをまとめると下記の表２のようになる。

次に、ヘテログループの反応パスを図６に示す。０５のヘテロは０２を経て００のヘテロに至る。１１のヘテロは０２のヘテロ又は０７のヘテロを経て００のヘテロに至る。同様に、１０のヘテロも０８のヘテロを経て００のヘテロに至るか、或いは０７のヘテロを経て００のヘテロに至る。さらに、１２のヘテロも０８のヘテロを経て００のへテロに至るか或いは０２のヘテロを経て０５や１１のヘテロと同様に００のヘテロに至る。また０４や０９のヘテロは００のヘテロに至る。類似した反応パスを有するヘテロ毎に分類した結果を表３に示す。

表２及び表３に示した基本コアとヘテログループの組合せにより、コアのランプを定義することができるようになる。さらにこれらの結果と、図５及び図６に示した反応ネットワークとを組み合わせることにより、ランプに基づいた反応ネットワークを構築することができる。上記のようにして構築された反応ネットワークを図７に示す。図７において、矢印はランプ間の反応の関係を示している。なお、図７では、個々のランプの名称をＡＸ、ＡＹ、ＡＶ・・・と表記したが、これは、便宜的に付した基本コアＡ、Ｂ、Ｃ・・・の名称と、ヘテログループＶ、Ｘ、Ｙ・・・とを組み合わせて表記したものである。例えば、図７のランプＦＸは、基本コアＸ及びヘテログループＸから構成されている。また、各ランプに含まれるコアの数は、下記の表４のようになる。

（３）ランプに基づいたコアの反応の重要度判定
上記のとおり、ランプによる反応ネットワークは、ランプした原料コアの反応を組み合わせたものである。しかしながら、原料のランプ数及び種類によっては反応が多数存在し、反応ネットワークはかなり複雑である。実際、ランプの反応ネットワークは２７のコア及び５０の反応から構成され複雑であり、反応速度定数の算出は膨大な計算が必要であり現実的ではないといえる。そこで、本発明においては、各ランプの反応速度モデルを構築する上で、重要ではない反応を除外し、重要である反応を保持することにより、構築した反応ネットワークを簡略化する。このようにして上記（２）のステップにおいて得られたランピングした反応ネットワークに含まれる反応をある程度削減することができる。反応の重要度の判定は以下のようにして行われる。

先ず、ある原料のコアについて、そのコアが含まれる全反応における個々の反応の影響度を下記式（I）及び（II）により算出する。

式（I）において、
Ｃ_ｉ，ｔ：は、反応時間ｔにおいて、反応ｉに存在する原料のモル濃度であり、
ｋ_ｉは、反応ｉの反応速度定数であり、
Ｎは、反応ネットワークに存在する反応の総数、
である。反応時間に応じてＧＬＯ_ｔの値は異なる。各反応のＧＬＯ_ｔを比較し、ＧＬＯ_ｔが最も小さい反応を除外した。下記の条件を満たさなくなるまでこの計算を繰り返す。

式（II）において、
ｋ_ｒは、反応ネットワークから除外された反応速度定数であり、
Ｃ_ｒは、反応ネットワークから除外された原料コアの濃度であり、
αは、０～１である。
式（II）は、反応ネットワークから除外した”重要ではない”反応の数を調節するために行う。αの値を高くするほど、”重要ではない”反応の数を削減することができる。例えば、後記する具体的なモデリングでは、αを０．９９５として、重要ではない反応の数を削減した。式（Ｉ）に示されるように、ＧＬＯ_ｔの値は主に反応速度定数及び原料モル濃度によって決まる値であることがわかる。

次に、ＧＬＯ_ｔによる判定後、残った反応は”重要”な反応と認識され、除外された反応は、下記式（III）により再度評価される。即ち、反応時間ｔ_ｑ及びｔ_ｑ－１における各原料のモル流量の差分を算出し、前述の判定において除外された反応の重要度を再度評価する。反応時間t_q-1におけるコアの濃度をＣ_{ｉ，ｔｑ－１}、反応時間t_qにおけるコアの濃度をＣ_ｉ，ｔｑとして、下記式（III）により除外された反応の重要度を再度評価する。ここでの閾値βは任意の値を用いる。算出されたＭＯＬの値がβよりも大きい場合は、重要なコアと判断する。

このようにして式（III）により重要度の再評価を行い、”重要ではない”コアを除外し、”重要”なコアを保持する。”重要”なコアは原料又は生成物として多くの反応に存在しているが、一部の反応は無視しても問題がないことから、保持されたコアの反応の重要度を下記式（IV）で示されるＬＯＣの値により評価を行う。

なお、上記式（IV）において、
Ｃ_ｉは、コアｉのモル濃度であり、
ｋ_ｉは、原料としてのコアｉの反応速度定数であり、
ｋ’_ｉは、 生成物としてのコアｉの反応速度定数であり、
ｍ及びｎは、原料及び生成物としてのコアｉを含む反応数であり、
χは、任意の値である。

上記（IV）により、あるコアについて、そのコアが含まれる全反応における個々の反応の影響度を算出する。ＬＯＣがχより大きい場合、その反応は”重要”であると判定される。

反応の重要度判定の手法を示すと、図８のようなフロー図になる。上記した式（I）～（III）の計算は、マサチューセッツ州ナティックのＭａｔｈＷｏｒｋｓ社によって開発されたＭＡＴＬＡＢ（登録商標）ソフトウェアによって提供される最適化ツールボックス等を用いて実施することができる。ＭＡＴＬＡＢ（登録商標）に組み込んだ場合、重要な反応又は重要ではない反応は１又は０のコードとして出力される。この手法に基づいて、重要ではないコアは原料及び生成物から除外される。また、重要ではない反応も除外される。

例えば、ＲＤＳには、種々の触媒が使用されており触媒の種類によって存在するコアや反応は変わらないものの、反応速度定数は異なるものと考えられる。したがって、本発明においても、触媒によって、上記した式（I）～（IV）の計算結果も異なるものと考えられる。以下、ＲＤＳに使用されている三種の触媒（ＨＤＭ＿１０、ＨＤＭ＿１１、ＨＤＳ３）のそれぞれについて行った計算について説明する。

図９は、３種の触媒が充填されたＲＤＳ反応器を示したものである。また、下記表５は、各触媒における主な反応条件を示したものである。

これらの触媒において生じる反応は同様であるが、反応速度定数は触媒により異なる。それぞれの触媒における反応は、２１０７の反応パス及び１２３３のコアから構成されている（但し、Ｈ_２、Ｈ_２Ｓ及びＮＨ_３のようなガス生成物は除く）。の条件において生じる主な反応は、核水添反応（ＳＡＴ）、水素化脱窒素反応（ＨＤＮ）及び水素化脱硫反応（ＨＤＳ）であり、シミュレーションにおいて、これら反応の速度パラメーター（頻度因子及び活性化エネルギー）及びコアのモル流量をＪＫＭＴにインプットする。

下記表６は、α、β及びχの値を、それぞれ０．９９５、１Ｅ－３、１Ｅ－３として、式（I）～（IV）の計算をＭＡＴＬＡＢ（登録商標）に組み込んで行ったものである。なお、α、β及びχの値は、反応にほとんど影響しないコア及び反応パスを削減するために設定した。

表６からも明らかなように、ＨＤＭ＿１０において”重要”だと判定され反応ネットワークに維持されたコア又は反応が、ＨＤＭ＿１１及びＨＤＳ＿３では”重要ではない”と判定されることもある。これより、３種の触媒に適用可能なランプによる反応ネットワークを構築する必要があるといえる。図８に示した応の重要度判定の手法を示すフローに従って３種の触媒における反応の重要度を判定した結果を下記の表７に示す。ここでの判定では、１つの触媒でも１となれば、最終判定では”重要”な反応と判定し、ランプした反応ネットワークの構築に使用した。全触媒で０となれば、その反応は除外した。

“重要ではない”反応の除外により、反応ネットワークを２７のコア及び３０の反応に簡略化することができる。さらに、反応先が同一となる類似した反応パスを組み合わせることにより、代わりのランプモデルを構築することができる。例えば、図７に示した反応ネットワークは、簡略化されて、図１０に示すように、２７のコア及び３０の反応から構築したランプ（簡略ランプ１）することができる。また、ランプＡ及びＥを組み合わせることにより簡略ランプ１を更に簡略化して図１１に示すように、２３のコア及び２５の反応から構築したランプ（簡略ランプ２）とすることができる。また、ランプＡ、Ｅ及びＧを１つのランプとすることにより、図１２に示すように１９のコア及び２１の反応から構築したランプ（簡略ランプ３）とすることができる。更に、簡略ランプ３のヘテログループＹ及びＶを組み合わせて簡略化することにより、図１３に示すように、１４のコア及び１５の反応から構築したランプ（簡略ランプ４）とすることもできる。

（４）速度モデルの構築
上記のようにしてランピングした反応ネットワークにより、各ランプの速度変化を表す速度モデルを構築する。反応は全て原料に対して擬一次反応と仮定する。

（５）反応速度定数の算出
生成物の収率及び組成推定のために、ランプに基づいた反応速度定数を推定する。推定にはＭＡＴＬＡＢ（登録商標）にて最小二乗法を用いた算出を行うことができる。各ランプの反応速度Ｌｉは、下記式（Ｖ）で表すことができる。なお、式（Ｖ）において、ｋ_ｌｐ，ｉは、他のランプへの反応における反応速度定数であり、ｋ_ｌｐ，ｊは、ランプからの反応による反応速度定数である。

ランプした反応速度モデルに基づき最小二乗法によりランプの反応速度定数を推定する。最小二乗法による反応速度定数の推定はＭＡＴＬＡＢ（登録商標）のｄｓｏｌｖｅ関数を使用することが好ましい。これより、反応速度定数及び反応時間を関数とした濃度が算出される。最小二乗法によりランプの反応速度定数を推定するために、反応速度定数の初期値を仮設定し、式（V）に代入する。その結果、種々の反応時間ｔにおけるランプの濃度Ｌ_{ｉ＿ｃａｌ，ｔ}が算出される。後記する最小二乗法によるランプの反応速度定数の算出にて、ここで算出した値と参考値（ＪＫＭＴより求めた濃度は）を用いて計算を行う。

最小二乗法によるランプの反応速度定数の算出は、下記式（VI）により行うことができる。なお、式（VI）において、Ｌ_{ｉ＿ｒｅｆ，ｔ}は反応時間ｔにおいて、ＪＫＭＴにより算出されたランプＬ_ｉに含まれる各コアの濃度の総和であり、Ｑはデータ数である。

上記式（ＶＩ）の計算は、ＭＡＴＬＡＢにおいて、Ｓｉｍｐｌｅｘ法により最初の計算の後に新たな反応速度定数が得られ、Ｆが最小になるまで計算を繰り返す。最終的にランプの反応速度定数が求められ、さらに各ランプの濃度が計算される。次いで、ランピングした反応の頻度因子Ａ_ｌｕｍｐ及び活性化エネルギーＥ_ｌｕｍｐを算出するためにアレニウスの式を使用する。下記式（VII）に示すように、アレニウスプロットはｌｎ（ｋ_ｌｕｍｐ）及び１／Ｔで表される線形式である。アレニウスプロットを使用することにより、種々の温度における反応速度定数を算出することができる。従って、種々の温度におけるランプの濃度を算出するための速度モデルを推定することができる。

例えば、図１３に示した簡略ランプ４を用いてランプの反応速度モデルを構築し、３種の触媒における反応挙動の描写及び生成物の収率を予測した。ＭＡＴＬＡＢ（登録商標）によるプログラムの構築において、各ランプの名前を設定した。各反応の反応速度定数とともに図１４に示す。図１４に示したＬ１～Ｌ１４の各反応速度は、以下の式で示される。

上記の式（V-1）～（V-14）をＭＡＴＬＡＢ（登録商標）のｄｓｏｌｖｅ関数を用いて積分すると、下記式（V-1’）～（V-14’）が得られる。

上記記式（V-1’）～（V-14’）で示されるように、各ランプのモル流量は反応時間、反応速度定数及び関連する反応のモル流量の関数である。式（V-1’）、 式（V-4’）、式（V-7’）、式（V-10’）を除き、多くの反応速度定数が式に含まれているため計算が複雑である。これらを解くため、ＭＡＴＬＡＢ（登録商標）にてＩｓｑｎｏｎｌｉｎ関数を使用して計算することができる。

ＪＫＭＴから得られたデータを、ランプした反応速度定数の推定及び構築したモデルの精度を比較するための参考値として用いた。ＪＫＭＴは種々の反応時間ｔにおけるコアのモル流量を得るために用い、これに基づき、種々の反応時間におけるランプのモル流量Ｌｉ＿ｒｅｆ（ｔ）を計算する。ｋｉの初期値を仮定し、Ｌｉ＿ｃａｌ（ｔ）を計算するために式（V-1’）～（V-14’）に導入する。上記式（VI）を用いてＬｉ＿ｒｅｆ（ｔ）とＬｉ＿ｃａｌ（ｔ）の誤差を算出する。Ｓｉｍｐｌｅｘ法によって新たにｋｉが算出され、Ｆが最小となるまで計算が繰り返される。Ｆが最小となった時のｋｉをランプによる反応の反応速度定数とすることができる。

各触媒において、種々の反応温度における反応速度定数を算出する。これを、上記式（VII）に導入して頻度因子Ａ_ｌｕｍｐ及び活性化エネルギーＥ_ｌｕｍｐを算出する。これらを他の温度域におけるランプの反応速度定数の算出に用い、これより構築したモデルの推定精度を確認することができる。このようにして算出した各触媒におけるランプに基づく反応速度定数を下記表に示す。表８は、ＨＤＭ－１０におけるランプに基づく各温度での反応速度定数であり、表９は、ＨＤＭ－１１におけるランプに基づく各温度での反応速度定数であり、表１０は、ＨＤＳ＿３におけるランプに基づく各温度での反応速度定数である。なお、表８中の６５１Ｋの温度箇所はモデル構築において検討した条件であり、表９中の６６４Ｋの温度箇所はモデル構築において検討した条件であり、表１０中の６７９Ｋの温度箇所はモデル構築において検討した条件である。

表８～１０からも明らかなように、ＨＤＭ＿１０のｋ３を除き、反応速度定数は温度の上昇に伴い増加することがわかる。３種の触媒において、ヘテログループ１１及び１２（ｋ１、ｋ５、ｋ９及びｋ１２）で示されるＨＤＮ及びＨＤＳの双方を含む反応は、ヘテログループ０２、０４、０５、０７、０８、０９及び１０（ｋ２、ｋ６及びｋ１０）よりも高い反応速度定数を示すことがわかった。また、核水添反応について、ヘテログループ（ｋ３、ｋ７、ｋ１１及びｋ１４）の反応速度定数は、ヘテロを含まないもの（ｋ４、ｋ８及びｋ１５）よりも高いことがわかる。ＨＤＮ及びＨＤＳ反応の反応速度定数は核水添反応よりも高い（ｋ１及びｋ２＞ｋ３及びｋ４、ｋ５及びｋ＞ｋ７及びｋ８）。３種の触媒において、ＨＤＳ＿３はヘテログループ０２、０４、０５、０７、０８、０９及び１０のＨＤＳ及びＨＤＮ反応の反応速度定数が高いことがわかる。ＨＤＭ＿１０はヘテログループ１１及び１２の反応に関して高い反応速度定数を示している。

また、図１５～図１７にモル流量の参考値及び推定値の比較結果を示す。広い温度域において、ＨＤＭ＿１１及びＨＤＳ＿３では良好な結果を示していた（Ｒ２＞０．９９９）。ＨＤＭ＿１０は、ＨＤＭ＿１１及びＨＤＳ＿３よりも精度が低かった（Ｒ２の平均＞０．９９）。ＨＤＭ＿１１及びＨＤＳ＿３では、モデル構築に用いた温度より高い温度では（ＨＤＭ１１＞６６４Ｋ、ＨＤＳ３＞６７９Ｋ）、精度が僅かに減少していることがわかる。ＨＤＭ＿１０では、６５１Ｋより高い温度において精度が急激に減少していた。

種々の反応時間において、構築したモデルの精度を検証した。図１８～図２２に反応器の反応方向に沿ったモル流量の変化を示す。ＨＤＭ＿１１及びＨＤＳ＿３において、反応器の各位置において精度は良好である。ＲＤＳの主な目的は油中の硫黄及び窒素を低減することであり、ランプのプロファイルに基づき最適な運転条件及び触媒を提案することが可能である。例えば、高温域においては、窒素及び硫黄を含むＬ１、Ｌ２、Ｌ４、Ｌ５のランプのモル流量は反応方向に沿って大きく減少していることがわかる。このプロファイルに基づき、最適な温度及び反応時間を提案することができる。ＨＤＭ＿１０において、ランプＬ８及びＬ１３では参考値と推定値の差異が極めて大きい。ランプによる算出値は参考値と比較して平均で２０及び１０％低い。これらの差異は高温域（＞６５１Ｋ）において極めて大きくなった。

＜石油の分子反応モデリング装置及び該装置を実行させるためのコンピュータプログラム＞
次に、本発明の石油の分子反応モデリング装置の実施形態を説明する。コンピュータに本発明のプログラムを実行させることにより、コンピュータが石油の分子反応モデリング装置として機能する。

本装置は、演算装置を備えている。演算装置は、１つのＣＰＵで構成してもよいし、通信回線を介して互いに接続された複数の演算装置で構成されてもよい。また、本装置は、記憶部を備えていてもよく、その場合、演算装置に内蔵されていてもよいし、外部装置であってもよいし、通信回線を介して接続された記憶装置であってもよい。

本演算装置は、モデリング成分に含まれる基本コア、ヘテログループ及びナフテン環数から成るコア因子情報を取得する成分情報取得部と、成分情報に基づき、基本コアの構造を比較した後に、ヘテログループの構造を比較することにより、前記構造に基づくランプへ、コアのランピングを行うランプ部と、ランプに含まれるコアの反応を追跡することにより、ランプに基づいたコアの反応ネットワークを構築する、反応ネットワーク構築部と、ランプに基づいたコアの反応の重要度を判定する重要度判定部と、反応ネットワークから、前記各ランプの速度変化を表す速度モデルを構築する、速度モデル構築部と各ランプの速度モデルに基づく反応速度定数を最小二乗法を用いて算出する、反応速度定数算出部とを備えている。

本発明において、上記したコアのランピング及びランピングした反応ネットワークの定義、ランピングした反応ネットワークの簡略化、並びに反応速度モデルの構築及び反応速度定数の算出は、ハードウェア又はソフトウェア、又はこれらを複合した構成によって実行することができる。ソフトウェアによる処理を実行する場合には、処理シーケンスを記録したプログラムを、専用のハードウェアに組み込まれたコンピュータ内のメモリにインストールして実行させるか、各種処理が実行可能な汎用コンピュータにプログラムをインストールして実行させることができる。

例えば、プログラムは、記録媒体としてのハードディスクやＲＯＭに予め記録しておくことができる。また、プログラムは、フレキシブルディスク、ＣＤ－ＲＯＭ、ＭＯディスク、ＤＶＤ、磁気ディスク、半導体メモリなどのリムーバブル記録媒体に、一時的又は永続的に格納（記録）しておくことができる。

なお、プログラムは、上述したようなリムーバブル記録媒体からコンピュータにインストールする他に、ダウンロードサイトから、コンピュータに無線転送したり、ＬＡＮ、インターネットといったネットワークを介して、コンピュータに有線で転送したりでき、コンピュータでは、そのようにして転送されてくるプログラムを受信し、内蔵するハードディスクなどの記録媒体にインストールすることができる。

また、本明細書に記載された各種の処理は、記載に従って時系列に実行されるだけではなく、処理を実行する装置の処理能力や必要に応じて並列的に又は個別に実行されてもよい。また、本明細書において、システムとは、複数の装置の論理的集合構成であり、各構成の装置が同一筐体内にあるものに限定されるものではない。

以上、本発明の実施形態を説明したが、本発明は、上述した実施形態に限定されるものではなく、本発明の範囲で種々の変更実施が可能である。

本発明によれば、コアのランピング及びランピングした反応ネットワークを定義し、コンピュータを用いて定義された反応ネットワークを簡略化し、それらの結果に基づいて分子の反応を描写し、構造に基づき生成物を推定する速度モデルを効率よく構築する新たな手法が提供され、石油精製設備の運転の安定性及び運転効率を飛躍的に向上させることに寄与するものである。

Claims (10)

1. コンピュータによる、石油の分子反応モデリング方法であって、
   （１）基本コア、ヘテログループ及びナフテン環数から成るコア因子に基づいて、基本コアの構造を比較した後に、ヘテログループの構造を比較することにより、前記構造に基づくランプへ、コアのランピングを行うステップと、
   （２）前記ステップ（１）で得られたランプに含まれるコアの反応を追跡することにより、前記ランプに基づいたコアの反応ネットワークを構築するステップと、
   （３）前記ステップ（２）で得られたランプに基づいたコアの反応の重要度の判定を、
   （Ａ）分子反応モデリングに供する原料のコアについて、該コアが含まれる全反応における個々の反応の影響度を、下式（I）：
   （式中、
   Ｃ_ｉ，ｔ：は、反応時間ｔにおいて、反応ｉに存在する原料のモル濃度であり、
   ｋ_ｉは、反応ｉの反応速度定数であり、
   Ｎは、反応ネットワークに存在する反応の総数、である）
   により算出し、各反応のＧＬＯ_ｔを比較し、ＧＬＯ_ｔが最も小さい反応を除外し、
   （Ｂ）前記ステップ（Ａ）の影響度の算出を、下式（II）：
   （式中、
   ｋ_ｒは、反応ネットワークから除外された反応速度定数であり、
   Ｃ_ｒは、反応ネットワークから除外された原料コアの濃度であり、
   αは、０～１である）
   により表される条件を満たさなくなるまで、繰り返し、
   （Ｃ）反応時間ｔ_ｑ及びｔ_ｑ－１における各原料のモル流量の差分を算出し、前記ステップ（Ａ）において除外された反応の重要度を、下式（III）：
   （式中、
   Ｃ_{ｉ，ｔｑ－１}は、反応時間t_q-1において、反応ｉに存在する原料のモル濃度であり、
   Ｃ_ｉ，ｔｑは、反応時間t_qにおいて、反応ｉに存在する原料のモル濃度であり、
   βは、任意の値である）
   により算出したＭＯＬがβよりも大きい場合には、重要なコアとして保持し、
   （Ｄ）前記ステップ（Ｃ）において得られた重要なコアを、下式（IV）：
   （式中、
   Ｃ_ｉは、コアｉのモル濃度であり、
   ｋ_ｉは、原料としてのコアｉの反応速度定数であり、
   ｋ’_ｉは、 生成物としてのコアｉの反応速度定数であり、
   ｍ及びｎは、原料及び生成物としてのコアｉを含む反応数
   χは、任意の値である）
   により算出したＬＯＣがχよりも大きい場合には、前記コアの含まれる反応は、重要であると判定すること、
   により行う、ステップと、
   （４）前記ステップ（２）で得られた反応ネットワークから、前記各ランプの速度変化を表す速度モデルを構築するステップと、
   （５）前記ステップ（４）で得られた各ランプの速度モデルにもとづいた反応速度定数を、最小二乗法を用いて算出するステップと、
   を含んでなることを特徴とする、石油の反応モデリング方法。
2. 石油が重質油である、請求項１に記載の方法。
3. 前記石油の分子反応が常圧残油水素化脱硫（ＲＤＳ）である、請求項１又は２に記載の方法。
4. 前記ステップ（５）における最小二乗法が、ｄｓｏｌｖｅ関数を用いて行われる、請求項１～３のいずれか一項に記載の方法。
5. 前記ステップ（５）における速度モデルの反応速度は、下式（V）：
   （式中、
   ｋ_ｌｐ，ｉは、他のランプへの反応における反応速度定数であり、
   ｋ_ｌｐ，ｊは、ランプからの反応による反応速度定数である）
   により算出される、請求項１～４のいずれか一項に記載の方法。
6. 前記式（Ｖ）をｄｓｏｌｖｅ関数を用いて積分することにより、各ランプの反応速度定数及び反応時間を関数とした濃度が算出され、
   前記反応速度定数の算出において、反応速度定数の初期値を仮設定し、前記式（V）に代入することにより、反応時間ｔにおけるランプ濃度Ｌ_{ｉ＿ｃａｌ，ｔ}が算出され、
   前記算出したランプ濃度Ｌ_{ｉ＿ｃａｌ，ｔ}から算出したランプのモル流量（Ｌ_{ｉ＿ｃａｌ，ｔ}）と、ＪＫＭＴにより算出した参考値（Ｌ_{ｉ＿ｒｅｆ，ｔ}）とを用いて、下式（VI）：
   （式中、
   Ｌ_{ｉ＿ｒｅｆ，ｔ}は反応時間ｔにおいて、ＪＫＭＴにより算出されたランプＬ_ｉに含まれる各コアの濃度の総和であり、
   Ｑはデータ数である）
   により算出されたＦを、Ｓｉｍｐｌｅｘ法によりＦ値が最小になるまで計算を繰り返し、
   前記Ｆ値が最小になるまで計算を繰り返すことにより、前記ランプの反応速度定数及び各ランプの濃度が算出され、
   前記ランプした反応の頻度因子Ａ_ｌｕｍｐ及び活性化エネルギーを算出するためにアレニウスの式（VII）：
   を使用し、該式により示されるアレニウスプロットをしようすることにより、種々の温度における反応速度定数を算出することにより、反応速度モデルを推定する、
   請求項５に記載の方法。
7. 分子反応モデリングは、ＪＫＭＴにより行われる、請求項１～６のいずれか一項に記載の方法。
8. 前記αは０．９９５であり、前記βは１Ｅ^－３であり、前記χは１Ｅ^－３である、請求項１～７のいずれか一項に記載の方法。
9. 石油の分子反応モデリング装置であって、
   モデリング成分に含まれる基本コア、ヘテログループ及びナフテン環数から成るコア因子情報を取得する成分情報取得部と、
   取得した成分情報に基づき、基本コアの構造を比較した後に、ヘテログループの構造を比較することにより、前記構造に基づくランプへ、コアのランピングを行うランプ部と、
   前記ランプに含まれるコアの反応を追跡することにより、前記ランプに基づいたコアの反応ネットワークを構築する、反応ネットワーク構築部と、
   前記ランプに基づいたコアの反応の重要度を、
   （Ａ）分子反応モデリングに供する原料のコアについて、該コアが含まれる全反応における個々の反応の影響度を、下式（I）：
   （式中、
   Ｃ_ｉ，ｔ：は、反応時間ｔにおいて、反応ｉに存在する原料のモル濃度であり、
   ｋ_ｉは、反応ｉの反応速度定数であり、
   Ｎは、反応ネットワークに存在する反応の総数、である）
   により算出し、各反応のＧＬＯ_ｔを比較し、ＧＬＯ_ｔが最も小さい反応を除外し、
   （Ｂ）前記ステップ（Ａ）の影響度の算出を、下式（II）：
   （式中、
   ｋ_ｒは、反応ネットワークから除外された反応速度定数であり、
   Ｃ_ｒは、反応ネットワークから除外された原料コアの濃度であり、
   αは、０～１である）
   により表される条件を満たさなくなるまで、繰り返し、
   （Ｃ）反応時間ｔ_ｑ及びｔ_ｑ－１における各原料のモル流量の差分を算出し、前記ステップ（Ａ）において除外された反応の重要度を、下式（III）：
   （式中、
   Ｃ_{ｉ，ｔｑ－１}は、反応時間t_q-1において、反応ｉに存在する原料のモル濃度であり、
   Ｃ_ｉ，ｔｑは、反応時間t_qにおいて、反応ｉに存在する原料のモル濃度であり、
   βは、任意の値である）
   により算出したＭＯＬがβよりも大きい場合には、重要なコアとして保持し、
   （Ｄ）前記ステップ（Ｃ）において得られた重要なコアを、下式（IV）：
   （式中、
   Ｃ_ｉは、コアｉのモル濃度であり、
   ｋ_ｉは、原料としてのコアｉの反応速度定数であり、
   ｋ’_ｉは、生成物としてのコアｉの反応速度定数であり、
   ｍ及びｎは、原料及び生成物としてのコアｉを含む反応数
   χは、任意の値である）
   により算出したＬＯＣがχよりも大きい場合には、前記コアの含まれる反応は、重要であると判定すること、
   により行う、重要度判定部と、
   前記反応ネットワークから、前記各ランプの速度変化を表す速度モデルを構築する、速度モデル構築部と、
   前記各ランプの速度モデルに基づく反応速度定数を最小二乗法を用いて算出する、反応速度定数算出部と、
   を備えることを特徴とする、石油の分子反応モデリング装置。
10. 請求項１～８のいずれか一項に記載の方法をコンピュータに実行させるためのコンピュータプログラム。