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