1.本发明涉及地球物理勘探
技术领域：
：，尤其涉及一种基于拉格朗日算子动态变化的地震反演方法及系统。
背景技术：
：：2.目前，在地震反演
技术领域：
：，通常采用由constable和degroot-hedlin在1990上提出的occam反演，其基本思想是在考虑数据拟合的前提下，引入自适应正则化因子和模型粗糙度，求取最光滑模型。其基本思想是将反演的目标函数定义为由数据项和模型项组成，并给模型项加入一个拉格朗日算子λ，进而得到目标函数。3.根据目标函数得到迭代公式之后，反演过程分为两个阶段。第一个阶段为计算一系列值下的模型更新量，同时比较数据拟合差的大小，并使数据拟合差达到期望值(例如1.0)；一旦达到了期望的拟合差，则通过改变不同的λ值来寻找最光滑的模型，这是第二阶段。因此，拉格朗日算子λ的给定对于地震反演具有重大影响，目前大多根据经验给定，没有明确和科学的给定方式。技术实现要素：4.本发明的一个目的在于提供一种基于拉格朗日算子动态变化的地震反演方法，通过在反演过程中采用动态的拉格朗日算子，有利于目标函数的快速减小，保证反演结果的速度和精度。本发明的另一个目的在于提供一种基于拉格朗日算子动态变化的地震反演系统。本发明的再一个目的在于提供一种计算机设备。本发明的还一个目的在于提供一种可读介质。5.为了达到以上目的，本发明一方面公开了一种基于拉格朗日算子动态变化的地震反演方法，包括：6.根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量，n为正整数；7.将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数；8.确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。9.优选的，所述根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量具体包括：10.根据预设先验速度模型和初始速度模型计算正演响应，并计算得到数据协方差矩阵和模型协方差矩阵；11.根据正演响应、数据协方差矩阵和模型协方差矩阵得到目标函数数据项和模型项；12.根据目标函数数据项和模型项得到参考系数，并根据所述参考系数计算得到2n 1个拉格朗日算子组成的列向量。13.优选的，所述数据协方差矩阵为n×n维的对角矩阵，所述模型协方差矩阵为m×m维的满秩矩阵，其中，m为地下介质剖分总单元个数，n为数据的总个数。14.优选的，计算得到模型协方差矩阵具体包括：15.通过有限差分算法计算得到模型协方差矩阵。16.优选的，所述有限差分算法为中心差分算法。17.优选的，所述n为大于5的正整数。18.优选的，所述将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型具体包括：19.将所述列向量代入预设的反演迭代公式；20.通过共轭梯度法求解代入列向量的反演迭代公式得到对应数量的反演速度模型。21.优选的，所述目标函数预设值为1。22.本发明还公开了一种基于拉格朗日算子动态变化的地震反演系统，包括：23.数据处理模块，用于根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量，n为正整数；24.目标函数确定模块，用于将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数；25.迭代求解模块，用于确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。26.本发明还公开了一种计算机设备，包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序，27.所述处理器执行所述程序时实现如上所述方法。28.本发明还公开了一种计算机可读介质，其上存储有计算机程序，29.该程序被处理器执行时实现如上所述方法。30.本发明通过在反演迭代过程中，利用数据项和模型项的关系式求取参考系数，根据参考系数给定每次迭代过程中的多个拉格朗日算子，带入到迭代公式中，并取目标函数最小的模型作为当次迭代的反演结果，可以兼顾数据项和模型项的权重，有利于目标函数的快速减小，保证反演结果的速度和精度，一改传统occam反演需要分为两个阶段进行的传统，在每次迭代时，均能同时有效的减小数据项和模型项，不需要先搜索减小数据项的模型，再进一步缩小模型项的模型，反演的效率大大提升。附图说明31.为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。32.图1示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体实施例的流程图；33.图2示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体实施例s100的流程图；34.图3示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体实施例s300的流程图；35.图4示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子的流程图；36.图5示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子数据项的变化曲线图；37.图6示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子模型项的变化曲线图；38.图7为示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子的目标函数与现有技术目标函数的对比图；39.图8示出本发明基于拉格朗日算子动态变化的地震反演系统一个具体实施例的结构图；40.图9示出适于用来实现本发明实施例的计算机设备的结构示意图。具体实施方式41.下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。42.1990年，由constable和degroot-hedlin提出了occam反演，其基本思想是在考虑数据拟合的前提下，引入自适应正则化因子和模型粗糙度，求取最光滑模型。其基本思想是将反演的目标函数定义为由数据项和模型项组成，并给模型项加入一个拉格朗日算子λ。令地下介质剖分总单元个数为m，数据的总个数为n，其目标函数定义如下：[0043][0044]其中，d为观测数据向量，f[m]为模型正演响应，cd为数据协方差矩阵，λ为拉格朗日算子，m为模型向量(可以是速度模型、电阻率模型、密度模型等)，m0为先验模型向量，cm为模型协方差矩阵，为模型粗糙度，求解最光滑模型，就是要求这一项尽量的小。[0045]对于第k 1次的模型响应做泰勒展开，并忽略高阶项，可得：[0046]f[mk 1]≈f[mk δm]＝f[mk] jk(mk 1-mk)ꢀꢀꢀꢀꢀꢀꢀꢀ(2)[0047]其中，k代表迭代次数，jk为第k次的雅克比矩阵，将公式(2)带入目标函数，并对模型向量求取偏导并令其等于零，可得：[0048][0049]则可以得到经典occam的迭代公式：[0050][0051]其中，xk＝d-f[mk] jk(mk-m0)。可以看到求取模型改变量时，需要解一个m×m维的线性方程组，其右端项为系数矩阵为[0052]得到迭代公式之后，反演过程分为两个阶段。第一个阶段为计算一系列λ值下的模型更新量，同时比较数据拟合差的大小，并使数据拟合差达到期望值(例如1.0)；一旦达到了期望的拟合差，则通过改变不同的λ值来寻找最光滑的模型，这是第二阶段。[0053]但是目前，没有相关理论研究拉格朗日算子λ如何给定。也有学者将λ直接翻译为惩罚因子(penaltyfactor)，是为了平衡数据项和模型项在目标函数中的权重。由公式(1)可以看出，如果λ的值给定的越大，则会导致目标函数中，模型项的权重越大，则可能会使得在迭代中，数据项的“贡献”较少，有可能会导致数据拟合差越来越大，导致反演不收敛。相反，如果λ的值给定的越小，则模型项的权重越小，导致反演迭代过程中，模型项的“贡献”较小，无法达到平滑模型的作用，反演出的模型与实际情况相差甚远。[0054]另一方面，在迭代过程中，数据项和模型项的大小在不断变化，必须动态调整λ的值，才能保证数据项和模型项在整个目标函数中的权重相当，更佳有利于反演算法的效率和精度的提升。如何在动态调整λ的值，目前还是根据经验给定，没有明确的公式。[0055]另外，传统的occam反演，首先考虑数据项的快速减小，第二阶段才能考虑模型项的减小，这虽然能够保证最终模型的平滑程度，但是把反演过程分成两个阶段，导致迭代次数往往较多，计算效率不够高。[0056]为了解决现有技术的问题，根据本发明的一个方面，本实施例公开了一种基于拉格朗日算子动态变化的地震反演方法。如图1所示，本实施例中，所述方法包括：[0057]s100：根据预设先验速度模型m0和初始速度模型m1得到2n 1个拉格朗日算子组成的列向量，n为正整数。[0058]s200：将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数。[0059]s300：确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。[0060]本发明通过在反演迭代过程中，利用数据项和模型项的关系式求取参考系数，根据参考系数给定每次迭代过程中的多个拉格朗日算子，带入到迭代公式中，并取目标函数最小的模型作为当次迭代的反演结果，可以兼顾数据项和模型项的权重，有利于目标函数的快速减小，保证反演结果的速度和精度，一改传统occam反演需要分为两个阶段进行的传统，在每次迭代时，均能同时有效的减小数据项和模型项，不需要先搜索减小数据项的模型，再进一步缩小模型项的模型，反演的效率大大提升。[0061]在优选的实施方式中，如图2所示，所述s100根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量具体可包括：[0062]s110：根据预设先验速度模型m0和初始速度模型m1计算正演响应f[m1]，并计算得到数据协方差矩阵cd和模型协方差矩阵cm。[0063]s120：根据正演响应f[m1]、数据协方差矩阵cd和模型协方差矩阵cm得到目标函数数据项wd和模型项wm。[0064]s130：根据目标函数数据项和模型项得到参考系数ε，并根据所述参考系数计算得到2n 1个拉格朗日算子组成的列向量λ。[0065]其中，优选的，所述数据协方差矩阵cd为n×n维的对角矩阵，所述模型协方差矩阵cm为m×m维的满秩矩阵，其中m为地下介质剖分总单元个数，n为数据的总个数。所述模型协方差矩阵cm可通过有限差分算法计算，更优选的，可通过中心差分算法计算。可以理解的是，在实际应用中，本领域技术人员也可采用其他可行的方法得到以上计算结果，本发明在此并不作限定。[0066]进一步的，根据正演响应f[m1]、数据协方差矩阵cd和模型协方差矩阵cm可根据目标函数数据项和模型项和模型项的表达式计算得到目标函数数据项wd和模型项wm。[0067]然后，根据数据项wd和模型项wm，通过参考系数ε的表达式可计算得到参考系数ε，其中，括号[]表示取整。根据参考系数ε计算2n 1个拉格朗日算子组成的列向量λ(2n 1)＝10ε-n,10ε-n 1,10ε-n 2…10ε n，其中，n为大于1的正整数，n优选的为大于5的正整数。[0068]在优选的实施方式中，如图3所示，所述s300将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型具体可包括：[0069]s310：将所述列向量代入预设的反演迭代公式。[0070]s320：通过共轭梯度法求解代入列向量的反演迭代公式得到对应数量的反演速度模型。[0071]具体的，可将上述得到的列向量λ(2n 1)带入迭代公式可使用共轭梯度法求解迭代公式。计算得到2n 1个拉格朗日算子λ对应的2n 1个反演速度模型，并计算对应的2n 1个目标函数。计算2n 1个拉格朗日算子λ对应的2n 1个反演速度模型，并计算对应的2n 1个目标函数。[0072]下面通过一个具体例子来对本发明作进一步的说明。该具体例子中，数据个数n取值为100，地下介质剖分总单元个数n取值为10000，n取5。数据协方差矩阵cd的对角线元素为100个数据的随机误差，随机误差为均值为0的白噪声。cm通过每个网格剖分之间的物性的有限差分算法计算，本实施例中使用中心差分计算cm。每次迭代过程中，cd保持不变，但cm随着每次迭代过程中物性的改变而需要重新计算。[0073]具体的，如图4所示，该具体例子中，基于拉格朗日算子动态变化的地震反演方法包括以下步骤：[0074]s1：给定先验速度模型m0和初始速度模型m1，计算初始速度模型的正演响应f[m1]。计算数据协方差矩阵cd和模型协方差矩阵cm。所述数据协方差矩阵cd为n×n维的对角矩阵，所述模型协方差矩阵cm为m×m维的满秩矩阵。[0075]s2：根据公式计算目标函数数据项和模型项本实施例中，首次迭代计算获得wd＝1579.64，wm＝17.87；[0076]s3：计算参考系数其中，该公式中的中括号表示取整。[0077]s4：根据参考系数ε计算2n 1＝11个拉格朗日算子组成的列向量λ(2n 1)＝10-3,10-2,10-1…107。[0078]s5：带入迭代公式使用共轭梯度法求解迭代公式。计算11个拉格朗日算子λ对应的11个反演速度模型，并计算对应的11个目标函数值。[0079]s6：取目标函数值最小的反演模型作为第一次迭代的反演速度模型。[0080]s7：重复步骤s2-s6，直到第k次迭代的n 1个反演速度模型所对应的n 1个目标函数中的最小值达到目标函数预设值，目标函数预设值为1。[0081]表1为本实施例中的迭代过程中拉格朗日算子的变化过程，可得到第1、2、3、4、5、6、7和13次迭代时对应的目标函数最小值分别为103、106、101、10-1、100、101、10-1和10-1。[0082]表1迭代过程拉格朗日算子动态搜索表[0083]k＝1，k＝2k＝3k＝4k＝5k＝6k＝7……k＝1310-310-410-510-510-510-510-5……10-510-210-310-410-410-410-410-4……10-410-110-210-310-310-310-310-3……10-310010-110-210-210-210-210-2……10-210110010-110-110-110-110-1……10-1102101100100100100100……100103102101101101101101……101104103102102102102102……102105104103103103103103……103106105104104104104104……104107106105105105105105……105[0084]由表中我们可以看到，如果在每次迭代过程中，都选择固定的拉格朗日算子，则不能搜索到目标函数下降最快的拉格朗日算子。图5和图6为本实施例中，数据项和模型项随着迭代次数的增加的曲线情况。由于数量级变化较大，因此对纵轴取了对数。由图中可以看出，数据项一直在单调递减，模型项经历了一个先增后减的过程，这是因为最开始给定的模型较为简单，模型项初始值较小，在迭代搜索过程中，模型趋于复杂。当数据项减小到一定程度后，搜索的方向又使得模型趋于平滑，因而又开始减小。图7为本发明的目标函数与现有技术目标函数的对比图，现有技术中拉格朗日系数取103保持不变，由此可见，从第二次迭代开始，现有技术的目标函数的减小趋势要慢于本发明，而且本发明在迭代13次之后就收敛，达到预定的值，而现有技术经历了17次迭代，由此可见本发明在缩短计算时间方面的效果良好。[0085]本发明通过在反演迭代过程中，求取参考系数，根据参考系数给定每次迭代过程中的多个拉格朗日算子，带入到迭代公式中，可以取得如下技术效果：[0086](1)本发明可以兼顾数据项和模型项的权重，使得数据项和模型项在目标函数中的贡献大小朝着最有利于目标函数的快速减小的方向发展。[0087](2)本发明可以在每次迭代过程中，并非生搬硬套的使用一系列的拉格朗日算子，而是根据每次数据项和模型项的相对大小，给出一系列随着迭代次数动态变化的系数，这样子在每次迭代过程中，都可以使得目标函数快速下降，保证反演结果的速度和精度。[0088](3)本发明一改传统occam反演需要分为两个阶段进行的传统，在每次迭代时，均能同时有效的减小数据项和模型项，不需要先搜索减小数据项的模型，再进一步缩小模型项的模型，反演的效率大大提升。[0089]基于相同原理，本实施例还公开了一种基于拉格朗日算子动态变化的地震反演系统。如图8所示，所述系统包括数据处理模块11、目标函数确定模块12和迭代求解模块13。[0090]其中，数据处理模块11用于根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量，n为正整数。[0091]目标函数确定模块12用于将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数。[0092]迭代求解模块13用于确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。[0093]由于该系统解决问题的原理与以上方法类似，因此本系统的实施可以参见方法的实施，在此不再赘述。[0094]上述实施例阐明的系统、装置、模块或单元，具体可以由计算机芯片或实体实现，或者由具有某种功能的产品来实现。一种典型的实现设备为计算机设备，具体的，计算机设备例如可以为个人计算机、膝上型计算机、蜂窝电话、相机电话、智能电话、个人数字助理、媒体播放器、导航设备、电子邮件设备、游戏控制台、平板计算机、可穿戴设备或者这些设备中的任何设备的组合。[0095]在一个典型的实例中计算机设备具体包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现如上所述方法。[0096]下面参考图9，其示出了适于用来实现本技术实施例的计算机设备600的结构示意图。[0097]如图9所示，计算机设备600包括中央处理单元(cpu)601，其可以根据存储在只读存储器(rom)602中的程序或者从存储部分608加载到随机访问存储器(ram))603中的程序而执行各种适当的工作和处理。在ram603中，还存储有系统600操作所需的各种程序和数据。cpu601、rom602、以及ram603通过总线604彼此相连。输入/输出(i/o)接口605也连接至总线604。[0098]以下部件连接至i/o接口605：包括键盘、鼠标等的输入部分606；包括诸如阴极射线管(crt)、液晶反馈器(lcd)等以及扬声器等的输出部分607；包括硬盘等的存储部分608；以及包括诸如lan卡，调制解调器等的网络接口卡的通信部分609。通信部分609经由诸如因特网的网络执行通信处理。驱动器610也根据需要连接至i/o接口605。可拆卸介质611，诸如磁盘、光盘、磁光盘、半导体存储器等等，根据需要安装在驱动器610上，以便于从其上读出的计算机程序根据需要被安装如存储部分608。[0099]特别地，根据本发明的实施例，上文参考流程图描述的过程可以被实现为计算机软件程序。例如，本发明的实施例包括一种计算机程序产品，其包括有形地包含在机器可读介质上的计算机程序，所述计算机程序包括用于执行流程图所示的方法的程序代码。在这样的实施例中，该计算机程序可以通过通信部分609从网络上被下载和安装，和/或从可拆卸介质611被安装。[0100]计算机可读介质包括永久性和非永久性、可移动和非可移动媒体可以由任何方法或技术来实现信息存储。信息可以是计算机可读指令、数据结构、程序的模块或其他数据。计算机的存储介质的例子包括，但不限于相变内存(pram)、静态随机存取存储器(sram)、动态随机存取存储器(dram)、其他类型的随机存取存储器(ram)、只读存储器(rom)、电可擦除可编程只读存储器(eeprom)、快闪记忆体或其他内存技术、只读光盘只读存储器(cd-rom)、数字多功能光盘(dvd)或其他光学存储、磁盒式磁带，磁带磁磁盘存储或其他磁性存储设备或任何其他非传输介质，可用于存储可以被计算设备访问的信息。按照本文中的界定，计算机可读介质不包括暂存电脑可读媒体(transitorymedia)，如调制的数据信号和载波。[0101]为了描述的方便，描述以上装置时以功能分为各种单元分别描述。当然，在实施本技术时可以把各单元的功能在同一个或2n 1个软件和/或硬件中实现。[0102]本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或2n 1个流程和/或方框图一个方框或2n 1个方框中指定的功能的装置。[0103]这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或2n 1个流程和/或方框图一个方框或2n 1个方框中指定的功能。[0104]这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或2n 1个流程和/或方框图一个方框或2n 1个方框中指定的功能的步骤。[0105]还需要说明的是，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、商品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、商品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、商品或者设备中还存在另外的相同要素。[0106]本领域技术人员应明白，本技术的实施例可提供为方法、系统或计算机程序产品。因此，本技术可采用完全硬件实施例、完全软件实施例或结合软件和硬件方面的实施例的形式。而且，本技术可采用在一个或2n 1个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd-rom、光学存储器等)上实施的计算机程序产品的形式。[0107]本技术可以在由计算机执行的计算机可执行指令的一般上下文中描述，例如程序模块。一般地，程序模块包括执行特定任务或实现特定抽象数据类型的例程、程序、对象、组件、数据结构等等。也可以在分布式计算环境中实践本技术，在这些分布式计算环境中，由通过通信网络而被连接的远程处理设备来执行任务。在分布式计算环境中，程序模块可以位于包括存储设备在内的本地和远程计算机存储介质中。[0108]本说明书中的各个实施例均采用递进的方式描述，各个实施例之间相同相似的部分互相参见即可，每个实施例重点说明的都是与其他实施例的不同之处。尤其，对于系统实施例而言，由于其基本相似于方法实施例，所以描述的比较简单，相关之处参见方法实施例的部分说明即可。[0109]以上所述仅为本技术的实施例而已，并不用于限制本技术。对于本领域技术人员来说，本技术可以有各种更改和变化。凡在本技术的精神和原理之内所作的任何修改、等同替换、改进等，均应包含在本技术的权利要求范围之内。当前第1页12当前第1页12
技术领域：
：，尤其涉及一种基于拉格朗日算子动态变化的地震反演方法及系统。
背景技术：
：：2.目前，在地震反演
技术领域：
：，通常采用由constable和degroot-hedlin在1990上提出的occam反演，其基本思想是在考虑数据拟合的前提下，引入自适应正则化因子和模型粗糙度，求取最光滑模型。其基本思想是将反演的目标函数定义为由数据项和模型项组成，并给模型项加入一个拉格朗日算子λ，进而得到目标函数。3.根据目标函数得到迭代公式之后，反演过程分为两个阶段。第一个阶段为计算一系列值下的模型更新量，同时比较数据拟合差的大小，并使数据拟合差达到期望值(例如1.0)；一旦达到了期望的拟合差，则通过改变不同的λ值来寻找最光滑的模型，这是第二阶段。因此，拉格朗日算子λ的给定对于地震反演具有重大影响，目前大多根据经验给定，没有明确和科学的给定方式。技术实现要素：4.本发明的一个目的在于提供一种基于拉格朗日算子动态变化的地震反演方法，通过在反演过程中采用动态的拉格朗日算子，有利于目标函数的快速减小，保证反演结果的速度和精度。本发明的另一个目的在于提供一种基于拉格朗日算子动态变化的地震反演系统。本发明的再一个目的在于提供一种计算机设备。本发明的还一个目的在于提供一种可读介质。5.为了达到以上目的，本发明一方面公开了一种基于拉格朗日算子动态变化的地震反演方法，包括：6.根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量，n为正整数；7.将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数；8.确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。9.优选的，所述根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量具体包括：10.根据预设先验速度模型和初始速度模型计算正演响应，并计算得到数据协方差矩阵和模型协方差矩阵；11.根据正演响应、数据协方差矩阵和模型协方差矩阵得到目标函数数据项和模型项；12.根据目标函数数据项和模型项得到参考系数，并根据所述参考系数计算得到2n 1个拉格朗日算子组成的列向量。13.优选的，所述数据协方差矩阵为n×n维的对角矩阵，所述模型协方差矩阵为m×m维的满秩矩阵，其中，m为地下介质剖分总单元个数，n为数据的总个数。14.优选的，计算得到模型协方差矩阵具体包括：15.通过有限差分算法计算得到模型协方差矩阵。16.优选的，所述有限差分算法为中心差分算法。17.优选的，所述n为大于5的正整数。18.优选的，所述将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型具体包括：19.将所述列向量代入预设的反演迭代公式；20.通过共轭梯度法求解代入列向量的反演迭代公式得到对应数量的反演速度模型。21.优选的，所述目标函数预设值为1。22.本发明还公开了一种基于拉格朗日算子动态变化的地震反演系统，包括：23.数据处理模块，用于根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量，n为正整数；24.目标函数确定模块，用于将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数；25.迭代求解模块，用于确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。26.本发明还公开了一种计算机设备，包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序，27.所述处理器执行所述程序时实现如上所述方法。28.本发明还公开了一种计算机可读介质，其上存储有计算机程序，29.该程序被处理器执行时实现如上所述方法。30.本发明通过在反演迭代过程中，利用数据项和模型项的关系式求取参考系数，根据参考系数给定每次迭代过程中的多个拉格朗日算子，带入到迭代公式中，并取目标函数最小的模型作为当次迭代的反演结果，可以兼顾数据项和模型项的权重，有利于目标函数的快速减小，保证反演结果的速度和精度，一改传统occam反演需要分为两个阶段进行的传统，在每次迭代时，均能同时有效的减小数据项和模型项，不需要先搜索减小数据项的模型，再进一步缩小模型项的模型，反演的效率大大提升。附图说明31.为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。32.图1示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体实施例的流程图；33.图2示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体实施例s100的流程图；34.图3示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体实施例s300的流程图；35.图4示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子的流程图；36.图5示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子数据项的变化曲线图；37.图6示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子模型项的变化曲线图；38.图7为示出本发明基于拉格朗日算子动态变化的地震反演方法一个具体例子的目标函数与现有技术目标函数的对比图；39.图8示出本发明基于拉格朗日算子动态变化的地震反演系统一个具体实施例的结构图；40.图9示出适于用来实现本发明实施例的计算机设备的结构示意图。具体实施方式41.下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。42.1990年，由constable和degroot-hedlin提出了occam反演，其基本思想是在考虑数据拟合的前提下，引入自适应正则化因子和模型粗糙度，求取最光滑模型。其基本思想是将反演的目标函数定义为由数据项和模型项组成，并给模型项加入一个拉格朗日算子λ。令地下介质剖分总单元个数为m，数据的总个数为n，其目标函数定义如下：[0043][0044]其中，d为观测数据向量，f[m]为模型正演响应，cd为数据协方差矩阵，λ为拉格朗日算子，m为模型向量(可以是速度模型、电阻率模型、密度模型等)，m0为先验模型向量，cm为模型协方差矩阵，为模型粗糙度，求解最光滑模型，就是要求这一项尽量的小。[0045]对于第k 1次的模型响应做泰勒展开，并忽略高阶项，可得：[0046]f[mk 1]≈f[mk δm]＝f[mk] jk(mk 1-mk)ꢀꢀꢀꢀꢀꢀꢀꢀ(2)[0047]其中，k代表迭代次数，jk为第k次的雅克比矩阵，将公式(2)带入目标函数，并对模型向量求取偏导并令其等于零，可得：[0048][0049]则可以得到经典occam的迭代公式：[0050][0051]其中，xk＝d-f[mk] jk(mk-m0)。可以看到求取模型改变量时，需要解一个m×m维的线性方程组，其右端项为系数矩阵为[0052]得到迭代公式之后，反演过程分为两个阶段。第一个阶段为计算一系列λ值下的模型更新量，同时比较数据拟合差的大小，并使数据拟合差达到期望值(例如1.0)；一旦达到了期望的拟合差，则通过改变不同的λ值来寻找最光滑的模型，这是第二阶段。[0053]但是目前，没有相关理论研究拉格朗日算子λ如何给定。也有学者将λ直接翻译为惩罚因子(penaltyfactor)，是为了平衡数据项和模型项在目标函数中的权重。由公式(1)可以看出，如果λ的值给定的越大，则会导致目标函数中，模型项的权重越大，则可能会使得在迭代中，数据项的“贡献”较少，有可能会导致数据拟合差越来越大，导致反演不收敛。相反，如果λ的值给定的越小，则模型项的权重越小，导致反演迭代过程中，模型项的“贡献”较小，无法达到平滑模型的作用，反演出的模型与实际情况相差甚远。[0054]另一方面，在迭代过程中，数据项和模型项的大小在不断变化，必须动态调整λ的值，才能保证数据项和模型项在整个目标函数中的权重相当，更佳有利于反演算法的效率和精度的提升。如何在动态调整λ的值，目前还是根据经验给定，没有明确的公式。[0055]另外，传统的occam反演，首先考虑数据项的快速减小，第二阶段才能考虑模型项的减小，这虽然能够保证最终模型的平滑程度，但是把反演过程分成两个阶段，导致迭代次数往往较多，计算效率不够高。[0056]为了解决现有技术的问题，根据本发明的一个方面，本实施例公开了一种基于拉格朗日算子动态变化的地震反演方法。如图1所示，本实施例中，所述方法包括：[0057]s100：根据预设先验速度模型m0和初始速度模型m1得到2n 1个拉格朗日算子组成的列向量，n为正整数。[0058]s200：将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数。[0059]s300：确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。[0060]本发明通过在反演迭代过程中，利用数据项和模型项的关系式求取参考系数，根据参考系数给定每次迭代过程中的多个拉格朗日算子，带入到迭代公式中，并取目标函数最小的模型作为当次迭代的反演结果，可以兼顾数据项和模型项的权重，有利于目标函数的快速减小，保证反演结果的速度和精度，一改传统occam反演需要分为两个阶段进行的传统，在每次迭代时，均能同时有效的减小数据项和模型项，不需要先搜索减小数据项的模型，再进一步缩小模型项的模型，反演的效率大大提升。[0061]在优选的实施方式中，如图2所示，所述s100根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量具体可包括：[0062]s110：根据预设先验速度模型m0和初始速度模型m1计算正演响应f[m1]，并计算得到数据协方差矩阵cd和模型协方差矩阵cm。[0063]s120：根据正演响应f[m1]、数据协方差矩阵cd和模型协方差矩阵cm得到目标函数数据项wd和模型项wm。[0064]s130：根据目标函数数据项和模型项得到参考系数ε，并根据所述参考系数计算得到2n 1个拉格朗日算子组成的列向量λ。[0065]其中，优选的，所述数据协方差矩阵cd为n×n维的对角矩阵，所述模型协方差矩阵cm为m×m维的满秩矩阵，其中m为地下介质剖分总单元个数，n为数据的总个数。所述模型协方差矩阵cm可通过有限差分算法计算，更优选的，可通过中心差分算法计算。可以理解的是，在实际应用中，本领域技术人员也可采用其他可行的方法得到以上计算结果，本发明在此并不作限定。[0066]进一步的，根据正演响应f[m1]、数据协方差矩阵cd和模型协方差矩阵cm可根据目标函数数据项和模型项和模型项的表达式计算得到目标函数数据项wd和模型项wm。[0067]然后，根据数据项wd和模型项wm，通过参考系数ε的表达式可计算得到参考系数ε，其中，括号[]表示取整。根据参考系数ε计算2n 1个拉格朗日算子组成的列向量λ(2n 1)＝10ε-n,10ε-n 1,10ε-n 2…10ε n，其中，n为大于1的正整数，n优选的为大于5的正整数。[0068]在优选的实施方式中，如图3所示，所述s300将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型具体可包括：[0069]s310：将所述列向量代入预设的反演迭代公式。[0070]s320：通过共轭梯度法求解代入列向量的反演迭代公式得到对应数量的反演速度模型。[0071]具体的，可将上述得到的列向量λ(2n 1)带入迭代公式可使用共轭梯度法求解迭代公式。计算得到2n 1个拉格朗日算子λ对应的2n 1个反演速度模型，并计算对应的2n 1个目标函数。计算2n 1个拉格朗日算子λ对应的2n 1个反演速度模型，并计算对应的2n 1个目标函数。[0072]下面通过一个具体例子来对本发明作进一步的说明。该具体例子中，数据个数n取值为100，地下介质剖分总单元个数n取值为10000，n取5。数据协方差矩阵cd的对角线元素为100个数据的随机误差，随机误差为均值为0的白噪声。cm通过每个网格剖分之间的物性的有限差分算法计算，本实施例中使用中心差分计算cm。每次迭代过程中，cd保持不变，但cm随着每次迭代过程中物性的改变而需要重新计算。[0073]具体的，如图4所示，该具体例子中，基于拉格朗日算子动态变化的地震反演方法包括以下步骤：[0074]s1：给定先验速度模型m0和初始速度模型m1，计算初始速度模型的正演响应f[m1]。计算数据协方差矩阵cd和模型协方差矩阵cm。所述数据协方差矩阵cd为n×n维的对角矩阵，所述模型协方差矩阵cm为m×m维的满秩矩阵。[0075]s2：根据公式计算目标函数数据项和模型项本实施例中，首次迭代计算获得wd＝1579.64，wm＝17.87；[0076]s3：计算参考系数其中，该公式中的中括号表示取整。[0077]s4：根据参考系数ε计算2n 1＝11个拉格朗日算子组成的列向量λ(2n 1)＝10-3,10-2,10-1…107。[0078]s5：带入迭代公式使用共轭梯度法求解迭代公式。计算11个拉格朗日算子λ对应的11个反演速度模型，并计算对应的11个目标函数值。[0079]s6：取目标函数值最小的反演模型作为第一次迭代的反演速度模型。[0080]s7：重复步骤s2-s6，直到第k次迭代的n 1个反演速度模型所对应的n 1个目标函数中的最小值达到目标函数预设值，目标函数预设值为1。[0081]表1为本实施例中的迭代过程中拉格朗日算子的变化过程，可得到第1、2、3、4、5、6、7和13次迭代时对应的目标函数最小值分别为103、106、101、10-1、100、101、10-1和10-1。[0082]表1迭代过程拉格朗日算子动态搜索表[0083]k＝1，k＝2k＝3k＝4k＝5k＝6k＝7……k＝1310-310-410-510-510-510-510-5……10-510-210-310-410-410-410-410-4……10-410-110-210-310-310-310-310-3……10-310010-110-210-210-210-210-2……10-210110010-110-110-110-110-1……10-1102101100100100100100……100103102101101101101101……101104103102102102102102……102105104103103103103103……103106105104104104104104……104107106105105105105105……105[0084]由表中我们可以看到，如果在每次迭代过程中，都选择固定的拉格朗日算子，则不能搜索到目标函数下降最快的拉格朗日算子。图5和图6为本实施例中，数据项和模型项随着迭代次数的增加的曲线情况。由于数量级变化较大，因此对纵轴取了对数。由图中可以看出，数据项一直在单调递减，模型项经历了一个先增后减的过程，这是因为最开始给定的模型较为简单，模型项初始值较小，在迭代搜索过程中，模型趋于复杂。当数据项减小到一定程度后，搜索的方向又使得模型趋于平滑，因而又开始减小。图7为本发明的目标函数与现有技术目标函数的对比图，现有技术中拉格朗日系数取103保持不变，由此可见，从第二次迭代开始，现有技术的目标函数的减小趋势要慢于本发明，而且本发明在迭代13次之后就收敛，达到预定的值，而现有技术经历了17次迭代，由此可见本发明在缩短计算时间方面的效果良好。[0085]本发明通过在反演迭代过程中，求取参考系数，根据参考系数给定每次迭代过程中的多个拉格朗日算子，带入到迭代公式中，可以取得如下技术效果：[0086](1)本发明可以兼顾数据项和模型项的权重，使得数据项和模型项在目标函数中的贡献大小朝着最有利于目标函数的快速减小的方向发展。[0087](2)本发明可以在每次迭代过程中，并非生搬硬套的使用一系列的拉格朗日算子，而是根据每次数据项和模型项的相对大小，给出一系列随着迭代次数动态变化的系数，这样子在每次迭代过程中，都可以使得目标函数快速下降，保证反演结果的速度和精度。[0088](3)本发明一改传统occam反演需要分为两个阶段进行的传统，在每次迭代时，均能同时有效的减小数据项和模型项，不需要先搜索减小数据项的模型，再进一步缩小模型项的模型，反演的效率大大提升。[0089]基于相同原理，本实施例还公开了一种基于拉格朗日算子动态变化的地震反演系统。如图8所示，所述系统包括数据处理模块11、目标函数确定模块12和迭代求解模块13。[0090]其中，数据处理模块11用于根据预设先验速度模型和初始速度模型得到2n 1个拉格朗日算子组成的列向量，n为正整数。[0091]目标函数确定模块12用于将所述列向量代入预设的反演迭代公式得到对应数量的反演速度模型，并根据所述反演速度模型得到目标函数。[0092]迭代求解模块13用于确定所述目标函数的函数值最小的反演速度模型为迭代模型，重复确定所述迭代模型直至目标函数的最小值达到目标函数预设值。[0093]由于该系统解决问题的原理与以上方法类似，因此本系统的实施可以参见方法的实施，在此不再赘述。[0094]上述实施例阐明的系统、装置、模块或单元，具体可以由计算机芯片或实体实现，或者由具有某种功能的产品来实现。一种典型的实现设备为计算机设备，具体的，计算机设备例如可以为个人计算机、膝上型计算机、蜂窝电话、相机电话、智能电话、个人数字助理、媒体播放器、导航设备、电子邮件设备、游戏控制台、平板计算机、可穿戴设备或者这些设备中的任何设备的组合。[0095]在一个典型的实例中计算机设备具体包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现如上所述方法。[0096]下面参考图9，其示出了适于用来实现本技术实施例的计算机设备600的结构示意图。[0097]如图9所示，计算机设备600包括中央处理单元(cpu)601，其可以根据存储在只读存储器(rom)602中的程序或者从存储部分608加载到随机访问存储器(ram))603中的程序而执行各种适当的工作和处理。在ram603中，还存储有系统600操作所需的各种程序和数据。cpu601、rom602、以及ram603通过总线604彼此相连。输入/输出(i/o)接口605也连接至总线604。[0098]以下部件连接至i/o接口605：包括键盘、鼠标等的输入部分606；包括诸如阴极射线管(crt)、液晶反馈器(lcd)等以及扬声器等的输出部分607；包括硬盘等的存储部分608；以及包括诸如lan卡，调制解调器等的网络接口卡的通信部分609。通信部分609经由诸如因特网的网络执行通信处理。驱动器610也根据需要连接至i/o接口605。可拆卸介质611，诸如磁盘、光盘、磁光盘、半导体存储器等等，根据需要安装在驱动器610上，以便于从其上读出的计算机程序根据需要被安装如存储部分608。[0099]特别地，根据本发明的实施例，上文参考流程图描述的过程可以被实现为计算机软件程序。例如，本发明的实施例包括一种计算机程序产品，其包括有形地包含在机器可读介质上的计算机程序，所述计算机程序包括用于执行流程图所示的方法的程序代码。在这样的实施例中，该计算机程序可以通过通信部分609从网络上被下载和安装，和/或从可拆卸介质611被安装。[0100]计算机可读介质包括永久性和非永久性、可移动和非可移动媒体可以由任何方法或技术来实现信息存储。信息可以是计算机可读指令、数据结构、程序的模块或其他数据。计算机的存储介质的例子包括，但不限于相变内存(pram)、静态随机存取存储器(sram)、动态随机存取存储器(dram)、其他类型的随机存取存储器(ram)、只读存储器(rom)、电可擦除可编程只读存储器(eeprom)、快闪记忆体或其他内存技术、只读光盘只读存储器(cd-rom)、数字多功能光盘(dvd)或其他光学存储、磁盒式磁带，磁带磁磁盘存储或其他磁性存储设备或任何其他非传输介质，可用于存储可以被计算设备访问的信息。按照本文中的界定，计算机可读介质不包括暂存电脑可读媒体(transitorymedia)，如调制的数据信号和载波。[0101]为了描述的方便，描述以上装置时以功能分为各种单元分别描述。当然，在实施本技术时可以把各单元的功能在同一个或2n 1个软件和/或硬件中实现。[0102]本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或2n 1个流程和/或方框图一个方框或2n 1个方框中指定的功能的装置。[0103]这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或2n 1个流程和/或方框图一个方框或2n 1个方框中指定的功能。[0104]这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或2n 1个流程和/或方框图一个方框或2n 1个方框中指定的功能的步骤。[0105]还需要说明的是，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、商品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、商品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、商品或者设备中还存在另外的相同要素。[0106]本领域技术人员应明白，本技术的实施例可提供为方法、系统或计算机程序产品。因此，本技术可采用完全硬件实施例、完全软件实施例或结合软件和硬件方面的实施例的形式。而且，本技术可采用在一个或2n 1个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd-rom、光学存储器等)上实施的计算机程序产品的形式。[0107]本技术可以在由计算机执行的计算机可执行指令的一般上下文中描述，例如程序模块。一般地，程序模块包括执行特定任务或实现特定抽象数据类型的例程、程序、对象、组件、数据结构等等。也可以在分布式计算环境中实践本技术，在这些分布式计算环境中，由通过通信网络而被连接的远程处理设备来执行任务。在分布式计算环境中，程序模块可以位于包括存储设备在内的本地和远程计算机存储介质中。[0108]本说明书中的各个实施例均采用递进的方式描述，各个实施例之间相同相似的部分互相参见即可，每个实施例重点说明的都是与其他实施例的不同之处。尤其，对于系统实施例而言，由于其基本相似于方法实施例，所以描述的比较简单，相关之处参见方法实施例的部分说明即可。[0109]以上所述仅为本技术的实施例而已，并不用于限制本技术。对于本领域技术人员来说，本技术可以有各种更改和变化。凡在本技术的精神和原理之内所作的任何修改、等同替换、改进等，均应包含在本技术的权利要求范围之内。当前第1页12当前第1页12
再多了解一些
本文用于企业家、创业者技术爱好者查询，结果仅供参考。
