一种计算孔洞体积的方法
刘文辉，张新明，唐建国，周卓平，李芳芳
(中南大学 材料科学与工程学院，湖南 长沙，410083)
摘 要：
摘 要：将孔洞看成一种弹性模量和泊松比都非常小的材料，在有限元计算中，通过输出该材料的体积来计算孔洞的体积。研究结果表明：采用该方法能很好地模拟孔洞的体积变化和材料的变形，避免了积分法的大量编程运算，又可以解决给定位移法中由于体积变化恒定，不能考虑晶体取向对孔洞长大影响的难题，而且在计算中能直接求解出单胞的有效力学性能，因而，该方法能较好地用于计算与分析含孔洞材料的孔洞体积和有效力学性能。
关键词：
中图分类号：TG146 文献标识码：A 文章编号：1672-7207(2007)05-0815-05
A method for calculating volume of void
LIU Wen-hui, ZHANG Xin-ming, TANG Jian-guo, ZHOU Zhuo-ping, LI Fang-fang
(School of Materials Science and Engineering, Central South University, Changsha 410083, China)
Abstract: The void is treated as a material with very small Young’s modulus and Poisson’s ratio, and the volume of void can be obtained by outputting the material’s volume in calculation. The results show that this method can be used to predict the change of void volume and shape during deformation. It can avoid compiling programs and be used to calculate the volume of void in the simulation of the effect of crystallographic orientation on void growth by using crystal plasticity finite element method compared with integral method and prescribed displacement method, and it can be used to predict the effective mechanical properties, so it is very useful to calculate and analyze the volume of voids and the effective mechanical properties of materials with voids.
Key words: crystal plasticity; finite element method; void; mechanical properties; crystallographic orientation
材料的韧性断裂一般包含孔洞的形核、长大、聚合，是材料失效的主要形式之一。一般地，有关孔洞长大和聚合的有限元计算都采用1个包含1个或多个孔洞的单胞进行分析，通常把基体当成各向同性的材料或采用改进的本构模型来研究孔洞的长大和聚合。由于晶体有限元可以考虑材料的各向异性和各种微观结构，如晶体滑移、晶界以及晶粒之间的相互作用的情况，最近，不少研究者采用晶体有限元研究孔洞的长大和聚合。关于孔洞体积的计算，一般采用给定位移边界条件法，近似法或编制相关程序得到孔洞的体积[1]。J.Y. Shu[2]采用给定位移边界条件，运用粘弹塑性应变梯度晶体塑性理论研究了含孔洞的单晶的变形情况，发现小孔洞长大的趋势比大孔洞的小。T. Schacht[3]采用给定位移边界条件研究了含孔洞的单晶的变形情况，研究表明，孔洞的长大和变形情况与晶体的取向有关。G.P. Potirniche等[4]采用晶体有限元研究了单轴和双轴载荷情况下晶体取向对fcc单晶中孔洞的长大和聚合的影响，计算结果表明，在双轴加载的情况下，不同晶体取向的孔洞体积变化相差较小。实际上，在给定位移边界条件下，单胞的体积增长给定，由于基体体积变化很小，因此，晶体取向对孔洞的长大速率的影响也很小，难以分析在不同应力三轴度下孔洞体积的变化。为了研究不同加载条件下，孔洞体积的变化情况，本文作者将孔洞当成一种弹性模量和泊松比都非常小的材料，通过输出该材料的体积来计算变形中孔洞体积的变化。
1 率相关晶体塑性本构关系
R. Hill等[5-9]认为可以将位移梯度张量进行如下的乘积分解：
。 (1)
式中：Fe表示晶格畸变和刚性转动所产生的变形梯度；Fp则表示晶体沿滑移方向的均匀剪切所对应的变形梯度；速度梯度L也分解成弹性和塑性2部分：
。 (2)
式中：
；
。
速度梯度张量为所有滑移系的速度梯度张量的叠加：
； (3)
。 (4)
式中：a为第a个滑移系；
为 a滑移系中的滑移率；Dp和Wp 分别为变形率张量和旋转张量。
； (5)
。 (6)
式中：m代表垂直于滑移面的单位向量；s代表相应滑移面的滑移方向的单位向量。在面心立方体系中，滑移系为{111}<110>。
剪切应变率采用方程(7)的率相关的硬化方程来描述：
当
时，
的初值为
，它的演变由式(8)确定：
。 (8)
式中：
为硬化系数，通常当
时称为潜在硬化系数，本文采用了广泛使用的一种描述方式，即：
； (9)
。 (10)
式中：q为潜在硬化率，对于共面滑移系，q＝1.0；对于非共面滑移系，q＝1.0~1.4；h0和a为常数，ts为饱和应力。
2 结果分析
2.1 一般塑性有限元计算
当采用一般塑性有限元计算时，基体的形状变得不规则，难以计算孔洞的体积。这里把孔洞当成一种弹性模量和泊松比都非常小的材料，在计算中输出该材料的体积，并与采用给定位移边界条件法计算的孔洞体积进行对比。
在计算中，铝基体为各向同性的双线性弹塑性体，铝基体的弹性模量和泊松比分别为：Em=70 GPa，um=0.33；铝基体屈服应力和切向斜率分别为：Y0=89 MPa，E′=13.6 GPa；初始孔洞体积含量f=0.1，在平面应力条件下，单胞沿X方向位移为0.13，沿Y方向位移为-0.065。图1所示为孔洞和把孔洞当成一种弹性模量和泊松比都非常小材料的Von Mises等效应力图。从图1可看出，采用给定位移法和把空洞当成材料处理计算的孔洞形状和等效应力的差别很小。
(a) 孔洞；(b) 把孔洞当成材料处理
图1 等效应力分布图
Fig.1 Distribution of equivalent stress
图2和3所示分别为采用给定位移法和把空洞当成材料处理所得的A点Y方向的位移和B点X方向位移。
图2 采用给定位移法和把空洞当成材料处理所得的A点位移
Fig.2 Displacements of point A obtained by prescribed displacement methods and treated void as a material
图3 采用给定位移法和把空洞当成材料处理所得的B点位移
Fig.3 Displacement of point B obtained by prescribed displacement methods and treated void as a material
图4所示为采用给定位移法和把空洞当成材料处理的方法计算的孔洞体积分数的对比情况。从图1~4可以看出，用这2种方法计算所得结果差别很小，这说明把孔洞看成一种弹性模量和泊松比都非常小的材料的处理方法是可行的，因此，可采用该方法计算单晶在变形时孔洞的体积变化。
图4 采用给定位移法和把空洞当成材料处理所得的孔洞体积
Fig.4 Volumes of void obtained by prescribed displacement methods and treated void as a material
2.2 晶体有限元计算
在采用晶体有限元的计算中，晶体的本构关系通过有限元软件ANSYS的用户子程序实现，子程序采用Fortran语言编写[10-11]。计算材料为纯铝，其材料参数为：弹性常数C11＝170 GPa，C22＝124 GPa，C44＝75 GPa，参考剪切应变速率
=0.001，应变速率敏感因子m＝0.02；方程(9)和(10)中的硬化参数为h0＝250 MPa，t0＝16 MPa，ts＝190 MPa，a＝2.5[12]。
本文把孔洞看成一种弹性模量和泊松比都非常小的各向同性材料，来研究晶体取向对孔洞长大的影响。采用含球形孔洞的三维单胞来分析晶体取向对孔洞长大的影响(如图5所示)，图中X轴平行于样品轧向，Z轴平行于板材法线方向，拉伸方向平行于轧制方向(X轴)，孔洞初始体积分数 f = 1%。
图5 含球形孔洞的三维单胞结构
Fig.5 3D unit cell including a void
单胞的等效应力可通过对各单元的应力进行体积平均得到。图6所示为含孔洞和不含孔洞的Cu取向铝单晶的等效应力—应变曲线。
图6 Cu取向等效应力—应变曲线
Fig.6 Curves of equivalent stress vs strain for Cu orientation
图7所示分别为计算结果经MATLAB程序处理后的含孔洞和不含孔洞的Cu取向铝单晶在ex=0.2时的{111}极图。
(a)含孔洞；(b) 单晶
图7 Cu取向在ex=0.2的极图
Fig.7 Pole figure of Cu orientation at ex=0.2
图8所示为含孔洞和不含孔洞的G取向铝单晶的等效应力—应变曲线。由于单向拉伸时，应力三轴度较小，孔洞长大较慢，因此，含孔洞单胞的有效应力下降较小。
图8 G取向等效应力—应变曲线
Fig.8 Curves of equivalent stress vs strain for G orientation
图9所示分别为含孔洞和不含孔洞的G取向铝单晶在ex=0.2时的极图。从图7和图9可以看出，对于不含孔洞的情况，偏离初始取向较小；而对于含孔洞的情况，由于孔洞周围单元约束较少，变形较为自由，偏离初始取向较大。
(a) 含孔洞; (b) 单晶
图9 G取向在ex=0.2的极图
Fig.9 Pole figure of G orientation at ex=0.2
图10所示为不同取向的铝单晶在变形过程中孔洞体积分数的变化情况。从图10可以看出，在相同的轧向应变，Cube取向的孔洞长大速度最慢，Cu取向长大速度最快，G取向次之。这表明孔洞的长大速率与晶体取向密切相关，与G. P. Potirniche等[4]的研究结果一致。
图10 不同取向的铝单晶孔洞体积
Fig.10 Volume of void for different orientations of Al single crystal
3 结 论
把孔洞看成一种弹性模量和泊松比都非常小的材料，研究变形中材料的变形和孔洞的体积的变化。计算结果表明，采用该方法能较好地模拟材料的变形和孔洞体积的变化，确定在各种加载条件下的孔洞体积，为晶体有限元中计算孔洞的体积提供了一种有效的分析方法。
参考文献：
[1] Thompson C I A, Worswick M J, Pilkey A K, et al. Modeling void nucleation and growth within periodic clusters of particles[J]. J Mech Phys Solids, 1999, 47(1): 1-26.
[2] Shu J Y. Scale dependent deformation of porous single crystals[J]. International Journal of Plasticity, 1998, 14(10/11): 1085-1107.
[3] Schacht T, Untermann N, Steck E. The influence of crystallographic orientation on the deformation behavior of single crystals containing microvoids[J]. International Journal of Plasticity, 2003, 19(10): 1605-1626.
[4] Potirniche G P, Hearndon J L, Horstemeyer M F, et al. Lattice orientation effects on void growth and coalescence in fcc single crystals[J]. International Journal of Plasticity, 2006, 22(5): 921-942.
[5] Hill R. Generalized constitutive relations for incremental deformation of metals crystals by multislip[J]. Journal of the Mechanics and Physics of Solids, 1966, 14(2): 95-102.
[6] Hill R, Rice J R. Constitutive analysis of elastic-plastic crystals at arbitrary strain[J]. Journal of the Mechanics and Physics of Solids, 1972, 20(6): 401-413.
[7] Asaro R J. Crystal plasticity[J]. Transactions of the ASME, Journal of Applied Mechanics, 1983, 50(4B): 921-934.
[8] Asaro R J. Geometrical effects in the inhomogeneous deformation of ductile single crystals[J]. Acta Metallurgica, 1979, 27(3): 445-453.
[9] Asaro R J. Micromechanics of crystals and polycrystals[J]. Advances in Applied Mechanics, 1983, 23: 1-115.
[10] Tang J G, Zhang X M, Chen Z Y, et al. Simulation of rolling textures of FCC metals with crystal plasticity finite element model[J], Mater Sci Technol, 2006, 22(10): 1171-1176.
[11] Tang J G, Zhang X M, Chen Z Y, et al. Simulation influence of inhomogeneous deformation on the development of rolling textures with crystal plasticity finite element[J]. Journal of Central South University, 2006, 13(2): 117-121.
[12] 唐建国. FCC金属不均匀滑移变形的织构多晶体有限元模拟[D]. 长沙: 中南大学材料科学与工程学院, 2006.
TANG Jian-guo. Simulation of inhomogeneous slip deformation in FCC metals with textured polycrystal plasticity finite element[D]. Changsha: School of Material Science and Engineering, Central Sourth University, 2006.
收稿日期：2007-01-18；修回日期：2007-02-26
基金项目：国家重点基础研究发展规划项目(2005CB623706); 国际科技合作重点项目(2004DFA00200)
作者简介：刘文辉(1978-)，男，湖南攸县人，博士研究生，从事材料科学计算、力学及高强铝合金研究
通信作者：张新明，男，教授, 博士生导师; 电话：0731-8830265；E-mail: xmzhang @mail.csu.edu.cn
