石のシンク

大成建設技術センター報
第38号(2005)
地下水－岩石反応を考慮した水理・物質移行解析
−解析コードTOUGHREACTとCO2地下貯留への適用−
山本 肇*1
Keywords : numerical simulation, groundwater, geochemical calculation, mass transport, CO2 sequestration
数値解析，地下水，地球化学計算，物質移行，CO2地下貯留
1. はじめに
うことができるのが大きな特徴である。また,
移流・
拡散解析と地球化学計算をカップリングし、多相（水相、
高レベル放射性廃棄物地層処分や炭酸ガスの地下貯留
においては、数万年以上の超長期にわたる貯蔵物質の移
気相、固相）間の化学反応による物質移動を考慮した移
行解析が可能である。
行特性を評価することが求められている。そのような超
物理・化学条件（温度、圧力、飽和度、イオン強度な
長期の物質移行を実験的に検討することは困難であり、
ど）は広い範囲に対応しており、温度：0℃～300℃,た
数値解析によるアプローチが有効と考えられる。
だし熱力学データベースに依存,
圧 力 ： 1bar ～ 数
この数値解析においては、地下水流動に伴う移流や分
100bar、イオン強度：塩水で最大 6 mol/kg H2O である。
散などの物理的機構とともに、物質の岩石表面への収着
化学種は気相、液相、固相に渡って存在でき、使用す
や鉱物生成など、地下水－岩石間における化学反応を同
る熱力学データベースに登録されていれば化学種数の制
時に考慮する必要がある。従来は地下水流動に伴う物質
限は無い。現在は、EQ3/6V7.2b2)の熱力学データベース
移行と化学反応は別途に評価されてきたが、最近は、計
に追加・改良を加えたデータベースを用いている。化学
算機能力の向上により、化学反応と物質移行の連成解析
反応モデルとしては、水溶性錯体、ガスの溶解／脱ガス、
コードが実用化されている。１次元問題であれば
イオン交換などの平衡反応の他、鉱物溶解／沈殿につい
1)
2)
PHREEQC やEQ3/6 などの標準的な地球化学計算コー
ては反応速度も考慮できる。その他、放射崩壊や線形収
3)
着モデル、表面錯体モデル（拡散二重層モデル）もサポ
ドで可能となっており, ３次元問題でもTOUGHREACT ,
5)
6)
7)
8)
NUFT-C 、MIN3P , PHAST 、CHEMTOUGH などの解
ートされている。
TOUGHREACT の特徴は次の通りである。
析コードが近年開発されている。
本報では、当社がLBNL（ローレンスバークレー国立
1) ３次元における熱－流体－化学の連成解析が可能で
ある。
研究所）から導入した化学反応を考慮できる水理・物質
移行解析コードTOUGHREACTの概要と、炭酸ガスの地
2) 水、ガスなどの複数の異なる流体からなる多相流体
11）
下貯留シミュレーションへの適用事例
を取り扱うことができる。
を紹介する。
3) 積分差分法（IFDM）を採用しているため、解析グ
2. 解析コードTOUGHREACT
リッド形状の自由度が高く、複雑な地形・地質構造
モデルや、亀裂・断層による不均質場の取り扱いが
3）
TOUGHREACT（Xu and Pruess ）は、非等温・多成
可能である。
分・多相系流体シミュレーターとして実績のある３次元
2.1 支配方程式 4)
積分差分法コード TOUGH2V2（Pruess 他, 1999）をベー
2.1.1 熱・流体移動
TOUGHREACT で用いられている熱、流動、物質移
スに開発された、化学反応を考慮できる物質移行解析コ
ードである。
行の支配方程式は、次式の熱・質量保存則により一般化
地下水だけでなくガスや石油などの多相流体を取り扱
＊１
される。
技術センター土木技術研究所地盤・岩盤研究室
48-1
地下水－岩石反応を考慮した水理・物質移行解析技術
∂M κ
= −∇Fκ + q κ
∂t
度が定まれば一意に定まる。そこで、移流・拡散による
(1)
物質移行を成分毎に考える。拡散係数が化学種によらず
3
ここに、Mκ：成分κの質量(kg/m )またはエネルギー密度
3
2
(J/m )、Fκ：成分κのフラックス(kg/m s)、qκ：成分κのシンク／
3
3
ソース項（kg/m s または J/m s）である。
水：
M w = φ( Sl ρl X wl + S g ρ g X wg )
一定と仮定すると、水相内での物質移行は、化学成分
（j=1,2...,Nc）ごとに、次の移流・拡散方程式により記述さ
れる。
M j = φS l C jl
Fw = X wl ρ l u l + X wg ρ g u g
F j = u l C jl − ( τφS l Dl )∇C jl
q w = q wl + q wg
ガス：
q j = q jl + q js + q jg
(2)
M c = φ( S l ρ l X cl + S g ρ g X cg )
ここに、Cjl：水相での成分 j の濃度、Dl：水相中での拡
Fc = X cl ρ l u l + X cg ρ g u g
散係数(m2/s)、τ：屈曲度、qjβ ：β相の化学反応により生
qc = qcl + qcg + qcr
熱：
じる成分 j のシンク／ソースである。化学反応と物質移
M h = φ( S l ρ lU l + S g ρ g U g ) + ( 1 − φ )ρ sU s
Fh =
∑h ρ u
β = l ,g
β
β
β
(5)
行は、このシンク／ソース項でカップリングされている。
− λ ∇T
2.1.3 化学反応
ここに、uβ：β相のダルシー流速(m/s)、Sβ：β相の飽和
(1) 水溶性錯体
3
度、ρβ：β相の密度（kg/m ）、φ：空隙率、Xκβ：β相にお
局所平衡（local equilibrium）の仮定のもと、 (4)式に質量
ける成分κの質量分率、λ：熱伝導率(W/m℃)、T：温度
作用則を適用すると、水溶性錯体(aqueous complexation)の
(℃)である。Uβ：β相の内部エネルギー(J/kg)、 w：水、
濃度は、成分濃度の関数として次式で求められる。
c：ガス、h：熱、l：液体、g：気体、s：固体、r：反応
Nc
生成である。
j =1
ここに、ci：i 番目の水溶性錯体のモル濃度(mol/kg)、cj：j 番
β相の流速は、次式のダルシー則により計算する。
uβ = −k
k rβ
μβ
( ∇P
β
− ρβ g
)
目の成分のモル濃度(mol/kg)、γj、γi：活動度係数、Ki：平衡
(3)
β = l,g
定数である。
ここに、Pβ：β相の圧力(MPa)、μβ: β相の粘性係数(kg/ms)、
k：絶対浸透率(m2)、krβ：β相の相対浸透率である。
2.1.2
(6)
ν
-1 -1
c i = K i γ i ∏ c νj ij γ j ij
(2) 酸・塩基反応、酸化還元反応
酸・塩基反応と酸化・還元反応の取り扱いは水溶性錯体の
場合と同様である。酸化・還元反応の定式化では、電子を溶
存化学種に加えることで、Fe2+⇔Fe3+ + e-のような半反応式
多成分反応系における物質移行
地球化学システムを構成する全ての化学種(chemical
を取り扱う。電子には仮想的な濃度が与えられるが、電子の
species)は、いくつかの独立した成分（component）の組み合
活動度係数は通常 1 とする。電子を溶存酸素で置き換えて
わせにより記述できる。この独立成分を基本化学種（basis
扱う場合もある。
species または primary species）と呼び、基本化学種の組み
(3) 鉱物の溶解／沈殿
一般に、水相中における鉱物 m の溶解／沈殿は、イ
合 わ せ で 表 現 さ れ る 化 学 種 を 二 次 化 学 種 （ secondary
オン活動度積と平衡定数 Km の比で定義される飽和指数
species）と呼ぶ。
例えば、HCO3-,
CO2(aq),
+
などの化学種は、Ca, H ,
CaHCO3+,
HCO3-の
-
CaCO3, CaOH , CH4
3 成分の組み合わせで
SIm によって判定できる（平衡：SIm=0、過飽和（沈
殿）
：SIm>0、不飽和（溶解）
：SIm<0）
。
NC
表現できる。これを一般化すると次式のようになる。
NC
S i = ∑ ν ij S j
j =1
∏ cν
j
i = 1,..., N R
SI m = log Ω m = log
(4)
ここに、Si：二次化学種、Sｊ：基本化学種、Nc：基本化学種
（成分）の数、NR:化学反応の数、νij：i 番目の化学反応にお
mj
γ νj
mj
j=1
X mλm K m
m = 1... N P
(7)
ここで、Xm：鉱物 m のモル分率、λm:活動度係数、NP:
鉱物数である。 通常、いずれも 1 とする。
鉱物の溶解／沈殿は、鉱物種や水質条件によるが、水
ける成分 j の化学量論係数(stoichiometric coefficient)である。
用いる基本化学種は、考慮する化学種や鉱物種を踏まえて
溶性化学種の反応よりもはるかにゆっくりと進行する。
解析者が選択する。
したがって、鉱物反応を平衡反応として扱うのではなく、
全ての水溶化学種が化学平衡状態にある場合、水相中
鉱物の反応速度を考慮することが必要となる。通常、反
の化学種濃度は、全ての成分（または基本化学種）の濃
応速度は実験的に求められ、それらは次式の形に集約さ
48-2
大成建設技術センター報
第38号(2005)
れる（Lasaga ら, 1994）
。
n
θ
r n = f ( c1 , c 2 ,..., c N C ) = ± k n An a H + Ω n − 1
η
(8)
n = 1,2 ,..., N q
ここで、Nq：反応速度モデルの鉱物数、rn が正であれば
溶解速度、負であれば沈殿速度を表す。kn は反応速度定
数（mol/m2/s）
、An は比反応表面積（m2/kg H2O）である。
θとηは実験的に決まるパラメータであり、多くの場合、
いずれも 1 とする。aH+は H+の活動度、n は H+の触媒作
用に関する実験的な定数である。反応速度定数の温度依
存性はアレニウスの式に従う。
⎡− E ⎛ 1
1 ⎞⎤
(9)
k = k 25 exp ⎢ a ⎜ −
⎟
298
R
T
.15 ⎠⎥⎦
⎝
⎣
、k25：温度 25℃
ここに、Ea：活性化エネルギー（J/mol）
での反応速度定数（mol/m2/s）、T：絶対温度（K）、R：
気体定数（=8.314 J/mol/K）である。
(4) ガスの溶解／脱ガス
水相中あるいは気相中におけるガスの反応は平衡反応
として扱う。質量作用則より、
ν
P f Γ f K f = ∏ Nj=c 1 c νj fj γ j fj
(10)
図-1 TOUGHREACT の計算フローチャート
（文献 4 の Figure4.2 を和訳、一部加筆）
Flow chart of TOUGHREACT program4)
ここで、添字 f は気体を表し、p は気体分圧（bar）
、Γは
気体のフガシティー係数である。低圧力条件下では理想
気体として扱えて、Γ =1 と仮定して良い。しかし地熱
地帯などの高温・高圧下では、その条件に応じたΓの補
正が必要となる 9)。
2.2 化学反応の数値解法
上記のうち、化学反応の数値解法について少し述べる。
化学成分ごとの物質収支を考え、時間 t=0 における成分
(5) 陽イオン交換
陽イオン交換は、交換性陽イオンとイオン交換基の間
j の全濃度が既知とすると、それは次式で与えられる。
換係数と呼ばれる。溶液のイオン強度により、イオン交
換係数が変化することが良く知られている。イオン交換
反応は、Gaines-Thomas convention に従うと、次のよう
に一般化して表される 10)。
(
)
(
)
1
1
1
1
Si +
Xνj − S j ⇔
X νi − S i + S j
νi
νj
νi
νj
電荷分率（equivalent fraction）βi で近似すると、成分 j の水
相濃度は次式により求められる。
j =1
⎛ β ⎞
c j γ j ⎜⎜ i ⎟⎟
⎝ ci γ i ⎠
n=1
j = 1... N C
(13)
度、c は化学種濃度、添字の j, k, m, n は、基本化学種、
錯化学種、鉱物（平衡モデル）、鉱物（反応速度モデ
論係数である。時間ステップをΔtr 進めると、
NX
Np
Nq
k =1
m= 1
n= 1
(
T j = c j + ∑ ν kj ck + ∑ ν mj cm + ∑ ν nj cno − rn Δtr
)
j = 1... N C (14)
ここで、rn は溶解速度（mol/kgH2O･s）である。質量保
存則より、
(15)
j = 1... N C
これら Nc 個の方程式を Newton/Raphson 法により解く。
(12)
残差項を F と書けば、
ここに、Kij*：イオン選択係数、Nw：層間陽イオンの数
である。
m=1
T j0 = T j
ν j / νi
=1
k =1
鉱物の数であり、νkj、νmj、νnj は同様に対応する化学量
は層間陽イオン交換サイトを表す。層間陽イオンの活動度を
* −ν j
ij
Nq
ここで、右上添字の 0 は時間を表す。Tj は成分 j の全濃
（電荷数）であり、Si、Sj は交換性陽イオン、Xνi－Si と Xνj－Sj
Nw
Np
ル）を各々表す。Nc、Nx、Np、Nq はそれらの化学種や
(11)
ここに、νi とνj は、各々、溶存、層間陽イオンの化学量論数
∑ (K )
NX
0
0
0
o
0
T j = c j + ∑ ν kj c k + ∑ ν mj cm + ∑ ν nj c n
の平衡反応によって扱われる。その平衡定数はイオン交
Fjc = T j − T j0 = 0
j = 1... N C
(16)式に(13)式と(14)式を代入すると、
48-3
(16)
地下水－岩石反応を考慮した水理・物質移行解析技術
NX
Np
Nq
k=1
m=1
n=1
0
c
F j = c j + ∑ ν kj ck + ∑ ν mj cm − ∑ ν nj rn Δtr − U j j = 1... N C
NX
Np
k=1
m=1
0
0
0
0
U j = c j + ∑ ν kj ck + ∑ ν mj cm
(17)
水溶性錯体の濃度 ck は(6)式の質量作用式により、cj
（j=1...Nc）の関数として求められる。同様に溶解速度
rn も(8)式に示したように cj の関数である。鉱物（平衡
モデル）の濃度 cm を求めるために、次式を追加する。
Fjp = SI m m = 1... N p
(18)
(17)式と(18)式の（Nc+Np）個の方程式を解くことにより、
（Nc+Np）の未知数（c1, c2, ..., cNc; c1, c2, ..., cNp）を求め
ることができる。
図-2 地下深部への CO2 貯留の概念
Concept of CO2 sequestration in deep underground
2.3 解析手順
(http://www.esru.strath.ac.uk/より)
TOUGHREACT の計算フローチャートを図-1 に示す。
TOUGH2V2 で計算された流速、温度を地球化学計算コ
ードに受け渡し、化学成分毎に物質移行（移流・拡散）
解析を行なう。ここで求めた各成分濃度を地球化学計算
コードに引き渡し、グリッドブロック毎に、水溶化学種
の濃度、鉱物生成・溶解量、ガス溶解・溶離量、表面収
着量などを計算する。所定の収束条件に達するまで物質
移行解析と地球化学計算を繰り返す。
CO2 が地中でトラップされる機構は、流体移動トラッ
プ、溶解トラップ、鉱物トラップの３つに分けられる。
流体移動トラップは、ガスが低浸透性の背斜構造などに
トラップされることであり、短期的には重要な機構であ
る。溶解トラップは地下水中への炭酸ガスの溶解を意味
する。炭酸塩鉱物化による固定化を意味する鉱物トラッ
プは最も長期的に安定性が高いトラップであり、固定化
に重要なメカニズムである。一般的に炭酸塩鉱物の生成
3. 適用事例 －CO2地下貯留シミュレーション－
速度は他の鉱物と比べて比較的早いものの、炭酸塩鉱物
TOUGHREACT は、米国の高レベル放射性廃棄物処
表-1 初期鉱物組成と想定した二次鉱物
分場に決定したヤッカマウンテンのサイト評価や CO2
地下貯留に関する解析、浅成銅鉱床の形成過程の解析な
ど、多くの適用が進められている。
Mineral
Initial mineral volume fractions and possible
secondary mineral phases
Chemical composition
Vol.%
of medium
以下に紹介するのは、米国テキサス州の深度 2000m
の砂岩帯水層を対象とし、超臨界 CO2 を注入した際に
生じる鉱物生成を評価することを目的して実施した筆者
らの数値解析事例である 11), 12)。
3.1 背景・目的
近年、地球温暖化ガスの排出を低減するための方策と
して、大規模発生源から回収した CO2 を地下深部の帯
水層に注入して貯留・隔離する方法が即効的で有効な手
段として有望視されている。すでに北欧、米国を中心に、
天然ガスから分離した CO2 の帯水層注入や、 EOR（石
油増進回収）を目的とした油層への CO2 圧入が行なわ
れている。国内でも、新潟県長岡にて実証実験が行なわ
れるなど、実現可能性の調査・研究が活発化しつつある。
地中貯留の対象となる地層としては、深部塩水の帯水層
やすでに枯渇した石油・ガス田、採掘不能な石炭層など
が考えられている。
48-4
Primary:
quartz
kaolinite
calcite
illite
kerogen-OS
oligoclase
K-feldspar
Na-smectite
chlorite
hematite
Secondary:
anhydrite
magnesite
low-albite
dolomite
siderite
Ca-smectite
pyrite
ankerite
dawsonite
alunite
SiO2
Al2Si2O5(OH)4
CaCO3
K0.6Mg0.25Al1.8(Al0.5Si3.5O10)(OH)2
C64H102O40S10
CaNa4Al6Si14O40
KAlSi3O8
Na0.290Mg0.26Al1.77Si3.97O10(OH)2
Mg2.5 Fe2.5Al2Si3O10(OH)8
Fe2O3
CaSO4
MgCO3
NaAlSi3O8
CaMg(CO3)2
FeCO3
Ca0.145Mg0.26Al1.77Si3.97O10(OH)
FeS2
CaMg0.3Fe0.7(CO3)2
NaAlCO3 (OH) 2
KAl3 (OH) 6(SO4) 2
40.522
1.4105
1.3503
0.6678
0.0
13.856
5.7253
2.7279
3.1892
0.3479
大成建設技術センター報
の生成に必要な Ca、Fe などの元素の供給に長石類や粘
第38号(2005)
表-2 鉱物の反応速度パラメータ（上段: 生成，下段: 溶解）
Kinetic rate parameters for minerals considered in
the simulations
土鉱物類などの溶解が必要な場合には、それらの鉱物溶
解が律速過程になる。長石類や粘土鉱物の生成・溶解は
Mineral
数万年以上の時間オーダーでの現象であるため、その炭
酸塩鉱物の生成挙動を室内実験で検証することは困難で
あり、数値解析的な予測に頼らざるを得ない。
3.2 解析条件
CO2 と水を 2：1 の質量比で混合したものである。注入
を解析した。砂岩の浸透性や造岩鉱物は、上記地点で取
率は 30%、絶対浸透率は 1×10
2
m 、岩石圧縮率は 1×
10-9 Pa-1 とした。相対浸透率および毛管圧力の関数は、
van
初期鉱物組成（表-1）は、対象帯水層を構成する新第
Equilibrium
1.26x10-14
87.5
0
1)
98
K-feldspar
1.00x10-12
57.8
0
2)
98
3.55x10-10
51.9
0.4
2)
98
1.00x10-13
62.8
0
3)
1516
4.37x10-12
62.8
0.17
3)
1516
4.47x10-10
62.8
0
4)
98
4.37x10-5
19.8
1.0
4)
98
1.26x10-9
62.8
0
5)
98
1.02x10-3
20.9
0.9
6)
98
1.26x10-9
62.8
0
A)
98
1.02x10-3
20.9
0.9
A)
98
1.26x10-9
62.8
0
A)
98
1.02x10-3
20.9
0.9
A)
98
1.26x10-9
62.8
0
A)
98
1.02x10-3
20.9
0.9
A)
98
kerogen-OS
1.00x10-13
0.0
0.0
Assumed
100
oligoclase
1.00x10-12
57.8
0
B)
100
3.55x10-10
51.8
0.4
B)
100
1.00x10-12
67.8
0
2)
98
2.04x10-10
59.8
0.5
2)
98
1.00x10-13
62.8
0
C)
1516
4.37x10-12
62.8
0.17
C)
1516
1.00x10-13
62.8
0
C)
1516
4.37x10-12
62.8
0.17
C)
1516
1.00x10-13
62.8
0
C)
1516
4.37x10-12
62.8
0.17
C)
1516
pyrite
4.00x10-11
62.8
0
7)
129
hematite
4.00x10-11
62.8
0
D)
129
chlorite
2.51x10-12
62.8
0
8)
98
kaolinite
magnesite
siderite
dolomite
ankerite
dawsonite
albite-low
三系砂岩（Frio 砂岩）の薄片観察や X 線分析結果から
決定した。石英、カオリナイト、方解石、イライト、斜
長石、カリ長石、Na-スメクタイト、緑泥石、ヘマタイ
ト、ケロジェン質の有機物を選定した。同表中には CO2
の注入後の生成を想定した二次鉱物をあわせて示した。
二次的な炭酸塩鉱物としてドロマイト、シデライト、マ
グネサイト、アンケライト、ドーソナイトを考えた。
鉱物の生成・溶解は前述の反応速度モデルに従うもの
A
cm2/g
quartz
Genuchten の式（ただし、気相の相対浸透率は
Corey 式)を用いた。
Reference
Equilibrium
得された岩石サンプルの試験結果により設定した。空隙
-13
n
calcite
圧力は 260 bar とした。圧入を 100 年間継続した後、そ
の後の 1 万年間にわたる鉱物溶解・生成などの化学反応
Ea
KJ/mol
anhydrite
米国テキサス州の地下深部 2,000m の砂岩帯水層（厚
さ 10m）への CO2 圧入を対象とした。圧入する流体は
k25
mol m-2s-1
Na-smectite
Ca-smectite
Illite
とした。ただし、今回の時間スケールに対して反応が十
1) Tester et al. (1994), 2) Blum and Stillings (1995), 3) Nagy (1995),
分速いと考えられる方解石と石膏は平衡モデルに従うと
4) Chou et al. (1989), 5) Steefel (2001), 6) Gautelier et al. (1999),
した。解析に用いた反応速度定数、反応表面積ならびに
7) Ague and Brimhall (1989), 8) Malmstrom et al. (1996)
活性化エネルギーなどのパラメータを表-2 に示した。
多くは表中の文献の実験値によるものであるが、適当な
実験値が無い場合には、性質の近い他の鉱物の値を仮定
した。例えば、ドロマイト、アンケライト、ドーソナイ
トなどの炭酸塩鉱物にはシデライトの値を使用した。
上記の鉱物組成を有する岩石と 1M の塩水のバッチ反
応を、CO2 分圧 1×10-2 bar、温度 75℃の条件で計算し、
準平衡状態に至った反応開始 10 年後の水質を解析に用
いる初期水質とした（表-3）
。
3.3 解析結果
地下水の飽和度 Sw の軸方向分布を図-3 に示す。CO2
ガス飽和度は Sg=1－Sw である。CO2 ガスのフロントは
圧入開始から 100 年後に約 2,000m まで達するのが分か
る。注入完了後、ガス飽和度は緩やかに減少している。
pH 変化を図-4 に示す。CO2 の注入初期には注入孔か
48-5
A) Set to siderite, B) Set to K-feldspar, C) Set to kaolinte, D) Set to pyrite
表-3 初期水質（化学成分ごとの全濃度）
Initial total dissolved component concentrations
Component
Ca2+
Mg2+
Na+
K+
Iron
SiO2(aq)
Carbon
S
Al+3
ClpH
Temperature
Conc. mol/kg
3.232 x10-3
1.532 x10-7
0.9889
7.517 x10-3
2.421 x10-5
7.262 x10-4
4.320 x10-2
1.323 x10-9
2.657 x10-8
1.000
7.343
75oC
地下水－岩石反応を考慮した水理・物質移行解析技術
ら約 20m の範囲に低 pH 領域を形成し、方解石が溶解
ドーソナイトの生成中にカオリナイトが不安定となる点
する。0.4 年ほど経つと全ての方解石が溶解し、注入孔
で解析結果の傾向と一致する。
Watson ら
周辺の pH は 3.4 まで低下する。しかしその後、カオリ
14)
は、南オーストラリアの高圧 CO2 ガス貯
ナイトなどの鉱物生成により pH は緩やかに上昇、注入
留層の鉱物を調査した。この貯留層は還元的であり、シ
が完了した 100 年以降は、pH=6 程度まで回復する。
デライトやアンケライトの安定性に必要な Fe2+が豊富で
方解石の存在量（全鉱物に対する体積比）を図-5a に
あるので、より今回の解析条件に近い。そして、緑泥石、
示す。方解石は注入後、非常に早い時間で溶解する。カ
方解石、長石類の溶解と、シデライト、アンケライト、
リ長石、斜長石、緑泥石なども注入孔周辺の低い pH 領
石英の生成などシミュレーションと同様の傾向が認めら
域で溶解する傾向にある（図-5）
。
れている。ただし、地下水がドーソナイトに関して過飽
それらの初期鉱物が溶解する一方、石英、イライト、
和状態にも関らず、ドーソナイトの存在は認められてい
カオリナイト、スメクタイトなどの二次鉱物が CO2 プ
ない。その理由は明らかでないが、CO2 圧力や Cl や Na
リューム全域に渡って生成する（図-5）
。
の濃度が今回の解析条件よりも低いためと考えられる。
CO2 を 固 定 化 す る 炭 酸 塩 鉱 物 は 、 マ グ ネ サ イ ト
（ MgCO3 ）、 シ デ ラ イ ト （ FeCO3 ）、 ア ン ケ ラ イ ト
1.00
（CaMg0.3Fe0.7(CO3)2）、ドーソナイト（NaAlCO3(OH)2）
である。ただし、マグネサイトの生成は注入孔近傍の酸
今回の計算で CO2 の固定化に最も寄与している炭酸
塩鉱物は、アンケライト、ドーソナイトと若干のシデラ
イトである。方解石はむしろ溶解しており、固定化には
ほとんど寄与していない。アンケライトの生成に必要な
Fe2+は、主に緑泥石の溶解により供給されている。また、
ドーソナイトの生成に必要な Na+は斜長石の溶解により
供給されている。従って、それらの鉱物の存在量が鉱物
トラップの容量を左右していると言える。今回のシミュ
Water saturation
性ゾーンに限られる。ドロマイトは生成していない。
10 yr
100 yr
1,000 yr
10,000 yr
0.95
0.90
0.85
0.80
0.75
-1
10
10
0
レーションでの 1 万年後の CO2 の固相への固定化量は、
注入孔から 1,000m の位置において 76kg/m3 に達し（図-
3
10
4
Water saturations alone radial distance
at different times
3.4 フィールド観察との比較
自然界での水－鉱物の反応現象は室内実験よりも複雑
8
であり、実験室内の理想的な条件下で求められた鉱物の
7
る。例えば、複数の鉱物からなる自然系での反応表面積
6
は、他の鉱物のコーティングの影響などにより、純粋な
5
実験の場合と比べて小さくなると考えられる。また、フ
ィールドでの鉱物観察から評価した反応速度が室内実験
pH
反応速度定数や反応表面積の適用性には多くの問題があ
4
3
結果より数オーダー遅いというような報告もある。
高圧 CO2 存在下における砂岩の変質と鉱物について
記載された例は非常に少ない。Moore ら
2
図-3 水相飽和度分布の時間変化
6, 7）
、帯水層（厚さ 10m）の全体では約 8Mt となる。
13)
1
10
10
10
Radial Distance, m
は、Arizona
州と New Mexico 州の境界付近の Springerville-St.Johns
のペルム系シルト岩において、斜長石、カリ長石の溶解
箇所でのドーソナイトの生成を確認している。すなわち、
ドーソナイトの生成に長石の溶解が必要となっている点
で今回のシミュレーションと整合的である。また、ドー
ソナイトの沈殿後にカオリナイトが生成していることは、
48-6
10 yr
100 yr
1,000 yr
10,000 yr
2
1
0
-1
10
10
0
1
2
10
10
10
Radial Distance, m
3
10
図-4 pH 分布の時間変化
pH distribution at different times
4
0.0E+00
-5.0E-03
10 yr
100 yr
1,000 yr
10,000 yr
100
101
102
103
Radial Distance, m
104
8.0E-03
6.0E-03
4.0E-03
2.0E-03
0.0E+00
10-1
Ankerite
4.0E-02
10 yr
100 yr
1,000 yr
10,000 yr
3.0E-02
2.0E-02
1.0E-02
0.0E+00
-1
10
10
0
1
2
10
10
10
Radial Distance, m
3
10
4
6.0E-02
10 yr
100 yr
1,000 yr
10,000 yr
-1.5E-01
-1
10
0
1
2
10
10
10
Radial Distance, m
3
10
4
Change of abundance (volume fraction)
Change of abundance (volume fraction)
-5.0E-02
10
10 yr
100 yr
1,000 yr
10,000 yr
1.5E-06
1.0E-06
5.0E-07
0.0E+00
0
1
2
10
10
10
Radial Distance, m
3
(j) Ca-スメクタイト
10
4
Change of abundance (volume fraction)
Change of abundance (volume fraction)
Ca-Smectite
10
8.0E-03
6.0E-03
4.0E-03
2.0E-03
0.0E+00
10-1
3.0E-02
2.0E-02
1.0E-02
0.0E+00
-1
10
10
0
1
2
10
10
10
Radial Distance, m
3
10
4
1.0E-02
100
1.0E-02
0.0E+00
-1.0E-02
-2.0E-02
-3.0E-02
10 yr
100 yr
1,000 yr
10,000 yr
-5.0E-02
-6.0E-02
-1
10
10
0
1
2
10
10
10
Radial Distance, m
1.4E-01
3
Chlorite
-1.0E-02
10 yr
100 yr
1,000 yr
10,000 yr
-2.0E-02
-3.0E-02
10-1
100
10
4
10 yr
100 yr
1,000 yr
10,000 yr
Na-Smectite
1.0E-01
8.0E-02
6.0E-02
4.0E-02
2.0E-02
10
0
1
2
10
10
10
Radial Distance, m
101
102
103
Radial Distance, m
3
(k) Na-スメクタイト
2.5E-02
10 yr
100 yr
1,000 yr
10,000 yr
Kaolinite
2.0E-02
1.5E-02
1.0E-02
5.0E-03
0.0E+00
10-1
100
101
102
103
Radial Distance, m
10
4
7.0E-02
104
10 yr
100 yr
1,000 yr
10,000 yr
Illite
6.0E-02
5.0E-02
4.0E-02
3.0E-02
2.0E-02
1.0E-02
0.0E+00
10-1
100
101
102
103
Radial Distance, m
(l) イライト
図-5 CO2 圧入後の鉱物の生成・溶解量（初期状態からの体積比変化、＋：生成、－：溶解）
Changes of mineral abundances due to dissolution /precipitation after CO2 injection
48-7
104
(i) カオリナイト
1.2E-01
0.0E+00
-1
10
104
(f) 緑泥石
K-feldspar
-4.0E-02
101
102
103
Radial Distance, m
0.0E+00
(h) カリ長石
2.0E-06
-5.0E-07
-1
10
1.0E-02
(c) シデライト
10 yr
100 yr
1,000 yr
10,000 yr
4.0E-02
(g) 斜長石
2.5E-06
104
10 yr
100 yr
1,000 yr
10,000 yr
1.2E-02
(e) ドーソナイト
Oligoclase
-1.0E-01
101
102
103
Radial Distance, m
Dawsonite
5.0E-02
(d) アンケライト
0.0E+00
100
Siderite
1.4E-02
(b) マグネサイト
Change of abundance (volume fraction)
Change of abundance (volume fraction)
(a) 方解石
5.0E-02
10 yr
100 yr
1,000 yr
10,000 yr
Change of abundance (volume fraction)
-1.5E-02
10-1
Magnesite
Change of abundance (volume fraction)
-1.0E-02
1.0E-02
Change of abundance (volume fraction)
Calcite
第38号(2005)
Change of abundance (volume fraction)
5.0E-03
Change of abundance (volume fraction)
Change of abundance (volume fraction)
大成建設技術センター報
104
地下水－岩石反応を考慮した水理・物質移行解析技術
謝辞
解析コードの導入や研究の遂行にご助力をいただいた
CO2 sequestererd (kg/m3 medium)
100
10 yr
100 yr
1,000 yr
10,000 yr
80
技術センター土木技術研究所の下茂道人氏に感謝する。
本報のうち,CO2 の地下貯留シミュレーションは,著者が
米国ローレンスバークレー国立研究所に客員研究員とし
60
て在籍中（2003.9～2004.8）に行った研究の一部である。
40
同研究所の Karsten Pruess 博士, Tianfu Xu 博士には終始
懇切な指導を頂いた。ここに記して謝意を表する。
20
0
参考文献
-20
-40
-1
10
0
10
1
2
10
10
Radial Distance, m
3
10
図-6 炭酸塩鉱物による CO2 固定化量分布の時間変化
CO2 sequestered (kg/m**3 medium)
Cumulative CO2 sequestered by carbonate mineral
precipitation at different times
80
70
60
50
40
30
20
10
0
0
2000
4000
6000
8000 10000
Time (yr)
図-7 CO2 固定化量（注入孔より 1000m 地点）
Cumulative CO2 sequestration with time at 1000 m
distance
4. おわりに
化学反応を考慮できる水理・物質移行解析コード
TOUGHREACT と CO2 地下貯留シミュレーションへの
適用事例を紹介した。筆者らは現在,北海道幌延町の核
燃料サイクル開発機構の地下研究施設周辺の広域地下水
水質形成シミュレーション
密性解析
15)
や LPG 地下備蓄基地の気
16)
などへ, 本コードの適用を進めている。今後
は,地下水・土壌汚染問題やメタンハイドレート開発な
ど幅広い環境・エネルギー分野に活用していく所存であ
る。
1) Parkhurst, D.L., and C.A.J. Appelo, 1999, User’s guide to
PHREEQC(version 2)--A computer program for speciation, batch-reaction,
one dimensional transport, and inverse geochemical calculations: U.S.
Geological Survey Water-Resources Investigations Report 99-4259, 321 p.
2) Wolery, T. J., 1992. EQ3/6: Software package for geochemical modeling of
aqueous systems: Package overview and installation guide (version 8.0),
Lawrence Livermore National Laboratory Report UCRL-MA-110662 PT I,
Livermore, California.
3) Xu, T., and K. Pruess, 2001. Modeling multiphase fluid and reactive
geochemical transport in variably saturated fractured rocks: 1.
Methodology, Am. J. Sci., v. 301, p.16-33.
4) Xu, T., E. Sonnenthal, N. Spycher, K. Pruess, 2004. TOUGHREACT
User’s Guide: A simulation program for non-isothermal multiphase
reactive geochemical transport in variably saturated geologic media,
Lawrence Berkeley National Laboratory Report LBNL-55460, Berkeley,
California.
5) Nitao, J.J. 1998. User's manual for the USNT module of the NUFT code.
Version 2.0. Rep. UCRL-MA-130653. Lawrence Livermore National
Laboratory, Livermore, CA.
6) Mayer, K.U., 1999. A numerical model for multicomponent reactive
transport in variably saturated porous media. Ph.D. thesis, Department of
Earth Sciences, University of Waterloo.
7) Parkhurst, D.L., Kipp, K.L., Engesgaard, Peter, and Charlton, S.R., 2004.
PHAST--A program for simulating ground-water flow, solute transport,
and multicomponent geochemical reactions: U.S. Geological Survey
Techniques and Methods 6-A8, 154 p.
8) White, S.P. 1995. Multiphase nonisothermal transport of systems of
reacting chemicals, Water resources research, Vol.31(7), pp.1761-1772.
9) Spycher, N. F., and Reed, M. H. 1988. Fugacity coefficients of H2, CO2,
CH4, H2O and of H2O-CO2-CH4 mixtures: A virial equation treatment for
moderate pressures and temperatures applicable to calculations of
hydrothermal boiling, Geochim. Cosmochim. Acta, v. 52, p. 739-749.
10) Appelo, C. A. J., and Postma, D., 1993, Geochemistry, groundwater and
pollution, Rotterdam, The Netherlands, Balkema, 536 pp.
11) Xu, T., J. Apps, K. Pruess, H. Yamamoto, 2004. Injection of CO2 with H2S
and SO2 and Subsequent Mineral Trapping in Sandstone-Shale Formation,
Lawrence Berkeley National Laboratory Report LBNL-57426, Berkeley,
California.
12) Xu, T., J. Apps, K. Pruess, H. Yamamoto, 2004. Mineral Sequestration of
CO2 Mixed with H2S and SO2 in Sandstone-Shale Formation, Eos Trans.
AGU, 85(47), Fall Meet. Suppl., Abstract GC51B-1058
13) Moore, J., Adams, M., Allis, R., Lutz, S., Rauzi, S., 2004. Mineralogical
and geochemical consequences of the long term presence of CO2 in natural
reservoirs: An example from the Springerville-St. Johns field, Arizona and
New Mexico, USA. Submitted to Chemical Geology.
14) Watson,M.N., Zwingmann, N., Lemon, N.M., 2002. The Ladbroke GroveKatnook carbon dioxide natural laboratory: A recent CO2 accumulation in a
lithic sandstone reservoir. Extended Abstract E1-5, Sixth International
Conference on Greenhouse Gas Control Technologies(GHCT-6), Kyoto,
Japan, 6p.
15) Yamamoto, H., T. Kunimaru, M. Shimo, Y. Fujiwara, M. Ono, T. Xu, M.
Laaksoharju, 2005. Reactive chemical transport simulations for evaluation
of ambient chemical conditions at Horonobe URL site, Hokkaido, Japan.
The 2nd Japan-Korea Joint Workshop on Radioactive Waste Disposal 2005,
Tokyo, pp.159-172.
16) Yamamoto, H., K. Pruess, 2004. Numerical simulation of leakage from
underground LPG storage caverns, LBNL Report LBNL-56175, Berkeley,
California.
48-8