CN106339547A - 一种爆破的数值模拟方法 - Google Patents

Abstract

一种爆破的数值模拟方法，包括:(1)在二维DDA方法中添加圆形单元，用以模拟炸药颗粒；(2)在DDA计算模型中设置爆破孔，根据实际情况中炸药点火后的体积膨胀过程，设置DDA炸药颗粒单元的体积随时间膨胀的函数；(3)在数值模拟过程中，对体积膨胀后的炸药颗粒单元与药孔孔壁单元的接触关系进行判定，根据二者的接触量计算接触反力，从而使炸药对孔壁产生挤压，使孔开裂,膨胀后的炸药颗粒单元始终以碰撞挤压的方式对周围块体产生作用，当炸药颗粒单元膨胀结束后，该作用力消失，被爆物体按照运动方程四散运动；(4)炸药颗粒与被爆物间、被爆物与被爆物间的相互作用，以及被爆物体的应力、变形、运动的力学行为均按照DDA方法进行模拟。

Description

一种爆破的数值模拟方法

技术领域

本发明涉及数值模拟的技术领域，尤其涉及一种爆破的数值模拟方法。

背景技术

爆破的数值模拟，是优化爆破参数，预测被爆物破碎程度及破碎块体运行轨迹，保证爆破安全的一种重要手段。目前，对爆破进行数值模拟的主要方法是，先通过实验测量出爆炸过程中炸药作用于装药孔壁的压力波形，然后将该压力波以力的形式施加在数值计算模型的孔壁上，从而模拟爆破过程中炸药爆炸对孔壁及被爆物体的作用。然而，由于爆破作用力波形的测量误差较大且爆破力波形与被爆物体的力学参数密切相关，模拟中所施加的荷载与真实荷载往往有较大偏差，因此这种方法的模拟精度难以保证，模拟的结果偏离实际。

20世纪80年代，石根华博士提出了一种基于非连续介质力学的数值方法，非连续变形分析(Discontinuous Deformation Analysis，简称DDA)。该方法以被天然节理裂隙切割成的任意块体为基本单元，以单元的刚体位移和变形为基本未知量，基于一套高效的接触搜索方法，能得到块体系统中任意可能的接触形式。同时，该方法充分考虑了块体间的相互作用，基于最小势能原理建立总体平衡方程，并采用隐式解法进行求解。非连续变形分析方法在模拟块体系统的接触问题时具有先天性的优势。目前，DDA方法已经广泛应用于各种岩体工程的仿真计算分析。

发明内容

为克服现有技术的缺陷，本发明要解决的技术问题是提供了一种爆破的数值模拟方法,其能够更好地模拟爆破过程。

本发明的技术方案是：这种爆破的数值模拟方法，该方法包括以下步骤:

(1)在二维非连续变形分析DDA方法中添加圆形单元，用以模拟炸药颗粒,对圆形单元的形状、力学参数进行定义，并确定圆形单元与圆形单元、圆形单元与多边形单元之间的接触形式；

(2)在非连续变形分析DDA的计算模型中设置爆破孔，在孔内设置炸药颗粒单元，根据实际情况中炸药点火后的体积膨胀过程，设置DDA炸药颗粒单元的体积随时间膨胀的函数；

(3)在数值模拟过程中，对体积膨胀后的炸药颗粒单元与药孔孔壁单元的接触关系进行判定，根据二者的接触量计算接触反力，从而使炸药对孔壁产生挤压，使孔开裂,膨胀后的炸药颗粒单元始终以碰撞挤压的方式对周围块体产生作用，当炸药颗粒单元膨胀结束后，该作用力消失，被爆物体按照运动方程四散运动，在体积膨胀过程中炸药颗粒单元始终保持质量守恒；

(4)炸药颗粒与被爆物间、被爆物与被爆物间的相互作用，以及被爆物体的应力、变形、运动的力学行为均按照DDA方法进行模拟。

本发明依据爆破的物理过程，设定炸药颗粒的体积随时间膨胀的函数，使膨胀后的炸药颗粒对药孔孔壁产生剧烈挤压，从而获得了更加接近真实的爆破数值模拟结果，规避了上述的现有方法中的问题，因此能够更好地模拟爆破过程。

附图说明

图1为炸药颗粒的体积膨胀函数。

图2为不同接触形式时的接触距离示意图。其中，图2(a)为圆与圆接触；图2(b)为圆与多边形的边接触；图2(c)为圆与多边形的角接触。

图3为不同接触形式时的接触反力示意图。其中，图3(a)为圆与圆接触；图3(b)为圆与多边形的边接触；图3(c)为圆与多边形的角接触。

图4为炸药颗粒膨胀法模拟流程框图。

具体实施方式

真实的爆破物理过程是，点火后炸药起爆，随后体积迅速膨胀急剧气化，对装药孔孔壁产生强烈的挤压作用，导致被爆物爆裂。被爆物炸开后，爆孔体积增大，膨胀气体对孔壁的作用力逐渐衰减，被爆物的碎块则在力的作用下作四散运动。因此，直接模拟爆炸过程中炸药颗粒体积的膨胀过程，从而对药孔孔壁四周物质产生挤压冲击作用，是正确模拟爆破的关键。

本发明的爆破的数值模拟方法，该方法包括以下步骤:

(1)在二维非连续变形分析DDA方法中添加圆形单元，用以模拟炸药颗粒,对圆形单元的形状、力学参数进行定义，并确定圆形单元与圆形单元、圆形单元与多边形单元之间的接触形式；

(2)在非连续变形分析DDA的计算模型中设置爆破孔，在孔内设置炸药颗粒单元，根据实际情况中炸药点火后的体积膨胀过程，设置DDA炸药颗粒单元的体积随时间膨胀的函数；

(3)在数值模拟过程中，对体积膨胀后的炸药颗粒单元与药孔孔壁单元的接触关系进行判定，根据二者的接触量计算接触反力，从而使炸药对孔壁产生挤压，使孔开裂,膨胀后的炸药颗粒单元始终以碰撞挤压的方式对周围块体产生作用，当炸药颗粒单元膨胀结束后，该作用力消失，被爆物体按照运动方程四散运动，在体积膨胀过程中炸药颗粒单元始终保持质量守恒；

(4)炸药颗粒与被爆物间、被爆物与被爆物间的相互作用，以及被爆物体的应力、变形、运动的力学行为均按照DDA方法进行模拟。

本发明依据爆破的物理过程，设定炸药颗粒的体积随时间膨胀的函数，使膨胀后的炸药颗粒对药孔孔壁产生剧烈挤压，从而获得了更加接近真实的爆破数值模拟结果，规避了上述的现有方法中的问题，因此能够更好地模拟爆破过程。

优选地，所述步骤(2)中，在DDA模型中设置不同的爆孔空间排布，不同的爆孔几何尺寸，不同的膨胀函数来模拟不同的爆破形式。

优选地，所述步骤(2)中，定义炸药颗粒的形状、力学及爆炸参数如下：

$x_i^0,y_i^0,{\mathrm{\γ}}_i^0,E_i^0,i=1,2,......N---\left(1\right)$

其中：为第i号炸药颗粒的初始坐标与半径，为炸药颗粒的初始弹性模量，N为炸药颗粒总数；

r_i(t)＝f_i(t)r_i ⁰ (2)

其中：r_i(t)为t时刻第i号炸药颗粒的半径，f_i(t)为t时刻第i号炸药颗粒体积膨胀倍数。

优选地，所述步骤(3)包括以下分步骤：

(3.1)在时步迭代之初，令炸药颗粒单元膨胀，设第k个计算时步的初始时间为t_k，则此时炸药的颗粒半径为r_i(t_k)，i＝1,2,……N在t_k+1＝t_k+Δt_k时刻，给定炸药颗粒的半径为r_i(t_k+1)，i＝1,2,……N；

(3.2)在接触搜索中，以新的颗粒单元半径r_i(t_k+1)进行圆与圆，圆与多边形的接触搜索，找出颗粒半径膨胀后的新的接触关系；两个圆形颗粒单元接触的判别公式为：

D(m,n)≤r_m(t_k+1)+r_n(t_k+1) (3)其中：D(m,n)为两个圆形单元m、n的形心间的距离，当两个圆形单元形心之间的距离小于两个单元的半径之和时，两个圆形单元发生接触；圆形单元与多边形单元的某一条边接触的判别公式为：

D(m,n_j)≤r_m(t_k+1) (4)

其中：D(m,n_j)为圆形单元m的形心到多边形单元n的第j条边的距离，当圆形单元m的形心到多边形单元n的第j条边的距离小于圆形单元m在t_k+1时刻的半径时，圆形单元m与多边形单元n的第j条边接触；

圆形单元与多边形单元的角点接触的判别公式为：

D(m,n_p)≤r_m(t_k+1) (5)

其中：D(m,n_p)为圆形单元m的形心与多边形单元n的第p个顶点之间的距离，当圆形单元m的形心与第n多边形的第p个顶点之间的距离小于或等于圆形单元m在t_k+1时刻的半径时，圆形单元m与多边形单元n的第p个顶点接触；

根据公式(3)、(4)，(5)搜索出炸药颗粒单元体积膨胀后各单元之间的接触关系；

(3.3)计算炸药颗粒单元膨胀后的力：首先计算相互接触的单元之间的接触量；

两个圆形颗粒单元接触时的接触量为：

d(m,n)＝D(m,n)-r_m(t_k+1)-r_n(t_k+1) (6)

圆形单元与多边形单元的某一条边接触时的接触量为：

d(m,n_j)＝D(m,n_j)-r_m(t_k+1) (7)

圆形单元与多边形单元的角点接触时的接触量为：

d(m,n_p)＝D(m,n_p)-r_m(t_k+1) (8)

利用公式(6)、(7)、(8)计算得到不同接触形式的接触力：

P＝Kd (9)

其中：P为接触力，K为接触处弹簧的刚度，d为接触量，对于不同的接触形式，接触力d可分别由公式(6)～(8)求得；

(3.4)当炸药颗粒单元由于体积膨胀与其他单元产生互相侵入时，在两个单元相互侵入的部位，施加一对大小相等，方向相反的力，力的作用点作用于两单元接触处，位于圆形单元的径向连线上，方向为沿圆形单元的半径方向指向圆心；力的大小按公式(9)计算；

(3.5)炸药颗粒单元保持质量守恒，在爆炸膨胀过程中虽然体积不断增大，但总质量保持不变，质量密度不断变小，单元的质量变化按下式：

$m\left(t\right)=\frac{m_0}{\stackrel{2}{f\left(t\right)}}---\left(10\right)$

其中，m₀为炸药颗粒的初始质量，f(t)为炸药颗粒的体积膨胀倍数，m(t)为t时刻炸药颗粒的质量。

其他如弹性刚度矩阵，接触刚度矩阵，力矩阵，开闭迭代等的计算方法仍沿用DDA方法。

以上所述，仅是本发明的较佳实施例，并非对本发明作任何形式上的限制，凡是依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属本发明技术方案的保护范围。

Claims (4)

1.一种爆破的数值模拟方法，其特征在于：该方法包括以下步骤:

(1)在二维非连续变形分析DDA方法中添加圆形单元，用以模拟炸药颗粒,对圆形单元的形状、力学参数进行定义，并确定圆形单元与圆形单元、圆形单元与多边形单元之间的接触形式；

(2)在非连续变形分析DDA的计算模型中设置爆破孔，在孔内设置炸药颗粒单元，根据实际情况中炸药点火后的体积膨胀过程，设置DDA炸药颗粒单元的体积随时间膨胀的函数；

(3)在数值模拟过程中，对体积膨胀后的炸药颗粒单元与药孔孔壁单元的接触关系进行判定，根据二者的接触量计算接触反力，从而使炸药对孔壁产生挤压，使孔开裂,膨胀后的炸药颗粒单元始终以碰撞挤压的方式对周围块体产生作用，当炸药颗粒单元膨胀结束后，该作用力消失，被爆物体按照运动方程四散运动，在体积膨胀过程中炸药颗粒单元始终保持质量守恒；

(4)炸药颗粒与被爆物间、被爆物与被爆物间的相互作用，以及被爆物体的应力、变形、运动的力学行为均按照DDA方法进行模拟。

2.根据权利要求1所述的爆破的数值模拟方法，其特征在于：

所述步骤(2)中，在DDA模型中设置不同的爆孔空间排布，不同的爆孔几何尺寸，不同的膨胀函数来模拟不同的爆破形式。

3.根据权利要求2所述的爆破的数值模拟方法，其特征在于：所述步骤(2)中，定义炸药颗粒的形状、力学及爆炸参数如下：

$x_i^0,y_i^0,{\mathrm{\γ}}_i^0,E_i^0,i=1,2,......N---\left(1\right)$

其中：为第i号炸药颗粒的初始坐标与半径，为第i号炸药颗粒的初始弹性模量，N为炸药颗粒总数；

r_i(t)＝f_i(t)r_i ⁰ (2)

其中：r_i(t)为t时刻第i号炸药颗粒的半径，f_i(t)为t时刻第i号炸药颗粒体积膨胀倍数。

4.根据权利要求3所述的爆破的数值模拟方法，其特征在于：所述步骤(3)包括以下分步骤：

(3.1)在计算时步迭代之初，令炸药颗粒单元膨胀，设第k个计算时步的初始时刻为t_k，则此时炸药颗粒的半径为r_i(t_k)，i＝1,2,……N在t_k+1＝t_k+Δt_k时刻，给定炸药颗粒的半径为r_i(t_k+1)，i＝1,2,……N；

(3.2)在接触搜索中，以新的颗粒单元半径r_i(t_k+1)进行圆与圆，圆与多边形的接触搜索，确定炸药颗粒膨胀后的新的接触关系；

两个圆形颗粒单元接触的判别公式为：

D(m,n)≤r_m(t_k+1)+r_n(t_k+1) (3)

其中：D(m,n)为两个圆形单元m、n的形心间的距离，当两个圆形单元形心之间的距离小于两个单元的半径之和时，两个圆形单元发生接触；圆形单元与多边形单元的某一条边接触的判别公式为：

D(m,n_j)≤r_m(t_k+1) (4)

其中：D(m,n_j)为圆形单元m的形心到多边形单元n的第j条边的距离，当圆形单元m的形心到多边形单元n的第j条边的距离小于圆形单元m在t_k+1时刻的半径时，圆形单元m与多边形单元n的第j条边接触；

圆形单元与多边形单元的角点接触的判别公式为：

D(m,n_p)≤r_m(t_k+1) (5)

其中：D(m,n_p)为圆形单元m的形心与多边形单元n的第p个顶点之间的距离，当圆形单元m的形心与第n多边形的第p个顶点之间的距离小于或等于圆形单元m在t_k+1时刻的半径时，圆形单元m与多边形单元n的第p个顶点接触；

根据公式(3)、(4)，(5)搜索出炸药颗粒单元体积膨胀后各单元之间的接触关系；

(3.3)计算炸药颗粒单元膨胀后的力：首先计算相互接触的单元之间的接触量；

两个圆形颗粒单元接触时的接触量为：

d(m,n)＝D(m,n)-r_m(t_k+1)-r_n(t_k+1) (6)

圆形单元与多边形单元的某一条边接触时的接触量为：

d(m,n_j)＝D(m,n_j)-r_m(t_k+1) (7)

圆形单元与多边形单元的角点接触时的接触量为：

d(m,n_p)＝D(m,n_p)-r_m(t_k+1) (8)

利用公式(6)、(7)、(8)计算得到不同接触形式的接触力：

P＝Kd (9)

其中：P为接触力，K为接触处弹簧的刚度，d为接触量，对于不同的接触形式，接触量d分别由公式(6)～(8)求得；

(3.4)当炸药颗粒单元由于体积膨胀与其他单元产生互相侵入时，在两个单元相互侵入的部位，施加一对大小相等，方向相反的力，力的作用点作用于两单元接触处，位于圆形单元的径向连线上，方向为沿圆形单元的半径方向指向圆心；力的大小按公式(9)计算；

(3.5)炸药颗粒单元保持质量守恒，在爆炸膨胀过程中虽然体积不断增大，但总质量保持不变，质量密度不断变小，单元的质量变化按下式：

$m\left(t\right)=\frac{m_0}{\stackrel{2}{f\left(t\right)}}---\left(10\right)$

其中，m₀为炸药颗粒的初始质量，f(t)为炸药颗粒的体积膨胀倍数，m(t)为t时刻炸药颗粒的质量。