CN109033705B - 一种空间褶皱薄膜弹性应变能和主应力有限元分析方法 - Google Patents

Abstract

本发明属于薄膜技术领域，公开了一种空间褶皱薄膜弹性应变能和主应力有限元分析一种空间薄膜结构褶皱的分析方法。本发明首先对薄膜结构褶皱研究的三种数值方法进行了比较，并与试验结果进行了对比，验证了方法的正确性，得出动力显示法是一种较为合适的褶皱分析方法的结论；采用动力显示法分析了在不同膜材厚度、不同网格密度和不同剪切位移下Kapton矩形剪切薄膜褶皱的主应力和弹性应变能分布规律及原因，得出以下结论：单元网格大小仅仅影响计算的精度和计算量，而影响薄膜褶皱的两个重要参数是膜材厚度和剪切位移，对薄膜褶皱的主应力大小、主应力波形曲线振幅、弹性应变能密度以及薄膜总弹性应变能影响都较大，其规律对工程实际具有重要参考价值。

Description

一种空间褶皱薄膜弹性应变能和主应力有限元分析方法

技术领域

本发明属于薄膜技术领域，尤其涉及一种空间褶皱薄膜弹性应变能和主应力有限元分析一种空间薄膜结构褶皱的分析方法。

背景技术

目前，业内常用的现有技术是这样的：

由于薄膜材料抗弯刚度极小，受到某种荷载作用致使内部产生压应力时，很容易发生平面外屈曲，产生褶皱。对于薄膜结构的褶皱研究，国内外学者主要集中在褶皱数量、幅度、半波长和剪切角度的影响因素以及应力的分布、褶皱对薄膜动态特性的影响等方面。主要研究理论有张力场理论、分叉屈曲理论和有限质点理论，对于褶皱的研究主要采用有限元分析方法，常见的有限元分析方法包括非线性屈曲法、直接扰动法和动力显示法等。

YWesleyWang和SergioPellegrino通过有限元软件ABAQUS结合实验研究了矩形薄膜受剪切和方形薄膜对角张拉两种情况的褶皱分布和应力分布；王震开发了向量式有限元分析方法并研究了褶皱扩展变形全过程；殷亮通过有限元软件ABAQUS结合实验研究了矩形薄膜受剪切的褶皱幅值、褶皱半波长、褶皱角度和褶皱数量的影响因素；李云良通过有限元软件ANSYS研究了矩形薄膜受剪切的二次屈曲行为、薄膜屈曲缺陷敏感性；马瑞通过有限元软件ANSYS/LS-YNA研究了褶皱对薄膜动力特性的影响。

薄膜褶皱形成过程中弹性应变能和主应力分布分析，有助于实际工程中薄膜结构规避高弹性应变能带来的破坏力，对于褶皱的形成机理的探究以及褶皱的预测与控制有着重要意义，能为薄膜结构的设计提供参考。但是国内外文献对褶皱产生时弹性应变能和应力的分布和影响报道很少。

综上所述，现有技术存在的问题是：

从文献资料来看，国内外对于薄膜结构的褶皱研究，主要集中在褶皱数量、幅度、半波长和剪切角度的影响因素以及应力的分布、褶皱对薄膜动态特性的影响等方面，缺乏对薄膜褶皱的形成机理的探究，因此，本发明试图从能量和应力的角度对褶皱的形成机理进行深入研究，为褶皱的预测与控制提供重要参考；

现有薄膜因其厚度薄、抗弯刚度小，当外载荷产生的压应力超过临界压应力时，就会出现局部屈曲现象而形成褶皱。褶皱的出现不仅会影响薄膜结构的表面精度，还会改变结构的应力分布状态，进而影响结构的静力和动力特性。

解决上述技术问题的难度和意义：

空间结构中常常使用的Kapton和Mylar薄膜厚度都非常薄，在太空环境或施工过程中，都会产生压应力，一旦外载荷产生的压应力超过临界应力，薄膜便会产生局部屈曲的现象，即褶皱。褶皱的存在，会使得表面张力减小、刚度降低，从而影响用薄膜制成的空间结构如反射面、太阳帆以及充气天线等结构的表面精度，直接影响收发信号能力的强弱，因此如何消除褶皱的形成和产生，是一项艰巨而重要的任务，而褶皱的形成离不开能量和应力。

本发明试图从能量和应力的角度研究褶皱的形成机理，对于褶皱的预测和防范具有直接而重要的意义。

发明内容

针对现有技术存在的问题，本发明提供了一种空间褶皱薄膜弹性应变能和主应力有限元分析一种空间薄膜结构褶皱的分析方法。

本发明是这样实现的，

一种空间褶皱薄膜弹性应变能和主应力有限元分析方法，包括：

以非线性有限元软件ABAQUS为平台，运用非线性屈曲法、直接干扰法以及动力显示法分别对薄膜褶皱的数量、面外变形及褶皱幅值进行分析，并与试验结果进行对比，验证数值分析方法的正确性；

以非线性有限元软件ABAQUS为平台，运用非线性屈曲法、直接干扰法以及动力显示法分别对薄膜褶皱的数量、面外变形及褶皱幅值进行分析，包括：

首先建模，建模过程参考图2的流程进行，具体建模过程见表1。

表1三种方法的建模过程

利用ABAQUS为平台进行计算，再对结果进行比较分析，尽管都是利用ABABQUS为平台，但由于三种方法的原理不同，采用算法和添加面外扰动的方式不同，因此其计算结果不同，计算结果见表2。

表2褶皱信息表

	褶皱数量(个)	中部褶皱幅值(mm)	中部褶皱波长(mm)	中部褶皱剪切角度(°)
					实验结果	13	0.481	25.7
非线性屈曲法	13	0.489	27.1	45.5
					直接扰动法	13	0.472	24.7	45.5
动力显示法	13	0.476	26.4	45.5

并与试验结果进行对比，包括：

将计算得到的结果(如表2所示)与文献(李云良，谭惠丰，王晓华.矩形薄膜和充气管的屈曲及后屈曲行为分析[J].航空学报，2008，29(4)：886-892)的试验结果，如褶皱数量、中部褶皱幅值、中部褶皱波长以及中部褶皱剪切角度等数值进行比较，看吻合的程度怎样。比较的结果是：褶皱数量与试验相同，都是13，但采用非线性屈曲法、直接扰动法和动力显示法计算的1-1截面中部褶皱幅值分别为0.489mm、0.472mm和0.476mm，与试验值0.481mm相比，分别增加1.5％、减小1.9％和减小1％。因此从比较的结果来看，动力显示法最接近试验结果，且动力显示法没有采用整体切线刚度矩阵，不需要平衡迭代和收敛控制。此方法分析速度快，建模和分析过程简单，容易收敛。因此动力显示法是一种最合适的分析褶皱的方法。

在非线性屈曲法、直接干扰法和动力显示法中选取最合适的方法，在不同网格密度、不同膜材厚度和不同剪切位移下对矩形薄膜受剪切褶皱的能量和应力分布进行进一步分析，包括：

第一步，分别在不同网格密度、不同膜材厚度和不同剪切位移等三种不同工况进行建模，同样参考图2的流程进行。

第二步，计算三种不同工况下薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量。

第三步，比较三种工况下计算得到的薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量与网格密度、膜材厚度和剪切位移等参数之间的关系，总结规律，得出有益结论，为薄膜褶皱的预测和防范提供工程实践参考。

进一步，所述非线性屈曲法包括：

首先通过模态屈曲分析得到薄膜结构的屈曲模态，再把模态的线性组合作为几何初始缺陷添加到薄膜结构中，触发褶皱的形成；

在采用Newton-Raphson算法，初始几何缺陷通过面外变形的方式添加到薄膜上时，组成方式包括：

式中Δz为面外位移值；ω_i为第i阶模态；

为折减系数。

进一步，所述直接扰动法包括：

在满足公式(2)的情况下撤去面外扰动力最为合适；

式中δ₁为总荷载；δ₂为预载；δ′₁为撤去扰动力时的荷载，δ″₁为撤去扰动力后的荷载。

进一步，所述动力显示法包括：

采用动力显示分析步和显式中心差分法对运动方程进行求解，方程如下：

式中u为位移；

为速度；ü为加速度；M为质量矩阵；P为荷载向量；I为内力向量；其中Δt^(n+1/2)＝t⁽ⁿ⁺¹⁾-t⁽ⁿ⁾，Δt^(n-1/2)＝t⁽ⁿ⁾-t^(n-1)。

本发明另一目的在于提供一种计算机程序，所述计算机程序运行所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法。

本发明另一目的在于提供一种终端，所述终端至少搭载实现所述空间褶皱薄膜弹性应变能和主应力有限元分析方法的控制器。

本发明另一目的在于提供一种计算机可读存储介质，包括指令，当其在计算机上运行时，使得计算机执行所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法。

本发明另一目的在于提供一种实现所述空间褶皱薄膜弹性应变能和主应力有限元分析方法的空间褶皱薄膜弹性应变能和主应力有限元分析系统。

本发明的优点及积极效果为：

本发明首先对薄膜结构褶皱研究的三种数值方法进行了比较，分析了各自的特点，并与试验结果进行了对比，验证了其正确性，得出动力显示法是一种较为合适的褶皱分析方法的结论。其次，采用动力显示法分析了在不同膜材厚度、不同网格密度和不同剪切位移下Kapton矩形剪切薄膜褶皱的主应力和弹性应变能分布规律及原因，得出以下结论：

(1)单元网格大小仅仅影响计算的精度和计算量，对薄膜褶皱的主应力大小、主应力波形曲线振幅、弹性应变能密度以及薄膜总弹性应变能影响都较小，其规律为：随着单元网格尺寸的增大，其值略微减小。

(2)膜材厚度对最大主应力影响较小，对最小主应力、主应力波形曲线振幅影响较大，其规律为：膜材厚度增大，最小主应力和主应力波形曲线振幅都增大；膜材厚度对弹性应变能密度和薄膜总弹性应变能影响也较大：随着膜材厚度增大，弹性应变能密度和薄膜总弹性应变能明显增大。

(3)剪切位移是研究的三个参数中影响最大的一个，随着剪切位移的增大，最大和最小主应力增大，主应力波形曲线振幅也增大；同时，剪切位移对弹性应变能的影响也更大，随着剪切位移增大，弹性应变能密度和薄膜总弹性应变能都明显增大。

(4)基于上述规律，膜材厚度和剪切位移是影响薄膜褶皱的两个重要参数，因此在实际工程中，为规避高弹性应变能和高主应力带来的破坏，可考虑采用较薄的膜材厚度，且形状控制方面尽量避免剪切位移所带来的不良影响。

附图说明

图1是本发明实施提供的空间褶皱薄膜弹性应变能和主应力有限元分析一种空间薄膜结构褶皱的分析方法流程图。

图2是本发明实施提供的三种方法流程图。

图3是本发明实施提供的矩形薄膜简图。

图4是本发明实施提供的三种方法与实验的面外变形图。

图中：(a)、非线性屈曲法；(b)、直接干扰法；(c)、动力显示法。

图5是本发明实施提供的不同单元网格尺寸下总弹性应变能图。

图6是本发明实施提供的不同膜材厚度下最大与最小主应力图。

图7是本发明实施提供的不同膜材厚度下总弹性应变能图。

图8是本发明实施提供的不同剪切位移下最大与最小主应力图。

图9是本发明实施提供的不同剪切位移下总弹性应变能图。

图10是本发明实施提供的1-1截面面外变形与主应力对比图。

图11是本发明实施提供的1-1截面面外变形与弹性应变能密度对比图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

下面结合附图对本发明的应用原理作进一步描述。

如图1所示，本发明实施例提供的空间褶皱薄膜弹性应变能和主应力有限元分析方法，包括以下步骤：

S101，以非线性有限元软件ABAQUS为平台，运用非线性屈曲法、直接干扰法以及动力显示法分别对薄膜褶皱的数量、面外变形及褶皱幅值进行了分析，并与试验结果进行了对比，验证了数值分析方法的正确性；

S102，在非线性屈曲法、直接干扰法和动力显示法中选取最合适的方法，在不同网格密度、不同膜材厚度和不同剪切位移下对矩形薄膜受剪切褶皱的能量和应力分布进行进一步分析，以分析薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量情况。

以非线性有限元软件ABAQUS为平台，运用非线性屈曲法、直接干扰法以及动力显示法分别对薄膜褶皱的数量、面外变形及褶皱幅值进行分析，包括：

首先建模，建模过程参考图2的流程进行，具体建模过程见表1。

表1三种方法的建模过程

利用ABAQUS为平台进行计算，再对结果进行比较分析，尽管都是利用ABABQUS为平台，但由于三种方法的原理不同，采用算法和添加面外扰动的方式不同，因此其计算结果不同，计算结果见表2。

表2褶皱信息表

并与试验结果进行对比，包括：

将计算得到的结果(如表2所示)与文献(李云良，谭惠丰，王晓华.矩形薄膜和充气管的屈曲及后屈曲行为分析[J].航空学报，2008，29(4)：886-892)的试验结果，如褶皱数量、中部褶皱幅值、中部褶皱波长以及中部褶皱剪切角度等数值进行比较，看吻合的程度怎样。比较的结果是：褶皱数量与试验相同，都是13，但采用非线性屈曲法、直接扰动法和动力显示法计算的1-1截面中部褶皱幅值分别为0.489mm、0.472mm和0.476mm，与试验值0.481mm相比，分别增加1.5％、减小1.9％和减小1％。因此从比较的结果来看，动力显示法最接近试验结果，且动力显示法没有采用整体切线刚度矩阵，不需要平衡迭代和收敛控制。此方法分析速度快，建模和分析过程简单，容易收敛。因此动力显示法是一种最合适的分析褶皱的方法。

在非线性屈曲法、直接干扰法和动力显示法中选取最合适的方法，在不同网格密度、不同膜材厚度和不同剪切位移下对矩形薄膜受剪切褶皱的能量和应力分布进行进一步分析，包括：

第一步，分别在不同网格密度、不同膜材厚度和不同剪切位移等三种不同工况进行建模，同样参考图2的流程进行。

第二步，计算三种不同工况下薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量。

第三步，比较三种工况下计算得到的薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量与网格密度、膜材厚度和剪切位移等参数之间的关系，总结规律，得出有益结论，为薄膜褶皱的预测和防范提供工程实践参考。

本发明实施例提供的非线性屈曲法如下：

非线性屈曲法是首先通过模态屈曲分析得到薄膜结构的屈曲模态，再把模态的线性组合作为几何初始缺陷添加到薄膜结构中，从而触发褶皱的形成。在采用Newton-Raphson算法，初始几何缺陷通过面外变形的方式添加到薄膜上时，其组成方式如下：

式中Δz为面外位移值；ω_i为第i阶模态；

为折减系数。

非线性屈曲法特点是建模过程复杂，分析过程繁琐，而且在分析过程中容易出现不收敛的情况，需反复调整参数和多次重启动分析。

本发明提供的直接扰动法如下：

直接扰动法是直接在薄膜结构上施加一组垂直于薄膜平面、方向彼此相反、随机分布的小荷载，使薄膜产生一定的面外变形，以此作为初始几何缺陷，触发褶皱的产生。同样采用Newton-Raphson算法，偶然荷载的取值应较小，与膜材厚度为同一数量级最为合适，正反两面偶然荷载的数量应相等，刚好抵消。研究发现：在满足公式(2)的情况下撤去面外扰动力最为合适。

式中δ₁为总荷载；δ₂为预载；δ′₁为撤去扰动力时的荷载，δ″₁为撤去扰动力后的荷载。

直接扰动法的特点是建模过程较为简单，但计算量大，分析过程繁琐，收敛困难，也需多次进行参数调整和重启动分析才能达到计算要求。

本发明提供的动力显示法如下：

与前两种方法不同，动力显示法不需要采用特殊的手段触发褶皱的形成，只需要利用特征模态或者扰动力扰动网格，使得薄膜结构的变形更光滑，更符合实际。动力显示法在后屈曲分析阶段，是采用动力显示分析步和显式中心差分法对运动方程进行求解，其方程如下：

式中u为位移；

为速度；ü为加速度；M为质量矩阵；P为荷载向量；I为内力向量；其中Δt^(n+1/2)＝t⁽ⁿ⁺¹⁾-t⁽ⁿ⁾，Δt^(n-1/2)＝t⁽ⁿ⁾-t^(n-1)；

动力显示法特点是没有采用整体切线刚度矩阵，不需要平衡迭代和收敛控制。此方法分析速度快，建模和分析过程简单，容易收敛。

本发明提供一种实现所述空间褶皱薄膜弹性应变能和主应力有限元分析方法的空间褶皱薄膜弹性应变能和主应力有限元分析系统。

下面结合具体分析对本发明作进一步描述。

1、三种褶皱方法的比较

非线性屈曲法、直接扰动法和动力显示法是研究薄膜褶皱最常用的三种有限元分析方法，其分析流程如图1所示。非线性屈曲法与直接扰动法算法采用的都是Newton-Raphson算法，区别在于二者添加面外扰动的方式不同；而动力显示法与前两种方法相比，不仅面外扰动的方式不同，而且算法也不同，详细内容介绍如下。

图2为三种方法流程图。

1.1非线性屈曲法

非线性屈曲法是首先通过模态屈曲分析得到薄膜结构的屈曲模态，再把模态的线性组合作为几何初始缺陷添加到薄膜结构中，从而触发褶皱的形成。在采用Newton-Raphson算法，初始几何缺陷通过面外变形的方式添加到薄膜上时，其组成方式如下：

式中Δz为面外位移值；ω_i为第i阶模态；

为折减系数。

非线性屈曲法特点是建模过程复杂，分析过程繁琐，而且在分析过程中容易出现不收敛的情况，需反复调整参数和多次重启动分析。

1.2直接扰动法

直接扰动法是直接在薄膜结构上施加一组垂直于薄膜平面、方向彼此相反、随机分布的小荷载，使薄膜产生一定的面外变形，以此作为初始几何缺陷，触发褶皱的产生。同样采用Newton-Raphson算法，偶然荷载的取值应较小，与膜材厚度为同一数量级最为合适，正反两面偶然荷载的数量应相等，刚好抵消。研究发现：在满足公式(2)的情况下撤去面外扰动力最为合适。

式中δ₁为总荷载；δ₂为预载；δ′₁为撤去扰动力时的荷载，δ″₁为撤去扰动力后的荷载。

直接扰动法的特点是建模过程较为简单，但计算量大，分析过程繁琐，收敛困难，也需多次进行参数调整和重启动分析才能达到计算要求。

1.3动力显示法

与前两种方法不同，动力显示法不需要采用特殊的手段触发褶皱的形成，只需要利用特征模态或者扰动力扰动网格，使得薄膜结构的变形更光滑，更符合实际。动力显示法在后屈曲分析阶段，是采用动力显示分析步和显式中心差分法对运动方程进行求解，其方程如下：

式中u为位移；

为速度；ü为加速度；M为质量矩阵；P为荷载向量；I为内力向量；其中Δt^(n+1/2)＝t⁽ⁿ⁺¹⁾-t⁽ⁿ⁾，Δt^(n-1/2)＝t⁽ⁿ⁾-t^(n-1)；

动力显示法特点是没有采用整体切线刚度矩阵，不需要平衡迭代和收敛控制。此方法分析速度快，建模和分析过程简单，容易收敛。

1.4算例分析

为了有效说明这三种方法的分析流程和特点，现通过一具体算例进行详细比较分析。

算例一：如图3所示，一矩形Kapton薄膜，长380mm，宽128mm，厚0.05mm，弹性模量和泊松比分别为3530MPa和0.34，薄膜密度为1400Kg/mm³。边界条件为固定AB、CD两长边，AC、BD两短边仅固定旋转自由度，通过张拉薄膜AB边产生预应力，张拉距离为δ₂＝0.1mm。加载条件为沿AB边水平向右移动δ₁＝2mm距离，以形成剪切变形。

为便于比较分析，三种方法采用的单元网格尺寸均为

建模过程参考图2的流程进行。

图4为采用三种方法得到的1-1截面的面外变形与实验结果的对比。图中：(a)、非线性屈曲法；(b)、直接扰动法；(c)、动力显示法。

从图4可以看出，除边缘部分存在差异外，1-1截面的面外变形与试验结果吻合较好，三种方法计算得到的褶皱数量与试验相同，都是13个，但采用非线性屈曲法、直接扰动法和动力显示法计算的1-1截面中部褶皱幅值分别为0.489mm、0.472mm和0.476mm，与试验值0.481mm相比，分别增加1.5％、减小1.9％和减小1％。

为了得到计算结果，非线性屈曲法和直接扰动法经过了多次参数调整和重启动分析才使得结果收敛，而动力显示法分析过程没有遇到不收敛问题。比较而言，动力显示法建模简单，计算速度快，比其他两种方法更容易收敛。因此，动力显示法在计算量较大，计算次数较多的情况下模拟分析过程中更具优势。

2、下面结合矩形受剪切薄膜褶皱弹性应变能与主应力分析对本发明作进一步描述。

2.1主应力和薄膜应变能

根据张力场理论，薄膜褶皱状态基于主应力的判断准则为：

①当σ₂≤0时，薄膜处于张紧状态；

②当σ₁＞0，σ2≤0时，薄膜处于褶皱状态；

③当σ₁≤0时，薄膜处于松弛状态。

其中，σ₁、σ₂分别为最大、最小主应力。基于能量变分原理，结构褶皱变形过程中总能量包括三部分：薄膜张拉应变能、薄膜弯曲应变能以及外载荷功。假设所有外载荷作用在薄膜边界上，考虑到在边界上，薄膜虚位移为零，因此基于薄膜褶皱特征，从能量的角度出发，薄膜褶皱能量看作薄膜面内张拉应变能与弯曲应变能之和。其大小分别由公式(4)和(5)确定：

其中：

式中：h为膜材厚度，E、υ分别为弹性模量和泊松比，u、v、w分别为x、y、z方向的位移，q为失稳波数。

本发明将对不同网格密度、不同膜材厚度、不同剪切位移三种工况下Kapton薄膜受剪切产生的褶皱弹性应变能和主应力分布情况进行分析，其结论对于褶皱的形成机理的探究以及实际工程对褶皱的预测和控制具有十分重要的参考价值。

2.2计算模型与建模方法

算例二：基本模型与算例一相似(如图3)，但采用的Kapton薄膜参数不同，如表3所示，以便与文献进行对比分析。从前述三种有限元分析方法的比较可知，动力显示法在薄膜褶皱的复杂分析中具有一定的优势，因此这里采用动力显示法进行分析，其建模过程参考图2的流程图。

表3膜材参数表

Tab.3Membrane material parameters

2.3不同网格密度下褶皱主应力与弹性应变能分析

工况一：假设薄膜厚度为0.05mm，剪切位移为3mm，分别对

五种不同网格密度下薄膜褶皱的弹性应变能和主应力分布进行分析。

不同网格密度下薄膜褶皱的最大和最小主应力分布情况中，从整体对比情况，截面最大和最小主应力分布情况。

可以看出，截面最大主应力波形曲线的平衡位置大约在44MPa附近，如以单元网格边长3mm的对应值为基数，可以得出规律：单元网格尺寸每增大1mm，其振幅约减小6％。同样地最小主应力波形曲线平衡位置大约在零附近，单元网格尺寸每增大1mm，其振幅约减小5％。

图5为网格密度与薄膜总弹性应变能的对应关系，从图5可以看出，不同单元网格密度的薄膜总弹性应变能非常接近，大致范围在0.550～0.566J之间，且单元网格边长每增加1mm，薄膜总弹性应变能仅仅增加0.3％左右。究其原因是因为有限元网格划分越细，其结果越趋近于理论值。

2.4不同膜材厚度下褶皱主应力与弹性应变能分析

工况二：假设薄膜单元网格密度为

剪切位移为3mm，分别对0.025mm、0.05mm、0.075mm、0.1mm、0.125mm五种不同膜材厚度下薄膜褶皱的弹性应变能和主应力分布进行分析。

图6为不同膜材厚度下薄膜褶皱的主应力大小及分布，从图中可以看出最大主应力都在44MPa处上下波动，这是因为最大主应力主要受拉应力控制，而张拉应力并没有因膜材厚度而变化。如以厚度为0.025mm的褶皱振幅为基数，则厚度每增加0.05mm，振幅大概增大1倍。最小主应力在零值上下波动，其变化规律一样，即厚度每增加0.05mm，振幅也增大1倍，这是因为厚度增大，弯曲变形增大，则褶皱振幅增大，从而导致最小主应力增大，又由于最小主应力受弯曲变形控制，但薄膜剪切变形中弯曲变形占的比例非常小，因此最小主应力变化不大。

图7为不同膜材厚度下薄膜总弹性应变能与厚度的关系。从图7可以看出，膜材厚度对总弹性应变能影响较大，厚度每增大1倍，总弹性应变能增大1倍，近似于线性关系，这可以从公式(4)、(5)体现的膜材厚度与总弹性应变能关系得到相同结论。

2.5不同剪切位移下褶皱主应力与弹性应变能分析

工况三：假设薄膜厚度为0.05mm，单元网格密度为

分别对1mm、1.5mm、2mm、2.5mm、3mm五种不同剪切位移下薄膜褶皱的弹性应变能和主应力进行分析。

图8为不同剪切位移下薄膜的最大和最小主应力分布情况，其规律为：以剪切位移1mm的对应值为基数，则剪切位移每增加1mm，最大主应力约增大14MPa，其波形曲线振幅约增大1倍，这是因为剪切位移的增大引起张拉变形的增大。最小主应力在零值附近上下波动，剪切位移每增加1mm，其振幅约增大1倍，这是因为剪切位移增大，弯曲变形增大，因此褶皱振幅增大，导致最小主应力增大，但变化很微小。图9为不同剪切位移下薄膜总弹性应变能与剪切位移的关系。从图9可以看出，剪切位移对总弹性应变能的影响较大，近似呈二次曲线关系，即剪切位移由1mm增大到2mm时，总弹性应变能增加0.183J，由2mm增大到3mm时，总弹性应变能增加0.196J。这是因为剪切位移增大导致公式(4)、(5)中的面内位移增大。

2.6褶皱局部位置弹性应变能与主应力分布

工况四：在算例二中假设单元网格密度为

膜材厚度为0.5mm，剪切位移为3mm，试分析褶皱波峰、波谷和平衡位置处主应力和弹性应变能情况。

正因为弹性应变能和主应力分布对于褶皱的形成有着举足轻重的作用，因此，有必要分析褶皱局部位置，如波峰、波谷和平衡位置处弹性应变能和主应力分布特点，这有助于弄清褶皱的形成机理，以规避高弹性应变能带来的破坏，为薄膜结构的设计提供参考。

图10、11分别为1-1截面面外变形与主应力、弹性应变能密度的对比。从图中可以看出，最大主应力在沿褶皱方向的对角区域最大，达到163.6MPa，而在垂直褶皱方向对角区域最小，只有-2.384MPa，这是因为沿褶皱方向膜材受拉，垂直褶皱方向膜材受压。最小主应力一般在褶皱波峰和波谷区域较大，能达到±7MPa，而在平衡位置为零，这是因为平衡位置膜材几乎没有发生弯曲，而波峰、波谷处弯曲程度最大。

同样可以看出，沿褶皱方向对角区域的弹性应变能密度最大，达到0.022J/mm³，而垂直褶皱方向对角区域应变能密度最小，只有0.0022J/mm³，其他区域能量分布比较均匀，大致范围为

而且波峰、波谷处弹性应变能变化很小，说明弹性应变能的主要来源是薄膜的受拉应变。

在上述实施例中，可以全部或部分地通过软件、硬件、固件或者其任意组合来实现。当使用全部或部分地以计算机程序产品的形式实现，所述计算机程序产品包括一个或多个计算机指令。在计算机上加载或执行所述计算机程序指令时，全部或部分地产生按照本发明实施例所述的流程或功能。所述计算机可以是通用计算机、专用计算机、计算机网络、或者其他可编程装置。所述计算机指令可以存储在计算机可读存储介质中，或者从一个计算机可读存储介质向另一个计算机可读存储介质传输，例如，所述计算机指令可以从一个网站站点、计算机、服务器或数据中心通过有线(例如同轴电缆、光纤、数字用户线(DSL)或无线(例如红外、无线、微波等)方式向另一个网站站点、计算机、服务器或数据中心进行传输)。所述计算机可读取存储介质可以是计算机能够存取的任何可用介质或者是包含一个或多个可用介质集成的服务器、数据中心等数据存储设备。所述可用介质可以是磁性介质，(例如，软盘、硬盘、磁带)、光介质(例如，DVD)、或者半导体介质(例如固态硬盘SolidState Disk(SSD))等。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (9)

1.一种空间褶皱薄膜弹性应变能和主应力有限元分析方法，其特征在于，所述空间褶皱薄膜弹性应变能和主应力有限元分析方法包括：

以非线性有限元软件ABAQUS为平台，运用非线性屈曲法、直接干扰法以及动力显示法分别对薄膜褶皱的数量、面外变形及褶皱幅值进行分析，并与试验结果进行对比，验证数值分析方法的正确性；

在非线性屈曲法、直接干扰法和动力显示法中选取最合适的方法，在不同网格密度、不同膜材厚度和不同剪切位移下对矩形薄膜受剪切褶皱的能量和应力分布进行进一步分析，分析薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量。

2.如权利要求1所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法，其特征在于，在非线性屈曲法、直接干扰法和动力显示法中选取最合适的方法，在不同网格密度、不同膜材厚度和不同剪切位移下对矩形薄膜受剪切褶皱的能量和应力分布进行进一步分析，具体包括：

第一步，分别在不同网格密度、不同膜材厚度和不同剪切位移三种不同工况进行建模；

第二步，计算三种不同工况下薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量；

第三步，比较三种工况下计算得到的薄膜弹性应变能密度分布、主应力分布以及整个薄膜所具有的能量与网格密度、膜材厚度和剪切位移参数之间的关系。

3.如权利要求1所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法，其特征在于，所述非线性屈曲法包括：

首先通过模态屈曲分析得到薄膜结构的屈曲模态，再把模态的线性组合作为几何初始缺陷添加到薄膜结构中，触发褶皱的形成；

在采用Newton-Raphson算法，初始几何缺陷通过面外变形的方式添加到薄膜上时，组成方式包括：

式中Δz为面外位移值；ω_i为第i阶模态；

为折减系数。

4.如权利要求1所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法，其特征在于，所述直接扰动法包括：

在满足公式(2)的情况下撤去面外扰动力最为合适；

式中δ₁为总荷载；δ₂为预载；δ₁′为撤去扰动力时的荷载，δ₁″为撤去扰动力后的荷载。

5.如权利要求1所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法，其特征在于，所述动力显示法包括：

采用动力显示分析步和显式中心差分法对运动方程进行求解，方程如下：

式中u为位移；

为速度；

为加速度；M为质量矩阵；P为荷载向量；

I为内力向量；其中Δt^(n+1/2)＝t⁽ⁿ⁺¹⁾-t⁽ⁿ⁾，Δt^(n-1/2)＝t⁽ⁿ⁾-t^(n-1)。

6.一种计算机程序，其特征在于，所述计算机程序运行权利要求1～5任意一项所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法。

7.一种终端，其特征在于，所述终端至少搭载实现权利要求1～5任意一项所述空间褶皱薄膜弹性应变能和主应力有限元分析方法的控制器。

8.一种计算机可读存储介质，包括指令，当其在计算机上运行时，使得计算机执行如权利要求1-5任意一项所述的空间褶皱薄膜弹性应变能和主应力有限元分析方法。

9.一种实现权利要求1所述空间褶皱薄膜弹性应变能和主应力有限元分析方法的空间褶皱薄膜弹性应变能和主应力有限元分析系统。

Images