板壳结构加筋布局的仿生脉序生长算法

薛开，李永欣

(哈尔滨工程大学 机电工程学院，黑龙江 哈尔滨 150001)

板壳结构由于轻量化的优势在航空航天、船舶海洋、汽车民用产品中被广泛采用，其加强筋的分布优化问题也一直受到国内外学者的重视.传统的加强筋未经过系统的优化设计，为保守的规则分布，随着优化设计方法尤其是以材料分布为目标的拓扑优化理论的逐步完善，加强筋的布局优化也得到了快速发展，变密度［1］、均匀化［2］和渐进结构［3］等拓扑优化技术都被应用到了加强筋的布局优化中.此外，一些学者尝试用智能算法来求解加筋布局问题，Kallassy等［4］用遗传算法对薄壳的加筋方式进行了进一步的探索，并用神经网络来描述优化目标.

加强筋的布局设计仍处在发展的阶段，用连续拓扑优化理论得到的是结构的加强区域，要转化为设计方案还需要设计师的处理.此外，离散拓扑优化由于可适用的算法少，应用受限，如果能够从自然界生物体获得规律，就可以为加强筋的离散拓扑优化研究提供一种新的思路.丁晓红等［5-6］基于根系等分枝系统采用自适应生长技术对加筋的布局进行了优化.作为生物体中特殊的分枝系统，植物叶脉与加筋薄板具有明显的生物相似性［7］，本文探索从叶脉脉序形态规律中提炼结构优化的思想，并将其用于板壳结构的加筋布局优化问题.

1 不规则加强筋的布局优化模型

工程中广泛应用的加强筋受工艺的限制，多为规则分布(等间距、正交或垂直排列)，目前数字化制造技术的快速发展使不规则布局的筋条成型成为可能.这里对不规则分布的加强筋进行建模，图1可以表示加强筋从规则到不规则的衍变过程，衍变过程中只改变相对的角度、高度和宽度，仍遵循如下原则:

1)筋条只分布于基板一侧;

2)筋条末端不应在内部停止，需与其他筋条或结构外壁连接;

3)筋条的最大宽度不应大于基板的厚度，最大高度不应大于3倍的基板厚度.

图1 加强筋的分布衍变Fig.1 Rib distribution evolution

本文基于有限单元技术对加强筋进行结构分析和优化，将所有的有限元节点作为可能生成加强筋的点集，每个节点只和其邻域内的节点建立联系，如图2所示.由此得到的加强筋间相对角度θ和筋条长度l都为离散集.基于有限元技术的加强筋分布优化，将不规则分布局限在一个有限离散集中，在一定程度上简化了加工工艺.

图2 板梁组合结构的有限元模型Fig.2 FE model of stiffened plate structure

加强筋布局优化的对象是加强筋的分布和筋条的尺寸参数，为拓扑优化和尺寸优化的综合.这里，取优化参数为节点的有效性γ(γ=1表示节点有效，γ=0表示无效)、筋条的宽度di和高度hi.将节点所组成的点集空间U作为布局优化的优化空间.优化的目标是得到一个顺序集合P(P⊂U)和尺寸的对应离散点集H、P中节点组成的几何拓扑形式即为加强筋的布局分布.薄板加强筋的布局优化问题可以描述为:在一定的力学性能要求和材料数量的限制下，找到满足要求的顺序集合P和尺寸参数di、hi.那么用最少的材料达到某许用位移的薄壳加强筋的布局优化模型为

2 双子叶植物叶脉形态规律

植物叶子可以看作分布有加强筋(叶脉)的薄板(叶片)结构，且在不同的环境中具有自适应性(不同生长环境下的叶脉布局略有不同)，是加筋薄板理想的生物模型.与单子叶植物的平行脉序或弧状脉序不同，双子叶植物叶脉脉序以等级和闭环为特点.植物叶脉为植株承担着物质传输的功能，从力学角度衡量，植物叶子也具有很明显的增强作用［8］.

2.1 双子叶植物叶脉的形态拓扑规律

叶脉的宏观构型规律:主脉是直径最粗且最早生长的，一级次脉在主脉两侧分布，二级次脉在一级次脉两侧分布;随着脉序级数的增加，直径逐渐减小;脉络趋向于形成闭环和网状结构;脉络之间以适当的角度连接，一级叶脉和主脉之间的角度多为30°～50°，高级脉之间的角度则接近90°，如图3(a)所示.由于叶脉是特殊的分枝循环系统，满足广义Murray法则［9］，更详细的描述为Bohn等［10］提出的叶脉分岔处的矢量平衡准则，如图3(b)所示.

图3 双子叶植物叶脉的拓扑形态规律Fig.3 Topology regularity of dicotyledonous venation

2.2 双子叶植物叶脉形态发生的结构力说

对于植物叶脉，尤其是具有等级和闭环的双子叶植物叶脉的形态形成学说，主要有激素说和结构力说.结构力说是由法国学者Y.Couder等［11］针对双子叶植物脉序提出的.认为不均匀生长造成的结构力在叶脉的形成初期，尤其是局部高次脉序的闭环形成过程中有重要作用.随后Laguna、Corson等［12-13］对结构力说进行了完善.结构力说突出了形态发生的力学因素，成功模拟了叶脉闭环的产生，弥补了生长素说只能形成树状结构的缺陷.其中，结构力说的主要观点有:

1)叶片分为表皮和叶肉，由于表皮的生长速度比叶肉慢，导致表皮受拉力而叶肉受压缩，形成了局部的应力场，叶脉就是在此应力场中形成的.双子叶植物叶脉的形态形成过程分为分化和重塑2个阶段.

2)在分化阶段，当叶肉细胞承受的压缩力超过了某临界值，就会分化成形成层—叶脉的雏形，且新形成的叶脉与主压缩力垂直，沿释放剪应力的方向生长，此时的形成层还是较凌乱的分布在叶肉里.

3)在重塑阶段，变形细胞垂直于形成层方向的粘性增大，产生了沿着叶脉方向的切向力，切向力使叶脉细胞微调以满足矢量平衡方程，切向力F和细胞变形能Γc可以表示为

式中:n是形成层方向脉细胞的个数，di和Sci分别为脉细胞i的宽度和面积，vpi和Ci分别为脉细胞i的粘度和切向系数，Δli为细胞壁的变形量差，Ti是脉细胞的内部肿胀力，Pi为外界均匀压力.

4)在叶脉形成过程中，由于生物体变形速度远远小于力信号的传播速度，因此在力学层面，局部满足超静定原则，也就是说在整个形态发生过程中保证局部弹性变形能为最小，即

由叶脉形态发生的结构力说可以看出，弹性变形能和剪应力在叶脉形态形成中为决定因素.而且剪应力的论述和Mattheck［14］关于树木纤维的分布结论是一致的，都沿着释放剪应力的方向延伸.因此，以变形能和剪应力为生长规则的形态力学规律可以指导仿生叶脉脉序的生长.

3 薄板加强件布局的脉序生长算法

3.1 VGA算法的概念定义和理论基础

脉序生长算法(venation growing algorithm，VGA)，是受自然界植物叶脉脉序结构力说及生物自适应生长观点的启发，基于有限元技术用于薄板/壳结构上加强脉布局优化的仿生算法.VGA为等级优化方法，生成的脉络由应变能决定的主脉和剪应力决定的次脉组成，并用矢量平衡准则进行加强筋几何尺寸的计算.在算法的实现过程中，采用设计灵敏度来衡量设计参数对性能参数的影响程度.

3.1.1 以应变能灵敏度为基础的主脉的生长

对于组合结构，结构的变形能由弹性变形能和耦合能组成，表示为

本文板梁组合结构的数学模型采用全离散板梁理论，用板壳单元和梁单元来代表基板和加强筋模型(如图2所示)，整体刚度矩阵由基板的刚度矩阵和考虑偏心的等效到板单元坐标系的梁单元刚度矩阵组成，即

式中:Kb'为梁单元的等效刚度矩阵，需要考虑偏心的影响.由于全离散模型，因此耦合能量为零.

由于［(δT)'kδ］= δTkδ'，且为常值，单元矩阵 ke为对称阵，因此应变能对体积灵敏度为

这里，刚度改变量等于离散板梁模型中将梁单元刚度等效到板单元坐标的刚度:

其中:

为具有偏心距e的板梁组合结构的转换矩阵，e=(δ+h)/2.

单元应变能对于体积的灵敏度集合:

式中:n为生长点集合中节点的个数，取灵敏度最小的点作为下一个主脉生长点.

3.1.2 以剪应力灵敏度为基础的次脉的生长

本文将一级次脉简称为次脉，次脉沿着剪应力的梯度方向ΔQ(x)生长.由于在薄板弯曲中，弯曲应力δx、δy和扭应力τxy为主要应力，横向剪应力和挤压应力为次要应力，因此此处的剪应力以扭应力τxy作为参考，取其对横坐标和纵坐标的灵敏度加权和为剪应力灵敏度:

取剪应力灵敏度最大的点作为下一个生长点.

3.1.3 优化加强筋高度

本文认为加强筋的宽度为一定，只优化其高度.并且由矢量平衡准则可以推得:加强筋不能在中间停止，分岔处需要满足矢量平衡方程.然后按下式优化截面尺寸:

3.2 VGA算法的实现流程和判定规则解释

基于有限元技术的VGA算法的实现，可以看作以初始有限元结构分析为输入变量而进行的多等级的寻优方法.三级寻优的实现过程的流程如图4所示.

图4 VGA的实现流程Fig.4 VGA realization flow

值得指出的是，在优化过程中，由于需要计算增加加强脉后的力学参数，因此每次迭代计算都要进行有限元分析.

4 基于VGA的典型薄板结构的加强筋分布优化设计

本文采用文献［6］关于自适应生长种子选取的结论，种子的位置选在载荷作用处和支撑处，次脉起点则选在主脉集合的剪应力极值处.

4.1 单角简支薄板加筋分布设计

单角固支薄板模型如图1所示.其初始几何参数:l=0.3m，δ=0.000 3m，F=1N，材料的弹性模量E=2×1011MPa，泊松比为 0.3，加筋的初始高度为0.000 5m，种子位于支撑点.得到加强脉的分布形态如图6所示.可见单角固支、对角加载的薄板模型的加强脉分布结果和植物叶子的分布结果相似.主脉从种子到对角延伸，而次脉从主脉分支到边界.

图5 单角固支薄板模型Fig.5 Single-point fixed plate model

图6 加强脉分布结果Fig.6 Rib layout result

4.2 四边简支薄板

四边简支薄板的几何参数和材料参数和上例相同，改变约束形式和载荷作用位置.将板的四边进行简支，载荷作用于板中心，如图7所示.种子的位置位于薄板中心，得到不同迭代过程中加强脉的生长形态结果如图8所示，最终加强筋使用的材料体积为3.62×10-8m3.将图8中的最终结果只取主脉作为简化结果，取相同体积材料，与传统型式、文献［3］，文献［5］进行对比，结果如图9所示，横坐标D为书点到板中心的距离.

图7 四边简支薄板模型和分析轨迹Fig.7 Simply supported square plate model

图8 不同迭代过程中加强脉的生长形态结果Fig.8 Rib growth result in iteration process

由图8可以看出，随着生长的进行，加强脉从板中心先以密集回旋，然后逐步生长，延伸至四角，次脉的生长则几乎与主脉垂直，连接四边的中点.

由图9可以看出，沿着2种轨迹得到的结果是基本一致的，3种加筋型式都对基板有不同程度的增强，文献［3］的加筋型式，比传统型式略好，但不如文献［5］和本文优化结果，而本文采用的VGA算法的改善效果是最显著的.

图9 优化结果对比Fig.9 Optimization results comparison

5 结论

通过对植物脉序规律的阐述和基于VGA的算例设计，可以得出如下规律:

1)双子叶植物叶脉，整体具有等级和闭环的特点，局部分岔满足矢量平衡条件，在形态发生过程中，变形能和剪应力起着重要作用;

2)以有限元理论和离散拓扑优化为基础，结合双子叶植物叶脉形态规律的VGA算法，适合薄板加强筋的分布设计;

3)VGA以单元应变能灵敏度和剪应变为准则，进行加强脉的等级优化，得到的加强筋主次分明，适合加工，且增强效果显著;

4)本算法为自适应生长算法，可以深入探索在结构动力优化方面的用途，并可以探索在复合材料加强纤维的铺设角度方面的应用.

［1］LAM Y C，SANTHIKUMAR S.Automated rib location and optimization for plate structures［J］.Structural Multidisciplinary Optimization，2003(25):35-45.

［2］CHUNG J，LEE K.Optimal design of rib structures using the topology optimization technique［J］.Proceedings of the Institution of Mechanical Engineers，1997(211)，Part C:425-437.

［3］章胜冬，杨军刚，张卫红.薄板结构加筋布局设计的渐进结构优化方法［J］.现代制造工程，2009(4):5-9.

ZHANG Shengdong，YANG Jungang，ZHANG Weihong.Evolutionary structural optimization for stiffener layout problem of plate structures［J］.Modern Manufacturing Engineering，2009(4):5-9.

［4］KALLASSY A，MARCELIN J L.Optimization of stiffened plates by genetic search［J］.Structural Optimization，1997(13):134-141.

［5］DING Xiaohong，YAMAZAKI K.Adaptive growth technique of stiffener layout pattern for plate and shell structures to achieve minimum compliance［J］.Engineering Optimization，2005，37(3):259-276.

［6］丁晓红，林建中，山崎光悦.利用植物根系形态形成机理的加筋薄壳结构拓扑优化设计［J］.机械工程学报，2008，44(4):201-205.

DING Xiaohong，LIN Jianzhong，YAMAZAKI K.Topology design optimization of stiffened thin-wall shell structures based on growth mechanism of root system［J］.Chinese Journal of Mechanical Engineering，2008，4(44):201-205.

［7］LI Yongxin，XUE Kai.Mechanics in leaf venation morphogenesis and their biomimetic inspiration to construct a 2-Dimensional reinforcement layout model［J］.Journal of Biomimetics， Biomaterials， and Tissue Engineering， 2011(10):81-93.

［8］KOBAYASHI H，DAIMARUYA M，KURIBAYASHI K K.Venation pattern of butterbur leaf and its mechanical contribution［J］.Journal of the Society of Materials Science，2000，49(12):1318-1323.

［9］DURAND M.Architecture of optimal transport networks［J］.Physical Review E，2006，73(016116):1-6.

［10］BOHN S，ANDREOTTI B，DOUADY S，MUNZINGER J，COUDER Y.Constitutive property of the local organization of leaf venation networks［J］.Physical Review，2002，65(061914):1-12.

［11］COUDER Y，PAUCHARD L，ALLAIN C，ADDA-BEDIAM，DOUADY S.The leaf venation as formed in a tensorial field［J］.The European Physical Journal B，2002，28:135-138.

［12］LAGUNA M F，BOHN S，JAGLA E A.The role of elastic stresses on leaf venation morphogenesis［J］.PLoS Computational Biology，2008，4(4):e1000055.

［13］CORSON F，ADDA-BEDIA M，BOUDAOUD A.In silico leaf venation networks:growth and reorganization driven by mechanical forces［J］.Journal of Theoretical Biology，2009(259):440-448.

［14］MATTHECK C.Designing in nature:learning from trees［M］.Berlin:Springer-Verlag，1998:129～140.