CN111881613B - 一种正应力和剪应力不同权重三维应力场反演方法及系统 - Google Patents

一种正应力和剪应力不同权重三维应力场反演方法及系统 Download PDF

Info

Publication number
CN111881613B
CN111881613B CN202010775092.6A CN202010775092A CN111881613B CN 111881613 B CN111881613 B CN 111881613B CN 202010775092 A CN202010775092 A CN 202010775092A CN 111881613 B CN111881613 B CN 111881613B
Authority
CN
China
Prior art keywords
stress
value
inversion
measuring point
normal
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
Application number
CN202010775092.6A
Other languages
English (en)
Other versions
CN111881613A (zh
Inventor
裴启涛
和礼红
熊永华
彭定新
陈昆鹏
陈渊
汤勇
付佩
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Wuhan Municipal Engineering Design and Research Institute Co Ltd
Original Assignee
Wuhan Municipal Engineering Design and Research Institute Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Wuhan Municipal Engineering Design and Research Institute Co Ltd filed Critical Wuhan Municipal Engineering Design and Research Institute Co Ltd
Priority to CN202010775092.6A priority Critical patent/CN111881613B/zh
Publication of CN111881613A publication Critical patent/CN111881613A/zh
Application granted granted Critical
Publication of CN111881613B publication Critical patent/CN111881613B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/10—Geometric CAD
    • G06F30/13—Architectural design, e.g. computer-aided architectural design [CAAD] related to design of buildings, bridges, landscapes, production plants or roads
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Geometry (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • General Physics & Mathematics (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Architecture (AREA)
  • Civil Engineering (AREA)
  • Structural Engineering (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Investigating Strength Of Materials By Application Of Mechanical Stress (AREA)

Abstract

本发明涉及一种正应力和剪应力不同权重三维应力场反演方法及系统，其方法包括获取不同测点的应力矢量及岩体弹性模量；根据工程地质信息建立有限元数值分析模型，并分别施加不同方向的荷载，分别确定对应荷载作用下测点处的应力计算值；建立地应力场反演目标函数模型并进行反演分析，获得不同权重下测点处的应力计算值和应力回归系数；建立地应力场反演效果评价的误差公式并确定最优的目标地应力场。本发明解决了三维应力场中剪应力与正应力间由于数值量级差异引起的剪应力易被忽略的问题，反演获得的地应力场精度大大提高，可准确、可靠、高效地获取地应力场分布特征，为科学识别开挖后的硐室围岩变形破坏特征提供技术支撑。

Description

一种正应力和剪应力不同权重三维应力场反演方法及系统
技术领域
本发明涉及岩土工程和数值分析技术领域，尤其涉及一种正应力和剪应力不同权重三维应力场反演方法及系统。
背景技术
随着社会经济的快速发展及西部大开发战略的实施，我国在水利水电工程、交通、国防等领域的建设正密集展开，与此相关的地下工程越来越多，并呈现出“深埋、超长、特大”的发展趋势。工程实践表明，地应力是引起地下工程围岩变形和破坏的根本作用力，伴随着工程安全问题的日益重视，预先对工程区地应力分布特征进行研究是一项必备的工作，已成为当今岩石力学界研究的热点和重要课题之一。
由于地应力成因的复杂性和不确定性，人们很难准确地认识其形成机理，只能通过外部的观测信息、现场有限的实测资料及地质构造条件，采用有效的分析方法进行岩体应力场反演。目前国内外模拟岩体应力场的方法主要包括位移反分析方法和应力反演分析方法，对于地质条件较为复杂、地应力场非线性特征较为突出区域，常常采用应力场反演分析法。然而，在采用应力反演法进行分析时，现有的研究成果在地应力场反演目标函数的建立方面，较常采用离差平方和来度量两个应力张量之间的接近程度。由于应力张量中的剪应力较之正应力要小得多，尤其对于深埋工程，甚至存在数量级的差异，若仍然采用该方法势必会导致剪应力拟合处于次要甚至被忽略不计的地位，从而容易产生对硐室围岩破坏特征的误判，对工程安全建设危害极大。然而，该方面现有成果相对较少，尚未建立起合理反映正应力和剪应力权重影响的地应力场反演方法。可见，如何构建正应力与剪应力不等权重下的应力场反演目标函数，并结合简易的室内试验、有限的现场测试成果及三维地质数值模型，准确、高效地获取地下工程应力场分布特征，为科学识别开挖后的硐室围岩变形破坏特征，以及永久支护设计的优化提供必备的基础资料，是急需解决的关键技术问题。
发明内容
本发明所要解决的技术问题是针对上述现有技术的不足，提供一种正应力和剪应力不同权重三维应力场反演方法及系统。
本发明解决上述技术问题的技术方案如下：一种正应力和剪应力不同权重三维应力场反演方法，包括如下步骤：
通过原位应力测试和承压板载荷试验，获取三维空间不同空间部位测点的力学参数测量信息；
根据工程地质信息建立有限元数值分析模型，并在三维空间内不同方向上分别施加荷载，分别确定对应荷载作用下测点处的应力计算值；
建立正应力和剪应力不同权重的三维地应力场反演目标函数模型：
Figure GDA0002944230280000021
其中，S残'为残余误差值，m为测点的总数,
Figure GDA0002944230280000022
为第k个测点、第i个正应力实测值，
Figure GDA0002944230280000023
为第k个测点、第j个剪应力实测值，w为正应力权重，w的范围为[0,1]；
Figure GDA0002944230280000024
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure GDA0002944230280000025
为第k个测点在对应荷载作用下的第j个剪应力计算值；
分别采用正应力和剪应力不同权重下的三维地应力场反演目标函数对目标地应力场进行反演分析，获得不同权重下测点处的应力计算值和应力回归系数；
建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场。
本发明的有益效果是：本发明的正应力和剪应力不同权重三维应力场反演方法，通过在有限元数值分析模型上施加不用方向的荷载，获取测点处的应力计算值，并进行反演分析，得到不同权重下测点处的应力计算值和应力回归系数，结合试验测得的应力矢量进行反演评价，从而准确确定最优目标地应力场，理论上解决了剪应力与正应力之间由于数值量级差异引起的剪应力在反演计算中处于次要甚至被忽略不计的问题，能够反映正应力和剪应力不等权重对地应力分布特征的影响，反演获得的应力矢量精度大大提高。通过采用本文的反演目标函数模型及反演方法对工程区三维初始应力场进行反演分析，可以准确、可靠、高效地获取地应力场分布特征，为科学识别开挖后的硐室围岩变形破坏特征，以及永久支护设计的优化提供必备的基础资料和技术支撑。
在上述技术方案的基础上，本发明还可以做如下改进：
进一步：所述通过原位应力测试和承压板载荷试验，获取三维空间不同空间部位测点的力学参数测量信息包括如下步骤：
所述通过原位应力测试和承压板载荷试验，获取不同空间部位测点的应力测量参数信息具体包括如下步骤：
在施工现场勘探平硐内进行原位应力测试，获取分析断面处的不同部位m个测点的应力测量值
Figure GDA0002944230280000032
和应力方向；
在施工现场勘探平硐内进行承压板载荷试验，采用逐级一次循环法分级施加荷载，获取现场试验实测的压力-变形(P-W)关系曲线，利用半无限空间弹性体公式求得岩体弹性模量E，具体公式如下：
Figure GDA0002944230280000031
其中，当l＝1,2,3时，
Figure GDA0002944230280000041
为第k个测点、第l个正应力测量值，当l＝4,5,6时，
Figure GDA0002944230280000042
为为第k个测点、第l个剪应力测量值，W为岩体弹性变形，单位为cm；P为按承压板面积计算的压力，单位为MPa；D为承压板直径单位为cm；ν为泊松比。
上述进一步方案的有益效果是：通过进行原位应力测试，可以准确获取不同测点的应力测量值
Figure GDA0002944230280000043
和岩体弹性模量，这样便于后续根据应力测量值以及反演得到的应力计算值评价反演精度，同时，通过所述岩体弹性模量便于对所述有限元数值分析模型进行初始化。
进一步：所述根据工程地质信息建立有限元数值分析模型，并在三维空间内不同方向上分别施加荷载，分别确定对应荷载作用下测点处的应力计算值具体包括如下步骤：
根据工程地质信息，建立有限元数值分析模型，并对所述有限元数值分析模型进行有限元网格划分，并根据所述岩体弹性模量E初始化所述有限元数值分析模型的参数信息；
沿所述有限元数值分析模型边界分别施加重力荷载、水平x方向单位挤压荷载、水平y方向单位挤压荷载及水平面内单位剪切荷载，通过有限元数值处理，并分别提取对应载荷作用下测点处的应力计算值。
上述进一步方案的有益效果是：通过建立有限元数值分析模型，并进行网格划分、初始化，这样在不同荷载的作用下可以分别提取对应的应力计算值，便于后续进行反演运算。
进一步：所述建立正应力和剪应力不同权重的地应力场三维反演目标函数模型具体包括如下步骤：
以测点处的所述应力测量值与应力计算值之差的平方和达到最小为优化目标，建立m个测点下的正应力与剪应力相互独立的三维反演目标函数模型：
Figure GDA0002944230280000051
其中：m为测点的总数；当l＝1,2,3时，
Figure GDA0002944230280000052
分别为第k个测点、第l个正应力测量值和计算值；当l＝4,5,6时，
Figure GDA0002944230280000053
分别为第k个测点、第l个剪应力测量值和计算值，
Figure GDA0002944230280000054
分别为第k个测点、第i个正应力测量值和计算值，
Figure GDA0002944230280000055
分别为第k个测点剪应力测量值和计算值；
根据上述公式(2)建立正应力权重为w时的三维地应力场反演目标函数模型，具体如下：
假设第k个测点、第i个正应力σik权重为w，相应的第k个测点、第j个剪切应力τjk权重为1-w，代入公式(2)中，可得：
Figure GDA0002944230280000056
其中，w为0～1.0之间的任意值；
根据上述公式(3)对应力测量值和计算值进行形式变换，建立一般形式下的三维地应力场反演目标函数模型，可得：
Figure GDA0002944230280000057
其中，
Figure GDA0002944230280000058
为第k个测点、第i个正应力实测值，
Figure GDA0002944230280000059
为第k个测点、第j个剪应力实测值，，w为正应力权重，w的范围为[0,1]；
Figure GDA00029442302800000510
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure GDA00029442302800000511
为第k个测点在对应荷载作用下的第j个剪应力计算值。
上述进一步方案的有益效果是：通过建立测点的正应力与剪应力相互独立的三维反演目标函数模型，并结合应力测量值和计算值进行形式变换，得到一般形式下的三维地应力场反演目标函数模型，便于后续针对应力测量值和应力计算值进行反演分析以及反演评价，得到准确的最优目标地应力场。
进一步：所述获得不同权重下测点处的应力计算值和应力回归系数具体包括如下步骤：
根据对应载荷下不同权重正应力的应力计算值和应力回归系数计算回归应力值，计算公式如下：
Figure GDA0002944230280000061
其中，当l＝1,2,3时，
Figure GDA0002944230280000062
为第k个测点、第j个正应力回归计算值；当l＝4,5,6时，
Figure GDA0002944230280000063
为第k个测点、第j个剪应力回归计算值，Li为相应于自变量的多元回归系数，当l＝1,2,3时，
Figure GDA0002944230280000064
为第k个测点、第j个正应力计算值的单列矩阵；当l＝4,5,6时，
Figure GDA0002944230280000065
为第k个测点、第j个剪应力计算值的单列矩阵；
根据多元回归分析理论，将上式(5)代入式(4)的公式中，可得：
Figure GDA0002944230280000066
根据最小二乘法原理，S残’为最小值的法方程式为：
Figure GDA0002944230280000067
采用线性回归分析法对上述法方程式(7)进行求解，得到正应力权重为w时的三个应力计算值和对应的4个应力回归系数L＝(L1,L2,...,Ln)T；
根据所述应力回归系数计算计算域内任一点P的回归应力值，具体公式如下：
Figure GDA0002944230280000068
其中，当l＝1,2,3时，σjp’为测点P处、第j个正应力计算值；当l＝4,5,6时，σjp’为测点P处、第j个剪应力计算值，j＝1,2,…,6对应初始应力六个分量。
上述进一步方案的有益效果是：通过将对应载荷下不同权重正应力的应力计算值和应力回归系数进行线性叠加，可以得到回归应力值，并结合最小二乘法求得回归应力值的最小值对应的法方程式，进而可以求解得到四个应力计算值和对应的应力回归系数。
进一步：所述建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场具体包括如下步骤:
采用2-范数定量评价地应力场反演分析后的应力计算值与应力测量值之间的误差值Δ，根据所述误差值Δ确定应力回归系数的最优值Li*，所述误差值的计算公式如下：
Figure GDA0002944230280000071
式中，||·||表示2-范数，当l＝1,2,3时，
Figure GDA0002944230280000072
分别为第k(k＝1,2，...，m)个测点、第l(l＝1,2,...,6)个正应力实测值和计算值；当l＝4,5,6时，
Figure GDA0002944230280000073
分别为第k个测点、第l个剪应力实测值和计算值；
将所述应力回归系数的最优值Li*与对应荷载作用下单位荷载的应力计算值相乘，并进行叠加，再迭代计算生成最终的最优目标地应力场，计算公式如下：
Figure GDA0002944230280000074
其中，当l＝1,2,3时，
Figure GDA0002944230280000075
为第i个荷载作用下计算域内任一测点P处的第j个正应力计算值；当l＝4,5,6时，
Figure GDA0002944230280000076
为第i个荷载作用下计算域内任一测点P处的第j个剪应力计算值，j＝1,2,…,6，σjP为计算域内任一测点P处的第j个回归应力值。
上述进一步方案的有益效果是：通过2-范数可以准确定量地评价地应力场反演分析后的应力计算值与应力测量值之间的误差，从而确定误差最小时对应的最优应力回归系数，再将应力回归系数的最优值与对应荷载作用下单位荷载的应力计算值进行线性叠加，并迭代，最终得到准确的最优目标地应力场。
进一步：所述根据所述误差值Δ确定应力回归系数的最优值Li*具体包括如下步骤：
根据所述误差值Δ判定是否满足反演精度，具体为：设定正应力权重为wq时的总应力误差Δq、正应力误差Δσq和剪应力误差Δτq；
当Δq>Δ0时，反演精度不足，本次地应力场反演计算无效；
当Δq≤Δ0时，反演精度满足工程要求；
对于满足工程要求正应力权重和剪应力权重，令Δστq＝max(Δτq,Δσq)，当且仅当正应力权重为wq时，存在Δστq*＝min(Δστq)，则此时的权重wq为正应力最优权重w*，即w*＝wq；
当且仅当正应力权重为最优权重w*时，应力回归系数为最优值，即Li*＝Li；
其中，Δσq、Δτq、Δq分别代表正应力、剪应力和总应力误差；Δ0代表工程反演精度；wq为正应力权重；w*为最优的剪应力权重；Δστq为正应力权重为w时的正应力误差Δσi、剪应力误差Δτi二者之间的最大值；min(Δστq)为Δστq的最小值；Δστq*为Δστq的最小值。
上述进一步方案的有益效果是：通过将所述地应力场反演分析后的应力计算值与应力测量值之间的总应力误差值与预设反应精度进行比较，可以确定反应精度是否满足预设工程精度要求，并在满足预设反演工程精度要求时，准确确定正应力最优权重，进而确定应力回归系数的最优值。
本发明还提供了一种正应力和剪应力不同权重三维应力场反演系统，包括：
测试与试验模块，用于通过原位应力测试和承压板载荷试验，获取不同空间部位测点的力学参数测量信息；
单位荷载作用及分析模块，用于根据工程地质信息建立有限元数值分析模型，并分别施加重力荷载、水平x方向单位挤压荷载、水平y方向单位挤压荷载及单位剪切荷载，分别确定对应荷载作用下测点处的应力计算值；
反演目标函数模型模块，用于建立正应力和剪应力不同权重的三维地应力场反演目标函数模型：
Figure GDA0002944230280000091
其中，S残'为残余误差值，
Figure GDA0002944230280000092
为第k个测点、第i个正应力实测值，
Figure GDA0002944230280000093
为第k个测点、第j个剪应力实测值，，w为正应力权重，w的范围为[0,1]；
Figure GDA0002944230280000094
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure GDA0002944230280000095
为第k个测点在对应荷载作用下的第j个剪应力计算值；
反演分析模块，用于分别采用正应力和剪应力不同权重下的三维地应力场反演目标函数对目标地应力场进行反演分析，获得不同权重下测点处的应力计算值和应力回归系数；
反演评价模块，用于建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场。
本发明的正应力和剪应力不同权重三维应力场反演系统，通过在有限元数值分析模型上施加不用方向的荷载，获取测点处的应力计算值，并进行反演分析，得到不同权重下测点处的应力计算值和应力回归系数，结合试验测得的应力矢量进行反演评价，从而准确确定最优目标地应力场，理论上解决了剪应力与正应力之间由于数值量级差异引起的剪应力在反演计算中处于次要甚至被忽略不计的问题，能够反映正应力和剪应力不等权重对地应力分布特征的影响，反演获得的应力矢量精度大大提高。通过采用本文的反演目标函数模型及反演方法对工程区三维初始应力场进行反演分析，可以准确、可靠、高效地获取地应力场分布特征，为科学识别开挖后的硐室围岩变形破坏特征，以及永久支护设计的优化提供必备的基础资料和技术支撑。
本发明还提供了一种计算机可读存储介质，其上存储有计算机程序，所述计算机程序被处理器执行时，实现所述的方法。
本发明还提供了一种正应力和剪应力不同权重三维应力场反演设备,其特征在于，包括所述的存储介质和处理器，所述处理器执行所述存储介质上的计算机程序时实现所述方法的步骤。
附图说明
图1为本发明一实施例的正应力和剪应力不同权重三维应力场反演方法流程示意图；
图2为本发明一实施例的有限差分数值模型示意图；
图3为本发明一实施例的自重荷载作用下对应的应力场形成示意图；
图4为本发明一实施例的水平x方向荷载作用下对应的应力场形成示意图；
图5为本发明一实施例的水平y方向荷载作用下对应的应力场形成示意图；
图6为本发明一实施例的剪切荷载作用下对应的应力场形成示意图；
图7a为本发明一实施例的正应力权重与总应力误差的变化关系示意图；
图7b为本发明另一实施例的正应力权重与正应力误差的变化关系示意图；
图7c为本发明又一实施例的正应力权重与剪应力误差的变化关系示意图；
图8为本发明一实施例的待开挖区域沿坝址区横剖面第一主应力等值线示意图；
图9为本发明一实施例的待开挖区域沿坝址区横剖面第二主应力等值线示意图；
图10为本发明一实施例的待开挖区域沿坝址区横剖面第三主应力等值线示意图；
图11为本发明一实施例的正应力和剪应力不同权重三维应力场反演系统结构框图。
具体实施方式
以下结合附图对本发明的原理和特征进行描述，所举实例只用于解释本发明，并非用于限定本发明的范围。
某在建的大型水电站地下厂房洞室群位于金沙江下游河段，该电站属山原峡谷地貌类型，坝址区峡谷两岸地形陡峻，河谷呈狭窄的“V”型，两岸地形不对称。金沙江以N20°E向流经坝址区，河谷狭窄，岸坡陡立，坡角一般为60°～75°，硐室水平埋深100～500m，竖向埋深180～430m。坝址区出露的地层岩性主要有因民组(Pt2y)中厚～厚层白云岩、灰质白云岩，落雪组(Pt2l1+2)中厚～厚层白云岩、灰质白云岩和薄层灰岩，落雪组(Pt2l3-1～Pt2l3-5)中厚层夹薄层灰岩、中厚层大理岩及厚层灰质白云岩，落雪组(Pt2l4-10)薄层～中厚层灰岩及灰白色薄层～厚层大理岩。坝址区地质构造总体以断层为主，规模较大的断层从上游至下游分别为F3、F14、F6及F7，且与河谷基本垂直。
如图1所示，一种正应力和剪应力不同权重三维应力场反演方法，包括如下步骤：
S1：通过原位应力测试和承压板载荷试验，获取三维空间不同空间部位测点的力学参数测量信息。
S11：在施工现场勘探平硐内进行原位地应力测试，获取三维空间空间不同部位6个测点的应力张量
Figure GDA0002944230280000111
如表1所示，计算坐标系下的应力分量见表2。
表1待反演的应力测量值
Figure GDA0002944230280000121
其中，α1～α3分别为第一、二、三主应力的水平面倾角，正为上倾，负为下倾；β1～β3分别为第一、二、三主应力的方位角，以正北顺时针转为正。
表2计算坐标系下实测点处的应力分量
应力值 测点1 测点2 测点3 测点4 测点5 测点6
σ<sub>xx</sub>/MPa 3.379 3.729 3.91 5.028 6.146 8.866
σ<sub>yy</sub>/MPa 4.077 4.141 4.327 5.959 4.519 6.746
σ<sub>zz</sub>/MPa 5.644 5.73 6.162 6.413 8.436 12.888
τ<sub>xy</sub>/MPa -0.186 -0.043 1.001 0.468 0.841 0.763
τ<sub>yz</sub>/MPa -0.125 -0.259 0.106 1.091 1.491 0.844
τ<sub>xz</sub>/MPa 0.316 0.243 1.081 -0.295 -0.711 -0.073
S12：在施工现场勘探平硐内进行承压板载荷试验，采用逐级一次循环法分级施加荷载，获取现场试验实测的压力-变形(P-W)关系曲线，利用半无限空间弹性体公式求得岩体弹性模量E，获得的岩体弹性模量表3，其余参数可通过现场工程地质勘察报告获取，岩体弹性模量计算具体公式如下：
Figure GDA0002944230280000122
其中，当l＝1,2,3时，
Figure GDA0002944230280000123
为第k个测点、第l个正应力测量值，当l＝4,5,6时，
Figure GDA0002944230280000124
为为第k个测点、第l个剪应力测量值，W为岩体弹性变形，单位为cm；P为按承压板面积计算的压力，单位为MPa；D为承压板直径单位为cm；ν为泊松比。
通过进行原位应力测试，可以准确获取不同测点的应力测量值
Figure GDA0002944230280000133
和岩体弹性模量，这样便于后续根据应力测量值以及反演得到的应力计算值评价反演精度，同时，通过所述岩体弹性模量便于对所述有限元数值分析模型进行初始化。
获得的岩体弹性模量等岩体物理力学参数如下表3所示。
表3岩体物理力学参数
Figure GDA0002944230280000131
S2：根据工程地质信息建立有限元数值分析模型，并在三维空间内不同方向上分别施加荷载，分别确定对应荷载作用下测点处的应力计算值。
S21：根据工程地质信息，建立有限元数值分析模型，并对所述有限元数值分析模型进行有限元网格划分，并输入相应的岩体弹性模量等地层岩性、地质构造等参数，初始化所述有限元数值分析模型的参数信息，如图2所示；
S22：沿所述有限元数值分析模型边界分别施加重力荷载(图3)、水平x方向单位挤压荷载(图4)、水平y方向单位挤压荷载(图5)及水平面内单位剪切荷载(图6)，通过有限元数值处理，并分别提取对应载荷作用下测点处的应力计算值，如表4-表7所示。
表4自重荷载作用下测点处的应力计算值
Figure GDA0002944230280000132
Figure GDA0002944230280000141
表5水平x方向单位荷载作用下测点处的应力计算值
应力值 测点1 测点2 测点3 测点4 测点5 测点6
σ<sub>xx</sub>/MPa -1.813 -1.777 -1.928 -2.025 -1.572 -1.631
σ<sub>yy</sub>/MPa -0.466 -0.436 -0.535 -0.567 -0.412 -0.431
σ<sub>zz</sub>/MPa -0.555 -0.466 -0.533 -0.528 -0.341 -0.301
τ<sub>xy</sub>/MPa 0.095 0.093 0.072 0.075 0.051 0.051
τ<sub>yz</sub>/MPa 0.089 0.085 0.070 0.070 0.076 0.069
τ<sub>xz</sub>/MPa -0.127 -0.126 -0.104 -0.101 -0.078 -0.069
表6水平y方向单位荷载作用下测点处的应力计算值
应力值 测点1 测点2 测点3 测点4 测点5 测点6
σ<sub>xx</sub>/MPa -0.204 -0.190 -0.258 -0.283 -0.219 -0.241
σ<sub>yy</sub>/MPa -0.723 -0.709 -0.869 -0.896 -0.982 -1.009
σ<sub>zz</sub>/MPa 0.064 0.097 0.050 0.048 0.067 0.075
τ<sub>xy</sub>/MPa 0.078 0.078 0.046 0.047 0.030 0.028
τ<sub>yz</sub>/MPa 0.053 0.051 0.074 0.071 0.075 0.068
τ<sub>xz</sub>/MPa -0.052 -0.053 -0.055 -0.054 -0.040 -0.035
表7水平面内单位剪切荷载作用下测点处的应力计算值
应力值 测点1 测点2 测点3 测点4 测点5 测点6
σ<sub>xx</sub>/MPa 0.171 0.190 0.002 0.014 0.052 0.056
σ<sub>yy</sub>/MPa -0.036 -0.026 0.054 0.054 0.040 0.028
σ<sub>zz</sub>/MPa -0.064 -0.068 0.106 0.109 0.049 0.050
τ<sub>xy</sub>/MPa -0.148 -0.148 -0.161 -0.164 -0.137 -0.140
τ<sub>yz</sub>/MPa -0.309 -0.307 -0.239 -0.228 -0.164 -0.151
τ<sub>xz</sub>/MPa 0.052 0.051 0.060 0.058 0.046 0.040
通过建立有限元数值分析模型，并进行网格划分、初始化，这样在不同荷载的作用下可以分别提取对应的应力计算值，便于后续进行反演运算。
S3：建立正应力和剪应力不同权重的三维地应力场反演目标函数模型：
Figure GDA0002944230280000151
其中，S残'为残余误差值，
Figure GDA0002944230280000152
为第k个测点、第i个正应力实测值，
Figure GDA0002944230280000153
为第k个测点、第j个剪应力实测值，，w为正应力权重，w的范围为[0,1]；
Figure GDA0002944230280000154
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure GDA0002944230280000155
为第k个测点在对应荷载作用下的第j个剪应力计算值。
S31：以测点处的所述应力测量值与应力计算值之差的平方和达到最小为优化目标，建立m个测点下的正应力与剪应力相互独立的三维反演目标函数模型：
Figure GDA0002944230280000156
其中：m为测点的总数；当l＝1,2,3时，
Figure GDA0002944230280000157
分别为第k个测点、第l个正应力测量值和计算值；当l＝4,5,6时，
Figure GDA0002944230280000158
分别为第k个测点、第l个剪应力测量值和计算值，
Figure GDA0002944230280000159
分别为第k个测点、第i个正应力测量值和计算值，
Figure GDA00029442302800001510
分别为第k个测点剪应力测量值和计算值；
S32：根据上述公式(2)建立正应力权重为w时的三维地应力场反演目标函数模型，具体如下：
假设第k个测点、第i个正应力σik权重为w，相应的第k个测点、第j个剪切应力τjk权重为1-w，代入公式(2)中，可得：
Figure GDA00029442302800001511
其中，w为0～1.0之间的任意值；
S33：根据上述公式(3)对应力测量值和计算值进行形式变换，建立一般形式下的三维地应力场反演目标函数模型，可得：
Figure GDA00029442302800001512
其中，
Figure GDA0002944230280000161
为第k个测点、第i个正应力实测值，
Figure GDA0002944230280000162
为第k个测点、第j个剪应力实测值，，w为正应力权重，w的范围为[0,1]；
Figure GDA0002944230280000163
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure GDA0002944230280000164
为第k个测点在对应荷载作用下的第j个剪应力计算值。
通过建立测点的正应力与剪应力相互独立的三维反演目标函数模型，并结合应力测量值和计算值进行形式变换，得到一般形式下的三维地应力场反演目标函数模型，便于后续针对应力测量值和应力计算值进行反演分析以及反演评价，得到准确的最优目标地应力场。
假设w＝0.1，对表4-表7中的单位荷载作用下测点处的应力计算值进行权重变换，见表8-11所示。
表8自重荷载作用下测点处的应力计算值
应力值 测点1 测点2 测点3 测点4 测点5 测点6
σ<sub>xx</sub>/MPa -4.278 -4.405 -5.157 -5.416 -5.24 -5.822
σ<sub>yy</sub>/MPa -2.066 -2.235 -2.554 -2.641 -3.407 -3.751
σ<sub>zz</sub>/MPa -5.982 -6.227 -8.318 -8.524 -10.931 -11.941
τ<sub>xy</sub>/MPa 0.345 0.333 0.387 0.375 0.405 0.387
τ<sub>yz</sub>/MPa 0.192 0.174 0.783 0.78 1.107 1.128
τ<sub>xz</sub>/MPa -0.474 -0.495 -0.459 -0.459 -0.408 -0.411
表9水平x方向单位荷载作用下测点处的应力计算值
应力值 测点1 测点2 测点3 测点4 测点5 测点6
σ<sub>xx</sub>/MPa -1.813 -1.777 -1.928 -2.025 -1.572 -1.631
σ<sub>yy</sub>/MPa -0.466 -0.436 -0.535 -0.567 -0.412 -0.431
σ<sub>zz</sub>/MPa -0.555 -0.466 -0.533 -0.528 -0.341 -0.301
τ<sub>xy</sub>/MPa 0.285 0.279 0.216 0.225 0.153 0.153
τ<sub>yz</sub>/MPa 0.267 0.255 0.21 0.21 0.228 0.207
τ<sub>xz</sub>/MPa -0.381 -0.378 -0.312 -0.303 -0.234 -0.207
表10水平y方向单位荷载作用下测点处的应力计算值
Figure GDA0002944230280000165
Figure GDA0002944230280000171
表11水平面内单位剪切荷载作用下实测点处的应力计算值
应力值 测点1 测点2 测点3 测点4 测点5 测点6
σ<sub>xx</sub>/MPa 0.171 0.19 0.002 0.014 0.052 0.056
σ<sub>yy</sub>/MPa -0.036 -0.026 0.054 0.054 0.04 0.028
σ<sub>zz</sub>/MPa -0.064 -0.068 0.106 0.109 0.049 0.05
τ<sub>xy</sub>/MPa -0.444 -0.444 -0.483 -0.492 -0.411 -0.42
τ<sub>yz</sub>/MPa -0.927 -0.921 -0.717 -0.684 -0.492 -0.453
τ<sub>xz</sub>/MPa 0.156 0.153 0.18 0.174 0.138 0.12
S4：分别采用正应力和剪应力不同权重下的三维地应力场反演目标函数对目标地应力场进行反演分析，获得不同权重下测点处的应力计算值和应力回归系数。
S41：根据对应载荷下不同权重正应力的应力计算值和应力回归系数计算回归应力值，计算公式如下：
Figure GDA0002944230280000172
其中，当l＝1,2,3时，
Figure GDA0002944230280000173
为第k个测点、第j个正应力回归计算值；当l＝4,5,6时，
Figure GDA0002944230280000174
为第k个测点、第j个剪应力回归计算值，Li为相应于自变量的多元回归系数，当l＝1,2,3时，
Figure GDA0002944230280000175
为第k个测点、第j个正应力计算值的单列矩阵；当l＝4,5,6时，
Figure GDA0002944230280000176
为第k个测点、第j个剪应力计算值的单列矩阵；
S42：根据多元回归分析理论，将上式(5)代入式(4)的公式中，可得：
Figure GDA0002944230280000181
根据最小二乘法原理，S残’为最小值的法方程式为：
Figure GDA0002944230280000182
S43：：采用线性回归分析法对上述法方程式(7)进行求解，得到正应力权重为w时的三个应力计算值和对应的4个应力回归系数L＝(L1,L2,...,Ln)T；
当w＝0.1时，计算得到的应力回归系数L＝(0.928,0.077,2.829,2.606)T。
同理，按照上述步骤，可得到正应力权重w＝0.2～0.9时的回归应力系数，见表12。
表12不同正应力权重下的应力回归系数
应力回归系数 w＝0.1 w＝0.2 w＝0.3 w＝0.4 w＝0.5 w＝0.6 w＝0.7 w＝0.8 w＝0.9
L<sub>1</sub> 0.928 0.925 0.930 0.931 0.933 0.934 0.936 0.938 0.941
L<sub>2</sub> 0.077 0.090 0.067 0.058 0.047 0.035 0.021 0.002 -0.020
L<sub>3</sub> 2.829 2.753 2.855 2.869 2.876 2.880 2.883 2.884 2.884
L<sub>4</sub> 2.606 2.235 2.938 3.276 3.641 4.059 4.535 5.136 5.837
S44：根据所述应力回归系数计算计算域内任一点P的回归应力值，具体公式如下：
Figure GDA0002944230280000183
其中，当l＝1,2,3时，σjp’为测点P处、第j个正应力计算值；当l＝4,5,6时，σjp’为测点P处、第j个剪应力计算值，计算域内任一点P的回归应力值如表10所示：
表13正应力不同权重下的应力回归值
Figure GDA0002944230280000184
Figure GDA0002944230280000191
Figure GDA0002944230280000201
通过将对应载荷下不同权重正应力的应力计算值和应力回归系数进行线性叠加，可以得到回归应力值，并结合最小二乘法求得回归应力值的最小值对应的法方程式，进而可以求解得到四个应力计算值和对应的应力回归系数。
S5：建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场。
S51：采用2-范数定量评价地应力场反演分析后的应力计算值与应力测量值之间的误差值Δ，根据所述误差值Δ确定应力回归系数的最优值Li*，所述误差值的计算公式如下：
Figure GDA0002944230280000211
式中，||·||表示2-范数，当l＝1,2,3时，
Figure GDA0002944230280000212
分别为第k(k＝1,2，...，m)个测点、第l(l＝1,2,...,6)个正应力实测值和计算值；当l＝4,5,6时，
Figure GDA0002944230280000213
分别为第k个测点、第l个剪应力实测值和计算值；
根据上述公式，对正应力和剪应力不等权重组合下的反演计算结果进行误差分析，获得正应力权重为wq时的总应力误差Δq、正应力误差Δσq和剪应力误差Δτq，见表14：
表14反演计算结果误差
正应力权重 正应力误差Δ<sub>σq</sub> 剪应力误差Δ<sub>τq</sub> 总应力误差Δ<sub>q</sub>
w＝0.1 0.06％ 15.56％ 2.59％
w＝0.2 0.20％ 4.05％ 2.89％
w＝0.3 0.16％ 20.6％ 2.48％
w＝0.4 0.22％ 23.45％ 2.41％
w＝0.5 0.26％ 25.33％ 2.37％
w＝0.6 0.29％ 27.00％ 2.33％
w＝0.7 0.32％ 29.02％ 2.30％
w＝0.8 0.36％ 32.34％ 2.28％
w＝0.9 0.40％ 39.66％ 2.27％
为了综合评价反演结果的实用性和可靠性，本发明采用了下述最优正应力和剪应力组合权重的判定准则。
可选地，在本发明的一个或多个实施例中，所述根据所述误差值Δ确定应力回归系数的最优值Li*具体包括如下步骤：
S51a：根据所述误差值Δ判定是否满足反演精度，具体为：设定正应力权重为wq时的总应力误差Δq、正应力误差Δσq和剪应力误差Δτq；
当Δq>Δ0时，反演精度不足，本次地应力场反演计算无效；
当Δq≤Δ0时，反演精度满足工程要求；
S51b：对于满足工程要求正应力权重和剪应力权重，令Δστq＝max(Δτq,Δσq)，当且仅当正应力权重为wq时，存在Δστq*＝min(Δστq)，则此时的权重wq为正应力最优权重w*，即w*＝wq；
当且仅当正应力权重为最优权重w*时，应力回归系数为最优值，即Li*＝Li；
其中，Δσq、Δτq、Δq分别代表正应力、剪应力和总应力误差；Δ0代表工程反演精度；wq为正应力权重；w*为最优的剪应力权重；Δστq为正应力权重为w时的正应力误差Δσi、剪应力误差Δτi二者之间的最大值；min(Δστq)为Δστq的最小值；Δστq*为Δστq的最小值。
通过将所述地应力场反演分析后的应力计算值与应力测量值之间的总应力误差值与预设反应精度进行比较，可以确定反应精度是否满足预设工程精度要求，并在满足预设反演工程精度要求时，准确确定正应力最优权重，进而确定应力回归系数的最优值。
根据工程经验，总应力误差精度一般取Δ0＝5％～10％，对于重要工程应取小值。根据表14，不同权重下的Δq均小于Δ0，表明反演精度均满足工程要求，但仍需进一步确定最优的权重及回归应力系数。令Δστq＝max(Δτq,Δσq)，则w＝0.1～0.9时的相对误差Δστq分别为15.56％、4.05％、20.6％、23.45％、25.33％、27.00％、29.02％、32.34％、39.66％。可见，当w＝0.2时，Δστq*＝min(Δστq)＝4.05％，则最优权重w*＝0.2，相应的最优回归应力系数L*＝(0.925,0.090,2.753,2.235)T。不同权重下的应力场反演获得的总应力误差、正应力误差、剪应力误差见图7。
S52：将所述应力回归系数的最优值Li*与对应荷载作用下单位荷载的应力计算值相乘，并进行叠加，再迭代计算生成最终的最优目标地应力场，计算公式如下：
Figure GDA0002944230280000221
其中，当l＝1,2,3时，
Figure GDA0002944230280000231
为第i个荷载作用下计算域内任一测点P处的第j个正应力计算值；当l＝4,5,6时，
Figure GDA0002944230280000232
为第i个荷载作用下计算域内任一测点P处的第j个剪应力计算值，j＝1,2,…,6，σjP为计算域内任一测点P处的第j个回归应力值。
最优权重时反演获得的工程区应力场等值线分布见图8-图10。
本发明的正应力和剪应力不同权重三维应力场反演方法，通过在有限元数值分析模型上施加不用方向的荷载，获取测点处的应力计算值，并进行反演分析，得到不同权重下测点处的应力计算值和应力回归系数，结合试验测得的应力矢量进行反演评价，从而准确确定最优目标地应力场，理论上解决了剪应力与正应力之间由于数值量级差异引起的剪应力在反演计算中处于次要甚至被忽略不计的问题，能够反映正应力和剪应力不等权重对地应力分布特征的影响，反演获得的应力矢量精度大大提高。通过采用本文的反演目标函数模型及反演方法对工程区三维初始应力场进行反演分析，可以准确、可靠、高效地获取地应力场分布特征，为科学识别开挖后的硐室围岩变形破坏特征，以及永久支护设计的优化提供必备的基础资料和技术支撑。
如图11所示，一种正应力和剪应力不同权重三维应力场反演系统，包括：
测试与试验模块，用于通过原位应力测试和承压板载荷试验，获取不同空间部位测点的力学参数测量信息；
单位荷载作用及分析模块，用于根据工程地质信息建立有限元数值分析模型，并分别施加重力荷载、水平x方向单位挤压荷载、水平y方向单位挤压荷载及单位剪切荷载，分别确定对应荷载作用下测点处的应力计算值；
反演目标函数模型模块，用于建立正应力和剪应力不同权重的三维地应力场反演目标函数模型：
Figure GDA0002944230280000233
其中，S残'为残余误差值，
Figure GDA0002944230280000241
为第k个测点、第i个正应力实测值，
Figure GDA0002944230280000242
为第k个测点、第j个剪应力实测值，，w为正应力权重，w的范围为[0,1]；
Figure GDA0002944230280000243
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure GDA0002944230280000244
为第k个测点在对应荷载作用下的第j个剪应力计算值；
反演分析模块，用于分别采用正应力和剪应力不同权重下的三维地应力场反演目标函数对目标地应力场进行反演分析，获得不同权重下测点处的应力计算值和应力回归系数；
反演评价模块，用于建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场。
本发明的正应力和剪应力不同权重三维应力场反演系统，通过在有限元数值分析模型上施加不用方向的荷载，获取测点处的应力计算值，并进行反演分析，得到不同权重下测点处的应力计算值和应力回归系数，结合试验测得的应力矢量进行反演评价，从而准确确定最优目标地应力场，理论上解决了剪应力与正应力之间由于数值量级差异引起的剪应力在反演计算中处于次要甚至被忽略不计的问题，能够反映正应力和剪应力不等权重对地应力分布特征的影响，反演获得的应力矢量精度大大提高。通过采用本文的反演目标函数模型及反演方法对工程区三维初始应力场进行反演分析，可以准确、可靠、高效地获取地应力场分布特征，为科学识别开挖后的硐室围岩变形破坏特征，以及永久支护设计的优化提供必备的基础资料和技术支撑。
本发明还提供了一种计算机可读存储介质，其上存储有计算机程序，所述计算机程序被处理器执行时，实现所述的方法。
本发明还提供了一种正应力和剪应力不同权重三维应力场反演设备,其特征在于，包括所述的存储介质和处理器，所述处理器执行所述存储介质上的计算机程序时实现所述方法的步骤。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种正应力和剪应力不同权重三维应力场反演方法,其特征在于，包括如下步骤：
通过原位应力测试和承压板载荷试验，获取三维空间不同空间部位测点的力学参数测量信息；
根据工程地质信息建立有限元数值分析模型，并在三维空间内不同方向上分别施加荷载，分别确定对应荷载作用下测点处的应力计算值；
建立正应力和剪应力不同权重的三维地应力场反演目标函数模型：
Figure FDA0002944230270000011
其中，S残'为残余误差值，m为测点的总数,
Figure FDA0002944230270000012
Figure FDA0002944230270000013
为第k个测点、第i个正应力实测值，
Figure FDA0002944230270000014
为第k个测点、第j个剪应力实测值，w为正应力权重，w的范围为[0,1]；
Figure FDA0002944230270000015
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure FDA0002944230270000016
为第k个测点在对应荷载作用下的第j个剪应力计算值；
分别采用正应力和剪应力不同权重下的三维地应力场反演目标函数对目标地应力场进行反演分析，获得不同权重下测点处的应力计算值和应力回归系数；
建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场。
2.根据权利要求1所述正应力和剪应力不同权重三维应力场反演方法,其特征在于，所述通过原位应力测试和承压板载荷试验，获取三维空间不同空间部位测点的力学参数测量信息包括如下步骤：
所述通过原位应力测试和承压板载荷试验，获取不同空间部位测点的应力测量参数信息具体包括如下步骤：
在施工现场勘探平硐内进行原位应力测试，获取分析断面处的不同部位m个测点的应力测量值
Figure FDA0002944230270000021
(k＝1,2，...，m；l＝1,2,...，6)和应力方向；
在施工现场勘探平硐内进行承压板载荷试验，采用逐级一次循环法分级施加荷载，获取现场试验实测的压力-变形(P-W)关系曲线，利用半无限空间弹性体公式求得岩体弹性模量E，具体公式如下：
Figure FDA0002944230270000022
其中，当l＝1,2,3时，
Figure FDA0002944230270000023
为第k个测点、第l个正应力测量值，当l＝4,5,6时，
Figure FDA0002944230270000024
为为第k个测点、第l个剪应力测量值，W为岩体弹性变形，单位为cm；P为按承压板面积计算的压力，单位为MPa；D为承压板直径单位为cm；ν为泊松比。
3.根据权利要求2所述正应力和剪应力不同权重三维应力场反演方法,其特征在于，所述根据工程地质信息建立有限元数值分析模型，并在三维空间内不同方向上分别施加荷载，分别确定对应荷载作用下测点处的应力计算值具体包括如下步骤：
根据工程地质信息，建立有限元数值分析模型，并对所述有限元数值分析模型进行有限元网格划分，并根据所述岩体弹性模量E初始化所述有限元数值分析模型的参数信息；
沿所述有限元数值分析模型边界分别施加重力荷载、水平x方向单位挤压荷载、水平y方向单位挤压荷载及水平面内单位剪切荷载，通过有限元数值处理，并分别提取对应载荷作用下测点处的应力计算值。
4.根据权利要求3所述正应力和剪应力不同权重三维应力场反演方法,其特征在于，所述建立正应力和剪应力不同权重的地应力场三维反演目标函数模型具体包括如下步骤：
以测点处的所述应力测量值与应力计算值之差的平方和达到最小为优化目标，建立m个测点下的正应力与剪应力相互独立的三维反演目标函数模型：
Figure FDA0002944230270000031
其中：m为测点的总数；当l＝1,2,3时，
Figure FDA0002944230270000032
分别为第k个测点、第l个正应力测量值和计算值；当l＝4,5,6时，
Figure FDA0002944230270000033
分别为第k个测点、第l个剪应力测量值和计算值，
Figure FDA0002944230270000034
分别为第k个测点、第i个正应力测量值和计算值，
Figure FDA0002944230270000035
分别为第k个测点剪应力测量值和计算值；
根据上述公式(2)建立正应力权重为w时的三维地应力场反演目标函数模型，具体如下：
假设第k个测点、第i个正应力σik权重为w，相应的第k个测点、第j个剪切应力τjk权重为1-w，代入公式(2)中，可得：
Figure FDA0002944230270000036
其中，w为0～1.0之间的任意值；
根据上述公式(3)对应力测量值和计算值进行形式变换，建立一般形式下的三维地应力场反演目标函数模型，可得：
Figure FDA0002944230270000037
其中，
Figure FDA0002944230270000038
Figure FDA0002944230270000039
为第k个测点、第i个正应力实测值，
Figure FDA00029442302700000310
为第k个测点、第j个剪应力实测值，w为正应力权重，w的范围为[0,1]；
Figure FDA00029442302700000311
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure FDA00029442302700000312
为第k个测点在对应荷载作用下的第j个剪应力计算值。
5.根据权利要求4所述的正应力和剪应力不同权重三维应力场反演方法,其特征在于，所述获得不同权重下测点处的应力计算值和应力回归系数具体包括如下步骤：
根据对应载荷下不同权重正应力的应力计算值和应力回归系数计算回归应力值，计算公式如下：
Figure FDA0002944230270000041
其中，当l＝1,2,3时，
Figure FDA0002944230270000042
为第k个测点、第j个正应力回归计算值；当l＝4,5,6时，
Figure FDA0002944230270000043
为第k个测点、第j个剪应力回归计算值，Li为相应于自变量的多元回归系数，当l＝1,2,3时，
Figure FDA0002944230270000044
为第k个测点、第j个正应力计算值的单列矩阵；当l＝4,5,6时，
Figure FDA0002944230270000045
为第k个测点、第j个剪应力计算值的单列矩阵；
根据多元回归分析理论，将上式(5)代入式(4)的公式中，可得：
Figure FDA0002944230270000046
根据最小二乘法原理，S残’为最小值的法方程式为：
Figure FDA0002944230270000047
采用线性回归分析法对上述法方程式(7)进行求解，得到正应力权重为w时的三个应力计算值和对应的4个应力回归系数L＝(L1,L2,...,Ln)T；
根据所述应力回归系数计算计算域内任一点P的回归应力值，具体公式如下：
Figure FDA0002944230270000048
其中，当l＝1,2,3时，σjp’为测点P处、第j个正应力计算值；当l＝4,5,6时，σjp’为测点P处、第j个剪应力计算值，j＝1,2,…,6对应初始应力六个分量。
6.根据权利要求5所述的正应力和剪应力不同权重三维应力场反演方法,其特征在于，所述建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场具体包括如下步骤:
采用2-范数定量评价地应力场反演分析后的应力计算值与应力测量值之间的误差值Δ，根据所述误差值Δ确定应力回归系数的最优值Li*，所述误差值的计算公式如下：
Figure FDA0002944230270000051
式中，||·||表示2-范数，当l＝1,2,3时，
Figure FDA0002944230270000052
分别为第k(k＝1,2，...，m)个测点、第l(l＝1,2,...,6)个正应力实测值和计算值；当l＝4,5,6时，
Figure FDA0002944230270000053
分别为第k个测点、第l个剪应力实测值和计算值；
将所述应力回归系数的最优值Li*与对应荷载作用下单位荷载的应力计算值相乘，并进行叠加，再迭代计算生成最终的最优目标地应力场，计算公式如下：
Figure FDA0002944230270000054
其中，当l＝1,2,3时，
Figure FDA0002944230270000055
为第i个荷载作用下计算域内任一测点P处的第j个正应力计算值；当l＝4,5,6时，
Figure FDA0002944230270000056
为第i个荷载作用下计算域内任一测点P处的第j个剪应力计算值，j＝1,2,…,6，σjP为计算域内任一测点P处的第j个回归应力值。
7.根据权利要求6所述的正应力和剪应力不同权重三维应力场反演方法,其特征在于：所述根据所述误差值Δ确定应力回归系数的最优值Li*具体包括如下步骤：
根据所述误差值Δ判定是否满足反演精度，具体为：设定正应力权重为wq时的总应力误差Δq、正应力误差Δσq和剪应力误差Δτq；
当Δq>Δ0时，反演精度不足，本次地应力场反演计算无效；
当Δq≤Δ0时，反演精度满足工程要求；
对于满足工程要求正应力权重和剪应力权重，令Δστq＝max(Δτq,Δσq)，当且仅当正应力权重为wq时，存在Δστq*＝min(Δστq)，则此时的权重wq为正应力最优权重w*，即w*＝wq；
当且仅当正应力权重为最优权重w*时，应力回归系数为最优值，即Li*＝Li；
其中，Δσq、Δτq、Δq分别代表正应力、剪应力和总应力误差；Δ0代表工程反演精度；wq为正应力权重；w*为最优的剪应力权重；Δστq为正应力权重为w时的正应力误差Δσi、剪应力误差Δτi二者之间的最大值；min(Δστq)为Δστq的最小值；Δστq*为Δστq的最小值。
8.一种正应力和剪应力不同权重三维应力场反演系统，其特征在于，包括：
测试与试验模块，用于通过原位应力测试和承压板载荷试验，获取三维空间不同空间部位测点的力学参数测量信息；
单位荷载作用及分析模块，用于根据工程地质信息建立有限元数值分析模型，并分别施加重力荷载、水平x方向单位挤压荷载、水平y方向单位挤压荷载及单位剪切荷载，分别确定对应荷载作用下测点处的应力计算值；
反演目标函数模型模块，用于建立正应力和剪应力不同权重的三维地应力场反演目标函数模型：
Figure FDA0002944230270000061
其中，S残'为残余误差值，
Figure FDA0002944230270000062
Figure FDA0002944230270000063
为第k个测点、第i个正应力实测值，
Figure FDA0002944230270000064
为第k个测点、第j个剪应力实测值，w为正应力权重，w的范围为[0,1]；
Figure FDA0002944230270000065
为第k个测点在对应荷载作用下的第i个正应力计算值，
Figure FDA0002944230270000066
为第k个测点在对应荷载作用下的第j个剪应力计算值；
反演分析模块，用于分别采用正应力和剪应力不同权重下的三维地应力场反演目标函数对目标地应力场进行反演分析，获得不同权重下测点处的应力计算值和应力回归系数；
反演评价模块，用于建立地应力场反演效果评价的误差计算公式，并根据所述应力参数测量信息以及不同权重下测点处的所述应力计算值和应力回归系数确定最优目标地应力场。
9.一种计算机可读存储介质，其特征在于，其上存储有计算机程序，所述计算机程序被处理器执行时，实现权利要求1至7任一项所述的方法。
10.一种正应力和剪应力不同权重三维应力场反演设备,其特征在于，包括权利要求9所述的存储介质和处理器，所述处理器执行所述存储介质上的计算机程序时实现如权利要求1至7任一项所述方法的步骤。
CN202010775092.6A 2020-08-05 2020-08-05 一种正应力和剪应力不同权重三维应力场反演方法及系统 Active CN111881613B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010775092.6A CN111881613B (zh) 2020-08-05 2020-08-05 一种正应力和剪应力不同权重三维应力场反演方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010775092.6A CN111881613B (zh) 2020-08-05 2020-08-05 一种正应力和剪应力不同权重三维应力场反演方法及系统

Publications (2)

Publication Number Publication Date
CN111881613A CN111881613A (zh) 2020-11-03
CN111881613B true CN111881613B (zh) 2021-04-06

Family

ID=73211633

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010775092.6A Active CN111881613B (zh) 2020-08-05 2020-08-05 一种正应力和剪应力不同权重三维应力场反演方法及系统

Country Status (1)

Country Link
CN (1) CN111881613B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112906162B (zh) * 2021-03-19 2024-03-29 大连科迈尔海洋科技有限公司 一种船用起重机的支撑柱的弯矩计算方法及装置
CN116796606B (zh) * 2023-06-26 2024-01-30 武汉市政工程设计研究院有限责任公司 一种考虑地应力分区特征的三维应力场反演方法及系统

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102659322A (zh) * 2012-04-06 2012-09-12 南昌航空大学 一种玻璃/pmma微纳界面结构层合材料制备方法
CN109543258A (zh) * 2018-11-09 2019-03-29 东南大学 一种基于模态叠加法的风电机组塔架应力状态计算方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20160121045A (ko) * 2015-04-09 2016-10-19 주식회사 엘지화학 중공 성형용 고밀도 폴리에틸렌 공중합체
CN105354370B (zh) * 2015-10-21 2017-12-19 河海大学 一种多层铁路路基结构安定性分析的三维有限元计算方法
CN106649925B (zh) * 2016-09-19 2019-10-18 华南理工大学 基于细宏观动态复杂应力监测的混凝土疲劳损伤分析方法
CN107700336A (zh) * 2017-10-23 2018-02-16 沈阳建筑大学 一种混凝土斜拉桥主梁施工阶段剪力滞的确定方法
CN109764986B (zh) * 2019-01-08 2020-11-27 哈尔滨工业大学（深圳） 一种基于超声横波相位谱的钢构件平面应力检测方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102659322A (zh) * 2012-04-06 2012-09-12 南昌航空大学 一种玻璃/pmma微纳界面结构层合材料制备方法
CN109543258A (zh) * 2018-11-09 2019-03-29 东南大学 一种基于模态叠加法的风电机组塔架应力状态计算方法

Also Published As

Publication number Publication date
CN111881613A (zh) 2020-11-03

Similar Documents

Publication Publication Date Title
CN111881612B (zh) 一种正应力和剪应力不同权重二维应力场反演方法及系统
Tovar-Valencia et al. Effect of surface roughness on the shaft resistance of displacement model piles in sand
Yuan et al. Investigation of deflection of a laterally loaded pile and soil deformation using the PIV technique
Guo et al. Simulation of three-dimensional tectonic stress fields and quantitative prediction of tectonic fracture within the Damintun Depression, Liaohe Basin, northeast China
Sun et al. Bayesian updating for progressive excavation of high rock slopes using multi-type monitoring data
Jeong et al. Proposed nonlinear 3-D analytical method for piled raft foundations
CN111881613B (zh) 一种正应力和剪应力不同权重三维应力场反演方法及系统
Wang et al. Bayesian stochastic soil modeling framework using Gaussian Markov random fields
He et al. Numerical modelling of transverse cracking in embankment dams
Wang et al. Probabilistic seismic displacement hazard assessment of earth slopes incorporating spatially random soil parameters
CN103605900A (zh) 跨尺度复杂地质体地应力场识别方法及装置
Huang et al. Time-dependent fragility functions for circular tunnels in soft soils
Wang et al. Lateral dynamic load tests of offshore piles based using the m-method
Mu et al. Defining the soil parameters for computing deformations caused by braced excavation.
Liu et al. Multipoint hybrid model for RCC arch dam displacement health monitoring considering construction interface and its seepage
Wang et al. Permanent accumulated rotation of an offshore monopile wind turbine in sand during a storm
CN114662360A (zh) 一种利用多软件联合构建复杂地质三维地应力场反演方法
Taşçi Analysis of dam deformation measurements with the robust and non-robust methods
Su et al. Nonprobabilistic reliability evaluation for in-service gravity dam undergoing structural reinforcement
Chen et al. Probabilistic back analysis for geotechnical engineering based on Bayesian and support vector machine
Yue et al. Soil deposit stochastic settlement simulation using an improved autocorrelation model
Platt et al. The mechanics of continental transforms: An alternative approach with applications to the San Andreas system and the tectonics of California
Hashash et al. The integration of numerical modeling and physical measurements through inverse analysis in geotechnical engineering
Konkol et al. Numerical modeling of cone penetration test in slightly overconsolidated clay with arbitrary Lagrangian-Eulerian formulation
CN116796606B (zh) 一种考虑地应力分区特征的三维应力场反演方法及系统

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