膜结构全过程计算方法
作者 向阳
摘要： 针对膜结构计算设计的全过程，进行了理论方法的研究，并编制了相应的计算设计程序。工程实践表明了本文所述的理论方法及计算程序的正确性及广泛适用性。 关键词： 膜结构、计算方法
一、引 言 膜结构是大跨空间结构的主要形式之一，新颖独特的建筑造型、优良的物理力学性能、简洁高效的结构体系，使其必将成为我国二十一世纪空间结构发展的主要方向。膜结构作为一种柔性张力体系，与传统的刚性结构在设计计算上有诸多不同之处，其设计计算的一般过程是：初始形态设计～静力性能分析～风振响应分析～裁剪设计计算，其中每一环节都是不可或缺的。 经验告诉我们，一种新兴结构的推广和发展，与其计算理论方法、相应的设计软件的成熟有相当大的关系。目前国内已经出现了诸如上海八万人体育场等少数应用膜结构的建筑，但遗憾的是其膜部分的设计、计算皆是国外公司所为，我国自行设计的大型膜结构还很少。这当然是国外控制其关键技术，也是国内尚无成熟的计算理论方法及相应的计算设计程序的缘由。正因为如此，本文作者一方面进行膜结构的计算理论方法的研究，一方面编制了相应的计算设计程序，旨在能为膜结构在我国的应用发展尽微薄之力。 二、计算理论方法 （一）非线性有限元基本方程 有限单元法是对工程结构进行数值分析的最有效方法，特别是在计算机应用越来越普及的今天。膜结构实际设计中，索及桁架等加强、边缘构件的应用是必不可少的，因此本文程序中包含了膜单元，索单元，杆单元，梁单元。其中空间膜单元定义为三结点的三角形等参元，考虑节点的ｘｙｚ三个方向的位移，但只计及面内的正应力σx、σy和剪应力Z xy。应用Ｕ•Ｌ法列式，可以得结构有限元基本迭代方程为： （1）其中{R}为外荷载向量；{F}为ｔ时刻单元应力节点等效力向量；[kL]为线性应变增量刚度矩阵；[KNL]为非线性应变增量刚度矩阵，非线性方程组的求解，采用增量形式的 Full Newton-Raphson 法。具体推导见文献[1]。
（二）初始形态设计 本文这里使用了“形态”这个词汇，所谓的“形”就是几何意义上的形状，所谓的“态”就是结构的内力分布状态。一种“形”对应一种“态”；反之亦然，一种“态”必然有一种“形”与它对应。膜结构这类柔性结构有一个显著的特点，即结构材料本身不具有刚度，由这些材料组成的构造体系还只是机构，只有当对其施加了预张力，它才具有了抵抗外荷载的结构刚度。这里应强调的是，不同的张力分布对应着不同的平衡形状。这就是本文所说的“形”和“态”。当然，这是一个动态跟踪过程，对这个过程的研究就是—初始形态设计。这是膜结构与传统刚性结构计算的一个显著区别，也膜结构设计计算中的一个关键问题。 膜结构的设计中，在找到初始形态之前，并不能准确确定膜结构的初始形状和与之对应的预张力分布状态，也就是说这时有两个未知数：一个是初始形状，一个是预张力分布状态。这时我们会给定其中的一个而求解另一个，从而产生两种思路。第一种思路产生比较自然，是将初始形状作为已知数，把初始预张力作为外荷载施加到结构上，求解达到平衡时的状态。此法求解比较直接方便，普通的静力计算程序就可以适用，并可以得到与设计者给定的曲面形状相近的结果。但是它最大的缺点是最终得到的预张力分布不再是设计给定的，当膜结构曲面比较复杂不规则时，预张力的分布将会非常不均匀，造成施工安装的困难和对受力性能的不利。 第二种思路亦是本文将采用的方法，是将初始预张力的分布状态作为已知数，而把与之对应的平衡的形状作为未知数来求解。此时的目的就是要得到具有给定初始预张力分布的一个平衡的初始形态，为了保证最终得到的预张力分布即为初始假定的预张力，我们将舍弃变形协调条件和材料本构关系。也就是设计给定控制点的位置，以及预张力的分布状态，然后寻求在此条件下的与之相对应的平衡的曲面形状。此种方法的优点是最终得到的初始形态的预张力分布即为初始假定的预张力，这会给膜结构的受力性能和施工安装带来极大的好处。前面的分析可知本文的初始形态设计过程中，将忽略变形协调条件和材料本构关系，方程中与变形协调条件及材料本构关系有关的项均将被舍弃。因此，线性刚度在成形计算过程中将不起作用，结构的刚度仅由非线性刚度构成，如果忽略外荷载作用，则膜结构的初始形态分析的有限元基本方程即为： [KNL].△{u}=-{F} （2） 求解膜结构初始形态的已知条件是，初始预张力（包括索的初始预拉力）的大小及分布、结构的边界控制点位置，要求解的是结构上其余节点的坐标。具体的理论推导和计算步骤详见文献[1]。
（三）静力性能分析 在推导得到含空间膜单元的非线性有限元基本方程后，膜结构的静力性能分析就变得简单了，直接应用式（1）即可，但需要指出的是，由于膜材为非抗压性材料，为此本文引入了由单元主应力判断褶皱的方法。设单元主应力为σ1、σ2，且σ1>σ2，则： σ1>0，σ2 >0，单元正常工作。σ1 <0，σ2 <0，单元退出工作。 σ1 >0，σ2 <0，单元为单向受拉，以转换的应力和刚度代入方程重新迭代。 （四）风振响应分析 首先将膜结构离散化，根据风的概率统计特性，将风速模拟成时间的函数，即人工仿真生成风速时程。然后利用风的空间相关特性，将在结构各个结点处生成的不相关的风速转换成空间相关的风速场。应用Morison公式将风速转变为风压，合理地确定风荷载，并作用于相应的单元结点上，利用有限元法在时间域内直接求解运动微分方程并求得结构的响应。在推导有限元基本方程的迭代公式时，在每一时间步长中考虑了结构的非线性特性、风与结构相互耦合作用的影响、风压分布系数（体型系数）等等因素。由响应值中求得所要的统计信息，如结构振动的位移、速度、加速度的均值和均方差，以及相应的功率谱，并从中获得膜结构的风振响应特性。膜结构质量很轻且刚度较弱，因此非线性风振响应分析将是其设计计算中的又一个关键问题、难点问题。本文应用随机模拟法的思路，发展建立了一种可以适用于大型膜结构非线性风振响应分析的方法。该方法的思路是， 对工程结构设计计算来说，风作用的大小一般以风压来表示，本文采用Morison公式来计算风荷载，它可更合理地反映物体的表面压力。当不考虑风与结构耦合作用时的风荷载公式为： （3）式中第一项与第二项之和是加速流中静止的弹性结构所受到的总惯性力，总称为惯性项。式中 是空气质量密度，A为面积，H为膜厚度， 为风压分布系数， 为风速，包含了平均风速和脉动风速。可以看出，Morison公式比一般公式增加了开头的两项，如果忽略前两项，则即为我们通常所说的风荷载公式。如果考虑风与结构的耦合作用，则上式变为： （4） 分别为结构振动的位移、速度和加速度，而结构的运动方程为： （5）文中质量矩阵M采用集中质量矩阵，阻尼矩阵C采用Rayleigh阻尼，非线性的刚度矩阵K如前文所述，而荷载项P采用式（4）。显然，这是一个复杂的非线性动力方程组，不仅刚度矩阵是随位移变化而变化的，而且荷载项里还包含了结构的加速度和速度项。为此本文利用Newmark法和Newton-Raphson迭代法的思想，推导了膜结构在风荷载作用下，考虑风与结构耦合作用的非线性动力增量平衡方程，最终整理得： （6）其中：
方程中个参数的具体含义及取法可参见文献[1]。
（五）裁剪设计计算 膜结构的曲面经初始形态设计得到，一般为复杂不规则、不可展的空间曲面形式，而且是由有限元的离散节点构成。因此，这里存在一个如何将平面的膜材料拼接成空间曲面的裁剪设计问题。裁剪设计是膜结构设计中的一个关键问题。裁剪下料图的准确与否直接关系到施工安装后的平整度，也就是膜结构的静、动力计算的前提—初始形态设计得到的初始状态与实际施工安装后是否吻合，进而影响静动力的计算结果。 本文发展建立了一种应用广义泛函变分取极值，得到膜结构曲面上的测地线，然后依据测地线计算裁剪线从而生成裁剪图的方法。具体方法详见文献[1]，其过程一般如下： 1、由初始形态设计程序求得膜结构的空间曲面离散点坐标。 2、根据建筑和结构上的要求确定测地线控制点的位置，计算出测地线轨迹坐标。这里建筑上的要求指的一是测地线的布置要美观；二是相邻两条测地线间距离要控制在膜材幅宽的范围内。而结构上的要求是要尽量避免结构在受荷时拉力大的方向与测地线方向垂直。 3、依据测地线求出裁剪线。4、生成裁剪图。当然，最终的施工下料图还要考虑绳边索套，焊缝宽度及初始预张力、温度应力、徐变等因素的影响。
结 语
膜结构是我国空间结构发展的一个主要方向，而比较成熟的计算理论方法和设计程序又是推动膜结构应用的一个主要力量。本文发展建立了一套从初始形态设计～静力性能分析～风振响应分析～裁剪设计计算的膜结构的计算方法和设计程序，并已应用于大型结构工程的计算设计。
参考文献 （1） 向阳，薄膜结构的初始形态设计、风振响应分析及风洞实验研究（博士学位论文），哈尔滨建筑大学，1998 （2） 武岳，向阳，沈世钊，威海体育场挑蓬膜结构风洞测压实验研究，国家自然科学基金重大项目专题年度研究报告（大连），1999