基于变分渐近法的双周期梯形波纹板等效模型

邓兵，钟轶峰，席森彪，禹辽

(重庆大学 土木工程学院；山地城镇建设与新技术教育部重点实验室，重庆 400045)

双周期梯形波纹板是沿两个平面方向均呈周期性梯形波纹变化的新型轻质结构，具有较强吸能质量比，较好抗弯曲性能。与蜂窝等结构相似，双周期梯形波纹结构也可以作为夹层结构的芯层使用，由于纵横波周期分布，双周期梯形波纹板具有很大的惯性矩和抗弯刚度，有效提高了整个结构的强度和承载能力，并大大降低自重而不影响其他使用性能。该结构既满足了强度和刚度的力学性能要求，又能保证材料对吸能和隔热等功能要求[1-5]。

由于波纹板复杂的形貌特征，宏观刚度特性难以通过弹塑性理论方法得到[6]。虽然有限元法可用于评估结构的刚度，但梯形波纹结构形状复杂，数值模拟需要精细网格划分，导致较高的计算成本。由于这类板壳结构的构造和材料分布在面内具有周期性，且周期的尺寸相比整个板的宏观尺寸小很多，所以在宏观尺寸上可以看作是均匀分布的，可使用等效模型以简化复杂波纹板结构的建模[7-9]。

学者们为获得复杂波纹板结构在各种工况载荷下力学行为的等效模型进行了大量的研究。Briassoulis[10]研究了矩形波纹板的等效抗弯刚度。Samanta等[11]在考虑拉伸和弯曲刚度的基础上，对梯形波纹板进行了静态和动态分析。Yokozeki等[12]通过实验和分析研究了由碳环氧复合材料制成的波纹层压板的性能。Peng等[13]通过无网格迦辽金法研究了正弦波形和梯形波纹板的等效弹性性质。Liew等[14]使用该方法进行波纹板的几何非线性分析。张勇等[15]提出了一种将波纹板等效为正交异性板的方法，根据等刚度公式可以将板件的材料参数进行等效。吴存利等[16]将复合材料层合板理论拓展到复合材料波纹板，得到拉伸、耦合和弯曲刚度的解析表达式。高轩能等[17]对槽型波纹板局部屈曲进行了研究，得出整体结构的极限承载力取决于其在半跨荷载作用下槽型波纹截面板的局部屈曲承载力的结论。

现阶段对波纹板的研究主要集中在单周期结构(仅一个方向有周期波纹变化)上，双周期梯形波纹板的力学性能分析缺少必要的研究。笔者利用波纹形状的周期远小于整个板结构尺寸的特点，以凸起为中心的典型代表单元——单胞为研究对象，利用变分渐近法[18-20]建立单胞等效刚度的数值计算模型。再通过均匀化技术将双周期梯形波纹板等效为具有相同刚度特性的正交异性板进行分析。基于该方法分别计算研究不同结构参数下的等效刚度和变形响应，并与三维有限元模拟结果进行对比，以验证等效刚度的精确性。

本文利用结构最小构建单元——单胞作为黑箱，通过渐近均匀化过程得到等效板模型的等效刚度，大大降低了三维有限元求解的繁琐程度，便于工程应用，尤其在求解大规模双周期梯形波纹板刚度中，具有计算简便的优点。

1 双周期梯形波纹板降维分析

对于构造上正交各向异性板的等效建模一般采用理论推导的方法，利用弹性力学理论得到该简化板在各方向的刚度，进而得到等效材料参数[13-14]。本文从能量概念角度建立双周期梯形波纹板的等效模型。由于双周期梯形波纹板沿平面两个坐标轴方向均具有异质性，其典型单胞为三维结构(图1(b))。板的参考面可视为二维连续介质，该连续介质的每个材料点都有一个相关的单胞作为其微观结构。

图1 双周期梯形波纹板降维分析及坐标系Fig.1 Dimensional reduction analysis and coordinate system of bi-periodic trapezoidal corrugated plate

将典型三维单胞分割为l1、l2、l3、l4、l55部分(其中l4=l2，l5=l1)，波纹板高度为h，板厚度为t，如图1(b)所示。引入两组坐标系：全局坐标系xi和局部坐标系yi。对于降维后的板结构，场变量是定义在参考面上的(x1,x2)的函数，x3消失(如图1所示)。全局和局部坐标的关系可表示为yi=xi/η(η是小参数，下标i,j=1,2,3；α=1,2)。在降维模型中，原三维结构的场变量可用剩余宏观坐标xα和局部坐标yj表示，其偏导数为

(1)

为使用变分渐近法进行双周期梯形波纹板分析，首先需用二维应变量表示原结构的三维位移场

(2)

式中：ui、vi分别表示原三维板和二维板的位移；wi是待求波动函数，下划线项为参考面(中面)产生的位移，需满足条件

hv3(xα)=[u3],hvα(xα)=[uα]+[ηy3]v3,2

(3)

式中：[·]表示对单胞体积域积分。

式(3)意味着对波动函数存在约束

〈wi〉=0

(4)

三维线弹性理论的应变场可表示为

(5)

将式(5)代入式(2)，并根据变分渐近法忽略高阶小项，得到三维应变场为

Γ11=ε11+ηy3κ11+w1,1

2Γ12=2ε12+2ηy3κ12+w1,2+w2,1

Γ22=ε22+ηy3κ22+w2,2

2Γ13=w1,3+w3,1

2Γ23=w2,3+w3,2

Γ33=w3,3

(6)

式中：二维面内拉伸应变εαβ与弯曲曲率καβ可分别表示为

εαβ(x1,x2)=(vα,β+vβ,α)/2,καβ(x1,x2)=-v3,αβ

(7)

三维应变场可表示为矩阵形式

Γt=Γ33=w3,3

(8)

控制梯形波纹板线弹性行为的应变能可表示为

(9)

式中：De、Des、Det、Ds、Dst、Dt是三维6×6阶材料矩阵的对应分块矩阵。

外载所做的虚功可表示为

(10)

(11)

(12)

δΠ=δU-δW*=

(13)

式中：仅wi是变化的。

2 梯形波纹板的变分渐近降维分析

为使用变分渐近法求解式(13)中的未知翘曲函数wi，首先需评估各项的阶数。由于梯形波纹板的厚度与边长之比很小(η=h/L≪1)，可利用该小参数进行降维分析。式(13)中各项的阶数评估为

(14)

式中：h与L分别表示板的厚度和边长；n为最小应变量的阶数；μ为材料属性的阶数。

2.1 零阶近似

总势能密度的显式表达式为

(15)

式中：下划线项为h/L阶或更高阶项，在零阶近似时可忽略；双下划线项为常量，不影响求解未知波动函数，也可忽略。

引入拉格朗日乘子λi以考虑对波动函数的约束，

δ(Π+λi〈wi〉)=0

(16)

求解波动函数的零阶近似变分表达式为

(17)

对式(17)部分积分，得到相应的Euler-Lagrange方程

(18)

根据自由表面条件，式(18)方括号内的表达式在板顶、底面应为零，可定义板顶、底面边界条件为

(19)

式中：上标“+/-”表示作用在板顶、底面上的量。

将上述条件代入式(18)，可求解w||,w3

(20)

式中：

(21)

将求解的波动函数代入式(17)，得到

(22)

式中：A为3×3阶拉伸刚度子矩阵;D为3×3阶弯曲刚度子矩阵;B为3×3阶拉弯耦合刚度子矩阵，其具体表达式为

(23)

零阶近似三维应变场可重构为

(24)

应力场重构为

(25)

2.2 一阶近似

由式(25)可知，零阶近似只能得到平面应力，为得到分析破坏至关重要的面外应力，需进行下阶近似。为此，可将零阶翘曲函数摄动为

(26)

将式(26)代入式(13)，得到一阶近似主导项为

(27)

式(27)的Euler-Lagrange方程为

(28)

(29)

式中：

其中：( )±=( )++( )-,( )∓=( )--( )+。

至此，修正到一阶的近似能量为

(30)

式中：

(31)

2.3 三维场重构

降维模型的可靠性取决于其对原三维结构位移、应力和应变场重构的准确性。为此，有必要提供重构关系以完善降维模型，即用二维变量和x3表征三维位移、应力和应变场。

由式(4)，可重构三维位移场为

Ui=ui+x3(C3i-δ3i)+Cjiwj

(32)

式中：Ui、ui分别为三维位移和二维位移。

由式(12)重构三维应变场为

(33)

三维应力场可重构为

(34)

由上述推导过程可知，通过提取出双周期波纹板结构的单胞，分析其单胞结构上的受力特点并结合周期性边界条件，即可得到等效刚度特性(见图2)；再通过均匀化技术将双周期梯形波纹板转换为具有相同刚度特性的正交异性板进行分析；最后，将等效板模型分析得到的全局部位移、应变代入所求位置的三维单胞内进行重构，可得到局部应力、位移场分布。

图2 双周期波纹板等效刚度特性求解过程Fig.2 The solution process of equivalent stiffness of bi-periodic trapezoidal corrugated plate

3 模型验证与讨论

为讨论和分析典型双周期梯形波纹板的刚度特性，同时，对构建的模型进行验证，选取结构参数为(见图3)：高度h=10 mm，厚度t=0.5 mm，周期长度l=30 mm；结构总长度为20个周期600 mm×600 mm；单个周期l1=5 mm，l2=5 mm，l3=10 mm。材料属性为各向同性：弹性模量E=206 GPa，泊松比v=0.3。

图3 结构的单胞及网格划分Fig.3 Unit cell and meshing of the unit cells

考虑周期结构在x1和x2方向上均具有20个周期，采用ABAQUS有限元软件建立三维模型和等效模型，并赋予相应的结构参数。

3.1 模型验证

由于波纹方向与平面坐标一致，A13、A23、D13、D23消失；对于正交各向异性板，x1和x2方向的扭转刚度是不同的，D22取平均值；双周期波纹形状沿x1和x2方向相对坐标原点对称，则拉弯耦合刚度B11=B13=B22=B33=0；因为两个方向周期间距相同，A11=A22，D11=D22。基于构建的等效模型计算得到的等效刚度列于表1，并与普通平板的刚度进行比较。由表1可看出：双周期梯形波纹板的拉伸刚度A较平板分别下降至37%、15%、37%、10%；弯曲刚度D分别增加了2.96、2.23、2.96、15.84倍。这主要是因为双周期梯形波纹板的平面内应变过程需先将梯形拉直，再产生类似平板的拉伸，使得其面内拉伸刚度有所降低；而等效模型为符合Kirchhoff假设，其相对厚度较平板提高很多，等效弯曲刚度较平板有较大增加。在保持总质量不变的前提下，弯曲刚度有较大增加，可实现结构的轻型化。

表1 双周期梯形波纹板与普通平板的刚度比较Table 1 Stiffness comparison between bi-periodic trapezoidal corrugated plates and flat Plates

分别建立结构的实体模型和以等效刚度为基础的等效板模型，分别使用411 804个实体单元(C3D20R)和1 575个壳单元(S8R5)。图4对比了两种模型在不同载荷(跨中施加集中载荷50 N和均匀载荷500 Pa)与边界条件(四边固结CCCC和两边自由两边固结FFCC)下的位移分布。可看出等效模型与三维模型在相同载荷下的位移分布吻合较好，四边固结(CCCC)施加集中力和均布载荷的最大误差分别为1%和1.2%；两边自由两边固结(FFCC)

图4 双周期梯形波纹板在不同载荷和边界条件下的位移比较Fig.4 Displacement comparison bi-periodic trapezoidal corrugated plate under different loads and boundary condtions

施加集中力和均布载荷的最大误差分别为3.17%和1.6%。验算比较最大误差在5%以内，等效模型具有较好的吻合性。因此，基于变分渐近法得到的等效刚度是准确的，可用于评估宏观刚度特性。

基于得到的等效刚度，同时对梯形波纹板和等效平板建模，通过线性屈曲分析求解结构的临界特征值列于表2。

表2 三维实体结构和等效板全局屈曲模态和临界荷载(N)比较(前5阶)Table 2 Comparison of first five buckling loads (N)and mode shapes between three-dimensional solid structures and equivalent plates

续表2

由表2可以看出，等效板模型的全局屈曲模态和三维梯形波纹板的有限元结果相吻合，屈曲临界荷载误差整体较小；随着模态阶数的增加，误差整体出现增大的趋势，原因在于，等效模型的屈曲模态受屈曲波长的影响，但最大误差(第5阶3.51%)仍低于5%。因此，等效模型具有较高的精确性和有效性。

3.2 局部场重构

基于2.3节推导公式，将四边固结CCCC和集中载荷50 N下等效板模型分析得到板中点处(图1(c)所示)的全局位移、应变代入对应三维单胞内进行重构，得到局部应力、位移分布(图5、图6)。值得注意的是，目前大多数简化模型无法预测出准确的局部场分布。

图5 板中点处单胞内局部应力分量分布Fig.5 Distribution of local stress components within the unit cell at the center of the plate

图6 板中点处单胞内局部位移分量分布Fig.6 Distribution of local displacements within the unit cell at the center of the plate

由图5应力分布云图得知，在发生应变时，梯形板的应力分布不均匀，凸起没有完全承受荷载，在x=5 mm和x=25 mm处(平板和凸起的相交处)出现最值，应力在平板和凸起的相交处变化比较明显。σ11和σ22在凸起位置有明显应力增大现象且均为负值。σ33沿平板基本保持不变，但在交点处出现应力突变。这就合理解释了在交点处应力先破坏的现象。

根据图6位移云图，在梯形凸起的顶点处出现极值，单胞突起的对称点呈现相反的位移变化。经过凸起位置，位移会发生突变，随后会保持平板的位移变化特征趋势。重构的U3最大位移与整体模型相一致，最大误差仅为1%。

3.3 几何参数影响分析

双周期梯形波纹板结构参数主要包括l1、l3、h和t。当材料属性不变时，利用图3(b)所示的三维单胞模型讨论双周期梯形波纹板结构参数对等效刚度特性的影响规律，保持单胞凸起部分长度不变，改变参数l1、l3、h和t，研究等效刚度的变化规律。

图7为其他参数不变，l1在6～18 mm范围内变化时结构等效刚度变化情况。由图7可知，随着梯形板的结构参数l1即凸起间距的增大，拉伸刚度逐渐增大，弯曲刚度逐渐降低，且都呈现出非线性变化。主要原因在于，随着凸起间距的增大，板中凸起结构越来越稀疏，结合应力云图可知单胞凸起的组合特性愈发不明显，结构更趋近于普通平板。

图7 凸起间隔长度对结构刚度特性的影响Fig.7 The influence of interval length on the structural stiffness

图8为其他参数不变，高度h在6～18 mm范围内变化时结构等效刚度变化情况。由图8可知，随着h增大，拉伸刚度逐渐减小，弯曲刚度逐渐增大，且弯曲刚度较拉伸刚度变化更加明显。主要原因在于随着梯形凸起高度的增加，结构的薄板特性越不明显，构造异性特征更加突出。

图8 板高度h对结构刚度特性的影响Fig.8 The influence of height h on the structural stiffness

图9为其他参数不变，厚度t在1.0～3.0 mm范围内变化时结构的等效刚度变化情况。由图9可知，结构随着几何参数板厚t的增加，拉伸刚度和弯曲刚度都逐渐增大。拉伸刚度A随着板厚度t呈现出线性变化，且A22增长较快；弯曲刚度D呈现出非线性增长，增长速率逐渐加快。分析其原因是随着板的厚度增加，梯形板的构造正交异性特性逐渐模糊，呈现实体的特性。

图9 板厚t对结构刚度特性的影响Fig.9 The influence of thickness t on the structural stiffness

腰部长度受到l2和l3的影响，为方便分析，只需要改变结构参数l3即可。图10为l1=12 mm、

图10 结构参数l3对刚度特性的影响Fig.10 The influence of l3 on the structural stiffness

h=12 mm、t=2mm、l3在4～12 mm范围内变化时梯形波纹板的等效刚度变化规律。由图10可知，随着梯形板的结构参数l3增大即顶部变宽时，拉伸刚度逐渐减小，弯曲刚度逐渐增大，分析原因是，在顶部变宽的过程中，凸起部分承受水平荷载能力逐渐减小。l3对结构等效刚度的影响比l1对等效刚度影响小。

4 结论

1)基于变分渐近法建立了双周期梯形波纹板的等效模型，得到各刚度系数解析式和局部场重构关系，适用于波纹的周期远小于结构尺寸的情况。由于等效板模型的近似能量与原三维波纹板能量尽可能接近(通过对能量泛函变分主导项的渐近扩展分析加以保证)，可用于不同材料波纹板的计算等效弯矩和等效最大拉、压应力等。

2)将构建的刚度计算方法以及建模方法应用于双周期梯型波纹钢板全局性能分析中，计算得到的全局位移和屈曲模态与Abaqus有限元结果相吻合，但计算效率大为提高。同时，将等效板模型分析得到的全局部位移、应变代入所求位置的三维单胞内进行重构，得到局部应力、位移分布。

3)利用三维单胞模型得到了双周期梯形波纹板结构参数(板高度h、板厚度t和凸起间隔长度l1和腰部长度)对等效刚度的影响规律。随着凸起间距的增大和板高度的减小，板凸起的组合特性愈发不明显，更趋近于普通平板，导致拉伸刚度逐渐增大，弯曲刚度逐渐降低；而随着板厚度和腰部长度的增加，梯形波纹板的构造正交异性特性逐渐降低，弯曲刚度呈上升的趋势。未来可结合隔热、减振等性能，对双周期梯形波纹板的实际应用提供理论依据，使其在土木建筑、航空航天等领域得到更广泛的应用。