模态法辨识针刺C C 复合材料弹性模量
收稿日期:2013-12-24 作者简介:严博燕,1988年出生,助理工程师,主要从事发动机喷管设计。 E-mail: boyyan0115@163. com 模态法辨识针刺C/ C复合材料弹性模量 严博燕 侯 晓 陈 慧 刘 芹 (西安航天动力技术研究所,西安 710025) 文 摘 针对传统方法测量C/ C复合材料弹性模量精度差问题,提出了一种利用模态试验辨识材料弹性 模量的新方法。数值仿真结果表明板件材料振动频率与弹性模量之间存在模糊对应关系,利用此关系采用二 分法及二阶响应面近似模型对材料弹性模量进行辨识,并分别利用各向同性薄壁板件及正交各向异性薄壁板 件对此方法进行验证,计算结果误差最大为16.67%,最小仅为0.5%,表明该方法准确有效。最后将其应用于 针刺C/ C复合材料,获得了其弹性模量。仿真及试验结果均表明该方法无损、高效且准确度高。 关键词 针刺C/ C复合材料,弹性模量,模态试验,有限元数值模拟,二阶响应面模型 中图分类号:V45 DOI:10.3969/ j. issn.1007-2330.2014.03.005 Inverse Method Based on Modal Test for Characterizing C/ C Material Elastic Constants YAN Boyan HOU Xiao CHEN Hui LIU Qin (Xi’an Institute of Aerospace Solid Propulsion Technology, Xi’an 710025) Abstract In order to solve the problem of poor precision for measuring C/ C composite elastic constant by tradi- tional test method, the present paper provided a new method for obtaining the main elastic constants from modal test. According to the modal test data and the relationship between vibration frequencies and elastic constant, the scope of the elastic constants were estimated preliminarily. Then the quadratic response surface approximate model about elas- tic constants and vibration frequencies was established to optimize the elastic constants. The isotropic and orthotropic thin plates were used to verify this method. The results show that the maximum error is 16.67%, the minimum is only 0. 5%, so the method is accurate and effective. At last, the elastic constant of needled C/ C composite is obtained through the present method. The simulation and experimental results proved that the method is nondestructive,effi- cient and accurate. Key words Needled C/ C composite,Elastic constants, Modal Test, Finite element numerical simulation, Quadratic response surface approximate model 0 引言 工程上获得材料弹性模量的传统方法多为为静 态拉伸与扭转试验。传统方法测试复合材料性能参 数时,由于复合材料多为各向异性或正交各向异性, 材料均匀性较差,故需多个试验样本且试验结果往往 为局部特性而与整体性能参数差距较大。有必要通 过寻求其他新方法获得更为理想的结果。 材料的弹性模量在结构的振动行为上起着重要 作用。在平板上作用一诱导力使其产生振动,其振动 特性将只与平板的几何尺寸、密度、边界条件和弹性 模量有关。若已知平板的几何尺寸、密度和边界条 件,其振动特性将由弹性模量唯一确定。这使得利用 动力学法反推材料弹性模量成为可能[1]。 F. Forster 最先利用模态法反推得到了横梁结构的弹性模 量[2],随后Marco Alfano利用Warburton方程求得了 各向同性薄板材料弹性模量[3],Sol等人提出了正交 各向异性材料弹性模量的获得方法[4]。以上方法均 基于振动基本方程,利用横梁或薄壁板件振动方程中 频率与弹性模量的关系进行推导计算,进而获得材料 弹性模量。然而以上方法共同缺点是所选模型假设 条件较多且关系复杂,需作大量近似,误差相对较大。 随着计算机技术的发展,有限元数值模拟与现代优化 算法在该领域得到了广泛应用,Cugnoni、林玉华等 人[5-6]先后利用基因遗传算法对多种材料弹性模量 确定进行了优化计算,获得了较好的结果。该方法对 材料要求较低,适用性好,不仅局限于薄壁板件,对厚 —12—宇航材料工艺 http:/ / www. yhclgy. com 2014年 第3期 板材料亦适用。但缺点是可控性较差,计算量偏大。 本文结合以上两种思想的优点,利用复合材料平 板结构,通过模态试验获得全自由边界条件下试件振 型及其对应频率,利用低阶模态对应的频率与部分弹 性模量之间的关系对弹性模量进行初步优化以确定 弹性模量范围,在此基础上建立其与对应频率的二阶 响应面近似模型,对弹性模量进行精细优化,并将其 应用于针刺C/ C复合材料辨识其弹性模量。 1 理论与建模方法 1.1 振动模型的建立 若板件由正交各向异性材料构成,且弹性主方向 沿着板件的坐标轴x,y,z,则其振动基本方程如式 (1)所示。 Dx ∂ 4W ∂x4 + 2 Dxνyx + 2Dh( ) ∂4W ∂x2∂y2 + Dy ∂4W ∂y4 + ρh ∂2W ∂t2 = q (1) 式中: Dx = Exh 3 12 1 - νxyνyx( ) ; Dy = Eyh3 12 1 - νxyνyx( ) ; Dxνyx = Dyνxy ; Dh = Gxyh3 12 。 利用梁函数组合法可求取边界条件的正交异性 矩形板的固有频率。现取一项梁函数组合的振型表 达式,如式(2)所示。 fmn x,y( ) = AmnXm x( ) Yn y( ) (2) 代入正交异性板的变分方程,可求得适用于各种边界 条件正交异性矩形板各阶频率的统一解析表达式,如 式(3)所示。 fmn = 1 b2 ρh AmDx baæè öø 4 + 2 baæè öø 2 KmnDxνyx + 2LmnDh( ) + BnDy (3) 整理可得: 12b4ρf2mn h2 = Am baæè öø 4 + 2 baæè öø 2 Kmnνyx 1 - νxyνyx Ex + Bn 1 - νxyνyxEy + 4 b a æ è ö ø 2 LmnGxy (4) 由式(4)可知正交各向异性板件与弹性模量之 间的关系。然而这一公式较为近似,仅依靠此式不足 以精确求出弹性模量。 对矩形板进行有限元模拟,逐个改变弹性模量 值,观察频率值的变化,发现Gxy对f11有显著影响; Ex对f20有显著影响;Ey对f02有显著影响。即L11≫ B1、A1+K11;A2+K20≫ B0、L20;B2≫ A0+K02、L 02 故可以利用式(5)及f11、f20、f02,采用二分法缩小 弹性模量的范围。 fi - fi0 fi0 × 100% ≤ 5% (i = 1,2,3) (5) 式中,fi为数值模拟获得频率值,fi0为试验获得频率 值,i取1、2、3分别代表模态(1,1)、(2,0)、(0,2)。 在此基础上建立弹性模量与对应频率的二阶响 应面模型,对弹性模量在小范围内进行精细优化。 1.2 优化模型的建立 在自变量的某一区域内,采用低阶的多项式进行 响应面拟合。二阶多项式构造响应面模型[7]如式 (6)所示 Y~ = β - 0 + ∑ m i = 1 β - ixi + ∑ m i = 1 β - iix 2 i + ∑ m 1≤ i j≤ m β - ijxixj (6) 式中, β - 0、β - i、β - ii、β - ij为待定系数,m为设计变量个数。 其多元线性回归处理的形式为: Y~ = β0 + ∑ p i = 1 βivi (7) 式中,p为式(7)中除常数项外的所有项数,p与设计 变量个数m的关系为: p = 2m + m(m - 1) /2 (8) 利用n (n ≥ p + 1)组方案的试验数据 (Y~;ν1,ν2,… ,νp)k(k = 1,2,… ,n) ,求出待定系数 βi i = 1,2,… ,p( ) 。当设计矩阵V的秩≥ p+1时, VTV为非奇异阵,β的最小二乘估计为: β^ = VTV( ) -1VTY (9) 式中, V = 1 ν11 ν12 … ν1p 1 ν21 ν22 … ν2p … … … … … 1 νn1 νn2 … νnp é ë ù û Y = y1 y2 … yn é ë ù û β^ = β1 β2 … βp é ë ù û 若每组设计变量对应多个状态变量,模型扩展 为: β^ = β11 β12 … β1q β21 β22 … β2q … … … … βp1 βp2 … βpq é ë ù û 式中,q为每组设计变量对应的状态变量数。 文中采用正交试验设计理论对初始撒点的设计 变量进行设计。利用式(9)编程求出β的最小二乘 估计解,结合式(10)求出弹性模量的最优解。 —22—宇航材料工艺 http:/ / www. yhclgy. com 2014年 第3期 e = ∑ 3 i = 1 f - i - fi0 fi0 æ è ö ø 2 (10) 式中, f - i为由二阶响应面模型计算得到的频率值。 e 为目标值,取最小。 以上模型已编制成优化子程序,用来寻求所需常 数最优解。具体算法如下: (1)给出任意一组初始弹性模量值; (2)初步采用二分法缩小弹性模量的范围, fi - fi0 fi0 × 100% ≤ 5% (i = 1,2,3) 若弹性模量满足,转入步骤(3),否则继续采用二分 法缩小范围; (3)建立二阶响应面模型并利用多组数值模拟 值反算出系数矩阵β^ 。 Y~ = β - 0 + ∑ m i = 1 β - ixi + ∑ m i = 1 β - iix 2 i + ∑ m 1≤ i j≤ m β - ijxixj (4)计算弹性模量:利用响应面模型找出使e最 小的一组弹性模量值,输出计算结果。 e = ∑ 3 i = 1 f - i - fi0 fi0 æ è ö ø 2 2 验证算例及讨论 2.1 各向同性薄板材料 采用已知材料性能参数的玻璃/环氧树脂薄板材 料进行计算,试验所得模态及对应频率值(除去刚体 位移)如表1所示。 表1 玻璃/环氧树脂薄板试验模态值 Tab.1 Results of modal test of glass/epoxy 阶数频率值/ Hz对应模态阶数频率值/ Hz对应模态 f1 151.00 (1,1) f4 352.25高阶模态 f2 174.25 (2,0) f5 431.00高阶模态 f3 317.75 (0,2) f6 497.25高阶模态 对于各向同性薄板材料,有两个未知弹性模量 [Ey =Ex,νxy =Ex / (2Gxy)-1]。利用上述方法进行二 分逼近。在二分法计算结果的基础上,对Ex和Gxy 在±5 GPa范围内建立两变量二阶响应面模型。由于 需要至少6组初始撒点值。为体现设计变量的均匀 性,对其取三水平,遍取9组值。然后利用该响应面 模型在±5 GPa范围内找出最优解,二分法及最终优 化结果如表2所示。 表2 各向同性薄板计算结果 Tab.2 Results of isotropic thin plate 数据来源Ex =Ey / GPa Gxy / GPa νxy =Ex / (2Gxy)-1 文献参考值68.99 28.7 0.217 二分法计算值70 32.5 0.076 最终优化值70 28.9 0.211 误差/ % 1.464 0.697 2.73 由表2可以看出,弹性模量的最大误差仅为2. 73%,说明此方法适用于各向同性材料弹性模量的计 算。 2.2 正交各向异性薄板材料 利用该方法对正交各向异性薄板材料的弹性模量 进行计算。采用碳纤维薄板,尺寸为240 mm×180 mm ×2 mm,Ex =120 GPa,Ey =10 GPa,Gxy =4. 9 GPa,νxy = 0.3,ρ=1 510 kg/ m3,试验模态值如表3所示。 表3 正交各向异性薄板模态试验值 Tab.3 Results of modal test of anisotropic thin plate 阶数频率值/ Hz对应模态 f1 88.582 (1,1) f2 163.23 (0,2) f4 318.89 (2,0) 在二分法基础上建立二阶响应面模型,对于正交 各向异性材料,有四个未知弹性模量,其响应面模型有 15组待定系数,故需要至少15组数据,此处采用正交 试验设计理论对初始撒点设计变量进行设计[8],采用 四因素四水平正交试验设计方法。其中Ex、Ey、和Gxy 的取值范围为±6 GPa(Gxy最小取1),νxy的取值范围 为0.1 ~0.4。二分法及最终优化结果如表4所示。 表4 正交各向异性薄板计算结果 Tab.4 Results of anisotropic thin plate 数据来源Ex / GPa Ey / GPa Gxy / GPa νxy 文献参考值120 10 4.9 0.3 二分法计算值125 10.625 5 - 最终计算值120.6 10.3 5.2 0.25 最终误差/ % 0.5 3 6.12 16.67 由表4可以看出,利用此方法所得到的弹性模量 的计算值与真实值相对误差最小仅为0. 5%,最大为 16.67%。泊松比相对误差稍微偏大,这是因为其变 化范围很小,在该范围内其对频率的影响相对而言很 不明显,故相对误差较大。 结合以上各向同性薄板与正交各向异性薄板的 计算结果,可以看出此方法准确有效。 3 模态法确定针刺C/C复合材料弹性模量 3.1 模态试验设计 试样设计:试件材料为针刺C/ C复合材料,在真 实扩张段产品大端取样机加,除去边界约束后的有效 尺寸为140 mm×105 mm×3 mm。试件及材料由西安 航天复合材料研究所提供,密度为1 650 kg/ m3。 测试系统:采用LMS数据采集分析系统对针刺 C/ C复合材料薄板进行模态试验。加载方式采用橡 皮绳悬挂方式模拟全自由边界条件。激励方式采用 冲击信号激励,在试件表面均布9个激励点,利用力 —32—宇航材料工艺 http:/ / www. yhclgy. com 2014年 第3期 锤逐次激励各点[9]。采用LMS Test. Lab模态分析系 统,完成对数据的采集与处理分析。 3.2 结果及分析 表5为针刺C/ C材料薄板模态试验值,图1为 试验过程薄板模态振型图。 表5 针刺C/C薄板模态试验结果 Tab.5 Results of modal test of C/C plate 对应模态试件1# / Hz试件2# / Hz (1,1) 520.338 518.449 (2,0) 668.071 662.582 (0,2) 1243.768 1240.888 (1,1) (2,0) (0,2) 图1 针刺C/ C薄板模态振型图 Fig.1 Modes of needled C/ C composites 由以上结果可以看出,两块试样的模态试验结果 相差极小,表明试样取样材料均匀,试验结果一致性 良好,可以进行常数确定。 对试验数据进行平均处理并利用本文方法进行 计算,获得结果如表6所示。 表6 针刺C/C薄板材料常数计算结果 Tab.6 Final result of C/C plate Ex / GPa Ey / GPa Gxy / GPa νxy 29.78 29.05 10 0.36 将此计算值再次带入有限元程序进行模拟计算, 并与试验值进行对比,结果如表7所示。 表7 数值模拟与试验对比结果 Tab.7 Comparison of test results and numerical simulation results 对应模态数模结果/ Hz试验值/ Hz相对误差/ % (1,1) 521.09 519.39 0.33 (2,0) 665.96 665.32 0.10 (0,2) 1244 1242.33 0.13 由表7可以看出,利用本文方法得到的弹性模量 值计算出的频率值与试验得到的频率值非常接近,最 大误差仅为0. 33%,表明本文方法获得的弹性模量 值可认为为近似真实材料性能常数。此外,从文中计 算结果可以看出Ex≈ Ey,Gxy≈ Ex /2(1+νxy),即该材 料可认为是宏观各向同性材料,这与其制造工艺所得 结论是相吻合的。 4 结论 (1)采用模态法辨识材料弹性模量方法在各向 同性薄壁板件和正交各向异性薄壁板件验证算例中, 计算结果误差最大为16.67%,最小仅为0. 5%,表明 该方法准确有效; (2)该方法获得的针刺C/ C材料弹性模量值计算 出的频率值与试验得到的频率值非常接近,最大误差仅 为0.33%,表明该方法获得的材料常数是近似真实值。 需要指出的是,本文采用的二阶响应面模型对设计 变量较少的模型能较好的满足其精度要求,若设计变量 较多,其精度将大大降低。文中只对四变量以下的板件 进行了计算,而对于设计变量更多的材料,例如正交各 向异性厚壁板件,此方法将较难得到理想结果。 参考文献 [1] Deobald L R, Gibson R F. Determination of elastic con- stants of orthotropic plates by a modal analysis/ rayleigh-ritz tech- nique [J]. Journal of Sound and Vibration,1988,124(2):269- 283 [2] Forster F. Ein neues messverfahren zur bestimmung des elastizitats-moduls und der dampfung[J]. Zeitschrift fuer Met- allkunde,1937,29:109-115 [3] Alfano Marco, Pagnotta Leonardo. Determining the elas- tic constants of isotropic materials by modal vibration testing of rectangular thin plate[J]. Journal of Sound and Vibration, 2006, 293(2):426-439 [4] Sol H. Identification of anisotropic plate rigidities using free vibration data[D]. Vrije Universiteit Brussel, Belgium,1986 [5] Cugnoni Joel,Gmur Thomal,Schorderet Alain. Inverse method based on modal analysisfor characterizing the constitutive properties of thick composite plates [J]. Computer and Struc- tures, 2007,56(3):1310-1320 [6]林玉华.以振动测试及最佳化方法反算材料之弹性 模量[D].国立云林科技大学,2010 [7]曾声宇.基于相应面近似的反舰导弹并行子空间优 化设计技术[D].西安:西北工业大学,2004 [8]王世磊.实验设计与数据处理[D].郑州:郑州大学, 2009 [9]张力.模态分析与实验[M].北京:清华大学出版社, 2011 (编辑 李洪泉) —42—宇航材料工艺 http:/ / www. yhclgy. com 2014年 第3期