等应力膜在膜结构设计中占有重要位置, 其膜面的找形分析[1]等价于求解极小曲面[2]问题。膜结构的极小曲面找形分析是一个非线性问题, 求解时需要反复迭代, 并需要一个合理的初始解以保证迭代的收敛。常见的找形方法 (如力密度法[3―4]、动力松弛法[5―6]和有限元法[7]) 中, 初始解通常由边界条件进行双线性插值[8]获得, 虽然该法比较简单, 但构造出的初始解曲面往往与真实解曲面的偏差较大, 用于找形分析中非线性迭代的初始解, 常造成迭代收敛困难, 甚至不收敛的情况。因此, 如何获得一个高质量的初始解, 成为困扰工程分析人员的一个难题。

本文提出一种线性化近似方法, 可给出膜结构极小曲面高质量的近似解。该近似解既可用于初步设计阶段了解膜面的大概形状, 亦可作为精细的找形分析中非线性迭代的初始解。其基本思路是利用积分中值定理和归一化手段对曲面面积的表达式进行简化和近似, 将原非线性问题转化为一线性问题。

有限元线法 (简称FEMOL) 作为一种基于常微分方程 (ODE) 求解的半解析方法[9], 用于膜结构的找形分析具有独到的潜在优势:解曲面在单元内光滑而且高阶连续, 利于后续的剪裁分析;常微分方程组的求解由ODE求解器COLSYS[10]自动完成, 无需人为控制求解过程;由于是半离散, 求解精度高。文[11]基于极小曲面找形的思想, 将该法用于膜结构的找形分析, 用有限元线法求解极小曲面这一非线性问题。但其初始解亦是人为设定的, 亦难以确保初始解的质量。

本文首先介绍所提出的线性化方法, 然后用有限元线法[9]求解该线性问题, 即可获得极小曲面的近似解。数值算例显示, 本法效果颇佳。

图1所示为一个一般的膜曲面单元示意图, 在参数坐标系 (ξ, η) 下, 膜曲面Ae的面积[12]为:

其中:g1、g2分别为曲面上任意一点P (x, y, z) 的协变基向量, 表达式为:

其中和以下 () ξ=∂ () /∂ξ, () η=∂ () /∂η;e1、e2、e3分别为x方向、y方向、z方向的单位矢量。被积函数可表示为:

其中:

由于ab>0、sinφ>0, 由积分中值定理有:

其中 (sinφ) e为积分的中值, 该值事先未知, 但当单元划分合理, 单元形状接近平行四边形时, φ在单元上变化不大, (sinφ) e可由单元角点上的值平均近似。记φie (i=1, 2, 3, 4) 为单元e四个角点上的φ值, 且记:

则有:

当单元有边未知待求时, 可用连接该边两端点的直线代替该边, 求得角点上近似的φ值。又, 记分别为单元上a、b的平均值, 则都是在1附近变化的量, 因此有:

从而:

进而有:

其中:

为反映单元两个方向长度比的参数。单元形状越接近平行四边形, 上述近似引起的误差就越小。单元上的值事先未知, 可取a、b在边界上的平均值代替。如图1中,

式中S代表弧长。当单元有边未知待求时, 可用连接单元四个端点的直边长度代替, 据此求得λ。记:

则有:

可见, A0e是对Ae的一个很好的简化近似, 且当单元越接近平行四边形, 上述近似引起的误差就越小, 但A0e是正定二次泛函, 而Ae为非线性泛函。记:

则有:

等应力膜曲面的找形分析等价于数学上极小曲面问题[2]:给定曲面边界坐标, 寻求一个满足边界条件并使曲面面积A最小的曲面, 即要求解:

由式 (17) 可以导出求解膜结构极小曲面找形的控制微分方程, 它是一组非线性方程, 求解时须提供有效的初始解。由于A0是对A很好的近似, 可考虑用使A0取得极小值的曲面来近似原极小曲面。即求解:

由于0A为正定二次泛函, 将常规的离散手段引入上式, 导出的控制方程是一组线性方程, 易于求解且无需迭代。本文将有限元线法这一半离散化方法引入该正定二次泛函极小值问题的求解, 导出的控制方程是一线性常微分方程组, 其解即为极小曲面的近似解。

图2所示为有限元线法的膜单元映射, 其中参数坐标ξ和η的定义域均为[-1, 1]。单元结线上的坐标向量为:

其中p为单元的次数。膜单元内任意一点的坐标 (x, y, z) 采用Lagrange插值得到:

其中:

由式 (20) 不难得到坐标x、坐标y和坐标z对ξ、η的各阶偏导数。

对于FEMOL膜单元, 式 (18) 的变分式可写为:

由式 (4) 可以得到:

将式 (20) 代入式 (23) 进而代入式 (22) 中, 形式上完成ξ方向的积分, 并在η方向上作必要的分部积分, 注意到{δx}e、{δy}e和{δz}e在η=±1的边界上均为0, 有:

其中:

由δA0=0可以导出一组二阶线性常微分方程:

其中, {和为给定的边界值。将边界上给定的结线解引入式 (26) 并整理, 方程右端出现已知边界结线生成的荷载项, 即得到一组标准的线性二阶常微分方程组:

注意, 式 (27) 中的矩阵和向量是引入边界结线解之后的量, 与式 (26) 中的有所不同;再有, 式 (27) 中三组方程是不耦合的, 因此可以分别求解, 这使得问题的求解更加高效。

如前文所述, 本文的目标是为极小曲面找形问题提供高质量的近似解。本节给出若干典型的数值算例用以展示本文方法的可靠性和有效性。本文使用ODE求解器COLSYS求解本文的线性常微方程体系。算例中, COLSYS的误差限tol统一取为10-5, 全部采用6次FEMOL单元。

例1. Scherk曲面。

Scherk曲面是为数不多的具有解析表达式的极小曲面之一, 其表达式为[13]:

本例中取k=1.0, c=5.0, , 四条边界曲线给定。本文用一个6次元即获得了很好的近似解, 如图3、图4所示, 可见精度很高。

例2. 悬链面。

旋转悬链面亦是有解析表达式的极小曲面, 本例求解下述悬链面[11]:

在a≤x, y≤5a, 0≤z≤2.29243167a部分的曲面, 并取参数a=10, 四条边界给定。本文采用一个6次FEMOL单元计算, 求出的近似解曲面与精确的悬链面如图5所示, x=y面的截线比较如图6所示。易见, 虽然本例单元与平行四边形的形状相差较大 (长边与短边之比高达5倍) , 本文方法仍然获得了逼近程度很高的近似解。将本文的近似解用作文[11]FEMOL找形分析的初始解, 迭代迅速收敛, 得到与精确解一致的极小曲面。

例3. 双曲抛物面。

结构外形投影为正方形, 投影后的边长距离为1.0, 对角点高度为0.25, 给定四个顶点坐标为1 (0, 0, 0.25) 、2 (1, 0, 0) 、3 (1, 1, 0.25) 和4 (0, 1, 0) , 边界条件为四条给定的直线边界。本文用一个6次元求解, 得到的解曲面与用ANSYS作找形分析得到的解曲面均列于图7中, 再次显示出线性解的极高逼近度。将本文的近似解用作文[11]FEMOL找形分析的初始解, 迭代亦迅速收敛, 得到与ANSYS解一致的极小曲面。

例4.鞍形曲面。

图7所示为一非规则鞍形曲面, 6条边界给定, 其中四条直线边界如图8所示, 两条曲边界分别为y=-30和面内的两条抛物线, 其方程均为z=10 (1-x2/900) 。本文用两个单元求解:底面为正方形的区域和底面为平行四边形的区域各划分为一个6次元。得到的解曲面如图9所示。用该解作文[11]方法的初始解, 得到的极小曲面亦示于图9中。值得一提的是, 对于本例, ANSYS未能给出收敛的结果。

例5.锥形曲面。

本例旨在检验本文方法对三角形区域的适用性。考察如下锥形曲面:顶点坐标是 (0, 0, 10) , 底边为半径为10的半圆。本文用两个6次元求解 (单元对称分布) , 得到的近似解曲面如图10 (a) 所示。用该解作文[11]方法的初始解, 得到的极小曲面亦列于图10 (b) 中。值得一提的是, 对于本例, ANSYS亦未能给出收敛的结果。

本例中单元面已退化为曲面三角形, 形状与平行四边形相距甚远, 但本文方法仍然给出了较为满意的曲面形状, 足见本文方法有很大的适用空间。

本文对膜结构极小曲面找形分析提出一线性化近似简化算法。通过对曲面面积公式的近似简化, 将原非线性问题近似为一线性问题, 用有限元线法求解该线性问题, 无需迭代, 一步即可得到逼近度很高的曲面解。该解既可用于初步设计阶段大致了解曲面的形状, 亦可用于精细的找形分析中的初始解。虽然本文方法是从接近于平行四边形的单元出发推导的, 但数值算例表明本文方法具有更广的适用范围, 且给出的解答质量很高, 用作精细的找形分析的初始解, 可确保找形分析的成功。本文方法对其他全离散或半离散极小曲面找形分析方法亦有借鉴意义和应用价值。

本文方法完善了膜结构有限元线法极小曲面的找形分析, 将本文方法与文[11]结合, 可实现对极小曲面完整的找形分析, 并使整个求解过程更加自动、规范、通用, 是一种简单、实用、高效、可靠的初始解曲面确定方法。