发明名称 一种用于复合材料结构失效有限元模拟中单元损伤耗散能量的估计方法
摘要 本发明涉及一种用于复合材料结构失效有限元模拟中单元损伤耗散能量的估计方法,首先建立了平面单元特征长度的计算模型,在此基础上建立了三维平行六面体单元的特征长度计算方法,再通过对一般形状单元和楔形单元进行等体积转换,将其特征长度计算问题转化为平行六面体单元的特征长度计算问题。本发明用于复合材料结构失效有限元模拟中单元损伤耗散能量的估计方法考虑了复合材料结构有限元模型中常用的典型单元构型,包括平行四边形单元,三角形单元,平行六面体单元,一般形状单元及楔形单元,建立了单元特征长度计算模型,提供了准确的单元特征长度计算公式,能够有效地计算单元特征长度,从而计算单元损伤的耗散能量,提高了复合材料结构失效有限元模拟的客观性和准确性。
申请公布号 CN103593567A 申请公布日期 2014.02.19
申请号 CN201310571786.8 申请日期 2013.11.13
申请人 北京航空航天大学 发明人 赵丽滨;秦田亮;山美娟;张建宇
分类号 G06F19/00(2011.01)I 主分类号 G06F19/00(2011.01)I
代理机构 北京科迪生专利代理有限责任公司 11251 代理人 杨学明;顾炜
主权项 1.一种用于复合材料结构失效有限元模拟中单元损伤耗散能量的估计方法,其特征在于,该估计方法的步骤如下:步骤一、裂纹带模型的本质是将材料损伤过程中耗散的总能量在整个单元体积即微元体内均匀化,使得有限元模拟过程中材料损伤耗散的能量恒定,而单元特征长度是裂纹带模型中的重要参数,它不仅与单元的几何形状有关,而且与材料的属性和裂纹的方向有关,准确地计算单元特征长度是保证准确地计算耗散能量的前提条件,复合材料结构失效有限元模拟中单元特征长度计算方法分为平面单元特征长度计算和三维单元特征长度计算,其中:平面单元特征长度计算步骤如下:设裂纹与单元一边的夹角为θ<sup>crc</sup>,裂纹长度为l<sub>crc</sub>,单元相邻两边的夹角为<img file="FDA0000414421720000011.GIF" wi="66" he="59" />边长分别为a和b,与AB边对应的距离为<img file="FDA0000414421720000012.GIF" wi="242" he="63" />步骤A1,首先假设一条裂纹穿过一层由n个平行四边形单元或三角形单元组成的区域;步骤B1,根据单元特征长度的定义,可得:<maths num="0001"><![CDATA[<math><mrow><msup><mi>l</mi><mo>*</mo></msup><mo>=</mo><mfrac><mi>nah</mi><msub><mi>l</mi><mi>crc</mi></msub></mfrac><mo>=</mo><mfrac><mi>nah</mi><mrow><mi>h</mi><mo>/</mo><mi>sin</mi><mrow><mo>(</mo><msup><mi>&theta;</mi><mi>crc</mi></msup><mo>)</mo></mrow></mrow></mfrac><mo>=</mo><mi>na</mi><mi>sin</mi><mrow><mo>(</mo><msup><mi>&theta;</mi><mi>crc</mi></msup><mo>)</mo></mrow><mo>;</mo></mrow></math>]]></maths>步骤C1,当<img file="FDA0000414421720000014.GIF" wi="509" he="65" />时,近似地认为n个平行四边形单元的面积为直角三角形AFE的2倍,<img file="FDA0000414421720000015.GIF" wi="194" he="71" />则上式可简化为:<img file="FDA0000414421720000016.GIF" wi="633" he="71" />步骤D1,当<img file="FDA0000414421720000017.GIF" wi="392" he="65" />时,平行四边形单元或三角形单元的特征长度为:<img file="FDA0000414421720000018.GIF" wi="403" he="71" />三维单元特征长度计算步骤如下:步骤A2,首先将六面体与纤维平行的平面按照平行四边形特征长度等效转化方法转换为两边与纤维平行的矩形E'F'G'H',然后将六面体的侧面四条边在平面E'F'G'H'的垂线上进行投影,得到一立方体;步骤B2,上述立方体A"B"C"D"E'F'G'H'各边的长度为:<maths num="0002"><![CDATA[<math><mrow><mover><mrow><msup><mi>E</mi><mo>&prime;</mo></msup><msup><mi>H</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover><mo>=</mo><mfenced open='{' close=''><mtable><mtr><mtd><mover><mi>EH</mi><mo>&OverBar;</mo></mover><mi>sin</mi><mrow><mo>(</mo><mo>&angle;</mo><mi>HEF</mi><mo>)</mo></mrow><mi>cos</mi><mrow><mo>(</mo><msup><mi>&theta;</mi><mi>f</mi></msup><mo>)</mo></mrow></mtd><mtd><mi>if</mi></mtd><mtd><mi>tan</mi><mrow><mo>(</mo><mo>&angle;</mo><mi>HEF</mi><mo>)</mo></mrow><mo>&le;</mo><mover><mi>EH</mi><mo>&OverBar;</mo></mover><mi>sin</mi><mrow><mo>(</mo><mo>&angle;</mo><mi>HEF</mi><mo>)</mo></mrow><mo>/</mo><mover><mi>EF</mi><mo>&OverBar;</mo></mover></mtd></mtr><mtr><mtd><mover><mi>EF</mi><mo>&OverBar;</mo></mover><mi>sin</mi><mrow><mo>(</mo><mo>&angle;</mo><mi>HEF</mi><mo>)</mo></mrow><mi>sin</mi><mrow><mo>(</mo><msup><mi>&theta;</mi><mi>f</mi></msup><mo>)</mo></mrow></mtd><mtd><mi>if</mi></mtd><mtd><mi>tan</mi><mrow><mo>(</mo><mo>&angle;</mo><mi>HEF</mi><mo>)</mo></mrow><mo>></mo><mover><mi>EH</mi><mo>&OverBar;</mo></mover><mi>sin</mi><mrow><mo>(</mo><mo>&angle;</mo><mi>HEF</mi><mo>)</mo></mrow><mo>/</mo><mover><mi>EF</mi><mo>&OverBar;</mo></mover></mtd></mtr></mtable></mfenced></mrow></math>]]></maths><maths num="0003"><![CDATA[<math><mrow><mover><mi>EF</mi><mo>&OverBar;</mo></mover><mo>=</mo><msub><mi>S</mi><mi>EFGH</mi></msub><mo>/</mo><mover><mrow><msup><mi>E</mi><mo>&prime;</mo></msup><msup><mi>H</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover></mrow></math>]]></maths><maths num="0004"><![CDATA[<math><mrow><mover><mrow><msup><mi>A</mi><mrow><mo>&prime;</mo><mo>&prime;</mo></mrow></msup><msup><mi>E</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover><mo>=</mo><mo>|</mo><mover><mi>AE</mi><mo>&RightArrow;</mo></mover><mo>&CenterDot;</mo><msub><mover><mi>e</mi><mo>&RightArrow;</mo></mover><mrow><mo>(</mo><mi>E</mi><mover><mi>F</mi><mo>&OverBar;</mo></mover><mo>&times;</mo><mi>E</mi><mover><mi>H</mi><mo>&OverBar;</mo></mover><mo>)</mo></mrow></msub><mo>|</mo></mrow></math>]]></maths>式中,<img file="FDA0000414421720000021.GIF" wi="379" he="64" />分别代表相应线段的长度,∠HEF、∠HEF代表角度,<img file="FDA0000414421720000022.GIF" wi="79" he="63" />代表用AE表示的向量,<img file="FDA0000414421720000023.GIF" wi="39" he="64" />代表单位向量;步骤C2,由立方体A"B"C"D"E'F'G'H'易知,对应于纤维“裂纹”的单元特征长度为:<maths num="0005"><![CDATA[<math><mrow><msubsup><mi>l</mi><mi>fbr</mi><mo>*</mo></msubsup><mo>=</mo><mover><mrow><msup><mi>E</mi><mo>&prime;</mo></msup><msup><mi>F</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover><mo>;</mo></mrow></math>]]></maths>步骤D2,设基体裂纹面与等效立方体单元的E'F'G'H'面所成夹角为<img file="FDA0000414421720000027.GIF" wi="97" he="57" />则对应于基体裂纹的单元特征长度为:<maths num="0006"><![CDATA[<math><mrow><msubsup><mi>l</mi><mi>mtr</mi><mo>*</mo></msubsup><mo>=</mo><mfenced open='{' close=''><mtable><mtr><mtd><mover><mrow><msup><mi>A</mi><mo>&prime;</mo></msup><msup><mi>E</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover><mi>cos</mi><mrow><mo>(</mo><msubsup><mi>&theta;</mi><mi>m</mi><mi>crc</mi></msubsup><mo>)</mo></mrow></mtd><mtd><mi>if</mi></mtd><mtd><mi>tan</mi><mrow><mo>(</mo><msubsup><mi>&theta;</mi><mi>m</mi><mi>crc</mi></msubsup><mo>)</mo></mrow><mo>&le;</mo><mover><mrow><msup><mi>A</mi><mrow><mo>&prime;</mo><mo>&prime;</mo></mrow></msup><msup><mi>E</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover><mo>/</mo><mover><mrow><msup><mi>E</mi><mo>&prime;</mo></msup><msup><mi>H</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover></mtd></mtr><mtr><mtd><mover><mrow><msup><mi>E</mi><mo>&prime;</mo></msup><msup><mi>H</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover><mi>sin</mi><mrow><mo>(</mo><msubsup><mi>&theta;</mi><mi>m</mi><mi>crc</mi></msubsup><mo>)</mo></mrow></mtd><mtd><mi>if</mi></mtd><mtd><mi>tan</mi><mrow><mo>(</mo><msubsup><mi>&theta;</mi><mi>m</mi><mi>crc</mi></msubsup><mo>)</mo></mrow><mover><mrow><msup><mi>A</mi><mrow><mo>&prime;</mo><mo>&prime;</mo></mrow></msup><msup><mi>E</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover><mo>/</mo><mover><mrow><msup><mi>E</mi><mo>&prime;</mo></msup><msup><mi>H</mi><mo>&prime;</mo></msup></mrow><mo>&OverBar;</mo></mover></mtd></mtr></mtable></mfenced></mrow></math>]]></maths>其中基体裂纹面与立方体单元边的相对位置可以通过裂纹面在模型整体坐标系中的方向坐标和平面E′F′H′G′的法向坐标求得;步骤E2,对于形状更为一般的六面体单元,可以首先将其等体积转化为平行六面体单元,其中平行六面体单元的底面由四边形EFGH等面积转化来确定,进而采用平行六面体的计算公式求得单元的特征长度;步骤F2,对于楔形单元,可以借鉴平面三角形单元特征长度的求解方法,将两个相同楔形单元组合、扩展为平行六面体后得到其单元特征长度;步骤二、复合材料结构失效有限元模拟中单元耗散能量计算步骤如下:步骤A3,在损伤生成阶段,微裂纹分散在整个单元体积内,随着载荷的增大,微裂纹生长、贯通,最终形成宏观裂纹。因此,单元在损伤即微裂纹形成、扩展过程中耗散的总能量W<sub>D</sub>就是细观裂纹形成消耗的断裂能W<sub>C</sub>,即:W<sub>D</sub>=W<sub>C</sub>;步骤B3,假设材料的临界应变能释放率G<sub>C</sub>恒定,则在单元内材料生成面积为A的裂纹消耗的断裂能W<sub>C</sub>为:W<sub>C</sub>=G<sub>C</sub>*A;步骤C3,基于上述单元特征长度计算方法得到的特征长度,由W<sub>D</sub>=W<sub>C</sub>,得单元的损伤耗散能密度为:<maths num="0007"><![CDATA[<math><mrow><mi>&Gamma;</mi><mo>=</mo><mfrac><msub><mi>G</mi><mi>c</mi></msub><mrow><mi>V</mi><mo>/</mo><mi>A</mi></mrow></mfrac><mo>=</mo><mfrac><msub><mi>G</mi><mi>c</mi></msub><msup><mi>l</mi><mo>*</mo></msup></mfrac></mrow></math>]]></maths>其中,G<sub>C</sub>为材料的临界应变能释放率,V为单元的体积,A为与宏观裂纹平行的截面积,l<sup>*</sup>为单元的特征长度,由上述单元特征长度计算方法得到;步骤D3,基于上面得到的单元损伤耗散能密度,最终得到从微裂纹产生、扩展到宏观裂纹形成过程中,单元损伤所耗散的能量W<sub>D</sub>为:W<sub>D</sub>=Γ*V。
地址 100191 北京市海淀区学院路37号