CN107766670B - 周期性手征蜂窝结构材料等效弹性模量预测方法 - Google Patents

Abstract

本发明公开了一种周期性手征蜂窝结构材料等效弹性模量预测方法，用于解决现有手征蜂窝结构材料等效弹性模量计算方法实用性差的技术问题。技术方案是取手征性蜂窝结构一个单胞为研究对象，把单胞简化为一系列梁和刚性节圆组成的结构，模拟结构单轴拉伸加载，运用梁的大变形理论计算单胞每个梁的变形，再列出结构整体力和力矩的平衡方程、结构的变形要求方程组，用牛顿迭代法求解方程组可得到结构受力及变形状态，进而得到结构的等效弹性模量。此方法基于结构大变形条件，可以预测周期性多孔材料在发生大变形下的非线性等效等弹性模量。

Description

周期性手征蜂窝结构材料等效弹性模量预测方法

技术领域

本发明涉及一种手征蜂窝结构材料等效弹性模量计算方法，特别涉及一种周期性手征蜂窝结构材料等效弹性模量预测方法。

背景技术

应用于自适应机翼设计中周期性手征蜂窝结构材料韧带长度与厚度比大，具有柔性特征，在分析计算柔性手征蜂窝结构的等效弹性常数时需考虑韧带的几何非线性变形。

文献“A.Spadoni,M.Ruzzene.Elasto-static micropolar behavior of achiral auxetic lattice[J].Journal of the Mechanics and Physics of Solids,60(2012),p156-171”公开了一种六韧带手性蜂窝结构的等效弹性矩阵的计算方法。该方法把结构简化为梁计算出结构的应变能，运用微极理论推导出等效本构方程，进而得出结构的各弹性常数。但是，此方法是基于微极理论，只适用于线弹性结构，无法反应结构在大变形情况下的非线性等效弹性模量等。

发明内容

为了克服现有手征蜂窝结构材料等效弹性模量计算方法实用性差的不足，本发明提供一种周期性手征蜂窝结构材料等效弹性模量预测方法。该方法取手征性蜂窝结构一个单胞为研究对象，把单胞简化为一系列梁和刚性节圆组成的结构，模拟结构单轴拉伸加载，运用梁的大变形理论计算单胞每个梁的变形，再列出结构整体力和力矩的平衡方程、结构的变形要求方程组，用牛顿迭代法求解方程组可得到结构受力及变形状态，进而得到结构的等效弹性模量。此方法基于结构大变形条件，可以预测周期性多孔材料在发生大变形下的非线性等效等弹性模量。

本发明解决其技术问题所采用的技术方案：一种周期性手征蜂窝结构材料等效弹性模量预测方法，其特点是包括以下步骤：

步骤一、手征蜂窝简化为半径为r的圆环及与其相切的长为L的韧带组成，韧带的有效梁变形部分长度为l：

l＝L-2l₁

其中α＝0.9。式中l₁表示韧带与节圆相切部分长度，α为系数。

步骤二、取手征性蜂窝结构单胞，由圆环和六条长为L/2的韧带组成。圆环的圆心距为R，圆环半径为r，圆环和韧带的壁厚为t，结构垂直于纸面的厚度为b，韧带与圆心连线夹角为β，材料本身杨氏模量E_s＝1.6GPa，泊松比v＝0.3。

步骤三、对周期性手征蜂窝结构进行单轴加载模拟，根据结构的变形机理，假设单胞第i条韧带的受力状态：受力大小F_i，力的方向角α_i，并设为未知量，i取值范围为1、2、3。

步骤四、根据结构的变形机理，设计未知量第i条韧带的转角γ_i，计算每条韧带的变形，得到每条韧带变形后的形状及位置；

(a)AB边变形分析。

A₀B₀为变形前位置，曲线AB为变形后状态，A’B为B点处切线。为节圆转角，α₁为A点受力F₁与A₀B₀的夹角，γ为AB上任一点处的转角。

根据梁的弹性弯曲理论，AB的变形微分方程为：

式中，I表示梁的惯性矩，s表示梁上任一点的位置。

化简得到：

式中，C为常数。已知当时，γ＝γ₁转角达到最大，将这一条件代入上式得将C代入上式整理得：

化简得微分表达式：

式中，S＝s/l(0≤S≤0.5)，分别表示无量纲位置和无量纲力。

引入中间量η，令得：

两边积分得到无量纲力表达式：

式中，参数F(m₁)表示第一类完全椭圆积分，F(η₁,m₁)表示第一类不完全椭圆积分。

然后求AB的投影长度。建立局部坐标系XBY，以力F₁方向为X轴，其垂直方向为Y轴。沿力方向及垂直方向投影无量纲长度分别为：

并进一步推导计算，得：

式中，E(m₁)表示第二类完全椭圆积分，E(η₁,m₁)表示第二类不完全椭圆积分。

AB在全局坐标系下的投影长度为：

x₁＝[X₁cos(α₁-β)+Y₁cos(α₁-β-π/2)]·l

y₁＝[X₁sin(α₁-β)+Y₁sin(α₁-β-π/2)]·l

A、B点坐标为：

x_A＝x_B+x₁

y_A＝y_B+y₁

(b)CD边变形分析。

C₀D₀为变形前位置，曲线CD为变形后状态，C’D为D点处切线。为节圆转角，α₂为C点受力F₂与初始位置C₀D₀的夹角，γ为CD上任一点处的转角。

坐标系xOy为全局坐标系。以力F₂方向为X轴，垂直方向为Y轴，建立局部坐标系XDY。

根据梁的弹性弯曲理论，CD的变形微分方程为：

化简得到无量纲力表达式：

式中，参数

求CD在局部坐标系下投影长度。沿力方向及垂直方向投影无量纲长度分别为：

CD在全局坐标系下的投影长度为：

x₂＝[X₂cos(α₂-β+π/3)+Y₂cos(α₂-β-π/6)]·l

y₂＝[X₂sin(α₂-β+π/3)+Y₂sin(α₂-β-π/6)]·l

求C、D点坐标为：

x_C＝x_D+x₂

y_C＝y_D+y₂

(c)EF边变形分析。

E₀F₀为变形前位置，曲线EF为变形后状态，E’F为F点处切线。为节圆转角，α₃为E点受力F₃与初始位置E₀F₀的夹角，γ为EF上任一点处的转角。

坐标系xOy为全局坐标系。以F为原点，以力F₃方向为X轴，垂直方向为Y轴，建立局部坐标系XFY。

根据梁的弹性弯曲理论，EF的变形微分方程为：

化简得无量纲力表达式：

式中，参数

求EF在局部坐标系下投影长度。沿力方向及垂直方向投影无量纲长度分别为：

EF在全局坐标系下的投影长度为：

x₃＝[X₃cos(α₃-β-π/3)+Y₃cos(α₃-β-5π/6)]·l

y₃＝[X₃sin(α₃-β-π/3)+Y₃sin(α₃-β-5π/6)]·l

求E、F点坐标为：

x_E＝x_F+x₃

y_E＝y_F+y₃

步骤五、根据单胞处于平衡状态及单胞变形为中心对称形式，列出力F_i及韧带位置满足的平衡方程组，其未知量为α_i、γ_i。

结构变形后C、E两点关于横向轴线上下对称，变形后仍然满足：

x_C-x_E＝0

y_C+y_E＝0

变形前横向节圆连线水平，变形后已然满足，即：

y_A＝0

变形前C、E点横坐标为A点横坐标的1/2，变形后仍然满足，即：

x_A-2x_E＝0

横向单轴拉伸时，手性蜂窝结构横向排与排之间的等效作用力为零，得：

F₃sin(α₃-β-π/3)-F₂sin(α₂-β+π/3)＝0

单轴拉伸时，假设单胞上中心对称的两条韧带受力大小相等方向相反，因此单胞自然满足静力平衡，同时单胞满足力矩平衡，即：

2F₁x_Asinω₁+2F₂(x_Csinω₂-y_Ccosω₂)+2F₃(x_Esinω₃-y_Ecosω₃)＝0

式中，定义参数为：ω₁＝α₁-β

ω₂＝α₂-β+π/3

ω₃＝α₃-β-π/3

根据前文中的分析，以上两式改写为：

对于手性蜂窝单轴拉伸时，给定一个转角值，6个基本变量：α₁、γ₁、α₂、γ₂、α₃、γ₃；并联立以上六个方程组成非线性方程组。

步骤六、利用牛顿迭代法求解非线性方程组，得出α_i、γ_i各量；

非线性方程组通过Newton-Raphson迭代法求得。为了保证迭代的收敛性，为分步加载。在迭代过程中，第一步迭代各变量的初值取结构线性变形结果的近似值，第k步迭代的初值取第k-1步的结果。至此，得到各个载荷步下的6个基本变量的解。

步骤七、计算结构等效应力、等效应变，得到等效杨氏模量。

求得变形后A点的横坐标后，A点的位移即是手性蜂窝结构的变形，等效应变为：

单轴拉伸下AB、CD、EF三条韧带合力在横向投影即为单胞的等效外力，等效应力为：

改写为：

最后得等效拉伸杨氏模量为：

由此得到六韧带手征蜂窝结构等效杨氏模量计算表达式。

本发明的有益效果是：该方法取手征性蜂窝结构一个单胞为研究对象，把单胞简化为一系列梁和刚性节圆组成的结构，模拟结构单轴拉伸加载，运用梁的大变形理论计算单胞每个梁的变形，再列出结构整体力和力矩的平衡方程、结构的变形要求方程组，用牛顿迭代法求解方程组可得到结构受力及变形状态，进而得到结构的等效弹性模量。此方法基于结构大变形条件，可以预测周期性多孔材料在发生大变形下的非线性等效等弹性模量。

考虑周期性手征蜂窝结构材料的柔性特征，采用本发明方法预测其等效弹性模量时考虑了韧带的几何非线性变形，在手征性蜂窝结构等效应变0～0.2范围内，计算出的等效应力应变曲线与仿真结果很吻合，等效弹性模量与仿真结果误差在10％以内，预测模型在更接近实际变形情况，预测结果更可靠，计算精度更高。

下面结合附图和具体实施方式对本发明作详细说明。

附图说明

图1是本发明周期性手征蜂窝结构材料等效弹性模量预测方法中手征蜂窝示意图。

图2是本发明方法实施例中手征蜂窝单胞受力图。

图3是本发明方法实施例中手征蜂窝AB边受拉变形图。

图4是本发明方法实施例中手征蜂窝CD边受拉变形图。

图5是本发明方法实施例中手征蜂窝EF边受拉变形图。

具体实施方式

参照图1-5。本发明周期性手征蜂窝结构材料等效弹性模量预测方法具体步骤如下：

步骤一、手征蜂窝可简化为半径为r的圆环及与其相切的长为L的韧带组成，韧带的有效梁变形部分长度为l：

l＝L-2l₁

其中α＝0.9。式中l₁表示韧带与节圆相切部分长度α为系数。

步骤二、取手征性蜂窝结构单胞，由圆环和六条长为L/2的韧带组成。圆环的圆心距为R，圆环半径为r，圆环和韧带的壁厚为t，结构垂直于纸面的厚度为b，韧带与圆心连线夹角为β，材料本身杨氏模量E_s＝1.6GPa，泊松比v＝0.3。

步骤三、对周期性手征蜂窝结构进行单轴加载模拟，根据结构的变形机理，假设单胞第i条韧带的受力状态：受力大小F_i，力的方向角α_i，并设为未知量，i取值范围为1、2、3。

步骤四、根据结构的变形机理，设计未知量第i条韧带的转角γ_i，计算每条韧带的变形，得到每条韧带变形后的形状及位置；

(a)AB边变形分析。

A₀B₀为变形前位置，曲线AB为变形后状态，A’B为B点处切线。为节圆转角，α₁为A点受力F₁与A₀B₀的夹角，γ为AB上任一点处的转角。

根据梁的弹性弯曲理论，AB的变形微分方程为：

式中，I表示梁的惯性矩，s表示梁上任一点的位置。

化简得到：

式中，C为一常数。已知当时，γ＝γ₁转角达到最大，将这一条件代入上式得将C代入上式整理得：

化简得微分表达式：

式中，S＝s/l(0≤S≤0.5)，分别表示无量纲位置和无量纲力。

引入中间量η，令可得：

两边积分可得到无量纲力表达式：

式中，参数F(m₁)表示第一类完全椭圆积分，F(η₁,m₁)表示第一类不完全椭圆积分。

然后求AB的投影长度。建立局部坐标系XBY，以力F₁方向为X轴，其垂直方向为Y轴。沿力方向及垂直方向投影无量纲长度分别为：

并进一步推导计算，可得：

式中，E(m₁)表示第二类完全椭圆积分，E(η₁,m₁)表示第二类不完全椭圆积分。

AB在全局坐标系下的投影长度为：

x₁＝[X₁cos(α₁-β)+Y₁cos(α₁-β-π/2)]·l

y₁＝[X₁sin(α₁-β)+Y₁sin(α₁-β-π/2)]·l

A、B点坐标为：

x_A＝x_B+x₁

y_A＝y_B+y₁

(b)CD边变形分析。

C₀D₀为变形前位置，曲线CD为变形后状态，C’D为D点处切线。为节圆转角，α₂为C点受力F₂与初始位置C₀D₀的夹角，γ为CD上任一点处的转角。

坐标系xOy为全局坐标系。以力F₂方向为X轴，垂直方向为Y轴，建立局部坐标系XDY。

根据梁的弹性弯曲理论，CD的变形微分方程为：

化简可得到无量纲力表达式：

式中，参数

求CD在局部坐标系下投影长度。沿力方向及垂直方向投影无量纲长度分别为：

CD在全局坐标系下的投影长度为：

x₂＝[X₂cos(α₂-β+π/3)+Y₂cos(α₂-β-π/6)]·l

y₂＝[X₂sin(α₂-β+π/3)+Y₂sin(α₂-β-π/6)]·l

求C、D点坐标为：

x_C＝x_D+x₂

y_C＝y_D+y₂

(c)EF边变形分析。

E₀F₀为变形前位置，曲线EF为变形后状态，E’F为F点处切线。为节圆转角，α₃为E点受力F₃与初始位置E₀F₀的夹角，γ为EF上任一点处的转角。

坐标系xOy为全局坐标系。以F为原点，以力F₃方向为X轴，垂直方向为Y轴，建立局部坐标系XFY。

根据梁的弹性弯曲理论，EF的变形微分方程为：

化简可得无量纲力表达式：

式中，参数

求EF在局部坐标系下投影长度。沿力方向及垂直方向投影无量纲长度分别为：

EF在全局坐标系下的投影长度为：

x₃＝[X₃cos(α₃-β-π/3)+Y₃cos(α₃-β-5π/6)]·l

y₃＝[X₃sin(α₃-β-π/3)+Y₃sin(α₃-β-5π/6)]·l

求E、F点坐标为：

x_E＝x_F+x₃

y_E＝y_F+y₃

步骤五、根据单胞处于平衡状态及单胞变形为中心对称形式，可列出力F_i及韧带位置满足的平衡方程组，其未知量为α_i、γ_i。

结构变形后C、E两点关于横向轴线上下对称，变形后仍然满足：

x_C-x_E＝0

y_C+y_E＝0

变形前横向节圆连线水平，变形后已然满足，即：

y_A＝0

变形前C、E点横坐标为A点横坐标的1/2，变形后仍然满足，即：

x_A-2x_E＝0

横向单轴拉伸时，手性蜂窝结构横向排与排之间的等效作用力为零，可得：

F₃sin(α₃-β-π/3)-F₂sin(α₂-β+π/3)＝0

单轴拉伸时，假设单胞上中心对称的两条韧带受力大小相等方向相反，因此单胞自然满足静力平衡，同时单胞也应该满足力矩平衡，即：

2F₁x_Asinω₁+2F₂(x_Csinω₂-y_Ccosω₂)+2F₃(x_Esinω₃-y_Ecosω₃)＝0

式中，定义参数为：ω₁＝α₁-β

ω₂＝α₂-β+π/3

ω₃＝α₃-β-π/3

根据前文中的分析，以上两式可改写为：

对于手性蜂窝单轴拉伸时，给定一个转角值，6个基本变量：α₁、γ₁、α₂、γ₂、α₃、γ₃；并联立以上六个方程组成非线性方程组。

步骤六、利用牛顿迭代法求解非线性方程组，得出α_i、γ_i各量；

非线性方程组通过Newton-Raphson迭代法求得。为了保证迭代的收敛性，应该分步加载。在迭代过程中，第一步迭代各变量的初值可取结构线性变形结果的近似值，第k步迭代的初值取第k-1步的结果。至此，可以得到各个载荷步下的6个基本变量的解。

步骤七、计算结构等效应力、等效应变，得到等效杨氏模量。

求得变形后A点的横坐标后，A点的位移即是手性蜂窝结构的变形，等效应变为：

单轴拉伸下AB、CD、EF三条韧带合力在横向投影即为单胞的等效外力，等效应力为：

可改写为：

最后可得等效拉伸杨氏模量为：

由此得到六韧带手征蜂窝结构等效杨氏模量计算表达式，该式考虑了韧带的几何非线性变形，体现了周期性手征手蜂窝结构材料的柔性特征。

Claims (1)

1.一种周期性手征蜂窝结构材料等效弹性模量预测方法，其特征在于包括以下步骤：

步骤一、手征蜂窝简化为半径为r的圆环及与其相切的长为L的韧带组成，韧带的有效梁变形部分长度为l：

l＝L-2l₁

其中α＝0.9；式中l₁表示韧带与节圆相切部分长度，α为系数；

步骤二、取手征性蜂窝结构单胞，由圆环和六条长为L/2的韧带组成；圆环的圆心距为R，圆环半径为r，圆环和韧带的壁厚为t，结构垂直于纸面的厚度为b，韧带与圆心连线夹角为β，材料本身杨氏模量E_s＝1.6GPa，泊松比v＝0.3；

步骤三、对周期性手征蜂窝结构进行单轴加载模拟，根据结构的变形机理，假设单胞第i条韧带的受力状态：受力大小F_i，力的方向角α_i，并设为未知量，i取值范围为1、2、3；

步骤四、根据结构的变形机理，设计未知量第i条韧带的转角γ_i，计算每条韧带的变形，得到每条韧带变形后的形状及位置；

(a)AB边变形分析；

A₀B₀为变形前位置，曲线AB为变形后状态，A’B为B点处切线；为节圆转角，α₁为A点受力F₁与A₀B₀的夹角，γ为AB上任一点处的转角；

根据梁的弹性弯曲理论，AB的变形微分方程为：

式中，I表示梁的惯性矩，s表示梁上任一点的位置；

化简得到：

式中，C为常数；已知当时，γ＝γ₁转角达到最大，将这一条件代入上式得将C代入上式整理得：

化简得微分表达式：

式中，S＝s/l(0≤S≤0.5)，分别表示无量纲位置和无量纲力；

引入中间量η，令得：

两边积分得到无量纲力表达式：

式中，参数F(m₁)表示第一类完全椭圆积分，F(η₁,m₁)表示第一类不完全椭圆积分；

然后求AB的投影长度；建立局部坐标系XBY，以力F₁方向为X轴，其垂直方向为Y轴；沿力方向及垂直方向投影无量纲长度分别为：

并进一步推导计算，得：

式中，E(m₁)表示第二类完全椭圆积分，E(η₁,m₁)表示第二类不完全椭圆积分；

AB在全局坐标系下的投影长度为：

x₁＝[X₁cos(α₁-β)+Y₁cos(α₁-β-π/2)]·l

y₁＝[X₁sin(α₁-β)+Y₁sin(α₁-β-π/2)]·l

A、B点坐标为：

x_A＝x_B+x₁

y_A＝y_B+y₁

(b)CD边变形分析；

C₀D₀为变形前位置，曲线CD为变形后状态，C’D为D点处切线；为节圆转角，α₂为C点受力F₂与初始位置C₀D₀的夹角，γ为CD上任一点处的转角；

坐标系xOy为全局坐标系；以力F₂方向为X轴，垂直方向为Y轴，建立局部坐标系XDY；

根据梁的弹性弯曲理论，CD的变形微分方程为：

化简得到无量纲力表达式：

式中，参数

求CD在局部坐标系下投影长度；沿力方向及垂直方向投影无量纲长度分别为：

CD在全局坐标系下的投影长度为：

x₂＝[X₂cos(α₂-β+π/3)+Y₂cos(α₂-β-π/6)]·l

y₂＝[X₂sin(α₂-β+π/3)+Y₂sin(α₂-β-π/6)]·l

求C、D点坐标为：

x_C＝x_D+x₂

y_C＝y_D+y₂

(c)EF边变形分析；

E₀F₀为变形前位置，曲线EF为变形后状态，E’F为F点处切线；为节圆转角，α₃为E点受力F₃与初始位置E₀F₀的夹角，γ为EF上任一点处的转角；

坐标系xOy为全局坐标系；以F为原点，以力F₃方向为X轴，垂直方向为Y轴，建立局部坐标系XFY；

根据梁的弹性弯曲理论，EF的变形微分方程为：

化简得无量纲力表达式：

式中，参数

求EF在局部坐标系下投影长度；沿力方向及垂直方向投影无量纲长度分别为：

EF在全局坐标系下的投影长度为：

x₃＝[X₃cos(α₃-β-π/3)+Y₃cos(α₃-β-5π/6)]·l

y₃＝[X₃sin(α₃-β-π/3)+Y₃sin(α₃-β-5π/6)]·l

求E、F点坐标为：

x_E＝x_F+x₃

y_E＝y_F+y₃

步骤五、根据单胞处于平衡状态及单胞变形为中心对称形式，列出力F_i及韧带位置满足的平衡方程组，其未知量为α_i、γ_i；

结构变形后C、E两点关于横向轴线上下对称，变形后仍然满足：

x_C-x_E＝0

y_C+y_E＝0

变形前横向节圆连线水平，变形后已然满足，即：

y_A＝0

变形前C、E点横坐标为A点横坐标的1/2，变形后仍然满足，即：

x_A-2x_E＝0

横向单轴拉伸时，手性蜂窝结构横向排与排之间的等效作用力为零，得：

F₃sin(α₃-β-π/3)-F₂sin(α₂-β+π/3)＝0

单轴拉伸时，假设单胞上中心对称的两条韧带受力大小相等方向相反，因此单胞自然满足静力平衡，同时单胞满足力矩平衡，即：

2F₁x_Asinω₁+2F₂(x_Csinω₂-y_Ccosω₂)+2F₃(x_Esinω₃-y_Ecosω₃)＝0

式中，定义参数为：

根据前文中的分析，以上两式改写为：

对于手性蜂窝单轴拉伸时，给定一个转角值，6个基本变量：α₁、γ₁、α₂、γ₂、α₃、γ₃；并联立以上六个方程组成非线性方程组；

,

其中

步骤六、利用牛顿迭代法求解非线性方程组，得出α_i、γ_i各量；

非线性方程组通过Newton-Raphson迭代法求得；为了保证迭代的收敛性，为分步加载；在迭代过程中，第一步迭代各变量的初值取结构线性变形结果的近似值，第k步迭代的初值取第k-1步的结果；至此，得到各个载荷步下的6个基本变量的解；

步骤七、计算结构等效应力、等效应变，得到等效杨氏模量；

求得变形后A点的横坐标后，A点的位移即是手性蜂窝结构的变形，等效应变为：

单轴拉伸下AB、CD、EF三条韧带合力在横向投影即为单胞的等效外力，等效应力为：

改写为：

最后得等效拉伸杨氏模量为：

由此得到六韧带手征蜂窝结构等效杨氏模量计算表达式。