计算力学课堂教学课件第2章-3.ppt
计算力学课堂教学课计算力学课堂教学课件第件第2章章-32.3 2.3 广义坐标有限元法的一般格式广义坐标有限元法的一般格式常见的单元类型：常见的单元类型：2.3.1 2.3.1 选择单元位移函数的一般原则选择单元位移函数的一般原则（1 1）位移模式中的待定系数（广义坐标）个数，）位移模式中的待定系数（广义坐标）个数，应与单元的结点位移数相等；应与单元的结点位移数相等；（2.2.12.2.1）（2 2）位移模式中，常数项）位移模式中，常数项与一次项必须完备；与一次项必须完备；（3 3）多项式选取由低次到高次，且尽可能选取）多项式选取由低次到高次，且尽可能选取完全的多项式，以提高单元的精度；如：完全的多项式，以提高单元的精度；如：一次一次二次二次三次三次四次四次0 次次Pascal 三角形三角形（4 4）当由于项数限制，不能选取取完全多项式时，应考当由于项数限制，不能选取取完全多项式时，应考虑多项式应具有对称性，如：虑多项式应具有对称性，如：4节点节点矩形单元矩形单元2.3.2 2.3.2 广义坐标有限单元的一般格式广义坐标有限单元的一般格式（1 1）以广义坐标以广义坐标 为待定参数，给出单元内位移为待定参数，给出单元内位移分布分布 u u；（2.3.12.3.1）对于二维问题：对于二维问题：对于对于3 3结点三角形单元：结点三角形单元：3 3节点节点三角形单元三角形单元1 12 23 3n n2用单元结点位移 表示广义坐标 惯用的单元结点位移排列是 为便于求解广义坐标，可采用另一表示方法，如（231）式中代入单元结点坐标得到 二维问题 用（232）式解出 n n3以单元结点位移ae 表示单元位移函数u，得到单元插值函数矩阵N N将（233）式代入（23l）式 二维问题 将结点位移 改为一般排列顺序ae，则有 n n4以单元结点位移ae 表示单元应变和应力应变：应力：由弹性变形产生的应力 当有初应力和初应变时，应力的一般式是 n n5用最小位能原理建立离散体系的结点平衡方程n n系统总位能的离散形式 将（235）（237）式代入上式并将单元结点位移ae用结构结点位移a表示，即ae=Ga，（238）式即为 n n总位能的变分 得到有限元求解方程 其中 式（2312）是单元刚度矩阵的普遍公式（239）式中的是作用在连续体边界上的力，包括作用在有关单元边界上的分布力和作用在结点上的集中力两部分。为了方便起见可以将这两部分外载分开，将T T作为分布面力，结点集中力用F F表示，则载荷列阵P P可以写作P PF F 结点集中力列阵（2314）式是计算单元等效结点载荷列阵的普遍公式。n n6引入强制边界条件n n7解方程得到结点位移n n8进行需要的辅助计算n n如利用（236）、（237）式计算单元应变和应力，也可按需要计算其他项目。n n由上面过程可以看到：n n13是假定位移模式、求解广义坐标，最后得到单元插值函数。这三步是广义坐标有限元的特征。n n45是利用变分原理建立有限元格式的一般方法。这里用的是最小位能原理，建立以位移为基本场变量的有限元求解方程，求解平衡问题。n n68是建立有限元方程后的一般解法和计算步骤。n n广义坐标有限元可能产生的困难是：n n当位移函数选择不恰当时，可能不存在A-1而使求解广义坐标 成为不可能。n n同时，当单元结点较多时，解广义坐标的过程显得繁琐，因此也可以利用自然坐标直接构造单元的插值函数，这样就可以避免求解广义坐标的过程，建立有限元的方程和求解只需从第4步开始。n n本章第5节将结合矩形单元和高精度三角形单元讨论直接建立单元插值函数的方法，关于建立单元插值函数更系统方法将在下一章中给出。2.4 2.4 有限单元解的性质与收敛性有限单元解的性质与收敛性一.Ritz法的收敛准则对连续介质问题，有泛函：其中:要求试函数必须满足：完备的函数系列（完备性）应满足连续性要求（协调性）1）2）Ritz 法的收敛条件：法的收敛条件：（1）近似函数）近似函数 u 具有具有完备性完备性（2）试探函数）试探函数 u 具有具有连续性连续性（取完全多项式）（取完全多项式）（C0 类连续）类连续）FEM 法与法与Ritz 法的区别：法的区别：Ritz 法在法在全域全域上假设近似函数上假设近似函数 u，FEM 法在法在单元单元上假设近似函数上假设近似函数 u，近似函数近似函数 u可以有多种类型。可以有多种类型。近似函数近似函数 u一般都为一般都为简单多项式简单多项式。问题：问题：在什么条件下，当单元尺寸趋于零时，在什么条件下，当单元尺寸趋于零时，FEM 法的解趋于精确解？法的解趋于精确解？引例：引例：设一标量场：设一标量场：存在标量泛函：存在标量泛函：假设泛函假设泛函 中含有：中含有：为非零的，为非零的，则近似函数则近似函数 至少为至少为 m 次多项式，即次多项式，即（2.4.3）显然，仅当显然，仅当 pm 时，时，各项都包含常数项，各项都包含常数项，意味着，当单元尺寸趋意味着，当单元尺寸趋各项趋于常数，即，各项趋于常数，即，有限元法的近似解有限元法的近似解收敛于收敛于精确解。精确解。于零时，于零时，收敛条件收敛条件1：必须是必须是m次以上的完备多项式；次以上的完备多项式；收敛条件收敛条件2：仅可能在相邻单元边界上连续。仅可能在相邻单元边界上连续。弹性力学问题中，有限元法的收敛准则：弹性力学问题中，有限元法的收敛准则：准则准则1：完备性要求完备性要求 若泛函若泛函p中含有未知函数的最高阶导数为中含有未知函数的最高阶导数为 m 阶，则未知函数至少是阶，则未知函数至少是 m 次次完全多项式。完全多项式。p中含：中含：要求：要求：为为 x、y、z 一次完全多项式一次完全多项式对于平面问题：对于平面问题：对于空间问题：对于空间问题：对于梁的弯曲问题？对于梁的弯曲问题？最高阶导数次数最高阶导数次数 m=1弹性力学问题：弹性力学问题：准则准则2：协调性要求协调性要求 若泛函若泛函p中含有未知函数的最高阶导数中含有未知函数的最高阶导数为为 m 阶，则要求未知函数在相邻单元交界面阶，则要求未知函数在相邻单元交界面上须满足上须满足 Cm-1 类连续性类连续性，即保证交界面上未即保证交界面上未知函数知函数 m1 阶的连续可导。阶的连续可导。当单元的插值函数满足上述要求时，称这当单元的插值函数满足上述要求时，称这种种单元是协调的单元是协调的。对对弹性力学问题弹性力学问题，泛函泛函p中含有未知函数中含有未知函数故要求近似位移函数故要求近似位移函数 u对于对于3结点三角形单元：结点三角形单元：满足协调条件满足协调条件的最高阶导数为的最高阶导数为 1 阶，阶，在相邻单元交界面上须满足在相邻单元交界面上须满足 C0 类连续性类连续性。总结：总结：同时满足同时满足完备性完备性、协调性协调性条件的单元条件的单元称为称为协调单元协调单元（3 3结点三角形单元结点三角形单元为完备、协调单元）为完备、协调单元）完备、协调单元的解一定收敛于精确解。完备、协调单元的解一定收敛于精确解。不满足不满足协调性协调性条件的单元条件的单元称为称为非非协调单元协调单元如：板壳问题中某些单元。如：板壳问题中某些单元。2.4.2 2.4.2 收敛准则的物理意义收敛准则的物理意义1.1.完备性准则完备性准则要求位移函数为一次完全多项式要求位移函数为一次完全多项式常数项：常数项：一次项：一次项：反映刚体位移反映刚体位移 反映常应变反映常应变2.2.协调性准则协调性准则 要求位移函数在相邻单元边界上连续，要求位移函数在相邻单元边界上连续，避免在单元交界面产生无限大的应变能。避免在单元交界面产生无限大的应变能。2.4.3 2.4.3 收敛速度与精度估计收敛速度与精度估计三.离散误差估计精确解可以在域内一点i 的某一邻域内泰勒展开：如果：单元(特征尺寸h)内位移插值函数是p次完全多项式,则可局部拟合泰勒展开p阶.x,y为h 量级，位移的误差约为应力和应变（为位移的一阶导数）则误差约为，精度低于位移。例：3结点三角形单元位移误差应力误差（应变）提高精度的方法：(1)单元尺寸变小,(2)插值函数，完备的多项式次数提高。由其他误差:计算误差，包括截断误差，舍入误差.提高精度的方法:(1)增长字长(双精度)(2)选取有效的计算方法和合理的程序结构。2.4.4 2.4.4 位移解的下限性质位移解的下限性质位移有限元法位移有限元法基于最小位能原理：基于最小位能原理：由第由第1 1章讨论，可知：章讨论，可知：变形位能（应变能）变形位能（应变能）设有限元解的总位能、应变能、刚度矩阵、设有限元解的总位能、应变能、刚度矩阵、结点位移向量为：结点位移向量为：由最小位能原理，有由最小位能原理，有即：即：由由即，精确解取近似解的下限。即，精确解取近似解的下限。位移解下限性质解释：位移解下限性质解释：单元原是连续体的一部分，具有无限个单元原是连续体的一部分，具有无限个自由度。假定单元的位移函数后，其自度自由度。假定单元的位移函数后，其自度数限制为只有以结点位移表示的有限个自数限制为只有以结点位移表示的有限个自由度，相当于由度，相当于对单元的变形施加了约束和对单元的变形施加了约束和限制限制，使，使单元的刚度较实际情况加强了单元的刚度较实际情况加强了，随之，整个连续体的刚度增加了，故求得随之，整个连续体的刚度增加了，故求得的位移总体上小精确解。的位移总体上小精确解。2.5 2.5 矩形单元和高精度三角形单矩形单元和高精度三角形单元元2.5.1 4 2.5.1 4 结点矩形单元结点矩形单元 2314aabb设:边长分别为2a 和2b。取坐标原点在单元形心上，x、y轴分别平行两对对边。2.5.1 4 2.5.1 4 结点矩形单元结点矩形单元 1.1.结点编号与单元结点位移向量结点编号与单元结点位移向量结点编号：结点编号：逆时针转向为正逆时针转向为正单元结点位移向量：单元结点位移向量：8 8个自由度个自由度2314aabb（2.5.22.5.2）2.2.单元位移模式单元位移模式 双线性位移模式双线性位移模式将结点的将结点的 x 方向位移与坐标代入方向位移与坐标代入一次一次二次二次三次三次四次四次0 次次Pascal 三角形三角形位移多项式：位移模式是具有完全一次式，非完全二次式。因为它在x,y方向呈线形变化，所以称之为双线形位移模式。2314aabb由此可解出：由此可解出：其中：其中：单元插值单元插值函数函数或或形函数形函数（2.5.42.5.4）2314aabb令：令：则有则有（2.5.5）局部坐标局部坐标表示的表示的单元插值函数单元插值函数或或形函数形函数 自然坐标自然坐标或或局部坐标局部坐标3.3.自然坐标与插值函数自然坐标与插值函数（1 1）自然（局部）自然（局部）坐标坐标2314aabb由：由：有：有：（1 1）（2 2）矩形单元四边的方程为：矩形单元四边的方程为：（3 3）角点值（坐标）角点值（坐标）2314aabbxy一般情形，将自然坐一般情形，将自然坐标原点取在单元的重心标原点取在单元的重心上，整体坐标与自然坐上，整体坐标与自然坐标的关系为：标的关系为：（2.5.6）（2 2）插值函数（形函）插值函数（形函数）数）Ni 的性质的性质（1 1）（2 2）（3 3）角点值：角点值：或：或：令：令：其中：其中：为它们在角点为它们在角点 i 处的值处的值则插值函数可表示为：则插值函数可表示为：（2.5.72.5.7）2314aabb4.4.单元位移单元位移 ae 的矩阵表示的矩阵表示2314aabb（2.5.8）其中：其中：（2.5.92.5.9）形函数矩阵形函数矩阵2314aabb5.5.单元应变与应力单元应变与应力（1 1）单元应变）单元应变 应变矩阵应变矩阵（2.5.102.5.10）2314aabb由由得得（2.5.12）（2 2）单元应力）单元应力2314aabb（2.5.13）应力矩阵应力矩阵其中：其中：说明：说明：矩形单元的应变、应力关于矩形单元的应变、应力关于 x、y 线线性分布性分布6.6.单元刚度矩阵单元刚度矩阵按结点的分块矩阵形式：按结点的分块矩阵形式：（2.5.162.5.16）其中：其中：（2.5.172.5.17）注意：注意：平面应变问题：平面应变问题：（P67P67）n n对于平面应力问题，应变矩阵、应力矩阵和单元刚度矩阵的显式如（2518）式（2520）。7.7.单元等效结点载荷矩阵单元等效结点载荷矩阵（1 1）体积力引起的：）体积力引起的：体力向量体力向量（2 2）边界面力引起的：）边界面力引起的：面力向量面力向量（3 3）初始应力引起的：）初始应力引起的：初应力向量初应力向量（4 4）初始应变引起的：）初始应变引起的：初应变向量初应变向量8.8.矩形单元插值函数构造矩形单元插值函数构造2314（1 1）插值函数构造的原则：）插值函数构造的原则：或：或：（2 2）（1 1）角点值：角点值：（2 2）插值函数构造的方法）插值函数构造的方法待定系数待定系数将将1 1点的坐标：点的坐标：代入代入由此可解得：由此可解得：代入，得代入，得又如：又如：将将2 2点的坐标：点的坐标：代入代入由此可解得：由此可解得：2314同理：同理：将将3 3点的坐标：点的坐标：代入代入由此可解得：由此可解得：将将4 4点的坐标：点的坐标：代入代入由此可解得：由此可解得：2314综合，得：综合，得：23148567例：例：8 8结点矩形单元插结点矩形单元插值函数的构造。值函数的构造。（1 1）4 4个角点：个角点：将将1 1点的坐标：点的坐标：代入代入由此可解得：由此可解得：同理，可求得：同理，可求得：23148567（2 2）4 4个中点：个中点：将将5 5点的坐标：点的坐标：代入代入由此可解得：由此可解得：23148567将将6 6点的坐标：点的坐标：代入代入由此可解得：由此可解得：2314856723148567 8 8结点矩形结点矩形单元插值函数单元插值函数或形函数或形函数（三次多项式）（三次多项式）23148567单元位移分布函数为：单元位移分布函数为：说明：说明：矩形单元的缺点：矩形单元的缺点：对边界形状的适应差。对边界形状的适应差。矩形单元的优点：矩形单元的优点：（1 1）插值函数（形函数）容易构造；）插值函数（形函数）容易构造；（2 2）单元矩阵）单元矩阵 k ke、P Pe 积分求解方便。积分求解方便。2.5.2 2.5.2 高精度三角形单元高精度三角形单元 6 6节点节点三角形单元三角形单元1 12 23 34 45 56 6Pascal Pascal 三角三角形形1.1.二次单元：二次单元：6 6节点三角形单元节点三角形单元其中：其中：4 4，5 5，6 6节点为三角形节点为三角形边的中点。边的中点。单元位移模式：单元位移模式：应变与应力向量：应变与应力向量：应变、应力应变、应力 线性分布线性分布 协调单元协调单元满足完备性；满足完备性；满足协调性。满足协调性。6 6节点三角形单元的精度较节点三角形单元的精度较3 3结点单元高。结点单元高。1 12 23 34 45 56 67 78 89 910101010节点节点三角形单元三角形单元2.2.三次单元：三次单元：10 10 节点三角形单元节点三角形单元 其中：其中：4-9 4-9 节点为三角形三边的三节点为三角形三边的三分点，分点，1010结点为三角形的中点。结点为三角形的中点。Pascal Pascal 三角三角形形单元位移模式：单元位移模式：应变、应力应变、应力 抛物线分布抛物线分布 协调单元协调单元满足完备性；满足完备性；满足协调性。满足协调性。其精度高于其精度高于6 6节点的三角形单元。节点的三角形单元。2.5.3 2.5.3 面积坐标为自然坐标时三角形单元的面积坐标为自然坐标时三角形单元的插值函数及单元矩阵的计算插值函数及单元矩阵的计算 1.1.面积坐标面积坐标i jm P(x,y)AiAmAj(xi,yi)(xj,yj)(xm,ym)设设 P（x,y）为单元内为单元内任一点，该点与单元三角任一点，该点与单元三角点的连线确定三个三角形，点的连线确定三个三角形，其面积分别用：其面积分别用：Ai、Aj、Am 表示，表示，(1,0,0)(0,1,0)(0,0,1)令：令：（2.5.232.5.23）P点的位置完全由此确定，即点的位置完全由此确定，即 称为称为面积坐标面积坐标其中：其中：A 为三角形单元为三角形单元的面积，有的面积，有显然，有显然，有i jm P(x,y)AiAmAj(xi,yi)(xj,yj)(xm,ym)(1,0,0)(0,1,0)(0,0,1)面积坐标的特点：面积坐标的特点：（1 1）角点值：角点值：（2 2）单元三边的方程：）单元三边的方程：jm边：边：mi边：边：ij边：边：i jm P(x,y)AiAmAj(xi,yi)(xj,yj)(xm,ym)(1,0,0)(0,1,0)(0,0,1)di（3 3）平行于单元边界直线上的面积坐标值：）平行于单元边界直线上的面积坐标值：l如：直线如：直线 l 平行于边平行于边 jm则直线则直线 l 任一点任一点 Li 值，有值，有（4 4）不独立不独立 与单元的与单元的具体形状具体形状、整体、整体坐标坐标x、y无关。无关。（5 5）为三角单元的为三角单元的局部坐标局部坐标i jm P(x,y)AiAmAj(xi,yi)(xj,yj)(xm,ym)(1,0,0)(0,1,0)(0,0,1)i jm P(x,y)AidiAmAj(1,0,0)(0,1,0)(0,0,1)l2.2.面积坐标与直角坐标的关系面积坐标与直角坐标的关系Ai 的值：的值：（2.5.26）（2.5.27）（2.5.28）用整体坐标用整体坐标 x、y 表示面表示面积坐标积坐标 Li、Lj、Lm。用矩阵表示，有用矩阵表示，有（2.5.292.5.29）将将分别乘以分别乘以xi、xj、xm，并相加，有，并相加，有用面积坐标用面积坐标 Li、Lj、Lm表示整体坐标表示整体坐标 x、y。3.3.面积坐标微积分运算面积坐标微积分运算由由有：有：微分运算微分运算（2.5.30）积分运算积分运算（1 1）面积坐标）面积坐标 Li、Lj、Lm的幂函数的幂函数在三角形全面积上的积分：在三角形全面积上的积分：（2.5.312.5.31）（2 2）面积坐标）面积坐标 Li、Lj、Lm的幂函数在三角形某的幂函数在三角形某一边（一边（i-j）上的积分：）上的积分：（2.5.32）如：如：4.4.面积坐标给出的三角形单元的插值函数面积坐标给出的三角形单元的插值函数（1 1）线性单元线性单元3 3节点三角形单元节点三角形单元1 23(1,0,0)(0,1,0)(0,0,1)（2.5.332.5.33）（2 2）二次单元二次单元6 6节点三角形单元节点三角形单元1 23(1,0,0)(0,1,0)(0,0,1)4563 3个角点：个角点：将将1 1节点的面积坐标：节点的面积坐标：L1=1代入代入求得：求得：所求所求 N1 为为同理：同理：1 23(1,0,0)(0,1,0)(0,0,1)4563 3个中点：个中点：将将4 4节点的坐标：节点的坐标：L1=L2=1/2 代入代入求得：求得：所求所求 N4 为为将将4 4节点的坐标：节点的坐标：L2=L3=1/2 代入代入求得：求得：所求所求 N5 为为同理：同理：综合得：综合得：（2.5.35）1 23(1,0,0)(0,1,0)(0,0,1)456（3 3）三次单元三次单元1010节点三角形单元节点三角形单元12345678910123456789105.5.采用面积坐标时，单元矩阵的计算采用面积坐标时，单元矩阵的计算单元刚度矩阵单元刚度矩阵 体力向量体力向量单元等效结点载荷矩阵单元等效结点载荷矩阵体力引起：体力引起：面力向量面力向量边界面力引起：边界面力引起：（2.5.31）（2.5.32）例：例：（1 1）均质等厚单元的自重）均质等厚单元的自重 单位体积的重量（容单位体积的重量（容重），沿重），沿 y 负方向。负方向。体积力向量体积力向量（2.2.44）（2 2）均布侧压力）均布侧压力q（2.2.45）（2.2.47）（3 3）x方向三角形分布力方向三角形分布力xyxy（2.2.49）（4 4）x方向三角形分布力方向三角形分布力（6 6结点三角形单元）结点三角形单元）0000所求单元等效结点载荷为：所求单元等效结点载荷为：求：图示求：图示 6 6 结点三角形单元的等效结点载荷结点三角形单元的等效结点载荷（1 1）自重作用）自重作用（2 2）均布侧压力）均布侧压力（3 3）边界）边界 x 方向方向均匀分布力均匀分布力思考题：思考题：
