直流电阻率2.5维有限元正演模拟程序

资 源 简 介

本程序旨在实现直流电阻率法的2.5维正演数值模拟，是地球物理电法勘探领域重要的理论研究工具。2.5维正演的含义是在地质结构具有二维特征（即沿走向方向无限延伸）的前提下，计算三维点源产生的电场分布。 程序的核心逻辑是利用傅里叶变换将三维电位控制方程转换到波数域，从而将三维问题简化为一系列准二维问题。具体实现步骤包括：首先对研究区域进行离散化网格划分（如采用三角形或矩形单元），利用有限单元法（FEM）构建不同波数下的刚度矩阵与线性方程组；接着选取一组最佳的离散波数，在波数域内分别求解代数方程组以获得各个波数对

详 情 说 明

直流电阻率2.5维正演模拟系统 (MATLAB)

项目介绍

本项目是一款基于 MATLAB 开发的直流电阻率法 2.5 维正演数值模拟程序。在地球物理电法勘探中，2.5 维正演是指在地质构造沿走向（通常为Y方向）无限延伸的前提下，计算三维点源产生的电场分布。本程序通过傅里叶变换将三维电位方程简化到波数域进行求解，利用有限单元法（FEM）实现高精度的电场数值模拟，是理解电法勘探物理机制、验证反演算法的重要工具。

功能特性

• 2.5维数值模拟：采用波数域离散化技术，将复杂的三维点源问题转化为多个准二维问题求解。
• 有限单元法（FEM）引擎：基于矩形单元和双线性插值函数构建刚度矩阵，能够精确处理非均匀介质。
• 灵活的模型构建：支持自定义背景电阻率及异常体形状、位置和电阻率参数。
• 温纳排列仿真：内置温纳（Wenner）电极排列扫描逻辑，自动计算极距变化下的观测响应。
• 自动可视化：程序运行后自动生成地质电阻率模型图及对应的视电阻率拟断面图。
• 优化波数选取：内置经过优化的离散波数及其对应的权重，利用 5 个关键波数即可实现高精度的反傅里叶变换还原。

系统要求

• 软件环境：MATLAB R2016a 或更高版本。
• 硬件要求：建议 8GB RAM 以上，以支持有限单元法中大规模稀疏矩阵的运算。

核心实现逻辑

程序的主体逻辑严格遵循地球物理数值模拟的标准流程：

1. 网格化与模型初始化 程序首先在 X-Z 平面上建立矩形网格，通过定义的网格数（nx, nz）和间距（dx, dz）生成坐标系统。电阻率模型以矩阵形式存储，默认设置为均匀背景，并允许通过索引操作在特定区域嵌入低阻或高阻异常体。

2. 波数域方程转换 为了处理三维点源，程序引入了 5 个优化的离散波数（$lambda$）。在每个波数下，三维泊松方程转化为波数域的偏微分方程：$nabla cdot (sigma nabla V) - lambda^2 sigma V = -I delta$。

3. 有限单元矩阵组装 对于每一个选定的波数，程序会组装全局刚度矩阵：

• 微分贡献（Ke）：计算电导率在空间导数上的贡献，反映电流在二维平面的扩散。
• 波数贡献（Me）：计算由波数项引起的电位衰减，这体现了点源三维特性的映射。
• 边界处理：采用 Dirichlet 边界条件，通过置 1 法强制将模型边界（左、右、底、顶边界）的电位设为零，模拟无穷远处电位趋于零的物理现象。

4. 观测系统扫描与求解 程序模拟温纳排列的移动过程：

• 根据给定的极距 $a$ 确定 $A, M, N, B$ 四个电极的位置。
• 在右端项（RHS）中为 A 电极注入正电流，B 电极注入负电流。
• 分别求解 5 个波数下的线性方程组，获取 $M, N$ 点在波数域的电位值。

5. 反傅里叶变换与视电阻率计算 通过加权组合算法将 5 个波数域的解进行叠加，恢复至空间域（$y=0$ 处）的真实电位：$V = frac{2}{pi} sum V_k cdot w_k$。随后利用温纳排列的几何因子 $K = 2pi a$ 计算视电阻率 $rho_a$。

关键函数与算法说明

• 刚度矩阵组装算法：

该算法遍历所有矩形单元，计算其 4x4 的局部局部贡献。局部矩阵结合了拉普拉斯项（Ke）和质量矩阵项（Me，受 $lambda^2$ 缩放）。这种处理方式确保了模型在波数域的稳定性。

• 稀疏矩阵存储：

考虑到有限单元法生成的矩阵具有高度稀疏性，程序采用 MATLAB 的 sparse 函数构建全局矩阵，极大地优化了内存占用并提升了求解线性方程组（K rhs）的速度。

• 结果插值平滑：

由于观测点是离散的，程序在绘图阶段利用 griddata 对拟断面数据进行线性插值，从而生成平滑的视电阻率彩色云图，便于分析异常体的响应特征。

使用方法

1. 打开 MATLAB 并将工作目录切换至程序所在文件夹。
2. 在命令行窗口输入主程序名并回车。
3. 程序将依次执行以下流程：

* 构建 $61 times 31$ 的地质网格。 * 在模型中部生成一个 $10,Omegacdottext{m}$ 的低阻异常体。 * 循环计算不同电极位置和极距下的电位响应。 * 自动弹出绘图窗口，展示“地质电阻率模型”与“视电阻率温纳排列拟断面图”。

1. 用户可以通过修改程序开头的参数（如 res_mod 的范围或 dx, dz）来模拟不同的地质场景。