一维定长悬索问题有限元法求解①
一维定长悬索问题的有限元法求解①
摘要:本文讨论一维定长的悬索问题,即长度为l的悬索在均匀载荷力和弹性力作用下产生位移变化,并将其有限元解与真解进行比较,验证有限元解的有效性和精确性。
关键词:有限元方法悬索问题刚度矩阵
中图分类号:u448 文献标识码:a 文章编
号:1674-098x(2011)06(c)-0037-01
1 引言
有限元方法[1-3]是求解各种复杂数学物理方程的重要方法,近
半个世纪来在工程计算中的作用越来越重要。利用该方法可以获取复杂工程问题的各种信息,直接对工程设计进行具体评价,还能对
工程中出现的问题进行具体详细的技术分析。如今有限元方法已经成为科学计算、工程设计、产品优化不可缺少的工具。
2 一维定长悬索问题
许多工程问题可以转化为与之等价的控制方程,通常由其基本方程和平衡方程给出。
2.1 解析解
设模型问题的基本模量为k,作用在悬索上垂直向上的外力为f,悬索沿y轴的正向向上偏移,设荷载f(x)和竖向位移v(x)均为正值,由位移理论知,悬索上的张力t为常数,将y方向上的所有力求和,有(1)
位移理论表明和,可得
一维有限元法
实习三、一维问题的有限元方法 一）实习问题： 设 ''1 4(0,1) (0)0,(1)x u u xe x u u e e -?-+=-∈??==-??， ~ 1()u x e e u -=--令 将原问题的边界条件齐次化 ''~~ 1 ~~4()(0,1) (0)0,(1)0 x xe x e e x u u u u -?-+=---∈??? ?==， 二）算法描述： 1 ()21,1 01101 ()1,011 101 (),1 011 10(),111 [ ()()()]1[()()()()]1[()()()()]1 [ ()(i i i i i i i i i i i i i i i i i i i i i i i i i i i i i i i i i i i i a p h h q h N d h a p h h q h N N d h a p h h q h N N d h a p h h q h h x x x x x x x x ξξξξξξξξξξξξξξξ------------=+++=- +++=- +++=+++???1 2010 )()()]N N d ξξξξ?1 ()2 1,1 0110 1[ ()()()]i i i i i i i i i a p h h q h N d h x x ξξξξ----=+++? 1 ()1 0101 ()110 ()()()()i i i i i i i i i i b h f h N d b h f h N d x x ξξξ ξξξ ---=+=+?? 1，单元剖分 (1,2,,)i i n e =L 2，i=1 ~ ~ 00A b == 3,计算数值积分：()()()()()()1,11,,1,1,,,,,i i i i i i i i i i i i i i i i a a a a b b -----即得单元上的i i A b 4，将i i A b 迭加到总的~~ A b 中 5，若i<=n,则i=i+1并转到底三步；否则继续下一步 6，根据边界条件调整~ ~ A b （掐头去尾），即得 A 和b 7，解线性方程组Au=b,得u 从而的h u
深入浅出的讲清楚有限元法
“有限元法基础及应用”补充讲义（一） 顾克秋 （2005年3月） 一、引子——弹簧单元与弹簧系统 目标：掌握离散结构直接刚度法分析的原理和形式。了解有限元位移法列式的形 式和基本概念。 1、典型弹簧单元分析 写成矩阵符号形式： ? =k F j i i j j ku ku u u k F f +-=-==)(? ?? ?????????--=??????j i j i u u k k k k f f 写成矩阵形式： kd f =（1-1） （1-2） （1-3） 1-2
式（1-2）、（1-3）为弹簧单元的刚度方程，反映了单元特性：节点力与节点位移之间的关系。式中： （注意：单元节点力是节点对单元的作用力） f d k ——单元节点力列阵 ——单元节点位移列阵 ——弹簧单元的刚度矩阵 弹簧单元刚度方程讨论： 1） 有何特点？ 对称、奇异、主对角元素恒正 2） 中元素代表什么含义？ 刚度系数大小等于弹簧刚度；每列元素代表一端固定、另一端产生单位 位移时加在弹簧单元上的节点力。 3）上面单元刚度方程可以求解吗？为什么？ 不可以。刚度方程仅仅表征一个典型单元的弹性特性，单元水平上无法确定单元节点位移。只有把系统中所有单元特性集成后，在系统水平上才可能求出所有未知位移和反力。单元水平上，若已知单元的节点位移，可由刚 度方程求出所有单元节点力分量。若节点力已知，单元节点位移不能确定， 单元可作刚体运动（小位移） 。这也是单元刚度矩阵奇异性的物理解释。 k k 2、弹簧系统整体分析原理
以右图的一个弹簧系统为例，研究如何由单元特性集成系统特性并建立对系统进行求解的控制方程。 由前面得到的弹簧单元的刚度方程公式（1-2），分别写出2个弹簧单元的特性方程如下： 单元1 单元2 （注：右端节点力分量的下标1，2为单元节点的局部编号，上标是单元号） 下面按两个方法完成系统特性的装配和控制方程的建立。并在特定条件下求解。 1）由节点平衡方程导出： 系统处于平衡时，考虑各节点（1，2，3节点）的平衡条件： 节点受到的外载荷与节点受到与其连接的所有单元对其作用力（单元节点力 的反作用力）之和等于零。因此有下列（节点）平衡方程（组）： 把单元特性（1-4），（1-5）代入（1-6）得到： 写成矩阵形式： 2 2 32 11 2211 1f F f f F f F =+==（1-4） （1-5） （1-6） 3 2223322211122 1111)(u k u k F u k u k k u k F u k u k F +-=-++-=-=（1-7） （1-8） 图 1-3
有限元方法讲义
第1讲抛物问题有限元方法 1、椭圆问题有限元方法 考虑椭圆问题边值问题： （1） 问题（1）的变分形式：求使满足 （2） 的性质，广义解的正则性结果。 区域的剖分，矩形剖分，三角剖分，剖分规则，正则剖分条件，拟一致剖分条件。 剖分区域上分片次多项式构成的有限元空间。 的逼近性质，逆性质： 这里，为的插值逼近。 问题（2）的有限元近似：求使满足 （3） （3）的解唯一存在，且满足。 （3）的解所满足的矩阵方程（离散方程组）形式： （4） 刚度矩阵的由单元刚度矩阵组装而成。 模误差分析：由（2）－（3）可得 （5） 由（5）可首先得到 则得到 （6） －模误差分析 设满足 用与此方程做内积，由（5）式和插值逼近性质得到 再利用模误差估计结果，得到 （7） 最优阶误差估计和超收敛估计概念。 当与时间相关时（为抛物问题准备），由（5）式得 （8） 利用（7），类似分析可得 （9） 2、抛物问题半离散有限元方法 考虑抛物型方程初边值问题：
（10） （10）的变分形式：求使满足 （11） （11）的半离散有限元近似：求使满足 （12） 令，代入（12），依次取可导出常微分方程组： （13） 其中为质量矩阵，K为刚度矩阵。。 求解常微分方程组（13），得到代回的表达式，即得半离散有限元解。 定理1.问题（12）的解唯一存在且满足稳定性估计： （14） 证明：在（12）中取得到 整理为（注意是正定的） 对此式积分，证毕。 误差分析。引进解的椭圆投影逼近：满足 （15） 根据椭圆问题的有限元结果可知 （16） 分解误差： 的估计由（16）式给出，只须估计。 由（11），（12）和（15）知，满足 取，类似稳定性论证可得 （17） 可取为的投影，插值逼近等。 由（17）式，三角不等式和（16），得到 （18） 3、抛物问题全离散有限元近似 剖分时间区间：。 引进差分算子： 规定，当为连续函数时，，则有 由此得到 （19） （20） 定义问题（11）的全离散向后Euler有限元近似：求，使满足 （21） 将代入（21）可导出全离散方程组 （22）
有限元一维杆问题解法及程序
解： 第一步——离散 对于一维杆问题，我们先离散成单元，对每个单元作如下计算 [][ ][][ ][]00,,,,0 0)(22 ,,,,2222222=?? ????-??????-?????????=?? ????-?? ? ?????????-?????????=-????-???=-??=-???? ? ?????ΓΓ Γdx x N N u u K x u N N dx x N N u u dx u u N N N N u u x u N N u u dx ux dx x u x u x u u dx ux dx x u u dx x x u u j i x x j i j i j i x x j i j i j i j i x x j i j i j i x x x x x x x x x x j i j i j i j i j i j i j i j i δδδδδδδδδδδδ 其中杆被平均离散为e 个单元（有限元不一定要均分），于是有node=e+1个结点，每个单元长度len=1/e,于是第n 个单元的左端点坐标len n x i )1(-=，右端点坐标nlen x j =； 第二步——刚度矩阵 线性插值有每个单元 ) () (j j e j j j e i x x e len x x N x x e len x x N -=-= -=-=???? ???-=e e B e ?? ? ????--===?e e e e len B B dx B B K T e e T e e x x e j i 对于整体叠加 ????? ? ? ?? ???--=e e e e e K 00020 （K 为一个node×node 阶矩阵） 程序用for 循环给K 赋值 K=zeros(node,node); K1=zeros(node,node); for n=1:(node-1); K1(n:n+1,n:n+1)=[e,-e;-e,e]; K=K1+K; K1=zeros(node,node);
一维有限元法解常微分方程
()()()''2 1 1 '00Fu=-u +u=sin x 1+=f x y(0)0,(1)0F Fu=0Fu vdx=0v ,v(0)=v(1)=0v (x y C ππ?∈ ?? ==?? ??∈? 一维问题的有限元法一.算法构思 考虑下面的两点边值问题 ，（0，1) 设是一个微分算子，则 ，即，且）连续 则把问题中的微分方程化为积分方程，得 ()1 1 11 1100 10011u 'v '+uvdx=fvdx ,u'v'+uvdx, (f,v)=fvdx u C a(u v)=(f,v), v u 1=x ,h=. n n n u v C C x b δ-=∈?∈∈<=???? 令 a 则问题就是求，使得，对于一般的，其范围很广泛，但样条函数理论给我们提供了解决问题的有力工具。 对[0，1]进行等分： 0 一维梁的M A T L A B有限元法分析 问题如下： 梁A B在A和B两端固定，中间点表示为C，在中间区域承受均匀分布的载荷q，如图所示。梁A B的抗弯刚度为E I。 1，使用R i t z法确定点C处的位移和弯矩，并讨论随着包含更多基本函数的准确性。提示：根据偏转曲线的形状，可以选择基函数的形式，其中系数，并且应该由点A或B处的边界条件确定。 2，采用一维有限元法解决问题，并讨论网格越细时的准确性。提示：使用1-D梁单元。 1R i t z法 一维欧拉-伯努利梁的势能如下： 设选择基函数，容易看出基函数满足边界条件 设， 代入势能表达式得到 由于三角函数是正交函数系，所以得到 令q=10N/m m，E=200000M P a，I=10000，L=200m m 在M A T L A B中计算A k前十项 得到 A= -0.008400754770396-0.000320811945459-0.000049922673823 -0.000020050746591-0.000009258470170-0.000003960641302 -0.000001943426801-0.000001253171662-0.000000837688763 -0.000000513299113 计算C点位移,使用1-10个试函数结果如下： 0.0084007547703960.008400754770396 0.008500600118041 0.0085006001180410.008519117058380 0.008519117058380 0.0085230039119830.008523003911983 0.0085246792895100.008524679289510 计算C点弯矩,,使用1-10个试函数结果如下： 0.8291212625436840.7024697829907620.746814316705690 0.7151514468174600.7379958063005830.723923419683592 0.7333220380018130.7254063205297550.732103122461494 0.727037063279377 可以看到，位移收敛是很快的，弯矩收敛速度慢于位移。 M A T L A B代码如下 %量纲为mm，N，k g %定义弹性模量泊松比,定义截面惯性矩 E=200000;v v v v v=0.3;I=10000; %定义梁的长度 一、里兹法与迦辽金法（摘自电磁场有限元方法 金建铭） 1. 里兹法 里兹法是一种变分方法，其中边值问题用变分表达式（也称泛函）表示，泛函的极小值对应于给定边界条件下的控制微分方程。通过求泛函相对于其变量的极小值可得到近似解。 2. 伽辽金法 伽辽金法属于残数加权方法类型，它通过对微分方程的残数求加权的方法得到方程的解。 若u 是方程的近似解，将u 代入方程可得到非零的残数： r Lu f =- u 的最佳近似应能使残数r 在Ω内所有点上有最小值。残数加权方法要求： 0i i R rd ωΩ =Ω=? 这里i R 表示残数的加权积分，i ω是所选的加权函数。 在伽辽金法中，加权函数与近似解展开中所用的函数相同。通常，这样可得到最精确的 解。 二、有限元方法 里兹法和伽辽金法中，在整个解域内找出能表示或至少近似表示问题真实解的试探函数是非常重要的。然而对于许多问题，这个步骤是十分困难的，对二维和三维问题尤其如此。为此，我们可将整个区域划分成小子域，并应用定义在每个子域上的试探函数。因为子域是小区域，因而在每一子域内函数的变化不大，所以定义在子域上的试探函数通常比较简单。这正是有限元法的基本思想。应用里兹法的过程通常称为里兹有限元法或变分有限元法，而应用伽辽金方法的过程通常称为伽辽金有限元方法。 有限元法与经典里兹法和伽辽金法的不同之处是在试探函数的公式上。在经典里兹法和伽辽金法中，试探函数由定义在全域上的一组基函数组成。这种组合必须能够（至少近似）表示真实解，也必须满足适当的边界条件。在有限元法中，试探函数是由定义在组成全域的子域上的一组基函数构成。因为子域很小，所以定义在子域上的基函数能够十分简单。 三、关于形函数（摘自有限元法在电磁计算中的应用 张榴晨） 对于一个待求的微分方程，用一组线性独立的尝试函数i ψ和待定系数i C 来表示方程的近似解，并用加权余数法（迦辽金法）来求解这些待定系数。求解待定系数的代数方程组为： 1 []1,2,,n i j i j i d C q d j n ψψψΩ Ω =??Ω=Ω =∑? ? 这里j ψ为所选择的加权函数，应用迦辽金法时，所选取的加权函数即为尝试函数。 有限元中应用的尝试函数代表了单元上近似解的一种插值关系，它决定了近似解在单元上的形状。因此尝试函数在有限元法中又称为形函数。对于一维有限元来说，形函数为一个直线段；对一维高阶有限元来说，形函数为一个曲线段；对二维一阶有限元来说，形函数为一个平面；对二维高阶有限元来说，形函数为一个曲面；三维有限元来说，形函数为多维平面或曲面。选择形函数时可以使一个任意元上的函数只与该元所对应的节点势函数值有关，而与其它各点的值无关。 1. 一维有限元 电磁学中有限元法简介 摘要：本文简单介绍了有限元法的历史、优点、基本原理及一维的有限元法，并用使用HFSS软件设计了一个各项参数都符合要求的3dB 的功率分配器。 1.有限元法的历史 有限元法起源于航空力学，最早思想是由Courant在1943年提出，但真正确定有限元的学科和命名的则是Clough 在1960年给出。我国著名学者冯康也对有限元法做了开创性贡献。20世纪70年代开始，开始在电磁领域移植有限元法，由于其本身的优点，逐渐成为了电磁场数值分析的一个主要分支。有限元法，简称FEM。有限元法有着扎实的理论基础，所给出的结果是变分稳定的。2.有限元法的优点 ?有限元法采用物理上离散与分片多项式插值，因此具有对材料、边界、激励的广泛适应性。 ?有限元法基于变分原理，将数理方程求解变成代数方程组的求解，因此非常简易 ?有限元法采用矩阵形式和单元组装方法，其各个环节易于标准化，程序通用性强，且有较高的计算精度，便于编制程序和维护，适宜于制作商业软件。?国际学术界对有限元法的理论、计算技术、以及各方面的应用做了大量的工作，许多问题均有现成的程序，可用的商业软件资源相对较多。 3.泛函和变分 最速降线问题 问题的提出：设点A与原点重合，点B的坐标是（a,b）,重物从A点下落，求出A，B之间的使重物从A运动到B的时间最短的路径。 可以求出A 到B 的总时间： T 中只含有y 和'y 。 变分命题的描述 ： T 是(),'y x y 的泛函； y 满足边界条件(0)0,()y y a b == ； 求一个函数y 使得泛函T 有最小值，也就是有极值。泛函T 取极值的时候其变分为0，记为：0T δ= 4. 一维有限元法 有限元法就是变分问题的数值解法，其基本思想是把场方程转化为能量积分的变分问题（能量最小）。 静电场中，所有满足相同边界条件的位函数u 中，真实的u 将保持能量最小。 在这里，通过研究平板电容器的电位分布问题来介绍有限元法的基本思想。 设电位()u x 满足： 220d u dx = ，00u = ，0N u =Φ 。 不妨取（0,1）区间为例来进行分析，这样也具有普遍意义。 T =?一维梁的MATLAB有限元法分析
有限元法的基础理论
有限元法简介