作業(yè)：一1用偽代碼方式寫出有限元求解步驟解：有限元分析流程圖：形成K消除K奇異性求解Ka=PI形成K消除K奇異性求解Ka=PNGNG:結(jié)構(gòu)的結(jié)點總數(shù)NE,MC,NX,NB.ND,EO,VO,TNE:結(jié)構(gòu)單元總數(shù)MC:計算控制類型參數(shù)=0平面應(yīng)力=1平面應(yīng)變NX:作用荷載組數(shù)NB：給定位移的個數(shù)NB:給定位移的個數(shù)ND：結(jié)構(gòu)剛度陣的半帶寬EO:彈性模量VO：泊松比T：單元的厚度3.NWA.NWE,NWKNVP,NWDNWA：單元參數(shù)的輸出控制參數(shù)NWE:單元剛度矩陣的輸出控制參數(shù)NWK：結(jié)構(gòu)剛度矩陣的輸出控制參數(shù)NWP：荷載向量的輸出控制參數(shù)NWD：結(jié)點位移的輸出控制參數(shù)輸出控制參數(shù)=1輸出=0不輸出4.IJM(3,NE):單元結(jié)點編碼數(shù)組UM(1,I),UM(2,I)JJM(3,I),第I個三角形單元的節(jié)點號，按單元編號順序填寫。XY(2,NG)：結(jié)點坐標(biāo)數(shù)組XY(1,I):第I個結(jié)點的坐標(biāo)，XY(2,I):第I個結(jié)點的Y坐標(biāo)按結(jié)點編號順序填寫。MB(2,NB).ZB(NB):給定位移約束的信息數(shù)組與值數(shù)組MB(1,I):第I個給定位移所在的結(jié)點號NB(2,1)=1：給定X方向位移=2：給定Y方向位移ZB(NB):給定位移值(以坐標(biāo)正向為正)NF,NPNF：作用于結(jié)點上的集中荷載的個數(shù)NP：作用于均布側(cè)壓的單元邊數(shù)若NF>0,填8MF(2,NF),ZF(NF):作用于結(jié)點上集中荷載的信息組與值數(shù)組MF(1,I):第I個集中荷載作用的結(jié)點號MF(2,1)=1:作用于x方向的集中力=2:作用y方向的集中力ZF(NF):作用的集中力值若NP>0,填9MP(2,NP),ZP(NP):作用于單元邊上均布荷載的信息數(shù)組與值數(shù)組MF(LI):第I個均布荷載作用邊的起始結(jié)點號MF(2,I):第I個均布荷載作用邊的終止結(jié)點號，逆時針排列ZP(NP)：第I個均布荷載值輸入數(shù)據(jù)格式，建立數(shù)據(jù)文件，文件名小于12個字符：NGNE,MC,NX,NB.ND,E,P,T3.NWA,NWE,NWKNWP,NWD4.IJM(3,NE)XY(2,NG)MB(2.NB),ZB(NB)NF,NP當(dāng)NF>0測填8MF(2,NF),ZF(NF)當(dāng)NP>0,則填9MP(2?NP),ZP(NP)若NX>1,即多組荷載情況，重復(fù)7?9，若計算多個結(jié)構(gòu)則重復(fù)1?9。結(jié)束(2)平面問題主程序框圖~——HNGNE.MC.NX.NB.NDEP.NVANWE.NWKNWRNWDL輸入子程序INPUT偽代碼DIMENSIONGUM(3,NE),XY(2,NG)?MB(2,NB),ZB(NB)READ(5,*)((IJM(I,L),I=1,3),L=I,NE)READ(5,*)((XY(I,J),I=1,2),J=LNG)READ(5,*)(MB(I,L),I=1,2),L=1,NB),(ZB(L),L=1,NB)WRITE(6,10)10FORMAT(/10X,4ELEMENTCODEBLOCKSUM(2,NE)'/)WRITE(6,20)((UM(MJ)M=l,3),1=1,NE)20FORMAT(IX,4(314,3X),3I4)WRITE(6,30)30FORMATE。」OX,4COORDINATESOFNODESXY(2,NG)/)WRITE(6,40)((XY(M,I),M=l,2),1=1,NG)40FORMATE(1X,6E12.5)WRITE(6,50)50FORMATE(/1OX.INRMATIONANDVALUESOFGIVENDISPLACEMENTS/20X,'MBL15X,^ZB\/)DO701=1.NBWRITE(6,60)(MB(K,I),K=l,2),ZB(I)60FORMATE(1OX,2I5,5X.E14.6)70CONGTINUERETURNEND計算單元參數(shù)子程序偽代碼1、X(2,5)當(dāng)前計算單元的結(jié)點坐標(biāo)2、B(7)當(dāng)前計算單元的單元參數(shù)3、BCA(7,NE)將B(7)按單元存儲于該數(shù)組中DIMENSIONIJM(3,NE),XY(2,NG),BCA(7,NE)X(2,5),B(7)IF(NEWA?EQ』)WRITE(6,5)5FORMAT(/1OX,TARAMETERSOFELEMENTSBCA(7,NE)7)DO801=1.NEDO10K=l,3K1=UM(KJ)DO10J=l,2X(J,K)=XY(J?K1)10CONTINUEDO20J=l,2X(J,4)=X(J,1)X(J,5)=X(J,2)20CONTINUEDO30K=l,3B(K)=X(2,K+l).X(2,K+2)B(K+3)=X(LK+2).X(LK+1)30CONTINUEB(7)=(B(1)*B(5).B(4)*B(2))*0.5IF(NWA.GT.O)WRITE(6,40)I,B40FORMAT。X；NE=',13,/3X,7E10.4)IF(B(7).LE.0.0)GOTO60DO50J=l,7BCA(J,I)=B(J)50CONTINUEDOTO8060WIRTE(6,70)I,(UM(J,I),J=l,3)70FORMAT(/5X,ELEMENT\I5,5X,5AREAISNONPOSITIVE\5X.1JM\315)STOP80CONTINUERETURNEND計算單元剛度矩陣程序KE偽代碼1、BCA(7,NE)單元參數(shù)組2、EK(6,6)單元剛度矩陣DIMENSIONB(7).BCA(7,NE),EK(6.6)DO101=1,7B(I)=BCA(I,IO)10CONTINUEA=A1/B(7)*TDO201=1,3DO20J=L311=2*1J1=2*JEK(Il.l,Jl-l)=A*(B(I)*B(J)+A2*B(I+3)*B(J+3))EK(Il.l,Jl)=A*(V*B(I)*B(J+3)+A2*B(I+3*B(J))EK(Il,Jl-l)=A*(V*B(I+3)*B(J)+A2*B(I)*B(J+3)EK(ILJl)=A*(B(I+3)*B(J+3)+A2*B(I)*B(J))20CONTINEDO301=3,6DO30J=1,IEK(IJ)=EK(J,I)30CONTINUEIF(NEW.EQ.O)GOTO60WRITE(6,40)IO40FORMATflX,4EKNE=\I5)WRITE(6,50)EK50FORMAT(1X,6E1L4)60RETURNEND4、總剛集成子程序偽代碼IO按單元循環(huán)的單元號EK(6,6)單元剛度矩陣IJM(3.NE)單元結(jié)點編碼數(shù)組SK(NT,ND)按二維等寬存儲的總剛度矩陣DOMENSIONIJ(3),SK(NT,ND),IJM(3.NE),EK(6,6)DO101=1,3IJ(I)=UM(I,IO)10CONTINUEDO201=1,3DO20J=l,3IF(IJ(I).IT(J))GOTO20M=2*IJ(I)-1N=2*IJ(J).U(I)+1MO=2*I-1NO=2*J-1SK(M,N)=SK(M,N)+EK(MO.NO)SK(M,N+1)=SK(M,N+1)+EK(MO,NO+1)SK(M+1,N+1)=SK(M+1,N+1TEK(MO+1,NO+1)IF(IJ(I).EQ?IJ(J))GOTO20SK(M+1,N.1)=SK(M+1,N.1)+EK(MO+1,NO)20CONTINUERETURNEND5、集中荷載子程序PF偽代碼參數(shù)及數(shù)組說明：NF集中荷載的個數(shù)MF(2,NF)荷載信息組數(shù)ZF(NF)集中荷載數(shù)組F(NT)等效結(jié)點荷載列陣DEMENSIONMF(2,NF),ZF(NF),F(NT)READ(5,*)(MF(I,L),I=1,2),L=1,NF),(ZF(L),L.1,NF)WRITE(6,10)10FORMAT^1OX/CONCEENTRATEDLOADSMF(2,NF)ZF(NF)7)DO30L=LNFWRITE(6,20)(MF(I,L),I=1,2),ZF(L))20FORMAT(IX,2110,10X,E12.6)30CONTNUEDO401=1,NFN=2*MF(1,I).MF(2,I)F(N)=F(N)+ZF(I)40CONTINUEIF(NP.GT.O)GOTO70IF(NWP.EQ?0)GOTO70WRITE(6,50)50FORMAT^1OX.TOADVECTOR?)WRITE(6,60)F60FORMAT(1X,5E14.6)70RETURNEND6、均布側(cè)壓等效結(jié)點荷載子程序PP偽代碼參數(shù)及數(shù)組說明NP作用側(cè)壓的單元邊數(shù)MP(2,NP)作用側(cè)壓信息數(shù)組ZP(NP)均布側(cè)壓值數(shù)組XY(2,NG)結(jié)點坐標(biāo)數(shù)組F(NT)等效結(jié)點荷載陣DIMENSIONMP(2,NP),ZP(NP),XY(2,NG),F(NT)READ(5,*)((MF(I,L),I=I,2),L=1,NF),(ZF(L),L=1,NF)WRITE(6,10)10FORMAT^1OX.^UNIFORMPRESSUREMP(2,NP)ZP(NP)7)DO30L=l,NPWRITE(6,20)(MP(I,L),I=1,2),ZP(L)20FORMAT(IX,2110,10X.E12.6)30CONTINUEDO401=1,NPNl=MP(bI)N2=MP(2,I)PX=XY(2,N1).XY(2,N2)PY=XY(1,N2).XY(LNl)PX=5*ZP(I)*PXPY=5*ZP(I)*PYF(2*N1.1)=F(2*N1.1)+PXF(2*N1)=F(2*N1)+PYF(2*N2-1)=F(2*N2.1)+PXF(2*N2)=F(2*N2)+PY40CONTINUEIF(NWP.EQ.O)GOTO70WRITE(6,50)50FORMAT(/,10X,IOADVECTOR?)WRITE(6,60)F60FORMAT(1X,5E14.6)70RETURNEND7、定位移子程序DBC偽代碼參數(shù)說明：A(NT,ND)系數(shù)矩陣B(NT)等效結(jié)點荷載矩陣NB給定位移個數(shù)MB(2,NB)給定位移信息數(shù)組ZB(NB)給定位移數(shù)組DIMENSIONMB(2,NB),ZB(NB).A(NT,ND),B(NT)DO601=1,NBN=2*MB(1,I).MB(2,I)Z=ZB(I)IF(ABS(Z).LT.1E-1O)GOTO20IF(NX.NE.NXl)GOTO10A(N,1)=A(N,1)*1E+1510B(N)=A((N,1)*ZGOTO6020IF(NX.NE.NX1)GOTO50A(N,l)=1.0DO30J=2,NDA(NJ)=0.030CONTINUEDO40K=2,NDIF(N.LT.K)GOTO50M=N-K+]A(M,K)=0.040CONTINUE50B(N)=0.060CONTINUERETURNEND單元應(yīng)力子程序stress偽代碼參數(shù)說明：1.A1,A2,V2.IJNI(3,NE)單元結(jié)點編碼數(shù)組BCA(7,NE)單元參數(shù)數(shù)組F(NT)結(jié)點位移DIMENSIONIJM(3,NE),BC(7,NE).F(NT),B(7),R(6)WRITE(6,5)5FORMAT^1OX,ELEMENTSTRESSES?)WRITE(6,10)10FORMAT(2X；NE\5X/S-X\8X；S-XY,,10X；S1\10X；S2\8X/AG7)DO601=1.NESl=0.S2=0.S3=0.DO20J=l,7B(J戶BCA(J,I)20CONTINUEA=2*AUB(7)DO30J=l,3N=UM(J,I)*2R(2*J-1)=F(N.1)R(2*J)=F(N)30CONTINUEDO40J=l,3K=2*JS1=S1+A*(B(J)*R(K.1)+V*B(J+3)*R(K))S2=S2+A*(V*B(J)*R(K.1)+B(J+3)*R(K))S3=S3+A*A2*(B(J+3)*R(K-1)+B(J)*R(K))40CONTINUEP=.5*(S1+S2)Q=.5*(S1?S2)X1=P+SQRT(Q*Q+S3*S3)X2=2*P-X1CTA=0.IF(S3.GT.O)CTA=ATAN((X1-S1)/S3)WRITE(6,50)I,SI,S2,S3,XLX2,CTA50FORMAT(1X,I3,1X,E10.4,2X,E10.4,2X,E10.4,2X,E10.4,2X,E10.4,,2X,F8.4)60CONTNUERETURNEND2求集中力P=20作用點的位移，己知各豎直邊邊長為1,水平邊長為2,E=2,u=0oft?：(i)先求單元剛度矩陣［KY<1000-10)(200、[B]=l/200020-2[D]=020k0120-2-1；ob<1000-10、[S]=[D][B]=00020-200.510-1-0.5r1000-10、00.510-1-0.5[K]e=[BY[D][B]tA=[BY[S]tA=^20120-2-100040-4-1-1-2031(