program tinkerinp character*80 s80,t80 parameter (nat0=1000) character*2 as(nat0) character*1 ok dimension r(nat0,3),ihot(nat0,3),ity(nat0) write(6,600) 600 format(' Makes Gaussian inputs from Tinker geometries',/, 1' listed in STRUCT.LST') write(6,*)' Maximum number of waters:' read(5,*)nwmax write(6,*)' Hydrogen bond waters only or all (B,A):' read(5,'(a)')ok if(ok.eq.'a')ok='A' if(ok.eq.'A')then write(6,*)'All waters' else write(6,*)'Hydrogen bond waters only' endif write(6,*)' Cut-off distance:' read(5,*)rmax is=0 inwmax=0 ix=0 open(3,file='STRUCT.LST',status='old') open(36,file='go') 1 read(3,300,end=999,err=999)s80 300 format(a80) write(6,300)s80 is=is+1 do ipoint=1,len(s80) if(s80(ipoint:ipoint).eq.'.')goto 5 enddo 5 do istart=1,len(s80) if(s80(istart:istart).ne.' ')goto 10 enddo 10 do iend=len(s80),1,-1 if(s80(iend:iend).ne.' ')goto 20 enddo 20 open(33,file=s80(istart:iend),status='old') open(34,file=s80(ipoint+1:iend)//'.inp') open(35,file=s80(ipoint+1:iend)//'.sp.inp') write(36,*)'g03 '//s80(ipoint+1:iend)//'.inp '// 1s80(ipoint+1:iend)//'.out' write(36,*)'g03 '//s80(ipoint+1:iend)//'.sp.inp '// 1'G.OUT' write(36,*)'stripg G.OUT '//s80(ipoint+1:iend)//'.out.sp' write(36,*)'rm '//s80(ipoint+1:iend)//'.chk' write(36,*)'rm Gau* ' write(36,*)'rm G.OUT' t80(1:iend-ipoint+9)='%chk='//s80(ipoint+1:iend)//'.chk' write(34,3410) 3410 format('%mem=124000000',/,'%nproc=1') write(35,3410) write(34,3420)(t80(i:i),i=1,iend-ipoint+9) 3420 format(80a1) write(35,3420)(t80(i:i),i=1,iend-ipoint+9) write(34,3430) 3430 format('#becke3lyp/6-311++G** td=(singlets,nstates=20)',/, 1'iop(2/21=1) geom=(nodistance,noangle) nosymm 6d 10f',/) write(35,3435) 3435 format('#becke3lyp/6-311++G** sp GFInput iop(5/33=3) ',/, 1'iop(2/21=1) geom=(nodistance,noangle) nosymm 6d 10f',/, 2'guess=checkpoint',/) write(34,300)s80 write(35,300)s80 write(34,*) write(35,*) write(34,*)'0 1' write(35,*)'0 1' read(33,*)nat if(nat.gt.nat0)then write(6,*)' Too many atoms' stop endif do 2 i=1,nat 2 read(33,333)as(i),(r(i,k),k=1,3) 333 format(8x,a2,1x,3f12.6) close(33) c c numbers of H(-N) and O(=C): ino=0 inh=0 do 3 i=1,nat xi=r(i,1) yi=r(i,2) zi=r(i,3) c c how many covalently bond Os,Ns,Hs and any atoms: io=0 in=0 ih=0 it=0 do 4 j=1,nat xj=r(j,1) yj=r(j,2) zj=r(j,3) dij=sqrt((xi-xj)**2+(yi-yj)**2+(zi-zj)**2) if(dij.gt.0.6.and.dij.lt.1.5.and.as(j)(1:1).eq.'O')io=io+1 if(dij.gt.0.6.and.dij.lt.1.5.and.as(j)(1:1).eq.'N')in=in+1 if(dij.gt.0.6.and.dij.lt.1.5.and.as(j)(1:1).eq.'H')ih=ih+1 4 if(dij.gt.0.6.and.dij.lt.1.5 )it=it+1 if(as(i)(1:1).eq.'H'.and.io.eq.0.and.in.eq.1)inh=i if(as(i)(1:1).eq.'O'.and.ih.eq.0.and.it.eq.1)ino=i c c ity - whether it is water (1) or not (0): if((as(i)(1:1).eq.'O'.and.ih.eq.2).or. 1 (as(i)(1:1).eq.'H'.and.io.eq.1))then ity(i)=1 else ity(i)=0 endif ihot(i,1)=ih ihot(i,2)=io 3 ihot(i,3)=it xh=r(inh,1) yh=r(inh,2) zh=r(inh,3) xo=r(ino,1) yo=r(ino,2) zo=r(ino,3) nw=0 i=0 31 i=i+1 ih=ihot(i,1) io=ihot(i,2) it=ihot(i,3) xi=r(i,1) yi=r(i,2) zi=r(i,3) c atoms that are not in water: if(ity(i).eq.0)then write(34,3401)as(i)(1:1),(r(i,k),k=1,3) write(35,3401)as(i)(1:1),(r(i,k),k=1,3) 3401 format(a1,3f12.6) else c water: c suppose that the other two water atoms follow i: x2=r(i+1,1) y2=r(i+1,2) z2=r(i+1,3) x3=r(i+2,1) y3=r(i+2,2) z3=r(i+2,3) icc=0 if(ok.eq.'A')then do 6 ii=1,nat if(ity(ii).eq.0)then xii=r(ii,1) yii=r(ii,2) zii=r(ii,3) diii=sqrt((xii-xi)**2+(yii-yi)**2+(zii-zi)**2) dii2=sqrt((xii-x2)**2+(yii-y2)**2+(zii-z2)**2) dii3=sqrt((xii-x3)**2+(yii-y3)**2+(zii-z3)**2) if(diii.lt.rmax.and.dii2.lt.rmax.and.dii3.lt.rmax)then icc=1 goto 9898 endif endif 6 continue 9898 continue else rhd=sqrt((xi-xh)**2+(yi-yh)**2+(zi-zh)**2) rod=sqrt((xi-xo)**2+(yi-yo)**2+(zi-zo)**2) rh2=sqrt((x2-xh)**2+(y2-yh)**2+(z2-zh)**2) ro2=sqrt((x2-xo)**2+(y2-yo)**2+(z2-zo)**2) rh3=sqrt((x3-xh)**2+(y3-yh)**2+(z3-zh)**2) ro3=sqrt((x3-xo)**2+(y3-yo)**2+(z3-zo)**2) if((rhd.lt.rmax.or.rod.lt.rmax).and. 1 (rh2.lt.rmax.or.ro2.lt.rmax).and. 2 (rh3.lt.rmax.or.ro3.lt.rmax))icc=1 endif if(icc.eq.1)then if(nw.eq.nwmax)then write(6,*)' Maximum number of waters overflow' else c all three atoms close, write them: write(34,3401)as(i)(1:1),(r(i,k),k=1,3) write(34,3401)as(i+1)(1:1),(r(i+1,k),k=1,3) write(34,3401)as(i+2)(1:1),(r(i+2,k),k=1,3) write(35,3401)as(i)(1:1),(r(i,k),k=1,3) write(35,3401)as(i+1)(1:1),(r(i+1,k),k=1,3) write(35,3401)as(i+2)(1:1),(r(i+2,k),k=1,3) nw=nw+1 endif i=i+2 endif endif if(i.lt.nat)goto 31 write(34,*) write(35,*) close(34) close(35) write(6,*)nw,' waters' if(nw.gt.inwmax)then inwmax=nw ix=is endif goto 1 999 close(3) close(36) write(6,*)' End of structure list' write(6,*)is,' files' write(6,*)'Maximum of ',inwmax,' waters in structure number ',ix stop end