大型龍門(mén)銑床主軸滑枕有限元分析說(shuō)明書(shū)

大型龍門(mén)銑床主軸滑枕有限元分析說(shuō)明書(shū),大型,龍門(mén),銑床,主軸,有限元分析,說(shuō)明書(shū),仿單 基于逆冪法的解決本地修改系統(tǒng)本征解的一種數(shù)值方法K. Krishnapillai, R.Jones摘要：快速再分析修改系統(tǒng)本征解的問(wèn)題具有相當(dāng)大的實(shí)際意義，現(xiàn)在已經(jīng)開(kāi)發(fā)出了用于計(jì)算特征值修改系統(tǒng)的幾種方法。本文提出的正是這種基于逆功率的方法，僅修改部分使用自由度。這種做法不論規(guī)模的修改和使用的濃縮程度的修改部分，都能使精確解的特征值盡快找到。提出的方法的優(yōu)勢(shì)是比較研究解決方案的逆功率方法的若干數(shù)值例子。這一做法將有助于改進(jìn)系統(tǒng)的自由度大和修改系統(tǒng)小的問(wèn)題。關(guān)鍵詞：局部修改系統(tǒng) 逆冪法 特征值 特征向量 再分析 1 引言振動(dòng)工程包括從機(jī)械震動(dòng)到電氣振蕩和擺動(dòng)的大結(jié)構(gòu)的廣泛的各種領(lǐng)域。動(dòng)態(tài)分析的特征值和特征向量在自由振動(dòng)的這種系統(tǒng)（動(dòng)力學(xué)）是一項(xiàng)基本的分支工程，而且大量的數(shù)值計(jì)算方法已經(jīng)被提出了有效地執(zhí)行這種分析。最近取得的進(jìn)展有限元法（ FEM ）和大型計(jì)算機(jī)，使很多分析家治療更大程度的自由（自由度） ，而且越來(lái)越多的軟件正在開(kāi)發(fā)，以適應(yīng)這些因素。但是，如果任何參數(shù)的變化（例如，形狀，材料，初始條件或環(huán)境條件）在先前分析系統(tǒng)，整個(gè)系統(tǒng)都必須再分析。特征值分析是一個(gè)更艱巨的任務(wù)，比較相應(yīng)的靜態(tài)分析，后者的計(jì)算時(shí)間一般是前者的幾倍甚至幾十倍。如果一個(gè)系統(tǒng)進(jìn)行多次修改，因?yàn)榉治鲂抻喓蟮南到y(tǒng)一樣需要時(shí)間，勞動(dòng)力和金錢(qián)，就像分析原始系統(tǒng)一樣，這可以使這一進(jìn)程相當(dāng)費(fèi)時(shí)?！霸俜治觥笔怯脕?lái)執(zhí)行部分分析（而不是全面分析），以便由于設(shè)計(jì)修改或其他一些變化使部分矩陣被更改時(shí)而獲得有效數(shù)據(jù)。如果可以找到采用先前計(jì)算特征值和特征向量表演再分析的方便方法，這將是非常有用的解決涉及部分改變矩陣的特征值問(wèn)題的。再分析方法在許多文件中得到高度重視，這些文件涉及變化的特征值在本地修改系統(tǒng)中的問(wèn)題1-9。許多這些文件都使用微擾理論1，4，7。Fox and Kapoor 1和Rogers 7應(yīng)用一階微分方程的特征值和特征向量的設(shè)計(jì)變量。這種擾動(dòng)的方法是使設(shè)計(jì)變量的初始值假定不變而有相對(duì)較小的調(diào)整時(shí)的是有效解決方案。但是，如果改變?cè)O(shè)計(jì)變量大甚至是增加，這些方法一般是提供低精度的擾動(dòng)。 除Hirai等人2證實(shí)了只使用程序的修改部分以獲取準(zhǔn)確的特征值和特征向量證明以外，還有Parazzola等人 6，10 采用此方法做為特征值和特征向量的阻尼系統(tǒng)的理論解決方案。這些方法在確切的解決辦法的基礎(chǔ)上使用特征值和特征向量的修改制度來(lái)確定一個(gè)基本公式具有相同程度的矩陣代表后修改的制度。他們不論規(guī)模的修改都是有效的。如果特征值可再分析問(wèn)題使用這種凝聚方程，這將意味著，矩陣可以簡(jiǎn)化程度較低，這是一個(gè)減少計(jì)算時(shí)間非常有效的方法。這基本方程是非線性的，它可以用來(lái)解決一個(gè)矩陣方程的組成是否合理的問(wèn)題。但解決這類問(wèn)題的試驗(yàn)和錯(cuò)誤是非常耗時(shí)的，而且通過(guò)這種過(guò)程無(wú)法顯示順序特征。此外，很少有前景的實(shí)現(xiàn)穩(wěn)定的計(jì)算過(guò)程中使用NewtonRaphson法或其他方法，其中的差別系數(shù)為當(dāng)?shù)亟鉀Q辦法和初步估計(jì)而設(shè)置為適當(dāng)?shù)闹?，其解決辦法必須合理地接近實(shí)際的解決辦法（如果不是合理地接近實(shí)際的解決辦法則其進(jìn)程不可預(yù)測(cè)）。即使大致范圍的解決方案已經(jīng)確定，很難找到穩(wěn)定的解決方案，利用逐次逼近的方法，如反向線性插值或多項(xiàng)式逼近。使用這些方法可能提供一些解決辦法，但眾所周知，數(shù)值計(jì)算可以提供解決方案的一個(gè)隨機(jī)秩序。當(dāng)一個(gè)解決方案已經(jīng)被發(fā)現(xiàn)使用，例如，NewtonRaphson法，通貨緊縮是建立和重復(fù)的程序，以確定下一個(gè)特征。因此，傳統(tǒng)的方法解決這些非線性方程組有一些根本性的弱點(diǎn)。Kashiwagi等人對(duì)Hirai等人1116出版的方程壓縮版本進(jìn)行了系統(tǒng)的研究 11-16日 。他們已經(jīng)解決了在當(dāng)?shù)亟?jīng)修改的制度中找到所有特征值的這個(gè)問(wèn)題，提出了找到低度系統(tǒng)為特征的組合DurandKerner 12 法與Newton法 16 和合理的功能 13 的可靠的解決方案。他們還表明，Sturm序列對(duì)某些壓縮的強(qiáng)烈非線性方程組 14,15 是有用的 ，討論了如何確定該地區(qū)的Sturm 序列的特征值是否存在問(wèn)題，并說(shuō)明了使用Sturm序列一分為二的解決方案。Sturm序列方法是計(jì)算特征值常用的傳統(tǒng)計(jì)算方法，一般，它最適合已轉(zhuǎn)化為三對(duì)角矩陣的矩陣 9 。轉(zhuǎn)換的過(guò)程中，需要進(jìn)行大量的計(jì)算，一般以三對(duì)角矩陣形式和代表半數(shù)以上的特征值來(lái)計(jì)算確定。由Kashiwagi等人制定的Sturm序列法的版本需要查明所有的特征值和特征向量的修改制度，但并不需要轉(zhuǎn)化成三對(duì)角矩陣的形式，所以它是這一領(lǐng)域獨(dú)特的貢獻(xiàn)。許多特征解決問(wèn)題，尋求公正的幾個(gè)特征值及其相應(yīng)的來(lái)自低度系統(tǒng)的載體。逆功率職能是這種系統(tǒng)用迭代方法來(lái)獲得解決方案的一個(gè)例子 17 。逆冪法是一個(gè)獲得解決方案的特征值分析的基本方法;許多分析都是根據(jù)它得到的。例如，子空間迭代 18,19 是一個(gè)確定本征值非常常用的方法。因此，我們必須早日完成工具箱的方法來(lái)確定特征值的本地修改的制度，但目前還沒(méi)有一種已提議的基于逆功率的職能。本文提出的正是這種僅修改部分使用自由度的基于逆功率的方法。這種做法不論使用的濃縮程度的修改部分的規(guī)模都能盡快找到精確解的特征值。當(dāng)程度的矩陣低，計(jì)算時(shí)間很短，特別是當(dāng)逆冪函數(shù)相結(jié)合轉(zhuǎn)變的起源，所有的特征值和特征向量的修改系統(tǒng)從最小的本征解開(kāi)始對(duì)前幾個(gè)本征解（ 10個(gè)或更少的實(shí)際系統(tǒng)）相對(duì)較小的系統(tǒng)是否有效必須被了解。以下各節(jié)中描述的理論和算法證明了數(shù)值求解典型特征值問(wèn)題這一有效的做法。2 本地修改系統(tǒng)的理論與算法的逆功率方法 本節(jié)描述修改系統(tǒng)的逆功率方法的理論和什么是本地修改系統(tǒng)轉(zhuǎn)向逆功率方法的理論。這些理論使只是濃縮版修改部分的本征解完全確定。2.1 本地修改系統(tǒng)的逆功率方法 一般特征值問(wèn)題如下，假設(shè)一個(gè)nn實(shí)對(duì)稱矩陣A和一個(gè)積極的實(shí)對(duì)稱矩陣B組： （1） 從最小的數(shù)值和是相應(yīng)的特征值特征向量開(kāi)始，是隨著那里的特征值。A和B有以下幾種形式： （2） （3）如果是模式矩陣的特征向量包含的問(wèn)題，則用等式（1）表示 （4）然后用下列關(guān)系表示： （5） （6）這里，是n n的矩陣，由等式 （ 5 ）可得 （7）通常由有限元分析的自由結(jié)構(gòu)振動(dòng)和屈曲有兩個(gè)特征值。在自由振動(dòng)問(wèn)題中，A是剛度矩陣K，B是質(zhì)量矩陣M。在屈曲中， A是剛度矩陣K，B是幾何剛度矩陣KG。根據(jù)當(dāng)?shù)氐男薷模?A是由A+A替代，B是由B+B替代的。讓我們假設(shè)A和B如下： （8） （9）如果零要素從A變到B ，我們會(huì)得到和： （10） （11）這是化簡(jiǎn)后的部分（mm矩陣，在等式（8）和（9）中假設(shè)m=2） 。我們考慮一個(gè)m n的布爾矩陣（使等式（12）中m = 2） ： （12）和由下式給出 （13） （14）一般等式（1）表示的特征值問(wèn)題可以用來(lái)表示上述化簡(jiǎn)的當(dāng)?shù)鼐仃嚕?（15）等式（15）可重寫(xiě)布爾矩陣形式： （16）讓我們重寫(xiě)上述矩陣，分別用P代表n n矩陣，Q代表n m矩陣，R代表m n矩陣，代表mm對(duì)稱矩陣。然后，我們可以寫(xiě)成 （17）的逆矩陣是 （18）這里，是修改部分元素對(duì)應(yīng)的mm矩陣 。用逆功率的方法，如果后面的特征值試驗(yàn)值為然后等式（16）和（18）可得 （19）我們又可以寫(xiě)成 （20）是下面等式（25）中的；是前迭代中獲得的近似特征向量?？紤]到式中的 （等式（9）中）可得 （21）這里，和是未修改部分的矢量 （22） （23）然后，由下式給出 （24）這里 （25）因此，我們可以用數(shù)值本身的表示來(lái)尋找該特征向量。同樣在等式（24）中我們只需要只用簡(jiǎn)明修改部分的自由度的化簡(jiǎn)特征向量： （26） （27）一旦等式（26）中的被找出和保存，則由等式（24）所得出的自由度就可以計(jì)算出近似的特征向量。2.2 本地修改系統(tǒng)的逆功率轉(zhuǎn)換法 讓我們指定當(dāng)前特征值和下一個(gè)特征值為。冪函數(shù)的收斂速度由|/|給出。當(dāng)原特征值為時(shí)，逆冪函數(shù)的收斂速度隨原值的改變由|/|給出，其值是一個(gè)比未修改的逆冪函數(shù)較低的值。 因此，選擇適當(dāng)?shù)闹翟谠俜治龅倪^(guò)程中收斂可能會(huì)被大大提高。在本節(jié)中，我們得出一個(gè)本地修改系統(tǒng)的轉(zhuǎn)換逆功率方法。等式（28）的基本方程具有跟一個(gè)完整的系統(tǒng)逆功率轉(zhuǎn)換方法相同的格式，因此預(yù)計(jì)將對(duì)本地修改系統(tǒng)是有效地。 此處重述等式（1）代表本地修改系統(tǒng)逆功率轉(zhuǎn)換方法作為一般特征值問(wèn)題： （28）等式（7）也可被重新寫(xiě)為： （29）我們按照上一節(jié)中給出的相同的方法： (30) (31)這確保我們說(shuō)明地點(diǎn)時(shí)尋找到近似的特征向量。而且，由等式（30）可得到濃縮修改部分基礎(chǔ)上的近似特征向量： (32) (33)一旦用等式（32）算出并保存，近似特征向量可以用等式（30）給予的自由度計(jì)算得到。2.3 初始近似特征向量和初始近似特征值 這種重新分析的方法需要以往修改系統(tǒng)的特征向量和特征值的知識(shí)。我們討論如何將這些數(shù)據(jù)可用于生成隨著i的初始近似特征值和特征向量，我們由下式開(kāi)始： (34)考慮到正常化的條件為：我們發(fā)現(xiàn)： (35)表明如果我們考慮到，可得： (36)初始的近似特征向量只占簡(jiǎn)明修改部分： (37)此外，考慮到初始近似特征值為：并考慮本地化的和，我們發(fā)現(xiàn) (38)2.4 特征向量正?；吞卣髦?我們現(xiàn)在描述如何隨著i近似特征向量正常化以及如何計(jì)算近似特征值。我們指定 (39)和條件由等式（24）、（30）和的本地化，我們可以得出 （40）然后 （41）并且初始近似特征向量根據(jù)濃縮修改部分 （42）我們指定初始近似特征值為并且由等式（24）、（30）和和本地化，我們可以得到為： （43）這操作要求用代替2.5 提取已取得的特征值 當(dāng)發(fā)現(xiàn)高階特征值，有必要提取使用逆冪函數(shù)已經(jīng)發(fā)現(xiàn)的本征解。革蘭氏Schmidt正交是用來(lái)配合等式（6）、（24）、（30）和的本地化： （44）并且： （45）方程僅由簡(jiǎn)明修改部分得到。如果我們定義的組成部分為，那么我們可以重寫(xiě)如下： （46）2.6 算法 本節(jié)說(shuō)明一個(gè)根據(jù)擬議法創(chuàng)建算法的例子：（1）K0（初始值）。和版本的特征值，特征向量和mm矩陣（m：一定程度的修改部分）可以得出，并且mq（數(shù)目的本征解）初始近似特征向量和特征值可以使用等式（37）和（38）算出?？梢杂玫仁剑?5）算出并保存。得出的特征值按升序重新排列。（2）KK1（計(jì)算K的特征值）。1、當(dāng)上一步計(jì)算的近似特征值和本步驟中的近似特征值之差的絕對(duì)值除以本步驟近似特征值所除得的值大于一些指定的值（在實(shí)驗(yàn)中用計(jì)算）則用逆功率近似法；當(dāng)小于時(shí)則用轉(zhuǎn)向逆功率近似法。在目前的數(shù)值試驗(yàn)中，轉(zhuǎn)向的起源（是上一步驟中的近似特征值）。2、簡(jiǎn)明特征向量和值可由等式（32）和（33）計(jì)算得出。3、（K1）組的特征向量可由等式（45）計(jì)算得出。4、簡(jiǎn)明近似特征向量的規(guī)范化使用等式（42）和（43），并可以計(jì)算出近似值。5、收斂測(cè)試；如果計(jì)算值被收斂，則用等式（30）所得的自由度計(jì)算特征值，否則，該程序返回到步驟1。（3）如果沒(méi)有得到指定的一些特征值，則該程序返回到步驟（2）。3 數(shù)值實(shí)驗(yàn)圖1. 平面框架模型（ A型）本實(shí)驗(yàn)使用的是如圖1顯示的利用有限元法求解包含關(guān)節(jié)強(qiáng)度的結(jié)構(gòu)框架的一般特征值問(wèn)題。當(dāng)上一步的近似特征值與目前步驟中的近似特征值所確定的差的絕對(duì)值除以目前步驟中的近似特征值小于時(shí)，每一個(gè)特征值則被決定。所有的數(shù)值計(jì)算都用雙精度。計(jì)算機(jī)使用的是一個(gè)富士通FMV（賽揚(yáng)處理器2.40GHz，RAM768MB，Windows XP系統(tǒng)，F(xiàn)ujitsu Fortran和C Academic Package Ver.3）。由鋼筋混凝土制成的柱子和梁都具有楊氏模量和密度。柱子和梁的橫截面分別是800mm800mm和400mm800mm（寬高）。柱子的長(zhǎng)度都是3000mm（落地式地板高度），梁的跨度都是6000mm。結(jié)構(gòu)中的每個(gè)節(jié)點(diǎn)有三個(gè)自由度（X方向的位移u ， Y方向的位移V和傾斜角度）。剛度矩陣A和質(zhì)量系數(shù)矩陣B中的元素值與變量的角度傾斜是對(duì)應(yīng)的，并允許是不同的層面。矩陣B表示分布式質(zhì)量，而不是集中質(zhì)量。因?yàn)橥ǔＴ诠こ虇?wèn)題中的需要，前18個(gè)特征值的計(jì)算可以用兩個(gè)建議的方法和逆冪函數(shù)的方法。表1給出了數(shù)值實(shí)驗(yàn)中的參數(shù)。為計(jì)算五個(gè)級(jí)別的自由度n，如圖1所顯示的修改了地點(diǎn)的和；這些列被刪除。當(dāng)只有列被移動(dòng)時(shí)標(biāo)明“m=3”，當(dāng)列都被移動(dòng)時(shí)則標(biāo)明“m=6”。表1 例舉參數(shù)類型mrnrnA88216B1212468C1616816D20201260E24241800F28282436mr:存儲(chǔ)數(shù)目，nr:跨度，n:自由度 表2 所得特征值的舉例（類型 A,m=3）序號(hào)提出的方法逆冪法10.2728886410530.27288864105322.7128022794292.71280227943034.1321810033104.13218100331048.9498442315268.949844231525517.11031834738117.110318347409618.45520191047218.455201911020718.69479049731018.694790496682820.12465734931120.124657349413921.27296814458621.2729681444881022.39193961816622.3919396184381123.16126469718323.1612646969281227.05263110457827.0526311046711328.73437771279828.7343777129461429.78342270647429.7834227062241534.87794693842834.8779469413301635.08136401958535.0813640166701740.49536753658540.4953675368771843.00477098790743.004770987659圖2.計(jì)算時(shí)間為最低特征值18（m=3）圖3. 迭代次數(shù)（類型E，m= 3 ）圖4. 計(jì)算時(shí)間為18最低特征值（m=6） 圖5.迭代次數(shù)（類型E，m=3） 表2列舉了一個(gè)所取得結(jié)果的例子。很顯然，本文提出的這種方法能正確計(jì)算出特征值，因?yàn)樗鼈兣c其他計(jì)算方法的預(yù)測(cè)是一致的。圖 2和圖4提供了這些結(jié)果和逆冪函數(shù)方法所預(yù)測(cè)的比較（帶系數(shù)矩陣）和CPU所需要的時(shí)間。在m=3時(shí)，用逆功率方法解決本地系統(tǒng)修改部分所需要的時(shí)間大約為該逆功率方法解決整個(gè)系統(tǒng)所需時(shí)間的1/25，而且轉(zhuǎn)向逆功率方法解決本地系統(tǒng)修改部分所需要的時(shí)間大約為逆功率方法解決整個(gè)系統(tǒng)所需時(shí)間的1/100。在m = 6時(shí) ，用逆功率方法解決本地系統(tǒng)修改部分所需要的時(shí)間大約為逆功率方法解決整個(gè)系統(tǒng)所需時(shí)間的一半，而且轉(zhuǎn)向逆功率方法解決本地系統(tǒng)修改部分所需要的時(shí)間大約為逆功率方法解決整個(gè)系統(tǒng)所需要時(shí)間的1/30。 圖3和圖5比較了三種方法所需要的迭代次數(shù)。在m=3時(shí)，用逆功率方法解決本地系統(tǒng)修改部分的收斂迭代次數(shù)遠(yuǎn)小于用逆功率方法解決整個(gè)系統(tǒng)的收斂迭代次數(shù)，而在m=6時(shí)，幾乎沒(méi)有區(qū)別這兩種方法。本地修改部分的逆功率轉(zhuǎn)向方法需要比其他兩種方法更少的迭代次數(shù)，保證使用這種方法時(shí)有較低的計(jì)算。4 總結(jié) 本文提出的理論研究，其目的是利用以前獲得的數(shù)據(jù)用于未修改系統(tǒng)，使修改系統(tǒng)再分析獲得其新的特征值和特征向量，是一種穩(wěn)定且有效的分析方法。這種方法的主要特點(diǎn)是經(jīng)過(guò)修改可以得到一個(gè)壓縮的低階矩陣，而不論修改的規(guī)模，這矩陣為修改后的特征值和特征向量提供了一個(gè)準(zhǔn)確的解決方案。這個(gè)較低階矩陣為大大縮短計(jì)算時(shí)間提供了解決方案。 在一個(gè)部分修改系統(tǒng)特征值分析的基本方程基礎(chǔ)上產(chǎn)生了逆功率法，并且這種方法的有效性和穩(wěn)定性證明了數(shù)值分析。在這一分析，也有人認(rèn)為，修改部分的轉(zhuǎn)向逆功率分析特別有效。它比修改部分和整個(gè)系統(tǒng)的未轉(zhuǎn)向逆功率分析更急劇地減少迭代而收斂。然而，數(shù)量的計(jì)算主要是確定原系統(tǒng)的自由度，特征值和特征向量，而且代表一定修改程度的矩陣中的數(shù)據(jù)必須在采用這種方法之前就已獲得。預(yù)計(jì)該方法也可與斯特姆序列或二等分法結(jié)合進(jìn)一步提高穩(wěn)定性。未來(lái)的研究需要提供這種算法的進(jìn)一步驗(yàn)證和處理理論方面的裁斷。參考文獻(xiàn)1 R.L. Fox, M.P. Kapoor, Rates of eigenvalues and eigenvectors, AIAA J. 6 (1968) 24262429.2 I. Hirai, T. Yoshimura, K. Takamura, On a direct eigensolution of locally modified structure, Int. J. Numer. Methods Eng. 6 (1973) 441442.3 S.A. Jasbir, Survey of structural reanalysis technique, in: ASCE ST4 Proceedings Paper, vol. 12056, 1976, pp. 783802.4 U. Kirsch, Approximate structural reanalysis for optimization along a line, Int. J. Numer. Methods Eng. 18 (1982) 635651. 5 L. Kitis, W.D. Pilky, I. Hirai, Reanalysis, in: Finite Element Handbook, McGraw- Hill, New York, 1987. 6 A.B. Parazzola, Vibration of locally modified mechanical and structural systems, Doctoral Thesis, University of Virginia, 1981. 7 L.C. Rogers, Derivatives of eigenvalues and eigenvectors, AIAA J. 8 (1970) 943945. 8 B.P. Wang, W.D. Pilkey, Efficient reanalysis of locally modified structures, in: Proceedings of the First Chautaqua on Finite Element Modeling, Schaeffer Analysis, 1980.9 J.H. Wilkinson, The Algebraic Eigenvalue Problem, Oxford University Press, Oxford, 1965. 10 B.P. Wang, A.B. Parazzola, W.D. Pilkey, Reanalysis Modal Synthesis and Dynamic Design in State-of-the-Art Surveys of Finite Element Methods, ASME, 1983. 11 I. Hirai, M. Kashiwagi, Derivatives of eigenvectors of locally modified structure, Int. J. Numer. Methods Eng. 11 (1977) 17691773. 12 M. Kashiwagi, A numerical method for eigensolution of locally modified systems by DurandKerner like method, J. Struct. Eng. 417 (1990) 7177. 13 M. Kashiwagi, I. Hirai, T. Katayama, A numerical method for eigensolution of locally modified systems by a rational function approximation, J. Struct. Eng. 37A (1991) 271278. 14 M. Kashiwagi, I. Hirai, S. Ohwaki, W.D. Pilkey, Stable eigensolution of locally modified systems based on the Sturm sequence property, Finite Elem. Ana. Des. 9 (1991) 133139.15 M. Kashiwagi, I. Hirai, T. Katayama, W.D. Pilkey, Sturm sequence recurrence formula of locally modified systems, Finite Elem. Anal. Des. 12 (1993) 3136. 16 T. Katayama, S. Miyamura, M. Kashiwagi, I. Hirai, Eigensolution of locally modified systems by Newton-like method, J. Struct. Eng. 38A (1992) 303310. 17 M. Kashiwagi, An eigensolution method of large sparse symmetric matrices by the CG method, Trans. Jpn. Soc Ind. Appl. Math. 15 (2005) 2943. 18 K.J. Bathe, Finite Element Procedures in Engineering Analysis, Prentice-Hall, Englewood Cliffs, NJ, 1982. 19 B. Nour-omid, B.N. Parlett, R.L. Taylor, Lanczos versus subspace iteration for solution of eigenvalue problems, Int. J. Numer. Method Eng. 19 (1983) 859971.