【文章內(nèi)容簡(jiǎn)介】 元 等效結(jié)點(diǎn)恢復(fù)力的計(jì)算式為： ? ? ? ? ? ? ? ?? ?11nn Te ijijF H H B t J???? ?? 其中， ??B 為應(yīng)變矩陣， J 為雅克比行列式。 ??B 和 J 的值與高斯積分點(diǎn)相對(duì)應(yīng)。 針對(duì)不同的單元，一般通過(guò)經(jīng)驗(yàn)和試算來(lái)確定最佳的積分點(diǎn)數(shù)。對(duì)于四邊形單元而言， 22? 的高斯積分點(diǎn)就可以較好地滿足精度要求 。表 中所示的是各個(gè)高斯積分點(diǎn)的坐標(biāo)和加權(quán)系數(shù)，圖 所示的是 22? 高斯積分點(diǎn)的分布情況。 表 高斯積分中的 i? 和 iH Table i? and iH in Gauss Integration 點(diǎn)數(shù)n 積分點(diǎn)坐標(biāo) i? 加權(quán)系數(shù) iH 點(diǎn)數(shù)n 積分點(diǎn)坐標(biāo) i? 加權(quán)系數(shù) iH 1 1? ＝ 0 1H ＝ 2 4 1? ＝ 2? ＝ 3? ＝ 2?? 4? ＝ 1?? 1H ＝ 2H ＝ 3H ＝ 2H 4H ＝ 1H 2 1? ＝ 1? ＝ 2?? 1H ＝ 1 2H ＝ 1 3 1? ＝ 2? ＝ 0 3? ＝ 1?? 1H ＝ 2H ＝ 3H ＝ 1H 5 1? ＝ 2? ＝ 3? ＝ 0 4? ＝ 2?? 5? ＝ 1?? 1H ＝ 2H ＝ 3H ＝ 4H ＝ 2H 5H ＝ 1H ζη 圖 高斯點(diǎn)的位置 Fig The location of Gauss Points NewtonRaphson 迭代法 基于切線剛度的 NewtonRapson 迭代法如圖 。 圖 基于切線剛度的 NewtonRaphson 算法 NewtonRaphson algorithm based on tangent stiffness NewtonRaphson 是一種變剛度 的迭代法。 首先 ， 取 單元的 初始剛度 為 ? ?0K ，根據(jù)所施加的外荷載等效結(jié)點(diǎn)力 ??P 可以求出 位移的第一次近似值 為 ： ? ? ? ? ? ?110KP? ?? () 由 得到的 初始位移可 以計(jì)算 單元 的 應(yīng)變 值，再結(jié)合本構(gòu)模型求到單元應(yīng)力，進(jìn)而得到相應(yīng)的結(jié)點(diǎn)恢復(fù)力 ??1P 。 然后，用外荷載等效結(jié)點(diǎn)力 ??P 和結(jié)構(gòu)恢復(fù)力 ??1P 求出殘余力 ? ?1P? ，并連同 相應(yīng) 于初始位移 ??1? 的即時(shí)切線 剛度 ? ?1K ， 求出 位移增量? ?2?? ，即 ? ? ? ? ? ?11P P P? ? ? () ? ? ? ? ? ?12 1 1KP? ?? ? ? () 得到 位移的第二次近似值 為： ? ? ? ? ? ?2 1 2? ? ?? ? ? () 以此類(lèi)推，第 k 步的計(jì)算過(guò)程為： P 0 1?0 2?0 3?0 1P?0 1P 2P?0 2P 1??0 2??0 3??0 ?0 ? ? ? ? ? ?kkP P P? ? ? () ? ? ? ? ? ?11k k kKP? ??? ? ? () ? ? ? ? ? ?11k k k? ? ???? ? ? () 當(dāng) ? ?1k?? 與 ??k? 的誤差達(dá)到所規(guī)定的限值 ，或 者 當(dāng) ? ?kP? 足夠小的時(shí)候，迭代結(jié)束 。 弧長(zhǎng)控制法 弧長(zhǎng)控制法的基本原理是，通過(guò) 在結(jié)構(gòu)分析的一般方程 ? ?iiF x f? ?? 中 添入一個(gè) 新的 變量 ? ， 來(lái)實(shí)現(xiàn)對(duì)整個(gè)計(jì)算過(guò)程的 自動(dòng)控制 。 加入變量 ? 后， 原方程組 將 變?yōu)?： ? ?iiF x f?? ? ? （ ） 上述 方程組 一 共有 n+1 個(gè)未知量： nxR?? 、 ? ， 但是 方程 卻 只有 n 個(gè) 。所以，需要補(bǔ)充如下條件來(lái)完成方程組的求解： 2222f x S? ? ? ? ? （ ） S 稱(chēng) 為 “ 弧長(zhǎng) ”，是 附加的控制變量 。這時(shí)，方程組變?yōu)?? ? ? ?11nnFR? ? ??? ?；￠L(zhǎng)法的迭代過(guò)程如圖 所示。 圖 弧長(zhǎng)法迭代過(guò)程示意圖 Fig The Arclength method 收斂準(zhǔn)則 基于切線剛度的 NewtonRaphson迭代法的收斂準(zhǔn)則有三種：力準(zhǔn)則、位移準(zhǔn)則以及能量準(zhǔn)則。選取合適的收斂準(zhǔn)則，有利于提高迭代法的計(jì)算效率，是增量求解策略中一個(gè)重要的組成部分。 當(dāng)?shù)Y(jié)束時(shí)，結(jié)構(gòu)應(yīng)當(dāng)處于一個(gè)穩(wěn)定的力平衡狀態(tài)，力準(zhǔn)則就是要求結(jié)構(gòu)的恢復(fù)力和外力 之差等于零。但是在數(shù)值計(jì)算中，不可能實(shí)現(xiàn)不平衡力完全為零，這就需要設(shè)定一個(gè)允許值。力準(zhǔn)則的表達(dá)式為： ? ? ? ? ? ? ? ?1 1 1 1 22im m m mFR F R F?? ? ? ?? ? ? 其中， F? 就是給殘余力設(shè)置的一個(gè)允許值。 在進(jìn)行位移分析時(shí)，計(jì)算位移應(yīng)該與真實(shí)值接近。位移準(zhǔn)則的表達(dá)式為： ? ? ? ? ? ? ? ? ? ?1122iimmDU U U? ??? ? ? 其中，2是 矢量的歐幾里德范數(shù)， ? ???iU? 是第 i 步迭代的位移增量， D? 為 設(shè)置的允許值。 能量準(zhǔn)則是用來(lái)度量力和位移同它們的平衡值之間的接近程度的，其表達(dá)式為： ? ? ? ? ? ? ? ?? ? ? ? ? ? ? ? ? ?? ?1 1 1i T i i Tm m m mEU R F U R F?? ? ?? ? ? ? ? （ ) 其中，不等式 的 左邊 表示殘余力 在位移增量上所做的功 ，不等式的 右邊是 殘余力所 做 的 功的初始值 ， E? 是 為能量設(shè)置 的 允許值。 鄧興龍和蔣曉華等編制的 FEAPpv程序采用的是能 量收斂準(zhǔn)則。 基于 MCFT 理論的 割線剛度迭代法 F. J. Vecchio 于 1990 年在 ACI 發(fā)表的論文中采用了基于 MCFT 理論的割線剛度迭代法來(lái)進(jìn)行鋼筋混凝土的有限元分析。通過(guò)比較割線剛度迭代法與基于切線剛度的 NewtonRaphson 法，可以得出以下幾點(diǎn)結(jié)論： 1. 在結(jié)構(gòu)的應(yīng)力 應(yīng)變曲線進(jìn)入到下降段后，對(duì)于切線剛度法，剛度會(huì)出現(xiàn)負(fù)值，從而給計(jì)算處理帶來(lái)難度。而割線剛度迭代法則沒(méi)有這種負(fù)剛度問(wèn)題，但是割線剛度迭代法會(huì)在卸載曲線穿越坐標(biāo)軸時(shí)出現(xiàn)數(shù)值問(wèn)題，這是由于此時(shí)割線剛度極小造成的。 2. 以往采用切線剛度法的程序，在更新應(yīng)變時(shí)，對(duì)于一個(gè)單元采用了 4 個(gè)高斯積分點(diǎn)，每次迭代得到了一個(gè)新的結(jié)構(gòu)位移，就需要在 4 個(gè)積分點(diǎn)上分別求出它們的應(yīng)變，再以此計(jì)算出 4 個(gè)新的割線剛度，用以判斷收斂與否。而 F. J. Vecchio 采用的割線剛度迭代法只涉及了 1 個(gè)積分點(diǎn)，即中心積分點(diǎn)。 根據(jù)附錄 1 可知，單元應(yīng)變和位移的關(guān)系為： jpimx i j m pNNNNu u u u ux x x x x? ?????? ? ? ? ?? ? ? ? ? jpimy i j m pNNNNv v v v vy y y y y? ?????? ? ? ? ?? ? ? ? ? j p j pi m i mx y i j m p i j m pN N N NN N N Nuv u u u u v v v vy x y y y y x x x x? ? ? ? ?? ? ? ???? ? ? ? ? ? ? ? ? ?? ? ? ? ? ? ? ? ? ?其中： 1 (1 )(1 )4i xyN ab? ? ? 1 (1 )(1 )4j xyN ab? ? ? 1 (1 )(1 )4m xyN ab? ? ? 1 (1 )(1 )4p xyN ab? ? ? 令 2ax? ， 2by? ，則 可以得到割線剛度迭代法中心積分點(diǎn)的應(yīng)變計(jì)算式為： 2j m i px u u u ua? ? ? ?? 2m p i jy v v v vb? ? ? ?? 22m p i j j m i pxy u u u u v v v vba? ? ? ? ? ? ??? 每次循環(huán)迭代出的新位移，僅需要換算成中心積分點(diǎn)上的應(yīng)變，通過(guò)一次計(jì)算便可以得出割線剛度用于收斂判斷。相比于切線剛度法，割線剛度迭代法簡(jiǎn)化了計(jì)算過(guò)程。 3. 切線剛度法的過(guò)程，需要通過(guò)外荷載節(jié)點(diǎn)力和結(jié)構(gòu)恢復(fù)力求解出一個(gè)殘余力，再以此計(jì)算出增量位移來(lái)判斷是否滿足收斂條件，若不收斂，便更新結(jié)構(gòu)位移并重新計(jì)算結(jié)構(gòu)的切線剛度矩陣。對(duì)于割線剛度迭代法，并不需要求解殘余力，而是直接利用外荷載節(jié)點(diǎn)力和結(jié)構(gòu)割線剛度來(lái)求解出一個(gè)全量位移，之后通過(guò)計(jì)算出的應(yīng)變來(lái)重新計(jì)算結(jié)構(gòu)的割線剛度矩陣。二者 相比起來(lái)，割線剛度迭代法不需要啟用一個(gè)新的數(shù)組來(lái)儲(chǔ)存殘余力值，降低了對(duì)程序的存儲(chǔ)空間要求并提高了程序的運(yùn)算速度。 從上面三點(diǎn)可以看出，切線剛度法和割線剛度迭代法各有特點(diǎn)和不足，在某些方面，割線剛度迭代法還體現(xiàn)出了一定的優(yōu)勢(shì)： 4. 割線剛度迭代法相比于切線剛度法的收斂準(zhǔn)則更為簡(jiǎn)單。 5. 割線剛度迭代法相比于切線剛度法的有限元模型更為簡(jiǎn)單。 6. 割線剛度迭代法的穩(wěn)定性較好。 所以對(duì)割線剛度迭代法進(jìn)行研究佐證是十分有必要的。 算法流程 以圖 所示的單元為例， Vecchio 在其論文中所述的割線剛度迭代法的具體計(jì)算 過(guò)程如下： 圖 單元尺寸 Element dimension 第 1 步：輸入結(jié)構(gòu)和材料的屬性數(shù)據(jù) 第 2 步：計(jì)算外荷載節(jié)點(diǎn)力 ??F 第 3 步：計(jì)算材料的剛度矩陣 ： 1) 確定初始應(yīng)變： x? ， y? ， xy? ； 2) 計(jì)算主方向應(yīng)變 1c? ， 2c? ， sx? ， sy? ，以及主拉應(yīng)變角 c? ； 3) 根據(jù)混凝土的本構(gòu)關(guān)系，確定混凝土主拉和主壓應(yīng)力 1cf ， 2cf 4) 根據(jù)鋼筋的本構(gòu)關(guān)系，確定鋼筋的主拉和主壓應(yīng)力 sxf ， syf ； 5) 計(jì)算混凝土和鋼筋的割線剛度元素值： 222cccfE ?? 111cccfE ?? 12ccccEEG EE? ? sysysyfE ?? sxsxsxfE ?? 第 4 步： 計(jì)算混凝土和鋼筋的材料割線剛度矩陣 1) 計(jì)算混凝土主方向上的材料割線剛度矩陣 ? ?39。cD ： ? ? 239。1000000cccEDEG??????? (32) 2) 計(jì)算轉(zhuǎn)換矩陣 ??cT ? ?2222c o s sin c o s sinsin c o s c o s sin2 c o s sin 2 c o s sin c o s sinT? ? ? ?? ? ? ?? ? ? ? ? ??????????? 3) 計(jì)算整體坐標(biāo)系（ XY 坐標(biāo)）下的混凝土割線剛度矩陣 ? ? ? ? ? ? ? ?39。TccD T D T? (33) 4) 計(jì)算整體坐標(biāo)系 (XY 坐標(biāo)系 )下的鋼筋割線剛度矩陣： ? ? 00000 0 0x sxs y syEDE????????? (34) 第 5 步：計(jì)算組合材料剛度 ??D ： ? ? ? ? ? ?csD D D?? (35) 第 6 步：計(jì)算單元?jiǎng)偠染仃???k ： ? ? ? ? ? ? ? ?Tk B D B d x d y d t? ??? (36) 第 7 步 ： 計(jì)算單元預(yù)應(yīng)變產(chǎn)生的位移 00,csrr? ? ? ?? ? ? ? 第 8 步：計(jì)算單元預(yù)應(yīng)變產(chǎn)生的節(jié)點(diǎn)力 *F???? 第 9 步：計(jì)算單元節(jié)點(diǎn)力： 39。*F F F? ? ? ? ? ???? ? ? ? ? ? 第 10 步：計(jì)算結(jié)構(gòu)總剛矩陣 ??K ： ? ? ? ?Kk?? (37) 第 11 步：計(jì)算結(jié)構(gòu)總剛矩陣的逆矩陣 ? ?1K? ： 第 12 步：計(jì)算節(jié)點(diǎn)位移 ? ? ? ? 1 39。r K F? ??? ?? (38) 第 13