平面三角形3節點有限元程序_第1頁
平面三角形3節點有限元程序_第2頁
平面三角形3節點有限元程序_第3頁
平面三角形3節點有限元程序_第4頁
平面三角形3節點有限元程序_第5頁
已閱讀5頁,還剩12頁未讀 繼續免費閱讀

下載本文檔

版權說明:本文檔由用戶提供并上傳,收益歸屬內容提供方,若內容存在侵權,請進行舉報或認領

文檔簡介

1、平面三角形結點有限元程序1、程序名:FEMT3.FOR,FEMT3.EXE2、程序功能本程序能計算彈性力學的平面應力問題和平面應變問題;能考慮自重和結點集中力兩種荷載的作用,在計算自重時y軸取垂直向上為正;能處理非零已知位移,如支座沉降的作用。主要輸出的內容包括:結點位移、單元應力、主應力、第一主應力與x軸的夾角以及約束結點的支座反力。程序采用Visual Fortran 5.0編制而成,輸入數據全部采用自由格式。3、程序流程及框圖輸入數據信息形成K形成R計算,輸出應力解KR,輸出啟 動停 機 圖-1 程序流程圖圖-2 程序框圖其中,各子程序的功能如下:INPUT輸入結點坐標、單元信息和材料參

2、數;MR形成結點自由度序號矩陣;FORMMA形成指標矩陣MA(N)并調用其他功能子程序,相當于主控程序;DIV取出單元的3個結點號碼和該單元的材料號并計算單元的bi,ci等;MGK形成整體勁度矩陣并按一維存放在SK(NH)中;LOAD形成整體結點荷載列陣F;OUTPUT輸出結點位移或結點荷載;TREAT由于有非零已知位移,對K和F進行處理;DECOMP整體勁度矩陣的分解運算;FOBA前代、回代求出未知結點位移;ERFAC計算約束結點的支座反力;KRS計算單元勁度矩陣中的子塊Krs。4、輸入數據及變量說明當程序開始運行時,按屏幕提示,鍵入數據文件的名字。在運行程序之前,必須根據程序中輸入要求建立

3、一個存放原始數據的文件,這個文件的名字由少于8個字符或數字組成。數據文件包括如下內容:總控信息,共一條,9個數據NP,NE,NM,NR,NI,NL,NG,ND,NCNP結點總數;NE單元總數;NM材料類型總數;NR約束結點總數;NI問題類型標識,0為平面應力問題,1為平面應變問題;NL受荷載作用的結點的數目;NG考慮自重作用為1,不計自重為0;ND非零已知位移結點的數目;NC要計算支座約束反力的結點數目。材料信息,共NM條,每條依次輸入EO,VO,W,tEO彈性模量(kN/m2);VO泊松比;W材料容重(kN/m3);t單元厚度(m)。這些信息都存放在數組AE(4,NM)中。坐標信息,共NP條

4、,每條依次輸入IP,X,YIP結點號;X,Y分別為結點的x坐標和y坐標。坐標信息存放在數組X(2,NP)中。單元信息,共NE條,每條依次輸入JE,L ,Io,Jo,MoJE單元號;L為該單元的材料類型號。Io,Jo,Mo分別為該單元i,j,m的整體編碼。單元信息存放在數組MEO(4,NE)中。約束信息共NR條,每條依次輸入一個數 ××IP IP結點號;Ix,Iy分別為該結點的約束情況,如果方向受約束時填0,如果自由則填1。荷載信息,共NL條,每條依次輸入yIP,Fx,FyIP結點號;Fx,Fy分別為該結點的x,y方向的荷載分量(kN)。結點號存放在數組NF(NL)中,結點荷

5、載分量存放在數組FV(2,NL)中。若ND0,輸入非零已知位移信息,共ND條,每條依次輸入IP,ux,uyIP結點號;ux,uy分別為該結點x,y方向已知位移分量(m),若其中某方向為自由,則其相應分量為0。結點號存放在數NDI(ND)中,已知位移分量存放在數組DV(2,ND)中。支座反力信息,共NC條,每條依次輸入IP,M1,M2,M3,M4IP 支座結點號;M1,M2,M3,M4 為與該支座結點相關的單元號,若不足4個,則用0補充。支座結點號存放在數組NCI(NC)中,相關單元號存放在數組NCE(4,NC)中。以上數據須按如上順序存放在數據文件中。除此之外,程序中還用到其他一些主要變量和數

6、組,說明如下:N結構自由度總數;NH 按一維存貯的整體勁度矩陣的總容量;NX 最大半帶寬;SK(10000)維存貯的勁度矩陣;R(1000)開始存放等效結點荷載,求解方程以后,用來存放結點位移;B(6)存放單元應力x,y,xy,1,2,;MA(1000)主元素序號指標矩陣;JR(2,500)結點自由度序號矩陣;ME(3)存放單元結點i,j,m的整體編碼;NN(6)單元結點自由度序號;BI(3),CI(3)單元勁度矩陣計算公式中的bi,bj,bm和ci,cj,cm;S三角形單元的面積;H11,H12,H21,H22單元勁度矩陣中子塊Krs的4個元素。5、算例一個正方形彈性體,厚度為1m,四邊受單

7、位均布法向力作用,由于對稱性,取其1/4進行計算,其有限元網格如圖2-3所示,設,不考慮自重。該問題的精確解應力為1,1,0。圖 -3 有限元網格 (1)輸入文件數據6 4 1 5 0 3 0 0 52000.0 0.0 0.0 1.01 0.0 2.02 0.0 1.03 1.0 1.04 0.0 0.05 1.0 0.06 2.0 0.01 1 3 1 2 2 1 2 4 5 3 1 3 2 5 4 1 5 6 3 1 0 1 2 0 1 4 0 0 5 1 0 6 1 01 -0.5 -0.53 -1.0 -1.06 -0.5 -0.51 1 0 0 02 1 2 3 04 2 0 0

8、05 2 3 4 06 4 0 0 0(2)輸出文件結果 NODAL DISPLACEMENTS NODE X-COMP. Y-COMP. 1 0.00000E+00 -0.10000E-02 2 0.00000E+00 -0.50000E-03 3 -0.50000E-03 -0.50000E-03 4 0.00000E+00 0.00000E+00 5 -0.50000E-03 0.00000E+00 6 -0.10000E-02 0.00000E+00 ELEMENT STRESSESELEMENT X-STRESS Y-STRESS XY-STRESS MAX-STRESS MIN-S

9、TRESS ANGLE 1 -1.000 -1.000 0.000 -1.000 -1.000 90.000 2 -1.000 -1.000 0.000 -1.000 -1.000 90.000 3 -1.000 -1.000 0.000 -1.000 -1.000 90.000 4 -1.000 -1.000 0.000 -1.000 -1.000 90.000 NODE STRESSESNODE X-STRESS Y-STRESS XY-STRESS MAX-STRESS MIN-STRESS ANGLE 1 -1.000 -1.000 0.000 -1.000 -1.000 90.000

10、 2 -1.000 -1.000 0.000 -1.000 -1.000 90.000 3 -1.000 -1.000 0.000 -1.000 -1.000 90.000 4 -1.000 -1.000 0.000 -1.000 -1.000 90.000 5 -1.000 -1.000 0.000 -1.000 -1.000 90.000 6 -1.000 -1.000 0.000 -1.000 -1.000 90.000 NODAL REACTIONS NODE X-COMP Y-COMP 1 0.000 0.000 2 1.000 0.000 4 0.500 0.500 5 0.000

11、 1.000 6 0.000 0.0006、源程序C FINITE ELEMENT PROGRAM FOR TWO DIMENSIONALC TRIANGLE ELEMENTC DIMENSION K(800000),COOR(2,3000),AE(4,11),* MEL(5,2000),MA(6000) CHARACTER*32 dat COMMON /CA/ NP,NE,NM,NR,NI,NL,NG,ND,NC WRITE(*,300) 300 FORMAT(/' ' * ':*'/'+PLEASE INPUT FILE NAME OF DATA&#

12、39;) READ(*,*)data OPEN (4,FILE=data,STATUS='OLD') OPEN (7,FILE='OUT',STATUS='UNKNOWN') READ (4,*) NP,NE,NM,NR,NI,NL,NG,ND,NC C WRITE (*,400) NP,NE,NM,NR,NI,NL,NG,ND,NC C WRITE (7,400) NP,NE,NM,NR,NI,NL,NG,ND,NC CALL INPUT (JR,COOR,MEL,AE) CALL CBAND (MA,JR,MEL) CALL SK0(SK,R

13、,COOR,MEL,MA,JR,AE) CALL LOAD (COOR,MEL,R,JR,AE) IF (ND.GT.0) CALL TREAT (SK,MA,JR,R) CALL DECOMP (SK,MA) CALL FOBA (SK,MA,R) WRITE(*,650) WRITE(7,650) CALL OUTPUT(JR,R) WRITE(*,700) WRITE(7,700) CALL CES (COOR,MEL,JR,R,AE) IF(NC.GT.0) call ERFAC (COOR,MEL,JR,R,AE) 400 FORMAT (/2X,'NP=',I3,2

14、X,'NE=',I3,2X,'NM=' *,I3,2X,'NR=',I3,2X,'NI='I3,2X,'NL=',I3,2X, *'NG=',I3,2X,'ND=',I3,2X,'NC=',I3)500 FORMAT(1X,'TOTAL DEGREES OF FREEDOM N=', *I4,1X,'TOTAL-STORAGE ','NH=',I5,1X, *'MAX-SEMI-BANDWIDTH MX='

15、;,I3)550 FORMAT(/20X,'TOTAL STORAGE IS *GREATER THAN 50000')600 FORMAT(30X,'NODAL FORCES'/8X,'NODE', *11X,'X-COMP.',14X,'Y-COMP.')650 FORMAT(/30X,'NODAL DISPLACEMENTS'/8X, *'NODE',13X,'X-COMP.',12X,'Y-COMP.')700 FORMAT(/30X,'

16、;ELEMENT STRESSES'/5X, *'ELEMENT',5X,'X-STRESS',3X,'Y-STRESS', *2X,'XY-STRESS',1X,'MAX-STRESS',1X, *'MIN-STRESS',6X,'ANGLE'/) STOP END C * SUBROUTINE KRS (BR,BS,CR,CS) COMMON /CB/ EO,VO,W,T,A,H11,H12,H21,H22 *,ME(3),BI(3),CI(3) ET=EO*T/(1.0

17、-VO*VO)/A/4.0 V=(1.0-VO)/2.0 H11=ET*(BR*BS+V*CR*CS) H12=ET*(VO*BR*CS+V*BS*CR) H21=ET*(VO*CR*BS+V*BR*CS) H22=ET*(CR*CS+V*BR*BS) RETURN END C * SUBROUTINE INPUT (JR,COOR,MEL,AE) DIMENSION JR(2,*),COOR(2,*),AE(4,*),MEL(3,*) COMMON /CA/ NP,NE,NM,NR COMMON /CC/ N,MX,NH DO 70 I=1,NP READ(4,*) IP,X,Y COOR(

18、1,IP)=X COOR(2,IP)=Y 70 CONTINUE DO 11 J=1,NE READ(4,*)NEE,NME,(MEL(I,NEE),I=1,3) MEL(3,NEE)=NME11 CONTINUE DO 10 I=1,NP DO 10 J=1,2 10 JR(J,I)=1 DO 20 I=1,NR READ(4,*) IP,IX,IY JR(1,IP)=IX JR(2,IP)=IY 20 CONTINUE N=0 DO 30 I=1,NP DO 30 J=1,2 IF (JR(J,I) 30,30,25 25 N=N+1 JR(J,I)=N 30 CONTINUE DO 55

19、 J=1,NM READ (4,*)jj,(AE(I,jj),I=1,4)C WRITE(*,910)jj,(AE(I,jj),I=1,4)IF(NI.eq.1)then AE(1,jj)=AE(1,jj)/(1.0-AE(2,jj)*AE(2,jj) AE(2,jj)=AE(2,jj)/(1.0-AE(2,jj) endif55 CONTINUE910 FORMAT (/20X,'MATERIAL PROPERTIES'/(3X,I5,4(1x,E8.3) RETURN ENDC * SUBROUTINE CBAND (MA,JR,MEL) DIMENSION MA(*),J

20、R(2,*),MEL(3,*),NN(6) COMMON /CA/ NP,NE,NM,NR COMMON /CC/ N,MX,NH DO 65 I=1,N65 MA(I)=0 DO 90 IE=1,NE DO 75 K=1,3 IEK=MEL(K,IE) DO 95 M=1,2 JJ=2*(K-1)+M NN(JJ)=JR(M,IEK)95 CONTINUE75 CONTINUE L=N DO 80 I=1,6 NNI=NN(I) IF(NNI.EQ.0) GO TO 80 IF(NNI.LT.L) L=NNI 80 CONTINUE DO 85 M=1,6 JP=NN(M) IF(JP.EQ

21、.0) GO TO 85 JPL=JP-L+1 IF(JPL.GT.MA(JP) MA(JP)=JPL 85 CONTINUE 90 CONTINUE MX=0 MA(1)=1 DO 10 I=2,N IF(MA(I).GT.MX) MX=MA(I) MA(I)=MA(I)+MA(I-1) 10 CONTINUE NH=MA(N) WRITE (*,500) N,MX,NH WRITE (7,500) N,MX,NH500 FORMAT (/5X,'FREEDOM N=' *,I5,3X,'SEMI-BANDWI. MX=',I5,3X, * 'STOR

22、AGE NH=',I7) RETURN END C* SUBROUTINE SK0(SK,R,COOR,MEL,MA,JR,AE) DIMENSION AE(4,*),COOR(2,*),MEL(3,*),JR(2,*),R(*), *MA(*),SK(*),SKE(6,6),NN(6) COMMON /CA/ NP,NE,NM,NR,NI,NL,NG,ND,NC COMMON /CB/ EO,VO,W,T,A,H11,H12,H21,H22, *ME(3),BI(3),CI(3) COMMON /CC/ N,NH DO 10 I=1,NH 10 SK(I)=0.0 DO 70 IE=

23、1,NE CALL DIV (IE,COOR,MEL,AE) DO 30 I=1,3 DO 30 J=1,3 CALL KRS (BI(I),BI(J),CI(I),CI(J) SKE(2*I-1,2*J-1)=H11 SKE(2*I-1,2*J)=H12 SKE(2*I,2*J-1)=H21 SKE(2*I,2*J)=H22 30 CONTINUE DO 40 I=1,3 J2=ME(I) DO 40 J=1,2 J3=2*(I-1)+J NN(J3)=JR(J,J2) 40 CONTINUE DO 60 I=1,6 DO 60 J=1,6 IF(NN(J).EQ.0.OR.NN(I).LT

24、.NN(J) GO TO 60 JJ=NN(I) JK=NN(J) JL=MA(JJ) JM=JJ-JK JN=JL-JM SK(JN)=SK(JN)+SKE(I,J) 60 CONTINUE 70 CONTINUE c WRITE(0,500) (SK(I),I=1,20) 500 FORMAT(/10X,'SK='/(6F12.5) RETURN END C * SUBROUTINE LOAD (COOR,MEL,R,JR,AE) DIMENSION AE(4,*),COOR(2,*),MEL(3,*),R(*),JR(2,*), * NF(50),FV(2,50) COM

25、MON /CA/ NP,NE,NM,NR,NI,NL,NG,ND,NC COMMON /CB/ EO,VO,W,T,A,A1(4),ME(3),BB(6) COMMON /CC/ N,NH DO 10 I=1,N 10 R(I)=0.0 IF(NG) 70,70,30 30 DO 60 IE=1,NE CALL DIV (IE,COOR,MEL,AE) DO 50 I=1,3 J2=ME(I) J3=JR(2,J2) IF(J3) 50,50,40 40 R(J3)=R(J3)-T*W*A/3.0 50 CONTINUE 60 CONTINUE 70 IF(NL) 110,110,80 80

26、READ(4,*) (NF(I),I=1,NL) READ(4,*) (FV(I,J),I=1,2),J=1,NL) C WRITE(*,500) (NF(I),I=1,NL) C WRITE(7,500) (NF(I),I=1,NL)C WRITE(*,600) (FV(I,J),I=1,2),J=1,NL) C WRITE(7,600) (FV(I,J),I=1,2),J=1,NL) DO 100 I=1,NL JJ=NF(I) J=JR(1,JJ) M=JR(2,JJ) IF (J.GT.0) R(J)=R(J)+FV(1,I) IF (M.GT.0) R(M)=R(M)+FV(2,I)

27、 100 CONTINUE 110 RETURN 500 FORMAT(/20X,'NODES OF APPLIED LOAD*NF=' */(1X,10I8)600 FORMAT(/30X,'LUMPED-LOADS*FV=' */(5X,5F15.3) END C * SUBROUTINE TREAT (SK,MA,JR,R) DIMENSION SK(*),MA(*),NDI(75),DV(2,75),JR(2,*),R(*) COMMON /CA/ NP,NE,NM,NR,NI,NL,NG,ND,NC COMMON /CC/ N,NH READ(4,*) (NDI(J),J=1,ND) READ(4,*) (DV(I,J),I=1,2),J=1,ND) C WRITE(*,500) (NDI(J),J=1,ND) C WRITE(7,500) (NDI(J),J=1,ND)C WRITE(*,550) (DV(I,J),I=1,2),J=1,ND) C WRITE(7,550) (DV(I,J),I=1,2),J=1,ND) DO 20 I=1,ND DO 20 J=1,2 IF(DV(J,I) 10,20,10 10 JJ=NDI(I) L=JR(J,JJ) J

溫馨提示

  • 1. 本站所有資源如無特殊說明,都需要本地電腦安裝OFFICE2007和PDF閱讀器。圖紙軟件為CAD,CAXA,PROE,UG,SolidWorks等.壓縮文件請下載最新的WinRAR軟件解壓。
  • 2. 本站的文檔不包含任何第三方提供的附件圖紙等,如果需要附件,請聯系上傳者。文件的所有權益歸上傳用戶所有。
  • 3. 本站RAR壓縮包中若帶圖紙,網頁內容里面會有圖紙預覽,若沒有圖紙預覽就沒有圖紙。
  • 4. 未經權益所有人同意不得將文件中的內容挪作商業或盈利用途。
  • 5. 人人文庫網僅提供信息存儲空間,僅對用戶上傳內容的表現方式做保護處理,對用戶上傳分享的文檔內容本身不做任何修改或編輯,并不能對任何下載內容負責。
  • 6. 下載文件中如有侵權或不適當內容,請與我們聯系,我們立即糾正。
  • 7. 本站不保證下載資源的準確性、安全性和完整性, 同時也不承擔用戶因使用這些下載資源對自己和他人造成任何形式的傷害或損失。

評論

0/150

提交評論