gprMax中多相随机介质探地雷达三维建模与模拟

一种基于gprMax的多相随机介质探地雷达三维建模与模拟

实际地下介质是非均匀介质，但数值模拟时常常把介质当做均匀介质，难以对实际介质产生准确认识。常规gprMax建模都是均匀介质建模、规则形状建模，本博文给出一种gprmax中多相随机介质建模的方法。

1、随机介质建模

利用matlab编程生成随机三维数组，在三维随机数组中填充材料，即可通过gprmax模拟出随机介质的探地雷达结果。理论上，三维随机介质建模可以对任意多种随机介质建模。首先看一下建模的效果，如下图所示。

模型分为三层，上下两层为均匀介质，中间层为三相随机介质，包含水、空气、泥三种介质。

2、paraview显示二维模型

将运行结果的*.vti文件导入到paraview软件中成图结果如下。

3、paraview显示三维模型

三维显示。

4、模拟结果

在matlab上显示探地雷达模拟的结果，如下。


从图中看，与均匀介质还是有较大区别，随机介质中有较多不均匀反射同相轴，比较符合实际情况，证明此建模方法是有效的。

三维情况下，一条测线的雷达剖面。

5、in文件

gprMax中的输入文件如下。
二维in文件。

#title: rand2D
#domain: 4.000 2.000 0.002
#dx_dy_dz: 0.002 0.002 0.002
#time_window: 40e-9
#waveform: ricker 1.0 400e6 my_ricker
#hertzian_dipole: z 0.05 1.95 0 my_ricker
#rx: 0.15 1.95 0
#src_steps: 0.050 0.000 0.000
#rx_steps: 0.050 0.000 0.000
#geometry_objects_read: 0 0 0 D:\Learnfile/stestbig01.h5 D:\Learnfile/materials.txt
#geometry_view: 0 0 0 4.000 2.000 0.002 0.002 0.002 0.002 Tuokong2Drand n

三维in文件

#title: rand3D

#domain: 4.000 2.000 2.000
#dx_dy_dz: 0.01 0.01 0.01
#time_window: 40e-9

#waveform: ricker 1.0 400e6 my_ricker
#hertzian_dipole: z 0.10 1.85 1.00 my_ricker

#rx: 0.20 1.85 1.00
#src_steps: 0.050 0.000 0.000
#rx_steps: 0.050 0.000 0.000

#geometry_objects_read: 0 0 0 D:\Learnfile/stestbig02.h5 D:\Learnfile/materials.txt

geometry_view: 0 0 0 4.000 2.000 2.000 0.01 0.01 0.01 rand3D n

6、matlab代码

close all
clear
clc

% 此程序是写gprmax随机介质建模的几何文件
% 创建任意随机介质
% 创建地下层状多相随机介质

% 给定网格单元边长
dx = 0.002;
dy = 0.002;
dz = 0.002;

% 网格单元个数，注意x,y,x的方向
nex = 1;
ney = 1000;
nez = 2000;

% 设置保存的文件名
file_name = 'stestbig01.h5';

% 创建一个HDF5文件
% h5create(file_name,'/data',[nex ney nez]);

% 创建一个充满随机数的大矩阵
data = zeros(nex,ney,nez);

% 第一层填充0
% 第一层默认为0
data(:,1:250,:) = 0;

% 第二层填充随机介质
% 创建一个随机的三维数组
data_rand = rand(1,500,2000);
data_rand(data_rand <= 0.3) = 0;
data_rand(data_rand >= 0.7) = 5;
data_rand(data_rand >= 0.3&data_rand <= 0.7) = 0.5;
data_rand(data_rand == 0) = 1;
data_rand(data_rand == 0.5) = 2;
data_rand(data_rand == 5) = 3;

data(:,251:750,:) = data_rand;
data(:,751:1000,:) = 4;

data = flip(data,2);

% 三维转二维，绘图
data1D = data(:);
data2D = reshape(data1D,ney,nez);
% imagesc(data2D);colorbar
pcolor(data2D);colormap(jet),colorbar
shading interp

具体的文件我打包成文件夹，上传至CSDN资源，欢迎下载。