dinput.f

关于网格自动生成delaunay算法的

**********************************************************************
*     Read Input data for EFGM with plastic-elastic condition
*
**********************************************************************
      subroutine gdata1(ng,mpoin,kbd,meb,ntype,nalgo,props,nbc,
     $     nbcbm,dbcbm,x,nex,nin,ibex,ibin,ibno,nob,nib,npoin,
     $     nszf,ncrit,hdpara,nincs,neb,nodb)
      implicit real*8(a-h,o-z)
      dimension  props(7),ncrit(3),hdpara(20)
      dimension  x(ng,mpoin+3)
      dimension  nbcbm(3,meb), dbcbm(4,meb)
      dimension  ibex(kbd),ibin(kbd),ibno(kbd,mpoin)
      dimension  nodb(2,meb)
      character  title*80
      read(5,*) ndyn
      read(5,10,err=1000) title
      write(6,200) title
      read(5,20,err=1000) ntype,nalgo,nincs,nbc,npoin,neb
      write(6,210) ntype,nbc,neb,npoin
      write(6,220) nalgo,nincs
*
*---- check dimension size
*
      if(ntype.gt.3) then
         write(7777,*) 'Error Message from subroutine gdata 1 !      '
         write(7777,*) '  ntype is illegal number !  ', ntype
         stop
      endif
      if(nalgo.gt.4) then
         write(7777,*) 'Error Message from subroutine gdata 2 !      '
         write(7777,*) '  nalgo is illegal number !  ', nalgo
         stop
      endif
      if(npoin.gt.mpoin) then
         write(7777,*) 'Error Message from subroutine gdata 3 !      '
         write(7777,*) '  npoin is greater than mpoin !  ', npoin, mpoin
         stop
      endif
      if(neb.gt.meb) then
         write(7777,*) 'Error Message from subroutine gdata 4 !      '
         write(7777,*) '  neb is greater than meb !  ', neb, meb
         stop
      endif
*
*---- read property
*
      read(5,20,err=1000) (ncrit(i),i=1,3)
      read(5,*,err=1000) (props(i),i=1,6)
*
*----- read stress strain curve
*
      if(ncrit(2).eq.1) then
         read(5,*,err=1000) hdpara(1)
      elseif(ncrit(2).eq.2) then
         nline = ncrit(3)
         do iline= 1,nline
            read(5,31) idmy,hdpara(iline*2-1),hdpara(iline*2)
         enddo
      elseif(ncrit(2).eq.3.or.ncrit(2).eq.4) then
         read(5,32) alpha,ann,ystres
         hdpara(1)= ystres
         hdpara(2)= ystres/props(1)
         hdpara(3)= alpha
         hdpara(4)= ann
      endif
*
*---- write property and stress strain curve
*
      write(6,460)
      write(6,470) (props(i),i=1,6)
      write(6,480) ncrit(1)
      if(ncrit(2).eq.1) then
         write(6,490)  hdpara(1)
      else if(ncrit(2).eq.2) then
         write(6,500) ncrit(3)
         write(6,510) 
         write(6,520) (j,hdpara(2*j-1),hdpara(2*j),
     $        j= 1,ncrit(3))
      else if(ncrit(2).eq.3) then
         write(6,530) (hdpara(j),j=1,4)
      else if(ncrit(2).eq.4) then
         write(6,540) (hdpara(j),j=1,4)
      endif
*
*---- Read delaunay information
*
      write(6,350)
      read(5,5100,err=1000) nex,nin
      write(6,5100) nex,nin
      if(nex.ne.0) then
         read(5,5100,err=1000) (ibex(i),i=1,nex)
         write(6,5100) (ibex(i),i=1,nex)
      endif
      if(nin.ne.0) then
         read(5,5100,err=1000) (ibin(i),i=1,nin)
         write(6,5100) (ibin(i),i=1,nin)
      else
         read(5,5100,err=1000) ibin(1)
         write(6,5100) ibin(1)
      endif
      do i=1,nex
         read(5,5100,err=1000) (ibno(i,j),j=1,ibex(i))
         write(6,5100) (ibno(i,j),j=1,ibex(i))
      enddo
      read(5,5100,err=1000) nob,nib
      write(6,5100) nob,nib
      write(6,240)
      read(5,5200,err=1000) (ii,x(1,i),x(2,i),i=1,nob+nib)
      write(6,250) (i,x(1,i),x(2,i),i=1,nob+nib)
*
*---- Boundary information
*
      read(5,40) (ieno,(nodb(i,j),i=1,2),j=1,neb)
      write(6,280)
      write(6,290)  (j,(nodb(i,j),i=1,2),j=1,neb)
*
*---- Boundary condition information
*
      write(6,300)
      read(5,50,err=1000) ((nbcbm(i,j),i=1,3),
     $     (dbcbm(i,j),i=1,4),j=1,nbc)
      write(6,310) ((nbcbm(i,j),i=1,3),
     $     (dbcbm(i,j),i=1,4),j=1,nbc)
      nszf=2*npoin
 10   format(a80)
 20   format(10i5)
 30   format(f15.0)
 31   format(i5,2f15.0)
 32   format(3f15.0)
 40   format(3i5)
 50   format(3i5,4f15.0)
 200  format(1h1//5x,70('*')/5x,a80/5x,70('*'))
 210  format(/5x,'Analysis type            =',i5,
     $     10x,'No. of total points      =',i5,
     $     /5x,'No. of boundary elements =',i5,
     $     /5x,'No. of elements preset the boudary condition =',i5)
 220  format(/5x,'Algorithm of conv.  =',i5/
     $     5x,'No. of Increment    =',i5)
 240  format(/5x,'** Nodal No. and Nodal Coordinates **'//
     $     5x,'  No.     X-coord.     Y-coord.')
 250  format(5x,i5,2d13.5,2x,i5,2d13.5)
 280  format(//5x,'** Information about Boundary Elements **'//5x,
     $     '  no.  n01  n02',10x,'  no.  n01  n02')
 290  format(5x,3i5,10x,3i5)
 300  format(//5x,'** Information of Boundary Condition **'
     $       //5x,'   no. x-dr y-dr  node 1 of x    node 2 of x',
     $     '    node 1 of y    node 2 of y'/)
 310  format(  5x,3i5,4d15.7)
 460  format(/5x,'** material property **'//
     $     5x,'Young modu.  Pois R  Thickness',
     $     '    Density   Yield Strs   Frct.')
 470  format(3x,d13.5,f7.3,2d12.4,d13.5,f7.3)
 480  format(/5x,'ncrit1= ',i5/)
 490  format(5x,'This is bi-linear curve !  Hard= ',d13.4)
 500  format(5x,'This is multi-linear curve !   Line No.= ',i5)
 510  format(/5x,'Line No.       Strain   R. of Hard')
 520  format(5x,i8,2d13.4)
 530  format(5x,'This is Lamberg-Osgood equation'/
     $     5x,'sigmy= ',d13.4,'  strainy= ',d13.4,
     $     '  alpha= ',d13.4,'  n= ',d13.4)
 540  format(5x,'This is Modified Lamberg-Osgood equation'/
     $     5x,'sigmy= ',d13.4,'  strainy= ',d13.4,
     $     '  alpha= ',d13.4,'  n= ',d13.4)
 350  format(/5x,'**Infomation of Delaunay Triangulation**'//)
 5100 format(10i5)
 5200 format(i5,2e15.7)
      return
 1000 write(7777,*) 'Error Message from subroutine gdata 5 !      '
      write(7777,*) 'read unexpected character!!!!!!!!!'
      stop
      end