CN116859478B - 一种基于瞬变电磁法成像的地下水模拟方法及系统 - Google Patents
一种基于瞬变电磁法成像的地下水模拟方法及系统 Download PDFInfo
- Publication number
- CN116859478B CN116859478B CN202311137070.7A CN202311137070A CN116859478B CN 116859478 B CN116859478 B CN 116859478B CN 202311137070 A CN202311137070 A CN 202311137070A CN 116859478 B CN116859478 B CN 116859478B
- Authority
- CN
- China
- Prior art keywords
- resistivity
- data
- model
- volume fraction
- clay
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 100
- 230000001052 transient effect Effects 0.000 title claims abstract description 69
- 239000003673 groundwater Substances 0.000 title claims abstract description 60
- 238000004088 simulation Methods 0.000 title claims abstract description 29
- 238000003384 imaging method Methods 0.000 title claims abstract description 18
- 238000009826 distribution Methods 0.000 claims abstract description 84
- 239000004576 sand Substances 0.000 claims abstract description 50
- 230000004044 response Effects 0.000 claims abstract description 41
- 238000011160 research Methods 0.000 claims abstract description 39
- 239000007787 solid Substances 0.000 claims abstract description 35
- 238000004364 calculation method Methods 0.000 claims abstract description 32
- 239000010410 layer Substances 0.000 claims description 103
- 239000004927 clay Substances 0.000 claims description 85
- 238000012549 training Methods 0.000 claims description 61
- 230000008569 process Effects 0.000 claims description 25
- 230000035945 sensitivity Effects 0.000 claims description 25
- 239000000463 material Substances 0.000 claims description 22
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims description 17
- 238000012360 testing method Methods 0.000 claims description 12
- 238000012795 verification Methods 0.000 claims description 12
- 238000005457 optimization Methods 0.000 claims description 10
- 238000004422 calculation algorithm Methods 0.000 claims description 9
- 239000011435 rock Substances 0.000 claims description 9
- 238000009933 burial Methods 0.000 claims description 8
- 239000002344 surface layer Substances 0.000 claims description 6
- 238000006243 chemical reaction Methods 0.000 claims description 5
- 238000004140 cleaning Methods 0.000 claims description 5
- 230000007246 mechanism Effects 0.000 claims description 5
- 238000012216 screening Methods 0.000 claims description 5
- 238000010276 construction Methods 0.000 claims description 4
- 238000000638 solvent extraction Methods 0.000 claims description 4
- 238000013473 artificial intelligence Methods 0.000 abstract description 2
- 230000006870 function Effects 0.000 description 35
- 238000013528 artificial neural network Methods 0.000 description 10
- 210000004027 cell Anatomy 0.000 description 6
- 230000008859 change Effects 0.000 description 6
- 238000011835 investigation Methods 0.000 description 6
- 230000035699 permeability Effects 0.000 description 6
- 241000282414 Homo sapiens Species 0.000 description 5
- 230000004913 activation Effects 0.000 description 4
- 238000012937 correction Methods 0.000 description 4
- 238000005070 sampling Methods 0.000 description 4
- 239000013049 sediment Substances 0.000 description 4
- 230000002159 abnormal effect Effects 0.000 description 3
- 238000013527 convolutional neural network Methods 0.000 description 3
- 238000010586 diagram Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 238000011156 evaluation Methods 0.000 description 3
- 230000008676 import Effects 0.000 description 3
- 239000011159 matrix material Substances 0.000 description 3
- 238000012545 processing Methods 0.000 description 3
- 239000002356 single layer Substances 0.000 description 3
- 230000008901 benefit Effects 0.000 description 2
- 230000015572 biosynthetic process Effects 0.000 description 2
- 238000001514 detection method Methods 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 239000003344 environmental pollutant Substances 0.000 description 2
- 239000000284 extract Substances 0.000 description 2
- 238000005755 formation reaction Methods 0.000 description 2
- 230000002779 inactivation Effects 0.000 description 2
- 238000013507 mapping Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000003062 neural network model Methods 0.000 description 2
- 210000002569 neuron Anatomy 0.000 description 2
- 231100000719 pollutant Toxicity 0.000 description 2
- 239000011148 porous material Substances 0.000 description 2
- 230000001105 regulatory effect Effects 0.000 description 2
- 239000002689 soil Substances 0.000 description 2
- 230000007480 spreading Effects 0.000 description 2
- 238000003892 spreading Methods 0.000 description 2
- 238000012876 topography Methods 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- ORILYTVJVMAKLC-UHFFFAOYSA-N Adamantane Natural products C1C(C2)CC3CC1CC2C3 ORILYTVJVMAKLC-UHFFFAOYSA-N 0.000 description 1
- 238000012935 Averaging Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000033228 biological regulation Effects 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000015556 catabolic process Effects 0.000 description 1
- 238000001311 chemical methods and process Methods 0.000 description 1
- 230000007423 decrease Effects 0.000 description 1
- 238000013135 deep learning Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000006731 degradation reaction Methods 0.000 description 1
- 238000000151 deposition Methods 0.000 description 1
- 230000008021 deposition Effects 0.000 description 1
- 238000005137 deposition process Methods 0.000 description 1
- 238000009792 diffusion process Methods 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 239000013505 freshwater Substances 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000010365 information processing Effects 0.000 description 1
- 238000012423 maintenance Methods 0.000 description 1
- 238000013508 migration Methods 0.000 description 1
- 230000005012 migration Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 238000012856 packing Methods 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 230000001737 promoting effect Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000000630 rising effect Effects 0.000 description 1
- 239000000126 substance Substances 0.000 description 1
- 239000002352 surface water Substances 0.000 description 1
- 230000004083 survival effect Effects 0.000 description 1
- 238000010200 validation analysis Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
- G01V3/38—Processing data, e.g. for analysis, for interpretation, for correction
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
- G01V3/08—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation operating with magnetic or electric fields produced or modified by objects or geological structures or by detecting devices
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
Landscapes
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Electromagnetism (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明提供了一种基于瞬变电磁法成像的地下水模拟方法及系统，属于水文地球物理和人工智能技术领域，方法包括：采用实际瞬变电磁响应数据，利用训练完毕的CNN‑LSTM网络模型预测地下的电阻率；利用研究区域内的测井资料，结合反演的地下电阻率，采用最优化插值方法，计算区域内不同测点不同深度层对应的砂/砂砾的含量分布；将不同测点不同深度对应的砂/砂砾的含量分布数据进行简化，并将简化的数据以测井数据形式导入GMS构建solid模型；将solid模型导入MODFLOW模块，结合有限的地质信息设置计算参数，计算区域地下水水位分布。本发明相较于传统反演方法，可以大幅提升反演速度，同时保证反演结果的稳定性与可靠性。
Description
技术领域
本发明属于水文地球物理和人工智能技术领域，更具体地，涉及一种基于瞬变电磁法成像的地下水模拟方法及系统。
背景技术
浅层的地下水是水循环系统中重要的组成部分，它所分布的近地表拥有完善的孔隙网格。浅层地下水通过土层孔隙网格进行物质输运和迁移，在改善土壤水力条件、促进陆地水、维系区域淡水生态系统等方面发挥着重要作用。由于气候变化和人类活动的影响，浅层地下水质量退化以及含量减少对人类生存发展和区域生态系统维系产生了直接的负面影响。特别是地下水的污染物调查和监测问题，一直都是人们关注的焦点，污染物的扩散以及分布关系到人类对于浅层地下水的正常利用。因此，为了平衡浅层地下水利用与人类发展的动态关系，建立科学的地下水管理计划和评估方法极为必要。
地下水建模的目的是帮助决策者评估地下水的赋存状态，进而实现对地下水的高效利用。通常，有价值的地下水模型应当包含对近地表水文构造以及有关重要物理、化学过程等关键信息。结合地球物理电磁法探测进行地下水建模已成为目前评估地下水资源赋存状态的有效方法。
地球物理电磁法可以无侵入、高效地获取研究区域地下介质的电阻率参数，进而获取地下水信息并可以进行大尺度区域的水文地质建模。这已经成为地下水调查领域一个重要的研究方向。需要注意的是，对于大尺度水文建模研究，为保障地质结构调查以及水文建模的精确性，数据采样密度应足够大。地面拖拽式瞬变电磁（tTEM）系统的开发满足了近地表地下结构高分辨率，密集数据采样的成像需求。为了应对密集、大规模的数据，有必要开发一种基于深度学习架构的大规模数据处理方法，提高地球物理数据的反演解释能力。
发明内容
针对现有技术的缺陷，本发明的目的在于提供一种基于瞬变电磁法成像的地下水模拟方法及系统，旨在解决现有的大规模瞬变电磁数据反演与大尺度水文建模研究效率较低的问题。
为实现上述目的，本发明提供了一种基于瞬变电磁法的地下水模拟方法，包括以下步骤：
步骤一：采用实际瞬变电磁响应数据，利用训练完毕的CNN-LSTM网络模型预测地下的电阻率；其中，CNN为卷积神经网络，LSTM为长短时记忆神经网络；
步骤二：利用研究区域内的测井资料，结合反演的地下电阻率，采用最优化插值方法，依次计算区域内不同测点不同深度层对应的砂/砂砾的含量分布；
步骤三：将研究区域内不同测点不同深度对应的砂/砂砾的含量分布数据进行简化，并将简化的数据以测井数据形式导入GMS构建solid模型；
步骤四：将solid模型导入MODFLOW模块，结合有限的地质信息设置计算参数，计算区域地下水水位分布；
其中，训练CNN-LSTM网络模型的方法为：
基于连续层状电阻率模型，结合已知的发射装置参数，通过瞬变电磁正演方法获取理论的瞬变电磁响应，形成用于训练CNN-LSTM网络的初始样本数据集，训练时CNN-LSTM网络的损失函数采用均方误差，同时在损失函数中添加由电阻率灵敏度约束构建的正则化项，输出为电阻率参数。
进一步优选地，训练CNN-LSTM网络模型的具体方法，包括以下步骤：
基于构建的固定层数的、固定层厚的连续层状电阻率模型，结合已知的发射装置参数，通过瞬变电磁正演方法计算对应时刻的瞬变电磁响应，形成初始样本数据集；
利用数据清洗方法清除初始样本数据集中的冗余数据，按一定比例划分为训练集、验证集和测试集；
构建包含电阻率灵敏度约束的损失函数，采用训练集和验证集训练CNN-LSTM网络，基于早停机制和测试集测试选择最优化CNN-LSTM网络模型。
进一步优选地，包含电阻率灵敏度约束的损失函数为：
其中，是正则化参数，MSE表示均方误差；/>为归一化后的灵敏度约束权重；
其中，是样本数量，/>是第/>个样本的真实电阻率参数，/>是对第/>个样本的预测电阻率参数；
其中，为数据集标签中理论观测瞬变电磁响应数据点上的噪声；/>
其中，为第i个数据集标签中理论观测瞬变电磁响应；/>表示第j个层状电阻率模型参数，/>为对应电阻率灵敏度计算值。
进一步优选地，最优化反演方法，包括以下步骤：
统计研究区域的岩性测井资料，以及测井及地球物理测点的空间位置，同时在研究区域表层构建均匀、离散正交的控制网格，并筛选出包围有地球物理测点的有效网格节点，在纵向上结合研究区域地质条件，构建均匀或具有不同层厚的纵向网格分布；
对岩性测井资料进行解译得到对应层的粘土体积分数及层厚度，并通过置信度算法对解译的粘土体积分布及层厚度进行调整，基于层厚度加权方式，计算岩性测井资料对应层的粘土体积分数作为最优化反演插值中的观测数据，同时设置有效网格节点上的初始电阻率阈值；
对于相同埋深的深度层，每次反演时采用双线性插值方法，由分布在单个正方形单元周围四个网格节点的初始电阻率阈值，计算包围电阻率模型的电阻率阈值，进而计算单个电阻率模型对应层的粘土体积分数；
以测井位置点为圆心，设置搜索半径，在搜索半径内对电阻率模型转换得到的粘土体积分数进行克里金插值，得到测井位置的粘土体积分数插值结果；
基于测井位置的粘土体积分数插值结果采用高斯-牛顿最优化反演方法实现最优化插值过程，依次计算不同深度水平层对应的粘土体积分数分布，基于粘土/砂总含量为1的假设可以得到砂/砂砾的含量分布。
进一步优选地，电阻率转换为粘土体积分数采用的转换函数为：
其中，和/>分别为电阻率的最高阈值和电阻率的最低阈值；erfc为误差函数；/>为电阻率。
进一步优选地，构建最优化插值过程中的反演目标函数为：
其中，是观测测井数目，/>是数据拟合项，/>是约束对数目，/>和分别表示第i个测井对应的插值粘土体积分数和观测粘土体积分数，/>表示数据方差，/>和/>表示相邻节点的粘土体积分数，/>为调控变化的测度常数。
进一步优选地，步骤三具体包括以下步骤：
对研究区域进行有效分区，并选择区域中心点作为区域所有结果的表征，从而对数据进行简化，其中，中心点数值为区块内所有对应深度测点砂/砂砾体积分数的平均；
根据每个区域的砂岩含量，设置三个砂/砾石体积分数分布截断值将区域材料分为四种不同属性的填充材料，对电阻率模型的不同层厚进行材料标记；
将简化后的砂/砂砾的含量分布数据以及材料标记以导入测井资料的形式导入GMS，进行Horizons ID智能识别，并判断识别的准确性，若识别有误则采用自行编码，并对相应的Horizons ID进行更正，在局部建模中增加对对应层的水平分布约束，从而构建GMSsolid实体模型。
进一步优选地，计算区域地下水分布的方法，包括以下步骤：
将solid模型导入新建的MODFLOW模块；
利用测井资料及相关信息定义填充材料属性，其中，填充材料属性为岩块的各种水力参数，包括垂直渗透率和水平渗透率；
在GMS内新建图层，定义研究区域、流动边界以及各种水文属性，并设置有效的计算条件；
对MODFLOW中的输入参数、设定的边界条件进行检验，进行计算并得到研究区域水位分布。
另一方面，本发明提供了一种基于瞬变电磁法成像的地下水模拟系统，包括：
地下电阻率预测模块，用于采用实际瞬变电磁响应数据，利用训练完毕的CNN-LSTM网络模型预测地下的电阻率；
砂/砂砾的含量分布计算模块，用于利用研究区域内的测井资料，结合反演的地下电阻率，采用最优化插值方法，依次计算区域内不同测点不同深度层对应的砂/砂砾的含量分布；
solid模型的构建模块，用于将研究区域内不同测点不同深度对应的砂/砂砾的含量分布数据进行简化，并将简化的数据以测井数据形式导入GMS构建solid模型；
地下水水位分布计算模块，用于将solid模型导入MODFLOW模块，结合有限的地质信息设置计算参数，计算区域地下水水位分布；
其中，训练CNN-LSTM网络模型的方法为：
基于连续层状电阻率模型，结合已知的发射装置参数，利用瞬变电磁正演方法获取理论的瞬变电磁响应，形成用于训练CNN-LSTM网络的初始样本数据集，训练时CNN-LSTM网络的损失函数采用均方误差，同时在损失函数中添加由电阻率灵敏度约束构建的正则化项，输出为电阻率参数。
进一步优选地，包含电阻率灵敏度约束的损失函数为：
其中，是正则化参数，MSE表示均方误差；/>为归一化后的灵敏度约束权重；
其中，是样本数量，/>是第/>个样本的真实电阻率参数，/>是对第/>个样本的预测电阻率参数；
其中，为数据集标签中理论观测瞬变电磁响应数据点上的噪声；
其中，为第i个数据集标签中理论观测瞬变电磁响应；/>表示第j个层状电阻率模型参数，/>为对应电阻率灵敏度计算值。
进一步优选地，砂/砂砾的含量分布计算模块包括：
网格构建单元，用于统计研究区域的岩性测井资料，以及测井及地球物理测点的空间位置，同时在研究区域表层构建均匀、离散正交的控制网格，并筛选出包围有地球物理测点的有效网格节点，在纵向上结合研究区域地质条件，构建均匀或具有不同层厚的纵向网格分布；
初始电阻率阈值的设置单元，用于设置有效网格节点上的初始电阻率阈值；
观测粘土体积分数的计算单元，用于对岩性测井资料进行解译得到对应层的粘土体积分数及层厚度，并通过置信度算法对解译的粘土体积分布及层厚度进行调整，基于层厚度加权方式，计算岩性测井资料对应层的粘土体积分数作为最优化插值中的观测数据；
相同埋深粘土体积分数的计算单元，用于基于纵向网格分布中具有相同埋深的深度层，每次反演时采用双线性插值方法，由分布在单个正方形单元周围四个网格节点的初始电阻率阈值，计算包围电阻率模型的电阻率阈值，进而计算单个电阻率模型对应层的粘土体积分数；
搜索半径内粘土体积分数的插值单元，用于以测井位置点为圆心，设置搜索半径，在搜索半径内对电阻率模型转换得到的粘土体积分数进行克里金插值，得到测井位置的粘土体积分数插值结果；
不同埋深粘土体积分数的计算单元，用于基于测井位置的粘土体积分数插值结果，采用高斯-牛顿最优化反演方法优化插值过程，依次计算不同深度水平层对应的粘土体积分数分布，基于粘土/砂总含量为1的假设可以得到砂/砂砾的含量分布。
进一步优选地，电阻率转换为粘土体积分数采用的转换函数为：
其中，和/>分别为电阻率的最高阈值和电阻率的最低阈值；erfc为误差函数；/>为电阻率。
进一步优选地，solid模型的构建模块包括：
砂/砂砾含量分布的简化单元，用于对研究区域进行有效分区，并选择区域中心点作为区域所有结果的表征，从而对不同测点不同深度层对应的砂/砂砾的含量分布进行简化，其中，中心点数值为区域内所有对应深度测点砂/砂砾体积分数的平均；
区域填充单元，用于根据每个区域的砂岩含量，设置三个砂/砾石体积分数分布截断值将区域材料分为四种不同属性的填充材料，对电阻率模型的不同层厚进行材料标记；
数据导入处理单元，用于将简化后的砂/砂砾的含量分布数据以及材料标记以导入测井资料的形式导入GMS，进行Horizons ID智能识别，并判断识别的准确性，若识别有误则采用自行编码，并对相应的Horizons ID进行更正，在局部建模中增加对对应层的水平分布约束，从而构建GMS solid实体模型。
总体而言，通过本发明所构思的以上技术方案与现有技术相比，具有以下有益效果：
本发明提供了一种基于瞬变电磁法成像的地下水模拟方法及系统，其中，所构建的反演流程（采用实际瞬变电磁响应数据，利用训练完毕的CNN-LSTM网络模型预测地下的电阻率）相较于传统反演方法，可以大幅提升反演速度，在网络完成训练后，可以实现大规模瞬变电磁数据的快速反演与成像；结合大量前期产生的电磁数据样本，可以降低处理大规模瞬变电磁数据的计算成本，同时保证反演结果的稳定性与可靠性。
本发明提供了一种基于瞬变电磁法成像的地下水模拟方法及系统，其中，本发明中采用CNN-LSTM网络模型预测地下的电阻率，通过卷积神经网络提取数据特征，长短期记忆神经网络学习特征中的时序性质，网络训练的损失函数采用均方误差，同时添加了具有地球物理意义的约束，使得网络收敛更稳定，增强了网络对浅层电阻率的反演识别能力，令网络在训练迭代中更稳定，更符合传统反演规律的假设。
本发明提供了一种基于瞬变电磁法成像的地下水模拟方法及系统，其中，采用三次样条插值构建连续层状电阻率模型描述地下电阻率的变化，生成数据集标签作为网络模型的数据集，该数据集可以使网络输出垂向更光滑的电阻率反演结果，同时能较好地反映电阻率模型的突变界面。在实际数据应用上，可以结合地面拖曳式瞬变电磁测线方向上密集采样的优势，反演结果能较好地描述实际地层在垂向及横向上的变化。
本发明提供了一种基于瞬变电磁法成像的地下水模拟方法及系统，其中，采用最优化插值技术可以将获取的一维电阻率模型转换为砂/粘土分布。因为砂/砂砾层的渗透性要强于淤泥/粘土层，所以这对于研究地下水的流动极为重要。包括通过统计区域的岩性测井资料，以及测井及地球物理测点的空间位置，同时在水平层面构建均匀、离散正交的控制网格，可以建立一个相对稀疏的结构化网格嵌套所有观测的地球物理测点；每次反演时，对于不同深度水平层，基于单个正方形单元的四个网格节点阈值，采用双线性插值技术更新每个电阻率模型的阈值，计算出对应粘土体积分数；以测井位置点为圆心，设置合适的搜索半径，在搜索半径内对电阻率模型转换的粘土体积分数进行克里金插值，从而得到测井位置的插值结果；结合高斯-牛顿最优化反演方法获取最终的插值结果，并依次计算不同深度层对应的粘土体积分数分布，基于粘土/砂总含量为1的假设可以得到砂/砂砾的含量分布。
本发明提供了一种基于瞬变电磁法成像的地下水模拟方法及系统，其中，采用GMS对地下水的建模过程可以严格计算出研究区域的地下水分布。包括对测点的观测数据进行了有效分区，并选择区块中心点作为该区域所有结果的表征，中心点数值为该区块对应深度数值的平均；考虑到了每个单元的砂岩含量，设置4种不同属性的填充材料，可以将得到的砂/砂砾分布的岩石物理模型简化得到对应的岩性模型；为了满足将数据以测井资料的形式导入GMS构建solid模型时Horizons ID的对应性，采用综合考虑输入数据的层分布特性，并自行编码获取相对准确的Horizons ID的方法，同时在局部建模中增加对对应层的水平分布约束，以获得准确的GMS solid实体模型；利用构建的GMS solid模型，在GMS创建MODFLOW 计算模块，计算过程选择的版本为MODFLOW 2000，并采用了均匀网格覆盖所创建的实体模型，结合有限的地质信息设置计算参数，计算区域地下水水位分布。
附图说明
图1是本发明实施例提供的基于瞬变电磁法成像的地下水模拟方法流程图；
图2是本发明实施例提供的连续层状电阻率模型示意图；
图3是本发明实施例提供的地面拖拽式瞬变电磁双发射机采集序列的正演响应示意图；
图4是本发明实施例提供的CNN-LSTM网络框架示意图；
图5（a）是本发明实施例提供的理论模型反演结果示意图；
图5（b）是本发明实施例提供的神经网络预测结果验证图；
图6是本发明实施例提供的研究区域数据拟合情况图；
图7是本发明实施例提供的最优化反演插值过程中的离散网格示意图；
图8（a）是本发明实施例提供的基于MODFLOW计算水位时测井和测点位置示意图；
图8（b）是本发明实施例提供的基于MODFLOW计算水位与真实水位结果对比图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。
一方面，本发明提供了一种基于瞬变电磁法成像的地下水模拟方法：
将采集地面的实际观测瞬变电磁响应数据输入至训练完毕的CNN-LSTM网络模型中，获取反演电阻率模型实现瞬变电磁成像；基于最优化反演插值技术，将反演的电阻率模型转换为地下介质对应粘土含量的空间分布，从而构建粘土/砂砾分布的岩性模型；基于GMS构建输入的solid模型，并结合有限离散的地层水力参数信息（渗透率），利用MODFLOW模拟浅层地下水位分布。
依托实施例，如图1所示，具体包括以下步骤：
S101：根据研究区域的地质背景生成随机的电阻率模型，并通过瞬变电磁正演获得理论瞬变电磁响应数据；
地面拖拽式瞬变电磁系统的有效探测深度范围为0-80m的近地表，在该范围内沉积物主要为碎屑岩，由于沉积过程通常是稳定的，因而沉积层的电阻率是连续渐变的；为了模拟真实地层，本发明采用三次样条插值方法使得生成模型在纵向上地层的电阻率是连续分布的，理论瞬变电磁响应获得方法包括：
由于该系统有效勘探深度通常在100m以内，因此设定最大反演深度为130m，并设置6个随机分布点控制电阻率模型曲线的展布形态，其中相邻控制点间隔应不小于10m；
控制点x坐标为电阻率，y坐标为深度，其中首尾控制点深度分别固定在0m与130m，其余4个控制点深度随机分布，所有控制点电阻率值根据测区地质信息在1-1000Ω·m范围内进行超立方采样；
确定控制点信息后，将整个电阻率模型剖分为30层，每层具有固定层厚，第一层层厚为1m，剩余层厚按深度呈指数递增，且递增系数约为100.0336，相邻控制点间的各层电阻率参数由三次样条插值得到，则最终将产生包含30层且具有固定层厚的电阻率模型，示例模型如图2所示；
根据所获得的电阻率模型，基于瞬变电磁正演方法求取理论瞬变电磁响应；实际研究中发射系统包括2种模式：低磁矩和高磁矩，分别对应早期和晚期的数据采集，两者仅在接收时间上存在差异，其中理论瞬变电磁响应如图3所示；其中低磁矩模式对应接收时间道为5个，高磁矩模式对应接收时间道为15个，每个模型对应输入时间道基本一致。
S201：对S101步骤中产生的样本数据集进行数据清洗获取最终的样本训练集，并按比例划分为训练集、验证集和测试集；
按S101步骤得到理论瞬变电磁响应数据后，还需要对数据进行清洗以确保得到足够的样本训练集；由于地球物理反演的多解性，一组输入数据可能对应多个数据标签，对于构建输入数据和监督模型的端对端非线性映射关系产生了严重的制约，因此需要通过数据清洗筛除不符合研究背景的训练模型；本发明采用了逐变量的数据清洗方法，选定一组电阻率模型及瞬变电磁响应作为参照，若另一组电阻率模型与该电阻率模型的相对误差较大，但理论响应间的平均相对误差低于2%，则删除该模型，直到获取100000组训练样本；在该过程中，规定了模型阈值和理论响应阈值，实际分别取0.02和0.45，通过设定两个阈值的评价标准可以迅速筛选训练样本；另外，在进行CNN-LSTM网络训练前，按8:1:1的比例划分训练集、验证集和测试集。
S301：数据归一化处理；
本发明把两次正演的电磁响应拼接为维度20×1的向量输入神经网络，输出为30层的电阻率模型参数；数据集整体输入为一个2维数组，大小为电磁响应数据输入个数乘以单个数据的时间道数目；由于瞬变电磁数据会跨越数量级分布，无法直接输入网络进行训练，因而需要对数据进行标准化处理，一方面使其符合神经网络的输入条件，另一方面使其训练权重分布合理；通常激活函数的最优输入值为0~1之间，过大的输入值会使得梯度消失，从而让神经网络训练失去作用；其次，这种数据差异也会对网络性能泛化性产生严重影响，特别是早期较大的瞬变电磁响应会在损失函数中占据较大权重，从而使得网络无法有效学习晚期数据的变化特征，丧失对深层电阻率的识别能力；在实施例中，本发明对同时间道上的瞬变电磁响应进行归一化处理，采用下述计算公式：
其中，x(t)表示对应时刻的电磁响应值，为变换后的值，x max与x min分别为对应时间道上电磁响应的最大值与最小值；变换后的数据将尽可能分布在0~1之间，填充坐标空间，增大不同曲线之间的差异，从而使神经网络可以更好地区别不同时间道上的数据特征；
S401：构建CNN-LSTM网络框架；
本发明采用的网络框架如图4所示，网络在输入电磁响应数据后，利用2个卷积层提取数据特征，卷积核大小为（2×1），卷积核移动步长为1，每个卷积层的卷积核个数均为256个；经过2个卷积层抽取数据特征后，再将特征输入到LSTM层中进行时序信息处理，其中LSTM层共计2层，并设置随机失活（dropout）比例为0.2；最后通过DNN（深度或全连接神经网络的代称）层进行输出调控；上述中，每个神经网络层后面都会附加激活函数层，选用的激活函数为ReLU，其中DNN层后附加ReLU激活函数层可以确保输出的对数电阻率不会小于0；另外，随机失活是训练中常用的方法，通过随机使选中神经元的输出权重为0，减弱了神经元间的依赖性，一定程度上能防止网络训练过拟合；
S501：在CNN-LSTM网络训练的损失函数项中加入电阻率灵敏度约束；
常见的损失函数形式包括1范数和2范数，2范数即是指均方误差（MSE），损失函数可以使每层的电阻率值获得相同权重，通过不断的误差反向传播实现严格的监督学习；事实上，瞬变电磁法的勘探分辨率随深度增加而减弱，特别是当深度达到某个值以后，产生的晚期电磁响应几乎一致，也意味着反演方法丧失了对于该深度以下地层电阻率的分辨能力；通过向损失函数中对电阻率模型在不同深度加入不同的灵敏度约束，根据瞬变电磁数据对CNN-LSTM网络的灵敏度为电阻率模型分配权重，使网络学习的反演机制更符合瞬变电磁响应的衰减规律，在实施例中，构建相应约束项的过程如下：
其中，是归一化权重，它表示各层灵敏度权重值与和的比值，系数/>为数据集标签中理论观测瞬变电磁响应数据点上的噪声，其中，单层灵敏度计算如下：/>
其中，为第i个数据集标签中理论观测瞬变电磁响应；/>表示第j个层状电阻率模型参数，/>为对应电阻率灵敏度计算值；则最终损失函数表示为：
其中，是样本数量，/>是第/>个样本的真实电阻率参数，/>对第/>个样本的预测电阻率参数，/>是正则化参数，通常在电磁反演中是调节正则化项在反演中的比重，在本发明中取固定值为0.05，实际数值仍然可以根据测试进行调整；优化算法采用Adam优化器来最小化损失函数，网络训练的超参数初始化如下：起始学习率为0.001，数据批尺寸为128，训练最大迭代次数为2000；训练策略采用早停机制，若验证集误差在连续迭代100次后仍未下降，则停止网络训练，防止网络对训练集数据过拟合；
S601：理论模型测试与网络训练及评价；
在训练过程中，利用平均相对误差来评价预测结果，它反映了预测值与真实值之间的偏差；在训练过程中，验证集和训练集需要同时导入，其中，训练集用于训练网络模型，验证集是用于实时评价训练网络的有效性；其中，网络训练是按批次导入数据，每次导入大小为一个数据批尺寸，即是128，直至所有训练集数据导入训练，该过程中样本导入是随机的，其中一个循环也就是遍历所有训练样本表示训练1次；网络评价是将训练集数据输入由训练集构建的临时网络模型，可以预测相应电阻率参数，从而计算验证集的损失函数值；若验证集的损失函数值在训练过程中整体趋势没有回升现象，并和训练集的损失函数都呈下降趋势，则表明网络没有过拟合；
其次，网络训练停止标准也是获取优良网络模型的关键，在该过程中选择了早停机制，即在损失函数值极低的情况下，若训练集的误差在100次训练后都没有下降，则停止训练过程，并寻找最优神经网络模型；本发明选择网络在第1150次迭代后的模型对测试数据进行反演；
为了测试网络的泛化性，构建了一个三维模型，黑色测点及白色线框异常体模型分布如图5（a）和5（b）所示，其中，沿X方向，点间距为15米，沿Y方向点间距为20米；异常体分布位置如图5（a）的Y=0方向切片所示，带有垂直黑线的区块即异常体单元的电阻率为10欧姆·米，背景即白色区域的电阻率为50欧姆·米；接着对每个测点的结果进行反演，然后拼接为三维模型，设定电阻率小于15欧姆·米的反演结果显示；图5（b）进一步给出了神经网络的反演结果，结果表明白色线框和包围的实体形态接近，位置上比较重合，由此反映了所开发的神经网络模型的有效性；
S701：对实测数据进行处理，获取反演电阻率模型；
为了进一步验证算法对实测数据的适用性，采用了瑞士Aare Valley所采集的地面拖拽式瞬变电磁数据（tTEM20AAR）进行了验证；该数据中，测线距离一般在10-25m之间，测线内的点距一般在3-10 m之间，共计包括37957个测点；图6展示了对所有测点的拟合误差统计，结果表明超过97%的结果平均相对拟合误差低于10%，绝大部分低于6%，表明获取的结果是可靠的；
S801：构建区域离散正交网格；
为了将上述反演的电阻率模型转化为对区域的地质解释，首先需要构建离散正交网格；在实施例中，统计了所有点墨卡托投影坐标的最值，这里用xy坐标代称；根据xy坐标的最小值，对其向下取整，然后按照网格单元边长300米，构建覆盖所有测点的平面正交网格；其次，在垂向上，考虑到系统的有效勘探深度80米，因此在实际岩性解释中只考虑了距离表面100米以内的介质；对应该过程，设置z方向单元长度4米，将纵向剖分为25段，则最大剖分深度即是100米，由此构建了1个三维正交网格；事实上，纵向剖分长度可根据实际研究变化，也可以是不均匀分布；由于最优化插值是针对单层水平网格，这里不考虑高程或地形因素，例如仅从网格角度来分析，单层是指都距离地下4米或其它固定剖分深度下的曲面；若不考虑表层地形因素，则表层的测点以及离散网格将如图7所示；当正方形单元网格区域中含有地球物理测点，则正方形单元上的4个节点即为有效节点，基于该准则，可以筛选出所有有效节点，相应的，其它剩余网格节点则为无效节点，它们将不参与最优化插值过程；
S901：对测井数据进行解译，形成最优化插值的观测数据；
尽管拥有超过3000口测井的资料，但是一些测井资料的时间、位置以及深度都不符合预期，最终筛选了66口测井的资料，作为进行岩性转换解释的依据，图7显示了筛选测井的位置分布情况；在实施中，研究区域近地表主要为第四系的碎屑岩沉积物，通常是粘土、砂等；由于电阻率以及测井等难以详细地区分岩性的具体构成，在水文研究中会假设沉积物的构成分为粘土和砂2类，其中粘土的电阻率较低，砂的电阻率较高；这个假设可以用下式描述：
在实施例中，先解译钻井剖面，得到对应的层厚度以及粘土含量等属性数据，并记录相关的地下水统计资料；基于S901构建的正交网格，需要对同一深度水平层的粘土体积分数进行计算，其过程遵循如下准则：
其中，是计算地层厚度对应的粘土综合体积分数，/>是计算地层厚度，/>是不同岩层对应厚度，/>则是该层对应的粘土体积分数。依照上述计算准则，可以得到不同深度段的粘土体积分数。接着采用置信度算法对其进行矫正，消除人为主观因素影响，其中对置信度/>的定义如下：
其中M表示距离置信度，N表示差异置信度。通常浅层的岩性是渐变的，和沉积作用相关，因为沉积过程通常是稳定的，即同一时间的水力环境应该是一致的；基于这个前提，假设500、1000、1500米内的粘土含量存在变化阈值，该值在实施例中规定为0.05、0.1、0.15，其中0.05表示对应500米内同一深度层的粘土含量差异不会大于5%，以上假设可以根据实际需求调整；
根据假设，在同一深度层上，首先确定一个目标测井，然后搜索最大定义半径内的其它测井，本发明实施例中其值为1500米；还需要确定排列测井与目标测井的体积分数差异是否满足假设，然后计算并列举所有不满足条件的测井，并按距离排列。选择距离最近的测井记为修正测井，如果不存在这样的情况，将跳过执行下一个测井；当筛选出目标测井和修正测井后，需要计算其对应的置信度，需要重复上述步骤，需再找出修正测井周围的在最大定义半径内的其它排列测井；
本发明对置信度解释如下，对于排列测井，不再考虑修正测井因素，同样对于修正测井而言，计算其置信度时也不考虑目标测井而只考虑其它排列测井，排列测井是按照最大定义半径搜索得到，意味着目标测井和修正测井的排列测井是不相同的；在计算置信度时，以目标测井为例，距离越近意味其粘土体积分数信息越可信，则置信度M值也就越高；另外按照假设，若排列测井对应的粘土体积分数与目标测井的粘土体积分数差异符合假设，则其置信度N的值将越大；通过计算所有排列测井的置信度，将其累加即得到目标测井的置信度L1，同样按照上述方法也可以得到修正测井的置信度L2；接着按照目标测井与修正测井的距离L0，先判断它们位于哪级，若在500米内，则按照0.05/500的系数乘以两者距离差L0变化因素修正，修正值为0.05*L0/500，同样其它情况下可依据假设依次计算，得到其修正值并记为X0，若两者之间的粘土体积分数差记为N0，则修正公式如下：
其中，和/>分别表示修改前的目标测井和修正测井的粘土体积分数含量，/>和/>表示修改后的值，/>表示若/>，则/>对应式子取/>，/>则取/>，相反则符号相应交换；重复上述过程直到所有测井满足假设即遍历完成；
S1001：利用高斯-牛顿方法实现最优化插值，构建砂/砂砾岩石物理模型；
本发明将电阻率转换为粘土体积分数，采用了下述转换函数：
/>
其中，和/>对应电阻率阈值并返回0.025和0.975的对应权重，erfc是误差函数，具体表示如下：
另外，对于最优化插值的理解包括：事实上，观测粘土体积分数和实际插值的粘土体积分数满足如下关系：
其中，是核矩阵，实际上是计算粘土体积分数/>和模型参数/>之间的矩阵映射，/>表示/>和真实粘土体积分数/>之间的差异；本实施例选用高斯-牛顿反演算法去优化上述过程，使得插值结果与测井解译结果趋于一致，从而获取最终模型参数/>，这里的模型参数/>包含电阻率阈值/>和/>，指的是有效节点电阻率阈值；其中高斯-牛顿反演算法的灵敏度矩阵求解采用扰动法，即通过微小扰动求解插值粘土体积分数对模型参数的差商，用其代替导数；
对于该过程插值的说明包括：首先通过给定有效节点的初始电阻率阈值，即是给定和/>，在后续将用电阻率阈值对其进行代称。利用有效节点的电阻率阈值进行双线性插值得到地球物理测点的电阻率阈值，有效节点参照S801的定义说明，其中双线性插值基于正方形单元，有效节点即为正方形单元节点；计算电阻率模型的粘土体积分数，计算方式为S901的给定的计算等式；以筛选的测井位置为圆心，设置1000米的搜索半径，并统计区域内的地球物理测点，搜索半径可适当调整，建议不超过S901的最大定义半径；在搜索半径内，利用克里金插值方法，基于统计的地球物理测点的粘土体积分数插值得到测井位置的粘土体积分数，以上过程即是一次正演过程；
对优化的目标函数说明如下：
其中，是观测测井数目，/>是数据拟合项，/>是约束对数目，/>和分别表示第i个测井对应的插值粘土体积分数和观测粘土体积分数，/>表示数据方差，/>和/>表示相邻节点的粘土体积分数，/>为其调控变化的测度常数；
选择优化和优化/>是一致的，通过优化目标函数，可以得到研究区域的粘土含量分布，然后基于粘土和砂总量为1的假设，可以直接得到区域的砂/砂砾分布，从而构建砂/砂砾岩石物理模型；
S1101：构建solid模型；
GMS允许多种数据导入形式， solid建模要求数据导入形式为测井。事实上，经过转换后的数据对应测点等效37957口测井，数据不方便直接导入GMS，因此需要简化；基于岩性渐变的假设，再次将区域剖分为单元长度为100米的正交网格，剖分方式与S801所述一致；依次统计所有正方形单元内的地球物理测点，并求其砂/砂砾含量数值的平均，同时把平均值赋予正方形单元中心位置，视为它的砂/砂砾含量；由此本实施例成功将37957个测点数据简化为1212个中心点所表征的数据；另外，考虑到每个单元的砂岩含量，本实施例设定了4种不同属性的填充材料，以匹配每个单元格的砂/砾石含量；在0~1的分布区间，设定了0.75、0.5、0.25等三个截断值，因此4种材料将分别对应不同的砂/砾石体积分数分布；应当注意的是，构建solid模型时，可以采用GMS自动识别功能，但需要严格检验Horizons ID的对应性，尤其是复杂数据构成的Horizons ID并不完全准确；通常GMS实体建模的依据是成层分布，意味着局部建模可能没有精准匹配对应层岩性的分布；为解决以上问题，综合考虑了输入数据的层分布特性，本实施例自行编码获取了相对准确的Horizons ID，也可选择手动调整，同时在局部建模中增加水平分布约束，约束局部岩层的展布形态从而获得solid模型；
S1201：基于构建solid模型，利用MODFLOW计算地下水位分布；
利用构建的GMS solid模型，在GMS创建MODFLOW 计算模块；根据实际地质特征，在GMS新建图层进行定义；本实施例中需定义河流属性，并计算实体模型对应的水位线分布；本发明首先圈定了研究区域，采用了均匀网格覆盖所创建的实体模型，实际模型在水平面的空间大小分布为8500m12200 m；为了使初始网格尽可能囊括实体模型，分别在水平方向设置了95/>130个数量的单元格，每个网格单元水平面大小为100m/>100 m；同时在z轴方向，设置起始剖分深度为405m，区域计算高度为150m，以对应导入数据的底层最低高程以及顶部最大高程；实际实体模型表明应至少划分12层才可以很好地与各岩层对应，因此在实际模型中设置了15层，每个垂向网格单元长度为10m。当将solid模型导入MODFLOW模块中，实际上每层厚度为7m；另外，为保证薄层计算的润湿性，对各固体层的剖分设置了2m最小离散层厚；同时在MODFLOW的Conceptual Model模块中的图层定义了层计算属性，一个重要的设置是将流经研究区域的Aare河设为固定水头；此外还在LPF Package导入材料属性，Material1~4的垂直渗透系数分别定义为30、15、5、1 m/d，水平渗透系数依靠相应的比值设定，设置初始Vert anisotropy 为3，孔隙度则默认为0.3；上述参数将根据估计的表层水位分布进行细微调整；
其中，计算水位与实际水位对比如图8（a）和图8（b）所示，图8（a）给出了标记测井位置，获取的地下水位统计信息来自2011年后的调查资料，拖曳式瞬变电磁数据的调查时间是2020年，尽管在时间上不统一，严格而言地下水位基本是接近的；图8（b）给出了对比结果，水位差异最大在3米以内，经统计发现这类时间通常较早；另外计算结果总体是一致的，表明了计算方法的可行性，由于缺少更详细的测井信息，采用模糊计算方法造成的偏差总体是可以接受的；由此给出了基于地球物理电磁数据到计算水位的一般方法。
本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。
Claims (10)
1.一种基于瞬变电磁法成像的地下水模拟方法，其特征在于，包括以下步骤：
步骤一：采用实际瞬变电磁响应数据，利用训练完毕的CNN-LSTM网络模型预测地下的电阻率；
步骤二：利用研究区域内的测井资料，结合反演的地下电阻率，采用最优化插值方法，依次计算区域内不同测点不同深度层对应的砂和砂砾的含量分布；
步骤三：将研究区域内不同测点不同深度对应的砂和砂砾的含量分布数据进行简化，并将简化的数据以测井数据形式导入GMS构建solid模型；
步骤四：将solid模型导入MODFLOW模块，结合有限的地质信息设置计算参数，计算区域地下水水位分布；
其中，训练CNN-LSTM网络模型的方法为：
基于连续层状电阻率模型，结合已知的发射装置参数，通过瞬变电磁正演方法获取理论的瞬变电磁响应，形成用于训练CNN-LSTM网络的初始样本数据集，训练时CNN-LSTM网络的损失函数采用均方误差，同时在损失函数中添加由电阻率灵敏度约束构建的正则化项，输出为电阻率参数。
2.根据权利要求1所述的地下水模拟方法，其特征在于，训练CNN-LSTM网络模型的具体方法，包括以下步骤：
基于构建的固定层数的、固定层厚的连续层状电阻率模型，结合已知的发射装置参数，通过瞬变电磁正演方法计算对应时刻的瞬变电磁响应，形成初始样本数据集；
利用数据清洗方法清除初始样本数据集中的冗余数据，按一定比例划分为训练集、验证集和测试集；
构建包含电阻率灵敏度约束的损失函数，采用训练集和验证集训练CNN-LSTM网络模型；
基于早停机制和测试集测试选择最优化CNN-LSTM网络模型；
其中，包含电阻率灵敏度约束的损失函数为：
其中，是正则化参数，MSE表示均方误差；/>为归一化后的灵敏度约束权重；
其中，是样本数量，/>是第/>个样本的真实电阻率参数，/>是对第/>个样本的预测电阻率参数；
其中，为数据集标签中理论观测瞬变电磁响应数据点上的噪声；
其中，为第i个数据集标签中理论观测瞬变电磁响应；/>表示第j个层状电阻率模型参数，/>为对应电阻率灵敏度计算值。
3.根据权利要求1或2所述的地下水模拟方法，其特征在于，最优化反演方法，包括以下步骤：
统计研究区域的岩性测井资料，以及测井及地球物理测点的空间位置，同时在研究区域表层构建均匀、离散正交的控制网格，并筛选出包围有地球物理测点的有效网格节点，在纵向上结合研究区域地质条件，构建均匀或具有不同层厚的纵向网格分布；
对岩性测井资料进行解译得到对应层的粘土体积分数及层厚度，并通过置信度算法对解译的粘土体积分布及层厚度进行调整，基于层厚度加权方式，计算岩性测井资料对应层的粘土体积分数作为最优化插值中的观测数据，同时设置有效网格节点上的初始电阻率阈值；
对于相同埋深的深度层，每次反演时采用双线性插值方法，由分布在单个正方形单元周围四个网格节点的初始电阻率阈值，计算包围电阻率模型的电阻率阈值，进而计算单个电阻率模型对应层的粘土体积分数；
以测井位置点为圆心，设置搜索半径，在搜索半径内对电阻率模型转换得到的粘土体积分数进行克里金插值，得到测井位置的粘土体积分数插值结果；
基于测井位置的粘土体积分数插值结果，采用高斯-牛顿最优化反演方法优化插值过程，依次计算不同深度水平层对应的粘土体积分数分布，基于粘土和砂总含量为1的假设可以得到砂和砂砾的含量分布。
4.根据权利要求3所述的地下水模拟方法，其特征在于，电阻率转换为粘土体积分数采用的转换函数为：
其中，和/>分别为电阻率的最高阈值和电阻率的最低阈值；erfc为误差函数；/>为电阻率。
5.根据权利要求1所述的地下水模拟方法，其特征在于，步骤三具体包括以下步骤：
对研究区域进行有效分区，并选择区域中心点作为区域所有结果的表征，从而对不同测点不同深度层对应的砂和砂砾的含量分布进行简化，其中，中心点数值为区域内所有对应深度测点砂和砂砾体积分数的平均；
根据每个区域的砂岩含量，设置三个砂和砾石体积分数分布截断值将区域材料分为四种不同属性的填充材料，对电阻率模型的不同层厚进行材料标记；
将简化后的砂和砂砾的含量分布数据以及材料标记以导入测井资料的形式导入GMS，进行Horizons ID智能识别，并判断识别的准确性，若识别有误则采用自行编码，并对相应的Horizons ID进行更正，在局部建模中增加对对应层的水平分布约束，从而构建GMSsolid实体模型。
6.根据权利要求1所述的地下水模拟方法，其特征在于计算区域地下水分布的方法，包括以下步骤：
将solid模型导入新建的MODFLOW模块，并创建离散网格；
利用测井资料及相关信息定义填充材料属性，其中，填充材料属性为岩块的各种水力参数；
在GMS内新建图层，定义研究区域、流动边界以及各种水文属性，并设置有效计算条件；
对MODFLOW中的输入参数、设定的边界条件进行检验，进行计算并得到研究区域水位分布。
7.一种基于瞬变电磁法成像的地下水模拟系统，其特征在于，包括：
地下电阻率预测模块，用于采用实际瞬变电磁响应数据，利用训练完毕的CNN-LSTM网络模型预测地下的电阻率；
砂和砂砾的含量分布计算模块，用于利用研究区域内的测井资料，结合反演的地下电阻率，采用最优化插值方法，依次计算区域内不同测点不同深度层对应的砂和砂砾的含量分布；
solid模型的构建模块，用于将研究区域内不同测点不同深度对应的砂和砂砾的含量分布数据进行简化，并将简化的数据以测井数据形式导入GMS构建solid模型；
地下水水位分布计算模块，用于将solid模型导入MODFLOW模块，结合有限的地质信息设置计算参数，计算区域地下水水位分布；
其中，训练CNN-LSTM网络模型的方法为：
基于连续层状电阻率模型，结合已知的发射装置参数，利用瞬变电磁正演方法获取理论的瞬变电磁响应，形成用于训练CNN-LSTM网络的初始样本数据集，训练时CNN-LSTM网络的损失函数采用均方误差，同时在损失函数中添加由电阻率灵敏度约束构建的正则化项，输出为电阻率参数。
8.根据权利要求7所述的地下水模拟系统，其特征在于，包含电阻率灵敏度约束的损失函数为：
其中，是正则化参数，MSE表示均方误差；/>为归一化后的灵敏度约束权重；
其中，是样本数量，/>是第/>个样本的真实电阻率参数，/>是对第/>个样本的预测电阻率参数；
其中，为数据集标签中理论观测瞬变电磁响应数据点上的噪声；
其中，为第i个数据集标签中理论观测瞬变电磁响应；/>表示第j个层状电阻率模型参数，/>为对应电阻率灵敏度计算值。
9.根据权利要求7或8所述的地下水模拟系统，其特征在于，所述砂和砂砾的含量分布计算模块包括：
网格构建单元，用于统计研究区域的岩性测井资料，以及测井及地球物理测点的空间位置，同时在研究区域表层构建均匀、离散正交的控制网格，并筛选出包围有地球物理测点的有效网格节点，在纵向上结合研究区域地质条件，构建均匀或具有不同层厚的纵向网格分布；
初始电阻率阈值的设置单元，用于设置有效网格节点上的初始电阻率阈值；
观测粘土体积分数的计算单元，用于对岩性测井资料进行解译得到对应层的粘土体积分数及层厚度，并通过置信度算法对解译的粘土体积分布及层厚度进行调整，基于层厚度加权方式，计算岩性测井资料对应层的粘土体积分数作为最优化插值中的观测数据；
相同埋深粘土体积分数的计算单元，用于基于纵向网格分布中具有相同埋深的深度层，每次反演时采用双线性插值方法，由分布在单个正方形单元周围四个网格节点的初始电阻率阈值，计算包围电阻率模型的电阻率阈值，进而计算单个电阻率模型对应层的粘土体积分数；
搜索半径内粘土体积分数的插值单元，用于以测井位置点为圆心，设置搜索半径，在搜索半径内对电阻率模型转换得到的粘土体积分数进行克里金插值，得到测井位置的粘土体积分数插值结果；
不同埋深粘土体积分数的计算单元，用于基于测井位置的粘土体积分数插值结果，采用高斯-牛顿最优化反演方法优化插值过程，依次计算不同深度水平层对应的粘土体积分数分布，基于粘土和砂总含量为1的假设可以得到砂和砂砾的含量分布。
10.根据权利要求9所述的地下水模拟系统，其特征在于，电阻率转换为粘土体积分数采用的转换函数为：
其中，和/>分别为电阻率的最高阈值和电阻率的最低阈值；erfc为误差函数；/>为电阻率。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311137070.7A CN116859478B (zh) | 2023-09-05 | 2023-09-05 | 一种基于瞬变电磁法成像的地下水模拟方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311137070.7A CN116859478B (zh) | 2023-09-05 | 2023-09-05 | 一种基于瞬变电磁法成像的地下水模拟方法及系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN116859478A CN116859478A (zh) | 2023-10-10 |
| CN116859478B true CN116859478B (zh) | 2023-11-28 |
Family
ID=88225408
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202311137070.7A Active CN116859478B (zh) | 2023-09-05 | 2023-09-05 | 一种基于瞬变电磁法成像的地下水模拟方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116859478B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117192628B (zh) * | 2023-11-03 | 2024-01-26 | 北京科技大学 | 深部破碎含水地层分布识别方法 |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CA2650105A1 (en) * | 2006-05-04 | 2007-11-15 | Exxonmobil Upstream Research Company | Time lapse analysis with electromagnetic data |
| RU2717740C1 (ru) * | 2019-05-28 | 2020-03-25 | Федеральное государственное бюджетное учреждение науки Институт физики Земли им. О.Ю. Шмидта Российской академии наук | Способ прогноза открытой пористости в пространстве между скважинами |
| CN113568056A (zh) * | 2021-07-16 | 2021-10-29 | 北京大学 | 一种基于卷积神经网络的航空瞬变电磁数据反演方法 |
| CN113568055A (zh) * | 2021-07-16 | 2021-10-29 | 北京大学 | 一种基于lstm网络的航空瞬变电磁数据反演方法 |
| CN114386464A (zh) * | 2022-01-12 | 2022-04-22 | 中国科学院地质与地球物理研究所 | 一种瞬变电磁激电信息的深度学习提取方法 |
| CN115508900A (zh) * | 2022-11-23 | 2022-12-23 | 中国地质大学（武汉） | 一种地面拖拽式瞬变电磁成像方法及系统 |
| CN115826059A (zh) * | 2022-11-30 | 2023-03-21 | 电子科技大学长三角研究院(湖州) | 一种瞬变电磁快速超前预测方法 |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US10877178B2 (en) * | 2019-05-03 | 2020-12-29 | Institute Of Geology And Geophysics, Chinese Academy Of Sciences | Method for suppressing airborne transient electromagnetic in-band vibration noise |
| US11782183B2 (en) * | 2021-04-26 | 2023-10-10 | Institute Of Geology And Geophysics, Chinese Academy Of Sciences | Magnetotelluric inversion method based on fully convolutional neural network |
-
2023
- 2023-09-05 CN CN202311137070.7A patent/CN116859478B/zh active Active
Patent Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CA2650105A1 (en) * | 2006-05-04 | 2007-11-15 | Exxonmobil Upstream Research Company | Time lapse analysis with electromagnetic data |
| RU2717740C1 (ru) * | 2019-05-28 | 2020-03-25 | Федеральное государственное бюджетное учреждение науки Институт физики Земли им. О.Ю. Шмидта Российской академии наук | Способ прогноза открытой пористости в пространстве между скважинами |
| CN113568056A (zh) * | 2021-07-16 | 2021-10-29 | 北京大学 | 一种基于卷积神经网络的航空瞬变电磁数据反演方法 |
| CN113568055A (zh) * | 2021-07-16 | 2021-10-29 | 北京大学 | 一种基于lstm网络的航空瞬变电磁数据反演方法 |
| CN114386464A (zh) * | 2022-01-12 | 2022-04-22 | 中国科学院地质与地球物理研究所 | 一种瞬变电磁激电信息的深度学习提取方法 |
| CN115508900A (zh) * | 2022-11-23 | 2022-12-23 | 中国地质大学（武汉） | 一种地面拖拽式瞬变电磁成像方法及系统 |
| CN115826059A (zh) * | 2022-11-30 | 2023-03-21 | 电子科技大学长三角研究院(湖州) | 一种瞬变电磁快速超前预测方法 |
Non-Patent Citations (3)
| Title |
|---|
| Target Electromagnetic Detection Method in Underground Environment: A Review;Xiaofen Wang 等;《IEEE SENSORS JOURNAL》;第22卷(第14期);13835-13852 * |
| 基于 Attention 机制的 CNN-BiLSTM 瞬变电磁实时反演方法;古瑶 等;《煤田地质与勘探》;第51卷(第8期);1−11 * |
| 基于深度学习的地面拖曳式瞬变电磁快速成像方法;冼锦炽 等;《工程地球物理学报》;第19卷(第4期);536-545 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN116859478A (zh) | 2023-10-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US10519766B2 (en) | Reservoir modelling with multiple point statistics from a non-stationary training image | |
| US20170011149A1 (en) | Geobody continuity in geological models based on multiple point statistics | |
| CN106605158B (zh) | 利用用于深度剖绘的参数化模板的沉积物运移模拟 | |
| Koch et al. | Challenges in conditioning a stochastic geological model of a heterogeneous glacial aquifer to a comprehensive soft data set | |
| CN111596978A (zh) | 用人工智能进行岩相分类的网页显示方法、模块和系统 | |
| CN116859478B (zh) | 一种基于瞬变电磁法成像的地下水模拟方法及系统 | |
| Song et al. | GANSim‐3D for conditional geomodeling: Theory and field application | |
| Ouadfeul et al. | Lithofacies classification using the multilayer perceptron and the self-organizing neural networks | |
| CN114723095A (zh) | 缺失测井曲线预测方法及装置 | |
| CN115659848B (zh) | 一种基于深度学习网络快速预测二维盆地基底界面的方法 | |
| Vasco et al. | Nonuniqueness in traveltime tomography: Ensemble inference and cluster analysis | |
| CN115508900A (zh) | 一种地面拖拽式瞬变电磁成像方法及系统 | |
| Straubhaar et al. | Conditioning Multiple‐Point Statistics Simulation to Inequality Data | |
| CN112100906B (zh) | 数据驱动的大尺度密度建模方法、计算设备及存储介质 | |
| Yanshu et al. | A three-dimensional model of deep-water turbidity channel in Plutonio oilfield, Angola: From training image generation, optimization to multi-point geostatistical modelling | |
| Baouche et al. | Intelligent methods for predicting nuclear magnetic resonance of porosity and permeability by conventional well-logs: a case study of Saharan field | |
| CN103543478A (zh) | 地质形态插值的km方法 | |
| Comunian et al. | Hybrid inversion method to estimate hydraulic transmissivity by combining multiple-point statistics and a direct inversion method | |
| Zare et al. | Reservoir facies and porosity modeling using seismic data and well logs by geostatistical simulation in an oil field | |
| Pan et al. | Efficient subsurface modeling with sequential patch generative adversarial neural networks | |
| CN114707597A (zh) | 一种河流相致密砂岩储层复杂岩相智能化识别方法及系统 | |
| CN115880455A (zh) | 基于深度学习的三维智能插值方法 | |
| CN113253342A (zh) | 联合测井数据构建复杂三维地震模型标签方法 | |
| Jo et al. | Conditioning stratigraphic, rule-Based models with generative adversarial networks: A deepwater lobe, deep learning example | |
| Neven et al. | Stochastic multi-fidelity joint hydrogeophysical inversion of consistent geological models |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |