• 四川郎酒股份有限公司获第十二届人民企业社会责任奖年度环保奖 2019-05-13
  • 银保监会新规剑指大企业多头融资和过度融资 2019-05-12
  • 韩国再提4国联合申办世界杯 中国网友无视:我们自己来 2019-05-11
  • 中国人为什么一定要买房? 2019-05-11
  • 十九大精神进校园:风正扬帆当有为 勇做时代弄潮儿 2019-05-10
  • 粽叶飘香幸福邻里——廊坊市举办“我们的节日·端午”主题活动 2019-05-09
  • 太原设禁鸣路段 设备在测试中 2019-05-09
  • 拜耳医药保健有限公司获第十二届人民企业社会责任奖年度企业奖 2019-05-08
  • “港独”没出路!“梁天琦们”该醒醒了 2019-05-07
  • 陈卫平:中国文化内涵包含三方面 文化复兴表现在其中 2019-05-06
  • 人民日报客户端辟谣:“合成军装照”产品请放心使用 2019-05-05
  • 【十九大·理论新视野】为什么要“建设现代化经济体系”?   2019-05-04
  • 聚焦2017年乌鲁木齐市老城区改造提升工程 2019-05-04
  • 【专家谈】上合组织——构建区域命运共同体的有力实践者 2019-05-03
  • 【华商侃车NO.192】 亲!楼市火爆,别忘了买车位啊! 2019-05-03
    • / 23
    • 下载费用:30 金币  

    重庆时时彩后一杀号法: 基于非线性优化的时空域交错网格有限差分方法和装置.pdf

    关 键 词:
    基于 非线性 优化 时空 交错 网格 有限 方法 装置
      专利查询网所有资源均是用户自行上传分享,仅供网友学习交流,未经上传用户书面授权,请勿作他用。
    摘要
    申请专利号:

    CN201310660960.6

    申请日:

    2013.12.09

    公开号:

    CN103630933A

    公开日:

    2014.03.12

    当前法律状态:

    授权

    有效性:

    有权

    法律详情: 授权|||实质审查的生效IPC(主分类):G01V 1/28申请日:20131209|||公开
    IPC分类号: G01V1/28 主分类号: G01V1/28
    申请人: 中国石油天然气集团公司; 中国石油大学(北京)
    发明人: 任志明; 刘洋
    地址: 100007 北京市东城区东直门北大街9号
    优先权:
    专利代理机构: 北京三友知识产权代理有限公司 11127 代理人: 王天尧
    PDF完整版下载: PDF下载
    法律状态
    申请(专利)号:

    CN201310660960.6

    授权公告号:

    ||||||

    法律状态公告日:

    2017.01.18|||2014.04.09|||2014.03.12

    法律状态类型:

    授权|||实质审查的生效|||公开

    摘要

    本发明提供了一种基于非线性优化的时空域交错网格有限差分方法和装置,其中,该方法包括:确定有限差分系数;基于时空域频散关系和非线性反演算法对有限差分系数进行优化;利用优化后的有限差分系数进行弹性波正演模拟。本发明解决了现有技术中采用泰勒级数展开和空间域频散关系的有限差分法获得有限差分系数进行弹性波正演模拟而导致的中高频段频散较大,模拟精度较低的技术问题,达到了减小中高频段的频散,提高模拟精度的技术效果。

    权利要求书

    权利要求书
    1.  一种基于非线性优化的时空域交错网格有限差分方法,其特征在于,包括:
    确定有限差分系数;
    基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    利用优化后的有限差分系数进行弹性波正演模拟。

    2.  如权利要求1所述的方法,其特征在于,所述确定有限差分系数,包括:
    按照以下公式确定有限差分系数:
    am=(-1)m+12m-1Π1nM,n≠m|(2n-1)2(2n-1)2-(2m-1)2|]]>
    其中,am为有限差分系数,M为算子长度,m为有限差分系数的序号,1≤m≤M,n为连乘的变量。

    3.  如权利要求1所述的方法,其特征在于,基于时空域频散关系和非线性反演算法对有限差分系数进行优化,包括:
    将有限差分系数作为初值确定P波和S波的频散大??;
    根据确定的P波和S波的频散计算共轭梯度矢量;
    根据所述共轭梯度矢量迭代对所述有限差分系数进行优化。

    4.  如权利要求3所述的方法,其特征在于:
    将有限差分系数作为初值确定的P波和S波的频散大小为:


    其中,A=(Σm=1Mamsin((m-0.5)κcosθ))2,B=(Σm=1Mamsin((m-0.5)κsinθ))2,]]>为P波的频散大小,为S波的频散大小,M为算子长度,κ=kh,k为波数,h为空间采样间隔,rP=vPΔt/h,rS=vSΔt/h,vP为P波的传播速度,vS为S波的传播速度,Δt为时间采样间隔,θ为波传播方向角度,am为有限差分系数;
    根据P波和S波的频散大小计算共轭梯度矢量包括:
    确定目标函数:

    其中,b1=max(2πfmaxvPh,π),b2=max(2πfmaxvSh,π),]]>fmax为最大频率;
    根据所述目标函数计算共轭梯度矢量:
    pk+1=-▿Ek+1+▿Ek+1T▿Ek+1▿EkT▿Ekpk,p0=-▿E0]]>
    其中,pk表示当前时刻的共轭梯度矢量,▿E=∂E∂a1∂E∂a2···∂E∂aMT,]]>k表示当前时刻,k+1表示当前时刻的下一时刻;
    按照以下公式,根据所述共轭梯度矢量迭代对有限差分系数进行优化:
    ak+1=ak+αkpk+1,a=[a1a2…aM]T,其中,a为不同时刻优化后的有限差分系数向量,a1,a2…aM为优化后的有限差分系数,αk为迭代步长。

    5.  如权利要求4所述的方法,其特征在于,基于时空域频散关系和非线性反演算法对有限差分系数进行优化,包括:
    基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    对优化后的有限差分系数进行校验;
    如果校验结果不满足约束条件,则改变算子长度,再根据改变后的算子长度对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件;
    或者,如果校验结果不满足约束条件,则改变最大频率值,再根据改变后的最大频率值对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件。

    6.  如权利要求1至5中任一项所述的方法,利用优化后的有限差分系数进行弹性波正演模拟,包括:
    将优化后的有限差分系数代入二维弹性介质速度应力方程以实现弹性波正演模拟。

    7.  如权利要求6所述的方法,其特征在于,所述二维弹性介质速度应力方程为:
    ρvxi,jn+1-vxi,jnΔt=1hΣm=1Mam(τxxi+m-1/2,jn+1/2-τxxi-m+1/2,jn+1/2)+(τxzi,j+m-1/2n+1/2-τxzi,j-m+1/2n+1/2)]]>
    ρvzi+1/2,j+1/2n+1-vzi+1/2,j+1/2nΔt=1hΣm=1Mam(τxzi+m,j+1/2n+1/2-τxzi-m+1,j+1/2n+1/2)+(τzzi+1/2,j+mn+1/2-τzzi+1/2,j-m+1n+1/2)]]>
    τxxi+1/2,jn+3/2-τxxi+1/2,jn+1/2Δt=1hΣm=1Mam(λ+2μ)(vxi+m,jn+1-vxi-m+1,jn+1)+λ(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τzzi+1/2,jn+3/2-τzzi+1/2,jn+1/2Δt=1hΣm=1Mamλ(vxi+m,jn+1-vxi-m+1,jn+1)+(λ+2μ)(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τxzi,j+1/2n+3/2-τxzi,j+1/2n+1/2Δt=1hΣm=1Mamμ(vzi+m-1/2,jn+1-vzi-m+1/2,jn+1)+μ(vxi,j+mn+1-vxi,j-m+1n+1)]]>
    其中,(vx,vz)为速度矢量,(τxx,τzz,τxz)为应力的三个分量,λ(x,z)和μ(x,z)为拉梅系数,ρ(x,z)为密度,Δt为时间采样间隔,h为空间采样间隔,am为优化后的有限差分系数,M为算子长度,x为沿水平方向的坐标,z为沿铅垂方向的坐标,i为沿水平方向的采样序号,j为沿铅垂方向的采样序号,n为沿时间方向的采样序号,m为有限差分系数的序号,1≤m≤M。

    8.  一种基于非线性优化的时空域交错网格有限差分装置,其特征在于,包括:
    确定???,用于确定有限差分系数;
    优化???,用于基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    模拟???,用于利用优化后的有限差分系数进行弹性波正演模拟。

    9.  如权利要求8所述的装置,其特征在于,所述确定??榫咛逵糜诎凑找韵鹿饺范ㄓ邢薏罘窒凳?BR>am=(-1)m+12m-1Π1nM,n≠m|(2n-1)2(2n-1)2-(2m-1)2|]]>
    其中,am为有限差分系数,M为算子长度,m为有限差分系数的序号,1≤m≤M,n为连乘的变量。

    10.  如权利要求8所述的装置,其特征在于,所述优化??榘ǎ?BR>确定单元,用于将有限差分系数作为初值确定P波和S波的频散大??;
    计算单元,用于根据确定的P波和S波的频散计算共轭梯度矢量;
    优化单元,用于根据所述共轭梯度矢量迭代对所述有限差分系数进行优化。

    11.  如权利要求10所述的装置,其特征在于:
    所述确定单元具体用于按照以下公式确定P波和S波的频散大?。?BR>

    其中,A=(Σm=1Mamsin((m-0.5)κcosθ))2,B=(Σm=1Mamsin((m-0.5)κsinθ))2,]]>为P波的频散大小,为S波的频散大小,M为算子长度,κ=kh,k为波数,h为空间采样间隔,rP=vPΔt/h,rS=vSΔt/h,vP为P波的传播速度,vS为S波的传播速度,Δt为时间采样间隔,θ为波传播方向角度,am为有限差分系数;
    所述计算单元具体按照以下方式计算共轭梯度矢量:
    确定目标函数:

    其中,b1=max(2πfmaxvPh,π),b2=max(2πfmaxvSh,π),]]>fmax为最大频率;
    根据所述目标函数计算共轭梯度矢量:
    pk+1=-▿Ek+1+▿Ek+1T▿Ek+1▿EkT▿Ekpk,p0=-▿E0]]>
    其中,pk表示当前时刻的共轭梯度矢量,▿E=∂E∂a1∂E∂a2···∂E∂aMT,]]>k表示当前时刻,k+1表示当前时刻的下一时刻;
    所述优化单元具体用于按照以下公式对有限差分系数进行优化:
    ak+1=ak+αkpk+1,a=[a1a2…aM]T,其中,a为不同时刻优化后的有限差分系数向量,a1,a2…aM为优化后的有限差分系数,αk为迭代步长。

    12.  如权利要求11所述的装置,其特征在于,所述优化??榘ǎ?BR>系数优化单元,用于基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    校验单元,用于对优化后的有限差分系数进行校验,如果校验结果不满足约束条件,则改变算子长度,再根据改变后的算子长度对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件;或者,如果校验结果不满足约束条件,则改变最大频率值,再根据改变后的最大频率值对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件。

    13.  如权利要求8至12中任一项所述的装置,所述模拟??榫咛逵糜诮呕蟮挠邢薏罘窒凳攵越橹仕俣扔αΨ匠桃允迪值圆ㄕ菽D?。

    14.  如权利要求13所述的装置,其特征在于,所述二维弹性介质速度应力方程为:
    ρvxi,jn+1-vxi,jnΔt=1hΣm=1Mam(τxxi+m-1/2,jn+1/2-τxxi-m+1/2,jn+1/2)+(τxzi,j+m-1/2n+1/2-τxzi,j-m+1/2n+1/2)]]>
    ρvzi+1/2,j+1/2n+1-vzi+1/2,j+1/2nΔt=1hΣm=1Mam(τxzi+m,j+1/2n+1/2-τxzi-m+1,j+1/2n+1/2)+(τzzi+1/2,j+mn+1/2-τzzi+1/2,j-m+1n+1/2)]]>
    τxxi+1/2,jn+3/2-τxxi+1/2,jn+1/2Δt=1hΣm=1Mam(λ+2μ)(vxi+m,jn+1-vxi-m+1,jn+1)+λ(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τzzi+1/2,jn+3/2-τzzi+1/2,jn+1/2Δt=1hΣm=1Mamλ(vxi+m,jn+1-vxi-m+1,jn+1)+(λ+2μ)(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τxzi,j+1/2n+3/2-τxzi,j+1/2n+1/2Δt=1hΣm=1Mamμ(vzi+m-1/2,jn+1-vzi-m+1/2,jn+1)+μ(vxi,j+mn+1-vxi,j-m+1n+1)]]>
    其中,(vx,vz)为速度矢量,(τxx,τzz,τxz)为应力的三个分量,λ(x,z)和μ(x,z)为拉梅系数,ρ(x,z)为密度,Δt为时间采样间隔,h为空间采样间隔,am为优化后的有限差分系数,M为算子长度,x为沿水平方向的坐标,z为沿铅垂方向的坐标,i为沿水平方向的采样序号,j为沿铅垂方向的采样序号,n为沿时间方向的采样序号,m为有限差分系数的序号,1≤m≤M。

    说明书

    说明书基于非线性优化的时空域交错网格有限差分方法和装置
    技术领域
    本发明涉及正演模拟技术领域,特别涉及一种基于非线性优化的时空域交错网格有限差分方法和装置。
    背景技术
    地震数值模拟技术就是对特定的地质、地球物理问题作适当的简化,以形成简化的数学模型,然后采用数值计算方法获取地震响应的过程。地震数值模拟技术是理解地震波在地下传播特点,帮助解释观测数据的有效手段。地震数值模拟还可以为新技术的提出、可行性分析和应用试验提供高质量的模拟数据;帮助地球物理工作者测试新的算法和处理技术,为地震反演问题提供思路和有效的验证数据。近年来,波动方程数值模拟方法被广泛应用于逆时偏移和全波形反演中。
    常用的数值模拟方法主要包括:有限元素法、有限差分法和伪谱法等。其中,有限差分法是偏微分方程的主要数值解法之一,也是最早出现的数值模拟方法,其主要优点是物理意义直观,易于实现,能够较精确地模拟任意非均匀介质中的地震波场。有限差分法根据不同的标准可以分为:显式有限差分、隐式有限差分、规则网格有限差分、交错网格有限差分和旋转交错网格有限差分。在有限差分法中,差分系数可以通过泰勒级数展开或最优化方法求得,分别对应以泰勒级数展开为基础的有限差分和以最优化为基础的有限差分。
    然而,使用基于泰勒级数展开和空间域频散关系的有限差分法存在如下问题:在低频段,频散接近于零,然而在中高频段,频散较大,从而导致模拟精度较低。
    发明内容
    本发明实施例提供了一种基于非线性优化的时空域交错网格有限差分方法,以达到减小中高频段的频散,提高模拟精度的目的,该方法包括:
    确定有限差分系数;
    基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    利用优化后的有限差分系数进行弹性波正演模拟。
    在一个实施例中,所述确定有限差分系数,包括:
    按照以下公式确定有限差分系数:
    am=(-1)m+12m-1Π1nM,n≠m|(2n-1)2(2n-1)2-(2m-1)2|]]>
    其中,am为有限差分系数,M为算子长度,m为有限差分系数的序号,1≤m≤M,n为连乘的变量。
    在一个实施例中,基于时空域频散关系和非线性反演算法对有限差分系数进行优化,包括:
    将有限差分系数作为初值确定P波和S波的频散大??;
    根据确定的P波和S波的频散计算共轭梯度矢量;
    根据所述共轭梯度矢量迭代对所述有限差分系数进行优化。
    在一个实施例中,将有限差分系数作为初值确定的P波和S波的频散大小为:


    其中,A=(Σm=1Mamsin((m-0.5)κcosθ))2,B=(Σm=1Mamsin((m-0.5)κsinθ))2,]]>为P波的频散大小,为S波的频散大小,M为算子长度,κ=kh,k为波数,h为空间采样间隔,rP=vPΔt/h,rS=vSΔt/h,vP为P波的传播速度,vS为S波的传播速度,Δt为时间采样间隔,θ为波传播方向角度,am为有限差分系数;
    根据P波和S波的频散大小计算共轭梯度矢量包括:
    确定目标函数:

    其中,b1=max(2πfmaxvPh,π),b2=max(2πfmaxvSh,π),]]>fmax为最大频率;
    根据所述目标函数计算共轭梯度矢量:
    pk+1=-▿Ek+1+▿Ek+1T▿Ek+1▿EkT▿Ekpk,p0=-▿E0]]>
    其中,pk表示当前时刻的共轭梯度矢量,▿E=∂E∂a1∂E∂a2···∂E∂aMT,]]>k表示当前时刻,k+1表示当前时刻的下一时刻;
    按照以下公式,根据所述共轭梯度矢量迭代对有限差分系数进行优化:
    ak+1=ak+αkpk+1,a=[a1a2…aM]T,其中,a为不同时刻优化后的有限差分系数向量,a1,a2…aM为优化后的有限差分系数,αk为迭代步长。
    在一个实施例中,基于时空域频散关系和非线性反演算法对有限差分系数进行优化,包括:
    基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    对优化后的有限差分系数进行校验;
    如果校验结果不满足约束条件,则改变算子长度,再根据改变后的算子长度对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件;
    或者,如果校验结果不满足约束条件,则改变最大频率值,再根据改变后的最大频率值对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件。
    在一个实施例中,利用优化后的有限差分系数进行弹性波正演模拟,包括:
    将优化后的有限差分系数代入二维弹性介质速度应力方程以实现弹性波正演模拟。
    在一个实施例中,所述二维弹性介质速度应力方程为:
    ρvxi,jn+1-vxi,jnΔt=1hΣm=1Mam(τxxi+m-1/2,jn+1/2-τxxi-m+1/2,jn+1/2)+(τxzi,j+m-1/2n+1/2-τxzi,j-m+1/2n+1/2)]]>
    ρvzi+1/2,j+1/2n+1-vzi+1/2,j+1/2nΔt=1hΣm=1Mam(τxzi+m,j+1/2n+1/2-τxzi-m+1,j+1/2n+1/2)+(τzzi+1/2,j+mn+1/2-τzzi+1/2,j-m+1n+1/2)]]>
    τxxi+1/2,jn+3/2-τxxi+1/2,jn+1/2Δt=1hΣm=1Mam(λ+2μ)(vxi+m,jn+1-vxi-m+1,jn+1)+λ(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τzzi+1/2,jn+3/2-τzzi+1/2,jn+1/2Δt=1hΣm=1Mamλ(vxi+m,jn+1-vxi-m+1,jn+1)+(λ+2μ)(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τxzi,j+1/2n+3/2-τxzi,j+1/2n+1/2Δt=1hΣm=1Mamμ(vzi+m-1/2,jn+1-vzi-m+1/2,jn+1)+μ(vxi,j+mn+1-vxi,j-m+1n+1)]]>
    其中,(vx,vz)为速度矢量,(τxx,τzz,τxz)为应力的三个分量,λ(x,z)和μ(x,z)为拉梅系数,ρ(x,z)为密度,Δt为时间采样间隔,h为空间采样间隔,am为优化后的有限差分系数,M为算子长度,x为沿水平方向的坐标,z为沿铅垂方向的坐标,i为沿水平方向的采样序号,j为沿铅垂方向的采样序号,n为沿时间方向的采样序号,m为有限差分系数的序号,1≤m≤M。
    本发明实施例还提供了一种基于非线性优化的时空域交错网格有限差分装置,以达到减小中高频段的频散,提高模拟精度的目的,该装置包括:
    确定???,用于确定有限差分系数;
    优化???,用于基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    模拟???,用于利用优化后的有限差分系数进行弹性波正演模拟。
    在一个实施例中,所述确定??榫咛逵糜诎凑找韵鹿饺范ㄓ邢薏罘窒凳?
    am=(-1)m+12m-1Π1nM,n≠m|(2n-1)2(2n-1)2-(2m-1)2|]]>
    其中,am为有限差分系数,M为算子长度,m为有限差分系数的序号,1≤m≤M,n为连乘的变量。
    在一个实施例中,所述优化??榘ǎ?
    确定单元,用于将有限差分系数作为初值确定P波和S波的频散大??;
    计算单元,用于根据确定的P波和S波的频散计算共轭梯度矢量;
    优化单元,用于根据所述共轭梯度矢量迭代对所述有限差分系数进行优化。
    在一个实施例中,所述确定单元具体用于按照以下公式确定P波和S波的频散大 ?。?


    其中,A=(Σm=1Mamsin((m-0.5)κcosθ))2,B=(Σm=1Mamsin((m-0.5)κsinθ))2,]]>为P波的频散大小,为S波的频散大小,M为算子长度,κ=kh,k为波数,h为空间采样间隔,rP=vPΔt/h,rS=vSΔt/h,vP为P波的传播速度,vS为S波的传播速度,Δt为时间采样间隔,θ为波传播方向角度,am为有限差分系数;
    所述计算单元具体按照以下方式计算共轭梯度矢量:
    确定目标函数:

    其中,b1=max(2πfmaxvPh,π),b2=max(2πfmaxvSh,π),]]>fmax为最大频率;
    根据所述目标函数计算共轭梯度矢量:
    pk+1=-▿Ek+1+▿Ek+1T▿Ek+1▿EkT▿Ekpk,p0=-▿E0]]>
    其中,pk表示当前时刻的共轭梯度矢量,▿E=∂E∂a1∂E∂a2···∂E∂aMT,]]>k表示当前时刻,k+1表示当前时刻的下一时刻;
    所述优化单元具体用于按照以下公式对有限差分系数进行优化:
    ak+1=ak+αkpk+1,a=[a1a2…aM]T,其中,a为不同时刻优化后的有限差分系数向量,a1,a2…aM为优化后的有限差分系数,αk为迭代步长。
    在一个实施例中,优化??榘ǎ?
    系数优化单元,用于基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    校验单元,用于对优化后的有限差分系数进行校验,如果校验结果不满足约束条 件,则改变算子长度,再根据改变后的算子长度对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件;或者,如果校验结果不满足约束条件,则改变最大频率值,再根据改变后的最大频率值对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件。
    在一个实施例中,所述模拟??榫咛逵糜诮呕蟮挠邢薏罘窒凳攵越橹仕俣扔αΨ匠桃允迪值圆ㄕ菽D?。
    在一个实施例中,所述二维弹性介质速度应力方程为:
    ρvxi,jn+1-vxi,jnΔt=1hΣm=1Mam(τxxi+m-1/2,jn+1/2-τxxi-m+1/2,jn+1/2)+(τxzi,j+m-1/2n+1/2-τxzi,j-m+1/2n+1/2)]]>
    ρvzi+1/2,j+1/2n+1-vzi+1/2,j+1/2nΔt=1hΣm=1Mam(τxzi+m,j+1/2n+1/2-τxzi-m+1,j+1/2n+1/2)+(τzzi+1/2,j+mn+1/2-τzzi+1/2,j-m+1n+1/2)]]>
    τxxi+1/2,jn+3/2-τxxi+1/2,jn+1/2Δt=1hΣm=1Mam(λ+2μ)(vxi+m,jn+1-vxi-m+1,jn+1)+λ(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τzzi+1/2,jn+3/2-τzzi+1/2,jn+1/2Δt=1hΣm=1Mamλ(vxi+m,jn+1-vxi-m+1,jn+1)+(λ+2μ)(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τxzi,j+1/2n+3/2-τxzi,j+1/2n+1/2Δt=1hΣm=1Mamμ(vzi+m-1/2,jn+1-vzi-m+1/2,jn+1)+μ(vxi,j+mn+1-vxi,j-m+1n+1)]]>
    其中,(vx,vz)为速度矢量,(τxx,τzz,τxz)为应力的三个分量,λ(x,z)和μ(x,z)为拉梅系数,ρ(x,z)为密度,Δt为时间采样间隔,h为空间采样间隔,am为优化后的有限差分系数,M为算子长度,x为沿水平方向的坐标,z为沿铅垂方向的坐标,i为沿水平方向的采样序号,j为沿铅垂方向的采样序号,n为沿时间方向的采样序号,m为有限差分系数的序号,1≤m≤M。
    在本发明实施例中,通过时空域频散关系和非线性反演算法对有限差分系数进行优化,并利用优化后的有限差分系数进行弹性波正演模拟,从而解决了现有技术中采用泰勒级数展开和空间域频散关系的有限差分法获得有限差分系数进行弹性波正演模拟而导致的中高频段频散较大,模拟精度较低的技术问题,达到了减小中高频段的频散,提高模拟精度的技术效果。
    附图说明
    此处所说明的附图用来提供对本发明的进一步理解,构成本申请的一部分,并不构成对本发明的限定。在附图中:
    图1是本发明实施例的基于非线性优化的时空域交错网格有限差分方法的流程图;
    图2是本发明实施例的进行有限差分系数优化的方法流程图;
    图3是本发明实施例的基于非线性优化的时空域交错网格有限差分装置的结构框图;
    图4是本发明实施例的基于非线性优化的时空域交错网格有限差分方法的具体流程图;
    图5是常规方法确定的不同传播角度时δ随频率的变化曲线示意图;
    图6是本发明实施例确定的不同传播角度时δ随频率的变化曲线示意图;
    图7是常规方法确定的不同算子长度时δ随频率的变化曲线;
    图8本发明实施例确定的不同算子长度时δ随频率的变化曲线。
    具体实施方式
    发明人发现,在常规的有限差分方法中,差分系数一般是通过极小化空间域的频散关系得到的,这种方式在低频段,频散误差接近于零,然而在中高频段,频散较大,不能很好地描述地震波在时空域传播的规律。为此,在本实施例中提出了一种基于非线性优化的时空域交错网格有限差分方法来进行弹性波正演模拟,该方法通过极小化时间域和空间域的频散关系,以及非线性反演算法来求取最优的差分系数,以达到减小在中高频段的频散,提高模拟精度的目的。
    在本发明实施例中,提出了一种基于非线性优化的时空域交错网格有限差分方法,如图1所示,包括以下步骤:
    步骤101:确定有限差分系数;
    步骤102:基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    步骤103:利用优化后的有限差分系数进行弹性波正演模拟。
    在上述实施例中,通过时空域频散关系和非线性反演算法对有限差分系数进行优化,并利用优化后的有限差分系数进行弹性波正演模拟,从而解决了现有技术中采用 泰勒级数展开和空间域频散关系的有限差分法获得有限差分系数进行弹性波正演模拟而导致的中高频段频散较大,模拟精度较低的技术问题,达到了减小中高频段的频散,提高模拟精度的技术效果。
    在上述步骤101中,可以采用常规的有限差分系数求取方法获得有限差分系数,例如可以按照以下公式计算有限差分系数:
    am=(-1)m+12m-1Π1nM,n≠m|(2n-1)2(2n-1)2-(2m-1)2|]]>
    其中,am为有限差分系数,M为算子长度,m为有限差分系数的序号,1≤m≤M,n为连乘的变量。
    具体的,上述步骤102基于时空域频散关系和非线性反演算法对有限差分系数进行优化可以如图2所示,包括以下步骤:
    步骤201:将有限差分系数作为初值确定P波和S波的频散大小,其中,P波是指纵波(Primary/Compressional wave),S波是指横波(Second/Shear wave):


    其中,A=(Σm=1Mamsin((m-0.5)κcosθ))2,B=(Σm=1Mamsin((m-0.5)κsinθ))2,]]>为P波的频散大小,为S波的频散大小,M为算子长度,κ=kh,k为波数,h为空间采样间隔,rP=vPΔt/h,rS=vSΔt/h,vP为P波的传播速度,vS为S波的传播速度,Δt为时间采样间隔,θ为波传播方向角度,am表示初始有限差分系数;
    步骤202:根据确定的P波和S波的频散计算共轭梯度矢量;
    先确定目标函数:

    其中,b1=max(2πfmaxvPh,π),b2=max(2πfmaxvSh,π),]]>fmax为最大频率;
    然后根据所述目标函数计算共轭梯度矢量:
    pk+1=-▿Ek+1+▿Ek+1T▿Ek+1▿EkT▿Ekpk,p0=-▿E0]]>
    其中,pk表示当前时刻的共轭梯度矢量,▿E=∂E∂a1∂E∂a2···∂E∂aMT,]]>k表示当前时刻,k+1表示当前时刻的下一时刻;
    步骤203:按照以下公式,根据所述共轭梯度矢量迭代对有限差分系数进行优化:
    ak+1=ak+αkpk+1,a=[a1a2…aM]T,其中,a表示不同时刻优化后的有限差分系数向量,a1,a2…aM为优化后的有限差分系数,αk表示迭代步长。
    考虑到算子长度M和最大频率在对有限差分系数的优化过程中起着重要的作用,可以通过对优化后的有限差分系数进行校验,判断其是否满足约束条件,如果不满足,可以更改算子长度M的值或者改变最大频率fmax的值,然后重新确定最佳有限差分系数。
    具体的可以包括:
    步骤1:基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    步骤2:对优化后的有限差分系数进行校验;
    步骤3:如果校验结果不满足约束条件,则改变算子长度,再根据改变后的算子长度对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件;或者,如果校验结果不满足约束条件,则改变最大频率值,再根据改变后的最大频率值对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件。
    在上述各个实施例中,上述步骤103利用优化后的有限差分系数进行弹性波正演模拟,可以包括:将优化后的有限差分系数代入二维弹性介质速度应力方程以实现弹性波正演模拟。
    其中,该二维弹性介质速度应力方程为:
    ρvxi,jn+1-vxi,jnΔt=1hΣm=1Mam(τxxi+m-1/2,jn+1/2-τxxi-m+1/2,jn+1/2)+(τxzi,j+m-1/2n+1/2-τxzi,j-m+1/2n+1/2)]]>
    ρvzi+1/2,j+1/2n+1-vzi+1/2,j+1/2nΔt=1hΣm=1Mam(τxzi+m,j+1/2n+1/2-τxzi-m+1,j+1/2n+1/2)+(τzzi+1/2,j+mn+1/2-τzzi+1/2,j-m+1n+1/2)]]>
    τxxi+1/2,jn+3/2-τxxi+1/2,jn+1/2Δt=1hΣm=1Mam(λ+2μ)(vxi+m,jn+1-vxi-m+1,jn+1)+λ(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τzzi+1/2,jn+3/2-τzzi+1/2,jn+1/2Δt=1hΣm=1Mamλ(vxi+m,jn+1-vxi-m+1,jn+1)+(λ+2μ)(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τxzi,j+1/2n+3/2-τxzi,j+1/2n+1/2Δt=1hΣm=1Mamμ(vzi+m-1/2,jn+1-vzi-m+1/2,jn+1)+μ(vxi,j+mn+1-vxi,j-m+1n+1)]]>
    其中,(vx,vz)为速度矢量,(τxx,τzz,τxz)为应力的三个分量,λ(x,z)和μ(x,z)为拉梅系数,ρ(x,z)为密度,Δt为时间采样间隔,h为空间采样间隔,am为优化的有限差分系数,M为算子长度,x为沿水平方向的坐标,z为沿铅垂方向的坐标,i为沿水平方向的采样序号,j为沿铅垂方向的采样序号,n为沿时间方向的采样序号,m为有限差分系数的序号,1≤m≤M。
    基于同一发明构思,本发明实施例中还提供了一种基于非线性优化的时空域交错网格有限差分装置,如下面的实施例所述。由于基于非线性优化的时空域交错网格有限差分装置解决问题的原理与基于非线性优化的时空域交错网格有限差分方法相似,因此基于非线性优化的时空域交错网格有限差分装置的实施可以参见基于非线性优化的时空域交错网格有限差分方法的实施,重复之处不再赘述。以下所使用的,术语“单元”或者“??椤笨梢允迪衷ざüδ艿娜砑?或硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现,但是硬件,或者软件和硬件的组合的实现也是可能并被构想的。图3是本发明实施例的基于非线性优化的时空域交错网格有限差分装置的一种结构框图,如图3所示,包括:确定???01、优化???02和模拟???03,下面对该结构进行具体说明。
    确定???01,用于确定有限差分系数;
    优化???02,用于基于时空域频散关系和非线性反演算法对有限差分系数进行优化;
    模拟???03,用于利用优化后的有限差分系数进行弹性波正演模拟。
    在一个实施例中,确定???01具体用于按照以下公式确定有限差分系数:
    am=(-1)m+12m-1Π1nM,n≠m|(2n-1)2(2n-1)2-(2m-1)2|]]>
    其中,am为有限差分系数,M为算子长度,m为有限差分系数的序号,1≤m≤M,n为连乘的变量。
    在一个实施例中,优化???02包括:确定单元,用于将有限差分系数作为初值确定P波和S波的频散大??;计算单元,用于根据确定的P波和S波的频散计算共轭梯度矢量;优化单元,用于根据所述共轭梯度矢量迭代对所述有限差分系数进行优化。
    在一个实施例中,确定单元具体用于按照以下公式确定P波和S波的频散大?。?


    其中,A=(Σm=1Mamsin((m-0.5)κcosθ))2,B=(Σm=1Mamsin((m-0.5)κsinθ))2,]]>为P波的频散大小,为S波的频散大小,M为算子长度,κ=kh,k为波数,h为空间采样间隔,rP=vPΔt/h,rS=vSΔt/h,vP为P波的传播速度,vS为S波的传播速度,Δt为时间采样间隔,θ为波传播方向角度,am为有限差分系数;
    计算单元具体按照以下方式计算共轭梯度矢量:确定目标函数:

    其中,b1=max(2πfmaxvPh,π),b2=max(2πfmaxvSh,π),]]>fmax为最大频率;
    根据所述目标函数计算共轭梯度矢量:
    pk+1=-▿Ek+1+▿Ek+1T▿Ek+1▿EkT▿Ekpk,p0=-▿E0]]>
    其中,pk表示当前时刻的共轭梯度矢量,▿E=∂E∂a1∂E∂a2···∂E∂aMT,]]>k表示当前时刻,k+1表示当前时刻的下一时刻;
    优化单元具体用于按照以下公式对有限差分系数进行优化:
    ak+1=ak+αkpk+1,a=[a1a2…aM]T,其中,a为不同时刻优化后的有限差分 系数向量,a1,a2…aM为优化后的有限差分系数,αk为迭代步长。
    在一个实施例中,优化???02包括:系数优化单元,用于基于时空域频散关系和非线性反演算法对有限差分系数进行优化;校验单元,用于对优化后的有限差分系数进行校验;如果校验结果不满足约束条件,则改变算子长度,再根据改变后的算子长度对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件;或者,如果校验结果不满足约束条件,则改变最大频率值,再根据改变后的最大频率值对有限差分系数进行优化,直至优化后的有限差分系数的校验结果可以满足所述约束条件。
    在一个实施例中,模拟???03具体用于将优化后的有限差分系数代入二维弹性介质速度应力方程以实现弹性波正演模拟。
    在一个实施例中,二维弹性介质速度应力方程为:
    ρvxi,jn+1-vxi,jnΔt=1hΣm=1Mam(τxxi+m-1/2,jn+1/2-τxxi-m+1/2,jn+1/2)+(τxzi,j+m-1/2n+1/2-τxzi,j-m+1/2n+1/2)]]>
    ρvzi+1/2,j+1/2n+1-vzi+1/2,j+1/2nΔt=1hΣm=1Mam(τxzi+m,j+1/2n+1/2-τxzi-m+1,j+1/2n+1/2)+(τzzi+1/2,j+mn+1/2-τzzi+1/2,j-m+1n+1/2)]]>
    τxxi+1/2,jn+3/2-τxxi+1/2,jn+1/2Δt=1hΣm=1Mam(λ+2μ)(vxi+m,jn+1-vxi-m+1,jn+1)+λ(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τzzi+1/2,jn+3/2-τzzi+1/2,jn+1/2Δt=1hΣm=1Mamλ(vxi+m,jn+1-vxi-m+1,jn+1)+(λ+2μ)(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τxzi,j+1/2n+3/2-τxzi,j+1/2n+1/2Δt=1hΣm=1Mamμ(vzi+m-1/2,jn+1-vzi-m+1/2,jn+1)+μ(vxi,j+mn+1-vxi,j-m+1n+1)]]>
    其中,(vx,vz)为速度矢量,(τxx,τzz,τxz)为应力的三个分量,λ(x,z)和μ(x,z)为拉梅系数,ρ(x,z)为密度,Δt为时间采样间隔,h为空间采样间隔,am为优化后的有限差分系数,M为算子长度,x为沿水平方向的坐标,z为沿铅垂方向的坐标,i为沿水平方向的采样序号,j为沿铅垂方向的采样序号,n为沿时间方向的采样序号,m为有限差分系数的序号,1≤m≤M。
    为了对上述弹性波正演模拟方法进行更为清楚的解释,下面结合一个具体的实施 例来进行说明,然而值得注意的是该实施例仅是为了更好地说明本发明,并不构成对本发明不当的限定。
    本发明在分析了常规有限差分法的不足的情况下、以时空域频散关系和非线性反演算法作为基础上提出了一种基于非线性优化的时空域交错网格有限差分法,以减小中高频段数值频散的大小、提高弹性波方程正演模拟的精度,该方法适用于弹性波方程数值模拟。
    如图4所示,是本发明实施例的基于非线性优化的时空域交错网格有限差分方法流程图,包括以下步骤:
    步骤401:根据有关参数对计算区域进行网格剖分;
    步骤402:求取常规有限差分法的差分系数:
    am=(-1)m+12m-1Π1nM,n≠m|(2n-1)2(2n-1)2-(2m-1)2|]]>
    步骤403:以常规有限差分系数为初值,基于时空域频散关系和非线性反演算法求取最佳有限差分系数(即优化后的有限差分系数);
    步骤404:将得到的差分系数代入二维弹性介质速度-应力方程,以实现弹性波正演模拟。
    其中,二维弹性介质速度-应力方程表示如下:
    ρ∂vx∂t=∂τxx∂x+∂τxz∂z]]>
    ρ∂vz∂t=∂τxz∂x+∂τzz∂z]]>
    ∂τxx∂t=(λ+2μ)∂vx∂x+λ∂vz∂z]]>
    ∂τzz∂t=λ∂vx∂x+(λ+2μ)∂vz∂z]]>
    ∂τxz∂t=μ(∂vz∂x+∂vx∂z)]]>
    其中,(vx,vz)为速度矢量,(τxx,τzz,τxz)为应力的三个分量,λ(x,z)和μ(x,z)为拉梅系数,ρ(x,z)为密度。
    ρvxi,jn+1-vxi,jnΔt=1hΣm=1Mam(τxxi+m-1/2,jn+1/2-τxxi-m+1/2,jn+1/2)+(τxzi,j+m-1/2n+1/2-τxzi,j-m+1/2n+1/2)]]>
    ρvzi+1/2,j+1/2n+1-vzi+1/2,j+1/2nΔt=1hΣm=1Mam(τxzi+m,j+1/2n+1/2-τxzi-m+1,j+1/2n+1/2)+(τzzi+1/2,j+mn+1/2-τzzi+1/2,j-m+1n+1/2)]]>
    τxxi+1/2,jn+3/2-τxxi+1/2,jn+1/2Δt=1hΣm=1Mam(λ+2μ)(vxi+m,jn+1-vxi-m+1,jn+1)+λ(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τzzi+1/2,jn+3/2-τzzi+1/2,jn+1/2Δt=1hΣm=1Mamλ(vxi+m,jn+1-vxi-m+1,jn+1)+(λ+2μ)(vzi+1/2,j+m-1/2n+1-vzi+1/2,j-m+1/2n+1)]]>
    τxzi,j+1/2n+3/2-τxzi,j+1/2n+1/2Δt=1hΣm=1Mamμ(vzi+m-1/2,jn+1-vzi-m+1/2,jn+1)+μ(vxi,j+mn+1-vxi,j-m+1n+1)]]>
    其中,Δt是时间采样间隔,h为空间采样间隔,am为优化后的最佳有限差分系数,M为算子长度。
    下面对图4所示的弹性波方程数值模拟流程进行具体描述。
    弹性波方程交错网格模拟,时空域频散关系如下:
    A+B≈rP-2sin2(0.5κrP)]]>
    A+B≈rS-2sin2(0.5κrS)]]>
    其中:
    A=(Σm=1Mamsin((m-0.5)κcosθ))2]]>
    B=(Σm=1Mamsin((m-0.5)κsinθ))2]]>
    其中,κ=kh,k为波数,rP=vPΔt/h,rS=vSΔt/h,vP和vS分别为P波和S波的传播速度。
    定义两个新函数:


    其中,和为上述公式两边的相对误差,用于分别表示P波和S波的频散大小,如果和都为0,则表示没有频散,通过极小化和可以求得最佳的交错网格的有限差分系数。
    定义目标函数为:

    对于给定的频率范围(0≤f≤fmax),P波和S波的κ是不同的,计算公式如下:
    b1=max(2πfmaxvPh,π),b2=max(2πfmaxvSh,π),]]>
    可以采用非线性共轭梯度法来求解该优化问题,共轭梯度矢量为:
    pk+1=-▿Ek+1+▿Ek+1T▿Ek+1▿EkT▿Ekpk,p0=-▿E0]]>
    其中:
    ▿E=∂E∂a1∂E∂a2···∂E∂aMT,]]>k和k+1分别表示当前时刻和下一时刻。
    迭代公式为:
    ak+1=ak+αkpk+1,a=[a1a2…aM]T,其中,αk为迭代步长。
    考虑到算子长度M在弹性波正演模拟方法中的最佳有限差分系数的确定过程中起着重要的作用,可以通过对求取的最佳有限差分系数进行校验,使其满足约束条件,如果不满足,可以更改算子长度M的值,然后重新确定最佳有限差分系数。具体的可以采用下面的公式来度量数值频散的大?。?
    δP(θ)=vPFDvP=2rPκarcsin(rPA+B)]]>
    δS(θ)=vSFDvS=2rPκarcsin(rSA+B)]]>
    其中,δP(θ)表示P波的相速度误差,δS(θ)表示S波的相速度误差,δP和δS越接近于1表示频散越小。
    在本例中,定义以下公式来平衡纵波和横波的频散误差:
    δ(f,θ,M)=0.5(δP-1)+0.5(δS-1)
    该公式可以近似表示纵波和横波相速度的相对误差,δ越接近于0表示频散越小。
    最大频散误差δmax如下:
    δmax(M,fmax)=maxf∈[0,fmax],θ∈[0,2π]|δ(f,θ,M)|]]>
    为了保证本文方法的有效性,在前面的最优化过程中加入约束条件:
    δmax<η,其中,η为最大允许误差。
    其中,δmax只与fmax和M有关。当M固定时,δmax决定于fmax,δmax随着fmax增加而变大。因此,若初始的fmax不满足上述约束条件,则通过逐渐减小fmax直至δmax<η来获得最佳fmax。另外,当fmax固定时,δmax只与M有关,δmax随着M增加而变小,因此如果初始的M不满足上述约束条件,则通过逐渐增大M直到δmax<η来获得最佳M。
    最终将确定的最佳有限差分系数代入到上述二维弹性介质速度-应力方程中就实现了弹性波数值模拟。
    以一个均匀弹性介质为例来说明本发明实施例的优点,在本例中,纵波速度为2800m/s,横波速度1500m/s,时间采样间隔为1ms,空间采用间隔为h=20m,算子长度M=14。
    图5和图6表示在不同差分算法被使用时数值频散随传播方向的变化规律,图7和图8表示在不同差分算法被使用时数值频散随算子长度的变化规律。由图5至图8可见,与常规交错网格有限差分相比,基于非线性最优化的时空域交错网格有限差分法具有更宽的有效频带及更小的数值频散。通过当算子长度相同时,基于非线性最优化的时空域交错网格有限差分法模拟的精度更高。
    本实施例中减小了中高频段的数值频散,以时空域频散关系为基础更符合实际,模拟精度更高,以常规有限差分法的差分系数为初值,结合非线性共轭梯度算法可以在较少的迭代次数下得到满意的结果。该方法可以应用到所有涉及到数值模拟的地球物理研究过程中,例如:逆时偏移、全波形反演,可以大大减少中高频段数值的频散,提高弹性波方程正演模拟的精度。
    在另外一个实施例中,还提供了一种软件,该软件用于执行上述实施例及优选实施方式中描述的技术方案。
    在另外一个实施例中,还提供了一种存储介质,该存储介质中存储有上述软件,该存储介质包括但不限于:光盘、软盘、硬盘、可擦写存储器等。
    从以上的描述中,可以看出,本发明实施例实现了如下技术效果:通过时空域频散关系和非线性反演算法对有限差分系数进行优化,并利用优化后的有限差分系数进行弹性波正演模拟,从而解决了现有技术中采用泰勒级数展开和空间域频散关系的有限差分法获得有限差分系数进行弹性波正演模拟而导致的中高频段频散较大,模拟精度较低的技术问题,达到了减小中高频段的频散,提高模拟精度的技术效果。
    显然,本领域的技术人员应该明白,上述的本发明实施例的各??榛蚋鞑街杩梢杂猛ㄓ玫募扑阕爸美词迪?,它们可以集中在单个的计算装置上,或者分布在多个计算装置所组成的网络上,可选地,它们可以用计算装置可执行的程序代码来实现,从而,可以将它们存储在存储装置中由计算装置来执行,并且在某些情况下,可以以不同于此处的顺序执行所示出或描述的步骤,或者将它们分别制作成各个集成电路???,或者将它们中的多个??榛虿街柚谱鞒傻ジ黾傻缏纺?槔词迪?。这样,本发明实施例不限制于任何特定的硬件和软件结合。
    以上所述仅为本发明的优选实施例而已,并不用于限制本发明,对于本领域的技术人员来说,本发明实施例可以有各种更改和变化。凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的?;し段е??!  ∧谌堇醋宰ɡ鴚ww.www.4mum.com.cn转载请标明出处

    关于本文
    本文标题:基于非线性优化的时空域交错网格有限差分方法和装置.pdf
    链接地址://www.4mum.com.cn/p-5779458.html
    关于我们 - 网站声明 - 网站地图 - 资源地图 - 友情链接 - 网站客服 - 联系我们

    [email protected] 2017-2018 www.4mum.com.cn网站版权所有
    经营许可证编号:粤ICP备17046363号-1 
     


    收起
    展开
  • 四川郎酒股份有限公司获第十二届人民企业社会责任奖年度环保奖 2019-05-13
  • 银保监会新规剑指大企业多头融资和过度融资 2019-05-12
  • 韩国再提4国联合申办世界杯 中国网友无视:我们自己来 2019-05-11
  • 中国人为什么一定要买房? 2019-05-11
  • 十九大精神进校园:风正扬帆当有为 勇做时代弄潮儿 2019-05-10
  • 粽叶飘香幸福邻里——廊坊市举办“我们的节日·端午”主题活动 2019-05-09
  • 太原设禁鸣路段 设备在测试中 2019-05-09
  • 拜耳医药保健有限公司获第十二届人民企业社会责任奖年度企业奖 2019-05-08
  • “港独”没出路!“梁天琦们”该醒醒了 2019-05-07
  • 陈卫平:中国文化内涵包含三方面 文化复兴表现在其中 2019-05-06
  • 人民日报客户端辟谣:“合成军装照”产品请放心使用 2019-05-05
  • 【十九大·理论新视野】为什么要“建设现代化经济体系”?   2019-05-04
  • 聚焦2017年乌鲁木齐市老城区改造提升工程 2019-05-04
  • 【专家谈】上合组织——构建区域命运共同体的有力实践者 2019-05-03
  • 【华商侃车NO.192】 亲!楼市火爆,别忘了买车位啊! 2019-05-03
  • 彩票平台快捷支付 手机支付宝赚钱软件 齐鲁风采30选5 中国体彩网官方网站22选5 360足彩胜负彩 福建快3N期有无号 排列五推荐 云南十一选五开奖今天 英超联赛利物浦对切尔西 钱柜彩票网址 辽宁快乐12选5预测 中国一重股票分析 宁夏十一选五游戏规则 福彩3d乐彩论坛 甘肃快三今日开奖结果 山西泳坛夺金开奖号码