你好，請問一下桁架結構帶未知參數(shù)的剛度矩陣如何在matlab中實現(xiàn)？

你提出的問題我之前剛好做過,使用有限元方法來進行桁架結構分析.

目前創(chuàng)新互聯(lián)建站已為上1000家的企業(yè)提供了網(wǎng)站建設、域名、網(wǎng)站空間、綿陽服務器托管、企業(yè)網(wǎng)站設計、豐臺網(wǎng)站維護等服務，公司將堅持客戶導向、應用為本的策略，正道將秉承"和諧、參與、激情"的文化，與客戶和合作伙伴齊心協(xié)力一起成長，共同發(fā)展。

Matlab編程實現(xiàn)平面桿單元分析

首先,明確Matlab程序要實現(xiàn)的5個重要模塊分別為：單元剛度矩陣的求解、單元組裝、節(jié)點位移的求解、單元應力的求解、節(jié)點力的求解.下面給出這5個模塊的實現(xiàn).

1．\x09單元剛度矩陣求解

定義函數(shù)Bar2D2Node_Stiffness,該函數(shù)計算單元的剛度矩陣,輸入彈性模量E,橫截面積A,兩個節(jié)點坐標輸出單元剛度矩陣k(4X4).具體代碼如下：

function k=Bar2D2Node_Stiffness(E,A,x1,y1,x2,y2)

L=sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));

x=acos((x2-x1)/L);

C=cos(x);

S=sin(x);

k=E*A/L*[C*C C*S -C*C -C*S; C*S S*S -C*S -S*S;

-C*C -C*S C*C C*S; -C*S -S*S C*S S*S];

2．\x09單元組裝

定義函數(shù)Bar2D2Node_Assembly,該函數(shù)進行單元剛度矩陣的組裝,輸入單元剛度矩陣k,單元的節(jié)點編號i、j.輸出整體剛度矩陣KK,具體代碼如下：

function z = Bar2D2Node_Assembly(KK,k,i,j)

DOF(1)=2*i-1;

DOF(2)=2*i;

DOF(3)=2*j-1;

DOF(4)=2*j;

for n1=1:4

for n2=1:4

KK(DOF(n1),DOF(n2))= KK(DOF(n1),DOF(n2))+k(n1,n2);

end

end

z=KK;

3．\x09節(jié)點位移的求解

定義函數(shù)Bar2D2Node_Disp(KK,num,p),該函數(shù)輸入KK為總體剛度矩陣；num為活動自由度編號數(shù)組；p為活動自由度方向上的節(jié)點力；輸出節(jié)點位移列陣.具體代碼如下：

function u = Bar2D2Node_Disp(KK,num,p)

k=KK(num,num)

u=k\p

4．\x09單元應力的求解

定義函數(shù)函數(shù)Bar2D2Node_Stress(E,x1,y1,x2,y2,u),該函數(shù)計算單元的應力輸入彈性模量E,第一個節(jié)點坐標（x1,y1）,第二個節(jié)點坐標（x2,y2）單元節(jié)點位移矢量u,返回單元應力標量stress .具體代碼如下：

L=sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));

x=acos((x2-x1)/L);

C=cos(x);

S=sin(x);

stress=E/L*[-C -S C S]*u;

5．\x09計算節(jié)點力

定義函數(shù)Bar2D2Node_Forces(KK,q),該函數(shù)用于計算節(jié)點力,KK為剛度矩陣,q為節(jié)點位移陣列

function P= Bar2D2Node_Forces(KK,q)

q=zeros(8,1);

q(num)=u;

P=KK*q;

至此,基于Matlab的桿單元有限元分析的程序設計已經完成,遇到實際問題時可以直接調用這些函數(shù)就可以解決問.

經典算例

如圖所示的結構,各個桿的彈性模量和橫截面積都為 , .試基于MATLAB平臺求解該結構的節(jié)點位移、單元應力以及支反力.

四桿桁架結構

對該問題進行有限元分析的過程如下

（1） 結構的離散化與編號

\x09對該結構進行自然離散,節(jié)點編號和單元編號如上圖所示

（2）計算各單元的剛度矩陣(基于國際標準單位)

\x09\x09輸入彈性模量E、橫截面積A,各點坐標.然后分別針對單元1,2,3和4,調用4次Bar2D2Node_Stiffness,就可以得到單元的剛度矩陣.

對應的主程序中代碼：

E=2.95e11;A=0.0001;x1=0;y1=0;x2=0.4;y2=0;x3=0.4;y3=0.3;x4=0;y4=0.3;

k1=Bar2D2Node_Stiffness (E,A,x1,y1,x2,y2)

k2=Bar2D2Node_Stiffness (E,A,x2,y2,x3,y3)

k3=Bar2D2Node_Stiffness (E,A,x1,y1,x3,y3)

k4=Bar2D2Node_Stiffness (E,A,x4,y4,x3,y3)

（3） 建立整體剛度方程

由于該結構共有4個節(jié)點,因此,設置結構總的剛度矩陣為KK(8×8),先對KK清零,然后四次調用函數(shù)Bar2D2Node _Assembly進行剛度矩陣的組裝.相關主程序代碼為：

KK=zeros(8,8);

KK=Bar2D2Node_Assembly (KK,k1,1,2);

KK=Bar2D2Node_Assembly (KK,k2,2,3);

KK=Bar2D2Node_Assembly (KK,k3,1,3);

KK=Bar2D2Node_Assembly (KK,k4,4,3)

（4）邊界條件的處理及剛度方程的求解

由圖可以看出,被約束的自由度有：節(jié)點1的x,y方向自由度,節(jié)點2的y方向自由度,4節(jié)點的x、y方向兩個自由度.則活動自由度編號為3,5,6.活動自由度對應的節(jié)點載荷F3=20000N,F5=0N,F6=25000N,采用高斯消去法進行求解,對應的代碼為：

num=[3,5,6];%可活動的自由度編號

p=[20000;0;-25000];

u=Bar2D2Node_Disp(KK,num,p)

（5）支反力的計算

在得到整個結構的節(jié)點位移后,由原整體剛度方程就可以計算出對應的支反力.這部分對應的主程序的代碼如下：

q=zeros(8,1);

q(num)=u;%節(jié)點位移陣列

P=Bar2D2Node_Forces(KK,q)

（6）單元應力的計算

先從整體位移列陣q中提取出單元的位移列陣,然后,調用計算單元應力的函數(shù)Bar2D2Node_ElementStress,就可以得到各個單元的應力分量.

u1=[q(1);q(2);q(3);q(4)]

stress1=Bar2D2Node_Stress(E,x1,y1,x2,y2,u1)

u2=[q(3);q(4);q(5);q(6)]

stress2=Bar2D2Node_Stress(E,x2,y2,x3,y3,u2)

u3=[q(1);q(2);q(5);q(6)]

stress3=Bar2D2Node_Stress(E,x1,y1,x3,y3,u3)

u4=[q(7);q(8);q(5);q(6)]

stress4=Bar2D2Node_Stress(E,x4,y4,x3,y3,u4)

（7）計算結果的整理

通過主程序的運行得計算結果.

主程序

%計算各單元的剛度矩陣(以國際標準單位)

E=2.95e11;

A=0.0001;

x1=0;

y1=0;

x2=0.4;

y2=0;

x3=0.4;

y3=0.3;

x4=0;

y4=0.3;

k1=Bar2D2Node_Stiffness (E,A,x1,y1,x2,y2)

k2=Bar2D2Node_Stiffness (E,A,x2,y2,x3,y3)

k3=Bar2D2Node_Stiffness (E,A,x1,y1,x3,y3)

k4=Bar2D2Node_Stiffness (E,A,x4,y4,x3,y3)

%建立整體剛度方程

%由于該結構共有4個節(jié)點,因此,結構總的剛度矩陣為KK(8×8),先對K清零,然后四次調用函數(shù)Bar2D2Node _Assembly進行剛度矩陣的組裝.

KK=zeros(8,8);

KK=Bar2D2Node_Assembly (KK,k1,1,2);

KK=Bar2D2Node_Assembly (KK,k2,2,3);

KK=Bar2D2Node_Assembly (KK,k3,1,3);

KK=Bar2D2Node_Assembly (KK,k4,4,3)

%邊界條件的處理及剛度方程求解

num=[3,5,6];%可活動的自由度編號

p=[20000;0;-25000];

u=Bar2D2Node_Disp(KK,num,p)

%支反力的計算

q=zeros(8,1);

q(num)=u;%節(jié)點位移陣列

P=Bar2D2Node_Forces(KK,q)

%各單元的應力計算

u1=[q(1);q(2);q(3);q(4)];

stress1=Bar2D2Node_Stress(E,x1,y1,x2,y2,u1)

u2=[q(3);q(4);q(5);q(6)];

stress2=Bar2D2Node_Stress(E,x2,y2,x3,y3,u2)

u3=[q(1);q(2);q(5);q(6)];

stress3=Bar2D2Node_Stress(E,x1,y1,x3,y3,u3)

u4=[q(7);q(8);q(5);q(6)];

stress4=Bar2D2Node_Stress(E,x4,y4,x3,y3,u4)

說的可能有些羅嗦,注意其中有5個function,和最后一個主程序,計算的時候直接運行主程序就可以了.希望能幫助到你.

對于一般層合板，耦合剛度矩陣的存在，表示什么 要素之間的耦合關系

給你說了也白搭，你大概知道就行了，這個需要代碼量。。。簡單地說就是，沒有多態(tài)，那么等號左邊是啥右邊就得是啥，這就叫耦合，有了多態(tài)，左邊是父類（或者接口），右邊是子類（或實現(xiàn)類），我只管調用接口里面的方法就是了，至于你實現(xiàn)類怎么去實現(xiàn)，那是你的事，你要修改一下實現(xiàn)，你只管去把實現(xiàn)類換掉，我這邊一行代碼都不用變，這就解耦了

知道剛度矩陣和質量矩陣后用matlab怎么求固有頻率和固有振型

特征向量都求出來了，用哪一種歸一化還不跟玩兒一樣嗎？okok已經回答你了，這里再貼一下，主要因為已經消失幾天了，出來透口氣：div class="blockcode"復制內容到剪貼板h5代碼:/h5code id="code0"function ziyouzhendong1(k,m)br /% k=600*[1 -1 0;-1 3 -2;0 -2 5];m=diag([1 1.5 2]);br /clcbr /[C,B]=eig(inv(m)*k);br /n=size(m,2);br /for i=1:nbr / Ct=C(:,i);br / zhenxing(:,i)=Ct/Ct(1);br /endbr /w=sqrt(diag(B));br /B1=diag(B)';br /Eigen=inv(m)*k;br /for i=1:nbr / t1(:,i)=B1(i)*zhenxing(:,i);br / t2(:,i)=Eigen*zhenxing(:,i);br /endbr /B,zhenxing,t1,t2,w/code/div其實完全可以不用循環(huán)，這是上學初學MATLAB一個月左右時寫的程序，目的是幫助同班美女從繁重的結構動力學作業(yè)中解脫出來并有充分的時間買化妝品and keep us amused...，當時還沒人學MATLAB這玩意兒，全拿計算器算動力學矩陣的問題，顯然即使最多只有4×4的K矩陣也足夠讓人受折磨了...br /

br /

轉載

patran/nastran如何輸出質量和剛度矩陣

在BDF文件中SOL之下加入以下語句，可以輸出剛度矩陣

compile semg

alter ’kjjz.*stiffness’ $

matprn kjjz// $

本文標題：剛度矩陣java代碼 單元剛度矩陣編程
路徑分享：http://jinyejixie.com/article24/hejgje.html

成都網(wǎng)站建設公司_創(chuàng)新互聯(lián)，為您提供響應式網(wǎng)站、動態(tài)網(wǎng)站、商城網(wǎng)站、外貿建站、移動網(wǎng)站建設、建站公司

廣告

聲明：本網(wǎng)站發(fā)布的內容（圖片、視頻和文字）以用戶投稿、用戶轉載內容為主，如果涉及侵權請盡快告知，我們將會在第一時間刪除。文章觀點不代表本網(wǎng)站立場，如需處理請聯(lián)系客服。電話：028-86922220；郵箱：631063699@qq.com。內容未經允許不得轉載，或轉載時需注明來源： 創(chuàng)新互聯(lián)