program sinp implicit none integer io,ia,nat,ix,n,ny,is,j,k,nu,ib,iy,iargc,sr,ic parameter (ny=11) real a(3,3,ny),tol,xa logical lex real ,allocatable::x(:),v(:) integer ,allocatable::ig(:),ii(:,:),si(:,:),am(:,:),lm(:) character*5 sy(ny) character*7 l1 character*7 ,allocatable::lb(:,:),lu(:) character*80 t80,s80 data sy/'C2x ','C2y ','C2z ','sxy ','sxz ','syz ', 1 'C4x ','C4y ','C4z ','i ','C1 '/ data a/ 1.0, 0.0, 0.0, 0.0,-1.0, 0.0, 0.0, 0.0,-1.0, 1 -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0,-1.0, 1 -1.0, 0.0, 0.0, 0.0,-1.0, 0.0, 0.0, 0.0, 1.0, 1 1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0,-1.0, 1 1.0, 0.0, 0.0, 0.0,-1.0, 0.0, 0.0, 0.0, 1.0, 1 -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0, 1 1.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0,-1.0, 0.0, 1 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, -1.0, 0.0, 0.0, 1 0.0, 1.0, 0.0, -1.0, 0.0, 0.0, 0.0, 0.0, 1.0, 1 -1.0, 0.0, 0.0, 0.0,-1.0, 0.0, 0.0, 0.0,-1.0, 1 1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0/ tol=0.2 if(iargc().gt.0)then call getarg(1,s80) read(s80,*)tol endif write(6,600)tol 600 format(' FILE.X -> FILE.INP (G.TXT,CM.TXT as par)',/, 1' Symmetric Cartesian Z-matrix',/, 2' Tolerance =',f10.4,' A') io=0 open(7,file='FILE.X') open(8,file='FILE.INP') 1 read(7,80,end=999,err=999)s80 80 format(a80) io=io+1 if(io.gt.1)write(8,803) write(6,601)io 601 format(/,' Geometry',i6,':') 803 format('--link1--') read(7,*)nat n=3*nat allocate(x(n),ig(nat),lb(nat,3),lu(n),si(nat,3),ii(nat,3),v(n), 1am(nat,nat),lm(nat)) open(71,file='G.TXT') 2 read(71,80,end=888,err=888)t80 write(8,80)t80 goto 2 888 close(71) write(8,*) write(8,80)s80 write(8,*) inquire(file='CM.TXT',exist=lex) if(lex)then open(71,file='CM.TXT') read(71,80)t80 close(71) write(8,80)t80 else write(8,*)'0 1' endif do 3 ia=1,nat write(lb(ia,1),500)' X',ia write(lb(ia,2),500)' Y',ia write(lb(ia,3),500)' Z',ia 500 format(A2,i5) 3 read(7,*)ig(ia),(x(ix+3*(ia-1)),ix=1,3) am=0 lm=0 do 5 is=1,ny 5 call tr(is,a,nat,x,ig,tol,sy(is),am,lm) do ia=1,nat write(6,*)ia,(am(ia,ix),ix=1,lm(ia)) enddo c number of unique variables: nu=1 lu(1)=' 0.0' l1='XXXXXXX' v(1)=0.0 do 43 ia=1,nat do 43 ix=1,3 if(abs(x(ix+3*(ia-1))).lt.tol)then ii(ia,ix)=1 si(ia,ix)=1 lb(ia,ix)=l1 endif 43 continue do 41 ia=1,nat do 41 ix=1,3 xa=x(ix+3*(ia-1)) if(lb(ia,ix).ne.l1)then nu=nu+1 lu(nu)=lb(ia,ix) v(nu)=xa lb(ia,ix)=l1 si(ia,ix)=1 ii(ia,ix)=nu do 42 ib=ia,nat is=1 if(ib.eq.ia)is=ix+1 do 42 iy=is,3 sr=0 do 501 ic=1,lm(ia) 501 if(am(ia,ic).eq.ib)sr=sr+1 if(sr.ne.0)then if(lb(ib,iy).ne.l1)then if(abs(x(iy+3*(ib-1))-xa).lt.tol)then lb(ib,iy)=l1 ii(ib,iy)=nu si(ib,iy)=1 endif if(abs(x(iy+3*(ib-1))+xa).lt.tol)then lb(ib,iy)=l1 ii(ib,iy)=nu si(ib,iy)=-1 endif endif endif 42 continue endif 41 continue do 47 ia=1,nat write(8,801)ig(ia) 801 format(i4,' 0 ',$) do 4 ix=1,3 l1=lu(ii(ia,ix)) if(si(ia,ix).eq.-1)l1(1:1)='-' k=0 do 441 j=1,len(l1) 441 if(l1(j:j).eq.' ')k=k+1 do 45 j=1,k 45 write(8,802)' ' do 44 j=1,len(l1) 44 if(l1(j:j).ne.' ')write(8,802)l1(j:j) 802 format(a1,$) 4 continue 47 write(8,*) write(8,*) do 6 ia=2,nu l1=lu(ia) k=0 do 442 j=1,len(l1) 442 if(l1(j:j).eq.' ')k=k+1 do 451 j=1,k 451 write(8,802)' ' do 62 j=1,len(l1) 62 if(l1(j:j).ne.' ') write(8,802)l1(j:j) 6 write(8,804)v(ia) 804 format(f12.6) write(8,*) inquire(file='AFTER.TXT',exist=lex) if(lex)then open(71,file='AFTER.TXT') 891 read(71,80,end=89,err=89)t80 write(8,80)t80 goto 891 89 close(71) write(8,*) endif deallocate(x,ig,lb,lu,si,ii,v,am,lm) goto 1 999 close(7) close(8) write(6,*)io,' structures' end subroutine asa(x,y,z,r,ia) implicit none integer ia real x,y,z,r(*) x=r(1+3*(ia-1)) y=r(2+3*(ia-1)) z=r(3+3*(ia-1)) return end subroutine tr(is,aa,nat,r,ig,tol,sy,am,lm) implicit none integer nat,ia,sr,ib,is,i,j,ny,ig(nat),am(nat,nat),lm(nat),ic parameter (ny=11) real r(3*nat),a(3,3),tol,x,y,z,xb,yb,zb,d, 1aa(3,3,ny) character*5 sy real,allocatable:: rs(:) allocate(rs(3*nat)) do 1 i=1,3 do 1 j=1,3 1 a(i,j)=aa(i,j,is) do 2 ia=1,nat call asa(x,y,z,r,ia) rs(1+3*(ia-1))=a(1,1)*x+a(1,2)*y+a(1,3)*z rs(2+3*(ia-1))=a(2,1)*x+a(2,2)*y+a(2,3)*z 2 rs(3+3*(ia-1))=a(3,1)*x+a(3,2)*y+a(3,3)*z do 3 ia=1,nat call asa(x,y,z,r,ia) sr=0 do 31 ib=1,nat call asa(xb,yb,zb,rs,ib) d=sqrt((x-xb)**2+(y-yb)**2+(z-zb)**2) if(ig(ia).eq.ig(ib).and.d.lt.tol)then do 5 ic=1,lm(ib) 5 if(am(ib,ic).eq.ia)goto 88 c add ia to symmetry-related atoms of ib lm(ib)=lm(ib)+1 am(ib,lm(ib))=ia 88 do 4 ic=1,lm(ia) 4 if(am(ia,ic).eq.ib)goto 99 c add ib to symmetry-related atoms of ia lm(ia)=lm(ia)+1 am(ia,lm(ia))=ib 99 sr=ib goto 3 endif 31 continue c this coordinate is not symmetry related to anything: if(sr.eq.0)return 3 continue write(6,600)sy 600 format(' Symmetry ',a5,' found') return end