matlab連續(xù)梁程序的編制與使用

1、mat lab連續(xù)梁程序的 編制與使用 -CAL-FENGHAI-(2020YEAR-YICAI)INGBIAN 第三章 連續(xù)梁程序的編制與使用 結構力學中的矩陣位移法是隨著電子汁算機進 入結構力學領域中而產生的一種方法，而Matlab 語言正是進行矩陣運算的強大工具，因此，用 Matlab語言編寫結構力學程序有更大的優(yōu)越性。 本章將詳細介紹如何利用Matlab語言編制連續(xù)梁 結構的計算程序。 矩陣位移法的解題思路是將結構離散為單元 （桿件），建立單元桿端力與桿端位移之間的關 系一單元剛度方程；再將各單元集成為原結構， 在滿足變形連續(xù)條件和平衡條件時，建立整體剛 度方程；在邊界條件處理完畢后，

2、由整體剛度方 程解出節(jié)點位移，進而求出結構內力。 用矩陣位移法計算連續(xù)梁的步驟如下： 1）整理原始數據，如材料性質、荷載 條件、約束條件等，并進行編碼：單元編 碼、結點編碼、結點位移編碼、選取坐標 系。 2) 陣。 建立局部坐標系下的單元剛度矩 3) 陣。 建立整體坐標系下的單元剛度矩 4) 集成總剛。 5) 建立整體結構的等效節(jié)點荷載和總 荷載矩陣 6) 邊界條件處理。 7) 解方程，求出節(jié)點位移。 8) 求出各單元的桿端內力。 實際上，上述步驟也是編制Matlab程序的基本步驟，在求出計算結果后, 還可以利用Matlab的繪圖功能繪制結構圖、內力圖、變形圖等等。 圖3-1程療:流程圖 程序

3、說明 %矩陣位移法解連續(xù)梁主程序 %*4（4（*4（）|（*床*穴末床床床末床沐沐沐沐沐審水承璟璟才才廉床床床甲卡床床床床床末床沐床沐沐沐審水承璟璟才才廉床末床床床床床床床沐 功能：運用矩陣位移法解連續(xù)梁的基本原理編制的計算主程序。 基本思想：結點（結點位移）編碼默認為從左至右，從1開始順序進行; 桿端彎矩的方向默認為逆時針。 荷載類型：可計算結點荷載，每單元作用的跨中集中力和均布荷載。 說明：主程序的作用是通過賦值語句、讀取和寫入文件、函數調用等完成 算法的全過程，即實現程序流程圖的程序表達。 % 1程序準備 format short e%設定輸出類型 clear all%清除所有已定 義變

4、量 clc%清屏 說明： format short e 一設定計算過程中顯示在屏幕上的數字類型為短格式、科 學計數法； clear all一清除所有已定義變量，目的是在本程序的運行過程中， 不會發(fā)生變量名相同等可能使計算出錯的情況； clc一清屏，使屏幕在本程序運行開始時 % FPl=fopen(；卍)； FP2=fopen(l,/wtl); 說明： 2打開文件 %打開輸入數據文件存放初始數據 %打開輸出數據文件存放訃算結果 FPl=fopen（；rf）; 一打開已存在的輸入數據文件，且設置其為只讀格式，使 程序在執(zhí)行過程中不能改變輸入文件中的數值，并用文件句柄FP1來 FP2=fopen（/

5、wt）;打開輸出數據文件，該文件不存在時，通過此命 令創(chuàng)建新文件，該文件存在時則將原有內容全部刪除。該文件設置為可寫格 式，可在程序執(zhí)行過程中向輸出文件寫入數據。 3讀入程序控制信息 NELEM =fscanf（FPl；%dzl）; %單元個數（單元編碼總數） NPOIN =fscanf（FPl；%dzl）; %結點個數（結點編碼總數） NVFIX =fscanf（FPl；%d,l）; %約束個數（零位移總數） NFPOIN=fscanf（FPl/%d/l）; %結點荷載個數（作用在結點上集中力偶總 數） NFPRES=fscanf（FPl；%d,/l）； %非結點荷載數（作用在單元上分布荷載

6、總 數） YOUNG=fscanf（FPl/%f/l）; %彈性模量 說明： 從輸入文件FP1中讀入單元個數，結點個數，約束個數，結點荷載個數， 非結點荷載個數，彈性模量； 程序中彈性模量僅輸入了一個值，表明本程序僅能求解一種材料構成的結 構，如：鋼筋混凝土結構、鋼結構，不能求解鋼筋混凝土一鋼組合結構。 采用了命令fscanf,其中：吆出表示讀入整數格式；1表示讀取1個數。 fprintf（FP2八n結構初始數據nrV）; fprintf（FP2,*單元總數=%d結點總數=%d約束總數=%d n,NELEM,NPOIN,NVFIX）; fprintf（FP2,結點荷載個數=%d非結點荷載個數=

7、%d彈性模量=%12.5g n,NFPOIN,NFPRES,YOUNG）; 說明： 在輸出文件FP2中顯示輸入的控制信息，便于程序執(zhí)行完畢后，從輸出文 件中查找輸入錯誤。 采用了命令fprintf,其中：n為換行標志，表示換一行；%12.5g表示輸出 12位且有5位小數的實數。 括在引號中的信息，單元總數二d結點總數=%d約束總數二dH為輸岀 到FP2文件中時的格式，其后的變量表NELEM,NPOIN,NVFIX依次將其代表的數 值輸岀到4所對應的位置。 4調用讀取初始數據函數，生成結構信息 LNODS,COORD,FPOINFPRES, FIXED二 eleNITIALDATA (FPlFP

8、NELEMPOINFPOINFPRENVFIX); 5調用總剛計算函數，生成單剛并集成總剛 HK = ele_HK(NPOIN,NELEM 丄 NODS,COORD,YOUNG); 6調用荷載計算函數，山結點荷載與非結點荷載生成總荷載向量列陣 FORCE=. ele_FORCE(NPOIN,NFPOIN,FPOIN,NFPRES/FPRES/LNODS/COORD); 7調用邊界條件處理函數，總剛、總荷載進行邊界條件處理 HK,FORCE = ele_BOUNDARY(NVF%FIXED,HK,FORCE); 8方程求解 DISP=zeros(NPOIN/l)； %建立并初始化結構的結點位移列

9、矩陣 DISP=HKFORCE; %計算出結構所有的結點位移 9調用單元桿端力計算函數，求單元桿端力 FE = ele_MOMENTS(FP2/NPOIN,DISP,NELEM/LNODS/. COORD,YOUNG,NFPRES,FPRES); 10關閉輸出數據文件 fclose(FP2); %讀取初始數據函數ele_lNITIALDATA % 入口參數：FPbFP2,NELEM,NPOIN,NFPOIN,NFPRES,NVFIX % 出口參數：LNODS,COORD,FPOIN,FPRES,FIXED % function LNODS,COORD,FPOIN,FPRES,FIXED=. e

10、le_INITIALDATA(FPl/FP2/NELEM,NPOIN/NFPOIN/NFPRES,NVFIX) 說明： FPl文件句柄，指示輸入數據文件； FP2文件句柄，指示輸出數據文件； %讀取結構信息 LNODS=fscanf(FPl；%f,3,NELEM)*; 說明： 建立LNODS矩陣，該矩陣指出了每一單元的連接信息和慣性矩。 矩陣的每一行針對每一單元，共計NELEM；每一列相應為單元左結點號 (編碼)、(編碼)、慣性矩。 命令中，3,NELEM表示讀取NELEM行3列數據賦值給LNODS矩陣。 顯然，LNODS(U：3)依次表示i單元的左結點號、右結點號、慣性矩。 從這種定義可見，

11、每一單元的慣性矩均為常數，均為等截面直桿。 % fprintf(FP2/單元號左結點號右結點號慣性矩n); for i=l:NELEM fprintf(FP2/ %3d%3d %3d% nzi,LNODS(i,:); end 說明： 在輸出文件FP2中顯示輸入的單元連接等信息。 從1單元到NELEM單元進行循環(huán)，并依次輸出每一單元的連接等信息。 COORD=fscanf(FPl,%f,NPOIN)*; % 坐標：x 坐標(共計 NPOIN 組) 說明： 建立COORD矩陣，該矩陣用來存儲各結點x方向的坐標值。 連續(xù)梁結構各結點均分布在x軸，以1結點為起始點順序編碼，因此表征 各結點位置時僅需存

12、儲其x方向的坐標即可。 從FP1文件中讀取全部結點個數NPOIN的坐標值。 COORD(i)表示第i個結點的x坐標。 fprintf(FP2/結點定義數據:結點號x坐標n); for i=l:NPOIN fprintf(FP 乙，%3d%6.2fn,/i,COORD(i); end 說明： 在輸出文件FP2中顯示輸入的結點坐標信息。 從1結點到NPOIN結點進行循環(huán)，并依次輸出每一結點的坐標信息。 FPOIN=fscanf(FPl,%f,2,NFPOIN)； 說明： 建立FPOIN矩陣，該矩陣用來存儲直接作用在結點上的荷載信息。 從FP1文件讀取NFPOIN行2列數據，賦值給FPOIN矩陣。

13、每一行的第一列表示荷載作用的結點號；第二列表示荷載的數值大小。連 續(xù)梁結構每一結點僅有1個轉角位移，相應地直接作用在結點上的荷載為力 偶，以順時針轉動為正。 FPOIN(i,l)表示第i個直接結點荷載作用的結點號，FPOIN(i,2)表示笫i個直 接結點荷載的數值大小。 若控制數據NFPOIN等于零，則此行命令不執(zhí)行。 訐 NFPOIN0 fprintf(FP2j結點荷載數據：結點號M力偶n); for i=l:NFPOIN fprintf(FP2, %3d%3d%6.2fn/i, FPOIN(i,:); end end 說明： 首先判斷控制數據NFPOIN是否為零，若為零則不需輸出數據；若不

14、為 零，表明本算例有直接作用在結點上的荷載，需要在輸岀文件FP2中顯示輸入 的結點荷載信息。 從第1到第NFPOIN個直接結點荷載進行循環(huán)，并依次輸出每個直接結點 荷載的順序號、結點號和荷載大小。 FPRES=fscanf(FPl,%f,/3,NFPRES); 說明： 建立FPRES矩陣，該矩陣用來存儲非結點荷載信息。 從FP1文件讀取NFPRES行3列數據，賦值給FPRES矩陣。 每一行的第一列表示非結點荷載作用的單元號；第二列表示荷載類型；第 三列表示荷載的數值大小。連續(xù)梁結構的局部坐標系x軸以向右為正、y軸以 向上為正、轉動方向以順時針為正。 FPRES(U)表示第i個非結點荷載作用的單

15、元號，FPRES(i,2)表示第i個非結點 荷載的類型(類型定義在函數ele_FORCE中)，FPRES(i,3)表示第i個非結點荷 載的數值大小。 若控制數據NFPRES等于零，則此行命令不執(zhí)行。 訐 NFPRES0 fprintf(FP2/非結點荷載數據：荷載號單元號荷載類型荷載大小卅)； for i=l:NFPRES fprintf(FP2/ %3d %3d %6.2f %6.2fn,i/FPRES(i,:); end end 說明： 首先判斷控制數據NFPRES是否為零，若為零則不需輸出數據；若不為 零，表明本算例有非結點荷載，需要在輸出文件FP2中顯示輸入的非結點荷載 信息。 從第1

16、到第NFPRES個非結點荷載進行循環(huán)，并依次輸出每個非結點荷載 作用的荷載序號、單元號、荷載類型和荷載大小。 FIXED=fscanf(FPl/,%f/NVFIX); 說明： 建立FIXED矩陣，該矩陣用來存儲零位移對應的位移編碼信息(即有限制 轉動的約束使結點的轉角位移為零，對應的結點位移編碼。對于連續(xù)梁結構即 固定支座對應的結點編碼)。 從FP1文件讀取NVFIX行數據，賦值給FIXED列矩陣。 每一行表示零位移對應的結點位移編碼。 FIXED(i)表示笫i個零位移對應的結點位移編碼。 若控制數據NFPRES等于零，則此行命令不執(zhí)行。 % if NVFIX0 fprintf(FP2/零位移

17、約束數據：位移編號位移編號rf); fprintf(FP2； %3d%3dn,FIXED(:); end 說明： 首先判斷控制數據NVFIX是否為零，若為零則不需輸出數據；若不為零, 表明本算例有零位移，需要在輸出文件FP2中顯示零位移約束數據。 從笫1到第NVFIX個零位移進行循環(huán)，依次輸出每個零位移對應的結點 return 說明： return表明主程序終I匕。 計算總剛矩陣函數ele HK %入口參數：結點數、單元數、單元信息數組、結點坐標、彈性模量 %出口參數：整體剛度矩陣 function HK = ele_HK(NPOIN,NELEM丄NODS,COORD,YOUNG) HK=ze

18、ros(NPOIN,NPOIN); %生成總剛矩陣并清零 %生成單剛并組成總剛 for i=l:NELEM %對單元個數循環(huán) %調用生成局部單剛(局部坐標)函數 EK=ele_EK(i,LNODS,COORD,YOUNG); HK(i:i+l/i:i+l) = HK(i:i+l,i:i+l)+EK; end 說明： 對每一單元順序循環(huán)，依次計算每一單元的局部單元剛度矩陣。EK為一個 臨時變量，計算得到i單元的單剛EK后，通過下一條語句集成入總剛矩陣中。 在下一次循環(huán)中，則計算i+1單元的單剛，并賦值給EK,再進行集成。顯然， 隨著循環(huán)的不斷進行，EK對應著不同單元的單剛，值不斷發(fā)生改變。 單剛

19、到總剛的集成是利用了連續(xù)梁結構的特點，以分塊的形式放入總剛矩 陣中。例如：對于第i個單元，其四個單剛元素分別放入總剛矩陣的第(i,i)、 (ij+l) (i+l,i)、 (i+lj+l)四個位置,因此通過HK的下標(i:i+lfi:i+l)來表 示相應的元素位置。山于相鄰單元的單剛集成時，有元素疊加的情況，因此相 應的總剛矩陣元素的疊加采用了 HKfiu+lJJ+l) = HKfiii+lJJ+ll+EK的表達形式。 return 訃算單元剛度矩陣函數ele EK %入口參數：單元號、單元信息、結點坐標、彈性模量 %出口參數：局部單元剛度矩陣 說明： 在ele_HK函數中調用本函數時，彈性模量

20、的入口參數符號為YOUNG,而 在本函數的輸入參數表中，符號為E。函數調用時可以采用此種用法。 % function EK=ele_EK(i,LNODS/COORD,E) % NL=LNODS(i,l)； %左結點號 NR=LNODS(i/2);%右結點號 L=COORD(NR)-COORD(NL); %單元長度 l=LNODS(i?3)；%慣性矩 i=E*l/L; 說明： 利用已知信息，計算每一單元的線剛度。 LNODS(U)為i單元左結點的結點號賦值給NL： LNODS(i,2)為i單元右結點的 結點號賦值給NR。 COORD(NR)為第NR個結點的x坐標；COORD(NL)為第NL個結點

21、的x坐 標；二者的差為單元的長度。 LNODS(U)為i單元的慣性矩。 %生成單剛(局部坐標)右手坐標系 EK = 4*i 2*i; 2*i 4*i; return 力璟*末溶卡床末床農床床床末沐沐沐水沐審沐* %生成總荷載向量函數ele.FORCE O)K)|c)K*4c*4c*)c)ic)ic)ic)|c4cHc4c|cHc*|c4c * %入口參數：結點數，結點荷載個數，結點荷載信息，非結點荷載個數, 非結點荷載信息，單元信息，結點坐標 %出口參數：單元固端力一左右兩端的桿端彎矩 function FORCE二 ele_FORCE(NPOIN,NFPOIN,FPOIN,NFPRES FP

22、RES,LNODS,COORD) % FORCE=zeros(NPOIN,l); %生成總荷載向量并清零 for i=l:NFPOIN%對結點荷載個數進行循環(huán) FORCE(FPOIN(i,l)=FORCE(FPOIN(i/l)+FPOIN(i,2); % 取出結點荷載 end 說明： 建立并初始化總荷載向量FORCE,其行數與總剛矩陣、總位移列矩陣相 同，即與總結點數NPOIN相同(對于連續(xù)梁結構)。 對結點荷載個數進行循環(huán)，將每個結點荷載的數值放入所在結點號對應的 總荷載向量FORCE對應的行處 FPOIN(iJ)為第i個結點荷載所作用的結點號；FPOIN(i,2)為第i個結點荷載 的大?。?/p>

23、FORCE(FPOIN(i,l)為總荷向量對應結點號處的荷載。 % for i=l:NFPRES %對非結點荷載個數進行循環(huán) FO=ele_FPRES(i,FPRES丄NODS,COORD); % 計算單元固端力 %對單元局部桿端力要進行坐標轉換 ele=FPRES(i,l); %取荷載所在的單元號 NL=LNODS(ele,l)； NR=LNODS(ele,2); %單元的左右結點號 FORCE(NL)=FORCE(NL)-FO(1); %將固端力變成等效結點荷載 FORCE(NR)=FORCE(NR)-FO(2); %固端力與等效荷載符號相反 end 說明： 訃算山非結點荷載引起的等效結點

24、荷載，并集成入已有的總荷向量中。 對非結點荷載個數進行循環(huán)，即有兒個非結點荷載循環(huán)兒次，若沒有非結 點荷載，則本循環(huán)語句不執(zhí)行。 調用計算非結點荷載引起的單元固端力函數，訃算出非結點荷載所在單元 的兩端固端力。 FPRES(U)為第i個非結點荷載所在的單元號賦值給ele； LNODS(eleJ)為ele 單元左結點的結點號賦值給NL； LNODS(ele,2)為ele單元右結點的結點號賦值給 NRo FORCE(NL). FORCE(NR)分別為單元的左、右桿端的桿端力；-F0、-F0 為非結點荷載引起的左、右端固端力的負值，即等效結點荷載。 FORCE(NL)=FORCE(NL)-FO為直接

25、作用在單元左端結點上的荷載與左端等效 結點荷載的疊加；FORCE(NR)=FORCE(NR)-FO(2)為直接作用在單元右端結點上的 荷載與右端等效結點荷載的疊加。 return %計算非結點荷載引起的單元固端力函數ele_FPRES % 符號規(guī)定：正方向為：X向右Y向上M順時針 %* * * *璟床4(末末床*床4(床穴* * * *沐沐沐沐* He沐* * *璟* *床末床*床璟*末末末床*床床*末床* * *沐沐沐沐審*璟承末*床床末末末 %入口參數：荷載序號，非結點荷載信息，單元信息，結點坐標 %出口參數：單元固端力(左右兩端固端彎矩) function FO=ele_FPRES(i,

26、FPRES/LNODSzCOORD) ele=FPRES(izl)； %取荷載所在的單元號 G=FPRES(i,3); %單元荷載大小 NL=LNODS(ele,l); NR=LNODS(ele,2); % 單元的左右結點號 L=COORD(NR)-COORD(NL); % 單元長度 說明： 將已知非結點荷載信息分別賦值給相應的變量。 FPRES(U)為第i個非結點荷載作用的單元號；G=FPRES(i,3)為第i個非結點 荷載的大??； LNODS(ele,l)為ele單元左結點的結點號賦值給NL； LNODS(ele,2)為ele單元 右結點的結點號賦值給NR。 COORD(NR)為第NR個結

27、點的x坐標；COORD(NL)為第NL個結點的x坐 標；二者的差為單元的長度。 F0=0;0; %單元固端彎矩清零 switch FPRES(i,2) case 1%橫向均布荷載 FO(1)=-G*LA2/; F0(2)=G*LA2/; case 2%橫向集中力 F0(l)=-G*L/8; F0(2)=G*L/8; end 說明： 根據荷載類型，求解單元桿端的桿端彎矩。 FPRES(i,2)為第i個非結點荷載的類型。該類型是在下面的switch判斷模塊 中定義的。當FPRES(i,2) = 1時，為橫向均布荷載(橫向指垂直于桿件軸線)： 當FPRES(i,2) = 2時，為橫向集中力。即荷載類

28、型分別為1、2,這兩個數值需要 在輸入非結點荷載信息時針對荷載的類型進行輸入。 每種非結點荷載作用下的單元固端力的訃算公式是按照結構力學教材得到 的。其中，FO(1)表示單元左端的固端彎矩；F0(2)表示單元右端的固端彎矩。 return %* * * *璟床4(末末床*床4(床穴*床* *沐沐沐沐* He沐* * *璟* *床末床*床* *末床* * *沐*沐沐* He沐* * *璟* *穴* * *床* *末* * %總剛、總荷載進行邊界條件處理函數ele_BOUNDARY %入口參數：零位移個數，約束信息，總剛矩陣，總荷載矩陣 %出口參數：總剛矩陣，總荷載矩陣 % function HK

29、,FORCE = ele_BOUNDARY(NVFIX,FIXED,HK,FORCE) % forj=l:NVFIX N1=FIXED0); %對約束個數進行循環(huán) %找出約束所在的結點位移編碼 HK(N1,N1)=HK(N1,N1)*1E1O; %乘大數法 end 說明： 采用乘大數法處理邊界條件，僅需對總剛元素進行處理即可。 FIXED(j)為第j個約束所在的結點位移編碼(連續(xù)梁即結點號)賦值給N1。 對約束個數進行循環(huán)，對每j個約束處理時，需對總剛矩陣中的第j行笫j 列的主對角線元素HK(N1,N1)乘以一個大數，如10。 若采用置0置1法處理邊界條件，其相應程序為 forj=l:NVFI

30、X Nl=FIXED(j); %找出約束所在的結點位移編碼 HK(:1)=0; %將總剛矩陣HK(:,N1)的第N1列所有元素賦值為0o HK(N1,:)=O; %將總剛矩陣HK(N1,:)的第N1行所有元素賦值為0。 HK(N1ZN1)=1; %將總剛矩陣HK(N1ZN1莊對角線元素賦值為1。 FORCE(N1)=0; %將總荷載矩陣FORCE(Nl)的笫N1行元素賦值為0 end % return %計算單元桿端力函數ele_MOMENTS %入口參數：輸出文件句柄，結點數，結構位移列陣，單元數，單元信息, 結點坐標，彈性模量，非結點荷載個數，非結點荷載信息 %出口參數：單元桿端力一左右兩

31、端的桿端彎矩 function FE = ele_MOMENTS(FP2,NPOIN,DISP/NELEM, LNODS,. COORD, YOUNG, NFPRES, FPRES) fprintf(FP2/n計算結果nn); fprintf(FP2z*輸出結構結點位移結點號轉角心; for i=l:NPOIN fprintf(FP2；%3d% n；iQISP(i); end 說明： 在輸出文件FP2中顯示提示信息：計算結果。 通過循環(huán)按照結點位移編碼(連續(xù)梁結構即結點編碼)的順序，在文件 FP2中輸出所有結點位移的數值。 EDISP=zeros(2,l); %單元位移列向量清零 fprint

32、f(FP2/輸出單元桿端內力，以順時針為正rf); for i=l:NELEM%對單元個數進行循環(huán)，依次計算每一單元的桿端力 forj=l:2 Nl=LNODS(izj); %單元桿端結點號 EDISP(j)=DISP(Nl); %取N1端的單元位移列向量 end 說明：通過循環(huán)，依次從結構位移列陣中對號，賦值給笫i個單元 的單元位移向量EDISP； j=l時N1為左端，j=2時N1為右端。 % EK=ele_EK(i,LNODS,COORD,YOUNG); % 生成局部單剛(局部坐標) FE=EK*EDISP; %計算局部坐標系下的曲結點位移產生的桿端力 說明：計算第i個單元的單剛EK： 由

33、單剛乘以單元位移得到單元僅曲結點位移產生的桿端力。 forj=l:NFPRES if FPRES(j,l) = i %成立時，當前單元上有非結點荷載 FO=ele_FPRES(j,FPRES,LNODS,COORD); %計算單元固端力 FE=FE+FO; %考慮山非結點荷載引起的桿端力，直接疊加 end end 說明：對非結點個數進行循環(huán)。 FPRESO,1)為第j個非結點荷載作用的單元號，若等于當前循環(huán)的單元號 i,則表明當前單元上有非結點荷載。 調用計算單元固端力函數，得到非結點荷載引起的單元固端力。 將桿端位移引起的桿端力和非結點荷載引起的單元固端力疊加，得到當 前單元的全部桿端力。

34、fprintf(FP2/單元(d)桿端力彎矩W； fprintf(FP2/ 左端3d 結點%7.2fn,/LNODS(i,l),FE(l); fprintf(FP2；右端3d 結點7.2fn,LNODS(i,2),FE(2); 說明：輸出單元的結點號，對應的桿端彎矩。 end %對單元循環(huán)結束。 % return 程序應用舉例 【例1】利用連續(xù)梁程序計算圖29所示的連續(xù)梁結構。EI = 輸入數據文件為: 34112 12 23 34 2 12 3 1 說明： 第一行：讀入程序控制信息 NELEM=fscanf(FPl；%d,l); %單元個數(單元編碼總數) NPOIN=fscanf(FPl；

35、%dzl); %結點個數(結點編碼總數) NVFIX=fscanf(FPl；%d/l); %約束個數(零位移總數) NFPOIN=fscanf(FPl,%d,l); %結點荷載個數(作用在結點上集中力偶總數) NFPRES=fscanf(FPl；%dzl)； %非結點荷載數(作用在單元上分布荷載總數) YOUNG=fscanf(FPl/%f/l); %彈性模量 第二、三四行：讀入單元連接信息： LN0DS=fscanf(FPl；%f,/3/NELEM);%單元號，左結點號,右結點號,慣性矩 第五、六、七、八行：讀入結點坐標 COORD=fscanf(FPl/,%f,/NPOIN)*; %結點坐

36、標值 笫九行：讀入結點荷載信息 FP0IN=fscanf(FPl,%f,2,NFP0IN); %結點號,結點荷載數值(以順時針為正) 第十、十一行：讀入非結點荷載信息 FPRES=fscanf(FPl,%f,/3/NFPRES);%荷載號，單元號，荷載類型,荷載大小 %荷載類型：1為均布荷載，2為集中力 第十二行：讀入零位移信息 FIXED=fscanf(FPl/%f,NVFIX); %零位移約束所在的結點位移編碼 輸出數據文件 結構初始數據 單元總數=3結點總數=4約束總數 結點荷載個數“非結點荷載個數=2彈性模量=1 單元定義數據：單元號左結點號右結點號慣性矩 112 +000 223+0

37、00 +000 結點定義數據：結點號X坐標 4 結點荷載數據：結點號M力偶 1 2 非結點荷載數據：荷載號單元號荷載類型荷載大小 23 零位移約束數據：位移編號 位移編號 4 計算結果 輸岀結構結點位移 結點號 轉角 +001 +001 +001 4 輸出單元桿端彎矩 單元桿端力 左端1結點 右端2結點 單元桿端力 左端2結點 右端3結點 單元桿端力 左端3結點 右端4結點 % 根據單元桿端力計算結果，可得結構的彎矩圖為: （順時針為正） 彎矩 彎矩 彎矩 M圖(kN m) 程序的改進要點 上述連續(xù)梁程序反映了矩陣位移法的基本思路，可以訃算簡單的連續(xù)梁問 題。在熟練掌握了程序的編制與使用后，可

38、在以下兒方面對程序進行改進，以 加深對矩陣位移法及MATLAB語言編程的理解： 1、本程序的單元信息輸入、結點坐標輸入均分別隱含要求從1單元、1結 點按編碼的序號順序輸入，考慮如何改進使相關信息的輸入可以更加靈活。 2、本程序的彈性模量僅能輸入一個數值，意味著程序僅能計算山同種材料 構成的結構。考慮如何改進使程序可計算由不同材料構成的組合結構。 3、本程序僅能計算兩種非結點荷載類型，考慮如何編制跨中作用集中力 偶、單元任意位置作用橫向集中力、單元部分范圍內作用均布荷載等類型的程 序。 4、本程序的單元桿端力僅能計算桿端彎矩，考慮如何編制計算桿端剪力的 程序 5、連續(xù)梁結構山于局部坐標系與整體坐標系相同，因此未涉及到坐標轉換 的問題，考慮如何編制坐標轉換的程序。 6、連續(xù)梁結構由于每個結點僅有一個轉角位移，所以程序中單剛、總剛、 荷載矩陣等多處均采用結點號來替代結點位移號，進行下表索引和循環(huán)控制。 但對于其他結構如平面剛架，每個結點一般有三個結點位移，則程序的編制過 程中需要嚴格分清結點號和結點位移號?？紤]如何編制平面剛架的訃算程序。 綜上所述，本章的連續(xù)梁程序體現了如何將矩陣位移法的計算方法和過程 用MATLAB程序語言表達出來，重點