2. 西北工业大学航天学院, 西安, 710072
2. School of Astronautics, Northwestern Polytechnical University, Xi′an, 710072, China
为了减重和提升飞机性能, 近年来, 在先进飞机设计和制造过程中采用了大量具有比强度高、比刚度大、抗疲劳和可设计性等诸多优异性能的复合材料。复合材料加筋壁板作为飞机复合材料结构的典型形式之一, 被广泛应用于飞机机翼、尾翼、机身和平尾等部位结构上, 如Boeing737平尾壁板、A320平尾及垂尾壁板等[1]。由于复合材料加筋壁板在压缩和剪切等载荷作用下易出现局部或整体失稳现象, 为了保证飞机使用安全, 在设计时必须对结构进行稳定性分析, 准确预测结构屈曲载荷。
国内外针对复合材料加筋壁板的稳定性分析进行了大量研究。穆朋刚等[2]分别采用特征值法和弧长法对“J”型复合材料加筋壁板稳定性进行了分析, 特征值法可确定加筋壁板的失稳临界载荷、位置及波形, 弧长法可判断加筋壁板的失稳临界载荷和位置。霍世慧等[3]利用工程及有限元方法分别分析了加筋壁板整体和局部稳定性问题, 利用有限元软件进行计算, 得到一种较合适的模型和边界条件。何吕龙等[4]通过剪切试验得到了不同筋条间距和腹板高度加筋壁板的屈曲载荷和破坏载荷并进行了分析。李乐坤等[5]分别采用特征值分析法和弧长法对加筋壁板进行了压缩屈曲分析, 并进行了试验验证。张国凡等[6]基于三维Hashin准则和平方应力准则建立了考虑渐进损伤的剪切破坏分析模型, 分析结果与试验结果符合良好。Zimmermann等[7]对复合材料加筋结构的屈曲分析准则进行研究, 并对计算结果与试验数据进行系统的总结分析, 探索一种快速、可靠的基于稳定性准则的设计方法; Kong等[8]通过数值和试验研究分析了复合材料加筋壁板的后屈曲特性, 讨论了铺层顺序和筋条形状对后屈曲承载能力的影响。Benedikt等[9]开发了一种估算复合材料加筋板轴压稳定性能的半解析方法。Riccio等[10]对复合材料加筋壁板在轴压载荷作用下的蒙皮-筋条脱粘过程进行了数值分析模拟, 计算结果与试验数据符合良好。Margarita等[11]采用有限元分析方法对大型复合材料加筋壁板的渐进破坏过程进行了模拟, 并与试验结果进行对比。陈金睿等[12]考虑筋条的扭转弹性支持作用, 采用里兹能量法建立了轴压复合材料加筋壁板蒙皮局部屈曲问题的理论模型。实验结果表明, 该方法可显著提高蒙皮局部屈曲载荷计算结果的精度。
在对复合材料加筋壁板的稳定性进行有限元分析时, 国内外学者几乎都采用常规有限元方法。常规有限元方法的计算单元位移形函数插值多项式最高阶固定在很低的值, 一般p=1或2, 求解精度通过改变网格粗细程度来控制, 计算收敛性通过减小最大单元尺寸(一般用h表示)实现, 因此被称为h型有限元。近年来才逐渐广泛应用的p型有限元技术是一种分析效率较高、误差可控的有限元分析技术。p型有限元在保持网格密度不变的情况下, 通过提高单元阶次p来提高计算精度。相较于h型有限元, p型有限元允许采用粗网格, 且单元长/细比最大可达200:1, 特别适合于复合材料层合结构这种铺层和胶层的厚度远小于面内尺寸的薄层结构力学性能分析。在复合材料层合结构有限元分析过程中, 为了得到更加精确的预测结果, 并且获取结构细节的变化过程, 常常需要对结构进行精细化建模, 即对结构的每个铺层、层间胶层、倒角区和填充区等进行建模。h型有限元不允许单元的长/细比过大, 在建立精细化模型时, 由于铺层和胶层的厚度尺寸限制, 单元尺寸与铺层和胶层的厚度应属于同一个数量级, 导致模型的单元数量巨大, 分析效率大大降低。p型有限元单元剖分时仅需要根据结构几何特征进行, 模型网格可以最大限度稀疏, 模型分析精度可通过提升单元阶次来保证, 而无需改变单元剖分形式及单元个数, 因此极大地降低了复合材料层合结构分析模型的单元数量和模型规模。
目前, 国内民用飞机复合材料用量较小、可用于借鉴的数据和资料较少, 新型飞机的研制和设计过程过多依赖试验验证技术。本文采用p型有限元软件建立了复合材料平尾加筋壁板的剪切稳定性分析模型, 结合初步技术方案选取四种典型复合材料平尾加筋壁板截面积, 研究不同筋条与蒙皮剖面面积比例对复合材料平尾加筋壁板剪切稳定性的影响。通过试验数据验证有限元分析方法的可行性和有限元模型的准确性。最后, 结合设计载荷情况, 为设计部门确定合适的复合材料平尾加筋壁板截面的结构型式和参数, 进行结构优化设计提供必要的建议。
1 试验件复合材料平尾加筋壁板剪切试验件由蒙皮和“T”型筋条通过胶结固化成型。蒙皮和筋条均采用M21/34%/UD194/IMA-12K碳带制造, 单铺层厚度为0.184 mm, 单向带力学性能参数如表 1所示。试验件共4组, 每组4件。第1组和第2组试验件蒙皮尺寸相同, 为等厚度蒙皮试验件; 第3组和第4组试验件蒙皮尺寸相同, 为变厚度蒙皮试验件。试验件尺寸参见图 1和表 2。试验件铺层参数见表 3。筋条与蒙皮剖面面积比例如表 4所示。为了便于载荷传递, 试验件筋条的缘条延伸直到试验件两端, 且在蒙皮加载区域上粘贴有铝加强片, 铝加强片厚度与缘条厚度相同。试验件筋条的腹板端部设计有45°倒角用于消除端部腹板缘条连接处的应力集中。根据复合材料平尾加筋壁板初步设计, 该复合材料平尾加筋壁板的剪切设计载荷为43.58 kN, 限制载荷(67%的设计载荷)为29.20 kN。
![]() |
表 1 单向带力学性能参数 Table 1 Mechanical properties of the prepreg |
![]() |
图 1 复合材料平尾加筋壁板试验件尺寸参数 Figure 1 Configuration of stiffened composite tail panels |
![]() |
表 2 复合材料平尾加筋壁板试验件几何尺寸 Table 2 Dimensions of stiffened composite tail panels |
![]() |
表 3 复合材料平尾加筋壁板铺层参数 Table 3 Stacking sequence of skin and stringers |
![]() |
表 4 复合材料平尾加筋壁板试验件筋条与蒙皮剖面面积比例 Table 4 Ratio of cross section area between stringers and skin |
2 有限元分析
本文采用p型商用有限元软件StressCheck[13]进行有限元分析。P型有限元技术是一种分析效率较高、误差可控的有限元分析技术。常规有限元方法的计算单元位移形函数插值多项式最高阶固定在很低的值, 一般p=1或2, 求解精度通过改变网格粗细程度来控制, 计算收敛性通过减小最大单元尺寸(一般用h表示)实现, 因此被称为h型有限元。P型有限元在保持网格密度不变的情况下, 通过提高单元阶次p来提高计算精度。相较于h型有限元, p型有限元允许采用粗网格, 且单元长/细比最大可达200:1, 特别适合于复合材料层合结构这种铺层和胶层的厚度远小于面内尺寸的薄层结构力学性能分析。
根据表 2所示的复合材料平尾加筋壁板尺寸参数和表 3所示的铺层参数, 分别建立与试验件对应的4类参数化模型。模型不仅包含了铺层和胶层细节, 而且包含了实际结构的众多关键细节, 如:壁板筋条根部三角填充区细节、下缘条与腹板之间的R区细节、蒙皮厚度过渡段细节等。在网格剖分时, 根据结构特点, 在下缘条与腹板之间的R区和腹板端部45°倒角区等关键部位采用加密的网格, 在其他区域采用稀疏的网格。两类等厚度蒙皮的复合材料平尾加筋壁板模型单元数目相同, 均包含8 287个单元。两类变厚度蒙皮的复合材料平尾加筋壁板模型单元数目相同, 均包含10 512个单元, 如图 2所示。所有单元均采用实体单元。
![]() |
图 2 复合材料平尾加筋壁板有限元模型及边界条件 Figure 2 Finite element models and the constraints |
为了模拟试验过程中试验件的约束条件, 在模型4个宽50 mm的侧边平面上施加对称约束, 限制模型的4个端面沿厚度方向(Y方向)的位移和面外转动。通过在4个侧边端面作用均匀分布的剪流来实现剪切载荷的施加。由于施加在4个侧边端面的剪流相等, 所以壁板两个长边的剪切载荷比短边剪切载荷大。在判断结构是否满足剪切稳定性要求时, 应取短边的剪切载荷F2与限制载荷进行对比。复合材料平尾加筋壁板剪切载荷示意图如图 3所示, F1和F2为剪切载荷, Ftotal为总载荷, 即试验机施加的载荷。在进行有限元计算时, 假设施加在短边上的集中载荷值为43.58 kN, 即复合材料平尾加筋壁板的设计载荷。根据短边的剪切载荷计算模型施加的剪流。通过计算可以得到模型的一阶模态及其对应的屈曲系数。模型的屈曲载荷通过式(1, 2) 计算得到。
![]() |
图 3 复合材料平尾加筋壁板剪切载荷示意图 Figure 3 Shear load of the stiffened panels |
$\begin{matrix} {{F}_{总屈曲}}= \\ {{F}_{设计}}\times (\frac{\sqrt{({{L}_{长边}}-100){{~}^{2}}+({{L}_{短边}}-100){{~}^{2}}}}{({{L}_{短边}}-100)}\text{ }) \\ \end{matrix}$ | (1) |
$F_{\rm 短边屈曲}=F_{\rm 设计}×屈曲系数$ | (2) |
式中:F总屈曲为复合材料平尾加筋壁板模型的屈曲载荷, F短边屈曲为屈曲时短边的剪切载荷, F设计为复合材料平尾加筋壁板的设计载荷, L长边为复合材料平尾加筋壁板长边的长度, L短边为复合材料平尾加筋壁板短边的长度。剪流施加在宽度为50 mm固定的区域上, 在计算剪切载荷时需将长边和短边的长度中减去固定区域的宽度。
P型有限元计算时单元的阶次一般可设在p=1~8之间, 随着单元阶次的增加, 模型计算量增大, 计算精度提高。对比试验数据发现, 当p=5时, 计算误差仅为0.13%。因此, 本文在进行模型分析时设定单元阶次为5。
计算得到复合材料平尾加筋壁板剪切一阶屈曲模态(如图 4所示)和屈曲系数。由计算结果可知, p型有限元模型预测的4组复合材料平尾加筋壁板的剪切屈曲模态基本一致。屈曲时, 壁板纵向均有2个半波(即中间2个加强筋之间蒙皮屈曲半波个数)。可见, 不同筋条和蒙皮剖面面积比例对复合材料平尾加筋壁板的剪切屈曲模态几乎没有影响。采用式(1, 2) 计算得到的复合材料平尾加筋壁板剪切屈曲载荷如表 5所示。
![]() |
图 4 复合材料平尾加筋壁板剪切屈曲模态 Figure 4 Predicted shear buckling modes of the stiffened panels |
![]() |
表 5 有限元预测的复合材料平尾加筋壁板屈曲载荷 Table 5 Predicted buckling load of the stiffened panels |
3 试验测试
复合材料平尾加筋壁板试验件蒙皮正面(有筋条的一侧)和反面上的应变花粘贴位置结合屈曲模态预测结果中波峰和波谷的位置确定, 如图 5所示。试验夹具和试验件安装如图 6所示。由于筋条的存在和变厚度蒙皮试验件蒙皮厚度不一致, 导致试验件剪心与蒙皮形心不在同一平面内。在试验件安装前, 计算4类试验件的剪心位置, 试验件安装过程中, 通过改变夹头豁口处的0.01 mm薄垫片数量来调整试验件的位置, 使试验件的剪心与试验加载轴线重合, 如图 7所示。剪切试验在电液伺服万能试验机上进行, 采用位移控制进行分级加载, 加载速率为1 mm/min。试验分两步进行:第(1) 步:预加载。按5%Pjsqq(计算屈曲载荷)级差逐级加载至40%Pjsqq, 然后按10%Pjsqq级差逐级卸载。对所有加载和测量设备进行检查, 并对应变花测量结果进行计算, 如正面和背面对应位置的应变花值偏差小于10%, 进行试验第(2) 步, 否则, 检查试验件安装并进行调整。第(2) 步:破坏加载。按5%Pjsqq级差逐级加载至65%Pjsqq, 按2%Pjsqq级差加载至67%Pjsqq, 按3%Pjsqq级差加载至70%Pjsqq, 再按2%Pjsqq级差加载至100%Pjsqq, 保载3 s, 再按2%Pjsqq级差加载至破坏。加载期间逐级对载荷和应变进行测量。试验过程中, 对试验件进行录像, 观察试验件的变形和屈曲情况。试验环境为室温干态。试验结束后, 对试验件破坏模式进行拍照记录。
![]() |
图 5 复合材料平尾加筋壁板剪切试验件应变花粘贴位置示意图 Figure 5 Location of the strain gauges |
![]() |
图 6 试验件安装示意图 Figure 6 Loading state of the specimen |
![]() |
图 7 试验件剪心位置调整示意图 Figure 7 Adjustment of location of the shear center |
预加载过程中, 夹具和试验件均未出现异常情况, 试验件也没有出现响声。破坏加载过程中, 当载荷接近屈曲载荷时, 试验件发出连续的噼啪声, 并可观测到明显的蒙皮鼓包现象, 即蒙皮出现屈曲。之后, 蒙皮与筋条脱粘, 试验件突然完全破坏, 失去承载能力。卸载后, 试验件蒙皮鼓包消失, 从试验件背面观察发现蒙皮无明显破坏。
复合材料平尾加筋壁板试验件剪切试验典型破坏模式如图 8所示。由图可知, 复合材料平尾加筋壁板在剪切载荷作用下蒙皮上可观测到鼓包现象, 结构典型破坏模式为蒙皮剪切失稳所致的蒙皮与筋条脱粘。剪应变可根据直角应变花沿0°, 45°, 90°三个方向的应变测量值ε0, ε45, ε90计算, 表达式为
![]() |
图 8 试验件典型破坏模式 Figure 8 Typical fracture modes of the specimens |
${{\gamma }_{xy}}={{\varepsilon }_{0}}+{{\varepsilon }_{90}}-2{{\varepsilon }_{45}}$ | (3) |
复合材料平尾加筋壁板剪切试验件加载轴线上蒙皮正面和背面对应位置剪应变-载荷曲线如图 9所示。由图可知, 在试验件弹性变形阶段, 剪应变-载荷曲线线性变化, 当载荷相同时, 所有位置的剪应变值几乎一致。试验件蒙皮正面和背面对应位置的剪应变值一致说明试验件在加载过程中没有出现弯曲变形; 试验件蒙皮加载轴线上不同位置的剪应变值一致说明试验机的集中载荷通过夹具均匀的施加在试验件上, 试验件远离边界区域承受纯剪切载荷。可见, 在加载过程中, 试验件剪心与加载轴线重合, 试验夹具满足要求, 加载路径合理。当试验件蒙皮发生屈曲时, 试验件蒙皮出现了鼓包, 导致应变值发生突变, 曲线出现明显的分叉, 曲线分叉拐点对应的载荷即为试验件的剪切屈曲载荷。复合材料平尾加筋壁板剪切屈曲试验结果见表 6。由试验结果可知, 试验件的屈曲载荷和破坏载荷测试结果离散系数最大为6%, 试验一致性非常好。
![]() |
表 6 复合材料平尾加筋壁板剪切试验结果 Table 6 Test results of the fracture load and the buckling load of the specimen |
![]() |
图 9 复合材料平尾加筋壁板剪切试验件典型剪应变-载荷曲线 Figure 9 Typical experimental local strains vs. load curves of the specimens |
分析图 9可知, 等厚度蒙皮试验件的剪应变-载荷曲线一直保持较好的线性度, 直到试验件发生屈曲, 且试验件破坏与屈曲几乎同时发生, 屈曲载荷与破坏载荷非常接近; 变厚度蒙皮试验件的剪应变-载荷曲线出现了整体分叉现象, 对应于试验件发生了整体剪切屈曲, 试验件在发生剪切屈曲后, 仍能够继续承受更大的载荷, 其破坏载荷大约比屈曲载荷高30%~40%。
4 分析与讨论 4.1 有限元分析结果与试验结果对比复合材料平尾加筋壁板剪切屈曲载荷有限元分析结果与试验测试结果如表 7所示。由表可知, 有限元分析结果与试验测试结果一致, 偏差小于18%。可见, 在计算过程中采用的建模方法和分析方法得当, 所建分析模型有效。有限元分析结果小于试验测试结果的主要原因为:有限元和试验施加的载荷不完全一致, 相对于有限元计算分析过程, 试验过程中夹具对试验件的约束过强, 不是理想简支约束, 导致试验件较晚出现失稳。
![]() |
表 7 复合材料平尾加筋壁板剪切屈曲载荷试验件有限元分析结果与试验结果比较 Table 7 Comparison of the buckling load between the test results and the numerical results |
在组号1模型的计算结果中提取载荷为240 kN时1105号应变花位置处的剪应变值(6 929微应变), 绘制剪应变-载荷曲线与试验结果进行比较, 如图 10所示。由图可知, 试验测得的剪应变-载荷曲线和有限元预测的结果吻合较好, 当载荷为240 kN时, 试验测得和有限元预测的剪应变分别为7 232微应变和6 929微应变, 相对误差为-4.1%。提取组号1等厚度蒙皮和组号3变厚度蒙皮模型在屈曲前相同载荷(66 kN)下的应力云图如图 11所示。由图可知, 在屈曲前蒙皮上的应力分布一致, 与图 9显示的结果一致。但是由于变厚度蒙皮模型加强筋之间的蒙皮变薄, 在相同的载荷作用下, 变厚度蒙皮模型的剪应变比等厚度蒙皮的模型的剪应变高约44%。
![]() |
图 10 有限元预测和试验得到的剪应变-载荷曲线对比 Figure 10 Comparison of shear strain vs. load curves between test results and numerical results |
![]() |
图 11 有限元预测应力分布云图 Figure 11 Predicted stress nephogra |
4.2 筋条和蒙皮剖面面积对试验件性能的影响
将复合材料平尾加筋壁板筋条和蒙皮的剖面面积与试验屈曲载荷、破坏载荷进行列表比较, 如表 8和表 9所示。由表 8可知, 在蒙皮剖面面积相同情况下, 增加筋条剖面面积不能有效提升试验件的屈曲载荷和破坏载荷, 筋条重量增加约15%, 试验件性能提升只有约5%。由表 9可以看出, 在筋条剖面面积相近的情况下, 虽然蒙皮截面积只降低了17%, 但是屈曲载荷下降了约39%。可见, 在筋条剖面面积/蒙皮剖面面积近似相同情况下, 增加蒙皮截面积和筋条剖面面积均可提高复合材料平尾加筋壁板的剪切屈曲性能和破坏性能, 但增加蒙皮截面积对剪切对屈曲性能和破坏性能影响程度会更大些。
![]() |
表 8 复合材料平尾加筋壁板试验件筋条和蒙皮剖面面积与试验件性能对比 Table 8 Effect of ratio between cross section of stringers and skin on properties of specimens |
![]() |
表 9 复合材料平尾加筋壁板试验件筋条和蒙皮剖面面积与试验件性能对比 Table 9 Effect of ratio between cross section of stringers and skin on properties of specimens |
5 结束语
本文采用p型有限元法预测了复合材料平尾加筋壁板剪切屈曲载荷, 然后通过试验获得复合材料平尾加筋壁板剪切屈曲载荷, 对有限元结果进行了验证。验证结果表明, p型有限元分析方法可以高效的完成复合材料平尾加筋壁板剪切稳定性分析。在剪切载荷作用下, 复合材料平尾加筋壁板典型破坏模式为蒙皮剪切失稳所致的蒙皮和筋条脱粘; 其中, 等厚度蒙皮试验件屈曲载荷与破坏载荷非常接近; 变厚度蒙皮试验件破坏载荷大约比屈曲载荷高30%~40%。在筋条剖面面积/蒙皮剖面面积近似相同情况下, 复合材料平尾加筋壁板的蒙皮截面积决定其承受剪切载荷的性能, 但对其屈曲模态几乎没有影响。该复合材料平尾加筋壁板的设计载荷仅为44 kN。可见, 综合考虑结构重量和承载能力, 第3组试验件的构型为较优构型。
[1] | 中国航空研究院. 复合材料结构设计手册[M]. 北京: 航空工业出版社, 2001: 93. |
[2] |
穆朋刚, 万小朋, 赵美英.
复合材料平尾加筋壁板稳定性分析研究[J]. 机械科学与技术, 2009, 28(9): 1190–1193.
MU Penggang, Wan Xiaopeng, ZHAO Meiying. A study of the stability of composite stiffened plates[J]. Mechanical Science and Technology for Aerospace Engineering, 2009, 28(9): 1190–1193. |
[3] |
霍世慧, 王富生, 王佩艳, 等.
复合材料机翼加筋壁板稳定性分析[J]. 应用力学学报, 2010, 27(2): 423–427.
HUO Shihui, WANG Fusheng, WANG Peiyan, et al. Stability analysis on the ribbed panel of the composite wing[J]. Chinese Journal of Applied Mechanics, 2010, 27(2): 423–427. |
[4] |
何吕龙, 尚柏林, 常飞, 等.
筋条参数对复合材料加筋壁板剪切承载能力的影响[J]. 机械工程材料, 2015, 39(9): 49–52.
DOI:10.11973/jxgccl201509011 HE Lvlong, SHANG Bolin, CHANG Fei, et al. Effects of stiffener parameters on shear carrying capaciy of composite stiffened wall panel[J]. Materials for Mechanical Engineering, 2015, 39(9): 49–52. DOI:10.11973/jxgccl201509011 |
[5] |
李乐坤, 李曙林, 常飞, 等.
复合材料加筋壁板压缩屈曲与后屈曲分析[J]. 南京航空航天大学学报, 2016, 48(4): 563–568.
LI Lekun, LI Shulin, CHANG Fei, et al. Buckling and post-buckling of composite stiffened panel under compresion[J]. Journal Nanjing University of Aeronautics & Astronautics, 2016, 48(4): 563–568. |
[6] |
张国凡, 孙侠生, 孙中雷.
复合材料加筋壁板剪切破坏试验与后屈曲分析[J]. 机械科学与技术, 2016, 35(8): 1280–1285.
SUN Guofan, SUN Xiasheng, SUN Zhonglei. Failure test and post-buckling analysis of composite stiffened panels under shear load[J]. Mechanical Science and Technology for Aerospace Engineering, 2016, 35(8): 1280–1285. |
[7] | ZIMMERMANN R, KLEIN H, KLING A. Buckling and post-buckling of stringer stiffened fiber composite curved panels-tests and computations[J]. Composite Structures, 2006, 73(2): 150–161. DOI:10.1016/j.compstruct.2005.11.050 |
[8] | KONG C W, LEE I C, KIM C G, et al. Postbuckling and failure of stiffened composite panels under axial compression[J]. Composite Structures, 1998, 42(1): 13–21. DOI:10.1016/S0263-8223(98)00044-0 |
[9] | BENEDIKT K, EELCO L J, RAIMUND R. Semi-analytic probabilistic analysis of axially compressed stiffenedcomposite panels[J]. Composite Structures, 2012, 94: 654–663. DOI:10.1016/j.compstruct.2011.08.033 |
[10] | RICCIO A, RAIMONDO A, SCARAMUZZINO F. A robust numerical approach for the simulation of skin-stringerdebonding growth in stiffened composite panels under compression[J]. Composites: Part B, 2015, 71: 131–142. DOI:10.1016/j.compositesb.2014.11.007 |
[11] | MARGARITA A, EELCO J, SINA H, et al. Efficient progressive failure analysis of multi-stringer stiffenedcomposite panels through a two-way loose coupling global-localapproach[J]. Composite Structures, 2017. |
[12] |
陈金睿, 陈普会, 孔斌, 等.
考虑筋条扭转弹性支持的轴压复合材料加筋板局部屈曲分析方法[J]. 南京航空航天大学学报, 2017, 41(1): 76–82.
CHEN Jinrui, CHEN Puhui, SUN Bin, et al. Local buckling analysis of axially compressed stiffened laminated panels considering rotational restraint of stiffeners[J]. Journal Nanjing University of Aeronautics & Astronautics, 2017, 41(1): 76–82. |
[13] | SZABO B, BABUSKA I. Stress check release 9.0 master guide[M]. [S.l.]: Engineering Software Research & Development Inc., 2009. |