生産シミュレータの開発
専用シミュレータの開発
2003年 5月 20日
研究実施機関：東京大学
日本オイルエンジニアリング株式会社

専用シミュレータの開発体制

東大，日本オイルエンジニ
アリング，産総研の共同研
究1フィールドスケールシミュレータによる予測
専用シミュレータ（ 2次元円筒・3次元直交座標）

現行シミュレータ
機能強化（2次元円筒座標）

フルフィード
定温閉境界

ガス

定温閉境界

減圧生産

コア実験

メタンハイドレート

水平坑井

定圧境界

ハイドレート層

熱単一坑井テスト

上部層

ガス生産井

熱
ハイドレート層
フリーガス層シミュレーションの例－ハイドレート層の直下のガス層からの減圧法によるガス生産# 1 free gas model/day)8000# 8 gas only7000

gas produciton rate (m9000

310000

6000
5000
4000

ハイドレート層＋ガス層のみ

フリーガス層

3000
2000
1000
002468101214161820time (years)

2シミュレーションの例－ハイドレート分解挙動の予測
垂直生産井

1年後２年後３年後平成14年度の研究項目
1. シミュレータ相平衡計算機能の強化
2. メタンハイドレート分解挙動予測モデルの改良

3

1. シミュレータ相平衡計算機能の強化

研究目的・内容及び必要性• 研究目的
6 専用シミュレータをフィールド産出試験の予測に適用させるように，相平衡計算プログラムの改良と機能追加を行う。

• 研究内容¾ 相平衡計算プログラムの改良作業
 インヒビター効果の予測機能の追加
 氷の生成に関する計算機能の追加
 エンタルピー計算機能の改良¾ 改良プログラムによる計算結果のシミュレータへの組み込み¾ 欧米研究機関訪問調査（JOEと共同）• 研究の必要性••

減圧法，熱刺激法，インヒビター注入法等の各種のフィールド産出試験を専用シミュレータによる計算で検討する場合に，各相の物性を適切に計算する必要がある。
専用シミュレータにおける各相の物性計算を，本研究の相平衡計算プログ
ラムが担う。

4

研究結果（１）
インヒビター効果の予測機能の追加
CH3OH，NaCl等のインヒビターが存在する場合の3相平衡圧力，各相平衡状態における相物性（相密度，インヒビター濃度，粘度等）を計算出力できるようにコードを改良した。
メタン―水系にメタノールを添加する場合の３相平衡圧力曲線の計算例

専用シミュレータへは簡便式を引渡し

研究結果（２）
インヒビター効果の予測機能の追加
メタン－水－メタノールの３成分系において，気相中のメタノールの微量成分濃度を計算するための気液平衡比KMeOHを与えた。

KMeOH = yMeOH / xMeOH

圧力P，温度T，液相組成xMeOHを入力して，KMeOHを計算

専用シミュレータへ関数として引渡し

5研究結果（３）氷の生成に関する計算機能の追加
シミュレータの計算過程では，ハイドレート分解により地層の温度が０℃以下に低下する場合がある。この場合の挙動予測に対応するために，氷が生成する場合の相平衡計算機能を追加した。
改良したプログラムによるI-V-Hの３相平衡線の予測

メタン―水系のP-T図

エタン―水系のP-T図•４重点温度をインヒビター濃度の関数として表して，簡便的な予測に対応•氷の生成速度をどのように表現するかが今後の課題

まとめ• 成果の総括¾ インヒビター効果の予測機能の追加
 CH3OH，NaCl等のインヒビターが存在する場合の3相平衡圧力，各相平衡状態における相物性（相密度，インヒビター濃度，粘度等）を計算出力できるようにコードを改良した。
 CH4-CH3OH-H2Oの３成分系で，メタノールの気液平衡比を与える関数を作成した。

¾ 氷の生成に関する計算機能の追加

 シミュレータの計算過程では，ハイドレート分解により地層の温度が０℃以下に低下する場合がある。この場合の挙動予測に対応するために，氷が生成する場合の相平衡計算機能を追加した。

¾ 以上の機能強化したプログラムで計算出力される各種の物性の推算式を専用シミュレータに組み込んだ。

• 今後の対応¾ 平成１５年度に相平衡計算プログラムを完成させて，その結果を専用シミュレー
タに組み込む。改良作業は，以下の通り。
 エンタルピー計算に関するコード改良とMH分解時のエネルギー保存の検証
 水相にイオン種が異なる濃度で含まれている系の計算に対応したコード改良
 相平衡計算プログラムを用いた各相物性の計算結果のまとめ（精度の把握） 計算結果を専用シミュレータへ組み込むにあたって，シミュレータの計算時間を短縮できるような各物性の推算式の導出

6

2. メタンハイドレート分解挙動予測モデルの改良研究目的・内容及び必要性• 研究目的
6 MH分解･流動予測の基礎となる物理現象のより正確な記述，数値誤差の減少，計算効率の改善により，プロトタイプシミュレータ機能の向上を図る。

• 研究内容¾ プロトタイプシミュレータの改良（要素モジュールの改良，メタノール成分の追加，計算速度の向上，貯留層パラメータのスケール依存性の検討）¾ 欧米研究機関訪問調査¾ プログラムの実践的使用に向けた改良（入出力画面の整備，マニュアルの整備）• 研究の必要性••

詳細なMH分解･採収機構の把握，最適なMH分解･採収法の模索のためには，数値シミュレーションは不可欠である。
より厳密にMH分解･採収機構を記述すると共に，より実践的かつ広範囲の分解･採収法の適用によるMH挙動予測に対応することは重要である。

7

（１）プロトタイプシミュレータの改良① 要素モジュールの改良（1/2）
 絶対浸透率（毛細管束，Kozeny，NBNLモデル）
k (Sh ) = k 0 (1 − Sh )N

従来モデル

1

有効浸透率比, k(Sh)/k0Kozenyモデル，
MHが孔隙中心から生成(1 − Sh )n + 2 k
k (Sh ) =02(1 + S )0.5hLBNLモデル
 1 − Sh − Sr k (Sh ) =  1 − Sr 0.5 1    1 − Sh − Sr  m 
1 − 1 −  1 − S r m

0.10.01
毛細管束モデル （3.1.2.2-7)式
毛細管束モデル （3.1.2.2-11）式
Kozeny粒子モデル (3.1.2.2-20）式
0.001

Kozeny粒子モデル （3.1.2.2-22）式

2LBNLモデル （3.1.2.2-23）式

 k0
0.0001
00.2

0.4

0.6

0.8

1ハイドレート飽和率（１）プロトタイプシミュレータの改良① 要素モジュールの改良（2/2）
相対浸透率Sgr

S h ≥ αS wc

S w* = S wc + S w(1 − S wc )Sw + S g

S h < αS wc

S w* = S wc + (S w − A)

(1 − S wc )(S + S − A)w

A = S wc −ShSH

Sgr

αSwc

S h ≤ 1 − S wc

SH-αSwc

S w* = S wc + (S w − S wc )wSgg“1-Swc”

α

不動水が優先的にハイドレートになると仮定

wc

Sg

g

S h > 1 − S wc

S w* = S wc −SwSH(1 − S wc )(S − S + S )( S h + S wc − 1)α

SwSw-Swc“1-Swc”

Swc

不動水はハイドレートにならないと仮定熱導率
K = (1 − φ)K r + φSg K g + φSw K w + φSh K hφS

K = K1r−φ K g g K φhSh K φwSw

従来モデル
３つの代替モデル（例）8（１）プロトタイプシミュレータの改良② メタノール成分の追加（1/2）
メタノールを独立した成分として取り扱う3成分系モデルへと変更メタノール平衡定数の導入
流動方程式に成分モル分率，相成分モル分率を導入

k Dx k rg
k k∇ ⋅ ρ g 1 − y CH 4
∇Φ g  + 1 − y CH 4 q̂ g + ∇ ⋅  ρ w x CH 3OH Dw rw ∇Φ w  + x CH 3OH q̂ wµgµw
=)[((()∂φ ρ w Sw 1 − y CH 4 + ρ gSg x CH 3OH∂t))]200

気・液フラッシュ計算ルーチンを付加
F(G ) =

3() − z H O + z CH = 0z CH 3OH K CH 3OH − 1∑ (yl − x l ) = 1 + G (K CH OH − 1)l =13241− GG

三相平衡圧力, MPa(150100
メタノール濃度0%
メタノール濃度X%

500270三相平衡圧力曲線を補正

280290300310320

温度, KJacobian行列要素計算ルーチン，線形方程式解法ルーチンを修正
流動方程式混合則を適用してメタノールを含む気･液相特性を推定（１）プロトタイプシミュレータの改良
Injector

Injector② メタノール成分の追加（2/2）
Producer

Injector

ガス相中のメタノール濃度予測結果
Injector

Injector

Producer

Cumulative Gas Production, Sm320000

Injector

15000

10000

メタノール無し（今年度）
メタノール5%（今年度）
メタノール無し（昨年度）
メタノール5%（昨年度）

5000

00

Injector

Injector

5101520Time, day

メタノール圧入に伴う生産ガス量の予測結果
水相中のメタノール濃度予測結果

9

（１）プロトタイプシミュレータの改良③ 計算速度の向上
配列の各要素を3重ポインターを使用した3次元配列から，1次元配列に変更（アクセス速度の向上）
 複数のIndexの使用によるJacobian行列構成ルーチンからの判定文の削除


計算速度の10％向上

例： 流出・入量計算時のIndex

7 7 8 8 9 97894564564561 1 2 2 3 3312i_ind

j_ind111212323441542653774875986z 従来if(block_num%3 == 1)
Qim = 0else
Qim = Tm*(p[i]-p[i-1])
z 改良後

Qim = Tm*(p[i]-p[i_ind[i]])（１）プロトタイプシミュレータの改良④ 分解速度定数のスケール依存性の検討（1/5）
貯留層パラメータの一部をグリッドのスケールに応じて変化させること（upscaling）により，ある程度サイズの大きなグリッドブロックを用いた場合にも，微小グリッドを用いた計算と同じ挙動を予測することを試みた。
分解機構
圧力減退

温度上昇

主な計算結果影響因子
初期浸透率（N, ハイドレート飽和率）
平衡圧力との乖離度
z 分解速度定数z ブロックサイズ
z 1次元･2次元モデル

アップスケール
対象パラメータＮ
平衡圧力
z 分解速度定数zzzz

初期浸透率（N, ハイドレート飽和率）
平衡温度との乖離
z ブロックサイズ
z 1次元モデルzz

Ｎz10

ハイドレート存在グリッドのアップスケール（1次元モデル）
閉境界

定圧境界

閉境界

定圧境界

微細グリッド

粗グリッド
810000

圧力, MPa100

分解可能ブロック（大）

圧力
浸透率

41

浸透率回復
（大）

20.01

6圧力, MPa

平衡圧力610000

平衡圧力

FLOW
浸透率比, Kd/Kd0FLOW

4100
圧力
浸透率

分解可能ブロック（小）120.01

浸透率比, Kd/Kd08

浸透率回復（小）000.0001

024コア長さ, m0.0001
062コア長さ, m

46Cumulative Gas Production, Sm31600

10分割
100分割
1000分割5000分割

1200800

4000036Time, day

912（１）プロトタイプシミュレータの改良④ 分解速度定数のスケール依存性の検討（2/5）
Cumulative Gas Production, Sm3250

N値の調整による粗グリッド
モデル挙動予測結果と微細グリッドモデル挙動予測結果のマッチング例

200

150

cum(Base_Case N8.0)
cum(N6.5)
cum(N6.7)
cum(N6.9)cum(N7.1)
cum(N7.3)
cum(N7.5)
cum(N8.0)

100

15ΔP/Pe = 100%（Case A）

012Time, day

34

N*0

グリッドサイズと適用すべきN値の関係（ N値のアップスケーリング）例

12

浸透率変化比, (1-Sh) /(1-Sh)

N

501.00m（10分割）
0.20m（50分割）
0.10m（100分割）

90.02m（500分割）
0.01m（1000分割）
6300.001

0.010.1

110100

初期浸透率, md11（１）プロトタイプシミュレータの改良④ 分解速度定数のスケール依存性の検討（3/5）
3000.10m0.10m (Pe + 1MPa)
0.01m

100

Cumulative Gas Production, Sm32005000分割50分割50分割（分解速度定数1000倍）

150

100

50

00003060Time, days

9010203040Time, day120

分解速度定数の調整による粗グリッドモデル挙動予測結果と微細グリッドモデル挙動予測結果

平衡圧力の調整による粗グ
リッドモデル挙動予測結果と微細グリッドモデル挙動予測結果のマッチング例（１）プロトタイプシミュレータの改良④ 分解速度定数のスケール依存性の検討（4/5）
700

閉境界・断熱2次元モデル減圧

坑井

Cumulative Gas Production, Sm3

Cumulative Gas Production, Sm

3200600

500

400
201×201（N6.00）
21×21(N6.00)
21×21（N5.31)
300

200

10000

2000000

4000000

6000000

8000000

10000000

Time, day

2次元モデル

N値の調整による粗グリッド
モデル挙動予測結果と微細グリッドモデル挙動予測結果

12

（１）プロトタイプシミュレータの改良④ 分解速度定数のスケール依存性の検討（5/5）
30

20

N*

浸透率変化比, (1-Sh) /(1-Sh)

N

ΔT=40K（Case C）

50分割
100分割
500分割
1000分割

10

00.000000001

0.0000001

0.00001

0.001

0.1

10

初期浸透率, md

N値の調整による粗グリッドモデル挙動予測結果と微細グリッド
モデル挙動予測結果 （温度支配ケース）（１）プロトタイプシミュレータの改良⑤ シミュレータ機能追加に伴う修正
幾つかの機能を付加
時間と共に変化する境界条件（圧力・温度）を設定可能→コア実験のより厳密なシミュレーション
氷点以下の温度におけるハイドレートの分解による吸熱及び浸透率を擬似的に計算できるようにした

13

（２）欧米研究機関訪問調査① シミュレータ開発現状調査（1/4）
訪問先訪問日

主な調査内容ClarksonUniversity

2003 年 3 月 17 日 MH 数値シミュレータ
MH 分解実験
油層工学関連実験
The University of 2003 年 3 月 19 日 MH 数値シミュレータ
Calgary
分解速度論
氷･MH 生成問題
Computer
2003 年 3 月 20 日 流体･岩石特性
Modelling Group
数値計算法
Ltd.
市販シミュレータを利用した MH シミュレーション
APA Petroleum
2003 年 3 月 20 日 MH フィールド産出試験
Engineering（２）欧米研究機関訪問調査① シミュレータ開発現状調査（2/4）
Clarkson University, Department of Mechanical and AeronauticalEngineering

Ahmadi教授の研究室はDOEの資金を得て，学生2人によりハイドレート関連の実験及び数値シミュレータの開発を実施。



実験では，長さ10 in，幅3 in，高さ2.5 inの直方体のアクリル容器の中にプロパンハイドレートを作成し，その1端を減圧することでこれを分解した際の分解フロントの動きを温度・圧力センサーにより検知。

数値シミュレータに関しては，現時点ではハイドレート相とガス相の2相に対する単純な温度・圧力の移流・拡散方程式を解くモデルを完成。現在は水相を含めた3相・2成分の，つまり日本側モデルと同じ機能を付加すべく，シミュレータを改良中。



今後予定している多孔質媒体中のハイドレートの分解実験，3次元・3相・多成分のシミュレータの開発を含め，実践的に利用可能なデータ・結果を提供するレベルに達するには，まだかなりの時間を要する。

14

（２）欧米研究機関訪問調査① シミュレータ開発現状調査（3/4）
University of Calgary, Department of Chemical & PetroleumEngineeringハイドレート分解速度論のKim-Bishnoiの式で有名なBishnoi教授が主としてハイドレートの分解実験を，PooladiDarvish助教授がハイドレートの分解・生産に関する数値シミュレータの作成を担当している。



日本側シミュレータを参考にして2次元円筒座標系のシミュレータを完成し，テストランによる日本側シミュレータの計算結果との比較，ケースランを行っている。日本側が計画している機能強化，特に氷とハイドレートの生成について意見を交換。



数値シミュレーションの大家であり，学科長であるSettari教授は不在で，地層変形シミュレーションに対する提言は得られなかった。



ハイドレートの分解実験，それに伴う物理現象に対する知見は深い。数値シミュレータに関しては，3～4年前の日本側シミュレータのコピーを作成したところで，現時点では特に注目すべき研究は行っていないが，今後の交流・情報交換を続ける意義は大きい。

（２）欧米研究機関訪問調査① シミュレータ開発現状調査（4/4）
Computer Modelling Group Ltd. (CMG)

CMG社は石油・ガスの貯留層シミュレータを開発・販売。ブラックオイル型のIMEX，成分型のGEM，熱を取扱うSTARSがその主力製品。これらのシミュレータは国際的に評価が高く，CMG技術者はこれらのシミュレータの開発を通して，熱力学特性を含めた流体特性，岩石特性，数値計算手法に多大な知見を有している。

潜熱を含む流体のエンタルピー，メタノールと水の混合物の特性，固相の存在下での相対浸透率，グリッドサイズに依存する問題の解法について主として議論。いずれも石油・ガス層で頻繁に取扱う問題ではない→すぐに明快な回答を出すことは不可能。

CMG社ではEOSによる流体混合則の適用やグリッドの局所分割等の技術を活かし，必要であれば更なる調査を実施することに対して積極的に参加したい旨を伝えていた。

STARSは熱，化学反応を取扱える為，ハイドレート特性を石油の特性として取扱うことで，ある程度ハイドレートの分解挙動を予測可能と考えられる。CMG社は日本側シミュレータの計算結果とSTARSによる（簡易）計算結果の比較も希望している。

15

（２）欧米研究機関訪問調査② 課題の抽出 計算速度の向上
 不確実要素の解明（ハイドレートの生成過程，氷の生成過程，等） 岩石特性モデルの確立（絶対浸透率，相対浸透率，熱伝導率，等） 相特性モデルの確立（密度，粘度，エンタルピー，内部エネルギー，等）
 坑井内流動機能の導入（圧力，温度） 地層変形計算機能の導入（３）プログラムの実践的使用に向けた改良 ① 入出力画面の整備（1/2）
MH Navigator
Main Menu

MENU BAR

Pull DownMenu

1
Data Entry

2
Run

3
Output View

11
Data Create

21
Execute31
Report

12
Save

22
View

32
Picture14
Data Revise

23
Kill Process

15
Exit

4
Help41
Contents

321
Line Graph

42
About

322
2D View

323
3D View

グラフィック入出力
メイン画面メニュー

Notes
1

FutureExpansion

16

（３）プログラムの実践的使用に向けた改良 ① 入出力画面の整備（2/2）

グラフィック出力例（ライングラフ）

グラフィック出力例（パラメータ分布図）（３）プログラムの実践的使用に向けた改良 ② マニュアルの整備プログラムの入出力
データの解説，理論的背景の説明等を記載した（英文）オンラインマニュアルを作成

17

まとめ• 成果の総括¾ 要素モジュール改良，メタノール成分の追加により， より複雑なMHの分解･流動機構の記述および広範囲の分解･採収法の予測への対応を可能とした。
¾ 行列計算ルーチンを改善することにより，計算速度を約10%向上させた。
¾ 貯留層パラメータの一部をグリッドのスケールに応じて変化させることにより，サイズの大きなグリッドブロックを用いた場合にも，微小グリッドを用いた計算と同じ挙動予測結果が得られる可能性があることを示した。
¾ 境界条件を時間と共に変化させられるように改良し，コア実験のより厳密なシミュレーションを可能とした。
¾ 入出力画面の整備，マニュアルの整備により，シミュレータのより広範囲な利用が期待できるようになった。

• 今後の対応¾ 計算誤差の減少，計算速度の向上に向けたシミュレータの改良は恒常的に継続する。
¾ より複雑なMHの分解･流動機構の記述を可能とすべく，更なるシミュレータの改良を行う

18