基于ANSYS的复杂多孔材料数值模拟方法研究.pdf

发布时间：2022-05-30 发布人：admin 分类：说明书资料大小：1.55M 资料格式：pdf 举报版权申诉

weixin_38640072-12129807-4744300845407620911.pdf-第1页.png

第1页 / 共12页

weixin_38640072-12129807-4744300845407620911.pdf-第2页.png

第2页 / 共12页

weixin_38640072-12129807-4744300845407620911.pdf-第3页.png

第3页 / 共12页

weixin_38640072-12129807-4744300845407620911.pdf-第4页.png

第4页 / 共12页

weixin_38640072-12129807-4744300845407620911.pdf-第5页.png

第5页 / 共12页

weixin_38640072-12129807-4744300845407620911.pdf-第6页.png

第6页 / 共12页

weixin_38640072-12129807-4744300845407620911.pdf-第7页.png

第7页 / 共12页

weixin_38640072-12129807-4744300845407620911.pdf-第8页.png

第8页 / 共12页

文本预览

中国科技论文在线 http://www.paper.edu.cn 基于 ANSYS 的复杂多孔材料数值模拟方法研究1 彭瑞东 1,2，凌天龙 2，杨永明 1,2，郭宏伟 2，孔令峰 2 1 中国矿业大学煤炭资源与安全开采国家重点实验室，北京 (100083) 2 中国矿业大学（北京）力学与建筑工程学院，北京 (100083) E-mail:prd@cumtb.edu.cn 摘要：运用数值模拟的方法来研究复杂多孔材料的力学行为是一种比较新颖的方法。本文结合 ANSYS 有限元软件，对复杂多孔材料数值模拟所涉及的几何建模、单元划分和材料选择等共性问题进行了研究探讨。利用程序生成了具有某种统计分布特征的随机数序列，并将其作为孔隙半径及位置坐标建立了复杂多孔材料的几何模型；提出了一种针对复杂多孔材料几何模型的网格优化措施，实现了整体映射网格与局部自由网格相结合的划分方法，提高了生成单元的质量；选择合适的材料本构模型，借助单元生死技术，模拟了多孔材料的损伤破坏过程。研究指出多孔材料的力学行为不仅与孔隙率、孔隙位置等几何参数有关，也和基体的弹性模量、屈服强度等材料属性有关。关键词：多孔材料；数值模拟；随机分布；网格划分；单元生死中图分类号：O346.5; O242.21 1 引言多孔材料是存在于自然界的一大类固体材料。地质体中的岩石、生物体中的骨骼以及木材等天然材料都是多孔材料，此外还有用金属、陶瓷、高聚物等制成的各种蜂窝材料和泡沫材料都是人造多孔材料[1]。这些多孔材料的使用遍及国民经济生产的各部门，深入了解其力学性能特点具有重要意义。对多孔材料力学性能的研究可在宏观、细观、微观三个不同尺度层次下展开。宏观研究采用的是连续介质理论，将多孔材料视为某种平均化后的等效的连续介质，不关心孔隙结构的具体形态。由于多孔材料往往具有非线性弹性的特点，本构方程中需采用增量描述形式。而且多孔材料的弹性范围通常很小，超出这个范围，材料就会发生屈服、损伤或断裂。宏观研究忽略了孔隙结构的拓扑特性，因此无法解释多孔材料的变形和破坏机制，只能用于粗略的估算和大致的预测。微观研究是采用物理力学及物理化学的方法研究孔隙结构的形成过程，着重解决孔隙的形核和长大问题，为材料学研究人员所重视。但由于材料结构的复杂性，目前还很难将微观尺度上的物理量直接与宏观材料性能关联起来，尤其是二者之间的定量关系极难确定。介于宏微观之间的细观研究作为固体力学中振奋人心的新领域受到国际学术界及工程界的广泛注意，得到了生气蓬勃的发展。细观力学的发展使得固体材料变形破坏问题的研究突破了传统的均匀连续性假设，从而在空间尺度上更为真实地描述了研究对象。按照细观力学的研究思路，多孔材料的孔隙结构可用由单孔、孔壁、孔棱及结点构成的代表性体积单元（RVE）表示，进而通过统计方法或建立离散几何模型方法加以研究。研究表明，孔隙材料的力学特性往往取决于孔壁与孔棱的弯曲、塑性屈服、失稳坍塌和断裂等局部变形，这也就是说宏观变形对细观几何形貌极为敏感。Gibson 和 Ashby 在单一孔单元的基础上建立了模型，获得了蜂窝、开孔和闭孔泡沫的很多力学参数结果，如弹性变形、弹性屈曲、塑性塌陷、脆性断裂、黏弹性变形、蠕变及蠕变屈曲等[2]。Grenestedt、Christensen 等人给出了规则周期性多孔弹性体性能的分析模型[3,4]。Zhu、Kyriakides 等人研究了十四面体单元结 1 本课题得到教育部博士学科点新教师基金项目"复杂孔隙材料损伤破坏的细观机理研究"(20070290011)、国家自然科学基金项目“岩石类材料中孔隙结构损伤演化的实时观测与数值模拟研究”(10802092)的资助。 -1-

中国科技论文在线 http://www.paper.edu.cn 构的不规则性对材料性能的影响[5-7]。Roberts 分别研究了开孔和闭孔泡沫的的弹性性能[8,9]。 Gao 和 Li 对孔隙形状及孔壁厚度的影响进行了研究[10-11]。Forest 对代表性体积单元(RVE) 的尺寸大小进行了研究[12]。Gurson 模型在韧性损伤断裂分析中有着广泛应用，在此基础上许多改进的模型被用于孔隙材料的损伤分析[13-17]。从现在的研究中发现，大部分研究都采用了一定的数值模拟方法，这得益于数值分析理论和现代计算机技术的发展[18]。在多孔材料的力学行为研究方面，有限元的方法也被引入进来。Chung 和 Waas 用有限元对聚碳酸醋酯蜂窝材料在面内受压缩的情况进行了数值模拟 [19]。Chen 等人利用 ABAQUS 梁单元模拟了金属蜂窝材料的胞壁并对其力学响应进行了预测，借助有限元方法研究了刚性夹杂和孔洞对规则六边形蜂窝材料的弹性模量和屈服强度的影响[20,21]。Guo 和 Gibson 用有限元方法分析了胞元缺失对蜂窝结构的弹性模量、弹性屈曲、塑性屈服等的影响[22]。Andrews 和 Gibson 用有限元方法分析了缺陷尺寸和胞元尺寸对铝蜂窝含圆孔、中心裂纹、缺口裂纹试件的影响 [23]。Lu 和 Ruan 等人借助 ABAQUS 用有限元方法模拟了蜂窝铝的动力学行为[24]。这些研究有力的推动了对多孔材料变形破坏行为及机理的研究。但目前的有限元模拟研究还主要局限于对蜂窝材料等规则多孔材料的研究，因此有必要进一步考虑不规则多孔材料的数值模拟研究。本文将借助 ANSYS 有限元软件，探讨如何建立具有随机分布孔隙的多孔材料的数值计算模型，并对不同材料属性下的计算结果进行分析研究，从而为进一步系统研究具有复杂孔隙结构的多孔材料在数值模拟方面提供借鉴。本文研究主要解决复杂多孔材料数值模拟方面的一些共性问题，包括几何建模、单元划分和材料选择等，研究结果有助于揭示孔隙分布对多孔材料整体力学行为的影响。 2 复杂多孔材料的几何建模多孔材料是由固相和通过固相形成的孔隙所组成的复合体，它区别于普通密实固体材料的最显著特点是具有孔隙，因此孔隙率大小以及孔隙分布特点将对材料的物理力学性能产生显著影响。就复杂多孔材料而言，这些孔隙的形貌、尺寸以及分布都表现出很大的随机性。通过大量实验研究表明[25,26]，复杂多孔材料中孔隙率、孔径分布、孔隙位置分布等参数的变化是服从一定统计分布规律的。不失一般性，可以将孔隙视为圆球，其半径大小通常服从正态分布或指数分布，其球心坐标通常服从均匀分布。因此建立复杂多孔材料几何模型的关键是按照一定的概率分布规律确定孔隙的大小及位置，即求各种已知分布下的随机数。从数学方法上看，只要有了一种分布规律的随机数，可以通过各种数学变换或抽样方法，产生具有任意分布的随机数。为了快速、高精度的产生随机数。通常分为两步进行：首先，产生最简单的[0，1]区间上均匀分布的随机数。然后，在此基础上将其转换成所需各种分布的随机数。随着对随机数的不断研究和应用，已提出了各种数学方法来生成均匀分布的随机数，如取中法、加同余法、乘同余法、混合同余法等。在计算机上可以方便的产生随机变量 X 的 nX }为[0，1]区间上均匀分布的随机变量 X 的随机数。需要注意抽样序列{ 的是，随机数不是单指一个数而言，而是就一组数值序列而言的，因此为了确保每次程序运 nX }，通常称{ 行时得到的随机序列是不同的，需要初始化生成随机序列的种子，通常可用系统时间来生成该种子。几乎在所有的程序设计语言和软件脚本语言中都提供了产生均匀随机数的函数，初始化种子函数后利用该函数可产生[0，1]区间均匀分布的随机数。在此基础上，对于在［a，b］区间上均匀分布的随机数 Y，可通过变换 Y＝a+(b-a)X 得 -2-

中国科技论文在线 http://www.paper.edu.cn 到。而满足其它分布形式的随机数也可由在［0，1］区间上均匀分布的随机数进行相应的变换而得到，这通常用反函数法来实现。反函数法是产生任意分布随机数的最常用的方法之一，即将所需概率分布的随机数分布函数 F(x)进行反变换，得其反函数 )(1 yF− ，若 y 是在区间 [0，1]上均匀分布的随机数，则布的随机分布函数为： )(1 yF− 就是概率分布函数为 F(x)的伪随机数。例如指数分 (1) 则 F(x)是单调下降函数，而且其区间在 0 和 l 之间，因此其反函数也必定是单调函数， xF 1)( −= ≥ 0 x e ax , − 为： )(1 yFx = − −= 1 α 1ln( − y 0，) ≤ y ≤ 1 (2) 当 y 是在区间[0，1]上均匀分布的随机数时，那么产生的随机数序列{ nx }为服从指数分 /1 α 。于是利用公式(2)就可实现指数分布的随机数发布的随机数，其均值为 α/1 ，方差为 2 生器，在 MATLAB 语言中只要写入语句“ x =Rand(*,*)； y =-log( x )”就能得到指数分布的随机数 y，因为 1- x 与 x 服从同样的分布。至于生成正态分布等其它分布的方法与此类似。通过利用 FORTRAN 或 MATLAB 等进行简单的程序设计，即可按特定的概率统计规律生成所需的孔隙大小及位置坐标，进而利用 APDL 语言可以在 ANSYS 中将其作为球心坐标及半径建立一系列圆或球，再通过布尔运算将这些圆或球从矩形区域中减去即可得到复杂多孔材料的二维或三维几何模型。表 1 具有 25 个孔的复杂多孔材料几何模型孔径分布情况 Tab.1 The distribution of the radius of 25 pores in a geometrical model of complex porous materials 1mm 孔径大小孔隙数目 14 6mm 5mm 4mm 3mm 2mm 1 1 6 1 2 本文将以某二维模型为例进行讨论，模型为 200mm*100mm 矩形区域，孔隙位置服从均匀分布，孔径分布服从指数分布，根据实验室的实验情况取 6mm，5mm，4mm，3mm， 2mm，1mm 六种半径尺寸。为考虑孔隙率的影响，分别建立了 5 个多孔模型，其中孔隙数目分别为 10、15、20、25、30 个。根据生成的孔径分布，可以计算出 5 个模型的孔隙率，分别为 0.66%、0.93%、1.46%、2.09%和 2.68%。以具有 25 个孔隙的模型为例，首先利用 MATLAB 生成孔径 r 的随机数，各种孔径的孔隙数目如表 1 所示；然后接着利用 MATLAB 生成孔隙圆心坐标(x,y)的随机数，需要注意的是要避免各个孔隙相交，可在生成随机坐标位置后进行判断，如相交就再取一个随机数；根据孔隙圆心坐标(x,y)和孔径 r 的就可利用 APDL 语言在 ANSYS 中建立二维复杂多孔材料几何模型，如图 1 所示。图 1 孔隙率为 2.09%的复杂多孔材料几何模型 Fig1 the geometrical model of complex porous materials which porosity is 2.09% -3-

中国科技论文在线 3 复杂多孔材料模型的网格划分 http://www.paper.edu.cn 在建立了复杂多孔材料的几何模型之后，需要选择合适的方法将几何模型分割为一系列有限元单元的组合体，即进行网格划分，这包括单元类型的选择以及网格密度的确定。对于有限元分析来说，网格划分的好坏直接影响到求解的精度和速度。在单元类型的选择上，要注意剪力自锁、体积自锁、沙漏等问题，合理选择一次单元或二次单元，完全积分或缩减积分。在网格密度的确定上，要综合考虑求解范围的大小、应力集中的部位以及计算机存储量的大小和计算时间的长短等因素。网格越密，单元和节点的总数就越多，所需计算机的内存也就越大，相应的求解速度就越慢。反之，网格过疏，单元和节点的总数太少，尽管求解速度很快，但求解精度就较差。一般来说，网格数目增加，计算精度会有所提高，但同时计算的规模也会增加。通常应在不同部位采用不同疏密的网格划分，在应力变化比较剧烈、应力梯度较大的区域划分得细密一些，而在应力变化比较平缓的区域网格适当划分得稀疏一些，这样比较有利于减少计算的时间又不至于对精度有大的影响。另外在可能的情况下，应尽量选择规则的四边形或六面体单元，这有利于提高求解速度，改善收敛性。在 ANSYS 中提供了自由网格划分和映射网格划分两种方法。自由网格划分是自动化程度最高的网格划分技术之一，它可在任意平面或曲面上自动生成三角形或四边形网格，在任意体上自动生成四面体网格。通常情况下，可利用 ANSYS 的智能尺寸控制技术来自动控制网格的大小和疏密分布，也可进行人工设置网格的大小并控制疏密分布以及选择分网算法等。对于复杂几何模型而言，这种分网方法省时省力，但缺点是单元数量通常会很大，计算效率降低。映射网格划分是对规整模型的一种规整网格划分方法，对划分区域的拓扑结构有特定要求。映射面网格生成比较规则的四边形或三角形单元，映射体网格生成比较规则的六面体单元。对于由面经过拖拉、旋转、偏移等方式生成的复杂三维实体而言，可先在原始面上生成面网格，然后在生成体的同时自动形成三维实体网格，但别忘了在体网格划分完毕后清除面网格，也可用专门用于辅助网格划分的虚拟单元类型（MESH200）来划分面网格，之后就不需清除虚拟面网格了。对于已经形成好了的三维复杂实体，如果其在某个方向上的拓扑形式始终保持一致，则可用扫略网格划分功能来划分网格；这两种方式形成的单元几乎都是六面体单元。对于复杂几何体，通常难以全部划分为规则的四边形或六面体单元，这时可采用混合网格划分，即根据各部位的特点分别采用自由、映射、扫略等多种网格划分方式，以形成综合效果尽量好的有限元模型。混合网格划分方式要在计算精度、计算时间、建模工作量等方面进行综合考虑。通常，为了提高计算精度和减少计算时间，应首先考虑对适合于扫略和映射网格划分的区域先划分六面体网格，这种网格既可以是线性的（无中节点）、也可以是二次的（有中节点）；然后对必须用四面体自由网格划分的区域，采用带中节点的六面体单元进行自由分网，部分网格会自动退化成适合于自由划分形式的四面体单元，此时在该区域与已进行扫略或映射网格划分的区域的交界面上，会自动形成金字塔过渡单元。最后还应采用 TCHG 命令来将退化形式的四面体单元自动转换成非退化的四面体单元以提高求解效率。对于复杂多孔材料模型，由于各种孔隙的随机分布，使得几何模型十分复杂，因此难以进行映射网格划分，直接划分网格只能采用自由网格划分，这就特别容易产生尖锐或边长比很大的不良单元。图 2 所示为孔隙率为 0.93%和 2.09%的两种二维复杂多孔材料模型的自由网格划分结果。可以明显看出，由于结构中孔隙的存在（尤其是小径孔隙的存在），将导致孔隙之间的网格划分极不规则，而且随着孔隙率的增大，网格划分将越来越不规则。这样自 -4-

中国科技论文在线 http://www.paper.edu.cn 由划分的网格，不仅单元的数目多，而且单元质量差，尤其是在几何形状突变处易形成过于细密和畸形的网格，这不仅会增加计算时间，还会造成模型刚度和质量矩阵奇异，对整个模型的求解非常不利，尤其是在非线性计算中甚至会导致计算不收敛。因此有必要寻找一种改进的网格划分方法。 (a)孔隙率为 0.93% (b)孔隙率为 2.09% 图 2 复杂多孔材料模型的自由网格划分 Fig2 the free meshing for geometrical models of complex porous materials 考虑到复杂多孔材料的应力应变主要在孔隙附近发生较大变化，我们只希望在孔隙结构附近区域划分比较细的网格，而在孔隙结构以外的区域能有相对粗疏一些的网格来表示，这样就能降低计算规模。因此可建立比孔隙直径稍大的一个矩形或正六面体区域将孔隙包裹起来，该区域内孔隙的边上划分较密的网格，区域外围的边与材料基体相连接，划分为较疏的网格，该区域也正好反映了孔隙附近的应力应变集中区域的大小。将各个孔隙周围的这些区域去除后，基体就成为一个规则的区域，可以进行映射网格划分，生成质量较好的四边形或六面体单元。而这些区域内部采用自由网格划分，生成由疏到密逐渐过渡的三角形或四面体单元。以前述二维复杂多孔材料几何模型为例，模型的基体是一个矩形结构（如图 3a 所示），在建立孔隙结构之前，在有孔隙存在的位置建立边长比圆孔直径大的正方形，用布尔运算中 DIVIDE 或 SUBTRACT 命令将矩形结构划分为几个形状规则的面的组合，再利用 GLUE 命令将其粘接成为一个实体（如图 3b 所示），并在此模型基础上建立孔隙结构（如图 3c 所示）。这样，可通过定义每个边的份数实现对网格的控制，得到形状较好的四边形单元（如图 3d 所示），并且能在应力集中区域（孔隙周围）通过优化网格划分得到比较精确的计算结果。 (a) (b) (c) (d) 图 3 复杂多孔材料模型的网格划分方法 Fig3 the optimized method of meshing for geometrical models of complex porous materials -5-

中国科技论文在线 http://www.paper.edu.cn 上面针对二维多孔材料这种复杂的几何模型的网格划分方法提出了一个方便、可行的改进方案，该方法的思想亦可用于三维复杂多孔材料的网格划分。图 4 为采用该方法生成的孔隙率为 0.93%和 2.09%的两种二维复杂多孔材料模型的网格划分结果，对比图 2 可以看到，与前面将模型作为一个整体自动划分的的网格情况相比，改进后的方法划分的网格规整的多，而且单元质量比较好。这样不仅可以减少计算时间，也大大提高了最后计算结果的精度。通过加载计算对两种模型进行比较，前者计算时间较长，而且不能在弹塑性范围内进行求解，而后者则在弹性和弹塑性范围内都能够求解，而且求解速度较快，求解精度也比较高。这就说明这种改进的方案是可行的。 (a)孔隙率为 0.93% (b)孔隙率为 2.09% 图 4 复杂多孔材料模型的优化网格划分 Fig4 the optimized meshing for geometrical models of complex porous materials 4 复杂多孔材料模型的变形破坏特点本文旨在探讨复杂多孔材料的数值模拟方法，除了前述几何模型的建立、网格划分的改进，材料属性的选择也是一个重要问题。下面结合前述二维复杂多孔材料模型，从弹性、弹塑性乃至断裂三个方面对材料属性的选择以及多孔材料的变形破坏特点加以简单分析。本文研究的多孔材料模型为岩土材料，其弹性模量为 20GPa，泊松比为 0.27。为了模拟模型在单轴载荷作用下的力学响应，对模型左端施加对称位移约束，模型右端施加位移载荷。 4.1 弹性变形行为模拟与分析考虑材料为线弹性材料时，只需在 ANSYS 中设置材料的弹性模量和泊松比即可。右端施加载荷为-1mm，即压缩载荷。求解结果如图 5 所示。由弹性理论可知，在单轴压缩条件下，孔洞的左、右两端将出现一拉应力区，从图 5 所示第一主应力云图中可以看到这些拉应力区。另外可见，孔隙对基体绝大部分区域的影响不大，其影响主要体现为在孔隙附近出现了不同程度的应力集中情况。对于理想线弹性材料，尤其在外加载荷较小的情况下，这一应力集中不会影响到分析结果，但实际的材料总是会出现屈服或破坏的，屈服之后的变形或破坏将可能影响到材料的力学响应，因此在实际材料中是不可能产生这么高的应力的。下面将对材料发生塑性或损伤破坏后的情况加以分析。 -6-

中国科技论文在线 http://www.paper.edu.cn (a)无孔隙 (b)孔隙率为 0.66% (c)孔隙率为 0.93% (d)孔隙率为 1.46% (e)孔隙率为 2.09% (f)孔隙率为 2.68% 图 5 线弹性复杂多孔材料模型的第一主应力云图 Fig5 the first principal stress contour of linear elastic models of complex porous materials 4.2 弹塑性变形行为模拟与分析材料在塑性阶段的本构关系有理想塑性和塑性强化两种情况，其中塑性强化又分为等向强化和随动强化。下面考虑理想塑性和具有线性变化规律的等向强化（如图 6 所示）两种情况，因此设置材料属性屈服强度为 100MPa，强化后模量为 2GPa。右端施加载荷仍为-1mm，即压缩载荷。 σ σ ε (e) 理想弹塑性 (f) 弹塑性线性强化图 6 复杂多孔材料模型的塑性本构关系示意图 Fig6 the plastic constitutive curve of complex porous materials ε -7-

中国科技论文在线 http://www.paper.edu.cn 经过前述的网格优化措施，可以顺利完成弹塑性变形的计算。图 7 所示为 0.93%和 2.09% 两种孔隙率的模型在线弹性、线性等向塑性强化、理想线弹性材料属性下求解得到的 Mises 等效应力云图。对比图 7 所示结果可以看出，由于塑性变形的发生，孔隙周围的应力集中程度有一定缓解，线性等向塑性强化模型的应力集中程度低于线弹性模型，而理想弹塑性模型的应力集中程度又低于线性等向塑性强化模型。但计算结果也表明孔隙周围出现了应变集中的情况，有关详细讨论将另文发表。这也就说明，在多孔材料中，其力学性能不仅与孔隙率、孔隙位置等几何参数有关，也和弹性模量、屈服强度等材料的属性有关。对于不同的材料属性，孔隙集中的区域会出现应力集中或应变集中的现象。 (a) 线弹性，孔隙率为 0.93% (b) 线弹性，孔隙率为 2.09% (c) 线性等向塑性强化，孔隙率为 0.93% (d) 线性等向塑性强化，孔隙率为 2.09% (e) 理想弹塑性，孔隙率为 0.93% (f) 理想弹塑性，孔隙率为 2.09% 图 7 复杂多孔材料模型的 Mises 等效应力云图 Fig7 the Mises equivalent stress contour of models of complex porous materials 4.3 损伤断裂行为模拟与分析对于脆性材料而言，在出现应力集中后一般不会发生塑性变形而产生应变集中现象，也没有塑性强化现象，而是会出现损伤软化现象，即表现为在高应力区出现空洞等损伤，从而使材料弹模、强度等下降。为了模拟材料的损伤，在 ANSYS 有限元分析中需借助单元生死的技术。所谓“单元生死”是在模型中加入或删除材料属性，使模型中相应的单元或“生” 或“死”。要达到“单元死”的效果，并非将“杀死”的单元从模型中删除，而是将其刚度矩阵乘以一个很小的因子，为了减少求解的方程数和避免病态条件，还要将“死”单元的自 -8-

分享到：

赞收藏

资料库

基于ANSYS的复杂多孔材料数值模拟方法研究.pdf

相关推荐

开发技术

热门标签

最新资料