平面梁单元有限元ortran程序.docVIP

program beam2d !character*64 fname1,fname2 dimension mea(100,4),jz(50,2),aai(10,2),aeu(10,2),qq(100) dimension r(6,6),rt(6,6),ra(6,6),pop(6),po(6) dimension akep(6,6),ake(6,6),p(300),x(100),y(100) dimension ld(300),is(6),a(50000),upe(6),ppe(6),ue(6) open(1,file=fname1_b.txt) open(2,file=fname2_b.txt) read(1,*) nn,ne,nm,na,nc write(*,*) NN,NE,NM,NA,NC nf=3 nd=2 nfd=nf*nd n=nn*nf do 50 i=1,nn 50 read(1,*) k,x(i),y(i),(p(nf*(i-1)+j),j=1,nf) read(1,*) ((jz(i,j),j=1,2),i=1,nc) do 100 i=1,ne 100 read(1,*) ie,(mea(i,j),j=1,4),qq(i) read(1,*) ((aai(i,j),j=1,2),i=1,na) read(1,*) ((aeu(i,j),j=1,2),i=1,nm) close(1) write(2,460) nn,ne,nm,na,nc write(2,465) (i,x(i),y(i),p(3*i-2),p(3*i-1),p(3*i),i=1,nn) write(2,470) ((jz(i,j),j=1,2),i=1,nc) write(2,475) (i,(mea(i,j),j=1,4),i=1,ne) write(2,480) ((aai(i,j),j=1,2),i=1,na) write(2,485) ((aeu(i,j),j=1,2),i=1,nm) 460 format(/5x,the input of nn,ne,nm,na,nc//(5x,5i5)) 465 format(/5x,the input of x,y,p//(5x,i5,5f10.2)) 470 format(/5x,the input of jz//(5x,i5,5x,i5)) 475 format(/5x,the input of mea,qq//(5x,i5,5x,4i5,f10.2)) 480 format(/5x,the input of aai//(5x,i5,5x,2e15.6)) 485 format(/5x,the input of aeu//(5x,i5,5x,2e15.6)) !c structure stiffness statement call fld(nn,ne,mea,nf,nd,n,nt,ld) do 500 i=1,nt 500 a(i)=0 do 600 ie=1,ne call kep(nn,ie,mea,aeu,aai,x,y,akep) call cr(nn,ie,mea,x,y,r) call tran(6,6,r,rt) call dot(6,6,6,rt,akep,ra) call dot(6,6,6,ra,r,ake) call fis(ie,mea,nf,nd,nfd,is) do 560 i=1,nfd do 560 j=1,nfd if(is(i)-is(j)) 560,520,520 520 ni=is(i) ij=ld(ni)-ni+is(j) a(ij)=a(ij)+ake(i,j) 560 continue 600 continue ! c equivalent node force do 700 ie=1,ne call fixf(nn,ne,ie,mea,x,y,qq,pop) call cr(nn,ie,mea,x,y,r) call tran(6,6,r,rt) call dot(6,6,1,rt,pop,po) call fis(ie,mea,nf,nd,nfd,is) do 650 i=1,nfd ni=is(i) p(ni)=p(ni)-po(i) 650 continue 700 continue call fcc(nc,n,nt,nf,jz,ld,a) call decom(n,nt,a,p,ld) write(2,850) (i,p(3*i-2),p(3*i-1),p(3*i),i=1,nn) w