program molecule_together
      real,allocatable::r(:)
      integer,allocatable::nt(:),bt(:,:),ty(:),ind(:),jnd(:),
     1im(:,:)
      logical,allocatable::used(:)
      character*3,allocatable::sy(:)
      character*80 fn,fnn
      write(6,*)' Renumbers atoms in Tinker files'
      write(6,*)' so that atoms in molecules are together'
      write(6,*)' Old filename:'
      read(5,'(a)')fn
      write(6,*)' New filename:'
      read(5,'(a)')fnn

      open(8,file=fn)
      read(8,*)nat
      allocate(r(3*nat),nt(nat),bt(nat,4),ty(nat),used(nat),
     1ind(nat),sy(nat),jnd(nat),im(nat,2))
      do 1 i=1,nat
      read(8,801)sy(i),(r(3*(i-1)+ix),ix=1,3),ty(i),(bt(i,ix),ix=1,4)
801   format(8x,a3,3f12.6,5i6)
      nt(i)=4
      used(i)=.false.
      do 11 ix=4,1,-1
11    if(bt(i,ix).eq.0)nt(i)=ix-1
1     continue
      close(8)
      
      nmol=1
      im(nmol,1)=1
      im(nmol,2)=1
      ind(1)=1
      jnd(1)=1
      ii=1
      used(1)=.true.
44    ic=0
      do 3 i=1,nat
      if(.not.used(i))then
c      is atom i bond to actual molecule?:
       do 5 ix=1,nt(i)
       do 5 jj=im(nmol,1),im(nmol,2)
       if(ind(jj).eq.bt(i,ix))then
        ic=ic+1
        ii=ii+1
        jnd(i)=ii
        im(nmol,2)=ii
        ind(ii)=i
        used(i)=.true.
c       write(6,*)nmol,':',i,'(',ii,') added'
        goto 501
       endif
5      continue
501    continue
      endif
3     continue
      if(ii.lt.nat)then
c      some atoms were added to nmol,investigate once more:
       if(ic.ne.0)goto 44
c      nothing else was found, try to initiate a new molecule:
       do 6 i=1,nat
       if(.not.used(i))then
        used(i)=.true.
        ii=ii+1
        ind(ii)=i
        jnd(i)=ii
        nmol=nmol+1
        im(nmol,1)=ii
        im(nmol,2)=ii
c       write(6,*)nmol,':',i,'(',ii,') started'
c       read(5,*)i5
        goto 44
       endif
6      continue
      endif

      open(9,file=fnn)
      write(9,901)nat
901   format(i6)
      do 7 imol=1,nmol
      write(6,6000)imol,im(imol,1),im(imol,2)
6000  format(i6,':',i6,' - ',i6)
      do 7 ii=im(imol,1),im(imol,2)
      i=ind(ii)
7     write(9,902)ii,sy(i),(r(3*(i-1)+ix),ix=1,3),ty(i),
     1(jnd(bt(i,ix)),ix=1,nt(i))
902   format(i6,2x,a3,3f12.6,5i6)
      close(9)

      end