program regrid_m

c    December 2002
c    This program is wriiten to work with EasyMesh
c    Its purpose is to transfor EasyMesh output
c    into input files for drain_m

      logical add
      character*80 aa
      integer nbln
      integer i,j,n,ntri,k
      integer inod,jtri,int,ii
      integer nside
      integer Tri(3000,0:3)
      integer Sup(3000,0:10)
      integer Bound(0:150)
      integer Bmark(150)
      integer nmat
      integer subbase  ! material index of subbase
      integer ntime

      real phead,qflux,tstep
      real lam1,sat,res,hd,ksat,lam2
      real stime

      real X(3000),Y(3000)
      real Blen(3000),Bside
      real Bval(150),Bcoef(150)
      real pinit
      real dum

      print *, 'Please Wait While Program Runs'

       
      Bound(0)=0
      open(3,file='ex3.e')
      open(4,file='ex3.n')
      open(8,file='ex3.s')
      open(5,file='drainin.dat')
      open(6,file='tri.dat')
      open(7,file='extra.dat')
    
  
c    number of tris and nodes
      read(3,*) ntri
      read(4,*) n

      do j=1,ntri
      read(3,*) aa,(Tri(J,I),I=1,3),(dum,i=1,8),Tri(J,0)
c     add 1 to each node number  
      do I=1,3
	Tri(j,I)=Tri(j,I)+1
      enddo
      enddo

      do j=1,n
      read(4,*) aa,X(j),Y(j),nbln
      if(nbln.eq.3.or.nbln.eq.4.or.nbln.eq.5) then
        Bound(0)=Bound(0)+1
        Bound(Bound(0))=j
	Bmark(Bound(0))=nbln
        BLen(j)=0.0
      endif
        
      enddo



      do inod=1,n
      sup(inod,0)=0
      do jtri=1,ntri
      do int=1,3
      if(tri(jtri,int).eq.inod) then
	if(int.eq.1) then
	 j=2
	 k=3
	endif
	if(int.eq.2) then
	 j=3
	 k=1
	endif
	if(int.eq.3) then
	 j=1
	 k=2
	endif

	add=.true.
	do ii=1,sup(inod,0)
	if(tri(jtri,j).eq.sup(inod,ii)) add=.false.
	enddo
	if(add) then
	sup(inod,0)=sup(inod,0)+1
	sup(inod,sup(inod,0))=tri(jtri,j)
	endif

	add=.true.
	do ii=1,sup(inod,0)
	if(tri(jtri,k).eq.sup(inod,ii)) add=.false.
	enddo
	if(add) then
	sup(inod,0)=sup(inod,0)+1
	sup(inod,sup(inod,0))=tri(jtri,k)
	endif

      endif
      enddo
      enddo
      enddo

c    Boundary lengths
c    Only count length in x-direction
      read(8,*) nside
      do j=1,nside
      read(8,*) aa,p1,p2,dum,dum,nbln
      if (nbln.eq.5) then
	Bside=sqrt((X(p1+1)-X(p2+1))**2)  !+(Y(p1+1)-Y(p2+1))**2)
        BLen(P1+1)=Blen(P1+1)+Bside/2.0
        Blen(P2+1)=Blen(P2+1)+Bside/2.0
      endif
      enddo

c    set up boundaries
      read(7,*) phead
      read(7,*) qflux
      read(7,*) tstep
      do j=1,Bound(0)
      if(Bmark(j).eq.3) then  !fixed value on pipe
      Bcoef(j)=1e20
      Bval(j)=1e20*phead
      endif
      if(Bmark(j).eq.4) then  !fixed value on base
      Bcoef(j)=5e18             !smaller value used to id pipe in drain
      Bval(j)=5e18*phead
      endif
      if(Bmark(j).eq.5) then  !fixed flux
      Bcoef(j)=0
      Bval(j)=tstep*qflux*Blen(Bound(j))
      endif
      enddo

c    print out
   
      write(5,*) n,ntri


      do i=1,n
      write(5,*) sup(i,0)
      write(5,*) (sup(i,j), j=1,sup(i,0))
      enddo

      do i=1,ntri
      write(5,*) tri(i,0),tri(i,1),tri(i,2),tri(i,3)

      write(6,*) x(tri(i,1)),y(tri(i,1)),Tri(i,0)
      write(6,*) x(tri(i,2)),y(tri(i,2)),Tri(i,0)
      write(6,*) x(tri(i,3)),y(tri(i,3)),Tri(i,0)
      write(6,*) x(tri(i,1)),y(tri(i,1)),Tri(i,0)

      enddo

 
      do i=1,n
       write(5,*) x(i),y(i)
      enddo
   
      write(5,*) Bound(0)
      write(5,*) (Bound(i),i=1,Bound(0))
      write(5,*) (Bcoef(i),i=1,Bound(0))
      write(5,*) (Bval(i),i=1,Bound(0))

      read(7,*) nmat
      write(5,*) nmat

      do i=1,nmat
      read(7,*) lam1,sat,res,hd,ksat,lam2
      write(5,*)lam1,sat,res,hd,ksat,lam2
      enddo

      read(7,*) subbase
      write(5,*) subbase

      read(7,*) stime
      stime=stime/tstep
      ntime=IFIX(stime)
      write(5,*) ntime
      write(5,*) tstep

      read(7,*) pinit
      write(5,*) pinit

      close(3)
      close(4)
      close(5)
      close(6)
      close(7)
      close(8)
      stop
      end