0 引　言

对于舵叶、机翼、螺旋桨、风机叶片等翼型结构，在工作载荷作用下，其主要变形模式为沿翼展方向发生弯曲变形，因此剖面惯性矩对结构的强度刚度性能极为重要。为了提高剖面惯性矩，往往对翼蒙皮离中性轴距离较远的中间区段进行局部加厚，形成与风机叶片类似的梁帽结构^{[1 – 2]}，以获得较高的材料利用率。如何精确确定中间加厚段的范围，以最少的材料获得较优的刚度和强度，是一个需要研究的问题。以往的研究主要借助有限元软件结合优化算法对翼型结构进行优化设计^{[3 – 6]}，其过程复杂，计算成本较高，若能提出一种简单的工程优化方法，会大大节约设计初期成本。

本文首先对变厚度剖面的比惯性矩（剖面惯性矩与面积的比值）进行理论推导，分析了影响比惯性矩的主要因素；然后针对复杂翼型剖面进行离散处理，采用数值分析方法编程计算，分析了中间加强区段的厚度及范围对比惯性矩的影响规律，根据最大比惯性矩提出了最优加厚区间概念，得到了最优分点位置坐标；在此基础上，针对常见的船用舵翼型NACA0012～NACA0025^[7]，总结了最优分点位置并绘制成图谱，以便为舵翼变厚度设计提供参考。本文所提出的方法简单有效，便于工程应用，尤其适用于复合材料舵翼蒙皮初期设计阶段，对于含加强筋的蒙皮，亦可通过将加强筋等效为蒙皮厚度^{[8 – 9]}的方式进行设计。

1 比惯性矩及其影响因素

典型的翼剖面型线如图1所示，因对称只取其一半，BC为剖面中线，A为剖面最高点。假设EAF为中间加厚段，其厚度为t₁，BE和CF为两边未加厚段，其厚度为t₂。为了计算方便，取翼型线最高点A作为原点，并以该点为分界点将翼型线分为AB与AC两段曲线，2段曲线分别采用共原点的2个独立坐标系。

以AB段为例（AC段计算方法相同），假设AB段的型线函数 $y = {f_1}(x),x \in [0,{x_B}]$ ，计算剖面对于中线BC的比惯性矩如下：

AE加厚段长度、面积及惯性矩分别为

${l_{AE}} = \int_0^{x_1} {\sqrt {1 + y{'^2}} } {\rm{d}}x\text{，}$

${s_{AE}} = {l_1} \cdot {t_1}\text{，}$

${j_{AE}} = \int_0^{{x_1}} {\sqrt {1 + y{'^2}} } \cdot {t_1} \cdot {({y_B} - y)^2}{\rm{d}}x\text{。}$

EB未加厚段长度、面积及惯性矩分别为

${l_{EB}} = \int_{{x_1}}^{{x_B}} {\sqrt {1 + y{'^2}} } {\rm{d}}x\text{，}$

${s_{EB}} = {l_2} \cdot {t_2}\text{，}$

${j_{EB}} = \int_{{x_1}}^{{x_B}} {\sqrt {1 + y{'^2}} } \cdot {t_2} \cdot {({y_B} - y)^2}{\rm{d}}x\text{。}$

则AB段比惯性矩为

$\begin{split}&j/s = \frac{{{j_1} + {j_2}}}{{{s_1} + {s_2}}} =\\&\frac{{{t_1}/{t_2} \cdot \int_0^{{x_1}} {\sqrt {1 + y{'^2}} } \cdot {{({y_B} - y)}^2}{\rm{d}}x + \int_{{x_1}}^{{x_B}} {\sqrt {1 + y{'^2}} } \cdot {{({y_B} - y)}^2}{\rm{d}}x}}{{{t_1}/{t_2} \cdot \int_0^{{x_1}} {\sqrt {1 + y{'^2}} } {\rm{d}}x + \int_{{x_1}}^{{x_B}} {\sqrt {1 + y{'^2}} } {\rm{d}}x}}\end{split}\text{。}$

从式中可以看出，比惯性矩只和 ${t_1}/{t_2}$ 及 ${x_1}$ 有关。然而对于复杂翼面曲线，通常很难给出具体函数表达式，即便拟合成函数，也是极为复杂的函数，从而使比惯性矩的最值求解变得困难。因此，可采用数值分析思想，通过几个已知型线点进行插值，将翼型线离散化，通过分析比惯性矩随 ${t_1}/{t_2}$ 及 ${x_1}$ 的变化规律，得到最优比惯性矩及其对应的分点坐标。

假设将翼型线AB段离散成N段，节点坐标依次为 $({x_i},{y_i}),0 \leqslant i \leqslant N$ ，1~M段为加厚段，则每小段长度为：

${d_i} = \sqrt {{{({x_i} - {x_{i - 1}})}^2} + {{({y_i} - {y_{i - 1}})}^2}} \text{，}$

AE段面积与惯性矩为

${s_A} = \sum\limits_{i = 1}^M {{d_i}} \cdot {t_1}\text{，}$

${j_A} = \sum\limits_{i = 1}^M {{d_i} \cdot {t_1} \cdot {{\Bigr({y_{\max }} - \frac{{{y_i} + {y_{i - 1}}}}{2}\Bigr)}^2}} \text{。}$

EB段面积与惯性矩为

${s_B} = \sum\limits_{i = M + 1}^N {{d_i}} \cdot {t_2}\text{，}$

${j_B} = \sum\limits_{i = M + 1}^N {{d_i} \cdot {t_2} \cdot {{\Bigr({y_{\max }} - \frac{{{y_i} + {y_{i - 1}}}}{2}\Bigr)}^2}} \text{。}$

则AB段比惯性矩为

$j/s = \frac{{{t_1}/{t_2} \cdot \sum\limits_{i = 1}^M {{d_i} \cdot {{\Bigr({y_B} - \frac{{{y_i} + {y_{i - 1}}}}{2}\Bigr)}^2}} + \sum\limits_{i = M + 1}^N {{d_i} \cdot {{\Bigr({y_B} - \frac{{{y_i} + {y_{i - 1}}}}{2}\Bigr)}^2}} }}{{{t_1}/{t_2} \cdot \sum\limits_{i = 1}^M {{d_i}} + \sum\limits_{i = M + 1}^N {{d_i}} }}\text{。}$

2 比惯性矩变化规律研究

以复合材料船用舵为研究对象，对翼剖面比惯性矩随 ${t_1}/{t_2}$ ， ${x_1}$ 及 ${x_2}$ 的变化规律进行研究。船用舵常见的翼型为NACA系列，其剖面型值如表1，其中b为弦长，e为剖面最大厚度。

首先对特定翼型特定厚度比的情况进行研究，采用第1节中的数值公式编程计算，代表翼型为NACA0025，计算时取弦长为1 m，注意将型值转换为图1坐标系中对应坐标。经过计算，绘制翼型NACA0025在特定厚度比（ ${t_1}/{t_2} = 4$ ）时，比惯性矩随不同分点E，F位置的变化曲线如图2所示，图中x₁为AB段的分点E的横坐标，x₂为AC段的分点F的横坐标，由图可知：随着分点E或F横坐标的增大，比惯性矩的变化趋势均为先增大后减小，其间存在的极大值点为最优的比惯性矩，其对应的分点位置为最佳的分点位置。

当厚度比 ${t_1}/{t_2}$ 变化时，最优分点位置及最优比惯性矩会发生改变，将不同厚度比情况下，最优比惯性矩的变化规律绘制成曲线如图3所示，图中 ${t_1}/{t_2}$ 变化范围为[2，20]，由图可知：随着 ${t_1}/{t_2}$ 的增大，AB段与AC的最优比惯性矩变化规律类似，均为逐渐增大，但增速变缓。这说明加厚段越厚效果越好，但实际中，由于工艺或其他因素对最大厚度有所限制，而且厚度比太大也会带来厚薄蒙皮的过渡等问题。

3 蒙皮厚度优化方法与分点图谱

在应用该方法确定加厚区间时，可先根据局部强度要求及最小工艺厚度要求确定最小厚度t₂；由以上分析可知，厚度比 ${t_1}/{t_2}$ 与最优比惯性矩成一一对应关系，而且当t₂确定以后，最优比惯性矩对应的惯性矩随即确定，其与 ${t_1}/{t_2}$ 的关系如图4所示，根据结构所需的最小惯性矩要求由图4可确定 ${t_1}/{t_2}$ ；继而可根据图1确定该厚度比下满足最大比惯性矩的最优分点坐标x₁和x₂。

实际上，船用舵以NACA0012-NACA0025翼型使用居多，因为在该范围内，舵叶的升力系数较大，超出该范围，升力系数较小^[7]。此外，对于复合材料舵，蒙皮最小工艺厚度约5 mm，最大厚度受工艺与声学性能限制一般不超过30 mm，因此蒙皮厚度比范围一般为2～6。因此，针对该范围内NACA0012-NACA0025翼型的最优比惯性矩进行研究，并将对应的最优分点位置x₁和x₂绘制成图谱，以供设计参考，如图5所示。分析可得出以下结论：1）对于弦长1 m的常见船用舵翼型，加厚段的最优分点位置坐标x₁的分布范围为[130 mm,174 mm]，x₂的分布范围为[220 mm, 297 mm]；2）随着厚度比 ${t_1}/{t_2}$ 的增大，最佳分点坐标x₁和x₂逐渐变小，即分点往中间移动，加厚段范围变小；3）随着翼型系数e/b的增大，x₁逐渐变大，x₂略有增大，即分点往两边移动，加厚段范围变大；4）随着翼型系数e/b的增大，同一厚度比的分点坐标x₁和x₂基本呈现线性变化，这说明对于变截面变翼型（厚度比不变）的船用舵，根据分点坐标的线性变化，只需确定两端面的最优分点位置，即可线性确定整个舵剖面的分点，极大地扩展了该方法的实用性和适用范围。

4 算例分析 4.1 算例描述

某船用舵为变截面变翼型舵，其展长为3 200 mm，大端面翼型NACA0018，弦长2 400 mm，小端面翼型NACA0016，弦长1 600 mm。采用玻璃钢复合材料设计，根据惯性矩要求确定优化前的蒙皮平均厚度为12.8 mm。在保证蒙皮质量不变的前提下，根据上述加厚区间优化设计方法，首先由局部强度限制两边未加厚段的厚度为6 mm，再由 ${t_1}/{t_2}$ 与惯性矩的曲线图，结合剖面惯性矩要求可确定加厚段厚度为24 mm，因此厚度比为4，此时分点坐标可查图5确定，并注意按弦长比例进行缩放，得到分点坐标为：

大端面分点坐标

x₁=146.0×2.4=350.4（mm），x₂=247.5×2.4=594.0（mm）；

小端面分点坐标

x₁=144.6×1.6=231.4（mm），x₂=247.2×1.6=395.5（mm）。

4.2 结果分析

通过Abaqus有限元软件，对优化前后的舵叶进行建模，并在舵轴端施加固支边界约束，在舵叶表面施加操舵水动力载荷，如图6所示。计算所得位移和蒙皮应力分布如图7和图8所示。优化前最大位移为38.95 mm，优化后最大位移为34.25 mm，减小了12.1%；优化前最大应力为51.74MPa，优化后最大应力为41.40MPa，减小了20.0%。从优化结果来看，优化后舵叶的强度与刚度都得到了较大提高，说明该方法的优化效果比较明显。

5 结　语

本文基于比惯性矩最大原则，提出了一种确定翼型结构最优加厚区间的工程优化方法，在此基础上将常用舵翼的最优加厚区间绘制成图谱，基于图谱分析了最优加厚区间变化规律，并通过一算例阐述了图谱使用方法，验证了该方法的有效性。经过研究，本文得出的主要结论如下：

1）在翼型与加厚比 ${t_1}/{t_2}$ 一定的情况下，随着分点E、F位置坐标的增大，比惯性矩先增大或减小，其间存在最优分点位置使比惯性矩最大。

2）随着厚度比 ${t_1}/{t_2}$ 的增大，最佳分点坐标x₁和x₂逐渐变小，即分点往中间移动，加厚段范围变小；随着翼型系数e/b的增大，x₁和x₂逐渐变大，即分点往两边移动，加厚段范围变大。

3）随着翼型系数e/b的增大，同一厚度比的分点坐标x₁，x₂呈现线性变化，这说明对于变截面变翼型（厚度比不变）的船用舵，只需确定两端面的最优分点位置，即可线性确定其他舵剖面的分点，极大地扩展了该方法的实用性和适用范围。