三维地质建模的数学模型与显示方法

2
R2
ln
rij
2
R2
+1−
rij
2
R2
)
式 中 rij2=(xi － xj)2+(yi － yj)2 。 对 N 个 实 测 数 据 点 求 解 此 线 性 方 程 组 可 得 待 定 系 数
A j (j=1,2, … ,n)，将这些值代回插值函数式，即为通过各实测数据点且处处连续光滑。
2.1.4 曲面样条插值法 该方法的优点是原始数据点可以任意排列，拟合出的光滑曲面通过已知数据点。其数学 表达式[3,21]为：
T
" r1,n −12 ln(r1, n −12 + ε ) r1, n 2 ln(r1, n 2 + ε ) 1 x1 2 2 " r2,n −1 ln(r2, n −1 + ε ) r2, n 2 ln(r2, n 2 + ε ) 1 x2 " # # # # " 0 rn −1, n 2 ln(rn −1, n 2 + ε ) 1 xn −1 " rn −1,n 2 ln(rn −1,n 2 + ε ) " 1 " xn −1 " yn −1
j =1
n
求解此线性方程组可获得待定系数 aj (j=1,2,…,n)，将这些值代回插值函数式，即为通过 各实测数据点且处处连续光滑的曲面方程。 在数据点数量不太大的情况下，Multiquadric 法计算不太麻烦。这一方法提出的近 20 年 间，在水文测量、大地测量、地质及采矿、地球物理等领域得到广泛应用，效果良好。 2.1.3 基于平面弹性理论的曲面样条插值法 对地表地形或地下水位等离散数据进行函数拟合可借助平面弹性理论，把计算域视为无 限延伸平面，在平面 (x 1 , y1 ) 点上给定一垂直位移 A 1 时，平面各处将均会随着产生位移，若 则由 A 1 引起距离 (x 1 , y1 ) 点为 r(r<R) 处的点 限制距离 (x 1 , y1 ) 点半径为 R 以外位移为零时， (x,y)的位移 z A (r) 为：
n ⎧ W a a x a y Fij rij 2 ln(rij 2 + ε ) = + + + ∑ 0 1 j 2 j ⎪ j i =1 ⎪ n n n ⎪ Fij = ∑ xi Fi = ∑ yi Fi ⎨ ∑ i =1 i =1 i =1 ⎪ 2 2 ⎪ rij = ( x j − xi ) + ( y j − yi ) 2 ⎪ j = 1, 2,..., n ⎩
［1］
在地质层面模拟
［2 ］
中的应用受到极大限制。另外，现今热门的通过 Voronoi 图和 Dulauny 三角剖分
在空间构
造不规则三角网（Triangular Irregular Net, TIN）方法构造层面，一方面由于原始数据点一般 相距较远，常需要进行三角网的插值加密，另一方面，插值曲面不光滑，无法求出层面上某 点沿坐标轴的坡向、坡度、曲率和产状（TIN 方法无数学曲面的解析表达式，就无法求得对 于 x 或 y 的一阶和二阶偏导，而地质层面的走向、倾向和倾角又与曲面方程的偏导数有一定 的关系
［3 ］
） ，用于地质界面的模拟也是不合适的。离散数据集的单值曲面拟合插值法是较为
成熟、使用频率高的方法。与距离成反比的加权方法、径向基函数插值方法、样条曲面方法 等，拟合过程大都是全局性的，当增加、删除和修改数据点时都要重新计算，有些方法还必 须要求解维数很大的线性方程组，这就使得这些方法的效率和稳定性大大降低，寻求局部插 值拟合方法就成为今后地质界面数学模拟的发展目标。对于用来模拟复杂褶皱的多值曲面拟
N
F ( x, y ) =
∑
i =1
Ai (
ri r r ln i 2 + 1 − i 2 ) 2 R R R
2
2
2
为使 (x 1 , y1 ) 点处 F(x,y)等于实测数据点值 f 1 ，需建立联立方程组：
4
工程地质计算机应用
2006 年 第 3 期
总 43 期
fi =
∑A
j =1
n
j
(
rij
工程地质计算机应用
※
2006 年 第 3 期 ※
总 43 期
1
专题论述
问题探讨
三维地质建模的数学模型与显示方法
曾钱帮 1 何小萍 2
（1 中国科学院地质与地球物理研究所 北京 100029） (2 北京软通动力科技有限公司 北京 100027)
【摘要】基于离散数据集的曲面插值拟合方法，精确通过工程勘察数据点，获得光滑连续的地质 界面的数学模型，可以用于表达地形、地下水位面、岩层面、构造面等各种地质界面和岩土体物 理力学参数的空间分布。单值界面的数学模型中的插值型滑动最小二乘法是局部插值方法，可避 免全局插值方法的缺憾。编制程序生成 AutoCAD 脚本文件的地质界面计算机显示方法的优点是对 编程技术要求不高，简单实用，可充分利用 AutoCAD 软件强大的图形显示功能。 【关键词】工程地质 三维建模 勘测数据 单值曲面 拟合函数 地层曲面 计算机显示
(x - x k ) 2 + (y - y k ) 2 表示由(x,y) 点到(xk, yk)点的距离， 则插值函数 F(x , y)
F ( x, y ) =
∑ [d
k =1 n
n
k
∑ [d
k =1
n ( x, y ) ] = ∑ f k ⋅ Wk ( x , y ) 1 k =1 2
fk
k
( x, y ) ]
F ( x, y ) = ∑ a j ( x − x j ) 2 + ( y − y j ) 2 + c 2
j =1
n
[
]
1
1
2
, j = 1,2,..., n
其中c为常数，一般取1。将 n 个点 (x i , y i ) 的实测值 f i 代入上式建立联立方程：
f i = ∑ a j [( xi − x j ) 2 + ( y i − y j ) 2 + c 2 ] 2 , i = 1,2,..., n
2
工程地质计算机应用
2006 年 第 3 期
总 43 期
合插值函数也是一个公认的难题。 地质体通常是不规则形体， 由多个各种成因类型形状各异的结构面
［4］
围限而成， 有一定
［5 ］
的物质组成、结构和赋存于一定的地质环境中，遭受过多期次的变形和破坏
，其工程地质
条件复杂、千变万化。对于单个单值层面和地质条件简单的多层连续地层的模拟和计算机图 形显示比较简单，但是，当把多个这样的单值层面在空间叠加，考虑出现断层错断岩层、地 层不整合和结构面的组合、相互穿插等地质现象，组成地质体的三维模型，并能进行切剖， 其复杂度和建模技术难度会爆炸性地增加，一些理论问题（如三维拓扑结构分析，曲面求交 运算，形体相互遮挡的消隐等）需亟待解决。
ቤተ መጻሕፍቲ ባይዱ
2
工程地质计算机应用
式中权函数
2006 年 第 3 期
总 43 期
3
Wk ( x, y ) =
∏ [d
n j≠k n n k =1 j ≠ k
j
( x, y )
j
]
2
∑∏ [d
( x, y )
]
2
有如下性质：① Wk (x, y) 非负；② Wk (x, y) 是C0连续 的；③ Wk (x, y) = δ kj ,当k=j时， δ kj =1，否则 δ kj =0； ④Σ Wk (x, y) =1，具有加权性质。
W ( x, y ) = a0 + a1 x + a2 y + ∑ Fi ri 2 ln(ri 2 + ε )
i =1
n
式中， a0，a1，a2，f i (i = 1, 2,..., n) 为待定系数；ri = (x-xi)2+(y-yi)2； ε 为调节曲率大小的经
2
验参数，一般取 10-2~10-6；n 为原始数据点个数。 n+3 个待定系数可通过下列方程组求得：
T
0 1 xn yn
1 0 0 0
xn 0 0 0
y1 ⎤ ⎥ y2 ⎥ # ⎥ ⎥， yn −1 ⎥ yn ⎥ ⎥ 0 ⎥ 0 ⎥ ⎥ ⎥ 0 ⎦
改写成矩阵形式: AX = B 其中，
⎡ 0 r12 2 ln(r12 2 + ε ) ⎢ 2 2 0 ⎢ r12 ln(r12 + ε ) ⎢ # # ⎢ 2 2 2 2 r ln( r + ε ) r ln( r 1, n −1 2, n −1 2, n −1 + ε ) A = ⎢ 1,n −12 ⎢ r ln(r 2 + ε ) r2, n 2 ln(r2, n 2 + ε ) 1, n ⎢ 1, n 1 1 ⎢ ⎢ x1 x2 ⎢ ⎢ y1 y2 ⎣
1 前言
工程地质三维建模与可视化是应用计算机图形学和图像处理技术，是将工程地质勘测数 据和工程地质岩土体力学数值模拟分析的计算结果转换为图形图像在计算机屏幕上显示出 来并进行交互处理的理论、方法和技术。复杂地质体中的各类地质信息都可以被看作是三维 空间的函数，利用各种野外勘测数据分别建立相应的曲面拟合函数，进而利用计算机建立三 维地质模型，逼真反映地质结构全貌，达到直观地表达地质信息的分布规律、提高对于地质 规律的认识、指导地质工程项目的勘测施工的目的。因此，工程地质三维建模与可视化研究 有其重要的理论和现实意义。 工程地质三维建模与可视化研究中，地质界面的数学模拟是基础。由于地质界面必须精 确通过控制点（工程勘测数据点） ，通过离散的工程勘测数据点建立三维地质模型及其计算 机图形可视化显示属于重构问题，与计算机图形学中对于机械设计特别有效的形体的构造问 题是有根本区别的， 所以针对机械设计发展起来的自由曲线曲面造型技术