新保存結(jié)構(gòu)力學(xué)課設(shè)

目錄結(jié)構(gòu)力學(xué)課程設(shè)計任務(wù)書 2第一章．緒論 7（1）matlab簡介 7（2）matlab優(yōu)勢特點 7（3）矩陣位移法簡介 10（4）矩陣位移法的特點 10（5）矩陣位移法的基本原理 10（6）矩陣位移法計算桿端力的步驟 11第二章．Matlab解法 12步驟一： 13步驟二： 13步驟三： 13步驟四： 14步驟五： 15步驟六： 16結(jié)論 19附錄程序 19程序 19子程序： 23第三章．組內(nèi)成員其他解法展示與對比 27（1）Ansys有限元解法： 27（2）矩陣位移法手算解法： 33第四章．用Matlab完成所選矩陣方程題目計算 38第五章．課程設(shè)計的心得體會。 40參考文獻 41結(jié)構(gòu)力學(xué)課程設(shè)計任務(wù)書設(shè)計題目結(jié)構(gòu)力學(xué)課程設(shè)計設(shè)計內(nèi)容查閱資料并學(xué)習(xí)相關(guān)軟件（Fortran或Matlab），對矩陣計算方程計算有所掌握，并完成所選矩陣方程題目計算；利用Fortran或Matlab軟件，完成對所選結(jié)構(gòu)力學(xué)題目的節(jié)點力及節(jié)點內(nèi)力計算工作；查閱資料并學(xué)習(xí)相關(guān)軟件（Ansys或Abaqus，Patran），對相關(guān)有限元軟件的基本分析有所掌握，并完成所選結(jié)構(gòu)力學(xué)題目的有限元分析計算；利用所學(xué)過的結(jié)構(gòu)力學(xué)知識，完成對所選結(jié)構(gòu)力學(xué)題目的手工計算工作，并完成相關(guān)內(nèi)力圖的繪制；完成對2）、3），4）三種計算結(jié)果進行對比分析。設(shè)計要求給出Fortran或Matlab矩陣計算程序的代碼，并給出詳細注釋，并把計算結(jié)果抓圖給出；給出Fortran或Matlab結(jié)構(gòu)力學(xué)計算程序的代碼，并給出詳細注釋，并把計算結(jié)果抓圖給出；給出有限元軟件計算結(jié)果（Mises應(yīng)力、應(yīng)變、位移、節(jié)點力）云圖結(jié)果；完成題目的手工計算，并輸入電腦，相關(guān)力學(xué)圖的繪制可用Office的Visio組件或Autocad繪制；給出三種計算結(jié)果的對比（Excel曲線圖或Origin）曲線圖分析。注意事項時間：2016學(xué)年第2學(xué)期校歷第17~18周（2016年6月地點：所有任務(wù)均可自行選擇完成地點；答疑地點為08A109；指導(dǎo)教師王毅、代君、王曉璐、趙輝、莫振偉第17周周一周二周三周四周五教師王毅代君王曉璐趙輝莫振偉第18周周一周二周三周四周五教師王毅代君王曉璐趙輝莫振偉題目選擇每位學(xué)生需要在矩陣方程計算題目及結(jié)構(gòu)力學(xué)計算題目中各選一道題目來完成；學(xué)生需要提交的課程設(shè)計報告應(yīng)包括：前附《課程設(shè)計任務(wù)書》、設(shè)計要求所涉及的五項內(nèi)容，參考文獻；成績評定及答辯由于工作量較大，最終提交任務(wù)完成的時間節(jié)點為2016學(xué)年第2學(xué)期校歷第19周周四（即2016年7月7日）。答辯時間暫定為2016年7月8日，成績分為優(yōu)秀、良好、中等、及格、不及格五等，具體時間另行通知。課程設(shè)計題目一）矩陣方程1.利用全選主元的高斯約當(dāng)（Gauss-Joadan）消去法求解如下方程組，并給出詳細的程序注解和說明：利用追趕法求解如下方程組，并給出詳細的程序注解和說明。利用全選主元的高斯約當(dāng)（Gauss-Joadan）消去法如下求解大型稀疏矩陣的大型方程組，并給出詳細注解及說明。結(jié)構(gòu)力學(xué)試求解圖示平面桁架各桿之軸力圖，已知各材料性能及截面面積相同，。（注：在有限元分析中，桁架桿的模擬只能選擇Ansys的Link單元）。試求解圖示平面剛架內(nèi)力圖（軸力圖、剪力圖和彎矩圖），已知各材料性能及截面面積相同，，泊松比。3.試求解圖示平面剛架內(nèi)力圖（軸力圖、剪力圖和彎矩圖），已知各桿。工作計劃所在小組成員及分工情況姓名張世秋，楊國，朱春暉顧凱強，袁介周劉永孫子楊，申正偉李虎，楊尚輝任務(wù)Matlab有限元分析Ansys有限元分析算法對比，誤差分析矩陣位移法手算Matlab列主元高斯消元法解矩陣方程我負責(zé)的是Matlab有限元分析，我們組分工的情況，如上所示。我一開始先查相關(guān)的資料，了解matlab工作原理，查閱相關(guān)算例，仿照書上的算例，編制相關(guān)的程序，嘗試計算。但是程序編制并不是一蹴而就，在編制的過程中，出了很多的相關(guān)問題，找不到問題的所在，最終我聽從老師的意見，然后終于查到問題所在，順利解決問題。課程負責(zé)人簽名：指導(dǎo)教師簽名：年月日年月日第一章．緒論（1）matlab簡介MATLAB[1]是美國MathWorks公司出品的商業(yè)數(shù)學(xué)軟件，用于算法開發(fā)、數(shù)據(jù)可視化、數(shù)據(jù)分析以及數(shù)值計算的高級技術(shù)計算語言和交互式環(huán)境，主要包括MATLAB和Simulink兩大部分。MATLAB是matrix&laboratory兩個詞的組合，意為矩陣工廠（矩陣實驗室）。是由美國mathworks公司發(fā)布的主要面對科學(xué)計算、可視化以及交互式程序設(shè)計的高科技計算環(huán)境。它將數(shù)值分析、矩陣計算、科學(xué)數(shù)據(jù)可視化以及非線性動態(tài)系統(tǒng)的建模和仿真等諸多強大功能集成在一個易于使用的視窗環(huán)境中，為科學(xué)研究、工程設(shè)計以及必須進行有效數(shù)值計算的眾多科學(xué)領(lǐng)域提供了一種全面的解決方案，并在很大程度上擺脫了傳統(tǒng)非交互式程序設(shè)計語言（如C、Fortran）的編輯模式，代表了當(dāng)今國際科學(xué)計算軟件的先進水平。MATLAB和Mathematica、Maple并稱為三大數(shù)學(xué)軟件。它在數(shù)學(xué)類科技應(yīng)用軟件中在數(shù)值計算方面首屈一指。MATLAB可以進行矩陣運算、繪制函數(shù)和數(shù)據(jù)、實現(xiàn)算法、創(chuàng)建用戶界面、連接其他編程語言的程序等，主要應(yīng)用于工程計算、控制設(shè)計、信號處理與通訊、圖像處理、信號檢測、金融建模設(shè)計與分析等領(lǐng)域。MATLAB的基本數(shù)據(jù)單位是矩陣，它的指令表達式與數(shù)學(xué)、工程中常用的形式十分相似，故用MATLAB來解算問題要比用C，F(xiàn)ORTRAN等語言完成相同的事情簡捷得多，并且MATLAB也吸收了像Maple等軟件的優(yōu)點，使MATLAB成為一個強大的數(shù)學(xué)軟件。在新的版本中也加入了對C，F(xiàn)ORTRAN，C++，JAVA的支持。（2）matlab優(yōu)勢特點1)高效的數(shù)值計算及符號計算功能，能使用戶從繁雜的數(shù)學(xué)運算分析中解脫出來；2)具有完備的圖形處理功能，實現(xiàn)計算結(jié)果和編程的可視化；3)友好的用戶界面及接近數(shù)學(xué)表達式的自然化語言，使學(xué)者易于學(xué)習(xí)和掌握；4)功能豐富的應(yīng)用工具箱(如信號處理工具箱、通信工具箱等)，為用戶提供了大量方便實用的處理工具。編程環(huán)境MATLAB由一系列工具組成。這些工具方便用戶使用MATLAB的函數(shù)和文件，其中許多工具采用的是圖形用戶界面。包括MATLAB桌面和命令窗口、歷史命令窗口、編輯器和調(diào)試器、路徑搜索和用于用戶瀏覽幫助、工作空間、文件的瀏覽器。隨著MATLAB的商業(yè)化以及軟件本身的不斷升級，MATLAB的用戶界面也越來越精致，更加接近Windows的標(biāo)準界面，人機交互性更強，操作更簡單。而且新版本的MATLAB提供了完整的聯(lián)機查詢、幫助系統(tǒng)，極大的方便了用戶的使用。簡單的編程環(huán)境提供了比較完備的調(diào)試系統(tǒng)，程序不必經(jīng)過編譯就可以直接運行，而且能夠及時地報告出現(xiàn)的錯誤及進行出錯原因分析。簡單易用Matlab是一個高級的矩陣/陣列語言，它包含控制語句、函數(shù)、數(shù)據(jù)結(jié)構(gòu)、輸入和輸出和面向?qū)ο缶幊烫攸c。用戶可以在命令窗口中將輸入語句與執(zhí)行命令同步，也可以先編寫好一個較大的復(fù)雜的應(yīng)用程序（M文件）后再一起運行。新版本的MATLAB語言是基于最為流行的C++語言基礎(chǔ)上的，因此語法特征與C++語言極為相似，而且更加簡單，更加符合科技人員對數(shù)學(xué)表達式的書寫格式。使之更利于非計算機專業(yè)的科技人員使用。而且這種語言可移植性好、可拓展性極強，這也是MATLAB能夠深入到科學(xué)研究及工程計算各個領(lǐng)域的重要原因。強處理能力MATLAB是一個包含大量計算算法的集合。其擁有600多個工程中要用到的數(shù)學(xué)運算函數(shù)，可以方便的實現(xiàn)用戶所需的各種計算功能。函數(shù)中所使用的算法都是科研和工程計算中的最新研究成果，而且經(jīng)過了各種優(yōu)化和容錯處理。在通常情況下，可以用它來代替底層編程語言，如C和C++。在計算要求相同的情況下，使用MATLAB的編程工作量會大大減少。MATLAB的這些函數(shù)集包括從最簡單最基本的函數(shù)到諸如矩陣，特征向量、快速傅立葉變換的復(fù)雜函數(shù)。函數(shù)所能解決的問題其大致包括矩陣運算和線性方程組的求解、微分方程及偏微分方程的組的求解、符號運算、傅立葉變換和數(shù)據(jù)的統(tǒng)計分析、工程中的優(yōu)化問題、稀疏矩陣運算、復(fù)數(shù)的各種運算、三角函數(shù)和其他初等數(shù)學(xué)運算、多維數(shù)組操作以及建模動態(tài)仿真等。圖形處理MATLAB自產(chǎn)生之日起就具有方便的數(shù)據(jù)可視化功能，以將向量和矩陣用圖形表現(xiàn)出來，并且可以對圖形進行標(biāo)注和打印。高層次的作圖包括二維和三維的可視化、圖象處理、動畫和表達式作圖。可用于科學(xué)計算和工程繪圖。新版本的MATLAB對整個圖形處理功能作了很大的改進和完善，使它不僅在一般數(shù)據(jù)可視化軟件都具有的功能（例如二維曲線和三維曲面的繪制和處理等）方面更加完善，而且對于一些其他軟件所沒有的功能（例如圖形的光照處理、色度處理以及四維數(shù)據(jù)的表現(xiàn)等），MATLAB同樣表現(xiàn)了出色的處理能力。同時對一些特殊的可視化要求，例如圖形對話等，MATLAB也有相應(yīng)的功能函數(shù)，保證了用戶不同層次的要求。另外新版本的MATLAB還著重在圖形用戶界面（GUI）的制作上作了很大的改善，對這方面有特殊要求的用戶也可以得到滿足。MATLAB對許多專門的領(lǐng)域都開發(fā)了功能強大的模塊集和工具箱。一般來說，它們都是由特定領(lǐng)域的專家開發(fā)的，用戶可以直接使用工具箱學(xué)習(xí)、應(yīng)用和評估不同的方法而不需要自己編寫代碼。領(lǐng)域，諸如數(shù)據(jù)采集、數(shù)據(jù)庫接口、概率統(tǒng)計、樣條擬合、優(yōu)化算法、偏微分方程求解、神經(jīng)網(wǎng)絡(luò)、小波分析、信號處理、圖像處理、系統(tǒng)辨識、控制系統(tǒng)設(shè)計、LMI控制、魯棒控制、模型預(yù)測、模糊邏輯、金融分析、地圖工具、非線性控制設(shè)計、實時快速原型及半物理仿真、嵌入式系統(tǒng)開發(fā)、定點仿真、DSP與通訊、電力系統(tǒng)仿真等，都在工具箱（Toolbox）家族中有了自己的一席之地。程序接口新版本的MATLAB可以利用MATLAB編譯器和C/C++數(shù)學(xué)庫和圖形庫，將自己的MATLAB程序自動轉(zhuǎn)換為獨立于MATLAB運行的C和C++代碼。允許用戶編寫可以和MATLAB進行交互的C或C++語言程序。另外，MATLAB網(wǎng)頁服務(wù)程序還容許在Web應(yīng)用中使用自己的MATLAB數(shù)學(xué)和圖形程序。MATLAB的一個重要特色就是具有一套程序擴展系統(tǒng)和一組稱之為工具箱的特殊應(yīng)用子程序。工具箱是MATLAB函數(shù)的子程序庫，每一個工具箱都是為某一類學(xué)科專業(yè)和應(yīng)用而定制的，主要包括信號處理、控制系統(tǒng)、神經(jīng)網(wǎng)絡(luò)、模糊邏輯、小波分析和系統(tǒng)仿真等方面的應(yīng)用。應(yīng)用軟件開發(fā)在開發(fā)環(huán)境中，使用戶更方便地控制多個文件和圖形窗口；在編程方面支持了函數(shù)嵌套，有條件中斷等；在圖形化方面，有了更強大的圖形標(biāo)注和處理功能，包括對性對起連接注釋等；在輸入輸出方面，可以直接向Excel和HDF5進行連接。（3）矩陣位移法簡介在結(jié)構(gòu)力學(xué)的計算中，通過采用對結(jié)點位移作為基本未知量，進而通過矩陣的形式對各基本參數(shù)進行組織，編排，求出未知量的方法，稱為矩陣位移法。（4）矩陣位移法的特點矩陣數(shù)學(xué)表達力強矩陣數(shù)學(xué)表達力強，運算簡潔方便并且適于計算機組織運算，是用計算機進行結(jié)構(gòu)數(shù)值分析的最強有力的數(shù)學(xué)工具。矩陣位移法與結(jié)構(gòu)力學(xué)的力法和位移法相對應(yīng)，也就是結(jié)構(gòu)的矩陣分析方法。矩陣位移法方便編制程序矩陣位移法便于編制程序，因而在工程界得到廣泛應(yīng)用。矩陣位移法并不因采用矩陣數(shù)學(xué)的描述手段，而改變位移法的基本原理。它與位移法的區(qū)別僅僅在于表達形式不同。（5）矩陣位移法的基本原理按位移法的基本原理運用矩陣計算內(nèi)力和位移的方法。是結(jié)構(gòu)矩陣分析方法中的一種，其基本未知數(shù)是結(jié)點位移，由于矩陣位移法較矩陣力法更適宜編制通用的計算程序，因而得到了更為廣泛的應(yīng)用。結(jié)構(gòu)矩陣分析方法首先把結(jié)構(gòu)離散成有限數(shù)目的單元，然后再合成為原結(jié)構(gòu)，因而也屬于有限元法。矩陣位移法常用的單元形式為一直桿。對于曲桿，如拱結(jié)構(gòu)，雖然也可取曲桿作為單元,但單元分析較煩,為簡化起見，可將它化成折線來處理，每一直線段作為一單元。當(dāng)單元承受非結(jié)點荷載時，可用等效結(jié)點荷載代替。其方法是將單元間的分界結(jié)點作為固端求出固端反力，然后反其向作用在結(jié)點上。（6）矩陣位移法計算桿端力的步驟矩陣位移法計算桿端力的步驟為：①劃分單元，求出等效結(jié)點荷載；②求單元剛度矩k(i),并轉(zhuǎn)換為整體坐標(biāo)的單元剛度矩陣;③由剛度法求出整體剛度矩陣K；④求出Krr；⑤求出結(jié)點位ur，再求出桿端s，實際桿端力應(yīng)再疊sf，即確定。第二章．Matlab解法為了便于計算和書寫對圖形做如下處理。其中p=20KN,L=2m。因為結(jié)構(gòu)是對稱的，載荷也是對稱的，為了簡化計算，可以利用對稱性只取結(jié)構(gòu)的一半進行計算，如下圖所示，因為節(jié)點2,4，沒有水平位移，所以各加了一根水平鉸鏈，又因為4,2，桿在對稱軸上，所以它的橫截面積取A/2，各節(jié)點位移的編號如圖所示，各基本數(shù)據(jù)的計算如下表所示。步驟一：離散化域，域分為8個單元5個節(jié)點，matlab中采用的單位是kN和m，下表給出了本題的單元連通性。步驟二：寫出單元剛度矩陣，通過matlab函數(shù)PlaneTrussElementStiffness，得到8個單元的剛度矩陣k1,k2,k3,……,k8。每個矩陣都是4*4矩陣。步驟三：集成整體剛度矩陣。由于該結(jié)構(gòu)一共有5個節(jié)點，所以整體剛度矩陣是10*10矩陣。因此為了得到整體剛度矩陣K,我們首先要生成一個10*10的零矩陣。由于該結(jié)構(gòu)有8個單元，8次調(diào)用matlab的PlaneTrussAssemble函數(shù)就可以得到整體剛度矩陣K。每次對該函數(shù)調(diào)用都集成一個單元。步驟四：引入邊界條件。用上一步得到的整體剛度矩陣，可以得到該結(jié)構(gòu)的矩陣方程。F1xF1yF1xF1yF2xF2yF3xF3yF4xF4yF5xF5yU1xU1yU2xU2yU3xU3yU4xU4yU5xU5y[K]* =0-20F2x000F4xF4yF5x0-20F2x000F4xF4yF5xF5yU1xU1y0U2yU3xU3y0000[K]* =步驟五：用分解（手動），和高斯消元法（應(yīng)用matlab）求解方程組。首先對方程組進行分解，提取整體剛度矩陣K的第1，2，4，5，6行及列的交集組成子矩陣[K1]。如下圖所示：從而得到：0-200-20000U1xU1yU2yU3xU3y[K]* =上述方程組的解用如matlab命令可以得到。解得u=-0.0000-0.0024-0.00080.0002-0.0020現(xiàn)在可以得知各節(jié)點位移。步驟六：在這一步中，首先我們建立結(jié)構(gòu)節(jié)點位移矢量UU=[0;-0.0024;0;-0.0008;0.0002;-0.0020;0;0;0;0];然后我們建立單元節(jié)點位移矢量u1,u2,u3,u4,u5,u6,u7,u8。u1=[U(1);U(2);U(3);U(4)];u2=[U(3);U(4);U(7);U(8)];u3=[U(7);U(8);U(5);U(6)];u4=[U(5);U(6);U(9);U(10)];u5=[U(9);U(10);U(1);U(2)];u6=[U(1);U(2);U(5);U(6)];u7=[U(3);U(4);U(5);U(6)];u8=[U(1);U(2);U(7);U(8)];然后調(diào)用matlab的PlaneTrussElementStress函數(shù)計算出單元應(yīng)力sigma1,sigma2,sigma3,sigma4,sigma5,sigma6,sigma7,sigma8。sigma1=0sigma2=84000sigma3=21000sigma4=-21000sigma5=-1.2600e+05sigma6=4.2000e+04sigma7=-5.2500e+04sigma8=1.2600e+05結(jié)論于是得到所有桿件的軸向應(yīng)力：附錄程序所有程序及其子程序如下：程序E=210*10^6;A=90.7*10^-6;A2=0.5*A;L1=2;theta1=0;L2=2;theta2=90;L3=2;theta3=0;L4=2;theta4=0;L5=PlaneTrussElementLength(0,0,2,2);theta5=45;L6=2;theta6=90;L7=L5;theta7=45;L8=L5;theta8=135;k1=PlaneTrussElementStiffness(E,A,L1,theta1);k2=PlaneTrussElementStiffness(E,A2,L2,theta2);k3=PlaneTrussElementStiffness(E,A,L3,theta3);k4=PlaneTrussElementStiffness(E,A,L4,theta4);k5=PlaneTrussElementStiffness(E,A,L5,theta5);k6=PlaneTrussElementStiffness(E,A,L6,theta6);k7=PlaneTrussElementStiffness(E,A,L7,theta7);k8=PlaneTrussElementStiffness(E,A,L8,theta8);K=zeros(10,10);K=PlaneTrussAssemble(K,k1,1,2);K=PlaneTrussAssemble(K,k2,4,2);K=PlaneTrussAssemble(K,k3,3,4);K=PlaneTrussAssemble(K,k4,5,3);K=PlaneTrussAssemble(K,k5,5,1);K=PlaneTrussAssemble(K,k6,3,1);K=PlaneTrussAssemble(K,k7,3,2);K=PlaneTrussAssemble(K,k8,4,1);k=K([1,2,4,5,6],[1,2,4,5,6])f=[0;-20;0;0;0];u=k\fu=[0;-0.0024;-0.0008;0.0002;-0.0020];U=[0;-0.0024;0;-0.0008;0.0002;-0.0020;0;0;0;0];u1=[U(1);U(2);U(3);U(4)];u2=[U(3);U(4);U(7);U(8)];u3=[U(7);U(8);U(5);U(6)];u4=[U(5);U(6);U(9);U(10)];u5=[U(9);U(10);U(1);U(2)];u6=[U(1);U(2);U(5);U(6)];u7=[U(3);U(4);U(5);U(6)];u8=[U(1);U(2);U(7);U(8)];sigma1=PlaneTrussElementStress(E,L1,theta1,u1)sigma2=PlaneTrussElementStress(E,L2,theta2,u2)sigma3=PlaneTrussElementStress(E,L3,theta3,u3)sigma4=PlaneTrussElementStress(E,L4,theta4,u4)sigma5=PlaneTrussElementStress(E,L5,theta5,u5)sigma6=PlaneTrussElementStress(E,L6,theta6,u6)sigma7=PlaneTrussElementStress(E,L7,theta7,u7)sigma8=PlaneTrussElementStress(E,L8,theta8,u8)子程序：functiony=PlaneTrussElementStress(E,L,theta,u)%該函數(shù)根據(jù)彈性模量E,長度L,角度theta(單位是度)以及單位節(jié)點位移矢量u計算單元應(yīng)力。它返回單位應(yīng)力大小，返回值是一個標(biāo)量而不是矢量。x=theta*pi/180;C=cos(x);S=sin(x);y=E/L*[-C-SCS]*u;functiony=PlaneTrussInclinedSupport(T,i,alpha)%該函數(shù)根據(jù)傾斜支柱的節(jié)點號i,以及傾斜角alpha計算得到傾斜支柱的變換矩陣。它返回2n*2n的變換矩陣。x=alpha*pi/180;T(2i-1,2i-1)=cos(x);T(2i-1,2i)=sin(x);T(2i-1,2i-1)=-sin(x);T(2i,2i)=cos(x);y=T;functiony=PlaneTrussElementStiffness(E,A,L,theta)%該函數(shù)根據(jù)每個平面桁架元的彈性模量E，橫截面積A,長度L,%以及角度theta(單位是度)計算得到單元剛度矩陣，它返回4*4的單元剛度矩陣。x=theta*pi/180;C=cos(x);S=sin(x);y=E*A/L*[C*CC*S-C*C-C*S;C*SS*S-C*S-S*S-C*C-C*SC*CC*S;-C*S-S*SC*SS*S];functiony=PlaneTrussElementLength(x1,y1,x2,y2)%該函數(shù)根據(jù)給出的第一個坐標(biāo)（x1,y1）和第二個坐標(biāo)（x2,y2）計算返回單元長度。y=sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));functiony=PlaneTrussElementForce(E,A,L,theta,u)%該函數(shù)根據(jù)彈性模量E,長度L,角度theta(單位是度)以及單元節(jié)點位移矢量u計算單元節(jié)點力。它以標(biāo)量形式返回單元節(jié)點力。x=theat*pi/180;C=cos(x);S=sin(x);y=E*A/L*[-C-SCS]*u;functiony=PlaneTrussAssemble(K,k,i,j)%該函數(shù)將連接i和j的平面桁架元的剛度矩陣k集成到整體剛度矩陣K。%每集成一個單元，該函數(shù)都會返回2n*2n的整體剛度矩陣K.K(2*i-1,2*i-1)=K(2*i-1,2*i-1)+k(1,1);K(2*i-1,2*i)=K(2*i-1,2*i)+k(1,2);K(2*i-1,2*j-1)=K(2*i-1,2*j-1)+k(1,3);K(2*i-1,2*j)=K(2*i-1,2*j)+k(1,4);K(2*i,2*i-1)=K(2*i,2*i-1)+k(2,1);K(2*i,2*i)=K(2*i,2*i)+k(2,2);K(2*i,2*j-1)=K(2*i,2*j-1)+k(2,3);K(2*i,2*j)=K(2*i,2*j)+k(2,4);K(2*j-1,2*i-1)=K(2*j-1,2*i-1)+k(3,1);K(2*j-1,2*i)=K(2*j-1,2*i)+k(3,2);K(2*j-1,2*j-1)=K(2*j-1,2*j-1)+k(3,3);K(2*j-1,2*j)=K(2*j-1,2*j)+k(3,4);K(2*j,2*i-1)=K(2*j,2*i-1)+k(4,1);K(2*j,2*i)=K(2*j,2*i)+k(4,2);K(2*j,2*j-1)=K(2*j,2*j-1)+k(4,3);K(2*j,2*j)=K(2*j,2*j)+k(4,4);y=K;第三章．組內(nèi)成員其他解法展示與對比（1）Ansys有限元解法：第一步：進入ANSYS(設(shè)定工作目錄和工作文件)進入ANSYS菜單路徑“程序>ANSYS>ANSYS12.0”設(shè)置工作文件名菜單路徑“file>ChangeJobname”，彈出“ChangeJobname”對話框，輸入“hengjia”，單擊OK確定并關(guān)閉對話框。設(shè)置計算類型菜單路徑“ANSYSMainMenu:Preferences”，在彈出的對話框中選擇“Structural”，單擊OK確定并關(guān)閉對話框。選擇單元類型菜單路徑“ANSYSMainMenu:Preprocessor>ElementType>Add/Edit/Delete…”，在彈出“LibraryofElementTypes”對話框中按照參數(shù)選擇，單擊OK確定并關(guān)閉對話框。第二步：定義實常數(shù)菜單路徑“ANSYSMainMenu:Preprocessor>RealConstants…>Add/Edit/Delete”，在彈出的對話框中單擊“Add>OK”，彈出“RealConstant…”對話框，參數(shù)設(shè)置“AREA0.00325”，單擊OK確定并關(guān)閉對話框。第三步：定義材料參數(shù)菜單路徑“ANSYSMainMenu:Preprocessor>MaterialProps>MaterialModels”，在彈出的菜單中打開“Structural>Linear>Elastic>Isotropic”，彈出“LinearIsotropicMaterial…”對話框，并設(shè)置參數(shù)。第四步：定義梁的截面菜單路徑“ANSYSMainMenu:Preprocessor>Sections>Beam>CommonSections”，彈出“BeamTool”對話框，并設(shè)置，單擊OK確定關(guān)閉對話框。生成幾何模型，結(jié)果如圖所示。第五步：劃分網(wǎng)格，生成有限元模型菜單路徑“ANSYSMainMenu:Preprocessor>Meshing>MeshTool”，“MeshTool”對話框。單擊其中Lines旁邊的Set按鈕，打開“ElementSizeonPickedLines”選擇對話框，然后輸入?yún)?shù)：第六步：在mesh框中單擊PickAll按鈕，生成有限元模型如下第七步:模型施加約束和外載在1點處施加自由度為零的約束菜單路徑“ANSYSMainMenu:Solution>DefineLoads>Apply>Structural>Displayment>OnKeypoints”，彈出的對話框，選擇關(guān)鍵點1，單擊Apply彈出“ApplyU，ROTonKPs”對話框，選擇“ALLDOF”，單擊OK確定關(guān)閉對話框。第八步：在7點處施加Y方向的約束方法同上。設(shè)置重力加速度g第九步：分析計算菜單路徑“ANSYSMainMenu:Solution>Solve>CurrentLS”，在彈出的對話框中單擊OK，并關(guān)閉文字窗口。第十步：顯示變形圖菜單路徑“ANSYSMainMenu:GeneralPostproc>PlotResults>DeformedShape…”，彈出“PlotDeformedShape”對話框，設(shè)置“Def+Undeformed”單擊OK確定關(guān)閉對話框。第十一步：顯示位移云圖菜單路徑“ANSYSMainMenu:GeneralPostproc>PlotResults>NodalSolu”，彈出“ContourNodalSolutionData”對話框，并設(shè)置，單擊OK查看位移云圖。結(jié)果如下圖所示。第十二步：顯示軸力云圖節(jié)點力圖最終輸出軸力的大小如下所示：PRINTELEMENTTABLEITEMSPERELEMENT*****POST1ELEMENTTABLELISTING*****STATCURRENTELEMSMIS111826.52-11559.3-3653.145166.25-1826.56-11559.7-7306.18-1826.59-11559.101826.5115166.212-3653.113-11559.140.0000150.0000MINIMUMVALUESELEM2VALUE-11559.MAXIMUMVALUESELEM4VALUE5166.2將其轉(zhuǎn)化為如下表所示：（2）矩陣位移法手算解法：1.試求解圖示平面桁架各桿之軸力圖，已知各材料性能及截面面積相同，。為了便于計算和書寫對圖形做如下處理。其中p=20KN,L=2m。因為結(jié)構(gòu)是對稱的，載荷也是對稱的，為了簡化計算，可以利用對稱性只取結(jié)構(gòu)的一半進行計算，如下圖所示，因為節(jié)點2,4，沒有水平位移，所以各加了一根水平鉸鏈，又因為4,2，桿在對稱軸上，所以它的橫截面積取A/2，各節(jié)點位移的編號如圖所示，各基本數(shù)據(jù)的計算如下表所示。應(yīng)用式：UUiviujvjXijXijYijXjiYjiK計算各單元勁度矩陣如下：KKKKK用對號入座的方法組成結(jié)構(gòu)勁度矩陣。K列出結(jié)構(gòu)勁度方程如下：EF/L1.7071001.7071由上列方程解得未知結(jié)點位移：u1v1v2根據(jù)式Nij=EF/L[COSθ(uj-ui)+sinθ(vj-vi)]列表計算出如下表所示的桁架各桿的內(nèi)力。由于上圖所示的桁架是由原來的桁架截取的一半，所以桿42的真實內(nèi)力應(yīng)將表中的數(shù)值乘以2.這樣如下表所示桁架的各桿的內(nèi)力。（3）三種解法的對比與分析：三種解法的數(shù)據(jù)對比顯示于下列表格和折線圖中：從圖表中很容易得出結(jié)論，整體來看，三種計算方法得出的結(jié)果是一致的，但使用ansys有限元方法進行分析，得出的計算結(jié)果更加接近于手算得出的結(jié)果。而且在實際的練習(xí)中發(fā)現(xiàn)使用ansys有限元方法比較生動形象，更加易于操作和求解，省去了很多編程計算所不必要的麻煩。而使用matlab求解此次出現(xiàn)的誤差較大主要是在編程的過程中，為了減少不必要的麻煩和出錯的可能性，選擇的計算方法略為粗糙，所以得到的計算結(jié)果不甚精確。第四章．用Matlab完成所選矩陣方程題目計算3.利用全選主元的高斯約當(dāng)（Gauss-Joadan）消去法如下求解大型稀疏矩陣的大型方程組，并給出詳細注解及說明。Matlab計算程序如下所示：clearA=input('輸入系數(shù)矩陣A：');b=input('輸入b向量（按行向量）：');B=[Ab'];n=length(b);RA=rank(A);RB=rank(B);zhica=RB-RA;ifzhica>0,disp('請注意：因為RA~=RB，所以此方程組無解.\n')returnendifRA==RBifRA==nfprintf('請注意：因為RA=RB=%d，所以此方程組有唯一解.\n',n)X=zeros(n,1);forp=1:n-1t=find(abs(B(p:end,p))==max(abs(B(p:end,p))))+p-1;ifabs(B(t,p))~=abs(B(p,p))l=B(t,:);B(t,:)=B(p,:);B(p,:)=l;end%列主元判斷fork=p+1:nm=B(k,p)/B(p,p);B(k,p:n+1)=B(k,p:n+1)-m*B(p,p:n+1);endend%把方程組系數(shù)矩陣A化為同解的上三角矩陣b=B(1:n,n+1);A=B(1:n,1:n);X(n)=b(n)/A(n,n);forq=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);end%從xn至x1逐個求解上三角方程組elsedisp('請注意：因為RA=RB<n，所以此方程組有無窮多解.')returnendenddisp('方程組的解為：');X求解過程及運行結(jié)果如下所示：gao_si_quan_zhu_yuan_xiao_yuan_fa輸入系數(shù)矩陣A：[0,0,-1,0,0,0,2,0;0,6,0,0,0,-6,0,0;0,0,0,2,0,0,0,-4;3,0,0,0,-2,0,1,0;0,0,6,0,0,0,5,0;1,0,0,0,-3,0,0,2;0,4,0,-1,0,0,0,0;0,0,1,0,-1,0,0,-2]輸入b向量（按行向量）：[4,6,-8,-2,27,-9,2,-4]請注意：因為RA=RB=8，所以此方程組有唯一解.方程組的解為：X=1.000002.0000-2.00004.0000-1.00003.00001.000第五章．課程設(shè)計的心得體會。課程設(shè)計，一個展示自己通過三年的大學(xué)學(xué)習(xí)而獨立思考問題、解決問題的平臺