绝对节点坐标法下斜率不连续问题处理方法讨论

张君茹，程耿东; Zhang Junru; Cheng Gengdong

引 en

绝对节点坐标法下斜率不连续问题处理方法讨论^*

张君茹
，程耿东

1. 大连理工大学工程力学系, 大连 116024；

DISCUSSION ON SLOPE DISCONTINUITY USING THE ABSOLUTE NODAL COORDINATE FORMULATION ^∗

Zhang Junru
， Cheng Gengdong

1. Department of Engineering Mechanics, Dalian University of Technology, Dalian 116024 , China；

DOI:10.6052/1672⁃6553⁃2019⁃060

全文
评论
参考文献
出版信息

参考文献 1

Shabana A A.Dynamics of multibody systems.New York:Cambridge University Press,2013

查找原文

参考文献 2

田强,刘铖,李培,等.多柔体系统动力学研究进展与挑战.动力学与控制学报,2017,15(5):385 ~ 405(Tian Q,Liu C,Li P,et al.Advances and challenges in dynamics of flexible multibody systems.Journal of Dynam⁃ ics and Control,2017,15(5):385~ 405(in Chinese))

查找原文

参考文献 3

Shabana A A.An absolute nodal coordinate formulation for the large rotation and deformation analysis of flexible bodies.Technical Report.no.MBS96⁃1⁃UIC,1996

查找原文

参考文献 4

Gerstmayr J,Sugiyama H,Mikkola A.Review on the ab⁃ solute nodal coordinate formulation for large deformation a⁃nalysis of multibody systems.Journal of Computational & Nonlinear Dynamics,2013,8(3):031016

查找原文

参考文献 5

Shabana A A,Mikkola A M.Use of the finite element ab⁃ solute nodal coordinate formulation in modeling slope dis⁃ continuities.Multibody System Dynamics,2003,125(2):357~ 387

查找原文

参考文献 6

Shabana A A,Maqueda L G.Slope discontinuities in the finite element absolute nodal coordinate formulation:gra⁃ dient deficient elements.Multibody System Dynamics,2008,20(3):239~ 249

查找原文

参考文献 7

Maqueda L G,Shabana A A.Numerical investigation of the slope discontinuities in large deformation finite element formulations.Nonlinear Dynamics,2009,58(1 ~ 2):23 ~ 37

查找原文

参考文献 8

Shabana A A.General method for modeling slope disconti⁃ nuities and T⁃sections using ANCF gradient deficient finite elements.Journal of Computational & Nonlinear Dynam⁃ ics,2011,6(2):267~ 274

查找原文

参考文献 9

Bayoumy A,Nada A,Megahed S.Methods of modeling slope discontinuities in large size wind turbine blades u⁃ sing absolute nodal coordinate formulation.Proceedings of the Institution of Mechanical Engineers Part K Journal of Multi⁃body Dynamics,2014,228(3):314~ 329

查找原文

参考文献 10

Romero I,Arribas J J.A simple method to impose rota⁃ tions and concentrated moments on ANC beams.Multibody System Dynamics,2009,21(4):307~ 323

查找原文

参考文献 11

Berzeri M,Shabana A A.Study of the centrifugal stiffe⁃ ning effect using the finite element absolute nodal coordi⁃ nate formulation.Multibody System Dynamics,2002,7(4):357~ 387

查找原文

参考文献 12

Yu L,Zhao Z,Tang J,et al.Integration of absolute nodal elements into multibody system.Nonlinear Dynamics,2010,62(4):931~ 943

查找原文

参考文献 13

章孝顺,章定国,陈思佳,等.基于绝对节点坐标法的大变形柔性梁几种动力学模型研究.物理学报,2016,65(9):148~ 157(Zhang X S,Zhang D G,Chen S J,et al.Several dynamic models of a large deformation flexible beam based on the absolute nodal coordinate formulation.Acta Physica Sinica,2016,65(9):148 ~ 157(in Chi⁃ nese))

查找原文

目录contents

摘要Abstract
关键词Keywords
0 引言
1 斜率不连续问题描述
1.1 绝对节点法缩减梁单元
1.2 折梁斜率不连续问题及非线性约束方程
1.3 直梁斜率不连续问题及非线性约束方程
2 文献中的处理办法
2.1 新坐标系的选择
2.2 新坐标系的应用
2.3 新坐标系带来的影响
3 算例
3.1 算例一
3.2 算例二
4 结论
参考文献

摘要

Shabana 提出的绝对节点坐标法,引入节点斜率坐标作为节点自由度描述转动.对于由梁板壳及块体组成的组合结构,在结构节点处相交单元的节点斜率自由度不连续,这给组合结构的建模和分析带来特殊的困难.本文讨论了文献中研究斜率不连续问题时的处理办法.在简要介绍绝对节点坐标法后,详细地讨论了经典折梁算例和截面呈阶梯变化的直梁算例中斜率不连续问题.对这两个算例,本文采用约束函数法和现有文献中的转换坐标方法,计算了在结构节点处相交杆件的轴向应变,对比这些数值结果,本文指出现有文献中的转换坐标办法,忽视了斜率自由度和转角自由度的差别,从而不能正确给出斜率不连续处相交杆件的轴向应变,需要进一步研究.

Abstract

The absolute nodal coordinate formulation(ANCF) proposed by Shabana describes the large rotation by introducing slope coordinates as nodal degrees of freedom. For composite structures composed of beams, plates,shells or blocks,the slope coordinates at the nodes of the elements intersecting at the structural joint points are discontinuous, which causes the difficulties in the modeling and analysis of composite structures. In this paper the methods of dealing with the problem of slope discontinuity in the literature were discussed. Firstly, the abso⁃ lute nodal coordinate method was briefly introduced. Secondly, the slope discontinuity problems in the classical example of folded beam and in the straight beam with stepped section were discussed in detail. Then by numeri⁃ cally calculating the axial strain on the neutral layer of two examples, the constrained function method was com⁃ pared with the coordinate transformation method in the existing literature. It was pointed out that the coordinate transformation method in current literatures neglects the difference between the slope degree of freedom and the rotation degree of freedom,and thus the correct axial strains cannot be achieved at the joint node with slope dis⁃ continuity, which needs a further study.

关键词

绝对节点坐标法；节点斜率不连续；应变抹平

Keywords

absolute nodal coordinate formulation ； slope discontinuity ； strain smoothing

0 引言
基于浮动坐标或绝对节点坐标的多体动力学方法是柔性多体动力学两类主流方法^{^[1,2]}.绝对节点坐标法( Absolute Nodal Coordinate Formulation, ANCF)是Shabana等人1996 年^{^[3]} 为了处理柔性体大转动和大变形问题而提出的,该方法引入节点斜率坐标作为节点自由度描述转动,不需要将柔性体的全局位置矢量分解为大范围刚体运动矢量和局部坐标下的弹性变形矢量,所建立的多体动力学方程的广义质量矩阵为常数矩阵,不包含科氏和离心惯性项,而且该列式能精确地描述结构刚体运动.
对于由梁板壳及块体组成的组合结构,采用基于绝对节点法进行有限元建模时,斜率不连续是常见的问题^{^[4]}.以杆系结构中的刚架结构为例,由于节点是刚节点,与同一个节点连结的杆件除了各杆端点的位移自由度相等,它们截面转角也相等.但是,斜率自由度并不是转角,转角相同并不等同于斜率自由度相等,一般地说,这些杆件的斜率自由度不相等.在结构力学的有限元法中,利用节点处各杆端点的位移和转角相等,可以进行单元刚度矩阵的组装,大幅度减少独立自由度数量,提高求解效率.但是,采用绝对节点坐标法时,如何由它们的节点自由度表示出转角相等的连接条件,并应用于多体系统的建模和分析是个困难问题^{^[10]}.该问题在文献中简称为斜率不连续( slope discontinuity).处理这样的斜率不连续问题的一个方法,是由绝对节点坐标法中各杆的节点自由度(包括位移及斜率自由度)表示出各杆端面的转角,再将刚架节点处各杆断面的转角相等的单元连接条件处理成约束函数,采用拉格朗日乘子引入动力学方程中,使得原方程变成带非线性约束的多体动力学方程.当多体系统中的斜率不连续很多时,例如大量杆件组成的刚架,动力学方程求解速度和收敛效率会降低.
Shabana和Mikkola等人^{^[5]} 认识到斜率不连续问题对于绝对节点坐标法的应用带来的困难,所以,他们提出坐标转换的思想,引入中间单元坐标系,将相连结杆件的不连续的斜率自由度转换成中间单元坐标系中的自由度,采用的坐标转换矩阵是常数矩阵^{^[6,7]},这样做不仅仅保证了质量矩阵是常数矩阵的性质,而且将单元之间非线性的连接条件转换成关于节点自由度的线性函数,实现了将斜率不连续结构带来的非线性约束消去,避免了斜率不连续问题带来的非线性约束.Shabana等人^{^[5]} 对全参数绝对节点法单元引入中间单元坐标系,但是此时的中间单元坐标系是全局坐标系,不适用于缩减单元,然后,Shabana和Maqueda等人^{^[6]} 提出一种新的中间单元坐标系,此时引入的中间单元坐标系是初始时刻未变形的单元坐标系,该坐标转换不仅仅适用于全参数单元,也适用于缩减单元.Bayoumy等人^{^[9]}将Shabana等人提出的全参数单元的中间坐标转换矩阵,用于风电叶片中大量斜率不连续结构建模中.
本文针对文献中斜率不连续问题的解决办法进行讨论,为了简化,本文只分析了绝对节点法中经典缩减梁单元建模下的折梁和直梁算例.
1 斜率不连续问题描述
折梁结构中斜率不连续问题是多篇文献中研究的经典算例^{^[7,8]}.在这些文献中,折梁采用三维梁单元或者缩减梁单元进行建模.为了简化起见,本节将以绝对节点法二维缩减梁单元为对象,研究折梁结构中斜率不连续问题,本文还给出了截面呈阶梯形变化的直梁结构中存在的斜率不连续算例.
1.1 绝对节点法缩减梁单元
对于绝对节点法缩减梁单元^{^[13]},采用平截面假定,不考虑剪切变形,采用经典欧拉⁃伯努利梁单元的形函数,其在长度方向是三次函数,高度方向的位移可用中性层内的位移线性表达.梁单元节点自由度,与结构力学中梁单元自由度不同,没有转角自由度,而是采用斜率自由度表示单元的转动.图1 给出了梁单元 $j$ ,该梁单元上任意一点 $P$ 的全局坐标系位置矢量为:

r_{P} = r_{x} i + r_{y} j

(1)

采用梁单元的形函数 $S^{j}$ 及单元自由度 $e^{j}$ ,位置矢量可以表示为:

r_{P} = S^{j} (x_{P}) e^{j}

(2)

其中,x_P 表示点P在梁单元轴线上到原点的距离, 记梁单元的左右节点为A､B, $e^{j}$ 为该单元的单元自由度,表示为:

e^{j} = {[\begin{matrix} e_{A}^{j} & e_{B}^{j} \end{matrix}]}^{T}

(3)

其中, $e_{A}^{j}$ 和 $e_{B}^{j}$ 分别表示第 $j$ 号梁单元的左右节点A､B的节点自由度,可以进一步写成:

e_{N}^{j} = {[r_{N}^{j} \frac{\partial r_{N}^{j}}{\partial x_{j}}]}^{T} N \in {A, B}

(4)

其中, $r_{N}^{j}$ 和 $\frac{\partial r_{N}^{j}}{\partial x_{j}}$ 分别为第 $j$ 号梁单元端点A或B的节点位移自由度及节点斜率自由度.
图1 编号为 $j$ 的平面缩减梁单元
Fig.1 Plane reduced beam element numbered $j$
由图1 可知,梁轴线AB在端点B的方向可以角度 $α_{B}^{j}$ 定义,其方向余弦为 ${[c o s α_{B}^{j}, s i n α_{B}^{j}]}^{T}$ ,而该点的斜率自由度为:

\frac{\partial r_{B}^{j}}{\partial x_{j}} = f_{B}^{j} {[c o s α_{B}^{j}, s i n α_{B}^{j}]}^{T}

(5)

其中, $f_{B}^{j}$ =| $\frac{\partial r_{B}^{j}}{\partial x_{j}}$ |, $f_{B}^{j}$ 与B点截面平均轴向应变 $ε_{B}^{j}$ 之间的关系如下:

ε_{B}^{j} = \frac{1}{2} (\frac{\partial r_{B}^{j}}{\partial x_{j}} \cdot \frac{\partial r_{B}^{j}}{\partial x_{j}} - 1) = \frac{1}{2} ({(f_{B}^{j})}^{2} - 1)

(6)

1.2 折梁斜率不连续问题及非线性约束方程
折梁是指方向不一样的两根梁刚接在一起的结构.简化模型如下图2 所示.
图2 折梁结构的简化模型
Fig.2 Simplified model of the folded beam
将二个梁单元分别编号为 $j$ ､k,x_j 和x_k 分别表示两个梁单元的坐标轴,其夹角为 θ.
第 $j$ 号单元的左右节点编号为A､B,第k号单元的左右节点编号为B､C.由公式(4) 知道,在B点,第 $j$ 号和k号梁的B节点位移自由度满足下列关系:

r_{B}^{j} = r_{B}^{k}

(7)

但是,一般地说,

\frac{\partial r_{B}^{j}}{\partial x_{j}} \neq \frac{\partial r_{B}^{k}}{\partial x_{k}}

(8)

这就是说B节点处两侧的梁单元的斜率自由度不相等,折梁结构为斜率不连续结构.下面给出折梁结构斜率不连续带来的非线性约束方程.
由于两个梁单元的夹角为 θ,则

α_{B}^{k} - α_{B}^{j} = θ

(9)

利用式(5),第 $j$ 号和第k号梁单元的方向余弦可以表示为

{[c o s α_{B}^{i}, s i n α_{B}^{i}]}^{T} = \frac{1}{f_{B}^{i}} \frac{\partial r_{B}^{i}}{\partial x_{i}} i \in {j, k}

(10)

可得B节点处两根梁的斜率自由度之间的非线性关系:

\frac{\partial r_{B}^{k}}{\partial x_{k}} = \frac{f_{B}^{k}}{f_{B}^{j}} A_{j}^{k} \frac{\partial r_{B}^{j}}{\partial x_{j}}

(11)

其中, $A_{j}^{k}$ 是坐标变换矩阵.

A_{j}^{k} = [\begin{matrix} c o s θ & - s i n θ \\ s i n θ & c o s θ \end{matrix}]

(12)

公式(11)可以进一步化简为两根杆的斜率自由度满足的非线性的约束方程:

C_{B} (q) = \frac{\partial r_{B}^{k}}{\partial x_{k}} \times A_{j}^{k} \frac{\partial r_{B}^{j}}{\partial x_{j}} = 0

(13)

因此,此时要求解的多体动力学方程为

{\begin{matrix} M \ddot{q} + C_{, q}^{T} λ = F (q) \\ C (q) = 0 \end{matrix}

(14)

其中,M为结构的常数质量矩阵,λ 为拉格朗日乘子,F(q)为非线性的弹性力向量,C(q)=0 为约束函数,在斜率不连续结构中,该约束函数为非线性函数.公式(14)中第一个方程是斜率不连续结构的多体动力学方程,公式(14) 中第二个方程是多体动力学方程的非线性约束,求解带非线性约束的多体动力学方程非常复杂,而且斜率不连续结构中非线性的约束定义也非常复杂.当系统中存在许多斜率不连续的结构时,非线性约束方程数目增加,算法的效率和收敛速度会降低.
1.3 直梁斜率不连续问题及非线性约束方程
下面我们讨论方向一样但是面积不一样或者弹性模量不一样的两根梁,刚接在一起形成的直梁.简化模型如图3 所示,采用如图1 相同的节点和单元编号.
图3 由2 根短梁组成的直梁结构的简化模型
Fig.3 Simplified model of the straight beam structures consisting of two short beams
由于x_j 与x_k 的方向一致,θ=0,由公式(9)得

α_{B}^{j} = α_{B}^{k}

(15)

但是,由于在B点两个梁的面积不一样或者弹性模量不一样,一般地说它们的应力不同,所以,B点两个梁的轴向应变不一样:

ε_{B}^{j} \neq ε_{B}^{k}

(16)

由公式(6),可以得到:

f_{B}^{j} \neq f_{B}^{k}

(17)

在B节点两根梁的斜率自由度不相等:

\frac{\partial r_{B}^{j}}{\partial x_{j}} \neq \frac{\partial r_{B}^{k}}{\partial x_{k}}

(18)

与折梁斜率不连续问题的非线性约束的推导过程类似,直梁斜率不连续问题的非线性约束表示如(13)所示,由于公式(15)成立, $θ = α_{B}^{k} - α_{B}^{j} = 0, A_{j}^{k}$ 变成单位矩阵.此时,直梁B点节点处的非线性约束为

C_{B} (q) = \frac{\partial r_{B}^{k}}{\partial x_{k}} \times \frac{\partial r_{B}^{j}}{\partial x_{j}} = 0

(19)

2 文献中的处理办法
针对绝对节点坐标法在斜率不连续结构的建模过程中带来的建模和数值求解困难,本文对文献 ^[6]提出的转换坐标系的方法进行研究和讨论
2.1 新坐标系的选择
文献^[6] 提出的坐标转换方法不仅应用在由全参数单元组成的结构,也应用在缩减单元组成的结构.以三维多体系统为例,选取系统中第i个构件上任意一个单元,单元编号为 $j$ ,对单元进行坐标转换时,需要用到4 套坐标系,依次介绍如下.
图4 中的OXYZ坐标系是系统坐标系,也称全局坐标系, $O_{i} X_{i} Y_{i} Z_{i}$ 是第 $i$ 个体的参考坐标系,称为体坐标系,该体坐标系与全局坐标系不一定平行, $O_{i}^{j} x_{i}^{j} y_{i}^{j} z_{i}^{j}$ 是第 $i$ 个体第 $j$ 号单元的参考坐标系.Shabana等人提出中间单元坐标系 $O_{i}^{j} x_{i}^{j} y_{i}^{j} z_{i}^{j}$ 也是第 $i$ 个体上第 $j$ 号单元的单元坐标系,初始时刻该坐标系坐标轴方向与单元坐标系 $O_{i}^{j} x_{i}^{j} y_{i}^{j} z_{i}^{j}$ 重合, 在运动过程中,该坐标系与体坐标系 $O_{i} X_{i} Y_{i} Z_{i}$ 方向总是固定在一起.
图4 三维单元的4 套坐标系示意图
Fig.4 Four sets of coordinate systems for three dimensional elements
2.2 新坐标系的应用
对斜率不连续结构,Shabana通过常数的坐标转换,将相连结杆件的斜率自由度转换到中间单元坐标系中的斜率自由度,转换关系如下:

\frac{\partial r_{B}^{j}}{\partial x_{j}} = A^{j} g_{B}, \frac{\partial r_{B}^{k}}{\partial x_{k}} = A^{k} g_{B}

(20)

其中, $g_{B}$ 为中间单元坐标系中的斜率自由度, $A^{j}$ 和 $A^{k}$ 为第 $j$ 号和第 $k$ 号单元的单元坐标系相对体坐标系的旋转矩阵,对于二维问题,文献^[8] 给出的表达式如下:

A^{i} = [\begin{matrix} c o s α_{B}^{i} & - s i n α_{B}^{i} \\ s i n α_{B}^{i} & c o s α_{B}^{i} \end{matrix}], i \in {j, k}

(21)

由于新单元坐标系与中间单元坐标系是固定在一起的, $A^{j}$ ､ $A^{k}$ 是常数矩阵.通过常数的中间坐标转换矩阵,能够保证绝对节点坐标法下常数质量矩阵的主要优势.转换后的第 $j$ 号和第 $k$ 号梁单元B节点自由度变成:

{\tilde{e}}_{B}^{j} = {[\begin{matrix} r_{B}^{j} & g_{B} \end{matrix}]}^{T}, {\tilde{e}}_{B}^{k} = {[\begin{matrix} r_{B}^{k} & g_{B} \end{matrix}]}^{T}

(22)

其中, ${\tilde{e}}_{B}^{j}$ 和 ${\tilde{e}}_{B}^{k}$ 是转换坐标系后的节点自由度,它们和 $e_{B}^{j}$ ､ $e_{B}^{j}$ 满足以下的关系:

{\tilde{e}}_{B}^{j} = T^{j} e_{B}^{j}, {\tilde{e}}_{B}^{k} = T^{k} e_{B}^{j}

(23)

T^{j} = [\begin{matrix} I & 0 \\ 0 & A^{j} \end{matrix}], T^{k} = [\begin{matrix} I & 0 \\ θ & A^{k} \end{matrix}]

(24)

并且,转换后的第 $j$ 号和第 $k$ 号梁单元B节点自由度满足下列关系:

r_{B}^{j} = r_{B}^{k}, {\tilde{e}}_{B}^{j} = {\tilde{e}}_{B}^{k}

(25)

通过上述的坐标转换,如公式(24) 所示,第 $j$ 号和第 $k$ 号梁单元B节点斜率自由度变得连续,这样能够带来的好处是,能够消除多体动力学方程(14)中的约束函数,即多体动力学方程变成:

\tilde{M} \tilde{\ddot{q}} = \tilde{F} (\tilde{q})

(26)

其中, $\tilde{M}$ 是转换坐标系后的常数质量矩阵, $\tilde{q}$ 是转换坐标系后的结构自由度, $\tilde{F}$ 是转换坐标系后的非线性的弹性力向量,关于转换后的质量矩阵请见文献^[8].文献^[6] 等认为新坐标系的引入有很大的优势,一方面是能够保证坐标转换后的结构质量矩阵仍然是常数矩阵,另外一方面是能够消去多体动力学方程中的约束函数,这为动力学方程的建模和求解带来很大的便利
2.3 新坐标系带来的影响
尽管新坐标系带来了动力学方程的建模和求解的便利,但是新坐标系的引入带来的影响需要进一步的讨论.
由公式(21)得 ${(A^{j})}^{T} A^{j} = {(A^{k})}^{T} A^{k} = I$ ,进而,

\begin{matrix} {(\frac{\partial r_{B}^{j}}{\partial x_{j}})}^{T} \frac{\partial r_{B}^{j}}{\partial x_{j}} = g_{B}^{T} {(A^{j})}^{T} A^{j} g_{B} = g_{B}^{T} g_{B} \\ {(\frac{\partial r_{B}^{k}}{\partial x_{k}})}^{T} \frac{\partial r_{B}^{k}}{\partial x_{k}} = g_{B}^{T} {(A^{k})}^{T} A^{k} g_{B} = g_{B}^{T} g_{B} \end{matrix}

(27)

又由公式(7)得到:

{(\frac{\partial r_{B}^{j}}{\partial x_{j}})}^{T} \frac{\partial r_{B}^{j}}{\partial x_{j}} = {(f_{B}^{j})}^{2}, {(\frac{\partial r_{B}^{k}}{\partial x_{k}})}^{T} \frac{\partial r_{B}^{k}}{\partial x_{k}} = {(f_{B}^{k})}^{2}

(28)

观察(27)和(28),发现必须满足:

f_{B}^{j} = f_{B}^{k}

(29)

所以,文献^[6]的方法实际上假设了在B点,第 $j$ 号和第 $k$ 号梁单元的轴向应变相等.但是,在斜率不连续的结构中,B点处不同梁单元上的轴向应变不一定相等.下面通过算例来观察文献^[6]中的方法对应变计算带来的影响.
3 算例
本文取了两种斜率不连续的代表性算例,分别是两根梁轴线一致,但面积不相同的直梁算例和轴线方向不相同的折梁算例.为了验证文献^[6]的方法中假设,我们对两种方法计算得到的轴向应变进行比较.第一个方法是采用带约束函数的多体动力学方程,见公式(14),计算斜率不连续处的应变.第二个方法是使用文献^[6]的方法,见公式(26),计算斜率不连续处的应变,比较二者应变.
3.1 算例一
考虑如下图5 所示的直梁,从水平位置开始自由下落,A端铰接,B点是刚接点,C端自由,取重力加速度为9.81m/s ²,其他参数见表1.
表1 直梁参数
Table 1 Data of the straight beam
图5 面积不一样的直梁简化模型
Fig.5 Simplified model of the straight beam with different area
每根梁划分4､8､32 个单元,时间积分步长为0.01s,积分总时长为2.5s.采用2 种方法算出在0.6s时,截面面积不连续B点处左右单元的中性层的轴向应变.图6(a)､( b)分别采用文献中的坐标转换法及约束函数法.随着单元数量的增加,这两个图中B点处中性层的轴向应变结果均趋于收敛.但是,图6(a)说明,采用文献中的坐标转换法,要精确算出应变分布结果需划分大量单元.图6( b) 采用约束函数法,即使采用4 个单元,也能够准确地描述轴向应变的跳跃.
为了进一步比较两种方法,取出每根梁划分32 个梁单元的应变结果,第二根梁考虑了三个不同的高度,其结果如图7 所示,图中放大斜率不连续处应变,缩小比较区间,只关注梁长度方向0.8~1m区间内中性层轴向应变变化.
图6 梁长度方向中性层上的轴向应变分布
Fig.6 Axial strain distribution on the neutral layer along the beam length direction
图7 第二根梁的高度变化时梁长度方向0.8~1m中性层上的轴向应变分布情况
Fig.7 Axial strain distribution on the 0.8~1m neutral layer along the beam length direction when the height of beam changes
在斜率不连续处,因为采用了32 个单元剖分梁,每个单元的长度很小,单元上分布的重力和惯性力转换到B点上的节点力很小,B点两侧梁的内力近似相等,弹性模量相等时,二根梁的轴向应变比应该近似反比于其面积比.对于两根杆件的面积比分别是2､4､8 的情况,使用带约束函数的方法时,算出B点左右截面的轴向应变比近似为1 ∶2､ 1 ∶4､1 ∶8,见表2,但是使用文献^[6]的方法时,这些轴向应变比始终是1 ∶1,见表3.
所以,在该直梁算例中,以上的结果证明了文献^[6]的方法不能准确给出斜率不连续处的相邻杆件的轴向应变.
表2 带约束方法下B点左右截面的轴向应变
Table 2 Axial strain of left and right section in B point under constraint method
表3 文献^[6]方法下B点左右截面的轴向应变
Table 3 Axial strain of left and right section in B point under literature ^[6] method
3.2 算例二
考虑如图8 所示的折梁,从水平位置开始自由下落,A端铰接,B点是刚接点,C端自由,取重力加速度为9.81m/s ²,其他参数见表4.
图8 方向不一样的折梁的简化模型
Fig.8 Simplified model of the folded beam with different orientation
表4 折梁参数
Table 4 Data of the folded beam
在有限元计算中,每根梁划分8 个单元,整个梁结构划分16 个梁单元,时间积分步长为0.01s, 积分总时长为0.6s.采用2 种方法算出在0.1s时, 斜率不连续处左右单元的B点处中性层的轴向应变.
在B点,因为两根梁的方向不一样,所以,第二根梁的轴向应变与第一根梁的轴向应变一般不相等.为了比较两种方法的应变,放大该B点附近的应变,缩小比较的区间范围,图9 中,横坐标是沿着梁的轴线的距离(在0.25m处梁的轴线发生了方向变化),给出了第一根梁0.2~0.25m处及第二根梁0~0.05m处中性层上的轴向应变分布情况.
图9 AB梁长度方向0.2~0.25m和BC梁长度方向0~0.05m中性层上的轴向应变分布情况
Fig.9 Axial strain distribution on the 0.2~0.25m neutral layer along AB length direction and on the 0~0.05m neutral layer along BC length direction
表5 两种方法下B点左右截面的轴向应变
Table 5 Axial strain of left and right section in B point under two methods
由图9 和表5 得到,使用带约束函数的方法时,算出B点左右截面的轴向应变不相等,但是使用文献^[6]的方法时,算出B点左右截面的轴向应变相等.所以,在折梁算例中,以上的结果也证明了文献^[6]的方法确实有斜率不连续处轴向应变相等的假设.
4 结论
本文研究了文献^[6] 中斜率不连续问题的处理方法.通过比较算例中斜率不连续处左右截面的轴向应变,发现文献^[6]的方法不能准确给出在斜率不连续处连结的杆件的轴向应变.值得注意的是,在采用绝对节点坐标法研究旋转梁的刚化效应时^{^[11]},由于绝对节点坐标法缺乏转角自由度,此文献也通过坐标转换的方式来实现转角边界条件的施加,但是,此文坐标转换矩阵是时变矩阵,转换后的质量矩阵非常矩阵,出现了科氏和离心惯性项, 影响了绝对节点坐标法的优势.这样的问题也存在于梁板壳及其组合结构,也存在于将绝对节点坐标法的单元和其他方法单元如何连接的问题中^{^[12]}.因此,绝对节点坐标法下斜率不连续问题的解决办法需要进一步研究,例如研究对于轴向应变足够小的问题,文献中的坐标转换法是否可以作为一种近似方法使用.
参考文献
- [1] Shabana A A.Dynamics of multibody systems.New York:Cambridge University Press,2013
- [2] 田强,刘铖,李培,等.多柔体系统动力学研究进展与挑战.动力学与控制学报,2017,15(5):385 ~ 405(Tian Q,Liu C,Li P,et al.Advances and challenges in dynamics of flexible multibody systems.Journal of Dynam⁃ ics and Control,2017,15(5):385~ 405(in Chinese))
- [3] Shabana A A.An absolute nodal coordinate formulation for the large rotation and deformation analysis of flexible bodies.Technical Report.no.MBS96⁃1⁃UIC,1996
- [4] Gerstmayr J,Sugiyama H,Mikkola A.Review on the ab⁃ solute nodal coordinate formulation for large deformation a⁃nalysis of multibody systems.Journal of Computational & Nonlinear Dynamics,2013,8(3):031016
- [5] Shabana A A,Mikkola A M.Use of the finite element ab⁃ solute nodal coordinate formulation in modeling slope dis⁃ continuities.Multibody System Dynamics,2003,125(2):357~ 387
- [6] Shabana A A,Maqueda L G.Slope discontinuities in the finite element absolute nodal coordinate formulation:gra⁃ dient deficient elements.Multibody System Dynamics,2008,20(3):239~ 249
- [7] Maqueda L G,Shabana A A.Numerical investigation of the slope discontinuities in large deformation finite element formulations.Nonlinear Dynamics,2009,58(1 ~ 2):23 ~ 37
- [8] Shabana A A.General method for modeling slope disconti⁃ nuities and T⁃sections using ANCF gradient deficient finite elements.Journal of Computational & Nonlinear Dynam⁃ ics,2011,6(2):267~ 274
- [9] Bayoumy A,Nada A,Megahed S.Methods of modeling slope discontinuities in large size wind turbine blades u⁃ sing absolute nodal coordinate formulation.Proceedings of the Institution of Mechanical Engineers Part K Journal of Multi⁃body Dynamics,2014,228(3):314~ 329
- [10] Romero I,Arribas J J.A simple method to impose rota⁃ tions and concentrated moments on ANC beams.Multibody System Dynamics,2009,21(4):307~ 323
- [11] Berzeri M,Shabana A A.Study of the centrifugal stiffe⁃ ning effect using the finite element absolute nodal coordi⁃ nate formulation.Multibody System Dynamics,2002,7(4):357~ 387
- [12] Yu L,Zhao Z,Tang J,et al.Integration of absolute nodal elements into multibody system.Nonlinear Dynamics,2010,62(4):931~ 943
- [13] 章孝顺,章定国,陈思佳,等.基于绝对节点坐标法的大变形柔性梁几种动力学模型研究.物理学报,2016,65(9):148~ 157(Zhang X S,Zhang D G,Chen S J,et al.Several dynamic models of a large deformation flexible beam based on the absolute nodal coordinate formulation.Acta Physica Sinica,2016,65(9):148 ~ 157(in Chi⁃ nese))

基本信息

DOI: 10.6052/1672⁃6553⁃2019⁃060

基金信息

国家自然科学基金资助项目(13332004)；
The project supported by the National Natural Science Foundation of China(13332004)；

引用信息

稿件历史

参考文献

[1] Shabana A A.Dynamics of multibody systems.New York:Cambridge University Press,2013
[2] 田强,刘铖,李培,等.多柔体系统动力学研究进展与挑战.动力学与控制学报,2017,15(5):385 ~ 405(Tian Q,Liu C,Li P,et al.Advances and challenges in dynamics of flexible multibody systems.Journal of Dynam⁃ ics and Control,2017,15(5):385~ 405(in Chinese))
[3] Shabana A A.An absolute nodal coordinate formulation for the large rotation and deformation analysis of flexible bodies.Technical Report.no.MBS96⁃1⁃UIC,1996
[4] Gerstmayr J,Sugiyama H,Mikkola A.Review on the ab⁃ solute nodal coordinate formulation for large deformation a⁃nalysis of multibody systems.Journal of Computational & Nonlinear Dynamics,2013,8(3):031016
[5] Shabana A A,Mikkola A M.Use of the finite element ab⁃ solute nodal coordinate formulation in modeling slope dis⁃ continuities.Multibody System Dynamics,2003,125(2):357~ 387
[6] Shabana A A,Maqueda L G.Slope discontinuities in the finite element absolute nodal coordinate formulation:gra⁃ dient deficient elements.Multibody System Dynamics,2008,20(3):239~ 249
[7] Maqueda L G,Shabana A A.Numerical investigation of the slope discontinuities in large deformation finite element formulations.Nonlinear Dynamics,2009,58(1 ~ 2):23 ~ 37
[8] Shabana A A.General method for modeling slope disconti⁃ nuities and T⁃sections using ANCF gradient deficient finite elements.Journal of Computational & Nonlinear Dynam⁃ ics,2011,6(2):267~ 274
[9] Bayoumy A,Nada A,Megahed S.Methods of modeling slope discontinuities in large size wind turbine blades u⁃ sing absolute nodal coordinate formulation.Proceedings of the Institution of Mechanical Engineers Part K Journal of Multi⁃body Dynamics,2014,228(3):314~ 329
[10] Romero I,Arribas J J.A simple method to impose rota⁃ tions and concentrated moments on ANC beams.Multibody System Dynamics,2009,21(4):307~ 323
[11] Berzeri M,Shabana A A.Study of the centrifugal stiffe⁃ ning effect using the finite element absolute nodal coordi⁃ nate formulation.Multibody System Dynamics,2002,7(4):357~ 387
[12] Yu L,Zhao Z,Tang J,et al.Integration of absolute nodal elements into multibody system.Nonlinear Dynamics,2010,62(4):931~ 943
[13] 章孝顺,章定国,陈思佳,等.基于绝对节点坐标法的大变形柔性梁几种动力学模型研究.物理学报,2016,65(9):148~ 157(Zhang X S,Zhang D G,Chen S J,et al.Several dynamic models of a large deformation flexible beam based on the absolute nodal coordinate formulation.Acta Physica Sinica,2016,65(9):148 ~ 157(in Chi⁃ nese))

分享给微信好友或者朋友圈

使用微信“扫一扫”功能。

绝对节点坐标法下斜率不连续问题处理方法讨论*

DISCUSSION ON SLOPE DISCONTINUITY USING THE ABSOLUTE NODAL COORDINATE FORMULATION ∗

摘要

Abstract

关键词

Keywords

0 引言

1 斜率不连续问题描述

1.1 绝对节点法缩减梁单元

(1)

(2)

(3)

(4)

(5)

(6)

1.2 折梁斜率不连续问题及非线性约束方程

(7)

(8)

(9)

(10)

(11)

(12)

(13)

(14)

1.3 直梁斜率不连续问题及非线性约束方程

(15)

(16)

(17)

(18)

(19)

2 文献中的处理办法

2.1 新坐标系的选择

2.2 新坐标系的应用

(20)

(21)

(22)

(23)

(24)

(25)

(26)

2.3 新坐标系带来的影响

(27)

(28)

(29)

3 算例

3.1 算例一

3.2 算例二

4 结论

参考文献

基本信息

基金信息

引用信息

稿件历史

参考文献

微信公众号二维码

手机版网站二维码

绝对节点坐标法下斜率不连续问题处理方法讨论^*

DISCUSSION ON SLOPE DISCONTINUITY USING THE ABSOLUTE NODAL COORDINATE FORMULATION ^∗