适用于高压、高应变率下的混凝土损伤本构模型构建方法与流程

1.本发明所属于防灾减灾工程及防护工程领域，具体是ls
‑
dyna商用有限元软件中混凝土本构模型的二次开发，应用于混凝土材料冲击爆炸等强动载作用下的数值仿真。


背景技术：

2.混凝土材料广泛应用于承受冲击和爆炸等强动载作用的军用与民用防护工程中。与静载、地震和风荷载等不同，强动载作用下混凝土材料受高静水压力、高应变率作用并产生大变形、开坑、震塌等损伤破坏，且在弹体及爆炸近区存在高围压区域。随着计算能力的迅速发展，数值模拟在混凝土结构在强动荷载作用下的响应分析中发挥着越来越重要的作用，适当的混凝土材料模型对于准确预测结构动力响应和失效是至关重要的。现有有限元计算软件中常用混凝土材料都是国外学者开发，如holmquist
‑
johnson
‑
cook(hjc)，riedel
‑
hiermaier
‑
thoma(rht)以及karagozian&case(k&c)模型等，模型的使用与开发长期受制于国外，为打破这一技术瓶颈，开发了能准确反映强动载下混凝土动态力学性能的新型本构模型。
3.强动载下的常用混凝土材料模型存在的问题与缺点：(1)hjc模型虽然已广泛应用于弹体侵彻的数值模拟中，可较好描述混凝土的动态压缩行为，但无法较好描述动态拉伸行为；此外hjc模型只采用了两个应力不变量来确定当前加载面，不能描述偏平面形状由低压时三角形向高压时圆形的过渡；并且hjc模型采用了j2流动法则，无法描述混凝土材料的剪胀行为；将偏量引起的微裂纹扩展和球量引起的孔隙坍塌完全分离，未考虑两种损伤机制之间的关联。
4.(2)k&c模型最初开发用来预测爆炸载荷作用下混凝土结构的动态响应，现已广泛应用于强动载下混凝土结构动态响应中。k&c模型采用三个独立的强度面，即初始屈服强度面，最大强度表面和残余强度面，当前破坏面由损伤因子内插屈服强度面和最大的强度面以及最大强度面和残余强度面得到。损伤因子先由零增大到1，对应于应变硬化阶段；然后再从1减小到零，对应于应变软化。k&c模型较为全面描述了冲击爆炸荷载作用下混凝土材料的动态力学行为，然而该模型仍存在如下几个缺点：k&c模型提供了参数自动生成算法并被广泛采用，自动生成的参数适用于较低静水压力情况(如远距离爆炸)，但不适用于高压力情况(如弹体侵彻问题等)；k&c模型预测的断裂应变和断裂能均随应变率增大而增大，与实验现象不符；采用了关联或部分关联的流动法则，仅考虑了偏量行为对球量行为的影响(剪切膨胀)，从而忽略了球量对偏量的影响，即孔隙坍塌引起的静水压缩损伤会降低混凝土单轴抗拉和抗压强度。(3)rht模型可以认为是hjc模型的改进模型，并广泛应用于强动载下混凝土结构的响应和破坏的数值模拟中，采用三个独立的强度面，即初始屈服强度面、最大强度面和残余强度面，当前加载面通过等效塑性应变内插屈服强度面和最大强度面确定应变硬化效应，并通过损伤因子内插最大强度面和残余强度面确定应变软化效应。由于残余强度面中缺少第三应力不变量，因此在特定应力状态下模型预测结果会出现问题；此外rht模型采用拉伸线性软化模型，与已有实验数据不符；与hjc模型一致，未考虑剪切损伤与
静水压缩之间的竞争机制。
5.可以看出，已有混凝土材料模型均侧重于强动载下响应的特定方面，存在其自身的缺陷和局限性，很难真实完整描述混凝土材料动态力学性能，其数值计算结果对工程设计可靠性及参考性有限，故急需一种适用于高压、高应变率下的新型混凝土损伤本构模型。


技术实现要素：

6.本发明旨在针对现有混凝土本构模型缺乏考虑高围压下的屈服面形式、忽略静水压缩损伤及未考虑各损伤之间相互影响关系和现有的应变率效应不能反映应变率突变时应力迟滞效应的缺陷，提供一种适用于高压、高应变率下的混凝土损伤本构模型构建方法。
7.为实现该技术目的，本发明采用以下方案。
8.一种适用于高压、高应变率下的混凝土损伤本构模型构建方法，包括引入修正系数来修正高围压段屈服强度面形式；通过引入拉伸损伤指标、剪切损伤指标和静水压缩损伤指标以及剪切损伤对拉伸损伤影响因子和静水压缩损伤对其他应力状态影响因子，确定等效塑性应变；引入时间松弛因子确定当前时间步的结束时刻的等效应变率。
9.进一步地，修正后的屈服强度面形式表示如下：
[0010][0011]
其中σ
y
是屈服强度面，r
p
＝p/dif，p是压力，dif为动态增强因子，t是混凝土单轴抗拉强度，f
c
是混凝土单轴抗压强度，f
yc
是混凝土抗压屈服强度，a
1y
、a
2y
是材料常数，n是修正系数。
[0012]
进一步地，所述等效塑性应变表示为：
[0013]
λ＝1
‑
(1
‑
δ
s
λ
t
)(1
‑
λ
s
)(1
‑
δ
h
λ
h
)
[0014]
其中λ是等效塑性应变，λ
t
是拉伸损伤指标、λ
s
是剪切损伤指标、λ
h
是静水压缩损伤指标，δ
s
是剪切损伤对拉伸损伤影响因子，δ
h
是静水压缩损伤对其他应力状态影响因子。通过等效塑性应变的表达式能够确定拉伸损伤与剪切损伤之间的相互影响和静水压缩损伤对其他应力状态的影响。
[0015]
再进一步地，拉伸损伤指标表示如下：
[0016][0017]
其中是等效塑性应变增量，b1是损伤常数。
[0018]
再进一步地，剪切损伤指标表示如下：
[0019][0020]
其中b2是损伤常数。
[0021]
再进一步地，静水压缩损伤指标表示如下:
[0022][0023]
其中dμ
p
是塑性体积应变增量，b3损伤常数。
[0024]
进一步地，当前时间步的结束时刻的等效应变率表示如下：
[0025]
其中是瞬时应变率，是当前时间步开始时刻的等效应变率，是当前时间步结束时刻的等效应变率，η0是时间松弛因子，δt是时间步，n和n 1是一个时间步的开始时刻与结束时刻，λ是等效塑性应变。
[0026]
进一步地，所述方法还包括：
[0027]
确定最大强度面和残余强度面，当等效塑性应变小于等于峰值应力时等效塑性应变时，通过插值屈服强度面和最大强度面描述混凝土材料应变硬化现象；
[0028]
当等效塑性应变大于峰值应力时等效塑性应变时，通过插值最大强度面和残余强度面描述混凝土材料应变软化或损伤现象。
[0029]
本发明所取得的有益技术效果：(1)本发明引入修正系数来修正高围压段屈服强度面形式，解决了现有屈服强度面无法准确描述高围压下混凝土材料力学性能，即试验结果无法找到混凝土材料的明显屈服点；
[0030]
(2)通过引入拉伸损伤指标、剪切损伤指标和静水压缩损伤指标以及剪切损伤对拉伸损伤影响因子和静水压缩损伤对其他应力状态影响因子，确定等效塑性应变，进而能够确定拉伸损伤与剪切损伤之间的相互影响和静水压缩损伤对其他应力状态的影响，解决了现有损伤机制无法考虑静水压缩损伤以及各损伤之间相互影响关系的问题；
[0031]
(3)引入了等效应变率效应，通过引入时间松弛因子达到描述应变率突变时应力增长滞后现象，而现有应变率效应无法考虑应变率突变时应力迟滞效应。
附图说明
[0032]
图1是本发明具体实施例中强度面示意图；
[0033]
图2是本发明具体实施例损伤机制示意图；
[0034]
图3是本发明具体实施例考虑应变率迟滞效应与未考虑应变率突变效应的对比示意图；
[0035]
图4是本发明具体实施例的流程示意图；
[0036]
图5是本发明具体实施例的应力更新过程流程图；
[0037]
图6是本发明具体实施实施例的技术方案示意图；
[0038]
图7是本发明具体实施例的单单元三轴压缩计算结果；
[0039]
图8是本发明具体实施例的单单元复杂应力状态计算结果；
[0040]
图9是本发明具体实施例的单单元动态应变率计算结果；
[0041]
图10(1)是利用本发明计算的有限厚度混凝土靶体侵彻结果；
[0042]
图10(2)是利用本发明计算的侵彻过程中靶体的拉伸损伤；
[0043]
图10(3)是利用本发明计算的侵彻过程中靶体的剪切损伤；
[0044]
图10(4)是利用本发明计算的侵彻过程中靶体的静水压缩损伤；
[0045]
图11(1)是利用本发明针对混凝土薄板的锥形装药接触爆炸计算结果；
[0046]
图11(2)是混凝土薄板的锥形装药接触爆炸试验结果；
[0047]
图11(3)是利用本发明针对混凝土厚板的锥形装药接触爆炸计算结果；
[0048]
图11(4)是混凝土厚板的锥形装药接触爆炸试验结果。
具体实施方式
[0049]
以下结合说明书附图和具体实施例对本发明做进一步详细说明。
[0050]
实施例：适用于高压、高应变率下的混凝土损伤本构模型构建方法，详细如下描述。
[0051]
图1是具体实施例中强度面示意图。图中1是屈服强度面，2是最大强度面，3是残余强度面。当前加载面为分段函数，当λ≤λ
m
时，通过插值屈服强度面和最大强度面描述混凝土材料应变硬化现象，当λ＞λ
m
时，通过插值最大强度面和残余强度面描述混凝土材料应变软化(损伤)现象。即
[0052][0053]
其中，y(σ
ij
,η)为当前加载面，σ
ij
是应力张量，η是强度面比例因子，λ是等效塑性应变，λ
m
是峰值应力时等效塑性应变，j2是第二应力偏量不变量，dif为动态增强因子，σ
y
是屈服强度面，σ
m
是最大强度面，σ
r
是残余强度面，r'是当前子午线与压缩子午线比值，即
[0054][0055]
其中，θ是lode角，即
[0056][0057]
其中，ψ是拉压子午线比值，s
ij
是偏应力张量。
[0058]
引入修正系数来修正高围压段屈服强度面形式，屈服强度面由多段函数组成，即
[0059][0060]
其中，r
p
＝p/t，t和f
c
分别是混凝土单轴抗拉强度和单轴抗压强度，f
yc
是混凝土抗压屈服强度，a
1y
、a
2y
分别是材料常数，a
1y
＝0.9088，a
2y
＝0.075/f
c
，n为修正系数，n的建议取值0.45～2，当n取1时与实验值吻合较好。
[0061]
当r
p
≤0时，峰值压力可达
‑
t，通过插值三轴拉伸点(
‑
t,0)和在零压力时的0.45倍最大强度面(0,1.35t)；
[0062]
当时，由于缺乏有效的试验数据，对初始屈服点(f
yc
/3,f
yc
)和在零压力时的0.45倍最大强度面(0,1.35t)线性插值，其中f
yc
＝0.45f
c
；
[0063]
当时，根据三轴压缩试验数据定义，σ
y
＝0.45σ
m
；
[0064]
当r
p
＞nf
c
时，屈服强度面与压力无关。
[0065]
最大强度面由多段函数构成，即
[0066][0067]
其中，a1和a2是材料常数，由三轴压缩试验数据确定，a1＝0.5857，a2＝0.025/f
c
。是拉压子午线比值，即
[0068][0069]
拉压子午线比值由离散点插值决定。由于最大强度面为多段函数定义，当r
p
＞f
c
时，可根据无侧限单轴压缩和三轴压缩实验数据定义，当0＜r
p
≤f
c
/3时，没有有效的实验数据，由于压力低于单轴抗压强度，故采用拉伸子午线除以拉压子午线比值的方法。
[0070]
对残余强度面进行了修正，弥补了混凝土靶体损伤后材料偏软的缺陷。即
[0071][0072]
其中，a3是材料常数，a3＝0.5。
[0073][0074]
其中，n、m、k是材料参数，其中，n、m、k是材料参数，m＝0.29，k＝1.86；当λ≤λ
m
时，η从0到1变化，此时由屈服强度面向最大强度面变化，捕捉应变硬化现象；当λ＞λ
m
时，η从1到0变化，此时由最大强度面向残余强度面变化，捕捉应变软化现象。
[0075]
图2是本发明实施例损伤机制示意图，本实施例通过引入拉伸损伤指标、剪切损伤指标和静水压缩损伤指标以及剪切损伤对拉伸损伤影响因子和静水压缩损伤对其他应力状态影响因子，确定等效塑性应变。本发明引入了由混凝土内部孔隙坍塌引起的静水压缩损伤，并引入了与应力状态相关的变量，考虑各损伤之间相互影响关系。即
[0076]
λ＝1
‑
(1
‑
δ
s
λ
t
)(1
‑
λ
s
)(1
‑
δ
h
λ
h
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)
[0077]
其中，λ
t
、λ
s
、λ
h
分别是拉伸损伤、剪切损伤和静水压缩损伤指标，δ
s
是剪切损伤对拉伸损伤影响因子，δ
h
是静水压缩损伤对其他应力状态影响因子。即
[0078][0079]
拉伸损伤指标、剪切损伤指标和静水压缩损伤指标如下所示，
[0080][0081][0082][0083]
其中，是等效塑性应变增量，dμ
p
是塑性体积应变增量，b1、b2、b3是损伤常数，b1＝1.6，b2＝1.0，b3＝1.5。
[0084]
通过4个简单工况对本发明中损伤机制做出解释。
[0085]
(1)当混凝土处在受拉状态，r
p
＜0，只有λ
t
累积，当处于压缩状态时，λ
s
和λ
h
同时存在，λ
h
是由孔隙坍塌引起的塑性体积应变累积而成，当混凝土处于高围压状态，δ
h
＝0，λ
h
对外没有影响；当混凝土处于静水压力状态，无偏量产生，λ
s
＝0，孔隙坍塌后仅存在体积硬化，与现有实验结果一致。
[0086]
(2)当混凝土应力状态由三轴压缩向非三轴压缩变化时，δ
h
＝1，λ
h
对λ有贡献，随着三向围压增大，λ
h
越大，与现有实验结果一致。
[0087]
(3)当混凝土从受拉状态变为受压状态，假设之前没有剪切损伤和静水压缩损伤，λ
s
＝0，δ
h
＝0，所以，λ＝1
‑
(1
‑
δ
s
λ
t
)，拉伸时，r
p
＜0，λ
s
＝1，λ＝λ
t
，当受压时，r
p
≥0，δ
s
＝0，λ＝λ
s
，显然，拉伸损伤对后续压缩行为没有影响。
[0088]
(4)当混凝土从受压状态变为受拉状态，假定之前没有拉伸损伤累积，λ
t
＝0，所以只有剪切损伤和静水压缩损伤存在，λ＝1
‑
(1
‑
λ
s
)(1
‑
λ
h
)，因此，前期的剪切损伤和静水压缩损伤对后续的拉伸性能有很大影响。当应力状态由静水压缩状态变为拉伸时，λ＝λ
h
反应了孔隙坍塌对混凝土抗拉强度的折减。
[0089]
为了便于结果输出和后处理，统一了输出的各损伤分量，即
[0090][0091]
其中，d
i
是输出的各损伤分量，λ
i
是不同应力状态下累积的等效塑性应变，t表示拉伸应力状态，s表示剪切应力状态，h表示静水压缩应力状态，ρ
c
为材料损伤参数，
[0092]
图3是本发明考虑应变率迟滞效应与未考虑应变率突变效应的对比示意图，本实施例引入时间松弛因子确定当前时间步的结束时刻的等效应变率。1是当在t时刻应变率突变时，应力突然增加；2是应力随应变率的突增而增加缓慢。混凝土是一种率敏感材料，对准确预测强动载下材料的破坏具有重要意义。本发明通过利用动态增强因子dif径向放大强度面实现应变率效应的考虑，现有应变率效应无法捕捉应变率突变实验下观测到的应力迟滞响应，故引入了时间松弛因子，解决这一问题。即
[0093][0094]
其中，是瞬时应变率，是是当前时间步开始时刻的等效应变率，η0是时间松弛因子，η0＝0.1ms，δt是时间步，n和n 1是一个时间步的开始与结束。
[0095]
图4是本发明实施例具体实施的流程图，包括：
[0096]
a1：基于流体弹塑性力学进行应力更新。
[0097]
用户可根据弹塑性力学中屈服准则、流动法则和一致性条件进行应力更新推导。
[0098]
图5是实施例中应力更新过程流程图。具体的当前加载面为
[0099][0100]
流动势函数为
[0101][0102]
其中，ω为塑性流动关联因子，当ω＝0时为非关联流动(j2流动)，当ω＝1时为完全关联，当ω介于0和1之间为部分关联。
[0103]
将式(16)带入式(17)，得
[0104][0105]
流动法则为
[0106][0107]
其中，du为一致性参数。
[0108]
将式(18)代入式(19)，得
[0109][0110]
其中，y
,p
＝dy(σ
ij
,η)/dp。
[0111]
等效塑性应变增量为
[0112][0113]
将式(20)代入式(21)，得
[0114][0115]
等效塑性应变为
[0116]
λ＝1
‑
(1
‑
δ
s
λ
t
)(1
‑
λ
s
)(1
‑
δ
h
λ
h
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(23)
[0117]
建立修正的等效塑性应变和等效塑性应变增量之间的关系：
[0118][0119]
将式(23)代入式(24)，得(考虑应变率)
[0120][0121]
假定
[0122][0123]
将式(25)代入式(26)，得
[0124][0125]
应变增量分解为弹性和塑性两部分：
[0126][0127]
将式(20)代入式(28)，得
[0128][0129]
式(29)两侧同时乘以弹性逆矩阵c
‑1[0130][0131]
式(30)两侧同时乘以屈服面对应力的梯度
[0132][0133]
应满足一致性条件(kuhn
‑
tucker准则)
[0134][0135]
其中，y,λ＝dy(σ
ij
,η)/dλ。
[0136]
将式(26)、(27)、(32)代入式(31)，得
[0137][0138]
注意到：
[0139][0140][0141][0142][0143]
[0144]
则：
[0145][0146]
将式(34)代入式(33)，得
[0147][0148]
式(35)提取公因式，移项，得
[0149][0150]
注意到，在时间步开始的时候泰勒展开另外，则一致性参数：
[0151][0152]
当前强度面更新
[0153]
y
n 1
＝y
*
y
,p
dp
p
y
,λ
dλ,dp
p
＝kωy
,p
du
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(38)
[0154]
更新应力及内变量
[0155][0156]
a2：基于fortran语言和visual studio进行编程并编译生成ls
‑
dyna有限元软件所需求解器。
[0157]
将上述理论推导过程通过fortran语言表达并编译成可适用于r5版本的ls
‑
dyna有限元软件的ls971.exe求解器，并将其拷贝至ls
‑
dyna软件安装目录下的program文件夹下。
[0158]
a3：基于本发明进行用户自定义工况数值计算并进行结果后处理。
[0159]
本发明适用于评估侵彻爆炸强动载作用下混凝土材料的损伤破坏，同样适用于混凝土结构的静载力学性能分析。首先，用户要采用国际单位制(m
‑
kg
‑
s)建立模型关键字文件(.k)，然后去更改混凝土材料模型及其参数，最后选择求解器并递交模型关键字文件进行求解即可，结果后处理不变。
[0160]
下表为本发明实施例中关键字参数输入格式，本发明实施方式已内嵌模量(包括弹性模量、体积模量和剪切模量)、拉伸强度计算公式和状态方程10组参数，本领域技术人员也可以根据经验输入。下表中前两行(mid～ig,ivect～ieos)参数选取已固定，其中，mid为材料模型编号，ro为材料密度，mt为用户自定义材料类型，lmc为材料参数数量，nhv为历史变量数量，iortho为定义材料各向同性和各向异性属性的参数，ibulk为体积模量参数位置，ig为剪切模量参数位置，ivect为用户提供向量化子程序的参数，ifail为控制壳单元和实体单元材料失效准则的参数，itherm为控制计算单元温度的参数，ihyper为控制计算变形梯度的参数，ieos为控制状态方程的参数。
[0161]
表1实施例中关键字参数输入格式
[0162][0163]
其中，f
c
：单轴抗压强度；e：弹性模量，k：体积模量，k＝e/3/(1
‑
2μ)；g：剪切模量，g＝e/2/(1 μ)；t：单轴抗拉强度，ω：塑性流动关联因子，0～1，建议0.5，可根据具体工况标定；frac_t：控制单元拉伸删除阈值，根据实验标定；frac_c：控制单元压缩删除阈值，本领域技术人员可根据实验标定；k1～k
10
：体积模量；μ1～μ
10
：体积应变；p1～p
10
：压力；rate：应变率效应开关，rate＜0时，不考虑应变率效应，rate≥0时，考虑应变率。
[0164]
(2)本发明的关键点
[0165]
1.引入三个强度面即屈服强度面、最大强度面和残余强度面，通过线性插值屈服强度面和最大强度面以达到描述混凝土材料应变硬化的目的，通过插值最大强度面和残余强度面以到达描述应变软化(损伤)的目的；
[0166]
2.根据已有三轴压缩试验数据显示屈服强度面在高围压下与压力关联不大，将高围压下的屈服强度面拉平，用以描述高围压下混凝土材料无明显屈服点的现象；
[0167]
3.引入基于微裂纹发展引起的拉伸损伤、剪切损伤和基于孔隙坍塌与压实引起的静水压缩损伤，健全了混凝土损伤机制的表达，并引入了拉伸损伤和剪切损伤之间的关联因子，同时在静水压缩损伤引入与应力状态相关的关联因子，以达到表征剪切损伤对拉伸损伤的影响和静水压缩损伤对其他应力状态影响的目的；
[0168]
4.为解决应力随应变率突变而突变的问题，考虑应变率时引入时间松弛因子，以达到反映应变率突变时应力迟滞效应的目的。
[0169]
图6是本发明技术方案，引入修正系数来修正高围压段屈服强度面形式；引入拉伸损伤指标、剪切损伤指标和静水压缩损伤指标以及剪切损伤对拉伸损伤影响因子和静水压缩损伤对其他应力状态影响因子，确定等效塑性应变；引入时间松弛因子确定当前时间步的结束时刻的等效应变率。
[0170]
图7是利用本发明进行的单单元三轴压缩计算结果，本实施例引入修正系数来修正高围压段屈服强度面形式，引入拉伸损伤指标、剪切损伤指标和静水压缩损伤指标以及剪切损伤对拉伸损伤影响因子和静水压缩损伤对其他应力状态影响因子，确定等效塑性应变，通过计算所体现的混凝土高围压下的力学性能。
[0171]
1是三轴压缩实验结果，2是利用本发明计算的混凝土三轴压缩下的应力应变曲线，3是轴向应力应变曲线，4是应力与径向应变曲线，5是应力与体积应变曲线。可以看出，计算结果与实验结果吻合较好，本发明能较好的考虑混凝土材料三轴力学性能。
[0172]
图8是利用本发明进行的单单元复杂应力状态计算结果，对应本实施例中引入拉伸损伤指标、剪切损伤指标和静水压缩损伤指标以及剪切损伤对拉伸损伤影响因子和静水压缩损伤对其他应力状态影响因子，通过计算所体现的混凝土高围压下静水压缩损伤累积对混凝土单轴受力状态的影响。复杂应力状态，即混凝土先承受静水压缩，然后释放静水压缩施加单轴压缩。1是低强度混凝土先静水压缩后单轴压缩的实验值，2是利用本发明计算的低强度混凝土先静水压缩后单轴压缩后与未受静水压缩后单轴压缩强度的比值，2是高强度混凝土先静水压缩后单轴压缩的实验值，4是利用本发明计算的高强度混凝土先静水压缩后单轴压缩后与未受静水压缩后单轴压缩强度的比值。可以看出，计算结果与实验值吻合较好，本发明能较好考虑静水压缩损伤对后续应力状态的影响，预测混凝土复杂应力状态性能。
[0173]
图9是利用本发明进行的单单元动态应变率计算结果，对应本实施例引入时间松弛因子确定当前时间步的结束时刻的等效应变率通过计算所体现的混凝土应变率突变引起的应力迟滞效应。1是应变率为5s
‑1单轴拉伸应力应变曲线，2是应变率为10s
‑1单轴拉伸应力应变曲线，3是应变率为20s
‑1单轴拉伸应力应变曲线，4是应变率由10s
‑1突变至5s
‑1单轴拉伸应力应变曲线，5是应变率由10s
‑1突变至20s
‑1单轴拉伸应力应变曲线。可以看出，本发明在计算应变率突变时，应力不发生突变，而是缓慢变化，能较好的反应实验过程中应变率突变引起的应力迟滞效应。
[0174]
图10是利用本发明进行的混凝土靶体弹体侵彻计算结果，本实施例采用以上技术方案，通过混凝土试件层次的侵彻计算，展示动态拉伸损伤、动态剪切损伤和动态静水压缩损伤的演化过程。图10(1)是利用本发明计算的有限厚度混凝土靶体侵彻结果，并注明了数
值计算和试验结果(1中括号内)对比，可以看出，本发明在计算混凝土靶侵彻过程中能很好地预测开坑和震塌现象；图10(2)是利用本发明计算的侵彻过程中靶体的拉伸损伤；图10(3)是利用本发明计算的侵彻过程中靶体的剪切损伤；图10(4)是利用本发明计算的侵彻过程中靶体的静水压缩损伤。可以看出，本发明能较好的预测侵彻过程靶体的损伤破坏情况。
[0175]
图11是利用本发明进行的混凝土靶体接触爆炸计算结果，本实施例采用以上技术方案，通过混凝土试件层次的爆炸计算，展示利用本发明计算混凝土遭受爆炸荷载作用下的损伤破坏情况。图11(1)是利用本发明针对混凝土薄板的锥形装药接触爆炸计算结果；图11(2)是混凝土薄板的锥形装药接触爆炸试验结果；图11(3)是利用本发明针对混凝土厚板的锥形装药接触爆炸计算结果；图11(4)是混凝土厚板的锥形装药接触爆炸试验结果。可以看出，本发明能较好的预测爆炸作用下混凝土毁伤破坏效应。
[0176]
以上结合附图对本发明的实施例进行了描述，但是本发明并不局限于上述的具体实施方式，上述的具体实施方式仅仅是示意性的，而不是限制性的，本领域的普通技术人员在本发明的启示下，在不脱离本发明宗旨和权利要求所保护的范围情况下，还可做出很多形式，这些均属于本发明的保护之内。