program v3atm IMPLICIT REAL*8 (A-H,O-Z) IMPLICIT INTEGER*4 (I-N) parameter (N0=1000) dimension r(3),ra(3),rb(3),rc(3),p(3) DIMENSION APT(N0*3,3),AAT(N0*3,3) DIMENSION ALPHA(N0*3,3,3),G(N0*3,3,3),A(N0*3,3,3,3) character*50 ias,ibs,ics write(6,*)'rewrites tensor derivatives into molecular frame' write(6,*)'defined by 3 atoms, vectors a-> b a-> c' Z0=0.0d0 if(iargc().lt.3)call report('Specify 3 integers in line!') call getarg(1,ias) call getarg(2,ibs) call getarg(3,ics) read(ias,*)ia read(ibs,*)ib read(ics,*)ic open(9,file='FILE.X') read(9,*) read(9,*)nat do 1 i=1,nat read(9,*)idum,(r(ix),ix=1,3) if(ia.eq.i)then do 4 ix=1,3 4 ra(ix)=r(ix) endif if(ib.eq.i)then do 2 ix=1,3 2 rb(ix)=r(ix) endif if(ic.eq.i)then do 3 ix=1,3 3 rc(ix)=r(ix) endif 1 continue close(9) do 5 ix=1,3 rb(ix)=rb(ix)-ra(ix) 5 rc(ix)=rc(ix)-ra(ix) call vp(rb,rc,ra) call vp(ra,rb,rc) call norm(ra) call norm(rb) call norm(rc) write(6,*)'FILE.Q.TEN' OPEN(15,FILE='FILE.Q.TEN') read(15,*) NAT if(NAT.gt.N0)call report('too many atoms') DO 10 L=1,NAT*3 read(15,1501) (p(ix),ix=1,3) APT(L,1)=0.0d0 APT(L,2)=0.0d0 APT(L,3)=0.0d0 do 10 ix=1,3 APT(L,1)=APT(L,1)+p(ix)*ra(ix) APT(L,2)=APT(L,2)+p(ix)*rb(ix) 10 APT(L,3)=APT(L,3)+p(ix)*rc(ix) DO 11 L=1,NAT*3 read(15,1501) (p(ix),ix=1,3) AAT(L,1)=0.0d0 AAT(L,2)=0.0d0 AAT(L,3)=0.0d0 do 11 ix=1,3 AAT(L,1)=AAT(L,1)+p(ix)*ra(ix) AAT(L,2)=AAT(L,2)+p(ix)*rb(ix) 11 AAT(L,3)=AAT(L,3)+p(ix)*rc(ix) close(15) OPEN(15,FILE='FILE.Q.A.TEN') WRITE(15,1500) NAT,NAT-6,0,ia,ib,ic 1500 FORMAT(3I5,' normal mode derivatives, frame:',3I5) DO 12 L=1,NAT*3 12 WRITE(15,1501) (APT(L,I),I=1,3),L 1501 FORMAT(3F14.8,I5) DO 221 L=1,NAT*3 221 WRITE(15,1501) (AAT(L,I),I=1,3),L DO 230 L=1,NAT*3 230 WRITE(15,1501) (Z0,I=1,3),L DO 100 L=1,NAT*3 100 WRITE(15,1501) (APT(L,I),I=1,3),L CLOSE(15) write(6,*)'FILE.Q.TTT' OPEN(2,FILE='FILE.Q.TTT',STATUS='OLD') READ(2,*) READ(2,*)NAT if(NAT.gt.N0)call report('too many atoms') READ(2,*) READ(2,*) DO 15 I=1,3 READ(2,*) DO 15 IND=1,NAT*3 15 READ(2,*)LDUM,IXDUM,(ALPHA(IND,I,J),J=1,3) READ(2,*) READ(2,*) DO 13 I=1,3 READ(2,*) DO 13 IND=1,NAT*3 13 READ(2,*)LDUM,IXDUM,(G(IND,I,J),J=1,3) READ(2,*) READ(2,*) DO 14 I=1,3 DO 14 J=1,3 READ(2,*) DO 14 IND=1,NAT*3 14 READ(2,*)LDUM,IXDUM,(A(IND,I,J,K),K=1,3) CLOSE(2) DO 141 IND=1,NAT*3 DO 141 I=1,3 ax=ALPHA(IND,I,1) ay=ALPHA(IND,I,2) az=ALPHA(IND,I,3) ALPHA(IND,I,1)=ax*ra(1)+ay*ra(2)+az*ra(3) ALPHA(IND,I,2)=ax*rb(1)+ay*rb(2)+az*rb(3) ALPHA(IND,I,3)=ax*rc(1)+ay*rc(2)+az*rc(3) ax=G(IND,I,1) ay=G(IND,I,2) az=G(IND,I,3) G(IND,I,1)=ax*ra(1)+ay*ra(2)+az*ra(3) G(IND,I,2)=ax*rb(1)+ay*rb(2)+az*rb(3) G(IND,I,3)=ax*rc(1)+ay*rc(2)+az*rc(3) DO 141 J=1,3 ax=A(IND,I,J,1) ay=A(IND,I,J,2) az=A(IND,I,J,3) A(IND,I,J,1)=ax*ra(1)+ay*ra(2)+az*ra(3) A(IND,I,J,2)=ax*rb(1)+ay*rb(2)+az*rb(3) 141 A(IND,I,J,3)=ax*rc(1)+ay*rc(2)+az*rc(3) DO 142 IND=1,NAT*3 DO 142 I=1,3 ax=ALPHA(IND,1,I) ay=ALPHA(IND,2,I) az=ALPHA(IND,3,I) ALPHA(IND,1,I)=ax*ra(1)+ay*ra(2)+az*ra(3) ALPHA(IND,2,I)=ax*rb(1)+ay*rb(2)+az*rb(3) ALPHA(IND,3,I)=ax*rc(1)+ay*rc(2)+az*rc(3) ax=G(IND,1,I) ay=G(IND,2,I) az=G(IND,3,I) G(IND,1,I)=ax*ra(1)+ay*ra(2)+az*ra(3) G(IND,2,I)=ax*rb(1)+ay*rb(2)+az*rb(3) G(IND,3,I)=ax*rc(1)+ay*rc(2)+az*rc(3) DO 142 J=1,3 ax=A(IND,I,1,J) ay=A(IND,I,2,J) az=A(IND,I,3,J) A(IND,I,1,J)=ax*ra(1)+ay*ra(2)+az*ra(3) A(IND,I,2,J)=ax*rb(1)+ay*rb(2)+az*rb(3) 142 A(IND,I,3,J)=ax*rc(1)+ay*rc(2)+az*rc(3) DO 143 IND=1,NAT*3 DO 143 I=1,3 DO 143 J=1,3 ax=A(IND,1,I,J) ay=A(IND,2,I,J) az=A(IND,3,I,J) A(IND,1,I,J)=ax*ra(1)+ay*ra(2)+az*ra(3) A(IND,2,I,J)=ax*rb(1)+ay*rb(2)+az*rb(3) 143 A(IND,3,I,J)=ax*rc(1)+ay*rc(2)+az*rc(3) open(22,file='FILE.Q.A.TTT') WRITE(22,2001)NAT,ia,ib,ic 2001 FORMAT(' ROA tensors, normal modes derivatives',/, 1I4,' atoms',3i5,/, 1' The electric-dipolar electric-dipolar polarizability:',/, 2' Atom/x jx jy jz') CALL WRITETTT(N0*3,NAT,ALPHA,A,G) close(22) stop end SUBROUTINE norm(a) REAL*8 a(3),an an=dsqrt(a(1)**2+a(2)**2+a(3)**2) a(1)=a(1)/an a(2)=a(2)/an a(3)=a(3)/an return end SUBROUTINE vp(a,b,c) REAL*8 a(3),b(3),c(3) c(1)=a(2)*b(3)-a(3)*b(2) c(2)=a(3)*b(1)-a(1)*b(3) c(3)=a(1)*b(2)-a(2)*b(1) return end SUBROUTINE REPORT(RS) IMPLICIT INTEGER*4 (I-N) IMPLICIT REAL*8 (A-H,O-Z) CHARACTER*(*) RS WRITE(3,3000) WRITE(6,3000) 3000 FORMAT(/,80(1H*)) WRITE(3,*)RS WRITE(6,*)RS WRITE(3,3001) WRITE(6,3001) 3001 FORMAT(80(1H*),/,/,'PROGRAM STOPPED') CLOSE(3) CLOSE(2) STOP END SUBROUTINE WRITETTT(N30,NAT,ALPHA,A,G) IMPLICIT REAL*8 (A-H,O-Z) IMPLICIT INTEGER*4 (I-N) DIMENSION ALPHA(N30,3,3),G(N30,3,3),A(N30,3,3,3) c DO 1 I=1,3 WRITE(22,2002)I 2002 FORMAT(' Alpha(',I1,',J):') DO 1 L=1,NAT DO 1 IX=1,3 II=3*(L-1)+IX 1 WRITE(22,2001)L,IX,(ALPHA(II,I,J),J=1,3) 2001 FORMAT(I5,1H ,I1,3F13.7) WRITE(22,2004) 2004 FORMAT(' The electric dipole magnetic dipole polarizability:',/, 1 ' Atom/x jx(Bx) jy(By) jz(Bz)') DO 2 I=1,3 WRITE(22,2003)I 2003 FORMAT(' G(',I1,',J):') DO 2 L=1,NAT DO 2 IX=1,3 II=3*(L-1)+IX c last inner index magnetic 2 WRITE(22,2001)L,IX,(G(II,I,J),J=1,3) WRITE(22,2005) 2005 FORMAT(' The electric dipole electric quadrupole polarizability:', 2/, ' Atom/x kx ky kz') DO 3 I=1,3 DO 3 J=1,3 WRITE(22,2006)I,J 2006 FORMAT(' A(',I1,'(u),',I1,',K):') DO 3 L=1,NAT DO 3 IX=1,3 II=3*(L-1)+IX 3 WRITE(22,2001)L,IX,(A(II,I,J,K),K=1,3) WRITE(22,2009) 2009 FORMAT('The velocity form of alpha',/, 1' Atom/x vx vy vz') DO 11 I=1,3 WRITE(22,2012)I 2012 FORMAT(' AlphaV(',I1,',J):') DO 11 L=1,NAT DO 11 IX=1,3 II=3*(L-1)+IX 11 WRITE(22,2001)L,IX,(ALPHA(II,I,J),J=1,3) RETURN END