CN117270072B - 基于改进型差分进化算法的重磁位场成像反演方法及系统 - Google Patents
基于改进型差分进化算法的重磁位场成像反演方法及系统 Download PDFInfo
- Publication number
- CN117270072B CN117270072B CN202311209659.3A CN202311209659A CN117270072B CN 117270072 B CN117270072 B CN 117270072B CN 202311209659 A CN202311209659 A CN 202311209659A CN 117270072 B CN117270072 B CN 117270072B
- Authority
- CN
- China
- Prior art keywords
- model
- inversion
- objective function
- data
- vector
- 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 71
- 238000003384 imaging method Methods 0.000 title claims abstract description 69
- 230000005484 gravity Effects 0.000 title claims abstract description 37
- 239000013598 vector Substances 0.000 claims description 85
- 239000011159 matrix material Substances 0.000 claims description 58
- 230000000704 physical effect Effects 0.000 claims description 9
- 238000012360 testing method Methods 0.000 claims description 8
- 230000035772 mutation Effects 0.000 claims description 7
- 238000010276 construction Methods 0.000 claims description 5
- 238000009933 burial Methods 0.000 claims description 4
- 238000005457 optimization Methods 0.000 abstract description 19
- 230000007547 defect Effects 0.000 abstract description 11
- 230000008569 process Effects 0.000 abstract description 8
- 230000000694 effects Effects 0.000 abstract description 7
- 238000013473 artificial intelligence Methods 0.000 abstract description 6
- 230000007246 mechanism Effects 0.000 abstract description 5
- 230000008901 benefit Effects 0.000 description 6
- 239000010977 jade Substances 0.000 description 5
- 238000009499 grossing Methods 0.000 description 4
- 230000009471 action Effects 0.000 description 3
- 238000013461 design Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 238000004364 calculation method Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000002159 abnormal effect Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000013401 experimental design Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 230000002068 genetic effect Effects 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 229910052500 inorganic mineral Inorganic materials 0.000 description 1
- 230000005389 magnetism Effects 0.000 description 1
- 239000011707 mineral Substances 0.000 description 1
- 239000002245 particle Substances 0.000 description 1
- 230000001575 pathological effect Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/12—Computing arrangements based on biological models using genetic models
- G06N3/126—Evolutionary algorithms, e.g. genetic algorithms or genetic programming
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V11/00—Prospecting or detecting by methods combining techniques covered by two or more of main groups G01V1/00 - G01V9/00
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- 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
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Life Sciences & Earth Sciences (AREA)
- Pure & Applied Mathematics (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Biology (AREA)
- Biophysics (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Health & Medical Sciences (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Computing Systems (AREA)
- Physiology (AREA)
- Molecular Biology (AREA)
- General Health & Medical Sciences (AREA)
- Evolutionary Computation (AREA)
- Computational Linguistics (AREA)
- Biomedical Technology (AREA)
- Artificial Intelligence (AREA)
- Genetics & Genomics (AREA)
- Operations Research (AREA)
- Probability & Statistics with Applications (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本申请提供一种基于改进型差分进化算法的重磁位场成像反演方法及系统,通过对重磁位场数据去噪、分离区域场,得到剩余场数据,生成地下介质的初始模型和物性参数范围约束,读入剩余场数据和初始模型,生成反演网格参数和先验信息矩阵,建立成像反演的反演目标函数,利用改进型差分进化算法优化反演目标函数,最终输出数据目标函数最小的模型向量作为最终的反演结果。以人工智能领域的差分进化方法为基础,建立一套适合差分进化的重磁位场的成像反演机制,能有效改进差分进化的反演效果,提升该方法的收敛速度,弥补全局方法成像反演的不足(成像结果包含大量噪声、搜索过程耗时较长等问题),增强全局反演的实用性,弥补局部优化反演的不足。
Description
技术领域
本申请涉及勘探地球物理领域,具体而言,涉及一种基于改进型差分进化算法的重磁位场成像反演方法及系统。
背景技术
重力和磁勘探作为两种经典的地球物理方法,已在区域地质调查、油气、矿产资源勘查等方面取得了广泛的应用。重磁位场的成像反演问题具有病态性,即存在多个不同的模型与观测数据吻合。因此,在解决重磁的成像反演问题的过程中,需要引入正则化约束技术。此时,反演的目标函数由两个部分构成:数据目标函数和模型目标函数。通过施加正则因子,数据目标函数和模型目标函数形成反演的目标函数。在重磁位场的反演目标函数中,数据与模型之间的关系是强非线性的,一般的局部优化方法依赖于初始模型的选取,容易陷入局部最优模型。
目前,重磁位场的成像反演主要通过局部优化方法来解决,如拟牛顿、共轭梯度、最速下降等。而基于人工智能思想的全局优化方法(如遗传算法、粒子群算法)主要用于重磁位场的场源参数的反演(如埋深、倾向、界面等),极少用于重磁位场的成像反演。
差分进化方法是一种基于种群的随机搜索方法,其是由Rainer Storn和KennethPrice共同提出。差分进化算法具有强的全局探索性,能在一定程度上保证所获得的模型是局部最优。差分进化具有结构简单、收敛速度快、鲁棒性强等优点。
然而,重磁位场的成像反演问题的维度较高,常在100维以上。如果要利用差分进化算法对这类高维复杂问题进行优化,最主要的问题是收敛速度慢和最终解的精度不足。此外,差分进化算法主要用于解决单目标优化问题,并不能直接用于求解重磁成像反演问题。
发明内容
本申请实施例的目的在于提供一种基于改进型差分进化算法的重磁位场成像反演方法及系统,以人工智能领域的差分进化方法为基础,建立一套适合差分进化的重磁位场的成像反演机制,能有效改进差分进化的反演效果,提升该方法的收敛速度,获得有效的地下物理介质的分布情况。
为了实现上述目的,本申请的实施例通过如下方式实现:
第一方面,本申请实施例提供一种基于改进型差分进化算法的重磁位场成像反演方法,包括:
S1:获取地下介质的重磁位场数据,并对重磁位场数据去噪、分离区域场,得到剩余场数据;
S2:生成地下介质的初始模型和物性参数范围约束;
S3:读入剩余场数据和初始模型;
S4:基于剩余场数据和初始模型,生成反演网格参数和先验信息矩阵;
S5:生成平滑矩阵S,并建立成像反演的反演目标函数,其中,反演目标函数包含用于衡量观测数据和反演预测数据的差异的数据目标函数和用于约束反演的模型目标函数;
S6:利用改进型差分进化算法优化反演目标函数,最终输出数据目标函数最小的模型向量作为最终的反演结果,反演结果揭示地下介质的密度或磁化率分布。
结合第一方面,在第一方面的第一种可能的实现方式中,先验信息矩阵包括:
深度加权矩阵Wz:
其中,i为单元编号,M为待反演的物性参数个数,zi和Vi分别表示单元i的中心埋深和体积大小,β为深度衰减因子,对重力数据,β取1,对磁数据,β取2;
数据权重矩阵Wd:
Wd=diag(Wd,1,Wd,2,…,Wd,i,…,Wd,N),
其中,N为观测数据个数;
Wd,i满足:
其中,为单元i的观测数据向量,std表示标准差,dobs为观测数据向量。
结合第一方面的第一种可能的实现方式,在第一方面的第二种可能的实现方式中,反演目标函数为:
Φ(m)=Φd(m)+λΦm(m),
其中,Φ(m)为反演目标函数,m为由M个待优化模型参数形成的向量,记为模型向量m,m=(m1,m2,…,mM)T,λ为正则因子,Φd(m)为数据目标函数,满足:
其中,F(m)为重磁位场正演函数;
Φm(m)为模型目标函数,满足:
其中,Wm为模型参数加权矩阵。
结合第一方面的第二种可能的实现方式,在第一方面的第三种可能的实现方式中,S6中,利用改进型差分进化算法优化反演目标函数,包括:
S61:迭代次数G=0时,设置初始交叉概率均值尺度因子的初始位置参数/>基于初始模型,生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵/>并初始化数据目标函数/>模型目标函数/>正则因子λ(0),计算数据目标函数的初始均值Mean(Φd)(0),其中,由NP个模型向量形成的矩阵记为种群P,P=(m1,m2,…,mNP);
S62:针对种群P中的每一组模型向量:生成相应的交叉概率CR和尺度系数F,并计算反演目标函数值表示模型向量mj在第G次迭代时的反演目标函数值;
S63:使用改进型变异策略生成变异向量v(G),并通过交叉生成新试验模型u(G);
S64:通过选择操作,将符合要求的模型存入新的种群P(G+1),对于模型向量mj:
j∈[1,NP],
对于种群P(G)中被淘汰的模型向量,将该模型向量及相应的目标函数值存入档案A中,档案A的大小为NP×M,同时,分别存储相应的尺度系数F、交叉概率CR至集合S F和集合SCR,并存储至集合SΦ;
S65:记录数据目标函数均值Mean(Φd)(G+1),更新差分进化的控制参数和
sk∈SΦ,k∈[1,|SΦ|],
其中,|SF|、|SCR|和|SΦ|分别表示集合SF、SCR和SΦ的大小;
S66:迭代次数G=G+1,更新权重系数正则因子λ(G+1),计算种群P(G+1)中所有模型向量的模型目标函数Φm和反演目标函数Φ;
S67:重复步骤S62到步骤S66,直至满足终止条件。
结合第一方面的第三种可能的实现方式,在第一方面的第四种可能的实现方式中,S61中,基于初始模型,生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵并初始化数据目标函数/>模型目标函数/>正则因子λ(0),计算数据目标函数的初始均值Mean(Φd)(0),包括:
以初始模型为基础,随机扰动生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵
j∈[1,NP],i∈[1,M],
对于种群P(0)中的任意一个模型向量计算数据目标函数值/>和模型目标函数值/>同时初始化正则因子λ(0):
计算数据目标函数的初始均值Mean(Φd)(0)并记录。
结合第一方面的第四种可能的实现方式,在第一方面的第五种可能的实现方式中,S62中,针对种群P中的每一组模型向量:生成相应的交叉概率CR和尺度系数F,并计算反演目标函数值包括:
针对种群P中的每一组模型向量:
根据每个模型向量对应的反演目标函数的zscore值计算相应的交叉概率,即,对于第j个模型向量mj,在第G次迭代的交叉概率为:
其中,为交叉概率均值,zscore根据下式计算:
MAD=median|Φ-median(Φ)|,
其中,Φ为所有模型的反演目标函数值构成的向量,定义为:
对给定的模型向量采用以下公式计算反演目标函数值/>
结合第一方面的第四种可能的实现方式,在第一方面的第六种可能的实现方式中,S63中,使用改进型变异策略生成变异向量v(G),并通过交叉生成新试验模型u(G),包括:
所使用的变异策略为:
式中,r1,r2互不相同且不等于j,且r2从种群P和档案A的并集中选择,尺度系数F由位置系数为的柯西分布生成,即:
S为平滑矩阵,是一个稀疏矩阵,对于给定单元t,设相邻单元的编号为t-nz-1,t-nz,t-nz+1,t-1,t+1,t+nz-1,t+nz,t+nz+1,则权重St,t=0.2,St,k(t≠k)=0.1;
由以下公式计算:
其中,pbest表示种群P(G)中反演目标函数值排名前pb×100%的模型,而pb满足:
Φ定义为种群P(G)中所有模型向量对应的反演目标函数值形成的向量,即:
结合第一方面的第四种可能的实现方式,在第一方面的第七种可能的实现方式中,S66中,更新权重系数的方式为:
若档案A非空,则第G+1次迭代的模型参数加权矩阵为:
其中,|A|表示档案A的大小。
结合第一方面的第四种可能的实现方式,在第一方面的第八种可能的实现方式中,S66中,更新正则因子λ(G+1)的方式为:
λ(G+1)=min(λ(G),λmax),
第二方面,本申请实施例提供一种基于改进型差分进化算法的重磁位场成像反演系统,包括:
数据获取单元,用于获取地下介质的重磁位场数据,并对重磁位场数据去噪、分离区域场,得到剩余场数据;
模型生成单元,用于生成地下介质的初始模型和物性参数范围约束;
数据读入单元,用于读入剩余场数据和初始模型;
参数生成单元,用于基于剩余场数据和初始模型,生成反演网格参数和先验信息矩阵;
函数构建单元,用于生成平滑矩阵S,并建立成像反演的反演目标函数,其中,反演目标函数包含用于衡量观测数据和反演预测数据的差异的数据目标函数和用于约束反演的模型目标函数;
模型反演单元,用于利用改进型差分进化算法优化反演目标函数,最终输出数据目标函数最小的模型向量作为最终的反演结果,反演结果揭示地下介质的密度或磁化率分布。
有益效果:
1.本发明提供了一种基于改进型差分进化算法的重磁位场成像反演方法及系统,以人工智能领域的差分进化方法为基础,建立一套适合差分进化的重磁位场的成像反演机制,能有效改进差分进化的反演效果,提升该方法的收敛速度,获得有效的地下物理介质的分布情况。基于改进型的差分进化,可以有效避免全局方法获得的成像结果包含大量噪声、搜索过程耗时较长的问题,有助于弥补全局方法成像反演的不足,增强全局反演的实用性,弥补局部优化反演的不足。
2.本方案克服了传统差分进化方法在高维复杂问题进行优化时收敛速度慢、难于产生有效反演结果的缺陷,提升了全局优化方法的实用性。此外,本方案引入了自适应正则化和模型参数动态重建技术,使得本方案的成像反演在保留了全局优化的优势的同时具有局部优化的部分特性:例如,不依赖初始模型、收敛速度快、能凸显和保持成像结果的轮廓信息。
为使本申请的上述目的、特征和优点能更明显易懂,下文特举较佳实施例,并配合所附附图,作详细说明如下。
附图说明
为了更清楚地说明本申请实施例的技术方案,下面将对本申请实施例中所需要使用的附图作简单地介绍,应当理解,以下附图仅示出了本申请的某些实施例,因此不应被看作是对范围的限定,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他相关的附图。
图1为本申请实施例提供的一种基于改进型差分进化算法的重磁位场成像反演方法的流程图。
图2为利用改进型差分进化算法进行成像反演的过程。
图3为实验设计的模型示意图。
图4为所设计的模型a和模型b的重磁数据曲线。
图5为本方案与现有JADE在重力成像反演方面的收敛对比。
图6为本方案与现有JADE在磁成像反演方面的收敛对比。
图7为模型b反演时是否结合Wf,i(i∈[1,M])的重力成像反演对比。
图8为模型b反演时是否结合Wf,i(i∈[1,M])的磁成像反演对比。
图9为本申请实施例提供一种基于改进型差分进化算法的重磁位场成像反演系统的结构框图。
图标:10-基于改进型差分进化算法的重磁位场成像反演系统;11-数据获取单元;12-模型生成单元;13-数据读入单元;14-参数生成单元;15-函数构建单元;16-模型反演单元。
具体实施方式
下面将结合本申请实施例中的附图,对本申请实施例中的技术方案进行描述。
为了便于对本方案的理解,此处先对本方案涉及的一些定义进行介绍:
定义1:模型向量,由M个待优化模型参数形成的向量,记为模型向量m,即:
m=(m1,m2,…,mM)T, (1)
定义2:种群,由NP个模型向量形成的矩阵,记为种群P,即:
P=(m1,m2,…,mNP), (2)
定义3:变异向量,由差分进化的变异策略产生的向量v。
定义4:试验向量,由差分进化的交叉策略产生的向量u。
请参阅图1,图1为本申请实施例提供的一种基于改进型差分进化算法的重磁位场成像反演方法的流程图。在本实施例中,基于改进型差分进化算法的重磁位场成像反演方法可以通过电子设备(例如计算机或者服务器)运行,可以包括步骤S1、步骤S2、步骤S3、步骤S4、步骤S5和步骤S6。
首先,电子设备运行步骤S1。
步骤S1:获取地下介质的重磁位场数据,并对重磁位场数据去噪、分离区域场,得到剩余场数据。
在本实施例中,电子设备可以获取地下介质的重磁位场数据,并对重磁位场数据进行去噪、分离区域场等处理,得到剩余场数据。
与此同时,电子设备可以运行步骤S2。
步骤S2:生成地下介质的初始模型和物性参数范围约束。
在本实施例中,生成地下介质的初始模型,可以是一个预先设定的初始模型;也可以是收集研究区地质及地球物理资料后,根据经验而设置的初始模型,此处不作限定。本方案的基于改进型差分进化算法的重磁位场成像反演方法可以不依赖初始模型,因此,以预先设定的初始模型为例进行说明,不作限定。而物性参数范围约束,可以是研究人员基于研究区地质及地球物理资料输入的经验值。
之后,电子设备可以运行步骤S3。
步骤S3:读入剩余场数据和初始模型。
此时,电子设备可以读入剩余场数据和初始模型。
读取到剩余场数据和初始模型后,电子设备可以运行步骤S4。
步骤S4:基于剩余场数据和初始模型,生成反演网格参数和先验信息矩阵。
在本实施例中,电子设备可以基于剩余场数据和初始模型,生成反演网格参数和先验信息矩阵。反演网格参数可以为设定大小的网格参数,而先验信息矩阵包括深度加权矩阵和数据权重矩阵:
深度加权矩阵Wz:
其中,i为单元编号,M为待反演的物性参数个数,即待优化模型参数的个数,zi和Vi分别表示单元i的中心埋深和体积大小,β为深度衰减因子,对重力数据,β取1,对磁数据,β取2。
数据权重矩阵Wd:
Wd=diag(Wd,1,Wd,2,…,Wd,i,…,Wd,N), (4)
其中,N为观测数据个数。
而Wd,i满足:
其中,为单元i的观测数据向量,std表示标准差,dobs为观测数据向量。
需要说明的是,步骤S1、步骤S2、步骤S3、步骤S4,与现有技术并无多大差异,因此简略介绍。
生成反演网格参数和先验信息矩阵后,电子设备可以运行步骤S5。
步骤S5:生成平滑矩阵S,并建立成像反演的反演目标函数,其中,反演目标函数包含用于衡量观测数据和反演预测数据的差异的数据目标函数和用于约束反演的模型目标函数。
在本实施例中,电子设备可以生成平滑矩阵S。平滑矩阵S是一个稀疏矩阵。
以及,电子设备可以建立成像反演的反演目标函数。反演目标函数为:
Φ(m)=Φd(m)+λΦm(m), (6)
其中,Φ(m)为反演目标函数,m为由M个待优化模型参数形成的向量,记为模型向量m,m=(m1,m2,…,mM)T,λ为正则因子。
Φd(m)为数据目标函数,满足:
其中,F(m)为重磁位场正演函数。
Φm(m)为模型目标函数,满足:
其中,Wm为模型参数加权矩阵。
建立好反演目标函数后,电子设备即可运行步骤S6。
步骤S6:利用改进型差分进化算法优化反演目标函数,最终输出数据目标函数最小的模型向量作为最终的反演结果,反演结果揭示地下介质的密度或磁化率分布。
在本实施例中,电子设备可以利用改进型差分进化算法优化反演目标函数,如图2所示(利用改进型差分进化算法进行成像反演的过程),具体包括以下步骤(S61、S62、S63、S64、S65、S66、S67):
S61:迭代次数G=0时,设置初始交叉概率均值尺度因子的初始位置参数/>基于初始模型,生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵/>并初始化数据目标函数/>模型目标函数/>正则因子λ(0),计算数据目标函数的初始均值Mean(Φd)(0),其中,由NP个模型向量形成的矩阵记为种群P,P=(m1,m2,…,mNP)。
具体的,电子设备可以以初始模型为基础,随机扰动生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵
j∈[1,NP],i∈[1,M],(13)
对于种群P(0)中的任意一个模型向量电子设备可以计算数据目标函数值和模型目标函数值/>同时初始化正则因子λ(0):
然后计算数据目标函数的初始均值Mean(Φd)(0)并记录。初始均值Mean(Φd)(0)表示种群P(0)中所有模型向量的数据目标函数的均值,记为Mean(Φd)(0),计算方式为:
由此,电子设备可以完成各项参数的初始化。
S62:针对种群P中的每一组模型向量:生成相应的交叉概率CR和尺度系数F,并计算反演目标函数值表示模型向量mj在第G次迭代时的反演目标函数值。
具体的,针对种群P中的每一组模型向量:
电子设备可以根据每个模型向量对应的反演目标函数的zscore值计算相应的交叉概率,即,对于第j个模型向量mj,在第G次迭代的交叉概率为:
其中,为交叉概率均值,而zscore根据下式计算:
MAD=median|Φ-median(Φ)|, (21)
由上式,优秀的模型向量被赋予较小的交叉概率,因此,使用该方案可以保留优秀模型向量中的有效信息。而Φ为所有模型的反演目标函数值构成的向量:
对给定的模型向量采用以下公式计算反演目标函数值/>
S63:使用改进型变异策略生成变异向量v(G),并通过交叉生成新试验模型u(G)。
具体的,在变异策略中加入平均模型的影响,所使用的变异策略为:
式中,r1,r2互不相同且不等于j,且r2从种群P和档案A的并集中选择,尺度系数F由位置系数为的柯西分布生成,即:
/>
其中,S为平滑矩阵,是一个稀疏矩阵,对于给定单元t,设相邻单元的编号为t-nz-1,t-nz,t-nz+1,t-1,t+1,t+nz-1,t+nz,t+nz+1,则权重St,t=0.2,St,k(t≠k)=0.1。
而由以下公式计算:
其中,pbest表示种群P(G)中反演目标函数值排名前pb×100%的模型,因此,本方案中差分进化的搜索方向由均值模型引导,pb由数据目标函数的收敛情况决定。而pb满足:
Φ则为种群P(G)中所有模型向量对应的反演目标函数值形成的向量,即:
而交叉的具体方式如下:
i∈[1,M],j∈[1,NP], (30)
其中,rand(0,1)为0和1之间的均匀分布随机数,jrand为均匀分布的随机整数,取值区间为[1,M]。
由此,可以完成变异和交叉操作。
S64:通过选择操作,将符合要求的模型存入新的种群P(G+1),对于模型向量mj:
j∈[1,NP], (32)
对于种群P(G)中被淘汰的模型向量,将该模型向量及相应的目标函数值存入档案A中,档案A的大小为NP×M,同时,分别存储相应的尺度系数F、交叉概率CR至集合SF和集合SCR,并存储至集合SΦ。
示例性的,还可以进一步改进选择操作的方式:
选择后,计算Mean(Φd)(G+1)和Mean(Φd)(G),若Mean(Φd)(G+1)<Mean(Φd)(G),则进行S65;否则,令正则因子λ(G)=0.1λ(G),重新计算当前种群P(G)中模型向量的反演目标函数值,重新计算试验模型u(G),再进行S65。
S65:记录数据目标函数均值Mean(Φd)(G+1),更新差分进化的控制参数和/>
sk∈SΦ,k∈[1,|SΦ|],(35)
其中,|SF|、|SCR|和|SΦ|分别表示集合SF、SCR和SΦ的大小。
S66:迭代次数G=G+1,更新权重系数正则因子λ(G+1),计算种群P(G+1)中所有模型向量的模型目标函数Φm和反演目标函数Φ。
具体的,采用以下方式更新权重系数
若档案A非空,则第G+1次迭代的模型参数加权矩阵为:
其中,|A|表示档案A的大小。更新后的模型参数加权矩阵将用于计算种群P(G+1)和档案A中模型向量的模型目标函数值。
以及,采用以下方式更新正则因子λ(G+1):
λ(G+1)=min(λ(G),λmax), (40)
正则因子截断为λmax,其作用在于避免正则因子过大,导致反演过程过度拟合模型目标部分。
以及,电子设备可以计算种群P(G+1)中所有模型向量的模型目标函数Φm和反演目标函数Φ。模型目标函数Φm和反演目标函数Φ的计算请参阅前文。
S67:重复步骤S62到步骤S66,直至满足终止条件。
满足终止条件后,电子设备可以最终输出数据目标函数最小的模型向量作为最终的反演结果,反演结果揭示地下介质的物理参数分布(如密度或磁化率分布)。
因此,基于改进型差分进化算法的重磁位场成像反演方法,以人工智能领域的差分进化方法为基础,建立一套适合差分进化的重磁位场的成像反演机制,能有效改进差分进化的反演效果,提升该方法的收敛速度,获得有效的地下物理介质的分布情况。基于改进型的差分进化,可以有效避免全局方法获得的成像结果包含大量噪声、搜索过程耗时较长的问题,有助于弥补全局方法成像反演的不足,增强全局反演的实用性,弥补局部优化反演的不足。并且,本方案克服了传统差分进化方法在高维复杂问题进行优化时收敛速度慢、难于产生有效反演结果的缺陷,提升了全局优化方法的实用性。此外,本方案引入了自适应正则化和模型参数动态重建技术,使得本方案的成像反演在保留了全局优化的优势的同时具有局部优化的部分特性:例如,不依赖初始模型、收敛速度快、能凸显和保持成像结果的轮廓信息。
为了对本方案提供的基于改进型差分进化算法的重磁位场成像反演方法的效果进行验证,发明人设计了以下实验:
如图3所示(左为模型a,右为模型b),设计模型(块体部分)的密度和磁化率值分别是1g/cm3和0.1SI,设计模型的重力场值和磁场值如图4所示(左为重力异常曲线,右为磁异常曲线)。
针对设计的模型,利用本方案的基于改进型差分进化算法的重磁位场成像反演方法进行反演,涉及以下方面测试。
(1)收敛性:本方案提出了改进型的差分进化方法,目的在于改善差分进化的收敛速度。和Zhang的JADE相比,本发明具有十分明显的优势。设置密度搜索范围为[-1,2]g/cm3,密度搜索范围为[0,0.2]SI,最大迭代次数300,种群大小NP为100,重复计算10次,本发明(MADE)和JADE的数据目标函数(拟合误差)的收敛情况如图5(左为模型a使用本方案和JADE在重力成像反演方面的收敛对比,右为模型b使用本方案和JADE在重力成像反演方面的收敛对比)和图6(左为模型a使用本方案和JADE在磁成像反演方面的收敛对比,右为模型b使用本方案和JADE在磁成像反演方面的收敛对比)所示。
由图5和图6可知,本方案的改进型差分进化优化方法相较于Zhang和Sanderson的JADE具有更快的收敛速度。
(2)模型权重有效性:本发明提出了物性权重矩阵确定方式,即模型加权矩阵中存在Wf,i(i∈[1,M])项。根据基本原理,物性权重矩阵将使反演的物性参数聚焦,缩小异常区域的分布范围。以模型b的重磁位场为例,设置密度搜索范围为[-1,2]g/cm3,密度搜索范围为[0,0.2]SI,最大迭代次数300,种群大小NP为100,重复计算10次,利用本发明的改进型差分进化结合含Wf,i(i∈[1,M])的模型权重矩阵进行反演,和无Wf,i(i∈[1,M])对比,效果如图7和图8所示。图7中左为无Wf,i(i∈[1,M])时模型b的重力成像反演,右为带Wf,i(i∈[1,M])时模型b的重力成像反演;图8中左为无Wf,i(i∈[1,M])时模型b的磁成像反演,右为带Wf,i(i∈[1,M])时模型b的磁成像反演。可见,模型b的重磁位场成像结果表明本发明的模型权重矩阵构建方式能明显提升成像结果的分辨率,有效凸显模型轮廓。
请参阅图9,基于同一发明构思,本申请实施例提供一种基于改进型差分进化算法的重磁位场成像反演系统10,包括:
数据获取单元11,用于获取地下介质的重磁位场数据,并对重磁位场数据去噪、分离区域场,得到剩余场数据。
模型生成单元12,用于生成地下介质的初始模型和物性参数范围约束。
数据读入单元13,用于读入剩余场数据和初始模型。
参数生成单元14,用于基于剩余场数据和初始模型,生成反演网格参数和先验信息矩阵。
函数构建单元15,用于生成平滑矩阵S,并建立成像反演的反演目标函数,其中,反演目标函数包含用于衡量观测数据和反演预测数据的差异的数据目标函数和用于约束反演的模型目标函数。
模型反演单元16,用于利用改进型差分进化算法优化反演目标函数,最终输出数据目标函数最小的模型向量作为最终的反演结果,反演结果揭示地下介质的密度或磁化率分布。
在本实施例中,模型反演单元16,利用改进型差分进化算法优化反演目标函数的过程如前文所述,此处不再赘述。
综上所述,本发明提供一种基于改进型差分进化算法的重磁位场成像反演方法及系统,以人工智能领域的差分进化方法为基础,建立一套适合差分进化的重磁位场的成像反演机制,能有效改进差分进化的反演效果,提升该方法的收敛速度,获得有效的地下物理介质的分布情况。基于改进型的差分进化,可以有效避免全局方法获得的成像结果包含大量噪声、搜索过程耗时较长的问题,有助于弥补全局方法成像反演的不足,增强全局反演的实用性,弥补局部优化反演的不足。本方案克服了传统差分进化方法在高维复杂问题进行优化时收敛速度慢、难于产生有效反演结果的缺陷,提升了全局优化方法的实用性。此外,本方案引入了自适应正则化和模型参数动态重建技术,使得本方案的成像反演在保留了全局优化的优势的同时具有局部优化的部分特性:例如,不依赖初始模型、收敛速度快、能凸显和保持成像结果的轮廓信息。
在本文中,诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来,而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。
以上所述仅为本申请的实施例而已,并不用于限制本申请的保护范围,对于本领域的技术人员来说,本申请可以有各种更改和变化。凡在本申请的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本申请的保护范围之内。
Claims (10)
1.一种基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,包括:
S1:获取地下介质的重磁位场数据,并对重磁位场数据去噪、分离区域场,得到剩余场数据;
S2:生成地下介质的初始模型和物性参数范围约束;
S3:读入剩余场数据和初始模型;
S4:基于剩余场数据和初始模型,生成反演网格参数和先验信息矩阵;
S5:生成平滑矩阵S,并建立成像反演的反演目标函数,其中,反演目标函数包含用于衡量观测数据和反演预测数据的差异的数据目标函数和用于约束反演的模型目标函数;
S6:利用改进型差分进化算法优化反演目标函数,最终输出数据目标函数最小的模型向量作为最终的反演结果,反演结果揭示地下介质的密度或磁化率分布,其中,改进型差分进化算法对交叉概率CR的生成方式进行了改进,以及,对变异策略进行了改进。
2.根据权利要求1所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,先验信息矩阵包括:
深度加权矩阵Wz:
其中,i为单元编号,M为待反演的物性参数个数,zi和Vi分别表示单元i的中心埋深和体积大小,β为深度衰减因子,对重力数据,β取1,对磁数据,β取2;
数据权重矩阵Wd:
Wd=diag(Wd,1,Wd,2,…,Wd,i,…,Wd,N),
其中,N为观测数据个数;
Wd,i满足:
其中,为单元i的观测数据向量,std表示标准差,dobs为观测数据向量。
3.根据权利要求2所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,反演目标函数为:
Φ(m)=Φd(m)+λΦm(m),
其中,Φ(m)为反演目标函数,m为由M个待优化模型参数形成的向量,记为模型向量m,m=(m1,m2,…,mM)T,λ为正则因子,Φd(m)为数据目标函数,满足:
其中,F(m)为重磁位场正演函数;
Φm(m)为模型目标函数,满足:
其中,Wm为模型参数加权矩阵。
4.根据权利要求3所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,S6中,利用改进型差分进化算法优化反演目标函数,包括:
S61:迭代次数G=0时,设置初始交叉概率均值尺度因子的初始位置参数基于初始模型,生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵/>并初始化数据目标函数/>模型目标函数/>正则因子λ(0),计算数据目标函数的初始均值Mean(Φd)(0),其中,由NP个模型向量形成的矩阵记为种群P,P=(m1,m2,…,mNP);
S62:针对种群P中的每一组模型向量:生成相应的交叉概率CR和尺度系数F,并计算反演目标函数值 表示模型向量mj在第G次迭代时的反演目标函数值;
S63:使用改进型变异策略生成变异向量v(G),并通过交叉生成新试验模型u(G);
S64:通过选择操作,将符合要求的模型存入新的种群P(G+1),对于模型向量mj:
对于种群P(G)中被淘汰的模型向量,将该模型向量及相应的目标函数值存入档案A中,档案A的大小为NP×M,同时,分别存储相应的尺度系数F、交叉概率CR至集合SF和集合SCR,并存储至集合SΦ;
S65:记录数据目标函数均值Mean(Φd)(G+1),更新差分进化的控制参数和
sk∈SΦ,k∈[1,|SΦ|],
其中,|SF|、|SCR|和|SΦ|分别表示集合SF、SCR和SΦ的大小;
S66:迭代次数G=G+1,更新权重系数正则因子λ(G+1),计算种群P(G+1)中所有模型向量的模型目标函数Φm和反演目标函数Φ;
S67:重复步骤S62到步骤S66,直至满足终止条件。
5.根据权利要求4所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,S61中,基于初始模型,生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵并初始化数据目标函数/>模型目标函数/>正则因子λ(0),计算数据目标函数的初始均值Mean(Φd)(0),包括:
以初始模型为基础,随机扰动生成预设数量为NP的初始化种群P(0)和初始化模型参数加权矩阵
对于种群P(0)中的任意一个模型向量计算数据目标函数值/>和模型目标函数值/>同时初始化正则因子λ(0):
计算数据目标函数的初始均值Mean(Φd)(0)并记录。
6.根据权利要求5所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,S62中,针对种群P中的每一组模型向量:生成相应的交叉概率CR和尺度系数F,并计算反演目标函数值包括:
针对种群P中的每一组模型向量:
根据每个模型向量对应的反演目标函数的zscore值计算相应的交叉概率,即,对于第j个模型向量mj,在第G次迭代的交叉概率为:
其中,为交叉概率均值,zscore根据下式计算:
MAD=median|Φ-median(Φ)|,
其中,Φ为所有模型的反演目标函数值构成的向量,定义为:
对给定的模型向量采用以下公式计算反演目标函数值/>
7.根据权利要求5所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,S63中,使用改进型变异策略生成变异向量v(G),并通过交叉生成新试验模型u(G),包括:
所使用的变异策略为:
式中,r1,r2互不相同且不等于j,且r2从种群P和档案A的并集中选择,尺度系数F由位置系数为的柯西分布生成,即:
S为平滑矩阵,是一个稀疏矩阵,对于给定单元t,设相邻单元的编号为t-nz-1,t-nz,t-nz+1,t-1,t+1,t+nz-1,t+nz,t+nz+1,则权重St,t=0.2,St,k(t≠k)=0.1;
由以下公式计算:
其中,pbest表示种群P(G)中反演目标函数值排名前pb×100%的模型,而pb满足:
Φ定义为种群P(G)中所有模型向量对应的反演目标函数值形成的向量,即:
8.根据权利要求5所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,S66中,更新权重系数的方式为:
若档案A非空,则第G+1次迭代的模型参数加权矩阵为:
其中,|A|表示档案A的大小。
9.根据权利要求5所述的基于改进型差分进化算法的重磁位场成像反演方法,其特征在于,S66中,更新正则因子λ(G+1)的方式为:
λ(G+1)=min(λ(G),λmax),
10.一种基于改进型差分进化算法的重磁位场成像反演系统,其特征在于,包括:
数据获取单元,用于获取地下介质的重磁位场数据,并对重磁位场数据去噪、分离区域场,得到剩余场数据;
模型生成单元,用于生成地下介质的初始模型和物性参数范围约束;
数据读入单元,用于读入剩余场数据和初始模型;
参数生成单元,用于基于剩余场数据和初始模型,生成反演网格参数和先验信息矩阵;
函数构建单元,用于生成平滑矩阵S,并建立成像反演的反演目标函数,其中,反演目标函数包含用于衡量观测数据和反演预测数据的差异的数据目标函数和用于约束反演的模型目标函数;
模型反演单元,用于利用改进型差分进化算法优化反演目标函数,最终输出数据目标函数最小的模型向量作为最终的反演结果,反演结果揭示地下介质的密度或磁化率分布,其中,改进型差分进化算法对交叉概率CR的生成方式进行了改进,以及,对变异策略进行了改进。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311209659.3A CN117270072B (zh) | 2023-09-19 | 2023-09-19 | 基于改进型差分进化算法的重磁位场成像反演方法及系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311209659.3A CN117270072B (zh) | 2023-09-19 | 2023-09-19 | 基于改进型差分进化算法的重磁位场成像反演方法及系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN117270072A CN117270072A (zh) | 2023-12-22 |
CN117270072B true CN117270072B (zh) | 2024-04-19 |
Family
ID=89203773
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202311209659.3A Active CN117270072B (zh) | 2023-09-19 | 2023-09-19 | 基于改进型差分进化算法的重磁位场成像反演方法及系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN117270072B (zh) |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6502037B1 (en) * | 1999-04-02 | 2002-12-31 | Conoco Inc. | Method for gravity and magnetic data inversion using vector and tensor data with seismic imaging and geopressure prediction for oil, gas and mineral exploration and production |
CN108802834A (zh) * | 2018-02-13 | 2018-11-13 | 中国科学院电子学研究所 | 一种基于联合反演的地下目标识别方法 |
CN114047554A (zh) * | 2021-11-05 | 2022-02-15 | 中国南方电网有限责任公司超高压输电公司检修试验中心 | 大地电阻率模型建模方法、装置、计算机设备和存储介质 |
CN115128700A (zh) * | 2022-07-21 | 2022-09-30 | 中国地质大学(武汉) | 一种基于gramian约束的重磁三维联合反演方法及系统 |
CN115407423A (zh) * | 2021-05-27 | 2022-11-29 | 中国自然资源航空物探遥感中心 | 重、磁测数据三维反演方法及装置 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
AU2012304174B2 (en) * | 2011-09-01 | 2014-11-06 | Seequent Solutions Canada Ii Ltd. | Methods and systems for the inversion of magnetic data from remnant and induced sources in geophysical exploration |
-
2023
- 2023-09-19 CN CN202311209659.3A patent/CN117270072B/zh active Active
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6502037B1 (en) * | 1999-04-02 | 2002-12-31 | Conoco Inc. | Method for gravity and magnetic data inversion using vector and tensor data with seismic imaging and geopressure prediction for oil, gas and mineral exploration and production |
CN108802834A (zh) * | 2018-02-13 | 2018-11-13 | 中国科学院电子学研究所 | 一种基于联合反演的地下目标识别方法 |
CN115407423A (zh) * | 2021-05-27 | 2022-11-29 | 中国自然资源航空物探遥感中心 | 重、磁测数据三维反演方法及装置 |
CN114047554A (zh) * | 2021-11-05 | 2022-02-15 | 中国南方电网有限责任公司超高压输电公司检修试验中心 | 大地电阻率模型建模方法、装置、计算机设备和存储介质 |
CN115128700A (zh) * | 2022-07-21 | 2022-09-30 | 中国地质大学(武汉) | 一种基于gramian约束的重磁三维联合反演方法及系统 |
Non-Patent Citations (1)
Title |
---|
印兴耀 ; 孔栓栓 ; 张繁昌 ; 张秀颀 ; .基于差分进化算法的叠前AVO反演.石油地球物理勘探.2013,(第04期),全文. * |
Also Published As
Publication number | Publication date |
---|---|
CN117270072A (zh) | 2023-12-22 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Colombo et al. | Physics-driven deep-learning inversion with application to transient electromagnetics | |
Chen et al. | A method of DEM construction and related error analysis | |
Jiang et al. | Electrical resistivity imaging inversion: An ISFLA trained kernel principal component wavelet neural network approach | |
Xiong et al. | Multiobjective particle swarm inversion algorithm for two-dimensional magnetic data | |
CN111143984A (zh) | 基于遗传算法优化神经网络的大地电磁二维反演方法 | |
Roy et al. | Gravity inversion for heterogeneous sedimentary basin with b-spline polynomial approximation using differential evolution algorithm | |
Gajda-Zagórska et al. | A hybrid method for inversion of 3D DC resistivity logging measurements | |
CN110927806B (zh) | 基于监督下降方法的大地电磁反演方法及系统 | |
Taufik et al. | Upwind, no more: Flexible traveltime solutions using physics-informed neural networks | |
Li et al. | Transient electromagnetic inversion based on particle swarm optimization and differential evolution algorithm | |
CN112926251B (zh) | 一种基于机器学习的滑坡位移高精度预测方法 | |
CN117270072B (zh) | 基于改进型差分进化算法的重磁位场成像反演方法及系统 | |
JIANG et al. | Ultra-high density resistivity nonlinear inversion based on principal component-regularized ELM | |
Yang et al. | Comparative analyses of covariance matrix adaptation and iterative ensemble smoother on high-dimensional inverse problems in high-resolution groundwater modeling | |
Dai et al. | Nonlinear inversion for electrical resistivity tomography based on chaotic DE-BP algorithm | |
WO2021002761A1 (en) | Improved inversions of geophysical data | |
CN115755199A (zh) | 一种实用的非结构网格三维电磁反演平滑正则化方法 | |
Xiao et al. | Distributed gauss-newton optimization with smooth local parameterization for large-scale history-matching problems | |
Xie et al. | 2D magnetotelluric inversion based on ResNet | |
Rozhenko | Comparison of radial basis functions | |
Zhang et al. | Petrophysical Regression regarding Porosity, Permeability, and Water Saturation Driven by Logging‐Based Ensemble and Transfer Learnings: A Case Study of Sandy‐Mud Reservoirs | |
CN113204905A (zh) | 一种接触式激发极化法有限单元数值模拟方法 | |
CN113970789A (zh) | 全波形反演方法、装置、存储介质及电子设备 | |
Zhdanov et al. | Inversion of gravity and gravity gradiometry data for density contrast surfaces using Cauchy-type integrals | |
Menzel | Constrained indicator data resampling—A parameter constrained irregular resampling method for scattered point data |
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 |