rw_symfile.f 5.18 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
      MODULE m_rwsymfile
      use m_juDFT
!----------------------------------------------------------------------!
!     writes spacegroup operations                                     ! 
!     and                                                              |  
!     rec. lattice vectors (for external k-point generator)            !
!----------------------------------------------------------------------!
      CONTAINS
      SUBROUTINE rw_symfile(
     >                      rw,symfh,symfn,nopd,bmat,
     X                      mrot,tau,nop,nop2,symor)

      IMPLICIT NONE

!===> Arguments
      CHARACTER(len=1), INTENT (IN)    :: rw
      CHARACTER(len=7), INTENT (IN)    :: symfn
      INTEGER,          INTENT (IN)    :: nopd,symfh
      REAL,             INTENT (IN)    :: bmat(3,3)
      INTEGER,          INTENT (INOUT) :: nop,nop2
      INTEGER,          INTENT (INOUT) :: mrot(3,3,nopd)
      REAL,             INTENT (INOUT) :: tau(3,nopd)
      LOGICAL,          INTENT (INOUT) :: symor

!===> Variables
26
      INTEGER i,j,n,ios,no2,no3,gen1,isrt(nop)
27
      REAL    t,d
28
      LOGICAL ex,op,l_exist
29 30 31 32 33 34 35 36 37 38 39
      CHARACTER(len=3) :: type
      CHARACTER(len=7) :: sym2fn

      sym2fn = 'sym.out'

      IF ( SCAN(rw,'wW') > 0 ) THEN

!===> write symfile

        OPEN (symfh, file=sym2fn, status='unknown', err=911, iostat=ios)
        WRITE (symfh,*) nop,nop2,symor,'    ! nop,nop2,symor '
40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55
        
        IF (nop == 2*nop2) THEN  ! film-calculation
          i = 1 ; j = nop2 + 1
          DO n = 1, nop
            IF (mrot(3,3,n) == 1) THEN
              isrt(n) = i ; i = i + 1
            ELSE
              isrt(n) = j ; j = j + 1
            ENDIF
          ENDDO
        ELSE
          DO n = 1, nop
            isrt(n) = n
          ENDDO
        ENDIF

56 57
        DO n = 1, nop
           WRITE (symfh,'(a1,i3)') '!', n
58 59
           WRITE (symfh,'(3i5,5x,f10.5)')
     &           ((mrot(i,j,isrt(n)),j=1,3),tau(i,isrt(n)),i=1,3)
60
        ENDDO
61

62 63 64 65 66 67
!        WRITE (symfh,*) '! reciprocal lattice vectors'
!        WRITE (symfh,'(3f25.15)') ((bmat(i,j),j=1,3),i=1,3)

      ELSEIF ( SCAN(rw,'rR') > 0 ) THEN

!===> read symfile
68 69 70 71 72
        INQUIRE(FILE=TRIM(ADJUSTL(symfn)),EXIST=l_exist)
        IF(.NOT.l_exist) THEN
           CALL juDFT_error("File "//TRIM(ADJUSTL(symfn))//
     +                      " is missing.",calledby="rw_symfile")
        END IF
73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118
        OPEN (symfh, file=trim(symfn),status='old',err=911,iostat=ios)
        READ (symfh,*) nop,nop2,symor
        IF (symfn.EQ.'sym.out') THEN
          gen1 = 0
        ELSEIF (trim(symfn).EQ.'sym') THEN
          gen1 = 1
        ELSE
           CALL juDFT_error("symfn should be sym or sym.out",calledby
     +          ="rw_symfile")
        ENDIF
        DO n = 1 + gen1, nop + gen1
          READ (symfh,*)
          READ (symfh,*) 
     &         ((mrot(i,j,n),j=1,3),tau(i,n),i=1,3)
        ENDDO
        IF (symor) THEN
          DO n=1,nop
            t= tau(1,n)**2 + tau(2,n)**2 + tau(3,n)**2
            IF (t > 1.e-8)  CALL juDFT_error("not symmorphic",calledby
     +           ="rw_symfile")
          ENDDO
        ELSE
          DO n=1,nop
            DO i = 1,3
             IF (ABS(tau(i,n)-0.33333) < 0.00001) THEN
               tau(i,n) = 1./3.
             ENDIF
             IF (ABS(tau(i,n)+0.33333) < 0.00001) THEN
               tau(i,n) = -1./3.
             ENDIF
             IF (ABS(tau(i,n)-0.66667) < 0.00001) THEN
               tau(i,n) = 2./3.
             ENDIF
             IF (ABS(tau(i,n)+0.66667) < 0.00001) THEN
               tau(i,n) = -2./3.
             ENDIF
             IF (ABS(tau(i,n)) > 0.00001) THEN
             ENDIF
            ENDDO
          ENDDO
        ENDIF

        DO n = 1,nop
!
! Determine the kind of symmetry operation we have here
!
119
          d = det(mrot(:,:,n))
120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157
          t =  mrot(1,1,n) + mrot(2,2,n) + mrot(3,3,n)

          IF (d.EQ.-1) THEN
            type = 'm  '
            IF (t.EQ.-3) type = 'I  '
          ELSEIF (d.EQ.1) THEN
            IF (t.EQ.-1) type = 'c_2'
            IF (t.EQ. 0) type = 'c_3'
            IF (t.EQ. 1) type = 'c_4'
            IF (t.EQ. 2) type = 'c_6'
            IF (t.EQ. 3) type = 'E  '
          ELSE
             CALL juDFT_error("determinant =/= +/- 1",calledby
     +            ="rw_symfile")
          ENDIF
 
          WRITE (6,FMT=8020) n, type
 8020     FORMAT (/,1x,i3,' : ',a3)
          DO i = 1,3
             WRITE (6,FMT=8030) (mrot(i,j,n),j=1,3),tau(i,n)
          ENDDO
 8030     FORMAT (5x,3i3,3x,f4.1)
        ENDDO

      ELSE
         CALL juDFT_error("ERROR! rw_symfile #1",calledby="rw_symfile")
      ENDIF
      CLOSE (symfh)
      RETURN

! === errors
 911  CONTINUE
      WRITE(*,*) 'Error in inquire. IOS=',ios
 912  CONTINUE
      WRITE(*,*) 'Error in open. IOS=',ios
       CALL juDFT_error("i/o ERROR",calledby="rw_symfile")

      END SUBROUTINE rw_symfile
158 159 160 161 162 163 164 165 166 167 168

!--------------------------------------------------------------------
      INTEGER FUNCTION det(m)
        INTEGER m(3,3)
        det = m(1,1)*m(2,2)*m(3,3) + m(1,2)*m(2,3)*m(3,1) +
     +        m(2,1)*m(3,2)*m(1,3) - m(1,3)*m(2,2)*m(3,1) -
     +        m(2,3)*m(3,2)*m(1,1) - m(2,1)*m(1,2)*m(3,3)
      END FUNCTION det
!--------------------------------------------------------------------


169
      END MODULE m_rwsymfile