【正文】
BLC C1=C/BL C2=C1*C1 C 集中力作用下的等效荷截列陣 P0(2) 20 F0(1)=*C*D*D/BL/BL F0(2)=*D*C2 END C ****************************************** c 計算第 NE0個單元的單元剛度矩陣 DK(2,2) C ****************************************** SUBROUTINE DG(NE0,NE,E0,GC,GX,DK) DIMENSION GC(NE),GX(NE),DK(2,2) DO 15 I=1,2 DO 10 J=1,2 DK(I,J)= 10 CONTINUE 15 CONTINUE C DL—桿長 ,DI慣性矩 ,S線剛度 DL=GC(NE0) DI=GX(NE0) S=E0*DI/DL DK(1,1)=*S DK(1,2)=*S DK(2,1)=*S DK(2,2)=*S END c ********************************************************* C 集成總體剛度矩陣 ,ZKl,ZK2分別存主對角元素和付對角元素 c ********************************************************* SUBROUTINE JCZK(NE,NJ,E0,ZK1,ZK2,DK,GC,GX) DIMENSION ZK1(NJ),ZK2(NJ),DK(2,2),GC(NE),GX(NE) DO 10 I=1,NJ ZK1(I)= ZK2(I)= 10 CONTINUE DO 100 I=1,NE CALL DG(I,NE,E0,GC,GX,DK) ZK1(I)=ZK1(I)+DK(1,1) ZK2(I)=ZK2(I)+DK(1,2) ZK1(I+1)=ZK1(I+1)+DK(2,2) 100 CONTINUE END C ***************************** C 形成總荷裁矩陣 C ***************************** SUBROUTINE XCP(NJ,NE,P,PF,F0,GC) DIMENSION P(NJ),PF(1,2),F0(2),GC(NE) DO 10 I=1,NJ P(I)= 10 CONTINUE CALL DJH(NE,PF,F0,GC) ND=INT(PF(1,2)+) P(ND)=P(ND)+F0(1) P(ND+1)=P(ND+1)+F0(2) END c ********************************** c 進行支承條件處理 c ********************************** SUBROUTINE ZCCL(NZ,NJ,ZK1,ZK2,P) DIMENSION ZK1(NJ),ZK2(NJ),P(NJ) GOTO(10,20,30,20),NZ 10 GOTO 100 20 ZK1(1)= P(1)= ZK2(1)= IF() GOTO 30 GOTO 100 30 ZK1(NJ)= ZK2(NJ1)= P(NJ)= 100 CONTINUE END c ****************************** c 解方程求節(jié)點位移 P c ****************************** SUBROUTINE FCQJ(NJ,ZK1,ZK2,P) DIMENSION ZK1(NJ),ZK2(NJ),P(NJ) DO 10 I=1,NJ1 ZK1(I+1)=ZK1(I+1)ZK2(I)*ZK2(I)/ZK1(I) P(I+1)=P(I+1)ZK2(I)*P(I)/ZK1(I) 10 CONTINUE P(NJ)=P(NJ)/ZK1(NJ) DO 20 I=1,NJ1 P(NJI)=(P(NJI)ZK2(NJI)*P(NJI+1))/ZK1(NJI) 20 CONTINUE END c ******************************* c 輸出位移 c ****************************** c 計算單元桿端內(nèi)力 c ***************** SUBROUTINE GDL(NE,NJ,E0,P,WY,F,F0,DK,PF,GC,GX,gd) DIMENSION P(NJ),WY(2),F(2),F0(2),DK(2,2),PF(1,2) DIMENSION GC(NE),GX(NE),gd(ne,2) DO 200 NE0=1,NE CALL DG (NE0,NE,E0,GC,GX,DK) WY(1)=P(NE0) WY(2)=P(NE0+1) DO 35 I=1,2 F(I)= DO 30 J=1,2 F(I)=F(I)+DK(I,J)*WY(J) 30 CONTINUE 35 CONTINUE IF(INT(PF(1,2)+).) THEN CALL DJH(NE,PF,F0,GC) DO 40 J=1,2 F(J)=F(J)F0(J) 40 CONTINUE END IF 50 CONTINUE gd(ne0,1)=f(1) gd(ne0,2)=f(2) 200 CONTINUE end c ********************************************* c 計算任意截面的彎矩和剪力 c ************************ subroutine jmnl(gd,sqjm,njm,qdnl,gc,ne,pf) dimension gd(ne,2),sqjm(njm,2),gc(ne),qdnl1(1,njm,2),pf(1,2),qdnl(, 1njm,2) c qdnl(njm,1)=彎矩 qdnl(njm,2)=剪力 c 取左側(cè)梁為隔離體,彎矩逆時針為正,剪力向上為正, do njm0=1,njm dyh=sqjm(njm0,1) qdnl(njm0,2)=(gd(dyh,2)+gd(dyh,1))/gc(dyh) qdnl(njm0,1)=(gd(dyh,1)+qdnl(njm0,2)*sqjm(njm0,2)) c=pf(1,1) d=gc(dyh)c if((abs(pf(1,2)dyh)).and.(pf(1,1)=sqjm(njm0,2))) then qdnl1(1,njm0,2)=*d/gc(dyh) qdnl1(1,njm0,1)=qdnl1(1,njm0,2)*sqjm(njm0,2) else if((abs(pf(1,2)dyh)).and.(pf(1,1)sqjm(njm0,2))) then if(d)then qdnl1(1,njm0,2)=0 qdnl1(1,njm0,1)=*d*c/gc(dyh)qdnl1(1,njm0,2)*(sqjm(njm0,2)c) else qdnl1(1,njm0,2)=*c/gc(dyh) qdnl1(1,njm0,1)=*d*c/gc(dyh)qdnl1(1,njm0,2)*(sqjm(njm0,2)c) end if else qdnl1(1,njm0,2)=0 qdnl1(1,njm0,1)=0 end if qdnl(njm0,1)=qdnl(njm0,1)+qdnl1(1,njm0,1) qdnl(njm0,2)=qdnl(njm0,2)+qdnl1(1,njm0,2) write(9,30) njm0,qdnl(njm0,1),qdnl(njm0,2) 30 format(1x,39。截面 39。,i3,2x,39。彎矩 =39。,10x,39。剪力 =39。,) end do end