lagrange方程中的系數(shù)矩陣

lagrange方程中的系數(shù)矩陣

1lagrange方程的求解等式束縛變換問題的兩種一般方法：拉格蘭乘數(shù)法和懲罰函數(shù)法。罰函數(shù)法的主要困難在于罰參數(shù)的選取，參數(shù)過大會(huì)導(dǎo)致方程組的病態(tài)和收斂上的困難，參數(shù)過小會(huì)導(dǎo)致解的失真。但是因?yàn)榱P函數(shù)法導(dǎo)致的方程組是正定的，便于求解，所以目前所流行的一些有限元或DDA程序等大都還是采用罰函數(shù)法。Lagrange乘子法能精確滿足約束方程，但它所導(dǎo)致的方程組(本文稱之為L(zhǎng)agrange方程組)的系數(shù)矩陣(本文稱之為L(zhǎng)agrange矩陣)是一對(duì)稱不定矩陣，求解起來會(huì)存在一些困難。A.F.Saleeb等目前關(guān)于Lagrange方程組的算法中，無論是迭代法，如Uzawa法Lagrange矩陣含帶狀稀疏的剛度陣和約束矩陣，因而具有所謂的“對(duì)角帶雙邊”(doublyborderedbanddiagonal)這樣一個(gè)特殊的稀疏結(jié)構(gòu)本文所建議的Lagrange方程組的直接解法，既適應(yīng)于滿秩的剛度矩陣，又適用于虧秩的剛度矩陣。同時(shí)本文還討論了Lagrange方程組的解存在唯一的充要條件以及在無單元Galerkin法中的應(yīng)用。2廢墟中的材料2.1根據(jù)順序ldl眾所周知，當(dāng)方程組為對(duì)稱正定時(shí)，基于LDL設(shè)矩陣K為n階對(duì)稱正定矩陣，則有如下形式的分解：式中：D=diag(l式中：2.2–m矩陣uv假定已知n階方陣A的逆A更進(jìn)一步，如果用一秩–m矩陣UV因?yàn)槿萘烤仃嘔一般認(rèn)為式(3)最早是由Sherman和Morrison所建議的，而后Woodburg將其推廣為式(4)3lagrange矩陣如果將Lagrange乘子法用于求解等式約束的變分問題，則可得到如下形式的Lagrange方程組：式中：K為n階剛度矩陣，是對(duì)稱并至少是半正定的；u為n階節(jié)點(diǎn)位移向量；f因?yàn)橐话闱闆r下m<<n且K是帶狀稀疏的，因此式(5)中的Lagrange矩陣A具有所謂的“對(duì)角帶雙邊”(doublyborderedbanddiagonal)這樣一個(gè)特殊的稀疏結(jié)構(gòu)本文始終假定C為滿秩?？梢宰C明Lagrange矩陣A非奇異的充分必要條件為N(K)∩IN(C證明：來考察齊次形式的Lagrange方程組：將u由于式(6)的第2個(gè)方程會(huì)使得式(7a)的第2項(xiàng)消失，所以有利用半正定矩陣的平方根分解因?yàn)閙<<n，通過解齊次方程組C4基于剛性矩陣分解的解算方法本文將根據(jù)K滿秩和K虧秩這兩種情況來討論式(5)的求解。4.1．求解式的確定假定K=LDL其中，式中：L分解式(8)，通過依次求解如下3個(gè)方程組來得到式(5)的解：從計(jì)算角度可將M寫成：求解式(9a)和(9c)時(shí)，需要計(jì)算L所以一旦求得了L值得注意的是，若矩陣C非列滿秩，由式(10)可知M必然奇異，從而也必將導(dǎo)致M的LDL4.2e為m+3階的矩陣這是Lagrange乘子法應(yīng)用中一種更常見的情況。以EFGM為例，K的虧秩數(shù)通常為3，對(duì)應(yīng)于N(K)中的3個(gè)剛體位移向量。為了克服EFGM中的剛體位移問題，龐作會(huì)等式中：K′為K的劃掉了最后3列和3行所對(duì)應(yīng)的主子陣；C′為塊K′和C組成，而K′是由K的最后3列并去掉其最后3行后組成的矩陣；E為(m+3)階的矩陣，其中除了K′的余子陣(三階)外，其他元素都為0。然而，如果研究一個(gè)多體接觸問題，其剛度矩陣K可能具有如圖1所示的存儲(chǔ)結(jié)構(gòu)如果仍然采用LDL假定在分解完K其中，式中：e其中，利用Sherman-Morrison公式并將-b將式(17)右乘以向量b可更簡(jiǎn)捷地表示成其中，這里及以后x由(18a)可知，為了求得x在隨后對(duì)K其中，亦即再次利用Sherman-Morrison公式，有將式(22)代入式(19)，可得再次記x其中，依此類推，假定在分解完K其中，今后，將統(tǒng)一地稱b仍然記x其中，注意式(26a)中所涉及的運(yùn)算屬于LINPACK定義的“saxpy”型運(yùn)算現(xiàn)在假定在隨后對(duì)K利用前述K滿秩情況下的Laplace方程組的算法，通過連續(xù)求解方程組：來得到向量x一旦求得x現(xiàn)在來考慮式(13)有多個(gè)給定向量b剩下的問題就是g由于本文堅(jiān)持了對(duì)小主元的擾動(dòng)策略，因此上述分解過程的數(shù)值穩(wěn)定性是得到保障的。5防止約束的形成基于上述分析，可設(shè)計(jì)出以下算法：(1)按照式(2)來分解K。注意每分解完一行，如第k(2)由式(10)來形成M，由式(2)來分解M。如果M的分解失敗，表明約束是線性相關(guān)的，此時(shí)應(yīng)終止程序的運(yùn)行，否則繼續(xù)向下進(jìn)行。(3)依次以p+1個(gè)向量b(4)按照迭代格式(26)和圖2依次將第j個(gè)向量層x6在admm中的應(yīng)用可以毫不夸張地說，無網(wǎng)格法是過去的10a間最為活躍的研究領(lǐng)域之一，迄今為止，已經(jīng)產(chǎn)生了約20種方法，幾本新近出版的專著在眾多的無網(wǎng)格法中，利用移動(dòng)最小插值的無單元Galerkin法(EFGM)既然本文所建議的方法能夠克服在求解Lagrange方程組時(shí)所導(dǎo)致的困難，作者建議Lagrange乘子法應(yīng)作為求解帶等式約束的變分問題的首選方法。若采用Lagrange乘子法來引入本質(zhì)邊界條件，則EFGM最終都將導(dǎo)致形如式(5)的Lagrange方程組，詳細(xì)推導(dǎo)過程可見相關(guān)研究圖3顯示了一個(gè)在其右端受均布面力p=1Pa的懸臂梁。其長(zhǎng)、寬、高分別為8,1,1m。彈性模量E=100kPa，泊松比ν=0.25。將用EFGM和FEM來分析本算例。網(wǎng)格含10×4個(gè)完全相同的矩形單元，它既作為EFGM的背景網(wǎng)格又作為FEM的網(wǎng)格，節(jié)點(diǎn)也同時(shí)為EFGM和FEM所共享。在進(jìn)行單元積分時(shí)，EFGM的積分階為4×4,FEM的積分階為2×2。EFGM的插值基取為線性基。在進(jìn)行有限元分析時(shí)，采用了兩種單元類型：等參元和Wilson元。因?yàn)樗械膯卧紴榫匦?，Wilson元自然滿足分片檢驗(yàn)圖4所示為對(duì)應(yīng)于解析解、EFGM、等參元和Wilson元的撓度。其中Wilson元的結(jié)果幾乎與精確解一致，所以此圖看起來似乎僅顯示了3種結(jié)果。由圖4可知，在相同的節(jié)點(diǎn)數(shù)下的EFGM的精度要遠(yuǎn)高于等參元，但比Wilson元的結(jié)果還是略差一些。順便指出，已用本文所建議的算法求解了許多算例，并且與MATLAB中的處理稀疏矩陣的全主元消元法的有關(guān)函數(shù)做了對(duì)比