论文部分内容阅读
【摘要】从岩体的热传导微分方程出发，考虑了岩体温度与岩体裂隙中水流温度存在温度差时，给出了裂隙岩体温度场的控制方程与定解条件，并推导了基于三角形单元的有限元法计算格式，编写相应程序计算岩体温度场中各节点的温度值，最后通过工程算例进行验证分析。研究表明：对于裂隙岩体而言，考虑裂隙水流温度与裂隙岩体温度之间的热量交换，更加符合实际工程的情况；裂隙岩体与裂隙水流之间的热量交换从二者的接触面开始，岩体内部的温度变化存在一定的滞后效应，且二者的初始温差越大，热量交换越快；若假定研究区域与外界隔热，且其中裂隙水流温度始终保持不变，则无论初始裂隙岩体的温度分布如何，当最终达到稳定时，裂隙岩体的温度分布基本一致。
【关键词】裂隙岩体；温度场；有限单元法；数值分析
1、引言
在大多数基础设施建设和对深层地下水资源的开采利用中，已越来越多地涉及到裂隙岩体[1～3]。因此在对裂隙岩体的研究过程中，渗流场与温度场的耦合分析已成为一项非常重要的研究内容，尤其是在岩体地下水开采、石油开发、地热开采和核废料处理等领域的作用则更为明显[4～7]。目前，对于裂隙岩体大都采用等效连续介质[8～9]来进行处理，因为它有非常好的理论基础和实践经验。然而对于一些断层或大的稀疏裂隙，采用等效连续介质处理就不能完全反映岩体的实际情况，因此采用裂隙网络法来进行处理就成为了一种必然趋势。于是假设水流仅在岩体裂隙中流动时，如何求解裂隙中水温变化情况下的岩体温度分布就成为一个核心问题[1，10]。
本文就是对裂隙岩体采用裂隙网络法来进行处理，假设水流仅在裂隙中流动，从岩体的热传导方程出发，借助有限单元法，并将裂隙水流的温度分布作为岩体温度场计算的边界条件，对裂隙岩体二维非稳定温度场的分布状况进行一些初步探讨，以便探索岩体裂隙网络中温度场的分布，从而为裂隙岩体渗流场与温度场的耦合分析奠定一定的基础。
2、裂隙岩体热传导微分方程和定解条件
热传导微分方程建立了温度与时间、空间的内在联系，由热传导理论得岩体的二维热传导方程为[11]
（1）
式中：T为岩体温度，℃；a为导温系数，可用来表示，其中λ为导热系数，；c为岩体比热，KJ/（kg·℃）；p为岩体容重，kg/m3；t为时间，s；w为单位时间内单位体积放出的热量。
满足热传导方程的解有无限多个，为了求解满足工程需要的温度场，还必须确定出相应的初始条件和边界条件。初始条件，即在初始瞬间时刻固体内部的温度场分布规律，可用下式来表示。当t=0时：T=T0（x，y）边界条件，为固体表面与周围介质相互作用的规律，可分为以下三类：第一类边界条件（给定温度边界条件）：
；
第二类边界条件（给定热流量边界条件）：
；
第三类边界条件：固体表面与流体接触时，通过固体表面的热流密度和固体表面温度T与流体表面温度T?之差成正比，即在边界Γ3上：
式中：β为表面放热系数，；lx、ly为边界表面向外法线的方向余弦。
3、裂隙岩体二维非稳定温度场求解
对（1）式假定w=0，根据变分原理，温度场的求解可化为泛函的极值问题[12]
（2）
式中：Γ3为满足第三类边界条件的边界；R为温度场研究区域。显然只有在与边界Γ3有接触的单元才有上式中的第二项。推导得其有限元求解格式。
（3）
式中：[Q]中的元素qij只有当有限元网格的节点i，j都在Γ3上时才不为零；[p]为已知常数项。对时间项取隐式有限差分，则
（4）
上式即为最终需要求解的线性代数方程组，已知前一时刻的温度便可求出下一时刻的温度。根据上式编写相应程序计算研究区域内各节点的温度值。
4、工程算例
4.1 计算模型与定解条件
以某工程坝址区的30m×30m岩体区域作为研究对象，假设共存在水平、竖直四条裂隙，如图1所示，故研究区域被剖分为9个子区域。假设岩体不透水，地下水流仅存在于裂隙当中。采用有限元法中的三角形单元对研究区域进行离散求解。对每个子区域均剖分为8个三角形单元，共计49个节点，72个三角形单元。共设计了5种计算工况，如表1所示。采用本文中所编写的程序代码对温度场的分布规律进行计算。用到的相关参数取工程中的常用值，见表2。
表1 非稳定温度场计算工况
工况 初始条件 边界条件
1 岩体温度自上边界（20℃）到下边界（2℃）线性分布 裂隙中的水流温度为2℃ 假设研究区域的四个边界均为绝热边界，且随时间不发生变化
2 裂隙中的水流温度为10℃
3 裂隙中的水流温度为20℃
4 岩体温度为20℃，裂隙中的水流温度10℃
5 岩体温度为2 ℃，裂隙中的水流温度20℃
4.2 非稳定温度场的计算
分别对表1中5种工况的非稳定温度场进行计算，采用surfer软件绘制了每种工况下各个时刻的温度场，由于篇幅所限，仅给出了每种工况下第10小时、第1000小时的温度场等值线分布图，见图2～图11所示。
4.3结果分析
（1）通过工况2和4（或3和5）比较分析可以看出，无论岩体的初始温度如何，由于水流温度相同且不发生变化，故在其二者热量交换达到稳定时，岩体温度场的分布基本一致，且均趋于水流温度；
（2）由图2、4、6可以看出，岩体与水流之间的初始温差越大，热量交换愈快，温度等值线图变化愈剧烈；
（3）从图8、10可以看出，岩块与水流之间的热量交换从二者的接触面开始，岩体内部存在一定的滞后效应。
5、结论
本文主要针对岩体裂隙网络中岩体温度与水流温度存在温度差时，对裂隙岩体二维非稳定温度场进行了分析。通过算例分析可以得出如下结论：
（1）裂隙岩体与裂隙水流之间的热量交换从二者的接触面开始，岩体内部的温度变化存在一定的滞后效应，且二者的初始温差越大，热量交换越快。
（2）若假设研究区域与外界隔热，且其中裂隙水流温度始终保持不变，则无论初始裂隙岩体的温度分布如何，当最终达到稳定时，裂隙岩体的温度分布则基本趋于一致。
（3）对于裂隙岩体而言，考虑裂隙水流温度与裂隙岩体温度之间的热量交换，更加符合实际工程的情况。
参考文献
[1]仵彦卿，张倬元.岩体水力学导论[M].成都：西南交通大学出版，1995.
[2]周志芳，王锦国.裂隙介质水动力学[M].北京：中国水利水电出版社，2004.
[3]柴军瑞.大坝工程渗流力学[M].拉萨：西藏人民出版社，2001.
[4]卓宁，孙家庆.工程对流换热[M].北京：机械工业出版社，1982.
[5]J.A.亚当斯，D.F.罗杰斯.传热学计算机分析[M].章靖武，蒋章焰译.北京：科学出版社，1982.
[6]井兰如，冯夏庭.放射性废物地质处置中主要岩石力学问题[J].岩石力学与工程学报，2006，25（4）：833～841.
[7]赵阳升.矿山岩石流体力学[M].北京：煤炭工业出版社，2004.
[8]许增光，李康宏，柴军瑞等.考虑渗流热学效应的大坝稳定温度场有限元数值分析[J].红水河，2006，25（2）：112～115.
[9]朱伯芳.大体积混凝土温度应力与温度控制[M].北京：中国电力出版社，1999.
[10]丁红瑞，王如宾.基于单裂隙岩体二维温度场数值模拟分析[J].灾害与防治工程，2005，59（2）：17～22.
[11]戴锅生.传热学（第二版）[M].北京：高等教育出版社，1999.
[12]朱伯芳.有限单元法原理与应用（第二版）[M].北京：中国水利水电出版社，1998.
【关键词】裂隙岩体；温度场；有限单元法；数值分析
1、引言
在大多数基础设施建设和对深层地下水资源的开采利用中，已越来越多地涉及到裂隙岩体[1～3]。因此在对裂隙岩体的研究过程中，渗流场与温度场的耦合分析已成为一项非常重要的研究内容，尤其是在岩体地下水开采、石油开发、地热开采和核废料处理等领域的作用则更为明显[4～7]。目前，对于裂隙岩体大都采用等效连续介质[8～9]来进行处理，因为它有非常好的理论基础和实践经验。然而对于一些断层或大的稀疏裂隙，采用等效连续介质处理就不能完全反映岩体的实际情况，因此采用裂隙网络法来进行处理就成为了一种必然趋势。于是假设水流仅在岩体裂隙中流动时，如何求解裂隙中水温变化情况下的岩体温度分布就成为一个核心问题[1，10]。
本文就是对裂隙岩体采用裂隙网络法来进行处理，假设水流仅在裂隙中流动，从岩体的热传导方程出发，借助有限单元法，并将裂隙水流的温度分布作为岩体温度场计算的边界条件，对裂隙岩体二维非稳定温度场的分布状况进行一些初步探讨，以便探索岩体裂隙网络中温度场的分布，从而为裂隙岩体渗流场与温度场的耦合分析奠定一定的基础。
2、裂隙岩体热传导微分方程和定解条件
热传导微分方程建立了温度与时间、空间的内在联系，由热传导理论得岩体的二维热传导方程为[11]
（1）
式中：T为岩体温度，℃；a为导温系数，可用来表示，其中λ为导热系数，；c为岩体比热，KJ/（kg·℃）；p为岩体容重，kg/m3；t为时间，s；w为单位时间内单位体积放出的热量。
满足热传导方程的解有无限多个，为了求解满足工程需要的温度场，还必须确定出相应的初始条件和边界条件。初始条件，即在初始瞬间时刻固体内部的温度场分布规律，可用下式来表示。当t=0时：T=T0（x，y）边界条件，为固体表面与周围介质相互作用的规律，可分为以下三类：第一类边界条件（给定温度边界条件）：
；
第二类边界条件（给定热流量边界条件）：
；
第三类边界条件：固体表面与流体接触时，通过固体表面的热流密度和固体表面温度T与流体表面温度T?之差成正比，即在边界Γ3上：
式中：β为表面放热系数，；lx、ly为边界表面向外法线的方向余弦。
3、裂隙岩体二维非稳定温度场求解
对（1）式假定w=0，根据变分原理，温度场的求解可化为泛函的极值问题[12]
（2）
式中：Γ3为满足第三类边界条件的边界；R为温度场研究区域。显然只有在与边界Γ3有接触的单元才有上式中的第二项。推导得其有限元求解格式。
（3）
式中：[Q]中的元素qij只有当有限元网格的节点i，j都在Γ3上时才不为零；[p]为已知常数项。对时间项取隐式有限差分，则
（4）
上式即为最终需要求解的线性代数方程组，已知前一时刻的温度便可求出下一时刻的温度。根据上式编写相应程序计算研究区域内各节点的温度值。
4、工程算例
4.1 计算模型与定解条件
以某工程坝址区的30m×30m岩体区域作为研究对象，假设共存在水平、竖直四条裂隙，如图1所示，故研究区域被剖分为9个子区域。假设岩体不透水，地下水流仅存在于裂隙当中。采用有限元法中的三角形单元对研究区域进行离散求解。对每个子区域均剖分为8个三角形单元，共计49个节点，72个三角形单元。共设计了5种计算工况，如表1所示。采用本文中所编写的程序代码对温度场的分布规律进行计算。用到的相关参数取工程中的常用值，见表2。
表1 非稳定温度场计算工况
工况 初始条件 边界条件
1 岩体温度自上边界（20℃）到下边界（2℃）线性分布 裂隙中的水流温度为2℃ 假设研究区域的四个边界均为绝热边界，且随时间不发生变化
2 裂隙中的水流温度为10℃
3 裂隙中的水流温度为20℃
4 岩体温度为20℃，裂隙中的水流温度10℃
5 岩体温度为2 ℃，裂隙中的水流温度20℃
4.2 非稳定温度场的计算
分别对表1中5种工况的非稳定温度场进行计算，采用surfer软件绘制了每种工况下各个时刻的温度场，由于篇幅所限，仅给出了每种工况下第10小时、第1000小时的温度场等值线分布图，见图2～图11所示。
4.3结果分析
（1）通过工况2和4（或3和5）比较分析可以看出，无论岩体的初始温度如何，由于水流温度相同且不发生变化，故在其二者热量交换达到稳定时，岩体温度场的分布基本一致，且均趋于水流温度；
（2）由图2、4、6可以看出，岩体与水流之间的初始温差越大，热量交换愈快，温度等值线图变化愈剧烈；
（3）从图8、10可以看出，岩块与水流之间的热量交换从二者的接触面开始，岩体内部存在一定的滞后效应。
5、结论
本文主要针对岩体裂隙网络中岩体温度与水流温度存在温度差时，对裂隙岩体二维非稳定温度场进行了分析。通过算例分析可以得出如下结论：
（1）裂隙岩体与裂隙水流之间的热量交换从二者的接触面开始，岩体内部的温度变化存在一定的滞后效应，且二者的初始温差越大，热量交换越快。
（2）若假设研究区域与外界隔热，且其中裂隙水流温度始终保持不变，则无论初始裂隙岩体的温度分布如何，当最终达到稳定时，裂隙岩体的温度分布则基本趋于一致。
（3）对于裂隙岩体而言，考虑裂隙水流温度与裂隙岩体温度之间的热量交换，更加符合实际工程的情况。
参考文献
[1]仵彦卿，张倬元.岩体水力学导论[M].成都：西南交通大学出版，1995.
[2]周志芳，王锦国.裂隙介质水动力学[M].北京：中国水利水电出版社，2004.
[3]柴军瑞.大坝工程渗流力学[M].拉萨：西藏人民出版社，2001.
[4]卓宁，孙家庆.工程对流换热[M].北京：机械工业出版社，1982.
[5]J.A.亚当斯，D.F.罗杰斯.传热学计算机分析[M].章靖武，蒋章焰译.北京：科学出版社，1982.
[6]井兰如，冯夏庭.放射性废物地质处置中主要岩石力学问题[J].岩石力学与工程学报，2006，25（4）：833～841.
[7]赵阳升.矿山岩石流体力学[M].北京：煤炭工业出版社，2004.
[8]许增光，李康宏，柴军瑞等.考虑渗流热学效应的大坝稳定温度场有限元数值分析[J].红水河，2006，25（2）：112～115.
[9]朱伯芳.大体积混凝土温度应力与温度控制[M].北京：中国电力出版社，1999.
[10]丁红瑞，王如宾.基于单裂隙岩体二维温度场数值模拟分析[J].灾害与防治工程，2005，59（2）：17～22.
[11]戴锅生.传热学（第二版）[M].北京：高等教育出版社，1999.
[12]朱伯芳.有限单元法原理与应用（第二版）[M].北京：中国水利水电出版社，1998.