本发明属于复合材料的有限元分析技术领域，尤其涉及一种变刚度复合材料的建模方法。
背景技术：
由于对轻质结构的需求，飞行器、汽车等结构设计的研究和发展一直与采用性能优越的新材料密切相关。先进复合材料具有比强度和比模量高、性能可设计和易于整体成形等许多优异特性，将其用于飞行器等结构上，可明显减轻结构重量、改善气动弹性特性，提高结构性能。先进复合材料在飞行器上应用的部位和用量的多少，现已成为衡量飞机结构先进性的重要指标之一。例如，新型客机波音787与空客A350的复合材料用量达到50％以上。
复合材料结构设计是复合材料在结构上科学合理应用和降低全寿命周期费用的关键，重点在于充分利用复合材料性能的方向性、结构性能的可设计性和大型构件整体成形的良好工艺性，实现结构、效率、性能、功能与成本的综合优化。
带有连续变化角度的纤维铺层的复合材料被称为变刚度复合材料，与传统的复合材料层合板相比，其重量和成本可以减少10％-30％，因为不断变化的纤维取向导致了每层的刚度在不同的空间位置各不相同。工程设计人员可据此调整内部载荷分布，使结构性能得到显著提高。变刚度复合材料是在已有主要材料体系基础上开发的先进复合材料，突破了传统的复合材料铺层设计理念，对制品的适应性更强，可以实现凸面、凹面和双曲率构件等复杂型面结构的铺叠，极大减少连接件的使用，减少结构重量。同时在使用连续变角度纤维铺层时可预留开孔结构，减少后期开孔，既可避免应力集中现象，提高结构性能，降低成本，又可以进行优化设计以及实现设计制造一体化，提高结构设计与生产效率，从而实现降低复合材料构件的全生命周期成本。
对于传统复合材料结构设计，一般采用有限元软件(如ABAQUS)中的复合材料模块进行铺层方向、铺层顺序、铺层层数的设计，该种方法中各个铺层中的纤维方向角为常数。然而，由于变刚度复合材料中各个铺层中的纤维方向角为连续变化的函数，而不是一个常数，所以该种方法不适用于变刚度复合材料。目前对变刚度复合材料的研究局限于平板、圆柱、圆锥面等规则的几何体。
因此，如何合理地将变刚度复合材料应用到不规则的复杂结构中，并对其力学性能和稳定性进行研究，十分迫切且具有重要意义。
技术实现要素：
有鉴于此，本发明提出了一种变刚度复合材料的建模方法，解决了将变刚度复合材料应用到机翼等复杂或不规则结构中的难题，具有广泛的推广意义。
本发明是通过以下技术方案来实现的：一种变刚度复合材料的建模方法，包括以下步骤：
S1、先对常刚度复合材料的结构模型进行优化设计；
S2、对步骤S1中的模型进行截面简化；
S3、根据步骤S2中的简化模型进行蒙皮铺层建模；
S4、对变刚度复合材料结构模型建立整体有限元分析模型；
进一步，该方法还包括步骤：
S5、对步骤S4中建立的整体有限元分析模型进行优化设计。
进一步，步骤S1中是通过求解数学模型的方式对所述的常刚度复合材料的结构模型进行优化设计，数学模型表达式如下：
min W(X)＝∑ρj*Vj(X)；
式中优化设计的目标是变刚度复合材料结构模型的质量W(X)最小；ρj为第j部分的材料密度；Vj(X)为第j部分的体积。
进一步，步骤2中，还包括对简化模型与原模型的力学性能进行分析，验证简化模型的有效性。
进一步，步骤S3中进行蒙皮铺层建模的方法包括：通过蒙皮中每个单元网格的几何数据，进行对比计算，判断出这个单元在哪一层，再使用相应的纤维路径函数给与该单元纤维方向的数据信息。
进一步，步骤S3中，通过判断蒙皮内所有点与芯体表面的距离确定每层蒙皮所在单层的位置。
进一步，步骤S4中，对变刚度复合材料结构模型建立整体有限元分析模型的方法，包括：蒙皮与芯体使用8节点六面体非协调单元，并在上蒙皮的截面相交区域进行了加密处理。
进一步，步骤S4中，还包括：对模型根弦截面施加固定约束以及在上下蒙皮施加与最优点相同的气动载荷。
进一步，步骤S5中，对步骤S4中建立的整体有限元分析建模进行优化设计，包括对变刚度复合材料模型的屈曲特性进行优化设计，并将纤维方向角函数的斜率与初始值作为设计变量。
进一步，步骤S5中，用基于梯度的序列二次规划法，并选取多个初始点进行优化计算。
本发明所提供的一种变刚度复合材料的建模方法，通过组合分步优化设计的方法，在第一步优化设计中先得到最优几何构型与材料参数，在最优几何构型上进行合理简化，实现了复杂结构的分层建模，然后在优化设计中通过包含复合判据的蒙皮铺层算法将纤维路径信息输入到有限元分析模型中，以变刚度复合材料的纤维铺设路径参数为设计变量，实现了变刚度复合材料结构模型的优化设计。解决了将变刚度复合材料应用到机翼等复杂结构中的难题。
附图说明
图1为一种机翼模型结构示意图；
图2为图1所示变刚度复合材料机翼分层建模的示意图；
图3为图1所示机翼模型的简化翼盒模型示意图；
图4为纤维铺层算法的一种具体实施方式的流程图；
图5为一种最优变刚度复合材料机翼蒙皮各个单层的纤维方向角分布示意图。
其中，图中标记表示为：1、机翼模型，11、芯体，12、蒙皮。
具体实施方式
下面将参考附图并结合实施例来详细说明本发明。以下结合附图对本发明的原理和特征进行描述，需要说明的是，在不冲突的情况下，本申请中的实施例及实施例中的特征可以相互组合。所举实例只用于解释本发明，并非用于限定本发明的范围。
以一种机翼模型1为例，请参考附图1，为一种机翼模型结构示意图。
先对常刚度复合材料的机翼结构模型进行优化设计。
复合材料机翼优化设计的数学模型表达如下：
min W(X)＝∑ρj*Vj(X)
s.t.0.5≤X1≤0.99,8≤X2≤15,0≤X5≤90,0.01≤X6≤0.9
0.01≤X7≤0.95,0.1≤X8≤30,0.1≤X9≤30
X3＝1,2,3,4,5,6,X4＝1,2,3,4,5
|σ1|≤X',|σ2|≤Y',|τ12|≤S
MaxDis≤10mm,MaxRot≤1°
式中优化设计的目标是机翼的质量W(X)最小；ρj为第j部分的材料密度；Vj(X)为第j部分的体积；设计变量X＝{X1,X2,X3,X4,X5,X6,X7,X8,X9}，如表1所示；σ1、σ2为机翼单层出现的最大正应力，τ12为机翼出现的最大剪切应力；X'、Y'、S为对应的材料强度；MaxDis为机翼最大变形量；MaxRot为机翼最大扭转角。
表1优化设计变量
表2复合材料机翼优化结果
根据上述优化的结构，对机翼模型1截面进行简化，并在相同载荷条件下，对简化模型与原模型的力学性能进行分析，验证简化模型的有效性，在后续的变刚度复合材料机翼的设计中，采用该简化模型，如图3所示。
由于已将翼型简化为由7条直线组成的平面，因此芯体截面可由A-J共10个控制点进行表征。将坐标原点O设为翼盒固定约束端的左下侧，建立直角坐标系Oxyz。其中Oxz平面与翼盒固定约束端重合，Oz轴与重合，坐标系遵守右手定则。
蒙皮12铺层角度取决于蒙皮12所在的单层的位置，其结构示意图如图2所示。由于单层蒙皮12厚度相同，因此可通过判断蒙皮12内所有点与芯体表面的距离作为依据进行判断。对于该模型而言，截面AOG对截面HIJ的放大倍数为1.5，这表明不能简单地通过x坐标以及z坐标来进行蒙皮12所在层位置的判定，而必须通过复合判据进行控制。以面ABH与面BCH附近的点P为例进行阐述，求解思想如图4中的流程图所示。
该算法的原理就是在给每一个单元赋予材料特性的时候，先通过该单元的几何数据进行对比计算，判断出这个单元在哪一层，再使用相应的纤维路径函数赋予该单元纤维方向的相关信息。
机翼的蒙皮12与芯体11使用8节点六面体(hex)非协调单元，并在上蒙皮12的截面相交区域进行了加密处理。蒙皮12与芯体11材料分别选用与上述优化设计得到的最优解相同材料。在机翼根部施加固定约束，在机翼上下表面施加与最优点相同的气动载荷。通过上节内容所研究的复合判据，判断蒙皮12内每个单层所在的位置，再通过用户子程序给每个单层赋予曲线纤维的方向角信息。
在得到的最优机翼模型1的基础上，可以对变刚度复合材料机翼的屈曲特性进行优化。在优化设计过程中，主要研究纤维方向角函数对变刚度机翼屈曲性能的影响。稍根比、展弦比、碳布的铺层层数与蒙皮材料等参数均直接使用上一步最优解对应的值，不再作为设计变量，只将纤维方向角函数的斜率与初始值作为设计变量。
变刚度复合材料的曲线纤维的方向角θ为沿机翼弦向方向线性变化的函数：
其中特征长度d为机翼根弦长度，T0为机翼根部最左端(x＝0)的纤维方向角，T1为机翼根部最右端(x＝d)的纤维方向角。为了方便后续表述，将上述公式写为θ(x)＝kx+b，其中k＝2(T1-T0)/d，b＝T0，则纤维参考路径可定义为：
因为纤维方向角函数的斜率与初始值都由T0和T1决定，所以直接把T0和T1作为优化设计变量。变刚度复合材料机翼优化设计的数学模型表达如下：
max f(X)
s.t.-90≤X1≤90
-90≤X2≤90
其中X＝{X1,X2}为设计变量，X1为机翼根部最左端的纤维方向角T0，X2为机翼根部最右端的纤维方向角T1，优化目标函数f(X)为变刚度复合材料机翼的屈曲状态的特征值，约束条件为T0和T1的取值范围。
在对变刚度复合材料机翼的研究中，选用与常刚度复合材料机翼相同的材料特性与几何参数，施加相同的载荷与边界条件，上下表面蒙皮选择三层正交铺层结构[<T0|T1>/<(T0+90)|(T1+90)>/<T0|T1>]。在保持重量不变的前提下，得到临界屈曲载荷最大的变刚度复合材料机翼。由于优化设计变量T0和T1都为连续变量，所以直接选用收敛速度快、精度高的基于梯度的序列二次规划法(SQP)，并选取多个初始点分别进行优化计算，通过比较验证，得到使变刚度复合材料机翼临界屈曲载荷最大的最优设计。
由于给变刚度复合材料机翼的上下表面蒙皮12的曲线纤维赋予方向特性时，需要使用用户子程序来判断纤维所在单层的位置，并通过相应的函数来定义纤维方向，因此需要把用户子程序实现参数化。在优化过程中的每一次迭代计算中，先将设计变量传输到用户子程序，在用户子程序中判断出纤维所在的单层后，按照该设计变量定义的纤维角函数来给曲线纤维赋予方向特性，然后完成变刚度复合材料机翼结构的屈曲分析。
最优设计点为X*<21|-34>，最优解对应的纤维方向角函数如下：
X*：θ(x)＝kx+b＝-0.4074x+21
图5为最优变刚度复合材料机翼蒙皮各个单层的纤维方向角分布，从图中可看出在机翼的任一截面的纤维方向角都沿x方向线性连续变化。常刚度复合材料机翼与最优变刚度复合材料机翼的一阶到四阶临界屈曲特征值如表3所示，变刚度复合材料机翼的一阶临界屈曲特征值为7.1953，与常刚度复合材料机翼的一阶临界屈曲特征值(4.8009)相比，增大了50％，优化效果非常明显。
表3复合材料机翼与变刚度机翼最优解的各阶临界屈曲特征值
综上所述，在本发明所提供的一种变刚度复合材料的建模方法的具体实施例中，由于机翼模型1的蒙皮12为变刚度复合材料层合板，需要对其进行分层建模(如图2所示为变刚度复合材料机翼分层建模的示意图)，然后通过用户子程序ORIENT给每个单层赋予纤维方向特性。而上述的最优机翼模型1为曲面结构，对其进行分层建模难度非常大，需要对机翼模型1的几何截面进行合理的简化，再对机翼的变刚度复合材料蒙皮进行铺层建模，然后建立变刚度复合材料机翼的整体有限元分析模型，最后，再对建立的整体有限元分析建模进行优化设计，对变刚度复合材料机翼的纤维铺放路径进行优化，使其得到最大屈曲承载能力。在调用用户子程序时，需要判断所铺设纤维所在的单层位置，而单层的位置是通过单层中的点到翼盒外表面的距离来确定的。
尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。