空間桁架結(jié)構(gòu)程序設(shè)計(Fortran)

空間桁架結(jié)構(gòu)程序設(shè)計(Fortran)空間桁架結(jié)構(gòu)程序設(shè)計(Fortran)空間桁架結(jié)構(gòu)程序設(shè)計(Fortran)xxx公司空間桁架結(jié)構(gòu)程序設(shè)計(Fortran)文件編號：文件日期：修訂次數(shù)：第1.0次更改批準審核制定方案設(shè)計，管理制度空間桁架靜力分析程序及算例1、變量及數(shù)組說明輸入數(shù)據(jù)控制數(shù)據(jù)NF單個節(jié)點的自由度數(shù)NP結(jié)構(gòu)離散節(jié)點的總數(shù)NE結(jié)構(gòu)離散單元的總數(shù)NM結(jié)構(gòu)中單元不同的特征數(shù)類的總數(shù)NR結(jié)構(gòu)受約束節(jié)點的總數(shù)NCF結(jié)構(gòu)受外荷載作用的節(jié)點總數(shù)ND一個單元的節(jié)點總數(shù)幾何數(shù)據(jù)X(NP)節(jié)點X坐標數(shù)組Y(NP)節(jié)點Y坐標數(shù)組Z(NP)節(jié)點Z坐標數(shù)組ME(ND,NE)單元節(jié)點信息存儲矩陣ME(1,NE)存儲桿件始端節(jié)點號ME(2,NE)儲存桿件末端節(jié)點號RR(2,NR)結(jié)構(gòu)約束信息矩陣RR(1,NR)存放受有約束的節(jié)點號RR(2,NR)存放節(jié)點位移約束情況單元特征數(shù)據(jù)AE(2,IN)單元特征數(shù)類數(shù)組AE(1,IN)單元的彈性模量AE(2,IN)單元的橫截面面積NAE(NE)單元特征類信息存儲數(shù)組荷載數(shù)據(jù)PF(4,NCF)外荷載信息數(shù)組PF(1,NCF)存放外荷載作用的節(jié)點號PF(2,NCF)存放X方向的外荷載PF(3,NCF)存放Y方向的外荷載PF(4,NCF)存放Z方向的外荷載輸出數(shù)據(jù)位移DIST(NPF)節(jié)點位移數(shù)組DIST(NF*I-2)存放I節(jié)點X方向的位移DIST(NF*I-1)存放I節(jié)點Y方向的位移DIST(NF*I)存放I節(jié)點Z方向的位移力SG(NE)單元內(nèi)力數(shù)組SM(NE)單元截面應(yīng)力數(shù)組FL(NF*NR)支座反力數(shù)組FL(NF*I-2)存放受約束的I節(jié)點X方向的反力FL(NF*I-1)存放受約束的I節(jié)點Y方向的反力FL(NF*I)存放受約束的I節(jié)點Z方向的反力中間變量NPF=NF*NP二維總剛度矩陣的最大行數(shù)NDF=ND*NF一個單元的自由度總數(shù)（2*3=6）IN單元特征類總數(shù)AKE(2,2)單元在局部坐標系中的剛度局矩陣BL桿件單元長度T(2,6)坐標轉(zhuǎn)換矩陣TAK(6,6)單元在總體坐標系中的剛度矩陣IT(NF,NP)節(jié)點聯(lián)系數(shù)組LMT(NDF,NE)單元聯(lián)系數(shù)組MAXA(NPF)結(jié)構(gòu)二維總剛度矩陣主對角元地址數(shù)組NWK結(jié)構(gòu)一維總剛度矩陣的總?cè)萘緾KK(NWK)結(jié)構(gòu)一維總剛度矩陣NN結(jié)構(gòu)矩陣方程的方程總數(shù)（去掉約束）NNMNNM=NN+1V(NN)已知節(jié)點荷載列陣數(shù)組，回代完成后為存放結(jié)構(gòu)位移PP(NPF)所有節(jié)點荷載列陣數(shù)組2、空間桁架結(jié)構(gòu)有限元分析程序源代碼!主程序（讀入文件，調(diào)用總計算程序，輸出結(jié)果） CHARACTERIDFUT*20,OUTFUT*20 WRITE(*,*)'InputDataFilename:' READ(*,*)IDFUT OPEN(11,FILE=IDFUT,STATUS='OLD') WRITE(*,*)'OutputFilename:' READ(*,*)OUTFUT OPEN(12,FILE=OUTFUT,STATUS='UNKNOWN') WRITE(12,*)'*****************************************' WRITE(12,*)'*ProgramforAnalysisofSpaceTrusses*' WRITE(12,*)'*SchoolofCivilEngineeringCSU*' WRITE(12,*)'*DesignedByMuZhaoxiang*' WRITE(12,*)'*****************************************' WRITE(12,*)''WRITE(12,*)'*****************************************'WRITE(12,*)'*************TheInputData****************'WRITE(12,*)'*****************************************' WRITE(12,100) READ(11,*)NF,NP,NE,NM,NR,NCF,ND WRITE(12,110)NF,NP,NE,NM,NR,NCF,ND100FORMAT(6X,'TheGeneralInformation'/2X,'NF',5X,'NP',5X,'NE',5X,'NM',5X,'NR',&5X,'NCF',5X,'ND')110FORMAT(2X,I2,6I7) NPF=NF*NP NDF=ND*NF CALLANALYSE(NF,NP,NE,NM,NR,NCF,ND,NPF,NDF) END!********************************************************************!總計算程序SUBROUTINEANALYSE(NF,NP,NE,NM,NR,NCF,ND,NPF,NDF)DIMENSIONX(NP),Y(NP),Z(NP),MM(NE),ME(ND,NE),IT(NF,NP),RR(ND,NR),NAE(NE),&AE(1,2),PF(4,NCF),LMT(NDF,NE),MAXA(NPF),CKK(1000),V(NPF),DIST(NPF),&PP(NPF),FF(NPF),SG(NE),SM(NE) READ(11,*)(X(I),Y(I),Z(I),I=1,NP) READ(11,*)(MM(I),ME(1,I),ME(2,I),NAE(I),I=1,NE) READ(11,*)(RR(1,J),RR(2,J),J=1,NR) READ(11,*)(AE(1,J),J=1,2) WRITE(12,120) WRITE(12,121)(I,X(I),Y(I),Z(I),I=1,NP) WRITE(12,130) WRITE(12,131)(MM(I),ME(1,I),ME(2,I),NAE(I),I=1,NE) WRITE(12,140) WRITE(12,141)(INT(RR(1,J)),RR(2,J),J=1,NR) WRITE(12,150) WRITE(12,151)(AE(1,J),J=1,2) IF(NCF/=0)THEN READ(11,*)((PF(I,J),I=1,4),J=1,NCF) WRITE(12,160) WRITE(12,161)(INT(PF(1,J)),PF(2,J),PF(3,J),PF(4,J),J=1,NCF) ENDIF120FORMAT(/6X,'TheInformationofJoints'/2x,'Joint',5X,'X',5X,'Y',5X,'Z')121FORMAT(1X,I4,130FORMAT(/6X,'TheInformationofMembers'/2x,'Member',2X,'START',4X,'END',6X,'NAE')131FORMAT(1X,I4,3I8)140FORMAT(/6X,'TheInformationofSUPPORTS'/2x,'Joint',5X,'S')141FORMAT(1X,I4,150FORMAT(/6X,'TheInformationofSections'/4x,'E0',8X,'A0')151FORMAT(1X,,160FORMAT(/6X,'TheLoadingatJoints'/2x,'Joint',5X,'FX',5X,'FY',7X,'FZ')161FORMAT(1X,I4, CALLFLMT(NP,NE,NN,NNM,NR,RR,ND,NF,NDF,ME,IT,LMT) CALLFMAXA(NNM,NE,LMT,MAXA,NWK,NPF,NDF) CALLLP(V,PP,IT,PF,NN,NCF,NF,NP,NPF) CALLCONKB(NP,NE,NM,NWK,ME,X,Y,Z,AE,NAE,LMT,MAXA,CKK,NNM) ISH=1 CALLLDLT(CKK,MAXA,NN,ISH,IOUT,NWK,NNM) CALLREBACK(CKK,V,MAXA,NN,NWK,NNM) CALLDISPLS(NP,NE,NPF,NM,NN,IT,V,DIST,AE,NAE,X,Y,Z,PP,FF,SG,SM,ME,NR,RR,NF)END!********************************************************************!矩陣轉(zhuǎn)置子程序SUBROUTINEMAT(M,N,A,B)DIMENSIONA(M,N),B(N,M)DOI=1,MDOJ=1,NB(J,I)=A(I,J)ENDDOENDDORETURNEND!單元剛度矩陣的形成SUBROUTINEFKE(NP,NE,NM,IE,X,Y,Z,ME,NAE,AE,AKE)DIMENSIONX(NP),Y(NP),Z(NP),ME(2,NE),NAE(NE),AE(2,NM),AKE(2,2)N1=ME(1,IE)N2=ME(2,IE)X1=X(N1);Y1=Y(N1);Z1=Z(N1)X2=X(N2);Y2=Y(N2);Z2=Z(N2)BL=SQRT((X2-X1)**2+(Y2-Y1)**2+(Z2-Z1)**2)NMI=NAE(IE)E0=AE(1,NMI);A0=AE(2,NMI)C=E0*A0/BLAKE(1,1)=CAKE(1,2)=-CAKE(2,1)=-CAKE(2,2)=CRETURNEND!單元坐標轉(zhuǎn)換矩陣SUBROUTINEFT(IE,NP,NE,X,Y,Z,ME,T)DIMENSIONX(NP),Y(NP),Z(NP),ME(2,NE),T(2,6)T=0N1=ME(1,IE);N2=ME(2,IE)X1=X(N1);Y1=Y(N1);Z1=Z(N1)X2=X(N2);Y2=Y(N2);Z2=Z(N2)BL=SQRT((X2-X1)**2+(Y2-Y1)**2+(Z2-Z1)**2)CX=(X2-X1)/BLCY=(Y2-Y1)/BLCZ=(Z2-Z1)/BLT(1,1)=CX;T(2,4)=CXT(1,2)=CY;T(2,5)=CYT(1,3)=CZ;T(2,6)=CZRETURNEND!生成單元聯(lián)系數(shù)組LMTSUBROUTINEFLMT(NP,NE,NN,NNM,NR,RR,ND,NF,NDF,ME,IT,LMT)DIMENSIONIT(NF,NP),LMT(NDF,NE),ME(ND,NE),RR(2,NR)NN=0;NNM=0;IT=0;LMT=0N=0DOI=1,NPC=0DOK=1,NRKR=RR(1,K)IFC=RR(2,K)ENDDONC=C !NC=0,提取了整數(shù)部分C=C-NC!C=0.***,例如C=DOJ=1,NFC=C*!例如C=L=C+!提取C整數(shù)部分,例如L=1,即提取了約束RR(2,K)十分位!上的數(shù)字,這里"+"是為了防止四舍五入是出現(xiàn)錯誤C=C-LIFN=N+1IT(J,I)=NELSEIT(J,I)=0ENDIFENDDOENDDONN=NNNM=NN+1DOIE=1,NEDOI=1,NDNI=ME(I,IE)DOJ=1,NFLMT((I-1)*NF+J,IE)=IT(J,NI)ENDDOENDDOENDDORETURNEND!二維總剛中對角線元地址數(shù)組SUBROUTINEFMAXA(NNM,NE,LMT,MAXA,NWK,NPF,NDF)DIMENSIONMAXA(NPF),LMT(NDF,NE)MAXA=0;NWK=0MAXA(1)=1DOI=2,NNMIP=I-1IG=IPDOIE=1,NEDOJ=1,NDFIF(LMT(J,IE).THENDOK=1,NDFIF(LMT(K,IE).IG=LMT(K,IE)ENDDOENDIFENDDOENDDOMAXA(I)=MAXA(I-1)+IP-IG+1ENDDONWK=MAXA(NNM)-1RETURNEND!生成一維存儲結(jié)構(gòu)總剛度矩陣SUBROUTINECONKB(NP,NE,NM,NWK,ME,X,Y,Z,AE,NAE,LMT,MAXA,CKK,NNM)DIMENSIONCKK(NWK),X(NP),Y(NP),Z(NP),AE(2,NM),NAE(NE),LMT(6,NE),ME(2,NE),&MAXA(NNM),AK(6,2),AKE(2,2),T(2,6),TT(6,2),TAK(6,6)CKK=0DO10IE=1,NETAK=0CALLFKE(NP,NE,NM,IE,X,Y,Z,ME,NAE,AE,AKE)CALLFT(IE,NP,NE,X,Y,Z,ME,T)CALLMAT(2,6,T,TT)AK=MATMUL(TT,AKE)TAK=MATMUL(AK,T)!總體坐標系下的單元剛度矩陣DO220I=1,6DO220J=1,6NI=LMT(I,IE)NJ=LMT(J,IE)IF((NJ-NI).THENIJ=MAXA(NJ)+NJ-NICKK(IJ)=CKK(IJ)+TAK(I,J)ENDIF220CONTINUE10CONTINUERETURNEND!生成荷載矩陣SUBROUTINELP(V,PP,IT,PF,NN,NCF,NF,NP,NPF)DIMENSIONV(NN),PP(NPF),IT(NF,NP),PF(4,NCF)V=0PP=0DOI=1,NFDOJ=1,NP DOK=1,NCF IF V(IT(I,J))=PF(I+1,K) ENDIF ENDDO ENDDOENDDODOK=1,NCFDOI=1,NPIF PP(NF*(I-1)+1)=PF(2,K)PP(NF*(I-1)+2)=PF(3,K) PP(NF*(I-1)+3)=PF(4,K)ENDIFENDDOENDDORETURNEND!對一維結(jié)構(gòu)總剛度矩陣進行矩陣分解（LDLT）SUBROUTINELDLT(A,MAXA,NN,ISH,IOUT,NWK,NNM)DIMENSIONA(NWK),MAXA(NNM)IFRETURNDO200N=1,NNKN=MAXA(N)KL=KN+1KU=MAXA(N+1)-1KH=KU-KLIF(KH)304,240,210210K=N-KHIC=0KLT=KUDO260J=1,KHKLT=KLT-1IC=IC+1KI=MAXA(K)ND=MAXA(K+1)-KI-1IF(ND)260,260,270270KK=MIN0(IC,ND)C=DO280L=1,KK280C=C+A(KI+L)*A(KLT+L)A(KLT)=A(KLT)-C260K=K+1240K=NB=DO300KK=KL,KUK=K-1KI=MAXA(K)C=A(KK)/A(KI)IF(ABS(C).GOTO290WRITE(IOUT,2010)N,CSTOP290B=B+C*A(KK)300A(KK)=CA(KN)=A(KN)-B304IF(A(KN))310,310,200310IFGOTO320IF(A(KN).A(KN)=GOTO200320WRITE(IOUT,2000)N,A(KN)STOP200CONTINUERETURN2000 FORMATxt31324 16120-1050-1025-10-250-10-50-10-25-1025-2-2 25-20-2 -2-25-200 0 01 1 7 12 1 8 13 2 8 14 2 9 15 3 9 16 3 10 17 4 10 18 4 11 19 5 11 110 5 12 111 6 12 112 6 7 113 7 8 114 8 9 115 9 10 116 10 11 117 11 12 118 7 12 119 7 13 120 8 13 121 9 13 122 10 13 123 11 13 124 12 13 11 2 3 4 5 6 210E6 13 0 -500 0②輸出文件.txt******************************************ProgramforAnalysisofSpaceTrusses**SchoolofCivilEngineeringCSU**DesignedByMuZhaoxiang************************************************************************************************TheInputData*********************************************************TheGeneralInformationNFNPNENMNRNCFND313241612TheInformationofJointsJointXYZ1.0234.056789.0101112.013.0.0.0TheInformationofMembersMemberSTARTENDNAE117121813281429153916310174101841119511110512111612112671137811489115910116101111711121187121197131208131219131221013123111312412131TheInformationofSUPPORTSJointS1.1112.1113.1114.1115.1116.111TheInformationofSectionsE0A0+08.0400TheLoadingatJointsJointFXFYFZ13.00.00*************************************************TheResultsofCalculation*********************************************