boxdim.f 4.29 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 35 36 37 38 39 40 41 42 43 44 45 46 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 82 83 84
      MODULE m_boxdim
      CONTAINS
      SUBROUTINE boxdim(
     >                  bmat,
     <                  arltv1,arltv2,arltv3)
c*********************************************************************
c      This subroutine determines the maximum number L, M and N
c      nrgvl, nrgvm, nrgvn, of reciprocal lattice vectors G(L,M,N) 
c      along the directions G(1), G(2), G(3), respectively, for which 
c
c                     | G(L,M,N) | <= GMAX.
c
c      This equation defines a "sphere" and the nrgvl,m,n define
c      the DIMension of the BOX in which the sphere is placed.
c
c      In reality the G(i)'s, do not form a carteasian coordinate 
c      system. Therefore the "sphere" is not a sphere, but an 
c      ellipsoid. For this ellipsoid the largest L, M and N component 
c      is determined as arltv1, arltv2, arltv3 to construct the boxes and
c
c                    nrgvl  = int( gmax/arltv1 ) + 1
c                    nrgvm  = int( gmax/arltv2 ) + 1
c                    nrgvn  = int( gmax/arltv3 ) + 1
c
c      G(i,xyz) is stored in bmat(i,xyz)
c
c      routine by s.bluegel from carpar-program
c 
c                         S. Bl"ugel, IFF, 13. Nov. 97    
c               tested by S. Heinze , IFF, 
c*********************************************************************
      use m_juDFT
c
C     .. Parameters ..
      IMPLICIT NONE
C
C     .. Scalar Arguments ..
      REAL,    INTENT (OUT) :: arltv1,arltv2,arltv3
C     ..
C     .. Array Arguments ..
      REAL,    INTENT (IN)  :: bmat(3,3)
C     ..
C     .. Local Scalars ..
      INTEGER ixyz,j,k
      REAL    denom,eps,one,zero
C     ..
C     .. Local Arrays ..
      REAL det(3,3),rr(3,3)
c     ..
      DATA one,eps,zero/1.0,1e-10,0.0/      
c
c--->  build up quadratic form for ellipsoid
c
      DO j = 1 , 3
         DO k = 1 , 3
            rr(k,j) = zero
            DO ixyz = 1 , 3
               rr(k,j) = rr(k,j) + bmat(k,ixyz)*bmat(j,ixyz)
            ENDDO
         ENDDO
      ENDDO
c
c---> build determinants for Cramer's rule
c
      det(1,1) = rr(2,2)*rr(3,3) - rr(3,2)*rr(2,3)
      det(1,2) = rr(2,1)*rr(3,3) - rr(3,1)*rr(2,3)
      det(1,3) = rr(2,1)*rr(3,2) - rr(3,1)*rr(2,2)
      det(2,1) = rr(1,2)*rr(3,3) - rr(3,2)*rr(1,3)
      det(2,2) = rr(1,1)*rr(3,3) - rr(3,1)*rr(1,3)
      det(2,3) = rr(1,1)*rr(3,2) - rr(3,1)*rr(1,2)
      det(3,1) = rr(1,2)*rr(2,3) - rr(2,2)*rr(1,3)
      det(3,2) = rr(1,1)*rr(2,3) - rr(2,1)*rr(1,3)
      det(3,3) = rr(1,1)*rr(2,2) - rr(2,1)*rr(1,2)
c
c---> check on the zeros of some determinants
c
      DO j = 1 , 3
         IF ( det(j,j) .lt. eps ) THEN
85
            WRITE (6,
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 119 120 121
     +                '('' problem with det('',i1,'','',i1,'')'')') j,j
            CALL juDFT_error(" boxdim: determinant",calledby ="boxdim")
         END IF
      ENDDO    
c
c---> scale determinants
c
      DO k = 1 , 3
         denom = one / det(k,k)
         DO j = 1 , 3
            det(k,j) = denom * det(k,j)
         ENDDO
         det(k,k) = one / denom
      ENDDO
c
c---> calculate the maximum l, m, n components of the ellipsoid
c
      arltv1 = sqrt( rr(1,1) + rr(2,2)*det(1,2)*det(1,2)
     >                       + rr(3,3)*det(1,3)*det(1,3)
     >                       - (rr(1,2)+rr(1,2))*det(1,2)
     >                       + (rr(1,3)+rr(1,3))*det(1,3)
     >                       - (rr(2,3)+rr(2,3))*det(1,2)*det(1,3))
      arltv2 = sqrt( rr(2,2) + rr(1,1)*det(2,1)*det(2,1)
     >                       + rr(3,3)*det(2,3)*det(2,3)
     >                       - (rr(1,2)+rr(1,2))*det(2,1)
     >                       + (rr(1,3)+rr(1,3))*det(2,1)*det(2,3)
     >                       - (rr(2,3)+rr(2,3))*det(2,3))
      arltv3 = sqrt( rr(3,3) + rr(1,1)*det(3,1)*det(3,1)
     >                       + rr(2,2)*det(3,2)*det(3,2)
     >                       - (rr(1,2)+rr(1,2))*det(3,1)*det(3,2)
     >                       + (rr(1,3)+rr(1,3))*det(3,1)
     >                       - (rr(2,3)+rr(2,3))*det(3,2))

      RETURN
      END SUBROUTINE
      END