CN116110520B - 基于近场动力学岩石类材料本构模型的建立与仿真方法 - Google Patents

基于近场动力学岩石类材料本构模型的建立与仿真方法 Download PDF

Info

Publication number
CN116110520B
CN116110520B CN202310059094.9A CN202310059094A CN116110520B CN 116110520 B CN116110520 B CN 116110520B CN 202310059094 A CN202310059094 A CN 202310059094A CN 116110520 B CN116110520 B CN 116110520B
Authority
CN
China
Prior art keywords
point
force
elongation
key
function
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
CN202310059094.9A
Other languages
English (en)
Other versions
CN116110520A (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.)
Shenyang University of Technology
Original Assignee
Shenyang University of Technology
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 Shenyang University of Technology filed Critical Shenyang University of Technology
Priority to CN202310059094.9A priority Critical patent/CN116110520B/zh
Publication of CN116110520A publication Critical patent/CN116110520A/zh
Application granted granted Critical
Publication of CN116110520B publication Critical patent/CN116110520B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16C—COMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C60/00—Computational materials science, i.e. ICT specially adapted for investigating the physical or chemical properties of materials or phenomena associated with their design, synthesis, processing, characterisation or utilisation
    • 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
    • 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
    • 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
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Computing Systems (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Investigating Strength Of Materials By Application Of Mechanical Stress (AREA)

Abstract

本发明涉及一种基于近场动力学岩石类材料本构模型的建立及仿真方法，建立方法的步骤为：根据岩石类材料变形特征，确定本构模型中点对力函数的表达形式；确定点对力函数的系数；引入键的损伤函数，将点对力函数表达成和损伤函数相关的形式；基于点对力函数，推导出键的微势能；基于键的微势能，推导出物质点的应变能密度；令物质点的应变能密度等于连续介质力学的应变能密度，求出微模量的表达式；将求得的微模量代入点对力函数，完成岩石类材料本构模型的构建。本发明提出的本构模型能反映岩石类材料应力随着应变先增大再减小最后破坏这一特性，解决现有方法不能有效模拟岩石类材料在不同荷载作用下裂纹起裂、扩展、贯穿全过程的问题。

Description

基于近场动力学岩石类材料本构模型的建立与仿真方法
技术领域
本发明涉及岩土工程领域，具体涉及一种基于近场动力学岩石类材料本构模型的建立与仿真方法。
背景技术
岩石作为一种天然介质体，其内部含有大量方向各异的微裂纹，这些微裂纹在不同荷载条件下萌生、扩展、相互作用，继而贯通最终导致岩体的失稳破坏。因此，研究裂纹的萌发、扩展、贯通机理对岩土工程建设具有重要意义。
近场动力学(peridynamic)是一种新兴的基于非局部作用无网格方法，它以积分形式重构了固体力学的运动方程，解决了传统方法在求解不连续处时没有定义的问题，无需外部定义的准则，就可以模拟裂纹自发扩展的过程，广泛的应用于材料断裂行为的研究。
近场动力学可分为“键”基、常规“态”基和非常规“态”基三种。其中键基近场动力学的微观弹脆性模型是应用最广泛的一种，该模型认为物质点是用键连接的，物质点之间的相互作用力用点对力函数来表示，点对力随着键的伸长率线性变化，当键的伸长超过临界伸长率，物质点之间的键断裂，物质点之间的相互作用力消失，这不适合模拟岩石类材料变形破坏特征。一些学者通过采用指数函数或对数函数的形式来表征岩石类材料渐进损伤的过程，而采用指数函数或对数函数的形式不容易得出一般的解析解。
发明内容
针对现有技术的不足，本发明提供一种基于近场动力学岩石类材料本构模型的建立与仿真方法，其目的是解决现有方法不能有效模拟岩石类材料在不同荷载作用下裂纹起裂、扩展、贯穿全过程的问题。
为了实现上述目的，本发明采用如下技术方案：
一种近场动力学岩石类材料本构模型的构建方法，包括以下步骤：
S11根据岩石类材料变形特征，确定本构模型中点对力函数的表达形式；
S12确定点对力函数的系数；
S13引入键的损伤函数，将点对力函数表达成和损伤函数相关的形式；
S14基于步骤S11中的点对力函数，推导出键的微势能；
S15基于步骤S14中键的微势能，推导出物质点的应变能密度；
S16令物质点的应变能密度等于连续介质力学的应变能密度，求出微模量c的表达式；
S17将步骤S16中求得的微模量代入步骤S13中的点对力函数，完成岩石类材料本构模型的构建。
进一步的，步骤S11中所述的点对力函数的表达形式为：
其中，sc表示压缩临界伸长率；s2c表示压缩强化临界伸长率，此时压缩近场力最大；s1c表示压缩弹性伸长率；st表示拉伸临界伸长率；s2t表示拉伸强化临界伸长率，此时拉伸近场力最大；s1t表示拉伸弹性伸长率。
进一步的，在步骤S11所述的本构力函数中，当s1c≤s≤s1t时，键处于线弹性变形阶段，点对力与伸长率的关系为线性关系，键无损伤；当s1t＜s＜st或s1c＜s＜sc时，键处于非线性变形阶段。对于非线性变形阶段又可分为非线性强化变形阶段和非线性软化变形阶段，当s1t＜s＜s2t或s1c＜s＜s2c时，键处于非线性强化变形阶段，键开始损伤，微裂纹开始产生，点对力随着伸长率增加非线性增加；当s2t＜s＜st或s2c＜s＜sc时，键处于非线性软化阶段，点对力随着伸长率增加非线性减小；当s＞st或s＜sc时，键的损伤为1，键发生断裂，物质点之间不存在相互作用。
进一步的，所述步骤12包括以下步骤：
S12.1以键的压缩为例，设键非线性变形阶段点对力标量函数为：
f＝A′s2+B′s+C′ (2)
S12.2确定步骤S12.1所述的点对力标量函数所过点O、P的坐标；
S12.3将步骤S12.2所述的点O，P两点坐标代入步骤S12.1所述的点对力标量函数中，得到点对力和键伸长率的关系式；
S12.4根据二次函数的性质，确定A′和B′的关系式；
S12.5根据步骤S12.3所述的点对力和键伸长率的关系式、步骤S12.4所述的A′和B′的关系式，求出系数A′、B′和C′的表达式。
进一步的，所述步骤S12.2中O，P的坐标分别为(s1c，cs1c)，(sc，0)。
进一步的，所述步骤S12.3中点对力和键伸长率的关系式为：
进一步的，步骤S12.5中所述的A′、B′、C′的表达式为：
其中，sc表示压缩临界伸长率；s2c表示压缩强化临界伸长率，此时压缩近场力最大；s1c表示压缩弹性伸长率。
同理可得
其中，st表示拉伸临界伸长率；s2t表示拉伸强化临界伸长率，此时拉伸近场力最大；s1t表示拉伸弹性伸长率。
进一步的，考虑键的损伤，步骤S11所述的点对力函数可以表示为：
其中，α(s)是表征键渐进损伤的标量函数，可表示成
其中，
进一步的，步骤S14所述的微势能通过推导可得
对于二维问题，步骤S15所述的应变能密度通过推导可得
对于三维问题，步骤S15所述的应变能密度通过推导可得
进一步的，对于平面应力问题，步骤S16所述的键的微模量c通过推导可以表示为：
对于三维问题，步骤S16所述的键的微模量c通过推导可以表示为：
进一步的，对于平面应力情况下步骤S17所述的点对力函数可表示成：
其中， h为厚度；E为弹性模量。
对于三维问题步骤S17所述的点对力函数可表示成：
本发明还提供一种基于近场动力学岩石类材料裂纹扩展仿真方法，包括以下步骤：
A1初始化求解域，输入材料参数，近场尺寸；
A2离散求解域，生成物质点坐标；
A3确定每个物质点领域内的其它物质点；
A4初始化所有物质点在其近场领域内的键；
A5判断是否需要预制裂纹；
A6施加边界条件；
A7如果采用自适应动态松弛法，则时间步长Δt＝1，如果不采用自适应动态松弛法，则确定时间步长；
A8计算任意物质点总的近场力；
A9如果伸长率大于临界伸长率，则物质点间的键断裂。
A10判断是否采用动态松弛算法，如果采用动态松弛算法，则利用动态松弛算法进行时域积分，如果不采用动态松弛法利用显式向前和向后差分公式进行时域积分，求出位移和速度。
所述采用动态松弛算法的差分格式可表示为：
其中，Δt表示时间步长，Λ表示虚拟对角密度矩阵；c为阻尼系数；L表示PD相互作用力的体积力密度，b表示外力的体积力密度。
所述显式向前和向后差分格式可表示为
A11输出计算结果；
A12绘制位移，应变能密度和损伤云图。
与现有技术相比，本发明的有益效果是：
(1)本发明提出的岩石类材料本构模型能反映岩石类材料应力随着应变先增大再减小最后破坏这一特性，弥补了传统键基近场动力学不能反映岩石类材料应力-应变变化这一特性的缺陷。
(2)采用二次函数的形式来描述键的非线性变形阶段，相比于指数形式和对数形式更容易求得解析解。
(3)对所提出本构模型中的应变能密度进行推导，给出了微模量的具体表达形式，解决了键在非线性变化过程中键的微模量为常数的情况，与岩石类材料在损伤阶段弹性模量随着应变增加而变化的情况更相符。
附图说明
图1为近场动力学物质点运动示意图；
图2为岩石类材料点对力函数模型图；
图3为岩石类材料本构模型仿真方法流程图；
图4为本实施例2计算模型图；
图5为单一非直裂隙试件裂纹扩展过程图；
图6为单一非直裂隙试件应力应变曲线图。
具体实施方式
以下结合本发明实施例中的附图来对本发明的技术方案进行完整清晰的描述。
实施例1
如图1所示，在近场动力学中，将宏观连续体在空间域Ω离散成大量带有体积Vx和质量密度ρ的物质点。物质点x只和其近场领域内H＝{x′∈R,||x′-x||＜δ}的物质点x′存在相互作用，与该范围以外的点没有相互作用。它们之间的相互作用力用点对力函数f来表示。根据牛顿第二定律，物质点在t时刻的运动方程可表示为
其中，x′为x近场范围内的物质点；u(x′,t)和u(x,t)分别表示物质点x′和x在t时刻的位移；ρ(x)和üuü(x,t)分别表示物质点的密度和加速度；dVx′为物质点x′的体积微元；H＝{x′∈R,||x′-x||＜δ}为物质点x近场范围内物质点x′的集合；f表示与材料本身属性相关的本构力函数，b(x,t)为体力密度。
键基近场动力学中的弹脆性材料模型点对力函数f可表示为：
其中，s0表示临界伸长率。
传统的键基近场动力学中的弹脆性模型认为当物质点对的伸长率超过临界伸长率s0，物质点间的键断裂，此时物质点间没有本构力，这种模型不适用于岩石类材料的破坏情况。据此，本发明提出了一种近场动力学岩石类材料本构模型的构建方法，如图2所示，一种近场动力学岩石类材料本构模型的构建方法包括以下步骤：
S11根据岩石类材料变形特征，确定本构模型中点对力函数的表达形式；
S12确定点对力函数的系数；
S13引入键的损伤函数，将点对力函数表达成和损伤函数相关的形式；
S14对步骤S11中的点对力函数进行积分，求出键的微势能；
S15对步骤S14中键的微势能进行积分，求出物质点的应变能密度；
S16令物质点的应变能密度等于连续介质力学的应变能密度，求出微模量c；
S17将步骤S16中求得的微模量代入步骤S13中的点对力函数，完成岩石类材料本构模型的构建；
步骤S11中所述点对力函数的表达形式为：
其中，sc表示压缩临界伸长率；s2c表示压缩强化临界伸长率，此时压缩近场力最大；s1c表示压缩弹性伸长率；st表示拉伸临界伸长率；s2t表示拉伸强化临界伸长率，此时拉伸近场力最大；s1t表示拉伸弹性伸长率。
在步骤S11所述的本构力函数中，当s1c≤s≤s1t时，键处于线弹性变形阶段，点对力与伸长率的关系为线性关系，键无损伤；当s1t＜s＜st或s1c＜s＜sc时，键处于非线性变形阶段。对于非线性变形阶段又可分为非线性强化变形阶段和非线性软化变形阶段，当s1t＜s＜s2t或s1c＜s＜s2c时，键处于非线性强化变形阶段，键开始损伤，微裂纹开始产生，点对力随着伸长率增加非线性增加；当s2t＜s＜st或s2c＜s＜sc时，键处于非线性软化阶段，点对力随着伸长率增加非线性减小；当s＞st或s＜sc时，键的损伤为1，键发生断裂，物质点之间不存在相互作用。
所述步骤12包括以下步骤：
S12.1以键的压缩为例，设键非线性变形阶段点对力标量函数为：
f＝A′s2+B′s+C′ (21)
S12.2确定步骤S12.1所述的点对力标量函数所过点O、P的坐标；
S12.3将步骤S12.2所述的点O，P两点坐标代入步骤S12.1所述的点对力标量函数中，得到点对力和键伸长率的关系式；
S12.4根据二次函数的性质，确定系数A′和B′的关系式；
S12.5根据步骤S12.3所述的点对力和键伸长率的关系式、步骤S12.4所述的A′和B′的关系式，求出系数A′、B′和C′的表达式。
所述步骤S12.2中O，P的坐标分别为(s1c，cs1c)，(sc，0)。
所述步骤S12.3中点对力和键伸长率的关系式为：
步骤S12.5中所述的系数A′、B′、C′的表达式为：
其中，sc表示压缩临界伸长率；s2c表示压缩强化临界伸长率，此时压缩近场力最大；s1c表示压缩弹性伸长率。
同理可得
其中，st表示拉伸临界伸长率；s2t表示拉伸强化临界伸长率，此时拉伸近场力最大；s1t表示拉伸弹性伸长率。
考虑键的损伤，步骤S11所述的点对力函数可以表示为：
其中，α(s)是表征键渐进损伤的标量函数，可表示成
其中， 当α(s)＝0键没有损伤，α(s)＝1键的损伤值为1，键断裂，物质点之间没有相互作用。
步骤S14所述的微势能可通过以下步骤求得：
S14.1微势能表示连接两个物质点之间键的势能，可通过对点对力求积分获得，对于键非线性变形阶段，微势能则包括线弹性变形阶段储存的能量及非线性变形阶段所储存的能量，则微势能函数表示为
S14.2将步骤S11中的点对力函数代入，则微势能可表示为：
S14.3推导前提为各向同性扩张，所以相对位移可表示为η＝sξ，则微势能函数又可表示为
S14.4对步骤S14.3中的微势能函数积分求得微势能的表达式为：
步骤S15所述的应变能密度通过对步骤S14中的微势能求积分获得，表示为：
其中，Vj表示物质点j的积分域，dVj为物质点j的体积微元。由于二维问题和三维问题的积分域不同，所以二维问题和三维问题的应变能密度也不同。
对于二维问题，应变能密度表示为：
对于三维问题，应变能密度可表示为：
步骤S16中所述的微模量c可通过物质点应变能密度和连续介质力学一点处的应变能密度相等求得，连续介质力学一点处的应变能密度可表示为：
对于平面应力问题，步骤S16所述的键的微模量c通过推导可以表示为：
对于三维问题，步骤S16所述的键的微模量c通过推导可以表示为：
对于平面应力情况下步骤S17所述的点对力函数可表示成：
其中， h为厚度；E为弹性模量。
对于三维问题步骤S17所述的点对力函数可表示成：
实施例2
如图3所示，本实施例还提供一种基于近场动力学岩石类材料仿真方法，包括以下步骤：
A1初始化求解域，输入材料参数，几何尺寸。
以二维问题为例，根据仿真需要，确定模型的长度L，宽度W，弹性模量E，密度ρ，泊松比μ，压缩临界伸长率sc，压缩强化临界伸长率s2c，压缩弹性伸长率s1c，表示拉伸临界伸长率st，拉伸强化临界伸长率s2t，拉伸弹性伸长率s1t，领域尺寸δ。
A2离散求解域，生成物质点坐标。
规定长度方向为x轴方向，宽度方向为y轴方向。将模型沿长度方向划分m个晶格，宽度方向划分n个晶格。每个晶格的边长每个物质点的坐标位于晶格的中心。模型物质点总个数为m×n。
按从左到右排序，每个物质点沿x方向的坐标可表示为：
其中，i＝1，2，…，m，表示沿x方向排列晶格的个数。
按从下到上排序，每个物质点沿y方向的坐标可表示为；
其中，j＝1，2，…，n，表示沿y方向排列晶格的个数。
A3确定每个物质点领域内的其它物质点；
检索每个物质点在其领域范围内H＝{x′∈R,||x′-x||＜δ}的其他物质点并编号。
A4初始化所有物质点在其近场领域内的键；
构建物质点领域范围内物质点之间的键，初始化标量函数μ(x,t,ξ)，μ(x,t,ξ)＝1，/>
其中，μ(x,t,ξ)是一个标量函数，用来判断键是否断裂，u(x,t,ξ)＝1，键完好，u(x,t,ξ)＝0，键断裂。代表物质点损伤程度，用其领域内消失的相互作用力在原总的相互作用力中所占的比值来表示：
其中，表示物质点没有损伤，/>表示物质点完全损伤。
A5判断是否需要预制裂纹；
如果需要预制裂纹，则需要把裂纹所在中心线的键断开，即裂纹中心线上的键满足u(x,t,ξ)＝0。裂纹所在中心线的表达式为：
其中，(x0,y0)为裂纹中心点坐标。
A6施加边界条件；
根据模型实际情况，选择合适的边界条件。对于速度和位移边界条件，需要在虚拟边界层施加速度和位移条件，所述虚拟边界层厚度为三倍的物质点间距。对于外荷载边界条件，将外部荷载转换成体力密度施加到边界层上，所述边界层的厚度为物质点间距。
A7确定时间步长。如果采用自适应动态松弛法，则时间步长Δt＝1s。如果不采用自适应动态松弛法，则通过时间步长模型确定时间步长，所述时间步长模型为：
其中，ρ为密度，c为微模量。
A8计算任意物质点总的近场力。计算每个物质点领域范围内物质点之间所有键的键力，并把所有的键力求和即可得到任意物质点总的近场力。
以平面应力问题为例，t时刻所述物质点之间的键力可表示为：
任意物质点总的近场力可表示为：
其中，k表示物质点i领域范围内物质点的总个数。
A9如果伸长率大于临界伸长率，则物质点间的键断裂。判断任意时刻物质点之间键的伸长率s是否满足s＞st或s＜sc，如果满足，则物质点之间的键力为0；如果不满足，按步骤A8计算物质点之间的键力。
所述键的伸长率可表示为：
其中，η＝u′-u为相对位移；ξ＝x′-x表示相对位置；η|表示键的初始长度；|ξ+η|表示键变形之后的长度。
A10如果采用动态松弛法则利用对系统中所有物质点引入虚拟惯性和局部阻尼利用动态松弛算法进行时域积分，如果不采用动态松弛法利用显式向前和向后差分公式进行时域积分，求出位移和速度。
所述采用动态松弛算法的差分格式可表示为：
其中，Δt表示时间步长，Λ表示虚拟对角密度矩阵；c为阻尼系数；L表示PD相互作用力的体积力密度，b表示外力的体积力密度。
其中显式向前和向后差分格式可表示为：
A11输出计算结果。将计算结果导出为txt、dat文件。计算结果是一个txt文件，包括速度、位移、应变能密度以及损伤。
A12根据A11的计算结果绘制位移云图、应变能密度和损伤云图绘制位移，应变能密度和损伤云图。将步骤A11中的txt、dat文件导入到paraview中，通过后处理软件paraview绘制位移、应变能密度和损伤云图。
图5为应用本发明仿真含单一非直裂隙岩石试件单轴压缩的计算效果图，其特点为基于本发明所构建的近场动力学岩石类材料仿真方法，构建了单轴压缩下含单一非直裂隙试件的数值模型。应用本发明所构建的近场动力学岩石类材料本构模型可以模拟裂纹起裂、扩展、贯穿的全过程。翼型拉伸裂纹首先从预制裂隙的尖端起裂并沿着与预制裂隙垂直的方向扩展，之后翼型裂纹扩展到一定长度时，逐渐沿着加载方向扩展。接着在裂纹尖端萌生出次生剪切裂纹，最后剪切裂纹扩展贯穿整个试件。试件最终的破坏形式为拉剪混合破坏。
图6为应用本发明仿真单轴压缩条件下，含单一非直裂隙试件得到的应力应变曲线，结合图5可知，本实施例能够有效仿真单轴压缩下含单一非直裂隙岩石类材料试件的裂纹扩展过程及应力-应变变化规律。
进一步的，所述计算效果的具体实施步骤包括：
A1模型离散化，生成物质点坐标，确定模型参数。
在本实施例中，如图4所示，模型的长80mm，宽160mm，弹性模量E＝13.94GPa，密度ρ＝2449kg/m3，泊松比预制单一非直裂隙长度为L＝30mm，宽度b＝2mm。裂隙长度方向与水平线的夹角α＝45°。
在本实例中，将模型沿长度方向划分160个物质点，沿宽度方向划分320个物质点，物质点间距Δx＝0.5mm，领域尺寸δ＝2.0075mm。模型物质点总数为51200。拉伸强化临界伸长率s2t＝1.42×10-3，拉伸弹性伸长率s1t＝0.3s2t，拉伸临界伸长率st＝5s2t，压缩强化临界伸长率s2c＝-25s2t，压缩弹性伸长率s1c＝0.8s2c，压缩临界伸长率sc＝3s2c。
本实施例中，按从左到右排序，每个物质点沿x方向的坐标可表示为：
x＝-40+0.25+0.5(i-1) (50)
其中，i＝1，2，…，160，表示沿x方向排列晶格的个数。
按从下到上排序，每个物质点沿y方向的坐标可表示为；
y＝-80+0.25+0.5(j-1) (51)
其中，j＝1，2，…，320，表示沿y方向排列晶格的个数。
A2确定每个物质点领域内的其它物质点并初始化领域内的键。
搜索每个物质点在其领域范围内H＝{x′∈R,||x′-x||＜2.0075mm}的其他物质点并编号。对领域内所有的键进行初始化，将标量函数μ(x,t,ξ)、分别设为1和0。
A3预制单一非直裂隙。
把裂纹所在中心线的键断开，即裂纹中心线上的键满足u(x,t,ξ)＝0。裂纹所在中心线的表达式为：
其中，L表示预制裂纹的长度，α表示预制裂纹的夹角，b表示裂纹的宽度。
A4施加边界条件。
在模型上下两端分别设定厚度为1.5mm的虚拟边界层，对虚拟边界层以加载速率±5.0×10-8m/s施加位移荷载。
A5计算物质点的近场力，求解物质点的速度，位移。
计算每个物质点近场领域范围内所有键的键力并求和，即为单个物质点的近场力。利用动态松弛算法求解物质点的速度，位移，时间步长为1s。
A6计算物质点的损伤。
如果键的伸长率s大于7.1×10-3或者s小于-1.065×10-1，则物质点间的键断裂，即u(x,t,ξ)＝0。物质点损伤可用其领域内断裂的键在原总的键中所占的比值来表示：
传统的键基近场动力学认为当键的伸长超过临界伸长率，物质点之间的键断裂，物质点之间的相互作用力消失，这不适合模拟岩石类材料变形破坏特征。本发明提出一种新的岩石类材料本构模型构建方法，考虑键的线性变形阶段和非线性变形阶段，物质点之间的点对力随着伸长率先增加再减小最后消失，更适合模型岩石类材料的破坏。弥补了传统键基近场动力学不能反映岩石类材料随着应变增加，应力先增大再减小最后破坏这一特性的缺陷。

Claims (3)

1.一种基于近场动力学岩石类材料本构模型的建立方法，其特征在于：包括以下步骤：
S11根据岩石类材料变形特征，确定本构模型中点对力函数的表达形式；
所述点对力函数的表达形式为：
其中，sc表示压缩临界伸长率；s2c表示压缩强化临界伸长率，此时压缩近场力最大；s1c表示压缩弹性伸长率；st表示拉伸临界伸长率；s2t表示拉伸强化临界伸长率，此时拉伸近场力最大；s1t表示拉伸弹性伸长率；
S12确定点对力函数的系数；
确定点对力函数系数的步骤为：
S12.1以键的压缩为例，设键非线性变形阶段点对力标量函数为：
f＝A′s2+B′s+C′ (2)
S12.2确定步骤S12.1所述的点对力标量函数所过点O、P的坐标；点O，P的坐标分别为(s1c，cs1c)，(sc，0)；
S12.3将步骤S12.2所述的点O，P两点坐标代入步骤S12.1所述的点对力标量函数中，得到点对力和键伸长率的关系式；
点对力和键伸长率的关系式为：
S12.4根据二次函数的性质，确定A′和B′的关系式；
S12.5根据步骤S12.3所述的点对力和键伸长率的关系式、步骤S12.4所述的A′和B′的关系式，求出系数A′、B′和C′的表达式；
系数A′、B′、C′的表达式为：
其中，sc表示压缩临界伸长率；s2c表示压缩强化临界伸长率，此时压缩近场力最大；s1c表示压缩弹性伸长率；
S13引入键的损伤函数，将点对力函数表达成和损伤函数相关的形式；
所述点对力函数表达成和损伤函数相关的形式为：
其中，c表示键的微模量常数，α(s)是表征键渐进损伤的标量函数，可表示成：
其中，
S14基于步骤S11中的点对力函数，推导出键的微势能；
所述微势能为：
其中，c表示键的微模量常数，s1c表示压缩弹性伸长率，s1t表示拉伸弹性伸长率，
S15基于步骤S14中键的微势能，推导出物质点的应变能密度；
对于二维问题，所述应变能密度为：
其中，h为厚度，c为微模量常数，δ为领域半径，s1c表示压缩弹性伸长率，s1t表示拉伸弹性伸长率，
对于三维问题，所述应变能密度为：
其中，c为微模量常数，δ为领域半径，s1c表示压缩弹性伸长率，s1t表示拉伸弹性伸长率，
S16令物质点的应变能密度等于连续介质力学的应变能密度，求出微模量c的表达式；
对于平面应力问题，所述微模量c的表达式为：
其中，E为弹性模量，h为厚度，δ为领域半径，s1c表示压缩弹性伸长率，s1t表示拉伸弹性伸长率，
对于三维问题，所述微模量c的表达式为：
其中，E为弹性模量，δ为领域半径，s1c表示压缩弹性伸长率，s1t表示拉伸弹性伸长率，
S17将步骤S16中求得的微模量c代入步骤S13中的点对力函数，完成岩石类材料本构模型的构建；
对于平面应力情况下所述点对力函数可表示成：
其中， h为厚度；E为弹性模量，δ为领域半径，s1c表示压缩弹性伸长率，s1t表示拉伸弹性伸长率；
对于三维问题所述点对力函数可表示成：
其中， E为弹性模量，s1c表示压缩弹性伸长率，s1t表示拉伸弹性伸长率，δ为领域半径。
2.根据权利要求1所述的一种基于近场动力学岩石类材料本构模型的建立方法，其特征在于：在步骤S11所述的点对力函数中，当s1c≤s≤s1t时，键处于线弹性变形阶段，点对力与伸长率的关系为线性关系，键无损伤；当s1t＜s＜st或s1c＜s＜sc时，键处于非线性变形阶段；对于非线性变形阶段又可分为非线性强化变形阶段和非线性软化变形阶段，当s1t＜s＜s2t或s1c＜s＜s2c时，键处于非线性强化变形阶段，键开始损伤，微裂纹开始产生，点对力随着伸长率增加非线性增加；当s2t＜s＜st或s2c＜s＜sc时，键处于非线性软化阶段，点对力随着伸长率增加非线性减小；当s＞st或s＜sc时，键的损伤为1，键发生断裂，物质点之间不存在相互作用。
3.一种如权利要求1所述的一种基于近场动力学岩石类材料本构模型的仿真方法，其特征在于：包括以下步骤：
A1初始化求解域，输入材料参数，近场尺寸；
A2离散求解域，根据步骤A1的尺寸生成物质点坐标；
A3确定每个物质点领域内的其它物质点；
A4初始化所有物质点在其近场领域内的键；
A5根据模型实际情况，判断模型是否需要预制裂纹；
A6根据模型实际情况，对模型施加边界条件；
A7如果采用自适应动态松弛法，则时间步长Δt＝1，如果不采用自适应动态松弛法，则确定时间步长；
A8计算每个物质点近场领域范围内所有键的键力并求和；
A9判断任意时刻物质点之间键的伸长率s是否满足s＞st或s＜sc，如果满足，则物质点之间的键力为0；如果不满足，按步骤A8计算物质点之间的键力；
A10如果采用动态松弛算法，则利用动态松弛算法进行时域积分；如果不采用动态松弛法利用显式向前和向后差分公式进行时域积分；最后求出位移和速度；
所述动态松弛算法的差分格式可表示为：
其中，Δt表示时间步长，Λ表示虚拟对角密度矩阵；c为阻尼系数；L表示PD相互作用力的体积力密度，b表示外力的体积力密度；
所述显式向前和向后差分格式可表示为
A11输出步骤A10的计算结果；
A12根据A11的计算结果绘制位移云图、应变能密度和损伤云图。
CN202310059094.9A 2023-01-16 2023-01-16 基于近场动力学岩石类材料本构模型的建立与仿真方法 Active CN116110520B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310059094.9A CN116110520B (zh) 2023-01-16 2023-01-16 基于近场动力学岩石类材料本构模型的建立与仿真方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310059094.9A CN116110520B (zh) 2023-01-16 2023-01-16 基于近场动力学岩石类材料本构模型的建立与仿真方法

Publications (2)

Publication Number Publication Date
CN116110520A CN116110520A (zh) 2023-05-12
CN116110520B true CN116110520B (zh) 2023-07-21

Family

ID=86257629

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310059094.9A Active CN116110520B (zh) 2023-01-16 2023-01-16 基于近场动力学岩石类材料本构模型的建立与仿真方法

Country Status (1)

Country Link
CN (1) CN116110520B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116895351A (zh) * 2023-09-11 2023-10-17 上海交通大学三亚崖州湾深海科技研究院 用于准脆性材料破坏预测的bpd有限元耦合方法及系统

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107727390B (zh) * 2017-08-21 2019-07-16 江铃汽车股份有限公司 汽车动力传动系弯曲特性测试方法
CN111767631B (zh) * 2019-08-02 2022-10-18 中国石油大学(华东) 一种基于多相数字岩心模拟岩石裂纹扩展的方法及系统
CN111274731B (zh) * 2020-02-18 2021-05-18 西南石油大学 一种裂缝性地层压裂裂缝延伸轨迹预测方法
CN112231965B (zh) * 2020-08-06 2024-03-05 沈阳工业大学 一种基于匹配优度的盾构机主控室显控界面综合测试分析方法
CN112131709B (zh) * 2020-08-25 2024-04-19 山东大学 基于近场动力学本构模型的节理岩体力学仿真方法及系统

Also Published As

Publication number Publication date
CN116110520A (zh) 2023-05-12

Similar Documents

Publication Publication Date Title
Bhat et al. A micromechanics based constitutive model for brittle failure at high strain rates
Zheng et al. Reformulation of dynamic crack propagation using the numerical manifold method
Zaccariotto et al. Examples of applications of the peridynamic theory to the solution of static equilibrium problems
Wang et al. Numerical simulation of propagation and coalescence of flaws in rock materials under compressive loads using the extended non-ordinary state-based peridynamics
Wu et al. Investigating the effects of micro-defects on the dynamic properties of rock using numerical manifold method
Munjiza et al. Structural applications of the combined finite–discrete element method
Chen et al. Simulating the breakage of glass under hard body impact using the combined finite-discrete element method
CN116110520B (zh) 基于近场动力学岩石类材料本构模型的建立与仿真方法
CN112131709A (zh) 基于近场动力学本构模型的节理岩体力学仿真方法及系统
Liu et al. Peridynamic modelling of impact damage in three-point bending beam with offset notch
Xu et al. An efficient solid-shell cohesive zone model for impact fracture analysis of laminated glass
CN113569442A (zh) 一种基于rkpm-pd耦合算法的岩石裂纹扩展预测方法
Chen et al. A multiscale method coupling peridynamic and boundary element models for dynamic problems
Pekau et al. Seismic fracture analysis of concrete gravity dams
CN111666699A (zh) 基于rev全区域覆盖的岩体工程跨尺度模拟计算方法
Gorgogianni et al. Mechanism-based energy regularization in computational modeling of quasibrittle fracture
Shahzamanian Implementation of a rate dependent tensile failure model for brittle materials in ABAQUS
Daneshyar et al. Seismic analysis of arch dams using anisotropic damage-plastic model for concrete with coupled adhesive-frictional joints response
Mazars et al. Dynamic behavior of concrete and seismic engineering
Hu et al. Numerical simulations of arbitrary evolving cracks in geotechnical structures using the nonlinear augmented finite element method (N-AFEM)
Li et al. On deformation and fracture of PBX simulant employing modified three-body potential peridynamics model with deformation-based failure criteria
CN111523217A (zh) 夹层玻璃抗冲击性能预测和结构优化方法
Long et al. Numerical verification of energy release rate and J-Integral in large strain formulation
Gaugele et al. Simulation of material tests using meshfree Lagrangian particle methods
Huang et al. A hybrid polymer–water peridynamics model for ballistic penetration damage of soft materials

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