od_mapatom.F90 4.46 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34
      MODULE m_od_mapatom
      use m_juDFT
      CONTAINS 
      SUBROUTINE od_mapatom(oneD,atoms,sym,cell)

!      written by Y.Mokrousov in order to generate the arrays connected
!      to the operations, transforming atoms into each other,
!      for more details look in mapatom.F.    year 2004 

      USE m_types
      IMPLICIT NONE
      TYPE(t_oneD),INTENT(IN)    :: oneD
      TYPE(t_atoms),INTENT(INOUT):: atoms
      TYPE(t_sym),INTENT(INOUT)  :: sym
      TYPE(t_cell),INTENT(IN)    :: cell


      REAL ij,pps(3),norm,aamat(3,3) 
      INTEGER i,j,n1,k,n,n2,np1,na,ix,iy,iz,nat1,nat2,na2
      REAL mt(3,3),sum_tau_lat(3),sum_taual(3)
      REAL, PARAMETER :: del = 1.0e-4

      aamat=matmul(transpose(cell%amat),cell%amat)
     
      n1 = 1
      DO n = 1,atoms%ntype
        n2 = n1 + atoms%neq(n) - 1
        IF (atoms%neq(n).EQ.1) THEN
35 36
           sym%ngopr(n2) = 1
           WRITE (6,FMT=8010) n2,n2,sym%ngopr(n2)
37 38 39 40 41 42 43 44
           n1 = n1 + atoms%neq(n)
           CYCLE
        END IF
        DO  na = n1,n2
           DO np1 = 1,oneD%odd%nop
                  pps =matmul(sym%mrot(:,:,np1),atoms%taual(:,n1))
                  pps(3) = pps(3)+sym%tau(3,np1)/cell%amat(3,3)
                  IF (all(abs(atoms%taual(:,na)-pps(:)).LE.1.e-4 )) THEN
45 46
                      sym%ngopr(na) = np1
                      WRITE (6,FMT=8010) na,n1,sym%ngopr(na)
47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81
 8010                 FORMAT (5x,'atom',i3,' can be mapped into atom', i3,' through group  operation',i4)
                      CYCLE
                   END IF
           END DO
      ENDDO
        n1 = n1 + atoms%neq(n)
    ENDDO

!---> defining inverse operations for the Hamiltonian and forces
!     where we do not need to consider the translational part

      DO n1 = 1,oneD%odd%nop
         n2loop:DO n2 = 1,oneD%odd%nop
            mt=matmul(sym%mrot(:,:,n1),sym%mrot(:,:,n2))
            DO n = 1,oneD%odd%nop
               if (all(abs(mt(:,:) - sym%mrot(:,:,n)).LE.1.e-06)) THEN
                     sym%multab(n1,n2) = n
                     IF (n.EQ.1) sym%invtab(n1) = n2
                     cycle n2loop
               endif
            enddo
            WRITE (6,FMT=8050) n1,n2
 8050       FORMAT (' error - n1,n2=',2i3)
            CALL juDFT_error("mult",calledby ="od_mapatom")
        ENDDO n2loop
    ENDDO

 8060 FORMAT (1x,24i3)

      WRITE (6,FMT='(//," inverse operations",//)')

      DO n1 = 1,oneD%odd%nop
         WRITE (6,FMT=8060) n1,sym%invtab(n1)
      END DO

Daniel Wortmann's avatar
Daniel Wortmann committed
82
      DO na = 1,atoms%nat
83
         sym%invsat(na) = 0
84 85 86 87 88 89 90 91 92
         sym%invsatnr(na) = 0
      END DO

      IF (oneD%odd%invs) THEN
         WRITE (6,FMT=*)
         nat1 = 1
         DO n = 1,atoms%ntype
            nat2 = nat1 + atoms%neq(n) - 1
            DO na = nat1,nat2 - 1
93
               IF (sym%invsat(na).EQ.0) THEN
94 95 96 97 98 99 100 101 102 103 104 105
                  naloop:DO na2 = na + 1,nat2
                     DO i = 1,3
                        sum_taual(i) = atoms%taual(i,na) + atoms%taual(i,na2)
                     END DO
                     DO ix = -2,2
                       sum_tau_lat(1) = sum_taual(1) + real(ix)
                       DO iy = -2,2
                         sum_tau_lat(2) = sum_taual(2) + real(iy)
                         DO iz = -2,2
                           sum_tau_lat(3) = sum_taual(3) + real(iz)
                           norm = sqrt(dot_product(matmul(sum_tau_lat,aamat),sum_tau_lat))
                           IF (norm.LT.del) THEN
106 107
                              sym%invsat(na) = 1
                              sym%invsat(na2) = 2
108 109 110 111 112 113 114 115 116 117 118 119 120 121
                              sym%invsatnr(na)  = na2
                              sym%invsatnr(na2) = na
                              WRITE (6,FMT=9000) n,na,na2
                              CYCLE naloop
                           END IF
                        END DO
                      END DO
                    END DO
                  ENDDO naloop
               END IF
            END DO
            nat1 = nat1 + atoms%neq(n)
         END DO
      END IF
122
      WRITE (6,FMT=*) sym%invsat
123 124 125 126
 9000 FORMAT ('atom type',i3,': atom',i3,' can be mapped into atom',i3, ' via 3d inversion')

      END SUBROUTINE od_mapatom
      END MODULE m_od_mapatom