CN112100762B - 一种陶瓷基复合材料超单元结构数值仿真计算方法 - Google Patents
一种陶瓷基复合材料超单元结构数值仿真计算方法 Download PDFInfo
- Publication number
- CN112100762B CN112100762B CN202010875363.5A CN202010875363A CN112100762B CN 112100762 B CN112100762 B CN 112100762B CN 202010875363 A CN202010875363 A CN 202010875363A CN 112100762 B CN112100762 B CN 112100762B
- Authority
- CN
- China
- Prior art keywords
- superunit
- axis
- strain
- calculating
- matrix
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000004364 calculation method Methods 0.000 title claims abstract description 40
- 239000000463 material Substances 0.000 title claims abstract description 36
- 238000004088 simulation Methods 0.000 title claims abstract description 34
- 239000011153 ceramic matrix composite Substances 0.000 title claims abstract description 26
- 239000011159 matrix material Substances 0.000 claims abstract description 49
- 238000000034 method Methods 0.000 claims abstract description 30
- 239000002131 composite material Substances 0.000 claims abstract description 12
- 230000000750 progressive effect Effects 0.000 claims abstract description 6
- 238000006073 displacement reaction Methods 0.000 claims description 61
- 239000013598 vector Substances 0.000 description 7
- 238000004458 analytical method Methods 0.000 description 3
- 239000000835 fiber Substances 0.000 description 3
- 239000011204 carbon fibre-reinforced silicon carbide Substances 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 230000006978 adaptation Effects 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 239000000919 ceramic Substances 0.000 description 1
- 239000003086 colorant Substances 0.000 description 1
- 238000010835 comparative analysis Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000001360 synchronised effect Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/17—Mechanical parametric or variational design
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/26—Composites
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Geometry (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Evolutionary Computation (AREA)
- Computer Hardware Design (AREA)
- General Engineering & Computer Science (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- Computational Mathematics (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
一种陶瓷基复合材料超单元结构数值仿真计算方法,包括如下步骤:步骤一、对结构进行超单元网格划分;步骤二、根据超单元的应变场计算超单元内部不同方向各个单胞的应变场;步骤三、根据单胞模型的材料参数计算出在此应变状态下的单胞的刚度矩阵;步骤四、根据复合材料刚度矩阵转角公式计算出超单元内部所有单胞在超单元坐标系下的刚度矩阵,再按体积平均法计算出超单元的刚度矩阵;步骤五、基于渐进损伤方法进行有限元数值仿真计算。本发明考虑了结构单胞尺度的细观结构,可以计算陶瓷基复合材料铺层不平行的情况,突破了现有方法无法对此类情况进行计算的技术瓶颈。
Description
技术领域
本发明属于复合材料有限元数值仿真计算领域,具体涉及一种用超单元划分网格并进行数值仿真计算的方法。
背景技术
陶瓷基复合材料(CMCs)先由纤维丝组成纤维束,纤维束编织成预制体,最后经基体沉积,表现为典型的各向异性非线性材料。当今对陶瓷基复合材料结构进行有限元数值仿真计算时,均简化为宏观均匀的各向异性材料,不考虑结构的细观结构。
一些学者通过力学试验来获取陶瓷基复合材料的弹性参数,将其作为结构在有限元模型中的单元材料属性,如文献[王姣,王成华,杨阳,等.C/SiC复合材料轴套渐进损伤分析[J].强度与环境,2016,43(05):30-37][卢子兴,廖强,杨振宇,等.C/SiC复合材料螺栓螺牙承载能力[J].复合材料学报,2015,32(01):182-187]。这种做法认为CMCs结构是宏观均匀的,即各处材料性能相同,没有考虑结构的细观结构。
另一部分学者建立CMCs的单胞模型,基于细观力学方法计算出单胞模型的弹性参数,在进行结构件的有限元仿真计算时将单胞的弹性参数作为单元的材料属性,如文献[张盛.编织陶瓷基复合材料力学行为的多尺度分析[D].南京:南京航空航天大学,2018][杨福树.2.5维编织陶瓷基复合材料疲劳行为研究[D].南京:南京航空航天大学,2011]。该方法只是通过细观力学方法得到材料弹性参数,实际进行结构的有限元计算时仍简化为宏观均匀的各向异性材料,没有考虑结构实际的细观结构。这两种方法在有限元仿真计算时都是将材料简化为宏观均匀的各向异性材料,区别在于得到有限元网格材料参数的途径不同。
当今的技术水平可以通过高精度仪器和先进测量技术绘制出结构件的细观尺度模型,例如基于同步辐射光源和微焦点XCT三维重构技术建立三维模型,如文献[于国强.陶瓷基纤维束复合材料各向异性力学模型与结构失效模拟[D].南京:南京航空航天大学,2020],但目前还没有考虑细观结构进行有限元仿真计算的研究。现有的方法都将材料简化为宏观均匀的各向异性材料,没有考虑结构的细观模型,必然会存在一定的误差,考虑细观结构进行有限元计算是必然的发展趋势。
现有方法只能计算结构件铺层平行规整的情况,只有当各个铺层平行,不同铺层的材料主方向才是相同的,才能简化为宏观均匀的各向异性材料。现有方法无法对不同铺层间铺层方向不同的结构进行数值仿真计算,由于实际结构复杂或根据使用需求,导致铺层间纤维束(纱线)并不平行,对此类结构进行有限元仿真计算已成为CMCs结构设计和力学分析的技术瓶颈。
综上:当今对陶瓷基复合材料进行数值仿真计算均是将其简化为宏观均匀的各向异性材料,没有考虑结构件的细观结构;且只能计算铺层平行规整的结构件,有必要对不同铺层间铺层方向不同的情况进行研究。
发明内容
本发明针对现有技术中的不足,为实现对不同铺层间铺层方向不同的结构进行仿真运算,提供一种陶瓷基复合材料超单元结构数值仿真计算方法,用超单元划分结构并进行有限元数值仿真。超单元尺度大于结构的单胞尺度,所以超单元内部包含多个单胞,且由于结构铺层不平行,超单元包含材料主方向不同的单胞。本发明基于超单元网格内部的单胞模型计算超单元的刚度矩阵,然后进行有限元数值仿真,本方法仿真计算过程中使用了实际结构的单胞模型,考虑了单胞尺度的细观结构。
为实现上述目的,本发明采用以下技术方案:
一种陶瓷基复合材料超单元结构数值仿真计算方法,其特征在于,包括如下步骤:
步骤一:对陶瓷基复合材料结构进行超单元网格划分;
步骤二:根据超单元的应变场计算超单元内部不同方向各个单胞的应变场;
步骤三:根据步骤二计算的各个单胞的应变场,计算出单胞模型的材料参数,进而计算出在此应变状态下的单胞刚度矩阵;
步骤四:基于步骤三计算的单胞刚度矩阵,结合复合材料刚度矩阵转角公式计算出超单元内部所有单胞在超单元坐标系下的刚度矩阵,再按体积平均法计算出超单元刚度矩阵;
步骤五:根据步骤四计算的超单元刚度矩阵,进行有限元数值仿真计算。
为优化上述技术方案,采取的具体措施还包括:
进一步地,所述步骤二中,计算各个单胞的应变场的具体步骤如下:
设超单元坐标系为Oxyz,超单元坐标系为有限元计算中的单元坐标系,超单元的应变场为εx、εy、εz、γxy、γxz、γyz,εx为x坐标轴方向的线应变,εy为y坐标轴方向的线应变,εz为z坐标轴方向的线应变,γxy为xy坐标轴构成的平面方向的切应变,γxz为xz坐标轴构成的平面方向的切应变,γyz为yz坐标轴构成的平面方向的切应变;
内部单胞坐标系为Ox'y'z',单胞坐标系坐标轴的方向为单胞材料主方向,内部单胞的应变场为εx'、εy'、εz'、γx'y'、γx'z'、γy'z',εx'为x'坐标轴方向的线应变,εy'为y'坐标轴方向的线应变,εz'为z'坐标轴方向的线应变,γx'y'为x'y'坐标轴构成的平面方向的切应变,γx'z'为x'z'坐标轴构成的平面方向的切应变,γy'z'为y'z'坐标轴构成的平面方向的切应变;
两单元坐标系夹角余弦值如表1所示。其中超单元坐标系为有限元计算中的单元坐标系,单胞坐标系坐标轴的方向为单胞材料主方向;
表1超单元坐标系与单胞坐标系的夹角余弦值
x | y | z | |
x' | l1 | m1 | n1 |
y' | l2 | m2 | n2 |
z' | l3 | m3 | n3 |
采用如下公式计算单胞应变场:
εx'=εxl1 2+εym1 2+εzn1 2+2γxyl1m1+2γxzl1n1+2γyzm1n1
εy'=εxl2 2+εym2 2+εzn2 2+2γxyl2m2+2γxzl2n2+2γyzm2n2
εz'=εxl3 2+εym3 2+εzn3 2+2γxyl3m3+2γxzl3n3+2γyzm3n3
γx'y'=εxl1l2+εym1m2+εzn1n2+γxy(l1m2+m1l2)+γxz(l1n2+n1l2)+γyz(m1n2+n1m2)
γx'z'=εxl1l3+εym1m3+εzn1n3+γxy(l1m3+m1l3)+γxz(l1n3+n1l3)+γyz(m1n3+n1m3)
γy'z'=εxl2l3+εym2m3+εzn2n3+γxy(l2m3+m2l3)+γxz(l2n3+n2l3)+γyz(m2n3+n2m3)
式中,l1是x'轴和x轴之间夹角的余弦值,m1是x'轴和y轴之间夹角的余弦值,n1是x'轴和z轴之间夹角的余弦值,l2是y'轴和x轴之间夹角的余弦值,m2是y'轴和y轴之间夹角的余弦值,n2是y'轴和z轴之间夹角的余弦值,l3是z'轴和x轴之间夹角的余弦值,m3是z'轴和y轴之间夹角的余弦值,n3是z'轴和z轴之间夹角的余弦值。
进一步地,所述步骤三中,根据各个单胞的应变场计算单胞刚度矩阵的具体步骤如下:
陶瓷基复合材料的应力应变曲线为非线性,根据应力应变曲线计算出在此应变状态下材料的割线模量;在不同应变状态下材料的弹性参数是不一样的,根据εx'、εy'、εz'、εx'y'、εx'z'、εy'z'计算出弹性模量Ex、Ey、Ez,Ex为x方向的弹性模量,Ey为y方向的弹性模量,Ez为z方向的弹性模量;
然后采用下式计算单胞刚度矩阵[C]:
μyx=μxyEy/Ex
μzx=μxzEz/Ex
μzy=μyzEz/Ey
Δ=(1-μxyμyx-μxzμzx-μyzμzy-2μxyμyzμzx)/(ExEyEz)
C11=(1-μyzμzy)/(EyEzΔ)
C22=(1-μxzμzx)/(ExEzΔ)
C33=(1-μxyμyx)/(ExEyΔ)
C12=C21=(μyx+μzxμyz)/(EyEzΔ)
C13=C31=(μxz+μxyμyz)/(ExEyΔ)
C23=C32=(μzy+μxyμzx)/(ExEzΔ)
C44=Gxy
C55=Gxz
C66=Gyz
式中,μxy为x对y方向的次泊松比(在单轴作用下,y方向的单位拉或压应变所引起的x方向的压或拉应变),μyx为y对x方向的次泊松比,μzx为z对x方向的次泊松比,μxz为x对z方向的次泊松比,μzy为z对y方向的次泊松比,μyz为y对z方向的次泊松比,Gxy为xy方向的剪切模量,Gxz为xz方向的剪切模量,Gyz为yz方向的剪切模量,Cij,i=1...6,j=1...6是刚度矩阵[C]中的元素。
进一步地,所述步骤四中,计算超单元刚度矩阵的具体步骤如下:
单胞刚度矩阵在超单元坐标系Oxyz和单胞坐标系Ox'y'z'分别为[C]xyz和[C]x'y'z',单胞刚度矩阵在超单元坐标系下通过如下公式进行计算:
[C]xyz=[T][C]x'y'z'[T]T
式中,[C]x'y'z'等于[C],矩阵[T]是根据两坐标系各轴线之间的夹角余弦值计算所得:
计算出超单元内部所有单胞在超单元坐标系下的刚度矩阵,然后按照各个单胞占超单元体积比的方法计算出超单元的刚度矩阵[K]:
式中,n为超单元内部的单胞数量,V单胞i为第i个单胞的体积。
进一步地,所述步骤五中,陶瓷基复合材料为非线性材料,在进行有限元仿真计算时需用渐进损伤方法刚度折减进行迭代计算;基于渐进损伤方法的有限元仿真计算用于计算应力应变曲线为非线性的材料,前文计算的超单元刚度矩阵是得到单元网格的刚度矩阵,是进行有限元仿真计算的必要条件。在进行有限元仿真计算时使用渐进损伤方法刚度折减进行迭代计算,并输出计算仿真结果进行应力应变位移云图显示。
本发明的有益效果是:
1、考虑了结构单胞尺度的细观结构;
2、可以计算陶瓷基复合材料铺层不平行的情况,突破了现有方法无法对此类情况进行计算的技术瓶颈。
附图说明
图1是超单元有限元计算程序的流程示意图。
图2是超单元和单胞的关系图。
图3是有限元软件和超单元有限元计算软件的总体位移云图对照图。
图4是有限元软件和超单元有限元计算软件的x方向位移云图对照图。
图5是有限元软件和超单元有限元计算软件的y方向位移云图对照图。
图6是有限元软件和超单元有限元计算软件的z方向位移云图对照图。
具体实施方式
现在结合附图对本发明作进一步详细的说明。
如图1所示的一种陶瓷基复合材料超单元结构数值仿真计算方法,包括以下内容:
一、对结构进行超单元网格划分,可在成熟有限元商业软件里进行,将划分好的几何模型和超单元模型导入到超单元有限元计算程序。其中超单元和单胞的关系如图2所示,两个不同的颜色为两个材料主方向不平行的铺层,超单元包含方向不同的单胞和不完整单胞。
二、材料为各向异性非线性,有限元仿真计算过程中需进行迭代计算。在初次迭代时,应力应变场为0,此时超单元的刚度矩阵根据初始弹性模量计算。
三、判断超单元内部含有哪些单胞。在计算某个超单元刚度矩阵时先选定这个超单元网格,遍历所有单胞得到此超单元具体包含的单胞。
判断某个单胞是否在超单元内部的算法如下。由超单元模型可以得到平面内两个不共线的向量,利用向量叉乘公式可以得到该平面的法向量,值得注意的是两向量的叉乘顺序,两向量应按照右手法则叉乘使得到的法向向量指向超单元外侧。得到平面的法向向量即得到平面方程a,b,c的系数,代入平面任一点得到空间平面方程,按照此种方法可求得超单元所有平面的空间平面方程。将需要判断的空间一点代入所有的平面方程,如果此点在超单元内部,那么代入求得的所有空间平面方程的值都将小于0。
四、按照前文所提及的应变场转轴公式基于超单元应变场计算超单元内部所有单胞的应变场。
五、计算各个单胞在当前应变状态下的割线模量,然后按照前文提及的复合材料力学基本公式计算单胞的刚度矩阵。
六、按照前文提及的复合材料刚度矩阵的转轴公式计算出各个单胞在超单元坐标系下的刚度矩阵。
七、计算所有单胞的体积,然后按照体积平均法计算超单元的刚度矩阵。
八、按照上述步骤可以计算出所有超单元的刚度矩阵,然后进行有限元求解计算,对比此次求解结果和上次结果的误差,迭代直至误差在精度要求内。
九、输出计算仿真结果到有限元商业软件里进行应力应变位移云图显示。
下面以一个12层平纹编织铺层的带孔板为例,其中奇数层和偶数层铺层方向有30°夹角。该模型分别在超单元有限元计算程序和有限元商业软件里进行数值仿真计算,验证本发明的可行性。
该带孔板模型在超单元有限元计算程序的具体流程包括以下步骤:
一、在有限元商业软件里对模型进行超单元有限元网格划分,将几何模型和超单元模型导入到超单元有限元计算程序。
二、遍历所有单胞判断超单元具体包含哪些单胞,具体方法前文中有提及。
三、根据应变场转轴公式基于超单元应变场计算超单元内部所有单胞的应变场。
四、计算各个单胞在当前应变状态下的割线模量,然后按照复合材料力学基本公式计算单胞的刚度矩阵。
五、根据复合材料刚度矩阵的转轴公式计算出各个单胞在超单元坐标系下的刚度矩阵。
六、计算所有单胞的体积,然后按照体积平均法计算超单元的刚度矩阵。
七、进行有限元求解计算,迭代计算直至能误差在精度要求范围内。
八、输出计算仿真结果到有限元商业软件里进行位移云图显示。
两者计算出的结果显示为位移云图如图3、图4、图5、图6所示
1、有限元仿真软件总体位移云图最大位移为0.002769mm,最小位移为0;超单元有限元计算程序总体位移云图最大位移为0.002737mm,最小位移为0。总体位移的最大位移误差为1.156%,总体位移的最小位移无误差。两者总体位移场基本一致。
2、有限元仿真软件计算的x方向位移云图最大位移为0.002769mm,最小位移为0;超单元有限元计算程序x方向位移云图最大位移为0.002737mm,最小位移为0。X方向最大位移误差为1.156%,X方向最小位移无误差。两者X方向位移场一致。
3、有限元仿真软件计算的y方向位移云图最大位移为0.000404mm,最小位移为-0.000404mm;超单元有限元计算程序y方向位移云图最大位移为0.000385mm,最小位移为-0.000385mm。Y方向最大位移和最小位移误差均为0.686%。两者Y方向位移场一致。
4、有限元仿真软件计算的z方向位移云图最大位移为0.0000898mm,最小位移为-0.0000898mm;超单元有限元计算程序y方向位移云图最大位移为0.0000780mm,最小位移为-0.0000780mm。Z方向最大位移和最小位移误差均为0.426%。
值得注意的是,两者z方向位移场看着不一致,这是因为带孔板的位移主要在x方向。z方向的位移实际位移差很小,只有x方向位移差的37%,只是z方向最大位移较小才导致位移云图看起来误差很大。
经过对总体位移场、X方向位移场、Y方向位移场和Z方向位移场的对比分析,发现总体位移场、X方向位移场、Y方向位移和Z方向位移场的误差很小。由此可见,本发明是可行的。
以上仅是本发明的优选实施方式,本发明的保护范围并不仅局限于上述实施例,凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出,对于本技术领域的普通技术人员来说,在不脱离本发明原理前提下的若干改进和润饰,应视为本发明的保护范围。
Claims (4)
1.一种陶瓷基复合材料超单元结构数值仿真计算方法,其特征在于,包括如下步骤:
步骤一:对陶瓷基复合材料结构进行超单元网格划分;
步骤二:根据超单元的应变场计算超单元内部不同方向各个单胞的应变场;
步骤三:根据步骤二计算的各个单胞的应变场,计算出单胞模型的材料参数,进而计算出在此应变状态下的单胞刚度矩阵;
步骤四:基于步骤三计算的单胞刚度矩阵,结合复合材料刚度矩阵转角公式计算出超单元内部所有单胞在超单元坐标系下的刚度矩阵,再按体积平均法计算出超单元刚度矩阵;
步骤五:根据步骤四计算的超单元刚度矩阵,进行有限元数值仿真计算;
所述步骤二中,计算各个单胞的应变场的具体步骤如下:
设超单元坐标系为Oxyz,超单元坐标系为有限元计算中的单元坐标系,超单元的应变场为εx、εy、εz、γxy、γxz、γyz,εx为x坐标轴方向的线应变,εy为y坐标轴方向的线应变,εz为z坐标轴方向的线应变,γxy为xy坐标轴构成的平面方向的切应变,γxz为xz坐标轴构成的平面方向的切应变,γyz为yz坐标轴构成的平面方向的切应变;
内部单胞坐标系为Ox'y'z',单胞坐标系坐标轴的方向为单胞材料主方向,内部单胞的应变场为εx'、εy'、εz'、γx'y'、γx'z'、γy'z',εx'为x'坐标轴方向的线应变,εy'为y'坐标轴方向的线应变,εz'为z'坐标轴方向的线应变,γx'y'为x'y'坐标轴构成的平面方向的切应变,γx'z'为x'z'坐标轴构成的平面方向的切应变,γy'z'为y'z'坐标轴构成的平面方向的切应变;
采用如下公式计算单胞应变场:
εx'=εxl1 2+εym1 2+εzn1 2+2γxyl1m1+2γxzl1n1+2γyzm1n1
εy'=εxl2 2+εym2 2+εzn2 2+2γxyl2m2+2γxzl2n2+2γyzm2n2
εz'=εxl3 2+εym3 2+εzn3 2+2γxyl3m3+2γxzl3n3+2γyzm3n3
γx'y'=εxl1l2+εym1m2+εzn1n2+γxy(l1m2+m1l2)+γxz(l1n2+n1l2)+γyz(m1n2+n1m2)
γx'z'=εxl1l3+εym1m3+εzn1n3+γxy(l1m3+m1l3)+γxz(l1n3+n1l3)+γyz(m1n3+n1m3)
γy'z'=εxl2l3+εym2m3+εzn2n3+γxy(l2m3+m2l3)+γxz(l2n3+n2l3)+γyz(m2n3+n2m3)
式中,l1是x'轴和x轴之间夹角的余弦值,m1是x'轴和y轴之间夹角的余弦值,n1是x'轴和z轴之间夹角的余弦值,l2是y'轴和x轴之间夹角的余弦值,m2是y'轴和y轴之间夹角的余弦值,n2是y'轴和z轴之间夹角的余弦值,l3是z'轴和x轴之间夹角的余弦值,m3是z'轴和y轴之间夹角的余弦值,n3是z'轴和z轴之间夹角的余弦值。
2.如权利要求1所述的一种陶瓷基复合材料超单元结构数值仿真计算方法,其特征在于:所述步骤三中,根据各个单胞的应变场计算单胞刚度矩阵的具体步骤如下:
根据εx'、εy'、εz'、εx'y'、εx'z'、εy'z'计算出弹性模量Ex、Ey、Ez,Ex为x方向的弹性模量,Ey为y方向的弹性模量,Ez为z方向的弹性模量;
然后采用下式计算单胞刚度矩阵[C]:
μyx=μxyEy/Ex
μzx=μxzEz/Ex
μzy=μyzEz/Ey
Δ=(1-μxyμyx-μxzμzx-μyzμzy-2μxyμyzμzx)/(ExEyEz)
C11=(1-μyzμzy)/(EyEzΔ)
C22=(1-μxzμzx)/(ExEzΔ)
C33=(1-μxyμyx)/(ExEyΔ)
C12=C21=(μyx+μzxμyz)/(EyEzΔ)
C13=C31=(μxz+μxyμyz)/(ExEyΔ)
C23=C32=(μzy+μxyμzx)/(ExEzΔ)
C44=Gxy
C55=Gxz
C66=Gyz
式中,μxy为x对y方向的次泊松比,μyx为y对x方向的次泊松比,μzx为z对x方向的次泊松比,μxz为x对z方向的次泊松比,μzy为z对y方向的次泊松比,μyz为y对z方向的次泊松比,Gxy为xy方向的剪切模量,Gxz为xz方向的剪切模量,Gyz为yz方向的剪切模量,Cij,i=1...6,j=1...6是刚度矩阵[C]中的元素。
3.如权利要求2所述的一种陶瓷基复合材料超单元结构数值仿真计算方法,其特征在于:所述步骤四中,计算超单元刚度矩阵的具体步骤如下:
单胞刚度矩阵在超单元坐标系Oxyz和单胞坐标系Ox'y'z'分别为[C]xyz和[C]x'y'z',单胞刚度矩阵在超单元坐标系下通过如下公式进行计算:
[C]xyz=[T][C]x'y'z'[T]T
式中,[C]x'y'z'等于[C],矩阵[T]是根据两坐标系各轴线之间的夹角余弦值计算所得:
计算出超单元内部所有单胞在超单元坐标系下的刚度矩阵,然后按照各个单胞占超单元体积比的方法计算出超单元的刚度矩阵[K]:
式中,n为超单元内部的单胞数量,V单胞i为第i个单胞的体积。
4.如权利要求1所述的一种陶瓷基复合材料超单元结构数值仿真计算方法,其特征在于:所述步骤五中,在进行有限元仿真计算时使用渐进损伤方法刚度折减进行迭代计算,并输出计算仿真结果进行应力应变位移云图显示。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010875363.5A CN112100762B (zh) | 2020-08-27 | 2020-08-27 | 一种陶瓷基复合材料超单元结构数值仿真计算方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010875363.5A CN112100762B (zh) | 2020-08-27 | 2020-08-27 | 一种陶瓷基复合材料超单元结构数值仿真计算方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112100762A CN112100762A (zh) | 2020-12-18 |
CN112100762B true CN112100762B (zh) | 2024-03-08 |
Family
ID=73757958
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202010875363.5A Active CN112100762B (zh) | 2020-08-27 | 2020-08-27 | 一种陶瓷基复合材料超单元结构数值仿真计算方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112100762B (zh) |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2016004543A (ja) * | 2014-06-19 | 2016-01-12 | マツダ株式会社 | 有限要素解析装置、該方法及び該プログラム |
JP2018055509A (ja) * | 2016-09-29 | 2018-04-05 | ファイフィット株式会社 | 複合有限要素のプリ処理方法、複合材料の解析方法、解析サービスシステムおよびコンピュータ読み取り可能な記録媒体 |
CN109920495A (zh) * | 2019-03-28 | 2019-06-21 | 南京航空航天大学 | 一种编织陶瓷基复合材料强度的多尺度预测方法 |
CN110348165A (zh) * | 2019-07-18 | 2019-10-18 | 南京航空航天大学 | 基于结构网格的陶瓷基复合材料细观建模与力学计算方法 |
CN110795873A (zh) * | 2019-09-30 | 2020-02-14 | 北京擎靖天启科技服务有限公司 | 一种考虑尺寸控制的跨尺度拓扑优化方法 |
-
2020
- 2020-08-27 CN CN202010875363.5A patent/CN112100762B/zh active Active
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2016004543A (ja) * | 2014-06-19 | 2016-01-12 | マツダ株式会社 | 有限要素解析装置、該方法及び該プログラム |
JP2018055509A (ja) * | 2016-09-29 | 2018-04-05 | ファイフィット株式会社 | 複合有限要素のプリ処理方法、複合材料の解析方法、解析サービスシステムおよびコンピュータ読み取り可能な記録媒体 |
CN109920495A (zh) * | 2019-03-28 | 2019-06-21 | 南京航空航天大学 | 一种编织陶瓷基复合材料强度的多尺度预测方法 |
CN110348165A (zh) * | 2019-07-18 | 2019-10-18 | 南京航空航天大学 | 基于结构网格的陶瓷基复合材料细观建模与力学计算方法 |
CN110795873A (zh) * | 2019-09-30 | 2020-02-14 | 北京擎靖天启科技服务有限公司 | 一种考虑尺寸控制的跨尺度拓扑优化方法 |
Non-Patent Citations (3)
Title |
---|
Application of the Finite Element Method in the Analysis of Composite Materials: A Review;Sarah David Müzel 等;Polymers 2020;第12卷(第4期);第1-59页 * |
复合材料高精度宏-细观统一本构模型及其应用研究;孙志刚;中国优秀博硕士学位论文全文数据库 (博士) 工程科技Ⅱ辑;第C031-8页 * |
陶瓷基复合材料损伤耦合的宏细观统一本构模型研究;高希光;中国博士学位论文全文数据库 工程科技Ⅰ辑;第B020-20页 * |
Also Published As
Publication number | Publication date |
---|---|
CN112100762A (zh) | 2020-12-18 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN109920495B (zh) | 一种编织陶瓷基复合材料强度的多尺度预测方法 | |
Zhu et al. | Evaluation of failure criteria for fiber composites using finite element micromechanics | |
Römelt et al. | A multi-scale finite element approach for modelling damage progression in woven composite structures | |
Santos et al. | A hybrid-mixed finite element formulation for the geometrically exact analysis of three-dimensional framed structures | |
CN112632819B (zh) | 一种连续纤维增强复合材料基本力学性能参数预测方法 | |
CN107273566B (zh) | 一种构建复杂形体引力梯度场的计算方法 | |
CN103366085A (zh) | 编织复合材料力学性能的多尺度预测方法 | |
CN105808893A (zh) | 2.5维机织复合材料刚度预测方法 | |
Patel et al. | Damage and failure modelling of hybrid three-dimensional textile composites: a mesh objective multi-scale approach | |
CN115879346A (zh) | 基于改进型四节点逆有限元理论的结构应变场反演方法 | |
CN112100762B (zh) | 一种陶瓷基复合材料超单元结构数值仿真计算方法 | |
CN114117839A (zh) | 一种陶瓷基复合材料耦合损伤力学性能预测方法 | |
Cater et al. | Experimental and numerical analysis of triaxially braided composites utilizing a modified subcell modeling approach | |
Gowayed et al. | Modification and application of a unit cell continuum model to predict the elastic properties of textile composites | |
Fang et al. | Improved unit cells to predict anisotropic thermal conductivity of three-dimensional four-directional braided composites by Monte-Carlo method | |
Li et al. | Topology optimization of the microstructure of solid oxide fuel cell cathodes | |
Meyer et al. | Tensile specimen design proposal for truss-based lattice structures | |
CN109948253B (zh) | 薄板无网格Galerkin结构模态分析的GPU加速方法 | |
Benedikt et al. | On elastic interactions between spherical inclusions by the equivalent inclusion method | |
CN104392032A (zh) | 基于有限元法的纱线材料参数识别方法 | |
CN114741744B (zh) | 一种适用于针刺复合材料的细观建模及多尺度分析方法 | |
Sorini et al. | Development of a subcell based modeling approach for modeling the architecturally dependent impact response of triaxially braided polymer matrix composites | |
Oakeshott et al. | Development of a representative unit cell model for bi-axial NCF composites | |
CN118036406A (zh) | 一种树脂基复合材料的钻削载荷仿真方法 | |
Haasemann et al. | A new modelling approach based on Binary Model and X-FEM to investigate the mechanical behaviour of textile reinforced composites |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |