实验四空间桁架程序设计.doc

实验四 空间桁架程序设计开 始标题及数组说明组集整体刚度矩阵计算并打印各杆轴力解方程并打印杆端位移形成结点载荷列阵结 束子程序READ子程序CSL子程序STIFFN输入并打印原始数据子程序MULV6子程序CALM子程序TARNS子程序SOLV一 空间桁架内力和位移计算总框图二 平面桁架静力分析源程序（STSAP.FOR）C ANALYSIS PROGRAM FOR SPACE TRUSS DIMENSION JE(2,100),JEA(50),JN(3,100),M(6),EA(2,50), & JPJ(100),PJ(100),P(200),FE(6),T(6,6),AKE(6,6),F(6), & AK(200,200),X(100),Y(100),Z(100),FF(100) OPEN(6,FILE=STSAP.IN,STATUS=OLD) OPEN(8,FILE=STSAP.OUT,STATUS=NEW) CALL READ(NJ,N,NNE,NMT,NPJ,JN,X,Y,Z,JE,JEA,EA,JPJ,PJ)10 DO 20 J=1,N DO 20 I=1,N20 AK(I,J)=0. DO 25 IE=1,NNE CALL CSL(CX,CY,CZ,AL,JE,IE,X,Y,Z) CALL STIFFN(AKE,JEA(IE),CX,CY,CZ,AL,EA) CALL CALM(M,IE,JE,JN) DO 25 I=1,6 IF(M(I).EQ.0)GOTO 25 DO 23 J=I,6 IF(M(J).EQ.0)GOTO 23 I1=MIN0(M(I),M(J) J1=MAX0(M(I),M(J) AK(I1,J1)=AK(I1,J1)+AKE(I,J)23 CONTINUE25CONTINUE DO 34 I=1,N34 P(I)=0. IF(NPJ.EQ.0)GOTO 39 DO 35 I=1,NPJ I1=JPJ(I)35 P(I1)=P(I1)+PJ(I)39 CALL SLOV(AK,N,P) WRITE(8,40)40 FORMAT(/10X,10(1H*), RESULTS OF CALCULATION ,10(1H*), & /5X,NO.N,4X,X-DISPLACEMENT,2X, & Y-DISPLACEMENT,3X,Z-DISPLACEMENT) DO 53 K=1,NJ DO 52 I=1,3 F(I)=0. I1=JN(I,K)52 IF(I1.GT.0)F(I)=P(I1)53 WRITE(8,60)K,F(1),F(2),F(3)60 FORMAT(I8,2X,3G16.5) WRITE(8,62)62 FORMAT(/5X,NO.E,5X,AXIAL FORCE,10X,NO.E,5X,AXIAL FORCE) DO 70 IE=1,NNE CALL CSL(CX,CY,CZ,AL,JE,IE,X,Y,Z) CALL CALM(M,IE,JE,JN) CALL STIFFN(AKE,JEA(IE),CX,CY,CZ,AL,EA) DO 65 I=1,6 L=M(I) F(I)=0.65 IF(L.GT.0)F(I)=P(L) CALL MULV6(FE,AKE,F) CALL TRANS(T,CX,CY,CZ) CALL MULV6(F,T,FE)70 FF(IE)=-F(1) WRITE(8,80)(I,FF(I),I=1,NNE)80 FORMAT(2(3X,I5,7X,G12.5,2X) STOP END SUBROUTINE READ(NJ,N,NNE,NMT,NPJ,JN,X,Y,Z,JE,JEA,EA,JPJ,PJ) DIMENSION JE(2,100),JEA(50),JN(3,100),EA(2,50), & JPJ(100),PJ(100),X(100),Y(100),Z(100),TL(20) READ(6,10) (TL(I),I=1,20)10 FORMAT(20A4) READ(6,20)NJ,N,NNE,NMT,NPJ20 FORMAT(5I5) WRITE(8,10)TL WRITE(8,22)NJ,N,NNE,NMT,NPJ22 FORMAT(/10X,NUMBER OF NODES =,I5, & /10X,NUMBER OF DEGREE OF FREEDOM =,I5, & /10X,NUMBER OF ELEMENTS =,I5, & /10X,NUMBER OF MATERIALS =,I5, & /10X,NUMBER OF NODAL LOADS =,I5) READ(6,30) (JN(J,I),J=1,3),X(I),Y(I),Z(I),I=1,NJ) WRITE(8,40)(I,(JN(J,I),J=1,3),X(I),Y(I),Z(I),I=1,NJ)30 FORMAT(2(3I5,3F10.0)40 FORMAT(/3X,NO.N,2X,NO.DISP.(X,Y,Z.),8X, & X-COORDINATE Y-COORDINATE Z-COORDINATE/(4I6,9X,3G14.6) READ(6,50)(JE(J,I),J=1,2),JEA(I),I=1,NNE)50 FORMAT(15I5) WRITE(8,60)(I,(JE(J,I),J=1,2),JEA(I),I=1,NNE)60 FORMAT(/8X,NO.E,5X,1(NODE),3X,2(NODE),3X & ,NO.MAT/4(I10) READ(6,70)(EA(I,J),I=1,2),J=1,NMT)70 FORMAT(6f10.0) WRITE(8,80)(J,(EA(I,J),I=1,2),J=1,NMT)80 FORMAT(/3X,NO.MAT,6X,ELASTIC MODULUS,8X, & AREA/(I6,9X,2G16.6) IF(NPJ.EQ.0)GOTO 100 READ(6,90) (JPJ(I),PJ(I),I=1,NPJ)90 FORMAT(5(I5,F10.0) WRITE(8,95)(JPJ(I),PJ(I),I=1,NPJ)95 FORMAT(/20X,NODAL LOAD(NO.DISP.,VALUE)/4(I6,G12.5)100 RETURN END SUBROUTINE CSL(CX,CY,CZ,AL,JE,IE,X,Y,Z) DIMENSION JE(2,1), X(1),Y(1),Z(1) I=JE(1,IE) J=JE(2,IE) S1=X(J)-X(I) S2=Y(J)-Y(I) S3=Z(J)-Z(I) AL=SQRT(S1*S1+S2*S2+S3*S3) CX=S1/AL CY=S2/AL CZ=S3/AL RETURN END SUBROUTINE STIFFN (AKE,I1,CX,CY,CZ,AL,EA) DIMENSION EA(2,1),AKE(6,6) CX2=CX*CX CY2=CY*CY CZ2=CZ*CZ B1=EA(1,I1)*EA(2,I1)/AL AKE(1,1)=CX2*B1 AKE(1,2)=Cx*CY*B1 AKE(1,3)=CX*CZ*B1 AKE(2,1)=AKE(1,2) AKE(2,2)=CY2*B1 AKE(2,3)=CY*CZ*B1 AKE(3,1)=AKE(1,3) AKE(3,2)=AKE(2,3) AKE(3,3)=CZ2*B1 AKE(4,1)=-AKE(1,1) AKE(4,2)=-AKE(1,2) AKE(4,3)=-AKE(1,3) AKE(5,1)=-AKE(2,1) AKE(5,2)=-AKE(2,2) AKE(5,3)=-AKE(2,3) AKE(6,1)=-AKE(3,1) AKE(6,2)=-AKE(3,2) AKE(6,3)=-AKE(3,3) DO 10 II=1,3 DO 10 JJ=1,3 AKE(II+3,JJ+3)=AKE(II,JJ)10 AKE(II,JJ+3)=AKE(JJ+3,II) RETURN END SUBROUTINE CALM(M,IE,JE,JN) DIMENSION M(6),JE(2,1),JN(3,2) DO 10 I=1,3 M(I)=JN(I,JE(1,IE)10 M(I+3)=JN(I,JE(2,IE) RETURN END SUBROUTINE SLOV(AK,N,P) DIMENSION P(200),AK(200,200) DO 40 K=1,N-1 DO 40 I=K+1,N CX=-AK(K,I)/AK(K,K) P(I)=P(I)+CX*P(K) DO 40 J=I,N40 AK(I,J)=AK(I,J)+CX*AK(K,J) P(N)=P(N)/AK(N,N) DO 50 I=N-1,1,-1 DO 48 J=I+1,N48 P(I)=P(I)-AK(I,J)*P(J)50 P(I)=P(I)/AK(I,I) RETURN END SUBROUTINE MULV6 (U,A,V) DIMENSION U(6),A(6,6),V(6) DO 10 I=1,6 U(I)=0. DO 10 J=1,610 U(I)=U(I)+A(I,J)*V(J) RETURN END SUBROUTINE TRANS (T,CX,CY,CZ) DIMENSION T(6,6) DO 10 J=1,6 DO 10 I=1,610 T(I,J)=0. T(1,1)=CX T(1,2)=CY T(1,3)=CZ T(4,4)=CX T(4,5)=CY T(4,6)=CZ RETURN END三 课后习题（1） 准备原始数据确定结点、划分单元，建立整体坐标系与局部坐标系如图所示。4m设本题标题为EXAMPLE(6.1)，根据题目所给的数据填写下表：空间桁架的输入数据表标题EXAMPLE(6.1)基本数据结点总数自由度单元总数单元类型数结点载荷数43311结点数据结点号位移号坐标结点号位移号坐标10, 0, 00.，0.，2.20, 0, 00.，0.，-2.30, 0, 00.，4.，0.41,2,34.，0.，0.单元数据单元号始点号终结点类型单元号始点号终结点类型141241341单元类型数据类型号弹性模量横截面积类型号弹性模量横截面积11.5E84.0E-4结点载荷数据位移号数值位移号数值位移号数值位移号数值位移号数值2-1.0（2） 建立数据文件STSAP.IN并输入以下数据：EXAMPLE(6.1)4,3,3,1,10,0,0,0.,0.,2.,0,0,0,0.,0.,-2.0,0,0,0.,4.,0.,1,2,3,4.,0.,0.1,4,1,2,4,1,3,4,11.5E8,4.E-42,-1.（3） 运行程序。从文件STSAP.OUT中得到如下结果： EXAMPLE(6.1) NUMBER OF NODES = 4 NUMBER OF DEGREE OF FREEDOM = 3 NUMBER OF ELEMENTS = 3 NUMBER OF MATERIALS = 1 NUMBER OF NODAL LOADS = 1 NO.N NO.DISP.(X,Y,Z.) X-COORDINATE Y-COORDINATE Z-COORDINATE 1 0 0 0 .000000 .000000 2.00000 2 0 0 0 .000000 .000000 -2.00000 3 0 0 0 .000000 4.00000 .000000 4 1 2 3 4.00000 .000000 .000000 NO.E 1(NODE) 2(NODE) NO.MAT 1 1 4 1 2 2 4 1 3 3 4 1 NO.MAT ELASTIC MODULUS AREA 1 .150000E+09 .400000E-03 NODAL LOAD(NO.DI