大地电磁法(MT)二维有限元正演模拟系统
项目介绍
本系统是一款用于地球物理电磁探测领域的二维正演模拟工具。它基于有限元法(FEM)实现,能够计算在各种地下电阻率结构下的横电波(TE)和横磁波(TM)电磁响应。该系统旨在模拟复杂地质条件(如高阻体、低阻体的嵌入)对大地电磁场的影响,并计算出生产实践中常用的视电阻率和相位参数,为地球物理资料解释和反演研究提供理论数据支撑。
功能特性
- 双极化模式求解:支持TE模式(电场沿走向)和TM模式(磁场沿走向)的独立数值模拟。
- 非均匀网格剖分:具备构建核心区域细分、四周边界扩展的非轴对称矩形网格能力,有效平衡计算精度与边界效应抑制。
- 复杂模型构建:支持在背景介质中灵活插入任意形状和电阻率的异常体。
- 高性能稀疏矩阵计算:采用稀疏矩阵装配与求解技术,优化了大型线性方程组的计算效率。
- 一维边界约束:通过一维解析解外推方法处理模型边界,提高了数值模拟的准确性。
- 多维度结果可视化:自动生成地电模型图、单站频率响应曲线以及双模式视电阻率剖面拟断面图。
系统要求
- 运行环境:MATLAB R2016b 或更高版本。
- 核心组件:基础MATLAB组件(需支持 sparse、interp1 等数值计算函数)。
- 硬件建议:8GB RAM 或以上,以支持大规模网格的矩阵运算。
逻辑流程与实现详情
1. 初始化与网格准备
程序首先定义频率序列(100Hz至0.01Hz对数分布)及物理常数。网格系统采用不均匀设计:
- 核心区:在水平和深度方向定义精细步长,用于捕捉异常体附近的电磁场变化。
- 扩展区:在核心区外围增加填充层,防止人工边界反射干扰计算区域。
- 测点布置:在地表特定范围均匀布置观测站。
2. 地电模型定义
系统初始化均匀背景(如100 $Omegacdot m$),并基于坐标索引向模型中注入特定异常:
- 高阻体:模拟地下坚硬岩石或空洞特征(5000 $Omegacdot m$)。
- 低阻体:模拟含水层、金属矿体等高导特征(10 $Omegacdot m$)。
3. 有限元数值求解核心
核心求解器通过以下步骤处理电磁波动方程:
- 单元刚度矩阵计算:针对矩形单元,通过变分原理将二阶偏微分方程转化为代数形式。根据TE或TM模式切换耦合系数($tau$ 和 $gamma$),其中包括了基于线性插值函数的梯度矩阵和质量矩阵的组合。
- 全局矩阵装配:利用稀疏格式构建大型刚度矩阵,将局部单元贡献合并至全局坐标系统。
- 边界条件加载:
-
顶部:设为Dirichlet边界(场强为1),模拟入射平面波。
-
底部:设为0,模拟电磁场在深部的衰减。
-
两侧:根据背景电阻率的一维解析解进行动态填充。
- 线性方程组求解:直接调用稀疏解法器求解复数线性方程组。
4. 电磁场分量提取与参数计算
- 辅助场计算:利用中心差分法,根据主场(E或H)计算其空间导数,从而获得辅助场分量。
- 插值处理:使用样条插值(spline)将网格节点上的场值精确映射到观测站坐标上。
- 视电阻率与相位:基于卡尼亚公式提取复阻抗 Z,并计算各频率下的视电阻率和阻抗相位。
关键算法分析
- 有限元离散化:本系统采用了基于矩形单元的双线性插值。刚度矩阵 A 的构建体现了电磁波在不同电阻率介质中的传播衰减特性。
- 矩阵稀疏化技术:由于二维有限元生成的矩阵具有极高的稀疏性,代码通过预分配索引数组并使用 sparse 函数,极大地节省了内存。
- TE/TM 模式切换逻辑:
-
TE模式:求解电场 $E_y$,通过其法向导数获取磁场 $H_x$。
-
TM模式:求解磁场 $H_y$,通过其法向导数乘以电阻率获取电场 $E_x$。
- 结果表现:可视化模块将抽象的复数解转化为工程直观的对数周期曲线和空间断面,能够清晰展示高/低阻异常体对地表观测点场值的影响。