******************************************************************************* program ddcalc ******************************************************************************* c c Calculate the dipole derivatives based on Mulliken charges c c To compile non-graphics programs use : c f77 $(OPT) ddcalc.f -o ddcalc c c character*72 atomform character*80 aform character*80 record, rec character*24 mulform character*80 freqform character*40 fname,gname character*40 hname character*4 rname character*8 mjunk character*1 tchar1,tempchar1 character*2 tchar2,tempchar2 character*2 AtomType(300) real*8 xr(20,300) real*8 yr(20,300) real*8 zr(20,300) real*8 xr0(300) real*8 yr0(300) real*8 zr0(300) real*8 charge(20,300) real*8 charge0(300) real*8 dd(300) real*8 xd(300) real*8 yd(300) real*8 zd(300) real*8 xdd(300) real*8 ydd(300) real*8 zdd(300) real*8 freq(300) real*8 junk1(300) real*8 junk2(300) real*8 spin(300) real*8 at(300) real*8 dev real*8 xd0,yd0,zd0 integer iatom,im,i,j,ierror integer ifreq,ii,il,ir integer iunit integer inum,icheck,itype integer fflag,mflag integer iAtomType(300) c fflag=0 mflag=0 c atomform = "(a2,11x,f16.12,4x,f16.12,4x,f16.12)" 30 format(a2,11x,f16.12,4x,f16.12,4x,f16.12) c mulform = "(a8,x,f7.3,2x,f6.3)" c freqform = "(3x,i2,4x,f8.6,5x,f7.1,6x,f6.2,39x)" c iunit=3 print *, 'Enter outmol filename (inside single quotes): ' read(5,*) fname open(unit=3,name=fname,err=999) print *,'Enter number of atoms: ' read(5,*) iatom ifreq=3*iatom-6 c c Read in Cartesian coordinates c do 900 j=1,30000 read(3,'(a80)',end=999,err=999) record if(record.eq.'$coordinates') then do i=1,iatom read(3,'(a80)',end=999,err=999) record read (record,atomform,err=999) AtomType(i),xr0(i), 1 yr0(i),zr0(i) if(AtomType(i).eq.'H') iAtomType(i)=1 if(AtomType(i).eq.'C') iAtomType(i)=6 if(AtomType(i).eq.'N') iAtomType(i)=7 if(AtomType(i).eq.'O') iAtomType(i)=8 if(AtomType(i).eq.'Fe') iAtomType(i)=26 if(AtomType(i).eq.'H') at(i)=1. if(AtomType(i).eq.'C') at(i)=6. if(AtomType(i).eq.'N') at(i)=7. if(AtomType(i).eq.'O') at(i)=8. if(AtomType(i).eq.'Fe') at(i)=26. write(6,55) AtomType(i),xr0(i),yr0(i),zr0(i) enddo endif c if(record.eq.' mode au_amu cm-1 km/mol') then do i=1,ifreq read(3,'(a80)',end=999,err=999) record read (record,freqform,err=999) im,junk1(i),freq(i),junk2(i) print *,freq(i) enddo if(mflag.eq.1) goto 1500 fflag=1 endif c if(record.eq.' Mulliken atomic charges:') then read(3,'(a80)',end=999,err=999) record do 1050 i=1,iatom read(3,'(a80)',end=999,err=999) record read(record,mulform,err=999) mjunk,charge0(i),spin(i) print *,charge0(i) 1050 continue if(fflag.eq.1) goto 1500 mflag=1 endif 900 continue 1500 continue close(3) c print *,'Enter displacement in reduced (rms) units:' read(5,*) dev c do 2700 k=1,ifreq print *, 'Enter shifted outmol filename (inside single quotes): ' read(5,*) fname open(unit=3,name=fname,err=999) c c c Read in Cartesian coordinates c do 1900 j=1,30000 read(3,'(a80)',end=999,err=999) record c print *, record if(record.eq.'$coordinates') then do i=1,iatom read(3,'(a80)',end=999,err=999) record read (record,atomform,err=999) AtomType(i),xr(k,i), 1 yr(k,i),zr(k,i) if(AtomType(i).eq.'H') iAtomType(i)=1 if(AtomType(i).eq.'C') iAtomType(i)=6 if(AtomType(i).eq.'N') iAtomType(i)=7 if(AtomType(i).eq.'O') iAtomType(i)=8 if(AtomType(i).eq.'Fe') iAtomType(i)=26 if(AtomType(i).eq.'H') at(i)=1. if(AtomType(i).eq.'C') at(i)=6. if(AtomType(i).eq.'N') at(i)=7. if(AtomType(i).eq.'O') at(i)=8. if(AtomType(i).eq.'Fe') at(i)=26. write(6,55) AtomType(i),xr(k,i),yr(k,i),zr(k,i) enddo endif c if(record.eq.' Mulliken atomic charges:') then read(3,'(a80)',end=999,err=999) record do 2050 i=1,iatom read(3,'(a80)',end=999,err=999) record read(record,mulform,err=999) mjunk,charge(k,i),spin(i) print *,charge(k,i) 2050 continue goto 2500 endif 1900 continue 2500 continue close(3) 2700 continue c do 1600 i=1,ifreq xd(i)=0.0 yd(i)=0.0 zd(i)=0.0 1600 continue xd0=0.0 yd0=0.0 zd0=0.0 do 1700 i=1,iatom xd0=xd0+xr0(i)*charge0(i)*sqrt(at(i)) yd0=yd0+yr0(i)*charge0(i)*sqrt(at(i)) zd0=zd0+zr0(i)*charge0(i)*sqrt(at(i)) do 1650 k=1,ifreq xd(k)=xd(k)+xr(k,i)*charge(k,i)*sqrt(at(i)) yd(k)=yd(k)+yr(k,i)*charge(k,i)*sqrt(at(i)) zd(k)=zd(k)+zr(k,i)*charge(k,i)*sqrt(at(i)) 1650 continue 1700 continue do 1750 k=1,ifreq xdd(k)=xd(k)-xd0 ydd(k)=yd(k)-yd0 zdd(k)=zd(k)-zd0 xdd(k)=xdd(k)/(dev*0.0165672*sqrt(freq(k))) ydd(k)=ydd(k)/(dev*0.0165672*sqrt(freq(k))) zdd(k)=zdd(k)/(dev*0.0165672*sqrt(freq(k))) 1750 continue do 1850 k=1,ifreq dd(k)=(xdd(k)*xdd(k)+ydd(k)*ydd(k)+zdd(k)*zdd(k))/3 1850 continue print *,'Enter output file name: ' read(5,*) gname junit=7 open(unit=7,file=gname,form='formatted',err=999) do 100 i=1,ifreq write (junit,20) freq(i),dd(i) 100 continue close(7) c 20 format(f10.5,2x,f10.5) 55 format(a7,f13.9,2x,f13.9,2x,f13.9,a30) goto 3001 995 continue print *,'Error in format' print *, record goto 999 996 continue print *,'Error parsing format' print *, record goto 999 997 continue print *,'End of file detected' 999 continue print *,'Error in file' goto 3002 3001 continue print *,'PROGRAM COMPLETE' print *,'DMOL INPUT FILES HAVE BEEN GENERATED' 3002 continue end *************************************************************************** Subroutine int2char(i,tempchar1,tempchar2) *************************************************************************** ! To convert the integer data to character data Implicit None Integer i Character*1 tempchar1 Character*2 tempchar2 Character*4 record if(i.lt.10) then Write (record,'(i1)') i Read (record,'(a1)') tempchar1 else Write (record,'(i2)') i Read (record,'(a2)') tempchar2 endif Return End