哈尔滨工业大学学报  2021, Vol. 53 Issue (12): 153-163  DOI: 10.11918/202004043
0

引用本文 

张源, 李范春, 贾德君. 常规立方晶格点阵轮的设计与静力学分析[J]. 哈尔滨工业大学学报, 2021, 53(12): 153-163. DOI: 10.11918/202004043.
ZHANG Yuan, LI Fanchun, JIA Dejun. Design and static analysis of conventional cubic lattice impeller[J]. Journal of Harbin Institute of Technology, 2021, 53(12): 153-163. DOI: 10.11918/202004043.

基金项目

教育部“双一流”学科建设基金(SSCXXM030)

作者简介

张源(1996—),男,硕士研究生;
李范春(1960—),男,教授,博士生导师

通信作者

李范春,lee_fc@126.com

文章历史

收稿日期: 2020-04-08
常规立方晶格点阵轮的设计与静力学分析
张源, 李范春, 贾德君     
大连海事大学 船舶与海洋工程学院,辽宁 大连 116026
摘要: 为实现压气机叶轮在刚度和强度具有可调性的同时降低叶轮质量,基于常规立方晶格结构设计出一种新型轻量化压气机叶轮。通过3种模型悬臂梁比较计算,证明渐进均匀化(asymptotic homogenization,AH)方法在大规模点阵计算问题上具有较好的计算精度。在此基础上,采用AH方法对不同单胞填充率的点阵轮在仅考虑80 000 r/min的惯性荷载时的静力学性能进行了数值计算。研究结果表明,对于所研究路径,点阵轮的变形和应力介于无填充轮和实心轮之间,大叶片外缘的周向变形最大与最小变形差值(极差)在所有研究工况下均小于无填充轮和实心轮。对于填充率为0.4的情况,其周向变形极差值比实心轮约低23.25%,比无填充轮约低55.46%,其质量相对于实心轮可降低17.08%。这意味着点阵轮相对于传统压气机叶轮除了能较大降低叶轮质量,还将具有更好的周向抗畸变能力和更高的工作效率。同时,点阵轮相对于无填充轮和实心轮叶片边缘具有更小的轴向应力,因此可为设计更高转速的压气机叶轮提供强度保证,同时也为叶轮结构的轻量化设计提供了一种新的设计思路。
关键词: 机械设计    叶轮    点阵    轻量化设计    AH方法    3D打印    
Design and static analysis of conventional cubic lattice impeller
ZHANG Yuan, LI Fanchun, JIA Dejun     
Naval Architecture and Ocean Engineering College, Dalian Maritime University, Dalian 116026, Liaoning, China
Abstract: In order to realize the adjustable stiffness and strength of compressor impeller and reduce the mass of the impeller, a new lightweight compressor impeller was designed based on the conventional cubic lattice structure. Based on the comparative calculation of three models of cantilever beam, it has been proved that the asymptotic homogenization (AH) method has good accuracy in large-scale lattice calculation. On this basis, the AH method was used to calculate the static performance of lattice impeller with different cell filling rates under the inertial load of 80 000 r/min. Results show that for the studied path, the deformation and stress of the lattice impeller were between those of the unfilled impeller and the solid impeller, and the difference (range) between the maximum and minimum circumferential deformation of the outer edge of the large blade was smaller than those of the unfilled impeller and the solid impeller. When the filling rate was 0.4, the circumferential deformation range of the lattice impeller was 23.25% lower than that of the solid impeller, 55.46% lower than that of the unfilled impeller, and the mass of the lattice impeller could be reduced by 17.08% compared with solid impeller. It indicates that in comparison with the traditional compressor impeller, lattice impeller can not only greatly reduce the impeller mass, but also has better circumferential distortion resistance and higher work efficiency. Meanwhile, the blade edge of lattice impeller has smaller axial stress than those of unfilled impeller and solid impeller, so it can provide strength guarantee for the design of compressor impellers with higher speed. The design of this paper also provides a new design idea for the lightweight design of impeller structures.
Keywords: machine design    impeller    lattice    lightweight design    asymptotic homogenization (AH) method    3D printing    

近年来,增材制造(additive manufacturing, AM)技术取得了前所未有的发展,结构轻量化设计也伴随着AM技术的发展使设计师拥有更灵活的设计空间,尤其是在航空航天[1-2],汽车[3-5],船舶[6-7],仿生学[8-9]和生物医学[8, 10-12]应用中。AM技术的快速发展使CAD工程师在自由形式的结构设计中创造了更多的可能性,点阵结构也伴随着AM技术的发展而飞速发展为目前最热门的研究领域之一。应用点阵结构可以实现飞行器结构的轻量化设计,压气机叶轮是飞行器上最关键且最难以轻量化设计的部件之一,刚度条件和气动性能的限制使得压气机叶轮的轻量化设计更加困难,点阵结构的出现为压气机叶轮的轻量化设计提供了一种新的可能。

AM技术的发展使得包括点阵结构在内的一系列复杂、轻质、高性能结构得以制作生产。点阵结构由于其高强重比,优异的散热性能,出色的减震和降噪功能成为设计工程师在结构轻量化设计时的重要考虑方案之一,点阵结构被大量的应用在多功能、工况复杂的设计中。Abdi等[5]基于拓扑优化技术和点阵技术,采用Iso-XFEM方法,该方法能够利用结构性能准则的等值线/等值面和扩展有限元法(extended finite element method, XFEM)生成高分辨率拓扑优化解,设计出一种轻质且能够增加制动稳定性汽车制动踏板。Robbins等[13]将具有变密度的点阵六面体晶胞应用到某一悬臂梁拓扑优化结构中。基于由均质化方法确定的共形点阵填充结果来创建具有外部拓扑的结构。为了保证具有相同单位密度的结构的相同质量,相应地调整了外部拓扑,从而使悬臂梁点阵密度较低的同时具有更小的挠度。Daynes等[14]提出了一种新的方法来产生优化的功能梯度晶格结构。首先通过拓扑优化得到最佳堆芯密度分布,从而在质量约束下最小化结构的柔度,然后根据局部主应力构造一系列等静线,生成晶格结构,该晶格结构在空间上随晶格单元大小、纵横比和取向而分级。为了验证这一新方法的有效性,通过对3点弯曲的夹层梁结构进行了优化,实验测试证明,与具有相同密度的均匀单元尺寸的基准芯相比,通过对晶格单元进行空间分级,芯的刚度和强度性能得到了极大改善。Lynch等[15]基于拓扑优化和点阵结构设计出一种独特的轻量化组件和增加多功能性的方法,用于根据工作载荷调整局部单元密度。为了验证所提方法的合理性,他们将其应用于一种壳体类套管试件的优化、分析、制造和力学试验验证。使用基于应力的均匀化拓扑优化方法对试件进行优化,与具有相同形状因子的实心、全致密套管相比,质量减轻了53%。通过高密度元素分析研究了优化的几何结构,然后进行了增材制造。通过力学试验证明用于优化的均质有限元模型、高密度有限元模型和实验结果之间的良好相关性。他们的研究结果验证了针对特定用途和负载情况的优化方法,并开始将该方法视为一种可接受的方法。除了基于点阵结构实现某种结构的设计,也有团队将点阵结构应用于3D打印的支撑设计中,如Cheng等[16]为探讨用拓扑优化方法设计支撑结构以减轻残余应力引起的建筑物破坏的可行性,采用固有应变法对AM结构中的残余应力进行了快速预测。由于周期性网架结构具有开孔、自支撑的特点,采用分级网架结构优化设计支撑结构,为了验证该方法的可行性,分别设计了双悬臂梁和髋关节种植体的支撑结构。他们的研究证明优化后的支撑结构可实现约60%的质量折减,并且优化后的支撑结构构件在AM实现设计后不再出现应力开裂现象,证明了该方法的有效性。

目前也有部分团队将3D打印技术应用于叶轮结构的3D打印制作,由于结构的复杂性,目前的应用主要是用于树脂材料完成叶轮结构加工或者用于叶轮结构铸造型芯的制备。例如,Meli等[17]提出了一种基于拓扑优化和增材制造的离心式压缩机叶轮工程化在设计域制造的总体框架,同时对一家意大利大型石油天然气公司的叶轮进行了完整的重新设计,证明了所提方法的可行性。Kim等[18]分别采用粉末床粘合剂喷射打印(binder jetting printing, BJP)和选择性激光烧结(selective laser sintering, SLS)技术完成了叶轮叶片型芯的制作,并对两种方式制作的芯样进行了比较,结果发现采用SLS技术制备的芯样具有更高的强度,并且具有更高的制造精度。由于复合材料制作的3D打印叶片具有较好的制作精度,对于部分受载较小的离心式叶片,也有团队采用聚氨酯(necuron)和ABS(3D打印)热塑性复材实现叶轮结构的增材制造[19]。然而,由于热力学过程和加工过程的复杂性,金属叶轮的3D打印相关研究还较少,本文在之前的研究工作中[20],已经利用SLM打印机验证了Ti6Al4V合金叶轮的可加工性和加工精度。

点阵结构的众多优势意味着将点阵结构应用于压气机叶轮的轻量化设计中也将使叶轮具备许多原始设计所不具备的优势。在文献[20]中首次提出了点阵压气机叶轮的设计并利用金属3D打印仿真技术对具有点阵填充的压气机叶轮的加工性能进行了研究分析,然而缺乏对其力学性能的研究,本文将进一步讨论点阵压气机叶轮在使用时的受力与变形情况。

1 方法选取

点阵结构具有复杂的几何形貌,对于具有点阵填充的复杂零件,直接进行有限元计算意味着极大的计算资源的消耗,甚至是不能求解,因此常用简化方法实现点阵结构的数值计算。目前针对点阵填充零件的常用计算方法主要有两种,一种是将点阵部分的实体单元简化为梁单元进行计算[21],梁单元节点处刚接,另一种则是基于多尺度算法的渐进均匀化方法(AH)[22-24]。为比较两种简化方法在不同填充率下相对于实体单元计算结果的误差大小,本文以某一点阵悬臂梁为例,对两种简化方法在某一路径上的变形与应力进行有限元计算。通过与实体单元计算出的点阵悬臂梁计算结果进行比较,选择出误差相对较小的简化方法,进一步为下文求解点阵压气机叶轮的力学性能奠定基础。悬臂梁的尺寸如图 1所示,A为模型设计域。悬臂梁上下表面均留取0.4 mm厚度的薄板,以便于施加载荷和选取研究路径,悬臂梁尺寸为4.0 mm×4.8 mm×40.0 mm,设计域尺寸为4 mm×4 mm×40 mm。分别采用梁单元立方点阵胞元,实体单元立方点阵胞元和均匀化等效立方点阵胞元材料对设计域A进行填充。

图 1 横梁尺寸及设计域 Fig. 1 Cross-beam size and design area

选择悬臂梁上表面中部的一条路径作为研究路径,如图 2所示,研究路径的方向由1端指向2端,即1处为路径上的坐标原点。

图 2 悬臂梁研究路径选择 Fig. 2 Research path selection of cantilever beam

悬臂梁模型的边界条件如图 3所示,B面施加0.625 MPa的均布荷载f,设计域的受载通过连接设置由B面传递过来,A面为固定端,采用固定约束。3种填充模型具有相同的边界条件,梁模型约束A面所有梁单元节点。

图 3 边界条件 Fig. 3 Boundary condition

材料设置为金属3D打印中常用材料之一的TiAl6V4合金。设计域A内采用具有不同填充率的立方点阵结构进行填充,单个立方点阵胞元的填充率分别设置为0.2, 0.4, 0.6和0.8。3种模型在经过网格无关性验证过后,对比不同尺寸下的计算误差,发现相对误差值小于1%。在保证各种模型的相网格质量满足计算要求前提下,对具有不同填充率的3种点阵梁模型进行有限元分析。

提取不同填充率下的点阵梁在图 2所示路径中的变形与应力分布。单个胞元填充率分别为0.2、0.4、0.6和0.8时,点阵梁在路径上的变形与应力值分别如图 45所示。对比上述4种填充率下的AH模型和梁单元模型在选择路径上的应力和变形分布,发现无论是变形或者应力,AH模型在各个填充率下的变形与应力值均与实体单元模型更为接近。这意味着采用AH方法对点阵模型进行简化比采用梁单元对模型进行简化能得到更好的精度,并且随着单个点阵填充率的增大,这种优势更加明显。基于以上分析,渐进均匀化方法(AH)在本文中被用于计算点阵材料的等效性能,进而实现点阵结构的简化计算。

图 4 3种计算模型在不同填充率下的变形分布 Fig. 4 Deformation distribution of three calculation models under different filling rates
图 5 3种计算模型在不同填充率下的应力分布 Fig. 5 Stress distribution of three calculation models under different filling rates
2 均匀化方法的晶格材料等效弹性特性[25]

图 6示意性地描绘出了AH方法的基本原理。假定设计域Βγ被点阵结构填充,牵引力t服从于边界Γt,位移m作用于边界Γu,而体力f在整个区域内均匀分布。采用AH方法,将晶格结构域视为具有弹性和屈服强度等效性质的连续实体域。因此,对于具有显式表示的晶格结构的全尺寸模拟,可以显著降低昂贵的计算成本。

图 6 均匀化方法的基本原理 Fig. 6 Basic principle of AH method

AH方法假定具有周期性微观结构的材料取决于两个尺度:一个在宏观尺度上平稳变化,另一个在微观尺度上周期性变化。通过引入放大因子γ(0 < γ≪1),微观坐标y=[yi](i=1, 2, 3)与宏观坐标x通过$y=\frac{x}{\gamma}$相关,这意味着微观水平上的位移等场量变化比宏观水平快1/γ倍。因此,任何场量都可以表示为两个尺度的渐近展开式:

$ u^{\gamma}(x, y)=u_{0}(x, y)+\gamma u_{1}(x, y)+\gamma^{2} u_{2}(x, y)+\cdots $ (1)

式中: uγ为场量的精确值,u0为宏观或平均值,u1u2为由于微观水平下的微观结构而引起的扰动。对于线弹性,u0为宏观尺度上的平均位移,而u1u2,…, 等为微观尺度上的微扰位移。计算式(1)关于坐标x的导数并在小变形假设下代入应变位移方程,应变张量为

$ \begin{aligned} \boldsymbol{\varepsilon}_{i j}=& \frac{1}{2}\left(\frac{\partial \boldsymbol{u}_{i}^{\gamma}}{\partial \boldsymbol{x}_{j}}+\frac{\partial \boldsymbol{u}_{i}^{\gamma}}{\partial \boldsymbol{x}_{i}}\right)=\\ & \frac{1}{2}\left(\frac{\partial \boldsymbol{u}_{0 i}}{\partial \boldsymbol{x}_{j}}+\frac{\partial \boldsymbol{u}_{0 j}}{\partial \boldsymbol{x}_{i}}\right)+\frac{1}{2}\left(\frac{\partial \boldsymbol{u}_{1 i}}{\partial \boldsymbol{x}_{j}}+\frac{\partial \boldsymbol{u}_{1 j}}{\partial \boldsymbol{x}_{i}}\right)+O(\gamma) \end{aligned} $ (2)

式中O(γ)为带有γ的项。通过忽略O(γ)和高阶项,可以将点阵结构域Βγ中的应变表示为平均应变和波动应变的组合为:

$ \boldsymbol{\varepsilon}_{i j}=\overline{\boldsymbol{\varepsilon}}_{i j}+\boldsymbol{\varepsilon}_{i j}^{*} $ (3)
$ \overline{\boldsymbol{\varepsilon}}_{i j}=\frac{1}{2}\left(\frac{\partial \boldsymbol{u}_{0 i}}{\partial \boldsymbol{x}_{j}}+\frac{\partial \boldsymbol{u}_{0 j}}{\partial \boldsymbol{x}_{i}}\right), \boldsymbol{\varepsilon}_{i j}^{*}=\frac{1}{2}\left(\frac{\partial \boldsymbol{u}_{1 i}}{\partial \boldsymbol{x}_{j}}+\frac{\partial \boldsymbol{u}_{1 j}}{\partial \boldsymbol{x}_{i}}\right) $ (4)

式中:εij为平均应变应变张量或宏观应变张量,εij*为假设波动应变是周期性变化的。通过将式(3)代入弱形式平衡方程,并假设虚拟位移在宏观水平上是一个常数,并且仅在微观水平上发生变化,因此可以得到以下单位胞元平衡方程:

$ \begin{aligned} &\int_{V_{\mathrm{RVE}}} C_{i j p q} \boldsymbol{\varepsilon}_{i j}^{1}(v) \boldsymbol{\varepsilon}_{p q}^{*}(u) \mathrm{d} V_{\mathrm{RVE}}= \\ &\ \ -\int_{V_{\mathrm{RVE}}} C_{i j p q} \boldsymbol{\varepsilon}_{i j}^{1}(v) \bar{\boldsymbol{\varepsilon}}_{p q}(u) \mathrm{d} V_{\mathrm{RVE}}, \quad \forall v \in U_{a d} \end{aligned} $ (5)

式中: VRVE为可以在晶格结构中重复的最小代表性晶格单位的体积,其中RVE为重复的最小代表性晶格单元; v为虚拟位移,Uad为允许的位移场,εij1(v) 为虚拟微观应变,在微观水平上周期性变化。用有限元法离散化并求解方程中的单元平衡方程(5)在RVE上具有周期性边界条件。为了求解该方程,波动应变用平均应变表示为

$ \boldsymbol{\varepsilon}_{i j}^{*}=-\boldsymbol{\varepsilon}_{i j}^{* k l} \overline{\boldsymbol{\varepsilon}}_{k l} $ (6)

式中εij*kl为应变场的周期性,通过求解式(5)获得。单位平均应变在等式(5)的右侧。等式(3)中RVE的局部应变可以使用结构张量表示为

$ \boldsymbol{\varepsilon}_{i j}=\boldsymbol{M}_{i j k l} \overline{\boldsymbol{\varepsilon}}_{k l} $ (7)

式中Mijkl为局部结构张量,可表示为

$ \boldsymbol{M}_{i j kl}=\frac{1}{2}\left(\boldsymbol{\delta}_{i k} \boldsymbol{\delta}_{j l}+\boldsymbol{\delta}_{i l} \boldsymbol{\delta}_{j k}\right)-\boldsymbol{\varepsilon}_{i j}^{* kl} $ (8)

代入式(7)转化为弱形式的平衡方程并在RVE上积分,可获得等效的晶格结构弹性张量为

$ \overline{\boldsymbol{C}}_{i j k l}=\frac{1}{\left|V_{\mathrm{RVE}}\right|} \int_{V_{\mathrm{RVE}}} \boldsymbol{C}_{i j m n} \boldsymbol{M}_{m n k l} \mathrm{~d} V_{\mathrm{RVE}} $ (9)

式中: Cijmn为RVE材料成分的弹性张量,而Cijkl为通过AH方法计算的等效弹性张量。本文以立方晶格材料为目标晶格结构。由于立方晶格结构的对称性,有3个独立的弹性常数。因此,立方晶格结构的弹性响应可以用以下本构方程表示:

$ \left[\begin{array}{c} \overline{\boldsymbol{\sigma}}_{11} \\ \overline{\boldsymbol{\sigma}}_{22} \\ \overline{\boldsymbol{\sigma}}_{33} \\ \overline{\boldsymbol{\sigma}}_{12} \\ \overline{\boldsymbol{\sigma}}_{23} \\ \overline{\boldsymbol{\sigma}}_{31} \end{array}\right]=\left[\begin{array}{cccccc} \bar{C}_{11} & \bar{C}_{12} & \bar{C}_{12} & 0 & 0 & 0 \\ \bar{C}_{12} & \bar{C}_{11} & \bar{C}_{12} & 0 & 0 & 0 \\ \bar{C}_{12} & \bar{C}_{12} & \bar{C}_{11} & 0 & 0 & 0 \\ 0 & 0 & 0 & \bar{C}_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & \bar{C}_{44} & 0 \\ 0 & 0 & 0 & 0 & 0 & \bar{C}_{44} \end{array}\right]\left[\begin{array}{l} \overline{\boldsymbol{\varepsilon}}_{11} \\ \overline{\boldsymbol{\varepsilon}}_{22} \\ \overline{\boldsymbol{\varepsilon}}_{33} \\ 2 \overline{\boldsymbol{\varepsilon}}_{12} \\ 2 \overline{\boldsymbol{\varepsilon}}_{23} \\ 2 \overline{\boldsymbol{\varepsilon}}_{31} \end{array}\right] $ (10)

式中: σijεij分别为宏观应力张量和宏观应变张量。C11C12C44为3个独立的等效弹性常数,由等式(9)计算得出。可以看出,等效特性是由材料组成和微观结构的几何形状决定的。对于这项工作的研究,晶格结构的几何是预先定义的(即立方晶格结构由RVE中的12个1/4圆柱体构成)。等效性质由材料组成(单位晶胞的体积分数或相对密度)独立确定。因此,假定等效性质是相对密度的函数。这里采用四阶多项式来描述等效关系。在不失一般性的前提下,等效弹性特性通过基础材料的特性进行正交化,可以这样写为

$ \boldsymbol{C}_{i j}=a_{0}+a_{1} \rho+a_{2} \rho^{2}+\cdots $ (11)

式中: Cij=Cij/Cij*为正交化等效弹性常数,Cij为等效弹性常数,而Cij*为固体材料组分的弹性常数。ρ为单位晶格单元的相对密度,ai(i=0, 1, 2, …)为通过拟合仿真或实验数据获得的模型参数。

3 基于AH的点阵压气机叶轮静力学分析

本文的点阵压气机叶轮是在实心压气机叶轮基础上,采用常规立方点阵结构叶轮轮毂部分进行局部填充得到的,所得点阵轮模型如图 7所示。在本文此前的研究工作中[17],已经利用有限元方法对点阵轮3D打印加工中的材料堆叠过程以及打印参数对打印结果的影响进行了系统分析,研究表明,相同加工工况下,点阵轮无论是从加工所得残余变形还是残余应力幅值上均是具有一定优势的。在本文的研究工作中,将利用AH方法,进一步对点阵轮在理想工况下运行时的力学性能进行分析。叶轮的几何参数见表 1

图 7 点阵压气机叶轮模型 Fig. 7 Lattice compressor impeller model
表 1 叶轮的几何参数 Tab. 1 Geometric parameters of impeller
3.1 基于AH方法的点阵轮计算模型简化

AH方法需要点阵结构的最小代表性晶格单元满足宏观尺度上的平稳性和微观尺度的周期性。本文在微观尺度上定义4种填充率的最小周期点阵胞元,单个胞元的填充率分别取0.2、0.4、0.6和0.8,胞元尺寸设置为2 mm。点阵轮的加工需要采用3D打印技术,本文点阵轮的材料设置为TiAl6V4合金。在单个晶格填充率分别为0.2、0.4、0.6和0.8下,点阵轮相对于实心轮质量可分别降低23.69%、17.08%、12.16%和6.53%。常规立方晶格在3个方向上具有相同的几何形状,采用AH方法能计算出3个方向上具有相同的弹性常数。

3.2 点阵轮的静力学计算

叶轮实际工作时,受载主要来自于流动的气流施加给叶片和轮毂的面压和由于自身质量所产生的“惯性力”,当考虑到气流扰动等不稳定因素时叶轮的实际工况则更为复杂。在本文的研究中,将只考虑均匀高转速时,叶轮在自身的惯性荷载下叶轮叶片和轮毂的受力及变形情况。叶轮的边界条件定义如图 8(a)所示,对轮毂小轴位置施加柱支撑,绕z轴施加80 000 r/min的匀转速。定义x轴为径向,y轴为周向,x轴为轴向。研究路径选择如图 8(b)所示。其中,大叶片直接影响叶轮的气动性能,因此路径1选择为叶轮的大叶片边缘。轮毂小轴孔位置与轮轴直接相连,应力较大,因此路径2选择为小轴孔表面。本文的所有仿真计算结果采用ANSYS完成。

图 8 叶轮研究路径选择及边界条件 Fig. 8 Boundary conditions and path selection of impeller

采用有限单元法实现叶轮在给定工况下的应力及变形分析。叶轮的网格划分情况如图 9所示。为减小计算误差,叶片和轮毂的网格单元均采用20节点六面体单元,其中叶片部分采用扫掠法实现网格划分。在采用不同尺寸的网格单元对叶轮的网格无关性进行验证后,发现同网格单元间计算出的相对误差值小于1%。经过计算最终确定轮毂部位网格尺寸为0.5 mm,叶片部位网格尺寸为0.4 mm,扫掠法中划分数目设置为20,此时叶轮具有良好的计算精度。

图 9 点阵轮的网格划分 Fig. 9 Grid division of lattice impeller

本文所涉及到的变形云图单位为mm,应力云图单位为MPa。对于单个点阵胞元填充率为0.4的点阵轮,其轮毂和叶片的变形及应力分布如图 10所示。其中,图 10(a)为采用AH方法计算的该填充率下的点阵轮变形云图,对于轮毂和叶片,此时的最大变形幅值为0.048 9 mm,变形主要分布在大叶片边缘和轮毂外边缘处。图 10(b)为该填充率下的应力分布,叶片和轮毂在该填充率下的应力幅值为309.27 MPa(注意此时叶轮的应力并未包含叶轮内部点阵区域的应力)。

图 10 点阵轮轮毂与叶片的变形与应力分布 Fig. 10 Deformation and stress distribution of lattice impeller hub and blades

为研究采用均匀化方法所得仿真结果的应力分布,在大规模的点阵计算问题上,采用有限元分析的子模型法分析点阵轮内部的应力分布。子模型方法的基本原理是在多尺度计算的基础上,以宏观变形或应力作为边界切割插值,求解局部结构的应力或变形情况。对于叶轮内部的点阵区域,由于边界上周期性条件的破坏,会导致均匀化中的“边界层现象”[26],这导致了边界上的应力预测不准确,但是这不影响AH区域外的应力预测结果,如图 5结果所示,AH模型的应力预测结果与全尺寸模型的预测结果仍然有非常好的吻合性。而对于离边界较远的区域,AH模型相对于全尺寸模型的预测误差则小于5%[22]。为了降低AH区域内由于“边界层现象”引入的误差,本文对设计域去除1 mm厚度的边界层。由于采用AH方法是将点阵材料模型等效为均匀的等效模型,因此通过AH方法不能直接观测到内部点阵结构的变形及应力分布情况。为研究点阵轮内部点阵的受力情况,对去除边界层的AH模型采用子模型的方法,将通过AH方法算得的变形施加至点阵部位的各个面上作为计算内部点阵结构的力学性能的边界条件。通过导入的边界条件,可以求得内部点阵结构的变形及应力大致分布。

由于点阵轮内部的点阵结构具有一定的对称性,选取点阵域的代表单元进行受力分析即可阐明整个点阵域的受力情况。本文选取整个点阵域的1/4作为点阵域的代表单元,此时已去除边界层。

图 11(a)11(b)分别表示所选点阵域代表单元的变形和应力分布。对于单个胞元填充率为0.4的点阵叶轮,在不考虑边界应力和变形时,其内部点阵的变形幅值为0.024 2 mm,变形值沿径向逐渐增大,靠近轮缘边缘处的变形最大。内部点阵结构的应力幅值为240.75 MPa,对于代表单元,该种类型的点阵填充形式应力较大处主要分布在45°位置,对于内部点阵单元的其他部位,有相同的应力分布。对于3D打印Ti6Al4V粉末材料,在铺粉层厚为60 μm,激光功率为400 W并且在未进行过热处理或热等静压处理的条件下,利用SLM280型金属3D打印机进行3D打印加工的打印后材料性能见表 2。其屈服强度远远大于由于惯性荷载在点阵轮内部产生的最大应力240.75 MPa。

图 11 点阵域代表单元选取及代表单元的应力及变形分布 Fig. 11 Selection of representative elements and corresponding stress and deformation distributions
表 2 3D打印后Ti6Al4V材料的物理性能 Tab. 2 Physical parameters of TiAl6V4 after 3D printing

作为比较组,实心轮和无填充轮的变形及应力分布分别如图 1213所示。从变形分布上来看,无论是对于实心轮、无填充轮还是点阵轮,三者的变形分布趋势大致是相同的,均是变形沿着径向变形幅值逐渐增大,最大变形发生在大叶片和轮毂边缘处。而应力分布则是点阵轮与无填充轮更为接近,在应力幅值上(此时不考虑点阵域边界层),点阵轮小于无填充轮和实心轮。

图 12 实心轮变形及应力分布 Fig. 12 Deformation and stress distribution of solid impeller
图 13 无填充轮变形及应力分布 Fig. 13 Deformation and stress distribution of unfilled impeller

本文能直观地从云图中判断出叶轮各个位置的变形和应力的大致范围以及幅值分布。但是对于本文所示的点阵轮结构,光从分布上难以直观看出点阵轮和实心轮以及无填充轮的力学性能差异,然而叶轮叶片或轮毂的微小差异都会造成实际叶轮气动性能与设计气动性能的不同。为进一步判断不同点阵填充率下点阵轮的静力学性能与工程上常用的实心轮和无填充轮的差异,需要对路径上的数据进行提取比较。因此取变形较大的路径1上3个方向上的变形和应力,取应力较大的路径2上3个方向上的应力,对3种叶轮进行比较分析。不同叶轮在路径1上的变形分布如图 14所示。图 14(a)为路径1上的径向变形,图 14(b)14(c)分别为路径1上的周向变形和轴向变形。

图 14 路径1变形 Fig. 14 Deformation of path 1

图 14(a)所示,叶轮的径向变形受单个晶格填充率影响明显。点阵轮在路径1上的径向变形位于实心叶轮和无填充叶轮之间,随着点阵填充率的减小,叶轮路径1上的径向变形呈非线性增大。这种变化趋势理论上是合理的,因为叶轮刚度受叶轮填充率的影响,在“离心力”驱动下,叶轮轮毂沿径向向外膨胀,当填充率较大时,叶轮刚度大,变形则较小,反之当点阵填充率较小时,叶轮的径向变形则较大。图 14(a)的分析结果能有效证明有限元计算的合理性。点阵轮在路径1上的周向变形也位于实心叶轮和无填充叶轮之间,如图 14(b)所示,叶轮在路径1不同位置的周向变形波动较大,与径向变形相同的是路径1上的周向变形值也随着点阵轮填充率的减小,由实心叶轮逐渐向无填充叶轮靠近,当单个晶格填充率为0.2时,点阵轮在路径1的周向变形曲线与无填充轮基本一致。点阵轮的周向变形和径向变形介于实心轮与无填充轮之间意味着存在这种可能性,即可以通过改变点阵轮的点阵填充率,在降低叶轮质量的同时使叶轮的气动性能实现在实心轮和无填充轮之间的调控,从而满足某些特定工况下的需要。

大叶片的形变情况直接影响叶轮结构的气动性能。在CAD工程师对叶轮结构进行模型设计时,通常在进行气动性能分析时,工程师是将压气机叶轮作为刚体来处理,因此设计性能是一种理想化的结果。而实际工况下叶片形状往往受叶轮刚度、转速等因素的影响,从而使得实际性能与设计性能发生偏差,实际模型相较于设计模型偏差越大,造成的实际气动性能与理想设计性能的不符合程度越大,这是工程中不希望看到的。气动性能受叶轮叶片的径向变形和周向变形影响较大,而周向变形对气动性能的影响又占主导[27]。路径1上各位置变形值的极差反映了叶片任意两位置变形的最大差值,当极差值越大时,意味着叶片的畸变越大。图 15显示了路径1上径向变形和周向变形的极差值,其中无填充轮的填充率为0,实心轮的填充率为1。如图 15(a)所示,路径1径向变形极差值随填充率的增大逐渐降低,这意味着填充率为0的无填充轮的畸变最大,而随着填充率的增大,畸变呈非线性减小,填充率为1的实心轮的畸变最小,这在理论上是符合实际情况的。对于路径1的周向变形极差值,则并非是随着填充率的增大,极差值降低,而是点阵轮的极差值均小于无填充叶轮和实心叶轮。这意味着在叶轮实际运行时,点阵轮周向变形产生的畸变反而会小于实心轮和无填充轮。在本文工况下,当点阵填充率为0.4时,叶轮的周向变形极差为0.019 31 mm,此时其径向变形极差为0.013 64 mm,此时其周向变形极差值比实心轮约低23.25%,比无填充轮约低55.46%。当点阵填充率为0.6时,叶轮的周向变形极差为0.021 17 mm,此时其径向变形极差为0.011 58 mm,此时其周向变形极差值比实心轮约低15.86%,比无填充轮约低51.16%,这意味着点阵轮实际运行时将比实心轮和无填充轮具有更符合理论设计的气动性能。

图 15 路径1径向和周向变形极差 Fig. 15 Range between radial deformation and circumferential deformation of path 1

图 16显示了无填充轮、实心轮和晶格填充率为0.4的点阵轮大叶片的周向变形和径向变形云图。为便于比较,本文将各个叶轮的仿真结果标尺调为一致。其中,图 16(a)16(b)16(c) 分别表示无填充轮、点阵轮和实心轮大叶片的径向变形,图 16(d)16(e)16(f)分别表示无填充轮、点阵轮和实心轮的周向变形。对比图 16中各图,可以发现点阵轮大叶片发生最大径向变形和周向变形的位置与实心轮近似,二者的最大径向变形均位于大叶片底缘,而最大周向变形均位于叶片外缘处。无填充轮的最大径向变形和最大周向变形与前二者相比有所差异,无填充轮径向变形较大的区域在底缘分布更广,而产生的最大径向变形分布也更靠近大叶片上缘叶梢的位置。

图 16 无填充轮、实心轮和单晶格填充率为0.4的点阵轮大叶片周向和径向变形 Fig. 16 Circumferential deformation and radial deformation of unfilled impeller, solid impeller, and lattice impeller with filling rate of 0.4

图 17描述了不同叶轮的von-Mises应力在路径1上的分布情况,其中图 17(a)为路径1的径向应力分布,图 17(b)为路径1的周向应力分布,图 17(c)为路径1的轴向应力分布,其中拉应力为正,压应力为负。图 17(c)上的应力分布,代表了叶片上缘受拉或者受压的受力状态。对于轴向,无填充轮和填充率为0.2的点阵轮叶片在20~30 mm段受拉应力,在25 mm处无填充轮拉应力最大为135.32 MPa,填充率为0.2的点阵轮在相同位置拉应力比无填充轮约低56.94%。实心轮和填充率分别为0.8,0.6和0.4的点阵轮在相同段处于压应力状态,实心轮压应力最大,最大值为71.137 MPa,而填充率为0.8,0.6和0.4的点阵轮压应力最大时相较于实心轮的压力最大值仍比其低7.50%,24.69%和46.29%。对于薄叶片而言,过大的压应力容易使叶片产生屈曲失稳,所以点阵填充率在0.2~0.4时叶片的受力状态更有利,优异的轴向受力性能可为设计更高转速的压气机叶轮提供强度保证。对比图 17中路径1各方向的von-Mises应力分布情况,发现点阵轮的应力分布与形变分布有着相似的分布规律:点阵轮各个方向上的应力均位于无填充叶轮和实心叶轮之间,并且随着点阵填充率由0.2~0.8的等差增大,各个方向的应力分布曲线逐渐由无填充叶轮向实心轮趋近,并且这种趋近是非线性的,随着点阵填充率更加接近于实心轮,趋近的速度逐渐变小。这种分布规律恰好能验证点阵轮的刚度确实是位于无填充轮和实心轮之间的。保证叶轮在合理的刚度条件下,最大化的实现轻量化设计,在工程上是需要的。

图 17 路径1应力 Fig. 17 Stress of path 1

路径2位于叶轮轮毂小轴孔处,由于实际转动时叶轮不产生周向和径向位移,但在“惯性力”作用下,该位置承受了较大荷载,因此,只研究路径2的von-Mises应力分布情况。图 18描述了各种叶轮在路径2上的von-Mises应力分布。其中,图 18(a)18(b)18(c)分别表示路径2上的径向、周向和轴向应力。对比图 18中各图,发现在路径2上,径向应力占主导。3个方向上的应力,仍然具有与路径1相似的规律:点阵轮的应力分布介于无填充叶轮和实心叶轮之间,随着点阵填充率的增大,叶轮在路径2上的应力分布曲线逐渐向实心轮靠近。不同的是,这种趋近相对于路径1更接近于线性。

图 18 路径2应力 Fig. 18 Stress of path 2
4 结论

1) 在单胞填充率为0.2、0.4、0.6和0.8,单胞尺寸为2 mm下,点阵轮的质量相对于实心轮可分别降低23.69%、17.08%、12.16%和6.53%。在对3种压气机叶轮进行有限元分析后,发现无论是点阵轮、无填充轮还是实心轮,三者的最大变形均位于叶轮大叶片外缘处,且变形沿径向向外逐渐增大,最大应力则位于叶轮轮毂小轴孔处。对于单胞填充率为0.4的点阵轮,其内部点阵域的变形也满足沿径向向外逐渐增大的规律,点阵域的最大变形幅值为0.024 2 mm,点阵域的最大应力幅值为240.75 MPa。

2) 路径1上3个方向的变形或应力分布以及路径2上的应力分布,点阵轮均位于无填充轮与实心轮之间,随着点阵填充率增大,点阵轮的应力或变形分布曲线逐渐由趋近无填充轮向实心轮趋近。这意味着点阵轮的刚度与变形可以通过改变点阵填充率在实心轮和无填充轮之间实现自由调控,即存在这种可能性,通过改变点阵填充率,使叶轮性能满足特定工况下的需要。

3) 各填充率下点阵轮的周向变形极差值均小于无填充轮和实心轮。对于填充率为0.4的情况,其周向变形极差值比实心轮约低23.25%,比无填充轮约低55.46%。这意味着点阵轮将具有更好的周向抗畸变能力和更高的工作效率。同时,对于轴向应力,无填充轮和填充率为0.2的点阵轮叶片在20~30 mm段受拉应力,填充率为0.2的点阵轮在相同位置拉应力比无填充轮约低56.94%。实心轮和填充率分别为0.8,0.6和0.4的点阵轮在相同段处于压应力状态,点阵轮压应力最大值比实心轮压力最大值分别低约7.50%,24.69%和46.29%。对于薄叶片而言,过大的压应力容易使叶片产生屈曲失稳,所以点阵填充率在0.2~0.4时叶片的受力状态更有利。

参考文献
[1]
SHAPIRO A A, BORGONIA J P, CHEN Q N, et al. Additive manufacturing for aerospace flight applications[J]. Journal of Spacecraft & Rockets, 2016, 53(5): 1. DOI:10.2514/1.A33544
[2]
HUO Fali, YANG Deqing, ZHAO Yinzhi. Vibration reduction design with hybrid structures and topology optimization[J]. Polish Maritime Research, 2016, 23(S1): 10. DOI:10.1515/pomr-2016-0040
[3]
WALTON D, MOZTARZADEH H. Design and development of an additive manufactured component by topology optimisation[J]. Procedia CIRP, 2017, 60: 205. DOI:10.1016/j.procir.2017.03.027
[4]
ANDREA M, ANDREA G, ELENA B, et al. Weight reduction by topology optimization of an engine subframe mount, designed for additive manufacturing production[J]. Materialstoday: Proceedings, 2019, 19(3): 1014. DOI:10.1016/j.matpr.2019.08.015
[5]
ABDI M, ASHCROFT I, WILDMAN R D. Design optimization for an additively manufactured automotivecomponent[J]. International Journal of Powertrains, 2018, 7(1/2/3): 142. DOI:10.1504/IJPT.2018.090371
[6]
JIA Dejun, LI Fanchun, ZHANG Cong, et al. Design and simulation analysis of trimaran bulkhead based on topological optimization[J]. Ocean Engineering, 2019, 191: 106304. DOI:10.1016/j.oceaneng.2019.106304
[7]
JIA Dejun, LI Fanchun, ZHANG Cong, et al. Design of bulkhead reinforcement of trimaran based on topological optimization[J]. Ocean Engineering, 2019, 191: 106498. DOI:10.1016/j.oceaneng.2019.106498
[8]
DU PLESSIS A, BROECKHOVEN C, YADROITSAVA I, et al. Beautiful and functional: A review of biomimetic design in additive manufacturing[J]. Additive Manufacturing, 2019, 27: 408. DOI:10.1016/j.addma.2019.03.033
[9]
JUNK S, KLERCH B, NASDALA L, et al. Topology optimization for additive manufacturing using a component of a humanoid robot[J]. Procedia CIRP, 2018(70): 102. DOI:10.1016/j.procir.2018.03.270
[10]
WALSH W R, PELLETIER M H, WANG Tian, et al. Does implantation site influence bone ingrowth into 3D printed porous implants?[J]. The Spine Journal, 2019, 19(11): 1885. DOI:10.1016/j.spinee.2019.06.020
[11]
WANG X J, XU S Q, ZHOU S W, et al. Topological design and additive manufacturing of porous metals for bone scaffolds and orthopaedic implants: A review[J]. Biomaterials, 2016, 83(1): 127. DOI:10.1016/j.biomaterials.2016.01.012
[12]
JIA Dejun, LI Fanchun, ZHANG Cong, et al. Design and simulation analysis of Lattice bone plate based on finite element method[J]. Mechanics of Advanced Materials & Structures, 2019, 28(13): 1311. DOI:10.1080/15376494.2019.1665759
[13]
ROBBINS J, OWEN S J, CLARK B W, et al. An efficient and scalable approach for generating topologically optimized cellular structures for additive manufacturing[J]. Additive Manufacturing, 2016(12): 296. DOI:10.1016/j.addma.2016.06.013
[14]
DAYNES S, FEIH S, LU Wenfeng, et al. Optimisation of functionally graded lattice structures using isostatic lines[J]. Materials & Design, 2017(127): 215. DOI:10.1016/j.matdes.2017.04.082
[15]
LYNCH M E, MORDASKY M, CHENG Lin, et al. Design, testing, and mechanical behavior of additively manufactured casing with optimized lattice structure[J]. Additive Manufacturing, 2018(22): 462. DOI:10.1016/j.addma.2018.05.021
[16]
CHENG Lin, LIANG Xuan, BAI Jiaxi, et al. On utilizing topology optimization to design support structure to prevent residual stress induced build failure in laser powder bed metal additive manufacturing[J]. Additive Manufacturing, 2019(27): 290. DOI:10.1016/j.addma.2019.03.001
[17]
MELI E, FURFERI R, RINDI A, et al. A general framework for designing 3D impellers using topology optimization and additive manufacturing[J]. IEEE Access, 2020(8): 60259. DOI:10.1109/ACCESS.2020.2982841
[18]
KIM E H, CHOI H H, JUNG Y G. Fabrication of a ceramic core for an impeller blade using a 3D printing technique and inorganic binder[J]. Journal of Manufacturing Processes, 2020(53): 43. DOI:10.1016/j.jmapro.2020.01.055
[19]
DRAGHICI S, VINTILA I S, MIHALACHE R, et al. Design and fabrication of thermoplastic moulds for manufacturing CFRP composite impeller blades[J]. Materiale Plastice, 2020, 57(1): 290. DOI:10.37358/MP.20.1.5338
[20]
JIA Dejun, LI Fanchun, ZHANG Yuan. 3D-printing process design of lattice compressor impeller based on residual stress and deformation[J]. Scientific Reports, 2020(10): 600. DOI:10.1038/s41598-019-57131-1
[21]
GUO Honghu, AIKIHIRO T, MASANORI H, et al. Finite element simulation of the compressive response of additively manufactured lattice structures with large diameters[J]. Computational Materials Science, 2020(175): 109610. DOI:10.1016/j.commatsci.2020.109610
[22]
HASSANI B, HINTON E. A review of homogenization and topology optimization I—homogenization theory for media with periodic structure[J]. Computers & Structures, 1998(69): 707. DOI:10.1016/S0045-7949(98)00131-X
[23]
SONG W, KRISHNASWAMY V, PUCHA R V, et al. Computational homogenization in RVE models with material periodic conditions for CNT polymer composites[J]. Composite Structures, 2016(137): 9. DOI:10.1016/j.compstruct.2015.11.013
[24]
CHENG Gengdong, CAI Yuanwu, XU Liang. Novel implementation of homogenization method to predict effective properties of periodic materials[J]. Acta Mechanica Sinica, 2013, 29(4): 550. DOI:10.1007/s10409-013-0043-0
[25]
CHENG Lin, BAI Jiaxi, ALBERT C T. Functionally graded lattice structure topology optimization for the design of additive manufactured components with stress constraints[J]. Computer Methods in Applied Mechanics and Engineering, 2019(344): 334. DOI:10.1016/j.cma.2018.10.010
[26]
KALAMKAROV A L, ANDRIANOV I V, DANISHEVS'K Y Y, et al. Asymptotic homogenization of composite materials and structures[J]. Applied Mechanics Reviews, 2009(62): 030802-1. DOI:10.1115/1.3090830;
[27]
郑赟, 田晓, 杨慧. 跨声速风扇叶片变形对气动性能的影响[J]. 航空动力学报, 2011, 26(7): 1621.
ZHENG Yun, TIAN Xiao, YANG Hui. Influence of transonic fan blade deformation on aerodynamic performance[J]. Journal of Aerospace Power, 2011, 26(7): 1621. DOI:10.13224/j.cnki.jasp.2011.07.028