📄 plastic_beam_fencheng_module_内部文件.for

📁 欧文的程序
💻 FOR
📖 第 1 页 / 共 2 页
字号:
12 下一页

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!	Timeshenko分层弹塑性梁单元材料非线性计算						  !
!	 材料采用弹塑性，混凝土材料同样采用改种材料						  !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!



$DEBUG
!	NPOIN(结点总数),NELEM(单元总数),NBOUN(边界点总数),NLAYR(分层数),NLOAD(),
!				NPROP(完全确定材料需要的参数),NNODE(每个单元的结点数),IINCS(),
!     IITER,KRESL,NCHEK,TOLER,NALGO(用以识别解法的指示变量),NSVAB,NDOFN(每个结点的自由度数),
!     NINCS(施加荷载增量的数目),NEVAB,NITER,NOUTP,FACTO,PVALU

!     IEVAB(结点自由度变量)


!	PROPS(5,25)(材料数组),COORD(26)(结点坐标),LNODS(25,2)(单元拓扑信息),IFPRE(52),
!     FIXED(52),TLOAD(25,4),RLOAD(25,4),ELOAD(25,4),
!     MATNO(25)(每个单元对应的材料号),STRES(25,2),PLAST(250),XDISP(52),
!     TDISP(26,2),TREAC(26,2),ASTIF(52,52),ASLOD(52),
!     REACT(52),FRESV(1352),PEFIX(52),ESTIF(4,4),
!     STRSL(250,2)
!*============================================================*!　
!*												 　			*!
!*============================================================*! 	

	MODULE CONFIG_ARRAY
	IMPLICIT NONE
	INTEGER::NPOIN,NELEM,NBOUN,NLAYR,NLOAD,NPROP,NNODE,IINCS,
     $			IITER,KRESL,NCHEK,NALGO,NSVAB,NDOFN,
     $			NINCS,NEVAB,NITER,NOUTP,NMATS
      DOUBLEPRECISION::FACTO,PVALU,TOLER


	DOUBLEPRECISION,ALLOCATABLE::PROPS(:,:),COORD(:),TLOAD(:,:) 
	DOUBLEPRECISION,ALLOCATABLE::RLOAD(:,:),ELOAD(:,:)
	DOUBLEPRECISION,ALLOCATABLE::STRES(:,:),PLAST(:),XDISP(:)
	DOUBLEPRECISION,ALLOCATABLE::TDISP(:,:),TREAC(:,:)
	DOUBLEPRECISION,ALLOCATABLE::ASTIF(:,:),ASLOD(:),REACT(:)
	DOUBLEPRECISION,ALLOCATABLE::FRESV(:),PEFIX(:)
	DOUBLEPRECISION,ALLOCATABLE::STRSL(:,:)
	DOUBLEPRECISION::ESTIF(4,4)
	INTEGER,ALLOCATABLE::LNODS(:,:),IFPRE(:),FIXED(:),MATNO(:)

	END MODULE CONFIG_ARRAY

!*============================================================*!　
!*												 　			*!
!*============================================================*! 	


	PROGRAM BEAM_FENGCHENG
	USE CONFIG_ARRAY
	IMPLICIT NONE
	CHARACTER PN*40,FN_1*12,FN_5*12,FN_6*12,FN_7*12
	WRITE(*,'(/A)')' 钢筋混凝土平面Timeshenko分层梁单元'
      write(*,'(A)')'	         承载能力计算程序'
	WRITE(*,'(/A)') ' 输入计算问题名(PN):'
	READ(*,'(A)')PN

	CALL FNAME(PN,'.DAT',FN_5)
	WRITE(*,'(/2A)') '  输入数据文件名为 ',FN_5
	CALL FNAME(PN,'.OUT',FN_6)
	WRITE(*,'(/2A)') '  输入数据文件名为 ',FN_6
	CALL FNAME(PN,'.PFX',FN_7)
	WRITE(*,'(/2A)') '  荷载挠度关系输出数据文件名为 ',FN_7
 	OPEN(1,STATUS='NEW',FORM='UNFORMATTED')
	OPEN(5,FILE=FN_5)
	OPEN(6,FILE=FN_6)
	OPEN(7,FILE=FN_7)
	CALL DATAIN
	CALL INITAL
	DO IINCS=1,NINCS
		CALL INCLOD
		DO IITER=1,NITER
			WRITE(*,*)'荷载步:' ,IINCS,'	迭代步:', IITER
			CALL NONAL
			IF (KRESL.EQ.1) CALL STIFBL
			CALL ASSEMB
			IF(KRESL.EQ.1)CALL GREDUC
			IF(KRESL.EQ.2) CALL RESOLV
			CALL BAKSUB
			CALL REFORBL
			CALL CONUND
			IF (NCHEK.EQ.0) GOTO 20
			IF(IITER.EQ.1.AND.NOUTP.EQ.1)CALL RESULT
			IF(NOUTP.EQ.2)CALL RESULT
		ENDDO
		WRITE(6,900)
900		FORMAT(1HO,5X,'SOLUTION NOT COVERGED')
		STOP
20		CALL RESULT
	ENDDO

	CLOSE(5)
	CLOSE(6)
	CLOSE(7)
	DEALLOCATE(PROPS,COORD,TLOAD )
	DEALLOCATE(RLOAD,ELOAD)
	DEALLOCATE(STRES,PLAST,XDISP)
	DEALLOCATE(TDISP,TREAC)
	DEALLOCATE(ASTIF,ASLOD,REACT)
	DEALLOCATE(FRESV,PEFIX)
	DEALLOCATE(STRSL)
	DEALLOCATE(LNODS,IFPRE,FIXED,MATNO)

	STOP
	END	PROGRAM BEAM_FENGCHENG

!*============================================================*!　
!*												 　			*!
!*============================================================*! 	
	SUBROUTINE DATAIN
	USE CONFIG_ARRAY
	IMPLICIT NONE
	INTEGER::IMATS,JMATS,IPROP,IELEM,JELEM,INODE
	INTEGER::IPOIN,JPOIN,ISVAB,IBOUN,NODFX,IDOFN,NPOSN,IEVAB
	INTEGER::ICODE(2)
	DOUBLEPRECISION::VALUE(2) 
!	CHARACTER TITLE*140 
!
!	READ(5,*) TITLE
!	WRITE(6,*) TITLE
	READ(5,*) NPOIN,NELEM,NBOUN,NMATS,NPROP,NNODE,NINCS,
     $NALGO,NDOFN,NLAYR
	WRITE(6,905)NPOIN,NELEM,NBOUN,NMATS,NPROP,NNODE,
     $NINCS,NALGO,NDOFN,NLAYR
905	FORMAT(//1X,'NPOIN=',I5,3X,'NELEM=',I5,3X,'NBOUN=',I5,3X,
     $		    'NMATS=',I5//1X,'NPROP=',I5,3X,'NNODE=',I5,3X,
     $			'NINCS=',I5,3X,'NALGO=',I5//1X,'NDOFN=',I5,3X,
     $			'NLAYR=',I5)
	NEVAB=NDOFN*NNODE
	NSVAB=NDOFN*NPOIN

	ALLOCATE (PROPS(NMATS,NPROP))
	ALLOCATE(LNODS(NELEM,NNODE),MATNO(NELEM))
	ALLOCATE(COORD(NPOIN))
	ALLOCATE(IFPRE(NSVAB),PEFIX(NSVAB),FIXED(NSVAB))
	ALLOCATE(RLOAD(NELEM,NEVAB))
	ALLOCATE(PLAST(NLAYR*NELEM),STRSL(NLAYR*NELEM,2))
	ALLOCATE(STRES(NELEM,NDOFN))
	ALLOCATE(ELOAD(NELEM,NEVAB),TLOAD(NELEM,NEVAB))
	ALLOCATE(TDISP(NPOIN,NDOFN),TREAC(NPOIN,NDOFN))
	ALLOCATE(ASLOD(NSVAB),ASTIF(NSVAB,NSVAB))
	ALLOCATE(FRESV((NSVAB-1)*NSVAB/2))
	ALLOCATE(REACT(NSVAB),XDISP(NSVAB)) 



	WRITE(6,910)
910	FORMAT(5X,'MATERIAL PROPERTIES')
 
	DO IMATS=1,NMATS
		READ(5,*) JMATS,(PROPS(JMATS,IPROP),IPROP=1,NPROP)
		WRITE(6,915)JMATS,(PROPS(JMATS,IPROP),IPROP=1,NPROP)
	ENDDO
915	FORMAT(I10,4F15.5/10X,4F15.5/10X,4F15.5/10X,4F15.5/10X,4F15.5)
	WRITE(6,920)
920	FORMAT(3X,'EL NODFS MAT.')

	DO IELEM=1,NELEM
	 READ(5,*) JELEM,(LNODS(JELEM,INODE),INODE=1,NNODE),MATNO(JELEM)
	 WRITE(6,925)JELEM,(LNODS(JELEM,INODE),INODE=1,NNODE),MATNO(JELEM)
	ENDDO
925	FORMAT(4I5)
	WRITE(6,930)
930	FORMAT(5X,'NODE',5X,'COORD.')

	DO IPOIN=1,NPOIN
		READ(5,*)JPOIN,COORD(JPOIN)
		WRITE(6,935)JPOIN,COORD(JPOIN)
	ENDDO
935	FORMAT(I10,F15.5)
	
	DO ISVAB=1,NSVAB
		IFPRE(ISVAB)=0.
		PEFIX(ISVAB)=0.0
	ENDDO

	IF(NDOFN.EQ.1)WRITE(6,940)
940	FORMAT(1X,'RES.NODE',2X,'CODE',3X,'PRES.VALUES')
	IF(NDOFN.EQ.2)WRITE(6,945)
945	FORMAT(1X,'RES.NODE',2X,'CODE',3X,'PRES.VALUES',2X,
     $		'CODE',3X,'PRES.VALUES')

	DO IBOUN=1,NBOUN
		READ(5,*) NODFX,(ICODE(IDOFN),VALUE(IDOFN),IDOFN=1,NDOFN)
		WRITE(6,950)NODFX,(ICODE(IDOFN),VALUE(IDOFN),IDOFN=1,NDOFN)
950	FORMAT(I10,2(I5,F15.5))
		NPOSN=(NODFX-1)*NDOFN
		DO IDOFN=1,NDOFN
			NPOSN=NPOSN+1
			IFPRE(NPOSN)=ICODE(IDOFN)
 			PEFIX(NPOSN)=VALUE(IDOFN)
		ENDDO
	ENDDO
	WRITE(6,955)
955	FORMAT(3X,'ELEMENT',10X,'NODAL LOADS')
	DO IELEM=1,NELEM
		DO IEVAB=1,NEVAB
			RLOAD(IELEM,IEVAB)=0.0
		ENDDO
	ENDDO
	
	DO IELEM=1,NELEM
 		READ(5,*)JELEM,(RLOAD(JELEM,IEVAB),IEVAB=1,NEVAB)
	ENDDO

	DO IELEM=1,NELEM
		WRITE(6,960)IELEM,(RLOAD(IELEM,IEVAB),IEVAB=1,NEVAB)
	ENDDO
960	FORMAT(I10,5F15.5)

	RETURN
	END
!*============================================================*!　
!*												 　			*!
!*============================================================*! 	
	SUBROUTINE INITAL
	USE CONFIG_ARRAY
	IMPLICIT NONE
	INTEGER::IXiaBiao,IELEM,IDOFN,IEVAB,IPOIN
	DO IXiaBiao=1,NLAYR*NELEM
		PLAST(IXiaBiao)=0.0
		STRSL(IXiaBiao,1)=0.0
		STRSL(IXiaBiao,2)=0.0
	ENDDO
	DO IELEM=1,NELEM
!		PLAST(IELEM)=0.0
		DO IDOFN=1,NDOFN
			STRES(IELEM,IDOFN)=0.0
		ENDDO
		DO IEVAB=1,NEVAB
			ELOAD(IELEM,IEVAB)=0.0
			TLOAD(IELEM,IEVAB)=0.0
		ENDDO
	ENDDO
	DO IPOIN=1,NPOIN
		DO IDOFN=1,NDOFN
			TDISP(IPOIN,IDOFN)=0.0
			TREAC(IPOIN,IDOFN)=0.0
		ENDDO
	ENDDO
	RETURN
	END
!*============================================================*!　
!*												 　			*!
!*============================================================*! 	
	SUBROUTINE NONAL
	USE CONFIG_ARRAY
	IMPLICIT NONE

	INTEGER::ISVAB
	KRESL=2
	IF(NALGO.EQ.1)KRESL=1
	IF(NALGO.EQ.2)KRESL=1
	IF(NALGO.EQ.3.AND.IINCS.EQ.1.AND.IITER.EQ.1)KRESL=1
	IF(NALGO.EQ.4.AND.IITER.EQ.1)KRESL=1
	IF(NALGO.EQ.5.AND.IITER.EQ.2)KRESL=1
	IF(IITER.EQ.1.OR.NALGO.EQ.1)GOTO 20
	DO ISVAB=1,NSVAB
		FIXED(ISVAB)=0.0
	ENDDO
	RETURN
20	DO ISVAB=1,NSVAB
		FIXED(ISVAB)=PEFIX(ISVAB)*FACTO
	ENDDO
	RETURN
	END	

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!     荷载增量
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	SUBROUTINE INCLOD
	USE CONFIG_ARRAY
	IMPLICIT NONE
	INTEGER::IELEM,IEVAB
	READ(5,*)NITER,NOUTP,FACTO,TOLER
	WRITE(6,905)IINCS,NITER,NOUTP,FACTO,TOLER
905	FORMAT(5X,'IINCS=',I5,3X,'NITER= ',I5,3X,' NOUTP=',I5,	
     $	3X,'FACTO=',E14.6,3X,'TOLER=',E14.6)
	DO IELEM=1,NELEM
		DO IEVAB=1,NEVAB
			ELOAD(IELEM,IEVAB)=ELOAD(IELEM,IEVAB)+
     $					RLOAD(IELEM,IEVAB)*FACTO
			TLOAD(IELEM,IEVAB)=TLOAD(IELEM,IEVAB)+
     $					RLOAD(IELEM,IEVAB)*FACTO
		ENDDO
	ENDDO
	RETURN
	END

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!     刚度矩阵的组集
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	SUBROUTINE ASSEMB 
	USE CONFIG_ARRAY
	IMPLICIT NONE
	INTEGER::ISVAB,JSVAB,IELEM,INODE,NODEI,NODEJ,NROWS,NROWE
	INTEGER::JNODE,NCOLS,NCOLE,I,J,IDOFN,JDOFN
 	REWIND 1
	DO ISVAB=1,NSVAB
		ASLOD(ISVAB)=0.0
	ENDDO

	IF(KRESL.EQ.2) GOTO 30
	DO ISVAB=1,NSVAB
		DO JSVAB=1,NSVAB
			ASTIF(ISVAB,JSVAB)=0.0
		ENDDO
	ENDDO
30	CONTINUE
C
C *** ASSEMBLE THE ELEMNT LAODS
C
	DO IELEM=1,NELEM
		DO I=1,4
			DO J=1,4
				ESTIF(I,J)=0.0	
			ENDDO
		ENDDO

			READ(1)ESTIF

 		DO INODE=1,NNODE
 			NODEI=LNODS(IELEM,INODE)
 			DO IDOFN=1,NDOFN
 				NROWS=(NODEI-1)*NDOFN+IDOFN
 				NROWE=(INODE-1)*NDOFN+IDOFN
 				ASLOD(NROWS)=ASLOD(NROWS)+ELOAD(IELEM,NROWE)
C
C *** ASSEMBLE THE ELEMNT STIFFNESS MATRICES 
C
 				IF(KRESL.EQ.2)GOTO 40
 				DO JNODE=1,NNODE
 					NODEJ=LNODS(IELEM,JNODE)
 					DO  JDOFN=1,NDOFN
 						NCOLS=(NODEJ-1)*NDOFN+JDOFN
 						NCOLE=(JNODE-1)*NDOFN+JDOFN
 						ASTIF(NROWS,NCOLS)=ASTIF(NROWS,NCOLS)
     $					+ESTIF(NROWE,NCOLE)
					ENDDO
				ENDDO
			ENDDO
40		ENDDO
	ENDDO
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!调试语句
!	OPEN(2,FILE='zg.dat')
!	DO I=1,NROWS
!			WRITE(2,111)(ASTIF(I,J),J=1,NCOLS)
!	ENDDO
!111	FORMAT(5X,22F25.5)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!调试语句
	RETURN
	END
!*============================================================*!　
!*												 　			*!
!*============================================================*! 	
	SUBROUTINE STIFBL
	USE CONFIG_ARRAY
	IMPLICIT NONE
	INTEGER::IELEM,LPROP,NODE1,NODE2,ISTIF,JSTIF
	DOUBLEPRECISION::EIVAL,SVALU,VALU1,VALU2,VALU3,VALU4,HARDS
	DOUBLEPRECISION::ELENG
 	REWIND 1
	DO IELEM=1,NELEM
		LPROP=MATNO(IELEM)
		CALL LAYER(IELEM,EIVAL,SVALU)
	WRITE(*,*)EIVAL,SAVLU
		HARDS=PROPS(LPROP,4)
		NODE1=LNODS(IELEM,1)
		NODE2=LNODS(IELEM,2)
		ELENG=ABS(COORD(NODE2)-COORD(NODE1))
!		IF(PLAST(IELEM).NE.0) EIVAL=EIVAL*(1.0-EIVAL/(EIVAL+HARDS))
		VALU1=0.5*SVALU
		VALU2=SVALU/ELENG
		VALU3=EIVAL/ELENG
		VALU4=0.25*SVALU*ELENG
		ESTIF(1,1)=VALU2
		ESTIF(1,2)=VALU1
		ESTIF(1,3)=-VALU2
		ESTIF(1,4)=VALU1
		ESTIF(2,2)=VALU3+VALU4
12 下一页
💿 文件大小 6 K
👤 上传用户 kyo
📂 所属分类 *行业应用
🏷️ 相关标签

#程序
⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -