桁架结构向量式有限元程序计算实现

同济大学软件学院 吴 瀚

1 绪论

计算机的诞生改变了世界的各行各业，计算机强大的计算能力无疑大量减少了人工的劳动力。在建筑行业中，随着力学理论的发展，很多的结构体的形变都可以通过数学公式求解出来，而计算机程序相比人工计算提高了极大的效率，并能将数据通过图形学展现在人们的眼前，使之更加的直观。向量式有限元（VFIFE）分析理论是近些年由美国普渡大学丁承先教授等人提出的一个创新概念，它是一个为解决大形变、大变位、碰撞、倒塌等力学行为而诞生的理论，在处理这些特殊问题上，不像传统理论的千变万化，向量式有限元的控制方程就是基于牛顿第二定律，对于一些拥有简单力学知识的人来说，就可以学会对这样的理论加以运用[1]。本文主要研究了空间桁架结构模型的向量式有限元分析，试图得到解决这类问题的一般性步骤以及程序计算的代码，为以后软件的行为分析模块奠定基础。

2 桁架结构模型数据分析

空间桁架结构（t r uss st r uct ur e）是建筑工程领域常见的一种格构化梁式结构。常用于大跨度的厂房、展览馆、体育馆和桥梁等公共建筑中。在空间桁架这类的结构体模型中，部件都是一根根刚性杆件或者柔性杆件，它们由铰接点或者刚接点连接组成。结点信息主要包括它们在空间坐标系中的坐标以及它们的类型（铰接点还是刚接点），杆件信息包括杆件两端的坐标以及材料性质，由于柔性杆件是可以划分成一个个刚性杆件组成，所以这里做一般性考虑就只做刚性杆件的定义。除了模型本身的数据之外，还有加载在模型上力的信息需要定义，力的定义包括力的三要素，即力的大小、方向和作用点。

对于结点信息，文件中的定义包括编号i d，结点类型t ype，位置信息posit ion.x，posit ion.y，posit ion.z。对于杆件信息，在文件中的定义包括编号i d，第一个端点的结点编号nodei d1，第二个端点的结点编号nodeid2，杨氏模量Young，截面积ar ea，密度densit y。对于力的信息，在文件中的定义包括编号id，作用点信息posi t i on.x、posi t i on.y、posi t i on.z，方向信息dir ect ion.x，dir ect ion.y，dir ect ion.z，力直接作用的杆件编号bar id。

由于力可以是不断变化的，所以需要在文件中定义力的变化信息，这里用另外一个文件去单独描述这样的信息。其中在val ue中定义的编号为3的力，在0.0-10.0时间段中，力的函数表达式就是t/10.0*10.0，在10.0之后的时间段中大小都是10.0，-1表示无穷大，编号为1的力就是一个大小不变的力，大小一直是1.0。dir ect ion中的定义也是一样的道理，最后三项是三个坐标的表达式。posi t ion中定义的最后一项代表力作用点距离力所作用杆件第二个端点的距离的随时间变化函数表达式。

3 程序计算的实现

程序计算主要依据就是模型中的数学公式，所以在编写程序的过程中，程序的通用化是需要解决的问题，程序采用面向对象的编程思想，将杆件、节点、力都抽象成一个类，由于程序目前以实现功能为主，所以并没有过多的继承层次，后面会展开说明程序将来的类图结构，在读入模型文件的数据之后，各个杆件、节点、力的数据实例化之后分别保存在一个数组中。

它主要是通过不断地循环迭代计算出各个时间点的帧数据，所以程序的框架也就是一个大循环，在循环之前需要求出n=0的时候点的坐标。循环的间隔就是时间步长，时间步长可以通过时间步长计算程序计算得到。由于这个公式适用于铰接点，符合牛顿第二定律的，所以在使用这个公式的时候需要判断节点类型，有的节点还有额外的约束，比如有的节点规定y坐标不变，那么这个公式就不能适用于它的y坐标计算。由公式可以得知，循环中的主要工作就是求出内力和外力的合力，在下一个迭代中带入差分公式中求解。在每个时刻的时候都需要将求出的点位置保存下来[3]～[5]。

程序中主要根据公式分为两大步，第一步是求解出节点数据初始值，然后在大循环中迭代求出每个时刻的值，其中h是根据步长公式求解得到。大循环中主要工作是求解每个节点内力和节点外力，然后代入差分公式得到结果。内力求解需要根据当前杆件的变形量得到杆件内力，再通过将内力沿杆件方向和垂直杆件方向分解，得到两端内力的向量值，将每个节点的向量值累加，最终得到当前时刻的内力合力。外力是通过解析加载公式得到，找到当前时刻对应分段函数中的某一段函数，然后将时间代入函数，解析表达式，从而得到当前时刻加载力的三要素。在计算出外力合力和内力合力之后，将每个节点对应的值代入公式，计算出当前时刻的位置数据。并且更新位置数据，在保存数据的时候只需保存nowpos数组即可。

4 计算结果分析

后处理软件是项目组已经完成的一个可实际应用的产品，本软件开发是为土木工程设计工作者进行结构力学行为分析后处理之用，对向量式结构力学前处理生成的海量数据进行解析，结合OpenGL技术对三维模型进行动态再现，主要用于对结构和行为进行准确的分析和预测，如高楼和桥梁受外力后的空间运动，大变位、大变形等行为进行准确预测。

将本文程序计算得到的node文件和el ement文件导入后处理软件中，得到的图像结果进行比较可以看出，程序计算的结果符合实际情况，向量式有限元理论的程序代码实现解决了结构体模型在传统理论上很棘手的大形变问题，本文的程序计算还停留在试验阶段，对于一些参数的设置还需要更加合理化，但是计算过程是没有问题的。本文通过对向量式有限元的程序计算实现，说明了这个理论的实施步骤，同时验证了理论的可实施性，在对解决空间杆件形变问题上找出了一般的程序步骤，将抽象的物理模型具体化到文件中，实现了从理论到实践的一个转换。

5 结论

本文通过简单模型的建立，说明了向量式有限元的计算步骤，程序可以遵循这样一个步骤去实现，最后将输出结果作为后处理的输入，通过后处理软件展示计算结果生成的模拟动画，验证了程序的可靠性。

向量式有限元理论目前还没有实现该理论的软件，所以需要我们不断地探索，这将是一个很有意义的课题，它能够提供给结构工程师们一个模拟仿真的平台。本文作为探索性地实现了刚性杆件的程序计算，这只是一个开始，向量式有限元中还有针对更多的模型进行探讨。

[1]丁承先,段元峰,吴东岳.向量式结构力学[M].上海:科学出版社,2012,09.

[2]卢哲刚,姚谏.向量式有限元——一种新型的数值方法[J].空间结构,2012,18(001):85-89.

[3]Ting E C,Shih C,Wang Y K.Fundamentals of a vector form intrinsic finite element:PartI. basic procedure and a plane frame element[J].Journal of Mechanics,2004,20(2):113-122.

[4]Ting E C,Shih C,Wang Y K.Fundamentals of a vector form intrinsic finite element:Part II.plane solid elements[J].Journal of Mechanics,2004,20(02):123-132.

[5]Shih C,Wang Y K,Ting E C.Fundamentals of a vector form intrinsic finite element:Part III.Convected material frame and examples[J].Journal of Mechanics,2004, 20(02):133-143.