变分原理的直接方法

1、第7章变分问题的直接方法为了书写方便,我们在这里先引入内积空间和线性算子两个概念。内积空间H是复线性空间，H上定义一个两元函数兀尹：H,HTC,x,尹H满足对称性，=尹兀ax+bx,y=ax,y+bx,ya,bC,x,xH,双线性，1212，12正定性，x,x-0,而且只有当x=0时等号才成立那么我们称该两元函数定义了线性空间H上的一个内积。定义了内积的线性空间称为内积空间。=Jb=Xxyiii=1例7.2:H是定义在a,b上连续函数所组成的线性空间Ca,b，(x),(x)H，那么下面定义的就是内积a通常，在连续函数空间中按如下来定义内积=Jbw(x)(x)卩(x)dxa其中w(x)0是个权函

2、数。线性算子A是定义在内积空间H上的一个线性映射A:peH如H；如果存在另一线性映射A*，使得H八AA叩，则A*称为A的伴映射。当H是函数空间时，A称为(线性)算子，A*称为共轭算子；特别当A*=A，A称为对称算子。如果当对称算子A满足A0，并且等号仅当=0时成立，则A称为对称正定算子。7.1里兹方法(Ritz)由线性对称正定算子A及函数f所确定的一个线性泛函为n(u)=-2f,u该泛函的变分为n，2Au,u2=2Auf,u泛函极值问题(也就是变分问题)其所对应的Euler方程为Au-f=0(Au=f)(712)不失一般性,我们假设泛函的边界条件是齐次的，否则我们总是可以通过函数变换来实现齐次

3、的边界条件:u，u+u0u其中0非齐次的边界条件，那么u满足齐次的边界条件。现选定一组满足泛函齐次边界条件的函数序列uiu2un,那么由该函数序列所张成的子空间为U，span(u,uu)，u|u，工au,aeR12nIiiiIIi，1丿该子空间上的每个函数都满足齐次边界条件.里兹法的核心思想就是用上述函数序列所张成的一个线性空间U，笃un)来近似地替代原泛函的定义域空间，然后在线性空间U中找到一个使得泛函口最小的一个u,u,.,u函数，该函数就是原问题的一个近似解。显然满足齐次边界条件的函数序列12n不是u,u,.,u唯一的，如果我们选择了比较合适的函数序列12n，而且该序列的个数足够多时(当

4、然函数序列的个数越多,其张成的子空间就越逼近原来的定义域空间)，那么里兹法所得到的近似解就能很好地逼近原问题的解。具体地讲，由Su2un线性组合成的一个函数为nu，厶auiii，1a,a,.,a其中系数1*2,n为待定的常数，对应的泛函为n(u)，2，2TOC o 1-5 h ziiiiiii，1i，1i，1由于A是线性正定对称算子，那么n(u)=2iiiiiii，1i，1i，1aa2乙ijijii2工fajiii，1i，1j，1i，1，工工Aaaiji，1j，1f,ii式中A，,ijijaa.aaa.a这是一个关于a12n的一个二次型。选择的a12an要使得该函数取到最小值，也1就是说，0,

5、s，1,2,.,nz-丿丿丿as从而有Aa-f，0,s，1,2,.nisisi，1a,a,.,a这是关于一个线性代数方程组。解此代数方程组后得到12n，由此得到原泛函极值问题(或者微分方程边值问题)的近似解为nu，auiii，1如果用向量的形式来表示p，(u,u，,u)T,a=(a,a，,a)T12n12n那么u，atpn(a)，n(u)=-2，aTKa-2qtF其中K，K,K，ijijijF,(f,f，,f,f=12niiK是nxn的矩阵，F是nx1向量。要使口(“)取到最小值，必须口(a),0a这也就是说Ka，F该方程的解为a，K-1F由于A是对称正定算子，可以证明K是对称正定矩阵，上述解

6、必定存在。所以说，通过里兹法，我们可以把一个泛函的极值问题转化成一个函数的极值问题，求解该函数极值问题所对应的代数方程组，就可以得到原问题的近似解。里兹法的关键在于函数序列的选择,如果选择合适的函数序列是该算法最核心之处。Jy，n(兀2y2+xy)dx,y(0)，y(1)，0例7.3求变分问题0的近似解。取(x),a1x(1-x)那么(x),a1(1-2x)Ja，f1a2(x2-4x3+4x4)+a(x2-x3)dx011J，0令da1得到a（兀2一4x+4x4）+（x2一x3）dx，001由此可以求得a，一5116一阶近似解为y(X)，一5X(1一X)16更进一步，可以取近似解为y(x)，a

7、x(1一x)i1i，1例7.4设Jy(X)，n(y2y22xy)dx,y(0)，y(1)，00求变分问题的近似解。该变分问题的精确解为y（x）,sinx一x(x)，现取sinlax(1一x)ki，11，如果n，1,也就是说y(x)，a1x(1x),上面求法一样,到a1=158,也就是说一阶近似解为y(x)，5x(1x)18TOC o 1-5 h zJJc，0如果n,2,巩x),x(1一x)(a1+a2x),代入泛函表达式，并令6a1aa2，得到3a+3a，1101202123a+13a，1102105220由此可以求得a，71,a，71369241也就是说两阶近似解为y(x)，x(1一x)(7

8、1+7x)36941与精确解相比，两阶近似解误差已经非常小。例7.5长度为1，抗弯刚度为EI的简支梁，受均布载荷q的作用。图7.1例7.5图取位移(挠度)的试函数为x,x,2x,3x,+a+a+a1丿21丿31丿41丿4w(x)=a1为了满足两端位移的简支边界条件，取a=-(a+a+a)1234x,2x,3x,4那么w(x)=-(a+a+a234梁内的应变能为d2w,2,dx作用在梁上的外力势能为U=-1qwdx20把挠度的试函数表达式代入总势能表达式中U=U+U12U=11EI12o26a+a1221231丿12x,2,+a2dx=-一2f,u得到一个关于ai(J)a2(I役(Xn)的新泛函

9、*(a(x),a(x).,a(x)1n2nkn于是问题就变为求函数iSn)ak(J)，使得新泛函能取到极小值。这是关于多个一元函数的变分问题，相应的Euler方程一般为常微分方程组(而原来变分问题得到的Euler方程一般为偏微分方程)，求解该常微分方程的边值问题就得到了原变分问题的近似解。和里兹法相比，康托罗维奇法稍显麻烦，因为里兹法最终得到的是代数方程，而康托罗维奇法最终得到的是常微分方程组。但是由于里兹法中的试函数一般都不满Euler方程，而康托罗维奇法中有一部分函数是通过求Euler方程的边值问题得到，所以康托罗维奇法的精度一般要比里兹法来得高。在实x际应用中，一般把变化比较复杂而较难选

10、择试函数的那个变量选作为n。和里兹法一样，要提高康托罗维奇法的精度，一般要增加试函数中的项数，但是相应的Euler方程个数也会有所增加(微分方程组的维数升高)，在实际计算中有很多麻烦。另外一种较为合理的方法是变量轮换法，也就是说交替轮换ak(xn)中的变量来提高计算精度。例7.6矩形截面的柱体扭转问题用应力函数表示的应变余能泛函为取一阶近似解为中(x,y)=(b2，y2)u(x)这里v(y)=(2,y2)满足y=-b边界上的边界条件中(xy)=0新的泛函为*(u)=Ja16b5u+$b3u2一16b3uh一a1533该泛函的Euler方程为u(x)，命u(x)=，忌相应的边界条件为u(a)=0

11、从而得到一阶近似解为中(x,y)=(b2，y2)u(x)这里u(x)=ch-ikachkx,1,k=1也可以直接用里兹法进行求解,取中(x,y)=(b2一y2)(a2一x2)Axmynmn对于矩形截面,取一项中(x,y)=(b2一y2)(a2一x2)A00代入到ababIax丿IT丿-4Pdxdy丄0dA00得到相应的解。7.3伽辽金法(Galerkin)伽辽金法的基础是虚功原理(虚位移原理)：一个平衡力系在任何虚位移中，外力的虚功等于虚应变能“-t”dQ=“pT”udB+“ft”udQ其中f为体积力，p为表面力，为与外力系所平衡的可能应力。”u为虚位移,虚位移对应的虚应变。如果在前面的恒等式

12、(421)中取-=-，u=”u，得到“t”dQ=“E(n)t”udB+“_E(V)t”udQ(7.3.1)”为与从而有“E(V),ft”udQ+“p一E(n)t”udB=0B如果选择的位移试函数不仅满足位移边界条件，同时还满足力的边界条件E(n)-=pB上2上Bi上(7.3.2)即那么就有“E(V)-+ft”udQ=0(7.3.3)和里兹法一样，我们取位移的试算函数为u=auiii=1那么虚位移为u=,auiii=1a因为i是独立的变分，因而有E(V)g+ftudQ=0(7.3.4)iQ这里Q是用位移来表示的应力。更一般地讲，如果微分方程为L(u)=f取试算函数为nu=auiii=1u其中/满

13、足所有的边界条件。那么伽辽金法的积分为(L(U)-f)udx=0,i=1,2,.,ni与里兹法不同的是，在伽辽金法中位移试算函数除了要满足位移边界条件外，一般还需要满足应力边界条件(也就是说要满足所有边界条件)。如果试算函数不满足应力边界条件，那么近似计算的结果可能不是很理想，甚至可能完全是错误的。有时我们也可以根据问题的需要，在伽辽金积分式中添加一个权函数W(x)0以保证收敛速度(L(U)-f)Wudx=0,i=1,2,.,ni该方法也称为加权残值法。例7.7yy=2x,y(0)=y=0取试函数为y(x)=(1-x)xk,k=1,2,k两阶近似为(x)=a(1-x)x+a(1-x)x212那

14、么(x)+(x)-2x=(一2+x-x2)a+(2-6x+x2一x3)a-2x12代入到伽辽金方程组中+x-x2)a+(2-6x+x2-x3)a-2x4(l-x)dx=012+x-x2)a+(2-6x+x2-x3)a-2(1-x)dx=0012积分并求解得到a=142,a=141369241于是近似解为y(x)=-出(1-x)x-廿(1-x)x236941例7.8长度为1，抗弯刚度为EI的悬臂梁，受均布载荷q的作用。图7.2例7.8图梁的平衡方程为EIw一1q(l一x)2=02或者EId4Wdx4取试算函数为w(x)1.Kx=a1一sin2/那么该试函数满足右端的边界条件(M=Q=)EIw(l

15、)=0,EIw(/)=0把位移的表达式积分，再根据左边的边界条件(W=W=)得到w(x)=a/21、2.兀x2l1x2sin一x22l冗这样的试函数可以满足所有的力和位移边界条件。把它代入到伽辽金积分式中EIw一iq(l一x)2w(x)dx=0从中可以求得a，从而得到相应的近似解。7.4有限元法前面讲到，基于最小势能原理的里兹法要求试算函数在整个区域内满足位移边界条件,这对于一些形状比较复杂或者边界条件比较复杂的问题就很难处理，有限元法正好能弥补里兹法的这种缺陷。下面我们用最小势能原理所对应的位移有限元方法加以说明。假如我们在整个区域内取一些点(称为节点)，用这些节点把整个区域划分成一个个子区

16、域(称为单元)，单元和单元之间通过节点连接。图7.3单元和节点假设这些节点的位移已经知道，我们记为d，它是所有节点位移分量所组成的一个向量。每个单元(子区域内)包含其中某几个节点，单元内节点的位移向量为de，显然它是总体位移向量d的一部分。在单元内我们对单元的位移we进行插值，也就是说单元内的位移we可以表示为weNede(7.4.1)其中Ne称为单元形函数矩阵，它需要满足一些特别的性质(见有限元的相关著作)。当然，这样得到的位移在单元内是连续的。因此整个区域内的位移W也可以表示为总体节点位移d的一个插值，也就是说(7.4.2)wNd其中N为形函数矩阵。对应的应变为ET()w=ET()N,d根

17、据本构关系可以得到应力应变关系为D(743)其中D为弹性矩阵，是对称正定。整个区域的总势能为niffftDsdV-fffutfdV-JJUTpdS-fffutG6(x-x)dV2iii这里，f,P,Gi分别为体积力，表面力和作用在位置巴的集中力，6(X_Xi)是脉冲函数。把位移表达式代入总势能中可以得到(7.4.4)口idTKd-dTF2其中K1斤Et()NTDEt()NdV称为总体刚度矩阵，它是对称。而FJJJNTfdV+JJNtpdS+JJJNtG6(x-x)dViii称为等效节点力。根据最小势能原理，弹性力学的解应使总势能取极小值，现在总势能是关于节点位移d的函数，因此有6口(6d)TK

18、d-(6d)tF0也就是说KdF(7.4.)这就是有限元的平衡方程。在代入边界条件(也就是边界上节点位移)后可以得到一个非奇异的线性代数方程组，求解该方程组后得到所有节点的位移，进而求得每个单元内的应变和应力。当然在实际求解的时候，往往是先在单元内进行插值wNede，再求得单元的刚度矩阵Ke和等效节点力向量Fe，然后再把单元刚度矩阵组装成总体刚度矩阵K，把单元等效节点力向量Fe组装成总体等效节点力向量F。因此说，有限元法可以看成是里兹法的一个推广。里兹法在整个区域内采用同一个插值函数，而有限元方法则把整个区域分成一个个子区域，然后在每个区域内对位移分别进行插值。由于有限元插值的区域较小(每个单

19、元)，因此插值函数的形式相对可以比较简单，不需要很高的精度。以上我们讨论的是根据最小势能原理得到的位移有限元方法，同样可以得到基于最小余能原理的应力有限元法。在实际应用中还可以得到基于广义变分原理的混合有限元法，他们的特点是包含两种或者两种以上的场变量。同样可以放松有限元在边界上的连续条件，通过引进交界面上的场变量来建立修正的广义变分原理，得到相应的杂交有限元(见7.6节)。关于有限元方法的详细讨论可以参考相关的著作(如辛可维奇的有限元著作)。7.5有限元方法的收敛性先看一个简单的例子。细绳在拉力N作用下的(小)挠度u可以用二阶微分方程来描述d2uNp(x),u=0atx=0,l(7.5.1)

20、(7.5.2)(7.5.3)的近似解。dx2这里p(x)为横向分布载荷。上述问题也可以化为下列泛函的极小问题来描述，(u)n1n(半)2p()udx02dx假定极小问题的解为ur，,(ur)min,(u)uH01H0是下列函数的集合：在0,1上存在平方可积的一阶导数，且u(0)u(l)0在广义解的意义下，(7.5.1)和(7.5.3)是等价的。现在用有限元方法求(7.5.3)首先，我们指出，要使(752)有意义的一个充分条件是u在1上连续、且存在分段连续的一阶导数。符合这一条件的有限元，称为协调元。将0,1分为n个线段x,x,i1,2,ni1i这里x0 x(n)1。每个线段上用线性插值函数xx

21、xx.u(x)_iu+u,i1,2,，n(7.5.4)xxi-1xxiii1ii1这样得到的函数满足上述的协调要求。u取上式中的/为(7.5.3)中极小问题的真解在节点上的值uur(x)ur,iiixxxxu(x)iUr+i1Ur,I=1,2,nXxi-1xx1ii1ii1u(Z)ur(x)u(x)Q(xx)(xx)=0(x-x)2,xxx2!i-1iii-1i-1i(7.5.5)代入(7.5.2),(u)=,(Ur)+0(d2),d=maxx-xii-1i(7.5.6)将(7.5.4)代入(7.5.3),则把一个泛函求极值问题化为函数求极值问题，记此时的解为U=Ue,ii兀-兀兀-兀.Ue(X)=iUe+i_Ue,I=1,2,