一种飞机薄壁件自动钻铆多姿态变形建模方法

摘要

公开了一种飞机薄壁件自动钻铆多姿态变形建模方法。该方法包括步骤：将飞机薄壁件和托板之间的定位夹紧关系简化为连续梁受多个支承座支承的模型；计算各托板的支持力，在所有中间托板位置将薄壁件截开并将托板支承简化为简支支承；求解出各个托板的弯矩后，进一步利用两个托板的受力平衡模型，计算自动钻铆过程中承载飞机薄壁件的各个托板的支持力；将薄壁件按照托板位置划分成相应的小区域，针对每个小区域求解薄壁件变形，其中，每个小区域薄壁件受力变形模型简化为托板位置薄壁件简支支承，薄壁件表面仅受均布载荷的变形计算；最后计算出薄壁件中面沿方向的变形量。

说明书

技术领域

本发明涉及一种飞机薄壁件自动钻铆变形建模方法，主要是针对 飞机薄壁件在自动钻铆过程中多个姿态下求解薄壁件变形的建模方 法。

背景技术

自动钻铆技术是目前飞机装配线上的先进技术之一，可以实现飞 机薄壁件的自动钻铆。自动钻铆系统主要由自动钻铆机和相配的托架 组成。自动钻铆过程中要求钻头或铆钉与钻铆点的法线重合，因此需 要计算托架和薄壁件的变形量并进行相应的补偿。自动钻铆系统中围 框质量大且跨度较长，造成飞机薄壁件在还没有上架之前就因自重而 产生较大的变形挠度，给钻铆点的法向找寻带来一定的困难。同时， 由于托板的数量及分布间距受钻铆空间和钻铆效益的限制，相邻托板 之间必须存在一定的间距。实际钻铆过程中飞机壁板多为大型薄壁 件，薄壁件安装在托板上会因为自重而造成相应的变形，同时自动钻 铆过程中托架的多自由度转动让薄壁件的变形也具有了动态多样性。

目前，较多的研究主要集中在自动钻铆系统中托架围框的截面优 化、托架自重在自动钻铆状态下的变形量计算及补偿方面。文献“自 动钻铆托架围框横梁的优化设计”(《机械制造》，2010年08期)对 自动钻铆托架的围框横梁进行了优化设计，使托架在满足刚度和强度 的要求下，尽可能轻巧，为设计、制造围框提供了较好的依据。文献 “基于CATIA V5的自动钻铆机托架变形研究”(《航空制造技术》， 2008年16期)，研发了基于CATIA V5的托架变形补偿技术及变形 模拟系统，进一步研究了托架变形。目前，对于薄壁件上架定位、夹 紧后因自重引起的变形研究较少，给自动钻铆系统的钻铆精度带来一 定的影响。

发明内容

为了计算薄壁件上架定位夹紧后因自重引起的变形量并验证托 板布局的合理性，本发明公开了一种飞机薄壁件自动钻铆多姿态变形 建模方法。该方法包括如下步骤：

1)将托板支撑壁板简化为连续梁受多个支承座支承的模型， 计算各托板的支撑力和弯矩；

2)利用各个托板位置确定的边界条件，结合壳体理论计算壁 板的变形量。

本发明的有益效果是：利用理论模型求解飞机薄壁件在自动钻铆 系统中定位、夹紧后各个托板的受力。为托板的强度设计、数目以及 在托架上的分布提供一定的依据。同时，计算托板间距已知的情况下 薄壁件的变形量，通过控制最大的变形量来进一步验证托板数量和分 布是否合理。

下面结合附图和实施例对本发明作详细说明。

附图说明

图1为根据本发明的具体实施方式的自动钻铆系统示意图；

图2为根据本发明的具体实施方式的自动钻铆系统简化示意图；

图3为根据本发明的具体实施方式的薄壁件上架后薄壁件与托板 相对位置及薄壁件分区图；

图4为根据本发明的具体实施方式的薄壁件与托板简化模型图；

图5为根据本发明的具体实施方式的中间托板及相邻两跨薄壁件 受力图；

图6为根据本发明的具体实施方式的中壁板分区及受力图；以及

图7为根据本发明的具体实施方式的飞机薄壁件的变形量建模方 法的流程框图。

具体实施方式

在以下的实施方式的详细描述中，参照构成该描述的一部分的附 图进行说明。附图以示例的方式展示出特定的实施方式，本发明被实 现在这些实施方式中。所示出的实施方式不是为了穷尽根据本发明的 所有实施方式。可以理解，其他的实施方式可以被利用，结构性或逻 辑性的改变能够在不脱离本发明的范围的前提下被做出。

为了研究薄壁件上架后在自动钻铆系统多个姿态下的变形，本发 明提供了一种飞机薄壁件自动钻铆多姿态变形建模方法。前提假设托 架的变形量在水平状态已经通过调节托板的位置得到补偿。考虑薄壁 件上架后在自动钻铆多姿态下其自重和托板分布引起薄壁件的动态 变形，建立自动钻铆多姿态下各个托板的受力以及薄壁件变形的建模 方法。

如图1所示，自动钻铆系统10主要由自动钻铆机100和相配的 托架200组成。其中，托架200由端梁202与横梁204焊接而成的围 框206、定位夹紧薄壁件的托板208和固定连接在围框206的两根横 梁204上用于为托板208及附件提供支承作用的角板210组成。自动 钻铆过程中要求钻头或铆钉与钻铆点的法线重合，因此需要计算托架 200和薄壁件的变形量并进行相应的补偿。如前所述，自动钻铆系统 中围框206质量大且跨度较长，造成飞机薄壁件在还没有上架之前就 因自重而产生较大的变形挠度，给钻铆点的法向找寻带来一定的困 难。

下面，详细地叙述本发明的建模方法。

(1)如图2所示，自动钻铆系统中，飞机薄壁件上架后将薄壁 件和托板208之间的定位、夹紧关系简化为连续梁受多个支承座支承 的模型(多姿态依据图2中的α角和β角的改变体现，实际钻铆过程 中β角很小(-4°～+4°)对薄壁件重力在托板208上的分布影响较小， 本文不予以考虑)。其中，参见图2和图3，图2所示的在X′方向延 伸的横梁204与图3所示的飞机薄壁件S的长度对应(图2中)，支 承座支承端梁202位置和托板支承薄壁件位置对应。同时，将端梁 202最左边支承座假定为固定支承，其余支承座假设为简支支承，简 化后模型如图4所示。在模型中各跨受力用均布载荷q_w表示， q_w＝G_i/L_i，其中，G_i表示第i跨薄壁件的重力，L_i表示第i跨的长度。

(2)计算各个托板208的支持力，即，求解各个支承座的支持 力。在所有中间托板208对应位置将薄壁件S截开，并将托板208支 承简化为简支支承。每两个支承之间为一个跨距，每个跨距上的薄壁 件S仅承受直接作用于该跨的载荷以及两端的支点弯矩(一个或两个 多余未知力)。如图5所示，可以很方便地求出薄壁件S在托板208 处的转角，并根据中间支座处相连两截面转角相同的条件建立补充方 程，从而确定全部托板208位置的弯矩：

$\frac{M_{i-1}{l}_i}{I_i{E}_i}+2M_i(\frac{l_i}{I_i{E}_i}+\frac{l_{i+1}}{I_{i+1}{E}_{i+1}})+\frac{M_{i+1}{l}_{i+1}}{I_{i+1}{E}_{i+1}}=-6(\frac{{\omega }_i{a}_i}{I_i{l}_i{E}_i}+\frac{{\omega }_{i+1}{b}_{i+1}}{I_{i+1}{l}_{i+1}{E}_{i+1}})$  (公式1-1)

式中：

l_i表示i-1到i托板208间的跨距m；

M_i表示i托板208支承位置截面上的弯矩KN·m；

a_i表示ω_i的形心c_i到i托板208的距离m；

ω_i表示i-1到i托板208间载荷弯矩图的面积m²；

b_i+1表示ω_i+1的形心c_i+1到i+1支承的距离m；

E_i表示i-1到i托板208间薄壁件材料的弹性模量GPa；

I_i表示i-1到i托板208间薄壁件简化为连续梁的截面惯性矩。

实际中，自动钻时薄壁件S的材料相同，所以薄壁件S依据托板 208划分区域后各区域薄板件的弹性模量E_i相同。

各个托板208处的弯矩满足：

$\frac{M_{i-1}{l}_i}{I_i}+2M_i(\frac{l_i}{I_i}+\frac{l_{i+1}}{I_{i+1}})+\frac{M_{i+1}{l}_{i+1}}{I_{i+1}}=-6(\frac{{\omega }_i{a}_i}{I_i{l}_i}+\frac{{\omega }_{i+1}{b}_{i+1}}{I_{i+1}{l}_{i+1}})$  (公式1-2)

式中：l_i、I_i、ω_i、a_i、b_i+1(i＝0，1，2，3，…，n)为已知条件。n+1个托 板208中有n-1个中间托板，可以列出n-1个方程，外加方程M₀＝0和 M_n＝0，可以求解出n-1个中间托板208所承受的弯矩M_i。

(3)求解出各个托板208的弯矩之后，进一步利用相邻两跨之 间薄壁件S的受力平衡模型，可以计算出各个托板208的支持力：

$F_{i-1}=0.5q_w{l}_i+\frac{M_i-M_{i-1}}{l_i}$(最前端托板i＝0)  (公式1-3)

$F_i=0.5q_w(l_i+l_{i+1})-\frac{M_i-M_{i-1}}{l_i}+\frac{M_{i+1}-M_i}{l_{i+1}}$(中间托板i＝1，2，…，n-1)(公式1-4)

$F_{i-1}=0.5q_w{l}_{i+1}-\frac{M_{i+1}-M_i}{l_{i+1}}$(最后端托板i＝n)  (公式1-5)

根据计算的托板支持力可以适当调整托板208数目、分布位置和 托板结构。

(4)自动钻铆过程中，为了防止薄壁件S在托架200转动过程 中与托板208产生相对运动，在围框上托板208对应位置配以拉紧带 拉紧薄壁件，保证托板208与薄壁件S的紧密贴合。计算薄壁件S上 架之后在自动钻铆各个姿态的变形量，即，求解薄壁件S上架后在多 位姿状态下因自重(计算时简化为均布力q_n)在薄壁件S法向的分量 不同而引起的动态变形量。实际钻铆过程中，飞机薄壁件S多为带曲 率的曲面，可以看作薄的壳体处理。计算薄壳体受复杂力的变形，目 前还较为困难。为了使计算简便再一次简化薄壁件S模型。将薄壁件 S依据托板208位置分割为小的区域，然后，针对每个小区域求解薄 壁件变形。最后，依据计算出的变形量进一步验证托板208的分布是 否合理。简化后每个小区域薄壁件除托架200位置受力外，薄壁件S 仅受自身重力。从而，各个区域薄壁件受力变形计算模型简化为托板 位置薄壁件简支支承，Y′方向为自由边界，薄壁件表面仅受均布载荷 (将薄壁件重力在自动钻铆过程的各个姿态下看成不同大小的均布 载荷)的变形计算如图6所示。其中，薄壁件简支支承边界条件满足：

T₁₁＝0、S＝ε₁₂Eh/2(1+μ)＝0、ω＝0、M₁₁＝0、M₁₂＝M₂₁＝0；

薄壁件自由边界条件满足：

μ₁＝0，μ₂＝0，ω＝0，M₁₁＝M⁰₁₁＝0，φ₂＝0。

建立满足以上边界条件的薄壳体变形计算，设：

$F=\underset{n=1}{\overset{\infty }{\Sigma }}{F}_n\left(y\right)\sin \left(\frac{\mathrm{n\pi x}}{L}\right)$  (公式1-6)

$f=\underset{n=1}{\overset{\infty }{\Sigma }}{f}_n\left(y\right)\cos \left(\frac{\mathrm{n\pi x}}{L}\right)$  (公式1-7)

将公式1-6和公式1-7代入薄壁件的静力方程并将均布载荷展开 成傅里叶级数得到公式1-8和公式1-9：

$\left(\begin{array}{c}{F^{\left(8\right)}}_n-4{\left(\frac{\mathrm{n\pi }}{L}\right)}^2{F^{\left(6\right)}}_n+6{\left(\frac{\mathrm{n\pi }}{L}\right)}^4{F^{\left(4\right)}}_n+[4{\left(\frac{\mathrm{n\pi }}{L}\right)}^6-\frac{1-{\mu }^2}{k{R}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^4]{F^{\left(2\right)}}_n+\\ [{\left(\frac{\mathrm{n\pi }}{L}\right)}^8-\frac{1-{\mu }^2}{k{R}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^6+\frac{1-{\mu }^2}{k{R}^2}{k}_0{\left(\frac{\mathrm{n\pi }}{L}\right)}^4]F_n=\frac{2q_n\cos \alpha \cos \beta [1-\cos \left(\mathrm{n\pi }\right)]}{\mathrm{Kk}(1-\mu )\mathrm{n\pi }}\end{array}\right)$

         (公式1-8)

${f^2}_n-[{\left(\frac{\mathrm{n\pi }}{L}\right)}^2+\frac{{2k}_0}{1-\mu }]f_n=0$  (公式1-9)

q_n^cosαcosβ中的α，β角在自动钻铆过程中因钻铆位置的不同而 改变是已知条件，公式1-8的特解为：

$F_p^n=\frac{{2q}_w[1-\cos \left(\mathrm{n\pi }\right)]}{\mathrm{Kk}(1-\mu )\mathrm{n\pi }[{\left(\frac{\mathrm{n\pi }}{L}\right)}^8-\frac{1-{\mu }^2}{{\mathrm{kR}}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^6+\frac{(1-{\mu }^2)k_0}{{\mathrm{kR}}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^4]}$  (公式1-10)

求解公式1-8的齐次解，设：

$F_n=e^{{\lambda }_{n}y}$  (公式1-11)

上式代入公式1-8得：

$\left(\begin{array}{c}{{\lambda }^8}_n-4{\left(\frac{\mathrm{n\pi }}{L}\right)}^2{{\lambda }^6}_n+6{\left(\frac{\mathrm{n\pi }}{L}\right)}^4{{\lambda }^4}_n+[4{\left(\frac{\mathrm{n\pi }}{L}\right)}^6-\frac{1-{\mu }^2}{k{R}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^4]{{\lambda }^2}_n+\\ {\left(\frac{\mathrm{n\pi }}{L}\right)}^8-\frac{1-{\mu }^2}{k{R}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^6+\frac{1-{\mu }^2}{k{R}^2}{k}_0{\left(\frac{\mathrm{n\pi }}{L}\right)}^4=0\end{array}\right)$(公式1-12)

进一步简化为：

${\lambda }_n^8+{\mathrm{a\lambda }}_n^6+{\mathrm{b\lambda }}_n^4+{\mathrm{c\lambda }}_n^2+d=0$  (公式1-13)

式中，$a=-4{\left(\frac{\mathrm{n\pi }}{L}\right)}^2,b=6{\left(\frac{\mathrm{n\pi }}{L}\right)}^4,c=4{\left(\frac{\mathrm{n\pi }}{L}\right)}^6-\frac{1-{\mu }^2}{k{R}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^4,$$d={\left(\frac{\mathrm{n\pi }}{L}\right)}^8-\frac{1-{\mu }^2}{k{R}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^6+\frac{1-{\mu }^2}{k{R}^2}{k}_0{\left(\frac{\mathrm{n\pi }}{L}\right)}^4.$

根据方程求出特征值λ_n(n＝1，2，…8)，则方程的齐次解为：

$F_n^h=A_1{e}^{{\lambda }_{1}y}+A_2{e}^{{\lambda }_{2}y}+A_3{e}^{{\lambda }_{3}y}+A_4{e}^{{\lambda }_{4}y}+A_5{e}^{{\lambda }_{5}y}+A_6{e}^{{\lambda }_{6}y}+A_7{e}^{{\lambda }_{7}y}+A_8{e}^{{\lambda }_{8}y}$

全解为：

$F_n=F_n^h+\frac{{2q}_n[1-\cos \left(\mathrm{n\pi }\right)]}{\mathrm{Kk}(1-\mu )\mathrm{n\pi }[{\left(\frac{\mathrm{n\pi }}{L}\right)}^8-\frac{{1-\mu }^2}{{\mathrm{kR}}^2}{\left(\frac{\mathrm{n\pi }}{L}\right)}^6+\frac{{1-\mu }^2}{{\mathrm{kR}}^2}{k}_0{\left(\frac{\mathrm{n\pi }}{L}\right)}^4]}$

对于公式1-9，设：

$f_n=e^{{\lambda }_{m}y}$  (公式1-14)

将上式代入公式1-9得到特征方程：

求得特征根为：

${\lambda }_m=\pm \sqrt{{\left(\frac{\mathrm{n\pi }}{L}\right)}^2+\frac{{2k}_0}{1-\mu }},(m=\mathrm{9,10})$

于是：$f_n=A_9{e}^{{\lambda }_{9}y}+A_9{e}^{{\lambda }_{10}y}$

由F_n、f_n可以得出F、f(在实际计算中n不可能取无穷多项， 在n取5时已经达到较高精度)，又由公式：

$F_1=({▿}^2-k_0)F$  (公式1-15)

$F_2=k_0\frac{1-\mu }{2}{▿}^{4}F$  (公式1-16)

计算出F₁、F₂。同时，

${ϵ}_{11}=\frac{{\partial }_{u1}}{\partial x};{ϵ}_{22}=\frac{{\partial }_{u2}}{\partial y}+\frac{\omega }{R};\omega =\frac{{\partial }_{u1}}{\partial y}+\frac{{\partial }_{u2}}{\partial x};{\kappa }_1=\frac{{\partial \phi }_1}{\partial x};{\kappa }_2=\frac{{\partial \phi }_2}{\partial y}$

${\tau }_{12}+{\tau }_{21}=\frac{{\partial \phi }_2}{\partial x}+\frac{{\partial \phi }_1}{\partial y};{ϵ}_{13}=\frac{\partial \omega }{\partial x}+{\phi }_1;{ϵ}_{23}=\frac{\partial \omega }{\partial y}+{\phi }_2;M_{11}=D({\kappa }_1+u{\kappa }_2);M_{22}=D({\kappa }_2+u{\kappa }_1)$$D=\frac{{\mathrm{Eh}}^3}{12(1-u^2)};M_{21}=M_{12}=D\frac{1-u}{2}({\tau }_{12}+{\tau }_{12}).$薄壁件中面位移u₁、u₂、ω、φ₁、 φ₂：

$u_1=-\frac{1-\mu }{2a}(\mu \frac{{\partial }^3}{{\partial x}^3}-\frac{{\partial }^3}{\partial x\partial y^2})F_1$  (公式1-17)

$u_2=-\frac{1-\mu }{2a}((2+\mu )\frac{{\partial }^3}{{\partial x}^2\partial y}-\frac{{\partial }^3}{{\partial y}^3})F_1$  (公式1-18)

$\omega =(-\frac{1-\mu }{2a}{▿}^2{▿}^2))F_1$  (公式1-19)

${\phi }_1=\frac{\partial F_2}{\partial x}-\frac{\partial f}{\partial y}$  (公式1-20)

${\phi }_2=\frac{\partial F_2}{\partial y}+\frac{\partial f}{\partial x}$  (公式1-21)

将已经计算出的F₁、F₂、f分别代入公式1-17到公式1-21，从 而计算出各个薄壁件区域在多姿态下因自重引起的变形，式中： ${▿}^2=\frac{{\partial }^2}{{\partial x}^2}+\frac{{\partial }^2}{{\partial y}^2}.$

图7示出了根据本发明的飞机薄壁件的变形量建模方法的流程框 图，包括：简化飞机薄壁件和托板之间的定位夹紧关系为连续梁支承 模型；计算各托板的支持力；按照托板位置划分飞机薄壁件为多个相 应的壁板区域；简化各壁板区域的受力变形模型为托板位置薄壁件简 支支承下的变形；计算各壁板区域在多姿态下因自重引起的变形。该 建模方法结合上述公式来实现。

实施例：

飞机薄壁件的材料为铝合金(2024-T3)，弹性模量E＝70GPa；泊松 比μ＝0.33，密度ρ＝2.6g/cm³，壁板厚h＝2mm。选择机身壁板，壁板内 侧半径R＝3000mm，θ＝30°；壁板总长L＝5128mm；壁板总重量m＝42kg； 壁板依靠7个托板支承，其中，壁板前后两端各一个托板，其余五个 均匀分布，即托板间距L＝855mm，将整体壁板分成6个区域，各个 小区域薄壁件的重力为68.6N，转化为均布载荷q_n＝51N/mm²。壁板自 重在壁板法向的分布用均布力q_n^cosαcosβ代替(q_n^cosαcosβ中的α，β角度 随着自动钻铆过程的不同而改变，为已知条件)，求解壁板变形量。 计算托板受力及壁板变形的步骤如下：

本实例中机身壁板托板间距相同，壁板厚度也相同，在简化为连 续梁时相邻托板间载荷弯矩图的面积也相同，而且都为 ω_i＝ω_i+1＝L²q_w/12，同时a_i＝b_i＝L/2，L₁＝L₂＝…＝L₆＝855mm，I₁＝I₂＝…I₆； 代入弯矩方程：

$\frac{M_{i-1}{l}_i}{I_i}+2M_i(\frac{l_i}{I_i}+\frac{l_{i+1}}{I_{i+1}})+\frac{M_{i+1}{l}_{i+1}}{I_{i+1}}=-6(\frac{{\omega }_i{a}_i}{I_i{l}_i}+\frac{{\omega }_{i+1}{b}_{i+1}}{I_{i+1}{l}_{i+1}})$

得到：$M_{i-1}+4M_i+M_{i+1}=-\frac{q_w{L}^2}{2}$

依照上面的公式得到各个托板处的弯矩方程：

托板0：M₀＝0；

托板1：M₀+4M₁+M₂＝-q_wL²/2；

托板2：M₁+4M₂+M₃＝-q_wL²/2；

托板3：M₂+4M₃+M₄＝-q_wL²/2；

托板4：M₃+4M₄+M₅＝-q_wL²/2；

托板5：M₄+4M₅+M₆＝-q_wL²/2；

托板6：M₅+4M₆+M₇＝-q_wL²/2；

托板7：M₇＝0；

求解得到：

M₇＝0；

M₆＝-0.42255q_wL²；

M₅＝-0.30985q_wL²；

M₄＝-0.33805q_wL²；

M₃＝-0.33805q_wL²；

M₂＝-0.30985q_wL²；

M₁＝-0.42255q_wL²；

M₀＝0。

并依据公式1-3、公式1-4、公式1-5求解托板的支持力，计算结 果为：

F₀＝0.6056q_wL；

F₁＝0.866q_wL；

F₂＝1.035q_wL；

F₃＝0.9925q_wL；

F₄＝0.9925q_wL；

F₅＝1.035q_wL；

F₆＝0.866q_wL；

F₇＝0.6056q_wL。

表1至表4列出自动钻铆过程中四种典型位姿下薄壁件相关点的 变形量。

表1：α＝0°，β＝0°时各区域薄壁件沿法向的最大变形量及相应位 置坐标

表2：α＝10°，β＝2°时各区域薄壁件沿法向的最大变形量及相应位置 坐标

表3：α＝35°，β＝0°时各区域薄壁件沿法向的最大变形量及相应位置 坐标

表4：α＝35°，β＝4°时各区域薄壁件沿法向的最大变形量及相应位置 坐标：

本领域技术人员应能理解，上述实施例均是示例性而非限制性 的。在不同实施例中出现的不同技术特征可以进行组合，以取得有益 效果。本领域技术人员在研究附图、说明书及权利要求书的基础上， 应能理解并实现所揭示的实施例的其他变化的实施例。