JP5334221B1 - 熱応答試験および揚水試験の解析方法および解析プログラム - Google Patents

Abstract

【課題】短時間で精度良く、熱応答試験結果に基づき地層の熱伝導率を求めることができる熱応答試験の解析方法を提案すること。
【解決手段】熱応答試験の解析方法では、コンピュータに熱応答試験による測定値を入力する。コンピュータにはパラメータ同定プログラムがインストールされており、当該同定プログラムにより、パウエルの共役傾斜法に基づく非線形最適化法を用いて、同定対象のパラメータに任意の値を設定して、指数関数の積分形で与えられる熱応答試験用の解析解を第３０項まで計算して地層の温度上昇量を算出する温度上昇量の算出し、算出値と熱応答試験によって得られた測定値との誤差を算出し、この誤差が最小となるようにパラメータの値を変更する処理を繰り返すことでパラメータを同定する。同定後のパラメータを用いて地層の熱伝導率を算出できる。
【選択図】図４

Description

本発明は熱応答試験および揚水試験の解析方法および解析プログラムに関する。さらに詳しくは、共役傾斜法を適用してパラメータ同定による逆解析を行って、短時間で精度良く地盤の熱交換特性および地下水流動特性を評価できる解析方法および解析プログラムに関する。

近年、原子力や火力に代わるクリーンなエネルギーとして、自然エネルギー(再生可能エネルギー）の利用が見直されてきている。その１つとして、地中熱の利用に関心が高まってきている。季節によって変動する外気温に対して、地中の温度は年間を通して約１５℃で安定している。このような特性を持つ地中の熱を有効利用したものとして、地中熱ヒート・ポンプ（Ｇｅｏｔｈｅｒｍａｌ Ｈｅａｔ Ｐｕｍｐ、以下ＧｅｏＨＰと略す）システムがある。これは、地中熱を利用した冷暖房・給湯システムのことである。安定した熱の供給が行なえる地中熱の特性に着目し、夏場には冷房として地中に熱を放出し、冬場には暖房の熱源として利用する。

しかしながら、日本ではＧｅｏＨＰシステムの普及は、欧米に比べ、遅れているのが実情である。その理由としては、認知度の低さ、井戸の掘削等における初期費用の問題、および、日本では欧米とは違い、その土地の風土によって熱物性値が異なるため、経済的なシステムづくりが依然として確立されていない、などが挙げられる。そのため、ＧｅｏＨＰシステムの普及には初期費用を抑えるためのシステムの効率化が必要不可欠である。

システムの効率化を図るための地盤調査試験として、熱応答試験が知られている。熱応答試験は主に地層の熱伝導率と熱交換井の熱抵抗を評価するための試験であり、地中熱交換井の熱交換挙動を予想する上で不可欠な情報である。熱応答試験を実施して、適正な井戸の本数・長さを決定することは、ＧｅｏＨＰシステムの初期費用の削減に極めて重要である。

地盤調査のための熱応答試験は、主にクローズド型方式における垂直型地中熱交換井を対象としている。このための熱応答試験装置は、地上部の水タンク内で電気ヒーターによって熱媒体を加熱し、循環ポンプにより配管内を循環させる。加熱された熱媒体は地中部を通過する際に放熱するため、Ｕ字管の入口・出口の温度差を測定し、その地層の熱交換特性を評価する。このような熱応答試験装置は、特許文献１に開示されており、熱応答試験装置を用いた熱応答試験方法は非特許文献１に開示されている。

特開２００４−３０１７５０号公報

「ボアホール型地中熱交換機に対する加熱法による熱応答試験の標準試験方法ｖｅｒ．２」、長野克則著、地下熱利用とヒートポンプシステム研究会編、財団法人ヒート・ポンプ・蓄熱センター

ここで、熱応答試験の代表的な解析法としてはＫｅｌｖｉｎの線源関数や円筒型熱源関
数に基づく解析的手法、および有限差分法や有限要素法を用いた数値モデルによる解析法がある。熱応答試験の解析的手法において最も一般的に用いられている手法としては、データ観測の簡便さや解析の容易さより、Ｋｅｌｖｉｎの線源関数に基づいた作図法が挙げられる。作図法では循環時における熱媒体の温度データの経時変化を用いる方法と循環停止後の地中温度の回復データを用いる方法があるが、循環時のデータを用いる解析が一般的に行われている。以下にＫｅｌｖｉｎの線源関数に基づいた熱応答試験の解析理論を説明する。

地下水流れのない均質、等方性、等圧で無限の広がりをもつ地層へ完全貫入している井戸へ一定量Ｑ_Ｈで熱負荷を与えているものとする。非定常の熱伝導方程式は円筒座標系に変換すると、式（１）で表される。

一方、井戸の半径をｒ_ｗとすると、初期条件および境界条件は次式（２）〜（４）で与えられる。

地層の温度上昇量をＴ_ｒ＝Ｔ−Ｔ_０として、式（１）と初期条件および境界条件を書き直すと次式（５）〜（８）のようになる。

以上の初期条件、境界条件のもとで式（５）を解くと次式（９）〜（１１）のようになる。

式（９）をＫｅｌｖｉｎの線源関数といい、Ｅ（Ｘ）は指数積分である。式（１１）におけるＥ（Ｘ）はＸが十分に小さい場合、一般にはＸ＜０．０５で誤差が２．５％以下の場合に、次式（１２）により近似される。

したがって、式（９）は次式（１３）により近似される。

式（１３）において、熱交換量一定条件の下で時間とともに変化する変数は時間ｔのみなので、時間に依存する項としない項に整理すると、式（１３）は次式（１４）のように
なる。

式（１４）において、ｌｎ（ｘ）＝２．３ｌｏｇ（ｘ）とすると、自然対数が常用対数に変換され次式（１５）となる。

ここで、式（１５）の傾きｍを次式（１６）で表す。

この場合、地層の見かけ熱伝導率λは次式（１７）で表される。

循環時データを用いた作図法では式（１６）に基づいて、熱媒体平均温度の経時変化を示す片対数グラフの傾きよりｍ値を決定し、式（１７）より地層の熱伝導率λを求める。また、近似解を用いて解析を行っていること、初期段階は地中熱交換器や充填剤の影響を受けることを考慮して、試験開始から十分な時間が経過した区間のデータよりｍ値を決定する。しかしながら、解析を行うデータ区間の選定に明確な基準は無く、解析を行なう人によって判断される。

このように、一般的に行われている熱応答試験の作図法による解析では、熱媒体平均温度の経時変化を示す片対数グラフの傾きよりｍ値を決定し、上記の式（１７）より地層の見かけ熱伝導率λを求める。ここでいうグラフの傾きとは、試験開始から数十時間以上経過した後に見られる直線部のデータである。しかしながら、作図法による解析には以下のような問題点がある。
（１）誤差が許容範囲とされている２．５％となる時間ｔか不明のため、傾きｍを決定する直線区間の選定に明確な基準は無く、どの区間のデータを直線部と定めるかは解析を行う人の主観によるため、解析結果が人によって異なる。
（２）初期段階のデータは全て捨ててしまうため、試験期間を長くとる必要があり、コストアップにつながる。
（３）作図法に用いられている解析解は、地下水流動、地温勾配を考慮していない。そのため、地下水流れがある場合、地温勾配がある場合には、正確な熱伝導率を求めることができない。したがって、熱応答試験の解析により求まる熱伝導率は、しばしば誤差を含んだ見かけの熱伝導率となる。

本発明の課題は、このような点に鑑みて、地盤の熱交換特性評価の時間短縮および精度向上を図るために、共役傾斜法を適用したパラメータ同定による熱応答試験の解析法および解析プログラムを提案することにある。

また、本発明の課題は、透水量係数、貯留係数等の地下水流動パラメータを求めるための原位置試験として一般的に行われる揚水試験の解析解と熱応答試験の解析解との同一性に着目し、地盤の地下水流動特性評価の時間短縮および精度向上を図るために、共役傾斜法を適用したパラメータによる揚水試験の解析法および解析プログラムを提案することにある。

上記の課題を解決するために、本発明の熱応答試験の解析方法および解析プログラムでは、地盤の温度上昇量を規定するＫｅｌｖｉｎの線源関数において、その指数積分Ｅ（Ｘ）のべき級数展開式を少なくとも第１６項、好ましくは第３０項まで考慮すると、熱応答試験の初期段階において実測値に近い計算値（準解析解）が得られることに着目し、熱応答試験の初期段階において得られる実測値に対して、パウエルの共役傾斜法による逆解析を用いた非線形最適化手法を適用して、温度上昇量を規定するＫｅｌｖｉｎの線源関数の未知のパラメータを同定して、地盤の熱伝導率を算出している。すなわち、パウエルの共役傾斜法に基づく非線形最適化法を用いて、未知のパラメータを所定の値に設定した場合に得られる温度上昇量の計算値と、熱応答試験によって測定された温度上昇量の実測値との誤差が最小となるように、パラメータの逆解析を行っている。

また、本発明の揚水試験の解析方法および解析プログラムでは、熱応答試験解析解と揚水試験解析解の同一性に着目すると共に、地盤の地下水位低下量を規定するタイスの井戸関数において、その指数積分Ｗ（ｕ）のべき級数展開式を少なくとも第１６項、好ましくは第３０項まで考慮すると、揚水試験の初期段階において実測値に近い計算値（準解析解）が得られることに着目して、揚水試験の初期段階において得られる実測値に対して、パウエルの共役傾斜法による逆解析を用いた非線形最適化手法を適用して、地下水位低下量を規定するタイスの井戸関数の未知のパラメータを同定して、地盤の透水量係数を算出している。

（パウエルの共役傾斜法）
本発明においては熱応答試験あるいは揚水試験の解析方法では、最適化手法による逆解析を利用している。先に述べたように、地下水モデルによる正解析は透水量係数や熱伝導率あるいは分散係数などを既知パラメータとして、地下水位、地下水温、物質の濃度などを計算する構造になっている。しかし、こうしたパラメータを直接精度よく測定することは難しく、むしろ解として得られるべき濃度や地下水温の方が容易に、精度よく観測・測定が可能である。本発明では、地下水位や地下水温などの実測値に計算値が一致するように非線形最適化手法を用いて、各種パラメータを同定する逆解析を採用している。逆解析では、まず目的関数を設定し、それらを最小化することによって各パラメータが得られる。

最適化手法によって最小化される目的関数Ｆ（ｘ）には、次式（２０）で示す実測値と計算値の２乗和などが用いられる。

本発明では、非線形最適化手法として、パウエル（Ｐｏｗｅｌｌ）の共役傾斜法を採用する。パウエル（Ｐｏｗｅｌｌ）の共役傾斜法は、他の多くの最適化手法で必要とされる目的関数の各パラメータに関する微係数を必要とせず、目的関数のみを用いて、その最適パラメータを求めることができ、高い汎用性を持つ。

パウエル（Ｐｏｗｅｌｌ）の共役傾斜法では、２つの離れた点からｖの方向に正値２次形式の最小点を見つけ出すとき、各々で求めた最小点ｘ１，ｘ２を結ぶベクトルは、ｖと共役であるという性質を利用している。この性質を利用して、ｎ回直接探索を繰り返すことにより、最小点を見つけ出すことができる。図１．１に共役方向の概略図を示し、図１．２にパウエル（Ｐｏｗｅｌｌ）の共役傾斜法の基本アルゴリズムを示す。

（熱応答試験の準解析解）
熱応答試験の解析解は以下の式（２１）〜（２３）で表せる。

従来における熱応答試験の解析法では、Ｋｅｌｖｉｎの線源関数の漸近解を第２項まで近似する近似解が用いられている。Ｋｅｌｖｉｎの線源関数の漸近解の項を増加させると正確解に近づく。図２．１は、縦軸をＥ（Ｘ）、横軸を１／Ｘとして、漸近解の項を第２項から徐々に増加させたグラフを示す。グラフから、漸近解の項を増加させると、近似解と正確解が同じ曲線を示す期間が長くなり、正確解に近づいていくことが分かる。このグ
ラフでは漸近解の項を第３０項まで増加させたグラフを描かせているが、第１６項以降は正確解とほとんど一致している。

図２．２は、漸近解を第３０項まで近似した近似解に対する相対誤差の一覧表である。全体として、１／Ｘの値が小さくなるほど誤差は大きくなる傾向にあるが、近似する項数によって誤差の大小が異なり、項数が増加するにつれて誤差が小さくなることが分かる。なお、「−」の表示がある部分は、誤差の値が大きくなりすぎるため省略してある。

以上のことから、漸近解の項数を増加させることで、従来の近似解よりも正確解に近い準解析解が求められることが分かる。本発明においては、Ｋｅｌｖｉｎの線源関数の漸近解を少なくとも第１６項まで、好ましくは、第３０項まで近似した漸近解（準解析解）を用いている。

（熱応答試験の解析における同定対象のパラメータ）
次に、熱応答試験においてΔＴ、Ｑ_Ｈ、ｂは観測データとして既知であり、ｒ、ρｃ、λが未知の状態である。本発明では、求めるパラメータを次のようにα１、α２として、共役傾斜法を組み込んだ同定プログラムにより逆解析を行う。

ここで、α_１、α_２を用いて式（２１）を書き直すと、以下の式（２４）〜（２６）のように表せる。

また、評価尺度としては、次式（２７）と式（２８）を用いることができる。式（２７）は、計算値と実測値の2乗誤差の総和を１／２乗したユークリッドノルムである。式（２８）は、計算値と実測値の絶対値誤差の総和を総データ数で除しており、データ１つあたりの誤差を表す。

（揚水試験解析への適用：揚水試験解析解と熱応答試験解析解の同一性）
次に、熱応答試験と揚水試験は、両者とも熱伝導理論を基として、解析解を導いている。そのため、本発明の熱応答試験の解析方法および解析プログラムは揚水試験の解析にも適応することが可能である。以下に、揚水試験解析解と熱応答試験解析解の同一性を説明する。

揚水試験は現場においては、透水量係数や貯留係数といった地下水流動系パラメータを求めるための原位置試験として最も一般的に行われている。ここでは、被圧帯水層にスクリーンを持つ井戸から揚水を行った時の非定常の地下水流動方程式に対する解析解として現在も広く利用されているタイス（Ｔｈｅｉｓ， １９３５）の正確解を誘導し、熱応答試験解析解との同一性を確認する。

均質、等方性、等圧で無限の広がりをもつ地層へ完全貫入している井戸から一定量Ｑ_Ｗ地下水が揚水されているものとする。均質・等方性の被圧帯水層における水平２次元の地下水の流れを円筒座標系に変換すると次式（３１）で表される。

一方、井戸の半径をｒ_Ｗとすると、初期条件および境界条件は次式（３２）〜（３４）で与えられる。

地下水位低下量をｓ＝ｈ_０−ｈとして式（３１）と初期条件および境界条件を書き直すと次式（３５）〜（３８）のようになる。

ここで、
とおくと、Ｔｈｅｉｓ（１９３５）は式（３４）を上記の初期および境界条件の下で解き
次式を得た。以上の初期条件、境界条件のもとで式（３１）を解くと次式（３９）〜（４１）のようになる。

式（３９）はタイスの井戸関数と呼ばれている。また式（３９）におけるＷ（ｕ）は指数積分であり、ｕが十分に小さい範囲において（一般ｕ＞０．０５で誤差が２．５％以下）、次式（４２）により近似される。

したがって、式（３９）は次式（４３）により近似される。

式（４３）はヤコブの近似解と呼ばれている。

以上のように熱応答試験解析解と揚水試験解析解は、元を辿れば熱伝導理論によるものであり、熱応答試験解析解のｂ（ρｃ）、ｂλが、揚水試験解析解ではＳ_ｃ、Ｔとなるが、式の形自体は同じであり、両者の同一性を確認することができた。これにより、本発明による熱応答試験の解析法を揚水試験にも同様に適応できることが分かる。

本発明によれば、従来の作図法に比べて、短時間で（すなわち、熱移動が安定して、温度上昇量の経時変化を表す曲線の傾きが一定になるまで待たずに）、客観的に、しかも精度良く、地盤の熱伝導率を求めることができる。また、本発明によれば、揚水試験による透水量係数の算出においても同様な作用効果を得ることができる。具体的には、本発明の解析方法および解析プログラムによれば次の効果を得ることができる。
（１）Ｐｏｗｅｌｌの共役傾斜法を適用した逆解析を用いることで、従来の作図法とは異
なり、客観的かつ精度良く、熱応答試験および揚水試験におけるパラメータを評価することができる。
（２）熱応答試験の解析においては、熱伝導率に加えて、従来の解析では求めることのできなかった体積熱容量も求めることが可能である。
（３）準解析解を用いることで従来の近似解では評価できなかった初期段階の試験結果を再現できるので、従来の作図法に比べて、試験期間の大幅な短縮を図ることができ、これにより、コスト縮減も図ることができる。
（４）本発明の解析プログラムを用いることにより、従来手法に比べてパーソナルコンピュータを用いて簡便にパラメータを同定できる。
（５）地温勾配、地下水流動によるｆ／Ｎ値、熱伝導率への影響を確認できる。
（６）試験の初期段階のデータを用いることで、流動や地温勾配に影響されない熱伝導率を求めることができる。
（７）本発明の熱応答試験の解析方法および解析プログラムによって得られる解析結果を用いて、準解析解により熱移動分布を描かせることで、適正な井戸配置を決定できる。

パウエルの共役傾斜法における共役方向を示す概念図である。 パウエルの共役傾斜法の基本アルゴリズムを示す説明図である。 Ｋｅｌｖｉｎの線源関数と近似解を比較したグラフである。 Ｋｅｌｖｉｎの線源関数の漸近解を第３０項まで近似した近似解に対する相対誤差を示す一覧表である。 本発明による同定プログラムの妥当性を検証するためのテストケースの条件を示す一覧表である。 テストケース１におけるＦ値とＤＩＦ値の変化を示すグラフである。 テストケース１におけるｂλの値とｂｒ^２（ρｃ）の値の変化を示すグラフである。 テストケース１における計算値と設定値の比較結果を示すグラフである。 テストケース１の解析結果の一覧表である。 テストケース２におけるＦ値とＤＩＦ値の変化を示すグラフである。 テストケース２におけるｂλの値とｂｒ^２（ρｃ）の値の変化を示すグラフである。 テストケース２における計算値と設定値の比較結果を示すグラフである。 テストケース２の解析結果の一覧表である。 熱応答試験装置の一例を示す説明図である。 第１地点周辺の水理地質断面図である。 第１地点の第１回目の熱応答試験条件の一覧表である。 第１地点の第２回目の熱応答試験条件の一覧表である。 第１地点の熱応答試験の測定機器と試験装置の諸元の一覧表である。 第１地点の第１回目の熱応答試験結果を示すグラフである。 第１地点の第１回目の熱応答試験の逆解析による同定結果の一覧表である。 第１地点の第１回目の熱応答試験の計算値と実測値の比較を示すグラフである。 第１地点の第１回目の熱応答試験の作図法による解析結果を示すグラフである。 第１地点の第２回目の熱応答試験結果を示すグラフである。 第１地点の第２回目の熱応答試験の逆解析による同定結果の一覧表である。 第１地点の第２回目の熱応答試験の計算値と実測値の比較を示すグラフである。 第１地点の第２回目の熱応答試験の作図法による解析結果を示すグラフである。 第２地点の温度検層結果を示すグラフである。 第２地点の熱応答試験条件の一覧表である。 第２地点の熱応答試験結果を示すグラフである。 第２地点の熱応答試験の逆解析による同定結果の一覧表である。 第２地点の熱応答試験の計算値と実測値の比較を示すグラフである。 第２地点の熱応答試験の作図法による解析結果を示すグラフである。 第３地点の熱応答試験条件の一覧表である。 第３地点の熱応答試験結果を示すグラフである。 第３地点の熱応答試験の逆解析による同定結果の一覧表である。 第３地点の熱応答試験の計算値と実測値の比較を示すグラフである。 第１地点での第１回目と第２回目の熱応答試験の解析結果の一覧表である。 第１地点の第２回目と第２地点の熱応答試験の解析結果の一覧表である。 第１地点の第２回目と第２地点の温度検層結果を示すグラフである。 第１地点の第２回目と第３地点の熱応答試験の解析結果の一覧表である。 第１地点の第２回目の熱応答試験における解析データ区間を変化させた場合の解析結果を示すグラフである。 第３地点の熱応答試験における解析データ区間を変化させた場合の解析結果を示すグラフである。 第１地点の第２回目において初期段階データと後半部データを用いた場合の解析結果を示すグラフである。 第２地点において初期段階データと後半部データを用いた場合の解析結果を示すグラフである。 第３地点において初期段階データと後半部データを用いた場合の解析結果を示すグラフである。 第１地点の第１回目における試験期間に対する熱伝導率の変化を示すグラフである。 第１地点の第１回目における試験期間に対するｆ／Ｎ値の変化を示すグラフである。 第１地点の第１回目における各試験期間の解析結果の一覧表である。 第１地点の第１回目における試験期間に対する解析結果を示すグラフである。 第１地点の第２回目における試験期間に対する熱伝導率の変化を示すグラフである。 第１地点の第２回目における試験期間に対するｆ／Ｎ値の変化を示すグラフである。 第１地点の第２回目における各試験期間の解析結果の一覧表である。 第１地点の第２回目における試験期間に対する解析結果を示すグラフである。 第２地点における試験期間に対する熱伝導率の変化を示すグラフである。 第２地点における試験期間に対するｆ／Ｎ値の変化を示すグラフである。 第２地点における各試験期間の解析結果の一覧表である。 第２地点における試験期間に対する解析結果を示すグラフである。 第３地点における試験期間に対する熱伝導率の変化を示すグラフである。 第３地点における試験期間に対するｆ／Ｎ値の変化を示すグラフである。 第３地点における各試験期間の解析結果の一覧表である。 第３地点における試験期間に対する解析結果を示すグラフである。 第１帯水層の揚水試験条件の一覧表である。 第２帯水層の揚水試験条件の一覧表である。 第１帯水層の揚水試験結果を示すグラフである。 第１帯水層についての逆解析による同定結果の一覧表である。 第１帯水層についての計算値と実測値の比較を示すグラフである。 第２帯水層の揚水試験結果を示すグラフである。 第２帯水層についての逆解析による同定結果の一覧表である。 第２帯水層についての計算値と実測値の比較を示すグラフである。

以下に、図面を参照して本発明の実施の形態を説明する。

［共役傾斜法を用いたパラメータの同定プログラムの検証］
本発明者等は、熱応答試験の解析法として採用した共役傾斜法によるパラメータ同定の手法を実行するために開発した同定プログラムをテストケースに適用し、当該同定プログラムの妥当性を検証した。

（検証手順）
開発した共役傾斜法を用いたパラメータの同定プログラムの妥当性の検証は以下の手順で実施した。
１）図３．１に一覧表として示す２つのテストケースにおいて、２つのパラメータα１、α２に対して予め任意の値を与え、計算ソフト「Ｅｘｃｅｌ」で準解析解（式（２４）〜（２６）、第３０項まで考慮）を用いて計算する。時間当たりの温度変化量を検証データに用いる。
２）同定パラメータとするα１（＝ｂｒ^２（ρｃ））とα２（＝ｂλ）について、予め与えた任意の値以外の値を初期値として与える。
３）逆解析によりパラメータα１、α２を同定する。
４）ＤＩＦ≦０．００００００１を計算終了条件とする。
（ＤＩＦ＝Ｆ_ｉ−Ｆ_ｉ−１：ｉ番目のＦ値とｉ−１番目のＦ値の差分）

（パラメータの評価）
＜テストケース１＞
図３．２は探索回数におけるＦ値、ＤＩＦ値を示すグラフであり、図３．３は探索回数におけるｂλ、ｂｒ^２（ρｃ）を示すグラフであり、図３．４は逆解析より求めたパラメータα１、α２を使った計算値と設定値との比較を示すグラフである。また、図３．５は解析結果の一覧表である。

Ｆ値が最小値となるように探索を繰り返し、これ以上値が変化しなくなったところで探索を終了する。またＤＩＦ値は１つ前のＦ値との差分を表しており、Ｆ値の値が変わらなくなったところで限りなく零に近い値を示す。

テストケース１では、設定値に対して初期推定値を２倍の大きさにして解析を行った。図３．３より、探索を繰り返して、パラメータの値が設定値に近づいていくことが分かり、図３．４より、計算値と設定値は非常に良く一致していることが分かる。また、図３．５に示すようにｆ／Ｎ値(データ１つ当たりの誤差)もほぼ零に等しく、精度の高い結果を得ることが出来たことが分かる。

＜テストケース２＞
図３．６は探索回数におけるＦ値、ＤＩＦ値を示すグラフであり、図３．７は探索回数におけるα２（＝ｂλ）、α１（＝ｂｒ^２（ρｃ））を示すグラフであり、図３．８は逆解析より求めたパラメータα１、α２を使った計算値と設定値との比較を示すグラフである。また、図３．９は解析結果の一覧表である。

テストケース２では設定値に対する初期推定値を１／２の大きさにして解析を行った。図３．７より、探索を繰り返して、パラメータα１、α２の値が設定値に近づいていくことが分かる。図３．８より、テストケース２においてもテストケース１と同様に設定値と計算値は非常に良く一致していることが分かる。また、図３．９に示すように、ｆ／Ｎ値(データ１つ当たりの誤差)もほとんど零に等しく、精度の高い結果が得られたことが分かる。

設定値に対する初期設定の値や有効深度などを変えた２つのテストケースによって検証したが、双方のテストケース共に非常に高い精度の結果が得られ、本発明の同定プログラムの妥当性を検証することが出来た。また、準解析解を用いることで、初期段階の結果も再現できることが確認できた。

１．原位置熱応答試験
本発明者等は、３地点で実施した熱応答試験について、本発明による共役傾斜法を用いたパラメータの同定プログラムを用いて解析を行った。以下に、使用した熱応答試験装置の構成および試験の実施手順を説明し、次に、熱応答試験を行った３箇所での現場サイトの地質条件、試験条件および解析結果を説明する。

［熱応答試験装置の構成］
図４は、熱応答試験の実測値から地盤の熱伝導率を算出するための熱応答試験装置の一例を示す全体構成図である。熱応答試験装置１は、主にクローズド型方式における垂直型地中熱交換井Ａを対象とした公知の構成のものであり、熱媒体を加熱し、循環を行うための循環系２と、熱媒体温度・消費電力・循環流量などを計測する計測系３より構成される。循環系２の主な構成要素は循環ポンプ４、電気ヒーター５、水タンク６などで、計測系３は温度センサー７、８、９、流量計１０、データロガー１１、パーソナルコンピュータ１２などである。

この構成の熱応答試験装置１では、地上部の水タンク６内で電気ヒーター５によって熱媒体を加熱し、循環ポンプ４により熱媒体を配管１３内に循環させる。加熱された熱媒体は、試験対象の地盤１４に形成されているボアホール孔１５を通過する際に放熱するため、Ｕ字管１６の入口・出口の温度差を温度センサー８、９でそれぞれ測定し、その地層の熱交換特性を評価する。

計測計３のパーソナルコンピュータ１２には、キーボード、マウスなどの入力装置１７と、液晶表示器などのモニター、プリンターなどの出力装置１８が接続されている。パーソナルコンピュータ１２には、本発明による熱応答試験の解析プログラム（以下の説明では、これを「同定プログラム」と呼ぶ場合もある。）がインストールされている。解析プログラムを実行することにより、データロガー１１から取り込んだ計測データに基づき、試験対象の地盤の熱交換特性の評価（地盤の熱伝導率の算出処理等）を行う。

ここで、熱負荷には、電圧の安定した電源に接続した電気ヒーター５を使用する。電気ヒーター５は、ＧｅｏＨＰシステム実使用時と同等程度の熱負荷を掛けることのできる容量を有する必要がある。また、電気ヒーター５には安全装置としてサーモスタットと温度ヒューズを設け、加熱事故を防止する。循環ポンプ４は実使用時と同等の流量を循環可能な能力を有し、流量を調整するバルブかインバータを設置する。流量計１０は±０．５％以下の精度を持つ電磁流量計が望ましい。また、熱交換井出入口温度を測定する温度センサー８、９には白金側温抵抗体を用いる。精度は想定される温度範囲において±０．１℃以下が望ましい。

地上部の配管１３の断熱が不十分だと負荷が安定しないため、地上部の配管１３の長さは最小限にし、保温材を巻いて十分な断熱と防水を施す。配管頂部には空気抜き弁１９を設置し、熱媒体の熱膨張による漏水を防ぐために、膨張タンク２０を設ける。さらに、気液分離器２１により、循環サイクル内で入り混じった気体と液体を分離する。データロガー１１に、熱交換井出入口温度、外気温度、流量のデータを記録させる。

［熱応答試験の実施手順］
熱応答試験装置１を用いた、垂直型Ｕ字管式地中熱交換井における熱応答試験の実施手順を以下に示す。
（１）熱交換井掘削・仕上げ後は最低３日間熱交換井を放置し、地温を自然状態へ回復させる。試験を行う熱交換井周辺において、他の掘削工事が行なわれていないことを確認し、近隣での地下水のくみ上げが行われている場合にはそれを休止しておく。
（２）熱応答試験装置１を設計する前に、±０．１℃も感度を持つ温度計を用いて坑位内温度（自然地温）を地表から坑底まで測定する。測定深度感覚は１０ｍ以下が望ましい。（３）断熱を施したホース等で熱応答試験装置１と熱交換井Ａを接続し、配管３内を熱媒体で満たした後、エア抜きをする。熱媒体としては、安価で取り扱いが容易な水が適当であるが、システム完成後に実際に使用する熱媒体を使用してもよい。
（４）温度センサー７〜９、流量計１０をデータロガー１１に接続する。熱媒体の熱交換井出入口温度の測定は熱交換井の出入口にできるだけ近い位置（可能ならば地表面以下）で行うのが望ましい。
（５）熱媒体を充填して数時間放置し、熱媒体温度が地層温度と等しくなった後、ヒーター５への通電、熱媒体の循環を開始する。流速はシステム実使用時と同等の流速とするが、熱媒体の流速が乱流域（Ｕ字管１６の直径長さを代表長さとして計算したレイノルズ数が２，３００以上）となること、温度測定の測定誤差（±０．１℃）を考慮して可能であれば熱交換井出入口での温度差が４℃以上生じるよう設定することが望ましい。
（６）ヒーターによる熱負荷は実際に設置予定のＧｅｏＨＰシステムの負荷に近い大きさとし、加熱期間は標準で６０時間、最短で４８時間が推奨されている。
（７）循環終了後は熱応答試験装置１と地上配管１３の水抜きをし、熱応答試験装置１を撤収する。

［実施例１：第１地点での熱応答試験とパラメータ同定結果］
(１)地質条件
図５．１は、熱応答試験を行った第１地点における深井戸の掘削データを元にした断面図である。第１地点付近の地質は主に砂礫層から構成されており、所々に粘土層・粘土混り砂層砂礫層を挟んでいる。このうち砂礫層が帯水層となり、粘土層・粘土混り砂層砂礫層が不（難）透水層となっている。

熱応答試験を実施したボーリング孔は深度７６ｍまでＵチューブが挿入されており、ほぼ第１帯水層および第２帯水層に相当する。図５．１に示す電気検層結果より、同じ帯水層内にも異なる比抵抗値を示す部分があり、地盤の透水性に不均一性が見られる。

(２)熱応答試験条件
第１地点における熱応答試験は２回実施し、試験実施期間は、１回目が２０１１年３月２５〜２８日であり、２回目が２０１１年１１月１７〜２０日である。図５．２および図５．３に、それぞれの加熱期間、平均加熱量、平均流量、地中熱交換井有効深度、地中熱交換井初期温度、充填剤、調査孔直径、使用したＵチューブを示す。また、図５．４に、使用した測定機器と熱応答試験装置の諸元を示す。第１回目、第２回目とも同様の測定機器を用いて試験を実施している。

(３）測定結果と同定結果
＜第１地点での第１回目（２０１１年３月２５日〜２８日）＞
図５．５は、熱応答試験により得られた測定データ（往き温度、還り温度、平均温度、流量、加熱量、外気温）を示すグラフである。加熱期間は２０１１年３月２５日〜２８日の３日間で、平均加熱量は３５５７（Ｗ）、平均流量は２０．９（Ｌ／ｍｉｎ）、有効深度は７６．０（ｍ）であった。加熱前の地層温度は約１４．１（℃）で安定していたが、試験開始７２時間後には往き温度約２６．５（℃）、還り温度約２３．５（℃）となった。

図５．６は逆解析による同定結果を示す一覧表であり、図５．７は同定結果を用いた準解析解の計算値と現地実測値との比較を示すグラフである。逆解析より得られた熱伝導率の値は２．３３（Ｗ／ｍ／Ｋ）で、この値を用いた計算値と実測値のグラフは非常によく一致している。また逆解析時に使用したデータ数は４３２０（個）で、ｆ／Ｎ値（データ１つ当たりの誤差）は０．０７５（℃）であった。

図５．８は、従来から用いられている作図法によって解析を行った結果を示すグラフである。作図法より得られた熱伝導率の値は２．４０（Ｗ／ｍ／Ｋ）で、逆解析により求めた結果と近い値を示した。解析に用いたデータ区間は１０．３（ｈ）〜７２．０（ｈ）である。

＜第１地点での第２回目（２０１１年１１月１７日〜２０日）＞
図５．９は、熱応答試験により得られた計測データ(往き温度、還り温度、平均温度、流量、加熱量、外気温)を示すグラフである。加熱期間は２０１１年１１月１７日〜２０日の３日間で、平均加熱量は３，９６４（Ｗ）、平均流量は２０．２（Ｌ／ｍｉｎ）、有効深度は７６．０ｍであった。加熱前の地層温度は約１４．１℃で安定していたが、試験開始７２時間後には往き温度約２７．０℃、還り温度約２４．０℃であった。

図５．１０は逆解析による同定結果を示す一覧表であり、図５．１１は同定結果を用いた準解析解の計算値と現地実測値との比較を示すグラフである。第１回目の結果と同様に計算値と実測値は非常によく一致しており、逆解析より得られた熱伝導率の値は２．５２（Ｗ／ｍ／Ｋ）で、第１回目の結果よりも少し大きい値となった。ちなみに逆解析時に使用したデータ数は４３２０（個）で、データ１つ当たりに対する誤差は０．０５６（℃）であった。

図５．１２に作図法を用いて解析を行った結果を示す。作図法より得られた熱伝導率の値は２．５９（Ｗ／ｍ／Ｋ）で、逆解析により求めた結果と近い値を示した。なお、解析に用いたデータ区間は２４．０（ｈ）〜７２．０（ｈ）である。

［実施例２：第２地点での熱応答試験とパラメータ同定結果］
第２地点における熱応答試験は２０１０年１２月２４日〜２９日に実施されものである。以下に、現場サイトの地質条件、熱応答試験条件、同定結果について説明する。

（１）地質条件
熱応答試験を行った第２地点ではコアボーリングによる地質確認は実施されていないため、採掘時の状況および既存資料などから、地質が想定されている。今回の掘削により想定される地質は、粘土や玉石などを混入する砂礫層を主体として構成されているものとみられるが、図６．１に示すように、深度８０ｍ付近以降で温度検層による温度上昇勾配が高くなる現象が認められることから、基盤岩とみられる新第三系の岩盤に変化している可能性が高い。なおＵチューブ設置深度は１４３．７５ｍである。

図６．２は、試験時の加熱期間、平均加熱量、平均流量、地中熱交換井有効深度、地中
熱交換井初期温度、充填剤、調査孔直径、使用したＵチューブを示す一覧表である。なお、使用した測定機器と熱応答試験装置の諸元は、第１地点における場合と同様である（図５．４参照）。

（２）測定結果と同定結果
図６．３は、熱応答試験により得られた計測データ（往き温度、還り温度、平均温度、流量、加熱量、外気温）を示すグラフである。加熱期間は２０１０年１２月２４日〜２７日の３日間で、平均加熱量は５，７１６（Ｗ）、平均流量は１０．４（Ｌ／ｍｉｎ）、有効深度は１４３．７５ｍであった。加熱前の地層温度は約１６．８℃で安定していたが、試験開始７２時間後には往き温度約３４．０℃、還り温度約２６．０℃となった。

図６．４は逆解析による同定結果を示す一覧表であり、図６．５は同定結果を用いた準解析解の計算値と現地実測値との比較を示すグラフである。計算値と実測値はよく一致している。逆解析より得られた熱伝導率の値は１．７９（Ｗ／ｍ／Ｋ）であり、第１地点での２回の試験における同定結果（２．３３ｋＷ、２．５１ｋＷ）と比較すると、小さい値となった。また、ｆ／Ｎ値（データ１つ当たりに対する誤差）は０．１５８（℃）であり、第１地点での同定結果に比べ、大きい結果を示した。これは、熱伝導理論で考慮されていない地下条件もしくは熱移動現象が存在するためだと考えられる。なお、逆解析時に使用したデータ数は４，３２０（個）である。

図６．６は作図法を用いて解析を行った結果を示すグラフである。作図法より得られた熱伝導率の値は１．５３（Ｗ／ｍ／Ｋ）であった。ちなみに解析に用いたデータ区間は１７．２（ｈ）〜７２．０（ｈ）である。

［実施例３：第３地点での熱応答試験とパラメータ同定結果］
第３地点における熱応答試験は、２０１１年１０月１３日〜１６日に実施された。以下に、地質条件、熱応答試験条件、同定結果について説明する。

（１）地質
第３地点は地下水が豊富であること、また地下水流れが存在することが予測される地点である。今回熱応答試験を実施した井戸は深度９９ｍまでＵチューブが挿入されている。また、地質構造については詳細な資料が得られなかったため、ここでは記述しない。

（２）熱応答試験条件
図７．１は、加熱期間、平均加熱量、平均流量、地中熱交換井有効深度、地中熱交換井初期温度、充填剤、調査孔直径、使用したＵチューブを示す一覧表である。なお、使用した測定機器と熱応答試験装置の諸元は、第１地点における場合と同様である（図５．４参照）。

（３）測定結果と同定結果
図７．２は、熱応答試験により得られた計測データ(往き温度、還り温度、平均温度、流量、加熱量、外気温)を示すグラフである。加熱期間は２０１１年１０月１３日〜１６日の３日間で、平均加熱量は４，６６７（Ｗ）平均流量は１０．４（Ｌ／ｍｉｎ）であった。有効深度は９９．０ｍである。加熱前の地層温度は約１７．８℃で安定していたが、試験開始７２時間後には往き温度約２８．８℃、還り温度約２５．８℃となった。

図７．３は逆解析による同定結果を示す一覧表であり、逆解析より得られた熱伝導率の値は３．０１（Ｗ／ｍ／Ｋ）であった。図７．４は同定結果を用いた準解析解の計算値と現地実測値との比較を示すグラフである。計算値と実測値のグラフは、一部において一致していない箇所もあるが、全体としてはほぼ一致している。また、ｆ／Ｎ値（データ１つ
当たりに対する誤差）は０．１８７（℃）で、第１地点での結果と比較すると、大きい結果を示した。これは、熱伝導理論で考慮されていない地下条件もしくは熱移動現象が存在するためだと考えられる。なお、今回の逆解析時に使用したデータ数は４，３３８（個）である。

［試験結果についての考察］
（１）解析精度および体積熱容量の算出
第１地点における２回の熱応答試験と、その結果の逆解析法による解析結果からは、非常に誤差が小さく、精度の高い解析結果が得られていることが分かる。また、従来は求めることのできなかった体積熱容量も、ｒ（井戸からの距離）が分かれば求められる。

（２）地下水位が熱伝導率に与える影響
図８．１は、第１地点における２回の熱応答試験の解析結果を示す一覧表である。それぞれの結果を比較すると、第１回目では熱伝導率が２．３３（Ｗ／ｍ／Ｋ）であるのに対して、第２回目では熱伝導率が２．５１（Ｗ／ｍ／Ｋ）とやや大きくなっている。２回共に同井戸で試験を行っているため、その違いについて考察する。両試験日の地下水位は季節変動により、第１回目の方が第２回目の場合よりも１．５（ｍ）の地下水位が深かったことが確かめられている。

文献値によると、水の熱伝導率と空気の熱伝導率はそれぞれ０．５９４（Ｗ／ｍ／Ｋ）、０．０２４１（Ｗ／ｍ／Ｋ）で、水の熱伝導率の方が高いため、地下水位が高くなると熱伝導率は大きくなる。このため、飽和部が厚い第２回目の方が第１回目より大きい熱伝導率となったと考えられる。

（３）地温勾配が熱伝導率に与える影響
図８．２は、第１地点における第２回目の解析結果と第２地点の解析結果を示す一覧表である。両者を比較すると、第２地点の熱伝導率は１．７９（Ｗ／ｍ／Ｋ）、第１地点の第２回目の熱伝導率は２．５１（Ｗ／ｍ／Ｋ）で、第２地点の熱伝導率の方が小さい値となった。また、ｆ／Ｎ値（データ１つ当たりの誤差）は第２地点が０．１５８（℃）、第１地点の第２回目は０．０５６（℃）で、第２地点の方がデータ１つ当たり約３倍の誤差となった。

そこで、両者の試験開始前の地温分布を調べ、図８．３にその結果を示す。図８．３より、第１地点の地温分布は深度によらずほぼ一定であるのに対して、第２地点の地温分布は不均一で、深度６０ｍ付近から深度が増すと温度が上昇している。深度２０ｍと１４０ｍの地温を比べると、２０ｍ地点の約１５℃に対して１４０ｍ地点では約２０℃となり、５℃の温度差がある。そのため、深度が増すごとに地中熱交換量が小さくなる。また、地温勾配は０．０３（℃／ｍ）が一般的であるが、第２地点の地温勾配は０．０５３（℃／ｍ）で少し大きい傾向にある。

第２地点のように地温勾配がある場合は、深度が増すごとに地温が上昇し、地中熱交換量が小さくなる。そのため、不均一な地温分布を示す第２地点の方が第１地点より小さい熱伝導率を示したと考えられる。線源理論の深度方向に対して地中熱交換量は一定であるという前提条件に反しているため、ｆ／Ｎ値（データ１つ当たりの誤差）が第１地点の結果と比べて、大きくなったことが考えられる。

（４）地下水流動が熱伝導率に与える影響
図８．４は、第１地点の第２回目の解析結果と第３地点の解析結果を示す一覧表である。両者を比較すると、第３地点の熱伝導率は３．０１（Ｗ／ｍ／Ｋ）、第１地点の熱伝導率は２．５１（Ｗ／ｍ／Ｋ）で、第３地点の熱伝導率の方が大きい値を示した。また、ｆ
／Ｎ値（データ１つ当たりの誤差）は、第３地点が０．１８６（℃）、第１地点が０．０５６（℃）で、第３地点の方が理論解と実測値の乖離大きくなった。

第３地点の地形条件を見ると、周辺に２つの川の合流地点があり、豊富な地下水と地下水流動が存在することが考えられる。

図８．５と図８．６は、第１地点と第３地点において、解析に用いるデータ区間を、全区間、後半部、初期段階の３段階に変化させた場合の解析結果を示すグラフである。

第１地点の場合は、どのデータ区間で解析を行っても、熱伝導率とｆ／Ｎ値（データ１つ当たりの誤差）にあまり違いは見られない。それに対して第３地点の場合は、解析を行うデータ区間によって、熱伝導率とｆ／Ｎ値（データ１つ当たりの誤差）に明確な差が表れている。初期段階のデータのみを用いた解析結果と、後半部分のデータを用いた解析結果を比較すると、熱伝導率に１．３１（Ｗ／ｍ／Ｋ）の差があり、ｆ／Ｎ値（データ１つ当たりの誤差）は、それぞれ０．０９（℃）と０．１４４（℃）で、初期段階のデータを用いた場合の方がｆ／Ｎ値（データ１つ当たりの誤差）は小さく、熱伝導理論により合致していることが分かる。

図８．６のようにデータ区間によって異なる結果が得られた理由として、地下水流動による影響が考えられる。熱伝導理論では、熱の移動は全て熱伝導によって行われると仮定されているため、地下水流動に伴う移流の影響が大きくなるほど、実測値と計算値の差は大きくなってしまう。第３地点では、地形的条件から地下水流動が存在することが考えられるが、熱伝導が支配的な初期段階のデータを用いた場合はｆ／Ｎ値（データ１つ当たりの誤差）が小さく、地下水流動の影響が大きくなる後半部分のデータを用いた場合はｆ／Ｎ値（データ１つ当たりの誤差）が大きくなることが分かる。

（５）解析データ区間による熱伝導率への影響
考察（４）で述べたように、今回の解析では、熱の移動は全て熱伝導によって行われると仮定しているため、熱伝導以外の移動がある場合、計算値と実測値の差は大きくなってしまう。しかしながら、熱伝導率以外の移動がある場合でも、温度が立ち上がり始めた初期段階では温度勾配が大きいため、熱伝導による移動が支配的である。従って、初期段階のデータを用いて解析を行うことで、真の熱伝導率を求めることができる。

図８．７、図８．８および図８．９に、それぞれ、第１地点の第２回目、第２地点および第３地点において、それぞれ初期段階のデータを用いた解析結果と後半部のデータを用いた解析結果を示す。

第１地点の第２回目の値は、解析データ区間が変化しても熱伝導率とｆ／Ｎ値（データ１つ当たりの誤差）に対する影響は見られない。一方、第２地点および第３地点は、後半部のデータを用いて解析を行った場合、それぞれ地温勾配と地下水流動による影響を受けている。かつ、初期段階のデータを用いた解析結果では、３ケースとも近い値の熱伝導率を示しており、ｆ／Ｎ値（データ１つ当たりの誤差）も小さい結果が得られた。

以上の結果より、地温勾配や地下水流動に影響されない真の熱伝導率を求めたい場合には、熱伝導が支配的な初期段階のデータを用いる方が良いことが分かる。

［試験期間を短縮した場合の解析結果の精度の検討］
本発明の逆解析による解析法では、準解析解を用いることで、従来の近似解では評価できなかった初期段階の試験結果が活用できる。このため、従来の方法よりも短い試験期間で良好な解析結果を得られる可能性がある。ここでは、逆解析法で使用していたデータ数
を段階的に減らして解析を行った場合の解析結果の精度の検討を行う。

（１）第１地点の第１回目
試験期間を０．５日刻みで、３．０日，２．５日，２．０日，１．５日，１．０日と段階的に減らした場合の熱伝導率とｆ／Ｎ値（データ１つ当たりの誤差）をそれぞれ図９．１と図９．２に示す。図９．３は解析結果のまとめである。図９．３より、試験期間を短縮した場合でも熱伝導率λ（Ｗ／ｍ／Ｋ）に大きな違いは見られず、ｆ／Ｎ値（データ１つ当たりの誤差）も０．０７７〜０．０７１（℃）と小さく、精度の高い解析結果が得られた。

図９．４は同定結果を用いた準解析解の計算値と現地実測値との比較である。計算値と現地実測値は非常に良く一致しており、準解析を用いることで従来の近似解では評価できなかった初期段階の試験結果を活用できた。

（２）第１地点の第２回目
試験期間を０．５日刻みで、３．０日，２．５日，２．０日，１．５日，１．０日と段階的に減らした場合の熱伝導率とｆ／Ｎ値（データ１つ当たりの誤差）をそれぞれ図９．５と図９．６に示す。図９．７は解析結果のまとめである。

図９．７より、試験期間を減らしても得られた熱伝導率λ（Ｗ／ｍ／Ｋ）に大きな違いは見られない。試験期間３．０日間と１．０日間の熱伝導率を比較すると、それぞれ２．５１（Ｗ／ｍ／Ｋ）、２．５８（Ｗ／ｍ／Ｋ）で、その差は１．０２８％と小さい。また、ｆ／Ｎ値（データ１つ当たりの誤差）は０．０６０（℃）〜０．０５３（℃）と小さく、試験期間を短縮した場合でも、精度の高い解析結果が得られた。

図９．８は同定結果を用いた準解析解の計算値と現地実測値との比較である。第１地点での第１回目の結果と同様に計算値と現地実測値は非常に良く一致しており、準解析を用いることで従来の近似解では評価できなかった初期段階の試験結果を活用できた。

（３）第２地点（２０１０年１２月２４〜２７日）
第２地点で実施された別応答試験において、試験期間を０．５日刻みで３．０日，２．５日，２．０日，１．５日，１．０日と段階的に減らした場合の熱伝導率とｆ／Ｎ値（データ１つ当たりの誤差）をそれぞれ図９．９と図９．１０に示す。図９．１１は解析結果のまとめである。

図９．１１より、同定された熱伝導率の値は１．７９（Ｗ／ｍ／Ｋ）〜２．０９（Ｗ／ｍ／Ｋ）で、試験期間を長くするにつれて熱伝導率の値が小さくなっている。試験が実施された地層では温度検層の結果より、地温勾配が確認されている。そのため、ボーリング孔下部で熱が供給され、見かけ熱伝導率が小さくなると考えられる。

図９．９より、熱伝導率の値は安定していないことが分かる。そのため、試験期間を延長すると、見かけ熱伝導率はさらに小さくなると予測される。また、ｆ／Ｎ値（データ１つ当たりの誤差）は０．１１２（℃）〜０．１６８（℃）であった。これは、試験結果についての考察（４）で述べたように、初期段階のデータでは熱伝導が支配的であることによる。

図９．１２は同定結果を用いた準解析解の計算値と現地実測値との比較である。初期段階のデータを用いた方がｆ／Ｎ値（データ１つ当たりの誤差）は小さくなり、計算値と現地実測値は良く一致している。準解析を用いることで従来の近似解では評価できなかった初期段階の試験結果を活用できた。

（４）第３地点（２０１１年１０月１３〜１６日）
第３地点で実施された熱応答試験において、試験期間を０．５日刻みで３．０日，２．５日，２．０日，１．５日，１．０日と段階的に減らした場合の熱伝導率とｆ／Ｎ値（データ１つ当たりの誤差）をそれぞれ図９．１３と図９．１４に示す。図９．１５は解析結果のまとめである。

図９．１５より、熱伝導率の値は２．４２（Ｗ／ｍ／Ｋ）〜３．０１（Ｗ／ｍ／Ｋ）で、試験期間を長くするにつれて熱伝導率の値が大きくなっている。ｆ／Ｎ値（データ１つ当たりの誤差）は０．１０７（℃）〜０．１８７（℃）で、試験期間が長いほど大きくなっている。これは、地下水流動の影響だと考えられる。温度が立ち上がり始めた初期段階では温度勾配が大きいため熱伝導が支配的になるが、後半部では移流が支配的となり、ｆ／Ｎ値（データ１つ当たりの誤差）は大きく、見かけ熱伝導率の値も大きくなったと考えられる。

図９．１３より、熱伝導率の値はまだ安定しておらず、試験期間を延長すると、見かけの熱伝導率はさらに大きくなると予測され、地下水流動がある場合は試験期間が多いほど真の熱伝導率と乖離していく。

図９．１６は同定結果を用いた準解析解の計算値と現地実測値との比較である。初期段階のデータを用いた方がｆ／Ｎ値（データ１つ当たりの誤差）は小さくなり、計算値と現地実測値は良く一致している。準解析を用いることで従来の近似解では評価できなかった初期段階の試験結果を活用できた。

２．揚水試験解析
次に、熱応答試験が行われた第１地点において、２０１１年１月１８日に行われた揚水試験に本発明の同定プログラムを適応した場合の結果について考察する。

（１）試験概要
第１帯水層、第２帯水層、それぞれの井戸能力を把握するため、層別揚水試験を行った。実施日は２０１１年１月１８日である。図１０．１と図１０．２に、各帯水層における揚水時間、揚水量、帯水層厚、揚水井と観測井との距離について示す。

（２）揚水試験結果と同定結果
ａ）第１帯水層
第１帯水層での揚水試験結果を図１０．３に示す。揚水井から３．８ｍの距離にある観測井の地下水位の変動を表している。地下水位は１ｍｉｎ後から減少を続け、６５ｍｉｎで、０．２４６ｍの水位降下が確認されている。６５ｍｉｎ以降は水位の変化がほとんど見られなくなったため、１２０ｍｉｎで試験を終了している。

図１０．４は逆解析による同定結果である。図１０．５は同定結果を用いた準解析解の計算値と現地実測値との比較を示している。この図より計算値と現地実測値は一致して見えるが、図１０．４よりｆ／Ｎ値（データ１つ当たりに対する誤差）は０．０４５（ｍ）で少し大きい値を示した。逆解析で得られた透水量係数は１．３９（ｍ^２／ｍｉｎ）となった。また、透水量係数１．３９（ｍ^２／ｍｉｎ）を帯水層厚３２．７（ｍ）で除して、透水係数は０．０４２（ｍ／ｍｉｎ）となった。この値は一般的な透水係数の値よりも大きい。これらから、第１帯水層は漏水していることが考えられる。帯水層が漏水している場合は、Ｔｈｅｉｓの式ではなく、Ｈａｎｔｕｓｈ−Ｊａｃｏｂの式が一般的に用いられているため、ｆ／Ｎ値（データ１つ当たりに対する誤差）も大きくなったと考えられる。

なお、今回の揚水試験では、揚水井と観測井の距離が分っていたので、貯留係数も求めることができた。逆解析時に使用したデータ数は２９（個）である。

ｂ）第２帯水層
第２帯水層での揚水試験結果を図１０．６に示す。揚水井から３．８ｍの距離にある観測井の地下水位の変動を表している。地下水位は１ｍｉｎ後から減少を続け、９０ｍｉｎで、０．２４０ｍの水位降下が確認されている。９０ｍｉｎ以降は水位の変化がほとんど見られなくなったため、１６０ｍｉｎで試験を終了している。

図１０．７は逆解析による同定結果であり、図１０．８は同定結果を用いた準解析解の計算値と現地実測値との比較を示している。図１０．８より計算値と現地実測値は一致して見えるが、第１帯水層と同様にｆ／Ｎ値（データ１つ当たりに対する誤差）は０．０６９（ｍ）で大きい値を示した。逆解析で得られた透水量係数は０．９２（ｍ^２／ｍｉｎ）となった。また、透水量係数０．９２（ｍ^２／ｍｉｎ）を帯水層厚２１．２（ｍ）で除して、透水係数は０．０４３（ｍ／ｍｉｎ）となった。この値は一般的な透水係数の値よりも大きい。これらから、第２帯水層は漏水していることが考えられる。帯水層が漏水している場合は、Ｔｈｅｉｓの式ではなく、Ｈａｎｔｕｓｈ−Ｊａｃｏｂの式が一般的に用いられているため、ｆ／Ｎ値（データ１つ当たりに対する誤差）も大きくなったと考えられる。なお逆解析時に使用したデータ数は３４（個）である。

１ 熱応答試験装置
２ 循環系
３ 計測系
４ 循環ポンプ
５ 電気ヒーター（加熱用ヒーター）
６ 水タンク
７〜９ 温度センサー
１０ 流量計
１１ データロガー
１２ パーソナルコンピュータ

Claims (12)

1. 熱応答試験によって測定された調査対象の地盤における温度上昇量の経時変化に基づき前記地盤の熱交換特性を評価する熱応答試験の解析方法であって、
   地盤の温度上昇量を式（Ａ）で表されるＫｅｌｖｉｎの線源関数で規定し、当該線源関数の近似式として、指数積分Ｅ（Ｘ）のべき級数展開式における少なくとも第１６項まで採用した近似式を使用し、
   パウエルの共役傾斜法に基づく非線形最適化法を用いて、前記近似式に含まれる未知のパラメータを所定の値に設定した場合に得られる温度上昇量の計算値Ｔｒと、前記熱応答試験によって測定された温度上昇量の実測値Ｔｒｉとの誤差が最小となるように、前記パラメータの逆解析を行い、
   当該逆解析においては、同定対象の前記パラメータを、前記線源関数における未知量である半径ｒ、密度ρ、比熱ｃおよび熱伝導率λを用いて、次式（Ｂ）、（Ｃ）で表されるα１およびα２とし、
   前記逆解析によって同定された前記パラメータα１およびα２を用いて、前記地盤の熱伝導率λおよび体積熱容量ρｃのうち、少なくとも熱伝導率λを算出することを特徴とする熱応答試験の解析方法。
2. 請求項１において、
   Ｋｅｌｖｉｎの線源関数の近似式として、指数積分Ｅ（Ｘ）のべき級数展開式における第３０項まで採用した近似式を用いることを特徴とする熱応答試験の解析方法。
3. 請求項１または２において、
   同一地点について異なる時点において行われた少なくとも２回の熱応答試験の実測値を用いて、それぞれの時点における同一地点における地層の熱伝導率を算出し、
   それぞれの時点における地下水位の変化の有無を判断するための指標として、算出された前記熱伝導率の差を求めることを特徴とする熱応答試験の解析方法。
4. 請求項１ないし３のうちのいずれかの項において、
   前記熱応答試験によって得られた実測値を、試験開始時点から所定時間経過した時点までの初期段階データと、当該時点以後の後半部データに分け、
   前記初期段階データを用いて熱伝導率および、式（Ｄ）により表されるデータ１つ当りの誤差であるｆ／Ｎ値を算出すると共に、前記後半部データを用いて熱伝導率および、前記ｆ／Ｎ値を算出し、
   地層深度方向の地熱勾配の有無、あるいは、地下水流動の有無を判断するための指標として、算出した２つの熱伝導率の差、および算出した２つのｆ／Ｎ値の差をそれぞれ求めることを特徴とする熱応答試験の解析方法。
5. 熱応答試験によって測定された調査対象の地盤における温度上昇量の経時変化に基づき前記地盤の熱交換特性を評価するために用いる熱応答試験の解析プログラムであって、
   コンピュータを、地盤の温度上昇量を、式（Ａ）で表されるＫｅｌｖｉｎの線源関数の近似式として、指数積分Ｅ（Ｘ）のべき級数展開式における少なくとも第１６項まで採用した近似式を用いて演算を行う演算手段として機能させ、
   コンピュータを、パウエルの共役傾斜法に基づく非線形最適化法を用いて、前記近似式に含まれる未知のパラメータを所定の値に設定した場合に得られる温度上昇量の計算値Ｔｒと、前記熱応答試験によって測定された温度上昇量の実測値Ｔｒｉとの誤差が最小となるように、前記パラメータの逆解析を行う逆解析手段として機能させ、
   当該逆解析手段においては、同定対象の前記パラメータを、前記線源関数における未知量である半径ｒ、密度ρ、比熱ｃおよび熱伝導率λを用いて、式（Ｂ）、（Ｃ）で表されるα１およびα２とし、
   さらに、コンピュータを、前記逆解析によって同定された前記パラメータα１およびα２を用いて、前記地盤の熱伝導率λおよび体積熱容量ρｃのうち、少なくとも熱伝導率λを算出する算出手段として機能させることを特徴とする熱応答試験の解析プログラム。
6. 請求項５において、
   Ｋｅｌｖｉｎの線源関数の前記近似式として、指数積分Ｅ（Ｘ）のべき級数展開式における第３０項まで採用した近似式を用いることを特徴とする熱応答試験の解析プログラム
   。
7. 請求項５または６において、
   前記逆解析手段および前記算出手段によって、同一地点について異なる時点において行われた少なくとも２回の熱応答試験の実測値を用いて、それぞれの時点における同一地点における地層の熱伝導率を算出し、
   さらに、コンピュータを、それぞれの時点における地下水位の変化の有無を判断するために、算出された熱伝導率の差を算出する熱伝導率差算出手段として機能させることを特徴とする熱応答試験の解析プログラム。
8. 請求項５ないし７のうちのいずれかの項において、
   前記逆解析手段は、前記熱応答試験によって得られた実測値を、試験開始時点から所定時間経過した時点までの初期段階データと、当該時点以後の後半部データに分け、前記初期段階データを用いて熱伝導率および、式（Ｄ）により表されるデータ１つ当りの誤差であるｆ／Ｎ値を算出すると共に、前記後半部データを用いて熱伝導率および、前記ｆ／Ｎ値を算出し、
   さらに、コンピュータを、地層深度方向の地熱勾配の有無、あるいは、地下水流動の有無を判断するために、算出した２つの熱伝導率の差、および算出した２つのｆ／Ｎ値の差をそれぞれ求める差分算出手段として機能させることを特徴とする熱応答試験の解析プログラム。
9. 揚水試験によって測定された調査対象の地盤における地下水位低下量の経時変化に基づき前記地盤の地下水流動特性を評価する揚水試験の解析方法であって、
   地盤の地下水位低下量を式（Ｅ）で表されるタイスの井戸関数で規定し、当該井戸関数の近似式として、指数積分Ｗ（ｕ）のべき級数展開式における少なくとも第１６項まで採用した近似式を使用し、
   パウエルの共役傾斜法に基づく非線形最適化法を用いて、前記近似式に含まれる未知のパラメータを所定の値にそれぞれ設定した場合に得られる地下水位低下量の計算値と、前記揚水試験によって測定された地下水位低下量の実測値との誤差が最小となるように、前記パラメータの逆解析を行い、
   当該逆解析においては、同定対象の前記パラメータを、前記井戸関数における未知量である半径ｒ、透水量係数Ｔおよび貯留係数Ｓｃを用いて、次式で表されるα１およびα２とし、
   α１＝ｒ^２Ｓｃ
   α２＝Ｔ
   前記逆解析によって同定された前記パラメータα１およびα２から、前記地盤の透水量係数Ｔおよび貯留係数Ｓｃのうち、少なくとも透水量係数Ｔを算出することを特徴とする揚水試験の解析方法。
10. 請求項９において、
    タイスの井戸関数の前記近似式として、指数積分Ｗ（ｕ）のべき級数展開式における第３０項まで採用した近似式を用いることを特徴とする揚水試験の解析方法。
11. 揚水試験によって測定された調査対象の地盤における地下水位低下量の経時変化に基づき前記地盤の地下水流動特性を評価するために用いる揚水試験の解析プログラムであって、
    コンピュータを、地盤の地下水位低下量を式（Ｅ）で表されるタイスの井戸関数の前記近似式として、指数積分Ｗ（ｕ）のべき級数展開式における少なくとも第１６項まで採用した近似式を用いて演算を行う演算手段として機能させ、
    コンピュータを、パウエルの共役傾斜法に基づく非線形最適化法を用いて、前記近似式に含まれる未知のパラメータを所定の値にそれぞれ設定した場合に得られる地下水位低下量の計算値と、前記揚水試験によって測定された地下水位低下量の実測値との誤差が最小となるように、前記パラメータの逆解析を行う逆解析手段として機能させ、
    当該逆解析手段においては、同定対象の前記パラメータを、前記井戸関数における未知量である半径ｒ、透水量係数Ｔおよび貯留係数Ｓｃを用いて、次式で表されるα１およびα２とし、
    α１＝ｒ^２Ｓｃ
    α２＝Ｔ
    さらに、コンピュータを、前記逆解析によって同定された前記パラメータα１およびα２から、前記地盤の透水量係数Ｔおよび貯留係数Ｓｃのうち、少なくとも透水量係数Ｔを算出する算出手段として機能させることを特徴とする熱応答試験の解析プログラム。
12. 請求項１１において、
    タイスの井戸関数の前記近似式として、指数積分Ｗ（ｕ）のべき級数展開式における第３０項まで採用した近似式を用いることを特徴とする揚水試験の解析プログラム。